Stiffness and mass matrices

The formula

K·u = f

The static problem. K is the global stiffness matrix, u is the vector of nodal displacements, f is the vector of applied nodal loads. Every linear static FEA run on the planet is a factorization of this system.

M·ü + K·u = 0   →   K·φ = ω²·M·φ

The undamped dynamics and the generalized eigenproblem it produces. Each eigenpair is a natural frequency ω and a mode shape φ. This is what a modal solution (NASTRAN SOL 103 and its cousins) actually computes.

k_axial = (E·A/L) · [ 1  −1 ]
                    [ −1  1 ]

The simplest element stiffness there is: a two-node axial rod. E·A/L is just a spring rate — a steel rod with E = 200 GPa, A = 100 mm², L = 1 m works out to 20 MN/m. Every fancier element (beam, shell, brick) is the same idea with more rows.

m_consistent = (ρ·A·L/6) · [ 2  1 ]        m_lumped = (ρ·A·L/2) · [ 1  0 ]
                           [ 1  2 ]                               [ 0  1 ]

Two ways to spread the same rod's mass onto its nodes. Consistent uses the element's own shape functions; lumped just splits the mass in half and drops it on the ends. Both conserve total mass; they disagree about the off-diagonal coupling, and that disagreement shows up in your frequencies.

K = Σₑ Lₑᵀ·kₑ·Lₑ

Assembly. Each element matrix gets scattered into the global matrix at the rows and columns of the degrees of freedom it touches, and overlapping entries add. Shared nodes are where element contributions stack — that summation is the model.

Where you meet it

Pre-test at a vibration lab. Before the shaker ever spins up, someone runs K·φ = ω²·M·φ on the test-article model to place accelerometers where the modes have motion and to predict which frequencies the fixture will pollute. If the model's first mode is 10% off the tap test, the argument that follows is about M and K, not about the lab.
A stalled solve at 4:45 p.m. The job that ran fine yesterday now dies with "matrix is singular" or crawls with pivot warnings. That is K talking. A free-free model, a forgotten constraint, an unmerged coincident node, or a mechanism in the joints all show up as zero (or near-zero) pivots during factorization.
Model correlation at a review board. Flight rationale for a launch vehicle bracket rests on a model whose mass matrix was tuned so total weight and CG match the scale and the as-built drawing. The board doesn't ask about shape functions; it asks why the model weighs 3 lb more than the hardware. That delta lives in M.
Component mode synthesis across contractors. When one company builds the avionics shelf model and another builds the structure it bolts to, what actually crosses the interface is a reduced K and M (Craig–Bampton matrices). Two teams who have never met integrate their hardware by mailing each other matrices.

How it works

Both matrices come from the same recipe. Pick shape functions that interpolate displacement between an element's nodes. Integrate strain energy against them and you get kₑ; integrate kinetic energy against them and you get mₑ. Assemble both the same way. Because they descend from energy integrals, K and M are symmetric, M is positive definite, and K is positive semi-definite — the "semi" being the rigid-body modes.

The zero-energy modes are a feature, then a bug. An unconstrained 3D structure can translate three ways and rotate three ways without straining, so a free-free K has exactly six zero eigenvalues. A spacecraft model should have exactly six; a fixed-base test stand model should have zero. Seven or more zero (or nearly zero) eigenvalues in a free-free check run means a mechanism — a loose fastener DOF, an unattached mass element, a pin joint nobody meant to model. The free-free eigenvalue check is the cheapest model-quality gate in the business, and skipping it is how mechanisms reach the review board.

Conditioning is a mesh property and a modeling property. For a second-order problem the condition number of K grows like 1/h² as you refine: halve the element size and cond(K) roughly quadruples. A uniform 1D bar shows it cleanly — 10, 20, 40, 80 elements give condition numbers of about 40, 161, 648, 2593. Refinement alone rarely kills a solve; stiffness contrast does. Bolt a 1e10 penalty spring onto a structure whose members run around 1e6 and you have donated four orders of magnitude of condition number for free. The rule of thumb: a direct solver loses roughly log₁₀(cond(K)) digits, and double precision only has about sixteen. Rigid elements exist so you don't have to fake rigidity with big numbers — use them.

