空间薄膜天线变形/姿态耦合动力学与控制

Deformation/Attitude Coupling Dynamics and Control of Space Membrane Antennas

摘要

针对空间薄膜天线自旋过程中柔性体多尺度耦合动力学建模与控制难题，提出了一种薄膜天线混合建模与姿态解耦控制策略。采用绝对节点坐标法的薄板单元与圆截面梁单元建立了薄膜-空心撑杆混合单元动力学模型。引入薄膜非线性本构描述薄膜天线作动过程中张紧-松弛状态的变化。利用绝对节点坐标法可以同时描述柔性体变形与刚体转动的特性，对选定节点的位置向量梯度矩阵进行极化分解，解算出薄膜天线姿态信息并在此基础上施加比例-微分力矩控制。仿真结果表明，提出的动力学建模与控制方法可以有效地实现薄膜天线姿态跟踪误差渐进收敛。

Abstract

Aiming to address the challenging issues of multiscale coupled dynamics modeling and control for flexible bodies during the spin deployment of space membrane antennas， a hybrid modeling approach and attitude decoupling control strategy was proposed herein. A hybrid element dynamic model integrating membrane and hollow support rods was developed using thin plate elements and circular cross-section beam elements based on ANCF. Nonlinear constitutive relations were introduced to characterize the tensioned-relaxed state transitions of space membrane antennas during deployment. By leveraging ANCF's capability to describe flexible body deformation and rigid body rotation simultaneously， polar decomposition was performed on the gradient matrices of selected nodal position vectors to extract attitude information. A proportional-derivative （PD） torque control scheme was subsequently implemented based on the derived posture parameters. Simulation results demonstrate that the proposed dynamics modeling and control methodologies may achieve asymptotic convergence of attitude tracking errors effectively.

Graphical abstract

关键词

薄膜天线 / 变形/姿态耦合 / 绝对节点坐标法 / 柔性多体动力学

Key words

space membrane antenna / deformation/attitude coupling / absolute nodal coordinate formulation （ANCF） / dynamics of flexible multibody system

0 引言

近年来，可展开空间结构技术快速发展，尤其在大功率空间发电、天基遥感等重大应用需求牵引下，空间可展开结构正呈现出大型化、轻量化、高精度等发展趋势^［1］。薄膜结构作为一种新型大型可展开结构，具有轻质、易折展等明显优点，被广泛应用于合成孔径雷达、可展开天线、太阳帆等各种空间装备与结构中。然而，薄膜结构本身弯曲刚度很小，易产生折痕和褶皱问题，从而影响薄膜的形面精度，进一步影响其结构和功能特性^［2］。通过深入研究薄膜天线的材料、结构及性能优化，可以进一步提高其效率和可靠性，推动航空通信技术的创新与发展，为未来智能社会的构建提供强有力的技术支持。

薄膜天线结构近年来成为空间可展开领域的研究热点。在平面薄膜天线研究中，美国喷气式推进实验室（Jet Propulsion Lab， JPL）团队^［3-4］在2002年提出了具有矩形框架的“电影屏幕”式膜阵列天线的概念，而后于2008年开发的2.2 m×2.2 m天线^［5］虽实现了0.17 mm表面平整度，但没有考虑极薄的膜在张力下易发生屈曲的几何不稳定性与非线性屈曲模态问题。加拿大航天局（Canadian Space Agency， CSA）^［6-7］在2007年开发了三层膜合成孔径雷达（synthetic aperture radar， SAR）天线。2008年德国航空航天中心^［8-9］开发的基于可卷曲的碳纤维增强聚合物桁架的可展开结构虽解决了宏观展开问题，但未建立膜面与支撑结构的刚柔耦合动力学模型。国内，肖薇薇等^［10］通过ABAQUS壳单元分析褶皱特性时，受限于线性应变假设，难以捕捉大变形区域的非线性应变分布。胥小勇等^［11］设计了薄膜在收换卷时的张力控制系统。操安博等^［12］针对空间薄膜结构高精度设计中粘接缝二维模型缺失的问题，开展了空间薄膜结构粘接缝模型的建立与验证工作。荣吉利等^［13］采用隐式非线性求解器SAMCEF分析太阳翼，使用膜单元模拟翼面，对其展开过程进行动力学分析。陈金铎等^［14］描述了薄膜结构大变形、耦合光压力与热载荷，构建10 μm量级薄膜非线性动力学主控方程，并对其变形特征进行分析。上述针对薄膜的研究大都会面临柔性体大变形的几何问题以及材料的弯曲或拉伸过程中刚度矩阵的变化等非线性问题。

1996年SHABANA^［15］提出的绝对节点坐标法（absolute nodal coordinate formulation， ANCF）可统一描述柔性体的刚体运动与弹性变形，特别适用于薄膜天线的大变形、大转动及非线性动力学问题。国内学者应用ANCF方法对薄膜类结构进行了深入研究。赵将等^［16］对简化的“IKAROS”太阳帆薄膜采用ANCF描述的三角形薄板单元进行建模并展开动力学分析；LIU等^［17］提出了ANCF梯度缺失薄壳单元和ANCF膜单元；薛鹏聪等^［18］整合张力场和不同模量弹性理论，提出了一种考虑不同模量的ANCF薄膜单元计算方法，对大型帆膜结构的设计有理论指导意义。

