with_dc.m

matlab仿真模糊控制器

%Fuzzy Controller
clear all;close all;

%**************************模糊控制单元*********************************%
a=newfis('fuzz_ls');                                    % 新建模糊推理系统

f1=1.0;
a=addvar(a,'input','e',[-3*f1,3*f1]);                     % 添加输入变量  e
a=addmf(a,'input',1,'NB','zmf',[-3*f1,-1*f1]);
a=addmf(a,'input',1,'NM','trimf',[-3*f1,-2*f1,0]);
a=addmf(a,'input',1,'NS','trimf',[-3*f1,-1*f1,1*f1]);
a=addmf(a,'input',1,'Z','trimf',[-2*f1,0,2*f1]);
a=addmf(a,'input',1,'PS','trimf',[-1*f1,1*f1,3*f1]);
a=addmf(a,'input',1,'PM','trimf',[0,2*f1,3*f1]);
a=addmf(a,'input',1,'PB','smf',[1*f1,3*f1]);

f2=1.0;
a=addvar(a,'input','ec',[-3*f2,3*f2]);                     % 添加输入变量 ec
a=addmf(a,'input',2,'NB','zmf',[-3*f2,-1*f2]);
a=addmf(a,'input',2,'NM','trimf',[-3*f2,-2*f2,0]);
a=addmf(a,'input',2,'NS','trimf',[-3*f2,-1*f2,1*f2]);
a=addmf(a,'input',2,'Z','trimf',[-2*f2,0,2*f2]);
a=addmf(a,'input',2,'PS','trimf',[-1*f2,1*f2,3*f2]);
a=addmf(a,'input',2,'PM','trimf',[0,2*f2,3*f2]);
a=addmf(a,'input',2,'PB','smf',[1*f2,3*f2]);

f3=1.0;
a=addvar(a,'output','u',[-3*f3,3*f3]);                     % 添加输出变量 u
a=addmf(a,'output',1,'NB','zmf',[-3*f3,-1*f3]);
a=addmf(a,'output',1,'NM','trimf',[-3*f3,-2*f3,0]);
a=addmf(a,'output',1,'NS','trimf',[-3*f3,-1*f3,1*f3]);
a=addmf(a,'output',1,'Z','trimf',[-2*f3,0,2*f3]);
a=addmf(a,'output',1,'PS','trimf',[-1*f3,1*f3,3*f3]);
a=addmf(a,'output',1,'PM','trimf',[0,2*f3,3*f3]);
a=addmf(a,'output',1,'PB','smf',[1*f3,3*f3]);

%Each rule is a PD rule: error=rin-yout (nagative feedback)
rulelist=[1 1 7 1 1;                                        % 编辑模糊规则
          1 2 7 1 1;
          1 3 6 1 1;
          1 4 6 1 1;
          1 5 5 1 1;
          1 6 5 1 1;
          1 7 4 1 1;
   
          2 1 7 1 1;
          2 2 6 1 1;
          2 3 6 1 1;
          2 4 5 1 1;
          2 5 5 1 1;
          2 6 4 1 1;
          2 7 3 1 1;
          
          3 1 6 1 1;
          3 2 6 1 1;
          3 3 5 1 1;
          3 4 5 1 1;
          3 5 4 1 1;
          3 6 3 1 1;
          3 7 3 1 1;
          
          4 1 6 1 1;
          4 2 5 1 1;
          4 3 5 1 1;
          4 4 4 1 1;
          4 5 3 1 1;
          4 6 3 1 1;
          4 7 2 1 1;
          
          5 1 5 1 1;
          5 2 5 1 1;
          5 3 4 1 1;
          5 4 3 1 1;
          5 5 3 1 1;
          5 6 2 1 1;
          5 7 2 1 1;
          
          6 1 5 1 1;
          6 2 5 1 1;
          6 3 4 1 1;
          6 4 3 1 1;
          6 5 2 1 1;
          6 6 2 1 1;
          6 7 1 1 1;
       
          7 1 4 1 1;
          7 2 3 1 1;
          7 3 3 1 1;
          7 4 2 1 1;
          7 5 2 1 1;
          7 6 1 1 1;
          7 7 1 1 1];
          
a=addrule(a,rulelist);                                        %添加模糊规则
%showrule(a)                                                  %显示模糊规则

