[发明专利]一种航天器轨道控制策略的鲁棒性分析方法

权利要求书

1.一种航天器轨道控制策略的鲁棒性分析方法，其特征在于：包括如下步骤：

步骤一：通过Gauss轨道要素摄动方程和平均轨道要素建立航天器轨道运动模型；

采用Gauss轨道要素摄动方程建立的航天器六个轨道要素半长轴a，偏心率e，倾角i，近地点角距ω，升交点赤经Ω和真近点角θ的轨道运动模型如下所示：

其中，

t表示时间，f_r是定义在轨道坐标系S_o的卫星摄动加速度的径向分量；f_u是定义在轨道坐标系S_o的卫星摄动加速度的横向分量，即在轨道平面内垂直于矢经，沿纬度幅角u增大的方向的分量；f_h是定义在轨道坐标系S_o的卫星摄动加速度的副法向分量，即垂直于轨道平面，与动量矩H同向的分量；p是半通径，μ是地球引力常数，r是地心与航天器之间的距离；

航天器轨道要素平均化方法如下：

公式(1)的向量表达形式如下：

其中

σ＝(a,e,i,Ω,ω,M)^T (3)

这里ε＝O(J₂)为通常所取的小参数；等式右端的f₁和f₂为非球形引力加速度的部分，它们分别对应J₂项和其他非球形引力项，并有

f₁/f₀＝O(ε),f₂/f₀＝O(ε²) (5)

公式(2)是一个非常典型的小参数方程；本步骤中采用古在由秀提出的根据非线性力学中平均法的思想，针对主要带谐项摄动的平均根数法；平均根数法将摄动变化项Δσ⁽¹⁾(t),Δσ⁽²⁾(t)，……中不同性质的项，即长期项、长周期项和短周期项区分开，更重要的是在此前提下将参考解取为只包含长期变化项的平均根数；

对Gauss轨道要素摄动方程中积分得出的瞬时轨道要素，利用古在由秀方法计算出其平均轨道要素，从而用于后续的控制和鲁棒性分析；

步骤二：对低轨卫星的地球非球形摄动以及大气阻力摄动进行分析

步骤2.1：地球非球形摄动；

地球非球形摄动中，由于地球的旋转，在运行时间超过一天时，扇形谐函数和田形谐函数的影响或多或少地被“平均掉”，因此对于大多数应用情况，只需考虑带形谐函数的影响；其形式如下：

其中，J_n是n阶带形谐函数的系数，P_n和P′_n是n阶勒让德多项式及其导数，x，y，z是惯性坐标系S_i下的航天器三维坐标，且R_E是地球的平均赤道半径；因此由式(6)得到的惯性坐标系下的地球非球形摄动加速度需要经过坐标转换，以符合式(1)中对于摄动加速度的定义环境；

由惯性坐标系到轨道坐标系的坐标转换分为三步：

其中，u是纬度幅角，u＝ω+θ；

转换矩阵具体形式如下：

因此，满足公式(1)的地球非球形摄动加速度在轨道坐标系的分量形式如下：

步骤3.2：大气阻力摄动；

由空气动力学可知，在大气中运动的航天器受到空气阻力为

其中，C_D为阻力系数，ρ为大气密度，S为横截面积，v_a是航天器质心相对于大气的速度；为简化计算，假设大气是球对称的；大气不随地球旋转；航天器的迎风面积不变；则大气阻力只产生沿轨道切线方向的大气阻力摄动加速度，所述大气阻力摄动加速度写在速度坐标系中如下：

f_n，f_t，f_h分别是大气阻力摄动加速度在速度坐标系的三个分量，其中f_t沿速度方向，f_h垂直于轨道平面指向角动量方向，f_n与其余两个向量垂直，满足右手定则；将速度坐标系S_v中的大气阻力摄动加速度转换到轨道坐标系S_o中,即该坐标系的转换只需要绕速度坐标系的Z轴进行一次欧拉转换，转换矩阵为：

