哑铃型航天器刚-柔耦合动力学建模与仿真分析*

王博洋 刘铸永† 郑鹏飞

（1.上海交通大学船舶海洋与建筑工程学院，上海 200240）（2.上海宇航系统工程研究所，上海 201109）

引言

随着科技发展，我们的目光投向了更深更广阔的太空.由于深空探测要求航天器长时间在轨工作，这给能源技术带来了很大的挑战［1］.核能源具有高能量密度、使用寿命长等特点，是解决深空探测能源问题的一个有效途径.然而核动力航天器的一个重要问题就是核辐射.在设计核动力航天器时，必须将核反应堆和航天器主体分离开，中间由桁架连接，这就会导致航天器整体构型出现质量集中在两端，中间由桁架连接的“哑铃型”.目前已有的哑铃型航天器概念设计航天器有Project Gerkule的Space Tug飞船，Project Promethe的JIMO和TPM飞船等［2］.这种类型的航天器主要由主星体、连接桁架和末端反应堆组成.

传统航天器一般由主星体和太阳电池阵构成［3］，有时根据任务需要还会附带网状天线等［4］.这些航天器主星体质量占主导地位，在建模时常采用“中心刚体+柔性附件”模型［5］.由于哑铃型航天器的构型与质量分布同传统航天器有很大差异，故传统的建模方法将不再适用于这种类型的航天器.目前，国内外对于这种新构型的哑铃型航天器轨道—姿态动力学和控制研究相对较少.文献［6，7］将系统两端简化为质点，中间由为无质量刚性杆连接.文献［8，9］将两端物体简化为质点，连接桁架简化为无质量的弹簧.文献［10］将两端物体简化为刚体，中间由无质量的弹簧连接.文献［11-13］将末端反应堆简化为质点，连接桁架简化为一根柔性梁，通过有限元方法将航天器整体离散，利用牛顿-欧拉方程和结构动力学方程建立了航天器的动力学模型.文献［14］两端物体简化为刚体，将连接桁架简化为一根梁，利用拉格朗日方程建立了平面问题的转动和振动的耦合姿态动力学方程.

图1 现有哑铃型航天器设计［2］Fig.1 Dumbbell shape spacecraft

但是现有的动力学模型也存在一定问题，例如：文献［6-9］中采用的刚体模型对于一些定性分析是可行的，但是显然无法反映柔性变形；文献［10］中采用的刚体+弹簧模型相比刚体模型，更贴近现实情况，但在模型准确性方面来说仍有不足；文献［11-13］中使用的整体离散模型，当航天器构型上发生改变时，便需要重新建立有限元模型计算整体模态，这给设计、优化分析等方面带来很大不便.此外，文献［11-13］中将末端反应堆简化为质点，没有考虑其转动惯量，这种简化的合理性尚值得商榷.

本文针对哑铃型航天器构型，基于多体系统动力学单向递推组集方法，建立了该系统准确高效的刚-柔耦合动力学模型，并与现有的动力学模型进行对比，探索哑铃型航天器动力学模型的有效性.

1 哑铃型航天器刚柔耦合动力学建模

由于哑铃型航天器中，通常两端物体的结构刚度较大，在建模时可以简化为刚体；中间的大跨度连接桁架具有较大的柔性，在建模时必须按照柔性体处理.虽然在运行过程中，桁架由于柔性较大会产生一定的变形，但是与桁架的整体长度相比，变形仍然可以算是小量，故此时仍然满足小变形假设.此时，系统为两端刚体外加柔性连接桁架的刚柔耦合系统.设主星体为刚体B1（Main body），连接桁架为柔性体B2（Truss structure），末端反应堆为刚体B3（Power unit）.

1.1 运动学递推关系

哑铃型航天器结构如图2所示.e1和e3为物体B1和B3的连体基，e2为物体B2的浮动基.P2和P3为桁架两端面中心点.对于柔性桁架而言，采用集中质量有限元离散，其上任意一节点k的绝对位置可以表示为［15］

图2 哑铃型航天器示意图Fig.2 Schematic diagram of dumbbell spacecraft

式中，r2为B2的浮动基基点的位置坐标阵，表示节点k变形后在惯性基的位置坐标阵，表示节点k未变形时在浮动基下的位置坐标阵，表示节点k的平移模态矩阵在浮动基下的位置坐标阵，计算时保留到s阶模态，a2表示柔性桁架的s阶模态坐标，A02表示连体基e2关于惯性基e0的方向余弦阵，右上标“’”为该坐标阵是在浮动基上的表示.

分别对式（1）求时间t一次和二次导数，并加以整理可得到物体B2上任意一节点k与其绝对速度和加速度的矩阵表达式：

式中，ω2为浮动基e2的角速度在惯性基下的坐标阵，上标“～”表示该坐标阵的反对称坐标方阵.类似地，定义刚体B1和B3的绝对速度矩阵：

