How we numerically verify a physics engine.
A finite-element engine that gives an agent a wrong number is worse than no engine. The agent trusts the number. Thus in physicsbase the test suite is not an afterthought. It is the specification.
Most software tests check that code does what the author intended. Numerical-physics tests add comparisons with analytical solutions and governing equations. A cantilever tip deflects by \(PL^3/3EI\). A column buckles at \(\pi^2 EI/L^2\). A restrained bar with a temperature rise \(\Delta T\) carries \(-E\alpha\,\Delta T\). Failure to reproduce those assumptions reveals an implementation or modeling defect. Agreement is numerical verification—not evidence that the model agrees with nature.
The test is the specification
Each element and each analysis in physicsbase has at least one test. The test builds a problem with an exact known answer. It solves the problem. Then it checks the result against the answer within the tolerance. We do not trust a result because it "looks correct." Here is a real test. It is the Euler buckling of a pinned column:
def test_buckling_pinned(): m, p = _column("pinned") # pin-ended steel column res = buckling(m, num_modes=1) expected = math.pi**2 * p["E"] * p["I"] / p["L"]**2 # Euler's Pcr assert res.critical_loads[0] == pytest.approx(expected, rel=1e-2)
There are no golden files. There is no "close enough by eye." The critical load from the engine must match the Euler formula to 1%. If it does not, the build fails.
From 53 checks to an evidence program
The current numerical guide records 987 passing comparisons across 181 case definitions. A separate seven-study benchmark guide records four passing solution benchmarks and three passing experimental benchmarks. The measurement comparisons do not meet the stricter uncertainty and preregistration gate, so physical validation remains NOT_ESTABLISHED. The representative families include:
| Area | Checked against |
|---|---|
| Bars, trusses, beams | axial \(PL/AE\), tip \(PL^3/3EI\), rotation \(ML/EI\) |
| 3D frames | biaxial bending and torsion \(TL/GJ\) |
| Continuum & solids | uniform-strain patch tests for every element (CST, Q4, LST, Q8, tet4, hex8, axisymmetric) |
| Dynamics | cantilever & simply-supported natural frequencies; step-response 2× overshoot |
| Stability | Euler buckling — pinned, cantilever, fixed-fixed |
| Thermal & loads | \(-E\alpha\,\Delta T\), self-weight, distributed loads, tractions |
| Multiphysics | conduction profiles, convection, 2D/3D field patch tests, convection-diffusion, coupled thermal stress |
| Nonlinear | the damage law's stress-strain response |
The uniform-strain patch test is important. It is the cleanest idea in the file. Apply a displacement field for a constant strain to the nodes of an element. Ask the element for its stress. Check that it equals the exact \(\mathbf{D}\!:\!\boldsymbol{\varepsilon}\). If an element passes the patch test, its strain-displacement matrix is correct. This is the most important property of a new element. Each solid element in physicsbase passes it.
pytest -q and python -m verification.suite. Then inspect the Numerical Verification Guide, the Benchmark Guide, and the NOT_ESTABLISHED physical-validation status. No one suite is an unconditional answer to “Can I trust this number?”