CN103838145A - 基于级联观测器的垂直起降飞机鲁棒容错控制系统及方法 - Google Patents

Abstract

本发明公开了基于级联观测器的垂直起降飞机鲁棒容错控制系统及方法，属于飞行器控制的技术领域，建立考虑了外部干扰以及执行器故障的垂直起降机动力学模型，转换动力学模型为标准形式的非线性模型，利用级联的干扰观测器、故障诊断观测器实时估计干扰、故障以及状态量，将非线性模型解耦为两个子系统并确定各子系统的控制率，再由两个子系统的控制率、故障估计值确定容错控制率。本发明考虑了外部干扰和故障作用的系统在容错控制器的作用下，实时补偿故障和干扰，迅速地减少了故障以及干扰的影响，提高了控制系统的鲁棒性。

Description

基于级联观测器的垂直起降飞机鲁棒容错控制系统及方法

技术领域

本发明公开了基于级联观测器的垂直起降飞机鲁棒容错控制系统及方法，属于飞行器控制的技术领域。

背景技术

垂直起降飞机是一种既具有直升飞机垂直起落和悬停的优点，又具有固定翼飞机高速平飞和机动的优点的一种新型飞行器。它是一种具有三个自由度、两个控制输入的典型的欠驱动非线性系统。近年来，由于其巨大的应用价值，尤其是军事应用价值，受到了许多国家相关学者的关注和研究。

由于垂直起降飞机具有欠驱动、强耦合和非线性的特性，在实际的飞行中，飞机的姿态不可避免地会受到周围环境中各种不确定干扰的影响，这一类扰动通常会与姿态角产生耦合进而影响到飞行的稳定性；同时，飞机的主要升力来源是推力发动机，发动机和相应执行机构的结构较为发杂，容易出现故障，导致飞机的飞行性能受损甚至出现不可预想的事故。因此，对干扰和故障进行实时检测，利用相应的冗余机制设计飞机的容错控制系统，具有现实意义和重要的工程应用价值。

近年来，滑模控制技术因其所具有的优良特性而受到越来越多的重视，该方法对参数变化和扰动不敏感，结构简单，适用于垂直起降飞机的飞行控制。针对飞行中存在的干扰影响，有学者设计了相应的干扰观测器抵消干扰的影响，改善滑模控制的抖振问题，实现高精度的轨迹跟踪控制和姿态控制，但是都没有考虑到同时存在的执行器故障问题。

发明内容

本发明所要解决的技术问题是针对上述背景技术的不足，提供了基于级联观测器的垂直起降飞机鲁棒容错控制系统及方法。

本发明为实现上述发明目的采用如下技术方案：

基于级联观测器的垂直起降飞机鲁棒容错控制系统，是串接在垂直起降飞机模型输入端、输出端之间闭环控制系统，具体包括：执行机构、观测器单元、容错控制器，

所述执行机构，在容错控制信号u_f以及故障f_a的作用下得到垂直起降飞机的操作指令；

所述垂直起降飞机模型在外部干扰f_d以及执行机构操作指令的作用下动作，以飞机质心位置(x,y)以及滚转角θ构成所述控制系统的输出向量Y=(x,y,θ)^T；

所述观测器单元包括级联的干扰观测器、故障诊断观测器，干扰观测器根据滚转角θ来得到干扰估计值

以及相应的状态变量观测值，故障诊断观测器根据干扰估计值

飞机质心位置(x,y)、容错控制信号u_f得到故障估计值

以及相应的状态变量观测值，x、y分别为飞机质心水平坐标和垂直坐标，所述状态变量观测值包括：滚转角和滚转速率的观测值

和飞机质心位置观测值以及相应的速率观测值

所述容错控制器，以跟踪误差信号、滚转角和滚转速率的观测值、故障估计值干扰估计值

为输入，得到容错控制信号修正值，所述跟踪误差信号为由飞机质心位置以及相应速率的观测值组成的观测信号与给定参考值

的差值，x_d(t),y_d(t)为飞机质心位移，

为飞机质心移动速率。

基于级联观测器的垂直起降飞机鲁棒容错控制方法，构建所述的控制系统，具体包括如下步骤：

步骤1，建立考虑了外部干扰以及执行器故障的垂直起降机动力学模型：