航天器姿态控制一直是一个需要不断深入研究的课题。支敬德等^［19］对原有的滑模控制律进行改进，设计了一种模糊滑模控制律，采用双曲正切函数代替原本的符号函数，并引入模糊控制对切换增益进行实时调整，可以有效抑制抖振。尹一蓁等^［20］针对薄膜航天器姿态控制问题，考虑大气扰动、飞轮摩擦力矩、转动惯量变化等因素，提出了一种基于非线性变结构滑模控制器的姿态机动控制方法。殷春武^［21］对现有的线性和非线性控制器设计方法进行了详细的综述，综合分析了现有姿态控制器的特点和应用情况。上述控制方法中的动力学模型大部分基于传统的小变形理论，目前基于ANCF建模的大变形结构的控制方法研究较少。LIU等^［22］针对基于ANCF梁单元建模的柔性机械臂控制，提出了两步优化策略来规划机械手的运动。

区别于现有的只考虑单一柔性梁或板结构，本文基于ANCF建模方法，整合薄膜单元与三维圆截面梁单元，提出了薄膜-撑杆混合单元建模方法，以实现空间薄膜天线结构的动力学模型的构建。此外，区别于现有的含中心刚体的刚柔耦合航天器，以中心刚体的位姿作为被控对象，本文提出了不含中心刚体的基于极化分解的柔性体结构位姿求解方法。最后采用简单成熟的PD控制方法验证求解薄膜姿态的动力学位姿的正确性和控制的有效性。

1 运动学描述

1.1 薄板单元

全参数的绝对节点坐标有限单元都可视作实体，只是根据其不同方向描述变形能力要求采用不同阶次的插值多项式。如，MIKKOLA等^［23］提出的含有48自由度的厚板单元，就是在长、宽方向采用三次插值，而在厚度方向采用一次插值。而DUFVA等^［24］指出，当建模对象的厚度很小时，采用全参数板单元会带来截面闭锁和数值效率下降的问题，并提出了一种基于Kirchhoff-Love板壳理论的36自由度绝对节点坐标薄板单元，其每个节点只保留位置向量\mathit{r}和长、宽两个方向的梯度向量{\mathit{r}}_x\mathrm{、}{\mathit{r}}_y作为节点坐标，如图1所示。薄板单元上任意一点的全局位置矢量在全局坐标系下定义为

\mathit{r}=\mathit{S}(x,y)\mathit{e}（1）

式中：\mathit{r}为板上任意一点的全局位置矢量；(x,y)为单元物质坐标；\mathit{S}为形函数矩阵；\mathit{e}为节点坐标矢量。

可以看出，单元形函数已不再含有厚度方向物质坐标。单元采用中间面的4个角点作为单元节点。每个节点j的坐标{\mathit{e}}_j为

{\mathit{e}}_j=({\mathit{r}}_j^{\mathrm{T}}{,\frac{\partial {\mathit{r}}_j^{\mathrm{T}}}{\partial x},\frac{\partial {\mathit{r}}_j^{\mathrm{T}}}{\partial y})}^{\mathrm{T}}j=\mathrm{1,2},\mathrm{3,4}（2）

式中：{\mathit{r}}_j为节点j的全局位置向量；\partial {\mathit{r}}_j/\partial x\mathrm{、}\partial {\mathit{r}}_j/\partial y为节点j的位置向量梯度。

板单元的节点坐标向量如下：

\mathit{e}=({\mathit{e}}_{\mathrm{1}}^{\mathrm{T}},{\mathit{e}}_{\mathrm{2}}^{\mathrm{T}},{\mathit{e}}_{\mathrm{3}}^{\mathrm{T}},{\mathit{e}}_{\mathrm{4}}^{\mathrm{T}}{,)}^{\mathrm{T}}（3）

单元形函数矩阵\mathit{S}表示如下：

\mathit{S}=\left[S_{\mathrm{1}}\mathit{I}{S}_{\mathrm{2}}\mathit{I}\dots S_{\mathrm{11}}\mathit{I}{S}_{\mathrm{12}}\mathit{I}\right]（4）

式中：\mathit{I}为维数为3×3的单位矩阵；形函数的具体形式见文献［24］。

薄板单元的质量矩阵表达式为

\mathit{M}={\int }_V^{}\rho {\mathit{S}}^T\mathit{S}\mathrm{d}V（5）

其中，\mathit{M}为薄板单元的常值质量阵；\rho为单元密度；V为薄板单元体积。应用薄板变形模式假设，体积分可以转化为单元中面上的积分，即\mathrm{d}V=t_{\mathrm{1}}\mathrm{d}S，t₁为单元厚度，S为单元中面上面积。

