螺旋管內(nèi)單相流動(dòng)周向非均勻傳熱現(xiàn)象的數(shù)值模擬

王 瑞， 肖 瑤， 顧漢洋， 葉亞楠

(上海交通大學(xué) 核科學(xué)與工程學(xué)院， 上海 200240)

符號說明

a—螺旋管內(nèi)徑， mm

c—螺旋管螺旋直徑，mm

Cp—比定壓熱容， kJ/(kg·℃)

F—合力，N

g—重力加速度，kg/(m·s2)

h—換熱系數(shù)，W/(m2·℃)

k—導(dǎo)熱系數(shù)， W/(m·℃)

Nu—努塞爾數(shù)

Nuθ—θ處的局部努塞爾數(shù)

q—熱流密度，kW/m2

Re—雷諾數(shù)

T—局部內(nèi)壁溫度， ℃

Tb—平均流體溫度，℃

Tw—平均內(nèi)壁溫度，℃

u—流速， m/s

α—螺旋管螺旋升角，(°)

β—離心力與重力合力與豎直方向的夾角，(°)

θ—截面周向角度，(°)

Θ—無量綱溫度

ρ—密度， kg/m3

φ—重力加速度與離心力加速度之比

螺旋管由于其結(jié)構(gòu)緊湊、易于制造、換熱效率高等優(yōu)點(diǎn)被廣泛地應(yīng)用于食品工業(yè)、核工業(yè)、廢熱回收、制冷、航空航天和許多其他工業(yè)場景[1-3].在核工業(yè)中，螺旋管是近年來廣受關(guān)注的小型模塊化反應(yīng)堆蒸汽發(fā)生器的理想選型.在設(shè)計(jì)螺旋管式蒸汽發(fā)生器時(shí)，充分了解螺旋管的傳熱特性才能保證蒸汽發(fā)生器換熱能力與堆芯功率匹配使反應(yīng)堆安全運(yùn)行.因此，對螺旋管內(nèi)流動(dòng)傳熱的研究具有重大意義.許多研究指出螺旋管內(nèi)流體存在復(fù)雜的流型[4].螺旋管內(nèi)流體在離心力的作用下在橫截面上產(chǎn)生了二次流.二次流的形態(tài)隨邊界條件的變化而變化，使得流體的速度分布不均，強(qiáng)化了螺旋管內(nèi)的傳熱[5-7].由于實(shí)驗(yàn)方法不能夠詳盡地描述螺旋管內(nèi)的局部單相流動(dòng)換熱特性，所以許多相關(guān)研究采用了數(shù)值模擬方法.這些研究中最常使用的湍流模型是k-ε模型和Reynolds應(yīng)力模型.馬越等[8]使用Reynolds應(yīng)力模型模擬了高溫氣冷堆中氦氣橫掠外壁面時(shí)螺旋管截面上的壁面溫度分布，發(fā)現(xiàn)螺旋管壁面最高溫度位于截面頂部與內(nèi)側(cè)之間.史建新等[9]使用重整化群(RNG)k-ε模型模擬了特定結(jié)構(gòu)的螺旋管內(nèi)單相流體的周向流動(dòng)與傳熱分布，發(fā)現(xiàn)螺旋管外側(cè)流速大、內(nèi)壁溫低、換熱系數(shù)高.Jayakumar等[10]使用了Realizablek-ε模型清晰地描述了螺旋管分別沿軸向和徑向的局部努塞爾數(shù)波動(dòng)，同時(shí)利用計(jì)算結(jié)果研究了螺旋管結(jié)構(gòu)參數(shù)對傳熱的影響，并驗(yàn)證了已有的傳熱經(jīng)驗(yàn)關(guān)系式的有效性.

綜上所述，對螺旋管截面周向的傳熱分布研究比較初步，如文獻(xiàn)[8]發(fā)現(xiàn)了螺旋管最高壁溫位于頂部及內(nèi)側(cè)之間，卻沒有指出其具體位置及變化規(guī)律；如文獻(xiàn)[10]給出了螺旋管周向傳熱分布的特征卻沒有從受力角度去解釋成因.本文對螺旋管的周向局部特性做了更進(jìn)一步的探討，系統(tǒng)地研究了螺旋管傳熱周向分布的特征，總結(jié)了管壁周向溫度分布隨加速度比、螺旋直徑、螺旋升角、螺旋管水力學(xué)直徑等參數(shù)的變化規(guī)律及成因.

1 數(shù)值模擬