$\left\{\begin{array}{c}\stackrel{\·\·}{x}=-(u_{1f}-f_a\left(t\right))\sin \mathrm{\θ}+{\mathrm{\ϵu}}_{2f}\cos \mathrm{\θ}+d_1\left(\mathrm{\θ}\right)f_d\left(t\right)\\ \stackrel{\·\·}{y}=(u_{1f}-f_a\left(t\right))\cos \mathrm{\θ}+{\mathrm{\ϵu}}_{2f}\sin \mathrm{\θ}+d_2\left(\mathrm{\θ}\right)f_d\left(t\right)-g\\ \stackrel{\·\·}{\mathrm{\θ}}={\mathrm{\λu}}_{2f}+d_3\left(\mathrm{\θ}\right)f_d\left(t\right)\end{array}\right.---\left(1\right),$

式（1）中，

分别为飞机在水平方向、垂直方向的加速度，是滚转角加速度，u_1f、u_2f分别为控制系统的推力容错控制输入、滚转力矩容错控制输入，ε为滚转力矩与水平方向加速度的耦合系数，ε≠0，f_d(t)、f_a(t)分别为外部干扰、故障的时间函数，d_i(θ),i=1,2,3为外部干扰和滚转角的耦合程度，g为重力加速度，归一化后取g=1，λ为滚转力矩的传输效率，0<λ≤1；

步骤2，选取状态变量向量x，x=(x₁,x₂,x₃,x₄,x₅,x₆)， $x_1=x,x_2=\stackrel{\&CenterDot;}{x},x_3=y,x_4=\stackrel{\&CenterDot;}{y},x_5=\mathrm{\&theta;},x_6=\stackrel{\&CenterDot;}{\mathrm{\&theta;}},$ 利用变量替换将步骤1建立的动力学模型转换为标准非线性模型：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_1=x_2\\ {\stackrel{\&CenterDot;}{x}}_2={-u}_{1f}\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+d_1\left(x_5\right)f_d+\sin x_5{f}_a\\ {\stackrel{\&CenterDot;}{x}}_3=x_4\\ {\stackrel{\&CenterDot;}{x}}_4=u_{1f}\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+d_2\left(x_5\right)f_d-\cos x_5{f}_a-1\\ {\stackrel{\&CenterDot;}{x}}_5=x_6\\ {\stackrel{\&CenterDot;}{x}}_6={\mathrm{\&lambda;u}}_{2f}+d_3\left(x_5\right)f_d\end{array}\right.---\left(2\right),$

选取第一至第三外部干扰变量f_1d=d₁(x₅)f_d,f_2d=d₂(x₅)f_d,f_3d=d₃(x₅)f_d以及第一、第二故障变量f_1a=sinx₅f_a,f_2a=-cosx₅f_a，确定故障系统的状态方程为：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_1=x_2\\ {\stackrel{\&CenterDot;}{x}}_2={-u}_{1f}\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+f_{1d}+f_{1a}\\ {\stackrel{\&CenterDot;}{x}}_3=x_4\\ {\stackrel{\&CenterDot;}{x}}_4=u_{1f}\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+f_{2d}-f_{2a}-1\\ {\stackrel{\&CenterDot;}{x}}_5=x_6\\ {\stackrel{\&CenterDot;}{x}}_6={\mathrm{\&lambda;u}}_{2f}+f_{3d}\end{array}\right.---\left(3\right),$

式（2）、（3）中：输出向量为Y=(y₁,y₂,y₃)^T=(x₁,x₃,x₅)^T，用于轨迹跟踪的输出为y₁和y₂；

步骤3，利用干扰观测器：

$\left\{\begin{array}{c}{\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_5={\hat{x}}_6+{\mathrm{\&lambda;}}_1\mathrm{sgn}(x_5-{\hat{x}}_5)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_6={\mathrm{\&lambda;u}}_{2f}+{\mathrm{\&lambda;}}_2{\mathrm{sgn}}_1({\tilde{x}}_6-{\hat{x}}_6)\end{array}\right.---\left(4\right),$

得到第三干扰变量估计值

进而变换得到干扰估计值式（4）中：

分别为滚转角、滚转速率的观测值，λ₁和λ₂是切换项的增益，

是滚转速率x₆的替代变量，sgn(·)、sgn₁(·)为符号函数，

${\tilde{x}}_6={\hat{x}}_6+{\left({\mathrm{\&lambda;}}_1\mathrm{sgn}(x_5-{\hat{x}}_5)\right)}_{\mathrm{eq}}---\left(5\right),$

$\mathrm{sgn}(\&CenterDot;)=\left\{\begin{array}{cc}1& x_5-{\hat{x}}_5>0\\ 0& x_5-{\hat{x}}_5=0\\ -1& x_5-{\hat{x}}_5<0\end{array}\right.---\left(6\right),$

${\mathrm{sgn}}_1(\&CenterDot;)=\left\{\begin{array}{cc}\mathrm{sgn}(\&CenterDot;)& x_5-{\hat{x}}_5=0\\ 0& x_5-{\hat{x}}_5\&NotEqual;0\end{array}\right.---\left(7\right),$

式（5）中，(·)_eq为等效的误差信号；

利用故障诊断观测器

$\left\{\begin{array}{c}{\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_1={\hat{x}}_2-l_1({\hat{x}}_1-x_1)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_2={-u}_{1f}\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+{\hat{f}}_{1a}+{\hat{f}}_{1d}-({\hat{x}}_1-x_1)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_3={\hat{x}}_4-l_2({\hat{x}}_3-x_3)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_4=u_{1f}\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+{\hat{f}}_{2a}+{\hat{f}}_{2d}-({\hat{x}}_3-x_3)-1\end{array}\right.---\left(8\right),$

得到故障估计值

式（8）中：

为飞机质心水平坐标、垂直坐标、水平方向运动速率、垂直方向运动速率的观测值，为第一、第二故障变量的估计值，

为第一、第二外部干扰变量的估计值，l₁和l₂为故障诊断观测器的增益；

步骤4，取六个误差变量e₁,e₂,e₃,e₄,η₁,η₂：e₁=x₁-x_d，

e₃=x₃-y_d， $e_4=x_4-{\stackrel{\&CenterDot;}{y}}_d,$ η₁=x₅， ${\mathrm{\&eta;}}_2=\frac{\mathrm{\&epsiv;}}{\mathrm{\&lambda;}}{x}_6-e_2\cos x_5-e_4\sin x_5,$ 利用坐标变换将控制系统坐标变为跟踪系统：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{e}}_1=e_2\\ {\stackrel{\&CenterDot;}{e}}_2=w_{s1}+f_{1d}\\ {\stackrel{\&CenterDot;}{e}}_3=e_4\\ {\stackrel{\&CenterDot;}{e}}_4=w_{s2}+f_{2d}\\ {\stackrel{\&CenterDot;}{\mathrm{\&eta;}}}_1=\frac{\mathrm{\&lambda;}}{\mathrm{\&epsiv;}}({\mathrm{\&eta;}}_2+e_2\cos {\mathrm{\&eta;}}_1+e_4\sin {\mathrm{\&eta;}}_1)\\ {\stackrel{\&CenterDot;}{\mathrm{\&eta;}}}_2=\frac{\mathrm{\&lambda;}}{\mathrm{\&epsiv;}}({\mathrm{\&eta;}}_2+e_2\cos {\mathrm{\&eta;}}_1+e_4\sin {\mathrm{\&eta;}}_1)(e_2\sin {\mathrm{\&eta;}}_1-e_4\cos {\mathrm{\&eta;}}_1)+{\stackrel{\&CenterDot;\&CenterDot;}{y}}_d\sin {\mathrm{\&eta;}}_1+\frac{\mathrm{\&lambda;}}{\mathrm{\&epsiv;}}{f}_{3d}\\ +{\stackrel{\&CenterDot;\&CenterDot;}{x}}_d\cos {\mathrm{\&eta;}}_1+\sin {\mathrm{\&eta;}}_1-f_{1d}\cos {\mathrm{\&eta;}}_1-f_{2d}\sin {\mathrm{\&eta;}}_1\end{array}\right.---\left(9\right),$

控制输入的可逆变换为：

$\left(\begin{array}{c}{w}_{s1}\\ w_{s2}\end{array}\right)=\mathrm{\&beta;}\left(x_5\right)\left(\begin{array}{c}{u}_{1f}-f_a\\ u_{2f}\end{array}\right)-\left(\begin{array}{c}{\stackrel{\&CenterDot;\&CenterDot;}{x}}_d\\ {\stackrel{\&CenterDot;\&CenterDot;}{y}}_d+1\end{array}\right)---\left(10\right),$

其中： $\mathrm{\&beta;}\left(x_5\right)=\left(\begin{array}{cc}-\sin x_5& \mathrm{\&epsiv;}\cos x_5\\ \cos x_5& \mathrm{\&epsiv;}\sin x_5\end{array}\right)$ 是系统在平衡点x₅=0附近的非奇异解耦矩阵。

利用线性相关法将跟踪系统解耦为最小相位子系统：

以及非最小相位子系统：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{e}}_1=e_2\\ {\stackrel{\&CenterDot;}{e}}_2=w_{s1}+f_{1d}\\ \stackrel{\&CenterDot;}{\mathrm{\&eta;}}=\mathrm{\&Gamma;}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)\end{array}\right.---\left(12\right),$

式（12）中：η=(η₁,η₂)^T，e=(e₁,e₂,e₃,e₄)^T，

为位置参考轨迹的二阶导数。

选择新的误差变量：z₁=(e₁,η_e ^T)^T，z₂=e₂，记 $\mathrm{\&psi;}({\mathrm{\&eta;}}_e,e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)=\left(\begin{array}{c}{e}_2\\ \mathrm{\&Gamma;}({\mathrm{\&eta;}}_e,e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)\end{array}\right),$ 则非最小相位子系统变为：

式（13）中， $A_2=\frac{\&PartialD;\mathrm{\&psi;}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)}{\&PartialD;(e_1,{\mathrm{\&eta;}}^T)}{|}_o,$ $B_2=\frac{\&PartialD;\mathrm{\&psi;}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)}{\&PartialD;\left(e_2\right)}{|}_o,$ $\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)=\mathrm{\&psi;}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)-A_2{z}_1-B_2{z}_2$ 是线性化剩余的高阶项，且有 $\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)=\left(\begin{array}{c}0\\ \frac{\mathrm{\&lambda;}}{\mathrm{\&epsiv;}}(e_2\cos {\mathrm{\&eta;}}_1+e_4\sin {\mathrm{\&eta;}}_1-e_2)\\ \frac{\mathrm{\&lambda;}}{\mathrm{\&epsiv;}}({\mathrm{\&eta;}}_2+e_2\cos {\mathrm{\&eta;}}_1+e_4\sin {\mathrm{\&eta;}}_1)(e_2\sin {\mathrm{\&eta;}}_1-e_4\cos {\mathrm{\&eta;}}_1)-{\mathrm{\&eta;}}_1\\ +{\stackrel{\&CenterDot;\&CenterDot;}{y}}_d\sin {\mathrm{\&eta;}}_1+{\stackrel{\&CenterDot;\&CenterDot;}{x}}_d\cos {\mathrm{\&eta;}}_1+\sin {\mathrm{\&eta;}}_1+\frac{\mathrm{\&lambda;}}{\mathrm{\&epsiv;}}{f}_{3d}-f_{1d}\cos {\mathrm{\&eta;}}_1-f_{2d}\sin {\mathrm{\&eta;}}_1\end{array}\right).$

