吳志勇 ，高 陽 ，麻宏強(qiáng) ，劉 洋 ，姜益強(qiáng)

（1.遼寧石油化工大學(xué) 石油天然氣工程學(xué)院，遼寧 撫順 113001；2.蘭州理工大學(xué) 土木工程學(xué)院，甘肅蘭州 730050；3.哈爾濱工業(yè)大學(xué) 市政與環(huán)境工程學(xué)院，黑龍江 哈爾濱 150090）

1 繞管式換熱器研究背景

繞管式換熱器是天然氣液化工藝中的關(guān)鍵設(shè)備，其特點是結(jié)構(gòu)緊湊、耐低溫、抗高壓、適合多種介質(zhì)換熱。這種換熱器是由許多小管徑換熱管以螺旋形狀分層纏繞在換熱器芯筒上而制成的，相鄰兩層換熱管的旋向相反，層間采用隔條保持間距[1,2]，繞管式換熱器內(nèi)部結(jié)構(gòu)見圖1[3]。對于天然氣液化工藝，殼側(cè)烷烴制冷劑首先進(jìn)入管側(cè)進(jìn)行過冷，天然氣與需要過冷的制冷劑在不同的換熱管內(nèi)形成管側(cè)多股流，并一同向上流動。過冷后的管側(cè)制冷劑在換熱器外部經(jīng)節(jié)流后導(dǎo)入到殼側(cè)，然后向下流動并汽化，對管側(cè)天然氣和需要過冷的制冷劑實現(xiàn)制冷。

在繞管式換熱器殼側(cè)研究方面，已公開發(fā)表的有關(guān)文獻(xiàn)比較少見。文獻(xiàn)[1]介紹了繞管式換熱器結(jié)構(gòu)參數(shù)間的耦合關(guān)系，給出了殼側(cè)流通面積的計算方法。Aunan[4]建立了液化天然氣用繞管式換熱器殼側(cè)實驗臺，對乙烷、丙烷在殼側(cè)的壓降與換熱特性進(jìn)行了實驗研究。Moawed[5]以實驗方式研究了定熱流密度下纏繞管曲率比、扭曲比與殼側(cè)換熱系數(shù)之間的關(guān)系，指出當(dāng)纏繞管曲率比、扭曲比變大時殼側(cè)單相換熱系數(shù)將會明顯增大。賈金才[6]建立了繞管式換熱器三維計算模型，數(shù)值研究了換熱管管徑、纏繞角、徑向比、軸向比對殼側(cè)換熱系數(shù)的影響，但所有模擬工況均為氣相流動，與殼側(cè)兩相流動為主的實際工況存在差別。Lu等[7]針對繞管式換熱器殼側(cè)氣態(tài)換熱過程展開了實驗與計算機(jī)仿真研究工作，建立了具有三層繞管的實驗裝置以及CFD仿真計算模型，實驗證實了所用CFD模型具有較高的計算精度。

圖1 繞管式換熱器內(nèi)部結(jié)構(gòu)[3]Fig.1 Internal structure of spiral wound heat exchanger[3]

通過上述文獻(xiàn)可以發(fā)現(xiàn)，繞管式換熱器殼側(cè)研究僅局限于壓降與換熱性能范疇，而對于空泡率這一兩相流動力學(xué)重要參數(shù)尚未涉及?？张萋适菤庖簝上嗔鲃恿W(xué)領(lǐng)域的重要研究內(nèi)容，其重要性在于它與流型識別、相間作用力等理論研究密不可分。目前，換熱器殼側(cè)空泡率研究集中于直列或錯列水平管束[8-14]，尚無關(guān)于螺旋管束的空泡率關(guān)聯(lián)式，因此本文通過數(shù)值模擬方式對現(xiàn)有空泡率關(guān)聯(lián)式進(jìn)行篩選，找出適用于繞管式換熱器殼側(cè)的關(guān)聯(lián)式，為殼側(cè)流型識別研究奠定基礎(chǔ)。

2 繞管式換熱器殼側(cè)模型的建立

2.1 殼側(cè)幾何模型

液化天然氣用繞管式換熱器整體尺寸巨大，其外部直徑約3～5m，高度在10～50m，加之內(nèi)部結(jié)構(gòu)復(fù)雜，所以極不利于開展數(shù)值模擬研究。本文在研究過程中是以文獻(xiàn)[4]中的實驗數(shù)據(jù)作為模擬計算檢驗標(biāo)準(zhǔn)，因此所建殼側(cè)幾何模型與文獻(xiàn)[4]中的實驗設(shè)備相一致。

