An earlier post discusses Mathematica code for Euler’s method. Here is some updated code.
A couple preliminary notes:
- I like to clear all variables at the beginning.
- For systems it is fun to superimpose the plot on top of the stream plot of the corresponding vector field
Example 1
For the IVP ,
with
we can use the code
Clear["Global`*"] f[t_, y_] := t - y deltat = 0.1; eulerValues = RecurrenceTable[{t[k + 1] == t[k] + deltat, y[k + 1] == y[k] + deltat*f[t[k],y[k]], t[0] == 0., y[0] == 2.}, {t, y}, {k, 0, 10}]; Grid[eulerValues] eulerPlot=ListLinePlot[eulerValues, PlotMarkers -> Automatic] realPlot = Plot[3 Exp[-t] + t - 1, {t, 0, 1}, PlotStyle -> Red]; Show[eulerPlot, realPlot]
Example 2
Consider ,
,
We use
Clear["Global`*"] deltat = .05; f[y_] := Exp[y]; values = RecurrenceTable[{ t[k + 1] == t[k] + deltat, y[k + 1] == y[k] + deltat*f[y[k]], t[0] == 0, y[0] == 1 }, {t, y}, {k, 0, 10}]; Grid[values] ListLinePlot[values, PlotMarkers -> Automatic]
Example 3
Here we consider the system
with and initial conditions
. The following code processes the first 100 time steps.
Clear["Global`*"] f[x_, y_] := x (1 - x) - x*y g[x_, y_] := -y + 2 x*y streamPlot = StreamPlot[ {f[x, y], g[x, y]}, {x, -1, 2}, {y, -1, 2}] deltat = 0.1; steps = 100; eulerValues = RecurrenceTable[{ t[0] == 0, x[0] == 1.5, y[0] == 1, t[k + 1] == t[k] + deltat, x[k + 1] == x[k] + deltat*f[x[k], y[k]], y[k + 1] == y[k] + deltat*g[x[k], y[k]]}, {t, x, y}, {k, 0, steps}]; eulerPlotValues = Table[{ eulerValues[[k, 2]], eulerValues[[k, 3]]}, {k, 1, steps + 1}]; eulerPlot = ListLinePlot[ eulerPlotValues, PlotMarkers -> Automatic, PlotStyle -> Red] Show[streamPlot, eulerPlot]