易繼軍 李艷梅 榮見(jiàn)華 曾 韜

1.中南大學(xué)，長(zhǎng)沙，410083 2.長(zhǎng)沙理工大學(xué)，長(zhǎng)沙，410004 3.湖南水利水電職業(yè)技術(shù)學(xué)院，長(zhǎng)沙，410131

0 引言

結(jié)構(gòu)布局優(yōu)化設(shè)計(jì)包括結(jié)構(gòu)拓?fù)鋬?yōu)化、結(jié)構(gòu)形狀優(yōu)化和結(jié)構(gòu)尺寸優(yōu)化［1］。多數(shù)的優(yōu)化方法是在概念設(shè)計(jì)階段進(jìn)行拓?fù)鋬?yōu)化，在詳細(xì)設(shè)計(jì)階段進(jìn)行尺寸和形狀優(yōu)化。此類(lèi)優(yōu)化算法由于設(shè)計(jì)空間的隔離，可能無(wú)法獲得最優(yōu)解［2］。連續(xù)體結(jié)構(gòu)尺寸和拓?fù)浣M合優(yōu)化是在優(yōu)化運(yùn)算中同時(shí)考慮結(jié)構(gòu)尺寸和拓?fù)渥兞?，存在著變量?lèi)型和尺寸度量上的差異，直接對(duì)它們進(jìn)行同時(shí)處理有一定困難?，F(xiàn)有的組合優(yōu)化主要針對(duì)桁架結(jié)構(gòu)。Steven等［3］應(yīng)用漸進(jìn)優(yōu)化方法（ESO）對(duì)離散結(jié)構(gòu)進(jìn)行了拓?fù)浜统叽绲慕M合優(yōu)化。劉濤等［4］提出了一種桁架拓?fù)洹⑿螤詈统叽缃M合優(yōu)化的漸進(jìn)優(yōu)化方法，該方法采用拓?fù)鋬?yōu)化和形狀尺寸優(yōu)化兩個(gè)子問(wèn)題分層處理。Rajan等［5］利用遺傳算法（GA）分別進(jìn)行了6節(jié)點(diǎn)桁架和14節(jié)點(diǎn)桁架的優(yōu)化。Chan等［6］采用準(zhǔn)則法和遺傳算法的混合算法對(duì)鋼框架結(jié)構(gòu)進(jìn)行了拓?fù)浜统叽鐑?yōu)化。在桁架形狀與尺寸組合優(yōu)化方面，常用的方法可以分為數(shù)學(xué)規(guī)劃法［7］和漸進(jìn)優(yōu)化方法［4］兩類(lèi)。在連續(xù)體組合優(yōu)化方面，Afonso等［8］提出了基于準(zhǔn)則優(yōu)化法的盤(pán)類(lèi)結(jié)構(gòu)拓?fù)浜统叽绲姆旨?jí)優(yōu)化策略。王偉等［9］為解決機(jī)翼結(jié)構(gòu)布局優(yōu)化設(shè)計(jì)問(wèn)題，提出了一種兩級(jí)三層優(yōu)化設(shè)計(jì)方法。榮見(jiàn)華等［10］基于ESO方法，提出了板梁組合的連續(xù)體結(jié)構(gòu)優(yōu)化方法。Zhou等［2］提出了一種拓?fù)?、尺寸和形狀集成的?yōu)化算法，由于沒(méi)有考慮變量尺寸度量上的差異，提出的算法存在不足。

本文提出一種連續(xù)體結(jié)構(gòu)組合優(yōu)化方法，同時(shí)考慮了設(shè)計(jì)變量類(lèi)型和尺度上的差異，實(shí)現(xiàn)了拓?fù)浜统叽缱兞康膮f(xié)同優(yōu)化。

1 組合結(jié)構(gòu)拓?fù)渥兞考皟?yōu)化數(shù)學(xué)模型

