蔡登安,周光明,曹 然,黃 翔
(南京航空航天大學(xué)機(jī)械結(jié)構(gòu)力學(xué)及控制國家重點(diǎn)實(shí)驗(yàn)室,南京 210016)
纖維增強(qiáng)復(fù)合材料以其輕質(zhì)、高強(qiáng)、抗疲勞、減振、耐高溫及可設(shè)計(jì)等一系列優(yōu)點(diǎn),在航空航天、能源、交通、機(jī)械和建筑等部門日益獲得了廣泛應(yīng)用。一般而言,構(gòu)件中復(fù)合材料多受雙向或多向應(yīng)力狀態(tài)[1]。復(fù)合材料的測試技術(shù)已能準(zhǔn)確表征纖維增強(qiáng)復(fù)合材料的單向力學(xué)行為。然而,同樣的試驗(yàn)方法并不能完全適用于復(fù)合材料多向甚至只是雙向響應(yīng)問題[2]。對于受復(fù)雜應(yīng)力狀態(tài)的復(fù)合材料,僅僅研究其單向載荷下的強(qiáng)度是遠(yuǎn)遠(yuǎn)不夠的。因此,不管是對強(qiáng)度準(zhǔn)則進(jìn)行試驗(yàn)驗(yàn)證,還是從損傷角度探討復(fù)合材料在復(fù)雜應(yīng)力狀態(tài)下的破壞,都需要對復(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),邊界元法針對一些問題的分析有著獨(dú)特的優(yōu)點(diǎn)。Lachat和Watson[7]為邊界元法作為一種有效數(shù)值方法的研究作出了貢獻(xiàn):他們提出了一種類似于有限元法的等參方法,并證明邊界元法可作為解決具有復(fù)雜構(gòu)型問題的一種有效工具。借助邊界元方法,Pineda León等[8]研究了在不同溫度和應(yīng)力狀態(tài)下材料的蠕變分析模型;劉云忠等[9-10]提出了正交各向異性復(fù)合材料裂紋問題的邊界元分析方法;商欣萍等[11]用邊界元法從理論上分析非織造土工織物的拉伸性能,與實(shí)驗(yàn)結(jié)果有較好的一致性;Araújo等[12]針對碳納米管增強(qiáng)復(fù)合材料的微觀結(jié)構(gòu)分析,編寫了邊界元并行計(jì)算程序,并據(jù)此提出了一種可擴(kuò)展的復(fù)合材料邊界元并行計(jì)算代碼;基于線性夾雜模型,Wang等[13]將邊界元方法應(yīng)用于纖維增強(qiáng)復(fù)合材料的3D大尺度熱分析問題的研究中。
針對十字型試樣雙向拉伸實(shí)驗(yàn)研究,本文采用邊界元法,從理論上研究了交織玻璃纖維增強(qiáng)復(fù)合材料的雙向拉伸強(qiáng)度。建立復(fù)合材料雙向拉伸邊界元分析模型,并與有限元結(jié)果及試驗(yàn)結(jié)果進(jìn)行比較,以期為復(fù)合材料雙向拉伸強(qiáng)度研究提供一些理論依據(jù)。
與長、寬尺寸相比,交織玻璃纖維增強(qiáng)復(fù)合材料十字型試樣的厚度很小,故可將其看作是平面正交各向異性材料。平面正交各向異性問題的基本方程可描述如下:
平衡方程:
幾何方程:
物理方程:
或
式中 [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)力邊界條件
位移邊界條件
邊界元法的基礎(chǔ)是用Fredholm[14]積分方程的位勢積分來表示問題的解,它是將支配物理現(xiàn)象的微分方程轉(zhuǎn)化成邊界上的積分方程,然后再進(jìn)行離散求解的數(shù)值計(jì)算方法。
圖1 平面區(qū)域V和邊界SFig.1 Planar region V and boundary S
對于平面正交各向異性材料,如圖1所示,由邊界S所包圍的平面區(qū)域上的任意點(diǎn)x的位移分量ui(x),在不考慮體積力時(shí),可由邊界上各點(diǎn)y的位移分量和面力分量通過以下邊界積分方程求得[5,11]:
式中 ui和ti分別為位移和面力;Uij(x,y)和Tij(x,y)分別為平面正交各向異性材料的位移基本解分量和面力基本解分量。
系數(shù)Cij(x)的值與x點(diǎn)處邊界形狀有關(guān),邊界S光滑時(shí),有
對于平面正交各向異性材料,位移和面力基本解可表示為[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分別為邊界外法線的方向余弦。
通過邊界積分方程,可求出邊界上所有未知的位移和表面力分量。然而,通常情況下很難用解析方法求解積分方程,只能采用數(shù)值解法。將區(qū)域的邊界劃分成N(N=N1+N2)個(gè)簡單單元。其中,N1個(gè)邊界單元屬于表面力已知,N2個(gè)邊界單元屬于位移已知,以N個(gè)邊界單元上的積分之和表示整個(gè)邊界上的積分。各邊界單元上的位移和表面力,可通過插值函數(shù)與邊界單元上有限個(gè)節(jié)點(diǎn)的位移和表面力以多項(xiàng)式近似。本文采用常單元將邊界積分方程離散,單元節(jié)點(diǎn)取為邊界單元中點(diǎn),邊界單元的位移和表面力等于單元節(jié)點(diǎn)的位移和表面力。
設(shè) yi(i=1,2,…,N)為邊界 Sj上的任意一點(diǎn),對任意節(jié)點(diǎn)xi,邊界積分方程可離散為
式中 下標(biāo)i,j分別為第i、第j個(gè)單元。
設(shè):
則對任意節(jié)點(diǎn),邊界積分方程可離散為
對于N個(gè)節(jié)點(diǎn),根據(jù)各節(jié)點(diǎn)的離散化方程,可得到包含2N個(gè)方程的一次方程組,用矩陣形式可表示為
式中 [H]、[G]均為2N×2N階的位移和面力系數(shù)矩陣;{U}、{T}分別為邊界單元節(jié)點(diǎn)的位移分量和表面力分量。
根據(jù)以上理論分析,利用FORTRAN語言編寫計(jì)算程序,來分析交織玻璃纖維增強(qiáng)復(fù)合材料雙向拉伸強(qiáng)度,程序的計(jì)算流程見圖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)過程中,試樣的AB端和CD端由夾頭夾持在拉力作用下移動,EF端和GH端由夾頭固定。與試驗(yàn)不同的是邊界元理論計(jì)算時(shí),為避免附加彎矩的影響,需釋放EF端在1方向的位移約束和GH端在2方向的位移約束,如圖4所示。若分別給AB端和CD端一個(gè)2方向和1方向上的伸長量,通過邊界元法計(jì)算程序,便可得出十字型試樣各自由邊上邊界單元的位移、各內(nèi)點(diǎn)的位移和應(yīng)力,以及達(dá)到該伸長量所施加的雙向拉伸載荷。
對十字型試樣進(jìn)行單元劃分,由于試件完全對稱,為減少計(jì)算量,只需對1/4試樣劃分單元,如圖5所示。
圖3 十字型試樣幾何形狀及尺寸Fig.3 Geometry and dimensions of cruciform specimen
圖4 十字型試樣的載荷與約束Fig.4 Loads and constraints of cruciform specimen
圖5 十字型試樣邊界單元劃分Fig.5 Boundary element of cruciform specimen
試樣的OI邊單元約束1方向的位移,OJ邊單元約束2方向的位移;于IB與JC端作用原十字型試樣拉伸載荷的一半。圖中標(biāo)記的各點(diǎn)為指定的內(nèi)點(diǎn)。通過計(jì)算各指定內(nèi)點(diǎn)的應(yīng)力,可了解十字型試樣在雙向拉伸過程中的內(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)則方程可簡化為
式中 Xt、Xc、Yt、Yc分別為材料縱向的拉壓強(qiáng)度和橫向的拉壓強(qiáng)度,可由單向拉伸(壓縮)試驗(yàn)測得;強(qiáng)度參數(shù)F12由試驗(yàn)獲得雙向強(qiáng)度后計(jì)算得到。
邊界元程序計(jì)算過程中,若中心試驗(yàn)區(qū)的縱向和橫向平均應(yīng)力使準(zhǔn)則方程左邊大于等于1,則判定十字型試樣斷裂。
在雙向拉伸試驗(yàn)之前,對交織玻璃纖維增強(qiáng)復(fù)合材料進(jìn)行基本力學(xué)性能測試,從而獲得交織玻璃纖維增強(qiáng)復(fù)合材料的工程彈性常數(shù)和縱、橫向拉壓強(qiáng)度,所得材料特性見表1。本文研究的復(fù)合材料在縱、橫向具有相同的力學(xué)性能,增強(qiáng)相為交織玻璃纖維織物,基體為WSR618型環(huán)氧樹脂。
表1 材料特性Table 1 Material properties
雙向拉伸試驗(yàn)在SDS100電液伺服動靜試驗(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í),對十字型試樣的1/4模型,建立了有限元分析模型,獲得其有限元分析結(jié)果。表2給出了邊界元模型與有限元模型的比較。
表2 邊界元與有限元分析模型的比較Table 2 Comparison for the models of BEM and FEM
由表2可見,針對同一個(gè)分析模型,邊界元法在單元數(shù)目的劃分上要遠(yuǎn)少于有限元法。因此,其模型的自由度數(shù)目也存在較大懸殊,邊界元法求解方程數(shù)目明顯較少,從而節(jié)約占機(jī)內(nèi)存及計(jì)算時(shí)間。由此可見,邊界元法在計(jì)算效率上比有限元法更具優(yōu)勢。
交織玻璃纖維增強(qiáng)復(fù)合材料的雙向拉伸縱向應(yīng)力-應(yīng)變曲線如圖6所示。從圖6中可看出,不同載荷比下,材料的雙向拉伸應(yīng)力、應(yīng)變保持著線性關(guān)系。邊界元法、有限元法計(jì)算結(jié)果與試驗(yàn)結(jié)果趨勢一致,且較為接近。說明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é)果比較見表3。分析表3中數(shù)據(jù)可知,十字型試樣邊界元模型與有限元模型計(jì)算所得的雙向拉伸強(qiáng)度結(jié)果較接近;與試驗(yàn)值比較,不同載荷比下2種計(jì)算方法的誤差范圍分別為5.9% ~6.8%和6.4% ~7.3%,且均為正誤差??烧J(rèn)為2種數(shù)值方法在計(jì)算精度上相當(dāng),2種方法的計(jì)算值均比試驗(yàn)值偏大。誤差來源可能有2種:一是實(shí)際試樣存在一定的缺陷和細(xì)微損傷,造成拉伸強(qiáng)度的降低;二是理論分析模型與實(shí)際情況不完全吻合。
同時(shí),由表3中結(jié)果可見,交織玻璃纖維增強(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%。由此可見,縱、橫向載荷比越小,該材料的雙向弱化效應(yīng)越明顯,雙向加載為等雙拉時(shí),雙向拉伸強(qiáng)度最小,材料抵抗雙向拉伸破壞的能力最弱。
表3 拉伸強(qiáng)度結(jié)果及比較Table 3 Results and comparison of tensile strength
十字型試樣雙向拉伸的破壞形式如圖7所示,破壞屬于脆性斷裂。
由圖7可見,不同載荷比下,交織玻璃纖維增強(qiáng)復(fù)合材料十字型試樣的破壞形式不同∶載荷比為1∶1時(shí),試樣沿試驗(yàn)區(qū)對角線斷裂;載荷比為2∶1和3∶1時(shí),試樣裂紋產(chǎn)生及斷裂方向與試驗(yàn)區(qū)對角線約成10°~15°夾角;載荷比為4∶1時(shí),試樣出現(xiàn)明顯的分層破壞,且分層方向偏向于載荷大的加載方向。
(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ù)合材料中裂紋問題的邊界元法分析[J].固體火箭技術(shù),1995,18(2):52-61.
[10] 劉云忠.邊界元法計(jì)算厚板表面裂紋問題的應(yīng)力強(qiáng)度因子[J].固體火箭技術(shù),1996,19(3):64-70.
[11] 商欣萍,儲才元.用邊界元法分析試樣寬度對非織造土工織物拉伸性能的影響[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.