目前国内外新一代大载重返回式载人飞船通常采用群伞技术实现回收着陆^[1-3]，包括美国多用途乘员飞行器(multi-purpose crew vehicle，MPCV)^[4-5]、乘员航天运输飞船^[6-8]、龙飞船(Dragon)^[9-10]等空间运输系统。我国新一代载人飞船采用全新气动外型返回舱，可全面兼顾近地轨道自主飞行、空间站运营和载人深空探测等多项任务目标，具备人员往返、货物运输等多种用途^[11]。

1 国内外研究现状

目前对于主伞包拉出过程，国内外学者通常采用结合经验公式的简化动力学模型进行分析，该方法带来的问题在于无法准确描述伞包拉出过程中柔性吊带的动态拉力载荷以及伞包与舱体的接触摩擦影响。由于影响该过程的力学因素较多，涉及柔性大变形复杂动力学运动过程，采用经典公式无法精准描述一些细节，因此需应用精细化动力学建模手段开展相关分析研究。

包进进等^[12]基于Newton经典力学通过建立减速伞释放后拉主伞包过程的动力学方程，计算了释放减速伞拉主伞包过程的减速伞伞绳、吊带和主伞包拖带的拉力随时间的变化情况，并在此基础上研究了拉力峰值的影响因素；林斌^[13]提出了一种典型的倒拉式回收降落伞伞包的受力情况分析和计算方法，得到了伞包拉力的计算公式，计算结果与实际情况相符；鲁媛媛等^[14]建立了在火星环境下降落伞拉直过程的三维动力学模型，模型将伞包视为变质量六自由度(six degrees of freedom, 6DoF)刚体，伞绳和伞衣采用质量阻尼弹簧模型模拟，并用动力学模型对“海盗号”第一次气球发射试验进行了仿真，从物伞系统的空间位形、伞绳张力曲线、拉直时间等方面进行对比，验证了模型的正确性；宋旭民等^[15]建立了包括降落伞系统拉直、充气和充满阶段的仿真模型，给出了主要的模型假设、分析方法和数学模型，可对飞船回收全过程进行集成分析；文[16-18]建立了伞舱盖弹射分离后拉伞包过程的动力学模型，分析在此过程中伞舱盖及伞包的动力学特性，以及主要设计参数对其产生的影响。

以上方法虽然为主伞包拉出过程动力学分析提供了崭新的解决思路，但引入了较多简化假设，如减速伞和主伞包在每一瞬时的运动都遵循同一条轨迹，且运动方向与重力方向相同，不考虑主伞包和主伞舱间的摩擦阻力、主伞包受力过程中的柔性变形和减速伞运动过程中的扭矩。以上简化假设不利于深入分析主伞包拉出过程的动力学特性，故需要一种能够对主伞包拉出过程进行精确建模和计算的方法，以兼顾各部分组件的实际运动状态和物理特性。为解决以上问题，本文提出了一种有限元分析与舱伞系统动力学计算相结合的仿真分析方法，开发了降落伞由气动力向动力学模型实时传递的专用算法，以完成主伞包拉出全过程的分析计算。

2 飞船返回舱主伞包拉出物理过程

如图 1所示，目前新一代载人飞船采用降落伞进行气动减速和回收，即首先利用2具减速伞进行初步减速，然后采用由3具主伞组成的群伞进行减速，最后采用缓冲气囊实现着陆缓冲。

飞船在初入大气阶段利用自身气动外形进行减速，到达8 km高空时，启动回收工作程序，其工作过程如图 2所示。1) 通过2个弹伞筒实现减速伞的弹射拉直，并进行第1次减速；2) 减速伞工作一段时间后与飞船分离，分离的同时拉出主伞舱盖和主伞包，当减速伞再次拉直时，将主伞从主伞包中拉出；3) 拉出的主伞充气展开进行第2次减速。在上述过程中，减速伞将主伞包从伞舱内拉出的过程涉及瞬间拉动伞舱盖和伞包2个大质量块，减速伞吊带受力剧烈，是整个工作过程中需要重点考虑的环节。