在初始設(shè)計(jì)域中，整個(gè)結(jié)構(gòu)被分成兩大類(lèi)區(qū)域。在一類(lèi)子區(qū)域中，單元的尺寸（如厚度）作為設(shè)計(jì)變量，即尺寸設(shè)計(jì)變量；在另外一類(lèi)子區(qū)域中，進(jìn)行拓?fù)湓O(shè)計(jì)，拓?fù)湓O(shè)計(jì)單元的拓?fù)渥兞孔鳛樵O(shè)計(jì)變量。在尺寸設(shè)計(jì)區(qū)域，不能進(jìn)行拓?fù)湓O(shè)計(jì)的單元用p（p＝1，2，…，P）編號(hào)，尺寸設(shè)計(jì)變量用z（z＝1，2，…，H）編號(hào)，Pz和Sz分別表示第z個(gè)子區(qū)域的單元數(shù)量和面積，尺寸設(shè)計(jì)變量用hz表示。在拓?fù)湓O(shè)計(jì)區(qū)域，Q表示拓?fù)湓O(shè)計(jì)變量的個(gè)數(shù)，拓?fù)湓O(shè)計(jì)單元用q（q＝1，2，…，Q）編號(hào)，Sq表示其單元面積，它們的拓?fù)渥兞坑忙裶（q＝1，2，…，Q）表示。

1.1 拓?fù)渥兞亢统叽缱兞?/h3>

類(lèi)似于SIMP（solid isotropic microstructure with penalization）方法，采用材料屬性的有理近似模型（rational approximation of material properties，RAMP）的過(guò)濾方法，將離散的拓?fù)渥兞哭D(zhuǎn)換成［0，1］區(qū)間的連續(xù)變量。通過(guò)懲罰因子對(duì)設(shè)計(jì)變量在（0，1）之間的中間密度值進(jìn)行懲罰，使連續(xù)變量的拓?fù)鋬?yōu)化模型能很好地逼近傳統(tǒng)的0/1離散變量的拓?fù)鋬?yōu)化模型。使用RAMP過(guò)濾函數(shù)fk（ρiq）識(shí)別單元?jiǎng)偠?，用過(guò)濾函數(shù)fv（ρiq）識(shí)別單元體積，單元特性參數(shù)識(shí)別采用如下公式：

對(duì)于兩類(lèi)變量在尺寸度量上的差異，采用歸一化處理方法將尺寸變量轉(zhuǎn)化成［0，1］區(qū)間的連續(xù)變量?jī)?yōu)化問(wèn)題。即，用hz/代替厚度變量hz，這樣厚度變量的約束區(qū)間變?yōu)椋踙z/，1］。

1.2 優(yōu)化數(shù)學(xué)模型

位移約束的最小體積的布局優(yōu)化問(wèn)題可以表示為

2 近似優(yōu)化模型和求解方法

2.1 變位移約束限的序列近似優(yōu)化模型

為了使式（2）中的位移近似函數(shù)在每次迭代循環(huán)中保持近似成立，建立了變位移約束限的近似系列模型。這些變化的位移約束限起到了設(shè)計(jì)變量的置信區(qū)間的作用。通過(guò)這個(gè)方法，能夠保持設(shè)計(jì)變量的小量變化，從而確保迭代過(guò)程求解的有效性［9-10］。模型如下：

2.2 近似優(yōu)化模型的求解方法

定義

這里ju（i）表示在第j組載荷作用下，第i號(hào)單元的位移矢量。

對(duì)于式（5），省略常數(shù)項(xiàng)，則式（3）轉(zhuǎn)換成：

這里

將目標(biāo)函數(shù)用二階泰勒展開(kāi)近似代替，那么式（6）的求解能夠轉(zhuǎn)換成對(duì)下式的求解：

其中

根據(jù)對(duì)偶理論，類(lèi)似于文獻(xiàn)［11-12］，式（7）轉(zhuǎn)換成下面的對(duì)偶規(guī)劃問(wèn)題：

這里

這里采用交互式的循環(huán)迭代策略來(lái)近似求解式（7）和式（8）的最優(yōu)解（i＝1，2，…，Q＋H）和拉格朗日乘子λ。為了求解拉格朗日乘子，將目標(biāo)函數(shù)φ（λ）≈L（x＊，λ）進(jìn)行近似二階泰勒展開(kāi)。對(duì)L（x＊，λ）求關(guān)于λ的一階和二階導(dǎo)數(shù)：

對(duì)式（8）應(yīng)用K－T條件，有

并且，由式（11）的第一式和第三式，有

根據(jù)式（9）～式（15），容易得到下面的方程：

