一种飞艇三维航迹跟踪的反步神经网络控制方法

摘要

本发明涉及一种飞艇三维航迹跟踪的反步神经网络控制方法。针对飞艇的航迹跟踪控制问题，本发明建立了飞艇的非线性动力学模型；以此为受控对象，将非线性动力学模型分解为两个子系统，采用反步法为每个子系统设计Lyapunov函数和中间虚拟控制量，通过确定适当的虚拟反馈，使得系统的前面状态达到渐近稳定，一直“反向推演”至整个系统，从而实现整个系统的渐近稳定；针对飞艇动力学模型不确定问题，采用神经网络逼近器精确估计未知的飞艇动力学模型，以提高控制精度和系统性能。由该方法控制的闭环系统能够高精度跟踪任意给定的参数化指令航迹，且具有良好的稳定性、适应性、鲁棒性和动态性能，为飞艇航迹控制的工程实现提供了有效方案。

说明书

技术领域

本发明涉及一种航天航空领域的飞行控制方法，它为飞艇航迹跟踪提供一种反 步神经网络控制方法，属于自动控制技术领域。

背景技术

飞艇是指一种依靠轻于空气的气体(如氦气、氢气等)提供静浮力升空，依靠 自动飞行控制系统实现定点驻留和低速机动的飞行器，具有滞空时间长、能耗低、 效费比高及定点驻留等优点，广泛应用于侦察监视、对地观测、环境监测、应急救 灾、科学探测等领域，具有重要应用价值和广阔的应用前景，当前已成为航空领域 的研究热点。

航迹跟踪是指飞艇从给定的初始状态出发并跟踪给定的惯性坐标系下的指令航 迹。飞艇的空间运动具有非线性、通道耦合、不确定、易受外界扰动等特点，因此， 航迹控制成为飞艇飞行控制的关键技术之一。众多研究人员针对飞艇的航迹跟踪问 题，提出了PID控制、反馈控制、滑模控制、鲁棒控制等方法，为飞艇航迹跟踪提 供了可供参考借鉴的技术方案。但是上述航迹控制方法尚未有效解决以下两类问题： 一是飞艇动力学模型不确定，存在建模误差及未建模动态；二是飞艇航迹控制系统 为一个复杂的多变量非线性系统，飞行包线内闭环控制系统的稳定性难以保证。

发明内容

针对现有技术存在的不足，本发明的目的是提供一种飞艇三维航迹跟踪的反步神 经网络控制方法。

本发明针对飞艇三维航迹跟踪问题，建立了飞艇的非线性动力学模型；以此为 受控对象，将非线性动力学模型分解为两个子系统，采用反步法为每个子系统设计 李雅普诺夫(Lyapunov)函数和中间虚拟控制量，通过确定适当的虚拟反馈，使得 系统的前面状态达到渐近稳定，一直“反向推演”至整个系统，从而实现整个系统 的渐近稳定；针对飞艇动力学模型不确定问题，采用神经网络精确逼近未知的飞艇 动力学模型，以提高控制精度和系统性能。本发明的优点表现在：①采用反步法设 计使得李雅普诺夫(Lyapunov)函数和控制律的设计过程系统化、结构化，确保了 系统的稳定性；②采用神经网络精确逼近飞艇的不确定模型，使得航迹跟踪控制系 统具有强适应性和强鲁棒性。

本发明的技术方案是：首先由给定的指令航迹和实际航迹计算航迹控制误差量， 然后采用反步方法设计航迹控制律，计算航迹控制量；为解决飞艇动力学模型不确 定问题，采用神经网络精确逼近未知的不确定模型。实际应用中，飞艇航迹由组合 导航系统测量得到，将由该方法计算得到的控制量传输至执行机构即可实现航迹控 制功能。

具体地，一种飞艇三维航迹跟踪的反步神经网络控制方法，包括以下步骤：

步骤一：给定指令航迹：η_d＝[x_d,y_d,z_d,θ_d,ψ_d,φ_d]^T；

其中：所述的指令航迹为广义坐标η_d＝[x_d,y_d,z_d,θ_d,ψ_d,φ_d]^T，x_d、y_d、z_d、θ_d、 ψ_d和φ_d分别为指令x坐标、指令y坐标、指令z坐标、指令俯仰角、指令偏航角和 指令滚转角，上标T表示向量或矩阵的转置。

步骤二：航迹控制误差量计算：计算指令航迹与实际航迹之间的航迹控制误差 量e；

所述航迹控制误差量e的计算方法为：

e＝η_d-η＝[x_d-x,y_d-y,z_d-z,θ_d-θ,ψ_d-ψ,φ_d-φ]^T       (1)

其中：η＝[x,y,z,θ,ψ,φ]^T为实际航迹，x、y、z、θ、ψ、φ分别为实际航迹 的x坐标、y坐标、z坐标、俯仰角、偏航角和滚转角。

步骤三：航迹控制律设计：选取Lypaunov函数和中间虚拟控制量，采用反步法 设计航迹控制律，计算航迹控制量u，具体包括以下步骤：

1)建立飞艇的动力学模型

飞艇空间运动的坐标系及运动参数定义如下：采用地面坐标系o_ex_ey_ez_e和体坐标 系o_bx_by_bz_b对飞艇的空间运动进行描述，CV为浮心，CG为重心，浮心到重心的矢量 为r_G＝[x_G,y_G,z_G]^T。运动参数定义：位置P＝[x,y,z]^T，x、y、z分别为轴向、侧向和 竖直方向的位移；姿态角Ω＝[θ,ψ,φ]^T，θ、ψ、φ分别为俯仰角、偏航角和滚转角； 速度v＝[u,v,w]^T，u、v、w分别为体坐标系中轴向、侧向和垂直方向的速度；角速度 ω＝[p,q,r]^T，p、q、r分别为滚转、俯仰和偏航角速度。记广义坐标 η＝[x,y,z,θ,ψ,φ]^T，广义速度为V＝[u,v,w,p,q,r]^T。

