rfnnic_yu_5.m

利用遞迴式模糊類神經網路(recerrent neural networks system identification, RFNN)進行系統辨識, 自己寫的請多包含

%Adaptive PID control based on RBF Identification
clear all;
close all;
clc;
format('long')


xite=0.08;
xiteC=0.06;
x=[0; 0];
xC=[0; 0];

ci=[-1, -.5, 0, .5, 1;...
    -1, -.5, 0, .5, 1];
%ci=ones(2,5);
bi=[ones(1,5); ones(1,5)];
w=0.001*ones(25,1); 
wi=0.1*ones(2,5);

ciC=[-1.5, -.5, 0, .5, 1.5;...
    -1.5, -.5, 0, .5, 1.5];
%ci=ones(2,5);
biC=[ones(1,5); ones(1,5)];
wC=0.01*ones(25,1); 
wiC=0.1*ones(2,5);

% initial feedback output
O12_1 = zeros(2,5)
O12C_1 = zeros(2,5)

u_1=0;
y_1=0;
error_1 = 0;
error_2 = 0;

tf=2;
ts=0.001;
n=floor(tf/ts);

for k=1:n
   time(k)=k*ts;
   rin(k)=1.0*sin(sin(2*pi*k*ts));
   
   yout(k)=(-0.1*y_1+u_1)/(1+y_1^2);  %Nonlinear plant
   
% RFNNI start --------------------------------------------- 
  % Input layer
        for i1 = 1:2
           for j1=1:5
              O12(i1,j1) = O12_1(i1,j1)*wi(i1,j1) + x(i1,end);
           end
        end
   % Fuzzification
  for i = 1:2 
   for j=1:5
      f2(i,j)=exp(-norm(O12(i,j)-ci(i,j))^2/(bi(i,j)^2));
   end
  end
   % Rules
 for j = 1:5
    m1(j)=f2(1,j);
    m2(j)=f2(2,j);
 end
  
   for i = 1:5
      for j = 1:5
         ff3(i,j) = m2(i)*m1(j); 
      end
   end
  
   f3 = [ff3(1,:), ff3(2,:), ff3(3,:), ff3(4,:), ff3(5,:)];
   % Defuzzification
   ymout(k)=(w'*f3');         
   
   wof3=w'*f3';
   e(k) = rin(k)-ymout(k);
% Learning 
   d_w=0*w;
   for j=1:25
      d_w(j)=xite*e(k)*f3(1,j);
   end
   w=w+d_w;
   
   d_bi=0*bi;
  for i = 1:2 
   for j=1:5
      d_bi(i,j)=xite*e(k)*wof3*(bi(i,j)^-3)*norm(O12(i,j)-ci(i,j))^2;
   end
  end
   bi=bi+ d_bi;
   %pause
   
   
     for i=1:2
         for j=1:5
      d_ci(i,j)=xite*e(k)*wof3*(O12(i,j)-ci(i,j))*(bi(i,j)^-2);
         end
     end

   ci=ci+d_ci;
 
   d_wi=0*wi;
  for i = 1:2 
   for j=1:5
      d_wi(j)=xite*e(k)*wof3*(O12(i,j)-ci(i,j))*(bi(i,j)^-2)*O12_1(i,j);
   end
  end
   wi=wi+d_wi;
% RFNNI end ---------------------------------------------  

%%%%%%%%%%%%%%%%%%%%%%Jacobian%%%%%%%%%%%%%%%%%%%%%%%
yu=0;
yuu = 0;
for i = 1:2
  for  j=1:5
      yu=yu+wi(i,j)*(-O12(i,j)+ci(i,j))/bi(i,j)^2; 
  end
  yuu=yuu+yu;
end
  dyout(k)=yuu*wof3;

% RFNNC start --------------------------------------------- 
% Input layer
        for i1 = 1:2
           for j1=1:5
              O12C(i1,j1) = O12C_1(i1,j1)*wiC(i1,j1) + xC(i1,end);
           end
        end
   % Fuzzification
  for i = 1:2 
   for j=1:5
      f2C(i,j)=exp(-norm(O12C(i,j)-ciC(i,j))^2/(biC(i,j)^2));
   end
  end
   % Rules
 for j = 1:5
    m1C(j)=f2C(1,j);
    m2C(j)=f2C(2,j);
 end
  
   for i = 1:5
      for j = 1:5
         ff3C(i,j) = m2C(i)*m1C(j); 
      end
   end
  
   f3C = [ff3C(1,:), ff3C(2,:), ff3C(3,:), ff3C(4,:), ff3C(5,:)];
   % Defuzzification
   u(k)=(wC'*f3C');         
   
   wof3C=wC'*f3C';
   eC(k) = rin(k)-yout(k);
% Learning 
   d_wC=0*wC;
   for j=1:25
      d_wC(j)=xiteC*eC(k)*f3C(1,j);
   end
   wC=wC+d_wC;
   
   d_biC=0*biC;
  for i = 1:2 
   for j=1:5
      d_biC(i,j)=xiteC*eC(k)*wof3C*(biC(i,j)^-3)*norm(O12C(i,j)-ciC(i,j))^2;
   end
  end
   biC=biC+ d_biC*dyout(k);
   %pause
   
     for i=1:2
         for j=1:5
      d_ciC(i,j)=xiteC*eC(k)*wof3C*(O12C(i,j)-ciC(i,j))*(biC(i,j)^-2);
         end
     end
   ciC=ciC+d_ciC*dyout(k);
 
   d_wiC=0*wiC;
  for i = 1:2 
   for j=1:5
      d_wiC(j)=xiteC*eC(k)*wof3C*(O12C(i,j)-ciC(i,j))*(biC(i,j)^-2)*O12C_1(i,j);
   end
  end
   wiC=wiC+d_wiC*dyout(k);
   
% RFNNC end --------------------------------------------- 

   % Return params
   % RFNNI
   O12_1 = O12;
   x(1,k+1) = u(k);
   x(2,k+1) = y_1;
   u_1=u(k);
   y_1=yout(k);
   
   % RFNNC
   O12C_1 = O12C;
   xC(1,k+1) = eC(k);
   xC(2,k+1) = eC(k)-2*error_1-error_2;
   error_2 = error_1;   
   error_1 = eC(k);

   k
  % pause
end
figure(1);
plot(time,rin,'r',time,ymout,'g');
xlabel('time(s)');ylabel('rin,ymout');
legend('reference', 'NNout');

figure(2);
plot(time,rin,'r',time,yout,'b');
xlabel('time(s)');ylabel('rin,yout');
legend('reference', 'Plant');

figure(3);
plot(time,u,'r');
xlabel('time(s)');ylabel('u');
legend('control');

figure(4);
plot(time,dyout,'r');
xlabel('time(s)');ylabel('dyout');
legend('Jacobian');