CN101364085B - 非线性微分黄金分割自适应控制方法 - Google Patents

Abstract

非线性微分黄金分割自适应控制方法，包括下列步骤：(1)针对单输入单输出线性时变系统建立特征模型；(2)针对所述的特征模型，构造非线性微分黄金分割自适应控制律；(3)对该控制律作用于特征模型组成的闭环系统进行稳定性分析，确定闭环系统的稳定性条件。本发明克服现有技术的不足，提供一种针对线性时变系统进行特征建模，以及采用非线性微分黄金分割的自适应控制方法，这种自适应控制方法能够实现对快变信号和具有突变斜率信号的跟踪。

Description

非线性微分黄金分割自适应控制方法

技术领域

本发明涉及一种自适应控制方法，特别是基于特征模型的非线性微分黄金分割自适应控制方法。

背景技术

为了对高阶复杂对象设计低阶控制器，文献“特征建模与挠性结构的控制”(中国科学E辑，2001，31(2)：P137-149，吴宏鑫等)针对线性定常系统提出了特征建模的方法，即结合对象的动力学特征和控制性能要求进行建模。该文献指出，对于高阶定常系统的低阶特征模型，其参数是慢时变的。在采样周期满足一定条件下，对于参数未知、线性定常或慢变的二阶对象，文献“黄金分割在自适应鲁棒控制器设计中的应用”(自动化学报，1992，18(2)：P177-185，解永春，吴宏鑫)已证明线性黄金分割自适应控制器可以保证控制系统在起动过程中的稳定性和良好的过渡特性，但对动态过程中的不同性能要求，线性黄金分割自适应控制器无法满足，为此，文献“非线性黄金分割自适应控制”(宇航学报，2002，23(6)：P1-8.吴宏鑫，王颖，解永春)针对空间交会对接过程中的常值跟踪问题，基于特征模型设计了一种比例系数为非线性的黄金分割自适应控制器。当跟踪信号不是常值，该文中的非线性比例黄金分割控制器不能满足要求。中国专利00132495.0，“一种基于对象特征模型描述的黄金分割智能控制方法”中介绍了一种基于特征模型的控制方法，该方法中的智能控制器中也包含线性黄金分割自适应控制律，并且这种方法复杂，不适用于跟踪快变信号和具有突变斜率的信号。宇航学报，1998，19(1)：P8-12，“载人飞船全系数自适应再入升力控制”一文中主要介绍了一种基于特征模型的全系数自适应控制方法，该方法只是上面所述专利00132495.0中的一种具体应用实例。需要指出的是，上述所有文献的特征建模都是针对线性定常系统的，对于线性实变系统来讲，其特征建模要比定常系统复杂得多。

发明内容

本发明的技术解决问题是：克服现有技术的不足，提供一种针对线性时变系统进行特征建模、采用非线性微分黄金分割的自适应控制方法，该方法能够实现对快变信号和具有突变斜率信号的跟踪。

本发明的技术解决方案是：非线性微分黄金分割自适应控制方法，包括下列步骤：

(1)针对单输入单输出线性时变系统建立特征模型；

(2)针对所述的特征模型，构造非线性微分黄金分割自适应控制律；

(3)对该控制律作用于特征模型组成的闭环系统进行稳定性分析，确定闭环系统的稳定性条件。

所述步骤(1)中单输入单输出线性时变系统为：

$y^{\left(n\right)}=a_0\left(t\right)y+a_1\left(t\right)\stackrel{.}{y}+...+a_{n-1}\left(t\right)y^{(n-1)}+b\left(t\right)u\left(t\right)$

其对应的特征模型为：

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)+β₁(k)u(k-1)

其对应的特征模型在工程上为：

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)

其中，

t为连续时间变量；

k为连续时间变量t的第k步离散值；

y为单输入单输出线性时变系统的输出；

为y的一阶导数；

y⁽ⁿ⁾为y的第n阶导数；

y(k)为y的第k步离散值；

α₁(k)＝2+(a₀(k)-a₀(k-1))ΔT+a₀(k)ΔT-ΔT+a₀(k)ΔT²；

α₂(k)＝-1-a₀(k)ΔT+ΔT；

β(k)＝(b(k)-b(k-1))ΔT+b(k)ΔT+b(k)ΔT²；

β₁(k)＝-b(k)ΔT；

a₀(t)，b(t)为单输入单输出线性时变系统方程系数；

a₀(k)，b(k)分别为a₀(t)，b(t)的第k步离散值；

u(t)为单输入单输出线性时变系统输入变量；

u(k)为u(t)的第k步离散值；

ΔT为采样周期。

所述步骤(2)中的非线性微分黄金分割自适应控制律为：

$u\left(k\right)=-[(L_1{\widehat{\mathrm{\α}}}_1\left(k\right)+L_2{\widehat{\mathrm{\α}}}_2\left(k\right))y\left(k\right)-L_2{\widehat{\mathrm{\α}}}_2\left(k\right)({\mathrm{\η}}_1{\left|y\left(k\right)\right|}^{\mathrm{\μ}}+{\mathrm{\η}}_2)(y\left(k\right)-y(k-1))]/\widehat{\mathrm{\β}}\left(k\right)$

所述步骤(3)中的闭环系统形式如下：

y(k+1)+f₁(k)y(k)+f₂(k)y(k-1)＝0

