連續體結構的變密度拓撲優化方法研究

王景良，朱天成，朱龍彪，許飛云

（1.江蘇海事職業技術學院輪機電氣與智能工程學院，江蘇 南京 211100；2.中國聯合網絡通信集團有限公司江蘇省分公司，江蘇 南京 211100；3.南通大學機械工程學院，江蘇 南通 226019；4.東南大學機械工程學院，江蘇 南京 211100）

拓撲優化是實現結構輕量化的一種有效的設計方法，蘊含如現代力學、現代數學、計算機科學等多學科知識。通過拓撲優化可產生材料分布合理、性能各異、質量最小的結構設計方案，故其已成為結構優化設計的熱點研究方向[1]。于20世紀90年代提出的備受關注的變密度法作為連續體結構拓撲優化設計的一種主流方法，因自身所具備的技術優勢而在航空航天、建筑、機械等領域的結構優化設計中得到廣泛應用[2-5]。如：為了實現高速列車底架的輕量化設計，陳秉智等[6]采用基于OptiStruct軟件的變密度法對高速列車的底架進行拓撲優化設計，獲得了可同時滿足強度、剛度及模態頻率等設計要求且應力分布得到改善的優質輕質底架；為了得到輕質和性能完善的汽車用擋位互換機構，王銘昭等[7]也采用基于OptiStruct軟件的變密度法對擋位互換機構的殼體進行了拓撲優化設計；周圍等[8]借助基于ANSYS Workbench軟件的變密度法研究了巖心鉆機送桿機構運送支架的拓撲優化設計問題，得到了新型輕質運送支架。

盡管變密度法有諸多優點，但所得優化結構普遍存在如棋盤格、網格依賴及灰度單元等數值不穩定問題[9]。為了解決此類問題，學者們提出了如密度過濾法、敏度過濾法、灰度單元抑制算子等多種方法[10-13]。Biyikli等[14]采用一種簡單有效的改進密度過濾法解決了采用變密度法引起的結構棋盤格問題。張國鋒等[15]為了消除采用變密度法引起的棋盤格、網格依賴和邊界擴散等問題，在保留原卷積因子的基礎上，通過引入新卷積因子，提出了一種改進敏度過濾方法。高翔等[16]為了解決采用變密度法引起的灰度單元和現有Heaviside函數須頻繁調用密度矩陣的問題，結合優化準則法和Heaviside函數，設計了基于改進Heaviside函數的拓撲優化方法。廉睿超等[17]從變密度法的材料插值方案出發，提出了一種灰度單元分層雙重懲罰方法，在有效抑制灰度單元產生的同時提高了拓撲優化的收斂速度。但是，目前所提出的策略或方法尚不能很好地兼顧解決多類型的數值不穩定問題與改進拓撲優化方法。

為了解決體積約束和柔順度最小的結構拓撲優化問題及采用變密度法引起的數值不穩定問題，以兼顧改進拓撲優化方法為前提，引入改進的固體各向同性材料懲罰（solid isotropic material with penalization，SIMP）法和基于高斯權重函數的敏度過濾法，并設計新的灰度單元抑制算子，對由變密度法和優化準則法形成的經典拓撲優化方法進行改進，提出了一種新的拓撲優化方法。

1 連續體結構拓撲優化

采用新的拓撲優方法研究連續體結構的優化設計問題，需要確定材料插值方案、建立結構拓撲優化模型、設計靈敏度、確定敏度過濾法、確定更新準則及給出評價指標等幾個步驟。

1.1 基于改進SIMP法的材料插值方案

采用能在0與1之間連續取值的基于SIMP法的材料插值方案可以有效降低結構優化設計的難度。參考Andreassen等[18]的研究成果，目前應用較為普遍的材料插值方案為：

式中：E0為實體材料的彈性模量；Emin為無效材料的彈性模量；p為無效材料的懲罰系數；xi為拓撲設計變量，i=1，2，…，n；Ei為當量彈性模量。

通過進一步研究發現，式（1）所示的插值方案對改進拓撲優化方法的效果有限。為了進一步改進拓撲優化方法，采用Zhu等[19]提出的改進SIMP法構建材料插值方案，具體如下：

式中：α為調整系數。

須說明的是，文獻[19]提出的材料插值方案中并沒有Emin，此處加入Emin的目的是規避奇異剛度矩陣的產生和確保xi的下限值可取0。

1.2 結構拓撲優化模型

建立優化模型是開展結構拓撲優化設計的前提。在式（2）的基礎上，參考文獻[9，18]，則使結構柔順度最小的拓撲優化模型為：

式中：x為結構拓撲設計變量的矢量；c為結構柔順度；K為結構的剛度矩陣；U為載荷矢量F作用在結構上而產生的位移矢量；V0和V分別為結構設計域的材料體積和結構優化設計后的材料體積；Vf為優化結構的許用材料體積分數。

1.3 靈敏度

