石油与天然气化工  2020, Vol. 49 Issue (1): 49-56
绕管式换热器壳侧汽化模拟时的液滴粒径研究
吴志勇1 , 张安昊1 , 国慧敏1 , 建伟伟1 , 蔡伟华2     
1. 辽宁石油化工大学石油天然气工程学院;
2. 哈尔滨工业大学能源科学与工程学院
摘要:使用ANSYS CFX软件下的欧拉模型对绕管式换热器壳侧乙烷剪切汽化流动进行了数值模拟,对控制方程组中影响相间作用的离散相液滴粒径展开了研究。从能量守恒角度推导出离散液相平均粒子直径公式,通过引入热相变比参数得以简化计算液相平均粒子直径,实现对壳侧剪切流汽化换热过程准确模拟,换热模拟偏差在±20%范围内。壳侧冷剂入口干度和工作压力是影响热相变比的主要因素,以入口干度和汽相密度为自变量对热相变比进行公式拟合,给出了适用于变工况下的热相变比计算式。
关键词绕管式换热器    剪切流    汽化换热    平均粒子直径    热相变比    
Study on liquid particle diameter during simulating vaporization of spiral wound heat exchanger shell-side
Wu Zhiyong1 , Zhang Anhao1 , Guo Huimin1 , Jian Weiwei1 , Cai Weihua2     
1. School of Petroleum and Natural Gas Engineering, Liaoning Shihua University, Fushun, Liaoning, China;
2. School of Energy Science and Engineering, Harbin Institute of Technology, Harbin, Heilongjiang, China
Abstract: The Euler model that attaches to the ANSYS CFX software is used to simulate the ethane shear and vaporizing flow on the shell-side of spiral wound heat exchanger (SWHE). The diameter of the disperse liquid particle that influences the phase interactions in the governing equations is investigated. Based on the energy conservation, formula for the average diameter of the disperse liquid particle is deduced. The heat ratio of phase change is proposed to compute simply the liquid particle average diameter so that the shear and vaporizing flow on the shell-side of SWHE could be simulated accurately and within the deviation of ±20%. The heat ratio of phase change is mainly influenced by the shell-side operating parameters including inlet quality and pressure, hence a formula of the heat ratio of phase change is fitted by both shell-side inlet quality and vapor density variables. The formula of the heat ratio of phase change is also applicable under different conditions such as mass velocity and heat flux variation.
Key words: spiral wound heat exchanger    shear flow    vaporization heat transfer    average diameter of liquid particle    heat ratio of phase change    

绕管式换热器(spiral wound heat exchanger,SWHE)作为一种结构复杂的换热器广泛用于天然气液化领域,具有抗高压、耐低温、多种介质可同时换热的特点。绕管式换热器内部构造见图 1[1-2],直径较小的换热管以螺旋方式分层缠绕在中央芯筒上,相邻两层换热管的缠绕方向不同,层间用金属隔条控制间距[3]

图 1     绕管式换热器构造

绕管式换热器用于天然气液化工艺时,天然气以多股流方式在换热管内自下向上地流动;而在换热器壳侧,烷烃制冷剂则自上而下地流过。管、壳两侧流动都存在相变换热,其中烷烃制冷剂以液相进入壳侧吸热汽化,先后经历降膜流、剪切流,并以过热气形式流出换热器。

在绕管式换热器壳侧换热性能研究方面,可借鉴的相关资料较少。Aunan[4]搭建了LNG绕管式换热器壳侧实验装置,对烷烃制冷剂在壳侧的流动与换热性能展开了研究,指出了干度、热流密度、压力、质量流率对换热的影响程度。Wu等[5]对绕管式换热器壳侧降膜流沸腾工况做出了数值模拟,基于Ansys Fluent软件下的VOF模型开展了传质时间松弛系数的研究,给出了乙烷、丙烷冷剂下的传质时间松弛系数的合理取值。李剑锐等[6]利用Fluent模拟出绕管式换热器壳侧降膜流下的沸腾换热特性,指出冷剂干度对换热系数影响大。季鹏等[7]对多个绕管式换热器壳侧换热关联式进行了对比分析,指出修正后的Abadzic关联式在天然气预冷段最为准确。Wu等[8]以数值模拟方法研究了绕管式换热器壳侧降膜沸腾时的换热特性,分析了利用VOF模型模拟剪切沸腾的不适用性。Ren等[9]对滚动工况下的绕管式换热器壳侧降膜流沸腾换热进行了数值模拟研究,认为滚动振幅与周期对换热的影响大于质量流率和压力的作用。

