电导增量法matlab,用s函数编写最大功率跟踪（电导增量法）

如何用S函数编写光伏电池的最大功率跟踪(电导增量法)仿真没有错误就是结果不对，求助各位大侠。

输入为V,I输出为Dref，3个状态变量。

function [sys,x0,str,ts] = s_MPPT2(t,x,u,flag)

%   The general form of an M-File S-function syntax is:

%       [SYS,X0,STR,TS] = SFUNC(T,X,U,FLAG,P1,...,Pn)

A=[0 0 0;0 0 0;0 0 1];

B=[1 0;0 1;0 0];

C=[0 0 1];

D=[0 0];

switch flag,

%%%%%%%%%%%%%%%%%%

% Initialization %

%%%%%%%%%%%%%%%%%%

case 0,

[sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D);

%%%%%%%%%%%%%%%

% Derivatives %

%%%%%%%%%%%%%%%

case 1,

sys=mdlDerivatives(t,x,u,A,B,C,D);

%%%%%%%%%%

% Update %

%%%%%%%%%%

case 2,

sys=mdlUpdate(t,x,u,A,B,C,D);

%%%%%%%%%%%

% Outputs %

%%%%%%%%%%%

case 3,

sys=mdlOutputs(t,x,u,A,B,C,D);

%%%%%%%%%%%%%%%%%%%%%%%

% GetTimeOfNextVarHit %

%%%%%%%%%%%%%%%%%%%%%%%

case 4,

sys=mdlGetTimeOfNextVarHit(t,x,u,A,B,C,D);

%%%%%%%%%%%%%

% Terminate %

%%%%%%%%%%%%%

case 9,

sys=mdlTerminate(t,x,u,A,B,C,D);

%%%%%%%%%%%%%%%%%%%%

% Unexpected flags %

%%%%%%%%%%%%%%%%%%%%

otherwise

error(['Unhandled flag = ',num2str(flag)]);

end

% end sfuntmpl

%

%=============================================================================

% mdlInitializeSizes

% Return the sizes, initial conditions, and sample times for the S-function.

%=============================================================================

%

function [sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D)

%

sizes = simsizes;

sizes.NumContStates  = 0;

sizes.NumDiscStates  = size(A,1);

sizes.NumOutputs     = size(D,1);

sizes.NumInputs      = size(D,2);

sizes.DirFeedthrough = 0;

sizes.NumSampleTimes = 1;   % at least one sample time is needed

sys = simsizes(sizes);

%

% initialize the initial conditions

%

x0  =[0 0 0.6];

%

% str is always an empty matrix

%

str = [];

%

% initialize the array of sample times

%

ts  = [1e-4 0];

% end mdlInitializeSizes

%=============================================================================

%

function sys=mdlDerivatives(t,x,u,A,B,C,D)

sys = [];

% end mdlDerivatives

%=============================================================================

%

function sys=mdlUpdate(t,x,u,A,B,C,D)

V=u(1);I=u(2);Vk=x(1);Ik=x(2);Dref=x(3);e=0.1;e0=0.005;De1=0.00025;De2=0.0001;De=0;

Ie=I-Ik;

Ve=V-Vk;

Pe=V*I-Vk*Ik;

if abs(Pe)

sys = A*x+B*u;

else

if abs(Pe)>e

De=De1;

else

De=De2;

end

end

if Ve==0

if Ie==0

Dref=Dref;

else if Ie>0

Dref=Dref+De;

else if Ie<0

Dref=Dref-De;

end

end

end

else if Ve~=0

if Ie/Ve==-I/V

Dref=Dref;

else if Ie/Ve>-I/V

Dref=Dref+De;

else if Ie/Ve

Dref=Dref-De;

end

end

end

end

end

x(3)=Dref;

sys = A*x+B*u;

% end mdlUpdate

%=============================================================================

%

function sys=mdlOutputs(t,x,u,A,B,C,D)

sys =x(3);

% end mdlOutputs

%=============================================================================

%

function sys=mdlGetTimeOfNextVarHit(t,x,u,A,B,C,D)

sys = [];

%=============================================================================

%

function sys=mdlTerminate(t,x,u,A,B,C,D)

sys = [];

% end mdlTerminate