防雷/防熱涂層的雷擊熱解損傷的跨尺度模擬

董 琪，趙致遠(yuǎn)，李 婷，賈玉璽，張如良

（1.山東科技大學(xué)，青島，266590；2.航天材料及工藝研究所，北京，100076；3.山東大學(xué)，濟(jì)南，250061）

0 引言

運(yùn)載火箭、導(dǎo)彈等航天器飛行將面臨超高溫?zé)崃鳝h(huán)境，其鋁合金或碳纖維復(fù)合材料殼體不能直接裸露在外，需涂覆硅橡膠類燒蝕防熱涂層［1］，以阻止高溫破壞殼體并向內(nèi)傳遞。而目前防熱涂層通常采用絕緣材料，在雷電環(huán)境中飛行會(huì)面臨雷擊破壞，致使航天器在后續(xù)飛行中失去防熱保護(hù)，導(dǎo)致任務(wù)失敗。因此，開(kāi)展防熱涂層雷擊損傷研究十分必要。

雷電作用下防雷/防熱涂層產(chǎn)生的直接效應(yīng)具有顯著的電-熱-化學(xué)-力多物理場(chǎng)耦合特征，涉及涂層材料的絕緣擊穿、電阻生熱、電弧燒蝕、樹(shù)脂熱解及陶瓷化以及力學(xué)破壞等。然而，雷擊試驗(yàn)的成本高、危險(xiǎn)性高、復(fù)雜瞬態(tài)過(guò)程難以表征，因此，使用計(jì)算機(jī)進(jìn)行模擬計(jì)算成為目前雷擊損傷研究的重要手段。其中，以有限元模擬為代表，能夠考慮雷擊損傷復(fù)雜的電-熱-化學(xué)-力多物理場(chǎng)耦合特征，量化不同耦合場(chǎng)作用下不同類型的損傷情況［2］，可以有效闡明防雷/防熱涂層雷擊損傷的機(jī)理，為雷擊損傷預(yù)測(cè)提供理論基礎(chǔ)。Dong 等［3］建立了碳纖維增強(qiáng)復(fù)合材料（Carbon Fiber Reinforced Plastic，CFRP）雷擊損傷過(guò)程的電-熱-化學(xué)分析模型并加以模擬分析；Alexandros 等［4］建立了CFRP 層合板雷擊損傷過(guò)程的電-熱耦合模型，對(duì)比分析了層合板由螺栓連接和粘貼連接時(shí)的雷擊損傷情況；Zhu 等［5］則利用CFRP 復(fù)合材料雷電損傷有限元模型討論了其幾何形狀、纖維排布順序等因素對(duì)損傷程度的影響。

但以有限元為代表的多物理場(chǎng)耦合技術(shù)由于時(shí)間尺度的限制，在模擬涂層在雷擊熱效應(yīng)下發(fā)生熱解這一過(guò)程時(shí)尚存在不足，采用的電-熱-化學(xué)耦合有限元模型計(jì)算熱解損傷時(shí)缺乏可靠的熱解反應(yīng)動(dòng)力學(xué)方程。這是由于熱解反應(yīng)動(dòng)力學(xué)方程的獲取通常來(lái)源于對(duì)材料熱重試驗(yàn)曲線的熱分析，而熱重試驗(yàn)條件遠(yuǎn)達(dá)不到實(shí)際雷擊所造成的極端升溫速率和極端高溫情況，所得的動(dòng)力學(xué)方程不能很好地與實(shí)際吻合。因此，為更合理地研究防雷/防熱涂層雷擊損傷機(jī)理，進(jìn)而準(zhǔn)確預(yù)測(cè)其雷擊燒蝕損傷，探明涂層材料在雷擊極端溫度條件下的熱解動(dòng)力學(xué)機(jī)制及行為具有重要意義。

目前針對(duì)此類問(wèn)題的研究較為少見(jiàn)，然而，針對(duì)這一問(wèn)題，分子模擬手段能夠從原子-分子的微觀尺度模擬聚合物分子在極端升溫速率和極端高溫下的熱解行為。Zhang等［6］利用分子模擬方法對(duì)PVC和巴塞木進(jìn)行了熱解模擬，并分析了模擬結(jié)果與雷擊損傷的關(guān)系；Paajanen等［7］模擬了高溫下纖維素的熱解并分析了其反應(yīng)機(jī)理與動(dòng)力學(xué)。

