Ra=104時方腔內熱圓管自然對流的位置優化

石宏巖，高 亮，張 偉，李鑫陽

（1.中國科學院大學工程科學學院，北京 100049；2.赤峰學院 建筑與機械工程學院，內蒙古 赤峰 024000）

Ra=104時方腔內熱圓管自然對流的位置優化

石宏巖1,2，高 亮2，張 偉2，李鑫陽2

（1.中國科學院大學工程科學學院，北京 100049；2.赤峰學院 建筑與機械工程學院，內蒙古 赤峰 024000）

本文采用計算流體力學（CFD）的模擬方法，計算了單個熱圓管在方腔內自然對流的換熱過程.在瑞利數Ra=104，普朗特數Pr=0.7070時，對于圓管在豎直軸（Y軸）不同的位置分別得到了穩態溫度場及流線場，同時進行了相關的物理分析.準確地給出了熱圓管的平均努塞爾數Nu隨著Y軸位置變化的經驗公式，進而可以直接的計算出自然對流換熱系數，在X/D=a-7時出現最大的Nu，可以知道方腔內熱圓管換熱的最佳位置，同時也可以為工業應用提供優化的參考價值.

自然對流；瑞利數；努塞爾數；數值模擬

自然對流是一種在熱能轉換與傳遞過程中常見的對流換熱方式，是由流體中的浮力引起的，沒有受迫速度[1].浮力是與密度成比例的物體力作用在有密度梯度的流體上而產生的，流體的密度與溫度有關，通常會隨著溫度的升高而降低，當重力作用在有溫度梯度的流體上時，就會發生自然對流換熱[2].自然對流獨特的物理起因和換熱特性，再加上其流動過程中不需要外力做功，相比于強迫對流更具有經濟性，使自然對流換熱長期以來一直是傳熱學研究的重點、熱點之一[2,3].大部分對流換熱需要解決的問題是找到對流換熱系數，在自然對流中Nu可以直接的反映出對流換熱系數，故而得到準確Nu的關聯公式是自然對流需要解決的最基本問題，也是本文的核心點.

在中國的北方供暖是非常普遍的現象，在高溫的暖通管道通入到房間后，此時將房間看為封閉腔體，而將暖通管道看為單個的熱圓管，熱圓管在房間內的位置直接的影響房間內的換熱情形，進而影響人的舒適程度.在確定房間內的最適宜供熱溫度以后，熱圓管在房間中的什么位置換熱最大就直接的影響到供暖所需熱量，從而決定熱量的供給.故而尋找房間內換熱的最佳位置是需要解決的生活實際問題.本文將房間的實際三維模型在水平方向看為無限長，其長度不影響其他的兩個方向.將三維簡化為二維的方腔，對單個熱圓管的位置進行優化，進而得到最佳的換熱位置.

1 數理模型及計算方法

本文的模擬結構如圖1所示，方腔的邊界長度L是熱圓管直徑D的11倍，熱圓管的位置沿著Y軸變化，而變量x是自圓管上表面到方腔上邊界的距離.在模型中熱圓管表面是定溫的，即TH，周圍靜止空氣的溫度也設為一定值，即，由于流動邊界和圓管附近的溫度不同TH＞TC，圓管附近的空氣與圓管表面存在一定的溫度梯度，在氣體的熱膨脹性下，圓管附近被加熱的空氣密度將會降低，導致熱空氣與遠處的冷空氣之間存在一個密度梯度，進而在空氣中產生浮力，空氣圍繞圓管附近發生對流傳熱[4,5].其工質環境可以被認為是二維、層流、穩態、不可壓縮，流體的熱物性參數等是不依賴于溫度的獨立常量.另外，流體的熱物性參數都用膜溫TM()的值，且除Y方向動量方程中的密度項外，其余變量都假設與溫度無關[2,3].在溫差較小的范圍內，它們之間的關系可由Boussinesq近似假設[1]得到：

式中ρ是工質的密度，ρ∞是無窮遠處的密度，β是體積熱膨脹系數，在常溫情況下它的定義是：

由于流動的溫差極小，β可以按照理想氣體方程，將其設定為β=1/TM.

圖1 模擬模型的示意圖

假設工質是不可壓縮的，并且忽略能量方程中的粘性耗散和輻射.控制方程組的數學表達式如下所示[1-3]，

連續性方程：

X軸的動量方程：

Y軸的動量方程：

能量方程：

