鈦鎳合金固體火箭發(fā)動機一維非穩(wěn)態(tài)導(dǎo)熱及熱應(yīng)力研究

（南昌航空大學(xué)飛行器工程學(xué)院，南昌，330063） 薛志遠 胡曉安

鈦鎳合金（TiNi）是近年來比較熱門的金屬材料，正被廣泛運用于電子、機械、醫(yī)療、能源、航空航天等領(lǐng)域。鈦鎳合金最為顯著的優(yōu)點是比強度高和耐腐蝕性好，具有優(yōu)良的形狀記憶和超彈性，具有比較好的機械性能，晶粒度很細小，可以被加工成很細小的絲材和薄板[1]。由于鈦合金具有優(yōu)異的熱強度以及較高的比強度，耐高溫，高溫環(huán)境下化學(xué)穩(wěn)定性好，得以在發(fā)動機上被逐漸使用。未來的航空航天飛行器及其推力系統(tǒng)需要擁有更高強度的、更高工作溫度、更高彈性模量、密度更高且更加經(jīng)濟的材料，鈦鎳高溫合金便成了主要的材料選擇對象之一。

基于厚壁圓筒溫度分布以及熱應(yīng)力的研究國內(nèi)外學(xué)者已經(jīng)做了諸多探索。美國學(xué)者Paton，Koppel[2]等人在上世紀八十年代對電感加熱器以及核燃料元件裝置建立厚壁圓筒模型，給出了具有內(nèi)熱源的圓筒熱應(yīng)力計算通式。日本學(xué)者平修二[3]在上世紀八十年代末，為了強化當(dāng)時工程技術(shù)人員對于設(shè)計制造高溫高壓下使用的動力機械或化學(xué)反應(yīng)裝置的熱疲勞以及高溫強度的理解，依據(jù)日本材料協(xié)會高溫強度委員會有關(guān)熱應(yīng)力和熱疲勞方面的研究成果，在對高溫下機械設(shè)備破損原因進行大量分析之后，出版了《熱應(yīng)力與熱疲勞》一書，系統(tǒng)闡述了熱應(yīng)力分析方法，并廣泛涉及了基礎(chǔ)理論和設(shè)計應(yīng)用；丸山武志[4]等人完善了關(guān)于多層圓通導(dǎo)熱系數(shù)的測定方法以及不穩(wěn)定溫度分布理論，該方法適用于高壓下測定各種多層化學(xué)反應(yīng)容器的導(dǎo)熱系數(shù)的測定，同時對容器壁厚沒有限制，具有較高的經(jīng)濟性和適用性，不穩(wěn)定溫度分布理論可以幫助測定容器焊接接頭附近由于溫度差引起的熱應(yīng)力，導(dǎo)熱系數(shù)的方法，為后來開展圓筒導(dǎo)熱研究奠定了基礎(chǔ)。在國內(nèi)，國防科技大學(xué)教授孟祥貴[5]等人為了對厚壁圓筒在受到溫度場的作用時在內(nèi)部產(chǎn)生的熱應(yīng)力進行精確分析，提出關(guān)于圓筒熱應(yīng)力問題的傳遞函數(shù)方法，該傳遞函數(shù)方法具有較大的靈活性和適應(yīng)性；青海師范大學(xué)楊能彪[3]等人針對瞬態(tài)導(dǎo)熱過程中，由于時間因素相關(guān)的非穩(wěn)態(tài)項是對時間的導(dǎo)數(shù)問題，給出了一維瞬態(tài)導(dǎo)熱問題的數(shù)值解答。

由于目前針對固體發(fā)動機各典型截面在非設(shè)計狀態(tài)（啟動過程）中，經(jīng)受瞬態(tài)溫度（物體內(nèi)部各點溫度隨時間改變叫做瞬態(tài)溫度，這種變化的溫度的熱量傳遞過程稱為非穩(wěn)態(tài)熱傳導(dǎo)）載荷條件下的瞬態(tài)溫度場以及熱應(yīng)力分析還有待進一步的研究，利用等厚壁圓筒簡化模型，通過研究以鈦鎳基合金為材料的某型固體發(fā)動機瞬態(tài)溫度和熱應(yīng)力分布，為固體發(fā)動機實驗時一維非穩(wěn)態(tài)導(dǎo)熱和熱應(yīng)力計算提供基本思路，同時也加深對解鈦鎳合金的認識，推動鈦鎳合金在發(fā)動機領(lǐng)域的運用。

1 傳熱分析模型

1.1 物理模型

