劉振國 胡 杰 胡 龍

(北京航空航天大學 航空科學與工程學院，北京100191)

由于復合材料的比強度、比剛度較高，近年來應用十分廣泛.通過改變鋪層形式和厚度來改變剛度、強度、屈曲和振動等性能是復合材料設計的一個優(yōu)勢.遺傳算法作為一種高效的離散空間搜索算法，在復合材料鋪層的優(yōu)化設計中適用性較好，國內外學者對其進行的研究較多［1］.文獻［2-5］對層合板的遺傳算法本身的改進作了研究，提高了收斂率，避免了早熟現(xiàn)象的產生.文獻［6］使用了神經(jīng)網(wǎng)絡和遺傳算法分別對層合板的鋪層比例和鋪層順序進行了優(yōu)化，但其將兩者的優(yōu)化看成是兩個獨立的過程，并不能達到最優(yōu)化的結果.另外該方法不對整體鋪層數(shù)量進行優(yōu)化.上述研究中在模擬結構力學性能方面均采用了經(jīng)典層合板理論，優(yōu)點是較為可靠、計算量較小，但是難以處理結構形狀復雜和邊界條件復雜的問題.目前在處理實際優(yōu)化問題上，有限元方法逐漸成為主流.文獻［7］借用有限元分析結合特殊設計的遺傳算子，實現(xiàn)了鋪層材料、鋪層順序和鋪層厚度的優(yōu)化.文獻［8］提出了一種改進的遺傳算法用以實現(xiàn)復合材料的多學科優(yōu)化.文獻［9］引入了層合板彎曲剛度等中間參數(shù)，降低了優(yōu)化過程的有限元計算量.但他們的研究僅限于對層合板的整體優(yōu)化.

由于在諸如飛機機翼翼面、螺旋槳葉片等大面積應用復合材料的部件上受力復雜，應力變化比較劇烈，采用相同的鋪層會降低材料的效率，而獨立分區(qū)優(yōu)化未考慮相互間的聯(lián)系，難以達到整體目標函數(shù)的最優(yōu)解，所以有必要研究在分區(qū)設計的基礎上進行整體的優(yōu)化.目前有關分區(qū)設計方面的研究成果較少.分區(qū)設計的缺點一方面是設計變量大大增加，另一方面區(qū)域與區(qū)域間差異較大的鋪層過渡也會對工藝性提出挑戰(zhàn).文獻［10］提出了一種復合材料全局優(yōu)化方法，通過對已經(jīng)分區(qū)的層合板進行補強的方式，對復合材料快艇艇身進行了優(yōu)化，得到了工藝性較好的結果.文獻［11］通過多島遺傳和序列二次規(guī)劃法相結合的形式對機翼翼盒段進行了整體優(yōu)化，翼盒間相互獨立，其內部保持鋪層一致，未考慮工藝.

本文結合整體和分區(qū)域設計的優(yōu)勢，提出了公共鋪層和局部鋪層共存的分級優(yōu)化模式，一方面保存了對局部區(qū)域特殊承載進行特殊設計的能力，另一方面全局鋪層的存在也保證了更好的工藝性，減少了設計參數(shù)可節(jié)約計算時間.另外，通過變異操作實現(xiàn)了鋪層的增減也降低了鋪層厚度變更的復雜程度.通過Python語言進行遺傳算法的流程的實現(xiàn)，驅動ABAQUS進行模型的分析計算，實現(xiàn)層合板分級遺傳、鋪層順序和鋪層厚度的優(yōu)化.

1 單元模型

通過ABAQUS的用戶自定義子程序，將經(jīng)典層合板理論計算的剛度矩陣和失效形式植入單元計算.使得優(yōu)化過程僅需輸出鋪層形式即可進行力學分析和失效判定.失效判定過程對每個單層逐一進行校核，單層失效即判定為失效.

1.1 剛度矩陣計算

由經(jīng)典層合板理論［12］得層合板合成軸力和彎矩為

其中，ε(k)為第 k層層合板產生的應變，ε(k)={εx(k)，εy(k)，εz(k)};Q-(k)為第k層層合板的剛度矩陣;zk為第k層中心沿厚度方向坐標.

由此，對于給定鋪層順序、厚度及鋪層角的層合板，即可求得層合板剛度.

1.2 破壞判定準則

判定層合板破壞存在多種準則，本文取平面應力狀態(tài)Tsai-Wu張量準則作為破壞判定條件［2］:

分析時，令式(3)左側計算值為強度因子STW，在每個單元內對鋪層進行逐層校核，判斷是否失效.

2 遺傳算法

