CN103412491B - 一种挠性航天器特征轴姿态机动指数时变滑模控制方法 - Google Patents

Abstract

本发明涉及一种挠性航天器特征轴姿态机动指数时变滑模控制方法，属于航天器控制技术领域。首先在航天器本体系下建立系统等效动力学模型、运动学模型和挠性振动模型，然后计算带有指数时变滑模控制律的闭环系统振动频率和阻尼比参数，根据单轴输入成型器设计方法，设计以特征轴为旋转轴的单轴多模态滤波输入成型器，抑制三轴运动中的挠性振动。同时，设计状态观测器实时估计挠性模态信息，构成输出反馈指数时变滑模控制方法。最后进行控制力矩的饱和性分析，以满足控制力矩的物理饱和约束。本发明扩大了现有输入成型的应用范围，将输入成型技术由单轴机动扩展到了三轴机动过程当中，增强了滤波输入成型自身鲁棒性，实现了航天器的姿态机动路径最短。

Description

一种挠性航天器特征轴姿态机动指数时变滑模控制方法

技术领域

本发明涉及一种挠性航天器特征轴姿态机动指数时变滑模控制方法，特别涉及基于滤波输入成型的挠性航天器特征轴姿态机动指数时变滑模控制方法，属于航天器控制技术领域。

背景技术

为了节约发射成本，现代航天器一般配有太阳能电池帆板、天线等其他轻质挠性结构附件。这种刚柔耦合的结构设计会导致航天器在快速姿态机动时出现强烈的挠性振动，继而影响姿态机动的控制精度，导致系统性能下降。长期的振动还会造成结构的疲劳损坏。因此，抑制挠性附件的残余振动就显得尤为重要。

目前，输入成型技术作为一种前馈控制方法，在挠性振动控制方面受到各国学者的青睐。Singhose[Singhose W E，Derezinski S，Singer N C.Extra-insensitive inputshapers for controlling flexible spacecraft[J].Journal of Guidance，Control，and Dynamics，1996，2：1122-1130]等人利用输入成形技术来抑制航天器大型挠性附件带来的残余振动，取得了很好的效果。孔宪仁等人[孔宪仁，杨正贤，叶东等.基于输入成形的柔性航天器振动闭环抑制方法研究[J]，振动与冲击，2010，29(3)：72-76]提出了一种将输入成形与PD反馈控制相结合的控制策略，在保证航天器完成姿态机动的同时有效地抑制了挠性附件的振动。

可是，输入成型器(IS)的设计需要精确的模型信息，模型误差会影响挠性振动的抑制效果。在一般情况下，挠性航天器的模型参数比如惯量阵、挠性模态频率和阻尼比等信息是部分未知的，外界干扰和参数不确定性也会严重影响输入成型器的振动抑制效果。因此，在外界干扰和参数不确定性的影响下，寻找一种鲁棒抑制挠性振动的方法已经变得十分重要。滑模控制是一种解决外界干扰和参数不确定性问题十分有效的控制方法，由于它的鲁棒性和简单性，滑模控制已经大量应用在航天控制领域上。胡庆雷[Hu Q L，Wang Z D，Gao H J.Sliding mode and shaped input vibration control of flexible systems[J].IEEE Transactions on Aerospace and Electronic Systems，2008，44(2)：503-519]提出了一种滑模控制与输入成型相结合的控制策略，解决了在参数不确定性和外部干扰影响下的挠性航天器大角度姿态机动和振动抑制问题。可是，因为聚合不确定的上界未知，滑模切换增益一般很难准确得到。针对这个问题，文[苗双全，丛炳龙，刘向东.基于输入成形的挠性航天器自适应滑模控制[J].航空学报，2013]提出了一种自适应滑模与输入成型相结合的控制方法，在抑制掉挠性振动的基础上，有效减小了因为切换增益取值过大而带来的滑模抖动。

但是，目前只是将输入成型器应用在单轴航天器模型上，而实际的航天器模型是三维的，执行的任务很多情况下都是三轴机动，所以，仅考虑单轴机动将给实际的航天器控制带来很多的局限性。三轴挠性航天器姿态机动和振动抑制的难点在于，三轴之间以及三轴和挠性振动模态之间的复杂耦合关系，使得很难用解析的方法单独对每一个轴进行输入成型的设计。因为上述原因，挠性航天器三轴机动的振动抑制问题很难用输入成形来解决，造成输入成形应用的局限性。

发明内容

本发明的目的是为抑制挠性航天器三轴姿态机动过程中因挠性附件产生的残余振动，提出一种基于滤波输入成型的挠性航天器特征轴姿态机动指数时变滑模控制方法，利用滑模控制理论提高输入成型器的鲁棒性。

本发明所采用的技术方案为：首先在航天器本体系下建立系统等效动力学模型、运动学模型和挠性振动模型，然后利用指数时变滑模控制，一方面在外界干扰和参数不确定性的影响下实现高精度的期望姿态机动，另一方面保证航天器为特征轴姿态运动。根据欧拉旋转定理，任一姿态变化都可以由绕特征轴旋转一个角度实现，因此采用单轴输入成型器设计方法，计算出带有指数时变滑模控制律的闭环系统振动频率和阻尼比参数，设计多模态滤波输入成型器，抑制三轴运动中的挠性振动。同时，为了解决挠性模态不可测的问题，设计状态观测器实时估计挠性模态信息，构成输出反馈指数时变滑模控制方法。最后进行控制力矩的饱和性分析，以满足控制力矩的物理饱和约束。

