基于COMSOL的熱壓罐成型模具溫度場(chǎng)分析

呂佳鎂，馬貴春，郭靖宇，岳 光，段連成，李云浩

(1.中北大學(xué) 機(jī)電工程學(xué)院，山西 太原 030051；2.莫斯科鮑曼國(guó)立技術(shù)大學(xué) 特種制造系，莫斯科 105005)

0 引 言

復(fù)合材料誕生于20世紀(jì)30年代，具有輕質(zhì)高強(qiáng)、耐腐蝕、絕緣等優(yōu)良性能，在我國(guó)的研究始于20世紀(jì)50年代末期，為滿足性能高、重量輕、成本低的新時(shí)代飛行條件的要求，復(fù)合材料成型技術(shù)被廣泛應(yīng)用于飛機(jī)零部件制造，并逐漸應(yīng)用于承力結(jié)構(gòu)，復(fù)合材料用量可作為衡量飛行器先進(jìn)性的重要指標(biāo)[1].同時(shí)，復(fù)合材料的廣泛應(yīng)用促進(jìn)了其成型技術(shù)的發(fā)展.其中，熱壓罐成型技術(shù)可用于大型復(fù)雜構(gòu)件成型，且成型質(zhì)量穩(wěn)定，是目前最為常用的成型方法[2].模具作為復(fù)合材料構(gòu)件成型的載體，其表面不同的溫度值決定了復(fù)合材料不同的熱收縮率與化學(xué)收縮率，并使構(gòu)件內(nèi)部積累應(yīng)力梯度，從而發(fā)生固化變形，導(dǎo)致構(gòu)件力學(xué)性能降低.可見，模具溫度場(chǎng)直接決定復(fù)合材料構(gòu)件的成型質(zhì)量.

通過(guò)實(shí)驗(yàn)研究大型復(fù)合材料構(gòu)件熱壓罐成型過(guò)程的成本較高，所以如何準(zhǔn)確、高效地模擬模具溫度場(chǎng)是當(dāng)下的研究熱點(diǎn)，工藝仿真模擬成為了重要支撐技術(shù).Park等[3]建立了一種通用的二維固化模擬模型，用于快速模擬任意形狀的復(fù)合材料結(jié)構(gòu).Hoon等[4]將復(fù)合材料固化的三維有限元公式編寫成相應(yīng)的程序用于模擬，所得結(jié)果與實(shí)驗(yàn)結(jié)果吻合.仿真軟件方面，COMPRO，ABAQUS等軟件常用于復(fù)合材料構(gòu)件制造過(guò)程中變形的模擬[5-6].李艷霞、徐強(qiáng)等[7-8]綜述了CFX，F(xiàn)LUENT等軟件在復(fù)材構(gòu)件溫度場(chǎng)方面的應(yīng)用，并對(duì)亟待解決的技術(shù)壁壘提出自己的思考.此外，研究中不斷引入新的仿真模擬軟件，如張晨群等[9]采用以格子玻爾茲曼法劃分網(wǎng)格的XFLOW軟件仿真模具溫度場(chǎng)，縮短了前處理時(shí)間.李彩林[10]基于PAM-AUTOCLAVE軟件平臺(tái)對(duì)成型工藝進(jìn)行模擬，分析了罐內(nèi)溫度場(chǎng)、流場(chǎng)的均勻性.

本文以用于某型飛機(jī)壁板的框架式成型模具為研究對(duì)象，利用COMSOL軟件模擬其溫度場(chǎng)分布情況，探究常見工藝參數(shù)對(duì)模具溫度分布的影響情況，并對(duì)未來(lái)模具溫度場(chǎng)的優(yōu)化提出幾點(diǎn)思考.

1 熱壓罐工作原理

如圖1 所示，風(fēng)扇驅(qū)動(dòng)罐內(nèi)氣體流經(jīng)熱源升溫后循環(huán)給模具加熱，再循環(huán)冷水吸收罐內(nèi)熱量從而降溫，其主要換熱方式為熱對(duì)流和熱傳導(dǎo).罐內(nèi)固化溫度變化如圖2 所示.

圖1 熱壓罐工作原理Fig.1 Working principle of autoclave

圖2 固化工藝溫度曲線Fig.2 Curing process temperature curve

2 數(shù)值模擬分析

2.1 COMSOL Multiphysics

目前計(jì)算流體力學(xué)(Computational Fluid Dynamics，CFD)軟件頗多，各有優(yōu)勢(shì).在處理高速、高精度和三維大空間等高難度問題時(shí)，對(duì)CFD軟件有很高的要求.但對(duì)于工程中常見的流體問題，主流的仿真軟件皆能保證計(jì)算精度.因此，在保證結(jié)果準(zhǔn)確性的前提下，運(yùn)用高效便捷、易于國(guó)人操作的工具解決實(shí)際問題是當(dāng)今研究不可忽視的需求.

