Finite element method

The formula

The whole method collapses to one linear system:

K · u = f

Read it as: the global stiffness matrix K times the vector of unknown nodal displacements u equals the vector of applied nodal loads f. Solve for u and you know how every node in the model moved. Everything else — stress, strain, temperature, mode shapes — is post-processing off u.

K is built by adding up element stiffness matrices. For a single axial bar element of length L, cross-section A, modulus E:

k_element = (A·E / L) · [ [ 1, -1 ],
                          [ -1,  1 ] ]

Reading: pull the two ends apart and the element resists with stiffness A·E/L — the same spring constant you'd write for a bar by hand, now packaged as a 2×2 matrix so it can be assembled with its neighbors.

The theory underneath comes from the weak (variational) form. Instead of forcing the governing PDE to hold at every point, you require it to hold on average against a set of test functions v:

a(u, v) = ℓ(v)   for all admissible v

Reading: a(·,·) is the strain-energy bilinear form (the stiffness), ℓ(·) is the work done by the loads. Pick approximate u and v from the same finite set of simple shape functions N_i(x) — one per node, equal to 1 at its own node and 0 at all others — and the weak form turns directly into K·u = f, with:

K_ij = a(N_i, N_j)        f_i = ℓ(N_i)

That is the Galerkin recipe, and it is the engine inside every solver you submit jobs to.

Where you meet it

The structural review board. Someone puts a von Mises contour plot of a bracket on the screen with a red hot-spot at a fillet. That plot is u solved, differentiated into strain, multiplied by the material matrix into stress. The first honest question in the room is "how fine is the mesh at that fillet," because the red number moves as the elements shrink.
The vibration bench at Redstone. Before a box goes on the shaker table, someone ran a modal FEA to predict natural frequencies. Same math, different right-hand side: instead of K·u = f you solve the eigenproblem K·φ = ω²·M·φ for mode shapes φ and frequencies ω. The test either confirms the model or exposes a boundary condition you faked.
Thermal analysis on an avionics chassis. Swap displacement for temperature, stiffness for conductivity, force for heat load. Identical K·u = f. FEM does not care whether u is inches or degrees.
The fatigue sign-off. The stress life calc downstream is only as good as the FEA stress it's fed. If the analyst reported peak stress at an unconverged singularity, the fatigue number is fiction with three decimal places.

How it works

The trade is exact-equation-you-can't-solve for approximate-equation-you-can. A real part obeys a partial differential equation over a continuous body — infinite unknowns. FEM replaces that with a finite number of nodal unknowns connected by simple element shapes, and the answer converges to the true solution as the elements get smaller. That's h-refinement: same element type, finer mesh. The alternative is p-refinement: keep the mesh, raise the polynomial order inside each element. Do both and you have the hp-method.

The mistakes are almost never in the solver. They're in what the analyst told it.

Singularities. A sharp re-entrant corner, a point load, a knife-edge boundary condition — the true stress there is infinite, so the FEA stress climbs every time you refine and never converges. Engineers chase that red spot with smaller elements and report a bigger number each pass, not realizing they're refining toward infinity. The fix is a fillet in the model or an understanding that the singularity isn't real.
No convergence study. One mesh is a guess. You need at least two, ideally three, showing the quantity of interest settling down. If halving the element size still moves your peak stress 20%, you are not converged and the number is not yours to report.
Element quality. Long, skinny, or badly warped elements make K ill-conditioned and the answer garbage. Solvers flag high aspect ratios and Jacobian warnings for a reason.
Boundary conditions. Over-constrain a model and you invent stiffness that carries load the real part can't. Fully fixing a bolt hole that's actually a slip-fit is the classic. The BCs are where the physics lives; the mesh is just bookkeeping.
Linear vs. nonlinear. K·u = f assumes stiffness doesn't change as the part deflects. Snap-through, contact, large rotation, or plasticity break that assumption. Then K depends on u and the solver iterates — and it can converge to the wrong branch or not at all.

The load path insight is worth stating plainly: FEM is exact for the model you built and only as good as your judgment about what the model leaves out. It computes; it does not understand.

History

The method has two family trees that didn't know they were related for fifteen years.

The math branch starts with Richard Courant, who in a 1941 lecture to the American Mathematical Society — published in 1943 — solved a torsion problem by splitting the cross-section into triangles and applying a Rayleigh–Ritz variational method with trial functions on each triangle ^[1]^[2]. That is the finite element method in embryo, though nobody built machines to run it. The same year, 1941, the Russian-Canadian structural engineer Alexander Hrennikoff published a lattice-framework model that replaced a continuous membrane or plate with a mesh of bars — the other early instance of discretizing a continuum into pieces ^[1]^[3].

The engineering branch grew out of the aircraft industry a decade later, when swept wings and pressurized fuselages outran the reach of hand calculation. At Boeing, M. J. Turner assembled element stiffness matrices for wing panels; joined by Ray Clough of UC Berkeley and Harold Martin, the group published the landmark Turner–Clough–Martin–Topp paper in 1956 that laid out the direct stiffness method engineers still use ^[1]^[2]. Around the same time in Britain, John Argyris connected the energy methods of structural analysis to the matrix formulation, giving the approach its theoretical spine ^[1]^[2]. It was Clough who coined the name "finite element" in a 1960 paper ^[1]^[2].

Olgierd Zienkiewicz then carried the method out of aircraft structures into heat, fluids, and everything else, and his textbook made it the standard tool of the profession ^[1]^[2]. The rigorous mathematical footing — proving when and how fast the approximate answer converges to the true one — came in 1973 with Gilbert Strang and George Fix's analysis ^[1]. Engineering had the method working for two decades before the mathematicians finished proving why it worked.

Related tools

/tools/beam-deflection
/tools/stress-strain
/tools/euler-buckling
/tools/vibration-natural-freq
/tools/moment-of-inertia-shapes
/tools/section-modulus-rect
/tools/thermal-expansion

Sources

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