Method # 1 - Euler (a custom-defined procedure for our uses)

DO NOT CHANGE ANYTHING IN THE FOLLOWING INPUT REGION - YOU HAVE BEEN WARNED !!!

> Euler := proc( MODEL, IC, DOMAIN, N )

> local a, b, dt, dy, f, i, tt, yy, LISTpts, vars;

> vars := (op(2,lhs(MODEL) ), op([1,0],lhs(MODEL) ) );

> f := unapply( rhs(MODEL), vars );

> a := op(1,rhs(DOMAIN)); b := op(2,rhs(DOMAIN));

> dt := (b-a)/N;

> tt := a;

> yy := rhs(IC);

> LISTpts := [ tt, yy ];

> for i from 1 to N do

> dy := evalf( f(tt,yy) ) * dt;

> yy := yy + dy;

> tt := tt + dt;

> LISTpts := LISTpts, [ tt, yy ];

> od; RETURN( [ LISTpts ] );

> end:

DO NOT CHANGE ANYTHING IN THE PRECEDING INPUT REGION - YOU HAVE BEEN WARNED !!!

Let's recall what we are working with.

> MODEL;

> IC;

> VAR;

> DOMAIN;

The number of steps to be used in Euler's method is

> Nsteps := ?;

>

The points generated by Euler's method are

> ptsE := Euler( MODEL, IC, DOMAIN, Nsteps );

>

These data points can be plotted using the command

> plotE := plot( ptsE, color = GREEN, thickness = 2 ):

> plotE;

>

Another way to view the approximate solution is to overlay the computed solution on the slope field..

> display( [ plotSLOPE, plotE ] );

>