遺傳算法是一種離散的空間搜索算法.通過將可能的解編碼成一向量——染色體，通過一代代的復制、交叉、變異操作，選出適應度最好的染色體，即獲得最優(yōu)解.

2.1 基因編碼

雖然層合板中鋪層順序不會改變材料的拉壓性能，但是會對承彎性能造成影響，所以仍需對鋪層進行排序.規(guī)定每層厚度T，用一定順序的數(shù)字串表示按此順序排布的鋪層.采用對稱形式鋪設，僅將一半用基因表示，鋪層基因編碼從左向右對應鋪層從內向外.為保證鋪層對稱性，±45°鋪層成對出現(xiàn)，兩層厚度值均設定為T/2.如下所示，其中 1，2，3 分別對應鋪層角 0°，±45°和90°，0 表示對應層不鋪:

對于給定鋪層而言，在何位置插入0，對鋪層本身不產生影響，但是破壞了基因-鋪層之間的一一對應關系.為防止這種情況的發(fā)生，0只存在于隊列末尾.如變異或交叉產生編碼形式不符合條件，自動傳遞至末尾，其后鋪層向前靠一位.

為實現(xiàn)多區(qū)域鋪層的全局優(yōu)化，對于每一個區(qū)域提供不同的鋪層建立獨立的基因編碼Nr，從而建立對于全局的公共鋪層Nw.計算過程中，全局鋪層也按對稱形式，鋪設于區(qū)域層表面，全局鋪層基因設置于基因鏈末尾，該鋪層也參與優(yōu)化過程，以獲取最優(yōu)結果.為了防止區(qū)域間厚度變化過大影響工藝性，需要對鋪層厚度變化范圍給出限定.限定 Nr和 Nw的最大值分別為 Nr，max和 Nw，max，最小值分別為 Nr，min和 Nw，min.如圖1所示 3 個區(qū)域，Nr，max=6，Nw，max=7 時的鋪層示意圖.如該基因所表示的鋪層1對應的區(qū)域1的實際鋪層(包含區(qū)域與全局鋪層，并且為對稱排布)應為

圖1 分級鋪層方式示意

區(qū)域鋪層與區(qū)域鋪層之間如果對應層鋪角不同，則采用錯位搭接方式，厚度變化不連續(xù)則采用逐層遞減方式過渡.由于過度區(qū)域強度增加，為降低分析復雜程度，減少計算時間，本文不再單獨對該區(qū)域進行力學分析.另外，過渡區(qū)域在質量上的增加可忽略不計.

2.2 遺傳算子

選擇與交叉:選擇過程采用蒙特卡羅方法，每個個體的選擇概率和其評估值成正比.

交叉:以一定交叉率pc進行交叉.為提高個體之間的交叉頻率，遵循多點交叉原則.父本與母本的對應區(qū)域基因相互交叉產生子代個體，不同區(qū)域之間不交叉.交點和選取的基因段均隨機確定.如圖2例舉了父代一種可能的交叉生成子代方式.

圖2 交叉過程

變異:以一定變異概率pm隨機抽取一個個體，選擇其中一個基因進行突變.根據(jù)給定的變異率進行變異操作，變異區(qū)域有相同概率增加或者減少一層.增加鋪層時，隨機選取其中一層，置于之前，鋪層角度隨機產生.若該區(qū)域鋪層數(shù)已達最大值，則刪去最后一層.減少鋪層時，隨機選取其中一層將該層刪去，若該區(qū)域鋪層數(shù)已減至最小值，則不再進行刪去操作.圖3顯示了當區(qū)域變異后層數(shù)超過限定層數(shù)時的操作過程.

圖3 變異過程

2.3 精英算法改型［13］

為了保存每代最優(yōu)的個體，使用的是保護精英的方法.交叉變異操作生成的子代與親代混合排名，取前Nelite個作為精英，直接傳遞到下一代，剩下的名額由之前產生的子代排名靠前的補滿.第一代個體基因隨機產生，第二代開始采用精英保護保存優(yōu)良個體.

3 評估函數(shù)設計

3.1 適應度函數(shù)

對應多目標優(yōu)化，適應度函數(shù)取權重函數(shù)，其形式為

其中，Dik為對應第i個個體的第k個權重參數(shù)的無量綱值;ωk為該參數(shù)的權重系數(shù)，總和為1;NF為適應度函數(shù)權重參數(shù)個數(shù).權重參數(shù)取設計過程中主要關注的變量，如結構形變、重量、整體剛度(應變能)等.Dik=dik/Sk，Sk根據(jù)初始設計對應計算結果給出.若無初始設計方案，則按給定最大層數(shù)下π/4板(0/±45/90)2計算.

