誘導輪非對稱空化流動數值計算和試驗研究

李龍賢，溫永鑫

（1.北京航天動力研究所；北京 100076；2.北京電子工程總體研究所；北京 100854）

0 引言

在航天推進系統中，對作為中大推力液體火箭發動機關鍵增壓裝置的渦輪泵的工作穩定性和可靠性等有嚴格的要求。誘導輪作為渦輪泵的前置增壓部件，對改善渦輪泵的抽吸性能、降低推進劑貯箱內部壓力及質量具有重要作用。高性能誘導輪轉速較高，通常允許流場內存在一定程度的空化，但空化的產生會使流場趨于不穩定，甚至對發動機的工作穩定性和可靠性產生影響。

誘導輪內的空化存在多種形式，非對稱空化是已知的一種存在于渦輪泵揚程下降過程中的復雜流動現象。對于附著型非對稱空化，誘導輪內各通道的空化區長度表現出明顯的不同。在此工況下，須單獨考慮和分析每個葉片，但這將使得空化模型變得非常復雜。Horiguchi等[1]利用奇點方法給出了預測交替型空化的近似模型，為非對稱型空化的預測提供了重要的借鑒。本文基于奇點法理論，并結合空化波動理論對某低溫介質三葉片誘導輪的非對稱空化特征進行研究。

1 研究方法

1.1 研究對象

以某低溫介質三葉片變螺距誘導輪為研究對象。為了計算及試驗方便，將離心輪、擴壓器、誘導輪等組成渦輪泵整體，模型如圖1所示，誘導輪實物如圖2所示。誘導輪的實際工作介質為低溫液氧，但在本文研究中由于采用低溫介質進行試驗難度和風險較大，故采用密度和黏度相近的常溫水代替。

圖1 液氧渦輪泵幾何模型Fig.1 Liquid oxygen turbopump geometry model

圖2 誘導輪實物圖Fig.2 Physical view of inducer of liquid oxygen turbopump

1.2 控制方程

控制方程主要由連續方程、動量方程、能量方程構成。其中連續方程為混合相質量守恒方程，在連續方程中增加一個氣態相或液態相質量守恒方程并帶入質量輸運源相，便可導入空化模型[2]。

（1）連續方程

式中：t為時間；ρ為氣液兩相流體密度，ρ=αlρl+αvρv；U為速度向量；?(ρU)為速度散度。αl為液體體積分數；ρl為液體密度；αv為氣體體積分數；ρv為氣體密度。

（2）動量方程

式中：U×U為向量積；τ為表面力，?τ為表面力τ的散度；p為流場壓力；SM為由體積力牽引的動量源項。

（3）能量方程

式中：e為工質內能；U為速度標量值；T為工質溫度；q?為工質與外界的熱交換量；SΦ為耗散函數。

（4）空化模型及空泡動力學方程

空化模型為基于均質多相質量輸運模型的空化模型（又稱1方程模型），并假設相間無滑移。數值計算過程中，包含多相流及空化的計算模型適用于常規的黏性流體力學計算方法，常規計算單相流的湍流模型對于多相流流場計算依然適用[3]。

式中：Re、Rc分別表示氣泡產生和氣泡潰滅的質量輸運源項。

空泡動力學方程采用Rayleigh-Plesset方程[4]，即：

式中：RB為空泡半徑；pv為飽和壓力；p∞為遠場壓力；S為液體表面張力系數；νl為液體運動黏度。

采用上述組合模型計算出的結果能較好地反映出空化區隨流場壓力的變化情況?；赗ayleigh-Plesset空泡動力學理論推導出的空化模型為均勻化模型，目前很難用質量輸運空化模型準確預測附著型非對稱流動[5]流場內空化區波動以及不對稱特性。因此，本文引入奇點法表示空化波動方程，即：

式中：Ue-iαz表示主流流動勢；N表示葉片數；wn為流場受第n個葉片的擾動，其物理意義是受源項及葉片上漩渦影響的流動勢，單位為J/g。

采用空化數σ定量描述空化的劇烈程度[6]，空化數為無量綱參數，其數學表達式為：

