宋佳峻，王雨晴，王鑫垚，徐琴琴，周丹，銀建中

大連理工大學(xué) 化工機(jī)械與安全學(xué)院，遼寧 大連 116024

天然氣作為一種清潔能源，其全球貿(mào)易量與日俱增。為了便于儲存運(yùn)輸，常將天然氣凈化處理后常壓冷卻至-162 ℃，轉(zhuǎn)變?yōu)橐簯B(tài)天然氣(liquid natural gas, LNG)。之后，LNG接收站將LNG加壓到外輸壓力，送至氣化器進(jìn)行氣化。氣化過程中，管內(nèi)的高壓天然氣在升溫時(shí)發(fā)生了跨臨界轉(zhuǎn)變，其傳熱過程屬于超臨界流體傳熱過程。

超臨界流體(supercritical fluid, SCF)是指溫度、壓力均高于臨界值(Tc,Pc)時(shí)物質(zhì)的一種狀態(tài)，但目前也將壓縮流體(即P>Pc,T

Schmidt等通過研究發(fā)現(xiàn)近擬臨界點(diǎn)附近SCF的自然對流強(qiáng)度非常大。Pioro等[1]則將超臨界水用于蒸汽發(fā)生器，以提高火電廠的效率。Yu等[2]通過實(shí)驗(yàn)，研究了管內(nèi)超臨界水傳熱過程中浮升力效應(yīng)的影響，并分析了質(zhì)量流量、熱流密度、壓力等操作條件的影響。Gu等[3]對超臨界深冷甲烷在毛細(xì)管內(nèi)的流動傳熱進(jìn)行實(shí)驗(yàn)研究，同樣探討了質(zhì)量流量、熱流密度、壓力等操作條件對傳熱的影響。粘權(quán)鑫[4]通過對換熱管進(jìn)行等距分段的方法，得到了沿程傳熱系數(shù)的變化。韓冬艷[5]通過理論計(jì)算，采用溫度分段的方法比較了現(xiàn)有描述管內(nèi)傳熱的三類9種模型，認(rèn)為Petukhov模型更適用于超臨界天然氣在管內(nèi)的流動傳熱。靳書武[6]利用標(biāo)準(zhǔn)k-ε模型和定壁溫條件，對水平圓管內(nèi)的超臨界甲烷流動傳熱進(jìn)行研究，得到了壓力、水浴溫度等對管內(nèi)表面?zhèn)鳠嵯禂?shù)影響。張康、董文平等[7-8]利用RNGk-ε模型和定熱流密度條件，分析了超臨界天然氣的傳熱特性。

本文以浸沒燃燒式氣化器蛇形盤管內(nèi)超臨界天然氣的流動傳熱為背景，利用Fluent 18.0進(jìn)行數(shù)值模擬分析，分析超臨界天然氣在定外壁溫?fù)Q熱管內(nèi)的流動傳熱特性，研究不同操作條件對流動傳熱的影響，并與現(xiàn)有超臨界流體管內(nèi)傳熱的關(guān)聯(lián)式進(jìn)行比較。

1 模型的建立及求解

1.1 問題描述與物性計(jì)算

以某在役浸沒燃燒式氣化器為模型，對其管內(nèi)超臨界天然氣的流動傳熱進(jìn)行模擬研究。其相關(guān)操作參數(shù)如表1所示。

表1 浸沒燃燒式氣化器操作參數(shù)

為了對超臨界天然氣在蛇形盤管內(nèi)的流動傳熱進(jìn)行模擬研究，做出以下假設(shè)：

1)由于LNG內(nèi)甲烷含量占85%以上，為簡化模型，管內(nèi)流體按甲烷單質(zhì)處理；

2)在浸沒燃燒式氣化器工作的過程中，其殼程水浴環(huán)境由于受燃燒煙氣的擾動而使各個(gè)位置的溫度分布均勻，因此設(shè)定其蛇形盤管的外壁溫度恒定；

3)流體在管內(nèi)流動時(shí)沿程壓降遠(yuǎn)小于操作壓力，故管內(nèi)按恒定壓力處理；