鉴于Wu等[8]指出了使用VOF模型模拟绕管式换热器壳侧剪切流沸腾换热的不适用性,所以本研究建议使用Ansys CFX软件下的欧拉模型对绕管式换热器壳侧剪切汽化流动进行数值模拟研究,并对欧拉模型下影响流动与换热强度的液相平均粒子直径这一关键参数的取值问题做深入研究,以便能够准确模拟换热器壳侧剪切流沸腾换热过程,进而掌握其流动与换热特性。

1 绕管式换热器壳侧模型的建立
1.1 壳侧几何模型

实际应用的绕管式换热器体积巨大,而内部结构却细致紧凑,在对其进行数值模拟研究时必将导致几何模型网格单元数量巨大,计算耗时长且难以收敛。因此,本研究按照文献[4]中实验装置建立壳侧几何模型,并通过文献[4]中的实验数据对模拟结果予以验证。

文献[4]中的实验装置是对绕管式换热器壳侧的简化,换热管分三层螺旋缠绕,中间层采用断面完整的换热管,而内、外两层使用断面为半圆的换热管,其结构见图 2。中间缠绕层下部并行的4根换热管壁面为加热壁面,由管内电加热丝控制热流密度,其他壁面均是绝热壁面。烷烃冷剂由换热器顶部30个分流孔进入壳侧,经过换热器内部多排换热管缓冲后进入到底部4根换热管测量段。测量段的加热壁面温度和流道断面温度由热电偶多点实测获得。该实验装置的几何参数见表 1

图 2     绕管式换热器壳侧模型

表 1    壳侧模型几何参数

鉴于该实验装置具有轴对称性,因此本研究在几何模型建立之后,将模型沿轴向切割36°作为仿真计算研究对象,其目的是降低网格数量,减少计算耗时。切割下来的壳侧36°几何模型见图 3

图 3     壳侧36°几何模型C

1.2 控制方程

在绕管式换热器壳侧,当冷剂入口干度x>0.2时,由于气、液相密度比很小,所以汽相占据壳侧绝大部分空间,液相以离散相液滴形式存在于汽相之中以及出现在换热管壁面上。在汽相吹扫作用下,换热管壁面会出现暂时的干涸现象。此时,汽、液间作用力明显,相间滑速比增大,在这种状态下适宜采用欧拉模型进行模拟计算。

欧拉模型控制方程组如下[10]

连续方程:

$ \frac{{\partial \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}} \right)}}{{\partial \tau }} + \nabla \cdot \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{{\vec v}_{\rm{g}}}} \right) = {{\dot m}_{{\rm{lg}}}} - {{\dot m}_{{\rm{gl}}}} $ (1)
$ \frac{{\partial \left( {{\beta _1}{\rho _1}} \right)}}{{\partial \tau }} + \nabla \cdot \left( {{\beta _1}{\rho _1}{{\vec v}_1}} \right) = {{\dot m}_{{\rm{gl}}}} - {{\dot m}_{{\rm{lg}}}} $ (2)

动量方程:

$ \begin{array}{*{20}{c}} {\frac{\partial }{{\partial \tau }}\left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{{\vec v}_{\rm{g}}}} \right) + \nabla \cdot \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{{\vec v}_{\rm{g}}}{{\vec v}_{\rm{g}}}} \right) = }\\ { - {\beta _{\rm{g}}}\nabla {p_{\rm{g}}} + \nabla \cdot \left[ {{\beta _{\rm{g}}}{\mu _{\rm{g}}}\left( {\nabla {{\vec v}_{\rm{g}}} + \nabla \vec v_{\rm{g}}^T} \right)} \right] + }\\ {{\beta _{\rm{g}}}{\rho _{\rm{g}}}\vec g + {{\vec F}_{1{\rm{g}}}} + {{\dot m}_{1{\rm{g}}}}{{\vec v}_1} - {{\dot m}_{{\rm{g1}}}}{{\vec v}_{\rm{g}}}} \end{array} $ (3)
$ \begin{array}{*{20}{c}} {\frac{\partial }{{\partial \tau }}\left( {{\beta _1}{\rho _1}{{\vec v}_1}} \right) + \nabla \cdot \langle {\beta _1}{\rho _1}{{\vec v}_1}{{\vec v}_1}) = }\\ { - {\beta _1}\nabla {p_1} + \nabla \cdot \left[ {{\beta _1}{\mu _1}\left( {\nabla {{\vec v}_1} + \nabla \vec v_1^{\rm{T}}} \right)} \right] + }\\ {{\beta _1}{\rho _1}\vec g + {{\vec F}_{{\rm{gl}}}} + {{\dot m}_{{\rm{gl}}}}{{\vec v}_{\rm{g}}} - {{\dot m}_{{\rm{lg}}}}{{\vec v}_1}} \end{array} $ (4)

