王 辰, 劉義暢, 陸宇帆, 賴章龍, 周明東
(上海交通大學 機械與動力工程學院,上海 200240)
隨著增材制造技術取得顯著進展,更多性能優異但傳統難以制造的結構進入人們的視野,其中多孔填充結構由于具有出色的抗屈曲性[1]、抵抗材料缺陷和外力變化的魯棒性[2]等性能,引起了國內外研究學者的廣泛關注.填充結構具有承受載荷和支撐外部殼體的重要作用,除結構剛度與輕量化指標外,填充結構設計需要兼顧結構強度,以避免出現應力集中現象所導致的結構斷裂失效問題[3-5].
實際工程應用領域中,填充結構設計首先確定多孔填充幾何形式,如直線晶格、蜂窩[6]等結構,再基于這些結構按照規定體積分數進行均勻填充,得到具有一定結構強度的輕質多孔填充結構.但該設計方法多孔填充幾何形式與填充結構體積分數的選取均依賴工程經驗,并且內部多孔填充形式在優化過程中不能根據載荷工況進行調整,所得填充結構力學性能有限.與先選型再填充的設計方法相比,結構拓撲優化是在給定設計空間內優化結構材料布局的設計方法,具有更大的設計自由度,能夠根據載荷逐點確定材料的分布,得到性能更加優良的填充結構.Wu等[2]提出了一種考慮局部體積約束的給定殼體填充結構拓撲優化設計方法,能夠得到具有優異魯棒性的輕質多孔填充結構.Liu等[7]基于Wu等[2]的工作提出了一種基于密度的自支撐輕量化填充結構拓撲優化設計方法,能夠得到剛度優異且滿足增材制造支撐條件的輕質填充結構.Qiu等[8]提出了一種殼-填充結構設計方法,該方法首先在假定填充結構由低強度材料組成的情況下開展殼體優化設計,然后基于優化所得殼體結構,考慮最小尺寸約束實現內填充結構優化設計.Groen等[9]提出了一種基于正交各向異性填充結構的殼-填充結構拓撲優化設計方法,通過將基于均質化的設計結果投影至更高分辨率,能夠在更低的計算成本下得到剛度優化的輕質殼-填充結構.Wadbro等[10]提出了一種多尺度殼-填充結構拓撲優化方法,能夠同時優化宏觀殼體與均勻填充的拓撲結構.Wu等[11]提出了一種包含兩組設計變量的參數化方案,能夠實現外部殼體和非均勻填充的協同拓撲優化設計,優化結果與相同材料體積下均勻填充結構相比具有更高的結構剛度.Zhou等[12]基于Wu等[11]的工作提出一種面向增材制造的自支撐殼-填充結構協同拓撲優化方法,該方法采用一種全新的增材制造自支撐約束,并通過優化結果后處理,能夠得到剛度優異且滿足增材制造工藝要求的殼-填充結構.上述填充結構拓撲優化方法能夠設計具有優異比剛度的填充結構,但在優化過程中均未考慮結構強度,優化后的結構易出現應力集中從而導致結構斷裂失效.
自Duysinx等[13]提出考慮應力的連續體結構拓撲優化研究以來,應力相關拓撲優化問題備受關注,考慮應力的結構拓撲優化主要存在兩項關鍵難點.首先是應力分布奇異性問題[14],當拓撲優化設計變量趨于某一臨界值時,部分密度較低的單元仍具有較高的應力值,導致優化問題無法穩定收斂.目前已有多種應力約束松弛方法能夠解決應力分布奇異性問題,如ε-松弛方法[13, 15]、二次規劃(QP)松弛方法[16]和應力懲罰方法[14]等.另一項關鍵難點是應力高度非線性特性問題[14],應力的大小和分布受填充結構密度變化影響較為明顯,這種現象在應力梯度較大的區域如高曲率結構邊界處更為嚴重,采用適當的優化公式和求解算法[17]能夠有效提高拓撲優化應力響應的收斂穩定性[4, 18].
目前針對實體結構的應力相關優化方法已較為成熟,但針對多孔填充結構的應力相關優化方法研究相對較少.Yu等[19]提出了一種考慮應力約束的殼-填充結構拓撲優化設計方法,通過對殼體與填充結構分別施加不同的強度準則,實現了最大應力可控的殼體與均勻多孔填充結構協同拓撲優化.然而該方法在給定多孔填充形式的基礎上優化填充密度,優化過程中填充結構形式無法根據載荷條件進行調整,力學性能存在提升空間.與給定多孔填充形式的優化結果相比,采用非均勻多孔填充的優化結果具有更加優異的力學性能[11],但目前尚未開展考慮非均勻多孔填充結構強度的拓撲優化方法研究.
本文針對非均勻多孔填充形式,提出了一種考慮增材制造填充結構強度的拓撲優化設計方法.該方法采用p范數應力近似表示結構最大應力值,并通過最小化該值提升填充結構強度;使用局部體積約束獲得輕質多孔填充構型,通過提出局部體積約束變上限優化策略,提高了填充結構優化過程中應力響應的收斂穩定性;采用自支撐約束保證填充結構滿足自支撐條件并能夠為給定薄壁外殼提供良好支撐;引入兩場模型控制優化結果的最小尺寸,避免出現無法增材制造的單節點連接結構.本文所提方法能夠得到滿足增材制造工藝要求且有效控制最大應力的內填充結構,所得優化結果與最小化柔度的填充結構優化結果相比,在相同體積分數情況下能夠顯著提升結構強度.
本文所采用的參數化模型如圖1所示.設計空間共分為3部分:薄壁外殼外部區域Ω0、給定薄壁外殼Ω1以及薄壁外殼內部區域Ω.其中Ω0為非設計域,在優化過程中始終為空材料;Ω1同樣為非設計域,在優化過程中始終為實體結構;Ω為填充結構設計域.采用平面單元進行優化,但優化后的結構代表平面單元沿垂直于平面方向拉伸所得結構,其邊緣部分由平面內的結構邊界拉伸所組成.

