靜動載約束桁架結(jié)構(gòu)拓撲優(yōu)化的微分演化算法

胡長遠，傅長榮，季光耀

（1.麗水學院工學院，浙江 麗水 323000；2.浙西南地質(zhì)災(zāi)害研究所，浙江 麗水 323000）

結(jié)構(gòu)拓撲優(yōu)化是指以連續(xù)體開孔的數(shù)量及位置，或以離散體結(jié)構(gòu)的節(jié)點布局、節(jié)點間的桿件連接關(guān)系為研究對象，通過改變結(jié)構(gòu)的受力體系尋求最優(yōu)的布局方案。拓撲尋優(yōu)雖然可以最大程度優(yōu)化結(jié)構(gòu)，但求解空間和可行域在優(yōu)化過程中都在不斷更迭，這給結(jié)構(gòu)重分析帶來困難。較之尺寸優(yōu)化和形狀優(yōu)化，拓撲優(yōu)化具有更大的設(shè)計空間、經(jīng)濟效益和復(fù)雜程度，故逐步成為優(yōu)化設(shè)計領(lǐng)域中的焦點和難點，受到更廣泛的關(guān)注［1-3］。

目前拓撲優(yōu)化方法多用于靜荷載優(yōu)化，而工程實際中結(jié)構(gòu)常受到動荷載的作用［4］。固有頻率是結(jié)構(gòu)的重要特征，共振狀態(tài)時固有頻率的影響變得異常重要。為了保證結(jié)構(gòu)的安全，須把固有頻率控制在一定范圍內(nèi)［5］。拓撲優(yōu)化設(shè)計靜載約束有桿件應(yīng)力、節(jié)點位移和歐拉屈曲等，動載約束有固有頻率，靜動荷載約束的疊加會使計算過程更為復(fù)雜、計算效率更低，而且難以收斂，無法得到較為理想的拓撲構(gòu)形［6］。桁架結(jié)構(gòu)拓撲優(yōu)化目前主要采用基結(jié)構(gòu)法，該方法以桿件截面面積為設(shè)計變量，通過逐步刪除截面面積為足夠小值的桿件來實現(xiàn)拓撲優(yōu)化，求解簡單方便，但存在奇異最優(yōu)解問題［7］。如果用很小的斷面積求解來代替零斷面求解將受到各種約束的桁架拓撲優(yōu)化問題的干擾，不僅會丟失可能存在的奇異最優(yōu)解，而且會增加不必要的分析。以往的解決辦法是二階段尺寸和拓撲優(yōu)化方法。該方法對尺寸和拓撲變量分別進行優(yōu)化，求得最優(yōu)拓撲后再進行尺寸優(yōu)化，求解相對簡單，拓撲構(gòu)形容易實現(xiàn)。然而，由于尺寸和拓撲變量在數(shù)學上不是獨立的，二階段優(yōu)化方法將真實搜索空間轉(zhuǎn)化為近似域，可能會剔除包含全局最優(yōu)的初始搜索域，大大增加了過早收斂的風險［8］。因此發(fā)展一些更加通用有效，能同時考慮尺寸和拓撲變量的全局優(yōu)化算法來解決拓撲優(yōu)化問題成為一種迫切的需要。

近年來，啟發(fā)式算法逐漸被引入到工程優(yōu)化領(lǐng)域，如遺傳算法（GA）、微粒群算法（PSO）、和聲搜索（HS）、教與學算法（TLBO）等［9］。它們不需要梯度信息，對目標函數(shù)和約束函數(shù)性態(tài)的要求較為寬松，且能兼顧靜動態(tài)約束，因而較傳統(tǒng)算法而言有更廣泛的應(yīng)用范圍。計算中可以處理設(shè)計變量離散、可行域狹小的問題，能有效地得到高度復(fù)雜非線性問題的全局最優(yōu)解，一些啟發(fā)式優(yōu)化方法已被應(yīng)用于拓撲優(yōu)化。Tejani 等［10］利用改進的共生生物搜索算法解決在多工況下的桁架拓撲優(yōu)化問題。Dehghani 等［11］在帝國競爭算法中采用新的同化策略求解考慮尺寸、形狀、拓撲變量的桁架優(yōu)化問題。Dang等［12］將改進差分進化算法用于多種變量組合的桁架結(jié)構(gòu)拓撲優(yōu)化。Wang 等［13］用基于黃金分割系數(shù)的人工蜂鳥算法求解桁架結(jié)構(gòu)的拓撲優(yōu)化問題。