能量方程:

$ \begin{array}{*{20}{c}} {\frac{\partial }{{\partial \tau }}\left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{h_{\rm{g}}}} \right) + \nabla \cdot \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{{\vec v}_{\rm{g}}}{h_{\rm{g}}}} \right) = }\\ {\nabla \cdot \left( {{\beta _{\rm{g}}}{\lambda _{\rm{g}}}\nabla {T_{\rm{g}}}} \right) + {Q_{{\rm{lg}}}} + {{\dot m}_{{\rm{lg}}}}{h_1} - {{\dot m}_{\rm{g}}}{h_{\rm{g}}}} \end{array} $ (5)
$ \begin{array}{l} \frac{\partial }{{\partial \tau }}\left( {{\beta _1}{\rho _1}{h_1}} \right) + \nabla \cdot \left( {{\beta _1}{\rho _1}{{\vec v}_1}{h_1}} \right) = \\ \;\;\;\;\;\;\;\;\;\;\;\;\nabla \cdot \left( {{\beta _1}{\lambda _1}\nabla {T_1}} \right) + {Q_{{\rm{gl}}}} + {{\dot m}_{{\rm{gl}}}}{h_{\rm{g}}} - {{\dot m}_{{\rm{lg}}}}{h_1} \end{array} $ (6)

式中:${\dot m_{{\rm{lg}}}}$为单位控制单元内液相向汽相的传质量,kg/(m3·s);${\dot m_{{\rm{lg}}}}$为单位控制单元内汽相向液相的传质量,kg/(m3·s);${\vec F_{{\rm{lg}}}}$为单位控制单元内液相对汽相的作用力,N/m3;${\vec F_{{\rm{gl}}}}$为单位控制单元内汽相对液相的作用力,N/m3,且有${\vec F_{{\rm{gl}}}}$=-${\vec F_{{\rm{lg}}}}$;Qlg为单位控制单元内液相向汽相的传热量,W/ m3Qgl为单位控制单元内汽相向液相的传热量,W/m3,且有Qgl=-Qlg。以上各变量通过下列模型进行求解。

(1) 相间传热模型

单位控制单元内液相向汽相的传热量可通过粒子模型进行求解,将液相定为离散相粒子公式如下[10]