固體火箭發(fā)動機燃燒室、發(fā)動機喉部以及發(fā)動機喉部出口某截面，認為溫度與周向和軸向無關(guān)，只是徑向坐標的函數(shù)，從而簡化為等厚壁圓筒。參照某型固體火箭發(fā)動機結(jié)構(gòu)尺寸，該等厚壁圓筒內(nèi)徑ra=1.8m，外徑rb=2m。圓筒結(jié)構(gòu)如圖1所。圓筒內(nèi)壁為高溫燃氣，由于發(fā)動機燃燒室、發(fā)動機喉部以及發(fā)動機喉部有不同的工作溫度，因而所選取的典型截面處具有不同的瞬態(tài)溫度值。內(nèi)壁有高溫燃氣充當(dāng)內(nèi)熱源，給定瞬態(tài)溫度邊界值。外壁與外界大氣相接觸，存在對流換熱。

圖1 圓筒結(jié)構(gòu)平面圖

1.2 物理模型的簡化

該等厚壁圓通為軸對稱旋轉(zhuǎn)體結(jié)構(gòu)，各典型截面所在圓筒溫度在軸向上無變化，溫度分布只與半徑和時間有關(guān)，故可以把圓筒投影到平面上。于是，該固體發(fā)動機導(dǎo)熱問題便簡化為一維下的非穩(wěn)態(tài)導(dǎo)熱。簡化后的圓筒一維模型在極坐標中進行分析，如圖2所示。

圖2 簡化后的圓筒一維模型在極坐標中分析圖

1.3 鈦鎳合金材料熱物性參數(shù)

通常在實際工程材料中，對于金屬材料而言，其物理性質(zhì)大多和溫度以及物質(zhì)的種類有關(guān)。由于發(fā)動機啟動過程溫度變化范圍很大，采用線性近似關(guān)系來確定相關(guān)熱物性。一般而言，熱膨脹系數(shù)主體上隨溫度升高從而線性增大，彈性模量隨著溫度的增加而大幅下降，泊松比隨溫度變化不明顯，導(dǎo)熱系數(shù)隨著溫度上升趨于降低，定壓比熱在溫度上升時也相應(yīng)增加[6]。這些參數(shù)均可作為溫度的線性函數(shù)。在非穩(wěn)態(tài)過程導(dǎo)熱中，當(dāng)給定初始條件和邊界條件時，導(dǎo)熱過程受材料熱物性中的導(dǎo)熱系數(shù)的影響十分明顯，特別是當(dāng)導(dǎo)熱系數(shù)是溫度的函數(shù)時，導(dǎo)熱過程遭受的影響更加顯著[6]。在溫度變化較大的發(fā)動機啟動工況下，故彈性模量、導(dǎo)熱系數(shù)、熱膨脹系數(shù)、定壓比熱以及換熱系數(shù)分別表示為E、λ、α、Cp、h，它們均為溫度的函數(shù)查材料參數(shù)手冊，相應(yīng)溫度下鈦鎳合金的熱物性參數(shù)如下表1所示，其余時刻對應(yīng)溫度下的材料熱物性采用線性方法進行材料參數(shù)的內(nèi)插和外推得到。鈦鎳合金密度為ρ=4.5*103kg/m3，密度為定值，不隨溫度變化。

1.4 基本假設(shè)

該等厚壁圓筒屬于軸對稱結(jié)構(gòu)，圓筒軸向上橫截面積為定值，各典型截面所在圓筒上的溫度在軸向上不發(fā)生變化，因而可取某個截面進行分析。溫度隨時間不斷變化的非穩(wěn)態(tài)過程；該圓筒內(nèi)壁均勻受熱，為純導(dǎo)熱作用沒有對流，按第一類邊界條件處理；外表面與外界干空氣接觸，涉及對流換熱，采用第三類邊界條件。圓環(huán)初始溫度場為外界干空氣溫度，取為T0=Ta0=Tb0=288K。

2 導(dǎo)熱過程的計算分析

當(dāng)相位角?=π/2時，圓環(huán)非穩(wěn)態(tài)導(dǎo)熱方程為

表1 材料熱物性參數(shù)

給定初始條件：當(dāng)t=0時，T=T0

邊界條件：

其中，T∞為外部大氣溫度

3 數(shù)值離散

對時間坐標進行考慮，因而一維瞬態(tài)導(dǎo)熱就變成了二維問題，擁有時間坐標t和空間坐標r兩個變量。但是時間坐標t是單向的，也就是前一時刻的狀態(tài)只能對后一時刻的狀態(tài)產(chǎn)生影響，反之則不產(chǎn)生影響。圖3顯示出了以r和t為坐標的環(huán)形區(qū)域的離散，時間t從0時刻開始，增量步為20s，經(jīng)過一個個時層增加到k、k+1時層[7]。

圖3 r-t環(huán)形區(qū)域的離散

3.1 微分方程的離散

對導(dǎo)熱方程（1）進行數(shù)值離散，則可得到離散的控制方程：

3.2 邊界條件的離散

