先進復(fù)合材料液體成型工藝多尺度樹脂流動數(shù)值模擬研究

李 晨，秦曉陳，陳 程，高麗敏，徐吉峰，2

（1. 中國商飛北京民用飛機技術(shù)研究中心民用飛機結(jié)構(gòu)與復(fù)合材料北京市重點實驗室，北京 102211；2. 太原理工大學(xué)，太原 030024）

相較于傳統(tǒng)結(jié)構(gòu)材料，先進復(fù)合材料因具有較高的比模量、比強度，抗疲勞性能和耐腐蝕性好，便于一體化成型等諸多優(yōu)點，被廣泛應(yīng)用于航空航天、交通設(shè)備制造等領(lǐng)域[1]。當(dāng)前，復(fù)合材料在民機結(jié)構(gòu)上的應(yīng)用比例在一定程度上已成為衡量民機先進性的重要指標(biāo)。民用航空領(lǐng)域?qū)?fù)合材料的經(jīng)濟性和安全性提出了更高要求，復(fù)合材料制品的成本主要源于原材料和成型工藝兩部分，其中制造方面的成本通常約占到總成本的70%[2]。

在諸多制造工藝中，復(fù)合材料熱壓罐成型是目前應(yīng)用最廣泛、用量最大的復(fù)合材料與結(jié)構(gòu)制備技術(shù)。作為一種具有代表性的先進復(fù)合材料低成本制造技術(shù)，復(fù)合材料液體成型（Liquid composites molding，LCM）工藝可以避免使用造價高昂的熱壓罐設(shè)備，而選擇使用成本較低的烘箱或者固化爐等，同時該工藝具有制件表面質(zhì)量好、精度高、制造效率高等諸多優(yōu)點而得到廣泛關(guān)注。如俄羅斯MC–21客機和空中客車公司A220系列飛機等都已經(jīng)將自動鋪絲液體成型技術(shù)用于飛機重要承力結(jié)構(gòu)的制造中[3]。真空輔助樹脂轉(zhuǎn)移成型 （Vacuum assisted resin transfer molding，VARTM）工藝作為液體成型工藝中較常見的方法，近年來也已被廣泛應(yīng)用于復(fù)合材料制件成型。典型的VARTM工藝樹脂流動浸漬的過程如圖1所示[4]。樹脂浸漬過程是該方法制備復(fù)合材料是否產(chǎn)生缺陷的關(guān)鍵環(huán)節(jié)。

圖1 VARTM工藝樹脂浸漬過程示意圖[4]Fig.1 Sketch of VARTM during resin infusion process[4]

傳統(tǒng)上采用液體成型工藝進行復(fù)合材料件制備和工藝優(yōu)化前，往往需要依靠經(jīng)驗設(shè)計樹脂注口位置并通過大量試驗試制對樹脂浸漬情況進行判定，不僅影響研制周期，還會帶來極大的試制成本。采用數(shù)值模擬對模具內(nèi)樹脂流動情況進行工藝仿真分析可以有效降低工藝設(shè)計成本，有利于優(yōu)化工藝參數(shù)，并提高制造效率。國外許多學(xué)者開發(fā)了LIMS程序，從試驗和仿真方面對液體成型工藝進行了較為系統(tǒng)的研究[5–8]。國內(nèi)學(xué)者也針對液體成型工藝仿真開展了廣泛研究[9–11]。戴福洪等[12]采用控制體積法模擬邊緣效應(yīng)，然而控制體積法計算精度較低，一定程度上影響了預(yù)報結(jié)果。Li等[13–14]采用有限元方法和相場法對樹脂在微觀和細觀層面的流動進行了研究。目前，國內(nèi)外對液體成型樹脂流動的數(shù)值模擬研究主要集中在宏觀尺度，然而工藝缺陷等往往是微觀和介觀尺度問題在宏觀層面的集中表現(xiàn)。本文通過有限元方法，基于Comsol Multiphysics多物理場仿真軟件研究了纖維預(yù)成型體從微觀到介觀尺度的滲透率，在宏觀尺度上采用理查德方程作為流動控制方程，實現(xiàn)了對充模過程樹脂瞬時流動前鋒和流場壓強進行預(yù)報，并在各個尺度層級上與現(xiàn)有理論或試驗進行了對比，驗證本文方法的正確性。