圖1 給定薄壁外殼填充結構參數化模型Fig.1 Parametrization of infill structures with a given shell

(1)

(2)

(3)
(4)
式中:ps為物理變量插值的冪指數,取ps=0.5.考慮到最大值函數不可微,采用所有單元插值應力σi的p范數σp近似表示全局最大應力:
(5)
式中:N為單元總數;p為應力范數的冪指數,用于控制p范數與實際最大值之間的近似程度,當p→∞時,σp趨近于實際的最大應力值,但p值過高會導致優化算法求解困難,而p值過小將導致p范數無法準確近似結構的最大應力,根據文獻[14]取p=8.
采用一種基于增材過濾器(AM filter)的自支撐約束[22],保證內填充結構滿足自支撐條件并能夠從結構內部支撐給定薄壁外殼.給定薄壁外殼的外部是否滿足支撐條件并非本文研究內容,因此假設給定薄壁外殼的外部已由支撐結構提供良好支撐.
增材過濾器的基本原理如圖2所示[12,23],定義單元(i,j)下方相鄰3個坐標為(i-1,j-1)、(i-1,j)、(i-1,j+1)的單元為該單元的支撐區域.如果單元(i,j)支撐區域的3個單元中至少一個單元為實體結構,單元(i,j)即滿足自支撐條件.如果單元(i,j)不滿足支撐條件,那么增材過濾器將移除該單元材料,通過自下而上逐層使用增材過濾器,最終能夠實現滿足自支撐條件的結構設計.增材過濾器的公式如下:
θij=
(6)
smin(a,b)=
(7)
(8)
式中:μij為增材過濾之前的單元;θij為增材過濾之后的單元;根據文獻[22],增材過濾器參數取為ε=10-4,P=40,Q=P+lg 3/lg 0.5.

(9)
變量場μ經過增材過濾后刪除所有不滿足自支撐條件的單元,得到變量場θ.變量場μ與θ的差值μ-θ表示不滿足自支撐條件的單元,通過限制這部分單元的體積,即可實現自支撐約束:
max{(μT-θT),0T}2IVe≤εr
(10)
式中:I為單位向量;Ve為單元體積;εr=2×10-5ITIVe[12].

圖3 自支撐約束示意圖Fig.3 Illustration of the self-supporting constraint
多孔填充結構與實體填充相比,在抵抗材料缺陷和外力變化時具有優異的魯棒性,本文通過使用局部體積約束實現多孔形式的填充結構設計[2],單元i的局部體積分數定義為
(11)

