基于有限元方法的液體火箭發(fā)動(dòng)機(jī)主動(dòng)冷卻技術(shù)研究

安國(guó)琛，李仁俊，臧月進(jìn)，李凌黎

(上海機(jī)電工程研究所，上海，201109)

0　引　言

液體火箭發(fā)動(dòng)機(jī)因其可以重復(fù)使用、比沖高等特點(diǎn)在助推級(jí)、姿軌控發(fā)動(dòng)機(jī)等多領(lǐng)域有越來(lái)越廣泛的應(yīng)用，對(duì)發(fā)動(dòng)機(jī)工作時(shí)長(zhǎng)、比沖等方面提出了更高的要求，進(jìn)而在熱防護(hù)設(shè)計(jì)方面也帶來(lái)了更大的挑戰(zhàn)[1-2]。主動(dòng)冷卻因其穩(wěn)定的冷卻性能在液體火箭發(fā)動(dòng)機(jī)上得到了廣泛的應(yīng)用[3]。本文以一種采用主動(dòng)冷卻的氣氧/煤油火箭發(fā)動(dòng)機(jī)推力室為例，通過(guò)有限元方法研究其內(nèi)流場(chǎng)結(jié)構(gòu)、室壁溫度及熱流密度分布、力學(xué)性能等方面，提出了一種具有普適意義的主動(dòng)冷卻研究方法。該方法能在同一算例內(nèi)考慮推進(jìn)劑霧化燃燒、流固耦合的影響，縮短主動(dòng)冷卻的研究周期。

1　數(shù)學(xué)模型及算例校驗(yàn)

1.1　熱交換過(guò)程數(shù)學(xué)模型

主動(dòng)熱防護(hù)整個(gè)換熱過(guò)程包括燃?xì)庀蚴冶诘膶?duì)流和輻射換熱、通過(guò)室壁的熱傳導(dǎo)以及室壁向冷卻劑的對(duì)流傳熱三個(gè)換熱環(huán)節(jié)[4-5]。當(dāng)傳熱處于穩(wěn)定時(shí)，傳熱過(guò)程如圖1。由能量守恒得傳熱過(guò)程計(jì)算公式(1)[6]。

圖1　傳熱過(guò)程Fig. 1　Heat transfer process

(1)

其中，q為總熱流密度，W/m2；q0為燃?xì)鈱?duì)流換熱傳遞給室壁的熱流密度，W/m2；qr為輻射熱流密度，W/m2；hg為燃?xì)馀c室壁間的總對(duì)流傳熱系數(shù)，W/(m2·K)；Taw為燃?xì)鉁囟?，K；Twg為燃?xì)鈧?cè)室壁溫度，K；λc為室壁的導(dǎo)熱系數(shù)，W/(m·K)；δc為室壁厚度，m；h1為冷卻劑與室壁的對(duì)流傳熱系數(shù)，W/(m2·K)；Twc為冷卻劑側(cè)室壁溫度，K；Tco為冷卻劑溫度，K。

主動(dòng)冷卻主要是通過(guò)冷卻劑的流量、流速和流道結(jié)構(gòu)形式的選擇來(lái)保證主動(dòng)冷卻結(jié)構(gòu)滿足以下條件：

1) 冷卻劑在吸收熱量后，其溫度不應(yīng)超過(guò)該冷卻劑的容許溫度。對(duì)于冷卻水而言，不超過(guò)水在該壓強(qiáng)下的沸點(diǎn)，或者說(shuō)防止壁面附近發(fā)生膜沸騰。

2) 燃?xì)鈧?cè)推力室壁面溫度不得超過(guò)室壁材料熱穩(wěn)定性條件所容許的溫度。

1.2　耦合傳熱及化學(xué)動(dòng)力模型

耦合傳熱方面，模擬發(fā)動(dòng)機(jī)中主動(dòng)冷卻過(guò)程的目的在于得到熱平衡狀態(tài)的流固耦合場(chǎng)，從而分析發(fā)動(dòng)機(jī)的長(zhǎng)時(shí)間工作特性。按穩(wěn)態(tài)傳熱過(guò)程進(jìn)行計(jì)算，求解耦合傳熱問(wèn)題的關(guān)鍵在于選取邊界條件。在燃?xì)饬鲃?dòng)過(guò)程中，燃燒室內(nèi)壁面的溫度與熱流分布均未知，需要通過(guò)燃?xì)鈱?duì)燃燒室壁面的強(qiáng)迫對(duì)流、燃燒室壁面的熱傳導(dǎo)以及冷卻液與冷卻通道的對(duì)流換熱耦合求解得到。流動(dòng)與傳熱之間的耦合發(fā)生在燃?xì)鈧?cè)壁面和冷卻劑側(cè)壁面。滿足耦合交界面條件為：溫度連續(xù)，熱流密度連續(xù)。

