mdarray-linalg-faer 0.1.2

Faer backend for mdarray-linalg
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

use dyn_stack::{MemBuffer, MemStack};
use faer_traits::ComplexField;
use mdarray::{DSlice, Layout};
use mdarray_linalg::svd::SVDError;
use num_complex::ComplexFloat;

use crate::{into_faer, into_faer_diag_mut, into_faer_mut, into_faer_mut_transpose};

pub fn svd_faer<
    T: ComplexFloat + ComplexField + Default + 'static,
    La: Layout,
    Ls: Layout,
    Lu: Layout,
    Lvt: Layout,
>(
    a: &DSlice<T, 2, La>,
    s_mda: &mut DSlice<T, 2, Ls>,
    u_mda: Option<&mut DSlice<T, 2, Lu>>,
    vt_mda: Option<&mut DSlice<T, 2, Lvt>>,
) -> Result<(), SVDError> {
    let (m, n) = *a.shape();
    let a_faer = into_faer(a);
    let par = faer::get_global_parallelism();
    // let par = faer::Par::Seq; // Faster for small matrices

    match (u_mda, vt_mda) {
        (Some(x), Some(y)) => {
            let mut s_faer = into_faer_diag_mut(s_mda);
            let u_faer = into_faer_mut(x);
            let vt_faer = into_faer_mut_transpose(y);

            let ret = faer::linalg::svd::svd(
                a_faer,
                s_faer.as_mut(),
                Some(u_faer),
                Some(vt_faer),
                par,
                MemStack::new(&mut MemBuffer::new(faer::linalg::svd::svd_scratch::<T>(
                    m,
                    n,
                    faer::linalg::svd::ComputeSvdVectors::Full,
                    faer::linalg::svd::ComputeSvdVectors::Full,
                    par,
                    faer::prelude::default(),
                ))),
                faer::prelude::default(),
            );
            match ret {
                Ok(()) => Ok(()),
                Err(_) => Err(SVDError::BackendDidNotConverge {
                    superdiagonals: (0),
                }),
            }
        }
        (None, None) => {
            let mut s_faer = into_faer_diag_mut(s_mda);
            let ret = faer::linalg::svd::svd(
                a_faer,
                s_faer.as_mut(),
                None,
                None,
                par,
                MemStack::new(&mut MemBuffer::new(faer::linalg::svd::svd_scratch::<T>(
                    m,
                    n,
                    faer::linalg::svd::ComputeSvdVectors::No,
                    faer::linalg::svd::ComputeSvdVectors::No,
                    par,
                    faer::prelude::default(),
                ))),
                faer::prelude::default(),
            );
            match ret {
                Ok(()) => Ok(()),
                Err(_) => Err(SVDError::BackendDidNotConverge {
                    superdiagonals: (0),
                }),
            }
        }
        _ => Err(SVDError::InconsistentUV),
    }
}