1.2 圆截面梁单元

本节阐述的以绝对节点坐标法为基础的圆截面梁单元是在三维矩形截面梁单元^［25］的基础上改进得到的。其单元的位移场表示可参考三维矩形梁单元，用与物质坐标(x,y,z)有关的多项式来表示，在x方向是三次的，在y\mathrm{、}z方向是线性的。圆截面梁单元的位移场离散方式表示如下：

\mathit{r}=\left[\begin{array}{c}{r}_{\mathrm{1}}\\ r_{\mathrm{2}}\\ r_{\mathrm{3}}\end{array}\right]=\left[\begin{array}{l}{a}_{\mathrm{0}}+a_{\mathrm{1}}x+a_{\mathrm{2}}y+a_{\mathrm{3}}z+a_{\mathrm{4}}xy+\\ a_{\mathrm{5}}xz+a_{\mathrm{6}}{x}^{\mathrm{2}}+a_{\mathrm{7}}{x}^{\mathrm{3}}\\ b_{\mathrm{0}}+b_{\mathrm{1}}x+b_{\mathrm{2}}y+b_{\mathrm{3}}z+b_{\mathrm{4}}xy+\\ b_{\mathrm{5}}xz+b_{\mathrm{6}}{x}^{\mathrm{2}}+b_{\mathrm{7}}{x}^{\mathrm{3}}\\ c_{\mathrm{0}}+c_{\mathrm{1}}x+c_{\mathrm{2}}y+c_{\mathrm{3}}z+c_{\mathrm{4}}xy+\\ c_{\mathrm{5}}xz+c_{\mathrm{6}}{x}^{\mathrm{2}}+c_{\mathrm{7}}{x}^{\mathrm{3}}\end{array}\right]（6）

式中，a_{\mathrm{0}}\mathrm{、}{a}_{\mathrm{1}}\mathrm{、}\dots \mathrm{、}{c}_{\mathrm{6}}\mathrm{、}{c}_{\mathrm{7}}均为位移插值函数的多项式系数。

图2所示为本节所研究的二节点圆截面梁单元，点A、B是单元的两个节点，每个节点由12个广义坐标定义，每个单元由24个广义坐标定义，节点和单元的广义坐标用矢量形式表示如下：

\left.\begin{array}{l}{\mathit{e}}_i=({\mathit{r}}^{\mathrm{T}},{\mathit{r}}_{,x}^{\mathrm{T}},{\mathit{r}}_{,y}^{\mathrm{T}},{\mathit{r}}_{,z}^{\mathrm{T}}{)}_i^{\mathrm{T}}i=A,B\\ \mathit{e}=({\mathit{e}}_A^{\mathrm{T}},{\mathit{e}}_B^{\mathrm{T}}{)}^{\mathrm{T}}\end{array}\right\}（7）

其中，{\mathit{r}}_{,x},{\mathit{r}}_{,y},{\mathit{r}}_{,z}分别为对梁轴线方向和横截面的两个方向的梯度向量。单元内任意一点的位置可以用绝对节点法表示为

\mathit{r}=\mathit{S}(x,y,z)\mathit{e}（8）

单元的形函数矩阵\mathit{S}可由插值多项式和节点坐标推导：

\mathit{S}=\left[\begin{array}{llllllll}{S}_{\mathrm{1}}\mathit{I}& S_{\mathrm{2}}\mathit{I}& S_{\mathrm{3}}\mathit{I}& S_{\mathrm{4}}\mathit{I}& S_{\mathrm{5}}\mathit{I}& S_{\mathrm{6}}\mathit{I}& S_{\mathrm{7}}\mathit{I}& S_{\mathrm{8}}\mathit{I}\end{array}\right]（9）

式中：\mathit{I}为维度是3×3的单位矩阵。

本节所阐述的圆截面梁单元的质量阵如式（5）所示。圆截面梁单元与三维矩形截面梁单元的不同之处在于对其体积积分时横截面的积分方式构造为高斯-切比雪夫积分^［26］。圆截面梁单元的轴向积分独立于横截面积分，体积积分方式如下：

I={\int }_{m_{\mathrm{1}}}^{n_{\mathrm{1}}}{\int }_A^{}f(\xi ,\eta ,\zeta )\mathrm{d}A\mathrm{d}x\approx \sum _{i}w(i){\int }_A^{}f(\xi (i),\eta ,\zeta )\mathrm{d}A\cdot

l=\sum _{i}w(i){\int }_a^b{\int }_{c(\eta )}^{d(\eta )}h(\eta ,\xi )\mathrm{d}\eta \mathrm{d}\zeta（10）

式中：w(i)为高斯积分的权重；(\eta ,\xi ,\zeta )为量纲一物质坐标；l为圆截面梁的长度；m_{\mathrm{1}}\mathrm{、}{n}_{\mathrm{1}}为轴线方向积分上下限；a\mathrm{、}b\mathrm{、}c(\eta )\mathrm{、}d(\eta )分别为横截面的范围；A为横截面的面积。