3.2 適應度函數(shù)改進

對于計算過程中需要獲取最小適應度值的個體，需要對適應度作適當調整.另外，使用統(tǒng)一的適應度算式存在算法前期適應度很高的個體優(yōu)勢特別突出的問題，以致產生超級個體，造成早熟的缺陷，并且在后期個體適應度差異縮小時也會降低遺傳效率.綜合考慮以上兩方面，對適應度值作以下調整:

其中，b為Fi理想最小值;a為Fi平均值與最小值之差，保證平均值的適應度為0.5時Fi與b的距離.取β=2，為了起到調節(jié)作用，對α早期取0.5，中期取1，后期取2.

由于計算過程中，每一個算例的理想值較難估算，過大過小都會延緩遺傳進化速率，為此，本文取b為群體本代的最小值，使得每一代均能且僅能產生一個適應度為1的個體.這么做的缺點是使得每一代調整后的適應度值F*i與本代最大適應度值Fi，max相關，不能逐代比較F*i.另外，采用精英算法時，考慮了父子兩代的適應度，所以，將b修正為父子兩代群體F*i最小值.

3.3 懲罰函數(shù)

對于實際遺傳計算過程中，不滿足約束條件的個體，如果直接對其拋棄，會對種群的多樣性產生影響，這種情況在最優(yōu)解在邊界附近時尤為明顯.為此，采用懲罰不可行解的方法將約束問題轉化為無約束問題.對于給出的約束條件，選取乘法形式的評估函數(shù)

其中Pi為罰函數(shù).采用文獻［14］提出的懲罰函數(shù)，對于約束gk≤bk，計算式如下:

其中，Δbk=max{0，gk-bk};Δbmaxk=max{ε，Δbk}，為避免分母為0給出的小正數(shù);Np為給出的約束條件總數(shù).同樣需要注意，采用精英算法時，將bmaxk修正為父子兩代的違反量最大值.

3.4 其他計算考慮

隨著遺傳過程的進行，產生相同基因個體的概率越來越大.由于采用ABAQUS計算一個實例的時間較長，重復計算將對資源造成極大的浪費.為防止此類狀況的發(fā)生，對歷代每個已計算的個體基因和計算結果進行保存，每次進行模型計算前，對待解基因與已存基因進行逐個比對，若已進行過計算，則將計算過程跳過，直接輸出之前的結果，否則進行計算.為方便比較，將基因序列看成4進制數(shù)值，轉化成10進制值Ngen保存，如鋪層為［1，1，2，2，3，2，0］時，Ngen=5 816.若將整個基因鏈整體轉化保存，Ngen會隨著區(qū)域數(shù)量增加呈幾何倍數(shù)上升，不利于比較，因此需進行逐區(qū)域保存，逐一比較兩個個體對應區(qū)域的Ngen，若相同則判定為同一個個體.

4 算例

參照文獻［15］算例1并根據(jù)對稱原則截取一半，設計算例如下:一塊長寬為40 mm×15 mm的層合板，受力形式如圖4所示，兩長邊自由，一短邊AB限制x向位移和繞y軸旋轉，并受不同方向均布線載 f1和 f2，取 f1=20 N/mm，f2=50 N/mm.考慮到剪力在平板內基本不變，板受彎沿x方向逐漸變化，將其沿x方向分為3個區(qū)域進行優(yōu)化.

圖4 模型受力及分區(qū)示意圖(單位:mm)

設定單層鋪層厚度T=0.2 mm，ρ=1.56 g/mm3，設定適應度權重參數(shù)為結構重量和層合板A，B點間平均位移.對等厚度(6.4 mm)π/4板計算后，取去量綱因子 S1=5×10-6kg，S2=0.5 mm，取ω1=0.5，ω2=0.5.約束為各區(qū)域每個單元每個單層Tsai-Wu準則下不失效，安全因子取1.遺傳算法方面，設定最高遺傳100代，每代種群個數(shù)為100.由于使用了精英算法保證了最優(yōu)結果遺傳的完整性，所以取稍大的交叉率和遺傳率:pc=0.9，pm=0.05.另外，鋪層方面取 Nr，max=7，Nw，max=9，Nr，min=3，Nw，min=4.

分級鋪層優(yōu)化遺傳過程見圖5.層合板最優(yōu)結果在第 41 代出現(xiàn)，最優(yōu)鋪層為［(2，2，2，1，1，1，1)，(2，2，2，2，1，1，0)，(2，2，2，0，0，0，0)，(1，1，1，1，1，1，1，1，0)］，對應質量 m=4.786 579 ×10-3kg，位移，w=0.190 642 3 mm，最優(yōu)適應度 F=0.6693002，STW，max=0.2123366.

