Skip to content

Commit 56973bc

Browse files
acgetchellacgetchell
authored
Merge pull request #91 from acgetchell/chose/67-bump-msrv-1.95
chore: bump MSRV to Rust 1.95 and adopt new stable features
2 parents 7f3d2e1 + f763b11 commit 56973bc

6 files changed

Lines changed: 96 additions & 9 deletions

File tree

‎Cargo.toml‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -2,7 +2,7 @@
22
name = "la-stack"
33
version = "0.4.0"
44
edition = "2024"
5-
rust-version = "1.94"
5+
rust-version = "1.95"
66
license = "BSD-3-Clause"
77
description = "Fast, stack-allocated linear algebra for fixed dimensions"
88
readme = "README.md"

‎rust-toolchain.toml‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
[toolchain]
22
# Pin to MSRV as specified in Cargo.toml
3-
channel = "1.94.0"
3+
channel = "1.95.0"
44

55
# Essential components for development
66
components = [

‎src/exact.rs‎

Lines changed: 20 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -43,6 +43,7 @@
4343
//! \[10\] for background on floating-point representation and exact
4444
//! rational reconstruction. Reference numbers refer to `REFERENCES.md`.
4545
46+
use core::hint::cold_path;
4647
use std::array::from_fn;
4748

4849
use num_bigint::{BigInt, Sign};
@@ -61,6 +62,7 @@ fn validate_finite<const D: usize>(m: &Matrix<D>) -> Result<(), LaError> {
6162
for r in 0..D {
6263
for c in 0..D {
6364
if !m.rows[r][c].is_finite() {
65+
cold_path();
6466
return Err(LaError::NonFinite {
6567
row: Some(r),
6668
col: c,
@@ -78,6 +80,7 @@ fn validate_finite<const D: usize>(m: &Matrix<D>) -> Result<(), LaError> {
7880
fn validate_finite_vec<const D: usize>(v: &Vector<D>) -> Result<(), LaError> {
7981
for (i, &x) in v.data.iter().enumerate() {
8082
if !x.is_finite() {
83+
cold_path();
8184
return Err(LaError::NonFinite { row: None, col: i });
8285
}
8386
}
@@ -324,6 +327,7 @@ fn gauss_solve<const D: usize>(m: &Matrix<D>, b: &Vector<D>) -> Result<[BigRatio
324327
mat.swap(k, swap_row);
325328
rhs.swap(k, swap_row);
326329
} else {
330+
cold_path();
327331
return Err(LaError::Singular { pivot_col: k });
328332
}
329333
}
@@ -419,6 +423,7 @@ impl<const D: usize> Matrix<D> {
419423
if val.is_finite() {
420424
Ok(val)
421425
} else {
426+
cold_path();
422427
Err(LaError::Overflow { index: None })
423428
}
424429
}
@@ -490,6 +495,7 @@ impl<const D: usize> Matrix<D> {
490495
for (i, val) in exact.iter().enumerate() {
491496
let f = val.to_f64().unwrap_or(f64::INFINITY);
492497
if !f.is_finite() {
498+
cold_path();
493499
return Err(LaError::Overflow { index: Some(i) });
494500
}
495501
result[i] = f;
@@ -537,23 +543,31 @@ impl<const D: usize> Matrix<D> {
537543
validate_finite(self)?;
538544

539545
// Stage 1: f64 fast filter for D ≤ 4.
540-
if let (Some(det_f64), Some(err)) = (self.det_direct(), self.det_errbound()) {
541-
// When entries are large (e.g. near f64::MAX) the determinant can
542-
// overflow to infinity even though every individual entry is finite.
543-
// In that case the fast filter is inconclusive; fall through to the
544-
// exact Bareiss path.
545-
if det_f64.is_finite() {
546+
//
547+
// When entries are large (e.g. near f64::MAX) the determinant can
548+
// overflow to infinity even though every individual entry is finite.
549+
// In that case the fast filter is inconclusive; fall through to the
550+
// exact Bareiss path.
551+
match self.det_direct() {
552+
Some(det_f64)
553+
if let Some(err) = self.det_errbound()
554+
&& det_f64.is_finite() =>
555+
{
546556
if det_f64 > err {
547557
return Ok(1);
548558
}
549559
if det_f64 < -err {
550560
return Ok(-1);
551561
}
552562
}
563+
_ => {}
553564
}
554565

555566
// Stage 2: integer Bareiss fallback — the 2^(D×e_min) scale factor
556567
// is always positive, so det_int.sign() == det(A).sign().
568+
// This is the cold path: the fast filter resolves the vast majority of
569+
// well-conditioned calls without allocating.
570+
cold_path();
557571
let (det_int, _) = bareiss_det_int(self);
558572
Ok(match det_int.sign() {
559573
Sign::Plus => 1,

‎src/ldlt.rs‎

Lines changed: 26 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -4,6 +4,8 @@
44
//! symmetric positive definite (SPD) and positive semi-definite (PSD) matrices (e.g. Gram
55
//! matrices) without pivoting.
66
7+
use core::hint::cold_path;
8+
79
use crate::LaError;
810
use crate::matrix::Matrix;
911
use crate::vector::Vector;
@@ -39,19 +41,22 @@ impl<const D: usize> Ldlt<D> {
3941
for j in 0..D {
4042
let d = f.rows[j][j];
4143
if !d.is_finite() {
44+
cold_path();
4245
return Err(LaError::NonFinite {
4346
row: Some(j),
4447
col: j,
4548
});
4649
}
4750
if d <= tol {
51+
cold_path();
4852
return Err(LaError::Singular { pivot_col: j });
4953
}
5054

5155
// Compute L multipliers below the diagonal in column j.
5256
for i in (j + 1)..D {
5357
let l = f.rows[i][j] / d;
5458
if !l.is_finite() {
59+
cold_path();
5560
return Err(LaError::NonFinite {
5661
row: Some(i),
5762
col: j,
@@ -69,6 +74,7 @@ impl<const D: usize> Ldlt<D> {
6974
let l_k = f.rows[k][j];
7075
let new_val = (-l_i_d).mul_add(l_k, f.rows[i][k]);
7176
if !new_val.is_finite() {
77+
cold_path();
7278
return Err(LaError::NonFinite {
7379
row: Some(i),
7480
col: k,
@@ -141,6 +147,7 @@ impl<const D: usize> Ldlt<D> {
141147
sum = (-row[j]).mul_add(*x_j, sum);
142148
}
143149
if !sum.is_finite() {
150+
cold_path();
144151
return Err(LaError::NonFinite { row: None, col: i });
145152
}
146153
x[i] = sum;
@@ -150,14 +157,17 @@ impl<const D: usize> Ldlt<D> {
150157
for (i, x_i) in x.iter_mut().enumerate().take(D) {
151158
let diag = self.factors.rows[i][i];
152159
if !diag.is_finite() {
160+
cold_path();
153161
return Err(LaError::NonFinite { row: None, col: i });
154162
}
155163
if diag <= self.tol {
164+
cold_path();
156165
return Err(LaError::Singular { pivot_col: i });
157166
}
158167

159168
let v = *x_i / diag;
160169
if !v.is_finite() {
170+
cold_path();
161171
return Err(LaError::NonFinite { row: None, col: i });
162172
}
163173
*x_i = v;
@@ -171,6 +181,7 @@ impl<const D: usize> Ldlt<D> {
171181
sum = (-self.factors.rows[j][i]).mul_add(*x_j, sum);
172182
}
173183
if !sum.is_finite() {
184+
cold_path();
174185
return Err(LaError::NonFinite { row: None, col: i });
175186
}
176187
x[i] = sum;
@@ -407,4 +418,19 @@ mod tests {
407418
let err = ldlt.solve_vec(b).unwrap_err();
408419
assert_eq!(err, LaError::NonFinite { row: None, col: 1 });
409420
}
421+
422+
#[test]
423+
fn nonfinite_solve_vec_diagonal_solve_overflow() {
424+
// Diagonal SPD matrix with a tiny diagonal entry just above the
425+
// singularity tolerance. Forward substitution passes through the
426+
// large RHS unchanged, then the diagonal solve z[1] = y[1] / D[1]
427+
// = 1e300 / 1e-11 = 1e311 overflows f64, exercising the
428+
// `!v.is_finite()` branch of the diagonal solve.
429+
let a = Matrix::<2>::from_rows([[1.0, 0.0], [0.0, 1.0e-11]]);
430+
let ldlt = a.ldlt(DEFAULT_SINGULAR_TOL).unwrap();
431+
432+
let b = Vector::<2>::new([0.0, 1.0e300]);
433+
let err = ldlt.solve_vec(b).unwrap_err();
434+
assert_eq!(err, LaError::NonFinite { row: None, col: 1 });
435+
}
410436
}

‎src/lu.rs‎

Lines changed: 27 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,7 @@
11
//! LU decomposition and solves.
22
3+
use core::hint::cold_path;
4+
35
use crate::LaError;
46
use crate::matrix::Matrix;
57
use crate::vector::Vector;
@@ -31,6 +33,7 @@ impl<const D: usize> Lu<D> {
3133
let mut pivot_row = k;
3234
let mut pivot_abs = lu.rows[k][k].abs();
3335
if !pivot_abs.is_finite() {
36+
cold_path();
3437
return Err(LaError::NonFinite {
3538
row: Some(k),
3639
col: k,
@@ -40,6 +43,7 @@ impl<const D: usize> Lu<D> {
4043
for r in (k + 1)..D {
4144
let v = lu.rows[r][k].abs();
4245
if !v.is_finite() {
46+
cold_path();
4347
return Err(LaError::NonFinite {
4448
row: Some(r),
4549
col: k,
@@ -52,6 +56,7 @@ impl<const D: usize> Lu<D> {
5256
}
5357

5458
if pivot_abs <= tol {
59+
cold_path();
5560
return Err(LaError::Singular { pivot_col: k });
5661
}
5762

@@ -63,6 +68,7 @@ impl<const D: usize> Lu<D> {
6368

6469
let pivot = lu.rows[k][k];
6570
if !pivot.is_finite() {
71+
cold_path();
6672
return Err(LaError::NonFinite {
6773
row: Some(k),
6874
col: k,
@@ -73,6 +79,7 @@ impl<const D: usize> Lu<D> {
7379
for r in (k + 1)..D {
7480
let mult = lu.rows[r][k] / pivot;
7581
if !mult.is_finite() {
82+
cold_path();
7683
return Err(LaError::NonFinite {
7784
row: Some(r),
7885
col: k,
@@ -132,6 +139,7 @@ impl<const D: usize> Lu<D> {
132139
sum = (-row[j]).mul_add(*x_j, sum);
133140
}
134141
if !sum.is_finite() {
142+
cold_path();
135143
return Err(LaError::NonFinite { row: None, col: i });
136144
}
137145
x[i] = sum;
@@ -148,14 +156,17 @@ impl<const D: usize> Lu<D> {
148156

149157
let diag = row[i];
150158
if !diag.is_finite() || !sum.is_finite() {
159+
cold_path();
151160
return Err(LaError::NonFinite { row: None, col: i });
152161
}
153162
if diag.abs() <= self.tol {
163+
cold_path();
154164
return Err(LaError::Singular { pivot_col: i });
155165
}
156166

157167
let q = sum / diag;
158168
if !q.is_finite() {
169+
cold_path();
159170
return Err(LaError::NonFinite { row: None, col: i });
160171
}
161172
x[i] = q;
@@ -474,4 +485,20 @@ mod tests {
474485
let err = lu.solve_vec(b).unwrap_err();
475486
assert_eq!(err, LaError::NonFinite { row: None, col: 1 });
476487
}
488+
489+
#[test]
490+
fn solve_vec_nonfinite_back_substitution_sum_overflow() {
491+
// Upper-triangular U with a very large off-diagonal in row 1 and a
492+
// very large x[2] produced by the RHS. The back-substitution
493+
// accumulator `sum = (-row[j]).mul_add(x[j], sum)` overflows while
494+
// reducing row 1, so the failure is detected via the `!sum.is_finite()`
495+
// branch of the combined diag/sum check (distinct from the
496+
// `q = sum / diag` overflow path covered above).
497+
let a = Matrix::<3>::from_rows([[1.0, 0.0, 0.0], [0.0, 1.0, 1.0e200], [0.0, 0.0, 1.0]]);
498+
let lu = a.lu(DEFAULT_PIVOT_TOL).unwrap();
499+
500+
let b = Vector::<3>::new([0.0, 0.0, 1.0e200]);
501+
let err = lu.solve_vec(b).unwrap_err();
502+
assert_eq!(err, LaError::NonFinite { row: None, col: 1 });
503+
}
477504
}

‎src/matrix.rs‎

Lines changed: 21 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,5 +1,7 @@
11
//! Fixed-size, stack-allocated square matrices.
22
3+
use core::hint::cold_path;
4+
35
use crate::LaError;
46
use crate::ldlt::Ldlt;
57
use crate::lu::Lu;
@@ -266,7 +268,11 @@ impl<const D: usize> Matrix<D> {
266268
(-r[0][1]).mul_add(c01, r[0][2].mul_add(c02, -(r[0][3] * c03))),
267269
))
268270
}
269-
_ => None,
271+
_ => {
272+
// Cold in the common D ≤ 4 case; callers fall back to LU for D ≥ 5.
273+
cold_path();
274+
None
275+
}
270276
}
271277
}
272278

@@ -296,6 +302,7 @@ impl<const D: usize> Matrix<D> {
296302
return if d.is_finite() {
297303
Ok(d)
298304
} else {
305+
cold_path();
299306
// Scan for the first non-finite entry to preserve coordinates.
300307
for r in 0..D {
301308
for c in 0..D {
@@ -703,6 +710,19 @@ mod tests {
703710
);
704711
}
705712

713+
#[test]
714+
fn det_returns_nonfinite_error_for_overflow_with_finite_entries() {
715+
// det_direct produces an overflowing f64 (1e300 * 1e300 = ∞) even
716+
// though every matrix entry is finite. The entry scan in `det`
717+
// falls through and returns NonFinite { row: None, col: 0 } to signal
718+
// a computed overflow rather than a NaN/∞ input.
719+
let m = Matrix::<2>::from_rows([[1e300, 0.0], [0.0, 1e300]]);
720+
assert_eq!(
721+
m.det(DEFAULT_PIVOT_TOL),
722+
Err(LaError::NonFinite { row: None, col: 0 })
723+
);
724+
}
725+
706726
#[test]
707727
fn det_direct_is_const_evaluable_d2() {
708728
// Const evaluation proves the function is truly const fn.

0 commit comments

Comments
 (0)