$ {Q_{{\rm{lg}}}} = {\alpha _{{\rm{lg}}}}{A_{{\rm{lg}}}}\left( {{T_1} - {T_{\rm{g}}}} \right) $ (7)
$ {\alpha _{{\rm{lg}}}} = \frac{{{\lambda _{\rm{g}}}N{u_{{\rm{lg}}}}}}{{{d_1}}} $ (8)
$ N{u_{{\rm{lg}}}} = \left\{ {\begin{array}{*{20}{l}} {2 + 0.6R{e^{0.5}}P{r^{0.33}}}&{0 \le Re < 776.06\quad 0 \le Pr < 250}\\ {2 + 0.27R{e^{0.62}}P{r^{0.33}}}&{776.06 \le Re\quad 0 \le Pr < 250} \end{array}} \right. $ (9)
$ Re = \frac{{{\rho _{\rm{g}}}\left| {{{\vec v}_1} - {{\vec v}_{\rm{g}}}} \right|{d_1}}}{{{\mu _{\rm{g}}}}} $ (10)
$ Pr = \frac{{{\mu _{\rm{g}}}C{P_{\rm{g}}}}}{{{\lambda _{\rm{g}}}}} $ (11)
$ {A_{{\rm{lg}}}} = \frac{{6{\beta _1}}}{{{d_1}}} $ (12)

式中:αlg是汽、液相间换热系数,W/(m2·K);Alg是汽、液相界面面积浓度,m2/m3dl是液相粒子直径,m。

(2) 相间传质模型

单位控制单元内汽、液相间传质量${\dot m_{{\rm{lg}}}}$和${\dot m_{{\rm{lg}}}}$可通过热相变模型进行求解,计算公式为[10]

$ \dot m = \frac{{{Q_{{\rm{sg}}}} + {Q_{{\rm{sl}}}}}}{{{h_{{\rm{ls}}}} - {h_{{\rm{gs}}}}}} $ (13)
$ {Q_{{\rm{sg}}}} = {\alpha _{\rm{g}}}{A_{{\rm{lg}}}}\left( {{T_{\rm{s}}} - {T_{\rm{g}}}} \right) $ (14)
$ {Q_{s{\rm{l}}}} = {\alpha _1}{A_{1{\rm{g}}}}\left( {{T_{\rm{s}}} - {T_1}} \right) $ (15)
$ {{\dot m}_{{\rm{lg}}}} = \dot m,{{\dot m}_{{\rm{gl}}}} = 0,{h_{{\rm{gs}}}} = {h_{{\rm{g}},{\rm{sat}}}},{h_{1{\rm{s}}}} = {h_1}\quad (\dot m > 0) $ (16)
$ {{\dot m}_{{\rm{lg}}}} = 0,{{\dot m}_{{\rm{g}}1}} = - \dot m,{h_{{\rm{gs}}}} = {h_{\rm{g}}},{h_{1{\rm{s}}}} = {h_{1,{\rm{ sat }}}}\quad (\dot m < 0) $ (17)

式中:QsgQsl分别为单位控制单元内相界面向汽相和液相的传热量,W/m3αgαl分别为汽相和液相的换热系数,W/(m2·K);Ts为相界面温度,K。当液相粒子采用零热阻换热模型时,有Ts=Tl

(3) 相间作用力模型

在研究文中,汽、液相间作用力仅考虑曳力,根据Schiller Naumann的曳力模型,此时液相对汽相的曳力Flg[10]

$ \left. {{{\bar F}_{{\rm{lg}}}} = \frac{3}{4}\frac{{{C_{\rm{D}}}}}{{{d_1}}}{\beta _1}{\rho _{\rm{g}}}\left| {{{\vec v}_1} - {{\vec v}_{\rm{g}}}} \right|\overrightarrow {\left( {{v_1}} \right.} - {{\vec v}_{\rm{g}}}} \right) $ (18)
$ {C_{\rm{D}}} = \max \left\{ {24\left( {1 + 0.15R{e^{0.687}}} \right)/Re,0.44} \right\} $ (19)

此外,欧拉模型在湍流计算时,汽相作为连续相使用标准k-ε湍流模型,液相以离散粒子形式使用零方程模型进行简化处理,此时对应的湍流方程分别为[10]

汽相湍动能方程:

$ \begin{array}{*{20}{c}} {\frac{{\partial \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{k_{\rm{g}}}} \right)}}{{\partial \tau }} + \nabla \cdot \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{{\vec v}_{\rm{g}}}{k_{\rm{B}}}} \right) = }\\ {\nabla \cdot \left[ {{\beta _{\rm{g}}}\left( {{\mu _{\rm{g}}} + \frac{{{\mu _{{\rm{t}},{\rm{B}}}}}}{{{\sigma _{\rm{k}}}}}\nabla {k_{\rm{g}}}} \right] + {\beta _{\rm{g}}}\left( {{R_{{\rm{k}},{\rm{g}}}} - {\rho _{\rm{g}}}{\varepsilon _{\rm{g}}}} \right) + } \right.}\\ {{{\dot m}_{{\rm{lg}}}}{k_{\rm{l}}} - {{\dot m}_{{\rm{gl}}}}{k_{\rm{g}}}} \end{array} $ (20)

汽相耗散率方程:

$ \begin{array}{l} \frac{{\partial \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{\varepsilon _{\rm{g}}}} \right)}}{{\partial \tau }} + \nabla \cdot \left( {{\beta _{\rm{g}}}{\rho _{\rm{g}}}{{\vec v}_{\rm{g}}}{\varepsilon _{\rm{g}}}} \right) = \\ \nabla \cdot \left[ {{\beta _g}\left( {{\mu _g} + \frac{{{\mu _{t,g}}}}{{{\sigma _\varepsilon }}}} \right)\nabla {\varepsilon _g}} \right] + \\ {\beta _{\rm{g}}}\left( {{c_{1\varepsilon }}{R_{{\rm{k}} \cdot {\rm{g}}}} - {c_{2\varepsilon }}{\rho _{\rm{g}}}{\varepsilon _{\rm{g}}}} \right)\frac{{{\varepsilon _{\rm{g}}}}}{{{k_{\rm{g}}}}} + {{\dot m}_{{\rm{lg}}}}{\varepsilon _1} - {{\dot m}_{{\rm{g}}1}}{\varepsilon _{\rm{g}}} \end{array} $ (21)

液相湍流黏度:

