Advanced Mathematics

Dynamical Systems: Flows, Chaos, and Bifurcations

Dynamical systems is the mathematical study of how systems evolve over time. From the pendulum and predator-prey models to the butterfly-shaped Lorenz attractor and the infinite complexity of the Mandelbrot set, dynamical systems unifies geometry, analysis, and physics into one of mathematics' most powerful frameworks.

Learning Objectives

After working through this page, you will be able to:

  • Sketch phase portraits for two-dimensional continuous systems and classify equilibria by their eigenvalues
  • Construct Lyapunov functions to prove stability or instability of equilibria
  • Apply the Hartman-Grobman theorem and center manifold theorem to analyze nonlinear systems near equilibria
  • Identify and classify saddle-node, pitchfork, and Hopf bifurcations from a parametric vector field
  • Use the Poincare-Bendixson theorem and Bendixson's criterion to rule out or confirm the existence of limit cycles
  • Compute Lyapunov exponents and understand their role in characterizing chaos
  • Analyze the logistic map, period-doubling cascade, and Feigenbaum universality
  • State Hamilton's equations and apply them to conservative mechanical systems

1. Continuous Dynamical Systems: Flows and Phase Portraits

A continuous dynamical system is specified by an autonomous differential equation x' = f(x), where x lives in a state space (usually R^n or a smooth manifold) and f is a vector field. Autonomous means f does not depend explicitly on time, so the behavior is determined entirely by the current state.

Flows

The flow of a vector field f is the family of maps phi(t, x0) giving the position at time t of a trajectory starting at x0. Formally, phi satisfies:

  • d/dt phi(t, x0) = f(phi(t, x0))
  • phi(0, x0) = x0
  • phi(t + s, x0) = phi(t, phi(s, x0)) (group property)

The group property encodes the fact that evolving for time s then time t is the same as evolving for time t + s. Existence and uniqueness of flows follow from the Picard-Lindelof theorem when f is Lipschitz.

Phase Portraits

The phase portrait is a picture of the state space showing representative trajectories. For two-dimensional systems, the phase portrait is a planar vector field together with the curves (called orbits) traced out by trajectories. Key features to identify:

Equilibria (Fixed Points)

Points where f(x*) = 0. The system does not move from an equilibrium. Classified by the eigenvalues of the Jacobian Df(x*): stable node, unstable node, saddle, stable spiral, unstable spiral, or center.

Nullclines

Curves where one component of f vanishes. The x-nullcline is where x' = 0 (vertical tangents), and the y-nullcline is where y' = 0 (horizontal tangents). Equilibria occur at their intersections.

Classification of Planar Equilibria

For the linear system x' = Ax, the character of the equilibrium at the origin is determined by the eigenvalues lambda1, lambda2 of A. Let T = trace(A) and D = det(A).

Eigenvalue conditionTypeStability
D < 0SaddleUnstable
D > 0, T < 0, T^2 - 4D > 0Stable nodeAsymptotically stable
D > 0, T > 0, T^2 - 4D > 0Unstable nodeUnstable
D > 0, T < 0, T^2 - 4D < 0Stable spiralAsymptotically stable
D > 0, T > 0, T^2 - 4D < 0Unstable spiralUnstable
D > 0, T = 0CenterLyapunov stable (not asymptotic)

Example: Van der Pol Oscillator

x'' - mu(1 - x^2) x' + x = 0, mu > 0

Written as a system: x' = y, y' = mu(1 - x^2)y - x. The origin is the only equilibrium. For mu small, the Jacobian has eigenvalues mu/2 plus-or-minus i*sqrt(1 - mu^2/4), giving an unstable spiral. By Poincare-Bendixson, the system has a stable limit cycle surrounding the origin.

2. Stability Theory and Lyapunov Functions

Lyapunov's method, developed by Aleksandr Lyapunov in 1892, is the most general tool for analyzing stability without solving the differential equation explicitly. The central idea: if a scalar function decreases along trajectories, the system must converge to a lower-energy state.

Lyapunov Stability Definitions

Lyapunov Stable

The equilibrium x* is Lyapunov stable if for every epsilon greater than 0, there exists delta greater than 0 such that |x(0) - x*| less than delta implies |x(t) - x*| less than epsilon for all t greater than or equal to 0.

Asymptotically Stable

The equilibrium x* is asymptotically stable if it is Lyapunov stable AND there exists delta > 0 such that |x(0) - x*| less than delta implies x(t) approaches x* as t approaches infinity.

Globally Asymptotically Stable

Asymptotically stable from every initial condition in the state space. This is the strongest form of stability and is often the goal in control engineering.

Lyapunov's Direct Method

