李家春 張 萍

貴州大學(xué)，貴陽，550003

0 引言

工程結(jié)構(gòu)的設(shè)計常常受制于材料屬性、幾何參數(shù)、作用載荷及其他參量的不確定性，處理不確定性問題的傳統(tǒng)方法是通過安全因子來進(jìn)行的。近年來，基于區(qū)間分析方法的工程結(jié)構(gòu)優(yōu)化設(shè)計逐漸得到重視并付諸于工程實踐［1－2］。

區(qū)間分析方法是在區(qū)間數(shù)學(xué)的基礎(chǔ)上發(fā)展起來的一種新的解決工程結(jié)構(gòu)不確定性問題的分析方法，它基于較少的數(shù)據(jù)（不確定參數(shù)的上下界）解決一系列不確定性問題，能得到令人滿意的結(jié)果。該方法作為隨機(jī)不確定模型的有益補(bǔ)充，具有較好的發(fā)展前景［3－5］。

在實際工程結(jié)構(gòu)設(shè)計中，需要合理地定量處理影響工程結(jié)構(gòu)性能的各種不確定性問題，將區(qū)間分析方法引入連續(xù)體拓?fù)鋬?yōu)化設(shè)計領(lǐng)域是解決此類問題的方法之一。邱志平等［6］以區(qū)間數(shù)學(xué)理論為基礎(chǔ)，提出了用區(qū)間方法進(jìn)行非線性結(jié)構(gòu)靈敏度分析的方法。本文將這種區(qū)間敏度分析思想引入到連續(xù)體結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計之中，建立能有效求解結(jié)構(gòu)最優(yōu)拓?fù)涞淖訁^(qū)間基移動算法。

1 基于雙向插值方法的區(qū)間攝動有限單元分析方法

盡管區(qū)間優(yōu)化問題的求解比確定性優(yōu)化問題的求解復(fù)雜，但區(qū)間分析方法比確定性分析方法得到的信息更多［7－8］。在連續(xù)體拓?fù)鋬?yōu)化中，敏度分析可采用區(qū)間分析方法進(jìn)行求解?；诿芏确椒ǖ倪B續(xù)體拓?fù)鋬?yōu)化敏度表示的是一種有限單元間的密度變化對目標(biāo)函數(shù)產(chǎn)生的不同影響。由于這種影響是一種相對影響，故采用區(qū)間分析方法來求解目標(biāo)函數(shù)或約束函數(shù)的相對敏度。這種相對敏度表征了在更廣的設(shè)計域上的目標(biāo)/約束函數(shù)的響應(yīng)。而基于區(qū)間的相對敏度是實數(shù)，故可利用成熟的規(guī)劃方法求解優(yōu)化問題。

1.1 基于雙向插值函數(shù)的區(qū)間優(yōu)化模型

在雙向插值函數(shù)模型中，假設(shè)材料的彈性張量是各向同性的。泊松比為常量，且與密度無關(guān)，而彈性模量隨相對密度ρ的變化而變化，為直觀起見，用變量xi代替i單元的相對密度ρ。

雙向插值模式的剛度矩陣Ki、柔度函數(shù)C可表示為

式中，φ（xi）為雙向插值函數(shù)；Ki0為第i單元的初始剛度矩陣；u為結(jié)構(gòu)位移；F為載荷。

φ（xi）的表達(dá)式為

式中，q為懲罰參數(shù)。

總剛度矩陣K為

彈性模量Ei為

式中，Ei0為上一迭代步的第i單元彈性模量。

雙向插值模型中結(jié)構(gòu)單元彈性模量的控制參數(shù)是xi和懲罰參數(shù)q。q取不同值時，中間密度單元xi的彈性模量參數(shù)有逼近0或1的趨勢。

設(shè)計變量x＝ （x1，x2，…，xn）T（0≤xi≤1），寫為區(qū)間向量形式：

在實際優(yōu)化分析中，不希望每次敏度分析的所用區(qū)間都是設(shè)計變量的原始區(qū)間，為此可以構(gòu)建既能體現(xiàn)局部信息，又能表征優(yōu)化問題全局性的子區(qū)間可變帶寬基移動優(yōu)化策略。其基本思想是：在初始計算分析時，在設(shè)計變量區(qū)間內(nèi)選一初始子區(qū)間：

根據(jù)區(qū)間分析方法進(jìn)行敏度分析計算，用相應(yīng)的標(biāo)準(zhǔn)優(yōu)化方法求設(shè)計區(qū)間中心（基）的移動xk＋1，c＝xk，c＋dk，k為迭代次數(shù)，dk為步長。

1.2 基于雙向插值模式的敏度分析攝動近似計算方法

以拓?fù)鋬?yōu)化問題中典型的柔度函數(shù)為目標(biāo)函數(shù)，以體積約束為約束條件，分析基于區(qū)間的相對拓?fù)涿舳取?/p>

體積函數(shù)表達(dá)式為

式中，Vi為結(jié)構(gòu)第i單元的體積；V＊為初始結(jié)構(gòu)體積；f為體積分?jǐn)?shù)。