$ {\mu _{{\rm{t}},1}} = \frac{{{\rho _1}}}{{{\rho _{\rm{g}}}}}\frac{{{\mu _{{\rm{t}},{\rm{g}}}}}}{\sigma } $ (22)

式中:σk=1.0,σε=1.3,c=1.44,c=1.92,σ=0.9。

1.3 边界条件

对于绕管式换热器壳侧36°几何模型,冷剂从上部3个分流孔进入,在底部端面流出,两个轴向端面设为对称边界,中间缠绕层下部4根换热管的外表面是加热壁面,其余壁面属于绝热壁面,边界条件设置见表 2图 3

表 2    边界条件

2 离散液相平均粒子直径推导

绕管式换热器壳侧剪切汽化流动模拟时,汽、液相分别设置成连续相和离散粒子,以便描述汽、液相间动量、能量交换,其中式(8)和式(18)的计算需要引入液相平均粒子直径,但该参数尚无计算方法,给剪切流汽化模拟带来困难,因此,需对液滴平均粒子直径作理论研究。

壳侧剪切汽化流动时,汽、液相温度并不相同,汽相温度高于液相温度,流场内弥散状液滴通过相界面从汽相获取热量实现汽化。液滴平均粒子直径可以通过微元控制体能量平衡推导出来。流场内部微元控制体见图 4,图中蓝色部分表示液相,其余部分表示汽相。

图 4     剪切流动下微元控制体

假设沿xyz方向流入微元控制体的单位体积能量为hV, xhV, yhV, z,相应的流速为uvw。则x方向流入能量为:

$ {H_x} = {h_{V,x}}udydz $ (23)

同理,yz方向流入能量为:

$ {H_y} = {h_{V,y}}vdxdz $ (24)
$ {H_z} = {h_{V,z}}wdxdy $ (25)

x方向流出能量为:

$ \begin{array}{*{20}{l}} {{H_{x + dx}} = \left( {u + \frac{{\partial u}}{{\partial x}}dx} \right)\left( {{h_{V,x}} + \frac{{\partial {h_{V \cdot x}}}}{{\partial x}}dx} \right)dydz}\\ { \approx \left( {{h_{V,x}}u + {h_{V,x}}\frac{{\partial u}}{{\partial x}}dx + u\frac{{\partial {h_{V,x}}}}{{\partial x}}dx} \right)dydz} \end{array} $ (26)

同理,yz方向流出能量为:

$ \begin{array}{*{20}{l}} {{H_{y + dy}} = \left( {v + \frac{{\partial v}}{{\partial y}}dy} \right)\left( {{h_{V,y}} + \frac{{\partial {h_{V,y}}}}{{\partial y}}dy} \right)dxdz}\\ {\quad \approx \left( {{h_{V,y}}v + {h_{V,y}}\frac{{\partial v}}{{\partial y}}dy + v\frac{{\partial {h_{V,y}}}}{{\partial y}}dy} \right)dxdz} \end{array} $ (27)
$ \begin{array}{*{20}{l}} {{H_{z + dz}} = \left( {w + \frac{{\partial w}}{{\partial z}}dz} \right)\left( {{h_{V,z}} + \frac{{\partial {h_{V,z}}}}{{\partial z}}dz} \right)dxdy}\\ {\quad \approx \left( {{h_{V,z}}w + {h_{V,z}}\frac{{\partial w}}{{\partial z}}dz + w\frac{{\partial {h_{V,z}}}}{{\partial z}}dz} \right)dxdy} \end{array} $ (28)

当汽、液相间换热采用两热阻模型且液相侧使用零热阻时,有αlg=αgTs=Tl。稳态时,微元控制体内汽相向液相传递的热量即是外界应向微元控制体输入的汽化热量,也即是微元体的净导入热量,该热量为:

$ {Q_{{\rm{lg}}}} = {\alpha _{\rm{g}}}{A_{{\rm{lg}}}}\left( {{T_{\rm{g}}} - {T_1}} \right)dxdydz $ (29)

稳定流动下根据微元控制体能量平衡可得到:

$ {H_x} + {H_y} + {H_z} + {Q_{{\rm{lg}}}} = {H_{x + dx}} + {H_{y + dy}} + {H_{z + dz}} $ (30)

dxdydz均趋于0,所以有:

$ {h_{V,x}} = {h_{V,y}} = {h_{V,z}} = {H_V} $ (31)

将式(23)~式(29)及式(31)代入至式(30),简化后得到:

