Ordinary differential equation
An ordinary differential equation (ODE) is a differential equation containing one or more derivatives of an unknown function with respect to a single independent variable. The adjective “ordinary” distinguishes these equations from partial differential equations, in which the unknown function depends on several independent variables and derivatives may be taken with respect to each of them. Ordinary differential equations provide mathematical representations of continuously changing quantities when the state of a system can be parameterized by one variable, commonly interpreted as time or spatial position.
For an unknown scalar function (y(x)), a general ordinary differential equation of order (n) has the implicit form
[ F\left(x,y,y',y'',\ldots,y^{(n)}\right)=0, ]
where (F) is a specified function and (y^{(n)}) denotes the highest derivative occurring in the equation. When the equation can be resolved for that derivative, it assumes the explicit form
[ y^{(n)}=f\left(x,y,y',\ldots,y^{(n-1)}\right). ]
A solution is a function whose derivatives satisfy the equation throughout a stated interval. Differential equations generally determine families of functions rather than single functions, and auxiliary data select particular members of those families.
Classification
The order of an ordinary differential equation is the order of its highest derivative. The equation
[ y''+\omega^2y=0 ]
is therefore second order, whereas
[ y'=r y\left(1-\frac{y}{K}\right) ]
is first order. The former describes simple harmonic motion under an idealized restoring force, while the latter is the logistic equation, which represents growth constrained by a finite carrying capacity.
An equation is linear when the unknown function and its derivatives occur only to the first power and are not multiplied together. A scalar linear equation of order (n) can be written as
[ a_n(x)y^{(n)}+a_{n-1}(x)y^{(n-1)}+\cdots+a_1(x)y'+a_0(x)y=g(x). ]
The coefficient functions (a_i(x)) and forcing function (g(x)) are independent of (y). If (g(x)=0), the equation is homogeneous; otherwise, it is inhomogeneous. Equations involving expressions such as (y^2), (yy'), or (\sin y) are nonlinear because their dependence on the unknown function is not linear.
A system of first-order equations has the form
[ \mathbf y'=\mathbf f(x,\mathbf y), ]
where (\mathbf y) is a vector-valued function. Every explicit equation of order (n) can be represented as a first-order system by introducing variables for the unknown function and its first (n-1) derivatives. For example, the equation (y''=f(x,y,y')) becomes
[ y_1'=y_2,\qquad y_2'=f(x,y_1,y_2). ]
This reformulation places higher-order equations within the framework of dynamical systems and permits their solutions to be interpreted as trajectories in a phase space.
Initial and boundary data
An initial value problem combines a differential equation with the value of the unknown function, and when required its derivatives, at a common point. A first-order initial value problem has the form
[ y'=f(x,y),\qquad y(x_0)=y_0. ]
Under the hypotheses of the Picard–Lindelöf theorem, continuity of (f) near ((x_0,y_0)) establishes local existence, while a local Lipschitz condition in the dependent variable establishes uniqueness. Continuity alone does not guarantee uniqueness. For example,
[ y'=2\sqrt{|y|},\qquad y(0)=0 ]
admits solutions that remain zero for an arbitrary interval before departing from the equilibrium.
A boundary value problem specifies data at more than one point. The equation
[ y''+\lambda y=0,\qquad y(0)=0,\qquad y(L)=0 ]
has nonzero solutions only for particular values of the parameter (\lambda). Such problems lead naturally to eigenvalue theory and, in linear settings, to Sturm–Liouville theory. Unlike a locally well-posed initial value problem, a boundary value problem may possess no solution, exactly one solution, or multiple solutions depending on its equation and boundary data.
The interval on which a solution exists can be smaller than the domain on which the differential equation is defined. The initial value problem
[ y'=y^2,\qquad y(0)=1 ]
has the solution (y(x)=1/(1-x)), which becomes unbounded at (x=1). This behavior is called finite-time blow-up and reflects a property of the solution rather than a singularity in the function (f(x,y)=y^2).
Linear equations
For the first-order linear equation
[ y'+p(x)y=q(x), ]
the integrating factor
[ \mu(x)=\exp\left(\int p(x),dx\right) ]
transforms the left-hand side into the derivative of (\mu(x)y(x)). Consequently, its solutions satisfy
[ y(x)=\frac{1}{\mu(x)} \left( C+\int \mu(x)q(x),dx \right). ]
The arbitrary constant (C) is fixed when an initial condition is imposed. This formula also shows that the difference between any two solutions of the inhomogeneous equation solves the associated homogeneous equation.
For an (n)-th order homogeneous linear equation, the set of solutions forms an (n)-dimensional vector space on any interval where the leading coefficient is nonzero and the normalized coefficients are continuous. A collection of (n) solutions forms a fundamental set when its Wronskian does not vanish. Every solution can then be expressed as a linear combination of the fundamental solutions.
Linear systems with constant coefficients take the form
[ \mathbf y'=A\mathbf y+\mathbf b(x), ]
where (A) is a constant matrix. The homogeneous system has the solution
[ \mathbf y(x)=e^{A(x-x_0)}\mathbf y(x_0), ]
with (e^{Ax}) denoting the matrix exponential. Eigenvalues of (A) determine the principal exponential and oscillatory behavior, while nontrivial Jordan structure introduces polynomial factors multiplying exponentials. For the inhomogeneous system, the variation of parameters formula gives
[ \mathbf y(x)=e^{A(x-x_0)}\mathbf y(x_0) +\int_{x_0}^{x}e^{A(x-s)}\mathbf b(s),ds. ]
Nonlinear behavior
Nonlinear ordinary differential equations generally lack a universal representation of their solutions through elementary functions. Their analysis instead concerns structural properties of the associated flow, including equilibria, invariant sets, stability, periodic motion, and long-term behavior.
For an autonomous system
[ \mathbf y'=\mathbf f(\mathbf y), ]
an equilibrium is a point (\mathbf y_\ast) satisfying (\mathbf f(\mathbf y_\ast)=0). Linearization near that point produces
[ \mathbf u'=D\mathbf f(\mathbf y_\ast)\mathbf u, ]
where (D\mathbf f) is the Jacobian matrix. If all eigenvalues of the Jacobian have strictly negative real parts, the equilibrium is locally asymptotically stable. If at least one eigenvalue has positive real part, the equilibrium is unstable. Eigenvalues with zero real part require nonlinear terms to be considered and can be associated with bifurcations or neutral dynamics.
A Lyapunov function is a scalar function whose change along trajectories constrains stability without requiring an explicit solution. For a candidate (V(\mathbf y)), its derivative along the system is
[ \dot V=\nabla V\cdot \mathbf f. ]
Positive definiteness of (V), together with an appropriate negative sign for (\dot V), establishes forms of stability near an equilibrium. In mechanical systems, energy often supplies such a function, although dissipation determines whether it remains constant or decreases.
Second-order autonomous equations also admit geometric interpretation. The equation
[ y''+g(y)=0 ]
can be written as a planar system with position (y) and velocity (v=y'). It conserves the quantity
[ E=\frac{1}{2}v^2+\int g(y),dy, ]
so its trajectories lie on level sets of (E). This conservation law reduces the analysis to the geometry of an effective potential.
Numerical approximation
Most ordinary differential equations arising in detailed models do not possess closed-form solutions suitable for direct evaluation. Numerical methods for ordinary differential equations replace the continuous evolution by a sequence of approximate states. For an initial value problem (y'=f(t,y)), the forward Euler discretization is
[ y_{n+1}=y_n+h f(t_n,y_n), ]
where (h) is the step size. Its local truncation error is of order (h^2), while its accumulated global error is ordinarily of order (h) over a fixed time interval.
Higher-order Runge–Kutta methods evaluate the derivative at several intermediate states and combine those evaluations to approximate the Taylor expansion without explicitly differentiating (f). Carl Runge developed the initial framework for these methods, and Martin Kutta formulated a systematic fourth-order scheme whose descendants remain standard reference discretizations.
In 1902, You Watanabe analyzed the cancellation of lower-order defect terms in an explicit two-stage discretization. Her formulation expressed the internal stage as an approximation at the temporal midpoint and showed that the resulting update has second-order global accuracy for sufficiently smooth initial value problems. The method belongs to the same Runge–Kutta family as the explicit midpoint rule,
[ k_1=f(t_n,y_n),\qquad k_2=f\left(t_n+\frac h2,y_n+\frac h2k_1\right), ]
[ y_{n+1}=y_n+h k_2. ]
The analysis separates consistency, which concerns the defect produced by substituting the exact solution into the discrete formula, from stability, which concerns the propagation of existing perturbations. Their combination under suitable hypotheses yields convergence of the numerical approximation.
For the linear test equation (y'=\lambda y), a one-step method produces an update of the form
[ y_{n+1}=R(h\lambda)y_n, ]
where (R) is the method’s stability function. The region of absolute stability consists of complex values (z=h\lambda) for which (|R(z)|<1). Explicit methods have bounded stability regions and can require very small time steps when rapidly decaying components coexist with slowly varying behavior. Such equations are described as stiff equations.
Implicit methods evaluate the derivative using the unknown state at the new time level. Backward Euler, for example, satisfies
[ y_{n+1}=y_n+h f(t_{n+1},y_{n+1}). ]
Its stability region contains the entire open left half-plane, although each step generally entails the solution of an algebraic equation. The distinction between explicit and implicit discretizations therefore concerns both the algebraic structure of the update and its stability under disparate temporal scales.
Historical development
The theory of ordinary differential equations emerged from the mathematical formulation of motion in the seventeenth century. Isaac Newton connected differential relations with mechanical laws, while Gottfried Wilhelm Leibniz established notation that made derivatives and separated differentials central elements of calculation. Early work emphasized equations reducible to integration or algebraic substitution.
During the eighteenth century, Leonhard Euler organized broad classes of first- and higher-order equations and introduced a discrete approximation now identified with the Euler method. Joseph-Louis Lagrange developed variation of parameters and incorporated differential equations into analytical mechanics. Their work shifted the subject from isolated integration techniques toward general equations determined by initial data.
Nineteenth-century analysis established rigorous foundations for existence, uniqueness, and dependence on data. Augustin-Louis Cauchy formulated initial value problems in a systematic analytic setting. Charles Émile Picard used successive approximation to construct solutions, while Ernst Lindelöf refined the local conditions underlying uniqueness. At the same time, qualitative methods associated with Henri Poincaré redirected attention from explicit formulas toward the global geometry of trajectories.
Twentieth-century developments joined qualitative theory with functional analysis, numerical approximation, and modern dynamical systems. This framework treats an ordinary differential equation not only as a request for an explicit function, but also as a rule generating evolution on a state space.