復合材料開口有限元建模方法研究

熊美蓉，陳琳，劉傳軍

(1.中國商飛上海飛機制造有限公司，復合材料中心，上海 200123；2.中國商飛北京民用飛機技術研究中心，民用飛機結構與復合材料北京市重點實驗室，北京 102211)

0 引言

復合材料開口產生的局部應力集中會引起結構強度下降，應力集中與層壓板鋪層比例、載荷情況、缺口尺寸和形狀、缺口約束情況有關。

在比較主流的分析方法中，對于單軸載荷下復合材料開口結構，可采用基于斷裂力學模型Mar-Lin準則[1]進行層合板剩余強度分析，將開口定義為損傷，通過獲得不同鋪層下的復合材料斷裂韌性、奇異性參數、構型因子等，擬合得到不同開口尺寸下的層合板剩余強度；另一種方法是采用基于應力失效模型的Withney-Nuismer方法[2]，獲得圓孔層合板在圓孔中心線上垂直于遠場載荷方向的應力分布，再通過點應力或平均應力方法開展失效分析。對于復雜載荷下，可采用Lehknitskii理論[3]中提出的復變函數計算含橢圓孔(圓孔)的正交各向異性無限大板，在復雜遠場載荷作用下的孔邊應力分布，再通過點應力或平均應力方法開展失效分析。

常見的飛機復合材料開口結構和其載荷狀況都較為復雜，一般采用精細有限元建模后提取孔邊最大應力或應變開展分析。而開口區的復雜應力分布對孔邊網格要求較高，網格尺寸、數量和質量等都會直接影響數值計算的結果精確性，同時網格尺寸過小也影響計算效率。本文主要研究復合材料開口結構在單軸和雙軸兩種載荷狀態下復合材料開口結構孔邊網格尺寸對有限元計算精度的影響，為開展復合材料開口有限元建模和分析提供指導。

1 單軸載荷分析

本文選用T800級碳纖維復合材料開口層合板，其鋪層為[45 °/135 °/0 °/90 °/90 °/0 °/135 °/45 °]，層合板為100 mm× 100 mm 方板，圓孔直徑為25.4 mm，主要分析層合板單軸載荷作用下，孔邊的應力應變分布情況。

1.1 分析方法

圖1 在單向拉伸作用下、含圓形孔開口區應力分布示意圖

對有限寬板，沿y軸方向的應力σx(y)表達式：

(1)

(2)

其中β為有限寬度修正系數，修正系數根據材料(各向同性和各向異性)和開口類型(圓孔或橢圓孔)變化；KT為孔的應力集中系數。

1.2 孔邊應力計算

采用1.1節中Withney-Nuismer方法計算在100 MPa沿X軸作用下，層合板孔邊沿Y軸的應力分布。

圖2 孔邊Quad單元與節點示意圖

當采用四節點的Quad單元劃分有限元網格時，對模型求解得到的孔邊最大應變，即為孔邊一圈所有單元中的最大應變。如圖2所示，可通過節點位置處的應變估算得到001單元的應變值，計算圖中1、2、3、4節點位置處應變的平均值，作為單元的中心點的應變，這樣即可預測不同有限元網格尺寸下孔邊的最大應力和應變。

取開口直徑為25.4 mm的復合材料層合板，計算網格大小為0.000 1 mm、0.001 mm至10 mm時，孔邊單元對應的最大應力、最大應變及相對于孔邊最大值精確解的誤差(表1)。從圖3中也可以看出，開口直徑為25.4 mm時，當網格尺寸在0.000 1 mm到0.01 mm時，孔邊的最大應變趨于收斂，精度不超過0.1%。隨著網格尺寸增大，孔邊最大應變略微下降，計算精度也逐漸降低。

表1 單軸載荷網格尺寸對計算精度的計算

圖3 不同開口區網格尺寸對孔邊最大應變的影響

1.3 開口尺寸及鋪層影響

從Withney-Nuismer方法的解析公式可以看出，孔邊應力分布主要是開口尺寸與層合板鋪層的函數，根據1.1中對有限元模型中孔邊最大應力應變的估算方法，可以分析不同開口直徑下不同網格尺寸的精度大小。如圖4所示，圖中曲線為針對準各向同性層合板開口直徑分別為6.35 mm、12.7 mm、25.4 mm、50 mm、100 mm及200 mm下不同網格尺寸下的精度大小。可以看出，開口直徑越小時對孔邊網格大小越為敏感，對較小直徑如6.35 mm的開口，孔邊網格為0.1 mm時，最大應變精度為3.54%，而較大直徑下200 mm的開口，孔邊網格為0.1 mm時，最大應變精度依然較高為0.58%。當直徑大于100 mm后，孔邊網格小于1 mm可獲得5%以上的精度。

