基于子結構的宏微結構協(xié)同優(yōu)化方法

吳紫俊 肖人彬

1.武漢紡織大學數(shù)字化紡織裝備湖北省重點實驗室，武漢，430020

2.華中科技大學人工智能與自動化學院，武漢，430074

0 引言

基于子結構方法的宏微結構一體化優(yōu)化方法中，利用子結構自由度縮減實現(xiàn)了結構宏微觀的一體化設計，但需要預先建立不同子結構構型的樣本[1]，宏觀結構優(yōu)化時只能在預定的微觀構型中進行插值，即限制了結構宏微觀的優(yōu)化空間[2]。

宏微結構的協(xié)同優(yōu)化設計是從宏微多個尺度尋找材料的最佳分布位置，組合成不同構型的微觀構型和宏觀結構，使得材料潛力充分發(fā)揮[3]。由于結構設計引入了微觀尺度，因而擴大了結構設計空間[4]，賦予了結構更多的功能，如最小頻率響應的結構[5-6]、具有良好可制造性的宏觀結構等[7]。

目前，周期性子結構的優(yōu)化方法主要有兩類[8-9]：基于均勻化方法的周期性結構設計[10-11]與基于宏觀結構周期性約束條件的設計方法。在基于均勻化方法[12-13]的周期性結構設計方法中，周期性位移邊界條件施加在單胞上[14]，并假設結構的應力和應變在宏觀結構上是周期性變化的[15-16]。而基于周期性約束的方法中，假設宏觀結構中任意兩個周期性結構中的單胞具有相同的單元密度和相等的靈敏度值。這兩類方法相比，盡管都能獲得相近的優(yōu)化構型，但基于周期性約束條件的設計方法的計算效率要低于均勻化方法。

作為結構拓撲優(yōu)化的熱點分支之一，周期性結構設計得到了許多關注。LI等[17-18]提出基于維度縮減的代理模型來設計周期性結構，降低了結構頻率帶隙計算的計算資源消耗。ZHENG等[19]基于不確定性負載的高斯分布提出了周期結構魯棒拓撲優(yōu)化框架，研究了周期結構設計空間小的結構對不確定性載荷的敏感度。HE等[20]針對不規(guī)則有限元網(wǎng)格提出了一種周期性結構優(yōu)化方法，用非結構設計點來清晰地描述周期結構的拓撲。RIVA等[21]提出了一種基于波在細胞界面反射和傳輸?shù)募兘馕龇椒ǎ膳袛嗨岢龅慕鉀Q方案的最優(yōu)性。CHEN等[22]提出了一種用于結構和蜂窩材料拓撲并行設計的移動等值面閾值(moving iso-surface threshold，MIST)方法，利用均勻化理論推導出了宏觀結構和微觀周期性單元兩個維度的物理響應函數(shù)。CHENG等[23]研究了周期性均質多孔材料板和均勻加筋實心板的最大面外屈曲載荷設計，并基于漸進均勻化方法和CAI方法建立了結構宏觀和微觀設計變量與靈敏度分析的關系。DAHLBERG等[24]引入周期性結構單元以獲得具有更高對稱性的結構，提出利用高對稱性設計導頻結構，大大減小了低頻電波傳輸時的漏波現(xiàn)象。

趙清海等[25]基于周期性約束提出了多材料結構穩(wěn)態(tài)熱傳導拓撲優(yōu)化設計方法，并結合Ordered-RAMP材料插值方法獲得了周期性熱傳導拓撲構型。焦洪宇等[26]利用固體各向同性材料懲罰(solid isotropic material with penalization，SIMP)方法，通過增加周期性約束條件，提出了各向同性微結構材料的周期性拓撲優(yōu)化方法，并在其隨后的研究中[27]，對循環(huán)對稱結構、梯形結構等進行了周期性布局優(yōu)化研究。杜義賢等[28]利用能量均勻化方法建立了基于宏觀力學性能的微觀點陣結構的優(yōu)化模型，以周期性單胞為研究對象，將材料用量和力學方程等作為約束條件，獲得了邊界清晰的周期性點陣結構；并在其后的研究中[10]，結合負泊松比表征吸能特性，提出了兼具吸能和承載特性的周期性宏細觀結構多尺度優(yōu)化方法。

在周期性結構的多尺度設計中，為在設計域中設計微結構構型的同時兼顧微結構在宏觀結構中的最優(yōu)分布[9，29]，本文提出了基于子結構法的周期性結構的宏微觀設計方法，并利用子結構法中的自由度凝聚來實現(xiàn)周期性單胞的矩陣規(guī)模縮減，提高其計算效率的同時，也可利用子結構的反求來完成周期性單胞結構與宏觀構型的協(xié)同設計。

