馬兵善,王 燁,趙興杰
(1. 蘭州交通大學 環境與市政工程學院,蘭州 730070;2. 蘭州理工大學 土木工程學院,蘭州 730050;3.蘭州交通大學 鐵道車輛熱工教育部重點實驗室,蘭州 730070)
在實際工業工程中,如建筑熱環境、電子元器件的冷卻、太陽能集熱器等方面均存在湍流自然對流現象.為深入理解其傳熱機理,許多學者進行了大量的研究[1-3].早期的研究經歷了從純自然對流換熱到考慮壁面輻射效應對換熱特性的影響過程[4-5],后來發展到同時考慮壁面輻射和內熱源條件下的封閉腔內湍流自然對流換熱研究[6-7].Kuznetsov等[8]采用類似直接數值計算的方法,對局部輻射板加熱的封閉腔內復雜傳熱(湍流自然對流、導熱和表面熱輻射)機制進行了研究.Miroshnichenko等[9]對水平放置、上下底面絕熱、含有內熱源的封閉腔內的復雜傳熱的研究發現:熱源表面的對流努塞爾數是表面發射率的減函數,壁面輻射有利于腔內的傳熱過程.Miroshnichenko等[10]對水平放置含內熱源的封閉腔內湍流自然對流——壁面熱輻射進行了數值研究,發現壁面發射率使壁面平均總努塞爾數增大,熱源位置引起了腔體內流動方式和熱羽流的顯著改變.Miroshnichenko和Sheremet[11]對內熱源位于底面中心的傾斜方腔內存在壁面熱輻射的湍流自然對流進行了數值分析,發現熱源表面輻射努塞爾數隨著腔體傾角的增大而減小.Sivaraj等[12]對中央放置方形等熱流源的傾斜方腔內壁面輻射——自然對流耦合傳熱進行了參數研究,發現熱源表面的平均對流努塞爾數是關于傾角的減函數,而且,不同壁面發射率所得熱源表面平均對流努塞爾數隨傾角的變化趨勢一致.但這些研究的局限性或者是其單一的邊界條件或者是熱源位置固定,與工程實際存在較大差異.
實際工程中對封閉空間內散熱設備冷卻問題進行優化設計時往往要考慮多種邊界條件作用下的耦合傳熱過程,其中,冷卻效果與散熱設備的位置及安裝角度密切相關.同時考慮腔體傾角、內熱源位置以及輻射效應對腔內對流換熱特性影響的文獻,目前還未見報道.基于這一研究現狀,為了獲得含內熱源封閉腔內的傳熱特性和最佳對流冷卻效果,本文對不同傾角下的封閉腔內自然對流換熱過程進行了數值分析,旨在揭示復雜條件下封閉腔內自然對流換熱機制,為改善工程實際中的對流冷卻效果提供理論參考.


圖1 物理模型Fig.1 Physical model
用于求解流動與傳熱的控制方程如下:
連續性方程:
·V=0,
(1)
動量方程:
(2)
(3)
能量方程:
(4)
湍流動能k方程:
(5)
湍流動能耗散率ε方程:
(6)
上式中,t為時間s;T為溫度K;T0為參考溫度,K;p為壓力Pa;u,v分別為x、y方向的速度m/s;ρ為密度kg/m3;λ為導熱系數W/(m·K);cp為定壓比熱容J/(kg·K);μeff為有效湍流黏度;β為熱膨脹系數K-1;φ為腔體的傾角°.
DO輻射模型
(7)

邊界條件:
熱壁面溫度恒為Th=323.15 K,冷壁面溫度恒為Tc=283.15 K,內熱源表面溫度恒為Ts=343.15 K,外界環境溫度T0=(Th+Tc)/2,外界與腔體外壁面間的對流換熱系數為h=7 W/(m2·K)[13],所有氣固交界面均為速度無滑移邊界條件.
腔體中熱壁面的局部對流努塞爾數Nuconv和局部輻射努塞爾數Nurad分別定義如下:
(8)
(9)
其中:qconv為壁面的對流熱流密度,W/m2;qrad為離開壁面的凈輻射熱流密度,W/m2;ΔT為計算溫差;H為特征長度.
幾何平均努塞爾數:
(10)
(11)
總平均努塞爾數:
(12)
利用SIMPLE算法對所研究問題的控制方程進行求解,采用非均分交錯網格.對流項采用QUICK格式,擴散項采用中心差分格式進行離散.為了驗證算法的正確性,對文獻[14]所研究的方腔內湍流自然對流問題進行了數值計算.計算結果如圖2所示,腔體Y=0.5水平線上無量綱垂直速度和無量綱溫度與文獻[14]中的實驗數據符合較好,因此,該方法可用于后續計算.對于本文所研究的問題,為了獲得網格無關解,分別采用網格數為150×150、200×200、250×250、300×300的四套網格進行了數值試驗,得到了如圖3所示計算時長為500 s時腔體X=0.5截面上的水平方向上的速度分布.可以看出,不同網格數所得曲線吻合較好,認為本文模型可以獲得網格無關解.考慮計算的經濟性,后續計算中取網格數150×150,時間步長τ=0.05 s.

