林元載

（東南大學(xué)機(jī)械工程學(xué)院，南京211189）

?

基于流體動(dòng)力學(xué)的低溫液化過(guò)程數(shù)值模擬研究

林元載

（東南大學(xué)機(jī)械工程學(xué)院，南京211189）

摘要：本文通過(guò)計(jì)算流體計(jì)算動(dòng)力學(xué)（CFD）方法同時(shí)模擬換熱器內(nèi)天然氣的主要成分甲烷深冷液化和多相流動(dòng)的過(guò)程，將氣相和液相甲烷看成混合物，通過(guò)混合物模型建立流體的動(dòng)量守恒方程、連續(xù)性方程、體積守恒方程以及相間熱傳遞方程，采用針對(duì)流體自由面運(yùn)動(dòng)的立體體積函數(shù)（VOF）對(duì)相間邊界層進(jìn)行追蹤。

關(guān)鍵詞：換熱器液化相變混合物模型VOF

引言

隨著各國(guó)環(huán)境保護(hù)意識(shí)的增強(qiáng)和世界能源結(jié)構(gòu)在逐漸變化，天然氣成為最受歡迎的能源之一。換熱器中預(yù)冷和深冷液化液化天然氣重要的一個(gè)環(huán)節(jié)。這是一個(gè)涉及能量、熱量和質(zhì)量傳遞的相變過(guò)程以及氣液兩相間界面的追蹤的過(guò)程。1981年，Hirts和Nichols［1］提出了VOF （Volume of Fluid）方法，并使用該方法對(duì)潰壩問題進(jìn)行數(shù)值模擬，論證方法的可行性。該方法的相界面構(gòu)造基本思想為運(yùn)動(dòng)界面追蹤問題的數(shù)值模擬起到了開創(chuàng)性的作用；1998年，Boris Halasz［2］從熱力學(xué)角度，以能量、動(dòng)量、質(zhì)量平衡為基礎(chǔ)，通過(guò)分析蒸發(fā)式冷凝器內(nèi)傳熱傳質(zhì)及流動(dòng)阻力，總結(jié)出了當(dāng)時(shí)所有類型的蒸發(fā)式冷卻裝置通用的數(shù)學(xué)模型。1988年，Osher和Sethian［3］提出Level Set（水平集）方法，較精確計(jì)算相界面曲率及相關(guān)的物理量。2015年，Li S［4］等人通過(guò)引入熱平衡模型，對(duì)亞音速蒸汽注入過(guò)冷池直接接觸冷凝進(jìn)行數(shù)值模擬，結(jié)果表明在管出口軸向溫度隨軸向速度降低而升高和壓力震蕩主要受蒸汽流速、蒸汽冷凝和低溫冷卻水的靜壓力影響。

本文的主要目的對(duì)天然氣在換熱器管程中進(jìn)行深冷液化的問題進(jìn)行研究，分析在低溫環(huán)境下甲烷深冷液化的主要影響因素。本文采用CFX中多相流混合模型，模擬低溫環(huán)境下天然氣深冷液化的相變傳熱傳質(zhì)過(guò)程，分析了不同因素對(duì)液化的影響及液化后的流型。

1 數(shù)值方法

1.1 幾何模型

由于在換熱器中管程通常為彎管結(jié)構(gòu)，本文幾何模型采用“S”型彎管（如圖1所示）。模擬計(jì)算過(guò)程中，入口通入高溫甲烷氣體；氣體與管道壁面換熱，使管道中甲烷遇冷液化；出口處則為甲烷的氣液混合物。

圖1 幾何模型

1.2 控制方程

混合物模型［5］（Mixture）是一種簡(jiǎn)化的多相流模型，它用平均速度的概念來(lái)模擬多相流中各相具有不同速度的情況。該模型能夠求解混合相的能量、動(dòng)量和連續(xù)性方程，以及各相的體積分?jǐn)?shù)、溫度、壓力等物理量。

