第5章
CHAPTER 5


RBF神经网络滑模控制








早在20世纪50年代，苏联学者Emelyanov等就提出了滑模控制(Sliding Mode Control，SMC)方法。在以后几十年中，滑模控制设计引起了国内外学者的广泛关注。 

滑模控制为含有不确定性的非线性系统鲁棒控制提供了有效的控制设计方法。滑模变结构控制的原理是根据系统所期望的动态特性来设计系统的切换超平面，通过滑动模态控制器使系统状态从超平面之外向切换超平面运动。系统一旦到达切换超平面，控制作用将保证系统状态沿切换超平面到达系统原点，这一沿切换超平面向原点的滑动过程称为滑模运动。由于系统的特性和参数只取决于设计的切换超平面而与外界干扰没有关系，所以滑模变结构控制具有很强的鲁棒性。

近年来，一些学者［12］ 将滑模控制结合神经网络用于非线性系统的控制中。滑模控制存在稳定性分析难，达到条件难以满足及抖振等问题［3］。如果系统的数学模型已知，滑模控制器可以使系统输出直接跟踪期望指令。但较大的外界扰动需要较大的切换增益，这就造成抖振，抖振是滑模控制中难以避免的问题。采用神经网络对滑模控制进行补偿，为解决这一问题提供了有效的途径。

5.1经典滑模控制器设计

针对线性系统
x·=Ax+bu,x∈Rn,u∈R(5.1)
滑模面设计为
s(x)= cTx=∑ni=1cixi=∑n－1i=1cixi+xn(5.2)
其中,x为状态向量，c=c1…cn－11T。

在滑模控制中，参数c1,c2,…,cn－1应满足多项式pn－1+cn－1pn－2+…+c2p+c1为Hurwitz，其中p为拉氏算子。

例如，n=2时，假设x·1=x2，则可设计s(x)=c1x1+x2，为了保证多项式p+c1为Hurwitz，需要多项式p+c1=0的特征根实数部分为负，即c1>0。例如，取c1=10，可得滑模面为s(x)=10x1+x2。

又例如，当n=3时，设x·1=x2，x·2=x3，则可设计s(x)=c1x1+c2x2+x3，针对n=2，为了保证多项式p2+c2p+c1为Hurwitz，需要多项式p2+c2p+c1=0特征根实数部分为负。可取p2+2λp+λ2=0，即p+λ2=0，取λ>0，即可满足多项式p2+2λp+λ2=0的特征根实数部分为负，从而可得c2=2λ，c1=λ2，例如，令λ=5，可得c1=25，c2=10，则s(x)=25x1+10x2+x3。

对于一个二阶系统，滑模控制设计可分为两步： 首先设计滑模面，以便使系统能够按照预定的“滑动模态”轨迹运动，然后设计控制器，使系统的状态沿滑模面运动。上述设计的基础是构建Lyapunov函数。

考虑如下被控对象： 
Jθ¨t=ut+dt(5.3)
其中，J为转动惯量； θt为角度； ut为控制输入； dt为外界扰动且满足dt≤D。

设计滑模函数为
st=cet+e·t(5.4)
其中，c必须满足Hurwitz条件，即c>0。

定义跟踪误差，并求导可得
et=θt－θdt，e·t=θ·t－θ·dt
其中，θdt为理想跟踪指令。

设计Lyapunov函数为
V=12s2
则
s·t=ce·t+e¨t=ce·t+θ¨t－θ¨dt=ce·t+1Ju+dt－θ¨dt(5.5)
且
ss·=sce·+1Ju+dt－θ¨d
为了保证ss·<0，设计滑模控制律为
u(t)=J［－ce·+θ¨d－ηsgns］－Dsgns(5.6)
则
ss·=sce·+1Ju+dt－θ¨d
ss·=－ηs－DJs<0
从而
V·≤0V·=0，当s=0
从控制律的表达式可知，滑模控制器具有很好的鲁棒性。然而，当干扰dt较大时，为了保证鲁棒性，必须保证足够大的干扰上界D，而较大的干扰上界D会导致切换增益过大，从而造成抖振。

另外，在控制律式(5.6)中，建模信息J必须是精确已知的，这在实际工程中是难以实现的，采用RBF神经网络逼近方法可有效地解决这一难题。

5.2基于RBF神经网络的二阶SISO系统的滑模控制
5.2.1系统描述

考虑如下二阶被控对象:
θ¨=fθ,θ·+gθ,θ·u+dt(5.7)
其中，f(·)和g(·)为非线性函数； u∈R和y∈R分别为控制输入和系统输出； dt为外界干扰，且满足dt≤D。

