有翼再入飛行器氣動外形集成設(shè)計優(yōu)化

李正洲，賀元元，高昌，張小慶，王琪

1. 中國空氣動力研究與發(fā)展中心 高超聲速沖壓發(fā)動機技術(shù)重點實驗室，綿陽 621000 2. 中國空氣動力研究與發(fā)展中心，綿陽 621000

氣動布局決定了飛行器的氣動性能，繼而影響飛行器的飛行性能、操穩(wěn)性能、防熱性能和有效載荷等，因此氣動布局設(shè)計在飛行器設(shè)計中具有舉足輕重的作用，良好的氣動布局是成功研制飛行器的基礎(chǔ)和前提條件[1-2]。有翼再入飛行器(Winged Re-entry Vehicles, RV-W)具有較優(yōu)的氣動性能，能夠更好地滿足航天運輸中的機動飛行/載荷投送等任務(wù),并具有潛在的可重復(fù)使用能力，已成為當(dāng)前各航天大國天地往返運輸系統(tǒng)中的研究熱點[3-4]。

再入環(huán)境決定了飛行器的氣動外形設(shè)計[5]。在有翼再入飛行器的再入飛行過程中，高度從100 km降至0 km，馬赫數(shù)從28降至0，雷諾數(shù)的量級在104～108范圍內(nèi)變化，迎角約從40°降至0°，如此劇烈的參數(shù)變化會導(dǎo)致氣動力、氣動熱環(huán)境的劇烈變化[6]，飛行器氣動外形必須適應(yīng)嚴(yán)酷的氣動力、熱環(huán)境[7]。

相對于彈道式再入航天器、半彈道式航天器，有翼再入飛行器氣動外形設(shè)計的特殊性主要體現(xiàn)在3方面[8]：① 再入減速,利用氣動外形的阻力特性，通過返回軌道設(shè)計，實現(xiàn)將極高的進入速度減小到能安全著陸的速度；② 飛行可控,通過氣動升力及對升力在軌道面分量的控制，實現(xiàn)再入飛行的軌跡、最大過載和落點精度控制；③熱防護安全,通過氣動外形避免產(chǎn)生過大的局部熱流設(shè)計，配置熱防護結(jié)構(gòu)保證氣動加熱安全。

為使機動飛行可控，良好的氣動特性成為有翼再入飛行器氣動外形的首先設(shè)計目標(biāo)。氣動外形設(shè)計要求在滿足內(nèi)/外尺寸約束下，機身盡量扁平且細長以減小阻力，機翼盡量增大展弦比和機翼面積以而增大升力系數(shù)。然而在高速再入過程中，飛行器將受到嚴(yán)重的氣動加熱，氣動熱環(huán)境除了與再入軌跡有關(guān)外，還與飛行器的氣動外形密切相關(guān)。從防熱的目的來看，有翼再入飛行器通常設(shè)計為具有鈍前緣機身與鈍前緣、大后掠角機翼等外形特征，這與氣動特性需要的外形布局矛盾。此外，有翼再入飛行器在再入軌跡上長時間以大迎角飛行，為使其底部迎風(fēng)面高效地進行飛行器表面輻射冷卻，其下表面通常設(shè)計為扁平狀，從而進一步增加了氣動阻力。再入過程中為確保宇航員生命安全和飛行器的儀器安全，要求有翼再入飛行器的熱防護系統(tǒng)(Thermal Protection System, TPS)隔熱性能要好，TPS既要能承受峰值加熱，又要能承受總加熱，還要最大限度地減少自身的重量，以增加有效載荷。

從上述分析可以看出，有翼再入飛行器的飛行速域、空域、時域跨度大，氣動力、氣動熱、熱防護等多學(xué)科緊密耦合，氣動布局綜合設(shè)計問題突出，氣動外形設(shè)計的復(fù)雜程度和難度都非常大。因此，必須在總體設(shè)計初期即采用氣動力/氣動熱/熱防護的集成設(shè)計優(yōu)化，解決力、熱設(shè)計矛盾，才能實現(xiàn)飛行器總體性能最優(yōu)，切實降低后期出現(xiàn)方案顛覆性問題的可能性[9]。

本文針對有翼再入飛行器氣動外形設(shè)計問題，首先探討了氣動外形設(shè)計規(guī)律與準(zhǔn)則；在此基礎(chǔ)上以類X-37B飛行器為研究對象，對其氣動布局進行集成設(shè)計優(yōu)化；發(fā)展了一種熱防護系統(tǒng)輕量化設(shè)計優(yōu)化方法，并采用該方法對類X-37B優(yōu)化構(gòu)型的熱防護系統(tǒng)進行了輕量化設(shè)計優(yōu)化。

1 氣動外形設(shè)計規(guī)律與準(zhǔn)則

對于再入飛行器，從軌跡優(yōu)化來看，主要的再入性能氣動參數(shù)有升阻比K(或用L/D表示)、升力系數(shù)CL以及彈道系數(shù)P=m/(CDA),其中m為飛行器質(zhì)量,CD為阻力系數(shù),A為截面面積。

