段治健 馬欣榮

（1.咸陽師范學(xué)院數(shù)學(xué)與統(tǒng)計(jì)學(xué)院 咸陽 712000）（2.西北工業(yè)大學(xué)機(jī)電學(xué)院 西安 710072）

1 引言

流體流動(dòng)系統(tǒng)和傳熱裝置廣泛應(yīng)用于電子設(shè)備、熱交換器網(wǎng)絡(luò)、汽車工業(yè)、發(fā)動(dòng)機(jī)、發(fā)電機(jī)、航天飛機(jī)、潛艇、火箭和衛(wèi)星設(shè)計(jì)中。尤其是發(fā)動(dòng)機(jī)、電子芯片散熱設(shè)計(jì)中冷卻流道設(shè)計(jì)不當(dāng)會導(dǎo)致零件嚴(yán)重過熱，從而限制使用性能和耐久性，甚至導(dǎo)致結(jié)構(gòu)損壞。因此，為了獲得更好的性能和更低的成本，必須設(shè)計(jì)高效的冷卻系統(tǒng)。

拓?fù)鋬?yōu)化設(shè)計(jì)為此提供了新穎思路，該方法最早由Bens?e 和Kikuchi［1］建立并成功應(yīng)用于固體力學(xué)問題的結(jié)構(gòu)優(yōu)化。優(yōu)點(diǎn)在于：允許在優(yōu)化過程中改變結(jié)構(gòu)拓?fù)?；最終設(shè)計(jì)幾乎不依賴于初始設(shè)計(jì)。Gersborg-Hansen 等［2］首次將有限體積法應(yīng)用于二維穩(wěn)態(tài)熱傳導(dǎo)問題。Borvall 和Peterson［3］的工作開創(chuàng)了流體傳熱系統(tǒng)的拓?fù)鋬?yōu)化設(shè)計(jì)研究。Yoon［4］首次在二維結(jié)構(gòu)網(wǎng)格上，忽略流體區(qū)域內(nèi)浮力和粘性耗散的影響，采用有限元方法研究了強(qiáng)迫對流溫度場的共軛拓?fù)鋬?yōu)化問題，并發(fā)現(xiàn)對流-熱傳導(dǎo)之間的平衡對傳熱系統(tǒng)中散熱結(jié)構(gòu)設(shè)計(jì)的影響較大。Sun等［5］提出了基于密度變量的三維熱流冷卻系統(tǒng)的大規(guī)模拓?fù)鋬?yōu)化方法。Joo 等［6］為了研究自然對流中的形狀依賴效應(yīng)，考慮傳熱系數(shù)在計(jì)算域內(nèi)變化，提出了一種新的替代模型，采用弱耦合拓?fù)鋬?yōu)化方法對自然對流散熱器進(jìn)行了傳熱優(yōu)化。Sato 等［7］以換熱最優(yōu)和壓降最小為目標(biāo)，提出了微通道熱沉的基于水平集方法的拓?fù)鋬?yōu)化策略。裴元帥等［8］采用變密度方法設(shè)計(jì)了風(fēng)冷熱沉結(jié)構(gòu)。Haertel 和Nellis［9］將基于密度變量的拓?fù)鋬?yōu)化方法應(yīng)用于干冷式電廠冷凝器氣側(cè)表面的設(shè)計(jì)，建立了穩(wěn)態(tài)熱流耦合充分發(fā)展的內(nèi)流拓?fù)鋬?yōu)化模型。Koga 等［10］首次開展了拓?fù)鋬?yōu)化微通道熱沉的實(shí)驗(yàn)研究，驗(yàn)證了拓?fù)錁?gòu)型數(shù)值結(jié)果的正確性。文獻(xiàn)［11～12］進(jìn)一步進(jìn)行了實(shí)驗(yàn)研究，驗(yàn)證了數(shù)值結(jié)果的有效性。Kontoleontos 等［13］將拓?fù)鋬?yōu)化方法推廣應(yīng)用于S-A湍流模型對流傳熱問題。此外，針對不同流態(tài)、離散求解方法、伴隨格式以及優(yōu)化求解方法得到了廣泛應(yīng)用和研究［14～20］。

表1 對流傳熱拓?fù)鋬?yōu)化方法對比

