再生冷却推力室热机疲劳寿命预测研究

黄峻峰, 贺尔铭, 易金翔, 杜大华, 王红建

(1.西北工业大学 航空学院, 陕西 西安 710072; 2.液体火箭发动机技术重点实验室, 陕西 西安 710100)

可重复使用火箭发动机是发展绿色航天工程的重点研究工作之一[1-2]。美国私立太空企业SpaceX公司近年成功研制的“梅林”火箭发动机使得“猎鹰9号”和“重型猎鹰”运载火箭具有了一定的重复发射能力。英国开展能重复使用超过100次的“佩刀”发动机的研制工作，并计划将其安装在“云霄塔”太空飞机[3]。欧洲宇航院(ESA)的“未来发射器”计划，基于先进3D打印技术研制的“普罗米修斯”低成本发动机及可重复使用飞行器样机“西弥斯女神号”都在进行紧锣密鼓的试验[2]。俄罗斯针对10～500 kg小载荷商业发射市场，正研制可重复使用的“飞翼-SV”轻型运载火箭及相关发动机[4]。我国已成功研制出50吨级氢氧发动机YF-77和120吨级液氧煤油发动机YF-100，并作为新一代运载火箭的动力系统推动近些年我国航天事业的巨大发展[5-6]。与此同时，我国已展开了可重复使用液氧煤油发动机、液氧液氢发动机、液氧甲烷发动机关键技术攻关，对可重复使用运载火箭的研制制定了“三步走”的发展战略[7]。

推力室作为液体火箭发动机的核心结构，承担着推进剂混合燃烧、释放热量、产生推力及完成能量转化的工作，自然成为了可重复使用关键技术研发和验证的重点结构。在推力室工作过程中，其内部的燃气压力、温度及流体物性变化十分剧烈，燃气压力最大可高达200个大气压，燃气温度最大可达到3 600 K，最大速度可超过6马赫[8]。这些复杂的载荷循环会使推力室薄壁夹层结构产生低周疲劳和高温蠕变效应，而且推力室在循环加载下受到的载荷往往是非对称应力控制的循环加载，这会导致结构无法恢复的塑性应变逐渐累积产生棘轮效应，大大缩短推力室的使用寿命。

低周疲劳、高温蠕变效应和棘轮效应这3种行为的共同作用会造成推力室薄壁夹层结构失效，大大降低了推力室的使用寿命。国内外专家学者先后对推力室结构的失效机理及寿命预测方法进行了大量有效工作。

羽中豪等[9]使用准二维传热计算方法、二维热-固耦合方法和局部应变法研究了室压、推力燃料混合比对推力室棘轮效应的影响。孙冰等[10]使用二维热-固耦合分析方法对推力室内壁进行非线性平面应变分析，并预测了内壁结构的使用寿命。孙冰和宋佳文[11-13]使用整场耦合方法对液氧甲烷发动机推力室内部的燃烧和冷却剂流动换热特性进行三维仿真，获得了推力室在不同试车阶段的应变分布，并比较了瞬态加载和稳态加载下的最大残余应变。Di等[14]研究了喷油器二维喷流导致推力室壁面周向温度分布不均的现象，并发现棘轮效应是导致结构失效的主要原因。Pizzarelli[15]用代数模型对推力室使用寿命进行预测分析，预测准确率超过80%。Riccius等[16]发现由燃烧造成的推力室铜合金内壁粗糙度的增加能降低推力室28%的疲劳寿命。

本文针对某液体火箭发动机再生冷却推力室的薄弱环节结构——流道喉部薄壁夹层结构，建立其局部二维平面应变有限元模型，通过传热仿真和非线性结构分析获得循环载荷下的温度分布和应力应变响应变化规律，研究了该结构产生的棘轮效应，并揭示了棘轮、蠕变及疲劳损伤对推力室结构使用寿命的影响。

1 物理模型与计算方法

1.1 物理模型

再生冷却推力室采用了大量如图1所示的带铣槽的钎焊薄壁夹层结构，图中从内至外依次为内壁、再生冷却通道及外壁，相邻冷却通道之间由肋片隔开。在设计阶段，为了使火箭发动机具有高推重比，在保证强度/刚度足够的条件下尽可能使结构轻量化，导致发动机的薄壁夹层结构设计几乎达到了材料的承载极限。

图1 推力室薄壁夹层结构及“狗窝”失效示意图

