李 皓
(中國直升機(jī)設(shè)計(jì)研究所,江西 景德鎮(zhèn) 333000)
隨著計(jì)算機(jī)性能的飛速發(fā)展以及計(jì)算流體力學(xué)(Computational Fuild Dynamics,CFD)技術(shù)的不斷進(jìn)步,數(shù)值模擬方法獲得了蓬勃發(fā)展。鑒于其經(jīng)濟(jì)快捷的特性,在設(shè)計(jì)初期使用數(shù)值模擬來指導(dǎo)設(shè)計(jì)過程能大大降低設(shè)計(jì)成本,減少設(shè)計(jì)周期。在飛機(jī)防除冰系統(tǒng)的設(shè)計(jì)過程中,首先需要通過數(shù)值模擬來確定防護(hù)范圍、防冰所需熱流等,進(jìn)而才能開始整個(gè)防除冰系統(tǒng)的初步設(shè)計(jì)[1]。
由于飛機(jī)飛行速度快,對(duì)流換熱是熱量傳遞過程中的主要熱流項(xiàng),對(duì)結(jié)冰過程、加熱除冰過程有至關(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 湍流模型等,具有較高的求解精度。該文通過計(jì)算不同底層網(wǎng)格厚度、不同疏密網(wǎng)格、不同參考溫度、不同換熱溫差下的對(duì)流換熱系數(shù),分析參數(shù)變化對(duì)計(jì)算結(jié)果的影響,確定了更精確的計(jì)算方法。最后將采用S-A 湍流模型和Transition SST 湍流模型的計(jì)算結(jié)果與國外文獻(xiàn)中的試驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比,驗(yàn)證了該數(shù)值計(jì)算結(jié)果的準(zhǔn)確性,可為數(shù)值仿真工作人員提供一定參考。
通過牛頓冷卻公式可以得到對(duì)流換熱系數(shù)h,如公式(1)所示。
式中:q為壁面總熱流;Tw為壁面溫度;Tref為參考溫度。
根據(jù)公式(1)可知,選取不同的參考溫度會(huì)獲得不同的對(duì)流換熱系數(shù)。為盡可能獲得準(zhǔn)確值,應(yīng)該選取更加貼近實(shí)際的計(jì)算方法。因此該文分別選取來流靜溫T∞、來流總溫Tt和絕熱壁溫Tad作為參考溫度,來計(jì)算對(duì)流換熱系數(shù),如公式(2)~公式(4)所示。
式中:r為恢復(fù)系數(shù),其值為1.4;Ma為馬赫數(shù)。
采用FLUENT軟件計(jì)算對(duì)流換熱系數(shù)實(shí)質(zhì)上是計(jì)算空氣流場(chǎng),并求解Navier-Stokes 方程。空氣流場(chǎng)的數(shù)值計(jì)算受到計(jì)算域的空間離散條件的影響,在貼近壁面附近需要足夠細(xì)密的網(wǎng)格來解析黏性底層。為了獲得網(wǎng)格無關(guān)解,需要進(jìn)行網(wǎng)格無關(guān)性驗(yàn)證。
首先,該文選擇半圓型的試驗(yàn)件作為仿真模型,試驗(yàn)件半圓直徑為44.5 mm,并在后面連接一個(gè)錐狀的后體以保持連續(xù)的氣動(dòng)外形,易于計(jì)算收斂。其次,利用Gridgen 軟件來劃分結(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)格的無量綱壁面距離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。最后,通過壁面總熱流來分析網(wǎng)格無關(guān)性,通過計(jì)算得到的壁面總熱流如圖2 所示。
根據(jù)圖2 可知,基于不同底層網(wǎng)格厚度計(jì)算出來的壁面總熱流存在較大的差異,可見底層網(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)格無關(guān)性。
圖2 不同底層網(wǎng)格厚度下的壁面總熱流
除了底層網(wǎng)格厚度對(duì)計(jì)算結(jié)果有較大影響之外,網(wǎng)格的疏密也會(huì)影響計(jì)算的準(zhǔn)確性。因此為了獲得網(wǎng)格無關(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ì),熱流曲線呈折線狀,過渡不夠平滑。精細(xì)網(wǎng)格的網(wǎng)格數(shù)量較多,反而出現(xiàn)了鋸齒形的波動(dòng),可能是由網(wǎng)格劃分的質(zhì)量引起的??梢娭械葦?shù)量網(wǎng)格的計(jì)算結(jié)果的連續(xù)性最好。
圖3 不同疏密網(wǎng)格下的壁面總熱流
該文選擇來流靜溫作為參考溫度來計(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ì)算值。
顯然,簡單地將來流靜溫作為參考溫度來計(jì)算對(duì)流換熱系數(shù)會(huì)存在較大誤差,為了消除氣動(dòng)加熱的影響,考慮將來流總溫作為參考溫度,在相同條件下計(jì)算得到的對(duì)流換熱系數(shù)曲線如圖5 所示。
根據(jù)圖5 可知,將來流總溫作為參考溫度進(jìn)行計(jì)算后,消除了駐點(diǎn)附近的計(jì)算誤差,基于不同壁面溫度計(jì)算得到的駐點(diǎn)附近的對(duì)流換熱系數(shù)完全相同。因?yàn)樵隈v點(diǎn)附近空氣速度接近于0 m/s,可以認(rèn)為實(shí)際溫度就是來流總溫。而在遠(yuǎn)離駐點(diǎn)的下游位置,將來流總溫作為參考溫度進(jìn)行計(jì)算則不夠準(zhǔn)確。根據(jù)圖5 也可知,遠(yuǎn)離駐點(diǎn)位置在換熱溫差較小時(shí)存在較大相對(duì)誤差??偨Y(jié)圖4 和圖5 的方法可知,采用來流靜溫作為壁面溫度會(huì)導(dǎo)致駐點(diǎn)附近的溫度與實(shí)際不符,而采用來流總溫作為壁面溫度則會(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)為在壁面附近的空氣溫度與壁面溫度相同,再將該溫度作為參考溫度來計(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ù)的簡化算法,無須空氣流場(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)無量綱距離弧長/弦長
式中:q為壁面總熱流;Tad為絕熱條件下的壁面溫度。
為了驗(yàn)證該文對(duì)流換熱系數(shù)計(jì)算的準(zhǔn)確性,該文選取了美國國家航空航天局(NASA)在Lewis 冰風(fēng)洞中對(duì)NACA0012翼型所做的對(duì)流換熱試驗(yàn)作為對(duì)比。該翼型弦長為0.533 m,取其中2 個(gè)算例作為驗(yà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)格無關(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é)果
該文采用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ì)量。過于密集的網(wǎng)格可能造成網(wǎng)格質(zhì)量的下降,過于稀疏的網(wǎng)格則會(huì)引入較大的計(jì)算誤差。2) 底層網(wǎng)格厚度對(duì)對(duì)流換熱系數(shù)計(jì)算準(zhǔn)確性影響很大,應(yīng)盡可能保證網(wǎng)格的無量綱壁面距離y+接近或小于1。3) 在仿真時(shí)盡可能增大換熱溫差,以縮小相對(duì)誤差。4) 選擇不同的參考溫度會(huì)影響對(duì)流換熱系數(shù)計(jì)算結(jié)果。應(yīng)采用絕熱壁面計(jì)算氣動(dòng)加熱作用下的壁面溫度,用該溫度作為參考溫度來計(jì)算對(duì)流換熱系數(shù),這樣能消除由于氣動(dòng)加熱引起的誤差。5) 不同湍流模型對(duì)計(jì)算結(jié)果也有影響。在該文的算例下,Transition SST 模型的計(jì)算結(jié)果表現(xiàn)出了更好的符合性。