張鈞波，張 敏

(1.南京師范大學(xué)泰州學(xué)院，江蘇泰州225300)

(2.南京理工大學(xué)能源與動(dòng)力工程學(xué)院，江蘇南京210094)

非牛頓流體在偏心圓環(huán)管中的流動(dòng)和傳熱在石油化工、生物、機(jī)械、食品等領(lǐng)域中有著廣泛的應(yīng)用，但因非牛頓流體剪切應(yīng)力與剪切速率的非線性關(guān)系，導(dǎo)致其數(shù)值計(jì)算比牛頓流體更加困難.近年來國(guó)內(nèi)外學(xué)者對(duì)非牛頓流體在圓環(huán)管道中流動(dòng)和傳熱或牛頓流體在偏心圓環(huán)中的流動(dòng)和傳熱進(jìn)行了很多研究［1－7］，但是對(duì)非牛頓流體在偏心圓環(huán)中的流動(dòng)和傳熱考慮的卻很少.文中在數(shù)值計(jì)算過程中采用基元中心法對(duì)非牛頓流體的非線性粘性系數(shù)進(jìn)行離散，解決了冪律流體模型迭代收斂困難的問題，并在內(nèi)壁恒溫，外壁絕熱的邊界條件下，對(duì)非牛頓流體在偏心圓環(huán)管中充分發(fā)展層流的流動(dòng)和傳熱進(jìn)行數(shù)值計(jì)算.分析了偏心圓環(huán)管的偏心率和非牛頓流體的冪律因子對(duì)圓環(huán)通道內(nèi)速度和溫度分布的影響，具有實(shí)際的理論意義和工程應(yīng)用價(jià)值.

1 控制方程和邊界條件

1.1 流動(dòng)控制方程和邊界條件

偏心圓環(huán)管橫截面如圖1，圖中ri和ro分別為內(nèi)外圓半徑，ε為偏心距，即內(nèi)外圓圓心距.文中以偏心圓環(huán)管內(nèi)充分發(fā)展層流為研究對(duì)象.

圖1 計(jì)算區(qū)域幾何圖形Fig.1 Geometric drawing of computational domain

對(duì)于充分發(fā)展層流問題，沿圓環(huán)管軸向的流動(dòng)控制方程為:

式中:w為流體軸向的流動(dòng)速度;dp/dz為圓環(huán)管軸向的壓力降;μ為流體動(dòng)力粘度;考慮偏心圓環(huán)管橫截面的對(duì)稱性，取偏心圓環(huán)管橫截面的右半部分為數(shù)值計(jì)算幾何模型并生成網(wǎng)格，計(jì)算網(wǎng)格如圖2.

圖2 計(jì)算區(qū)域網(wǎng)格Fig.2 Grid of computation domain

邊界條件為:

為簡(jiǎn)化計(jì)算和增加程序的通用性，將各計(jì)算參數(shù)無量綱化，取外圓半徑ro為長(zhǎng)度量綱，則無量綱坐標(biāo)(x*，y*)為:

偏心率和半徑比分別為:

定義無量綱速度w*:

式中:K和n是流體的無量綱常數(shù).將各參數(shù)的無量綱形式代入控制方程中:

同樣可得邊界條件:

1.2 傳熱控制方程和邊界條件

對(duì)于充分發(fā)展層流的穩(wěn)態(tài)導(dǎo)熱問題，溫度只與流體的軸向速度有關(guān)，因此傳熱控制方程可表示為:

式中:ρ為密度;Cp為定壓比熱;k為熱傳導(dǎo)系數(shù);ST為外熱源或熱匯為溫度梯度.

式中:Tw為壁面溫度;Tb為平均體溫度;Tm為平均溫度.

假設(shè)無外熱源，則控制方程(12)可簡(jiǎn)化為

取偏心圓環(huán)管的傳熱邊界條件為內(nèi)壁恒壁溫外壁絕熱，即

