黃 飛，程曉麗，俞繼軍，姜貴慶

(中國(guó)航天空氣動(dòng)力技術(shù)研究院，北京 100074)

0 引言

新型高超聲速飛行器要經(jīng)歷長(zhǎng)時(shí)間的氣動(dòng)加熱環(huán)境，在材料非燒蝕的要求下，對(duì)飛行器的隔熱性能模擬提出了更高的要求。目前，隨著多相復(fù)合材料傳熱機(jī)制研究的不斷深入，材料的隔熱性能也相應(yīng)的得到大幅度提高。多相隔熱材料是利用空氣分子在低壓下碰撞頻率的降低來(lái)實(shí)現(xiàn)材料整體的良好隔熱性能。如當(dāng)空氣足夠稀薄時(shí)，可以近似認(rèn)為是不存在傳導(dǎo)導(dǎo)熱的，因此在材料中引入多空隙結(jié)構(gòu)，并降低組成纖維的尺寸大小，可以有效地提高材料的隔熱性能?，F(xiàn)有研究結(jié)果表明具有多孔結(jié)構(gòu)、密度較低、由纖維或顆粒組成的多相材料在飛行器的隔熱設(shè)計(jì)中具有良好的應(yīng)用前景，可以有效發(fā)揮阻斷或者延緩熱量傳遞的功能。

多孔隔熱材料存在明顯的微納米高孔隙結(jié)構(gòu)特征(圖1)，其熱導(dǎo)率對(duì)壓力的響應(yīng)非常敏感［1］，因此了解材料內(nèi)部壓力的響應(yīng)過(guò)程和響應(yīng)時(shí)間，對(duì)確定多孔隔熱材料的熱導(dǎo)率就顯得尤為重要，從現(xiàn)有研究手段來(lái)看多以實(shí)驗(yàn)為主，在某一給定壓力環(huán)境下測(cè)出材料的等效熱導(dǎo)率，由于壓力響應(yīng)時(shí)間的差異，很難提取響應(yīng)過(guò)程中其隨環(huán)境壓力變化的經(jīng)驗(yàn)關(guān)系式，國(guó)內(nèi)外許多學(xué)者開(kāi)展了相關(guān)方面的研究工作［2－5］，然而將DSMC方法應(yīng)用于材料內(nèi)部壓力響應(yīng)方面的研究工作還不多見(jiàn)。

基于此，文中在文獻(xiàn)［1］的研究基礎(chǔ)上進(jìn)一步對(duì)孔隙率、纖維尺寸、較大的宏觀尺寸等影響因素進(jìn)行了深入分析，研究了纖維尺寸 1.25μm、2.5μm、10μm，孔隙率85%、90%、95%，內(nèi)外壓比，以及模型的宏觀尺寸對(duì)內(nèi)部壓力響應(yīng)時(shí)間的影響，為材料的隔熱機(jī)理研究提供理論基礎(chǔ)。

圖1 多孔隔熱材料的微觀結(jié)構(gòu)Fig.1 Micromechanism of multiphase thermal insulation material

1 計(jì)算方法－DSMC方法

1.1 計(jì)算方法介紹:

DSMC［6］方法的基本思想是:用有限個(gè)仿真分子代替真實(shí)氣體分子，并在計(jì)算機(jī)中存儲(chǔ)仿真分子的位置坐標(biāo)、速度分量以及內(nèi)能，其值隨仿真分子的運(yùn)動(dòng)、與邊界的作用以及仿真分子之間的碰撞而改變，最后通過(guò)統(tǒng)計(jì)網(wǎng)格內(nèi)仿真分子的運(yùn)動(dòng)狀態(tài)實(shí)現(xiàn)對(duì)真實(shí)氣體流動(dòng)問(wèn)題的模擬。計(jì)算時(shí)除考慮平動(dòng)能外，還考慮了內(nèi)部能量。作為近似，分子的內(nèi)能只考慮轉(zhuǎn)動(dòng)能，未考慮振動(dòng)能。分子間相互作用模型采用VHS模型，物面邊界條件采用完全漫反射條件，平動(dòng)能與內(nèi)能的能量交換采用L-B模型［7－8］，按照Bird的能量按自由度分配原則采用取舍法進(jìn)行抽樣分配［8］。

本文根據(jù)材料特征對(duì)材料的?；W(wǎng)格進(jìn)行隨機(jī)地選擇標(biāo)注來(lái)區(qū)分纖維和孔隙，運(yùn)用改進(jìn)后的DSMC程序?qū)?biāo)記為實(shí)心的網(wǎng)格面進(jìn)行物面的漫反射處理，其它的處理方法與傳統(tǒng)DSMC方法一致。主要模擬步驟如下:

