碳基復(fù)合材料熱化學(xué)燒蝕模型

碳基復(fù)合材料熱化學(xué)燒蝕模型

在模型的超聲速和高聲速再入層的過程中，由于粘性散散效應(yīng)和強光束激發(fā)的強烈壓縮，巨大的動能損失中的一部分轉(zhuǎn)化為激發(fā)層中的氣體內(nèi)部。空氣在穿過頭激波時被加熱,溫度突然升高,使空氣分子產(chǎn)生離解和電離,空氣比熱比不再是常數(shù)。隨著高度的進一步降低,雷諾數(shù)增加,在表面附近形成逐漸變薄的邊界層,由于粘性作用,溫度進一步升高,一般可達到7000～8000K。因此研究熱防護系統(tǒng)燒蝕問題多年來一直是航天技術(shù)的一個關(guān)鍵問題。早期研究對于馬赫數(shù)不高的情況下,大多采用無粘性理論求解邊值條件,然后再確定熱邊界層的摩擦系數(shù),利用雷諾相似定律求出熱交換系數(shù),如層流的Lees方法、湍流的參考焓法。隨著再入馬赫數(shù)的提高,此時必須考慮物面氣流與材料物理化學(xué)性質(zhì)的相互作用,但是限于當(dāng)時的計算條件,只按固定邊界的假設(shè),用化學(xué)邊界層方程求出燒蝕質(zhì)量和壁面的對流熱流,再由壁面能量平衡條件結(jié)合熱傳導(dǎo)方程求解材料內(nèi)部的熱響應(yīng)。如今隨著航天技術(shù)的發(fā)展,結(jié)構(gòu)輕質(zhì)化、小型化和飛行、再入速度的大幅度提高,作為熱防護的碳基復(fù)合材料具有耐高溫、耐燒蝕、抗沖擊等特點,并且通過自身燒蝕引起質(zhì)量損失,吸收并帶走大量的熱量,阻止外部熱量向結(jié)構(gòu)內(nèi)部傳遞,從而保護內(nèi)部結(jié)構(gòu)在一定溫度范圍內(nèi)正常工作。燒蝕引起結(jié)構(gòu)宏觀和材料的細觀形貌變化,直接影響宏觀、細觀流場特性;材料的燒蝕過程伴隨著復(fù)雜的物理化學(xué)變化,同時也是材料與熱環(huán)境的耦合作用過程,文獻等通過細觀結(jié)構(gòu)流場的熱環(huán)境分析,結(jié)合材料表面的細觀材料組成,材料表面的細觀燒蝕結(jié)構(gòu)及材料的熱化學(xué)燒蝕性質(zhì),研究常壓高焓條件下材料在細觀尺度上的燒蝕機理,并給出材料燒蝕后退率的預(yù)測方法;結(jié)構(gòu)內(nèi)部溫度劇烈變化必然產(chǎn)生復(fù)雜的熱應(yīng)力,出現(xiàn)熱裂破壞、熱損毀,因此熱力荷載引起的瞬態(tài)復(fù)雜熱應(yīng)力,對結(jié)構(gòu)的強度和安全性至關(guān)重要。本文作者主要研究燒蝕碳化學(xué)反應(yīng)機理、燒蝕模型的建立,基于移動邊界條件下的瞬態(tài)熱邊界、物性隨溫度變化的物體內(nèi)非線性傳熱特性,對瞬態(tài)的熱力耦合進行有限元求解計算;并對熱應(yīng)力的演化規(guī)律進行分析。通過數(shù)值模擬表征了材料在高溫?zé)g條件下的動態(tài)熱力學(xué)行為;計算結(jié)果吻合電弧加熱器燒蝕試驗對比結(jié)果。因此,可以實現(xiàn)高溫、高過載等復(fù)雜因素耦合作用下材料燒蝕的實時仿真。1氧化燒蝕速率碳在較低溫度下首先氧化,一般在低于1000K時,氧化過程是速率控制,燒蝕速率由表面化學(xué)動力學(xué)決定。當(dāng)溫度超過1700K后,氧化速率急劇增加,燒蝕速率完全受氧化擴散控制,表面氧氣被完全燃盡。在擴散控制的條件下,材料表面的燒蝕組元及空氣組元滿足表面的質(zhì)量守恒原理和表面元素的熱化學(xué)平衡原理。1.1化學(xué)元素k的濃度和質(zhì)量擴散流率對于在燒蝕條件下材料表面具有k種元素的質(zhì)量交換系統(tǒng),元素的質(zhì)量守恒如圖1所示。其表面的質(zhì)量守恒形式為～Jkw+(ρv)w～ckw=˙ms～cks(1)J～kw+(ρv)wc～kw=m˙sc～ks(1)式中:～JkwJ～kw為化學(xué)元素k在壁面的擴散流率;ρ為壁面氣體的密度;v為壁面引射氣體的速度;～ckwc～kw為元素k在壁面的質(zhì)量濃度;˙msm˙s為壁面固體的質(zhì)量損失率;～cks為壁面固體的質(zhì)量損失中元素k的質(zhì)量濃度。對k種元素進行求和,則得∑k～Jkw+(ρv)w=˙ms(2)質(zhì)量擴散流率Jkw通常通過傳導(dǎo)勢的方法求得,元素k在壁面擴散流率為Jkw≈ρeueCΜ(Ζ*kw-Ζ*ke)(3)其中Z*是質(zhì)量比數(shù)和摩爾比數(shù)的擴散系數(shù)加權(quán)平均量。在等擴散系數(shù)的條件下～ck與Z*k是等效的,則有(ρv)w=˙ms(4)1.2化學(xué)平衡狀態(tài)對于具有i種組元Ni,k種元素的熱化學(xué)平衡系統(tǒng),一般可選取k種基本組元Nk,則j=i-k種非基本組元Nj的生成反應(yīng)可表示為∑kνkjΝk→Νj(5)這里νkj是基本組元的化學(xué)計量系數(shù)。由于系統(tǒng)處于化學(xué)平衡狀態(tài),由熱力學(xué)第二定律熵增為零條件可得如下的依賴于溫度的化學(xué)平衡常數(shù)的關(guān)系式:Κpj(Τ)=pj∏kp-νkjk(6)氣體組元滿足道爾頓分壓定律,即p=∑ipi(7)同時根據(jù)元素當(dāng)量濃度及平均分子量的定義,可得系統(tǒng)內(nèi)元素的質(zhì)量濃度及平均分子量的關(guān)系式～ck=ΜkpˉΜ∑ixkipi(8)ˉΜ=1p∑iΜiΡi(9)式中:p為混合氣體的總壓力;Mk為k元素的分子量;pi為i氣體的分壓;xki為i組元中元素k的原子數(shù)目。對于封閉的熱化學(xué)平衡系統(tǒng),上述關(guān)系式在已知元素比數(shù)的條件下是完備可解的,但在燒蝕的過程中,由于邊界層內(nèi)的質(zhì)量交換,在邊界層內(nèi)的熱化學(xué)平衡系統(tǒng)中應(yīng)滿足表面的質(zhì)量守恒原理。1.3whcs+aqr-4w10的熱價計算材料表面的能量守恒如圖2所示。材料表面的守恒關(guān)系為qΝ=qw-[˙mwhw-˙mwhCs]+aqr-εσΤ4w(10)式中:qN為進入材料內(nèi)部的凈熱流;qw為外界對材料壁面的對流加熱;˙mwhw為壁面氣體帶走的熱量;˙mwhCs為固壁質(zhì)量損失帶走的熱量;aqr為外界對壁面的輻射加熱量;εσT4w為壁面對外界的氣體輻射量。2熱平衡方程的建立根據(jù)進入壁面的凈熱流即可確定物體溫度。燒蝕移動邊界傳熱問題的主要特點:控制方程為非穩(wěn)態(tài)導(dǎo)熱方程,區(qū)域內(nèi)存在著一個隨時間移動的兩相界面,在該界面上吸收潛熱,傳熱和相變一起發(fā)生。C/C材料在高溫下的材料性能與溫度有關(guān),熱傳導(dǎo)方程是非線性的,又由于相變的存在,邊界條件也是非線性的。變熱物性實際上把求解溫度場的問題與求解物性場問題聯(lián)系在一起,互相耦合,將溫度場與物性場進行耦合迭代求解。先假設(shè)一個物性場,在該物性場下求解溫度場,由求得的溫度場再求解物性場,循環(huán)交替進行,直至溫度場和物性場不再變化為止。物體的初始溫度分布為T(x,y,z,t),求解該物體的溫度場便歸結(jié)為尋找滿足下列微分方程和定解條件的溫度分布函數(shù),熱傳導(dǎo)微分方程為?Τ?x(Κx?Τ?x)+?Τ?y(Κy?Τ?y)+?Τ?z(Κz?Τ?z)=ρc(?Τ/?t)(11)式中:c=c(T);K=K(T);ρ為材料的密度;c為比熱容;K為熱傳導(dǎo)系數(shù)且具有方向性;t為時間;T為材料表面的溫度。材料的非燒蝕表面為絕熱邊界條件,燒蝕表面為對流加熱條件。根據(jù)變分原理,引入虛溫度(δT),并對方程(11)進行體積積分,則有∫V[?Τ?x(Κx?Τ?x)+?Τ?y(Κy?Τ?y)+?Τ?z(Κz?Τ?z)]δΤdV=∫Vρc(?Τ/?t)δΤdV利用有限元方法離散方程,單元內(nèi)任意點溫度可由節(jié)點溫度表示為Τ=m∑iΝiΤi式中:Ni和Ti分別為單元形函數(shù)和第i節(jié)點的溫度;m為單元所包含的節(jié)點個數(shù),經(jīng)過一系列推導(dǎo),最后的離散方程可表示為[C]{˙Τ}+[Κ]{Τ}={F}式中:[C]是比熱矩陣;{˙Τ}是結(jié)點溫度對時間的導(dǎo)數(shù)向量;{T}是結(jié)點溫度向量;[K]是熱導(dǎo)率矩陣;{F}是熱荷載向量。對有限元數(shù)值求解時,必須考慮當(dāng)材料發(fā)生燒蝕時,物體表面不斷向后退縮。因此,考慮燒蝕情況下的溫度場、熱力耦合應(yīng)力場計算實際上是求解移動邊界的溫度場、熱力耦合應(yīng)力場問題。以表面一點為例,節(jié)點m的溫度變化率可以寫成?Τm?t=?Τm?t+?Τm?y?y?t(12)與固定邊界條件相比,方程右邊第二項表示有節(jié)點m的坐標移動所引起的溫度變化。如何處理此項,是計算移動邊界條件問題的關(guān)鍵所在。以下用節(jié)點熱平衡方程來說明節(jié)點坐標移動的物理意義及處理方法。如果在Δt開始時刻的初始溫度分布為T=f(y,t),在Δt終了時刻,由于材料燒蝕的結(jié)果,外表面移動至坐標ˉy處,在分割數(shù)N不變的情況下,此時的分割長度為Δy′=(δ-ˉy)/Ν由內(nèi)部節(jié)點m的熱平衡方程FΔtλΤm-1-ΤmΔy=FΔtλΤm-Τm+1Δy+FΔtρcp(Τm-Τ0m)得出λΤm-1-ΤmΔy′=λΤm-Τm+1Δy′+ρcp(Τm-Τ0m′)(13)式(13)中T0m′表示與時間間隔Δt終了時節(jié)點m的位置相對應(yīng)處在Δt開始時刻的溫度。根據(jù)開始時刻各個節(jié)點的溫度及坐標位置,用插值方法可以求得T0m′的值。即用有限元處理移動邊界時,根據(jù)每一瞬時的燒蝕外形,計算出各個節(jié)點的坐標,然后用插值的方法求出各個節(jié)點坐標相對應(yīng)Δt開始時刻的初始溫度。3應(yīng)力-應(yīng)變關(guān)系有限元求解方程由于熱變形產(chǎn)生的應(yīng)變可以看作為物體的初應(yīng)變。其中ε0為溫度變化引起溫度應(yīng)變,它現(xiàn)在是作為初應(yīng)變出現(xiàn)在應(yīng)力-應(yīng)變關(guān)系式中。通過虛位移原理可以得到包含溫度應(yīng)變在內(nèi)、用以求解熱應(yīng)力問題的最小位能原理,泛函表達式如下:Πp=∫Ω(12εΤDε-εΤDε0-uΤf)dΩ-∫ΓσuΤˉΤdΓ(14)將求解域Ω進行有限元離散,從δΠp=0可以得到有限元求解方程。有限元求解方程的荷載向量中包括由溫度應(yīng)變引起的溫度荷載,用虛功原理的應(yīng)力-應(yīng)變關(guān)系式代替一般的應(yīng)力-應(yīng)變關(guān)系有[k]e[δ]e={pf}e+{pΤ}e+{pε0}e(15)式中:[k]e為單元剛度矩陣;[δ]e為單元位移列陣;{pf}e、{pT}e為單元體積荷載和表面荷載引起的等效荷載列陣;{pε0}e為考慮溫度變化而增加的荷載列陣,相當(dāng)于考慮荷載變化而施加于結(jié)點的假想等效荷載,稱為溫度荷載。其有限元解法與前述溫度場解法相似;考慮材料物性隨溫度的變化,材料為橫觀各向同性,在應(yīng)力、應(yīng)變與溫度的關(guān)系式成為非線性的情況下,把求解熱力耦合應(yīng)力場的問題與求解物性場問題聯(lián)系在一起,互相耦合,將應(yīng)力場與物性場進行耦合迭代求解。4c/c高溫?zé)g材料模擬系統(tǒng)的試驗結(jié)果本文中依據(jù)以上原理,對高溫、超音速氣流的流場和防熱材料燒蝕行為進行數(shù)值模擬。首先,對流場計算一個時間步長,計算此時刻的邊界熱流,然后根據(jù)邊界熱流對結(jié)構(gòu)燒蝕計算一步,得到一個新的后退邊界條件下的物體溫度場分布及物面溫度,再由新的物體溫度場分布及結(jié)構(gòu)邊界條件計算應(yīng)力場,這就是一步完整的耦合計算。試件為C/C平頭圓柱體,為了表示方便,將燒蝕面中心至邊緣的半徑上各點作為表示流-熱參數(shù)的位置,如圖3所示?？紤]材料的方向性,利用實驗測得的C/C材料變物性性能如下:比熱(J·g-1℃-1)C=0.997+2.930×10-5Τ熱導(dǎo)率(Cal·cm-1·s-1·℃-1)Κx=Κy=172.528-0.236Τ+1.994×10-5Τ2-7.792×10-9Τ3+1.140×10-11Τ4Κz=82.785-0.103Τ+9.846×10-6Τ2-4.372×10-9Τ3+7.315×10-13Τ4膨脹系數(shù)(10-6℃-1)αx=αy=-0.4315+0.00105Τ-3.06818×10-7Τ2αz=-0.37867+4.7×10-4Τ+7.57574×10-8Τ2彈性模量(GPa)Ex=Ey=5.28961-0.031681Τ+3.57617×10-5Τ2-1.08245×10-8Τ3+7.87254×10-13Τ4Ez=3.03906-0.01496Τ-3.58779×10-5Τ2+2.69612×10-8Τ3-5.98254×10-12Τ4剪切模量(GPa)Gxz=Gyz=6.39,Gxy=8.878泊松比γxz=γyz=0.28,γxy=0.31由計算獲得燒蝕各個時刻(間隔3s)的模擬結(jié)果。流場加熱熱流見圖4;進入材料內(nèi)部凈熱流見圖5;物面溫度見圖6;燒蝕表面后退位置見圖7。電弧加熱器防熱材料模擬系統(tǒng)可以提供的連續(xù)可控超高溫、超音速等離子流,能夠基本再現(xiàn)飛行熱環(huán)境中材料的復(fù)雜反應(yīng),并通過先進的檢測技術(shù)和手段,對熱化學(xué)燒蝕產(chǎn)物的成分、形態(tài)、含量,燒蝕試件溫度分布以及燒蝕面的質(zhì)量損失過程等進行實時檢測、分析。以此解決燒蝕過程數(shù)值仿真所需參數(shù)及依據(jù),揭示防熱材料的高溫?zé)g機理,實踐與理論預(yù)報相互校驗等問題。在射流燒蝕環(huán)境下,通過測試C/C復(fù)合材料試件的溫度變化以及試件燒蝕表面的形態(tài)變化等過程,給出了狀態(tài)參數(shù)和材料的質(zhì)量燒蝕率及表面后退率(線燒蝕率)見表1。C/C材料燒蝕速度的測量結(jié)果與數(shù)值仿真結(jié)果比較表明,二者是相符的。實驗中分別采用輻射式高溫計及熱電偶,測量試件燒蝕面和試件沿軸向各點隨時間變化的溫度。紅外測溫儀測量范圍為900～3000℃,試件軸向的溫度采用Ni/Cr-Ni/Al熱電偶,最高測量溫度為1370℃。將熱電偶用高溫陶瓷(BN)保護套封裝,沿試件軸向植入試件如圖8所示,通過數(shù)據(jù)采集模塊與計算機相連,保存溫度信號。通過A、B、C測點的數(shù)值計算與試驗測量的比對,表明計算結(jié)果與試驗測量結(jié)果相符,如圖9所示。在氣流壓力與熱共同作用下引起內(nèi)應(yīng)力。如在表1中對于狀態(tài)2,圖10給出燒蝕5s、10s時x方向應(yīng)力場計算結(jié)果;應(yīng)力場的應(yīng)力分布并不均勻,存在多個大小不等的應(yīng)力區(qū)即應(yīng)力核;仍以x方向為例,通過兩個時刻5s與10s應(yīng)力分布對比,可以看出,應(yīng)力核的大小和其內(nèi)部應(yīng)力方向是隨時間變化的。這一現(xiàn)象的形成是因為在模型上下方向存在溫度梯度,而水平剖面上溫度差異較小,同層各個點自由變形應(yīng)該相同,但核心材料受外部材料的約束,失去了自由變形的機會,形成應(yīng)力核;此外,變形協(xié)調(diào)性影響又使試件分成上下多個核心區(qū)。因此,動態(tài)的熱邊界條件導(dǎo)致一個動態(tài)的溫度場,進而導(dǎo)致一個動態(tài)的熱應(yīng)力場。又由于溫度