與此同時，文獻[1]和文獻[2]中提出的基于自定義微結構構型的宏微結構一體化設計方法ARSP(approximation of reduced substructure with penalization)，由于限定了微結構構型的設計范圍，其宏觀結構與自定義微結構構型之間存在兼容匹配問題。本文為突破預定義微結構構型的設計限制，利用子結構的自由度縮減與反求來實現(xiàn)兩個尺度上結構的協(xié)同設計，有助于實現(xiàn)大規(guī)模個性化設計[30]。

1 子結構與優(yōu)化插值模型

1.1 子結構凝聚與反求

在優(yōu)化子結構中，單個子結構的本構方程可表示為[31]

(1)

由式(1)的第二式可得子結構內(nèi)部節(jié)點位移向量：

(2)

將式(2)代入式(1)的第一式，可得方程：

(3)

式(3)可簡單表示為

(4)

(5)

本文中，將利用式(5)進行子結構凝聚組建宏觀本構方程的過程稱為子結構凝聚，通過凝聚后的本構方程進行結構設計稱為宏觀結構設計；利用式(2)計算子結構內(nèi)部自由度的位移的過程稱為子結構反求，通過子結構反求獲得對應的周期性單胞的過程稱為微觀結構設計[32]。

1.2 優(yōu)化插值模型的構建

在周期性宏觀結構設計中，假設宏觀結構由若干個子結構組成，如圖1所示。

圖1 周期性子結構與宏觀結構

根據(jù)子結構可知，宏觀結構設計域Ω的計算可直接利用凝聚后的剛度矩陣表示：

KsUs=Fs

(6)

其中，Ks為子結構凝聚后所形成的整體剛度矩陣；Us、Fs為凝聚后子結構交界面上節(jié)點的位移向量及載荷列陣。

針對宏觀結構設計，設計變量設置為子結構所含的材料體積，由于宏觀結構由周期性結構組成，故只需要判斷子結構材料的含量。此時宏觀子結構材料含量設計變量的取值為(0，1]，對應的剛度矩陣可表示為

(7)

其中，ρ為設計變量；p為ρ的懲罰因子。

子結構密度決定了周期性子結構構型的材料含量。微觀子結構構型設計中，設計密度為子結構中各個單元的材料含量，利用式(2)可獲得子結構的內(nèi)部節(jié)點值，通過優(yōu)化變量更新即可獲得微觀結構的構型。假設子結構內(nèi)部單元的剛度矩陣為K0，它與子結構優(yōu)化密度變量之間的關系可表示為

Kb=(ρb)qK0

(8)

其中，ρb為該子結構內(nèi)部單元的材料相對含量；q為ρb的懲罰因子。

2 基于周期性微觀單胞的宏微觀結構優(yōu)化

2.1 宏微觀優(yōu)化模型定義

以結構柔度最小為優(yōu)化目標，實現(xiàn)宏微觀結構的協(xié)同優(yōu)化：

(9)

式中，X為宏觀結構中所有子結構的密度集合；c為宏觀結構柔度設計目標函數(shù)；V為宏觀結構的材料含量總量；K、U、F分別為整體剛度矩陣、位移場和邊界約束；υ為體積約束；ρi為各子結構設計變量；ρmin是為避免奇異問題設定的最小密度，ρmin=0.001。

基于有限元子結構方法，優(yōu)化問題(式(9))在宏觀層面的優(yōu)化可表示為

(10)

微觀結構是構型相同的周期性結構，優(yōu)化模型可表示為

(11)

式中，cb為子結構柔度優(yōu)化目標；K0、Uj，i分別為第i個子結構中單元的單元剛度矩陣和位移場；Fbb，i為當前子結構凝聚后的邊界自由度對應的位移場；m為當前子結構中單元總數(shù)；ρb為微結構構型的設計變量；Vb為當前子結構的體積約束，此時該約束是變動的，與子結構在宏觀結構的分布數(shù)量相關。

通過上述宏微觀的優(yōu)化模型可知，結合SIMP方法計算，宏觀結構的拓撲構型由宏觀子結構的材料含量確定。而在微觀結構設計中，微觀周期性結構設計的材料約束不是固定值，而是根據(jù)宏觀結構的子結構的材料密度分布進行更新。當宏觀結構中周期性子結構分布在整個設計域時，微觀結構的材料含量約束與宏觀材料約束相同；當周期性結構分布在宏觀設計域的子域時，微觀結構的材料含量約束將增大。

2.2 宏微觀兩尺度的靈敏度分析

宏觀結構由構型相同的周期性結構組成，假設對式(10)的目標函數(shù)計算密度導數(shù)：

(12)

微觀結構的靈敏度可表示為

(13)