控制方程的無量綱式為:

wm為平均速度，無量綱坐標(biāo)為:

Dh為水力直徑，無量綱溫度為:

2 數(shù)值計(jì)算方法

對(duì)于流動(dòng)控制方程(8)和傳熱控制方程(19)的離散，采用二階精度的有限體積法，與文獻(xiàn)［8－9］的計(jì)算方法相同.以下重點(diǎn)介紹采用基元中心法對(duì)冪律非牛頓流體表觀粘性系數(shù)的離散.粘性系數(shù)可用下式表示:式中:K是流體物性系數(shù);n為冪律因子;γ為應(yīng)變速度.對(duì)于充分發(fā)展層流，w為流體軸向流動(dòng)速度，運(yùn)用基元中心法對(duì) ?w/?x和 ?w/?y進(jìn)行離散，假設(shè)兩控制體中心點(diǎn)連線上速度梯度呈線性變化關(guān)系，則此連線中點(diǎn)處的速度梯度值即為相鄰控制體交界面處的速度梯度［10］.

式中:i為相鄰控制體編號(hào)

則速度梯度在x和y方向上的分量分別為:

式中:θi為控制體中心連線與x軸正向的夾角.

采用基元中心法，求得:

m為相鄰控制體個(gè)數(shù)，得到粘性系數(shù)的離散式為:

3 計(jì)算結(jié)果及分析

文中選用的非牛頓流體模型剪切應(yīng)變力與剪切應(yīng)變存在冪函數(shù)關(guān)系，即Ostwald-de Waele關(guān)系式:

式中:K是物性系數(shù)，n為冪律因子.選取n=0.6，0.8，1.0的剪切稀化流體，在半徑比r*=0.5，偏心率分別為ε*=0.2，0.5，0.8的偏心圓環(huán)管中進(jìn)行數(shù)值計(jì)算分析.

3.1 偏心圓環(huán)中非牛頓流體流動(dòng)數(shù)值計(jì)算

圖3為n=0.6的剪切稀化流在ε*=0.2，0.5，0.8，r*=0.5的圓環(huán)管中的等速線分布圖.從圖中對(duì)比可以看出偏心率造成流體速度分布不均，且隨著偏心率的逐漸增大，位于流道截面較窄區(qū)域的流體速度逐漸減少為零，而流道截面較寬區(qū)域的流體流動(dòng)速度逐漸增大.

圖3 冪律因子n=0.6的非牛頓流體在不同偏心圓環(huán)管中流動(dòng)等速線分布圖Fig.3 Distribution of isovels in different eccentric annular channels(n=0.6)

圖4為不同冪律因子的非牛頓流體在r*=0.5，ε*=0.8的偏心圓環(huán)管較寬區(qū)域?qū)ΨQ軸上的速度分布圖.圖中橫坐標(biāo)通過公式(33)由偏心圓環(huán)管較寬區(qū)域?qū)ΨQ軸上的縱坐標(biāo)轉(zhuǎn)換得到的.圖中流道較窄區(qū)域(窄域)的速度趨近為零.而流體速度在較寬區(qū)域(寬域)的二分之一半徑處達(dá)到最大值，且速度分布曲線的斜率隨著冪律因子減小而減小，這也反應(yīng)了非牛頓流體粘性系數(shù)隨著速度梯度的增大而減小.

圖5為ε*=0.2的偏心圓環(huán)管窄域和寬域?qū)ΨQ軸線上的速度分布曲線.從圖中可以看出不同冪律因子的非牛頓流體在流體通道窄域與寬域的速度分布曲線相似，最大速度都是隨著冪律因子n減小而減小，速度分布斜率也是隨著冪律因子的減小而減小.分析比較圖4，5可得，隨著偏心率的減小，流體通道寬域的速度逐漸減小，而窄域的速度則逐漸增大，最終會(huì)回歸均勻分布，原因是隨著偏心率逐漸減小，流體逐漸由寬域回流到窄域.

