一种基于多执行机构的飞行器的姿态控制分配方法

摘要

本发明公开一种基于多执行机构的飞行器的姿态控制分配方法，包括：S1、建立飞行器的姿态运动学方程；S2、通过姿态转换矩阵确定飞行器的姿态角；S3、建立飞行器的姿态动力学方程；S4、根据控制律计算飞行器的期望控制力矩；S5、利用控制分配算法获得当前时刻飞行器的各执行机构的状态组合。本发明适合飞行器大角度姿态快速机动与高精度稳定的机动任务需求。

说明书

技术领域

本发明涉及飞行器智能控制领域。更具体地，涉及一种基于多执行机构的飞行器的姿态控制分配方法。

背景技术

基于多执行机构的飞行器是一种在智能材料、智能自适应非线性控制分配技术、微推进系统的基础上，发挥各自的优势，组合在一起的多用途、多功能的新型片状材料组合体式飞行器。该飞行器具有可变蒙皮，具备变轨机动性能大幅提升，实现最佳气动性能，提高可重复使用能力，可重构能力大幅提升等优点。

该型飞行器需要在姿态机动、高精度稳定等控制模式之间快速切换.快速机动为攻击敌方目标、规避敌方攻击争取更多的时间，而机动后的快速高精度稳定则是为了保证飞行平稳性，使其在高精度飞行中提高可重复使用能力。因此大角度姿态快速机动与高精度稳定成为基于多执行机构飞行器的关键技术之一。

执行机构的选择决定飞行器控制系统的机动能力，选择不同的执行机构也对系统的机动性能起决定作用。针对上述复杂任务，单一的执行机构已无法满足需求，多个执行机构联合控制方式能够发挥各个机构组合的长处，实现优势互补，以最小代价实现控制目标。

因此，需要提供一种基于多执行机构的飞行器的姿态控制分配方法。

发明内容

本发明的目的在于提供一种基于多执行机构的飞行器的姿态控制分配方法，解决传统姿态控制分配方法假设条件理想化的问题。

为达到上述目的，本发明采用下述技术方案：

一种基于多执行机构的飞行器的姿态控制分配方法，该方法包括如下步骤：

S1、建立飞行器的姿态运动学方程；

S2、通过姿态转换矩阵确定飞行器的姿态角；

S3、建立飞行器的姿态动力学方程；

S4、根据控制律计算飞行器的期望控制力矩；

S5、利用控制分配算法获得当前时刻飞行器的各执行机构的状态组合。

优选地，步骤S1中飞行器的姿态运动学方程为：

其中，q₀、q₁、q₂和q₃为当前时刻飞行器姿态的四元数；和为当前时刻飞行器姿态的四元数导数；w_x为飞行器绕X轴转动角速度；w_y为飞行器绕Y轴转动角速度；w_z为飞行器绕Z轴转动角速度。

优选地，步骤S2中姿态转换矩阵为：

优选地，步骤S3中飞行器的姿态动力学方程为：

其中，I_x为飞行器绕X轴转动惯量；I_y为飞行器绕Y轴转动惯量；I_z为飞行器绕Z轴转动惯量；dw_x为飞行器绕X轴转动角速度的导数；dw_y为飞行器绕Y轴转动角速度的导数；dw_z为飞行器绕Z轴转动角速度的导数；T_x为作用于飞行器上的所有外力对质心之力矩在X轴上的分量；T_y为作用于飞行器上的所有外力对质心之力矩在Y轴上的分量；T_z为作用于飞行器上的所有外力对质心之力矩在Z轴上的分量。

优选地，步骤S4进一步包括如下子步骤：

S4.1、计算当前时刻飞行器姿态的四元数微分：

dq₀＝q_c0×q₀+q_c1×q₁+q_c2×q₂+q_c3×q₃

dq₁＝-q_c1×q₀+q_c1×q₁+q_c2×q₂+q_c3×q₃

dq₂＝-q_c3×q₁-q_c4×q₂+q_c2×q₂+q_c3×q₃

dq₃＝-q_c3×q₀+q_c2×q₁-q_c1×q₂+q_c0×q₃

其中，dq₀、dq₁、dq₂和dq₃为当前时刻飞行器姿态的四元数微分；q_c0、q_c1、q_c2和q_c3为飞行器期望姿态的四元数；

S4.2、计算飞行器的期望控制力矩的非线性补偿项T_R：

T_R＝J×dw_c+w×J×w

其中，w为飞行器当前姿态角速度；dw_c为飞行器期望角加速度；J为飞行器的转动惯量矩阵，

S4.3、计算飞行器的期望控制力矩的比例项T_p：

其中，K_p为比例系数矩阵；

S4.4、计算飞行器的期望控制力矩的微分项T_d：

T_d＝-K_d×J×(w-w_c)

其中，w_c为飞行器期望角速度；