圖4 單軸載荷不同開口直徑下網格尺寸對計算精度的影響

復合材料層合板鋪層主要影響孔邊應力集中系數，從而影響孔邊應力應變分布。圖5主要分析了四種鋪層比例下不同網格尺寸的精度大小。對于開口直徑為25.4 mm的層合板，當AML(AML=±45 °鋪層百分比-0 °鋪層百分比)越小，即0 °鋪層比例越高時對孔邊網格大小更為敏感。對絕大多數鋪層，在孔邊網格尺寸在0.01 mm時精度可以控制在0.21%以下。相對于開口尺寸，鋪層的影響要較小一些，因此可以主要根據開口直徑和精度要求估算出孔邊網格大小，再針對鋪層作簡單分析來確定有限元模型網格劃分，以滿足強度計算要求。

2.2.5 水煎煮提取工藝正交試驗 根據上述的藥效篩選結果，對水煎煮提取工藝進行正交試驗，以君藥黃芪中黃芪甲苷含量（Y1）和固形物質量（Y2）為評價指標，并分別賦予權重系數0.6和0.4，計算綜合評分值[Y，Yi＝(X1i/X1max×0.6＋X2i/X2max×0.4)×100]，對影響煎煮的加水量（A）、煎煮時間（B）和煎煮次數（C）進行優選。因素與水平見表3。

圖5 單軸載荷不同鋪層比例下網格尺寸對計算精度的影響

2 雙軸載荷分析

本文選用T800級碳纖維復合材料開口層合板，其鋪層為[45 °/135 °/0 °/90 °/90 °/0 °/135 °/45 °]，層合板為100 mm× 100 mm 方板，圓孔直徑為25.4 mm，主要分析在雙軸應力100 MPa作用下，孔邊的應力應變分布情況。

2.1 分析方法

雙軸載荷下開口分析采用Lehknitskii方法，Lehknitskii模型如圖 6所示。分析假設如下：

(1)材料主軸方向12與橢圓孔主軸方向xy一致。

(2)平板尺寸與橢圓孔或者圓孔尺寸相比可認為是無限大。

(4)考慮平面應變或者廣義平面應力情形。

圖6 在遠場力作用下、含開口的無限大平板

孔緣的應力由遠場應力和擾動應力構成，形式如下，

(3)

(4)

式中：Re——復數的實部；

Zk——復數變量，形式為Zk=x+yμk(k=1,2)；

Φk(Zk)為復應力函數， 由橢圓孔存在引起的擾動應力對應的應力函數為：

(5)

2.2 孔邊應力計算

根據2.1中對有限元模型中孔邊最大應力應變的估算方法，取開口直徑為25.4 mm的復合材料層合板，采用Lehknitskii方法，計算復合材料開口層合板在計算在X軸遠場載荷100 MPa與Y軸遠場載荷100 MPa作用下，網格大小為0.000 1 mm、0.001 mm至10 mm時，孔邊單元對應的最大應力、最大應變及相對于孔邊最大值精確解的誤差。從表2中也可以看出，對于開口直徑為25.4 mm時，當網格尺寸在0.000 1～0.01 mm時，孔邊的最大應變趨于收斂，精度不超過0.1%，與單軸載荷趨勢一致。

表2 雙軸載荷網格尺寸對計算精度的計算

2.3 開口尺寸及鋪層影響

在雙軸載荷下，復合材料層合板開口孔邊的應力應變分布受鋪層影響較大。對于25 °/50 °/25 °準各向同性鋪層的層合板，在100 MPa、100 MPa雙軸載荷作用下，其0～180 °孔邊的應力分布如圖7所示，可知最大主應力出現在孔邊0～360 °所有方向上。而對于純±45 °鋪層層合板，在100 MPa、100 MPa雙軸載荷作用下，其0～180 °孔邊的應力分布如圖8所示，可知最大主應力出現在孔邊45 °、135 °、225 °、275 °方向。此外，應力應變分布還與x軸、y軸載荷比例有關，本文暫不做具體分析。

