含中心孔復(fù)合材料變剛度板孔邊應(yīng)力解析法分析

李 陽，牛雪娟，2，潘文峰

(1.天津工業(yè)大學(xué) 機械工程學(xué)院，天津 300387；2.天津市現(xiàn)代機電裝備技術(shù)重點實驗室，天津 300387)

0 引言

炭纖維復(fù)合材料具有優(yōu)異的力學(xué)性能及可設(shè)計性等優(yōu)點，被廣泛用于船舶、建筑、醫(yī)學(xué)等領(lǐng)域，日益受到人們的關(guān)注。在實際應(yīng)用中，由于結(jié)構(gòu)或功能方面的需求，常需要在復(fù)合材料層合板上開孔或槽，這將導(dǎo)致在孔或槽周圍產(chǎn)生應(yīng)力集中，降低層合板的承載能力，形成潛在威脅。在發(fā)生應(yīng)力集中的情況下，孔邊應(yīng)力大小是研究者非常關(guān)心的內(nèi)容。

Nageswara等[1]根據(jù)薩文對帶孔各向異性板的研究，提出一種能夠計算具有對稱結(jié)構(gòu)的層合板及各向同性板孔邊應(yīng)力分布的方法，該方法還能夠基于一些準(zhǔn)則確定失效強度。Kaltackle等[2]及Sharma等[3]利用不同方法研究了單層板的不同鋪層角度、不同材料特性對孔邊應(yīng)力的影響。William[4]研究了復(fù)合材料層合板在不同載荷方向下孔邊的應(yīng)力集中情況。Russo等[5]求出在單軸拉伸載荷作用下，復(fù)合材料層合板孔邊的應(yīng)力分布，這種計算方法包含一定的修正功能。Michele等[6]基于已有的正交各向異性彈性理論，提出能夠計算槽邊應(yīng)力集中系數(shù)的工程公式，并利用有限元方法及已有文獻驗證其準(zhǔn)確性。國內(nèi)，韓小平等[7]推導(dǎo)出能夠計算層合板孔邊應(yīng)力及應(yīng)力集中系數(shù)的公式，根據(jù)公式分別討論了材料彈性參數(shù)、鋪層結(jié)構(gòu)等因素對孔邊應(yīng)力集中的影響。李成等[8]用積分方法對含孔層合板進行應(yīng)力分析，針對不同載荷狀況對孔邊的應(yīng)力集中情況進行探討。謝慈航等[9]采用各向彈性理論及保角變換法，總結(jié)出一個包含孔形參數(shù)的表達式，利用該表達式求出復(fù)合材料層合板在一定受力情況下的孔邊應(yīng)力，并對不同曲率孔形的應(yīng)力集中系數(shù)進行了討論。

上述國內(nèi)外學(xué)者所做解析法研究以及仿真研究和實驗研究，都是利用傳統(tǒng)的直線鋪放含孔復(fù)合材料層合板進行的，隨著復(fù)合材料層合板制備工藝的提高，為更加充分地利用炭纖維的可設(shè)計性及其他特點，對于改變纖維角度的這種變剛度[10]層合板，國內(nèi)外學(xué)者進行了相關(guān)研究。Gürdal等[11-12]采用角度沿幾何軸線性變化的曲線來定義鋪放參考軌跡，并采用先進的鋪放技術(shù)制備了變剛度層合板；然后利用實驗證明了變剛度復(fù)合材料層合板其承載能力高于直線鋪放的層合板。ZHU Yingdan等[13]采用分段函數(shù)擬合最大主應(yīng)力方向，通過優(yōu)化目標(biāo)函數(shù)獲得構(gòu)造參數(shù)確定鋪放軌跡函數(shù)，之后進行了實驗和仿真驗證，但他們對變剛度復(fù)合材料層合板進行研究時，采用的都是實驗和仿真的方法，沒有采用解析法。

本文提出一種對變剛度板孔邊應(yīng)力進行解析分析的方法，根據(jù)牛雪娟等[14]提出的利用流場函數(shù)構(gòu)造變剛度鋪放軌跡的理論，采用合適的流場函數(shù)擬合最大主應(yīng)力方向，通過優(yōu)化目標(biāo)函數(shù)確定構(gòu)造參數(shù)，得到適用于整個變剛度板的連續(xù)的鋪放軌跡函數(shù)即優(yōu)化函數(shù)，基于獲得的優(yōu)化函數(shù)，對變剛度板孔邊應(yīng)力進行解析分析。

1 規(guī)劃變剛度鋪放軌跡

1.1 確定最大主應(yīng)力方向

