板翅式换热器微通道内汽化相变过程的数值模型

李剑锐，王皓显，武春林，陈慧，胡海涛*，丁国良，邢占洋，王丽春

（1-上海交通大学制冷与低温工程研究所，上海 200240；2-山西汾西机电有限公司，山西太原 030027；3-山西汾西重工有限责任公司，山西太原 030027）

0 引言

我国主要依赖煤炭、石油的能源结构已经造成越来越大的环境污染压力。逐步提升天然气这一清洁能源在我国能源供应结构中的比重是解决上述问题的重要措施[1]。板翅式换热器作为一种高效、紧凑、轻巧的换热设备，在天然气液化领域得到广泛的应用，其优点为传热效率高、结构紧凑、体积小、重量轻、牢固、适用范围大，但制造工艺复杂、开发费用高、流动通道压降大、易堵塞、易腐蚀、加工工艺复杂[2-3]。为了减少热损失、降低压降损耗、增加LNG液化效率，针对板翅式换热器的设计优化和性能改善是十分必要的。

在板翅式换热器中换热过程都伴随着相变，冷流体发生蒸发相变，相对于单相区对流换热，相变区域的潜热换热过程传热性能更强，是换热器中的关键的换热发生区域[4]。板翅式换热器性能提升的关键在于正确掌握两相相变过程中的热值传递特性。

针对板翅式换热器的两相研究主要集中于不同结构分配器的分配性能分析以及两相均布器的性能优化。袁培等[5]以空气-水为介质，通过实验模拟系统研究了入口两相分配器的分配特性；张哲等[6]通过实验系统研究了气液分配器的特性，结果表明高雷诺数下分配不均现象更显著；林彬彬等[7]利用ASPEN PLUS软件，模拟了不同两相流分配的情况下对传热性能的影响；潘良高等[8]考虑表面张力和重力的影响，建立了预测微通道气液两相流型的数值模型；焦安军等[9]研究了封头结构对两相分配的影响，通过改变封头结构的设计改善换热器的物流分配；万智华等[10]通过CFD研究了实际板翅式换热器冷箱并联管路的两相流均布特性，结果表明当管道水平布置时，优先考虑水平上进下出式；当管道垂直布置时，优先考虑垂直向下流入式。

针对板翅式换热器的数值模拟研究主要集中于不同结构间的单相传热以及压降性能比较。祝银海等[11]通过CFD对平直形和锯齿翅片形通道流体的流动与传热特性进行了研究，结果表明在锯齿形翅片的相邻2个锯齿的交错面上，流体的局部传热系数和压力存在突变，流体的边界层厚度要薄于在平直形翅片中的厚度；王臣等[12]通过CFD模拟了耦合导热情况下的不同结构平直形翅片形通道的流动传热特性，并提出了优化结构设计；魏进家[13]通过模拟单相工况下的流动通道内温度、速度分布，拟合得到了Nu、α*及f与Re的关联式；董其伍等[14]通过数值模拟方法，对比不同结构下的Nu和Re变化规律；曲乐等[15]通过分段物性曲线的方法进行了板翅式换热器内部传质过程的初步模拟，并与MUSE结果进行比较，计算结果有一定合理性。

而目前关于板翅式换热器内的两相区域传热过程机理的研究则较少。因此，需要建立适用的传热传质模型，以反映板翅式换热器的汽化相变过程。

1 数值模型

1.1 模型对象描述

通常板翅式内部结构如图1所示，由流道的入口流入后，通过封条未堵住的通道开口，流入对应流道的换热翅片层，经过一定长度的均流分布段后，均匀流入包含翅片的换热通道，并与相邻层间不同流道的介质发生换热，在经过均流段汇集至出口对应的封头区域，最后从出口流出。

流体在流动通道中发生的传热过程根据流体相态可以分为3个阶段：过冷段、两相段、过热段。对于冷流体在整个过程中存在汽化沸腾现象，其中带来的不同相间的热值过程对传热的影响十分显著。则正确描述两相段的汽化相变过程在板翅式换热器的传热性能研究中十分关键。

板翅式换热器中，冷流体在低干度下的汽化相变过程属于流动沸腾，此时主要的传质发生区域是聚集在换热表面。对于产生的气泡，热容相对小，局部温度会上升，使得周围液相汽化，气泡尺寸增长，然后脱离壁面。对于中干度情况下，环状流的气液发生分层，换热表面接触的主要相仍为液相，因此传质过程与低干度类似，但在气泡成长后会发生突入中间气相部分的现象，在液膜较薄的情况下，可能会发生液膜断裂，气相与换热壁面直接接触，传热性能开始发生下降。在高干度雾状流，换热壁面主要接触气相，为显热传热，气相热容小，温度上升快，传热性能急剧下降，传质区域主要发生在两相的相界面。

换热通道内介质两相间的相互作用主要有以下3种：1）气相成核、脱离过程；2）气泡成长突破液膜的过程；3）气相流速增大，撕裂液膜、夹带液滴的过程。在整个汽化相变过程中，起到决定性作用的力有5种，分别是气相粘性力、液相粘性力、气液相间表面张力、气液相间的剪切作用力、重力，如图2所示。

