張國棟,羅宇翔,李龍飛,唐桂華
(1. 西安交通大學能源與動力工程學院,710049,西安; 2. 西安交通大學熱流科學與工程教育部重點實驗室, 710049,西安; 3. 西安航天動力研究所液體火箭發動機科學技術實驗室,710100,西安)
對于液體火箭發動機,有效的冷卻技術能延長發動機的壽命,使其可重復使用[1-3]。冷卻技術包括燒蝕冷卻、輻射冷卻、排放冷卻、液膜冷卻以及再生冷卻等[4]。由于小推力液體火箭推進劑流量小、燃燒室體積小,主要采用液膜冷卻,其原理是噴注器將液膜冷卻劑噴射在燃燒室內壁,冷卻劑在內壁表面流動過程中受熱蒸發與內壁表面形成冷卻屏障,從而降低發動機表面溫度[5-6]。
在液膜冷卻模擬研究方面,楊薇等[7]和孫冰等[8]對液膜再生冷卻問題建立液膜傳熱一維分析計算模型,推導了液膜長度和厚度的工程計算方法。Amato等[9]建立了液膜冷卻的二維軸對稱CFD模型,研究了不同工況的冷卻效果。Strokach等[10]采用基于拉格朗日方法的離散相模型模擬了煤油的冷卻,該模型使用拉格朗日液滴在冷卻表面形成一層薄膜,并假設每一個撞擊表面的液滴都會轉化為液滴顆粒,通過動量守恒、傳熱傳質方程可以計算出液膜的速度、厚度和溫度等變量[11]。符鵬飛和韓兆鵬等[12-14]采用同樣的方法模擬得到了推進劑液膜在推力室內壁上的鋪展厚度,雖然離散相模型能夠獲得液膜厚度,但不能展現液膜流動形態。隨著數值模擬技術的發展,流體體積法(volume of fluid,VOF)在液膜冷卻方面得到了初步應用,林慶國[15]構建了基于VOF方法的射流模型,研究了單股液柱噴射在冷態壁面上的流動規律,同時獲得了液膜的厚度。在熱態壁面上利用VOF方法研究液膜冷卻涉及較少,需要進一步論證。
實驗研究方面,初期主要關注液膜在冷態條件下的流動規律。Good和Nollet[16]、唐亮等[17]利用高速攝相機獲得透明平板上液膜的分布情況,研究射流參數對液膜外形的影響規律。Sako等[18]利用可視化研究方法探究了單股水膜射流在不銹鋼平板上射流流量、噴嘴孔徑、射流角對液膜鋪展寬度的影響,但實驗中未能獲取液膜厚度。隨著實驗技術的發展出現了接觸式和非接觸式的液膜厚度測量方法,其中接觸式測量方法主要有電導法和探針法。Lee等[5]和Yuan等[19]采用電導法進行了液膜厚度測定。林慶國[15]采用接觸式探針法獲得液膜在金屬壁面上的鋪展厚度。非接觸式的液膜厚度測量方法主要是利用先進的光學測量儀器,根據不同位置液膜波的反射情況來確定液膜的厚度。龍黃祥[20]采用一套基于光學原理的白光共焦位移傳感器測量了降膜冷卻系統中的液膜厚度數值,且具備一定的精度要求。Stumpf等[21]采用一套共焦位移傳感器測量噴霧液滴濺落在透明平板上的液膜厚度,測試系統可獲得液滴沖擊平板前后不同時刻下的薄膜厚度,并與液滴動力學計算模型吻合較好。此外,平面激光誘導熒光技術(planar laser-induced fluorescence,PLIF)也是一種利用光學法測量液膜厚度的非接觸式測量技術[22]。王博[23]采用PLIF技術對波形板干燥器壁面上的液膜厚度進行測量,獲得了液膜厚度隨時間的變化規律。
綜上所述,對于實驗研究,目前很少關注液膜在熱態金屬壁面上的流動鋪展規律,而在熱態壁面上研究液膜冷卻更具有實際意義。在熱態壁面條件下,液膜鋪展形態沒有明顯的濺射、下游破裂及逆流現象,鋪展形態較為穩定,且金屬壁面的熱性能會影響液膜的傳熱特性[24]。在液膜厚度測量方面,相比于接觸式的測量方法,非接觸測量方法操作方便、測量精度高,因此本文采用白光共焦位移傳感器測量液膜厚度。在數值模擬方面,本文構建了基于VOF方法的熱態條件下液膜射流冷卻模型計算冷卻劑的蒸發吸熱、流動鋪展過程,驗證了VOF方法對射流撞壁形成液膜模擬的可行性。本文的研究結果可為液體火箭發動機液膜冷卻技術提供一定的參考。
液膜冷卻熱態不銹鋼曲壁實驗系統如圖1所示,其原理為:利用大功率電加熱裝置持續輸出最大500 A的恒定電流將不銹鋼曲壁加熱至150℃~180℃,制造熱態壁面環境;液體儲箱中的冷卻水通過水泵沿管路流動,通過控制水泵的轉速來實現射流流量的調節,通過玻璃轉子流量計獲得流量大小,測量誤差為示值的±0.5%;此后,冷卻劑通過噴嘴