1.1 數(shù)學(xué)物理模型

所模擬的螺旋管物理結(jié)構(gòu)參數(shù)如表1所示.應(yīng)用計(jì)算流體力學(xué)(CFD)分析軟件CFX16.0進(jìn)行模擬計(jì)算，湍流模型采用Reynolds應(yīng)力模型.文獻(xiàn)[8]指出，Reynolds應(yīng)力模型更多地考慮了旋轉(zhuǎn)流動(dòng)的特性，更適合于計(jì)算旋轉(zhuǎn)流動(dòng)中的流體；而常用的k-ε模型則采用各向同性的湍動(dòng)黏度計(jì)算湍流應(yīng)力，在計(jì)算旋轉(zhuǎn)流動(dòng)方面偏差較大.同時(shí)，文獻(xiàn)[8]還驗(yàn)證了Reynolds應(yīng)力模型的計(jì)算結(jié)果與經(jīng)驗(yàn)關(guān)系式最為吻合，故本文采用Reynolds應(yīng)力模型進(jìn)行計(jì)算.

在CFX計(jì)算中advection scheme和turbulence numerical選項(xiàng)均設(shè)置為high resolution，即動(dòng)量方程對流項(xiàng)和湍流輸運(yùn)方程對流項(xiàng)的離散格式均采用比較精細(xì)的離散格式.當(dāng)質(zhì)量、動(dòng)量、能量等相關(guān)量的均方根(RMS)殘差為1×10-7時(shí)，判斷計(jì)算收斂.邊界條件設(shè)置如下：入口為恒定質(zhì)量流速入口，并設(shè)定恒定入口溫度；出口為定壓出口；壁面為無滑移的光滑壁面，壁面上設(shè)置恒定均勻熱流密度加熱；設(shè)重力加速度大小為9.8 kg/(m·s2)，方向?yàn)樨Q直向下.以水為工作介質(zhì)，水物性采用IAPWS IF97的變物性程序包進(jìn)行計(jì)算，計(jì)算區(qū)域水保持單相.計(jì)算涉及到的各工況邊界條件如表2所示.

表2 工況邊界條件Tab.2 Boundary conditions

1.2 網(wǎng)格的劃分

分別對不同結(jié)構(gòu)的螺旋管使用ICEM CFD軟件創(chuàng)建結(jié)構(gòu)化六面體網(wǎng)格.為了保證計(jì)算精度，進(jìn)行網(wǎng)格的獨(dú)立性實(shí)驗(yàn).在工況1的邊界條件下由不同數(shù)量的網(wǎng)格計(jì)算出的某截面(截面Tb=141.2 ℃)周向最高內(nèi)壁溫，如圖1所示.其中：N為體積單元網(wǎng)格數(shù)量；Tm為最高內(nèi)壁溫度.由圖1可見，隨著網(wǎng)格數(shù)量的上升，最高內(nèi)壁溫度逐漸增大.當(dāng)體積單元網(wǎng)格數(shù)量達(dá)到7×106后最高內(nèi)壁溫度基本保持不變，故認(rèn)為此時(shí)的網(wǎng)格符合獨(dú)立性要求.經(jīng)過獨(dú)立性驗(yàn)證后的各螺旋管體積單元網(wǎng)格數(shù)量如表3所示，

圖1 螺旋管1的網(wǎng)格獨(dú)立性實(shí)驗(yàn)Fig.1 The grid independence experiment of Helical Pipe 1

表3 各螺旋管的網(wǎng)格數(shù)量Tab.3 Grid numbers of different helical pipes

優(yōu)化后的螺旋管1網(wǎng)格如圖2所示.螺旋管1橫截面上總節(jié)點(diǎn)數(shù)為 1 464，第一層網(wǎng)格質(zhì)心到壁面的無量綱距離Y+接近30，符合Reynolds應(yīng)力模型的要求.

圖2 螺旋管1的網(wǎng)格結(jié)構(gòu)Fig.2 Grid structure of Helical Pipe 1

2 計(jì)算結(jié)果的實(shí)驗(yàn)驗(yàn)證

在計(jì)算換熱系數(shù)時(shí)需要用到截面平均流體溫度.由于流體的物性與溫度有關(guān)，采用文獻(xiàn)[10]的方式計(jì)算截面平均流體溫度：

(1)

式中：A為平均流體軸向截面.

截面換熱系數(shù)如下式：

式中：Tθ為θ處的局部內(nèi)壁溫；θ1，θ2為任意兩個(gè)截面的軸向角度位置.

