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
60 changes: 34 additions & 26 deletions benches/exact.rs
Original file line number Diff line number Diff line change
Expand Up @@ -19,8 +19,17 @@
use criterion::{BenchmarkGroup, Criterion, measurement::WallTime};
use la_stack::{Matrix, Vector};
use pastey::paste;
use std::fmt::Display;
use std::hint::black_box;

/// Return a successful benchmark operation result or panic with the named operation.
fn require_ok<T, E: Display>(result: Result<T, E>, operation: &str) -> T {
match result {
Ok(value) => value,
Err(err) => panic!("{operation} failed: {err}"),
}
}

#[inline]
#[allow(clippy::cast_precision_loss)]
const fn matrix_entry<const D: usize>(r: usize, c: usize) -> f64 {
Expand Down Expand Up @@ -134,34 +143,34 @@ fn bench_extreme_group<const D: usize>(
) {
group.bench_function("det_sign_exact", |bencher| {
bencher.iter(|| {
let sign = black_box(m)
.det_sign_exact()
.expect("finite matrix entries");
let sign = require_ok(black_box(m).det_sign_exact(), "exact determinant sign");
black_box(sign);
});
});

group.bench_function("det_exact", |bencher| {
bencher.iter(|| {
let det = black_box(m).det_exact().expect("finite matrix entries");
let det = require_ok(black_box(m).det_exact(), "exact determinant");
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 x = require_ok(
black_box(m).solve_exact(black_box(rhs)),
"exact linear solve",
);
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 x = require_ok(
black_box(m).solve_exact_f64(black_box(rhs)),
"exact linear solve converted to f64",
);
let _ = black_box(x);
});
});
Expand All @@ -178,9 +187,7 @@ macro_rules! gen_exact_benches_for_dim {
// === f64 baselines ===
[<group_d $d>].bench_function("det", |bencher| {
bencher.iter(|| {
let det = black_box(a)
.det()
.expect("diagonally dominant matrix is non-singular");
let det = require_ok(black_box(a).det(), "f64 determinant");
black_box(det);
});
});
Expand All @@ -195,47 +202,48 @@ 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("finite matrix entries");
let det = require_ok(black_box(a).det_exact(), "exact determinant");
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("det representable in f64");
let det = require_ok(
black_box(a).det_exact_f64(),
"exact determinant converted to 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("finite matrix entries");
let sign = require_ok(black_box(a).det_sign_exact(), "exact determinant sign");
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("diagonally dominant matrix is non-singular");
let x = require_ok(
black_box(a).solve_exact(black_box(rhs)),
"exact linear solve",
);
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("solution representable in f64");
let x = require_ok(
black_box(a).solve_exact_f64(black_box(rhs)),
"exact linear solve converted to f64",
);
black_box(x);
});
});
Expand Down
82 changes: 46 additions & 36 deletions benches/vs_linalg.rs
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,22 @@ use faer::linalg::solvers::{PartialPivLu, Solve};
use faer::perm::PermRef;
use la_stack::{DEFAULT_PIVOT_TOL, Matrix, Vector};
use pastey::paste;
use std::fmt::Display;
use std::hint::black_box;

/// Return a successful benchmark operation result or panic with the named operation.
fn require_ok<T, E: Display>(result: Result<T, E>, operation: &str) -> T {
match result {
Ok(value) => value,
Err(err) => panic!("{operation} failed: {err}"),
}
}

/// Return a present third-party benchmark result or panic with the named operation.
fn require_some<T>(value: Option<T>, operation: &str) -> T {
value.unwrap_or_else(|| panic!("{operation} returned no result"))
}