Consistent vs. lumped mass is a real choice, not a checkbox. Archer proved that consistent-mass frequencies bound the exact answer from above ^[4]; lumped mass usually errs low. A one-element axial rod, fixed-free, tells the story: the exact first-frequency coefficient is π/2 ≈ 1.571, consistent mass gives √3 ≈ 1.732 (10% high), lumped gives √2 ≈ 1.414 (10% low). A one-element cantilever beam with consistent mass gets the exact 3.516 coefficient to within half a percent (3.533). Refine the mesh and both converge — but explicit dynamics codes lump anyway, because a diagonal M inverts for free every time step, and that trade is the entire reason crash and pyroshock codes run at all.

The mistakes people actually make.

Diagnosing a singular K as a solver problem. The solver is telling you your structure can move without resisting. Fix the constraints or the connectivity, not the tolerance settings.
Grounding a free structure with stiff springs "to make it solve." Soft seismic springs (a few orders below the structure's stiffness) are legitimate for free-free checks; springs stiff enough to look like constraints quietly shift your low modes and wreck conditioning at the same time.
Trusting a mass matrix nobody audited. Nonstructural mass, paint, harness, and fastener mass get smeared in by hand, and units errors in density (kg/mm³ vs kg/m³ is a factor of 10⁹) sail through assembly without complaint. Check total mass, CG, and inertias against hardware before you believe a single frequency.
Assuming symmetry survives everything. Follower forces, spinning-rotor gyroscopics, and some contact formulations produce unsymmetric terms. The moment K stops being symmetric, half your solver guarantees and most of your intuition about real eigenvalues go with it.

History

The matrix view of structures came out of the aircraft industry, because swept jet wings in the early 1950s had outrun every hand method. Between October 1954 and May 1955, John Argyris serialized "Energy Theorems and Structural Analysis" in the journal Aircraft Engineering, recasting structural analysis as matrix operations organized around energy principles — the force method and displacement method as two faces of one algebra ^[2]^[3]. In September 1956, M. J. Turner of Boeing, with Ray Clough, H. C. Martin, and L. J. Topp, published "Stiffness and Deflection Analysis of Complex Structures" in the Journal of the Aeronautical Sciences — the paper that laid out the direct stiffness method and derived the stiffness matrix for a triangular plate element, so that a continuous wing skin could be assembled from pieces ^[1]^[2]. Clough supplied the name in 1960, calling the building blocks "finite elements" in a plane-stress paper at the ASCE electronic computation conference in Pittsburgh ^[2]^[8].

The mass matrix got its modern form in 1963, when John Archer's "Consistent Mass Matrix for Distributed Mass Systems" in the ASCE Journal of the Structural Division showed how to build M from the same shape functions as K — and proved the resulting frequencies were upper bounds on the true ones ^[4]^[5]. The whole apparatus then became infrastructure: after a 1964 NASA review found its centers each writing their own one-off structural codes, the agency commissioned a single general-purpose program, demonstrated it at Goddard in May 1967, and released NASTRAN in 1968 ^[6]^[7]. Six decades on, the matrices those aircraft engineers defined are still what the solver factorizes.

Related tools

/tools/vibration-natural-freq — the single-DOF ω = √(k/m) that K and M generalize
/tools/spring-rate — every stiffness matrix entry is ultimately a spring rate
/tools/beam-deflection — closed-form answers the beam element stiffness reproduces
/tools/shaft-speed-critical — a rotor eigenproblem in miniature
/tools/moment-of-inertia-shapes — the I inside every beam element's EI/L³
/tools/section-modulus-rect — section properties feeding element stiffness
/tools/euler-buckling — the geometric-stiffness eigenproblem next door
/tools/stress-strain — the E that scales the whole stiffness matrix

Sources

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