孫 郡， 宋 杰， 李清廉*， 崔 朋， 陳蘭偉

(1.國防科技大學(xué)空天科學(xué)學(xué)院， 長沙 410073；2.國防科技大學(xué)高超聲速?zèng)_壓發(fā)動(dòng)機(jī)技術(shù)重點(diǎn)實(shí)驗(yàn)室， 長沙 410073)

1 引言

近年來，憑借比沖性能高、結(jié)焦積炭性能好且成本低等優(yōu)勢，液氧/甲烷推進(jìn)劑組合在液體火箭發(fā)動(dòng)機(jī)中得到廣泛應(yīng)用。 為減少燃燒室內(nèi)的高溫環(huán)境對(duì)結(jié)構(gòu)及系統(tǒng)部件造成的損傷，液氧甲烷發(fā)動(dòng)機(jī)多采用再生冷卻的方式降低壁溫，達(dá)到熱防護(hù)的目的。

甲烷在冷卻通道入口溫度通常較低，而其臨界溫度又較高(190 K)，因此在冷卻通道內(nèi)甲烷通常會(huì)經(jīng)歷跨臨界狀態(tài)。 然而在臨界溫度附近，甲烷物性變化劇烈，冷卻效果下降，引起壁溫的急劇變化。 因此開展甲烷的相變傳熱過程研究及再生冷卻通道的優(yōu)化設(shè)計(jì)具有很大的理論及工程價(jià)值。

陳尊敬等研究了超臨界條件下甲烷在矩形通道中的湍流換熱過程，并探討了不同熱流密度及通道截面形狀對(duì)流固耦合的超臨界壓力傳熱過程的影響；孫冰等針對(duì)推力室中跨臨界甲烷三維耦合傳熱進(jìn)行了數(shù)值模擬，得出燃?xì)鈧?cè)壁面溫度的最大值出現(xiàn)在喉部上游，且由于突擴(kuò)突縮處冷卻劑湍流強(qiáng)度增強(qiáng)，冷卻劑側(cè)表面?zhèn)鳠嵯禂?shù)提高的結(jié)論；Votta 等通過實(shí)驗(yàn)?zāi)M發(fā)動(dòng)機(jī)冷卻通道的工況，研究了跨臨界甲烷的傳熱特性，得到了非對(duì)稱加熱通道的換熱關(guān)聯(lián)式；Urbano 等研究了圓管中跨臨界甲烷的對(duì)流換熱，總結(jié)得出了跨臨界甲烷發(fā)生傳熱惡化現(xiàn)象所需的條件。

目前國內(nèi)外針對(duì)甲烷再生冷卻主要研究的是超臨界和跨臨界傳熱特性，對(duì)于亞臨界甲烷的傳熱特性研究比較有限。 因此，本文基于微小通道低溫工質(zhì)相變傳熱模型，開展亞臨界甲烷的再生冷卻通道設(shè)計(jì)和傳熱特性分析，旨在補(bǔ)充亞臨界甲烷再生冷卻研究的不足，為變推力液氧/甲烷發(fā)動(dòng)機(jī)再生冷卻設(shè)計(jì)提供借鑒和參考。

2 輸入?yún)?shù)及計(jì)算模型

2.1 推力室及計(jì)算工況參數(shù)

為滿足設(shè)計(jì)需求，發(fā)動(dòng)機(jī)推力室設(shè)計(jì)參數(shù)如表1 所示。

表1 推力室設(shè)計(jì)參數(shù)[10]Table 1 Design parameters of the thrust chamber[10]

由于推力室再生冷卻方案設(shè)計(jì)的難點(diǎn)在于滿足低工況、亞臨界條件下的熱防護(hù)要求，考慮到甲烷的臨界壓力為4.64 MPa 及最大室壓為4.5 MPa，選取50%、30%、10% 3 種工況進(jìn)行甲烷再生冷卻方案設(shè)計(jì)及傳熱特性的分析驗(yàn)證，由RPA 軟件計(jì)算得到對(duì)應(yīng)工況參數(shù)如表2 所示。

表2 推力室變工況參數(shù)Table 2 Variable working condition parameters of thrust chamber

2.2 再生冷卻傳熱計(jì)算模型

再生冷卻是指利用推進(jìn)劑中的1 種或2 種組元，在噴入燃燒室之前先使其通過推力室壁上的冷卻通道對(duì)發(fā)動(dòng)機(jī)推力室壁面進(jìn)行冷卻，以達(dá)到熱防護(hù)的目的，其傳熱過程如圖1 所示。