通过在模型的连续性方程中加入气液两相质量传递源项，在动量方程中加入表面张力源项和剪切力源项，在能量方程中加入潜热传热源项，从而将汽化相变过程需要考虑的各个因素反映到控制方程中。

图1 板翅换热器结构示意图

图2 气泡和相界面的受力分析

1.2 控制方程

根据N-S方程，对于图示的3种类型的控制单元，针对汽化相变模型可建立如下的基本控制方程。

气相和液相的连续性方程如下：

动量和能量方程如下：

式中，αl和α1分别代表控制单元内部气相和液相的体积分数；右侧为传质质量源项Sm，反映传质过程中两相间的传质质量；Fσ为表面张力项；∇⋅(k∇T)为显热源项，Q为潜热源项；Fσ表面张力项，可通过连续表面张力（CSF）模型求取：

其中ˆn为相界面函数、θ为接触角。

传质源项Sm，分为蒸发和冷凝2个部分进行考虑。对于传热介质为单工质时，不会存在蒸发与冷凝同时发生的情况，其中一项的值应为0。

在板翅式换热器的换热通道中，气相相间传质主要在2种情况下发生：高干度下的相界面传质、换热通道壁面的液相汽化传质，分别介绍如下。

1）高干度下的相界面传质

由于高干度工况下气相温度高，此时的传质过程主要为蒸发过程（Tl＞Tsat）：

2）换热通道壁面的液相汽化传质

在低干度工况下，换热表面的流体为液相，所有获得的热量都用于汽化相变。根据壁面当地的热流密度，蒸发传质质量计算如下（α1＞ 0，Tl＞Tsat）：

潜热传热源项Q，根据质量传质源项Sm进行计算，如下所示：

基于对以上各个模型进行综合，建立板翅式换热器通道内部的汽化相变数值模型，从而实现对冷流体侧汽化相变过程进行完整的数学描述。

2 求解方法

本文数值模拟使用的模型如图3所示，选用了平直形板翅式换热器的换热通道为例，基于商用软件FLUENT进行流动模拟。采用VOF模型作为模拟两相流模型，连续表面张力模型（CSF）作为模拟表面张力模型，实现两相分布、流型转变过程的模拟；同时通过FLUENT的用户自定义方程(UDFs)分别建立针对壁面及非壁面的传质模型，以实现通过对网格类型的判断以采用不同方式计算传质质量，模拟实现气泡形成、突破液膜的过程；采用VOF-CSF模型作为表面张力模型。

模拟使用六面体结构化网格进行网格划分，对通道表面的边界层区域以及可能出现气液交界面的区域进行了加密。网格独立性验证结果如下，在最大网格尺寸不超过0.05 mm时且壁面附近的网格单元尺寸不超过0.0025 mm时，换热系数及传质质量计算误差低于2%。因此采用0.05 mm作为基准网格尺寸及0.0025 mm作为边界层加密尺寸，保证模拟结果精度。

图3 几何模型及网格划分

3 结果分析

3.1 模型验证

选取文献[16]圆管中天然气两相传热工况进行两相传热数值模拟，并与实验数据[16]进行对比，以验证本文中模型的有效性。模型验证结果如图4所示。

通过与实验结果对比，两相工况下传热趋势与实验结果趋势一致，最大的误差为10.57%，出现在低干度区域。平均传热系数的偏差为6.6%。

图4 模型的实验验证

3.1 不同干度下流型模拟

分别针对质流密度对应的工况进行模拟，对应的近饱和流的相变情况如图5所示。

在低干度工况下，气相占据空间较少，气泡以离散形式存在，受到液相拖曳而流动，此时气液两相的流速是相同的，不存在滑移速度。

在中干度工况下，由于气相增加，占据空间也随之增加，气泡之间形成连续的气流，从而将液相阻断，形成气液交界面。在换热器的流动通道为水平且流动速度较低的情况下，由于重力的作用，使得液相会向通道一侧聚集，此时形成的流型为分层流；在换热器流动通道方向为竖直方向的情况下，液相不会发生单层聚集的情况，此时由于液相粘性力及气液相间表面张力的作用下，液相会趋向于附着在通道表面，而气相在通道中心流动，气相的流动速度会高于液相，产生滑移速度，此时的流型为环状流；在通道水平且高质流密度的情况下，也会形成环状流，但上下两侧的液膜厚度会有所差异。

在高干度工况下，由于气相已经占据了绝大部分的通道空间。由于液相体积不足无法形成液膜，液相只能以小液滴的形式存在，受到高速气流的卷袭而随之流动，而部分液相会因重力作用趋于一侧的现象依然存在。

图5 流型模拟结果（质流密度均为20 kg/m2s）

3.2 相变传热传质过程模拟

在引入了传热传质模型后，板翅换热器通道内传热传质过程的模拟结果如图6所示。