对于空心圆管，只需将横截面面积的积分范围更改为

A=\{(\eta ,\zeta ),(R_{\mathrm{1}}/R_{\mathrm{2}}{)}^{\mathrm{2}}\le {\eta }^{\mathrm{2}}+{\zeta }^{\mathrm{2}}\le \mathrm{1}\}（11）

式中：R_{\mathrm{1}}为圆环的内径；R_{\mathrm{2}}为圆环的外径。

空心圆环的高斯-切比雪夫积分方式详见文献［26］。

2 薄膜本构

通常，薄板单元的应变能^［24］可分为两部分：

U=U^{\mathrm{m}\mathrm{i}\mathrm{d}}+U^{\mathrm{b}}=\frac{\mathrm{1}}{\mathrm{2}}{\int }_V^{}({\mathit{\epsilon }}^{\mathrm{T}}{\mathit{E}}_{\mathrm{\epsilon }}\mathit{\epsilon }+{\mathit{\kappa }}^{\mathrm{T}}{\mathit{E}}_{\kappa }\mathit{\kappa })\mathrm{d}V（12）

式中：U^{\mathrm{m}\mathrm{i}\mathrm{d}}为壳单元面内变形应变能；U^{\mathrm{b}}为面外弯曲变形应变能；\mathit{\epsilon }为描述壳单元面内变形的应变张量；{\mathit{E}}_{\epsilon }为与面内变形有关的弹性系数矩阵；\mathit{\kappa }为描述壳单元面外弯曲的应变张量；{\mathit{E}}_{\kappa }为与面外弯曲相关的弹性系数矩阵。

2006年，MIYAZAKI^［27］提出了各向同性膜的褶皱等状态判定准则及修正后的主应力和弹性矩阵。本节通过整合薄膜褶皱模型与ANCF薄板单元来对薄膜进行离散。由于薄膜抗弯刚度很小，因此弯曲刚度可以忽略不计。在各向同性材料的本构关系中，第二类皮奥拉-基尔霍夫应力张量对应的主应力\mathit{\tau }与格林-拉格朗日应变张量对应的主应变\mathit{\gamma }之间的关系可表示为

\mathit{\tau }=\mathit{\Gamma }\mathit{\gamma }（13）

\mathit{\tau }=\left[\begin{array}{c}{\sigma }_{\mathrm{1}}\\ {\sigma }_{\mathrm{2}}\end{array}\right]\mathit{\gamma }=\left[\begin{array}{c}{\gamma }_{\mathrm{1}}\\ {\gamma }_{\mathrm{2}}\end{array}\right]（14）

式中：\mathit{\Gamma }为弹性矩阵^［27］；{\sigma }_{\mathrm{1}}为最大主应力；{\sigma }_{\mathrm{2}}为最小主应力；{\gamma }_{\mathrm{1}}\mathrm{、}{\gamma }_{\mathrm{2}}为最大最小主应变。

薄膜单元中的格林-拉格朗日应变张量可以改造为ANCF薄板中的应变张量\mathit{\epsilon }：

\left[\begin{array}{c}{\gamma }_{\mathrm{1}}\\ {\gamma }_{\mathrm{2}}\end{array}\right]=\frac{{\epsilon }_{\mathrm{11}}+{\epsilon }_{\mathrm{22}}}{\mathrm{2}}\pm \sqrt[]{{(\frac{{\epsilon }_{\mathrm{11}}-{\epsilon }_{\mathrm{22}}}{\mathrm{2}})}^{\mathrm{2}}+{(\frac{\mathrm{2}{\epsilon }_{\mathrm{12}}}{\mathrm{2}})}^{\mathrm{2}}}（15）

\left.\begin{array}{l}{\epsilon }_{\mathrm{11}}=\frac{\mathrm{1}}{\mathrm{2}}({{\mathit{r}}_x}^{\mathrm{T}}{\mathit{r}}_x-\mathrm{1})=\frac{\mathrm{1}}{\mathrm{2}}[({\mathit{S}}_x{\mathit{e})}^{\mathrm{T}}({\mathit{S}}_x\mathit{e})-\mathrm{1}]\\ {\epsilon }_{\mathrm{22}}=\frac{\mathrm{1}}{\mathrm{2}}({{\mathit{r}}_y}^{\mathrm{T}}{\mathit{r}}_y-\mathrm{1})=\frac{\mathrm{1}}{\mathrm{2}}[({\mathit{S}}_y{\mathit{e})}^{\mathrm{T}}({\mathit{S}}_y\mathit{e})-\mathrm{1}]\\ {\epsilon }_{\mathrm{12}}=\frac{\mathrm{1}}{\mathrm{2}}{{\mathit{r}}_x}^{\mathrm{T}}{\mathit{r}}_y\end{array}\right\}（16）

修正后的主应力矢量和弹性矩阵如下：

