导读:6月25日,嫦娥六号一路狂飙,跑了38万公里,带着2公斤月球土特产,顶着3000度的高温,以11.2公里每秒的速度冲进大气层,玩了31次“大气层跳水”,14点07分最后在内蒙古四子王旗精准落地,可喜可贺。
从视频画面来看,嫦娥六号的返回舱主伞是单点悬挂式,在空中下落的时候就不是竖直姿态,而是侧倾着左右摇晃。无独有偶,4月30日,我国的神舟十七号载人飞船,通过和空间站对接,在从空间站组合体撤离后,也是经过一系列的制动减速、再入大气层、打开降落伞,然后降落到着陆点。
但是,很多朋友也看到过美国载人龙飞船返航时的画面,在着陆的最后阶段,载人龙飞船会打开4个降落伞,缓缓下降后降落到大海。因此,有人就问为什么美国的载人飞船要使用3-4个伞呢?群伞仿真告诉你真相。
下面我们来看一篇由南京航空航天大学航空宇航学院连亮同学公开发表的一篇《基于 ALE 方法的群伞稳降阶段的数值模拟》学术论文,论文采用任意拉格朗日—欧拉法对美国“阿波罗”飞船降落伞回收系统的稳降阶段进行了数值模拟。通过建立结构动力 学与计算流体力学的耦合模型,研究了 3 具环帆伞的流固耦合过程,并模拟了降落伞的外形变化和伞衣 内外流场的变化。流固耦合计算结果表明,各伞最终稳定在 20°左右,环帆伞的群伞效率因数约为 0.93, 比 C-9 伞大。文章的分析结果对中国探月工程回收系统的主伞设计具有参考意义。
1、写在文前
在航天回收和重型空投等领域,单具伞已经快达到实际使用极限。群伞可以替代大尺寸单伞,且比 单伞更有优势。当群伞中的一具伞失效时,在剩余几具伞的作用下回收系统的着陆速度要比单个主伞失 效后的自由落体安全得多。
早在 20 世纪 60~70 年代,美国就在“阿波罗”探月计划中对群伞的有关技术做了大量研究,但这 些研究大多是空投试验和风洞试验。之后随着计算机技术的迅猛发展,对降落伞工作过程的仿真也取得 了不少成果。文献[1]使用 LS-DYNA 对 1 具伞的充气过程进行了模拟。文献[2]采用有限质量法对 1 具伞 的开伞过程进行了仿真。文献[3]使用 LS-DYNA 对美国先进战术降落伞系统(advanced tactical parachute system, ATPS)备选伞进行了数值分析。文献[4]采用流固耦合方法对超声速来流条件下降落伞的流场进行 了分析。文献[5]采用任意拉格朗日—欧拉法(ArbitraryLagrangian-Eulerian method, ALE)模拟了无限质 量情况下不同风速的开伞情况。文献[6]采用流固耦合的计算方法对带牵顶伞的主伞的充气过程进行了数 值模拟。以上这些文献拓宽了研究者的视野,然而都只针对某种类型的单伞,至今国内外对群伞的数值 模拟还很少。文献[7]首次用 CFD 方法对 3 具 C-9 伞组成的群伞在充满后的降落阶段进行了数值计算。但是该文献发表于 1997 年,方法已显陈旧:该模型所使用的是刚性外形,与实际中降落伞是柔性织物不 符合;该文只是针对 3 具 C-9 伞在几个典型角度时做了 CFD 分析,并没有真正模拟群伞相对竖直轴的摆 动过程。
本文以“阿波罗”主伞系统为研究对象,尝试采用 ALE 方法模拟出 3 具环帆伞在稳定下降阶段的工 作过程,计算出群伞系统在稳定降落阶段的阻力,分析出组成群伞的各单伞之间的流场干扰情况,判断 出群伞最终稳定角等。
2、数学建模
模拟群伞的工作过程难度很大,涉及到空气动力学、纺织材料结构力等学诸多复杂的力学问题。流 场取决于伞衣外形,且流场又影响伞衣外形。必须将降落伞的结构受力和周围的气动力相耦合才有可能 实现准确模拟。众所周知,Lagrange 方法的优势在于有限元节点即为物质点,物质不会在单元间发生流 动,材料变形后的自由表面可被网格边界捕捉到,但对于大变形问题经常因网格过度扭曲而增加数值计 算的难度。该方法多用于固体结构的应力应变分析。Euler 方法的优势在于网格始终是固定不动的,有限 元节点即为空间点,材料在固定网格中流动,可用于分析变形程度较大的问题,但为了精确捕捉到固体 材料的变形,需要将网格划分得非常精细,这大大增加了计算成本。本文对 Lagrange 方法和 Euler 方法 两者取长补短,采用 ALE 算法对群伞稳降阶段进行流固耦合计算,其控制方程主要包括流场控制方程和 网格控制方程,通过对这些控制方程的耦合计算来实现流固耦合。
(1)流场控制方程
本文采用 CAD 三维建模技术做出 3 具环帆伞的初始外形,该外形符合“阿波罗”飞船所用环帆伞的实际制造尺寸,近似于充满状态。“阿波罗”飞船降落伞回收系统如图 1 所示。其主要参数如下:
Fig. 1 The parachute recovery system ofApollo spacecraft
每具主伞的名义直径为 25.5m,名义面积为 508.8m2,有68个伞衣幅。每个伞衣幅上有4个环、8个帆。群伞载值量5897kg,着陆速度 9.57m/s[8]。充满状态单伞的阻力特征为390m2 [9]。伞绳长度为 37.74m,吊带长度为 2.37m,伞衣质量为38.73kg,伞绳质量为 13.24kg,吊带加连接绳质量为 5.67kg[10]。限于篇幅,环帆伞的具体尺寸等数据不再多做说明。
物伞系统坐标系和单伞坐标系,如图 2 所示。单伞坐标系的原点在伞衣底边中心,Ys 轴即伞衣对称轴(指向伞衣顶),Xs 轴在伞衣底边平面内,Zs 与 Xs、Ys 组成右手坐标系。物伞系统坐标系的坐标原点 Od 在 3 伞连接带的汇交点,Yd 轴垂直向上,Zd 与 1 号伞(伞的标号见图 3)的 Zs轴同向,Xd 与 Yd、Zd 组成右手坐标系,OdXdYd 平面与 1 号伞的 OsXsYs 平面共面,伞 2 和伞 3 关于 OdXdYd 平面对称。3 具伞的 Xs 轴在物伞系统坐标系的 OdXdZd 平面上的投影与 Xd 轴分别成 0°、120°和–120°。3 伞相接触时的位置示意图,如图3。图2中的 角即3伞相接触时 Ys 轴与 Yd 轴的夹角。
本文基于两点假设。
1)降落伞属于拖曳式气动减速装置,它的气动特性必然受前方物体尾流的影响, 降落伞在强湍流气流中的动压比自由来流中的动压要低。对于单伞系统来说,降落伞总是处于尾流之中。但是对于群伞系统来说,由于前体体积较小且伞绳较长,各伞相对于竖直轴有一定夹角,大多数时候降落伞都不在前体尾流中而在自由来流中。此外,亚声速来流条件下,前体尾流的影响并不明显,不会像 超声速来流那样前体产生脱体激波。所以,本文暂不考虑前体尾流对降落伞流场的影响,只计算 3 具伞 在自由来流中的流场。
2)降落伞伞衣的厚度很小,远小于伞衣的几何尺度。因此本文忽略了伞衣厚度的 影响,在 CAD 建模时将伞衣作为片体处理。
稳降过程是伞衣充满后到系统合力为零的过程。稳降过程不是恒速过程,也不是恒加速度过程,是非平衡态到平衡态的过程。在伞衣充满瞬间,由于惯性作用,系统还来不及减速到稳降速度,此时系统阻力大于系统重力,从而产生减速度,群伞的阻力也随着速度的减小而减小。降落伞持续减速,因为阻尼作用降落伞最终会在气动力和重力对连接带汇交点合力矩为零的那个角度停止摆动,从而整个回收系统以稳降速度匀速下降。
稳降速度是指气动力合力在垂直方向上的分力与物伞系统重力相等时的物伞系统速度。本文模拟的是“无限质量”条件下群伞的充气特性,充气时间为15s。流场入口为速度入口边界,来流速度取着陆速度9.57m/s,其余流场边界采用无反射边界。仿真模型的网格划分:伞绳和加强带有3´104 个梁单元(beam),伞衣有4´104 个壳单元(shell),流场有7.5´105 个六
面体固体单元(solid)。有限元模型如图4 所示。
图 4 群伞的有限元模型Fig.
5、计算结果与分析
将上述方法运用于群伞的仿真研究。入口边界条件同样是标准大气压下着陆时的动压56Pa,3 具伞与对称轴的初始夹角为30°。采用ALE 算法计算出3 具环帆伞与对称轴夹角随时间变化如图7,3 具伞相对竖直轴来回摆动,3 伞的摆动幅度大致相同,且最终都会稳定在20°左右。阻力曲线如图8 所示,气动阻力最终稳定在70 000N 左右。群伞效率因数定义为群伞阻力系数与单伞阻力系数乘以伞的数目的比率。因为尾流的相互作用,群伞效率因数通常小于1.0。文献[6]中3 具C-9 伞组成的群伞在25°至40°范围内的效率因数介于0.82 到0.92 之间,且群伞效率因数随着夹角增大而减少,在稳定角33°时群伞效率因数约为0.86。本文算得在稳定角20°时环帆伞组成的群伞效率因数约为0.93,要比C-9 伞的更大。其原因一方面可能是伞型不同,另一方面在于环帆伞的稳定角要小于C-9 伞。
对比图7 和图8 可以发现,3 具环帆伞在5°(三者相碰)到30°之间摆动。在大约8s 时,3 伞相距最远,此时群伞阻力也相应处于波谷。在大约11s 时,3 伞再次相互靠拢,此时群伞阻力也相应处于波峰。
图9 和图10 是用CFD 方法计算出来的C-9 伞和环帆伞在25°倾角时的流场图谱。从这两图可以较为清晰地看出2 种不同伞型所组成的群伞之间的差别。环帆伞比平面圆形伞多了很多结构开孔。因为环帆伞每块帆片的前缘都比同一幅中与它相邻的下一块帆片的后缘要长,在降落伞内外压差的作用下,上下相邻的帆片之间会形成月牙形缝隙,气流从这些月牙形缝隙流出,且方向斜向下。这数百股向下气流有效地阻止了已经分离出去的气流再附到伞衣外表面上。同时环形缝中流出的气流加剧了尾流区的气流掺混,使得环帆伞的尾流区要比平面圆形伞的更大。这在某种程度上解释了环帆伞的群伞效率因数更高的原因。
图11 显示了第4s 时OdXdYd 平面的流场图,即在降落伞的气动力矩和重力矩共同作用下,3 具伞由30°摆动到相互靠拢时的流场图。图12 显示了第8s 时,3 具伞又分离时的流场图。图13 显示了第11s时,3 具伞再次相互靠拢时的流场图。可以看出,当3 伞相互靠近时它们的尾流有强烈的相互干扰和掺混现象,当3 伞相互分离时尾流的相互干扰现象会减弱。对比这3 个图,可见降落伞的外形有微小变化,
并非固定不变,这是由于不同状态时降落伞内外压差不同导致的;且每具伞的外形并非关于Ys 轴对称,这是因为有倾角时伞衣各处内外压差各不相同;而这些都是单纯的CFD 方法所不能模拟的。此外,由图13 还可以看出,3 具伞相对于连接带汇交点发生了转动,2 号伞和3 号伞已经不关于OdXdYd 平面对称了。这也印证了实际工作过程中在广义坐标内整个群伞系统有15 个自由度,其中单具伞相对汇交点有3 个转动自由度。
群伞回收系统在着陆时受到的气动阻力直接决定了该回收系统的负载质量。“阿波罗”飞船加上3 具环帆伞的总质量约为60000N。降落伞的工作过程是伞衣结构变形与伞衣周围流场变化相互耦合的复杂过程,因而数值模拟难免会有误差。加之本文所采用的网格密度还不够细,所以算出来的70 000N 阻力要大于负载质量。群伞的计算误差约为16%,与单伞的计算误差是一致的。用本文算出群伞效率因数0.93 乘以单具环帆伞的阻力特征390m2 再乘以动压56Pa 与伞的个数3,算得阻力为60 933N,与实际负载非常接近。说明本文采用ALE 方法模拟群伞在稳降阶段的工作过程是可行的。
本文首次采用ALE 方法描述了3具全尺寸的环帆伞在稳降阶段的工作特性。主要结论如下:1)预测该群伞最终的稳定角度约为20°; 2)在稳定角附近,该群伞系统的效率因数约为0.93。 本文首次将流固耦合的方法引入到群伞的研究中,虽然没能模拟群伞的整个开伞过程,但对下一步群伞的仿真研究具有启发意义,并对中国探月工程回收系统的主伞选型具有指导意义。
我们的神舟载人飞船采用单个伞降落,并不是说我们的技术不如美国,我们其实也在发展群伞降落的技术,也取得了成果。不过现在不会使用到神舟飞船上,因为如果采用新的伞降方式,需要对飞船进行大改,未来新一代载人飞船会使用到群伞降落的方式。
最后,向读者朋友推荐一个视频教程《LS-Dyna流固耦合仿真技术应用15讲-详解水射流破岩、降落伞开伞、无人机投伞、重物入水和水箱晃动分析》欢迎读者朋友试看和批评指正。
[10] Robert B W. ApolloExperience Report Earth Landing System[R]. NASA TN D-7437, 1973.