确定镇定最小相位子系统的最优控制律w_s2为：

$w_{s2}={-\mathrm{\&beta;}}_1{\hat{e}}_3-{\mathrm{\&beta;}}_2{\hat{e}}_4-{\hat{f}}_{2d}---\left(14\right),$

式（14）中，

是f_2d的估计值，可由干扰观测器获得的

进行变换获得，

是最优控制增益，R₁和P₁是满足相应Riccati方程的正定矩阵，在最优控制律的作用下，最小相位子系统显然能够指数渐进稳定；

针对线性化后的非最小相位子系统，选择滑模平面为：

s=e₂+c₁e₁+c₂η₁+c₃η₂                 （15），

式（15）中，C=(-c₁ -c₂ -c₃)为使得A₂+B₂C为Hurwitz矩阵的向量，是f_1d的估计值，

是滑模平面s的估计值，σ是滑动切换项的增益，饱和函数sat(i)用来代替符号函数，消除滑模控制律的抖振现象，

表达式为：

$\mathrm{sat}\left(\hat{s}\right)=\left\{\begin{array}{cc}\mathrm{sgn}\left(\hat{s}\right)& \left|\hat{s}\right|>\mathrm{\&gamma;}\\ \frac{\hat{s}}{\mathrm{\&gamma;}}& \left|\hat{s}\right|\&le;\mathrm{\&gamma;}\end{array}\right.---\left(16\right),$

式（16）中，γ为消除抖振的因子，同时也能保证系统快速收敛，一般取在(0.02,0.08)范围内。

相应的滑模控制律w_s1为

$w_{s1}=C(A_2{\hat{z}}_1+B_2{\hat{z}}_2+\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\widehat{\mathrm{\&eta;}},\hat{e},{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d))-{\hat{f}}_{1d}-\mathrm{\&sigma;sat}\left(\hat{s}\right)---\left(17\right),$

步骤5，根据步骤4确定的最优控制律以及滑模控制律，结合故障估计值和可逆的控制变换，确定整个容错控制器的控制律：

$u_f={\mathrm{\&beta;}}^{-1}\left(x_5\right)\left(\begin{array}{c}{w}_{s1}+{\stackrel{\&CenterDot;\&CenterDot;}{x}}_d\\ w_{s2}+{\stackrel{\&CenterDot;\&CenterDot;}{y}}_d+1\end{array}\right)+\hat{F}---\left(18\right),$

式（18）中：β^-1(x₅)是可逆解耦矩阵β(x₅)的逆矩阵，

为重构故障向量的估计值。

本发明采用上述技术方案，具有以下有益效果：考虑了外部干扰和故障作用的系统在容错控制器的作用下，实时补偿故障和干扰，迅速地减少了故障以及干扰的影响，提高了控制系统的鲁棒性。

附图说明

图1(a)是垂直起降飞机模型的示意图。

图1(b)是垂直起降飞机容错控制系统的结构框图。

图2(a)、图2(b)是干扰观测器的干扰估计和状态估计曲线。

图3(a)、图3(b)、图3(c)是故障诊断观测器的故障估计和状态估计曲线。

图4是容错控制器的曲线。

图5(a)、图5(b)是容错控制器作用下的系统输出跟踪和滚转姿态（零动态）响应曲线。

具体实施方式

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

飞机的受力分析如图1(a)所示，考虑推力输入通道执行器故障和外部干扰的存在，为了使得垂直起降飞机能够跟踪指定的参考指令且保持姿态稳定，通过建立级联的观测器系统对故障和干扰信息的实时估计，进而利用估计的干扰、故障和状态信息设计容错控制器，使得受影响的飞机尽可能地完成轨迹跟踪任务并且保持姿态平稳。在具体实施中，所提容错控制方法的仿真和检验都借助于MATLAB中的Simulink工具箱来实现。

本发明提出了如图1(b)所示的基于级联观测器的垂直起降飞机鲁棒容错控制系统，包括：执行机构、观测器单元、容错控制器。执行机构，在容错控制信号u_f以及故障f_a的作用下得到垂直起降飞机动的操作指令。垂直起降飞机在外部干扰f_d以及执行机构操作指令的作用下动作，以飞机质心位置(x,y)以及滚转角θ构成所述控制系统的输出向量Y=(x,y,θ)^T。观测器单元包括级联的干扰观测器、故障诊断观测器，干扰观测器根据飞机滚转角θ得到干扰估计值以及相应状态变量观测值，故障诊断观测器根据干扰估计值

飞机质心位置(x,y)、容错控制信号u_f得到故障估计值以及相应状态变量观测值，x、y分别为飞机质心在惯性坐标系中的水平坐标和垂直坐标，状态变量观测值包括：滚转角及角速率观测值

和

飞机质心位置

及相应的速率观测值和

容错控制器，以跟踪误差信号、故障估计值

干扰估计值为输入，得到容错控制信号修正值，跟踪误差信号为由飞机质心位置以及速率的观测值组成的观测信号与给定参考值 $(x_d\left(t\right),y_d\left(t\right),{\stackrel{\&CenterDot;}{x}}_d\left(t\right),{\stackrel{\&CenterDot;}{y}}_d\left(t\right))$ 的差值。

基于级联观测器的垂直起降飞机鲁棒容错控制方法，构建如图1(b)所示的基于级联观测器的垂直起降飞机鲁棒容错控制系统，具体包括如下步骤。

步骤1，建立外部干扰和滚转角的耦合程度为d₁(θ)=cosθ，d₂(θ)=(sinθ+θcosθ)，d₃(θ)=(cosθ+2)，第一、第二故障变量为f_1a=sinθf_a，f_2a=-cosθf_a时的垂直起降机动力学模型，表达式（1）变换为：

$\left\{\begin{array}{c}\stackrel{\&CenterDot;\&CenterDot;}{x}=-(u_{1f}-f_a\left(t\right))\sin \mathrm{\&theta;}+{\mathrm{\&epsiv;u}}_{2f}\cos \mathrm{\&theta;}+{\cos \mathrm{\&theta;f}}_d\left(t\right)\\ \stackrel{\&CenterDot;\&CenterDot;}{y}=(u_{1f}-f_a\left(t\right))\cos \mathrm{\&theta;}+{\mathrm{\&epsiv;u}}_{2f}\sin \mathrm{\&theta;}+{(\sin \mathrm{\&theta;}+\mathrm{\&theta;}\cos \mathrm{\&theta;})f}_d\left(t\right)-1\\ \stackrel{\&CenterDot;\&CenterDot;}{\mathrm{\&theta;}}={\mathrm{\&lambda;u}}_{2f}+(\cos \mathrm{\&theta;}+2)f_d\left(t\right)\end{array}\right..$

步骤2、选取状态变量向量x，x=(x₁,x₂,x₃,x₄,x₅,x₆)， $x_1=x,x_2=\stackrel{\&CenterDot;}{x},x_3=y,x_4=\stackrel{\&CenterDot;}{y},x_5=\mathrm{\&theta;},x_6=\stackrel{\&CenterDot;}{\mathrm{\&theta;}},$ 利用变量替换将步骤1建立的动力学模型转换为标准形式的非线性模型，表达式（2）变换为：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_1=x_2\\ {\stackrel{\&CenterDot;}{x}}_2=-(u_{1f}-f_a\left(t\right))\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+\cos x_5{f}_d\left(t\right)\\ {\stackrel{\&CenterDot;}{x}}_3=x_4\\ {\stackrel{\&CenterDot;}{x}}_4=(u_{1f}-f_a\left(t\right))\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+({\sin x}_5+x_5\cos x_5)f_d\left(t\right)-1\\ {\stackrel{\&CenterDot;}{x}}_5=x_6\\ {\stackrel{\&CenterDot;}{x}}_6={\mathrm{\&lambda;u}}_{2f}+(\cos x_5+2)f_d\left(t\right)\end{array}\right.$

输出向量为Y=[y₁,y₂,y₃]^T=[x₁,x₃,x₅]^T，用于轨迹跟踪的外部输出为y₁和y₂，故系统的相对阶为r=4，

选取第一至第三外部干扰变量f_1d=cosx₅f_d，f_2d=(sinx₅+x₅cosx₅)f_d，f_3d=(cosx₅+2)f_d，以及第一、第二故障变量f_1a=sinx₅f_a，f_2a=-cosx₅f_a，确定表达式（3）所示的故障系统的状态方程。

步骤3、由表达式（4）所示的干扰观测器以及表达式（8）所示的故障诊断观测器对干扰、故障在线估计，得到如图2（a）、图2（b）所示的干扰估计值

以及相应的状态估计值，如图3（a）。图3（b）、图3（c）所示的故障估计值

以及相应的状态估计值。

步骤4、利用坐标变换将非线性模型变换为表达式（9）所示的跟踪系统，再由表达式（10）所示的控制输入的可逆变换将跟踪系统解耦为表达式（11）所示的最小相位子系统以及表达式（13）所示的非最小相位子系统，进而得到表达式（14）所示的镇定最小相位子系统的最优控制律以及表达式（17）所示的镇定非最小相位子系统的滑模控制律。

步骤5，根据步骤4确定的最优控制律以及滑模控制律，结合故障估计值和可逆的控制变换，确定如表达式（18）所示的容错控制器的控制律。

下面分析非最小相位位子系统在上述滑模控制律作用下的稳定性

取李雅普诺夫函数对其求导，则有：

$\begin{array}{c}\stackrel{\&CenterDot;}{V}=s\stackrel{\&CenterDot;}{s}=s({\stackrel{\&CenterDot;}{z}}_2-C{\stackrel{\&CenterDot;}{z}}_1)\\ =s\left\{C\right[A_2{\hat{z}}_1+B_2{z}_2+\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\widehat{\mathrm{\&eta;}},\hat{e},{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)]-{\hat{f}}_{1d}+f_{1d}-\mathrm{\&sigma;sat}\left(\hat{s}\right)-C[A_2{z}_1+B_2{z}_2+\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)\}\\ \&le;\left|\hat{s}\right|\left(\right|\left|{\mathrm{CA}}_2\right|\left|\right||{\hat{z}}_1-z_1||+|{\mathrm{CB}}_2\left|\right|{\hat{z}}_2-z_2|+|f_{1d}-{\hat{f}}_{1d}|+|\left|C\right|\left|\right||\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)-\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)||-\mathrm{\&sigma;})\\ +|s-\hat{s}|\left(\right|\left|{\mathrm{CA}}_2\right|\left|\right||{\hat{z}}_1-z_1||+|{\mathrm{CB}}_2\left|\right|{\hat{z}}_2-z_2|+|f_{1d}-{\hat{f}}_{1d}|-\mathrm{\&sigma;}+|\left|C\right|\left|\right||\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)-\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)|\left|\right)\end{array}$

