直升機工作熱載荷條件下耐熱功能材料變熱物性模型

黃 波，高志武，李 黔

(中國直升機設計研究所，江西 景德鎮 333001)

0 引言

在耐熱功能材料燒蝕模型研究方面，國內外研究學者針對耐熱材料的配方、燒蝕機理、燒蝕模型等展開了深入研究。Resch針對硅橡膠和氟橡膠兩種隔熱材料創建了考慮受熱化學反應的熱解燒蝕模型，并通過試驗驗證了模型。Blackwell等提出了一種利用朗道變換解決變面積熱防護材料燒蝕問題的方法，該方法基于應用積分形式的守恒方程來控制體積，求解控制體的熱傳導和對流通量。Chen Y等開發了用于預測熱防護材料炭化燒蝕和形狀變化的計算程序，其控制方程包括能量守恒和三層熱解層模型，程序能夠模擬大長徑比超音速飛行器的燒蝕。張杰等運用熱解動力學建立耐熱材料的一維燒蝕模型，通過有限體積法離散計算，分析反映熱解材料的壁面燒蝕率隨壁面溫度變化，材料因熱解的內部空隙率為梯度分布。徐義華等結合多孔介質傳質傳熱理論針對具有多孔介質特性的三元乙丙橡膠(EPDM，ethylene propylene diene monomer)耐熱材料炭化層建立了以炭化層空隙率為參數的力學模型，參考所測定的炭化層結構確定力學強度系數，進而導出破壞準則，為耐熱材料的燒蝕數值仿真提供了相關力學參數。Maurizio Natali等重點考慮了EPDM的可逆及不可逆膨脹過程和熱解氣體經過多孔炭化層的焦化粘結致密過程，建立了一維燒蝕模型。該模型的建立過程為熱防護材料的能量傳遞研究提供了良好的理論基礎。劉驍等基于化學反應動力學模型對炭化材料熱解過程進行描述，使用彈簧松弛法對耐熱材料的表面燒蝕后退展開模擬，且考慮材料熱解、熱解氣體與材料的熱交換和熱解氣體流動。耐熱材料三維溫度場的模擬結果與試驗數據對比驗證了計算程序的準確性。Yang Liu等通過對不同燒蝕環境下獲取的炭化層的特征分析，開發了一種模擬絕緣材料整體燒蝕過程的數值計算方法，針對炭化層的多孔介質特性，建立了一種全面的燒蝕模型，其中著重考慮了耐熱材料熱化學燒蝕和顆粒侵蝕。

然而，以上針對耐熱功能材料建立的傳熱燒蝕模型大都基于常物性或分段常物性假設，即材料的密度、熱導率、比熱容等物性參數保持不變，忽略了材料熱物性參數在燒蝕過程中的變化，影響燒蝕模型仿真的準確性。實際上，耐熱材料的熱物性與外界能量作用時間和耐熱材料內部溫度分布均相關，因此，考慮耐熱材料的變熱物性能夠更加準確地描述耐熱材料的傳熱燒蝕過程。

1 變熱物性模型

1.1 熱解氣體變熱物性模型

芳綸/EPDM耐熱材料在熱解期間產生多種氣體，且熱解氣體在流動過程中會吸收和釋放能量，因此充分了解熱解氣體的組成和性質是建立準確熱解氣體模型的前提。根據Jia Xiao-long等人的研究，芳綸/EPDM耐熱材料熱解氣體的主要成分為：2-丁烯(CH)，2-戊烯(CH)，3-己烯(CH)，苯蒸氣(CH)，庚烷(CH)，2-甲基-1-己烯(CH)。表1給出了熱解氣體主要產物的質量分數分布。

表1 熱解氣體主要產物質量分數分布

a

T

+

a

T

+

a

T

(1)

(2)

式中：

R

為通用氣體常數；系數

a

-

a

為熱力學數據表(NIST-JANAF Thermochemical Tables，National Institute of Standards and Technology-Joint Army Navy Air Force Thermochemical Tables)中的擬合系數。熱解氣體主要組分的擬合系數可查表獲取。

(3)

式中：

n

為總摩爾數；

n

為第

i

種組分的摩爾數；

X

為組分

i

