基于渦動力學的離心泵內部流動診斷

王維軍，王 洋 ，黎義斌，陸 榮，劉 敏

(1. 江蘇大學鎮(zhèn)江流體工程裝備技術研究院，江蘇 鎮(zhèn)江 212009；2.江蘇大學流體機械工程技術研究中心，江蘇 鎮(zhèn)江 212013； 3.蘭州理工大學能源與動力工程學院，蘭州 730010)

離心泵是工業(yè)上應用最為廣泛的動力機械之一，特別是電力系統(tǒng)、冶金機械、鋼鐵行業(yè)、建筑行業(yè)、有色、石油化工領域、水處置等領域。離心泵中流體是有渦運動，流動形態(tài)非常復雜，是非定常的和非線性的，泵進、出口管路中存在明顯的有渦流動存在脫流、軸向渦旋、旋轉失速等二次流[1,2]，直接關系著泵運行過程中效率及振動噪聲的高低。近年來，隨著計算流體動力學和計算機軟、硬件的快速發(fā)展，為數(shù)值求解研究離心泵提供了有效的手段。但僅靠單一地提取流場中壓力、速度云圖已經(jīng)不能滿足流場分析，很難定量精確地捕捉到流場中物理量的變化，只能定性描述。

離心泵中渦旋運動和分離流動是常見到的流動現(xiàn)象?？諝鈩恿W家柯齊曼曾說過：“渦是流體運動的肌腱”，北航陸士嘉教授指明“流體的本質就是渦”[3]。渦動力學[3]是計算流體動力學的重要分支，目前發(fā)展已經(jīng)較為成熟。國內、外研究學者采用渦動力學診斷各種三維高度復雜的湍流流動取得了很大的進步[5-8]。在流體力學中，大量的特征都是有精確定義的，然而渦旋至今在數(shù)學上沒有準確的定義，沒有一種數(shù)學方法能精確描述流體渦旋的結構，加之不同流場中渦旋的特性也有不同，仍然沒有一種通用的方法描述渦流結構。在大量實驗和數(shù)值計算復雜流動中研究者逐步形成了不依賴坐標的選擇和旋轉等變化的渦判據(jù)方法：正則化螺旋度、最大渦量法、Q準則、Δ準則、λ2準則、螺旋度法、渦流參數(shù)法等[9]。張翔[10]采用正則化螺旋度法對不銹鋼沖壓焊接離心泵流向渦運動結構和葉輪內的流動分離進行了診斷，結果顯示葉輪流道內的渦旋和流動分離不強烈，泵內流動呈現(xiàn)一種“有序”的結構。小流量工況下，泵內的“有序”的結構將不復存在。曹璞鈺[11]采用渦旋判別Q準則通過研究管道式離心泵吸入式流動不穩(wěn)定發(fā)現(xiàn)，吸入室中的彎管流動和消旋板繞流構成了進口畸變流，存在雙龍卷風式分離渦，該渦旋形成的渦空化引起了離心泵運行噪聲。王洋[12]采用渦方法發(fā)現(xiàn)多級離心泵大流量工況下?lián)P程急劇下降的原因是導葉前緣產(chǎn)生較大的負沖角導致壓力面產(chǎn)生非穩(wěn)定渦旋，加大了導葉內部的損失。

本文針對蝸殼式離心泵運行中壓力脈動、振動過大等現(xiàn)象，采用渦動力學診斷內部流場渦結構和流動分離，研究這些現(xiàn)象的形成機理。

1 數(shù)值計算

1.1 計算模型

離心泵采用懸臂式(端吸)直聯(lián)結構設計，電機與泵同軸，中間取消聯(lián)軸器，減小了泵的軸向尺寸，結構更為緊湊。泵參數(shù)如表1所示，計算模型水力參數(shù)[13]如表2所示。

表1 泵參數(shù)Tab.1 Pump parameters

表2 計算模型水力參數(shù)Tab.2 Calculation model of hydraulic parameters

模型泵流道的三維造型采用Pro/E三維造型軟件根據(jù)葉輪木模圖、蝸殼水力圖造型，確定計算模型。將葉輪進口段和蝸殼出口段分別進行了適當延長，葉輪和蝸殼三維如圖1所示。

圖1 計算模型和子午面Fig.1 Pump model and meridian plane

1.2 湍流模型

