Decide whether inf or NaN is legitimate#

The problem#

An implementation returns inf or NaN for some elements. That may be a real overflow that any implementation would hit, an overflow caused by a poor order of additions, NaN from uninitialized memory, or NaN that the inputs force. torch.isfinite and allclose(..., equal_nan=True) cannot tell these apart.

Why the usual approach falls short#

Whether inf is acceptable depends on the exact value and on how close it is to the top of the format's range. Whether NaN is acceptable depends on the special values in the inputs. Neither is visible from the output.

The solution#

covers() judges infinities and NaN too: an inf or NaN output passes only where may_be() says a correct implementation can produce it. A sum near the top of float32's range may overflow in a correct implementation:

Code
import torch
import exactref as er

# A float32 sum whose first value is float32's largest finite value
top = er.sum(torch.tensor([3.4028234e38, 2e31, -1e31]))
# The exact sum, rounded once to float32
print("rounded:", top.rounded.tolist())
# Where a correct implementation may output +inf
print("may be +inf:", top.may_be("+inf").tolist())

# Check an output of +inf
print("+inf covered:", top.covers(torch.tensor(float("inf"))).all())
Output
rounded: 3.4028234663852886e+38
may be +inf: True
+inf covered: True

A float16 loop that overflows in its running sum, on inputs whose sum fits, is not correct by either standard. Reordering, or a float32 accumulator, avoids the overflow:

Code
# Three float16 values whose exact sum, 60000, fits in float16
x = torch.tensor([60000.0, 60000.0, -60000.0], dtype=torch.float16)
# A float16 loop: the first addition overflows to inf
loop = (x[0] + x[1]) + x[2]
print("loop:", loop.item())

# ExactRef's reference for the same sum
r = er.sum(x)
# The exact sum rounded once to float16, and whether +inf is allowed
print("rounded:", r.rounded.tolist())
print("may be +inf:", r.may_be("+inf").tolist())
# The second argument is the compute format; float32 is the default here
print("covered, float32 accumulation:", r.covers(loop).all())
print("covered, float16 accumulation:", r.covers(loop, "float16").all())

# Normal accumulation in float16, in index order: where adding the next
# value would overflow, it starts a new partial sum instead
print("normal:", r.normal("float16").tolist())
Output
loop: inf
rounded: 60000.0
may be +inf: False
covered, float32 accumulation: False
covered, float16 accumulation: False
normal: 60000.0

NaN passes where the inputs force it, and fails anywhere else:

Code
# A sum of +inf and -inf: the inputs force NaN
both_infinities = er.sum(torch.tensor([float("inf"), float("-inf")]))
print("NaN covered:", both_infinities.covers(torch.tensor(float("nan"))).all())

# The float16 sum from above, whose finite inputs force no NaN
nan16 = torch.tensor(float("nan"), dtype=torch.float16)
print("float16 sum, NaN covered:", r.covers(nan16).all())
Output
NaN covered: True
float16 sum, NaN covered: False

For complex outputs, each part is judged on its own, so a finite imaginary part is still checked next to an infinite real part.