對(duì)流換熱數(shù)值模擬

李 皓

（中國(guó)直升機(jī)設(shè)計(jì)研究所，江西 景德鎮(zhèn) 333000）

0 引言

隨著計(jì)算機(jī)性能的飛速發(fā)展以及計(jì)算流體力學(xué)（Computational Fuild Dynamics，CFD）技術(shù)的不斷進(jìn)步，數(shù)值模擬方法獲得了蓬勃發(fā)展。鑒于其經(jīng)濟(jì)快捷的特性，在設(shè)計(jì)初期使用數(shù)值模擬來(lái)指導(dǎo)設(shè)計(jì)過(guò)程能大大降低設(shè)計(jì)成本，減少設(shè)計(jì)周期。在飛機(jī)防除冰系統(tǒng)的設(shè)計(jì)過(guò)程中，首先需要通過(guò)數(shù)值模擬來(lái)確定防護(hù)范圍、防冰所需熱流等，進(jìn)而才能開(kāi)始整個(gè)防除冰系統(tǒng)的初步設(shè)計(jì)[1]。

由于飛機(jī)飛行速度快，對(duì)流換熱是熱量傳遞過(guò)程中的主要熱流項(xiàng)，對(duì)結(jié)冰過(guò)程、加熱除冰過(guò)程有至關(guān)重要的影響。能否獲得準(zhǔn)確的對(duì)流換熱系數(shù)會(huì)直接關(guān)系到防除冰系統(tǒng)設(shè)計(jì)的準(zhǔn)確性，因此進(jìn)行對(duì)流換熱系數(shù)計(jì)算準(zhǔn)確性的研究是非常必要的[2]。

該文的數(shù)值計(jì)算工作基于FLUENT 軟件。作為如今CFD 領(lǐng)域最全面的軟件之一，F(xiàn)LUENT 提供了多種求解方法和湍流模型，包括Spalart-Allmaras（S-A）湍流模型和Transition SST 湍流模型等，具有較高的求解精度。該文通過(guò)計(jì)算不同底層網(wǎng)格厚度、不同疏密網(wǎng)格、不同參考溫度、不同換熱溫差下的對(duì)流換熱系數(shù)，分析參數(shù)變化對(duì)計(jì)算結(jié)果的影響，確定了更精確的計(jì)算方法。最后將采用S-A 湍流模型和Transition SST 湍流模型的計(jì)算結(jié)果與國(guó)外文獻(xiàn)中的試驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比，驗(yàn)證了該數(shù)值計(jì)算結(jié)果的準(zhǔn)確性，可為數(shù)值仿真工作人員提供一定參考。

1 對(duì)流換熱系數(shù)的計(jì)算

通過(guò)牛頓冷卻公式可以得到對(duì)流換熱系數(shù)h，如公式（1）所示。

式中：q為壁面總熱流；Tw為壁面溫度；Tref為參考溫度。

根據(jù)公式（1）可知，選取不同的參考溫度會(huì)獲得不同的對(duì)流換熱系數(shù)。為盡可能獲得準(zhǔn)確值，應(yīng)該選取更加貼近實(shí)際的計(jì)算方法。因此該文分別選取來(lái)流靜溫T∞、來(lái)流總溫Tt和絕熱壁溫Tad作為參考溫度，來(lái)計(jì)算對(duì)流換熱系數(shù)，如公式（2）～公式（4）所示。

式中：r為恢復(fù)系數(shù)，其值為1.4；Ma為馬赫數(shù)。

2 網(wǎng)格無(wú)關(guān)性驗(yàn)證

采用FLUENT軟件計(jì)算對(duì)流換熱系數(shù)實(shí)質(zhì)上是計(jì)算空氣流場(chǎng)，并求解Navier-Stokes 方程?？諝饬鲌?chǎng)的數(shù)值計(jì)算受到計(jì)算域的空間離散條件的影響，在貼近壁面附近需要足夠細(xì)密的網(wǎng)格來(lái)解析黏性底層。為了獲得網(wǎng)格無(wú)關(guān)解，需要進(jìn)行網(wǎng)格無(wú)關(guān)性驗(yàn)證。