1.3 網(wǎng)格劃分與邊界條件設置

網(wǎng)格質量是影響數(shù)值結果精度的重要因素之一，隨著網(wǎng)格數(shù)量增加，數(shù)值預測精度將有所提高。但隨之求解時間也會大大加長，增大了運算成本，綜合考量后本文對葉輪和蝸殼采用非結構四面體網(wǎng)格，泵進、出口延伸段采用六面體結構網(wǎng)格。圖2為設計工況下網(wǎng)格無關性曲線。從圖2可以看出，最小網(wǎng)格為4.0 mm、網(wǎng)格數(shù)484 429時，求解精度最差，耗時3.2 h；最小網(wǎng)格數(shù)為2.2，網(wǎng)格數(shù)1 326 179時，求解耗時5.8 h；隨著網(wǎng)格數(shù)增加揚程逐步下降，當最小網(wǎng)格為1.2 mm、網(wǎng)格數(shù)2 233 499時，耗時7.3 h。最終確定劃分的最小網(wǎng)格數(shù)為1.2 mm，網(wǎng)格總數(shù)為2 233 499，其中葉輪網(wǎng)格數(shù)為883 423，蝸殼網(wǎng)格數(shù)900 476，進口段網(wǎng)格數(shù)為163 600，出口段網(wǎng)格數(shù)為286 000。

圖2 網(wǎng)格無關性假設曲線Fig.2 Grid Independence check curve

本文設置進口邊界條件為壓力進口一個大氣壓(1atm)；出口邊界條件設置為質量流量；固壁條件采用無滑移壁面；旋轉坐標系采用右手定則(旋轉方向為-Z)；葉輪與進口、蝸殼間交互面采用凍結轉子法(Frozen Rotor)；本文采用的是ANSYS-CFX，其基于有限元法的有限體積法對N-S方程進行求解的，在保證了有限體積法的守恒特性的基礎上，吸收了有限元法的數(shù)值精確性；離散過程中的對流項和湍流數(shù)值項均采用高分辨率格式(High Resolution Scheme)；收斂精度設為1×10-5。

2 試驗驗證

圖3為江蘇大學流體中心離心泵閉式實驗臺示意圖，實驗臺精度滿足GB/T 3216-2016《回轉動力泵水力性能驗收試驗1級和2級》中1級精度要求。圖4為0.6Qd、0.8Qd、1.0Qd、1.2Qd和1.4Qd5個工況下的數(shù)值計算與實驗測得的外特性對比值。從圖4可以看出，5個工況的數(shù)值計算的效率值高于實驗結果，誤差依次為：10.86%、3.87%、3.38%、3.81%和3.12%，除了0.6Qd工況下誤差較高之外，其余工況數(shù)值計算的效率較精確；揚程誤差分別為：1.81%、0.06%、3.21%、3.24%和1.30%；功率誤差較大，這是因為在數(shù)值模擬過程忽略了軸承、摩擦副等引起的機械損失的原因。數(shù)值計算外特性變化趨勢和實驗值是一致的，最高效率點流量也大致相當，CFD計算的結果是基本準確。

圖3 閉式實驗臺Fig.3 Closed testing bed

圖4 數(shù)值計算和實驗的外特性曲線Fig.4 Results of numerical simulation and test

3 流動診斷

3.1 常規(guī)診斷

常規(guī)診斷是數(shù)值計算后處理中提取壓力、速度云圖來描述流場變化，定性描述流體運動狀態(tài)，觀察葉輪做功。圖5為設計工況下不同截面上的壓力云圖，可以看出，隨著半徑的增大，葉輪內壓力逐步增大，葉輪中的壓力呈軸對稱分布，葉輪進口為最低壓力區(qū)，出口存在明顯的射流尾跡區(qū)；蝸殼內部的壓力值范圍在315～430 kPa，靠近隔舌位置的壓力最大，并延伸至蝸殼出口；從blade to blade 截面可以看出，同一半徑處葉片壓力面的壓力值高于吸力面的，葉輪低壓區(qū)位于葉片頭部吸力面處。

圖5 設計工況下不同截面上的壓力云圖Fig.5 Pressure clouds on different sections at design conditions