1)選取纖維的平均尺寸作為網(wǎng)格的尺寸，以孔隙率為標(biāo)準(zhǔn)對(duì)網(wǎng)格進(jìn)行隨機(jī)地選擇標(biāo)注，其中標(biāo)注網(wǎng)格1、0來(lái)區(qū)分纖維與孔隙，并用數(shù)組對(duì)每一網(wǎng)格的1、0標(biāo)注信息進(jìn)行記錄。

2)對(duì)于標(biāo)注為0的孔隙網(wǎng)格布置粒子，采用DSMC方法進(jìn)行分子的碰撞運(yùn)動(dòng)模擬，其中分子每次運(yùn)動(dòng)后對(duì)于標(biāo)注為纖維的網(wǎng)格4個(gè)面將視為固體物面，判斷分子是否與其碰撞，若分子與其作用將采用物面漫反射處理，反射后將以反射速度繼續(xù)運(yùn)動(dòng)剩余時(shí)間。

3)對(duì)材料中心點(diǎn)附近的幾個(gè)網(wǎng)格進(jìn)行采樣統(tǒng)計(jì)，可得出壓力隨時(shí)間的響應(yīng)規(guī)律。

1.2 主要關(guān)鍵技術(shù)點(diǎn):

1)分子與物面碰撞點(diǎn)的求法:

分子與物面的碰撞點(diǎn)就是分子與網(wǎng)格面的碰撞點(diǎn)，碰撞點(diǎn)的坐標(biāo)按照直線(即分子運(yùn)動(dòng)軌跡)與平面的交點(diǎn)求解(如圖2)。

分子運(yùn)動(dòng)前的坐標(biāo)為(XI，YI，ZI)，運(yùn)動(dòng)后的坐標(biāo)為(XI+DX，YI+DY，ZI+DZ)，設(shè)碰撞點(diǎn)的坐標(biāo)為(XC，YC，ZC)，P為物面網(wǎng)格的基本點(diǎn)，n為物面的法向矢量，n=RI×RJ，RI、RJ是由網(wǎng)格面的頂點(diǎn)所決定的基于P頂點(diǎn)的兩個(gè)矢量。設(shè)物面方程為:

分子運(yùn)動(dòng)軌跡的直線方程為:

由以上兩方程可得到:

從而可得到碰撞點(diǎn)的坐標(biāo)為:

本文為二維問(wèn)題，只需去掉Z方向的相關(guān)量即可。

2)分子與纖維碰撞的判斷方法:

由于分子頻繁地與周圍網(wǎng)格面作用，這對(duì)于分子所在網(wǎng)格的判斷帶來(lái)困難，尤其是當(dāng)分子越過(guò)兩個(gè)網(wǎng)格時(shí)的判斷更為復(fù)雜，在較小時(shí)間步長(zhǎng)下越過(guò)兩個(gè)以上網(wǎng)格的分子數(shù)極少，因此，作為近似，我們對(duì)越過(guò)兩個(gè)網(wǎng)格以上的分子進(jìn)行減速處理不會(huì)帶來(lái)很大誤差，而對(duì)于越過(guò)一個(gè)網(wǎng)格的分子進(jìn)行判斷相對(duì)容易。將分子運(yùn)動(dòng)前所在的網(wǎng)格記為II、JJ，運(yùn)動(dòng)后的網(wǎng)格記為ICN、JCN;此時(shí)我們可以將分子的判斷分為9種情況進(jìn)行區(qū)分判斷:

以上所有的情況中只要與網(wǎng)格面相交，求出交點(diǎn)后就需要判斷相鄰的網(wǎng)格是否標(biāo)注為纖維(即數(shù)組Ibiaozhu為1)，是纖維的進(jìn)行物面反射，反射后將以反射速度繼續(xù)運(yùn)動(dòng)剩余時(shí)間，不是纖維時(shí)分子將以原速度繼續(xù)運(yùn)動(dòng)剩余時(shí)間。

本文程序的算例驗(yàn)證見(jiàn)文獻(xiàn)［9］。

圖2 分子碰撞示意圖Fig.2 Molecular collision with wall

2 計(jì)算結(jié)果與分析

2.1 計(jì)算狀態(tài):

計(jì)算的網(wǎng)格示意圖如圖3，其中網(wǎng)格尺寸采用纖維的近似尺寸，分別選取了 1.25μm、2.5μm、10μm，其中孔隙率選取了85%、90%、95%，在網(wǎng)格中隨機(jī)地選取相應(yīng)比例的實(shí)心網(wǎng)格作為纖維，選取后的隔熱層簡(jiǎn)化模型如圖4，其中黑色的點(diǎn)代表場(chǎng)中的纖維，其余為孔隙結(jié)構(gòu)。宏觀尺寸選取了1mm ×1mm、2.5mm×2.5mm、5mm×5mm、10mm ×10mm、20mm ×20mm 等，材料內(nèi)部的壓力 Pin為1Pa、10Pa、100Pa等，外部環(huán)境的壓力 Pout為 0.1Pa、1Pa、10Pa、0.2Pa、2Pa、20Pa 等。

