張 章, 吳 杰, 侯安平, 王立武, 竺梅芳

(1. 北京空間機(jī)電研究所 中國(guó)空間技術(shù)研究院航天器無損著陸技術(shù)核心專業(yè)實(shí)驗(yàn)室, 北京 100094；2. 北京航空航天大學(xué) 能源與動(dòng)力工程學(xué)院 航空發(fā)動(dòng)機(jī)氣動(dòng)熱力國(guó)家級(jí)重點(diǎn)實(shí)驗(yàn)室, 北京 100191)

0 引 言

隨著航天事業(yè)的發(fā)展，空天之間的物資來往越來越頻繁，在當(dāng)今的航天器返回過程中，普遍使用剛性飛行器回收太空物資，但這一技術(shù)正迅速接近其有效載荷極限。剛性飛行器一般被安裝在運(yùn)載火箭內(nèi)送入太空，因而結(jié)構(gòu)的最大直徑尺寸受限，從而影響其承載能力[1-2]?？臻g再入充氣結(jié)構(gòu)是美國(guó)宇航局研發(fā)的一種新興的再入與返回式航天器，利用充氣薄膜的可壓縮性，在發(fā)射過程中折疊以占有盡可能小的體積，返回時(shí)通過充氣使得具有防熱性能的柔性結(jié)構(gòu)形成氣動(dòng)外形，并通過這種充氣外形得到熱防護(hù)和氣動(dòng)減速功能[3-4]。具有質(zhì)量輕、發(fā)射體積小、結(jié)構(gòu)簡(jiǎn)單、靈活機(jī)動(dòng)、有效載荷比大等優(yōu)點(diǎn)，可以有效降低返回成本，并能夠應(yīng)用于多種復(fù)雜的再入返回任務(wù)中[5]?？臻g再入充氣結(jié)構(gòu)在返回過程中面臨嚴(yán)峻的飛行環(huán)境，在進(jìn)入大氣層后需承受超高聲速的氣流沖擊和巨大的氣動(dòng)熱載荷[6]。高超聲速氣動(dòng)力會(huì)導(dǎo)致阻力型面的變形，產(chǎn)生較大的結(jié)構(gòu)應(yīng)力，影響結(jié)構(gòu)強(qiáng)度和固有振動(dòng)特性[7-8]；氣動(dòng)熱載荷會(huì)使得充氣結(jié)構(gòu)駐點(diǎn)溫度上升、內(nèi)充壓氣體受熱劇烈膨脹，極易造成柔性充氣結(jié)構(gòu)破壞，影響再入返回任務(wù)的順利進(jìn)行[9-11]。因此，考慮氣動(dòng)力載荷下的非線性結(jié)構(gòu)動(dòng)力學(xué)研究成為空間再入充氣結(jié)構(gòu)設(shè)計(jì)階段不可忽視的環(huán)節(jié)[12]。

針對(duì)空間再入充氣結(jié)構(gòu)的結(jié)構(gòu)動(dòng)力學(xué)特性，國(guó)內(nèi)外學(xué)者已開展了相應(yīng)的研究工作。Lindell采用有限元方法對(duì)空間再入充氣結(jié)構(gòu)進(jìn)行了結(jié)構(gòu)分析，研究了在給定內(nèi)充壓情形下的應(yīng)力分布和模態(tài)振型[13]；Kinney在有限元數(shù)值建模中將柔性充氣薄膜簡(jiǎn)化為拉伸彈簧和彎曲彈簧兩種單元類型，開展了穩(wěn)態(tài)、瞬態(tài)結(jié)構(gòu)動(dòng)力學(xué)分析[14]；許家裕在獲得高超聲速定常氣動(dòng)力結(jié)果后進(jìn)行瞬態(tài)結(jié)構(gòu)動(dòng)力學(xué)計(jì)算，分析了IRVE系統(tǒng)的瞬態(tài)響應(yīng)情況[15]，但其僅將最大壓力加載到迎風(fēng)面頭部，未完整地考慮高超聲速流場(chǎng)真實(shí)氣動(dòng)力、氣動(dòng)加熱對(duì)空間再入充氣結(jié)構(gòu)的結(jié)構(gòu)動(dòng)力學(xué)特性的影響。

