閆 波，張會(huì)強(qiáng)，王 兵

(1.清華大學(xué) 航天航空學(xué)院，北京 100084；2.中國(guó)運(yùn)載火箭技術(shù)研究院 研究發(fā)展中心，北京 100076)

0 引言

空間探測(cè)是航天活動(dòng)的一個(gè)熱點(diǎn)方向[1-2]，而空間探測(cè)離不開小推力發(fā)動(dòng)機(jī)。小推力發(fā)動(dòng)機(jī)的工作環(huán)境通常是太空，所以需要進(jìn)一步考慮太空星體對(duì)發(fā)動(dòng)機(jī)的影響，特別是太陽(yáng)輻射對(duì)模型的非均勻性熱影響。當(dāng)飛行器在空間預(yù)定軌道飛行時(shí)，由于飛行器相對(duì)于太陽(yáng)和地球的位置和方向發(fā)生著不斷的變化，其外熱流和輻射交換也發(fā)生著改變，從而使得發(fā)動(dòng)機(jī)遭受高低溫變化的影響。因此，對(duì)小推力液體火箭發(fā)動(dòng)機(jī)進(jìn)行空間在軌熱分析是十分有意義的。

本文利用ANSYS Workbench 有限元分析軟件，對(duì)小推力發(fā)動(dòng)機(jī)推力室在地球同步軌道(GEO)的結(jié)構(gòu)熱特性開展了分析研究，研究成果為小推力液體火箭發(fā)動(dòng)機(jī)結(jié)構(gòu)熱設(shè)計(jì)提供了依據(jù)。

1 GEO軌道外熱流計(jì)算方法

空間飛行器繞地球軌道運(yùn)行時(shí)，受到太陽(yáng)直接輻射、地球反照和地球的紅外輻射，這三部分通常稱為空間外熱流(簡(jiǎn)稱外熱流)，屬于空間飛行器的外部熱環(huán)境[3-5]。

當(dāng)空間飛行器在約200 km以上高度的低地球軌道(LEO)上飛行時(shí)，其接受的空間外熱流主要是太陽(yáng)輻射、地球反照輻射和地球紅外輻射三種，對(duì)于GEO軌道，地球反射和紅外輻射可以忽略不計(jì)，空間飛行器接受的空間外熱流主要是太陽(yáng)輻射[5]。以下只針對(duì)太陽(yáng)輻射進(jìn)行分析。

從LEO至GEO軌道的高度上，太陽(yáng)光被認(rèn)為是均勻的平行光束，太陽(yáng)輻射強(qiáng)度為太陽(yáng)常數(shù)S，目前規(guī)定S=(1 353±21)W/m2。

空間飛行器外表面任一微元面積dA上(如圖1)受到的太陽(yáng)輻射外熱流為

dq1=ScosθdA

(1)

令

φ1=dq1/SdA

則

φ1=cosθ

式中φ1為太陽(yáng)輻射角系數(shù)。任一面元的太陽(yáng)輻射外熱流為

dq1=αsSφ1dA

(2)

式中αs為面元表面的太陽(yáng)吸收率，太陽(yáng)輻射角系數(shù)需要分析空間飛行器在空間的運(yùn)動(dòng)規(guī)律，以確定微元表面與地球、太陽(yáng)的相對(duì)關(guān)系[5]。

圖1 太陽(yáng)輻射熱流圖Fig.1 Heat flow diagram of solar radiation

綜上所述，本文對(duì)熱源的基本假定如下：

1)太陽(yáng)輻射計(jì)算時(shí)采用的太陽(yáng)常數(shù)為1 353 W/m2(熱流密度)；

2)只考慮太陽(yáng)對(duì)發(fā)動(dòng)機(jī)的直射和斜射，不考慮太陽(yáng)的散射，太陽(yáng)光為平行光；

3)不考慮地球的紅外輻射，地球反照；

4)忽略空間飛行器對(duì)發(fā)動(dòng)機(jī)的輻射和導(dǎo)熱，不計(jì)其他行星的熱輻射；

5)外層空間是絕對(duì)黑體；

6)只考慮推力室自身遮擋光線，不考慮其他物體遮擋[6]。

2 有限元數(shù)值仿真模型