其中，L₁＝0.382，L₂＝0.618为黄金分割系数；

y为单输入单输出线性时变系统输出变量；

y(k)为y的第k步离散值；

η₁，η₂为非负常数，μ为常数，其取值满足步骤(3)中的稳定性条件；

和

为特征模型对应系数α₁(k)、α₂(k)、β(k)的估计值。

$f_1\left(k\right)=-{\mathrm{\α}}_1\left(k\right)+[L_1{\widehat{\mathrm{\α}}}_1\left(k\right)+L_2{\widehat{\mathrm{\α}}}_2\left(k\right)-L_2{\widehat{\mathrm{\α}}}_2\left(k\right)({\mathrm{\η}}_1{\left|y\left(k\right)\right|}^{\mathrm{\μ}}+{\mathrm{\η}}_2)]\mathrm{\β}\left(k\right)/\widehat{\mathrm{\β}}\left(k\right);$

$f_2\left(k\right)=-{\mathrm{\α}}_2\left(k\right)\left[L_2{\widehat{\mathrm{\α}}}_2\left(k\right)({\mathrm{\η}}_1{\left|y\left(k\right)\right|}^{\mathrm{\μ}}+{\mathrm{\η}}_2)\right]\mathrm{\β}\left(k\right)/\widehat{\mathrm{\β}}\left(k\right);$

L₁＝0.382，L₂＝0.618为黄金分割系数。

所述步骤(3)中闭环系统的稳定性条件为：

(a)f₁(k)的变化率满足：

$\frac{N_1\left(k\right)}{2}-\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}<-{\mathrm{\&epsiv;}}_1\mathrm{\&Delta;}{f}_1\left(k\right)<\frac{N_1\left(k\right)}{2}+\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}$

(b)f₁ ²(k)和f₂ ²(k)的变化率满足：

