Numerical ordinary differential equations is the part of numerical analysis which studies the numerical solution of ordinary differential equations (ODEs). This field is also known under the name numerical integration, but some people reserve this term for the computation of integrals. Numerical analysis is the study of approximate methods for the problems of continuous mathematics (as distinguished from discrete mathematics). ... A simulation of airflow into a duct using the Navier-Stokes equations A differential equation is a mathematical equation for an unknown function of one or several variables which relates the values of the function itself and of its derivatives of various orders. ... Numerical Integration with the Monte Carlo method: Nodes are random equally distributed. ... In calculus, the integral of a function is an extension of the concept of a sum. ...

Many differential equations cannot be solved analytically, in which case we have to satisfy ourselves with an approximation to the solution. The algorithms studied here can be used to compute such an approximation. An alternative method is to use techniques from calculus to obtain a series expansion of the solution. In mathematics, computing, linguistics, and related disciplines, an algorithm is a finite set of well-defined instructions for accomplishing some task which, given an initial state, will terminate in a defined end-state. ... Calculus (from Latin, pebble or little stone) is a major area in mathematics where infinitesimal data yields global information. ... In mathematics, a series is often represented as the sum of a sequence of terms. ...

Ordinary differential equations occur in many scientific disciplines, for instance in mechanics, chemistry, biology, and economics. In addition, some methods in numerical partial differential equations convert the partial differential equation into an ordinary differential equation, which must then be solved. Mechanics (Greek ) is the branch of physics concerned with the behaviour of physical bodies when subjected to forces or displacements, and the subsequent effect of the bodies on their environment. ... Chemistry - the study of atoms, made of nuclei (conglomeration of center particles) and electrons (outer particles), and the structures they form. ... This article or section does not cite any references or sources. ... Face-to-face trading interactions on the New York Stock Exchange trading floor. ... Numerical partial differential equations is the branch of numerical analysis that studies the numerical solution of partial differential equations. ... In mathematics, a partial differential equation (PDE) is a relation involving an unknown function of several independent variables and its partial derivatives with respect to those variables. ...

The problem

We want to approximate the solution of the differential equation

where f is a function that maps [t₀,∞) × R^d to R^d, and the initial condition y₀ ∈ R^d is a given vector.

The above formulation is called an initial value problem (IVP). The Picard-Lindelöf theorem states that there is a unique solution, if f is Lipschitz continuous. In contrast, boundary value problems (BVPs) specify (components of) the solution y at more than one point. Different methods need to be used to solve BVPs, for example the shooting method, multiple shooting or global methods like finite differences or collocation methods. In mathematics, an initial value problem is a statement of a differential equation together with specified value of the unknown function at a given point in the domain of the solution. ... In mathematics, the Picard-Lindelöf theorem on existence and uniqueness of solutions of differential equations (Picard 1890, Lindelöf 1894) states that an initial value problem has exactly one solution if f is Lipschitz continuous in , continuous in as long as stays bounded. ... In mathematics, a function f : M → N between metric spaces M and N is called Lipschitz continuous (or is said to satisfy a Lipschitz condition) if there exists a constant K > 0 such that d(f(x), f(y)) ≤ K d(x, y) for all x and y... Shows a region where a differential equation is valid and the associated boundary values In mathematics, in the field of differential equations, a boundary value problem is a differential equation together with a set of additional restraints, called the boundary conditions. ... In numerical analysis, the shooting method is a method for solving a boundary value problem by reducing it to the solution of an initial value problem. ... A finite difference is a mathematical expression of the form f(x + b) − f(x +a). ... In mathematics, a collocation method is a method for the numerical solution of ordinary differential equation and partial differential equations and integral equations. ...

Note that we restrict ourselves to first-order differential equations (meaning that only the first derivative of y appears in the equation, and no higher derivatives). However, a higher-order equation can easily be converted to a first-order equation by introducing extra variables. For example, the second-order equation y'' = −y can be rewritten as two first-order equations: y' = z and z' = −y.

Methods

Two elementary methods are discussed to give the reader a feeling for the subject. After that, pointers are provided to other methods (which are generally more accurate and efficient). The methods mentioned here are analysed in the next section.

The Euler method

Starting with the differential equation (1), we replace the derivative y' by the finite difference approximation In mathematics and computational science, Euler integration is the most basic kind of numerical integration for calculating trajectories from forces at discrete timesteps. ... A finite difference is a mathematical expression of the form f(x + b) − f(x +a). ...

which yields the following formula