飞艇的动力学模型描述如下：

$\stackrel{·}{\eta }=J\left(\eta \right)=\left(\begin{array}{cc}{J}_1& 0_{3\times 3}\\ 0_{3\times 3}& J_2\end{array}\right)V---\left(2\right)$

$M\stackrel{·}{V}=\bar{N}+\bar{G}+\tau ---\left(3\right)$

式中

$J_1=\left(\begin{array}{ccc}\cos \psi \cos \theta & \cos \psi \sin \theta \sin \phi -\sin \psi \cos \phi & \cos \psi \sin \theta \cos \phi +\sin \psi \sin \phi \\ \sin \psi \cos \theta & \sin \psi \sin \theta \sin \phi +\cos \psi \cos \phi & \sin \psi \sin \theta \cos \phi -\cos \psi \sin \phi \\ -\sin \theta & \cos \theta \sin \phi & \cos \theta \cos \phi \end{array}\right)---\left(4\right)$

$J_2=\left(\begin{array}{ccc}0& \cos \phi & -\sin \phi \\ 0& \sec \theta \sin \phi & \sec \theta \cos \phi \\ 1& \tan \theta \sin \phi & \tan \theta \cos \phi \end{array}\right)---\left(5\right)$

$M=\left(\begin{array}{cccccc}m+m_{11}& 0& 0& 0& {\mathrm{mz}}_G& -{\mathrm{my}}_G\\ 0& m+m_{22}& 0& -{\mathrm{mz}}_G& 0& {\mathrm{mx}}_G\\ 0& 0& m+m_{33}& {\mathrm{my}}_G& -{\mathrm{mx}}_G& 0\\ 0& -{\mathrm{mz}}_G& {\mathrm{my}}_G& I_x+I_{11}& 0& 0\\ {\mathrm{mz}}_G& 0& -{\mathrm{mx}}_G& 0& I_x+I_{22}& 0\\ 0& {\mathrm{mx}}_G& 0& -I_{\mathrm{xz}}& 0& I_x+I_{33}\end{array}\right)---\left(6\right)$

$\bar{G}=\left(\begin{array}{c}(B-G)\sin \theta \\ (G-B)\cos \theta \sin \phi \\ (G-B)\cos \theta \cos \phi \\ y_{G}G\cos \theta \cos \phi -z_{G}G\cos \theta \sin \phi \\ -x_{G}G\cos \theta \cos \phi -z_{G}G\sin \theta \\ x_{G}G\cos \theta \sin \phi +y_{G}G\sin \theta \end{array}\right)---\left(7\right)$

$\tau =\left(\begin{array}{c}{\tau }_u\\ {\tau }_v\\ {\tau }_w\\ {\tau }_l\\ {\tau }_m\\ {\tau }_n\end{array}\right)=\left(\begin{array}{c}T\cos \mu \cos \upsilon \\ T\sin \mu \\ T\cos \mu \sin \upsilon \\ T\cos \mu \sin \upsilon l_y-T\sin {\mathrm{\mu l}}_z\\ T\cos \mu \cos \upsilon l_z-T\cos \mu \sin \upsilon l_x\\ T\sin {\mathrm{\mu l}}_x-T\cos \mu \cos \upsilon l_y\end{array}\right)---\left(8\right)$

$\bar{N}={[N_u,N_v,N_w,N_p,N_q,N_r]}^T---\left(9\right)$

其中

N_u＝(m+m₂₂)vr-(m+m₃₃)wq+m[x_G(p²+r²)-y_Gpq-z_Gpr]

                                         (10)

+QV^2/3(-C_X cosαcosβ+C_Y cosαsinβ+C_Z sinα)

N_v＝(m+m₃₃)wp-(m+m₁₁)ur-m[x_Gpq-y_G(p²+r²)+z_Gqr]

                                              (11)

+QV^2/3(C_X sinβ+C_Y cosβ)

N_w＝(m+m₂₂)vp-(m+m₁₁)uq-m[x_Gpr+y_Gqr-z_G(p²+q²)]

                                         (12)

+QV^2/3(-C_X sinαsinβ+C_Y sinαcosβ-C_Z cosα)

N_p＝[(I_y+m₅₅)-(I_z+I₆₆)]qr+I_xzpq-I_xypr-I_yz(r²-q²)+

                                          (13)

[mz_G(ur-wp)+y_G(uq-vp)]+QVC_l

N_q＝[(I_z+m₆₆)-(I_x+I₄₄)]pr+I_xyqr-I_yzpq-I_xz(p²-r²)

                                             (14)

+m[x_G(vp-uq)-z_G(wp-vr)]+QVC_m

N_r＝[(I_y+m₅₅)-(I_x+I₄₄)]pq-I_xzqr-I_xy(q²-p²)+I_yzpr

                                             (15)

+m[y_G(wq-vr)-x_G(ur-wp)]+QVC_n