具体包括以下步骤：

步骤1，在航天器姿态运动的构型空间内定义相对姿态变量和挠性振动变量，在航天器本体系下建立挠性航天器动力学方程、姿态运动学方程和振动方程。具体方法为：

挠性航天器动力学方程为：

振动方程为：

式中为系统实际的正定对称转动惯量矩阵，为系统名义惯量阵，ΔJ为由系统质量变化引起的惯量阵误差。ω＝[ω₁ ω₂ ω₃]^T为在本体系下的刚体角速度表示，ω^×为ω的斜对称矩阵。T_c＝[T_c1 T_c2 T_c3]^T和T_d＝[T_d1 T_d2 T_d3]^T分别表示控制力矩和外界干扰力矩。η为挠性结构弹性形变的广义坐标，ζ和Λ分别为挠性附件的模态阻尼矩阵和模态频率矩阵，C为挠性附件与星体的刚柔耦合矩阵。

姿态误差σ_e和角速度误差ω_e为：

${\mathrm{\σ}}_e=\mathrm{\σ}\⊗{{\mathrm{\σ}}_d}^{-1}=\frac{({\mathrm{\σ}}^T\mathrm{\σ}-1){\mathrm{\σ}}_d+(1-{{\mathrm{\σ}}_d}^T{\mathrm{\σ}}_d)\mathrm{\σ}+2{{\mathrm{\σ}}_d}^{\×}\mathrm{\σ}}{1+\left({\mathrm{\σ}}^T\mathrm{\σ}\right)\left({{\mathrm{\σ}}_d}^T{\mathrm{\σ}}_d\right)+2\left({\mathrm{\σ}}^T{\mathrm{\σ}}_d\right)}---\left(3\right)$

${\mathrm{\ω}}_e=\mathrm{\ω}-{\mathrm{\ω}}_d^b---\left(4\right)$

式中σ为当前姿态角，σ_d为期望姿态角，σ_d ^×为σ_d的斜对称矩阵，为航天器在本体系下的期望角速度，ω_d为航天器在惯性系下的期望角速度，为从惯性系到本体系的转移矩阵，σ^×为σ的斜对称矩阵。

姿态运动学方程在本体系下表示为

$\stackrel{\·}{\mathrm{\σ}}=M\left(\mathrm{\σ}\right)\mathrm{\ω}---\left(5\right)$

式中满足条件M^TM＝m_aI_3×3，其中I_3×3为3×3的单位矩阵。

设φ＝[η^T ψ^T]^T为航天器挠性振动变量，其中令D＝2ζΛ，K＝Λ²。对式(1)、(2)和(5)进行整理，得到等效的数学模型：

式中E＝[K D]，表示由惯量阵不确定性和外界干扰引起的聚合扰动，其中I表示单位阵，d约束于未知上界d_max。

步骤2，针对步骤1建立的等效数学模型，设计状态反馈指数时变滑模控制律，一方面保证系统在外界干扰和参数不确定性的影响下，实现精确的姿态机动，另一方面，希望能够实现特征轴姿态机动，以便于后续的进一步控制。具体方法为：

设计滑模面函数为

$S({\mathrm{\ω}}_e,{\mathrm{\σ}}_e,t)={\mathrm{\ω}}_e+l\frac{M^T}{m_a}{\mathrm{\σ}}_e+\frac{M^T}{m_a}{\mathrm{Qe}}^{-\mathrm{at}}---\left(7\right)$

式中l和a为正的标量。当t→∞时，S(ω_e，σ_e，t)＝0。根据时变滑模控制理论，系统初始状态在滑模面上，即S(ω_e，σ_e，0)＝0，则有Q＝-lσ_e(0)，σ_e(0)即σ_e在0时刻的状态值。

设计状态反馈指数时变滑模控制律如下

$T_c^{\′}=-\mathrm{CE\φ}+{\mathrm{CDC}}^T{\mathrm{\ω}}_e+{\mathrm{\ω}}^{\×}\hat{J}\mathrm{\ω}+{\hat{J}}_m{\stackrel{\·}{\mathrm{\ω}}}_d^b-l{\hat{J}}_m\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)-{\hat{J}}_m\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-\mathrm{at}}\right)-\mathrm{\γsgn}\left(S\right)---\left(8\right)$

γ为滑模的切换增益，γ＞0，E＝[K D]；sgn(S)为符号函数，当S大于0时，sgn(S)为1，当S的等于0时，sgn(S)为0，当S小于0时，sgn(S)为-1。

步骤3，利用步骤2实现的特征轴运动特性，将系统三轴运动转变为绕特征轴旋转的单轴运动，根据现有的单轴输入成型器设计方法，计算出带有步骤2中指数时变滑模控制律的闭环系统振动频率和阻尼比参数，最后根据这些参数设计多模态滤波输入成型器，抑制三轴运动中的挠性振动。

多模态滤波输入成型器的具体设计方法为：

设状态变量为X＝[σ^T ω^T η^T ψ^T]^T，对闭环系统建立状态空间方程，求出系数矩阵的特征值λ_sys，再通过下式求出闭环系统的振动频率和阻尼比信息：

${\mathrm{\λ}}_{\mathrm{sys}}=-{\mathrm{\ζ}}_{\mathrm{sys},i}{\mathrm{\ω}}_{\mathrm{sys},i}\±{\mathrm{j\ω}}_{\mathrm{sys},i}\sqrt{1-{{\mathrm{\ζ}}_{\mathrm{sys},i}}^2},i=\mathrm{1,2},\·\·\·,n---\left(9\right)$