1 復(fù)合材料預(yù)成型體多尺度結(jié)構(gòu)與滲透率

通常，纖維單絲 （μm級別）層面被認為是不可被滲透的最小尺度，即微觀尺度；而纖維束內(nèi)和束間 （μm～mm級別）介于微觀和宏觀制件尺度間則被稱為介觀或細觀。圖2展示了典型的纖維單絲微觀尺度、纖維束內(nèi)和束間介觀尺度及宏觀尺度形貌特征，其中纖維束外輪廓一般呈橢圓形結(jié)構(gòu)。

在微觀尺度上，纖維單絲間的隨機排布形成樹脂流道，符合層流自由流的特征；而在介觀層面，當(dāng)纖維束被看作一個整體時，該區(qū)域符合多孔介質(zhì)滲流特點，纖維束外仍是層流自由流；在宏觀尺度上，則可從整體上認為是一個多孔介質(zhì)滲流區(qū)域。

滲透率是衡量多孔介質(zhì)結(jié)構(gòu)允許樹脂流動滲透能量的表征量，是材料固有屬性，與所滲透的液體和外施條件無關(guān)，滲透率也是數(shù)值模擬中重要的參數(shù)。滲透率常用K表示，其定義由Darcy定理引出，即

三維滲透率張量可表示為

宏觀尺度由微觀和介觀經(jīng)過尺度升級獲得，逐級將滲透率數(shù)據(jù)代入代表性單胞結(jié)構(gòu)有限元模型中，多尺度的幾何結(jié)構(gòu)模型成為整個多尺度樹脂流動數(shù)值模擬研究的基礎(chǔ)和載體，而多尺度的滲透率則成為連接不同尺度進行數(shù)值模擬尺度升級的主線，不同尺度滲透率之間的聯(lián)系如圖3所示。

圖3 不同尺度滲透率的聯(lián)系Fig.3 Relationship of multi-scale permeability

2 基于有限元方法的多尺度樹脂流動數(shù)值模擬研究流程

本文從微觀尺度結(jié)構(gòu)建模和滲透率預(yù)報研究入手，進而擴展到介觀尺度并逐級對研究尺度進行升級，最后得到宏觀尺度滲透率，并基于此開展全局樹脂流動數(shù)值模擬研究與試驗結(jié)果進行比對。整個研究流程如圖4所示。

圖4 多尺度樹脂流動數(shù)值模擬研究流程Fig.4 Flowchart of numerical simulation of resin flow at multi-scale

對樹脂浸漬流動進行數(shù)值模擬可以歸屬為計算流體力學(xué)的研究范疇，控制方程是數(shù)值模擬計算過程的核心，是對物理過程的公式化表述；有限元法則通過對計算域進行網(wǎng)格劃分、控制方程和初始條件的離散，進行求解以達到迭代計算的結(jié)果滿足收斂誤差條件。

3 多尺度樹脂流動數(shù)值模擬理論

本文為滿足控制方程的使用條件，進行如下假設(shè)：

（1）認為纖維預(yù)成型體為符合多孔介質(zhì)定義的連續(xù)結(jié)構(gòu)，且在樹脂流動浸漬過程中對結(jié)構(gòu)的形變影響忽略不計；

（2）樹脂在纖維預(yù)成型體中的流動為多孔介質(zhì)滲流并滿足達西流體流動特性；

（3）溫度恒定且樹脂流動不受溫度變化影響，即樹脂流動時黏度保持一定，且不考慮樹脂流動時的相變和物理化學(xué)反應(yīng)等。

3.1 微觀尺度樹脂流動理論

Navier–Stokes方程是描述黏性不可壓縮流體動量守恒運動的經(jīng)典方程，通過求解式（2）和質(zhì)量守恒方程，即連續(xù)性方程（式 （3）），可以得到速度場u，代入Darcy定理計算得到滲透率。

式中，ρ為樹脂的密度；u為流動速度；p為壓強；μ為動態(tài)黏度；F為體積力；I為單位向量；ε為纖維預(yù)成型體的孔隙率。在穩(wěn)態(tài)流動下，流體不可被壓縮，ρ為常數(shù)，式 （3）可化簡為式（4），即

式中，v為速度向量。

3.2 介觀尺度層流和滲流耦合流動方法研究