圖3 工況1管道沿程計(jì)算和實(shí)驗(yàn)結(jié)果的對比Fig.3 Comparison of calculated and experimental results along the flow direction for boundary Condition 1

圖4 螺旋管截面上T和T′的對比Fig.4 Comparison of T and T′ on helical pipe cross section

圖5 換熱系數(shù)計(jì)算值和實(shí)驗(yàn)值比較Fig.5 Comparison of calculated and experimental heat transfer coefficient

3 結(jié)果分析與討論

截面上周向努塞爾數(shù)的計(jì)算式為

(4)

截面努塞爾數(shù)計(jì)算式為

(5)

采用歸一化努塞爾數(shù)Nuθ/Nu討論局部換熱情況.螺旋管截面的周向傳熱分布受重力和離心力的共同影響，故提出以下無量綱參數(shù)：

(6)

即φ可用于表征重力與離心力合力的方向.此外將壁溫也進(jìn)行了無量綱化，無量綱的溫度形式為

(7)

以下將就Θ與Nuθ/Nu的周向分布隨φ和螺旋管幾何結(jié)構(gòu)參數(shù)的變化情況分別進(jìn)行討論.

3.1 φ對截面周向傳熱分布的影響

首先研究了φ變化較小的情況下周向傳熱分布的變化情況，選取工況1管道沿程的4個(gè)截面進(jìn)行比較.工況1螺旋管截面?zhèn)鳠岱植记闆r如圖6所示.由圖6可知，由于這4個(gè)截面沿流動(dòng)方向的流體溫度遞增、密度遞減、黏度遞減，質(zhì)量流速恒定時(shí)u和Re均遞增，φ由12.66遞減至11.79.此時(shí)，Θ和h在周向上的變化曲線形狀基本不變：螺旋管內(nèi)壁溫從θ=0°開始緩慢下降，在90°左右下降到最低，之后又開始緩慢上升，θ=225°之后溫度開始急劇上升，在θ=285°～300°之間溫度達(dá)到最高.總體而言，低溫區(qū)溫度變化比較平緩，高溫區(qū)溫度變化劇烈；與Θ變化趨勢相反，h在高換熱區(qū)域，即θ=45°～135°區(qū)域變化較為平緩，在θ=80°～100°之間達(dá)到最大，換熱性能最強(qiáng)，之后又開始緩慢下降，θ=225°之后h急劇下降，在θ=290°～300°之間達(dá)到最低.隨著Re的上升，離心力略有增大，重力和離心力的合力方向發(fā)生改變，使得最高壁溫位置由θ=295°向θ=285°方向作小角度的移動(dòng)，整體而言變化十分微小.因此可以認(rèn)為，φ變化較小時(shí)截面的周向傳熱分布特性基本不變.此外，Re越大，盡管h的最低值基本不變，但h的最高值及平均值均有所上升，即Re越大，整體換熱性能越強(qiáng).

圖6 φ微小變化時(shí)的無量綱溫度分布和換熱系數(shù)分布Fig.6 Dimensionless temperature distribution and heat transfer coefficient distribution with small changes in φ

φ變化較大的情況下，周向傳熱分布的變化情況如圖7所示.分別選取工況1、2、3、4的某橫截面進(jìn)行研究.這4個(gè)工況的物理模型均為螺旋管1，但因流速差異較大，其φ值有明顯的差異.不同φ對應(yīng)的最高壁溫點(diǎn)位置等分布特征如表4所示.綜合表4和圖7可以看出，無重力作用時(shí)，最高壁溫點(diǎn)在θ=270° 左右即最內(nèi)側(cè)， 最低壁溫點(diǎn)在θ=90° 左右即最外側(cè).此時(shí)，螺旋管內(nèi)流體在截面上主要受到離心力和壓強(qiáng)梯度力的作用.離心力指向管外側(cè)，壓強(qiáng)梯度力指向管內(nèi)側(cè).離心力的大小與流體軸向速度的平方呈正比，在螺旋管截面中心區(qū)的流體流速較高，且離心力大于壓強(qiáng)梯度力，管內(nèi)側(cè)壁面處的高溫流體經(jīng)過管道中心區(qū)向管外側(cè)流動(dòng)，并與主流的低溫流體混合后降溫，而貼近管壁處的流體由于黏性力的影響，軸向速度變小，離心力小于壓強(qiáng)梯度力，到達(dá)外側(cè)壁面的低溫流體在壓強(qiáng)梯度力驅(qū)動(dòng)下分兩支沿兩側(cè)管壁流回管內(nèi)側(cè)并由兩側(cè)管壁加熱后升溫，如此循環(huán)往復(fù)形成了截面上下兩個(gè)Dean渦，即圖8中的二次流，也造成了θ=90°位置流體溫度及壁溫低，而θ=270°位置流體溫度及壁溫高的現(xiàn)象[11-13].而當(dāng)重力相對作用逐漸增強(qiáng)時(shí)，由圖7及表4可知，最高壁溫點(diǎn)由內(nèi)側(cè)開始向頂部移動(dòng).當(dāng)φ=63.96時(shí)，最高壁溫點(diǎn)甚至移動(dòng)到了θ=338°處，同時(shí)高換熱區(qū)域也逐漸從外側(cè)向管底部移動(dòng)，并且該區(qū)域逐漸變寬，溫度變化更加緩和.同時(shí)，Nuθ/Nu變化幅度也逐漸變小，這是因?yàn)棣赵龃髸r(shí)，螺旋管內(nèi)流速降低、離心力變小、二次流強(qiáng)度變?nèi)?，截面周向傳熱不均勻性程度降?