用K+1時刻相應(yīng)節(jié)點的差分，替換式（2）（3）中的微分，得到邊界節(jié)點的差分方程顯式解：

4 各典型截面的溫度場的計算

參照某型固體發(fā)動機實際工況，計算等厚壁固體發(fā)動機圓筒各典型截面各節(jié)點的溫度值，并和ABAQUS計算得到的溫度值進行對比，以論證數(shù)值計算結(jié)果的準確性。

4.1 該型固體發(fā)動機地面試車設(shè)計點飛行條件

Ma=0，H=0，T0=288K。穩(wěn)定工作時燃燒室出口溫度為T4=1182K。

地面試車時，啟動工況是發(fā)動機從停止狀態(tài)到穩(wěn)定運行狀態(tài)的過渡過程，此過渡過程每一瞬時的燃氣溫度、負荷都不相同，圓筒內(nèi)側(cè)邊界溫度以及外側(cè)邊界換熱條件都和時間有關(guān)，因而這一過渡階段是非穩(wěn)態(tài)的，本論文對這一階段進行研究。

4.2 各典型截參考節(jié)點的溫度值

將表1數(shù)據(jù)代入式（4）（5）中得各截面各節(jié)點的溫度分布，其中，i=1,2,…,6，K=0,1,2…,13。利用MATLAB繪制出各典型截面上不同位置節(jié)點處溫度隨時間變化圖像如圖4所示：

圖4 燃燒室出口、噴管喉部、噴管出口截面各參考節(jié)點溫度曲線圖

4.3 應(yīng)用ABAQUS對各典型截面進行溫度場分析

如圖5、圖6、圖7各典型截面溫度場分布云圖所示?？梢姼鹘孛嫠矐B(tài)溫度場分布和實際數(shù)值計算所得溫度分布大致相同。

圖5 發(fā)動機燃燒室截面溫度場

圖6 噴管喉部截面溫度場

5 基于MATLAB的等厚壁圓桶熱應(yīng)力計算分析

物體由于熱變形使受約束物體內(nèi)產(chǎn)生的應(yīng)力稱為熱應(yīng)力。熱應(yīng)力的產(chǎn)生有兩種原因：結(jié)構(gòu)內(nèi)外表面受熱不均，產(chǎn)生溫度梯度；受到外部約束，這樣熱變形就要受到內(nèi)部各部分的相互制約和外界的限制，從而產(chǎn)生熱應(yīng)力。等厚壁圓筒形結(jié)構(gòu)，內(nèi)部存在高溫氣體作為熱源，內(nèi)外表面形成溫度梯度，因而產(chǎn)生了較大的熱應(yīng)力[6]。

不考慮軸向應(yīng)力，圓筒徑向、周向應(yīng)力關(guān)于溫度的函數(shù)如下：

徑向應(yīng)力：

周向應(yīng)力：

其中，a為熱膨脹系數(shù)，ra為圓筒內(nèi)徑，rb為圓桶外徑，E為彈性模量。利用MATLAB編制徑向與周向熱應(yīng)力計算程序，即把由（4）（5）式得到的各節(jié)點，運用曲線擬合法，得到同一時刻徑向上各個節(jié)點溫度值關(guān)于半徑r的函數(shù)關(guān)系，再代入（6）（7）式中即可得到相應(yīng)的徑向熱應(yīng)力和周向應(yīng)力，熱應(yīng)力隨時間變化的關(guān)系曲線如下圖8、圖9、圖10所示。

圖8 燃燒室出口截面半徑方向上各節(jié)點徑向、周向熱應(yīng)力

圖9 噴管喉部截面半徑方向上各節(jié)點徑向、周向熱應(yīng)力

圖10 噴管出口截面半徑方向上各節(jié)點徑向、周向熱應(yīng)力

6 討論

基于傳熱學(xué)理論，參考單層無限長圓筒的溫度分布一般情況，對一維圓筒模型的非穩(wěn)態(tài)導(dǎo)熱方程用差分法進行數(shù)值離散，借助于MATLAB計算出圓筒半徑方向上所選取的節(jié)點各時間節(jié)點處的瞬態(tài)溫度解析計算值，同時用ABAQUS有限元分析軟件對圓筒導(dǎo)熱過程進行模擬得到瞬態(tài)溫度場的有限元解法，與解析計算結(jié)果進行比較。再采用曲線擬合的方法，得到發(fā)動機啟動時各個考察時間點上，半徑方向上各節(jié)點溫度值T關(guān)于圓筒半徑r的函數(shù)關(guān)系，再代入圓筒熱應(yīng)力計算公式（6）、（7）中，利用MATLAB編寫程序，計算得到徑向熱應(yīng)力和周向熱應(yīng)力。