本文以美國(guó)第三代充氣式返回飛行器IRVE為研究對(duì)象，計(jì)算了一定初始再入條件下的二維平面彈道方程。以此為出發(fā)點(diǎn)，通過CFD數(shù)值模擬研究了再入過程中的阻力系數(shù)變化曲線，分析了半錐角、初始再入角的改變對(duì)飛行過載和結(jié)構(gòu)頂端熱流密度的影響關(guān)系。同時(shí)基于帶預(yù)緊力的薄膜振動(dòng)方程和非線性有限元理論，建立了空間再入充氣結(jié)構(gòu)的結(jié)構(gòu)動(dòng)力學(xué)模型，研究其在給定充氣壓力下的結(jié)構(gòu)靜力學(xué)特性和固有模態(tài)特征，并比較了不同充氣壓力和薄膜厚度對(duì)靜力學(xué)和模態(tài)特征的影響。最后本文運(yùn)用流固與熱固單向耦合的方法，建立相應(yīng)的耦合分析模型，以CFD數(shù)值模擬得到的真實(shí)氣動(dòng)力作為結(jié)構(gòu)有限元分析的輸入，研究了再入過程中各高度處的氣動(dòng)力、氣動(dòng)熱載荷作用下的結(jié)構(gòu)熱應(yīng)力及熱模態(tài)變化，極大程度上描述了真實(shí)氣動(dòng)力對(duì)空間再入充氣結(jié)構(gòu)材料非線性結(jié)構(gòu)動(dòng)力學(xué)行為的影響。

1 數(shù)值方法和計(jì)算模型

1.1 技術(shù)路線

本文針對(duì)空間再入充氣結(jié)構(gòu)的彈道方程推導(dǎo)、流場(chǎng)數(shù)值模擬及結(jié)構(gòu)有限元計(jì)算等具體研究?jī)?nèi)容，采用如圖1所示的技術(shù)路線。具體方法是：基于空間再入充氣結(jié)構(gòu)在返回過程中的運(yùn)動(dòng)方程自主編寫飛行彈道計(jì)算程序，求解飛行軌跡中高度、速度及過載等運(yùn)動(dòng)參數(shù)隨時(shí)間的變化關(guān)系；以彈道方程所得結(jié)果為邊界條件，結(jié)合國(guó)際標(biāo)準(zhǔn)大氣表，對(duì)空間再入充氣結(jié)構(gòu)所處流體域進(jìn)行CFD數(shù)值模擬，得到飛行過程中阻力系數(shù)、頂端熱流密度的變化，以及結(jié)構(gòu)表面氣動(dòng)壓力及溫度的分布特征；通過提取CFD計(jì)算得到的空間再入充氣結(jié)構(gòu)表面氣動(dòng)力及氣動(dòng)熱載荷分布，分別加載到有限元模型上進(jìn)行結(jié)構(gòu)動(dòng)力學(xué)分析，得到考慮真實(shí)流場(chǎng)作用下的柔性充氣結(jié)構(gòu)熱應(yīng)力分布與熱模態(tài)振動(dòng)特征。

圖1 計(jì)算流程圖Fig.1 Flow chart of analysis

1.2 薄膜靜力方程與振動(dòng)方程

根據(jù)薄膜靜力學(xué)方程計(jì)算求解預(yù)應(yīng)力和開展靜力學(xué)分析。由于再入充氣結(jié)構(gòu)需要內(nèi)充壓氣體產(chǎn)生的張力才能維持一定的型面和剛度，因此一切結(jié)構(gòu)特性計(jì)算的前提是具有一定壓力的內(nèi)充壓氣體。該結(jié)構(gòu)是一種具有預(yù)應(yīng)力的結(jié)構(gòu)，對(duì)這一類結(jié)構(gòu)進(jìn)行模態(tài)分析前需進(jìn)行靜力學(xué)分析。薄膜靜力學(xué)方程可以表示如下形式:

(1)

其中，h為薄膜厚度，p為充氣壓力，E為彈性模量，z為橫向位移，z是x、y的函數(shù)。通過薄膜自由振動(dòng)方程及系統(tǒng)的振動(dòng)微分方程，求解固有頻率。微分形式的薄膜自由振動(dòng)方程可以表示為如下形式:

(2)

式中，ms為單位面積上的質(zhì)量，p為充氣壓力，Tx和Ty為由變形引起的內(nèi)部張力。

