基于雙曲正弦函數的拓撲優化在熱傳導結構設計中的應用

2020-06-27 04:35:46薛紅軍陶才勇黨思娜

航空工程進展 2020年3期

薛紅軍，陶才勇，黨思娜

(西北工業大學航空學院，西安 710072)

0 引言

飛機電子設備不斷增多，其散熱結構是航空工程結構設計過程中必須考慮的重要問題之一[1]，飛機中的傳熱結構設計既要滿足一定的強度和剛度要求以達到特定的功能，又要滿足散熱要求，目的是保證關鍵設備不會因為溫度過高而出現故障。通過熱力學計算設計的傳統散熱裝置，已無法滿足對熱防護結構具有高性能要求的實際工程需求[2]，因此，尋求一種新的散熱裝置設計方法成為解決熱傳導問題的關鍵。

針對熱傳導問題，H.Rodrigues等[3-4]和C.Jog[5]將拓撲優化應用在熱彈性結構中，初次實現了拓撲優化在熱傳導中的應用；A.Bejan[6]在熱傳導的背景下介紹了熱傳導通道構造理論，建立了所謂的“體-點(VP)問題”或“通道問題”，并討論了在發熱域(例如CPU)中布置固定數量的材料的需要；Li Q等[7]明確地將ESO用于熱傳導問題。

現階段應用于熱傳導結構的拓撲優化方法主要有：漸進結構優化方法(ESO)、雙向漸進結構優化方法(BESO)[8]、變密度方法[9]、水平集方法[10]、拓撲導數方法[11]。就目前的研究成果而言，應用比較廣泛的拓撲優化方法是變密度法。

基于變密度方法的散熱結構的拓撲優化設計主要使用：帶懲罰的各向同性固體材料模型[12](Solid Isotropic Material with Penalization，簡稱SIMP)和材料屬性的理性近似模型[13](Rational Approximation of Material Properties，簡稱RAMP)。F.H.Burger等[14]利用SIMP和移動漸近線(MMAs)的方法，探討了體-點問題的3D解決方案；T.Van Oevelen等[15]采用SIMP模型演示了用雙層簡化模型表示全三維解的共軛傳熱問題的解；G.Marck等[16]使用SIMP方法進行了多目標優化；Li Jiachun等[17]使用了材料特性的理性近似模型(RAMP)，進行了熱傳導結構優化設計的探索。兩種模型對于導熱結構的拓撲優化設計都能得到比較好的布局，但是對于“體-點”問題[18]，SIMP方法會產生很多樹狀末梢結構，這樣的結果有一定的缺點，如可制造性差，并且結果并不是最靠近最優值的結果[19]。因此有學者通過RAMP方法得到了一定的改進，能過獲取相對更加良好的拓撲結果[20]。但是，這兩種方法有其局限性，在解決二維、三維均勻熱源的散熱問題時，SIMP方法不夠精確，RAMP方法迭代步長過多。

本文在變密度法的基礎上，建立基于sinh函數的插值模型，與SIMP和RAMP模型的計算結果進行對比，并將此模型應用于機載LRM模塊導熱拓撲優化設計中進行驗證。

1 基于sinh函數插值模型的拓撲優化在熱傳導結構設計中的模型建立

基于變密度法的拓撲優化的關鍵是選擇合適的插值函數和懲罰技術，將問題的物理量表示為連續設計變量的函數。本文采用基于sinh函數的插值模型[19]進行求解計算。

1.1 基于sinh函數的插值模型

基于sinh函數的插值模型，其實質是在插值模型中應用了雙曲正弦函數(sinh函數)：

(1)

可以簡化為

(2)

sinh函數插值模型采用雙曲正弦函數來對中間密度進行懲罰，其函數關系如圖1所示。

圖1 sinh插值模型

sinh函數插值模型的微分形式為

(3)

當ρ=0時，

當ρ=1時，

1.2 sinh插值模型靈敏度分析

計算導數最好的辦法是使用伴隨方法，在此也采用伴隨方法推導sinh函數插值模型的目標函數的靈敏度。

本文以熱柔度為目標函數，基于sinh函數插值模型的目標函數為

(4)

式中：K0為單元導熱系數矩陣。

對于上面給出的熱柔度目標函數，通過添加零函數，將其重寫為

C=QTT-λT(KT-Q)

(5)

式中：λ為拉格朗日乘數，是一個任意的、但是固定的實數向量。

然后將其對設計變量求導得到：

(6)

當λ滿足伴隨方程

QT-λTK=0

(7)

目標函數靈敏度可寫為

(8)

根據伴隨方程，可以直接得到λ=T，因此目標函數靈敏度為

(9)

其中，對于傳熱系數矩陣對設計變量的導數為

(10)

因此，當傳熱系數矩陣采用sinh函數插值模型對其進行插值，則目標函數的靈敏度為

(11)