a1=setfis(a,'DefuzzMethod','centroid');                       %重心法反模糊化
writefis(a1,'ls');                                            %保存
a2=readfis('ls');                                             %读取
%**********************************************************************%

%*******************u2有阶跃的隶属度函数和模糊规则设计********************%
dca=newfis('fuzz_ls1');

% f1=1.0;
% dca=addvar(dca,'input','au1',[-3*f1,3*f1]);               % Parameter au1
% dca=addmf(dca,'input',1,'NB','zmf',[-2.5*f1,-1.5*f1]);
% dca=addmf(dca,'input',1,'NS','trimf',[-2.5*f1,-1.5*f1,0*f1]);
% dca=addmf(dca,'input',1,'Z','trimf',[-1.5*f1,0*f1,1.5*f1]);
% dca=addmf(dca,'input',1,'PS','trimf',[0*f1,1.5*f1,2.5*f1]);
% dca=addmf(dca,'input',1,'PB','smf',[1.5*f1,2.5*f1]);
% 
% f2=1.0;
% dca=addvar(dca,'input','ae2',[-3*f1,3*f1]);               % Parameter ae2
% dca=addmf(dca,'input',2,'NB','zmf',[-2.5*f1,-1.5*f1]);
% dca=addmf(dca,'input',2,'NS','trimf',[-2.5*f1,-1.5*f1,0*f1]);
% dca=addmf(dca,'input',2,'Z','trimf',[-1.5*f1,0*f1,1.5*f1]);
% dca=addmf(dca,'input',2,'PS','trimf',[0*f1,1.5*f1,2.5*f1]);
% dca=addmf(dca,'input',2,'PB','smf',[1.5*f1,2.5*f1]);
% 
% 
% %f3=1.4;
% f3=1.0;
% dca=addvar(dca,'output','du12',[-3*f1,3*f1]);             % Parameter du12
% dca=addmf(dca,'output',1,'NB','zmf',[-2.5*f1,-1.5*f1]);
% dca=addmf(dca,'output',1,'NS','trimf',[-2.5*f1,-1.5*f1,0*f1]);
% dca=addmf(dca,'output',1,'Z','trimf',[-1.5*f1,0*f1,1.5*f1]);
% dca=addmf(dca,'output',1,'PS','trimf',[0.5*f1,1.5*f1,2.5*f1]);
% dca=addmf(dca,'output',1,'PB','smf',[1.5*f1,2.5*f1]);

f1=1.0;
dca=addvar(dca,'input','au1',[-3*f1,3*f1]);                 % Parameter au1
dca=addmf(dca,'input',1,'NB','zmf',[-2*f1,-1*f1]);
dca=addmf(dca,'input',1,'NS','trimf',[-2.0*f1,-1*f1,0*f1]);
dca=addmf(dca,'input',1,'Z','trimf',[-1.0*f1,0*f1,1.0*f1]);
dca=addmf(dca,'input',1,'PS','trimf',[0*f1,1*f1,2.0*f1]);
dca=addmf(dca,'input',1,'PB','smf',[1*f1,2*f1]);

f2=1.0;
dca=addvar(dca,'input','ae2',[-3*f1,3*f1]);                 % Parameter ae2
dca=addmf(dca,'input',2,'NB','zmf',[-2*f1,-1*f1]);
dca=addmf(dca,'input',2,'NS','trimf',[-2*f1,-1*f1,0*f1]);
dca=addmf(dca,'input',2,'Z','trimf',[-1*f1,0*f1,1*f1]);
dca=addmf(dca,'input',2,'PS','trimf',[0*f1,1*f1,2*f1]);
dca=addmf(dca,'input',2,'PB','smf',[1*f1,2*f1]);


%f3=1.4;
f3=0;
dca=addvar(dca,'output','du12',[-3*f1,3*f1]);              % Parameter du12
dca=addmf(dca,'output',1,'NB','zmf',[-2*f1,-1*f1]);
dca=addmf(dca,'output',1,'NS','trimf',[-2.0*f1,-1*f1,0*f1]);
dca=addmf(dca,'output',1,'Z','trimf',[-1.0*f1,0*f1,1.0*f1]);
dca=addmf(dca,'output',1,'PS','trimf',[0.5*f1,1*f1,2.0*f1]);
dca=addmf(dca,'output',1,'PB','smf',[1*f1,2*f1]);

