劉宇寧，王秀勇，劉志遠，白小榜

（1.蘭州理工大學 能源與動力工程學院，蘭州 730050；2.重慶水泵廠有限責任公司國家級企業(yè)技術中心，重慶 400030）

0 引言

多級離心泵在工農業(yè)生產中的應用極為廣泛，由于多級泵水力模型的開發(fā)制造成本較高，因而在設計過程中通常借助于CFD技術對其水力性能進行預測，并對泵內的流場結構進行分析研究，為水力性能的進一步優(yōu)化提供有效的參考信息［1］。雖然CFD技術在流體機械行業(yè)中已被廣泛應用，但數值計算精度問題依然是當前研究的重點內容。由于目前關于多級離心泵的數值模擬研究多集中于泵內流動特征的分析，而關于提高水力性能預測精度的研究還很少，因而在既有理論基礎上，有必要對影響多級離心泵數值計算精度的因素進行研究，從而為提高水力性能的預測精度提供可靠而又可行的數值模擬方法。

研究結果表明，影響數值計算精度的主要因素有流動狀態(tài)是否定常、湍流模型的選取以及近壁面網格尺度等［2-5］。在工程應用中，關于流動控制方程的求解通常采用雷諾時均化（RANS）方法，該方法具有較強的計算經濟性，對計算機資源要求不高，在離心泵數值模擬中的應用最為廣泛［6-7］，但由于時均化過程會導致部分湍流脈動信息丟失，流場保真性較差；大渦模擬（LES）方法在反映湍流脈動特征方面具有RANS無可比擬的優(yōu)勢，常用于捕捉泵內部的渦核分布［8］，計算精度較高，但對計算資源的要求也很高；近年來，結合了LES與RANS各自優(yōu)點的分離渦模擬（DES）方法在流體機械的數值模擬研究中取得了良好的效果［9-10］，但計算過程中可能出現(xiàn)模化應力不足問題，而延遲分離渦模擬（DDES）方法可以對該問題進行修正，并且在離心泵數值模擬中的應用效果比較明顯［11］，有望成為今后離心泵數值模擬所選用湍流模型的主流。

由于目前探討離心泵數值計算精度影響因素的研究對象主要集中于單級單吸結構，為了系統(tǒng)地研究各因素對多級離心泵數值計算精度的影響，本文以某型號的多級離心泵為研究對象，采用ANSYS FLUENT計算軟件，首先對基于RNGk-ε湍流模型的定常和非定常計算的結果進行對比；然后比較分析RNGk-ε和DDES兩種湍流模型在多級離心泵數值模擬中的應用效果；最后探討葉輪近壁面網格尺度的不同對數值計算結果的影響。通過上述研究，為多級離心泵的數值模擬方法提供有價值的參考。

1 數值計算方法

1.1 計算模型

本文以并聯(lián)結構的節(jié)段式多級泵為研究對象，結合試驗泵的實際情況，級數取為兩級，其設計參數為：流量Q=2 400 m3/h，揚程H=190 m，轉速n=1 480r/min，首級和次級葉輪的葉片數均為6枚，導葉的葉片數為8枚。由于并聯(lián)式結構關于末級壓出室中心面兩側對稱，為減少數值計算的工作量，計算域在此取其中一側，包括吸入室、首級葉輪、次級葉輪、導葉、壓出室、葉輪前后腔以及口環(huán)間隙，整體幾何模型如圖1所示。

應用ICEM CFD進行網格劃分。由于泵體結構復雜，若進行六面體結構化網格劃分，則工作量巨大，耗時很長，技術要求極高，難以達到在一般工程應用中進行普遍推廣的目的，因而對吸入室、導葉和壓出室采用容易實施的四面體非結構網格劃分方案。葉輪是影響泵性能好壞的核心部件，葉輪內流場結構的模擬精度是決定泵整體性能計算精度的關鍵，為了研究葉輪近壁面處的網格尺度對數值計算精度的影響，對葉輪采用兩種網格劃分方案，一種是進行四面體非結構化網格劃分，另一種是進行六面體結構化網格劃分并加邊界層。葉輪的網格劃分情況如圖2所示。