式中，m为飞艇质量，m₁₁、m₂₂、m₃₃为附加质量,I₁₁、I₂₂、I₃₃为附加惯量， V为飞艇体积；；Q为动压，α为迎角,β为侧滑角，C_X、C_Y、C_Z、C_l、C_m、C_n为气动系数；I_x、I_y、I_z分别为绕o_bx_b、o_by_b、o_bz_b的主惯量；I_xy、I_xz、I_yz分别 为关于平面o_bx_by_b、o_bx_bz_b、o_by_bz_b的惯量积；T为推力大小，μ为推力矢量与o_bx_bz_b面之间的夹角，规定其在o_bx_bz_b面之左为正，υ为推力矢量在o_bx_bz_b面的投影与o_bx_b轴 之间的夹角，规定其投影在o_bx_b轴之下为正；l_x、l_y、l_z表示推力作用点距原点o_b的 距离。

式(3)为关于广义速度V的表达式，需要将其变换为关于广义坐标η的表达式。

由式(1)可得：

$V=J^{-1}\left(\eta \right)\stackrel{·}{\eta }=R\left(\eta \right)\stackrel{·}{\eta }=\left(\begin{array}{cc}A& 0_{3\times 3}\\ 0_{3\times 3}& B\end{array}\right)\stackrel{·}{\eta }---\left(16\right)$

式中，J^-1(η)为J(η)的逆矩阵。

$A=\left(\begin{array}{ccc}\cos \psi \cos \theta & \sin \psi \cos \theta & -\sin \theta \\ \cos \psi \sin \theta \sin \phi -\sin \psi \cos \phi & \sin \psi \sin \theta \sin \phi +\cos \psi \cos \phi & \cos \theta \sin \phi \\ \cos \psi \sin \theta \cos \phi +\sin \psi \sin \phi & \sin \psi \sin \theta \cos \phi -\cos \psi \sin \phi & \cos \theta \cos \phi \end{array}\right)---\left(17\right)$

$B=\left(\begin{array}{ccc}0& -\sin \theta & 1\\ \cos \phi & \cos \theta \sin \phi & 0\\ -\sin \phi & \cos \theta \cos \phi & 0\end{array}\right)---\left(18\right)$

对式(16)微分，可得

$\stackrel{·}{V}=\stackrel{·}{R}\stackrel{·}{\eta }+R\stackrel{··}{\eta }---\left(19\right)$

式中

$\stackrel{·}{R}=\left(\begin{array}{cc}\stackrel{·}{A}& 0_{3\times 3}\\ 0_{3\times 3}& \stackrel{·}{B}\end{array}\right)---\left(20\right)$

式(19)左乘可得

$R^{T}M\stackrel{·}{V}=R^{T}M\stackrel{·}{R}\stackrel{·}{\eta }+R^T\mathrm{MR}\stackrel{··}{\eta }---\left(21\right)$

综合式(3)、式(19)以及式(21)可得:

$M_{\eta }\left(\eta \right)\stackrel{··}{\eta }+N_{\eta }(\eta ,\stackrel{·}{\eta })\stackrel{·}{\eta }+G_{\eta }\left(\eta \right)=u---\left(22\right)$

式中

M_η(η)＝R^TMR       (23)

$N_{\eta }(\eta ,\stackrel{·}{\eta })=R^{T}M\stackrel{·}{R}---\left(24\right)$

$G_{\eta }\left(\eta \right)=-R^T(\bar{N}+\bar{G})---\left(25\right)$

u＝R^Tτ      (26)

令x₁＝η，则动力学模型式(22)可写为如下形式：

$\left(\begin{array}{c}{\stackrel{·}{x}}_1=x_2\\ {\stackrel{·}{x}}_2=M_{\eta }^{-1}u-M_{\eta }^{-1}{N}_{\eta }\stackrel{·}{\eta }-G_{\eta }\left(\eta \right)\end{array}\right)---\left(27\right)$

式中，表示矩阵M_η的逆矩阵；

以式(27)所描述的数学模型为被控对象，采用反步方法设计航迹控制律。

2)设计航迹控制律

根据指令航迹与实际航迹之间的航迹控制误差量e，定义如下虚拟量：

                                           $a_1=-k_{1}e+\stackrel{·}{\eta }d$

                                                  (28)

式中，α₁为虚拟量，k₁为可调的控制参数。

定义虚拟量α₁与x₂之间的误差e′：

e′＝x₂-α₁         (29)

式(29)对时间微分并将式(27)代入，可得：

${\stackrel{·}{e}}^{\prime }={\stackrel{·}{x}}_2-{\stackrel{·}{a}}_1=M_{\eta }^{-1}u-M_{\eta }^{-1}{N}_{\eta }\left(\eta \right)\stackrel{·}{\eta }-G_{\eta }\left(\eta \right)+k_1\stackrel{·}{e}-{\stackrel{··}{\eta }}_d---\left(30\right)$

令

$f\left(x\right)=-M_{\eta }^{-1}{N}_{\eta }\left(\eta \right)\stackrel{·}{\eta }-G_{\eta }\left(\eta \right)---\left(31\right)$

则式(31)可表示为：

${\stackrel{·}{e}}^{\prime }=f\left(x\right)+M_{\eta }^{-1}u+k_1\stackrel{·}{e}-{\stackrel{··}{\eta }}_d---\left(32\right)$

选取李雅普诺夫函数(Lyapunov函数)V₁

$V_1=\frac{1}{2}{e}^{T}e+\frac{1}{2}{e}^{\prime T}{e}^{\prime }---\left(33\right)$

式(33)对时间微分，并将式(32)代入，可得：

${\stackrel{·}{V}}_1=-k_1{e}^{T}e+e^{\prime T}(f\left(x\right)+M_{\eta }^{-1}u+k_1\stackrel{·}{e}-{\stackrel{··}{\eta }}_d)---\left(34\right)$