本文選用功能靈活的COMSOL展開研究.COMSOL Multiphysics衍生于MATLAB中的Toolbox，它與MATLAB有完整接口，非常便于后續(xù)更為深入的研究.同時(shí)，其具備完全開放的構(gòu)架、友好的操作界面和語(yǔ)言環(huán)境，用戶可對(duì)多個(gè)物理場(chǎng)自由耦合，根據(jù)自己的需求定義偏微分方程，再采用自適應(yīng)劃分網(wǎng)格與求解器進(jìn)行準(zhǔn)確的計(jì)算、仿真.綜上，COMSOL具有良好的可擴(kuò)展性和靈活性.

2.2 控制方程

2.2.1 對(duì)流傳熱

根據(jù)連續(xù)介質(zhì)假設(shè)，罐內(nèi)所發(fā)生的流體流動(dòng)和傳熱現(xiàn)象均遵守牛頓力學(xué)中的三大守恒定律，即質(zhì)量、動(dòng)量、能量守恒，可將其作為仿真過(guò)程中的基本控制方程，在三維笛卡兒坐標(biāo)系下，方程的具體表現(xiàn)形式如下：

1)質(zhì)量守恒方程

(1)

2)動(dòng)量守恒(牛頓運(yùn)動(dòng)定律)方程

(2)

(3)

(4)

3)能量守恒(熱力學(xué)第一定律)方程

-ρdiv(U)+div(λ·grad(T))+SA+φ,

(5)

式中：p為流體壓力；U為流體流速；T為流體溫度；ρ為流體密度；η為動(dòng)力粘度；λ為流體導(dǎo)熱率；φ為流體耗散函數(shù)；h由流體壓強(qiáng)和流體溫度決定；Su,Sv,Sw為動(dòng)量方程的廣義源項(xiàng)；SA為流體內(nèi)熱源；u,v,w為速度矢量在x,y,z方向上的分量.

2.2.2 傳導(dǎo)傳熱

(6)

式中：ρS為Q235密度；CS為Q235比熱容；TS為Q235溫度；λS為Q235導(dǎo)熱率；QT為熱源項(xiàng).

2.2.3 熱輻射

通常情況下，在熱壓罐工作溫度范圍內(nèi)，罐壁與框架式模具間的熱輻射影響極弱，并且國(guó)內(nèi)外有研究證明是否考慮熱輻射對(duì)仿真結(jié)果影響不大，所以，本文不考慮熱輻射的影響[11-12].

2.3 模型的建立

本文研究材料為Q235鋼的框架式模具，尺寸為2 040 mm×1 140 mm×500 mm，型板厚度為12 mm，支撐結(jié)構(gòu)為7×3的厚度為10 mm的隔板，并帶有一定數(shù)量的散熱孔與通風(fēng)孔結(jié)構(gòu).選用工作尺寸為Φ3 mm×6 mm的熱壓罐，并將其簡(jiǎn)化為圓柱體結(jié)構(gòu).該模型為軸對(duì)稱結(jié)構(gòu)，可對(duì)一半結(jié)構(gòu)進(jìn)行分析以提高效率，如圖3 所示.

圖3 成型模具簡(jiǎn)化建模Fig.3 Simplified modeling of forming die

2.4 網(wǎng)格劃分

本文采用變密度網(wǎng)格.在COMSOL網(wǎng)格界面下選擇用戶控制模式，設(shè)定邊界層網(wǎng)格拉伸因子、厚度調(diào)節(jié)因子及第一層厚度等參數(shù)，其他部分生成自由四面體網(wǎng)格.手動(dòng)在需要細(xì)致反映數(shù)據(jù)變化處加密網(wǎng)格，如模具型面的過(guò)渡處等；而在數(shù)據(jù)變化平緩處疏化網(wǎng)格，如外部熱壓罐區(qū)域等.劃分結(jié)果見圖4.此方法可在保證計(jì)算結(jié)果準(zhǔn)確度的基礎(chǔ)上，減少計(jì)算成本，提高效率.

(a)整體網(wǎng)格劃分

2.5 參數(shù)設(shè)置與邊界條件

根據(jù)雷諾數(shù)所在范圍，流體的流動(dòng)狀態(tài)可分為：

1)層流：Re≤2 300；

2)層流與湍流間：2 300≤Re≤8 000；

3)湍流：Re≥8 000～12 000.

雷諾數(shù)計(jì)算公式為

(7)

式中：u為風(fēng)速；ρ為流體密度；μ為流體動(dòng)力粘度；d為熱壓罐工作直徑.

根據(jù)表1 中空氣的相關(guān)參數(shù)，計(jì)算本文雷諾數(shù)Re=310 977，所以罐內(nèi)流體的流動(dòng)狀態(tài)為湍流.