圖1 幾何模型

圖2 葉輪網格劃分

根據網格無關性驗證結果，當網格總數增加到1 100×104后，揚程和效率的計算值趨于穩(wěn)定，此時吸入室、導葉和壓出室的網格尺度設定為12 mm；葉輪采用四面體網格劃分時設定葉片表面網格尺度為6 mm，全局網格尺度為8 mm；葉輪采用六面體網格劃分時，設定近壁面第一層網格的高度為0.05 mm，網格尺度增長比率為1.2。

1.2 湍流模型

RNGk-ε湍流模型［12］采用重正化群（RNG）的統(tǒng)計技術，它會根據不同強度的旋流來適當地改變湍流黏度，進而對旋流流動進行修正，同時它也兼顧了平均流動中的旋轉及旋流情況。

DDES湍流模型［13］是在DES的基礎上，將延遲函數引入到方程當中，混合長度尺度被重新構建，從而解決了DES模型中由網格誘導產生的分離流動（Grid Induced Separation，GIS）問題，進而減弱了DES對網格的依賴性。SST-DDES對長度尺度的修正公式為：

式中l(wèi)DDES——DDES混合長度尺度，m；

lRANS——RANS混合長度尺度，m；

fd——DDES模型的延遲函數；

CDES——DES函數校準常數，CDES=0.61；

Δmax——本地網格單元最大尺度，m；

Cd1，Cd2——經驗常數，Cd1=20，Cd2=3；

rd——當地湍流尺度與到壁面距離的比值；

v——運動黏度，m2/s；

vt——旋渦運動黏度，m2/s；

k——卡門通用常數，k=0.412；

d——網格單元到壁面距離，m；

S——應變率張量，s-1；

Ω——渦量張量，s-1。

1.3 計算方法

采用有限體積法對控制方程進行離散。壓力和速度的耦合采用SIMPLEC算法。定常計算時，采用動參考系模型；非定常計算時，采用滑移網格模型，并使用前面的定常計算結果作為初始流場。壁面采用無滑移邊界條件；入口給定速度邊界條件；出口給定自由出流邊界條件；靜止域與旋轉域之間采用Interface邊界條件進行數據傳遞。將葉輪每轉過3°圓周角所需要的時間定義為時間步長，其值為3.378 38×10-4s。定常計算時，當泵出口總壓趨于穩(wěn)定時，可視為計算結果收斂；非定常計算時，當泵出口總壓呈現(xiàn)周期性波動時，取物理量在最后一個旋轉周期內的時均值進行性能計算。

2 計算結果分析

2.1 定常與非定常計算結果對比

圖3示出全計算域采用四面體網格時，在全流量工況范圍內，RNGk-ε湍流模型分別采用定常與非定常2種計算方法預測得到的揚程、軸功率和效率的計算結果。

圖3 RNG k-ε湍流模型定常與非定常計算結果比較

由圖3（a）所示的揚程計算結果來看，定常與非定常計算在全流量工況范圍內的結果有明顯的不同，其中非定常計算的揚程在各個工況點均高于試驗值，但揚程計算曲線隨流量增加時的變化趨勢與試驗曲線吻合良好；而定常計算的揚程在（0～0.3）Qd以及（1.1～1.2）Qd工況范圍內高于試驗值，在（0.3～1.1）Qd工況范圍內低于試驗值，揚程的計算曲線與試驗曲線存在兩個交點，導致?lián)P程計算曲線隨流量增加時的變化趨勢與試驗曲線不完全一致，在（0～0.8）Qd工況范圍內下降速度較快，在（0.8～1.2）Qd工況范圍內下降速度又相對緩慢。從揚程相對計算誤差大小的角度來看，在（0.2～0.5）Qd以及（0.9～1.2）Qd工況范圍內，定常計算的揚程計算誤差小于非定常計算，原因是這兩個工況范圍處于定常計算曲線與試驗曲線的交點位置附近，當偏離這兩個工況范圍時，定常計算的相對計算誤差大于非定常計算；定常計算的揚程相對計算誤差最大值在關死工況點，為12.23%，其次是0.8Qd工況點，為9.10%，其他工況點的計算誤差位于2.54%～6.61%之間；非定常計算的揚程相對計算誤差最大值在1.2Qd工況點，為10.47%，其他工況點的計算誤差相差不大，均位于4.13%～7.32%之間。