計(jì)算的來(lái)流溫度T∞為300K，壁面溫度Tw為300K，邊界AB和BD采用對(duì)稱邊界，AC和CD為自由來(lái)流邊界條件，即分子可以在四周邊界自由進(jìn)出。纖維物面邊界以壁溫Tw進(jìn)行漫反射處理。

2.2 壓力的時(shí)間響應(yīng)

圖5為隔熱層尺寸在1mm×1mm(網(wǎng)格100×100)時(shí)，隔熱層的內(nèi)壓為1Pa時(shí)的壓力初場(chǎng)，其中的實(shí)心點(diǎn)代表場(chǎng)中的纖維，場(chǎng)中其余位置布置空氣分子，其壓力值約為1Pa。圖 6中的(a)圖和(b)圖分別為 t=3000Δt、t=30000Δt時(shí)刻圖5中所示位置處的壓力響應(yīng)及速度矢量的局部放大圖。從t=3000Δt中可以看出，由于內(nèi)壓pin=1Pa大于外壓pout=0.1Pa，在這種壓力驅(qū)動(dòng)下，氣體不斷地向周圍環(huán)境擴(kuò)散，在隔熱層內(nèi)形成了明顯地壓力梯度，當(dāng)t=30000Δt時(shí)，內(nèi)壓達(dá)到與環(huán)境壓力基本相等，這種壓力梯度隨之消失，氣體表現(xiàn)為圖中隨機(jī)的擴(kuò)散運(yùn)動(dòng)。

孔隙率效應(yīng):

圖7給出了在纖維尺寸10μm、宏觀尺寸10mm×10mm下，不同孔隙率時(shí)的壓力響應(yīng)歷程，從圖7(a)中可以看出，在內(nèi)外壓力比pin/pout=10時(shí)，孔隙率降低，分子運(yùn)動(dòng)受到固體纖維的阻擋增強(qiáng)使得響應(yīng)時(shí)間有所增大，在此模擬條件下，孔隙率每降低5%，響應(yīng)時(shí)間約增加了0.002s，響應(yīng)時(shí)間的相對(duì)增量較大，從而可知，孔隙率對(duì)響應(yīng)時(shí)間的影響較大。圖7(b)為內(nèi)外壓比等于5時(shí)的響應(yīng)曲線，其中孔隙率對(duì)響應(yīng)時(shí)間的影響規(guī)律與圖7(a)的規(guī)律相似。

圖3 網(wǎng)格示意圖Fig.3 Grid in the computation region

圖4 隔熱層二維簡(jiǎn)化模型Fig.4 Simplified model of insulating material

2.2.2 纖維尺寸效應(yīng)

圖8給出了孔隙率90%時(shí)不同宏觀尺寸下，纖維尺寸對(duì)響應(yīng)時(shí)間的影響規(guī)律。從圖8(a)中可以發(fā)現(xiàn)，在相同的宏觀尺寸及孔隙率下，隨著纖維尺寸的減小，響應(yīng)時(shí)間增大。這主要是因?yàn)椋嗤紫堵氏?，纖維的尺寸越小，代表纖維的網(wǎng)格數(shù)越多，纖維面的總面積越大，氣體受纖維壁面的阻擋作用越強(qiáng)烈。圖中還可看出，在此條件下，響應(yīng)時(shí)間的相對(duì)增量較大，響應(yīng)時(shí)間對(duì)纖維的尺寸較為敏感。圖8(b)對(duì)應(yīng)的宏觀尺寸狀態(tài)下具有與圖8(a)中相似的變化規(guī)律，僅不同之處在于由于宏觀尺寸的增加導(dǎo)致每條曲線的響應(yīng)時(shí)間相對(duì)更長(zhǎng)。

2.2.3 內(nèi)外壓比效應(yīng)

圖9為纖維尺寸10μm，孔隙率90%，宏觀尺寸10mm時(shí)不同內(nèi)外壓下的響應(yīng)曲線。經(jīng)過(guò)厘米量級(jí)尺寸的計(jì)算，從圖9(a)中可看出，即使在較大尺寸下，當(dāng)尺寸不變時(shí)，內(nèi)外壓比相等，響應(yīng)時(shí)間也相等的結(jié)論仍然成立，這一結(jié)論與文獻(xiàn)［5］的結(jié)論相同。因此，為了節(jié)省計(jì)算時(shí)間，同樣尺寸下，對(duì)高壓狀態(tài)下的模擬可以采用低壓狀態(tài)來(lái)進(jìn)行近似模擬。圖9(b)中給出了不同壓比下的響應(yīng)曲線。從中可以發(fā)現(xiàn)，內(nèi)外壓比增加，壓力達(dá)到穩(wěn)態(tài)時(shí)需要的響應(yīng)時(shí)間將更長(zhǎng)，然而，由于壓比增加時(shí)，氣體擴(kuò)散的速度將有所增加，這在一定程度上又會(huì)使降低響應(yīng)時(shí)間有所降低，兩種相反的作用使得最終響應(yīng)時(shí)間受內(nèi)外壓比的影響較弱。宏觀尺寸效應(yīng):

