基于有限元壓縮方法的復合材料RVE 創建1)

田文龍 齊樂華 晁許江

* (西北工業大學機電學院,西安 710072)

? (西北工業大學深圳研究院,廣東深圳 518057)

引言

復合材料具有低密度、優異的力學性能(高比強度、高比剛度、低熱膨脹系數和良好的耐磨性能等)及物理性能(優良的導電、導熱和減振性能等),在航空、航天、國防和汽車等領域得到廣泛的應用[1-5],例如在著名的雙發寬體遠程客機空客A350和波音Boeing787 中使用的復合材料的體積和質量占比分別達80.0%和50.0%.因此,準確表征復合材料的熱-力學性能,建立其微觀結構與熱-力學性能的映射關系,具有非常重要的科學研究意義和工程應用價值,既可以為其工程應用提供理論指導,還有利于其微觀結構優化及制備成形.

實驗測試法具有設備昂貴復雜、測試周期較長和不能得到復合材料微觀變形場等缺點[6],而傳統分析方法不能充分考慮復合材料復雜的微觀結構[7-8],因此數值方法常被用來預測復合材料的熱-力學性能.數值方法既可以充分地考慮復合材料復雜的微觀結構、準確地預測復合材料的熱-力學性能,還能夠提供復合材料各組分的微觀變形信息.當采用數值方法預測復合材料熱-力學性能時,首先需要創建能夠準確表征復合材料微觀結構的代表性體胞單元(RVE),而本文的側重點就在于復合材料RVE 的創建方法.

目前,復合材料RVE 的創建方法主要包括以下3 類: (1)增強體非接觸類方法,包括隨機順序吸附(RSA)法[9-10]和基于分子動力學(MD)的方法[11-12]等;(2)增強體相交移除類方法,主要是增強體遷移法[13-14];(3)基于圖像三維重構技術的方法[15-16].基于圖像三維重構技術的方法通過采用X 射線斷層攝影技術和三維圖像分析技術獲取復合材料微觀結構的層析圖像,然后使用專業的軟件重構出復合材料RVE.該類方法創建的復合材料RVE 與真實的復合材料微觀結構非常接近,然而卻需要昂貴、復雜的實驗設備和專業的軟件,且重構過程非常復雜,導致其應用范圍較小.由于原理及執行算法的簡單,RSA 法是前述提到的復合材料RVE 創建方法中應用最廣泛的方法.然而,該類方法生成的復合材料RVE 中增強體體積分數存在特定的極限值(jamming limit)[17],同時創建增強體數量較多的復合材料RVE 時執行效率較低,且無法控制復合材料RVE 中增強體的取向分布.基于MD 法的復合材料RVE 創建方法可以創建較高增強體體積分數的復合材料RVE 及提高復合材料RVE 創建效率.然而,該類方法的數值執行算法非常復雜,同時也無法控制復合材料RVE 中增強體的取向分布.增強體遷移法能夠創建具有高增強體體積分數的復合材料RVE、提高復合材料RVE 的創建效率及控制復合材料RVE 中增強體的取向分布,然而目前文獻中報道的該類方法只能創建復合材料非周期性RVE,且創建增強體數量較多的復合材料RVE 時執行效率也較低.

近年來,相關學者提出了一種新型增強體非接觸類復合材料RVE 創建方法,即有限元壓縮方法[18-19]:首先創建具有較低增強體體積分數的復合材料RVE,然后采用有限元方法對其進行壓縮,獲取具有較高增強體體積分數的復合材料RVE.有限元壓縮方法簡單、高效,并且具有創建較高增強體體積分數復合材料RVE 的能力.Islam 等[18]將RSA 算法和有限元壓縮方法相結合,創建了體積分數高達50.0%和大長徑比的隨機纖維增強復合材料的RVE(其他RVE 創建方法很難達到).為創建平面纖維增強復合材料的RVE,Zhang 等[19]結合虛擬填裝算法和有限元壓縮方法,創建了纖維體積分數高達45.0%復合材料的RVE,且在纖維動態壓縮過程中基本不改變其取向分布.然而,有限元壓縮方法目前還無法創建復合材料周期性RVE.