圖5 層合板優(yōu)化計算結果

總體上，由于剪切剛度對鋪層厚度不敏感，而外層的承彎效率更高，因此優(yōu)化結果中內層由±45°鋪層占據(jù)，外層均為0°，符合實際情況.剪切力沿板x方向基本不變，所以±45°層厚度基本相同.而區(qū)域1由于要承擔更大的彎曲應力，所以較區(qū)域2的0°鋪層所占的比例更大，區(qū)域3所承擔的彎曲應力最小，由于有全局鋪層的存在，所以不再出現(xiàn)0°鋪層，而層數(shù)減小至3層.由于算例中不存在沿y方向的外力，所以優(yōu)化結果中不含有90°鋪層.

為證明該鋪層形式為設計條件下的最優(yōu)解，現(xiàn)通過有條件的窮舉法證明.由上述分析可知，區(qū)域最優(yōu)形式鋪層分布為內層為±45°，外層為0°，0°在區(qū)域1中所占的比例最大，在區(qū)域3最小.另外，全局鋪層在區(qū)域鋪層最外層均為0°的情況下，均為0°.對窮舉滿足條件的105717種情況進行計算統(tǒng)計，證明該鋪層方案是最優(yōu)解.

若采用不分級僅分區(qū)優(yōu)化方式，由于將所有整體鋪層都拆分成了區(qū)域鋪層，相較于分級鋪層方式，增加了設計變量，因此在計算時間上增加了37%.計算結果是最優(yōu)鋪層在區(qū)域1中增加了1層0°鋪層，其他區(qū)域中與分級鋪層方式優(yōu)化結果相同，適應度值減少1%，差別較小，因此付出的計算效比不高.

對于同一個算例，其他參數(shù)保持不變，僅取ω1在0.1～0.9之間變化，得到分區(qū)優(yōu)化結果如表1所示.由表1可見，隨著質量權重的增加，各區(qū)域和全局鋪層數(shù)量逐漸減少，直至ω1=0.9時全體減至最低鋪層數(shù).ω1=0.7時區(qū)域2鋪層的增加是建立在全局鋪層減少2層的基礎之上的，綜合而言，該區(qū)域鋪層總數(shù)是減小的.另外，鋪層中0°鋪角比例逐漸增加，以便在有限質量內保證一定的彎曲剛度.所列最優(yōu)鋪層結果均未破壞，STW，max=0.3851623，出現(xiàn)在 ω1=0.9 處.

表1 分級鋪層優(yōu)化結果隨質量權重因素ω1變化

本文也對不分區(qū)的整體優(yōu)化方式進行了計算(即令區(qū)域數(shù)為1).為保證設計空間的一致性，設最大鋪層數(shù)為上例 Nr，max與 Nw，max之和，最小鋪層數(shù)為 Nr，min與 Nw，min之和，其他條件不變.圖6表示了分級優(yōu)化比整體優(yōu)化最優(yōu)適應度增加的百分比隨ω1的變化趨勢.整體上隨著質量權重系數(shù)的增加(ω1值增加)，兩者差距逐漸增大.由此可見，對質量要求越苛刻，則分級設計的優(yōu)勢越明顯.ω1=0.9處時百分比產生較大下降的原因是本文設置的最低層數(shù)影響了質量的進一步降低、區(qū)域間差異擴大的實現(xiàn)，事實上，在ω1=0.4處時，最低層數(shù)的限制已經(jīng)開始影響到了區(qū)域3的鋪層改變.

圖6 整體比分級優(yōu)化最優(yōu)適應度增量百分比隨ω1的變化

5 結論

本文建立了層合板剛度模型，實現(xiàn)了對層合板鋪層角和厚度的遺傳優(yōu)化.一方面，通過對分區(qū)鋪層進行單獨優(yōu)化和提出整體結構目標函數(shù)，考慮各個局部應力，能夠比整體鋪層優(yōu)化得到受力形式更加細致合理的鋪層分布.另一方面，由于一般情況下應力在結構上是連續(xù)分布的，全局鋪層的存在，體現(xiàn)了針對整體應力分布的鋪層設計，保護了該部分鋪層的完整性，也在一定程度上節(jié)約了計算資源.本文算例表明，分級鋪層遺傳優(yōu)化過程收斂較快，優(yōu)化結果合理可信，在應力分布較為復雜的層合板中更能發(fā)揮優(yōu)勢，更能達到降低結構質量、提高結構剛度等目的.另外本文優(yōu)化分析過程結合了有限元計算方法，使得即使應用復雜形狀的層板模型，其分析也變得簡單易行.