(a)實驗系統示意圖

(b)實驗系統實物圖
形成射流液柱,撞擊到高溫實驗件壁面形成穩定的液膜形態,利用高速攝像機拍攝液膜流動形態;同時在液膜上方固定采集光筆,連接白光共焦位移傳感器獲得壁面上的液膜厚度;在不銹鋼曲壁實驗件的外表面焊接熱電偶探頭,不銹鋼壁面加熱過程中的溫升及液膜冷卻后的溫降過程由溫度采集系統獲得;液膜形態攝像、液膜厚度采集、壁面溫度采集的信號傳輸可由同一臺計算機控制。為了讓冷卻水從不銹鋼曲壁上順利流入排水槽中,將不銹鋼曲壁的夾具工作臺沿著液膜流動方向傾斜5°~10°。
在實驗過程中,將電加熱裝置的電源兩極連接到特制的銅板夾具上,讓電流通過不銹鋼曲壁實驗工件,利用電流的熱效應進行加熱,電加熱裝置及實驗件夾具如圖2所示。
液膜射流冷卻實驗系統中,對于射流角的調節是通過將噴嘴固定在一個旋轉微調機構上實現的,微調機構的精度為1°,具體射流角調節裝置如圖3所示。

圖3 射流角調節裝置Fig.3 Jet angle adjustment device
采用K型熱電偶對不銹鋼外壁的溫度進行監測,在不銹鋼曲壁的外表面布置8個溫度測點,取8個測點的平均溫度作為不銹鋼曲壁的溫度,溫度測點布置示意圖如圖4所示。具體的實驗所用儀器規格見表1。

圖4 溫度測點布置示意圖Fig.4 Schematic of temperature measurement point layout

表1 實驗所用儀器規格
在定量分析之前,需要對液膜厚度測量實驗系統進行驗證。將本文實驗與林慶國[15]的研究進行了比較,液膜厚度對比結果如圖5所示。文獻[15]中采用接觸式探針測量方法,而本文采用非接觸式的光學測量方法。圖5(a)對應文獻[15]中的實驗工況:入射角為20°、噴注壓力為0.3 MPa、射流孔徑為0.34 mm。圖5(b)對應本文的實驗工況:入射角為25°、射流流量為200 mL·min-1、射流孔徑為1 mm。雖然兩者在測量液膜厚度方法、射流工況上有所差別,但曲壁中軸線上的液膜厚度分布規律具有共同的特點:起初液膜在撞擊點上的厚度較大,隨著液膜不斷向下游流動鋪展,液膜厚度逐漸減小,達到最小值后又逐漸增大。這是因為射流液膜剛開始撞擊曲壁后,沖擊動量較大,導致中心軸線上的厚度減小,隨著液膜向下游流動鋪展逐漸穩定,在液體表面張力的作用下,液膜逐漸向中心軸線上收縮,在下游處匯聚,因此中心軸線上液膜厚度會逐漸增大,驗證了本文實驗中液膜厚度測量系統的正確性。