由图中可看出，液膜在壁面与高温的管壁发生换热，由于达到饱和压力，发生汽化成核，并随着传质质量的增加而长大，然后跟随液相的流动离开。随着干度增加后，部分较大气泡相互接触，形成更大的气泡，并受到液相的拖曳作用而拉长，由于重力的作用下，气泡开始附着到了通道的一侧，形成了分层流。在低质流密度的工况下，由于流速较慢，气泡受到的拖曳作用也越小，更趋于形成气液分层的情况，此时气相与壁面直接换热，传热系数下降。而对于高质流密度的工况而言，由于气泡难以在低干度下聚集变大，因此壁面上主要还是液相，从而传热系数相对于低质流密度要高。

图6 不同质流密度下的流型

3.3 不同传热温差工况传热性能对比

质流密度为20 kg/m2kg，不同热流密度的工况的传热系数如图7所示，传热系数随着热流密度增大而增大，并在中低干度下增幅明显，而在高干度下传热系数的增幅下降。而传热系数随干度上升呈现先上升后下降的趋势，不同热流密度下传热系数的下降节点都在0.4干度附近。

热流密度为6,000 W/m2，不同质流密度的工况的传热系数如图8所示，传热系数随着质流密度增大而增大，并呈现正比上升的趋势。不同质流密度下随着质流密度上升，传热系数最高点出现时的干度上升，质流密度60 kg/m2s时，最高点则出现在0.4附近，而质流密度60 kg/m2s时出现在0.6附近。

图7 不同热流密度下的传热系数

图8 不同质流密度下的传热系数

4 结论

1）基于VOF模型、连续表面张力模型、接触角模型，建立了板翅式换热器汽化相变过程的数值模拟模型。

2）针对相变过程中发生的传质现象，在连续性方程引入质量传递项，以预测在通道内部发生的气泡产生过程；在能量方程加入潜热项，以预测伴随着传质过程的介质温度变化。

3）对不同干度下的饱和状流、泡状流、环状流、雾状流等流型进行了模拟，结果显示流型不仅与干度工况有关，还与质流密度有关。

4）对不同质流密度工况下的全液相的传质模型进行模拟，高质流密度工况下不易形成气液分层现象，有利于传热。

5）传热系数随着热流密度增大而增大，并在中低干度下增幅明显；传热系数随干度上升呈现先上升后下降的趋势；传热系数随着质流密度增大而增大，传热系数最高点出现时的干度上升。

[1]曹寅, 姚斌. 天然气热电联产系统在上海地区的应用[J]. 制冷技术, 2003, 23(4): 26-29.

[2]凌祥, 涂善东, 陆卫权. 板翅式换热器的研究与应用进展[J]. 石油机械, 2000, 28(5): 54-58.

[3]贾承造, 张永峰, 赵霞. 中国天然气工业发展前景与挑战[J]. 天然气工业, 2014, 34(2): 8-18.

[4]文键, 李亚梅, 王斯民, 等. 板翅式换热器平直翅片表面流动及传热特性[J]. 化学工程, 2012, 40(10): 25-28, 39.

[5]袁培, 姜国宝, 张菲妮, 等. 板翅式换热器两相流分配器[J]. 化工学报, 2011, 62(S1): 31-36.

[6]张哲, 田津津, 厉彦忠, 等. 板翅式换热器两相流流体分配特性的实验研究[J]. 低温工程, 2013(6): 21-25.

[7]林彬彬, 韦小雄, 周寒秋, 等. 板翅式换热器两相流混合性能模拟分析[J]. 石化技术, 2016, 23(4): 10-12.

[8]潘良高, 徐琛, 柏祥华, 等. 微通道内气液两相流型的数值模拟[J]. 制冷技术, 2014, 34(4): 8-12.

[9]焦安军, 厉彦忠, 张瑞, 等. 封头结构对板翅式换热器物流分配不均匀性的影响[J]. 化工学报, 2003, 54(7): 907-912.

[10]万智华, 朱蒙生, 蔡伟华, 等. 并联管路气液两相流量均布特性数值模拟研究[J]. 化学工业与工程技术,2015, 36(2): 52-57.

[11]祝银海, 厉彦忠. 板翅式换热器翅片通道中流体流动与传热的计算流体力学模拟[J]. 化工学报, 2006, 57(5):1102-1106.

[12]王臣, 钱奕枝, 陈斌. 板翅式换热器的数值模拟研究[J]. 建筑节能, 2015(5): 41-44.

[13]魏进家, 刘海燕, 龙延. 板翅式换热器流动和换热性能研究[J]. 工程热物理学报, 2012, 33(10): 1785-1788.

[14]董其伍, 王丹, 刘敏珊, 等. 板翅式换热器数值模拟研究[J]. 化工设备与管道, 2008, 45(2): 25-27.

[15]曲乐, 贾林祥. 相变换热混合工质板翅式换热器流动与传热数值模拟[J]. 低温与超导, 2008, 36(4): 23-28.

[16]陈东升. 液化天然气管内两相流动与传热特性研究[D].上海: 上海交通大学, 2014.