首先，該文選擇半圓型的試驗(yàn)件作為仿真模型，試驗(yàn)件半圓直徑為44.5 mm，并在后面連接一個(gè)錐狀的后體以保持連續(xù)的氣動(dòng)外形，易于計(jì)算收斂。其次，利用Gridgen 軟件來(lái)劃分結(jié)構(gòu)網(wǎng)格，如圖1 所示。計(jì)算域尺寸約為7 m，只有范圍足夠大，才能保證遠(yuǎn)場(chǎng)邊界條件成立，不會(huì)影響到壁面附近的空氣流場(chǎng)。并選取0.001 mm、0.005 mm、0.01 mm、0.05 mm 和0.1 mm 共5 種不同的底層網(wǎng)格厚度進(jìn)行計(jì)算，不同網(wǎng)格的無(wú)量綱壁面距離y+在0.1～45.0。

圖1 Gridgen 劃分C 字型結(jié)構(gòu)網(wǎng)格

再次，選用Spalart-Allmaras 湍流模型計(jì)算空氣流場(chǎng)?？諝獠捎美硐霘怏w，比熱cp=1 006.43 J/（kg·K），空氣導(dǎo)熱系數(shù)λ=0.023 56 W/（m·K），動(dòng)力黏性系數(shù)μ=1.789 4×10-5kg/（m·s）。并設(shè)置壓力遠(yuǎn)場(chǎng)和等溫壁面，環(huán)境壓力P=34 500 Pa，馬赫數(shù)Ma=0.25，壁面溫度Tw=273.15 K，環(huán)境溫度Tinf=253.15 K。最后，通過(guò)壁面總熱流來(lái)分析網(wǎng)格無(wú)關(guān)性，通過(guò)計(jì)算得到的壁面總熱流如圖2 所示。

根據(jù)圖2 可知，基于不同底層網(wǎng)格厚度計(jì)算出來(lái)的壁面總熱流存在較大的差異，可見(jiàn)底層網(wǎng)格厚度能顯著影響計(jì)算結(jié)果。而隨著底層網(wǎng)格厚度減少，計(jì)算結(jié)果也逐漸趨于穩(wěn)定，厚度為0.001 mm 和0.005 mm 網(wǎng)格的計(jì)算結(jié)果幾乎完全重合。由于0.005 mm 厚度的網(wǎng)格的y+≈1，因此可以認(rèn)為0.005 mm 厚度的網(wǎng)格已經(jīng)足夠細(xì)密，計(jì)算結(jié)果即表現(xiàn)出網(wǎng)格無(wú)關(guān)性。

圖2 不同底層網(wǎng)格厚度下的壁面總熱流

除了底層網(wǎng)格厚度對(duì)計(jì)算結(jié)果有較大影響之外，網(wǎng)格的疏密也會(huì)影響計(jì)算的準(zhǔn)確性。因此為了獲得網(wǎng)格無(wú)關(guān)解，還需要驗(yàn)證網(wǎng)格疏密對(duì)計(jì)算結(jié)果的影響。該文選取3套不同疏密的網(wǎng)格，底層網(wǎng)格厚度均為0.005 mm，網(wǎng)格數(shù)量分別為70 850 個(gè)、58 000 個(gè)和30 240 個(gè)。計(jì)算得到的壁面總熱流結(jié)果如圖3 所示。

根據(jù)圖3 可知，基于3 套不同疏密的網(wǎng)格計(jì)算得到的壁面總熱流在總體趨勢(shì)和范圍上差距不大，但細(xì)節(jié)上存在較大的不同。粗糙網(wǎng)格的網(wǎng)格數(shù)量少，劃分不夠精細(xì)，熱流曲線呈折線狀，過(guò)渡不夠平滑。精細(xì)網(wǎng)格的網(wǎng)格數(shù)量較多，反而出現(xiàn)了鋸齒形的波動(dòng)，可能是由網(wǎng)格劃分的質(zhì)量引起的?？梢?jiàn)中等數(shù)量網(wǎng)格的計(jì)算結(jié)果的連續(xù)性最好。

圖3 不同疏密網(wǎng)格下的壁面總熱流

3 參考溫度對(duì)對(duì)流換熱系數(shù)的影響

該文選擇來(lái)流靜溫作為參考溫度來(lái)計(jì)算對(duì)流換熱系數(shù)h，計(jì)算條件模型和網(wǎng)格同第2 節(jié)，底層網(wǎng)格厚度0.005 mm，網(wǎng)格數(shù)量適中。采用Transition SST 湍流模型，空氣采用理想氣體，并設(shè)置壓力遠(yuǎn)場(chǎng)和等溫壁面，環(huán)境壓力P=34 500 Pa，馬赫數(shù)Ma=0.25，攻角為0°，環(huán)境溫度T∞=253.15 K。只改變壁面溫度，分別取Tw為263.15 K、273.15 K、283.15 K、293.15 K、303.15 K 和313.15 K。計(jì)算結(jié)果如圖4 所示。