(a)文獻[15]結果

(b)本文實驗結果
實驗的誤差來源主要有3個,實驗儀器誤差δX1、測量誤差δX2以及數據處理誤差δX3,總誤差的計算公式[25]如下

(1)
式中:Rtot為總誤差;δXi為每個來源的誤差;n為誤差來源總數,其中?R/?Xi為間接測量物理量對直接測量物理量的偏導。實驗中的測量參數有液膜鋪展長度L、液膜鋪展寬度W、液膜鋪展擴張角α、不銹鋼曲壁溫度Tw、液膜厚度δ,每個測量參數對應不同來源誤差的數值見表2,最終5個測量參數總誤差計算結果分別是1.06%、1.70%、3.03%、8.70、2.01%。

表2 測量參數誤差計算
在本研究中,噴嘴的孔徑為1 mm,噴嘴距離壁面的高度為15 mm,射流角的范圍為25°~45°,液膜射流流量范圍為200~400 mL·min-1,設計了15組實驗工況,見表3。

表3 15組實驗工況統計表
在討論射流沖擊液膜厚度的變化規律前,可以定性觀察液膜在熱態不銹鋼曲壁上的流動鋪展形態。以工況13為例,圖6展示了典型的液膜流動鋪展形態,可見射流液膜撞擊曲壁后,液滴顆粒具備一定的動量,其中軸向動量轉化為水膜局部壓力,促使液膜向曲壁下游流動;而徑向分量迫使水膜以撞擊點為中心不斷向外圍流動,形成液膜沖擊區,然后匯入兩側邊緣的水躍區,水躍區為液膜厚度達到峰值的窄條形區域,唐亮等[17]將水躍區的邊界分為內邊界和外邊界。在液膜沖擊區,隨著液膜不斷向下游流動,液膜的慣性力不斷克服表面張力,從而使液膜在沖擊區域上出現最大鋪展寬度,之后液膜的動量減小,在表面張力作用下,液膜鋪展寬度逐漸減小,在沖擊區域的下游匯集。為了方便測量,主要研究了液膜撞擊點到匯集點附近曲壁中軸線上的液膜厚度變化規律。此外,本文建立的液膜流動形態參數示意圖如圖7所示。

圖6 典型液膜鋪展圖像Fig.6 Typical liquid film spreading image

圖7 液膜流動形態參數示意圖Fig.7 Schematic of film flow shape parameters
圖7中,o為射流撞擊點;x方向為液膜曲壁中軸線;Lmax為液膜沖擊區域最大鋪展長度;Wmax為液膜沖擊區域最大鋪展寬度;α為液膜圍繞撞擊點展開的擴張角。圖8顯示了不同工況的液膜鋪展形態。

圖8 不同工況下的液膜鋪展形態Fig.8 Liquid film spreading patterns under various operating condition
由圖8可以看出,在相同的液膜流量的工況下,隨著射流角的增大,射流撞擊曲面的撞擊點位置發生了明顯的偏移,但是測量液膜長度時都是以液膜撞擊點為起始點,匯集點為終點,且壁面傾斜角較小,產生的射流動量分量可忽略不計,因此射流角變化引起的撞擊點偏移不會對液膜長度測量產生影響。同時,液膜沖擊區域的鋪展寬度、擴張角增加,但是鋪展長度減小,對比工況1與工況5、工況6與工況10、工況11與工況15的3組液膜形態圖可以發現:隨著入射角的增大,下游處的匯集點與撞擊點的間距減小;而在射流角相同的工況下,隨著液膜流量的增加,射流沖擊區域液膜鋪展的長度、寬度和擴張角均有所增加。
為了定量考核不同工況下的液膜鋪展區域范圍,通過高速攝相機的Phantom Pcc軟件測量圖像中射流沖擊區域的最大長度Lmax、寬度W、擴張角α,對每個參數進行3次測量,取平均值及標準差,不同工況下液膜鋪展長度以及寬度和擴張角的結果分別如圖9、圖10所示。