自Storn 和Price 于1996 年提出微分演化算法（differential evolution，DE）之后，該方法在解決復(fù)雜的優(yōu)化問題上得到了很多關(guān)注［14］。微分演化算法是一種混合算法，它結(jié)合了遺傳算法的更大種群概念和進化算法的自適應(yīng)變異以及采用了貪婪選擇策略。這些特征使微分演化算法相比進化算法和遺傳算法具有更好的魯棒性和更快的收斂速度。該算法已被廣泛應(yīng)用于工程優(yōu)化、可靠度分析、光譜分析及無線電系統(tǒng)感應(yīng)等領(lǐng)域［15-18］。

在結(jié)構(gòu)工程中，微分演化算法應(yīng)用于桁架拓撲優(yōu)化方面的研究分析還很少。筆者將DE 算法應(yīng)用于考慮尺寸和拓撲變量耦合的具有多工況靜動荷載響應(yīng)約束的桁架結(jié)構(gòu)拓撲優(yōu)化，通過對典型算例的求解并與相關(guān)文獻的結(jié)果比較來評價筆者所提方法的有效性。

1 微分演化算法（DE）

一個包括t個參數(shù)的優(yōu)化問題可以用一個t維的向量來描述，該向量可以表示為：xi=（xi1，xi2，…，xit）T∈S，i=1，2，3，…，NP。其中S∈Rn為優(yōu)化問題的搜索空間；DE算法利用NP作為向量每一代的個體數(shù)。類似于遺傳算法，DE算法通過變異、交叉和選擇過程實現(xiàn)種群的更新進化。具體過程如下［15-16］。

1.1 變異過程

變異是為了保證種群的多樣性，同時用合適的參數(shù)變化來指導已有的目標向量在合適的時間內(nèi)達到一個更好的結(jié)果，從而保證了搜索的魯棒性。

1.2 交叉過程

與GA 算法相似，DE 算法中的個體經(jīng)過變異后也進行交叉操作。對于群體中第G+1代經(jīng)過變異過程后的向量個體按照式（2）進行交叉，將產(chǎn)生新的個體：

式中：j=1，2，…，t；rand（j）∈［0，1］，是t個0～1 相互獨立的隨機數(shù)中的第j個；randn（i）是隨機從集合{1，2，…，t} 中取得個體向量維度的序號；CR為交叉因子，將決定個體之間交叉的概率。

1.3 選擇過程

DE 算法采用和GA 算法不同的貪婪準則：通過比較由變異和交叉產(chǎn)生的子代個體和父代個體，選擇適應(yīng)值好的變量，即如果父代個體適應(yīng)值更優(yōu)將繼續(xù)保留在種群中，否則保留子代個體。選擇過程如式（3）所示：

2 桁架結(jié)構(gòu)拓撲優(yōu)化

2.1 桁架結(jié)構(gòu)拓撲優(yōu)化的數(shù)學模型

結(jié)合拓撲優(yōu)化的設(shè)計思想，構(gòu)建以結(jié)構(gòu)構(gòu)件尺寸和拓撲變量為設(shè)計變量，以結(jié)構(gòu)重量極小化為目標函數(shù)，受多工況荷載作用，考慮固有頻率、桿件應(yīng)力、節(jié)點位移、歐拉屈曲以及動態(tài)失穩(wěn)等約束的離散型結(jié)構(gòu)拓撲優(yōu)化的數(shù)學模型，表達式為：