A Lyapunov function V: U to R is a continuously differentiable function on a neighborhood U of x* satisfying:

  • V(x*) = 0
  • V(x) > 0 for x not equal to x* in U (positive definite)
  • dV/dt = grad(V) cdot f(x) leq 0 in U (non-increasing along orbits)

If dV/dt is strictly negative (negative definite), then x* is asymptotically stable. Think of V as a generalized energy function that the trajectory cannot increase.

LaSalle's Invariance Principle

LaSalle's principle extends Lyapunov's method to cases where dV/dt is only negative semi-definite (zero on some set, not just at the equilibrium). Let E = (x : dV/dt = 0) and let M be the largest invariant set contained in E. Then every trajectory starting in a compact sublevel set (V(x) less than c) converges to M as t approaches infinity.

Example: Damped Pendulum via LaSalle

theta'' + b*theta' + sin(theta) = 0, b > 0

Take V = (1/2)*theta'^2 + (1 - cos theta) (total mechanical energy). Then dV/dt = -b*(theta')^2, which is zero when theta' = 0. The set E = (theta' = 0). On E, the equation gives sin(theta) = 0, so theta must be a multiple of pi. The only invariant set in E near the origin is theta = 0, theta' = 0. By LaSalle, the equilibrium is asymptotically stable.

Instability: Chetaev's Theorem

Chetaev's theorem gives sufficient conditions for instability. If there exists a function V and a region U near x* where V takes positive values, and dV/dt is strictly positive wherever V is positive in U, with V = 0 on the boundary of U, then x* is unstable. This provides a counterpart to Lyapunov's stability theorem.

3. Linear Systems: Matrix Exponential and Invariant Manifolds

The linear system x' = Ax (A an n-by-n matrix) has an explicit global solution: x(t) = exp(At) * x(0). Understanding the matrix exponential and the spectral decomposition of A is the gateway to both the linear theory and the local theory of nonlinear systems.

The Matrix Exponential

The matrix exponential is defined by the convergent power series:

exp(At) = I + At + (At)^2/2! + (At)^3/3! + ...

This series converges for all matrices A and all times t. Key properties:

  • d/dt exp(At) = A * exp(At)
  • exp(A*0) = I
  • exp(A(s+t)) = exp(As) * exp(At) (group property)
  • If A = PJP^(-1), then exp(At) = P * exp(Jt) * P^(-1)

Jordan Normal Form

Every matrix A over C is similar to a Jordan matrix J, a block-diagonal matrix with Jordan blocks. A Jordan block for eigenvalue lambda is an upper bidiagonal matrix with lambda on the diagonal and 1s on the superdiagonal. The matrix exponential of a Jordan block J_k(lambda) (size k) is:

exp(J_k(lambda) * t) = e^(lambda*t) * [1, t, t^2/2!, ..., t^(k-1)/(k-1)!; 0, 1, t, ...; ...]

The polynomial factor t^j multiplying the exponential explains why repeated eigenvalues can cause polynomial growth before exponential decay dominates.

Stable, Unstable, and Center Subspaces

For x' = Ax, the state space decomposes into three invariant subspaces determined by the eigenvalues of A:

Stable Subspace E^s

Spanned by (generalized) eigenvectors for eigenvalues with negative real part. Trajectories on E^s decay to 0 exponentially.

Unstable Subspace E^u

Spanned by eigenvectors for eigenvalues with positive real part. Trajectories on E^u grow exponentially.

Center Subspace E^c

Spanned by eigenvectors for purely imaginary (or zero) eigenvalues. Trajectories on E^c neither grow nor decay exponentially.

For a hyperbolic linear system (no eigenvalues on the imaginary axis), E^c = 0 and the state space is the direct sum of E^s and E^u. The general solution x(t) = exp(At)x(0) decays to 0 if and only if all eigenvalues of A have negative real part.

4. Nonlinear Systems: Hartman-Grobman and Center Manifolds

For nonlinear systems, local analysis near an equilibrium reduces to the linear theory via two fundamental theorems. The Hartman-Grobman theorem handles hyperbolic equilibria; the center manifold theorem handles the degenerate case where the linearization has eigenvalues on the imaginary axis.

The Hartman-Grobman Theorem

Let x* be a hyperbolic equilibrium of x' = f(x) (meaning no eigenvalue of Df(x*) has zero real part). Then there is a homeomorphism h defined in a neighborhood of x* such that h maps trajectories of the nonlinear system to trajectories of the linearization x' = Df(x*)(x - x*). The homeomorphism h preserves the direction of time but is generally not differentiable.

Consequence

Near a hyperbolic equilibrium, the qualitative behavior (node, saddle, spiral) is completely determined by the eigenvalues of the Jacobian Df(x*). In particular, if all eigenvalues have negative real part, x* is asymptotically stable; if any eigenvalue has positive real part, x* is unstable. The nonlinear terms do not affect the qualitative local picture at a hyperbolic equilibrium.