$\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)<-f_2^2(k+1)+(1+2{\mathrm{\&epsiv;}}_1-{\mathrm{\&delta;}}_2)f_1^2(k+1)-f_1^4(k+1)+{\mathrm{\&delta;}}_2-{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{22}$

$\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)<(1-{\mathrm{\&delta;}}_2)f_2^2(k+1)-f_1^2(k+1)f_2^2(k+1)+{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{11}$

其中，

Δf₁(k)＝f₁(k+1)-f₁(k)，Δf₂(k)＝f₂(k+1)-f₂(k)，

N₁(k)＝-2p₁₂(k+1)-2[p₁₂(k+1)f₂(k+1)-p₂₂(k+1)f₁(k+1)f₂(k+1)]

$\left\{\begin{array}{c}{p}_{11}\left(k\right)=f_2^2\left(k\right)+{\mathrm{\&delta;}}_1\\ p_{12}\left(k\right)={\mathrm{\&epsiv;}}_1{f}_1\left(k\right)\\ p_{22}\left(k\right)=f_1^2\left(k\right)+{\mathrm{\&delta;}}_2\end{array}\right.$

$\left\{\begin{array}{c}{q}_{11}\left(k\right)=p_{11}\left(k\right)-p_{22}(k+1)f_2^2(k+1)\\ q_{12}\left(k\right)=p_{12}\left(k\right)+p_{12}(k+1)f_2(k+1)-p_{22}(k+1)f_1(k+1)f_2(k+1)\\ q_{22}\left(k\right)=p_{22}\left(k\right)-p_{11}(k+1)+{2p}_{12}(k+1)f_1(k+1)-p_{22}(k+1)f_1^2(k+1)\end{array}\right.$

0＜δ₁＜δ₂，0＜δ₁₁，0＜δ₂₂，

M₁为|f₁₍k)|的上界。

本发明与现有技术相比有益效果为：

(1)本发明针对线性实变系统进行特征建模、采用非线性微分黄金分割自适应控制方法，克服了反馈线性化方法要求精确获得对象模型的缺点，并且与现有的线性黄金分割和非线性比例黄金分割控制方法相比，本发明能够实现对具有突变斜率信号的跟踪，见图2～6。

(2)本发明基于特征模型的非线性微分黄金分割自适应控制方法比线性黄金分割自适应控制方法和的非线性比例黄金分割自适应控制方法具有更快的跟踪速度和更小的超调量，见图4～6。

附图说明

图1为本发明方法应用过程流程图；

图2为本发明实施例1中对参考阻力加速度的跟踪曲线；

图3为与图2相对应的跟踪误差曲线；

图4为本发明对方波的跟踪曲线；

图5为线性黄金分割自适应控制器对方波的跟踪曲线；

图6为非线性比例黄金分割自适应控制器对方波的跟踪曲线。

具体实施方式

本发明非线性微分黄金分割自适应控制方法，具体步骤如下：

(1)针对单输入单输出线性时变系统建立特征模型；

(2)针对所述的特征模型，构造非线性微分黄金分割自适应控制律；

(3)对该控制律作用于特征模型组成的闭环系统进行稳定性分析，确定闭环系统的稳定性条件。

下面具体介绍上述步骤中各步的实现及推导过程：

1、单输入单输出线性时变系统的特征模型

步骤(1)中的特征建模过程如下：

对单输入单输出线性时变系统(方程1)进行特征建模：

$y^{\left(n\right)}=a_0\left(t\right)y+a_1\left(t\right)\stackrel{.}{y}+...+a_{n-1}\left(t\right)y^{(n-1)}+b\left(t\right)u\left(t\right)---\left(1\right)$

其中，u(t)为输入变量，y为输出变量。假设方程(1)中的系数a_i(t)，b(t)均有界，y及其各阶导数有界。由方程(1)可得

$\stackrel{.}{y}=a_0\left(t\right)y+b\left(t\right)u\left(t\right)+F(\stackrel{.}{y},...,y^{\left(n\right)})---\left(2\right)$

式中， $F(\stackrel{.}{y},...,y^{\left(n\right)})=a_1\left(t\right)\stackrel{.}{y}+...+a_{n-1}\left(t\right)y^{(n-1)}-y^{\left(n\right)}+\stackrel{.}{y}.$

方程(2)两边对时间求导得

$\stackrel{..}{y}={\stackrel{.}{a}}_0\left(t\right)y+a_0\left(t\right)\stackrel{.}{y}+\stackrel{.}{b}\left(t\right)u\left(t\right)+b\left(t\right)\stackrel{.}{u}\left(t\right)+\frac{\mathrm{dF}(\stackrel{.}{y},...,y^{\left(n\right)})}{\mathrm{dt}}---\left(3\right)$

对方程(2)两边进行差分得

$\frac{y\left(k\right)-y(k-1)}{\mathrm{\&Delta;T}}=a_0\left(k\right)y\left(k\right)+b\left(k\right)u\left(k\right)+F\left(k\right)---\left(4\right)$

对方程(3)两边进行差分得

$\frac{y(k+1)-2y\left(k\right)+y(k-1)}{{\mathrm{\&Delta;T}}^2}=\frac{a_0\left(k\right)-a_0(k-1)}{\mathrm{\&Delta;T}}y\left(k\right)+a_0\left(k\right)\frac{y\left(k\right)-y(k-1)}{\mathrm{\&Delta;T}}$

$+\frac{b\left(k\right)-b(k-1)}{\mathrm{\&Delta;T}}u\left(k\right)+b\left(k\right)\frac{u\left(k\right)-u(k-1)}{\mathrm{\&Delta;T}}+\frac{F\left(k\right)-F(k-1)}{\mathrm{\&Delta;T}}---\left(5\right)$

将方程(4)和(5)两边相加整理得

y(k+1)＝[2+(a₀(k)-a₀(k-1))ΔT+a₀(k)ΔT-ΔT+a₀(k)ΔT²]y(k)

+[-1-a₀(k)ΔT+ΔT]y(k-1)+[(b(k)-b(k-1))ΔT+b(k)ΔT+b(k)ΔT²]u(k)

+[-b(k)ΔT]u(k-1)+F(k)ΔT²+[F(k)-F(k-1)]ΔT    (6)

令

α₁(k)＝2+(a₀(k)-a₀(k-1))ΔT+a₀(k)ΔT-ΔT+a₀(k)ΔT²；

α₂(k)＝-1-a₀(k)ΔT+ΔT；

β(k)＝(b(k)-b(k-1))ΔT+b(k)ΔT+b(k)ΔT²；

β₁(k)＝-b(k)ΔT；

W(k)＝F(k)ΔT²+[F(k)-F(k-1)]ΔT

若要实现恒值控制，在稳态时，由于F(k)是y的各阶导数的线性组合，则F(k)＝F(k-1)＝0，因此W(k)＝0。在动态过程中，由于方程(1)中的各系数和y的各阶导数均有界，则F(k)和F(k-1)有界，于是，可设|F(k)|＜K(K为正常数)，故|W(k)|≤2KΔT+KΔT²，当ΔT→0时，W(k)→0，因此W(k)可看作动态过程中的建模误差。由此可得到方程(1)的特征模型为

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)+β₁(k)u(k-1)    (7)

工程上一般可取特征模型为

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)    (8)

由W(k)的表达式可知，采样周期越小，则动态过程中的建模误差越小。给定ε＞0，取0＜ΔT＜δ，δ＝min{1，ε/(3K)}，则

|W(k)|≤2KΔT+KΔT²＝(2K+KΔT)ΔT＜(2K+K)ε/(3K)＝ε

从以上推导过程可知，当采样周期满足0＜ΔT＜δ时，特征模型的建模误差小于ε，因此，特征模型与原系统是等价的，即在相同的输入作用下，特征模型的输出与实际对象的输出在稳态情况下相等，在动态情况下误差保持在允许范围内，故特征模型反映了系统的特性。

2、控制律设计

步骤(2)中的非线性微分黄金分割自适应控制律为：

$u\left(k\right)=-[(L_1{\widehat{\mathrm{\&alpha;}}}_1\left(k\right)+L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right))y\left(k\right)-L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)(y\left(k\right)-y(k-1))]/\widehat{\mathrm{\&beta;}}\left(k\right)---\left(9\right)$

其中，L₁＝0.382，L₂＝0.618为黄金分割系数，η₁，η₂为非负常数，μ为常数，其取值满足步骤(3)中的稳定性条件，

和

为步骤(1)中特征模型对应系数α₁(k)、α₂(k)、β(k)的估计值，可采用任何一种估计算法确定，如梯度算法，或递推最小二乘法。(参见《最优状态估计与系统辨识》，王志贤编著，西北工业大学出版社，2004)

3、闭环系统稳定的充分条件