本文采用反應(yīng)分子動(dòng)力學(xué)模擬方法研究甲基苯基硅橡膠在雷電極端升溫速率下的熱解反應(yīng)機(jī)制，進(jìn)而通過(guò)雷擊過(guò)程電-熱-化學(xué)耦合有限元模擬，跨尺度計(jì)算甲基苯基硅橡膠基防雷/防熱涂層材料體系的雷擊熱解損傷，并將模擬結(jié)果與雷擊試驗(yàn)損傷結(jié)果作對(duì)比，旨在構(gòu)建能更準(zhǔn)確反映涂層雷擊實(shí)際過(guò)程的跨尺度模擬方法，為涂層雷擊熱解損傷預(yù)測(cè)以及防熱層防雷材料開(kāi)發(fā)提供理論支持。

1 基礎(chǔ)理論與模擬方法

1.1 涂層樹(shù)脂熱解反應(yīng)分子動(dòng)力學(xué)模擬

選取甲基苯基硅橡膠的交聯(lián)網(wǎng)絡(luò)狀結(jié)構(gòu)作為研究體系，其分子的通用結(jié)構(gòu)式如圖1所示。

苯基為硅橡膠提供卓越的耐熱、耐低溫及耐輻照性能，本文以30%苯基含量構(gòu)建甲基苯基硅橡膠分子鏈模型作為交聯(lián)前的單體模型，如圖2所示。

圖2 甲基苯基硅橡膠的單體分子模型Fig.2 Molecular model of methyl phenyl silicone rubber monomer

將上述16 個(gè)硅橡膠單體分子鏈裝入大小為30.4 ?×30.4 ?×30.4 ? 的具有用周期性邊界條件的立方盒子中，構(gòu)建交聯(lián)體系，初始建模密度為1 g/cm3。根據(jù)氧化脫氫的硫化反應(yīng)法，使用交聯(lián)腳本程序［8］代替催化劑的作用以簡(jiǎn)化交聯(lián)過(guò)程，將硅橡膠鏈上甲基的C 原子標(biāo)記為反應(yīng)位點(diǎn)直接進(jìn)行交聯(lián)，形成C—C 鍵的同時(shí)刪除甲基上的一個(gè)氫原子。模擬交聯(lián)過(guò)程采用COMPASS Ⅱ力場(chǎng)，目標(biāo)交聯(lián)度為100%。為增加分子的移動(dòng)性，設(shè)置溫度為500 K，最終得到交聯(lián)度為95.1%的硅橡膠交聯(lián)結(jié)構(gòu)模型，如圖3所示。

圖3 甲基苯基硅橡膠交聯(lián)分子模型Fig.3 Molecular model of cross-linked methyl phenyl silicone rubber

基于上述硅橡膠分子模型，本文通過(guò)Lammps 2020 軟件，采用Adri van Duin 等［9］開(kāi)發(fā)的ReaxFF 反應(yīng)力場(chǎng)模擬硅橡膠的熱解過(guò)程。

在實(shí)際的雷擊過(guò)程中，復(fù)合材料雷擊作用點(diǎn)處升溫速率可達(dá)1011K/min以上，極短時(shí)間內(nèi)可達(dá)3 000 K以上的極高溫度［10］。為使熱解反應(yīng)模擬的條件盡可能接近雷擊的實(shí)際情況，同時(shí)兼顧計(jì)算成本和計(jì)算機(jī)運(yùn)行效能，本文進(jìn)行了升溫速率位于1012～1014K/s 區(qū)間范圍內(nèi)的多次熱解反應(yīng)分子模擬，且將模擬的終止溫度設(shè)置在7 000 K以上，以保證分子模型完全熱解，具體算例的溫度設(shè)置條件如表1所示。此外，模擬在NPT系綜下進(jìn)行，設(shè)置步長(zhǎng)為0.25 fs。

