发动机喷管内流场对流换热系数影响因素的数值分析

王明华，陈劲松

（北京航天发射技术研究所，北京100076）

0 引言

目前，结构件在燃气流场中的烧蚀及热防护受到广泛关注。燃气流热力特性可为燃气流烧蚀影响区域耐烧蚀材料选用及热防护措施选择提供直接依据。国内外学者和技术人员对飞行器飞行状态气动热CFD计算进行了研究，研究的内容主要集中在格式效应和网格效应两个方面[1-4]，对燃气流热力特性做系统研究的比较少。

工程上，很多结构件在燃气流场中遭受冲刷，燃气流对结构件表面的对流传热类似于喷管内部燃气流对于管壁的对流传热。本文通过对某液氧煤油发动机喷管内流场结构及对流换热系数进行研究，分析网格雷诺数以及恒温壁壁温对对流换热系数的影响。

1 计算模型及边界条件

1.1 计算模型

在计算喷管内流场时，假设流动是轴对称的，运用Roe-FDS格式所对应的密度基耦合求解器求解。FDS不仅具有很好的激波等间断分辨率，而且粘性分辨率也相当优秀，这是一般格式所不具备的。Roe格式是目前采用较广泛的一种格式。湍流模型采用标准k-ε模型。考虑可压缩条件下同时模拟燃烧反映动力过程、传热过程难度较大，同时难以分析主要因素对计算结果的影响，本文暂不考虑燃烧过程，计算过程按一维定常等熵流基本关系计及燃气物性参数、气体状态方程。

1.2 边界条件

控制方程的边界条件为：在推力室进口截面设定压力进口边界条件，喷口出口设定压力出口边界条件，壁面无滑移、非穿透。由对流传热热流密度计算公式可知[5-6]，热流密度与对流换热系数及恢复温度有关，同时本文不考虑热传导和辐射传热，为了简化计算取壁面为恒温壁。开发的喷管内流场网格模型如图1所示：ab为喷管入口边界，给定燃气的总温和总压，流动方向与入口边界相垂直；bod为喷管轴线，取轴线边界条件；cd为喷管出口边界，取压力出口边界条件。

图1 某液氧煤油发动机喷管网格模型Fig.1 Grid model of nozzle in a liquid oxygen-kerosene rocket engine

根据研究需要，基于边界层理论及网格剖分方法，采用四种网格划分方法。边界层理论中控制对流传热计算精度主要受边界层厚度、边界层网格缩放规则、邻近壁面网格法向厚度等三个方面因素影响。本文采用的边界层厚度预估公式如下：

在来流气流参数一定的情况下，邻近壁面网格法向厚度n在气动领域常以另一个术语形式体现，即网格雷诺数。其定义见式 (2)[1]：

公式 (1)和 (2)中：ρ∞，v∞和μ∞分别为燃气流密度、轴向速度以及动力粘性系数；R，Cf和Tw分别为气体常数、当地摩擦系数和壁面温度；y+是无量纲参数，用于反映并控制近壁面网格质量，本文计算控制y+不超过10。

为综合考虑网格雷诺数效应，邻近壁面网格法向厚度n分别取n=0.1 mm，n=0.05 mm，n= 0.01 mm，n=0.005 mm，来流速度利用进口压力边界条件控制，进口压力分别取6 MPa，9 MPa，12 MPa和18 MPa。为有效控制网格精度及数量，边界层法向网格按线性比例1.2缩放，边界层法向及纵向网格比例控制在1：100以内。

2 计算基本结果

2.1 喷管内流场结构分布规律

基于上述边界层网格模型开发及边界条件设置方法，完成所有工况计算。图2显示了四种边界层网格模型条件下沿轴线静压、静温、马赫数以及壁面y+值的变化曲线。

新兴的工作重塑理论还有很多未解之谜，为今后的研究留下了广阔的空间。对本土学者而言，只有扎根中国文化和管理实践，结合质性与定量研究等多种方法，对各类行业和职业劳动者的工作重塑行为及其影响因素进行深入探索和细致分析，才能为构建中国本土的工作重塑理论做出贡献。

