推力变向无人机的姿态控制系统及控制方法

摘要

本发明公开了一种推力变向无人机的姿态控制系统及控制方法。本发明首次建立了推力变向无人机的非线性全量数学模型和线性化纵横向数学模型，并利用融合估计和最优控制的方法，提出了本发明的基于姿态回路信息融合的推力变向无人机的姿态控制系统及控制方法。对于推力变向无人机的设计开发具有极为重要的现实意义。本发明在控制性能方面优于传统单变量控制方法，且计算简单，完全满足实际飞行控制系统的实时性要求。

说明书

技术领域

本发明涉及无人机的姿态控制系统及控制方法，尤其涉及一种推力变向无人机的姿态控制系统及控制方法，属于飞行控制技术领域。

背景技术

侦察攻击一体化无人机提高了部队的战场空间感知能力、对高风险目标的突防能力、电子战能力和对敌防空系统的压制能力，实现了对敌机动目标精确、快速攻击，是一种适应信息战特点和需求的全新武器系统，使无人机发展的重要趋势之一。

推力矢量控制技术是目前比较先进的飞行控制手段，在现代飞行器中得到成功应用。推力矢量技术在有人机上已被验证能够改善飞机的隐身性、机动能力和敏捷性等，因此推力矢量技术的应用已成为无人机察打一体化发展的重要技术之一。所谓推力矢量技术，是指发动机推力通过喷管或尾喷流的偏转产生的推力分量来替代原飞机的操纵面或增强飞机的操纵功能，对飞机的飞行进行实时控制的技术。与一般的推力矢量技术相比较而言，变推力轴线技术既不特指喷气发动机，也不特指燃气舵，它可以通过增加额外的辅助装置，实现对螺旋桨飞机或喷气飞机的推力偏转，从而使飞机获得额外的控制能力。变推力轴线技术近年来刚开始有人研究，仍然处于理论研究阶段，理论研究和数值仿真结果表明它具有与推力矢量技术相似的作用。此处的推力变向技术，是涡喷式无人机的推力矢量技术与螺旋桨无人机的变推力轴线技术的统称。所谓推力变向技术，是指通过控制发动机使其推力产生纵向或横向的偏转，从而对飞机产生附加的推力分量和附加推力力矩，实现对原有气动操纵面的增强或替代的一种技术。

关于推力矢量无人机、变推力轴线无人机的非线性数学模型已有公开文献，但是，这些文献中的无人机飞行控制系统主要实现姿态控制、高度控制和航路控制等飞行模态。飞行航路控制系统的任务是使无人机能自动地按预定航线飞行，对飞行器航向和偏离参考航路线距离的控制。

为了适应无人机察打一体化的发展趋势，对推力变向无人机的姿态控制提出更高的要求。传统无人机飞行控制策略主要基于经典的单变量控制回路频域或根轨迹方法。然而，无人机飞行控制系统是一个典型的多输入多输出系统，传统方法难于协调处理。无人机飞行滚转角速度引起偏航力矩，偏航角速度引起滚转力矩，而侧滑角会引起滚转和偏航力矩，同时滚转和航向操纵也常常是交联的。现代控制理论的优势就在于解决多变量系统的控制问题，与经典设计方法相比，现代控制技术不再是单回路逐个设计，矩阵运算使得所有的控制增益同时计算，也即所有的控制回路同时闭合，并且能够利用更多的系统状态信息。利用现代控制技术设计飞行控制系统的试验和验证一直在不断地进行，并取得了一些成功的经验。最优控制理论是现代控制理论的重要组成部分，传统最优控制方法主要有极小值原理、动态规划方法、Lyapunov函数方法、矩阵配方法等。传统方法求解最优控制问题都是基于严格的数学原理，推导过程比较复杂，且不适合解决非线性、时滞性、状态相关或控制受约束等系统的最优控制问题。融合估计和最优控制都属于二次型最优化问题，融合控制的思想启发于融合估计，是一种崭新的最优控制方法。利用融合估计的方法来求解最优控制问题，推导过程因为具有物理含义而显得容易理解，并且有利于解决上述特定系统的最优控制问题。

发明内容

本发明所要解决的技术问题在于解决现有无人机姿态控制技术所存在的不足，提供一种基于姿态回路信息融合的推力变向无人机的姿态控制系统及控制方法。

本发明采用以下技术方案：

一种推力变向无人机的姿态控制系统，所述推力变向无人机采用单发动机，发动机安装位置处于飞机纵向对称面上，所述姿态控制系统包括根据给定姿态指令信号解算控制律，得到姿态控制信号并将姿态控制信号输出给执行机构的姿态控制器，所述姿态控制器包括纵向姿态信息融合控制器和横向姿态信息融合控制器；

