CN103116280A - 一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法 - Google Patents

Abstract

本发明公开了一种存在网络随机延迟的微小型无人飞行器纵向控制方法，属于飞行控制技术领域，该方法首先建立无人飞行器纵向系统模型并设计不存在网络延迟的飞行控制系统；然后在不同网络随机延迟条件下，进行飞行控制系统随机鲁棒性分析与设计；确定飞行控制系统的增益调度策略；最后闭环六自由度非线性蒙特卡洛仿真验证。本发明缓解了网络延迟导致的飞行控制系统品质恶化，能够增强分布式系统的稳定性；利用随机鲁棒分析与设计方法以及基于线性插值的增益调度策略，得到的控制律具有简单便于工程实现的优点，能够在大网络时延变化的情况下较好的维持原控制系统的控制品质。

Description

一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法

技术领域

本发明涉及存在变分布网络随机延迟的微小型无人飞行器纵向控制，属于飞行控制技术领域，具体涉及一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法。

背景技术

当前已经有多达三十多个国家投入大量人力和财力从事无人机的研究和生产。经过二十年的发展，该项技术已经比较成熟，在军民各个领域发挥着作用，尽管如此，单架无人机在遂行任务时存在着一些问题，例如单架无人机可能受到传感器的数量限制，不能从多角度全方位的对目标区域进行观察，面临大面积搜索任务时，不能有效的覆盖整个搜索区域；如果执行的是救援任务，单架无人机在载荷方面受到限制，往往影响整个救援的效能，带来更大损失，另外，一旦单架无人机出现故障，必须立即中断任务返回，可能会延误救援时机。

针对单架无人机的上述缺点，近些年提出了编队飞行控制的概念并取得了一定得研究成果，其中网络飞行控制系统设计是编队飞行控制是否成功的决定性因素，目前无人机群编队控制主要采用集中式控制结构：将编队成员分成一个长机和若干僚机，长机将自身的制导回路信息通过无线通信网络发送给僚机，僚机接收到长机的信息后与自身的制导回路信息比较，按照某种队形结构计算自身的制导指令进而与长机以一定队形结构飞行。由于无人机制导回路的惯性大，带宽较窄，而按照最优随机控制的观点，网络延迟往往可以等效为频率相对较高的白噪声，因此集中式控制体制下的网络延迟对编队队形控制的影响不大。近些年来，随着分布式编队控制的发展，要求编队成员不仅要在制导回路进行协调，甚至是控制回路也要进行协调，这就对网络飞行控制系统设计提出了较高的要求，由于控制回路的带宽相对制导回路宽，很有可能将网络延迟近似的噪声纳入到带宽之内，故网络延迟可能导致飞控系统控制品质恶化。

网络飞行控制系统设计的概念是：针对无人机编队飞行过程中的无线通信网络延迟、丢包和时序错乱情况，设计飞行控制系统，满足系统对超调量和调节时间等性能指标的要求。大量的网络测试研究表明网络状况是随时间变化而变化的，在不同的网络状态下，采用与之相适应的状态反馈增益能提高网络控制系统的控制品质，因此依据不同的网络状况调度相适应的状态反馈增益是非常有必要的。增益调度的方法分为传统增益调度、基于LPV（线性变参数）的增益调度和基于回路成形的鲁棒增益调度三种，其中后两种方法有严格的数学证明和可量化的控制品质评价标准，但是增益调度结构复杂，不利于实现和在线调整；增益调度设计过程复杂，不易为工程设计人员掌握。

发明内容

本发明的目的是为了解决上述问题，针对无人飞行器编队飞行过程中不可避免的无线自组织网络随机延迟问题，引入了随机鲁棒分析与设计方法，设计了无人机的控制律，解决了传统的线性二次型调节器（LQR）控制对通信质量要求较高和鲁棒性较差的缺点，降低了对编队无线数据链更新率的要求；根据不同的网络随机延迟，提出了基于线性插值的增益调度策略，拓宽了线性控制系统的控制范围，使飞行器在存在较大网络随机延迟变化的情况下控制品质仍能满足要求。控制策略简单不含复杂数学运算，便于有实时性要求的情况下进行工程实现。采用这种控制策略可以降低对无人机编队无线数据链更新率的要求，进而能够增强无人机群编队的鲁棒性。

本发明提供的一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法，具体包括如下步骤：

步骤一：建立无人飞行器纵向系统模型并设计不存在网络延迟的飞行控制系统；

具体为：

（1）根据风洞吹风获得无人飞行器动力学参数和物理参数；

（2）建立飞行器纵向非线性动力学与运动学方程；

（3）解耦线性化，得到无人飞行器纵向线性状态方程；

（4）LQR设计控制律；

步骤二：在不同网络随机延迟条件下，进行飞行控制系统随机鲁棒性分析与设计；

步骤三：确定飞行控制系统的增益调度策略；

采用线性插值进行网络飞行控制系统的增益调度，假设网络延迟在τ₁和τ_i之间，网络随机延迟在τ₁～N(μ₁,σ₁)，τ₂～N(μ₂,σ₂)......τ_i～N(μ_i,σ_i)时分别利用随机鲁棒分析与设计方法设计的控制系统参数

控制器当前测量的网络延迟为τ_N，则实际采用的控制器反馈系数K_N为：

$K_N=\left\{\begin{array}{c}\frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_1}{{\mathrm{\&tau;}}_2-{\mathrm{\&tau;}}_1}\&CenterDot;(K_2^N-K_1^N)+K_1^N,{\mathrm{\&tau;}}_1\&le;{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_2\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_2}{{\mathrm{\&tau;}}_3-{\mathrm{\&tau;}}_2}\&CenterDot;(K_3^N-K_2^N)+K_2^N,{\mathrm{\&tau;}}_2<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_3\\ \&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_{i-1}}{{\mathrm{\&tau;}}_i-{\mathrm{\&tau;}}_{i-1}}\&CenterDot;(K_i^N-K_{i-1}^N)+K_{i-1}^N,{\mathrm{\&tau;}}_{i-1}<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_i\end{array}\right.$

步骤四：闭环六自由度非线性蒙特卡洛仿真验证。

本发明的优点在于：

（1）本发明在一定程度上缓解了网络延迟导致的飞行控制系统品质恶化，能够增强分布式系统的稳定性；其次控制系统的随机鲁棒分析与设计针对的网络传输模式是事件驱动的，这种传输模式有以下两个优点，一是降低了数据链设计的难度与成本，二是能较好的与现有民用网络协议兼容；

（2）利用随机鲁棒分析与设计方法设计了基于单一分布网络延迟的控制系统后，采用基于线性插值的增益调度策略设计了适用于大网络时延变化的控制系统，控制律具有简单便于工程实现的优点，能够在大网络时延变化的情况下较好的维持原控制系统的控制品质；