依據(jù)一定的原則把分布質(zhì)量換算為集中的質(zhì)量，接著像處理集中質(zhì)量一樣建立質(zhì)量矩陣。其中換算的方法有很多，如果質(zhì)量的分布是均勻的，那么只需要把質(zhì)量平均分配到所有節(jié)點(diǎn)。至于質(zhì)量分布不均勻的情況，可以預(yù)先分配給各個(gè)節(jié)點(diǎn)所要負(fù)責(zé)的區(qū)域，接著將各個(gè)區(qū)域的質(zhì)量分配給這些節(jié)點(diǎn)，這樣建立起來的質(zhì)量矩陣具有除對(duì)角線元素外其他元素為0的特點(diǎn)。利用此方法，建立再入充氣結(jié)構(gòu)的集中質(zhì)量矩陣，根據(jù)多自由度系統(tǒng)自由振動(dòng)時(shí)候的頻率行列式

[KA]-ω2[MA]=0

(3)

求解系統(tǒng)各階固有頻率。其中KA為系統(tǒng)剛度矩陣，MA為系統(tǒng)質(zhì)量矩陣，而ω是所求的固有頻率。

1.3 CFD計(jì)算方法

根據(jù)國(guó)際標(biāo)準(zhǔn)大氣表以及空間再入充氣結(jié)構(gòu)的真實(shí)飛行彈道方程，對(duì)流場(chǎng)計(jì)算設(shè)置邊界條件，選取合適的湍流模型和求解方程。本文中,針對(duì)不同的進(jìn)口條件相應(yīng)地選擇亞聲速或超聲速邊界條件:對(duì)于亞聲速飛行，進(jìn)口給定總溫、速度，出口給定靜壓；對(duì)于超聲速飛行，進(jìn)口給定總溫、總壓和速度，出口不給定邊界條件。高超聲速氣動(dòng)力和氣動(dòng)熱計(jì)算采用基于有限體積法的雷諾平均N-S方程進(jìn)行求解，湍流模型采用SST模型，空間離散采用迎風(fēng)格式，時(shí)間離散采用二階歐拉后差格式來保證時(shí)間方向上的守恒。

1.4 結(jié)構(gòu)熱傳導(dǎo)及內(nèi)充壓氣體熱膨脹方程

對(duì)具有一定溫度邊界的結(jié)構(gòu)，加載溫度邊界條件至有限元節(jié)點(diǎn)上，根據(jù)材料物理特性確定熱膨脹系數(shù)和導(dǎo)熱系數(shù)，根據(jù)熱傳導(dǎo)方程計(jì)算穩(wěn)態(tài)溫度分布，并采用理想氣體狀態(tài)方程修正內(nèi)充壓氣體參數(shù)變化，在傳熱計(jì)算與內(nèi)充壓氣體參數(shù)修正的基礎(chǔ)上計(jì)算求解熱應(yīng)力。其中，熱傳導(dǎo)微分方程可以表示為如下形式：

(4)

式中，T為溫度，c為比熱容，λ為導(dǎo)熱系數(shù)。

理想氣體狀態(tài)方程可以表示為如下形式：

pM=ρRT

(5)

式中，p為充氣壓力，M為充壓氣體摩爾質(zhì)量，R為氣體常數(shù)8.314 J/(mol·K)，T是傳熱計(jì)算后得到的內(nèi)充壓氣體溫度。

熱應(yīng)力計(jì)算由熱變形方程和胡克定律求解得到:

(6)

式中，α為熱膨脹系數(shù)，E為彈性模量，σ為應(yīng)力，ε為應(yīng)變。

1.5 計(jì)算模型

根據(jù)空間再入充氣結(jié)構(gòu)外形參數(shù)，采用建模軟件建立幾何模型，并涵蓋環(huán)形氣囊、隔層、蒙皮等真實(shí)幾何特征。本文中所參考的結(jié)構(gòu)實(shí)物圖為美國(guó)IRVE飛行器，如圖2所示。

圖2 IRVE結(jié)構(gòu)示意圖Fig.2 Schematic diagram of IRVE system

圖3展示了空間再入充氣結(jié)構(gòu)的結(jié)構(gòu)有限元模型，采用Shell 181四節(jié)點(diǎn)面單元，網(wǎng)格數(shù)約為10萬。由于飛行器具有軸對(duì)稱的外形，所以取飛行器的一半進(jìn)行計(jì)算，并把對(duì)稱面的邊界條件設(shè)置為對(duì)稱邊界。選取計(jì)算域?yàn)橹睆綖?0.0 m、高為42.0 m的圓柱形，圓柱里又有一個(gè)直徑為3.75 m、高為10.5 m的圓柱形區(qū)域，該區(qū)域的建立是為了進(jìn)行近壁加密。圖4是取一半后的計(jì)算流場(chǎng)，采用混合網(wǎng)格進(jìn)行劃分，總網(wǎng)格數(shù)約為800萬。

