Computational Physics - Department of Physics

(Axel Boer) #1

8.5 Adaptive Runge-Kutta and multistep methods 253


8.5 Adaptive Runge-Kutta and multistep methods


In case the function to integrate varies slowly or fast in different integration domains, adap-
tive methods are normally used. One strategy is always to decrease the step size. As we have
seen earlier, this leads to more computations and may eventually even lead to the loss of
numerical precision. An alternative is to use higher-orderRunge-Kutta methods for example.
However, this leads again to more cycles, furthermore, there is no guarantee that higher-
order leads to an improved error, see for example the discussions in Ref. [42]
Assume the exact result isy ̃and that we are using a Runge-Kutta method of orderM.
Suppose we run two calculations, one with a step lengthh(which we will labely 1 ) and one
with step lengthh/ 2 (labelledy 2 ). The exact solution in terms ofy 1 is


y ̃=y 1 +ChM+^1 +O(hM+^2 ),

whereCis some constant and


y ̃=y 2 + 2 C(h/ 2 )M+^1 +O(hM+^2 ).

Note that we need to perform two calculations in the last equation, one for each interval
defined byh/ 2 .calculate two halves in the last equation. The difference between the two
solutions is then
|y 1 −y 2 |=ChM+^1 ( 1 −


1

2 M),

from which we can define the constantCas


C= |y^1 −y^2 |
( 1 − 2 −M)hM+^1

. (8.14)

We rewrite then the exact solution in terms of a quantityε


y ̃=y 2 +ε+O((h)M+^2 ),

with


ε=
|y 1 −y 2 |
2 M− 1.
If we employ our fourth-order Runge-Kutta scheme, we have


y ̃=y 2 +ε+O(h^6 ),

with
ε=|y^1 −y^2 |
15


.

The estimate is one order higher than the original Runge-Kutta method to fourth order. But
this method is normally rather inefficient since it requiresa lot of computations. We solve
typically the equation three times at each time step. However, we can compare the estimate
εwith some by us given accuracyξsay for exampleξ= 10 −^8. We can then ask the following
question: what is, with a givenyjandtj, the largest possible step sizeh ̃that leads to an error
belowξ? We want
Ch ̃M+^1 ≤ξ,


which leads to, using Eq. (8.14),


( ̃
h
h

)M+ 1

|y 1 −y 2 |
( 1 − 2 −M)
≤ξ,
Free download pdf