从图2可以看出：静压、静温、马赫数的分布规律与理论计算值相比，有很好的一致性，四种边界层网格划分情况下，静温的仿真结果与理论计算结果相比最大误差不超过14%；四种计算模型壁面y+的变化趋势基本相同，y+数值随着近壁面法向第一层网格的高度增大而增大。当n取0.005 mm时，整个喷管壁面y+的值在10以下。数值仿真结果与理论计算结果之间的差异主要在喉部及出口处，原因是喷管型面是轴对称的，在进行理论计算时采用了一维等熵流关系式，与实际的二维流动计算结果有所差别。

图2 沿轴线方向的静压、静温、马赫数、壁面y+分布曲线Fig.2 Distribution curves of static pressure,static temperature,Mach number and y+of wall surface in axial direction

2.2 喷管内部燃气流对壁面的对流传热规律

基于n取0.005 mm，总压18 MPa条件说明对流传热基本规律。

图3给出了喷管壁面对流换热系数的分布曲线。从图中可以看出，在喷管入口段对流换热系数变化缓慢，靠近喉部区域时迅速上升，在喉部对流换热系数达到最大值。此后对流换热系数急剧下降，最后趋于平缓。

从图3可以看出，除了在喷管入口段以外，其余位置数值模拟结果与Bartz公式预示结果一致性较好。入口段，数值模拟结果比Bartz公式预示结果大。原因是数值模拟过程中，通过对入口湍流强度的设置，限制了边界层的厚度，而边界层越薄，对流换热系数越大。同时在喷管入口段边界层并未完全发展起来。Bartz公式本身则是在大量数据基础上拟合得到的经验公式[7]，考虑了燃烧非平衡性等因素，但Bartz公式本身与试验数据之间也存在差异[8]。因此数值模拟和Bartz公式两者的计算结果在喷管入口段有较大差异。

图3 喷管壁面对流换热系数分布曲线Fig.3 Convective heat transfer coefficient of wall surface in nozzle

2.3 网格雷诺数对对流换热系数的影响

由网格雷诺数的定义可知，影响网格雷诺数大小的因素主要有近壁面第一层网格的法向高度以及来流速度。因此，分别从近壁面第一层网格的法向高度和入口总压两个方面研究网格雷诺数对计算对流换热系数的影响。

1）近壁面第一层网格的法向高度对对流换热系数的影响

由图4可以看出，近壁面第一层网格的法向高度对对流换热系数有重要影响。对于整个喷管而言，n值的变化并不影响壁面对流换热系数的整体变化趋势，但是对于喷管壁面的某一位置来说，对流换热系数的值随着n的减小而增大。

2）来流雷诺数对对流换热系数的影响

为分析来流雷诺数对对流换热系数的影响，对入口压强进行了不同设置。基于同样的网格模型，分别计算入口总压分别为6 MPa，9 MPa，12 MPa和18 MPa时喷管内部燃气流对壁面对流换热系数的影响。

图4 n=0.1 mm,n=0.05 mm,n=0.01 mm,n=0.005 mm时壁面对流换热系数Fig.4 Convective heat transfer coefficient of wall surface when n=0.1 mm,0.05 mm,0.01 mm and 0.005 mm

图5 不同入口压强条件下对壁面的对流换热系数Fig.5 Convective heat transfer coefficient of wall surface under conditions of different inlet pressure

图5给出了不同入口压强条件下，燃气流对壁面的对流换热系数的分布曲线。可以看出，对流换热系数的大小与入口压强有关。入口段和喉部附近的区域其壁面的对流换热系数随着入口压强的增大有着明显的增大趋势；出口段壁面的对流换热系数在不同入口压强条件下相差不大。就沿喷管轴向的某一截面而言，入口压强越大，该截面处燃气流对壁面的对流换热系数的值也越大，但其值的大小并不随压强的线性增大而呈现出线性比例关系。由对流传热热流密度公式可知，热流密度的数量级为MW。

2.4 壁温对对流换热系数的影响

恒温壁的壁温也是影响燃气对壁面对流传热量的因素之一。根据上述分析，选择n=0.005 mm的模型，将壁温分别设为300 K和800 K，研究壁温对对流换热系数的影响。两种壁温条件下，静温等值线图、喷管壁面y+的分布曲线以及壁面对流换热系数的分布规律如图6所示。