上述的控制方程中，u、v分別是x、y方向的速度，p是相對于無窮遠處標準大氣壓 p∞的表壓力，▽2=?2/?x2+?2/?y2,g、?、α是重力加速度、工質的運動粘度和熱擴散系數.

綜合上述數理模型及無滑移條件，設置的邊界條件[4,5]如下：

圓管表面：無滑移和定溫邊界條件：

方腔四周：無滑移和定溫邊界條件：

引入變量圓管直徑D，對上述控制方程及邊界條件中的所有變量進行無量綱化：

帶入方程中，得到無量綱數：

瑞利數（Ra）:

其中,Gr是格拉曉夫數，代表著作用在流體上的浮力與黏性力之比，定義式為 Gr=gβΔTD3/?2；Pr=?/α.Ra是衡量自然對流傳熱強度不可避免的參數，通常其值越大，表明圓管周圍的流動越劇烈，本文選擇Ra=104.

努塞爾數（Nu）：衡量對流熱傳遞速率的無量綱參數.圓管表面某個點的Nu稱為局部努塞爾數，其定義為：

其中，h、k分別是對流換熱系數和熱導率，r表示熱圓管表面向外法線方向的單位法向量，θ則表示在圓管表面上的點.將圓管表面上所有點的Nuθ積分求和再平均得到的就是平均努塞爾數，其定義為：

上述控制方程在CFD里通過有限體積法被離散化[5].壓力-速度耦合格式采用壓力耦合方程組的半隱式方法SIMPLE聯立求解，動量和能量方程中的對流項采用二階迎風格式進行離散化[5,6].為了獲得更加準確的結果，連續性、動量和能量方程的收斂標準都設定為10-6.

2 計算結果及分析

2.1 結果驗證

在自然對流的數值研究中，網格是影響計算結果的一個重要因素，對其合理的劃分可以得到準確的結果，減小網格數量，提高計算效率.圖2是方腔內熱圓管計算的網格示意圖，在熱圓管附近的區域，網格采用加密結構網格，有利于精確捕獲邊界層內的速度和溫度的變化，這是由于圓管附近的對流換熱邊界層變化劇烈.其自熱圓管表面向外輻射稀疏的網格，采用的網格膨脹比（即后一個網格的尺寸與前一個網格尺寸的比值）為1.01.需要說明的是：本文熱圓管附近兩個網格之間的距離δ與圓管的直徑D的比值，即為δ/D≤0.0025, 這樣就能夠保證計算結果的準確性，滿足網格質量不產生影響的基本要求[7].

圖2 方腔內熱圓管的網格示意圖

表1 網格無關性驗證

網格的數量多少是對計算結果產生直徑影響的因素，本文在變量a=5的時候對其網格的無關性進行了驗證，如表1 所示，可以看到在最小網格間距δ/D≤0.0025時對結果的影響非常微小了，其主要變化參數Nu和單位面積上的熱流密度q的相對誤差全部在0.3%以內，這充分的說明網格的密度已經達到了模擬網格無關性驗證的基本需求.

圖3 不同位置時溫度場的分布

2.2 溫度場分布

在分析傳熱問題的時候，借助于溫度場分布可以很方便的展示要研究的傳熱問題.圖3是a不同位置時溫度場的分布，在a的值發生變化的過程中，方腔內部的溫度場也會產生較為明顯的變化，圖示中原點位置始終在方腔的中心，邊界的尺寸和溫度都進行了無量綱處理.總體而言，在重力作用下，熱流體繞著圓管表面向下游流動，故而在圓管的后方形成尾流.尾流會受到方腔上方邊界的影響而產生向兩邊的分流，進而溫度場產生了蘑菇的形狀，其分流的影響區域隨著a的增加而增大，蘑菇扇區邊界尺寸也隨之增加.在a=1的時候，扇區尺寸邊界為3.2，在a=5的時候，扇區尺寸邊界為-1.3，而在a=8.5的時候，扇區尺寸邊界達到了-4.8.溫度的分布沿著熱圓管表面向外法線方向變厚，在熱圓管附近的區域受到熱流的影響較強，故而熱邊界層的梯度較大，邊界層區域較薄，而在相對于熱圓管遠處的區域，熱圓管的影響較小，故而梯度小邊界層較厚.

2.3 流線場分布

