SM212
Euler's method and direction fields

The examples below are taken from the board problems in class but they are not necessarily in order.

sm212_euler.mws, wdj, 9-2-98

with(plots):
with(DEtools):

eulers_method:=proc(f,x0,y0,h,n::integer)
local a,i,x,y;
x:=x0;
y:=y0;
a:=[[x,y]];
for i from 1 to n do
y:=y+h*f(x,y);
x:=x+h;
a:=[op(a),[x,y]];
od;
RETURN(a);
end;

[Maple Math]
[Maple Math]
[Maple Math]
[Maple Math]

improved_euler:=proc(f,x0,y0,h,n::integer)
local a,i,x,y;
x:=x0;
y:=y0;
a:=[[x,y]];
for i from 1 to n do
y:=y+(h/2)*(f(x,y)+f(x+h,y+h*f(x,y)));
x:=x+h;
a:=[op(a),[x,y]];
od;
RETURN(a);
end;

[Maple Math]
[Maple Math]
[Maple Math]
[Maple Math]
[Maple Math]
[Maple Math]
[Maple Math]
[Maple Math]

#5 y' = xy, y(0)=1

f0:=(x,y)->x*y;

[Maple Math]

em:=eulers_method(f0,0,1,1/3,3);
iem:=improved_euler(f0,0,1,1/3,3);

[Maple Math]

[Maple Math]

evalf(em); evalf(iem);

[Maple Math]

[Maple Math]

A0:=listplot(iem0,color=red):
A1:=listplot(iem1,color=blue):
display([A0,A1]);

[Maple Plot]

#4

DEplot(diff(y(x),x)=y(x)*x,y(x),
x=(-0.1)..(1.1),[[y(0)=1]],y=0..2,stepsize=.05);

[Maple Plot]

soln:=dsolve({diff(y(x),x)=y(x)*x,y(0)=1},y(x)):
y0:=t->subs(x=t,rhs(soln)):
y0(x);
evalf(y0(1));

[Maple Math]

[Maple Math]

#2 y' = y - x, y(0)=0

f0:=(x,y)->y-x;

[Maple Math]

em:=eulers_method(f0,0,0,1/3,3);
iem:=improved_euler(f0,0,0,1/3,3);

[Maple Math]

[Maple Math]

evalf(em);evalf(iem);

[Maple Math]

[Maple Math]

A0:=listplot(iem0,color=red):
A1:=listplot(iem1,color=blue):
display([A0,A1]);

[Maple Plot]

#1 y' = y-x, y(0)=0

DEplot(diff(y(x),x)=y(x)-x,y(x),
x=(-0.1)..(1.1),[[y(0)=0]],y=-2..2,stepsize=.05);

[Maple Plot]

soln:=dsolve({diff(y(x),x)=y(x)-x,y(0)=0},y(x)):
y0:=t->subs(x=t,rhs(soln)):
y0(x);
evalf(y0(1));

[Maple Math]

[Maple Math]

#7 y' = y^2 + x, y(0)=0

f0:=(x,y)->y^2+x;

[Maple Math]

em:=eulers_method(f0,0,0,1/2,2);
iem:=improved_euler(f0,0,0,1/2,2);

[Maple Math]

[Maple Math]

evalf(em);evalf(iem);

[Maple Math]

[Maple Math]

A0:=listplot(iem0,color=red):
A1:=listplot(iem1,color=blue):
display([A0,A1]);

[Maple Plot]

# 6

DEplot(diff(y(x),x)=y(x)^2+x,y(x),
x=(-0.1)..(1.1),[[y(0)=0]],y=0..1,stepsize=.01);

[Maple Plot]

#3 y' = -xy, y(0)=1, y(0)=-1, y(-1)=0

A1:=DEplot(diff(y(x),x)=-y(x)*x,y(x),
x=-2..2,[[y(0)=1]],y=-5..5,stepsize=.05,thickness=3):
A2:=DEplot(diff(y(x),x)=-y(x)*x,y(x),
x=-2..2,[[y(0)=-1]],y=-5..5,stepsize=.05,thickness=3):
A3:=DEplot(diff(y(x),x)=-y(x)*x,y(x),
x=-2..2,[[y(-1)=0]],y=-5..5,stepsize=.05,thickness=3):
display([A1,A2,A3]);

[Maple Plot]

DEplot(diff(y(x),x)=-y(x)*x,y(x),
x=(-0.1)..(1.1),[[y(0)=1]],y=0..2,stepsize=.05);

[Maple Plot]