式中ω1和ω3分别为连体基e1和e3的角速度在惯性基下的坐标阵.

此外，在小变形假设下，由于变形引起的物体B1-B2，B2-B3之间的方向余弦阵可以近似表达为

其中，I3为3阶单位阵，分别表示点P2点和P3点的转动模态矩阵在惯性基下的位置坐标阵.基于多体系统动力学单向递推组集方法［16］，可以得到系统中各物体绝对速度矩阵以及绝对加速度矩阵之间的递推关系式

式中，03为3阶零矩阵，0s×3以及0s×1为其下标所示阶数的零矩阵.

航天器在运行过程中，系统质心时刻都在变化.如图2所示，定义连体基e与物体B1固结，其原点C为桁架未变形时系统的质心.当桁架变形时，系统质心在C点附近运动.定义r为原点C在惯性基的坐标阵.为了表示C点在空间中的运动轨迹，利用以下转换关系：

1.2 动力学方程

基于速度变分原理［17］，可得到整个系统的动力学方程：

其中mi和Ji为物体Bi的质量和转动惯量张量阵.将运动学递推关系式（10）、式（11）和式（17）代入式（21）得

由于采用了单向递推组集的方法，非独立的变分项在递推过程中由运动学关系消去，此时变分独立，系统的动力学方程为

其中广义质量阵的各项为

式中，Mtt为系统质量阵；Mtr为系统平移-转动耦合矩阵，在系统未变型时为零矩阵；Mta为系统平移-振动耦合矩阵；Mrr为系统相对于O点的转动惯量阵；Mra为转动-振动耦合矩阵；Maa为振动质量阵.

式中，mi为物体Bi的质量，Ji为物体Bi的转动惯量张量在惯性基上的坐标阵，Ca和Ka分别为柔性桁架的模态阻尼阵和模态刚度阵.此外还有一些为表示方便而定义的矩阵：

该动力学模型保留了所有变形产生的高次耦合项，为全耦合的“刚体-桁架-刚体”模型.此外，按照文献［11-13］，将末端物体简化为质点，即J'3=0时，该模型称为“刚体-桁架-质点”模型.

2 数值仿真与分析

利用推导的动力学方程，编写动力学计算程序，设计工况对其进行动力学仿真.由于本文着重讨论的是具有哑铃状构型特点的航天器的动力学建模过程，并不针对某一特定型号进行建模，故在仿真计算时可用相对理想的参数模型，系统具体参数如表1所示.其中转动惯量矩阵分别为物体相对各自质心连体坐标系（B2为相对质心浮动坐标系）.利用有限元软件Abaqus对30m长的柔性桁架进行有限元离散，单元类型采用B33欧拉梁单元，节点数为1292，单元数为1612.

表1 物体质量与惯量Table 1 Mass and inertia of bodies

取自由—自由边界条件，动力学仿真时保留除刚体模态外的前11阶模态.各阶模态频率如表2所示，模态振型如图3所示.当中包含8阶弯曲模态，如第1-6阶和第8-9阶；两阶扭转模态，如第7阶和第10阶；一阶纵向振动模态，如第11阶.

图3 桁架构型与各阶模态振型图Fig.3 Modal shape of truss structure

表2 桁架物理参数Table 2 Parameters of truss

由于本文考虑的是空间问题，故桁架的弯曲模态中会出现频率相同、但振型正交的两组模态.

2.1 动力学仿真

分别采用本文推导的“刚体-桁架-刚体”模型、文献［11-13］提出的“刚体-桁架-质点”模型以及文献［10］提出的“刚体-弹簧-刚体”模型三种模型进行数值仿真对比.仿真时，“刚体-桁架-刚体”模型和“刚体-桁架-质点”模型动力学仿真时保留前11阶模态.仿真时采用滚转、偏航、俯仰角描述姿态.

在实际运行过程中，执行机构集中在主星体B1处，故仿真时主动力只施加在主星体B1上.本文采用的仿真工况为航天器三个方向均受外力矩作用，模拟航天器作姿态机动时的情况，具体数值大小见表3.由于航天器构型的原因，其滚转轴方向的转动惯量要远远小于偏航和俯仰方向的转动惯量，故该方向施加的外力也相应较小.

表3 桁架模态频率Table 3 Modal frequencies of truss

表4 受力情况Table 4 Acting forces

动力学仿真时长为100s，仿真结果如图4-6所示，对比了三种模型的姿态坐标、桁架最大变形和末端物体相对姿态的时间历程.

图4 主星体姿态角俯仰角对比Fig.4 Comparison of attitude angles of main body

图5 桁架最大变形量对比Fig.5 Comparison of truss maximum deformation

图6 末端物体相对姿态Fig.6 Relative attitude angles of power unit

2.2 动力学模型对比分析

