朱尊红, 戈新生

(北京信息科技大学 机电工程学院，北京 100192)

自1958年“探索者1号”发射成功以来，研究人员就不断探讨航天器的功能，他们成功地将复杂的航天器引向更遥远的星空[1]。然而，未来大型空间结构的物理复杂性是需要迫切考虑的，且建模问题就成了重中之重[2]。以前航天器的功能比较简单，挠性附件所承担的任务较少，但随着航天领域的发展，现有航天器的功能性大大加强，规模越来越大，挠性附件安装的数量越来越多，就会使挠性附件在整个航天器之中所占的比重越来越大，与此同时，挠性附件的挠性问题就成为一个值得关注的问题[3]。由于外部扰动引起挠性附件振动，会影响航天器姿态稳定和姿态机动(刚柔耦合特性)，这种新特点对传统的航天器姿态运动柔性附件振动耦合动力学提出了严峻挑战。

挠性航天器中心体的运动，会引起大型挠性结构的变形，它们之间会相互耦合[4]，即航天器挠性附件的振动可能造成航天器运动失稳。现有航天器包括刚性中心体和挠性附件，对于此类航天器的振动分析，我们一般采用两种模态：约束模态和非约束模态。约束模态：刚性中心体被固定，而附件振动，也可称为悬臂模态或附件模态。非约束模态：不是简单的附件振动，而是整个航天器都在振动，并且航天器各个部分产生相互作用，也可称为载体模态。在非约束模态中，把航天器的运动分为两部分：① 航天器的整体刚体运动，包括平动和转动，此时把整个航天器看做刚体；② 航天器刚体和挠性附件的小尺度振动，此时其振动是同步的，且相互耦合。Hughes等[5-6]针对多体挠性航天器，创建了一般动力学方程，并讨论了模态恒等式，即惯性完备性准则，描述了耦合系数与频率之间的关系。徐小胜等[7]阐述了包括柔性核、柔性耦合系数(包括航天器的平动和转动)、两种模态的概念，并利用模态恒等式进行惯性完备性降阶。以上研究只是阐述两种模态的模型，并未对特征值问题进行求解。吕旺等[8]探讨了一种挠性航天器模态计算方法，通过先求导约束模态下的系统方程，然后再求解航天器非约束模态。Barbieri[9]针对柔性回转杆，导出两种模态下的运动方程及边界条件，并求出相关的频率方程和振型方程，比较了约束模态和非约束模态的差异。袁秋帆等[10]采用一种全局模态，推导了全局模态下的动力学方程，并推导了非约束模态的特征值问题求解方法。

以上在求解特征值问题时，把挠性附件当做二维梁处理，并未按照实际情况当做帆板处理。李东旭[11]将挠性航天器等效为“航天器刚性主体+悬臂板模型”，进行了航天器-板式挠性结构系统刚柔耦合动力学分析。Hablani[12-13]针对大挠性结构板在约束和非约束模态下建立模型，并推导了动力学方程。程亮亮等[14]利用ANSYS软件和试验的方法对悬臂板进行模态分析，对相关特征值问题进行求解。杜圆等[15]提出一种基于改进傅里叶级数的方法，对任意边界条件下的矩形板进行振动分析。Leissat[16]详细地描述了各种特征方程解的存在性，给出全面、准确的矩形板自由振动分析结果。Bhat等[17-19]利用瑞利-里兹方法中的特征多项式或多项式计算矩形板及椭圆板的振动固有频率。曹志远[20]提出一种适合于求解各种边界条件矩形板的梁函数组合法，采用双向梁函数组合级数逼近方法对一般性边界条件的矩形板进行固有振动分析。

