垂直着陆过程推进剂流动行为特性及影响分析

尕永婧，王浩苏，张青松，徐珊姝，吴义田

（北京宇航系统工程研究所，北京，100076）

引 言

运载火箭可回收并多次使用是降低发射成本的重要途径，然而，返回着陆过程的实现充满了困难与挑战。垂直返回是指火箭子级在完成任务后，通过自身携带的控制系统和动力装置，按照设定的轨迹自主飞回着陆场，并以垂直的箭体姿态缓慢稳定地降落到指定着陆场。运载火箭以垂直起降方式实现重复使用的需求对动力系统、姿控系统、载荷与热环境等专业提出了新的挑战。垂直起降过程不仅要求动力系统能够提供上升段推力，还需实现垂直着陆的平稳减速和控制。为实现精准垂直着陆，发动机需具有深度推力调节及多次启动的能力，同时为满足返回过程中发动机正常启动与工作要求，需对贮箱内的推进剂进行有效管理。如安装防晃档板和消漩防塌装置抑制推进剂晃幅和减小不可用量、选用推进剂管理装置（PMD）提高推进剂的蓄流能力、采用脉冲式沉底和连续沉底等方案保证推进剂可靠沉底，从而防止进入发动机的推进剂管路夹气，以及由推进剂晃动加速换热带来的气枕压力大幅下降无法满足发动机最低泵入口压力的情况发生。

本文在调研国内外运载火箭贮箱内推进剂流动行为特性的基础上，采用Flow 3D数值模拟垂直着陆过程贮箱内推进剂流动过程，分析垂直着陆过程中箱内推进剂的流动行为特性及相关影响因素。

1 国内外研究现状及发展趋势

当前实现运载火箭子级回收的主要技术途径有伞降回收、带翼飞回以及垂直着陆回收3种方式。伞降回收方式技术具有成熟度高、可靠性高和运载能力损失少的优点，已成功在我国的“神舟”载人飞船返回舱、美国航天飞机固体助推器等多项工程任务中应用[1～2]。垂直着陆回收方式相较伞降回收方式落点精度更高，相较有翼水平着陆方式，对火箭外形和总体布局影响更小，因此近年来，垂直着陆回收方式获得了广泛的关注和大量的研究。同时，在线轨迹规划技术的发展应用也可显著提升火箭在返回着陆过程的自适应调整和在线实时优化能力[3]。最为典型的代表就是蓝色起源公司（Blue Origin）的“新谢帕德号”（New Shepard）火箭以及SpaceX公司利用“猎鹰9号”（Falcon 9）系列运载火箭执行的多次一级可靠垂直着陆回收试验[4]。

“猎鹰9号”可重复使用火箭子级返回时靠主发动机点火反推实现减速，垂直着陆时由安装在子级下部的着陆支架支撑。SpaceX公司在进行子级回收试验前，曾进行了8次低空（ ≤ 3.5 km）低速的蚱蜢（Grasshopper）试验（见图1）和4次高空（3.5～9.1 km）中速的“猎鹰9R Dev”验证机试验。2013年9月29日，“猎鹰9V1.1”火箭在范登堡空军基地首次执行发射任务，并开展一子级海面软着陆试验。但在回收过程中火箭旋转产生的离心力使燃料无法进入燃料管，导致发动机因燃料供应不足而提前关机，无法有效控制火箭飞行速度和姿态，致使一子级砸在水面上，实验失败（见图2）[4]。SpaceX公司后续通过携带更多姿控推进剂，保证姿控系统可靠工作，成功消除了箭体返回过程的滚转，实现了有效的推进剂管理。

图1 蚱蜢验证机[4]Fig. 1 Grasshopper verification machine[4]

图2 一级着陆装置展开[4]Fig. 2 Deployment of stage I landing gear[4]

2015年12月，SpaceX公司实现了火箭一子级陆地发射场成功回收；2016年4—5月成功完成了一子级3次海上平台垂直回收（见图3）。

图3 海上着陆Fig. 3 Landing at sea