根據(jù)圖4 可知，基于不同溫差（Tw-Tinf）計(jì)算得到的對(duì)流換熱系數(shù)相差極大，即隨著壁面溫度的升高，駐點(diǎn)附近的對(duì)流換熱系數(shù)逐漸增大，而遠(yuǎn)離駐點(diǎn)的下游處對(duì)流換熱系數(shù)逐漸減小。并且隨著溫差的增大，計(jì)算結(jié)果的相對(duì)誤差逐漸減小。產(chǎn)生該結(jié)果的原因主要在于未考慮氣動(dòng)加熱的影響，由于駐點(diǎn)附近的空氣速度顯著降低，動(dòng)能轉(zhuǎn)化為熱能，導(dǎo)致溫度有所升高，因此實(shí)際的換熱溫差小于計(jì)算值。

顯然，簡(jiǎn)單地將來(lái)流靜溫作為參考溫度來(lái)計(jì)算對(duì)流換熱系數(shù)會(huì)存在較大誤差，為了消除氣動(dòng)加熱的影響，考慮將來(lái)流總溫作為參考溫度，在相同條件下計(jì)算得到的對(duì)流換熱系數(shù)曲線如圖5 所示。

根據(jù)圖5 可知，將來(lái)流總溫作為參考溫度進(jìn)行計(jì)算后，消除了駐點(diǎn)附近的計(jì)算誤差，基于不同壁面溫度計(jì)算得到的駐點(diǎn)附近的對(duì)流換熱系數(shù)完全相同。因?yàn)樵隈v點(diǎn)附近空氣速度接近于0 m/s，可以認(rèn)為實(shí)際溫度就是來(lái)流總溫。而在遠(yuǎn)離駐點(diǎn)的下游位置，將來(lái)流總溫作為參考溫度進(jìn)行計(jì)算則不夠準(zhǔn)確。根據(jù)圖5 也可知，遠(yuǎn)離駐點(diǎn)位置在換熱溫差較小時(shí)存在較大相對(duì)誤差?？偨Y(jié)圖4 和圖5 的方法可知，采用來(lái)流靜溫作為壁面溫度會(huì)導(dǎo)致駐點(diǎn)附近的溫度與實(shí)際不符，而采用來(lái)流總溫作為壁面溫度則會(huì)導(dǎo)致遠(yuǎn)離駐點(diǎn)區(qū)域的溫度與實(shí)際不符。也即這2 種計(jì)算方法都存在一定不足，均會(huì)引入較大的相對(duì)誤差，因此采用等溫壁面進(jìn)行計(jì)算是不合適的。

圖4 不同換熱溫差下對(duì)流換熱系數(shù)計(jì)算結(jié)果（Transition SST 模型，壓力遠(yuǎn)場(chǎng)和等溫壁面邊界條件）

為了減少誤差，如果采用絕熱壁面計(jì)算出在空氣流場(chǎng)作用下壁面的實(shí)際溫度，并認(rèn)為在壁面附近的空氣溫度與壁面溫度相同，再將該溫度作為參考溫度來(lái)計(jì)算壁面的對(duì)流換熱系數(shù)，就可以消除氣動(dòng)加熱引起的誤差，改善結(jié)果的準(zhǔn)確性。計(jì)算結(jié)果如圖6 所示。

根據(jù)圖6 可知，與圖4、圖5 的2 種方法相比，圖6 的計(jì)算結(jié)果有極大改善。駐點(diǎn)附近計(jì)算結(jié)果完全相同，遠(yuǎn)離駐點(diǎn)的下游區(qū)域的相對(duì)誤差也顯著減小。因此采用該種方法作為對(duì)流換熱系數(shù)的簡(jiǎn)化算法，無(wú)須空氣流場(chǎng)的反復(fù)迭代，誤差又能控制在可接受范圍內(nèi)。其對(duì)流換熱系數(shù)h的計(jì)算如公式（5）所示。

圖5 不同換熱溫差下對(duì)流換熱系數(shù)計(jì)算結(jié)果（Transition SST 模型，壓力遠(yuǎn)場(chǎng)和等溫壁面邊界條件，考慮滯止溫升）

圖6 不同換熱溫差下對(duì)流換熱系數(shù)計(jì)算結(jié)果（Transition SST 模型，壓力遠(yuǎn)場(chǎng)和絕熱壁面邊界條件）

駐點(diǎn)無(wú)量綱距離弧長(zhǎng)/弦長(zhǎng)