式中ζ_sys，i和ω_sys，i分别为第i阶闭环系统的阻尼比和振动频率。

ZVD即输入成型器形式为

$\begin{array}{ccc}{A}_1=\frac{1}{1+2K^{\′}+K^{\′2}}& A_2=\frac{2K^{\′}}{1+2K^{\′}+K^{\′2}}& A_3=\frac{K^{\′2}}{1+2K^{\′}+K^{\′2}}\\ T_1=0& T_2=\frac{\mathrm{\π}}{{\mathrm{\ω}}_{\mathrm{sys},i}\sqrt{1+{{\mathrm{\ζ}}_{\mathrm{sys},i}}^2}}& T_3=\frac{2\mathrm{\π}}{{\mathrm{\ω}}_{\mathrm{sys},i}\sqrt{1-{{\mathrm{\ζ}}_{\mathrm{sys},i}}^2}}\end{array}---\left(10\right)$

式中A₁，A₂，A₃和T₁，T₂，T₃分别为脉冲的幅值和时间。

对于多模态系统的振动抑制问题，输入成型器的脉冲序列通过各个单模态的脉冲序列相卷积得到，公式如下：

A_mult＝A_mult1*A_mult2*…A_mulj…*A_multn

式中A_multi代表第j阶单模态的脉冲序列，*代表卷积运算，j＝1，2，…，n；

将设计好的输入成型器与一个一阶惯性环节相连接(即相乘)，构成多模态滤波输入成型器(FIS)。FIS不仅可以提高抑制残余振动的能力，而且可以极大改善控制力矩的输出特性，避免因输入成型器的脉冲响应造成控制力矩的阶跃跳变，使得控制力矩变得更加平滑。所述一阶惯性环节的形式如下

$H\left(s\right)=\frac{1}{T_{\mathrm{ci}}s+1}---\left(11\right)$

式中T_ci为惯性时间常量。

步骤4，因为步骤2是在挠性状态可测的情况下进行设计的。但是一般情况下挠性模态不可测，因此引入挠性状态观测器，在线估计挠性模态信息，设计输出反馈指数时变滑模控制律，解决挠性模态不可测量的问题，对步骤2提出的姿态控制律进一步改进。对控制律进行饱和性分析，使控制律满足执行器饱和约束。具体方法为：

设计挠性状态观测器如下

$\stackrel{\·}{\widehat{\mathrm{\φ}}}=A\widehat{\mathrm{\φ}}-\mathrm{AB}{C}^T{\mathrm{\ω}}_e-B{C}^T{\stackrel{\·}{\mathrm{\ω}}}_d^b+P^{-1}{\left(S^T\mathrm{CE}\right)}^T---\left(12\right)$

其中为航天器挠性振动变量φ的估计值，P为观测器正定对称增益矩阵，且PA＜0。则输出反馈指数时变滑模控制律如下

${T_c}^{\′\′}-CE\widehat{\mathrm{\φ}}+{\mathrm{CDC}}^T{\mathrm{\ω}}_e+{\mathrm{\ω}}^{\×}\hat{J}\mathrm{\ω}+{\hat{J}}_m{\stackrel{\·}{\mathrm{\ω}}}_d^b-l{\hat{J}}_m\frac{d}{dt}(\frac{M^T}{m_a}-{\mathrm{\σ}}_e)-{\hat{J}}_m\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}\right)-\mathrm{\γ}\mathrm{sgn}\left(S\right)---\left(13\right)$

将输出反馈指数时变滑模控制律输入航天器模型，航天器模型在该控制律控制下得到σ，将σ输入挠性状态观测器，然后将在挠性状态观测器下得到的输入控制器，形成闭环系统。

所述闭环系统包括控制器、挠性状态观测器和航天器模型。

有益效果

本发明方法能够有效地抑制挠性航天器机动过程中的残余振动，与现有技术相比的优点在于：

1)本发明扩大了现有输入成型的应用范围，将单轴机动中振动抑制扩展到了三轴机动当中。设计者避免了三轴机动与挠性振动模态之间的复杂耦合关系，不需要针对每个轴设计独立的输入成型器，降低了计算难度，只需针对绕特征轴旋转的单轴机动设计单轴输入成型，就可以解决三轴旋转运动的挠性振动抑制问题。

2)输出反馈指数时变滑模控制与滤波输入成型的结合，不仅弥补了滤波输入成型自身鲁棒性不强的缺点，而且实现了航天器的姿态机动路径最短。引入的挠性状态观测器可以实时估计挠性模态信息，解决了挠性模态不可测量的问题。此外，控制力矩的饱和分析使得控制器输出满足执行机构的物理饱和约束。

附图说明

图1为本发明基于滤波输入成型(FIS)的挠性航天器特征轴机动输出反馈指数时变滑模控制律的设计流程图；

图2为具体实施中基于滤波输入成型的挠性航天器特征轴机动输出反馈指数时变滑模控制律的设计框图；

图3为具体实施中挠性航天器姿态机动路径图，其中：(a)为PD+IS控制下的系统姿态机动路径曲线；(b)为OFETVSMC+FIS控制下的系统姿态路径曲线。

图4为PD控制、PD+IS控制和OFETVSMC+FIS控制下的姿态角响应曲线。其中(a)为在三种控制律作用下姿态角分量σ₁曲线图，(b)为在三种控制律作用下姿态角分量σ₂曲线图，(c)为在三种控制律作用下姿态角分量σ₃曲线图；