的摩爾百分比。對于熱解氣體中的多原子分子組分，根據A.Eucken熱導率修正公式，同時結合已得到的各組分標準摩爾定壓比熱隨溫度變化的函數可推導出主要熱解產物的熱導率

k

(W/m·K)隨溫度的變化關系：

(4)

式中：

M

為各組分的氣體摩爾質量；

Ψ

為碰撞直徑；

Ω

為倫納德·瓊斯參數，稱為碰撞積分。熱解氣體各組分的碰撞直徑和碰撞積分可以查詢倫納德·瓊斯勢參數表獲取。根據倫納德·瓊斯輸運系數計算模型，已知組分

i

的熱導率為

k

，則熱解氣體的熱導率

k

為：

(5)

其中：

(6)

式中：

X

和

X

分別為組分

i

和

j

的摩爾百分比；

M

和

M

分別為組分

i

和

j

的摩爾質量。

1.2 材料密度模型

根據Arrhenius方程，熱解過程中材料密度的控制方程為：

T

>

T

且

ρ

>

ρ

(7)

式中：

A

為指前因子；

E

為熱解反應的活化能；

T

為耐熱材料熱解反應的臨界溫度；

ρ

為炭化層密度。根據“多組分模型”定義，芳綸/EPDM耐熱材料的炭化層密度為：

ρ

=

γVρ

+(1-

V

)

ρ

(8)

式中：

γ

為EPDM橡膠的熱解成碳率；

V

為EPDM橡膠在原始材料中的體積分數；

ρ

、

ρ

分別為EPDM橡膠和芳綸纖維的密度。

對材料密度的控制方程在當前時間步內積分，結果如下：

(9)

式中：

ρ

+為材料在當前時間步內的當地密度；

ρ

為材料在上一時間步內的密度。

1.3 材料比熱容和熱導率模型

熱解過程中，熱解層內高聚物裂解釋放出小分子熱解氣體并剩余殘碳，故而此狀態下熱解層主要包含未裂解的EPDM橡膠基體、芳綸纖維、基體裂解的剩余殘碳、大量孔隙和彌散在孔隙中的小分子熱解氣體，所以該狀態下材料密度定義為：

ρ

=(1-

V

)

ρ

+(1-

m

)

Vρ

+

m

γVρ

(10)

式中：m為已裂解的EPDM橡膠基體占當地初始基體的質量分數。由于芳綸/EPDM耐熱材料的原始層密度

ρ

為：

ρ

=

Vρ

+(1-

V

)

ρ

(11)

則

(12)

即在原始層中，

m

=0；在熱解層中，0<

m

<1；在炭化層中，0<

m

<1。由此可知熱解過程中芳綸纖維、EPDM橡膠和熱解殘碳所占的質量分數分別為：

(13)

然后再根據“多組分模型”可得熱解過程中芳綸/EPDM耐熱材料的比定壓熱容(TPS，thermal protection system)為：

c

=

m

c

+

m

c

+

m

c

(14)

式中：

c

、

c

和

c

分別為芳綸纖維、EPDM橡膠和熱解殘碳的比定壓熱容。

在耐熱材料的傳熱模型中，炭化層是傳質傳熱的主要場所，直接關系到芳綸/EPDM耐熱材料的抗燒蝕性能，所以對炭化層的結構表征至關重要。炭化層的微觀結構為連續相固體骨架分布于所占據的體積空間，其間的孔隙空間多數相互聯通，且充滿單相或多相介質，故而孔隙率是多孔介質結構的重要表征參數之一。所以，結合“有效熱導率法”定義熱解過程中芳綸/EPDM耐熱材料的熱導率為：

k

=(1-

V

)

k

+(1-

m

)

Vk

+(1-

P

)

m

Vk

+

Pm

Vk

(15)

式中：

k

、

k

和

k

分別為芳綸纖維、EPDM橡膠和熱解殘碳的熱導率。

1.4 控制方程

本文假設耐熱功能材料的熱導率為各向同性，其內部的二維非定常能量傳遞控制方程的積分形式為:

(16)

式中：

U

=

ρc

T

，

ρ

、

c

、

k

分別為材料當前密度、比定壓熱容和熱導率，這三者均為時間

t

和材料當前溫度