圖4 不同冪律因子非牛頓流體在偏心圓環(huán)管對(duì)稱軸上的速度分布曲線Fig.4 Distribution of flow velocity along the symmetry axis

圖5 ε*=0.2時(shí)，窄域和寬域?qū)ΨQ軸線上的速度分布曲線Fig.5 Distribution of flow velocity along the symmetry axis in narrow area and wide area

3.2 非牛頓流體傳熱數(shù)值計(jì)算

圖6為n=0.6的非牛頓流體在偏心率ε*=0.2，0.5，0.8，半徑比 r*=0.5 的圓環(huán)通道中的等溫線分布圖.從圖中可以看出隨著偏心率的增大，偏心圓環(huán)管流道橫截面寬域的溫度逐漸降低，且等溫線分布曲線相對(duì)稀疏.

圖7為n=0.6，0.8，1.0的冪律非牛頓流體在r*=0.5，ε*=0，0.2，0.5，0.8 的 4 個(gè)圓環(huán)管中的溫度分布比較圖.圖中橫軸負(fù)值區(qū)域?yàn)槠膱A環(huán)管橫截面窄域，正值區(qū)域?yàn)槠膱A環(huán)管寬域.從圖中分析可得，同心圓環(huán)管(ε*=0)的溫度分布均勻，而偏心圓環(huán)管中，隨著偏心率的增大流體通道窄域的溫度趨近于零，流體通道寬域的溫度也逐漸降低，溫度梯度也在逐漸增大，原因是偏心率的增大使得流體速度分布不均勻，從而使流體傳熱性能降低.

圖6 冪律因子n=0.6的非牛頓流體在偏心圓環(huán)管中等溫線分布Fig.6 Distribution of isotherm in different eccentric annular channels(n=0.6)

圖7 非牛頓流體在不同偏心率的圓環(huán)管中的溫度分布比較Fig.7 Distribution of temperature along the symmetry axis in different eccentric annular channels(n=0.6)

對(duì)于冪律因子對(duì)傳熱的影響，在同心圓環(huán)管中，流場(chǎng)溫度分布幾乎不受流體冪律因子影響，而在偏心圓環(huán)管中，冪律因子對(duì)流場(chǎng)的溫度分布產(chǎn)生比較大的影響，這種影響在偏心率大于0.8時(shí)，幾乎可以忽略，而在偏心率小于0.8時(shí)，流場(chǎng)的溫度梯度變化才隨冪律因子的增大逐漸顯現(xiàn)出來.

表1為不同偏心圓環(huán)管中的無量綱努塞爾特?cái)?shù)(Nu).從表中數(shù)據(jù)分析可得，同心圓環(huán)管中冪律因子的變化對(duì)努塞爾特?cái)?shù)幾乎不產(chǎn)生影響.而在半徑比和偏心率較大的偏心圓環(huán)管中，冪律因子對(duì)努塞爾特?cái)?shù)(Nu)的影響比較明顯，Nu隨n逐漸增大，這也驗(yàn)證了流體冪律因子對(duì)溫度分布的影響受到偏心率和半徑比的影響.從表中數(shù)據(jù)整體分析，努塞爾特?cái)?shù)Nu隨著偏心率和半徑比的增大而減小，這說明較大的偏心率和半徑比使流體在偏心圓環(huán)管中的傳熱性能降低.

表1 非牛頓流體在不同偏心圓環(huán)管中無量綱努塞爾特?cái)?shù)(Nu)值Table 1 Dimensionless nusselt number in different eccentric annular channels

4 結(jié)論

對(duì)偏心圓環(huán)管中冪律非牛頓流體的流動(dòng)和傳熱進(jìn)行了數(shù)值計(jì)算.得到以下結(jié)論:

