Generalized eigenvalue problem

The formula (the core equation(s), each with a one-line reading of what it says)

K·φ = λ·M·φ

Reads: find the special vectors φ where pushing against stiffness K and accelerating mass M produce the same shape — those shapes are the mode shapes, and each one comes with its own number λ.

λ = ω²,   ω = √λ  (rad/s),   f = √λ / (2π)  (Hz)

Reads: the eigenvalue is the natural frequency squared, in rad²/s² — two unit conversions stand between the solver output and the number you compare to the shaker sweep.

R(x) = (xᵀ·K·x) / (xᵀ·M·x)

Reads: the Rayleigh quotient — strain energy over kinetic energy for any trial shape x — and it never falls below the first eigenvalue, so any guessed shape gives an upper bound on the fundamental.

φᵢᵀ·M·φⱼ = 0  for i ≠ j,    φᵢᵀ·M·φᵢ = 1  (mass-normalized)

Reads: modes are orthogonal through the mass matrix, not to each other directly — this is what lets you split a coupled thousand-DOF model into independent single-DOF oscillators.

M = L·Lᵀ   ⇒   (L⁻¹·K·L⁻ᵀ)·y = λ·y,   y = Lᵀ·φ

Reads: when M is symmetric positive definite, a Cholesky factorization converts the pair (K, M) into an ordinary symmetric eigenproblem with the same eigenvalues — the theoretical bridge every textbook uses and production solvers mostly avoid.

K·ψ = λ·K_G·ψ

Reads: the same machinery with the geometric stiffness matrix in the M slot solves linear buckling — λ is now the load multiplier at which the structure loses stability, not a frequency.

Where you meet it (2-4 concrete engineering situations, specific: bench, test stand, review board)

The FEA modal run before a vibe test. Every NASTRAN SOL 103, every ANSYS or Abaqus modal solution, assembles K and M from the mesh and hands the pair to a generalized eigensolver. The frequencies and mode shapes that come back set the shaker profile, decide where the response accelerometers go, and tell you whether the fixture's first mode is going to sit inside the test band and contaminate the data.