圖3 有限元模型Fig.3 Finite element model

圖4 CFD模型Fig.4 CFD model

2 二維彈道方程計(jì)算

2.1 運(yùn)動(dòng)學(xué)方程推導(dǎo)

通過飛行動(dòng)力學(xué)計(jì)算可得到彈道方程、熱流密度、阻力系數(shù)以及溫度分布等參數(shù)的工程算法結(jié)果，其計(jì)算流程如圖5所示。具體方法是:建立空間再入充氣結(jié)構(gòu)在再入返回過程中的動(dòng)力學(xué)方程，開發(fā)飛行彈道計(jì)算程序，結(jié)合國(guó)際標(biāo)準(zhǔn)大氣表得到不同再入初始條件下的飛行軌跡及各高度狀態(tài)參數(shù)；并根據(jù)氣體分子疏密程度，選擇不同的駐點(diǎn)熱流計(jì)算公式，開發(fā)熱流密度工程算法；采用熱流密度工程算法對(duì)空間再入充氣結(jié)構(gòu)在自由分子流區(qū)、過渡流區(qū)、連續(xù)流區(qū)進(jìn)行外熱流數(shù)值估算，求出各流態(tài)剛性球頭駐點(diǎn)的熱流值，結(jié)合LEES熱流分布公式即可得到非駐點(diǎn)區(qū)域的熱流分布。

圖5 軌跡計(jì)算流程圖Fig.5 Trajectory calculation flow chart

對(duì)于彈道式的再入，因飛行時(shí)間短，可忽略地球自轉(zhuǎn)的影響，其運(yùn)動(dòng)更接近平面運(yùn)動(dòng)，因此彈道式再入的質(zhì)心運(yùn)動(dòng)描述，使用二維彈道方程更加簡(jiǎn)潔有效[16-17]。為了獲得簡(jiǎn)化后的二維彈道運(yùn)動(dòng)方程，本文忽略地球自轉(zhuǎn)對(duì)彈道的影響并假設(shè)飛行器以穩(wěn)定的姿態(tài)飛行，0°攻角，無側(cè)滑?；诖思僭O(shè)，飛行器的運(yùn)動(dòng)方向由再入初始條件確定。由再入點(diǎn)的位置矢量和速度矢量就可以唯一確定整個(gè)再入彈道所在的平面，在這個(gè)平面內(nèi)可以建立一個(gè)再入坐標(biāo)系，如圖6所示。

圖6 再入坐標(biāo)系Fig.6 Re-entry coordinate

圖6中oe-xeye就是由再入初始條件確定的再入坐標(biāo)系，原點(diǎn)oe是再入點(diǎn)在地面的投影。Θ是速度與當(dāng)?shù)厮骄€的夾角，這里稱為飛行路徑角。βe是r和ye之間的夾角，稱為射程角。經(jīng)過推導(dǎo)，可以得到二維運(yùn)動(dòng)軌跡方程：

(7)

2.2 半錐角和初始再入角對(duì)過載的影響

過載直接反映了飛行器除重力外所受合力的大小，由于結(jié)構(gòu)的強(qiáng)度存在限度，因此飛行過程中對(duì)過載峰值也有著明確的要求。對(duì)于空間再入充氣結(jié)構(gòu)，一般要求過載峰值不能超過10倍重力加速度[18]。本文探討了半錐角和初始再入角改變對(duì)空間再入充氣結(jié)構(gòu)飛行過載的影響。圖7是再入過程中不同半錐角的過載曲線，可以看出半錐角對(duì)再入段過載影響不大。其原因是，半錐角雖然會(huì)影響飛行過程中的阻力系數(shù)，但阻力的大小還和速度的平方成正比。因?yàn)楫?dāng)阻力系數(shù)升高的同時(shí)，加速度增大導(dǎo)致速度降低更多，從而阻力并沒有產(chǎn)生較大的差異，因此在圖中表現(xiàn)為不同半錐角的過載曲線差異較小。

圖7 不同半錐角下的過載曲線Fig.7 Drag acceleration curve at different half taper angles