流線分布可以展示要研究的流動問題，圖4是a不同位置時流場的分布.總體而言，流動繞著圓管表面到方腔上方邊界，受其影響而產生分流，在兩邊形成對稱的旋渦，其旋渦的尺寸隨著a值的變化而變化.在熱流體到達上邊界以后產生分流，其一分部沿著邊界流動，由于邊界的溫度要低于熱圓管的溫度，熱圓管就會卷吸這部分溫度低的氣體形成裹夾，進而產生回流，形成旋渦.在a＜7的時候，其熱圓管水平方向的上方只會形成兩個對稱的旋渦，在a≥7以后會形成多個旋渦，這是由于卷吸效應的影響.在a=7時，其上方和下方均形成旋渦，只是流形和尺寸不同，這也是在方腔內熱圓管最特別的位置，由于流形的差異和旋渦的數量都會影響熱圓管的換熱能力.在熱圓管附近流場的分布，由于密度梯度大，速度邊界層的梯度較大，邊界層區域較薄，而在相對于熱圓管遠處的區域，密度梯度較小，故而速度邊界層較小.

圖4 不同位置時流場的分布

2.4 平均努塞爾數

Nu是衡量幾何結構整體平均傳熱速率的一個重要參數，在實際工程的計算中通常是不可或缺的.對于工業應用，Nu可以直接的反映出熱圓管的整體換熱能力強弱，在知道圓管直徑，工質物性的時候，根據公式（11）、（12）可以得知，如果知道Nu的公式就可以得出對流換熱系數，而自然對流的核心問題就是計算出對流換熱系數.本文對Ra=104，Pr=0.7070時的不同a的值進行公式擬合，根據泰勒公式[8]，可以知道一個函數均可以用多項式的形式進行表達，本文的擬合公式就將Nu擬合為多項式的形式，如圖5所示：

擬合公式的可決系數為0.9993，其中擬合公式與模擬數據的最大相對誤差出現在a=7的時候，為0.08%，可以看出擬合的公式是相當精確的.在圖5中可以看到，a∈1-9的時候，Nu的數值先隨著a值的增大而增加，在到達最大值以后，Nu與a呈負變化關系，其最大值出現在7處，Nu=4.3538（按公式13最大值在6.9處，為4.3539，與7處的值非常接近，所以本文直接選擇7處的值為最大值）這也符合上文的流場分布的分析.

圖5 模擬的數據及擬合公式

3 結論

本文利用CFD數值模擬的方式，在，Ra=10-4，Pr=0.7070時，準確的給出Nu隨著不同a值的關聯式，進而可以得到最佳的位置和最大的對流換熱系數，為方腔內的優化參數提供基礎.結果表明：

1）尾流會受到方腔上方邊界的影響而產生分流，溫度場形成蘑菇的形狀，蘑菇扇區邊界尺寸與a的值呈正相關性.

2）在a＜7時，熱圓管水平方向的上方只會形成兩個對稱的旋渦，在a≥7以后會形成多個旋渦.

3）Nu的數值隨著a值的增大先增加后減小，在7時出現最大值，Nu=4.3537也就是熱圓管在方腔中的最佳位置，對流換熱系數最大.

〔1〕Incropera F P,Dewitt D P,Bergman T L,et al.Fundamentals of Heat and Mass Transfer [M].Seventh ed.Indiana:John Wiley&Sons,2011:595-630.

〔2〕楊世銘,陶文銓.傳熱學 [M].北京:高等教育出版社,2006.197-276.

〔3〕Holman J P.Heat Transfer[M].Tenth ed.New York:Mcgraw-Hill Education,2010:215-375.

〔4〕田麗亭,賈篤雨,閔春華.封閉方腔內叉排管間熔融鹽自然對流換熱[J].熱能動力工程，2017,32(4):86-92.

〔5〕Liu J,Liu H,Zhen Q,et al.Numerical investigation of the laminar natural convection heat transfer from two horizontally attached horizontal cylinders[J].International Journal of Heat and Mass Transfer,2017,104:517-532.

〔6〕陶文銓.數值傳熱學[M].西安:西安交通大學出版社,2001.

〔7〕Shyam R,Sasmal C,Chhabra R P,Natural convection heat transfer from two vertically aligned circular cylinders in power-law fluids [J].International Journal of Heat and Mass Transfer,2013,64:1127-1152.

〔8〕菲赫金哥爾茨.微積分學教程[M].北京:高等教育出版社,2006.205-220.

TK124

A

1673-260X（2017）11-0035-03

2017-08-12