燃气射流冲击传热特性的数值模拟*

刘小军,傅德彬,牛青林,李 霞

(北京理工大学宇航学院，北京 100081)



燃气射流冲击传热特性的数值模拟*

刘小军,傅德彬,牛青林,李 霞

(北京理工大学宇航学院，北京 100081)

针对射流传热问题，利用基于RNGk-ε湍流模型的数值方法模拟了射流垂直冲击平板的流动过程，并与实验数据比较，验证了模型的可行性。在此基础上，以火箭喷管入口参数为入口条件，建立了超音速燃气射流垂直冲击平板和冲击浸没平板的计算模型，分析了不同冲击条件下努塞尔数分布规律和温度分布规律, 论述了超音速射流传热的特性及影响传热特性的因素。得到了冲击距离为 (14～18)D的努塞尔数取值范围，并表明冲击距离和射流温度是影响传热效率的关键因素；冲击距离增加，传热效率降低，冲击平板表面的射流温度越高，传热效率越高。

流体力学；射流传热；RNGk-ε湍流模型;浸没平板；努塞尔数

燃气射流与工程结构之间的冲击传热广泛存在航天发射和航天器级间分离过程中，其传热过程和传热特性对结构的设计和防护有着重要影响。因此，针对燃气射流与平板间的传热特性的研究有着重要的理论意义和工程价值。

针对射流传热过程，利用理论分析、数值计算和实验测量开展了较多研究。在实验测量方面，E.Baydar等[1]通过实验得出入口雷诺数和冲击距离与喷口直径比H/D对努塞尔数Nu的影响。Y.C.Chen等[2]通过实验，测出入口雷诺数Re=600～1 200、冲击距离与喷口直径比H/D=2～20时努塞尔数的变化情况，结果表明雷诺数增大和冲击距离减小都可以使努塞尔数升高。Q.Liu等[3]通过实验和数值模拟，观察了压力和温度在雷诺数从10 000～60 000变化时的分布规律，得出努塞尔数出现第二次峰值的位置。V.Ramanujachari等[4]通过实验测量了在超音速燃气射流冲击下平板传热特性，得到在H/D变化时沿平板径向的努塞尔数比亚音速燃气射流冲击下的努塞尔数大了一个数量级。数值计算方面，M.Behnia等[5]采用垂直壁面湍流模型 (V2F模型)和标准k-ε模型分别模拟了H/D为2和6情况下的冲击射流，计算结果与实验数据对比表明标准k-ε模型在模拟湍流流动传热时误差较大。陈庆光等[6-7]采用标准k-ε模型并结合壁面函数法，对冲击射流场的平均速度和湍动能分布进行了数值计算，结果表明数值计算能够预测湍动能分布的趋势，但对滞止点附近湍动能的预测值过高；射流入口条件对湍动能分布具有较大的影响；要提高冲击射流场总体的预测效果，需要对k-ε模型进行修正。B.Merci等[8]应用低雷诺数的三次k-ε模型得出H/D为2、6时的努塞尔数与实验数据吻合良好，尤其滞止点努塞尔数的接近程度有很大改进。许坤梅等[9]采用标准k-ε模型、RNG (renormalization group)k-ε模型结合壁面函数法和低雷诺数k-ε模型, 对半封闭圆管冲击射流流场进行了数值计算, 并将计算与测量结果进行了比较，结果表明RNGk-ε模型的结果与其他两种模型相比更接近实验值。

目前，对于射流冲击传热的研究，大都是围绕垂直冲击传热进行研究，对浸没结构的传热问题研究较少。本文中，利用数值计算方法，采用RNGk-ε湍流模型对圆管冲击射流流场进行数值计算，通过与文献[8]中实验数据比较验证计算模型的可行性，进而建立超音速燃气射流传热计算模型，对垂直冲击平板和冲击浸没平板传热进行研究分析，得出超音速射流与平板间传热的特性和规律，拟为工程应用提供参考。

1 计算模型

1.1 数学模型