由于和

都渐近趋向零，选择切换项的增益满足如下条件：

$\mathrm{\&sigma;}>\left|\right|{\mathrm{CA}}_2\left|\right|\left|\right|{\hat{z}}_1-z_1\left|\right|+\left|{\mathrm{CB}}_2\right||{\hat{z}}_2-z_2|+|f_{1d}-{\hat{f}}_{1d}|+\left|\right|C\left|\right|\left|\right|\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)-\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)\left|\right|$

则有，

$\stackrel{\&CenterDot;}{V}\&le;\left|\hat{s}\right|\left(\right|\left|{\mathrm{CA}}_2\right|\left|\right||{\hat{z}}_1-z_1||+|{\mathrm{CB}}_2\left|\right|{\hat{z}}_2-z_2|+|f_{1d}-{\hat{f}}_{1d}|-\mathrm{\&sigma;}+|\left|C\right|\left|\right||\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\widehat{\mathrm{\&eta;}},\hat{e},{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)-\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)|\left|\right)\&le;0$

因此，当t→∞时，上述线性化的非最小相位子系统各个状态都能够从初始状态滑动到达平面s=0，即有z₂=Cz₁。在滑模控制律的作用下，有闭环系统：

${\stackrel{\&CenterDot;}{z}}_1=(A_2+B_{2}C)z_1+\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)$

