基于拟沸腾理论的超临界CO2 管内传热恶化量纲分析*

张海松 徐进良2)† 朱鑫杰

1) (华北电力大学低品位能源多相流与传热北京市重点实验室，北京 102206)

2) (华北电力大学电站能量传递转化与系统教育部重点实验室，北京 102206)

超临界流体广泛应用于工程技术领域, 其流动传热特性对工程设计具有重要意义, 但是, 由于超临界流体的物理微观和宏观行为的机理尚不清晰, 所以其异常的流动传热特性并未得到很好的解决.普遍认为超临界流体在分子尺度上可分为类气和类液两种不同的特性, 直到最近通过实验在宏观上监测到超临界水类液和类气之间的转变, 且这一过程与拟沸腾理论一致, 使得问题逐渐变得清晰.本文基于拟沸腾理论对超临界CO2 异常流动传热行为进行了研究, 在假设类液和类气转换过程不均匀的情况下, 从经典的量纲分析和亚临界过冷沸腾理论模型出发, 提出了一个适用于超临界流体拟沸腾换热过程的分析方法.通过引入表征类气膜生长速度与流体主流平均流速之比π = (qw·ρl)/(G·Δi·ρg)和表征近壁区类气膜温度梯度π13 = (qw·βpc·di)/λg 两个无量纲数, 来表征拟沸腾如何导致传热恶化, 解释了超临界CO2 竖直向上加热流动过程中的异常换热特性, 即较大的类气膜生长速度使近壁区快速聚集了较多的高温流体, 而较大的类气膜温度梯度使类气膜覆盖在壁面.当核心的冷类液不能充分润湿热壁面时, 传热恶化.新无量纲数较好的诠释了超临界流体拟沸腾诱导传热恶化机制, 为超临界拟沸腾传热研究提供了理论依据.

1 引 言

近年来, 超临界流体(supercutical fluid, SCF)广泛应用于食品加工、药物制取、材料合成、微纳米系统、污水处理及蒸汽发电等领域[1−2].在能源利用领域, 当流体的压力和温度升高时, 通常会提高系统热力学循环效率, 使系统结构更加紧凑.因此, 在基于热力学循环性能的工程系统中, SCF 已经为高压流体系统打开了一个新的领域.但是, 按照传统热力学思路, SCF 为绝对均匀的单相流体,在研究其强化传热及传热恶化现象过程中, 引入单相传热中的浮升力和流动加速效应, 而各研究者提出的计算公式适用参数范围窄, 不能统一其他研究者的实验数据[3], 这在工程应用中是不能接受的.按照单相传热思路, 难以理解SCF 异常的流动换热机理及较好地预测其传热特性, 这导致SCF 相关动力循环设计需要依赖大量的实验数据, 以及付出很大的人力和物力[4].