The Center Manifold Theorem

When Df(x*) has eigenvalues with zero real part, the center subspace E^c captures the slow dynamics. The center manifold theorem guarantees the existence (local, typically not global) of a center manifold W^c tangent to E^c at x*. Crucially, the stability of x* is determined by the dynamics on W^c alone, which is a lower-dimensional system.

Computation via Taylor Expansion

For a system written in coordinates (x, y) with x in E^c and y in E^s, write y = h(x) for the center manifold. Substitute into the equation for y' to get the invariance equation: Dh(x) * f_c(x, h(x)) = f_s(x, h(x)), where f_c, f_s are the projections onto the center and stable directions. Solve for h order by order in Taylor series.

Normal Forms

Normal form theory provides a systematic way to simplify a vector field near an equilibrium by near-identity coordinate changes. The goal is to eliminate as many nonlinear terms as possible while preserving the topological structure. Resonance conditions determine which terms cannot be eliminated: a monomial x^alpha in the k-th equation is resonant if the eigenvalue relation lambda_k = alpha cdot lambda holds, where lambda = (lambda_1, ..., lambda_n) are the eigenvalues.

Poincare Linearization Theorem

If no resonance relations hold among the eigenvalues, the system can be linearized by a formal power series transformation. When eigenvalues satisfy a Poincare condition (they all lie in a half-plane not containing the origin of the complex plane), the transformation converges, giving an actual smooth linearization. This is Poincare's theorem, which precedes the Hartman-Grobman result.

5. Bifurcation Theory

Bifurcation theory studies how the qualitative behavior of a parametric family of systems x' = f(x, mu) changes as the parameter mu varies. A bifurcation occurs at a parameter value where the system is not structurally stable — a small change in mu produces a qualitatively different phase portrait.

Saddle-Node Bifurcation

The simplest and most fundamental bifurcation: two equilibria (one stable, one unstable) collide and annihilate as the parameter crosses the bifurcation value. The normal form is:

x' = mu - x^2

  • mu > 0: two equilibria at x = plus-or-minus sqrt(mu); x = +sqrt(mu) stable, x = -sqrt(mu) unstable
  • mu = 0: one equilibrium at x = 0 (semi-stable); this is the bifurcation point
  • mu < 0: no real equilibria; trajectories drift to -infinity

In applications, saddle-node bifurcations model the sudden appearance or disappearance of steady states, such as the onset of bistability in chemical reactions or saddle-node bifurcations on invariant circles (SNIC), which are a route to oscillation in neural models.

Pitchfork Bifurcation

The pitchfork bifurcation occurs in systems with a symmetry x to -x. The normal forms are:

Supercritical (soft)

x' = mu*x - x^3

For mu less than 0: only x = 0 (stable). For mu greater than 0: x = 0 becomes unstable and two new stable equilibria appear at x = plus-or-minus sqrt(mu). Safe, smooth transition.

Subcritical (hard)

x' = mu*x + x^3

For mu less than 0: x = 0 is stable alongside two unstable equilibria at x = plus-or-minus sqrt(-mu). For mu greater than 0: only x = 0 remains and it is unstable. Trajectories jump discontinuously — hysteresis.

Hopf Bifurcation

The Hopf bifurcation is the birth (or death) of a periodic orbit from an equilibrium. It occurs when a complex conjugate pair of eigenvalues crosses the imaginary axis as the parameter varies.

Conditions (Poincare-Andronov-Hopf)

  • The Jacobian at mu = mu_0 has eigenvalues alpha(mu) plus-or-minus i*omega(mu) with alpha(mu_0) = 0 and omega(mu_0) not equal to 0
  • The eigenvalues cross the imaginary axis transversally: d(alpha)/d(mu) at mu_0 is not zero
  • A first Lyapunov coefficient l_1 determines the direction:
  • l_1 < 0: supercritical Hopf (stable limit cycle is born)
  • l_1 > 0: subcritical Hopf (unstable limit cycle collapses)

Example: Glycolytic Oscillator

In biochemistry, the Selkov model for glycolytic oscillations undergoes a supercritical Hopf bifurcation as the substrate injection rate passes a critical value. Below the bifurcation, the system reaches a steady state for ATP production. Above it, sustained oscillations in ATP and ADP concentrations occur — a limit cycle in the phase plane.

Bifurcation Diagrams

A bifurcation diagram plots the location and stability of equilibria (and periodic orbits) as a function of the parameter mu. Solid curves indicate stable states; dashed curves indicate unstable states. Fold points (turning points) on the bifurcation curve correspond to saddle-node bifurcations. The diagram encodes hysteresis: the system state can depend on the direction of parameter variation, not just its current value.

6. Periodic Orbits: Limit Cycles and the Poincare Map

