考慮漸進損傷的纖維纏繞復合材料圓筒鋪層順序優化設計①

顏 標，郭凱特，校金友

(1.中國電子科技集團公司第三十八研究所，合肥 230088；2.西北工業大學 航天學院，西安 710072)

0 引言

纖維纏繞復合材料圓筒由于其高比強度、比模量、輕質化和良好的抗腐蝕等優越性能，目前正被廣泛用于航空航天和石油化工等領域[1-7]。復合材料的設計和制造過程參數嚴重影響著最終結構件的力學性能，因此如何確定纖維纏繞復合材料圓筒的設計及制造參數(鋪層順序、纏繞角、纖維體積分數和張力梯度等參數)，最大化其剛度和強度，同時盡可能的減輕構件質量，一直是工程設計中不斷追求的目標。實驗和理論研究表明，鋪層順序和纏繞角對于纏繞殼體的質量和力學性能有著重要的影響[8-11]。對于復合材料纖維纏繞圓筒而言合理的鋪層順序，不僅能夠提高圓筒的強度和剛度，還能減小復合材料圓筒的質量。對于內壓載荷下的復合材料圓筒，網格理論表明最佳纏繞角度是±55°，然而網格理論未考慮基體的作用，同時該結論又是將鋪層形式先定位±θ下所獲得的結果，并不一定是最佳結果，因為實際每一層鋪層角均可以變化即鋪層形式并不確定。因此，在鋪層形式未確定的情況下，研究特定載荷工況下復合材料圓筒的最佳鋪層形式對提高圓筒承載性能具有重要意義。

近十年來，采用智能優化算法尋找復合材料層合板或者復合材料圓筒最佳鋪層順序，眾多學者進行了研究。Park[10]運用GA、第一剪切理論與Tsai-Hill失效準則優化了對稱鋪層層合板在不同外載和不同邊界條件下的鋪層序列 ，使其承載能力最強。Pelletier[11]運用多目標遺傳算法優化了承內壓圓筒的鋪層，使得爆破壓強與環向剛度均滿足設計要求。Rich與Haftka[12]以及Nagendra[13]等優化了復合材料平板的鋪層，使其屈曲強度最大。Soremekun[14]運用GA優化了彎-扭耦合復合材料梁的鋪層，使其強度最大。Martins[15]通過實驗與模擬均發現當復合材料結構失效(最終失效)強度要大于功能失效(初始損傷)強度，而且優化鋪層時所采用的強度計算方法對優化結果有很大的影響。Mertiny[16]通過實驗方法研究了雙軸載荷作用下鋪層順序對于復合材料圓筒強度的影響規律，結果表明鋪層順序對于結構最終失效模式及極限載荷有顯著影響。Tabakov[17]基于3D彈性力學理論，結合Tsai-Wu失效準則，以最大爆破壓強為目標，優化了復合材料圓筒的鋪層順序，結果表明優化的鋪層并非常規鋪層，沒有明顯規律。目前，雖然以強度為目標函數的復合材料鋪層優化研究很多，但大部分是基于首層失效，未考慮漸進損傷，而復合材料損傷是漸進式的，出現初始損傷之后，載荷會重新分配，結構還能繼續承載。因此，基于初始損傷與GA優化出的鋪層不一定意味著結構具有最大的結構失效強度，有必要在鋪層優化設計中考慮結構漸進損傷行為，研究鋪層參數對圓筒結構極限強度的影響規律。

采用有限元進行漸進損傷包含兩個步驟：損傷準則和損傷演化。損傷準則用來判斷損傷的起始，損傷演化來控制損傷的擴展。關于失效準則，首先提出來的是唯象失效理論，該類失效理論依賴于大量精確的實驗數據，運用多項式將不同失效模式耦合起來，例如Tsai-Wu失效模式、Tsai-Hill失效模式、Hollfman準則等[18]。Hashin和Rotem[19]首先提出了基于損傷機理的新一代損傷準則，分開考慮不同損傷模式，該準則將復合材料損傷分為基體損傷與纖維損傷，其中基體損傷考慮了剪切與拉壓的相互作用，采用了二次耦合關系。Puck[20]與Schurmann[21]將Mohr和Coulomb理論[11]引入Hashin準則，提出了失效平面與失效角的概念，并重新推導出了纖維與基體的損傷準則并提出了Puck準則，但是該準則十分復雜，包含許多需要實驗擬合的參數。后來，Davila[22]基于Puck準則，提出了LaCR準則，不包含任何擬合參數，所有參數都具有明確的物理意義。因此，本文復合材料圓筒損傷起始采用該準則。