分析開孔復(fù)合材料直線鋪放單層板在一定受力工況下的數(shù)據(jù)，獲得層合板上各個單元的最大主應(yīng)力方向。根據(jù)文獻[15]所述，在相同的邊界條件下，各向同性板與各向異性板經(jīng)有限元分析后其最大主應(yīng)力方向相同。分析一含中心孔復(fù)合材料0°單層板，材料為45#鋼，單層板長155 mm，寬36 mm，厚2 mm，中心圓孔的直徑為6 mm。利用abaqus 6.12建立模型，層合板右端在長度方向承受面內(nèi)拉伸載荷10 N，邊界條件為左端完全固定，原點設(shè)置在孔的中心，網(wǎng)格劃分單元類型為C3D8R，運行后從.rpt文件中提取每個單元的S11S22S12。其中，S11為沿著x軸方向的應(yīng)力，S22為沿著y軸方向的應(yīng)力，S12為作用在yoz平面上沿y軸方向的剪力。利用式(1)、式(2)獲得其最大主應(yīng)方向，如圖1所示。

(1)

可知：

(2)

1.2 設(shè)計變剛度鋪放軌跡

由圖1可知，在遠離孔的區(qū)域最大主應(yīng)力方向呈直線狀分布，在孔的周圍區(qū)域最大主應(yīng)力方向呈余弦狀分布，整個板的最大主應(yīng)力方向的分布類似于圓柱繞流流線簇的分布。根據(jù)牛雪娟提出的理論，本文采用直線流與偶極流疊加來擬合最大主應(yīng)力方向，通過建立目標(biāo)函數(shù)獲得優(yōu)化的構(gòu)造參數(shù)，也就獲得了鋪放軌跡函數(shù)。

簡化后的流函數(shù)表達式為

(3)

在此表達式中構(gòu)造參數(shù)Γ未知，可通過如下優(yōu)化目標(biāo)函數(shù)獲得：

(4)

此優(yōu)化目標(biāo)函數(shù)的中心思想為孔周圍區(qū)域單元中心點處流函數(shù)的速度矢(uxi，uyi)方向與最大主應(yīng)力方向(cosθi，sinθi)的夾角最小，同時考慮到這些單元中心點處的夾角在優(yōu)化目標(biāo)函數(shù)中的重要性其實是不同的。

式(4)中的σ1i為最大主應(yīng)力代表權(quán)重就說明了不同的夾角其在優(yōu)化目標(biāo)函數(shù)中地位的不同。若某一點的最大主應(yīng)力較大，說明該點在實際的受力情況下承受了較大載荷。因此，最好使該點的纖維方向即鋪放軌跡方向與最大主應(yīng)力方向一致，使纖維更有效的承受載荷，這一點在優(yōu)化準(zhǔn)則中就占有重要地位；若某一點的最大主應(yīng)力較小，那么該點的纖維方向與最大主應(yīng)力方向是否一致，就不是十分關(guān)心的部分，其在優(yōu)化目標(biāo)函數(shù)中占有地位較小。

(5)

(6)

速度矢(uxi，uyi)可通過式(6)獲得，式(6)中的(xi，yi)為孔邊單元的中心點坐標(biāo)，利用該坐標(biāo)也可獲得變剛度板在x-y總坐標(biāo)系下，以變剛度板中心孔的中心為原點的極角θ。最大主應(yīng)力σ1i可通過式(5)獲得，利用該優(yōu)化目標(biāo)函數(shù)能夠求得構(gòu)造參數(shù)Γ，最終獲得該含孔單層板鋪放軌跡函數(shù)。

2 含孔單合板孔邊解析法分析理論基礎(chǔ)

2.1 經(jīng)典層合板理論

圖2所示為單元體平衡，1-2坐標(biāo)表示材料的局部坐標(biāo)，x-y表示層合板總坐標(biāo)，α表示從x軸到1軸的角度，逆時針方向為正。若圖3單元體x軸方向平衡，可得

σx=σ1cos2α+σ2sin2α-2τ12cosαsinα

(7)

2.2 傳統(tǒng)直線鋪放含中心孔單層板孔邊應(yīng)力分析

含中心孔單層板承受沿x軸方向的載荷且單層板的纖維方向與載荷方向成α角，如圖3所示。

對于這樣的各向異性板，已有學(xué)者[2]給出求其孔邊環(huán)向應(yīng)力的解析式，如式(8)所示。

(8)

其中，γ1、γ2可由式(9)給出：

(9)

式中E1、E2、G12、ν12為材料的彈性常數(shù)。

接下來，利用求得的孔邊環(huán)向應(yīng)力σθα可求出在孔邊不同極角θ處纖維方向的應(yīng)力σ1，垂直于纖維方向的應(yīng)力σ2以及剪應(yīng)力τ12由式(10)[2]給出：