表1 不同極端升溫速率下硅橡膠熱解的分子模擬條件Tab.1 Molecular simulation conditions for the pyrolysis of silicone rubber under different extreme heating rates

根據(jù)數(shù)據(jù)處理需要，本文設(shè)置每100步輸出一次當(dāng)前時(shí)刻模擬系統(tǒng)內(nèi)產(chǎn)物信息，該信息包括當(dāng)前時(shí)刻下系統(tǒng)中所有存在的產(chǎn)物的分子式及其相應(yīng)數(shù)量。

1.2 樹(shù)脂熱解反應(yīng)分子模擬的動(dòng)力學(xué)方程

基于上述硅橡膠熱解反應(yīng)分子模擬輸出的不同時(shí)刻產(chǎn)物信息，本文根據(jù)式（1）計(jì)算不同時(shí)刻下固相產(chǎn)物所占的質(zhì)量比：

式中W為質(zhì)量保留率；x為當(dāng)前體系中存在的化合物種類總數(shù)量；N0為反應(yīng)開(kāi)始時(shí)體系的微粒數(shù)；M0為反應(yīng)開(kāi)始時(shí)體系的相對(duì)分子質(zhì)量；Nti為t時(shí)刻下第i種剩余固相產(chǎn)物的微粒數(shù)；Mti為t時(shí)刻下第i種剩余固相產(chǎn)物的相對(duì)分子質(zhì)量。

由此可得到硅橡膠熱解反應(yīng)的虛擬熱解度曲線。通過(guò)不同升溫速率下的虛擬熱解度曲線，進(jìn)一步計(jì)算動(dòng)力學(xué)參數(shù)（活化能Ea、指前因子A、反應(yīng)級(jí)數(shù)n）。本文采用的熱解動(dòng)力學(xué)參數(shù)計(jì)算方法如下：

a）Starink法。

Starink［11］對(duì)比了熱分析動(dòng)力學(xué)方法中經(jīng)典的Kissinger 方程、Ozawa 方程和Boswell 方程的計(jì)算精度，并提出了較高計(jì)算精度的Starink方程，其具體形式如下：

式中β為升溫速率；T為熱力學(xué)溫度；Rg為普適氣體常量；G(α)為反應(yīng)機(jī)理函數(shù)的積分式。

根據(jù)該方程進(jìn)行擬合可計(jì)算出活化能Ea，且在確定機(jī)理函數(shù)的前提下，可進(jìn)一步得到指前因子A。

b）Malek法。

硅橡膠的熱解過(guò)程是聚合物受熱解聚、無(wú)規(guī)則斷鏈的結(jié)構(gòu)轉(zhuǎn)變，是分子量變小的反應(yīng)，符合n級(jí)化學(xué)反應(yīng)機(jī)理，Malek法［12］能夠在得到活化能Ea的基礎(chǔ)上計(jì)算反應(yīng)級(jí)數(shù)n。對(duì)于n級(jí)化學(xué)反應(yīng)機(jī)理，Malek法計(jì)算反應(yīng)級(jí)數(shù)n的具體形式如下：

式中αp為熱解速率最大時(shí)的熱解度；up表示Ea/RgTp；π(up)為L(zhǎng)uke近似式［13］中引入的有理函數(shù)。

將αp、Tp、Ea、Rg代入式（3）中迭代可求出反應(yīng)級(jí)數(shù)n。

在得到活化能Ea、指前因子A、反應(yīng)級(jí)數(shù)n后，代入式（4）中，即得到熱解反應(yīng)分子模擬動(dòng)力學(xué)方程。

1.3 防雷/防熱涂層雷擊熱解損傷的有限元模擬

首先，構(gòu)建與雷擊試驗(yàn)中試樣一致的涂層平板幾何模型，如圖4所示。模型由三層材料構(gòu)成，表層為130 μm的防雷層、次表層為4 mm的防熱層、底層為3 mm 的鋁基板。為提高計(jì)算效率，構(gòu)建的幾何模型為1/4的平板模型，平面尺寸為250 mm×250 mm。

圖4 雷擊損傷電-熱-化學(xué)耦合有限元分析1/4涂層平板幾何模型Fig.4 Quarter geometric model for electrical-thermal-chemical coupling analysis of lightning damage of the coating