理想跟踪指令为θd，定义跟踪误差为 
e=θd－θ
设计滑模面为
s= e·+ce(5.8)
其中，c>0，则
s·=e¨+ce·=θ¨d－θ¨+ce·=θ¨d－f－gu－dt+ce·(5.9)
如果f和g是已知的，可设计控制律为
u=1g［－f+θ¨d+ce·+ηsgns］(5.10)
将式(5.10)代入式(5.9)，可得
s·= e¨+ce·=θ¨d－θ¨+ce·=θ¨d－f－gu－dt+ce·=－ηsgns－dt
如果选择η≥D，可得
ss·=－ηs－s·dt≤0
如果f(·)未知，可通过逼近f(·)来实现稳定控制设计。下面介绍RBF神经网络对未知项f(·)的逼近算法。

5.2.2基于RBF神经网络逼近f(·)的滑模控制 

采用RBF神经网络逼近f(·)，RBF神经网络算法为
hj=exp‖x－cj‖22b2j
f= W*Thx+ε
其中，x为网络的输入； i为网络的输入个数； j表示网络隐含层第j个节点； h=hjT为高斯函数的输出； W*为网络的理想权值； ε为网络的逼近误差； 且|ε|≤εN。

网络的输入取x=ee·T，则RBF神经网络的输出为
f⌒x= W^Thx(5.11)
其中，hx为RBF神经网络的高斯函数。

则控制输入式(5.10)可写为
u=1g－f^+θ¨d+ce·+ηsgns(5.12)
将控制律式(5.12)代入式(5.9)，可得
s·=θ¨d－f－gu－dt+ce·=θ¨d－f－－f^+θ¨d+ce·+ηsgns－dt+ce·
=－f+f^－ηsgns－dt=－f～－dt－ηsgns(5.13)
其中，
f～=f－f^=W*Thx+ε－W^Thx=W～Thx+ε(5.14)
并定义W～=W*－W^。

定义Lyapunov函数为
L=12s2+12γW～TW～
其中，γ>0。

对Lyapunov函数L求导，综合式(5.12)和式(5.13)，可得
L·=ss·+γW～TW=s［－f～－dt－ηsgns］－γW～TW
=s［－W～Thx－ε－dt－ηsgns］－γW～TW
=－W～T［shx+γW］－s［ε+dt+ηsgns］
设计自适应律为
W=－1γshx(5.15)
则
L·=－s［ε+dt+ηsgns］=－s［ε+dt］－ηs
由于逼近误差ε可以限制得足够小，取η≥εN+D，可得L·≤0。



存在η0>0，η≥η0+εN+D，使得


L·≤－η0s≤0


由于L≥0，L·≤0，从而s和W～有界。当L·≡0时，s=0，根据LaSalle不变性原理［8,9］，闭环系统为渐近稳定，当t→∞时，s→0，从而e→0，e·→0。

5.2.3仿真实例

考虑单级倒立摆动力学方程
x·1=x2
x·2=gsinx1－mlx22cosx1sinx1/mc+ml［4/3－mcos2x1/mc+m］+cosx1/mc+ml［4/3－mcos2x1/mc+m］u
其中，x1和x2分别为摆角和摆速； g=9.8m/s2为重力加速度； mc=1kg为小车质量；  m=0.1kg为摆的质量； l=0.5m为摆长的一半； u为控制输入。

取x1=θ，期望轨迹为θdt=0.1sint，系统的初始状态为π/60,0。采用控制律式(5.12)和自适应律式(5.15)，控制参数取c=15，η=0.1和自适应参数取γ=0.05。 

神经网络的结构取为251，ci和bi分别设置为-1.0-0.500.51.0和bj=0.50，网络的初始权值为0.10。仿真结果如图5.1~图5.3所示。



图5.1摆角跟踪




图5.2控制输入




图5.3fx和f^x


仿真主程序为chap5_1sim.mdl，详见附录。

5.3基于RBF逼近未知函数f(·)和g(·)的滑模控制
5.3.1引言

考虑二阶非线性系统式(5.7)，假设f(·)和g(·)都是未知的非线性函数，u∈R和y∈R分别是控制输入和系统的输出，dt为外界干扰，且满足dt≤D。

类似5.2节，设理想跟踪指令为θd，定义跟踪误差为e=θd－θ，设计滑模面为s=e·+ce，其中c>0。