圖1 再生冷卻傳熱示意圖Fig.1 Schematic diagram of heat transfer of regenerative cooling

考慮到熱流密度的量級(jí)及簡化計(jì)算的需求，認(rèn)為推力室中的傳熱過程主要包括3 部分:①高溫燃?xì)馀c推力室內(nèi)壁的對(duì)流換熱；②燃?xì)鈧?cè)壁面與冷卻劑側(cè)壁面之間的導(dǎo)熱；③冷卻劑側(cè)壁面與冷卻劑之間的對(duì)流換熱。 達(dá)到穩(wěn)態(tài)平衡后上述3 個(gè)過程的熱流密度相等。

2.2.1 燃?xì)鈧?cè)對(duì)流換熱系數(shù)

燃?xì)鈧?cè)對(duì)流換熱系數(shù)由半經(jīng)驗(yàn)巴茲公式求得，如式(1)所示。

式中，為定性溫度變換系數(shù)，如式(2)所示。

其中，為喉部直徑，為燃?xì)鉁箿囟认碌膭?dòng)力粘度，為燃?xì)鉁箿囟认碌亩▔罕葻?，為燃?xì)鉁箿囟认碌钠绽侍財(cái)?shù)，為燃燒室室壓，為燃燒室特征速度，為定性溫度變化系數(shù)，為喉部曲率半徑，，分別為喉部面積和分段處截面積，為燃燒室溫度，為燃?xì)獗葻岜取?/p>

2.2.2 冷卻劑側(cè)對(duì)流換熱系數(shù)

1)液相:當(dāng)2300 時(shí)，根據(jù)齊德－泰勒公式(3)計(jì)算液相層流區(qū)對(duì)流換熱系數(shù)。

2)氣相/液相:當(dāng)＞2300 時(shí)，根據(jù)Gnielinski 公式(4)、(5)計(jì)算液相湍流區(qū)或氣相對(duì)流換熱系數(shù)。

3)兩相:基于實(shí)驗(yàn)得到的微小通道中液氮流動(dòng)沸騰實(shí)驗(yàn)關(guān)系式(6)，計(jì)算兩相區(qū)對(duì)流換熱系數(shù)。

2.3 再生冷卻壓降計(jì)算模型

冷卻通道內(nèi)甲烷兩相流壓降Δ包括摩擦壓降Δ、加速壓降Δ及重力壓降Δ，摩擦壓降采用均相模型計(jì)算，如式(13)所示。

各壓降計(jì)算式如式(14)～(16)所示。

式中相關(guān)無量綱參數(shù)根據(jù)式(17)、(18)計(jì)算:

其中，為質(zhì)量流速，為通道長度，為干度，為液相密度，為氣相密度，為重力加速度，為水力直徑，為空泡份額，為均相下的兩相流密度。

3 通道方案設(shè)計(jì)及結(jié)果分析

3.1 再生冷卻通道方案設(shè)計(jì)

采用變截面冷卻通道的設(shè)計(jì)原則，通道為漸縮－漸擴(kuò)結(jié)構(gòu)，首先確定再生冷卻通道出、入口及喉部的槽寬，沿冷卻劑的流動(dòng)方向，設(shè)計(jì)入口處的槽寬為0.9 mm，并沿著噴管擴(kuò)張段逐漸減小，直到喉部減小至0.6 mm，經(jīng)過喉部后又沿噴管收縮段逐漸增加，到圓柱段處恢復(fù)至0.9 mm；根據(jù)設(shè)計(jì)經(jīng)驗(yàn)取肋高為1.5 mm， 壁厚為1 mm。冷卻通道結(jié)構(gòu)圖如圖2 所示。

圖2 冷卻通道結(jié)構(gòu)圖Fig.2 Structure diagram of cooling passage

針對(duì)5000 N 級(jí)變推力液氧/甲烷發(fā)動(dòng)機(jī)50%工況對(duì)應(yīng)參數(shù)設(shè)計(jì)再生冷卻方案，即燃燒室總壓2.25 MPa，再生冷卻通道入口壓力為3.45 MPa，入口溫度為111 K。 選取導(dǎo)熱系數(shù)較高的紫銅(398 W/(m·K))作為壁面材料，根據(jù)約束條件求得冷卻通道數(shù)量為40，對(duì)應(yīng)肋寬為1.8 mm。

