pub fn pseudo_inverse(ai: &mut Matrix, a: &Matrix) -> Result<(), StrError>
Expand description

Computes the pseudo-inverse matrix

Finds ai such that:

a⋅ai⋅a == a
  ai  :=  pseudo_inverse(a)
(n,m)          (m,n)

Output

  • ai – (n,m) pseudo-inverse matrix
  • If a is invertible, its pseudo-inverse equals its inverse

Input

  • a – (m,n) matrix, symmetric or not

Example

use russell_lab::{pseudo_inverse, Matrix, StrError};

fn main() -> Result<(), StrError> {
    // set matrix
    let mut a = Matrix::from(&[
        [1.0, 0.0],
        [0.0, 1.0],
        [0.0, 1.0],
    ]);
    let a_copy = a.clone();

    // compute pseudo-inverse matrix (because it's square)
    let mut ai = Matrix::new(2, 3);
    pseudo_inverse(&mut ai, &mut a)?;

    // compare with solution
    let ai_correct = "┌                ┐\n\
                      │ 1.00 0.00 0.00 │\n\
                      │ 0.00 0.50 0.50 │\n\
                      └                ┘";
    assert_eq!(format!("{:.2}", ai), ai_correct);

    // compute a⋅ai
    let (m, n) = a.dims();
    let mut a_ai = Matrix::new(m, m);
    for i in 0..m {
        for j in 0..m {
            for k in 0..n {
                a_ai[i][j] += a_copy[i][k] * ai[k][j];
            }
        }
    }

    // check if a⋅ai⋅a == a
    let mut a_ai_a = Matrix::new(m, n);
    for i in 0..m {
        for j in 0..n {
            for k in 0..m {
                a_ai_a[i][j] += a_ai[i][k] * a_copy[k][j];
            }
        }
    }
    let a_ai_a_correct = "┌           ┐\n\
                          │ 1.00 0.00 │\n\
                          │ 0.00 1.00 │\n\
                          │ 0.00 1.00 │\n\
                          └           ┘";
    assert_eq!(format!("{:.2}", a_ai_a), a_ai_a_correct);
    Ok(())
}