文獻(xiàn)[4]中的實驗設(shè)備為LNG繞管式換熱器的簡化模型，由三層換熱管以交替的纏繞方向纏繞而成，每層的換熱管并管纏繞數(shù)目由內(nèi)向外依次為3、4、5。由完整的換熱管纏繞而成的中間層及兩側(cè)流道作為實驗研究的核心部位；內(nèi)、外兩側(cè)纏繞層均由半剖換熱管纏繞而成，并作為內(nèi)、外邊壁使用；換熱器簡化模型見圖2。文獻(xiàn)[4]以給定熱流方式研究殼側(cè)換熱特性，其中內(nèi)、外纏繞層及中間纏繞層上部為絕熱壁面，中間纏繞層下部4排換熱管為受熱壁面，通過管內(nèi)電加熱絲實現(xiàn)定熱流加熱。烷烴制冷劑由頂部30個直徑為10mm的分流孔進(jìn)入殼側(cè)，經(jīng)過上方多排纏繞管束緩沖后進(jìn)入到壓降與換熱測試區(qū)。

圖2 換熱器簡化模型Fig.2 Simplified model of SWHE

按照文獻(xiàn)[4]中的換熱器幾何參數(shù)，本文利用UGNX 6.0建立了殼側(cè)模型，其幾何參數(shù)見表1。由于模擬相變換熱需要精細(xì)的邊界層網(wǎng)格，因此在保證計算精度的前提下，殼側(cè)模型網(wǎng)格數(shù)量將達(dá)到1000萬以上，由此將導(dǎo)致模擬計算耗時巨增?？紤]到該模型具有軸對稱性，因此將模型沿軸向切割36°作為模擬研究對象，使得網(wǎng)格數(shù)量降至100萬左右，此時網(wǎng)格數(shù)量與計算精度能夠得以兼顧。切割下來的殼側(cè)36°幾何模型見圖3。

表1 殼側(cè)幾何模型參數(shù)Table 1 Geometric parameters of model

圖3 殼側(cè)36°幾何模型Fig.3 Angel model with 36 degrees

2.2 控制方程

烷烴制冷劑流入殼側(cè)時，換熱管管束被液膜覆蓋，殼側(cè)形成降膜流動，適宜采用兩相流VOF模型進(jìn)行計算。ANSYSFLUENT下的VOF模型控制方程組如下：

連續(xù)方程

動量方程

能量方程

其中，

式中：βg、β1-汽、液相體積分?jǐn)?shù)；ρg、ρ1-汽、液相密度，kg/m3；μg、μ1-汽、液相動力粘度，kg/(m·s)；v→-汽、液相共享速度，m/s；E-比內(nèi)能，J/kg；Γm-連續(xù)方程源項，kg/(m3·s)；Γe-能量方程源項，W/m3。

控制方程組中的連續(xù)方程和能量方程源項通過下列模型進(jìn)行計算：

（1）相變傳質(zhì)模型

根據(jù)Lee[15]的傳質(zhì)模型，當(dāng)流體溫度高于汽化溫度即T≥Tsat時，單位控制單元內(nèi)由液相轉(zhuǎn)變?yōu)槠嗟乃矔r汽化量即為，并有：

反之，當(dāng)T＜Tsat時由汽相轉(zhuǎn)變?yōu)橐合嗟乃矔r冷凝量為：

式中，βm-相變傳質(zhì)過程中的時間松弛參數(shù)，s-1。

（2）能量傳輸模型

當(dāng)汽、液相間傳質(zhì)量Γm確定后，相應(yīng)的能量傳輸源項Γe可通過下式來表達(dá)：

式中，hLH-流體汽化潛熱，J/kg，負(fù)號表示需要外界向控制單元輸入熱量。

另外，VOF模型用于湍流時，本文使用RNG kε湍流應(yīng)力模型[6]作為求解封閉方程。

2.3 邊界條件

對于殼側(cè)36°幾何模型，流體從頂部3個分流孔流入，在底部端面流出，兩側(cè)剖面均近似處理為對稱邊界，中間纏繞層底部4排換熱管（處在流道截面1、2之間）為加熱壁面，其余壁面均為絕熱壁面，邊界條件設(shè)置見表2和圖3。

表2 邊界條件Table 2 Boundary conditions

3 殼側(cè)沸騰時空泡率關(guān)聯(lián)式的篩選

由于氣液兩相流動時空泡率難以測定，因此確定繞管式換熱器殼側(cè)適用的空泡率關(guān)聯(lián)式是通過空泡率數(shù)值模擬結(jié)果與關(guān)聯(lián)式計算結(jié)果相比對的方式來進(jìn)行的。具體來說即是，以反映繞管式換熱器殼側(cè)實際運(yùn)行狀態(tài)的工況作為模擬工況，通過數(shù)值模擬方式獲得流道截面2處的空泡率，然后與備選的空泡率關(guān)聯(lián)式計算結(jié)果作比較，進(jìn)而篩選出繞管式換熱器殼側(cè)適用的空泡率關(guān)聯(lián)式。