式中：X是設(shè)計向量，A1，A2，…，Am為桿件截面積，B1，B2，…，Bm為拓撲變量；Ai、Ei、Li、ρi、σi、σcri和Bi分別是第i類設(shè)計變量的截面面積、彈性模量、桿件長度、密度、桿件應(yīng)力、歐拉屈曲應(yīng)力和拓撲值。設(shè)Alim為桿件的臨界面積，當Ai小于Alim時，Bi等于0，表示桿件刪除；當Ai大于等于Alim時，Bi等于1，表示桿件保留。bj和δj分別表示節(jié)點j的質(zhì)量和位移；ki為歐拉屈曲系數(shù)；fr是第r階固有頻率。上標 “comp” 表示抗壓， “max” 和 “min” 分別表示上下限值。

數(shù)學模型中的約束函數(shù)還包括結(jié)構(gòu)動態(tài)穩(wěn)定性檢查g6和有效性檢查g7，前者是為了避免不穩(wěn)定結(jié)構(gòu)的出現(xiàn)，后者則是為了保證荷載作用點和支承節(jié)點等基本節(jié)點不被刪除。由于結(jié)構(gòu)的拓撲構(gòu)形是由DE算法迭代時隨機操作產(chǎn)生的，難免會出現(xiàn)一些不穩(wěn)定的拓撲結(jié)構(gòu)，因此有必要對新生成的結(jié)構(gòu)進行穩(wěn)定性檢查。結(jié)構(gòu)動態(tài)穩(wěn)定性檢查分兩步進行：1）檢查結(jié)構(gòu)自由度是否違反格魯布里爾準則；2）檢查結(jié)構(gòu)的總剛度矩陣是否為正定矩陣。

采用懲罰函數(shù)法處理上述約束，如果沒有違反約束，懲罰函數(shù)F（X）為零；否則，按照以下標準執(zhí)行處罰：

其中：

式（12）中：q為違反約束的個數(shù)；pi和分別表示違反的約束值及其限值；β1和β2為系數(shù)，本文均取1.5。

2.2 基于DE算法的桁架結(jié)構(gòu)拓撲優(yōu)化程序

為了提高DE算法的拓撲優(yōu)化效率，采用以下改進策略和技巧：1）對于離散型結(jié)構(gòu)拓撲優(yōu)化，取截面下限Amin=-Amax，將搜索空間轉(zhuǎn)為約兩倍的設(shè)計變量限值，這樣會豐富拓撲構(gòu)形，有利于算法快速尋找到全局最優(yōu)拓撲方案。2）對DE 初始隨機生成以及經(jīng)迭代交叉變異后生成的拓撲構(gòu)形進行動態(tài)穩(wěn)定性和有效性檢查。若不滿足要求，則更新個體，不再進入有限元計算，減少計算量。3）計算過程中，一旦桿件被刪除，接續(xù)刪除其連接，重新構(gòu)建有限元模型，消除應(yīng)力矩陣的奇異性，避免不必要的分析。

結(jié)合DE算法及前述拓撲優(yōu)化數(shù)學模型，可以給出圖1所示的離散型結(jié)構(gòu)拓撲優(yōu)化設(shè)計流程圖。由算法原理及優(yōu)化設(shè)計流程圖可知，DE算法原理簡單，能方便地應(yīng)用到結(jié)構(gòu)的拓撲優(yōu)化設(shè)計中。

圖1 基于DE算法的離散型結(jié)構(gòu)拓撲優(yōu)化設(shè)計流程圖

3 算例分析