為創建具有較高增強體體積分數的復合材料周期性RVE,本文提出一種改進型有限元壓縮方法,該方法主要包括3 個步驟: (1)利用RSA 算法生成具有較低體積分數的復合材料周期性RVE;(2)基于周期性邊界條件和有限元壓縮方法得到有限元網格格式、具有較高體積分數的復合材料周期性RVE;(3)后處理有限元模擬結果以得到復合材料周期性RVE 中所有增強體的位置及取向,進而創建CAD格式、具有較高體積分數的復合材料周期性RVE.采用提出的改進型有限元壓縮方法成功創建了球形增強體復合材料的周期性RVE,其中增強體體積分數可達50%;然后,利用第1,2 和3 階最近鄰距離的概率分布函數[20]和最近鄰取向角的累積概率分布函數[17]、Ripleys-K函數K(r)[17]和對關聯函數G(r)[21]來表征創建的復合材料周期性RVE 中增強體的分布;最后,基于創建的復合材料周期性RVE和有限元均質法預測復合材料的彈性性能,并與實驗測試和雙夾雜模型[22]預測的結果進行對比,驗證提出的改進型有限元壓縮方法的有效性.

1 復合材料RVE 創建方法

本文提出一種用于創建具有較高增強體體積分數的復合材料周期性RVE 的有限元壓縮方法,該方法的示意圖如圖1 所示,主要包括3 個步驟: (1)采用改進的RSA 算法創建具有較低增強體體積分數的復合材料周期性RVE;(2)在周期性邊界條件的約束下,采用有限元方法壓縮第1 步創建的較低增強體積分數的復合材料周期性RVE,以得到有限元網格格式、具有較高增強體體積分數的復合材料周期性RVE;(3)對第2 步的有限元模擬結果進行后處理,得到復合材料周期性RVE 中所有增強體的位置及取向,然后創建CAD 格式、具有較高增強體體積分數的復合材料周期性RVE.

圖1 基于有限元壓縮方法的復合材料周期性RVE 創建過程示意圖Fig.1 Schematic diagram of the generation process of periodic RVE of composites using the FE compression method

本文后續將以球形增強體復合材料RVE 為研究對象,詳細介紹提出的有限元壓縮方法的基本原理.

1.1 復合材料稀疏RVE 創建方法

由于RSA 算法具有原理簡單和執行效率高(當RVE 中增強體體積分數較低時)的優點,本文采用改進RSA 算法創建具有較低增強體體積分數的復合材料周期性RVE,其主要步驟如下:

(1) 創建尺寸為 [-1.5L,1.5L]×[-1.5L,1.5L]×[-1.5L,1.5L](L表示最后得到的復合材料周期性RVE 的邊長)的基體(注意: RVE 的中心點為全局坐標系xyz的原點o);

(2) 在基體范圍內隨機(位置隨機)生成一個球形增強體Ec;

(3) 檢查增強體Ec是否超出基體的邊界,若超出,則在基體的對應位置創建其周期性鏡像Ep;否則,進行下一步操作;

學校發展亦或教改中都存在很多實際問題，面臨種種實際困難。這些問題本身既是問題又是契機，我們必須以問題為導向，抓住學校發展的困難以及各學科獨特的困難，這樣才能精準發力。

(4) 檢查增強體Ec及其周期性鏡像Ep與基體中已存在的增強體Ei是否相交,若不相交,則在基體中保留增強體Ec及其周期性鏡像Ep;否則,刪除增強體Ec及其周期性鏡像Ep,并返回第2 步;

(5) 重復執行第(2)～(4)步,直至基體中增強體體積分數或數量達到預設值.

在創建上述復合材料周期性RVE 過程中,需要保證后續添加的增強體(及其周期性映像)與基體中已存在的增強體之間不能發生相交,即

式中,rc ,rp和ri分別表示后續添加的增強體Ec、其周期性鏡像Ep和基體中已存在的增強體Ei的中心,R表示增強體的半徑,ξ是預設置的增強體間最小分離距離[23].同時,若后續添加的增強體超出基體邊界,則在基體相應的位置創建特定數量的周期性鏡像.對于增強體Ec,其周期性鏡像Ep的數量Np及位置rp可以由下述方法確定.

圖2 復合材料RVE 基體表面標號示意圖Fig.2 Planar and vertex notations of an RVE of composites

在采用RSA算法創建復合材料RVE的過程中,最耗時的步驟是后續添加的增強體(及其周期性映像)與基體中已存在的增強體之間相交狀態的判斷.若基體中存在N個增強體,則后續添加的增強體(及其周期性映像)與基體中已存在的增強體需要進行o(N)次相交判斷.若在RSA算法中引入能夠減少增強體之間相交判斷次數的算法,則可以提高RSA算法創建復合材料RVE的效率.對于后續添加的增強體Ec(及其周期性映像Ep),它們僅可能與其附近的增強體發生相交,因此本文引入RVE分區算法[24],以排除基體中肯定不與增強體Ec(及其周期性映像Ep)發生相交的增強體,進而提高RSA算法創建復合材料周期性RVE的效率.