根据式(34)，设计如下航迹控制律：

$u=M_{\eta }(-e-k_1\stackrel{·}{e}-k_2{e}^{\prime }+{\stackrel{··}{\eta }}_d-f\left(x\right))---\left(35\right)$

3)稳定性分析

将航迹控制律式(35)代入式(34)，可得：

$\left(\begin{array}{c}{\stackrel{·}{V}}_1=-k_1{e}^{T}e+e^T{e}^{\prime }+e^{\prime T}\{f\left(x\right)+M_{\eta }^{-1}[M_{\eta }(-e-k_1\stackrel{·}{e}-k_2{e}^{\prime }+{\stackrel{··}{\eta }}_d-f\left(x\right))]+k_1\stackrel{·}{e}-k{\stackrel{··}{\eta }}_d\}\\ =-k_1{e}^{T}e+e^{\prime T}{e}^{\prime }+e^{\prime T}(-e-k_2{e}^{\prime })=-k_1{e}^{T}e-k_2{e}^{\prime T}{e}^{\prime }\end{array}\right)---\left(36\right)$

式(36)表明：采用航迹控制律(35)能够保证闭环系统的稳定性。

步骤四：神经网络逼近器设计：以航迹控制误差e及其变化率实际航迹η及 其变化率为神经网络的输入变量，以飞艇动力学模型的估计值为神经网络的 输出变量设计了神经网络逼近器，利用神经网络无限逼近功能估计未知的不确定模 型，以提供控制精度，具体步骤如下：

1)由于在实际飞行过程中难以对飞艇进行精确建模，f(x)为未知函数，难以根 据式(35)进行控制律解算，因此，须采用f(x)的估计值对航迹控制律式(35)进 行解算；采用神经网络逼近未知函数f(x),则有：

f(x)＝w^Th(x)+ε          (37)

式中，w为神经网络的权重向量，ε为逼近误差，h(x)＝[h_i(x)]^T，h_i(x)为高斯 基函数，下标i表示第i个高斯基函数；

2)选择输入输出变量

令航迹控制误差e及其变化率实际航迹η及其变化率为神经网络逼近器的 输入变量，令估计值为神经网络逼近器的输出变量。

3)设计神经网络结构

神经网络结构包括输入层、隐层和输出层。

输入层：选取神经网络的输入变量为$x={\left(\begin{array}{cccc}e& \stackrel{·}{e}& \eta & \stackrel{·}{\eta }\end{array}\right)}^T.$

隐层：选取高斯函数作为隐层节点的基函数

$h_i\left(x\right)=\exp \left(\frac{{\left|\right|x-c\left|\right|}^2}{{2\sigma }_i^2}\right)---\left(38\right)$

其中，c＝[c_i]^T，c_i为第i个高斯函数的中值，σ_i为第i个节点的基宽度参数，||·|| 表示欧几里德范数。

输出层：神经网络的输出为

$\hat{f}\left(x\right)={\hat{w}}^{T}h\left(x\right)---\left(39\right)$

其中，为w的估计值。

4)稳定性分析

定义与w差值：

$\tilde{w}=w-\hat{w}---\left(40\right)$

选取Lyapunov函数：

$V_2=\frac{1}{2}{\xi }^T\xi +\frac{1}{2}\mathrm{tr}\left({\tilde{W}}^T{Q}^{-1}\tilde{W}\right)---\left(41\right)$

式中，$\xi ={\left(\begin{array}{cc}{e}_1& e_2\end{array}\right)}^T,\tilde{W}=\left(\begin{array}{cc}0& 0\\ 0& \tilde{w}\end{array}\right),Q=\left(\begin{array}{cc}0& 0\\ 0& \Gamma \end{array}\right),$Γ为可调的正定矩阵，Q^-1表 示矩阵Q的逆矩阵。

对式(41)微分，可得：

$\left(\begin{array}{c}{\stackrel{·}{V}}_3={\xi }^T\xi +\mathrm{tr}\left({\tilde{W}}^T{Q}^{-1}\stackrel{·}{\stackrel{~}{W}}\right)\\ =-k_1{e}_1^T{e}_1-k_2{e}_2^T{e}_2+e_2^T({\hat{w}}^{T}h\left(x\right)+ϵ)+\mathrm{tr}\left({\tilde{W}}^T{Q}^{-1}\stackrel{·}{\stackrel{~}{W}}\right)\end{array}\right)---\left(42\right)$

定义

$k=\left(\begin{array}{cc}{k}_1& 0\\ 0& k_2\end{array}\right)---\left(43\right)$

则式(42)可写为：

$\left(\begin{array}{c}{\stackrel{·}{V}}_3={\xi }^T\stackrel{·}{\xi }+\mathrm{tr}\left({\tilde{w}}^T{Q}^{-1}\stackrel{·}{\stackrel{~}{w}}\right)\\ =-{\xi }^T\mathrm{k\xi }+{\xi }^Tϵ+{\xi }^T\tilde{W}\Psi +\mathrm{tr}\left({\tilde{W}}^T{Q}^{-1}\stackrel{·}{\stackrel{~}{W}}\right)\\ =-{\xi }^T\mathrm{k\xi }+{\xi }^Tϵ+\mathrm{tr}(\tilde{W}\Psi {\xi }^T+{\tilde{W}}^T{Q}^{-1}\stackrel{·}{\stackrel{~}{W}})\end{array}\right)---\left(44\right)$