在控制系统设计中，采用两个RBF神经网络分别逼近f(·)和g(·)，图5.4为基于神经网络的闭环自适应控制系统。




图5.4基于神经网络的自适应控制系统


采用RBF神经网络逼近f(x)，网络算法为
hj=exp‖x－cj‖22b2j

f(·)= W*Thfx+εf， g(·)= V*Thgx+εg
其中，x为网络的输入； i为网络的输入个数； j表示网络隐含层第j个节点； h=hjT为高斯函数的输出； W*和V*为网络的理想权值； εf和εg为网络的逼近误差，且εf≤εMf， εg≤εMg，f(·)和g(·)分别为理想RBF神经网络的输出。

定义网络的输入为x=x1x2T，则RBF神经网络的输出为
f^x=
W^Thfx，g^x= V^Thgx(5.16)

其中，hfx和hgx为RBF神经网络的高斯函数。

则控制律式(5.10)可写为
u=1g^x－f^x+θ¨d+ce·+ηsgns(5.17)
其中，η待设定。

将式(5.17)代入式(5.9)，可得
s·=e¨+ce·=θ¨d－θ¨+ce·=θ¨d－f－gu－dt+ce·
=θ¨d－f－g^u+g^－gu－dt+ce·
=θ¨d－f－g^1g^x－f^x+θ¨d+ce·+ηsgns+g^－gu－dt+ce·
=f^－f－ηsgns+g^－gu－dt=f～－ηsgns+g～u－dt
=W～Tφfx－εf－ηsgns+V～Tφgx－εgu－dt(5.18)
其中，W～=W*－W^， V～=V*－V^，且
f～=f^－f=W^Thfx－W*Thfx－εf=W～Thfx－εf
g～=g^－g=V^Thgx－V*Thgx－εg=V～Thgx－εg(5.19)
定义闭环系统Lyapunov函数为
L=12s2+12γ1W～TW～+12γ2V～TV～
其中，γ1>0，γ2>0。

对L求导，综合式(5.18)，可得
L·=ss·+1γ1W～TW+1γ2V～TV
=s{W～Thfx－εf－ηsgns
+［V～Thgx－εg］u－dt}－1γ1W～TW－1γ2V～TV
=W～Tshfx－1γ1W+V～Tshgxu－1γ2V+
s［－εf－ηsgns－εgu－dt］
设计自适应律为
W=－γ1shfx(5.20)
V=－γ2shgxu(5.21)
则
L·=s［－εf－ηsgns－εgu－dt］
=［－εf－εgu－dt］s－ηs
由于逼近误差εf和εg可以限制得足够小，如取η≥εf+εgu+dt，则可得L·≤0。



存在η0>0，η≥η0+εMf+εMgu+D，使得


L·≤－η0s≤0


由于L≥0，L·≤0，从而s、W～和V～有界。当L·≡0时，s=0，根据LaSalle不变性原理，闭环系统为渐近稳定，当t→∞时，s→0，从而e→0，e·→0。L·≤0只能保证W～和V～有界，但无法得到W～和V～收敛于零，故f～和g～无法收敛于零。


5.3.2仿真实例

被控对象为单级倒立摆，其动力学方程为 
x·1=x2
x·2=fx+gxu
其中，fx=gsinx1－mlx22cosx1sinx1/mc+ml［4/3－mcos2x1/mc+m］； gx=cosx1/(mc+m)l［4/3－mcos2x1/(mc+m)］，其中x1和x2分别为摆角和摆速； g=9.8m/s2为重力加速度； mc=1kg为小车质量；  m=0.1kg为摆的质量； l=0.5m为摆长的一半； u为控制输入。

取x1=θ，期望轨迹为θdt=0.1sint，系统的初始状态为π/60,0。网络输入取x1和x2，考虑到x1和x2的取值范围，高斯函数的参数取ci=［-1.0-0.500.51.0］和bi=5.0。网络的初始权值取0.10。控制律采用式(5.17)，自适应律采用式(5.20)和式(5.21)，参数选为γ1=10，γ2=1.0，c=5.0和η=0.01。 

仿真结果如图5.5~图5.7所示。



图5.5摆角和摆速跟踪





图5.6控制输入





图5.7f(·)和g(·)及逼近


仿真主程序为chap5_2sim.mdl，详见附录。


5.4基于神经网络最小参数学习法的自适应滑模控制

采用神经网络最小参数学习法［4,5］，取神经网络权值的上界估计值作为神经网络权值的估计值，通过设计参数估计自适应律代替神经网络权值的调整，自适应算法简单，便于实际工程应用。