微觀上靈敏度與宏觀子結構相關，結構優(yōu)化過程中，在微觀結構設計變量增加宏觀的懲罰因子項可加快其收斂速度。

在基于子結構的優(yōu)化方法中，兩個維度上的變量更新策略均采用傳統(tǒng)的優(yōu)化準則法(optimality criteria method，OC)[33]：

(14)

在上述更新策略中，ε是變量迭代步長。η=0.5，而優(yōu)化條件Di可由子結構的目標函數(shù)和體積導數(shù)定義：

(15)

式中，κ為拉格朗日因子；?V/?ρi=vi為每個子結構的體積約束。

2.3 周期性結構的宏微觀拓撲優(yōu)化算法

在該宏微觀結構協(xié)同優(yōu)化中，通過子結構自由度凝聚和反求計算，可以獲得設計域上所有子結構邊界節(jié)點和內(nèi)部節(jié)點的位移場。通過式(12)和式(13)計算對應的宏觀結構和微觀結構的靈敏度。由于微觀周期性結構的構型會影響宏觀子結構凝聚的剛度矩陣，因此在進行變量更新時，應先更新周期性結構的設計變量，然后更新宏觀結構設計變量。對應的優(yōu)化設計流程如圖2所示，關鍵算法如下所示。

DeginMacroStructure()∥宏觀結構設計 while(判斷是否收斂) { substrcutureMacro(網(wǎng)格劃分參數(shù));∥根據(jù)優(yōu)化參數(shù)獲得凝聚的子結構有限元模型 calcMacroModel(邊界條件);∥計算宏觀設計域的有限元模型 for 1:N iSensitivityMacro();∥計算宏觀結構每個子結構的設計變量靈敏度 updatingMacro();∥根據(jù)靈敏度值更新設計變量 DeginMicroStructure();∥微觀結構設計 }DeginMicroStructure()∥微觀結構設計

while(判斷是否收斂) { for 1:N { calcsubstrcutureMicro(子結構節(jié)點編號);∥反求子結構內(nèi)部節(jié)點位移場 iSensitivityMicro();∥計算宏觀結構每個子結構的設計變量靈敏度 } updatingMicro();∥根據(jù)靈敏度值更新微結構設計變量 }

圖2 結構設計流程

該協(xié)同優(yōu)化計算框架中，每一步迭代均會判斷宏微觀結構是否收斂，而微觀結構優(yōu)化過程中，其材料約束受宏觀子結構分布限制，因此，微觀結構優(yōu)化中需要從宏觀結構中更新微觀材料約束。

3 周期性結構優(yōu)化實例

3.1 懸臂梁結構的設計(實例1)

以懸臂梁的計算為例，如圖3所示，設計域大小為L1×L2=2×1，其左側固定，右上角承受F=1 N的垂直向下壓力。所選用材料的彈性模量為E0=1，泊松比υ0=0.3，體積約束為0.4。Nx與Ny分別是水平與垂直方向上的子結構數(shù)目。

圖3 優(yōu)化模型(實例1)

設計域劃分為320×160個四邊形平面單元，并選擇不同大小的子結構進行宏觀結構凝聚。宏觀結構設計中，懲罰因子p設置為2，過濾半徑設置為與單元邊長相等；微觀結構設計中，懲罰因子q設置為3，過濾半徑設置為1.1。具體參數(shù)設置如表2所示?；赟IMP方法，不考慮宏微觀協(xié)同設計時，該模型的優(yōu)化構型和柔度值c如圖4所示。

表2 宏微觀優(yōu)化設計參數(shù)

圖4 基于傳統(tǒng)SIMP方法的優(yōu)化設計結構

在不同子結構大小下的優(yōu)化結果如圖5所示。從優(yōu)化結果可看出，設計域劃分的子結構越多則結構桿系特征越明顯，隨著子結構數(shù)量的增多，其宏觀結構的柔度值先增大后減小。

(a)Nx×Ny=2×1, (b)Nx×Ny=8×4,c=79.72 c=90.03

微觀結構優(yōu)化過程中，隨著宏觀結構的桿系特征越來越明顯，材料逐漸向微結構中集中，使得微觀結構周期性子結構中的材料增多，微結構的桿系特征越來越不明顯。

隨著設計域中的子結構數(shù)量增多，通過式(11)可得出宏觀結構的材料會逐漸集中到各個子結構中，使得周期子結構單胞中的材料逐步增多。對比圖5c與圖5d，宏觀結構的桿系特征越明顯，材料在子結構單胞中的集中程度越大。由此可以推導出，在利用所提出的周期性子結構設計方法進行結構設計時，需要平衡宏觀結構劃分和周期子結構單胞大小。

3.2 固端梁的結構優(yōu)化設計(實例2)