式中，$\Psi =\left(\begin{array}{c}0\\ h\left(x\right)\end{array}\right).$

设计如下自适应律：

$\stackrel{·}{\hat{W}}=\mathrm{Q\Psi }{\xi }^T-\mathrm{\gamma Q}\left|\right|\xi \left|\right|\hat{W}---\left(45\right)$

式中，γ＞0为可调的参数，$\hat{W}=\left(\begin{array}{cc}0& 0\\ 0& \hat{w}\end{array}\right).$

将自适应律代入式(44)，可得：

$\left(\begin{array}{c}{\stackrel{·}{V}}_3={\xi }^T\stackrel{·}{\xi }+\mathrm{tr}\left({\tilde{w}}^T{Q}^{-1}\stackrel{·}{\stackrel{~}{w}}\right)\\ =-{\xi }^T\mathrm{k\xi }+{\xi }^Tϵ+\mathrm{tr}(\tilde{W}\Psi {\xi }^T-{\tilde{W}}^T{Q}^{-1}(\mathrm{Q\Psi }{\xi }^T-\mathrm{\gamma Q}|\left|\xi \right|\left|\hat{W}\right))\\ =-{\xi }^T\mathrm{k\xi }+{\xi }^Tϵ+\mathrm{tr}\left(\gamma {\tilde{W}}^T\right|\left|\xi \right|\left|\hat{W}\right)\\ =-{\xi }^T\mathrm{k\xi }+{\xi }^Tϵ+\gamma \left|\right|\xi \left|\right|\mathrm{tr}\left({\tilde{W}}^T(W-\tilde{W})\right)\end{array}\right)---\left(46\right)$

根据施瓦茨(Schwarz)不等式，有：

$\mathrm{TR}\left({\tilde{W}}^T(W-\tilde{W})\right)\le {\left|\right|\tilde{W}\left|\right|}_F{\left|\right|W\left|\right|}_F-{\left|\right|\tilde{W}\left|\right|}_F^2---\left(47\right)$

式中，||·||_F表示伏柔宾纽斯(Frobenius)范数。

将式(47)代入式(46)，可得：

$\left(\begin{array}{c}{\stackrel{·}{V}}_3\le -k_{\min }{\left|\right|\xi \left|\right|}^2+{ϵ}_N\left|\right|\xi \left|\right|+\gamma \left|\right|\xi \left|\right|\left[{\left|\right|\tilde{W}\left|\right|}_F({\left|\right|W\left|\right|}_F-{\left|\right|\tilde{W}\left|\right|}_F)\right]\\ \le -\left|\right|\xi \left|\right|\left[k_{\min }\right|\left|\xi \right||-{ϵ}_N+\gamma {\left|\right|\tilde{W}\left|\right|}_F({\left|\right|\tilde{W}\left|\right|}_F-W_M)]\end{array}\right)---\left(48\right)$

式中，ε_N为逼近误差的上界，k_min为可调参数矩阵k的最小特征值，W_M为权重 矩阵W的最大值元素。

考虑到如下等式：

$\gamma {\left|\right|\tilde{W}\left|\right|}_F({\left|\right|\tilde{W}\left|\right|}_F-W_M)=\gamma {({\left|\right|\tilde{W}\left|\right|}_F-\frac{1}{2}{W}_M)}^2-\frac{\gamma }{4}{W}_M^2---\left(49\right)$

若使得则须有以下不等式成立：

$\left|\right|\xi \left|\right|>\frac{{ϵ}_N+\frac{1}{4}\gamma W_M^2}{k_{\min }},$或${\left|\right|\tilde{W}\left|\right|}_F>\frac{1}{2}{W}_M+\sqrt{\frac{W_M^2}{4}+\frac{{ϵ}_N}{\gamma }}---\left(50\right)$

若则有||ξ||和一致最终有界，从||ξ||的收敛性可得：航迹跟踪精度 与神经网络逼近误差上界ε_N、可调参数矩阵k有关。

由此，通过上述的神经网络逼近器能够精确估计不确定的飞艇非线性动力学模 型。

本发明的有益技术效果：

1)该方法直接基于飞艇的非线性动力学模型设计航迹控制律，考虑了各项非线 性因素以及多变量耦合作用，克服了线性化模型仅适于平衡态的局限性，能够高精 度跟踪任意给定的参数化指令航迹。

2)该方法将复杂的非线性航迹控制系统分解成两个不超过系统阶数的子系统， 然后为每个子系统设计Lypaunov函数和中间虚拟控制量，通过确定适当的虚拟反 馈，使得系统的前面状态达到渐近稳定，一直“后推”至整个系统，确保了整个系 统的渐近稳定。

3)该方法不需要精确已知对飞艇动力学模型，采用神经网络逼近器估计未知的 飞艇动力学模型，提高了系统的适应性和控制精度。

控制工程师在应用过程中可以根据实际飞艇给定任意指令航迹，并将由该方法 得到的控制量传输至执行机构实现航迹控制功能。

附图说明

图1为本发明所述飞艇航迹控制系统结构图

图2为本发明所述飞艇三维航迹跟踪控制方法步骤流程图

图3为本发明所述飞艇坐标系及运动参数定义

图4为本发明所述神经网络结构图

图5为本发明所述飞艇三维航迹跟踪控制结果

图6为本发明所述飞艇三维航迹跟踪控制误差

图7为神经网络逼近结果

图中符号说明如下：

η η＝[x,y,z,θ,ψ,φ]^T为飞艇航迹，其中x、y、z、θ、ψ、φ分别为实际航迹的 x坐标、y坐标、z坐标、俯仰角、偏航角和滚转角；