采用經典優化準則法求解式（3）所示的優化模型時，須設計目標函數和約束關于拓撲設計變量xi的靈敏度。根據文獻[18]，基于式（2）和式（3），求解c和V關于xi的靈敏度：

式中：ui為單元的位移矩陣；k0為具有單位彈性模量的剛度矩陣。

1.4 基于高斯權重函數的敏度過濾法

為了有效解決采用變密度法引起的數值不穩定問題及改進拓撲優化方法，選用合適的敏度過濾法極為重要。Cheng等[20]發現，采用基于高斯權重函數的敏度過濾法可以得到輕質且拓撲構型好的優化結構。本文采用基于高斯權重函數的敏度過濾法對c關于xi的靈敏度進行過濾處理，具體為：

式中：γ為為了規避當xi取為0時式中分母為0而引入的一個較小正值參數；Ni為與單元i相鄰的單元j所組成的集合；Hi，j為敏度過濾的權重系數。

根據文獻[10，20]，得到：

式中：di，j為單元i、j之間的距離；σ為高斯參數，與過濾半徑R有關。

1.5 基于新灰度單元抑制算子的更新準則

基于優化準則法求解式（3）時，需要一種啟發式更新準則。根據文獻[18]，常用的更新準則為：

式中：m和η分別為正值移動步長和阻尼系數；xi與的結合物可以視為簡單的灰度單元抑制算子，用于更新xi。

中間變量Bi為：

式中：λ為拉格朗日乘子，可通過二分法計算得到。

為了有效抑制灰度單元的產生，受Groenwold等[21]提出的灰度單元抑制算子的啟發，在原灰度單元抑制算子的基礎上，引入正值比例參數β，設計了一種新的灰度單元抑制算子，為：

式中：q為灰度抑制參數。

為了提高拓撲優化方法的收斂速度，在綜合考慮每次迭代計算所得的結構柔順度的基礎上，設計了新的參數q的更新方案，具體為：

式中：q0和Δq分別為q的初值和增量；ξ和c*分別為強度因子和結構柔順度的臨界值。

須說明的是，通過測試發現，c*值設置合理與否會對優化結構的柔順度、拓撲構型及拓撲優化方法的收斂速度產生較大影響。

根據式（7）和式（9），基于新灰度單元抑制算子的更新準則為：

1.6 評價指標

選擇結構的柔順度c、程序循環執行次數I和灰度單元比例指數G作為評價拓撲優化方法的指標。其中：c用于反比表征獲取柔順度小的優化結構的能力；I用于反比表征收斂速度或收斂性能；G用于反比表征抑制灰度單元產生的能力。

參考文獻[9]，可得：

式中：n*為滿足約束條件（xi＜ 0.01或xi＞ 0.99）的單元數量。

2 算例分析與討論

采用2個經典的數值算例，分別測試改進SIMP法、基于高斯權重函數的敏度過濾法、新灰度單元抑制算子和新拓撲優化方法對其的優化效果，并與采用文獻[14，8]提出的拓撲優化方法的優化結果進行對比。設置與文獻[14，18]相同的基本參數：E0=1 MPa，Emin=1×10-9MPa，材料泊松比υ=0.3，p=3，m=0.2，η=1×10-3。采用雙線性矩形單元對結構進行網格劃分。另外，須特別說明的是，本文設計的拓撲優化方法的程序是在文獻[18]所提出方法代碼的基礎上編寫的，本文方法與文獻[18]方法的區別在于材料插值方案、敏度過濾法和灰度單元抑制算子這幾點上。

2.1 算例1

MBB（Messerschmitt-Bolkow-Blohm）梁結構如圖1所示。其幾何尺寸為280 mm×50 mm，頂端中部承受大小為2 N的外載荷F。對MBB梁結構進行拓撲優化，來測試采用改進SIMP法、基于高斯權重函數的敏度過濾法、新灰度單元抑制算子和新拓撲優化方法的優化效果。根據對稱性，采用具有1/2設計域的MBB梁為優化對象。對梁結構進行網格劃分，可以得到單元總數量為140×50個的有限元模型。設置R=3.1，Vf=0.4。

圖1 MBB梁結構Fig.1 MBB beam structure

不同SIMP法下MBB梁結構的拓撲構型和優化結果分別如圖2和表1所示；在本文SIMP法下，采用不同敏度過濾法得到的結構拓撲構型和優化結果分別如圖3和表2所示；在本文SIMP法和敏度過濾法下，采用不同灰度單元抑制算子得到的結構拓撲構型和優化結果分別如圖4和表3所示；在不同拓撲優化方法下結構拓撲構型和優化結果分別如圖5和表4所示。

圖2 不同SIMP法下MBB梁結構的拓撲構型Fig.2 Topological configuration of MBB beam structure under different SIMP methods

表1 不同SIMP法下MBB梁結構優化結果Table 1 Optimization results of MBB beam structure under different SIMP methods

圖3 不同敏度過濾法下MBB梁結構的拓撲構型Fig.3 Topological configuration of MBB beam structure under different sensitivity filtration methods

