Systems of ODEs

The formula

The core object is a first-order system written in matrix form:

x'(t) = A·x(t) + B·u(t)     (state equation)
y(t)  = C·x(t) + D·u(t)     (output equation)

Reading it: x is the state vector — every quantity the system needs to remember to know where it's going next. A is the system matrix that couples the states to each other. B injects the inputs u, C picks out what your sensors actually see, and D is the direct feed-through from input to output.

For the free (unforced) case u = 0, the whole solution is one clean expression:

x(t) = e^(A·t)·x(0)

Reading it: the state at any time is the matrix exponential of A·t applied to the starting state. The matrix exponential is the same power series as the scalar one, with A in place of a number:

e^(A·t) = I + A·t + (A·t)²/2! + (A·t)³/3! + ...

Reading it: identity matrix plus the matrix, plus half its square, and so on — the series always converges, so this is a real, computable object, not a formal trick.

The behavior hides in the eigenvalues of A. Solve:

det(A − λ·I) = 0

Each eigenvalue λ = σ ± j·ω is one mode. σ (the real part) sets growth or decay; ω (the imaginary part) sets oscillation frequency. The poles of the transfer function are a subset of the eigenvalues of A; they coincide only when the realization is minimal — both controllable and observable. Pole-zero cancellations can hide unstable internal modes from the transfer function entirely ^[3], so the eigenvalues of A are the ground truth, not the poles you see from the outside.

Where you meet it

Any multi-DOF vibration model on the test stand. A single mass-spring-damper m·x'' + c·x' + k·x = 0 is second order, but you rewrite it as two first-order states — position and velocity — and it becomes x' = A·x. Bolt three masses together with cross-coupling springs and you have a 6-state system whose six eigenvalues are your three mode shapes and their damping. This is what the modal-analysis software is doing under the hood when it hands you natural frequencies.
Control-loop design at the design review. When the controls lead puts a state-space block on the screen and starts talking about pole placement, that block is x' = A·x + B·u. Moving the closed-loop eigenvalues to the left in the complex plane is the entire game — it's how you make a loop faster or better damped without touching the hardware.
RF and analog filter chains. A ladder of L and C stages, or a multi-pole active filter, is a coupled linear system. The state vector is the inductor currents and capacitor voltages. The eigenvalues tell you where the poles sit and whether the thing rings.
Coupled thermal or fluid nodes. Lump a structure into thermal masses connected by conductances and you get T' = A·T + (heat in). Same math as vibration, slower time constants. Every eigenvalue is a thermal time constant.

How it works

The trick that makes the whole subject work: any higher-order ODE becomes a first-order system if you add states. A fourth-order equation becomes four coupled first-order equations by naming the derivatives as new variables. Once everything is first order and linear, one matrix A carries all the coupling, and the eigenvalues of that matrix tell you everything about the free response.

Run a real case. Take m = 1, k = 100, c = 2, so the undamped natural frequency is ωn = √(k/m) = 10 rad/s and the damping ratio is ζ = c/(2·√(k·m)) = 0.1. Writing states as position and velocity:

A = [  0     1 ]
    [ -100  -2 ]

Its eigenvalues come out to λ = −1 ± 9.9499·j. Read them straight off: the real part −1 equals −ζ·ωn = −0.1·10, and the imaginary part 9.9499 equals the damped frequency ωn·√(1−ζ²) = 10·√(0.99). The matrix handed you the exact same numbers the textbook formulas give — that's the point. The eigenvalues are the physics.

Now the gotchas, because this is where people get burned:

