萬(wàn)小朋, 楊光猛, 趙美英

(西北工業(yè)大學(xué) 航空學(xué)院, 陜西 西安 710072)

纖維增強(qiáng)復(fù)合材料具有比強(qiáng)度高、比剛度大、可設(shè)計(jì)強(qiáng)以及良好的抗疲勞性能等優(yōu)點(diǎn),在航空航天等眾多領(lǐng)域得到廣泛應(yīng)用。然而復(fù)合材料復(fù)雜的損傷失效過(guò)程以及破壞模式給力學(xué)性能仿真分析以及結(jié)構(gòu)應(yīng)用帶來(lái)巨大的難度[1-2]。

近年來(lái)基于代表性體積單元(representative volume element,RVE)的微觀力學(xué)損傷模型在預(yù)測(cè)復(fù)合材料損傷過(guò)程和強(qiáng)度性能中得到廣泛應(yīng)用[3-5]。微觀力學(xué)模型中對(duì)纖維分布形貌以及組分材料本構(gòu)關(guān)系的準(zhǔn)確描述是預(yù)測(cè)單向復(fù)合材料力學(xué)性能的關(guān)鍵。González等[6]對(duì)比了6組纖維隨機(jī)分布模型仿真得到的復(fù)合材料橫向壓縮應(yīng)力-應(yīng)變曲線(xiàn),結(jié)果表明纖維分布形式對(duì)損傷后的應(yīng)力-應(yīng)變曲線(xiàn)變化有一定影響。Totry等[7-8]進(jìn)一步研究表明不同隨機(jī)分布形式在預(yù)測(cè)復(fù)合材料C/epoxy和C/peek橫向壓縮強(qiáng)度的偏差分別為3%和7%。上述研究表明纖維分布形式會(huì)對(duì)單向復(fù)合材料力學(xué)性能的預(yù)測(cè)精度產(chǎn)生影響,需要在數(shù)值分析中加以考慮。

組分材料性能參數(shù)是影響單向復(fù)合材料力學(xué)性能的一個(gè)重要因素。一般來(lái)說(shuō)組分材料中纖維與基體的性能可以通過(guò)強(qiáng)度試驗(yàn)獲得,然而界面強(qiáng)度很難通過(guò)傳統(tǒng)試驗(yàn)方法獲得,并且對(duì)于同一材料,界面脫黏起始點(diǎn)判定、試驗(yàn)方法以及試樣幾何差異等因素的不同,同樣會(huì)導(dǎo)致測(cè)出的數(shù)據(jù)分散性很大[9]。當(dāng)前數(shù)值分析中界面強(qiáng)度的選取主要有3種方案：①采用Mohr-Coulomb(M-C)彈塑性準(zhǔn)則對(duì)基體力學(xué)行為進(jìn)行描述時(shí),選取M-C準(zhǔn)則中基體的黏聚力作為界面強(qiáng)度[6,10-11]；②采用基體強(qiáng)度作為界面強(qiáng)度[12]；③通過(guò)分析界面強(qiáng)度對(duì)仿真結(jié)果的影響規(guī)律,反演出與試驗(yàn)結(jié)果較為吻合的界面參數(shù)值[5,13]。然而3種方案中方案①和方案②只對(duì)于特定材料性能的仿真具有很好的精度，方案③需要進(jìn)行復(fù)合材料力學(xué)性能試驗(yàn),花費(fèi)較高且時(shí)間歷程過(guò)長(zhǎng)。

基于上述原因,本文通過(guò)C++編寫(xiě)纖維隨機(jī)分布算法,使用python在Abaqus中參數(shù)化生成包含纖維、基體和界面的復(fù)合材料代表性體積單元模型,采用擴(kuò)展的線(xiàn)性Drucker-Prager(D-P)彈塑性準(zhǔn)則[14]和雙線(xiàn)性?xún)?nèi)聚力模型分別描述基體和界面的本構(gòu)關(guān)系,最終建立了復(fù)合材料二維損傷微觀力學(xué)模型。在此基礎(chǔ)上研究分析纖維分布形式與界面強(qiáng)度對(duì)單向復(fù)合材料橫向壓縮強(qiáng)度的影響。

1 微觀力學(xué)數(shù)值分析模型

