demonstration (follow along)
For[t = List[0]; y = List[1];
i = 1, i <= 10, i++,
t = Append[t, t[[i]] + 1];
y = Append[y, y[[i]] + (t[[i]] - y[[i]])]];
t
y
coords = Table[{t[[i]], y[[i]]}, {i, Length[y]}]
p1 = ListLinePlot[coords]
This is Euler's method with step size 1 for approximating y' = t-y
starting from the initial condition y(0) = 1, and going out to t=10.
We are constructing lists of values of t_i and y_i so that
t_i = t_{i-1} + Delta(t), y_i = y_{i-1} + Delta(y),
with Delta(t) = 1 and Delta(y) = y' * Delta(t).
Notes on syntax:
- the List command makes a list with some initial elements.
- the Append command takes in a list and an object, and adds the object to the end
of the list
- the i-th element of a list
y is accessed by y[[i]] , where i=1, not i=0, corresponds to the first element of the list
- the For command executes a loop. statements before the i=1 command are initial
statements before the loop, and statements after the i++ (which means
increment i by 1 each time) are statements done repeatedly. be careful about
the differentiation between commas and semicolons here, which is not
quite intuitive. (in our formatting, all the commas go on the second line,
semicolons everywhere else.)
- when h is small, the
coords list will be long
and you may want to hide the output. you can do this by putting a semicolon
at the end of that line.
- in the above case, we just have
t[[i]] = i-1 and
below we will have t[[i]] = h*(i-1) . however, i constructed
a table of t_i's in the above way make it more clear to generalize this to
systems. to show you want i mean, here is a simplification of the above code
using that t[[i]] = i-1 :
For[y = List[1];
i = 1, i <= 10, i++,
y = Append[y, y[[i]] + ((i-1) - y[[i]])]];
t
y
coords = Table[{(i-1), y[[i]]}, {i, Length[y]}]
p1 = ListLinePlot[coords]
you can copy and paste to check it gives you the same output.
- this code is not particularly efficient. there are two things
that are inefficient here: using loops in mathematica and the repeated use of
the
Append command. the latter issue is that mathematica makes a
new list (array) each time this is called. this is okay for our purposes,
and i did not try to make it more efficient so as not to make the code
unnecessarily confusing (also, i am no mathematica expert), but if you want,
you can try to make it more efficient (say for your project).
Now we will see how to do this for a general step size h=Delta(t).
For[t = List[0]; y = List[1]; h = 0.5;
i = 1, i*h <= 10, i++,
t = Append[t, t[[i]] + h];
y = Append[y, y[[i]] + (t[[i]] - y[[i]])*h]];
coords = Table[{t[[i]], y[[i]]}, {i, Length[y]}]
p2 = ListLinePlot[coords]
Show[p1, p2]
Show[p1, p2, PlotRange -> {{0, 3}, {0, 3}}]
do on your own
- Continuing with the above example, try various step sizes h. About
how small do you need to take h to be confident in your approximation
for y for 0 <= t <= 10.
for homework
-
The Euler method works similarly for systems of first-order equations,
and thus indirectly for higher-equations.
Recall the system: x' = -sin y, y' = x
from lab 3 for the pendulum y'' = -sin y.
Use this to plot approximations to y(t) with varying step sizes to get
good approximations to solutions for 0 <= t<= 10 for the following ICs
(a) y(0)=1, y'(0) = 0.
y(0)=0, y'(0) = 1.
y(0)= 4, y'(0) = 0.
y(0)= 1, y'(0) = 3.
(compare the solutions with the phase portrait to see differences in types
of solutions)
-
Revisit the predator-prey example from lab 3.
That is, x' = 2x - 0.1xy, y' = 0.005xy - 0.3y, where x(t) is the rabbit
population and y(t) is the fox population. Plot approximations to both
x(t) and y(t) for 0 <= t <= 50 for the following sets of initial conditions:
(a) 100 rabbits, 20 foxes
(b) 50 rabbits, 20 foxes
(c) 10 rabbits, 10 foxes
(d) 100 rabbits, 50 foxes
(e) Among scenarios (a)--(d), are there ones where
the model doesn't make sense in terms of
populations? Which scenarios seem to be robust, in the sense that small
deviations from the model due to random events will be "corrected" over time?
(E.g., if a population gets
down to 1, or is all males, it's probably not going to bounce back up.)
Explain.
help
see the lab 1 and
lab 2 pages for reminders on how to do some things
in mathematica and pointers to documentation.
i also found this book by
shifrin helpful for getting deeper into mathematica.
|