本文讨论中心刚体-单翼大挠性结构的航天器模型，这里的挠性结构为板式模型，即中心刚体带单侧矩形薄板简化模型。着重讨论板式挠性航天器在非约束模态下的动力学建模及特性，首先利用哈密顿原理推导单侧帆板航天器的运动方程，方程表述了中心刚体和挠性附件之间的耦合，然后分别在两种模态下进行了模态离散化，分别研究相关的特征值问题，比较两种模态下频率和振型，求解单翼板式挠性航天器的特征值问题，比较了约束模态与非约束模态之间的差异。

1 动力学建模及方程的建立

将挠性航天器等效为“航天器刚性主体+单翼挠性帆板模型”，不考虑板间连接铰链的影响，且航天器主体考虑为不变形的刚体，模型如图1所示。

图1 单翼帆板航天器

为更好地描述航天器中心刚体与挠性附件的变形及刚柔耦合作用，从X-Z和Y-Z两个视图描述变形情况。挠性航天器沿Z轴的平动z0、绕Y轴的转动θy和绕X轴的转动θx，如图2、图3所示。图2中中心刚体上方的小长方形是由于挠性航天器绕X轴的转动θx才能看到，而图3中中心刚体下方的的小长方形是由于挠性航天器绕Y轴的转动θy才能看到。

图2 挠性航天器X-Z视图

图3 挠性航天器Y-Z视图

设挠性航天器由主刚体R和挠性帆板ε组成，刚体平台是边长为a的立方体。其中OXYZ为惯性坐标系，O点位于变形前系统的质心位置，X轴水平向右，Y轴垂直向里，Z轴垂直向上，Og为变形后系统质心，Or为中心刚体质心。定义坐标系O0X0Y0Z0为浮动坐标系，O0点位于中心刚体与挠性帆板连接处的中点，X0轴沿板伸展方向，Y0轴垂直向里。由于中心刚体和挠性板在沿X轴和Y轴方向上以及绕Z转动方向几乎不发生耦合，所以只考虑沿Z轴的运动位移z0和绕X轴的转动角θx以及绕Y轴的转动角θy。

假定航天器姿态角θx,θy很小，得到

sinθx=θx,sinθy=θy,cosθx=1,cosθy=1

(1)

则浮动坐标系O0X0Y0Z0转换到惯性坐标系OXYZ的矩阵C0p为

(2)

P点位移向量Rp

(3)

式中，u为P点沿Z轴的位移，帆板变形后P点的位置向量R0为

(4)

其中初始位置向量

(5)

变形后中心刚体相对于惯性坐标系的位置向量为

(6)

则点P在惯性坐标系中的向量R为

R=r+C0pR0=

(7)

刚体平台是边长为a的立方体，绕X轴转动惯量为Jrx，绕Y轴转动惯量为Jry,航天器系统动能为

(8)

式中，μ为太阳能帆板的面密度。

考虑太阳能帆板，则航天器系统势能为

(9)

通过哈密顿原理

(10)

得到航天器系统动力学方程

(11)

(12)

(13)

(14)

2 模态离散化

此过程中，采用两种模态即约束模态和非约束模态进行模态离散化。

2.1 约束模态

约束模态下相当于把挠性帆板当做悬臂板处理，此时z0=0,θx=0,θy=0,式(14)为

(15)

对于此模态的特征值问题求解，本文运用康特洛维奇法。

设u(x,y,t)=U(x,y)sin(ωt+φ)，代入式(15)，得

∇4U=α4U

(16)

其中

(17)

边界条件

x=0:U=0,Ux=0

(18)

x=l:Uxx+νUyy=0,Uxxx+(2-ν)Uxyy=0

(19)

y=0:Uxx+νUyy=0,Uxxx+(2-ν)Uxyy=0

(20)

y=b:Uxx+νUyy=0,Uxxx+(2-ν)Uxyy=0

(21)

得到相应的振型变分方程

(22)

设振型函数为在X方向上的未知函数X(x)乘以在Y方向上满足边界条件的已知函数Yn(y)

U(x,y)=X(x)Yn(y)

(23)

