svd3.rs - source
Description: Source of the Rust file `/Users/erlendbasso/.cargo/registry/src/github.com-1ecc6299db9ec823/nalgebra-0.32.1/src/linalg/svd3.rs`.
Keywords: rust, rustlang, rust-lang

    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

    use crate::{Matrix3, SVD, U3};
    use simba::scalar::RealField;

    // For the 3x3 case, on the GPU, it is much more efficient to compute the SVD
    // using an eigendecomposition followed by a QR decomposition.
    //
    // This is based on the paper "Computing the Singular Value Decomposition of 3 x 3 matrices with
    // minimal branching and elementary floating point operations" from McAdams, et al.
    pub fn svd_ordered3<T: RealField>(
        m: &Matrix3<T>,
        compute_u: bool,
        compute_v: bool,
        eps: T,
        niter: usize,
    ) -> Option<SVD<T, U3, U3>> {
        let s = m.tr_mul(&m);
        let mut v = s.try_symmetric_eigen(eps, niter)?.eigenvectors;
        let mut b = m * &v;

        // Sort singular values. This is a necessary step to ensure that
        // the QR decompositions R matrix ends up diagonal.
        let mut rho0 = b.column(0).norm_squared();
        let mut rho1 = b.column(1).norm_squared();
        let mut rho2 = b.column(2).norm_squared();

        if rho0 < rho1 {
            b.swap_columns(0, 1);
            b.column_mut(1).neg_mut();
            v.swap_columns(0, 1);
            v.column_mut(1).neg_mut();
            std::mem::swap(&mut rho0, &mut rho1);
        }
        if rho0 < rho2 {
            b.swap_columns(0, 2);
            b.column_mut(2).neg_mut();
            v.swap_columns(0, 2);
            v.column_mut(2).neg_mut();
            std::mem::swap(&mut rho0, &mut rho2);
        }
        if rho1 < rho2 {
            b.swap_columns(1, 2);
            b.column_mut(2).neg_mut();
            v.swap_columns(1, 2);
            v.column_mut(2).neg_mut();
            std::mem::swap(&mut rho0, &mut rho2);
        }

        let qr = b.qr();

        Some(SVD {
            u: if compute_u { Some(qr.q()) } else { None },
            singular_values: qr.diag_internal().map(|e| e.abs()),
            v_t: if compute_v { Some(v.transpose()) } else { None },
        })
    }