Computational Cost and Numerical Accuracy

Performance matters scientifically when it changes which validated calculations can be completed, repeated, refined, or independently reproduced. Raw floating-point rate, device occupancy, or one fast step does not establish that improvement. The comparison unit should be resources required to reach the same declared physical error for the same problem, including setup, communication, adaptation, I/O, failed steps, and verification. This page distinguishes diagnostic performance measures from evidence that an algorithm expands reliable scientific reach.

Helpful background. Formulations of General Relativity supplies the evolved physical systems, Ordinary Differential Equations supplies tolerance and stability concepts, Computational Relativity supplies evidence standards, and Probability and Statistical Inference supplies uncertainty summaries for timings and observables.

Define time to a validated result

Choose an observable $Q$, target uncertainty $\epsilon_Q$, physical duration, and benchmark input. For each admissible configuration $c$ of implementation $A$, define a resource vector rather than one number,

$$\mathcal R_A(c)= (T_{\rm wall},\;N_{\rm node\,h},\;E_{\rm energy},\;M_{\rm peak},\;V_{\rm output})$$

and collect the nondominated resource vectors at the requested accuracy,

$$\mathcal P_A(\epsilon_Q)=\operatorname{Pareto} \left\{\mathcal R_A(c):e_Q(c)\leq\epsilon_Q\right\}.$$

Report which components were measured and how. There is no universal scalar conversion between wall time, allocation, energy, and memory, so “least expensive” is undefined until a scalar objective or explicit resource constraints are chosen. If one is chosen, state it and select from $\mathcal P_A$ rather than concealing tradeoffs in an implicit ordering.

The target error includes spatial and temporal truncation, iterative tolerance, boundary and extraction errors, and any stochastic or precision component relevant to $Q$. Two codes at the same grid spacing need not have the same error. Conversely, a costlier step may win if its convergence order reaches the target with far fewer degrees of freedom.

Build a fair comparison

Pin physical equations, initial and boundary data, units, requested output, evolution interval, and failure criteria. Allow each method to choose a documented resolution sequence, then estimate its continuum error using the same observable norm. Include compilation or initialization only if it matters at the intended production cadence; always state the choice. Include I/O when production requires it.

Use a real baseline that runs efficiently and fits memory. Strong-scaling efficiency relative to one underfilled node can exaggerate gains. Report problem size, machine, compiler, software stack, process and device layout, repetitions, median and variation, and profiling categories. Kidder et al. 2017, §§5.1–5.3 and Figures 4–6 ties task-based performance tests to explicit numerical workloads.

Compare not only elapsed time but compute, communication, idle, load imbalance, adaptation, solver iterations, and output. These components explain why scaling saturates and whether a proposed improvement transfers to a production regime.

A useful first estimate follows from $e_Q\simeq Ch^p$. In $d$ dimensions, a fixed-duration explicit hyperbolic calculation has $O(h^{-d})$ cells and $O(h^{-1})$ time steps, giving work $W\simeq Kh^{-(d+1)}$. Eliminating $h$ gives

$$W(\epsilon_Q)\simeq K \left(\frac{C}{\epsilon_Q}\right)^{(d+1)/p}.$$

In three dimensions, reducing error by a factor sixteen multiplies idealized work by $256$ for a second-order method and by $16$ for a fourth-order method. The constants, communication, shocks, and memory can reverse a comparison at practical resolution, so this is an explanatory scaling estimate rather than a performance prediction. It shows why equal cell counts are an incomplete comparison.

Specified validated-performance record

Use the analytic wave

$$\partial_tu+\partial_xu=0, \qquad u_{\rm ex}(t,x)=\sin\!\bigl[2\pi(x-t)\bigr]$$