式中：p∞和u∞分別為液體流動參考壓力和參考速度，在葉輪機械流場中，p∞通常取入口總壓pin；u∞取葉尖圓周速度uT，因此空化數可以表示為[7]：

式中：Ω為角速度；rT為葉片圓周半徑。若不計比尺效應，理論上如果空化數σ相同，空化程度也相同[8]。

1.3 試驗研究

圖3為誘導輪空化試驗及測試系統示意圖。系統主要由三個分系統組成：泵組合體外特性參數試驗系統、流場脈動壓力測量系統和空化流場可視化試驗系統。由于直接用低溫工質進行試驗的難度大，成本較高，試驗風險較大，因此，采用與實際低溫工質液態氧黏度和密度等物性參數相近的常溫水進行試驗。

圖3 誘導輪空化試驗及測試系統Fig.3 Diagram of inducer cavitation test system

圖4為試驗系統現場圖。試驗過程中采用高速相機通過觀察窗口連續高速拍攝流場的瞬態圖像。

圖4 可視化試驗系統現場圖Fig.4 Visualization test system

誘導輪空化流場觀測窗口采用透明有機玻璃制成，圖5為觀測窗口實物圖。

圖5 誘導輪空化觀測窗口實物圖Fig.5 Inducer cavitation observation window

2 數值計算和試驗結果分析

2.1 數值計算結果分析

圖6為試驗和數值仿真計算結果的對比。

圖6 空化數-泵出口壓力系數曲線Fig.6 Cavitation number-pump outlet pressure coefficient curve

由圖6可以看出，由數值計算得出的泵出口壓力系數ψ隨空化數σ的變化關系曲線與試驗結果吻合較好，說明由數值計算得出的泵組合體外特性參數比較準確，最大誤差在5%以內。壓力系數：

式中：pout為泵出口壓力，壓力系數和空化數均為無量綱參數。

通過數值計算和試驗可以發現，本文研究的誘導輪，非對稱空化區開始出現的空化數σ≈0.03，空化數區間為0.01≤σ≤0.03，非對稱空化出現的空化數與斷裂揚程發生時空化數很接近（斷裂揚程指泵出口壓力系數急劇下降階段所對應的泵揚程）。

以空化數σ=0.022為例，介紹試驗過程中觀測到的非對稱空化區的特征。圖7為σ=0.022情況下誘導輪葉片和導葉葉片表面動態壓力分布云圖?？梢钥吹?，此時低壓區已覆蓋葉片前緣部分，壓力最低的區域分布在吸力面修圓的表面和葉尖處。

圖7 葉片表面壓力動態分布試驗云圖Fig.7 Test image of pressure dynamic distribution on blade surface

圖8為σ=0.022工況下誘導輪葉片通道內空化區隨時間變化的動態云圖，（a）中的0°（0 s）為初始角度（初始時刻）流場試驗云圖，（b）（c）（d）依次為誘導輪轉過90°、180°、270°時在同一觀測位置觀察到的誘導輪葉片流場云圖?？梢钥闯?，此時空化區主要分布在葉片前緣附近吸力面表面及葉片修圓吸力面上，葉片通道間的空化區表現出明顯的非對稱性，同時每個葉片上的附著空化區也在不斷地伸縮變化。

圖8 葉片通道內空化區動態分布試驗云圖Fig.8 Cavitation area dynamic distribution on stator and rotor blade surface

圖9為采用數值計算得到的葉片通道內空化區 的動態分布圖。除了附著在葉片修圓段外緣以外的空化區，靠近葉根的空化區以云狀氣團的形式 分布在葉片通道入口，對通道造成阻塞。

圖9 葉片通道內非對稱空化區數值計算云圖Fig.9 Nonsymmetrical cavitation area in flow passge obtained by numerical simulation