本文采用遺傳算法和復合材料結構漸進損傷理論，針對無約束條件以及纏繞工藝與環向剛度約束下三種工況，尋求內壓載荷下使得結構性能最佳的鋪層順序。遺傳算法基于Matlab平臺實現，漸進損傷采用Abaqus軟件實現，采用Python語言讀寫文件，建立兩個軟件平臺之間信息傳遞橋梁，搭建起一種內壓載荷作用下考慮漸進損傷的復合材料圓筒鋪層優化方法，并對結果進行了分析討論。

1 復合材料圓筒有限元模型

本文采用的有限元模型為8層纖維纏繞圓筒，由于該圓筒為纏繞制成，因此鋪層模式為均衡鋪層，即[±θ1±θ2…±θn]。為減小計算量，取1/8結構分析。其邊界條件如圖1所示，環向與圓筒右端為對稱邊界，左端為端部封閉邊界條件，即固定1，2方向的位移，沿筒身方向自由，且其端部載荷為

圖1 圓筒有限元模型Fig.1 Finite element model of cylinder

(1)

式中p為圓筒承受內壓；r為圓筒內半徑；r1為圓筒外徑；σA為圓筒端部表面壓強。

本文采取的圓筒模型內半徑r=50.8 mm，厚度t=2.75 mm，長度L=1100 mm，材料采用Glass/Epoxy復材，參數取自文獻[23]，如表1所示。模型劃分單元數目為環向15個，軸向80個。采用C3D8R單元，層間性能用基體性能近似替代。沙漏控制為Enhanced。圓筒的破壞通過監控點(圖1右上角)徑向位移-內壓曲線判斷。

表1 復合材料單向板材料性能Table 1 Mechanical properties of composite

2 復合材料圓筒損傷起始及損傷演化

在對復合材料圓筒進行漸進損傷分析時，首先應給出相應的復合材料起始損傷判據。由于LaCR準則不僅預測精度較高，不包含任何擬合參數，且所有參數都具有明確的物理意義。本文運用LaCR準則來預測復合材料的面內損傷，LaCR準則是基于損傷機理推導出的，損傷模式分為四種：纖維拉伸損傷、纖維壓縮損傷、基體拉伸損傷與基體壓縮損傷。下面針對這四種損傷模式做以敘述。

2.1 損傷起始準則

(1)基體拉伸失效

圖2給出了各種載荷作用下復合材料斷裂平面示意圖。其中，橫向拉伸和面內剪切產生的斷裂平面均與纖維方向平行，LaCR準則的基體拉伸失效考慮了拉伸與剪切的耦合，其表達式如下：

Fmt=(1-g)σ22/YT+g(σ22/YT)2+(τ22/SL)2

=1

(2)

圖2 不同載荷下斷裂平面Fig.2 Fracture plane under different loads

式中g=GⅠ/GⅡ，GⅠ為一型裂紋斷裂韌性，GⅡ為二型裂紋斷裂韌性；YT為橫向拉伸強度；SL為縱向剪切強度。

(2)基體壓縮失效

基體壓縮損傷較為復雜，基體損傷的失效平面的法向方向n與橫向壓縮載荷方向成一定的夾角α，如圖3所示。

圖3 壓縮損傷失效平面Fig.3 Failure plane of compression damage

將直角坐標系應力繞1軸旋轉α度得到：

(3)

σn對失效平面的壓縮會產生摩擦力，會削弱剪切應力對損傷的影響。因此，用等效切應力來描述失效平面的切應力：

(4)

其中，

(5)

式中YC為橫向壓縮強度；α0為失效平面法向方向與2方向的夾角，α0根據Puck實驗數據， 一般取53°。