(12)
式中:α為局部體積約束的上限值.

(13)
式中:P1為局部體積分數p范數的冪指數,取P1=8.
當優化過程中投影銳度參數β增加時,單元密度發生變化使得局部體積超出約束上限α,為滿足局部體積約束,優化算法將會刪除部分結構以降低局部體積分數,這將導致應力分布發生突變進而影響優化過程收斂穩定性.針對該問題,本文提出了一種優化過程中變局部體積約束上限α的優化策略,優化初始階段α取為較小值,在投影銳度參數β增加至8、16、32時,局部體積約束上限α分別增加0.02,最終達到初始給定的局部體積約束上限,通過在優化過程中局部體積超出上限時適當提高約束值,實現優化過程的穩定收斂.

(14)

基于兩場模型提出了一種面向增材制造的填充結構最小化應力拓撲優化模型:
(15)


(16)
物理場與腐蝕場均采用帶有懲罰的固體各向同性材料(SIMP)法對彈性模量進行插值:

(17)
(18)

拓撲優化流程如圖4所示,具體的優化求解過程如下.

圖4 拓撲優化流程Fig.4 Workflow of topology optimization
步驟1設置設計域Ω、給定薄壁外殼非設計域Ω1以及空白的非設計域Ω0,初始化設計變量ρ與投影銳度參數β.

步驟3求解有限元分析方程并計算設計響應.
步驟4計算各設計響應的靈敏度.
步驟5通過移動漸近線(MMA)算法[24]更新設計變量.
步驟6判斷投影銳度參數β的連續變化條件:① 投影銳度參數β<βmax;② 目標函數連續5次迭代變化量小于0.2%;③ 當前β已至少進行50次迭代;④ 所有約束函數均滿足.如果4個條件全部滿足,則投影銳度參數β增加并轉至步驟2,否則轉至步驟7.
步驟7檢查收斂性.如果滿足4個條件(① 目標函數連續5次迭代變化量小于0.2%;② 投影銳度參數β=βmax;③ 當前β已至少進行50次迭代;④ 所有約束函數均滿足),則轉至步驟8,否則重復步驟2~6.
步驟8優化結束.
通過二維平面的數值算例討論所提填充結構拓撲優化方法有效性,所有算例中的實體與空材料彈性模量分別取為 210 GPa和1×10-9MPa,泊松比均取為0.3.設計域為圖5所示的L型梁填充區域.該結構尺寸為長L=300 mm,高H=300 mm,設計域長度Ld=120 mm,設計域高度Hd=120 mm,L型梁薄壁外殼左上固定,對于右端呈正方形分布的相鄰9個節點,每點均施加豎直向下、大小F=1 N的載荷.投影銳度參數β初始取1,如果所有約束均滿足且在當前β條件下已至少迭代50步,則投影銳度參數β按如下規則遞增:當β<32時,β由1連續翻倍直至達到32;當β≥32時,β逐次加8直至達到最大值64.

圖5 L型梁設計域與邊界條件示意圖Fig.5 Illustration of the designable domain and the boundary condition for the L-bracket
取整體體積約束上限為0.55,局部體積約束上限為0.7,r=3,R=10,通過求解最小化應力的優化模型(式(15)),最小化應力優化結果如圖6所示,該優化結果的最大應力為1.899 MPa,柔度為0.091 mJ,整體體積分數為0.55.由應力云圖可知,在L型梁易出現應力集中的轉角處,應力分布較為平均,顯著緩解了優化結果的應力集中現象.

圖6 最小化應力優化結果和應力云圖Fig.6 Optimization result for minimizing stress and stress contour
優化過程中不同迭代步數對應的填充結構如圖7所示,最小化應力優化模型的設計響應收斂曲線如圖8所示.如圖8(a)所示,目標函數即應力響應在投影銳度參數變化時會發生突變并逐步恢復穩定,當達到第180迭代步之后應力開始平穩下降,收斂穩定且并未出現明顯波動.圖8(d)所示為物理變量場與腐蝕變量場不滿足自支撐條件的單元體積之和,這一數量在投影銳度參數β增加時發生突變,但是能夠逐漸下降直至低于約束值,最終優化結果的所有約束函數值均滿足約束條件.