4)模擬的流動過程已達(dá)到穩(wěn)定工況。

由于超臨界甲烷在升溫過程中物性變化劇烈，故采用美國國家標(biāo)準(zhǔn)與技術(shù)研究院(National Institute of Standards and Technology, NIST)開發(fā)的軟件REFPROP查詢8.88 MPa下甲烷的相關(guān)物性，再利用Fluent自帶的線性插值功能定義物性，插值結(jié)果如圖1所示。

(a) 比熱

(b) 密度

(d) 黏度圖1 8.88 MPa下甲烷的物性與插值結(jié)果

1.2 物理模型與網(wǎng)格劃分

圖2所示為該浸沒燃燒式氣化器的換熱管結(jié)構(gòu)尺寸，其中換熱管外徑D為25.4 mm，壁厚2 mm，兩端長直管段長L1為7 600 mm，中間短直管段長L2為7 200 mm，大彎管半徑R1為114 mm，小彎管半徑R2為56 mm。

圖2 模型尺寸

邊界條件的設(shè)置上，選用質(zhì)量入口、壓力出口為換熱管的進(jìn)出口條件，設(shè)定進(jìn)口質(zhì)量流量和溫度，出口壓力等；換熱管外壁選用定壁溫條件，設(shè)定為水浴溫度；內(nèi)壁設(shè)為耦合壁面。本文利用Gambit 2.4.6對計(jì)算模型進(jìn)行網(wǎng)格劃分，網(wǎng)格全部采用結(jié)構(gòu)化網(wǎng)格，同時(shí)為了更好地模擬彎頭處復(fù)雜的流動情況，對該位置進(jìn)行局部加密處理。網(wǎng)格結(jié)構(gòu)如圖3。

(a) 徑向網(wǎng)格

(b) 軸向網(wǎng)格圖3 網(wǎng)格劃分示意

由于計(jì)算時(shí)采用了標(biāo)準(zhǔn)壁面函數(shù)，故與壁面相鄰的第1層網(wǎng)格節(jié)點(diǎn)應(yīng)位于湍流的充分發(fā)展區(qū)，在黏性底層內(nèi)不布置任何節(jié)點(diǎn)。在劃分網(wǎng)格時(shí)，根據(jù)量綱為1的數(shù)y+對第1層網(wǎng)格高度進(jìn)行估算，y+定義如式(1)所示，保證計(jì)算完后y+值在30～300之間，且接近30。

(1)

式中：Uτ為壁面摩擦速度，m為壁面切應(yīng)力，MPa；y為第1層網(wǎng)格的高度，m。

1.3 數(shù)學(xué)模型與求解器設(shè)置

計(jì)算過程中，流體的流動傳熱滿足穩(wěn)態(tài)情況下的控制方程，如式(2)～(6)。

連續(xù)性方程：

(2)

式中：xi為i方向；ρ為流體的密度，kgm-3；ui為流體速度在方向上的分量，ms-1。

動量方程：

(3)

式中：p為流體靜壓力，Pa；τij為流體受到的切應(yīng)力，Pa；gi為重力加速度在i方向上的分量，ms-2。

能量方程：

(4)

式中：T為流體溫度，K；λ為流體的熱導(dǎo)率，Wm-1K-1；cp為流體的定壓比熱容，Jkg-1K-1。

湍動能k的控制方程：

(5)

式中：k為流體的湍動能，m2s-2；μ為流體的動力黏度，Pas；μt為流體湍流黏度，Pas，μt=cμρk2/ε，cμ為系數(shù)；σk是湍動能的普朗特?cái)?shù)；ε是湍流耗散率。

湍流耗散率ε的控制方程：

(6)

式中：c1、c2為系數(shù)；σε為湍流耗散率的普朗特?cái)?shù)。不同k-ε模型的區(qū)別主要體現(xiàn)在k和ε方程中的系數(shù)cμ、c1、c2、σk、σε的不同。本文使用對二次流模擬較好的RNGk-ε模型，其系數(shù)為理論分析得出而非實(shí)驗(yàn)數(shù)據(jù)，其取值可見文獻(xiàn)[9]，在此不再贅述。