圖6為設計工況下Z=0截面上的速度云圖和速度矢量圖，可以觀察到，葉片壓力面上存在明顯的低速區(qū)，該處的速度基本上在0 m/s，稱為流動滯止區(qū)，范圍占據(jù)葉片壓力面的一半，有研究表明通過縫隙引射流技術可以有效控制該處的流動狀態(tài)，從而提高葉輪做功能力[16,17]；同一半徑處吸力面上的速度值大于壓力面的，尤其在靠近葉片出口位置區(qū)；葉輪出口存在5處高速區(qū)，與葉片數(shù)吻合，該處的最大速度值為18.77 m/s。

圖6 設計工況下Z=0截面上的速度云圖Fig.6 Velocity clouds on Z=0 section under design conditions

3.2 渦動力學診斷

渦包括渦量場(Vorticity)和渦旋(vortex)，渦量場是指渦量的空間分布，渦旋是集中渦，渦是湍流的一種基本結構[3]。渦量的來源是流場中存在速度梯度，是描述有旋運動的一個運動學物理量，渦量場在離心泵中會演變成一個個離散渦量聚集的渦旋。在黏性流體力學中所述渦是造成液流能量損失的主要原因之一，研究離心泵中的黏性流動，特別是二次流動，就離不開研究渦旋的運動。一方面可以弄清楚渦的形成、發(fā)展和衰亡與渦旋之間的流動規(guī)律和相互作用機理，另一方面可以研究如何采用可行的控制手段對渦旋起主導作用的流動實施有效的控制。尤其對離心泵而言，弄清不同部件內部的流動結構尤為重要，可以從根本上控制流動分離。

3.2.1 |ω|判據(jù)

渦是渦量集中的區(qū)域，因此自然可以用渦量的模作為渦的判據(jù)，渦量極大的地方是渦的中心。根據(jù)渦量ω的數(shù)學定義，可以理解為流體微團繞其中心做剛性旋轉的角速度之兩倍[4]

ω=▽×V

(1)

式中：V為流體質點的速度矢量；ω為流體質點的渦量。

從圖7可以看出，在葉片進口吸力面位置的渦量模值|ω|最大，隨著半徑的增大，逐步減??；蝸型段內的值在1.72～214.46 s-1之間，無法有效捕捉渦的運動狀態(tài)。因此采用渦量模值|ω|無法判斷離心泵內渦結構、捕捉渦特征。主要體現(xiàn)在兩個方面：其一，渦量模值極大不等于有渦存在；其二，不具有廣義伽利略不變性[18]。

圖7 渦量的模云圖Fig.7 The value of

3.2.2 正則化螺旋度

螺旋度[19](Helicity)是度量湍流渦量場拓撲結構的一個重要的物理量。離心泵三維流場中流體運動的螺旋度的定義為：

(2)

被積函數(shù)稱為螺旋度密度。正則化螺旋度(Hn)根據(jù)速度矢量與渦量的夾角判斷渦核，這種方法可以捕捉到渦核位置，Hn定義為速度與渦量的點積除以速度的模與渦量的模的乘積，如下式：

(3)

在流場區(qū)域內，它是在除速度矢量V與渦量矢量ω為零的特殊點外均有定義的標量場，其值在[-1，1]區(qū)間。在渦核區(qū)域，速度矢量方向與渦量矢量方向近于平行，正則化螺旋度Hn趨于±1。正則化螺旋度Hn的符號表明渦旋轉的方向，以流動方向為正方向，Hn為正，渦旋轉方向為逆時針方向；Hn為負，渦旋轉方向為順時針方向。圖8為Z=0截面上的正則化螺旋度分布。從圖8可以看出，葉輪吸力面靠近頭部位置存在旋轉方向與流動相反的渦旋，Hn在[-1，0]；靠近葉輪出口位置吸力面上存在和流動方向相一致的渦旋，Hn在[0，1]，圖9為葉刪局部放大圖尤為明顯，該現(xiàn)象與文獻[10]獲得的結果一致，葉片吸力面上存在一對旋轉方向相反的渦旋，且隨著葉輪旋轉周期的生成與潰滅；葉輪和蝸殼的交界面上的Hn值近似對稱(見圖10)。

圖8 離心泵Z=0截面與葉片、交界面正則化螺旋度Fig.8 Hn of blades, interface and plane of Z=0

圖9 葉柵流道與局部放大圖Fig.9 Hn of cascade and partial enlarged drawing