纖維束內(nèi)部是由纖維單絲和它們之間的間隙構(gòu)成的多孔介質(zhì)結(jié)構(gòu)。樹脂流動穿過多孔介質(zhì)結(jié)構(gòu)為多孔介質(zhì)流動，樹脂在纖維束之間的空隙流動為自由流。當(dāng)忽略慣性項時，方程變?yōu)镾tokes方程。由于雷諾數(shù)遠小于1，且認為流體不可壓縮，自由流變?yōu)槿鋭恿?，并且在多孔介質(zhì)外部和內(nèi)部的流動是連續(xù)的。兩種流動在微觀尺度上可以用統(tǒng)一的控制方程Brinkman公式來描述。Brinkman公式在Stokes方程基礎(chǔ)上引入Darcy項，拓展了應(yīng)用范圍，所形成的Stokes–Brinkman方程恰當(dāng)?shù)孛枋隽藢恿骱蜐B流的耦合流動，即

在多孔介質(zhì)中，Darcy項“μK–1·u”起著關(guān)鍵作用，而在自由流域中，由于ε= 1，K趨于無窮大，Darcy項消失，方程可簡化為Stokes方程。

3.3 宏觀尺度樹脂流動模擬

宏觀尺度控制方程采用多孔介質(zhì)滲流理查德方程，即

式（7）為樹脂流動輪廓方程，建立了樹脂流動飽和度S和壓強p之間的關(guān)系，直觀反映了宏觀尺度上樹脂流動前鋒的位置并與實際樹脂浸漬自動干絲鋪放的預(yù)成型體試驗及理論解析解進行對比驗證，能夠?qū)渲髁?、用量和樹脂充模完成時間進行預(yù)報。

4 多尺度有限元建模及驗證

4.1 微觀尺度滲透率計算及驗證

根據(jù)采集的微觀尺度纖維束內(nèi)結(jié)構(gòu)形貌，基于Monte Carlo隨機算法建立纖維束內(nèi)正方形區(qū)域幾何模型 （孔隙率ε= 0.5，纖維單絲直徑rf= 0.0035 mm，纖維數(shù)量n= 122）。對于滲透率計算的穩(wěn)態(tài)數(shù)值模擬，本文采用P1–P1進行有限元離散化，用線性單元P1對速度u和壓力p進行離散化，網(wǎng)格劃分情況和速度場如圖5所示。

圖5 纖維束內(nèi)正方形區(qū)域網(wǎng)格劃分與速度場分布Fig.5 Fibers intra square domain meshing and velocity distribution

由圖6可知，基于Monte Carlo隨機排布結(jié)構(gòu)計算的滲透率值低于Kozeny–Carman方程和Gebart方程[15]等經(jīng)驗公式，但與其他采用隨機建模方法學(xué)者計算值接近。這是由于規(guī)則排布的纖維間容易形成較通暢的“流道”，而隨機排布的結(jié)構(gòu)纖維相互交錯容易造成“流道”的堵塞，進而導(dǎo)致流體難以滲透通過，并直觀表現(xiàn)為滲透率值的下降[16–17]。

圖6 隨機排布與經(jīng)驗公式滲透率結(jié)果對比Fig.6 Comparison of permeability results between random configuration and empirical formula

4.2 介觀尺度滲透率計算和驗證

本文基于平紋玻璃纖維編織結(jié)構(gòu)，采用TexGen軟件建立介觀尺度纖維束內(nèi)和束外單胞結(jié)構(gòu)，如圖7所示，前后兩側(cè)分別為入口和出口條件，其他4個面為周期性邊界條件。

圖7 介觀尺度可重復(fù)單胞結(jié)構(gòu)Fig.7 Repeatable unit cell structure at mesoscale

將微觀尺度計算獲得的滲透率數(shù)值代入纖維束區(qū)域，計算纖維束內(nèi)與束間單胞結(jié)構(gòu)的滲透率，所得的速度流場如圖8所示。

圖8 單胞結(jié)構(gòu)速度場Fig.8 Velocity field of unit cell

將基于本文模型方法的計算結(jié)果與文獻[18–19]改進的單胞滲透率計算值相比，誤差在5%以內(nèi)，如表1所示。

表1 滲透率對比[18–19]Table 1 Comparison of permeability[18–19]

4.3 宏觀尺度數(shù)值模擬與試驗驗證

