CN109782787B - 一种太阳光压辅助下欠驱动航天器姿态的双模mpc控制方法 - Google Patents

Abstract

本发明公开了一种太阳光压辅助下欠驱动航天器姿态的双模MPC控制方法，针对高轨卫星及深空探测器等采用飞轮作为执行机构、受到的环境力矩主要是太阳光压力矩的挠性欠驱动航天器姿态控制系统。首先，建立以飞轮为执行机构的挠性航天器动力学与运动学模型，并根据航天器构型，建立太阳光压力矩解析式；然后，根据状态方程建立离散预测模型并确定约束条件；最后，根据近似的线性模型，确定终端惩罚矩阵和终端域，设计双模MPC控制律。本发明考虑了帆板转角、转速的限制，以及帆板的刚性转动和振动对航天器姿态的干扰，通过控制帆板转角产生太阳光压力矩作为辅助力矩，实现了欠驱动挠性航天器三轴姿态稳定。

Description

一种太阳光压辅助下欠驱动航天器姿态的双模MPC控制方法

技术领域

本发明涉及一种太阳光压辅助下欠驱动航天器姿态的双模MPC控制方法，实现了航天器在只有两个飞轮的情况下，将太阳光压力矩作为辅助力矩，利用双模MPC控制算法完成三轴姿态稳定。可用于高轨及深空探测等欠驱动航天器的姿态稳定控制中。

背景技术

欠驱动航天器是指执行机构提供的独立控制力矩个数少于系统运动自由度的航天器。对于微小卫星和深空探测器等需长时间工作的卫星，研究欠驱动航天器的姿态动力学与控制问题，能够有效提高其姿态控制系统的可靠性，延长航天器的工作寿命。同时能够减少执行机构配置，减小航天器的质量和功耗。

目前已有的对使用角动量交换装置的欠驱动航天器的研究都是在忽略环境干扰力矩，系统总角动量守恒的条件下得到的。而航天器在实际运行中，环境力矩是不可忽略的，这导致目前已有的欠驱动控制研究成果很难应用于工程实际。环境力矩在全驱动航天器姿态控制系统设计中通常作为干扰力矩处理，而对于欠驱动航天器，若作为干扰力矩，势必会进一步增大欠驱动控制系统设计的难度。本发明借鉴早期的航天器常利用环境力矩做被动/半被动稳定控制的思想，在欠驱动航天器的控制问题研究中，将环境干扰力矩作为辅助力矩联合控制姿态，以有效改善姿态可控性并提升控制性能。

高轨卫星以及深空探测器等航天器常采用飞轮作为姿态控制执行机构，而长时间的在轨运行，飞轮有可能失效从而导致航天器成为欠驱动航天器。而此类航天器其主要受到的主要环境力矩是太阳光压力矩，因此本发明针对使用飞轮的欠驱动航天器，提出一种太阳光压力矩辅助下的联合姿态控制方法。

发明内容

本发明要解决的技术问题是：针对高轨卫星以及深空探测器等采用飞轮作为执行机构、受到的环境力矩主要是太阳光压力矩的挠性欠驱动航天器姿态稳定系统，提供一种利用太阳光压力矩辅助完成姿态控制的方法，提出将基于整体数学模型的非线性模型预测控制(Nonlinear Model Predictive Control，NMPC)用于太阳光压辅助的欠驱动航天器姿态控制。实现了航天器在只有两轴姿态控制力矩输出能力的情况下，完成三轴姿态控制的目的，可用于高轨及深空探测等欠驱动航天器的姿态稳定控制。

本发明解决上述技术问题采用的技术方案为：针对高轨卫星以及深空探测器等采用飞轮作为执行机构、受到的环境力矩主要是太阳光压力矩的挠性欠驱动航天器姿态稳定系统；首先，建立以飞轮为执行机构的挠性航天器动力学与运动学模型，并根据航天器构型，建立太阳光压力矩解析式；然后，根据状态方程建立离散预测模型并确定约束条件；最后，根据近似的线性模型，确定终端惩罚矩阵和终端域，设计双模MPC控制律；具体实施步骤如下：