\tilde{\mathit{\tau }}=\left[\begin{array}{c}{\tilde{\sigma }}_{\mathrm{1}}\\ {\tilde{\sigma }}_{\mathrm{2}}\end{array}\right]=\tilde{\mathit{\Gamma }}\mathit{\gamma }（17）

\tilde{\mathit{\Gamma }}=\frac{E}{\mathrm{1}-k_{\mathrm{1}}{k}_{\mathrm{2}}{\upsilon }^{\mathrm{2}}}\left[\begin{array}{cc}{k}_{\mathrm{1}}& k_{\mathrm{1}}{k}_{\mathrm{2}}\upsilon \\ k_{\mathrm{1}}{k}_{\mathrm{2}}\upsilon & k_{\mathrm{2}}\end{array}\right]（18）

式中：k_{\mathrm{1}}\mathrm{、}{k}_{\mathrm{2}}为松弛度和褶皱度；E为薄膜的弹性模量；\upsilon为泊松比。

定义两个参数f_{\mathrm{1}}\mathrm{、}{f}_{\mathrm{2}}来表示褶皱判据：

\left.\begin{array}{l}{f}_{\mathrm{1}}={\gamma }_{\mathrm{2}}+\upsilon {\gamma }_{\mathrm{1}}\\ f_{\mathrm{2}}={\gamma }_{\mathrm{1}}+u_{\mathrm{1}}\upsilon {\gamma }_{\mathrm{2}}\end{array}\right\}（19）

薄膜的张紧、褶皱、松弛三种应力状态由如下判定准则决定：

（20）

式中：参数u_{\mathrm{1}}\mathrm{、}{u}_{\mathrm{2}}为接近于零的特定数值，本文中u_{\mathrm{1}}=u_{\mathrm{2}}=\frac{\mathrm{1}}{\mathrm{35}}\times {\mathrm{10}}^{-\mathrm{8}}。

薄膜单元的应变能表示为

U^{\mathrm{m}\mathrm{e}\mathrm{m}}={\int }_V\frac{\partial [\begin{array}{cc}{\gamma }_{\mathrm{1}}& {\gamma }_{\mathrm{2}}\end{array}]}{\partial \mathit{e}}\tilde{\mathit{\Gamma }}\left[\begin{array}{c}{\gamma }_{\mathrm{1}}\\ {\gamma }_{\mathrm{2}}\end{array}\right]\mathrm{d}V（21）

弹性力可表示为

{\mathit{F}}^{\mathrm{m}\mathrm{e}\mathrm{m}}=\frac{\partial U^{\mathrm{m}\mathrm{e}\mathrm{m}}}{\partial \mathit{e}}（22）

弹性力雅可比公式表示为

{\mathit{J}}^{\mathrm{m}\mathrm{e}\mathrm{m}}=\frac{\partial {\mathit{F}}^{\mathrm{m}\mathrm{e}\mathrm{m}}}{\partial \mathit{e}}=\frac{{\partial }^{\mathrm{2}}{U}^{\mathrm{m}\mathrm{e}\mathrm{m}}}{\partial {\mathit{e}}^{\mathrm{2}}}=

{\int }_V\frac{{\partial }^{\mathrm{2}}[\begin{array}{cc}{\gamma }_{\mathrm{1}}& {\gamma }_{\mathrm{2}}\end{array}]}{\partial {\mathit{e}}^{\mathrm{2}}}\tilde{\mathit{\Gamma }}\left[\begin{array}{c}{\gamma }_{\mathrm{1}}\\ {\gamma }_{\mathrm{2}}\end{array}\right]\mathrm{d}V+{\int }_V\frac{\partial [\begin{array}{cc}{\gamma }_{\mathrm{1}}& {\gamma }_{\mathrm{2}}\end{array}]}{\partial \mathit{e}}\tilde{\mathit{\Gamma }}\frac{\partial \left[\begin{array}{c}{\gamma }_{\mathrm{1}}\\ {\gamma }_{\mathrm{2}}\end{array}\right]}{\partial \mathit{e}}\mathrm{d}V（23）

3 薄膜天线建模与结果讨论

将薄膜式天线结构简化为薄膜-撑杆混合模型，如图3所示。OXYZ为全局坐标系，Puvw为以模型中心P为原点的固定坐标系。

薄膜采用ANCF薄板单元建模，应力应变的本构关系采用第2节描述的薄膜本构关系构造；两根撑杆采用中空圆截面梁单元建模，交叉布置于方形薄膜的对角线处，对薄膜起到支撑作用。薄膜和撑杆之间在A、B、C、D这4个节点处及薄膜中心单元E、F、G、H这4个节点处的相对位置采用耦合约束，添加约束方程，即{\mathit{e}}_i^{\mathrm{B}}-{\mathit{e}}_j^{\mathrm{P}}=\mathrm{0}，{\mathit{e}}_i^{\mathrm{B}}代表梁单元节点i的位置向量，{\mathit{e}}_j^{\mathrm{P}}代表薄膜单元节点j的位置向量，薄膜单元的中心点P采用球铰约束。薄膜-撑杆混合模型的物理参数和材料参数如表1所示。