三维可压缩射流冲击传热涉及气体的流动和流体与固体之间的传热过程。在本文中，气体流动过程通过求解流体力学基本方程组得出，流体和固体之间的传热通过边界条件确定。

燃气射流具有高温、高速特性，其流动过程中的可压缩性和黏性不可忽略，三维可压缩、有黏性流动的控制方程可表示为[10]：

(1)

式中：Q为守恒变矢量；Fc、Gc、Hc分别为3个坐标方向的对流通量，Fv、Gv、Hv分别为3个坐标方向的黏性通量。他们分别可表示为：

其中，应力项为：

热流量与温度关系为：

为了求解上述控制方程，采用基于涡黏性假设的RNGk-ε湍流模型模拟雷诺应力和湍流标量输运项以使方程封闭。该模型考虑了流体可压缩性对湍流的影响，其湍动能k和湍流耗散率ε的输运方程表示为：

(2)

(3)

式中：Gk为平均速度梯度所引起的湍动能生成项；Gb为浮力所引起的湍动能生成项；YM为可压缩流动中脉动扩散所引起的耗散率；Sk、Sε代表源项。

为了模拟流体与固体壁面传热特性，壁面采用第一类边界条件[11]，即恒壁温边界条件。流体流过物体表面时，流体与物体表面的对流换热可表示为：

(4)

式中：q为通过壁面的热流密度；tf和tw分别是流体温度和壁面温度；h称为表面传热系数，又称对流换热系数[12]。为获得具有一般意义的传热特性，射流冲击的传热特性采用量纲一特征量努塞尔数表示：

(5)

1.2 几何模型与计算条件

射流垂直冲击平板的几何模型如图1所示。入口为压力入口，气体经过一段充分发展的圆管或喷管，获得一定的射流速度，然后垂直冲击到平板。平板设置为等温壁面，在喷口设置一块与冲击平板相同大小的绝热壁面阻挡板，射流出口为两板之间的间隙，设置为压力出口。喷口的直径为D，喷口与平板间距为H，平板径向长为R。

超音速燃气射流冲击浸没平板的对流传热模型采用三维模型，如图2所示。浸没平板位于喷管下方，平板厚度为5 mm，浸没平板采用等温壁面的边界条件，喷管入口为压力入口，外场为压力远场。

为提高壁面附近的模拟精度，壁面附近采用近壁面函数，并通过加密壁面附近网格使得y+<1。模型采用计算流体力学软件Fluent求解计算，收敛条件设置为归一化残差收敛到10-3。

图1 射流冲击垂直平板模型Fig.1 Model of vertical plate impinged by jet flow

图2 射流冲击浸没平板模型Fig.2 Model of submerged plate impinged by jet flow

1.3 模型验证

图3 亚音速射流冲击平板的努塞尔数分布Fig.3 Nusselt number distribution in plate impingedby subsonic jet flow

为验证上述计算模型和计算网格，对亚音速气体射流垂直冲击平板的模型进行计算分析。气流经过长为20D的长直管充分发展，以13.765 m/s的速度从管口喷出，垂直冲击到平板底部，喷口直径D=0.026 mm，喷口与平板间距H=2D，平板长R=15D。计算获得的努塞尔数，与文献[8]中实验数据和立方k-ε计算结果对比，如图3所示。

从图中可以看出，采用RNGk-ε模型模拟的数据与实验数据吻合较好，最大误差出现在滞止点附近，为约10.7%。在滞止点附近， RNGk-ε模型未能准确地模拟努塞尔数出现第二次峰值的位置。U.Heck等[13]也提到RNGk-ε模型模拟近壁面射流区的传热结果很接近实验值，但是在滞止点附近产生的误差会在10%左右，原因在于RNGk-ε模型对射流扩散率的预测值比实验值高2倍[14]。尽管采用RNGk-ε湍流模型模拟气体射流冲击平板传热时有一定的误差，但对于工程应用，该误差在可以接受范围内，模型可用于传热计算分析。