表1 材料相關(guān)參數(shù)Tab.1 Related parameters of materials

由于風(fēng)扇驅(qū)動(dòng)使得流體強(qiáng)制對(duì)流傳熱，在COMSOL物理場(chǎng)中選擇湍流(k-ε)接口，流體入口速度為1.5 m/s，出口壓力為罐內(nèi)工作壓力0.6 MPa，熱壓罐外壁邊界設(shè)置為壁條件.又由于熱壓罐內(nèi)共軛傳熱的特點(diǎn)，選擇固體和流體傳熱(ht)接口，入口條件為圖2所示的固化溫度函數(shù)，出口為默認(rèn)條件，熱壓罐外壁為熱絕緣，中間面設(shè)置為對(duì)稱.再將以上物理場(chǎng)耦合為非等溫流動(dòng)(nitf1)多物理場(chǎng)，設(shè)置溫度條件為來(lái)自傳熱接口，從而實(shí)現(xiàn)流體與傳熱的雙向耦合.參數(shù)設(shè)置如圖5 所示.

圖5 相關(guān)參數(shù)設(shè)置Fig.5 Related parameter setting

2.6 求解步驟

研究選擇穩(wěn)態(tài)、瞬態(tài)兩個(gè)求解步驟，將穩(wěn)態(tài)的求解結(jié)果作為瞬態(tài)求解的初始值，有利于計(jì)算穩(wěn)定.

2.7 模擬結(jié)果及分析

1)模擬結(jié)果的準(zhǔn)確性驗(yàn)證

根據(jù)參考文獻(xiàn)[13]，仿真時(shí)在模具型面上與實(shí)驗(yàn)監(jiān)測(cè)點(diǎn)對(duì)應(yīng)的位置加入域點(diǎn)探針，如圖6 所示.提取出5處具體溫度值，并與實(shí)驗(yàn)結(jié)果比較，相關(guān)結(jié)果見表2 和表3.對(duì)比表2 和表3 數(shù)據(jù)可知，與實(shí)驗(yàn)值相比，各監(jiān)測(cè)點(diǎn)溫度的仿真值誤差小于3%，所以模擬計(jì)算精度較高，具有可信度.

圖6 監(jiān)測(cè)點(diǎn)位置Fig.6 Location of monitoring points

表2 仿真各監(jiān)測(cè)點(diǎn)溫度Tab.2 Monitoring simulated temperature data

表3 實(shí)驗(yàn)各監(jiān)測(cè)點(diǎn)溫度值Tab.3 Monitoring test temperature data

2)模擬結(jié)果的分析

圖7 為模具型面在不同時(shí)間節(jié)點(diǎn)下的溫度分布云圖.從圖中可以得出：在升溫階段從左端迎風(fēng)區(qū)至右端背風(fēng)區(qū)型面溫度呈下降趨勢(shì)，模具背風(fēng)端中間部位存在一個(gè)低溫區(qū)域；保溫階段溫度分布有均勻化的趨勢(shì)；降溫階段從迎風(fēng)端至背風(fēng)端型面溫度呈上升趨勢(shì)，模具背風(fēng)端中間部位存在一個(gè)高溫區(qū).綜上可得，模具型板中后部存在一個(gè)“溫度滯后區(qū)”.該區(qū)域的溫度變化滯后于型面的其他部位，從而加大模具表面溫度分布的不均勻程度.

(a)t=1 800 s(第一升溫階段)

對(duì)于模具型面溫度分布的仿真模擬結(jié)果可做以下幾點(diǎn)解釋：

1)流體從迎風(fēng)端至背風(fēng)端以強(qiáng)制對(duì)流的方式給模具加熱，根據(jù)上文計(jì)算可知流體雷諾數(shù)較高，在流體與固體交界面處形成帶有脈沖波動(dòng)的湍流邊界層，見圖8.由于流體粘性作用損耗其動(dòng)能，導(dǎo)致從迎風(fēng)端至背風(fēng)端的邊界層流體速度逐漸減小，厚度逐漸增大；又因?yàn)闊嶙枧c厚度成正比，所以流體傳熱能力逐漸下降.

圖8 模具型面邊界層Fig.8 Boundary layer of mold plane

2)迎風(fēng)端的支撐結(jié)構(gòu)在一定程度上會(huì)阻礙流體流動(dòng)并消耗能量，導(dǎo)致流體在背風(fēng)端的沖擊換熱強(qiáng)度降低.

3)模具四周區(qū)域直接與流體接觸，受外界溫度變化影響大，而模具中部區(qū)域會(huì)受到四周結(jié)構(gòu)的阻礙作用.

以上三點(diǎn)充分解釋了模具“溫度滯后區(qū)”的溫度變化落后于其他部位的原因.

3 工藝參數(shù)對(duì)模具溫度場(chǎng)的影響規(guī)律

