## Mathematics & Science Learning Center Computer Laboratory

### Mathematica Implementation

(continued from last page...)

euler[f_,{x_,x0_,xn_},{y_,y0_},steps_]:=

Block[{ xold=x0,
yold=y0,
sollist={{x0,y0}},
x,y,h
},

h=N[(xn-x0)/steps];

Do[ xnew=xold+h;
ynew=yold+h*(f/.{x->xold,y->yold});
sollist=Append[sollist,{xnew,ynew}];
xold=xnew;
yold=ynew,
steps
];

Return[sollist]
]

The red-highlighted parts of the code are the only locations in the program that actually require any modification. If you sift through the rest of the code carefully you should be able to see that it has nothing to do with the algorithm being used, and is merely concerned with the nitty-gritty of constructing our list of solution points, etc., as discussed in great detail in the previous lab.

Obviously the name of the program should no longer be "euler," since we're now planning on using the Heun algorithm. A more appropriate name would be "heun," of course, so we'll make that modification.

More importantly, the method we're using to predict the next y-value is now much more sophisticated, namely:

yn+1 = yn + (h/2)   (f(xnyn) + f(xn + hyn +  h f(xnyn)))

We need to modify the line:

ynew=yold+h*(f/.{x->xold,y->yold});

accordingly, making appropriate changes to put the iteration formula we found for Heun's method into the language of Mathematica. This is something I'm going to leave you to do by yourself. (I warned you that by the end of the course I'd be expecting more independent thinking from you, and this is really not asking too much!) There is actually more than one way to modify the code and get the correct results from the program, some of which even utilize inserting more than one line of new code!