本文研究对象为某型液氧煤油火箭发动机推力室喉部结构，该结构的内壁和肋片采用NARloy-Z铜合金材料，外壁则采用1Cr21Ni5Ti不锈钢材料，如图2所示。

图2 推力室结构及其喉部局部细节

推力室喉部截面是循环对称结构，内壁和外壁之间均匀分布着270条再生冷却通道。根据其轴对称特性，选取如图3所示的二分之一冷却通道截面作为计算模型。

图3 二分之一冷却通道物理模型及有限元模型

使用八节点四边形单元对模型进行网格划分，共生成15 592个节点和4 901个单元。为了便于后续工作描述模型温度变化和预测寿命，在模型上依次选取7个参考点，其中A,B点分别为外壁面上、下端，C,D点为冷却通道肋片的上下顶点，E,F,G点位于靠近高温燃气通道的内壁面上。

1.2 控制方程

本文采用单向耦合方法，先通过传热分析计算整个工作循环中结构的温度场分布，再在其基础上进行瞬态结构分析得到结构的应力应变演化结果。

传热分析的控制方程为无内部热源的导热微分方程

(1)

式中:T为结构温度；λ为推力室材料导热率；ρ为材料密度；c为推力室材料比热容；为哈密顿算子。

在进行结构分析时,结构的总应变εtotal由热应变εth、弹性应变εel、塑性应变εpl和蠕变应变εcr组合而成

εtotal=εth+εel+εpl+εcr

(2)

在计算热应变时使用的公式为

εth=aΔT

(3)

式中：a为材料的热膨胀系数；ΔT为相对参考温度的温度变化量。

在计算弹性应变和塑性应变时,使用Von Mises屈服准则来判断材料是否达到屈服,其表达式为

(4)

式中：f为屈服函数；s为偏应力张量；α为背应力张量,即屈服面中心；σy为屈服应力。达到屈服面后塑性应变产生,塑性应变沿着屈服面梯度方向增加

(5)

式中：dλ为塑性乘子,即等效塑性应变增量；σ为应力张量。

本文使用Chaboche随动硬化模型来描述背应力张量的演化过程,在该本构关系中,总背应力由多级背应力分量叠加得到,如(6)式所示

式中:M为背应力分量的级数;Ci和γi为随动硬化参数;dp为累积塑性应变,其表达式为

(8)

由于金属材料力学性能受温度影响较大,本文将Esposito等[17]得到的NARloy-Z材料在高低温环境下的单轴拉伸试验数据,使用背应力级数为4的Chaboche模型进行拟合,得到不同温度下的随动硬化参数，如表1所示。试验曲线与拟合曲线如图4所示,可以看到拟合曲线与试验曲线[17]在弹性段和塑性段的吻合度都很高,说明随动硬化模型参数合理准确,能够准确描述该材料的弹塑性变形。

表1 NARloy-Z材料随动硬化参数

图4 NARloy-Z材料不同温度下的单轴拉伸曲线

在计算高温蠕变应变时,使用Norton蠕变模型来进行描述,其表达式为

(9)

表2 Norton蠕变模型参数

1.3 寿命预估模型

推力室喉部结构在每个循环下产生的总损伤Dtotal由低周疲劳损伤Df、棘轮损伤Dr和蠕变损伤Dc组成,可按照Miner线性累积损伤理论来计算每个循环下的总损伤[10]

Dtotal=4×(Df+Dr+Dc)

(10)

式中,4为安全因子。

在计算低周疲劳损伤时,以该次循环应力水平下的疲劳循环寿命Nf的倒数来定义低周疲劳损伤

(11)

本文使用基于平均应力的Morrow修正模型来计算每个循环下的低周疲劳损伤,其表达式为

(12)

每个循环下棘轮损伤可以用残余应变增量与材料断裂应变的比值来表示

(13)

式中：εend为当前循环结束时的残余应变；εbegin为当前循环开始时的残余应变；εu为材料单轴拉伸时的断裂应变。

蠕变损伤定义为当前应力水平和温度下的累计保载时间t与该应力下的断裂时间tc之比

(14)

1.4 边界条件设置

该型液体火箭发动机一个完整的循环加载分为启动-热试车-后冷-关机4个阶段,各阶段时长分别为5,300,20,600 s,环境温度为295.15 K。在不同的试验阶段下,对模型边界条件的设置是不同的。

