本文章以最简单的二阶系统为例，介绍其simulink仿真实现和m代码实现

案例中的二阶系统如下所示

 经典ADRC的基本结构如下：

 本案例中的simulink仿真整体结构（为便于理解，结构图与上述ADRC整体结构类似）

仿真参数初始化所需m文件程序（文章尾部附有本仿真模型及m，文件支持matlab2017b以上）

%-----------ADRC参数初始化------------%
%参数初始化
%跟踪微分器
r=100;%表示跟踪快慢
h0=5*h;%h0代表信号的平滑程度（滤波效果）
v1_last=0;
v2_last=0;
v0_last=0;
%扩张状态观测器
beta01=10;
beta02=200;
beta03=30;
alpha1=0.5;%文献里给定值
alpha2=0.25;%给定值
delta=0.0025;
b=1;
z1_last=0;
z2_last=0;
z3_last=0;
%非线性误差反馈
nlsef_alpha1=0.7;
nlsef_alpha2=1;
%被控对象初始化
temp_y=[0.5;0];
u_last=0;

TD微分跟踪器结构如下所示

 图中fst函数利用matlab_function搭建，内部代码如下所示

%fst函数
function fn=fst(x1,x2,r,h)
d=h*r;
d0=h*d;
y=x1+h*x2;
a0=sqrt(d^2+8*r*abs(y));
if abs(y)<=d0
    a=x2+y/h;
else
    a=x2+0.5*(a0-d)*sgn(y);
end
fn=-r*sat(a,d);
end
%符号函数
function y=sgn(x)
if x>0
    y=1;
elseif x<0
    y=-1;
else
    y=0
end
end
%sat函数
function y=sat(a,d)
if abs(a)<=d
    y=a/d;
else
    y=sgn(a);
end
end
%符号函数

ESO扩张状态观测器结构如下所示

 图中fal函数利用matlab_function搭建，内部代码如下所示

function y=fal(e,alpha,delta)
if abs(e)>delta
    y=sign(e)*abs(e)^alpha;
else
    y=e/(delta^(1-alpha));
end
end

BLSEF非线性误差反馈控制律结构如下所示

  图中fal函数利用matlab_function搭建，内部代码同上

典型被控系统结构如下所示

 阶跃信号下输入/输出/跟踪信号仿真结果

 系统输出的微分观测效果

 系统非线性项的观测效果

系统仿真的m文件实现（想要深入理解实现过程的朋友可以尝试复现该代码）

 m代码（ADRC控制器）

clc;clear all;close all;
%设定运行时间
time=10;
%设定仿真步长
h=0.01;
%时间定义
t=0.01:h:time;
%跟踪信号
v0=zeros(1,time/h);
for i=time/h/2+1:time/h;
    v0(i)=1;
end
rand_noise=0.05*randn(1,time/h);
%跟踪信号中加入随机噪声
vn=v0+rand_noise;
%-----------ADRC------------%
%参数初始化
%跟踪微分器
r=100;%表示跟踪快慢
h0=5*h;%h0代表信号的平滑程度（滤波效果）
v1_last=0;
v2_last=0;
v0_last=0;
%扩张状态观测器
beta01=10;
beta02=200;
beta03=30;
alpha1=0.5;%文献里给定值
alpha2=0.25;%给定值
delta=0.0025;
b=1;
z1_last=0;
z2_last=0;
z3_last=0;
%非线性误差反馈
nlsef_alpha1=0.7;
nlsef_alpha2=1;
%被控对象初始化
temp_y=[0.5;0];
u_last=0;
%----ADRC正式开始------%
for k=1:time/h
    %第一轮迭代处理
    %两个参数分别为控制量和当前时间
    parameter1=u_last;
    parameter2=k*h;
    tSpan=[0 h];
    %利用龙格库塔法求解微分方程
    [~,total_y]=ode45('PlantModel',tSpan,temp_y,[],parameter1,parameter2);
    %total_state里面的元素都是龙格塔库一点点计算的结果，直接使用最后一列，即计算结果即可
    temp_y=total_y(length(total_y),:);%寻访最后一行，全部列的元素
    %记录下输出和输出的微分
    y(k)=temp_y(1);
    dy(k)=temp_y(2);
    %---跟踪微分器TD----%
    v1(k)=v1_last+h*v2_last;
    v2(k)=v2_last+h*fst(v1_last-vn(k),v2_last,r,h0);
    x3(k)=-v1_last^2;
    v1_last=v1(k);
    v2_last=v2(k);
    v0_last=vn(k);
    %----扩张状态观测器--%
    e=z1_last-y(k);
    z1(k)=z1_last+h*(z2_last-beta01*e);
    z2(k)=z2_last+h*(z3_last-beta02*(fal(e,alpha1,delta))+b*u_last);
    z3(k)=z3_last-h*beta03*(fal(e,alpha2,delta));
    z1_last=z1(k);
    z2_last=z2(k);
    z3_last=z3(k);
    %---非线性误差反馈----%
    e1(k)=v1(k)-z1(k);
    e2(k)=v2(k)-z2(k);
    u0(k)=beta01*fal(e1(k),nlsef_alpha1,delta)+beta02*fal(e2(k),nlsef_alpha2,delta);
    u(k)=u0(k)-z3(k)/b;
    u_last=u(k);
end
figure(1);
plot(t,u,'r');
figure(2);
subplot(311);
plot(t,z1,'r',t,y,'k',t,vn,'b','linewidth',2);
xlabel('time(s)');ylabel('z1,y');
legend('目标输出信号','估计输出信号','实际输出信号');
subplot(312);
plot(t,z2,'r',t,dy,'k','linewidth',2);
xlabel('time(s)'),ylabel('z2,dy');
legend('估计输出微分信号','实际输出微分信号');
subplot(313);
plot(t,z3,'r',t,x3,'k','linewidth',2);
xlabel('time(s)'),ylabel('z3,x3');
legend('估计扰动','实际扰动');
%---------函数部分---------%
%sat函数
function y=sat(a,d)
if abs(a)<=d
    y=a/d;
else
    y=sgn(a);
end
end
%符号函数
function y=sgn(x)
if x>0
    y=1;
elseif x<0
    y=-1;
else
    y=0
end
end
%fst函数
function fn=fst(x1,x2,r,h)
d=h*r;
d0=h*d;
y=x1+h*x2;
a0=sqrt(d^2+8*r*abs(y));
if abs(y)<=d0
    a=x2+y/h;
else
    a=x2+0.5*(a0-d)*sgn(y);
end
fn=-r*sat(a,d);
end
%fal函数
function y=fal(e,alpha,delta)
if abs(e)>delta
    y=abs(e)^alpha*sign(e);
else
    y=e/(delta^(1-alpha));
end
end

系统模型

%系统方程
function dy=PlantModel(t,y,flag,p1,p2)
u=p1;
time=p2;
dy=zeros(2,1);
dy(1)=y(2);
dy(2)=-y(1)^2+u;
end

m文件实现的仿真结果

仿真源文件网盘链接

链接：https://pan.baidu.com/s/1bZ-teW6aMX2XjmLkBW1ozg 
提取码：4m5p