由于

可以看成是系统的高阶项，在控制律的作用下能够迅速地收敛到零，因此闭环系统指数渐近稳定。根据之前的可逆变换可得，对于任意二阶可微的期望输出x_d(t)，当t→∞时都有y₁(t)→x_d(t)和同时也有η₁(t)→0和η₂(t)→0，即飞机的滚转姿态也渐近稳定。执行器卡死故障的重构信号可由步骤3中的估计信息变换获得，为：

其中：M和μ均为足够小的正数，用于保证的连续性。

本发明在MATLAB7.0环境下对所设计的鲁棒容错控制方法进行仿真验证实验，具体如下：

(1)VTOL飞机仿真模型的初始状态选择为(0.1,0.1,0.1,0.1,0.1,0)^T，系统参数选择为ε=0.5，λ=1，选取参考轨迹为：

x_d(t)=5        0s≤t≤500s

$y_d\left(t\right)=\left\{\begin{array}{cc}0& t\&le;50s\\ 0.01(t-50)& 50s<t\&le;500s\end{array}\right.$

(2)外部干扰和推力发动机的执行器卡死故障选取为：

f_d(t)=0.05sin(0.1t)+0.04cos(0.05t)

$f_a\left(t\right)=\left\{\begin{array}{cc}0& 0s\&le;t<150s\\ 0.4& 150s\&le;t<250s\\ 0& 250s\&le;t<350s\\ 0.6& 350s\&le;t\&le;500s\end{array}\right.$

(3)干扰观测器的切换项增益为λ₁=3，λ₂=20，故障诊断观测器增益为l₁=3，l₂=5，故障更新参数选择为T₁=T₂=5。选取权值矩阵Q₁=diag(100,100)，R₁=50，经计算可得K₁=(1.4142 2.1974)。滑模面的增益参数C=(4 -10 -5)，饱和函数sat(·)的消除抖振参数选为γ=0.05，滑模控制律的切换项增益选取为σ=5。故障估计信号变换中的参数M=0.01，μ=0.05。

结果说明：

如图2(a)和图2(b)所示，可以看出，干扰观测器的干扰和状态估计曲线能够在很短的时间内迅速精确地收敛到其实际值；

如图3(a)、图3(b)和图3(c)所示，可看出，故障诊断观测器的故障和状态估计曲线能够迅速准确地收敛到真实值；

如图4所示的容错控制器曲线，控制器利用观测器系统估计到的信息，实时补偿故障和干扰，最大限度地减少飞机飞行性能的下降；

如图5(a)、图5(b)所示，可以看出，系统受到故障和干扰的影响出现了一定的性能下降，但在容错控制器的作用下，迅速地减少了影响，使得飞机飞行性能表现满足要求。

由以上可知：本发明针对垂直起降飞机同时存在外部干扰和推力通道执行器故障的情况，提出一种实时有效的鲁棒容错控制方案，能够较好地在线实时精确重构外部扰动和故障；并且所设计基于估计信息的容错控制器通过补偿能够达到很好地容错效果，完成飞行要求。

Claims (3)

1.基于级联观测器的垂直起降飞机鲁棒容错控制系统，是串接在垂直起降飞机模型输入端、输出端之间闭环控制系统，具体包括：执行机构、观测器单元、容错控制器，其特征在于：

所述执行机构，在容错控制信号以及故障的作用下得到垂直起降飞机的操作指令；

所述垂直起降飞机模型在外部干扰以及执行机构操作指令的作用下动作，以飞机质心位置以及滚转角构成所述控制系统的输出向量；

所述观测器单元包括级联的干扰观测器、故障诊断观测器，干扰观测器根据飞机滚转角和滚转速率来得到干扰估计值以及状态变量观测值，故障诊断观测器根据干扰估计值、飞机质心位置、容错控制信号得到故障估计值以及状态变量观测值，所述状态变量观测值包括：滚转角和滚转速率的观测值、飞机质心位置观测值以及速度观测值；

所述容错控制器，以跟踪误差信号、滚转角和滚转速率的观测值、故障估计值、干扰估计值为输入，得到容错控制信号修正值，所述跟踪误差信号为由飞机质心位置以及相应速率的观测值组成的观测信号与给定参考值的差值。

2.基于级联观测器的垂直起降飞机鲁棒容错控制方法，其特征在于构建如权利要求1所述的控制系统，具体包括如下步骤：

步骤1，建立考虑了外部干扰以及执行器故障的垂直起降机动力学模型：

$\left\{\begin{array}{c}\stackrel{\&CenterDot;\&CenterDot;}{x}=-(u_{1f}-f_a\left(t\right))\sin \mathrm{\&theta;}+{\mathrm{\&epsiv;u}}_{2f}\cos \mathrm{\&theta;}+d_1\left(\mathrm{\&theta;}\right)f_d\left(t\right)\\ \stackrel{\&CenterDot;\&CenterDot;}{y}=(u_{1f}-f_a\left(t\right))\cos \mathrm{\&theta;}+{\mathrm{\&epsiv;u}}_{2f}\sin \mathrm{\&theta;}+d_2\left(\mathrm{\&theta;}\right)f_d\left(t\right)-g\\ \stackrel{\&CenterDot;\&CenterDot;}{\mathrm{\&theta;}}={\mathrm{\&lambda;u}}_{2f}+d_3\left(\mathrm{\&theta;}\right)f_d\left(t\right)\end{array}\right.,$

其中，

分别为飞机在水平方向、垂直方向的加速度，

为滚转角加速度，u_1f、u_2f分别为控制系统的推力容错控制输入、滚转力矩容错控制输入，ε为滚转力矩与水平方向加速度的耦合系数，f_d(t)、f_a(t)分别为外部干扰、故障的时间函数，d_i(θ),i=1,2,3为外部干扰和滚转角的耦合程度，g为重力加速度，λ为滚转力矩的传输效率；

步骤2，选取状态变量向量x，x=(x₁,x₂,x₃,x₄,x₅,x₆)， $x_1=x,x_2=\stackrel{\&CenterDot;}{x},x_3=y,x_4=\stackrel{\&CenterDot;}{y},x_5=\mathrm{\&theta;},x_6=\stackrel{\&CenterDot;}{\mathrm{\&theta;}},$ 利用变量替换将步骤1建立的动力学模型转换为标准形式的非线性模型：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_1=x_2\\ {\stackrel{\&CenterDot;}{x}}_2={-u}_{1f}\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+d_1\left(x_5\right)f_d+\sin x_5{f}_a\\ {\stackrel{\&CenterDot;}{x}}_3=x_4\\ {\stackrel{\&CenterDot;}{x}}_4=u_{1f}\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+d_2\left(x_5\right)f_d-\cos x_5{f}_a-1\\ {\stackrel{\&CenterDot;}{x}}_5=x_6\\ {\stackrel{\&CenterDot;}{x}}_6={\mathrm{\&lambda;u}}_{2f}+d_3\left(x_5\right)f_d\end{array}\right.,$

选取第一至第三外部干扰变量f_1d=d₁(x₅)f_d,f_2d=d₂(x₅)f_d,f_3d=d₃(x₅)f_d以及第一、第二故障变量f_1a=sinx₅f_a,f_2a=-cosx₅f_a，确定故障系统的状态方程为：

$\left\{\begin{array}{c}{\stackrel{\&CenterDot;}{x}}_1=x_2\\ {\stackrel{\&CenterDot;}{x}}_2={-u}_{1f}\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+f_{1d}+f_{1a}\\ {\stackrel{\&CenterDot;}{x}}_3=x_4\\ {\stackrel{\&CenterDot;}{x}}_4=u_{1f}\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+f_{2d}-f_{2a}-1\\ {\stackrel{\&CenterDot;}{x}}_5=x_6\\ {\stackrel{\&CenterDot;}{x}}_6={\mathrm{\&lambda;u}}_{2f}+f_{3d}\end{array}\right.,$

其中：x、y分别为飞机质心在惯性坐标系中的水平坐标和垂直坐标，分别为飞机在水平方向、垂直方向的运动速率，

是滚转速率；

步骤3，利用干扰观测器 $\left\{\begin{array}{c}{\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_5={\hat{x}}_6+{\mathrm{\&lambda;}}_1\mathrm{sgn}(x_5-{\hat{x}}_5)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_6={\mathrm{\&lambda;u}}_{2f}+{\mathrm{\&lambda;}}_2{\mathrm{sgn}}_1({\tilde{x}}_6-{\hat{x}}_6)\end{array}\right.$ 得到干扰估计值

其中：

分别是滚转角、滚转速率的观测值，λ₁和λ₂是切换项的增益，

是滚转速率x₆的替代变量；

利用故障诊断观测器 $\left\{\begin{array}{c}{\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_1={\hat{x}}_2-l_1({\hat{x}}_1-x_1)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_2={-u}_{1f}\sin x_5+{\mathrm{\&epsiv;u}}_{2f}\cos x_5+{\hat{f}}_{1a}+{\hat{f}}_{1d}-({\hat{x}}_1-x_1)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_3={\hat{x}}_4-l_2({\hat{x}}_3-x_3)\\ {\stackrel{\stackrel{\&CenterDot;}{^}}{x}}_4=u_{1f}\cos x_5+{\mathrm{\&epsiv;u}}_{2f}\sin x_5+{\hat{f}}_{2a}+{\hat{f}}_{2d}-({\hat{x}}_3-x_3)-1\end{array}\right.$ 得到故障估计值

其中：

为飞机水平坐标、垂直坐标、水平速率、垂直速率的观测值，

为第一、第二故障变量的估计值，

为第一、第二外部干扰变量的估计值，l₁和l₂为故障诊断观测器的增益；

步骤4，选取控制输入可逆变换矩阵将系统输入输出线性化 $\left(\begin{array}{c}{w}_{s1}\\ w_{s2}\end{array}\right)=\mathrm{\&beta;}\left(x_5\right)\left(\begin{array}{c}{u}_{1f}-f_a\\ u_{2f}\end{array}\right)-\left(\begin{array}{c}{\stackrel{\&CenterDot;\&CenterDot;}{x}}_d\\ {\stackrel{\&CenterDot;\&CenterDot;}{y}}_d+1\end{array}\right),$ 其中 $\mathrm{\&beta;}\left(x_5\right)=\left(\begin{array}{cc}-\sin x_5& \mathrm{\&epsiv;}\cos x_5\\ \cos x_5& \mathrm{\&epsiv;}\sin x_5\end{array}\right)$ 是系统在平衡点x₅=0附近的非奇异解耦矩阵，再取六个误差变量e₁,e₂,e₃,e₄,η₁,η₂，e₁=x₁-x_d，

e₃=x₃-y_d， $e_4=x_4-{\stackrel{\&CenterDot;}{y}}_d,$ η₁=x₅， ${\mathrm{\&eta;}}_2=\frac{\mathrm{\&epsiv;}}{\mathrm{\&lambda;}}{x}_6-e_2\cos x_5-e_4\sin x_5,$ 将控制系统坐标变换为跟踪系统，将系统解耦为最小相位子系统：

及非最小相位子系统:

确定镇定最小相位子系统的最优控制律ws2为：

镇定非最小相位子系统的滑模控制律w_s1为： $w_{s1}=C[A_2{\hat{z}}_1+B_2{\hat{z}}_2+\stackrel{\widehat{\&OverBar;}}{\mathrm{\&psi;}}(\widehat{\mathrm{\&eta;}},\hat{e},{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)]-{\hat{f}}_{1d}-\mathrm{\&sigma;sat}\left(\hat{s}\right),$

其中：z₁=(e₁,η^T)^T，z₂=e₂，η=(η₁,η₂)^T，e=(e₁,e₂,e₃,e₄)^T，

$\stackrel{\&OverBar;}{\mathrm{\&psi;}}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)=\mathrm{\&psi;}(\mathrm{\&eta;},e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)-A_2{z}_1-B_2{z}_2,$ $\mathrm{\&psi;}({\mathrm{\&eta;}}_e,e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)=\left(\begin{array}{c}{e}_2\\ \mathrm{\&Gamma;}({\mathrm{\&eta;}}_e,e,{\stackrel{\&CenterDot;\&CenterDot;}{Y}}_d)\end{array}\right),$ $({\mathrm{\&beta;}}_1,{\mathrm{\&beta;}}_2)=K_1={R_1}^{-1}{B_1}^T{P}_1$ 是最优反馈增益，R₁和P₁是满足相应Riccati方程的正定矩阵，C为使得A₂+B₂C为Hurwitz矩阵的向量，

为误差变量e₁,e₂,e₃,e₄的估计值，

分别为水平方向、垂直方向参考加速度，σ为滑动切换项的增益，

是以滑模平面估计值

为变量的饱和函数；

步骤5，根据步骤4确定的最优控制律以及滑模控制律，结合故障估计值和可逆的控制变换，确定容错控制器的控制律：

$u_f={\mathrm{\&beta;}}^{-1}\left(x_5\right)\left(\begin{array}{c}{w}_{s1}+{\stackrel{\&CenterDot;\&CenterDot;}{x}}_d\\ w_{s2}+{\stackrel{\&CenterDot;\&CenterDot;}{y}}_d+1\end{array}\right)+\hat{F},$

其中：β^-1(x₅)是可逆解耦矩阵β(x₅)的逆矩阵，

为重构故障向量的估计值。

3.根据权利要求2所述的基于级联观测器的垂直起降飞机鲁棒容错控制方法，其特征在于：外部干扰和滚转角的耦合程度为d₁(θ)=cosθ,d₂(θ)=(sinθ+θ₅cosθ),d₃(θ)=(cosθ+2)，第一、第二故障变量为f_1a=sinθf_a,f_2a=-cosθf_a时，步骤1建立的动力学模型为：

$\left\{\begin{array}{c}\stackrel{\&CenterDot;\&CenterDot;}{x}=-(u_{1f}-f_a\left(t\right))\sin \mathrm{\&theta;}+{\mathrm{\&epsiv;u}}_{2f}\cos \mathrm{\&theta;}+{\cos \mathrm{\&theta;f}}_d\left(t\right)\\ \stackrel{\&CenterDot;\&CenterDot;}{y}=(u_{1f}-f_a\left(t\right))\cos \mathrm{\&theta;}+{\mathrm{\&epsiv;u}}_{2f}\sin \mathrm{\&theta;}+{(\sin \mathrm{\&theta;}+\mathrm{\&theta;}\cos \mathrm{\&theta;})f}_d\left(t\right)-1\\ \stackrel{\&CenterDot;\&CenterDot;}{\mathrm{\&theta;}}={\mathrm{\&lambda;u}}_{2f}+(\cos \mathrm{\&theta;}+2)f_d\left(t\right)\end{array}\right..$

Images