Yn(y)一般可取自由-自由单向板第n阶振型函数(梁函数)，即

Y1=1

(24)

(25)

Y3=(coshαny+cosαny)-an(sinhαny+

sinαny),n>2

(26)

表1 自由-自由边界单向板频率与振型系数

而未知函数X(x)将通过变分方程满足，将式(23)代入式(22)，然后通过变分运算，得到常微分方程

(27)

及边界条件

(28)

(29)

其中

(30)

常微分方程式(27)满足边界条件式(29)的一般解为

(31)

其中

(32)

将式(31)代入边界条件式(28)，可得关于系数A，B的二阶线代方程组。从而可得频率方程

(33)

上述两式中

(34)

将式(32)及式(17)代入式(33)，对应每个n值，即可求得一系列频率值ωmn(m=1,2,3，…)。

2.2 非约束模态

非约束模态下，此时涉及中心刚体运动，则航天器上任意一点P的振动位移为

w(x,y,t)=z0e(t)+yθxe(t)-xθye(t)+

(35)

式中：b为挠性附件的宽度；x∈[-a+r0,r0+l]，引入模态坐标，将z0e(t),θxe(t),θye(t),u(x,y,t)表示为

(36)

(37)

(38)

(39)

则航天器上任意一点P的振动位移w(x,y,t)为

(40)

非约束模态表示航天器在无外力作用下的自由振动，则非约束模态振型为

wm(x,y)=z0m+yθxm-xθym+

(41)

则采用非约束模态时，整个航天器的运动为

(42)

式中，z0r,θxr,θyr为将整个挠性航天器看为刚体时，航天器的平动和绕X轴、Y轴的转动位移。当式(15)中约束模态下的位移换作非约束模态下的振动位移时，得到非约束模态下特征值问题

(43)

3 非约束模态频率计算

本文选用瑞利里兹法，用特征正交多项式簇作为基振型，设模态基振型X方向为一端固定一端自由，Y方向为两端自由，挠性板长度为1。

首先进行X方向的正交多项式函数，设多项式簇第一项基函数表示为

X1(x)=a0+a1x+a2x2+a3x3+a4x4

(44)

边界条件为

(45)

将式(45)代入式(44)，得

X1(x)=a4(6x2-4x3+x4)

(46)

式中，a4可取任意值，对第一阶多项式函数进行模态归一化，设

(47)

则后续正交基可式(48)得出

(48)

每一阶正交基模态归一化，即

(49)

当挠性板长度为l时，正交基为

(50)

然后进行Y方向的正交多项式函数，设多项式簇第一项基函数表示为Y1=1，后续正交基由式(51)求出

(51)

每一阶正交基模态归一化，即

(52)

当挠性板宽度为b时，正交基为

(53)

选取第m阶非约束模态

w(x,y,t)=[z0m+yθxm-xθym+um(x,y)]ηm(t)

(54)

式中，ηm(t)=sin(ωmt)，根据式(8)和式(9)，得到航天器第m阶系统动能Km和航天器第m阶系统势能Pm

