View source for Initial value problems

==Discrete Variable Methods==

Discrete variable methods approximate <math>y(t)</math> on a mesh
<math>t_0 < t_1 < \ldots < t_f\ .</math>
They start with the initial value <math>y_0 = A</math> and on
reaching <math>t_n\ ,</math> step to <math>t_{n+1} = t_n + h_n</math>
by computing <math>y_{n+1} \approx y(t_{n+1})\ .</math>
Methods are often analyzed with the assumption that the step sizes
<math>h_n</math> are constant, but general-purpose solvers
vary the step size for reasons discussed below. A great many
methods have been proposed, but three kinds dominate: Runge-Kutta,
Adams, BDFs (backward differentiation formulas).
On reaching <math>t_n</math>
there are available values <math>y_n,y_{n-1},\ldots </math> and
<math>F_n = F(t_n,y_n), F_{n-1},\ldots </math> that might be exploited.
Methods with memory like Adams and BDFs use some of these values; one-step methods like
Runge-Kutta evaluate <math>F</math> only at points in <math>[t_n,t_{n+1}]\ .</math>

===Euler Methods===

Any smooth function <math>u(t)</math> can be expanded in a Taylor series
:<math>
 u(t_{n+1}) = u(t_n) + h_n u't_n) + h_n^2
 \frac{u^{(2)}(t_n)}{2} + \ldots
</math>
If <math>u(t)</math> is the local solution defined by
:<math>\label{lsol}
u' = F(t,u),\quad u(t_n) = y_n
</math>

this expansion suggests an approximation of <math>u(t_{n+1})</math> called 
'''Euler's method''' or the '''forward Euler method'''
:<math>
 y_{n+1} = y_n + h_n F(t_n,y_n)
</math>
This is an explicit formula that requires one evaluation of <math>F</math>
per step. The series shows that the ''local error''
:<math>\label{le}
le_n = u(t_{n+1}) - y_{n+1}
</math>

is <math>O(h_n^2)\ .</math> When <math>u(t)</math> in this definition is the solution
<math>y(t)</math> of the initial value problem,
the local error is called the ''truncation error'' or ''discretization error''.

The '''backward Euler method''' is
:<math>
 y_{n+1} = y_n + h_n F(t_{n+1},y_{n+1})
</math>
It is not obvious why we might be interested in a formula that
defines <math>y_{n+1}</math> implicitly as the solution of a system of
algebraic equations and has the same accuracy as an explicit
formula, but it turns out that this formula is effective for stiff
systems and the forward Euler method is not. These two one-step
methods can be derived in several ways with extensions that lead to the popular
kinds of methods.

===Backward Differentiation Formulas===

The backward Euler method can be derived from a difference
approximation to the derivative in equation \eqref{odes}
:<math>
 \frac{y_{n+1} - y_n}{h_n} = F(t_{n+1},y_{n+1})
</math>
Difference approximations that use more of the values
<math>y_n,y_{n-1},\ldots </math> result in other 
'''backward differentiation formulas''' (BDFs). Interpolation
theory or Taylor series expansion shows that each additional 
<math>y_{n-j}</math> used raises the order of the truncation error by one.
The backward Euler formula is BDF1, the lowest order member of the
family. These formulas are the most popular for solving stiff
systems.

===Adams Methods===

'''Adams methods''' arise from integrating \eqref{lsol} to get
:<math>
 u(t_{n+1}) = y_n + \int^{t_{n+1}}_{t_n} F(t,u(t))\,dt
</math>
If the function <math>F(t,u(t))</math> here is approximated by a polynomial
interpolating values <math>F_n, F_{n-1},\ldots \ ,</math> an explicit
'''Adams-Bashforth (AB) formula''' is obtained. Euler's method is AB1,
the lowest order member of the family. If the polynomial also
interpolates <math>F(t_{n+1},y_{n+1})\ ,</math> an implicit '''Adams-Moulton (AM) formula''' is obtained. As with the BDFs, each additional 
<math>F_{n-j}</math> used raises the order of the truncation error by one.
The backward Euler method is AM1. Linear interpolation results in the
second order AM2,
:<math>
 y_{n+1} = y_n + \frac{h_n}{2}\, [F(t_n,y_n) + F(t_{n+1},y_{n+1})]