在RVE 分區算法中,基體被分割成Nc×Nc×Nc個尺寸大于增強體直徑的子區域;然后,根據分割的子區域與基體已存在的增強體之間的相交關系,將基體已存在的增強體分配到不同的子區域中;接下來,采用類似的方法,將后續添加的增強體Ec(及其周期性鏡像Ep)分配到相應的子區域中;最后,只需判斷增強體Ec(及其周期性鏡像Ep)與包含增強體Ec(及其周期性鏡像Ep)的子區域中的增強體之間的相交關系.注意,單個增強體可以被分配到多個子區域中(最多可達8 個).

1.2 復合材料稠密RVE 創建方法

在上述的有限元模擬中,增強體被設置為離散剛體,而RVE 的表面被設置為解析剛體,則采用三維4 節點雙線性剛性四邊形單元(Abaqus/Explicit中的單元類型為R3D4,每個節點具有3 個平移自由度和3 個旋轉自由度)對增強體進行離散.為準確模擬增強體之間的接觸,需將增強體離散成足夠數量的單元,在本文中增強體的網格尺寸s選擇為s=R/5.為考慮增強體之間的接觸,本文采用面-面接觸算法[18],其中接觸發生在兩個方向上: 切向和法向(垂直于接觸表面).對于法向接觸,采用“硬”接觸算法來保證零穿透及增強體之間接觸壓力的傳遞;對于切向接觸,采用庫倫摩擦算法來保證增強體表面的相對滑動,摩擦系數選擇為0.3.此外,為了實現復合材料RVE 壓縮過程的準靜態模擬,分析(壓縮)步的時間分別選擇為 200 s,最小增量步長選擇為0.1 s,增強體的質量設置為 1.0 kg.

復合材料RVE 通常都是周期性的,在采用數值均質法預測復合材料熱-力學性能時,這有利于施加周期性邊界條件[25-26].為了生成復合材料周期性RVE,需將穿過基體邊界的增強體復制并平移到基體的相應位置.因此,本文開發了一種基于周期性邊界條件約束的有限元壓縮方法創建具有較高增強體體積分數的復合材料周期性 RVE,即對于增強體(Ei) 及其周期性鏡像 (Ep) 采用下述周期性邊界條件進行約束.

1.3 CAD 格式復合材料RVE 創建

提出的有限元壓縮方法創建的復合材料周期性RVE 中增強體以R3D4 單元的形式存在,需要將其轉換為CAD 格式的復合材料周期性RVE,以便用于預測復合材料的結構和功能性能.因此,本文通過開發兩個Abaqus Python API 腳本來自動檢測基于有限元壓縮方法得到的復合材料周期性RVE中增強體的中心點和創建CAD 格式的復合材料周期性RVE.注意: 在較低增強體體積分數的復合材料RVE 的壓縮過程中,偶爾會出現一定數量(常小于10 對)的增強體相交現象,此時需要手動刪除這些相交增強體.

2 結果與討論

2.1 復合材料周期性RVE 中增強體分布統計

基于提出的有限元壓縮方法,分別創建了增強體體積分數分別為 30.0%,40.0%和 50.0%的復合材料周期性RVE (如圖3 所示,增強體的數量分別為579,764 和955),其中基體的尺寸為L=100.0 μm、增強體的半徑為R=5.0 μm 和最小分離距離為 ξ=0.05R.在這一部分,本文將分別采用第1、第2 和第3 階最近鄰距離的概率分布函數[20]和最近鄰取向角的累積概率分布函數[17]、Ripleys-K函數K(r)[17]和對關聯函數G(r)[21]來表征創建的復合材料周期性RVE 中增強體的分布規律.

圖3 基于提出的有限元壓縮方法創建的增強體體積分數分別為30.0%,40.0%和50.0%的復合材料周期性RVEFig.3 Periodic RVEs of composites with inclusion volume fractions of 30.0%,40.0% and 50.0% generated using the proposed FE compression method