（3）本发明设计过程没有繁琐的数学推导，而是充分利用了计算机高速处理数据的能力，相对于以往的针对线性时滞系统依靠数学推导获得控制律的方法，该方法工程性与可操作性较强，加之采用的线性插值的增益调度策略，能够放宽线性控制系统的控制范围，对网络时延分布变化较大的系统有较好的控制品质。

附图说明

图1是本发明提供的存在变分布网络随机延迟的微小型无人飞行器纵向控制方法总体框图；

图2a～图2h是不存在网络随机延迟时LQR设计的控制律在飞行器纵向通道响应曲线；

图3a～图3c为τ₁～N(μ₁,σ₁ ²)的AID值和高度响应的蒙特卡洛仿真曲线；

图4a～图4c为τ₂～N(μ₂,σ₂ ²)的AID值和高度响应的蒙特卡洛仿真曲线；

图5a～5c为τ₃～N(μ₃,σ₃ ²)的AID值和高度响应的蒙特卡洛仿真曲线；

图6a～6d为τ₁～N(μ₁,σ₁ ²)的随机鲁棒设计结果；

图7a～7d为τ₂～N(μ₂,σ₂ ²)的随机鲁棒设计结果；

图8a～8d为τ₃～N(μ₃,σ₃ ²)的随机鲁棒设计结果；

图9为K₁在网络随机延迟服从τ₃～N(μ₃,σ₃ ²)分布时蒙特卡洛仿真曲线；

图10a～10b为正态分布闭环六自由度非线性蒙特卡洛仿真验证图；

图11a～11b为均匀分布闭环六自由度非线性蒙特卡洛仿真验证图。

具体实施方式

下面将结合附图和实施例对本发明作进一步的详细说明。

本发明是一种存在变分布网络随机延迟的微小型无人飞行器控制方法，流程如图1所示，包括以下几个步骤：

步骤一：建立无人飞行器纵向系统模型并设计不存在网络延迟的飞行控制系统；

具体为：

（1）根据风洞吹风获得无人飞行器动力学参数和物理参数；

无人飞行器动力学参数和物理参数可以根据实际风洞吹风得到，采用英美坐标系，具体为：

获取动力学参数包括纵向力和力矩系数：包括升力系数C_L0、C_La、C_Lq、

阻力系数C_D0、C_Da、C_Dq、

俯仰力矩系数C_m0、C_ma、C_mq、

其中，C_L0为攻角为0度时的升力系数，C_La为升力关于攻角的升力系数，C_Lq为升力关于俯仰角速度的升力系数，

为升力关于升降舵的升力系数，C_D0为攻角为0度时的阻力系数，C_Da为阻力关于攻角的阻力系数，C_Dq为阻力关于俯仰角速度的阻力系数，

为阻力关于升降舵的阻力系数，C_m0为攻角为0度时的俯仰力矩系数，C_ma为俯仰力矩关于攻角的力矩系数，C_mq为俯仰力矩关于俯仰角速度的力矩系数、

为俯仰力矩关于升降舵的力矩系数；

获取物理参数包括无人飞行器的质量m（kg）、平均几何弦长C_A（m）、翼展b（m）、飞行器参考面积S_w（m²）、X轴转动惯量I_X（kg·m²）、Y轴转动惯量I_Y（kg·m²）、Z轴转动惯量I_Z（kg·m²）、惯量积I_XZ（kg·m²）、舵机时间常数T_δ(s)、舵机放大系数K_δ、发动机时间常数T_t(s)、发动机放大系数K_t、怠速推力t_A（N）、最大推力t_max（N）和最小推力t_min（N）。

（2）建立飞行器纵向非线性动力学与运动学方程；

选取无人机飞行状态向量为：

$\stackrel{-}{X}\left(t\right)={\left[\begin{array}{ccccccccc}V& \mathrm{\&alpha;}& q& \mathrm{\&theta;}& x& h& m& {\mathrm{\&delta;}}_e& {\mathrm{\&delta;}}_t\end{array}\right]}^T$

其中，V表示速度、α表示攻角、q表示俯仰角速度、θ表示俯仰角、x表示东向位置、h表示高度、m表示质量、δ_e表示升降舵偏角、δ_t表示发动机推力。

建立无人飞行器非线性六自由度运动学与动力学方程如下：