2 超音速燃气射流垂直冲击平板传热问题

超音速燃气射流垂直冲击平板模型采用三维模型，气体以一定的压力和温度从直径D=0.024 3 m的喷管喷出垂直冲击到圆形平板，冲击距离H=14,16,18D，圆形平板径向R=15D。气流在喷管入口总压为4 MPa，总温为3 117 K，经过拉伐尔喷管的压缩和扩张段后，气流以超音速喷出，垂直冲击到平板上。

计算获得的努塞尔数及文献[13]实验数据如图4所示，可以看出射流冲击距离为16D的数值计算结果与实验数据吻合良好，而且证实了文献[4]中提到的，在超音速气体射流冲击下，沿平板径向的努塞尔数比亚音速气体射流冲击下的努塞尔数大一个数量级的结论。

从不同冲击高度得到的3条数值模拟数据曲线可以看出，随着冲击距离的增加，努塞尔数逐渐减小，说明了冲击距离增加，传热的效果逐渐减弱。为了便于分析，在距离冲击平板0.5 mm的截面上取沿着平板径向的气流温度曲线，如图5所示。从图5可以看出，沿着冲击平板径向的射流温度下降速度快，这主要是由于射流冲击到垂直平板上时，在轴线附近形成一个滞止区，在该区射流速度被滞止，压强和温度升高；当射流离开滞止区形成贴壁射流后，气体流速加快，压强和温度有所降低。对应到传热问题中，由于滞止区的气流温度高，平板表面气流边界层的量纲一过余温度梯度大，表征射流与平板之间传热强弱的努塞尔数大，传热效率较高；离开滞止区，量纲一过余温度梯度变小，传热效率相应下降。而且从图4～5可以看出，随着喷口与平板间距H的增加，贴壁射流温度曲线依次降低，努塞尔数相应减小，也表明冲击到平板上的气流温度越高，传热效果越明显。

图4 超音速射流垂直冲击平板的努塞尔数分布Fig.4 Nusselt number distribution in vertical plate impinged by supersonic jet flow

图5 在距垂直平板0.5 mm截面的温度分布Fig.5 Temperature distribution of the section which is 0.5 mm away from the impinged plate

3 超音速燃气射流冲击浸没平板传热问题

3.1 努塞尔数分布规律

采用总压为4 MPa、总温为3 117 K的超音速射流从直径D=0.024 3 m的喷管口喷出，冲刷在顶端面距离喷口分别为14D、16D、18D的浸没平板上，平板长L=2D，宽R=2D，厚度Z=0.005 m。

为了便于分析，设置正对着喷口的面为顶端面，顶端面与喷口之间距离为H，沿着喷管直径方向为径向，距离为R，平板上沿喷管轴线方向为其长度方向，距离为L，平板上与射流轴线平行的壁面在文中称为平行冲刷面。

图6依次给出浸没平板顶端面与喷口距离为14D、16D、18D时的努塞尔数等值线分布图。从图中可以看出，浸没条件下的3种状态中沿着径向R，在同距离点处的努塞尔数依次减小，说明随着冲击距离H增大，努塞尔数减小，传热效率减弱；冲击距离H的增加影响传热效率的因素，如密度、比热容、导热系数、黏度等，而这些因素都是温度的函数[15]。图7依次给出了贴近浸没平板顶端面的燃气射流层温度等值线分布图。从图中可以得到，随着冲击距离从14D到18D，流过浸没平板顶端面的气体温度依次减小，努塞尔数减小，传热减弱，进一步说明了气体射流的温度对传热效率影响较大，而且温度梯度越大，热传导方式的作用就越强烈，传热效率就越高。再与图4中垂直冲击平板的努塞尔数相比，在数值上减小了，原因是浸没平板正对着喷口的面积小，射流冲击形成的滞止区小，传热效率也就明显降低。

图6 浸没平板顶端面的努塞尔数分布Fig.6 Nusselt number distribution on top surface of submerged plate