圖8是再入過程中改變初始再入角時(shí)的過載曲線。與改變半錐角不同，初始再入角對(duì)過載曲線有顯著影響。當(dāng)初始再入角增大時(shí)，最大過載值迅速增加：初始再入角從-1°增加到-6°時(shí)，過載峰值從-3.03g迅速增加到-16.24g，此時(shí)結(jié)構(gòu)面臨嚴(yán)重的不穩(wěn)定性。

圖8 不同初始再入角下的過載曲線Fig.8 Drag acceleration curve at different initial re-entry angles

2.3 半錐角和初始再入角的最優(yōu)化設(shè)計(jì)

工程中希望在滿足一定的約束條件下，過載峰值能盡可能小，從而為空間再入充氣結(jié)構(gòu)的保形設(shè)計(jì)和結(jié)構(gòu)安全性設(shè)計(jì)提供有價(jià)值的參考，并為其他性能參數(shù)的設(shè)計(jì)選取提供較大的裕度[19]。為探尋能使過載峰值達(dá)到最小的再入角和半錐角選擇，本文采用遺傳算法對(duì)其進(jìn)行優(yōu)化。遺傳算法是進(jìn)化計(jì)算的一個(gè)分支，是一種模擬自然界生物進(jìn)化過程的隨機(jī)搜索算法，借鑒生物界自然選擇原理和自然遺傳機(jī)制而形成的一種迭代式自適應(yīng)概率性全局優(yōu)化搜索算法[20]。

本文問題可表述為：在再入過程的總時(shí)間內(nèi)，求解最優(yōu)控制變量——半錐角和初始再入角，使給定目標(biāo)過載峰值在滿足狀態(tài)方程、邊界條件、約束條件的前提下達(dá)到最小。流程圖如圖9所示。具體思路是通過一定的編碼方式，將半錐角和初始再入角編碼成一定的二進(jìn)制字符串，代表遺傳過程中每一代的染色體組成。由此得到個(gè)體在問題域中的適應(yīng)度值，結(jié)合從遺傳學(xué)中借鑒的再造方法進(jìn)行個(gè)體選擇，產(chǎn)生一個(gè)新的近似解。這個(gè)過程使種群中個(gè)體不斷進(jìn)化，得到新的個(gè)體比原個(gè)體具有更小的過載峰值，從而滿足研究要求。

圖9 遺傳算法示意圖Fig.9 Genetic algorithm schematic diagram

表1展示了優(yōu)化后的角度和性能參數(shù)，與沒有優(yōu)化相比，返回過程中的過載峰值降低到-1.96g，此時(shí)半錐角的取值為58.65°，初始再入角的取值為1.35°，此時(shí)再入總時(shí)間相比優(yōu)化前有所延長(zhǎng)。

表1 優(yōu)化前后參數(shù)表Table 1 Parameters table before and after optimization

3 高超聲速流場(chǎng)計(jì)算

3.1 氣動(dòng)壓力及表面溫度分布

將飛行彈道計(jì)算結(jié)果作為CFD計(jì)算的邊界條件，對(duì)空間再入充氣結(jié)構(gòu)不同高度下的流場(chǎng)進(jìn)行數(shù)值模擬。表2給出了再入過程中的馬赫數(shù)變化情況，其數(shù)值從23.6逐漸下降到0.7。

表2 馬赫數(shù)分布表Table 2 Mach number distribution table

圖10是IRVE系統(tǒng)表面最大氣動(dòng)壓力與最高溫度隨飛行高度的變化曲線。由圖可知，隨高度的增加，最大氣動(dòng)壓力有下降趨勢(shì)，而表面最大溫度有上升趨勢(shì)。

圖10 最大氣動(dòng)壓力和表面溫度隨高度變化曲線Fig.10 Maximum aerodynamic pressure and surface temperature curves vs. height

圖11和圖12展示了60 km高度處空間再入充氣結(jié)構(gòu)迎風(fēng)面和背風(fēng)面的氣動(dòng)壓力和溫度云圖。由圖可知，最大壓力位于迎風(fēng)面中心，大小為0.591 kPa。表面最大溫度同樣位于迎風(fēng)面中心，為640.5 K。

圖11 結(jié)構(gòu)迎風(fēng)面壓力分布Fig.11 Pressure distributions in windward side

圖12 結(jié)構(gòu)迎風(fēng)面溫度分布Fig.12 Temperature distributions in windward side