S4.5、计算飞行器的期望控制力矩的误差补偿项T_e：

其中，ww为误差补偿项归一化系数，K_ep为第一误差补偿比较项系数，K_ed为第二误差补偿比较项系数；

S4.6、计算飞行器的期望控制力矩T_c：

T_c＝T_R+T_p+T_e+T_d。

优选地，步骤S5进一步包括如下子步骤：

S5.1、确定控制指令所在的象限；

S5.2、根据控制指令的方向选择执行机构组合；

S5.3、确定执行机构组合中各执行机构的可能状态；

S5.4、计算执行机构组合中各执行机构的可能状态对应的力矩；

S5.5、比较得出执行机构组合中各执行机构的可能状态对应的力矩之和与期望控制力矩之间误差最小值对应的各执行机构的状态组合，作为当前时刻飞行器的各执行机构的状态组合。

本发明的有益效果如下：

本发明所述技术方案在执行机构存在冗余的条件下，可以实现某一优化准则或约束下的优化分配，提高系统对执行机构故障的容错能力，为冗余异构执行机构情况下的姿态快速机动与高精度稳定控制问题提供了一条可行的解决途径。具体的，本发明所述技术方案在建立执行机构控制方案基础上，提出了一种基于MEMS微执行机构的姿态控制方法，建立了姿态运动学和动力学方程，通过制定控制分配方案给出各个微执行机构的状态，即确定动态分配控制策略，仿真中设置三种典型仿真场景，从仿真结果中可以看出，当前的姿态角和姿态角速度与期望的姿态角和姿态角速度误差在0.5度范围之内，满足误差要求，并且开关机次数满足该飞行器执行任务的需要，实现了最大机动角速度以及目标函数约束下的最优机动.该方案适合飞行器大角度姿态快速机动与高精度稳定的机动任务需求。

附图说明

下面结合附图对本发明的具体实施方式作进一步详细的说明。

图1示出基于多执行机构的飞行器的姿态控制分配方法的流程图。

图2示出飞行器的当前姿态角和期望控制角度图。

图3示出飞行器的当前姿态角速度和期望控制姿态角速度图。

图4示出飞行器的3s后的姿态跟踪误差图。

图5示出飞行器的期望控制力矩图。

图6示出飞行器的实际输出力矩图。

图7示出飞行器的期望控制力矩示意图。

图8示出飞行器的微执行机构开关曲线图。

具体实施方式

为了更清楚地说明本发明，下面结合优选实施例和附图对本发明做进一步的说明。附图中相似的部件以相同的附图标记进行表示。本领域技术人员应当理解，下面所具体描述的内容是说明性的而非限制性的，不应以此限制本发明的保护范围。

如图1所示，本实施例提供的基于多执行机构的飞行器(以下简称飞行器)的姿态控制分配方法，包括如下步骤：

S1、建立飞行器的姿态运动学方程；

S2、通过姿态转换矩阵确定飞行器的姿态角；

S3、建立飞行器的姿态动力学方程；

S4、根据控制律计算飞行器的期望控制力矩；

S5、利用控制分配算法获得当前时刻飞行器的各执行机构的状态组合。

其中，

步骤S1中飞行器的姿态运动学方程为：

其中，q₀、q₁、q₂和q₃为当前时刻飞行器姿态(惯性坐标系到弹体坐标系)的四元数；和为当前时刻飞行器姿态的四元数导数；w_x为飞行器绕(惯性坐标系)X轴转动角速度；w_y为飞行器绕Y轴转动角速度；w_z为飞行器绕Z轴转动角速度。

步骤S2的具体过程为：

按照Z-Y-X顺序旋转惯性坐标系，将飞行器的姿态四元数转换为飞行器的欧拉角。姿态转换矩阵为：

步骤S3中飞行器的姿态动力学方程为：

其中，I_x为飞行器绕X轴转动惯量；I_y为飞行器绕Y轴转动惯量；I_z为飞行器绕Z轴转动惯量；dw_x为飞行器绕X轴转动角速度的导数；dw_y为飞行器绕Y轴转动角速度的导数；dw_z为飞行器绕Z轴转动角速度的导数；T_x为作用于飞行器上的所有外力(包括空气动力和推力)对质心之力矩在(弹体坐标系)X轴上的分量；T_y为作用于飞行器上的所有外力(包括空气动力和推力)对质心之力矩在Y轴上的分量；T_z为作用于飞行器上的所有外力(包括空气动力和推力)对质心之力矩在Z轴上的分量。

步骤S4“根据控制律计算飞行器的期望控制力矩”进一步包括如下子步骤：

S4.1、计算当前时刻飞行器姿态的四元数微分：

dq₀＝q_c0×q₀+q_c1×q₁+q_c2×q₂+q_c3×q₃

