针对冬季天然气门站调压时由于节流效应产生的管道冻堵问题，国内外使用最多的方式为预加热天然气弥补调压温降法[1-2]。浸没燃烧天然气加热装置作为一种新型天然气加热装置，相较于其他加热方式，具有结构紧凑、占地面积小、热效率高等优势[3]。
装置主要由燃烧系统、换热系统、控制系统等部分组成。换热系统包括烟气-水换热与水-天然气换热两部分，燃烧系统燃烧产生的高温烟气通过鼓泡孔鼓出进入水浴中，形成气液两相流冲刷天然气换热器以达到加热天然气的目的。因气液两相流环境比较复杂，导致在进行装置前期换热器传热设计计算时的关键参数水浴流速取值比较复杂，而设计时，水浴流速取值较小是造成管外对流换热系数较小的主要原因[4]。因此，获得准确的管外水浴最大流速对浸没燃烧天然气加热装置的换热器设计至关重要。
由于浸没燃烧天然气加热装置结构特殊，现有的测量方法均不能有效地确定换热器管外水浴最大流速。王浩等[4]基于首例装置的实际使用情况，利用装置实测数据，根据经典传热公式计算管外对流换热系数，然后通过实际管外对流换热系数反算得到水浴最大流速值，拟合得到功率-水浴流速的关系式，并进行了验证。
本研究在此基础上通过数值模拟的方法，对浸没燃烧天然气加热装置水箱内水浴的换热过程进行了传热特性分析，并计算验证了管外水浴最大流速值，提供了一种可参考的计算水浴流速的数值模拟思路。本研究对后续浸没燃烧天然气加热装置的换热器设计也具有指导意义。
浸没燃烧天然气加热装置换热原理如图 1所示，燃烧产生的高温烟气通过烟管输送，经由鼓泡孔鼓入水浴中形成气液两相流横掠天然气换热器，从而达到加热天然气的目的。
以实际运行中的加热装置为模型选取计算区域。加热装置的整体实物模型见图 2，其内部主要换热区由天然气换热器、烟管、鼓泡孔、挡板4部分组成。内部局部结构模型见图 3。由于实际模型尺寸相对较大，本研究以单元厚度对水箱截取三维几何模型，结果如图 4所示，沿着换热管束轴线方向，模型中换热区域结构具有较高的对称性，因此，在保证计算结果准确的前提下，最终建立了二维模型进行模拟计算。简化后的二维计算模型如图 5所示，计算域尺寸为0.80 m×1.10 m，换热管束直径D=34 mm，烟管直径D=159 mm，鼓泡孔出口直径D=20 mm，均为实际尺寸。
燃烧后的高温烟气在水浴中的流动与传热过程严格遵循物理守恒定律，基本的守恒定律主要包括质量守恒定律、动量守恒定律和能量守恒定律。基本控制方程如下：
(1) 质量守恒方程
式中：ux、uy、uz分别为x、y、z 3个方向的速度分量，m/s；t为时间，s；ρ为密度，kg/m3。
(2) 动量守恒方程
式中：p为微元体上所受到的压力，Pa；τxx、τxy、τxz为因分子黏性作用而产生的作用在微元体表面上的黏性应力τ的分量，Pa；Fx、Fy、Fz为微元体上的体积力，Pa。
(3) 能量守恒方程
式中：cp为流体的定压比热容，J/(kg·K)；T为流体的温度，K；K为流体的传热系数，W/(m2·K)；ST为流体由于黏性作用内部热源及机械能转化为热能的部分，W。
从湍流模型的适用场合和计算精度考虑，本研究中所选湍流模型为RNG k-ε湍流模型，具体系数取值可参考文献[5]。其湍动能、湍动能耗散率的控制方程分别见式(6)、式(7)。
式中：Gk、Gb分别为平均速度梯度和浮力影响引起的湍动能，m2/s2；YM为耗散率，受可压缩湍流脉动膨胀的影响；C1、C2均为模型常量；αk、αε分别为湍动能和耗散率对应的有效普朗特数。
利用前处理软件Gambit进行网格划分，考虑到计算模型中换热管束区域结构的不规则性，本研究使用结构网格和非结构网格相结合的混合网格。挡板内的换热管束区域及烟管部分均采用非结构网格进行划分并对其进行加密，其余区域使用计算质量较高的结构化网格进行划分。考虑计算结果精度和计算速度，最终选择数目为154 966的网格作为最优计算模型，换热管束区域网格的最小尺寸为10-3 m，局部网格分布如图 6所示。
将网格导入FLUENT 2D求解器中，选择基于压力的分离求解算法，计算采用瞬态求解器，速度形式选Absolute，重力加速度沿着y轴方向，设为-9.81 m/s2。模型设置里，选择VOF多相流模型、能量方程、RNG k-ε湍流模型。在Database中添加H2O，两相的物性参数均设置为常数，设置高温空气来代替烟气，因为高温烟气具有温度高、密度小、压力低的特点，可看做理想气体，而且在传热过程中，烟气在水浴中的少量溶解对传热影响较小，可忽略[6]。烟气物性参数根据入口状态下的烟气温度计算，烟气密度ρ、比热cp、动力黏度μ、运动黏度γ的计算方法参照文献[7]。导热系数λ、普朗特数Pr的计算方法参照文献[8]和文献[9]。
(1) 入口条件
鼓泡孔出口设为速度进口，速度值由烟气总体积流量(考虑温度修正)与烟气出口总面积换算，温度和烟气体积流量均由燃烧计算所得[10]，燃烧计算时过剩空气系数均为1.4，具体设置参数见表 1。
(2) 出口条件
由于烟气换热结束后从箱体顶部自然排出，故烟气出口设为压力出口，其值p=101 325 Pa。
(3) 固体边界条件
水箱的底部、侧壁、换热管束、烟管和挡板设置为壁面(wall)，其材料设置为不锈钢，壁面设置为300 K恒定温度。定义为不可渗透、无滑移的绝热边界条件；
(4) 区域条件
计算区域的主相(primary phase)设置为水，次相(second phase)设置为烟气，相间作用力(inter-phase forces)设置为0.072 5。
Solution Methods中，压力速度耦合方式采用PISO算法求解，压力项选择PRESTO！，其中的能量、动量方程均选二阶迎风离散格式，湍流动能、耗散率方程均选一阶迎风离散格式。在迭代过程中，各项亚松弛因子均保持默认数值。创建监测面，对烟气出口温度、水浴温度等参数进行监测，以其数值趋于稳定、残差收敛来判断计算达到收敛。在残差收敛标准中，对能量方程设定为10-6，其他变量收敛标准均设定为10-4。对全局初始化，时间步长设定为10-4 s。利用Patch功能设定水浴高度。初始化后的流场压力为101 325 Pa，水浴高度为0.8 m，水浴温度为305 K。初始化后的水浴相态、温度云图如图 7、图 8所示。相态分布图里，蓝色部分为水，红色部分为气体。
用模拟结果与加热装置功率为80 kW、128 kW、180 kW时的实际运行数据对比验证模型，结果如图 9所示。由图 9可以发现，模拟得到的3个功率下的换热量、水浴温度与实际运行数据基本吻合，换热量误差在5%以内，水浴温度误差在2%以内。因此，可认为上述数值模型选型可真实反映水浴与换热管束的传热过程。
图 10为加热装置功率80 kW、128 kW、180 kW稳定运行时的相态云图，蓝色部分为液相，红色部分为气相。模拟结果直观反映了加热装置运行时的气液分布情况。首先，从图 10中可以看出，各个功率燃烧产生的高温烟气从鼓泡孔鼓出，先形成较大的气泡，气泡上升到换热管束区域后，受换热管束的干涉和扰动作用，大气泡破碎形成体积较小、大小不一且数量较多的气泡，并且均匀分布在换热管束区域，功率越大，这种现象越明显。这表明装置换热区域结构布置合理，气泡破碎效果较好，增大了高温烟气与水浴的接触面积，进而强化了气液间传热。
其次，从图 10中还可看出，装置小功率运行时，换热一开始集中在换热管束区域内，随着功率的增大，烟气量增加，烟气入口速度升高，气泡增多，水浴中含气率增大，气泡对水的扰动效果明显增加，水浴液位抬升形成溢流，上升到水浴表面的大部分气泡没有及时排出便被水浴裹挟，重新进入到循环溢流中，增强了扰动。液位高度抬升后，挡板左右两侧区域的水浴中也分布了大量气泡，说明在挡板左右两侧形成了液体回流。
图 11是加热装置在功率为80 kW、128 kW、180 kW运行时的温度云图。从图 11中模拟结果可明显看到，鼓泡孔鼓出的高温烟气温度远大于水浴，但是在其上升与水浴混合的过程中，在换热管束区域下方两者温度就已经达到平衡，这是因为烟气的比热远小于液体，烟气喷射到水浴中与水换热时，其携带的热量迅速与水浴进行交换，使其温度骤降，在换热管束区域下方就达到了水浴平衡温度，说明在换热管壁处气泡几乎不参与换热。而且在整个模拟时间范围里，换热管束区域的水浴温度分布整体都比较均匀，温度几乎是稳定不变的，说明气液两相流是以一个相对恒定的温度冲刷天然气换热器的。
图 12为加热装置在功率为80 kW、128 kW、180 kW运行时的速度云图。数值模拟结果显示，在各功率运行状态下，受燃烧功率影响，烟气进入到水浴中的速度并不相同，但烟气离开鼓泡孔后，在换热管束下方区域均保持较高速状态。在换热器区域，由于换热管束的扰动作用，速度受到阻挡开始向两侧分散，管束中部区域速度大于两侧区域，最终速度分布呈现蝶形态势。在换热管器区域上方的局部位置存在速度增大的现象，这是因为气泡从水浴底部上升到自由界面时，具有一定的初始速度，同时气泡内部压强小于外部压强发生破裂，烟气迅速脱离水浴。由挡板左右两侧的速度分布可以看出，形成了水浴循环溢流。对比不同功率下水浴速度分布，可以看出功率增大时烟气对水浴的扰动作用较强，水浴整体流速增大，冲刷换热管束的过程加剧，更有利于换热。
在管外换热计算中，水浴流速取值为管束最窄处的最大值，如图 13所示，计算收敛后，以换热管束间隙方向为X轴，高度方向为Y轴，在管束中间定义了10个X方向坐标和8个Y方向坐标，分别对这80个位置点的瞬时速度大小分布进行了研究。从相态云图中可知，换热管束区域存在大量烟气形成的气泡，为了更加直观的看到气液两相的速度分布，本研究通过各监测点的密度对烟气与水浴的速度进行区分。读取数据时选用点的Magnitude Velocity, 其为标量，代表点的综合速度，单位为m/s。
各功率下模拟计算出的各监测点的速度分布如图 14所示，图中蓝色点表示烟气速度，红色点表示水浴速度。从图 14中可看出, 水浴和烟气形成的气液两相流受换热管束影响，水浴速度在换热区域高度方向分布杂乱无章，在宽度方向呈现中间高两边低的趋势，烟气整体速度大于水浴速度，这是因为本研究中加热装置换热管束下方只设置1根烟管，烟气携带的初始动能只能从管束中部向两侧传递。当直接取所有点即烟气与水的平均速度时，管外水浴流速值与理论计算所得偏差较大，而且从图 14中也可看出，其平均水浴速度曲线波动幅度较大。
根据2.2节中对水浴温度分布的分析可知，尽管换热介质为气液两相流，但是烟气几乎在到达换热器区域前就达到了水浴平衡温度，在换热管壁处时气泡几乎不参与换热，因此可忽略。依据上述分析，假定靠近换热管束壁面处的速度都是液相的速度，将烟气速度剔除，处理后的平均水浴速度分布曲线如图 14中所示，其速度值相对稳定，并且处理后水浴速度模拟结果与结合实测数据使用经典传热公式推导计算出来的数据很吻合[4]，如表 2所示，模拟计算最大相对误差为8.33%，最小相对误差为2.94%，表明水浴速度模拟计算结果的准确性，可为设计相同结构换热器时水浴最大流速取值提供参考。同时也证明，经典传热公式仍适用于浸没燃烧天然气加热装置换热器的设计计算。
结合理论与装置实测数据，使用商用流体力学软件FLUENT对浸没燃烧天然气加热装置水浴箱内水浴的换热过程进行数值模拟，得到以下结论：
(1) 加热装置的换热量与水浴温度的模拟结果与实测数据吻合较好，相对误差在5%以内，说明该模型能够较为准确地模拟装置的水浴换热过程。
(2) 模拟结果直观反映了水箱内部水浴流动及烟气与水浴间的换热过程，装置换热区域结构布置合理，气泡破碎效果很好，增大了高温烟气与水浴的接触面积，进而强化了气液间传热。
(3) 在加热装置水浴传热过程中，气液两相流是以一个相对恒定的温度冲刷换热管束的，烟气温度降低后形成的气泡在换热器管壁处几乎不参与换热。
(4) 通过在换热器区域创建监测点、分析各点的速度分布，最终验证了前者结合实测数据使用经典传热公式推导计算出的水浴流速值，最大相对误差为8.33%，为后续浸没燃烧天然气加热装置的研发和设计提供了一种可参考的计算水浴流速的数值模拟思路。同时也证明，经典传热公式仍适用于浸没燃烧天然气加热装置换热器的设计计算。
2020, Vol. 49 Issue (6): 58-65
2020, Vol. 49 Issue (6): 58-65