(1)建立以飞轮为执行机构的挠性航天器动力学与运动学模型；

建立以飞轮为执行机构的挠性航天器动力学与运动学模型：

基于欧拉角描述的航天器运动学方程可以写为：

其中，ω_b＝[ω_x ω_y ω_z]^T是本体坐标系相对惯性系的角速度在本体系下的分量列阵，

为欧拉角角速度列阵，

θ和ψ分别代表航天器的滚动角、俯仰角以及偏航角，ω₀是轨道角速度；

假定航天器的姿态角、姿态角速度都为小量，则运动学方程式可简化为：

假设航天器由中心本体和两个帆板组成，中心刚体和帆板均为分布均匀的六面体，且帆板具有一个转动自由度，中心刚体可视为质量均匀分布的六面体，飞轮安装构型为三正交加一斜装。建立以飞轮为执行机构的挠性航天器的动力学模型：

其中，I_t为加入执行机构后整个系统的转动惯量矩阵，C为飞轮安装矩阵，I_w为飞轮组转动惯量；Ω为飞轮转速列向量，

表示飞轮组提供的力矩，T_srp表示太阳光压力矩。Λ_ak为帆板的模态频率对角阵，ξ_ak为帆板的模态阻尼矩阵，I_ak为帆板的转动惯量在本体系下的分量列阵，T_ak为作用在太阳帆板上的外力矩，R_bak为帆板转动对中心刚体转动的刚性耦合系数矩阵，η_ak为归一化后的模态坐标，F_bak为帆板振动对中心刚体转动的柔性耦合系数矩阵，F_ak为帆板振动对自身转动的柔性耦合系数矩阵。

在控制器设计中将这些帆板挠性振动和刚性转动视为内部扰动项，并且通过施加约束可以将扰动项对中心刚体的影响降低，因此可以得到如下的刚体转动动力学方程：

其中，h为飞轮角动量。

(2)根据给出的航天器构型，建立太阳光压力矩解析式；

根据给定航天器构型，建立航天器受到的太阳光压力矩T_srp解析式：

对于带帆板的航天器，根据几何遮挡算法对中心本体和太阳帆板所受的太阳光压力矩分别进行分析。得到太阳帆板的光压力矩在本体轴上的分量为：

其中，β_i(i＝1,2)表示帆板绕安装轴Y_b轴旋转的角度。根据上式可知，当β₁＝β₂时，T_{s_srp}＝0，当β₁≠β₂，安装在本体系Y_b轴的帆板可以产生X_b和Z_b轴上的力矩。

(3)根据状态方程建立离散预测模型并确定约束条件；

根据动力学与运动学模型，建立仿射非线性系统状态方程：

其中，状态向量选为

控制向量

系统方程中的非线性项及控制系数矩阵表达式如下：

其中，

λ₂＝0，

根据系统状态方程，应用四阶龙格库塔法，得到系统离散的预测模型，其预测方程可以写成：

根据实际要求确定状态量以及控制量的约束：由于姿态角通过3-1-2转序的欧拉角描述，为避免奇异，欧拉角的范围约束为

θ∈(-π,π]，ψ∈(-π,π]；因为航天器姿态均处于小角速度范围内，取

由于帆板本身具有对称性，为在保证光压力矩的大小的同时，不至于造成帆板翻滚的情况，将帆板转角的约束为：

帆板的转速限制为

飞轮角加速度的约束为

综上，系统的状态量约束和控制量约束写为：

x(k)∈X,X＝{x∈Rⁿ|x∈[x_min,x_max]},k＞0

u(k)∈U,U＝{u∈R^m|u∈[u_min,u_max]},k＞0