圖10 葉輪一個周期內兩個時刻典型的渦核變化[10]Fig.10 Vetex core in different time under one cycle

圖11和圖12為蝸殼4個斷面上正則化螺旋度與截面流線分布?？梢钥闯鑫仛嗝嫔洗嬖谝粚πD方向相反的渦旋，在Plane1和Plane3上表現(xiàn)尤為明顯，基本呈對稱分布，在該對渦旋的結合部分有明顯的相互融合的趨勢；Plane2上流線和其他截面一樣較為對稱，但正則化螺旋度提取的渦核分布比較雜亂，呈現(xiàn)4個正向和1個反向渦核聚集區(qū)，渦旋有明顯的破裂；Plane1靠近蝸殼隔舌，該區(qū)域的渦核聚集不明顯。

圖11 蝸殼斷面上的正則化螺旋度與流線Fig.11 Hn and velocity curves of section in volute

圖12 放大圖Fig.12 Enlarged drawing

3.2.3Q判據(jù)

由二階張量特性[3]可知，離心泵不可壓縮流的局部速度梯度張量▽V的特征方程可寫為：

λ3+Qλ-R=0

(4)

如果λ1、λ2、λ3是它的3個根，它們之間存在3個獨立的不變量：

P=λ1+λ2+λ3=divV=0

(5)

(6)

(7)

式中：Ωij分別為應變速率張量和渦量張量。

(8)

渦旋有3種定義方式，分別是Q判據(jù)，Δ判據(jù)和 判據(jù)。Hunt等人提出把Q>0的區(qū)域定義為渦，即意味著‖Ω‖2>‖E‖2，即在離心泵渦旋的區(qū)域內流體的旋轉(渦量大小)起主導地位，而流體的應變率大小次之，這種方法稱之為Q判據(jù)。圖13為蝸殼內的渦旋Q等值面。當Q=140 626 s-2時，葉輪出口位置脫落的一對渦旋進入了蝸殼，且呈軸對稱分布，如圖13(a)所示；隨著Q值減小到58 533 s-2時，蝸殼蝸形段出現(xiàn)了一對長度較大的渦旋，可以預知隨著葉輪旋轉，該對渦旋將周期性的拍打隔舌，形成明顯的壓力脈動[20]，也就是說蝸殼式離心泵產(chǎn)生壓力脈動的根源是葉輪出口位置產(chǎn)生周期性的渦旋導致的，隔舌位置的渦旋強度最大。因此，有研究者通過改變隔舌來改善離心泵壓力脈動強度是符合渦動力學診斷規(guī)律，如賈程麗[21]研究隔舌倒圓，祝磊[22]采用階梯隔舌等。

圖13 蝸殼蝸形段內渦結構Fig.13 Votex structures of volute spiral contour

圖14為Q=21 072 s-2蝸殼擴散段渦結構。以間隔20 mm做不同的截面，通過渦動力學診斷發(fā)現(xiàn)，蝸殼擴散段上靠近隔舌位置有兩對渦旋強度較大的區(qū)域，主渦來自葉輪上脫落后進入蝸殼的旋轉方向相反的一對渦，這對旋轉方向相反的主渦比較穩(wěn)定，渦心一端在葉輪上，一端在蝸殼出口上；次渦是蝸殼隔舌上脫落的一對類卡門渦街的結構，沿著蝸殼出口位置次渦逐步和主渦融合； 次渦不穩(wěn)定，在各種外力的作用下被拉伸，渦有明顯的撕裂，渦心變成兩個區(qū)域；從圖14(c)可以看出，沿著蝸殼出口方向，渦量的模值由規(guī)則的對稱分布逐步呈現(xiàn)紊亂的結構，這是由于正反渦在相互作用下融合的結果，出口位置的渦強度明顯要若于靠近隔舌位置的。

圖14 蝸殼擴散段渦結構Fig.14 Votex structures of volute diffusion section

3.2.4 葉輪表面流動分離診斷

上文可知，葉輪和蝸殼內均出現(xiàn)了一對旋轉方向相反的主渦，葉輪出口位置脫落的主渦是離心泵壓力脈動的根源。在渦量動力學中邊界渦量[3](Boun- dary Vorticity Flux，BVF)是其核心的概念，描述單位時間內通過單位面積擴散的渦通量的物理量，即：