LaCR準則運用二次相互作用來描述基體的壓縮損傷：

(6)

其中，〈·〉為Cauley算子，定義為

(7)

(3) 纖維拉伸損傷

對于纖維拉伸損傷，最大應力準則便能給出滿意的結果：

(8)

(4)纖維壓縮損傷

纖維壓縮損傷較為復雜，其失效準則直到現在仍然在不斷的發展中。許多研究認為，纖維的損傷是由基體失效引起的纖維屈曲導致的[24]。

LaCR準則采取了相同的觀點推導纖維的壓縮損傷準則，假設纖維的局部彎曲角度為φ，如圖4所示?？梢运愕镁植孔鴺讼迪碌膽顟B為：由于基體的拉伸與壓縮損傷準則表達式不同，因此將分基體受拉與受壓來給出纖維的壓縮損傷準則。當基體受拉時，將式(9)代入式(2)，可得式(10)。

(9)

(10)

圖4 纖維局部屈曲失效模式Fig.4 Failure mode of local fiber buckling

當基體受拉時，其失效準則表達式如下：

(11)

其中，局部屈曲角φ表達式如下：

(12)

(13)

2.2 損傷演化

當材料發生損傷時對材料剛度矩陣進行折減。纖維與基體損傷分別用變量df與dm描述，那么當損傷發生后的剛度矩陣可表示為

(14)

其中，K=1-(1-df)(1-dm)ν12ν21，E11、E22、G12、ν12、ν21為材料彈性參數。當完全損傷發生后，若直接將剛度矩陣相應分量折減為0，會引起剛度矩陣奇異。因此，當完全損傷后df(m)一般不取1，可取一個十分接近1的值，本文取1-10-6。

圖5 雙線性損傷演化準則Fig.5 Bilinear damage evolution criterion

由于損傷的不可逆性，損傷變量只能單調增加，且其值不能大于1：

(15)

在有限元中運用應變軟化模型會導致嚴重的網格依賴性，損傷過程的能量耗散與網格大小相關，為緩解網格的依賴性，本文采用Bazant提出的裂紋帶模型(Crack band model)，認為裂紋產生所耗損的能量應該等于裂紋表面形成所需的能量即斷裂韌性：

(16)

式中Gc為斷裂韌性；Lc為單元特征長度。

該方法不僅有效緩解了有限元模擬過程的網格依賴性，同時將損傷過程的能量耗散與實際物理量斷裂韌性聯系起來了，更加符合實際。

材料本構模型存在應變軟化或性能衰退會導致收斂困難，采用粘性正則化可以改善收斂性，該方法在計算損傷演化時不使用當前損傷變量，而是使用正則化后的損傷變量：

(17)

3 考慮漸進損傷的圓筒鋪層優化方法

遺傳算法是模擬達爾文生物進化論的自然選擇和孟德爾遺傳學機理的生物進化過程的計算模型，是一種通用模擬自然進化過程搜索最優解的方法。求解的過程主要由問題的數學建模、基因編碼、適應度函數求解與運用遺傳算子生成下一代四部分組成。以上過程反復迭代求解直到求出滿足要求的解或迭代步數達到指定值，遺傳算法的主要過程可分為編碼、解碼、雜交、突變和個體適應度評價幾個步驟。

本文的遺傳算法優化模型在Matlab中實現，爆破內壓運用漸進損傷模型預測，在商業有限元軟件Abaqus中實現。因此，鋪層優化算法主要分為優化算法與有限元求解兩部分組成，如圖6所示。其中，黃色實線框內為GA算法求解過程，綠色虛線框內為有限元求解過程，GA算法模塊與有限元求解模塊的數據交流通過輸入文件Sequence.in與輸出文件Pressure.out進行橋接。

圖6 鋪層優化流程圖Fig.6 Flow chart of optimization process

