C7MW.html

Euler.mw

dsolve in Maple

> restart;
ff:=(x,y)->-y;
ode:=diff(y(x),x)=ff(x,y(x));

> d1:=dsolve({ode,y(0)=1},y(x));

> f1:=unapply(rhs(d1),x);

> with(plots):with(linalg):
p0:=plot(f1(x),x=0..2):
display(p0);

Warning, the name changecoords has been redefined

Warning, the protected names norm and trace have been redefined and unprotected

> d2:=dsolve({ode,y(0)=1},type=numeric,range=0..2);
p1:=odeplot(d2,color=black):
display(p1,p0);

Thermometer

Euler method

> a:=x0;
subs(x-a=h,series(f(x),x=a,2));
convert(%,polynom);

> x0:=0;
N:=10;
h:=2/N;
xx:=array(0..N,sparse);
for i from 0 to N do
xx[i]:=i*h;
end do:
convert(xx,list);

> yy:=array(0..N,sparse);yy[0]:=1;
for i from 0 to N-1 do
yy[i+1]:=yy[i]+ff(xx[i],yy[i])*h;
end do:

> v1:=transpose([convert(xx,list),convert(yy,list)]):
p2:=pointplot(v1):
display(p0,p2);

Predictor-Corrector method

Implicit Euler method

Runge-Kutta

> #x0:=0;
#N:=10;
#h:=2/N;
xx:=array(0..N,sparse);
for i from 0 to N do
xx[i]:=i*h;
end do:
convert(xx,list);

> yy:=array(0..N,sparse);yy[0]:=1;
for i from 0 to N-1 do
k1:=h*ff(xx[i],yy[i]);
k2:=h*ff(xx[i]+h/2,yy[i]+k1/2);
k3:=h*ff(xx[i]+h/2,yy[i]+k2/2);
k4:=h*ff(xx[i]+h,yy[i]+k3);
yy[i+1]:=yy[i]+(k1+2*k2+2*k3+k4)/6;
end do:

> v1:=transpose([convert(xx,list),convert(yy,list)]):
p3:=pointplot(v1):
display(p0,p2,p3);

Logistics

> restart;
func:=(y,x)->y-y^2;

> ode:=diff(y(x),x)=func(y(x),x);
ds:=dsolve({ode,y(0)=0.1},type=numeric,range=0..80);
with(plots):
odeplot(ds);

Warning, the name changecoords has been redefined

> g0:=0.1;h:=2.4;N:=80/h;x:=0;
g1:=g0+func(g0,x)*h;
x:=x+h;
total:=[[0,g0],[x,g1]];
for k from 2 to N do
g2:=g1+func(g1,x)*h;
x:=x+h;
total:=[op(total),[x,g2]];
g1:=g2;
end do:

> plot(total);

> z:=y->h/(1+h)*y;
g1:=g0+func(g0,x)*h:
x:=x+h:
logistic:=[[z(0.1),0],[z(0.1),z(g1)]]:
for k from 2 to N do
g2:=g1+func(g1,x)*h;
x:=x+h;
logistic:=[op(logistic),[z(g1),z(g1)],[z(g1),z(g2)]];
g1:=g2;
end do:

> plot([logistic,z,(1+h)*z*(1-z)],z=0..1);