1)偏心率是造成偏心圓環(huán)管內(nèi)流體速度分布不均的主要原因，偏心率的增大導(dǎo)致寬域和窄域的速度分布極不對(duì)稱，寬域速度增大，窄域速度趨于停滯.

2)非牛頓流體在流道橫截面較窄區(qū)域與較寬區(qū)域的速度分布曲線相似，最大流速和速度分布斜率隨冪律因子的減小而減小.

3)在偏心圓環(huán)管中，偏心率的增大使流場(chǎng)溫度逐漸降低，溫度梯度和努塞爾特?cái)?shù)Nu也逐漸減小.說明偏心率對(duì)偏心圓環(huán)管中的流體傳熱性能有著比較大的不利影響.

4)冪律因子對(duì)偏心圓環(huán)管內(nèi)溫度分布的影響受到偏心率的影響，且這種影響隨著偏心率的增大先增強(qiáng)后減弱，這個(gè)偏心率臨界值處于0.8附近.

References)

［1］ 袁世偉，賴煥新.冪律非牛頓流體的有限體積算法［J］.華東理工大學(xué)學(xué)報(bào):自然科學(xué)版，2013，39(3):364－369.Yuan Shiwei，Lai Huanxin.A finite volume method for calculating flows of power-law non-newton fluids［J］.Journal of East China University of Science and Technology:Natural Science Edition，2013，39(3):364 －369.(in Chinese)

［2］ 尹析明.環(huán)形套管內(nèi)非牛頓流體非定常旋轉(zhuǎn)流動(dòng)譜方法研究［J］.四川大學(xué)學(xué)報(bào):自然科學(xué)版，2009，46(4):1016－1020.Ying Ximing.Research on analytical approach to timedependent rotation flow of non-Newtonian fluid in the pap between cylinders［J］.Journal of Sichuan University:Natural Science Edition，2009，46(4):1016 －1020.(in Chinese)

［3］ 沈丹丹，朱宏武，丁礦，等.水下臥式采油樹內(nèi)部環(huán)空通道的流動(dòng)特性分析［J］.石油機(jī)械，2014，42(2):42－44.Shen Dandan，Zhu Hongwu，Ding Kuang，et al.Analysis of the flow characteristic of the annular channel in sub sea horizontal christmas tree［J］.China Petroleum Machinery，2014，42(2):42 －44.(in Chinese)

［4］ 姚壽廣，程清芳，王公利，等.以泡沫金屬為吸液芯的納米流體熱管傳熱性能試驗(yàn)研究［J］.江蘇科技大學(xué)學(xué)報(bào):自然科學(xué)版，2013，27(6):556 －560.Yao Shouguang，Cheng Qingfang，Wang Gongli，et al.Experimental investigation on heat transfer performance of foam metal wick heat pipe［J］.Journal of Jiangsu U-niversity of Science and Technology:Natural Science E-dition，2013，27(6):556 －560.(in Chinese)

［5］ Nobari M R H，Mehrabani M T.A numerical study of fluid flow and heat transfer in eccentric curved annuli［J］.International Journal of Thermal Sciences，2010，49:380－396.

［6］ Bartosik B.Application of rheological models in prediction of turbulent slurry flow［J］.Flow Turbulence Combust，2010，84:277 －279.

［7］ Ho C J，Lin Y H，Chen T C.A numerical study of natural convection in concentric and eccentric horizontal cylindrical annuli with mixed boundary conditions［J］.International Journal of Heat Fluid Flow，1989，10(1):40－47.

［8］ Patankar S V.Numerical and heat transfer fluid flow［M］.Washington D.C:Hemisphere，1980.

［9］ Zhang M.Modeling of radiative heat transfer and diffusion processes using unstructured grid［D］.USA:Tennessee Technological University，2000.

［10］ Fang P P，Manglik R M，Jog M A.Characteristics of laminar viscous shear-thinning fluid flows in eccentric annular channels［J］.Journal of Non-Newtonian Fluid Mechanics，1999，84:1 －17.