图5为PD控制、PD+IS控制和OFETVSMC+FIS控制下的前三阶挠性坐标曲线。其中(a)为在三种控制律作用下第一阶挠性模态坐标η₁曲线图，(b)为在三种控制律作用下第二阶挠性模态坐标η₂曲线图，(c)为在三种控制律作用下第三阶挠性模态坐标η₃曲线图；

图6为PD控制、PD+IS控制和OFETVSMC+FIS控制下的控制力矩曲线。其中(a)为在三种控制律作用下控制力矩分量T₁曲线图，(b)为在三种控制律作用下控制力矩分量T₂曲线图，(c)为在三种控制律作用下控制力矩分量T₃曲线图。

具体实施方式

为了更好地说明本发明的目的和优点，下面结合附图和实施例加以进一步阐述。

本发明的设计流程图如图1所示。依据本发明方法实施的基于滤波输入成型的挠性航天器特征轴机动输出反馈指数时变滑模控制律设计框图如图2所示，该控制律能够有效地抑制挠性航天器机动过程中的残余振动。

步骤1，在航天器姿态运动的构型空间内定义相对姿态变量和挠性振动变量，在航天器本体系下建立挠性航天器动力学方程、姿态运动学方程和振动方程。具体方法为：

挠性航天器动力学方程和振动方程如下所示：

$\stackrel{..}{\mathrm{\η}}+2\mathrm{\ζ\Λ}\stackrel{.}{\mathrm{\η}}+{\mathrm{\Λ}}^2\mathrm{\η}+C^T\stackrel{.}{\mathrm{\ω}}=0$

式中为系统实际的正定对称转动惯量矩阵，为系统名义惯量阵，ΔJ为由系统质量变化引起的惯量阵误差。ω＝[ω₁ ω₂ ω₃]^T为在本体系下的刚体角速度表示，ω^×为ω的斜对称矩阵。T_c＝[T_c1 T_c2 T_c3]^T和T_d＝[T_d1 T_d2 T_d3]^T分别表示控制力矩和外界干扰力矩。η为挠性结构弹性形变的广义坐标，ζ和Λ分别为挠性附件的模态阻尼矩阵和模态频率矩阵，C为挠性附件与星体的刚柔耦合矩阵。

姿态误差σ_e和角速度误差ω_e为：

${\mathrm{\σ}}_e=\mathrm{\σ}\⊗{{\mathrm{\σ}}_d}^{-1}=\frac{({\mathrm{\σ}}^T\mathrm{\σ}-1){\mathrm{\σ}}_d+(1-{{\mathrm{\σ}}_d}^T{\mathrm{\σ}}_d)\mathrm{\σ}+2{{\mathrm{\σ}}_d}^{\×}\mathrm{\σ}}{1+\left({\mathrm{\σ}}^T\mathrm{\σ}\right)\left({{\mathrm{\σ}}_d}^T{\mathrm{\σ}}_d\right)+2\left({\mathrm{\σ}}^T{\mathrm{\σ}}_d\right)}---\left(16\right)$

${\mathrm{\ω}}_e=\mathrm{\ω}-{\mathrm{\ω}}_d^b---\left(17\right)$

式中σ为当前姿态角，σ_d为期望姿态角，σ_d ^×为σ_d的斜对称矩阵，为航天器在本体系下的期望角速度，ω_d为航天器在惯性系下的期望角速度，为从惯性系到本体系的转移矩阵σ^×为σ的斜对称矩阵。

姿态运动学方程在本体系下表示为

$\stackrel{\·}{\mathrm{\σ}}=M\left(\mathrm{\σ}\right)\mathrm{\ω}---\left(18\right)$

式中满足条件M^TM＝m_aI_3×3，其中I_3×3为3×3的单位矩阵。

设φ＝[η^T ψ^T]^T为航天器挠性振动变量，其中令D＝2ζΛ，K＝Λ²。对式(14)、(15)和(18)进行整理，得到等效的数学模型：

式中E＝[K D]，表示由惯量阵不确定性和外界干扰引起的聚合扰动，其中I表示单位阵，不失一般性，设d约束于未知上界d_max。

步骤2，针对步骤1建立的等效数学模型，设计状态反馈指数时变滑模控制律，一方面保证系统在外界干扰和参数不确定性的影响下，实现精确的姿态机动，另一方面，希望能够实现特征轴姿态机动，以便于后续的进一步控制。具体方法为：

针对步骤1建立的等效动力学模型，首先设计滑模面函数为

$S({\mathrm{\ω}}_e,{\mathrm{\σ}}_e,t)={\mathrm{\ω}}_e+l\frac{M^T}{m_a}{\mathrm{\σ}}_e+\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}---\left(20\right)$

式中l和a为正的标量。显然，当t→∞时，S(ω_e，σ_e，t)＝0。根据时变滑模控制理论，系统初始状态在滑模面上，即S(ω_e，σ_e，0)＝0，则有Q＝-lσ_e(0)，σ_e(0)即σ_e在0时刻的状态值。

设计状态反馈指数时变滑模控制律如下

${T_c}^{\′}=-CE\mathrm{\φ}+{\mathrm{CDC}}^T{\mathrm{\ω}}_e+{\mathrm{\ω}}^{\×}\hat{J}\mathrm{\ω}+J_m{\stackrel{\·}{\mathrm{\ω}}}_d^b-l{\hat{J}}_m\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)-{\hat{J}}_m\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}\right)-\mathrm{\γ}\mathrm{sgn}\left(S\right)---\left(21\right)$