（1）試驗設(shè)計和有限元建模。為驗證本文方法的正確性，本文將宏觀尺度與樹脂浸漬試驗 （圖9）進行對比，試驗設(shè)備由樹脂注射系統(tǒng)和位于樹脂注射系統(tǒng)上方的監(jiān)測系統(tǒng)兩部分構(gòu)成，其中監(jiān)測系統(tǒng)負責(zé)記錄試驗中不同時刻樹脂流動前鋒位置。本文采用Textile Aero Tarn公司生產(chǎn)的平紋編織玻璃纖維材料作為預(yù)成型體進行單層鋪放，外覆真空袋。

圖9 樹脂流動浸漬試驗設(shè)計Fig.9 Experiment design of resin flow and infusion

在纖維預(yù)成型體上預(yù)設(shè)了兩個鏤空區(qū)域，在臨近下端的位置預(yù)設(shè)阻礙區(qū)域與整個矩形邊緣區(qū)域相聯(lián)通，忽略厚度方向樹脂浸漬影響，采用二維平面模型模擬樹脂流動。運用三角形單元對幾何模型進行有限元劃分后的結(jié)果如圖10所示，采用三角形網(wǎng)格，共劃分為5760個域單元。

圖10中除了兩個位于幾何模型中部的圓形樹脂注入口，其他邊界條件均設(shè)為不可滲透。初始狀態(tài)模具內(nèi)壓強與出口壓強差為0；樹脂注入初始階段，入口壓強 （Pin）和出口壓強 （Pout）壓差為0.1 MPa，試驗和數(shù)值模擬選用的樹脂黏度和纖維預(yù)成型體的滲透率值沿豎直方向K1和水平方向K2由前述章節(jié)中方法計算獲得。具體參數(shù)如表2所示。

表2 材料參數(shù)Table 2 Material property

圖10 網(wǎng)格劃分Fig.10 Meshing

（2）結(jié)果對比驗證。不同時刻的樹脂流動前鋒位置和壓力場的試驗結(jié)果和數(shù)值模擬結(jié)果如圖11所示。其中，圖11（a）和（e）為不同時刻試驗所得樹脂流動前鋒的記錄。圖11（c）和（g）為Gantois等[20]采用邊界元法得到的壓強場云圖。圖11（b）和（f）、圖11（d）和（h）分別為本文所采用方法獲得樹脂流動前鋒和壓強場的結(jié)果。

圖11 不同時刻的樹脂流動前鋒位置和壓力場的試驗結(jié)果和數(shù)值模擬結(jié)果對比Fig.11 Comparison of experimental results and numerical simulation results of resin flow front position and pressure field at different time

為量化試驗與數(shù)值模擬對比結(jié)果，本文提出了浸潤比的定義，即λ=S樹脂/S纖維，其中S樹脂代表瞬時投影面內(nèi)被樹脂流動浸漬區(qū)域面積；S纖維代表整個纖維預(yù)成型體鋪放區(qū)域面積。對比不同時刻數(shù)值模擬結(jié)果與試驗結(jié)果的浸潤比λ絕對值差|λ試驗–λ模擬|均小于10%；同時，將不同時刻壓強場結(jié)果與Gantois等[20]結(jié)果進行了對比，亦得到近似結(jié)果，證明了本文方法具有良好的預(yù)報準(zhǔn)確性。

5 結(jié)論

本文研究了從纖維束內(nèi)部微觀尺度到纖維束內(nèi)和束間介觀尺度滲透率預(yù)報方法，并進行對比驗證；基于微觀和介觀尺度滲透率進行尺度升級，在宏觀尺度實現(xiàn)了對樹脂浸漬纖維預(yù)成型體過程瞬時流動前鋒的數(shù)值模擬預(yù)測并通過提出浸潤比的概念對樹脂浸潤情況進行量化與試驗結(jié)果對比。本文初步建立了一套較為完整的跨越微觀、介觀和宏觀多尺度域集預(yù)成型體滲透率預(yù)報與樹脂滲流的數(shù)值模擬方法體系。

在后續(xù)數(shù)值模擬研究工作中，可以考慮溫度變化對流體流動特性的影響（如樹脂黏度的影響）。本文所建立的方法，目前僅對較為簡單的復(fù)合材料結(jié)構(gòu)的樹脂流動仿真進行了驗證，在后續(xù)研究中可以繼續(xù)將該方法擴展到更復(fù)雜結(jié)構(gòu)模型中進行探討研究。