增強體Ei的 第n階最近鄰距 離dnth定義為該增強體中心點到它的第n階最近鄰增強體中心點的距離,而增強體Ei的第n階最近鄰取向角 θnth和 φnth定義為該增強體中心點到它的第n階最近鄰增強體中心點的連線與z軸的夾角和該增強體中心點到它的第n階最近鄰增強體中心點的連線在xoy平面的投影與x軸的夾角.圖4 和圖5 給出了創建的復合材料周期性RVE 中增強體的第1、第2 和第3 階最近鄰距離的概率分布函數和最近鄰取向角的累積概率函數,可以發現: 第1 階最近鄰距離的概率分布函數曲線初始非常尖銳,然后快速下降,第2 和第3 階最近鄰距離的概率分布函數曲線相對平滑,而最近鄰取向角 θ和 φ的累積概率函數曲線與增強體空間隨機分布的理論累積概率函數曲線(Ψ(θ)=(1-cosθ)/2和 Ψ (φ)=φ/(2π))非常接近.因此,基于提出的有限元壓縮方法創建的復合材料周期性RVE 中增強體空間隨機(CSR)分布.

圖4 復合材料周期性RVE 中增強體的第1、第2 和第3 階最近鄰距離的概率分布函數Fig.4 Probability distribution functions of the 1st,2nd and 3rd nearest neighbor distances of the inclusions in the generated periodic RVEs of composites

圖5 復合材料周期性RVE 中增強體的最近鄰取向角的累積概率分布函數Fig.5 Cumulative probability distribution functions of the 1st,2nd and 3rd nearest neighbor orientation angles of the inclusions in the generated periodic RVEs of composites

Ripleys-K函數K(r) 計算半徑為r的搜索球中包含的增強體中心點的數量與創建的復合材料周期性RVE 中的增強體中心點密度的比值,即

式中,N表示創建的復合材料周期性RVE 中的增強體中心點的數量,V是創建的復合材料周期性RVE的體 積,κij表 示增強體Fi和Fj中心點之 間的距離.對于I(·),若括號中的條件為真,則I(·)的取值為1.0,否則I(·) 的取值為0.對于 ω (ri,rj),若中心點為ri,半徑為 |ri-rj|的球完全包含于創建的復合材料周期性RVE,則 ω (ri,rj)返回數值1.0,否則返回包含在創建的復合材料周期性RVE 中球的體積比例.當創建的復合材料周期性RVE 中增強體空間隨機分布時,Ripleys-K函數K(r) 的計算公式為Kν(r)=4πr3/3.對關聯函數G(r)計算在距離給定的增強體中心點為r范圍內找到另一個增強體中心點的概率,通常被視為Ripleys-K函數的空間導數[21],因此具有下述形式的表達式

式中,G(r)=1.0表示創建的復合材料周期性RVE 中增強體的分布符合CSR 分布.創建的復合材料周期性RVE 中增強體的Ripleys-K函數K(r)和對關聯函數G(r)見圖6 和圖7.可以發現: 隨著搜索半徑r的增加,Ripleys-K函數K(r) 和對關聯函數G(r)的曲線逐漸趨于曲線Kν(r)=4πr3/3 和G(r)=1.0.因此,可以得到下述結論: 基于提出的有限元壓縮方法創建的復合材料周期性RVE 中增強體空間隨機分布.

圖6 復合材料周期性RVE 中增強體的Ripleys-K 函數K(r)Fig.6 Ripleys-K function of the inclusions in the generated periodic RVEs of composites

圖7 復合材料周期性RVE 中增強體的對關聯函數G(r)Fig.7 Pair correlation function G(r) of the inclusions in the generated periodic RVEs of composites

2.2 基于RVE 的復合材料彈性性能預測

在這部分,基于前述創建的復合材料周期性RVE,采用有限元均質方法[27-28]和周期性邊界條件[25-26]對隨機分布球形增強體復合材料的彈性性能進行數值預測.復合材料周期性RVE 的尺寸為L/R=20.0,網格單元選擇4 節點四面體單元(即Abaqus/Standard中C3D4 單元),網格尺寸選為l=R/10.需要說明的是,上述給定的復合材料周期性RVE 和網格的尺寸可以保證得到收斂的復合材料彈性性能.

研究的第1 類SiC/Al 復合材料由2080 鋁合金基體和隨機分布的球形碳化硅(SiC)增強體組成[29],其中基體和增強體的各向同性彈性模量和泊松比分別為E0=74.0 GPa,ν0=0.33 和E1=410.0 GPa,ν1=0.19.研究的第2 類復合材料是球形氫氧化鋁顆粒增強PMMA,其中氫氧化鋁顆粒的體積分數、彈性模量和泊松比分別為v1=0.48,E1=70.0 GPa和 ν1=0.24,而PMMA 的彈性模量和泊松比分別為E0=3.5 GPa 和 ν0=0.31[30].研究的第3 類復合材料是球形玻璃顆粒增強樹脂基復合材料[31],其中玻璃顆粒和樹脂基體的彈性模量和泊松比分別為E0=74.0 GPa,ν0=0.33 和E1=410.0 GPa,ν1=0.19.