本文采用SIMPLEC算法。在對對流項(xiàng)的離散上，湍流脈動能方程和湍流耗散率方程采用二階迎風(fēng)格式離散；動量方程和能量方程采用QUICK格式離散。計(jì)算時(shí)對出口溫度和質(zhì)量流量進(jìn)行監(jiān)控，當(dāng)相關(guān)出口參數(shù)穩(wěn)定且殘差達(dá)到收斂標(biāo)準(zhǔn)后，即認(rèn)為計(jì)算收斂。

管程內(nèi)壁表面?zhèn)鳠嵯禂?shù)按照式(7)進(jìn)行計(jì)算：

(7)

式中：hi為內(nèi)壁表面?zhèn)鳠嵯禂?shù)，Wm-2K-1；q為局部壁面熱通量，Wm-2；Tw為截面管壁平均溫度，K；Tb為截面流體平均溫度，K。

1.4 模型驗(yàn)證

為了驗(yàn)證計(jì)算結(jié)果與網(wǎng)格的無關(guān)性，改變軸向和徑向網(wǎng)格尺寸，如表2所示，對計(jì)算結(jié)果進(jìn)行對比，從圖4可以發(fā)現(xiàn)case2與case1計(jì)算結(jié)果基本一致，而case3偏差較大。考慮到計(jì)算的準(zhǔn)確性和計(jì)算效率，本文選用case2進(jìn)行計(jì)算。

表2 網(wǎng)格無關(guān)性驗(yàn)證

圖4 不同Case計(jì)算結(jié)果比較

根據(jù)在役SCV運(yùn)行結(jié)果，對模型進(jìn)行驗(yàn)證。模型模擬出口溫度為15 ℃，滿足其外輸天然氣的溫度要求，證明了本模型采用的數(shù)值方法的準(zhǔn)確性，可利用其進(jìn)行進(jìn)一步研究。

2 結(jié)果與討論

2.1 流動與傳熱分析

圖5所示為流體流動過程中，流體主流溫度和內(nèi)壁溫度、相關(guān)物性、主流平均速度和管內(nèi)截面平均傳熱系數(shù)隨換熱管長的變化規(guī)律?？梢钥闯觯S流動換熱的進(jìn)行，壁溫和主流溫度逐漸趨于水浴溫度，且溫差逐漸變小。在彎頭處由于傳熱強(qiáng)化，使局部內(nèi)壁溫度偏低，但波動幅度逐漸減小。另外，流體主流溫度的變化在擬臨界點(diǎn)位置附近較平緩，偏離該區(qū)域后再次大幅上升。而從速度分布可以看出，由于流體密度隨溫度升高而降低，故管內(nèi)流速不斷增大。

(a) 流體溫度

(b) 物性

(c) 速度

(d) 傳熱系數(shù)圖5 管內(nèi)流體溫度、物性、速度和傳熱系數(shù)沿管長的變化規(guī)律

從圖5中可以看出，在入口段，管內(nèi)傳熱系數(shù)先下降后增大，這與入口段流體流動還未充分發(fā)展有關(guān)，由于入口段熱邊界層較薄，局部表面?zhèn)鳠嵯禂?shù)較充分發(fā)展段偏高，隨流體逐漸發(fā)展充分，在湍流擾動作用下表面?zhèn)鳠嵯禂?shù)沿主流方向先下降后增大。在充分發(fā)展區(qū)，管內(nèi)傳熱系數(shù)隨換熱管長先增大后減小，在擬臨界點(diǎn)位置出現(xiàn)峰值，這與該位置定壓比熱容較大有關(guān)。另外，在彎管位置可以看出，由于流動方向改變造成二次流，使傳熱得到顯著強(qiáng)化，并且彎頭曲率半徑越小，傳熱強(qiáng)化效果越明顯。

