1. Why systems show up
There are two roads that lead straight to systems of first-order equations, and they meet at the same destination.
Road one: higher-order ODEs in disguise. Any nth-order ODE can be flattened into $n$ first-order equations by introducing the lower derivatives as new variables. Take a damped spring:
$$ m\,x'' + c\,x' + k\,x = 0. $$Let $x_1 = x$ (position) and $x_2 = x'$ (velocity). Then
$$ \begin{aligned} x_1' &= x_2, \\ x_2' &= -\tfrac{k}{m}\,x_1 - \tfrac{c}{m}\,x_2. \end{aligned} $$The same content, expressed as two coupled first-order equations. This isn't a trick — it's the canonical form, and almost every numerical ODE solver wants its input in exactly this shape.
Road two: genuinely coupled processes. Many real-world systems are natively coupled: a chemical reaction with three concentrations changing in tandem, an ecosystem where predators eat prey while themselves being eaten, two tanks of brine swapping fluid. There's no scalar ODE waiting underneath — the system is the model.
Whichever road you arrive by, the same mathematical machinery applies: a state vector, a vector field, and trajectories that trace out the system's life in phase space.
2. Vector form: $\mathbf{x}' = \mathbf{f}(\mathbf{x})$
Collect the unknowns into a single column vector $\mathbf{x}(t) \in \mathbb{R}^n$. The dynamics are then a single equation:
$$ \mathbf{x}'(t) = \mathbf{f}(\mathbf{x}(t), t). $$If $\mathbf{f}$ doesn't depend on $t$ explicitly, the system is called autonomous, and the right side is a static vector field: at every point $\mathbf{x}$ in space, the rule prescribes a velocity $\mathbf{f}(\mathbf{x})$ that the state must obey.
$\mathbf{x}(t)$ is the state — everything you need to know about the system at time $t$. For the spring, that's position and velocity. The map $\mathbf{x} \mapsto \mathbf{f}(\mathbf{x})$ is the vector field — at each state, it tells the state how to move next.
The picture to hold in your head: imagine an arrow drawn at every point of the plane. A trajectory is what you get if you drop a particle at some starting point and let it ride the arrows. Solving the ODE means tracing that path.
3. Linear systems $\mathbf{x}' = A\mathbf{x}$
The cleanest systems are linear with constant coefficients:
$$ \mathbf{x}'(t) = A\,\mathbf{x}(t), \qquad A \in \mathbb{R}^{n \times n}. $$By analogy with the scalar equation $x' = a x$, whose solution is $x(t) = e^{at}\,x(0)$, the system's solution is
$$ \mathbf{x}(t) = e^{A t}\,\mathbf{x}(0), $$where $e^{At}$ is the matrix exponential, defined by the same Taylor series as the scalar one:
$$ e^{At} = I + At + \tfrac{1}{2!}(At)^2 + \tfrac{1}{3!}(At)^3 + \cdots $$That definition is conceptually satisfying but practically painful — you don't want to sum an infinite series of matrices. The shortcut is to diagonalize $A$. If $A\mathbf{v}_k = \lambda_k \mathbf{v}_k$ for $k = 1, \ldots, n$, then each eigenvector points along a direction in which the system behaves like a scalar exponential:
$$ \mathbf{x}(t) = c_1 e^{\lambda_1 t}\,\mathbf{v}_1 + c_2 e^{\lambda_2 t}\,\mathbf{v}_2 + \cdots + c_n e^{\lambda_n t}\,\mathbf{v}_n, $$with the constants $c_k$ determined by the initial condition $\mathbf{x}(0)$.
Eigenvectors are the directions along which the dynamics decouple into independent one-dimensional exponentials. Every other trajectory is a linear combination of those clean ones.
Complex eigenvalues — oscillation
If $A$ is real but $\lambda = \alpha + i\beta$ is complex, the eigenvector is too. The two complex solutions combine into two real ones via Euler's formula, giving an oscillatory factor:
$$ e^{(\alpha + i\beta)t} = e^{\alpha t}\,(\cos\beta t + i\sin\beta t). $$The real part $\alpha$ governs whether the oscillation grows ($\alpha > 0$), shrinks ($\alpha < 0$), or holds steady ($\alpha = 0$); the imaginary part $\beta$ is the angular frequency.
4. Eigenvalues, equilibria, and the zoo of 2D portraits
Every linear system $\mathbf{x}' = A\mathbf{x}$ has at least one equilibrium — the origin — where $\mathbf{x}' = \mathbf{0}$. The local behavior near the origin is determined entirely by the eigenvalues of $A$. For 2×2 systems, the classification is finite and worth memorizing.
| Eigenvalues of $A$ | Equilibrium type | Stability |
|---|---|---|
| Real, both negative | Stable node | Stable (all trajectories flow in) |
| Real, both positive | Unstable node | Unstable (all flow out) |
| Real, opposite signs | Saddle | Unstable (one direction in, one out) |
| Complex, $\operatorname{Re}\lambda < 0$ | Stable spiral | Stable (inward spiral) |
| Complex, $\operatorname{Re}\lambda > 0$ | Unstable spiral | Unstable (outward spiral) |
| Pure imaginary, $\operatorname{Re}\lambda = 0$ | Center | Marginally stable (closed orbits) |
You can read the eigenvalue regime straight off the trace and determinant of $A$. Since $\lambda^2 - (\operatorname{tr} A)\lambda + \det A = 0$, the discriminant $\Delta = (\operatorname{tr} A)^2 - 4 \det A$ splits the trace-determinant plane into the regions above.
5. The phase plane: trajectories, nullclines, direction fields
For a 2D system $(x_1, x_2)$, the phase plane is the $(x_1, x_2)$ plane viewed without any explicit time axis. A trajectory is the curve traced out as $t$ runs from $-\infty$ to $+\infty$. Time is implicit — you can recover it by watching which way the arrows point.
Three tools dominate the analysis.
Direction field
At each point $(x_1, x_2)$, draw a small arrow in the direction of $\mathbf{f}(x_1, x_2)$. A trajectory is any curve that's everywhere tangent to those arrows. You can sketch one by eye from a few sample arrows.
Nullclines
The $x_1$-nullcline is the set of points where $x_1' = 0$ — the state is momentarily not changing in the $x_1$ direction, so trajectories crossing it have purely vertical velocity. The $x_2$-nullcline is the analog for $x_2' = 0$. Equilibria are exactly the intersections of all the nullclines.
Trajectories
The curve $\mathbf{x}(t)$ as $t$ varies. Trajectories never cross each other in an autonomous system — uniqueness forbids it. If two paths met at a point, the vector field would have to pick a direction, but it can only pick one.
A time-series plot has $t$ on the horizontal axis and $x_1$ (or $x_2$) on the vertical. A phase-plane plot has $x_1$ horizontally and $x_2$ vertically — no $t$. Same dynamics, two completely different pictures. Both are useful; confusing them is a common early mistake.
6. Equilibria and linearization for nonlinear systems
Most interesting systems are nonlinear: $\mathbf{x}' = \mathbf{f}(\mathbf{x})$ with $\mathbf{f}$ not of the form $A\mathbf{x}$. There's usually no closed-form solution. But you can still understand the system qualitatively by finding its equilibria and zooming in.
An equilibrium $\mathbf{x}^*$ is any point where $\mathbf{f}(\mathbf{x}^*) = \mathbf{0}$. Near it, write $\mathbf{x} = \mathbf{x}^* + \mathbf{u}$ for small perturbation $\mathbf{u}$, and Taylor-expand:
$$ \mathbf{u}' \approx J(\mathbf{x}^*)\,\mathbf{u}, \quad \text{where } J = \frac{\partial \mathbf{f}}{\partial \mathbf{x}} $$is the Jacobian matrix evaluated at $\mathbf{x}^*$. Now apply the 2D classification table to $J(\mathbf{x}^*)$: the type of $J$'s eigenvalues tells you the type of the equilibrium locally. (One caveat: this fails at marginal cases where $J$ has eigenvalues with zero real part — nonlinear terms can tip a linear center into a stable or unstable spiral.)
The whole nonlinear story is approximately the union of linear stories told around each equilibrium.
7. Lotka-Volterra: predators, prey, and a closed orbit
A classical model from mathematical biology. Let $x$ be a prey population (rabbits) and $y$ a predator population (foxes):
$$ \begin{aligned} x' &= \alpha x - \beta xy, \\ y' &= -\gamma y + \delta xy, \end{aligned} $$with positive constants $\alpha, \beta, \gamma, \delta$. The story each term tells:
- $\alpha x$: rabbits reproduce in proportion to their number.
- $-\beta xy$: rabbits get eaten at a rate proportional to the chance of a rabbit-fox encounter.
- $-\gamma y$: foxes die off naturally.
- $+\delta xy$: foxes reproduce when they catch rabbits.
Equilibria: $\mathbf{f} = \mathbf{0}$ gives $(x, y) = (0, 0)$ and $(x, y) = (\gamma/\delta,\ \alpha/\beta)$. The Jacobian at the nontrivial equilibrium has purely imaginary eigenvalues — a center. So trajectories near that point are closed loops: the populations oscillate forever, predators lagging prey by a quarter cycle. It's not a coincidence that real ecological data shows this lag — Lotka and Volterra wrote down the model because Volterra's son-in-law (a marine biologist) had noticed Adriatic shark populations cycling out of phase with their fish prey.
8. Common pitfalls
A trajectory in the phase plane is the path $(x_1(t), x_2(t))$ traces out — it has no $t$ axis. A graph of $x_1(t)$ has $t$ on the horizontal axis. They're different objects. A spiral in the phase plane corresponds to an oscillation with decaying amplitude in the time series, not a spiral-shaped time series.
For nonlinear systems, the linearization at a marginal equilibrium (purely imaginary eigenvalues) can be misleading. A linear center can become a slowly drifting spiral once higher-order terms kick in. Always check whether the nonlinear terms preserve the structure.
If "eigenvalue" and "eigenvector" aren't fluent, this whole chapter is brick wall. The system stuff is genuinely just $\mathbf{x}' = A\mathbf{x}$ with linear-algebra glue. Spend an hour on the eigenstuff if you have to — it pays for itself many times over.
In an autonomous system, distinct trajectories never intersect. If your sketch shows two trajectories crossing, you've made an error — either the vector field is wrong or you've miscounted equilibria. (Trajectories can approach the same equilibrium asymptotically, but they never touch it in finite time and never cross each other elsewhere.)
9. Worked examples
Example 1 · Real distinct eigenvalues (saddle)
Solve $\mathbf{x}' = A\mathbf{x}$ with $A = \begin{pmatrix} 1 & 2 \\ 2 & 1 \end{pmatrix}$.
Step 1. Characteristic equation: $\det(A - \lambda I) = (1-\lambda)^2 - 4 = 0$, so $\lambda = 3$ or $\lambda = -1$.
Step 2. Eigenvectors. For $\lambda = 3$: $(A - 3I)\mathbf{v} = 0$ gives $\mathbf{v}_1 = (1, 1)^T$. For $\lambda = -1$: $\mathbf{v}_2 = (1, -1)^T$.
Step 3. General solution:
$$ \mathbf{x}(t) = c_1 e^{3t} \begin{pmatrix} 1 \\ 1 \end{pmatrix} + c_2 e^{-t} \begin{pmatrix} 1 \\ -1 \end{pmatrix}. $$One eigenvalue is positive and one is negative — the origin is a saddle. Trajectories race outward along the $(1,1)$ direction and decay inward along $(1,-1)$.
Example 2 · Complex eigenvalues (stable spiral)
Solve $\mathbf{x}' = A\mathbf{x}$ with $A = \begin{pmatrix} -1 & -2 \\ 2 & -1 \end{pmatrix}$.
Step 1. Characteristic equation: $(\lambda + 1)^2 + 4 = 0$, so $\lambda = -1 \pm 2i$.
Step 2. The real part is $-1$ (decay), the imaginary part is $2$ (angular frequency). Real general solution:
$$ \mathbf{x}(t) = e^{-t}\!\left[\,c_1\!\begin{pmatrix} \cos 2t \\ -\sin 2t \end{pmatrix} + c_2\!\begin{pmatrix} \sin 2t \\ \cos 2t \end{pmatrix}\right]\!. $$Trajectories spiral inward toward the origin, completing one full turn every $\pi$ time units while shrinking by factor $e^{-\pi}$.
Example 3 · Finding equilibria of a nonlinear system
Find the equilibria of $x' = x - x^2 - xy$, $y' = 2y - 2y^2 - xy$.
Step 1. Set both right-hand sides to zero. Factor:
$$ x(1 - x - y) = 0, \qquad y(2 - 2y - x) = 0. $$Step 2. Case analysis. Either $x = 0$ or $1 - x - y = 0$; either $y = 0$ or $2 - 2y - x = 0$. The four combinations:
- $x = 0, y = 0$: equilibrium $(0,0)$.
- $x = 0, 2 - 2y - x = 0$: gives $(0, 1)$.
- $1 - x - y = 0, y = 0$: gives $(1, 0)$.
- $1 - x - y = 0, 2 - 2y - x = 0$: solve simultaneously to get $(0, 1)$ again, so only three distinct equilibria.
Step 3. The Jacobian is $J = \begin{pmatrix} 1 - 2x - y & -x \\ -y & 2 - 4y - x \end{pmatrix}$. Evaluate at each equilibrium and classify by its eigenvalues.
Example 4 · Damped spring as a 2D system
Take $x'' + x' + x = 0$. Let $x_1 = x$, $x_2 = x'$. Then
$$ \begin{pmatrix} x_1' \\ x_2' \end{pmatrix} = \begin{pmatrix} 0 & 1 \\ -1 & -1 \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \end{pmatrix}. $$Eigenvalues of the matrix: $\lambda^2 + \lambda + 1 = 0$, so $\lambda = -\tfrac{1}{2} \pm \tfrac{\sqrt{3}}{2} i$. Real part negative, imaginary part nonzero — the origin is a stable spiral, exactly matching the under-damped oscillation you'd expect from a spring with friction.