如果φ（λ）用其二階近似表達(dá)，同時(shí)省略常數(shù)項(xiàng)，我們能夠得到下面的二次規(guī)劃方程：

這里，H為L(zhǎng)×J階向量；D為（L×J）×（L×J）矩陣。

通過(guò)二次規(guī)劃法求解式（18），得到拉格朗日乘子。然后通過(guò)式（12）得到變量（i＝1，2，…，Q＋H）的新的近似值。用新的（i＝1，2，…，Q＋H）更新系數(shù)矩陣D和H。再次求解式（18），得到新的（i＝1，2，…，Q＋H）。這個(gè)循環(huán)過(guò)程一直反復(fù)進(jìn)行，直到滿(mǎn)足收斂條件 ‖x（s＋1）－x（s）‖/‖x（s）‖ ≤ε3（s是求解二次規(guī)劃的循環(huán)迭代次數(shù)）。

3 優(yōu)化策略和設(shè)計(jì)空間調(diào)整

3.1 優(yōu)化策略

在本文研究的方法中，整個(gè)優(yōu)化過(guò)程被分成兩個(gè)優(yōu)化調(diào)整階段和一個(gè)相變轉(zhuǎn)換階段。在第一個(gè)階段獲得包括灰色區(qū)域和尺寸改變的拓?fù)浣Y(jié)構(gòu)，當(dāng)設(shè)計(jì)結(jié)構(gòu)到達(dá)優(yōu)化目標(biāo)的鄰域，如果滿(mǎn)足下式優(yōu)化進(jìn)入相變轉(zhuǎn)換階段：

在此階段，進(jìn)行設(shè)計(jì)空間調(diào)整，移除部分灰色區(qū)域，然后優(yōu)化進(jìn)程進(jìn)入第二個(gè)優(yōu)化調(diào)整階段，在此階段，進(jìn)一步移除灰色區(qū)域，在滿(mǎn)足位移約束的前提下，如果滿(mǎn)足收斂準(zhǔn)則，則迭代終止。收斂準(zhǔn)則為

這里k表示當(dāng)前循環(huán)迭代步，τVmax是最大允許誤差。N是一個(gè)整數(shù)，在本文中N＝3。這意味著在6個(gè)收斂迭代步中達(dá)到穩(wěn)定［13］。

優(yōu)化流程如圖1所示。

圖1 優(yōu)化流程圖

3.2 設(shè)計(jì)空間調(diào)整

為了移除灰色區(qū)域，同時(shí)為了易于求解具有較大規(guī)模有限元網(wǎng)格模型的優(yōu)化問(wèn)題，這里采用在一些優(yōu)化迭代步將拓?fù)渥兞亢苄〉膯卧獜慕Y(jié)構(gòu)上去掉，而其他單元拓?fù)渥兞坎蛔兊牟呗?。即設(shè)置，使得＋ε（i＝1，2，…，Q），ε為一小量。把每個(gè)可設(shè)計(jì)的保留單元的拓?fù)渥兞颗c閾值相比較來(lái)確定單元是否存在的狀態(tài)，基于式（22），挑選用于單元?jiǎng)h除的候選單元，即將這些單元看成將要消失的單元。

為了確保優(yōu)化迭代中的結(jié)構(gòu)非奇異，同時(shí)結(jié)構(gòu)優(yōu)化時(shí)具有設(shè)計(jì)空間的擴(kuò)展（增添單元）功能，這里采用在結(jié)構(gòu)邊界和孔洞周?chē)郊尤斯げ牧蠁卧霓k法［11，14］來(lái)確保結(jié)構(gòu)剛度矩陣非奇異。將人工材料單元的彈性模量設(shè)為占據(jù)該空間的原始材料單元的彈性模量，僅將人工材料單元的拓?fù)渥兞咳棣裶。人工材料單元參與結(jié)構(gòu)參數(shù)計(jì)算和優(yōu)化迭代。根據(jù)人工材料單元拓?fù)渥兞康淖兓闆r，進(jìn)行設(shè)計(jì)空間的擴(kuò)展。設(shè)增添單元的閾值為，基于式（23），挑選出用于單元增添的候選單元，即將這些單元看成將要增添的單元。

設(shè)計(jì)空間縮減和擴(kuò)展的準(zhǔn)則表為