式中：K0為具有單元熱傳導矩陣，當材料的導熱系數為k0時，對于二維的四節點矩形單元，單元的尺寸為單位長度，K0可寫為

(12)

對于三維八節點正六面體單元，單元的尺寸為單位長度，單元導熱系數矩陣K0可寫為

(13)

當懲罰因子p=3時，靈敏度隨著設計變量的變化曲線如圖2所示，可以看出：所有的靈敏度值

都是負值，也就意味著隨著導熱材料的增加，目標函數在不斷地減小，導熱性能更好。

圖2 靈敏度變化曲線

為了避免棋盤格現象和網格依賴現象，在計算中采用靈敏度過濾技術，將靈敏度改寫為

(14)

式中：Hei為單元i的一個集合；γ為一個比較小的常數，為避免分母為0，可取γ=10-3；單元i中心到單元e的距離Δ(e,i)小于最小濾波半徑rmin。

1.3 基于sinh函數插值模型的優化準則

本文以結構的散熱弱度最小為目標函數，體積分數、結構平衡方程和設計變量約束為約束條件，建立基于sinh函數插值模型的熱傳導結構的拓撲優化數學模型為

min:

(15)

式中：C為散熱弱度；K為結構的總體傳熱系數矩陣；F為溫度載荷向量；T為節點溫度向量；Ve為單元體積；V0為結構初始狀態下總的體積；f為體積分數，是優化后的結構總體積與初始狀態下總體積之比；ximin為設計變量取值的下限，通常取一個較小的數，如0.001；ximax為設計變量取值的上限，取為ximax=1。

利用拉格朗日乘數法將多為高次多約束問題轉換為無約束最優化問題，對應于拓撲優化問題的拉格朗日函數為

(16)

式中：λ1為拉格朗日乘子；λ2，λ3，λ4為拉格朗日乘子向量；bi和ci為松弛因子；x為xi組成的列向量。

(17)

(18)

考慮上式等于0的情況，并且代入C=TTK(X)T有：

(19)

利用導熱系數矩陣的對稱性得出：

(20)

式中：λ2為列矢量，其取值無限制。取λ2=-2T代入式(20)得出：

(21)

進而得出：

(22)

因此可以得到基于優化準則法的迭代公式為

(23)

式中：m為移動極限，通常取一個較小的正數；η為阻尼系數，引入η的目的是為了保證數值計算的穩定性和收斂性。

(24)

(25)

將式(21)代入式(25)得到：

(26)

通過求解上述方程便可以求解出λ1，通常采用二分法進行求解。

設計變量通過優化準則更新，其收斂性的檢查方法如式(27)所示：

|max (xk+1-xk)|<ε

(27)

式中:xk+1為更新的設計變量；xk為舊的設計變量；ε為可行性誤差，即每個單元更新后的密度值與前一次的密度值之差的絕對值不超過規定的誤差值。

2 模型算法實現以及算例分析

2.1 算法實現

本文在MATLAB上進行程序編譯，實現模型功能，整體的實現步驟如圖3所示。

圖3 基于sinh的熱傳導結構拓撲優化方法算法流程

2.2 算例分析

2.2.1 二維算例

(1) 體-點問題

算例1：整塊板均勻加熱，左邊界中心處給定溫度T=0 ℃，其余邊界絕熱，其結構示意圖如圖4所示。

圖4 導熱板均勻加熱，左側中心處冷卻模型示意圖

高熱導率材料的導熱系數為k0=1 W/m·℃的假想材料，低熱導率材料的導熱系數為kmin=0.001 W/m·℃的假想材料，模型離散網格為180×180。優化參數設置為：體積分數0.4；懲罰因子p=3；采用敏度濾波，濾波半徑r=2。分別采用SIMP模型，RAMP模型和基于sinh函數模型對其進行優化計算，其優化結果如圖5所示，對比結果如表1所示。

(a) SIMP (b) RAMP (c) sinh模型

圖5 導熱板均勻加熱，左側中心處冷卻模型優化結果

Fig.5 Optimization results of heat-conducting plate evenly heating and cooling at the left center

表1 導熱板均勻加熱，左側中心處冷卻模型優化數據對比

(2) 均勻熱源問題其他邊界條件

算例2：整塊板均勻加熱，左邊界給定溫度T=0 ℃，其余邊界絕熱，其結構示意圖如圖6所示。高熱導率材料的導熱系數為k0=1 W/m·℃的假想材料，低熱導率材料的導熱系數為kmin=0.001 W/m·℃的假想材料，優化參數設置為：體積分數0.5；懲罰因子p=3；采用敏度濾波，濾波半徑r=3。分別采用SIMP模型，RAMP模型和基于sinh函數模型對其進行優化計算，其優化結果如圖7所示，對比結果如表2所示。

圖6 導熱板均勻加熱，左側邊界冷卻模型示意圖

(a) SIMP (b) RAMP (c) sinh模型

圖7 導熱板均勻加熱，左側邊界冷卻模型優化結果