γ为滑模的切换增益，γ＞0，现在对该控制律进行稳定性分析。定义正定李雅普诺夫函数为

$V=\frac{1}{2}{S}^T{\hat{J}}_{m}S---\left(22\right)$

滑模面函数S的导数为

$\stackrel{\·}{S}={\stackrel{\·}{\mathrm{\ω}}}_e+l\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)+\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}\right)$

然后对李雅普诺夫函数V求导

$\begin{array}{c}\stackrel{\·}{V}=S^T{\hat{J}}_m\stackrel{\·}{S}\\ =S^T({\hat{J}}_m{\stackrel{\·}{\mathrm{\ω}}}_e+l{\hat{J}}_m\frac{d}{dt}(\frac{M^T}{m_a}{\mathrm{\σ}}_e)+{\hat{J}}_m\frac{d}{dt}(\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}\left)\right)\\ =S^T(T_c+CE\mathrm{\φ}-{\mathrm{CDC}}^T{\mathrm{\ω}}_e-{\hat{J}}_m{\stackrel{\·}{\mathrm{\ω}}}_d^b-{\mathrm{\ω}}^{\×}\hat{J}\mathrm{\ω}+d(t)+l{\hat{J}}_m\frac{d}{dt}(\frac{M^T}{m_a}{\mathrm{\σ}}_e)+{\hat{J}}_m\frac{d}{dt}(\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}\left)\right)\\ =S^T\left(d\right(t)-\mathrm{\γ}\mathrm{sgn}(S\left)\right)\\ =S^{T}d\left(t\right)-\mathrm{\γ}\left|\right|S|{|}_1\\ \≤(d_{\max }-\mathrm{\γ})\left|\right|S|{|}_1\end{array}$

式中为向量1范数。假如满足γ＞d_max，则李雅普诺夫函数的导数为负，即根据李雅普诺夫稳定性定理，闭环系统是渐进稳定的。进一步考虑当(ω_e，σ_e)→(∞，∞)时，V→∞。因此闭环系统是全局渐进稳定的。

从上述分析可以看出，又因为当t＝0，李雅普诺夫函数V(ω_e，σ_e，t)＝0，所以V(ω_e，σ_e，t)≤0。又知道李雅普诺夫函数为正定函数，因此可以有以下结论

S(ω_e，σ_e，t)＝0当t≥0 (23)

上式表明系统状态始终保持在滑模面上，保证该控制方法具有全局鲁棒性。

接下来分析系统特征轴机动，从式(20)和(23)(23)可以得出

${\mathrm{\ω}}_e=-(l\frac{M^T}{m_a}{\mathrm{\σ}}_e+\frac{M^T}{m_a}{\mathrm{Qe}}^{-at})---\left(24\right)$

将式(10)代入式(35)中，经过化简有

${\stackrel{\·}{\mathrm{\σ}}}_e+{\mathrm{l\σ}}_e+{\mathrm{Qe}}^{-at}=0---\left(25\right)$

上面的二阶微分方程的解为

σ_e(t)＝f(t)σ_e(0) (26)

其中是一个标量函数。将式(26)代入式(24)中，有

ω_e＝g(t)σ_e(0) (27)

式中是一个标量函数。因为特征轴为而且ω(t)＝ω_e(t)，所以有ω(t)×n≡0。因此，可以看出角速度方向与特征轴始终共线，即实现了特征轴机动。

步骤3，将步骤2设计的控制律代入系统的动力学控制模型当中，得到带有指数时变滑模控制律的闭环系统，进而设计滤波输入成型器。

首先ZVD输入成型器形式如下

$\begin{array}{ccc}{A}_1=\frac{1}{1+2K^{\′}+K^{\′2}}& A_2=\frac{2K^{\′}}{1+2K^{\′}+K^{\′2}}& A_3=\frac{K^{\′2}}{1+2K^{\′}+K^{\′2}}\\ T_1=0& T_2=\frac{\mathrm{\π}}{\mathrm{\ω}\sqrt{1-{\mathrm{\ζ}}^2}}& T_3=\frac{2\mathrm{\π}}{\mathrm{\ω}\sqrt{1-{\mathrm{\ζ}}^2}}\end{array}---\left(28\right)$

式中ω和ζ分别为闭环系统的振动频率和阻尼比。

对于多模态系统的振动抑制问题，ZVD输入成型器的脉冲序列通过各个单模态的脉冲序列相卷积得到，公式如下：

A_mult＝A_mult1*A_mult2*…A_mulj…*A_multn (29)

式中A_mulj代表第j阶闭环系统单模态的脉冲序列，*代表卷积运算，j＝1，2，…，n。

本发明将输入成型与闭环反馈相结合，来抑制挠性航天器机动过程中的残余振，因此设计输入成型器所需的挠性信息应该为代入控制律之后闭环系统的挠性模态频率和阻尼比。首先设状态变量为X＝[σ^T ω^T η^T ψ^T]^T，然后对闭环系统建立状态空间方程，求出系数矩阵的特征值λ_sys，就可以通过下式求出闭环系统的振动频率和阻尼比参数。

${\mathrm{\λ}}_{sys}=-{\mathrm{\ζ}}_{sys,i}{\mathrm{\ω}}_{sys,i}\±{\mathrm{j\ω}}_{sys,i}\sqrt{1-{{\mathrm{\ζ}}_{sys,i}}^2},i=1,2,...,n---\left(30\right)$

