蔡登安，周光明，曹 然，黃 翔

(南京航空航天大學(xué)機(jī)械結(jié)構(gòu)力學(xué)及控制國(guó)家重點(diǎn)實(shí)驗(yàn)室，南京 210016)

0 引言

纖維增強(qiáng)復(fù)合材料以其輕質(zhì)、高強(qiáng)、抗疲勞、減振、耐高溫及可設(shè)計(jì)等一系列優(yōu)點(diǎn)，在航空航天、能源、交通、機(jī)械和建筑等部門日益獲得了廣泛應(yīng)用。一般而言，構(gòu)件中復(fù)合材料多受雙向或多向應(yīng)力狀態(tài)［1］。復(fù)合材料的測(cè)試技術(shù)已能準(zhǔn)確表征纖維增強(qiáng)復(fù)合材料的單向力學(xué)行為。然而，同樣的試驗(yàn)方法并不能完全適用于復(fù)合材料多向甚至只是雙向響應(yīng)問(wèn)題［2］。對(duì)于受復(fù)雜應(yīng)力狀態(tài)的復(fù)合材料，僅僅研究其單向載荷下的強(qiáng)度是遠(yuǎn)遠(yuǎn)不夠的。因此，不管是對(duì)強(qiáng)度準(zhǔn)則進(jìn)行試驗(yàn)驗(yàn)證，還是從損傷角度探討復(fù)合材料在復(fù)雜應(yīng)力狀態(tài)下的破壞，都需要對(duì)復(fù)合材料試樣進(jìn)行雙向加載研究。獲得雙向應(yīng)力狀態(tài)的傳統(tǒng)方式是采用圓管形試件，目前更多的研究者把焦點(diǎn)放到了十字型試樣上［3－4］，其優(yōu)點(diǎn)在于十字型試樣雙向拉伸試驗(yàn)?zāi)苤庇^地反映板殼復(fù)合材料的雙向受力狀態(tài)，因而成為當(dāng)下最受重視的一種雙向拉伸試驗(yàn)方法。

邊界元法(BEM)是一種將研究的連續(xù)體區(qū)域的邊界劃分成有限個(gè)單元的數(shù)值方法，在某些領(lǐng)域已經(jīng)成為工程設(shè)計(jì)的重要工具［5－6］。相比有限元法(FEM)，邊界元法針對(duì)一些問(wèn)題的分析有著獨(dú)特的優(yōu)點(diǎn)。Lachat和Watson［7］為邊界元法作為一種有效數(shù)值方法的研究作出了貢獻(xiàn):他們提出了一種類似于有限元法的等參方法，并證明邊界元法可作為解決具有復(fù)雜構(gòu)型問(wèn)題的一種有效工具。借助邊界元方法，Pineda León等［8］研究了在不同溫度和應(yīng)力狀態(tài)下材料的蠕變分析模型;劉云忠等［9－10］提出了正交各向異性復(fù)合材料裂紋問(wèn)題的邊界元分析方法;商欣萍等［11］用邊界元法從理論上分析非織造土工織物的拉伸性能，與實(shí)驗(yàn)結(jié)果有較好的一致性;Araújo等［12］針對(duì)碳納米管增強(qiáng)復(fù)合材料的微觀結(jié)構(gòu)分析，編寫了邊界元并行計(jì)算程序，并據(jù)此提出了一種可擴(kuò)展的復(fù)合材料邊界元并行計(jì)算代碼;基于線性?shī)A雜模型，Wang等［13］將邊界元方法應(yīng)用于纖維增強(qiáng)復(fù)合材料的3D大尺度熱分析問(wèn)題的研究中。

針對(duì)十字型試樣雙向拉伸實(shí)驗(yàn)研究，本文采用邊界元法，從理論上研究了交織玻璃纖維增強(qiáng)復(fù)合材料的雙向拉伸強(qiáng)度。建立復(fù)合材料雙向拉伸邊界元分析模型，并與有限元結(jié)果及試驗(yàn)結(jié)果進(jìn)行比較，以期為復(fù)合材料雙向拉伸強(qiáng)度研究提供一些理論依據(jù)。

1 邊界元模型

1．1 基本方程

與長(zhǎng)、寬尺寸相比，交織玻璃纖維增強(qiáng)復(fù)合材料十字型試樣的厚度很小，故可將其看作是平面正交各向異性材料。平面正交各向異性問(wèn)題的基本方程可描述如下:

平衡方程:

幾何方程:

物理方程:

或