圖7 不同φ下橫截面上的無量綱溫度總體分布和Nuθ/Nu分布Fig.7 Dimensionless temperature distribution and Nuθ/Nu distribution on cross sections at different φ values

圖8 無重力作用下的螺旋管二次流Fig.8 Secondary flow of helical pipes without gravity

圖7對應(yīng)曲線的詳細(xì)溫度分布及二次流情況如圖9所示，其中Tbl為局部Tb.圖9(a)中φ=63.96，即重力占主導(dǎo)作用，這時(shí)截面中心的流體從θ=340°方向向大約θ=150°方向流動(dòng)，θ=150°一側(cè)流速最大、溫度最低，而θ=340°方向溫度最高，整個(gè)溫度分布呈現(xiàn)出明顯的分層.從圖9(a)～(d)離心力作用逐漸增強(qiáng)，截面中心二次流方向逐漸由傾斜變?yōu)樗剑瑴囟确植嫉姆謱忧闆r也和二次流流動(dòng)狀況相符合，與圖7的規(guī)律一致.此外由表4可知，當(dāng)有重力作用時(shí)，壁溫最高點(diǎn)與壁溫最低點(diǎn)的角度差在θ=185°～220°之間，超過θ=180°.結(jié)合圖9的二次流流動(dòng)情況可以發(fā)現(xiàn)，重力對內(nèi)外側(cè)流體二次流的影響不同，圖9(d)中僅有離心力作用，截面中心的二次流水平流動(dòng)，因而最高壁溫和最低壁溫位置相差接近180°；圖(b)和(c)中重力作用增強(qiáng)，重力的浮升力效應(yīng)使得內(nèi)側(cè)的熱流體上浮，該側(cè)的中心水平二次流變成了由θ=300°～320° 位置斜向下的流動(dòng)，而外側(cè)流體的流動(dòng)方向也由水平方向略微向下偏移，但由于內(nèi)側(cè)流體二次流流動(dòng)方向在重力作用下偏移更大，導(dǎo)致了最高壁溫位置與最低壁溫位置的角度差超過180°，這反映出重力對內(nèi)外側(cè)二次流不同的影響.

表4 圖7中不同φ下的周向傳熱特征

圖9 不同φ下橫截面上的溫度分布及二次流情況Fig.9 Temperature distribution and secondary flow on cross section at different φ values

3.2 螺旋管結(jié)構(gòu)參數(shù)對截面周向傳熱分布的影響

為探究a對截面溫度分布的影響，在工況7和工況8的計(jì)算結(jié)果中分別選取了φ=22.25，20.22兩個(gè)截面進(jìn)行比較.由表2可以看出，工況7和工況8的物理模型除a不同外，其他參數(shù)均相同.不同水力學(xué)直徑下橫截面的溫度及Nuθ/Nu的分布如圖10所示.由圖10可知，a對螺旋管周向傳熱分布沒有明顯的影響，相同φ下工況7與工況8的Nuθ/Nu分布曲線基本一致.

圖10 不同水力學(xué)直徑下橫截面的溫度及Nuθ/Nu分布Fig.10 Temperature and Nuθ/Nu distribution of cross sections at different hydraulic diameters

