基于熵產理論的超低揚程雙向臥式軸流泵裝置飛逸特性

許 哲，鄭 源,，闞 闞，黃佳程

（1. 河海大學水利水電學院，南京 210098；2. 河海大學能源與電氣學院，南京 210098）

0 引 言

臥式軸流泵具有結構緊湊、安裝及檢修方便等優點，多用于平原地區的低揚程泵站，以滿足灌溉、排水、防洪及調水等需求。近年來，隨著南水北調東線工程穩步推進，泵站運行的安全穩定性得到越來越多的關注。當泵站機組遭遇斷電事故，若出水流道側閘門拒動而不能截斷水流，水泵機組會進入水流倒流、葉片反轉的水輪機工況。當葉輪轉速增至最大轉速時，水泵處于穩定飛逸狀態，此時葉輪的最大轉速稱為飛逸轉速。在飛逸過渡過程中，機組流道內將出現流態復雜、大尺度漩渦流動等不穩定現象，極易引發流場內部劇烈的壓力脈動和周期性尾水渦帶，威脅到機組的穩定運行。為了保證泵站的穩定運行和電站設備的安全，需要對機組的瞬態過程開展研究。

目前，針對水泵系統瞬態過程的研究，主要是基于一維特征線法（Method of Characteristics，MOC）的數值解法，優點是計算迅速和外特性模擬準確。在理論創新方面，Afshar等[1]提出了隱式特征線法，對管路系統中的閥門、蓄水池和泵等設備進行單元式定義，并推導相應的控制方程，與顯示特征線法相比，計算得到的揚程與流量變化更加準確；Rohani等[2]基于隱式特征線法，提出了點隱式特征線法，進一步簡化了計算量。在工程應用方面，葛強等[3]研究了燈泡貫流泵站在不同電機過載系數及不同葉片角度下的起動過程；余國鋒等[4]研究了水輪發電機組導葉按線性規律開啟的過渡過程；黃晨等[5]研究了三機式抽水蓄能電站在不同調速器轉速整定值下的停泵過程。一維數值模擬能夠快速準確地得到水泵系統的外特性參數變化過程，但不能反映機組內部流態的瞬態特性。近年來，隨著現代數值方法的快速發展，計算流體動力學（Computational Fluid Dynamics, CFD）已成為水力機械過渡過程數值模擬的有力工具[6-8]。李琪飛等[9]對水泵水輪機飛逸過程進行模擬，發現無葉區高速水環是導致隨機脈動與能量損失的原因；Fortin等[10]對水輪機飛逸過程中葉片上的壓力進行模擬，并通過模型試驗對壓力結果進行驗證；羅興锜等[11]分析了燈泡貫流泵飛逸時葉片進口相對液流角與速度矩的演變規律；陳會向等[12]研究了槳葉調節方式對軸流式水輪機甩負荷特性的影響。

在分析機組運行過程中的水力損失時，常通過效率公式來間接評估總水力損失，不能直觀判斷不同部位能量損失的來源及具體分布。與此同時，當機組處于瞬態過程中時，各參數尚在變化之中，因此常規的計算方法不能得到準確的水力損失。如今越來越多的學者開始關注熵產與水力損失的聯系，即通過計算內能耗散來評估水力損失，但目前熵產理論在水力機械方面的應用仍然較少。Kock等[13]基于雷諾時均方程首次提出熵產率的計算公式；張翔等[14]提出了壁面摩擦損失的計算方法；Yu等[15]推導了空化熵產率的公式。隨著熵產理論的發展，越來越多的學者開始運用熵產理論來分析流動過程中的能量損失。在優化設計方面，王威等[16-18]運用熵產理論對翼型、離心風機及軸流式反向發電機組進行了優化設計，提高了能量的利用效率。在飛逸過程的模擬方面，Zhou等[19-20]基于熵產理論對抽水蓄能機組及離心泵的斷電飛逸特性進行分析。以上結果均表明，熵產理論可定量描述各過流部件內的能量損失。

以往水力機械飛逸過程的研究中，研究對象大多局限于水力機械裝置本身，未涉及上下游水域的演變過程；與此同時，研究方法也少有涉及熵產理論，少有從能量損失的角度來闡述機組的瞬態變化過程。超低揚程的雙向泵具有雙向輸水能力，同時也具有雙向飛逸的風險，而目前針對雙向飛逸過渡過程的研究較少。本文采用熵產理論分析機組飛逸過程的演變規律，并對比正反向飛逸工況下能量損失；并通過模型試驗來驗證三維模擬及熵產理論的準確性；最后，以正向飛逸工況下葉片展開面處漩渦與熵產率的分布為例，分析能量損失與流場漩渦的關聯，以期為機組安全穩定運行提供參考。