This formula is usually applied in the following way. We choose a step size h, and we construct the sequence t₀, t₁ = t₀ + h, t₂ = t₀ + 2h, … We denote by y_n a numerical estimate of the exact solution y(t_n). Motivated by (3), we compute these estimates by the following recursive scheme A visual form of recursion known as the Droste effect. ...

This is the Euler method, named after Leonhard Euler who described this method in 1768. Leonhard Euler (pronounced Oiler; IPA ) (April 15, 1707 – September 18 [O.S. September 7] 1783) was a pioneering Swiss mathematician and physicist, who spent most of his life in Russia and Germany. ... 1768 was a leap year starting on Friday (see link for calendar). ...

Euler's method is an example of an explicit method. This means that the new value y_n+1 is defined in terms of things that are already known, like y_n. In applied mathematics, explicit and implicit methods are approaches used in computer simulations of physical processes, or in other words, they are numerical methods for solving time-variable ordinary and partial differential equations. ...

The backward Euler method

If, instead of (2), we use the approximation

we get the backward Euler method:

The backward Euler method is an implicit method, meaning that we have to solve an equation to find y_n+1. One often uses fixed point iteration or (some modification of) the Newton-Raphson method to achieve this. Of course, it costs time to solve this equation; this cost must be taken into consideration when one selects the method to use. The advantage of implicit methods such as (6) is that they are usually more stable for solving a stiff equation, meaning that a larger step size h can be used. In applied mathematics, explicit and implicit methods are approaches used in computer simulations of physical processes, or in other words, they are numerical methods for solving time-variable ordinary and partial differential equations. ... In numerical analysis, fixed point iteration is a method of computing fixed points of functions. ... In numerical analysis, Newtons method (also known as the Newton–Raphson method or the Newton–Fourier method) is an efficient algorithm for finding approximations to the zeros (or roots) of a real-valued function. ... In mathematics, a stiff equation is a differential equation for which certain numerical methods for solving the equation are numerically unstable, unless the step size is taken to be extremely small. ...

Generalizations

The Euler method is often not accurate enough. In more precise terms, it only has order one (the concept of order is explained below). This caused mathematicians to look for higher-order methods.

One possibility is to use not only the previously computed value y_n to determine y_n+1, but to make the solution depend on more past values. This yields a so-called multistep method. Almost all practical multistep methods fall within the family of linear multistep methods, which have the form Multistep methods are methods used in numerical analysis to solve numerically ordinary differential equations. ...

Another possibility is to use more points in the interval [t_n,t_n+1]. This leads to the family of Runge-Kutta methods, named after Carle Runge and Martin Kutta. One of their fourth-order methods is especially popular. In numerical analysis, the Runge-Kutta methods are a family of techniques for the approximation of solutions of ordinary differential equations. ... Carle David Tolmé Runge (August 30, 1856 – January 3, 1927) was a German mathematician, physicist, and spectroscopist. ... Martin Wilhelm Kutta (November 3, 1867 - December 25, 1944) was a German mathematician. ...

Both ideas can also be combined. The resulting methods are called general linear methods.

Advanced features

A good implementation of one of these methods for solving an ODE entails more than the time-stepping formula.

It is often inefficient to use the same step size all the time, so variable step-size methods have been developed. Usually, the step size is chosen such that the (local) error per step is below some tolerance level. This means that the methods must also compute an error indicator, an estimate of the local error.

An extension of this idea is to choose dynamically between different methods of different orders (this is called a variable order method). Extrapolation methods are often used to construct various methods of different orders. In mathematics, extrapolation is the process of constructing new data points outside a discrete set of known data points. ...

Other desirable features include:

• dense output: cheap numerical approximations for the whole integration interval, and not only at the points t₀, t₁, t₂, ...
• event location: finding the times where, say, a particular function vanishes.
• support for parallel computing.

Parallel computing is the simultaneous execution of the same task (split up and specially adapted) on multiple processors in order to obtain results faster. ...

Alternative methods

Many methods do not fall within the framework discussed here. Some classes of alternative methods are:

• multiderivative methods, which use not only the function f but also its derivatives. This class includes Hermite-Obreschkoff methods and Fehlberg methods, as well as methods like the Parker–Sochacki method, which compute the coefficients of the Taylor series of the solution y recursively.
• methods for second order ODEs. We said that all higher-order ODEs can be transformed to first-order ODEs of the form (1). While this is certainly true, it may not be the best way to proceed. In particular, Nyström methods work directly with second-order equations.
• geometric integration methods are especially designed for special classes of ODEs (e.g., symplectic integrators for the solution of Hamiltonian equations). They take care that the numerical solution respects the underlying structure or geometry of these classes.