(4)根据近似的线性模型，确定终端惩罚矩阵和终端域，设计双模MPC控制律；

假设仿射非线性系统存在局部稳定的线性控制器，计算系统方程在平衡点(0,0)处的雅克比矩阵，得到

然后采用如下算法离线计算终端惩罚矩阵和终端域：

a.基于雅克比线性化模型

设计线性反馈增益K使得A+BK渐进稳定的；

b.令Q^*＝Q+K^TRK∈R^n×n，求解离散Lyapunov方程

κ＞1为常数，得到唯一正定对称解P；

c.寻找尽可能大的α₁＞0，使得对于所有的x∈Ω₁，有Ω₁∈X，Kx∈U成立，其中Ω₁＝{x∈Rⁿ|x^TPx≤α₁}；

d.定义V(x(k))＝x^T(k)Px(k)，寻找尽可能大的α∈(0,α₁]，使得在Ω内局部满足Hamilton-Jacobian-Bellman(HJB)不等式：V(x(k+1))-V(x(k))≤-x^T(k)Q^*x(k)；

得到终端惩罚矩阵和终端域后，可以求解满足动力学约束和时域约束的每个采样时刻k的优化问题：

其中，N为预测时域，Q∈R^n×n和R∈R^m×m是正定对称加权矩阵，P∈R^n×n为终端惩罚矩阵。求解优化问题，并将控制序列的第一项u^*(k|k)作为终端域外的控制量，因此双模MPC控制律可以写为：

本发明与现有技术相比的优点在于：传统欠驱动航天器控制大多是利用两个飞轮设计控制律实现三轴姿态稳定，控制律较为复杂，较难实现。本发明将航天器受到的扰动力矩(太阳光压力矩)作为辅助控制力矩，并建立了光压力矩的解析模型，通过控制帆板转角实现航天器三轴姿态稳定，较容易实现，可用于高轨及深空探测等欠驱动航天器的姿态稳定控制。

附图说明

图1为一种太阳光压辅助下欠驱动航天器姿态的双模MPC控制方法的流程框图；

图2为航天器在本体系下的示意图。航天器中心本体尺寸为50×25×20cm³，两个太阳帆板的几何尺寸为80×25×1cm³，两个太阳帆板安装点在中心体机械坐标系下的坐标分别为(0,±21,0)cm。

具体实施方式

下面结合附图及实施例对本发明进一步详细说明。

如图1所示，本发明一种太阳光压辅助下欠驱动航天器姿态的双模MPC控制方法，具体实现步骤如下：

第一步，建立以飞轮为执行机构的挠性航天器动力学与运动学模型，并根据给出的航天器构型，建立太阳光压力矩解析式；

建立以飞轮为执行机构的挠性航天器动力学与运动学模型：

基于欧拉角描述的航天器运动学方程可以写为：

其中，ω_b＝[ω_x ω_y ω_z]^T是本体坐标系相对惯性系的角速度在本体系下的分量列阵，

为欧拉角角速度列阵，

θ和ψ分别代表航天器的滚动角、俯仰角以及偏航角，ω₀是轨道角速度；

假定航天器的姿态角、姿态角速度都为小量，则运动学方程式可简化为：

假设航天器由中心本体和两个帆板组成，中心刚体和帆板均为分布均匀的六面体，且帆板具有一个转动自由度，中心刚体可视为质量均匀分布的六面体，飞轮安装构型为三正交加一斜装。建立以飞轮为执行机构的挠性航天器的动力学模型：

其中，I_t为加入执行机构后整个系统的转动惯量矩阵，C为飞轮安装矩阵，I_w为飞轮组转动惯量；Ω为飞轮转速列向量，