1) 热分析的边界条件设置

将图3中对称边界Γsym设置为绝热边界,即在该边界法线方向上的温度梯度为0

(15)

式中:n为边界的法线方向;hsym为对称边界的对流传热系数。

而内壁与外部环境、冷却通道和高温环境直接接触的边界为耦合边界,在该边界上存在如(16)式所示的对流关系

hcon(Tbulk-Twall)=n·λTwall, ∀x∈Γcon

(16)

式中:Γcon为耦合边界;hcon为耦合边界的对流传热系数;Tbulk为流体温度;Twall为固体壁面温度。在热分析中各耦合边界条件设置如表3所示,表中下标out、cool和hot分别表示外壁面、冷却通道壁面和内壁面,CFD代表由流-热耦合计算得到的壁面参数设置。

表3 热分析边界条件设置

2) 结构分析的边界条件设置

在进行结构分析时,将传热分析得到的每一步结构温度分布作为温度体载荷。由于推力室薄壁夹层结构对称特性,将模型两边设置为对称边界。将各阶段下耦合边界上的载荷分别施加在相应边界上,各阶段下边界载荷如表4所示。

表4 结构分析边界条件设置 MPa

2 计算结果与分析

2.1 热分析结果

图5为不同阶段下的结构温度分布,图6为结构上7个参考点温度随时间变化曲线。可以看到,

图5 各阶段结构温度分布

图6 参考点温度变化曲线

由于外壁材料的导热率小于内壁材料,内壁上的参考点温度会很快到达稳定,而外壁面上的参考点则会出现一定的延迟才能到达稳定。在整个工作过程中,下表面出现最高温度且该部位温度变化最大,温度最大点出现在热试车阶段,此时内壁下表面G,F点处的最高温度分别为862 K,851 K。

2.2 结构分析结果

图7和图8分别为各阶段下结构的应力和应变分布。

图7 各阶段结构应力分布 图8 各阶段结构应变分布

在开机阶段,结构总体温度仅提升到309 K左右,虽然结构整体开始膨胀,但由于冷却通道内煤油的压力载荷较低,整体结构的压应变和压应力都处于较低水平。

在热试车阶段,外壁和肋片上半部的温度普遍低于450 K,该区域膨胀时内部的压应变较小。同时由于外壁的热膨胀系数小于肋片的热膨胀系数,肋片阻碍了外壁的膨胀,所以外壁C点处出现了最大压应力434.23 MPa。高温燃气会在内壁下表面产生较大的压力载荷进而阻止附近结构的膨胀,越靠近燃气的区域受到的阻碍作用越大,因此内壁下表面结构内部的压应力要小于内壁上表面。此时内壁下表面的G点处出现了最大压应变2.307×10-2,这是因为该阶段下冷却通道内煤油作用在内壁上表面的压力大于高温燃气对内壁下表面的压力,使得内壁向燃气侧凸出并弯曲,由此造成G点附近区域受压严重。

在后冷阶段,内壁区域温度在短时间内降低了很多,其中下表面温度最高降幅达550 K,导致该区域出现了收缩趋势,因此该区域的应力和应变普遍为正。冷却通道内煤油对周围壁面仍有压力作用,而燃气作用域内壁下表面的压力在该阶段因燃烧终止而消失,导致内壁向燃烧室凸出的趋势更加明显,F点出现了最大拉应变4.467×10-3,拉应力值为205.13 MPa。

卢一平没有马上回驳她。这时候,回驳就是争论,就是口角,就是龃龉,就是搁浅。卢一平沉默了,沉默地盯视着郝桂芹。盯一会儿,有意无意地抬起手腕,看了眼手表。

在关机阶段,由于结构在之前的试验阶段已经出现了塑性应变,在卸掉热载荷和压力载荷后这些不可恢复的塑性应变作为残余应变仍存在结构中。关机阶段结束时的结构温度相较于后冷阶段结束时的温度只降低了约14 K,所以两时刻下的应力和应变分布云图相差不大。F点处的最大残余应变为4.437×10-3。

2.3 循环载荷下结构棘轮效应研究

在获得单次循环加载下结构的温度变化和应力应变演化规律后,对10个连续循环加载下结构的应力应变演化过程进行了仿真。图9为不同循环加载次数后的结构残余应变分布,可以看到残余应变集中出现在内壁区域,并且每个循环下残余应变的分布云图基本相同,只是残余应变的绝对值随着循环数的增加而逐渐增大。