本文以繞管式換熱器殼側(cè)常用的乙烷、丙烷單質(zhì)作為制冷劑，模擬工況取自文獻(xiàn)[4]。模擬工況分為兩類，一類用于初步篩選空泡率關(guān)聯(lián)式，另一類涵蓋多壓力、多熱流、多流率情況，用于二次篩選關(guān)聯(lián)式并檢驗其通用性。兩類工況下，殼側(cè)入口干度均為0，其余參數(shù)見表3～4。

為了證實數(shù)值模擬結(jié)果具有準(zhǔn)確性，本文將表4中23個檢驗工況下的換熱系數(shù)模擬結(jié)果和文獻(xiàn)[4]中的實驗結(jié)果進(jìn)行了比較。圖4和圖5分別展示了丙烷和乙烷在各檢驗工況下的換熱系數(shù)偏差，兩種工質(zhì)的換熱偏差均控制在-20%～10%范圍內(nèi)，說明在變壓力、變熱流、變流率情況下繞管式換熱器殼側(cè)數(shù)值模擬結(jié)果具有較高的準(zhǔn)確性。

表3 空泡率關(guān)聯(lián)式篩選工況Table 3 Selection conditions of void correlation

表4 空泡率關(guān)聯(lián)式檢驗工況Table 4 Verification conditions of void correlation

圖4 丙烷檢驗工況下沸騰換熱系數(shù)偏差Fig.4 Propane heat transfer deviation under verification conditions

圖5 乙烷檢驗工況下沸騰換熱系數(shù)偏差Fig.5 Ethane heat transfer deviation under verification conditions

根據(jù)以往的空泡率研究文獻(xiàn)，本文選取11個關(guān)聯(lián)式作為篩選對象，其表達(dá)式見表5（見下頁）。使用關(guān)聯(lián)式時，流道截面2處的干度x可根據(jù)理論汽化量而得到，計算公式如下：

式中，Qin-殼側(cè)輸入熱流量，W；m-殼側(cè)質(zhì)量流率，kg·s-1·m-2；A-殼側(cè)流道平均截面積，m2。

表6（見下頁）顯示的是P1、E1篩選工況下在流道截面2處由備選關(guān)聯(lián)式計算出來的空泡率與數(shù)值模擬出來的空泡率之間的對比結(jié)果，其中的空泡率偏差是以模擬結(jié)果作為基準(zhǔn)。從表中可以看出，關(guān)聯(lián)式VFC1～VFC4計算偏差比較小，因此將VFC1～VFC4關(guān)聯(lián)式作為初選結(jié)果，并且將在丙烷、乙烷檢驗工況下做進(jìn)一步考察。

圖6 丙烷檢驗工況下空泡率關(guān)聯(lián)式比較Fig.6 Comparison between void correlations under propane verification conditions

圖7 乙烷檢驗工況下空泡率關(guān)聯(lián)式比較Fig.7 Comparison between void correlations under ethane verification conditions

圖6和圖7分別顯示出空泡率關(guān)聯(lián)式VFC1～

VFC4在丙烷、乙烷檢驗工況下與模擬結(jié)果之間的比較情況。經(jīng)過殘差分析后發(fā)現(xiàn)由關(guān)聯(lián)式VFC1計算出的空泡率偏差達(dá)到最小，即關(guān)聯(lián)式VFC1優(yōu)于VFC2～ VFC4。

表5 備選的空泡率關(guān)聯(lián)式Table 5 Alternative void correlation

表6 空泡率關(guān)聯(lián)式計算結(jié)果與模擬結(jié)果對比Table 6 Comparison between void correlation and simulation

圖8 丙烷檢驗工況下VFC1計算偏差Fig.8 VFC1 predicting deviation under propane verification conditions

圖9 乙烷檢驗工況下VFC1計算偏差Fig. 9 VFC1 predicting deviation under ethane verification conditions

圖8和圖9給出了Chisholm空泡率關(guān)聯(lián)式（VFC1）在丙烷、乙烷檢驗工況下的計算偏差，相對偏差是基于數(shù)值模擬結(jié)果而得出的。從圖中可以看出，Chisholm空泡率關(guān)聯(lián)式計算偏差可以控制在-15%～0%范圍內(nèi)，并且適用于變熱流、變壓力、變流率情況，因此推薦使用該關(guān)聯(lián)式預(yù)測殼側(cè)空泡率。

4 結(jié)論

通過數(shù)值計算方式對繞管式換熱器殼側(cè)烷烴冷劑沸騰過程進(jìn)行了模擬，得到以下兩點結(jié)論：

