lax/
rcond.rs

1use crate::{error::*, layout::MatrixLayout, *};
2use cauchy::*;
3use num_traits::Zero;
4
5pub trait Rcond_: Scalar + Sized {
6    /// Estimates the the reciprocal of the condition number of the matrix in 1-norm.
7    ///
8    /// `anorm` should be the 1-norm of the matrix `a`.
9    fn rcond(l: MatrixLayout, a: &[Self], anorm: Self::Real) -> Result<Self::Real>;
10}
11
12macro_rules! impl_rcond_real {
13    ($scalar:ty, $gecon:path) => {
14        impl Rcond_ for $scalar {
15            fn rcond(l: MatrixLayout, a: &[Self], anorm: Self::Real) -> Result<Self::Real> {
16                let (n, _) = l.size();
17                let mut rcond = Self::Real::zero();
18                let mut info = 0;
19
20                let mut work = unsafe { vec_uninit(4 * n as usize) };
21                let mut iwork = unsafe { vec_uninit(n as usize) };
22                let norm_type = match l {
23                    MatrixLayout::C { .. } => NormType::Infinity,
24                    MatrixLayout::F { .. } => NormType::One,
25                } as u8;
26                unsafe {
27                    $gecon(
28                        norm_type,
29                        n,
30                        a,
31                        l.lda(),
32                        anorm,
33                        &mut rcond,
34                        &mut work,
35                        &mut iwork,
36                        &mut info,
37                    )
38                };
39                info.as_lapack_result()?;
40
41                Ok(rcond)
42            }
43        }
44    };
45}
46
47impl_rcond_real!(f32, lapack::sgecon);
48impl_rcond_real!(f64, lapack::dgecon);
49
50macro_rules! impl_rcond_complex {
51    ($scalar:ty, $gecon:path) => {
52        impl Rcond_ for $scalar {
53            fn rcond(l: MatrixLayout, a: &[Self], anorm: Self::Real) -> Result<Self::Real> {
54                let (n, _) = l.size();
55                let mut rcond = Self::Real::zero();
56                let mut info = 0;
57                let mut work = unsafe { vec_uninit(2 * n as usize) };
58                let mut rwork = unsafe { vec_uninit(2 * n as usize) };
59                let norm_type = match l {
60                    MatrixLayout::C { .. } => NormType::Infinity,
61                    MatrixLayout::F { .. } => NormType::One,
62                } as u8;
63                unsafe {
64                    $gecon(
65                        norm_type,
66                        n,
67                        a,
68                        l.lda(),
69                        anorm,
70                        &mut rcond,
71                        &mut work,
72                        &mut rwork,
73                        &mut info,
74                    )
75                };
76                info.as_lapack_result()?;
77
78                Ok(rcond)
79            }
80        }
81    };
82}
83
84impl_rcond_complex!(c32, lapack::cgecon);
85impl_rcond_complex!(c64, lapack::zgecon);