針對再入飛行的任務(wù)需求，本文首先對氣動特性參數(shù)如最大升阻比、最大升力系數(shù)和彈道系數(shù)對再入飛行性能的影響進行研究，為氣動外形設(shè)計優(yōu)化提出合理建議。

1.1 升阻比對再入飛行的影響

1) 升阻比與航程的關(guān)系

根據(jù)文獻[10]的研究，假設(shè)飛行器以平衡滑翔方式再入飛行，則升阻比K與航程L的關(guān)系式為

(1)

從式(1)中可以看出，飛行器升阻比越高，則可獲得的航程L越大，機動范圍也更大。因此高升阻比再入飛行器能更精確地控制再入段軌道和著陸段軌道，實現(xiàn)在跑道上水平著陸，從而為重復(fù)使用創(chuàng)造了條件。

2) 升阻比與過載的關(guān)系

再入飛行器中過載n與升阻比K有直接關(guān)系[11]：

(2)

式中：B為彈道系數(shù)的倒數(shù)，即質(zhì)量-阻力參數(shù);ρ0為大氣密度;v為速度。

從式(2)可以看出，再入過程中使用升阻比越大，飛行器承受的過載也越大。飛行器再入時必須考慮航天員、儀器能承受過載的峰值，也就是必須考慮再入軌跡中升阻比變化。

3) 升阻比與氣動熱的關(guān)系

圖1[12]為飛行器以最大升阻比姿態(tài)進入大氣層時總加熱量Qs與進入速度ventry的關(guān)系。圖2[12]為飛行器以最大升力系數(shù)姿態(tài)進入大氣層時總加熱量與進入速度的關(guān)系，總加熱量表示進入過程飛行器受到的累積氣動加熱量。以最大升阻比為2.0的飛行器為例(圖中紅線)：當(dāng)飛行器以最大升阻比姿態(tài)進入大氣層時，受到的總加熱量最高可達7×109J/m2，而以最大升力系數(shù)姿態(tài)進入大氣層受到的總加熱量約為2×109J/m2(以14 km/s 進入時)，說明飛行器以最大升阻比進入大氣層時受到的總加熱量要遠大于以最大升力系數(shù)進入。

4) 升阻比與結(jié)構(gòu)質(zhì)量的關(guān)系

圖3為再入飛行器結(jié)構(gòu)重量與最大升阻比的關(guān)系，其中縱軸表示假設(shè)有效載荷相同時，其他飛行器需要的結(jié)構(gòu)重量W與返回艙結(jié)構(gòu)重量W0之比；橫軸表示不同氣動外形飛行器的升阻比變化。從圖中可以看出，為獲得高升阻比，飛行器的變化趨勢為拉長機身，壓縮飛行器橫截面積，機頭由鈍頭體逐漸變尖銳，這樣將使飛行器空間部位的利用率的降低，同時熱防護系統(tǒng)質(zhì)量明顯增加[12]。

圖1 最大升阻比下進入時總加熱量與進入速度的關(guān)系[12]Fig.1 Integrated heat load as a function of entry velocity at maximum lift-to-drag ratio[12]

圖2 最大升力系數(shù)下進入時總加熱量與進入速度的關(guān)系[12]Fig.2 Integrated heat load as a function of entry velocity at maximum lift coefficient[12]

圖3 再入飛行器結(jié)構(gòu)重量與最大升阻比的關(guān)系Fig.3 Re-entry vehicle structural weight as a function of maximum lift-to-drag ratio

1.2 升力系數(shù)對再入飛行的影響

對于推力為零的再入滑翔機動飛行，一般選取升力系數(shù)和滾轉(zhuǎn)角做為控制參數(shù)。升力系數(shù)的取值范圍受最大升力系數(shù)的約束，因此最大升力系數(shù)對再入機動有著很大影響[12]。

圖4為軌道傾角if與最大升力系數(shù)CL,max的關(guān)系。從圖中可以看出隨著最大升力系數(shù)的增加，所能改變的最大軌道傾角也相應(yīng)增加，即機動飛行能力相應(yīng)增加；但當(dāng)最大升力系數(shù)大于一定值時，軌道傾角增加趨勢變緩。

圖4 軌道傾角隨最大升力系數(shù)變化Fig.4 Orbit inclination as a function of maximum lift coefficient

1.3 彈道系數(shù)對再入飛行的影響

彈道系數(shù)又叫質(zhì)阻比，是飛行器質(zhì)量和外形的組合參數(shù),用來衡量飛行器克服空氣阻力、維持飛行速度的能力。

圖5為彈道系數(shù)與升阻特性的變化關(guān)系示意圖。從圖中可以看出，對于再入飛行器，升阻比大于1.5之后，若要獲得更好的升阻特性，則必須減小大幅飛行器受到的氣動阻力，形成尖前緣帶翼構(gòu)型的氣動外形。但這種外形布局的容積利用率也隨之降低[12]。