基于前述創建的復合材料周期性RVE,采用有限元均質方法和周期性邊界條件預測的增強體體積分數為 50.0%SiC/Al 復合材料的剛度矩陣如下(單位GPa)

因此,可以得到該復合材料的彈性模量、剪切模量和泊松比,發現E1≈E2≈E3, ν12≈ν21≈ν13≈ν31≈ν23≈ν32和G12≈G13≈G23,則可以確定該復合材料的彈性性能是宏觀各向同性的,原因在于該復合材料中增強體的空間隨機分布.后續研究中,我們將使用宏觀各向同性的彈性模量E、剪切模量G和泊松比 μ表征隨機分布球形增強體復合材料的彈性性能.

基于創建的不同增強體體積分數的SiC/Al 復合材料周期性RVE、采用有限元均質方法預測該復合材料的彈性性能,并與實驗測試結果[29]和雙夾雜模型[22]的預測結果進行對比(見圖8),發現有限元均質方法預測的該復合材料的彈性性能與實驗測試和文獻給出的結果及雙夾雜模型的預測結果偏差很小.同時,基于有限元均質方法預測了的顆粒體積分數為0.48 的氫氧化鋁/PMMA 復合材料彈性模量為E=11.06 GPa,與實驗測試的該復合材料的彈性模量(E=10.4 GPa)[30]對比,可以發現兩者之間的相對偏差小于6.35%.對于不同顆粒體積分數的玻璃/樹脂基復合材料,基于有限元均質方法預測和實驗測試的復合材料的彈性性能如表1 所示,發現有限元均質方法預測的該復合材料的彈性性能與實驗測試結果吻合良好.

表1 基于有限元均質方法預測和實驗測試的玻璃顆粒增強樹脂復合材料的彈性性能[31]Table 1 Elastic properties of glass particles reinforced polymer composites predicted using the FE homogenization method and the experimental tests[31]

圖8 基于有限元均質方法、雙夾雜模型和實驗測試得到的不同增強體體積分數SiC/Al 復合材料的彈性性能Fig.8 Elastic properties of SiC/Al composites with different inclusion volume fractions using FE homogenization method,double-inclusion model and available experimental tests

因此,創建的復合材料周期性RVE 可以用于準確預測隨機分布球形增強體復合材料的彈性性能,則驗證了本研究提出的有限元壓縮方法創建復合材料周期性RVE 的有效性.

3 結論

為了簡便、高效地創建具有較高增強體體積分數的復合材料周期性RVE,本文提出了一種改進型有限元壓縮方法.基于提出的多步有限元壓縮方法,成功創建了不同增強體體積分數的隨機分布球形增強體復合材料的周期性RVE.然后,采用最近鄰距離的概率分布函數、最近鄰取向角的累積概率分布函數、Ripleys-K函數和對關聯函數對創建的復合材料周期性RVE 中增強體的分布進行統計.最后,基于創建的復合材料周期性RVE,采用有限元均質方法和周期性邊界條件預測了球形增強體復合材料的彈性性能,并與實驗測試、文獻給出和雙夾雜模型預測的結果進行對比,驗證創建的復合材料周期性RVE 及提出的有限元壓縮方法的有效性.本研究的結論如下.

(1)提出的改進型有限元壓縮方法的具體步驟如下: 生成具有較低增強體體積分數的復合材料周期性RVE;在周期性邊界條件約束下,采用有限元方法壓縮第1 步創建的復合材料周期性RVE,得到具有較高增強體體積分數的復合材料周期性RVE;通過后處理得到復合材料周期性RVE 中所有增強體的位置,進而創建CAD 格式的高增強體體積分數的復合材料周期性RVE.

(2)采用提出的改進型有限元壓縮方法,成功創建了體積分數達50.0%的球形增強體復合材料的周期性RVE;基于第1、第2 和第3 階最近鄰距離的概率分布函數和最近鄰取向角的累積概率分布函數、Ripleys-K函數和對關聯函數對創建的復合材料周期性RVE 中增強體的分布進行統計,發現創建的復合材料周期性RVE 中球形增強體空間隨機分布.

(3)基于有限元均質方法預測的復合材料彈性性能與實驗測試和文獻給出的結果及雙夾雜模型的預測結果偏差很小,說明創建的復合材料周期性RVE可以用于準確預測復合材料的彈性性能,驗證了提出的有限元壓縮方法創建復合材料周期性RVE 的有效性.