Test-analysis correlation at the review board. The modal survey measured 87 Hz for the first bending mode; the model predicted 94. Before anyone signs the verification matrix, someone has to explain the delta — and the argument happens in the language of this equation: is the stiffness wrong (a joint modeled rigid that isn't) or is the mass wrong (cable harnesses nobody put in the model)? K and M are the only two knobs.

Spacecraft and launch-vehicle coupled loads. The payload's K and M matrices get reduced (Craig–Bampton), shipped to the launch provider, and merged into the vehicle model. The combined generalized eigenproblem produces the system modes that drive the loads cycle. A frequency requirement like "payload first lateral mode above 10 Hz" is a floor on the smallest λ of this problem.

Rotating machinery. Critical speeds of a shaft are the eigenvalues of a stiffness-mass pair; the whirl map on the test cell wall is this problem solved repeatedly across a speed range.

How it works (the real substance — behavior, gotchas, limits of validity, the mistake people make)

Run the smallest honest example. Two masses in a chain, m₁ = 2 kg driven wall-side, m₂ = 1 kg on the free end, two springs of k = 1000 N/m:

K = [ 2000  -1000 ]      M = [ 2  0 ]
    [-1000   1000 ]          [ 0  1 ]

Solving K·φ = λ·M·φ gives λ₁ = 292.9 and λ₂ = 1707.1 (rad/s)², so f₁ = 2.72 Hz and f₂ = 6.58 Hz. Mode 1 is both masses moving the same direction (the outer mass swinging √2 times as far); mode 2 is the masses opposed. Check the orthogonality: φ₁ᵀ·M·φ₂ = 0 to machine precision, while the plain dot product φ₁ᵀ·φ₂ is not zero. That is the point people miss — with unequal masses the mode shapes are not perpendicular as ordinary vectors. Any hand-rolled modal superposition that projects response onto modes with a plain dot product instead of an M-weighted one quietly mixes modes and gives wrong participation factors.

The symmetry of K and M is what buys the good behavior: real eigenvalues, real mode shapes, M-orthogonality, the Rayleigh bound. Add damping, gyroscopic effects, or follower forces and the problem stops being this one — you get complex eigenvalues and the quadratic eigenproblem (λ²M + λC + K)·φ = 0. Undamped modal analysis stays honest because structural damping is small; a squeeze-film damper or a spinning rotor is where the assumption dies.

Gotchas that cost real test time:

λ is ω², not ω, not f. The classic embarrassment is a factor of 2π or a squared frequency in a hand check. From λ to Hz is a square root and a divide by 2π, in that order.
Rigid-body modes. An unconstrained model — a spacecraft, an aircraft in free-free test configuration — has K singular and six eigenvalues at zero. Solvers handle this with a shift; a free-free run that reports six modes at 0.0001–0.02 Hz is healthy, and the first flexible mode is mode seven.
Singular M. Lumped-mass models often carry zero mass on rotational DOFs, which pushes eigenvalues to infinity. The symmetric shift-invert machinery tolerates it; naive reduction to M⁻¹·K does not, because M⁻¹ doesn't exist. This is also why nobody actually forms M⁻¹·K on big models — it destroys both symmetry and sparsity.
Missed modes. Iterative solvers extract a subset of eigenvalues and can skip one when two modes sit close together. Production codes run a Sturm sequence check — factor K − σ·M at a shift σ and count negative pivots, which equals the number of eigenvalues below σ. If the count says 12 and the solver returned 11, a mode escaped. When a modal survey finds a mode the model "doesn't have," check the extraction before rebuilding the mesh.
Normalization. Mode shapes are shapes; their amplitude is arbitrary. FEA codes default to mass normalization (φᵀMφ = 1), test-derived modes often come out unit-peak. Comparing them raw makes the correlation look worse than it is — that's what the MAC computation and consistent scaling are for.
Trust the bottom of the list. A mesh converges its low modes first. The 200th eigenvalue of a 10,000-DOF model is a property of the discretization, not the hardware. The working rule: extract to 1.5–2× the top of the frequency band you care about and treat everything above that as unconverged.

History (who derived it and when, told as a short story with inline [n] citations)

The algebra was settled before anyone had a mass matrix to put in it. Karl Weierstrass, working in Berlin's small circle of bilinear-form theorists, published the canonical theory of a nonsingular pair of forms in 1868 — his elementary divisors, read off the determinant of P + s·Q, classify exactly the regular pencils that a stiffness-mass pair forms ^[1]^[2]. Leopold Kronecker closed the remaining case, singular pencils, in 1890, and the combined Weierstrass–Kronecker canonical form still underlies how modern solvers reason about degenerate K, M pairs ^[3]^[4].

The physics ran in parallel. Lord Rayleigh's Theory of Sound (1877) worked fundamental frequencies of strings, bars, and plates by equating peak kinetic and potential energy over an assumed shape — the Rayleigh quotient, and with it the guarantee that a guessed mode shape overestimates the first frequency ^[5]^[6]. Every "energy method" hand check on a fixture design is Rayleigh's move, 150 years on.

Computation is the modern chapter. Cornelius Lanczos published his iteration at the National Bureau of Standards in 1950 ^[7]^[8]; in shift-inverted form it became the workhorse eigensolver inside NASTRAN-class codes. Cleve Moler and G. W. Stewart published the QZ algorithm in 1973, the first method to attack A·x = λ·B·x directly — with particular care for singular B — instead of reducing it to a standard problem first ^[9]^[10]. The same year, Klaus-Jürgen Bathe and Edward Wilson at Berkeley published their survey of eigensolution methods for K·φ = ω²·M·φ in structural mechanics, including the subspace iteration method Bathe had built for exactly the large sparse matrices FEA produces ^[11]^[12]. Between Lanczos, QZ, and subspace iteration, the overnight modal run became a lunch-break modal run.

Related tools (bullet list of HE calculator slugs that use or neighbor this topic, as /tools/ links; omit section if none fit)

/tools/vibration-natural-freq — the 1-DOF case, ω = √(k/m), which is this problem with 1×1 matrices
/tools/shaft-speed-critical — critical speeds are the eigenvalues of a rotor's stiffness-mass pair
/tools/euler-buckling — the buckling form K·ψ = λ·K_G·ψ, collapsed to a single column
/tools/spring-rate — where the entries of K come from
/tools/beam-deflection — the stiffness side of the pair for the beams your fixture is made of

Sources

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