由圖3（b）所示的軸功率計算結果來看，定常與非定常計算的結果差異明顯，軸功率的非定常計算值在全流量工況點都高于試驗值，但其計算曲線隨流量增加時的變化趨勢與試驗曲線吻合良好，能夠準確反映出軸功率與流量之間的特性關系；而軸功率的定常計算值在0.6Qd工況點之后小于試驗值，之前大于試驗值，并且在關死工況點附近遠大于試驗值，導致定常計算曲線的變化趨勢僅在（0.6～1.2）Qd工況范圍內與試驗曲線基本一致，但在（0～0.6）Qd工況范圍內隨流量的增加呈現(xiàn)下降的趨勢，與試驗曲線的變化趨勢相反，不能正確反映軸功率在小流量工況區(qū)間內的變化特征。從軸功率相對計算誤差大小的角度來看，定常與非定常計算均在關死工況點存在最大計算誤差，其中定常計算的相對計算誤差為53.36%，非定常計算為22.80%；定常計算的最小誤差在0.6Qd工況點，僅有0.09%，該工況點剛好是軸功率的定常計算曲線與試驗曲線的交點；非定常計算的誤差隨著流量的增加逐漸減小，在1.2Qd工況點達到最小值，為5.06%。

由圖3（c）所示的效率計算結果來看，定常與非定常計算的效率計算值均在設計工況點之前小于試驗值，之后大于試驗值；在（0～0.8）Qd工況范圍內，效率的非定常計算曲線與試驗曲線的吻合度明顯高于定常計算；在（0.8～1.2）Qd工況范圍內，效率的定常與非定常計算結果比較接近，隨著流量的增加效率的計算值也在增大，而試驗值則變化平緩，效率的計算曲線與試驗曲線在設計工況附近有明顯的區(qū)別。從效率相對計算誤差大小的角度來看，定常與非定常計算的最大誤差均在0.2Qd工況點，其中定常計算的相對誤差大小為23.72%，非定常計算為7.36%；定常計算誤差的次高點在0.3Qd工況點，為14.89%，其他工況點相對計算誤差的大小位于1.95%～7.25%之間；非定常計算誤差的次高點在1.2Qd工況點，為6.22%，其他工況點相對計算誤差的大小位于1.98%～3.97%之間。

總體來看，在全流量工況范圍內，多級泵各性能參數的非定常計算曲線與試驗曲線的吻合度都遠高于定常計算，尤其是效率的非定常計算結果，無論是在計算精度方面還是在計算曲線隨流量的變化趨勢方面，都具有定常計算不可比擬的明顯優(yōu)勢。因此，對多級泵進行性能預測時，采用非定常計算的方法可以獲得更精確的計算結果。

2.2 RNG k-ε與DDES湍流模型計算結果對比

圖4示出全計算域采用四面體網格時，在全流量工況范圍內進行非定常計算的前提下，RNGk-ε與DDES 2種湍流模型關于揚程、軸功率和效率的計算結果。

圖4 RNG k-ε與DDES湍流模型非定常計算結果比較