自20 世纪60 年代拟沸腾概念提出以来, Knapp 和Sabersky[5]研究了超临界CO2池式传热, 采用0.254 mm 直径镍铬金属丝为加热元件, 在近临界或略大于临界压力下, 观察到了类气泡传热,qw~T 曲线斜率变陡, 具有明显的沸腾传热特征,还发现了类膜态沸腾传热, 并观察到气、液界面.Ackerman[6]认为在超临界水管内流动传热实验过程中, 存在类膜态沸腾现象, 内螺纹结构可抑制膜态沸腾发生.Stewart 等[7]在超临界水均匀加热实验条件下, 观测到了压力脉动和壁面温度的周期性振荡, 认为这是类似“气膜”破裂造成的.Ambrosini[8]研究结果表明, 在加热通道内的超临界流体, 可能会出现密度波振荡和Ledinegg 不稳定性, 这和亚临界压力下的沸腾表现出了非常相似的定性和定量行为.从公开的文献看, 超临压力下没有宏观可见的气、液界面[9], 但是, 是否存在与流体密度变化相关的相变一直存在争议.近年来, SCF 基础研究已初步揭示不均匀异质结构端倪, 2010 年, Simeoni等[10]采用非弹性X 射线散射及MD 模拟表明, 尽管流体穿过Widom 线(Widom line, WL)不属于一级相变, 仍可将SCF 区域划分为类液(liquidlike, LL)区和类气(gas-like, GL)区.Ha 等[11]引入机器学习方法对SCF 分子进行标记, 任一分子要么标记为LL, 要么标记为GL, 并发现存在Widom三角形, 三角形内是LL 和GL 的混合物, 即类沸腾区.至此, 可以确定SCF 在分子量级上确实可分为GL 和LL 两个不同的区域.2015 年, Banuti[12]采用理论分析方法, 证明了存在类沸腾现象, 并给出了临界点以上气液共存的拟沸腾曲线.最近, Maxim 等[13]在等压加热过程中, 通过中子成像技术成功地监测到了超临界水通过Widom 线时, 密度发生LL 和GL 间的转变, 并且这个过程和Banuti提出的拟沸腾理论[12]一致.本课题组采用超临界拟沸腾理论进行了一系列研究[14−16], Zhu 等[14]在内径为10 mm 的光滑圆管内, 对超临界CO2竖直向上流动传热过程进行了实验研究, 发现当超临界拟沸腾数SBO 超过5.126 × 10–4时, 发生传热恶化, Xu 等[15]进一步将这个结论推广到超临界H2O,R134 a 和R22, 不同的流体的临界值SBO 不一样.张海松等[16]通过理论分析了浮升力和加速度对SCF 传热的影响, 并通过实验数据验证了经典的衡量浮升力和加速度准则, 发现传统的经典湍流理论采用太多假设, 忽略了SCF 的真实气体效应,提出的准则与实验数据相差较大.虽然, 超临界拟沸腾或拟相变现象在数值、理论及实验方面均取得了较大进展, 但是, 还远远没有形成完整的理论体系, 可以用来解决工程实际问题, 因此, 需要进一步深入的研究.

本文假设超临界流体为不均匀的异质结构, 采用拟沸腾理论处理超临界流体的传热恶化现象, 通过量纲分析得出了两个具有明确物理意义的新无量纲数, 表征了拟沸腾如何诱导传热恶化发生.超临界CO2实验数据与理论分析结果吻合较好, 实验数据来自华北电力大学低品位能源多相流与传热北京市重点实验室.最后, 基于无量纲准则分析,提出了一个拟沸腾诱导超临界流体传热恶化机制,为超临界拟沸腾传热现象研究提供了理论依据.

2 超临界流体物性及拟沸腾理论

超临界流体的GL 和LL 转换发生在气液共存线的延长线上[12], 这条线被称为WL, 这个过程几乎是不连续的, 图1(a)展示了不同工质在超临界和亚临界压力下跨过饱和曲线和WL 时的相图.当流体跨过WL 或拟临界温度Tpc时, 流体状态发生急剧改变, 当温度低于Tpc时, 流体处于LL 态,当温度高于Tpc时, 流体处于GL 态.这两种状态在亚临界压力下分别对应液体和蒸气, 超临界压力下的拟沸腾温度或拟相变温度定义为Tpc, 相当于亚临界压力下的饱和温度Ts.WL 的位置由3 个热力学响应函数的极值确定, 即比热容cp、等温压缩系数KT和膨胀系数βp:

其中, ρ 和i 分别是流体的密度和焓值, 下标p 和T 分别代表为压力和温度.当发生LL 与GL 转变时, 比热容、等温压缩系数和膨胀系数在拟临界温度处达到最大值.当Ma 数足够低和加热速率较高时, 对整个加热系统有dp/p ≪ dT/T, 压力对密度的影响相对较小.拟临界点处的比热容和膨胀系数对能量和动量传递过程中的影响相对较大.不同压力下, 比热容、膨胀系数和等温压缩系数随温度变化趋势, 如图1(b)—(d)所示, 其中Tc为临界温度.当流体穿过WL 时, 比定压热容、膨胀系数和等温压缩系数出现峰值, 随着压力升高, 峰值逐渐减小, 当压力接近p=3pc(临界压力pc=7.38 MPa)时, 峰值逐渐变得平缓, 这表明压力升高物性畸变变小.

图1 不同压力下的SCF 穿越WL 的物性变化Fig.1.The physical properties of SCF crossing WL under different pressures.