3.2 阻力系數(shù)求解與反饋修正

為精確計(jì)算阻力系數(shù)的大小，本文提取CFD模擬后的表面壓力分布，分別在迎風(fēng)面和背風(fēng)面處進(jìn)行壓力積分，得到空間再入充氣結(jié)構(gòu)外表面所受氣動(dòng)阻力的大小，結(jié)合結(jié)構(gòu)尺寸及速度大小，計(jì)算阻力系數(shù)的數(shù)值，并反饋至彈道方程中進(jìn)行修正。圖13是不同半錐角構(gòu)型時(shí)，空間再入充氣結(jié)構(gòu)的阻力系數(shù)隨馬赫數(shù)變化曲線。由圖可知，在跨聲速區(qū)，隨著馬赫數(shù)的增加，各構(gòu)型的阻力系數(shù)先增大后減??；當(dāng)馬赫數(shù)繼續(xù)增加到高超聲速區(qū)時(shí)，50°半錐角構(gòu)型的阻力系數(shù)有明顯的下降，而其他構(gòu)型的阻力系數(shù)變化幅度較小；當(dāng)馬赫數(shù)大于15時(shí)，各構(gòu)型的阻力系數(shù)曲線差異逐漸明顯；最后當(dāng)馬赫數(shù)超過20時(shí)，阻力系數(shù)趨于不變。

圖13 阻力系數(shù)曲線Fig.13 Drag coefficient curve

3.3 熱流密度計(jì)算

通過對(duì)CFD計(jì)算結(jié)果進(jìn)行后處理，并改變半錐角和初始再入角，獲得飛行過程中空間再入充氣結(jié)構(gòu)頂端熱流密度的變化曲線。圖14是再入過程中不同半錐角對(duì)駐點(diǎn)熱流密度的影響曲線。結(jié)果表明，當(dāng)初始再入角不變時(shí)，隨著半錐角的增大，駐點(diǎn)的熱流密度呈減小趨勢(shì)，但駐點(diǎn)最大熱流所在高度不斷上升。

圖15是改變初始再入角時(shí)駐點(diǎn)熱流密度的變化曲線。結(jié)果表明，當(dāng)固定半錐角時(shí)，隨著初始再入角的增大，再入過程中駐點(diǎn)熱流密度隨之增大，且駐點(diǎn)最大熱流所在高度呈下降趨勢(shì)。

圖14 熱流密度隨半錐角變化曲線Fig.14 Heat flux curve at different half taper angles

圖15 熱流密度隨初始再入角變化曲線Fig.15 Heat flux curve at different initial angles

4 考慮內(nèi)充壓作用的靜力學(xué)和模態(tài)分析

4.1 靜力學(xué)分析

根據(jù)之前建立的有限元結(jié)構(gòu)，在其內(nèi)部氣囊的有限元節(jié)點(diǎn)上加載一定壓力的載荷來等效代替內(nèi)充壓氣體。具體方法是打開預(yù)應(yīng)力選項(xiàng)，輸入材料參數(shù)，進(jìn)行靜力計(jì)算，得到一定內(nèi)充壓作用下系統(tǒng)應(yīng)力分布。圖16展示了充氣壓力20 kPa、薄膜厚度0.5 mm條件下的結(jié)構(gòu)靜應(yīng)力分布，從云圖可以看出結(jié)構(gòu)最大應(yīng)力為15.464 MPa，位于最底層氣囊和隔層的連接處；最小應(yīng)力為61.084 kPa，位于蒙皮上。

圖16 結(jié)構(gòu)靜應(yīng)力分布Fig.16 Static stress distribution

4.2 模態(tài)分析

根據(jù)多自由度系統(tǒng)自由振動(dòng)時(shí)的頻率行列式求解系統(tǒng)各階固有頻率，并對(duì)求解得到的初始解進(jìn)行模態(tài)擴(kuò)展，然后在后處理器中查看振型。根據(jù)此方法對(duì)充氣壓力20 kPa、薄膜厚度0.5 mm條件下的結(jié)構(gòu)進(jìn)行模態(tài)分析，得到固有頻率及模態(tài)振型如表3和圖17所示。由圖可知，前兩階固有頻率一致，為106.89 Hz，對(duì)應(yīng)的是扭轉(zhuǎn)振型，分別沿著兩個(gè)不同的方向；三階、四階為頻率不同的拉伸振型。