图7 浸没平板顶端面的射流温度分布Fig.7 Temperature distribution on top surface of submerged plate

图8 H/D=14时浸没平板的努塞尔数分布Fig.8 Nusselt number distributionon submerged plate under H/D=14

图9 H/D=14时浸没平板的温度分布Fig.9 Temperature distribution on submerged plateunder H/D=14

图8～9分别给出H/D为14时平行于射流轴线的冲刷面上努塞尔数分布和温度分布等值线图。可以看出在平行冲刷面上，沿着径向R，努塞尔数明显减小，传热主要发生在一倍喷口直径的范围内。在努塞尔数分布图中，沿着长度L方向，在L/D=0.5～1范围内，出现一个低的传热区，而在L/D=1～1.5范围内出现一个传热高的区域。这是由于射流冲击到浸没平板的顶端面后，由于浸没平板的阻滞作用，射流绕过平板并在平行冲刷面附近形成局部绕流，在L/D=0.5～1范围内，贴近平行冲刷面的射流流速梯度较大，于是在该处形成一个低压、低温的区域，量纲一过余温度梯度小；绕过该区之后，燃气在L/D=1处开始沿喷管轴线又形成贴壁流动，由于黏性作用，射流速度减小，压强和温度相应升高，量纲一过余温度梯度变大。在浸没平板平行冲刷面表层的燃气射流温度分布图和平行冲刷面上努塞尔数分布图呈现相同的情形，说明了温度梯度的大小对努塞尔数分布影响很大，贴壁射流温度是影响传热效率的关键因素。

3.2 浸没平板传热统计特征

图10 浸没平板冲刷面的平均努塞尔数分布Fig.10 Average Nusselt number distribution on eroded surface of submerged plate

对比图6和图8，顶端面的努塞尔数比平行冲刷面上的高，原因在于受浸没平板顶端面对气流阻滞作用的影响，浸没平板上垂直于喷管轴线面的传热效率高于平行面。为便于理论分析和工程应用，对浸没平板上努塞尔数的特征进行了统计分析，取浸没平板平行冲刷面上沿着垂直于射流轴线的平均努塞尔数，如图10所示。从图中看出，平均努塞尔数在L/D=0时最大，然后随着冲击距离的增加而减小，而在L/D=1～1.5范围内出现了一个峰值，这与上面分析图6努塞尔数等值线的情况相吻合。

4 结 论

利用数值计算方法，对超音速燃气射流冲击平板传热特性进行研究分析，得出如下结论：

(1)采用基于RNGk-ε湍流模型的数值计算方法，求解气体射流与固体壁面之间的对流传热问题，具有较好的效果；

(2)超音速燃气射流垂直冲击平板的努塞尔数比亚音速射流冲击平板所得努塞尔数大一个数量级，滞止点附近燃气温度较高，量纲一过余温度梯度大，所以努塞尔数大，在滞止点附近射流的传热效率高；

(3)冲击距离对传热效率影响较大。对于H/D=14～18时，垂直冲击平板的努塞尔数在滞止点从1 200～1 600变化；而对于浸没平板的平行冲刷面,平均努塞尔数在70～170之间。

[1] Baydar E, Ozmen Y. An experimental and numerical investigation on a confined impinging air jet at high Reynolds numbers[J]. Applied Thermal Engineering, 2005,25(2):409-421.

[2] Chen Y C, Ma C F, Qin M, et al. Forced convective heat transfer with impinging slot jets of mesoscale[J]. International Journal of Heat and Mass Transfer, 2006,49(1):406-410.

[3] Liu Q, Sleiti A K, Kapat J S. Application of pressure and temperature sensitive paints for study of heat transfer to a circular impinging air jet[J]. International Journal of Thermal Sciences, 2008,47(6):749-757.

[4] Ramanujachari V , Vijaykant S, Roy R D, et al. Heat transfer due to supersonic flow impingement on a vertical plate[J]. International Journal of Heat and Mass Transfer, 2005,48(17):3707-3712.