1.1 微觀結(jié)構(gòu)中纖維分布的生成

在單向復(fù)合材料固化過(guò)程中,纖維單絲浸潤(rùn)在流動(dòng)性樹(shù)脂中,鋪設(shè)出的單向?qū)永w維呈隨機(jī)分布。為了研究纖維隨機(jī)分布對(duì)單向?qū)恿W(xué)性能的影響,采用隨機(jī)擴(kuò)張算法[15]生成纖維隨機(jī)分布模型,生成步驟如下所示:

1) 輸入纖維體積含量和纖維半徑,計(jì)算出目標(biāo)區(qū)域尺寸；

2) 利用隨機(jī)函數(shù)生成第一根纖維的圓心坐標(biāo)(x1,y1),設(shè)為當(dāng)前纖維；

3) 生成第二根纖維(x2,y2)。此點(diǎn)位置由(x2,y2)與當(dāng)前纖維(x1,y1)之間的距離l和夾角θ決定,l和夾角θ為隨機(jī)數(shù),滿(mǎn)足lmin+2r

4) 重復(fù)步驟3),直到在當(dāng)前纖維周?chē)荒茉偕尚碌睦w維。判斷新生成纖維與已生成纖維之間距離是否小于lmin,若小于則刪去此纖維嘗試重新生成。若嘗試次數(shù)達(dá)到規(guī)定,則進(jìn)行步驟5)；

5) 選取第二根纖維為當(dāng)前纖維,重復(fù)步驟3)和步驟4)；

6) 對(duì)于每一根纖維都重復(fù)步驟3)～步驟5),直至纖維體積含量滿(mǎn)足或目標(biāo)區(qū)域已填滿(mǎn),結(jié)束。

本文纖維隨機(jī)分布算法采用C++編程,再通過(guò)python語(yǔ)言在Abaqus中參數(shù)化生成包含纖維、基體和界面的復(fù)合材料RVE模型。

1.2 FEM有限元模型

為了準(zhǔn)確預(yù)測(cè)復(fù)合材料的宏觀力學(xué)性能,RVE模型必須包含足夠多的纖維。文獻(xiàn)[6]指出RVE模型包含至少30根纖維時(shí)足以表征材料的性能。本文采用上述算法生成5組不同的纖維隨機(jī)分布模型,模型包含36根纖維,纖維半徑為5 μm,體積含量為60%,RVE模型的尺寸為68.7 μm×68.7 μm。其中一組纖維隨機(jī)分布RVE模型如圖1所示。

圖1 復(fù)合材料纖維隨機(jī)分布RVE有限元模型

有限元模型采用平面應(yīng)變單元CPE3劃分基體和纖維網(wǎng)格,以及二維單元COH2D4劃分界面網(wǎng)格,每個(gè)模型平均包含單元數(shù)約19 000。為了克服非線(xiàn)性數(shù)值求解中的收斂性問(wèn)題,基于Abaqus/explicit進(jìn)行RVE模型的準(zhǔn)靜態(tài)分析,設(shè)定模型質(zhì)量放大系數(shù)為1 000,以保證準(zhǔn)靜態(tài)求解的同時(shí)減小計(jì)算時(shí)間。

RVE模型施加周期性邊界條件,在相對(duì)平行的2條邊以及角點(diǎn)上滿(mǎn)足

(1)

式中，u表示邊界上節(jié)點(diǎn)的位移向量,下標(biāo)L,R,T,B代表4條邊,LB,LT,RB,RT代表4個(gè)角點(diǎn),如圖1所示。在Abaqus中通過(guò)“equation constraint”來(lái)實(shí)現(xiàn)公式(1)。同時(shí)分別約束角點(diǎn)RB和LT的水平自由度和豎直自由度,固定角點(diǎn)LB,以限制RVE模型的剛體運(yùn)動(dòng)。角點(diǎn)RB處施加水平位移載荷,如圖1所示。

1.3 基體與界面本構(gòu)關(guān)系

樹(shù)脂基體在橫向壓縮載荷作用下表現(xiàn)出明顯的彈塑性力學(xué)行為,并且其屈服行為對(duì)靜水壓力敏感。結(jié)合擴(kuò)展的線(xiàn)性D-P屈服準(zhǔn)則和塑性準(zhǔn)則模擬基體材料在靜水壓力下的塑性以及損傷退化力學(xué)行為,D-P屈服準(zhǔn)則控制方程為