本文以風(fēng)速、升溫和降溫速率為代表探究工藝參數(shù)對(duì)模具表面溫度分布的影響規(guī)律.僅分析模具型面各點(diǎn)溫度值或溫度差值難以判斷整體溫度分布的均勻性，而方差能夠反映一組數(shù)據(jù)的波動(dòng)情況，所以下文以溫度方差作為模具溫度分布均勻性的衡量標(biāo)準(zhǔn)[14].

3.1 流體流速對(duì)模具型面溫度場(chǎng)的影響規(guī)律

采用第2節(jié)中的模型并保持其他工藝參數(shù)不變，分別設(shè)置流體流速為1.5 m/s，2.5 m/s，3.5 m/s，4.5 m/s，利用COMSOL得到不同流體流速下模具型面的溫度分布云圖，提取1 500個(gè)網(wǎng)格節(jié)點(diǎn)的溫度值并計(jì)算其方差.不同流體流速下，模具型面溫度方差隨時(shí)間變化的趨勢(shì)如圖9 所示.

由圖9 可知改變流體流速對(duì)模具型面溫度分布有顯著影響：1)相同流速下，升、降溫過(guò)程中，模具型面溫度方差不斷增大；保溫過(guò)程中，模具型面溫度方差逐漸減小，溫度分布逐漸趨于均勻.2)流速不同其他條件相同時(shí)，流體流速越快模具型面溫度方差越小，溫度分布均勻性越好，且在保溫階段不同流速的溫度均勻性差異表現(xiàn)更明顯.3)流體流速由 1.5 m/s 增大至 2.5 m/s，由 2.5 m/s 增大至3.5 m/s，由3.5 m/s增大至 4.5 m/s 的溫度方差變化量分別為ΔT1，ΔT2，ΔT3，有ΔT1>ΔT2>ΔT3，由此可知，在一定范圍內(nèi)增大流體流速對(duì)提高模具型面溫度分布的均勻性有非常明顯的效果.

圖9 不同流體流速下型面溫度方差趨勢(shì)Fig.9 The variance trend of plane temperature at different fluid velocities

3.2 升降溫速率對(duì)模具型面溫度場(chǎng)的影響規(guī)律

采用第2節(jié)中的模型并保持其他工藝參數(shù)不變，分別以升溫速率u=1.5 K/min，u=2.5 K/min，u=3.5 K/min，u=4.5 K/min完成熱壓罐第一階段升溫過(guò)程.計(jì)算4種情況下的溫度分布方差，如圖10 所示.

圖10 不同升溫速率下型面溫度方差趨勢(shì)Fig.10 Comparison of temperature variance at different heating rates

采用第2節(jié)中的模型，分別以降溫速率u=1.5 K/min，u=2.5 K/min，u=3.5 K/min，u=4.5 K/min從初始溫度453 K降溫至294 K.計(jì)算4種情況下的溫度分布方差，如圖11 所示.

圖11 不同降溫速率下型面溫度方差趨勢(shì)Fig.11 Comparison of temperature variance at different cooling rates

由仿真結(jié)果可知，改變升、降溫速率對(duì)模具型面溫度分布有顯著影響.升、降溫速率越慢，流體與模具間換熱越充分，溫度方差越小.升、降溫速率基數(shù)越大，增加相同數(shù)值后溫度方差減小得越慢.所以，應(yīng)在適當(dāng)?shù)姆秶鷥?nèi)，綜合考慮成本等因素減小升、降溫速率將有利于提高模具溫度分布的均勻性.

4 結(jié) 論

本文將COMSOL仿真引入到熱壓罐成型模具的溫度場(chǎng)模擬，相關(guān)的研究思想可用于以溫度為導(dǎo)向的多物理場(chǎng)耦合的工程研究中，對(duì)復(fù)合材料成型模具的優(yōu)化具有一定指導(dǎo)意義.

1)選用COMSOL Multiphysics軟件模擬模具型板溫度分布情況，誤差小于3%.

2)研究了流體流速、升降溫速率對(duì)模具溫度分布情況的影響.在適當(dāng)范圍內(nèi)提高流體速度、考慮成本情況下減小升降溫速率均可提高溫度分布的均勻性.

3)模具型板背風(fēng)端存在“溫度滯后區(qū)”，提高此區(qū)域溫度隨固化工藝溫度曲線的變化率，是解決模具溫度分布均勻性問題的關(guān)鍵.

另外，還可以從以下兩方面改進(jìn)熱壓罐成型工藝：從對(duì)流傳熱角度出發(fā)，可以周期性改變流體流動(dòng)方向以避免局部區(qū)域過(guò)熱；從傳導(dǎo)換熱的角度出發(fā)，可以適當(dāng)優(yōu)化模具結(jié)構(gòu)從而提高模具和流體間的傳熱能力.