on $[0,1]^3$ with periodic boundaries and $u(0,x)=\sin(2\pi x)$. Use fourth-order centered differences, RK4, CFL $0.2$, and $N^3$ grids with $N=32,64,128,256$ to $t=1$. The analytic solution is the reference. Compare using the volume-weighted relative $L_2$ error $E_2$ and the phase error $\Delta\phi$ of an independently computed fundamental Fourier coefficient. A configuration is scientifically admissible only if the $256^3$ run has $E_2\leq10^{-5}$ and $|\Delta\phi|\leq10^{-5}\,\mathrm{rad}$ and if the final two refinement ratios approach fourth order to within $20\%$.

For strong scaling set $P_0=1$ node and hold the $256^3$ problem, diagnostics, and output cadence fixed for $P=1,2,4,8$ nodes. For weak scaling tile $1,2,4,8$ unit cubes along $x$, retain the same exact wave with $128^3$ cells per node, and keep ranks, threads, and devices per node fixed. Perform one warm-up and seven measured repetitions, reporting median and interquartile range for wall time and every available component of $\mathcal R_A$. Include a run in the Pareto set only when its own $E_2$ and $\Delta\phi$ pass the same limits; report solution differences between layouts in relative $L_2$ norm and require them to be at most $10^{-10}$ on corresponding cubes.

The discrete Fourier calculation explains the per-tile resolution: $64^3$ cells already exceed the phase target before parallel overhead is considered.

Observed discrepancy: not available—specification only.

Scientific checks accompany timings

For an analytic wave, require convergence to the exact phase and amplitude. For a black hole, require constraints and horizon invariants. For GRMHD, include conserved quantities, primitive-recovery failures, and shock benchmarks. Performance runs that omit expensive verification or use fewer outputs are a different workflow unless the production claim also omits them.

Cross architecture and decomposition changes can alter floating-point order. Demonstrate that those changes remain within the physical error budget. Zhu et al. 2024, §§3–5 combines standard relativity tests with large accelerator scaling; its device-specific efficiencies illustrate why workload and baseline must remain attached to every number.

The strongest evidence is a Pareto curve of resources versus validated error. A method dominates only where it reaches equal or smaller error with no greater resource on every chosen axis. Otherwise report the tradeoff. One can legitimately conclude that a method is faster on a named machine while consuming more energy or memory.

Environment and failure criteria

Pin hardware, allocation, compiler and flags, libraries, code revision, workload checksum, numerical tolerances, resolution and refinement policy, output cadence, profiling tools, warm-up, repetition protocol, and queue or filesystem conditions. A clean interface could be validated-performance --config wave-target.yaml --output run/, emitting observables, error estimates, constraint checks, timings, resources, profiles, and environment metadata. It specifies a reproducible experiment and does not imply portability beyond the tested system.

Reject comparisons with different target errors, missing failures, cherry-picked fastest trials, synthetic baselines, changed physics, or excluded production costs. Other failure modes include measuring only kernel time while communication dominates, counting nominal rather than occupied resources, and extrapolating weak scaling to a different AMR pattern. A speedup whose uncertainty overlaps one, or whose observable exceeds tolerance, is not supporting evidence.

Diagnostic

Method A takes half the wall time of method B at the same cell count. A is second order and misses the waveform-phase tolerance; B is fourth order and meets it. At the resolution where A meets the tolerance, it costs twice as many node-hours as B. Which is faster scientifically?

Solution

For the declared phase target and allocation measure, B reaches the validated result with fewer node-hours. A’s same-cell-count wall-time advantage answers a different question. Report both performance curves, the common error estimator, wall time, node-hours, and memory; do not compress them into a context-free speedup.

References

• Kidder, L. E., Field, S. E., Foucart, F., et al. (2017). “SpECTRE: A task-based discontinuous Galerkin code for relativistic astrophysics.” Journal of Computational Physics 335, 84–114. arXiv.
• Zhu, H., Fields, J., Zappa, F., et al. (2024). “Performance-portable numerical relativity with AthenaK.” arXiv:2409.10383. arXiv.