Fig.7 Optimization results of heat-conducting plate evenly heating, and the left edge cooling

表2 導熱板均勻加熱，左側邊界冷卻模型優化數據對比

(3) 集中熱源問題優化分析

算例3：平板中心加熱，四個角點給定溫度T=0 ℃，其結構示意圖如圖8所示,其優化參數與算例2一致。分別采用SIMP模型，RAMP模型和基于sinh函數模型進行優化計算，優化結果如圖9所示，對比結果如表3所示。

圖8 導熱板的中心加熱，四個頂點冷卻模型示意圖

(a) SIMP (b) RAMP (c) sinh模型

圖9 導熱板的中心加熱，四個頂點冷卻模型優化結果

Fig.9 Optimization Results of heat-conducting plate center heated and four apexes cooled

表3 導熱板的中心加熱，四個頂點冷卻模型優化結果數據對比

算例4：板上有5個集中熱源，四周邊界給定溫度T=0 ℃，其結構示意圖10如圖所示,其優化參數與算例2一致。高熱導率材料的導熱系數為k0=1 W/m·℃的假想材料，低熱導率材料的導熱系數為kmin=0.001 W/m·℃的假想材料。優化參數設置為：體積分數0.5；懲罰因子p=3；采用敏度濾波，濾波半徑r=2.5。分別采用SIMP模型，RAMP模型和基于sinh函數模型進行優化計算，優化結果如圖11所示，對比結果如表4所示。

圖10 導熱板上5個集中熱源加熱，四周邊界散熱模型示意圖

(a) SIMP (b) RAMP (c) sinh模型

圖11 導熱板上5個集中熱源加熱，四周邊界散熱模型優化結果

Fig.11 Optimization results of the heating of 5 centralized heat sources on the thermal plate and the heat dissipation at the surrounding boundary

表4 導熱板上5個集中熱源加熱，四周邊界散熱模型優化結果數據對比

2.2.2 三維算例

(1) 均勻熱源問題

算例5：設計域下面均勻加熱，上平面四個頂點溫度固定T=0 ℃，模型示意圖如圖12所示。高熱導率材料的導熱系數為k0=1 W/m·℃的假想材料，低熱導率材料的導熱系數為kmin=0.001 W/m·℃的假想材料。優化參數設置為：體積分數0.3；懲罰因子p=3；采用密度濾波，濾波半徑r=1.4。分別采用SIMP模型，RAMP模型和基于sinh函數模型進行優化計算，優化結果如圖13所示，對比結果如表5所示。

圖12 下表面均勻加熱，上表面四個頂點冷卻模型示意圖

(a) SIMP (b) RAMP (c) sinh模型

圖13 下表面均勻加熱，上表面四個頂點冷卻模型優化結果

Fig.13 Optimization results of uniform heating of the lower surface and cooling of the four vertices of the upper surface

表5 下表面均勻加熱，上表面四個頂點冷卻模型優化結果數據對比

(2) 集中熱源問題

算例6:正方體中心有集中熱源，八個頂點溫度固定為T=0 ℃，模型示意圖如圖 14所示。高熱導率材料的導熱系數為k0=1 W/m·℃的假想材料，低熱導率材料的導熱系數為kmin=0.001 W/m·℃的假想材料。優化參數設置為：體積分數0.3；懲罰因子p=3；采用密度濾波，濾波半徑r=1.4。分別采用SIMP模型，RAMP模型和基于sinh函數模型進行優化計算，優化結果如圖15所示，對比結果如表6所示。

圖14 正方體中心加熱，八個頂點冷卻模型示意圖

(a) SIMP (b) RAMP (c) sinh模型

圖15 正方體中心加熱，八個頂點冷卻模型優化結果

Fig.15 Optimization results of cube center heating and eight vertices cooling

表6 正方體中心加熱，八個頂點冷卻模型優化結果數據對比

通過算例將SIMP模型和RAMP模型與基于sinh函數的插值模型進行對比，得出基于sinh函數的插值模型較SIMP插值模型精確，較RAMP插值模型的迭代步長少，計算效率高。

3 機載LRM模塊導熱拓撲優化設計

機載LRM(現場可更換模塊)模塊已成為飛機電子系統設備的硬件標準,其優點是標準程度高、模塊功能化強、維修性好、成本低。由于機載LRM模塊是由超大規模集成電路構成的，包含大量的VLSI和VHSIC芯片，而這些設備在工作時會產生大量熱量，使得LRM模塊溫度急劇上升，因此，機載LRM模塊輕量化的冷卻成為設計的首要問題。

由于LRM模塊空間的限制，傳導冷卻方法被優先考慮。此方法是通過將熱量從LRM模塊芯片傳至導熱體，再經由導熱體傳至機架冷板，最后通過通風冷卻或通液冷卻等散熱形式將熱量帶走。簡化后的LRM模塊導熱模型如圖16所示，導熱體材料為銅，導熱系數為377 W/m·℃。