opensrdk_linear_algebra/matrix/st/
evd.rs

1use super::SymmetricTridiagonalMatrix;
2use crate::{matrix::ge::Matrix, matrix::*};
3use lapack::dstevd;
4
5impl SymmetricTridiagonalMatrix {
6    /// # Eigen decomposition
7    /// return (lambda, pt)
8    pub fn stevd(self) -> Result<(Vec<f64>, Matrix), MatrixError> {
9        let (mut d, mut e) = self.eject();
10        let n = d.len();
11        let mut z = Matrix::new(n, n);
12        let lwork = 1.max(1 + 4 * n + n.pow(2));
13        let mut work = vec![0.0; lwork];
14        let liwork = 1.max(3 + 5 * n);
15        let mut iwork = vec![0; liwork];
16        let mut info = 0;
17
18        let n = n as i32;
19
20        unsafe {
21            dstevd(
22                'V' as u8,
23                n,
24                &mut d,
25                &mut e,
26                &mut z.elems_mut(),
27                n,
28                &mut work,
29                lwork as i32,
30                &mut iwork,
31                liwork as i32,
32                &mut info,
33            )
34        }
35
36        match info {
37            0 => Ok((d, z)),
38            _ => Err(MatrixError::LapackRoutineError {
39                routine: "dstev".to_owned(),
40                info,
41            }),
42        }
43    }
44}