圖9 液膜鋪展長度L變化曲線Fig.9 Variation of liquid film spreading length L
由圖9可知,不同射流冷卻工況下,工況11(θ=25°,Q=400 mL·min-1)的液膜長度達到最大值62.45 mm,而工況5(θ=45°,Q=200 mL·min-1)的液膜長度達到最小值26.43 mm。總體而言,液膜流動鋪展的區域的長度會隨著液膜流量的增大而增大,特別地當Q從300 mL·min-1增加至400 mL·min-1時,L存在明顯的遞增,ΔLmax出現在工況6(θ=25°,Q=300 mL·min-1)到工況11(θ=25°,Q=400 mL·min-1),ΔLmax為20.94 mm。這是因為當液膜流量增大時,射流動能增加,慣性力的作用加強,液膜越容易在壁面上鋪展,從而促進液膜向下游流動。此外,當Q=200 mL·min-1,θ由25°變化到30°時,L存在1.23 mm的增加量,這是由于液膜在射流流量較低、射流角較小的實驗工況下,液膜鋪展變化不明顯,實驗測量液膜長度存在一個小幅度的增加量,但是在Q大于200 mL·min-1射流流量的情況,隨著射流角θ的增加,液膜軸向動量減小,導致液膜向曲壁下游流動鋪展的長度減小。
液膜寬度和擴張角的變化曲線如圖10所示。當入射角θ增加時,液膜射流在壁面上的液滴出現一定的飛濺和反彈現象,使得大部分的液膜沖擊動量轉化為徑向分量,從而使得液膜徑向擴展寬度增加。例如工況15(θ=45°,Q=400 mL·min-1)最大液膜寬度為25.01 mm。類似地,液膜擴張角也隨著液膜入射角的增加而增加,例如工況15對應的α=120.20°。當θ超過35°時,液滴的沖擊行為減弱,W和α的增加趨勢不明顯,但射流液滴仍然具有一定的動量,并將液膜推向撞擊點外圍形成水躍區,增加了沖擊區的覆蓋范圍。因此,液膜流量和液膜入射角對于液膜冷卻規律的影響非常大。

圖10 液膜鋪展寬度W和擴張角α變化曲線Fig.10 Variation of liquid film spreading width W and spreading angle α
分析不銹鋼曲壁上液膜的流動鋪展規律后,需要進一步研究不同射流工況時高溫壁面的冷卻情況,結果如圖11所示。以θ=35°為例,由圖11(a)可知,不銹鋼曲壁被加熱900 s后外壁面平均溫度達172.67℃,900 s時開啟射流冷卻工況,壁面平均溫度迅速下降,當Q從200增加至400 mL·min-1時,冷卻至1 800 s,平均壁面溫度分別為40.78℃和35.40℃,溫度降低幅度分別為131.89℃和137.27℃。可見射流角一定時,增加液膜流量能夠有效降低壁面平均溫度。此外,液膜流量為300 mL·min-1,不同射流角θ的外壁面平均溫度變化曲線如圖11(b)所示。

(a)射流流量的影響