5.4.1问题描述

考虑如下二阶非线性系统： 
θ¨=fθ,θ·+gθ,θ·u+dt(5.22)
其中，f为未知非线性函数； g为已知非线性函数； u∈R和y=θ∈R分别为系统的输入和输出； dt为外加干扰，dt≤D。

设位置指令为θd，令e=θ－θd，设计切换函数为
s=e·+ce(5.23)

其中，c>0，则
s·=e¨+ce·=θ¨+ce·－θ¨d=f+gu+d－θ¨d+ce·(5.24)
在实际工程中，模型不确定项f为未知，为此，控制律无法设计，需要对f进行逼近。

5.4.2基于RBF神经网络逼近的自适应控制

采用RBF神经网络对不确定项f进行自适应逼近。RBF神经网络算法为
hj=exp－‖x－cj‖22b2j,j=1,2,…,m
f=WThx+ε
其中，x为网络的输入信号； j表示网络隐含层节点的个数； h=h1,h2,…,hmT为高斯函数的输出； W为理想神经网络权值； ε为神经网络逼近误差，ε≤εN。

采用RBF神经网络逼近f，根据f的表达式，网络输入取x=θθ·T,RBF神经网络的输出为
f^x= W^Thx(5.25)
采用神经网络最小参数学习法［5］，令=‖W‖2，为正常数，^为的估计，~=^－。

设计控制律为
u=1g－12s^hTh+θ¨d－ce·－ηsgns－μs(5.26)
其中,f^为RBF神经网络的估计值； η≥εN+D； μ>0。

将控制律式(5.26)代入式(5.24)，得
s·= WTh+ε－12s^hTh－ηsgns+d－μs(5.27)
定义Lyapunov函数： 
L=12s2+12γ~2
其中,γ>0。

对L求导，并将式(5.26)和式(5.27)代入，得
L·=ss·+1γ~=sWTh+ε－12s^hTh－ηsgns+d－μs+1γ~
≤12s2hTh+12－12s2^hTh+ε+ds－ηs+1γ~－μs2
=－12s2~hTh+12+ε+ds－ηs+1γ~－μs2
=~－12s2hTh+1γ+12+ε+ds－ηs－μs2
≤~－12s2hTh+1γ+12－μs2 
设计自适应律为
=γ2s2hTh－kγ^(5.28)
其中,k>0。

则
L·≤－k~^+12－μs2≤－k2~2－2+12－μs2=－k2~2－μs2+k22+12

则


L·≤－k~^+12－μs2
≤－k2~2－2+12－μs2
=－k2~2－μs2+k22+12
≤－μs2+k22+12

可见，当s2≥1μk22+12时，L·≤0，闭环系统稳定，则闭环系统的收敛结果为

t→∞时，|s|≤1μk22+12


根据引理4.1，可得

limt→∞e≤1c1μk22+12，limt→∞e·≤21μk22+12

可见，当取μ足够大，t→∞时，s→0，从而e→0，e·→0。


注： 推导中采用了以下两个结论。

(1) s2hTh+1=s2‖W‖2hTh+1=s2‖W‖2‖h‖2+1≥s2‖WTh‖2+1≥2sWTh，即
sWTh≤12s2hTh+12 
(2) 由于~+2≥0，则~2+2~+2≥0，~2+2~^－~+2≥0，即2~^≥~2－2。

采用神经网络最小参数学习法的不足之处为： 由于采用了神经网络权值的上界作为神经网络权值的估计值，而且又利用了不等式进行了放大，所设计的控制算法过于保守。

5.4.3仿真实例

被控对象取单级倒立摆，其动态方程如下： 
x·1=x2
x·2=gsinx1－mlx22cosx1sinx1/mc+ml［4/3－mcos2x1/mc+m］+cosx1/mc+ml［4/3－mcos2x1/mc+m］u
其中，f(·)=gsinx1－mlx22cosx1sinx1/mc+ml［4/3－mcos2x1/mc+m］； g(·)=cosx1/(mc+m)l［4/3－mcos2x1/(mc+m)］； x1和x2 分别为摆角和摆速； g=9.8m/s2； mc=1kg为小车质量； mc=1kg；  m为摆杆质量； m=0.1kg； l为摆长的一半； l=0.5m； u为控制输入。 