1 數值計算方法

1.1 控制方程與湍流模型

臥式軸流泵裝置內水流的流動受質量守恒定律和動量守恒定律的約束，為簡化計算，本文采用雷諾時均化的方程形式。上下游計算域自由液面附近的水和空氣并不混摻，故可用流體體積函數（Volume of Fluid，VOF）來確定自由液面的位置及計算水氣兩相的體積分數，控制方程[6,11]如下：

連續性方程：

動量方程：

VOF方程：

剪切壓力傳輸（Shear-Stress Transport，SST）k-ω[12]湍流模型修正了湍流黏度公式，可以更好地傳遞壁面處的剪切應力，有助于預測近壁區的流動，因此選擇SSTk-ω湍流模型來封閉控制方程。

1.2 熵產理論

根據熱力學第二定律，熵的微增量總是大于0。通過熵產理論可以對水泵飛逸過程中能量的耗散進行量化評估，如回流現象、壓力脈動、尾水渦帶、流動分離等會加劇能量的損耗。在機組飛逸的過程中，考慮到水的比容較大，可以認為機組飛逸過程中水的溫度保持不變，因此不考慮傳熱引起的熵產。基于雷諾時均的湍流運動，熵產率（Entropy Production Rate, EPR）主要包含3項：1）由時均速度引起的單位體積熵產率，稱為直接耗散項；2）由湍流脈動引起的單位體積熵產率，稱為湍流耗散項；3）由壁面摩擦損失引起的單位面積熵產率，稱為壁面耗散項。其中，單位體積總熵產率S˙D′[13]可表示為

式中tμ表示湍流動力黏度，Pa·s。

本文采用雷諾時均方法進行數值模擬，無法直接計算湍流耗散項。根據Kock等[13]提出的思想，當雷諾數趨向于無窮大時，可以將SSTk-ω湍流模型中的ω與脈動速度分量產生的熵產率關聯起來，即

式中經驗系數β通過直接數值模擬率定得到，β=0.09[13,21]；k表示湍動能，m2/s2；ω表示湍流渦黏頻率，s-1。

根據張翔等[14]提出的壁面摩擦損失的計算方法，壁面耗散熵產WS可通過積分得到：

式中τw表示壁面剪切應力，Pa；uw表示壁面區第一層網格中心速度矢量，m/s；A表示計算域表面積，m2。

計算區域內由時均速度與脈動速度引起的熵產亦可通過積分獲取，即

式中SD表示直接耗散熵產，W/K；S D′表示湍流耗散熵產，W/K；V表示計算域體積，m3。

各計算域內的總熵產S(W/K)為

1.3 算法實現

當水泵處于斷電飛逸狀態時，葉輪的轉速受自身水力扭矩的控制，因此使用Fluent軟件的自定義函數功能（User defined function, UDF）來控制葉輪的轉速，力矩平衡方程[11]如下：

式中M是葉輪的總力矩，N·m，可通過UDF中的Compute_Force_And_Moment語句進行實時獲取，機械摩擦力矩和轉子風阻力矩較小，在總力矩中暫不考慮；J是轉動慣量，kg·m2；n是轉速，r/min。任意時刻下的葉輪轉速離散化求解公式為

式中Mr表示第r時刻的總力矩，N·m；nr和nr+1分別表示第r和r+1時刻的轉速，r/min；Δt表示第r和r+1時刻之間的間隔，s。

2 數值計算模型

2.1 幾何模型

本文研究對象為超低揚程雙向臥式軸流泵全過流系統，如圖1所示。

泵段采用SZM35水力模型。葉輪直徑為1.6 m，葉片數為4，導葉數為5，葉片安放角為-4°。在實際工程運行中，不同調水任務下的設計參數不同，正向工況下額定轉速、設計流量及設計揚程分別為170 r/min、5.0 m3/s及0.91 m；反向工況下額定轉速、設計流量及揚程分別為170 r/min、4.5 m3/s及1.03 m。這些參數均由實際工程需要確定。數值模擬的計算域包括葉輪、導葉、進出水流道及上下游計算域。為了保證模擬的真實有效性，以實際工程中泵裝置的間隔距離、設計水位高度與開挖高度來確定上下游區域的寬度、自由液面與底面的位置等。

2.2 工況設置

