opensrdk-linear-algebra 0.8.5

Standard linear algebra library using blas and lapack for OpenSRDK toolchain.
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

use super::TridiagonalMatrix;
use crate::matrix::MatrixError;
use crate::number::*;
use lapack::dgttrf;
use lapack::zgttrf;

#[derive(Clone, Debug)]
pub struct GTTRF<T = f64>(pub Vec<T>, pub [Vec<T>; 3], pub Vec<i32>)
where
    T: Number;

impl TridiagonalMatrix {
    /// # LU decomposition
    /// for tridiagonal matrix
    pub fn gttrf(self) -> Result<GTTRF, MatrixError> {
        let (mut dl, mut d, mut du) = self.eject();
        let n = d.len();
        let mut du2 = vec![0.0; n.max(2) - 2];
        let mut ipiv = vec![0; n];
        let mut info = 0;

        let n = n as i32;

        unsafe { dgttrf(n, &mut dl, &mut d, &mut du, &mut du2, &mut ipiv, &mut info) }

        if info != 0 {
            return Err(MatrixError::LapackRoutineError {
                routine: "dgttrf".to_owned(),
                info,
            });
        }

        let u = [d, du, du2];

        Ok(GTTRF(dl, u, ipiv))
    }
}

impl TridiagonalMatrix<c64> {
    /// # LU decomposition
    /// for tridiagonal matrix
    pub fn gttrf(self) -> Result<GTTRF<c64>, MatrixError> {
        let (mut dl, mut d, mut du) = self.eject();
        let n = d.len();
        let mut du2 = vec![c64::default(); n.max(2) - 2];
        let mut ipiv = vec![0; n];
        let mut info = 0;

        let n = n as i32;

        unsafe { zgttrf(n, &mut dl, &mut d, &mut du, &mut du2, &mut ipiv, &mut info) }

        if info != 0 {
            return Err(MatrixError::LapackRoutineError {
                routine: "dgttrf".to_owned(),
                info,
            });
        }

        let u = [d, du, du2];

        Ok(GTTRF::<c64>(dl, u, ipiv))
    }
}