Behruzi等[5]采用Flow 3D对“阿丽亚娜5号”（Ariane 5） ESC-A 上面级分离过程中液氢箱推进剂的流动过程进行数值模拟，结果如图4和图5所示。分析表明，在有分离、晃动等外界干扰因素的情况下，推进剂晃动并向贮箱顶部运动，加剧了气枕与推进剂的接触面积和换热，致使贮箱内的气枕压力快速下降。由此可见，外界干扰引起的晃动可能是贮箱内气枕压力下降的极为重要的影响因素。

图4 Ariane 5火箭第一次飞行上面级氢箱压力变化及仿真计算[5]Fig. 4 Hydrogen tank Pressure variation of the Ariane 5 upper stage[5]

图5 分离过程中贮箱内推进剂形态[5]Fig. 5 Flow characteristics of propellant in tank during seperation[5]

Konopka[6]和Rose等[7]将GNC算法纳入推进剂晃动抑制过程中。诸桂敏[11]采用试验和Flow 3D数值模拟的方法，研究了低温上面级推进剂晃动的主动抑制技术，结果表明：在推进剂晃动过程中引入主动抑制技术，可有效抑制推进剂晃动，并能使箱内推进剂晃幅快速减小，详见图6和图7。

图6 基于控制模型、实物模型和数值仿真的实验系统[6]Fig. 6 Experimental system based on control model, physical model and numerical simulation[6]

图7 开环和闭环控制实验对推进剂晃动抑制效果（推进剂液位高度635 mm）[6]Fig. 7 Effect of open-loop and closed-loop control on sloshing suppression of propellant （propellant height is 635 mm）[6]

Ohashi[13]、Himeno[14]、Aoki[15]等对推进剂晃动问题进行了研究。Ohashi[13]采用试验和数值模拟的方法研究了带环形隔板贮箱中，低温推进剂晃动对贮箱压力的影响，结果表明：隔板的贮箱在某些情况下，可能会加剧压力的下降过程。 Himeno[14]和Aoki[15]采用试验和数值模拟的方法对垂直起降过程与可回收火箭贮箱内推进剂的流动过程进行了研究，结果表明：防晃板的贮箱可有效抑制箱内推进剂的晃动过程。

上述资料表明，回收过程中推进剂的可靠有效管理对回收成败有至关重要的影响。回收子级再入过程和垂直降落过程中，贮箱外部气动加热明显，虽有绝热层保护，但仍有一部分热量可通过贮箱侧壁与箱内推进剂进行换热，致使推进剂温度升高，进而发动机正常工作需要的增压气体压力进一步提高。然而，若贮箱内部推进剂处于微重力环境，稍有干扰，贮箱内的推进剂即会剧烈晃动，如推进剂晃动至贮箱顶部，可能造成贮箱异常排气或推进剂输送管路夹气，并且推进剂在箱内可能与气枕剧烈掺混，导致推进剂温度进一步升高而气枕压力降低，使箱内推进剂无法有效供应发动机，无法保证发动机的可靠启动和正常工作（图8及图9为有无防晃板两种情况下，箱内推进剂在垂直起降过程中的晃动情况）。

图8 不带防晃板箱内推进剂在垂直起降过程中的晃动情况（推进剂FC-77，横向加速度与轴向加速度比为0.8∶1.0）[15]Fig. 8 The sloshing of propellant without baffle during vertical take-off and landing（propellant FC-77, the ratio of lateral acceleration to axial acceleration is 0.8∶1.0）[15]

图9 带防晃板箱内推进剂的晃动情况（推进剂FC-77，横向加速度与轴向加速度比为0.8∶1.0） [15]Fig. 9 The sloshing of propellant wirh baffle （propellant FC-77, the ratio of lateral acceleration to axial acceleration is 0.8∶1.0）[15]

2 计算模型

本文以三维椭球底圆柱贮箱为例，进行垂直着落过程贮箱内低温推进剂流动过程的仿真计算，贮箱结构如图10所示，网格选用结构化网格。

图10 贮箱结构及网格化分Fig. 10 Tank structure and meshing

2.1 控制方程

假设流体不可压缩，控制方程由不可压连续方程、动量方程、体积输送方程组成，公式如下