式中，nd、na分別為集合Nd和Na的單元個(gè)數(shù)，上述和都是小的經(jīng)驗(yàn)參數(shù)值。

如果式（24）成立，則執(zhí)行以下操作：

當(dāng)執(zhí)行式（25）和式（26）的操作后，非零的所對(duì)應(yīng)的單元材料特性參數(shù)編號(hào)變?yōu)楸Ａ魡卧奶匦詤?shù)編號(hào)，為零的所對(duì)應(yīng)的單元材料特性參數(shù)編號(hào)變?yōu)闊o(wú)材料單元所對(duì)應(yīng)的特性參數(shù)編號(hào)。這一輪的拓?fù)鋬?yōu)化就完成了設(shè)計(jì)空間的減縮和擴(kuò)展，否則，該輪拓?fù)鋬?yōu)化不進(jìn)行設(shè)計(jì)空間的減縮和擴(kuò)展。對(duì)于相變轉(zhuǎn)換階段的設(shè)計(jì)空間調(diào)整參照式（25）和式（26）進(jìn)行。

當(dāng)拓?fù)湓O(shè)計(jì)區(qū)域發(fā)生改變時(shí)，如果不調(diào)整厚度，將會(huì)使兩類(lèi)變量變化不一致，因此有必要擬定厚度調(diào)整方案。根據(jù)迭代過(guò)程中拓?fù)渥兓瘏^(qū)域和厚度變化區(qū)域體積小量變化以及相關(guān)變量在前后迭代步的連續(xù)性原則，可以建立以下近似式：

4 算例

4.1 算例1—— 立柱

初始的結(jié)構(gòu)模型和最大的設(shè)計(jì)域如圖2所示。整個(gè)結(jié)構(gòu)被分成了8個(gè)子區(qū)域，其中子區(qū)域6～8為拓?fù)湓O(shè)計(jì)區(qū)域，其單元采用四節(jié)點(diǎn)平面應(yīng)力單元，子區(qū)域1～5為尺寸變量設(shè)計(jì)區(qū)域，四節(jié)點(diǎn)平面應(yīng)力單元厚度作為設(shè)計(jì)變量，且在一個(gè)子區(qū)域中單元厚度變量相同。兩組載荷被施加到底部固定的結(jié)構(gòu)上。一組載荷P1＝500N，方向沿水平向右，施加到結(jié)構(gòu)左上角；另一組載荷P2＝500N，方向水平向左，施加到結(jié)構(gòu)右上角。最大設(shè)計(jì)區(qū)域1.0m×3.0m×0.002m （Lx＝1.0m，Ly＝3.0m）離散成50×150共7500個(gè)四節(jié)點(diǎn)平面單元，它們的初始厚度是0.002m，允許的最小厚度是0.0002m。位移約束點(diǎn)在結(jié)構(gòu)頂部的中間位置，水平方向的初始位移是0.153mm，允許值是0.400mm。彈性模量E＝210GPa，泊松比ν＝0.3。在第一階段單元?jiǎng)h除閾值和增加閾值都為0.001，在第二階段單元?jiǎng)h除閾值和增加閾值都為0.01，在階段轉(zhuǎn)換迭代步刪除閾值為0.1。

圖2 模型初始結(jié)構(gòu)和最大的設(shè)計(jì)域（算例1）

結(jié)構(gòu)的拓?fù)鋬?yōu)化歷程如圖3所示。在第70步滿(mǎn)足收斂條件式（20），優(yōu)化迭代進(jìn)入相變轉(zhuǎn)換階段，移除部分拓?fù)鋬?yōu)化的灰色區(qū)域并調(diào)整厚度，圖3b為執(zhí)行相變轉(zhuǎn)換后的結(jié)構(gòu)拓?fù)?。在?0步（圖3c）滿(mǎn)足收斂條件式（21），優(yōu)化迭代結(jié)束，其約束點(diǎn)的位移是0.400mm。

圖3 結(jié)構(gòu)的拓?fù)鋬?yōu)化歷程（算例1）

結(jié)構(gòu)的體積和約束點(diǎn)優(yōu)化歷程如圖4、圖5所示。

圖4 體積優(yōu)化歷程（算例1）

圖5 位移優(yōu)化歷程（算例1）

