Euler's Method: Numerical and Graphical Solutions
The approximation method suggested above is known as Euler's method . It uses the tangent line to approximate successive points on a solution curve of a differential equation. Notice that in the implementation below exactly the same method is used to generate the sequence of points as was used for Newton's method. Namely, we define a Maple function Euler that implements one step of Euler's method, and we then iterate this function. This generates the desired list of points.
f:=x->24*x^3;
P0:=[0,-1];
deltax:=0.1;
Euler:=P->[P[1]+deltax,P[2]+deltax*f(P[1])];
Note that in defining the function Euler we gave it a point P = (
,
) as its independent variable. In Maple we represent a point as a list of two numbers, thus
= P[1] and
= P[2]. This explains the strange notation in the definition of the function Euler. Now we use Euler to compute a list of points on an approximate solution curve:
Euler(P0);
Euler(%);
Euler(%);
Euler(%);
Euler(%);
Euler(%);
Euler(%);
Euler(%);
Euler(%);
Euler(%);
The final touch is to generate the whole list of points by calling upon the function iterate , that you activated earlier in this worksheet, to apply Euler successively 10 times, starting with the initial point P0.
soln:=iterate(Euler,P0,10);
This is a list of points on the approximate solution curve. From this list we can generate either a table of values of the approximate solution or a plot of the approximate solution:
with(math3):
plot(soln,style=point);
plot(soln);
printtable(soln,`Numerical solution of y'=24*x^3`,[`x`,`y`]);
0.000000 -1.000000
.100000 -1.000000
.200000 -.997600
.300000 -.978400
.400000 -.913600
.500000 -.760000
.600000 -.460000
.700000 .058400
.800000 .881600
.900000 2.110400
1.000000 3.860000
Greater accuracy may be obtained by choosing smaller values of
. In Exercise 1 you will experiment with different value of
. Of course the full power of Euler's method is not yet apparent. The examples given above involved very simple differential equations that could be solved more easily by simply finding antiderivatives. But there is no reason why we need to restrict ourselves to such simple functions. EulerŐs method applies just as well to the initial-value problem
,
in which the slope function involves both variables
and
.
The next example carries out the very same Euler approximation in a more general setting. The first four lines just define the inputs to the problem: namely they define the function that is the right-hand side of the differential equation, specify the initial conditions by giving the point P0, give the end of the interval on which the plot is to be drawn, and specify how many points to plot. The remaining six lines are the ŇprogramÓ that carries out EulerŐs method. The idea is to use these ten lines as a program for solving many initial-value problems. Only the first four lines need be changed to consider different examples. You can cut and paste this program into many other settings. (For convenience this program is included later in the section entitled Program for Euler's Method )
f:=(x,y)->-x/y;
P0:=[0,1];
endOfInterval:=1;
numPts:=100;
Euler:=P->[P[1]+deltaI,P[2]+deltaI*f(P[1],P[2])]:
deltaI:=evalf((endOfInterval-P0[1])/numPts):
with(math3):
clear(i):
soln:=iterate(Euler,P0,numPts):
plot(soln);
As a final example let us revisit a previous CSC. In that problem the amount of water remaining in a cylindrical tank was modeled by the initial value problem
,
where the tank had a diameter of 20 feet (
), the hole in the bottom through which water was escaping had a diameter of 6 inches (
), and the initial depth of the water was 15 feet. (We took the value of
to be 32.16
.) In the earlier CSC you solved this differential equation explicitly and determined when the tank would be empty. Now we use Euler's method to approximate the water level in the tank and estimate the time when it will be empty. We simply copy the input cell containing the program for Euler's Method, above, and change the first four lines to the initial value problem at hand. The function
that defines the right side of the differential equation is (
,
)->-
, the initial depth of the water determines the initial point P0, and endOfInterval is the amount of time we wish to monitor the water level. We can experiment with the value of endOfInterval until it is large enough to observe the tank becoming empty. For example try values 500, 1000, 1200, 1400, 1500, 1550.
a:=Pi/16; A:=100*Pi; g:=32.16; H:=15;
f:=(t,h)->-a/A*sqrt(2*g*h);
P0:=[0,H];
endOfInterval:=1550;
numPts:=100;
Euler:=P->[P[1]+deltaI,P[2]+deltaI*f(P[1],P[2])]:
deltaI:=evalf((endOfInterval-P0[1])/numPts):
with(math3):
clear(i):
soln:=iterate(Euler,P0,numPts):
plot(soln);
Go back