表3 前六階固有頻率Table 3 Natural frequency of first six orders

對(duì)內(nèi)充氣壓力10～30 kPa、薄膜厚度在0.2～1 mm內(nèi)的不同結(jié)構(gòu)進(jìn)行了計(jì)算，得到結(jié)構(gòu)最大靜應(yīng)力和一階固頻隨內(nèi)充氣壓力和薄膜厚度的變化規(guī)律分別如圖18和圖19所示。由圖可知，最大靜應(yīng)力和一階固頻都隨內(nèi)充氣壓力的增加呈線性遞增關(guān)系；另一方面，隨薄膜厚度的增加，最大靜應(yīng)力呈反比下降規(guī)律，而一階固頻以近似線性規(guī)律下降。

圖17 前四階固有振型Fig.17 Modal shapes of first four orders

圖18 最大靜應(yīng)力和一階固頻隨內(nèi)充壓變化曲線Fig.18 Maximum stress and first order frequency curves vs. internal pressure

圖19 最大靜應(yīng)力和一階固頻隨薄膜厚度變化曲線Fig.19 Maximum stress and first order frequency curves vs. film thickness

5 高超聲速流場(chǎng)作用下的靜力學(xué)和模態(tài)分析

5.1 考慮氣動(dòng)壓力下的靜力學(xué)和模態(tài)分析

在流固單向耦合研究中，提取氣動(dòng)力載荷中表面壓力分布信息，通過節(jié)點(diǎn)插值的方法加載至有限元結(jié)構(gòu)上進(jìn)行靜力學(xué)分析，得到其整體質(zhì)量矩陣與剛度矩陣，從而求出在一定內(nèi)充壓氣體作用下空間再入充氣結(jié)構(gòu)的整體靜應(yīng)力分布，再以靜力學(xué)分析得到的預(yù)緊力為基礎(chǔ)，對(duì)其進(jìn)行模態(tài)分析，得到結(jié)構(gòu)固有頻率及各階振型等模態(tài)振動(dòng)特性，其具體加載示意圖如圖20所示。

圖21展示了考慮氣動(dòng)壓力時(shí)，結(jié)構(gòu)最大靜應(yīng)力和一階固頻隨飛行高度的變化曲線?？梢缘贸觯?dāng)飛行高度低于40 km時(shí)，最大靜應(yīng)力隨高度增加而降低，一階固頻隨高度增加而增加；當(dāng)飛行高度超過40 km時(shí)，隨高度增加，最大靜應(yīng)力變化不明顯，第一階頻率略有下降趨勢(shì)。

圖20 表面壓力加載示意圖Fig.20 Surface pressure loading diagram

圖21 考慮氣動(dòng)壓力時(shí)最大靜應(yīng)力和一階固頻隨高度變化曲線Fig.21 Maximum static stress and first order frequency curves vs. height considering aerodynamic pressure

5.2 考慮氣動(dòng)加熱下的靜力學(xué)和模態(tài)分析

5.2.1 結(jié)構(gòu)熱傳導(dǎo)

在熱固單向耦合研究中，基于CFD數(shù)值模擬得到的氣動(dòng)熱載荷，提取可展開式氣動(dòng)減速與熱防護(hù)結(jié)構(gòu)的表面溫度分布信息，通過節(jié)點(diǎn)插值的方法將表面溫度加載至有限元模型上，進(jìn)行熱傳導(dǎo)計(jì)算?？臻g再入充氣結(jié)構(gòu)在返回過程中會(huì)面臨巨大的氣動(dòng)熱載荷沖擊，因此其結(jié)構(gòu)表面通常配置熱防護(hù)系統(tǒng)(TPS)以阻擋熱流穿透。在熱防護(hù)結(jié)構(gòu)各功能層瞬態(tài)熱傳導(dǎo)分析中，需推導(dǎo)適用于多功能層的瞬態(tài)熱傳導(dǎo)方程，建立防熱層、絕熱層、氣密層的瞬態(tài)熱傳導(dǎo)模型，設(shè)定材料參數(shù)及溫度、熱流邊界條件，選用合適的熱傳導(dǎo)單元進(jìn)行瞬態(tài)熱傳導(dǎo)的有限元分析，進(jìn)而得到各功能層溫度分布云圖。圖22展示了返回過程中，TPS最外部防熱層和最內(nèi)部氣密層的平均溫度隨時(shí)間的變化，可以看出溫度的傳遞存在明顯的延時(shí)效應(yīng)。