a）有限元邊界條件與載荷。

在模型的兩個(gè)非對(duì)稱側(cè)面設(shè)置零電勢(shì)，并在模型上表面施加對(duì)流和輻射熱邊界；將雷電弧的電流密度和熱流密度定義為面載荷，施加在模型上表面。雷電流密度分布方程定義如下：

式中J為電流密度；r為徑向坐標(biāo)；c為常數(shù)；Q為熱流密度；R(t)為雷電通道半徑，見(jiàn)式（7）。

式中Ip為峰值電流強(qiáng)度；t為時(shí)間。

電流載荷遵循人工模擬雷電試驗(yàn)SAE APR5412［14］標(biāo)準(zhǔn)，該標(biāo)準(zhǔn)依據(jù)自然雷電現(xiàn)象中的電流變化定義了用于試驗(yàn)研究的電流波形。為使模擬中加載的雷電流波形與試驗(yàn)電流波形保持一致，提取文獻(xiàn)［1］中雷擊試驗(yàn)所用雷電流波形曲線上測(cè)得的數(shù)據(jù)點(diǎn)進(jìn)行擬合，見(jiàn)圖5。

圖5 雷電流波形Fig.5 Lightning current waveform

模型施加電流載荷分兩個(gè)階段，第1 階段施加500 μs 的雷電分量A，第2 階段不施加雷電流，散熱至60 s。其中，第1階段電流載荷為圖5b中擬合得到的紅色曲線，其具體函數(shù)形式如下：

b）控制方程。

涂層中雷電流傳導(dǎo)過(guò)程由如下方程控制：

式中V為控制單元面積；S為控制單元表面；rc為單位體積內(nèi)體積電流源；m表示S的法向量；σ為電導(dǎo)率矩陣；φ為電勢(shì)。

雷電流在涂層中傳導(dǎo)過(guò)程中電能會(huì)以焦耳熱的形式部分轉(zhuǎn)換為熱能，相應(yīng)的溫度場(chǎng)控制方程為

式中ρ為密度；cp為比熱；k為熱導(dǎo)率矩陣；θ為內(nèi)熱源。

其中，內(nèi)熱源又分為焦耳熱θe和熱解反應(yīng)熱Qp，兩者分別由式（11）與式（12）表示。

式中η為能量轉(zhuǎn)換系數(shù)；Hs為聚合物熱解反應(yīng)的焓。

式（12）的熱解反應(yīng)熱可由1.2 節(jié)中熱解反應(yīng)分子模擬動(dòng)力學(xué)方程式（4）中定義的熱解速率代入計(jì)算得到。

上述涂層樹(shù)脂熱解反應(yīng)在有限元熱-化學(xué)反應(yīng)模塊的定義是通過(guò)用戶自定義子程序HETVAL 和USDFLD實(shí)現(xiàn)的。

本文利用子程序USDFLD計(jì)算并定義樹(shù)脂材料的熱解反應(yīng)和熱解程度，從而控制雷擊過(guò)程中材料的電學(xué)、熱學(xué)參數(shù)隨溫度和熱分解反應(yīng)程度的演變；另一方面利用子程序HETVAL定義了材料的升華吸熱、熱分解吸熱等熱通量。在仿真計(jì)算的每一個(gè)增量步開(kāi)始時(shí)，Abaqus求解器會(huì)將當(dāng)前的溫度、時(shí)間和狀態(tài)變量等數(shù)據(jù)傳遞給子程序USDFLD，隨后USDFLD根據(jù)這些數(shù)據(jù)和定義的算法等計(jì)算增量步結(jié)束時(shí)狀態(tài)變量的值，隨后，這些變量的值將被傳遞給Abaqus 求解器和HETVAL，由HETVAL計(jì)算當(dāng)前增量步所產(chǎn)生的熱通量，并將計(jì)算結(jié)果返回給Abaqus 求解器。Abaqus求解器首先根據(jù)USDFLD返回的狀態(tài)變量確定當(dāng)前的熱導(dǎo)率、電導(dǎo)率、比熱容等材料參數(shù)，而后再結(jié)合HETVAL返回的熱通量計(jì)算增量步結(jié)束時(shí)的溫度。如此不斷迭代，直至仿真結(jié)束。