Periodic orbits — closed trajectories in phase space — are among the most important geometric structures in dynamical systems. They correspond to sustained oscillations and are studied via the Poincare return map, which reduces a continuous flow to a discrete map on a lower-dimensional section.

Poincare Return Map

Choose a codimension-1 surface S (the Poincare section) that is transverse to the flow. For a trajectory starting at a point x on S, define P(x) as the first return of the trajectory to S. A fixed point of P corresponds to a periodic orbit of the flow. The stability of the periodic orbit is determined by the eigenvalues of DP(x*) — the Floquet multipliers (also called characteristic multipliers).

Floquet Theory

For a periodic orbit of period T, the linearization of the flow around the orbit is a T-periodic linear system. By Floquet's theorem, the fundamental matrix solution has the form Phi(t) = P(t) * exp(Bt), where P(t) is T-periodic and B is constant. The eigenvalues of exp(BT) are the Floquet multipliers. One multiplier is always 1 (corresponding to perturbations along the orbit). The orbit is stable if all other multipliers lie inside the unit circle.

Bendixson's Criterion

Bendixson's criterion rules out the existence of closed orbits in simply connected planar regions. If div f = (partial f_1/partial x) + (partial f_2/partial y) does not change sign (and is not identically zero) on a simply connected region D, then the system x' = f(x) has no closed orbits lying entirely in D.

Proof idea via Green's theorem

If a closed orbit C existed, Green's theorem would give the integral of div f over the interior equal to the line integral of f cdot n around C, which is zero since f is tangent to C. But if div f does not change sign, the area integral is nonzero — contradiction.

Dulac's Theorem (Generalized Bendixson)

Dulac's theorem generalizes Bendixson by allowing a multiplier function. If there exists a smooth function rho(x, y) on a simply connected domain D such that div(rho*f) does not change sign on D, then there are no closed orbits in D. The choice of rho provides more flexibility to rule out limit cycles.

Poincare-Bendixson Theorem

The Poincare-Bendixson theorem is one of the most powerful results for planar systems. It states: if a trajectory remains in a compact (closed and bounded) region of the plane for all future time, and that region contains only finitely many equilibria, then the omega-limit set of the trajectory is either an equilibrium, a closed orbit, or a trajectory connecting equilibria (homoclinic or heteroclinic orbit).

Applying Poincare-Bendixson

To prove existence of a limit cycle: (1) Find an annular region that is positively invariant (flow points inward on both boundaries). (2) If the only equilibrium in the region has been ruled out as an omega-limit set (e.g., by showing it is an unstable spiral), then the omega-limit set must be a closed orbit — a limit cycle.

7. Chaos: Sensitive Dependence and Strange Attractors

Chaos is a precise mathematical phenomenon: deterministic systems exhibiting sensitive dependence on initial conditions, making long-term prediction practically impossible even though the system follows exact laws. Chaos was dramatically demonstrated by Edward Lorenz in 1963 when he discovered that a simplified model of atmospheric convection produced aperiodic, unpredictable behavior.

Sensitive Dependence on Initial Conditions

A system is sensitively dependent on initial conditions if there exists delta greater than 0 such that for every x and every neighborhood U of x, there exists y in U and time t such that the distance from phi(t, x) to phi(t, y) is greater than delta. Intuitively, nearby trajectories eventually diverge by a macroscopic amount, regardless of how close they started.

The Butterfly Effect

Lorenz illustrated this with his famous metaphor: a butterfly flapping its wings in Brazil could, in principle, set off a tornado in Texas. This is not because the butterfly causes the tornado, but because the atmosphere is a chaotic system where infinitesimal differences in initial conditions grow exponentially, making long-range weather forecasting fundamentally limited, not just computationally.

Lyapunov Exponents

Lyapunov exponents measure the average exponential rate of divergence of nearby trajectories. For a trajectory starting at x0, the maximal Lyapunov exponent is:

lambda = lim (t to infinity) (1/t) * ln |Dphi_t(x0) * v| / |v|

where v is a generic initial tangent vector. A positive maximal Lyapunov exponent is the defining quantitative signature of chaos. For an n-dimensional system, there are n Lyapunov exponents (the Lyapunov spectrum). For a dissipative system, the sum of all Lyapunov exponents equals the time-averaged divergence of f (equals trace of Df averaged over the attractor), which is negative for dissipative systems. A chaotic attractor must have at least one positive and one negative exponent.

The Lorenz System

The Lorenz equations are a three-dimensional autonomous system originally derived from a simplified model of Rayleigh-Benard convection:

x' = sigma(y - x)

y' = x(rho - z) - y

z' = xy - beta*z

Standard parameters: sigma = 10, rho = 28, beta = 8/3. The system has three equilibria: the origin and two symmetric points at plus-or-minus (sqrt(beta(rho-1)), sqrt(beta(rho-1)), rho-1). At the standard parameters all three are unstable.