為驗證筆者方法對桁架結(jié)構(gòu)拓撲優(yōu)化設(shè)計的有效性，以下對3 個典型桁架結(jié)構(gòu)進行分析，通過改變桿件尺寸以及拓撲構(gòu)形獲得最小重量，并與其他算法的結(jié)果進行比較。對所有問題，歐拉屈曲系數(shù)（ki，i=1，2，…，m）和節(jié)點質(zhì)量（bj，j=1，2，...，n）分別取為4.0和5 kg。DE算法的參數(shù)選擇如下：F1=0.75，F(xiàn)=0.6，CR=0.85，20 桿和24 桿平面桁架算例種群規(guī)模NP 均取為100，72 桿空間桁架算例則取150。

3.1 24桿平面桁架

圖2 所示為24 桿平面桁架的基結(jié)構(gòu)，材料的彈性模量E=6.9×1010Pa，密度ρ=2 740 kg/m3，節(jié)點3上有雙向集中質(zhì)量500 kg。頻率約束為f1≥30 Hz，應(yīng)力約束為±172.43 MPa，節(jié)點5、6 豎向位移約束為δmax=10 mm。設(shè)計變量為連續(xù)變量，其下限為-40 cm2，上限為40 cm2，臨界截面面積為1 cm2。兩工況下節(jié)點荷載見表1。

表1 24桿平面桿桁架荷載工況

圖2 24桿平面桁架的基結(jié)構(gòu)

24桿平面桁架拓撲優(yōu)化的目標函數(shù)收斂曲線見圖3，優(yōu)化后桁架最優(yōu)拓撲布局見圖4。為了驗證筆者方法的有效性，與文獻［8，12，19-22］相比較，對比結(jié)果見表2。

表2 24桿平面桁架拓撲優(yōu)化結(jié)果比較

圖3 24桿平面桁架拓撲優(yōu)化的目標函數(shù)收斂曲線

圖4 24桿平面桁架的最優(yōu)拓撲形狀

由圖3可知：考慮靜動荷載多種約束時，24桿平面桁架拓撲優(yōu)化在迭代約150次后已經(jīng)收斂，收斂速度快，不易陷入局部最優(yōu)解。在表2 中，第1階固有頻率和節(jié)點6的豎向位移均已達到限值，優(yōu)化后的桁架總重量為118.895 5 kg，優(yōu)于文獻［8，12，19-22］的結(jié)果。根據(jù)表2中的優(yōu)化結(jié)果繪制其最優(yōu)拓撲布局，如圖4 所示，最優(yōu)拓撲布局與文獻［8，12，20，22］的結(jié)果相一致，桿件1～6，10～11，14，17～21，23以及節(jié)點1、4被刪除，均為穩(wěn)定性結(jié)構(gòu)。

3.2 20桿平面桁架

本例的基結(jié)構(gòu)由9 個節(jié)點和20 根桿件組成，如圖5所示，其中節(jié)點1和9為支座，節(jié)點4上作用有荷載F1和F2。材料彈性模量E=6.9×1010Pa，密度ρ=2 740 kg/m3，兩工況下節(jié)點荷載見表3。頻率約束為f1≥60 Hz，f2≥100 Hz，應(yīng)力約束為±172.43 MPa，位移約束以節(jié)點4 的豎向位移等于10 mm 和60 mm兩種情形進行計算。設(shè)計變量為連續(xù)變量，其下限為-100 cm2，上限為100 cm2，臨界截面面積為1 cm2。

表3 20桿平面桿桁架荷載工況

圖5 20桿平面桁架的基結(jié)構(gòu)

考慮不同位移約束的20桿平面桁架拓撲優(yōu)化目標函數(shù)曲線見圖6，優(yōu)化的最終拓撲形狀見圖7。筆者計算結(jié)果與相關(guān)文獻的比較見表4和表5。

表4 δ4y=10 mm時20桿平面桁架拓撲優(yōu)化結(jié)果比較

表5 δ4y=60 mm時20桿平面桁架拓撲優(yōu)化結(jié)果比較

圖6 20桿平面桁架拓撲優(yōu)化的目標函數(shù)收斂曲線

圖7 20桿平面桁架的最優(yōu)拓撲形狀