所述纵向姿态信息融合控制器的输入为期望的飞行速度增量ΔV^*、期望的俯仰角增量Δθ^*、期望的迎角增量Δα^*、期望的俯仰角速度增量Δq^*、传感器反馈的俯仰角速度增量Δq、传感器反馈的俯仰角增量Δθ、传感器反馈的迎角增量Δα、传感器反馈的飞行速度增量ΔV、期望的升降舵偏角增量Δδ_e^*、油门开度增量Δδ_T^*、发动机推力纵向偏角增量Δα_T^*；输出为升降舵偏角增量Δδ_e、油门开度增量Δδ_T、发动机推力纵向偏角增量Δα_T，其信息融合最优控制律如下，

${\hat{U}}_{\mathrm{lon}}\left(k\right)={[R_{\mathrm{lon}}\left(k\right)+{R_{\mathrm{lon}}}^T\left(k\right){R_{\mathrm{lon}}}^{-1}(k+1)B_{\mathrm{lon}}\left(k\right)]}^{-1}\{{B_{\mathrm{lon}}}^T\left(k\right){P_{\mathrm{lon}}}^{-1}(k+1)$

$[{\hat{X}}_{\mathrm{lon}}(k+1)-A_{\mathrm{lon}}\left(k\right)X_{\mathrm{lon}}\left(k\right)]+R_{\mathrm{lon}}\left(k\right){U_{\mathrm{lon}}}^{*}\left(k\right)\}$

式中，为虚拟纵向状态向量X_lon的融合估计，P_lon^-1为融合估计的融合信息量，且融合估计及其信息量P_lon^-1由下面公式逆向递推计算获得，

$\left(\begin{array}{c}{\hat{X}}_{\mathrm{lon}}\left(j\right)=P_{\mathrm{lon}}\left(j\right)\left\{{A_{\mathrm{lon}}}^T\left(j\right){[P_{\mathrm{lon}}(j+1)+B_{\mathrm{lon}}\left(j\right){R_{\mathrm{lon}}}^{-1}\left(j\right){B_{\mathrm{lon}}}^T\left(j\right)]}^{-1}\right[{\hat{X}}_{\mathrm{lon}}(j+1)-B_{\mathrm{lon}}\left(j\right){U_{\mathrm{lon}}}^{*}\left(j\right)]+Q_{\mathrm{lon}}\left(j\right){X_{\mathrm{lon}}}^{*}\left(j\right)\}\\ {P_{\mathrm{lon}}}^{-1}\left(j\right)=Q_{\mathrm{lon}}\left(j\right)+{A_{\mathrm{lon}}}^T\left(j\right){[P_{\mathrm{lon}}(j+1)+B_{\mathrm{lon}}\left(j\right){R_{\mathrm{lon}}}^{-1}\left(j\right){B_{\mathrm{lon}}}^T\left(j\right)]}^{-1}{A}_{\mathrm{lon}}\left(j\right),{P_{\mathrm{lon}}}^{-1}\left(N\right)=Q_{\mathrm{lon}}\left(N\right)\\ j=k+1,k+2,L,N-1\end{array}\right),$

式中，X_lon＝[ΔV，Δα，Δθ，Δq]^T，U_lon＝[Δδ_T，Δδ_e，Δα_T]^T，A_lon(k)，B_lon(k)为连续时间纵向线性化运动方程经过离散化后的系统矩阵，R_lon(k)，Q_lon(k)为最优性能指标权重矩阵，X_lon^*＝[ΔV^*，Δα^*，Δθ^*，Δq^*]^T为期望的纵向状态向量，U_lon^*＝[Δδ_T^*，Δδ_e^*，Δα_T^*]^T为期望的纵向控制向量；

所述横向姿态信息融合控制器的输入为期望的侧滑角增量Δβ^*、期望的滚转角增量Δφ^*、期望的滚转角速度增量Δp^*、期望的偏航角速度增量Δr^*、传感器反馈的滚转角速度增量Δp、传感器反馈的滚转角增量Δφ、传感器反馈的偏航角速度增量Δr、传感器反馈的侧滑角增量Δβ、期望的副翼偏角增量Δδ_α^*、期望的方向舵偏角增量Δδ_r^*、期望的发动机推力横向偏角增量Δβ_T^*，输出为副翼偏角增量Δδ_α、方向舵偏角增量Δδ_r、发动机推力横向偏角增量Δβ_T，其信息融合最优控制律如下，