上述5个方程分别为流体连续性方程、动量守恒方程、体积函数方程、密度和黏性体积加权方程。其中：V、P、ρ、μ、f、fs、F分别为流场的速度、压力、密度、黏度、体积力、表面张力和网格单元内流体体积分数；下标L，v分别表示液相、气相，本文推进剂选用液氧，ρL＞ 10ρv，可以仅考虑液相的影响（ρv＝ 0，μv＝ 0）。

在飞行过程中，贮箱内推进剂和增压气体间的自由界面用VOF方法进行追踪和重构，计算过程不考虑推进剂与气枕换热影响，即不考虑推进剂相变。采用CSF（Continuous Surface Force）模型[8]考虑表面张力的影响，将表面张力视为界面单元所受的体积力，如式（8）所示

其中：k、 σ 分别为界面曲率和表面张力系数；n为界面法向量。表面张力较大时容易导致数值计算发散，数值求解时必须满足稳定性准则

2.2 初边值条件

假设贮箱壁面绝热，流体和壁面完全浸润，速度无滑移V=0。推进剂选用低温推进剂液氧，初始气体压力选取1个标准大气压，初始温度90 K，密度1 142 kg/m3，初始充液无量纲高度0.34，表面张力系数取0.013 2 kg/s2。

2.3 气液界面稳定条件

飞行器在轨道飞行过程中，不可避免会遇到各种干扰加速度的作用，如姿控系统工作和低地轨道的空气阻力作用等，均有可能破坏贮箱内液气界面的稳定。研究表明，微重力环境下，贮箱内的液体动力学特性与常重力相比有显著区别，微重力环境下，贮箱内液体的表面张力不可忽略，其对箱内液体的形状和位置起主导作用。早期的理论分析和试验研究表明，有加速度作用时液气界面的稳定条件可用无量纲Bo表示[9]

其中：邦德数Bo为表征表面张力影响的无量纲参数；ρ为液体密度；σ为表面张力；g为重力加速度；R为容器半径。

根据国内外飞行经验，在有效控制飞行器姿态及推进剂晃动的基础上，Bo不小于400～500时，可维持界面稳定，根据本文条件计算出维持界面稳定需轴向加速度需不小于1.14 × 10-4（无量纲）。

2.4 晃幅放大条件

“半人马座号”（Centaurus）、“宇宙神”（Atlas ）SLV-3C第4次飞行（AC-4）失利表明，在推进剂管理过程中，不能仅单纯考虑邦德数，还必须考虑施加给液体的能量，即对推进剂有干扰作用的外界影响因素。主发动机关机、发动机再起动、横向干扰带来的推进剂晃动、贮箱保险阀排气形成的不平衡推力等因素均有可能对贮箱内的推进剂添加能量，这些干扰因素均需在飞行过程中得到有效控制，避免发动机关机时由加速度突然减小引起的晃动放大。

假定关机前后的液体晃动幅值为hh，hl；相应加速度为ah，al；关机时晃动波的相位角ϕ；由关机前后的能量守恒可以导出[10]

如果关机时ϕ= 0 或180°，即平液面，那么放大系数K为

AC-4 主发动机关机时加速度2.44g，滑行段加速度3 × 10-4g，放大系数K= 90，关机前晃动幅值为0.12 m，经90倍放大，滑行段的最大晃幅达11 m，而关机时液面距箱顶仅4 m，造成液氢从排气口排出，最终导致任务的失败[11-12]。

本文在计算过程中，给定各工况下作用于流体上的3个方向加速度随时间的变化过程，并将其转化为体积力代入控制方程（2）中，根据计算结果分析加速度变化对推进剂晃幅的影响分析。

3 典型垂直着陆过程推进剂流动过程仿真

以飞行器在原地垂直起飞，飞行一段时间后部分发动机关机，然后利用自重与推力的平衡过程实现飞行器的垂直平稳着陆过程为研究对象。典型飞行器起降过程主要特征点参数和加速度曲线（无量纲化）详见图11和表1，其中–z方向为起飞方向。

图11 垂直起降过程加速度曲线Fig. 11 Acceleration curve of VTVL

表1 垂直起降过程主要特征点参数（无量纲化）Table 1 Characteristic parameters of VTVL(dimensionless)