式中ζ_sys，i和ω_sys，i分别为第i阶闭环系统振动阻尼和频率。利用求出的挠性振动信息，根据公式(28)和(29)，就可以求出多模态输入成型器。

将设计好的输入成型器与一个一阶惯性环节相连接，就可以构成多模态滤波输入成型器(FIS)。FIS不仅可以提高抑制残余振动的能力，而且可以极大改善控制力矩的输出特性，避免因输入成型器的脉冲响应造成控制力矩的阶跃跳变，使得控制力矩变得更加平滑。一阶惯性环节的形式如下

$H\left(s\right)=\frac{1}{T_{ci}s+1}---\left(31\right)$

式中T_ci为惯性时间常量。

步骤4，因为步骤2是在挠性状态可测的情况下进行设计的。但是一般情况下挠性模态不可测，因此引入挠性状态观测器，在线估计挠性模态信息，设计输出反馈指数时变滑模控制律，解决挠性模态不可测量的问题，对步骤2提出的姿态控制律进一步改进。并对控制律进行饱和性分析，使控制律满足执行器饱和约束。具体方法为：

针对步骤2设计的指数时变滑模控制律，设计状态观测器如下

$\stackrel{\·}{\widehat{\mathrm{\φ}}}=A\widehat{\mathrm{\φ}}-{\mathrm{ABC}}^T{\mathrm{\ω}}_e-{\mathrm{BC}}^T{\stackrel{\·}{\mathrm{\ω}}}_d^b+P^{-1}{\left(S^{T}CE\right)}^T---\left(32\right)$

其中为挠性振动变量φ的估计值，P为观测器正定对称增益矩阵，且PA＜0。矩阵A，B，C、E和滑模函数S(ω_e，σ_e，t)定义同上，则输出反馈指数时变滑模控制律设计如下

${T_c}^{\′\′}=-CE\widehat{\mathrm{\φ}}+{\mathrm{CDC}}^T{\mathrm{\ω}}_e+{\mathrm{\ω}}^{\×}\hat{J}\mathrm{\ω}+J_m{\stackrel{\·}{\mathrm{\ω}}}_d^b-l{\hat{J}}_m\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)-{\hat{J}}_m\frac{d}{dt}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-at}\right)-\mathrm{\γ}\mathrm{sgn}\left(S\right)---\left(33\right)$

式中参数定义同上。现对上述控制律进行稳定性分析。设挠性状态变量误差为李雅普诺夫函数定义如下

$V=\frac{1}{2}{S}^T{\hat{J}}_{m}S+\frac{1}{2}{e_{\mathrm{\φ}}}^T{\mathrm{Pe}}_{\mathrm{\φ}}---\left(34\right)$

对李雅普诺夫函数V求导可得

$\begin{array}{c}\stackrel{\·}{V}=S^T{\hat{J}}_m\stackrel{\·}{S}+{e_{\mathrm{\φ}}}^{T}P{\stackrel{\·}{e}}_{\mathrm{\φ}}\\ =S^T({\hat{J}}_m{\stackrel{\·}{\mathrm{\ω}}}_e+l{\hat{J}}_m\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)+{\hat{J}}_m\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-\mathrm{at}}\right))+{e_{\mathrm{\φ}}}^{T}P{\stackrel{\·}{e}}_{\mathrm{\φ}}\\ =S^T(T_c^{\′}+\mathrm{CE\φ}-{\mathrm{CDC}}^T{\mathrm{\ω}}_e-{\hat{J}}_m{\stackrel{\·}{\mathrm{\ω}}}_d^b-{\mathrm{\ω}}^{\×}\hat{J}\mathrm{\ω}+d+l{\hat{J}}_m\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)+{\hat{J}}_m\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{Qe}}^{-\mathrm{at}}\right))+{e_{\mathrm{\φ}}}^{T}P{\stackrel{\·}{e}}_{\mathrm{\φ}}\\ =S^T(\mathrm{CE}(\mathrm{\φ}-\widehat{\mathrm{\φ}})+d-\mathrm{\γsgn}\left(S\right))+{e_{\mathrm{\φ}}}^T(\mathrm{PA}(\mathrm{\φ}-\widehat{\mathrm{\φ}})-{\left(S^T\mathrm{CE}\right)}^T)\\ =S^T(d-\mathrm{\γsgn}\left(S\right))+{e_{\mathrm{\φ}}}^T\mathrm{PA}{e}_{\mathrm{\φ}}\\ \≤(d_{\max }-\mathrm{\γ}){\left|\right|S\left|\right|}_1+{e_{\mathrm{\φ}}}^T\mathrm{PA}{e}_{\mathrm{\φ}}\end{array}$

假如满足γ＞d_max，则李雅普诺夫函数导数为负，即此外，当S→∞，e_φ→∞时，有V→∞。根据李雅普诺夫稳定性理论，在控制律式(21)的作用下，该闭环系统为全局渐进稳定。特征轴证明与步骤3相同。上述理论分析表明所提出的输出反馈指数时变滑模控制不仅可以保证航天器完成期望姿态机动任务，而且实现了特征轴机动过程。

由于执行器输出的物理饱和限制，控制饱和已经成为工程上一个很普遍的问题。为了解决这个问题，本实例通过参数近似计算来整定控制器参数，解决了控制饱和问题，而且保证闭环系统渐进收敛。一般情况下，假设控制饱和要求为||T_c′||_∞≤T_max，其中T_max为执行器物理饱和量。因为角速度ω_e和挠性振动变量φ很小，所以多项式CEφ，CDC^Tω_e和可以忽略不计。