式中 ［S］和［Q］分別表示柔度矩陣和剛度矩陣，且［Q］=［S］－1。

平面正交各向異性材料有4個(gè)獨(dú)立的彈性常數(shù)，即 E1、E2、ν12和 G12，分別表示交織玻璃纖維增強(qiáng)復(fù)合材料2個(gè)方向拉伸時(shí)的彈性模量、泊松比和剪切模量。交織玻璃纖維增強(qiáng)復(fù)合材料的柔度矩陣［S］可表示為

邊界條件:

應(yīng)力邊界條件

位移邊界條件

1．2 邊界積分方程與基本解

邊界元法的基礎(chǔ)是用Fredholm［14］積分方程的位勢(shì)積分來(lái)表示問(wèn)題的解，它是將支配物理現(xiàn)象的微分方程轉(zhuǎn)化成邊界上的積分方程，然后再進(jìn)行離散求解的數(shù)值計(jì)算方法。

圖1 平面區(qū)域V和邊界SFig．1 Planar region V and boundary S

對(duì)于平面正交各向異性材料，如圖1所示，由邊界S所包圍的平面區(qū)域上的任意點(diǎn)x的位移分量ui(x)，在不考慮體積力時(shí)，可由邊界上各點(diǎn)y的位移分量和面力分量通過(guò)以下邊界積分方程求得［5，11］:

式中 ui和ti分別為位移和面力;Uij(x，y)和Tij(x，y)分別為平面正交各向異性材料的位移基本解分量和面力基本解分量。

系數(shù)Cij(x)的值與x點(diǎn)處邊界形狀有關(guān)，邊界S光滑時(shí)，有

對(duì)于平面正交各向異性材料，位移和面力基本解可表示為［15］

其中，Kα、Ai、αi、Mi(i=1，2)為與彈性常數(shù)有關(guān)的參數(shù);ri和θi分別為復(fù)平面中的矢徑和幅角，且有

式中 l1、l2分別為y點(diǎn)與x點(diǎn)橫、縱坐標(biāo)之差;n1、n2分別為邊界外法線的方向余弦。

1．3 積分方程離散化

通過(guò)邊界積分方程，可求出邊界上所有未知的位移和表面力分量。然而，通常情況下很難用解析方法求解積分方程，只能采用數(shù)值解法。將區(qū)域的邊界劃分成N(N=N1+N2)個(gè)簡(jiǎn)單單元。其中，N1個(gè)邊界單元屬于表面力已知，N2個(gè)邊界單元屬于位移已知，以N個(gè)邊界單元上的積分之和表示整個(gè)邊界上的積分。各邊界單元上的位移和表面力，可通過(guò)插值函數(shù)與邊界單元上有限個(gè)節(jié)點(diǎn)的位移和表面力以多項(xiàng)式近似。本文采用常單元將邊界積分方程離散，單元節(jié)點(diǎn)取為邊界單元中點(diǎn)，邊界單元的位移和表面力等于單元節(jié)點(diǎn)的位移和表面力。

設(shè) yi(i=1，2，…，N)為邊界 Sj上的任意一點(diǎn)，對(duì)任意節(jié)點(diǎn)xi，邊界積分方程可離散為

式中 下標(biāo)i，j分別為第i、第j個(gè)單元。

設(shè):

則對(duì)任意節(jié)點(diǎn)，邊界積分方程可離散為

對(duì)于N個(gè)節(jié)點(diǎn)，根據(jù)各節(jié)點(diǎn)的離散化方程，可得到包含2N個(gè)方程的一次方程組，用矩陣形式可表示為

式中 ［H］、［G］均為2N×2N階的位移和面力系數(shù)矩陣;{U}、{T}分別為邊界單元節(jié)點(diǎn)的位移分量和表面力分量。

2 計(jì)算方法

根據(jù)以上理論分析，利用FORTRAN語(yǔ)言編寫計(jì)算程序，來(lái)分析交織玻璃纖維增強(qiáng)復(fù)合材料雙向拉伸強(qiáng)度，程序的計(jì)算流程見(jiàn)圖2。計(jì)算前，要輸入試樣4個(gè)獨(dú)立的彈性常數(shù)(縱向和橫向的彈性模量、泊松比及剪切模量);同時(shí)，要輸入試樣的強(qiáng)度參數(shù)、幾何尺寸、各單元的坐標(biāo)以及邊界條件。

圖2 程序計(jì)算流程圖Fig．2 Flow-process diagram of program