图 3为飞船回收系统。飞船返回舱在回收着陆过程中存在减速伞分离拉出主伞包关键环节，为了量化分析上述过程中过载、速度和姿态变化等动力学特性，需开展相关的仿真分析计算工作。本文分析了返回舱在该过程中的动力学特性，即通过对返回舱气动减速过程中减速伞拉出主伞包过程建立有限元模型，量化了减速伞拉出主伞包过程中各部件之间的相对运动及受力情况，厘清了各种因素对该过程中舱体动力学特性的影响，从而有效指导了方案设计。

新一代载人飞船主伞组件安装布局在返回舱下部，而主伞的连接吊带连接点则设置在返回舱顶部。减速伞分离后首先需将主伞拉出吊带从返回舱侧壁的防热槽道内拖出，随后将解除约束的主伞舱盖和主伞包一起从返回舱下部拉出，并将连接吊点转换至顶部。以上各部组件的拉出过程连接关系如图 4所示。

以上过程中，减速伞在拉直后，首先拉动解除固定约束的主伞舱盖，该舱盖质量为50 kg；随后将长度为1.8 m的主伞包拖带拉直，并将质量为210 kg的主伞包从伞舱中拉出。以上拉动主伞舱盖和主伞包过程涉及大质量块瞬间运动的过程，产生作用在各级吊带的瞬态冲击载荷，对吊带设计提出了较高要求。由于伞舱为异型伞舱，在主伞包拉出过程中，伞包与主伞舱侧壁间存在复杂的相互作用力。拉出主伞包后主伞拉出吊带拉直，拉出主伞，主伞充气张满完成开伞过程。

3 飞船返回舱主伞包拉出过程动力学

为了精准描述返回舱主伞包从伞舱拉出过程的动力学特性，建立了返回舱、各级吊带、伞绳、主伞包、伞舱盖及防热结构的数值分析刚柔组合体非线性有限元模型，应用显式动力学方法对伞包拉出过程进行分析，并提出了降落伞气动载荷的实时施加方法，发展了基于ABAQUS外部接口的气动力传递与载荷转换技术，建立了气动力在动力学模型中实时引入的路径，分析了舱体速度、防热层拉脱力及舱盖结构质量等关键参数对拉出过程中部件受力及运动特性的影响。分析过程中提出将减速伞伞衣气动力实时加载到动力学模型中，从而实现减速伞由气动力向动力学计算受力边界条件的转换。在此基础上，分析了在多种影响因素下的主伞包拉出过程各项参数的动力学响应，获得了在单因素下不同量级的影响分析结果。

3.1 返回舱及主伞舱盖结构建模

返回舱主体由金属结构组成，在主伞包的拉出过程中可近似为刚体，采用6DoF模型分析其运动及姿态变化，在实际计算过程中将返回舱模型采用4节点壳单元进行离散，表面单元节点与质心节点采用刚性约束连接，刚体表面单元节点为从节点，质心虚拟节点为主节点，刚体质量和惯量施加于主节点；伞舱盖建模方法与返回舱类似。返回舱在动力学分析过程中受到的主要气动载荷是重力及降落伞通过吊带传递的牵引力、接触力。返回舱动力学控制方程表示如下：

$ \boldsymbol{F}=\boldsymbol{F}_{\mathrm{M}}+\sum\limits_{i=1}^N \boldsymbol{F}_i, $

(1)

$ \boldsymbol{M}=\boldsymbol{M}_{\mathrm{M}}+\sum\limits_{i=1}^N \boldsymbol{G}_i \times \boldsymbol{F}_i . $

(2)

其中：F和M分别为作用在返回舱及舱盖上的力和力矩；F_M为主节点上的力矢量；F_i为编号为i的从节点上的力矢量，i=1, 2, …, N；M_M和M_i分别为主节点和从节点上的力矩矢量；G_i为由从节点指向主节点的矢量。

分解后的运动加速度和角加速度分别为

$ \boldsymbol{a}=\boldsymbol{F} / m, $

(3)