在飞行器飞行至25.31 s时刻，关掉部分发动机，此时轴向过载（z向，下同）由–0.92下降至–0.49，降幅约50%，如同时引入横向干扰（x向，下同），则可能会放大该横向干扰，从而导致贮箱内推进剂晃动。本文在表1飞行特征点的基础上，对垂直起降过程中贮箱内推进剂流动特性进行分析，重点通过设计25.31 s部分发动机关机时刻的横向过载与轴向过载变化，研究轴向过载大幅变化过程中，横向干扰对箱内推进剂流动特性的影响。

3.1 横向干扰对垂直起降过程箱内推进剂流动特性的影响

横向干扰反映贮箱在飞行过程中所受的横向外力，进而引起箱内推进剂晃动。根据2.3节“气液界面稳定条件”计算出维持界面稳定所需轴向加速度应不小于1.14 × 10–4（无量纲），以此作为25.31 s时刻的轴向过载，同时在25.31 s时，分别向箱内推进剂施加无量纲12.13 × 10–4和12.13 × 10–4的横向过载，分析相同轴向过载变化情况下，不同幅度横向干扰对箱内推进剂流动特性的影响过程。不考虑Y方向加速度影响，25～27 s间横向过载曲线和轴向过载变化情况详见表2。

表2 横向干扰随时间的变化Table 2 Variation of lateral interference with time

计算结果详见图12～13，反映了在横向干扰作用下，贮箱内推进剂的流动特性。从图12可以看出，25.31 s引入12.13 × 10–4的横向过载，推进剂液面虽有小幅晃动，但液面基本保持稳定，与2.3节提到的静稳定条件（Bo不小于400～500）预测结果一致。

图12 横向干扰对箱内推进剂流动特性影响（工况1，无量纲ax = 12.13 × 10–4）Fig. 12 Influence of lateral overload on propellant(Case 1, dimensionless ax=12.13× 10–4)

然而，若在25.31 s将横向干扰放大至10倍（工况2），液面出现了较大幅度的晃动见图13。晃幅在25.31～27 s时轴向过载剧烈变化过程中，出现放大的趋势。25.31～30.8 s时间内，箱内液面晃动明显，随着27 s后轴向推力增大，晃动逐渐被抑制，此间推进剂一直沉底，可以满足发动机的启动要求。但是推进剂在晃动过程中，会加速与气枕的换热，由此会带来贮箱压力不满足发动机最低泵入口压力要求，若不加以控制，可能会导致发动机无法正常启动等问题，因此在垂直起降过程中应尽量保证飞行器姿态，避免出现大幅横向干扰等作用。

图13 横向干扰对箱内推进剂流动特性影响（工况2，无量纲ax = 12.13 × 10–3）Fig. 13 Influence of lateral overload on propellant(Case 2, dimensionless ax=12.13× 10–3)

推进剂液面晃幅详见图14～15，从图中可以看出，横向干扰放大10倍后，尽管此时轴向过载的施加可以克服表面张力的影响，但是液面晃幅明显增大，无法使推进剂液面保持静平衡。

图14 液面晃幅随位置和时间变化（无量纲ax = 12.13 × 10–4，X向）Fig. 14 Variation of sloshing amplitude with position and time(dimensionless ax=12.13 × 10–4, X direction)

图15 液面晃幅随位置和时间变化（无量纲ax = 12.13 × 10–3，X向）Fig. 15 Variation of sloshing amplitude with position and time(dimensionless ax=12.13× 10–3, X direction)

对此采用2.4节“晃幅放大条件”进行分析。弗劳德数Fr（Froude）数反映流体惯性力与重力效应之比，当Fr数大于1，表明惯性力对流动的影响大于重力影响。对比两种工况的Fr数（Fr=v2/gr）见图16。在25～27 s时，横向干扰增大后，Fr数相应增大，工况1的Fr数最大约0.4，惯性力略小于重力，量级基本相当；工况2的Fr数最大约2.6，惯性力大于重力。Fr数表明若外界带入推进剂更多的能量，且这些能量产生的惯性作用与重力效应相比无法忽略，则会致使推进剂不稳定。根据公式（8）可知，作用于推进剂的横向过载所引入的推进剂晃幅，在轴向过载突降过程中（1.0下降至1.14 × 10–4）会被放大，放大比例与变化前后轴向过载比值及相应晃动波相位的均方根成正比。