[5] Behniat M, Parneix S. Prediction of heat transfer in an axisymmetric turbulent jet impinging on a flat plate[J]. International Journal of Heat Mass Transfer, 1998,41(12):1845-1855.

[6] 陈庆光,徐忠,张永建.半封闭圆管湍流射流冲击平板的数值研究[J].西安交通大学学报,2001,35(11):1206-1208. Chen Qing-guang, Xu Zhong, Zhang Yong-jian. Numerical investigation of semi-confined round turbulent jet impinging on a flat plate[J]. Journal of Xi’an Jiaotong University, 2001,35(11):1206-1208.

[7] 陈庆光,徐忠,张永建.用改进的RNG模式数值模拟湍流冲击射流流动[J].西安交通大学学报,2002,36(9):916-920. Chen Qing-guang, Xu Zhong, Zhang Yong-jian. Numerical investigation of semi-confined round turbulent jet impinging on a flat plate[J]. Journal of Xi’an Jiaotong University, 2002,36(9):916-920.

[8] Merci B, Dick E. Heat transfer predictions with a cubick-εmodel for axisymmetric turbulent jets impinging onto a flat plate[J]. International Journal of Heat and Mass Transfer, 2003,46(3):469-480.

[9] 许坤梅,张平.半封闭圆管冲击射流湍流换热数值模拟[J].北京理工大学学报,2003,23(5):540-544. Xu Kun-mei, Zhang Ping. Numerical prediction of turbulent heat transfer in a semi-confined impinging jet[J].Transactions of Beijing Institute of Technology, 2003,23(5):540-544.

[10] 傅德彬.数值仿真及其在航天发射技术中的应用[M].北京：国防工业出版社，2011:28-34.

[11] 李鹏飞，徐敏义，王飞飞.精通CFD工程仿真与案例实战[M].北京：人民邮电出版社，2011:329-338.

[12] 赵承庆.气体射流动力学[M].北京：北京理工大学出版社，1998:103-106.

[13] Heck U, Fritsching K, Bauckhage K. Fluid flow and heat transfer in gas jet quenching of a cylinder[J]. International Journal of Numerical Methods for Heat & Fluid Flow, 2001,11(1):36-49.

[14] Wilcox D C. Turbulence modeling for CFD[M]. La Canada, California:DCW Industries,2002.

[15] 赵镇南.传热学[M].北京：高等教育出版社，2008:173-192.

(责任编辑 丁 峰)

Numerical simulation of heat transfer for exhausted gases jet impinging

Liu Xiao-jun, Fu De-bin, Niu Qing-lin, Li Xia

(SchoolofAerospaceEngineering,BeijingInstituteofTechnology,Beijing100081,China)

To the case of heat transfer, the flowing process of jet flow impacting on a plate vertically is simulated by employing RNGk-εturbulence model, which is compared with the experimental data, to verify the feasibility of the model. Based on the simulating results, the models of the impact of supersonic jet flow on plate vertically and on submerged plate are built respectively by considering the parameters of rocket nozzle entrance as the inlet conditions. In addition, the distributions of the Nusselt number and temperature are calculated under different impacting conditions. Moreover the characteristics and factors of supersonic jet flow heat transfer are analyzed. The results show the range of Nusselt number under different impacting distances are between 14Dand 18D, and reflect that the impacting distance and jet flow temperature are the key factors which influence the heat transfer rate. Furthermore, when the impacting distance increases, the heat transfer rate decreases. In contrast, the higher the temperature of jet flow on the plate is, the greater the efficiency of heat transfer is.

fluid mechanics; heat transfer for jet flow; RNGk-εturbulence model; submerged plate; Nusselt number

10.11883/1001-1455(2015)02-0229-07

2013-08-29；

2013-12-13

国家自然科学基金项目(51306019)

刘小军(1985— )，男，硕士研究生; 通讯作者： 傅德彬,fdb007@bit.edu.cn。

O358;V411.3 国标学科代码： 1302531

A