槽寬，肋寬，肋高和通道水力直徑的沿程分布情況如圖3 所示。

圖3 通道參數(shù)沿推力室軸向分布圖Fig.3 Distribution of channel parameters along the axial direction of thrust chamber

3.2 再生冷卻通道傳熱特性分析

根據(jù)3.1 方案得到的再生冷卻通道結(jié)構(gòu)的具體參數(shù)，采用2.2、2.3 的計(jì)算模型對(duì)再生冷卻通道的傳熱過程進(jìn)行計(jì)算并分析。

3.2.1 冷卻劑壓降和流速

如圖4、5 所示，甲烷壓力沿流動(dòng)方向逐漸減小，通道入口處甲烷壓力為3.45 MPa，出口時(shí)壓力降為2.96 MPa，總壓降為0.49 MPa。 流速則逐漸增加，入口流速為7.98 m/s，到喉部處加速至16.34 m/s，出口處流速升至57.35 m/s。

圖4 冷卻劑壓力沿推力室軸向分布圖Fig.4 Distribution of coolant pressure along the axial direction of thrust chamber

圖5 冷卻劑速度沿推力室軸向分布圖Fig.5 Distribution of coolant velocity along the axial direction of thrust chamber

3.2.2 冷卻劑溫度和干度

圖6 冷卻劑溫度沿推力室軸向分布圖Fig.6 Distribution of coolant temperature along the axial direction of thrust chamber

冷卻劑干度沿推力室軸向分布情況如圖7 所示，根據(jù)曲線走勢變化可將相變過程分為3 個(gè)不同的區(qū)域。

圖7 冷卻劑干度沿推力室軸向分布圖Fig.7 Distribution of coolant dryness along the axial direction of thrust chamber

1)區(qū)域I:0.182 m＜＜0.3 m，此時(shí)＝0，屬于單相液態(tài)區(qū)；

2) 區(qū)域II:0.004 m＜＜0.182 m，此時(shí)01，屬于兩相區(qū)；

3) 區(qū)域III:0＜＜0.004 m，此時(shí)＝1，屬于單相氣態(tài)區(qū)。

其中，區(qū)域II 可以再分為強(qiáng)化換熱階段(0＜＜0.6)和傳熱惡化階段(0.6＜＜1)。

3.2.3 燃?xì)鈧?cè)對(duì)流換熱系數(shù)和冷卻劑側(cè)對(duì)流換熱系數(shù)

圖8 給出了燃?xì)鈧?cè)對(duì)流換熱系數(shù)沿推力室軸向分布情況，燃?xì)鈧?cè)對(duì)流換熱系數(shù)總體呈先增后減的趨勢，在喉部偏上游處達(dá)到峰值，由入口處2.01 kW/(m·K)增至喉部處9.71 kW/(m·K)，后逐步降至出口處2.67 kW/(m·K)。

圖8 燃?xì)鈧?cè)對(duì)流換熱系數(shù)沿推力室軸向分布圖Fig.8 Distribution of convective heat transfer coefficient of gas side along the axial direction of thrust chamber

圖9 給出了冷卻劑側(cè)對(duì)流換熱系數(shù)沿推力室軸向分布情況，當(dāng)冷卻劑處于單相液態(tài)區(qū)(＝0)時(shí)，采用公式(3)、(4)進(jìn)行計(jì)算，冷卻劑對(duì)流換熱系數(shù)由入口處26.46 kW/(m·K)緩慢增加至相態(tài)區(qū)轉(zhuǎn)折點(diǎn)處的109.34 kW/(m·K)；經(jīng)過轉(zhuǎn)折點(diǎn)后，進(jìn)入兩相區(qū)(01) ，采用公式(6)計(jì)算，由轉(zhuǎn)折點(diǎn)處的109.34 kW/(m·K)陡然增加至峰值處的3984 kW/(m·K)，進(jìn)入強(qiáng)化傳熱階段，然后逐步減小至1676 kW/(m·K)，發(fā)生傳熱惡化，陡然降至5.22 kW/(m·K)；之后進(jìn)入單相氣態(tài)區(qū)(＝1) ，根據(jù)公式(5)計(jì)算，在兩相區(qū)至單相氣態(tài)區(qū)轉(zhuǎn)折點(diǎn)處，由5.7 kW/(m·K)陡然增至77 kW/(m·K)。 由計(jì)算結(jié)果可知，與單相液態(tài)區(qū)相比，兩相區(qū)的冷卻劑側(cè)對(duì)流換熱系數(shù)更高，而氣態(tài)區(qū)的冷卻劑側(cè)對(duì)流換熱系數(shù)則高于兩相區(qū)的傳熱惡化階段。