組分輸運(yùn)及化學(xué)動(dòng)力方面，燃燒過(guò)程是湍流和化學(xué)反應(yīng)相互作用的結(jié)果，常用的湍流/化學(xué)反應(yīng)相互作用模型有渦耗散模型、有限速率/渦耗散模型。渦耗散模型基于化學(xué)反應(yīng)速率快于流動(dòng)的傳遞過(guò)程的假設(shè)，用于說(shuō)明反應(yīng)速率受湍流控制的情況。對(duì)于計(jì)算區(qū)域內(nèi)既存在較強(qiáng)的湍流混合，又存在動(dòng)力學(xué)控制燃燒過(guò)程的情況，常采用有限速率/渦耗散模型。考慮到本文模擬對(duì)象的燃燒流動(dòng)過(guò)程較復(fù)雜，采用有限速率/渦耗散模型。

1.3　算例校驗(yàn)

主動(dòng)冷卻模擬過(guò)程主要考慮流固耦合的換熱過(guò)程。為了驗(yàn)證下文數(shù)值模擬的合理性，對(duì)換熱過(guò)程進(jìn)行算例校驗(yàn)。模擬文獻(xiàn)[7]中的管內(nèi)對(duì)流換熱過(guò)程，并與文獻(xiàn)中試驗(yàn)結(jié)果進(jìn)行對(duì)比。數(shù)值模擬區(qū)域包括圓管以及外場(chǎng)，環(huán)境溫度290.5 K，圓管與外場(chǎng)兩相交界面為耦合型WALL，圓管內(nèi)壁面設(shè)置30 W熱流。計(jì)算結(jié)果如圖2，結(jié)果表明壁面平均溫度為303.7 K與文獻(xiàn)中試驗(yàn)測(cè)量結(jié)果303 K一致，可認(rèn)為本文采用的方法可以準(zhǔn)確模擬流固耦合的換熱過(guò)程。

圖2　算例結(jié)果Fig. 2　Calculation results

2　模型設(shè)置與網(wǎng)格劃分

本文主動(dòng)冷卻設(shè)計(jì)主要針對(duì)推力室等直段進(jìn)行。等直段室壁采用銅合金內(nèi)襯加鋼外殼的結(jié)構(gòu)，通過(guò)在銅合金內(nèi)襯外側(cè)刻槽來(lái)形成冷卻流道。冷卻水由頭部流進(jìn)而后從推力室等直段末端流出。為了使冷卻水均勻的進(jìn)入冷卻通道，在冷卻水入口設(shè)計(jì)圓環(huán)結(jié)構(gòu)的集液腔，冷卻水先進(jìn)入集液腔外環(huán)而后溢入環(huán)內(nèi)并分配進(jìn)各冷卻通道。

結(jié)構(gòu)方面，推力室長(zhǎng)度為200 mm，內(nèi)徑62 mm，噴管喉徑12 mm，頭部為三個(gè)沿周向均布的噴嘴，推進(jìn)劑為氣氧/煤油。冷卻通道采用30個(gè)5 mm×5 mm正方形截面的等直流道。室壁內(nèi)側(cè)厚度為2 mm。冷卻通道結(jié)構(gòu)如圖3。推進(jìn)劑總流量400 g/s，氧氣/煤油流量比為1，冷卻水流量為5 kg/s，銅合金導(dǎo)熱系數(shù)為352 W/(m·k)。

計(jì)算區(qū)域方面，包括三部分：推力室內(nèi)流場(chǎng)、室壁以及冷卻水。計(jì)算過(guò)程中將推力室內(nèi)流場(chǎng)與室壁以及冷卻劑耦合在一起。網(wǎng)格劃分方面，采用非結(jié)構(gòu)網(wǎng)格，在兩相交界面進(jìn)行網(wǎng)格加密。網(wǎng)格總數(shù)65萬(wàn)，如圖4。

模型選取方面，計(jì)算過(guò)程考慮以下過(guò)程：燃料霧化燃燒，室壁與燃?xì)庖约袄鋮s劑之間熱量交換。采用DPM模型模擬從噴嘴射入的煤油液滴顆粒，采用有限速率/渦耗散模型模擬氧氣/煤油反應(yīng)[8]，氣/固以及固/液耦合面設(shè)置為Coupled類(lèi)型Wall。

(a) 推力室等直段(a) The straight section of chamber

(b) 冷卻通道橫截面(b) The cross-section of cooling channel

圖4　網(wǎng)格劃分Fig4　Grid division

3　數(shù)值模擬結(jié)果分析

3.1　流固耦合計(jì)算結(jié)果分析