rulelist1=[1 1 1 1 1;                                      % Edit rule base
          1 2 2 1 1;
          1 3 2 1 1;
          1 4 3 1 1;
          1 5 3 1 1;        
   
          2 1 2 1 1;
          2 2 2 1 1;
          2 3 3 1 1;
          2 4 3 1 1;
          2 5 4 1 1;          
          
          3 1 2 1 1;
          3 2 3 1 1;
          3 3 3 1 1;
          3 4 4 1 1;
          3 5 4 1 1;          
          
          4 1 3 1 1;
          4 2 3 1 1;
          4 3 4 1 1;
          4 4 4 1 1;
          4 5 5 1 1;          
          
          5 1 3 1 1;
          5 2 4 1 1;
          5 3 4 1 1;
          5 4 5 1 1;
          5 5 5 1 1];          
dca=addrule(dca,rulelist1);
%showrule(a)                                  % Show fuzzy rule base

dca1=setfis(dca,'DefuzzMethod','centroid');   % Defuzzy
writefis(dca1,'ls1');                         % save to fuzzy file "ljk.fis" 
                                              % simulated with fuzzy tool
dca2=readfis('ls1');

%**********************************************************************%
     
%****************u2有阶跃的隶属度函数和模糊规则设计***********************%

dcb=newfis('fuzz_ls2');

f1=1.0;
dcb=addvar(dcb,'input','au2',[-3*f1,3*f1]);                 % Parameter au2
dcb=addmf(dcb,'input',1,'NB','zmf',[-2.5*f1,-1.5*f1]);
dcb=addmf(dcb,'input',1,'NS','trimf',[-2.5*f1,-1.5*f1,0*f1]);
dcb=addmf(dcb,'input',1,'Z','trimf',[-1.5*f1,0*f1,1.5*f1]);
dcb=addmf(dcb,'input',1,'PS','trimf',[0*f1,1.5*f1,2.5*f1]);
dcb=addmf(dcb,'input',1,'PB','smf',[1.5*f1,2.5*f1]);

f2=1.0;
dcb=addvar(dcb,'input','ae1',[-3*f1,3*f1]);                 % Parameter ae1
dcb=addmf(dcb,'input',2,'NB','zmf',[-2.5*f1,-1.5*f1]);
dcb=addmf(dcb,'input',2,'NS','trimf',[-2.5*f1,-1.5*f1,0*f1]);
dcb=addmf(dcb,'input',2,'Z','trimf',[-1.5*f1,0*f1,1.5*f1]);
dcb=addmf(dcb,'input',2,'PS','trimf',[0*f1,1.5*f1,2.5*f1]);
dcb=addmf(dcb,'input',2,'PB','smf',[1.5*f1,2.5*f1]);


%f3=1.4;
f3=1.0;
dcb=addvar(dcb,'output','du12',[-3*f1,3*f1]);               % Parameter du21
dcb=addmf(dcb,'output',1,'NB','zmf',[-2.5*f1,-1.5*f1]);
dcb=addmf(dcb,'output',1,'NS','trimf',[-2.5*f1,-1.5*f1,0*f1]);
dcb=addmf(dcb,'output',1,'Z','trimf',[-1.5*f1,0*f1,1.5*f1]);
dcb=addmf(dcb,'output',1,'PS','trimf',[0.5*f1,1.5*f1,2.5*f1]);
dcb=addmf(dcb,'output',1,'PB','smf',[1.5*f1,2.5*f1]);

rulelist2=[1 1 5 1 1;                                      % Edit rule base
          1 2 4 1 1;
          1 3 4 1 1;
          1 4 3 1 1;
          1 5 3 1 1;        
   
          2 1 4 1 1;
          2 2 4 1 1;
          2 3 3 1 1;
          2 4 3 1 1;
          2 5 2 1 1;          
          
          3 1 4 1 1;
          3 2 3 1 1;
          3 3 3 1 1;
          3 4 2 1 1;
          3 5 2 1 1;          
          
          4 1 3 1 1;
          4 2 3 1 1;
          4 3 2 1 1;
          4 4 2 1 1;
          4 5 1 1 1;          
          
          5 1 3 1 1;
          5 2 2 1 1;
          5 3 2 1 1;
          5 4 1 1 1;
          5 5 1 1 1];

