Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 4 additions & 1 deletion .github/workflows/benchmarks.yml
Original file line number Diff line number Diff line change
Expand Up @@ -114,8 +114,11 @@ jobs:

if [ -d target/criterion/exact_d2/det/main ]; then
echo "::notice::Baseline found — comparing against main"
# --baseline-lenient rather than --baseline: benches added on the
# PR branch that don't yet exist in the main baseline get a
# "no baseline data" notice instead of aborting the whole run.
cargo bench --features bench,exact --bench exact \
-- --baseline main 2>&1 | tee bench-output.txt
-- --baseline-lenient main 2>&1 | tee bench-output.txt
else
echo "::notice::No baseline found — running without comparison"
cargo bench --features bench,exact --bench exact \
Expand Down
8 changes: 6 additions & 2 deletions AGENTS.md
Original file line number Diff line number Diff line change
Expand Up @@ -221,7 +221,10 @@ When creating or updating issues:

- `exact` — enables exact arithmetic methods via `BigRational`:
`det_exact()`, `det_exact_f64()`, `det_sign_exact()`, `solve_exact()`, and `solve_exact_f64()`.
Also re-exports `BigRational` from the crate root and prelude.
Re-exports `BigInt`, `BigRational`, and the commonly needed `num-traits`
items (`FromPrimitive`, `ToPrimitive`, and `Signed`) from the crate root and prelude
(so consumers get usable `from_f64` / `to_f64` / `is_positive` etc. without adding
`num-bigint` / `num-rational` / `num-traits` as their own deps).
Gates `src/exact.rs`, additional tests, and the `exact_det_3x3`/`exact_sign_3x3`/`exact_solve_3x3` examples.
Clippy, doc builds, and test commands have dedicated `--features exact` variants.