交織玻璃纖維增強(qiáng)復(fù)合材料十字型試樣的幾何尺寸(單位:mm)，如圖3所示。十字型試樣中心36 mm×36 mm的正方形區(qū)域?yàn)橹行脑囼?yàn)區(qū)，該區(qū)域滿足應(yīng)力分布均勻，且應(yīng)力水平較高，以保證初始破壞發(fā)生在中心試驗(yàn)區(qū);同時(shí)，試驗(yàn)區(qū)剪應(yīng)力遠(yuǎn)小于正應(yīng)力，以至可忽略不計(jì)。受雙向拉伸時(shí)，十字型試樣在實(shí)際拉伸試驗(yàn)過(guò)程中，試樣的AB端和CD端由夾頭夾持在拉力作用下移動(dòng)，EF端和GH端由夾頭固定。與試驗(yàn)不同的是邊界元理論計(jì)算時(shí)，為避免附加彎矩的影響，需釋放EF端在1方向的位移約束和GH端在2方向的位移約束，如圖4所示。若分別給AB端和CD端一個(gè)2方向和1方向上的伸長(zhǎng)量，通過(guò)邊界元法計(jì)算程序，便可得出十字型試樣各自由邊上邊界單元的位移、各內(nèi)點(diǎn)的位移和應(yīng)力，以及達(dá)到該伸長(zhǎng)量所施加的雙向拉伸載荷。

對(duì)十字型試樣進(jìn)行單元?jiǎng)澐?，由于試件完全?duì)稱，為減少計(jì)算量，只需對(duì)1/4試樣劃分單元，如圖5所示。

圖3 十字型試樣幾何形狀及尺寸Fig．3 Geometry and dimensions of cruciform specimen

圖4 十字型試樣的載荷與約束Fig．4 Loads and constraints of cruciform specimen

圖5 十字型試樣邊界單元?jiǎng)澐諪ig．5 Boundary element of cruciform specimen

試樣的OI邊單元約束1方向的位移，OJ邊單元約束2方向的位移;于IB與JC端作用原十字型試樣拉伸載荷的一半。圖中標(biāo)記的各點(diǎn)為指定的內(nèi)點(diǎn)。通過(guò)計(jì)算各指定內(nèi)點(diǎn)的應(yīng)力，可了解十字型試樣在雙向拉伸過(guò)程中的內(nèi)力分布情況。

要計(jì)算交織玻璃纖維增強(qiáng)復(fù)合材料雙向拉伸強(qiáng)度，需確定十字型試樣斷裂的條件，即試件在雙向載荷下，拉伸到何種程度達(dá)到斷裂。本文采用的強(qiáng)度準(zhǔn)則為Tsai-Wu張量準(zhǔn)則。平面應(yīng)力狀態(tài)下，Tsai-Wu張量準(zhǔn)則可表述為

十字型試樣中心試驗(yàn)區(qū)的剪應(yīng)力較小，忽略不計(jì)，因而準(zhǔn)則方程可簡(jiǎn)化為

式中 Xt、Xc、Yt、Yc分別為材料縱向的拉壓強(qiáng)度和橫向的拉壓強(qiáng)度，可由單向拉伸(壓縮)試驗(yàn)測(cè)得;強(qiáng)度參數(shù)F12由試驗(yàn)獲得雙向強(qiáng)度后計(jì)算得到。

邊界元程序計(jì)算過(guò)程中，若中心試驗(yàn)區(qū)的縱向和橫向平均應(yīng)力使準(zhǔn)則方程左邊大于等于1，則判定十字型試樣斷裂。

3 實(shí)驗(yàn)驗(yàn)證

3．1 基本力學(xué)性能測(cè)試

在雙向拉伸試驗(yàn)之前，對(duì)交織玻璃纖維增強(qiáng)復(fù)合材料進(jìn)行基本力學(xué)性能測(cè)試，從而獲得交織玻璃纖維增強(qiáng)復(fù)合材料的工程彈性常數(shù)和縱、橫向拉壓強(qiáng)度，所得材料特性見(jiàn)表1。本文研究的復(fù)合材料在縱、橫向具有相同的力學(xué)性能，增強(qiáng)相為交織玻璃纖維織物，基體為WSR618型環(huán)氧樹(shù)脂。

表1 材料特性Table 1 Material properties

3．2 試驗(yàn)結(jié)果與分析

雙向拉伸試驗(yàn)在SDS100電液伺服動(dòng)靜試驗(yàn)機(jī)上進(jìn)行。主要進(jìn)行4個(gè)載荷比(f=縱向載荷∶橫向載荷=1∶1，2∶1，3∶1，4∶1)的雙向拉伸試驗(yàn)，試驗(yàn)采用位移控制方式，按f×0．5mm/min的加載速率加載。