彈性靜力結(jié)構(gòu)控制方程為

由雙向插值方法可知，剛度矩陣K＝K（x），x∈xI，是設(shè)計變量的直接函數(shù)，而u為設(shè)計變量的間接函數(shù)，故位移響應(yīng)可表示為

如果在優(yōu)化過程中載荷F保持不變，則柔度函數(shù)對設(shè)計變量xi的相對區(qū)間敏度可通過對設(shè)計變量xi的合適區(qū)間擴(kuò)展得到，敏度求解可表示為

式中，rad（·）為求取區(qū)間半徑的函數(shù)。

上述設(shè)計變量中，只有是區(qū)間參數(shù)，其他設(shè)計變量x1c，x2c，…，x（i－1）c，x（i＋1）c，…，xnc皆為實參數(shù)。

要求解柔度函數(shù)的敏度信息，關(guān)鍵是如何正確高效地求解Δ。令設(shè)計向量x在第i設(shè)計變量上的區(qū)間擴(kuò)展為x（i＋1）c，…xnc），則K（x）u＝F的矩陣區(qū)間擴(kuò)展為區(qū)間參數(shù)方程為

根據(jù)結(jié)構(gòu)力學(xué)單元剛度組裝的原則，將每個單元剛度矩陣擴(kuò)展為維數(shù)與總剛度矩陣相同的子矩陣，其實現(xiàn)方法是將單元剛度矩陣的值對應(yīng)賦予擴(kuò)展單元剛度矩陣相對位置，其余元素為0。如此便得

式中，Kc為標(biāo)稱剛度矩陣，δKI為攝動區(qū)間剛度矩陣。

由結(jié)構(gòu)有限元靜力攝動理論，標(biāo)稱結(jié)構(gòu)有限元控制方程為

相應(yīng)的攝動控制方程為

分解式（15），得

由式（16）有

將式（16）代入式（10）中的柔度函數(shù)相對敏度計算公式得

從式（18）可以看出由于Kc、uc、FT只與標(biāo)稱結(jié)構(gòu)有關(guān)，而δKI、Δxi只與結(jié)構(gòu)第i單元區(qū)間參數(shù)有關(guān)，簡化高階項，可將式（18）簡化為

2 基于區(qū)間敏度分析方法的拓?fù)鋬?yōu)化算法

在構(gòu)建連續(xù)體拓?fù)鋬?yōu)化算法的過程中，采用了固定帶寬的敏度分析子區(qū)間，優(yōu)化的參數(shù)為區(qū)間的中心值。換言之，就是在優(yōu)化過程中，用一定帶寬的盒子對整個設(shè)計域進(jìn)行搜索，以確定最優(yōu)設(shè)計變量所在的區(qū)間，盒子的搜索方向即基移動方向由盒子上的區(qū)間敏度決定。合適的固定帶寬可使優(yōu)化問題得到收斂全局解，優(yōu)化過程中帶寬可變，并采用一定策略控制帶寬按一定尺度逐漸減小，構(gòu)建變帶寬基移動優(yōu)化方法。

2.1 優(yōu)化模型

以結(jié)構(gòu)總體柔度為拓?fù)鋬?yōu)化的目標(biāo)函數(shù)，將結(jié)構(gòu)的體積作為優(yōu)化的約束條件，在給定載荷和邊界條件基礎(chǔ)上，基于雙向插值的區(qū)間優(yōu)化模型為

式（20）中，為避免總剛度矩陣奇異，取單元最小相對密度xmin＝0.0015。

取優(yōu)化初始點為

初始設(shè)計變量帶寬取Δx10＝Δx20＝…Δxn0＝0.1。第k次迭代計算時目標(biāo)函數(shù)的敏度可根據(jù)式（19）計算。

2.2 優(yōu)化算法

基于雙向插值的區(qū)間敏度分析連續(xù)體拓?fù)鋬?yōu)化準(zhǔn)則設(shè)計變量的迭代策略如下：

其中，xi（k＋1）c為第k＋1次迭代第i區(qū)間變量中心值，λk為第k次迭代時體積約束的拉格朗日乘子，拉格朗日乘子λk按罰乘子法思想進(jìn)行迭代計算：

式中，V（k）為第k次迭代設(shè)計域體積；m為增長系數(shù)，取m＝2～4。計算時取初始值：λ0＝1，r（0）＝0.25。

3 基于區(qū)間敏度分析方法的拓?fù)鋬?yōu)化算法數(shù)值實現(xiàn)

3.1 懸臂梁優(yōu)化數(shù)值實現(xiàn)

在拓?fù)鋬?yōu)化鄰域，懸臂梁設(shè)計問題被廣泛用作優(yōu)化問題研究時的測試、檢驗問題。以柔度最小化為目標(biāo)函數(shù)，保留25%的體積。材料彈性模量E＝200GPa，泊松比υ＝0.3。梁的左端固支約束，右端受垂直向下作用力F，模型可近似為二維平面應(yīng)力問題。結(jié)構(gòu)離散為60×30四節(jié)點四邊形單元。圖1為結(jié)構(gòu)模型示意圖。圖2所示為優(yōu)化結(jié)果，結(jié)構(gòu)優(yōu)化邊界從模糊到清晰。圖3所示為目標(biāo)函數(shù)值收斂情況。