為探究c對截面溫度分布的影響，在工況5和工況7的計(jì)算結(jié)果中分別選取了φ=20.23，19.27兩個(gè)截面進(jìn)行比較.不同螺旋直徑下截面溫度及Nuθ/Nu的分布如圖11所示.由表2可以看出工況5和工況7的物理模型除c不同外，其他均相同.圖11的結(jié)果表明，當(dāng)φ一定時(shí)，c變大，周向壁溫峰值和谷值位置不變，而Nuθ/Nu波動(dòng)幅度明顯變小.結(jié)合上文所述的二次流的形成機(jī)理對此現(xiàn)象解釋如下：當(dāng)φ不變時(shí)，重力和離心力的合力方向不變，二次流動(dòng)方向不變，所以周向壁溫峰值和谷值位置不變；當(dāng)φ不變時(shí)，離心力大小不變，螺旋管截面中心區(qū)的二次流主要在離心力作用下產(chǎn)生，因此強(qiáng)度變化不大；當(dāng)c增大即螺旋管曲率減小時(shí)，截面上壓強(qiáng)梯度力變小，貼近壁面處的二次流由壓強(qiáng)梯度力克服離心力產(chǎn)生，這部分二次流強(qiáng)度明顯減弱，導(dǎo)致了二次流強(qiáng)度整體變?nèi)?，截面周向傳熱不均勻程度減小，Nuθ/Nu波動(dòng)幅度減小.

圖11 不同螺旋直徑下截面溫度及Nuθ/Nu分布Fig.11 Temperature and Nuθ/Nu distribution of cross sections at different hydraulic diameters

為探究α對截面溫度分布的影響，在工況7和工況6的計(jì)算結(jié)果中分別選取了φ=23.53，17.5的兩個(gè)截面進(jìn)行了比較.從表2可以看出工況7和工況6的物理模型除α不同外其他均相同.不同螺旋升角下橫截面的溫度及Nuθ/Nu分布如圖12所示.圖12結(jié)果表明，α對螺旋管周向傳熱并沒有明顯的影響，相同φ下工況7與工況6的Nuθ/Nu分布曲線基本一致.螺旋管受力分析如圖13所示，其中F為向心力與重力的合力.由圖13可知，重力與橫截面1-1的夾角為α，故重力加速度在橫截面1-1上的分量為gcosα，對于螺旋線運(yùn)動(dòng)，離心力方向?yàn)樗奖畴x螺旋線中心軸方向，又知離心力方向與螺旋管軸線切線和重力方向垂直，可以得出離心力與重力在橫截面1-1上的分量垂直，故離心力方向如圖12中的橫截面1-1所示，其加速度大小為2u2/c.因此可以得出F與離心力的夾角β滿足：

圖12 不同螺旋升角下橫截面的溫度及Nuθ/Nu分布Fig.12 Temperature and Nuθ/Nu distribution of cross sections at different spirally ascend angles

圖13 螺旋管受力分析Fig.13 Force diagram of helical pipes

(8)

工況7和工況6的α分別為12°和3.5°，而cos12°/cos3.5°=0.98，即12°與3.5°升角下β差異極小，即合力F方向變化極小，故α(在本文的升角范圍內(nèi))對螺旋管傳熱在周向的分布影響不大.

4 結(jié)論

本文采用CFD軟件CFX 16.0對螺旋管單相流動(dòng)換熱展開了數(shù)值模擬計(jì)算，并基于計(jì)算結(jié)果對螺旋管內(nèi)單相流動(dòng)周向非均勻傳熱現(xiàn)象做了分析，分析結(jié)論如下：

(1) 重力和離心力的共同作用決定了螺旋管周向傳熱分布，影響截面周向傳熱分布的主要因素是加速度之比φ.φ越大，離心力作用減弱，重力作用相對增強(qiáng)；最高壁溫點(diǎn)由管內(nèi)側(cè)向管頂部方向移動(dòng)，低溫高換熱區(qū)域也逐漸由外側(cè)向管底部移動(dòng)，并且該區(qū)域逐漸變寬；Θ與Nuθ/Nu在周向上的變化幅度也逐漸變小.

(2) 由于重力對橫截面內(nèi)外側(cè)二次流不同程度的扭曲，壁面周向最高內(nèi)壁溫位置與最低內(nèi)壁溫位置相距超過180°，位于185°～220°之間.

(3)a和α(本文的升角范圍內(nèi))對螺旋管周向的傳熱分布影響不明顯.c增大時(shí)，螺旋管周向壁溫峰值和谷值的位置基本不變，Nuθ/Nu波動(dòng)幅度變小.