圖9 冷卻劑側(cè)對(duì)流換熱系數(shù)沿推力室軸向分布圖Fig.9 Distribution of convective heat transfer coefficient of coolant side along the axial direction of thrust chamber

3.2.4 熱流密度

圖10 給出了熱流密度沿推力室軸向分布情況，當(dāng)0.182 m ＜＜0.3 m 時(shí)，冷卻劑處于單相液態(tài)區(qū)(＝0) ，熱流密度經(jīng)歷一個(gè)先增后減的過程， 并在喉部處靠近入口位置達(dá)到峰值(10.16 MW/m)；當(dāng)0.064 m ＜＜0.182 m，此時(shí)冷卻劑進(jìn)入強(qiáng)化傳熱階段( 00.6)，冷卻劑側(cè)對(duì)流換熱系數(shù)高于單相液態(tài)區(qū)，因此熱流密度也更大；當(dāng)0.004 m ＜＜0.064 m 時(shí)，冷卻劑處于兩相流傳熱惡化階段，此時(shí)冷卻劑側(cè)對(duì)流換熱系數(shù)急劇減小，故而導(dǎo)致熱流密度驟降至3.79 MW/m；當(dāng)冷卻劑進(jìn)入氣相區(qū)(＝1) ，即0＜0.006 m 時(shí)，冷卻劑側(cè)對(duì)流換熱系數(shù)減小，故熱流密度由3.90 MW/m增大至5.66 MW/m。

圖10 熱流密度沿推力室軸向分布圖Fig.10 Distribution of heat flux along the axial direction of thrust chamber

3.2.5 燃?xì)鈧?cè)壁溫和冷卻劑側(cè)壁溫

圖11 給出了燃?xì)鈧?cè)及冷卻劑側(cè)壁溫沿推力室軸向分布情況，可以看出，燃?xì)鈧?cè)壁溫與冷卻劑壁溫變化趨勢一致。 冷卻劑處于單相液態(tài)區(qū)時(shí)，燃?xì)鈧?cè)及冷卻劑側(cè)壁溫均經(jīng)歷一個(gè)先增后減的過程，在喉部偏上游處達(dá)到峰值；進(jìn)入兩相區(qū)后，甲烷先經(jīng)歷強(qiáng)化傳熱，冷卻劑側(cè)對(duì)流換熱系數(shù)增大，故兩側(cè)壁溫在這一階段均有減小；當(dāng)0.006 m ＜＜0.064 m 時(shí)，出現(xiàn)傳熱惡化現(xiàn)象，燃?xì)鈧?cè)及冷卻劑側(cè)壁溫均出現(xiàn)驟增，燃?xì)鈧?cè)壁溫達(dá)到峰值868 K，冷卻劑側(cè)壁溫達(dá)到峰值862 K；當(dāng)冷卻劑進(jìn)入氣相區(qū)后，同樣得益于對(duì)流換熱系數(shù)的增大，冷卻劑側(cè)壁溫減小至260 K，燃?xì)鈧?cè)壁溫減小至280 K。

圖11 壁溫沿推力室軸向分布圖Fig.11 Distribution of wall temperature along the axial direction of thrust chamber

4 傳熱特性影響因素分析

本節(jié)針對(duì)不同工況下燃燒室室壓及其對(duì)應(yīng)的冷卻劑流量對(duì)傳熱性能的影響展開分析，采取3.1 節(jié)的冷卻通道方案。 由于高工況下冷卻劑流量較大，冷卻效果較好，易滿足推力室冷卻要求，故選取50%、30%、10% 3 個(gè)工況進(jìn)行探討并驗(yàn)證方案的可行性，不同工況對(duì)應(yīng)參數(shù)如表2 所示。

4.1 冷卻劑干度x 和溫度Tf 沿程變化