為了研究太陽(yáng)輻射對(duì)推力室結(jié)構(gòu)熱特性的影響，本文采用有限元法對(duì)深空工作環(huán)境下推力室穩(wěn)態(tài)工作時(shí)和不工作時(shí)分別進(jìn)行了模擬計(jì)算[7-8]，兩種分析都采用三維穩(wěn)態(tài)熱分析模型。

2.1 推力室工作參數(shù)

本文主要研究對(duì)象是某小推力液體發(fā)動(dòng)機(jī)，該型發(fā)動(dòng)機(jī)相關(guān)參數(shù)如表1所示，發(fā)動(dòng)機(jī)額定推力為445 N，N2O4為氧化劑，混肼-50或一甲基肼(MMH)為燃料，采用液膜加輻射組合的冷卻方式。發(fā)動(dòng)機(jī)有較長(zhǎng)的工作壽命，穩(wěn)態(tài)工作時(shí)間大于500 s。該發(fā)動(dòng)機(jī)由推力室(噴注器、燃燒室、噴管)、各種閥門、調(diào)節(jié)器及機(jī)架等各部分組成[9]。

表1 發(fā)動(dòng)機(jī)綜合參數(shù)Tab.1 Integrated parameters of engine

2.2 定義材料屬性及單元類型

推力室三維實(shí)體模型在Pro/E中創(chuàng)建，通過(guò)CAD和ANSYS Workbench的無(wú)縫連接將模型導(dǎo)入。建立三維實(shí)體模型后，需要定義推力室結(jié)構(gòu)材料屬性，包括材料密度、比熱、導(dǎo)熱等參數(shù)。

單元類型選用ANSYS 軟件中的solid87，單元尺寸取為0.002 m，因?yàn)橥屏κ抑型瑫r(shí)存在熱對(duì)流和熱輻射換熱方式，為了避免施加載荷的覆蓋情況出現(xiàn)，這里要在內(nèi)表面和外表面各設(shè)定一個(gè)表面效應(yīng)單元SURF152。在內(nèi)表面把燃?xì)廨椛錈崃髅芏仁┘釉趯?shí)體單元上，把對(duì)流換熱施加在表面效應(yīng)單元上，在外表面把太陽(yáng)對(duì)推力室的輻射施加在表面效應(yīng)單元上。

2.3 邊界條件與加載過(guò)程

對(duì)于小推力液體火箭發(fā)動(dòng)機(jī)，由于受發(fā)動(dòng)機(jī)結(jié)構(gòu)本身及擠壓式供應(yīng)系統(tǒng)供應(yīng)壓力的限制，一般采用液膜冷卻結(jié)合輻射冷卻的方法。同時(shí)，其中包含的傳熱過(guò)程有燃?xì)鈱?duì)結(jié)構(gòu)及結(jié)構(gòu)和環(huán)境的輻射換熱、燃?xì)鈱?duì)結(jié)構(gòu)及推進(jìn)劑對(duì)結(jié)構(gòu)的對(duì)流換熱和結(jié)構(gòu)內(nèi)部的熱傳導(dǎo)。小推力發(fā)動(dòng)機(jī)的工作過(guò)程，涉及到三種基本傳熱方式，即熱傳導(dǎo)、熱對(duì)流和熱輻射。

加載載荷包括：模型初場(chǎng)溫度分布、工作環(huán)境溫度、燃?xì)廨椛錈崃髅芏群腿細(xì)饣謴?fù)溫度、流體與壁面間對(duì)流換熱系數(shù)等。本文全場(chǎng)發(fā)射率均設(shè)置為0.9，初場(chǎng)溫度分布均為298 K，工作環(huán)境溫度為4 K。利用巴茲公式和發(fā)動(dòng)機(jī)的熱力計(jì)算數(shù)據(jù)，得到了推力室內(nèi)的燃?xì)獗诿姹砻鎸?duì)流換熱系數(shù)、燃?xì)饣謴?fù)溫度和燃?xì)鉁囟鹊姆植?，如圖2和圖3中所示[10-11]。噴注面上的冷卻劑孔噴出冷卻液，在燃燒室前段形成冷卻液膜，液膜長(zhǎng)度約為燃燒室的一半，液膜溫度約為430 K。因?yàn)橐耗ず穸群苄?，且其透射率較高，故忽略燃?xì)夂捅诿鎸?duì)液膜的輻射傳熱過(guò)程。噴注器前端面和噴管末端面絕熱處理。

圖2 燃?xì)?液膜)沿軸向溫度分布Fig.2 Distribution of gas(liquid film)temperature along the axial

