吉成才, 關(guān)英杰, 鄭 皓, 張 明
(長沙礦冶研究院有限責(zé)任公司 深海礦產(chǎn)資源開發(fā)利用國家重點實驗室,湖南長沙410012)
隨著陸地礦產(chǎn)資源的衰竭,深海礦產(chǎn)資源開發(fā)技術(shù)越來越被世界各國所重視[1],目前最有商業(yè)前景的開采方案是管道提升法[2-3]。 礦漿泵作為采礦輸送系統(tǒng)中最核心的動力裝備,其水力結(jié)構(gòu)復(fù)雜,顆粒極易堵塞泵管流道。 國內(nèi)開發(fā)的礦漿泵過流和回流能力往往不足,導(dǎo)致流道堵塞[4-5]。 流道堵塞與顆粒及顆粒群的流態(tài)、顆粒流的時空分布和動態(tài)遷移直接關(guān)聯(lián),要解決上述工程問題,必須從泵內(nèi)顆粒兩相流的深層規(guī)律和堵塞機理上深入研究。
目前對于深海礦漿泵的研究,主要是基于傳統(tǒng)的固液兩相流模型[6-10]。 深海采礦礦漿泵帶有空間導(dǎo)葉,與普通的蝸殼式離心泵結(jié)構(gòu)有較大差別,輸送工況也具有大顆粒低濃度的特點,目前將離散元和流體動力學(xué)耦合方法應(yīng)用到礦漿泵流場數(shù)值計算中的研究還不多。 本文采用CFD-DEM 耦合算法對礦漿泵內(nèi)顆粒運動特性進行數(shù)值模擬,分析顆粒粒徑對礦漿泵內(nèi)流動特性的影響,并分析流道堵塞的機理,旨在為未來的深海采礦工程泵開發(fā)提供理論依據(jù)。
CFD-DEM 耦合算法的求解策略由CFD 求解器求解連續(xù)相、由DEM 求解器求解離散相,所以CFD-DEM求解分為3 個部分:連續(xù)相求解、離散相求解和兩相間的耦合。
流體連續(xù)相求解主要應(yīng)用連續(xù)性方程和動量方程(Navier-Stokes 方程)。
連續(xù)性方程又稱為質(zhì)量守恒方程,其表述為:
式中ρ 為流體密度;t 為時間;vi為速度矢量在笛卡爾直角坐標(biāo)系i 方向的分量。
動量守恒方程也稱作Navier-Stokes 方程,x、y、z 3 個方向的動量守恒方程為:
式中P 為靜壓;ρgi和Fi分別為i 方向上的重力體積力和外部體積力;τij為由流體運動引起的黏性應(yīng)力張量。
離散相的求解主要是通過接觸模型求解顆粒碰撞過程的受力,并運用牛頓第二定律計算出顆粒的加速度,然后更新顆粒的速度和位移。 對粒子的受力進行分析,其對應(yīng)的動量守恒和角動量守恒方程為:
式中mi、vi分別為顆粒的質(zhì)量和速度;Fn,ji為顆粒間的法向接觸力;Fτ,ji為切向接觸力; Ii為顆粒的轉(zhuǎn)動慣量; ωi為顆粒的角速度; ri為顆粒半徑,Ffp,i為來自液相對顆粒的作用力;Mi為滾動摩擦力矩。
CFD-DEM 的計算核心在于兩個物理場的數(shù)據(jù)交換,在每個時間步長中,首先CFD 求解器計算連續(xù)相流場,仿真迭代至收斂,然后將網(wǎng)格單元的流體條件傳遞給DEM 求解器,從而計算出作用在顆粒上的阻力,并將此力代入顆粒運動方程,求解離散相的位置、速度等信息。 DEM 迭代計算結(jié)束后,估計出計算單元內(nèi)的孔隙率并連同相間作用力傳遞回CFD 求解器。 CFD求解器利用這些數(shù)據(jù)求解連續(xù)相流場,更新流動區(qū)域,以此循環(huán)進入下一個時間步長。
深海采礦礦漿泵的結(jié)構(gòu)如圖1 所示。
圖1 深海采礦礦漿泵結(jié)構(gòu)示意
依據(jù)我國“十三五”深海采礦重點研發(fā)計劃中揚礦泵的工況參數(shù),按照相似原理,限定二者比轉(zhuǎn)速相同,設(shè)計模型泵進行計算和試驗研究。 模型基本參數(shù)為:流量Q=120 m3/h,揚程H =20 m,轉(zhuǎn)速n =1 450 r/min,輸送濃度5%~15%,輸送顆粒粒徑小于20 mm。 模型水力部件的主要幾何參數(shù)如表1 所示。
表1 礦漿模型泵主要幾何參數(shù)
流場計算域由進口延長段、葉輪、空間導(dǎo)葉、出口延長段組成,整泵豎直放置,進口端處于最下方,漿料垂直提升。 將構(gòu)造的流體計算域三維模型導(dǎo)入ICEM軟件中進行計算網(wǎng)格劃分,得到如圖2 所示的網(wǎng)格單元。 DEM 是一種無網(wǎng)格方法,需要消耗大量的計算資源,在CFD-DEM 耦合計算中70%左右的時間耗費在DEM 求解中。 網(wǎng)格數(shù)過多,每個時間步長內(nèi)CFD-DEM之間都需要進行大量的數(shù)據(jù)交換(如計算連續(xù)相各網(wǎng)格內(nèi)的空隙率和固相平均速度等),耦合模擬非常緩慢;網(wǎng)格單元過小,顆粒將占據(jù)網(wǎng)格單元大部分體積,此時較高的空隙度將導(dǎo)致DEM 求解不穩(wěn)定,所以網(wǎng)格選擇應(yīng)在CFD 應(yīng)用基礎(chǔ)上尋求較粗糙網(wǎng)格為宜。 本文計算模型網(wǎng)格經(jīng)歷了無關(guān)性分析,最終網(wǎng)格單元總數(shù)為54 163。
圖2 流動域三維造型
選取常溫清水作為連續(xù)相,泵入口固相體積濃度為5%。 設(shè)定顆粒為球形,選用5 mm、10 mm、15 mm 3 種粒徑顆粒,顆粒與顆粒、顆粒與泵體間的碰撞采用Hertz-Mindlin 無滑動接觸模型。 泵體材料、顆粒物料的有關(guān)參數(shù)見表2,顆粒與顆粒、顆粒與泵體間的相互影響系數(shù)見表3。
表2 物性參數(shù)
本文運用Fluent 和EDEM 耦合計算對模型泵進行非定常計算,總體時間5 s(約120Tn,Tn 為葉輪旋轉(zhuǎn)一周的時間),流體計算時間步長設(shè)置為0.000 1 s,約為葉輪旋轉(zhuǎn)1°的時間步長。 EDEM-Fluent 耦合計算中顆粒計算時間步長要求不能大于流體計算時間步長,參考瑞利時間步長大小,5 mm 粒徑工況時間步長設(shè)為0.000 05 s,10 mm 和15 mm 粒徑工況時間步長設(shè)為0.000 1 s,顆粒計算時間步長控制在瑞利時間步長的60%以內(nèi)。 湍流模型采用RNG κ-ε 模型。 葉輪設(shè)置于旋轉(zhuǎn)坐標(biāo)系中,導(dǎo)葉處于靜止坐標(biāo)系中。 重力加速度為9.81 m/s2,重力方向與泵進口方向相反。 進口邊界采用速度入口,假定進口處的流動是軸對稱和無預(yù)旋的,出口邊界采用自由出流邊界。 壁面處默認為無滑移邊界條件,近壁處采用標(biāo)準(zhǔn)壁面函數(shù)處理。
顆粒在泵內(nèi)的空間分布、運動特征對于整機的性能,包括正向流動的通暢性和過流壁面的磨損等特性影響很大。 對于給定的流道,相同濃度下,顆粒粒徑越大則顆粒數(shù)越少,反之則越多;顆粒越大,流體拖曳其遷移所需要的動能也將越大;局部濃度偏大將容易導(dǎo)致流道堵塞。
對3 種粒徑下礦漿泵內(nèi)的顆粒情況進行分析,泵內(nèi)的顆粒體積濃度隨運行時間變化情況如圖3 所示。礦漿泵初始工作的一段時間內(nèi),泵內(nèi)顆粒相體積濃度逐漸升高,在運行一段時間后,濃度不再繼續(xù)增加,結(jié)合數(shù)值計算的殘差收斂特征,當(dāng)數(shù)學(xué)收斂和物理收斂均達到要求時,認為計算穩(wěn)定。 由圖3 可知,不同粒徑下的顆粒體積濃度數(shù)值呈現(xiàn)一定的波動,且波動幅度隨顆粒粒徑增大而增大。
對3 種工況在顆粒體積濃度相對穩(wěn)定的時間段內(nèi)(1~5 s)進行分析,得到不同工況下的平均濃度如表4所示。 隨著顆粒粒徑增加,顆粒體積濃度也逐漸升高,但均高于入口處的體積濃度(5%),說明泵內(nèi)顆粒有局部聚集的情況,且隨粒徑增加,聚集更加明顯,但各部分濃度保持基本恒定,說明流道并未發(fā)生堵塞。 圖4 為泵流道內(nèi)各部分的顆粒相對體積濃度,可見相對體積濃度從大到小排序為:進口>葉輪>導(dǎo)葉>出口。
泵內(nèi)各部分的顆??臻g分布如圖5 所示,在進口、導(dǎo)葉及出口段,顆粒體積濃度均隨著顆粒粒徑增大而增高。 在進口延長段,顆粒與清水以相同速度射入,由于重力的作用,顆粒會發(fā)生“沉降”,顆粒越大,沉降速度越大。 顆粒由進口延長段入口到葉輪入口,顆粒相速度下降,且顆粒越大,下降幅度越大。 經(jīng)統(tǒng)計,5 mm,10 mm 和15 mm 粒徑顆粒在該區(qū)域的垂直方向平均速度分別為2.07 m/s,1.97 m/s 和1.89 m/s,均低于該區(qū)域清水流速(2.33 m/s)。 因此,該區(qū)域的顆粒相體積濃度大于入射狀態(tài)的5%,且粒徑越大,濃度越大。
圖5 不同粒徑顆粒工況礦漿泵內(nèi)顆粒分布
在總體體積濃度相同的條件下,小顆粒數(shù)量更多,大顆粒數(shù)量較少(單個顆粒的體積比值為1 ∶8 ∶27),圖6 為葉輪內(nèi)的顆粒分布,它展現(xiàn)出兩個顯著特征:一是葉片進口位置低速顆粒聚集,二是流道中部小顆粒貼近工作面移動,而大顆粒則相對均勻分布。 關(guān)于葉片頭部顆粒聚集效應(yīng),可認為顆粒低速進入流道時,先進入的顆粒與葉輪頭部和蓋板發(fā)生碰撞,反彈后會與后進入的顆粒再度碰撞,碰撞次數(shù)隨著顆粒數(shù)增加而增加,這種聚集和堵塞反過來又對流體起到了加速的作用。 關(guān)于流道中部顆粒的運動和分布,大顆粒在流道中的速度顯著低于小顆粒,表明其獲得的能量也更少,但是其空間分布相對于前者更加均衡。
圖6 葉輪內(nèi)顆粒的空間分布
圖7為導(dǎo)葉流道內(nèi)的顆??臻g分布。 大顆粒工況顆粒體積濃度大于小顆粒體積濃度。 經(jīng)統(tǒng)計,5 mm,10 mm 和15 mm 粒徑顆粒在導(dǎo)葉流道內(nèi)的平均速度分別為9.20 m/s,8.83 m/s 和7.12 m/s。 顆粒在導(dǎo)葉流道內(nèi)沿著導(dǎo)葉葉片移動,外側(cè)為海水介質(zhì),由于顆粒與導(dǎo)葉壁面碰撞摩擦,顆粒速度逐漸減小。 小顆粒呈現(xiàn)多層連續(xù)滑移,大顆粒則多以單層和離散滑移為主。
圖7 導(dǎo)葉內(nèi)顆粒的空間分布
顆粒由導(dǎo)葉以較高速度進入出口延長段,顆粒的速度由于重力和清水的“阻尼”作用,將會逐漸減小,5 mm,10 mm 和15 mm 粒徑顆粒在出口延長段的速度分別為3.22 m/s,3.12 m/s 和3.05 m/s,均高于出口延長段清水的速度(2.89 m/s)。
為了研究不同粒徑顆粒的運動對葉輪和導(dǎo)葉流道內(nèi)流場的影響,監(jiān)測葉輪和導(dǎo)葉流道中間位置處的流動參數(shù)進行分析,同時為詳盡闡述流動機理,特選擇5 mm 的流場,將葉輪按照葉柵周向展開,對比說明葉柵內(nèi)的流動特征。
2.3.1 葉輪內(nèi)的流場特性
圖8 為不同粒徑下的葉輪流道內(nèi)速度和壓力變化曲線,右側(cè)云圖為相對應(yīng)的葉柵內(nèi)流動特征(下同)。不同粒徑的顆粒在葉輪流道內(nèi)流體速度呈現(xiàn)相似的變化趨勢,存在兩個顯著的拐點,出口速度比之進口有一定增加。 實際上,離心泵的葉輪過流面積是逐漸增大的,因此理論上單相介質(zhì)在流道內(nèi)的流速是減小的,造成速度波動的主要原因是流道幾何和局部流動擴散。 圖中2 點位置為葉片進口,由于葉片排擠作用,流速增大,從2 點到5 點流速的下降主要是流道的擴散,5 點以后到出口的流速增加主要是局部的低速區(qū)和流動分離引起的,將圖8 的流線分布與圖6 比對,可知顆粒在工作面附近聚集,形成局部低速區(qū),流體則在外側(cè)加速通過。
圖9為3 種顆粒作用下的流體在葉輪流道內(nèi)的靜壓分布,總體而言其壓力特征為葉片核心區(qū)域壓力均勻上升,后段(5 點以后)出現(xiàn)下降和波折,壓力下降的原因主要是流速增加,而壓力拐點出現(xiàn)在6 點位置處,既與出口截面積擴散有關(guān),也與流道中后部的漩渦有關(guān),漩渦引起的不穩(wěn)定流一定程度上降低了流動的效率。 位置5 點以后小顆粒流體的壓力略高于大顆粒,表明小顆粒群從葉片獲取能量的能力略高于大顆粒,進一步被加速的顆粒又因速度滑移和流固黏性力的作用,將一部分動能轉(zhuǎn)化為流體的能量。
圖9 葉輪流道內(nèi)流線方向上的壓力變化
2.3.2 導(dǎo)葉內(nèi)的流場特性
圖10 為導(dǎo)葉內(nèi)的流場特性,流速整體下降,導(dǎo)葉流道內(nèi)出現(xiàn)了顯著的渦旋,從進口稍后位置一直延伸到中后部,不同流道流動特征存在較大差異。 相對于葉輪內(nèi)葉片的持續(xù)做功,導(dǎo)葉內(nèi)流體是一個動量矩守恒的動壓和靜壓轉(zhuǎn)化的過程,因此理論上單相介質(zhì)在導(dǎo)葉內(nèi)的能量變化特性應(yīng)該是降速增壓的過程。 實際由于設(shè)計階段顆粒過流的考慮,將導(dǎo)葉過流截面做了較大的修正,同時顆粒兩相的內(nèi)部耗散導(dǎo)致靜壓出現(xiàn)了一定的下降。 此外,粒徑的差異也引起了速度下降快慢的變化,10 mm 的顆粒流場流速下降相對均勻,15 mm 顆粒流場流速前段后段均下降較快,5 mm 顆粒流場流速前段平緩后段陡降,拐點位置出現(xiàn)在11~13 點處。
圖10 導(dǎo)葉內(nèi)流線方向上的速度變化
導(dǎo)葉內(nèi)流線方向上的壓力變化見圖11。 可見導(dǎo)葉內(nèi)的壓力變化特征與流速相類似,整體下降,10 mm顆粒的變化幅度小于5 mm 和15 mm 顆粒。 壓力云圖出現(xiàn)了一定的局部高壓區(qū),與拐點區(qū)域相對應(yīng),此外靜壓的周向梯度顯著,壓力面顯著高于流道中間和吸力面區(qū)域,證實顆粒兩相的流動分離。
圖11 導(dǎo)葉內(nèi)流線方向上的壓力變化
2.3.3 葉輪-導(dǎo)葉級間流場特性
深海采礦礦漿泵的單級水力包含動葉柵和靜葉柵,二者的動靜耦合往往伴隨著能量損失和流動不穩(wěn)定。 葉輪-導(dǎo)葉級間流場特性見圖12。 由圖12(a)可知,葉輪導(dǎo)葉銜接處為局部高速區(qū),無論是葉輪段還是導(dǎo)葉段,該處均存在顯著的周向速度梯度,吸力面一側(cè)高速,壓力面低速,導(dǎo)葉進口后端的壓力面一側(cè)漩渦明顯。 圖12(b)顯示,二者過渡區(qū)域總壓變化頻繁,導(dǎo)葉內(nèi)的總壓梯度并不是呈現(xiàn)流動方向,徑向的總壓梯度顯示了固液流動的分層和周向壓力面一側(cè)流動的不穩(wěn)定。 即顆粒對于流動的影響更多地體現(xiàn)在周向的分層(或者過流斷面流動特征),而非流線方向。
圖12 葉輪-導(dǎo)葉級間流場特性
為驗證數(shù)值計算的可靠性,搭建如圖13 所示的礦漿模型泵試驗系統(tǒng)。 試驗時通過調(diào)節(jié)閥門開度來控制泵出口壓力,獲得泵礦漿工作特性。 測試系統(tǒng)包含泵進出口的壓力變送器、管線流量計、轉(zhuǎn)速轉(zhuǎn)矩儀等傳感器、信號采集和數(shù)據(jù)處理系統(tǒng)等,所有的傳感器誤差均在3%以內(nèi),保證試驗數(shù)據(jù)的準(zhǔn)確。
圖13 漿料輸送系統(tǒng)示意
考慮到模型泵的兩級和泄漏效應(yīng),參照經(jīng)驗公式對數(shù)值計算的結(jié)果進行修正[11](修正系數(shù)0.96~0.98),最終數(shù)值計算的泵揚程和效率數(shù)值與試驗對比數(shù)據(jù)如表5 所示。 揚程和效率的數(shù)值模擬相對誤差均在5%以內(nèi),表明數(shù)值模擬結(jié)果較準(zhǔn)確。 礦漿泵輸送漿料時的揚程和效率均低于清水工況,且揚程和效率隨顆粒粒徑增大而下降,這與清水和顆粒之間速度差異而導(dǎo)致的水力損失有關(guān)。 清水和顆粒的密度差導(dǎo)致的顆粒“沉降”效應(yīng),顆粒與顆粒間的碰撞以及顆粒和清水在與泵內(nèi)壁面作用時表現(xiàn)出不同性能等原因均會形成兩者的速度差,且顆粒粒徑越大,速度差越大,由于液體“阻尼”,能量會有損耗,造成泵揚程和效率下降。
表5 數(shù)值計算與試驗結(jié)果對比
1) 基于離散元方法的CFD-DEM 耦合策略能較為準(zhǔn)確地計算大顆粒固液兩相介質(zhì)在離心泵內(nèi)的復(fù)雜流動問題,計算結(jié)果與經(jīng)驗公式計算結(jié)果吻合較好。
2) 顆粒在泵流道內(nèi)產(chǎn)生局部聚集效應(yīng),葉輪和導(dǎo)葉內(nèi)的顆粒體積濃度均高于管線平均值,且顆粒粒徑越大,局部體積濃度越高。 顆粒和流體在葉輪和導(dǎo)葉內(nèi)呈現(xiàn)顯著的分層運動,在葉輪內(nèi)貼近壓力面滑移,而在導(dǎo)葉內(nèi)則偏向吸力面。 小顆粒獲得的能量特性優(yōu)于大顆粒。
3) 顆粒的存在導(dǎo)致了流動分離和流道漩渦存在,過流面積的擴散和固液兩相的內(nèi)部損耗是流道內(nèi)能量損失的主要原因。 隨著粒徑增大,礦漿泵性能下降。
4) 本文研究結(jié)果可為大顆粒固液兩相流的計算和深海采礦礦漿泵的優(yōu)化設(shè)計提供理論指導(dǎo)。