當(dāng)?shù)谝粋€(gè)優(yōu)化調(diào)整階段結(jié)束的時(shí)候，得到的拓?fù)錁?gòu)型接近最優(yōu)結(jié)構(gòu)，通過(guò)階段轉(zhuǎn)換步的操作，結(jié)構(gòu)設(shè)計(jì)區(qū)域中的黑/白分布清晰，最優(yōu)結(jié)構(gòu)布局在第二階段獲得。最優(yōu)解的結(jié)構(gòu)體積是0.002 538 7m3，結(jié)構(gòu)位移是0.400mm，結(jié)構(gòu)的厚度分布如表1所示。結(jié)構(gòu)的厚度優(yōu)化歷程見(jiàn)圖6。

表1 優(yōu)化后的厚度值 mm

圖6 厚度迭代歷程（算例1）

4.2 算例2——簡(jiǎn)支梁

初始的結(jié)構(gòu)模型和最大的設(shè)計(jì)域如圖7所示。整個(gè)結(jié)構(gòu)被分成了9個(gè)子區(qū)域，其中子區(qū)域7～9為拓?fù)湓O(shè)計(jì)區(qū)域，其單元采用四節(jié)點(diǎn)平面應(yīng)力單元，子區(qū)域1～6為尺寸變量設(shè)計(jì)區(qū)域，四節(jié)點(diǎn)平面應(yīng)力單元厚度作為設(shè)計(jì)變量，且在一個(gè)子區(qū)域中單元厚度變量相同。一組載荷P＝－2000N被施加到結(jié)構(gòu)中部（如圖7a所示），方向向下。結(jié)構(gòu)底部左端固定，右端簡(jiǎn)支。最大設(shè)計(jì)區(qū)域1.0m×3.0m×0.002m （Lx＝1.0m，Ly＝3.0m）離散成50×150×1共7500個(gè)四節(jié)點(diǎn)平面單元，它們的初始厚度是2mm，允許的最小厚度是0.2mm。位移約束點(diǎn)在結(jié)構(gòu)頂部的中間，垂直向下初始位移是－0.0647mm，允許值是0.300mm。在第一階段單元?jiǎng)h除閾值和增加閾值都為0.001，在第二階段單元?jiǎng)h除閾值和增加閾值都為0.01，在階段轉(zhuǎn)換迭代步刪除閾值為0.05。

圖7 模型初始結(jié)構(gòu)和最大設(shè)計(jì)域（算例2）

結(jié)構(gòu)的拓?fù)鋬?yōu)化歷程如圖8所示。在第135步滿(mǎn)足收斂條件（式（20）），優(yōu)化迭代進(jìn)入相變轉(zhuǎn)換階段，移除部分拓?fù)鋬?yōu)化的灰色區(qū)域并調(diào)整厚度，圖8b為執(zhí)行相變轉(zhuǎn)換后的結(jié)構(gòu)拓?fù)洹Ｔ诘?96步（圖8c所示）滿(mǎn)足收斂條件（式（21）），優(yōu)化迭代結(jié)束，其約束點(diǎn)的位移是0.300mm。

圖8 結(jié)構(gòu)的拓?fù)鋬?yōu)化歷程（算例2）

結(jié)構(gòu)的體積優(yōu)化和約束點(diǎn)優(yōu)化歷程如圖9、圖10所示。

圖9 體積優(yōu)化歷程（算例2）

圖10 位移優(yōu)化歷程（算例2）

當(dāng)?shù)谝粋€(gè)優(yōu)化調(diào)整階段結(jié)束的時(shí)候，得到的拓?fù)錁?gòu)型接近最優(yōu)結(jié)構(gòu)，通過(guò)階段轉(zhuǎn)換步的操作，結(jié)構(gòu)設(shè)計(jì)區(qū)域中的黑/白分布清晰，最優(yōu)結(jié)構(gòu)布局在第二階段獲得。最優(yōu)解的結(jié)構(gòu)體積是0.002 538 7m3，結(jié)構(gòu)位移是－0.300mm，結(jié)構(gòu)的厚度分布如表2所示。結(jié)構(gòu)的厚度優(yōu)化歷程見(jiàn)圖11。

表2 優(yōu)化后的厚度值 mm

圖11 厚度迭代歷程（算例2）

5 結(jié)論