The Lorenz attractor is a strange attractor — a fractal subset of the phase space with Hausdorff dimension approximately 2.06. It has an intricate two-lobed butterfly shape. Rigorous mathematical proof that the Lorenz attractor is chaotic (not merely numerically observed) was achieved by Tucker in 2002 using interval arithmetic computations.

Strange Attractors and Fractal Dimension

An attractor is a compact invariant set that attracts nearby trajectories. A strange attractor is an attractor with fractal (non-integer) Hausdorff dimension. The fractal structure reflects the simultaneous stretching (sensitivity to initial conditions) and folding (confinement to a bounded region) that characterizes dissipative chaos. The Kaplan-Yorke formula estimates the information dimension from the Lyapunov spectrum:

D_KY = j + (lambda_1 + ... + lambda_j) / |lambda_(j+1)|

where j is the largest index such that the sum of the first j Lyapunov exponents is non-negative. For the Lorenz attractor with standard parameters, D_KY approximately equals 2.06.

8. Discrete Dynamical Systems: Iterated Maps and the Logistic Map

Discrete dynamical systems are defined by a map f: X to X and the iteration x(n+1) = f(x(n)). Despite their simple description, discrete systems can exhibit the full complexity of continuous dynamics — and their lower dimensionality makes them tractable for rigorous analysis.

Fixed Points and Periodic Orbits

A fixed point satisfies f(x*) = x*. A period-n orbit satisfies f^n(x) = x but f^k(x) not equal to x for k = 1, ..., n-1. The stability of a fixed point is determined by the derivative: if |f'(x*)| is less than 1, the fixed point is stable (attracting); if |f'(x*)| is greater than 1, it is unstable (repelling); if |f'(x*)| equals 1, the stability is inconclusive from the linear analysis alone.

The Logistic Map

The logistic map is the canonical example of deterministic chaos in a one-dimensional discrete system:

x(n+1) = r * x(n) * (1 - x(n)), x in [0,1], r in [0,4]

Parameter rangeBehavior
0 < r < 1Population goes extinct; x* = 0 is stable
1 < r < 3Stable fixed point x* = 1 - 1/r
r = 3Period-2 orbit born via period-doubling bifurcation
3 < r < ~3.57Cascade of period doublings: 2, 4, 8, 16, ...
r > ~3.57Chaotic behavior with windows of periodicity

Feigenbaum Universality

Let r_n be the parameter value at which the period-2^n orbit is born. The ratios (r_(n+1) - r_n) / (r_(n+2) - r_(n+1)) converge to Feigenbaum's constant delta approximately 4.669201... This universal constant is the same for any smooth unimodal map (one with a single maximum), regardless of the specific functional form. A second Feigenbaum constant alpha approximately 2.5029... governs the self-similar scaling of the bifurcation diagram.

Universality

Feigenbaum universality was one of the first discoveries in chaos theory with practical experimental consequences. The period-doubling cascade and the constants delta and alpha were subsequently observed in physical experiments: fluid convection, electronic circuits, and cardiac rhythms all show bifurcation cascades converging at Feigenbaum's ratio, regardless of the physical mechanism.

The Mandelbrot Set

The Mandelbrot set M is the set of complex parameters c for which the iteration z(n+1) = z(n)^2 + c, starting from z(0) = 0, does not escape to infinity. The boundary of M is the most complex known fractal — it has infinite topological complexity (every boundary point is a limit of both interior and exterior points) and Hausdorff dimension 2 (the full plane dimension, despite being a one-dimensional curve in many senses).

Julia Sets and the MLC Conjecture

For each parameter c, the corresponding Julia set J_c is the boundary of the set of initial conditions that do not escape. When c is in the Mandelbrot set, J_c is a connected set; when c is outside M, J_c is a Cantor set (totally disconnected). The Mandelbrot-Local Connectivity conjecture (MLC), still open, would imply that M is locally connected and completely determine the combinatorial structure of the Mandelbrot set.

9. Hamiltonian Systems: Symplectic Structure and Conservation

