254 8 Differential equations
meaning that we can define this optimal step length as
̃h=h
(
ξ
ε
) 1 /(M+ 1 )
.
Using this equation, we can design the following algorithm:
- If the two answers are close, use the current value for the step lengthh.
- Ifε>ξwe need to decrease the step size in the next time step.
- Ifε<ξwe need to increase the step size in the next time step.
At each step, two different approximations for the solutionare made and compared. If the
two answers are in close agreement, the approximation is accepted. If the two answers do
not agree to a specified accuracy, the step size is reduced. Ifthe answers agree to more
significant digits than required, the step size is increased. Even though this algorithm is
rather simple to implement, it requires unnecessarily manycomputations.
It is possible to reduce the number of operations by combining Runge-Kutta algorithms
of different orders. A much used algorithm is the so-called Runge-Kutta-Fehlberg algorithm
which uses a combination of fourth and fifth order Runge-Kutta methods, normally abbrevi-
ated to RKF45. Without going into much details, the philosophy of such methods consists in
evaluating the functionfsuch that the function values can be used for both the fourth order
and the fifth order method, avoiding thereby additional computations. The RKF45 method
requires at each step the computations of the following six values
k 1 =h f(tk,yk),
k 2 =h f(tk+
1
4
h,yk+
1
4
k 1 ),
k 3 =h f(tk+
3
8
h,yk+
3
32
k 1 +
9
32
k 2 ),
k 4 =h f(tk+
12
13 h,yk+
1932
2197 k^1 +
7200
2197 k^2 +
7296
2197 k^3 ),
k 5 =h f(tk+h,yk+
439
216
k 1 − 8 k 2 +
3680
513
k 3 +
845
4104
k 4 ),
and
k 6 =h f(tk+^1
2
h,yk−^8
27
k 1 + 2 k 2 −^3544
2565
k 2 +^1859
4104
k 4 −+^11
40
k 5 ).
Then an approximation to the solution of the ordinary differential equation is made using
a Runge-Kutta method of order four:
yk+ 1 =yk+
25
216
k 1 +
1408
2565
k 3 +
2197
4101
k 4 −
1
5
k 5 ,
where the four function valuesk 1 ,k 3 ,k 4 , andk 5 are used. Notice thatk 2 is not used here.
A better value for the solution is determined using a Runge-Kutta method of order five as
follows
zk+ 1 =yk+
16
135
k 1 +
6656
12825
k 3 +
28561
56430
k 4 −
9
50
k 5 +
2
55
k 6.
The optimal time stepαhis then determined by
α=
(
ξh
2 |zk+ 1 −yk+ 1 |