本文研究了拓?fù)渥兞亢统叽缱兞拷M合的結(jié)構(gòu)優(yōu)化設(shè)計(jì)方法。采用材料屬性的有理近似模型的過(guò)濾方法，將離散的拓?fù)渥兞哭D(zhuǎn)換成［0，1］區(qū)間的連續(xù)變量；通過(guò)歸一化的處理方法使兩類(lèi)變量變成同一尺度。將組合優(yōu)化問(wèn)題轉(zhuǎn)化為單一類(lèi)型變量?jī)?yōu)化問(wèn)題。采用變位移約束限和二次規(guī)劃法求解優(yōu)化問(wèn)題，通過(guò)設(shè)計(jì)空間調(diào)整的方法移除灰色區(qū)域，同時(shí)調(diào)整厚度，得到了清晰的拓?fù)錁?gòu)型和優(yōu)化的厚度值。給出的算例說(shuō)明了本文方法的正確性和有效性。

［1］Rozvany G I N，Zhou M.Advances in Overcoming Computational Pitfalls in Topology Optimization［C］//6th AIAA/USAF/NASA/ISSMO Symposium on Multidisciplinary Analysis and Optimization.Reston，1996：1122－1132.

［2］Zhou M，Pagaldipti N，Thomas H L，et al.An Integrated Approach to Topology，Sizing and Shape Optimization［J］.Structural and Multidisciplinary Optimization，2004，26（5）：308－317.

［3］Steven G，Querin O，Xie M.Evolutionary Structural Optimisation （ESO）for Combined Topology and Size Optimisation of Discrete Structures［J］.Computer Methods in Applied Mechanics and Engineering，2000，188（4）：743－754.

［4］劉濤，鄧子辰.桁架結(jié)構(gòu)尺寸和形狀、拓?fù)涞臐u進(jìn)優(yōu)化方法［J］.西北工業(yè)大學(xué)學(xué)報(bào)，2004，22（6）：739－743.

［5］Rajan S D.Sizing，Shape，and Topology Design Optimization of Trusses Using Genetic Algorithm［J］.Journal of Structural Engineering，1995，121（10）：1480－1487.

［6］Chan C－M，Wong K－M.Structural Topology and Element Sizing Design Optimisation of Tall Steel Frameworks Using a Hybrid OC– GA Method［J］.Structural and Multidisciplinary Optimization，2008，35（5）：473－488.

［7］劉軍偉，姜節(jié)勝.桁架動(dòng)力學(xué)形狀優(yōu)化的統(tǒng)一設(shè)計(jì)變量方法［J］.振動(dòng)工程學(xué)報(bào)，2000，13（1）：84－88.

［8］Afonso S B，Belblidia F，Hinton E，et al.A Combined Structural Topology and Sizing Optimization Procedure for Optimum Plates Design［C］//European Congress on Computational Methods in Applied Sciences and Engineering.Barcelona，2000：125－131.

［9］王偉，楊偉，趙美英.大展弦比飛翼結(jié)構(gòu)拓?fù)?、形狀與尺寸綜合優(yōu)化設(shè)計(jì)［J］.機(jī)械強(qiáng)度，2008，30（4）：596－600.

［10］Rong J H，Jiang J S，Xie Y M.Evolutionary Structural Topology Optimization for Continuum Structures with Structural Size and Topology Variables［J］.Advances in Structural Engineering，2007，10（6）：681－695.

［11］榮見(jiàn)華，邢曉娟，鄧果.一種變位移約束限的結(jié)構(gòu)拓?fù)鋬?yōu)化方法［J］.力學(xué)學(xué)報(bào)，2009，41（3）：431－439.

［12］榮見(jiàn)華，張強(qiáng)，葛森，等.基于設(shè)計(jì)空間調(diào)整的結(jié)構(gòu)拓?fù)鋬?yōu)化方法［J］.力學(xué)學(xué)報(bào)，2010，42（2）：256－267.

［13］Huang X，Xie Y M.Optimal Design of Periodic Structures Using Evolutionary Topology Optimization［J］.Structural and Multidisciplinary Optimization，2008，36（6）：597－606.

［14］Stolpe M，Svanberg K.An Alternative Interpolation Scheme for Minimum Compliance Topology Optimization［J］.Structural and Multidisciplinary Optimization，2001，22（2）：116－124.