由圖4（a）所示的揚程計算結果來看，在全流量工況范圍內，DDES湍流模型的揚程計算值明顯低于RNGk-ε湍流模型，其計算曲線與試驗曲線之間具有極高的吻合度，計算精度遠高于RNGk-ε湍流模型。DDES的相對計算誤差最大值在1.2Qd工況點，為2.48%，其次是0.4Qd工況點，為2.09%，在其他工況點的相對計算誤差值均小于1.30%。

由圖4（b）所示的軸功率計算結果來看，在全流量工況范圍內，DDES湍流模型的軸功率計算值明顯低于RNGk-ε湍流模型，并且也都低于試驗值，但二者的軸功率計算曲線隨流量增加時的變化趨勢仍然保持較高的相似度。從軸功率相對計算誤差大小的角度來看，DDES湍流模型在全流量工況范圍內的計算誤差都明顯低于RNGk-ε湍流模型，其相對計算誤差的最小值在關死工況點，為0.05%，其次是1.2Qd工況點，為1.95%，在（0.2～1.0）Qd工況范圍內的相對計算誤差值比較穩(wěn)定，基本維持在4.22%左右。

由圖4（c）所示的效率計算結果來看，在全流量工況范圍內，DDES湍流模型的效率計算值明顯高于RNGk-ε湍流模型，并且也都高于試驗值，但DDES湍流模型的效率計算曲線在設計工況點附近更為平坦一些，可以更準確地反映效率隨流量增加時的變化趨勢。從效率相對計算誤差大小的角度來看，在全流量工況范圍內，DDES湍流模型的計算誤差值均高于RNGk-ε湍流模型，其計算誤差的平均值比RNGk-ε高出1.7%左右。

總體來看，當計算域全局采用四面體網格時，DDES湍流模型在多級泵的揚程和軸功率的計算精度上比RNGk-ε湍流模型具有明顯的優(yōu)勢，雖然在預測效率計算曲線隨流量增加時的變化趨勢上也具有一定的優(yōu)勢，但在計算精度上與RNGk-ε湍流模型相比卻處于明顯的劣勢狀態(tài)。實際上，DDES湍流模型對近壁面的網格尺度具有較高的要求，而四面體網格劃分的方法遠達不到其要求，因而下面將進一步探討近壁面網格尺度的大小對DDES湍流模型計算結果的影響。

2.3 葉輪近壁面網格尺度對計算結果的影響

吸入室、導葉和壓出室均采用四面體網格，葉輪分別采用四面體（無邊界層）和六面體（有邊界層）網格進行數值計算。根據數值計算的統(tǒng)計結果，葉輪近壁面網格的y+值分布情況如圖5所示，從圖可以看出，網格無邊界層時的y+值集中分布在250～1 000的區(qū)間內，有邊界層時的y+值集中分布在10～30的區(qū)間內。雖然葉輪的兩種網格劃分方案的y+值均不滿足DDES湍流模型y+＜5的要求，但仍可以用來定性討論分析近壁面網格尺度的變化對計算結果的影響，更細致的研究將在后期展開。

圖5 不同網格劃分方案近壁面y+值分布情況

圖6示出采用非定常計算時，DDES湍流模型在葉輪有無邊界層網格的情況下關于揚程、軸功率和效率的計算結果。

圖6 DDES湍流模型有無邊界層計算結果比較

由圖6（a）所示的揚程計算結果來看，當在葉輪壁面附近加入邊界層使y+值大幅減小時，揚程的計算曲線與試驗曲線幾乎完全重合，與無邊界層的計算結果相比，加入邊界層后揚程的計算精度有了明顯提高。葉輪網格加入邊界層后，揚程的最大相對計算誤差在1.2Qd工況點，僅有1.43%，其次是設計工況點，為1.06%，其他工況點的相對計算誤差值均在0.50%左右。