图9 多次循环加载后的残余应变分布

图10为内壁上D,E,F,G点在10次循环加载中的应力-应变曲线,可见每个点都是应力控制的非对称循环加载,并且都产生了棘轮效应。而且在后续几个循环中参考点的应力-应变曲线形状趋势基本相同,呈现了一种渐进稳定性,这就是材料在非对称应力循环中的调整行为。

图10 各点应力-应变曲线

为了定量评估各点的棘轮效应大小,采用棘轮应变作为评估指标。使用每个循环中初始残余应变和循环结束时的残余应变之差来确定棘轮(残余)应变

εr=εend-εbegin

(17)

图11为各点的棘轮应变发展曲线,可以看到在第二次循环后各点的棘轮应变变化与循环次数基本上呈线性变化,即体现了非对称应力循环中的调整行为。随着循环数的增加,位于内壁上表面的D,E点棘轮应变大小变化不大,而位于下表面的F,G点的棘轮应变逐渐增大。其中F点在每个循环下残余拉应变率约为2.3×10-3/次,10次循环加载后残余拉应变值为2.529×10-2;G点的残余压应变率约为2.5×10-3/次,10次循环加载后的残余压应变分别为3.193×10-2。

图11 各点棘轮应变发展曲线

根据结构残余应变分布,综合各点的应力-应变曲线和棘轮应变发展曲线,发现内壁下表面的F点和G点会成为结构的潜在失效点。这是因为该两点的残余应变最大,而且该两点在每个循环中的应力和应变变化范围也较大导致其低周疲劳损伤也很大。

图12为10次循环后结构的蠕变应变分布及位移云图,可以看出高温蠕变效应集中出现在温度较高的内壁下半部,位于内壁下表面的F点和G点的蠕变效应明显,进一步加速了这两点所在区域的结构损伤。

图12 10次循环后结构蠕变应变及径向位移分布

图13为10个循环下内壁上下表面的E点和F点的径向位移发展曲线。随着循环数的增加,内壁会逐渐向燃烧室凸出并减薄,这就是推力室结构典型的“狗窝”失效模式。

图13 E点和F点的径向位移发展曲线

2.4 推力室结构寿命预测

F点和G点在10次循环下的损伤计算结果如图14所示。可以看到F点的损伤形式为棘轮损伤、蠕变损伤及低周疲劳损伤,其中棘轮损伤为最主要的损伤形式,而G点的损伤由蠕变损伤和低周疲劳损伤构成。并且由于调整行为,两点每个循环下的损伤基本保持不变,其中F点的单循环损伤约为2.95×10-2,总损伤为31.12×10-2,G点的单循环损伤约为2.32×10-2,总损伤为24.124×10-2。根据(10)式计算寿命,F点、G点的失效循环数分别为33次及42次。

图14 F,G两点损伤计算结果

综上，在经过33次循环加载后，推力室内壁下表面中心点会首先发生失效破坏。在导致该点的总损伤中，棘轮损伤、蠕变损伤和低周损伤分别占比52%,32%和16%，说明了棘轮效应是导致推力室结构失效的主要原因。

3 结 论

本文建立了某液氧煤油火箭发动机推力室喉部结构的周期对称模型，通过三维传热分析获得结构的循环载荷特征，基于各工况下结构应变及损伤分析初步揭示棘轮、蠕变及疲劳对推力室结构使用寿命的影响，得到以下研究结论：

1) 本文建立发动机推力室喉部结构的循环对称模型，进行推力室流-热-固载荷传递及结构损伤分析，文中寿命预测模型及方法计算效率高、具有重要工程参考价值；

2) 定量分析了推力室内壁结构“狗窝”失效的影响机理，在冷却通道总损伤中棘轮损伤、蠕变损伤和低周疲劳损伤分别占比52%,32%和16%，因此棘轮效应是导致其失效的主要原因;

3) 结构的蠕变效应主要集中在内壁的下半部分，这是因为在热试车阶段该区域的温度很高，进而引起应力松弛现象。

本文的研究方法及结果可为液体火箭发动机再生冷却推力室结构优化设计及快速寿命预测提供重要的工程参考。后续工作可进一步考虑推力室工作时序效应、各损伤模式之间的耦合作用及相关的试验验证，将会使推力室结构寿命预测结果更加准确。