c）材料性能參數(shù)。

防熱層、防雷層和鋁基板的材料屬性設(shè)置見(jiàn)表2至表4。其中，防雷層的電導(dǎo)率隨溫度變化的規(guī)律如圖6所示。

表2 防熱層材料性能參數(shù)Tab.2 Material performance parameters of the thermal insulation layer

圖6 防雷層電導(dǎo)率隨溫度變化的曲線Fig.6 Curve of electrical conductivity of the lightning protection layer as a function of temperature

表2中防熱層材料性能參數(shù)取自參考文獻(xiàn)［15］。由于防雷層是在與防熱層相同的樹(shù)脂基體中摻雜了金屬銅粉制備得到，因此表3中給出的防雷層材料性能是基于表2中硅橡膠的性能參數(shù)，根據(jù)銅、硅橡膠的比例，按照混合法則計(jì)算得到。圖6中防雷層電導(dǎo)率隨溫度變化的數(shù)據(jù)為試驗(yàn)測(cè)試得到。由于超過(guò)300 ℃后，防雷層材料的樹(shù)脂基體將開(kāi)始發(fā)生熱分解反應(yīng)，繼續(xù)升溫過(guò)程（高于300 ℃）中的電導(dǎo)率測(cè)試難以進(jìn)行，而在實(shí)際仿真中，雷擊區(qū)域的溫度將達(dá)到上千攝氏度，電導(dǎo)率的實(shí)測(cè)數(shù)據(jù)顯然不滿足仿真要求。因此，考慮到防雷層中樹(shù)脂基體的熱分解會(huì)使其中導(dǎo)電金屬粒子接觸電阻降低，本文通過(guò)唯象的方法，在防雷材料溫度超過(guò)300 ℃后，將其電導(dǎo)率設(shè)置為熱解度的線性函數(shù)，假設(shè)防雷材料完全熱解時(shí)的電導(dǎo)率為純銅的電導(dǎo)率（5.714×107S/m）；而在達(dá)到銅的沸點(diǎn)（2 562 ℃）時(shí)，防雷材料不再導(dǎo)電，電導(dǎo)率退化至1×10-14S/m。

表3 防雷層材料性能參數(shù)Tab.3 Material performance parameters of the lightning protective layer

由于雷電流并未擊穿防熱層到達(dá)鋁基板，即鋁基板并不會(huì)參與導(dǎo)電，所以其溫度變化不明顯，表4中鋁基板的材料參數(shù)無(wú)需考慮溫度的影響。

表4 鋁基板材料性能參數(shù)Tab.4 Material performance parameters of the aluminum substrate

2 結(jié)果與分析

2.1 涂層樹(shù)脂熱解反應(yīng)分子模擬的動(dòng)力學(xué)參數(shù)計(jì)算

涂層樹(shù)脂熱解過(guò)程分子模擬結(jié)果可視化處理如圖7所示，可以發(fā)現(xiàn)，在溫度上升的過(guò)程中交聯(lián)甲基苯基硅橡膠分子不斷解聚，熱分解產(chǎn)生的小分子體現(xiàn)出氣體擴(kuò)散現(xiàn)象，系統(tǒng)體積也隨之?dāng)U大。

圖7 1×1013 K/s升溫速率下硅橡膠熱解反應(yīng)分子模擬可視化結(jié)果Fig.7 Visualization results of molecular simulation of silicone rubber pyrolysis reaction at a heating rate of 1×1013 K/s

基于硅橡膠熱解過(guò)程產(chǎn)物數(shù)據(jù)，根據(jù)式（1）計(jì)算表1中所列算例的剩余固相產(chǎn)物，進(jìn)而得到硅橡膠在不同升溫速率下的熱解度隨溫度的變化曲線，部分升溫速率下所得的虛擬熱解度曲線如圖8所示。

圖8 在不同極端升溫速率下硅橡膠熱解反應(yīng)分子模擬所得虛擬熱解度曲線Fig.8 Virtual pyrolysis degree curves obtained from molecular pyrolysis simulations of silicone rubber at different extreme heating rates