針對超燃沖壓發(fā)動(dòng)機(jī)冷卻結(jié)構(gòu)的簡化模型，本文提出了流固耦合傳熱結(jié)構(gòu)的拓?fù)鋬?yōu)化方法。采用有理近似方法處理材料屬性，結(jié)合移動(dòng)漸近線方法，給出了約束條件下超燃沖壓發(fā)動(dòng)機(jī)冷卻結(jié)構(gòu)的最佳流道分布，為實(shí)際工程問題提供理論依據(jù)。

2 幾何模型及邊界條件

2.1 幾何模型

超燃沖壓發(fā)動(dòng)機(jī)二維模型設(shè)計(jì)尺寸：長度為1m，寬度為1m，其中水力直徑L 為0.2m，如圖1（a）所示。設(shè)計(jì)域由流體與固體組成，給定設(shè)計(jì)區(qū)域材料初始密度為0.5。

圖1 超燃沖壓發(fā)動(dòng)機(jī)二維幾何模型

2.2 邊界條件

控制方程為連續(xù)方程、動(dòng)量方程以及能量方程，基于非結(jié)構(gòu)網(wǎng)格，采用間斷Galerkin 有限元方法［21］離散求解偏微分方程組。為了簡化計(jì)算，僅采用半模進(jìn)行計(jì)算，對稱邊界條件設(shè)置如圖1（b）所示，固體外表面采用絕熱無滑移邊界條件，材料設(shè)計(jì)區(qū)域施加恒定熱源。設(shè)計(jì)域的初始溫度為293.15 K，設(shè)計(jì)域中的流體與固體材料屬性值如表2所示。

表2 流體、固體材料屬性參數(shù)

3 拓?fù)鋬?yōu)化方法

3.1 流動(dòng)傳熱方程

忽略重力和浮升力的影響，基于Boussinesq 假設(shè)條件下的穩(wěn)態(tài)不可壓Navier-Stokes方程為

式中：ρf為流體密度；μf為流體動(dòng)力黏度；u為流體速度；p為壓力，f為源項(xiàng)。

3.2 材料屬性插值方法

優(yōu)化結(jié)構(gòu)的形狀由材料分布決定，因此多孔材料屬性采用有理近似插值方法（Rational Approximation of Material Properties，RAMP）［14］計(jì)算。構(gòu)造材料的懲罰系數(shù)、熱導(dǎo)率、密度及比熱容為

其中：αsf為懲罰系數(shù)，ρsf為材料密度，ksf為材料熱傳導(dǎo)系數(shù)，Cpsf為材料比熱容，qγ=0.01。γ為設(shè)計(jì)變量，當(dāng)γ=1 時(shí)，αsf=αs，代表固體區(qū)域；當(dāng)γ=0時(shí)，αsf=αf，代表流體區(qū)域，因此拓?fù)鋬?yōu)化問題轉(zhuǎn)化為0-1規(guī)劃問題。

3.3 靈敏度分析

流動(dòng)傳熱耦合方程可以表達(dá)為殘差形式：

與直接微分法基本一致，但計(jì)算效率明顯提升。

3.4 拓?fù)鋬?yōu)化模型

在設(shè)計(jì)區(qū)域中給定材料初始密度、體積約束，通過溫度場控制流體物性參數(shù)，進(jìn)而由流場影響溫度場，實(shí)現(xiàn)雙向耦合，從而構(gòu)建對流傳熱雙向耦合拓?fù)淠Ｐ汀?/p>

滿足：

其中：γk為材料分布值，Sv=0.4為材料體積約束上限占比，采用GCMMA求解拓?fù)鋬?yōu)化問題［22］。

4 計(jì)算結(jié)果與分析

4.1 算法驗(yàn)證

考慮流動(dòng)問題［23］：3.5mm×5mm 的矩形設(shè)計(jì)區(qū)域，兩入口對應(yīng)兩出口，進(jìn)出口高度為1mm，長度為2mm，如圖2所示。

圖2 不同雷諾數(shù)Re的材料與速度流線分布

圖2（a）、（b）分別給出了Re 為20 與200 時(shí)，設(shè)計(jì)域最優(yōu)材料與速度流線分布。低雷諾數(shù)下，為了減少功率損耗，流體通道應(yīng)盡可能短且寬，因此在Re=20 時(shí)，形成了U 形通道；隨著雷諾數(shù)增加，當(dāng)Re=200 時(shí)，彎曲區(qū)域流體的功耗隨之增加，當(dāng)慣性占主導(dǎo)地位時(shí)，在彎道處存在較大的速度梯度，因此通道自然優(yōu)化為兩條平直通道。計(jì)算結(jié)果與文獻(xiàn)［23］保持一致，驗(yàn)證了本文算法的可靠性與準(zhǔn)確性。