(b)射流角的影響
由圖11(b)可知,壁面加熱900 s后平均溫度為180.71℃,在射流入射角為35°時,冷卻后外壁面平均溫度達到最低值38.90℃。當入射角從25°增加至45°時,外壁面平均溫度分別降低136.01℃、137.84℃、141.81℃、138.91℃。由此可推斷,入射角為35°時能夠有效降低壁面溫度,隨著入射角θ的繼續增加,射流液膜在壁面的反彈和飛濺加強,損失了一部分液膜冷卻質量,壁面溫度降低效果減弱。
本文研究過程中主要關注曲壁中軸線上射流撞擊點到匯集點的液膜厚度變化規律。以工況1為例,當入射角為25°、液膜流量為200 mL·min-1時,不同位置處的液膜厚度變化如圖12所示。由圖12(a)可以看出,不同位置處的瞬態液膜厚度是隨機的,瞬態液膜厚度會在平均液膜厚度附近進行不對稱波動。在撞擊點x=0 mm處液膜的厚度波動最劇烈,液膜撞擊壁面后初始動量較大,液膜表面波很不穩定。隨著液膜向下游繼續流動鋪展,在x=15 mm處的波動幅度最小,液膜厚度瞬態峰值的波動相對平緩,在匯集點x=35 mm附近,分股后的液膜相遇,液膜厚度增大。此外,當入射角為25°時,不同射流流量對應的曲壁中軸線穩態液膜厚度變化情況分別如圖12(b)所示。

(a)瞬態液膜厚度變化曲線

(b)穩態液膜厚度變化曲線
由圖12(b)可知,在液膜撞擊點處液膜厚度存在一個峰值,且液膜流量越大對應的峰值較高,例如在撞擊點x=0 mm處液膜的厚度就是這樣,當Q=200 mL·min-1時,d=363.64 mm;當Q=400 mL·min-1時,d=679.32 mm。隨著液膜不斷向下游流動,液膜軸向動量減小,厚度逐漸降低,在軸線上的某一位置出現一個厚度波谷。然后在射流沖擊區的下游,液膜在表面張力的作用下匯聚,液膜厚度開始增加。
本文建立與實驗對照的液膜射流冷卻壁面模型,其中壁面曲率為50 mm,射流噴嘴的直徑為1 mm,計算域的長度為40 mm,高度為3.5 mm。為了節省網格數與模擬計算量,采用一半的模型進行模擬計算。非結構化網格適用于VOF模型的計算,容易收斂,因此本文使用ANSYS Meshing劃分將對稱部分幾何模型進行非結構網格劃分,整體網格尺寸為0.3 mm,最終的網格數為24.37萬,平均網格質量為0.73。為了精確計算液膜厚度,在壁面處添加6層邊界層網格,第一層高度為0.01 mm,增長率為1.2,并在噴嘴處對網格進行局部加密,液膜射流仿真模型網格及邊界如圖13所示。

圖13 液膜射流仿真模型網格及邊界Fig.13 Grid and boundary of liquid film simulation model
邊界條件包括液膜入口處的速度入口,冷卻工質為水,射流速度分別為4.3、6.4、8.5 m/s,分別對應實驗射流流量200、300、400 mL·min-1。液膜冷卻劑的入口溫度為290 K;出口為壓力出口,壓力為0 Pa,回流溫度為300 K;設置底面的加熱熱流為 6 800 W/m2,對應實驗中加熱壁面最高溫度為180℃。其余壁面均為絕熱壁面。模擬中考慮液膜重力,射流液膜與底面的傳熱計算采用realizablek-ε湍流模型及增強壁面函數處理,且底面邊界層網格滿足y+≈1的要求。模擬計算時間步長為10-5s,計算2 000步直至底面平均溫度達到穩定。以射流流量為200 mL·min-1、入射角為45°為例,模擬獲得不同網格數下的液膜形態參數變化如圖14所示。

圖14 不同網格數下的液膜形態參數變化Fig.14 Changes in liquid film parameters under different grid numbers
由圖14可知,分別采用網格數為20.35萬、22.38萬、24.37萬、26.32萬的4組網格對仿真模型進行計算,可以發現當網格數超過24.37萬時,模擬獲得的液膜形態參數L、W和α的變化量不足1%,因此本文采用的仿真模型網格數為24.37萬。
在本文的數值模擬中,使用基于Euler-Euler法的VOF模型求解射流計算域內的液相與氣相,分別求解各項的連續性方程如下

(2)


(3)
然后在整個計算域內求解動量方程,得到的速度場在各相之間共用。動量方程取決于所有相體積分數計算的密度和黏度值,方程如下

