錢煒祺,周宇,何開鋒,邵元培
(1.空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,四川 綿陽 621000;2.中國(guó)空氣動(dòng)力研究與發(fā)展中心,四川 綿陽 621000)
考慮燒蝕情況下的表面熱流辨識(shí)
錢煒祺1,2,周宇1,2,何開鋒1,2,邵元培1,2
(1.空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,四川 綿陽 621000;2.中國(guó)空氣動(dòng)力研究與發(fā)展中心,四川 綿陽 621000)
針對(duì)燒蝕傳熱問題,在熱解面模型的基礎(chǔ)上通過伴隨方程推導(dǎo)建立了基于測(cè)點(diǎn)溫度辨識(shí)表面熱流的方法,并進(jìn)行了算例考核。結(jié)果表明:熱流辨識(shí)結(jié)果與真值符合較好,辨識(shí)結(jié)果與真值之間的偏差隨測(cè)量誤差的增加而增加;燒蝕后退量測(cè)量結(jié)果的誤差對(duì)辨識(shí)結(jié)果有較為顯著的影響。然后,將該辨識(shí)方法用于鈍頭型碳酚醛材料Narmco4028試件在陶瓷加熱風(fēng)洞中的燒蝕試驗(yàn)結(jié)果分析,結(jié)果表明辨識(shí)出的表面熱流與加熱功率基本符合,辨識(shí)方法是有效的,在工程實(shí)際中有較好的應(yīng)用前景。
燒蝕;表面熱流辨識(shí);熱解面模型
飛行器在再入大氣層過程中,高速氣流流過飛行器表面,由于氣體粘性的阻滯作用,會(huì)帶來氣動(dòng)加熱問題,使飛行器表面溫度明顯升高,以美國(guó)阿波羅飛船指令艙和前蘇聯(lián)聯(lián)盟號(hào)飛船返回艙為例,再入時(shí)的最高表面溫度分別達(dá)到了2950 K和1900K,這會(huì)給飛行安全帶來影響。因此,再入飛行器需進(jìn)行防熱系統(tǒng)設(shè)計(jì)。燒蝕防熱是目前主要采用的防熱系統(tǒng)設(shè)計(jì)方法,即通過防熱材料的熱解、燒蝕后退來吸收并帶走氣動(dòng)加熱的熱量,從而達(dá)到熱防護(hù)的目的[1-2]。
表面熱流的辨識(shí),或稱為表面熱流的反演,是指在飛行器表面溫度較高或表面有燒蝕后退的情況下,不適宜在表面直接安裝傳感器來進(jìn)行熱流和溫度測(cè)量,只能在飛行器的防熱層內(nèi)某些點(diǎn)上安裝溫度傳感器進(jìn)行溫度測(cè)量,然后根據(jù)這些點(diǎn)上的溫度測(cè)量歷程來反算飛行器的表面熱流。這類問題是一類典型的熱傳導(dǎo)逆問題IHCP(Inverse Heat Conduction Problem),是傳熱學(xué)研究的一個(gè)重要領(lǐng)域。目前,針對(duì)固定幾何域的熱傳導(dǎo)逆問題研究較多,但由于防熱層的燒蝕后退對(duì)應(yīng)著防熱層幾何外形的變化,因而此時(shí)防熱層內(nèi)的熱傳導(dǎo)逆問題是一個(gè)變幾何域的熱傳導(dǎo)逆問題,研究工作相對(duì)較少,尤其是在防熱層燒蝕后退過程中還伴有材料熱解反應(yīng)、熱解氣體流動(dòng)等復(fù)雜物理化學(xué)過程,使得此時(shí)的逆問題更為復(fù)雜。文獻(xiàn)[3]中,Kanevce利用地面燒蝕試驗(yàn)的溫度測(cè)量數(shù)據(jù)對(duì)酚醛材料的熱物性參數(shù)和材料熱解反應(yīng)系數(shù)進(jìn)行了辨識(shí);文獻(xiàn)[4-5]建立了考慮防熱材料燒蝕情況下的表面熱流反演辨識(shí)算法,但文獻(xiàn)[4]采用的是一層的簡(jiǎn)化燒蝕模型,并且通過給定燒蝕表面溫度來確定外表面位置,這樣處理會(huì)導(dǎo)致辨識(shí)結(jié)果出現(xiàn)非物理振蕩;文獻(xiàn)[5]采用了兩層燒蝕熱解面模型,建立了相應(yīng)的表面熱流辨識(shí)算法,但未用實(shí)測(cè)數(shù)據(jù)進(jìn)行驗(yàn)證。本文在前期變幾何域表面熱流辨識(shí)方法[6]研究的基礎(chǔ)上,借鑒文獻(xiàn)[5]的方法,采用兩層燒蝕熱解面模型,實(shí)現(xiàn)了表面熱流辨識(shí)算法,并用典型試驗(yàn)數(shù)據(jù)對(duì)辨識(shí)算法進(jìn)行了驗(yàn)證。
防熱材料的燒蝕是一個(gè)復(fù)雜的物理化學(xué)過程。本文針對(duì)目前工程上通常采用的碳化復(fù)合材料開展研究。碳化復(fù)合材料的燒蝕模型主要有熱解層模型、熱解面模型和簡(jiǎn)化模型,本項(xiàng)工作中采用的是熱解面模型[1],如圖1示,左邊界為材料受熱邊界,右邊界為防熱材料內(nèi)部邊界,即圖中的s3。熱解面模型是將防熱材料分為兩層:碳化層和原始材料層,熱解層則簡(jiǎn)化為一個(gè)熱解面進(jìn)行處理。圖1中給定碳化邊界和外邊界燒蝕后退規(guī)律,此時(shí)的傳熱區(qū)域包含兩個(gè)運(yùn)動(dòng)邊界:一個(gè)是表面燒蝕后退的邊界,即圖中的s1;一個(gè)是碳化邊界,即圖中的s2,對(duì)應(yīng)的數(shù)學(xué)模型為:
圖1 熱解面模型Fig.1 Pyrogeneration-plane ablation model
本文采用變幾何域的有限控制體積法(ECV:Einite Control Volume method)[7]來對(duì)式(1)進(jìn)行數(shù)值求解,式中的s1(t)和s2(t)可通過兩種方法來給出,一種是給定為測(cè)量結(jié)果;另一種是給定燒蝕碳化溫度后通過迭代計(jì)算來確定其位置。同時(shí),此時(shí)由于計(jì)算域幾何外形隨時(shí)間變化,計(jì)算網(wǎng)格相應(yīng)隨時(shí)間變化,因而在離散方程中需考慮由于控制體邊界運(yùn)動(dòng)而引起的控制體內(nèi)能變化項(xiàng)。該計(jì)算方法處理變幾何域傳熱問題的有效性在文獻(xiàn)[6]中已得到了有解析解的一維半無限區(qū)域燒蝕表面恒溫恒速后退算例的驗(yàn)證。
通常的燒蝕問題計(jì)算是在表面熱流Q(t)已知的情況下,求解式(1)得出材料表面的后退率和內(nèi)部溫度分布。而表面熱流辨識(shí)問題則是在表面熱流Q(t)未知的情況下,根據(jù)材料端面或內(nèi)部測(cè)點(diǎn)的溫度歷程來反演計(jì)算出表面熱流Q(t)的值。不失一般性,下面將測(cè)量點(diǎn)取在右端面處,建立由測(cè)量信息來辨識(shí)表面熱流Q(t)的方法。首先,利用拉格朗日(Lagrange)乘數(shù)法將辨識(shí)問題轉(zhuǎn)化為使選取表面熱流Q(t)值使如下目標(biāo)函數(shù)達(dá)極小的優(yōu)化問題:
式中的λ為伴隨變量。對(duì)此式做分部積分,代入式(1)中的相應(yīng)邊界條件后,再做變分,可知為使式(2)中的目標(biāo)函數(shù)達(dá)極小,伴隨變量需滿足如下伴隨方程:
而目標(biāo)函數(shù)對(duì)Q的梯度為:
在此基礎(chǔ)上即可采用梯度類優(yōu)化算法,如最速下降法和共軛梯度法[6,8],來對(duì)Q(t)進(jìn)行辨識(shí)計(jì)算,以最速下降法為例,優(yōu)化迭代計(jì)算公式為:
式中上標(biāo)“l(fā)”和“l(fā)+1”表示迭代步數(shù);下標(biāo)“i”表示時(shí)間方向上的離散;β為步長(zhǎng)。在優(yōu)化計(jì)算過程中,為克服不適定性的影響,優(yōu)化計(jì)算設(shè)置了如下停止準(zhǔn)則:J≤δ;δ=σ2tf;σ為溫度測(cè)量結(jié)果的標(biāo)準(zhǔn)差。
下面進(jìn)行算法驗(yàn)證:設(shè)燒蝕材料的原始材料層ρ=1300kg/m3,比熱Cp=1507J/kgK,kv=1W/(mK),初始厚度為L(zhǎng)=20mm;碳化后的ρ=406.3kg/m3,比熱Cp=2000J/kg K,kc=3W/(m K);裂解化學(xué)反應(yīng)吸收熱ΔHp=831500 J/kg;熱解面溫度為873K;裂解氣比熱Cg=1300J/kg K。給定表面熱流如圖2中的“Qexact”,給定燒蝕端面s1(t)的后退歷程如圖3中的“s1without noise”,采用第1節(jié)中的數(shù)值方法對(duì)溫度場(chǎng)進(jìn)行求解,s2(t)根據(jù)熱解面溫度通過迭代計(jì)算來確定,圖4示出了計(jì)算出的熱解面位置。
圖2 表面熱流給定值與辨識(shí)結(jié)果對(duì)比Fig.2 Comparison of exact and estimated value of heat flux
圖3 燒蝕后退量S1隨時(shí)間變化歷程Fig.3 Time history of ablated surface recession
圖4 熱解面時(shí)間變化歷程Fig.4 Time history of calculated pyrogeneration-plane
以此時(shí)計(jì)算出的內(nèi)壁測(cè)點(diǎn)溫度歷程為實(shí)測(cè)值,采用第2節(jié)介紹的方法對(duì)表面熱流進(jìn)行辨識(shí),圖2中的虛線“Qestimated(σ=0K)”示出了此時(shí)的辨識(shí)結(jié)果,與給定值符合較好,但在末時(shí)刻有差異,這是由于傳熱存在一定的滯后性,末段時(shí)刻的表面熱流信息在內(nèi)壁溫度測(cè)量結(jié)果中還沒有體現(xiàn)出來,所以不能有效辨識(shí)。
接下來分析測(cè)量噪聲的影響,在計(jì)算出內(nèi)壁溫度歷程上疊加標(biāo)準(zhǔn)差σ=6K(對(duì)應(yīng)測(cè)量的相對(duì)誤差約為1%)的測(cè)量白噪聲來作為測(cè)量結(jié)果(圖5中“exp.”),對(duì)表面熱流進(jìn)行辨識(shí),圖2中的虛線“Qestimated(σ=6K)”示出了此時(shí)的辨識(shí)結(jié)果,與給定值仍符合較好;圖5中的“Eitted”給出了利用辨識(shí)結(jié)果計(jì)算出的測(cè)點(diǎn)溫度歷程,和測(cè)量值符合較好。進(jìn)一步將疊加的測(cè)量白噪聲標(biāo)準(zhǔn)差由σ=6K增大到σ=12K(對(duì)應(yīng)2%的測(cè)量誤差),對(duì)表面熱流進(jìn)行辨識(shí),圖2中的虛線“Qestimated(σ=12K)”示出了此時(shí)的辨識(shí)結(jié)果,可以看到:在相對(duì)誤差2%的情況下,表面熱流的辨識(shí)結(jié)果與給定值進(jìn)一步偏離,但仍較好地反映出了表面熱流的變化趨勢(shì)與強(qiáng)度,辨識(shí)算法具有較好的魯棒性。
圖5 測(cè)點(diǎn)溫度對(duì)比Fig.5 Comparison of temperature measurement
在實(shí)際應(yīng)用中,燒蝕后退量的測(cè)量通常有較大的測(cè)量誤差,因而在圖3中s1(t)上疊加σ=0.11mm(相當(dāng)于燒蝕后退量最大值的10%)的白噪聲作為測(cè)量噪聲,得到的測(cè)量結(jié)果如圖3中的“s1with noise(σ=0.00011m)”示。利用不考慮s1(t)測(cè)量噪聲時(shí)的測(cè)點(diǎn)溫度計(jì)算值作為實(shí)測(cè)值(不疊加測(cè)量噪聲),對(duì)表面熱流進(jìn)行辨識(shí),得到的辨識(shí)結(jié)果與給定值的比較如圖6示,從圖中可以看到,由于燒蝕后退量測(cè)量結(jié)果存在明顯波動(dòng),表面熱流的辨識(shí)結(jié)果中也出現(xiàn)了明顯波動(dòng),二者之間存在較直接的關(guān)聯(lián)關(guān)系。
圖6 表面熱流給定值與辨識(shí)結(jié)果對(duì)比Fig.6 Comparison of exact and estimated value of heat flux
圖7 鈍頭型碳酚醛材料地面燒蝕試驗(yàn)Fig.7 Ablation experiment of blunt Carbon-phenonic material
圖8 燒蝕表面、碳化邊界計(jì)算結(jié)果和表面溫度試驗(yàn)結(jié)果Fig.8 Calculated front surface and interface of ablation and experimental surface temperature
圖9 表面熱流時(shí)間歷程辨識(shí)結(jié)果和加熱熱流值對(duì)比Fig.9 Comparison of estimated heat flux and heating power
下面通過一地面試驗(yàn)的實(shí)例來對(duì)辨識(shí)方法的有效性做進(jìn)一步分析驗(yàn)證。鈍頭型碳酚醛材料(Narmco4028)試件在Langley中心11英寸陶瓷加熱風(fēng)洞中的燒蝕試驗(yàn)示意圖[9]如圖7所示,材料的幾何參數(shù)在圖中標(biāo)注,來流中氧的質(zhì)量分?jǐn)?shù)為0.13,駐點(diǎn)壓力5.88atm,駐點(diǎn)焓值2.55MJ/kg。圖8為文獻(xiàn)[9]中給出的燒蝕表面、碳化邊界的后退歷程計(jì)算值(圖中分別記為"front surface"和"interface")和表面溫度測(cè)量結(jié)果(圖中記為"front surface temperature"),其中燒蝕表面和碳化邊界在試驗(yàn)結(jié)束時(shí)刻的計(jì)算結(jié)果與實(shí)測(cè)值吻合?;谶@組結(jié)果,利用前面建立的辨識(shí)算法,可辨識(shí)出表面加熱熱流的時(shí)間變化歷程,如圖9中“Qestimated”示??梢钥吹?,熱流辨識(shí)結(jié)果和陶瓷加熱風(fēng)洞的加熱熱流3.43MW/m2是基本符合的,偏差約20%。這一結(jié)果表明:本文所建立的考慮燒蝕情況下表面熱流辨識(shí)方法是基本有效的,經(jīng)過更多算例的檢驗(yàn)校核后有望在工程實(shí)際中得到更進(jìn)一步的應(yīng)用。
本文針對(duì)燒蝕傳熱問題,在熱解面模型的基礎(chǔ)上,推導(dǎo)了伴隨方程,建立了基于測(cè)點(diǎn)溫度辨識(shí)表面熱流的方法,并對(duì)算法的有效性進(jìn)行了算例考核。結(jié)果表明:在不考慮測(cè)量噪聲時(shí),熱流辨識(shí)結(jié)果與真值符合較好;隨著測(cè)量誤差的增大,辨識(shí)結(jié)果與真值之間的偏差增大;燒蝕后退量測(cè)量結(jié)果的誤差對(duì)辨識(shí)結(jié)果有較為顯著的影響。此后,將該辨識(shí)方法用于鈍頭型碳酚醛材料Narmco4028試件在陶瓷加熱風(fēng)洞中的燒蝕試驗(yàn)結(jié)果的分析,結(jié)果表明辨識(shí)出的表面熱流與加熱功率基本符合,辨識(shí)方法是有效的,在工程實(shí)際中有較好的應(yīng)用前景。
[1]姜貴慶,劉連元.高速氣流傳熱與燒蝕熱防護(hù)[M].北京:國(guó)防工業(yè)出版社,2003.
[2]王希季,林華寶,李頤黎,等.航天器進(jìn)入與返回技術(shù)(下)[M].北京:宇航出版社,1991.
[3]KANEVCE L P,KANEVCE G H,ANGELEVSKI Z Z.Comparison of two kinds of experiments for estimation of thermal properties of ablative composite[R].In:Inverse Problems in Engineering:Theory and Practice[C].3rd Int.Conference on Inverse Problems in Engineering.Port Ludlow,WA,USA,1999.
[4]OLIVEIRA A P D,ORLANDE H R B.Estimation of the heat flux at the surface of ablating materials by using temperature and surface position measurements[J].Inverse Problems in Science and Engineering,2004,12(5):563-577.
[5]HAKKAKI-EARD A,KOWSARY E.Heat flux estimation in a charring ablator[J].Numerical Heat Transfer,Part A,2008,53:543-560.
[6]邵元培,錢煒祺,周宇,等.變幾何域傳熱的表面熱流反演方法[J].計(jì)算力學(xué)學(xué)報(bào),2013,30(2):296-301.
[7]AMAR A J,BLACKWELL B E,EDWARDSJ R.One-dimensional ablation using a full Newton’s method and finite control volume procedure[J].Journal of Thermophysics and Heat transfer,2008,22(1):71-82.
[8]錢煒祺,周宇,何開鋒,等.非線性熱傳導(dǎo)逆問題的表面熱流辨識(shí)方法[J].空氣動(dòng)力學(xué)學(xué)報(bào),2012,30(2):145-150.
[9]SUTTON K.An experimental study of a carbon-phenolic ablation material[R].NASA TND-5930,1970.
Heat flux estimation for heat transfer problem with ablation
QIAN Weiqi1,2,ZHOU Yu1,2,HE Kaifeng1,2,SHAO Yuanpei1,2
(1.State Key Laboratory of Aerod ynamics,Mianyang Sichuan 621000,China;2.China Aerodynamics Research and Development Center,Mianyang Sichuan 621000,China)
Ablation is one of the main means of thermal protection for hypersonic flight vehicles.In this paper,based on pyrogeneration-plane ablation model and temperature measurement at some specific locations,the surface heat flux estimation method isdeveloped with adjoint equation analysis and validated by numerical examples.It is shown that the estimated heat flux agrees with the exact value well,and thedifference between the estimated heat flux and the exact value grows as the measurement noise increases.It is also revealed that the measured noise in the ablated surface recession value may have a significant influence on the estimated heat flux.Eurthermore,the estimation method is utilized to analyze the experimental ablation results of the blunt Carbon-phenonic material Narmco4028 in ceramic-heated tunnel,it is found that the estimated surface heat flux value is close to the heating power of the arc-heater.This result further verifies that the estimation method is effective and may have a bright prospect in engineering practice.
ablation;heat flux estimation;pyrogeneration-plane ablation model
V211.3
Adoi:10.7638/kqdlxxb2014.0038
0258-1825(2014)06-0772-05
2014-05-16;
2014-11-06
國(guó)家自然科學(xué)基金項(xiàng)目(11372338);空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室基金資助(JBKY11030903)
錢煒祺(1973-),男,江蘇無錫,研究員,主要從事空氣動(dòng)力學(xué)、飛行力學(xué)與數(shù)理反問題研究.
錢煒祺,周宇,何開鋒,等.考慮燒蝕情況下的表面熱流辨識(shí)[J].空氣動(dòng)力學(xué)學(xué)報(bào),2014,32(6):772-776.
10.7638/kqdlxxb-2014.0038 QIAN W Q,ZHOU Y,HE K E,et al.Heat flux estimation for heat transfer problem with ablation[J].ACTA Aerodynamica Sinica,2014,32(6):772-776.