dq₁＝-q_c1×q₀+q_c1×q₁+q_c2×q₂+q_c3×q₃

dq₂＝-q_c3×q₁-q_c4×q₂+q_c2×q₂+q_c3×q₃

dq₃＝-q_c3×q₀+q_c2×q₁-q_c1×q₂+q_c0×q₃

其中，dq₀、dq₁、dq₂和dq₃为当前时刻飞行器姿态的四元数微分；q_c0、q_c1、q_c2和q_c3为飞行器期望姿态的四元数；

S4.2、计算飞行器的期望控制力矩的非线性补偿项T_R：

T_R＝J×dw_c+w×J×w

其中，w为飞行器当前姿态角速度；dw_c为飞行器期望角加速度；J为飞行器的转动惯量矩阵，

S4.3、计算飞行器的期望控制力矩的比例项T_p：

其中，K_p为比例系数矩阵；

S4.4、计算飞行器的期望控制力矩的微分项T_d：

T_d＝-K_d×J×(w-w_c)

其中，w_c为飞行器期望角速度；

S4.5、计算飞行器的期望控制力矩的误差补偿项T_e：

其中，ww为误差补偿项归一化系数，K_ep为第一误差补偿比较项系数，K_ed为第二误差补偿比较项系数；

S4.6、计算飞行器的期望控制力矩T_c：

T_c＝T_R+T_p+T_e+T_d。

步骤S5“利用控制分配算法获得当前时刻飞行器的各个执行机构的状态组合”的具体过程为：

已知凸锥法向量：

利用扩展几何法确定出50度安装角构型的执行机构组合如表1-14所示，根据各执行机构产生力矩在三维空间的分布图像，可以进一步确定出预装订的数据：各执行机构组合围成多锥形体的锥面法向量、各象限中包含的执行机构组合、以及期望控制力矩属于各执行机构组合围成区域的内部时所应满足的点积的正负条件。计算的结果分别如表1-15以及表1-16所示：

表1-14基于扩展几何法制定执行机构组合表的结果

表1-15各执行机构组合围成多锥形体的锥面法向量

锥面法向量锥面法向量锥面法向量<T₇,T₁>_n1＝[0.8609,-0.5068,0.0443]^T<T₈,T₆>n₅＝[-0.8609,-0.5068,-0.0443]^T<T₁,T₈>n₉＝[0.2042,-0.6922,0.6922]^T<T₁,T₃>n₂＝[0.8609,-0.0443,-0.5068]^T<T₂,T₈>n₆＝[-0.8609,0.0443,-0.5068]^T<T₃,T₂>n₁₀＝[0.2042,-0.6922,-0.6922]^T<T₃,T₅>n₃＝[0.8609,0.5068,-0.0443]^T<T₄,T₂>n₇＝[-0.8609,0.5068,0.0443]^T<T₅,T₄>n₁₁＝[0.2042,0.6922,-0.6922]^T<T₅,T₇>n₄＝[0.8609,0.0443,0.5068]^T<T₆,T₄>n₈＝[-0.8609,-0.0443,0.5068]^T<T₇,T₆>n₁₂＝[0.2042,0.6922,0.6922]^T

表1-15中的<T_i,T_j>代表由第i个和第j个执行机构的力矩矢量构成的锥面，在实际使用过程中，需要注意法线的方向与表1-16中判断逻辑的统一。

表1-16各象限中包含的执行机构组合以及判断逻辑

表1-16中的(u_3d,n_i)代表期望控制力矩u_3d与面法向量n_i的点积，而符号“&”代表逻辑与。从表1-16可以看出对于50度推力构型来讲，在预装定数据时，需要存储12个锥面方向矢量即可，它们分别为n₁、n₂、n₃、n₄、n₅、n₆、n₇、n₈、n₉、n₁₀、n₁₁以及n₁₂。而每一个象限会包含三个可能的执行机构组合，在确定最优执行机构组合时，需要计算三个点积即可。

姿态角响应曲线如图2所示，姿态角速度响应曲线如图3所示，3s后的姿态跟踪误差曲线如图4所示，期望控制力矩曲线如图5所示，执行机构实际输出力矩的曲线如图6所示，5s后的期望控制力矩与执行机构实际输出力矩的曲线如图7所示，8个执行机构的开关状态曲线如图8所示。

最优执行机构组合：

步骤S5“利用控制分配算法获得当前时刻飞行器的各个执行机构的状态组合”进一步包括如下子步骤：

S5.1、确定控制指令所在的象限：

设定象限符号为num，飞行器的期望控制力矩的标志符flag，满足：

num＝num+flag×2^i-1(i＝1，2，3)

S5.2、根据控制指令的方向选择执行机构组合：