取x1=θ，摆的角度指令为θd=0.1sint。倒立摆初始状态为π/60,0，控制律取式(5.26)，取η=0.5,μ=40,自适应律取式(5.28)，自适应参数取γ=150,k=1.0。在滑模函数中，取c=15，采用RBF神经网络逼近函数f(x1,x2)，网络结构为251，根据网络输入的实际范围来设计高斯函数的参数，参数cj和bj的取值分别为－1－0.500.51－1－0.500.51和5.0。仿真结果如图5.8和图5.9所示。仿真主程序为chap5_3sim.mdl。



图5.8角度和角速度跟踪





图5.9控制输入信号




附录仿真程序
5.2节的程序
1. 仿真主程序: chap5_1sim.mdl








2. 控制律设计程序: chap5_1ctrl.m

function ［sys，x0，str，ts］=spacemodel(t，x，u，flag)

switch flag，

case 0，

［sys，x0，str，ts］=mdlInitializeSizes;

case 1，

sys=mdlDerivatives(t，x，u);

case 3，

sys=mdlOutputs(t，x，u);

case {2，4，9}

sys=［］;

otherwise

error(［'Unhandled flag = '，num2str(flag)］);

end

function ［sys，x0，str，ts］=mdlInitializeSizes

global cij bj c

sizes = simsizes;

sizes.NumContStates  = 5;

sizes.NumDiscStates  = 0;

sizes.NumOutputs     = 2;

sizes.NumInputs      = 2;

sizes.DirFeedthrough = 1;

sizes.NumSampleTimes = 0;

sys = simsizes(sizes);

x0  = 0*ones(1，5);

str = ［］;

ts  = ［］;

cij=0.10*［-1 -0.5 0 0.5 1;

-1 -0.5 0 0.5 1］;

bj=5.0;

c=15;

function sys=mdlDerivatives(t，x，u)

global cij bj c

e=u(1);

de=u(2);

s=c*e+de;

 

xi=［e;de］;

h=zeros(5，1);

for j=1:1:5

h(j)=exp(-norm(xi-cij(:，j))^2/(2*bj^2));

end

gama=0.015;

W=［x(1) x(2) x(3) x(4) x(5)］';

for i=1:1:5

sys(i)=-1/gama*s*h(i);

end

function sys=mdlOutputs(t，x，u)

global cij bj c

e=u(1);

de=u(2);

thd=0.1*sin(t);

dthd=0.1*cos(t);

ddthd=-0.1*sin(t);

x1=thd-e;

 

s=c*e+de;

W=［x(1) x(2) x(3) x(4) x(5)］';

xi=［e;de］;

h=zeros(5，1);

for j=1:1:5

h(j)=exp(-norm(xi-cij(:，j))^2/(2*bj^2));

end

fn=W'*h;

 

g=9.8;mc=1.0;m=0.1;l=0.5;

S=l*(4/3-m*(cos(x1))^2/(mc+m));

gx=cos(x1)/(mc+m);

gx=gx/S;

 

if t<=1.5

xite=1.0;

else

xite=0.10;

end

ut=1/gx*(-fn+ddthd+c*de+xite*sign(s));

sys(1)=ut;

sys(2)=fn;

3. 控制对象程序: chap5_1plant.m

function ［sys，x0，str，ts］=s_function(t，x，u，flag)

switch flag，

case 0，

［sys，x0，str，ts］=mdlInitializeSizes;

case 1，

sys=mdlDerivatives(t，x，u);

case 3，

sys=mdlOutputs(t，x，u);

case {2， 4， 9}

sys = ［］;

otherwise

error(［'Unhandled flag = '，num2str(flag)］);

end

function ［sys，x0，str，ts］=mdlInitializeSizes

sizes = simsizes;

sizes.NumContStates  = 2;

sizes.NumDiscStates  = 0;

sizes.NumOutputs     = 2;

sizes.NumInputs      = 1;

sizes.DirFeedthrough = 0;

sizes.NumSampleTimes = 0;

sys=simsizes(sizes);

x0=［pi/60 0］;

str=［］;

ts=［］;

function sys=mdlDerivatives(t，x，u)

g=9.8;mc=1.0;m=0.1;l=0.5;

S=l*(4/3-m*(cos(x(1)))^2/(mc+m));

fx=g*sin(x(1))-m*l*x(2)^2*cos(x(1))*sin(x(1))/(mc+m);

fx=fx/S;

gx=cos(x(1))/(mc+m);

gx=gx/S;

%%%%%%%%%

dt=0*10*sin(t);

%%%%%%%%%

 

sys(1)=x(2);

sys(2)=fx+gx*u+dt;

function sys=mdlOutputs(t，x，u)