Hamiltonian systems are the mathematical framework for classical mechanics. Unlike dissipative systems, Hamiltonian systems preserve volume in phase space (Liouville's theorem) and energy, and have a rich geometric structure based on symplectic geometry.

Hamilton's Equations

A Hamiltonian system is defined by a scalar function H(q, p, t) — the Hamiltonian — on the phase space with coordinates (q_1, ..., q_n, p_1, ..., p_n). The equations of motion are:

dq_i/dt = partial H / partial p_i

dp_i/dt = -(partial H / partial q_i)

For an autonomous system (H does not depend on t), H is a first integral: dH/dt = 0 along any trajectory. This is conservation of energy. The q_i are generalized coordinates and the p_i are generalized momenta (conjugate momenta).

Symplectic Structure

The phase space of a Hamiltonian system carries a natural symplectic form omega = dq_1 wedge dp_1 + ... + dq_n wedge dp_n. The symplectic form is a closed, non-degenerate 2-form. Hamilton's equations can be written compactly as z' = J * grad H, where z = (q, p) and J is the 2n-by-2n standard symplectic matrix.

J = [[0, I_n], [-I_n, 0]]

A transformation (q, p) to (Q, P) is symplectic (canonical) if it preserves the symplectic form: J_T = omega. Symplectic maps preserve volumes, areas of 2D projections, and all higher symplectic invariants (by Gromov's nonsqueezing theorem, the most subtle).

Liouville's Theorem

The Hamiltonian flow preserves volume in phase space: if a region R(0) evolves under the flow to R(t), then volume(R(t)) = volume(R(0)) for all t. This follows from the fact that the divergence of the Hamiltonian vector field is identically zero (the div of (partial H/partial p, -partial H/partial q) cancels by equality of mixed partials). Liouville's theorem has profound consequences: Hamiltonian systems cannot have attractors in the usual sense, and statistical mechanics relies on it for the ergodic hypothesis.

Integrable Systems and the KAM Theorem

A Hamiltonian system with n degrees of freedom is completely integrable (in the Liouville sense) if it has n independent first integrals F_1 = H, F_2, ..., F_n that are in involution (the Poisson bracket of any pair is zero). By the Liouville-Arnold theorem, the level sets of the first integrals are invariant tori, and the motion on each torus is quasi-periodic.

KAM Theorem (Kolmogorov-Arnold-Moser)

When an integrable Hamiltonian system is perturbed by a small perturbation epsilon*H_1, most invariant tori (those with sufficiently irrational frequency ratios, satisfying a Diophantine condition) persist and are only slightly deformed. The measure of destroyed tori is of order epsilon^(1/2). The remaining tori (KAM tori) form a Cantor-like set of positive measure. Between the KAM tori, destroyed rational tori leave behind chains of island and chaotic layers (the Birkhoff islands and Chirikov resonances).

Poisson Brackets and Canonical Transformations

The Poisson bracket of two functions F and G on phase space is:

(F, G) = sum_i (partial F/partial q_i * partial G/partial p_i - partial F/partial p_i * partial G/partial q_i)

A function F is a first integral if and only if (F, H) = 0. Hamilton's equations become dF/dt = (F, H) + partial F/partial t. A transformation is canonical if it preserves all Poisson brackets: (Q_i, Q_j) = 0, (P_i, P_j) = 0, (Q_i, P_j) = delta_ij. Canonical transformations generated by Hamilton's principal function S are the foundation of the Hamilton-Jacobi theory.

10. Applications: From Biology to Robotics

Dynamical systems provides the mathematical language for virtually every scientific domain where change over time matters. The following applications illustrate how the abstract theory connects directly to real phenomena.

Population Dynamics (Ecology)

The Lotka-Volterra predator-prey system:

x' = ax - bxy (prey)
y' = cxy - dy (predator)

The system has a nontrivial equilibrium (d/c, a/b) that is a center for the linearization. The system is Hamiltonian with H = c*x - d*ln x + b*y - a*ln y, giving closed orbits (periodic population cycles). Adding competition or carrying-capacity terms breaks the Hamiltonian structure and creates stable limit cycles via Hopf bifurcation.

Chemical Reactions (Oregonator)

The Belousov-Zhabotinsky reaction exhibits chemical oscillations. The Field-Koros-Noyes (Oregonator) model reduces the mechanism to three species. Near a Hopf bifurcation, the system transitions from a steady state (chemical equilibrium) to limit cycle oscillations visible as color-changing waves in a petri dish.

Climate Models (Energy Balance)

Simple energy balance climate models take the form C*dT/dt = Q*(1 - alpha(T)) - epsilon*sigma*T^4, where T is global mean temperature, alpha(T) is the ice-albedo feedback (nonlinear). This one-dimensional system can exhibit saddle-node bifurcations corresponding to the abrupt climate transitions between snowball-Earth and warm states — tipping points with hysteresis.

Neural Networks (Hodgkin-Huxley)

The Hodgkin-Huxley equations model action potential generation in a neuron as a four-dimensional dynamical system. The resting potential is a stable equilibrium. When the applied current exceeds a threshold, a Hopf bifurcation (or SNIC bifurcation, depending on the neuron type) produces repetitive spiking — a limit cycle. The FitzHugh-Nagumo model is a two-dimensional simplification amenable to phase-plane analysis.

Robotics and Control (Passive Walking)

Passive dynamic walkers — robots with no actuators — walk down a gentle slope through the interaction of gravity and their leg dynamics. The Poincare map for the step-to-step dynamics has a stable fixed point corresponding to steady walking. Period-doubling bifurcations as the slope increases lead to period-2 (limping) gait and eventually to chaotic falling.

Celestial Mechanics (Three-Body Problem)

The gravitational three-body problem is a Hamiltonian system that is not integrable in general. Near Lagrange points, the linearized system can be used to find families of periodic orbits (Lyapunov orbits, halo orbits) used by spacecraft. KAM theory explains why the solar system is stable over geological timescales despite perturbations: most orbits lie on KAM tori, though the system is technically chaotic on timescales of tens of millions of years.

11. Practice Problems with Full Solutions

Problem 1: Classifying an Equilibrium

Consider the system x' = y - x^3, y' = -x. Find all equilibria and classify them.

Show Solution

Step 1: Find equilibria. Set y - x^3 = 0 and -x = 0. From the second equation, x = 0. Substituting, y = 0. The only equilibrium is (0, 0).

Step 2: Compute the Jacobian.

Df = [[-3x^2, 1], [-1, 0]]
At (0,0): Df = [[0, 1], [-1, 0]]

Step 3: Compute eigenvalues. det(Df - lambda*I) = lambda^2 + 1 = 0, so lambda = plus-or-minus i. The linearization is a center.

Step 4: Nonlinear analysis. Since the eigenvalues are purely imaginary, Hartman-Grobman does not apply (non-hyperbolic). We need a Lyapunov function. Try V = x^2 + y^2. Then dV/dt = 2x*x' + 2y*y' = 2x(y - x^3) + 2y(-x) = 2xy - 2x^4 - 2xy = -2x^4 leq 0. Since dV/dt is negative semi-definite (zero only when x = 0), LaSalle's principle applies. On the set (x = 0), the trajectory must have y' = -x = 0 and y'' = -x' = -(y - x^3) = -y, giving y = 0. The origin is asymptotically stable, not just a center.

Problem 2: Lyapunov Function Construction

Prove that the origin is globally asymptotically stable for x' = -x^3 - y, y' = x - y^3.

Show Solution

Candidate: V(x, y) = x^2 + y^2. This is positive definite, radially unbounded (V(x,y) to infinity as |(x,y)| to infinity), and V(0,0) = 0.

Compute dV/dt:

dV/dt = 2x*(-x^3 - y) + 2y*(x - y^3)
= -2x^4 - 2xy + 2xy - 2y^4
= -2x^4 - 2y^4
leq 0

Since dV/dt = -2(x^4 + y^4) is strictly negative for (x,y) not equal to (0,0), V is a strict Lyapunov function. By Lyapunov's theorem, the origin is asymptotically stable. Since V is also radially unbounded (a proper function), the origin is globally asymptotically stable by Barbalat's lemma / La Salle's theorem applied globally. (Every sublevel set (V leq c) is compact and positively invariant.)

Problem 3: Ruling Out Limit Cycles

Show that the system x' = x - y - x(x^2 + 2y^2), y' = x + y - y(x^2 + 2y^2) has no limit cycles.

Show Solution

Strategy: We will show the system has a globally attracting invariant set (so no separate limit cycle can exist), or use Bendixson on a suitable region. Actually, let's compute in polar coordinates to identify the attractor.

Let x = r cos(theta), y = r sin(theta). Then:

r' = x*x'/r + y*y'/r
= x(x - y - x(x^2 + 2y^2))/r + y(x + y - y(x^2 + 2y^2))/r
= (x^2 + y^2)/r - (xy - xy)/r - (x^2(x^2 + 2y^2) + y^2(x^2 + 2y^2))/r
= r - r(x^2 + 2y^2)
= r(1 - (r^2 cos^2 + 2r^2 sin^2))
= r(1 - r^2(1 + sin^2(theta)))

Since 1 + sin^2(theta) varies between 1 and 2, there is no single closed curve where r' = 0 for all theta. The system spirals inward for large r (specifically for r > 1) and outward for small r (r < 1/sqrt(2)), so it approaches a region between r = 1/sqrt(2) and r = 1. The omega-limit set is not a circle, but by a more careful analysis with Dulac's theorem using rho = 1/(x^2 + y^2), one can show div(rho*f) is sign-definite, ruling out closed orbits in the annular region and the entire plane (except near origin where the system spirals outward).

Problem 4: Bifurcation Analysis

For the system x' = mu*x - x^2 (one-dimensional), find all equilibria, determine their stability, and identify the type of bifurcation at mu = 0.

Show Solution

Equilibria: Set mu*x - x^2 = x(mu - x) = 0, so x = 0 or x = mu.

Stability: f(x) = mu*x - x^2, so f'(x) = mu - 2x.

  • At x = 0: f'(0) = mu. Stable if mu < 0, unstable if mu > 0.
  • At x = mu: f'(mu) = mu - 2mu = -mu. Stable if mu > 0, unstable if mu < 0.

Bifurcation: At mu = 0 the two equilibria x = 0 and x = mu merge. This is a transcritical bifurcation — two equilibria exchange stability as they pass through each other. Unlike a saddle-node, both equilibria persist on both sides of mu = 0, but their stability swaps.

The transcritical bifurcation is common in models where x = 0 is always an equilibrium (e.g., disease-free equilibrium in epidemiology), and mu crossing 0 corresponds to R_0 = 1 (the basic reproduction number threshold).

Problem 5: Hamiltonian System

The simple pendulum has Hamiltonian H(theta, p) = p^2/2 - cos(theta). Write Hamilton's equations, find all equilibria, classify them, and describe the phase portrait.

Show Solution

Hamilton's equations:

d(theta)/dt = partial H / partial p = p
dp/dt = -(partial H / partial theta) = -sin(theta)

This is the undamped pendulum: theta'' = -sin(theta).

Equilibria: theta' = p = 0 and p' = -sin(theta) = 0 give theta = k*pi for integer k, p = 0. The stable equilibria are at theta = 2k*pi (pendulum hanging down); the unstable equilibria are at theta = (2k+1)*pi (pendulum pointing up).

Classification: Jacobian Df = [[0, 1], [-cos(theta), 0]]. At theta = 0 (mod 2pi): eigenvalues are plus-or-minus i, a center (Lyapunov stable but not asymptotic, consistent with Hamiltonian nature — no damping, no attractor). At theta = pi (mod 2pi): eigenvalues are plus-or-minus 1, a saddle (unstable).

Phase portrait: Level curves H = c separate into three regimes: oscillations (|c| < 1, closed orbits around centers), the separatrix (c = 1, the homoclinic orbit connecting the saddle to itself), and rotation (|c| > 1, unbounded orbits corresponding to the pendulum spinning over the top). The separatrix forms the boundary between oscillation and rotation — it connects the unstable equilibria at theta = plus-or-minus pi.

12. Exam Tips and Common Mistakes

Key Strategies

  • Always check hyperbolicity before applying Hartman-Grobman. If any eigenvalue has zero real part, you need center manifold theory or a Lyapunov function.
  • For Lyapunov functions, try V = x^2 + y^2 first, then weighted versions V = ax^2 + by^2, then energies from the physics of the problem. If the system has a conserved quantity, that is a natural candidate.
  • In two dimensions, use the trace-determinant plane to classify equilibria without computing eigenvalues explicitly: T = trace(Df), D = det(Df) determines the type.
  • For bifurcation problems, draw the bifurcation diagram carefully: plot equilibria vs. parameter, indicate stability with solid/dashed curves, and label the bifurcation type at each special point.
  • To apply Poincare-Bendixson: construct an annular trapping region explicitly and verify the flow points inward on both boundaries (check the sign of f cdot n on each boundary).

Common Mistakes

  • Concluding an equilibrium is a center from the linearization alone. Centers are fragile under nonlinear perturbations — they can become spirals of either sign. Always confirm with a Lyapunov analysis.
  • Forgetting that positive Lyapunov exponents imply chaos only if the attractor is bounded. An unbounded trajectory can have positive exponents without being chaotic in the usual sense.
  • Confusing Lyapunov stability with asymptotic stability. Lyapunov stable means trajectories stay close; asymptotically stable means they also converge. Centers are Lyapunov stable but not asymptotically stable.
  • Applying Bendixson's criterion to a non-simply-connected domain. The criterion requires the region to be simply connected (no holes). For an annulus, use Dulac's theorem or check both boundaries explicitly.
  • Forgetting that Hamiltonian systems cannot have asymptotically stable equilibria or limit cycles — Liouville's theorem forbids volume contraction in the full phase space.

Quick Reference: Theorem Checklist

Stability

  • Lyapunov stability: V > 0, dV/dt leq 0
  • Asymptotic stability: V > 0, dV/dt < 0
  • Global: add radial unboundedness of V
  • LaSalle: extends to dV/dt semi-definite
  • Instability: Chetaev's theorem

Local Analysis

  • Hyperbolic: Hartman-Grobman (topological conjugacy)
  • Non-hyperbolic: center manifold theorem
  • No resonance: Poincare linearization
  • Normal forms: eliminate non-resonant terms

Periodic Orbits

  • Poincare-Bendixson: existence in planar systems
  • Bendixson/Dulac: non-existence criteria
  • Poincare map: stability via Floquet multipliers

Chaos and Hamiltonians

  • Chaos: positive max Lyapunov exponent
  • Feigenbaum: period-doubling cascade
  • Hamiltonian: Liouville, KAM, Poisson brackets
  • Integrable: n first integrals in involution

Related Topics