ExactRef

Exact reference answers, with proven error bounds, for testing numerical code#

For an operation and its inputs, ExactRef gives you the exact answer rounded once to the output format, and a per-element tolerance: how far the output of a correct finite-precision implementation can be from that answer, on those inputs. Use it to check PyTorch operations, custom kernels and numerical libraries.

Get started GitHub

Install#

pip install exactref

ExactRef needs Python 3.10 or newer. Its one dependency, python-flint, is installed with it. NumPy and PyTorch are not needed, but ExactRef reads their arrays and tensors directly when you have them.

Quick start#

A plain Python loop adds ten copies of 0.1 and gets 0.9999999999999999. Is that a correct answer?

Code

import exactref as er

# Ten copies of 0.1
xs = [0.1] * 10

# A plain loop, rounding to float64 at each addition
total = 0.0
for x in xs:
    total += x

print("loop:", total)

# ExactRef's reference for the same sum
r = er.sum(xs)
# The exact sum, rounded once to float64
print("rounded:", r.rounded.tolist())
# How far a correct loop can land from 1.0
print("tolerance:", r.tolerance().tolist())
# Can a correct loop give this output?
print("covered:", r.covers(total).all())

Output

loop: 0.9999999999999999
rounded: 1.0
tolerance: 8.881784197001252e-16
covered: True

Yes: the loop is off by one unit in the last place, well inside what rounding at each addition allows. The same check works on a PyTorch matrix product:

Code

import torch
import exactref as er

# Seeded, so every run uses the same values
torch.manual_seed(0)
# A random 64x256 and 256x32 float32 matrix
A, B = torch.randn(64, 256), torch.randn(256, 32)

# The reference for A @ B
r = er.matmul(A, B)
# Every element must be one a correct implementation can produce
assert r.covers(A @ B).all()

What every call returns#

Every operation returns a Result:

r.rounded: the exact answer, rounded once to the output format.
r.exact: certified enclosures of the exact answer, for comparisons finer than the output format.
r.tolerance(): per element, how far a correct implementation's finite output can be from rounded.
r.covers(output): per element, whether a correct implementation can produce output.
r.normal(): what an ordinary implementation computes, adding in index order, emulated exactly.

Results describes each of them and their options.

How to read a tolerance#

An implementation's output, element by element:

Within tolerance(): no worse than normal accumulation can be. Normal accumulation is the way ordinary kernels add: in index order, or in a few SIMD lanes.
Beyond it, but within tolerance(standard="any_order"): still correct, but worse than normal accumulation can be, for example from a poor order of additions.
Beyond both: no implementation that rounds only where the operation requires, without overflowing, could produce it.

r.covers(output) makes the first judgment and r.covers(output, standard="any_order") the second, infinities and NaN included. How tolerance works shows each case in a few lines of Python.

What you can do with it#

Replace hand-picked test tolerances with bounds computed from your inputs.
Test a custom kernel, or a PyTorch operation across dtypes.
Check low-precision kernels and fused attention kernels.
Accept nondeterministic reductions that differ from run to run.

Solutions has these and more, each with code.

Operations#

Operation	Call	Inputs	Backward	Normal answer
Matrix product	`er.matmul`	real, complex, integer	yes	yes
Matrix product added to a matrix (gemm)	`er.gemm`	real, complex, integer	yes	yes
Elementwise sum and product	`er.add`, `er.multiply`	real, complex, integer	yes	yes
Sum over dimensions	`er.sum`	real, complex, integer, bool	yes	yes
Matrix exponential	`er.matrix_exp`	real, complex	yes	no
Its Frechet derivative	`er.matrix_exp_frechet`	real, complex	no	no
Scaled dot-product attention	`er.scaled_dot_product_attention`	real	yes	yes

Formats: float64, float32, float16, bfloat16, the float8 formats, complex32, complex64 and complex128, the signed and unsigned integer formats, and bool, with all five IEEE rounding modes.

Why trust it#

Every guarantee and every tolerance has a written proof, in the proofs folder on GitHub.
When it cannot certify an exact answer, ExactRef raises ReferenceUnavailableError instead of guessing.
Its answers are tested against independent references: MPFR and MPC, SymPy, NumPy, PyTorch and the Python standard library.