This function explores tiny values with nonzero low-order digits, and it looks at the differences between them. The loop has two layers, skipping over the inner layer as soon as the difference between tiny values fails the test
def tiny_values_and_difference(one_plus):
"""Compute the tiniest values from a starting value with nonzero
low-order bits. This function is parallel to tiny_powers_of_B(),
but looks for nonzero differences of tiny values, which are
evidence of gradual underflow.
If multiplication is too inaccurate, raise an error and return default
values, shown in parens under Returns.
Args:
one_plus: 1 plus enough ulps that they shouldn't be lost when
multiplied by a power of B
if (tiny_delta + tiny_delta) <= tiny_delta:
which says tiny_delta has taken a nonzero and non-𝛜 value.
This indicates that there isome loss of information due to underflow.
While too_tiny is in the range of normal numbers,
tinier scales back to tiny
by 1/H exactly.
The next diagram shows the situation in binary in the typical case
that multiplication is accurate enough that a single unit in the last place
will detect differences between two nearly similar values after scaling by
1/H.
Returns:
tiny: smallest well-behaved value (0) *** why not tiny_B?
tinier: first value to exhibit underflow (tiny_B)
too_tiny: underflows to 0 or an epsilon value (too_tiny_B)
u_threshold: smallest well-behaved value, if underflow degrades
into nonzero values, or 0 if underflow flushes to 0 (tiny_B)
tiny_delta: difference between two values caused by underflow (0)
save_tiny1: value of tiny where tiny_delta computed (0)
save_tiny2: recomputed tiny where tiny_delta computed (0)
Basic 4450-4520
"""
If too_tiny underflows to zero,
the loop ends with tiny_delta still zero.
But with any form of gradual underflow, the loss of precision is detected.
Depending on the design of the arithmetic and which mode of rounding is
selected, the low-order bit of tinier may be 0 or 1.
This explains the use of
the absolute value in computing tiny_delta.
d = C * one_plus
# 1 + SAFE_ULPS_OF_ONE should be enough ulps of 1 that d > C.
if d <= C:
# But if not, we punt with default values.
bad_cond(err_failure,
"multiplication gets too many last digits wrong.")
pause()
return 0, tiny_B, too_tiny_B, tiny_B, 0, 0, 0
The code sets the underflow threshold value to tiny the first time this behavior is discovered.
[TODO: Provide an example of an arithmetic with multiplication so bad that it falls back to default values.]
else:
tinier = d
too_tiny = tinier * H
u_threshold = ZERO
tiny_delta = ZERO # indicates no value computed yet
save_tiny1 = ZERO
save_tiny2 = ZERO
while True:
tiny = tinier
tinier = too_tiny
# This loop mirrors the computation in tiny_powers_of_B(),
# except that the launch point is C*s with nonzero low-order
# bits, versus C, which is just a power of the radix.
if (tiny_delta + tiny_delta) <= tiny_delta:
if ultra_verbose:
print("\tUnderflow loop: within if, tinier = {:.7e}"
.format(tinier))
save_tiny2 = tinier * ONE_OVER_H
tiny_delta = fabs(tiny - save_tiny2) # may underflow to 0
save_tiny1 = tiny
if (u_threshold == ZERO) and (tiny != save_tiny2):
# Test for 0 captures first time through, in case
# tiny_delta underflows to 0.
u_threshold = tiny
# Placing this statement after the (tiny*H)/H test above
# emphasizes that only nonzero values of tinier qualify for
# values degraded by underflow.
too_tiny = too_tiny * H
if ultra_verbose:
print("\tUnderflow loop: too_tiny = {:.7e}".format(too_tiny))
if (tinier <= too_tiny
or too_tiny + too_tiny <= too_tiny): break
return (tiny, tinier, too_tiny, u_threshold,
tiny_delta, save_tiny1, save_tiny2)