Euler.mw

dsolve in Maple

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

ode:=diff(y(x),x)=ff(x,y(x));

ff := proc (x, y) options operator, arrow; -y end proc

ode := diff(y(x), x) = -y(x)

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

d1 := y(x) = exp(-x)

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

f1 := proc (x) options operator, arrow; exp(-x) end proc

> 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

[Plot]

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

display(p1,p0);

d2 := proc (x_rkf45) local res, data, vars, solnproc, outpoint, ndsol, i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; _EnvDSNumericSaveDigits := Digits; Digits := 15; if _...

[Plot]

Thermometer

Euler method

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

convert(%,polynom);

a := x0

series(f(x0)+D(f)(x0)*h+O(h^2),h,2)

f(x0)+D(f)(x0)*h

> 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);

x0 := 0

N := 10

h := 1/5

xx := array(sparse, 0 .. 10, [])

[0, 1/5, 2/5, 3/5, 4/5, 1, 6/5, 7/5, 8/5, 9/5, 2]

> 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:

yy := array(sparse, 0 .. 10, [])

yy[0] := 1

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

display(p0,p2);

[Plot]

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);

xx := array(sparse, 0 .. 10, [])

[0, 1/5, 2/5, 3/5, 4/5, 1, 6/5, 7/5, 8/5, 9/5, 2]

> 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:

yy := array(sparse, 0 .. 10, [])

yy[0] := 1

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

display(p0,p2,p3);

[Plot]

>

Logistics

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

func := proc (y, x) options operator, arrow; y-y^2 end proc

> 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);

ode := diff(y(x), x) = y(x)-y(x)^2

ds := proc (x_rkf45) local res, data, vars, solnproc, outpoint, ndsol, i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; _EnvDSNumericSaveDigits := Digits; Digits := 15; if _...

Warning, the name changecoords has been redefined

[Plot]

> 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:

g0 := .1

h := 2.4

N := 33.33333334

x := 0

g1 := .316

x := 2.4

total := [[0, .1], [2.4, .316]]

> plot(total);

[Plot]

> 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:

z := proc (y) options operator, arrow; h*y/(1+h) end proc

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

[Plot]

>