对于超临界流体, 就区分LL 和GL 特性而言,最重要的是密度变化, 又称为高密度和低密度流体.不同压力下, CO2的密度随温度分布, 如图2所示, 对于亚临界压力, 随着温度升高, 密度逐渐降低, 在相应压力的饱和温度下, CO2由液态向气态转变的过程是不连续的.但是, 当压力超过临界点, 随着温度升高, 超临界CO2的LL 和GL 在WL处转化是连续的.虽然, SCF 跨过WL 不属于一级相变, 但是, 在一个有限的温度区间内, 剧烈的密度变化, 使其和亚临界流体跨过饱和温度的密度变化仍高度相似.因此, 在超临界压力下, LL 和GL的转换可认为发生在1 个有限的温度区间内, 这个过程和亚临界气、液相变过程类似.

图2 不同压力下的CO2 密度随温度变化分布Fig.2.Density of CO2 varies with temperature under different pressures.

对亚临界气液加热系统, 流体在饱和温度下的气、液界面上吸收汽化潜热, 由液态转化为气态, 这个过程是等温的.对于超临界压力下, 根据Banuti[12]的拟沸腾理论, 流体由LL 转变为GL 时, 这个过程不仅流体的微观结构发生变化, 而且必须伴随着流体的温度升高.因此, SCF 的LL 和GL 转变发生在1 个低于Tpc的温度T–和高于Tpc的温度T+区间内.

从热力学看, 可将SCF 按温度分为3 个区域,当 TbT+时, 流体为单一的GL 态; 当 T−≤Tb≤T+时, 流体发生了LL 态和GL 态间的转变, 这个温度区间是LL 和GL 共存状态.因此, SCF 相变焓Δi 定义为

对于压力为10 MPa 的超临界CO2, 拟沸腾相变焓Δi 和跨过的温度区间的确定方法[12], 如图3(a)和图3(b)所示, 类液态比热容cp,l为温度3Tpc/4 处对应的比热容, 如图3(a)中的点对应的比热容为类气态比热容, 由(5)式确定:

其中: 下标g 和l 分别表示类气和类液; γ 是比热容比(γ=1.15); R 为气体常数, R=Rg/M, Rg气体分子常数, 其值为8.31451 J/mol·K; M 为流体分子质量, 对于CO2, 其值为44 g/mol.那么, 确定类气流体焓ig直线方程为

确定类液流体焓的直线方程的方法(涉及的温度单位均是℃), 如图3(b)所示.点 B¯ 对应温度为3Tpc/4 处的焓值 iB¯, 在点 B¯ 处的曲线斜率为cp, l,与图3(a)中的点重合, 因此, 过点B 且斜率为cp, l的直线方程可以表达为

其中下标pc 表示拟临界.红色曲线为过点(Tpc,ipc)的切线方程(见图3(b)), 这个切线分别和确定类液流体焓和类气流体焓的直线方程存在1 个交点A 和B, 因此, 可分别确定T–和T+的值, 再根据4)式可以确定SCF 的拟相变焓Δi.根据这个方法可确定出超临界CO2在不同压力下的 T−和T+, 如图3(c)所示.

3 超临界流体管内竖直向上加热流动换热特性

3.1 不同参数对超临界CO2 壁温分布的影响

图4(a)—(d)给出了不同参数下, 内壁温Twi随焓值ib的变化分布情况.当压力为7.995 MPa、质量流速为520 kg/(m2·s)和管径为8 mm 时, 一旦传热恶化(heat transfer deterioration, HTD)发生, 壁温先是上升, 然后出现1 个明显的峰值, 然后再下降, 随着热流密度逐渐增大, 壁温飞升值越大, 并且恶化点向低焓值区间移动, 这个过程和亚临界膜态沸腾非常相似, 如图4(a)所示.如果给定其他参数(见图4 中蓝色字体对应的工况), 分别改变质量流速G 和压力p 及内径di, 内壁温Twi随焓值ib变化分布如图4(b)—(d)所示.增大G 或p 可以减小传热恶化程度, 甚至使恶化传热变为正常传热(normal heat transfer, NHT), 而管径对恶化起始点没有影响, 但是改变了传热恶化时的壁温飞升量级, 当传热恶化发生后, 管径减小, 恶化程度减小, 这和亚临界膜态沸腾规律类似[17].