同時(shí)，對(duì)十字型試樣的1/4模型，建立了有限元分析模型，獲得其有限元分析結(jié)果。表2給出了邊界元模型與有限元模型的比較。

表2 邊界元與有限元分析模型的比較Table 2 Comparison for the models of BEM and FEM

由表2可見(jiàn)，針對(duì)同一個(gè)分析模型，邊界元法在單元數(shù)目的劃分上要遠(yuǎn)少于有限元法。因此，其模型的自由度數(shù)目也存在較大懸殊，邊界元法求解方程數(shù)目明顯較少，從而節(jié)約占機(jī)內(nèi)存及計(jì)算時(shí)間。由此可見(jiàn)，邊界元法在計(jì)算效率上比有限元法更具優(yōu)勢(shì)。

交織玻璃纖維增強(qiáng)復(fù)合材料的雙向拉伸縱向應(yīng)力-應(yīng)變曲線如圖6所示。從圖6中可看出，不同載荷比下，材料的雙向拉伸應(yīng)力、應(yīng)變保持著線性關(guān)系。邊界元法、有限元法計(jì)算結(jié)果與試驗(yàn)結(jié)果趨勢(shì)一致，且較為接近。說(shuō)明2種數(shù)值方法的模擬結(jié)果都是比較可信的。因此，可利用邊界元法，從理論上分析復(fù)合材料雙向拉伸力學(xué)性能。

圖6 材料雙向拉伸應(yīng)力-應(yīng)變曲線(縱向)Fig．6 Stress-strain curve of the material under biaxial tensile(longitudinal direction)

交織玻璃纖維增強(qiáng)復(fù)合材料不同載荷比的拉伸強(qiáng)度結(jié)果比較見(jiàn)表3。分析表3中數(shù)據(jù)可知，十字型試樣邊界元模型與有限元模型計(jì)算所得的雙向拉伸強(qiáng)度結(jié)果較接近;與試驗(yàn)值比較，不同載荷比下2種計(jì)算方法的誤差范圍分別為5．9% ～6．8%和6．4% ～7．3%，且均為正誤差。可認(rèn)為2種數(shù)值方法在計(jì)算精度上相當(dāng)，2種方法的計(jì)算值均比試驗(yàn)值偏大。誤差來(lái)源可能有2種:一是實(shí)際試樣存在一定的缺陷和細(xì)微損傷，造成拉伸強(qiáng)度的降低;二是理論分析模型與實(shí)際情況不完全吻合。

同時(shí)，由表3中結(jié)果可見(jiàn)，交織玻璃纖維增強(qiáng)復(fù)合材料的雙向拉伸行為存在明顯的雙向弱化現(xiàn)象:在載荷比為1∶1的加載條件下，材料的雙向拉伸強(qiáng)度(縱向強(qiáng)度)僅為單向拉伸強(qiáng)度的60．5%;載荷比為2∶1時(shí)，雙向拉伸強(qiáng)度為單向拉伸強(qiáng)度的71．5%;載荷比為3∶1和4∶1時(shí)，雙向拉伸強(qiáng)度分別為單向拉伸強(qiáng)度的77．5%和81．2%。由此可見(jiàn)，縱、橫向載荷比越小，該材料的雙向弱化效應(yīng)越明顯，雙向加載為等雙拉時(shí)，雙向拉伸強(qiáng)度最小，材料抵抗雙向拉伸破壞的能力最弱。

表3 拉伸強(qiáng)度結(jié)果及比較Table 3 Results and comparison of tensile strength

十字型試樣雙向拉伸的破壞形式如圖7所示，破壞屬于脆性斷裂。

由圖7可見(jiàn)，不同載荷比下，交織玻璃纖維增強(qiáng)復(fù)合材料十字型試樣的破壞形式不同∶載荷比為1∶1時(shí)，試樣沿試驗(yàn)區(qū)對(duì)角線斷裂;載荷比為2∶1和3∶1時(shí)，試樣裂紋產(chǎn)生及斷裂方向與試驗(yàn)區(qū)對(duì)角線約成10°～15°夾角;載荷比為4∶1時(shí)，試樣出現(xiàn)明顯的分層破壞，且分層方向偏向于載荷大的加載方向。

4 結(jié)論