圖5 升阻比隨彈道系數(shù)變化Fig.5 Lift-to-drag ratios as a function of ballistic coefficients

此外，再入飛行器通常具有很高的初始動能和勢能。飛行器在再入過程中利用氣動阻力耗散大部分能量，才能夠?qū)崿F(xiàn)以較小速度安全返回地面。因此對于有翼再入飛行器，其阻力不能太小，也就是彈道系數(shù)不能過大。

1.4 有翼再入飛行器氣動外形設(shè)計建議

根據(jù)氣動特性參數(shù)對再入飛行性能影響的分析，對有翼再入飛行氣動布局設(shè)計優(yōu)化可提出如下建議：

1) 升阻比是氣動外形設(shè)計的重要參數(shù)，它綜合反映飛行器的氣動性能。但升阻比與過載、氣動加熱、結(jié)構(gòu)質(zhì)量占比、有效容積利用率等有密切關(guān)系。因此對有翼再入飛行器氣動外形設(shè)計優(yōu)化時，不能一味追求最大升阻比，必須對升阻比與其他因素進行綜合權(quán)衡。

2) 氣動布局設(shè)計應(yīng)使其升力系數(shù)保持一個較大的范圍，滿足升力系數(shù)最大準(zhǔn)則。但氣動外形設(shè)計優(yōu)化時并不能總是追求最大升力系數(shù)，而應(yīng)根據(jù)實際再入機動飛行要求設(shè)計合理的最大升力系數(shù)。

3) 彈道系數(shù)體現(xiàn)了飛行器的氣動控制能力，是氣動減速性能、有效容積利用率、防熱性能等的綜合體現(xiàn)。有翼再入飛行器氣動外形設(shè)計優(yōu)化時必須合理選擇飛行器彈道系數(shù)。

2 優(yōu)化策略與流程

本文以類X-37B飛行器為研究對象，通過集成氣動力、氣動熱、熱防護等學(xué)科的設(shè)計優(yōu)化，研究有翼再入飛行器氣動外形的設(shè)計優(yōu)化技術(shù)。

圖6 X-37B飛行器示意圖及CAD建模Fig.6 X-37B analog sketch and modeling by CAD software

X-37B飛行器是美國為了驗證可重復(fù)使用空間技術(shù)和在軌空間飛行任務(wù)而啟動的項目[13]，并計劃在X-37B飛行器的技術(shù)上繼續(xù)進行能夠投送6名宇航員進入太空的X-37C計劃[14]。圖6為X-37B飛行器示意圖及本文采用CAD軟件重建的三維幾何外形。

2.1 集成優(yōu)化模型

在對氣動力、熱進行優(yōu)化時，需要確定飛行條件才進行計算分析。本文選擇高度H=50 km、馬赫數(shù)Ma=15、迎角α=15°與側(cè)滑角β=0°。

通過氣動參數(shù)對飛行性能影響的分析可知，氣動特性參數(shù)和局部熱環(huán)境是再入飛行器氣動外形設(shè)計的關(guān)鍵因素。根據(jù)有翼再入飛行器氣動外形設(shè)計建議，綜合考慮氣動性能、過載、氣動熱、結(jié)構(gòu)質(zhì)量占比、有效容積利用率等要求，第1類優(yōu)化目標(biāo)可以總結(jié)為使飛行器獲得一定范圍內(nèi)的高升阻比與升力系數(shù)，即

(3)

高速再入飛行時，會遇到氣動加熱的問題。因此第2類優(yōu)化目標(biāo)為令駐點熱流密度Qws與熱載荷Φ(氣動加熱在單位在時間內(nèi)傳入飛行器的熱量，W)最小，即

(4)

對以上4個優(yōu)化目標(biāo)f1～f4，通過線性加權(quán)法建立一個綜合指標(biāo)f，這個值越大，表示飛行器氣動性能越優(yōu)異：

(5)

式中：下標(biāo)0表示飛行器初始外形的性能參數(shù)，該項作為基準(zhǔn)值；ki(i=1,2,3,4)為權(quán)重系數(shù)，表示各優(yōu)化目標(biāo)的相對重要程度。

對飛行器熱防護系統(tǒng)設(shè)計優(yōu)化時，需要有確定的再入軌跡。盡管再入軌跡對熱防護重量有顯著影響[15-16]，然而本文的目的是建立有翼再入飛行器的力、熱、熱防護多學(xué)科設(shè)計優(yōu)化體系。因此本文在有翼再入飛行器的再入走廊區(qū)間內(nèi)，選取典型再入軌跡對熱防護系統(tǒng)進行輕量化設(shè)計優(yōu)化，而不對再入軌跡進行優(yōu)化。圖7為在過載、動壓、氣動熱等約束下建立的速度-高度再入走廊[17]，圖8為本文TPS輕量化設(shè)計優(yōu)化采用的再入軌跡[18]，圖中t為時間。