In mathematics, the Runge–Kutta–Fehlberg method (or Fehlberg method) is a method for the numerical solution of ordinary differential equations developed by the German mathematician Erwin Fehlberg. ... The Parker-Sochacki method is a new algorithm for solving systems of differential equations, which has been developed by J. Edgar Parker and James Sochacki, of the James Madison University Mathematics Department. ... As the degree of the Taylor series rises, it approaches the correct function. ... In the mathematical field of numerical ODEs, a geometric integrator is a numerical method that preserves geometric properties of the exact flow of a differential equation. ... In mathematics, a symplectic integrator (SI) is a numerical integration scheme for a specific group of differential equations related to classical mechanics. ... Hamiltonian mechanics is a re-formulation of classical mechanics that was invented in 1833 by William Rowan Hamilton. ...

Analysis

Numerical analysis is not only the design of numerical methods, but also their analysis. Three central concepts in this analysis are: Numerical analysis is the study of approximate methods for the problems of continuous mathematics (as distinguished from discrete mathematics). ...

• convergence: whether the method approximates the solution,
• order: how well it approximates the solution, and
• stability: whether errors are damped out.

In the mathematical subfield of numerical analysis, numerical stability is a property of numerical algorithms. ...

Convergence

A numerical method is said to be convergent if the numerical solution approaches the exact solution as the step size h goes to 0. More precisely, we require that for every ODE (1) with a Lipschitz function f and every t^* > 0, In mathematics, a function f : M → N between metric spaces M and N is called Lipschitz continuous (or is said to satisfy a Lipschitz condition) if there exists a constant K > 0 such that d(f(x), f(y)) ≤ K d(x, y) for all x and y...

All the methods mentioned above are convergent. In fact, convergence is a condition sine qua non for any numerical scheme. Sine qua non or condicio sine qua non was originally a Latin legal term for without which it could not be (but for). It refers to an indispensable and essential action, condition, or ingredient. ...

Consistency and order

Suppose the numerical method is

The local error of the method is the error committed by one step of the method. That is, it is the difference between the result given by the method, assuming that no error was made in earlier steps, and the exact solution:

The method is said to be consistent if

The method has order p if

Hence a method is consistent if it has an order greater than 0. The (forward) Euler method (4) and the backward Euler method (6) introduced above both have order 1, so they are consistent. Most methods being used in practice attain higher order. Consistency is a necessary condition for convergence, but not sufficient; for a method to be convergent, it must be both consistent and zero-stable.

A related concept is the global error, the error sustained in all the steps one needs to reach a fixed time t. Explicitly, the global error at time t is y_N − y(t) where N = (t−t₀)/h. The global error of a pth order one-step method is O(h^p); in particular, such a method is convergent. This statement is not necessarily true for multi-step methods.

Stability and stiffness

Main article: Stiff equation