取圖8 中每條熱解度曲線上α為0.1、0.2、…、0.9時(shí)的數(shù)據(jù)，依據(jù)式（2）的Starink方程，分別計(jì)算ln（β/T1.92）和-1.000 8/RgT，繪制散點(diǎn)圖并擬合，如圖9 所示。圖9 中每組擬合直線的斜率為對(duì)應(yīng)不同熱解度下的活化能Ea，其數(shù)值與擬合相關(guān)系數(shù)Rc如表5所示。

表5 Starink法所得硅橡膠熱解反應(yīng)分子模擬對(duì)應(yīng)不同熱解度下的活化能Tab.5 Activation energies at different pyrolysis degrees obtained from molecular pyrolysis simulation of silicone rubber by Starink method

圖9 Starink法所得硅橡膠熱解分子模擬對(duì)應(yīng)不同熱解度的活化能擬合結(jié)果Fig.9 Fitting results of activation energies for different pyrolysis degrees obtained from molecular pyrolysis simulation of silicone rubber by Starink method

根據(jù)式（3）的Malek法，將表5計(jì)算得到的活化能Ea分別代入不同升溫速率下的αp，Tp，Ea，Rg，迭代計(jì)算得到反應(yīng)級(jí)數(shù)n如表6 所示，取算數(shù)平均值，得到的反應(yīng)級(jí)數(shù)n為1.041。

表6 Malek法所得硅橡膠熱解反應(yīng)分子模擬的反應(yīng)級(jí)數(shù)nTab.6 Reaction order n obtained from the molecular pyrolysis simulation of silicone rubber by Malek method

將活化能Ea、n級(jí)化學(xué)反應(yīng)機(jī)理函數(shù)的積分式G(α)代入式（2）可得到指前因子A。其中，n級(jí)化學(xué)反應(yīng)機(jī)理函數(shù)的積分式G(α)代入反應(yīng)級(jí)數(shù)n后形式如下：

將該式與活化能Ea代入式（2），得到指前因子A如表7所示。

表7 Starink法所得硅橡膠熱解反應(yīng)分子模擬對(duì)應(yīng)不同熱解度的指前因子ATab.7 Pre-exponential factors A at different pyrolysisi degrees obtained from molecular pyrolysis simulations of silicone rubber by Starink method

將活化能Ea、指前因子A、反應(yīng)級(jí)數(shù)n代入式（4），就得到了涂層熱解分子模擬的反應(yīng)動(dòng)力學(xué)方程。

2.2 有限元模擬結(jié)果及試驗(yàn)驗(yàn)證

根據(jù)上述分子模擬計(jì)算得到的反應(yīng)動(dòng)力學(xué)方程編寫(xiě)用戶自定義子程序HETVAL 和USDFLD，在Abaqus軟件中模擬防雷/防熱涂層的雷擊熱解損傷。

圖10 為雷電流分量A 結(jié)束時(shí)刻和雷擊結(jié)束后自然冷卻60 s時(shí)的溫度分布云圖。結(jié)果表明，雷擊結(jié)束時(shí)，雷擊中心處溫度接近2 000 ℃，而散熱60 s 后導(dǎo)電層與防熱層溫度均降至30 ℃左右，遠(yuǎn)低于硅橡膠的熱解溫度，因此認(rèn)為雷擊結(jié)束后60 s時(shí)防護(hù)涂層的熱解損傷過(guò)程已經(jīng)結(jié)束。

圖10 標(biāo)準(zhǔn)雷電波形A作用下1/4涂層平板模型溫度分布Fig.10 Temperature distribution of the quarter flat model under the standard lightning waveform A

圖11為標(biāo)準(zhǔn)雷電波形A作用下涂層結(jié)構(gòu)在不同時(shí)刻的熱解云圖。結(jié)果顯示，標(biāo)準(zhǔn)雷電波形A作用下涂層的熱解損傷形貌呈圓形。表面導(dǎo)電防護(hù)層熱解損傷面積在雷擊作用開(kāi)始的前100 μs 快速擴(kuò)大，100～500 μs變化速率減慢，而在雷擊作用結(jié)束后的散熱過(guò)程中幾乎不再擴(kuò)大。在厚度方向上，熱解損傷在雷擊作用階段及其結(jié)束后的前1 ms內(nèi)主要集中在防雷層，此后才開(kāi)始向防熱層擴(kuò)展。熱解損傷深度在雷擊結(jié)束后約5 ms時(shí)達(dá)到最大值。