圖7 有翼再入飛行器再入走廊示意圖[17]Fig.7 Re-entry corridor sketch of RV-W[17]

圖8 再入軌跡示意圖[18]Fig.8 Re-entry trajectory sketch[18]

2.2 集成優(yōu)化流程

飛行器氣動外形優(yōu)化主要有3個核心問題：首先要對飛行器的幾何外形進行參數(shù)化建模，參數(shù)化建模方法以能夠適用復(fù)雜飛行器的幾何外形、并且用盡可能少的參數(shù)來描述飛行器的特征為優(yōu)；其次是建立各個學(xué)科的性能評估模型，由于需要進行大規(guī)模的優(yōu)化迭代運算，性能評估模型的計算速度顯得尤為重要，要在精度和速度之間取得平衡，否則工程應(yīng)用價值將大打折扣；最后便是優(yōu)化算法，能夠進行多目標(biāo)、大規(guī)模、高效率的全局優(yōu)化算法也是非常重要的一環(huán)。

飛行器外形參數(shù)化建模是指用一組參數(shù)描述飛機幾何的數(shù)學(xué)模型，根據(jù)幾何參數(shù)化模型生成三維幾何圖形。本文采用基于CAD軟件二次開發(fā)的方法對類X-37B飛行器進行參數(shù)化建模，以飛機典型部件(如機翼、機翼、垂尾等)形狀為目標(biāo)，通過部件參數(shù)(如機翼展弦比、前緣后掠角、根弦長等)來定義出部件的幾何形狀，所采用的參數(shù)在工程上具有明確的意義。類X-37B飛行器參數(shù)化建模外形如圖6所示。

優(yōu)化過程中要用到數(shù)值優(yōu)化方法進行尋優(yōu)?？紤]本次氣動外形優(yōu)化可能存在多峰性、非線性、非連續(xù)性，傳統(tǒng)的梯度優(yōu)化算法和直接搜索法可能無法找到全局最優(yōu)解，因此本文采用多目標(biāo)遺傳算法(NSGA-II)進行數(shù)值優(yōu)化。

本文類X-37B飛行器氣動外形集成設(shè)計優(yōu)化流程如圖9所示。具體可描述為：首先以升阻比、升力系數(shù)、駐點熱流、總價熱量為目標(biāo)，對氣動外形進行優(yōu)化；待飛行器得到綜合性能較優(yōu)氣動布局之后，再對飛行器熱防護系統(tǒng)進行輕量化設(shè)計優(yōu)化。

圖9 氣動外形集成優(yōu)化流程Fig.9 Aeroshape integrated optimization process

3 學(xué)科方法

3.1 氣動力快速預(yù)測

對于氣動力快速預(yù)測，采用飛行器全速域氣動特性快速預(yù)測工具(Rapid Aerodynamics Prediction Tool, RAPT)[19-20]。該工具在流場計算方面采用3種方法：① 亞聲速情況,采用基于自適應(yīng)直角網(wǎng)格求解全速勢方程的有限體積法；② 跨 /超聲速情況,采用基于自適應(yīng)直角網(wǎng)格求解歐拉方程的有限體積法；③ 高超聲速情況,采用基于牛頓撞擊理論的工程面元法，并且上述3種方法均加入黏性修正。圖10為RAPT氣動特性分析流程。

本文氣動外形設(shè)計優(yōu)化主要針對再入飛行的高速段，該階段RAPT采用基于牛頓撞擊理論的工程面元法預(yù)測飛行器的氣動特性。工程面元法的思路是將飛行器表面劃分為迎風(fēng)面和背風(fēng)面，再分別采取不同的經(jīng)驗公式進行計算。迎風(fēng)面、背風(fēng)面的劃分是根據(jù)自由來流條件及面元法矢共同確定的。對于迎風(fēng)面，采用采用活塞理論與修正牛頓理論結(jié)合的方法進行計算[21]；對于背風(fēng)面，則采用普朗特-邁耶膨脹波方法。對于黏性阻力，采用參考溫度法進行修正[22]。此外，工程面元法加入了對“遮擋面元”的判斷，即采用光線投射算法(Ray-casting Algorithm)[23]判斷當(dāng)前面元是否被遮擋。若為遮擋面元，則將該面元作為背風(fēng)面處理。

選取NASA Langley中心對高升力體(HL-20)的風(fēng)洞試驗數(shù)據(jù)[24]考核本文氣動力快速預(yù)測方法的準(zhǔn)確性。該試驗提供了從超聲速到高超聲速(Ma=1.5～10)、從負迎角到正迎角(α=-5°～30°)等大量工況的試驗數(shù)據(jù)。本文選取高升力體在Ma=6的工況作為算例。

圖10 RAPT氣動特性分析流程Fig.10 Aerodynamic characteristic analysis process of RAPT