3.2 超临界和亚临界压力下传热恶化

对于亚临界压力, 随着液体不断被加热, 液体的主流过冷度降低, 壁面及靠近壁面的液体边界层温度逐渐升高, 故壁面上的汽化核心增多, 以至于气泡布满壁面, 如果气泡不从壁面跃离, 而是在近壁处形成气泡层, 这会阻碍主流液体润湿壁面, 而气泡下面的液体快速蒸干形成气膜, 最终导致传热恶化.对于超临界压力, 在高热流qw和低质量流速G 等特定条件下, 壁面温度会出现峰值, 如图5中的红色曲线, 黑色曲线为亚临界压力下偏离核态沸腾(departure from nucleate boiling, DNB), 这两个压力下的传热恶化时的壁温分布非常类似, 均出现了严重的壁温飞升, 随后再恢复的过程.这种现象在实际工程中应该避免, 但还没有可靠的理论来解释和预测超临界压力下的传热恶化.

图4 不同参数下CO2 内壁温随焓值变化分布Fig.4.The inner wall temperature distribution of S-CO2 with enthalpy under varies parameters.

图5 亚临界DNB (黑色)和超临界压力下传热恶化壁温(红色)分布对比Fig.5.Comparison of wall temperature distribution subcritical DNB (black) and supercritical heat transfer deterioration (red).

为进一步解释这两种现象的相似性, 给出了公开文献中R134a[18]、R22[19−20]和H2O[21−22]在超临界压力下的内壁温Twi、主流温度Tb及拟相变温度Tpc相对大小的关系, Twi/Tpc和Tb/Tpc随焓值的分布关系, 如图6(a)—(f)所示.当SBO 相对较小时, 不同的SCF 无量纲内壁温Twi/Tpc随焓值变化平稳, 壁面没有明显的峰值, 但当SBO 相对较大时, 发现3 种工质均在Twi/Tpc略微大于1,而Tb/Tpc小于1 时, 发生传热恶化, 甚至会出现多次恶化现象.这是因为当Twi/Tpc略微大于1 时, 拟临界点刚好出现在近壁区, 根据我们之前的研究, 当SBO 分别超过各自的临界值时[15], 类气膜(又称为低密度层)会在近壁区充分膨胀, 导致其厚度增加, 从而导致传热恶化, 而当SBO 相对较小时, 类气膜不能充分膨胀, 传热没有明显的恶化.故当Tb< Tpc< Twi时, LL 开始向GL 转化, 这个过程导致了超临界和亚临界传热恶化时的壁温分布具有相似性.

图6 不同超临界流体在正常传热和恶化传热过程中的Twi/Tpc 和Tb/Tpc 随焓值变化分布Fig.6.Distribution of Twi/Tpc and Tb/Tpc with enthalpy during normal heat transfer (NHT) and heat transfer deterioration(HTD)of different supercritical fluids.

4 量纲分析

如果做如下假设: 1) 假定超临界流体类气和类液转化是不均匀的, 这个过程和亚临界压力相变过程中的液体和气体转化类似, 超临界流体的类气和类液转化发生在拟相变温度Tpc处, 这样超临界类气和类液的物性参数都具有明确的定义; 2) 超临界流体的相变焓Δi 由Banuti[12]的拟沸腾理论确定, 当拟相变过程发生在1 个温度区间内时, 这个过程和亚临界沸腾类似.