fn faer_perm_sign(p: PermRef<'_, usize>) -> f64 {
// Sign(det(P)) for a permutation matrix P is +1 for even permutations, -1 for odd.
// Parity can be computed from the number of cycles:
Expand Down Expand Up @@ -145,9 +159,7 @@ macro_rules! gen_vs_linalg_benches_for_dim {
let fv2 = faer::Mat::<f64>::from_fn($d, 1, |i, _| vector_entry(i, 1.0));

// Precompute LU once for solve-only / det-only benchmarks.
let a_lu = a
.lu(DEFAULT_PIVOT_TOL)
.expect("matrix should be non-singular");
let a_lu = require_ok(a.lu(DEFAULT_PIVOT_TOL), "precomputed la_stack LU");
let na_lu = na.clone().lu();
let fa_lu = fa.partial_piv_lu();

Expand All @@ -156,13 +168,11 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Determinant via LU (factor + det) ===
[<group_d $d>].bench_function("la_stack_det_via_lu", |bencher| {
bencher.iter(|| {
let lu = black_box(a)
.lu(DEFAULT_PIVOT_TOL)
.expect("matrix should be non-singular");
let det = match lu.det() {
Ok(det) => det,
Err(err) => panic!("finite benchmark matrix determinant failed: {err}"),
};
let lu = require_ok(
black_box(a).lu(DEFAULT_PIVOT_TOL),
"la_stack LU factorization",
);
let det = require_ok(lu.det(), "la_stack LU determinant");
black_box(det);
});
});
Expand All @@ -186,17 +196,18 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Determinant via det() (closed-form for D≤4, LU for D≥5) ===
[<group_d $d>].bench_function("la_stack_det", |bencher| {
bencher.iter(|| {
let det = black_box(a).det().expect("matrix should be non-singular");
let det = require_ok(black_box(a).det(), "la_stack determinant");
black_box(det);
});
});

// === LU factorization ===
[<group_d $d>].bench_function("la_stack_lu", |bencher| {
bencher.iter(|| {
let lu = black_box(a)
.lu(DEFAULT_PIVOT_TOL)
.expect("matrix should be non-singular");
let lu = require_ok(
black_box(a).lu(DEFAULT_PIVOT_TOL),
"la_stack LU factorization",
);
let _ = black_box(lu);
});
});
Expand All @@ -218,22 +229,22 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === LU solve (factor + solve) ===
[<group_d $d>].bench_function("la_stack_lu_solve", |bencher| {
bencher.iter(|| {
let lu = black_box(a)
.lu(DEFAULT_PIVOT_TOL)
.expect("matrix should be non-singular");
let x = lu
.solve_vec(black_box(rhs))
.expect("solve should succeed");
let lu = require_ok(
black_box(a).lu(DEFAULT_PIVOT_TOL),
"la_stack LU factorization",
);
let x = require_ok(
lu.solve_vec(black_box(rhs)),
"la_stack LU solve",
);
let _ = black_box(x);
});
});

[<group_d $d>].bench_function("nalgebra_lu_solve", |bencher| {
bencher.iter(|| {
let lu = black_box(na.clone()).lu();
let x = lu
.solve(black_box(&nrhs))
.expect("solve should succeed");
let x = require_some(lu.solve(black_box(&nrhs)), "nalgebra LU solve");
black_box(x);
});
});
Expand All @@ -249,18 +260,20 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Solve using a precomputed LU ===
[<group_d $d>].bench_function("la_stack_solve_from_lu", |bencher| {
bencher.iter(|| {
let x = a_lu
.solve_vec(black_box(rhs))
.expect("solve should succeed");
let x = require_ok(
a_lu.solve_vec(black_box(rhs)),
"precomputed la_stack LU solve",
);
let _ = black_box(x);
});
});

[<group_d $d>].bench_function("nalgebra_solve_from_lu", |bencher| {
bencher.iter(|| {
let x = na_lu
.solve(black_box(&nrhs))
.expect("solve should succeed");
let x = require_some(
na_lu.solve(black_box(&nrhs)),
"precomputed nalgebra LU solve",
);
black_box(x);
});
});
Expand All @@ -275,10 +288,7 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Determinant from a precomputed LU ===
[<group_d $d>].bench_function("la_stack_det_from_lu", |bencher| {
bencher.iter(|| {
let det = match a_lu.det() {
Ok(det) => det,
Err(err) => panic!("finite benchmark matrix determinant failed: {err}"),
};
let det = require_ok(a_lu.det(), "precomputed la_stack LU determinant");
black_box(det);
});
});
Expand All @@ -300,7 +310,7 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Vector dot product ===
[<group_d $d>].bench_function("la_stack_dot", |bencher| {
bencher.iter(|| {
let result = black_box(v1).dot(black_box(v2)).unwrap();
let result = require_ok(black_box(v1).dot(black_box(v2)), "la_stack dot");
black_box(result);
});
});
Expand All @@ -327,7 +337,7 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Vector norm squared ===
[<group_d $d>].bench_function("la_stack_norm2_sq", |bencher| {
bencher.iter(|| {
let result = black_box(v1).norm2_sq().unwrap();
let result = require_ok(black_box(v1).norm2_sq(), "la_stack norm2_sq");
black_box(result);
});
});
Expand All @@ -354,7 +364,7 @@ macro_rules! gen_vs_linalg_benches_for_dim {
// === Matrix infinity norm (max absolute row sum) ===
[<group_d $d>].bench_function("la_stack_inf_norm", |bencher| {
bencher.iter(|| {
let result = black_box(a).inf_norm().unwrap();
let result = require_ok(black_box(a).inf_norm(), "la_stack inf_norm");
black_box(result);
});
});
Expand Down
6 changes: 4 additions & 2 deletions examples/const_det_4x4.rs
Original file line number Diff line number Diff line change
Expand Up @@ -20,18 +20,20 @@ const DET: f64 = match MAT.det_direct() {
Err(_) => panic!("matrix entries must be finite"),
};

fn main() {
fn main() -> Result<(), LaError> {
println!("4×4 matrix:");
for r in 0..4 {
print!(" [");
for c in 0..4 {
if c > 0 {
print!(", ");
}
print!("{:5.1}", MAT.get(r, c).unwrap());
print!("{:5.1}", MAT.get_checked(r, c)?);
}
println!("]");
}
println!();
println!("det (computed at compile time) = {DET}");

Ok(())
}
10 changes: 6 additions & 4 deletions examples/exact_sign_3x3.rs
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@

use la_stack::prelude::*;

fn main() {
fn main() -> Result<(), LaError> {
// Base matrix: rows in arithmetic progression → exactly singular (det = 0).
// [[1, 2, 3],
// [4, 5, 6],
Expand All @@ -23,8 +23,8 @@ fn main() {
[7.0, 8.0, 9.0],
]);

let sign = m.det_sign_exact().unwrap();
let det_f64 = m.det().unwrap();
let sign = m.det_sign_exact()?;
let det_f64 = m.det()?;

println!("Near-singular 3×3 matrix (perturbation = 2^-50 ≈ {perturbation:.2e}):");
for r in 0..3 {
Expand All @@ -33,7 +33,7 @@ fn main() {
if c > 0 {
print!(", ");
}
print!("{:22.18}", m.get(r, c).unwrap());
print!("{:22.18}", m.get_checked(r, c)?);
}
println!("]");
}
Expand All @@ -42,4 +42,6 @@ fn main() {
println!("det_sign_exact() = {sign}");
println!();
println!("The exact sign is −1 (negative), matching the analytical result.");

Ok(())
}
2 changes: 1 addition & 1 deletion examples/ldlt_solve_3x3.rs
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,7 @@ fn main() -> Result<(), LaError> {
if c > 0 {
print!(", ");
}
print!("{:5.1}", a.get(r, c).unwrap());
print!("{:5.1}", a.get_checked(r, c)?);
}
println!("]");
}
Expand Down
Loading
Loading