Code as Physics

Physics expressed and tested as code. The other pages in this series turn code into something to look at; this one turns a physical claim into something to rerun. Nothing here is asserted on authority: every claim is an experiment, with the code that tests it, the hardware it ran on, the numbers it produced and the result that would have proved it wrong.

Experiment log

1 published.

• 001 · Energy drift: explicit Euler against symplectic integrators
  Experiment 001published2026-10-26NVIDIA Warp 1.18.0

  Question

  Over 10,000 to 100,000 steps, does the choice of integrator decide whether a simulated system keeps its energy?

  Hypothesis

  Explicit Euler adds energy every step, so its relative energy error grows without bound. Semi-implicit (symplectic) Euler and velocity Verlet keep the error bounded and oscillating, with a size that scales as dt and dt squared respectively.

  Method

  Three integrators (explicit Euler, semi-implicit Euler, velocity Verlet) on two systems: a harmonic oscillator (m = k = 1, start x = 1, v = 0) and a Kepler two-body orbit reduced to the relative coordinate (GM = 1, eccentricity 0.5, start at perihelion r = 1). 100,000 steps of dt = 0.01 in float64, about 159 oscillator periods and 56 orbits. Each trajectory is one Warp thread that records |E - E0| / |E0| at every step. Three checks: a 64-trajectory ensemble per system with random amplitude, phase, eccentricity, perihelion and orientation; a rerun at dt/2 over the same physical time to measure the convergence order; and the same update rules in plain numpy as a second, independent instrument.

  Engine

  NVIDIA Warp 1.18.0, CPU backend (no GPU on this host); numpy 2.5.4 (cross-check), Python 3.14.8

  Hardware

  Apple M4 Pro Mac mini, macOS 27.0.1, CPU only

  Run

  seed 20261026, dt = 0.01, 100,000 steps, float64, ensemble of 64 per system. The two canonical runs are deterministic; the seed draws the 64-trajectory ensembles only.

  Result

  Supported. Explicit Euler on the oscillator multiplies the energy by exactly (1 + dt²) every step: relative error 1.72 at step 10,000 and 2.20 × 10⁴ at step 100,000, the same to three figures for all 64 starting conditions. On the Kepler orbit it keeps adding energy (error 0.752 at the end, so the orbit has widened from a semi-major axis of 2 to about 8), and halving dt cuts the worst error only by a factor of 1.16. Both symplectic methods stay bounded: their largest error in the second half of the run equals the first half to three figures. At the last step semi-implicit Euler reads 4.60 × 10⁻³ (oscillator) and 8.83 × 10⁻⁴ (orbit); velocity Verlet reads 1.72 × 10⁻⁵ and 3.35 × 10⁻⁵. Halving dt divides the largest error by 2.01 and 2.02 for semi-implicit Euler and by 4.00 and 4.00 for Verlet: first and second order, as predicted. The numpy reproduction agrees with Warp on every final error to within 1.08 × 10⁻⁶ relative.

  Would falsify it

  Any of these would have refuted it: a symplectic method whose largest error in the second half of the run is clearly larger than in the first (secular drift); a dt/2 ratio far from 2 for semi-implicit Euler or 4 for Verlet; an ensemble member that drifts while the canonical run does not; or Warp and numpy disagreeing beyond rounding.

  Limits

  This is a textbook result, chosen on purpose to calibrate the instrument before anything new is claimed. Bounded energy error does not mean an exact orbit: the symplectic runs still accumulate phase error, which this experiment does not measure.

  Source

  code-as-physics/experiments/001_energy_drift at 76b3df9 — lab repository, MIT, not yet published

  System	Integrator	Error at step 10,000	Error at final step	Max error, 2nd half / 1st half	Max error, dt / (dt/2)	Ensemble final, median / max
  harmonic oscillator	explicit Euler	1.72	2.20×104	149.00	149.00	2.20×104 / 2.20×104
  harmonic oscillator	semi-implicit (symplectic) Euler	4.38×10-3	4.60×10-3	1.00	2.01	5.09×10-3 / 8.27×10-3
  harmonic oscillator	velocity Verlet	6.40×10-6	1.72×10-5	1.00	4.00	1.64×10-5 / 2.07×10-5
  Kepler two-body orbit	explicit Euler	0.343	0.752	1.18	1.16	0.77 / 0.963
  Kepler two-body orbit	semi-implicit (symplectic) Euler	2.24×10-4	8.83×10-4	1.00	2.02	1.10×10-3 / 0.0108
  Kepler two-body orbit	velocity Verlet	3.34×10-5	3.35×10-5	1.00	4.00	2.23×10-5 / 2.85×10-4

  Relative energy error |E − E₀| / |E₀|. Every number is output of the run.

The first experiments reproduce known results on purpose: an instrument earns trust by getting textbook answers right before it is pointed at anything new. Records are written as drafts while an experiment runs and appear here only once they have been run, measured and given a falsifier.