$ \boldsymbol{I}_1 \boldsymbol{a}_1=\boldsymbol{M}_1-\left(\boldsymbol{I}_3-\boldsymbol{I}_2\right) \boldsymbol{\omega}_2 \boldsymbol{\omega}_3, $

(4)

$ \boldsymbol{I}_2 \boldsymbol{a}_2=\boldsymbol{M}_2-\left(\boldsymbol{I}_1-\boldsymbol{I}_3\right) \boldsymbol{\omega}_1 \boldsymbol{\omega}_3, $

(5)

$ \boldsymbol{I}_3 \boldsymbol{a}_3=\boldsymbol{M}_3-\left(\boldsymbol{I}_2-\boldsymbol{I}_1\right) \boldsymbol{\omega}_1 \boldsymbol{\omega}_2 . $

(6)

其中：a和m分别为舱体的加速度和质量，I₁、I₂、I₃分别为转动惯量张量矩阵第1、2、3列元素组成的列向量，a₁、a₂、a₃分别为刚体体轴局部坐标系内角加速度张量矩阵第1、2、3列元素组成的列向量，ω₁、ω₂、ω₃分别为刚体体轴局部坐标系内角速度张量矩阵第1、2、3列元素组成的列向量，M₁、M₂、M₃分别为刚体体轴局部坐标系内力矩张量矩阵第1、2、3列元素组成的列向量。

3.2 吊带及伞绳结构建模

降落伞的各类吊带、伞包拖带和伞绳结构均以一维线单元进行离散。根据其只受拉伸力的特性，以及应力应变的变化关系，定义了非线弹性索单元，以更加真实地模拟绳带部件的力学本构模型，单元拉伸力F_s可表示为

$ F_{\mathrm{s}}=f\left(l-l_0\right)+c \frac{\mathrm{d} l}{\mathrm{~d} t} . $

(7)

其中：f(l)为绳带拉伸非线性应变函数，由材料拉伸试验曲线通过多项式拟合方式得出；l为伞绳长，l₀为初始伞绳长；c为阻尼系数；$\frac{\mathrm{d} l}{\mathrm{~d} t}$为伸长变形率，t为时间。

3.3 飞船返回舱主伞包拉出过程动力学建模

首先根据飞船返回舱系统的外形及几何参数建立了全尺寸分析模型。如图 5所示，利用三维建模软件建立了不含伞衣的返回舱主伞包拉出全尺寸几何模型，随后建立了返回舱、伞包、吊带、拖带和伞绳的有限元模型，考虑吊带及拖带与舱体的接触力，系统的动力学控制方程可表示如下：

$ \boldsymbol{Q} \frac{\mathrm{d} \boldsymbol{v}}{\mathrm{d} t}=\left\{\boldsymbol{F}_{\text {ext }}\right\}-\left\{\boldsymbol{F}_{\text {int }}\right\}+\left\{\boldsymbol{F}_{\mathrm{hgr}}\right\}+\left\{\boldsymbol{F}_{\text {cont }}\right\} . $

(8)

其中：Q为质量矩阵，v为物质运动速度，F_ext、F_int、F_hgr和F_cont分别为外部力、内力、沙漏力以及接触力。为简化模型，将舱体的外表面设计为空壳，主伞包为实体，伞包拖带及连接方式按照实际连接形式进行简化，伞舱和防热层具体位置参照实际安装位置建立斜拉主伞包分析模型。

接触力的受力分析模型采用主从节点方式，舱体表面为主面，壳单元节点作为主节点，绳带单元节点作为从节点。

如图 6所示，在几何模型中，主伞包分3组，左右并列紧挨排布，伞包拖带从伞包顶面边缘处引出，每个伞包有16根拖带，拖带从伞舱盖方向伸出后向下对折放置，按照拖带实际长度，每个伞包的16根拖带再汇成一股拖带，继续对折排布，再通过伞舱盖上的固定点固定在伞舱盖上。

