[n,t,y,yp]=ddf(n,t,y,yp);

function [n,t,y,yp]=ddf(n,t,y,yp);
persistent alfa ; 

if isempty(alfa), alfa=0; end;
y_shape=size(y);y=reshape(y,1,[]);
yp_shape=size(yp);yp=reshape(yp,1,[]);
alfa = y(n+1);
yp(1) = 1.0d0 + alfa.*(y(2)-y(1)) - y(1).*y(3);
yp(2) = alfa.*(y(1)-y(2)) - y(2).*y(3);
yp(3) = 1.0d0 - y(3).*(y(1)+y(2));
y_shape=zeros(y_shape);y_shape(:)=y(1:numel(y_shape));y=y_shape;
yp_shape=zeros(yp_shape);yp_shape(:)=yp(1:numel(yp_shape));yp=yp_shape;
end %subroutine ddf