非對稱空化出現時伴隨著大量的空化區氣泡的脫落和破碎，按照空化的類別，此時的空化屬于片空化和云空化的混合態[9]。云空化是由片空化發展到一定程度氣泡在下游脫落破碎產生的，是由大量微泡群和液態工質混合而成的兩相流。云空化的形成與片空化的尾部回射流以及片空化的整體失穩密切相關。云空化的出現總是與流動失穩、旋渦和湍流脈動等復雜流動現象聯系在一起，其發生、發展和潰滅會帶來結構的震動、噪聲和空蝕[10]。此時泵在外特性上表現為揚程降低3%，工程上認為泵內開始發生空蝕。云空化的出現從微觀的角度解釋了空蝕出現的原因。

圖10為空化數σ=0.022時計算得到的誘導輪葉片修圓末端橫截面上的壓力動態分布?？梢钥闯?，此時截面大部分被低壓區覆蓋，壓力較高的區域只出現在葉片壓力面葉尖處及葉片表面小部分區域，壓力的等壓線隨時間變化非常明顯。

圖10 葉片修圓末端橫截面壓力動態分布數值計算云圖Fig.10 Pressure dynamic distribution on blade tip relief cross section by numerical calculation

圖11為σ=0.022時計算得到的該截面空化區隨時間的動態分布?？梢钥吹?，空化區主要分布在吸力面葉尖及葉中表面上。從圖11可以清晰地看到空化區在三個葉片的非對稱性分布。為便于分析，將誘導輪葉片編號，圓心處的箭頭方向為誘導輪旋轉方向，附著空化區面積最大的為葉片1，沿著與誘導輪旋轉方向相反的方向依次為葉片2和葉片3；葉片1、2、3下游的葉片通道分別為通道1、通道2和通道3。葉片1表面的空化區體積分數明顯高于葉片2和葉片3表面；從圖11可以較明顯地看到空化區氣泡脫落過程。空化區氣泡脫落發生在葉片吸力面，部分空化區氣泡脫落并逐漸消失（圖11紅色虛線圈內為空化區變化過程）?？栈瘏^氣泡脫落是云空化形成的誘導因素。

圖11 葉片修圓末端橫截面空化區動態分布數值計算云圖Fig.11 Cavitation area dynamic distribution on blade tip relief cross section by numerical calculation

圖12為空化數σ=0.022時誘導輪三個葉片表面平均空化區體積分數?？梢钥吹?，此時葉片上的空化區以大致相同的周期隨時間不斷波動，空化區波動周期約為葉輪旋轉周期的3倍；葉片3表面的空化體積分數平均值比葉片2高出約25%，三個葉片上的空化區呈現出明顯的非對稱性。

圖12 誘導輪葉片表面空化區體積分數波動曲線Fig.12 Fluctuation curve of inducer blade surface cavitation area volume fraction（σ=0.022）

2.2 試驗結果分析

圖13～15分別為誘導輪3個葉片依次經過觀察窗口的連續成像，試驗工況為轉速4 000 r/min，空化數σ=0.022，流量系數φ=0.059，在此階段發現了比較強烈的非對稱空化，通道1的空化區比其他兩個要劇烈得多。Tsuijimoto[9]也在試驗中發現了這種非對稱空化，并證實這種空化不穩定形式出現的空化數要小于旋轉空化。這種形式的非對稱流動流道中的空化區較穩定，空化區隨葉片旋轉而同步旋轉，這與旋轉空化截然不同，將這種非對稱空化稱之為附著型非對稱空化。這種流道中的非對稱空化區為云空化的一種，誘導輪中云空化的出現是揚程開始顯著下降和葉片空化破壞的誘導因素，云空化出現時對應的空化數與斷裂揚程對應的空化數非常接近。

圖13 通道1空化區連續成像Fig.13 Continuous imaging for the first blade channel cavitation area

圖16為從試驗結果中提取的空化區光學圖像，從該組圖像中可以更清楚地看到誘導輪3個葉片通道中的空化區差異。

圖14 通道2空化區連續成像Fig.14 Continuous imaging for the second blade channel cavitation area

圖16 空化區光學圖像Fig.16 Cavitatoin area optical image

圖15 通道3空化區連續成像Fig.15 Continuous imaging for the third blade channel cavitation area