圖2 數值方法驗證Fig.2 Validation of the numerical method
對含內熱源的傾斜方腔,為了研究腔體傾角φ、內熱源距腔體熱壁面的無量綱長度D(D=l/L)及腔體表面發射率對湍流自然對流的影響,計算中,Ra=1.58×109,Pr=0.71,傾角φ分別取0°、15°、30°、45°、60°、75°、90°,熱源位置的無量綱長度D分別取0.1、0.2、0.3、0.4、0.5、0.6、0.7,壁面發射率分別取0、0.3、0.6、0.9.下面給出計算時長為t=1 000 s時的結果.

圖3 網格獨立性驗證Fig.3 Validation of the mesh independence

圖4 不同φ時等流函數線(左)和等溫線圖(右)Fig.4 Streamlines(left) and isotherms(right) with φ 000 s)

圖5 熱壁面平均努塞爾數隨傾角的變化Fig.5 Variation of average Nusselt number on the hot wall with inclination angle
由3.1節分析可知,腔體傾角為75°時對流冷卻效果最佳.所以,下文只討論腔體傾角為75°時的情況.
圖6為φ=75°時不同無量綱長度D下的等流函數線圖、等溫線圖、等湍流動能線圖及等湍流動能耗散率線圖.由圖6(a)可以看出,腔體中心區域的溫度變化對傾角變化不是很敏感,即保持了穩定的熱層結構.當D=0.1時,在熱源的上部頂角出現了熱羽流,這主要是由于在浮升力的作用下,此處形成了兩個反向旋轉的渦,熱羽流產生的位置在兩渦的交界處.隨著D的增大,熱羽流向著冷壁面的方向擴展,是因為熱源與冷壁面間的距離減小引起冷面附近的溫度梯度增大,熱源與冷壁面間的換熱因此被強化.
圖6(b)為腔體內的等流函數線的分布.腔體呈雙渦結構,D=0.1時,左下側逆時針旋轉的渦相比于右上側渦強度較大,但所占空間較小,主要是因為熱源距離熱壁面較近,受到熱羽流的阻隔所導致.隨著D增大,左側渦所占腔體空間進一步擴展,成為主渦,D=0.2時取得最大流函數值,隨著D的增大最大流函數值逐漸減小,是由于熱源逐漸靠近冷壁面導致冷壁面附近流體在熱源輻射作用下溫度升高,浮升力減小,流動減弱.
圖6(c)為腔體內湍流動能的等值線.D=0.1時,腔體內的湍流區域主要集中在熱源的上部頂角和腔體的左側,湍流動能的核心區與產生熱羽流的位置一致.隨著D的增大,湍流動能的核心向冷壁面方向轉移,腔體中心區域的湍流動能變化較小,這與溫度場結構有關.圖6(d)為腔體中湍流動能耗散率的等值線圖.D=0.1時,在熱源上部頂角出現了一個渦結構,此處氣流動能耗散較大,也正是兩渦的交界面(6(b)所示);D>0.3時,腔體中熱源上方和熱、冷壁面附近的湍流動能耗散率的梯度減小,說明腔體中的湍流強度減弱,流動變緩,腔體中熱源位置對流動換熱的影響顯著.




圖6 (a)等溫線圖、(b)等流函數線、(c)等湍流動能線圖和(d)湍流動能耗散率Fig.6 Isotherms (a),streamlines (b),turbulent kinetic energy fields(c) and turbulent kinetic energy dissipation rate field(d) at
如圖7(b)所示,不同熱源位置對腔內速度場分布影響顯著.D=0.1時豎向速度在X=0.4處取得正的最大值,對應熱羽流的形成位置;D=0.5時豎向速度在X=0.1處取得負的最大值.熱源右移,熱壁面附近流體被加速,這與圖7(a)中熱源右移有利于熱壁面的冷卻散熱是一致的.0.2≤D≤0.6時熱源位置對豎向速度的影響很微弱,這與圖7(a)中溫度場的分布特征一致,這一點體現了自然對流換熱中溫度場與速度場間的耦合關系.

圖7 Y=0.5水平線上溫度和速度隨D的變化Fig.7 Variations of temperature and vertical velocity with D at




圖8 熱壁面隨D的變化Fig.8 Variation of average Nusselt number on the hot

圖9 熱壁面Nuconv及Nurad隨的變化(φ=75°,D=0.6)Fig.9 Variations of Nuconv and Nurad on the hot wall with at φ=75°,D=0.6
對含內熱源的傾斜封閉方腔內輻射—對流耦合傳熱特性進行了數值計算,得到了腔體傾角、熱源位置及壁面發射率對湍流自然對流傳熱的影響規律,具體如下:

2) 輻射效應對總換熱能力的貢獻受腔體傾角的影響微弱,熱源對熱壁面冷卻效果的貢獻隨其距熱壁面無量綱長度D增大而增大.
3) 壁面發射率的增大對熱壁面邊界層起始段的對流換熱以及輻射換熱激勵作用顯著,較大的壁面發射率使得Y=0.5附近傳熱能力產生了階躍性突增,但熱邊界層下游區域的輻射換熱被削弱.