$ \frac{{\partial \left( {{h_V}u} \right)}}{{\partial x}} + \frac{{\partial \left( {{h_V}v} \right)}}{{\partial y}} + \frac{{\partial \left( {{h_V}\omega } \right)}}{{\partial z}} = {\alpha _g}{A_{\lg }}\left( {{T_{\rm{g}}} - {T_1}} \right) $ (32)

将式(8)及式(12)代入至式(32)中,得到:

$ \frac{{\partial \left( {{h_V}u} \right)}}{{\partial x}} + \frac{{\partial \left( {{h_V}v} \right)}}{{\partial y}} + \frac{{\partial \left( {{h_V}w} \right)}}{{\partial z}} = \frac{{{\lambda _{\rm{g}}}N{u_{{\rm{g}}1}}}}{{{d_1}}}\frac{{6{\beta _1}}}{{{d_1}}}\left( {{T_{\rm{g}}} - {T_1}} \right) $ (33)

在流场范围Ω内对式(33)进行三重积分,有:

$ \begin{array}{*{20}{c}} {\iiint\limits_\Omega {\frac{\partial }{\partial }\left[ {\frac{{\partial \left( {{h_V}u} \right)}}{{\partial x}} + \frac{{\partial \left( {{h_V}v} \right)}}{{\partial y}} + \frac{{\partial \left( {{h_V}w} \right)}}{{\partial z}}} \right]dV = } }\\ {\iiint\limits_\Omega {\frac{{{\lambda _g}N{u_{g1}}}}{{{d_1}}}\frac{{6{\beta _1}}}{{{d_1}}}\left( {{T_{\rm{g}}} - {T_1}} \right)dV} } \end{array} $ (34)

式(34)左侧积分利用高斯定理可转化为:

$ \begin{array}{l} \;\;\;\;\iiint\limits_\Omega {\left[ {\frac{{\partial \left( {{h_V}u} \right)}}{{\partial x}} + \frac{{\partial \left( {{h_V}v} \right)}}{{\partial y}} + \frac{{\partial \left( {{h_V}w} \right)}}{{\partial z}}} \right]dV = } \\ \sum\limits_{i = 1}^n {\iint\limits_{{A_i}} {\left( {{h_V}u\cos {\vartheta _x} + {h_V}v\cos {\vartheta _y} + {h_V}w\cos {\vartheta _z}} \right)} } d{A_i} \end{array} $ (35)

式中:n为流场边界曲面数量,${\vartheta _x}, {\vartheta _y}, {\vartheta _z}$为边界曲面外法线与坐标轴夹角。将壳侧模型对称面视为壁面并且以无滑移处理壁面后,式(35)右侧可简化为:

$ \begin{array}{l} \sum\limits_{i = 1}^n {\iint\limits_{{A_i}} {\left( {{h_V}u\cos v{\vartheta _x} + {h_V}v\cos {\vartheta _y} + {h_V}w\cos {\vartheta _z}} \right)} } d{A_i} = \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;{H_{{\rm{out }}}} - {H_{{\rm{in }}}} = {Q_{{\rm{in }}}} = {A_{{\rm{heat }}}}q \end{array} $ (36)

式中:Aheat为换热管加热壁面的面积,m2q为换热管壁面热流密度,W/m2

式(34)右侧积分时利用中值定理可变为:

$ \begin{array}{l} \iiint\limits_\Omega {\frac{{{\lambda _{\rm{g}}}N{u_{{\rm{gl}}}}}}{{{d_1}}}} \frac{{6{\beta _1}}}{{{d_1}}}\left( {{T_{\rm{g}}} - {T_1}} \right)dV = \\ \frac{{6{\lambda _{\rm{g}}}\overline {N{u_{{\rm{gl}}}}} \cdot \overline {{\beta _1}} }}{{\overline {d_1^2} }}(\overline {{T_{\rm{g}}} - {T_1}} ){V_{{\rm{model}}}} \end{array} $ (37)

式中:Vmodel为模型壳侧体积,m3dl为液滴平均粒子直径,m;其余参量上方加横线表示为平均值。

对于球形粒子,Nugl可以近似取为2[10],式(37)变为:

$ \begin{array}{l} \iiint\limits_\Omega {\frac{{{\lambda _{\rm{g}}}N{u_{{\rm{gl}}}}}}{{{d_1}}}} \frac{{6{\beta _1}}}{{{d_1}}}\left( {{T_{\rm{g}}} - {T_1}} \right)dV = \\ \frac{{12{\lambda _{\rm{g}}}{{\overline \beta }_1}}}{{d_1^2}}(\overline {{T_{\rm{g}}} - {T_1}} ){V_{{\rm{model }}}} \end{array} $ (38)