圖11 標(biāo)準(zhǔn)雷電波形A作用下1/4涂層平板模型在不同時(shí)刻下的熱解度云圖Fig.11 Pyrolysis degree contours of quarter plat model at different times under the standard lightning waveform A

這是由于防雷層的熱解損傷主要由雷電流在材料中傳導(dǎo)產(chǎn)生的焦耳熱導(dǎo)致，雷擊過(guò)程的前100 μs雷電流較大，在材料中產(chǎn)生了大量的焦耳熱，導(dǎo)致熱解區(qū)域快速擴(kuò)張。在100 μs 后，隨著雷電流強(qiáng)度的減弱，熱解損傷的擴(kuò)展速率降低。而防熱層不導(dǎo)電，其升溫?zé)嵩礊榉罒釋优c防雷層之間的熱傳遞，傳遞效率受兩者間溫差及材料導(dǎo)熱效率的影響。而由于防熱層具有較好的隔熱性能，因此需要其與防雷層溫差較大時(shí)溫度才會(huì)顯著升高并導(dǎo)致材料的熱解損傷，具有一定的滯后性。

為了驗(yàn)證仿真模型與結(jié)果的合理性，對(duì)比雷擊熱解損傷有限元模擬結(jié)果與硅橡膠防護(hù)涂層人工雷擊試驗(yàn)的損傷形貌，如圖12 所示。由圖12 可以發(fā)現(xiàn)，試驗(yàn)結(jié)果與雷擊熱解損傷模擬結(jié)果均體現(xiàn)出圓形損傷形貌，但試驗(yàn)得到的損傷形貌由靠近中心的完全熱解區(qū)和外側(cè)的放射狀未完全熱解區(qū)組成。考慮到試驗(yàn)中試樣的接地狀態(tài)，雷電流的傳導(dǎo)難以形成未完全熱解區(qū)的放射狀損傷形貌，因此，推測(cè)放射狀的損傷由雷電電弧掃略并燒蝕形成。而本文所建立的有限元模型僅考慮了雷電流在模型中傳導(dǎo)產(chǎn)生的焦耳熱及其導(dǎo)致的燒蝕損傷，暫未考慮電弧掃略注入的電弧熱等，因此，僅取試驗(yàn)結(jié)果中的完全熱解區(qū)與仿真結(jié)果進(jìn)行對(duì)比。試驗(yàn)結(jié)果中完全熱解區(qū)的直徑約為220 mm，仿真結(jié)果得到的熱解區(qū)域直徑約為195 mm，兩者相對(duì)誤差約為11.4%。以上結(jié)果表明，仿真結(jié)果與試驗(yàn)結(jié)果較為吻合，所建立的有限元仿真模型較為可靠。

圖12 防雷/防熱涂層雷擊損傷模擬與試驗(yàn)結(jié)果對(duì)比Fig.12 Comparison of lightning damage simulation and experimental results for lightning protection/thermal insulation coatings

3 結(jié)束語(yǔ)

針對(duì)甲基苯基硅橡膠基防雷/防熱涂層，采用反應(yīng)分子動(dòng)力學(xué)模擬與有限元模擬兩種手段，建立了在雷擊極端溫度環(huán)境下涂層熱解損傷的跨尺度模擬方法。模擬所得的涂層完全熱解形貌及面積與人工模擬雷擊試驗(yàn)結(jié)果吻合較好。但由于有限元模型未考慮電弧掃略造成的燒蝕作用，仿真結(jié)果未能體現(xiàn)試驗(yàn)結(jié)果中涂層未完全熱解區(qū)的放射狀形貌。綜上，文中構(gòu)建的跨尺度模擬方法能較可靠地預(yù)測(cè)硅橡膠基防雷/防熱涂層的雷擊熱解損傷，具有較好的工程應(yīng)用價(jià)值。