η_d η_d＝[x_d,y_d,z_d,θ_d,ψ_d,φ_d]^T为指令航迹，其中x_d、y_d、z_d、θ_d、ψ_d和φ_d分别为 指令x坐标、指令y坐标、指令z坐标、指令俯仰角、指令偏航角和指令滚转角；

o_ex_ey_ez_e o_ex_ey_ez_e表示地面坐标系；

o_bx_by_bz_b o_bx_by_bz_b表示飞艇体坐标系；

e e＝[x_e,y_e,z_e,θ_e,ψ_e,φ_e]^T为航迹控制误差量，分别为航迹控制的x坐标误差、y坐 标误差、z坐标误差；

u u为系统控制量；

f(x) f(x)为飞艇的不确定动力学模型。

具体实施方式

本发明“一种飞艇三维航迹跟踪的反步神经网络控制方法”，其具体步骤如下：

步骤一：给定指令航迹

给定指令航迹为：

η_d＝[x_d,y_d,z_d,θ_d,ψ_d,φ_d]^T＝[(3t)m,(0.93t)m,10m,0rad,0.3rad,0rad]^T，x_d、y_d、z_d、 θ_d、ψ_d和φ_d分别为指令x坐标、指令y坐标、指令z坐标、指令俯仰角、指令偏航 角和指令滚转角；

步骤二：航迹控制误差量计算

计算指令航迹与实际航迹之间的航迹控制误差量：

e＝η_d-η＝[x_d-x,y_d-y,z_d-z,θ_d-θ,ψ_d-ψ,φ_d-φ]^T，

其中，η＝[x,y,z,θ,ψ,φ]^T为实际航迹，x、y、z、θ、ψ、φ分别为实际航迹 的x坐标、y坐标、z坐标、俯仰角、偏航角和滚转角，为连续变化值。

初始航迹为：

η₀＝[x₀,y₀,z₀,θ₀,ψ₀,φ₀]^T＝[100m,-200m,5m,0.02rad,0.02rad,0.1rad]^T。

初始速度：

V₀＝[u₀,v₀,w₀,p₀,q₀,r₀]^T＝[8m/s,0m/s,0m/s,0rad/s,0rad/s,0rad/s]^T

步骤三：设计航迹控制律：

1)建立飞艇动力学模型

飞艇空间运动的数学模型可表示为：

$\stackrel{·}{\eta }=J\left(\eta \right)=\left(\begin{array}{cc}{J}_1& 0_{3\times 3}\\ 0_{3\times 3}& J_2\end{array}\right)V---\left(51\right)$

$M\stackrel{·}{V}=\bar{N}+\bar{G}+\bar{\tau }---\left(52\right)$

式中

$J_1=\left(\begin{array}{ccc}\cos \psi \cos \theta & \cos \psi \sin \theta \sin \phi -\sin \psi \cos \phi & \cos \psi \sin \theta \cos \phi +\sin \psi \sin \phi \\ \sin \psi \cos \theta & \sin \psi \sin \theta \sin \phi +\cos \psi \cos \phi & \sin \psi \sin \theta \cos \phi -\cos \psi \sin \phi \\ -\sin \theta & \cos \theta \sin \phi & \cos \theta \cos \phi \end{array}\right)---\left(53\right)$

$J_2=\left(\begin{array}{ccc}0& \cos \phi & -\sin \phi \\ 0& \sec \theta \sin \phi & \sec \theta \cos \phi \\ 1& \tan \theta \sin \phi & \tan \theta \cos \phi \end{array}\right)---\left(54\right)$

$M=\left(\begin{array}{cccccc}m+m_{11}& 0& 0& 0& {\mathrm{mz}}_G& -{\mathrm{my}}_G\\ 0& m+m_{22}& 0& -{\mathrm{mz}}_G& 0& {\mathrm{mx}}_G\\ 0& 0& m+m_{33}& {\mathrm{my}}_G& -{\mathrm{mx}}_G& 0\\ 0& -{\mathrm{mz}}_G& {\mathrm{my}}_G& I_x+I_{11}& 0& 0\\ {\mathrm{mz}}_G& 0& -{\mathrm{mx}}_G& 0& I_x+I_{22}& 0\\ 0& {\mathrm{mx}}_G& 0& -I_{\mathrm{xz}}& 0& I_x+I_{33}\end{array}\right)---\left(55\right)$

$\bar{G}=\left(\begin{array}{c}(B-G)\sin \theta \\ (G-B)\cos \theta \sin \phi \\ (G-B)\cos \theta \cos \phi \\ y_{G}G\cos \theta \cos \phi -z_{G}G\cos \theta \sin \phi \\ -x_{G}G\cos \theta \cos \phi -z_{G}G\sin \theta \\ x_{G}G\cos \theta \sin \phi +y_{G}G\sin \theta \end{array}\right)---\left(56\right)$

$\tau =\left(\begin{array}{c}{\tau }_u\\ {\tau }_v\\ {\tau }_w\\ {\tau }_l\\ {\tau }_m\\ {\tau }_n\end{array}\right)=\left(\begin{array}{c}T\cos \mu \cos \upsilon \\ T\sin \mu \\ T\cos \mu \sin \upsilon \\ T\cos \mu \sin \upsilon l_y-T\sin {\mathrm{\mu l}}_z\\ T\cos \mu \cos \upsilon l_z-T\cos \mu \sin \upsilon l_x\\ T\sin {\mathrm{\mu l}}_x-T\cos \mu \cos \upsilon l_y\end{array}\right)---\left(57\right)$