(9)

式中：n為壁面外法線單位矢量；μ為動力黏性系數(shù)。

通常認為BVF的峰值會引起下游邊界渦量的峰值，三維黏性可壓縮有旋運動：邊界渦量σ包含固體壁面加速度產(chǎn)生的σa、體積力產(chǎn)生的σf、切向應力產(chǎn)生的στ和法向應力產(chǎn)生的σπ四部分。對于離心泵而言，葉片是勻速轉動的無滑移壁面，體積力通?？梢院雎圆挥?，且流動屬于完全湍流，σπ遠大于στ。因此，只需考慮由法向應力產(chǎn)生的邊界渦量。對于離心泵而言，法向應力等于壓力，即：

▽P

(10)

在離心泵中作用在葉片表面的力矩主要是由壓力梯度提供的，根據(jù)導數(shù)矩理論可得到：

(11)

式中：Mz為流體作用離心泵葉輪上的扭矩；S為面積；Z為旋轉軸方向；σpz為軸向壓力梯度引起的邊界渦量流。

式(11)右邊第一項表示BVF矩在葉片上的面積分，第二項表示壓力矩在邊界上的線積分。式(11)建立了離心泵葉輪做功與葉片表面上BVF之間的關系。根據(jù)BVF的分布狀況可以找到對離心泵葉輪做功起到負貢獻的區(qū)域。從圖15葉輪表面BVF分布可以看出，BVF極大值位于葉片頭部吸力面，該處是葉片與流體剛剛接觸位置，此處流動不穩(wěn)定，葉片進口安放角可能設計不合理導致BVF出現(xiàn)了峰值，葉輪出口位置同樣存在如此現(xiàn)象。

圖15 葉片表面BVF分布Fig.15 BVF on the blades

離心泵葉輪分離區(qū)特性[23]：不僅表面摩擦力線匯聚，而且渦量線有大的曲率。分離線判據(jù)：渦量線曲率達到局部最大值。分離先兆：BVF線折向基本沿表面摩擦力線方向，數(shù)學表達式為：

(12)

τ=μωn

(13)

從式(12)可以看出，雷諾數(shù)小則BVF線和摩擦力線夾角大。吳介之等基于渦動力學理論[23]，提出了任意曲面上三維流動分離和分離渦形成的新判據(jù)：BVF線和摩擦力線方向趨于一致就可以判定流動分離了，表面摩擦力線匯聚，并且BVF有峰值。以下根據(jù)一組正交的曲線簇(τ，ω)對葉片表面的流動分離進行診斷。

從圖16可以看出，葉片吸力面靠近蓋板位置摩擦力線發(fā)生了匯聚，有流動分離現(xiàn)象；葉片壓力面表面摩擦力線在葉片中間位置形成了匯聚，并且在摩擦力線開始匯聚時渦量線表現(xiàn)出大曲率，并延伸至葉片出口位置，這與前文采用正則化螺旋度診斷獲得的結論一致；因此可以判斷葉片表面黏性流動分離的起點位于葉片骨線1/2處，并延伸至葉片出口；渦動力學診斷方法可以準確捉到葉片表面的流動分離區(qū)，可以為優(yōu)化葉輪提供一種有效的方法。

圖16 葉片表面流動分離結構(紅色為τ線，黑線為ω)Fig.16 Flow separation structure on blades

4 結 語

(1)通過采用常規(guī)流動診斷發(fā)現(xiàn)，設計工況下葉輪出口位置存在射流尾跡區(qū)，葉片壓力面存在一個明顯的低速區(qū)，其尺度約為葉片骨線的一半。

(2)渦動力學診斷可以準確捕捉到離心泵葉輪、蝸殼中的渦結構的變化，蝸殼蝸形段內存在一對有葉輪表面分離出來的旋轉方向相反的渦旋，該渦旋是蝸殼式離心泵脈動的根源，數(shù)量和葉片數(shù)一致。

(3)蝸殼擴散段存在兩對渦，一對是蝸型段的主渦，一對是隔舌位置脫落的次渦。沿著蝸殼出口的方向，兩對渦逐步融合，強度逐步減弱。

(4)采用BVF邊界渦量法可以獲得葉片表面流動分離區(qū)，為采用渦動力學方法控制離心泵葉片表面流動分離提供了一條有效的思路。