圖11、圖12分別為本文方法預(yù)測的升力系數(shù)、阻力系數(shù)與風(fēng)洞試驗數(shù)據(jù)對比，圖中紅色曲線為本文方法預(yù)測結(jié)果、黑色曲線為風(fēng)洞試驗數(shù)據(jù)。從圖中可以看出，升力系數(shù)與試驗數(shù)據(jù)總體吻合較好，僅在小迎角時存在較明顯誤差(最大誤差不超過6%)；本文方法預(yù)測的阻力系數(shù)總體上大于風(fēng)洞試驗得出的數(shù)據(jù)(最大誤差發(fā)生在α=0°，誤差為12.8%)，但兩者趨勢一致。從本文方法氣動力預(yù)測結(jié)果精度來看，能夠滿足飛行器概念研究與初步設(shè)計階段的精度要求。

圖11 升力系數(shù)對比Fig.11 Comparison of lift coefficients

圖12 阻力系數(shù)對比Fig.12 Comparison of drag coefficients

3.2 氣動加熱快速預(yù)測

RAPT同樣用于飛行器氣動熱環(huán)境快速預(yù)測[19,25]。圖13為RAPT采用軸對稱比擬參考焓法快速預(yù)測飛行器氣動熱環(huán)境的流程。

校企深度融合，對企業(yè)來說，可從學(xué)校畢業(yè)生中選擇優(yōu)秀畢業(yè)生充實員工隊伍，學(xué)校成為企業(yè)人才的儲備庫；對學(xué)校來說，可以通過企業(yè)了解企業(yè)對人才的需要狀況，從而修訂人才培養(yǎng)方案。由于企業(yè)的參與，突出了學(xué)生專業(yè)技能的培養(yǎng)和素質(zhì)的全面提高，培養(yǎng)學(xué)生的良好職業(yè)道德和職業(yè)能力，實現(xiàn)學(xué)生的專業(yè)技能訓(xùn)練與崗位需求“零距離”對接，達到學(xué)校人才培養(yǎng)與社會試產(chǎn)該需求全面對接，真正實現(xiàn)校企雙贏。

圖13 氣動加熱快速預(yù)測流程Fig.13 Rapid prediction process of aerodynamic heating

軸對稱比擬法的原理是[26]：基于工程方法或無黏數(shù)值解求出無黏流場變量(飛行器邊界層外緣速度、密度、壓力等)，在流線坐標(biāo)系內(nèi)，給出三維邊界層方程，然后在冷壁和小橫向流假設(shè)下，三維邊界層方程可簡化為軸對稱形式的邊界層方程，軸對稱半徑由無黏表面流線的尺度因子代替；此時，沿某一根流線的熱流密度的確定，就相當(dāng)于求解某一零迎角軸對稱物體上的熱流密度。其中，駐點熱流采用Kemp-Riddell公式進行計算[27]，非駐點采用Zoby等建立的軸對稱物體氣動加熱熱流計算公式[28]:

(6)

(7)

式中:ρ、μ、h、v、q分別為密度、黏性系數(shù)、焓、速度、熱流密度；上標(biāo)*為基于參考焓得出的參數(shù)；下標(biāo)w、s、ws、e、aw和∞分別為壁面、駐點、駐點壁面、邊界層外緣、絕熱壁以及自由來流的相關(guān)參數(shù)；Reθ,e為動量厚度雷諾數(shù)；c1為試驗修正系數(shù)；vc為第一宇宙速度；RN為駐點曲率半徑。

采用文獻[29]空天飛機模型對本文氣動熱計算方法進行考核。圖14為空天飛機表面熱流密度分布云圖(qref為駐點熱流參考值)，圖15為空天飛機迎風(fēng)中心線上熱流密度分布對比，圖16為空天飛機軸向x=62%L站位上熱流密度分布對比。從圖中可以看出：本文方法對飛行器駐點、機身迎風(fēng)面等部位氣動熱環(huán)境預(yù)測精度較高；同時本文方法也能捕捉翼前緣等局部氣動加熱較嚴(yán)酷的部位。除了機翼前緣熱流密度比試驗數(shù)據(jù)低約20%以外，其余點都與風(fēng)洞試驗數(shù)據(jù)符合較好。

圖14 空天飛機模型熱流密度云圖Fig.14 Aerodynamic heating contour of aerospace plane

圖15 迎風(fēng)中心線上熱流密度分布Fig.15 Heat flux distribution at windward center line

圖16 空天飛機62%站位上熱流密度分布Fig.16 Heat flux distribution of aerospace plane at 62% position

3.3 熱防護系統(tǒng)輕量化設(shè)計與優(yōu)化

本文TPS輕量化設(shè)計的思想是：以最小質(zhì)量為目標(biāo)，TPS既能承受再入軌跡上最嚴(yán)酷點的瞬時表面氣動加熱、又能承受整個軌跡氣動加熱產(chǎn)生的內(nèi)部熱傳導(dǎo)，還要保證TPS分布的完整性。