表示飞轮组提供的力矩，T_srp表示太阳光压力矩。Λ_ak为帆板的模态频率对角阵，ξ_ak为帆板的模态阻尼矩阵，I_ak为帆板的转动惯量在本体系下的分量列阵，T_ak为作用在太阳帆板上的外力矩，R_bak为帆板转动对中心刚体转动的刚性耦合系数矩阵，η_ak为归一化后的模态坐标，F_bak为帆板振动对中心刚体转动的柔性耦合系数矩阵，F_ak为帆板振动对自身转动的柔性耦合系数矩阵。

在控制器设计中将这些帆板挠性振动和刚性转动视为内部扰动项，并且通过施加约束可以将扰动项对中心刚体的影响降低，因此可以得到如下的刚体转动动力学方程：

其中，h为飞轮角动量。

假设航天器位于距离太阳一个天文单位的日心轨道上，太阳方向矢量在该轨道系中始终为s_o＝L_bo[0 0 1]^T，航天器中心刚体转动惯量取I_b＝diag[0.02,0.058 0.0625]kggm²；执行机构为三正交加一斜装的飞轮组，初始转速为Ω＝[0 0 0 0]^Trad/s，假设只有安装在本体系Y轴和Z轴的飞轮正常工作，故障的飞轮转速视为0。

中心刚体对称安装两个太阳帆板，其相对系统质心的惯量矩阵分别为：

则考虑执行机构以及太阳帆板的航天器整体转动惯量矩阵I_t的计算公式为：

I_t＝I_b+I_bal+I_bar+CI_wC^T

其中，C表示飞轮组的安装矩阵：

第二步，根据图2给出的航天器构型，建立太阳光压力矩T_srp解析式；

对于带帆板的航天器，根据几何遮挡算法对中心本体和太阳帆板所受的太阳光压力矩分别进行分析。得到太阳帆板的光压力矩在本体轴上的分量为：

其中，β_i(i＝1,2)表示帆板绕安装轴Y_b轴旋转的角度。根据上式可知，当β₁＝β₂时，T_{s_srp}＝0，当β₁≠β₂，安装在本体系Y_b轴的帆板可以产生X_b和Z_b轴上的力矩。

由于在控制器设计中只考虑到了帆板的光压力矩，且帆板上下两个帆面为主要的光压力矩来源。因此帆板主要受照表面面积A＝0.2m²，且该表面中心距离航天器质心的矢量长度为r＝0.61m，帆板表面材料特性参数分别为ρ_a＝0.75，ρ_d＝0，ρ_s＝0.25。太阳光压系数P＝4.5598×10^-6N/m²，干扰力矩T_d的表达式为：

航天器初始的姿态参数为θ₀＝[-8 10 9]°，ω_b0＝[1.5 1.5 1.5]°/s。

第三步，根据状态方程建立离散预测模型并确定约束条件；

根据动力学与运动学模型，建立仿射非线性系统状态方程：

其中，状态向量选为

控制向量

系统方程中的非线性项及控制系数矩阵表达式如下：

其中，

根据系统状态方程，应用四阶龙格库塔法，得到系统离散的预测模型，其预测方程可以写成：

根据实际要求确定状态量以及控制量的约束：由于姿态角通过3-1-2转序的欧拉角描述，为避免奇异，欧拉角的范围约束为

θ∈(-π,π]，ψ∈(-π,π]；因为航天器姿态均处于小角速度范围内，取

由于帆板本身具有对称性，为在保证光压力矩的大小的同时，不至于造成帆板翻滚的情况，将帆板转角的约束为：

帆板的转速限制为

飞轮角加速度的约束为

综上，系统的状态量约束和控制量约束写为：

x(k)∈X,X＝{x∈Rⁿ|x∈[x_min,x_max]},k＞0

u(k)∈U,U＝{u∈R^m|u∈[u_min,u_max]},k＞0

第四步，根据近似的线性模型，确定终端惩罚矩阵和终端域，设计双模MPC确定的控制律。

假设仿射非线性系统存在局部稳定的线性控制器，计算系统方程在平衡点(0,0)处的雅克比矩阵，得到

