竖直螺旋管中超临界RP-3 航空煤油换热数值研究

王彦红，陆英楠，李洪伟，*，李素芬，东明

(1.东北电力大学 能源与动力工程学院，吉林 132012；2.大连理工大学 能源与动力学院，大连 116024)

空-油换热器可以实现航空发动机有效的热防护[1]。通道内超临界航空煤油出现变热物性换热问题。空-油换热器广泛采用螺旋管结构，离心力导致通道换热周向不均匀，加剧了换热的复杂度[2]。需要开展螺旋管中超临界航空煤油的换热特性研究。

国内外学者对螺旋管中超临界水和二氧化碳的换热开展了大量的研究[3-10]。Zhao 等[3]研究表明，离心力和浮升力共同引起水的二次流，螺旋管出现周向非均匀换热特征，拟临界区高的热导率削弱了非均匀换热问题。Li 和Bai[4]发现半周加热螺旋管内水的换热强于全周加热情况，因为半周加热二次流增强了流速和湍流能。Zhang 等[5]通过边界层理论和二次流特征分析了浮升力和离心力对螺旋管内二氧化碳换热的影响，建立了适用于螺旋管的浮升力参数和准则。Wang 等[6]进行了螺旋管内二氧化碳的换热研究，管截面流速和湍动能具有斜向对称轴，浮升力起到强化换热的作用。Xu 等[7]通过内侧-外侧和顶部-底部2 个径向位置温度、流速、湍动能、动力黏度的分布情况，阐述了周向的非均匀换热机制。Liu 等[8]考察了倾斜螺旋管中二氧化碳的换热特征，探讨了倾角变化对离心力和浮升力的合力作用及其对周向换热的影响。李洪瑞等[9]研究表明，竖直流动不同流向时换热差别可以忽略，浮升力引起水平流动换热的振荡问题。黄小锐等[10]探究了运行参数对二氧化碳换热的影响机制。

超临界航空煤油换热研究集中于圆截面直通道的情况。浮升力影响竖直通道内边界层流体微团交混作用，抑制上升流动的湍流强度，增强下降流动的湍流强度[11]。水平通道煤油密度产生非规则分层，出现二次流和周向非均匀换热问题[12]。Wen 等[13]实验研究了竖直螺旋管内RP-3 航空煤油的换热机制，离心力致使管外侧换热显著优于管内侧。因为实验研究不能揭示温度场和流场的详细信息，不能深入地阐述离心力的作用机制，需要通过数值方法解决这些问题。

本文对竖直螺旋圆管内超临界压力RP-3 航空煤油的换热开展了数值研究，探究了换热特性和换热机理，阐述了运行参数和结构参数对换热的影响机制，及其对二次流强度的影响。基于传统换热公式，通过误差分析得到了适用于航空煤油的换热关联式，研究成果可为空-油换热器设计提供一定的指导。

1 数值模型与数值方法

1.1 物理模型与网格方案

图1 为竖直螺旋管物理模型。圆管外径do为3 mm，内径di为2 mm。螺旋管结构参数包括绕径D、螺距P、加热段长度L。进、出口的绝热段长度均为150 mm，用以避免进、出口效应的影响。加热段外表面施加均匀的热流密度。流体进口边界给定质量流速和进口温度，出口给定静压，进口和出口的圆环为绝热壁面，流固之间为耦合面。g为重力加速度，完成上升流动工况。取管顶部（φ=0°）、管外侧（φ=90°）、管底部（φ=180°）、管内侧（φ=270°）这4 个周向位置开展换热分析。

图1 竖直螺旋管物理模型Fig.1 Physical model of vertical helical tube

不同通道结构的圆管截面尺寸相同，长度差别不大，制定相同的网格划分方案。固体域径向等间隔划分了10 层网格。流体域近壁进行网格加密，第1 层网格y+<1，前3 层网格y+≤5。通过网格无关性分析，3 200×1 200（管截面网格数量×流动方向网格数量）的网格方案满足网格无关性要求。网格划分情况如图2 所示。

图2 螺旋管网格Fig.2 Meshes of helical tube

1.2 控制方程

质量守恒方程：

式中：ρ为密度；u为流速。

动量守恒方程：

式中：µe为有效黏度；p为压力；i、j和m表示3 个坐标方向。

能量守恒方程：

式中：cp为定压比热容；T为温度；Ke为有效热导率；δij为克罗内克符号。

RNGk-ε湍流模型：