设定(u_3d,n₁)为(1)式，(u_3d,n₂)为(2)式，(u_3d,n₉)为(3)式，(u_3d,n₃)为(4)式，(u_3d,n₁₀)为(5)式，(u_3d,n₄)为(6)式，(u_3d,n₁₁)为(7)式，(u_3d,n₁₂)为(8)式，(u_3d,n₅)为(9)式，(u_3d,n₆)为(10)式，(u_3d,n₇)为(11)式，(u_3d,n₈)为(12)式。

如果num＝7:

若(1)式≥0且(2)式≥0，则执行机构组合为[1 3 5 7]

若(1)式<0且(3)式<0，则执行机构组合为[1 6 7 8]

若(2)式<0且(3)式≥0，则执行机构组合为[1 2 3 8]

如果num＝6:

若(2)式≥0且(4)式≥0，则执行机构组合为[1 3 5 7]

若(2)式<0且(5)式<0，则执行机构组合为[1 2 3 8]

若(4)式<0且(5)式≥0，则执行机构组合为[2 3 4 5]

如果num＝4:

若(4)式≥0且(6)式≥0，则执行机构组合为[1 3 5 7]

若(4)式<0且(7)式<0，则执行机构组合为[2 3 4 5]

若(6)式<0且(7)式≥0，则执行机构组合为[4 5 6 7]

如果num＝5:

若(6)式≥0且(1)式≥0，则执行机构组合为[1 3 5 7]

若(6)式<0且(8)式<0，则执行机构组合为[4 5 6 7]

若(1)式<0且(8)式≥0，则执行机构组合为[1 6 7 8]

如果num＝3:

若(9)式≥0且(10)式≥0，则执行机构组合为[2 4 6 8]

若(9)式<0且(3)式<0，则执行机构组合为[1 2 3 8]

若(10)式<0且(3)式≥0，则执行机构组合为[2 3 4 5]

如果num＝2:

若(10)式≥0且(11)式≥0，则执行机构组合为[2 4 6 8]

若(10)式<0且(5)式<0，则执行机构组合为[1 2 3 8]

若(11)式<0且(5)式≥0，则执行机构组合为[2 3 4 5]

如果num＝0：

若(11)式≥0且(12)式≥0，则执行机构组合为[2 4 6 8]

若(11)式<0且(7)式<0，则执行机构组合为[2 3 4 5]

若(12)式<0且(7)式≥0，则执行机构组合为[4 5 6 7]

如果num＝1：

若(12)式≥0且(9)式≥0，则执行机构组合为[2 4 6 8]

若(12)式<0且(8)式<0，则执行机构组合为[4 5 6 7]

若(9)式<0且(8)式≥0，则执行机构组合为[1 6 7 8]。

S5.3、确定执行机构组合中各执行机构的可能状态

计算在16种情况(8个推力器，每个推力器有开和关两种状态)下所有可能的开关状态，设定第j个执行机构的控制指令为U₀(j)，第j个执行机构的第i种开关组合UU(i,j)，将0-15种组合转化为四位二进制数：

abcd即为四位二进制数的控制指令(i＝0-15)。

S5.4、计算执行机构组合中各执行机构的可能状态对应的力矩：

若U₀(j)＝0且UU(i,j)＝0，则：

F(j)＝0

若U₀(j)＝1且UU(i,j)＝1，则：

F(j)＝f_max

其中，f_max为执行机构最大推力；

若U₀(j)≠1或UU(i,j)≠1

F(j)＝0.5×f_max

T(:,i)＝M×F'(j)

上式中，M为执行机构安装矩阵(代表控制力矩矩阵，将八个执行机构所受的力乘以力矩得到俯仰、偏航和滚转三个方向的控制力矩，可以表明8个执行机构)，F(j)为第j个执行机构产生的力矩，T(:,i)为16种开关状态下产生的控制力矩之和。

S5.5、比较得出执行机构组合中各执行机构的可能状态对应的力矩之和与期望控制力矩之间误差最小值对应的各执行机构的状态组合，作为当前时刻飞行器的各执行机构的状态组合：

设定error＝||T(:,1)-T_c||，

如果||T(:,i)-T_c||＜error

则error＝||T(:,i)-T_c||(i＝0-15)

通过遍历找到各执行机构产生的力矩之和与期望控制力矩之间误差，记录误差最小值对应的第i种开关状态下的组合序号。

本发明能够实现针对基于多执行机构的高精度姿态控制，其突出特点是发挥多个执行机构组合控制优势，从而实现高精度姿态控制，满足飞行器在姿态机动、高精度稳定等控制模式之间快速切换的需要。

显然，本发明的上述实施例仅仅是为清楚地说明本发明所作的举例，而并非是对本发明的实施方式的限定，对于所属领域的普通技术人员来说，在上述说明的基础上还可以做出其它不同形式的变化或变动，这里无法对所有的实施方式予以穷举，凡是属于本发明的技术方案所引伸出的显而易见的变化或变动仍处于本发明的保护范围之列。