(4)
式中:ρ為控制體密度;p為壓力;μ為黏度;F為相之間的表面張力。此外,能量方程在各相間共用,方程的具體形式如下

(5)
式中:keff為有效導熱率;hj,q為流體j在q相里的焓;Jj,p為流體j在q相的擴散通量;τeff為有效剪應力;Sh為體積熱源;E為能量。此外,熱態壁面上液膜射流冷卻模型需引入傳質計算模型,本文采用Lee模型描述液膜受熱蒸發過程。在Lee模型中,從液相到氣相的傳質過程由氣體傳輸方程計算,表達式如下

(6)

當入射角為45°時,不同射流流量下通過數值模擬得到的液膜流動鋪展形態如圖15所示。

(a)流量為200 mL·min-1

(b)流量為300 mL·min-1

(c)流量為400 mL·min-1
由圖15可知,隨著液膜射流流量的增加,可以定性地看到液膜的鋪展形態不斷向四周擴展,液膜的鋪展面積不斷增加。當射流流量為300 mL·min-1時,在液膜鋪展的下游區域出現了明顯的液膜分股現象,當液膜射流流量增加到400 mL·min-1時,分股后的液膜繼續向下流動,液膜的鋪展范圍進一步擴大。此外沿著垂直于液膜流動x方向,距離噴嘴分別為5 mm和15 mm獲取兩個截面,用于液膜厚度的定量分析。當射流流量為200 mL·min-1時,液膜鋪展截面上液相體積分數分布云圖如圖16所示。

(a)距離射流噴嘴下游5 mm截面

(b)距離射流噴嘴下游15 mm截面
利用圖像處理方法可從液相體積分數分布云圖中獲取中心區和水躍區的液膜厚度,然后分別研究不同射流流量下液膜在上游區域及下游區域截面上厚度變化規律,將液膜厚度模擬值與實驗值進行對比,結果見表4、表5。

表4 距離噴嘴下游5 mm截面液膜厚度對比

表5 距離噴嘴下游15 mm截面液膜厚度對比
根據表4及表5的數據可知,在x=5 mm截面處的中心區域,射流流量為300 mL·min-1時模擬結果與實驗結果的最大相對偏差為7.9%;射流流量為400 mL·min-1時,水躍區液膜厚度最大模擬值為0.510 mm,而在下游x=15 mm截面處的水躍區液膜厚度最大模擬值為0.507 mm,均與實驗結果相近,說明基于VOF方法的熱態壁面液膜冷卻模型,能夠較好地預測射流液膜的流動鋪展形態及液膜厚度,為后續的液膜冷卻仿真模擬研究提供了基礎。
本文開展了液膜冷卻熱態金屬曲壁的實驗及仿真研究,得到以下結論。
(1)在相同的液膜流量的工況下,隨著射流角的增大,射流撞擊點會向曲壁上端移動,鋪展長度減小,同時液膜沖擊區域的鋪展寬度、擴張角增加。在射流角相同的工況下,隨著液膜流量的增加,液膜鋪展的長度、寬度和擴張角都有所增加。
(2)射流角一定時,增加液膜流量能夠有效降低壁面平均溫度。而當液膜流量一定時,隨著入射角的繼續增加,射流液膜在壁面的反彈和飛濺加強,損失了一部分液膜冷卻質量,壁面溫度降低效果減弱。
(3)曲壁中軸線上不同位置處的液膜厚度變化是隨機的,瞬態液膜厚度會在平均液膜厚度附近呈現不對稱波動,隨著液膜向下游繼續流動鋪展,液膜厚度峰值的波動相對平緩,到了下游匯集點附近,分股后的液膜相遇,液膜厚度增大。
(4)基于VOF方法的熱態條件下液膜射流冷卻仿真模型能夠較好地預測射流液膜鋪展形態及液膜厚度,模擬結果與實驗結果最大相對偏差為7.9%,控制在工程應用允許的10%誤差范圍內。仿真模型能夠擴展工況范圍,為進一步揭示液膜射流冷卻機理及規律奠定基礎。