ULP (Units in the Last Place) distance helpers for f32 and f64 comparison.
A ULP is the gap between two adjacent representable float values at a given
magnitude — it grows as values get larger. Comparing by ULP distance rather
than absolute error gives a scale-invariant correctness criterion: a result
0 ULP away is bit-for-bit identical; 1 ULP means the nearest representable
neighbour; anything beyond that is a rounding error.
ULP distance is computed by reinterpreting both floats as (unsigned, for
f32; signed, for f64) integers of the same width and taking their absolute
difference. IEEE 754 bit patterns are ordered the same way as the real line
(for same-sign values), so this directly counts how many representable
values lie between a and b.
±0 guard — WGSL §15.7.1 permits implementations to return either +0
or -0 interchangeably. Their bit patterns differ (0x00000000 vs
0x80000000), which would give a huge ULP distance even though both
represent the same real number. Both ulpDiff and ulpDiff64 special-case
this: if both values have absolute value 0 the distance is 0.
Same-sign assumption — like the f32 functions, ulpDiff64/maxUlp64
only give a meaningful distance when a/b share a sign; every current
caller compares non-negative values (e.g. an asum-style result), so this
has never needed a general any-sign fixup.
When not to use ULP
For routines that involve multiply-add (saxpy, sdot, srot, srotm), the GPU
may execute the operation as a fused multiply-add (FMA) — a single rounding
instead of two. Near cancellation this makes the raw ULP difference between
a fused and an unfused result unbounded even though both are individually
correct. Use the forward error factor from tests/<routine>/helpers.js
instead. See tests/helpers/fixtures.js for which metric each routine uses.
For f64-emulated routines (e.g. dasum), a huge ULP64 count from a tiny
absolute error used to come from aux (the packed extra bits — see
src/util/f64pack.mjs) landing on a NaN/Infinity bit pattern: any ordinary
double, or any intermediate sum during a reduction, can produce this — not
just unusual inputs — and a real Float32Array silently canonicalizes
(quiets) a NaN bit pattern on any round trip, corrupting it. That's now
fixed at the source: aux is transported exclusively as raw u32 bits (never
an actual float value) throughout f64pack.mjs, GpuVector, and the WGSL f64
shaders, so this class of corruption can't happen at all — raw ULP64 is a
meaningful metric and is what dasum's own tests use.
ULP (Units in the Last Place) distance helpers for f32 and f64 comparison.
A ULP is the gap between two adjacent representable float values at a given magnitude — it grows as values get larger. Comparing by ULP distance rather than absolute error gives a scale-invariant correctness criterion: a result 0 ULP away is bit-for-bit identical; 1 ULP means the nearest representable neighbour; anything beyond that is a rounding error.
ULP distance is computed by reinterpreting both floats as (unsigned, for f32; signed, for f64) integers of the same width and taking their absolute difference. IEEE 754 bit patterns are ordered the same way as the real line (for same-sign values), so this directly counts how many representable values lie between
aandb.±0 guard — WGSL §15.7.1 permits implementations to return either
+0or-0interchangeably. Their bit patterns differ (0x00000000vs0x80000000), which would give a huge ULP distance even though both represent the same real number. BothulpDiffandulpDiff64special-case this: if both values have absolute value 0 the distance is 0.Same-sign assumption — like the f32 functions,
ulpDiff64/maxUlp64only give a meaningful distance whena/bshare a sign; every current caller compares non-negative values (e.g. anasum-style result), so this has never needed a general any-sign fixup.When not to use ULP
For routines that involve multiply-add (saxpy, sdot, srot, srotm), the GPU may execute the operation as a fused multiply-add (FMA) — a single rounding instead of two. Near cancellation this makes the raw ULP difference between a fused and an unfused result unbounded even though both are individually correct. Use the forward error factor from
tests/<routine>/helpers.jsinstead. Seetests/helpers/fixtures.jsfor which metric each routine uses.For f64-emulated routines (e.g. dasum), a huge ULP64 count from a tiny absolute error used to come from
aux(the packed extra bits — see src/util/f64pack.mjs) landing on a NaN/Infinity bit pattern: any ordinary double, or any intermediate sum during a reduction, can produce this — not just unusual inputs — and a real Float32Array silently canonicalizes (quiets) a NaN bit pattern on any round trip, corrupting it. That's now fixed at the source: aux is transported exclusively as raw u32 bits (never an actual float value) throughout f64pack.mjs, GpuVector, and the WGSL f64 shaders, so this class of corruption can't happen at all — raw ULP64 is a meaningful metric and is what dasum's own tests use.