然后采用如下算法离线计算终端惩罚矩阵和终端域。

a.基于雅克比线性化模型

设计线性反馈增益K使得A+BK渐进稳定的；

b.令Q^*＝Q+K^TRK∈R^n×n，求解离散Lyapunov方程

κ＞1为常数，得到唯一正定对称解P；

c.寻找尽可能大的α₁＞0，使得对于所有的x∈Ω₁，有Ω₁∈X，Kx∈U成立，其中Ω₁＝{x∈Rⁿ|x^TPx≤α₁}；

d.定义V(x(k))＝x^T(k)Px(k)，寻找尽可能大的α∈(0,α₁]，使得在Ω内局部满足Hamilton-Jacobian-Bellman(HJB)不等式：V(x(k+1))-V(x(k))≤-x^T(k)Q^*x(k)；

得到终端惩罚矩阵和终端域后，可以求解满足动力学约束和时域约束的每个采样时刻k的优化问题：

其中，N为预测时域，Q∈R^n×n和R∈R^m×m是正定对称加权矩阵，P∈R^n×n为终端惩罚矩阵。求解优化问题，并将控制序列的第一项u^*(k|k)作为终端域外的控制量，因此双模MPC控制律可以写为：

目标函数的加权矩阵Q和R取：Q＝diag(1 1 1 1×10⁴ 1 1 1×10^-6 1×10^-6)，R＝diag(0.001 0.001 0.001 0.001)。首先对仿射非线性模型进行雅克比线性化，求解线性二次型最优控制问题，得到线性的状态反馈控制增益u＝Kx；再根据第三步计算Lyapunov方程求得唯一正定对称解P。将该矩阵作为性能指标函数的终端惩罚矩阵；然后计算满足控制约束和状态约束的区域Ω₁＝{x∈Rⁿ|x^TPx≤α₁}，得α₁＝6.8859；最后找到一个Ω＜Ω₁使得d中的不等式成立，通过从α＝α₁开始不断减小α直到不等式成立，我们得到了一个不太保守的终端域Ω＝{x∈Rⁿ|x^TPx≤7.5×10^-6}。

系统仿真步长为0.1s，预测步长为1s，预测时域N＝5，控制时域M＝2。

本发明说明书中未作详细描述的内容属于本领域专业技术人员公知的现有技术。

Claims (1)

1.一种太阳光压辅助下欠驱动航天器姿态的双模MPC控制方法，其特征在于：包括以下步骤：

第一步，建立以飞轮为执行机构的挠性航天器动力学与运动学模型，并根据航天器构型，建立太阳光压力矩解析式；

第二步，根据状态方程建立离散预测模型并确定约束条件；

第三步，根据近似的线性模型，确定终端惩罚矩阵和终端域，设计双模MPC控制律；

其中，所述第一步具体实现如下：

建立以飞轮为执行机构的挠性航天器动力学与运动学模型：

基于欧拉角描述的航天器运动学方程可以写为：

其中，ω_b＝[ω_bx ω_by ω_bz]^T是本体坐标系相对惯性系的角速度在本体系下的分量列阵，

为欧拉角角速度列阵，

θ和ψ分别代表航天器的滚动角、俯仰角以及偏航角，ω₀是轨道角速度；

假定航天器的姿态角、姿态角速度都为小量，则运动学方程式可简化为：

假设航天器由中心本体和两个帆板组成，中心刚体和帆板均为分布均匀的六面体，且帆板具有一个转动自由度，中心刚体可视为质量均匀分布的六面体，飞轮安装构型为三正交加一斜装，建立以飞轮为执行机构的挠性航天器的动力学模型：

其中，I_t为加入执行机构后整个系统的转动惯量矩阵，C为飞轮安装矩阵，I_w为飞轮组转动惯量；Ω为飞轮转速列向量，