3.1 薄膜天线下摆算例

如图4所示，薄膜-撑杆混合模型一个角点采用球铰约束，对角点为监测点，在微重力g=1 m/s²作用下摆动。采用商业软件ANSYS与本文提出的方法进行验证，薄膜采用shell181单元离散，离散网格为42×42，撑杆采用beam188单元，每根撑杆单元个数为200，仿真时间为2 s。撑杆划分不同单元的仿真结果如图5~图7所示，其中，2×2-4代表薄膜结构采用2×2 ANCF薄板单元离散，每个撑杆采用4个ANCF梁单元离散，其他同理。随着ANCF薄板单元和梁单元的增加，位移逐渐趋于收敛，而且与商业软件结果相比，整体结果吻合较好，误差不超过0.02 m。

3.2 薄膜天线机动控制研究

在以往的文献中，研究人员引入绝对节点坐标参考节点描述航天器结构中的刚体，可以建立刚柔耦合模型并选取和定义刚体转动。通过前文的建模过程不难看出，薄膜天线是一个大型柔性结构，在对其进行主动控制时，如何选取表征其姿态的角度是一个难题。绝对节点坐标法的优势在于其描述的柔性体上任意一点的位置向量梯度矩阵同时包含刚体转动和柔性体变形信息，通过极化分解的方法即可实现解耦。根据极化分解定理，任一非奇异矩阵均可以分解为一个正交矩阵和一个对称矩阵相乘的形式，具体形式如下：

{\mathit{J}}^{\mathbf{R}}=\mathit{R}\mathit{U}（24）

式中： J^R为柔性体上任意点的位置向量梯度矩阵； U 为对称正定的拉伸矩阵； R 为正交的转动矩阵。

在本研究中，选择薄膜天线中心点处位置向量梯度矩阵作极化分解，用得到的旋转变化矩阵 R 表征薄膜天线整体转动。该过程可以视作模型由参考构型分别绕自身Puvw坐标系依次旋转合成：

\mathit{R}={\mathit{R}}_w(\gamma )\cdot {\mathit{R}}_v(\beta )\cdot {\mathit{R}}_u(\theta )（25）

式中：{\mathit{R}}_u(\theta )表示绕原始u轴旋转角度\theta；{\mathit{R}}_v(\beta )表示绕原始v轴旋转角度\beta；{\mathit{R}}_w(\gamma )表示绕原始w轴旋转角度\gamma。

PD控制方法是工程上最常用的简单且有效的控制方法，PD控制数学模型简单，具有响应快速、抗干扰能力强的特点，因此，将旋转矩阵转换成欧拉角后，利用目标角度和当前角度的差值，采用PD控制方法对模型姿态进行调整。

设置两种自旋工况验证模型及控制算法的有效性，即分别对模型中心点P施加使模型绕着固定轴w轴及u轴旋转的驱动力矩，力矩大小采用PD反馈控制实时调整，仿真时间为40 s，在Puvw坐标系下模型整体相对于三个轴的目标旋转角度如表2所示。薄膜的网格划分结构为3×3（长×宽），每根圆截面梁撑杆含有6个单元。

模型整体的旋转角度检测点为薄膜中心点P。采用基于广义\alpha的系统主动控制算法进行求解，算法流程如图8所示。