(10)

2.3 曲線鋪放含中心孔單層板孔邊應(yīng)力分析

通過圖1可知，在孔周圍區(qū)域及遠離孔區(qū)域最大主應(yīng)力的方向，可預(yù)計利用優(yōu)化函數(shù)獲得的鋪放軌跡大致形狀如圖4所示。下面對曲線鋪放含中心孔單層板孔邊應(yīng)力分析，將利用圖4進行說明。

圖4中，1-2坐標(biāo)的1軸方向為鋪放軌跡在該點的切線方向即纖維的鋪放方向，2軸方向為垂直于切線方向。α表示切線與x軸之間的夾角，θ表示該點孔邊極角。α可由式(11)求出：

(11)

因為載荷施加在x軸方向，所以本文主要考察載荷方向的應(yīng)力σx。先利用式(6)、式(8)、式(9)、式(11)求出孔邊不同極角θ處的應(yīng)力σθαi，再利用式(7)、式(10)求出σθαi沿x軸方向的分力σxi。

3 算例分析及有限元驗證

3.1 算例分析

本文所用算例材料為AS4/9973，彈性常數(shù)E1=130 MPa，E2=9.24 MPa，G12=5.24 MPa，ν12=0.36幾何模型，網(wǎng)格類型，載荷大小，邊界條件等均與本文開頭時所用各向同性模型相同。利用優(yōu)化目標(biāo)函數(shù)即式(4)，獲得該模型的流函數(shù)的構(gòu)造參數(shù)Γ=-5.87，離散獲得的鋪放軌跡函數(shù)，相鄰兩條等流線間的數(shù)值差為1.5，利用MatLab繪制鋪放軌跡，如圖5所示。利用式(6)、式(11)求得變剛度板孔邊不同極角θ處鋪放軌跡切線的角度，即纖維鋪放角度α，如圖6所示。從圖6可看出，在θ約為90°、270°處，α約為0°，與最大主應(yīng)力方向近乎相同，可預(yù)計其孔邊應(yīng)力σx將會下降。

利用在本文提出的方法，計算后可獲得變剛度板孔邊不同極角θ處沿x軸方向的應(yīng)力σx，如圖7所示。紅線表示0°板孔邊不同極角處的σx，黑線表示變剛度板孔邊不同極角θ處的σx，從圖7可看出，在θ約為90°、270°處，變剛度板的σx比0°板的σx下降11.90%，傳統(tǒng)直線鋪放單層板的應(yīng)力集中得到較大改善。為更準(zhǔn)確地說明改善情況，在表1中列出了變剛度板及0°板的極角θ在90°及270°附近角度處的σx。

極角θ/(°)變剛度板/MPa0°板/MPa76．513．223614．267585．521．614224．534894．521．368124．5348103．513．146514．2675256．513．223614．2675265．521．621324．5348274．521．334224．5348283．513．146514．2675

從圖7中也可看出，變剛度板孔邊環(huán)向其他角度θ處的σx要比0°板的大。分析其原因，在于變剛度板實現(xiàn)了應(yīng)變的轉(zhuǎn)移，使原本集中在垂直受力方向的應(yīng)變轉(zhuǎn)移到其他孔邊位置，同時使材料更有效地承受載荷。

3.2 有限元驗證

利用式(6)、式(11)獲得鋪放軌跡在每個單元中心點處切線的角度，然后利用MatLab編程，再用生成的文件更改在3.1節(jié)中由有限元模型獲得的.INP文件，賦予每個單元纖維的鋪放角度，即2.1節(jié)中所說的α，由此獲得變剛度模型。將得到的變剛度模型重新導(dǎo)入有限元軟件Abaqus 6.12，提取運行后變剛度模型孔邊不同極角處的σx，并將其與3.1節(jié)中利用解析法獲得的σx進行對比，如圖8所示。紅線表示利用有限元方法得到的σx，黑線表示利用解析法得到的σx。經(jīng)分析可知，該解析法的最大誤差是應(yīng)力集中處的誤差，為5.63%。

4 結(jié)論

(1)利用流場函數(shù)構(gòu)造變剛度板的鋪放軌跡，建立起變剛度板力學(xué)解析分析與變剛度鋪放軌跡設(shè)計間的橋梁，使得采用解析法分析變剛度板孔邊應(yīng)力能夠?qū)崿F(xiàn)。