本文考慮2種不同上下游水位工況下的飛逸過渡過程，正反向飛逸工況下的初始計算參數如表1所示。對于正向飛逸工況（Forward Runaway Condition, FRC），臥式軸流泵全過流系統的初始流場如圖2a所示，水流依次流過直錐型進水流道、葉輪、導葉及S型出水流道；對于反向飛逸工況（Backward Runaway Condition, BRC），初始流場如圖2b所示，水流依次流過S型進水流道、導葉、葉輪及直錐型出水流道。

表1 初始計算參數Table 1 Initial simulation parameters

在上述2種工況中，上下游計算域的進出口條件均為壓力邊界條件，并采用UDF功能定義壓力邊界條件，即壓力沿水深變化，而非定值。采用Pressure-linked equation-consistent（SIMPLEC）算法來求解壓力項和速度項，采用一階隱式格式離散時間項，采用二階迎風格式離散對流項和擴散項，時間步長設為0.001 s。

2.3 網格劃分

本文采用ICEM軟件對整個計算域進行結構化六面體網格劃分，利用O形網格劃分進出水流道，并對壁面處網格進行加密。y+值為第一層網格質心到壁面的無量綱距離[21]，以此可判別近壁區的節點分布是否合理。當壁面處y+值小于1.5時，可滿足SSTk-ω湍流模型中低雷諾數k-ω在近壁處的要求。為了驗證網格的獨立性，采用網格收斂指數[22-23]（Grid Convergence Index，GCI）來評估網格方案引起的數值誤差。本文選取3種由疏至密的網格方案，網格數量分別為547萬、770萬和984萬，依據表1中的流量與轉速等參數進行正反向工況的定常模擬，并選取葉輪域的總熵產值參與獨立性網格驗證。

網格的獨立性驗證結果如表2所示，其中網格細化因子1與2分別為密中及中疏網格數量的相對值的立方根，數值解1～3分別對應網格由密至疏時的葉輪域總熵產值，當網格數量為984萬時，正反向工況下細網格收斂指數分別為0.55%與0.21%，均小于3%，說明方案3（984萬網格）的網格精度滿足要求。圖3為網格劃分示意圖，總網格數為984萬。

表2 網格獨立性驗證Table 2 Grid independence verification

2.4 模型驗證試驗

模型試驗是驗證數值模擬結果可靠性的重要方式之一，河海大學高精度水力機械多功能試驗臺如圖4所示。進行雙向水泵工況試驗時，流量從大到小依次測試至少15個工況點，測量點應合理分布在整個能量曲線上，2個工況點之間應保持足夠的穩定時間，待試驗數據穩定后，測得流量、揚程、轉速和軸功率等試驗數據。

在正反向不同流量工況下對水泵進行數值模擬，并計算單位流量Q11與單位轉速n11。其中，單位流量Q11與單位轉速n11公式[9]如下：

3 結果與分析

3.1 模型驗證結果

將模擬的外特性曲線與試驗值比較（圖5a），正向工況下試驗值與模擬值的誤差小于2%，反向工況下的誤差則在3%以內，說明模擬方法準確可靠。與此同時，將熵產計算得到的揚程損失[21]與模型試驗的揚程損失進行對比（圖5b），正向工況下試驗值與模擬值的誤差小于2.5%，反向工況下的誤差則在3%以內，從而驗證了熵產計算的準確性。

3.2 外特性規律

在對水泵進行正反向飛逸工況的模擬之前，先進行30個周期的非定常計算，保留最后5 s的計算結果作為飛逸數值計算的初始場，并將外特性參數進行無量綱處理，如圖6所示。正向飛逸工況（FRC）的初始流量為5.86 m3/s，轉速為170 r/min，扭矩為7.24 kN·m。反向飛逸工況（BRC）的初始流量為5.08 m3/s，轉速為170 r/min，扭矩為7.93 kN·m。

根據機組流量、轉速及扭矩的變化情況，水泵在飛逸過渡過程中將經歷水泵狀態、制動狀態、水輪機狀態以及飛逸狀態（圖6）。在飛逸發生之前（0～5 s），機組保持在水泵工況運行，此時外特性參數分別與對應初始值保持一致，正反向飛逸工況的初始揚程均為0.91 m，而正向飛逸工況的初始流量更大且初始扭矩更小，可見該水泵在進入正向飛逸前作泵工況時水力損失更小，這與試驗結果相一致，即同一單位流量下正向飛逸工況的揚程損失更小（圖5b）；當機組剛斷電后（5～10.8 s），葉輪不再對水流提供動力矩，扭矩、轉速、流量等參數值持續減小；隨后機組進入制動狀態（10.8～11.7 s），水流開始倒流，葉輪轉速方向保持泵工況旋轉方向但轉速值逐漸趨向于0；之后機組進入水輪機狀態（11.7～27.5 s），葉輪在水流的扭矩作用下開始反向旋轉，且轉速隨流量的增大而增大，此時正反向飛逸工況下轉速的差值開始增大；當機組進入穩定飛逸狀態時（27.5～42.5 s），外特性參數值基本保持穩定，正向飛逸工況下的流量為10.30 m3/s，轉速為271.6 r/min，而反向飛逸工況下的流量為9.32 m3/s，轉速為220.5 r/min；可見正向飛逸工況的流量與轉速均高于反向飛逸工況，此時2種工況下的扭矩值均在零值附近波動，但正向飛逸工況的扭矩波動的幅值更大。