同时，因为成型器第一个脉冲幅值很小，导致初始角误差ω_e很小，使得Q很小。所以多项式也可以忽略不计。因此，控制律可以重新写为

$T_c^{\′\′}\≅l\tilde{J}\frac{d}{\mathrm{dt}}\left(\frac{M^T}{m_a}{\mathrm{\σ}}_e\right)-\mathrm{\γsgn}\left(S\right)---\left(35\right)$

将式(27)代入式(35)中，则

$T_c^{\′\′}=\frac{8l\stackrel{\·}{f}\hat{J}{\mathrm{\σ}}_e\left(0\right)}{{(1+f^2{{\mathrm{\σ}}_e}^T\left(0\right){\mathrm{\σ}}_e\left(0\right))}^2}(2M-f^2{{\mathrm{\σ}}_e}^T\left(0\right){\mathrm{\σ}}_e\left(0\right)I_{3\×3})-\mathrm{\γsgn}\left(S\right)---\left(36\right)$

式中f(t)定义同上，则有

$\left|\right|T_c^{\′}|{|}_{\mathrm{\∞}}\≤l|h\left(t\right)\left|\right|\left|\hat{J}{\mathrm{\σ}}_e\left(0\right)\right|{|}_{\mathrm{\∞}}+\left|\right|\mathrm{\γ}\mathrm{sgn}\left(S\right)|{|}_{\mathrm{\∞}}\≤l{\left|h\left(t\right)\right|}_{max}|\left|\hat{J}{\mathrm{\σ}}_e\left(0\right)\right|{|}_{\mathrm{\∞}}+\left|\mathrm{\γ}\right|---\left(37\right)$

式中

当参数a和k给定，则函数h(t)随着时间运行总||T_c′||_∞≤T_max会在一个确定的范围里，如果满足下式

$l{\left|h\left(t\right)\right|}_{max}\left|\right|\hat{J}{\mathrm{\σ}}_e\left(0\right)|{|}_{\mathrm{\∞}}+|\mathrm{\γ}|\≤T_{max}---\left(38\right)$

则满足饱和限制|T_c′||_∞≤T_max。因此，通过数学分析和计算机工具箱可以得到一个合适的γ的取值，使得饱和要求满足。γ取正值，则

$\mathrm{\γ}\≤T_{max}-l\left|h\left(t_1\right)\right|\left|\right|\hat{J}{\mathrm{\σ}}_e\left(0\right)|{|}_{\mathrm{\∞}}---\left(39\right)$

实施例

本发明在Matlab20011b环境下进行仿真验证。航天器的惯量阵为

$\hat{J}=\left[\begin{array}{ccc}3472& 0& 0\\ 0& 2280& 0\\ 0& 0& 2992\end{array}\right]kg\·m^2$

外界干扰力矩和惯量阵误差如下

$T_d=\left[\begin{array}{c}0.1sin\left(0.1t\right)\\ 0.1\sin (0.1t+\frac{\mathrm{\π}}{2})\\ 0.2\sin \left(0.1t\right)\end{array}\right],\mathrm{\Δ}J=0.2sin\left(t\right)\·\hat{J}$

执行器物理饱和限制为T_max＝1N。为了简单起见，这里只考虑前三阶振动模态，频率矩阵为Λ＝diag(0.31 0.83 1.33)rad/s，阻尼比矩阵为ζ＝diag(0.05 0.05 0.05)，刚柔耦合矩阵为

$C=\left[\begin{array}{ccc}-10.65& 0.14& -0.1\\ 0.05& -10.33& 0\\ 0.21& 0.1& -3.03\end{array}\right]{\mathrm{kg}}^{1/2}m$

姿态信息的初始值和期望值如下所示

σ(0)＝[0 0 0]rad，ω(0)＝[0 0 0]rad/s，η(0)＝[0 0 0]，

ω_d＝[0.2 0.2 -0.2]rad/s，ω_d＝[0 0 0]rad/s

为了满足物理饱和限制，通过计算可以得到γ≤0.88，在这里γ取值为0.88。滑模控制器的参数为k＝0.008，a＝0.5，边界层厚度ε＝0.0001。观测器增益矩阵P＝I₆，一阶惯性环节时间常数T_ci＝15。

为了体现本发明提出的控制方法的优越性，现对以下三种情况进行比较：比例为分(PD)控制系统响应；PD+IS控制系统响应；带有一阶惯性环节的输出反馈指数时变滑模控制(OFETVSMC)+FIS控制系统响应。

本发明带有一阶惯性环节的OFETVSMC+FIS控制律、PD控制律和PD+IS控制律的挠性模态坐标η响应曲线如图4所示。以第一阶模态坐标为例，三种控制律作用下的最大幅值分别为0.02m，0.12m和0.038m。可以看出，本发明带有一阶惯性环节的OFETVSMC+FIS控制律可以有效抑制挠性附件的残余振动。三种控制律的姿态角σ响应曲线如图3所示。以σ₁为例，三种控制律作用下的稳定误差分别为1×10^-4rad、0.02rad和0.02rad。仿真数据表明指数时变滑模控制作用下的稳定误差要比传统PD控制作用下的小很多，说明本发明不仅可以抑制挠性附件的残余振动，而且比传统控制方法具有很强的鲁棒性。