F=t-ptanβ-d=0

(2)

式中,p為靜水壓力,q為Mises等效應(yīng)力,r為偏應(yīng)力第三不變量,β和d分別是D-P準(zhǔn)則的內(nèi)摩擦角和黏聚力。k是三軸拉伸強(qiáng)度與三軸壓縮強(qiáng)度之比,引入了基體材料在拉伸和壓縮情況下不同的屈服行為。k為0.78,β為19.11°,d為83.1 MPa。

M-C屈服準(zhǔn)則同樣對(duì)基體的彈塑性行為進(jìn)行了描述,非相關(guān)聯(lián)流動(dòng)法則中M-C準(zhǔn)則參數(shù)可以通過(guò)與D-P準(zhǔn)則的相互轉(zhuǎn)換關(guān)系(3)式獲得

(3)

式中，φ,c分別為M-C準(zhǔn)則的內(nèi)摩擦角和黏聚力。

采用塑性準(zhǔn)則定義材料初始損傷的等效塑性應(yīng)變?chǔ)舙l,它表示為應(yīng)力三軸度η的函數(shù)。η刻畫(huà)了不同的三維應(yīng)力狀態(tài),取0,1/3,-1/3分別表示純剪切、單軸拉/壓等應(yīng)力狀態(tài),對(duì)應(yīng)的εpl分別為0.5,0.05和0.5[10]。初始損傷發(fā)生后,基體表現(xiàn)出脆性斷裂,取斷裂能Gm=0.5 kJ/m2。

纖維/基體的界面采用經(jīng)典的雙線(xiàn)性?xún)?nèi)聚力單元表征其損傷起始與退化的力學(xué)行為

(4)

(5)

(4)式為初始損傷控制方程,(5)式為損傷擴(kuò)展控制方程。T和G分別表示為界面應(yīng)力與能量釋放率,下標(biāo)n,t,s分別表示正向和2個(gè)切向方向。選取內(nèi)聚力單元的初始剛度K=108GPa[5]。界面正向強(qiáng)度Tn0=120 MPa,與基體壓縮強(qiáng)度相同,切向強(qiáng)度Ts0=180 MPa,正向斷裂能Gnc=3 kJ/m2,切向斷裂能Gsc=6 kJ/m2。

2 橫向壓縮性能計(jì)算及結(jié)果分析

2.1 RVE模型驗(yàn)證

為了驗(yàn)證所建立的RVE有限元模型在預(yù)測(cè)復(fù)合材料橫向壓縮強(qiáng)度的準(zhǔn)確性,選取試驗(yàn)中復(fù)合材料E-glass/MY750作為研究對(duì)象[16]。文獻(xiàn)[5]對(duì)該材料進(jìn)行了有限元數(shù)值仿真,與本文仿真結(jié)果進(jìn)行對(duì)比驗(yàn)證。纖維E-glass彈性模量Ef=74.0 GPa,泊松比vf=0.2,在橫向壓縮載荷作用下認(rèn)為不發(fā)生破壞[5,8]。基體彈性模量Em=3.35 GPa,泊松比vm=0.35,拉伸強(qiáng)度σmt和壓縮強(qiáng)度σmc分別為80 MPa和120 MPa。

圖2所示為RVE有限元模型預(yù)測(cè)的橫向壓縮應(yīng)力-應(yīng)變曲線(xiàn)與文獻(xiàn)[5]仿真結(jié)果的對(duì)比?？梢钥闯?RVE模型預(yù)測(cè)的復(fù)合材料橫向壓縮強(qiáng)度為138.8 MPa,與文獻(xiàn)[16]中試驗(yàn)值145 MPa的誤差在4.3%以?xún)?nèi),比文獻(xiàn)[5]仿真結(jié)果更接近試驗(yàn)值,說(shuō)明本文所建立的RVE模型可以很好地預(yù)測(cè)單向復(fù)合材料橫向壓縮強(qiáng)度。

圖2 復(fù)合材料橫向壓縮應(yīng)力-應(yīng)變曲線(xiàn)