（1）VOF模型能夠較好地模擬出繞管式換熱器殼側(cè)冷劑沸騰現(xiàn)象，在變壓力、變熱流、變流率情況下均具有較高的準(zhǔn)確性，其換熱偏差可控制在-20%～10%范圍內(nèi)。

（2）通過數(shù)值模擬方式可以得到繞管式換熱器殼側(cè)沸騰時的空泡率數(shù)據(jù)，并與以往文獻(xiàn)推薦的11個空泡率關(guān)聯(lián)式相比較，通過兩次篩選確定出Chisholm關(guān)聯(lián)式能夠很好地預(yù)測殼側(cè)沸騰時的空泡率，其預(yù)測偏差在-15%～0%范圍內(nèi)。

參考文獻(xiàn)

[1] 吳志勇,陳杰,浦暉,等.LNG繞管式換熱器結(jié)構(gòu)與流通參數(shù)計算方法[J].煤氣與熱力,2014,34(3)：A34-A39.

[2] 吳志勇,陳杰,浦暉,等.LNG繞管式換熱器殼側(cè)過熱流動的數(shù)值模擬[J].煤氣與熱力,2014,34(8)：B6-11.

[3] Steffen H.Near-optimal operation of LNG liquefaction processes by means of regulation [D].Berlin：Berlin Institute of Technology,2011.

[4] Aunan B.Shell-side heat transfer and pressure drop in coil-wound LNG heat exchanger [D]. Trondheim：Norwegian University of Science and Technology,2000.

[5] Moawed M.Experimental study of forced convection from helical coiled tubes with different parameters[J].Energy Convers Manage,2011,52：1150-1156.

[6] 賈金才.幾何參數(shù)對繞管式換熱器傳熱特性影響的數(shù)值研究[J].流體機(jī)械,2011,39(8)：33-37.

[7] Lu X,Du X P,Zeng M,et al.Experimental and numerical investigation on shell-side performance of multilayer spiral-wound heat exchangers[J].Chem Eng Trans,2013,35：445-450.

[8] Schrage D S,Hsu J T,Jensen M K.Two-phase pressure drop in vertical cross?ow across a horizontal tube bundle[J].AIChE J,1988,34：107-115.

[9] Thome.Engineering Data Book Ⅲ [M].Decatur：Wolverine Tube Incorporation,2010.

[10]Fair JR,klip A.Thermal design of horizontal reboilers[J].Chem Eng Prog,1983,79(3)：86-96.

[11]Dowlati R,Kawaji M,Chan A M C.Void fraction and friction pressure drop in two-phase flow across a horizontal tube bundle[C]//AIChE Symposium Series.New York：AIChE,1988,84(263)：126-132.

[12]Dowlati R,Kawaji M,Chan A M C.Pitch-to-diameter effect on two-phase flow across an in-line tube bundle[J].AIChE J,1990,36：765-772.

[13]Dowlati R,Chan A M C,Kawaji M.Hydrodynamics of two-phase flow across horizontal inline and staggered rod bundles[J].JFluids Eng,1992,114：450-456.

[14]Dowlati R,Kawaji M,Chisholm D,et al.Void fraction prediction in two-phase flow across a tube bundle[J].AIChE J,1992,38：619-622.

[15]Lee W H.A Pressure Iteration Scheme for Two-phase Flow Modeling [M].Washington,DC：Hemisphere Publishing,1980：407-431.

[16]Chisholm D.Pressure gradients due to friction during the flow of evaporating two-phase mixtures in smooth tubes and channels[J].Int J Heat Mass Transfer,1973,16(2)：347-358.

[17]Spedding P L,Spence D R.Prediction of holdup in two phase flow[J].Int J Eng Fluid Mechanics,1988,1(1)：67-82.

[18]Ishiraha K,Palen J W,Taborek J.Critical review of correlations for predicting two-phase flow pressure drop across tube banks[J].Heat Transfer Eng,1980,1(3)：23-32.

[19]Bankoff S G.A Variable density single-fluid model for two-phase flow with particular reference to steam-water flow[J].JHeat Transfer,1960,82(4)：265-272.

[20]Zivi SM.Estimation of steady-state steam void-fraction by means of the principle of minimum entropy production[J].JHeat Transfer,1964,86(2)：247-251.

[21]Wallis G B.One Dimensional Two-phase Flow[M].New York：McGraw Hill Incorporation,1969.

[22]Domanski P,Didion D.Computer modeling of the vapor compression cycle with constant flow area expansion device[J].JHigh Energy Phys,2008,26：742-746.

[23]Xu G P,Tso C P,Tou K W.Hydrodynamics of two-phase flow in vertical up and down-flow across a horizontal tube bundle[J].Int JMultiphase Flow,1998,24：1317-1324.