<Text-field style="Heading 2257" layout="Heading 2257">Euler method</Text-field>

We start with a Maple program for finding approximate solution of the following first-order differential equation LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYpLUkjbWlHRiQ2I1EhRictRiM2KUYrLUkmbWZyYWNHRiQ2KC1GIzYlLUYsNiVRI2R5RicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLyUrYmFja2dyb3VuZEdRLlsyNTUsMjU1LDI1NV1GJy9GPVEnbm9ybWFsRictRiM2JS1GLDYlUSNkeEYnRjlGPEY/RkIvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRk4vJSliZXZlbGxlZEdRJmZhbHNlRictSSNtb0dGJDYtUSI9RidGQi8lJmZlbmNlR0ZTLyUqc2VwYXJhdG9yR0ZTLyUpc3RyZXRjaHlHRlMvJSpzeW1tZXRyaWNHRlMvJShsYXJnZW9wR0ZTLyUubW92YWJsZWxpbWl0c0dGUy8lJ2FjY2VudEdGUy8lJ2xzcGFjZUdRLDAuMjc3Nzc3OGVtRicvJSdyc3BhY2VHRmJvLUYjNigtRiw2JVEiZkYnRjlGPC1GVTYtUTAmQXBwbHlGdW5jdGlvbjtGJ0ZCRlhGWkZmbkZobkZqbkZcb0Zeby9GYW9RJjAuMGVtRicvRmRvRl5wLUYjNiUtSShtZmVuY2VkR0YkNiQtRiM2J0YrLUYjNictRiw2JVEieEYnRjlGPC1GVTYtUSIsRidGQkZYL0ZlbkY7RmZuRmhuRmpuRlxvRl5vRl1wL0Zkb1EsMC4zMzMzMzMzZW1GJy1GLDYlUSJ5RidGOUY8Rj9GQkYrRj9GQkZCRj9GQkYrRj9GQkYrRj9GQkZccS1GIzYpLUZVNi1RIn5GJ0ZCRlhGWkZmbkZobkZqbkZcb0Zeb0ZdcEZfcEZicS1GY3A2JC1GIzYnRistRiM2JS1GLDYlUSJhRidGOUY8Rj9GQkYrRj9GQkZCRlQtRiw2JVEiQUYnRjlGPEY/RkJGK0Y/RkI=, on an interval [a,b] by using Euler method. Such program can be represented in Maple as a procedure which requires the knowledge of function f(x,y), numbers a, b and A and an integer number of points N for the number of steps (or the number of subintervals in a partition of the interval [a,b]).

We use command proc to define such procedure which we call eu .

eu:=proc(f,a,b,A,N) local n, xn, yn, h, k1, w; h:=evalf((b-a)/N); xn:=evalf(a); yn:=evalf(A); w[0]:=[xn,yn]; for n from 1 to N do k1:=evalf(f(xn,yn)); yn:=evalf(yn+h*k1); xn:=evalf(xn+h); w[n]:=[xn,yn]; end do; w; end:

COMMENT: the first line of this procedure tells Maple that our procedure should use function f, numbers a, b , A and an integer N. The second line tells that we use the following local variables: n to enumerate steps in Euler method, xn and yn for values of the xn and yn variables at the n-th step, h for stepsize, k1 for the falue of f(xn,yn), w for the array to keep results of our computation. On the third and fourth lines we compute a stepsize h, values xn and yn from the initial condition y(a)=A. Note that we use command evalf to tell Maple that we require numerical values (not symbolic ones). The fifth line starts with a loop command for ... do to repeat Euler method computations for every n-th step with n changing from 1 to N. Note that end do tells that the loop command is over. Before the last line command end (of procedure) we put the name of the array w to indicate that the output of our procedure is the array w which every n-th element w[n] consists of [xn,yn].

Now we define below a function f (x,y), observe how we can do it

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

Then we can call our procedure eu where a should be 0, b should be 1, A is 0 and a number of steps N, say, is 10 (then stepsize is (b-a)/N=0.1). The results of computations are kept in the next array v1

v1:=eu(f,0,1,0,10);

v1[5] ; The first number is the value of xn and the second number is yn for n=5

v1[5][2]; This is the yn value (we choose the second part of the element v1[5] of the array v1)

Below we print values of xn and yn in different forms (the third command line below contains also values of the exact solution y(xn) to compare with numerical approximations yn)

for n from 0 to 10 do lprint(v1[n][1],v1[n][2]) od;

for n from 0 to 10 do lprint(x||n=v1[n][1],y||n=v1[n][2]) od;

for n from 0 to 10 do lprint(x||n=v1[n][1],y||n=v1[n][2], y(x||n)=evalf(tan(v1[n][1]))) od;

<Text-field style="Heading 2" layout="Heading 2">Changing stepsize</Text-field>

Here we want to decrease the stepsize h in 10 times (which means we need to increase the number of steps N in 10 times) and compare new approximate solutions with the exact one. Please, pay attention how we print values yn corresponding to n which is multiple of 10 in for do loop command

v2:=eu(f,0,1,0,100);

for n from 0 to 100 by 10 do lprint(x||n=v2[n][1],y||n=v2[n][2], y(x||n)=evalf(tan(v2[n][1]))) od;

<Text-field style="Heading 2" layout="Heading 2">Visualizing approximate solutions</Text-field>

We can use use Maple graphing capabilities and visualize approximate solution using Maple package plots and

command seq . Command seq tells Maple to produce from the given array v1 or v2 an ordered sequences of points. Then we use command plot to prepare plotting data arrays p1, p2 from these sequences, and p3 for the known solution After that we use command display to plot simultaneously p1,p2, p3.

p1:=plot([seq(v1[k],k=0..10)],color=red,style=point,thickness=2): p2:=plot([seq(v2[k],k=0..100)],color=green,style=point,thickness=2): p3:=plot(tan(x),x=0..1,color=blue,thickness=2):

with(plots):

display(p1,p2,p3,view=[0..1,0..1.6]); The option view indicates to Maple the size of window for plot in the format [xmin . . xmax, ymin .. ymax].