g=9.8;

mc=1.0;

m=0.1;

l=0.5;

 

S=l*(4/3-m*(cos(x(1)))^2/(mc+m));

fx=g*sin(x(1))-m*l*x(2)^2*cos(x(1))*sin(x(1))/(mc+m);

fx=fx/S;

 

sys(1)=x(1);

sys(2)=fx;

4. 作图程序: chap5_1plot.m

close all;

 

figure(1);

plot(t，y(:，1)，'k'，t，y(:，2)，'r:'，'linewidth'，2);

xlabel('time(s)');ylabel('Position tracking');

legend('ideal signal'，'practical signal');

 

figure(2);

plot(t，u(:，1)，'k'，'linewidth'，2);

xlabel('time(s)');ylabel('Control input');

 

figure(3);

plot(t，fx(:，1)，'k'，t，fx(:，2)，'r:'，'linewidth'，2);

xlabel('time(s)');ylabel('fx and estiamted fx');

legend('fx'，'estiamted fx');

5.3节的程序
1. 仿真主程序: chap5_2sim.mdl





2. 控制律设计程序: chap5_2ctrl.m

function ［sys，x0，str，ts］=spacemodel(t，x，u，flag)

switch flag，

case 0，

［sys，x0，str，ts］=mdlInitializeSizes;

case 1，

sys=mdlDerivatives(t，x，u);

case 3，

sys=mdlOutputs(t，x，u);

case {2，4，9}

sys=［］;

otherwise

error(［'Unhandled flag = '，num2str(flag)］);

end

function ［sys，x0，str，ts］=mdlInitializeSizes

global xite cij bj h c

sizes = simsizes;

sizes.NumContStates  = 10;

sizes.NumDiscStates  = 0;

sizes.NumOutputs     = 3;

sizes.NumInputs      = 5;

sizes.DirFeedthrough = 1;

sizes.NumSampleTimes = 0;

sys = simsizes(sizes);

x0  = 0.1*ones(10，1);

str = ［］;

ts  = ［］;

cij=［-1 -0.5 0 0.5 1;

-1 -0.5 0 0.5 1］;

bj=5;

h=［0，0，0，0，0］;

c=5;

xite=0.01;

function sys=mdlDerivatives(t，x，u)

global xite cij bj h c

thd=u(1);

dthd=0.1*cos(t);

ddthd=-0.1*sin(t);

 

x1=u(2);

x2=u(3);

e=thd-x1;

de=dthd-x2;

 

s=c*e+de;

 

xi=［x1;x2］;

for j=1:1:5

h(j)=exp(-norm(xi-cij(:，j))^2/(2*bj^2));

end

 

for i=1:1:5

wf(i，1)=x(i);

end

for i=1:1:5

wg(i，1)=x(i+5);

end

fxn=wf'*h';

gxn=wg'*h'+0.01;

 

ut=1/gxn*(-fxn+ddthd+xite*sign(s)+c*de);

 

gama1=10;gama2=1.0;

S1=-gama1*s*h;

S2=-gama2*s*h*ut;

for i=1:1:5

sys(i)=S1(i);

end

for j=6:1:10

sys(j)=S2(j-5);

end

 

function sys=mdlOutputs(t，x，u)   

global xite cij bj h c

thd=u(1);

dthd=0.1*cos(t);

ddthd=-0.1*sin(t);

 

x1=u(2);

x2=u(3);

e=thd-x1;

de=dthd-x2;

 

s=c*e+de;

 

for i=1:1:5

wf(i，1)=x(i);

end

for i=1:1:5

wg(i，1)=x(i+5);

end

 

xi=［x1;x2］;

for j=1:1:5

h(j)=exp(-norm(xi-cij(:，j))^2/(2*bj^2));

end

 

fxn=wf'*h';

gxn=wg'*h'+0.01;

 

ut=1/gxn*(-fxn+ddthd+xite*sign(s)+c*de);



sys(1)=ut;

sys(2)=fxn;

sys(3)=gxn;

3. 控制对象程序: chap5_2plant.m

function ［sys，x0，str，ts］=s_function(t，x，u，flag)

switch flag，

case 0，

［sys，x0，str，ts］=mdlInitializeSizes;

case 1，

sys=mdlDerivatives(t，x，u);

case 3，

sys=mdlOutputs(t，x，u);

case {2， 4， 9}

sys = ［］;

otherwise

error(［'Unhandled flag = '，num2str(flag)］);

end

function ［sys，x0，str，ts］=mdlInitializeSizes

