Pendulum period vs. amplitude

2D rigid bodies + joints

A disk on a 1 m rod from a fixed pin, released from rest at 5° to 80°. Small swings take T₀ = 2π√(L/g) = 2.0066 s; larger ones take longer, exactly T = 4√(L/g)·K(sin(θ₀/2)) where K is the complete elliptic integral of the first kind. That curve is far from the small-angle formula at 80° (+13.7%), so matching it tests the full nonlinear motion and the joint solver's energy behavior.

Quantity	Expected	Measured
Period at 5°	2.008 s	2.007 s (-0.01%)	Pass
Period at 20°	2.022 s	2.022 s (-0.01%)	Pass
Period at 40°	2.069 s	2.069 s (-0.01%)	Pass
Period at 60°	2.153 s	2.153 s (-0.02%)	Pass
Period at 80°	2.283 s	2.282 s (-0.02%)	Pass

Pendulum energy over 30 s

2D rigid bodies + joints

The same pendulum swinging 60° for 30 s, with no friction: its energy should stay constant. The joint solver aims each step at the joint's predicted end-of-step position (a SHAKE-style scheme), so the error oscillates within a bounded band instead of steadily draining. (The textbook velocity-only constraint loses 23% of this pendulum's energy in 10 s.)

Quantity	Expected	Measured
Largest energy error	< 2%	1.24%	Pass
Drift (last 5 s vs. first 5 s)	< 0.5%	+0.088%	Pass

Rolling vs. sliding down a ramp

3D rigid bodies

A ball and a box released on a 20° slope. The ball rolls without slipping, so part of its energy goes into spin: it accelerates at (5/7) g sin θ. The box slides against Coulomb friction μ = 0.2: g (sin θ − μ cos θ). Both are distance-vs-time parabolas.

Quantity	Expected	Measured
Ball, rolling: distance after 1.5 s	2.695 m	2.709 m (+0.50%)	Pass
Box, sliding (μ = 0.2): distance after 1.5 s	1.700 m	1.709 m (+0.56%)	Pass

Bounce height vs. restitution

3D rigid bodies

A ball dropped 3 m onto the floor with coefficient of restitution e leaves at e times its impact speed, so it rises to e² of the drop height. The engine applies restitution to the velocity at the start of the contact step, which includes that step's gravity, so lively balls bounce about 2% higher than e²; the check allows ±0.03 h.

Quantity	Expected	Measured
Rebound height, e = 0.2	0.04000 × h	0.03299 × h (-0.007)	Pass
Rebound height, e = 0.4	0.1600 × h	0.1572 × h (-0.003)	Pass
Rebound height, e = 0.6	0.3600 × h	0.3640 × h (+0.004)	Pass
Rebound height, e = 0.8	0.6400 × h	0.6530 × h (+0.013)	Pass
Rebound height, e = 0.9	0.8100 × h	0.8283 × h (+0.018)	Pass

Mass on a spring

2D joints (spring)

Masses from 0.25 to 4 kg hanging on a 40 N/m spring joint, released at the spring's rest length. Each oscillates about the stretched equilibrium with period T = 2π√(m/k). The spring is integrated implicitly (a soft constraint), so it's stable at any stiffness.

Quantity	Expected	Measured
Period, m = 0.25 kg	0.4967 s	0.4986 s (+0.37%)	Pass
Period, m = 0.5 kg	0.7025 s	0.7038 s (+0.18%)	Pass
Period, m = 1 kg	0.9935 s	0.9944 s (+0.09%)	Pass
Period, m = 2 kg	1.405 s	1.406 s (+0.05%)	Pass
Period, m = 4 kg	1.987 s	1.987 s (+0.02%)	Pass

Soft-body stiffness is physical

2D soft bodies (XPBD)

Five 0.1 kg particles hanging from a pin by springs of stiffness k. Edge i carries the weight of everything below it, so the chain stretches m g n(n+1)/(2k) in total. Plain position-based dynamics gets stiffer the more iterations it runs; XPBD's compliance term makes stiffness a real N/m, so 10 and 40 substeps give the same, correct answer.

Quantity	Expected	Measured
Stretch at 200 N/m, 10 substeps	0.07355 m	0.07362 m (+0.09%)	Pass
Stretch at 200 N/m, 40 substeps	0.07355 m	0.07355 m (+0.01%)	Pass

Hydrostatic pressure in still water

2D SPH fluid

A 1 m column of water left to settle, then particle pressures averaged over 1 s in 0.1 m depth bands away from the walls. Still water's pressure grows linearly with depth, p = ρ₀ g h (ρ₀ = 10 kg/m², a 1 cm slab of water in 2D). The bands nearest the free surface and the floor read low and high: SPH particles there have one-sided neighborhoods.

Quantity	Expected	Measured
Pressure gradient dp/dh (interior fit)	98.07 N/m²	98.26 N/m² (+0.20%)	Pass

Archimedes' principle

3D SPH fluid + rigid bodies

0.3 m cubes of density 250, 500, and 750 kg/m³ floating in water (1000 kg/m³). A floating body displaces its own weight of water, so the fraction under water equals the density ratio. The water and the cubes push on each other with exactly opposite impulses.

Quantity	Expected	Measured
Submerged fraction, 250 kg/m³	0.2500	0.2552 (+2.10%)	Pass
Submerged fraction, 500 kg/m³	0.5000	0.4994 (-0.12%)	Pass
Submerged fraction, 750 kg/m³	0.7500	0.7448 (-0.69%)	Pass

Heat diffusion

2D smoke & heat grid