由圖6 可知：節(jié)點4 的豎向位移約束為10 mm和60 mm時，計算分別在迭代約40次和120次后完全收斂。表4 中，節(jié)點4 的豎向位移已達到限值。由表4和表5可知，兩種情況下筆者優(yōu)化后的桁架總質(zhì)量分別為315.913 0 kg和150.637 8 kg，略好于文獻［12］的結(jié)果，優(yōu)于文獻［4，19-22］的結(jié)果。對比表4和表5還可知，位移約束嚴苛時結(jié)構(gòu)的總重量會增加。圖7所示的最優(yōu)拓撲形狀中節(jié)點3、6、7被刪除，布局上具有對稱性，僅與文獻［4］不同。

3.3 72桿空間桁架

如圖8 所示的72 桿空間桁架基結(jié)構(gòu)，由20 個節(jié)點和72 根桿件組成，依據(jù)結(jié)構(gòu)對稱特性，桿件被分為16 組。材料彈性模量E=6.895×1010Pa，密度ρ=2 767.99 kg/m3，節(jié)點1～4 上有三向集中質(zhì)量2 270 kg。頻率約束為f1≥4 Hz，f3≥6 Hz，應(yīng)力約束為±172.375 MPa，節(jié)點1～4 的X和Y向位移限值為6.35 mm。設(shè)計變量為連續(xù)變量，其下限為-30 cm2，上限為30 cm2，臨界截面面積為1 cm2。兩工況下節(jié)點荷載見表6。

表6 72桿空間桿桁架荷載工況

圖8 72桿空間桁架的基結(jié)構(gòu)

靜動荷載約束下72桿空間桁架拓撲優(yōu)化目標函數(shù)收斂曲線見圖9，與相關(guān)文獻的比較列于表7，最優(yōu)拓撲形式見圖10。

表7 72桿空間桁架拓撲優(yōu)化結(jié)果比較

圖9 72桿空間桁架拓撲優(yōu)化的目標函數(shù)收斂曲線

圖10 72桿空間桁架的最優(yōu)拓撲形狀

圖9 結(jié)果顯示迭代約180 次時已找到最優(yōu)拓撲布局方案，評價函數(shù)趨于收斂。由表7 可知：筆者DE 算法得出的桁架總質(zhì)量為434.894 9 kg，優(yōu)于其他所有結(jié)果，最優(yōu)解受固有頻率約束控制。圖10最優(yōu)拓撲形狀保留了2層外圍系桿，4層水平系桿及全部豎桿和腹桿，豎桿截面積向下逐層增大，與其他文獻結(jié)果一致。

4 結(jié)論

DE 算法是一種啟發(fā)式搜索的混合進化算法，和其他進化算法相比，DE算法不易陷入局部最優(yōu)解，具有原理簡單、易于實現(xiàn)以及魯棒性強等優(yōu)點。以往用二階段優(yōu)化方法求解桁架拓撲優(yōu)化，雖然降低了求解難度與計算量，但割裂了解空間而難以保證得到全局最優(yōu)解。筆者在DE 算法基礎(chǔ)上，通過對每個桿件引入拓撲變量，便于實現(xiàn)同時考慮兩種不同性質(zhì)設(shè)計變量的拓撲優(yōu)化問題?？紤]到桁架結(jié)構(gòu)在實際使用時可能會因動載作用產(chǎn)生噪聲過大、共振等現(xiàn)象，需要在約束條件中引入固有頻率約束。雖然考慮靜動載約束耦合會造成求解難度的增加，但是更符合實際情況。此外，在優(yōu)化過程中對有限元模型進行修正，重新構(gòu)建剛度、質(zhì)量和荷載矩陣，解決了矩陣奇異和不必要分析給計算帶來的困難。數(shù)值分析結(jié)果表明：筆者方法很好地保持了所得最優(yōu)解的穩(wěn)定性、有效性以及較強的收斂性，可以成功地進行桁架結(jié)構(gòu)的拓撲優(yōu)化設(shè)計。