sizes = simsizes;

sizes.NumContStates  = 2;

sizes.NumDiscStates  = 0;

sizes.NumOutputs     = 4;

sizes.NumInputs      = 1;

sizes.DirFeedthrough = 0;

sizes.NumSampleTimes = 0;

sys=simsizes(sizes);

x0=［pi/60 0］;

str=［］;

ts=［］;

function sys=mdlDerivatives(t，x，u)

g=9.8;mc=1.0;m=0.1;l=0.5;

 

S=l*(4/3-m*(cos(x(1)))^2/(mc+m));

fx=g*sin(x(1))-m*l*x(2)^2*cos(x(1))*sin(x(1))/(mc+m);

fx=fx/S;

gx=cos(x(1))/(mc+m);

gx=gx/S;

 

sys(1)=x(2);

sys(2)=fx+gx*u;

function sys=mdlOutputs(t，x，u)

g=9.8;mc=1.0;m=0.1;l=0.5;

 

S=l*(4/3-m*(cos(x(1)))^2/(mc+m));

fx=g*sin(x(1))-m*l*x(2)^2*cos(x(1))*sin(x(1))/(mc+m);

fx=fx/S;

gx=cos(x(1))/(mc+m);

gx=gx/S;

 

sys(1)=x(1);

sys(2)=x(2);

sys(3)=fx;

sys(4)=gx;


4. 作图程序: chap5_2plot.m


close all;

 

figure(1);

subplot(211);

plot(t，x(:，1)，'r'，t，x(:，2)，'k:'，'linewidth'，2);

xlabel('time(s)');ylabel('Position tracking');

legend('Ideal position signal'，'Position signal tracking');

subplot(212);

plot(t，0.1*cos(t)，'r'，t，x(:，3)，'k:'，'linewidth'，2);

xlabel('time(s)');ylabel('Speed tracking');

legend('Ideal speed signal'，'Speed signal tracking');

 

figure(2);

plot(t，u(:，1)，'r'，'linewidth'，2);

xlabel('time(s)');ylabel('Control input');

 

figure(3);

subplot(211);

plot(t，f(:，1)，'r'，t，f(:，2)，'k:'，'linewidth'，2);

xlabel('time(s)');ylabel('f and estiamted f');

legend('True f'，'Estimation f');

subplot(212);

plot(t，g(:，1)，'r'，t，g(:，2)，'k:'，'linewidth'，2);

xlabel('time(s)');ylabel('g and estimated g');

legend('True g'，'Estimation g');

5.4节的程序
1. Simulink主程序： chap5_3sim.mdl




2. 控制器S函数： chap5_3ctrl.m

function ［sys,x0,str,ts］=spacemodel(t,x,u,flag)

switch flag,

case 0,

［sys,x0,str,ts］=mdlInitializeSizes;

case 1,

sys=mdlDerivatives(t,x,u);

case 3,

sys=mdlOutputs(t,x,u);

case {2,4,9}

sys=［］;

otherwise

error(［'Unhandled flag = ',num2str(flag)］);

end

function ［sys,x0,str,ts］=mdlInitializeSizes

global cc bb c miu

sizes = simsizes;

sizes.NumContStates  = 1;

sizes.NumDiscStates  = 0;

sizes.NumOutputs     = 1;

sizes.NumInputs      = 4;

sizes.DirFeedthrough = 1;

sizes.NumSampleTimes = 0;

sys = simsizes(sizes);

x0  = 0;

str = ［］;

ts  = ［］;

cc=［-2 -1 0 1 2;

-2 -1 0 1 2］;

bb=1;

c=200;

miu=30;

function sys=mdlDerivatives(t,x,u)

global cc bb c miu

x1d=u(1);

dx1d=u(2);

x1=u(3);

x2=u(4);

 

e=x1-x1d;

de=x2-dx1d;

s=c*e+de;

 

xi=［x1;x2］;

h=zeros(5,1);

for j=1:1:5

h(j)=exp(-norm(xi-cc(:,j))^2/(2*bb*bb));

end

gama=150;

k=1.0;

sys(1)=gama/2*s^2*h'*h-k*gama*x;

 

function sys=mdlOutputs(t,x,u)

global cc bb c miu

x1d=u(1);

dx1d=u(2);

x1=u(3);

x2=u(4);

 

e=x1-x1d;

de=x2-dx1d;

thd=0.1*sin(t);

dthd=0.1*cos(t);

ddthd=-0.1*sin(t);

s=c*e+de;

 

fi=x;

