From 5143a672113419fb2c57ee46b4917c7ebbe54877 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 13:30:39 +0200 Subject: [PATCH 01/12] 21b tie classification: all exact-path ties are four-node lattice-cocircular under QW2's conditions; QW2 admits 99.97 % of all incircle calls Measurement only; no production code. Instrumented d3ee2ce in a scratch worktree (scripts/instrument.patch, kept unapplied), 1 m benchmark, quarter circle and tile, tolerance 1, AC. Refine loop: 146,962 / 154,502 exact-path calls, every one four node corners, int64 lattice det on (col, -row) = 0, dx == dy, spread <= 2^14, frame exact. Only 38 % are axis-aligned rectangles. Of all refine-loop incircle calls, 99.973 % (quarter) and 100 % (tile) meet QW2's conditions; lattice sign agreed with the kernel on every four-node call. Output identical to 21a's acceptance meshes. Co-Authored-By: Claude Opus 5.5 --- docs/benchmarks/2026-09-27/21b-ties/README.md | 202 +++++++++++ .../2026-09-27/21b-ties/data/pmset-start.txt | 2 + .../2026-09-27/21b-ties/data/tables.md | 60 ++++ .../21b-ties/data/ties_quarter_t1.jsonl | 1 + .../21b-ties/data/ties_quarter_t8.jsonl | 1 + .../21b-ties/data/ties_tile_t1.jsonl | 1 + .../21b-ties/data/ties_tile_t8.jsonl | 1 + .../21b-ties/scripts/classifier_selftest.cpp | 21 ++ .../21b-ties/scripts/counter_selftest.cpp | 16 + .../21b-ties/scripts/instrument.patch | 315 ++++++++++++++++++ .../21b-ties/scripts/prof_driver.py | 51 +++ .../2026-09-27/21b-ties/scripts/tables.py | 102 ++++++ 12 files changed, 773 insertions(+) create mode 100644 docs/benchmarks/2026-09-27/21b-ties/README.md create mode 100644 docs/benchmarks/2026-09-27/21b-ties/data/pmset-start.txt create mode 100644 docs/benchmarks/2026-09-27/21b-ties/data/tables.md create mode 100644 docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t1.jsonl create mode 100644 docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t8.jsonl create mode 100644 docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t1.jsonl create mode 100644 docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t8.jsonl create mode 100644 docs/benchmarks/2026-09-27/21b-ties/scripts/classifier_selftest.cpp create mode 100644 docs/benchmarks/2026-09-27/21b-ties/scripts/counter_selftest.cpp create mode 100644 docs/benchmarks/2026-09-27/21b-ties/scripts/instrument.patch create mode 100644 docs/benchmarks/2026-09-27/21b-ties/scripts/prof_driver.py create mode 100644 docs/benchmarks/2026-09-27/21b-ties/scripts/tables.py diff --git a/docs/benchmarks/2026-09-27/21b-ties/README.md b/docs/benchmarks/2026-09-27/21b-ties/README.md new file mode 100644 index 00000000..a2cf951e --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/README.md @@ -0,0 +1,202 @@ +# 21b tie classification: which incircle tests could QW2's integer path answer? + +This is the measurement that `docs/increments/21-parallel-refine.md` §3 QW2 asks +for before 21b is designed in detail. It counts cases and does not time +anything. + +**Finding.** In the refine loop, every exact-path incircle test is one of these +ties: a quad with four DEM-node corners, a lattice determinant of exactly 0, and +all of QW2's conditions met. That is 146,962 of 146,962 on the quarter circle +and 154,502 of 154,502 on the tile. The 146,962 is the same count the +2026-09-27 profile gave. The ties are not mostly grid rectangles: on the quarter +circle 38 % are axis-aligned rectangles or squares, and 62 % are other +lattice-cocircular quads. Of **all** incircle tests in the refine loop, QW2's +conditions hold on 99.973 % (quarter circle) and 100 % (tile). + +## Method + +- **Commit** `d3ee2ce` (branch `increment21b-lattice-incircle`, 21a's head). + It was measured in a scratch worktree (`git worktree add --detach d3ee2ce`) + with `scripts/instrument.patch` applied. That is a local patch and was never + committed; it is kept here unapplied. `git apply --check + scripts/instrument.patch` passes on `d3ee2ce`. +- **What the patch records.** `FilteredKernel::incircle_ccw` sets a flag + saying whether it took the exact (`E::incircle_ccw`) path. + `detail::must_flip` passes the four `MeshVertex` corners of each incircle + call to `mesh::prof_ties::record` (new header + `include/terrain/mesh/prof_ties.hpp`), in the order they went to the kernel + (`a, b, c` counter-clockwise in the frame, then `d`), together with the + kernel's answer and that flag. `refine` tags the phase: `p0` is + `legalise_all`, `p1` the quality pass, `p2` the refine loop. At the end of + `refine` it appends one JSON line to `$RASPUTIN_TIES_OUT`. For each call it + records: + 1. the number of corners that are nodes (`MeshVertex::is_node`); + 2. for four nodes, the lattice determinant on `(col, -row)` translated to + `d`, computed in `__int128` so it is exact at any spread, with its sign + compared against the kernel's answer; + 3. QW2's conditions, per call: `dx == dy`; every coordinate difference from + `d` at most 2^14 nodes; the frame exact, meaning `col * dx` and + `row * dy` are exact for all four corners (checked with `fma(a, b, -a*b) + == 0`). Once per refine call it also records the sufficient condition + QW2 states: the significant bits of `dx` plus + `bit_width(max(rows, cols) - 1)` are at most 53; + 4. for four nodes and a zero determinant, the shape. A quad is an + axis-aligned square or rectangle when it has two distinct rows and two + distinct columns, and the size is recorded. It is a rotated square or + rectangle when its diagonals have equal midpoints and equal lengths. It + is an isosceles trapezoid when two disjoint chords are parallel: along + rows or columns, or in another direction. Anything else is "other + cyclic". A quad with three collinear corners is counted separately. +- **Probes that can fail.** `scripts/classifier_selftest.cpp` checks the shape + classifier on hand-made quads, one per class, including a circle of radius + 5. `scripts/counter_selftest.cpp` feeds `record` a case for each refusal or + mismatch counter: an inexact frame at `dx = 0.1`, `dx != dy`, a spread above + 2^14, a nonzero lattice determinant on the exact path, a lattice sign that + disagrees with the kernel, and three node corners. Every one of those + counters fired. To build either one: + `c++ -std=c++20 -Wall -Wextra -Wpedantic -Werror -I/include `. +- **The instrumented build produces the same output.** The quarter circle's + binary VTK hashes to `ff705683…`, the value the serial profile got. Hashed + from `POINTS` onward, the way `bench.py` does it, the ASCII meshes match the + 21a acceptance run (`docs/benchmarks/2026-09-27/21a/run.json`): quarter + circle `1e531976…`, tile `11741a81…`. Rounds, insertions and flips match as + well. +- **Build**: Release, configured the way `bench.py`'s `build()` does it + (`-DCMAKE_BUILD_TYPE=Release -DRASPUTIN_BUILD_PYTHON=ON + -DRASPUTIN_BUILD_TESTS=OFF`), in `/build-instr`. The `_core` target was + built (exit 0) and copied into a fresh `build-instr/pkg/tin_engine/` + (symlinks to `src_python/tin_engine/*` plus the `.so`). The `.so` sha256 is + `731a9609…`. +- **Runs.** The inputs are the 1 m benchmark: DEM + `tests/fixtures/dem_archive/7908_3_10m_z33.tif` (5051 × 5051 nodes, + `dx = dy = 10`), tolerance 1, all other CLI defaults. The domains are + `docs/benchmarks/2026-09-26/quarter.geojson` and the whole tile (no + `--domain`). The command runs through `scripts/prof_driver.py`, copied + verbatim from `serial-profile/scripts/`: + + ```bash + RASPUTIN_TIES_OUT=/ties__t.jsonl .venv/bin/python scripts/prof_driver.py \ + --pkg /build-instr/pkg --threads --repeat 1 -- \ + mesh --dem tests/fixtures/dem_archive/7908_3_10m_z33.tif --tolerance 1 \ + [--domain docs/benchmarks/2026-09-26/quarter.geojson] --out /_t.vtk --binary + python scripts/tables.py data/ties_quarter_t1.jsonl data/ties_tile_t1.jsonl > data/tables.md + ``` + + Each domain was run at 1 and 8 threads. The count files are byte-identical + between the two thread counts (`cmp`), as they should be, since the split + phase is serial. +- **Machine**: Apple M1 Max, macOS 27.0. **Power: AC**, 100 % + (`data/pmset-start.txt`). The counts do not depend on power; it is recorded + because the persona requires it. +- **Timings** printed by the driver are inflated by the instrumentation, which + does a `std::map` update per call. None is reported. + +## Results + +`data/tables.md` is produced by `scripts/tables.py` from `data/ties_*_t1.jsonl`. +The rows for the refine loop: + +| measure | quarter | tile | +|---|---:|---:| +| incircle calls | 1,615,895 | 1,639,699 | +| exact path | 146,962 (9.095 %) | 154,502 (9.423 %) | +| exact path, answer Cocircular | 146,962 | 154,502 | +| exact path, four node corners | 146,962 | 154,502 | +| exact path, lattice det = 0 | 146,962 | 154,502 | +| exact path, lattice det != 0 | 0 | 0 | +| exact path, all QW2 conditions hold | 146,962 | 154,502 | +| filtered path, four node corners | 1,468,493 | 1,485,197 | +| filtered path, all QW2 conditions hold | 1,468,493 | 1,485,197 | +| **all calls QW2 would answer** | **1,615,455 (99.973 %)** | **1,639,699 (100 %)** | +| calls with fewer than four node corners | 440 | 0 | +| lattice sign differs from the kernel's answer (four nodes) | 0 | 0 | +| max coordinate difference from `d`, four-node calls | 951 nodes | 88 nodes | +| max coordinate difference from `d`, exact-path calls | 146 nodes | 48 nodes | +| frame-exact sufficient condition | 3 + 13 = 16 ≤ 53 | 3 + 13 = 16 ≤ 53 | + +The per-call conditions never refused a four-node quad. Every four-node call +had `dx == dy`, a spread at most 2^14 and an exact frame, on both domains and +in every phase. + +The shapes of the exact-path quads in the refine loop: + +| shape | quarter | tile | +|---|---:|---:| +| other cyclic quad (no parallel sides) | 47,547 (32.4 %) | 47,269 (30.6 %) | +| axis-aligned square | 39,977 (27.2 %) | 39,911 (25.8 %) | +| rotated square | 17,100 (11.6 %) | 17,203 (11.1 %) | +| axis-aligned rectangle, not square | 16,261 (11.1 %) | 16,567 (10.7 %) | +| isosceles trapezoid, parallel sides on rows or columns | 12,697 (8.6 %) | 20,146 (13.0 %) | +| isosceles trapezoid, other direction | 11,169 (7.6 %) | 11,080 (7.2 %) | +| rotated rectangle, not square | 2,211 (1.5 %) | 2,326 (1.5 %) | +| three corners collinear | 0 | 0 | + +The most common axis-aligned sizes are 1×1 (37,887 on the quarter circle), 1×2 +(13,119), 2×2 (1,845), 1×3 and 2×3. `data/tables.md` has the per-size lists and +the other phases, which also go through `must_flip`: + +- On the quarter circle, `legalise_all` makes 533 calls, of which QW2 would + answer 122. The quality pass makes 7,677, of which QW2 would answer 4,622. + The remaining 411 and 3,055 calls have an off-node corner (the arc + vertices). +- On the tile, QW2 would answer all of `legalise_all`'s 48,133 calls, 16,129 + of them exact ties: 15,876 are 40×40 squares and 252 are 10×40 rectangles, + from the starting mesh. It would also answer all 4,017 of the quality pass's + calls. + +## What is measured and what is inferred + +Measured: + +- QW2's premise holds. All 146,962 (quarter circle) and 154,502 (tile) + exact-path ties in the refine loop have four node corners and a zero int64 + lattice determinant on `(col, -row)`, and each of QW2's per-call conditions + holds for each of them. +- On every four-node call, in every phase and on both domains (3.3 million + calls), the lattice sign agreed with `FilteredKernel`'s answer + on the frame doubles. This is QW2's bit-identity argument checked on this + input. It does not prove the argument. +- The share of all refine-loop incircle calls that QW2's conditions would + admit: 99.973 % and 100 %. +- The ties are lattice-cocircular quads of many shapes. Only 38 % (quarter + circle) and 37 % (tile) are axis-aligned rectangles. + +Not measured, and so inferred or open: + +- **The time QW2 saves.** Answering 99.97 % of calls says nothing about the + cost of an int64 determinant compared with the filtered double determinant + plus its permanent. Whether QW2 reaches past the 5.0 % exact-path cost into + the rest of the 13.2 % flip-test cost has to be measured on an + implementation. +- **Why the ties occur.** One conjecture is that small cocircular quads are + common because refinement inserts nodes near one another on a grid. That + has not been checked. +- **Other DEMs.** This covers one DEM, with `dx = dy = 10`. A DEM with + `dx != dy`, or with a `dx` like 0.1 whose frame is not exact, would get + nothing from QW2 by design. How common such DEMs are was not surveyed. + +## What bears on QW2's design + +1. **The general determinant is needed.** An axis-aligned-rectangle shortcut + would handle only 38 % of the ties, so QW2's general int64 determinant is + the right shape. The "points on a circle of radius 5" case in the planned + suite matches what the data shows. +2. **The spread bound never binds here.** The largest spread was 951 nodes + against a limit of 2^14. By arithmetic, not measurement: for a spread of at + most 2^12, every intermediate of the plain double determinant on the + integer `(col, -row)` differences is an integer below 2^53, so it would be + exact in doubles too. The choice between int64 and doubles can be decided + by timing, not by correctness. +3. **The conditions are per refine, not per call, on this data.** The per-call + checks for `dx == dy` and an exact frame never refused anything. The + once-per-refine decision QW2 already plans is enough, leaving the node check + and the spread check as the only per-call tests. +4. **`must_flip` is also used outside the refine loop.** `legalise_all` + (16,129 ties on the tile) and the quality pass also go through it, so + QW2's call at the top of `must_flip` covers them without extra work. + +## Not in the repository + +The meshes (15-17 MB binary, 22-23 MB ASCII) and the scratch worktree with its build +stayed in the session scratchpad and are gone with it. To regenerate them, +apply `scripts/instrument.patch` to `d3ee2ce` and follow the Method above. diff --git a/docs/benchmarks/2026-09-27/21b-ties/data/pmset-start.txt b/docs/benchmarks/2026-09-27/21b-ties/data/pmset-start.txt new file mode 100644 index 00000000..62504dbf --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/data/pmset-start.txt @@ -0,0 +1,2 @@ +Now drawing from 'AC Power' + -InternalBattery-0 (id=7929955) 100%; finishing charge; 0:13 remaining present: true diff --git a/docs/benchmarks/2026-09-27/21b-ties/data/tables.md b/docs/benchmarks/2026-09-27/21b-ties/data/tables.md new file mode 100644 index 00000000..8415df32 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/data/tables.md @@ -0,0 +1,60 @@ +| measure | quarter | tile | +|---|---:|---:| +| refine loop: incircle calls | 1,615,895 | 1,639,699 | +| refine loop: exact path | 146,962 (9.095 %) | 154,502 (9.423 %) | +| refine loop: exact path, Cocircular | 146,962 | 154,502 | +| refine loop: exact path, four node corners | 146,962 | 154,502 | +| refine loop: exact path, lattice det = 0 | 146,962 | 154,502 | +| refine loop: exact path, lattice det != 0 | 0 | 0 | +| refine loop: exact path, QW2 conditions hold | 146,962 | 154,502 | +| refine loop: filtered path, four node corners | 1,468,493 | 1,485,197 | +| refine loop: filtered path, QW2 conditions hold | 1,468,493 | 1,485,197 | +| refine loop: all calls QW2 would answer | 1,615,455 (99.973 %) | 1,639,699 (100.000 %) | +| refine loop: calls with fewer than four node corners | 440 | 0 | +| refine loop: lattice sign differs from kernel (four nodes) | 0 | 0 | +| legalise_all: incircle calls | 533 | 48,133 | +| legalise_all: exact path | 61 (11.445 %) | 16,129 (33.509 %) | +| legalise_all: exact path, Cocircular | 61 | 16,129 | +| legalise_all: exact path, four node corners | 61 | 16,129 | +| legalise_all: exact path, lattice det = 0 | 61 | 16,129 | +| legalise_all: exact path, lattice det != 0 | 0 | 0 | +| legalise_all: exact path, QW2 conditions hold | 61 | 16,129 | +| legalise_all: filtered path, four node corners | 61 | 32,004 | +| legalise_all: filtered path, QW2 conditions hold | 61 | 32,004 | +| legalise_all: all calls QW2 would answer | 122 (22.889 %) | 48,133 (100.000 %) | +| legalise_all: calls with fewer than four node corners | 411 | 0 | +| legalise_all: lattice sign differs from kernel (four nodes) | 0 | 0 | +| quality pass: incircle calls | 7,677 | 4,017 | +| quality pass: exact path | 22 (0.287 %) | 0 (0.000 %) | +| quality pass: exact path, Cocircular | 22 | 0 | +| quality pass: exact path, four node corners | 22 | 0 | +| quality pass: exact path, lattice det = 0 | 22 | 0 | +| quality pass: exact path, lattice det != 0 | 0 | 0 | +| quality pass: exact path, QW2 conditions hold | 22 | 0 | +| quality pass: filtered path, four node corners | 4,600 | 4,017 | +| quality pass: filtered path, QW2 conditions hold | 4,600 | 4,017 | +| quality pass: all calls QW2 would answer | 4,622 (60.206 %) | 4,017 (100.000 %) | +| quality pass: calls with fewer than four node corners | 3,055 | 0 | +| quality pass: lattice sign differs from kernel (four nodes) | 0 | 0 | +| max spread from d, four-node calls, refine loop (nodes) | 951 | 88 | +| max spread from d, exact-path calls, refine loop (nodes) | 146 | 48 | +| dx, dy; rows x cols | 10, 10; 5051 x 5051 | 10, 10; 5051 x 5051 | +| frame-exact sufficient condition (bits of dx + bit_width) <= 53 | 3 + 13 = 16: True | 3 + 13 = 16: True | +| rounds, inserted, flips | 41, 213,464, 445,657 | 53, 219,837, 445,675 | + +Shapes of the exact-path (tie) quads, refine loop: + +| shape | quarter | tile | +|---|---:|---:| +| other cyclic quad (no parallel sides) | 47,547 (32.353 %) | 47,269 (30.594 %) | +| axis-aligned square | 39,977 (27.202 %) | 39,911 (25.832 %) | +| rotated square | 17,100 (11.636 %) | 17,203 (11.134 %) | +| axis-aligned rectangle, not square | 16,261 (11.065 %) | 16,567 (10.723 %) | +| isosceles trapezoid, parallel sides on rows or columns | 12,697 (8.640 %) | 20,146 (13.039 %) | +| isosceles trapezoid, other direction | 11,169 (7.600 %) | 11,080 (7.171 %) | +| rotated rectangle, not square | 2,211 (1.504 %) | 2,326 (1.505 %) | + +Most frequent axis-aligned sizes (short x long side, nodes), refine loop exact path: + +- quarter: square_1x1 37,887, rectangle_1x2 13,119, square_2x2 1,845, rectangle_1x3 1,182, rectangle_2x3 1,107, rectangle_1x4 240 +- tile: square_1x1 37,866, rectangle_1x2 13,103, square_2x2 1,794, rectangle_2x3 1,133, rectangle_1x3 1,131, rectangle_1x4 207 diff --git a/docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t1.jsonl b/docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t1.jsonl new file mode 100644 index 00000000..939baf10 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t1.jsonl @@ -0,0 +1 @@ +{"dx": 10, "dy": 10, "rows": 5051, "cols": 5051, "dx_significant_bits": 3, "index_bit_width": 13, "global_frame_exact_condition": true, "rounds": 41, "inserted": 213464, "flips": 445657, "max_spread_all4_refine": 951, "max_spread_exact4_refine": 146, "counts": {"p0.calls": 533, "p0.exact.calls": 61, "p0.exact.n4.dx_eq_dy": 61, "p0.exact.n4.frame_exact": 61, "p0.exact.n4.lattice_det_zero": 61, "p0.exact.n4.lattice_sign_agrees": 61, "p0.exact.n4.qw2_qualifies": 61, "p0.exact.n4.shape.isosceles_trapezoid_other": 61, "p0.exact.n4.spread_le_2^14": 61, "p0.exact.nodes4": 61, "p0.exact.nodes4.result.cocircular": 61, "p0.exact.result.cocircular": 61, "p0.filtered.calls": 472, "p0.filtered.n4.dx_eq_dy": 61, "p0.filtered.n4.frame_exact": 61, "p0.filtered.n4.lattice_det_neg": 61, "p0.filtered.n4.lattice_sign_agrees": 61, "p0.filtered.n4.qw2_qualifies": 61, "p0.filtered.n4.spread_le_2^14": 61, "p0.filtered.nodes1": 58, "p0.filtered.nodes1.result.outside": 58, "p0.filtered.nodes2": 352, "p0.filtered.nodes2.result.outside": 352, "p0.filtered.nodes3": 1, "p0.filtered.nodes3.result.outside": 1, "p0.filtered.nodes4": 61, "p0.filtered.nodes4.result.outside": 61, "p0.filtered.result.outside": 472, "p1.calls": 7677, "p1.exact.calls": 22, "p1.exact.n4.dx_eq_dy": 22, "p1.exact.n4.frame_exact": 22, "p1.exact.n4.lattice_det_zero": 22, "p1.exact.n4.lattice_sign_agrees": 22, "p1.exact.n4.qw2_qualifies": 22, "p1.exact.n4.shape.isosceles_trapezoid_axis": 22, "p1.exact.n4.spread_le_2^14": 22, "p1.exact.nodes4": 22, "p1.exact.nodes4.result.cocircular": 22, "p1.exact.result.cocircular": 22, "p1.filtered.calls": 7655, "p1.filtered.n4.dx_eq_dy": 4600, "p1.filtered.n4.frame_exact": 4600, "p1.filtered.n4.lattice_det_neg": 1850, "p1.filtered.n4.lattice_det_pos": 2750, "p1.filtered.n4.lattice_sign_agrees": 4600, "p1.filtered.n4.qw2_qualifies": 4600, "p1.filtered.n4.spread_le_2^14": 4600, "p1.filtered.nodes2": 2637, "p1.filtered.nodes2.result.inside": 2093, "p1.filtered.nodes2.result.outside": 544, "p1.filtered.nodes3": 418, "p1.filtered.nodes3.result.inside": 368, "p1.filtered.nodes3.result.outside": 50, "p1.filtered.nodes4": 4600, "p1.filtered.nodes4.result.inside": 2750, "p1.filtered.nodes4.result.outside": 1850, "p1.filtered.result.inside": 5211, "p1.filtered.result.outside": 2444, "p2.calls": 1615895, "p2.exact.calls": 146962, "p2.exact.n4.dx_eq_dy": 146962, "p2.exact.n4.frame_exact": 146962, "p2.exact.n4.lattice_det_zero": 146962, "p2.exact.n4.lattice_sign_agrees": 146962, "p2.exact.n4.qw2_qualifies": 146962, "p2.exact.n4.shape.axis_rectangle_1x10": 1, "p2.exact.n4.shape.axis_rectangle_1x13": 1, "p2.exact.n4.shape.axis_rectangle_1x17": 1, "p2.exact.n4.shape.axis_rectangle_1x2": 13119, "p2.exact.n4.shape.axis_rectangle_1x3": 1182, "p2.exact.n4.shape.axis_rectangle_1x4": 240, "p2.exact.n4.shape.axis_rectangle_1x5": 61, "p2.exact.n4.shape.axis_rectangle_1x6": 21, "p2.exact.n4.shape.axis_rectangle_1x7": 5, "p2.exact.n4.shape.axis_rectangle_1x8": 5, "p2.exact.n4.shape.axis_rectangle_1x9": 2, "p2.exact.n4.shape.axis_rectangle_2x10": 1, "p2.exact.n4.shape.axis_rectangle_2x3": 1107, "p2.exact.n4.shape.axis_rectangle_2x4": 179, "p2.exact.n4.shape.axis_rectangle_2x5": 64, "p2.exact.n4.shape.axis_rectangle_2x6": 17, "p2.exact.n4.shape.axis_rectangle_2x7": 11, "p2.exact.n4.shape.axis_rectangle_2x8": 5, "p2.exact.n4.shape.axis_rectangle_2x9": 2, "p2.exact.n4.shape.axis_rectangle_3x10": 1, "p2.exact.n4.shape.axis_rectangle_3x4": 108, "p2.exact.n4.shape.axis_rectangle_3x5": 36, "p2.exact.n4.shape.axis_rectangle_3x6": 16, "p2.exact.n4.shape.axis_rectangle_3x7": 7, "p2.exact.n4.shape.axis_rectangle_3x8": 1, "p2.exact.n4.shape.axis_rectangle_3x9": 1, "p2.exact.n4.shape.axis_rectangle_4x5": 27, "p2.exact.n4.shape.axis_rectangle_4x6": 12, "p2.exact.n4.shape.axis_rectangle_4x7": 3, "p2.exact.n4.shape.axis_rectangle_4x8": 2, "p2.exact.n4.shape.axis_rectangle_4x9": 1, "p2.exact.n4.shape.axis_rectangle_5x10": 1, "p2.exact.n4.shape.axis_rectangle_5x6": 8, "p2.exact.n4.shape.axis_rectangle_5x7": 1, "p2.exact.n4.shape.axis_rectangle_5x8": 2, "p2.exact.n4.shape.axis_rectangle_5x9": 1, "p2.exact.n4.shape.axis_rectangle_6x10": 1, "p2.exact.n4.shape.axis_rectangle_6x12": 1, "p2.exact.n4.shape.axis_rectangle_6x8": 1, "p2.exact.n4.shape.axis_rectangle_6x9": 2, "p2.exact.n4.shape.axis_rectangle_7x8": 2, "p2.exact.n4.shape.axis_rectangle_8x11": 1, "p2.exact.n4.shape.axis_rectangle_9x16": 1, "p2.exact.n4.shape.axis_square_1x1": 37887, "p2.exact.n4.shape.axis_square_2x2": 1845, "p2.exact.n4.shape.axis_square_3x3": 193, "p2.exact.n4.shape.axis_square_4x4": 43, "p2.exact.n4.shape.axis_square_5x5": 6, "p2.exact.n4.shape.axis_square_6x6": 2, "p2.exact.n4.shape.axis_square_7x7": 1, "p2.exact.n4.shape.isosceles_trapezoid_axis": 12697, "p2.exact.n4.shape.isosceles_trapezoid_other": 11169, "p2.exact.n4.shape.other_cyclic": 47547, "p2.exact.n4.shape.rotated_rectangle": 2211, "p2.exact.n4.shape.rotated_square": 17100, "p2.exact.n4.spread_le_2^14": 146962, "p2.exact.nodes4": 146962, "p2.exact.nodes4.result.cocircular": 146962, "p2.exact.result.cocircular": 146962, "p2.filtered.calls": 1468933, "p2.filtered.n4.dx_eq_dy": 1468493, "p2.filtered.n4.frame_exact": 1468493, "p2.filtered.n4.lattice_det_neg": 1022985, "p2.filtered.n4.lattice_det_pos": 445508, "p2.filtered.n4.lattice_sign_agrees": 1468493, "p2.filtered.n4.qw2_qualifies": 1468493, "p2.filtered.n4.spread_le_2^14": 1468493, "p2.filtered.nodes1": 2, "p2.filtered.nodes1.result.outside": 2, "p2.filtered.nodes2": 127, "p2.filtered.nodes2.result.inside": 38, "p2.filtered.nodes2.result.outside": 89, "p2.filtered.nodes3": 311, "p2.filtered.nodes3.result.inside": 111, "p2.filtered.nodes3.result.outside": 200, "p2.filtered.nodes4": 1468493, "p2.filtered.nodes4.result.inside": 445508, "p2.filtered.nodes4.result.outside": 1022985, "p2.filtered.result.inside": 445657, "p2.filtered.result.outside": 1023276}} diff --git a/docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t8.jsonl b/docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t8.jsonl new file mode 100644 index 00000000..939baf10 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/data/ties_quarter_t8.jsonl @@ -0,0 +1 @@ +{"dx": 10, "dy": 10, "rows": 5051, "cols": 5051, "dx_significant_bits": 3, "index_bit_width": 13, "global_frame_exact_condition": true, "rounds": 41, "inserted": 213464, "flips": 445657, "max_spread_all4_refine": 951, "max_spread_exact4_refine": 146, "counts": {"p0.calls": 533, "p0.exact.calls": 61, "p0.exact.n4.dx_eq_dy": 61, "p0.exact.n4.frame_exact": 61, "p0.exact.n4.lattice_det_zero": 61, "p0.exact.n4.lattice_sign_agrees": 61, "p0.exact.n4.qw2_qualifies": 61, "p0.exact.n4.shape.isosceles_trapezoid_other": 61, "p0.exact.n4.spread_le_2^14": 61, "p0.exact.nodes4": 61, "p0.exact.nodes4.result.cocircular": 61, "p0.exact.result.cocircular": 61, "p0.filtered.calls": 472, "p0.filtered.n4.dx_eq_dy": 61, "p0.filtered.n4.frame_exact": 61, "p0.filtered.n4.lattice_det_neg": 61, "p0.filtered.n4.lattice_sign_agrees": 61, "p0.filtered.n4.qw2_qualifies": 61, "p0.filtered.n4.spread_le_2^14": 61, "p0.filtered.nodes1": 58, "p0.filtered.nodes1.result.outside": 58, "p0.filtered.nodes2": 352, "p0.filtered.nodes2.result.outside": 352, "p0.filtered.nodes3": 1, "p0.filtered.nodes3.result.outside": 1, "p0.filtered.nodes4": 61, "p0.filtered.nodes4.result.outside": 61, "p0.filtered.result.outside": 472, "p1.calls": 7677, "p1.exact.calls": 22, "p1.exact.n4.dx_eq_dy": 22, "p1.exact.n4.frame_exact": 22, "p1.exact.n4.lattice_det_zero": 22, "p1.exact.n4.lattice_sign_agrees": 22, "p1.exact.n4.qw2_qualifies": 22, "p1.exact.n4.shape.isosceles_trapezoid_axis": 22, "p1.exact.n4.spread_le_2^14": 22, "p1.exact.nodes4": 22, "p1.exact.nodes4.result.cocircular": 22, "p1.exact.result.cocircular": 22, "p1.filtered.calls": 7655, "p1.filtered.n4.dx_eq_dy": 4600, "p1.filtered.n4.frame_exact": 4600, "p1.filtered.n4.lattice_det_neg": 1850, "p1.filtered.n4.lattice_det_pos": 2750, "p1.filtered.n4.lattice_sign_agrees": 4600, "p1.filtered.n4.qw2_qualifies": 4600, "p1.filtered.n4.spread_le_2^14": 4600, "p1.filtered.nodes2": 2637, "p1.filtered.nodes2.result.inside": 2093, "p1.filtered.nodes2.result.outside": 544, "p1.filtered.nodes3": 418, "p1.filtered.nodes3.result.inside": 368, "p1.filtered.nodes3.result.outside": 50, "p1.filtered.nodes4": 4600, "p1.filtered.nodes4.result.inside": 2750, "p1.filtered.nodes4.result.outside": 1850, "p1.filtered.result.inside": 5211, "p1.filtered.result.outside": 2444, "p2.calls": 1615895, "p2.exact.calls": 146962, "p2.exact.n4.dx_eq_dy": 146962, "p2.exact.n4.frame_exact": 146962, "p2.exact.n4.lattice_det_zero": 146962, "p2.exact.n4.lattice_sign_agrees": 146962, "p2.exact.n4.qw2_qualifies": 146962, "p2.exact.n4.shape.axis_rectangle_1x10": 1, "p2.exact.n4.shape.axis_rectangle_1x13": 1, "p2.exact.n4.shape.axis_rectangle_1x17": 1, "p2.exact.n4.shape.axis_rectangle_1x2": 13119, "p2.exact.n4.shape.axis_rectangle_1x3": 1182, "p2.exact.n4.shape.axis_rectangle_1x4": 240, "p2.exact.n4.shape.axis_rectangle_1x5": 61, "p2.exact.n4.shape.axis_rectangle_1x6": 21, "p2.exact.n4.shape.axis_rectangle_1x7": 5, "p2.exact.n4.shape.axis_rectangle_1x8": 5, "p2.exact.n4.shape.axis_rectangle_1x9": 2, "p2.exact.n4.shape.axis_rectangle_2x10": 1, "p2.exact.n4.shape.axis_rectangle_2x3": 1107, "p2.exact.n4.shape.axis_rectangle_2x4": 179, "p2.exact.n4.shape.axis_rectangle_2x5": 64, "p2.exact.n4.shape.axis_rectangle_2x6": 17, "p2.exact.n4.shape.axis_rectangle_2x7": 11, "p2.exact.n4.shape.axis_rectangle_2x8": 5, "p2.exact.n4.shape.axis_rectangle_2x9": 2, "p2.exact.n4.shape.axis_rectangle_3x10": 1, "p2.exact.n4.shape.axis_rectangle_3x4": 108, "p2.exact.n4.shape.axis_rectangle_3x5": 36, "p2.exact.n4.shape.axis_rectangle_3x6": 16, "p2.exact.n4.shape.axis_rectangle_3x7": 7, "p2.exact.n4.shape.axis_rectangle_3x8": 1, "p2.exact.n4.shape.axis_rectangle_3x9": 1, "p2.exact.n4.shape.axis_rectangle_4x5": 27, "p2.exact.n4.shape.axis_rectangle_4x6": 12, "p2.exact.n4.shape.axis_rectangle_4x7": 3, "p2.exact.n4.shape.axis_rectangle_4x8": 2, "p2.exact.n4.shape.axis_rectangle_4x9": 1, "p2.exact.n4.shape.axis_rectangle_5x10": 1, "p2.exact.n4.shape.axis_rectangle_5x6": 8, "p2.exact.n4.shape.axis_rectangle_5x7": 1, "p2.exact.n4.shape.axis_rectangle_5x8": 2, "p2.exact.n4.shape.axis_rectangle_5x9": 1, "p2.exact.n4.shape.axis_rectangle_6x10": 1, "p2.exact.n4.shape.axis_rectangle_6x12": 1, "p2.exact.n4.shape.axis_rectangle_6x8": 1, "p2.exact.n4.shape.axis_rectangle_6x9": 2, "p2.exact.n4.shape.axis_rectangle_7x8": 2, "p2.exact.n4.shape.axis_rectangle_8x11": 1, "p2.exact.n4.shape.axis_rectangle_9x16": 1, "p2.exact.n4.shape.axis_square_1x1": 37887, "p2.exact.n4.shape.axis_square_2x2": 1845, "p2.exact.n4.shape.axis_square_3x3": 193, "p2.exact.n4.shape.axis_square_4x4": 43, "p2.exact.n4.shape.axis_square_5x5": 6, "p2.exact.n4.shape.axis_square_6x6": 2, "p2.exact.n4.shape.axis_square_7x7": 1, "p2.exact.n4.shape.isosceles_trapezoid_axis": 12697, "p2.exact.n4.shape.isosceles_trapezoid_other": 11169, "p2.exact.n4.shape.other_cyclic": 47547, "p2.exact.n4.shape.rotated_rectangle": 2211, "p2.exact.n4.shape.rotated_square": 17100, "p2.exact.n4.spread_le_2^14": 146962, "p2.exact.nodes4": 146962, "p2.exact.nodes4.result.cocircular": 146962, "p2.exact.result.cocircular": 146962, "p2.filtered.calls": 1468933, "p2.filtered.n4.dx_eq_dy": 1468493, "p2.filtered.n4.frame_exact": 1468493, "p2.filtered.n4.lattice_det_neg": 1022985, "p2.filtered.n4.lattice_det_pos": 445508, "p2.filtered.n4.lattice_sign_agrees": 1468493, "p2.filtered.n4.qw2_qualifies": 1468493, "p2.filtered.n4.spread_le_2^14": 1468493, "p2.filtered.nodes1": 2, "p2.filtered.nodes1.result.outside": 2, "p2.filtered.nodes2": 127, "p2.filtered.nodes2.result.inside": 38, "p2.filtered.nodes2.result.outside": 89, "p2.filtered.nodes3": 311, "p2.filtered.nodes3.result.inside": 111, "p2.filtered.nodes3.result.outside": 200, "p2.filtered.nodes4": 1468493, "p2.filtered.nodes4.result.inside": 445508, "p2.filtered.nodes4.result.outside": 1022985, "p2.filtered.result.inside": 445657, "p2.filtered.result.outside": 1023276}} diff --git a/docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t1.jsonl b/docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t1.jsonl new file mode 100644 index 00000000..5e227269 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t1.jsonl @@ -0,0 +1 @@ +{"dx": 10, "dy": 10, "rows": 5051, "cols": 5051, "dx_significant_bits": 3, "index_bit_width": 13, "global_frame_exact_condition": true, "rounds": 53, "inserted": 219837, "flips": 445675, "max_spread_all4_refine": 88, "max_spread_exact4_refine": 48, "counts": {"p0.calls": 48133, "p0.exact.calls": 16129, "p0.exact.n4.dx_eq_dy": 16129, "p0.exact.n4.frame_exact": 16129, "p0.exact.n4.lattice_det_zero": 16129, "p0.exact.n4.lattice_sign_agrees": 16129, "p0.exact.n4.qw2_qualifies": 16129, "p0.exact.n4.shape.axis_rectangle_10x40": 252, "p0.exact.n4.shape.axis_square_10x10": 1, "p0.exact.n4.shape.axis_square_40x40": 15876, "p0.exact.n4.spread_le_2^14": 16129, "p0.exact.nodes4": 16129, "p0.exact.nodes4.result.cocircular": 16129, "p0.exact.result.cocircular": 16129, "p0.filtered.calls": 32004, "p0.filtered.n4.dx_eq_dy": 32004, "p0.filtered.n4.frame_exact": 32004, "p0.filtered.n4.lattice_det_neg": 32004, "p0.filtered.n4.lattice_sign_agrees": 32004, "p0.filtered.n4.qw2_qualifies": 32004, "p0.filtered.n4.spread_le_2^14": 32004, "p0.filtered.nodes4": 32004, "p0.filtered.nodes4.result.outside": 32004, "p0.filtered.result.outside": 32004, "p1.calls": 4017, "p1.filtered.calls": 4017, "p1.filtered.n4.dx_eq_dy": 4017, "p1.filtered.n4.frame_exact": 4017, "p1.filtered.n4.lattice_det_neg": 2760, "p1.filtered.n4.lattice_det_pos": 1257, "p1.filtered.n4.lattice_sign_agrees": 4017, "p1.filtered.n4.qw2_qualifies": 4017, "p1.filtered.n4.spread_le_2^14": 4017, "p1.filtered.nodes4": 4017, "p1.filtered.nodes4.result.inside": 1257, "p1.filtered.nodes4.result.outside": 2760, "p1.filtered.result.inside": 1257, "p1.filtered.result.outside": 2760, "p2.calls": 1639699, "p2.exact.calls": 154502, "p2.exact.n4.dx_eq_dy": 154502, "p2.exact.n4.frame_exact": 154502, "p2.exact.n4.lattice_det_zero": 154502, "p2.exact.n4.lattice_sign_agrees": 154502, "p2.exact.n4.qw2_qualifies": 154502, "p2.exact.n4.shape.axis_rectangle_1x11": 2, "p2.exact.n4.shape.axis_rectangle_1x12": 1, "p2.exact.n4.shape.axis_rectangle_1x13": 1, "p2.exact.n4.shape.axis_rectangle_1x2": 13103, "p2.exact.n4.shape.axis_rectangle_1x3": 1131, "p2.exact.n4.shape.axis_rectangle_1x33": 1, "p2.exact.n4.shape.axis_rectangle_1x4": 207, "p2.exact.n4.shape.axis_rectangle_1x5": 55, "p2.exact.n4.shape.axis_rectangle_1x6": 28, "p2.exact.n4.shape.axis_rectangle_1x7": 5, "p2.exact.n4.shape.axis_rectangle_1x8": 3, "p2.exact.n4.shape.axis_rectangle_1x9": 1, "p2.exact.n4.shape.axis_rectangle_2x11": 1, "p2.exact.n4.shape.axis_rectangle_2x14": 1, "p2.exact.n4.shape.axis_rectangle_2x3": 1133, "p2.exact.n4.shape.axis_rectangle_2x4": 173, "p2.exact.n4.shape.axis_rectangle_2x5": 70, "p2.exact.n4.shape.axis_rectangle_2x6": 18, "p2.exact.n4.shape.axis_rectangle_2x7": 8, "p2.exact.n4.shape.axis_rectangle_2x8": 7, "p2.exact.n4.shape.axis_rectangle_2x9": 1, "p2.exact.n4.shape.axis_rectangle_34x35": 1, "p2.exact.n4.shape.axis_rectangle_34x40": 70, "p2.exact.n4.shape.axis_rectangle_35x40": 124, "p2.exact.n4.shape.axis_rectangle_3x10": 1, "p2.exact.n4.shape.axis_rectangle_3x11": 1, "p2.exact.n4.shape.axis_rectangle_3x4": 115, "p2.exact.n4.shape.axis_rectangle_3x5": 37, "p2.exact.n4.shape.axis_rectangle_3x6": 16, "p2.exact.n4.shape.axis_rectangle_3x7": 3, "p2.exact.n4.shape.axis_rectangle_3x8": 4, "p2.exact.n4.shape.axis_rectangle_3x9": 1, "p2.exact.n4.shape.axis_rectangle_4x11": 1, "p2.exact.n4.shape.axis_rectangle_4x12": 1, "p2.exact.n4.shape.axis_rectangle_4x20": 1, "p2.exact.n4.shape.axis_rectangle_4x5": 15, "p2.exact.n4.shape.axis_rectangle_4x6": 11, "p2.exact.n4.shape.axis_rectangle_4x7": 3, "p2.exact.n4.shape.axis_rectangle_4x8": 2, "p2.exact.n4.shape.axis_rectangle_4x9": 3, "p2.exact.n4.shape.axis_rectangle_5x10": 1, "p2.exact.n4.shape.axis_rectangle_5x40": 125, "p2.exact.n4.shape.axis_rectangle_5x6": 3, "p2.exact.n4.shape.axis_rectangle_5x7": 1, "p2.exact.n4.shape.axis_rectangle_6x40": 70, "p2.exact.n4.shape.axis_rectangle_6x7": 2, "p2.exact.n4.shape.axis_rectangle_6x8": 2, "p2.exact.n4.shape.axis_rectangle_8x9": 1, "p2.exact.n4.shape.axis_rectangle_9x11": 1, "p2.exact.n4.shape.axis_rectangle_9x13": 1, "p2.exact.n4.shape.axis_square_1x1": 37866, "p2.exact.n4.shape.axis_square_2x2": 1794, "p2.exact.n4.shape.axis_square_3x3": 201, "p2.exact.n4.shape.axis_square_4x4": 37, "p2.exact.n4.shape.axis_square_5x5": 12, "p2.exact.n4.shape.axis_square_7x7": 1, "p2.exact.n4.shape.isosceles_trapezoid_axis": 20146, "p2.exact.n4.shape.isosceles_trapezoid_other": 11080, "p2.exact.n4.shape.other_cyclic": 47269, "p2.exact.n4.shape.rotated_rectangle": 2326, "p2.exact.n4.shape.rotated_square": 17203, "p2.exact.n4.spread_le_2^14": 154502, "p2.exact.nodes4": 154502, "p2.exact.nodes4.result.cocircular": 154502, "p2.exact.result.cocircular": 154502, "p2.filtered.calls": 1485197, "p2.filtered.n4.dx_eq_dy": 1485197, "p2.filtered.n4.frame_exact": 1485197, "p2.filtered.n4.lattice_det_neg": 1039522, "p2.filtered.n4.lattice_det_pos": 445675, "p2.filtered.n4.lattice_sign_agrees": 1485197, "p2.filtered.n4.qw2_qualifies": 1485197, "p2.filtered.n4.spread_le_2^14": 1485197, "p2.filtered.nodes4": 1485197, "p2.filtered.nodes4.result.inside": 445675, "p2.filtered.nodes4.result.outside": 1039522, "p2.filtered.result.inside": 445675, "p2.filtered.result.outside": 1039522}} diff --git a/docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t8.jsonl b/docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t8.jsonl new file mode 100644 index 00000000..5e227269 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/data/ties_tile_t8.jsonl @@ -0,0 +1 @@ +{"dx": 10, "dy": 10, "rows": 5051, "cols": 5051, "dx_significant_bits": 3, "index_bit_width": 13, "global_frame_exact_condition": true, "rounds": 53, "inserted": 219837, "flips": 445675, "max_spread_all4_refine": 88, "max_spread_exact4_refine": 48, "counts": {"p0.calls": 48133, "p0.exact.calls": 16129, "p0.exact.n4.dx_eq_dy": 16129, "p0.exact.n4.frame_exact": 16129, "p0.exact.n4.lattice_det_zero": 16129, "p0.exact.n4.lattice_sign_agrees": 16129, "p0.exact.n4.qw2_qualifies": 16129, "p0.exact.n4.shape.axis_rectangle_10x40": 252, "p0.exact.n4.shape.axis_square_10x10": 1, "p0.exact.n4.shape.axis_square_40x40": 15876, "p0.exact.n4.spread_le_2^14": 16129, "p0.exact.nodes4": 16129, "p0.exact.nodes4.result.cocircular": 16129, "p0.exact.result.cocircular": 16129, "p0.filtered.calls": 32004, "p0.filtered.n4.dx_eq_dy": 32004, "p0.filtered.n4.frame_exact": 32004, "p0.filtered.n4.lattice_det_neg": 32004, "p0.filtered.n4.lattice_sign_agrees": 32004, "p0.filtered.n4.qw2_qualifies": 32004, "p0.filtered.n4.spread_le_2^14": 32004, "p0.filtered.nodes4": 32004, "p0.filtered.nodes4.result.outside": 32004, "p0.filtered.result.outside": 32004, "p1.calls": 4017, "p1.filtered.calls": 4017, "p1.filtered.n4.dx_eq_dy": 4017, "p1.filtered.n4.frame_exact": 4017, "p1.filtered.n4.lattice_det_neg": 2760, "p1.filtered.n4.lattice_det_pos": 1257, "p1.filtered.n4.lattice_sign_agrees": 4017, "p1.filtered.n4.qw2_qualifies": 4017, "p1.filtered.n4.spread_le_2^14": 4017, "p1.filtered.nodes4": 4017, "p1.filtered.nodes4.result.inside": 1257, "p1.filtered.nodes4.result.outside": 2760, "p1.filtered.result.inside": 1257, "p1.filtered.result.outside": 2760, "p2.calls": 1639699, "p2.exact.calls": 154502, "p2.exact.n4.dx_eq_dy": 154502, "p2.exact.n4.frame_exact": 154502, "p2.exact.n4.lattice_det_zero": 154502, "p2.exact.n4.lattice_sign_agrees": 154502, "p2.exact.n4.qw2_qualifies": 154502, "p2.exact.n4.shape.axis_rectangle_1x11": 2, "p2.exact.n4.shape.axis_rectangle_1x12": 1, "p2.exact.n4.shape.axis_rectangle_1x13": 1, "p2.exact.n4.shape.axis_rectangle_1x2": 13103, "p2.exact.n4.shape.axis_rectangle_1x3": 1131, "p2.exact.n4.shape.axis_rectangle_1x33": 1, "p2.exact.n4.shape.axis_rectangle_1x4": 207, "p2.exact.n4.shape.axis_rectangle_1x5": 55, "p2.exact.n4.shape.axis_rectangle_1x6": 28, "p2.exact.n4.shape.axis_rectangle_1x7": 5, "p2.exact.n4.shape.axis_rectangle_1x8": 3, "p2.exact.n4.shape.axis_rectangle_1x9": 1, "p2.exact.n4.shape.axis_rectangle_2x11": 1, "p2.exact.n4.shape.axis_rectangle_2x14": 1, "p2.exact.n4.shape.axis_rectangle_2x3": 1133, "p2.exact.n4.shape.axis_rectangle_2x4": 173, "p2.exact.n4.shape.axis_rectangle_2x5": 70, "p2.exact.n4.shape.axis_rectangle_2x6": 18, "p2.exact.n4.shape.axis_rectangle_2x7": 8, "p2.exact.n4.shape.axis_rectangle_2x8": 7, "p2.exact.n4.shape.axis_rectangle_2x9": 1, "p2.exact.n4.shape.axis_rectangle_34x35": 1, "p2.exact.n4.shape.axis_rectangle_34x40": 70, "p2.exact.n4.shape.axis_rectangle_35x40": 124, "p2.exact.n4.shape.axis_rectangle_3x10": 1, "p2.exact.n4.shape.axis_rectangle_3x11": 1, "p2.exact.n4.shape.axis_rectangle_3x4": 115, "p2.exact.n4.shape.axis_rectangle_3x5": 37, "p2.exact.n4.shape.axis_rectangle_3x6": 16, "p2.exact.n4.shape.axis_rectangle_3x7": 3, "p2.exact.n4.shape.axis_rectangle_3x8": 4, "p2.exact.n4.shape.axis_rectangle_3x9": 1, "p2.exact.n4.shape.axis_rectangle_4x11": 1, "p2.exact.n4.shape.axis_rectangle_4x12": 1, "p2.exact.n4.shape.axis_rectangle_4x20": 1, "p2.exact.n4.shape.axis_rectangle_4x5": 15, "p2.exact.n4.shape.axis_rectangle_4x6": 11, "p2.exact.n4.shape.axis_rectangle_4x7": 3, "p2.exact.n4.shape.axis_rectangle_4x8": 2, "p2.exact.n4.shape.axis_rectangle_4x9": 3, "p2.exact.n4.shape.axis_rectangle_5x10": 1, "p2.exact.n4.shape.axis_rectangle_5x40": 125, "p2.exact.n4.shape.axis_rectangle_5x6": 3, "p2.exact.n4.shape.axis_rectangle_5x7": 1, "p2.exact.n4.shape.axis_rectangle_6x40": 70, "p2.exact.n4.shape.axis_rectangle_6x7": 2, "p2.exact.n4.shape.axis_rectangle_6x8": 2, "p2.exact.n4.shape.axis_rectangle_8x9": 1, "p2.exact.n4.shape.axis_rectangle_9x11": 1, "p2.exact.n4.shape.axis_rectangle_9x13": 1, "p2.exact.n4.shape.axis_square_1x1": 37866, "p2.exact.n4.shape.axis_square_2x2": 1794, "p2.exact.n4.shape.axis_square_3x3": 201, "p2.exact.n4.shape.axis_square_4x4": 37, "p2.exact.n4.shape.axis_square_5x5": 12, "p2.exact.n4.shape.axis_square_7x7": 1, "p2.exact.n4.shape.isosceles_trapezoid_axis": 20146, "p2.exact.n4.shape.isosceles_trapezoid_other": 11080, "p2.exact.n4.shape.other_cyclic": 47269, "p2.exact.n4.shape.rotated_rectangle": 2326, "p2.exact.n4.shape.rotated_square": 17203, "p2.exact.n4.spread_le_2^14": 154502, "p2.exact.nodes4": 154502, "p2.exact.nodes4.result.cocircular": 154502, "p2.exact.result.cocircular": 154502, "p2.filtered.calls": 1485197, "p2.filtered.n4.dx_eq_dy": 1485197, "p2.filtered.n4.frame_exact": 1485197, "p2.filtered.n4.lattice_det_neg": 1039522, "p2.filtered.n4.lattice_det_pos": 445675, "p2.filtered.n4.lattice_sign_agrees": 1485197, "p2.filtered.n4.qw2_qualifies": 1485197, "p2.filtered.n4.spread_le_2^14": 1485197, "p2.filtered.nodes4": 1485197, "p2.filtered.nodes4.result.inside": 445675, "p2.filtered.nodes4.result.outside": 1039522, "p2.filtered.result.inside": 445675, "p2.filtered.result.outside": 1039522}} diff --git a/docs/benchmarks/2026-09-27/21b-ties/scripts/classifier_selftest.cpp b/docs/benchmarks/2026-09-27/21b-ties/scripts/classifier_selftest.cpp new file mode 100644 index 00000000..f471dea8 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/scripts/classifier_selftest.cpp @@ -0,0 +1,21 @@ +#include +#include +using namespace terrain::mesh; +using namespace terrain::mesh::prof_ties; +static MeshVertex V(double c, double r) { return MeshVertex{c, r}; } +int main() { + // shape on (col, row); the classifier reads them as a set + std::printf("%s\n", shape({V(0,0), V(1,0), V(1,1), V(0,1)}).c_str()); // axis_square_1x1 + std::printf("%s\n", shape({V(0,0), V(2,0), V(2,1), V(0,1)}).c_str()); // axis_rectangle_1x2 + std::printf("%s\n", shape({V(0,0), V(1,1), V(0,2), V(-1,1)}).c_str()); // rotated_square + std::printf("%s\n", shape({V(0,0), V(2,2), V(1,3), V(-1,1)}).c_str()); // rotated_rectangle + std::printf("%s\n", shape({V(0,0), V(3,0), V(2,1), V(1,1)}).c_str()); // trapezoid axis? (not cyclic in general, shape only) + std::printf("%s\n", shape({V(5,0), V(3,4), V(-4,3), V(0,-5)}).c_str()); // other_cyclic (r = 5) + std::printf("%s\n", shape({V(0,0), V(1,0), V(2,0), V(0,1)}).c_str()); // collinear_triple + // det: unit square cocircular -> 0; d at centre-ish inside -> positive + std::printf("det sq %d\n", (int)lattice_det({V(0,0), V(1,0), V(1,1), V(0,1)})); + // a,b,c CCW in (col,-row): (0,0),(2,0) at row 0; (1,-1) is row 1 -> (1,1) in world up? use row -> -y + std::printf("det in %d\n", (int)(lattice_det({V(0,0), V(2,0), V(1,-2), V(1,-1)}) > 0)); // world (0,0),(2,0),(1,2), query (1,1): inside + std::printf("fma 0.1*3 exact %d, 10*4999 exact %d\n", product_exact(3, 0.1), product_exact(4999, 10)); + std::printf("bits 10=%d 0.1=%d 1=%d\n", significant_bits(10), significant_bits(0.1), significant_bits(1)); +} diff --git a/docs/benchmarks/2026-09-27/21b-ties/scripts/counter_selftest.cpp b/docs/benchmarks/2026-09-27/21b-ties/scripts/counter_selftest.cpp new file mode 100644 index 00000000..03a7dd89 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/scripts/counter_selftest.cpp @@ -0,0 +1,16 @@ +#include +#include +using namespace terrain::mesh; +using namespace terrain::mesh::prof_ties; +static MeshVertex V(double c, double r) { return MeshVertex{c, r}; } +int main() { + std::printf("%s\n", shape({V(5,0), V(3,4), V(-5,0), V(0,-5)}).c_str()); // other_cyclic + phase = 2; + record(V(0,0), V(1,0), V(1,1), V(0,1), 0.1, 0.1, 0, true); // frame inexact? 0.1*1 exact; row 1 -> exact + record(V(0,3), V(7,0), V(7,3), V(0,0), 0.1, 0.1, 0, true); // 7*0.1, 3*0.1 inexact + record(V(0,0), V(1,0), V(1,1), V(0,1), 10, 20, 0, true); // dx != dy + record(V(0,0), V(20000,0), V(1,1), V(0,1), 10, 10, 0, true); // spread > 2^14, det != 0 + record(V(0,0), V(1,0), V(1,1), V(0,1), 10, 10, 1, false); // sign differs + record(V(0.5,0), V(1,0), V(1,1), V(0,1), 10, 10, 0, true); // 3 nodes + for (auto& [k, n] : counts) std::printf("%s %llu\n", k.c_str(), n); +} diff --git a/docs/benchmarks/2026-09-27/21b-ties/scripts/instrument.patch b/docs/benchmarks/2026-09-27/21b-ties/scripts/instrument.patch new file mode 100644 index 00000000..6b3652fc --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/scripts/instrument.patch @@ -0,0 +1,315 @@ +diff --git a/include/terrain/mesh/lawson.hpp b/include/terrain/mesh/lawson.hpp +index 2796055..c4dae58 100644 +--- a/include/terrain/mesh/lawson.hpp ++++ b/include/terrain/mesh/lawson.hpp +@@ -24,6 +24,8 @@ + #include + #include + ++#include ++ + #include + #include + #include +@@ -67,10 +69,18 @@ template + // the frame. Both branches flip exactly when the same exact determinant, + // on the frame points of the quad in mesh CCW order, is positive; a flip + // lowers the sum of signed lifted volumes, so the process terminates. +- if (K::orient2d(a, b, c) == pred::Orientation::CounterClockwise) +- return K::incircle(a, b, c, d) == pred::Incircle::Inside; +- return K::orient2d(b, a, d) == pred::Orientation::CounterClockwise +- && K::incircle(b, a, d, c) == pred::Incircle::Inside; ++ const MeshVertex va = v[tri[e]], vb = v[tri[(e + 1) % 3]], vc = v[tri[(e + 2) % 3]], ++ vd = v[m.triangles()[u][(j + 2) % 3]]; ++ if (K::orient2d(a, b, c) == pred::Orientation::CounterClockwise) { ++ const auto r = K::incircle(a, b, c, d); ++ prof_ties::record(va, vb, vc, vd, f.dx, f.dy, static_cast(r), pred::prof_ties::last_exact); ++ return r == pred::Incircle::Inside; ++ } ++ if (K::orient2d(b, a, d) != pred::Orientation::CounterClockwise) ++ return false; ++ const auto r = K::incircle(b, a, d, c); ++ prof_ties::record(vb, va, vd, vc, f.dx, f.dy, static_cast(r), pred::prof_ties::last_exact); ++ return r == pred::Incircle::Inside; + } + + } // namespace detail +diff --git a/include/terrain/mesh/prof_ties.hpp b/include/terrain/mesh/prof_ties.hpp +new file mode 100644 +index 0000000..5d516f4 +--- /dev/null ++++ b/include/terrain/mesh/prof_ties.hpp +@@ -0,0 +1,212 @@ ++#pragma once ++ ++// 21b tie classification (docs/increments/21-parallel-refine.md, section 3 QW2). ++// Local measurement instrumentation, never committed: kept as ++// docs/benchmarks/2026-09-27/21b-ties/scripts/instrument.patch. ++// ++// At every incircle call made by must_flip, record, on the four MeshVertex ++// corners as passed to the kernel (a, b, c counter-clockwise in the frame, ++// d the query point): ++// - how many are DEM nodes (MeshVertex::is_node); ++// - for four nodes: the exact lattice determinant on (col, -row) translated ++// to d (in __int128, so it is exact whatever the spread), QW2's ++// conditions (dx == dy; every |difference from d| <= 2^14; the frame ++// products col * dx and row * dy exact, checked per call with fma), the ++// agreement of the lattice sign with the kernel's answer, and the shape; ++// - whether the kernel took the exact path, and its answer. ++// Counters are keyed by phase (0 legalise_all, 1 quality, 2 refine loop). The ++// callers are serial (the split phase), so plain globals suffice. ++ ++#include ++ ++#include ++#include ++#include ++#include ++#include ++#include ++#include ++#include ++#include ++#include ++#include ++ ++namespace terrain::mesh::prof_ties { ++ ++inline int phase = -1; ++inline std::map counts; ++inline double g_dx = 0, g_dy = 0; ++inline std::size_t g_rows = 0, g_cols = 0; ++inline long long max_spread_all4 = 0, max_spread_exact4 = 0; ++ ++inline void reset(double dx, double dy, std::size_t rows, std::size_t cols) { ++ counts.clear(); ++ g_dx = dx; ++ g_dy = dy; ++ g_rows = rows; ++ g_cols = cols; ++ max_spread_all4 = max_spread_exact4 = 0; ++} ++ ++// Significant bits of a positive finite double: 53 less the trailing zeros of ++// its 53-bit integer significand. ++inline int significant_bits(double x) { ++ int e = 0; ++ const double m = std::frexp(x, &e); ++ const auto s = static_cast(std::ldexp(m, 53)); ++ return 53 - std::countr_zero(s); ++} ++ ++inline bool product_exact(double a, double b) { ++ const double p = a * b; ++ return std::fma(a, b, -p) == 0.0; ++} ++ ++__extension__ typedef __int128 I; ++ ++inline I lattice_det(const std::array& q) { ++ const auto X = [](const MeshVertex& v) { return static_cast(static_cast(v.col)); }; ++ const auto Y = [](const MeshVertex& v) { return -static_cast(static_cast(v.row)); }; ++ const auto& d = q[3]; ++ const I adx = X(q[0]) - X(d), ady = Y(q[0]) - Y(d); ++ const I bdx = X(q[1]) - X(d), bdy = Y(q[1]) - Y(d); ++ const I cdx = X(q[2]) - X(d), cdy = Y(q[2]) - Y(d); ++ const I al = adx * adx + ady * ady, bl = bdx * bdx + bdy * bdy, cl = cdx * cdx + cdy * cdy; ++ return al * (bdx * cdy - cdx * bdy) + bl * (cdx * ady - adx * cdy) + cl * (adx * bdy - bdx * ady); ++} ++ ++inline long long orient_ll(const MeshVertex& a, const MeshVertex& b, const MeshVertex& c) { ++ const auto ax = static_cast(a.col), ay = -static_cast(a.row); ++ const auto bx = static_cast(b.col), by = -static_cast(b.row); ++ const auto cx = static_cast(c.col), cy = -static_cast(c.row); ++ return (bx - ax) * (cy - ay) - (by - ay) * (cx - ax); ++} ++ ++// Shape of a four-node quad with a zero lattice determinant. ++inline std::string shape(const std::array& q) { ++ for (int i = 0; i < 4; ++i) ++ for (int j = i + 1; j < 4; ++j) ++ for (int k = j + 1; k < 4; ++k) ++ if (orient_ll(q[i], q[j], q[k]) == 0) ++ return "collinear_triple"; ++ std::set cols, rows; ++ for (const auto& v : q) { ++ cols.insert(static_cast(v.col)); ++ rows.insert(static_cast(v.row)); ++ } ++ const auto c = [&](int i) { return std::pair{static_cast(q[i].col), static_cast(q[i].row)}; }; ++ constexpr std::array, 3> pairings{{{0, 1, 2, 3}, {0, 2, 1, 3}, {0, 3, 1, 2}}}; ++ for (const auto& p : pairings) { ++ const auto [x0, y0] = c(p[0]); ++ const auto [x1, y1] = c(p[1]); ++ const auto [x2, y2] = c(p[2]); ++ const auto [x3, y3] = c(p[3]); ++ const bool same_mid = x0 + x1 == x2 + x3 && y0 + y1 == y2 + y3; ++ const bool same_len = (x1 - x0) * (x1 - x0) + (y1 - y0) * (y1 - y0) ++ == (x3 - x2) * (x3 - x2) + (y3 - y2) * (y3 - y2); ++ if (!same_mid || !same_len) ++ continue; ++ // Rectangle with diagonals (p0,p1), (p2,p3); sides p0-p2 and p0-p3. ++ const auto s1 = (x2 - x0) * (x2 - x0) + (y2 - y0) * (y2 - y0); ++ const auto s2 = (x3 - x0) * (x3 - x0) + (y3 - y0) * (y3 - y0); ++ const bool square = s1 == s2; ++ if (cols.size() == 2 && rows.size() == 2) { ++ const long long w = *cols.rbegin() - *cols.begin(), h = *rows.rbegin() - *rows.begin(); ++ return std::string(square ? "axis_square" : "axis_rectangle") + "_" + std::to_string(std::min(w, h)) ++ + "x" + std::to_string(std::max(w, h)); ++ } ++ return square ? "rotated_square" : "rotated_rectangle"; ++ } ++ if (cols.size() == 2 || rows.size() == 2) ++ return "isosceles_trapezoid_axis"; // two parallel sides on rows or columns ++ // A cyclic quad with a pair of parallel (disjoint) chords is an isosceles ++ // trapezoid; the diagonals cross, so testing all three pairings is safe. ++ for (const auto& p : pairings) { ++ const auto [x0, y0] = c(p[0]); ++ const auto [x1, y1] = c(p[1]); ++ const auto [x2, y2] = c(p[2]); ++ const auto [x3, y3] = c(p[3]); ++ if ((x1 - x0) * (y3 - y2) - (y1 - y0) * (x3 - x2) == 0) ++ return "isosceles_trapezoid_other"; ++ } ++ return "other_cyclic"; ++} ++ ++inline void add(const std::string& key, unsigned long long n = 1) { counts[key] += n; } ++ ++inline void record(MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d, double dx, double dy, int result, ++ bool exact) { ++ const std::string P = "p" + std::to_string(phase) + "."; ++ const std::string R = result > 0 ? "inside" : result < 0 ? "outside" : "cocircular"; ++ const std::array q{a, b, c, d}; ++ int nodes = 0; ++ for (const auto& v : q) ++ nodes += v.is_node() ? 1 : 0; ++ const std::string path = exact ? "exact" : "filtered"; ++ add(P + "calls"); ++ add(P + path + ".calls"); ++ add(P + path + ".result." + R); ++ add(P + path + ".nodes" + std::to_string(nodes)); ++ add(P + path + ".nodes" + std::to_string(nodes) + ".result." + R); ++ if (nodes != 4) ++ return; ++ long long spread = 0; ++ for (int i = 0; i < 3; ++i) { ++ spread = std::max(spread, std::llabs(static_cast(q[i].col) - static_cast(d.col))); ++ spread = std::max(spread, std::llabs(static_cast(q[i].row) - static_cast(d.row))); ++ } ++ if (phase == 2) { ++ max_spread_all4 = std::max(max_spread_all4, spread); ++ if (exact) ++ max_spread_exact4 = std::max(max_spread_exact4, spread); ++ } ++ bool frame_exact = true; ++ for (const auto& v : q) ++ frame_exact = frame_exact && product_exact(v.col, dx) && product_exact(v.row, dy); ++ const bool eq = dx == dy; ++ const bool small = spread <= (1LL << 14); ++ const bool qualifies = eq && small && frame_exact; ++ const I det = lattice_det(q); ++ const int s = det > 0 ? 1 : det < 0 ? -1 : 0; ++ const std::string N = P + path + ".n4."; ++ add(N + (eq ? "dx_eq_dy" : "dx_ne_dy")); ++ add(N + (small ? "spread_le_2^14" : "spread_gt_2^14")); ++ add(N + (frame_exact ? "frame_exact" : "frame_inexact")); ++ add(N + (qualifies ? "qw2_qualifies" : "qw2_refuses")); ++ add(N + "lattice_det_" + (s > 0 ? "pos" : s < 0 ? "neg" : "zero")); ++ add(N + (s == result ? "lattice_sign_agrees" : "lattice_sign_differs")); ++ if (qualifies && s != result) ++ add(N + "qw2_qualifies_and_differs"); ++ if (s == 0) ++ add(N + "shape." + shape(q)); ++ if (exact && s != 0) ++ add(N + "exact_path_lattice_nonzero"); ++} ++ ++inline void dump(std::size_t rounds, std::size_t inserted, std::size_t flips) { ++ const char* path = std::getenv("RASPUTIN_TIES_OUT"); ++ if (path == nullptr) ++ return; ++ std::FILE* f = std::fopen(path, "a"); ++ if (f == nullptr) ++ return; ++ const int bits = significant_bits(g_dx); ++ const int width = std::bit_width(std::max(g_rows, g_cols) - 1); ++ std::fprintf(f, ++ "{\"dx\": %.17g, \"dy\": %.17g, \"rows\": %zu, \"cols\": %zu, \"dx_significant_bits\": %d, " ++ "\"index_bit_width\": %d, \"global_frame_exact_condition\": %s, \"rounds\": %zu, " ++ "\"inserted\": %zu, \"flips\": %zu, \"max_spread_all4_refine\": %lld, " ++ "\"max_spread_exact4_refine\": %lld, \"counts\": {", ++ g_dx, g_dy, g_rows, g_cols, bits, width, ++ (g_dx == g_dy && bits + width <= 53) ? "true" : "false", rounds, inserted, flips, ++ max_spread_all4, max_spread_exact4); ++ bool first = true; ++ for (const auto& [k, n] : counts) { ++ std::fprintf(f, "%s\"%s\": %llu", first ? "" : ", ", k.c_str(), n); ++ first = false; ++ } ++ std::fprintf(f, "}}\n"); ++ std::fclose(f); ++} ++ ++} // namespace terrain::mesh::prof_ties +diff --git a/include/terrain/predicates/kernel.hpp b/include/terrain/predicates/kernel.hpp +index b794353..a4f7eff 100644 +--- a/include/terrain/predicates/kernel.hpp ++++ b/include/terrain/predicates/kernel.hpp +@@ -33,6 +33,12 @@ concept GeometryKernel = requires(Point2 a, Point2 b, Point2 c, Point2 d) { + { K::incircle(a, b, c, d) } -> std::same_as; + }; + ++// 21b tie classification (local instrumentation, never committed): whether the ++// last FilteredKernel incircle call took the exact path. ++namespace prof_ties { ++inline bool last_exact = false; ++} // namespace prof_ties ++ + namespace detail { + + // The naive 3x3 lifted determinant, translated so d is the origin. Positive iff +@@ -172,8 +178,10 @@ private: + const double permanent = detail::incircle_permanent(a, b, c, d); + + if (std::fabs(det) > incircle_bound_a * permanent) { ++ prof_ties::last_exact = false; + return incircle_of_sign(det); + } ++ prof_ties::last_exact = true; + return E::incircle_ccw(a, b, c, d); + } + }; +diff --git a/include/terrain/refinement/refine.hpp b/include/terrain/refinement/refine.hpp +index c0c6a25..9003583 100644 +--- a/include/terrain/refinement/refine.hpp ++++ b/include/terrain/refinement/refine.hpp +@@ -282,9 +282,12 @@ template + const auto since = [](clock::time_point t0) { + return std::chrono::duration(clock::now() - t0).count(); + }; ++ mesh::prof_ties::reset(g.delta_x(), g.delta_y(), g.rows(), g.cols()); ++ mesh::prof_ties::phase = 0; + auto t0 = clock::now(); + out.flips = mesh::legalise_all(m, frame, [](std::uint32_t) {}); + out.legalise_seconds = since(t0); ++ mesh::prof_ties::phase = 1; + if (options.min_angle_deg > 0.0) { + t0 = clock::now(); + const auto q = mesh::improve( +@@ -297,6 +300,7 @@ template + std::vector results; + std::set> footed; // (row, col), R2 step 5 + mesh::FlipStack flip_stack; // one buffer for every legalise_around ++ mesh::prof_ties::phase = 2; + std::vector active(m.triangle_count()); + for (std::uint32_t t = 0; t < active.size(); ++t) + active[t] = t; +@@ -378,6 +382,7 @@ template + detail::rebuild_active(touched, skipped, active); + } + ++ mesh::prof_ties::dump(out.rounds, out.inserted, out.flips); + // By the stopping rule a void triangle holds no valid node, so `uncovered` + // sums zeros unless that rule changes; it is reported so a change shows. + for (const ScanResult& r : results) { diff --git a/docs/benchmarks/2026-09-27/21b-ties/scripts/prof_driver.py b/docs/benchmarks/2026-09-27/21b-ties/scripts/prof_driver.py new file mode 100644 index 00000000..ca8eacd5 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/scripts/prof_driver.py @@ -0,0 +1,51 @@ +"""prof_driver.py --pkg DIR --threads N --repeat K [--pause S] -- + +bench.py's child technique (pkg first on sys.path, wrap cli.refine), but the +wrapped refine is called K times with the same arguments, each timed, and the +RefineOutcome phase seconds printed per call as one JSON line (PHASES ...). +--pause S sleeps S seconds before the first call and prints the pid, so +`sample` can attach to a steady window of refine calls. +""" +import json +import os +import sys +import time + +argv = sys.argv[1:] +split = argv.index("--") +own, rasputin = argv[:split], argv[split + 1 :] +pkg = own[own.index("--pkg") + 1] +threads = int(own[own.index("--threads") + 1]) +repeat = int(own[own.index("--repeat") + 1]) +pause = float(own[own.index("--pause") + 1]) if "--pause" in own else 0.0 +sys.meta_path[:] = [f for f in sys.meta_path if "ScikitBuild" not in type(f).__name__] +sys.path.insert(0, pkg) +import tin_engine._core as core # noqa: E402 +import tin_engine.cli as cli # noqa: E402 + +assert core.__file__.startswith(pkg), core.__file__ +real = vars(cli)["refine"] + + +def refine(*args, **kwargs): + kwargs["threads"] = threads + if pause: + print(f"PID {os.getpid()}", file=sys.stderr, flush=True) + time.sleep(pause) + out = None + for i in range(repeat): + t0 = time.perf_counter() + out = real(*args, **kwargs) + wall = time.perf_counter() - t0 + rec = dict(call=i, threads=threads, refine_s=wall, legalise_s=out.legalise_seconds, + quality_s=out.quality_seconds, scan_s=out.scan_seconds, + split_s=out.split_seconds, rounds=out.rounds, inserted=out.inserted, + flips=out.flips, max_error=out.max_error, + triangles=len(out.triangles), vertices=len(out.vertices)) + rec["rest_s"] = wall - rec["legalise_s"] - rec["quality_s"] - rec["scan_s"] - rec["split_s"] + print("PHASES " + json.dumps(rec), file=sys.stderr, flush=True) + return out + + +vars(cli)["refine"] = refine +cli.app(args=rasputin, prog_name="rasputin", standalone_mode=False) diff --git a/docs/benchmarks/2026-09-27/21b-ties/scripts/tables.py b/docs/benchmarks/2026-09-27/21b-ties/scripts/tables.py new file mode 100644 index 00000000..5eda6dba --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-ties/scripts/tables.py @@ -0,0 +1,102 @@ +"""tables.py data/ties_quarter_t1.jsonl data/ties_tile_t1.jsonl > data/tables.md + +Summarises instrument.patch's RASPUTIN_TIES_OUT records (one JSON line per +refine call) into the README's tables. Phase p2 is the refine loop, p0 +legalise_all, p1 the quality pass. +""" +import json +import re +import sys +from pathlib import Path + + +def load(path): + return json.loads(Path(path).read_text().splitlines()[-1]) + + +def shape_group(name): + if name.startswith("axis_square"): + return "axis-aligned square" + if name.startswith("axis_rectangle"): + return "axis-aligned rectangle, not square" + return { + "rotated_square": "rotated square", + "rotated_rectangle": "rotated rectangle, not square", + "isosceles_trapezoid_axis": "isosceles trapezoid, parallel sides on rows or columns", + "isosceles_trapezoid_other": "isosceles trapezoid, other direction", + "other_cyclic": "other cyclic quad (no parallel sides)", + "collinear_triple": "three corners collinear", + }[name] + + +def pct(n, d): + return f"{100.0 * n / d:.3f} %" if d else "-" + + +def main(paths): + recs = {Path(p).stem.split("_")[1]: load(p) for p in paths} + names = list(recs) + print("| measure | " + " | ".join(names) + " |") + print("|---|" + "---:|" * len(names)) + rows = [] + for ph, label in (("p2", "refine loop"), ("p0", "legalise_all"), ("p1", "quality pass")): + def g(r, k, ph=ph): + return r["counts"].get(f"{ph}.{k}", 0) + + rows += [ + (f"{label}: incircle calls", lambda r, g=g: f"{g(r, 'calls'):,}"), + (f"{label}: exact path", lambda r, g=g: f"{g(r, 'exact.calls'):,} ({pct(g(r, 'exact.calls'), g(r, 'calls'))})"), + (f"{label}: exact path, Cocircular", lambda r, g=g: f"{g(r, 'exact.result.cocircular'):,}"), + (f"{label}: exact path, four node corners", lambda r, g=g: f"{g(r, 'exact.nodes4'):,}"), + (f"{label}: exact path, lattice det = 0", lambda r, g=g: f"{g(r, 'exact.n4.lattice_det_zero'):,}"), + (f"{label}: exact path, lattice det != 0", lambda r, g=g: f"{g(r, 'exact.n4.exact_path_lattice_nonzero'):,}"), + (f"{label}: exact path, QW2 conditions hold", lambda r, g=g: f"{g(r, 'exact.n4.qw2_qualifies'):,}"), + (f"{label}: filtered path, four node corners", lambda r, g=g: f"{g(r, 'filtered.nodes4'):,}"), + (f"{label}: filtered path, QW2 conditions hold", lambda r, g=g: f"{g(r, 'filtered.n4.qw2_qualifies'):,}"), + (f"{label}: all calls QW2 would answer", lambda r, g=g: ( + lambda q: f"{q:,} ({pct(q, g(r, 'calls'))})")( + g(r, 'exact.n4.qw2_qualifies') + g(r, 'filtered.n4.qw2_qualifies'))), + (f"{label}: calls with fewer than four node corners", lambda r, g=g: f"{sum(g(r, f'{p}.nodes{k}') for p in ('exact', 'filtered') for k in range(4)):,}"), + (f"{label}: lattice sign differs from kernel (four nodes)", lambda r, g=g: f"{g(r, 'exact.n4.lattice_sign_differs') + g(r, 'filtered.n4.lattice_sign_differs'):,}"), + ] + rows += [ + ("max spread from d, four-node calls, refine loop (nodes)", lambda r: f"{r['max_spread_all4_refine']}"), + ("max spread from d, exact-path calls, refine loop (nodes)", lambda r: f"{r['max_spread_exact4_refine']}"), + ("dx, dy; rows x cols", lambda r: f"{r['dx']:g}, {r['dy']:g}; {r['rows']} x {r['cols']}"), + ("frame-exact sufficient condition (bits of dx + bit_width) <= 53", lambda r: f"{r['dx_significant_bits']} + {r['index_bit_width']} = {r['dx_significant_bits'] + r['index_bit_width']}: {r['global_frame_exact_condition']}"), + ("rounds, inserted, flips", lambda r: f"{r['rounds']}, {r['inserted']:,}, {r['flips']:,}"), + ] + for label, f in rows: + print(f"| {label} | " + " | ".join(f(recs[n]) for n in names) + " |") + + print() + print("Shapes of the exact-path (tie) quads, refine loop:") + print() + print("| shape | " + " | ".join(names) + " |") + print("|---|" + "---:|" * len(names)) + groups = {} + for n in names: + tot = recs[n]["counts"].get("p2.exact.calls", 0) + for k, v in recs[n]["counts"].items(): + m = re.fullmatch(r"p2\.exact\.n4\.shape\.(.+)", k) + if m: + groups.setdefault(shape_group(m.group(1)), {}).setdefault(n, 0) + groups[shape_group(m.group(1))][n] += v + order = sorted(groups, key=lambda s: -groups[s].get(names[0], 0)) + for s in order: + cells = [] + for n in names: + v = groups[s].get(n, 0) + cells.append(f"{v:,} ({pct(v, recs[n]['counts'].get('p2.exact.calls', 0))})") + print(f"| {s} | " + " | ".join(cells) + " |") + print() + print("Most frequent axis-aligned sizes (short x long side, nodes), refine loop exact path:") + print() + for n in names: + sizes = sorted(((v, k.split(".")[-1]) for k, v in recs[n]["counts"].items() + if k.startswith("p2.exact.n4.shape.axis_")), reverse=True)[:6] + print(f"- {n}: " + ", ".join(f"{s.removeprefix('axis_')} {v:,}" for v, s in sizes)) + + +if __name__ == "__main__": + main(sys.argv[1:]) From 93a8066fdffbbee25ca3e3463b5d01109b0c34c0 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 13:31:41 +0200 Subject: [PATCH 02/12] 21: QW2's premise measured (21b-ties) Co-Authored-By: Claude Opus 5.5 --- docs/increments/21-parallel-refine.md | 9 +++++++-- 1 file changed, 7 insertions(+), 2 deletions(-) diff --git a/docs/increments/21-parallel-refine.md b/docs/increments/21-parallel-refine.md index 7597c778..2f39f9b2 100644 --- a/docs/increments/21-parallel-refine.md +++ b/docs/increments/21-parallel-refine.md @@ -784,8 +784,13 @@ Each is answered or placed. cheaper exact decision for quads whose corners are all nodes? That is 5 % of refine and about 17 % of the serial phase."* **Answered: yes,** an int64 lattice incircle under four stated conditions, - bit-identical to today (QW2, 21b). The premise is still an inference; the - classification measurement in QW2 comes first. + bit-identical to today (QW2, 21b). **The premise is measured** + (`docs/benchmarks/2026-09-27/21b-ties/README.md`): every exact-path tie + has four node corners, lattice determinant 0 and meets QW2's conditions + (146,962 of 146,962 on the quarter circle, 154,502 of 154,502 on the tile); + QW2 would answer 99.97 % and 100 % of all refine-loop incircle calls. Only + 38 % of the ties are axis-aligned rectangles, so the general determinant is + needed. The time saved is not measured. 4. *"The scan loses 16.5 ms of 61.8 ms to imbalance at 8 threads. Is the chunking free to change, given that results are written per slot? And does From 042d3e441d46b337d8361d18e2b2f1651b1966e1 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 13:43:10 +0200 Subject: [PATCH 03/12] red: increment 21b (lattice incircle) 21b is bit-identical, so 14's T6 and 18's golden digests stay untouched as the output oracle. The new suite test_mesh_lattice_incircle is 21b's invariant-critical suite (section 7). It pins the interface QW2 left open; the pins are under "Pinned by the red suite (21b)" in 21-parallel-refine.md: lattice_frame(dx, dy, rows, cols) decides once per refine call whether the integer path may answer, and lattice_incircle(a, b, c, d, frame) returns an optional Incircle. A directly built LatticeFrame{dx, dy} never enables the path. Agreement with DetriaExact::incircle_ccw on the frame doubles, for every quad it answers: an exhaustive 5 x 5 neighbourhood at the origin and at the benchmark's far corner in four exact frames; the tie shapes 21b-ties lists, in every rotation; every quad on the radius-5, radius-sqrt(65) and radius-8085 lattice circles, and those quads with d moved one node; random uniform and near-circle quads up to the spread bound. It answers at a spread of exactly 2^14 and refuses at 2^14 + 1, for each corner, axis and direction. A named test fails on a (col, row) determinant. It refuses an off-node corner in each position; dx != dy; dx that is zero, negative, infinite or NaN; an inexact frame (0.1, with rows and cols tested separately, and a 53-bit dx); and a direct frame. It must answer where QW2's sufficient condition holds, including exactly 53 bits. A search finds a lattice tie that the 0.1 frame breaks, which shows why that refusal is needed. must_flip against a copy of today's version, on random lattice meshes with ties and off-node feet, in five frames: the same decisions edge by edge, and the same legalise_all flips, writes and mesh. A counting kernel shows that must_flip asks the kernel nothing when the integer path answers. Mutation round against a scratch implementation that is not committed. All 18 mutants were killed: sign inverted; (col, row); >= at the bound; a bound of 2^14 + 1; a bound of 2^15; the bound checked on one side only; the bound checked on columns only; no exact-frame check; the exact-frame test at < 53; the exact-frame test on columns only; on rows only; no node check; a node check that skips d; dx != dy accepted; dx <= 0 accepted; a direct frame enabling the path; no call in must_flip; a determinant in doubles. The scratch implementation also passes test_mesh_lawson, test_mesh_lawson_stack and test_mesh_quality. The suite is registered only once lawson.hpp names lattice_incircle, as in 21a, and the header is a configure dependency. That the guard fires was checked by building with the scratch header in place (the suite built and passed), then restoring the header. The suite starts no threads, so it is not in the TSan job. Co-Authored-By: Claude Opus 5.5 --- docs/increments/21-parallel-refine.md | 54 ++ tests/cpp/CMakeLists.txt | 21 + tests/cpp/unit/test_mesh_lattice_incircle.cpp | 838 ++++++++++++++++++ 3 files changed, 913 insertions(+) create mode 100644 tests/cpp/unit/test_mesh_lattice_incircle.cpp diff --git a/docs/increments/21-parallel-refine.md b/docs/increments/21-parallel-refine.md index 2f39f9b2..e1cc6e60 100644 --- a/docs/increments/21-parallel-refine.md +++ b/docs/increments/21-parallel-refine.md @@ -487,6 +487,60 @@ is no difference beyond noise (refine 515 vs 519 ms), as QW1 predicts. The mesh sha256 that `bench.py` prints for the tile and the quarter circle was the same before and after. +### Pinned by the red suite (21b) + +QW2 names `lattice_incircle(a, b, c, d) -> std::optional` but not +where `dx == dy` and the exact-frame decision come from. `@tester` chose the +following in the red step; `tests/cpp/unit/test_mesh_lattice_incircle.cpp` +holds it. Both functions are in `namespace terrain::mesh` and reachable +through `include/terrain/mesh/lawson.hpp`. Whether they live there or in a new +header beside it that `lawson.hpp` includes is `@developer`'s choice. + +- `[[nodiscard]] LatticeFrame lattice_frame(double dx, double dy, std::size_t rows, std::size_t cols) noexcept` + returns a frame with the same `dx` and `dy` that also records, once per + refine call, whether the integer path may answer on it. It may answer only + when `dx == dy`, `dx` is finite and positive, and `col * dx` and `row * dy` + are exact for every node with `col < cols` and `row < rows`. It must answer + wherever QW2's sufficient condition holds: the significant bits of `dx` plus + `bit_width(max(rows, cols) - 1)` are at most 53. The suite tests exactly 53. + Between that condition and exactness, for example `dx = 0.1` on a + 3 × 3 grid, which is exact but fails the condition, nothing is pinned. + `refine` builds its frame with this function instead of + `LatticeFrame{g.delta_x(), g.delta_y()}`. No test can see that, because the + output is bit-identical by design. The reviewer checks it by reading the + code, and `@perf`'s acceptance run shows it as time saved. +- **A `LatticeFrame` built directly never enables the integer path.** + `LatticeFrame{dx, dy}`, the way every caller builds one today, keeps compiling + and keeps today's kernel path. The existing Lawson and quality suites + therefore still exercise that path, and the 21b suite uses such a frame as + its "without". +- `[[nodiscard]] std::optional lattice_incircle(MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d, const LatticeFrame& f) noexcept`. + Precondition: `a, b, c` strictly counter-clockwise on `(col, -row)`, which + `must_flip`'s triangle is (the `LatticeMesh` invariant). It returns empty + unless `f` came from an enabling `lattice_frame`, all four corners are nodes, + and `|col_x - col_d|` and `|row_x - row_d|` are at most 2^14 for each `x` in + `a, b, c`. The bound is measured from `d`, so `a` and `b` may be 2^15 apart. + When it answers, the answer is the sign `DetriaExact::incircle_ccw` gives on + the frame points `(col * dx, -(row * dy))`. +- **`must_flip` calls it first.** When it answers, `must_flip` returns + `answer == Inside` and asks the kernel nothing, not even the frame + `orient2d`. The suite counts kernel calls to check this. When it refuses, + today's path runs unchanged. + +The suite is registered only once `lawson.hpp` names `lattice_incircle`, as +in 21a. It starts no threads, so it is not in the TSan job. + +**Mutation round** against a scratch implementation that is not committed. All +of these were killed: the sign inverted; the determinant on `(col, row)`; +refusal at exactly 2^14 (`>=`); a bound of 2^14 + 1; a bound of 2^15; the +bound checked on one side only; the bound checked on columns only; the +exact-frame check dropped; the exact-frame test at `< 53`; the exact-frame +test on columns only, and on rows only; the node check dropped; the node check +skipping `d`; `dx != dy` accepted; `dx <= 0` accepted; a directly built frame +enabling the path; no call in `must_flip`. The determinant computed in doubles +rather than `int64` was also killed, and only by the radius-8085 circle, where +spreads come close to 2^14. + ## 4. Determinism levels Today's contract (14 R5, 14b R1, tested by 14's T6 and 18's T3 golden diff --git a/tests/cpp/CMakeLists.txt b/tests/cpp/CMakeLists.txt index a21cc2c2..42d6d7b3 100644 --- a/tests/cpp/CMakeLists.txt +++ b/tests/cpp/CMakeLists.txt @@ -291,3 +291,24 @@ add_terrain_backend_test(test_refinement_active unit/test_refinement_active.cpp) target_link_libraries(test_refinement_active PRIVATE Threads::Threads) add_terrain_backend_test(test_mesh_lawson_stack unit/test_mesh_lawson_stack.cpp) target_link_libraries(test_mesh_lawson_stack PRIVATE Threads::Threads) + +# The integer incircle, increment 21b (docs/increments/21-parallel-refine.md, +# section 3 QW2 and "Pinned by the red suite (21b)"). Bit-identical: 14's T6 +# and 18's golden digests are the output oracle and stay as they are. +# test_mesh_lattice_incircle is the invariant-critical suite: agreement with +# DetriaExact on the frame doubles wherever lattice_incircle answers, refusal +# wherever QW2's conditions fail, and must_flip's decisions as today's. It +# starts no threads, so it is not in the TSan job. +# +# RED UNTIL THE INTERFACE EXISTS: registered only once lawson.hpp names +# lattice_incircle (defined there or in a header it includes), so the red does +# not break the rest of the build. The header is a configure dependency, so the +# green commit's edit registers the suite without a manual cmake. The review +# drops this guard. +set_property(DIRECTORY APPEND PROPERTY CMAKE_CONFIGURE_DEPENDS + "${PROJECT_SOURCE_DIR}/include/terrain/mesh/lawson.hpp") +file(READ "${PROJECT_SOURCE_DIR}/include/terrain/mesh/lawson.hpp" _lawson_21b) +string(FIND "${_lawson_21b}" "lattice_incircle" _lattice_incircle_at) +if(NOT _lattice_incircle_at EQUAL -1) + add_terrain_backend_test(test_mesh_lattice_incircle unit/test_mesh_lattice_incircle.cpp) +endif() diff --git a/tests/cpp/unit/test_mesh_lattice_incircle.cpp b/tests/cpp/unit/test_mesh_lattice_incircle.cpp new file mode 100644 index 00000000..bf22f061 --- /dev/null +++ b/tests/cpp/unit/test_mesh_lattice_incircle.cpp @@ -0,0 +1,838 @@ +// Increment 21b (docs/increments/21-parallel-refine.md, section 3 QW2, section +// 7's 21b row, and "Pinned by the red suite (21b)"): an integer incircle for +// quads whose four corners are DEM nodes, answered at the top of must_flip. +// This is 21b's invariant-critical suite. +// +// 21b is bit-identical. Where lattice_incircle answers, its answer must be the +// sign DetriaExact computes on the frame doubles (col * dx, -(row * dy)), +// because that is the number must_flip decides on today. Where any of QW2's +// conditions fails it must refuse, and must_flip must then run today's path. +// +// Interface pinned here (reachable through include/terrain/mesh/lawson.hpp): +// +// namespace terrain::mesh { +// // A frame for this grid that also records, once, whether QW2 may answer +// // on it: dx == dy, finite and > 0, and col * dx, row * dy exact for +// // every node with col < cols, row < rows. +// [[nodiscard]] LatticeFrame lattice_frame(double dx, double dy, +// std::size_t rows, std::size_t cols) noexcept; +// // Precondition: a, b, c strictly counter-clockwise on (col, -row). +// [[nodiscard]] std::optional lattice_incircle( +// MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d, const LatticeFrame& f) noexcept; +// } +// +// A LatticeFrame built directly, LatticeFrame{dx, dy} as every caller builds +// it today, never enables the integer path. +// +// The oracle is DetriaExact::incircle_ccw on the frame points, computed here +// from col * dx and -(row * dy), not through LatticeFrame::at. For today's +// must_flip the oracle is a copy of it (lawson.hpp at 93a8066). + +#include +#include +#include + +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +using terrain::Point2; +using terrain::TriangleIndices; +using terrain::mesh::kNoNeighbour; +using terrain::mesh::lattice_frame; +using terrain::mesh::lattice_incircle; +using terrain::mesh::LatticeFrame; +using terrain::mesh::LatticeMesh; +using terrain::mesh::LatticeVertex; +using terrain::mesh::MeshVertex; +using terrain::pred::DefaultKernel; +using terrain::pred::DetriaExact; +using terrain::pred::Incircle; +using terrain::pred::Orientation; + +namespace { + +constexpr std::int64_t kBound = std::int64_t{1} << 14; // QW2's spread bound, in nodes +constexpr std::size_t kBench = 5051; // the 1 m benchmark's rows and cols + +// A node at (col, row). MeshVertex is (col, row); LatticeVertex is {row, col}. +MeshVertex node(std::int64_t col, std::int64_t row) { + return MeshVertex{static_cast(col), static_cast(row)}; +} + +// The frame point, written out rather than taken from LatticeFrame::at. +Point2 fp(double dx, double dy, MeshVertex v) { return Point2{v.col * dx, -(v.row * dy)}; } + +// Twice the signed area on (col, -row), exact for nodes of the sizes used here. +std::int64_t iorient(MeshVertex a, MeshVertex b, MeshVertex c) { + const auto ax = static_cast(a.col), ay = -static_cast(a.row); + const auto bx = static_cast(b.col), by = -static_cast(b.row); + const auto cx = static_cast(c.col), cy = -static_cast(c.row); + return (bx - ax) * (cy - ay) - (by - ay) * (cx - ax); +} + +// The oracle: DetriaExact on the frame doubles. Empty if a, b, c is not +// counter-clockwise in the frame (DetriaExact's precondition), which under an +// exact frame and a counter-clockwise lattice triangle cannot happen. +std::optional oracle(double dx, double dy, MeshVertex a, MeshVertex b, MeshVertex c, + MeshVertex d) { + const Point2 pa = fp(dx, dy, a), pb = fp(dx, dy, b), pc = fp(dx, dy, c), pd = fp(dx, dy, d); + if (DetriaExact::orient2d(pa, pb, pc) != Orientation::CounterClockwise) + return std::nullopt; + return DetriaExact::incircle_ccw(pa, pb, pc, pd); +} + +std::string name(Incircle s) { + return s == Incircle::Inside ? "Inside" : s == Incircle::Outside ? "Outside" : "Cocircular"; +} +std::string describe(MeshVertex v) { + return "(col " + std::to_string(v.col) + ", row " + std::to_string(v.row) + ")"; +} + +// Agreement bookkeeping: every quad the caller hands in must be answered, with +// the oracle's sign. The first failure is kept for the report, and the three +// signs are counted so a run that saw only one of them is visible. +struct Agreement { + std::size_t quads = 0, refused = 0, disagreed = 0; + std::array by_sign{}; // Outside, Cocircular, Inside + std::string first; + + void check(const LatticeFrame& f, double dx, MeshVertex a, MeshVertex b, MeshVertex c, + MeshVertex d) { + ++quads; + const auto want = oracle(dx, dx, a, b, c, d); + const auto got = lattice_incircle(a, b, c, d, f); + if (!want || !got || *got != *want) { + (got ? disagreed : refused) += 1; + if (first.empty()) + first = describe(a) + " " + describe(b) + " " + describe(c) + " d " + describe(d) + + " dx " + std::to_string(dx) + ": got " + + (got ? name(*got) : std::string{"nullopt"}) + ", oracle " + + (want ? name(*want) : std::string{"not ccw in frame"}); + return; + } + ++by_sign[static_cast(static_cast(*want) + 1)]; + } + void require_all_agree() const { + CAPTURE(quads, refused, disagreed, first); + REQUIRE(quads > 0); + REQUIRE(refused == 0); + REQUIRE(disagreed == 0); + } +}; + +// The four points in cyclic order, turned counter-clockwise on (col, -row). +std::array ccw_cycle(std::array p) { + if (iorient(p[0], p[1], p[2]) < 0) + std::reverse(p.begin(), p.end()); + return p; +} + +// Every lattice point on the circle of squared radius r2 about (cx, cy), in +// angular order. +std::vector circle(std::int64_t cx, std::int64_t cy, std::int64_t r2) { + std::vector> pts; + const auto root = [](std::int64_t v) { // floor(sqrt(v)), exact for the sizes here + auto s = static_cast(std::sqrt(static_cast(v))); + while (s * s > v) + --s; + while ((s + 1) * (s + 1) <= v) + ++s; + return s; + }; + const auto r = root(r2); + for (std::int64_t x = -r; x <= r; ++x) { + const std::int64_t y = root(r2 - x * x); + if (x * x + y * y != r2) + continue; + for (const std::int64_t sy : {y, -y}) { + pts.emplace_back(std::atan2(static_cast(sy), static_cast(x)), node(cx + x, cy + sy)); + if (y == 0) + break; + } + } + std::sort(pts.begin(), pts.end(), [](const auto& l, const auto& r) { return l.first < r.first; }); + std::vector out; + for (const auto& [angle, v] : pts) + out.push_back(v); + return out; +} + +// Every choice of four points from a cyclically ordered set, as quads in +// cyclic order. +std::vector> quads_of(const std::vector& ring) { + std::vector> out; + const std::size_t n = ring.size(); + for (std::size_t i = 0; i < n; ++i) + for (std::size_t j = i + 1; j < n; ++j) + for (std::size_t k = j + 1; k < n; ++k) + for (std::size_t l = k + 1; l < n; ++l) + out.push_back(ccw_cycle({ring[i], ring[j], ring[k], ring[l]})); + return out; +} + +// A cocircular quad in each of its four rotations, (a, b, c) the triangle and +// d the fourth corner, as must_flip hands them over. +template +void for_each_rotation(const std::array& q, Fn&& fn) { + for (std::size_t k = 0; k < 4; ++k) + fn(q[k], q[(k + 1) % 4], q[(k + 2) % 4], q[(k + 3) % 4]); +} + +} // namespace + +// ---------------------------------------------------------------- agreement + +TEST_CASE("lattice_incircle agrees with DetriaExact on every quad of a 5 x 5 neighbourhood", + "[lattice_incircle][agreement]") { + // Every ordered (a, b, c, d) of distinct nodes with a, b, c strictly + // counter-clockwise: all ties, all near misses and all clear cases that + // fit. At the origin and at the far corner of the benchmark's grid, in + // four exact frames (dx = 10 is the benchmark's). + const double dx = GENERATE(1.0, 10.0, 0.5, 0.375); + const std::int64_t base = GENERATE(std::int64_t{0}, std::int64_t{5046}); + CAPTURE(dx, base); + const LatticeFrame f = lattice_frame(dx, dx, kBench, kBench); + std::vector nodes; + for (std::int64_t r = 0; r < 5; ++r) + for (std::int64_t c = 0; c < 5; ++c) + nodes.push_back(node(base + c, base + r)); + Agreement agree; + for (const auto& a : nodes) + for (const auto& b : nodes) + for (const auto& c : nodes) { + if (iorient(a, b, c) <= 0) + continue; + for (const auto& d : nodes) + if (d != a && d != b && d != c) + agree.check(f, dx, a, b, c, d); + } + agree.require_all_agree(); + // All three answers occur, ties included, so a constant answer cannot pass. + REQUIRE(agree.by_sign[0] > 1000); + REQUIRE(agree.by_sign[1] > 1000); + REQUIRE(agree.by_sign[2] > 1000); +} + +TEST_CASE("lattice_incircle answers Cocircular on the tie shapes 21b-ties found, in every rotation", + "[lattice_incircle][agreement][ties]") { + // docs/benchmarks/2026-09-27/21b-ties/README.md, "The shapes of the + // exact-path quads": axis-aligned squares and rectangles (1x1, 1x2, 2x2, + // 1x3, 2x3, and legalise_all's 40x40 and 10x40), rotated squares and + // rectangles, isosceles trapezoids with parallel sides on rows or columns + // and in another direction, and other cyclic quads. Each is a tie: the + // answer must be Cocircular, and that must be the oracle's answer too. + const double dx = GENERATE(1.0, 10.0); + CAPTURE(dx); + const LatticeFrame f = lattice_frame(dx, dx, kBench, kBench); + const auto rect = [](std::int64_t c, std::int64_t r, std::int64_t w, std::int64_t h) { + return ccw_cycle({node(c, r), node(c + w, r), node(c + w, r + h), node(c, r + h)}); + }; + const auto quad = [](std::array, 4> p) { + return ccw_cycle({node(p[0][0], p[0][1]), node(p[1][0], p[1][1]), node(p[2][0], p[2][1]), + node(p[3][0], p[3][1])}); + }; + const std::vector>> shapes{ + {"axis square 1x1", rect(7, 9, 1, 1)}, + {"axis square 2x2", rect(7, 9, 2, 2)}, + {"axis square 40x40", rect(40, 80, 40, 40)}, + {"axis rectangle 1x2", rect(7, 9, 1, 2)}, + {"axis rectangle 2x1", rect(7, 9, 2, 1)}, + {"axis rectangle 1x3", rect(7, 9, 1, 3)}, + {"axis rectangle 2x3", rect(7, 9, 2, 3)}, + {"axis rectangle 10x40", rect(40, 80, 10, 40)}, + {"rotated square", quad({{{10, 11}, {11, 10}, {12, 11}, {11, 12}}})}, + {"rotated square, side (3, 2)", quad({{{10, 12}, {13, 10}, {15, 13}, {12, 15}}})}, + {"rotated rectangle", quad({{{10, 11}, {11, 10}, {13, 12}, {12, 13}}})}, + {"isosceles trapezoid, parallel sides on rows", quad({{{10, 10}, {14, 10}, {13, 12}, {11, 12}}})}, + {"isosceles trapezoid, parallel sides on columns", quad({{{10, 10}, {12, 11}, {12, 13}, {10, 14}}})}, + // On the radius-5 circle about (20, 20): chords (5,0)-(0,5) and + // (4,-3)-(-3,4) are both along (-1, 1). + {"isosceles trapezoid, parallel sides along a diagonal", + quad({{{25, 20}, {20, 25}, {17, 24}, {24, 17}}})}, + // Radius 5 again: no two sides parallel. + {"other cyclic quad", quad({{{25, 20}, {23, 24}, {15, 20}, {16, 17}}})}, + }; + for (const auto& [label, q] : shapes) { + CAPTURE(std::string{label}); + for_each_rotation(q, [&](MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d) { + CAPTURE(describe(a), describe(b), describe(c), describe(d)); + REQUIRE(iorient(a, b, c) > 0); + REQUIRE(oracle(dx, dx, a, b, c, d) == Incircle::Cocircular); + REQUIRE(lattice_incircle(a, b, c, d, f) == Incircle::Cocircular); + }); + } +} + +TEST_CASE("lattice_incircle agrees with DetriaExact on every quad of lattice circles, and one node off them", + "[lattice_incircle][agreement][ties]") { + // Every four of the 12 lattice points on the radius-5 circle and of the 16 + // on the radius-sqrt(65) circle, in every rotation: all ties. Then d moved + // one node in each axis direction, which puts it strictly inside or outside. + // The scaled radius-5 circle (x 1617 = 3 * 7^2 * 11, radius 8085, still 12 + // lattice points) has spreads up to 16170 nodes, just under 2^14: there a + // one-node move of d changes the determinant by about 2^40 in terms near + // 2^56, beyond what plain double arithmetic on the differences resolves. + const double dx = GENERATE(1.0, 10.0); + const std::int64_t r2 = GENERATE(std::int64_t{25}, std::int64_t{65}, std::int64_t{25} * 1617 * 1617); + CAPTURE(dx, r2); + const std::int64_t centre = 9000; // keeps every point at col, row >= 0 + const auto ring = circle(centre, centre, r2); + REQUIRE(ring.size() == (r2 == 65 ? 16u : 12u)); + const LatticeFrame f = lattice_frame(dx, dx, 20000, 20000); + Agreement ties, moved; + std::size_t cocircular = 0; + for (const auto& q : quads_of(ring)) + for_each_rotation(q, [&](MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d) { + ties.check(f, dx, a, b, c, d); + cocircular += lattice_incircle(a, b, c, d, f) == Incircle::Cocircular; + for (const auto& [mc, mr] : {std::pair{1, 0}, std::pair{-1, 0}, std::pair{0, 1}, std::pair{0, -1}}) { + const MeshVertex e{d.col + mc, d.row + mr}; + if (e != a && e != b && e != c) + moved.check(f, dx, a, b, c, e); + } + }); + ties.require_all_agree(); + moved.require_all_agree(); + REQUIRE(cocircular == ties.quads); + REQUIRE(ties.by_sign[1] == ties.quads); + REQUIRE(moved.by_sign[0] > 0); + REQUIRE(moved.by_sign[2] > 0); +} + +TEST_CASE("lattice_incircle agrees with DetriaExact on random quads up to the spread bound", + "[lattice_incircle][agreement][random]") { + // Two families, spreads drawn log-uniformly from 1 to 2^14 nodes: + // - uniform: a, b, c anywhere within the spread of d; + // - near-circle: four points rounded from one circle, so the + // determinant is small against its terms, where rounding would show. + const std::uint32_t seed = GENERATE(range(1u, 5u)); + const double dx = GENERATE(1.0, 10.0, 0.25, 3.0); + CAPTURE(seed, dx); + const LatticeFrame f = lattice_frame(dx, dx, 1u << 16, 1u << 16); + std::mt19937_64 rng{seed}; + std::uniform_real_distribution unit{0.0, 1.0}; + const std::int64_t mid = std::int64_t{1} << 15; + Agreement uniform, near; + const auto within = [](MeshVertex p, MeshVertex d) { + return std::fabs(p.col - d.col) <= static_cast(kBound) + && std::fabs(p.row - d.row) <= static_cast(kBound); + }; + for (int i = 0; i < 4000; ++i) { + const auto s = static_cast(std::exp2(14.0 * unit(rng))); + const auto off = [&] { return static_cast(std::llround((2.0 * unit(rng) - 1.0) * s)); }; + // One draw per statement: argument evaluation order is unspecified. + std::array o{}; + for (auto& x : o) + x = off(); + const MeshVertex d = node(mid + o[0], mid + o[1]); + MeshVertex a = node(mid + o[0] + o[2], mid + o[1] + o[3]); + MeshVertex b = node(mid + o[0] + o[4], mid + o[1] + o[5]); + const MeshVertex c = node(mid + o[0] + o[6], mid + o[1] + o[7]); + if (iorient(a, b, c) == 0 || d == a || d == b || d == c) + continue; + if (iorient(a, b, c) < 0) + std::swap(a, b); + uniform.check(f, dx, a, b, c, d); + } + for (int i = 0; i < 4000; ++i) { + const double radius = std::exp2(1.0 + 12.9 * unit(rng)); // diameter <= 2^14 + const double cx = static_cast(mid) + 1000.0 * unit(rng); + const double cy = static_cast(mid) + 1000.0 * unit(rng); + std::array p{}; + for (auto& v : p) { + const double t = 2.0 * std::numbers::pi * unit(rng); + v = MeshVertex{std::round(cx + radius * std::cos(t)), std::round(cy + radius * std::sin(t))}; + } + auto [a, b, c, d] = p; + if (iorient(a, b, c) == 0 || d == a || d == b || d == c || !within(a, d) || !within(b, d) + || !within(c, d)) + continue; + if (iorient(a, b, c) < 0) + std::swap(a, b); + near.check(f, dx, a, b, c, d); + } + uniform.require_all_agree(); + near.require_all_agree(); + REQUIRE(uniform.quads > 3000); + REQUIRE(near.quads > 3000); + REQUIRE(near.by_sign[0] > 100); + REQUIRE(near.by_sign[2] > 100); +} + +// ---------------------------------------------------------------- the spread bound + +namespace { + +// a, b, c with one corner `far` from d along one axis, placed at `pos` in the +// counter-clockwise triple; the other two close to d. +std::array far_triple(MeshVertex d, std::int64_t dc, std::int64_t dr, unsigned pos) { + const auto dcol = static_cast(d.col), drow = static_cast(d.row); + std::array t{node(dcol + dc, drow + dr), node(dcol + 1, drow + 2), node(dcol - 2, drow + 1)}; + if (iorient(t[0], t[1], t[2]) < 0) + std::swap(t[1], t[2]); + std::rotate(t.begin(), t.begin() + (3 - pos) % 3, t.end()); + return t; +} + +} // namespace + +TEST_CASE("lattice_incircle answers at a spread of exactly 2^14 nodes and refuses at 2^14 + 1", + "[lattice_incircle][bound]") { + // QW2: every coordinate difference from d at most 2^14 nodes, in absolute + // value. For each corner position, each axis and each direction: at 2^14 + // it answers, with the oracle's sign; one node further it refuses. + const std::int64_t spread = GENERATE(kBound, kBound + 1); + const unsigned pos = GENERATE(0u, 1u, 2u); + const auto [dc, dr] = GENERATE(std::pair{1, 0}, std::pair{-1, 0}, std::pair{0, 1}, std::pair{0, -1}); + CAPTURE(spread, pos, dc, dr); + const double dx = 10.0; + const LatticeFrame f = lattice_frame(dx, dx, 1u << 16, 1u << 16); + const MeshVertex d = node(1 << 15, 1 << 15); + const auto [a, b, c] = far_triple(d, dc * spread, dr * spread, pos); + REQUIRE(iorient(a, b, c) > 0); + const auto got = lattice_incircle(a, b, c, d, f); + if (spread == kBound) { + REQUIRE(got.has_value()); + REQUIRE(got == oracle(dx, dx, a, b, c, d)); + } else { + REQUIRE_FALSE(got.has_value()); + } +} + +TEST_CASE("lattice_incircle measures the spread from d, not between the other corners", + "[lattice_incircle][bound]") { + // a and b are 2^14 from d on either side, so 2^15 apart; the bound is on + // differences from d (QW2), and this quad is within it. + const LatticeFrame f = lattice_frame(1.0, 1.0, 1u << 16, 1u << 16); + const MeshVertex d = node(1 << 15, 1 << 15); + MeshVertex a = node((1 << 15) + kBound, (1 << 15) + 3), b = node((1 << 15) - kBound, (1 << 15) + 1), + c = node(1 << 15, (1 << 15) - kBound); + if (iorient(a, b, c) < 0) + std::swap(a, b); + const auto got = lattice_incircle(a, b, c, d, f); + REQUIRE(got.has_value()); + REQUIRE(got == oracle(1.0, 1.0, a, b, c, d)); +} + +// ---------------------------------------------------------------- orientation + +TEST_CASE("lattice_incircle takes the determinant on (col, -row), as the frame is oriented", + "[lattice_incircle][orientation]") { + // In the frame (col, -row): a (0, -2), b (2, -2), c (1, 0) turn + // counter-clockwise; their circle has centre (1, -1.25) and radius 1.25, + // and d (1, -1) is 0.25 from the centre, strictly inside. On (col, row) + // the four points are the mirror image: a, b, c turn clockwise, so the + // same determinant formula has the opposite sign and says Outside. That is + // the answer an implementation taking (col, row) gives here. + const LatticeFrame f = lattice_frame(1.0, 1.0, 16, 16); + const MeshVertex a = node(0, 2), b = node(2, 2), c = node(1, 0), d = node(1, 1); + REQUIRE(iorient(a, b, c) > 0); + REQUIRE(oracle(1.0, 1.0, a, b, c, d) == Incircle::Inside); + REQUIRE(lattice_incircle(a, b, c, d, f) == Incircle::Inside); + // And a point clearly outside, which (col, row) would call Inside. + const MeshVertex e = node(5, 5); + REQUIRE(oracle(1.0, 1.0, a, b, c, e) == Incircle::Outside); + REQUIRE(lattice_incircle(a, b, c, e, f) == Incircle::Outside); +} + +// ---------------------------------------------------------------- refusals + +TEST_CASE("lattice_incircle refuses a quad with an off-node corner, in each position", + "[lattice_incircle][refuse]") { + // QW2 answers only when all four corners are nodes (MeshVertex::is_node). + // Off-node corners are the domain's arc vertices (16 R2) and 20b's feet. + const LatticeFrame f = lattice_frame(10.0, 10.0, kBench, kBench); + std::array q{node(10, 12), node(12, 12), node(11, 10), node(11, 11)}; + REQUIRE(iorient(q[0], q[1], q[2]) > 0); + REQUIRE(lattice_incircle(q[0], q[1], q[2], q[3], f).has_value()); // all nodes: answers + const std::size_t which = GENERATE(0u, 1u, 2u, 3u); + const auto [dc, dr] = GENERATE(std::pair{0.5, 0.0}, std::pair{0.0, 0.25}, + std::pair{std::ldexp(1.0, -40), 0.0}, std::pair{0.0, std::ldexp(1.0, -40)}); + CAPTURE(which, dc, dr); + q[which] = MeshVertex{q[which].col + dc, q[which].row + dr}; + REQUIRE_FALSE(q[which].is_node()); + REQUIRE(terrain::mesh::orient_sign(q[0], q[1], q[2]) > 0); // still counter-clockwise + REQUIRE_FALSE(lattice_incircle(q[0], q[1], q[2], q[3], f).has_value()); +} + +TEST_CASE("lattice_incircle refuses every frame QW2 excludes", "[lattice_incircle][refuse]") { + // The quad is a plain all-node quad the benchmark frame answers; only the + // frame changes. + const MeshVertex a = node(3, 7), b = node(6, 7), c = node(4, 3), d = node(4, 5); + REQUIRE(iorient(a, b, c) > 0); + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(10.0, 10.0, 100, 100)).has_value()); + + const double inf = std::numeric_limits::infinity(); + const double nan = std::numeric_limits::quiet_NaN(); + SECTION("dx != dy: the frame is not the lattice times one constant") { + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(1.0, 2.0, 100, 100)).has_value()); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(10.0, 5.0, 100, 100)).has_value()); + REQUIRE_FALSE( + lattice_incircle(a, b, c, d, lattice_frame(1.0, std::nextafter(1.0, 2.0), 100, 100)).has_value()); + } + SECTION("dx not finite and positive") { + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(0.0, 0.0, 100, 100)).has_value()); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(-1.0, -1.0, 100, 100)).has_value()); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(inf, inf, 100, 100)).has_value()); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(nan, nan, 100, 100)).has_value()); + } + SECTION("an inexact frame: some col * dx or row * dy of the grid rounds") { + // 3 * 0.1 != 0.3 exactly; the grid holds col and row 3. + REQUIRE(std::fma(3.0, 0.1, -(3.0 * 0.1)) != 0.0); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(0.1, 0.1, 100, 100)).has_value()); + // 0.1 has 52 significant bits, so a grid two nodes wide is exact on + // that axis; the other axis is long enough to round. Rows and cols + // must both count. The quads lie inside those grids. + const MeshVertex ta = node(0, 5), tb = node(1, 5), tc = node(0, 3), td = node(1, 4); // cols 0..1 + const MeshVertex wa = node(3, 1), wb = node(6, 1), wc = node(4, 0), wd = node(5, 0); // rows 0..1 + REQUIRE(iorient(ta, tb, tc) > 0); + REQUIRE(iorient(wa, wb, wc) > 0); + REQUIRE_FALSE(lattice_incircle(ta, tb, tc, td, lattice_frame(0.1, 0.1, 100, 2)).has_value()); + REQUIRE_FALSE(lattice_incircle(wa, wb, wc, wd, lattice_frame(0.1, 0.1, 2, 100)).has_value()); + // A 53-bit dx: 3 * dx needs 54 bits. + const double full = 1.0 + std::ldexp(1.0, -52); + REQUIRE(std::fma(3.0, full, -(3.0 * full)) != 0.0); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(full, full, 100, 100)).has_value()); + } + SECTION("a LatticeFrame built directly never enables the integer path") { + REQUIRE_FALSE(lattice_incircle(a, b, c, d, LatticeFrame{10.0, 10.0}).has_value()); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, LatticeFrame{1.0, 1.0}).has_value()); + } +} + +TEST_CASE("lattice_frame enables the integer path wherever QW2's sufficient condition holds", + "[lattice_incircle][frame]") { + // QW2: the frame is exact when the significant bits of dx plus + // bit_width(max(rows, cols) - 1) are at most 53. Inside that, it must + // answer (that is the speed-up, and must_flip's zero-kernel-call test + // below depends on it). + const MeshVertex a = node(0, 1), b = node(1, 1), c = node(0, 0), d = node(1, 0); // unit square + REQUIRE(iorient(a, b, c) > 0); + // The benchmark: 3 + 13 = 16. + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(10.0, 10.0, kBench, kBench)) == Incircle::Cocircular); + // Exactly 53: a 52-bit dx on a 2 x 2 grid (bit_width(1) = 1). + const double dx52 = 1.0 + std::ldexp(1.0, -51); + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(dx52, dx52, 2, 2)) == Incircle::Cocircular); + // A power of two on a 2^20 grid: 1 + 20. + const double tiny = std::ldexp(1.0, -30); + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(tiny, tiny, 1u << 20, 1u << 20)) == Incircle::Cocircular); + // 0.1 on a 2 x 2 grid (52 + 1 = 53): col and row are 0 or 1, exact. + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(0.1, 0.1, 2, 2)) == Incircle::Cocircular); +} + +TEST_CASE("an inexact frame changes the answer on some lattice tie, so its refusal is load-bearing", + "[lattice_incircle][refuse]") { + // QW2's determinism argument: without the exact-frame condition the + // integer answer is the true lattice answer, and the rounded frame's + // answer can differ on a tie. Find one: a cocircular lattice quad on which + // DetriaExact at dx = dy = 0.1 does not say Cocircular. lattice_incircle + // must refuse it (it would otherwise say Cocircular and change a flip). + const double dx = 0.1; + const LatticeFrame f = lattice_frame(dx, dx, 100, 100); + std::optional> found; + for (const std::int64_t r2 : {25, 50, 65, 85, 125}) + for (const auto& q : quads_of(circle(40, 40, r2))) + for_each_rotation(q, [&](MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d) { + const auto o = oracle(dx, dx, a, b, c, d); + if (!found && o && *o != Incircle::Cocircular) + found = std::array{a, b, c, d}; + }); + REQUIRE(found.has_value()); + const auto [a, b, c, d] = *found; + CAPTURE(describe(a), describe(b), describe(c), describe(d)); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, f).has_value()); +} + +// ---------------------------------------------------------------- must_flip + +namespace { + +// must_flip as of 93a8066, the oracle for decisions. +bool reference_must_flip(const LatticeMesh& m, std::uint32_t t, unsigned e, const LatticeFrame& f) { + const auto u = m.neighbours(t)[e]; + if (u == kNoNeighbour || m.is_constrained(t, e)) + return false; + const auto& tri = m.triangles()[t]; + unsigned j = 0; + while (m.triangles()[u][j] != tri[(e + 1) % 3]) + ++j; + const auto v = m.vertices(); + const Point2 a = fp(f.dx, f.dy, v[tri[e]]), b = fp(f.dx, f.dy, v[tri[(e + 1) % 3]]), + c = fp(f.dx, f.dy, v[tri[(e + 2) % 3]]), d = fp(f.dx, f.dy, v[m.triangles()[u][(j + 2) % 3]]); + if (DefaultKernel::orient2d(a, b, c) == Orientation::CounterClockwise) + return DefaultKernel::incircle(a, b, c, d) == Incircle::Inside; + return DefaultKernel::orient2d(b, a, d) == Orientation::CounterClockwise + && DefaultKernel::incircle(b, a, d, c) == Incircle::Inside; +} + +// DefaultKernel, counting its calls: must_flip asks it nothing when the +// integer path answers. +struct CountingKernel { + static inline std::size_t calls = 0; + static Orientation orient2d(const Point2& a, const Point2& b, const Point2& c) { + ++calls; + return DefaultKernel::orient2d(a, b, c); + } + static Incircle incircle(const Point2& a, const Point2& b, const Point2& c, const Point2& d) { + ++calls; + return DefaultKernel::incircle(a, b, c, d); + } +}; +static_assert(terrain::pred::GeometryKernel); + +// An (n+1) x (n+1) grid of nodes `step` apart from (base, base), two triangles +// per cell, the ring constrained. Mixed diagonals give ties both ways. +LatticeMesh grid(std::uint32_t n, std::uint32_t step, std::uint32_t base) { + std::vector v; + for (std::uint32_t r = 0; r <= n; ++r) + for (std::uint32_t c = 0; c <= n; ++c) + v.push_back(LatticeVertex{base + r * step, base + c * step}); + const auto at = [n](std::uint32_t r, std::uint32_t c) { return r * (n + 1) + c; }; + std::vector t; + for (std::uint32_t r = 0; r < n; ++r) + for (std::uint32_t c = 0; c < n; ++c) { + const auto tl = at(r, c), tr = at(r, c + 1), bl = at(r + 1, c), br = at(r + 1, c + 1); + if ((r + c) % 2 == 0) { + t.push_back({tl, bl, br}); + t.push_back({tl, br, tr}); + } else { + t.push_back({tl, bl, tr}); + t.push_back({bl, br, tr}); + } + } + std::vector bits; + std::vector> masks; + const std::uint32_t lo = base, hi = base + n * step; + for (const auto& tri : t) { + std::uint8_t b = 0; + for (unsigned k = 0; k < 3; ++k) { + const auto p = v[tri[k]], q = v[tri[(k + 1) % 3]]; + if ((p.row == q.row && (p.row == lo || p.row == hi)) || (p.col == q.col && (p.col == lo || p.col == hi))) + b = static_cast(b | (1u << k)); + } + bits.push_back(b); + masks.push_back({b & 1u, (b >> 1) & 1u, (b >> 2) & 1u}); + } + auto m = LatticeMesh::build(std::move(v), std::move(t), std::move(bits), std::move(masks)); + REQUIRE(m.has_value()); + return *m; +} + +// Insert node p strictly inside a triangle or on an edge; false when p is a +// vertex. Exact orientation on (col, -row), so off-node vertices are fine. +bool insert_node(LatticeMesh& m, LatticeVertex p) { + const MeshVertex x{p}; + for (std::uint32_t t = 0; t < m.triangle_count(); ++t) { + const auto& tri = m.triangles()[t]; + const auto v = m.vertices(); + std::array o{}; + for (unsigned k = 0; k < 3; ++k) + o[k] = terrain::mesh::orient_sign(v[tri[k]], v[tri[(k + 1) % 3]], x); + if (o[0] < 0 || o[1] < 0 || o[2] < 0) + continue; + const int zeros = (o[0] == 0) + (o[1] == 0) + (o[2] == 0); + if (zeros >= 2) + return false; + if (zeros == 0) + m.split_inside(t, p); + else + m.split_edge(t, o[0] == 0 ? 0u : o[1] == 0 ? 1u : 2u, x); + return true; + } + return false; +} + +// Split a constrained ring edge at an off-node point strictly inside it (a +// foot, as 20b inserts them). false if no ring edge is long enough. +bool insert_foot(LatticeMesh& m, std::mt19937& rng) { + const auto n = m.triangle_count(); + const auto start = static_cast(rng() % n); + for (std::uint32_t i = 0; i < n; ++i) { + const auto t = (start + i) % n; + for (unsigned k = 0; k < 3; ++k) { + if (!m.is_constrained(t, k)) + continue; + const auto p = m.vertices()[m.triangles()[t][k]], q = m.vertices()[m.triangles()[t][(k + 1) % 3]]; + const double len = std::fabs(p.col - q.col) + std::fabs(p.row - q.row); // axis-aligned + if (len < 2.0) + continue; + const double s = (std::floor(len / 2.0) + 0.5) / len; // off-node, strictly inside + m.split_edge(t, k, MeshVertex{p.col + s * (q.col - p.col), p.row + s * (q.row - p.row)}); + return true; + } + } + return false; +} + +struct Decisions { + std::size_t edges = 0, flips = 0, ties = 0, differ = 0; + std::string first; +}; + +// Every (t, e): must_flip under the lattice frame against today's must_flip. +void compare_decisions(const LatticeMesh& m, double dx, double dy, std::size_t grid_n, Decisions& out) { + const LatticeFrame with = lattice_frame(dx, dy, grid_n, grid_n); + const LatticeFrame without{dx, dy}; + for (std::uint32_t t = 0; t < m.triangle_count(); ++t) + for (unsigned e = 0; e < 3; ++e) { + if (m.neighbours(t)[e] == kNoNeighbour || m.is_constrained(t, e)) + continue; + ++out.edges; + const bool want = reference_must_flip(m, t, e, without); + const bool got = terrain::mesh::detail::must_flip(m, t, e, with); + out.flips += want; + const auto& tri = m.triangles()[t]; + const auto u = m.neighbours(t)[e]; + unsigned j = 0; + while (m.triangles()[u][j] != tri[(e + 1) % 3]) + ++j; + const auto v = m.vertices(); + out.ties += DefaultKernel::incircle(fp(dx, dy, v[tri[0]]), fp(dx, dy, v[tri[1]]), fp(dx, dy, v[tri[2]]), + fp(dx, dy, v[m.triangles()[u][(j + 2) % 3]])) + == Incircle::Cocircular; + if (got != want && out.differ++ == 0) + out.first = "t " + std::to_string(t) + " e " + std::to_string(e) + ": got " + + std::to_string(got) + ", today " + std::to_string(want); + } +} + +struct FrameCase { + double dx, dy; + bool integer; // whether lattice_frame enables the integer path on it +}; +const std::array kFrameCases{FrameCase{1.0, 1.0, true}, FrameCase{10.0, 10.0, true}, + FrameCase{0.5, 0.5, true}, FrameCase{0.1, 0.1, false}, + FrameCase{10.0, 5.0, false}}; + +} // namespace + +TEST_CASE("must_flip decides as today on random lattice meshes, with ties and off-node feet", + "[lattice_incircle][must_flip]") { + // Random nodes inserted into a mixed-diagonal grid without legalising, so + // many edges must flip and many are ties; every few insertions a foot + // splits a ring edge off the lattice. After every insertion, every + // interior unconstrained edge is decided by must_flip under + // lattice_frame(dx, dy, ...) and by a copy of today's must_flip under + // LatticeFrame{dx, dy}: the decisions must be equal. Near the origin and + // near the benchmark's far corner. + const std::uint32_t seed = GENERATE(range(1u, 7u)); + const std::size_t frame = GENERATE(range(std::size_t{0}, kFrameCases.size())); + const std::uint32_t base = GENERATE(0u, 5000u); + const auto [dx, dy, integer] = kFrameCases[frame]; + CAPTURE(seed, dx, dy, base); + auto m = grid(4, 6, base); // nodes base .. base + 24 + std::mt19937 rng{seed}; + Decisions dec; + std::size_t inserted = 0, feet = 0; + for (int attempt = 0; attempt < 120; ++attempt) { + if (attempt % 15 == 7 && insert_foot(m, rng)) + ++feet; + const LatticeVertex p{base + static_cast(rng() % 25), + base + static_cast(rng() % 25)}; + if (!insert_node(m, p)) + continue; + ++inserted; + compare_decisions(m, dx, dy, kBench, dec); + } + CAPTURE(dec.edges, dec.flips, dec.ties, dec.first); + REQUIRE(dec.differ == 0); + REQUIRE(inserted >= 40); + REQUIRE(feet >= 3); + REQUIRE(dec.flips >= 100); + REQUIRE(dec.ties >= 100); + REQUIRE(dec.edges - dec.flips - dec.ties >= 100); +} + +TEST_CASE("legalise_all under lattice_frame makes the same flips, writes and mesh as under LatticeFrame{dx, dy}", + "[lattice_incircle][must_flip]") { + // The mesh-level form of bit-identity: the whole Lawson run, not one + // decision at a time. LatticeFrame{dx, dy} never enables the integer path, + // so it is today's run. + const std::uint32_t seed = GENERATE(range(1u, 7u)); + const std::size_t frame = GENERATE(range(std::size_t{0}, kFrameCases.size())); + const auto [dx, dy, integer] = kFrameCases[frame]; + CAPTURE(seed, dx, dy); + auto m = grid(5, 5, 0); + std::mt19937 rng{seed}; + for (int attempt = 0; attempt < 150; ++attempt) { + if (attempt % 20 == 3) + insert_foot(m, rng); + insert_node(m, LatticeVertex{static_cast(rng() % 26), static_cast(rng() % 26)}); + } + LatticeMesh today = m; + std::vector w_today, w_new; + const auto f_today = terrain::mesh::legalise_all( + today, LatticeFrame{dx, dy}, [&](std::uint32_t s) { w_today.push_back(s); }); + const auto f_new = terrain::mesh::legalise_all( + m, lattice_frame(dx, dy, 26, 26), [&](std::uint32_t s) { w_new.push_back(s); }); + REQUIRE(f_today >= 20); + REQUIRE(f_new == f_today); + REQUIRE(w_new == w_today); + REQUIRE(m.triangle_count() == today.triangle_count()); + for (std::uint32_t t = 0; t < m.triangle_count(); ++t) { + CAPTURE(t); + REQUIRE(m.triangles()[t] == today.triangles()[t]); + REQUIRE(m.neighbours(t) == today.neighbours(t)); + } +} + +TEST_CASE("must_flip asks the kernel nothing when lattice_incircle answers, and today's questions when it refuses", + "[lattice_incircle][must_flip]") { + // QW2: lattice_incircle is called at the top of must_flip, and when it + // answers, must_flip also skips the frame orient2d (the mesh triangle is + // counter-clockwise in integers, and so in an exact frame). That is the + // whole speed-up, so it is pinned: zero kernel calls on an all-node quad + // under an enabling frame; kernel calls under a direct LatticeFrame, an + // excluded frame, or with an off-node corner. + auto m = grid(4, 3, 0); + std::mt19937 rng{7}; + for (int i = 0; i < 40; ++i) + insert_node(m, LatticeVertex{static_cast(rng() % 13), static_cast(rng() % 13)}); + const auto run = [&m](const LatticeFrame& f, bool node_quads_only, std::size_t min_edges) { + std::size_t asked = 0, edges = 0; + for (std::uint32_t t = 0; t < m.triangle_count(); ++t) + for (unsigned e = 0; e < 3; ++e) { + const auto u = m.neighbours(t)[e]; + if (u == kNoNeighbour || m.is_constrained(t, e)) + continue; + bool nodes = true; + for (const auto i : m.triangles()[t]) + nodes = nodes && m.vertices()[i].is_node(); + for (const auto i : m.triangles()[u]) + nodes = nodes && m.vertices()[i].is_node(); + if (nodes != node_quads_only) + continue; + ++edges; + CountingKernel::calls = 0; + const bool got = terrain::mesh::detail::must_flip(m, t, e, f); + REQUIRE(got == reference_must_flip(m, t, e, LatticeFrame{f.dx, f.dy})); + asked += CountingKernel::calls; + } + REQUIRE(edges >= min_edges); + return asked; + }; + REQUIRE(run(lattice_frame(10.0, 10.0, kBench, kBench), true, 20) == 0); + REQUIRE(run(LatticeFrame{10.0, 10.0}, true, 20) > 0); + REQUIRE(run(lattice_frame(0.1, 0.1, 100, 100), true, 20) > 0); + REQUIRE(run(lattice_frame(10.0, 5.0, kBench, kBench), true, 20) > 0); + for (int i = 0; i < 6; ++i) + REQUIRE(insert_foot(m, rng)); + REQUIRE(run(lattice_frame(10.0, 10.0, kBench, kBench), false, 3) > 0); +} From f7073224d1987d5c81ed340c6f131d747ce34a1a Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 18:30:09 +0200 Subject: [PATCH 04/12] green: increment 21b (lattice incircle) QW2 (docs/increments/21-parallel-refine.md, section 3, and "Pinned by the red suite (21b)"), in lawson.hpp: - LatticeFrame becomes a class with a private flag that only lattice_frame() sets, so LatticeFrame{dx, dy} keeps compiling and never enables the integer path. A constructor rather than a third aggregate member keeps GCC's -Wmissing-field-initializers out of every existing two-value brace init. - lattice_frame(dx, dy, rows, cols) enables it when dx == dy, dx is a positive normal double, the significant bits of dx plus bit_width(max(rows, cols) - 1) are at most 53 (QW2's sufficient condition, used as the test), and the largest col * dx is finite. - lattice_incircle(a, b, c, d, f): four node corners, every difference from d at most 2^14 nodes, then the int64 incircle determinant on (col, -row) translated to d (below 3 * 2^58). - must_flip calls it first; when it answers, the answer is used and the kernel is not asked, not even the frame orient2d. refine builds its frame with lattice_frame(g.delta_x(), g.delta_y(), g.rows(), g.cols()). No test file and not the CMake guard are touched. Bit-identical: T6 and test_refine_golden.py pass unchanged, and the 1 m benchmark's binary meshes hash as 21a's. Timed back to back against d3ee2ce on battery: split phase -26 to -27 %, refine -7 to -8 % at 1 thread and -17 % at 8 (recorded in the 21 file; @perf's AC acceptance is still owed). 54 production lines added, 3 removed (CLAUDE.md section 2's unit), against an estimate of ~50. Co-Authored-By: Claude Opus 5.5 --- docs/increments/21-parallel-refine.md | 24 +++++++++ include/terrain/mesh/lawson.hpp | 78 ++++++++++++++++++++++++++- include/terrain/refinement/refine.hpp | 2 +- 3 files changed, 101 insertions(+), 3 deletions(-) diff --git a/docs/increments/21-parallel-refine.md b/docs/increments/21-parallel-refine.md index e1cc6e60..5caf8952 100644 --- a/docs/increments/21-parallel-refine.md +++ b/docs/increments/21-parallel-refine.md @@ -541,6 +541,30 @@ enabling the path; no call in `must_flip`. The determinant computed in doubles rather than `int64` was also killed, and only by the radius-8085 circle, where spreads come close to 2^14. +### 21b: the integer path is faster (green, not acceptance) + +Measured by `@developer` at the green commit, back to back against 21a's head +`d3ee2ce`, **on battery** (98 %, before and after), so it compares with no AC +figure. Both trees were built by `bench.py`'s `build()` (Release; the base +`.so` hashes to the 21a acceptance run's `654a2442…`), and each sample was one +`rasputin mesh` process on the 1 m benchmark at tolerance 1, reading +`RefineOutcome`'s phase times. Medians of 5, interleaved base/green: + +| domain, threads | refine, 21a -> 21b | split phase, 21a -> 21b | +|---|---:|---:| +| quarter, 1 | 0.444 -> 0.414 s (-6.7 %) | 0.122 -> 0.090 s (-26 %) | +| quarter, 8 | 0.198 -> 0.164 s (-17 %) | 0.132 -> 0.098 s (-26 %) | +| tile, 1 | 0.503 -> 0.464 s (-7.8 %) | 0.124 -> 0.091 s (-27 %) | +| tile, 8 | 0.225 -> 0.187 s (-17 %) | 0.134 -> 0.097 s (-27 %) | + +The scan is unchanged, as it should be; the whole gain is in the split phase, +and it is larger than QW2's estimate (-4 % and -8 %), so the int64 determinant +with the frame `orient2d` it skips is cheaper than the filtered path too, not +only than the exact path (5.0 % of refine in the profile) it replaces. +The binary meshes hash the same as 21a's (quarter `ff705683…`, tile +`645919aa…`), with the same flip and insertion counts. The scratch driver and +its samples were not kept; `@perf`'s AC acceptance run is still owed. + ## 4. Determinism levels Today's contract (14 R5, 14b R1, tested by 14's T6 and 18's T3 golden diff --git a/include/terrain/mesh/lawson.hpp b/include/terrain/mesh/lawson.hpp index 27960558..2a2f8922 100644 --- a/include/terrain/mesh/lawson.hpp +++ b/include/terrain/mesh/lawson.hpp @@ -24,23 +24,92 @@ #include #include +#include +#include +#include +#include #include #include +#include #include #include #include namespace terrain::mesh { -struct LatticeFrame { +// A frame built directly, LatticeFrame{dx, dy}, never enables the integer +// incircle; only lattice_frame() below can (docs/increments/21-parallel-refine.md, +// "Pinned by the red suite (21b)"). +class LatticeFrame { +public: double dx = 1.0; double dy = 1.0; + constexpr LatticeFrame() noexcept = default; + constexpr LatticeFrame(double x, double y) noexcept : dx{x}, dy{y} {} + [[nodiscard]] Point2 at(MeshVertex v) const noexcept { return Point2{v.col * dx, -(v.row * dy)}; } + // Whether lattice_incircle may answer on this frame. + [[nodiscard]] bool integer() const noexcept { return integer_; } + +private: + bool integer_ = false; + friend LatticeFrame lattice_frame(double, double, std::size_t, std::size_t) noexcept; }; +// The frame of a rows x cols grid, with the integer path enabled when the frame +// is the lattice times one exact constant (QW2): dx == dy, dx a positive normal +// double, and col * dx, row * dy exact for every node. Exactness is decided by +// QW2's sufficient condition, which is also what it takes: the significant bits +// of dx plus bit_width(max(rows, cols) - 1) are at most 53, so every product +// fits the mantissa; the largest one must also be finite. +[[nodiscard]] inline LatticeFrame lattice_frame(double dx, double dy, std::size_t rows, + std::size_t cols) noexcept { + LatticeFrame f{dx, dy}; + if (dx != dy || !std::isnormal(dx) || dx < 0.0) + return f; + const std::uint64_t mantissa = + (std::bit_cast(dx) & ((std::uint64_t{1} << 52) - 1)) | (std::uint64_t{1} << 52); + const auto significant = 53 - std::countr_zero(mantissa); + const std::size_t top = std::max(rows, cols); + const auto span = top == 0 ? 0 : std::bit_width(top - 1); + f.integer_ = significant + span <= 53 + && std::isfinite(static_cast(top == 0 ? 0 : top - 1) * dx); + return f; +} + +// The incircle sign of the quad (a, b, c, d) from an int64 determinant on the +// lattice (col, -row), translated to d. Answers only on an enabling frame, four +// node corners and every difference from d at most 2^14 nodes; the determinant +// is then below 3 * 2^58 and exact. Under those conditions it is the sign +// DetriaExact gives on the frame points, since the frame is the lattice scaled +// by dx > 0 without rounding. Precondition: a, b, c counter-clockwise on +// (col, -row). Otherwise empty, and the caller asks the kernel. +[[nodiscard]] inline std::optional lattice_incircle(MeshVertex a, MeshVertex b, + MeshVertex c, MeshVertex d, + const LatticeFrame& f) noexcept { + if (!f.integer() || !a.is_node() || !b.is_node() || !c.is_node() || !d.is_node()) + return std::nullopt; + constexpr double bound = 1 << 14; + std::array, 3> p{}; // (col, -row) minus d's + const std::array abc{a, b, c}; + for (std::size_t i = 0; i < 3; ++i) { + const double x = abc[i].col - d.col, y = d.row - abc[i].row; // exact: integers < 2^53 + if (std::fabs(x) > bound || std::fabs(y) > bound) + return std::nullopt; + p[i] = {static_cast(x), static_cast(y)}; + } + const auto lift = [](const std::array& q) { return q[0] * q[0] + q[1] * q[1]; }; + const auto cross = [](const std::array& u, const std::array& v) { + return u[0] * v[1] - u[1] * v[0]; + }; + const std::int64_t det = + lift(p[0]) * cross(p[1], p[2]) + lift(p[1]) * cross(p[2], p[0]) + lift(p[2]) * cross(p[0], p[1]); + return det > 0 ? pred::Incircle::Inside : det < 0 ? pred::Incircle::Outside : pred::Incircle::Cocircular; +} + namespace detail { // Whether t's edge e must flip: interior, unconstrained, and the quad's @@ -57,8 +126,13 @@ template while (m.triangles()[u][j] != tri[(e + 1) % 3]) ++j; const auto v = m.vertices(); + const MeshVertex d_vertex = v[m.triangles()[u][(j + 2) % 3]]; + // The integer path first; the mesh triangle is counter-clockwise on + // (col, -row), and it answers only where the kernel would give the same sign. + if (const auto s = lattice_incircle(v[tri[e]], v[tri[(e + 1) % 3]], v[tri[(e + 2) % 3]], d_vertex, f)) + return *s == pred::Incircle::Inside; const Point2 a = f.at(v[tri[e]]), b = f.at(v[tri[(e + 1) % 3]]), c = f.at(v[tri[(e + 2) % 3]]), - d = f.at(v[m.triangles()[u][(j + 2) % 3]]); + d = f.at(d_vertex); // Orientation is exact on (col, -row), but the frame rounds col * dx and // row * dy, so a triangle counter-clockwise in the mesh can be collinear // or clockwise here. The kernel answers Cocircular for a collinear triple diff --git a/include/terrain/refinement/refine.hpp b/include/terrain/refinement/refine.hpp index c0c6a256..1a6a7823 100644 --- a/include/terrain/refinement/refine.hpp +++ b/include/terrain/refinement/refine.hpp @@ -276,7 +276,7 @@ template return std::move(*refused); auto& m = std::get(built); - const mesh::LatticeFrame frame{g.delta_x(), g.delta_y()}; + const auto frame = mesh::lattice_frame(g.delta_x(), g.delta_y(), g.rows(), g.cols()); RefineOutcome out; using clock = std::chrono::steady_clock; const auto since = [](clock::time_point t0) { From 7eb1220807b6add162b7157f8e7a553d449eceb2 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 18:31:19 +0200 Subject: [PATCH 05/12] 21b: remove the red-step scaffold, register the suite plainly Co-Authored-By: Claude Opus 5.5 --- docs/increments/21-parallel-refine.md | 5 +++-- tests/cpp/CMakeLists.txt | 15 ++------------- 2 files changed, 5 insertions(+), 15 deletions(-) diff --git a/docs/increments/21-parallel-refine.md b/docs/increments/21-parallel-refine.md index 5caf8952..db73d6e7 100644 --- a/docs/increments/21-parallel-refine.md +++ b/docs/increments/21-parallel-refine.md @@ -527,8 +527,9 @@ header beside it that `lawson.hpp` includes is `@developer`'s choice. `orient2d`. The suite counts kernel calls to check this. When it refuses, today's path runs unchanged. -The suite is registered only once `lawson.hpp` names `lattice_incircle`, as -in 21a. It starts no threads, so it is not in the TSan job. +At the red commit the suite was registered only once `lawson.hpp` named +`lattice_incircle`, as in 21a; the guard was removed after green. It starts no +threads, so it is not in the TSan job. **Mutation round** against a scratch implementation that is not committed. All of these were killed: the sign inverted; the determinant on `(col, row)`; diff --git a/tests/cpp/CMakeLists.txt b/tests/cpp/CMakeLists.txt index 42d6d7b3..f179e0a1 100644 --- a/tests/cpp/CMakeLists.txt +++ b/tests/cpp/CMakeLists.txt @@ -299,16 +299,5 @@ target_link_libraries(test_mesh_lawson_stack PRIVATE Threads::Threads) # DetriaExact on the frame doubles wherever lattice_incircle answers, refusal # wherever QW2's conditions fail, and must_flip's decisions as today's. It # starts no threads, so it is not in the TSan job. -# -# RED UNTIL THE INTERFACE EXISTS: registered only once lawson.hpp names -# lattice_incircle (defined there or in a header it includes), so the red does -# not break the rest of the build. The header is a configure dependency, so the -# green commit's edit registers the suite without a manual cmake. The review -# drops this guard. -set_property(DIRECTORY APPEND PROPERTY CMAKE_CONFIGURE_DEPENDS - "${PROJECT_SOURCE_DIR}/include/terrain/mesh/lawson.hpp") -file(READ "${PROJECT_SOURCE_DIR}/include/terrain/mesh/lawson.hpp" _lawson_21b) -string(FIND "${_lawson_21b}" "lattice_incircle" _lattice_incircle_at) -if(NOT _lattice_incircle_at EQUAL -1) - add_terrain_backend_test(test_mesh_lattice_incircle unit/test_mesh_lattice_incircle.cpp) -endif() + +add_terrain_backend_test(test_mesh_lattice_incircle unit/test_mesh_lattice_incircle.cpp) From 60fc65fc497b4296fac4556280f5d1856a42b5f5 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 18:47:54 +0200 Subject: [PATCH 06/12] 21b review fixes: stale 'every caller' and 'not measured', a rounding, the 10x10 tie Co-Authored-By: Claude Opus 5.5 --- docs/benchmarks/2026-09-27/21b-ties/README.md | 3 ++- docs/increments/21-parallel-refine.md | 7 ++++--- 2 files changed, 6 insertions(+), 4 deletions(-) diff --git a/docs/benchmarks/2026-09-27/21b-ties/README.md b/docs/benchmarks/2026-09-27/21b-ties/README.md index a2cf951e..7b27c1da 100644 --- a/docs/benchmarks/2026-09-27/21b-ties/README.md +++ b/docs/benchmarks/2026-09-27/21b-ties/README.md @@ -140,7 +140,8 @@ the other phases, which also go through `must_flip`: The remaining 411 and 3,055 calls have an off-node corner (the arc vertices). - On the tile, QW2 would answer all of `legalise_all`'s 48,133 calls, 16,129 - of them exact ties: 15,876 are 40×40 squares and 252 are 10×40 rectangles, + of them exact ties: 15,876 are 40×40 squares, 252 are 10×40 rectangles and + one is a 10×10 square, from the starting mesh. It would also answer all 4,017 of the quality pass's calls. diff --git a/docs/increments/21-parallel-refine.md b/docs/increments/21-parallel-refine.md index db73d6e7..4dd1167d 100644 --- a/docs/increments/21-parallel-refine.md +++ b/docs/increments/21-parallel-refine.md @@ -510,7 +510,7 @@ header beside it that `lawson.hpp` includes is `@developer`'s choice. output is bit-identical by design. The reviewer checks it by reading the code, and `@perf`'s acceptance run shows it as time saved. - **A `LatticeFrame` built directly never enables the integer path.** - `LatticeFrame{dx, dy}`, the way every caller builds one today, keeps compiling + `LatticeFrame{dx, dy}`, the way every caller except refine builds one, keeps compiling and keeps today's kernel path. The existing Lawson and quality suites therefore still exercise that path, and the 21b suite uses such a frame as its "without". @@ -556,7 +556,7 @@ figure. Both trees were built by `bench.py`'s `build()` (Release; the base | quarter, 1 | 0.444 -> 0.414 s (-6.7 %) | 0.122 -> 0.090 s (-26 %) | | quarter, 8 | 0.198 -> 0.164 s (-17 %) | 0.132 -> 0.098 s (-26 %) | | tile, 1 | 0.503 -> 0.464 s (-7.8 %) | 0.124 -> 0.091 s (-27 %) | -| tile, 8 | 0.225 -> 0.187 s (-17 %) | 0.134 -> 0.097 s (-27 %) | +| tile, 8 | 0.225 -> 0.187 s (-17 %) | 0.134 -> 0.097 s (-28 %) | The scan is unchanged, as it should be; the whole gain is in the split phase, and it is larger than QW2's estimate (-4 % and -8 %), so the int64 determinant @@ -869,7 +869,8 @@ Each is answered or placed. (146,962 of 146,962 on the quarter circle, 154,502 of 154,502 on the tile); QW2 would answer 99.97 % and 100 % of all refine-loop incircle calls. Only 38 % of the ties are axis-aligned rectangles, so the general determinant is - needed. The time saved is not measured. + needed. The time saved was measured on battery at green ("21b: the + integer path is faster"); `@perf`'s acceptance run is still owed. 4. *"The scan loses 16.5 ms of 61.8 ms to imbalance at 8 threads. Is the chunking free to change, given that results are written per slot? And does From 24c280197bed67542b324493a4a85b8023034f1d Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 18:51:48 +0200 Subject: [PATCH 07/12] 21b tests: both sides of the exact-frame bound, overflowing extent, non-finite corners Test amendment after green, in its own commit as docs/increments/README.md allows; no production code. @reviewer found two surviving mutants in lattice_frame and one unguarded input in lattice_incircle: - The limit taken one bit loose (bit_width(top - 1) - 1) survived: the 0.1 refusals sat only at 7-bit spans. Now pinned on both sides: 0.1 and 1.5 + 2^-51 (52 bits, 3 * dx rounds) must refuse on 4x4, 4x2 and 2x4, and must answer on 2x2; 1.5 + 2^-50 (51 bits) must answer on 4x4. A sweep of every ordered quad of the 4x4 grid requires any answer to match DetriaExact on the 0.1 frame, and the reviewer's quad (a lattice tie the 0.1 frame calls Inside) is pinned by name. - Dropping the isfinite extent clause survived. 2^1023 on a 3x3 grid (1 + 2 bits, extent inf) must refuse; on 2x2 it must answer. - is_node() is true for +-inf, and inf - inf = NaN passes `fabs > bound` into the int64 cast. lattice_incircle must refuse non-finite corners. This case is RED at this commit: with a, b, c and d on the same infinity it returns a value (4 failing assertions); the developer adds the check next. Also fixes the stale "as every caller builds it today" comment, and records the finite-corner and extent pins and the two mutants in the increment file. Co-Authored-By: Claude Opus 5.5 --- docs/increments/21-parallel-refine.md | 15 +- tests/cpp/unit/test_mesh_lattice_incircle.cpp | 143 +++++++++++++++++- 2 files changed, 153 insertions(+), 5 deletions(-) diff --git a/docs/increments/21-parallel-refine.md b/docs/increments/21-parallel-refine.md index 4dd1167d..a4fa0e8c 100644 --- a/docs/increments/21-parallel-refine.md +++ b/docs/increments/21-parallel-refine.md @@ -502,7 +502,10 @@ header beside it that `lawson.hpp` includes is `@developer`'s choice. when `dx == dy`, `dx` is finite and positive, and `col * dx` and `row * dy` are exact for every node with `col < cols` and `row < rows`. It must answer wherever QW2's sufficient condition holds: the significant bits of `dx` plus - `bit_width(max(rows, cols) - 1)` are at most 53. The suite tests exactly 53. + `bit_width(max(rows, cols) - 1)` are at most 53. The suite tests exactly 53, + and 54 with a `dx` whose `3 * dx` rounds (it must refuse), on each axis. A + largest product `(max(rows, cols) - 1) * dx` that overflows is not exact: + `2^1023` on a 3 × 3 grid must refuse. Between that condition and exactness, for example `dx = 0.1` on a 3 × 3 grid, which is exact but fails the condition, nothing is pinned. `refine` builds its frame with this function instead of @@ -517,8 +520,9 @@ header beside it that `lawson.hpp` includes is `@developer`'s choice. - `[[nodiscard]] std::optional lattice_incircle(MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d, const LatticeFrame& f) noexcept`. Precondition: `a, b, c` strictly counter-clockwise on `(col, -row)`, which `must_flip`'s triangle is (the `LatticeMesh` invariant). It returns empty - unless `f` came from an enabling `lattice_frame`, all four corners are nodes, - and `|col_x - col_d|` and `|row_x - row_d|` are at most 2^14 for each `x` in + unless `f` came from an enabling `lattice_frame`, all four corners are nodes + with finite coordinates (`is_node()` is true for ±inf, and `inf - inf` is a + NaN that the spread bound lets through to the `int64` cast), and `|col_x - col_d|` and `|row_x - row_d|` are at most 2^14 for each `x` in `a, b, c`. The bound is measured from `d`, so `a` and `b` may be 2^15 apart. When it answers, the answer is the sign `DetriaExact::incircle_ccw` gives on the frame points `(col * dx, -(row * dy))`. @@ -542,6 +546,11 @@ enabling the path; no call in `must_flip`. The determinant computed in doubles rather than `int64` was also killed, and only by the radius-8085 circle, where spreads come close to 2^14. +Two mutants survived that round and were found by `@reviewer` after green: the +exact-frame limit one bit loose (`bit_width(top - 1) - 1`) and the finite-extent +clause dropped. A test amendment after green kills both, and adds the +non-finite-corner refusal, which HEAD at the amendment did not have. + ### 21b: the integer path is faster (green, not acceptance) Measured by `@developer` at the green commit, back to back against 21a's head diff --git a/tests/cpp/unit/test_mesh_lattice_incircle.cpp b/tests/cpp/unit/test_mesh_lattice_incircle.cpp index bf22f061..a316bf5f 100644 --- a/tests/cpp/unit/test_mesh_lattice_incircle.cpp +++ b/tests/cpp/unit/test_mesh_lattice_incircle.cpp @@ -21,8 +21,9 @@ // MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d, const LatticeFrame& f) noexcept; // } // -// A LatticeFrame built directly, LatticeFrame{dx, dy} as every caller builds -// it today, never enables the integer path. +// A LatticeFrame built directly, LatticeFrame{dx, dy}, never enables the +// integer path; refine builds its frame with lattice_frame, every other caller +// directly. // // The oracle is DetriaExact::incircle_ccw on the frame points, computed here // from col * dx and -(row * dy), not through LatticeFrame::at. For today's @@ -537,6 +538,144 @@ TEST_CASE("lattice_frame enables the integer path wherever QW2's sufficient cond REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(0.1, 0.1, 2, 2)) == Incircle::Cocircular); } +TEST_CASE("lattice_frame refuses one bit past QW2's limit and answers at it, on both axes", + "[lattice_incircle][frame][limit]") { + // Both sides of "significant bits of dx + bit_width(max(rows, cols) - 1) + // <= 53". One bit past it the frame can be inexact, and is for these dx: + // a grid with col or row 3 and a 52-bit dx whose 3 * dx needs 54 bits. + // A limit taken one bit loose (e.g. bit_width(top - 1) - 1) accepts them. + const double dx51 = 1.5 + std::ldexp(1.0, -50); // 51 significant bits + // 0.1 and 1.5 + 2^-51 both have 52 significant bits. + const double dx = GENERATE(0.1, 1.5 + std::ldexp(1.0, -51)); + CAPTURE(dx); + REQUIRE(std::fma(3.0, dx, -(3.0 * dx)) != 0.0); // col or row 3 rounds + // 52 + bit_width(3) = 54: refused, whichever axis carries the 4. + REQUIRE_FALSE(lattice_frame(dx, dx, 4, 4).integer()); + REQUIRE_FALSE(lattice_frame(dx, dx, 4, 2).integer()); + REQUIRE_FALSE(lattice_frame(dx, dx, 2, 4).integer()); + // 52 + bit_width(1) = 53: at the limit, it must answer. + REQUIRE(lattice_frame(dx, dx, 2, 2).integer()); + // 51 + bit_width(3) = 53 on the 4 x 4 grid itself: must answer. + REQUIRE(std::fma(3.0, dx51, -(3.0 * dx51)) == 0.0); + REQUIRE(lattice_frame(dx51, dx51, 4, 4).integer()); + REQUIRE(lattice_frame(dx51, dx51, 4, 2).integer()); + REQUIRE(lattice_frame(dx51, dx51, 2, 4).integer()); +} + +TEST_CASE("lattice_incircle never disagrees with DetriaExact on a frame one bit past the limit", + "[lattice_incircle][frame][limit]") { + // Every ordered quad of the 4 x 4 grid, a, b, c counter-clockwise, under + // lattice_frame(dx, dx, 4, 4) with dx one bit past the limit. Wherever it + // answers, it must give the oracle's sign; and the grid holds quads whose + // lattice sign differs from the frame's, so an answer there is wrong. + const double dx = GENERATE(0.1, 1.5 + std::ldexp(1.0, -51)); + CAPTURE(dx); + const LatticeFrame f = lattice_frame(dx, dx, 4, 4); + const LatticeFrame exact = lattice_frame(1.0, 1.0, 4, 4); // the lattice sign + std::vector nodes; + for (std::int64_t r = 0; r < 4; ++r) + for (std::int64_t c = 0; c < 4; ++c) + nodes.push_back(node(c, r)); + std::size_t quads = 0, sign_differs = 0, wrong = 0; + std::string first; + for (const auto& a : nodes) + for (const auto& b : nodes) + for (const auto& c : nodes) { + if (iorient(a, b, c) <= 0) + continue; + for (const auto& d : nodes) { + if (d == a || d == b || d == c) + continue; + ++quads; + const auto want = oracle(dx, dx, a, b, c, d); + const auto lattice = lattice_incircle(a, b, c, d, exact); + REQUIRE(lattice.has_value()); + sign_differs += want != lattice; + const auto got = lattice_incircle(a, b, c, d, f); + if (got && got != want && wrong++ == 0) + first = describe(a) + " " + describe(b) + " " + describe(c) + " d " + describe(d); + } + } + CAPTURE(quads, sign_differs, wrong, first); + REQUIRE(sign_differs > 0); + REQUIRE(wrong == 0); +} + +TEST_CASE("lattice_incircle refuses the reviewer's quad on 0.1 at 4 x 4, a lattice tie the frame calls Inside", + "[lattice_incircle][frame][limit]") { + // In (col, -row): a (0, 0), b (0, -1), c (3, 0), d (1, -2). The circle + // through a, b, c has centre (1.5, -0.5) and squared radius 2.5, and d is + // on it; at dx = 0.1 the rounded frame puts d strictly inside. + const MeshVertex a = node(0, 0), b = node(0, 1), c = node(3, 0), d = node(1, 2); + REQUIRE(iorient(a, b, c) > 0); + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(1.0, 1.0, 4, 4)) == Incircle::Cocircular); + REQUIRE(oracle(0.1, 0.1, a, b, c, d) == Incircle::Inside); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(0.1, 0.1, 4, 4)).has_value()); +} + +TEST_CASE("lattice_frame refuses a frame whose extent overflows, and answers where it does not", + "[lattice_incircle][frame][extent]") { + // dx = 2^1023 has one significant bit, so the bit count passes on any + // small grid; but (top - 1) * dx is infinite once top - 1 >= 2, and an + // infinite frame coordinate is not exact. + const double huge = std::ldexp(1.0, 1023); + REQUIRE(std::isinf(2.0 * huge)); + REQUIRE(std::isfinite(1.0 * huge)); + REQUIRE_FALSE(lattice_frame(huge, huge, 3, 3).integer()); + REQUIRE_FALSE(lattice_frame(huge, huge, 3, 2).integer()); + REQUIRE_FALSE(lattice_frame(huge, huge, 2, 3).integer()); + const MeshVertex a = node(0, 1), b = node(1, 1), c = node(0, 0), d = node(1, 0); // unit square + REQUIRE(iorient(a, b, c) > 0); + REQUIRE_FALSE(lattice_incircle(a, b, c, d, lattice_frame(huge, huge, 3, 3)).has_value()); + // On a 2 x 2 grid every product is 0 or 2^1023, finite and exact. + REQUIRE(lattice_frame(huge, huge, 2, 2).integer()); + REQUIRE(lattice_incircle(a, b, c, d, lattice_frame(huge, huge, 2, 2)) == Incircle::Cocircular); +} + +TEST_CASE("lattice_incircle refuses a corner with a non-finite coordinate", "[lattice_incircle][refuse][nonfinite]") { + // MeshVertex::is_node is true for +-inf (floor(inf) == inf), so the node + // check does not stop them. A corner infinite alone gives an infinite + // difference from d, which the spread bound refuses; but a corner and d + // on the same infinity give inf - inf = NaN, which passes `fabs > bound` + // and reaches the int64 cast, undefined for NaN. When a, b, c and d all + // share it, no difference is refused and a value comes back. NaN is not a + // node (NaN != floor(NaN)). lattice_incircle is public, so it must refuse + // all of these. Unreachable from refine's grid mesh. The precondition + // (a, b, c counter-clockwise) has no meaning here; the refusal comes first. + const LatticeFrame f = lattice_frame(10.0, 10.0, kBench, kBench); + const std::array base{node(10, 12), node(12, 12), node(11, 10), node(11, 11)}; + REQUIRE(iorient(base[0], base[1], base[2]) > 0); + REQUIRE(lattice_incircle(base[0], base[1], base[2], base[3], f).has_value()); + const double bad = GENERATE(std::numeric_limits::infinity(), -std::numeric_limits::infinity(), + std::numeric_limits::quiet_NaN()); + const bool on_col = GENERATE(true, false); + const auto spoil = [&](MeshVertex v) { return on_col ? MeshVertex{bad, v.row} : MeshVertex{v.col, bad}; }; + SECTION("one corner non-finite, in each position") { + for (std::size_t which = 0; which < 4; ++which) { + auto q = base; + q[which] = spoil(q[which]); + CAPTURE(which, bad, on_col); + REQUIRE_FALSE(lattice_incircle(q[0], q[1], q[2], q[3], f).has_value()); + } + } + SECTION("a corner and d on the same non-finite value; the other two corners infinitely far") { + for (std::size_t which = 0; which < 3; ++which) { + auto q = base; + q[which] = spoil(q[which]); + q[3] = spoil(q[3]); + CAPTURE(which, bad, on_col); + REQUIRE_FALSE(lattice_incircle(q[0], q[1], q[2], q[3], f).has_value()); + } + } + SECTION("every corner on the same non-finite value, so every difference from d is NaN") { + auto q = base; + for (auto& v : q) + v = spoil(v); + CAPTURE(bad, on_col); + REQUIRE_FALSE(lattice_incircle(q[0], q[1], q[2], q[3], f).has_value()); + } +} + TEST_CASE("an inexact frame changes the answer on some lattice tie, so its refusal is load-bearing", "[lattice_incircle][refuse]") { // QW2's determinism argument: without the exact-frame condition the From a03ae65f06773d7191f7bcb66088a8ce087892d7 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 18:54:47 +0200 Subject: [PATCH 08/12] 21b: lattice_incircle refuses non-finite corners; grid and exact-kernel preconditions stated is_node() is true for +-inf, and a corner and d on the same infinity give inf - inf = NaN, which passed the spread bound and reached the int64 cast (UB; UBSan fired on [nonfinite]). Every corner is now checked finite before any difference is taken. Comments only, no behaviour change: the DetriaExact-sign claim holds for nodes inside the rows x cols grid the frame was built for; must_flip answering without K assumes K's incircle sign is exact. Co-Authored-By: Claude Opus 5.5 --- include/terrain/mesh/lawson.hpp | 18 +++++++++++++++--- 1 file changed, 15 insertions(+), 3 deletions(-) diff --git a/include/terrain/mesh/lawson.hpp b/include/terrain/mesh/lawson.hpp index 2a2f8922..a4ffca55 100644 --- a/include/terrain/mesh/lawson.hpp +++ b/include/terrain/mesh/lawson.hpp @@ -85,12 +85,21 @@ class LatticeFrame { // node corners and every difference from d at most 2^14 nodes; the determinant // is then below 3 * 2^58 and exact. Under those conditions it is the sign // DetriaExact gives on the frame points, since the frame is the lattice scaled -// by dx > 0 without rounding. Precondition: a, b, c counter-clockwise on -// (col, -row). Otherwise empty, and the caller asks the kernel. +// by dx > 0 without rounding. That holds only for nodes inside the rows x cols +// grid f was built for: lattice_frame checked exactness of col * dx and +// row * dy up to that extent and no further, so a node beyond it can have a +// rounded frame point whose DetriaExact sign differs. Precondition: a, b, c +// counter-clockwise on (col, -row). Otherwise empty, and the caller asks the +// kernel. A corner with a non-finite coordinate is refused before any +// difference is taken: is_node() is true for +-inf, and inf - inf is a NaN that +// would pass the spread bound and reach the int64 cast. [[nodiscard]] inline std::optional lattice_incircle(MeshVertex a, MeshVertex b, MeshVertex c, MeshVertex d, const LatticeFrame& f) noexcept { - if (!f.integer() || !a.is_node() || !b.is_node() || !c.is_node() || !d.is_node()) + const auto usable = [](MeshVertex v) { + return std::isfinite(v.col) && std::isfinite(v.row) && v.is_node(); + }; + if (!f.integer() || !usable(a) || !usable(b) || !usable(c) || !usable(d)) return std::nullopt; constexpr double bound = 1 << 14; std::array, 3> p{}; // (col, -row) minus d's @@ -129,6 +138,9 @@ template const MeshVertex d_vertex = v[m.triangles()[u][(j + 2) % 3]]; // The integer path first; the mesh triangle is counter-clockwise on // (col, -row), and it answers only where the kernel would give the same sign. + // Answering without consulting K assumes K's incircle sign is exact, as + // DefaultKernel's (FilteredKernel) is: for an inexact K the + // integer sign could differ from K's, and the flip sequence with it. if (const auto s = lattice_incircle(v[tri[e]], v[tri[(e + 1) % 3]], v[tri[(e + 2) % 3]], d_vertex, f)) return *s == pred::Incircle::Inside; const Point2 a = f.at(v[tri[e]]), b = f.at(v[tri[(e + 1) % 3]]), c = f.at(v[tri[(e + 2) % 3]]), From 560a99d9f74286c975b4d15d32b1a8e46e1b73a3 Mon Sep 17 00:00:00 2001 From: Ola Skavhaug Date: Sun, 27 Sep 2026 19:10:10 +0200 Subject: [PATCH 09/12] 21b acceptance: ACCEPTED on battery, refine -15 % at 8 threads, mesh identical Back to back against master 071e4f3 (21a merge), two pairs, bench.py blob 77765b1, threads 1..20, 5 repeats, quarter and tile, 1 m DEM, tolerance 1. Split phase -21..-24 %, scan unchanged; ceiling 2.45-2.55x. Head's split phase is 4-6 % slower than green f707322 (non-finite refusal in lattice_incircle, attributed by elimination). Co-Authored-By: Claude Opus 5.5 --- docs/benchmarks/2026-09-27/21b-acceptance.md | 205 ++ .../2026-09-27/21b-base-071e4f3-r2/README.md | 80 + .../2026-09-27/21b-base-071e4f3-r2/raw.tsv | 210 ++ .../2026-09-27/21b-base-071e4f3-r2/run.json | 2956 +++++++++++++++++ .../2026-09-27/21b-base-071e4f3/README.md | 79 + .../2026-09-27/21b-base-071e4f3/raw.tsv | 210 ++ .../2026-09-27/21b-base-071e4f3/run.json | 2955 ++++++++++++++++ docs/benchmarks/2026-09-27/21b-r2/README.md | 79 + docs/benchmarks/2026-09-27/21b-r2/raw.tsv | 210 ++ docs/benchmarks/2026-09-27/21b-r2/run.json | 2955 ++++++++++++++++ docs/benchmarks/2026-09-27/21b/README.md | 79 + docs/benchmarks/2026-09-27/21b/raw.tsv | 210 ++ docs/benchmarks/2026-09-27/21b/run.json | 2955 ++++++++++++++++ 13 files changed, 13183 insertions(+) create mode 100644 docs/benchmarks/2026-09-27/21b-acceptance.md create mode 100644 docs/benchmarks/2026-09-27/21b-base-071e4f3-r2/README.md create mode 100644 docs/benchmarks/2026-09-27/21b-base-071e4f3-r2/raw.tsv create mode 100644 docs/benchmarks/2026-09-27/21b-base-071e4f3-r2/run.json create mode 100644 docs/benchmarks/2026-09-27/21b-base-071e4f3/README.md create mode 100644 docs/benchmarks/2026-09-27/21b-base-071e4f3/raw.tsv create mode 100644 docs/benchmarks/2026-09-27/21b-base-071e4f3/run.json create mode 100644 docs/benchmarks/2026-09-27/21b-r2/README.md create mode 100644 docs/benchmarks/2026-09-27/21b-r2/raw.tsv create mode 100644 docs/benchmarks/2026-09-27/21b-r2/run.json create mode 100644 docs/benchmarks/2026-09-27/21b/README.md create mode 100644 docs/benchmarks/2026-09-27/21b/raw.tsv create mode 100644 docs/benchmarks/2026-09-27/21b/run.json diff --git a/docs/benchmarks/2026-09-27/21b-acceptance.md b/docs/benchmarks/2026-09-27/21b-acceptance.md new file mode 100644 index 00000000..0295a991 --- /dev/null +++ b/docs/benchmarks/2026-09-27/21b-acceptance.md @@ -0,0 +1,205 @@ +# Increment 21b acceptance run (@perf, 2026-09-27) + +**Verdict: ACCEPTED.** On battery, measured back to back against the previous +merge (`071e4f3`, 21a): refine at 8 threads is 14.5-14.8 % faster on both +domains, in both pairs. The mesh sha256 is identical on both domains. No +measure regressed. + +One thing to act on: the head `a03ae65` has a split phase 4-6 % slower than +21b's green commit `f707322`. The only production change between them is the +non-finite refusal in `lattice_incircle` (see "Green against head"). + +## Method + +- **Why a re-measured base.** The Mac was on battery. The stored `bench.py` + runs of 21a (`21a*/`) are AC. The only battery `bench.py` run, + `serial-profile-bench/`, predates 21a. So there was no comparable baseline. + Master `071e4f3` (the 21a merge, #102) was checked out in a temporary + `git worktree` and measured with `--tree`, back to back with 21b's head + `a03ae65`. The worktree was removed afterwards. +- **`--baseline` is explicit.** `071e4f3` is not an ancestor of `a03ae65`: the + branch forked from 21a's head, not from the merge + (`git merge-base --is-ancestor 071e4f3 a03ae65` exits 1). `find_baseline` + would therefore skip it. +- **Two pairs**, run in the order base, 21b, base, 21b between 18:56 and + 19:05 CEST. The verdicts can be reproduced from the stored evidence: + + ``` + .venv/bin/python tools/bench.py compare docs/benchmarks/2026-09-27/21b \ + --baseline docs/benchmarks/2026-09-27/21b-base-071e4f3 # ACCEPTED + .venv/bin/python tools/bench.py compare docs/benchmarks/2026-09-27/21b-r2 \ + --baseline docs/benchmarks/2026-09-27/21b-base-071e4f3-r2 # ACCEPTED + .venv/bin/python tools/bench.py compare docs/benchmarks/2026-09-27/21b-base-071e4f3-r2 \ + --baseline docs/benchmarks/2026-09-27/21b-base-071e4f3 # REGRESSION (same build; see Spread) + .venv/bin/python tools/bench.py compare docs/benchmarks/2026-09-27/21b-r2 \ + --baseline docs/benchmarks/2026-09-27/21b # ACCEPTED + ``` + + The two cross pairs (`21b-r2` against `21b-base-071e4f3`, and `21b` against + `21b-base-071e4f3-r2`) are also ACCEPTED. +- **Command.** All four runs used the main tree's `tools/bench.py`, blob + `77765b1`, the same as 21a's. Each run did its own Release build into + `/build-bench`. + + ``` + .venv/bin/python tools/bench.py run --label