e^(A+B) is not e^A·e^B unless A and B commute. For scalars you never worry about this. For matrices, A·B ≠ B·A in general, and the exponential product rule breaks. Verified numerically: two diagonal (commuting) matrices satisfy e^A·e^B = e^(A+B); two non-commuting ones do not, and the difference is large. Anytime you split a matrix exponential into a product to make life easier, first prove the pieces commute.
Don't compute the matrix exponential by literally truncating the series when A has widely spread eigenvalues. Fast modes and slow modes together make the problem stiff — the series terms swing huge before they settle, and floating point loses precision. Use a scaling-and-squaring routine (what expm in MATLAB, SciPy, and NumPy actually do) rather than a naive sum. This is the source of the classic paper title "Nineteen Dubious Ways to Compute the Exponential of a Matrix," and it is not academic — a stiff structural or thermal model will bite you.
Repeated eigenvalues break the clean eigenvector picture. When A has a repeated eigenvalue but not enough independent eigenvectors, you can't diagonalize it. The solution then carries t·e^(λ·t) terms — a mode that grows linearly before the exponential takes over. This is a defective matrix, and it shows up in critically damped systems and in some sensor-fusion models. The matrix exponential still works fine; the eigenvector shortcut doesn't.
Stability is a real-part question, not a magnitude question. For a continuous system x' = A·x, every eigenvalue must have a strictly negative real part for the response to decay. A simple (non-repeated) eigenvalue pair on the imaginary axis is marginal — bounded oscillation forever, or a constant offset if the eigenvalue is zero. A repeated defective eigenvalue on the axis is not marginal at all: the t·e^(λ·t) term from the previous bullet makes the response grow without bound — the free-free rigid-body mode, a double pole at the origin, is the classic case every structural test engineer meets. One eigenvalue in the right half-plane and the whole system diverges no matter how nice the others are. A single bad mode poisons everything. This is why controls people obsess over the rightmost pole.
Linearity is a validity limit. x' = A·x assumes A is constant and the physics is linear. Big deflections, saturating actuators, temperature-dependent stiffness — all of them make A a function of state, and the eigenvalue story only holds locally, around the operating point you linearized about. It's a superb local model and a liar far from that point.

History

The bones of this go back further than the control engineers who made it famous. In 1887 Giuseppe Peano — better known for his axioms of arithmetic — worked out a method for solving systems of linear differential equations by successive approximations, and published it the next year in Mathematische Annalen as "Intégration par séries des équations différentielles linéaires" ^[1]^[2]. Émile Picard hit on the same successive-approximation idea independently around the same time, which is why undergraduates meet it as Picard iteration. Peano's series is the direct ancestor of the matrix exponential: for constant A, his infinite series collapses to exactly e^(A·t).

Henry Frederick Baker carried the idea further in a 1905 paper in the Proceedings of the London Mathematical Society, tightening the theory of what is now called the Peano–Baker series and giving a fuller account of where it came from ^[6]. For decades this stayed a piece of pure analysis.

The engineering half of the story is Rudolf Kálmán. In a burst of papers around 1960 — the same year as his famous filtering paper, "A New Approach to Linear Filtering and Prediction Problems" — Kálmán reframed dynamical systems in terms of an internal state vector and the matrices A, B, C, D, and gave the field the notions of controllability and observability ^[3]^[4]. He didn't invent the matrix exponential; he took an old analytic tool and made it the native language of a machine age that finally had computers to run it. The reason every controls textbook now opens with x' = A·x + B·u is Kálmán deciding that was the right way to talk to a computer about a system.

Related tools

/tools/vibration-natural-freq — the single-DOF spring-mass case that becomes a 2-state system; the eigenvalues of A reproduce its natural frequency and damping
/tools/lc-resonance — an LC pair is a 2-state system whose eigenvalues sit on the imaginary axis (undamped oscillation)
/tools/frequency-period — the ω from an eigenvalue's imaginary part converts straight to the mode's period
/tools/heat-conduction — lump a conductive path into thermal nodes and you get T' = A·T, the same first-order form with slow eigenvalues

Sources

https://link.springer.com/article/10.1007/BF01443609
https://mathshistory.st-andrews.ac.uk/Biographies/Peano/
https://en.wikipedia.org/wiki/State-space_representation
https://en.wikipedia.org/wiki/Rudolf_E._K%C3%A1lm%C3%A1n
https://en.wikipedia.org/wiki/Matrix_exponential
https://arxiv.org/abs/1011.1775 — Baake & Schlägel, "The Peano-Baker series," Proc. Steklov Inst. Math. 275 (2011); documents Peano 1888 and H. F. Baker, "Note on the Integration of Linear Differential Equations," Proc. London Math. Soc. Ser. 2, Vol. 2 (1905), pp. 293–296

Written by HE in our own words from the cited sources — engineering judgment included, your stamp still required. All entries →