本文的計算過程中，區(qū)域是人為劃分的，全局鋪層、區(qū)域鋪層的最大層數(shù)和最小層數(shù)均是人為給定的，并沒有進行定量的規(guī)劃，影響了最優(yōu)結果性能的進一步提高.因此如何根據(jù)初步設計的應力分布選擇區(qū)域數(shù)量、進行自動區(qū)域劃分，如何根據(jù)應力變化幅度設定鋪層變化范圍，將成為進一步提高優(yōu)化能力的研究方向.

References)

［1］ Akbulut M，Sonmez F O.Design optimization of laminated composites using a new variant of simulated annealing［J］.Computers and Structures，2011，89(17/18):1712-1724

［2］晏飛，李為吉.基于自適應遺傳算法的復合材料層合板鋪層順序優(yōu)化方法［J］.西北工業(yè)大學學報，2001，19(1):156-159 Yan Fei，Li Weiji.A stacking sequence optimization method of composite laminated panel based on an adaptive genetic algorithm［J］.Journal of Northwestern Polytechnical University，2001，19(1):156-159(in Chinese)

［3］楊建國，李蓓智，俞蕾.基于免疫遺傳算法的優(yōu)化設計［J］.機械設計，2002，19(9)14-17 Yang Jianguo，Li Beizhi，Yu Lei.Optimization design based on genetic algorithm of Immunity［J］.Journal of Machine Design，2002，19(9):14-17(in Chinese)

［4］魯大偉，李書.應用免疫遺傳算法優(yōu)化設計層合板鋪層順序［J］.北京航空航天大學學報，2005，31(2):247-250 Lu Dawei，Li Shu.Stacking secuence optimization of composite laminate applying genetic algorithm of immunity［J］.Journal of Beijing University of Aeronautics and Astronautics，2005，31(2):247-250(in Chinese)

［5］任茶仙，張鐸.復合材料層合結構鋪層順序優(yōu)化設計的免疫遺傳算法［J］.強度與環(huán)境，2007，34(2):43-50 Ren Chaxian，Zhang Duo.Immunity genetic algorithm for stacking sequence optimization in composite laminated structures［J］.Structure ＆ Environment Engineering，2007，34(2):43-50(in Chinese)

［6］修英姝，崔德剛.復合材料層合板穩(wěn)定性的鋪層優(yōu)化設計［J］.工程力學，2005，22(6):212-216 Xiu Yingshu，Cui Degang.Plyoptimization design for stability of composite laminates［J］.Engineering Mechanics，2005，22(6):212-216(in Chinese)

［7］ Almeida F S，Awruch A M.Design optimization of composite laminated structures using genetic algorithms and finite element analysis［J］.Composite Structures，2009(88):443-454

［8］ Park C H，Lee W I，Han W S.Improved genetic algorithm for multidisciplinary optimization of composite laminates［J］.Comput Struct 2008，86(19/20):1894-1903

［9］喬巍，姚衛(wèi)星.復合材料加筋板鋪層優(yōu)化設計的等效彎曲剛度法［J］.計算力學學報，2011，28(2):158-162 Qiao Wei，Yao Weixing.Equivalent bending stiffness method for stacking sequence optimization of composite stiffened panel［J］.Chinese Journal of Computational Mechanics 2011，28(2):158-162(in Chinese)

［10］ Zehnder N，Ermanni P.A methodology for the global optimization of laminated composite structures［J］.Composite Structures，2006(72):311-320

［11］唐振杰，李玉龍.基于等效剛度矩陣的復合材料機翼盒段優(yōu)化設計［J］.南京航空航天大學學報，2011，43(4):521-525 Tang Zhenjie，Li Yulong.Optimization design for box of composite wing based on effective stiffness matrix［J］.Journal of Nanjing University of Aeronautics and Astronautics，2011，43(4):521-525(in Chinese)

［12］王耀先，復合材料結構設計［M］.北京:化學工業(yè)出版社，2001:83-87，34 Wang Yaoxian.Composite structure design［M］.Beijing:Chemical Industry Press，2001:83-87，34(in Chinese)

［13］ Almeida F S，Awruch A M.Design optimization of composite laminated structures using genetic algorithms and finite element analysis［J］.Composite Structures，2009(88):443-454

［14］ Gen M，Cheng R.Optimal design of system reliability using intervel programming and genetic algorithms［J］.Computers Ind Engng，1996，31(1/2):1-29

［15］ Huang Zhengming，Inelastic and failure analysis of laminate structures by ABAQUS incorporated with a general constitutive relationship ［J］.Journal of Reinforced Plastics and Composites，2007，26(11):1135-1181