基于以上思想，TPS輕量化設(shè)計優(yōu)化流程如圖17所示：對優(yōu)化后的飛行器幾何外形輸出網(wǎng)格文件，選取軌跡上氣動加熱最嚴(yán)酷的軌跡點計算飛行器表面熱平衡溫度；根據(jù)熱平衡溫度對飛行器表面進行分區(qū)、布置熱防護材料層；結(jié)合熱防護材料屬性和再入軌跡進行傳熱分析及熱防護厚度優(yōu)化；待優(yōu)化結(jié)束后輸出熱防護分布、厚度、重量等TPS數(shù)據(jù)。

圖17 TPS輕量化設(shè)計優(yōu)化流程Fig.17 TPS lightweight design optimization process

多層熱防護材料傳熱分析采用簡化傳熱模型。將TPS物理模型離散為一維模型(如圖18所示)進行傳熱分析,Tj為各層溫度，j表示物理模型的分層。

一維拋物型非定常熱傳導(dǎo)方程為

(8)

將多層材料的熱防護系統(tǒng)物理模型離散后，將一維熱傳導(dǎo)方程改寫為統(tǒng)一的標(biāo)準(zhǔn)形式：

(9)

式中：Aj=kj/(ρjcj)。

對上表面施加氣動加熱邊界條件，采用隱式向前差分格式離散；下表面施加絕熱壁邊界條件，采用隱式向后差分格式離散；內(nèi)部節(jié)點采用隱式中心差分格式離散。再對所有節(jié)點構(gòu)造非線性微分方程組，采用Newton-Raphson法求解[30]。圖19 為采用文獻[30]中算例對本文傳熱分析方法進行驗證，從圖中可以看出，內(nèi)、外表面?zhèn)鳠釟v程與文獻結(jié)果一致性較好。

圖18 TPS物理模型及離散模型Fig.18 TPS physical and discrete models

圖19 熱防護材料溫度響應(yīng)歷程Fig.19 Temperature responses history of TPS material

4 設(shè)計優(yōu)化結(jié)果與分析

4.1 幾何參數(shù)對氣動特性的敏感性

對類X-37B飛行器幾何外形進行參數(shù)化描述，如圖20所示。從圖中可以看出，本次參數(shù)化包括展弦比、根梢比、后掠角等共涉及32個參數(shù)。為了減少眾多參數(shù)對優(yōu)化帶來的“維數(shù)災(zāi)難”，首先進行單幾何參數(shù)對氣動特性影響的敏感性進行分析，目的是找出顯著影響性能的參數(shù)，將這些參數(shù)用于設(shè)計優(yōu)化中，從而實現(xiàn)可行、高效的多學(xué)科設(shè)計優(yōu)化。

進行單幾何參數(shù)對氣動特性的敏感性分析，將機身長度、機身高度、機頭半徑、機翼展弦比、機身襟翼長度等個10參數(shù)作為幾何外形參數(shù)對氣動特性影響的顯著因素傳遞給多學(xué)科優(yōu)化模型，同時固定其他參數(shù)，這樣大大減少了優(yōu)化變量，提高了優(yōu)化可行性與效率。

圖20 類X-37B飛行器幾何外形參數(shù)化示意圖Fig.20 X-37B analog geometry parametric modeling sketch

4.2 氣動特性優(yōu)化結(jié)果

圖21為優(yōu)化過程中升阻比收斂歷程，其他優(yōu)化目標(biāo)收斂歷程類似。表1為類X-37B飛行器優(yōu)化前后性能參數(shù)對比。從表中可以看出：升阻比提升最為明顯(12.82%)；15°迎角下升力系數(shù)提高了6.18%，說明飛行器是同時增升、減阻以提高升阻比的；駐點熱流略微有所增大，原因是飛行器為了減阻從而減小了機鼻半徑；熱載荷減少了5.22%說明單位時間內(nèi)傳入飛行器的氣動加熱量減小，這有利于減輕熱防護重量。

圖22為初始構(gòu)型與優(yōu)化構(gòu)型對比，其中粉色代表初始構(gòu)型、藍色代表優(yōu)化構(gòu)型。從圖22中可以看出，為了提高升阻特性，機身優(yōu)化的總體趨勢是采用較小前緣半徑的弧形鼻錐、上部傾斜的機身側(cè)面；機身最大橫截面積呈減小的趨勢并使機身橫截面積分布沿飛行器流向光滑過渡；低熱載荷的要求使飛行器下表面呈現(xiàn)扁平、鈍頭的趨勢。機翼優(yōu)化的總體趨勢是增大機翼面積以有效增大升力系數(shù)，同時為了減小機翼前緣氣動加熱，機翼前緣后掠角增加。

圖21 升阻比優(yōu)化收斂歷程Fig.21 History of lift-to-drag ratio optimization convergence

表1 優(yōu)化前后性能參數(shù)對比Table 1 Comparison of performance values before and after optimization

圖22 初始構(gòu)型與優(yōu)化構(gòu)型對比Fig.22 Comparison of initial configuration and optimal configuration

4.3 壓心后移對穩(wěn)定性的影響