${\hat{U}}_{\mathrm{lat}}\left(k\right)={[R_{\mathrm{lat}}\left(k\right)+{B_{\mathrm{lat}}}^T\left(k\right){P_{\mathrm{lat}}}^{-1}(k+1)B_{\mathrm{lat}}\left(k\right)]}^{-1}\{{B_{\mathrm{lat}}}^T\left(k\right){P_{\mathrm{lat}}}^{-1}(k+1)$

$[{\hat{X}}_{\mathrm{lat}}(k+1)-A_{\mathrm{lat}}\left(k\right)X_{\mathrm{lat}}\left(k\right)]+R_{\mathrm{lat}}\left(k\right){U_{\mathrm{lat}}}^{*}\left(k\right)\}\text{'}$

式中，为虚拟横向状态向量X_lat的融合估计，P_lat^-1为融合估计的融合信息量，且融合估计及其信息量P_lat^-1由下面公式逆向递推计算获得，

$\left(\begin{array}{c}{\hat{X}}_{\mathrm{lat}}\left(j\right)=P_{\mathrm{lat}}\left(j\right)\left\{{A_{\mathrm{lat}}}^T\left(j\right){[P_{\mathrm{lat}}(j+1)+B_{\mathrm{lat}}\left(j\right){R_{\mathrm{lat}}}^{-1}\left(j\right){B_{\mathrm{lat}}}^T\left(j\right)]}^{-1}\right[{\hat{X}}_{\mathrm{lat}}(j+1)-B_{\mathrm{lat}}\left(j\right){U_{\mathrm{lat}}}^{*}\left(j\right)]+Q_{\mathrm{lat}}\left(j\right){X_{\mathrm{lat}}}^{*}\left(j\right)\}\\ {P_{\mathrm{lat}}}^{-1}\left(j\right)=Q_{\mathrm{lat}}\left(j\right)+{A_{\mathrm{lat}}}^T\left(j\right){[P_{\mathrm{lat}}(j+1)+B_{\mathrm{lat}}\left(j\right){R_{\mathrm{lat}}}^{-1}\left(j\right){B_{\mathrm{lat}}}^T\left(j\right)]}^{-1}{A}_{\mathrm{lat}}\left(j\right),{P_{\mathrm{lat}}}^{-1}\left(N\right)=Q_{\mathrm{lat}}\left(N\right)\\ j=k+1,k+2,L,N-1\end{array}\right),$

式中，X_lat＝[Δβ，Δφ，Δp，Δr]^T，U_lat＝[Δδ_a，Δδ_r，Δβ_T]^T，A_lat(k)，B_lat(k)为连续时间横侧向线性化运动方程经过离散化后的系统矩阵，R_lat(k)，Q_lat(k)为最优性能指标权重矩阵，X_lat^*＝[Δβ^*，Δφ^*，Δp^*，Δr^*]^T为期望的横向状态向量，U_lat^*＝[Δδ_a^*，Δδ_r^*，Δβ_T^*]^T为期望的横向控制向量，其它参数为中间矩阵变量；

所述纵向线性化运动方程如下式所示，

式中，

X_lon＝[ΔV，Δα，Δθ，Δq]^T

U_lon＝[Δδ_T，Δδ_e，Δα_T]^T

所述横侧向线性化运动方程如下式所示，

式中，

X_lat＝[Δβ，Δφ，Δp，Δr]^T，

U_lat＝[Δδ_a，Δδ_r，Δβ_T]^T，

$A_{\mathrm{lat}}=\left(\begin{array}{cccc}\frac{Y_{\beta }}{{\mathrm{mV}}^{*}}& \frac{{g\cos \theta }^{*}}{V^{*}}& {\sin \alpha }^{*}+\frac{Y_p}{{\mathrm{mV}}^{*}}& \frac{Y_r}{{\mathrm{mV}}^{*}}-{\cos \alpha }^{*}\\ 0& 0& 1& {\tan \theta }^{*}\\ \frac{I_z{\bar{L}}_{\beta }+I_{\mathrm{zx}}{N}_{\beta }}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& 0& \frac{I_z{\bar{L}}_p+I_{\mathrm{zx}}{N}_p}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_z{\bar{L}}_r+I_{\mathrm{zx}}{N}_r}{I_x{I}_z-{I_{\mathrm{zx}}}^2}\\ \frac{I_{\mathrm{zx}}{\bar{L}}_{\beta }+I_x{N}_{\beta }}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& 0& \frac{I_{\mathrm{zx}}{\bar{L}}_p+I_x{N}_p}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_{\mathrm{zx}}{\bar{L}}_r+I_x{N}_r}{I_x{I}_z-{I_{\mathrm{zx}}}^2}\end{array}\right),$