上述鋪層優化算法計算過程可總結為：初始種群產生后，將其鋪層序列寫入文件Sequence.in中，然后Python程序Runcae.py從中讀入鋪層序列建立CAE模型，并調用漸進損傷子程序PD.for進行承內壓圓筒有限元分析，完成后生成結果文件，運用Getdata.py腳本讀出爆破內壓并寫出到文件Pressure.out中，最后將爆破內壓作為個體適應度返回到GA部分并生成下一代種群。在本次遺傳算法求解中，雜交概率取100%，突變概率取2%，種群數目取20，最大迭代次數取50。

4 基于遺傳算法的復合材料圓筒鋪層順序優化

4.1 復合材料圓筒有限元模型驗證

為進一步驗證上述損傷模型預測圓筒爆破壓強的準確性，運用該模型計算[±45]4與[±60]4鋪層的圓筒爆破內壓并與實驗對比，圖 7為圓筒內壓與監測點徑向變形關系圖，可看到[±45]4鋪層在壓強達到8.34 MPa時，徑向位移急劇增加，即圓筒破壞。

圖7 壓強-監測點位移關系曲線Fig.7 Curves between pressure and monitoring node displacement

對于[±60]4鋪層，當圓筒爆破內壓達11.7 MPa時失效，但徑向位移并未急劇增加，這主要是由于該鋪層下失效模式為圓筒中部軸向失效，此時軸向位移急劇增大后圓筒失效。損傷模型預測強度值與實驗值對比如表2所示，實驗值取自文獻[26]。兩種鋪層下誤差分別為1.65%與11.3%，驗證了上述損傷模型的準確性，接下來將運用基于上述損傷模型與GA的鋪層優化算法對承內壓復合材料圓筒進行鋪層設計。

表2 實驗與模擬爆破壓強對比Table 2 Comparison of burst pressure for experiment and simulation

4.2 無約束的鋪層設計