将式(36)和式(38)代入至式(34)中,变形后得到:

$ \overline {{d_1}} = \sqrt {\frac{{12{\lambda _{\rm{g}}}\overline {{\beta _1}} \cdot (\overline {{T_{\rm{g}}} - {T_1}} ){V_{{\rm{model}}}}}}{{{A_{{\rm{heat}}}}q}}} $ (39)

式(39)即为推导出的平均粒子直径公式。使用该式时,λgq已知,VmodelAheat可以由壳侧几何模型确定,βl可由文献[11]推荐的Chisholm空泡率关联式先算出汽相平均体积分数,进而求得液相平均体积分数,但因$\overline {{T_{\rm{g}}} - {T_1}} $难以确定,因此使用式(39)仍存在问题。

3 热相变比的应用

为了解决$\overline {{T_{\rm{g}}} - {T_1}} $难以确定的问题,本研究引入热相变比的概念,其定义为:

$ \zeta = 1 - \frac{{m\left( {{x_{{\rm{in}}}} + \Delta x} \right)C{p_{\rm{g}}}\left( {{T_{\rm{g}}} - {T_1}} \right)}}{{{A_{{\rm{heat}}}}q}} $ (40)

式中:Δx为壳侧沸腾时的干度增加量。从定义式可以看出,热相变比ζ表征了向壳侧输入热量中,用于液滴汽化的热量占总热量的比例。

将式(40)变形后代入至式(39),得到:

$ {\bar d_1} = \sqrt {\frac{{12{\lambda _{\rm{g}}}{{\bar \beta }_1}(1 - \zeta ){V_{{\rm{model}}}}}}{{m\left( {{x_{{\rm{in}}}} + \Delta x} \right)C{p_{\rm{g}}}}}} $ (41)

式中:Δx可根据壳侧热流量和质量流量经计算而得到,故平均粒子直径dl的求解完全在于如何确定热相变比ζ

本研究依据文献[4]的实验数据,研究了乙烷在变工况条件下做剪切汽化流动时热相变比的取值问题。通过设定不同的ζ值,将模拟工况下的换热结果与实验数据做比对来确定合理的热相变比。模拟工况共选取23个,涵盖了变压力、变干度、变热流、变流量的情况,经过百余次计算,确定出来的ζ值见表 3

表 3    乙烷液相平均粒子直径及其工况参数

通过表 3可以看出,壳侧冷剂入口干度(x)对热相变比的影响最为敏感,两者间呈同向变化,即干度增大时要求热相变比相应变大,从而减小液相平均粒子直径,以满足提升换热系数的要求。而当质量流率(G)和热流密度(g)变化时,两因素对热相变比的影响不明显,即无法改变冷剂干度对热相变比的影响走势,因此,热相变比的取值可以忽略质量流率和热流密度的影响。

图 5是乙烷23个模拟工况按定压条件绘制出的热相变比变化图。从图 5可以看出,热相变比受压力的影响较为明显。在壳侧入口干度相同的条件下,压力越小则要求热相变比越大。其原因在于,高压向低压转变时, 相对于其他物性参数而言,汽相密度变化最为明显,汽相密度变小将导致汽相流速增大,随之而来的相间动量、能量传输增强,必然要求液相平均粒子直径进一步变小,即热相变比需要相应增大。

图 5     乙烷变工况下的热相变比

鉴于变压力时,汽相密度对热相变比有较大的影响,因此本研究对热相变比进行公式拟合时自变量包含了干度x和汽相密度ρg两个物理量,根据模拟工况下的计算数据拟合出的热相变比公式见式(42):

$ \zeta = \frac{{{a_1}{x^2} + {a_2}x + {a_3} + {b_1}{{(\rho /7.36)}^2} + {b_2}(\rho /7.36)}}{{{a_4}{x^2} + {a_5}x + {a_6} + {b_3}{{(\rho /7.36)}^2} + {b_4}(\rho /7.36)}} $ (42)

a1=-0.155,a2=0.088,a3=7.606,a4=-0.207,a5=0.034,a6=7.564;

b1=8.122,b2=-14.247,b3=7.879,b4=-13.859

式中:ρ/7.36为变压力时汽相密度发生改变而做出的修正,以乙烷0.4 MPa下汽相密度7.36 kg/m3作为基准进行热相变比修正。上式拟合程度较好,拟合相关系数R=0.993 9。