(2)觀察利用解析法及有限元法得到的應(yīng)力圖，發(fā)現(xiàn)二者的趨勢基本相同，發(fā)生應(yīng)力集中的位置一樣，最大誤差是應(yīng)力集中處的誤差。解析法及有限元法的最大誤差為5.63%。通過對比變剛度板及0°板孔邊應(yīng)力峰值，可知變剛度板應(yīng)力σx的峰值比0°板的應(yīng)力峰值σx下降11.90%。

(3)通過本文分析可知，若邊界條件、受力工況、幾何模型等條件不同，將導(dǎo)致其最大主應(yīng)力方向不同，其中載荷施加的方式和位置是影響最大主應(yīng)力方向的主要因素。這時需根據(jù)實際情況采用合適的流場函數(shù)擬合最大主應(yīng)力方向，利用優(yōu)化目標(biāo)函數(shù)獲得優(yōu)化函數(shù)。

[1] Nageswara D K Rao,Ramesh M Babu,Raja K Narender reddy,et al.Stress around square and rectangular cutouts in symmetric laminates[J].Composite Structures,2010,92(12):2845-2859.

[2] Yasar M Kaltackle,Arslan H M.Stress concentions of symmetrically laminated composite plates containing circular hole[J].Iranian Journal of science & Technology Transaction B,Engineering,2006,30(B4):447-460.

[3] Dharmendra S.Sharma,Nirav P Patel,Reena R Trivedi.Opti-mum design of laminates containing an elliptical hole[J].International Journal of Mechanical Sciences,2014,85(8):76-87.

[4] William L Ko.Stress concentration around a small circular hole in the HIMAT composite plate[R].California,NASA TM 86038,1985:1-16.

[5] Russo A,Zuccarello B.An accurate method to predict the stress concentration laminates with a circular hole under tensile loading[J].Mechanic of Composite Materials,2007,43(4):59-376.

[6] Michele Zappalorto,Paolo Andrea Carraro.An engineering formula for the stress concentration factor of orthotropic composite plates[J].Composite Part B:Engineering,2015,68(6):51-58.

[7] 韓小平,郭章新,朱西平,等.含孔復(fù)合材料層合板孔邊的應(yīng)力集中[J].復(fù)合材料學(xué)報,2009,26(1):168-173.

HAN Xiaoping,GUO Zhangxin,ZHU Xiping,et al.Stress distribution near the opening and stress concentration for composite laminates with acircular opening[J].Acta Materiae Compositae Sinica,2009,26(1):168-173.

[8] 李成,鄭艷萍,王迎佳.積分方程求解復(fù)合材料開口板的應(yīng)力分布[J].玻璃鋼/復(fù)合材料,2007,34(1):9-12.

LI Cheng,ZHENG Yanping,WANG Yingjia.Stress analysis of composite plate with holes by integral tunction[J].Fiber Reinforced Plastics/Composites,2007,34(1):9-12.

[9] 謝慈航,薛璞.含孔復(fù)合材料層合板的孔邊應(yīng)力分析及孔形優(yōu)化[J].陜西師范大學(xué)學(xué)報(自然版),2014,42(4):27-32.

XIE Cihang,XUE Pu.Optimization and analysis of stress concentration in composite laminates with different cutouts[J].Journal of Shaanxi Normal University(Natural Scienc Edition),2014,42(4):27-32.

[10] Gürdal Z,Olmedo R.In-plane response of laminates with spatially varying fibre orientations:Variable stiffness concept[J].AIAA Journal,1993,31(4):751-758.

[11] Wu K,Gürdal Z.Thermal testing of tow-placed,variable stiffness panels[C]//Proceedings of the AIAA/ASME/ ASCE/AHS/ASC 42nd Structures,Structural Dynamics and Materials Conference,Seattle,WA,AIAA 2001:1190-1212.

[12] Wu K,Gürdal Z,Starnes Jr J.Structural response of compression-loaded,tow-placed,variable stiffness panels[C]//Proceedings of the AIAA/ASME/ASCE/ AHS/ASC 43rd Structures,Structural Dynamics and Materials Conference.Denver,CO,AIAA 2002:1512-1535.

[13] ZHU Yangdan,LIU Jiancai,LIU Dong,et al.Fiber path optimization based on a family of curves in composite laminate with a center hole[J].Composites Part B:Engineering,2017,111(2):91-102.

[14] 牛雪娟,楊濤,杜宇.基于流場函數(shù)的變剛度層合板鋪放設(shè)計[J].固體火箭技術(shù),2014,37(6):848-855.

Niu Xuejuan,Yang Tao,Du Yu.Path planning of variable stiffness laminates based on the flow field function[J].Journal of Solid Rocket Technology,2014,37(6):848-855.

[15] Tosh M W,Kelly D W.Fiber steering for a composite C-beam[J].Composite Structures,2001,53(2):133-141.