$B_{\mathrm{lat}}=\left(\begin{array}{ccc}\frac{Y_{\mathrm{\delta a}}}{{\mathrm{mV}}^{*}}& \frac{Y_{\mathrm{\delta r}}}{{\mathrm{mV}}^{*}}& \frac{T^{*}}{{\mathrm{mV}}^{*}}\\ 0& 0& 0\\ \frac{I_z{\bar{L}}_{\mathrm{\delta a}}+I_{\mathrm{zx}}{N}_{\mathrm{\delta a}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_z{\bar{L}}_{\mathrm{\delta r}}+I_{\mathrm{zx}}{N}_{\mathrm{\delta r}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_z{\bar{L}}_{\mathrm{\beta T}}+I_{\mathrm{zx}}{N}_{\mathrm{\beta T}}}{I_z{I}_x-{I_{\mathrm{zx}}}^2}\\ \frac{I_{\mathrm{zx}}{\bar{L}}_{\mathrm{\delta a}}+I_x{N}_{\mathrm{\delta a}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_{\mathrm{zx}}{\bar{L}}_{\mathrm{\delta r}}+I_x{N}_{\mathrm{\delta r}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_{\mathrm{zx}}{\bar{L}}_{\mathrm{\beta T}}+I_x{N}_{\mathrm{\beta T}}}{I_z{I}_x-{I_{\mathrm{zx}}}^2}\end{array}\right).$

一种推力变向无人机的姿态控制方法，所述推力变向无人机采用单发动机，发动机安装位置处于飞机纵向对称面上，所述推力变向无人机采用权利要求1所述姿态控制系统进行姿态控制，包括以下步骤：

步骤1、姿态控制系统接收给定姿态指令信号、传感器反馈姿态信号以及期望控制信号，所述给定姿态指令信号包括，纵向期望姿态信号X_lon^*＝[ΔV^*，Δα^*，Δθ^*，Δq^*]^T，横向期望姿态信号X_lat^*＝[Δβ^*，Δφ^*，Δp^*，Δr^*]^T，所述传感器反馈姿态信号包括，纵向实际姿态信号X_lon＝[ΔV，Δα，Δθ，Δq]^T，横向实际姿态信号X_lat＝[Δβ，Δφ，Δp，Δr]^T，所述期望控制信号包括，纵向期望控制信号U_lon^*＝[Δδ_T^*，Δδ_e^*，Δα_T^*]^T，横向期望控制信号U_lat^*＝[Δδ_a^*，Δδ_r^*，Δβ_T^*]^T；

步骤2、姿态控制系统根据接收到的给定姿态指令信号解算姿态控制律，得到纵向姿态控制信号U_lon＝[Δδ_T，Δδ_e，Δα_T]^T，横向姿态控制信号U_lat＝[Δδ_a，Δδ_r，Δβ_T]^T；

步骤3、姿态控制系统将姿态控制信号发送给执行机构，执行机构操纵无人机的气动舵面、油门开度、推力方向，实现控制无人机飞行姿态和飞行速度的目的。

本发明首次建立了推力变向无人机的非线性全量数学模型和线性化纵横向数学模型，并利用融合估计和最优控制的方法，提出了本发明的基于姿态回路信息融合的推力变向无人机的姿态控制系统及控制方法。对于推力变向无人机的设计开发具有极为重要的现实意义。

附图说明

图1为本发明所建立的推力变向无人机非线性全量模型的结构示意图；

图2为本发明的推力变向无人机的姿态控制系统控制原理示意图；

图3为采用传统单回路姿态控制方法所得到的仿真结果；

图4为采用本发明姿态控制方法所得到的仿真结果。

具体实施方式

下面结合附图对本发明的技术方案进行详细说明：

要得到推力变向无人机的姿态控制系统及方法，首先要建立推力变向无人机的非线性全量模型。该推力变向无人机采用单发动机，其安装位置处于飞机纵向对称面上，若是螺旋桨无人机，则推力变向机械装置安装于发动机上，可实现发动机的两自由度偏转运动。设发动机推力T在飞机对称面上的投影与x_b轴的夹角记为α_T，称之为推力纵向偏转角，投影在机体轴下为正；推力在x_b-y_b平面的投影与x_b轴的夹角记为β_T，称之为推力横向偏转角，投影在对称平面右侧为正。推力变向无人机的非线性全量模型由动力学方程组和运动学方程组组成。动力学方程组包含质心移动的动力学方程组(又称力方程组)和绕质心转动的动力学方程组(又称力矩方程组)。无人机的运动学方程组包括质心移动的运动学方程组(又称导航方程组)和绕质心转动的运动学方程组(又称运动方程组)。描述推力变向无人机的数学模型需要明确输入输出关系，为此建立如图1所示的模型结构框图，图中虚线框描述的就是推力变向无人机的动力学和运动学模型。推力变向无人机系统的控制输入信号包括升降舵偏转角δ_e、副翼舵偏转角δ_a、方向舵偏转角δ_r、油门大小δ_T、推力偏转角α_T、β_T，状态或输出信号包括飞行速度V、迎角α、侧滑角β、航迹倾斜角μ、航迹方位角空间位置坐标x，y，z、空间姿态角θ，φ，ψ及其角速率p，q，r，飞机的合外力由气动力A_x，A_y，A_z、重力G_x，G_y，G_z和推力T_x，T_y，T_z组成，合外力矩由气动力矩M_A，N_A和推力力矩M_T，N_T组成。

基于推力变向无人机的受力和力矩分析，容易获得其非线性动力学和运动学方程组。与常规无人机不同，变推力无人机的控制输入信号增加了发动机推力的两个偏转角，使飞机y_b轴和z_b轴向的推力分量不再是零，因此无人机的合外力发生变化，由于发动机的安装位置不在质心处，因此也引起无人机的合外力矩发生改变。

推力变向无人机的数学模型与常规无人机相比，具有以下特点：

(1)对于常规无人机，合外力和外力矩计算式中的推力偏转角为α_T＝β_T＝0，即推力在机体坐标系y_b轴和z_b轴向没有分量。由于推力变向无人机的合外力和外力矩与常规无人机不同，导致动力学特性也不同。在力方程组中，推力变向无人机增加了推力在y_b轴和z_b轴向的分量所引起的部分。在力矩方程组中，推力变向无人机增加了推力在y_b轴和z_b轴向的分量引起的力矩项。

(2)无人机质心移动的运动学方程描述了飞机在空间位置的变化规律，绕质心转动的运动学方程描述了飞机空间姿态的变化规律，它与飞机旋转角速度密切相关。推力变向无人机的运动学方程组的描述形式与常规无人机相同，但是两者的各个状态物理量的动力学特性不同，导致运动学特性也会发生改变。

下面根据引入推力变向技术后无人机合外力和外力矩的变化情况，建立推力变向无人机的动力学方程组。

根据推力变向无人机的发动机推力的两个偏转角在机体坐标系中的几何关系，可得推力T在机体坐标轴系的投影为

$\left(\begin{array}{c}{T}_x=\frac{T}{\sqrt{1+{\tan }^2{\alpha }_T+{\tan }^2{\beta }_T}}\\ T_y=T_x\tan {\beta }_T=\frac{T}{\sqrt{1+{\tan }^2{\alpha }_T+{\tan }^2{\beta }_T}}\tan {\beta }_T\\ T_z=T_x\tan {\alpha }_T=\frac{T}{\sqrt{1+{\tan }^2{\alpha }_T+{\tan }^2{\beta }_T}}\tan {\alpha }_T\end{array}\right)---\left(1\right)$

对于常规无人机，由于α_T＝β_T＝0，因此T_x＝T，T_y＝T_z＝0。

推力变向无人机的合外力由发动机推力、气动力和重力构成。根据机体坐标系中飞行器在合外力作用下的线运动方程，可以推导得到力方程组

式中，D、Y、L分别表示空气阻力、侧力和升力，m为飞机质量，G_xa，G_ya，G_za分别表示飞机重力在气流坐标系中的三轴分量。

在研究无人机飞行轨迹特性时，常采用航迹坐标系上投影的质心动力学方程，此时无人机所受的推力、气动力和重力均需投影到航迹坐标系中。根据飞行器在合外力作用下的线运动方程，可得力方程组为

式中，γ为航迹滚转角。

由于常规无人机的发动机推力始终沿着轴线方向保持不变，因此，推力变向无人机的质心移动方程组与常规无人机的区别之处在于，常规无人机的力方程组中推力沿机体坐标系的三个分量分别为T_x＝T，T_y＝T_z＝0，将其代入上述力方程组，即可获得常规无人机的力方程组。

设推力变向无人机的发动机推力作用点在机体坐标系的坐标为(l_x，l_y，l_z)，由于推力作用点位于飞机的对称面Oxz上，因此l_y＝0，于是发动机推力力矩的三轴分量表示为

$\left(\begin{array}{c}{\bar{L}}_T\\ M_T\\ N_T\end{array}\right)=\left(\begin{array}{c}-T_y{l}_z+T_z{l}_y\\ T_x{l}_z-T_z{l}_x\\ -T_x{l}_y+T_y{l}_x\end{array}\right)=\left(\begin{array}{c}-T_y{l}_z\\ T_x{l}_z-T_z{l}_x\\ T_y{l}_x\end{array}\right)---\left(4\right)$

式中，若l_z＜0，l_x＜0，当β_T＜0时，有T_y＜0，N_T＞0，即推力产生的滚转力矩为负，推力产生的偏航力矩为正，当α_T＜0时，有M_T＜0，即推力产生低头力矩。对于常规无人机，则有M_T＝Tl_z，

若由动量矩定理表示质心的转动运动，可得变向无人机绕质心转动的动力学方程组为

由于推力变向无人机与常规无人机的合外力力矩不同，使得绕质心转动的动力学特性发生改变。

飞行器运动方程组的线性化一般是在小扰动假设条件下进行的，在该条件下，纵向气动力和力矩对横侧参数在基准运动状态下的导数为零，横侧向气动力和力矩对纵向参数在基准运动下的导数也是零。由于推力变向技术能够补偿飞行操纵效率，因此在对推力变向无人机数学模型进行线性化处理时，有必要考虑发动机推力的两个偏转角。在定直平飞基准运动状态下，推力变向无人机的飞行状态满足：(3)α_T^*＝β_T^*＝0。本节中变量符号带上标“*”表示的是基准运动状态下的变量。

利用小扰动线性化方法，求出发动机推力的动力学偏导数、阻力、升力和侧力的偏导数、纵向力矩和横向力矩的偏导数，基于上述对推力、空气动力和合外力力矩的偏量分析，对推力变向无人机的动力学和运动学方程组进行线性化。求出推力变向无人机的阻力线性化方程、推力变向无人机的升力线性化方程、推力变向无人机的侧力线性化方程，求出机体坐标系下推力变向无人机的力矩线性化方程组。导航方程组和运动方程组的偏量表达式、线性化几何关系方程可参考有关文献，与常规无人机相同。

推力变向无人机的纵向耦合方程组包括ΔV，Δα，Δθ，Δq，Δδ_e，Δδ_T，Δα_T等七个变量。基于上述线性表达式，可得推力变向无人机纵向线性化运动方程为：

式中，

X_lon＝[ΔV，Δα，Δθ，Δq]^T

U_lon＝[Δδ_T，Δδ_e，Δα_T]^T

横侧向耦合的方程组包含Δβ，Δφ，Δp，Δr，Δδ_a，Δδ_r，Δβ_T等七个变量。基于上述线性表达式，可得推力变向无人机的横侧向线性化运动方程为：

式中，

X_lat＝[Δβ，Δφ，Δp，Δr]^T

U_lat＝[Δδ_a，Δδ_r，Δβ_T]^T

$A_{\mathrm{lat}}=\left(\begin{array}{cccc}\frac{Y_{\beta }}{{\mathrm{mV}}^{*}}& \frac{{g\cos \theta }^{*}}{V^{*}}& {\sin \alpha }^{*}+\frac{Y_p}{{\mathrm{mV}}^{*}}& \frac{Y_r}{{\mathrm{mV}}^{*}}-{\cos \alpha }^{*}\\ 0& 0& 1& {\tan \theta }^{*}\\ \frac{I_z{\bar{L}}_{\beta }+I_{\mathrm{zx}}{N}_{\beta }}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& 0& \frac{I_z{\bar{L}}_p+I_{\mathrm{zx}}{N}_p}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_z{\bar{L}}_r+I_{\mathrm{zx}}{N}_r}{I_x{I}_z-{I_{\mathrm{zx}}}^2}\\ \frac{I_{\mathrm{zx}}{\bar{L}}_{\beta }+I_x{N}_{\beta }}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& 0& \frac{I_{\mathrm{zx}}{\bar{L}}_p+I_x{N}_p}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_{\mathrm{zx}}{\bar{L}}_r+I_x{N}_r}{I_x{I}_z-{I_{\mathrm{zx}}}^2}\end{array}\right),$

$B_{\mathrm{lat}}=\left(\begin{array}{ccc}\frac{Y_{\mathrm{\delta a}}}{{\mathrm{mV}}^{*}}& \frac{Y_{\mathrm{\delta r}}}{{\mathrm{mV}}^{*}}& \frac{T^{*}}{{\mathrm{mV}}^{*}}\\ 0& 0& 0\\ \frac{I_z{\bar{L}}_{\mathrm{\delta a}}+I_{\mathrm{zx}}{N}_{\mathrm{\delta a}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_z{\bar{L}}_{\mathrm{\delta r}}+I_{\mathrm{zx}}{N}_{\mathrm{\delta r}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_z{\bar{L}}_{\mathrm{\beta T}}+I_{\mathrm{zx}}{N}_{\mathrm{\beta T}}}{I_z{I}_x-{I_{\mathrm{zx}}}^2}\\ \frac{I_{\mathrm{zx}}{\bar{L}}_{\mathrm{\delta a}}+I_x{N}_{\mathrm{\delta a}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_{\mathrm{zx}}{\bar{L}}_{\mathrm{\delta r}}+I_x{N}_{\mathrm{\delta r}}}{I_x{I}_z-{I_{\mathrm{zx}}}^2}& \frac{I_{\mathrm{zx}}{\bar{L}}_{\mathrm{\beta T}}+I_x{N}_{\mathrm{\beta T}}}{I_z{I}_x-{I_{\mathrm{zx}}}^2}\end{array}\right).$

为了便于公众了解本发明的技术方案，用下表1将本发明涉及的物理量进行说明：

表1

上述高阶线性常微分方程组描述了推力变向无人机偏离基准运动后纵横向运动的偏离情况，有关变量的量纲导数可参考文献。与常规无人机不同，推力变向无人机的纵向线性状态方程中增加了控制量Δα_T，横向线性状态方程中增加了控制量Δβ_T，相应的控制输入矩阵也增加了对应的分量。

推力变向无人机的姿态控制问题中，飞行纵向姿态的稳定和控制就是俯仰角的稳定和控制，通过控制升降舵和推力纵向偏角来实现，飞行横侧向姿态的稳定和控制就是要保证高精度的偏航角和滚转角的稳定和控制，通过方向舵、副翼舵以及推力横向偏转角操纵来实现。由于控制器设计方法是基于离散状态方程的，因此有必要对无人机纵横向线性化方程转化为离散状态方程。

对上述推力变向无人机的纵向、横侧向线性化模型，即式(6)、(7)所示，进行离散化，然后设计

首先设计纵向姿态的融合控制器，建立纵向姿态控制的最优性能指标模型，如下式所示：

$J={X_{\mathrm{lon}}}^T\left(N\right)Q_{\mathrm{lon}}\left(N\right)X_{\mathrm{lon}}\left(N\right)+$

$\underset{k=0}{\overset{N-1}{\Sigma }}\left[{[{X_{\mathrm{lon}}}^{*}\left(k\right)-X_{\mathrm{lon}}\left(k\right)]}^T{Q}_{\mathrm{lon}}\left(k\right)\right[{X_{\mathrm{lon}}}^{*}\left(k\right)-X_{\mathrm{lon}}\left(k\right)]+{[{U_{\mathrm{lon}}}^{*}\left(k\right)-U_{\mathrm{lon}}\left(k\right)]}^T{R}_{\mathrm{lon}}\left(k\right)[{U_{\mathrm{lon}}}^{*}\left(k\right)-U_{\mathrm{lon}}\left(k\right)\left]\right]---\left(7\right)$

式中，X_lon^*为期望姿态信号，U_lon^*(k)为期望控制信号，Q_lon＝Q_lon^T＞0、R_lon＝R_lon^T＞0为权重矩阵，k为时间序列，N为终端时刻。根据信息融合控制理论与方法，建立如下信息表达式

$\left(\begin{array}{c}{X_{\mathrm{lon}}}^{*}\left(k\right)=X_{\mathrm{lon}}\left(k\right)+m_{\mathrm{lon}}\left(k\right)\\ {U_{\mathrm{lon}}}^{*}\left(k\right)=U_{\mathrm{lon}}\left(k\right)+n_{\mathrm{lon}}\left(k\right)\\ {\hat{X}}_{\mathrm{lon}}(j+1)=X_{\mathrm{lon}}(j+1)+w_{\mathrm{lon}}(j+1),j=k+1,k+2,L,N\\ X_{\mathrm{lon}}(k+1)=A_{\mathrm{lon}}{X}_{\mathrm{lon}}\left(k\right)+B_{\mathrm{lon}}{U}_{\mathrm{lon}}\left(k\right)\end{array}\right)---\left(8\right)$

式中，m_lon和n_lon为白噪声序列，E[m_lon(k+1)]＝0，var[m_lon(k+1)]＝Q_lon^-1(k+1)，E[n_lon(k)]＝0，var[n_lon(k)]＝R_lon^-1(k)。基于信息融合估计定理，得到纵向姿态的融合控制律为

${\hat{U}}_{\mathrm{lon}}\left(k\right)={[R_{\mathrm{lon}}\left(k\right)+{R_{\mathrm{lon}}}^T\left(k\right){R_{\mathrm{lon}}}^{-1}(k+1)B_{\mathrm{lon}}\left(k\right)]}^{-1}\{{B_{\mathrm{lon}}}^T\left(k\right){P_{\mathrm{lon}}}^{-1}(k+1)$(9)

$[{\hat{X}}_{\mathrm{lon}}(k+1)-A_{\mathrm{lon}}\left(k\right)X_{\mathrm{lon}}\left(k\right)]+R_{\mathrm{lon}}\left(k\right){U_{\mathrm{lon}}}^{*}\left(k\right)\}$

式中，

$\left(\begin{array}{c}{\hat{X}}_{\mathrm{lon}}\left(j\right)=P_{\mathrm{lon}}\left(j\right)\left\{{A_{\mathrm{lon}}}^T\left(j\right){[P_{\mathrm{lon}}(j+1)+B_{\mathrm{lon}}\left(j\right){R_{\mathrm{lon}}}^{-1}\left(j\right){B_{\mathrm{lon}}}^T\left(j\right)]}^{-1}\right[{\hat{X}}_{\mathrm{lon}}(j+1)-B_{\mathrm{lon}}\left(j\right){U_{\mathrm{lon}}}^{*}\left(j\right)]+Q_{\mathrm{lon}}\left(j\right){X_{\mathrm{lon}}}^{*}\left(j\right)\}\\ {P_{\mathrm{lon}}}^{-1}\left(j\right)=Q_{\mathrm{lon}}\left(j\right)+{A_{\mathrm{lon}}}^T\left(j\right){[P_{\mathrm{lon}}(j+1)+B_{\mathrm{lon}}\left(j\right){R_{\mathrm{lon}}}^{-1}\left(j\right){B_{\mathrm{lon}}}^T\left(j\right)]}^{-1}{A}_{\mathrm{lon}}\left(j\right),{P_{\mathrm{lon}}}^{-1}\left(N\right)=Q_{\mathrm{lon}}\left(N\right)\\ j=k+1,k+2,L,N-1\end{array}\right)---\left(10\right)$

类似地，可得到推力变向无人机的横侧向姿态的融合控制律为

${\hat{U}}_{\mathrm{lat}}\left(k\right)={[R_{\mathrm{lat}}\left(k\right)+{B_{\mathrm{lat}}}^T\left(k\right){P_{\mathrm{lat}}}^{-1}(k+1)B_{\mathrm{lat}}\left(k\right)]}^{-1}\{{B_{\mathrm{lat}}}^T\left(k\right){P_{\mathrm{lat}}}^{-1}(k+1)$(11)

$[{\hat{X}}_{\mathrm{lat}}(k+1)-A_{\mathrm{lat}}\left(k\right)X_{\mathrm{lat}}\left(k\right)]+R_{\mathrm{lat}}\left(k\right){U_{\mathrm{lat}}}^{*}\left(k\right)\}$

式中，

$\left(\begin{array}{c}{\hat{X}}_{\mathrm{lat}}\left(j\right)=P_{\mathrm{lat}}\left(j\right)\left\{{A_{\mathrm{lat}}}^T\left(j\right){[P_{\mathrm{lat}}(j+1)+B_{\mathrm{lat}}\left(j\right){R_{\mathrm{lat}}}^{-1}\left(j\right){B_{\mathrm{lat}}}^T\left(j\right)]}^{-1}\right[{\hat{X}}_{\mathrm{lat}}(j+1)-B_{\mathrm{lat}}\left(j\right){U_{\mathrm{lat}}}^{*}\left(j\right)]+Q_{\mathrm{lat}}\left(j\right){X_{\mathrm{lat}}}^{*}\left(j\right)\}\\ {P_{\mathrm{lat}}}^{-1}\left(j\right)=Q_{\mathrm{lat}}\left(j\right)+{A_{\mathrm{lat}}}^T\left(j\right){[P_{\mathrm{lat}}(j+1)+B_{\mathrm{lat}}\left(j\right){R_{\mathrm{lat}}}^{-1}\left(j\right){B_{\mathrm{lat}}}^T\left(j\right)]}^{-1}{A}_{\mathrm{lat}}\left(j\right),{P_{\mathrm{lat}}}^{-1}\left(N\right)=Q_{\mathrm{lat}}\left(N\right)\\ j=k+1,k+2,L,N-1\end{array}\right)---\left(12\right)$

基于上述纵向、横侧向姿态的融合控制律，可得到本发明的推力变向无人机的姿态控制系统，如图2所示。

采用本发明的姿态控制系统进行推力变向无人机姿态控制，按照以下方法：

步骤1、姿态控制系统接收给定姿态指令信号、传感器反馈姿态信号以及期望控制信号，所述给定姿态指令信号包括，纵向期望姿态信号X_lon^*＝[ΔV^*，Δα^*，Δθ^*，Δq^*]^T，横向期望姿态信号X_lat^*＝[Δβ^*，Δφ^*，Δp^*，Δr^*]^T，所述传感器反馈姿态信号包括，纵向实际姿态信号X_lon＝[ΔV，Δα，Δθ，Δq]^T，横向实际姿态信号X_lat＝[Δβ，Δφ，Δp，Δr]^T，所述期望控制信号包括，纵向期望控制信号U_lon^*＝[Δδ_T^*，Δδ_e^*，Δα_T^*]^T，横向期望控制信号U_lat^*＝[Δδ_a^*，Δδ_r^*，Δβ_T^*]^T；

步骤2、姿态控制系统根据接收到的给定姿态指令信号解算姿态控制律，得到纵向姿态控制信号U_lon＝[Δδ_T，Δδ_e，Δα_T]^T，横向姿态控制信号U_lat＝[Δδ_a，Δδ_r，Δβ_T]^T；

步骤3、姿态控制系统将姿态控制信号发送给执行机构，执行机构操纵无人机的气动舵面、油门开度、推力方向，实现控制无人机飞行姿态和飞行速度的目的。

为了验证本发明在推力变向无人机姿态控制中的有效性和先进性，对受到气流扰动下的某型推力方向可变的螺旋桨无人机纵横向姿态的稳定过程进行仿真研究，控制系统的结构图如图2所示。受控对象是全量非线性方程组描述的推力变向无人机，其气动参数通过风洞实验获得，动态方程采用龙格-库塔法进行计算，仿真工具采用Matlab软件。考虑无人机在高度h＝500m、飞行速度V＝50m/s，处于定直平飞平衡状态。假设无人机从初始时刻开始4s内持续受到阵风干扰，气流方向为θ_aero＝45°，扰动气流速度为V_aero＝20m/s，说明无人机的俯仰、滚转和偏航姿态均受到强干扰。

采用常规单回路设计方法的姿态控制仿真结果如图3所示；采用本发明姿态控制系统的仿真结果如图4所示。其中虚线表示纯气动舵面控制方式下的姿态响应曲线，实线表示增加了推力变向控制的姿态响应曲线。

从仿真结果可以得出以下几点结论：

(1)推力变向技术能够提高飞机操纵效率，减小气动舵面的偏转量，一定程度上能够避免气动舵面饱和的不利情况，充分说明了推力变向控制对气动舵面有补偿作用；

(2)气动舵面加推力变向控制下的姿态稳定过程要比纯气动舵面控制下的姿态稳定过程更加快速，推力变向无人机受阵风扰动的影响要比常规无人机的更小。

(3)信息融合控制方法对受阵风干扰下的姿态稳定效果要优于线性单回路设计方法，从而验证了融合最优控制方法作为多变量控制方法的优越性。