圖7 優化過程中的填充結構Fig.7 Infill structures during the iterative optimization process

圖8 不同設計響應的收斂曲線Fig.8 Convergence curves for different design responses
為對比最小化應力與最小化柔度的優化結果,建立最小化柔度問題優化模型,其約束函數與最小化應力優化模型(式(15))完全相同,僅將式(15)所示優化模型的目標函數改為腐蝕變量場的柔度,即可得到給定薄壁外殼填充結構的最小化柔度拓撲優化模型:
(19)


圖9 最小化柔度優化結果和應力云圖Fig.9 Optimization result for minimizing compliance and stress contour

圖10 不同局部體積約束策略的目標函數和局部體積分數收斂曲線Fig.10 Convergence curves of the objective function and the local volume fraction for different local volume constraint strategies
此外,為對比局部體積約束固定上限與變上限策略對優化過程收斂穩定性的影響(見圖10),固定局部體積上限為0.7,其余優化參數與3.1節相同,優化結果如圖10(b)右側所示.由目標函數收斂曲線可知,從第377次迭代開始,固定局部體積約束上限的目標函數即應力響應突然增大,這是由于投影銳度參數β在第376步由32翻倍至64,結構發生明顯變化導致局部體積分數超出當前約束值,為滿足局部體積約束,部分結構消失從而導致整體應力響應發生振蕩,最終優化結果中出現斷開的結構,優化結果的最大應力為2.470 MPa.考慮變局部體積上限的優化結果如圖10(b)左側所示,優化參數均與3.1節相同,由目標函數收斂曲線可知,變局部體積上限的優化結果目標函數收斂穩定,應力響應并未隨β翻倍而出現明顯波動,與固定上限的優化結果相比,最大應力由2.470 MPa降至 2.010 MPa.通過圖10所示的優化結果對比可以看出,采用局部體積約束變上限的優化策略,可以避免最終優化結果出現斷開的結構,這是由于當投影銳度參數β增加時,適當提高局部體積約束值能夠避免β增加后的局部體積超出約束上限,進而避免優化算法為滿足局部體積約束而刪除部分結構,有效地提高應力優化問題的收斂穩定性.
根據3.1節的最小化應力優化結果可知,所提方法能夠有效地降低優化結果的最大應力,降低應力集中的影響,但上述優化結果僅考慮結構強度而并未考慮結構剛度,對此在式(15)最小化應力優化模型基礎上添加柔度約束:
cero≤c*
(20)
研究了不同柔度約束情況下優化結果的應力變化規律.為控制優化結果的最小尺寸,柔度約束施加于腐蝕變量場.
考慮柔度約束的最小化應力優化模型參數與3.1節相同,圖11所示為柔度約束c*分別為0.090、0.088、0.086 mJ的填充結構優化結果.根據3.1節的優化結果,不考慮柔度約束情況下最大應力值為1.899 MPa,柔度值為0.091 mJ.圖11所示結果表明,隨著柔度約束逐漸收緊,優化結果的最大應力值逐漸增大.此外通過圖11所示各結構柔度值對比可得,各優化結果柔度值均嚴格滿足柔度約束,表明該方法能夠有效實現給定柔度限制下的最小化應力設計,得到綜合考慮結構強度與結構剛度的輕量化內填充結構.

圖11 不同柔度約束c*的最小化應力優化結果Fig.11 Optimization results for minimizing stress with different compliance constraint c*
提出了一種考慮增材制造填充結構強度的拓撲優化設計方法,能夠設計滿足增材制造工藝要求的輕質高強填充結構.該方法提出局部體積約束變上限優化策略,避免優化過程中由于局部體積約束過緊出現應力響應振蕩問題,顯著改善了優化過程收斂穩定性.最小化應力優化結果表明,該方法所得優化結果與最小化柔度優化結果相比顯著降低了結構最大應力,有效地減輕了應力集中效應.在此基礎上,通過在最小化應力優化模型中添加柔度約束,能夠實現兼具剛度與強度的填充結構設計.在工程應用中,本文方法能夠針對考慮強度的給定薄壁外殼填充結構設計提供一定指導.但薄壁外殼分布可變的填充結構同樣具有優化潛力,其力學性能存在提升空間,考慮結構強度的殼體與填充結構協同拓撲優化方法有待開展進一步研究.