$\bar{N}={[N_u,N_v,N_w,N_p,N_q,N_r]}^T---\left(58\right)$

其中

N_u＝(m+m₂₂)vr-(m+m₃₃)wq+m[x_G(p²+r²)-y_Gpq-z_Gpr]

                                               (59)

+QV^2/3(-C_X cosαcosβ+C_Y cosαsinβ+C_Z sinα)

N_v＝(m+m₃₃)wp-(m+m₁₁)ur-m[x_Gpq-y_G(p²+r²)+z_Gqr]

                                        (60)

+QV^2/3(C_X sinβ+C_Y cosβ)

N_w＝(m+m₂₂)vp-(m+m₁₁)uq-m[x_Gpr+y_Gqr-z_G(p²+q²)]

                                (61)

+QV^2/3(-C_X sinαsinβ+C_Y sinαcosβ-C_Z cosα)

N_p＝[(I_y+m₅₅)-(I_z+I₆₆)]qr+I_xzpq-I_xypr-I_yz(r²-q²)+

                          (62)

[mz_G(ur-wp)+y_G(uq-vp)]+QVC_l

N_q＝[(I_z+m₆₆)-(I_x+I₄₄)]pr+I_xyqr-I_yzpq-I_xz(p²-r²)

                                           (63)

+m[x_G(vp-uq)-z_G(wp-vr)]+QVC_m

N_r＝[(I_y+m₅₅)-(I_x+I₄₄)]pq-I_xzqr-I_xy(q²-p²)+I_yzpr

                             (64)

+m[y_G(wq-vr)-x_G(ur-wp)]+QVC_n

式中，m为飞艇质量，m₁₁、m₂₂、m₃₃为附加质量,I₁₁、I₂₂、I₃₃为附加惯量；Q为 动压，α为迎角,β为侧滑角，C_X、C_Y、C_Z、C_l、C_m、C_n为气动系数；I_x、I_y、 I_z分别为绕o_bx_b、o_by_b、o_bz_b的主惯量；I_xy、I_xz、I_yz分别为关于平面o_bx_by_b、o_bx_bz_b、 o_by_bz_b的惯量积；T为推力大小，μ为推力矢量与o_bx_bz_b面之间的夹角，规定其在 o_bx_bz_b面之左为正，υ为推力矢量在o_bx_bz_b面的投影与o_bx_b轴之间的夹角，规定其投 影在o_bx_b轴之下为正；l_x、l_y、l_z表示推力作用点距原点o_b的距离。

式(52)为关于广义速度V的表达式，需要将其变换为关于广义坐标η的表达式。

由式(51)可得：

$V=J^{-1}\left(\eta \right)\stackrel{·}{\eta }=R\left(\eta \right)\stackrel{·}{\eta }=\left(\begin{array}{cc}A& 0_{3\times 3}\\ 0_{3\times 3}& B\end{array}\right)\stackrel{·}{\eta }---\left(65\right)$

式中，J^-1(η)为J(η)的逆矩阵，

$A=\left(\begin{array}{ccc}\cos \psi \cos \theta & \sin \psi \cos \theta & -\sin \theta \\ \cos \psi \sin \theta \sin \phi -\sin \psi \cos \phi & \sin \psi \sin \theta \sin \phi +\cos \psi \cos \phi & \cos \theta \sin \phi \\ \cos \psi \sin \theta \cos \phi +\sin \psi \sin \phi & \sin \psi \sin \theta \cos \phi -\cos \psi \sin \phi & \cos \theta \cos \phi \end{array}\right)---\left(66\right)$

$B=\left(\begin{array}{ccc}0& -\sin \theta & 1\\ \cos \phi & \cos \theta \sin \phi & 0\\ -\sin \phi & \cos \theta \cos \phi & 0\end{array}\right)---\left(67\right)$

对式(65)微分，可得

$\stackrel{·}{V}=\stackrel{·}{R}\stackrel{·}{\eta }+R\stackrel{··}{\eta }---\left(68\right)$

式中

$\stackrel{·}{R}=\left(\begin{array}{cc}\stackrel{·}{A}& 0_{3\times 3}\\ 0_{3\times 3}& \stackrel{·}{B}\end{array}\right)---\left(69\right)$

式(68)左乘可得

$R^{T}M\stackrel{·}{V}=R^{T}M\stackrel{·}{R}\stackrel{·}{\eta }+R^T\mathrm{MR}\stackrel{··}{\eta }---\left(70\right)$

综合式(52)、式(68)以及式(70)可得:

$M_{\eta }\left(\eta \right)\stackrel{··}{\eta }+N_{\eta }(\eta ,\stackrel{·}{\eta })\stackrel{·}{\eta }+G_{\eta }\left(\eta \right)=u---\left(71\right)$

式中

M_η(η)＝R^TMR             (72)

$N_{\eta }(\eta ,\stackrel{·}{\eta })=R^{T}M\stackrel{·}{R}---\left(73\right)$

$G_{\eta }\left(\eta \right)=-R^T(\bar{N}+\bar{G})---\left(74\right)$

u＝R^Tτ            (75)

令x₁＝η，则动力学模型(71)可写为如下形式：

$\left(\begin{array}{c}{\stackrel{·}{x}}_1=x_2\\ {\stackrel{·}{x}}_2=M_{\eta }^{-1}u-M_{\eta }^{-1}{N}_{\eta }\left(\eta \right)\stackrel{·}{\eta }-G_{\eta }\left(\eta \right)\end{array}\right)---\left(76\right)$