图16 Froude数随时间变化情况Fig. 16 Variation of Froude with time

3.2 轴向过载变化对垂直起降过程箱内推进剂流动特性的影响

轴向过载反映贮箱在飞行过程中所承受的轴向推力。在3.1节分析结果的基础上，研究相同横向干扰，不同幅度轴向过载变化对箱内推进剂流动特性的影响。在25.31 s时，向给箱内推进剂施加12.13 × 10–3（无量纲）的横向干扰，分析5种轴向过载工况下箱内推进剂的流动特性。不考虑Y方向加速度影响，25～27 s时间段内的横向过载和轴向过载变化情况详见表3。

表3 轴向过载随时间的变化Table 3 Variation of axial overload with time

计算结果详见图17～19。图17反映不同轴向过载工况下，箱内推进剂液面晃幅随横向位置（X）和轴向过载的变化情况。从图中可以看出，横向干扰一定的情况下，25.31～27 s时间段内随着轴向过载的增加，液面晃幅逐渐减小。其中，工况1～工况3晃幅大小基本一致，箱内推进剂的流动分布与3.1节中的图12较为一致。工况4～5的推进剂流动特性见图18和图19。从图中可以看出，随着25～27 s内轴向过载变化程度的减小，相同横向干扰引起的推进剂特性存在一定的差别，轴向过载变化越小，推进剂越稳定。然而，根据本文的计算结果，当变化后的轴向过载小于–1.14 × 10–1时，箱内推进剂的晃动特性并无较大差别，但是当轴向过载大于1.14 × 10–1时，箱内推进剂的晃幅大幅减小，从侧面表明轴向过载越大，箱内推进剂越不容易大幅晃动。

图17 液面晃幅随位置和轴向过载变化情况（t = 26.8 s）Fig. 17 Variation of sloshing amplitude with position and axial overload(t=26.8 s)

图18 轴向过载变化对箱内推进剂流动特性影响（工况4，无量纲ax = –12.13 × 10–3，az = –1.14 × 10–1）Fig. 18 Influence of axial overload on propellant(Case 4, dimensionless ax=12.13× 10–3， az = –1.14 × 10–1）

图19 轴向过载变化对箱内推进剂流动特性影响（工况5，无量纲ax = –12.13 × 10–3，az = –0.92）Fig. 19 Influence of axial overload on propellant(Case 5, dimensionless ax = –12.13 × 10–3， az = –0.92)

图20反映了不同轴向过载下，Fr数随时间的变化情况。从图中可以看出，在横向干扰一定的情况下，随着轴向过载的增加，Fr数逐渐减小。与图15对比分析可知，Fr数可以解释贮箱内出现大幅晃动的原因，并且在一定程度上其数值可以反映推进剂是否会出现大幅晃动，但无法预测晃幅的大小。根据本文所研究的表2和表3共6种工况算例可知，横向过载Fr数需小于0.01，才能维持箱内液面稳定。但这个数并不是绝对的，例如图16中在横向过载缩小10倍，Fr数约0.4即可维持液面平稳，仅有小幅晃动。

图20 Froude数随时间变化情况Fig. 20 Variation of Froude with time

4 结 论

本文研究运载器垂直着陆过程贮箱内推进剂流动行为特性，并开展部分发动机关机时刻的横向过载与轴向过载变化对箱内推进剂流体行为特性的影响分析。本文研究结果表明：轴向过载变化一定的情况下，箱内推进剂的晃动幅值与横向干扰幅值相关；较大的横向干扰，会给箱内推进剂施加较大的横向动能，若该动能引起的惯性力与轴向过载效应相比可以忽略，则箱内推进剂满足静稳定条件可保持液面稳定（小幅晃动，见表2中工况1），否则箱内推进剂会出现大幅晃动，在部分发动机关机或其它轴向力突然减小的情况下，推进剂晃幅会被大幅放大，从而加速推进剂与气枕的换热，使得气枕压力大幅下降，若不加以控制，可能会影响发动机的正常启动和工作。因此在垂直起降过程中应尽量保证飞行器姿态，避免出现大幅横向干扰。