S函数写的BP-PID控制程序

function [sys,x0,str,ts]=nnbp_pid(t,x,u,flag,T,nh,xite,alfa,kF1,kF2) %T是采样时间，nh是隐含层数量，xite是学习率，alfa是动量因子，kF1和kF2的值分别来选择隐含层和输出层的激活函数 %以上参数在双击main中BP_PID模块后可以设置 switch flag, case 0, [sys,x0,str,ts] = mdlInitializeSizes(T,nh); case 3, sys = mdlOutputs(t,x,u,T,nh,xite,alfa,kF1,kF2); case {1, 2, 4, 9},sys = []; otherwise, error(['Unhandled flag =',num2str(flag)]); end; % 初始化函数 function [sys,x0,str,ts] = mdlInitializeSizes(T,nh) sizes = simsizes; % 读入模板，得出默认的控制量 sizes.NumContStates = 0; sizes.NumDiscStates = 0; sizes.NumOutputs = 4+7*nh; %输出包括u,kp,ki,kd,wi(4*nh)个及wo(3*nh)个 sizes.NumInputs = 7+14*nh; %输入包括error(k),error(k-1),error(k-2),yout(k),yout(k-1)... %rin,u(k-1),wi_2(4*nh)个,wo_2(3*nh),wi_1(4*nh)个及wo_1(3*nh)个 sizes.DirFeedthrough = 1; sizes.NumSampleTimes = 1; sys = simsizes(sizes); x0 = []; str = []; ts = [T 0]; % 系统输出计算函数 function sys = mdlOutputs(t,x,u,T,nh,xite,alfa,kF1,kF2) wi_2=reshape(u(8:7+4*nh),nh,4); wo_2=reshape(u(8+4*nh:7+7*nh),3,nh); wi_1=reshape(u(8+7*nh:7+11*nh),nh,4); wo_1=reshape(u(8+11*nh:7+14*nh),3,nh); %权值更新，但有初始化问题（若不初始化，权值初值均为0，pid三个参数经过计算会一直保持为0.5） % if sum(wi_2)==0 %初始化 % wi_2=0.5*rands(nh,4); % wo_2=0.5*rands(3,nh); % wi_1=0.5*rands(nh,4); % wo_1=0.5*rands(3,nh); % end %以上注释的if---end部分为随机初始权值，效果不理想，采用以下代码对权值进行二次选择后作为初始权值，效果比较理想 %w_data.mat文件为权值随机初始化，输出稳定后确定的权值数据，可以load后作为网络初始权值，控制效果比较理想 if sum(wi_2)==0 %权值二次选择后初始化 load w_data; u=WW'; %WW里面前7个数据没有用，只是为了保持代码一致，随机给赋了7个值 wi_2=reshape(u(8:7+4*nh),nh,4); wo_2=reshape(u(8+4*nh:7+7*nh),3,nh); wi_1=reshape(u(8+7*nh:7+11*nh),nh,4); wo_1=reshape(u(8+11*nh:7+14*nh),3,nh); end xi=[u([6,4,1])', 1]; %xi=[rin,yout,error,1] xx=[u(1)-u(2);u(1);u(1)+u(3)-2*u(2)]; %增量式PID三个本量 I=xi*wi_1'; Oh=non_transfun(I,kF1); K=non_transfun(wo_1*Oh',kF2); uu=u(7)+K'*xx; dyu=sign((u(4)-u(5))/(uu-u(7)+0.0000001)); dK=non_transfun(K,3); %输出层激活函数的导数 delta3=u(1)*dyu*xx.*dK; wo=wo_1+xite*delta3*Oh+alfa*(wo_1-wo_2); dO=2*non_transfun(I,3); %隐含层激活函数的导数 wi=wi_1+xite*(dO.*(delta3'*wo))'*xi+alfa*(wi_1-wi_2); sys=[uu; K;wi(:);wo(:)]; %输出数量为4+7*nh % 激活函数及导数 function W1=non_transfun(W,key) switch key case 1, W1=(exp(W)-exp(-W))./(exp(W)+exp(-W)); %隐层激活函数 case 2, W1=exp(W)./(exp(W)+exp(-W)); %输出层激活函数 case 3, W1=2./(exp(W)+exp(-W)).^2; %输出层激活函数的导数 end