從優(yōu)化構(gòu)型與初始構(gòu)型的對比可以看出，優(yōu)化構(gòu)型的機翼位置后移，其結(jié)果將會導(dǎo)致全機壓心后移。由于本文是基于一種實際的外形進行優(yōu)化，因此需對壓心變化造成的靜、動穩(wěn)定性進行分析，才能判定優(yōu)化結(jié)果的實際有效性。

圖23為優(yōu)化構(gòu)型與初始構(gòu)型在高馬赫數(shù)下的俯仰力矩系數(shù)(Cmy)對比。從圖中可以看出，在實際再入飛行所需的大迎角下，優(yōu)化構(gòu)型的俯仰力矩系數(shù)絕對值小于初始構(gòu)型，也就是，優(yōu)化構(gòu)型在該階段所需的配平舵偏角更小。因此，優(yōu)化構(gòu)型的配平升阻比更高。

圖24為優(yōu)化構(gòu)型與初始構(gòu)型在高馬赫數(shù)、大迎角下的靜穩(wěn)定性對比，其中黑色線條代表俯仰靜導(dǎo)數(shù)(Cmα)，藍色線條代表壓心與質(zhì)心的相對位置(xp)。從圖中可以看出，根據(jù)Cmα<0、xp>0的俯仰靜穩(wěn)定性判據(jù)，初始構(gòu)型在實際再入迎角附近具有靜穩(wěn)定性，這與文獻[31]的結(jié)果一致；而優(yōu)化構(gòu)型的靜穩(wěn)定姿態(tài)出現(xiàn)在更早的迎角，這說明在不考慮實際配平舵偏需求下，優(yōu)化構(gòu)型機翼位置后移帶來的靜穩(wěn)定性變化能夠更好滿足再入飛行要求。

圖23 俯仰力矩系數(shù)對比Fig.23 Comparison of pitch moment coefficients

圖24 靜穩(wěn)定性對比Fig.24 Comparison of static stability

對于動態(tài)穩(wěn)定性分析，采用RAPT氣動工具對飛行器動導(dǎo)數(shù)進行快速預(yù)測[32]，原理為首先對飛行器施加簡諧運動，然后基于當(dāng)?shù)亓骰钊碚?，快速預(yù)測周期非定常氣動力，再對非定常氣動力提取、辨識飛行器動態(tài)穩(wěn)定性導(dǎo)數(shù)。為求解俯仰組合動導(dǎo)數(shù)，對飛行器施加低頻小幅周期強迫振蕩運動。在定軸轉(zhuǎn)動情況下，強迫運動方程為

α=α0+αmsin(ωt+φ)

(10)

式中:α0和αm分別為初始迎角和迎角振幅;ω為強迫振動減縮頻率;φ為初始相位角。

本文對飛行器施加強迫振動的振幅取1°，縮減頻率取0.1。圖25為類X-37B飛行器優(yōu)化構(gòu)型與初始構(gòu)型的俯仰組合動導(dǎo)數(shù)對比，并以文獻[33] 的CFD結(jié)果作為參考。從圖中可以看出，優(yōu)化前后飛行器的俯仰組合動導(dǎo)數(shù)變化很小。由于有翼再入飛行器的姿態(tài)基本不會劇烈變化(不會產(chǎn)生高速旋轉(zhuǎn)之類的運動)，根據(jù)實際工程經(jīng)驗，這些微小的動導(dǎo)數(shù)改變量基本不會對飛行器的軌道設(shè)計、控制系統(tǒng)設(shè)計以及飛行品質(zhì)產(chǎn)生影響。

圖25 動導(dǎo)數(shù)對比Fig.25 Comparison of dynamic derivatives

4.4 優(yōu)化構(gòu)型TPS輕量化設(shè)計優(yōu)化

在熱防護方案設(shè)計時，初始設(shè)定飛行器表面將分為高溫區(qū)、中溫區(qū)與低溫區(qū)，對不同分區(qū)布置不同材料：高溫區(qū)外層布置單片增韌抗氧化復(fù)合材料(TUFROC)[34]，內(nèi)層布置先進柔性隔熱氈(AFRSI)；中溫區(qū)外層布置高溫重復(fù)使用陶瓷防熱瓦(LI-2000)，內(nèi)層布置芳綸柔性防熱氈(Nomex FRSI)；低溫區(qū)外層布置低溫重復(fù)使用陶瓷防熱瓦(LI-900)，內(nèi)層布置芳綸柔性防熱氈(Nomex FRSI)。材料屬性見文獻[35]。

圖26為TPS質(zhì)量優(yōu)化收斂歷程。從圖中可以看出，經(jīng)過16次迭代優(yōu)化，最終得到TPS理想質(zhì)量約為427 kg。