对于竖直向上加热的管内SCF 流动传热过程,选择以下变量: 1) 描述系统的宏观参数, 即: 热流密度qw、质量流速G、超临界进口欠热度ipc–iin、重力加速度g、径向长度di、轴向长度L; 2) 影响系统流动与传热的物性参数, 即: 比热容cp, l和cp,g、密度ρl和ρg、导热系数λl和λg、黏度µl和µg、体膨胀系数βl和βg以及拟相变焓Δi.所选取的变量没有包含压力和温度, 这是因为压力和温度均是通过物性影响系统的流动与传热过程, 这些物性参数已经包含在所选取的变量中.如果忽略热能和机械能之间的转化, 那么, 合理地假设在SCF 流动传热过程中, 存在如下关系: f (qw, G, Δi, g, di, L, cp,l, cp,g,ρl, ρg, λl, λg, µl, µg, βl, βg) = 0.这是关于管壁与流体之间的热传递问题, 这里选择基本量纲为M、L、τ、T, 上述16 个有量纲数很容易构造出7 个无量纲数: L/di, cp, l/cp, g, ρl/ρg, λl/λg, µl/µg, βl/βg以及(ipc–iin)/Δi, 其中无量纲数(ipc–iin)/Δi 类似于亚临界中的过冷度数, 这个无量纲数表示在超临界压力下类液流体温度低于相应压力下拟临界温度差值.这里暂时剔除几何长度L 和物性参数cp,g,ρg, λg, µg, ipc–iin, 此时, 将上述物理量的量纲列成矩阵如下, 每1 个纵列表示1 个变量的量纲矢量:

将上述量纲矢量表示成矩阵的形式有:

这是1 个4 行10 列的矩阵, 运用几何方法很容易证明向量组是线性无关的, 并且这个矩阵的矩阵秩为4, 组成6 个无量纲数π1, π2, …, π6, 因此,在完整的集合中共得到13 个无量纲数, 根据π 定理, 表达成无量纲方程的形式如下:

进一步可得:

实际上, 从上述的物理量可得出的无量纲数远比表中的多, 但是Buckingham[23]研究表明, 合理的选择无量纲数, 剔除无用的无量纲数, 对实际问题分析更有益处, 这里得到的无量纲数组主要是基于以下考虑: 1) 无量纲数包含的独立变量容易被实验控制, 尽量避免同1 个变量出现在不同的无量纲数组中, 可最大限度的优化实验; 2) 保留了经典的无量纲数, 包括雷诺数、普朗特数及弗劳德-雷诺数等, 这些无量纲数在大量的实验数据中得到验证, 物理意义已经非常明确, 具有普遍的适用性;3) 这里对传热恶化分析是类比亚临界流动沸腾过程中的过冷沸腾, 认为SCF 发生传热恶化主要是由于类气膜或低密度层覆盖壁面引起的, 因此, 管内传热主要机理由类气膜导热和类液体对流组成;4) 所选择的无量纲数包含的物理量应受实验影响产生的误差尽量小.

在上述无量纲数中, π13这个无量纲数表征了近壁区的类气膜内的温度梯度大小的无量纲数, 当这个无量纲数较大时, 表明近壁区的温度梯度较大.从物性来说, 当发生拟相变过程中的温度梯度足够大时, SCF 密度分布更接近亚临界压力下的相变过程中的密度分布(见图2), 从而产生了类似于亚临界相变过程中的表面张力作用, 为了更好地表达近壁区附近类气膜的温度梯度对类气膜膨胀过程的影响, 将π13表达为

其中, 类气膜密度ρg的定性温度为(Tb+Twi)/2,类液密度ρl的定性温度为Tb, 如果类气膜的厚度很小, 类气膜的密度定性温度也可采用Tw估计,βpc为拟临界温度下的膨胀系数.在图7 中, 红色代表正常传热, 黑色代表恶化传热, 当热流密度为239.1 kW/m2时, 传热发生明显恶化时(见图7(a)),对应的类气膜内的温度梯度也同样存在1 个先增大后减小的趋势(见图7(b)), 这表明在恶化时, 近壁区物性剧烈变化, 较大的温度梯度使类气膜覆盖在壁面上, 热量集聚在近壁区.但是, 对于低热流密度的正常传热, 类气膜内的温度梯度相对较小, 而是随着焓值缓慢地增大, 传热没有明显的恶化现象.

图7 正常和恶化传热下类气膜内的温度梯度和内壁温随焓值分布Fig.7.Distribution of temperature gradients and inner wall temperature with enthalpy in gas-like film under normal heat transfer(NHT) and heat transfer deterioration(HTD).

如果对π1和π6简单的组合, 可以得到1 个新的无量数:

