94 Ordinary Differential Equations
his procedure, two Runge-Kutta methods of different order are run simultane-
ously. The user specifies the desired upper bound for the local discretization
error. At each step, the algorithm produces an estimate Yn+I of the exact
solution. A second and more accurate method uses the same function eval-
uations as the first method plus one or two additional function evaluations
to produce a more accurate estimate Zn+I of the exact solution. The dif-
ference J Zn+i - Yn+i I yields an estimate of the local discretization error. If
the prescribed error bound is satisfied, the step is taken and a new stepsize
is estimated for the next step. If the prescribed error bound is not satis-
fied, the stepsize is reduced and another attempt is made to satisfy the error
bound using the smaller stepsize. In 1974, L. F. Shampine and H. A. Watts
implemented the Runge-Kutta-Fehlberg 4(5) method, RKF45. This method
requires six function evaluations per each successful step. Four function val-
ues are combined to produce a fourth order Runge-Kutta estimate Yn+I and
all six function values are combined to produce a fifth order Runge-Kutta
estimate Zn+I · The RKF45 numerical approximation method is available in
several commercial computer software packages.
!Comments on Computer Software I Historically, the numerical solution of
initial value problems originated in the 1700s with the work of Euler and Tay-
lor. In the late 1800s and early 1900 s, Runge, Heun, and Kutta developed
their numerical integration techniques. All of these techniques were used to
generate numerical solutions of initial value problems by hand- that is, with
paper and pencil. There were no electronic calculators or computers available
to compute the solutions. The numerical results which appear in examples
1, 3, 4, and 5 and Table 2.2 were generated using relatively simple computer
programs. If you know a scientific computing language, you should be able to
write computer programs to produce similar output. The software which ac-
companies this text does not include programs to produce a numerical solution
using the Taylor series method, Euler's method, the improved Euler's method,
the modified Euler's method, or the fourth order Runge-Kutta method. None
of these methods is capable of producing an extremely accurate numerical
solut ion over a "large" interval of the independent variable, because none of
these methods is an adaptive method and because most of these methods
are of low order. The software included with this text contains a program
named SOLVEIVP which uses a sophisticated, variable order, variable step-
size, numerical integration scheme to generate a numerical solution to an
initial value problem. Computer algebra systems (CAS) such as Mathemat-
ica, MATLAB, and MAPLE contain commands which allow you to generate a
numerical solution using the methods listed above. For example, the following
six MAPLE statements enable you to generate Euler's approximation to the
IVP (7) y' = y + x; y(O) = 1 on the interval [O, l] with a constant stepsize
of h = .2. When rounded to six significant digits, the values output by this
program are the same values which appear in the columns labeled Xn and Yn
of the table of values in example 3.