xi=［x1;x2］;

h=zeros(5,1);

for j=1:1:5

h(j)=exp(-norm(xi-cc(:,j))^2/(2*bb*bb));

end

 

g=9.8;mc=1.0;m=0.1;l=0.5;

S=l*(4/3-m*(cos(x1))^2/(mc+m));

gx=cos(x1)/(mc+m);

gx=gx/S;

 

xite=0.5;

miu=40;

 

ut=1/gx*(-0.5*s*fi*h'*h+ddthd-c*de-xite*sign(s)-miu*s);

sys(1)=ut;


3. 被控对象S函数： chap5_3plant.m

function ［sys,x0,str,ts］=s_function(t,x,u,flag)

switch flag,

case 0,

［sys,x0,str,ts］=mdlInitializeSizes;

case 1,

sys=mdlDerivatives(t,x,u);

case 3,

sys=mdlOutputs(t,x,u);

case {2, 4, 9}

sys = ［］;

otherwise

error(［'Unhandled flag = ',num2str(flag)］);

end

function ［sys,x0,str,ts］=mdlInitializeSizes

sizes = simsizes;

sizes.NumContStates  = 2;

sizes.NumDiscStates  = 0;

sizes.NumOutputs     = 2;

sizes.NumInputs      = 1;

sizes.DirFeedthrough = 0;

sizes.NumSampleTimes = 0;

sys=simsizes(sizes);

x0=［pi/60 0］;

str=［］;

ts=［］;

function sys=mdlDerivatives(t,x,u)

g=9.8;mc=1.0;m=0.1;l=0.5;

S=l*(4/3-m*(cos(x(1)))^2/(mc+m));

fx=g*sin(x(1))-m*l*x(2)^2*cos(x(1))*sin(x(1))/(mc+m);

fx=fx/S;

gx=cos(x(1))/(mc+m);

gx=gx/S;

%%%%%%%%%

dt=0.1*10*sin(t);

%%%%%%%%%

 

sys(1)=x(2);

sys(2)=fx+gx*u+dt;

function sys=mdlOutputs(t,x,u)

g=9.8;

mc=1.0;

m=0.1;

l=0.5;

 

S=l*(4/3-m*(cos(x(1)))^2/(mc+m));

fx=g*sin(x(1))-m*l*x(2)^2*cos(x(1))*sin(x(1))/(mc+m);

fx=fx/S;

 

sys(1)=x(1);

sys(2)=x(2);


4. 作图程序： chap5_3plot.m

close all;

 

figure(1);

subplot(211);

plot(t,y(:,1),'r',t,y(:,2),'b:','linewidth',2);

xlabel('time(s)');ylabel('angle tracking');

legend('ideal angle','angle tracking');

subplot(212);

plot(t,0.1*cos(t),'r',t,y(:,3),'b:','linewidth',2);

xlabel('time(s)');ylabel('speed tracking');

legend('ideal speed','speed tracking');

 

figure(2);

plot(t,ut(:,1),'r','linewidth',2);

xlabel('time(s)');ylabel('Control input,u');


参考文献

［1］HUANG S,TAN K K,TONG H L,et al. Adaptive Control of Mechanical Systems Using Neural Networks［J］. IEEE Transactions on Systems,Man and Cybernetics Part C,2014,37(5): 897903.

［2］TSAI C H,CHUNG H Y,YU F M. Neurosliding mode control with its applications to seesaw systems［J］. IEEE Transactions Neural Networks,2004,15(1): 124134.

［3］EDWARDS C,SPURGEON S. Sliding mode control: theory and applications［M］.London: Taylor & Francis,1998.

［4］YANG Y S,REN J S. Adaptive fuzzy robust tracking controller design via small gain approach and its application［J］.IEEE Transactions on Fuzzy Systems,2003,11(6): 783795.

［5］CHEN B,LIU X P,LIU K F,et al. Direct adaptive fuzzy control of nonlinear strictfeedback systems［J］. Automatica,2009,45: 15301535.

［6］刘金琨.滑模变结构控制MATLAB仿真［M］.2版.北京: 清华大学出版社，2012.

［7］LIU J K,LU Y.Adaptive RBF neural network control of robot with actuator nonlinearities［J］.Journal of Control Theory and Applications,2010,8(2): 150156.

［8］LASALLE J,LEFSCHETZ S.Stability by Lyapunovs direct method［M］.New York: Academic Press,1961.

［9］HASSAN K H.Nonlinear systems［M］.3rd ed.New Jersey: Prentice Hall,2002.