From: Torsten Hennig on
> Here is react.m Hope you can help
> function dydt1 = react (t,y)
> y1=6;
> k=2;
> L=1;
> z=0:0.00011999:1;
> t=z/L;

You can not modify t here ; it's the actual time
of integration.

> c=0:0.0003:2.5;
> a=4*k;
> s=y1/a;
> m=1./(1+(s.*c).^2).^0.5;
> cons=((2.*k)./m);
> [sn0,cn0,dn0]=ellipj(cons,m);
> t(8335)=[];

??? t is the time of integration, not an array.

> A1=(c.*((dn0.*(t-1)+1)./(2.*dn0.*(t-1)))).^0.5;
> A2=(c.*((1-dn0.*(t-1))./(2.*dn0.*(t-1)))).^0.5;
> psi=0:0.0003:2.5;
> % dydt1 = zeros(size(y));
> % disp (size(dydt1))
> % size(A1)
> % size(A2)
>
> C=(A2./A1)+(A1./A2);
> R=(A1.^2+A2.^2);
> % A1=y(1);
> % A2=y(2);
> % psi=y(3);
> % dydt1(1)=k*A2.*sin(psi);
> % dydt1(2)=k*A1.*sin(psi);
> % dydt1(3)=k.*cos(psi).*C+y1.*R;
>
> dydt1=[(k*A2.*sin(psi)).';(k*A1.*sin(psi)).';(k.*cos(p
> si).*C+y1.*R).'];

You have to form a vector of size 3 here ;
this is not consistent with psi which seems to be an
array.

> disp(dydt1)

Best wishes
Torsten.