2z0 m[um+yθxm-(r0+x)θym]-2(r0+

(55)

(56)

(57)

(58)

4 仿真分析

本文采用中心刚体加挠性帆板结构进行仿真分析，中心刚体为边长为2 m的立方体，中心刚体密度为300 kg/m3。挠性板长度l分别为4 m,8 m,12 m，宽度b为2 m，厚度h为0.02 m，挠性板材料密度为2.8×103kg/m3，弹性模量E为6.89×1010Pa，泊松比ν为0.33。

4.1 频率求解

在此过程中取X、Y方向正交基个数分别为4，因为只包括一侧帆板，所以并不会有正对称和反对称变形，求得前5阶非约束模态频率，约束模态亦考虑前5阶模态。注意非约束模态下刚体模态频率均为0，所以不考虑。然后通过有限元软件计算非约束模态频率，与解析法计算求得的非约束模态比较，并与约束模态作比较分析，其中利用有限元软件计算过程中节点选择2 604，单元为433。C为约束模态的数值解，CFEM为约束模态的有限元软件解,UC为非约束模态的数值解，UCFEM为非约束模态的有限元软件解。频率对比如表2、表3、表4所示。

表2 频率对比(l=4 m)

表3 频率对比(l=8 m)

本文模型为中心刚体加大挠性单翼帆板，该模型的有限元解作为参考标准解，约束模态与非约束模态数值解可通过与UCFEM有限元解进行对比。由表2～表4可知，瑞利里兹法求得的非约束模态频率与有限元软件求得的非约束频率接近，而约束模态频率一般小于非约束模态频率。随着太阳能帆板长度的增加，即刚柔惯性比、质量比的减小，挠性附件在整个航天器所占的比例越来越大，非约束模态频率相比于非约束有限元之间差别较小，而约束模态频率相对于非约束有限元频率会出现偏差，因此，对于小中心刚体与大挠性结构的挠性航天器构型非约束模态动力学模型要更加精确。

表4 频率对比(l=12 m)

与双翼帆板不同的是，由于只有一侧帆板，本文中并不存在正对称和反对称变形，所以非约束模态频率一般均比约束模态频率高，这说明在任何模态下帆板振动与刚性中心体运动都产生耦合，并不会出现双翼帆板中一阶扭转反对称下约束模态与非约束模态频率几乎相同的情况。

4.2 振型求解

图4 一阶外弯曲(UCFEM)

图5 二阶外弯曲(UCFEM)

图6 一阶扭转(UCFEM)

图7 三阶外弯曲(UCFEM)

图8 二阶扭转(UCFEM)

图9 一阶外弯曲(UC)

图10 二阶外弯曲(UC)

图11 一阶扭转(UC)

图12 三阶外弯曲(UC)

图13 二阶扭转(UC)

图中振型为挠性附件与刚性中心体运动的耦合，其中一阶外弯曲、二阶外弯曲体现了在Z方向上平移运动的耦合和绕Y轴的姿态运动耦合，一阶扭转则体现了绕X轴的姿态运动耦合，值得注意的是，在以往双翼帆板振型求解中，耦合只会发生在单一情况中，即每阶振型所体现的耦合作用只能是Z方向上平移运动的耦合或者绕Y轴的姿态运动耦合，又或者绕X轴的姿态运动耦合。但是单翼帆板却不同，单翼帆板会同时出现沿Z方向上平移运动的耦合和绕Y轴的姿态运动耦合，这是由于单翼帆板并不存在正对称或者反对称，所以某些平移运动或者姿态运动无法抵消。

通过比较图中表示的利用有限元法和瑞利里兹法求得的非约束模态振型可知，瑞利里兹法求得的非约束模态振型与有限元软件求得的非约束振型接近，当中心刚体在整体中所占的比重比较大时，约束模态振型能够准确描述实际情况，但是随着挠性附件在整体中占的比重越来越大时，非约束模态与约束模态之间差异就会比较明显，与约束模态相比，非约束模态更贴近实际情况。

5 结 论

本文针对中心刚体-单翼大挠性结构的航天器，提出了航天器太阳能帆板的板式简化模型，建立了一种非约束模态动力学模型，探讨了板式挠性航天器非约束模态的动力学建模及动态特性，分别在约束模态和非约束模态下研究相关的特征值问题，对两种模态的频率和模态振型进行比较，求解了单翼板式挠性航天器的特征值问题，比较了约束模态与非约束模态之间的差异。结果说明在单翼帆板情况下，当中心刚体在整体中所占的比重比较大时，约束模态能够准确描述实际情况，但是随着挠性附件在整体中占的比重越来越大时，非约束模态与约束模态之间差异比较明显，与约束模态相比，非约束模态更贴近实际情况。