基于趋近律的滑模鲁棒控制simulink仿真

滑模控制刚入门菜鸟一枚，找个实例练练手。参考刘金琨老师的《滑模变结构控制MATLAB仿真》中的基于趋近律的滑模鲁棒控制，对部分地方做出了修改。

考虑如下的被控对象：，其中，b>0，d(t)为外部干扰。

设计滑模面为：，其中

求滑模面的一阶导数：

 采用指数趋近律：

 可以计算出系统所需的控制器（过程略）：

 但是由于干扰d未知，上述控制器无法实现，书中采用干扰的界来设计控制器。

设dc为干扰d的界相关的正实数，此时控制器为：

 利用Lyapunov第二法证明系统稳定性，

设计Lyapunov函数：，对时间进行微分，得到：

 要使系统稳定就得使负定，这就需要对dc进行合适地选取。

dc选取原则：

假设dL≤d≤dU，其中dL和dU分别为d的上下界，

①当s>0时，为了使<0，取dc=dL；

②当s<0时，为了使 >0，取dc=dU；

取d1=（du-dL）/2，d2=（du+dL）/2，那么可以设计满足上述两个条件的dc为：

这样可以始终保证 <0，从而保证系统的稳定。

但是刘金琨老师在书中使用的是符号函数sign(s)，这里我改为了饱和函数sat(s)，从而有效地减小了抖震。

改为饱和函数后，

控制器为：

dc为：

 系统各参数取值： f(θ,t)=25，b=133，d(t)=10sin(πt)，θd=sint，初始状态为[-0.15 -0.15]；

控制率中各参数取值：c=25（书中取的15）,ε=0.5，k=10。

控制器部分的代码(ctrl.m)为：

function [sys,x0,str,ts,simStateCompliance] = ctrl(t,x,u,flag,pa)

switch flag
case 0 % 初始化
[sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes;

case 1 % 连续时间导数
sys=mdlDerivatives(t,x,u);

case 2 % 更新离散状态量
sys=mdlUpdate(t,x,u);

case 3 % 计算输出
sys=mdlOutputs(t,x,u,pa);

case 4 % 计算下一步采样时刻
sys=mdlGetTimeOfNextVarHit(t,x,u);

case 9 % 结束仿真
sys=mdlTerminate(t,x,u);

otherwise % 未知flag值
DAStudio.error('Simulink:blocks:unhandledFlag', num2str(flag));
end % S函数主程序结束

function [sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes


sizes = simsizes;                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                  simsizes;
sizes.NumContStates = 0; %连续状态个数
sizes.NumDiscStates = 0; %离散状态个数
sizes.NumOutputs = 1; %输出量个数
sizes.NumInputs = 3; %输入量个数
sizes.DirFeedthrough = 1; %直接馈通标志
sizes.NumSampleTimes = 1; % 至少有一个采样时刻
sys = simsizes(sizes);

x0 = []; % 状态初始化
str = []; % str 始终为空
ts = [0 0];% 初始化采样时间


simStateCompliance = 'UnknownSimState';

function sys=mdlDerivatives(t,x,u)

sys = [];


function sys=mdlUpdate(t,x,u)

sys = [];
% 子函数 mdlUpdate 结束


function sys=mdlOutputs(t,x,u,pa)
thetad=u(1);          %θd
thetad_dot=cos(t);    %θd`
thetad_ddot=-sin(t);  %θd``
theta=u(2);           %θ
theta_dot=u(3);       %θ`

b=pa.b;
epsilon=pa.epsilon;
k=pa.k;
c=pa.c;
dU=pa.dU;
dL=pa.dL;

e=thetad-theta;
de=thetad_dot-theta_dot;
s=c*e+de;
delta=pa.delta;
if abs(s)>delta
    sat=sign(s);
else
    sat=s/delta;
end

fx=25*theta;
d1=(dU-dL)/2;
d2=(dU+dL)/2;
dc=d2-d1*sat;


ut=1/b*(epsilon*sat+k*s+c*(thetad_dot-thetad)+thetad_ddot+fx-dc);

sys(1) = ut;


function sys=mdlGetTimeOfNextVarHit(t,x,u)

sampleTime = 1; % 例子。设置下一个采样时刻为1s后。
sys = t + sampleTime;


function sys=mdlTerminate(t,x,u)

sys = [];

% 子函数 mdlTerminate结束

系统部分的代码(plant.m)为：

function [sys,x0,str,ts,simStateCompliance] = plant(t,x,u,flag,pa)

switch flag
case 0 % 初始化
[sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes;

case 1 % 连续时间导数
sys=mdlDerivatives(t,x,u,pa);

case 2 % 更新离散状态量
sys=mdlUpdate(t,x,u);

case 3 % 计算输出
sys=mdlOutputs(t,x,u);

case 4 % 计算下一步采样时刻
sys=mdlGetTimeOfNextVarHit(t,x,u);

case 9 % 结束仿真
sys=mdlTerminate(t,x,u);

otherwise % 未知flag值
DAStudio.error('Simulink:blocks:unhandledFlag', num2str(flag));
end % S函数主程序结束

function [sys,x0,str,ts,simStateCompliance]=mdlInitializeSizes


sizes = simsizes;                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                                               simsizes;
sizes.NumContStates = 2;  %连续状态个数
sizes.NumDiscStates = 0;  %离散状态个数
sizes.NumOutputs = 2;     %输出量个数
sizes.NumInputs = 1;      %输入量个数
sizes.DirFeedthrough = 0; %直接馈通标志
sizes.NumSampleTimes = 1; % 至少有一个采样时刻
sys = simsizes(sizes);

x0 = [-0.15 -0.15]; % 状态初始化
str = []; % str 始终为空
ts = [0 0];% 初始化采样时间


simStateCompliance = 'UnknownSimState';

function sys=mdlDerivatives(t,x,u,pa)
b=pa.b;
ut=u(1);

theta_ddot = -25*x(2)+b*ut+10*sin(pi*t);

sys = [x(2);theta_ddot];


function sys=mdlUpdate(t,x,u)
sys = [];

function sys=mdlOutputs(t,x,u)

sys = [x(1); x(2)];

function sys=mdlGetTimeOfNextVarHit(t,x,u)

sampleTime = 1; % 例子。设置下一个采样时刻为1s后。
sys = t + sampleTime;

function sys=mdlTerminate(t,x,u)

sys = [];

% 子函数 mdlTerminate结束

 参数部分（parameter.m）代码：

pa=struct('c',25,'b',133,'dL',-10,'dU',10,'epsilon',0.5,'k',10,'delta',0.05);

 simulink主程序：

 运行后得到控制器输出：

实际输出和理想输出的拟合情况还是比较理想的：