T

的函數；

Ψ

為軸對稱源項的開關函數，

Ψ

=1時為二維軸對稱控制方程，

Ψ

=0時為二維控制方程；

Q

表示材料單位時間、單位體積內的總熱量，即：

(17)

(18)

式中：

H

為熱解氣體顯焓，是熱解氣體當前溫度

T

的函數；

ε

為材料表面發射率；

σ

為斯忒藩·玻爾茲曼常數；

T

為外界環境參考溫度；

H

為耐熱材料的熱解潛熱。

2 模型驗證

為驗證本文所提出的變熱物性燒蝕模型與計算方法，將分別對固體熱傳導程序和變熱物性燒蝕模型進行考核。采取具有已知解析解的半無限大平板熱傳導算例和二維熱傳導算例對程序計算結果比對分析，充分驗證本文所發展程序與數值計算方法的精度。

2.1 半無限大平板熱傳導算例驗證

建立的變熱物性燒蝕模型與計算方法屬于非定常導熱問題，為驗證該程序，選取了半無限大平板瞬態熱傳導過程進行計算。假設，有一半無限大平板，其初始狀態溫度為

T

，現對半無限大平板設定兩種邊界條件進行瞬態熱傳導驗證。工況1：在半無限大平板表面加載等熱流密度邊界條件

q

。工況2：在半無限大平板表面加載等溫壁邊界條件

T

。工況1和工況2中半無限大平板內部的溫度-時間的解析解分別為：

(19)

(20)

式中：

erf

()表示誤差函數，

α

=

k/

(

ρc

)為組成半無限大平板材料的熱擴散系數，計算中組成半無限大平板的固體材料熱導率取

k

=49.8W/(m·K)，密度取

ρ

=7840kg/m，比熱容取

c

=454J/(kg·K)。圖1和圖2分別對兩種工況下的理論值與計算值進行了對比。

圖1 等熱流密度邊界條件下半無限大平板內部溫度變化

圖2 等溫壁邊界條件下半無限大平板內部溫度變化

從兩個工況的對比圖中可得，數值計算的結果與解析解吻合度極高。半無限大平板的導熱問題屬于典型的瞬態熱傳導問題，計算值與理論值的高度吻合驗證了本章所建立的模型與計算方法的準確性。

2.2 二維熱傳導算例驗證

上一節通過半無限大平板熱傳導對模型中的熱傳導計算方法進行了驗證。究其本質，半無限大平板熱傳導其實為一維導熱問題。故而接下來采用二維矩形物體的熱傳導，進一步對模型計算二維導熱問題進行驗證。如圖3所示，一個二維矩形物體，其邊長分別為

L

=60mm，

L

=40mm，該物體熱物性參數定義為常數。邊界條件設置：左、右及下三個邊界均為等溫壁，

T

1=300K；上壁面也為等溫壁，

T

(

x

)=

T

1+

T

2sin(

πx/L

)，其中

T

2=600K。顯然，該二維熱傳導問題即為二維穩態熱傳導問題，控制方程如下：

圖3 二維矩形物體及邊界條件

(21)

已知二維矩形物體的尺寸及邊界條件，則可得該物體內部溫度分布的解析解為：

(22)

對整個二維矩形物體進行結構化網格劃分求解，網格分布為30×20。采用所提出的模型對該熱傳導問題計算，得到的二維矩形物體內部溫度分布如圖4所示。由圖4可見，該物體上壁面溫度最高，通過熱傳導溫度不斷向下推進，溫度由上至下呈現下降趨勢。左、右及下邊界受等溫壁條件限制，溫度維持在300K。

圖4 二維矩形物體內部溫度云圖

為便于與解析解進行對比，圖5直觀地給出了二維矩形物體內部橫縱中心截面(即

y

= 0.5

L

處截面和

x

= 0.5

L

處截面)的理論值溫度與計算值溫度對比。由圖5可見，理論值與計算值差距極小，吻合度極高，展示了本文所建立的模型和計算方法對二維熱傳導問題求解的準確性與可靠性。

圖5 二維矩形物體內部橫縱中心截面理論值與計算值溫度對比圖

3 直升機工作熱載荷條件下熱響應模擬