图8中，t₀、t_p分别为动力学仿真的开始和结束时间；a_n 为中间类加速度变量；t_tol为算法收敛阈值。在动力学求解部分，{\mathit{G}}_{\mathrm{1}}是力向量的残差；\beta是广义\alpha求解参数；{\beta }^{\mathrm{\text{'}}}是另一个用来表示伪加速度和速度之间关系的广义\alpha的计算参数；h 是时间步长；\mathit{M}是系统质量矩阵​；\mathit{J}是力雅可比矩阵；\mathrm{\Delta }{\mathit{q}}_n^k对应节点坐标的增量；{\mathit{Q}}_{\mathrm{e}}为系统的弹性力；{\mathit{Q}}_{\mathrm{e}\mathrm{x}\mathrm{t}}为恒定外力；{\mathit{Q}}_{\mathrm{c}}(t)为引入的经过PD控制的时变外力；{\mathit{D}}_{\mathrm{0}}为模型整体在Puvw坐标系下相对于三个轴的初始旋转角度（0°，0°，0°）；{\mathit{D}}_n为经过PD控制后返回的模型实时旋转角度。

图9~图12所示为工况1的仿真结果。其中，图9为薄膜天线模型经PD控制绕自身w轴自旋稳定过程中不同时刻的构型图，薄膜的不同颜色变化表示薄膜张紧、褶皱、松弛三种应力状态的变化。图10为模型绕u、v、w三个轴旋转角度的变化曲线。图11为模型在w轴方向的控制力矩变化曲线。图12为模型绕自身w轴自旋过程中，点A在全局坐标系下沿着X轴方向及Y轴方向的位移变化曲线。工况1的仿真结果表明，薄膜天线整体做平面内转动时表现出了较好的刚度和控制效果，收敛后角度误差小于0.2%。

图13~图15所示为工况2的仿真结果。其中，图13为模型绕自身u轴自旋稳定过程中不同时刻的构型图，图14为模型绕u、v、w三个轴旋转角度的变化曲线，图15为模型绕u轴自旋过程中，点A在全局坐标系下的位移变化曲线。结果表明，当施加平面外的驱动力矩后，薄膜天线由于整体结构平面外抗弯刚度相对较弱，出现了一定程度的振动，但角度控制误差依旧不超过2.5%。两种不同工况的仿真结果均验证了本文提出的薄膜天线混合单元建模、变形/姿态耦合动力学和控制方法的有效性。

4 结束语

本文建立了一种基于绝对节点坐标法混合单元的空间薄膜天线柔性体变形/姿态耦合动力学与控制模型。其中，薄膜本体采用薄板单元离散并引入薄膜本构描述其张紧-松弛的状态分布，而交叉撑杆则采用空心圆截面梁单元进行建模。二者通过节点自由度耦合的方式进行组装。选取薄膜天线中心节点对其位置向量梯度矩阵进行极化分解，获取其解耦的柔性变形/刚体转动信息，从而解算出欧拉角描述的薄膜天线整体姿态，最终实现基于PD控制的天线角度调整。但是因为弹性力及其雅可比矩阵高度非线性，整体计算效率较低，不适用于高精度实时控制的场景。本文提出的方法可为以空间薄膜天线为代表的大型柔性航天器的主动控制提供参考。

参考文献

[1] 彭福军，谢超，张良俊. 面向空间应用的薄膜可展开结构研究进展及技术挑战［J］. 载人航天，2017， 23（4）： 427-439.

[2] PENG Fujun， XIE Chao， ZHANG Liangjun. Research Progress and Technological Challenges of Thin Film Deployable Structures for Space Applications ［J］. Journal of Astronautics， 2017， 23（4）： 427-439.

[3] 李鹏飞.薄膜结构褶皱分析及抑制方法研究［D］.西安：西安电子科技大学，2023.

[4] LI Pengfei. Study on Folding Analysis and Suppression Method of Thin Film Structure ［D］. Xi'an：Xidian University， 2023.

[5] FANG H， LOU M， HUANG J， et al. Development of a Three-meter Ka-band Reflectarray Antenna ［R］.Denver： AIAA Electronic Library， 2002.

[6] HUANG J， FANG H， RICHARD L， et al. The Development of Large Flat Inflatable Antenna for Deep-space Communications ［R］. San Diego： AIAA Electronic Library， 2004.

[7] FANG H， KNARR K， QUIJANO U， et al. In-space Deployable Deflectarray Antenna ［R］. Schaumburg： AIAA Electronic Library， 2008.

[8] HEALD J， POTVIN M， JIANG X. Experimental Investigations to Support a Multi-layer Deployable Membrane Structure for Space Antennae ［R］. Austin： AIAA Electronic Library， 2005.

[9] SHEN Y， MONTMITY S， ZHENG W， et al. Large SAR Membrane Antenna Deployable Structure Design and Dynamic Simulation ［R］. Honolulu： AIAA Electronic Library， 2007.

[10] JOACHIM B， MSINNAPIUS M. Deployable Composite Booms for Various Gossamer Space Structures ［R］. Denver： AIAA Electronic Library， 2011.

[11] STRAUBEL M， SEEFELDT P， SPIEYTZ P， et al. The Design and Test of the Gossamer-1 Boom Deployment Mechanisms Engineering Model ［R］. Kissimmee： AIAA Electronic Library， 2015.

[12] 肖薇薇，陈务军，付功义. 空间薄膜阵面结构褶皱分析［J］. 宇航学报， 2010，31 （11）：2604-2609.

[13] XIAO Weiwu， CHEN Wujun， FU Gongyi. Analysis of Folded Structures on the Surface of Space Film Membrane ［J］. Journal of Astronautics， 2010， 31（11）：2604-2609.

[14] 胥小勇，孙宇，蒋清海.薄膜张力控制系统的建模与设计［J］.中国机械工程，2013，24（18）：2452-2457.

[15] XU Xiaoyong， SUN Yu， JIANG Qinghai. Modeling and Design of Film Tension Control System ［J］. China Mechanical Engineering， 2013， 24（18）： 2452-2457.

[16] 操安博，林秋红，邱慧，等. 空间薄膜结构粘接缝模型的建立与验证［J］. 宇航学报， 2023， 44（11）：1777-1784.

[17] CAO Anbo， LIN Qiuhong， QIU Hui， et al. Establishment and Verification of the Bonding Seam Model for Space Film Structures ［J］. Journal of Astronautics， 2023， 44（11）： 1777-1784.

[18] 荣吉利，宋逸博，刘志超，等. 圆形薄膜太阳翼展开动力学分析与模态分析［J］. 宇航学报， 2020， 41（9）：1125-1131.

[19] RONG Jili， SONG Yibo， LIU Zhichao， et al. Dynamic Analysis and Modal Analysis of Circular Film Solar Wing Deployment ［J］. Acta Aeronautica et Astronautica Sinica， 2020， 41（9）： 1125-1131.

[20] 陈金铎，史爱明，樊世超. 考虑光压力和热载荷的太阳帆薄膜力学行为分析［J］. 宇航学报， 2024， 45（1）：74-80.

[21] CHEN Jinduo， SHI Aiming， FAN Shichao. Analysis of the Mechanical Behavior of Solar Sail Films Considering Light Pressure and Thermal Load ［J］. Journal of Astronautics， 2024， 45（1）： 74-80.

[22] SHABANA A A. Dynamics of Multibody Systems ［M］. Chicago： Cambridge University Press， 2005.

[23] 赵将，刘铖，田强，等. 黏弹性薄膜太阳帆自旋展开动力学分析［J］. 力学学报，2013，45（5）：746-752.

[24] ZHAO Jiang， LIU Cheng， TIAN Qiang，et al. Spin Deployment Dynamics Analysis of Viscoelastic Membrane Solar Sail［J］. Chinese Journal of Theoretical and Applied Mechanics，2013，45（5）：746-752.

[25] LIU Cheng， TIAN Qiang， HU Haiyan，et al. Dynamic Analysis of Membrane Systems Undergoing Overall Motions，Large Deformations and Wrinkles via Thin Shell Elements of ANCF［J］. Computer Methods in Applied Mechanics and Engineering，2013，258：81-95.

[26] 薛鹏聪，刘铖，水小平. 基于ANCF的薄膜太阳帆自旋展开动力学模拟［J］.深空探测学报（中英文），2022，9（2）：217-229.

[27] XUE Pengcong， LIU Cheng， SHUI Xiaoping. Kinematic Simulation of Spin Deployment of Thin Film Solar Sail Based on ANCF ［J］. Journal of Deep Space Exploration （Chinese and English）， 2022， 9（2）： 217-229.

[28] 支敬德，戈新生.基于模糊滑模控制的挠性航天器姿态机动及抖振抑制研究［J］.应用力学学报，2020，37（5）：1972-1979+2317-2318.

[29] ZHI Jingde， GE Xinsheng. Attitude Maneuver and Buffeting Suppression of Flexible Spacecraft based on Fuzzy Sliding Mode Control［J］. Chinese Journal of Applied Mechanics，2020，37（5）：1972-1979+2317-2318.

[30] 尹一蓁，康国华，武俊峰，等.基于非线性变结构滑模控制器的薄膜航天器姿态机动控制研究［J］.上海航天（中英文）， 2024， 41（6）：95-105.

[31] YIN Yizhen， KANG Guohua， WU Junfeng， et al. Research on Attitude Maneuver Control of Thin-Film Spacecraft Based on Nonlinear Variable Structure Sliding Mode Controller［J］. Aerospace Shanghai （Chinese and English）， 2024， 41（6）：95-105.

[32] 殷春武. 航天器姿态控制方法及其关键技术综述［J］.海军航空大学学报，2024，39（6）：660-670.

[33] YIN Chunwu. Review on Spacecraft Attitude Control Method and Its Key Technologies［J］. Journal of Naval Aviation University， 2024， 39（6）：660-670.

[34] LIU， Yuhang， LUO Kai， TIAN Qiang， et al. Nonlinear Dynamics Design for In-space Assembly Motion of Manipulators on Flexible Base Structures［J］. Nonlinear Dynamics 2025，113，9485-9507.

[35] MIKKOLA M A， SHABANA A A. A Non-Incremental Finite Element Procedure for the Analysis of Large Deformation of Plates and Shells in Mechanical System Applications ［J］. Multibody System Dynamics，2003，9（3）：283-309.

[36] DUFVA K， SHABANA A A. Analysis of Thin Plate Structures Using the Absolute Nodal Coordinate Formulation［J］. Proceedings of the Institution of Mechanical Engineers， Part K： Journal of Multi-body Dynamics，2005，219（4）：345-355.

[37] SHABANA A A， YAKOUB R Y. Three Dimensional Absolute Nodal Coordinate Formulation for Beam Elements Implementation and Application［J］. Journal of Mechanical Design， 2001， 123：614-621.

[38] LAN Peng， TIAN Qinglong， YU Zuqing. A New Absolute Nodal Coordinate Formulation Beam Element with Multilayer Circular Cross Section ［J］. Acta Mechanica Sinica， 2020， 36（1）： 82-96.

[39] MIYAZAKI Y. Wrinkle/Slack Model and Finite Element Dynamics of Membrane［J］. International Journal for Numerical Methods in Engineering， 2006， 66（7）： 1179-1209

基金资助

国家自然科学基金(12272123)

国家自然科学基金(12302047)

江苏省基础研究计划自然科学基金(BK20231185)

江苏省青年科技人才托举工程项目(JSTJ-2024-481)