将图 6中建好的几何模型导入分析软件，建立如图 7所示的返回舱主伞包拉出过程有限元分析模型，并分别赋予材料属性。在该分析模型中：返回舱舱体采用Shell单元并进行刚性约束，赋予集中质量；伞绳、主伞吊带和垂挂吊带等简化为一维单元，采用索单元模拟。以上单元组件的连接采用绑定模式模拟，建立部组件间的作用关系。防热层与返回舱之间采用失效单元连接，用以模拟可破坏的拉脱力，具体破坏力的数值通过设定材料的本构模型获得。对上述材料及结构赋予相应的材料属性，并为各部组件赋予分析所需要的各类初始条件，包括初始速度、重力加速度、姿态和位置等。坐标系的原点O为拉出主伞包初始时刻主伞包拉出吊带与舱体头罩中心的位置点，Z轴为返回舱轴线方向，Y轴为头罩平面内远离伞舱的方向，根据右手定则确定X轴。

3.4 飞船返回舱主伞包拉出过程气动力加载设计

在主伞包被拉出伞舱的过程中，上一级降落伞已处于完全张满状态，伞衣外形变化甚微，并且降落伞中心轴线与吊带方向近似重合，且在将降落伞气动力转换为动力学分析模型的受力边界时，忽略降落伞横向气动力，仅考虑轴向气动力，即气动阻力。在阻力计算过程中，由于伞绳呈周向均布式分布，垂直于伞中轴线方向的气动力合力为0，沿中轴线方向的气动力相较于伞衣承受的气动力量级较小，因此分析中仅考虑伞衣受到的气动力F_f。

在计算降落伞流固耦合过程中，网格数量是影响计算效率的关键因素，而网格以降落伞伞衣网格为主，由于本文研究重点在于整个主伞包拉出的动力学响应过程，而非伞衣的局部受力，因此采用等效模拟的方式对伞衣受力进行了模拟，即将伞衣受到的气动力通过气动载荷程序算法加载实现，分析中伞衣受到的气动力可表示为

$ \boldsymbol{F}_{\mathrm{f}}=\frac{1}{2} \rho C_{\mathrm{s}} \boldsymbol{v}^2 . $

(9)

其中：ρ为主伞包拉出时的空气密度, C_s为阻力面积。

气动力的大小会直接影响伞衣的速度，使其发生变化，而伞衣的速度又影响气动力的大小，二者之间的作用是相互耦合的。分析中采用ABAQUS自带的子程序接口VUAMP来实现伞衣气动力幅值的实时动态控制，即当伞衣下降时，速度方向沿轴线方向向下，伞衣受到气动力的作用；当伞衣加速下降时，除受气动力作用外，还受附加质量施加给伞衣的阻力，该阻力可认为是降落伞在理想流体中作变速运动时，受到的惯性阻力作用，该惯性阻力可表示为

$ \boldsymbol{F}_{\mathrm{m}}=k m_0 \boldsymbol{a}_{\mathrm{u}}. $

(10)

其中：k为修正系数，取值范围为0~1；m₀为附加质量；a_u为伞衣加速度。

图 8为气动力加载有限元分析模型，图中g为地球重力加速度，分析中通过子程序实现连接伞衣和舱体的伞绳、连接带、吊带等变形对系统受力分布的量化，以更真实地反映伞衣气动力和速度相互耦合的关系，通过对算法的修正使气动力的计算值更加精确。分析中首先设定有限元计算边界条件，如重力、质量和接触等，随后开始增量步前处理和计算，并通过接口子程序不断迭代，即首先提取伞的v，随后判断伞的速度方向，并计算F_f，基于a_u和伞衣轴向方向向量n判断是否叠加附加质量。伞衣载荷的计算主要利用伞衣阻力特征变化和附加质量的方法，建立物伞系统动力学方程，求解开伞充气过程的载荷变化，在整个主伞包拉出过程中，不断迭代计算至稳定状态^[19-21]。

3.5 飞船返回舱主伞包拉出过程计算

基于气动力加载有限元模型和设置开展了计算分析，包括标准工况下过载、速度和姿态变化等动力学响应特性，此外，围绕不同的初始速度、防热层拉脱力、舱盖质量等工况进行了对比分析，最大限度地获得了各项指标参数的影响，量化了返回舱各部分组件的动力学响应特性。

1) 在标准工况下计算。