步骤(3)中的闭环系统形式如下：

y(k+1)+f₁(k)y(k)+f₂(k)y(k-1)＝0    (10)

其中，

$f_1\left(k\right)=-{\mathrm{\&alpha;}}_1\left(k\right)[L_1{\widehat{\mathrm{\&alpha;}}}_1\left(k\right)+L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)-L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)]\mathrm{\&beta;}\left(k\right)/\widehat{\mathrm{\&beta;}}\left(k\right)$

$f_2\left(k\right)=-{\mathrm{\&alpha;}}_2\left(k\right)+\left[L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)\right]\mathrm{\&beta;}\left(k\right)/\widehat{\mathrm{\&beta;}}\left(k\right)$

将闭环系统(10)写为状态空间的形式：

Y(k+1)＝A(k+1)Y(k)    (11)

其中，

$A(k+1)=\left[\begin{array}{cc}0& 1\\ -f_2(k+1,y(k+1))& -f_1(k+1,y(k+1))\end{array}\right]$

Y(k)＝[y(k)y(k+1)]^T

首先选取李亚普诺夫函数为

V(k)＝Y^T(k)P(k)Y(k)

其中， $P\left(k\right)=\left[\begin{array}{cc}{p}_{11}\left(k\right)& p_{12}\left(k\right)\\ p_{12}\left(k\right)& p_{22}\left(k\right)\end{array}\right],$

$\left\{\begin{array}{c}{p}_{11}\left(k\right)=f_2^2\left(k\right)+{\mathrm{\&delta;}}_1\\ p_{12}\left(k\right)={\mathrm{\&epsiv;}}_1{f}_1\left(k\right)\\ p_{22}\left(k\right)=f_1^2\left(k\right)+{\mathrm{\&delta;}}_2\end{array}\right.,$

0＜δ₁＜δ₂，

$0<{\mathrm{\&epsiv;}}_1<\frac{{\mathrm{\&delta;}}_1}{\sqrt{3}{M}_1},$

M₁为|f₁(k)|的上界。

则

ΔV(k)＝-Y^T(k)[P(k)-A^T(k+1)P(k+1)A(k+1)]Y(k)    (12)

令

Q(k)＝P(k)-A^T(k+1)P(k+1)A(k+1)    (13)

则

$Q\left(k\right)=\left[\begin{array}{cc}{q}_{11}\left(k\right)& q_{12}\left(k\right)\\ q_{12}\left(k\right)& q_{22}\left(k\right)\end{array}\right]$

$\left\{\begin{array}{c}{q}_{11}\left(k\right)=p_{11}\left(k\right)-p_{22}(k+1)f_2^2(k+1)\\ q_{12}\left(k\right)=p_{12}\left(k\right)+p_{12}(k+1)f_2(k+1)-p_{22}(k+1)f_1(k+1)f_2(k+1)\\ q_{22}\left(k\right)=p_{22}\left(k\right)-p_{11}(k+1)+{2p}_{12}(k+1)f_1(k+1)-p_{22}(k+1)f_1^2(k+1)\end{array}\right.$

假设系数α₁(k)，α₂(k)，β₀(k)及估计值

属于有界闭凸集，假设

和|y(k)|^μ＜M，M＞0为常数，那么f₁(k)和f₂(k)有界。

由于f₁(k)和f₂(k)均有界，则P(k)的各元素均有界，Q(k)的各元素均有界。

其次，证明P(k)是一致有界且正定的矩阵。

P(k)的一阶顺序主子式：

$M_{P_1}=p_{11}\left(k\right)=f_2^2\left(k\right)+{\mathrm{\&delta;}}_1>{\mathrm{\&delta;}}_1---\left(14\right)$

P(k)的二阶顺序主子式：

$M_{P_2}=p_{11}\left(k\right)p_{22}\left(k\right)-p_{12}^2\left(k\right)=(f_2^2\left(k\right)+{\mathrm{\&delta;}}_1)(f_1^2\left(k\right)+{\mathrm{\&delta;}}_2)-{\mathrm{\&epsiv;}}_1^2{f}_1^2\left(k\right)>{\mathrm{\&delta;}}_1^2-\frac{1}{3}{\mathrm{\&delta;}}_1^2---\left(15\right)$

根据式(14)和式(15)，矩阵P(k)的一阶顺序主子式和二阶顺序主子式均大于正的小常数。由于P(k)的各元素均有界，则

和

有界。那么可以判定P(k)是一致有界且正定的矩阵。

第三步：证明Q(k)是一致有界且正定的矩阵。

定义：

Δf₁(k)＝f₁(k+1)-f₁(k)

Δf₂(k)＝f₂(k+1)-f₂(k)

$\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)=f_1^2(k+1)-f_1^2\left(k\right)$

$\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)=f_2^2(k+1)-f_2^2\left(k\right)$

0＜δ₁₁，0＜δ₂₂， $0<\stackrel{\&OverBar;}{\mathrm{\&delta;}}<{\mathrm{\&delta;}}_{11}{\mathrm{\&delta;}}_{22},$

则

$q_{11}\left(k\right)=-\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)+(1-{\mathrm{\&delta;}}_2)f_2^2(k+1)-f_1^2(k+1)f_2^2(k+1)+{\mathrm{\&delta;}}_1---\left(16\right)$