式中：q為壁面總熱流；Tad為絕熱條件下的壁面溫度。

4 算例驗(yàn)證

為了驗(yàn)證該文對(duì)流換熱系數(shù)計(jì)算的準(zhǔn)確性，該文選取了美國(guó)國(guó)家航空航天局（NASA）在Lewis 冰風(fēng)洞中對(duì)NACA0012翼型所做的對(duì)流換熱試驗(yàn)作為對(duì)比。該翼型弦長(zhǎng)為0.533 m，取其中2 個(gè)算例作為驗(yàn)證，見(jiàn)表1。計(jì)算同時(shí)采用了S-A 湍流模型和Transition SST 湍流模型，并采用遠(yuǎn)場(chǎng)邊界條件，調(diào)節(jié)馬赫數(shù)使自由流速度與試驗(yàn)數(shù)據(jù)吻合，空氣采用理想氣體，密度為1.326 kg/m3，動(dòng)力黏度為1.676 7×10-5kg/（m·s），比熱為1 006.43 J/（kg·K），空氣導(dǎo)熱系數(shù)為0.023 56 W/（m·K），環(huán)境溫度為266.45 K。為了減少相對(duì)誤差，盡可能增大換熱溫差，取壁面溫度比環(huán)境溫度高60 ℃。計(jì)算所用網(wǎng)格滿足第2 節(jié)網(wǎng)格無(wú)關(guān)性驗(yàn)證中的條件y+＜1，第一層網(wǎng)格厚度為0.005 mm。

表1 對(duì)流換熱系數(shù)計(jì)算條件

圖7、圖8 分別為算例1、算例2 對(duì)流換熱系數(shù)計(jì)算結(jié)果與文獻(xiàn)試驗(yàn)結(jié)果的比較。根據(jù)圖7、圖8 可知，不同的湍流模型計(jì)算結(jié)果存在較大差異。在駐點(diǎn)附近，S-A 湍流模型和Transition SST 湍流模型的計(jì)算結(jié)果基本相同，同時(shí)也與文獻(xiàn)試驗(yàn)結(jié)果基本吻合，很好地描述了實(shí)際的換熱情況。但隨著位置往后移動(dòng)，S-A 模型的轉(zhuǎn)捩很快出現(xiàn)，對(duì)流換熱系數(shù)不再下降，反而有微小回升，比文獻(xiàn)試驗(yàn)結(jié)果偏高。而Transition SST 模型則在整個(gè)試驗(yàn)測(cè)量范圍內(nèi)有良好的仿真效果。因此該文后面對(duì)對(duì)流換熱系數(shù)的計(jì)算全部采用Transition SST 湍流模型。

圖7 算例1 對(duì)流換熱系數(shù)計(jì)算結(jié)果

圖8 算例2 對(duì)流換熱系數(shù)計(jì)算結(jié)果

5 結(jié)論

該文采用FLUENT 軟件計(jì)算對(duì)流換熱系數(shù)，根據(jù)上述分析和比較可知，網(wǎng)格疏密、第一層網(wǎng)格的厚度、參考溫度、換熱溫差和湍流模型都會(huì)對(duì)對(duì)流換熱系數(shù)的計(jì)算結(jié)果產(chǎn)生影響。要想得到盡可能準(zhǔn)確的計(jì)算結(jié)果，應(yīng)遵循以下規(guī)律：1） 保證較好的網(wǎng)格質(zhì)量。過(guò)于密集的網(wǎng)格可能造成網(wǎng)格質(zhì)量的下降，過(guò)于稀疏的網(wǎng)格則會(huì)引入較大的計(jì)算誤差。2） 底層網(wǎng)格厚度對(duì)對(duì)流換熱系數(shù)計(jì)算準(zhǔn)確性影響很大，應(yīng)盡可能保證網(wǎng)格的無(wú)量綱壁面距離y+接近或小于1。3） 在仿真時(shí)盡可能增大換熱溫差，以縮小相對(duì)誤差。4） 選擇不同的參考溫度會(huì)影響對(duì)流換熱系數(shù)計(jì)算結(jié)果。應(yīng)采用絕熱壁面計(jì)算氣動(dòng)加熱作用下的壁面溫度，用該溫度作為參考溫度來(lái)計(jì)算對(duì)流換熱系數(shù)，這樣能消除由于氣動(dòng)加熱引起的誤差。5） 不同湍流模型對(duì)計(jì)算結(jié)果也有影響。在該文的算例下，Transition SST 模型的計(jì)算結(jié)果表現(xiàn)出了更好的符合性。