圖3 燃?xì)?液膜)沿壁面對(duì)流換熱系數(shù)分布Fig.3 Distribution of gas(liquid film)convection heat transfer coefficient along the wall

外壁對(duì)太空的輻射采用軟件自帶輻射模型施加。在太空中，由于不同時(shí)刻太陽(yáng)光線照射到推力室的位置不同，所以不同時(shí)刻施加到發(fā)動(dòng)機(jī)表面的熱量也不同，直射時(shí)通過(guò)單元熱流密度為q=S，S為太陽(yáng)常數(shù)，斜射時(shí)q=Scosθ，太陽(yáng)輻射施加的載荷類型為熱流密度[6]。太陽(yáng)光在坐標(biāo)系中方向余弦為(cosα，cosβ，cosγ)，在ANSYS Workbench 下插入APDL 語(yǔ)言加載燃?xì)廨椛?、?duì)流載荷和太陽(yáng)輻射，將燃?xì)廨椛錈崃髅芏?、燃?xì)饣謴?fù)溫度和對(duì)流換熱系數(shù)設(shè)置為表格數(shù)組加載，把太陽(yáng)光方向余弦設(shè)為參數(shù)變量，可根據(jù)發(fā)動(dòng)機(jī)在太空中所處的方位設(shè)計(jì)不同的太陽(yáng)輻射角度。本節(jié)假設(shè)太陽(yáng)沿x軸垂直照射，即方向余弦為(-1，0，0)，如圖4所示。

圖4 太陽(yáng)輻射方向Fig.4 Direction of solar radiation

噴注面有燃?xì)廨椛?，需要施加輻射邊界條件，本文計(jì)算模型中對(duì)噴注面輻射的燃?xì)鉁囟热? 000 K是合適的。需要對(duì)噴注器管道施加對(duì)流換熱邊界條件，需要得到管道對(duì)流換熱系數(shù)和流體溫度對(duì)流換熱系數(shù)，如表2所示。

表2 噴注器管道對(duì)流換熱系數(shù)Tab.2 Convective heat transfer coefficient of injector pipe

推進(jìn)劑入口溫度取290 K。但在施加邊界條件時(shí)，為簡(jiǎn)便起見，各段流體溫度設(shè)置為平均溫度，氧化劑主噴孔、燃料主噴孔、冷卻劑孔內(nèi)流體溫度取320 K[12]。噴注器的輻射和對(duì)流換熱采用軟件自帶模型施加。

推力室不工作工況只有推力室對(duì)太空的輻射和太陽(yáng)對(duì)推力室的輻射，不計(jì)地球反照輻射、地球紅外輻射和其他行星的熱輻射。

3 數(shù)值仿真結(jié)果分析

3.1 發(fā)動(dòng)機(jī)穩(wěn)態(tài)工作工況

經(jīng)模擬計(jì)算得到有太陽(yáng)輻射和沒(méi)有太陽(yáng)輻射推力室室壁溫度分布規(guī)律，有太陽(yáng)輻射推力室溫度分布云圖如圖5所示，圖6為沒(méi)有太陽(yáng)輻射推力室穩(wěn)態(tài)工作時(shí)溫度分布。

圖5 有太陽(yáng)輻射穩(wěn)態(tài)工作推力室壁溫分布Fig.5 Wall temperature distribution of thrust chamber with solar radiation in steady working

圖6 無(wú)太陽(yáng)輻射穩(wěn)態(tài)工作推力室壁溫分布Fig.6 Wall temperature distribution of thrust chamber without solar radiation in steady working

有太陽(yáng)輻射時(shí)選取如圖5中兩條曲線上的溫度進(jìn)行分析，經(jīng)過(guò)分析比較發(fā)現(xiàn)，曲線1溫度略高于曲線2溫度，最大溫差出現(xiàn)在噴管出口處，溫度相差10 K左右，其余位置溫度相差在3 K以下。對(duì)有太陽(yáng)輻射和沒(méi)有太陽(yáng)輻射推力室室壁溫度分布云圖比較可知，溫度分布和大小幾乎一致，噴管喉部溫差在1 K以下，噴管出口處相差最大，最大差值10 K左右。由此可以得出結(jié)論，太陽(yáng)輻射對(duì)推力室穩(wěn)態(tài)工作時(shí)溫度有一定影響，但影響不大。

3.2 發(fā)動(dòng)機(jī)不工作工況

