在计算模型中,最关键的是燃气侧的传热,主要问题在于推力室喉部区域几何和燃气参数变化过于激烈,另外推进剂的喷注和燃烧过程还影响燃气的流动状态,在近壁面的边界层具有很高的温度梯度,一个纯粹的理论模型根本无法描述真实的情况,因此工程上常采用经过实验数据验证后的半经验公式^[6]。

传统Bartz法采用的是下式计算对流换热系数:

h_g,Bartz=(0.026)/(D_t^0.2)·(μ_g^0.2·C_pg)/(P_rg^0.6)·((P_c)/(C^*))^0.8·((A_t)/(A⁾)^0.9·σ (1)

式中σ为定性温度变换系数,由下式计算

σ={[1/2(T_wg)/(T_c)·(1+(k-1)/2M^_a2)+1/2]^0.68[1+(k-1)/2M^2]0.12}-1(2)

燃气的动力粘度μ_g和普朗特数P_rg根据燃气的其他热力参数近似计算,在喉部附近换热系数可以采用喉部曲率半径修正。

为了考虑圆柱段和喷管边界层厚度变化带来的影响,增加了一个修正系数K_x:

K_x=(x/(x_th))^a(3)

式中x_th指喉部处横向坐标,圆柱段初始位置为原点。

因此Modified Bartz公式为:

h_g,Bartz,mod=h^g_,BartzK_x(4)

另外一种Pavli方程同样考虑了相关修正系数:

h_g,Pavli=0.023Re^-0.2Pr^-0.6ρ·u·C_pg·K_x(5)

为了更加准确的计算燃烧室的热流密度,需要考虑喷注面板附近喷雾现象,在实际过程中,因为这部分区域推进剂的不完全混合和反应,靠近喷注器面板附近的燃烧区域热流密度在减小。对于特定燃烧室而言,燃烧区域长度x_max是确定的,可以根据喷注器和燃烧室几何尺寸计算得到,在燃烧区域x₀xx_max,可以通过下式计算燃气对流传热系数^[3]:

(h^*_g(x))/(h_g(x_max))=1/4arctan[7(x/(x_max)-0.63)]+0.7(6)

最后再增加了一个流动加速相关的修正系数K_acc,因为实际测量的喷管喉部附近热流密度比采用传统巴兹法计算得到的值要小,这可能与喷管轮廓以及加速导致的边界层厚度增加相关,为了更加简单方便的计算这种因素导致的热流密度变化,定义了K_acc如下式:

K_acc=(1-|(dr)/(dx)|)^1/2(7)

式中:dr为燃烧室半径的变化; dx为横向坐标变化量,在圆柱段K_acc=1。

冷却剂对流换热公式采用NASA甲烷电热管试验数据得到的经验公式:

Nu=0.023·Re^0.8Pr^0.4((T_l)/(T_bl))^0.45(8)

再生冷却推力室有铣槽式冷却通道,通道间高热导率的锆铜肋条加强了内壁与冷却剂的换热,因此计算时必须考虑冷却肋的影响^[7]。在程序中,冷却剂侧热流密度需乘上肋条的散热系数。

根据热传导定律^[8],通过内壁的热流为:

q_w=(λ_w)/(δ_w)(T_wg-T_wl)(9)

在实际发动机沿轴向和径向内壁温度是变化的,因此对内壁材料锆铜的导热系数进行随温度的插值处理。

采用上述三种方法对该型号发动机开展额定工况下的再生冷却传热计算,分别得到了热流密度、气壁温、冷却剂温度及压力的轴向分布,计算结果如图2-图5所示。

图2 热流密度轴向分布
Fig.2 Axial distribution of heat flux density

图3 气壁温
Fig.3 Temperature of gas wall

图4 冷却剂温度曲线
Fig.4 Temperature curves of coolant

图5 冷却剂压力曲线
Fig.5 Pressure curves of coolant

从图2可以看出,考虑燃烧区域的影响后,采用Modified Bartz法和Pavli法得到的热流密度曲线,均能反映出推力室头部附近低热流的状况,在燃烧区域x₀xx_max热流密度不断增长。而Bartz法在整个圆柱段热流密度基本保持不变,与实际情况不符合。三种方法算出的最大热流密度均发生在喉部,分别为:25.468,25.024和19.009 MW/m²,其中Pavli法算出的最大热流密度与前两种方法计算差异较大,但是整体趋势与修正后的Bartz法较为一致。图3结果表明,推力室内壁气壁温度在圆柱段,头部以及喉部较高,尾部较低,这是由于喉部热流密度大,换热效果强,而圆柱段冷却剂温度较高,对内壁的冷却效果减弱。推力室身部内壁最高气壁温不超过800 K,可以满足身部热防护要求。

图4和图5分别为冷却剂的温度和压力曲线,在工程实践中,对于冷却夹套而言,最重要的参数就是冷却剂的温升及流阻,只有保证足够大的温升以及足够小的流阻,从夹套出来后的气态甲烷才有足够的能量带动涡轮的转动,从而满足发动机工作要求。从图4可以知道,冷却剂温度沿流动方向逐步升高,在出口处达到最大值,三者出口温度分别为467.5 K,441.3 K和393.7 K。冷却剂压力由于流阻的影响,不断减小,其中在喉部区域减小迅速,因为喉部区域冷却剂流通面积小,流速大,流阻也较大,在出口处压力达到最小值,出口压力分别为5.381 MPa,5.349 MPa和5.185 MPa。该型号推力室30 s挤压点火试验结果显示,冷却夹套入口与出口压力稳态平均值分别为7.30和5.25 MPa,入口与出口温度稳态平均值分别为143 K和428 K。

与试验值相比后,可以计算三种方法得到温升与流阻的误差,从而反映出三种方法在计算再生冷却温升与流阻方面的准确性。结果如表1所示。

表1 温升与流阻误差计算
Tab.1 Error calculation of temperature rise and pressure loss

从表中可以看出,Bartz法和Pavli法得到的温升误差均在10%以上,而Modified Bartz法误差只有4.67%,与实验结果非常符合,说明考虑燃烧区域、流动加速性的影响后极大地减小了计算的误差。在流阻方面,误差均较为接近,其中Pavli方法算得的误差最小。对比两种重要指标的计算误差后,可以发现Modified Bartz法计算值与试验值最为接近,可以较为准确的指导工程实践。

发动机在变工况工作状态下,室压和混合比会偏离额定值,燃烧室温度也会随之改变,Modified Bartz法能否在不同室压及混合比下良好的工作,直接影响了发动机动态仿真的结果。

图 6所示为Modified Bartz法在混合比均为2.92,室压分别是3 MPa,3.46 MPa和4 MPa下热流密度分布曲线。在三种室压下,热流密度分布趋势与额定压力下基本一致,可以满足该型号发动机在变工况情况下的使用。室压越高,燃气温度越大,会导致喉部的热流密度急剧增大,对换热的要求也会更高。图7所示为在室压均为3.46 MPa,混合比分别是2.6,2.92和3.2下热流密度分布曲线。在上述三种混合比下,圆柱段以后的热流密度差异很小,但是依然可以发现混合比越大,热流密度也越大,因为混合比越大,特征速度越小,在燃烧效率一定下,导致燃烧室流量增大,燃烧室温度也变大。因此采用Modified Bartz法能够较为准确地计算不同工况下推力室再生冷却结果。

图6 不同室压下热流密度
Fig.6 Heat flux density at different chamber pressures

图7 不同混合比下热流密度
Fig.7 Heat flux density at different mixing ratios