4.2 優(yōu)化結(jié)果對比

圖3（a）顯示了設(shè)計(jì)域中的初始溫度分布；圖3（b）給出了拓?fù)鋬?yōu)化后的溫度分布，其中流體占比為0.6。與超燃沖壓發(fā)動(dòng)機(jī)主動(dòng)冷卻系統(tǒng)的初始溫度場相比，優(yōu)化后的冷卻結(jié)構(gòu)溫度明顯降低，能量損失顯著改善，且溫度分布更加均勻，避免了高溫集中區(qū)的發(fā)生。初始流場的最高溫度與平均溫度分別為325.74 K、313.18 K，拓?fù)浣Y(jié)構(gòu)的最高溫度與平均溫度分別為293.88 K、293.29 K，拓?fù)浣Y(jié)構(gòu)最高溫度與平均溫度分別降低了9.8%和6.35%。因此，計(jì)算結(jié)果表明拓?fù)鋬?yōu)化結(jié)構(gòu)具有良好的散熱性能。

圖3 拓?fù)鋬?yōu)化前后溫度分布

圖4 給出了拓?fù)鋬?yōu)化結(jié)構(gòu)流線與流道材料分布，（a）為速度流線分布，（b）為流道結(jié)構(gòu)材料分布。通過在設(shè)計(jì)域內(nèi)施加均勻熱源，實(shí)現(xiàn)超燃沖壓發(fā)動(dòng)機(jī)的底部加熱效應(yīng)。進(jìn)出口壓差為2 Pa，拓?fù)鋬?yōu)化形成的通道結(jié)構(gòu)包括：從進(jìn)出口延伸的寬較的主通道，主通道內(nèi)流速較高，相鄰主通道之間形成多條細(xì)支通道，使冷卻通道覆蓋整個(gè)設(shè)計(jì)區(qū)域。

圖4 拓?fù)鋬?yōu)化速度流線與材料分布

由于入口與出口的壓降固定，整個(gè)通道結(jié)構(gòu)表現(xiàn)為主通道的寬度大于分支通道的寬度，減少了其他部分的流動(dòng)損失，通道之間的寬度差較小，流速分布更均勻，在設(shè)計(jì)域內(nèi)具有更好的冷卻效果。

圖5 給出了拓?fù)鋬?yōu)化計(jì)算中不同迭代步數(shù)（n=0，10，50，100）的材料和速度流線分布。顯然，隨著迭代次數(shù)的增加，主通道中的速度逐漸增加，流線分支的數(shù)量逐漸增加。當(dāng)?shù)綌?shù)n＞50 后，主流通道速度變化微小。流道結(jié)構(gòu)逐漸清晰，最終優(yōu)化后的結(jié)構(gòu)邊界輪廓清晰。拓?fù)鋬?yōu)化結(jié)構(gòu)明顯降低了設(shè)計(jì)域內(nèi)的平均溫度，與工程實(shí)際情況相符。

圖5 不同迭代步的材料分布和流線分布

4.3 傳熱性能對比

圖6 拓?fù)鋬?yōu)化模型收斂歷史

5 結(jié)語

本文提出了一種超燃沖壓發(fā)動(dòng)機(jī)主動(dòng)冷卻系統(tǒng)的拓?fù)鋬?yōu)化方法，并采用非結(jié)構(gòu)網(wǎng)格間斷有限元法對主動(dòng)冷卻系統(tǒng)的流動(dòng)與傳熱特性進(jìn)行了數(shù)值研究，得出以下結(jié)論：

1）結(jié)合GCMMA 和RAMP 插值方法，實(shí)現(xiàn)了主動(dòng)冷卻系統(tǒng)內(nèi)部冷卻通道的拓?fù)鋬?yōu)化結(jié)構(gòu)設(shè)計(jì)。通過經(jīng)典算例驗(yàn)證了算法的可靠性，與初始設(shè)計(jì)域相比，拓?fù)浣Y(jié)構(gòu)溫度明顯降低，分布更加均勻。

2）隨著拓?fù)鋬?yōu)化迭代次數(shù)的增加，拓?fù)浣Y(jié)構(gòu)呈現(xiàn)為主支和分支的混合結(jié)構(gòu)且分支逐漸增多，流道邊界輪廓逐漸清晰。拓?fù)鋬?yōu)化方法具有良好的收斂性，有效降低了設(shè)計(jì)區(qū)域的平均溫度，符合工程實(shí)際需求。