$q_{22}\left(k\right)=-\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)-f_2^2(k+1)+(1+2{\mathrm{\&epsiv;}}_1-{\mathrm{\&delta;}}_2)f_1^2(k+1)-f_1^4(k+1)+{\mathrm{\&delta;}}_2-{\mathrm{\&delta;}}_1---\left(17\right)$

当Δ(f₂ ²(k))满足(18)，Δ(f₁ ²(k))满足(19)时，则q₁₁(k)＞δ₁₁，q₂₂(k)＞δ₂₂。

$\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)<(1-{\mathrm{\&delta;}}_2)f_2^2(k+1)-f_1^2(k+1)f_2^2(k+1)+{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{11}---\left(18\right)$

$\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)<-f_2^2(k+1)+(1+2{\mathrm{\&epsiv;}}_1-{\mathrm{\&delta;}}_2)f_1^2(k+1)-f_1^4(k+1)+{\mathrm{\&delta;}}_2-{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{22}---\left(19\right)$

Q(k)的一阶顺序主子式：

$M_{Q_1}=q_{11}\left(k\right)>{\mathrm{\&delta;}}_{11}---\left(20\right)$

Q(k)的二阶顺序主子式：

$M_{Q_2}=q_{11}\left(k\right)q_{22}\left(k\right)-q_{12}^2\left(k\right)=-{[p_{12}\left(k\right)-p_{12}(k+1)]}^2+N_1\left(k\right)[p_{12}\left(k\right)-p_{12}(k+1)]+N_0\left(k\right)---\left(21\right)$

其中

N₁(k)＝-2p₁₂(k+1)-2[p₁₂(k+1)f₂(k+1)-p₂₂(k+1)f₁(k+1)f₂(k+1)]    (22)

$N_0\left(k\right)=p_{12}^2(k+1)+N_1\left(k\right)p_{12}(k+1)$

$-{[p_{12}(k+1)f_2(k+1)-p_{22}(k+1)f_1(k+1)f_2(k+1)]}^2+q_{11}\left(k\right)q_{22}\left(k\right)---\left(23\right)$

二次多项式-[p₁₂(k)-p₁₂(k+1)]²+N₁(k)[p₁₂(k)-p₁₂(k+1)]+N₀(k)的判别式为：

$\mathrm{\&Delta;}\left(k\right)=N_1^2\left(k\right)+4N_0\left(k\right)=4q_{11}\left(k\right)q_{22}\left(k\right)>{4\mathrm{\&delta;}}_{11}{\mathrm{\&delta;}}_{22}---\left(24\right)$

取 $\mathrm{\&epsiv;}=4\stackrel{\&OverBar;}{\mathrm{\&delta;}},$ $0<\stackrel{\&OverBar;}{\mathrm{\&delta;}}<{\mathrm{\&delta;}}_{11}{\mathrm{\&delta;}}_{22},$ 则Δ(k)＞ε。

$\frac{N_1\left(k\right)-\sqrt{\mathrm{\&Delta;}\left(k\right)-\mathrm{\&epsiv;}}}{2}<p_{12}\left(k\right)-p_{12}(k+1)<\frac{N_1\left(k\right)+\sqrt{\mathrm{\&Delta;}\left(k\right)-\mathrm{\&epsiv;}}}{2}---\left(25\right)$

$\frac{N_1\left(k\right)}{2}-\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}<-{\mathrm{\&epsiv;}}_1\mathrm{\&Delta;}{f}_1\left(k\right)<\frac{N_1\left(k\right)}{2}+\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}---\left(26\right)$

当式(25)成立，即式(26)成立时，由文献“四阶时变离散系统的一致渐近稳定性”(控制理论与应用，2006，23(6)：P845-852，孙多青，吴宏鑫)中的定理知

$M_{Q_2}>\frac{\mathrm{\&epsiv;}}{4M}=\stackrel{\&OverBar;}{\mathrm{\&delta;}}---\left(27\right)$

根据式(21)和式(28)，矩阵Q(k)的一阶顺序主子式

和二阶顺序主子式

均大于正的小常数。由于Q(k)的各元素均有界，则

和

有界，由此可以判定Q(k)是一致有界且正定的矩阵。

综上所述，可以得出闭环系统(10)在以原点为平衡位置处一致渐近稳定的充分条件为：

(a)f₁(k)的变化率满足：

$\frac{N_1\left(k\right)}{2}-\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}<-{\mathrm{\&epsiv;}}_1\mathrm{\&Delta;}{f}_1\left(k\right)<\frac{N_1\left(k\right)}{2}+\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}$

(b)f₁ ²(k)和f₂ ²(k)的变化率满足：

$\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)<-f_2^2(k+1)+(1+{2\mathrm{\&epsiv;}}_1-{\mathrm{\&delta;}}_2)f_1^2(k+1)-f_1^4(k+1)+{\mathrm{\&delta;}}_2-{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{22}$

$\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)<(1-{\mathrm{\&delta;}}_2)f_2^2(k+1)-f_1^2(k+1)f_2^2(k+1)+{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{11}$

其中，

0＜δ₁＜δ₂，0＜δ₁₁，0＜δ₂₂，

M₁为|F₁(k)|的上界，Δf₁(k)＝f₁(k+1)-f₁(k)，Δf₂(k)＝f₂(k+1)-f₂(k)， $\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)=f_1^2(k+1)-f_1^2\left(k\right),$ $\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)=f_2^2(k+1)-f_2^2\left(k\right)$