</math>
which is known as the '''trapezoidal rule'''.

How an implicit formula is evaluated is crucial to its use.
An explicit '''predictor''' formula is first used to compute a
tentative value <math>p_{n+1} \ .</math>  It is substituted for
<math>y_{n+1}</math> in the right hand side of the implicit
'''corrector''' formula and the formula evaluated to get a
better approximation to <math>y_{n+1}\ .</math>
For example, if the formulas AB1 and AM1 are used as a pair,
this process of ''simple iteration'' starts with the predicted value
:<math>
 y_{n+1}^{[0]} = y_n + h_n F(t_n, y_{n})
</math>
and successively corrects values
:<math>
 y_{n+1}^{[m+1]} = y_n + h_n F(t_{n+1}, y_{n+1}^{[m]})
</math>
This process is repeated until the implicit formula is satisfied well
enough or the step is deemed a failure.  If the step is a failure, <math>h_n</math> 
is reduced to improve the rate of convergence and the accuracy of the
prediction, and the program again tries to take a step.
Simple iteration is the standard way of evaluating implicit formulas for
non-stiff problems. If a predetermined, ''fixed'' number of iterations is done,
the resulting method is called a '''predictor-corrector''' formula.
If the predictor has a truncation error that is of the same order as
the implicit formula, a single iteration produces a result with a
truncation error that agrees to leading order with that of iterating
to completion. If the predictor is of order one lower, the result
has the same order as the corrector, but the leading terms in the
truncation error are different. There are two ways that Adams
formulas are implemented in popular solvers. Both predict with an
Adams-Bashforth formula and correct with an Adams-Moulton method.
One way is to iterate to completion, so that the integration is effectively done
with an implicit Adams-Moulton method. The other is to do a single
correction, which amounts to an explicit formula. For non-stiff
problems there is little difference in practice, but the two
implementations behave very differently when solving stiff problems.
As an example, if we predict with Euler's method (AB1) and correct
once with the trapezoidal rule (AM2), we obtain a formula
:<math>
 p_{n+1} = y_n + h_n F(t_n,y_n)
</math>
:<math>
 y_{n+1} = y_n + \frac{h_n}{2}\, [F(t_n,y_n) + F(t_{n+1},p_{n+1})]
</math>
that is known as ''Heun's'' method. It is an explicit method that
costs two evaluations of <math>F</math> per step and has a local error
that is <math>O(h_n^3)\ .</math>

===Runge-Kutta Methods===

Euler's method, the backward Euler method, the trapezoidal rule, and
Heun's method are all examples of '''Runge-Kutta (RK) methods'''.
RK methods are often derived by writing down the form of the
one-step method, expanding in Taylor series, and choosing
coefficients to match terms in a Taylor series expansion of the
local solution to as high order as possible. The general form of an explicit Runge-Kutta method involving two evaluations of <math>F\ ,</math> or as they are called in this context, ''stages'', is
:<math>
 y_{n,1} = y_n + \alpha_{1} h_n F(t_n,y_n)
</math>
:<math>
 y_{n+1} = y_n + h_n \, [ \beta_{2,0} F(t_n,y_n) + 
 \beta_{2,1} F(t_{n+1},y_{n,1})]
</math>
It turns out that it is not possible to get order three with just two stages, but the coefficients of Heun's method show that it is possible to get order two.  Another formula of order two that has been used in some solvers is
:<math>
 y_{n,1} = y_n + \frac{h_n}{2} F(t_n,y_n)
</math>
:<math>
 y_{n+1} = y_n + h_n\, F(t_{n+1},y_{n,1})
</math>
Explicit RK methods are very popular for solving non-stiff IVPs. Implicit RK methods are very popular for solving BVPs and also used for solving stiff IVPs.