$\left\{\begin{array}{c}\stackrel{\&CenterDot;}{V}=({\mathrm{\&delta;}}_t\&CenterDot;\cos \left(\mathrm{\&alpha;}\right)-D-m\&CenterDot;g\&CenterDot;\sin (\mathrm{\&theta;}-\mathrm{\&alpha;}))/m\\ \stackrel{\&CenterDot;}{\mathrm{\&alpha;}}=(-{\mathrm{\&delta;}}_t\&CenterDot;\sin \left(\mathrm{\&alpha;}\right)-L+m\&CenterDot;V\&CenterDot;q+m\&CenterDot;g\&CenterDot;\cos (\mathrm{\&theta;}-\mathrm{\&alpha;\&alpha;}))/(m\&CenterDot;V)\\ \stackrel{\&CenterDot;}{q}=M/I_Y\\ \stackrel{\&CenterDot;}{\mathrm{\&theta;}}=q\\ \stackrel{\&CenterDot;}{x}=V\&CenterDot;\cos \left(\mathrm{\&mu;}\right)\\ \stackrel{\&CenterDot;}{h}=V\&CenterDot;\sin \left(\mathrm{\&mu;}\right)\\ \stackrel{\&CenterDot;}{m}=K_m\&CenterDot;{\mathrm{\&delta;}}_T\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e=(-{\mathrm{\&delta;}}_e+K_{\mathrm{\&delta;}}\&CenterDot;{\mathrm{\&delta;}}_e^{*})/T_{\mathrm{\&delta;}}\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e=(-{\mathrm{\&delta;}}_t+K_t\&CenterDot;{\mathrm{\&delta;}}_t^{*})/T_t\\ \mathrm{\&mu;}=\mathrm{\&theta;}-\mathrm{\&alpha;}\end{array}\right.---\left(1\right)$

其中：

表示速度的导数、

表示攻角的导数、俯仰角速度的导数、

表示俯仰角的导数、

表示东向位置导数、

表示高度的导数、

表示质量消耗率、表示升降舵偏角的导数、表示发动机推力的导数、

表示升降舵偏角指令、表示油门指令、μ为航迹倾角、K_m为发动机耗油率、δ_T为发动机推力、T_δ为舵机时间常数、T_t为发动机时间常数、K_δ为舵机放大系数、K_t为发动机放大系数，I_Y为Y轴转动惯量，式（1）简记为：

其中

为飞行状态向量的导数，

为飞行控制向量。

$\left\{\begin{array}{c}{c}_L=C_{L0}+C_{\mathrm{L\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{\mathrm{Lq}}\&CenterDot;q+C_{L{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e\\ C_D=C_{D0}+C_{\mathrm{D\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{\mathrm{Dq}}\&CenterDot;q+C_{D{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e\\ C_m=C_{m0}+C_{\mathrm{m\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{\mathrm{mq}}\&CenterDot;q+C_{m{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e\\ L=\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_L\\ D=\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_D\\ M=\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_m\&CenterDot;C_A\end{array}\right.---\left(2\right)$

其中：C_L为总升力系数、C_D为总阻力系数、C_m为总俯仰力矩系数、ρ为空气密度、L为总升力、D为总阻力、M为总俯仰力矩、S_w为飞行器参考面积、C_A为飞行器平均几何弦长、b为翼展。

（3）解耦线性化；

假设无人飞行器处于定高稳定平飞状态，则必有升力-重力平衡、推力-阻力平衡以及俯仰力矩平衡，因此下面三个等式成立：

$\left\{\begin{array}{c}\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{L0}+C_{\mathrm{L\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{L{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e)=m\&CenterDot;g\\ \frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{D0}+C_{\mathrm{D\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{D{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e)=t_A\\ \frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{m0}+C_{\mathrm{m\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{m{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e)\&CenterDot;C_A=0\end{array}\right.---\left(3\right)$

解方程（3），得到无人飞行器的配平攻角α₀，配平升降舵偏角δ_e0，配平飞行速度V₀，在该状态点利用小扰动线性化原理可得无人飞行器纵向线性状态方程为：

$\left[\begin{array}{c}\stackrel{\&CenterDot;}{V}\\ \stackrel{\&CenterDot;}{\mathrm{\&alpha;}}\\ \stackrel{\&CenterDot;}{q}\\ \stackrel{\&CenterDot;}{\mathrm{\&theta;}}\\ \stackrel{\&CenterDot;}{h}\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_t\end{array}\right]=A\left[\begin{array}{c}V\\ \mathrm{\&alpha;}\\ q\\ \mathrm{\&theta;}\\ h\\ {\mathrm{\&delta;}}_e\\ {\mathrm{\&delta;}}_t\end{array}\right]+B\left[\begin{array}{c}{\mathrm{\&delta;}}_e^{*}\\ {\mathrm{\&delta;}}_t^{*}\end{array}\right]---\left(4\right)$

其中A为纵向线性系统状态矩阵，B为纵向线性系统控制矩阵；

（4）LQR设计控制律；

针对纵向：选取Q₁和R₁矩阵，采用线性二次型调节器LQR设计纵向控制律K，K为2×7的矩阵。用k_ij表示K中的第i行第j列元素，为使控制系统结构尽量简化，速度反馈只加入发动机油门控制（k₂₂＝k₂₃＝k₂₄＝k₂₅＝0.0）中，而不加入升降舵控制（k₁₁＝0.0）中，无人飞行器舵机一般采用电动舵机，在飞行过程中无法准确测量舵面偏角，因此k₁₆＝k₁₇＝k₂₆＝k₂₇＝0.0，因为速度反馈不参与升降舵控制，而发动机油门调节只依赖速度，故应增加其他状态量对升降舵控制的反馈力度，最终选择以下控制器结构：

$\left\{\begin{array}{c}{\mathrm{\&delta;}}_e^{*}={\mathrm{\&delta;}}_{e0}+k_{12}\&CenterDot;(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0)-k_{13}\&CenterDot;q-k_{14}\&CenterDot;(\mathrm{\&theta;}-{\mathrm{\&alpha;}}_0)-k_{15}\&CenterDot;(h-h^{*})\\ {\mathrm{\&delta;}}_t^{*}=t_A-k_{21}\&CenterDot;(V-V_0)\end{array}---\left(5\right)\right.$

其中：δ_e0为配平升降舵偏角，α₀为无人飞行器的配平攻角，h^*为高度指令，V₀为配平飞行速度。微小型飞行器由于结构较小无法安装精确的攻角传感器，故将攻角反馈换算为过载反馈，由于微小型飞行器在飞行过程中，机翼和水平安定面是主要升力面，升降舵面产生的升力相对较小，主要是提供俯仰力矩诱导主升力面发生变化，因此：

$\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{L0}+C_{\mathrm{L\&alpha;}}\&CenterDot;(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0))+{\mathrm{\&delta;}}_t\&CenterDot;\sin \left(\mathrm{\&alpha;}\right)=(N_Z-1)\&CenterDot;m\&CenterDot;g---\left(6\right)$

其中：N_Z是无人飞行器的法向过载，C_L0相对C_La较小，可以忽略，δ_t·sin(α)相对 $\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_{\mathrm{L\&alpha;}}\&CenterDot;(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0)$ 较小也可以忽略，代入数据得：

$(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0)=(N_Z-1)\&CenterDot;m\&CenterDot;g/(\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_{\mathrm{L\&alpha;}})---\left(7\right)$

联立式（5）、式（7）可得实际飞行中的控制器结构：

$\left\{\begin{array}{c}{\mathrm{\&delta;}}_e^{*}={\mathrm{\&delta;}}_{e0}-k_1\&CenterDot;(N_z-1)-k_2\&CenterDot;q-k_3\&CenterDot;(\mathrm{\&theta;}-{\mathrm{\&alpha;}}_0)-k_4\&CenterDot;(h-h^{*})\\ {\mathrm{\&delta;}}_t^{*}=t_A-k_5\&CenterDot;(V-V_0)\end{array}\right.---\left(8\right)$

其中：k₁k₂k₃k₄k₅表示按照（5）和（7）换算后的最终应用的控制器反馈系数。

步骤二：在不同网络随机延迟条件下，进行飞行控制系统随机鲁棒性分析与设计；

（1）飞行控制系统随机鲁棒性分析。

网络随机延迟的定义是：在传感器、执行器和控制器等多个节点通过网络交换数据时，由于网络带宽有限且网络中的数据流量变化不规则，不可避免的造成数据碰撞、多路传输、连接中断和网络拥塞等现象。本发明主要采用相互独立的随机延迟作为网络飞行控制系统的主要延迟类型，假设网络随机延迟τ_i服从均值μ_i，方差σ_i的正态分布，即：

τ_i∈N(μ_i,σ_i),i＝1,2,...,N    （9）

由于控制系统采用事件驱动的方式处理数据，故延迟发生的时刻并不均匀，无人飞行器非线性系统方程由原来的 $\stackrel{\&CenterDot;}{\stackrel{\&OverBar;}{X}}\left(t\right)=f(\stackrel{\&OverBar;}{X}\left(t\right),u\left(t\right))$ 变为了：

$\stackrel{\&CenterDot;}{\stackrel{-}{X}}\left(t\right)=f(\stackrel{-}{X}\left(t\right),u(t-{\mathrm{\&tau;}}_i))---\left(10\right)$

其中：

为飞行状态向量，

为飞行状态向量的导数，u(t-τ_i)为加入控制延迟的飞行控制向量。网络延迟本身的随机性为随机鲁棒分析的应用提供了条件。对服从（9）式的网络随机延迟进行蒙特卡洛仿真考察（8）式控制律在网络延迟情况下对无人飞行器控制的鲁棒性，鲁棒性的衡量标准取为一定时间内实际高度与指令高度的平均偏差。

（2）飞行控制系统随机鲁棒性设计：

随机鲁棒设计过程包括控制器结构设计和现代优化算法两部分，控制器结构采用传统的PID控制，在步骤一中已叙述，现代优化算法采用粒子群算法，优化指标是飞行控制系统随机鲁棒性分析中的鲁棒性衡量标准，因为控制器反馈系数k₁,k₂,k₃,k₄,k₅变化是连续的，粒子群算法在处理多连续变量寻优相对遗传算法和蚁群算法具有优势。

采用标准粒子群算法求解该优化问题，问题的解对应于搜索空间中的一个粒子，每个粒子都有自己的位置和速度及一个由被优化函数决定的代价函数。各个粒子记忆、追随当前的最优粒子，每次迭代中，粒子通过跟随两个极值来更新自己的位置和速度：一个是粒子本身获得的最优解一个是整个粒子群中全部粒子历代搜索获得的最优解

标准粒子群算法中，速度更新和位置更新如式(11)所示：

$\left\{\begin{array}{c}{v}_{\mathrm{id}}^{k+1}={\mathrm{\&omega;v}}_i^k+c_1{r}_1(p_{\mathrm{id}}^k-x_{\mathrm{id}}^k)+c_2{r}_2(g_{\mathrm{id}}^k-x_{\mathrm{id}}^k)\\ x_{\mathrm{id}}^{k+1}=x_{\mathrm{id}}^k+v_{\mathrm{id}}^{k+1}\end{array}---\left(11\right)\right.$

式中，i∈N(1，n_s),n_s为粒子群中粒子的个数，N代表整数；d∈N(1,n)，n为解向量的维数，k为迭代次数，c₁和c₂为学习因子，初始值均为0.6，随迭代次数增加降为0.2；ω为惯量权重，初始值为0.9，随迭代次数增加降为0.4；r₁和r₂为[0,1]之间的随机数，判断粒子群中粒子优劣的标准是蒙特卡洛仿真获得的代价函数W，代价函数越小，在下一代该粒子存在的概率越大。

设定粒子最大的迭代次数，判断迭代是否到达最大迭代次数，如果到达，输出最好粒子代表的网络延迟τ_i下的控制系统参数

其中i代表不同的网络状况，否则，继续进行迭代。

针对不同的网络随机延迟分布τ_i利用随机鲁棒性设计方法设计对应的控制律，保证该控制律在对应的网络随机延迟分布情况下随机鲁棒性分析的鲁棒性最优。

步骤三：确定飞行控制系统的增益调度策略；

随机鲁棒设计的控制器只能保证网络随机延迟在一定范围内的鲁棒性能满足一定要求，当网络随机延迟变化较大时，控制器的鲁棒性能可能不满足指标要求。采用线性插值进行网络飞行控制系统的增益调度。假设网络延迟在τ₁和τ_i之间，网络随机延迟在τ₁～N(μ₁,σ₁)，τ₂～N(μ₂,σ₂)......τ_i～N(μ_i,σ_i)时分别利用随机鲁棒分析与设计方法设计控制系统参数

控制器当前测量的网络延迟为τ_N，则实际采用的控制器反馈系数K_N为：

$K_N=\left\{\begin{array}{c}\frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_1}{{\mathrm{\&tau;}}_2-{\mathrm{\&tau;}}_1}\&CenterDot;(K_2^N-K_1^N)+K_1^N,{\mathrm{\&tau;}}_1\&le;{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_2\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_2}{{\mathrm{\&tau;}}_3-{\mathrm{\&tau;}}_2}\&CenterDot;(K_3^N-K_2^N)+K_2^N,{\mathrm{\&tau;}}_2<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_3\\ \&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_{i-1}}{{\mathrm{\&tau;}}_i-{\mathrm{\&tau;}}_{i-1}}\&CenterDot;(K_i^N-K_{i-1}^N)+K_{i-1}^N,{\mathrm{\&tau;}}_{i-1}<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_i\end{array}\right.$

步骤四：闭环六自由度非线性蒙特卡洛仿真验证；

在存在网络随机延迟τ的情况下，对步骤二中随机鲁棒设计出的控制系统参数

和步骤三中设计的控制增益调度策略进行蒙特卡洛仿真验证，结果显示，采用随机鲁棒分析与设计后的控制系统参数

较步骤一得出的最终应用的控制器反馈系数k₁k₂k₃k₄k₅抗网络随机延迟能力增强，最大超调量获得明显改善；在线性增益调度的调节下，放大了控制律保证控制品质的范围，可以在网络随机延迟存在较大变化的情况下仍能维持控制器的鲁棒稳定性和鲁棒性能。所以采用本发明方法获得的控制律调度策略，相对于LQR设计出的控制律抗网络延迟的能力加强，把本发明方法设计出的控制律应用到无人机自动驾驶仪的纵向通道中，在无人机自主飞行过程中，能够在存在较大网络延迟的情况下维持无人机的正常飞行。

图1给出了本专利方法的总体流程图。本发明能够解决一类具有大变化网络随机延迟的网络飞行控制律设计与调节问题，近些年来，随着分布式编队控制的发展，要求编队成员不仅要在制导回路进行协调，甚至是控制回路也要进行协调，这就对网络飞行控制系统设计提出了较高的要求，由于控制回路的带宽相对制导回路宽，很有可能将网络延迟近似的噪声纳入到带宽之内，故网络延迟可能导致飞控系统控制品质恶化，本发明在一定程度上缓解了网络延迟导致的飞行控制系统品质恶化，能够增强分布式系统的稳定性；其次控制系统的随机鲁棒分析与设计针对的网络传输模式是事件驱动的，这种传输模式有以下两个优点，一是不需要设计定时精确的定时器，降低了数据链设计的难度与成本，二是能较好的与现有民用网络协议兼容。

实施例：

步骤一：建立无人飞行器系统模型并确定不存在网络延迟的飞行控制系统的控制律；

无人飞行器动力学参数和物理参数可以根据实际风洞吹风得到，采用英美坐标系，具体数值见下表：

（2）建立非线性六自由度动力学与运动学方程；

选取无人机飞行状态向量为：

$\stackrel{-}{X}\left(t\right)={\left[\begin{array}{ccccccccc}V& \mathrm{\&alpha;}& q& \mathrm{\&theta;}& x& h& m& {\mathrm{\&delta;}}_e& {\mathrm{\&delta;}}_t\end{array}\right]}^T$

建立如公式（1）所示的无人机六自由度动力学与运动学方程；

（3）解耦线性化；

解公式（3），得到无人飞行器的配平攻角α₀＝4.4821°，配平升降舵偏角δ_e0＝-2.2736°，配平飞行速度V₀＝41.0214m/s，在该状态点利用小扰动线性化原理可得无人飞行器纵向线性系统状态方程为：

$\left[\begin{array}{c}\stackrel{\&CenterDot;}{V}\\ \stackrel{\&CenterDot;}{\mathrm{\&alpha;}}\\ \stackrel{\&CenterDot;}{q}\\ \stackrel{\&CenterDot;}{\mathrm{\&theta;}}\\ \stackrel{\&CenterDot;}{h}\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_t\end{array}\right]=\begin{array}{}\end{array}\left[\begin{array}{ccccccc}-0.129607& -26.4232& 0& -9.81& 0& 2.40894& 0.0483686\\ -0.0109844& -5.57054& 1.00102& 0& 0& -0.321125& -0.000067\\ -0.0000225& -70.3391& -511.921& 0& 0& -54.0864& 0\\ 0& 0& 1& 0& 0& 0& 0\\ 0& -42& 0& 42& 0& 0& 0\\ 0& 0& 0& 0& 0& -20.0& 0\\ 0& 0& 0& 0& 0& 0& -5.0\end{array}\right]\left[\begin{array}{c}V\\ \mathrm{\&alpha;}\\ q\\ \mathrm{\&theta;}\\ h\\ {\mathrm{\&delta;}}_e\\ {\mathrm{\&delta;}}_t\end{array}\right]$

$+\left[\begin{array}{cc}0& 0\\ 0& 0\\ 0& 0\\ 0& 0\\ 0& 0\\ 20.0& 0\\ 0& 5.0\end{array}\right]\left[\begin{array}{c}{\mathrm{\&delta;}}_e^{*}\\ {\mathrm{\&delta;}}_t^{*}\end{array}\right]---\left(12\right)$

（4）LQR设计控制律；

针对纵向：选取Q₁和R₁矩阵如下：

$Q=\left[\begin{array}{ccccccc}0.01& & & & & & \\ & 131.3& & & & 0& \\ & & 820.82& & & & \\ & & & 131.3& & & \\ & & & & 0.0025& & \\ & 0& & & & 131.3& \\ & & & & & & 0.01\end{array}\right],$ $R=\left[\begin{array}{cc}131.3& 0\\ 0& 0.01\end{array}\right]---\left(13\right)$

采用线性二次型调节器设计纵向控制律：

$K=\left[\begin{array}{ccccccc}-0.00435& 0.133262& -0.016633& -2.4003& -0.00417& 0.443827& -0.00003\\ 0.46118& -4.6277& 0.100893& 56.3062& 0.147224& -0.0995& 0.417404\end{array}\right]---\left(14\right)$

根据公式（7）（8）（9）（10）的运算最终可得

$\left\{\begin{array}{c}{\mathrm{\&delta;}}_e^{*}={\mathrm{\&delta;}}_{e0}-k_1\&CenterDot;(N_z-1)-k_2\&CenterDot;q-k_3\&CenterDot;(\mathrm{\&theta;}-{\mathrm{\&alpha;}}_0)-k_4\&CenterDot;(h-h^{*})\\ {\mathrm{\&delta;}}_t^{*}=t_A-k_5\&CenterDot;(V-V_0)\end{array}---\left(15\right)\right.$

无人机的升降舵不接入速度反馈，油门控制只接入速度控制，由于除速度以外的其他状态量不能通过油门调节，故加强升降舵对其他状态量的反馈，最终选择如下控制器反馈系数：

k₁＝-0.48,k₂＝0.02632,k₃＝2.4003,k₄＝0.008,k₅＝-0.46118；

假设初始无人机飞行高度在300m，这时给高度指令h^*＝250m，控制计算机频率为50Hz，系统响应如图2所示，分别包括速度、攻角、俯仰角速度、俯仰角、高度、升降舵偏角、油门、法向过载的相应图。

取仿真时间100s内的平均指令偏差（Average Instruction Deviation，AID）作为系统控制品质的度量：

$\mathrm{AID}=\frac{1}{T_0}{\&Integral;}_0^{T_0}|h-h^{*}|\mathrm{dt},$ 其中T₀＝100s

此时系统的AID为6.4394m。

步骤二：在不同网络随机延迟条件下，进行飞行控制系统随机鲁棒性分析与设计；

（1）飞行控制系统随机鲁棒性分析

假设无人机控制系统的网络延迟可以分为两个等级，一类网络延迟属于短时间延迟，记为τ₁，服从均值μ₁为700ms，方差σ₁为100ms的标准正态分布；一类网络延迟属于中时间延迟，记为τ₂，服从均值μ₂为900ms，方差σ₂为100ms的标准正态分布；一类网络延迟属于长时间延迟，记为τ₃，服从均值μ₃为1100ms，方差σ₃为100ms的标准正态分布，即：

τ₁～N(μ₁,σ₁ ²),τ₂～N(μ₂,σ₂ ²),τ₃～N(μ₃,σ₃ ²)

对网络延迟服从上述分布的网络飞行控制系统进行蒙特卡洛仿真，控制律采用步骤一中设计的无人机纵向控制律，图3是τ₁～N(μ₁,σ₁ ²)的AID值和高度响应的蒙特卡洛仿真曲线，此时最大超调量达到20%，AID最大值为6.87；图4是τ₂～N(μ₂,_σ2 ²)的AID值和高度响应的蒙特卡洛仿真曲线，此时最大超调量达到35%，AID最大值为7.48；图5是τ₃～N(μ₃,σ₃ ²)的AID值和高度响应的蒙特卡洛仿真曲线，此时最大超调量达到45%，AID最大值为8.35，所以三种分布的网络随机延迟造成的控制系统的控制品质较差，控制系统的控制律是针对不含网络随机延迟的连续系统，采用LQR方法设计的，设计的控制律对含有网络延迟的飞行器纵向通道控制品质不高，因此采用随机鲁棒性设计方法进一步提高控制品质。

（2）飞行控制系统随机鲁棒性设计

控制器结构采用传统的PID控制，利用粒子群算法优化控制器参数，粒子群中粒子的个数取为30，解向量的维数n为5，迭代次数k为20，对控制器反馈参数k₁…k₅进行寻优，网络随机延迟τ₁～N(μ₁,σ₁ ²)的随机鲁棒设计结果如图6所示，图6给出了标准粒子群算法的代价收敛曲线，获得AID的最小值为6.26927，蒙特卡洛仿真中获得的最大AID为6.33，获得的最终控制系统参数 $K_1^N=[k_{11}^N,k_{12}^N,k_{13}^N,k_{14}^N,k_{15}^N]:$

$k_{11}^N=-0.267848,$ $k_{12}^N=0.032116,$ $k_{13}^N=2.63188,$ $k_{14}^N=0.00724909,$ $k_{15}^N=-0.290237$

控制器结构采用传统的PID控制，利用粒子群算法优化控制器参数，粒子群中粒子的个数取为30，解向量的维数n为5，迭代次数k为20，对控制器反馈参数k₁…k₅进行寻优，网络随机延迟τ₂～N(μ₂,σ₂ ²)的随机鲁棒设计结果如图7所示，图7给出了标准粒子群算法的代价收敛曲线，获得AID的最小值为6.47923，蒙特卡洛仿真中获得的最大AID为6.7，获得的最终控制系统参数 $K_2^N=[k_{21}^N,k_{22}^N,k_{23}^N,k_{24}^N,k_{25}^N]:$

$k_{21}^N=-0.303185,$ $k_{22}^N=0.0232822,$ $k_{23}^N=2.5781,$ $k_{24}^N=0.00706601,$ $k_{25}^N=-0.269634$

控制器结构采用传统的PID控制，利用粒子群算法优化控制器参数，粒子群中粒子的个数m取为30，解向量的维数n为5，迭代次数k为20，对控制器反馈参数k₁…k₅进行寻优，网络随机延迟τ₃～N(μ₃,σ₃ ²)的随机鲁棒设计结果如图8所示，图8给出了标准粒子群算法的代价收敛曲线，获得AID的最小值为6.51899，蒙特卡洛仿真中获得的最大AID为6.82，获得的最终控制系统参数 $K_3^N=[k_{31}^N,k_{32}^N,k_{33}^N,k_{34}^N,k_{35}^N]:$

$k_{31}^N=-0.221949,$ $k_{32}^N=0.032926,$ $k_{33}^N=2.44389,$ $k_{34}^N=0.00662643,$ $k_{35}^N=-0.231925$

步骤三：确定飞行控制系统的增益调度策略；

由于随机鲁棒设计的控制器只能保证网络随机延迟在一定范围内的鲁棒性能满足一定要求，例如K₁在网络随机延迟服从τ₃～N(μ₃,σ₃ ²)分布时蒙特卡洛仿真曲线如图9所示，此时的鲁棒性能已经下降到不可接受的地步，所以网络控制系统增益调度策略是有必要的。采用线性插值进行网络飞行控制系统的增益调度。假设网络延迟在600ms和1200ms之间，控制器测量的网络延迟为τ_N，则实际采用的控制器反馈系数K_N为：

$K_N=\left\{\begin{array}{c}\frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_1}{{\mathrm{\&tau;}}_2-{\mathrm{\&tau;}}_1}\&CenterDot;(K_2^N-K_1^N)+K_1^N,{\mathrm{\&tau;}}_1\&le;{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_2\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_2}{{\mathrm{\&tau;}}_3-{\mathrm{\&tau;}}_2}\&CenterDot;(K_3^N-K_2^N)+K_2^N,{\mathrm{\&tau;}}_2<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_3\end{array}\right.$

步骤四：闭环六自由度非线性蒙特卡洛仿真验证；

图10给出了高度指令和网络随机延迟满足正态分布

$\left\{\begin{array}{c}{h}^{*}=250m,{\mathrm{\&tau;}}_1~N({\mathrm{\&mu;}}_1,{{\mathrm{\&sigma;}}_1}^2),t\&le;50s\\ h^{*}=200m,{\mathrm{\&tau;}}_2~N({\mathrm{\&mu;}}_2,{{\mathrm{\&sigma;}}_2}^2),50s<t\&le;100s\\ h^{*}=150m,{\mathrm{\&tau;}}_3~N({\mathrm{\&mu;}}_3,{{\mathrm{\&sigma;}}_3}^2),100s<t\&le;150s\end{array}\right.$

时的控制律增益调度策略50次蒙特卡洛仿真曲线，可见控制品质可以满足要求。

图11给出了高度指令和网络随机延迟满足均匀分布

$\left\{\begin{array}{c}{h}^{*}=250m,{\mathrm{\&tau;}}_1~U(a,b),t\&le;50s\\ h^{*}=200m,{\mathrm{\&tau;}}_2~U(a,b),50s<t\&le;100s\\ h^{*}=150m,{\mathrm{\&tau;}}_3~U(a,b),100s<t\&le;150s\end{array}\right.$

时的控制律增益调度策略50次蒙特卡洛仿真曲线，其中a=600ms，b=1200ms，这是控制器增益调度中是最差的情况，可见控制品质可以满足要求。

图10和图11的结果的前提是网络随机延迟可以准确的测量，这时可以根据测量值采用响应的线性调度策略，实际情况中应研究网络随机延迟的估计方法。

Claims (4)

1.一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法，其特征在于，包括如下步骤：

步骤一：建立无人飞行器纵向系统模型并设计不存在网络延迟的飞行控制系统，具体为：

（1）根据风洞吹风获得无人飞行器动力学参数和物理参数；

（2）建立飞行器纵向非线性动力学与运动学方程；

（3）解耦线性化；

假设无人飞行器处于定高稳定平飞状态，则必有升力-重力平衡、推力-阻力平衡以及俯仰力矩平衡，因此下面三个等式成立：

$\left\{\begin{array}{c}\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{L0}+C_{\mathrm{L\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{L{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e)=m\&CenterDot;g\\ \frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{D0}+C_{\mathrm{D\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{D{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e)=t_A\\ \frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{m0}+C_{\mathrm{m\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{m{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e)\&CenterDot;C_A=0\end{array}\right.---\left(3\right)$

解方程（3），得到无人飞行器的配平攻角α₀，配平升降舵偏角δ_e0，配平飞行速度V₀，在该状态点利用小扰动线性化原理得无人飞行器纵向线性状态方程为：

$\left[\begin{array}{c}\stackrel{\&CenterDot;}{V}\\ \stackrel{\&CenterDot;}{\mathrm{\&alpha;}}\\ \stackrel{\&CenterDot;}{q}\\ \stackrel{\&CenterDot;}{\mathrm{\&theta;}}\\ \stackrel{\&CenterDot;}{h}\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_t\end{array}\right]=A\left[\begin{array}{c}V\\ \mathrm{\&alpha;}\\ q\\ \mathrm{\&theta;}\\ h\\ {\mathrm{\&delta;}}_e\\ {\mathrm{\&delta;}}_t\end{array}\right]+B\left[\begin{array}{c}{\mathrm{\&delta;}}_e^{*}\\ {\mathrm{\&delta;}}_t^{*}\end{array}\right]---\left(4\right)$

其中A为纵向线性系统状态矩阵，B为纵向线性系统控制矩阵；

（4）LQR设计控制律；

步骤二：在不同网络随机延迟条件下，进行飞行控制系统随机鲁棒性分析与设计；

（1）飞行控制系统随机鲁棒性分析：

采用相互独立的随机延迟作为网络飞行控制系统的主要延迟类型，假设网络随机延迟τ_i服从均值μ_i，方差σ_i的正态分布，即：

τ_i∈N(μ_i,σ_i),i＝1,2,...,N    （9）

无人飞行器非线性系统方程为：

$\stackrel{\&CenterDot;}{\stackrel{\&OverBar;}{X}}\left(t\right)=f(\stackrel{-}{X}\left(t\right),u(t-{\mathrm{\&tau;}}_i))---\left(10\right)$

其中：

为飞行状态向量，

为飞行状态向量的导数，u(t-τ_i)为加入网络随机延迟的飞行控制向量；

（2）飞行控制系统随机鲁棒性设计：

随机鲁棒性设计过程包括控制器结构设计和现代优化算法两部分，控制器结构设计采用传统的PID控制，现代优化算法采用标准粒子群算法，优化指标是飞行控制系统随机鲁棒性分析中的鲁棒性衡量标准；

采用标准粒子群算法求解该优化问题，问题的解对应于搜索空间中的一个粒子，每个粒子都有自己的位置和速度及一个由被优化函数决定的代价函数；各个粒子记忆、追随当前的最优粒子，每次迭代中，粒子通过跟随两个极值来更新自己的位置和速度：一个是粒子本身获得的最优解

一个是整个粒子群中全部粒子历代搜索获得的最优解标准粒子群算法中，速度更新和位置更新如式(11)所示：

$\left\{\begin{array}{c}{v}_{\mathrm{id}}^{k+1}={\mathrm{\&omega;v}}_i^k+c_1{r}_1(p_{\mathrm{id}}^k-x_{\mathrm{id}}^k)+c_2{r}_2(g_{\mathrm{id}}^k-x_{\mathrm{id}}^k)\\ x_{\mathrm{id}}^{k+1}=x_{\mathrm{id}}^k+v_{\mathrm{id}}^{k+1}\end{array}---\left(11\right)\right.$

式中，i∈N(1,n_s),n_s为粒子群中粒子的个数，N代表整数；d∈N(1,n)，n为解向量的维数，k为迭代次数，c₁和c₂为学习因子；ω为惯量权重；r₁和r₂为[0,1]之间的随机数，判断粒子群中粒子优劣的标准是蒙特卡洛仿真获得的代价函数W，代价函数越小，在下一代该粒子存在的概率越大；

设定粒子最大的迭代次数，判断迭代是否到达最大迭代次数，如果到达，输出最好粒子代表的网络随机延迟τ_i下的控制系统参数

否则，继续进行迭代；

步骤三：确定飞行控制系统的增益调度策略；

采用线性插值进行网络飞行控制系统的增益调度，假设网络随机延迟在τ₁和τ_i之间，网络随机延迟在τ₁～N(μ₁,σ₁)，τ₂～N(μ₂,σ₂)......τ_i～N(μ_i,σ_i)时分别利用随机鲁棒分析与设计方法设计的控制系统参数

控制器当前测量的网络随机延迟为τ_N，则实际采用的控制器反馈系数K_N为：

$K_N=\left\{\begin{array}{c}\frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_1}{{\mathrm{\&tau;}}_2-{\mathrm{\&tau;}}_1}\&CenterDot;(K_2^N-K_1^N)+K_1^N,{\mathrm{\&tau;}}_1\&le;{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_2\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_2}{{\mathrm{\&tau;}}_3-{\mathrm{\&tau;}}_2}\&CenterDot;(K_3^N-K_2^N)+K_2^N,{\mathrm{\&tau;}}_2<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_3\\ \&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\&CenterDot;\\ \frac{{\mathrm{\&tau;}}_N-{\mathrm{\&tau;}}_{i-1}}{{\mathrm{\&tau;}}_i-{\mathrm{\&tau;}}_{i-1}}\&CenterDot;(K_i^N-K_{i-1}^N)+K_{i-1}^N,{\mathrm{\&tau;}}_{i-1}<{\mathrm{\&tau;}}_N\&le;{\mathrm{\&tau;}}_i\end{array}\right.$

步骤四：闭环六自由度非线性蒙特卡洛仿真验证。

2.根据权利要求1所述的一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法，其特征在于：所述的动力学参数包括纵向力和力矩系数：包括升力系数C_L0、C_La、C_Lq、阻力系数C_D0、C_Da、C_Dq、俯仰力矩系数C_m0、C_ma、C_mq、

其中，C_L0为攻角为0度时的升力系数，C_La为升力关于攻角的升力系数，C_Lq为升力关于俯仰角速度的升力系数，为升力关于升降舵的升力系数，C_D0为攻角为0度时的阻力系数，C_Da为阻力关于攻角的阻力系数，C_Dq为阻力关于俯仰角速度的阻力系数，

为阻力关于升降舵的阻力系数，C_m0为攻角为0度时的俯仰力矩系数，C_ma为俯仰力矩关于攻角的力矩系数，C_mq为俯仰力矩关于俯仰角速度的力矩系数、

为俯仰力矩关于升降舵的力矩系数；

所述的物理参数包括无人飞行器的质量m，单位kg；平均几何弦长C_A，单位m；翼展b，单位m；飞行器参考面积S_w单位m²；X轴转动惯量I_X，单位kg·m²；Y轴转动惯量I_Y，单位kg·m²；Z轴转动惯量I_Z，单位kg·m²；惯量积I_XZ，单位kg·m²；舵机时间常数T_δ，单位s；舵机放大系数K_δ；发动机时间常数T_t，单位s；发动机放大系数K_t、怠速推力t_A，单位N；最大推力t_max和最小推力t_min，单位N。

3.根据权利要求1所述的一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法，其特征在于：步骤一（2）中所述的飞行器纵向非线性动力学与运动学方程具体如下：

选取无人机飞行状态向量为：

$\stackrel{-}{X}\left(t\right)={\left[\begin{array}{ccccccccc}V& \mathrm{\&alpha;}& q& \mathrm{\&theta;}& x& h& m& {\mathrm{\&delta;}}_e& {\mathrm{\&delta;}}_t\end{array}\right]}^T$

其中，V表示速度、α表示攻角、q表示俯仰角速度、θ表示俯仰角、x表示东向位置、h表示高度、m表示质量、δ_e表示升降舵偏角、δ_t表示发动机推力；

建立无人飞行器非线性六自由度运动学与动力学方程如下：

$\left\{\begin{array}{c}\stackrel{\&CenterDot;}{V}=({\mathrm{\&delta;}}_t\&CenterDot;\cos \left(\mathrm{\&alpha;}\right)-D-m\&CenterDot;g\&CenterDot;\sin (\mathrm{\&theta;}-\mathrm{\&alpha;}))/m\\ \stackrel{\&CenterDot;}{\mathrm{\&alpha;}}=(-{\mathrm{\&delta;}}_t\&CenterDot;\sin \left(\mathrm{\&alpha;}\right)-L+m\&CenterDot;V\&CenterDot;q+m\&CenterDot;g\&CenterDot;\cos (\mathrm{\&theta;}-\mathrm{\&alpha;}))/(m\&CenterDot;V)\\ \stackrel{\&CenterDot;}{q}=M/I_Y\\ \stackrel{\&CenterDot;}{\mathrm{\&theta;}}=q\\ \stackrel{\&CenterDot;}{x}=V\&CenterDot;\cos \left(\mathrm{\&mu;}\right)\\ \stackrel{\&CenterDot;}{h}=V\&CenterDot;\sin \left(\mathrm{\&mu;}\right)\\ \stackrel{\&CenterDot;}{m}=K_m\&CenterDot;{\mathrm{\&delta;}}_T\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e=(-{\mathrm{\&delta;}}_e+K_{\mathrm{\&delta;}}\&CenterDot;{\mathrm{\&delta;}}_e^{*})/T_{\mathrm{\&delta;}}\\ {\stackrel{\&CenterDot;}{\mathrm{\&delta;}}}_e=(-{\mathrm{\&delta;}}_t+K_t\&CenterDot;{\mathrm{\&delta;}}_t^{*})/T_t\\ \mathrm{\&mu;}=\mathrm{\&theta;}-\mathrm{\&alpha;}\end{array}\right.---\left(1\right)$

其中：表示速度的导数、表示攻角的导数、

俯仰角速度的导数、

表示俯仰角的导数、

表示东向位置导数、

表示高度的导数、表示质量消耗率、

表示升降舵偏角的导数、

表示发动机推力的导数、

表示升降舵偏角指令、

表示油门指令、μ为航迹倾角、K_m为发动机耗油率、δ_T为发动机推力、T_δ为舵机时间常数、T_t为发动机时间常数、K_δ为舵机放大系数、K_t为发动机放大系数，I_Y为Y轴转动惯量，式（1）简记为：

其中

为飞行状态向量的导数，

为飞行控制向量；

$\left\{\begin{array}{c}{c}_L=C_{L0}+C_{\mathrm{L\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{\mathrm{Lq}}\&CenterDot;q+C_{L{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e\\ C_D=C_{D0}+C_{\mathrm{D\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{\mathrm{Dq}}\&CenterDot;q+C_{D{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e\\ C_m=C_{m0}+C_{\mathrm{m\&alpha;}}\&CenterDot;\mathrm{\&alpha;}+C_{\mathrm{mq}}\&CenterDot;q+C_{m{\mathrm{\&delta;}}_e}\&CenterDot;{\mathrm{\&delta;}}_e\\ L=\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_L\\ D=\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_D\\ M=\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_m\&CenterDot;C_A\end{array}\right.---\left(2\right)$

其中：C_L为总升力系数、C_D为总阻力系数、C_m为总俯仰力矩系数、ρ为空气密度、L为总升力、D为总阻力、M为总俯仰力矩、S_w为飞行器参考面积、C_A为飞行器平均几何弦长、b为翼展。

4.根据权利要求1所述的一种存在变分布网络随机延迟的微小型无人飞行器纵向控制方法，其特征在于：步骤一（4）中所述的LQR设计控制律，具体为得到最终应用的控制器反馈系数，针对纵向：

选取Q₁和R₁矩阵，采用线性二次型调节器LQR设计纵向控制律K，K为2×7的矩阵；用k_ij表示K中的第i行第j列元素，选择以下控制器结构：

$\left\{\begin{array}{c}{\mathrm{\&delta;}}_e^{*}={\mathrm{\&delta;}}_{e0}+k_{12}\&CenterDot;(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0)-k_{13}\&CenterDot;q-k_{14}\&CenterDot;(\mathrm{\&theta;}-{\mathrm{\&alpha;}}_0)-k_{15}\&CenterDot;(h-h^{*})\\ {\mathrm{\&delta;}}_t^{*}=t_A-k_{21}\&CenterDot;(V-V_0)\end{array}---\left(5\right)\right.$

其中：δ_e0为配平升降舵偏角，α₀为无人飞行器的配平攻角，h^*为高度指令，V₀为配平飞行速度；对于

$\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;(C_{L0}+C_{\mathrm{L\&alpha;}}\&CenterDot;(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0))+{\mathrm{\&delta;}}_t\&CenterDot;\sin \left(\mathrm{\&alpha;}\right)=(N_Z-1)\&CenterDot;m\&CenterDot;g---\left(6\right)$

其中：N_Z是无人飞行器的法向过载，进一步化简并代入数据得：

$(\mathrm{\&alpha;}-{\mathrm{\&alpha;}}_0)=(N_Z-1)\&CenterDot;m\&CenterDot;g/(\frac{1}{2}\mathrm{\&rho;}{V}^2\&CenterDot;S_w\&CenterDot;C_{\mathrm{L\&alpha;}})---\left(7\right)$

联立式（5）、式（7）得实际飞行中的控制器结构：

$\left\{\begin{array}{c}{\mathrm{\&delta;}}_e^{*}={\mathrm{\&delta;}}_{e0}-k_1\&CenterDot;(N_z-1)-k_2\&CenterDot;q-k_3\&CenterDot;(\mathrm{\&theta;}-{\mathrm{\&alpha;}}_0)-k_4\&CenterDot;(h-h^{*})\\ {\mathrm{\&delta;}}_t^{*}=t_A-k_5\&CenterDot;(V-V_0)\end{array}---\left(8\right)\right.$

其中：k₁k₂k₃k₄k₅表示按照式（5）和式（7）换算后的最终应用的控制器反馈系数。

Images