从图4中可以看出，由于连接桁架柔性较大，故航天器在受到外力矩时，会产生较大幅度的低频振动，从图5（c）中可以看出，桁架末端的最大变形幅度可达1.2m.这种振动对于航天器的正常运行来说是十分不利的.故在进行姿态机动时，应采取一定措施，控制桁架变形，减小振动幅度.

从“刚体-桁架-刚体”模型和“刚体-弹簧-刚体”模型的对比中可以看出，“刚体-桁架-刚体”动力学模型能够很好地反映出航天器由于连接桁架的弹性振动与航天器姿态之间的耦合.从图5可以看出，由于桁架纵向刚度远远大于横向刚度，桁架的变形主要为弯曲变形，而不是纵向的伸缩变形，导致“刚体-弹簧-刚体”模型的响应与刚体类似，无法体现出弹性振动与航天器姿态之间的耦合.由此说明对于哑铃型的航天器，文献［9］中的“刚体-弹簧-刚体”模型无法准确描述系统的动力学行为.

从图4中的“刚体-桁架-刚体”模型和“刚体-桁架-质点”模型的对比中可以看出，忽略末端物体的转动惯量后，系统的动力学响应会产生一定的差异.从图4（a）可以看出，两模型在主星体滚转角上的差异明显，且趋于发散.这是由于哑铃型航天器的构型特点造成的，其滚转轴上的转动惯量会远小于其他两个轴，这导致末端物体滚转轴方向上的转动惯量占整体的比例会远大于其他两个方向.此时如果忽略了末端物体的转动惯量，会造成较大的误差.两模型差异不仅体现在大范围刚体运动上，还体现在其扭转振动上.从图6（a）中相对姿态的对比还可以看出，忽略了末端物体的转动惯量后，“刚体-桁架-质点”模型的扭转变形会远小于“刚体-桁架-刚体”模型.从图4（b）、图4（c）、图6（b）、图6（c）可以看出，虽然两模型在偏航角和俯仰角上的差异小于滚转角，但是存在扩大的趋势.这是因为空间物体的三方向转动是耦合的，随着时间的增加，两模型在偏航角和俯仰角上的差异会因为滚转角差异的增加而增加.

“刚体-桁架-刚体”模型和“刚体-桁架-质点”模型差异不仅体现在响应幅值上，还体现在响应频率上.在图4（b）—图4（c）和图5（b）—图5（c）中，两模型响应幅度差异不大，但两模型响应频率有较明显的差异.从定性分析的角度，将整体简化为一个二阶线性系统［18］，由于忽略了末端物体的转动惯量，系统的广义质量减小，导致了整体频率的升高，其频率相对误差在本构型参数下为2%左右.

此外，如表5所示，表中成对出现的频率为同一阶弯曲模态在相互垂直的两个方向的频率.由于忽略了末端的转动惯量，“刚体-桁架-质点”模型的扭转频率远高于“刚体-桁架-刚体”模型，系统的振动性质发生了改变.系统在扭转方向的频率误差达519%，这个现象从图6（a）中可以清晰地看出，如果在建模时忽略了末端的转动惯量，则系统的扭转方向刚度会大大增加，与现实情况不符.

表5 两模型频率对比Table 5 Frequency comparison of two models

综上所述，由于将末端物体简化为质点，无论从时域还是从频域的角度，都会带来较大的误差，因此在动力学建模时不能将末端物体简化为质点，必须考虑其转动惯量.

3 结论

为了满足深空探测的需求，需要构建一种全新的哑铃型航天器.由于构型和质量分布的差异，传统“中心刚体+柔性附件”动力学模型将不再适用.为了满足哑铃型航天器总体设计的需要，必须针对该构型建立一种高效准确的动力学模型.

本文针对哑铃型构型，在小变形假设下，采用多体系统动力学单向递推组集的方法，基于速度变分原理建立该“刚体-桁架-刚体”系统的刚-柔耦合动力学模型.该动力学模型是保留了所有弹性变形和大范围刚体运动的耦合项.对比了本文所建立“刚体-桁架-刚体”模型和已有的“刚体-桁架-质点”模型及“刚体-弹簧-刚体”模型的动力学响应.研究结果表明：现有的“刚体-弹簧-刚体”不能很好地反映航天器姿态与桁架变形之间的耦合效应.“刚体-桁架-刚体”模型与“刚体-桁架-质点”模型在滚转角响应上有较大差异，且随着时间的推移，会导致偏航角和俯仰角也产生较大差异；将末端物体简化为质点还会造成系统整体频率的升高，说明哑铃型航天器建模时必须考虑末端物体的转动惯量.“刚体-桁架-刚体”很好地反映了哑铃型航天器的刚-柔耦合动力学特性.此外，由于哑铃型航天器质量分布的特点，以及整体尺寸和桁架柔性都较大，系统在运动过程中可能会产生较大的弹性变形，今后应采取适当措施抑制桁架振动.本文研究将为哑铃型航天器的总体设计，特别是控制规律的设计，提供重要的技术支持.