dcb=addrule(dcb,rulelist2);
%showrule(a)                                 % Show fuzzy rule base

dcb1=setfis(dcb,'DefuzzMethod','centroid');  % Defuzzy
writefis(dcb1,'ls2');                        % save to fuzzy file "ljk.fis" 
                                             % simulated with fuzzy tool
dcb2=readfis('ls2');
%*********************************************************************%

%******************* Using Fuzzy Controller****************************%
e1_1=0.0;
ec1_1=0.0;
e2_1=0.0;
ec2_1=0.0;

u1(1)=0.0; 
u2(1)=0.0;

X1(1)=0.0;
X2(1)=0.0;

% G=[0.9986 0;0 0.9972];
% H=[0.2776 0.2776;-5.5487 5.5487];
% C=[1 0;0 1];

G=[1 0;0 1];
H=[0.5 0.5;-4 2];
C=[1 0;0 1];


yout1(1)=0;
yout2(1)=0;
        
ts=0.5;
%********************Start of Control**********************************%

for k=2:1:600
     T(k)=k*ts;
   rin1(k)=0.5;     
   rin2(k)=0.5;

    if k>=200                                 %100s时rin1有阶跃变化
       rin1(k)=1;   
    end
    
    if k>=400                                 %200s时rin2有阶跃变化
       rin2(k)=1;   
    end
error1(k)=yout1(k-1)-rin1(k);                 %计算偏差
error2(k)=yout2(k-1)-rin2(k);

e1_1=error1(k);
ec1_1=error1(k)-error1(k-1);

e2_1=error2(k);
ec2_1=error2(k)-error2(k-1);

u1(k)=evalfis([e1_1 ec1_1],a2)/3;             %Using fuzzy inference
u2(k)=evalfis([e2_1 ec2_1],a2)/3;             %Using fuzzy inference

%**********************模糊控制补偿***************************************%
au1=u1(k)-u1(k-1);
ae2=error2(k)-error2(k-1);
du12(k)=evalfis([au1 ae2],dca2)/4; 

au2=u2(k)-u2(k-1);
ae1=error1(k)-error1(k-1);
du21(k)=evalfis([au2 ae1],dcb2)/4;

u1(k)=u1(k)+du21(k);                         %Using fuzzy inference
u2(k)=u2(k)+du12(k);                         %Using fuzzy inference
%************************************************************************%

%********************耦合单元********************************************%
X1(k)=G(1,1)*X1(k-1)+G(1,2)*X2(k-1)+H(1,1)*u1(k)+H(1,2)*u2(k);
X2(k)=G(2,1)*X1(k-1)+G(2,2)*X2(k-1)+H(2,1)*u1(k)+H(2,2)*u2(k);

yout1(k)=C(1,1)*X1(k); 
yout2(k)=C(2,2)*X2(k);
%************************************************************************%

end
%************************************************************************%
  %保存 调用 
save '\ls\with_dc11.txt' yout1 -ascii 
load '\ls\with_dc11.txt'
g1=with_dc11;

save '\ls\with_dc21.txt' yout2 -ascii
load '\ls\with_dc21.txt'
g2=with_dc21;

save '\ls\with_dc12.txt' u1 -ascii 
load '\ls\with_dc12.txt'
h1=with_dc12;

save '\ls\with_dc22.txt' u2 -ascii
load '\ls\with_dc22.txt'
h2=with_dc22;

load '\ls\without_dc11.txt'
g3=without_dc11;
load '\ls\without_dc12.txt'
h3=without_dc12;
load '\ls\without_dc21.txt'
g4=without_dc21;
load '\ls\without_dc22.txt'
h4=without_dc22;
   %  结束
%************************************************************************%

%************************************************************************%
%    输出图形
figure(1);
subplot 211;
plot(T,rin1,'b',T,g1,'r',T,g3,':black');
xlabel('Time(second)');ylabel('rin1,yout1');
subplot 212;
plot(T,rin2,'b',T,g2,'r',T,g4,':black');
xlabel('Time(second)');ylabel('rin2,yout2');

figure(2);
subplot 211;
plot(T,h1,'r',T,h3,':black');
xlabel('Time(second)');ylabel('u1');
subplot 212;
plot(T,h2,'r',T,h4,':black');
xlabel('Time(second)');ylabel('u2');
%    结束
%************************************************************************%