图3中给出的是PD+IS和OFETVSMC+FIS两种控制策略的姿态机动路径曲线图。从图中可以看出，两种控制都实现了系统的特征轴机动，姿态机动路径为一条沿着特征轴方向的直线，因此姿态机动路径最短。另一方面，在外界干扰和参数不确定性的影响下，PD+IS下的曲线偏离了特征轴方向，而OFETVSMC+FIS控制下的姿态机动路径曲线始终与特征轴方向重合，表明本发明相比PD控制具有更强的鲁棒性。

图6给出了三种控制律的控制力矩响应曲线图。使用IS控制的控制力矩最大幅值比不使用IS控制的要小很多，说明输入成型可以有效减小控制力矩的最大峰值，提高了控制力矩的响应性能。同时，在带有一阶惯性环节的OFETVSMC+FIS控制律作用下的控制力矩最大峰值都不超过1N，满足执行机构的物理饱和约束。控制力矩的饱和性分析优化了控制律参数，使得控制力矩既可以保证系统渐进稳定，又可以满足一定的饱和限制约束。

Claims (3)

1.一种挠性航天器特征轴姿态机动指数时变滑模控制方法，其特征在于：具体包括以下步骤：

步骤1，在航天器姿态运动的构型空间内定义相对姿态变量和挠性振动变量，在航天器本体系下建立挠性航天器动力学方程、姿态运动学方程和振动方程；

具体方法为：

挠性航天器动力学方程为：

振动方程为：

式中为系统实际的正定对称转动惯量矩阵，为系统名义惯量阵，ΔJ为由系统质量变化引起的惯量阵误差；ω＝[ω₁ ω₂ ω₃]^T为在本体系下的刚体角速度表示，ω^×为ω的斜对称矩阵；T_c＝[T_c1 T_c2 T_c3]^T和T_d＝[T_d1 T_d2 T_d3]^T分别表示控制力矩和外界干扰力矩；η为挠性结构弹性形变的广义坐标，ζ和Λ分别为挠性附件的模态阻尼矩阵和模态频率矩阵，C为挠性附件与星体的刚柔耦合矩阵；

姿态误差σ_e和角速度误差ω_e为：

式中σ为挠性航天器当前姿态角，σ_d为期望姿态角，σ_d ^×为σ_d的斜对称矩阵，为航天器在本体系下的期望角速度，ω_d为航天器在惯性系下的期望角速度，为从惯性系到本体系的转移矩阵，σ^×为σ的斜对称矩阵；

姿态运动学方程在本体系下表示为

式中满足条件M^TM＝m_aI_3×3，其中I_3×3为3×3的单位矩阵；

设φ＝[η^T ψ^T]^T为航天器挠性振动变量，其中令D＝2ζΛ，K＝Λ²；整理得到等效的数学模型：

式中E＝[K D]，表示由惯量阵不确定性和外界干扰引起的聚合扰动，I表示单位阵，d约束于未知上界d_max；

步骤2，针对步骤1建立的等效数学模型，设计状态反馈指数时变滑模控制律，具体方法为：

设计滑模面函数为

式中l和a为正的标量，Q＝-lσ_e(0)，σ_e(0)即σ_e在0时刻的状态值，e^-at以自然指数e为底数的指数函数，t为时间；

设计状态反馈指数时变滑模控制律如下

γ为滑模的切换增益，γ＞0,E＝[K D]；S即为滑模面，S＝S(ω_e,σ_e,t)，sgn(S)为符号函数，当S大于0时，sgn(S)为1，当S的等于0时，sgn(S)为0，当S小于0时，sgn(S)为-1；

步骤3，将系统三轴运动转变为绕特征轴旋转的单轴运动，计算出带有步骤2中指数时变滑模控制律的闭环系统振动频率和阻尼比参数，设计多模态滤波输入成型器；

多模态滤波输入成型器的具体设计方法为：

设状态变量为X＝[σ^T ω^T η^T ψ^T]^T，对闭环系统建立状态空间方程，求出 系数矩阵的特征值λ_sys，再通过下式求出闭环系统的振动频率和阻尼比参数：

式中ζ_sys,i和ω_sys,i分别为第i阶闭环系统的阻尼比和振动频率；

零振动输入成型器形式为

T₁＝0

式中A₁，A₂，A₃和T₁，T₂，T₃分别为脉冲的幅值和时间；

零振动输入成型器的脉冲序列通过各个单模态的脉冲序列相卷积得到，公式如下：

A_mult＝A_mult1*A_mult2*…A_multj…*A_multn

式中A_multj代表第j阶单模态的脉冲序列，*代表卷积运算，j＝1,2,…,n；

将设计好的零振动输入成型器与一个一阶惯性环节相连接，构成多模态滤波输入成型器；

步骤4，设计挠性状态观测器如下

其中为航天器挠性振动变量φ的估计值，P为观测器正定对称增益矩阵，且PA＜0；则输出反馈指数时变滑模控制律如下

将输出反馈指数时变滑模控制律输入航天器模型，航天器模型在该控制律 控制下得到σ，将σ输入挠性状态观测器，然后将在挠性状态观测器下得到的输入控制器，从而使航天器模型、挠性状态观测器、控制器形成一个闭环系统。

2.根据权利要求1所述的一种挠性航天器特征轴姿态机动指数时变滑模控制方法，其特征在于：所述一阶惯性环节的形式如下

式中T_ci为惯性时间常量。

3.根据权利要求1所述的一种挠性航天器特征轴姿态机动指数时变滑模控制方法，其特征在于：所述闭环系统包括控制器、挠性状态观测器和航天器模型。