圖22 TPS最外層及最內(nèi)層溫度變化Fig.22 Temperature variation in outer and inner layers of TPS

5.2.2 內(nèi)充壓氣體熱膨脹

通過提取最內(nèi)層氣密層節(jié)點(diǎn)溫度，進(jìn)行算術(shù)平均，近似模擬氣動(dòng)熱載荷作用下的內(nèi)充壓氣體溫度，并結(jié)合理想氣體狀態(tài)方程計(jì)算，得到溫度變化后內(nèi)充壓氣體的壓力變化情況，從而有效地模擬出內(nèi)充壓氣體受熱膨脹的物理現(xiàn)象。圖23展示了氣密層平均溫度及受熱膨脹后的新充氣壓力隨飛行高度的變化曲線，可見在45 km時(shí)，氣體受熱程度達(dá)到最大。

5.2.3 熱應(yīng)力及熱模態(tài)分析

由于各功能層溫度變化及內(nèi)充壓氣體的壓力變化都能引起結(jié)構(gòu)動(dòng)力學(xué)行為發(fā)生較大改變，因此將迎風(fēng)面和背風(fēng)面的溫度分布加載至結(jié)構(gòu)有限元模型中，并改變內(nèi)充壓氣體壓力數(shù)值，進(jìn)行結(jié)構(gòu)熱應(yīng)力計(jì)算，即可得到在氣動(dòng)熱載荷作用下，由于結(jié)構(gòu)溫度和內(nèi)充壓氣體壓力變化產(chǎn)生的結(jié)構(gòu)變形及熱應(yīng)力分布。圖24是氣動(dòng)熱載荷的加載示意圖。

圖23 熱膨脹下內(nèi)部氣體壓力變化Fig.23 Pressure of internal gas after thermal expansion

通過將不同高度氣動(dòng)熱載荷依次施加到結(jié)構(gòu)上，得到隨高度變化的結(jié)構(gòu)熱應(yīng)力及一階固頻曲線，如圖25所示。該曲線僅考慮了高超聲速流場(chǎng)中的氣動(dòng)熱載荷。與之前考慮氣動(dòng)力載荷相比，可以發(fā)現(xiàn)當(dāng)飛行高度小于40 km時(shí)，氣動(dòng)力對(duì)結(jié)構(gòu)靜應(yīng)力的影響更大；當(dāng)飛行高度大于40 km時(shí)，氣動(dòng)熱對(duì)結(jié)構(gòu)靜應(yīng)力的影響更大。原因是，隨高度的增加，氣動(dòng)壓力呈下降趨勢(shì)而表面溫度呈上升趨勢(shì)，因此在高度較低時(shí)氣動(dòng)壓力的作用更為顯著，而高度上升時(shí)氣動(dòng)熱載荷的作用明顯提升。

圖25 考慮氣動(dòng)加熱時(shí)最大靜應(yīng)力和一階固頻隨高度變化曲線Fig.25 Maximum static stress and first order frequency curves vs. height considering aeroheating

6 結(jié) 論

(1) 本文基于空間再入充氣結(jié)構(gòu)的二維彈道方程，分析了半錐角和初始再入角對(duì)飛行過載和結(jié)構(gòu)頂端熱流密度的影響規(guī)律，并通過遺傳算法對(duì)軌跡進(jìn)行了優(yōu)化，為該類再入飛行器的軌跡分析提供了一般性的思路。

(2) 建立空間再入充氣結(jié)構(gòu)的結(jié)構(gòu)動(dòng)力學(xué)模型，研究了不同充氣壓力和薄膜厚度對(duì)靜力學(xué)特性和模態(tài)特征的影響規(guī)律，得出最大應(yīng)力和固有頻率隨內(nèi)充壓的增加而增加，隨薄膜厚度的增加而減小的結(jié)論。

(3) 通過流固及熱固單向耦合的方法，研究了氣動(dòng)壓力和氣動(dòng)加熱對(duì)結(jié)構(gòu)靜應(yīng)力及模態(tài)的影響。研究成果有望準(zhǔn)確描述空間再入充氣結(jié)構(gòu)的材料非線性結(jié)構(gòu)動(dòng)力學(xué)行為，為空間充氣式再入飛行器的保形設(shè)計(jì)和結(jié)構(gòu)安全性設(shè)計(jì)提供有價(jià)值的參考。