選取如圖6所示的換熱管環(huán)向4個(gè)位置，分析換熱管環(huán)向局部傳熱系數(shù)沿管長的變化規(guī)律，結(jié)果如圖7所示?？梢钥闯鱿啾扔诮孛嫫骄鶄鳠嵯禂?shù)，在第一排管束，即劇烈升溫部分的換熱管頂部的傳熱系數(shù)偏低約32%，但彎頭處局部強(qiáng)化效果更明顯；管底部在劇烈升溫段系數(shù)偏大約22%，但在彎頭處卻發(fā)生了傳熱的惡化；管側(cè)邊局部傳熱系數(shù)與截面平均傳熱系數(shù)較為接近。由于流體在水平段流動時(shí)，在浮升力與重力作用下密度較小的熱流體積聚在管的頂部，密度較大的冷流體積聚在管的底部，而流體熱導(dǎo)率隨溫度的升高而降低，造成頂部傳熱熱阻較大，底部熱阻較小，從而使傳熱系數(shù)出現(xiàn)差異。而在彎頭處，由于流動方向改變，原來在管頂部的熱流體在經(jīng)過彎頭后未重新回到管頂部時(shí)積聚在新管段的底部，因此分別造成了頂部底部的傳熱強(qiáng)化與惡化。

圖6 換熱管環(huán)向4個(gè)位置分布

(a) 管內(nèi)截面

(b) 管頂部局部

(c) 管底部局部

(d) 管側(cè)邊局部圖7 換熱管環(huán)向局部傳熱系數(shù)沿管長的變化

2.2 操作參數(shù)的影響

1) 質(zhì)量流量

圖8為質(zhì)量流量對傳熱系數(shù)的影響。從圖8中可以看出，質(zhì)量流量越大，傳熱系數(shù)越高，且流體跨臨界位置相對延后。當(dāng)管內(nèi)流速增大時(shí)，流體的湍流程度增加，有利于削薄邊界層，對傳熱起到強(qiáng)化作用。但是流速增大使流體流過同一控制體時(shí)的時(shí)間減小，雖然傳熱系數(shù)增大，但傳熱不充分，故升溫較慢，跨臨界位置延后。

圖8 質(zhì)量流量對傳熱系數(shù)的影響

2) 入口溫度

如圖9所示為入口溫度對傳熱系數(shù)的影響。從圖9中可以看出，入口溫度增加，傳熱系數(shù)變化不大，但流體跨臨界位置稍稍提前，且第1個(gè)彎頭處傳熱系數(shù)偏小。這是由于入口溫度增加使流體能更快地達(dá)到擬臨界溫度，且第1個(gè)彎頭處溫度偏高(擬臨界溫度之上)，比熱容偏小，傳熱系數(shù)偏小。

圖9 入口溫度對傳熱系數(shù)的影響

3) 水浴溫度

如圖10所示為水浴溫度對管內(nèi)傳熱系數(shù)的影響?？梢钥闯觯囟仍龃笸瑯訒沽黧w跨臨界位置提前，且影響效果強(qiáng)于入口溫度。另外，傳熱系數(shù)峰值隨水浴溫度升高而降低，這是由于水浴溫度越高，相同位置處的熱導(dǎo)率越低，增大了導(dǎo)熱熱阻，使傳熱系數(shù)下降。

圖10 水浴溫度對傳熱系數(shù)的影響

2.3 模擬結(jié)果與傳熱關(guān)聯(lián)式的對比

目前，關(guān)于超臨界流體在管內(nèi)流動的傳熱特性已有諸多關(guān)聯(lián)式對其進(jìn)行描述，但尚不能統(tǒng)一。本文挑選其中較為經(jīng)典的4個(gè)關(guān)聯(lián)式，與模擬結(jié)果進(jìn)行對比。

Dittus與Boelter等[10]提出的管內(nèi)單相流體強(qiáng)制對流傳熱的關(guān)聯(lián)式是目前工程上應(yīng)用最為廣泛的，其公式形式為

Nub=0.023Reb0.8Prbn

(8)

式中：Nu是流體對流傳熱的平均努塞爾數(shù)，Nu=hid/λ；Re為流體的雷諾數(shù)，Re=ρud/μ；Pr為流體的普朗特?cái)?shù)，Pr=μcp/λ；下標(biāo)b指管橫截面流體的平均參數(shù)；上標(biāo)n，當(dāng)流體被加熱時(shí)取0.4，流體被冷卻時(shí)取0.3。