本文RVE模型與文獻(xiàn)[5]預(yù)測(cè)的失效應(yīng)變都略高于試驗(yàn)的失效應(yīng)變1.2%,這主要是由于數(shù)值仿真假設(shè)基體的硬化模量為零,而實(shí)際上在橫向壓縮載荷作用下,單向復(fù)合材料中基體會(huì)產(chǎn)生塑性硬化,延性減弱,進(jìn)而使得失效應(yīng)變下降。

2.2 纖維分布對(duì)復(fù)合材料橫向壓縮力學(xué)性能的影響

為了研究纖維隨機(jī)分布對(duì)單向復(fù)合材料橫向壓縮力學(xué)性能的影響,建立5組不同纖維分布的RVE模型進(jìn)行分析,計(jì)算獲得在橫向壓縮載荷作用下復(fù)合材料的應(yīng)力-應(yīng)變曲線(xiàn)如圖3所示。

圖3 不同纖維分布形式預(yù)測(cè)的復(fù)合材料橫向壓縮應(yīng)力-應(yīng)變曲線(xiàn)

由應(yīng)力-應(yīng)變曲線(xiàn)可以看出,5組曲線(xiàn)在彈性段幾乎重合,隨著初始損傷的發(fā)生,曲線(xiàn)開(kāi)始分離,且隨著應(yīng)變的增大差別越來(lái)越明顯,最終導(dǎo)致壓縮強(qiáng)度出現(xiàn)差異。這主要是因?yàn)楦鹘MRVE模型中纖維分布位置、尤其是纖維間最近距離分布不同,使得初始損傷產(chǎn)生區(qū)域以及最終破壞路徑出現(xiàn)不同,從而引起復(fù)合材料應(yīng)力-應(yīng)變曲線(xiàn)以及橫向壓縮強(qiáng)度的差別。

表1 不同纖維分布形式預(yù)測(cè)的復(fù)合材料橫向壓縮強(qiáng)度與試驗(yàn)值誤差

表1列出了5組纖維分布形式預(yù)測(cè)的復(fù)合材料橫向壓縮強(qiáng)度與試驗(yàn)值的誤差?？梢钥闯?組纖維分布形式的預(yù)測(cè)結(jié)果與試驗(yàn)值較為吻合,最大誤差為6.2%,說(shuō)明本文RVE模型尺寸以及界面強(qiáng)度的取值可以很好地預(yù)測(cè)單向復(fù)合材料的宏觀強(qiáng)度性能。但不同纖維分布之間預(yù)測(cè)的橫向壓縮強(qiáng)度存在偏差,模型RVE-2與RVE-4之間的偏差達(dá)到8.5%。由此可見(jiàn)為了更準(zhǔn)確地預(yù)測(cè)單向復(fù)合材料力學(xué)性能以及損傷破壞形貌,需要建立符合真實(shí)微觀纖維分布特征的RVE模型。

2.3 界面強(qiáng)度對(duì)復(fù)合材料橫向壓縮力學(xué)性能的影響

為了對(duì)比現(xiàn)有研究中界面強(qiáng)度的3種選取方案對(duì)復(fù)合材料橫向壓縮力學(xué)性能的影響,選取界面強(qiáng)度Tn0分別為:

a) M-C準(zhǔn)則的黏聚力

根據(jù)關(guān)系(3)式可以獲得為49 MPa;

b) 基體壓縮強(qiáng)度

基體壓縮強(qiáng)度為120 MPa;

c) 根據(jù)試驗(yàn)數(shù)據(jù)反演界面強(qiáng)度

通過(guò)選取不同界面強(qiáng)度,仿真獲得相應(yīng)的復(fù)合材料橫向壓縮強(qiáng)度。比較發(fā)現(xiàn)當(dāng)界面強(qiáng)度為140 MPa時(shí),復(fù)合材料橫向壓縮強(qiáng)度與試驗(yàn)值最為吻合,因此選擇140 MPa作為反演法得到的界面強(qiáng)度。

3種界面強(qiáng)度選取方案預(yù)測(cè)的復(fù)合材料橫向壓縮強(qiáng)度以及與試驗(yàn)值的誤差如表2所示。

表2 3種界面選取方案預(yù)測(cè)的復(fù)合材料橫向壓縮強(qiáng)度與試驗(yàn)誤差

