When the Jacobian $Jf_p$ of a smooth system $\dot{X} = F(X)$ at an equilibrium $p$ has a pair of purely imaginary eigenvalues $\pm i\omega_0$ with $\omega_0 > 0$ — and all other eigenvalues have nonzero real part — the linearisation produces exact rotation at frequency $\omega_0$. It neither contracts nor expands, making the linear theory inconclusive about stability. Resolving this ambiguity requires examining the nonlinear terms through **complexification**, which isolates the single scalar quantity — the First Lyapunov Coefficient — that determines stability. When a parameter is introduced and the eigenvalues cross the imaginary axis, the same coefficient governs the **Hopf bifurcation**.
## Definition
[definition: Equilibrium With Purely Imaginary Eigenvalues]
Let $F: U \subseteq \mathbb{R}^n \to \mathbb{R}^n$ be a smooth vector field with an equilibrium $p$ (i.e., $F(p) = 0$). The equilibrium has **purely imaginary eigenvalues** if the Jacobian $A := Jf_p$ has a simple pair of eigenvalues $\lambda_{1,2} = \pm i\omega_0$ with $\omega_0 > 0$, and every other eigenvalue $\lambda_j$ satisfies $\operatorname{Re}(\lambda_j) \neq 0$.
[/definition]
The condition $\operatorname{Re}(\lambda_{1,2}) = 0$ places the equilibrium exactly on the stability [boundary](/page/Boundary): the imaginary axis for [continuous](/page/Continuity) systems plays the role of the unit circle for maps (see [Fixed Points of Maps](/page/Fixed%20Points%20of%20Maps)). The linear flow on the center eigenspace is rigid rotation — all orbits are circles of constant radius. Whether the full nonlinear orbits spiral inward (stable) or outward (unstable) depends entirely on higher-order terms.
When $n = 2$, the entire phase plane is the center eigenspace and no reduction is needed. When $n > 2$, the [Center Manifold Theorem](/page/Center%20Manifold%20Theorem) provides a 2D invariant surface on which the stability is decided, and the remaining directions are exponentially contracting or expanding.
## Complexification
The nonlinear analysis begins by writing the system in complex coordinates that diagonalise the linear part. In real coordinates, the linearisation is a $2 \times 2$ rotation matrix that mixes the two variables at every step. In complex coordinates, rotation becomes multiplication by $e^{i\omega_0 t}$, cleanly separating the rotational dynamics from the amplitude dynamics.
[definition: Complexification Of A Planar System]
Consider a smooth system on $\mathbb{R}^2$ with an equilibrium at the origin whose linearisation has eigenvalues $\pm i\omega_0$. After a real linear change of coordinates (if needed) to put the linear part in rotation form $\begin{pmatrix} 0 & -\omega_0 \\ \omega_0 & 0 \end{pmatrix}$, define the **complex coordinate**:
\begin{align*}
z := x_1 - ix_2, \qquad \bar{z} = x_1 + ix_2,
\end{align*}
with inverse $x_1 = \frac{z + \bar{z}}{2}$, $x_2 = \frac{z - \bar{z}}{-2i}$. In these coordinates the system becomes:
\begin{align*}
\dot{z} = i\omega_0 z + \sum_{j+k \ge 2} \frac{1}{j!\,k!}\,g_{jk}\,z^j\bar{z}^k,
\end{align*}
where the Taylor coefficients $g_{jk} \in \mathbb{C}$ are determined by the nonlinear terms of the original system.
[/definition]
The computation proceeds mechanically: [set](/page/Set) $z = x_1 - ix_2$, compute $\dot{z} = \dot{x}_1 - i\dot{x}_2$ by substituting the original ODEs, then re-express everything in terms of $z$ and $\bar{z}$ using the inverse relations. The result is matched to the expansion $\dot{z} = i\omega_0 z + \frac{1}{2}g_{20}z^2 + g_{11}z\bar{z} + \frac{1}{2}g_{02}\bar{z}^2 + \frac{1}{2}g_{21}z^2\bar{z} + \cdots$ to read off the coefficients $g_{20}, g_{11}, g_{02}, g_{21}$.
[example: Complexification Of A Planar System]
Consider the system:
\begin{align*}
\dot{x} &= y, \\
\dot{y} &= -x - 6x^2 + 2xy.
\end{align*}
**Step 1: Eigenvalues.** The Jacobian at the origin is $A = \begin{pmatrix} 0 & 1 \\ -1 & 0 \end{pmatrix}$ with eigenvalues $\pm i$, so $\omega_0 = 1$. The linear part is already in rotation form (with the convention that gives $\dot{z} = iz$ for $z = x - iy$).
**Step 2: Complex coordinate.** Set $z = x - iy$, so $x = \frac{z+\bar{z}}{2}$ and $y = i\frac{z - \bar{z}}{2}$. Compute $\dot{z} = \dot{x} - i\dot{y}$:
\begin{align*}
\dot{z} = y + ix + 6ix^2 - 2ixy.
\end{align*}
The linear part gives $y + ix = i(x - iy) = iz$. For the nonlinear terms, substitute the inverse relations:
\begin{align*}
6ix^2 &= \frac{6i}{4}(z + \bar{z})^2 = \frac{3i}{2}z^2 + 3iz\bar{z} + \frac{3i}{2}\bar{z}^2, \\
-2ixy &= -2i \cdot \frac{z+\bar{z}}{2} \cdot i\frac{z - \bar{z}}{2} = \frac{1}{2}z^2 - \frac{1}{2}\bar{z}^2.
\end{align*}
Combining all terms:
\begin{align*}
\dot{z} = iz + \left(\tfrac{1}{2} + \tfrac{3i}{2}\right)z^2 + 3iz\bar{z} + \left(-\tfrac{1}{2} + \tfrac{3i}{2}\right)\bar{z}^2.
\end{align*}
**Step 3: Read off coefficients.** Matching to $\dot{z} = iz + \frac{1}{2}g_{20}z^2 + g_{11}z\bar{z} + \frac{1}{2}g_{02}\bar{z}^2 + \cdots$:
\begin{align*}
g_{20} &= 1 + 3i, \qquad g_{11} = 3i, \qquad g_{02} = -1 + 3i, \qquad g_{21} = 0.
\end{align*}
There is no $z^2\bar{z}$ term in the expansion, so $g_{21} = 0$.
[/example]
## The First Lyapunov Coefficient
With the system in complex form, the stability reduces to a single number: the coefficient of the leading nonlinear term that systematically changes the amplitude of oscillation. The following motivation explains why this is the $z^2\bar{z}$ term and how the formula arises.
[motivation]
### Polar Coordinates and the Amplitude Equation
Substituting $z = re^{i\theta}$ into $\dot{z} = i\omega_0 z + \sum c_{jk}z^j\bar{z}^k$ and using $z^j\bar{z}^k = r^{j+k}e^{i(j-k)\theta}$, the amplitude equation becomes:
\begin{align*}
\dot{r} = \sum_{j+k \ge 2} r^{j+k}\,\operatorname{Re}\!\left(c_{jk}\,e^{i(j-k-1)\theta}\right).
\end{align*}
Each monomial contributes a term to $\dot{r}$ that oscillates in $\theta$ at frequency $j - k - 1$. To determine whether $r$ systematically grows or shrinks, we compute the net change in $r$ over one full rotation by integrating over $\theta$ from $0$ to $2\pi$:
\begin{align*}
\Delta r = \frac{1}{\omega_0}\int_0^{2\pi} \dot{r}\,d\theta + O(r^5).
\end{align*}
The change of integration variable from $t$ to $\theta$ uses $d\theta = \dot{\theta}\,dt \approx \omega_0\,dt$, with corrections of order $r^2$ that only affect higher-order terms.
### Resonant vs Non-Resonant Terms
Each [integral](/page/Integral) evaluates to:
\begin{align*}
\int_0^{2\pi} e^{in\theta}\,d\theta = \begin{cases} 2\pi & \text{if } n = 0, \\ 0 & \text{if } n \neq 0. \end{cases}
\end{align*}
A monomial $z^j\bar{z}^k$ produces $n = j - k - 1$ in the amplitude equation. When $n \neq 0$, the contribution oscillates and integrates to exactly zero — the term pushes $r$ up for part of the rotation and pulls it down equally for the rest, with no net effect. When $n = 0$ (i.e., $j = k + 1$), the contribution is constant in $\theta$ and produces a nonzero cumulative change. These are the **resonant terms**: the only monomials that systematically affect the amplitude.
### Why Cubic Is the Leading Order
At quadratic order ($j + k = 2$), the resonance condition $j = k + 1$ gives $j = 3/2$, which is not an integer — no quadratic monomial is resonant. At cubic order ($j + k = 3$), it gives $j = 2, k = 1$: the monomial $z^2\bar{z}$. This is the leading resonant term, and its coefficient determines:
\begin{align*}
\Delta r = \frac{2\pi}{\omega_0}\operatorname{Re}(c_{21})\,r_0^3 + O(r_0^5).
\end{align*}
### The Indirect Contribution of Quadratic Terms
Although quadratic terms have zero direct contribution to $\Delta r$, they distort the trajectory within each rotation. A near-identity coordinate change $z \mapsto z + h_{20}z^2 + h_{11}z\bar{z} + h_{02}\bar{z}^2$ eliminates them, but generates new cubic terms through interaction with the linear part. The $L_1$ formula accounts for both the direct resonant cubic term ($g_{21}$) and the indirect cubic contributions generated by eliminating the quadratics ($ig_{20}g_{11}$).
[/motivation]
[definition: Resonant Monomial]
A monomial $z^j\bar{z}^k$ in the complex expansion of a planar system with linear part $\dot{z} = i\omega_0 z$ is **resonant** if $j = k + 1$. Equivalently, its contribution to the amplitude equation $\dot{r}$ is independent of $\theta$ and does not integrate to zero over one rotation. The resonant monomials at successive orders are $z^2\bar{z}$ (cubic), $z^3\bar{z}^2$ (quintic), $z^4\bar{z}^3$ (septic), with corresponding Lyapunov coefficients $L_1, L_2, L_3$.
[/definition]
After performing successive near-identity coordinate changes to eliminate all non-resonant monomials, the system reduces to its **normal form**, in which only resonant terms survive:
\begin{align*}
\dot{z} = i\omega_0 z + c_1 z^2\bar{z} + c_2 z^3\bar{z}^2 + c_3 z^4\bar{z}^3 + \cdots
\end{align*}
where $c_k \in \mathbb{C}$. In polar coordinates $z = re^{i\theta}$, each resonant monomial $z^{k+1}\bar{z}^k = r^{2k+1}e^{i\theta}$ contributes purely to the amplitude, giving:
\begin{align*}
\dot{r} = \operatorname{Re}(c_1)\,r^3 + \operatorname{Re}(c_2)\,r^5 + \operatorname{Re}(c_3)\,r^7 + \cdots
\end{align*}
[definition: Lyapunov Coefficients]
The **$k$-th Lyapunov coefficient** is
\begin{align*}
L_k := \operatorname{Re}(c_k),
\end{align*}
where $c_k$ is the coefficient of the $k$-th resonant monomial $z^{k+1}\bar{z}^k$ in the normal form. The amplitude equation becomes:
\begin{align*}
\dot{r} = L_1 r^3 + L_2 r^5 + L_3 r^7 + \cdots
\end{align*}
The stability of the equilibrium is determined by the **first nonzero** $L_k$: for $r$ sufficiently small, that term dominates all higher powers of $r$, so $\operatorname{sgn}(\dot{r}) = \operatorname{sgn}(L_k)$.
[/definition]
The sign of the leading nonzero Lyapunov coefficient determines stability. The following result makes this precise and provides the formula for $L_1$.
[definition: First Lyapunov Coefficient]
Consider the complex expansion at the critical equilibrium:
\begin{align*}
\dot{z} = i\omega_0 z + \tfrac{1}{2}g_{20}\,z^2 + g_{11}\,z\bar{z} + \tfrac{1}{2}g_{02}\,\bar{z}^2 + \tfrac{1}{2}g_{21}\,z^2\bar{z} + \cdots
\end{align*}
The **First Lyapunov Coefficient** is:
\begin{align*}
L_1 := \frac{1}{2\omega_0}\,\operatorname{Re}\!\left(ig_{20}g_{11} + \omega_0\,g_{21}\right).
\end{align*}
[/definition]
[quotetheorem:926]
The formula $L_1 = \frac{1}{2\omega_0}\operatorname{Re}(ig_{20}g_{11} + \omega_0 g_{21})$ has two terms: $g_{21}$ is the direct contribution from the resonant cubic monomial $z^2\bar{z}$, and $ig_{20}g_{11}$ is the indirect contribution generated when the non-resonant quadratic terms are eliminated by a near-identity coordinate change. When $g_{21} = 0$ (no cubic resonant term in the original expansion), the stability is determined entirely by the interaction between the quadratic terms.
[remark: Convention Dependence]
The formula for $L_1$ depends on the normalisation convention for the Taylor coefficients $g_{jk}$. In the convention used here and in the [Hopf Bifurcation Theorem](/theorems/234), the expansion uses $\frac{1}{j!k!}g_{jk}$, so the coefficient of $z^2$ is $\frac{1}{2}g_{20}$, the coefficient of $z\bar{z}$ is $g_{11}$, and the coefficient of $z^2\bar{z}$ is $\frac{1}{2}g_{21}$. Different textbooks and lecture courses may absorb these factorial factors differently, leading to different-looking formulas that give the same sign of $L_1$.
[/remark]
[example: Computing The First Lyapunov Coefficient]
Continuing the example from the complexification section, we have $\omega_0 = 1$, $g_{20} = 1 + 3i$, $g_{11} = 3i$, and $g_{21} = 0$. Applying the formula:
\begin{align*}
ig_{20}g_{11} &= i(1 + 3i)(3i) = i(-9 + 3i) = -3 - 9i, \\
L_1 &= \frac{1}{2}\operatorname{Re}(-3 - 9i + 0) = \frac{1}{2}(-3) = -\frac{3}{2}.
\end{align*}
Since $L_1 < 0$, the origin is asymptotically stable — all sufficiently small perturbations spiral back to the equilibrium.
[/example]
### Hamiltonian and Reversible Systems
If the system admits a conserved quantity $H$ with $\dot{H} = 0$ (Hamiltonian) or a time-reversal symmetry $(x, t) \mapsto (Sx, -t)$ (reversible), then every nearby orbit is periodic and the equilibrium is a true nonlinear center. In such cases, **all Lyapunov coefficients vanish**: $L_1 = L_2 = \cdots = 0$. This can be concluded without computation — if the system preserves energy or has a reversibility, no systematic amplitude change is possible. Recognising a Hamiltonian or reversible structure is a valuable shortcut on exams: it immediately gives $L_k = 0$ for all $k$.
## Center Manifold Reduction
When the system lives in $\mathbb{R}^n$ with $n > 2$, the purely imaginary pair $\pm i\omega_0$ spans a 2D center eigenspace, and the remaining eigenvalues (with nonzero real part) span the stable and unstable eigenspaces. The center manifold theorem guarantees that the stability analysis reduces to the 2D center manifold. The following result makes this reduction precise for the Hopf setting.
[quotetheorem:233]
This theorem says that no matter how high-dimensional the original system is, the birth of periodic orbits and the stability at purely imaginary eigenvalues is a **two-dimensional phenomenon**. The remaining directions contract or expand exponentially and do not participate in the bifurcation. All that matters are the Taylor coefficients $g_{jk}$ of the restricted system on the center manifold.
In practice, the reduction proceeds as follows. Decompose the state as $(w, \bar{w}, v)$ where $w, \bar{w}$ are complex center coordinates and $v$ collects the non-center variables. The center manifold is the graph $v = h(w, \bar{w})$ with $h(0) = 0$ and $Dh(0) = 0$. To find $h$, either solve the invariance equation (see [Center Manifold Theorem](/page/Center%20Manifold%20Theorem)) or use the **normal form transformation shortcut**: choose a coordinate change in the non-center variables that eliminates the center-variable terms from the non-center equations up to the required order, making the center manifold simply $v = 0$ in the new coordinates. After substituting $v = h(w, \bar{w})$ into the center equations, the complexification and $L_1$ computation proceed exactly as in the 2D case.
## The Hopf Bifurcation
Everything above applies to a **single system** at a fixed configuration where the eigenvalues are exactly $\pm i\omega_0$. The Lyapunov coefficient $L_1$ determines whether the equilibrium is stable or unstable at this critical moment.
A fundamentally richer question arises when a **parameter** $\mu$ is introduced and the eigenvalues become $\alpha(\mu) \pm i\omega(\mu)$ with $\alpha(0) = 0$. As $\mu$ varies, the real part $\alpha(\mu)$ crosses zero — the eigenvalues traverse the imaginary axis — and the equilibrium transitions between stable and unstable. The [Hopf bifurcation](/page/Hopf%20Bifurcation) describes what is born at this transition: a periodic orbit ([limit](/page/Limit) cycle) branches from the equilibrium, and its stability and direction are governed by the same $L_1$.
[quotetheorem:234]
The computation of $L_1$ is identical in both settings — the parameter $\mu$ only affects the **interpretation**:
| Setting | $L_1$ answers | $L_1 < 0$ means | $L_1 > 0$ means |
|---|---|---|---|
| No parameter (fixed system) | Is the equilibrium stable? | Asymptotically stable | Unstable |
| With parameter (Hopf) | What type of bifurcation? | Supercritical: stable limit cycle born | Subcritical: unstable limit cycle born |
The transversality condition $\alpha'(0) \neq 0$ ensures that the eigenvalues genuinely cross the imaginary axis rather than tangentially touching it. Combined with $L_1 \neq 0$ (non-degeneracy), it guarantees the existence and uniqueness of the bifurcating periodic orbit family.
This parallels the discrete-time story developed in [Fixed Points of Maps](/page/Fixed%20Points%20of%20Maps): for maps, the First Lyapunov Coefficient at the unit circle determines stability without parameters, and determines the type of [Neimark-Sacker bifurcation](/page/Neimark-Sacker%20Bifurcation) (supercritical vs subcritical) with parameters. The imaginary axis for ODEs corresponds to the unit circle for maps, and the limit cycle corresponds to the invariant closed curve.
## References
1. Y. A. Kuznetsov, *Elements of Applied Bifurcation Theory*, 3rd ed., Springer (2004). §3.5, §5.4.
2. J. Guckenheimer and P. Holmes, *Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields*, Springer (1983). §3.4.
3. S. H. Strogatz, *Nonlinear Dynamics and Chaos*, 2nd ed., Westview Press (2015). §8.2–8.3.