(1)建立了交織玻璃纖維增強(qiáng)復(fù)合材料雙向拉伸的邊界元理論分析模型，理論分析結(jié)果與試驗(yàn)結(jié)果具有較好的一致性。因此，可利用邊界元法，從理論上分析該類復(fù)合材料的雙向拉伸性能。

圖7 十字型試樣破壞形貌Fig．7 Fracture morphology of cruciform specimen

(2)在計(jì)算交織玻璃纖維增強(qiáng)復(fù)合材料雙向拉伸強(qiáng)度時(shí)，邊界元法在計(jì)算效率上優(yōu)于有限元法，計(jì)算精度上兩者相當(dāng);同時(shí)，2種方法的計(jì)算值均比試驗(yàn)值偏大。

(3)交織玻璃纖維增強(qiáng)復(fù)合材料的雙向拉伸強(qiáng)度，與其單向拉伸強(qiáng)度相比，具有明顯的雙向弱化效應(yīng)?？v、橫向載荷比越小，材料的雙向弱化現(xiàn)象越明顯，等雙拉時(shí)，材料的雙向拉伸強(qiáng)度(縱向強(qiáng)度)僅為單向拉伸強(qiáng)度的60．5%。

(4)不同載荷比下，交織玻璃纖維增強(qiáng)復(fù)合材料雙向拉伸的破壞形式有所不同，主要表現(xiàn)為試驗(yàn)區(qū)纖維斷裂方向的不同，載荷比為4∶1時(shí)出現(xiàn)了分層破壞。

［1］ Lin W P，Hu H T．Parametric study on failure of fiber-reinforced composite laminates under biaxial tensile load［J］．Journal of Composite Materials，2002，36(12):1481-1503．

［2］ Smits A，Van Hemelrijck D，Philippidis T P，et al．Design of a cruciform specimen for biaxial testing of fibre reinforced composite laminates［J］．Composites Science and Technology，2006，66(7-8):964-975．

［3］ Makris A，Vandenbergh T，Ramault C，et al．Shape optimisation of a biaxially loaded cruciform specimen［J］．Polymer Testing，2010，29(2):216-223．

［4］ Serna Moreno M C，López Cela J J．Failure envelope under biaxial tensile loading for chopped glass-reinforced polyester composites［J］．Composites Science and Technology，2011，72(1):91-96．

［5］ Liu Y J．A new fast multipole boundary element method for solving large-scale two-dimensionalelastostatic problems［J］．International Journal Numerical Methods in Engineering，2006，65(6):863-881．

［6］ 陳龍，朱興一．求解具有弱界面的復(fù)合材料有效性能的近似模型與邊界元法［J］．復(fù)合材料學(xué)報(bào)，2013，30(2):137-143．

［7］ Lachat J C，Watson J O．Effective numerical treatment of boundary integral equations:A formulation for three‐dimensional elastostatics［J］．International Journal for Numerical Methods in Engineering，1976，10(5):991-1005．

［8］ Pineda León E，Samayoa Ochoa D，Rodríguez-Castellanos A，et al．Creeping analysis with variable temperature applying the boundary element method［J］．Engineering Analysis with Boundary Elements，2012，36(12):1715-1720．

［9］ 劉云忠，雷海峰．正交各向異性復(fù)合材料中裂紋問(wèn)題的邊界元法分析［J］．固體火箭技術(shù)，1995，18(2):52-61．

［10］ 劉云忠．邊界元法計(jì)算厚板表面裂紋問(wèn)題的應(yīng)力強(qiáng)度因子［J］．固體火箭技術(shù)，1996，19(3):64-70．

［11］ 商欣萍，儲(chǔ)才元．用邊界元法分析試樣寬度對(duì)非織造土工織物拉伸性能的影響［J］．東華大學(xué)學(xué)報(bào)(自然科學(xué)版)，2004，30(1):1-5．

［12］ Araújo F C，d＇Azevedo E F，Gray L J．Boundary-element parallel-computing algorithm for the microstructural analysis of general composites［J］．Computers and Structures，2010，88(11):773-784．

［13］ Wang H T，Yao Z H．Large-scale thermal analysis of fiber composites using a line-inclusion model by the fast boundary element method［J］．Engineering Analysis with Boundary Elements，2013，37(2):319-326．

［14］ Fredholm I．Sur une classe d＇équations fonctionnelles［J］．Acta Mathematica，1903，27(1):365-390．

［15］ Rizzo F J，Shippy D J．A method for stress determination in plane anisotropic elastic bodies［J］．Journal of Composite Materials，1970，4(1):36-61．