首先，分析标准工况下的主伞包拉出过程，选用工况和参数为：2具减速伞同时工作，拉出时刻初始速度为60 m/s，防热层拉脱力为5 kN，舱盖质量为50 kg，主伞包质量为210 kg。图 9和10分别为飞船返回舱拉出主伞包过程中主吊带载荷F_main以及伞包拖带载荷F_pack随时间变化的曲线。以2条吊带初始相对夹角作矢量合成计算可得到主吊带总吊带载荷和伞包汇股后总体拖带载荷，拉出主伞包过程中各吊带和拖带受力较为接近。减速伞先后拉开防热层、伞舱盖和3个主伞包，在此过程中主吊带载荷随着拉动伞舱盖和伞包先迅速上升，主吊带和总吊带载荷与伞包拖带和汇股后总体拖带载荷峰值分别为105.05和291.62 kN，峰值出现时刻集中在0.13 s前后。随着伞包相继脱出，载荷立即下降，之后再次出现峰值，但主吊带载荷在峰值附近的波动相对平缓。此后，舱体与主伞包分离，载荷最后在较低水平波动。

与之相比，伞包拖带载荷变化更复杂，不同伞包拖带载荷峰值和出现时间较为随机，这与伞包拖带的布置形式密切相关。3个主伞包被拉出时，不仅受到减速伞的拉力作用，而且同时存在主伞包与舱体、主伞包与伞舱盖以及主伞包之间的相互作用。此外，由于伞包运动姿态较为复杂且随机，16根伞包拖带交替受力，因此载荷峰值出现时间随机，且不同伞包拖带载荷峰值大小有差异。

图 11和12分别为飞船返回舱拉出主伞包过程中返回舱速度V和过载a_o变化曲线，a_o用地球重力加速度g的倍数表示。2具减速伞和返回舱的初始速度均为60 m/s，由于主伞包质量与返回舱质量(6 730 kg)相差较大，且主伞包拉出后与返回舱体完全分离，返回舱与主伞包之间没有连接，因此在主伞包拉出过程中飞船返回舱做匀加速运动，速度变化较小，仅在主伞包与返回舱相互作用时产生过载；而2具减速伞受拉出主伞包过程的影响，速度波动较大，初始时由气动力的作用速度降低，当伞舱盖拉出时速度突然上升，接着在拉出主伞包瞬间速度再次突升。拉出主伞包过程中，返回舱过载曲线有2个主要峰值区域，其中，第1个过载峰值对应防热层拉脱及伞舱盖拉出；第2个过载峰值与主伞包拖带载荷峰值相对应，是3个主伞包之间及其与舱体相互作用的结果。

2) 在不同速度工况下计算。

返回舱在不同的海拔高度弹射减速伞，其初始速度有一定的差别，设置拉出主伞包过程的初始速度V₀分别为60、65和76 m/s共3种工况，分析不同初始速度对主伞包拉出过程的影响。图 13和14分别为不同初始速度下飞船返回舱拉出主伞包过程中主吊带载荷和伞包拖带载荷随时间的变化曲线。初始速度越大，主吊带载荷越大。当初始速度变化时，载荷峰值出现时刻的差异较大。初始速度为60、65和76 m/s时，主伞包完全被拉出所需的时间分别为0.20、0.19和0.17 s。3种初始速度下主伞包拉出过程相似，但伞包拉出速度随着初始速度的增加有所加快，导致伞包在脱出过程中相对运动位移增加。与65 m/s工况相比，76 m/s工况主吊带载荷存在明显上升现象，这主要是由于气动力随着速度的平方增加。此外，伞包拖带缠绕形式复杂，同一时间不同工况可能出现不同的运动状态，载荷峰值出现的时间和大小较为随机。

图 15和16分别为在不同初始速度下减速伞及返回舱速度、过载随时间的变化曲线。结果表明：当初始速度变化时，舱体速度变化趋势一致，均作匀加速运动；减速伞速度变化趋势相近，即先减速后加速，伴有多个峰值，不同初始速度工况的减速伞速度峰值大小有明显差异，但出现时刻大致相同。主伞包拉出过程中仅在主伞包与返回舱相互作用时产生过载，因此当初始速度不同时，主伞包运动状态不同，从而对伞包与返回舱相互作用产生影响。当初始速度较大，即76 m/s时，返回舱过载峰值明显大于其他工况，且峰值出现时刻及曲线走向存在明显差异。