Expand All @@ -245,7 +248,8 @@ When creating or updating issues:
- The public API re-exports these items from `src/lib.rs`.
- The `justfile` defines all dev workflows (see `just --list`).
- Dev-only benchmarks live in `benches/vs_linalg.rs` (Criterion + nalgebra/faer comparison)
and `benches/exact.rs` (exact arithmetic across D=2–5).
and `benches/exact.rs` (exact arithmetic across D=2–5, plus adversarial-input groups
`exact_near_singular_3x3`, `exact_large_entries_3x3`, `exact_hilbert_4x4`, `exact_hilbert_5x5`).
- Python scripts under `scripts/`:
- `bench_compare.py`: exact-arithmetic benchmark comparison across releases (generates `docs/PERFORMANCE.md`)
- `criterion_dim_plot.py`: benchmark plotting (CSV + SVG + README table update)
Expand Down
9 changes: 7 additions & 2 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -190,8 +190,13 @@ assert!((x[0] - 1.0).abs() <= f64::EPSILON);
assert!((x[1] - 2.0).abs() <= f64::EPSILON);
```

`BigRational` is re-exported from the crate root and prelude when the `exact`
feature is enabled, so consumers don't need to depend on `num-rational` directly.
With the `exact` feature enabled, `BigInt` and `BigRational` are re-exported
from the crate root and prelude, alongside the most commonly needed
`num-traits` items (`FromPrimitive`, `ToPrimitive`, `Signed`). This lets
consumers construct exact values (`BigRational::from_f64`, `from_i64`), query
sign (`is_positive` / `is_negative`), and convert back to `f64` (`to_f64`)
with a single `use la_stack::prelude::*;` — no need to add `num-bigint`,
`num-rational`, or `num-traits` to their own `Cargo.toml`.

For `det_sign_exact()`, D ≤ 4 matrices use a fast f64 filter (error-bounded
`det_direct()`) that resolves the sign without allocating. Only near-degenerate
Expand Down
191 changes: 166 additions & 25 deletions benches/exact.rs
Original file line number Diff line number Diff line change
@@ -1,10 +1,22 @@
//! Benchmarks for exact arithmetic operations.
//!
//! These benchmarks measure the performance of the `exact` feature's
//! arbitrary-precision methods across dimensions D=2..5 (the primary
//! target for geometric predicates).
//! arbitrary-precision methods. They are organised into two classes:
//!
//! 1. **General-case benches** (`exact_d{2..5}`) — a single
//! well-conditioned diagonally-dominant matrix per dimension. These
//! measure typical-case performance and track regressions against a
//! reproducible input.
//! 2. **Adversarial / extreme-input benches** — matrices chosen to
//! stress specific corners of the exact-arithmetic pipeline:
//! near-singularity (forces the Bareiss fallback), large f64 entries
//! (stresses intermediate `BigInt` growth), and Hilbert-style
//! ill-conditioning (wide range of `(mantissa, exponent)` pairs in
//! the `f64_decompose → BigInt` path). These measure tail behaviour
//! that fixed well-conditioned inputs miss and provide stronger
//! empirical evidence for `docs/PERFORMANCE.md`.

use criterion::Criterion;
use criterion::{BenchmarkGroup, Criterion, measurement::WallTime};
use la_stack::{Matrix, Vector};
use pastey::paste;
use std::hint::black_box;
Expand Down Expand Up @@ -47,8 +59,13 @@ fn make_vector_array<const D: usize>() -> [f64; D] {
}

/// Near-singular matrix: base singular matrix + tiny perturbation.
/// This forces the exact Bareiss fallback in `det_sign_exact` (the fast
/// f64 filter cannot resolve the sign).
///
/// The base `[[1,2,3],[4,5,6],[7,8,9]]` is exactly singular; adding
/// `2^-50` to entry (0,0) makes `det = -3 × 2^-50 ≠ 0`. The f64 filter
/// in `det_sign_exact` cannot resolve this sign, so Bareiss is forced;
/// `solve_exact` is the primary use case for near-degenerate inputs
/// (exact circumcenter etc.) and exercises the largest intermediate
/// `BigInt` values in the hybrid solve.
#[inline]
fn near_singular_3x3() -> Matrix<3> {
let perturbation = f64::from_bits(0x3CD0_0000_0000_0000); // 2^-50
Expand All @@ -59,6 +76,97 @@ fn near_singular_3x3() -> Matrix<3> {
])
}

/// Large-entry 3×3: strictly diagonally-dominant matrix with diagonal
/// entries near `f64::MAX / 2` and ones elsewhere.
///
/// Each big entry decomposes into a 53-bit mantissa with exponent `~970`;
/// the unit off-diagonals have exponent `0`, so the shared `e_min = 0`
/// shift in `component_to_bigint` produces `BigInt`s of `~1023` bits for
/// the diagonal and small integers elsewhere. Bareiss fraction-free
/// updates then multiply these together, stressing the big-integer
/// multiply and allocator along the full `O(D³)` elimination phase. The
/// matrix is non-singular (det ≈ `big³`) so both `det_*` and `solve_*`
/// paths complete.
#[inline]
fn large_entries_3x3() -> Matrix<3> {
let big = f64::MAX / 2.0;
Matrix::<3>::from_rows([[big, 1.0, 1.0], [1.0, big, 1.0], [1.0, 1.0, big]])
}

/// Hilbert matrix `H[i][j] = 1 / (i + j + 1)`.
///
/// Most entries (`1/3`, `1/5`, `1/6`, `1/7`, …) are non-terminating in
/// binary, so every cell has a distinct 53-bit mantissa and a small
/// negative exponent. `f64_decompose` therefore produces a wide mix of
/// `(mantissa, exponent)` pairs with no shared power-of-two factors,
/// and the scaling shift to the common `e_min` yields `BigInt` values
/// of varied bit-lengths — a different kind of adversarial input from
/// the large-entries case. Hilbert matrices are also classically
/// ill-conditioned (condition number grows exponentially with D), so
/// they are a realistic stand-in for the near-degenerate geometric
/// predicate inputs that motivate exact arithmetic.
#[inline]
#[allow(clippy::cast_precision_loss)]
fn hilbert<const D: usize>() -> Matrix<D> {
let mut rows = [[0.0; D]; D];
let mut r = 0;
while r < D {
let mut c = 0;
while c < D {
rows[r][c] = 1.0 / ((r + c + 1) as f64);
c += 1;
}
r += 1;
}
Matrix::<D>::from_rows(rows)
}

/// Populate a Criterion group with the four headline exact-arithmetic
/// benches on a single `(matrix, rhs)` pair: `det_sign_exact`,
/// `det_exact`, `solve_exact`, `solve_exact_f64`.
///
/// Used by every adversarial-input group so each one measures the same
/// operations, making the resulting tables directly comparable.
fn bench_extreme_group<const D: usize>(
group: &mut BenchmarkGroup<'_, WallTime>,
m: Matrix<D>,
rhs: Vector<D>,
) {
group.bench_function("det_sign_exact", |bencher| {
bencher.iter(|| {
let sign = black_box(m)
.det_sign_exact()
.expect("finite matrix entries");
black_box(sign);
});
});

group.bench_function("det_exact", |bencher| {
bencher.iter(|| {
let det = black_box(m).det_exact().expect("finite matrix entries");
black_box(det);
});
});

group.bench_function("solve_exact", |bencher| {
bencher.iter(|| {
let x = black_box(m)
.solve_exact(black_box(rhs))
.expect("non-singular matrix with finite entries");
let _ = black_box(x);
});
});

group.bench_function("solve_exact_f64", |bencher| {
bencher.iter(|| {
let x = black_box(m)
.solve_exact_f64(black_box(rhs))
.expect("solution representable in f64");
let _ = black_box(x);
});
});
}

macro_rules! gen_exact_benches_for_dim {
($c:expr, $d:literal) => {
paste! {{
Expand All @@ -72,7 +180,7 @@ macro_rules! gen_exact_benches_for_dim {
bencher.iter(|| {
let det = black_box(a)
.det(la_stack::DEFAULT_PIVOT_TOL)
.expect("should not fail");
.expect("diagonally dominant matrix is non-singular");
black_box(det);
});
});
Expand All @@ -87,39 +195,47 @@ macro_rules! gen_exact_benches_for_dim {
// === det_exact (BigRational result) ===
[<group_d $d>].bench_function("det_exact", |bencher| {
bencher.iter(|| {
let det = black_box(a).det_exact().expect("should not fail");
let det = black_box(a).det_exact().expect("finite matrix entries");
black_box(det);
});
});

// === det_exact_f64 (exact → f64) ===
[<group_d $d>].bench_function("det_exact_f64", |bencher| {
bencher.iter(|| {
let det = black_box(a).det_exact_f64().expect("should not fail");
let det = black_box(a)
.det_exact_f64()
.expect("det representable in f64");
black_box(det);
});
});

// === det_sign_exact (adaptive: fast filter + exact fallback) ===
[<group_d $d>].bench_function("det_sign_exact", |bencher| {
bencher.iter(|| {
let sign = black_box(a).det_sign_exact().expect("should not fail");
let sign = black_box(a)
.det_sign_exact()
.expect("finite matrix entries");
black_box(sign);
});
});

// === solve_exact (BigRational result) ===
[<group_d $d>].bench_function("solve_exact", |bencher| {
bencher.iter(|| {
let x = black_box(a).solve_exact(black_box(rhs)).expect("should not fail");
let x = black_box(a)
.solve_exact(black_box(rhs))
.expect("diagonally dominant matrix is non-singular");
black_box(x);
});
});

// === solve_exact_f64 (exact → f64) ===
[<group_d $d>].bench_function("solve_exact_f64", |bencher| {
bencher.iter(|| {
let x = black_box(a).solve_exact_f64(black_box(rhs)).expect("should not fail");
let x = black_box(a)
.solve_exact_f64(black_box(rhs))
.expect("solution representable in f64");
black_box(x);
});
});
Expand All @@ -140,25 +256,50 @@ fn main() {
gen_exact_benches_for_dim!(&mut c, 5);
}

// Near-singular 3×3: forces Bareiss fallback in det_sign_exact.
// === Adversarial / extreme-input groups ===
//
// Each group runs the same four exact-arithmetic benches
// (`det_sign_exact`, `det_exact`, `solve_exact`, `solve_exact_f64`)
// via `bench_extreme_group`, so the resulting tables are directly
// comparable across input classes.

// Near-singular 3×3: forces Bareiss fallback in det_sign_exact and
// exercises the largest intermediate BigInt values in solve_exact
// (the primary motivating use case for exact solve).
{
let m = near_singular_3x3();
let mut group = c.benchmark_group("exact_near_singular_3x3");
bench_extreme_group(
&mut group,
near_singular_3x3(),
Vector::<3>::new([1.0, 2.0, 3.0]),
);
group.finish();
}

group.bench_function("det_sign_exact", |bencher| {
bencher.iter(|| {
let sign = black_box(m).det_sign_exact().expect("should not fail");
black_box(sign);
});
});
// Large-entry 3×3: diagonal entries near `f64::MAX / 2` stress
// BigInt growth during Bareiss forward elimination.
{
let mut group = c.benchmark_group("exact_large_entries_3x3");
bench_extreme_group(
&mut group,
large_entries_3x3(),
Vector::<3>::new([1.0, 1.0, 1.0]),
);
group.finish();
}

group.bench_function("det_exact", |bencher| {
bencher.iter(|| {
let det = black_box(m).det_exact().expect("should not fail");
black_box(det);
});
});
// Hilbert 4×4 and 5×5: classically ill-conditioned matrices whose
// entries span many orders of magnitude in `(mantissa, exponent)`
// space, exercising the f64 → BigInt scaling path.
{
let mut group = c.benchmark_group("exact_hilbert_4x4");
bench_extreme_group(&mut group, hilbert::<4>(), Vector::<4>::new([1.0; 4]));
group.finish();
}

{
let mut group = c.benchmark_group("exact_hilbert_5x5");
bench_extreme_group(&mut group, hilbert::<5>(), Vector::<5>::new([1.0; 5]));
group.finish();
}

Expand Down
14 changes: 14 additions & 0 deletions docs/BENCHMARKING.md
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,20 @@ la-stack has two Criterion benchmark suites:
(`det_exact`, `solve_exact`, `det_sign_exact`, etc.) alongside f64
baselines (`det`, `det_direct`) across D=2–5. Use this to understand
the cost of exact arithmetic and track optimization progress.
In addition to the per-dimension groups (`exact_d{2..5}`), the suite
includes four adversarial-input groups designed to stress specific
corners of the pipeline:
- `exact_near_singular_3x3` — a 2^-50 perturbation of a singular base
matrix; forces the Bareiss fallback in `det_sign_exact` and
exercises the largest intermediate `BigInt` values in `solve_exact`.
- `exact_large_entries_3x3` — diagonal entries near `f64::MAX / 2`
stress `BigInt` growth during Bareiss forward elimination.
- `exact_hilbert_4x4` / `exact_hilbert_5x5` — classically
ill-conditioned matrices whose non-terminating-in-binary entries
stress the `f64_decompose → BigInt` scaling path.
Each adversarial group runs the same four benches (`det_sign_exact`,
`det_exact`, `solve_exact`, `solve_exact_f64`) so the resulting tables
are directly comparable across input classes.

## Quick reference

Expand Down
19 changes: 17 additions & 2 deletions scripts/bench_compare.py
Original file line number Diff line number Diff line change
Expand Up @@ -38,12 +38,22 @@
# ---------------------------------------------------------------------------

# Groups and the benchmarks within each group that we track.
#
# Mirrors the structure of `benches/exact.rs`: general-case per-dimension
# groups (`exact_d{2..5}`) plus adversarial/extreme-input groups that
# share a fixed four-bench layout (`det_sign_exact`, `det_exact`,
# `solve_exact`, `solve_exact_f64`).
_EXTREME_BENCHES: list[str] = ["det_sign_exact", "det_exact", "solve_exact", "solve_exact_f64"]

EXACT_GROUPS: dict[str, list[str]] = {
"exact_d2": ["det", "det_direct", "det_exact", "det_exact_f64", "det_sign_exact", "solve_exact", "solve_exact_f64"],
"exact_d3": ["det", "det_direct", "det_exact", "det_exact_f64", "det_sign_exact", "solve_exact", "solve_exact_f64"],
"exact_d4": ["det", "det_direct", "det_exact", "det_exact_f64", "det_sign_exact", "solve_exact", "solve_exact_f64"],
"exact_d5": ["det", "det_direct", "det_exact", "det_exact_f64", "det_sign_exact", "solve_exact", "solve_exact_f64"],
"exact_near_singular_3x3": ["det_sign_exact", "det_exact"],
"exact_near_singular_3x3": _EXTREME_BENCHES,
"exact_large_entries_3x3": _EXTREME_BENCHES,
"exact_hilbert_4x4": _EXTREME_BENCHES,
"exact_hilbert_5x5": _EXTREME_BENCHES,
}


Expand Down Expand Up @@ -197,11 +207,16 @@ def _group_by_group(items: list[_T]) -> dict[str, list[_T]]:

def _group_heading(group: str) -> str:
"""Turn a Criterion group name into a readable heading."""
# exact_d3 -> "D=3", exact_near_singular_3x3 -> "Near-singular 3x3"
# exact_d3 -> "D=3", exact_near_singular_3x3 -> "Near-singular 3x3",
# exact_hilbert_4x4 -> "Hilbert 4x4", etc.
if group.startswith("exact_d"):
return f"D={group.removeprefix('exact_d')}"
if group == "exact_near_singular_3x3":
return "Near-singular 3x3"
if group == "exact_large_entries_3x3":
return "Large entries 3x3"
if group.startswith("exact_hilbert_"):
return f"Hilbert {group.removeprefix('exact_hilbert_')}"
return group


Expand Down
9 changes: 9 additions & 0 deletions scripts/tests/test_bench_compare.py
Original file line number Diff line number Diff line change
Expand Up @@ -77,6 +77,15 @@ def test_dimension_group(self) -> None:
def test_near_singular(self) -> None:
assert bench_compare._group_heading("exact_near_singular_3x3") == "Near-singular 3x3"

def test_large_entries(self) -> None:
assert bench_compare._group_heading("exact_large_entries_3x3") == "Large entries 3x3"

def test_hilbert_4x4(self) -> None:
assert bench_compare._group_heading("exact_hilbert_4x4") == "Hilbert 4x4"

def test_hilbert_5x5(self) -> None:
assert bench_compare._group_heading("exact_hilbert_5x5") == "Hilbert 5x5"

def test_unknown_passthrough(self) -> None:
assert bench_compare._group_heading("something_else") == "something_else"

Expand Down
Loading
Loading