Dtypes, devices and numerics

Five element types, two backends, and a tolerance policy derived up front rather than tuned until the tests pass.

Five types, deliberately

ggml carries around 30 types because quantised inference needs them. Training does not, and every extra type multiplies the kernel matrix — quantisation is an explicit non-goal here. BF16 is absent because the target GPU does not support it, measured rather than assumed.

Only floating tensors can carry gradients, matching PyTorch. backward on an integer root raises DTypeError rather than producing zeros.

f16 is storage, never an accumulator

Values widen to float at the memory boundary and narrow once on the store. Both backends implement it the same way and the source says so — the CPU's widen and the shader's dtype-switched load are deliberately written to look alike.

The C++ side goes further: Half is a trivial storage wrapper with no arithmetic operators, on purpose. An implicit-conversion half type makes it far too easy to accumulate in 16 bits by accident, and the whole contract is that you cannot.

It matters most in matmul: an f16 accumulator over K = 784 would lose roughly three decimal digits, well outside the 1e-3 that the f16 tolerance allows.

The conversion itself

IEEE-754 binary16, handling subnormals, infinities and NaN correctly. The naive bit-shuffle most tutorials show silently flushes subnormals to zero, which would surface as a tolerance failure against PyTorch only for very small values and would be miserable to track down later. The branch-free approach follows Fabian Giesen's float_to_half_fast3, the same lineage ggml's implementation comes from.

On the GPU the narrowing is done in the integer domain rather than with float16_t(value), because SPIR-V leaves OpFConvert's rounding mode implementation-defined — see the shader page for that story.

Tolerance is a property of the operation

Individual tests do not choose their own tolerances. Twice during the project a test "failed" because its tolerance model was wrong rather than because the code was:

• matmul at K = 4096 missed a 1e-5 relative check while sitting 431× inside its backward-error bound — the dot product was ill-conditioned (Σ|aᵢbᵢ| = 4051 against a result of 3.17), so relative-to-result was meaningless.
• Vulkan exp() differed from glibc's by 1.14e-5 absolute, which is 3 ULP at magnitude 46 — inside the Vulkan specification's own allowance.

Both were the check being wrong, and both cost real debugging time. The tolerance for an operation is therefore derived from a citable source and stated once, in tests/python/tolerance.py, which carries 69 entries across four kinds — 39 exact, 16 relative, 10 ULP and 4 backward:

The sources are named: the Vulkan 1.3 specification's per-instruction ULP allowances (which are permissions granted to the driver, so a conforming implementation may legitimately differ from libm by that much); IEEE-754 for the correctly-rounded operations; and Higham for the backward-error bounds on summation.

Above that sits the class-level policy decided in advance — 1e-6 for element-wise f32, 1e-5 for reductions and matmul within their K bounds, 1e-3 for f16 storage with f32 accumulation. Any failure is investigated as a bug first, and a tolerance that genuinely needs to change must come with the error analysis that justifies it.

Determinism

Identical inputs give bit-identical outputs on the same device, and across drivers wherever the contract claims it. Two mechanisms carry it:

• Fixed reduction trees. Determined by shape, never by which workgroup finished first. An atomic-accumulation reduction would give a different answer each run — which is also why the absence of a global float atomicAdd on this hardware costs less than it appears to.
• Software f32→f16 narrowing, because the hardware instruction's rounding mode is implementation-defined.

Verified across a discrete RX 5600M and an integrated Renoir: MNIST and CIFAR-100 both produce the same accuracy and the same loss to the last digit on the two.

NaN follows PyTorch, on both backends

This was unwritten until it drifted: relu(nan) returned 0 while maximum(x, 0) and clamp_min(x, 0) — the same function spelled differently — returned NaN, and the Vulkan amax/amin reductions dropped NaN where the CPU propagated it.

A tolerance cannot express any of this. NaN is not far from a number, it is a different kind of answer.

The rule is that torch is the reference, because a user should get the same answer from vkML as from torch and the same answer from either backend. Where torch and NumPy disagree — sign(nan) is +0.0 in torch and NaN in NumPy — torch wins.

One mechanism explains every case: every comparison against NaN is false. So x > 0 ? x : 0 falls through to 0 and destroys a NaN, while x <= 0 ? 0 : x falls through to x and keeps it. The two are identical on numbers, so choosing between them is not a matter of style — and the same choice appears in relu's gradient, whose mask is x <= 0 for exactly this reason.

Comparison alone cannot make a min/max reduction propagate NaN at all, which is why those kernels test isnan explicitly at both fold stages.

One documented divergence: subnormals

Vulkan permits flush-to-zero for float32 denormals unless shaderDenormPreserveFloat32 is requested, and vkML does not request it. So relu(1e-45) is 0 on the GPU and 1e-45 on the CPU, and exp(-89) is 0 rather than the subnormal 2.227e-39.

This is pinned by tests rather than papered over, and the tolerance below which a disagreement carries no information is FLT_MIN — the smallest positive normal float.