對(duì)上述模型，通過(guò)建立動(dòng)量守恒方程、連續(xù)性方程、體積守恒方程以及相間熱傳遞方程。再進(jìn)行數(shù)值模型的穩(wěn)態(tài)求解，不考慮控制方程在時(shí)間上的連續(xù)性，控制方程如下：

（1）動(dòng)量守恒方程

（2）連續(xù)性方程

式中，SM描述用戶指定的質(zhì)量源，Γαβ指從β相到α相的單位體積質(zhì)量流量

（3）體積守恒方程

式中，rα表示α相體積比，Np是流體域內(nèi)所有相的數(shù)量。

（4）相間熱傳遞方程

式中，hα,Tα,λα分別表示α相靜焓、溫度、熱導(dǎo)率，SEα是外部熱源，Qα其他相傳遞到α相的內(nèi)能。

1.3 網(wǎng)格劃分

利用Hypermesh劃分模型網(wǎng)格，流體域采用四面體網(wǎng)格劃分，在壁面添加10層邊界層，第一層厚度0.01mm，增長(zhǎng)率1.2。保持壁面邊界層參數(shù)不變，改變網(wǎng)格最小單元尺寸。網(wǎng)格最小單元尺寸：0.5mm、0.7mm、1mm、1.5mm、2mm和2.5mm，網(wǎng)格數(shù)量分別為2605506、1286105、449675、203434、93534和55231。分別對(duì)以上5中網(wǎng)格進(jìn)行模擬計(jì)算，提取管道出口位置溫度的平均值與最大值。根據(jù)分析可以得出從網(wǎng)格數(shù)量449675以后，網(wǎng)格數(shù)量的增加對(duì)模擬結(jié)果的影響很小。本文選定最小單元尺寸0.7、網(wǎng)格數(shù)量1286105作為分析對(duì)象。

1.4 邊界條件

本節(jié)主要針對(duì)流體域動(dòng)力學(xué)參數(shù)和流體域熱力學(xué)參數(shù)進(jìn)行分析。由于液化后的液體受重力影響對(duì)流動(dòng)影響較大，需要考慮重力和浮力的作用；甲烷氣相和液相會(huì)形成交界面，因此選用CFX多相流自由液面模型；流體域入口控制參數(shù)設(shè)置為總壓、溫度及氣態(tài)甲烷的體積比；出口設(shè)置為靜壓；所有壁面均為無(wú)滑移壁面模型。

流體域內(nèi)的熱力學(xué)模型為總能量模型，相間設(shè)置熱傳遞系數(shù)1000W/m2k；管程和殼程之間存在熱量傳遞，為方便計(jì)算設(shè)定殼程流體溫度均為100K，故管外溫度設(shè)置為100K，管的傳熱系數(shù)為20000W/m2k。其入口壓強(qiáng)隨重力方向而改變。

2 數(shù)值結(jié)果與討論

（1）數(shù)值結(jié)果的驗(yàn)證。甲烷在換熱器管程中深冷液化涉及氣液相變及多相流動(dòng)問題。高溫甲烷氣體流經(jīng)壁面，遇冷液化。同時(shí)，甲烷氣體與液化后產(chǎn)生的液體同時(shí)在管道中流動(dòng)。氣體液化主要發(fā)生在管道表面，當(dāng)重力作用可忽略不計(jì)時(shí)，將會(huì)形成環(huán)狀的流行（如圖2所示）。圖3展現(xiàn)的為高溫甲烷氣體從管道入口到管道出口氣液百分比，顏色越深液體比例越高。由圖可以看出，越靠近出口，甲烷液體越多。與Osher S等實(shí)驗(yàn)結(jié)果基本符合，證明了仿真的正確性。

圖2 環(huán)狀流［6］

圖3 冷凝結(jié)果

（2）驅(qū)動(dòng)壓力的影響。保持管道壁面的傳熱系數(shù)不變，設(shè)置7組遞增驅(qū)動(dòng)壓力的計(jì)算組，分別為100Pa、200Pa、300Pa、400Pa、500Pa、600Pa和700Pa。通過(guò)仿真模擬計(jì)算，獲得管道出口處甲烷的冷凝量。圖5為不同驅(qū)動(dòng)壓力時(shí)，出口處甲烷冷凝量的變化曲線。由圖4可以看出，在100～300Pa階段，冷凝量隨著驅(qū)動(dòng)壓力的增大成直線增大；300～500Pa階段，冷凝量成穩(wěn)定狀態(tài)；500～700Pa，冷凝量又回到直線增長(zhǎng)。