因此，满足公式(1)的大气阻力摄动加速度在轨道坐标系的分量形式如下：

将地球非球形摄动加速度和大气阻力摄动加速度带入公式(1)中的摄动加速度中，即

则得到了带有摄动的航天器轨道运动模型；

步骤三：使用轨道保持策略对航天器轨道进行控制；

名义轨道的轨道要素记为σ^*，由步骤二中的轨道运动模型数值积分所得任意时刻航天器的密切轨道要素记为σ(t)，则轨道要素偏移量记为δσ(t)＝σ(t)-σ^*，记可容忍的最大轨道要素偏差区间为轨道保持采取偏差修正的控制策略，即当则施加脉冲进行轨道修正，其余轨道要素的控制同理可得；

步骤四：对步骤三中计算所需速度脉冲采用微分修正算法进行修正，从而完成精确的轨道机动；

步骤五：建立误差模型

定轨误差分别加在卫星的三个方向的位置和速度上，发动机执行误差则加在三个方向的发动机脉冲分量上；随机误差均服从正态分布；

步骤六：蒙特卡洛仿真分析

采用带有微分修正的轨控策略对带有摄动的航天器轨道进行控制，并加入步骤五中的误差模型，采用MATLAB中normrnd函数对具有同一控制误差标准差的航天器轨道重复多次仿真，消除随机性，模拟体现真实情况；

步骤七、计算模型

对控制误差的分析主要集中于所有未落在控制范围中的点的均值和分布情况作为控制误差的评价方法；首先计算点的均值和标准差，因此控制误差的计算方法如下所示：

其中，

σ'＝0,Δσ∈[c,d] (15)

Δσ是轨道要素相对于其初始值的变化量；最大容许误差范围是[c,d]，其中c和d分别是控制范围的下边界和上边界；N为所有数据点的总个数；在式(15)中，当轨道要素的变化量在最大容许误差范围内时，认为控制误差为零；当轨道要素的变化量超出了控制范围，σ'为超出控制范围的控制结果点到最大容许误差边界的距离；第二个参数是控制误差标准差，它表示了控制误差的分布情况，具体计算公式如下：

第三个参数是控制误差分布比例，表示超出控制范围的点占总数据点的比例；具体计算公式如下：

这三个参数可充分反映控制误差的各项分布；因此，可通过计算ξ，k的值对其在误差存在的情况下的轨道控制性能进行评价，并给出该控制策略鲁棒性的分析结果，为轨控策略的实际应用提供依据和保障。

2.如权利要求1所述的一种航天器轨道控制策略的鲁棒性分析方法，其特征在于：步骤四所述的微分修正方法为：

公式(1)能够简化为

其中x是描述航天器运动的状态量；

假设航天器的初始位置状态向量为x₀，在施加轨道机动ΔV的情况下到达位置为x^t，而其目标到达位置为x^d；

约束函数F(X)表示航天器目标位置与施加机动后的位置之差，如下所示：

F(X)＝x^t-x^d (20)

当F(X)＝0，则表示当前施加的轨道机动ΔV能精确地将航天器机动到目标位置；而为了寻找合适的ΔV，采用牛顿迭代法进行实现；根据Taylor多项式的展开法则，省略后面的高阶项，可以得到如下表达式：

F(X^j)+DF(X^j)(X^j+1-X^j)＝0 (21)

其中DF(X)表示约束函数对状态变量的偏导数，即状态变量的雅克比矩阵，形式如下：

该矩阵表示了轨道要素从初态转移到末态的状态转移矩阵，具体的数值大小可通过给定初值的函数积分得到；

给定一个合理的初值，将式(22)带入式(21)中，则可通过迭代计算得到新的X^j向量，理论上

||F(X^j+1)||＜||F(X^j)|| (23)

于是，通过数次迭代计算，可使F(X)函数的值趋向于0；当

||F(X^j+1)||＜η (24)

其中||F(X^j+1)||是F(X^j+1)的模；η是一给定的极小量；通过上述方法，最终找到合适的机动向量，从而完成精确轨道机动。