圖5(a)中噴嘴射流下游的煤油液滴蒸發(fā)燃燒產(chǎn)生局部高溫約3 000 K，主流區(qū)域溫度在2 000～2 500 K。近壁低溫區(qū)域溫度為1 000～1 500 K，燃?xì)馀c室壁對(duì)流傳熱，熱量通過(guò)室溫度壁傳遞給冷卻劑。存在主動(dòng)冷卻的情況下，推力室內(nèi)壁面附近流場(chǎng)形成明顯的溫度梯度。圖5(b)為推力室內(nèi)流場(chǎng)流線圖，可以看出推力室內(nèi)流場(chǎng)上游有明顯的回流區(qū)域，推力室具備較好的穩(wěn)定持續(xù)燃燒條件，同時(shí)回流也促進(jìn)了高溫燃?xì)馀c推力室內(nèi)壁的對(duì)流換熱過(guò)程。推力室內(nèi)大部分區(qū)域流速很低，噴嘴射流以及推力室收斂段流速較高，但均未達(dá)到音速，噴嘴射流為亞音速流。計(jì)算得推力室壓力為4.7 MPa。

圖6(a)為推力室內(nèi)壁面熱流密度云圖，其中高熱流密度區(qū)域位于與噴嘴射流相對(duì)應(yīng)的位置，與反應(yīng)放熱區(qū)域一致，熱流密度最高達(dá)到107W/m2。圖6(b)為推力室內(nèi)壁面溫度云圖，其平均溫度為483 K，最高溫度不超過(guò)600 K，高溫區(qū)域位于推力室中部區(qū)域。圖6(c)為冷卻水在冷卻通道內(nèi)的流速，出入口位置流速最高約20 m/s。30個(gè)5 mm×5 mm正方形通道內(nèi)冷卻水流速均勻，約5 m/s，冷卻結(jié)構(gòu)在冷卻水分配及流速控制上效果較好。圖6(d)為冷卻通道壁面溫度，冷卻通道內(nèi)為流動(dòng)的冷卻水，冷卻通道壁面平均壁溫為384 K，溫度最高約530 K。計(jì)算得冷卻通道內(nèi)壓力約2.5 MPa，此時(shí)冷卻水不易發(fā)生膜沸騰。

(a) 推力室內(nèi)流場(chǎng)溫度分布(a) Temperature distribution of flow field　　(b) 內(nèi)流場(chǎng)馬赫數(shù)及流線圖(b) The mach number and streamline圖5　推力室內(nèi)流場(chǎng)云圖Fig. 5　The picture of chamber flow field

(a) 推力室內(nèi)壁面熱流密度(a) Chamber wall heat flux density　　　　　　(b) 推力室內(nèi)壁面溫度分布(b) Temperature of chamber wall

(c) 冷卻通道內(nèi)冷卻水流速(c) Flow of cooling water　　　(d) 冷卻通道壁面溫度分布(d) Temperature of channel wall圖6　推力室壁面及冷卻通道云圖Fig. 6　The picture of chamber wall and cooling channel

3.2　冷卻結(jié)構(gòu)力學(xué)分析

銅合金內(nèi)襯為冷卻結(jié)構(gòu)的核心部件，具有壁薄且內(nèi)外承壓的特點(diǎn)，因此針對(duì)內(nèi)襯進(jìn)行靜力學(xué)分析。因?yàn)橛糜诠潭▋?nèi)襯的鋼外殼較厚，因而應(yīng)力分析時(shí)在銅合金內(nèi)襯兩端以及肋的上端面添加固定約束。銅襯套內(nèi)外兩側(cè)承受壓力分別為4.7 MPa和2.5 MPa。網(wǎng)格及等效應(yīng)力如圖7。圖7(b)中銅合金內(nèi)襯等效應(yīng)力最高為250 MPa，應(yīng)力分布均勻，沒(méi)有明顯的應(yīng)力集中。

(a) 冷卻通道內(nèi)襯網(wǎng)格(a) Grid of cooling channel　　　　(b) 力學(xué)分析云圖(b) The picture of stress圖7　冷卻通道力學(xué)分析Fig. 7　Mechanics analysis of cooling channel

4　結(jié)束語(yǔ)

采用考慮霧化燃燒過(guò)程的流固耦合數(shù)值模擬方法，以本文所述的液體火箭發(fā)動(dòng)機(jī)推力室為例，分析得到了內(nèi)流場(chǎng)結(jié)構(gòu)、室壁溫度及熱流密度分布以及力學(xué)特性等結(jié)果。文中推力室采用銅合金內(nèi)襯加鋼外殼的冷卻結(jié)構(gòu)，通過(guò)在銅合金內(nèi)襯外側(cè)刻槽來(lái)形成冷卻流道，分析得到設(shè)計(jì)工況下燃?xì)鈧?cè)室壁溫度不高于600 K，銅合金內(nèi)襯應(yīng)力不高于250 MPa，冷卻水在通道內(nèi)均勻流動(dòng)。本文采用的流固耦合換熱計(jì)算方法經(jīng)過(guò)算例驗(yàn)證，其結(jié)果有較高的可靠性，是主動(dòng)冷卻設(shè)計(jì)的一種可行數(shù)值模擬方法。該方法可以同時(shí)考慮推力室內(nèi)流場(chǎng)結(jié)構(gòu)、推進(jìn)劑霧化燃燒以及流固耦合、對(duì)流換熱的綜合影響，進(jìn)而節(jié)省試驗(yàn)時(shí)間、縮短研制周期。