圖1 懸臂梁有限元模型

圖2 雙向插值區(qū)間敏度拓?fù)鋬?yōu)化算法優(yōu)化結(jié)果

3.2 MBB梁優(yōu)化數(shù)值實現(xiàn)對比驗證

圖3 目標(biāo)函數(shù)收斂情況

MBB梁的拓?fù)鋬?yōu)化在歐洲空中客車底扳支架設(shè)計中得到應(yīng)用［9－11］，MBB梁有各種不同的變化形式。如圖4所示，一個邊緣固定（不可設(shè)計）的簡支梁，其內(nèi)部為可設(shè)計區(qū)域，在規(guī)定載荷下以結(jié)構(gòu)柔度最小化為目標(biāo)函數(shù)的優(yōu)化拓?fù)淙鐖D所示。應(yīng)用本節(jié)區(qū)間敏度分析方法的優(yōu)化結(jié)果如圖5所示，表示結(jié)構(gòu)優(yōu)化邊界由模糊逐漸變清晰過程。

圖4 MBB梁優(yōu)化模型及優(yōu)化結(jié)果

圖5 MBB梁拓?fù)鋬?yōu)化結(jié)果及目標(biāo)函數(shù)迭代過程

上述算例數(shù)值結(jié)果表明，將區(qū)間分析方法應(yīng)用到拓?fù)鋬?yōu)化領(lǐng)域敏度分析中完全可行，結(jié)合雙向插值模式的拓?fù)鋬?yōu)化具有無需過濾、全局收斂和邊界清晰的優(yōu)點。

4 結(jié)語

本文研究區(qū)間分析方法特性，分析了區(qū)間分析方法在包含界約束拓?fù)鋬?yōu)化問題中使用可行性；根據(jù)區(qū)間敏度分析優(yōu)點與傳統(tǒng)敏度分析的不足，將區(qū)間敏度分析思想引入結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計，建立了子區(qū)間基（中心）移動拓?fù)鋬?yōu)化方法，并用數(shù)值實驗驗證方法的正確性及實用性。本文結(jié)論為拓?fù)鋬?yōu)化研究提供了一種新的研究思路和方法。

對于載荷等設(shè)計變量具有不確定性的結(jié)構(gòu)，將區(qū)間分析方法與拓?fù)鋬?yōu)化方法相結(jié)合，研究結(jié)構(gòu)可靠性拓?fù)鋬?yōu)化、結(jié)構(gòu)失效模式拓?fù)鋬?yōu)化、智能結(jié)構(gòu)拓?fù)鋬?yōu)化將是很有前景的研究課題。

［1］李新蘭，姜潮，韓旭.基于區(qū)間的不確定多目標(biāo)優(yōu)化方法及應(yīng)用［J］.中國機(jī)械工程，2011，22（9）：1100－1106.

［2］蔡新，昊威，昊黎華，等.土石壩結(jié)構(gòu)區(qū)間優(yōu)化設(shè)計河［J］.河海大學(xué)學(xué)報，1998，26（3）：1－5.

［3］李方義，李光耀，鄭剛.基于區(qū)間的不確定多目標(biāo)優(yōu)化方法研究［J］.固體力學(xué)學(xué)報，2010，31（1）：86－93.

［4］吳杰.不確定參數(shù)結(jié)構(gòu)的區(qū)間優(yōu)化方法［J］.機(jī)械設(shè)計與制造，2008（1）：218－220.

［5］王海軍，曹德欣，李蘇北，等.非線性等式約束全局優(yōu)化問題的區(qū)間算法［J］.中國礦業(yè)大學(xué)學(xué)報，2003，32（2）：204－208.

［6］邱志平，王曉軍.結(jié)構(gòu)靈敏度分析的區(qū)間方法［J］.兵工學(xué)報，2005，26（6）：798－802.

［7］邱志平，顧元憲.有界不確定參數(shù)結(jié)構(gòu)位移范圍的區(qū)間攝動法［J］.應(yīng)用力學(xué)學(xué)報，1999，16（1）：1－10.

［8］邱志平，顧元憲.有界不確定性參數(shù)結(jié)構(gòu)靜力位移范圍的區(qū)間參數(shù)攝動法［J］.兵工學(xué)報，1998，19（3）：255－258.

［9］Diaz A，Sigmund O.Checkerboard Patterns in Layout Optimization［J］.Structural Optimization，1995，10：40－45.

［10］Jog C S，Haber R B.Stability of Finite Element Models for Distributed－parameter Optimization and Topology Design［J］.Computer Methods in Applied Mechanics and Engineering，1996，130：203－226.

［11］Eschenauer H，Olhoff N.Topology Optimization of Continuum Structures：a Review［J］.Appl.Mech.Rev.，2001，50（4）：331－390.