untitled.m

在MATLAB中

num=20;
den=[1.6 4.4 1];
[a b c d]=tf2ss(num,den);
x=[0;0];
T=0.01;h=T;
N=250;
R=1.5*ones(1,N);%参考输入
%定义输入输出变量与隶属度函数
a=newfis('Simple');
a=addvar(a,'input','e',[-6 6]);
a=addmf(a,'input',1,'NB','trapmf',[-6,-6,-5,-3]);
a=addmf(a,'input',1,'NS','trapmf',[-5,-3,-2,0]);
a=addmf(a,'input',1,'ZR','trimf',[-2,0,2]);
a=addmf(a,'input',1,'PS','trapmf',[0,2,3,5]);
a=addmf(a,'input',1,'PB','trapmf',[3,5,6,6]);
a=addvar(a,'input','de',[-6,6]);
a=addmf(a,'input',2,'NB','trapmf',[-6,-6,-5,-3]);
a=addmf(a,'input',2,'NS','trapmf',[-5,-3,-2,0]);
a=addmf(a,'input',2,'ZR','trimf',[-2,0,2]);
a=addmf(a,'input',2,'PS','trapmf',[0,2,3,5]);
a=addmf(a,'input',2,'PB','trapmf',[3,5,6,6]);
a=addvar(a,'output','u',[-3,3]);
a=addmf(a,'output',1,'NB','trapmf',[-3,-3,-2,-1]);
a=addmf(a,'output',1,'NS','trimf',[-2,-1,0]);
a=addmf(a,'output',1,'ZR','trimf',[-1,0,1]);
a=addmf(a,'output',1,'PS','trimf',[0,1,2]);
a=addmf(a,'output',1,'PB','trapmf',[1,2,3,3]);
%模糊规则矩阵
rr=[5 5 4 4 3
    5 4 4 3 3
    4 4 3 3 2
    4 3 3 2 2
    3 3 2 2 1];
rr=zeros(prod(size(rr)),3);
k=1;
for i=1:size(rr,1)
    for j=1:size(rr,2)
        r1(k,:)=[i,j,rr(i,j)];
        k=k+1;
    end
end
[r,s]=size(r1);
r2=ones(r,2);
rulelist=[r1,r2];
a=addrule(a,rulelist);
%模糊控制系统仿真
e=0;
de=0;
ke=30;
kd=20;
ku=1;
for k=1:N
    e1=ke*e;
    de1=kd*de;
    if e1>=6
        e1=6;
    elseif e1<=-6
        e1=-6;
    end
    if de1>=6
        de1=6;
    elseif de1<=-6
        de1=-6;
    end
%模糊推理，计算被控对象输入
in=[e1 de1];
u=ku*evalfis(in,a);
uu(1,k)=u;
%控制作用于被控系统，计算系统输出
k0=a1*x+b*u;
k1=a1*(x+h*k0/2)+b*u;
k2=a1*(x+h*k1/2)+b*u;
k3=a1*(x+h*k2)+b*u;
x=x+(k0+2*k1+2*k2+k3)*h/6;
y=c*x+d*u;
yy(1,k)=y;
%计算系统输出误差及误差变化率
e1=e;
e=y-R(1,k);
de=(e-e1)/T;
end
%模糊控制输出曲线
kk=[1:N]*T;
figure(1);
plot(kk,R,'k',kk,yy,'r');
grid on