圖7 準各向同性層合板孔邊應力應變分布(σx ∶σy=1 ∶1)

圖8 純±45 °鋪層層合板孔邊應力應變分布(σx ∶σy=1 ∶1)

根據2.1中對有限元模型中孔邊最大應變的估算方法，可以分析不同開口直徑下不同網格尺寸的精度大小。如圖9所示，圖中曲線為針對雙軸載荷下的準各向同性層合板，開口直徑分別為6.35 mm、12.7 mm、25.4 mm、50 mm及100 mm下不同網格尺寸下的精度大小。可以看出，開口直徑越小時對孔邊網格大小越為敏感，對較小直徑如6.35 mm的開口，孔邊網格為0.1 mm時，最大應變精度為3.54%，而較大直徑100 mm的開口，孔邊網格為0.1 mm時，最大應變精度依然較高為0.58%。當直徑大于100 mm后，孔邊網格小于1 mm可獲得5%以上的精度。相比于單軸，不同開口直徑層合板在雙軸載荷下精度略高，主要原因是層合板的泊松效應一定程度降低了孔邊的應力集中，應力分布受網格尺寸的影響也略為降低。

圖9 準各向同性層合板不同開口直徑開口網格尺寸對孔邊最大應力的影響(σx ∶σy=1 ∶1)

采用同樣方法，對雙軸載荷下的純±45 °鋪層層合板開展分析，取45 °方向最大應力，分析不同開口直徑下不同網格尺寸的精度大小如圖10所示。可以看出，與前文中分析結果一致，開口直徑越小時對孔邊網格大小越為敏感。當開口直徑為6.35 mm和12.7 mm層合板時，圖中可明顯看出，當網格尺寸>1 mm后曲線出現拐點，主要原因為±45 °鋪層層合板在x，y軸雙軸加載時，沿著45 °方向離孔中心超過5 mm后，應力有所升高，如圖11所示。

圖10 純±45 °鋪層層合板，不同開口直徑開口網格尺寸對孔邊最大應力的影響(σx ∶σy=1 ∶1)

圖11 純±45 °鋪層層合板，E11隨距離孔邊位置變化趨勢(σx：σy=1 ∶1)

通過對比兩種鋪層下不同開口的精度分布，可發現純±45 °鋪層相比準各向同性鋪層對孔邊網格尺寸更為敏感，且鋪層比例影響孔邊最大應力點的分布位置。因此對雙軸載荷的復合材料開口結構建模時，需按實際鋪層進行估算，以確定有限元模型網格劃分，以滿足強度計算要求。

3 有限元網格尺寸

為驗證上述分析方法，本文在ABAQUS CAE中建立了復合材料開口結構有限元模型，定義載荷和邊界條件，對模型網格劃分，并對孔周邊網格進行細化。圖12與圖13分別為準各向同性鋪層層合板在雙軸載荷下的孔邊最大主應力和最大應變的分布云圖。在最大主應力的分布云圖中可知，孔邊最大主應變為3 466ε，E11方向最大應變為3 458ε，與解析法計算的E11最大應變相差9ε，偏差在0.26%。

圖12 復合材料開口結構雙軸加載有限元模型

圖13 復合材料開口結構雙軸載荷的最大主應力和最大應變分布云圖

本文同時開展了單軸載荷下有限元方法、Withney-Nuismer方法和Lehknitskii方法計算結果的對比，如圖14所示。從圖14可看出，解析解和數值解的計算結果一致性較好，有效驗證了解析方法計算的合理性。

圖14 有限元方法和解析法計算結果對比

4 結論

本文通過研究復合材料開口層合板在不同載荷狀態、開口尺寸及鋪層比例下，有限元網格尺寸與數值計算精度的關系，通過分析，得到如下主要結論：

(1)在單軸載荷下，開口直徑越小時對孔邊網格大小越為敏感，0 °鋪層比例越高時對孔邊網格大小更為敏感。對開口直徑25.4 mm、大多數鋪層的層合板，在孔邊網格尺寸0.01 mm時精度在0.21%以上。

(2)在雙軸載荷下，開口直徑越小時對孔邊網格大小越為敏感;相對單軸載荷，存在一定泊松效應降低應力集中，降低網格敏感度。孔邊最大應力位置與鋪層比例和載荷狀態相關。