為了進一步分析葉輪扭矩的波動特性，采用短時傅立葉變換方法（Short-time Fourier Transform，SFFT）對無量綱扭矩值進行處理，窗函數設置為漢寧窗函數，以獲取更加準確的脈動頻率，如圖7所示。由于葉輪扭矩是其正反表面壓力繞軸的積分量，說明扭矩的波動將受葉片表面壓力脈動的影響，而葉輪室內壓力脈動現象與動靜干涉相關[24-25]，使得扭矩的波動頻率受葉頻（Blade Passing Frequency,fBPF）控制。當機組處于水輪機狀態或飛逸狀態時，正向飛逸工況的脈動頻率以葉頻及其高次諧波為主，反向飛逸工況的脈動主頻為葉頻，且正向飛逸工況扭矩波動幅值更高，可見正向飛逸工況下動靜干涉作用更強烈。

3.3 熵產演變規律

圖8 a所示為總計算域的3項熵產的變化曲線。其中直接耗散熵產占據主導地位，湍流耗散熵產SD′次之，壁面耗散熵產WS最少。根據熵產的計算公式可知，直接耗散熵產與速度梯度密切相關，簡化后計算的湍流耗散熵產SD′受湍動能控制，而壁面耗散熵產WS與近壁區水流流速相關。因此機組內流量增大時，流道內的速度梯度、湍動能以及近壁區的流速在宏觀上相應增大，使得總計算域3項熵產值與流量的變化規律皆為先減少后增加。對比正反向飛逸工況的熵產值，在斷電飛逸前，反向飛逸工況的各項熵產值較大，表明反向泵工況運行時的能量損失較大；在最終飛逸狀態下，正向飛逸工況的各項熵產值較大，這可能是因為流態的惡化導致更多的水力能量損失。

水泵在飛逸過渡過程中裝置內各過流部件及上下游計算域的總熵產值變化情況分別如圖8b與圖8c所示。當機組處于穩定飛逸狀態時，正反向飛逸工況中熵產最大的過流部件均為葉輪，其次為進水流道。這是因為高速旋轉的葉輪內速度梯度的量值較大，同時強烈的動靜干涉作用提高了湍動能的量值，且葉輪域具有強壁面效應，因此葉輪域的總熵產值最大；而高速旋流自葉輪室延伸至進水流道，繼而使得進水流道的總熵產值僅次于葉輪室。泵裝置正反向飛逸過程中，下游計算域由泵工況的穩定來流域轉變為飛逸工況的出流區域，出流域受飛逸工況復雜的出水流態影響，因此下游熵產值隨飛逸過程的發展逐漸增大；上游計算域同理，由泵工況的出流域轉變為飛逸工況的穩定均勻來流域，因此熵產值隨時間由大變小。同時飛逸前后對比可見，飛逸工況出水區域相比于泵工況出水區域有著明顯更大的熵產值，代表飛逸工況的復雜出流帶來了更大的能量損失。

圖9 所示分別展示了葉輪、導葉、進水流道、出水流道、下游及上游計算域在飛逸過渡過程中三項熵產值的變化情況。由圖9a可知，在機組飛逸后，葉輪域內直接耗散熵產占比最多，壁面耗散熵產SW略高于湍流耗散熵產SD′，可見葉輪域的強壁面效應產生了明顯的能量損耗；同時正向飛逸工況下湍流耗散熵產SD′大于反向飛逸工況，這源于正向飛逸工況下更為強烈的動靜干涉作用（圖7）。由圖9b可知，正向飛逸前機組尚處于水泵狀態時，導葉位于葉輪的出流方向，而反向飛逸前導葉位于葉輪的入流側，因此正向飛逸前處于水泵工況下的導葉內復雜的流態致使各項熵產值更大；當機組處于水輪機狀態或飛逸狀態時，正向飛逸工況下導葉位于葉輪的入流方向，導葉內相對平順的流態使得各項熵產值較小。由圖9c與9d可知，正反向飛逸工況下壁面耗散熵產WS遠小于其他2項熵產，即進出水流道中能量損失主要為直接耗散熵產和湍流耗散熵產，近壁區的能量損耗較小。