N₁(k)＝-2p₁₂(k+1)-2[p₁₂(k+1)f₂(k+1)-p₂₂(k+1)f₁(k+1)f₂(k+1)]

$\left\{\begin{array}{c}{p}_{11}\left(k\right)=f_2^2\left(k\right)+{\mathrm{\&delta;}}_1\\ p_{12}\left(k\right)={\mathrm{\&epsiv;}}_1{f}_1\left(k\right)\\ p_{22}\left(k\right)=f_1^2\left(k\right)+{\mathrm{\&delta;}}_2\end{array}\right.$

$\left\{\begin{array}{c}{q}_{11}\left(k\right)=p_{11}\left(k\right)-p_{22}(k+1)f_2^2(k+1)\\ q_{12}\left(k\right)=p_{12}\left(k\right)+p_{12}(k+1)f_2(k+1)-p_{22}(k+1)f_1(k+1)f_2(k+1)\\ q_{22}\left(k\right)=p_{22}\left(k\right)-p_{11}(k+1)+{2p}_{12}(k+1)f_1(k+1)-p_{22}(k+1)f_1^2(k+1)\end{array}\right.$

图1为本发明方法应用的流程图，其具体应用过程如下：

第一步，将整个应用系统的参考输入y_r(k)与系统输出y_out(k)作差，作为非线性微分黄金分割自适应控制器的输入y(k)，根据第k步的参数估计值

和

按照公式(9)得到非线性微分黄金分割自适应控制器的输出u(k+1)；

第二步，将该输出值u(k+1)作为被控对象的输入得到系统下一时刻的输出y_out(k+1)；

第三步，通过u(k+1)与y_out(k+1)进行参数估计，得到特征模型在k+1时刻的参数估计值

和

令k＝k+1，循环执行第一步到第三步，直至控制结束。

实施例1

下面针对再入飞行器结合上述说明详细介绍本发明的方法。

首先根据再入飞行器的阻力加速度方程确定阻力加速度的摄动方程。阻力加速度动力学方程：

${\stackrel{..}{a}}_D=\frac{{\stackrel{.}{a}}_D^2}{a_D}-2\frac{a_D{\stackrel{.}{a}}_D}{v}-\frac{{2a}_D^3}{v^2}+\frac{{2a}_D{g}^2}{v^2}-\frac{{2a}_{D}g}{r}-\frac{a_D{v}^2}{h_{s}r}+\frac{a_{D}g}{h_s}+\frac{a_D{\stackrel{..}{C}}_D}{C_D}-\frac{a_D{\stackrel{.}{C}}_D^2}{C_D^2}+(-\frac{{2\mathrm{ga}}_D^2}{v^2}-\frac{a_D^2}{h_s})u$

其中a_D＝ρSC_Dv²/(2m)为阻力加速度，v为速度，r为飞行器质心到地心的距离，g为重力加速度，h_s为一常数，C_D为阻力系数，ρ为大气密度，m为飞行器的质量，S为飞行器的参考面积，控制输入u为铅垂面内的升阻比。

定义偏离参考轨道的小偏差量为

下标0表示各量对应于参考轨道的取值。对阻力加速度动力学方程两边取变分获得摄动方程，公式如下：

$\mathrm{\&delta;}{\stackrel{..}{a}}_D=a_1\left(t\right)\mathrm{\&delta;}{\stackrel{.}{a}}_D+a_0\left(t\right)\mathrm{\&delta;}{a}_D+b\left(t\right)\mathrm{\&delta;u}$

其中，

$a_1\left(t\right)={\stackrel{.}{2a}}_{D_0}/a_{D_0}-{2a}_{D_0}/v_0$

$a_0\left(t\right)=-\frac{{\stackrel{.}{a}}_{D_0}^2}{a_{D_0}^2}-\frac{{2\stackrel{.}{a}}_{D_0}}{v_0}-\frac{{6a}_{D_0}^2}{v_0^2}+\frac{{2g}^2}{v_0^2}-\frac{2g}{r_0}-\frac{v_0^2}{h_s{r}_0}+\frac{g}{h_s}+\frac{{\stackrel{..}{C}}_{D_0}}{C_{D_0}}-\frac{{\stackrel{.}{C}}_{D_0}^2}{C_{D_0}^2}-\frac{{4\mathrm{ga}}_{D_0}}{v_0^2}{u}_0-\frac{{2a}_{D_0}}{h_s}{u}_0$

$b\left(t\right)=(-2g/v_0^2-1/h_s)a_{D_0}^2$

(1)针对所述的摄动方程建立特征模型；

以u(k)＝δu(k)作为输入，y(k)＝δa_D(k)作为输出，按照步骤(1)建立特征模型为：

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)

其中，

α₁(k)＝2+a₁(k)ΔT+a₀(k)ΔT²

α₂(k)＝-1-a₁(k)ΔT

β(k)＝b(k)ΔT²

ΔT为采样周期

(2)针对步骤(1)的特征模型，构造非线性微分黄金分割自适应控制律，具体如下：

$u\left(k\right)=-[(L_1{\widehat{\mathrm{\&alpha;}}}_1\left(k\right)+L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right))y\left(k\right)-L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)(y\left(k\right)-y(k-1))]/\widehat{\mathrm{\&beta;}}\left(k\right)$

