wyj
7c980c35aa
Generate plain long-double and experimental FMA variants from the production entry kernel. Compare exact-rational references, geometric validation, fallback counts and repeated timings with self-contained fixtures. Validate cached build dependencies and reject stale timing output. Count all unconfirmed entry outcomes independently of reference classification.
442 lines
17 KiB
Python
442 lines
17 KiB
Python
#!/usr/bin/env python3
|
|
"""Deterministic fixture generation for the quadratic precision benchmark.
|
|
|
|
All values are frozen as IEEE-754 hex literals in the generated C header and as
|
|
hex strings in cases.json, so the C probe and the Python oracle see bit-for-bit
|
|
identical inputs. No external catalog/observer/slab data is required.
|
|
|
|
Case families (kernel):
|
|
curated explicit adversarial / production fixtures
|
|
fixed fixed sphere, axis and oblique directions
|
|
moving constant sphere velocity (radial/transverse), inward photons
|
|
growing rr < 0 (radius grows along the past parameter), both inward
|
|
(toward-center) and initially-outward photons, axis + oblique
|
|
shrinking rr > 0, kernel only (radius would go negative on the far past)
|
|
cancellation translated/far-origin style input cancellation
|
|
|
|
Route cases use only rr == 0 (positive radius on the whole open past segment)
|
|
plus the two production grazing rows.
|
|
"""
|
|
from __future__ import annotations
|
|
|
|
import json
|
|
import math
|
|
from pathlib import Path
|
|
|
|
# Radii spanning subnormal-adjacent to huge scales.
|
|
RADII = [1e-100, 1e-10, 1.0, 1e10, 1e100]
|
|
DIST = [
|
|
math.nextafter(1.0, math.inf),
|
|
2.0,
|
|
10.0,
|
|
100.0,
|
|
512.0,
|
|
1024.0,
|
|
1e4,
|
|
1e8,
|
|
1e10,
|
|
]
|
|
IMPACTS = [
|
|
0.0,
|
|
0.5,
|
|
0.99,
|
|
1.0 - 1e-6,
|
|
math.nextafter(1.0, 0.0),
|
|
1.0,
|
|
math.nextafter(1.0, math.inf),
|
|
1.0 + 1e-6,
|
|
1.1,
|
|
]
|
|
VELS = [0.0, 0.1, 2.0, 10.0]
|
|
RR_GROW = [
|
|
-1.0,
|
|
math.nextafter(-1.0, -math.inf),
|
|
math.nextafter(-1.0, math.inf),
|
|
-0.999999,
|
|
]
|
|
|
|
|
|
def norm3(v):
|
|
n = math.sqrt(v[0] * v[0] + v[1] * v[1] + v[2] * v[2])
|
|
return [v[0] / n, v[1] / n, v[2] / n]
|
|
|
|
|
|
def cross(a, b):
|
|
return [a[1] * b[2] - a[2] * b[1], a[2] * b[0] - a[0] * b[2],
|
|
a[0] * b[1] - a[1] * b[0]]
|
|
|
|
|
|
# Fixed direction frames: (tow, perp). `tow` points from the sphere centre to
|
|
# the camera; `w = -tow` is the inward (toward-centre) past direction.
|
|
DIRS_AXIS = ([1.0, 0.0, 0.0], [0.0, 1.0, 0.0])
|
|
DIRS_OBLIQUE1 = (
|
|
norm3([1.0, 2.0, 3.0]),
|
|
norm3(cross(norm3([1.0, 2.0, 3.0]), [0.0, 0.0, 1.0])),
|
|
)
|
|
DIRS_OBLIQUE2 = (
|
|
norm3([-2.0, 1.0, 0.7]),
|
|
norm3(cross(norm3([-2.0, 1.0, 0.7]), [0.0, 1.0, 0.0])),
|
|
)
|
|
DIRS = [DIRS_AXIS, DIRS_OBLIQUE1, DIRS_OBLIQUE2]
|
|
|
|
|
|
def _case(cid, category, x, c, w, v, R0, rr):
|
|
return {
|
|
"id": cid,
|
|
"category": category,
|
|
"x": list(x),
|
|
"c": list(c),
|
|
"w": list(w),
|
|
"v": list(v),
|
|
"R0": R0,
|
|
"rr": rr,
|
|
}
|
|
|
|
|
|
def curated_kernel_cases():
|
|
"""Explicit adversarial and production-reproduction fixtures."""
|
|
cases = []
|
|
t63 = float.fromhex("0x1.fa8f5c28f5c29p+3")
|
|
t64 = float.fromhex("0x1.fb17e4b17e4b1p+3")
|
|
w63 = [-float.fromhex("0x1.d4afba4704cap-2"),
|
|
float.fromhex("0x1.c7378f8e872d1p-1"),
|
|
float.fromhex("0x1.15bad4e30e8ddp-8")]
|
|
w64 = [-float.fromhex("0x1.f98ae1a782104p-2"),
|
|
float.fromhex("0x1.bd3bb364ac492p-1"),
|
|
float.fromhex("0x1.102d2a1c6ac74p-7")]
|
|
# Production Alcubierre grazing rows: x=(0,-24,0), centre=2t, v=2, R=5.
|
|
cases.append(_case(0, "real63", [0.0, -24.0, 0.0],
|
|
[2.0 * t63, 0.0, 0.0], w63, [2.0, 0.0, 0.0], 5.0, 0.0))
|
|
cases.append(_case(1, "real64", [0.0, -24.0, 0.0],
|
|
[2.0 * t64, 0.0, 0.0], w64, [2.0, 0.0, 0.0], 5.0, 0.0))
|
|
|
|
D = 10.0
|
|
R = 5.0
|
|
cases += [
|
|
_case(2, "headon_hit", [D, 0.0, 0.0], [0.0] * 3, [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, R, 0.0),
|
|
_case(3, "headon_miss", [D, 0.0, 0.0], [0.0] * 3, [1.0, 0.0, 0.0],
|
|
[0.0] * 3, R, 0.0),
|
|
_case(4, "clear_miss", [D, 6.0, 0.0], [0.0] * 3, [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, R, 0.0),
|
|
_case(5, "grazing_in", [D, math.nextafter(R, 0.0), 0.0], [0.0] * 3,
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, R, 0.0),
|
|
_case(6, "exact_tangent", [D, R, 0.0], [0.0] * 3, [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, R, 0.0),
|
|
_case(7, "grazing_out", [D, math.nextafter(R, math.inf), 0.0],
|
|
[0.0] * 3, [-1.0, 0.0, 0.0], [0.0] * 3, R, 0.0),
|
|
_case(8, "near_boundary_hit",
|
|
[math.nextafter(R, math.inf), 0.0, 0.0], [0.0] * 3,
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, R, 0.0),
|
|
# 1e10 + 0.5, R=1: c loses the transverse term, b*b-4ac rounds to 0.
|
|
_case(9, "cancel_1e10_05", [1e10, 0.5, 0.0], [0.0] * 3,
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, 1.0, 0.0),
|
|
_case(10, "cancel_1e10_1", [1e10, 1.0, 0.0], [0.0] * 3,
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, 1.0, 0.0),
|
|
]
|
|
# Large-t positive reconstructed minimum (tests/test_asymptotic.c).
|
|
camera_x = float.fromhex("0x1.6bcc41e901908p+46")
|
|
t0 = 1e15
|
|
vx = 0.1
|
|
cases.append(_case(11, "large_t_01t", [camera_x, 0.0, 0.0],
|
|
[vx * t0, 0.0, 0.0], [-1.0, 0.0, 0.0], [vx, 0.0, 0.0],
|
|
10.0, 0.0))
|
|
# Exact linear growing sphere (a == 0), inward.
|
|
cases.append(_case(12, "linear_grow_in", [100.0, 0.0, 0.0], [0.0] * 3,
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, 10.0, -1.0))
|
|
# Exact linear, initially outward: constant gap, honest miss.
|
|
cases.append(_case(13, "linear_grow_out", [100.0, 0.0, 0.0], [0.0] * 3,
|
|
[1.0, 0.0, 0.0], [0.0] * 3, 10.0, -1.0))
|
|
# Near-linear oblique, inward vs initially outward, rr just around -1.
|
|
tow = DIRS_OBLIQUE1[0]
|
|
perp = DIRS_OBLIQUE1[1]
|
|
for rr in RR_GROW:
|
|
for sign, tag in ((-1.0, "in"), (1.0, "out")):
|
|
w = [sign * tow[i] for i in range(3)]
|
|
x = [50.0 * tow[i] + 0.25 * perp[i] for i in range(3)]
|
|
cases.append(_case(len(cases), f"nearlin_{tag}", x, [0.0] * 3, w,
|
|
[0.0] * 3, 4.0, rr))
|
|
# Translated-origin style cancellation: d = x - c with x = c + small.
|
|
base = 1e10
|
|
cbase = [base, -base, base * 0.5]
|
|
xb = [cbase[0] + 10.0, cbase[1] + 0.5, cbase[2] + 0.0]
|
|
cases.append(_case(len(cases), "translated_origin", xb, cbase,
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, 1.0, 0.0))
|
|
return cases
|
|
|
|
|
|
def _broad_fixed(cid):
|
|
for di, (tow, perp) in enumerate(DIRS):
|
|
for R in RADII:
|
|
for dr in DIST:
|
|
for ir in IMPACTS:
|
|
D = dr * R
|
|
b = ir * R
|
|
x = [D * tow[i] + b * perp[i] for i in range(3)]
|
|
w = [-tow[i] for i in range(3)]
|
|
yield _case(cid, f"fixed_d{di}", x, [0.0] * 3, w, [0.0] * 3,
|
|
R, 0.0)
|
|
cid += 1
|
|
return cid
|
|
|
|
|
|
def _broad_moving(cid):
|
|
for (tow, perp) in DIRS:
|
|
for R in (1e-10, 1.0, 1e10):
|
|
for dr in (2.0, 100.0, 1e4):
|
|
for ir in (0.0, 0.99, 1.0, 1.1):
|
|
for vel in VELS:
|
|
D = dr * R
|
|
b = ir * R
|
|
x = [D * tow[i] + b * perp[i] for i in range(3)]
|
|
w = [-tow[i] for i in range(3)]
|
|
v = [vel * perp[i] for i in range(3)]
|
|
yield _case(cid, "moving", x, [0.0] * 3, w, v, R, 0.0)
|
|
cid += 1
|
|
return cid
|
|
|
|
|
|
def _broad_growing(cid):
|
|
for di, (tow, perp) in enumerate(DIRS):
|
|
for sign, tag in ((-1.0, "in"), (1.0, "out")):
|
|
for rr in RR_GROW:
|
|
for dr in (2.0, 10.0, 100.0, 1e4, 1e8):
|
|
for ir in (0.0, 0.5, 0.99, 1.0, 1.1):
|
|
D = dr * 4.0
|
|
b = ir * 4.0
|
|
x = [D * tow[i] + b * perp[i] for i in range(3)]
|
|
w = [sign * tow[i] for i in range(3)]
|
|
yield _case(cid, f"growing_{tag}", x, [0.0] * 3, w,
|
|
[0.0] * 3, 4.0, rr)
|
|
cid += 1
|
|
return cid
|
|
|
|
|
|
def _broad_shrinking(cid):
|
|
for (tow, perp) in (DIRS_AXIS, DIRS_OBLIQUE1):
|
|
for dr in (2.0, 10.0, 100.0):
|
|
for ir in (0.0, 0.99, 1.0):
|
|
D = dr * 10.0
|
|
b = ir * 10.0
|
|
x = [D * tow[i] + b * perp[i] for i in range(3)]
|
|
w = [-tow[i] for i in range(3)]
|
|
yield _case(cid, "shrinking", x, [0.0] * 3, w, [0.0] * 3, 10.0,
|
|
0.1)
|
|
cid += 1
|
|
return cid
|
|
|
|
|
|
def _materialize(gen, cases):
|
|
for c in gen:
|
|
cases.append(c)
|
|
|
|
|
|
def _build_all():
|
|
cases = curated_kernel_cases()
|
|
cid = 1000
|
|
for gen in (_broad_fixed, _broad_moving, _broad_growing, _broad_shrinking):
|
|
gen_cases = []
|
|
_materialize(gen(cid), gen_cases)
|
|
if gen_cases:
|
|
cid = gen_cases[-1]["id"] + 1
|
|
cases.extend(gen_cases)
|
|
return cases
|
|
|
|
|
|
# ---------------------------------------------------------------------------
|
|
def _route(cid, category, t0, obs, direction, c0, v, R0, rr, model=0,
|
|
valid_t_min=-1e300):
|
|
return {
|
|
"id": cid,
|
|
"category": category,
|
|
"t0": t0,
|
|
"obs": list(obs),
|
|
"dir": list(direction),
|
|
"c0": list(c0),
|
|
"v": list(v),
|
|
"R0": R0,
|
|
"rr": rr,
|
|
"valid_t_min": valid_t_min,
|
|
"model": model,
|
|
}
|
|
|
|
|
|
def route_cases():
|
|
cases = []
|
|
t63 = float.fromhex("0x1.fa8f5c28f5c29p+3")
|
|
t64 = float.fromhex("0x1.fb17e4b17e4b1p+3")
|
|
w63 = [-float.fromhex("0x1.d4afba4704cap-2"),
|
|
float.fromhex("0x1.c7378f8e872d1p-1"),
|
|
float.fromhex("0x1.15bad4e30e8ddp-8")]
|
|
w64 = [-float.fromhex("0x1.f98ae1a782104p-2"),
|
|
float.fromhex("0x1.bd3bb364ac492p-1"),
|
|
float.fromhex("0x1.102d2a1c6ac74p-7")]
|
|
# Production rows reproduced through the Alcubierre-style callback (model 1).
|
|
cases.append(_route(0, "real63", t63, [0.0, -24.0, 0.0], w63,
|
|
[2.0 * t63, 0.0, 0.0], [2.0, 0.0, 0.0], 5.0, 0.0,
|
|
model=1))
|
|
cases.append(_route(1, "real64", t64, [0.0, -24.0, 0.0], w64,
|
|
[2.0 * t64, 0.0, 0.0], [2.0, 0.0, 0.0], 5.0, 0.0,
|
|
model=1))
|
|
# Same rows through the input-stable callback (model 0).
|
|
cases.append(_route(2, "real63_stable", t63, [0.0, -24.0, 0.0], w63,
|
|
[2.0 * t63, 0.0, 0.0], [2.0, 0.0, 0.0], 5.0, 0.0,
|
|
model=0))
|
|
cases.append(_route(3, "real64_stable", t64, [0.0, -24.0, 0.0], w64,
|
|
[2.0 * t64, 0.0, 0.0], [2.0, 0.0, 0.0], 5.0, 0.0,
|
|
model=0))
|
|
R = 5.0
|
|
cases += [
|
|
_route(10, "headon_hit", 0.0, [10.0, 0.0, 0.0], [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(11, "headon_miss", 0.0, [10.0, 0.0, 0.0], [1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(12, "clear_miss", 0.0, [10.0, 6.0, 0.0], [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(13, "grazing_in", 0.0, [10.0, math.nextafter(R, 0.0), 0.0],
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(14, "exact_tangent", 0.0, [10.0, R, 0.0], [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(15, "grazing_out", 0.0,
|
|
[10.0, math.nextafter(R, math.inf), 0.0], [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(16, "camera_inside", 0.0, [2.0, 0.0, 0.0], [1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
_route(17, "near_boundary_hit", 0.0,
|
|
[math.nextafter(R, math.inf), 0.0, 0.0], [-1.0, 0.0, 0.0],
|
|
[0.0] * 3, [0.0] * 3, R, 0.0),
|
|
]
|
|
# Broad fixed-sphere subset through the public route.
|
|
cid = 100
|
|
for (tow, perp) in (DIRS_AXIS, DIRS_OBLIQUE1):
|
|
for R0 in (1e-10, 1.0, 1e10):
|
|
for dr in (2.0, 10.0, 100.0, 512.0, 1024.0, 1e4):
|
|
for ir in (0.0, 0.5, 0.99, 1.0, 1.1):
|
|
D = dr * R0
|
|
b = ir * R0
|
|
obs = [D * tow[i] + b * perp[i] for i in range(3)]
|
|
direction = [-tow[i] for i in range(3)]
|
|
cases.append(_route(cid, f"route_fixed_d{dr:g}", 0.0, obs,
|
|
direction, [0.0] * 3, [0.0] * 3, R0,
|
|
0.0))
|
|
cid += 1
|
|
# Moving spheres (v = 2 tow), model 0.
|
|
for ir in (0.0, 0.99, 1.1):
|
|
cases.append(_route(cid, "route_moving", 0.0, [100.0, 0.0, 0.0],
|
|
[-1.0, 0.0, 0.0], [0.0] * 3, [2.0, 0.0, 0.0], 5.0,
|
|
0.0))
|
|
cid += 1
|
|
# Growing worldtubes (rr <= 0, positive radius on the whole past), inward
|
|
# and initially-outward photons, axis + oblique. These mirror the
|
|
# growing_out kernel false-MISS family and exercise the public route where
|
|
# a kernel MISS bypasses the fallback entirely. valid_t_min is far below
|
|
# any sampled time so no artificial history clip is introduced.
|
|
R0 = 4.0
|
|
for (tow, perp) in DIRS:
|
|
for sign, tag in ((-1.0, "in"), (1.0, "out")):
|
|
for rr in RR_GROW:
|
|
for dr in (2.0, 10.0, 100.0, 512.0, 1024.0, 1e4):
|
|
for ir in (0.0, 1.0, 1.1):
|
|
D = dr * R0
|
|
b = ir * R0
|
|
obs = [D * tow[i] + b * perp[i] for i in range(3)]
|
|
direction = [sign * tow[i] for i in range(3)]
|
|
cases.append(_route(cid, f"route_grow_{tag}_d{dr:g}",
|
|
0.0, obs, direction, [0.0] * 3,
|
|
[0.0] * 3, R0, rr))
|
|
cid += 1
|
|
return cases
|
|
|
|
|
|
# ---------------------------------------------------------------------------
|
|
def _hex(x):
|
|
return float(x).hex()
|
|
|
|
|
|
def _fmt(v):
|
|
return _hex(v)
|
|
|
|
|
|
def write_header(kernel, route, out_path: Path):
|
|
lines = []
|
|
lines.append("/* Generated by benchmarks/quadratic_precision/cases.py. */\n")
|
|
lines.append("#ifndef QUADRATIC_CASES_H\n#define QUADRATIC_CASES_H\n")
|
|
lines.append("typedef struct {\n int id;\n const char *category;\n"
|
|
" double x[3], c[3], w[3], v[3];\n double R0, rr;\n"
|
|
"} QuadKernelCase;\n\n")
|
|
lines.append("typedef struct {\n int id;\n const char *category;\n"
|
|
" double t0;\n double obs[3];\n double dir[3];\n"
|
|
" double c0[3];\n double v[3];\n double R0, rr;\n"
|
|
" double valid_t_min;\n int model;\n} QuadRouteCase;\n\n")
|
|
lines.append("static const QuadKernelCase quad_kernel_cases[] = {\n")
|
|
for c in kernel:
|
|
lines.append(
|
|
" {.id=%d,.category=\"%s\","
|
|
".x={%s,%s,%s},.c={%s,%s,%s},.w={%s,%s,%s},.v={%s,%s,%s},"
|
|
".R0=%s,.rr=%s},\n"
|
|
% (c["id"], c["category"], *[_fmt(z) for z in c["x"]],
|
|
*[_fmt(z) for z in c["c"]], *[_fmt(z) for z in c["w"]],
|
|
*[_fmt(z) for z in c["v"]], _fmt(c["R0"]), _fmt(c["rr"])))
|
|
lines.append("};\n")
|
|
lines.append("static const int quad_kernel_case_count = %d;\n\n"
|
|
% len(kernel))
|
|
lines.append("static const QuadRouteCase quad_route_cases[] = {\n")
|
|
for c in route:
|
|
lines.append(
|
|
" {.id=%d,.category=\"%s\",.t0=%s,"
|
|
".obs={%s,%s,%s},.dir={%s,%s,%s},.c0={%s,%s,%s},.v={%s,%s,%s},"
|
|
".R0=%s,.rr=%s,.valid_t_min=%s,.model=%d},\n"
|
|
% (c["id"], c["category"], _fmt(c["t0"]),
|
|
*[_fmt(z) for z in c["obs"]], *[_fmt(z) for z in c["dir"]],
|
|
*[_fmt(z) for z in c["c0"]], *[_fmt(z) for z in c["v"]],
|
|
_fmt(c["R0"]), _fmt(c["rr"]), _fmt(c["valid_t_min"]),
|
|
c["model"]))
|
|
lines.append("};\n")
|
|
lines.append("static const int quad_route_case_count = %d;\n"
|
|
% len(route))
|
|
lines.append("#endif\n")
|
|
out_path.write_text("".join(lines))
|
|
|
|
|
|
def write_json(kernel, route, out_path: Path):
|
|
def enc(c):
|
|
d = {}
|
|
for k, v in c.items():
|
|
if isinstance(v, float):
|
|
d[k] = v.hex()
|
|
elif isinstance(v, list):
|
|
d[k] = [z.hex() for z in v]
|
|
else:
|
|
d[k] = v
|
|
return d
|
|
|
|
out_path.write_text(
|
|
json.dumps(
|
|
{"deterministic": True, "kernel": [enc(c) for c in kernel],
|
|
"route": [enc(c) for c in route]},
|
|
indent=1,
|
|
)
|
|
+ "\n")
|
|
|
|
|
|
def build():
|
|
kernel = _build_all()
|
|
route = route_cases()
|
|
return kernel, route
|
|
|
|
|
|
if __name__ == "__main__":
|
|
k, r = build()
|
|
print(f"kernel cases: {len(k)}, route cases: {len(r)}")
|
|
cats = {}
|
|
for c in k:
|
|
cats[c["category"]] = cats.get(c["category"], 0) + 1
|
|
for key in sorted(cats):
|
|
print(f" {key}: {cats[key]}")
|
|
rcats = {}
|
|
for c in r:
|
|
rcats[c["category"]] = rcats.get(c["category"], 0) + 1
|
|
print("route categories:")
|
|
for key in sorted(rcats):
|
|
print(f" {key}: {rcats[key]}")
|