这个无量纲数表征了类气膜的径向生长速度qw/(ρg·Δi)和主流流体的平均速度G/ρl之比.如图8(a)所示, 给出了正常传热和恶化传热下的新无量纲数π、内壁温和换热系数随焓值分布, 对于正常传热, 类气膜生长速度相对于主流流体的速度较小, 热量被及时带走, 因此, 传热没有明显的恶化现象.但是, 当类气膜生长速度相对较大时, 主流流体没有及时通过对流带走近壁区的热量, 热量集聚在近壁面, 传热恶化发生, 如图8(b)所示, 在传热恶化时, 这个无量纲数同样出现1 个峰值.

图8 正常传热和恶化传热下的新无量纲数π、内壁温Twi 和换热系数h 随焓值分布Fig.8.Distribution of the new dimensionless numberπ, inner wall temperature Twi and heat transfer coefficient h with enthalpy under normal heat transfer (NHT) and deteriorated heat transfer (HTD).

5 新无量纲数描述超临界流体拟沸腾作用机制

采用新无量纲数描述拟沸腾诱导SCF 传热恶化过程, 如图9(a)所示, 当传热发生明显的恶化时, 管内类气膜膨胀变厚, 在近壁区形成1 个凸起.如图9(b)所示, 当传热恶化时, 较厚的热类气膜覆盖在壁面, 管中心为相对较冷的类液流体, 较冷的类液流体不能充分润湿壁面导致传热恶化, 当类气膜变薄时, 传热又发生恢复, 类气膜越厚, 传热恶化越严重.图9(c)给出了新无量纲数π 和π13在管内轴向距离的局部分布, 较大的π 表明较大的类气膜生长速度, 热量没有被及时带走, 而较大的类气膜温度梯度使类气膜覆盖在近壁区, 在峰值附近出现了1 个峰值, 这个值越大, 流动越接近亚临界过冷沸腾过程, 故传热恶化, 这两个无量纲数一个表征类气膜生长速度与主流平均速度大小, 一个表征类气膜内的温度梯度大小.

图9 新无量纲数作用超临界流体拟沸腾传热恶化机制Fig.9.New dimensionless action on number deterioration of pseudo-boiling heat transfer mechanism.

6 结 论

由于SCF 物性剧烈的非线性变化, 导致其流动传热特性异常复杂, 如何更好的处理其内部温度与速度场之间的关系方面, 许多研究者对其进行了诸多探索[24−26], 这极大推动了SCF 理论和实验方面的进展, 这些研究大都是从理想的单相流体概念出发, 忽略了SCF 的实际气体效应.目前, 关于SCF的异常传热特性研究已进入了瓶颈期, 尤其是对传热恶化的预测, 在没有明确其内部的物理过程的情况下, 难以有进一步的较大突破.因此, 如何利用壁面温度、主流温度及已知参数等少量信息最大限度地表征其内部物理过程, 仍然是目前的主要任务.

本文放弃了传统的理想单相流体概念, 依据拟沸腾理论, 从经典的量纲分析和亚临界过冷沸腾理论模型出发, 对SCF 传热恶化进行了研究, 得出了2 个新无量纲数, 即π = (qwρl)/(G Δi ρg)和π13= (qw·βpc·di)/λg.前者表征了超临界类气膜生长速度与流体主流流速之比, 这个比值越大, 类气膜在近壁区生长的越快, 大量的热量聚集在近壁区, 而后者表征了类气膜内的温度梯度(qw·βpc·di)/λg,这个数越大, 表明流动越接近亚临界过冷沸腾, 较大的温度梯度使类气膜覆盖在壁面, 最终导致传热恶化.实验结果表明, 超临界CO2传热恶化壁温和亚临界过冷沸腾恶化壁温相似, 与拟临界温度Tpc相关的拟沸腾现象能够较好的解释二者的相似性.对于超临界CO2竖直向上加热流动过程中, 增大压力或减小质量流速能够减小传热恶化程度, 甚至能够消除传热恶化.在目前的实验中, 发现管径对传热恶化起始点没有影响, 但是, 对传热恶化的壁温飞升程度有明显的影响, 管径越大, 恶化越严重.根据拟沸腾理论研究了其内部流动传热场的特性与已知参数和物性之间的联系和规律, 量纲分析得到的结果可用于不同流体间的相似理论分析, 这对推动SCF 传热理论研究具有一定意义.