表2為典型再入飛行器構(gòu)型(包括返回艙/探測器構(gòu)型、升力體構(gòu)型、翼身組合體構(gòu)型)的熱防護重量占比統(tǒng)計數(shù)據(jù)[36-39]，其中對于返回艙/探測器構(gòu)型，由于統(tǒng)計包括了進入其他行星大氣層的返回艙/探測器，因此熱防護重量占比區(qū)間較大；典型再入地球大氣層的Apollo返回艙的熱防護重量占比約為18%。本文對類X-37B飛行器熱防護系統(tǒng)設(shè)計優(yōu)化得到的重量占比約為8.7%，低于同類飛行器的熱防護重量占比。

圖26 TPS質(zhì)量優(yōu)化收斂歷程Fig.26 TPS mass convergence history

表2 典型再入飛行器熱防護重量占比Table 2 TPS mass fractions of typical re-entry vehicles

圖27為飛行器表面根據(jù)最嚴(yán)酷氣動加熱狀態(tài)進行的熱防護分區(qū)。從圖中可以看出，飛行器頭部、機身下迎風(fēng)面前半部、機翼前緣等位置都被分為高溫區(qū)；背風(fēng)面等位置被分為低溫區(qū)，這與飛行器大攻角再入氣動熱環(huán)境相符。

表3為TPS質(zhì)量分布統(tǒng)計。從表中可以看出，低溫區(qū)熱防護材料占比最大，原因是低溫區(qū)占飛行器表面積最大。

圖27 飛行器表面TPS分區(qū)Fig.27 Vehicle surface TPS partition

表3 TPS質(zhì)量統(tǒng)計Table 3 TPS mass statistics

圖28為TPS厚度分布。從圖中可以看出：① 在氣動加熱嚴(yán)酷的位置TPS厚度較厚，而在背風(fēng)面、機翼遮擋處，所需的TPS厚度較?。虎?熱防護材料厚度分布具有一定連續(xù)性，這既符合飛行器表面熱載荷分布規(guī)律，又能滿足整體式熱防護系統(tǒng)的設(shè)計要求[40]。

圖28 TPS厚度分布Fig.28 TPS thickness distribution

圖29為飛行器在軌跡最終時刻的內(nèi)表面(熱防護底層與艙內(nèi)結(jié)構(gòu)接觸的表面)溫度分布云圖。從圖中可以看出，飛行器內(nèi)表面溫度均在500 K附近(優(yōu)化初始設(shè)定內(nèi)艙結(jié)構(gòu)的許用溫度為500 K)，說明熱防護厚度優(yōu)化起到了良好的效果。

為了判斷熱防護材料演軌跡的溫度響應(yīng)是否超過材料本身的最高許用溫度，對飛行器典型部位的溫度響應(yīng)歷程進行監(jiān)測。選取的3個監(jiān)測點分別為飛行器駐點、下表面中心線上x=25%L處及x=50%L處。圖30為2個監(jiān)測點的外表面、中間層及內(nèi)表面溫度響應(yīng)歷程。從圖中可以看出：① 監(jiān)測點的所有材料層均為超過材料的最高許用溫度，例如駐點外表面在軌跡上的最高溫度為1 885 K，沒有超過駐點表層熱防護材料RCC的許用溫度2 100 K；② 監(jiān) 測點的內(nèi)表面溫度在軌跡重點基本都收斂到設(shè)定的500 K附近，這說明本文對熱防護材料厚度的優(yōu)化，既能滿足艙內(nèi)許用溫度的要求，又能使熱防護重量最輕；③ 內(nèi) 表面降溫比外表面有一定的滯后性，這是由于當(dāng)飛行器過了軌跡上最嚴(yán)酷的氣動加熱段后，表面熱流開始減小，外表面的溫度開始下降，但此時飛行器外表面溫度任然比內(nèi)表面溫度高，溫度繼續(xù)向內(nèi)傳遞造成的。

圖29 TPS最終時刻內(nèi)表面溫度分布Fig.29 Inner surface temperature distribution of TPS at the final moment

圖30 典型位置溫度響應(yīng)歷程Fig.30 Temperature responses history of typical points

5 結(jié) 論

對有翼再入飛行器氣動外形進行了集成設(shè)計優(yōu)化，相關(guān)結(jié)果可得到如下結(jié)論：

1) 再入飛行過程對升阻比、升力系數(shù)、彈道系數(shù)等氣動參數(shù)的要求，直接決定飛行器氣動外形，這些影響因素可以總結(jié)為有翼再入飛行器的氣動外形設(shè)計規(guī)律與準(zhǔn)則。

2) 基于氣動外形設(shè)計規(guī)律與準(zhǔn)則，集成幾何參數(shù)化建模、氣動力快速預(yù)測、氣動熱加熱快速預(yù)測等學(xué)科模塊，采用多學(xué)科設(shè)計優(yōu)化的方法，能夠得到較優(yōu)氣動性能的飛行器外形。

3) 對類X-37B飛行器熱防護輕量化設(shè)計優(yōu)化得到了優(yōu)于同類飛行器熱防護重量占比統(tǒng)計數(shù)據(jù)的結(jié)果，且能夠滿足再入飛行使用要求，表明了本文熱防護輕量化設(shè)計優(yōu)化方法的有效性。