含中心孔復合材料變剛度板孔邊應力解析法分析

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

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

0 引言

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

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

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

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

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

1.1 確定最大主應力方向

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

(1)

可知：

(2)

1.2 設計變剛度鋪放軌跡

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

簡化后的流函數表達式為

(3)

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

(4)

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

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

(5)

(6)

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

2 含孔單合板孔邊解析法分析理論基礎

2.1 經典層合板理論

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

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

(7)

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

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

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

(8)

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

(9)

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

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

(10)

2.3 曲線鋪放含中心孔單層板孔邊應力分析

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

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

(11)

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

3 算例分析及有限元驗證

3.1 算例分析

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

利用在本文提出的方法，計算后可獲得變剛度板孔邊不同極角θ處沿x軸方向的應力σx，如圖7所示。紅線表示0°板孔邊不同極角處的σx，黑線表示變剛度板孔邊不同極角θ處的σx，從圖7可看出，在θ約為90°、270°處，變剛度板的σx比0°板的σx下降11.90%，傳統(tǒng)直線鋪放單層板的應力集中得到較大改善。為更準確地說明改善情況，在表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°板的大。分析其原因，在于變剛度板實現了應變的轉移，使原本集中在垂直受力方向的應變轉移到其他孔邊位置，同時使材料更有效地承受載荷。

3.2 有限元驗證

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

4 結論

(1)利用流場函數構造變剛度板的鋪放軌跡，建立起變剛度板力學解析分析與變剛度鋪放軌跡設計間的橋梁，使得采用解析法分析變剛度板孔邊應力能夠實現。

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

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

[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] 韓小平,郭章新,朱西平,等.含孔復合材料層合板孔邊的應力集中[J].復合材料學報,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] 李成,鄭艷萍,王迎佳.積分方程求解復合材料開口板的應力分布[J].玻璃鋼/復合材料,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] 謝慈航,薛璞.含孔復合材料層合板的孔邊應力分析及孔形優(yōu)化[J].陜西師范大學學報(自然版),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] 牛雪娟,楊濤,杜宇.基于流場函數的變剛度層合板鋪放設計[J].固體火箭技術,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.