直升機工作過程中，排氣管持續噴出大量高溫尾氣。受下洗流作用，尾氣與附近直升機殼體主要通過熱傳導發生劇烈的熱交換現象，導致動力艙靠近排氣口位置處的溫度較高，如圖6所示。殼體溫度最高處達443K，且殼體較大面積溫度超414K。因此，直升機殼體選用耐熱功能材料，并對其工作過程中的溫度準確預測及監控尤為重要。

圖6 直升機排氣高溫區

復合材料因其耐熱氧老化性能良好、密度低、抗疲勞性強和層鋪可設計性佳等優異特性，作為直升機上的耐熱功能材料得到日益廣泛的應用。其中，復合型耐熱材料廣受青睞：直升機結構用的增強耐熱材料如芳綸、玻璃纖維和碳纖維及其織物，旋翼、尾槳常用的熱防護材料如短切纖維和預浸玻璃布，耐熱功能材料中的耐高溫抗燒蝕材料如超級纖維和碳纖類等等。

本章以芳綸/EPDM耐熱材料為例，運用變熱物性燒蝕模型，模擬直升機工作環境下耐熱材料的熱防護工作過程。芳綸/EPDM耐熱材料熱物性參數隨溫度和時間的數學模型，能夠保證耐熱材料原始層、熱解層和炭化層傳熱燒蝕控制方程的統一求解，不必對某一結構層次進行單獨求解，簡化了數學模型，提高了數值運算效率。本文通過熱平衡法保證能量守恒，采用中心差分法對控制方程進行離散，并將上述變熱物性傳熱模型和源項加入離散方程中，得到離散控制方程。計算所需的芳綸/EPDM耐熱材料相關參數如表2所示。

表2 芳綸/EPDM耐熱材料相關參數

依據上述所描述的傳熱模型和離散方程編寫計算程序，計算流程如圖7所示，圖中

t

表示計算時間，Δ

t

表示時間步長。通過比對試驗與仿真結果考核模型和計算方法的可靠性。

圖7 計算流程圖

以參考文獻[6]中的芳綸/EPDM耐熱材料燒蝕試驗對變熱物性數學模型和數值計算方法進行驗證。耐熱材料內部深度5mm和10mm處溫度的變化歷程如圖8所示，并與常熱物性計算結果進行了對比。

圖8 試驗結果與仿真結果對比

從5mm處溫度變化歷程可以看出，兩者相比，變熱物性模型計算結果與試驗值吻合較好，最大誤差僅為10.61K。常熱物性模型計算結果在前50s與變熱物性模型計算結果基本一致，但是在50s-90s之間，計算結果差異較為明顯，溫度最大值絕對誤差達到20.12K。

另外，從10mm處溫度變化歷程可以發現，變熱物性模型計算結果與試驗值吻合較好，最大誤差僅為4.22K。而10mm處的常熱物性模型計算結果在80s吻合較好，之后吻合較差，期間最大誤差達11.31K。綜上所述，變熱物性模型計算結果與試驗結果吻合較好，而常熱物性模型計算結果有一定差異。結果表明本文所發展的變熱物性模型以及數值計算方法對芳綸/EPDM耐熱材料熱解過程中能量傳遞的熱響應仿真有較高的準確性，對預測直升機上所用復合型耐熱功能材料的熱防護工作過程具有一定參考價值。

4 結論

針對耐熱功能材料傳熱燒蝕過程，提出了一種耐熱材料比熱容、熱導率等熱物性參數隨時間和溫度變化的數學模型，著重介紹、推導了熱解氣體變熱物性模型、材料密度模型和材料比熱容及熱導率模型的建立過程。

采取具有已知解析解的半無限大平板熱傳導算例和二維熱傳導算例對模型程序考核驗證，計算值與理論值高度吻合。通過與經典的燒蝕試驗結果和常熱物性計算結果進行對比，開發的變熱物性燒蝕程序計算結果與試驗值吻合較好，而常熱物性模型計算結果與試驗值有一定差異，表明變熱物性數學模型對耐熱材料的熱解過程仿真具有較高的準確性，充分展示了計算模型的可靠性和可信度，對預測直升機上所用復合型耐熱功能材料的熱防護工作過程具有一定參考價值。