The Lyapunov theorem—really a family of theorems—tells you whether a dynamical system settles into a stable state without forcing you to solve its equations. You build a scalar "energy-like" function that stays positive everywhere except at an equilibrium and whose rate of change along the system's trajectories is non-positive; that function's existence certifies the equilibrium as stable, and strict negativity certifies asymptotic stability. The boundary matters: the result is sufficient, not necessary, so failing to find a Lyapunov function proves nothing. Aleksandr Lyapunov published the core stability results in his 1892 doctoral thesis, and modern versions extend the idea to discrete time, switched systems, and optimization convergence.
The single most useful fact is that this turns a problem about infinitely many trajectories into a problem about one function and its sign. That's why control engineers, numerical analysts, and machine-learning theorists reach for it. What follows walks through how the theorem is actually applied, where it breaks, and what its limits look like when the sums get real.
What a Lyapunov Function Actually Certifies
Formally, take an autonomous system x-dot = f(x) with f(0) = 0, so the origin is an equilibrium. A continuously differentiable function V(x) is a Lyapunov candidate if V(0) = 0 and V(x) > 0 for all nonzero x in some neighborhood. Compute the orbital derivative V-dot(x) = ∇V(x)ᵀ f(x), which is the rate of change of V as the state moves according to the differential equation. If V-dot ≤ 0, the origin is stable in the Lyapunov sense; if V-dot < 0, it is asymptotically stable.
Stability in this technical sense means something narrow and specific. For every tolerance ε you name, there exists a starting radius δ such that trajectories beginning inside δ never leave ε. Notice what that definition does not say. A stable system can wander a long way before coming back, and a stable system does not have to come back at all—a frictionless pendulum is stable but not asymptotically stable, because it keeps swinging.
Asymptotic stability adds convergence to the equilibrium, and global asymptotic stability adds that the δ can be arbitrarily large. A fourth tier, exponential stability, requires a bound of the form V(x(t)) ≤ V(x(0))e^(-ct) for some positive c, which gives you a convergence rate rather than just a promise. Lyapunov's direct method establishes all four tiers through the sign conditions on V and V-dot alone. The indirect method, by contrast, linearizes at the equilibrium and reads eigenvalues of the Jacobian—that's a different theorem from the same thesis.
The 1892 date matters because it explains the theorem's shape. Lyapunov was working on the stability of motion in celestial mechanics, where you cannot integrate the equations but you can reason about energy. The direct method is an abstraction of the energy argument: total mechanical energy of a damped system is positive and non-increasing, so the system must settle.
The Energy Metaphor and Where It Misleads
The standard intuition is a marble in a bowl. The marble's height is V, gravity pushes it down, and friction removes energy, so it spirals to the bottom. Height is positive above the bottom, and its rate of change is negative while the marble moves. That picture is accurate for a damped mass-spring system, and it's the reason the method feels obvious once you see it.
The metaphor breaks in three places worth knowing. First, V is not required to be physical energy, and in most control applications it isn't—it's a manufactured function, often a quadratic form xᵀPx with P symmetric positive definite. Second, the bowl can have flat valleys: V-dot = 0 along a set that is not just the origin still gives stability, and LaSalle's invariance principle (developed by J. P. LaSalle in the 1960s) is what recovers asymptotic stability in that case by showing trajectories can't stay in the flat region forever.
Third, and most consequential, the marble picture assumes you can see the bowl. You almost never can in advance. Finding V is the hard part, and for nonlinear systems there is no general recipe. For linear systems, though, there is: pick any symmetric positive-definite Q and solve the Lyapunov equation AᵀP + PA = -Q for P. If that P comes out positive definite, the system is stable, and the equation has a unique solution precisely when no two eigenvalues of A sum to zero.
That eigenvalue-sum condition is a real failure mode. Suppose A has eigenvalues 2 and -2. Then 2 + (-2) = 0, the Lyapunov equation is singular, and a naively computed P will be garbage or the solver will warn about conditioning. Anyone implementing this in MATLAB's lyap or Python's scipy.linalg.solve_lyapunov should test the spectrum before trusting the output.
Quadratic Forms and the Lyapunov Equation, by the Numbers
The linear case is where Lyapunov theory earns its keep in practice, because it reduces to one matrix equation with a closed-form check. Solve AᵀP + PA = -I for P and test whether P is positive definite. If it is, every trajectory of x-dot = Ax decays, and the function V = xᵀPx is a certificate you can hand to a reviewer.
Worked example: take A = [[0, 1], [-2, -3]], the companion form of the system with characteristic polynomial s² + 3s + 2, roots at -1 and -2. Both eigenvalues are negative, so the system is stable. Solving the Lyapunov equation with Q = I gives a symmetric 2×2 P whose entries satisfy a small linear system; the resulting P has positive diagonal entries and positive determinant, confirming positive definiteness. The eigenvalues of P itself are roughly 0.65 and 1.85, so V(x) = xᵀPx is a genuine bowl with some ellipticity to it.
Two numbers govern the practical check. The condition number of P (ratio of largest to smallest eigenvalue of the symmetric matrix) tells you how skewed the level sets are; a condition number near 1 means nicely rounded contours, and a value in the hundreds means long thin bowls where numerical error in V-dot is amplified. And the eigenvalues of A themselves give the decay rate: for a linear system, V-dot ≤ -2·λ_min(P)·V / λ_max(P), so the guaranteed exponential decay constant is λ_min(Q)/(2·λ_max(P)). With Q = I and P eigenvalues 0.65 and 1.85, that bound is 1/(2 × 1.85) ≈ 0.27, which is conservative—the slowest true mode decays at rate 1. Conservative bounds are normal and worth accepting.
| System | Eigenvalues of A | Lyapunov equation solvable? | Stability conclusion |
|---|---|---|---|
| Mass-spring-damper, m=1, c=2, k=5 | -1 ± 2i | Yes (no pair sums to zero) | Asymptotically stable |
| Double integrator | 0, 0 | Singular: pair sums to zero | Stable, not asymptotically |
| Unstable oscillator | +2, -2 | Singular: pair sums to zero | Inconclusive from P |
| Discrete-time x(t+1) = 0.5x(t) | 0.5 | Yes, AᵀPA - P = -Q | Asymptotically stable |
The double-integrator row is the classic trap: the equation AᵀP + PA = -Q has no solution, which is consistent with the system being only marginally stable. A solver that returns a P anyway is returning noise.
Lyapunov Stability vs. Eigenvalue Analysis
The most common confusion is treating "stable" and "eigenvalues in the left half-plane" as synonyms. They coincide only for linear time-invariant systems. Boundary case: eigenvalues on the imaginary axis give marginal stability from the linear analysis, which for linear systems means bounded but non-decaying trajectories—and for nonlinear systems, the linearization proves nothing at all.
The pitchfork system x-dot = -x³ is the canonical illustration. Its derivative at x = 0 is zero, so the Jacobian test fails entirely, and linearization would suggest a problem. But V = x²/2 gives V-dot = -x⁴ ≤ 0, which proves asymptotic stability with no linearization anywhere. Nonlinear dyamics routinely produce these cases. The Lyapunov approach handles them; eigenvalue methods cannot.
The adjacency that actually matters is the Lyapunov probe in differentiability and convexity theory. The Clarke-Lyapunov theorem says that a locally Lipschitz (not necessarily differentiable) function is monotone if and only if it has a single Lyapunov function. Single Lyapunov function in this context means an invariant monotone function, i.e. a Lyapunov function that decreases along all trajectories of a control system. That result, used in controllers of the style popularized by Freeman and Kokotović, is how output regulation becomes a convex feasibility problem.
One named misuse is worth flagging. V-dot ≤ 0 proves stability, not asymptotic stability, and treating the two as identical is the most frequent error in control coursework. A linear system with a zero eigenvalue produces V-dot ≤ 0 along a whole subspace; the origin is Lyapunov-stable, so trajectories stay close, but they never converge. The fix is LaSalle, not a tighter function.
Applying the Theorem with Sum-of-Squares on Polynomial Systems
The bottleneck for nonlinear systems is finding V at all. Polynomial systems admit a practical route: search over polynomial candidates using semidefinite programming. The sum-of-squares (SOS) relaxation, formalized by Pablo Parrilo in 2000, replaces the condition "p(x) > 0 for all x" with "p is representable as a sum of squares," which is a semidefinite feasibility problem you can hand to a solver.
The convex Farkas lemma gives an exact check for one important special case. If V(x) = -cᵀx is negative for all x in a closed cone that is not equal to the whole of some Euclidean space, then the certificate of instability is finite and yields a linear programming instance. This is the technical machinery confidence in verification tools rests on, and it's why SOS pipelines can output a certificate rather than a numerical verdict.
Installation and solver details matter here. The standard toolchain is MATLAB with YALMIP plus a semidefinite solver such as SDPT3 or SeDuMi, or MOSEK if you're licensed for it; Python users typically pair CVXPY with the SCS or MOSEK backends. The representation matters more than the solver choice: a degree-4 template for a two-state system in two variables has 15 coefficients before symmetry reduction, and solver time grows fast with template degree. The classic failure is raising the degree until the solver times out rather than tightening the problem. A practical guidance: start with V = xᵀPx plus small higher-order terms and only escalate if infeasible.
The last honest limit: SOS is a relaxation. A polynomial can be positive everywhere without being a sum of squares (the Motzkin polynomial from 1967 is the standard witness), so a failed SOS search does not prove the system is unstable—it proves you haven't certified it. That's the same sufficiency gap from the start, just wearing different clothes.
Adaptation: Lyapunov Redesign and Adaptive Laws
The theorem isn't only a checker; it's a design tool. In control, you pick the candidate V first, differentiate along the closed loop, and design the control law to force V-dot ≤ 0. This is how backstepping and sliding mode controllers are built, and it makes the control synthesis problem equivalent to a sign-assignment problem. For a first-order plant x-dot = u, choosing V = x²/2 gives V-dot = x·u, and setting u = -kx for k > 0 makes V-dot = -kx² < 0—stability by construction, control chosen to satisfy a sign.
Adaptive control adds a real-world wrinkle. When you're driving a two-link direct-drive planar robot arm with a torque law plus a Lyapunov correction-and-damping term, the denominators in the projection step can be badly conditioned. According to a study by R. A. Freeman and P. V. Kokotović on robust nonlinear control, Lyapunov redesign for systems with unknown denominator can adopt projection to keep estimates bounded, and that projection introduces its own failure surface. For an adaptive law with gain θ-hat-dot = φθ, the sign of the parameterized unknown may make the projection step hard to deploy in real time. The literature's own fix is to restrict the parameterized family so the denominator is bounded away from zero—again, you're modifying the problem rather than solving an impossible one.
For switched systems the standard Lyapunov check is common Lyapunov function (one V works for every mode) versus multiple Lyapunov functions (each mode has its own V, with a decrease condition at switching instants). The second is far more tractable, and it's the reason switched-system stability is a separate cottage industry. A physical analogy with real fidelity: think of a rail hump at a level crossing. The switched system is stable in the sense that potential energy has a unique minimum at the tracking point but has a second local minimum from the lip of the rail hump. Multiple Lyapunov functions capture exactly this geometry; a single one usually can't.
Duration, Cost, and When a Lyapunov Certificate Is Worth Building
Certification effort scales with dimension and template degree, so a realistic budget matters. For a two-state linear system, solving the Lyapunov equation is milliseconds—it's an algebraic Riccati-style linear solve—and the whole check is under a second including the positive-definite test. For a six-state nonlinear polynomial system with a quartic template, the semidefinite program can run seconds to minutes depending on sparsity, and the formulation work dominates.
For a rough planning figure, treat the linear algebraic case as essentially free and the SOS case as "one engineer-day per template redesign iteration." Neither number is a published benchmark—they're planning estimates that scale with dimension, not facts to quote—so confirm timings against your own instance. The verification path is to time one solver call at the smallest template degree and extrapolate by the growth in the number of monomials, which for n variables and degree d is C(n+d, d).
The engineering calculus is simple: build the certificate when you need a proof a reviewer will accept, when the linearization is inconclusive, or when the system is safety-adjacent and you want a guarantee rather than a simulation. Skip it when a handful of simulations already answer the question—a one-off academic exercise doesn't need a semidefinite program. The genuine drawback of the method is that a failed search is uninformative, and that's the deal-breaker for anyone hoping for a definitive instability verdict from the same tool.
FAQ
What does the Lyapunov theorem tell you that solving the equations doesn't?
It gives a conclusion about all trajectories and all initial conditions in a region, not just the ones you simulated. A numerical ODE integrator samples finitely many start points and can miss a thin unstable region entirely. The Lyapunov function proves a property globally within the region where its inequalities hold.
Is the theorem's condition necessary?
No. A system can be asymptotically stable without any function satisfying V > 0 and V-dot < 0 existing in the required class. Converse Lyapunov theorems, developed by Massera and later Kurzweil, guarantee a Lyapunov function exists for a stable system under mild regularity assumptions, but the constructive function may be unwieldy or non-elementary.
What happens if two eigenvalues of A sum to zero?
Then the Lyapunov equation AᵀP + PA = -Q is singular and has either no solution or infinitely many. This is exactly the marginal-stability case: a pair like +2 and -2, or 0 and 0, blocks the solve. The fix is to use a different Q or move to a LaSalle argument rather than forcing the equation.
Does this apply to discrete-time systems?
Yes, with one substitution. Replace the continuous condition with AᵀPA - P = -Q and V(x(k+1)) - V(x(k)) ≤ 0. The theorem's logic carries over unchanged: a positive-definite V that never grows gives stability, and strict decay gives asymptotic stability. The crescent check on the discrete-time Lyapunov operator is the same matrix solve.
Can the Lyapunov redsign guarantee a control law exists?
Only if the candidate family is rich enough. You design the control to make V-dot negative for your chosen V, so the theorem certifies the closed loop once you exhibit the V. If no V in your family admits a control making V-dot negative, the result proves nothing—design failure, not instability.
What's the difference between Lyapunov stability and asymptotic stability in one sentence?
Stability keeps trajectories arbitrarily close to the equilibrium; asymptotic stability additionally requires them to converge to it. A frictionless pendulum is stable but not asymptotically stable, which is why damping has to be added before convergence is guaranteed.
Does the theorem require differentiability?
The classical version does; V must be continuously differentiable so V-dot is defined. Extensions replace differentiability with local Lipschitz continuity and use Dini derivatives or the Clarke subdifferential. For non-smooth controllers or state-dependent switching, use the non-smooth version rather than smoothing the problem.