圖17為經過圖像處理后得到的空化區等高線圖像。

圖17 不同通道的空化區等高線圖像Fig.17 Cavitation area contour line image

對圖17的圖像進行光譜亮度分析，得到空化區頻譜能量概率統計圖，如圖18所示，這是對空化圖像特征的定量分析[11]。可以看到，通道1的頻譜能量峰值均明顯大于其他兩個通道，表明通道1內的空泡數量、空泡密度均大于另外兩個通道，通道2和通道3內的空化區分布規律比較接近。

圖18 頻譜能量概率統計Fig.18 Spectrum energy probability statistics

利用圖像的灰度直方圖提取諸如均值、方差、能量及熵等來描述圖像特征[13]。設（x，y）為圖像中的一點，該點和離它只有微小距離的點（x+Δx，y+Δy）的灰度差值gΔ為：

式中：設灰度差值的所有可能取值有m級，點（x，y）在整個圖像上移動，統計出gΔ取各個數值的次數，由此做出gΔ的直方圖。由直方圖可以得到gΔ取值的概率p（i），p（i）表示第i級灰度出現的概率，相關的特征參數有[14]：灰度平均值灰度對比度，其中灰度平均值表示圖像的亮度，平均值越大，圖像平均亮度越亮，圖像中的空泡越多。對比度表示圖像的亮度差值，對比度越大，圖像中的空泡越密集，熵表示圖像的平整程度，熵越大，空化區分布越均勻。

對表1空化區圖像特征參數進行分析，可以看到，通道1的灰度平均值、對比度及熵均明顯大于其他兩個通道的空化區，表明通道1內的空泡數量、空泡密度均大于另外兩個通道。通道2和通道3內的空化區分布規律非常接近。

表1 空化圖像特征分析（σ=0.022）Tab.1 Characteristic analysis for cavitation image（σ=0.022）

圖19為經過平均化處理的不同葉片背面上的附著空化區長度及不同通道的流量隨空化數的變化關系。

圖19 通道空化區長度及流量的歸一化處理Fig.19 Averaging treatment for cavitation area length and flux in different channel

圖19曲線上標注的數字對應葉片編號，lc表示空化區長度，lm表示空化混合區長度，Q表示葉片流道內流量?？梢钥吹?，葉片上出現不對稱附著空化區的空化數范圍為0.02≤σ≤0.035，非對稱空化出現時通道內的流量也隨之不對稱，空化最嚴重的通道1流量最小，旋轉阻塞只是出現在空化程度最嚴重的通道1中[5]。隨著空化數繼續降低，空化區的不對稱性也逐漸降低，當空化數進一步降低，誘導輪揚程接近斷裂揚程，此時空化區已蔓延至葉片后緣，通道內的空化區不對稱性消失。

3 結論

誘導輪內的空化存在多種形式，非對稱空化是已知的存在于渦輪泵揚程下降過程中的復雜流動現象的一種。本文采用數值計算和試驗方法對某低溫介質誘導輪流場中的非對稱空化進行了研究，得出以下結論：

（1）數值計算和可視化試驗證明了非對稱空化的存在，數值計算結果與可視化試驗結果一致性較好，說明本文中所采用的奇點法空化波動方程數值計算方法預測的非對稱空化具有較好的準確性；

（2）揚程明顯下降時葉片通道出現明顯的非對稱性，對應誘導輪的空化數為σ=0.03；

（3）三葉片誘導輪內發生非對稱空化時，其中一個通道內的空化區明顯大于其他兩個通道，且每個通道的空化區大小以一個特定值為中心隨時間成周期性波動；

（4）非對稱空化出現后，通道內的流量也隨之變得不對稱，空化最嚴重的通道對應的流量最小。旋轉阻塞只是出現在空化程度最嚴重的通道中；

（5）非對稱空化的出現只發生在特定空化區域，對于本文研究的誘導輪0.02≤σ≤0.035。空化數進一步降低，誘導輪揚程接近斷裂揚程，此時空化區已蔓延至葉片后緣，葉片通道內的空化區不對稱性消失。