CN108891625B - 固体微推进器阵列与磁力矩器联合控制方法 - Google Patents

固体微推进器阵列与磁力矩器联合控制方法 Download PDF

Info

Publication number
CN108891625B
CN108891625B CN201810585104.1A CN201810585104A CN108891625B CN 108891625 B CN108891625 B CN 108891625B CN 201810585104 A CN201810585104 A CN 201810585104A CN 108891625 B CN108891625 B CN 108891625B
Authority
CN
China
Prior art keywords
array
control
satellite
solid micro
thruster
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810585104.1A
Other languages
English (en)
Other versions
CN108891625A (zh
Inventor
侯晓磊
周康博
刘勇
潘泉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Northwestern Polytechnical University
Original Assignee
Northwestern Polytechnical University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Northwestern Polytechnical University filed Critical Northwestern Polytechnical University
Priority to CN201810585104.1A priority Critical patent/CN108891625B/zh
Publication of CN108891625A publication Critical patent/CN108891625A/zh
Application granted granted Critical
Publication of CN108891625B publication Critical patent/CN108891625B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B64—AIRCRAFT; AVIATION; COSMONAUTICS
    • B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
    • B64G1/00—Cosmonautic vehicles
    • B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
    • B64G1/24—Guiding or controlling apparatus, e.g. for attitude control
    • B64G1/244—Spacecraft control systems
    • B64G1/245—Attitude control algorithms for spacecraft attitude control

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Chemical & Material Sciences (AREA)
  • Combustion & Propulsion (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Automation & Control Theory (AREA)
  • Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)

Abstract

本发明公开了一种固体微推进器阵列与磁力矩器联合控制方法，设计固体微推进器阵列优化点火模型，通过优化点火模型控制固体微推进器阵列对卫星进行大角度速度阻尼控制，直至卫星的角速度
Figure DDA0001689304650000011
通过磁力矩器和固体微推进器阵列对卫星进行小角度阻尼联合控制，直至卫星角速度小于等于
Figure DDA0001689304650000012
通过磁力矩器和固体微推进器阵列对卫星进行联合姿态捕获控制；通过设计MSPT点火模型设计，可以优化固体微推进器的点火，给出最优的点火组合，根据任务需求，将固体微推进器阵列和磁力矩器相结合进行联合控制，并分别通过理论推导和仿真实验证明了算法的有效性，可以降低卫星的姿态控制周期，提升控制精度，减少控制能耗。

Description

固体微推进器阵列与磁力矩器联合控制方法
【技术领域】
本发明属于微星姿态控制技术领域，尤其涉及固体微推进器阵列与磁力矩器联合控制方法。
【背景技术】
在单位企业里，微小卫星以其实现成本低，设计时间周期短等的快速响应能力应用在应急通信、空间组网等领域。在高校教育中，其低成本、高集成的特性使得皮纳卫星可以很好的融入创新与学习，对学生的发展很有帮助，所以被很多大学所重视。当下已有多所大学的小卫星发射成功。皮纳卫星的姿态确定与控制系统(Attitude Determina tion andControlSystem简写为ADCS)为小卫星提供姿态控制能力，它由姿态敏感器、姿态控制执行器与姿控计算机组成。姿态敏感器配合姿控计算机可以提供卫星姿态信息，通过姿控计算机对陀螺仪、卫星导航接收机(GNSS)和磁强计、太阳敏感器、星敏感器等姿态敏感器采集的信息进行处理，解算出当前卫星的姿态，为姿态控制提供保障。
微型固体微推进器阵列的控制精度取决单个推进器的推力大小，而整体推进能力取决于于推进器的推力及阵列规模。当单个推进器推力较大时，姿态控制精度将无法保证；而推进器推力较小时，在完成大角度机动或多次机动任务中，需要的推进器规模庞大，难以达到功耗、体积等方面的限制。因此单独采用固体微推进器阵列进行姿态控制无法满足如很多实际飞行任务的需求，为充分利用微型固体微推进器的优势，仍需探究更好的组合姿态控制方法。
【发明内容】
本发明的目的是提供一种固体微推进器阵列与磁力矩器联合控制方法，以磁力矩器作为主执行器，同时为充分发挥微型固体微推进器阵列的作用，提高控制效率，降低对卫星的姿态控制周期、提升控制精度、减少控制能耗。
本发明采用以下技术方案：固体微推进器阵列与磁力矩器联合控制方法，具体包括以下步骤：
步骤1、设计固体微推进器阵列优化点火模型，通过优化点火模型控制固体微推进器阵列对卫星进行大角度速度阻尼控制，直至卫星的角速度
Figure BDA0001689304630000021
步骤2、基于控制律
Figure BDA0001689304630000022
通过磁力矩器和固体微推进器阵列对卫星进行小角度阻尼联合控制，直至卫星角速度小于等于
Figure BDA0001689304630000023
其中，
Figure BDA0001689304630000024
mb为磁力矩器产生的输出力矩，K1、K2均为正定义系数矩阵，
Figure BDA0001689304630000025
Bb表示本体系地磁场矢量，
Figure BDA0001689304630000026
表示卫星在本体系内的真实角速度，tk表示固体微推进器阵列的工作时间点，k∈{1,2,…,M}，M为大于等于1的整数；Ak为固体微推进器阵列在小角度阻尼控制过程中每次工作产生的输出冲量；
步骤3、基于控制律
Figure BDA0001689304630000027
通过磁力矩器和固体微推进器阵列对卫星进行联合姿态捕获控制；
其中，uc(t)为磁力矩器在捕获控制时的输出力矩，udk为固体微推进器阵列在捕获控制时产生的输出冲量，kc为连续控制律增益，
Figure BDA0001689304630000028
表示本体系磁场的斜对称矩阵、J为转动惯量矩阵，P1、P2均为正定矩阵，θ(t)表示卫星在t时刻欧拉角，
Figure BDA0001689304630000031
的导数，kd为离散控制律增益，
Figure BDA0001689304630000032
为卫星在
Figure BDA00016893046300000319
时刻的欧拉角，
Figure BDA0001689304630000033
为
Figure BDA0001689304630000034
的导数。
进一步地，固体微推进器阵列优化点火模型具体为：
Figure BDA00016893046300000318
Figure BDA0001689304630000035
其中，Sk为固体微推进器阵列第k次点火时的目标函数，
Figure BDA0001689304630000036
分别表示向量
Figure BDA0001689304630000037
和
Figure BDA0001689304630000038
中元素的个数，nk为固体微推进器阵列在第k次点火的倍数参数，
Figure BDA0001689304630000039
表示固体微推进器阵列B阵列中第i行第j列的固体微推进器状态，
Figure BDA00016893046300000310
为固体微推进器阵列D阵列中第i行第j列的固体微推进器状态，B阵列
Figure BDA00016893046300000311
和D阵列
Figure BDA00016893046300000312
表示卫星不同的两个面的微推进器阵列；
Figure BDA00016893046300000313
分别为B阵列的状态矩阵和D阵列的状态矩阵中第i行第j列元素。
进一步地，步骤2中控制律设计方法为：
选取卫星相对于惯性系的旋转能量作为Lypunov函数V(x)：
Figure BDA00016893046300000314
V(x)相对于时间求导，得出：
Figure BDA00016893046300000315
其中，
Figure BDA00016893046300000316
为
Figure BDA00016893046300000317
的导数，Tb表示作用于卫星的总力矩，
则有：
Figure BDA0001689304630000041
其中，δ(t-tk)表示t-tk时刻的脉冲，teq表示燃烧周期，Atot为固体微推进器阵列在大角度速度阻尼控制阶段的总输出冲量；
令
Figure BDA0001689304630000042
且有正定义系数矩阵 K1＝diag{[k1x k1yk1z]}和K2＝diag{[k2x k2y k2z]}，则有：
Figure BDA0001689304630000043
其中，S为斜对称矩阵，
Figure BDA0001689304630000044
根据Lypunov稳定性定理，系统渐进稳定，磁场矢量在轨道平面内变化率接近于0，则有：
Figure BDA0001689304630000045
得出对卫星进行小角度阻尼控制的控制律：
Figure BDA0001689304630000046
进一步地，步骤3中控制律设计方法：
初始状态量x(t0)＝x0，得到在t≠tk的连续函数表达式为：
Figure BDA0001689304630000047
其中，Αc(t)、Bc(t)、Cc(t)、Dc(t)表示连续状态空间函数，
Figure BDA0001689304630000048
为x(t)的导数， yc(t)表示输出矢量，uc(t)表示磁力矩器在卫星捕获控制阶段的输出力矩，
在t＝tk的离散函数表达式为：
Figure BDA0001689304630000051
其中，
Figure BDA0001689304630000052
表示tk +时刻状态向量的导数，Αdk、Bdk、Cdk、Ddk分别表示离散状态空间函数，
Figure BDA0001689304630000053
表示tk -时刻状态向量,ydk表示状态矢量；
根据姿态运动学和动力学方程，得出混合系统的线性模型：
Figure BDA0001689304630000054
其中，
Figure BDA00016893046300000511
为θ(t)的一阶导数，
Figure BDA0001689304630000056
为θ(t)的二阶导数；I3×3表示3*3的单位矩阵；m(t)表示输出磁场；Αd、Bd分别表示矩阵
Figure BDA0001689304630000057
和
Figure BDA0001689304630000058
Αc表示连续状态空间函数；
Figure BDA0001689304630000059
表示磁场的斜对称矩阵；
则系统输入表示为：
Figure BDA00016893046300000510
本发明的有益效果是：通过设计MSPT点火模型设计，可以优化固体微推进器的点火，给出最优的点火组合，根据任务需求，将固体微推进器阵列和磁力矩器相结合进行联合控制，并分别通过理论推导和仿真实验证明了算法的有效性，可以降低卫星的姿态控制周期，提升控制精度，减少控制能耗。
【附图说明】
图1为本发明实施例中微分轨道元系统的能控图；
图2为本发明实施例中圆形近极点轨道的一轨的可控性函数的最小特征值示意图；
图3为本发明实施例中的角速度曲线图；
图4为本发明实施例中的输出磁矩曲线图；
图5为本发明实施例中小角度速度阻尼的角速度曲线图；
图6为本发明实施例中小角度速度阻尼的输出磁矩曲线图；
图7为本发明实施例中各种控制算法的角速度曲线图；
图8为本发明实施例中各种控制算法的姿态角曲线图；
图9为本发明实施例中各种控制算法的输出磁矩曲线图。
【具体实施方式】
下面结合附图和具体实施方式对本发明进行详细说明。
本发明公开了固体微推进器阵列与磁力矩器联合控制方法，具体包括以下步骤：
步骤1、设计MSPT(固体微推进器阵列)优化点火模型。
针对消旋或三轴稳定控制的皮纳卫星，研究固体脉冲式推进器阵列构成的直接力姿控系统在给定位置产生制定力矩的组合点火算法。这个问题最终可以转化为以所点火的推进器数量为目标，同时满足响应的力矩约束条件的最优组合选择问题。求解该组合选择问题的一种非常合适的方法是0-1规划问题。
MSPT阵列点火算法主要应用于纯MSPT姿态控制，提出了微推进器最少前提下，越靠近阵列中心的推进器优先级越高，以此保证剩余推进器提供更多的推进器组合。这种方法主要考虑思想前提是，靠近阵列对称中心的推进器可以产生能够产生的较少的力矩组合，但是在考虑微推进器推力方向的前提下，每行中的推进器产生的组合数是完全一致的，设单行中推进器数量为x，则该行中每个推进器可以产生的力矩组合数均为
Figure BDA0001689304630000071
此外，其提出的阵列规划方法并未全局考虑行选择对整体多样性的影响，采取逐行点火的方式，降低了固体微推进器的利用率。
在上述前提下，我们提出了MSPT优化点火模型设计：
a、使用的推进器组合数量最少；
b、同行不同列的微推进器的优先级相同；
c、同列不同行的微推进器的优先级取决于该行剩余微推进器数量，剩余数量越多，优先级越高。
下面以卫星X轴为例介绍MSPT优化点火算法。
定义{Xk,Yk}为第k次点火后控制X轴的两个推进器阵列的状态矩阵。其中， Xk为B阵列状态矩阵(控制卫星X轴的两个固体微推进器阵列分别为B阵列和 D阵列，其中B阵列状态矩阵指用N阶方阵表示的B阵列面点火状态，D阵列状态矩阵指用N阶方阵表示的D阵列面点火状态数值1表示未点火，数值0表示已点火)，Yk为D阵列状态矩阵。
Figure BDA0001689304630000072
为Xk中第i行的行向量，且有
Figure BDA0001689304630000073
Figure BDA0001689304630000074
为Xk中第i行第j列元素，i∈[1,N]，j∈[1,N]。其中：
Figure BDA0001689304630000075
定义X0为状态矩阵初始状态，式中所有元素均为1。同理，可以定义Yk中的元素为
Figure BDA0001689304630000076
第i行的行向量为
Figure BDA0001689304630000077
定义
Figure BDA0001689304630000078
为第k次根据控制律选择的微推进器集合，
Figure BDA0001689304630000079
对应B阵列，
Figure BDA00016893046300000710
对应D阵列。
Figure BDA0001689304630000081
为
Figure BDA0001689304630000082
中第i行的行向量，且有
Figure BDA0001689304630000083
Figure BDA0001689304630000084
为
Figure BDA0001689304630000085
中第i行第j列元素，其中：
Figure BDA0001689304630000086
同理，定义
Figure BDA0001689304630000087
中的元素为
Figure BDA0001689304630000088
第i行的行向量为
Figure BDA0001689304630000089
Figure BDA00016893046300000810
选择满足控制律设计，则：
Figure BDA00016893046300000811
根据上述定义，第k次点火后的状态矩阵表示为：
Figure BDA00016893046300000812
其中，Xk-1、Yk-1分别表示什么第k-1次点火后的状态矩阵，
第k次点火后第i行的行向量表示为：
Figure BDA00016893046300000813
定义Sk为第k次点火时的目标函数，满足上述约束条件，则MSPT优化点火模型目标函数为：
Figure BDA00016893046300000814
其中|x|表示向量x中元素的个数。
综上，MSPT的优化点火问题可以描述为：
Figure BDA0001689304630000091
Figure BDA0001689304630000092
其中，Sk为固体微推进器阵列第k次点火时的目标函数，
Figure BDA0001689304630000093
分别表示向量
Figure BDA0001689304630000094
和
Figure BDA0001689304630000095
中元素的个数，nk表示固体微推进器阵列在第k次点火的倍数参数，
Figure BDA0001689304630000096
表示B阵列
Figure BDA00016893046300000911
中第i行第j列的固体微推进器状态，
Figure BDA0001689304630000097
为D阵列
Figure BDA0001689304630000098
第i 行第j列的固体微推进器状态；B阵列和D阵列表示卫星不同的两个面的微推进器阵列；
Figure BDA0001689304630000099
分别为B阵列状态矩阵和D阵列状态矩阵中第i行第j列元素； B阵列状态矩阵和D阵列状态矩阵分别表示B阵列和D阵列点火情况的数学矩阵。
上述问题是典型的0/1规划算法，由于计算机性能要求，为简化计算，我们可以采用贪心算法进行模型求解，可得到具体的推进器组合。
固体微推进器阵列的执行时间问题：
可以根据磁力矩器控制能力的强弱来确定一组推进器执行的时间，即可以利用推进器来补偿磁力矩器的控制缺陷以提高指向精度。通过计算小时间窗(小时间窗是获得足够采用数据所需的特定一段时间)的宽度Δt跨越总旅行时间的可控性来评估系统的逐点的可控性。然后通过确定矩阵Wj的最小特征值到达每轨的局部最小值来选择脉冲时间。
Figure BDA00016893046300000910
式中：Wj为第j个时间窗的可控性，从
Figure BDA0001689304630000101
时刻开始，e表示自然底数，
Figure BDA0001689304630000102
为起始时刻，τ为积分变量。
Wj若是非奇异的，则卫星的控制系统是可控的。Wj的值主要根据轨道倾角以及轨道高度的不同产生变化。如图1所示，绘制了可控性函数(Wj)最小特征值的绝对值，并将其与第二小特征值进行了比较。对于不同的半长轴，轨道高度选择范围为6478km的LEO轨道到16378km的中轨，轨道倾角范围为0°到90°。结果表明，对于测试的整个轨道，在可控制性Gramian中始终具有一个接零的特征值，表明线性化微分动力学不能仅靠磁力矩器。
如图2所示，为轨道高度为682km，轨道倾角为97°的圆形近极点轨道的一轨的可控性函数的最小特征值，其中时间窗Δt＝0.05T。随着时间的推移，控制能力最差时刻约从
Figure BDA0001689304630000103
和
Figure BDA0001689304630000104
时刻开始，周期与轨道周期相同。这表明脉冲(指固体微推进器阵列的作用)应用于这些控制能力最差的时间窗口内。在此例中，脉冲应用于tk＝0.225T+(k-1)×.05T,k∈{1,2,…,N}。
通过所述MSPT优化点火模型结合固体微推进器阵列对卫星进行大角度速度阻尼控制，直至卫星的角速度
Figure BDA0001689304630000105
MSPT大角度速度阻尼目的在于快速将卫星的角速度降低到其最小精度以内。已知卫星初始时刻角速度为
Figure BDA0001689304630000106
目标角速度为
Figure BDA0001689304630000107
则需要机动的角速度大小为
Figure BDA0001689304630000108
此时，输出力矩全部由MSPT产生，根据角动量定理的积分形式， MSPT的总冲量输出为
Figure BDA0001689304630000109
倍数参数n满足
Figure BDA00016893046300001010
其中。 Amems为控制律求解得到的固体微推进器冲量输出理论值，Amin为最小值。
根据任务需求，当卫星需要进行大角速度阻尼时，由于MSPT最大可以产生 25°/s的速度增量，首先利用微推进器阵列大机动能力，使卫星速度状态迅速变化到既定阈值以内(即角速度
Figure BDA0001689304630000111
)，再利用磁力矩器结合固体微推进器阵列进行微调，在磁力矩器控制性能较低的点，利用微推进器进行补充，从而大幅度加快控制速度，缩短控制周期。
步骤2、基于控制律
Figure BDA0001689304630000112
通过磁力矩器和固体微推进器阵列对卫星进行小角度阻尼控制，直至卫星角速度≤
Figure BDA0001689304630000113
Figure BDA0001689304630000114
为预设阈值且
Figure BDA0001689304630000115
联合速度阻尼控制中，首先使用MSPT进行大角度脉冲，使得角速度稳定在 1°以内，然后采用联合控制策略。采用B-dot控制律驱动磁力矩器，同时，在磁力矩器控制能力较弱的时刻，采用MSPT控制，从而缩短控制周期。
其中，mb为磁力矩器产生的输出力矩，K1、K2均为正定义系数矩阵，
Figure BDA0001689304630000116
Bb本体系下的地磁场矢量；
Figure BDA0001689304630000117
表示卫星在本体系内的真实角速度， tk表示固体微推进器阵列的工作时间点，k∈{1,2,…,M}，M为大于等于1的整数； Ak为固体微推进器阵列在小角度阻尼控制过程中每次工作所产生的输出冲量。
控制律设计方法为：
假设MSPT作用时间为tk，k∈{1,2,…,M}。忽略其他干扰力矩，卫星的动力学方程为：
Figure BDA0001689304630000118
其中，J表示转动惯量；
Figure BDA0001689304630000119
表示卫星在本体系内的真实角速度；mb表示磁力矩器的输出力矩；Bb表示本体系下的地磁场矢量；δ(t-tk)表示t-tk时刻的脉冲； Ak表示固体微推进器阵列在大角度速度阻尼控制阶段的总输出冲量；teq表示燃烧周期。
用δ(t-tk)描述MSPT用以描述脉冲输出的特性。选取卫星相对于惯性系的旋转能量作为Lypunov函数V(x)：
Figure BDA0001689304630000121
其中，
Figure BDA0001689304630000122
为卫星在本体系内的真实角速度；J表示卫星转动惯量。
由Lypunov函数表达式可以看出，当
Figure BDA0001689304630000123
时，V(x)正定。忽略其他干扰力矩，根据能量公式，V(x)相对于时间求导，得出：
Figure BDA0001689304630000124
其中，
Figure BDA0001689304630000125
为
Figure BDA0001689304630000126
的导数；Tb表示作用于卫星的总力矩。
忽略其他干扰力矩，作用于卫星上的力矩全部由磁力矩器和推进器产生，则有：
Figure BDA0001689304630000127
其中，δ(t-tk)表示t-tk时刻的脉冲，Bb表示本体系地磁场矢量，Atot为固体微推进器阵列在大角度速度阻尼控制阶段的总输出冲量。
令
Figure BDA0001689304630000128
且有正定义系数矩阵
K1＝diag{[k1x k1y k1z]}和K2＝diag{[k2x k2y k2z]}，则：
Figure BDA0001689304630000129
其中，S斜对称矩阵，
Figure BDA0001689304630000131
当
Figure BDA0001689304630000132
时，
Figure BDA0001689304630000133
负定。根据Lypunov稳定性定理，系统渐进稳定，磁场矢量在轨道平面内变化率接近于0，则有：
Figure BDA0001689304630000134
综上，得出对卫星进行小角度阻尼控制：
Figure BDA0001689304630000135
倍数参数
Figure BDA0001689304630000136
步骤3、当卫星需要进行大角度姿态捕获时。由于微推进器阵列无法避免空间环境力矩干扰，且控制精度有限。因此，利用固体微推进器阵列作为辅助，在弱磁控的时刻，使用连续/脉冲结合的方案，采用PID控制的概念，增强控制性能，从而缩短姿态捕获周期，提高控制精度。
基于控制律
Figure BDA0001689304630000137
通过磁力矩器和固体微推进器阵列对卫星进行姿态捕获控制，其中，uc(t)为磁力矩器在捕获控制时的输出力矩， udk为固体微推进器阵列在捕获控制时产生的输出冲量，kc连续控制律增益，
Figure BDA0001689304630000138
J转动惯量矩阵，P1、P2均为正定矩阵，θ(t)为卫星在t时刻的欧拉角，
Figure BDA0001689304630000139
为θ(t) 的导数，kd离散控制律增益，
Figure BDA00016893046300001310
为卫星在
Figure BDA00016893046300001311
时刻的欧拉角，
Figure BDA00016893046300001312
为
Figure BDA00016893046300001313
的导数。
控制律具体通过以下方法得出：
对于一个连续/脉冲的混合系统，脉冲施加在tk，tk∈N时刻，可以表示为 y＝ζhu，其中，u＝{uc(t),udk}同时包含连续和离散的控制输入，同理，y＝{yc(t),ydk} 同时包含连续和离散系统的输出。假设uc和yc、udk和ydk分别具有相同的维数。
初始状态量x(t0)＝x0，得到在t≠tk的连续函数表达式为：
Figure BDA0001689304630000141
其中，Αc(t)、Bc(t)、Cc(t)、Dc(t)表示连续状态空间函数，Αc(t)和Αc含义相同，但意义略微不同，Αc(t)是Αc在t时刻的状态。
Figure BDA0001689304630000142
为x(t)的导数，yc(t)表示输出矢量，x(t)表示状态矢量，uc(t)表示磁力矩器在卫星捕获控制阶段的输出力矩。
在t＝tk的离散函数表达式为：
Figure BDA0001689304630000143
其中，
Figure BDA0001689304630000144
表示tk +时刻状态向量的导数，Αdk、Bdk、Cdk、Ddk分别表示离散状态空间函数，
Figure BDA0001689304630000145
表示tk -时刻状态向量,ydk表示状态矢量。
根据姿态运动学和动力学方程，得出混合系统的线性模型：
Figure BDA0001689304630000146
其中，
Figure BDA00016893046300001413
为θ(t)的一阶导数，
Figure BDA0001689304630000148
为θ(t)的二阶导数；I3×3表示3*3的单位矩阵；m(t)表示输出磁场；Αd、Bd分别表示矩阵
Figure BDA0001689304630000149
和
Figure BDA00016893046300001410
Αc表示连续状态空间函数；
Figure BDA00016893046300001411
表示磁场的斜对称矩阵。
此模型中，假设角速度和欧拉角均较小，因此θ≈2ε，
Figure BDA00016893046300001412
(ε、w表示两个任意小的常数。
系统输入表示为：
Figure BDA0001689304630000151
其中，uc(t)为磁力矩器在捕获控制时的输出力矩，udk为固体微推进器阵列在捕获控制时产生的输出冲量，kc为连续控制律增益，
Figure BDA0001689304630000152
表示本体系磁场的斜对称矩阵、J为转动惯量矩阵，P1、P2均为正定矩阵，θ(t)表示卫星在t时刻欧拉角，
Figure BDA0001689304630000153
为θ(t)的导数，kd为离散控制律增益，
Figure BDA0001689304630000154
为卫星在
Figure BDA0001689304630000155
时刻的欧拉角，
Figure BDA0001689304630000156
为
Figure BDA0001689304630000157
的导数。
令P＝[P1 P2]，将系统输入带入混合系统线性模型，得到：
Figure BDA0001689304630000158
为了稍微简化符号，令
Figure BDA0001689304630000159
由此得到连续/脉冲混合系统状态转移矩阵：
Figure BDA00016893046300001510
Figure BDA00016893046300001511
系统控制初始时间为
Figure BDA00016893046300001512
脉冲每轨道施加时间tk∈(0,T)，T为轨道周期， k∈{1,2,…,M}，且
Figure BDA00016893046300001513
由此可以将系统状态(一个周期的系统状态)分为k个区间，每个区间的大小为(tk-1,tk)。由于
Figure BDA00016893046300001514
且
Figure BDA00016893046300001515
则控制系统第一个轨道周期的状态可以描述为：
Figure BDA00016893046300001516
为了方便讨论，令
Figure BDA0001689304630000161
若矩阵Ψk对所有k＜M+1可逆，Ξ(t)表示系统的状态转移矩阵，即有 x(t)＝Ξ(t)x0。将Ξ(t)的定义扩展到多轨控制中，则Ξ(t+iT)＝Ξ(t)Ξi(T)，t∈[0,T)。式中，Αc和Bc是连续的，同样假设Bd是连续的。同时根据地磁场矢量的周期性变化， Bc(t)也是周期变化的。P是定常系数，若脉冲在每周期的施加时间一致，那么这个系统即可以成为混合线性周期系统，运用混合系统的Floquet方法进行稳定性证明，从而得到本系统是Lyapunov稳定的。
控制律校验：
为了验证本发明提出的联合控制律，根据初始角速度的大小，共设计了两组仿真对比实验。第一组仿真实验在分离角速度下，将本发明提出的联合控制方法与传统的纯磁控方法、纯MSPT控制方法以及本发明提出的MSPT辅助控制方法进行比较。第二组仿真实验主要对比小角速度情况下联合控制算法与传统的纯磁控方法及MSPT控制方法的控制效果。
不考虑其他干扰力矩的影响，选取初始姿态角为[60° 60° 60°]，磁控增益系数矩阵K＝diag(1.2×105,1.2×105,1.2×105)，MSPT控制律为 K2＝diag(0.91×10-2,1.08×10-2,0.62×10-2)。大角速度阻尼实验的初始角速度为 [5 5 5]°/s，小角速度阻尼实验的初始角速度为[1 1 1]°/s。仿真测试的结果如图 3、图4所示。
其中，MSPT与磁力矩器联合控制方法的稳定周期为1602s，1752s和1643s， 4000s后控制精度可以达到5×10-3°/s；0.017°/s和0.013°/s。X轴磁矩输出在22s 和135s左右时达到峰值，Y轴在60s时达到峰值，Z轴峰值时间为115s，平均幅值低于0.05Am2，三轴磁矩输出绝对值分别为6.72Am2，8.40Am2和8.64Am2。 Am2。表1给出了大角速度联合控制MSPT的点火位置。从表1中结果可以看出， MSPT共进行了4次点火，三轴消耗的推进器数量分别为13，15和8。
表1大角速度联合控制MSPT点火位置表
Figure BDA0001689304630000171
Figure BDA0001689304630000181
纯磁控方法三轴稳定周期需要3045s，3815s和3351s，4500s后，三轴稳定精度为0.015°/s，0.013°/s和8×10-3°/s。三轴磁矩输出曲线约从900s左右开始振荡幅度下降并逐渐进入稳定输出，磁矩输出绝对值分别为114.53Am2，156.10Am2和 162.56Am2。对比表明，联合控制方法相较与传统磁控方法在相同的控制精度的前提下，控制角速度周期缩短了48％，54％和51％，且大幅度降低了振荡幅度。同时降低了系统的能量输出。
纯MSPT控制所需要的推机器数量为9，12和7对，小于联合控制方法所用的推机器数量。这种方法虽然能够快速降低卫星的旋转，但是控制精度受限与本身推进器的性能以及敏感器误差。仿真结果显示其三轴控制精度仅有0.3325°/s， 0.2823°/s和0.3143°/s，无法满足实际任务需求和姿态捕获任务的初始精度要求，且可靠性无法得以保证。
同样对只采用MSPT辅助控制的方案进行验证，其控制精度与纯磁控方法以及联合控制方法相近，约为0.011°/s；0.012°/s和0.015°/s。控制时长处于两者之间，约为2103s，2245s和1986s。所需的磁矩同样处于两者之间：40.24Am2， 56.43Am2和57.06Am2。表2给出了大角速度辅助控制MSPT的点火位置。从表中结果可以看出，MSPT共进行了4次点火，三轴消耗的推进器数量分别为5，15和7对，略小于联合控制所需要的微推进器数量。在微推进器数量有限且时间宽松的情况下，可以选择使用辅助控制器。
表2大角速度辅助控制MSPT点火位置表
Figure 1
小角速度阻尼时，如图5和图6所示，纯磁控方法角速度控制周期平均约为 3102s，而联合控制方法周期平均仅有722s，相比纯磁控方法提高了76.72％。表给出了小角速度条件下辅助控制MSPT点火位置，表中结果表明，联合控制方法仅需进行两次脉冲，三轴共消耗4对微推进器阵列。因此在初始角速度较小时，在保证控制精度相同的前提下，联合控制律仅需消耗较少的微推进器组合，就能有效缩短控制周期，降低振荡幅度。
大角度速度辅助控制MSPT点火位置表表3：
Figure BDA0001689304630000201
综上，MSPT与磁力矩器联合速度阻尼控制方法在保证角速度的精度的前提下，能够有效缩短控制周期，同时降低系统的振荡幅度，从而验证了该算法的有效性。
姿态捕获联合控制律校验：为了验证本发明提出的联合控制律，将本发明提出的联合控制方法与传统的纯磁控方法以及纯MSPT控制方法进行比较。不考虑其他干扰力矩的影响，选取初始姿态角为[60° 60° 60°]，旋转角速度为 [0.05°/s 0.05°/s 0.05°/s]。纯磁控方法角速度增益系数矩阵为 K1＝diag(3.11×102,3.11×102,3.11×102)，欧拉角增益系数矩阵为K2＝K2/1000；MSPT控制角速度增益矩阵为Kmspt1＝diag(2×10-2,2×10-2,2×10-2)。欧拉角增益系数矩阵为 Kmspt2＝Kmspt1/1000；联合控制方法测试结果如图7、图8、图9所示，图7是四种控制算法的角速度曲线图，其中MSPT与磁力矩器联合姿态控制方法的姿态稳定周期为5466s，6165s和6435，20000s后欧拉角控制精度可以达到0.3°；0.2°和 0.4°，角速度控制精度可达0.8×10-3°/s，1.1×10-3°/s和2×10-3°/s；X轴磁矩输出曲线在762s时达到峰值，Y轴磁矩输出曲线从零时刻起逐步上升至754s时出现波峰，而后逐渐下降并趋于稳定。Z轴曲线从零时刻开始下降，762s达到波峰，而后逐渐上升直至稳定。三轴磁矩输出绝对值分别为6.72Am2，8.40Am2和8.64Am2。表4给出了姿态捕获联合控制MSPT的点火位置。从表中结果可以看出，MSPT 共进行了6次点火，三轴消耗的推进器数量分别为9，16和15对。
Figure BDA0001689304630000211
Figure BDA0001689304630000221
相比而言，纯磁控方法达到上述姿态稳定精度时需要16874s，17923s和 18567s，因此MSPT与磁力矩器联合捕获控制的效率相比于该算法三轴分别提高约208％，191％和188％，从而大大缩短捕获周期。另一方面，捕获时间超过20000s 后，联合控制方法平均角速度控制精度为1.3×10-3°/s，相比于纯磁控方法的 5×10-3°/s提高了约2.3倍。此外，在60000s左右时，联合控制算法角速度曲线和欧拉角曲线均已稳定，而纯磁控方法的两个姿态曲线图的最大幅值分别达到 0.025°/s和40°。此外纯磁控方法轴磁矩输出曲线约从8200s左右开始振荡幅度下降并逐渐进入稳定输出，磁矩输出绝对值分别为114.53Am2，156.10Am2和 162.56Am2。对比表明，联合控制方法相较与传统磁控方法在相同的控制精度的前提下，能有效缩短姿态控制周期，降低振荡幅度，节省系统输出能量。
纯MSPT控制所需要的推机器数量为70，35和100对，远远多于联合控制方法所用的推机器数量，这是由于MSPT方法最小控制精度取决于自身设计，控制精度较差，仿真结果显示其欧拉角角和角速度的稳定精度分别为8°和1.3°/s。另一方面，受限于自身的控制精度，Z轴角速度曲线一直处于持续振荡，直至所有推进器耗光，因而可靠性无法得以保证。因此这种方法虽然能够在30s内快速将姿态调整到一定范围内，但是无法满足对地指向和对地通信等实际任务需求。
综上所述，MSPT与磁力矩器联合姿态捕获控制方法在保证角速度与欧拉角精度的前提下，能够有效缩短控制周期，同时降低系统的振荡幅度，从而验证了该算法的有效性。
速度阻尼联合控制律校验：
Figure BDA0001689304630000231
实验结果表明，组合控制方法相较于传统的纯磁控方法，能够有效提高控制精度，大幅度缩短控制周期，且一定程度上减少了能量损耗；相对于轮控及轮控 /磁控联合控制，虽然控制精度较低，但控制器整体能耗仅为后者的20％，且重量轻。因此，针对于低能耗高精度的姿态控制系统，磁力矩器和MSPT组合姿态控制方法具有明显的优势。

Claims (4)

1.固体微推进器阵列与磁力矩器联合控制方法，其特征在于，具体包括以下步骤：
步骤1、设计固体微推进器阵列优化点火模型，通过所述优化点火模型控制固体微推进器阵列对卫星进行大角度速度阻尼控制，直至卫星的角速度
Figure FDA0002354527380000011
步骤2、基于控制律
Figure FDA0002354527380000012
通过磁力矩器和固体微推进器阵列对卫星进行小角度阻尼联合控制，直至卫星角速度小于等于
Figure FDA0002354527380000013
其中，
Figure FDA0002354527380000014
mb为磁力矩器产生的输出力矩，K1、K2均为正定义系数矩阵，
Figure FDA0002354527380000015
Bb表示本体系地磁场矢量，
Figure FDA0002354527380000016
表示卫星在本体系内的真实角速度，tk表示固体微推进器阵列的工作时间点，k∈{1,2,…,M}，M为大于等于1的整数；Ak为固体微推进器阵列在小角度阻尼控制过程中每次工作产生的输出冲量；
步骤3、基于控制律
Figure FDA0002354527380000017
通过磁力矩器和固体微推进器阵列对卫星进行联合姿态捕获控制；
其中，uc(t)为磁力矩器在捕获控制时的输出力矩，udk为固体微推进器阵列在捕获控制时产生的输出冲量，kc为连续控制律增益，
Figure FDA0002354527380000018
表示本体系磁场的斜对称矩阵、J为转动惯量矩阵，P1、P2均为正定矩阵，θ(t)表示卫星在t时刻欧拉角，
Figure FDA0002354527380000019
为θ(t)的导数，kd为离散控制律增益，
Figure FDA00023545273800000110
为卫星在
Figure FDA00023545273800000111
时刻的欧拉角，
Figure FDA00023545273800000112
为
Figure FDA00023545273800000113
的导数。
2.如权利要求1所述的固体微推进器阵列与磁力矩器联合控制方法，其特征在于，所述固体微推进器阵列优化点火模型具体为：
Figure FDA0002354527380000021
其中，Sk为固体微推进器阵列第k次点火时的目标函数，
Figure FDA0002354527380000022
分别表示向量
Figure FDA0002354527380000023
和
Figure FDA0002354527380000024
中元素的个数，nk为固体微推进器阵列在第k次点火的倍数参数，
Figure FDA0002354527380000025
表示固体微推进器阵列B阵列中第i行第j列的固体微推进器状态，
Figure FDA0002354527380000026
为固体微推进器阵列D阵列中第i行第j列的固体微推进器状态，B阵列
Figure FDA0002354527380000027
和D阵列
Figure FDA0002354527380000028
表示卫星不同的两个面的微推进器阵列；
Figure FDA0002354527380000029
分别为B阵列的状态矩阵和D阵列的状态矩阵中第i行第j列元素。
3.如权利要求1所述的固体微推进器阵列与磁力矩器联合控制方法，其特征在于，步骤2中控制律设计方法为：
选取卫星相对于惯性系的旋转能量作为Lypunov函数V(x)：
Figure FDA00023545273800000210
V(x)相对于时间求导，得出：
Figure FDA00023545273800000211
其中，
Figure FDA00023545273800000212
为
Figure FDA00023545273800000213
的导数，Tb表示作用于卫星的总力矩，
则有：
Figure FDA0002354527380000031
其中，δ(t-tk)表示t-tk时刻的脉冲，teq表示燃烧周期，Atot为固体微推进器阵列在大角度速度阻尼控制阶段的总输出冲量；
令
Figure FDA0002354527380000032
且有正定义系数矩阵K1＝diag{[k1x k1y k1z]}和K2＝diag{[k2x k2y k2z]}，则有：
Figure FDA0002354527380000033
其中，S为斜对称矩阵，
Figure FDA0002354527380000034
根据Lypunov稳定性定理，系统渐进稳定，磁场矢量在轨道平面内变化率接近于0，则有：
Figure FDA0002354527380000035
得出对卫星进行小角度阻尼控制的控制律：
Figure FDA0002354527380000036
4.如权利要求1-3任一所述的固体微推进器阵列与磁力矩器联合控制方法，其特征在于，步骤3中控制律设计方法：
初始状态量x(t0)＝x0，得到在t≠tk的连续函数表达式为：
Figure FDA0002354527380000037
其中，Αc(t)、Bc(t)、Cc(t)、Dc(t)表示连续状态空间函数，
Figure FDA0002354527380000038
为x(t)的导数，yc(t)表示输出矢量，uc(t)表示磁力矩器在卫星捕获控制阶段的输出力矩，
在t＝tk的离散函数表达式为：
Figure FDA0002354527380000041
其中，
Figure FDA0002354527380000042
表示tk +时刻状态向量的导数，Αdk、Bdk、Cdk、Ddk分别表示离散状态空间函数，
Figure FDA0002354527380000043
表示tk -时刻状态向量,ydk表示状态矢量；
根据姿态运动学和动力学方程，得出混合系统的线性模型：
Figure FDA0002354527380000044
其中，
Figure FDA0002354527380000045
为θ(t)的一阶导数，
Figure FDA0002354527380000046
为θ(t)的二阶导数；I3×3表示3*3的单位矩阵；m(t)表示输出磁场；Αd、Bd分别表示矩阵
Figure FDA0002354527380000047
和
Figure FDA0002354527380000048
Αc表示连续状态空间函数；
Figure FDA0002354527380000049
表示磁场的斜对称矩阵；
则系统输入表示为：
Figure FDA00023545273800000410
CN201810585104.1A 2018-06-08 2018-06-08 固体微推进器阵列与磁力矩器联合控制方法 Active CN108891625B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810585104.1A CN108891625B (zh) 2018-06-08 2018-06-08 固体微推进器阵列与磁力矩器联合控制方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810585104.1A CN108891625B (zh) 2018-06-08 2018-06-08 固体微推进器阵列与磁力矩器联合控制方法

Publications (2)

Publication Number Publication Date
CN108891625A CN108891625A (zh) 2018-11-27
CN108891625B true CN108891625B (zh) 2020-03-17

Family

ID=64344317

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810585104.1A Active CN108891625B (zh) 2018-06-08 2018-06-08 固体微推进器阵列与磁力矩器联合控制方法

Country Status (1)

Country Link
CN (1) CN108891625B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109901600A (zh) * 2019-03-08 2019-06-18 宁波天擎航天科技有限公司 一种航天器飞行控制方法、系统及装置
CN115057006A (zh) * 2022-06-15 2022-09-16 中国科学院软件研究所 一种基于强化学习的蒸馏策略评估的方法、装置及介质

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
DE102007033076A1 (de) * 2007-07-13 2009-01-15 Dr. Fritz Faulhaber Gmbh & Co. Kg Elektrischer Kleinstantrieb sowie Rückschlusselement und Verfahren zu dessen Herstellung
US20110274134A1 (en) * 2010-05-04 2011-11-10 Alexander Rubenchik Solar pumped laser microthruster
CN104176275B (zh) * 2014-07-16 2016-05-04 北京控制工程研究所 一种使用动量轮与磁力矩器联合的速率阻尼方法
CN107643688B (zh) * 2017-09-26 2020-08-28 北京航空航天大学 一种针对固体微推力器阵列的两步控制分配方法

Also Published As

Publication number Publication date
CN108891625A (zh) 2018-11-27

Similar Documents

Publication Publication Date Title
Li et al. A survey on moving mass control technology
Horri et al. Practical implementation of attitude-control algorithms for an underactuated satellite
CN104527994B (zh) 异面交叉快变轨道固定时间稳定姿态指向跟踪控制方法
Yeh Sliding-mode adaptive attitude controller design for spacecrafts with thrusters
JP3244717B2 (ja) モーメンタムバイアス宇宙船の傾斜軌道姿勢制御を行うための方法及び装置
CN109911249B (zh) 低推重比飞行器的星际转移有限推力入轨迭代制导方法
US20080315039A1 (en) System and methods for space vehicle torque balancing
Sutherland et al. Attitude control of a 2U cubesat by magnetic and air drag torques
US9428285B2 (en) System and method for managing momentum accumulation
CN112572835B (zh) 一种具有姿态切换的卫星在轨角动量管理及控制方法
CN108891625B (zh) 固体微推进器阵列与磁力矩器联合控制方法
Wu et al. Sliding-mode control for staring-mode spacecraft using a disturbance observer
CN111638643B (zh) 一种位移模式无拖曳控制动力学协调条件确定方法
Yoshimura Optimal formation reconfiguration of satellites under attitude constraints using only thrusters
CN103818564B (zh) 一种采用小推力的航天器轨道维持与对地定向姿态保持一体化控制方法
CN113778047A (zh) 考虑测量误差及综合故障的复杂航天器容错控制方法
CN110502025B (zh) 考虑反作用飞轮特性及功率受限的航天器姿态控制方法
Calaon et al. Momentum management of a spacecraft equipped with a dual-gimballed electric thruster
Бурнашев et al. Micro-satellite orientation control system with rotation wheels
Zhang The output torque estimation of MCMG for agile satellites
Ashok et al. Cmg configuration and steering approach for spacecraft rapid maneuvers
Hua et al. Cooperative game theory-based steering law design of a CMG system
CHUJO et al. Optimization for Solar Angle Free Attitude Motion Model by Solar Sail IKAROS.
Saiki et al. Development of a small-spin-axis controller and its application to a solar sail subpayload satellite
Lee et al. Non-certainty-equivalent adaptive satellite attitude control using solar radiation pressure

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant