ndarray_linalg/
svddc.rs

1//! Singular-value decomposition (SVD) by divide-and-conquer (?gesdd)
2
3use super::{convert::*, error::*, layout::*, types::*};
4use ndarray::*;
5
6pub use lax::JobSvd;
7
8/// Singular-value decomposition of matrix (copying) by divide-and-conquer
9pub trait SVDDC {
10    type U;
11    type VT;
12    type Sigma;
13    fn svddc(&self, uvt_flag: JobSvd) -> Result<(Option<Self::U>, Self::Sigma, Option<Self::VT>)>;
14}
15
16/// Singular-value decomposition of matrix by divide-and-conquer
17pub trait SVDDCInto {
18    type U;
19    type VT;
20    type Sigma;
21    fn svddc_into(
22        self,
23        uvt_flag: JobSvd,
24    ) -> Result<(Option<Self::U>, Self::Sigma, Option<Self::VT>)>;
25}
26
27/// Singular-value decomposition of matrix reference by divide-and-conquer
28pub trait SVDDCInplace {
29    type U;
30    type VT;
31    type Sigma;
32    fn svddc_inplace(
33        &mut self,
34        uvt_flag: JobSvd,
35    ) -> Result<(Option<Self::U>, Self::Sigma, Option<Self::VT>)>;
36}
37
38impl<A, S> SVDDC for ArrayBase<S, Ix2>
39where
40    A: Scalar + Lapack,
41    S: Data<Elem = A>,
42{
43    type U = Array2<A>;
44    type VT = Array2<A>;
45    type Sigma = Array1<A::Real>;
46
47    fn svddc(&self, uvt_flag: JobSvd) -> Result<(Option<Self::U>, Self::Sigma, Option<Self::VT>)> {
48        self.to_owned().svddc_into(uvt_flag)
49    }
50}
51
52impl<A, S> SVDDCInto for ArrayBase<S, Ix2>
53where
54    A: Scalar + Lapack,
55    S: DataMut<Elem = A>,
56{
57    type U = Array2<A>;
58    type VT = Array2<A>;
59    type Sigma = Array1<A::Real>;
60
61    fn svddc_into(
62        mut self,
63        uvt_flag: JobSvd,
64    ) -> Result<(Option<Self::U>, Self::Sigma, Option<Self::VT>)> {
65        self.svddc_inplace(uvt_flag)
66    }
67}
68
69impl<A, S> SVDDCInplace for ArrayBase<S, Ix2>
70where
71    A: Scalar + Lapack,
72    S: DataMut<Elem = A>,
73{
74    type U = Array2<A>;
75    type VT = Array2<A>;
76    type Sigma = Array1<A::Real>;
77
78    fn svddc_inplace(
79        &mut self,
80        uvt_flag: JobSvd,
81    ) -> Result<(Option<Self::U>, Self::Sigma, Option<Self::VT>)> {
82        let l = self.layout()?;
83        let svd_res = A::svddc(l, uvt_flag, self.as_allocated_mut()?)?;
84        let (m, n) = l.size();
85        let k = m.min(n);
86
87        let (u_col, vt_row) = match uvt_flag {
88            JobSvd::All => (m, n),
89            JobSvd::Some => (k, k),
90            JobSvd::None => (0, 0),
91        };
92
93        let u = svd_res
94            .u
95            .map(|u| into_matrix(l.resized(m, u_col), u).unwrap());
96
97        let vt = svd_res
98            .vt
99            .map(|vt| into_matrix(l.resized(vt_row, n), vt).unwrap());
100
101        let s = ArrayBase::from(svd_res.s);
102        Ok((u, s, vt))
103    }
104}