式中：k为湍动能；ε为耗散率；ak和aε为湍流普朗特数；C1ε=1.42；C2ε=1.68；C3ε=0.084 5；Gk为剪切力项；Gb为浮升力项；Rε为附加项。该湍流模型更适用于螺旋管中超临界流体换热模拟，模型中的附加项提高了螺旋管内快速应变流动预测的准确性[14]，结合增强壁面处理可进一步提高湍流换热的预测精度。

热传导方程：

式中：K为固壁热导率。

通过二阶迎风差分格式离散控制方程，进一步选择双精度分离求解器、SIMPLEC 算法、隐式Gauss-Seidel 处理，控制方程的收敛标准设定为10−6。

1.3 热物性参数

RP-3 航空煤油（临界压力为2.34 MPa[15]、临界温度为645 K[15]）在超临界条件下的密度、定压比热容、热导率和动力黏度通过分段线性函数加入Fluent 14.5，拟临界点附近增加线性分段数以提高热物性计算精度。热物性参数选取实验测量数据[11]和燃料3 组分替代模型计算数据[16]，即以温度800 K为界，低于该温度采用实验测量数据，高于该温度采用替代模型数据。图3 为2 种压力下定压比热容随温度的变化情况。

图3 RP-3 航空煤油比热容随温度的变化情况Fig.3 Specific heat capacity variation with temperature of RP-3 aviation kerosene

1.4 模型验证

基于竖直螺旋圆管中RP-3 航空煤油的换热实验数据[13]做数值模型验证，如图4 所示。螺旋管参数：do=2.2 mm，di=1.82 mm，D=20 mm，P=10 mm，L=1 200 mm，n=8。l为局部加热长度，管内壁温度Twi为周向平均值。质量流速G设定为1 572 kg/(m2·s)，进口温度Tin为473 K，运行压力p为5 MPa，壁面热流密度q范围为360～600 kW/m2，向上流动。由图4可知，数值计算管内壁温度Twi的变化特征与实验数据符合良好，相对偏差落在±6.25%的范围，说明数值模型有效且合理。随着压力降低，拟临界区热物性奇异性更强，不会影响湍流模型的可靠性和预测精度[12]。

图4 模型验证Fig.4 Model validation

2 数值结果与分析

2.1 运行压力的影响

探讨不同运行压力下的换热特性，运行参数：Tin=400 K，G=1 200 kg/(m2·s)，q=600 kW/m2。材料热导率为20 W/(m·K)。螺旋管结构参数：D=30 mm，P=10 mm，L=758 mm，共8 圈。

图5 为不同压力下管内壁温度Twi和换热系数h（均为周向平均值Twi,av、hav）沿流动方向的变化情况。由图5 可知，压力为3 MPa 时，管内壁温度沿流动方向不断增大，l/di=120 位置管壁温度出现局部峰值，说明该管段出现了传热恶化现象。因为边界层流体处于拟临界点（定压比热容峰值对应的温度，3 MPa 下约为672 K）附近时，煤油密度随温度剧烈改变且热导率为最小值，导热性能和流动性能差的流体层覆盖壁面，阻碍热量由壁面向主流传递，起到弱化换热的作用。换热系数出现局部谷值，验证了传热恶化问题。管下游换热系数出现峰值，这与定压比热容的变化规律相同，即定压比热容峰值是该强化换热的机制。压力增加管内壁温度沿流动方向近似线性增加，局部传热恶化没有出现，换热仅为管下游的强化换热机制。因为高压力下拟临界温度提高，密度增大且变化平缓，热导率增大，传热弱化被削弱。同时，压力增加时壁温上升斜率增大，运行压力提高致使前段换热增强，而后段换热减弱。

图5 不同压力下Twi,av 和hav 沿流动方向的分布情况Fig.5 Twi,av and hav distributions along flow direction at various pressures

图6 为不同压力下Re沿流动方向的变化情况。由图6 可知，Re峰值与换热系数峰值相对应，均出现在管道下游，说明高Re有利于强化换热。随着压力增加，热物性变缓，Re降低，致使换热系数下降。

图6 不同压力下Re 沿流动方向的分布情况Fig.6 Re distribution along flow direction at various pressures

图7 为2 种压力下管内壁温度和换热系数（l/di=160）沿管周向的变化情况。由图7 可知，两者沿管周向不是均匀分布，而是出现周向非均匀换热问题。管外侧（φ=90°）位置管壁温度最低，换热系数最大。管内侧（φ=270°）管壁温度最高，换热系数最小。管顶部（φ=0°）和管底部（φ=180°）的管壁温度和换热系数差别不大。这说明离心力导致了管内侧和管外侧的换热差别，浮升力对管顶部和管底部的换热影响可以忽略。高压力下周向的换热差别减弱。

图7 不同压力下Twi 和h 沿管周向的分布情况Fig.7 Twi and h distributions along circumferential direction at various pressures