可以看出,界面強(qiáng)度為M-C準(zhǔn)則黏聚力時(shí)對(duì)復(fù)合材料橫向壓縮強(qiáng)度的估值明顯偏低。隨著界面強(qiáng)度提高,復(fù)合材料橫向壓縮強(qiáng)度升高,并且逐漸逼近試驗(yàn)值。這說(shuō)明此時(shí)在橫向壓縮載荷作用下,界面損傷對(duì)復(fù)合材料橫向壓縮強(qiáng)度預(yù)測(cè)精度產(chǎn)生了影響。

方案c中界面強(qiáng)度140 MPa時(shí)預(yù)測(cè)的橫向壓縮強(qiáng)度顯然與試驗(yàn)值的誤差最小。此時(shí)復(fù)合材料壓縮強(qiáng)度基本不再隨界面強(qiáng)度增大而變化,這說(shuō)明實(shí)際上基體塑性損傷,而不是界面損傷在復(fù)合材料橫向壓縮破壞過(guò)程中起主導(dǎo)作用,這一結(jié)論也與文獻(xiàn)[6]得到的結(jié)果相符。

圖4進(jìn)一步給出了界面強(qiáng)度為140 MPa時(shí),復(fù)合材料在橫向壓縮載荷作用下的損傷斷裂形貌。可以看出在橫向壓縮載荷作用下,復(fù)合材料最終形成與加載垂直方向成48°的剪切斷裂帶。這與M-C彈塑性理論計(jì)算得到的基體斷裂角50.8°相吻合[5],說(shuō)明基體剪切塑性斷裂是引起復(fù)合材料橫向壓縮破壞的主要因素,而界面損傷只是在基體出現(xiàn)損傷的局部區(qū)域內(nèi)由于應(yīng)力集中引起的。

由上述分析可知,復(fù)合材料在實(shí)際壓縮破壞過(guò)程中,界面損傷不是導(dǎo)致材料破壞主要因素,此時(shí)繼續(xù)增大界面強(qiáng)度仿真得到的復(fù)合材料橫向壓縮強(qiáng)度基本不再升高。因此本文采用一種忽略界面單元的有限元模型,用此簡(jiǎn)化的無(wú)界面單元模型仿真得到的復(fù)合材料橫向壓縮應(yīng)力-應(yīng)變曲線(xiàn)以及損傷破壞形貌如圖5所示?？梢钥闯?無(wú)界面單元模型預(yù)測(cè)的復(fù)合材料橫向壓縮強(qiáng)度為142.9 MPa,與試驗(yàn)值最為吻合,誤差僅為1.4%。同時(shí)在橫向壓縮載荷作用下復(fù)合材料最終斷裂形貌與界面強(qiáng)度為140 MPa時(shí)的情況相近。這說(shuō)明了通過(guò)準(zhǔn)確描述基體的塑性損傷退化行為,無(wú)界面單元可以更好地預(yù)測(cè)單向復(fù)合材料的橫向壓縮力學(xué)性能,避免了數(shù)值仿真中界面參數(shù)選取的不確定性問(wèn)題。

圖5 無(wú)界面單元模型預(yù)測(cè)的復(fù)合材料橫向壓縮應(yīng)力-應(yīng)變曲線(xiàn)及破壞形貌

3 結(jié) 論

1) 本文建立的纖維隨機(jī)分布代表性體積單元微觀力學(xué)模型可以準(zhǔn)確預(yù)測(cè)單向復(fù)合材料的力學(xué)性能與損傷形貌。

2) 纖維隨機(jī)分布會(huì)引起復(fù)合材料損傷初始位置與斷裂路徑改變,是引起復(fù)合材料橫向壓縮強(qiáng)度差異性的一個(gè)因素。

3) 基體塑性損傷而不是界面損傷在復(fù)合材料橫向壓縮破壞過(guò)程中起主導(dǎo)作用。通過(guò)準(zhǔn)確描述基體的塑性損傷退化行為,無(wú)界面單元可以很好地預(yù)測(cè)單向復(fù)合材料的橫向壓縮力學(xué)性能,避免了數(shù)值仿真中界面參數(shù)選取的不確定性問(wèn)題。