其中，L₁＝0.382，L₂＝0.618为黄金分割系数；η₁＝1，η₂＝50，μ＝1。

和

为步骤(1)中特征模型对应系数的估计值，估计方法采用梯度投影算法。

(3)对该控制律作用于特征模型组成的闭环系统进行稳定性分析，得出的结论是在上述控制律的作用下，该闭环系统是稳定的。

图2、3为本发明应用到再入飞行器跟踪分段线性的参考阻力加速度的结果，从图中可以看出，本发明能够很好的跟踪具有突变斜率的信号，本例中对参考阻力加速度的最大跟踪误差为|Δa_D|＝0.52m/s²。

实施例2

范德堡方程： $m\stackrel{..}{y}\left(t\right)+2c(y^2\left(t\right)-1)\stackrel{.}{y}\left(t\right)+\mathrm{ky}\left(t\right)=\mathrm{bu}\left(t\right)$

其中，m＝2，c＝3，k＝4，b＝1，要求y(t)跟踪频率为100HZ的方波信号y_r(t)，采用本发明设计控制器实现对参考信号y_r(t)的跟踪。具体实现过程如下：

定义偏离参考曲线的小偏差量为

δu＝u(t)-u_r(t)。

下标r表示各量对应于参考曲线的取值。对范德堡方程两边取变分获得摄动方程，公式如下：

$\mathrm{\&delta;}\stackrel{..}{y}=a_1\left(t\right)\mathrm{\&delta;}\stackrel{.}{y}+a_0\left(t\right)\mathrm{\&delta;y}+b_0\left(t\right)\mathrm{\&delta;u}$

其中，

$a_1\left(t\right)=-\frac{2c[y_r^2\left(t\right)-1]}{m}$

$a_0\left(t\right)=-\frac{[{4\mathrm{cy}}_r\left(t\right){\stackrel{.}{y}}_r\left(t\right)+k]}{m}$

$b_0\left(t\right)=\frac{b}{m}$

(1)对上述摄动方程建立特征模型

以u(k)＝δu(k)作为输入，y(k)＝δy(k)作为输出，按照步骤(1)建立特征模型为：

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)

其中，

α₁(k)＝2+a₁(k)ΔT+a₀(k)ΔT²

α₂(k)＝-1-a₁(k)ΔT

β(k)＝b₀(k)ΔT²

ΔT为采样周期

(2)构造非线性微分黄金分割自适应控制律，过程同前对图1的介绍；其自适应控制律具体如下：

$u\left(k\right)=-[(L_1{\widehat{\mathrm{\&alpha;}}}_1\left(k\right)+L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right))y\left(k\right)-L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)(y\left(k\right)-y(k-1))]/\widehat{\mathrm{\&beta;}}\left(k\right)$

其中，L₁＝0.382，L₂＝0.618为黄金分割系数；η₁＝10，η₂＝0.1，μ＝-0.1。

和

为步骤(1)中特征模型对应系数的估计值，估计方法采用梯度投影算法。

(3)分析闭环系统的稳定性，通过验证证明该闭环系统是稳定的；

图4为本发明设计的控制器对方波的跟踪情况，图6为现有的线性黄金分割控制器的跟踪情况，图7为现有的非线性比例黄金分割控制器的跟踪情况，通过比较图4、5和6可知，本发明比线性黄金分割自适应控制方法和的非线性比例黄金分割自适应控制方法具有更快的跟踪速度和更小的超调量。

本发明未详细说明部分属本领域技术人员公知常识。

Claims (6)

1.非线性微分黄金分割自适应控制方法，其特征在于包括下列步骤：

(1)根据再入飞行器的阻力加速度方程确定阻力加速度的摄动方程；针对所述的摄动方程建立特征模型；

(2)针对步骤(1)的特征模型，构造非线性微分黄金分割自适应控制律；

(3)对该控制律作用于特征模型组成的闭环系统进行稳定性分析，确定闭环系统的稳定性条件；得到该特征模型的非线性微分黄金分割自适应控制器；

(4)将整个应用系统的参考阻力加速度y_r(k)与应用系统输出y_out(k)作差，作为非线性微分黄金分割自适应控制器的输入y(k)，根据第k步的参数估计值

和

得到非线性微分黄金分割自适应控制器的输出u(k+1)；所述的应用系统包括非线性微分黄金分割自适应控制器和再入飞行器；

(5)将该输出值u(k+1)作为再入飞行器的输入得到应用系统下一时刻的输出y_out(k+1)；

(6)通过u(k+1)与y_out(k+1)进行参数估计，得到特征模型在k+1时刻的参数估计值和

令k＝k+1，循环执行步骤(4)到(6)，实现对参考阻力加速度y_r(k)的跟踪，进而实现对再入飞行器的控制。

2.根据权利要求1所述的非线性微分黄金分割自适应控制方法，其特征在于：所述步骤(1)中特征模型为：

y(k+1)＝α₁(k)y(k)+α₂(k)y(k-1)+β(k)u(k)

其中，

α₁(k)＝2+a₁(k)ΔT+a₀(k)ΔT²；

α₂(k)＝-1-a₁(k)ΔT；

β(k)＝b(k)ΔT²；

ΔT为采样周期；

a₀(k)，b(k)，a₁(k)分别为a₀(t)，b(t)，a₁(t)的第k步离散值；