對于超臨界流體在管內(nèi)的混合對流傳熱，Jackson和Hall等[11]提出了新的經(jīng)驗(yàn)關(guān)聯(lián)式：

(9)

式中Tpc指流體在對應(yīng)壓力下的的擬臨界溫度，K。

Jackson和Fewster[12]也提出了對應(yīng)的關(guān)聯(lián)式：

(10)

Petukhov關(guān)聯(lián)式[13]：

(11)

圖11為模擬結(jié)果與式(9)～(11)4個(gè)傳熱關(guān)聯(lián)式計(jì)算結(jié)果的對比。從圖11中可以看出，Jackson& Hall關(guān)聯(lián)式與模擬結(jié)果基本吻合，平均相對誤差約為4.84%，但在流體跨臨界之后，該關(guān)聯(lián)式的計(jì)算結(jié)果偏低約5%；而Dittus& Boelter公式、Jackson& Fewster公式與Petukhov公式與模擬結(jié)果偏差較大，平均相對誤差分別為8.17%、5.03%和10.85%，且均在近臨界區(qū)域的計(jì)算結(jié)果偏高，超臨界區(qū)域結(jié)果偏低。另外，從圖中還可以看出4個(gè)公式對超臨界流體在彎頭處由于二次流引起的傳熱強(qiáng)化作用不能體現(xiàn)或體現(xiàn)不夠，偏差可達(dá)30%以上，但考慮到彎頭部分換熱管長占總管長的比例較小，對總平均傳熱系數(shù)影響較小，故可以對其進(jìn)行簡化處理。

圖11 模擬結(jié)果與已有傳熱關(guān)聯(lián)式的對比

3 結(jié)論

本文通過建立超臨界天然氣在浸沒燃燒式氣化器蛇形盤管中流動的相關(guān)模型，對管內(nèi)超臨界流體流動傳熱過程進(jìn)行了相關(guān)研究，得到以下結(jié)論：

1)換熱管內(nèi)壁溫和流體主流溫度逐漸趨于水浴溫度，但在彎頭處，由于傳熱強(qiáng)化使局部內(nèi)壁溫度偏低，波動幅度隨換熱進(jìn)行逐漸變小；由于流體密度隨溫度升高而降低，使管內(nèi)流速不斷增大；入口段未充分發(fā)展的流體流動使傳熱系數(shù)先下降后增大，而在充分發(fā)展區(qū)，管內(nèi)傳熱系數(shù)隨換熱管長先增大后減小，在擬臨界點(diǎn)位置出現(xiàn)峰值。另外，在彎管處二次流的影響下，使傳熱得到顯著強(qiáng)化，并且彎頭曲率半徑越小，傳熱強(qiáng)化效果越明顯。

2)比較換熱管環(huán)向局部傳熱系數(shù)沿管長的變化發(fā)現(xiàn)，管頂部的局部傳熱系數(shù)偏低，但彎頭處局部強(qiáng)化效果更明顯；管底部局部傳熱系數(shù)偏大，但在彎頭處卻發(fā)生了傳熱的惡化；管側(cè)邊局部傳熱系數(shù)與截面平均傳熱系數(shù)較為接近。

3)提高管內(nèi)質(zhì)量流量能夠顯著強(qiáng)化管內(nèi)傳熱，但流體跨臨界位置相對延后；提高管程入口溫度和水浴溫度均使流體跨臨界位置提前，且較高的水浴溫度會使傳熱系數(shù)峰值偏低。

4)通過比較已有關(guān)聯(lián)式與模擬結(jié)果發(fā)現(xiàn)，Jackson& Hall關(guān)聯(lián)式與模擬結(jié)果基本吻合，平均相對誤差約為4.84%，最為吻合；而Dittus& Boelter公式、Jackson& Fewster公式與Petukhov公式與模擬結(jié)果偏差均較大；各個(gè)關(guān)聯(lián)式對超臨界流體在彎頭處的傳熱強(qiáng)化預(yù)測不夠準(zhǔn)確，偏差最大可達(dá)30%以上。