圖4 不同驅(qū)動(dòng)壓力時(shí)冷凝量變化曲線

（3）管壁傳熱系數(shù)的影響。保持初始溫度和驅(qū)動(dòng)壓力不變，設(shè)置9組管道壁面?zhèn)鳠嵯禂?shù)遞增的計(jì)算組，分別為1、2、6、8、10、12、12.5、15、17.5和20kW/m2k。通過(guò)模擬計(jì)算，獲得管道出口處甲烷的冷凝量。圖3-4為不同管道壁面?zhèn)鳠嵯禂?shù)時(shí)，出口處甲烷冷凝量的變化曲線。由圖5可以看出，隨著管道壁面?zhèn)鳠嵯禂?shù)的增大，出口處冷凝量呈直線增加。管道壁面?zhèn)鳠嵯禂?shù)增大，管壁的熱阻減少，使氣體與管壁的傳熱效率提高，從而使冷凝的液體增多。

圖5 不同管壁傳熱系數(shù)時(shí)冷凝量變化曲線

3 結(jié)論

本文主要研究驅(qū)動(dòng)壓力和管壁傳熱系數(shù)分別對(duì)目標(biāo)參數(shù)（出口處甲烷冷凝量）的影響。通過(guò)簡(jiǎn)化管道模型，并采用流體動(dòng)力學(xué)混合物模型結(jié)合自由液面模型（VOF）對(duì)其進(jìn)行數(shù)值模擬，得到以下結(jié)論：

（1）數(shù)值模擬結(jié)果與文獻(xiàn)中的流型一致性較好，從而證明了本文所建模型的合理性。

（2）在優(yōu)化的范圍內(nèi)，若要提高甲烷的冷凝量，可通過(guò)提高驅(qū)動(dòng)壓力和使用傳熱系數(shù)較高的管道材料。

參考文獻(xiàn)

［1］Hirt CW,Nichols BD.Volume of Fluid（Vof）Method for the Dynamics of Free Boundaries.J Comput Phys.1981;39（1）:201-25. doi: 10.1016/0021-9991（81）90145-5.

［2］Halasz B.A general mathematical model of evaporative cooling devices. Rev Gen Therm.1998;37（4）:245-55.doi:10.1016/S0035-3159（98）80092-5.

［3］Osher S,Sethian JA.Fronts Propagating with Curvature-Dependent Speed-Algorithms Based on Hamilton-Jacobi Formulations.J Comput Phys.1988;79（1）:12-49. doi:10.1016/0021-9991（88）90002-2.

［4］Li SQ,Wang P,Lu T.Numerical simulation of direct contact condensation of subsonic steam injected in a water pool using VOF method and LES turbulence model.Prog Nucl Energ.2015;78: 201-15. doi: 10.1016/j.pnucene.2014.10.002.

［5］江帆，黃鵬.Fluent高級(jí)應(yīng)用.2008.

［6］Garimella S. Condensation flow mechanisms in microchannels: Basis for pressure drop and heat transfer models. Heat Transfer Eng. 2004;25（3）:104-16. doi: 10.1080/01457630490280489.

Numericalsimulation Offluid Dynamicsbased Oncryogenic Liquefactionprocess

LIN Yuanzai
（MechanicalEngineering,Southeast University,Nanjing 211189）

Abstract：This paper simulated methane gas cryogenic liquefaction process and multiphase flows at the same time in the heat exchanger by computing computational fluid dynamics（CFD）method. The paper established momentum conservation equation, continuity equation, volume conservation equation and interphase heat transfer equation of methane gas-liquid mixture by the Mixture Model, and tracked phase boundary by the Volume of Fluid（VOF）method..

Key Words:heat exchanger, liquefaction, phase change, Mixture Model, VOF