A Gaussian hot spot (σ₀ = 0.2 m, 10 K) in still air with thermal diffusivity κ = 0.01 m²/s. It stays Gaussian with variance σ² = σ₀² + 2κt, and in 2D its peak falls as σ₀²/σ² so the total heat is unchanged. The solver diffuses heat implicitly, so it's stable at any κ and step size, and its insulated walls conserve heat exactly.

Quantity	Expected	Measured
Peak at t = 0.5 s	8.000 K	8.025 K (+0.31%)	Pass
Peak at t = 2 s	5.000 K	5.025 K (+0.50%)	Pass
Peak at t = 5 s	2.857 K	2.869 K (+0.40%)	Pass

Warm air starts to rise

2D and 3D smoke & heat grids

A patch of air 10 K warmer than its surroundings is lighter (Boussinesq: g′ = g β ΔT with β = 1/T for an ideal gas, here 0.3345 m/s²), but it can't rise without pushing the surrounding air out of the way. That air's inertia acts as added mass: equal to the patch's own for a 2D cylinder, half of it for a sphere. So the patch starts out accelerating at g′/2 in 2D and 2g′/3 in 3D, a sharp test of the pressure solve.

Quantity	Expected	Measured
Initial acceleration, 2D cylinder	0.5000 × g′	0.4881 × g′ (-2.37%)	Pass
Initial acceleration, 3D sphere	0.6667 × g′	0.6581 × g′ (-1.29%)	Pass

Pressure through two liquids

2D SPH fluid, several liquids

0.3 m of mercury under 0.5 m of water, settled, then pressures averaged over 2 s in 5 cm bands. Pressure grows at ρg in each liquid: 13.5× faster in the mercury. This needs the multiphase formulation (each particle sums neighbors by number, not mass, so a water particle beside mercury doesn't count mercury's mass as its own density). Bands at the floor and at the interface read off: one-sided SPH neighborhoods there.

Quantity	Expected	Measured
dp/dh in the water	98.07 N/m²	104.7 N/m² (+6.73%)	Pass
dp/dh in the mercury	1328 N/m²	1347 N/m² (+1.43%)	Pass

Heat conduction in a liquid

2D SPH fluid, heat

Still mercury with a sinusoidal temperature pattern, 5 K amplitude over a 5 cm wavelength (20 particles per wavelength). The heat equation says each Fourier mode decays independently, A(t) = A₀ e^(−αk²t), with mercury's thermal diffusivity α = k/(ρc) = 4.39e-6 m²/s. The SPH conduction term (Cleary & Monaghan) uses each liquid's real conductivity and specific heat and conserves energy pair by pair. The box is closed and there's no gravity (so nothing convects): with an open top, surface tension, which rules in zero gravity, would pull the mercury into a drop and stir it.

Quantity	Expected	Measured
Decay rate αk²	0.06936 1/s	0.06591 1/s (-4.97%)	Pass

Heat from a hot plate into water

3D SPH fluid, heat

A 0.6 m square plate held at 40, 60, or 80 °C under 20 °C water. The thermal boundary layer on a heated wall is far thinner than a particle, so the engine carries heat across the first layer at the natural-convection rate q = hΔT, with h = 0.13 k (gβΔT/να)^(1/3) (turbulent, McAdams) from water's properties at the film temperature: the flux grows as ΔT^(4/3). This checks that the SPH wall term delivers that flux (plate area 0.36 m²; the plate's edges add a little).

Quantity	Expected	Measured
Heat in over 0.25 s, plate at 40 °C	1125 J	1151 J (+2.37%)	Pass
Heat in over 0.25 s, plate at 60 °C	3276 J	3354 J (+2.36%)	Pass
Heat in over 0.25 s, plate at 80 °C	6311 J	6459 J (+2.35%)	Pass

Water flash-boiling in vacuum

3D SPH fluid, boiling

A 0.4 m pool of 20 °C water with nothing above it. Water boils where the pressure on it is below its vapor pressure (2.34 kPa at 20 °C), and the only pressure here is the water's own weight, ρgh. So the top 24 cm boils, each layer cooling until it reaches its own boiling point, T_sat(ρgh) from the IAPWS-IF97 steam-table equation, while the water below stays at 20 °C. The heat that leaves is exactly the latent heat of the vapor that escapes. The profile is read after 0.25 s: the top band lags because a free surface can only lose molecules at the Hertz–Knudsen rate, and later the chilled water (denser than 20 °C water) sinks and mixes, as it would in reality. (Without freezing in the model, the surface stops at 0.01 °C instead of icing over.)

Quantity	Expected	Measured
Largest band deviation from T_sat(ρgh)	0 K	1.136 K (+1.136)	Pass
Water (liquid + vapor) conserved	64.00 kg	64.00 kg (-0.00%)	Pass
Heat lost = latent heat of the vapor	1757156 J	1757156 J (+0.00%)	Pass

Kepler's third law

N-body orbits

Planets placed on circular orbits from Mercury's distance to Neptune's, each timed over two orbits. Kepler's third law says the period grows as the distance to the 1.5 power: T = 2π√(a³/GM), one year at 1 AU. The integrator (leapfrog, at least 256 steps per orbit) runs orbits long by a fixed (2π/256)²/3 = 0.02%, and being symplectic, its energy error stays bounded instead of drifting: the full Solar System, with every planet and moon pulling on every other, is checked over 20 years.

Quantity	Expected	Measured
Period at 0.387 AU	0.2408 years	0.2408 years (+0.02%)	Pass
Period at 1 AU	1.000 years	1.000 years (+0.02%)	Pass
Period at 30.1 AU	165.1 years	165.2 years (+0.02%)	Pass
Solar System energy change, 20 years	< 10⁻⁶	1.7e-8	Pass
Solar System angular momentum change, 20 years	< 10⁻⁹	2.1e-14	Pass