如圖6所示的固端梁結構優(yōu)化問題中，設計域大小為L1×L2=2×1，其左右兩側固定，設計域中心受F=1 N的載荷。宏觀結構劃分為Nx×Ny=8×4個子結構，子結構大小m×m分別設置為10×10、20×20、30×30、40×40，體積約束為0.3。宏微觀優(yōu)化參數(shù)同上一實例。

圖6 優(yōu)化模型(實例2)

在不同子結構大小下，該優(yōu)化模型的設計結果如圖7所示。從優(yōu)化構型來看，宏觀結構基本相似，其結構柔度值隨著子結構的增大先增大后減小。從微觀構型來看，無論子結構尺寸的大小，均可得到相似的微觀構型，隨著子結構尺寸的增大，構型越來越清晰，其周期性子結構如圖8所示。

(a)m×m=10×10, (b)m×m=20×20,c=12.97 c=13.55

(a)m×m=10×10 (b)m×m=20×20

從優(yōu)化目標函數(shù)分析，本文所提的周期性結構設計方法所得結構柔度值均小于ARSP方法[1]所得結構柔度值。相比于ARSP方法，本文方法由于擴展了子結構內(nèi)部的設計空間，可以得到柔度更小的優(yōu)化構型。

在宏微觀協(xié)同優(yōu)化過程中，由于微觀結構優(yōu)化嵌入到了宏觀結構優(yōu)化中，為了得到清晰的微觀構型，微觀優(yōu)化變量迭代步長應小于宏觀優(yōu)化步長，避免迭代步長太大使得微觀構型變化劇烈，從而引起宏觀優(yōu)化時的收斂速度緩慢的問題。圖7d的優(yōu)化收斂過程如圖9所示。在迭代計算中，最后一步的柔度值有較大程度的減小，這是由于在該迭代步的前一步已經(jīng)得到了最優(yōu)構型，即該迭代步不屬于優(yōu)化過程的迭代步，而是為了獲得優(yōu)化構型的實際柔度值，通過有限元方法進行結構的直接計算，即懲罰因子設置為1時獲得的結構柔度值。通過SIMP方法，在整個優(yōu)化過程中宏觀材料含量不變，其材料在微觀周期性單胞構型和宏觀周期性結構分布之間進行匹配，以同時獲得宏微觀構型。

圖9 圖7d的收斂過程

對比上述兩個實例，為了避免由于宏觀子結構數(shù)增多而導致的周期性子結構大小減小，進而導致材料在周期性子結構中過度集中，難以獲得具有明顯幾何特征的周期性子結構，在利用所提方法進行周期性宏微觀結構設計時，盡量采用固定周期性子結構單胞大小的形式進行優(yōu)化計算，既可獲得清晰的宏觀構型，又可避免由于材料集中引起的周期子結構單胞結構特征不清晰的問題。

4 結論

本文提出了基于子結構法的周期性微觀構型的宏微觀協(xié)同設計方法，根據(jù)子結構凝聚和反求構建對應的宏微觀優(yōu)化模型?；赟IMP方法的優(yōu)化框架，利用子結構的自由度凝聚實現(xiàn)結構的宏觀構型設計，利用子結構內(nèi)部自由度反求來實現(xiàn)周期性微觀構型設計，同時給出了該優(yōu)化方法的靈敏度計算和變量更新迭代方法。通過兩個實例驗證了所提方法的合理性。

相比于傳統(tǒng)基于均勻化的宏微觀協(xié)同優(yōu)化方法，基于子結構的宏微觀協(xié)同優(yōu)化方法可利用子結構及其對應超單元實現(xiàn)宏觀結構的快速迭代和微觀結構的再設計。由于子結構是宏觀設計域的一部分，網(wǎng)格、節(jié)點等是一一對應的，故解決了均勻化方法的尺度分離問題。同時通過子結構自由度凝聚提高了宏觀結構的設計效率。

盡管本文基于子結構實現(xiàn)了周期性結構的宏微觀協(xié)同設計，但仍然存在以下幾個問題，需要進一步研究。一是優(yōu)化效率需要進一步研究。通過構建子優(yōu)化問題實現(xiàn)微觀周期性結構的設計，子結構劃分的網(wǎng)格數(shù)越多，則微結構構型越精細，但也制約著整個優(yōu)化問題的計算效率。二是宏微觀設計參數(shù)的匹配仍需要進一步研究。較大的迭代步長會加快結構設計的收斂速度，本文中宏觀結構與微結構通過兩個關聯(lián)的優(yōu)化問題分別設計，較小的子結構設計變量迭代步長將影響宏觀結構的收斂，較大的宏觀設計變量迭代步長可能會讓子結構得不到合適的微觀構型。