表示飞轮组提供的力矩，T_srp表示太阳光压力矩，Λ_ak为帆板的模态频率对角阵，ξ_ak为帆板的模态阻尼矩阵，I_ak为帆板的转动惯量在本体系下的分量列阵，T_ak为作用在太阳帆板上的外力矩，R_bak为帆板转动对中心刚体转动的刚性耦合系数矩阵，η_ak为归一化后的模态坐标，F_bak为帆板振动对中心刚体转动的柔性耦合系数矩阵，F_ak为帆板振动对自身转动的柔性耦合系数矩阵；

在控制器设计中将这些帆板挠性振动和刚性转动视为内部扰动项，并且通过施加约束可以将扰动项对中心刚体的影响降低，因此可以得到如下的刚体转动动力学方程：

其中，h为飞轮角动量；

根据航天器构型，建立航天器受到的太阳光压力矩T_srp解析式：

对于带帆板的航天器，根据几何遮挡算法对中心本体和太阳帆板所受的太阳光压力矩分别进行分析，得到太阳帆板的光压力矩在本体轴上的分量为：

其中，β_j表示帆板绕安装轴Y_b轴旋转的角度，其中j＝1,2，根据上式可知，当β₁＝β₂时，T_{s_srp}＝0，当β₁≠β₂，安装在本体系Y_b轴的帆板可以产生X_b和Z_b轴上的力矩；

所述第二步具体实现如下：

根据动力学与运动学模型，建立仿射非线性系统状态方程：

其中，状态向量选为

控制向量

系统方程中的非线性项及控制系数矩阵表达式如下：

其中，

λ₂＝0，

根据系统状态方程，应用四阶龙格库塔法，得到系统离散的预测模型，其预测方程可以写成：

根据实际要求确定状态量以及控制量的约束：由于姿态角通过3-1-2转序的欧拉角描述，为避免奇异，欧拉角的范围约束为

θ∈(-π,π]，ψ∈(-π,π]；因为航天器姿态均处于小角速度范围内，取

由于帆板本身具有对称性，为在保证光压力矩的大小的同时，不至于造成帆板翻滚的情况，将帆板转角的约束设为：

帆板的转速限制为

飞轮角加速度的约束为

综上，系统的状态量约束和控制量约束写为：

x(k)∈X,X＝{x∈Rⁿ|x∈[x_min,x_max]},k＞0

u(k)∈U,U＝{u∈R^m|u∈[u_min,u_max]},k＞0

所述第三步具体实现如下：

根据近似的线性模型，确定终端惩罚矩阵和终端域：

假设仿射非线性系统存在局部稳定的线性控制器，计算系统方程在平衡点(0,0)处的雅克比矩阵，得到

然后采用如下算法离线计算终端惩罚矩阵和终端域；

(1)基于雅克比线性化模型

设计线性反馈增益K使得A+BK渐进稳定的；

(2)令Q^*＝Q+K^TRK∈R^n×n，求解离散Lyapunov方程

κ＞1为常数，得到唯一正定对称解P；

(3)寻找尽可能大的α₁＞0，使得对于所有的x∈Ω₁₁，有Ω₁₁∈X，Kx∈U成立，其中Ω₁₁＝{x∈Rⁿ|x^TPx≤α₁}；

(4)定义V(x(k))＝x^T(k)Px(k)，寻找尽可能大的α∈(0,α₁]，使得在Ω内局部满足Hamilton-Jacobian-Bellman(HJB)不等式：V(x(k+1))-V(x(k))≤-x^T(k)Q^*x(k)；

得到终端惩罚矩阵和终端域后，可以求解每个采样时刻k的优化问题：

其中

其中，N为预测时域，Q∈R^n×n和R∈R^m×m是正定对称加权矩阵，P∈R^n×n为终端惩罚矩阵，并将控制序列的第一项u^*(k|k)作为终端域外的控制量，因此双模MPC控制律可以写为：

Images