推力室發(fā)動(dòng)機(jī)不工作工況只有推力室對(duì)太空的輻射和太陽(yáng)對(duì)推力室的輻射，模擬計(jì)算得到推力室溫度分布云圖如圖7所示。推力室最高溫度在太陽(yáng)直射面噴管尾部為338.81 K，由于有熱傳導(dǎo)作用，溫度最小值不是太空溫度4 K而是267.27 K，最小值在推力室背光面的噴管尾部，推力室溫度沿負(fù)x軸方向逐漸遞減。

如圖4所示在推力室身部?jī)?nèi)壁面兩條母線上的不同位置選取了6個(gè)采樣點(diǎn)，分別位于x1=-0.074 m，x2=-0.022 m，x3=0 m，x4=0.070 m，x5=0.123 m，x6=0.186 m處。各處溫度變化曲線如圖8所示。從圖8中可以看出曲線1上的溫度變化趨勢(shì)為溫度先上升后有小幅度下降，過(guò)了喉部再次上升，噴管出口處達(dá)到最大值；曲線2上的溫度一直保持下降變化趨勢(shì)，噴管出口處達(dá)到最小值。從兩曲線溫度差值可以看出，曲線1和曲線2溫度差值在燃燒室和噴注器連接處最小，只相差4.15 K，然后溫差開始逐漸增加直到燃燒室結(jié)束，而后溫度差又開始減少，喉部溫度差為8.95 K，過(guò)了喉部溫度差開始逐漸增加，到達(dá)噴管尾部達(dá)到最大溫度差值71.41 K。

圖7 有太陽(yáng)輻射發(fā)動(dòng)機(jī)不工作時(shí)推力室壁溫分布Fig.7 Wall temperature distribution of thrust chamber with solar radiation under not working condition

圖8 推力室壁溫沿軸向變化曲線Fig.8 Changing curve of thrust chamber wall temperature along the axial

由此可以得出結(jié)論，當(dāng)發(fā)動(dòng)機(jī)穩(wěn)態(tài)工作時(shí)太陽(yáng)輻射對(duì)推力室溫度有一定影響，但影響不大。而當(dāng)發(fā)動(dòng)機(jī)不工作時(shí)，由于小推力發(fā)動(dòng)機(jī)工作環(huán)境是太空，發(fā)動(dòng)機(jī)機(jī)體一半接受太陽(yáng)輻射，溫度較高，一半面對(duì)深空冷環(huán)境，溫度較低，太陽(yáng)輻射對(duì)模型的非均勻性影響較大。

4 結(jié)論

對(duì)考慮太陽(yáng)輻射的深空工作環(huán)境下推力室穩(wěn)態(tài)工作時(shí)和發(fā)動(dòng)機(jī)不工作時(shí)分別進(jìn)行了模擬計(jì)算，太陽(yáng)輻射沿x軸方向垂直照射推力室，兩種分析都采用三維穩(wěn)態(tài)熱分析模型。

1)當(dāng)發(fā)動(dòng)機(jī)穩(wěn)態(tài)工作時(shí)太陽(yáng)輻射對(duì)推力室溫度有一定影響，但影響不大。有太陽(yáng)輻射時(shí)推力室最大溫差出現(xiàn)在直射母線和背光母線噴管出口處，溫度相差10 K左右，其余位置溫度相差在3 K以下。有太陽(yáng)輻射和無(wú)輻射室壁溫度分布和大小幾乎一致，噴管喉部溫差在1 K以下，噴管出口處相差最大，最大差值10 K左右。

2)發(fā)動(dòng)機(jī)不工作工況只有推力室對(duì)太空的輻射和太陽(yáng)對(duì)推力室的輻射，最高溫度在太陽(yáng)直射面噴管尾部為338.81 K，由于有熱傳導(dǎo)作用，溫度最小值不是太空溫度4 K而是267.27 K，最小值在推力室背光面的噴管尾部，推力室溫度沿負(fù)x軸方向逐漸遞減。當(dāng)發(fā)動(dòng)機(jī)不工作時(shí)，由于小推力發(fā)動(dòng)機(jī)工作環(huán)境是太空，發(fā)動(dòng)機(jī)機(jī)體一半接受太陽(yáng)輻射，溫度較高，一半面對(duì)深空冷環(huán)境，溫度較低，太陽(yáng)輻射對(duì)模型溫度分布均勻性有一定影響。