768 lines
30 KiB
Python
768 lines
30 KiB
Python
#!/usr/bin/env python3
|
|
# genFpTables.py - expected-value tables for the soft-float smoke checks in
|
|
# scripts/smokeTest.sh ("float/double compares" and "float/double NaN and
|
|
# infinity"). Paste the output over the matching table in smokeTest.sh.
|
|
#
|
|
# python3 tests/fpNan/genFpTables.py cmp # kFBits / kDBits / kExp
|
|
# python3 tests/fpNan/genFpTables.py special # kCases
|
|
# python3 tests/fpNan/genFpTables.py edge # kEdge
|
|
#
|
|
# Comparisons come from Python's own float compares, which are IEEE 754
|
|
# (every ordered compare with a NaN is false). Float operands are widened
|
|
# to double exactly, so float and double share one truth table per value.
|
|
#
|
|
# The special-value table is a model of the documented runtime rules, not
|
|
# of the host FPU: a NaN operand comes back quieted with its payload kept
|
|
# (the first NaN operand when both are), an invalid operation gives the
|
|
# default NaN 0x7FC00000 / 0x7FF8000000000000 (x86 would give the negative
|
|
# one), and float/double to integer saturates with NaN converting to 0.
|
|
# Finite cases use operands whose results are exact, because softFloat.c
|
|
# truncates instead of rounding.
|
|
#
|
|
# The math.c cases follow C99 Annex F's special values (modelled here where
|
|
# Python raises instead of answering). The series-based functions (exp,
|
|
# pow with a non-integer exponent, sinh, cbrt, ...) are only ~1e-6
|
|
# accurate, so their ordinary-value cases carry an `approx` flag: the check
|
|
# accepts a result within 1e-6 relative of the host's value.
|
|
#
|
|
# The edge table is math.c's exact and bit-level functions, checked bit for
|
|
# bit against Python's (C99 Annex F) answers: fmod and remainder must be
|
|
# exact for any quotient, modf / modff keep the sign of zero and split the
|
|
# infinities, atan2 has every signed-zero and infinity case, fpclassify /
|
|
# isnormal tell subnormals apart in both widths, hypot scales instead of
|
|
# overflowing or underflowing, frexp / ldexp handle subnormals and huge
|
|
# exponents, and asinh / acosh / atanh keep signs and do not overflow.
|
|
|
|
import math
|
|
import struct
|
|
import sys
|
|
from fractions import Fraction
|
|
|
|
F_QUIET = 0x00400000
|
|
D_QUIET = 0x0008000000000000
|
|
F_DNAN = 0x7FC00000
|
|
D_DNAN = 0x7FF8000000000000
|
|
F_INF = 0x7F800000
|
|
D_INF = 0x7FF0000000000000
|
|
F_SIGN = 0x80000000
|
|
D_SIGN = 0x8000000000000000
|
|
|
|
# The compare operand set: signed infinities and zeros, normals of both
|
|
# signs, a quiet NaN, a negative quiet NaN and a signaling-pattern NaN.
|
|
CMP_F = [0xFF800000, 0xC0200000, 0xBF800000, 0x80000000, 0x00000000,
|
|
0x3F800000, 0x41200000, 0x7149F2CA, 0x7F800000, 0x7FC00000,
|
|
0xFFC00000, 0x7F800001]
|
|
# The same values as doubles (1e30f widened exactly, so both tables share
|
|
# one truth table).
|
|
CMP_D = [0xFFF0000000000000, 0xC004000000000000, 0xBFF0000000000000,
|
|
0x8000000000000000, 0x0000000000000000, 0x3FF0000000000000,
|
|
0x4024000000000000, 0x46293E5940000000, 0x7FF0000000000000,
|
|
0x7FF8000000000000, 0xFFF8000000000000, 0x7FF0000000000001]
|
|
|
|
# Predicate order matches ALL() in the check: the six C operators, their
|
|
# negations (the unordered LLVM predicates uge/ugt/ule/ult), then the C99
|
|
# isunordered / !isunordered / islessgreater / !islessgreater.
|
|
PREDS = [
|
|
lambda a, b: a < b,
|
|
lambda a, b: a <= b,
|
|
lambda a, b: a > b,
|
|
lambda a, b: a >= b,
|
|
lambda a, b: a == b,
|
|
lambda a, b: a != b,
|
|
lambda a, b: not (a < b),
|
|
lambda a, b: not (a <= b),
|
|
lambda a, b: not (a > b),
|
|
lambda a, b: not (a >= b),
|
|
lambda a, b: math.isnan(a) or math.isnan(b),
|
|
lambda a, b: not (math.isnan(a) or math.isnan(b)),
|
|
lambda a, b: a < b or a > b,
|
|
lambda a, b: not (a < b or a > b),
|
|
]
|
|
|
|
|
|
def fval(bits):
|
|
return struct.unpack('<f', struct.pack('<I', bits))[0]
|
|
|
|
|
|
def dval(bits):
|
|
return struct.unpack('<d', struct.pack('<Q', bits))[0]
|
|
|
|
|
|
def fbits(v):
|
|
return struct.unpack('<I', struct.pack('<f', v))[0]
|
|
|
|
|
|
def dbits(v):
|
|
return struct.unpack('<Q', struct.pack('<d', v))[0]
|
|
|
|
|
|
def emitCmp():
|
|
for f, d in zip(CMP_F, CMP_D):
|
|
assert (math.isnan(fval(f)) and math.isnan(dval(d))) or fval(f) == dval(d)
|
|
bits = []
|
|
for vals in ([fval(b) for b in CMP_F], [dval(b) for b in CMP_D]):
|
|
for p in PREDS:
|
|
for a in vals:
|
|
for b in vals:
|
|
bits.append(1 if p(a, b) else 0)
|
|
out = []
|
|
for i in range(0, len(bits), 8):
|
|
byte = 0
|
|
for k, bit in enumerate(bits[i:i + 8]):
|
|
byte |= bit << k
|
|
out.append(byte)
|
|
print('#define NV %d' % len(CMP_F))
|
|
print('#define NP %d' % len(PREDS))
|
|
print('static const u32 kFBits[NV] = {')
|
|
print(' ' + ', '.join('0x%08XUL' % b for b in CMP_F))
|
|
print('};')
|
|
print('static const u64 kDBits[NV] = {')
|
|
for i in range(0, len(CMP_D), 3):
|
|
print(' ' + ', '.join('0x%016XULL' % b for b in CMP_D[i:i + 3]) + ',')
|
|
print('};')
|
|
print('static const u8 kExp[%d] = {' % len(out))
|
|
for i in range(0, len(out), 16):
|
|
print(' ' + ', '.join('0x%02X' % b for b in out[i:i + 16]) + ',')
|
|
print('};')
|
|
|
|
|
|
# ---- Special-value model ------------------------------------------------
|
|
|
|
def fclass(b):
|
|
e = (b >> 23) & 0xFF
|
|
m = b & 0x7FFFFF
|
|
if e == 0:
|
|
return 'zero'
|
|
if e == 0xFF:
|
|
return 'inf' if m == 0 else 'nan'
|
|
return 'normal'
|
|
|
|
|
|
def dclass(b):
|
|
e = (b >> 52) & 0x7FF
|
|
m = b & 0xFFFFFFFFFFFFF
|
|
if e == 0:
|
|
return 'zero'
|
|
if e == 0x7FF:
|
|
return 'inf' if m == 0 else 'nan'
|
|
return 'normal'
|
|
|
|
|
|
def binop(isD, op, a, b):
|
|
cls = dclass if isD else fclass
|
|
quiet, dnan, inf, sign = ((D_QUIET, D_DNAN, D_INF, D_SIGN) if isD else
|
|
(F_QUIET, F_DNAN, F_INF, F_SIGN))
|
|
val, tobits = (dval, dbits) if isD else (fval, fbits)
|
|
if op == 'sub':
|
|
# __subXf3(a, b) is __addXf3(a, b with its sign flipped), as in
|
|
# compiler-rt, so a NaN b comes back with its sign flipped too.
|
|
b ^= sign
|
|
op = 'add'
|
|
ca, cb = cls(a), cls(b)
|
|
if ca == 'nan':
|
|
return a | quiet
|
|
if cb == 'nan':
|
|
return b | quiet
|
|
s = (a ^ b) & sign
|
|
if op == 'add':
|
|
if ca == 'inf':
|
|
return dnan if cb == 'inf' and s else a
|
|
if cb == 'inf':
|
|
return b
|
|
if ca == 'zero' and cb == 'zero':
|
|
return a & b & sign
|
|
elif op == 'mul':
|
|
if ca == 'inf' or cb == 'inf':
|
|
return dnan if 'zero' in (ca, cb) else s | inf
|
|
if ca == 'zero' or cb == 'zero':
|
|
return s
|
|
elif op == 'div':
|
|
if ca == 'inf':
|
|
return dnan if cb == 'inf' else s | inf
|
|
if cb == 'inf':
|
|
return s
|
|
if cb == 'zero':
|
|
return dnan if ca == 'zero' else s | inf
|
|
if ca == 'zero':
|
|
return s
|
|
x, y = val(a), val(b)
|
|
if op == 'add':
|
|
r = x + y
|
|
elif op == 'mul':
|
|
r = x * y
|
|
else:
|
|
r = x / y
|
|
if not isD and abs(r) > 3.4028234663852886e38:
|
|
r = math.copysign(math.inf, r)
|
|
got = tobits(r)
|
|
assert val(got) == r or math.isinf(r), (op, hex(a), hex(b), r)
|
|
return got
|
|
|
|
|
|
def f2d(a):
|
|
if fclass(a) == 'nan':
|
|
return ((a & F_SIGN) << 32) | D_INF | ((a & 0x7FFFFF) << 29) | D_QUIET
|
|
return dbits(fval(a))
|
|
|
|
|
|
def d2f(a):
|
|
if dclass(a) == 'nan':
|
|
return ((a >> 32) & F_SIGN) | F_DNAN | (((a & 0xFFFFFFFFFFFFF) >> 29) & 0x7FFFFF)
|
|
v = dval(a)
|
|
if abs(v) > 3.4028234663852886e38:
|
|
return fbits(math.copysign(math.inf, v))
|
|
got = fbits(v)
|
|
assert fval(got) == v
|
|
return got
|
|
|
|
|
|
def toInt(v, unsigned):
|
|
if math.isnan(v):
|
|
return 0
|
|
lo, hi = (0, 0xFFFFFFFF) if unsigned else (-0x80000000, 0x7FFFFFFF)
|
|
if math.isinf(v):
|
|
r = hi if v > 0 else lo
|
|
else:
|
|
r = max(lo, min(hi, int(v)))
|
|
return r & 0xFFFFFFFF
|
|
|
|
|
|
def classify(v, bitsSign):
|
|
nan = 1 if math.isnan(v) else 0
|
|
inf = 1 if math.isinf(v) else 0
|
|
fin = 1 if not (nan or inf) else 0
|
|
cls = 0 if nan else 1 if inf else 4 if v == 0 else 2
|
|
return nan | (inf << 1) | (fin << 2) | (bitsSign << 3) | (cls << 4)
|
|
|
|
|
|
# Operation codes, matching the switch in the check's run().
|
|
OPS = ['fadd', 'fsub', 'fmul', 'fdiv', 'dadd', 'dsub', 'dmul', 'ddiv',
|
|
'f2d', 'd2f', 'f2i', 'f2u', 'd2i', 'd2u', 'fcls', 'dcls',
|
|
'sqrt', 'log', 'asin', 'acos', 'atanh',
|
|
'floor', 'ceil', 'trunc', 'round', 'rint', 'nearbyint',
|
|
'floorf', 'ceilf', 'truncf', 'roundf', 'rintf',
|
|
'exp', 'exp2', 'expm1', 'sinh', 'cosh', 'tanh', 'log1p', 'cbrt',
|
|
'pow', 'fmod', 'fdim', 'hypot', 'remainder',
|
|
'nearbyintf', 'fmaxf', 'fminf']
|
|
|
|
FQ = 0x7FC00123 # quiet NaN with a payload
|
|
FS = 0x7F800005 # signaling-pattern NaN
|
|
FNQ = 0xFFC00000 # negative quiet NaN
|
|
DQ = 0x7FF8000000000123
|
|
DS = 0x7FF0000000000005
|
|
DNQ = 0xFFF8000000000000
|
|
F1, FM1, F15, F2, F225, FM2, F3E38 = (fbits(1.0), fbits(-1.0), fbits(1.5),
|
|
fbits(2.0), fbits(2.25), fbits(-2.0),
|
|
fbits(3e38))
|
|
D1, DM1, D15, D2, D225, DM2, D1E300 = (dbits(1.0), dbits(-1.0), dbits(1.5),
|
|
dbits(2.0), dbits(2.25), dbits(-2.0),
|
|
dbits(1e300))
|
|
FZ, FNZ, FI, FNI = 0, F_SIGN, F_INF, F_INF | F_SIGN
|
|
DZ, DNZ, DI, DNI = 0, D_SIGN, D_INF, D_INF | D_SIGN
|
|
|
|
|
|
def fpairs():
|
|
return [(FQ, F1), (F1, FS), (FS, FQ), (FNQ, F2), (FZ, FQ), (FI, F1),
|
|
(FI, FI), (FI, FNI), (FNI, F2), (F2, FNI), (F1, FM1), (FM1, F1),
|
|
(FNZ, FNZ), (FNZ, FZ), (FZ, FNZ), (FZ, FZ), (F225, F15),
|
|
(F3E38, F3E38), (FZ, FI), (FI, FZ), (FNI, FM2), (FNI, FNI),
|
|
(FNZ, F2), (FZ, FZ), (F1, FZ), (FM1, FZ), (F1, FNZ), (F1, FI),
|
|
(FM1, FI), (FI, F2), (F15, F2)]
|
|
|
|
|
|
def dpairs():
|
|
return [(DQ, D1), (D1, DS), (DS, DQ), (DNQ, D2), (DZ, DQ), (DI, D1),
|
|
(DI, DI), (DI, DNI), (DNI, D2), (D2, DNI), (D1, DM1), (DM1, D1),
|
|
(DNZ, DNZ), (DNZ, DZ), (DZ, DNZ), (DZ, DZ), (D225, D15),
|
|
(D1E300, D1E300), (DZ, DI), (DI, DZ), (DNI, DM2), (DNI, DNI),
|
|
(DNZ, D2), (DZ, DZ), (D1, DZ), (DM1, DZ), (D1, DNZ), (D1, DI),
|
|
(DM1, DI), (DI, D2), (D15, D2)]
|
|
|
|
|
|
def emitSpecial():
|
|
cases = []
|
|
for isD, prefix, pairs in ((False, 'f', fpairs()), (True, 'd', dpairs())):
|
|
for op in ('add', 'sub', 'mul', 'div'):
|
|
for a, b in pairs:
|
|
cases.append((OPS.index(prefix + op), a, b, binop(isD, op, a, b)))
|
|
for a in (FQ, FS, 0xFFC00001, FI, FNI, FNZ, FZ, F15, FM2):
|
|
cases.append((OPS.index('f2d'), a, 0, f2d(a)))
|
|
for a in (DQ, 0xFFF8000020000000, DS, DI, DNI, DNZ, DZ, D15, DM2, D1E300,
|
|
dbits(-1e300)):
|
|
cases.append((OPS.index('d2f'), a, 0, d2f(a)))
|
|
fconv = (FQ, FNQ, FS, FI, FNI, fbits(3e9), fbits(-3e9), fbits(5e9),
|
|
fbits(-2.5), fbits(3.5), FM1, FNZ)
|
|
dconv = (DQ, DNQ, DS, DI, DNI, dbits(3e9), dbits(-3e9), dbits(5e9),
|
|
dbits(-2.5), dbits(3.5), DM1, DNZ)
|
|
for a in fconv:
|
|
cases.append((OPS.index('f2i'), a, 0, toInt(fval(a), False)))
|
|
cases.append((OPS.index('f2u'), a, 0, toInt(fval(a), True)))
|
|
for a in dconv:
|
|
cases.append((OPS.index('d2i'), a, 0, toInt(dval(a), False)))
|
|
cases.append((OPS.index('d2u'), a, 0, toInt(dval(a), True)))
|
|
for a in (FQ, FNQ, FS, FI, FNI, FZ, FNZ, F15, FM2):
|
|
cases.append((OPS.index('fcls'), a, 0, classify(fval(a), a >> 31)))
|
|
for a in (DQ, DNQ, DS, DI, DNI, DZ, DNZ, D15, DM2):
|
|
cases.append((OPS.index('dcls'), a, 0, classify(dval(a), a >> 63)))
|
|
# math.c's domain edges. The check folds every NaN result to the
|
|
# default NaN, because how a NaN travels through a series (and so its
|
|
# sign) is not part of these functions' contract.
|
|
halfPi = dbits(math.pi / 2)
|
|
mathCases = [
|
|
('sqrt', dbits(-1.0), D_DNAN), ('sqrt', DI, DI), ('sqrt', DNI, D_DNAN),
|
|
('sqrt', DQ, D_DNAN), ('sqrt', dbits(4.0), dbits(2.0)),
|
|
('sqrt', DNZ, DNZ),
|
|
('log', DI, DI), ('log', DNI, D_DNAN), ('log', DQ, D_DNAN),
|
|
('log', DM1, D_DNAN), ('log', DZ, DNI), ('log', D1, DZ),
|
|
('asin', D2, D_DNAN), ('asin', DM2, D_DNAN), ('asin', D1, halfPi),
|
|
('asin', DM1, halfPi | D_SIGN), ('asin', DQ, D_DNAN),
|
|
('acos', D2, D_DNAN), ('acos', D1, DZ),
|
|
('atanh', D2, D_DNAN), ('atanh', DM2, D_DNAN), ('atanh', D1, DI),
|
|
('atanh', DM1, DNI), ('atanh', DQ, D_DNAN),
|
|
]
|
|
cases = [c + (0,) for c in cases]
|
|
for name, a, want in mathCases:
|
|
cases.append((OPS.index(name), a, 0, want, 0))
|
|
# Two tables, one per build of the check: both together (with the math
|
|
# functions they pull in) overflow a bank-0 test binary.
|
|
for part, table in ((1, cases), (2, annexFCases())):
|
|
print('#if FP_TABLE == 1' if part == 1 else '#else')
|
|
print('#define NC %d' % len(table))
|
|
print('static const CaseT kCases[NC] = {')
|
|
for op, a, b, want, approx in table:
|
|
print(' { %-12s 0x%016XULL, 0x%016XULL, 0x%016XULL, %d },' %
|
|
('O_' + OPS[op].upper() + ',', a, b, want, approx))
|
|
print('};')
|
|
print('#endif')
|
|
|
|
|
|
# ---- math.c Annex F cases ----------------------------------------------
|
|
|
|
INF = math.inf
|
|
NAN = math.nan
|
|
|
|
|
|
def signedZero(r, x):
|
|
# An integral result of zero keeps x's sign (floor(-0.0) is -0.0,
|
|
# ceil(-0.5) is -0.0).
|
|
return math.copysign(0.0, x) if r == 0 else r
|
|
|
|
|
|
def cfloor(x):
|
|
return x if not math.isfinite(x) else signedZero(float(math.floor(x)), x)
|
|
|
|
|
|
def cceil(x):
|
|
return x if not math.isfinite(x) else signedZero(float(math.ceil(x)), x)
|
|
|
|
|
|
def ctrunc(x):
|
|
return x if not math.isfinite(x) else signedZero(float(math.trunc(x)), x)
|
|
|
|
|
|
def cround(x):
|
|
# C round: half away from zero, computed exactly.
|
|
if not math.isfinite(x):
|
|
return x
|
|
r = float(math.floor(Fraction(abs(x)) + Fraction(1, 2)))
|
|
return math.copysign(r, x)
|
|
|
|
|
|
def crint(x):
|
|
# Round to nearest, ties to even (Python's round() on a float is exact).
|
|
if not math.isfinite(x):
|
|
return x
|
|
return math.copysign(float(round(x)), x)
|
|
|
|
|
|
def guard(f, *args):
|
|
# Python raises where C returns: a domain error is NaN, an overflow is
|
|
# an infinity of the true result's sign.
|
|
try:
|
|
return f(*args)
|
|
except ValueError:
|
|
return NAN
|
|
except OverflowError:
|
|
return INF
|
|
|
|
|
|
def isOddInt(y):
|
|
return math.isfinite(y) and y == math.floor(y) and abs(y) < 2.0 ** 53 and int(y) % 2 == 1
|
|
|
|
|
|
def cpow(x, y):
|
|
# C99 F.9.4.4.
|
|
if y == 0 or x == 1.0:
|
|
return 1.0
|
|
if math.isnan(x) or math.isnan(y):
|
|
return NAN
|
|
if math.isinf(y):
|
|
if abs(x) == 1:
|
|
return 1.0
|
|
return INF if (abs(x) > 1) == (y > 0) else 0.0
|
|
if x == 0:
|
|
r = INF if y < 0 else 0.0
|
|
return math.copysign(r, x) if isOddInt(y) else r
|
|
if math.isinf(x):
|
|
r = INF if y > 0 else 0.0
|
|
return -r if (x < 0 and isOddInt(y)) else r
|
|
if x < 0 and y != math.floor(y):
|
|
return NAN
|
|
r = guard(math.pow, x, y)
|
|
# guard() cannot know an overflow's sign: a negative x to an odd power.
|
|
return -r if (math.isinf(r) and x < 0 and isOddInt(y)) else r
|
|
|
|
|
|
def cfmod(x, y):
|
|
# C99 F.9.7.1.
|
|
if math.isnan(x) or math.isnan(y) or math.isinf(x) or y == 0:
|
|
return NAN
|
|
if math.isinf(y):
|
|
return x
|
|
return math.fmod(x, y)
|
|
|
|
|
|
def cfdim(x, y):
|
|
if math.isnan(x) or math.isnan(y):
|
|
return NAN
|
|
return x - y if x > y else 0.0
|
|
|
|
|
|
def cremainder(x, y):
|
|
# C99 F.9.7.2 (Python's math.remainder raises for the NaN cases).
|
|
if math.isnan(x) or math.isnan(y) or math.isinf(x) or y == 0:
|
|
return NAN
|
|
return math.remainder(x, y)
|
|
|
|
|
|
def ccbrt(x):
|
|
if not math.isfinite(x) or x == 0:
|
|
return x
|
|
return math.copysign(abs(x) ** (1.0 / 3.0), x)
|
|
|
|
|
|
def fround32(v):
|
|
return fval(fbits(v))
|
|
|
|
|
|
MATH1 = {
|
|
'floor': cfloor, 'ceil': cceil, 'trunc': ctrunc, 'round': cround,
|
|
'rint': crint, 'nearbyint': crint,
|
|
'exp': lambda x: guard(math.exp, x),
|
|
'exp2': lambda x: guard(math.pow, 2.0, x),
|
|
'expm1': lambda x: guard(math.expm1, x),
|
|
'sinh': lambda x: guard(math.sinh, x) if not (math.isinf(x) or abs(x) > 700) else math.copysign(INF, x),
|
|
'cosh': lambda x: guard(math.cosh, x),
|
|
'tanh': math.tanh,
|
|
'log1p': lambda x: -INF if x == -1 else guard(math.log1p, x),
|
|
'cbrt': ccbrt,
|
|
}
|
|
MATH2 = {'pow': cpow, 'fmod': cfmod, 'fdim': cfdim, 'hypot': math.hypot,
|
|
'remainder': cremainder}
|
|
FLOATOPS = {'floorf': cfloor, 'ceilf': cceil, 'truncf': ctrunc,
|
|
'roundf': cround, 'rintf': crint, 'nearbyintf': crint}
|
|
|
|
|
|
def cfmax(x, y):
|
|
# C99 F.9.9.2: a NaN operand is ignored in favour of the other.
|
|
if math.isnan(x):
|
|
return y
|
|
if math.isnan(y):
|
|
return x
|
|
return max(x, y)
|
|
|
|
|
|
def cfmin(x, y):
|
|
if math.isnan(x):
|
|
return y
|
|
if math.isnan(y):
|
|
return x
|
|
return min(x, y)
|
|
|
|
|
|
FLOAT2 = {'fmaxf': cfmax, 'fminf': cfmin}
|
|
|
|
# (op, x[, y], approx). Every function gets its special values and a few
|
|
# ordinary ones, so a normal-path regression shows too.
|
|
MATH_CASES = [
|
|
('floor', 2.7, 0), ('floor', -2.7, 0), ('floor', -0.0, 0), ('floor', 0.0, 0),
|
|
('floor', -0.3, 0), ('floor', 0.3, 0), ('floor', INF, 0), ('floor', -INF, 0),
|
|
('floor', NAN, 0), ('floor', 1e300, 0), ('floor', -0.5, 0),
|
|
('ceil', 2.3, 0), ('ceil', -2.3, 0), ('ceil', 0.0, 0), ('ceil', -0.0, 0),
|
|
('ceil', -0.5, 0), ('ceil', 0.5, 0), ('ceil', INF, 0), ('ceil', NAN, 0),
|
|
('trunc', -2.7, 0), ('trunc', 2.7, 0), ('trunc', -0.5, 0), ('trunc', -0.0, 0),
|
|
('trunc', INF, 0),
|
|
('round', 2.5, 0), ('round', -2.5, 0), ('round', 0.49999999999999994, 0),
|
|
('round', -0.49999999999999994, 0), ('round', 4503599627370497.0, 0),
|
|
('round', -0.0, 0), ('round', 0.5, 0), ('round', 1.5, 0), ('round', INF, 0),
|
|
('round', NAN, 0), ('round', -0.4, 0),
|
|
('rint', 2.5, 0), ('rint', 3.5, 0), ('rint', -2.5, 0), ('rint', -0.3, 0),
|
|
('rint', 0.3, 0), ('rint', 1.5, 0), ('rint', 4503599627370497.0, 0),
|
|
('rint', INF, 0), ('rint', NAN, 0), ('rint', -0.0, 0),
|
|
('nearbyint', 2.5, 0), ('nearbyint', -3.5, 0),
|
|
('floorf', -0.0, 0), ('floorf', 2.75, 0), ('floorf', -2.75, 0),
|
|
('ceilf', 0.0, 0), ('ceilf', -0.5, 0), ('truncf', -2.5, 0),
|
|
('roundf', 2.5, 0), ('roundf', -0.0, 0),
|
|
('rintf', 2.5, 0), ('rintf', 3.5, 0), ('rintf', -0.25, 0),
|
|
('nearbyintf', 2.5, 0), ('nearbyintf', -1.5, 0),
|
|
('fmaxf', 1.0, NAN, 0), ('fmaxf', NAN, -2.0, 0), ('fmaxf', -0.5, 2.0, 0),
|
|
('fminf', 1.0, NAN, 0), ('fminf', -0.5, 2.0, 0), ('fminf', 3.0, -INF, 0),
|
|
('exp', NAN, 0), ('exp', INF, 0), ('exp', -INF, 0), ('exp', 1e5, 0),
|
|
('exp', -1e5, 0), ('exp', 1024.0, 0), ('exp', -1024.0, 0), ('exp', 800.0, 0),
|
|
('exp', -800.0, 0), ('exp', 0.0, 0), ('exp', -0.0, 0), ('exp', 1.5, 1),
|
|
('exp', -2.0, 1),
|
|
('exp2', 3.0, 1), ('exp2', NAN, 0), ('exp2', INF, 0), ('exp2', -INF, 0),
|
|
('expm1', -INF, 0), ('expm1', INF, 0), ('expm1', NAN, 0), ('expm1', -0.0, 0),
|
|
('expm1', 0.5, 1),
|
|
('sinh', INF, 0), ('sinh', -INF, 0), ('sinh', NAN, 0), ('sinh', -0.0, 0),
|
|
('sinh', 800.0, 0), ('sinh', -800.0, 0), ('sinh', 1.0, 1), ('sinh', -1.0, 1),
|
|
('cosh', INF, 0), ('cosh', -INF, 0), ('cosh', NAN, 0), ('cosh', 0.0, 0),
|
|
('cosh', 1.0, 1),
|
|
('tanh', INF, 0), ('tanh', -INF, 0), ('tanh', NAN, 0), ('tanh', -0.0, 0),
|
|
('tanh', 40.0, 0), ('tanh', -40.0, 0), ('tanh', 400.0, 0), ('tanh', 0.5, 1),
|
|
('tanh', -0.5, 1),
|
|
('log1p', -0.0, 0), ('log1p', 0.0, 0), ('log1p', -1.0, 0), ('log1p', INF, 0),
|
|
('log1p', NAN, 0), ('log1p', 1.0, 1),
|
|
('cbrt', INF, 0), ('cbrt', -INF, 0), ('cbrt', NAN, 0), ('cbrt', -0.0, 0),
|
|
('cbrt', 27.0, 1),
|
|
('pow', NAN, 0.0, 0), ('pow', NAN, -0.0, 0), ('pow', 1.0, NAN, 0),
|
|
('pow', 1.0, INF, 0), ('pow', -1.0, INF, 0), ('pow', -1.0, -INF, 0),
|
|
('pow', 2.0, NAN, 0), ('pow', NAN, 2.0, 0), ('pow', 0.0, -3.0, 0),
|
|
('pow', -0.0, -3.0, 0), ('pow', -0.0, -2.0, 0), ('pow', -0.0, -INF, 0),
|
|
('pow', 0.0, -0.5, 0), ('pow', -0.0, 3.0, 0), ('pow', -0.0, 2.0, 0),
|
|
('pow', 0.0, 0.5, 0), ('pow', 0.5, -INF, 0), ('pow', 2.0, -INF, 0),
|
|
('pow', 0.5, INF, 0), ('pow', 2.0, INF, 0), ('pow', -INF, -3.0, 0),
|
|
('pow', -INF, -2.0, 0), ('pow', -INF, 3.0, 0), ('pow', -INF, 2.0, 0),
|
|
('pow', -INF, 0.5, 0), ('pow', INF, -1.0, 0), ('pow', INF, 2.0, 0),
|
|
('pow', -8.0, 1.0 / 3.0, 0), ('pow', -2.0, 0.5, 0), ('pow', 2.0, 10.0, 0),
|
|
('pow', 2.0, -3.0, 0), ('pow', -2.0, 3.0, 0), ('pow', -2.0, 2.0, 0),
|
|
('pow', 10.0, 20.0, 0), ('pow', -2.0, 1e20, 0), ('pow', -2.0, 40001.0, 0),
|
|
('pow', 2.0, 0.5, 1), ('pow', 2.0, 1.5, 1), ('pow', 10.0, 0.3, 1),
|
|
('pow', 1.0000001, 100000.0, 1),
|
|
('fmod', 7.5, 2.0, 0), ('fmod', -7.5, 2.0, 0), ('fmod', -4.0, 2.0, 0),
|
|
('fmod', 4.0, 2.0, 0), ('fmod', 5.0, 3.0, 0), ('fmod', 1.0, INF, 0),
|
|
('fmod', -1.0, -INF, 0), ('fmod', 1.0, 0.0, 0), ('fmod', 1.0, -0.0, 0),
|
|
('fmod', INF, 2.0, 0), ('fmod', -INF, 2.0, 0), ('fmod', NAN, 2.0, 0),
|
|
('fmod', 2.0, NAN, 0), ('fmod', 0.0, 5.0, 0), ('fmod', -0.0, 5.0, 0),
|
|
('fmod', 0.0, 0.0, 0), ('fmod', INF, INF, 0),
|
|
('fdim', NAN, 1.0, 0), ('fdim', 1.0, NAN, 0), ('fdim', 5.0, 3.0, 0),
|
|
('fdim', 3.0, 5.0, 0), ('fdim', INF, 1.0, 0), ('fdim', 1.0, INF, 0),
|
|
('fdim', INF, INF, 0),
|
|
('hypot', INF, NAN, 0), ('hypot', NAN, -INF, 0), ('hypot', -INF, 1.0, 0),
|
|
('hypot', NAN, 1.0, 0), ('hypot', 3.0, 4.0, 1), ('hypot', 0.0, 0.0, 0),
|
|
('remainder', 5.0, 2.0, 0), ('remainder', 7.0, 2.0, 0),
|
|
('remainder', 7.0, 4.0, 0), ('remainder', 1.0, INF, 0),
|
|
('remainder', -3.0, -INF, 0), ('remainder', 1.0, 0.0, 0),
|
|
('remainder', INF, 2.0, 0), ('remainder', NAN, 1.0, 0),
|
|
]
|
|
|
|
|
|
def nanFold(bits, isF):
|
|
# The check folds every NaN a math function returns to the default NaN.
|
|
v = fval(bits) if isF else dval(bits)
|
|
if math.isnan(v):
|
|
return F_DNAN if isF else D_DNAN
|
|
return bits
|
|
|
|
|
|
def annexFCases():
|
|
out = []
|
|
for c in MATH_CASES:
|
|
name = c[0]
|
|
if name in FLOAT2:
|
|
x, y, approx = fround32(c[1]), fround32(c[2]), c[3]
|
|
want = fround32(FLOAT2[name](x, y))
|
|
out.append((OPS.index(name), fbits(x), fbits(y), nanFold(fbits(want), True), approx))
|
|
continue
|
|
if name in MATH2:
|
|
x, y, approx = c[1], c[2], c[3]
|
|
want = MATH2[name](x, y)
|
|
a, b = dbits(x), dbits(y)
|
|
else:
|
|
x, approx = c[1], c[2]
|
|
y = 0.0
|
|
if name in FLOATOPS:
|
|
x = fround32(x)
|
|
want = fround32(FLOATOPS[name](x))
|
|
out.append((OPS.index(name), fbits(x), 0, nanFold(fbits(want), True), approx))
|
|
continue
|
|
want = MATH1[name](x)
|
|
a, b = dbits(x), 0
|
|
out.append((OPS.index(name), a, b, nanFold(dbits(want), False), approx))
|
|
return out
|
|
|
|
|
|
# ---- math.c edge table -------------------------------------------------
|
|
|
|
EDGE_OPS = ['modfFrac', 'modfInt', 'modffFrac', 'modffInt', 'atan2',
|
|
'clsD', 'clsF', 'hypot', 'fmod', 'remainder', 'frexpM', 'frexpE',
|
|
'ldexp', 'asinh', 'acosh', 'atanh']
|
|
|
|
FP_NAN_CLS, FP_INFINITE_CLS, FP_NORMAL_CLS, FP_SUBNORMAL_CLS, FP_ZERO_CLS = 0, 1, 2, 3, 4
|
|
D_MIN_SUB = 5e-324
|
|
F_MIN_SUB_BITS = 0x00000001
|
|
F_MAX_SUB_BITS = 0x007FFFFF
|
|
D_MAX_SUB = dval(0x000FFFFFFFFFFFFF)
|
|
|
|
|
|
def clsBits(bits, expBits, mantBits):
|
|
e = (bits >> mantBits) & ((1 << expBits) - 1)
|
|
m = bits & ((1 << mantBits) - 1)
|
|
if e == (1 << expBits) - 1:
|
|
c = FP_NAN_CLS if m else FP_INFINITE_CLS
|
|
elif e:
|
|
c = FP_NORMAL_CLS
|
|
else:
|
|
c = FP_SUBNORMAL_CLS if m else FP_ZERO_CLS
|
|
return c | ((1 if c == FP_NORMAL_CLS else 0) << 4)
|
|
|
|
|
|
def cmodf(x):
|
|
if math.isnan(x):
|
|
return (NAN, NAN)
|
|
return math.modf(x)
|
|
|
|
|
|
def cfmodExact(x, y):
|
|
if math.isnan(x) or math.isnan(y) or math.isinf(x) or y == 0:
|
|
return NAN
|
|
if math.isinf(y):
|
|
return x
|
|
return math.fmod(x, y)
|
|
|
|
|
|
def catan2(y, x):
|
|
return NAN if math.isnan(x) or math.isnan(y) else math.atan2(y, x)
|
|
|
|
|
|
def chypot(x, y):
|
|
if math.isinf(x) or math.isinf(y):
|
|
return INF
|
|
return math.hypot(x, y)
|
|
|
|
|
|
def cldexp(x, n):
|
|
try:
|
|
return math.ldexp(x, n)
|
|
except OverflowError:
|
|
return math.copysign(INF, x)
|
|
|
|
|
|
def cacosh(x):
|
|
return NAN if math.isnan(x) or x < 1 else math.acosh(x)
|
|
|
|
|
|
def catanh(x):
|
|
if math.isnan(x) or abs(x) > 1:
|
|
return NAN
|
|
if abs(x) == 1:
|
|
return math.copysign(INF, x)
|
|
return math.atanh(x)
|
|
|
|
|
|
# (op, a, b, approx); a and b are doubles, or floats for the float ops.
|
|
EDGE_CASES = [
|
|
('modfFrac', -2.0, 0, 0), ('modfInt', -2.0, 0, 0), ('modfFrac', -0.0, 0, 0),
|
|
('modfInt', -0.0, 0, 0), ('modfFrac', INF, 0, 0), ('modfInt', INF, 0, 0),
|
|
('modfFrac', -INF, 0, 0), ('modfInt', -INF, 0, 0), ('modfFrac', NAN, 0, 0),
|
|
('modfInt', NAN, 0, 0), ('modfFrac', -3.75, 0, 0), ('modfInt', -3.75, 0, 0),
|
|
('modfFrac', 0.5, 0, 0), ('modfInt', 0.5, 0, 0), ('modfFrac', 1e300, 0, 0),
|
|
('modfInt', 4503599627370497.5, 0, 0), ('modfFrac', 4503599627370497.5, 0, 0),
|
|
('modfFrac', D_MIN_SUB, 0, 0), ('modfInt', -D_MIN_SUB, 0, 0),
|
|
('modffFrac', -2.0, 0, 0), ('modffInt', -2.0, 0, 0), ('modffFrac', -0.0, 0, 0),
|
|
('modffFrac', INF, 0, 0), ('modffInt', -INF, 0, 0), ('modffFrac', NAN, 0, 0),
|
|
('modffInt', NAN, 0, 0), ('modffFrac', -3.75, 0, 0), ('modffInt', -3.75, 0, 0),
|
|
('modffFrac', 8388609.5, 0, 0), ('modffFrac', 1e-45, 0, 0),
|
|
('atan2', 0.0, 0.0, 0), ('atan2', 0.0, -0.0, 0), ('atan2', -0.0, 0.0, 0),
|
|
('atan2', -0.0, -0.0, 0), ('atan2', 0.0, -1.0, 0), ('atan2', -0.0, -1.0, 0),
|
|
('atan2', 0.0, 1.0, 0), ('atan2', -0.0, 1.0, 0), ('atan2', 1.0, 0.0, 0),
|
|
('atan2', -1.0, -0.0, 0), ('atan2', 1.0, -INF, 0), ('atan2', -1.0, -INF, 0),
|
|
('atan2', 1.0, INF, 0), ('atan2', -1.0, INF, 0), ('atan2', INF, 1.0, 0),
|
|
('atan2', -INF, -1.0, 0), ('atan2', INF, -INF, 0), ('atan2', -INF, -INF, 0),
|
|
('atan2', INF, INF, 0), ('atan2', -INF, INF, 0), ('atan2', NAN, 1.0, 0),
|
|
('atan2', 1.0, NAN, 0), ('atan2', 1.0, 1.0, 1), ('atan2', -1.0, -2.0, 1),
|
|
('atan2', 2.0, -1.0, 1), ('atan2', 1e300, 1e-300, 0), ('atan2', 1e-300, -1e300, 0),
|
|
('clsD', 1.0, 0, 0), ('clsD', -0.0, 0, 0), ('clsD', D_MIN_SUB, 0, 0),
|
|
('clsD', -D_MAX_SUB, 0, 0), ('clsD', 2.2250738585072014e-308, 0, 0),
|
|
('clsD', INF, 0, 0), ('clsD', NAN, 0, 0), ('clsD', -1e300, 0, 0),
|
|
('clsF', 1.0, 0, 0), ('clsF', -0.0, 0, 0), ('clsF', fval(F_MIN_SUB_BITS), 0, 0),
|
|
('clsF', -fval(F_MAX_SUB_BITS), 0, 0), ('clsF', 1.1754943508222875e-38, 0, 0),
|
|
('clsF', INF, 0, 0), ('clsF', NAN, 0, 0), ('clsF', -3e38, 0, 0),
|
|
('hypot', 1e200, 1e200, 1), ('hypot', 3e200, -4e200, 1), ('hypot', 1e-200, 1e-200, 1),
|
|
('hypot', -3e-300, 4e-300, 1), ('hypot', 1e308, 1e308, 0), ('hypot', 1e300, 1.0, 0),
|
|
('hypot', D_MIN_SUB, D_MIN_SUB, 0), ('hypot', -0.0, -0.0, 0), ('hypot', -0.0, -5.0, 0),
|
|
('hypot', -INF, NAN, 0), ('hypot', NAN, 2.0, 0), ('hypot', 3.0, 4.0, 1),
|
|
('fmod', 1e300, 3.0, 0), ('fmod', -1e300, 7.0, 0), ('fmod', 1e20, 0.1, 0),
|
|
('fmod', 5.5, -2.0, 0), ('fmod', -5.5, -2.0, 0), ('fmod', -6.0, 3.0, 0),
|
|
('fmod', 6.0, -3.0, 0), ('fmod', 1.0, 3.0, 0), ('fmod', 3 * D_MIN_SUB, 2 * D_MIN_SUB, 0),
|
|
('fmod', 1.5, D_MIN_SUB, 0), ('fmod', 1.7976931348623157e308, 1e-300, 0),
|
|
('fmod', 2.0 ** 60 + 2.0 ** 8, 2.0 ** -40, 0), ('fmod', 7.5, 2.0, 0),
|
|
('fmod', -0.0, 3.0, 0), ('fmod', 1.0, -INF, 0), ('fmod', INF, 1.0, 0),
|
|
('fmod', NAN, 1.0, 0), ('fmod', 1.0, 0.0, 0),
|
|
('remainder', 1e300, 3.0, 0), ('remainder', -1e300, 7.0, 0), ('remainder', 1e20, 0.1, 0),
|
|
('remainder', -4.0, 2.0, 0), ('remainder', 4.0, -2.0, 0), ('remainder', 3.0, 2.0, 0),
|
|
('remainder', -3.0, 2.0, 0), ('remainder', 5.0, -2.0, 0), ('remainder', 1.5, 1e-300, 0),
|
|
('remainder', 3 * D_MIN_SUB, 2 * D_MIN_SUB, 0), ('remainder', 1.7976931348623157e308, -1.7e308, 0),
|
|
('remainder', 1.0, INF, 0), ('remainder', INF, 1.0, 0), ('remainder', 1.0, -0.0, 0),
|
|
('frexpM', D_MIN_SUB, 0, 0), ('frexpE', D_MIN_SUB, 0, 0), ('frexpM', -D_MAX_SUB, 0, 0),
|
|
('frexpE', -D_MAX_SUB, 0, 0), ('frexpM', 1e300, 0, 0), ('frexpE', 1e300, 0, 0),
|
|
('frexpM', -0.0, 0, 0), ('frexpE', -0.0, 0, 0), ('frexpM', -INF, 0, 0),
|
|
('frexpM', 0.75, 0, 0), ('frexpE', 0.75, 0, 0),
|
|
('ldexp', 1.0, 1023, 0), ('ldexp', 1.0, 1024, 0), ('ldexp', -1.0, 32767, 0),
|
|
('ldexp', 1.0, -32768, 0), ('ldexp', 0.0, 30000, 0), ('ldexp', -0.0, -30000, 0),
|
|
('ldexp', 1.0, -1074, 0), ('ldexp', 1.5, -1074, 0), ('ldexp', 1.0, -1075, 0),
|
|
('ldexp', 3.0, -1075, 0), ('ldexp', D_MIN_SUB, 1074, 0), ('ldexp', -D_MAX_SUB, 1, 0),
|
|
('ldexp', 0.75, -1022, 0), ('ldexp', INF, -5, 0), ('ldexp', NAN, 5, 0),
|
|
('asinh', -0.0, 0, 0), ('asinh', 0.0, 0, 0), ('asinh', -INF, 0, 0), ('asinh', INF, 0, 0),
|
|
('asinh', NAN, 0, 0), ('asinh', -1e10, 0, 1), ('asinh', 1e200, 0, 1), ('asinh', -1e300, 0, 1),
|
|
('asinh', 0.5, 0, 1), ('asinh', -2.0, 0, 1),
|
|
('acosh', 1e200, 0, 1), ('acosh', INF, 0, 0), ('acosh', 1.0, 0, 0), ('acosh', 0.5, 0, 0),
|
|
('acosh', -INF, 0, 0), ('acosh', NAN, 0, 0), ('acosh', 2.0, 0, 1),
|
|
('atanh', -0.0, 0, 0), ('atanh', 0.0, 0, 0), ('atanh', -0.5, 0, 1), ('atanh', 1.0, 0, 0),
|
|
]
|
|
|
|
EDGE_MODEL = {
|
|
'modfFrac': lambda x, y: cmodf(x)[0], 'modfInt': lambda x, y: cmodf(x)[1],
|
|
'atan2': catan2, 'hypot': chypot, 'fmod': cfmodExact, 'remainder': cremainder,
|
|
'frexpM': lambda x, y: x if not math.isfinite(x) else math.frexp(x)[0],
|
|
'ldexp': lambda x, n: x if not math.isfinite(x) else cldexp(x, n),
|
|
'asinh': lambda x, y: math.asinh(x), 'acosh': lambda x, y: cacosh(x),
|
|
'atanh': lambda x, y: catanh(x),
|
|
}
|
|
|
|
|
|
def edgeCase(name, x, y, approx):
|
|
op = EDGE_OPS.index(name)
|
|
if name in ('modffFrac', 'modffInt'):
|
|
xf = fround32(x)
|
|
r = cmodf(xf)[0 if name == 'modffFrac' else 1]
|
|
return (op, fbits(xf), 0, nanFold(fbits(fround32(r)), True), approx)
|
|
if name == 'clsF':
|
|
b = fbits(fround32(x)) if not math.isnan(x) else F_DNAN
|
|
return (op, b, 0, clsBits(b, 8, 23), approx)
|
|
if name == 'clsD':
|
|
b = dbits(x) if not math.isnan(x) else D_DNAN
|
|
return (op, b, 0, clsBits(b, 11, 52), approx)
|
|
if name == 'frexpE':
|
|
return (op, dbits(x), 0, math.frexp(x)[1] & 0xFFFF, approx)
|
|
if name == 'ldexp':
|
|
return (op, dbits(x), y & 0xFFFF, nanFold(dbits(EDGE_MODEL[name](x, y)), False), approx)
|
|
want = EDGE_MODEL[name](x, y)
|
|
return (op, dbits(x), dbits(y), nanFold(dbits(want), False), approx)
|
|
|
|
|
|
def emitEdge():
|
|
table = [edgeCase(*c) for c in EDGE_CASES]
|
|
print('#define NE %d' % len(table))
|
|
print('static const CaseT kEdge[NE] = {')
|
|
for op, a, b, want, approx in table:
|
|
print(' { %-14s 0x%016XULL, 0x%016XULL, 0x%016XULL, %d },' %
|
|
('E_' + EDGE_OPS[op].upper() + ',', a, b, want, approx))
|
|
print('};')
|
|
|
|
|
|
if __name__ == '__main__':
|
|
modes = {'cmp': emitCmp, 'special': emitSpecial, 'edge': emitEdge}
|
|
if len(sys.argv) != 2 or sys.argv[1] not in modes:
|
|
sys.exit('usage: genFpTables.py cmp|special|edge')
|
|
modes[sys.argv[1]]()
|