3) 防热结构在不同拉脱力下计算。

主伞拉出吊带撕裂防热层结构并将主伞包从防热层沟槽中拉出，两者之间具有复杂的非线性作用关系，该过程可能会严重影响舱体和主伞包的运动特性，进而影响后续工作程序。因此，对返回舱拉出主伞包过程中防热结构在不同拉脱力下的工况进行了对比，本文分析了在5 kN和10 kN这2种不同防热层拉脱力下的工况。

图 17和18分别为在不同拉脱力下飞船返回舱拉出主伞包过程中主吊带和伞包拖带载荷随时间变化曲线。图 19和20分别为不同拉脱力下减速伞及返回舱速度和过载随时间的变化曲线。由计算可知，各吊带载荷及过载出现峰值的时刻相近。拉脱力为5 kN和10 kN对主吊带和伞包拖带载荷变化的影响较小。因此，拉脱力的变化对主伞包拉出过程的影响并不显著。

4) 在不同伞舱盖质量条件下计算。

本文对比了拉出主伞包过程中采用不同伞舱盖质量对整舱受力、速度和过载的影响。伞舱盖质量分别采用40和50 kg，系统初始速度采用60 m/s，防热层拉脱力工况采用5 kN。图 21和22为在不同伞舱盖质量下拉出主伞包过程中主吊带和伞包拖带载荷随时间变化曲线，图 23和24分别为在不同伞舱盖质量下减速伞及返回舱速度和过载随时间的变化曲线。由分析可知，伞舱盖质量的改变对吊带受力、减速伞与伞舱的运动速度、伞舱运动姿态等的影响较不明显，而对伞包拖带载荷曲线峰值影响较大。主要原因是伞舱盖质量的变化使伞包拖带被拉出过程中姿态有所变化，进而使各伞包拖带的受力分布差异明显，导致总载荷发生变化，这与不同拉脱力工况的对比结果相似。

在实际空投实验中，伞舱盖质量为50 kg，通过在返回舱模型质心位置安装加速度传感器对过载进行测量，并对数据进行50 Hz低通滤波处理，与实测舱盖过载时间变化曲线进行对比，峰值过载误差约为19.4%，整体变化趋势的吻合度较好。因此，采用动力学方法对伞包拉出过程进行数值模拟的方法可行，能够有效预测力学和运动参数变化。

4 结论

本文针对群伞系统减速伞分离拉出主伞包环节涉及的非线性、大变形、复杂柔性系统动力学仿真建模困难、细节载荷计算精度差的问题，提出了一种有限元分析与舱伞系统动力学计算相结合的仿真分析技术，开发了由降落伞气动载荷向动力学模型实时传递的专用算法，采用有限元方法，建立真实降落伞及绳带系统计算模型，对主伞包拉出过程中各部件的相对运动及载荷情况进行精确的量化计算，获得了不同连接部位的载荷分布，为系统方案设计和评估提供支撑。研究结论如下：

1) 基于有限元分析技术能够精准描述主伞包拉出过程各降落伞部组件间的耦合关系，简化伞衣气动力计算，引入动态迭代气动载荷子程序可降低计算量，为返回舱动力学分析提供了高效的分析手段。

2) 随着拉出主伞包过程初始速度增加，主伞包完全拉出所需时间缩短，各吊带载荷呈现增加趋势，整舱过载也显现出明显增加的趋势，因此需要严格控制初始拉出边界条件。

3) 防热层拉脱力由5 kN增加至10 kN的变化对主伞包拉出过程有一定影响，但总体而言，这种影响并不显著，可忽略。

4) 主伞舱盖质量的增大对伞包拖带总载荷曲线峰值影响显著，伞包拖带总载荷存在明显增大的趋势。因此在方案设计中要严格控制主伞舱盖的重量，从而降低伞包拖带的载荷。