圖8為鋪層角度無約束下的優化，大概30代后達到最大爆破壓強17.56 MPa，平均爆破壓強不斷增加直到22代后在穩定值附近波動，平均值沒有趨近于最大值是由于通過遺傳算子雜交、變異等產生的下一代中總可能包含非最優解。正是由于此，遺傳算法可跳出局部最優從全局尋求最優解，而不是收斂到局部最優。該種工況(#1)下優化出的最佳鋪層如表3所示。

圖8 無約束情況下鋪層設計Fig.8 Stacking sequence design under unconstrained condition

可看到，在考慮漸進損傷優化出的鋪層沒有規律，而實際生產中，一般采用均衡反對稱鋪層[±55]4。實際上圓筒承受內壓時，當某些層損傷后，圓筒原有的鋪層規律也將消失，Pelletier[11]與Irisarri[8]運用GA優化出的層合板鋪層序列也沒有規律。而工業上優化出的鋪層有規律是由于其優化前假設的鋪層序列就具有規律，如假設鋪層序列為[±θ]n。

表3 不同約束條件下優化的最佳鋪層序列Table 3 Optimized stacking sequence under different constrains

圖9對比了[±55]4鋪層與優化出的鋪層的壓強-檢測點位移關系曲線，可以看到優化出的鋪層爆破壓強為17.58 MPa，大于傳統采用的[±55]4鋪層。可以發現，兩種鋪層下基體損傷的產生均遠早于最終失效，當纖維失效后結構不再能承受更大壓強，纖維失效意味著最終失效說明鋪層較為合理，因為充分發揮了纖維的高強度特性。同時，可發現優化的鋪層具有更大的剛度，這主要是由于在承力方向上布置了更多的纖維，這也導致了更大的爆破壓強。此外，優化出的鋪層的失效應變更小，這主要是由于優化出的鋪層在承力方向上具有更多纖維，因此導致了比[±55]4鋪層提前出現最終失效(失效位移更小)。

圖9 不同條件下優化鋪層的壓強-檢測點位移曲線Fig.9 Curves between pressure and monitoring node for optimized stacking sequence under different constrains

4.3 制造工藝約束下的鋪層設計

由于封頭與摩擦系數的限制，纏繞工藝中纖維的纏繞角不能過小，否則無法在封頭處從+θ角過渡到-θ。本文根據圓筒尺寸與摩擦系數限制纏繞角不能小于20°。由表3看到無約束條件下優化出的鋪層為[±84/±0/±88/±60]，其中的0°鋪層無法纏繞。因此考慮工藝限制，上述鋪層在實際生產中無法實現，現在工藝限制下對鋪層進行優化。

圖10為考慮工藝約束下的鋪層優化結果，可以看到在第12代后到達最大壓強17.11 MPa。平均值與無約束條件下一樣，在不斷波動中趨于穩定值，在27代后便在最大值附近波動。此種情況下優化的鋪層角為[±44/±89/±40/±69]，所有角度均在要求內，而且優化出的角度也沒有明顯規律。將該鋪層角度下的位移壓強繪制在圖9中，如藍色曲線所示，其壓強略小于無約束下的鋪層。同時，可發現不同鋪層形式下基體損傷的時間差別不大，纖維的損傷的時間差別較大，爆破壓強越大纖維損傷的時間越早。

圖10 制造工藝約束情況下鋪層設計Fig.10 Stacking sequence design under manufacturing constrains

4.4 剛度約束下的鋪層設計

本節針對高剛度、高強度要求進行約束，其中強度同樣由漸進損傷模型給出。對于剛度，其可由損傷發生前的壓強-監測點位移曲線獲得，但是每次均調用Abaqus將帶來巨大的計算量，使得計算效率較低，因此在本次優化中剛度由層合板理論計算得到。當給出鋪層模式后，首先在Matlab中計算該鋪層模式下的環向剛度E_y，若環向剛度小于指定值E_0，則直接跳過調用Abaqus，并將適應度函數直接賦值為0，從而該鋪層序列會在后續迭代中迅速被淘汰。只有剛度滿足要求的鋪層才會調用Abaqus計算強度，從而減小計算量，本次模擬中E_0取20 GPa。

圖11為剛度約束下優化結果圖，20代后壓強到達穩定值15.09 MPa，平均值與上文兩種情況類似，在震蕩中慢慢趨于穩定。優化出的鋪層為[±67/±78/±20/±86]，可看到角度均大于等于20，符合工藝要求。將該鋪層下的爆破壓強-監測點位移曲線繪制在圖9中如紅色曲線所示，壓強明顯小于其余鋪層模式。同時，可看到與綠色曲線斜率(剛度)基本一樣，但綠色曲線(無約束下鋪層應力應變曲線)對應的鋪層不符合工藝約束，其存在0°層。而且該鋪層模式下圓筒破壞時參考點徑向位移沒有發生突變，這主要是由于該鋪層下失效模式為圓筒中部軸向失效，此時軸向位移急劇增大后圓筒失效。對比圖中四條爆破曲線，可看到[±55]4不僅強度較小，同時剛度也較小，并不是最合適的鋪層，相比于[±55]4鋪層，本文給出的鋪層具有明顯優勢。

圖11 剛度約束下的鋪層設計Fig.11 Stacking sequence design under constrained modulus

5 結論

(1)建立了內壓載荷作用下復合材料圓筒鋪層順序優化方法，該方法考慮了復合材料漸進損傷過程。優化算法基于Matlab平臺建立，漸進損傷分析基于Abaqus軟件平臺分析。采用Python語言編寫腳本進行兩個軟件平臺的數據交換，搭建了基于Matlab和Abaqus軟件平臺的復合材料圓筒鋪層順序的優化設計平臺。

(2)基于本文所建立的復合材料漸進損傷分析方法，對[±45]4與[±60]4鋪層的圓筒爆破內壓并與實驗對比。結果表明，誤差在13%以內，從而驗證了本文漸進損傷模型的正確性，為后續鋪層順序的優化奠定了基礎。

(3)基于搭建的復合材料圓筒鋪層順序的優化平臺，對三種工況條件下的鋪層優化做了分析計算。結果表明，無約束條件和考慮工藝約束條件下所得到的鋪層均無明顯規律，無約束條件下爆破壓強和考慮工藝約束下的爆破壓強均高于傳統鋪層[±55]4所對應的爆破壓強，分別高出6.5%和3.7%，說明該方法可在一定程度上提高復合材料圓筒的爆破壓強。