以式(76)所描述的数学模型为被控对象，采用反步方法设计航迹控制律。

本实施例中的飞艇参数见表1。

表1 飞艇参数

2)设计控制律

根据指令航迹与实际航迹之间的航迹控制误差量e，定义如下虚拟量：

                                                       $a_1=-k_{1}e+\stackrel{·}{\eta }d$

                                                                (77)

式中，α₁为虚拟量，k₁为可调的控制参数。

定义虚拟量α₁与x₂之间的误差e′：

e′＝x₂-α₁             (78)

式(78)对时间微分并将式(76)代入，可得：

$e^{\prime }={\stackrel{·}{x}}_2-{\stackrel{·}{a}}_1=M_{\eta }^{-1}u-M_{\eta }^{-1}{N}_{\eta }\left(\eta \right)\stackrel{·}{\eta }-G_{\eta }\left(\eta \right)+k_1\stackrel{·}{e}-{\stackrel{··}{\eta }}_d---\left(79\right)$

令

$f\left(x\right)=-M_{\eta }^{-1}{N}_{\eta }\left(\eta \right)\stackrel{·}{\eta }-G_{\eta }\left(\eta \right)---\left(80\right)$

则式(80)可表示为：

${\stackrel{·}{e}}^{\prime }=f\left(x\right)+M_{\eta }^{-1}u+k_1\stackrel{·}{e}-{\stackrel{··}{\eta }}_d---\left(81\right)$

其中，k₁取值为12。

选取Lyapunov函数V₁

$V_1=\frac{1}{2}{e}^{T}e+\frac{1}{2}{e}^{\prime T}{e}^{\prime }---\left(82\right)$

式(82)对时间微分并将式(81)代入，可得：

${\stackrel{·}{V}}_1=-k_1{e}^{T}e+e^{\prime T}(f\left(x\right)+M_{\eta }^{-1}u+k_1\stackrel{·}{e}-{\stackrel{··}{\eta }}_d)---\left(83\right)$

根据式(83)，设计如下控制律：

$u=M_{\eta }(-e-k_1\stackrel{·}{e}-k_2{e}^{\prime }+{\stackrel{··}{\eta }}_d-f\left(x\right))---\left(84\right)$

其中，k₂取值为10。

由于在实际飞行过程中难以对飞艇进行精确建模，f(x)为未知函数，难以根据 式(35)进行控制律解算，因此，须采用f(x)的估计值对航迹控制律式(35)进行 解算；

其中，在步骤四中所述的设计神经网络逼近器，其设计方法为：

1)采用神经网络逼近未知函数f(x),则有：

f(x)＝w^Th(x)+ε         (37)

式中，w为神经网络的权重向量，ε为逼近误差，h(x)＝[h_i(x)]^T，h_i(x)为高斯 基函数，下标i表示第i个高斯基函数；

2)选择输入输出变量

令航迹控制误差e及其变化率实际航迹η及其变化率为神经网络逼近器的 输入变量，令估计值为神经网络逼近器的输出变量。

3)设计神经网络结构

神经网络结构包括输入层、隐层和输出层，如图4所示。

输入层：选取神经网络的输入变量为$x={\left(\begin{array}{cccc}e& \stackrel{·}{e}& \eta & \stackrel{·}{\eta }\end{array}\right)}^T.$

隐层：选取高斯函数作为隐层节点的基函数

$h_i\left(x\right)=\exp \left(\frac{{\left|\right|x-c\left|\right|}^2}{{2\sigma }_i^2}\right)---\left(86\right)$

其中，

$c=\left(\begin{array}{ccccccc}-1.5& -1& -0.5& 0& 0.5& 1& 1.5\\ -1.5& -1& -0.5& 0& 0.5& 1& 1.5\\ -1.5& -1& -0.5& 0& 0.5& 1& 1.5\\ -1.5& -1& -0.5& 0& 0.5& 1& 1.5\end{array}\right),{\sigma }_i=2.$

输出层：神经网络的输出为

$\hat{f}\left(x\right)={\hat{w}}^{T}h\left(x\right)---\left(87\right)$

其中，的取值取为$0.15\times \left(\begin{array}{ccccccc}1& 1& 1& 1& 1& 1& 1\\ 1& 1& 1& 1& 1& 1& 1\\ 1& 1& 1& 1& 1& 1& 1\\ 1& 1& 1& 1& 1& 1& 1\\ 1& 1& 1& 1& 1& 1& 1\\ 1& 1& 1& 1& 1& 1& 1\end{array}\right).$

由此，通过上述的神经网络逼近器能够精确估计不确定的飞艇非线性动力学模 型。

实施例中的飞艇三维航迹跟踪控制结果如图5-图7所示。图5给出了飞艇三维 航迹跟踪控制结果，由图5可得：飞艇由初始位置出发，能够准确地跟踪指令航迹， 验证了本发明所提出的航迹跟踪控制方法的有效性；图6给出了航迹跟踪控制误差， 由图6可得，本发明所提出的航迹控制方法能够高精度地跟踪给定的指令航迹。图7 给出了神经网络逼近结果，由图7可得，本发明所设计的神经网络逼近器能够精确 地估计不确定的飞艇动力学模型。

以上包含了本发明优选实施例的说明，这是为了详细说明本发明的技术特征， 并不是想要将发明内容限制在实施例所描述的具体形式中，依据本发明内容主旨进 行的其他修改和变型也受本专利保护。本发明内容的主旨是由权利要求书所界定， 而非由实施例的具体描述所界定。