rustebra 0.5.0

A hybrid no_std/alloc linear algebra crate for Rust, scaling from embedded targets to dynamic Krylov subspace solvers.
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88

/// The `K x K` symmetric tridiagonal matrix produced by [`crate::krylov::lanczos`]: the
/// projection `Qᵗ * A * Q` of a symmetric matrix onto the orthonormal Krylov basis `Q`.
///
/// Only the diagonal (`K` entries) and the first sub/super-diagonal (`K - 1` entries, shared
/// by symmetry) are stored; every other entry is zero by construction. Its eigenvalues
/// approximate eigenvalues of the projected matrix — exactly its spectrum when `K` equals the
/// matrix dimension and no breakdown occurred.
///
/// # Examples
///
/// ```
/// use rustebra::krylov::lanczos;
/// use rustebra::storage::{Basis, StaticStorage};
///
/// // Already tridiagonal: [[2, 1], [1, 2]]. Starting from e1, Lanczos reproduces it.
/// let a = StaticStorage::new([2.0_f64, 1.0, 1.0, 2.0]);
/// let v0 = StaticStorage::new([1.0, 0.0]);
/// let mut buffer = [0.0; 4];
/// let mut basis = Basis::<f64, 2>::new(&mut buffer, 2).unwrap();
/// let mut scratch = [0.0; 2];
/// let t = lanczos(&a, 2, &v0, 1e-12, &mut basis, &mut scratch).unwrap();
///
/// assert!((t.diagonal()[0] - 2.0).abs() < 1e-12);
/// assert!((t.diagonal()[1] - 2.0).abs() < 1e-12);
/// assert_eq!(t.off_diagonal().len(), 1);
/// assert!((t.off_diagonal()[0] - 1.0).abs() < 1e-12);
/// ```
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct TridiagonalMatrix<T, const K: usize> {
    diagonal: [T; K],
    // `[T; K - 1]` isn't expressible on stable Rust, so the array carries one padding slot;
    // `off_diagonal()` hides it.
    off_diagonal: [T; K],
}

impl<T, const K: usize> TridiagonalMatrix<T, K> {
    /// Assembles the result of a Lanczos run; `off_diagonal[K - 1]` is padding, never read.
    pub(super) fn new(diagonal: [T; K], off_diagonal: [T; K]) -> Self {
        Self {
            diagonal,
            off_diagonal,
        }
    }

    /// The `K` diagonal entries `α_0 ..= α_{K-1}`.
    ///
    /// # Examples
    ///
    /// ```
    /// use rustebra::krylov::lanczos;
    /// use rustebra::storage::{Basis, StaticStorage};
    ///
    /// let a = StaticStorage::new([3.0_f64]);
    /// let v0 = StaticStorage::new([1.0]);
    /// let mut buffer = [0.0; 1];
    /// let mut basis = Basis::<f64, 1>::new(&mut buffer, 1).unwrap();
    /// let mut scratch = [0.0; 1];
    /// let t = lanczos(&a, 1, &v0, 1e-12, &mut basis, &mut scratch).unwrap();
    ///
    /// assert_eq!(t.diagonal(), &[3.0]);
    /// ```
    pub fn diagonal(&self) -> &[T; K] {
        &self.diagonal
    }

    /// The `K - 1` sub/super-diagonal entries `β_0 ..= β_{K-2}` (shared by symmetry):
    /// `off_diagonal()[j]` couples rows `j` and `j + 1`.
    ///
    /// # Examples
    ///
    /// ```
    /// use rustebra::krylov::lanczos;
    /// use rustebra::storage::{Basis, StaticStorage};
    ///
    /// let a = StaticStorage::new([3.0_f64]);
    /// let v0 = StaticStorage::new([1.0]);
    /// let mut buffer = [0.0; 1];
    /// let mut basis = Basis::<f64, 1>::new(&mut buffer, 1).unwrap();
    /// let mut scratch = [0.0; 1];
    /// let t = lanczos(&a, 1, &v0, 1e-12, &mut basis, &mut scratch).unwrap();
    ///
    /// // A 1x1 tridiagonal matrix has no off-diagonal at all.
    /// assert!(t.off_diagonal().is_empty());
    /// ```
    pub fn off_diagonal(&self) -> &[T] {
        &self.off_diagonal[..K.saturating_sub(1)]
    }
}