由圖6（b）所示的軸功率計算結果來看，在葉輪壁面附近加入邊界層時，軸功率的計算值比無邊界層時有所提高，但在0.2Qd工況點之后仍小于試驗值。在考慮了邊界層對計算結果的影響后，除了在關死點軸功率的相對計算誤差比無邊界層時增加了1.17%外，在其他工況點的相對計算誤差均明顯減小，最低可以減小0.55%，最高減小了2.87%，平均能夠減小1.20%。

由圖6（c）所示的效率計算結果來看，有邊界層時的效率計算曲線與試驗曲線的吻合度進一步提高，在全流量工況范圍內比無邊界層時的計算結果更精確，除了關死工況點的效率在理論上不存在計算誤差外，其他工況點有邊界層時的相對計算誤差均比無邊界層時減小了2.50%左右。

總體來看，當減小葉輪近壁面處的網格尺度時，DDES湍流模型關于多級泵各性能參數的計算精度均有明顯提高，并且可以推測，當葉輪尤其是所有過流部件近壁面網格尺度的y+值完全滿足DDES湍流模型的要求時，其計算精度會有進一步的提高。

2.4 各數值計算方法計算誤差綜合分析

為了綜合比較上述各數值計算方法的優(yōu)劣，在此對各性能參數在全流量工況范圍內相對計算誤差大小的平均值及相對計算誤差的均方差進行統(tǒng)計，如表1所示。

表1 全流量工況點計算誤差分析 （%）

相對計算誤差大小的平均值可以反映計算值與試驗值之間的偏離程度，而均方差則可以反映計算曲線與試驗曲線在曲線形狀上的相似度；平均值越小，則計算值越接近試驗值；均方差越小，則計算曲線隨流量增加時的變化趨勢越接近試驗曲線。

由表1所示各計算方法相對計算誤差的對比來看，在全流量工況點，定常計算方法的計算誤差最大，計算曲線與試驗曲線在曲線形狀的相似度上也最低，而非定常計算方法在上述兩個問題上都有明顯改善。計算域全局采用四面體網格并且進行非定常計算時，雖然DDES湍流模型關于效率的計算誤差比RNGk-ε湍流模型高出了1.72%，但關于揚程、軸功率的計算誤差卻都比RNGk-ε湍流模型低了4%，并且各性能參數相對計算誤差的均方差也都比RNGk-ε湍流模型低很多，也就是DDES湍流模型在預測性能曲線隨流量增加時的變化趨勢上比RNGk-ε湍流模型要精確很多，所以綜合來看，在相同的計算方法下，DDES湍流模型比RNGk-ε湍流模型更適合于預測多級泵的性能參數。當采用DDES湍流模型進行非定常計算時，在葉輪近壁面區(qū)域加入邊界層網格后，各性能參數的計算誤差進一步減小，計算曲線與試驗曲線的吻合度進一步提高。

綜上所述，對多級泵的性能參數進行預測時，采用DDES湍流模型和非定常計算方法，能夠獲得較高的計算精度，尤其是對葉輪的近壁面網格進行加密處理時，可以使計算精度進一步提高。在該數值計算方法中，網格劃分容易實現(xiàn)，網格數量容易控制，非定常計算所需要的時間也僅為定常計算的4～5倍，但計算精度卻比常規(guī)的采用RNGk-ε湍流模型和定常計算的方法高出很多，適合于一般工程應用。

3 結論

（1）在全流量工況點，非定常計算方法不僅計算誤差比定常計算方法小，并且性能計算曲線的形狀也更接近于試驗曲線；定常計算時在個別工況點的計算精度極高具有一定的偶然性，與該工況點距離計算曲線和試驗曲線交點位置的遠近有關。

（2）全計算域采用四面體網格和非定常計算的條件下，DDES湍流模型關于效率的計算精度稍低于RNGk-ε湍流模型，但關于揚程和軸功率的計算精度卻遠高于RNGk-ε湍流模型，并且各性能參數計算曲線的形狀更接近于試驗曲線。

（3）當對近壁面網格進行加密處理時，DDES湍流模型的計算精度可以進一步提高。