圖12 給出了冷卻劑干度的沿程變化情況，冷卻劑在燃燒室圓柱段依次進(jìn)入兩相區(qū)，并且室壓/流量越小，進(jìn)入兩相區(qū)的位置越靠近上游。 圖13給出了冷卻劑溫度的變化，比較室壓為2.25 MPa和1.35 MPa 的2 條曲線可以看出，隨著燃燒室室壓的降低，冷卻劑出口溫度也隨之降低；然而當(dāng)室壓降至0.45 MPa 時(shí)，由于低工況下流量的急劇減小，導(dǎo)致傳熱能力下降，出口溫度升至386 K，遠(yuǎn)高于另2 種工況的出口溫度。

圖12 冷卻劑干度沿推力室軸向分布圖Fig.12 Distribution of coolant dryness along the axial direction of thrust chamber

圖13 冷卻劑溫度沿推力室軸向分布圖Fig.13 Distribution of coolant temperature along the axial direction of thrust chamber

4.2 冷卻劑壓力p 和流速u 沿程變化

冷卻劑的壓力及流速變化情況如圖14、15 所示，可以看出，冷卻劑壓降隨著燃燒室室壓的減小而減小，當(dāng)室壓為0.45 MPa 時(shí)，冷卻劑壓降為0.12 MPa，當(dāng)室壓增至2.25 MPa 時(shí)，冷卻劑壓降為0.488 MPa；隨著室壓的降低，冷卻劑出口流速逐漸增大，當(dāng)室壓為0.45 MPa 時(shí)，由于甲烷流量的減小導(dǎo)致冷卻能力下降，甲烷在冷卻通道內(nèi)更早轉(zhuǎn)變?yōu)闅鈶B(tài)，因此出口速度較高。

圖14 冷卻劑壓力沿推力室軸向分布圖Fig.14 Distribution of coolant pressure along the axial direction of thrust chamber

圖15 冷卻劑速度沿推力室軸向分布圖Fig.15 Distribution of coolant velocity along the axial direction of thrust chamber

4.3 燃?xì)鈧?cè)壁溫Twg

圖16 給出了不同工況下燃?xì)鈧?cè)壁溫沿推力室軸向變化情況。 由圖可知，＝0.45 MPa 時(shí)的曲線明顯區(qū)別于其他2 條曲線，這同樣是因?yàn)榈凸r下甲烷流量的減小導(dǎo)致冷卻能力下降，冷卻劑更早進(jìn)入兩相區(qū)傳熱惡化階段和氣態(tài)區(qū)，引起燃?xì)鈧?cè)壁溫的急劇升高。 另2 種工況下，冷卻劑在燃燒室圓柱段先經(jīng)歷強(qiáng)化傳熱階段，此時(shí)燃?xì)鈧?cè)壁溫較低，約為200 K 左右，后進(jìn)入傳熱惡化階段，燃?xì)鈧?cè)壁溫陡然增加，并且隨著燃燒室室壓的減小，進(jìn)入傳熱惡化階段的位置提前，兩相區(qū)燃?xì)鈧?cè)壁溫峰值也隨之減小。 10%工況時(shí)，燃?xì)鈧?cè)壁溫的最大值高達(dá)877 K，略高于材料的極限溫度，因此在低工況時(shí)，應(yīng)當(dāng)考慮改變推進(jìn)劑混合比或采取膜冷卻、輻射冷卻等方式來達(dá)到熱防護(hù)的目的。

圖16 燃?xì)鈧?cè)壁溫沿推力室軸向分布圖Fig.16 Distribution of gas side wall temperature along axial direction of the thrust chamber

5 結(jié)論

1)基于實(shí)驗(yàn)獲得了再生冷卻通道內(nèi)甲烷相變傳熱模型，以甲烷干度0.6 為界限，0～0.6 為強(qiáng)化傳熱階段，0.6～1 為傳熱惡化階段；

2)設(shè)計(jì)了5000 N 級(jí)變推力液氧/甲烷推力室再生冷卻方案，并分析驗(yàn)證其在寬變推力范圍內(nèi)的可靠性，在30%～100%工況下，設(shè)計(jì)的冷卻通道能夠滿足發(fā)動(dòng)機(jī)的熱防護(hù)要求；

3)分析了燃燒室室壓及對(duì)應(yīng)的冷卻劑流量對(duì)傳熱性能的影響，流量過低時(shí)，甲烷冷卻能力急劇下降，需考慮改變推進(jìn)劑混合比或采取膜冷卻、輻射冷卻等方式來達(dá)到熱防護(hù)的目的。