220
8 Solving Ordinary Differential Equations
time point t n . The second idea is that we want to progress the solution from t n to
t n+1 as a straight line.
We know that the line must go through the solution at t n , i.e., the point (t n , u n ).
The differential equation tells us the slope of the line: u (t n ) = f (u n , t n ) = ru n .
That is, the differential equation gives a direct formula for the further direction of the
solution curve. We can say that the differential equation expresses how the system
(u) undergoes changes at a point.
There is a general formula for a straight line y = ax + b with slope a that goes
through the point (x 0 , y 0 ): y = a(x − x 0 ) + y 0 . Using this formula adapted to the
present case, and evaluating the formula for t n+1 , results in
u
n+1
= ru
n (t n+1 − t n ) + u
n
= u
n
+ Δt ru
n ,
which is nothing but the Forward Euler formula. You are now encouraged to do
Exercise 8.2 to become more familiar with the geometric interpretation of the
Forward Euler method.
8.2.5 Programming the FE Scheme; the General Case
Our previous program was just a simple main program without function definitions,
tailored to a special differential equation. When programming mathematics, it is
always good to consider a (large) class of problems and making a Python function
to solve any problem that fits into the class. More specifically, we will make software
for the class of differential equation problems of the form
u
(t) = f (u, t), u = U 0 , t ∈ [0, T ],
for some given function f , and numbers U 0 and T . We also take the opportunity
to illustrate what is commonly called a demo function. As the name implies,
the purpose of such a function is solely to demonstrate how the function works
(not to be confused with a test function, which does verification by use of
assert). The Python function calculating the solution must take f , U 0 , Δt, and
T as input, find the corresponding N t , compute the solution, and return an array
with u 0 , u 1 , . . . , u N t and an array with t 0 , t 1 , . . . , t N t . The Forward Euler scheme
reads
u
n+1
= u
n
+ Δtf (u
n , t n ), n = 0, . . . , N t − 1 .
The corresponding program may now take the form (file ode_FE.py):
import numpy as np
import matplotlib.pyplot as plt
def ode_FE(f, U_0, dt, T):
N_t = int(round(T/dt))
u = np.zeros(N_t+1)
t = np.linspace(0, N_t*dt, len(u))
u[0] = U_0
8 Solving Ordinary Differential Equations
time point t n . The second idea is that we want to progress the solution from t n to
t n+1 as a straight line.
We know that the line must go through the solution at t n , i.e., the point (t n , u n ).
The differential equation tells us the slope of the line: u (t n ) = f (u n , t n ) = ru n .
That is, the differential equation gives a direct formula for the further direction of the
solution curve. We can say that the differential equation expresses how the system
(u) undergoes changes at a point.
There is a general formula for a straight line y = ax + b with slope a that goes
through the point (x 0 , y 0 ): y = a(x − x 0 ) + y 0 . Using this formula adapted to the
present case, and evaluating the formula for t n+1 , results in
u
n+1
= ru
n (t n+1 − t n ) + u
n
= u
n
+ Δt ru
n ,
which is nothing but the Forward Euler formula. You are now encouraged to do
Exercise 8.2 to become more familiar with the geometric interpretation of the
Forward Euler method.
8.2.5 Programming the FE Scheme; the General Case
Our previous program was just a simple main program without function definitions,
tailored to a special differential equation. When programming mathematics, it is
always good to consider a (large) class of problems and making a Python function
to solve any problem that fits into the class. More specifically, we will make software
for the class of differential equation problems of the form
u
(t) = f (u, t), u = U 0 , t ∈ [0, T ],
for some given function f , and numbers U 0 and T . We also take the opportunity
to illustrate what is commonly called a demo function. As the name implies,
the purpose of such a function is solely to demonstrate how the function works
(not to be confused with a test function, which does verification by use of
assert). The Python function calculating the solution must take f , U 0 , Δt, and
T as input, find the corresponding N t , compute the solution, and return an array
with u 0 , u 1 , . . . , u N t and an array with t 0 , t 1 , . . . , t N t . The Forward Euler scheme
reads
u
n+1
= u
n
+ Δtf (u
n , t n ), n = 0, . . . , N t − 1 .
The corresponding program may now take the form (file ode_FE.py):
import numpy as np
import matplotlib.pyplot as plt
def ode_FE(f, U_0, dt, T):
N_t = int(round(T/dt))
u = np.zeros(N_t+1)
t = np.linspace(0, N_t*dt, len(u))
u[0] = U_0