a₀(t)，b(t)，a₁(t)为步骤(1)中阻力加速度摄动方程的系数。

3.根据权利要求1所述的非线性微分黄金分割自适应控制方法，其特征在于：所述步骤(2)中的非线性微分黄金分割自适应控制律为：

$u\left(k\right)=-[(L_1{\widehat{\mathrm{\&alpha;}}}_1\left(k\right)+L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right))y\left(k\right)-L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)(y\left(k\right)-y(k-1))]/\widehat{\mathrm{\&beta;}}\left(k\right)$

其中，L₁＝0.382，L₂＝0.618为黄金分割系数；

y为再入飞行器阻力加速度摄动方程的输出变量；

y(k)为y的第k步离散值；

η₁，η₂为非负常数，μ为常数，其取值满足步骤(3)中的稳定性条件；

和

为特征模型对应系数α₁(k)、α₂(k)、β(k)的估计值。

4.根据权利要求3所述的非线性微分黄金分割自适应控制方法，其特征在于：所述

和的估计值采用估计算法确定，如梯度算法，或递推最小二乘法。

5.根据权利要求1所述的非线性微分黄金分割自适应控制方法，其特征在于：所述步骤(3)中的闭环系统形式如下：

y(k+1)+f₁(k)y(k)+f₂(k)y(k-1)＝0

其中，

$f_1\left(k\right)=-{\mathrm{\&alpha;}}_1\left(k\right)+[L_1{\widehat{\mathrm{\&alpha;}}}_1\left(k\right)+L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)-L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)]\mathrm{\&beta;}\left(k\right)/\widehat{\mathrm{\&beta;}}\left(k\right);$

$f_2\left(k\right)=-{\mathrm{\&alpha;}}_2\left(k\right)+\left[L_2{\widehat{\mathrm{\&alpha;}}}_2\left(k\right)({\mathrm{\&eta;}}_1{\left|y\left(k\right)\right|}^{\mathrm{\&mu;}}+{\mathrm{\&eta;}}_2)\right]\mathrm{\&beta;}\left(k\right)/\widehat{\mathrm{\&beta;}}\left(k\right);$

L₁＝0.382，L₂＝0.618为黄金分割系数；

y为再入飞行器阻力加速度摄动方程的输出变量；

y(k)为y的第k步离散值；

η₁，η₂为非负常数，μ为常数，其取值满足步骤(3)中的稳定性条件；

和

为特征模型对应系数α₁(k)、α₂(k)、β(k)的估计值。

6.根据权利要求5所述的非线性微分黄金分割自适应控制方法，其特征在于：所述步骤(3)中闭环系统的稳定性条件为：

(a)f₁(k)的变化率满足：

$\frac{N_1\left(k\right)}{2}-\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}<-{\mathrm{\&epsiv;}}_1\mathrm{\&Delta;}{f}_1\left(k\right)<\frac{N_1\left(k\right)}{2}+\sqrt{q_{11}\left(k\right)q_{22}\left(k\right)-\stackrel{\&OverBar;}{\mathrm{\&delta;}}}$

(b)f₁ ²(k)和f₂ ²(k)的变化率满足：

$\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)<-f_2^2(k+1)+(1+2{\mathrm{\&epsiv;}}_1-{\mathrm{\&delta;}}_2)f_1^2(k+1)-f_1^4(k+1)+{\mathrm{\&delta;}}_2-{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{22}$

$\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)<(1-{\mathrm{\&delta;}}_2)f_2^2(k+1)-f_1^2(k+1)f_2^2(k+1)+{\mathrm{\&delta;}}_1-{\mathrm{\&delta;}}_{11}$

其中，

Δf₁(k)＝f₁(k+1)-f₁(k)，Δf₂(k)＝f₂(k+1)-f₂(k)， $\mathrm{\&Delta;}\left(f_1^2\left(k\right)\right)=f_1^2(k+1)-f_1^2\left(k\right),$

$\mathrm{\&Delta;}\left(f_2^2\left(k\right)\right)=f_2^2(k+1)-f_2^2\left(k\right)$

N₁(k)＝-2p₁₂(k+1)-2[p₁₂(k+1)f₂(k+1)-p₂₂(k+1)f₁(k+1)f₂(k+1)]

$\left\{\begin{array}{c}{p}_{11}\left(k\right)=f_2^2\left(k\right)+{\mathrm{\&delta;}}_1\\ p_{12}\left(k\right)={\mathrm{\&epsiv;}}_1{f}_1\left(k\right)\\ p_{22}\left(k\right)=f_1^2\left(k\right)+{\mathrm{\&delta;}}_2\end{array}\right.$

$\left\{\begin{array}{c}{q}_{11}\left(k\right)=p_{11}\left(k\right)-p_{22}(k+1)f_2^2(k+1)\\ q_{12}\left(k\right)=p_{12}\left(k\right)+p_{12}(k+1)f_2(k+1)-p_{22}(k+1)f_1(k+1)f_2(k+1)\\ q_{22}\left(k\right)=p_{22}\left(k\right)-p_{11}(k+1)+2p_{12}(k+1)f_1(k+1)-p_{22}(k+1)f_1^2(k+1)\end{array}\right.$

0＜δ₁＜δ₂，0＜δ₁₁，0＜δ₂₂，

M₁为|f₁(k)|的上界。

Images