表2 不同敏度過濾法下MBB梁結構優化結果Table 2 Optimization results of MBB beam structure under different sensitivity filtration methods

圖4 不同灰度單元抑制算子下MBB梁結構的拓撲構型Fig.4 Topological configuration of MBB beam structure under different gray-scale unit suppression operators

表3 不同灰度單元抑制算子下MBB梁結構優化結果Table 3 Optimization results of MBB beam structure under different gray-scale suppression operators

圖5 不同拓撲優化方法下MBB梁結構的拓撲構型Fig.5 Topology configuration of MBB beam structure under different topology optimization methods

表4 不同拓撲優化方法下MBB梁結構優化結果Table 4 Optimization results of MBB beam structure under different topology optimization methods

由圖2可知，采用不同SIMP法得到的MBB梁結構的拓撲構型存在一定差別：采用本文SIMP法所得優化結構在“區域A”和“區域B”處的顏色比較淺。意味著采用文獻[18]的SIMP法所得結構在這些區域的xi值比較大。

由表1可知，相比于文獻[18]的SIMP法，采用本文SIMP法所得結構的c和I分別減少了2.09% 和22.89% ，G增加了13.80% 。說明采用本文SIMP法可以使拓撲優化方法的收斂速度變得更快，而且可以增強獲取柔順度小的優化結構的能力，但也會削弱抑制灰度單元產生的能力。

由圖3可知，圖3（a）中“區域A”和“區域B”處的形狀不同于圖3（b），且在“區域A”處的區別更為明顯，說明使用不同的敏度過濾法會得到不同的結構拓撲構型。

由表2可知，相比于文獻[18]的敏度過濾法，采用本文敏度過濾法所得的c和G分別減少了1.30% 和10.25% ，I則相同。說明采用本文的敏度過濾法可以增強獲取柔順度小的優化結構和抑制灰度單元產生的能力。

由圖4可知，采用本文灰度單元抑制算子所得優化結構的拓撲構型更為清晰。

由表3可知，相比于文獻[18]的灰度單元抑制算子，采用本文灰度單元抑制算子所得的c、G和I分別減少了4.46% 、85.71% 和45.31% 。說明采用本文灰度單元抑制算子不僅可以增強獲取柔順度小的優化結構和抑制灰度單元產生的能力，而且可以提高收斂速度。

由圖5可知，采用文獻[14，18]方法所得結構的拓撲構型存在更多的灰度單元，而采用本文方法所得結構的輪廓非常清晰，且從宏觀上看不到灰度單元，所得結構比較好。

由表4可知，相比于文獻[14，18]，采用本文拓撲優化方法所得的c分別減少了7.14% 和7.68% ，G分別減少了83.97% 和85.41% ，I分別減少了51.39% 和57.83% 。說明本文方法在收斂速度、獲取柔順度小和拓撲構型好的優化結構及抑制灰度單元產生等方面具有明顯優勢。

2.2 算例2

懸臂梁結構如圖6所示。其幾何尺寸為160mm×70 mm，梁的左端被固定，右端中部承受大小為1 N的豎直向下的外載荷F。網格劃分后結構的單元總數量為160×70個。設置R=3.5，Vf=0.4。采用3種拓撲優化方法所得懸臂梁結構的拓撲構型和優化結果分別如圖7和表5所示。

圖6 懸臂梁結構Fig.6 Cantilever beam structure

圖7 不同拓撲優化方法下懸臂梁結構的拓撲構型Fig.7 Topological configuration of cantilever beam structure under different topology optimization methods

表5 不同拓撲優化方法下懸臂梁結構優化結果Table 5 Optimization results of cantilever beam structure under different topology optimization methods

由圖7可知，采用不同拓撲優化方法所得結構的拓撲構型存在較明顯差異。相比于文獻[14，18]，采用本文方法所得結構的輪廓更為清晰。由表5可知，采用本文方法所得的c、I和G分別為104.98 N·mm、38次和4.73% ，較文獻[14]方法分別減少了8.63% 、67.80% 和80.90% ，較文獻[18]方法分別減少了8.90% 、77.25% 和82.32% 。說明本文方法在獲取柔順度小且拓撲構型好的優化結構和抑制灰度單元產生等方面均明顯優于文獻[14，18]的方法，此外，收斂速度也明顯加快。

3 結 論

為了實現使連續體結構的體積約束和柔順度最小的拓撲優化及解決采用經典變密度法引起的數值不穩定問題，在兼顧改進拓撲優化方法的前提下，通過引入改進的SIMP法、基于高斯權重函數的敏度過濾法和新灰度單元抑制算子，并結合優化準則法，提出了新的拓撲優化方法。通過優化算例可知，新的拓撲優化方法能夠有效解決連續體結構的拓撲優化問題，具有收斂速度較快、可更好地獲取柔順度小且拓撲構型好的優化結構和抑制灰度單元產生等優勢。