Mercurial > repos > alg_tools / changeset

--- a/Cargo.toml	Tue Apr 08 13:30:08 2025 -0500
+++ b/Cargo.toml	Tue Apr 08 13:30:12 2025 -0500
@@ -8,7 +8,13 @@
 homepage = "https://tuomov.iki.fi/software/alg_tools/"
 repository = "https://tuomov.iki.fi/repos/alg_tools/"
 license-file = "LICENSE"
-keywords = ["iterative", "optimization", "bisection", "branch-and-bound", "numerical"]
+keywords = [
+    "iterative",
+    "optimization",
+    "bisection",
+    "branch-and-bound",
+    "numerical",
+]
 categories = ["mathematics", "data-structures"]
 
 [dependencies]
@@ -27,7 +33,7 @@
 anyhow = "1.0.95"
 
 [package.metadata.docs.rs]
-rustdoc-args = [ "--html-in-header", "katex-header.html" ]
+rustdoc-args = ["--html-in-header", "katex-header.html"]
 
 
 [profile.release]
@@ -36,5 +42,4 @@
 [features]
 default = ["nightly"]
 use_custom_thread_pool = []
-nightly = [] # enable for higher-performance implementations of some things
-
+nightly = []                # enable for higher-performance implementations of some things

--- a/src/linsolve.rs	Tue Apr 08 13:30:08 2025 -0500
+++ b/src/linsolve.rs	Tue Apr 08 13:30:12 2025 -0500
@@ -8,33 +8,33 @@
 
 /// Gaussian elimination for $AX=B$, where $A$ and $B$ are both stored in `ab`,
 /// $A \in \mathbb{R}^{M \times M}$ and $X, B \in \mathbb{R}^{M \times K}$.
-pub fn linsolve0<F : Float, const M : usize, const N : usize, const K : usize>(
-    mut ab : [[F; N]; M]
+pub fn linsolve0<F: Float, const M: usize, const N: usize, const K: usize>(
+    mut ab: [[F; N]; M],
 ) -> [[F; K]; M] {
-    assert_eq!(M + K , N);
+    assert_eq!(M + K, N);
 
     let mut k = 0;
 
     // Convert to row-echelon form
-    for h in 0..(M-1) {
+    for h in 0..(M - 1) {
         // Find pivotable column (has some non-zero entries in rows ≥ h)
         'find_pivot: while k < N {
             let (mut î, mut v) = (h, ab[h][k].abs());
             // Find row ≥ h of maximum absolute value in this column
-            for i in (h+1)..M {
+            for i in (h + 1)..M {
                 let ṽ = ab[i][k].abs();
-                if ṽ > v  {
+                if ṽ > v {
                     î = i;
                     v = ṽ;
                 }
             }
             if v > F::ZERO {
                 ab.swap(h, î);
-                for i in (h+1)..M {
+                for i in (h + 1)..M {
                     let f = ab[i][k] / ab[h][k];
                     ab[i][k] = F::ZERO;
-                    for j in (k+1)..N {
-                        ab[i][j] -= ab[h][j]*f;
+                    for j in (k + 1)..N {
+                        ab[i][j] -= ab[h][j] * f;
                     }
                 }
                 k += 1;
@@ -56,8 +56,8 @@
         //unsafe { std::mem::MaybeUninit::uninit().assume_init() };
         for i in (0..M).rev() {
             for 𝓁 in 0..K {
-                let mut tmp  = ab[i][M+𝓁];
-                for j in (i+1)..M {
+                let mut tmp = ab[i][M + 𝓁];
+                for j in (i + 1)..M {
                     tmp -= ab[i][j] * unsafe { *(x[j][𝓁].assume_init_ref()) };
                 }
                 tmp /= ab[i][i];
@@ -71,11 +71,11 @@
     }
     #[cfg(not(feature = "nightly"))]
     {
-        let mut x : [[F; K]; M] = [[F::ZERO; K]; M];
+        let mut x: [[F; K]; M] = [[F::ZERO; K]; M];
         for i in (0..M).rev() {
             for 𝓁 in 0..K {
-                let mut tmp  = ab[i][M+𝓁];
-                for j in (i+1)..M {
+                let mut tmp = ab[i][M + 𝓁];
+                for j in (i + 1)..M {
                     tmp -= ab[i][j] * x[j][𝓁];
                 }
                 tmp /= ab[i][i];
@@ -89,12 +89,11 @@
 /// Gaussian elimination for $Ax=b$, where $A$ and $b$ are both stored in `ab`,
 /// $A \in \mathbb{R}^{M \times M}$ and $x, b \in \mathbb{R}^M$.
 #[inline]
-pub fn linsolve<F : Float, const M : usize, const N : usize>(ab : [[F; N]; M]) -> [F; M] {
-    let x : [[F; 1]; M] = linsolve0(ab);
-    unsafe { *((&x as *const [F; 1]) as *const [F; M] ) }
+pub fn linsolve<F: Float, const M: usize, const N: usize>(ab: [[F; N]; M]) -> [F; M] {
+    let x: [[F; 1]; M] = linsolve0(ab);
+    unsafe { *((&x as *const [F; 1]) as *const [F; M]) }
 }
 
-
 #[cfg(test)]
 mod tests {
     use super::*;
@@ -103,8 +102,11 @@
     fn linsolve_test() {
         let ab1 = [[1.0, 2.0, 3.0], [2.0, 1.0, 6.0]];
         assert_eq!(linsolve(ab1), [3.0, 0.0]);
-        let ab2 = [[1.0, 2.0, 0.0, 1.0], [4.0, 0.0, 0.0, 0.0], [0.0, 0.0, 1.0, 0.0]];
+        let ab2 = [
+            [1.0, 2.0, 0.0, 1.0],
+            [4.0, 0.0, 0.0, 0.0],
+            [0.0, 0.0, 1.0, 0.0],
+        ];
         assert_eq!(linsolve(ab2), [0.0, 0.5, 0.0]);
     }
 }
-