圖10為纖維尺寸10μm，孔隙率90%時(shí)不同宏觀尺寸下壓力響應(yīng)的對(duì)數(shù)曲線。從圖中可以看出，當(dāng)材料長(zhǎng)度增大2倍即面積增大4倍時(shí)，響應(yīng)時(shí)間的對(duì)數(shù)值的增量近似相等。這一結(jié)果表明，即使宏觀尺寸進(jìn)一步增加，響應(yīng)時(shí)間的對(duì)數(shù)值近似呈現(xiàn)線性增加的結(jié)論仍然成立，這一結(jié)論與文獻(xiàn)［5］中小尺寸下的結(jié)論相同。

圖10 不同宏觀尺寸下的壓力響應(yīng)Fig.10 Pressure response at different material thickness

3 結(jié)論

基于多孔隔熱材料的特征，本文采用文獻(xiàn)［5］中建立的模擬手段，進(jìn)一步選取了 1.25μm、2.5μm、10μm 的纖維尺寸，85%、90%、95%的孔隙率，宏觀尺寸達(dá)到厘米量級(jí)的模擬環(huán)境，對(duì)影響壓力響應(yīng)的眾多因素進(jìn)行了系統(tǒng)性地分析。結(jié)果表明:

(1)孔隙率對(duì)響應(yīng)時(shí)間的影響較大，同一條件下，孔隙率減小，響應(yīng)時(shí)間增大。

(2)響應(yīng)時(shí)間對(duì)纖維尺寸較為敏感，同一條件下，纖維尺寸降低，響應(yīng)時(shí)間增大。

(3)響應(yīng)時(shí)間對(duì)內(nèi)外壓比的敏感性相對(duì)較弱，同一條件下，內(nèi)外壓比增加，響應(yīng)時(shí)間有所增加。

(4)即使在較大的宏觀尺寸下，當(dāng)材料長(zhǎng)度增大2倍即面積增大4倍時(shí)，響應(yīng)時(shí)間的對(duì)數(shù)值的增量也近似相等，即響應(yīng)時(shí)間的對(duì)數(shù)值近似呈現(xiàn)線性增加的結(jié)論仍然成立，這一結(jié)論與文獻(xiàn)［5］中小尺寸下的結(jié)論相同。

［1］黃飛，程曉麗，俞繼軍.多孔隙隔熱材料內(nèi)壓的時(shí)間響應(yīng)［J］.宇航學(xué)報(bào)，2010，31(1)∶233-238.(HUANG F，CHENG X L，YU J J.The study of pressure respond to time on multiaperture insulating material［J］.Journal of Astronautics，2010，31(1)∶233-238.)

［2］MARSCHALL J，MADDREN J.Internal radiation transport and effective thermal conductivity of fibrous ceramic insulation［R］.AIAA 2001-2822，2001.

［3］KAMRAN D.Thermal analysis and design of multi-layer insulation for re-entry aerodynamic［R］.AIAA 2001-2834，2001.

［4］KAMRAN D.Heat transfer in high-temperature fibrous insulation［R］.AIAA 2002-3332，2002.

［5］HANLON M A，MA H B.Evaporation heat transfer in sintered porous media［R］.AIAA 2002-3092，2002.

［6］BIRD G A.Molecular gas dynamics and direct simulation of gas flow clarendon［M］.Oxford，1994.

［7］BORGANOFF C，LARSEN P S.Statistical collision model for monte carlo simulation of polyatomic gas mixture［J］.Journal of Computational Physics，1975，18∶405-420.

［8］陳偉芳，吳其芬.高溫稀薄氣體熱化學(xué)非平衡流動(dòng)的DSMC方法［M］.長(zhǎng)沙:國(guó)防科技大學(xué)出版社，1999.(CHEN W F，WU Q F.DSMC method of nonequilibrium thermodynamics and chemical reaction in the high temperature rarefied gas flow［M］.Chang Sha:National University of Defense Technology Press，1999.)

［9］黃飛，程曉麗，沈清.高超聲速平板近空間氣動(dòng)特性的計(jì)算分析研究［J］.宇航學(xué)報(bào)，2009，30(3)∶900-907.(HUANG F，CHENG X L，SHEN Q.Numerical investigation of hypersonic aerodynamics of flat plate flying in near space［J］.Journal of Astronautics，2009，30(3)∶900-907.)