图8 为2 种压力下管截面的温度分布情况。由图8 可知，离心力导致流体温度横向分层，温度等值线呈月牙型，高温流体趋向于管内侧，低温流体趋向于管外侧。流体温度异常分层导致流体密度周向不均匀，周向密度梯度引起周向不平衡动能，进而出现自然对流问题。流体温度异常分层对固体域起到不同程度冷却作用。管内侧受到高温流体冷却，冷却作用差，固壁温度较高，而管外侧受到低温流体冷却，冷却作用好，固壁温度较低，固体域出现横向非规则分层。因此，固体域出现非均匀的热传导过程，内表面热流密度重新分配，高温区热流密度减小，更多热流由其他位置吸收利用。压力改变造成管截面温度分布差别，低压力下流体域温度异常分层更显著，高压力下固体域温度异常分层更突出。

图8 不同压力下管截面的温度分布情况Fig.8 Temperature distribution in tube cross section at various pressures

图9 为2 种压力下管截面二次流分布情况。二次流流速usec=((ux)2+(uy)2)0.5，取平行于xy面的管截面。可以看到，周向不平衡动能致使流体截面出现二次流，管顶部和管底部各有1 个涡，2 个涡不完全对称。二次流流速最大值位于管顶部和管底部的壁面附近，数值沿流动方向不断增大，说明二次流沿流动方向逐渐增强。二次流驱动热流体沿管周向从管外侧流向管内侧，再从中心线返回主流，周向流动差别引起周向不均匀换热问题。水平通道浮升力引起的二次流机制：管底部为低温高密度流体，管顶部为高温低密度流体，高密度流体沿管壁从管底部流向管顶部。而螺旋通道离心力引起的二次流机制：管内侧为高温低密度流体，管外侧为低温高密度流体，低密度流体沿管壁从管外侧流向管内侧。两者的作用机制迥异，即水平通道二次流增强管顶部和管底部的换热差别，而螺旋通道二次流增强管内侧和管外侧的换热差别。由于管道内径较小，浮升力作用可以忽略。因此，二次流没有出现浮升力和离心力共同作用的斜向分布情况。离心力加速度为2u2/D。高压力下密度提高，流速减小，离心力作用削弱。流体密度梯度减小，周向不平衡动能减弱，管截面具有减小的二次流速度。

图9 不同压力下管截面二次流分布情况Fig.9 Secondary flow distribution in tube cross section at various pressures

图10 为2 种压力下P1～P4位置流速和湍动能的径向分布情况。P1：l=0.189 m；P2：l=0.379 m；P3：l=0.569 m；P4：l=0.758 m。ri为管内半径，r为径向位置，r/ri=−1 为管内侧位置，r/ri=1 为管外侧位置。由图10(a)可以看到，管外侧流速高于管内侧。管内侧流体为低密度低流速，局部质量流速ρu较小，流动性能较差；管外侧流体高密度高流速，高质量流速有利于强化换热，故两侧出现换热差别。提高压力，密度增大，流速下降，两侧质量流速和换热差别减弱。图10(b)的结果表明，管外侧湍动能高于管内侧，离心力增强了管外侧湍流的交混作用，抑制了管内侧的湍流交混作用。同时，管内侧具有厚的热边界层，而管外侧具有薄的热边界层，这也是导致两侧换热差别的原因。

图10 流速和湍动能的径向分布情况Fig.10 Velocity and turbulent kinetic energy radial distribution

2.2 热质比的影响

图11 为不同热质比下平均管内壁温度和平均换热系数沿流动方向的变化情况。热质比定义为外壁热流密度与质量流速的比值。选取3 种热质比工况，考察壁面热流密度和质量流速的影响。p=3 MPa，Tin=400 K。q/G=0.4 J/g，q=600 kW/m2，G=1 500 kg/(m2·s)；q/G=0.5 J/g，q=600 kW/m2，G=1 200 kg/(m2·s)；q/G=0.6 J/g，q=720 kW/m2，G=1 200 kg/(m2·s)。螺旋管结构参数与2.1 节相同。可以看到，随着热质比提高，管壁温度整体增大，传热恶化起始位置前移。这是因为质量流速下降，流速减小，起到削弱换热的作用；而热流密度增大，冷却需求增加，对换热状况造成不利影响。

图11 不同热质比下Twi,av 和hav 沿流动方向的分布情况Fig.11 Twi,av and hav distributions along flow direction at various heat-mass ratios