图 6图 7分别为乙烷于换热器壳侧做剪切汽化流动时在轴向剖面、中间缠绕层壁面处的汽相体积分数模拟云图,其模拟工况参数p=0.4 MPa、G=59.36 kg/(m2·s)、x=0.91、q=7 839.13 W/m2

图 6     剪切流轴向剖面汽相体积分数分布

图 7     剪切流中间层壁面汽相体积分数分布

通过图 6图 7可以看出,壳侧冷剂高干度条件下的两相流动趋于汽态流动,液相是以极细的液滴粒子形式弥散于汽相当中,因此,流道内每个网格单元的汽相体积分数接近于1。在中间缠绕层壁面上,由于液相润湿效应而使壁面局部位置处的液相体积分数略有增加。在近似于汽态流动、换热的条件下,只有将热相变比增至0.992 1,液滴平均粒子直径减小到0.000 004 m,才能准确地模拟出换热系数高达12 153 W/(m2·K)的汽化工况,并使模拟偏差低至1.66%。

4 结论

(1) 两相流欧拉模型适用于绕管式换热器壳侧剪

切汽化流动时的数值模拟研究,通过理论推导以及引入热相变比参数,可以确定出离散液相的平均粒子直径,进而实现对壳侧剪切流汽化换热的准确模拟,其换热模拟偏差在±20%范围内。

(2) 以乙烷在绕管式换热器壳侧剪切汽化流动为例,通过换热模拟结果与实验数据相对比的方式,确定出变工况条件下的热相变比取值。经分析发现,壳侧冷剂入口干度对热相变比的影响最为强烈,冷剂工作压力对热相变比的影响次之,而壳侧质量流率、热流密度对热相变比的影响可以忽略不计。

(3) 在分析热相变比影响因素之后,选取壳侧冷剂入口干度和汽相密度作为自变量对热相变比进行公式拟合,给出了相应的计算式,其拟合程度较好,相关系数R=0.993 9。

(4) 数值模拟研究结果表明,通过热相变比确定离散相平均粒子直径的方法,可以使欧拉模型准确地模拟出两相流相变换热过程,解决了相变模拟失真这一突出矛盾。该方法对两相流仿真计算具有实用价值。

    符号表

    下角标

参考文献
[1]
HELDT S. Near-optimal operation of LNG liquefaction processes by means of regulation[D]. Berlin: Berlin Institute of Technology, 2011.
[2]
王楚琦, 马颖. 绕管式换热器在天然气处理中的应用浅析[J]. 石油与天然气化工, 2014, 43(6): 602-606. DOI:10.3969/j.issn.1007-3426.2014.06.005
[3]
吴志勇, 陈杰, 浦晖, 等. LNG绕管式换热器结构与流通参数计算[J]. 煤气与热力, 2014, 34(3): A34-A39.
[4]
AUNAN B. Shell-side heat transfer and pressure drop in coil-wound LNG heat exchangers, laboratory measurements and modeling[D]. Trondheim: Norwegian University of Science and Technology, 2000.
[5]
WU Z Y, CAI W H, QIU G D, et al. Prediction of mass transfer time relaxation parameter for boiling simulation on the shell-side of LNG spiral wound heat exchanger[J]. Advances in Mechanical Engineering, 2014, doi: 10.1155/2014/275708.
[6]
李剑锐, 陈杰, 浦晖, 等. 绕管式换热器壳侧降膜流动和相变传热的数值模拟[J]. 化工学报, 2015, 66(增刊2): 40-49.
[7]
季鹏, 李玉星, 朱建鲁, 等. LNG绕管式换热器壳侧单相传热模型的优化[J]. 制冷学报, 2015, 36(2): 21-26. DOI:10.3969/j.issn.0253-4339.2015.02.004
[8]
WU Z Y, WANG H, CAI W H, et al. Numerical investigation of boiling heat transfer on the shell-side of spiral wound heat exchanger[J]. Heat and Mass Transfer, 2016, 52(9): 1973-1982. DOI:10.1007/s00231-016-1867-5
[9]
REN Y, CAI W H, CHEN J, et al. The heat transfer characteristics of shell-side film flow in spiral wound heat exchanger under rolling working conditions[J]. Applied Thermal Engineering, 2018, 132: 233-244. DOI:10.1016/j.applthermaleng.2017.12.092
[10]
ANSYS Incorporation. ANSYS CFX-solver theory guide[M]. Canonsburg: ANSYS Inc., 2011.
[11]
吴志勇, 高阳, 麻宏强, 等. 绕管式换热器壳侧沸腾时空泡率关联式筛选[J]. 天然气化工-C1化学与化工, 2018, 43(2): 93-99.