For some differential equations, application of standard methods—such as the Euler method, explicit Runge-Kutta methods, or multistep methods (e.g., Adams-Bashforth methods)—exhibit instability in the solutions, though other methods may produce stable solutions. This "difficult behaviour" in the equation (which may not necessarily be complex itself) is described as stiffness, and is often caused by the presence of different time scales in the underlying problem. Stiff problems are ubiquitous in chemical kinetics, control theory, solid mechanics, weather prediction, biology, and electronics. In mathematics, a stiff equation is a differential equation for which certain numerical methods for solving the equation are numerically unstable, unless the step size is taken to be extremely small. ... In numerical analysis, the Runge-Kutta methods are an important family of implicit and explicit iterative methods for the approximation of solutions of ordinary differential equations. ... Multistep methods are methods used in numerical analysis to solve numerically ordinary differential equations. ... In physical chemistry, chemical kinetics or reaction kinetics study reaction rates in a chemical reaction. ... In engineering and mathematics, control theory deals with the behavior of dynamical systems. ... Solid mechanics is the branch of physics and mathematics that concern the behavior of solid matter under external actions (e. ... Weather is a term that encompasses phenomena in the atmosphere of a planet. ... This article or section does not cite any references or sources. ... Electronics is the study of the flow of charge through various materials and devices such as, semiconductors, resistors, inductors, capacitors, nano-structures, and vacuum tubes. ...

History

Below is a timeline of some important developments in this field. For the novel by Michael Crichton, see Timeline (novel). ...

• 1768 - Leonhard Euler publishes his method.
• 1824 - Augustin Louis Cauchy proves convergence of the Euler method. In this proof, Cauchy uses the implicit Euler method.
• 1855 - First mention of the multistep methods of John Couch Adams in a letter written by F. Bashforth.
• 1895 - Carle Runge publishes the first Runge-Kutta method.
• 1905 - Martin Kutta describes the popular fourth-order Runge-Kutta method.
• 1910 - Lewis Fry Richardson announces his extrapolation method.
• 1952 - Charles F. Curtiss and Joseph Oakland Hirschfelder coin the term stiff equations.

1768 was a leap year starting on Friday (see link for calendar). ... Leonhard Euler (pronounced Oiler; IPA ) (April 15, 1707 – September 18 [O.S. September 7] 1783) was a pioneering Swiss mathematician and physicist, who spent most of his life in Russia and Germany. ... 1824 was a leap year starting on Thursday (see link for calendar). ... Augustin Louis Cauchy (August 21, 1789 – May 23, 1857) was a French mathematician. ... 1855 was a common year starting on Monday (see link for calendar). ... Multistep methods are methods used in numerical analysis to solve numerically ordinary differential equations. ... John Couch Adams (June 5, 1819 – January 21, 1892), was a British mathematician and astronomer. ... 1895 (MDCCCXCV) was a common year starting on Tuesday (see link for calendar) of the Gregorian calendar (or a common year starting on Thursday of the 12-day-slower Julian calendar). ... Carle David Tolmé Runge (August 30, 1856 – January 3, 1927) was a German mathematician, physicist, and spectroscopist. ... In numerical analysis, the Runge-Kutta methods are a family of techniques for the approximation of solutions of ordinary differential equations. ... 1905 (MCMV) was a common year starting on Sunday (link will display the full calendar). ... Martin Wilhelm Kutta (November 3, 1867 - December 25, 1944) was a German mathematician. ... In numerical analysis, the Runge-Kutta methods are a family of techniques for the approximation of solutions of ordinary differential equations. ... 1910 (MCMX) was a common year starting on Saturday (link will display calendar) of the Gregorian calendar or a common year starting on Sunday of the 13-day slower Julian calendar. ... Lewis Fry Richardson (October 11, 1881 - September 30, 1953) was a mathematician, physicist and psychologist. ... In mathematics, extrapolation is the process of constructing new data points outside a discrete set of known data points. ... 1952 (MCMLII) was a Leap year starting on Tuesday (link will take you to calendar). ... In mathematics, a stiff equation is a differential equation for which certain numerical methods for solving the equation are numerically unstable, unless the step size is taken to be extremely small. ...

See also

• Verlet integration
• Trapezoidal rule
• Courant–Friedrichs–Lewy condition

Verlet integration is a method for calculating the trajectories of particles in molecular dynamics simulations. ... The function f(x) (in blue) is approximated by a linear function (in red). ... In mathematics, the Courant–Friedrichs–Lewy condition (CFL condition) is a condition for certain algorithms for solving partial differential equations to be convergent (not to be confused with numerically stable, though for an explicit it often is within a small constant factor of the stability condition). ...

References

• Ernst Hairer, Syvert Paul Nørsett and Gerhard Wanner, Solving ordinary differential equations I: Nonstiff problems, second edition, Springer Verlag, Berlin, 1993. ISBN 3-540-56670-8.
• Ernst Hairer and Gerhard Wanner, Solving ordinary differential equations II: Stiff and differential-algebraic problems, second edition, Springer Verlag, Berlin, 1996. ISBN 3-540-60452-9.
  (This two-volume monograph systematically covers all aspects of the field.)
• Arieh Iserles, A First Course in the Numerical Analysis of Differential Equations, Cambridge University Press, 1996. ISBN 0-521-55376-8 (hardback), ISBN 0-521-55655-4 (paperback).
  (Textbook, targeting advanced undergraduate and postgraduate students in mathematics, which also discusses numerical partial differential equations.)
• John Denholm Lambert, Numerical Methods for Ordinary Differential Systems, John Wiley & Sons, Chichester, 1991. ISBN 0-471-92990-5.
  (Textbook, slightly more demanding than the book by Iserles.)

Numerical partial differential equations is the branch of numerical analysis that studies the numerical solution of partial differential equations. ...

External links

• Joseph W. Rudmin, Application of the Parker–Sochacki Method to Celestial Mechanics, 1998.
• Dominique Tournès, L'intégration approchée des équations différentielles ordinaires (1671-1914), thèse de doctorat de l'université Paris 7 - Denis Diderot, juin 1996. Réimp. Villeneuve d'Ascq : Presses universitaires du Septentrion, 1997, 468 p. (Extensive online material on ODE numerical analysis history, for English-language material on the history of ODE numerical analysis, see eg the paper books by Chabert and Goldstine quoted by him.)