由图6可以看出，两个工况的沿轴向的静温、壁面y+的差值很小，说明在燃气流的温度与壁温相差较大的情况下，当采用同样的计算模型而只改变恒温壁壁温时，壁温对流场结构影响很小。由图6(c)可以看出，壁温800 K与壁温300 K时计算出的燃气流对壁面的对流换热系数整体趋势一致，壁温为800K时的对流换热系数较壁温为300K时有所减小，但差值很小，最大相对误差不超过5.5%。

图6 壁温为300 K和800 K时沿轴线方向的静温、壁面y+、壁面对流换热系数分布曲线Fig.6 Distribution curves of static temperature,wall surface y+and convective heat transfer coefficient in axial direction when wall temperature is 300 K and 800 K

2.5 仿真计算结果与Bartz公式计算结果比较

工程上常采用Bartz公式来估算喷管内部对流换热系数[9-10]，

式中：dt为喉部直径；μ为动力黏度；Pr为普朗特数；po为燃烧室压强；c*为燃气特征速度；Rt为喉部曲率半径；σ为修正系数；C为亚声速，取0.026，超声速取0.023。

图7比较了壁面对流换热系数的仿真计算结果与Bartz公式计算结果。从图中可以看出，仿真结果与Bartz公式计算出的对流换热系数变化趋势相同，数值上有差异。计算得到的某一位置的对流换热系数随着n的减小而变大。

从整体看，n=0.005 mm模型的计算结果与工程估算结果的拟合程度最好，但在入口段误差较大，其原因在前文已分析。其他三个模型对喷管入口段和喉部区域的计算结果与工程估算值相比较小，且差值较大。而在出口段，随着n的减小，仿真结果与Bartz公式计算值的差异越来越大。在近壁面第一层网格的法向高度为10-3mm数量级条件下，对流换热系数的数值仿真结果与Bartz公式估算结果一致性较好。

图7 壁面对流换热系数与Bartz公式计算结果的对比Fig.7 Comparison of convective heat transfer coefficients and results calculated with Bartz formula

3 结论

1）网格雷诺数影响喷管内流场对流传热。在适当控制网格纵向收缩因子以及纵横比的情况下，网格雷诺数效应集中体现在控制近壁面第一层网格的法向高度上。在来流参数一定的情况下，当高度为10-3mm数量级时，对流换热系数的数值模拟结果与Bartz公式估算结果一致性较好。

2）依据本文的研究工况，在网格雷诺数合理的条件下，喷管壁面对流换热系数随着来流速度的增大而增大。

3）在量热完全气体条件下，当燃气流温度远高于恒温壁壁温时，恒温壁壁温对流场结构及壁面对流换热系数的影响较小。

[1]阎超,禹建军,李君哲.热流CFD计算中格式和网格效应若干问题研究[J].空气动力学学报,2006,24(1)：127-132.

[2]SWAMINATHAN V,SETTY Speenivasa,RAMAMADHAV B S.A method for the evaluation of temperature in rocket motor nozzles with varying heat transfer coefficients,AIAA 75-1347[R].USA：AIAA,1975.

[3]李君哲,阎超,柯伦,等.气动热CFD计算的格式效应研究[J].北京航空航天大学学报,2003,29(11)：76-79.

[4]覃文洁,胡春光,郭良平,等.近壁面网格尺寸对湍流计算的影响[J].北京理工大学学报,2006,26(5)：16-20.

[5]杨世铭,陶文铨.传热学 (第四版)[M].高等教育出版社, 2006.

[6]B K科什金.航空与火箭技术中的传热原理[M].北京：科学出版社,1963.

[7]张宏伟,陶文铨,何雅玲,等.再生冷却推力室耦合传热数值模拟[J].航空动力学报,2006,21(5)：148-154.

[8]TOMIFUMI Godai,YOSHINORI Yuzawa,KATUYA Ito,et al.An experimental investigation of heat transfer in the nozzle of high-aluminized solid rocket,NASA TT F-12 [R].[S.l.]：NASA TT,1968.