由圖9e與9f可知，正反向飛逸工況下游計算域的湍流耗散熵產SD′占比均高于直接耗散熵產，可見下游計算域流速相對均勻，且水流的湍動效應占據主導；上下游計算域內壁面（計算域底部）耗散熵產WS均近乎為0，占比小于總熵產的0.1%，此時壁面耗散熵產的影響較小。

3.4 全流域流態分布

圖10 所示為正反向飛逸工況下全流域子午面處的流線圖。在t0時，正反向飛逸工況下進水流道內速度的跨度值分別為4.07與4.22 m/s，正反向飛逸工況下出水流道內速度的跨度值分別為4.47與5.47 m/s，可見反向飛逸工況下進出水流道內速度的跨度值均相應大于正向飛逸工況下的跨度值，這在一定程度上說明反向飛逸工況下的速度梯度較大，這也是反向飛逸工況下直接耗散熵產大于正向飛逸工況的原因（圖9c與9d）。在tQ=0時，進出水流道內流態嚴重失穩，存在強烈的漩渦和回流，但此時流速較低，因而速度梯度與湍動能的量值較小，使得正反向飛逸工況下進出水流道的總熵產值較小（圖8b）。在tn=0時，機組已處于倒流狀態，此時流道內的流態惡化的程度稍有緩解。在tM=0時，進出水流道內流速近似穩定，同時上游計算域內的流態均已恢復平順狀態，因此上游計算域內熵產值逐漸趨0（圖9f）。在t=42.5 s時，機組穩定在飛逸狀態，正反向飛逸工況下游域中速度的跨度值分別為2.05 m/s與3.20 m/s，可見反向飛逸工況下速度的跨度值較大，這可能是反向飛逸工況下游域的直接耗散熵產S

D高于正向飛逸工況的原因（圖9e）。

3.5 葉輪室渦核分布

圖11 所示為tQ=0時正向飛逸工況下渦核與熵產率（EPR）不同圓周截面分布的對比圖，并采用特征值Q作為渦識別準則。正向飛逸工況下渦核均聚集于泵工況的葉輪葉片進水邊，且靠近輪轂側葉片展開面處（span=0.05）脫落渦呈聚集態，而臨近外殼側葉輪展開面處（span=0.95）脫落渦呈離散態，可見渦核沿徑向呈現逐漸脫離葉片的趨勢。與此同時，渦核的出現意味著水流在葉片進水邊形成了強度較高的漩渦，較大的速度梯度則產生了較大熵產率，不同截面熵產率與渦核分布較高程度對應相似，這不僅表示渦核聚集之處產生了明顯的能量損失，也說明漩渦的產生與演變是熵產率較大和能量損失的重要原因。

4 結 論

本文以超低揚程下臥式軸流泵全過流系統為研究對象，對比分析了正反向飛逸工況下外特性和熵產的演變規律，并從流態與渦核分布著手分析能量損失的特點與原因。主要研究結論如下：

1）在高精度水力機械試驗臺上進行雙向泵的模型試驗，一方面對比模擬與試驗得到單位轉速值，另一方面對比熵產計算與試驗得到的揚程損失，結果說明模擬方法與熵產理論準確可靠。

2）在機組正反向飛逸過程中，葉輪扭矩的波動頻率受葉頻控制。正向飛逸工況的脈動頻率以葉頻及其高次諧波為主，反向飛逸工況的脈動主頻為葉頻，且正向飛逸工況下扭矩波動的幅值更高。

3）在機組正反向飛逸過程中，計算域的流量與總熵產值均呈現先減小后增大的規律。這是因為較大的流量將產生較大的速度梯度、湍動能與近壁區流速，繼而分別使得直接耗散熵產、湍流耗散熵產及壁面耗散熵產的增大。

4）在最終飛逸狀態下，高速旋轉的葉輪使得水流流速激增，因此葉輪域內的熵產值遠高于其他過流部件。與此同時，正向飛逸工況下動靜干涉作用比反向飛逸工況更加強烈，使得正向飛逸工況中扭矩的波動幅度更大。

5）在零流量時刻下，渦核聚集于泵工況的葉輪葉片進水邊，即水流在葉片進水邊形成強度較高的漩渦。而不同葉片展開面處渦核與熵產率的分布較為相似，因為渦核附近產生較大的速度梯度，繼而造成了較大的熵產率，表明漩渦是造成能量損失的原因。