图12 为不同热质比下管内壁温度和换热系数（l/di=160）沿管周向的变化情况。可以看到，热质比提高致使管内侧和管外侧的管内壁温度和换热系数差别增大。质量流速不变，壁面热流密度提高，热质比增大，主流密度减小增强，流速增加更加显著，离心力作用增强。壁面热流密度不变，质量流速提高，热质比减小，流速增大使离心力增强，主流高流速对离心力削弱作用也加剧，周向非均匀换热受到抑制。

图12 不同热质比下Twi 和h 沿管周向的分布情况Fig.12 Twi and h distributions along circumferential direction at various heat-mass ratios

2.3 管道尺寸的影响

讨论通道结构（D和P）对换热的影响，3 种通道方案：D=30 mm，P=10 mm，L=758 mm，共8 圈；D=50 mm，P=10 mm，L=787 mm，共5 圈；D=30 mm，P=20 mm，L=770 mm，共8 圈。运行参数选取：p=3 MPa，Tin=400 K，G=1 200 kg/(m2·s)，q=600 kW/m2。图13为不同通道参数下平均管内壁温度和平均换热系数沿流动方向的变化情况。由图13 可知，随着绕径和螺距增加，上游平均管内壁温度提高，下游平均换热系数下降，不同通道结构的传热恶化程度相近。

图13 不同通道参数下Twi,av 和hav 沿流动方向的分布情况Fig.13 Twi,av and hav distributions along flow direction at various channel parameters

图14 为不同通道参数下管内壁温度和换热系数（l/di=160）沿管周向的变化情况。可以看到，绕径和螺距提高均使管内侧和管外侧的管内壁温度和换热系数差别减小。这是因为绕径增大，离心力加速度减小，离心力作用减弱；而提高螺距，离心力分力减小，离心力作用减弱。

图14 不同通道参数下Twi 和h 沿沿管周向的分布情况Fig.14 Twi and h distributions along circumferential direction at various channel parameters

2.4 二次流强度

二次流强弱通过二次流强度Se[17]描述，形式如下：

式中：A为流体域管截面面积；下标b 表征主流数值。

图15 为二次流强度沿流动方向的变化情况。基 本 工 况：p=3 MPa，Tin=400 K，G=1 200 kg/(m2·s)，q=600 kW/m2。由图15 可知，螺旋管上游低温区热物性变化平缓，二次流较弱；下游高温区热物性变化显著，二次流较强。提高运行压力、降低热质比、增大绕径、增大螺距均使二次流强度减弱，提高热质比造成二次流增强。

图15 Se 沿流动方向的变化情况Fig.15 Se variation along flow direction

2.5 换热关联式

螺旋管换热关联式通常以管周向平均温度为基准，并考虑螺旋管结构参数的影响。螺旋管换热经验关联式为

Merkel 关联式[18]：

式中：Nu为努塞尔数；Pr为普朗特数。

式(8)～式(10)均是针对二氧化碳等非碳氢燃料提出的，通过偏差分析，选取适用于RP-3 航空煤油的换热关联式。

图16 为Nu沿流动方向的变化情况。工况参数与图5 相同。可以看到，3 个关联式计算的努塞尔数沿流动方向呈现了先增大后减小的变化特征。随着压力增大，关联式预测精度提高。相比其他2 个关联式，Merkel 换热式具有最优的预测性能，相对偏差处于±20%的范围。

图16 Nu 沿流动方向的变化情况Fig.16 Nu variation along flow direction

图17 为其他工况下Merkel 关联式的预测误差情况。选取图11 和图13 工况。由图17 可知，Merkel 关联式的预测偏差处于±20%的范围，证明该换热公式具有较好的精度。

图17 Nu 数值结果与预测数据的比较情况Fig.17 Comparison of Nu between numerical results and prediction data

3 结 论

1）低压力下存在局部传热恶化问题，管下游表现为强化换热机制。浮升力对管顶部和管底部的换热影响可忽略。离心力导致管内侧换热显著弱于管外侧，出现周向换热差别。

2）离心力造成管截面流体温度横向异常分层，出现周向不平衡动能，产生二次流。管内侧边界层厚、密度小、流速低、局部流量小、湍流强度弱，这是其换热弱于管外侧的原因。

3）提高压力，密度增大，流速减小，离心力减弱。提高热流密度，管内流体加速更显著，离心力增强。提高质量流速，高流速对离心力削弱增强。提高绕径，离心力加速度减小，离心力减弱。而提高螺距，离心力分力减小，离心力减弱。二次流强度随离心力作用增强而增大。

4）Merkel 关联式可以运用于螺旋管内超临界压力RP-3 航空煤油的换热预测，预测偏差处于±20%的范围。预测精度更好的换热关联式需要通过后续的实验研究加以补充。