考慮了第一類和第三類邊界條件，即規(guī)定圓筒內(nèi)側(cè)溫度隨時間變化，圓筒外側(cè)與外界干空氣表面?zhèn)鳠嵯禂?shù)h以及周圍空氣的溫度。圓筒內(nèi)部高溫燃氣充當(dāng)內(nèi)熱源。但是，在固體發(fā)動機實際工作過程中，熱傳遞將變得更為復(fù)雜，為了更好的對固體發(fā)動機筒壁進行溫度場和熱應(yīng)力測定分析，還應(yīng)當(dāng)考慮輻射換熱以及不同的升溫速率對溫度場和熱應(yīng)力的影響。因為篇幅限制，本論文也只考慮了啟動工作狀態(tài)下的非穩(wěn)態(tài)導(dǎo)熱傳熱過程，更好的改進方法是對固體發(fā)動機停止過程的非穩(wěn)態(tài)導(dǎo)熱過程進行熱應(yīng)力分析。同時可以在本文一維非穩(wěn)態(tài)基礎(chǔ)上進行多維非穩(wěn)態(tài)導(dǎo)熱及熱應(yīng)力的分析計算，由此才能夠更接近實際工況，獲得更為精確的瞬態(tài)溫度以及瞬態(tài)熱應(yīng)力，同時也能夠更好的論證某種材料在固體發(fā)動機上運用的可行性，以達到固體火箭發(fā)動機地面試車的實驗熱應(yīng)力計算分析的要求。

7 結(jié)論

（1）以固體發(fā)動機地面試車啟動工況為背景，將固體發(fā)動機噴管簡化為等厚壁圓筒，給出了一維非穩(wěn)態(tài)工作條件下瞬態(tài)溫度分布以及熱應(yīng)力的分析方法的基本思路。

（2）利用邊界條件以及差分法離散導(dǎo)熱偏微分方程求解，研究了等厚壁圓筒一維非穩(wěn)態(tài)導(dǎo)熱的數(shù)值解法，得到了一維非穩(wěn)態(tài)導(dǎo)熱方程離散形式的顯式數(shù)值解。

（3）利用了ABAQUS對瞬態(tài)導(dǎo)熱過程進行模擬，得到了相應(yīng)的瞬態(tài)溫度場云圖，與數(shù)值溫度計算結(jié)果進行了對比，以確保數(shù)值解的準確性。利用MATLAB編制程序進行熱應(yīng)力的求解，得到了各節(jié)點不同時刻熱應(yīng)力分布，提高了計算效率。

（4）由于固體發(fā)動機的工作溫度通常是超高溫的，溫度變化范圍大，故本文考慮了鈦鎳基合金的物性隨溫度的線性變化，同時考慮了實際工況下的對流換熱邊界條件，從而計算過程更接近于實際情況。

（5）本文推導(dǎo)出的溫度數(shù)值解以及熱應(yīng)力計算相對較為簡單，但是為固體發(fā)動機概念設(shè)計和參數(shù)敏感性分析提供方法支持，對實際工程應(yīng)用有一定的參考意義。

（6）本文基于鈦鎳基合金的熱物性，加深了對該合金的認識，推動了鈦鎳合金在固體火箭發(fā)動機領(lǐng)域的運用。

[1]曹繼敏，楊冠軍，鄧炬.鈦鎳合金的研究開發(fā)及市場狀況[J].鈦工業(yè)發(fā)展，2002（2）.

[2]Paton，Koppel.Thermalstressof thickcylinder[J].Machinedesign.1983（25）.

[3]平修二.熱應(yīng)力與熱疲勞[M].北京：國防工業(yè)出版社，1984.

[4]丸山武志，十河宏，麥春生.多層圓筒的熱傳導(dǎo)[J].化工與通用機械，1977(6).

[5]孟祥貴，陳棣湘，潘孟春.厚壁圓筒熱應(yīng)力問題的傳遞函數(shù)方法[J].強度與環(huán)境，2008（35）.

[6]蘇俠杰.無限長FGM圓筒穩(wěn)態(tài)熱應(yīng)力的有限元研究[J].河北：河北工程大學(xué)，2008.

[7]楊能彪.一維非穩(wěn)態(tài)導(dǎo)熱問題的數(shù)值計算[J].青海師范大學(xué)學(xué)報.2006（8）.

[8]盧學(xué)長.多層圓筒非定常穩(wěn)態(tài)傳熱的計算[J].遼寧：東北大學(xué)，2008.

[9]喬春珍，吳照云，項新耀.一維非穩(wěn)態(tài)導(dǎo)熱過程熱傳遞的規(guī)律及計算[J].熱科學(xué)與技術(shù)，2000（3）.

[10]封恒林.高溫下厚壁長圓筒熱應(yīng)力的有限元分析(一)[J].橡塑技術(shù)與裝備，1985（2）.