#matrix #vector #linspace

russell_lab

Matrix-vector laboratory including linear algebra tools

10 releases

new 0.2.4 Oct 22, 2021
0.2.3 Oct 13, 2021
0.2.1 Sep 19, 2021
0.1.4 Sep 9, 2021
0.1.0 Jun 22, 2021

#49 in Math

Download history 1/week @ 2021-07-20 2/week @ 2021-07-27 2/week @ 2021-08-03 1/week @ 2021-08-10 1/week @ 2021-08-17 16/week @ 2021-08-24 19/week @ 2021-08-31 36/week @ 2021-09-07 62/week @ 2021-09-14 9/week @ 2021-09-21 13/week @ 2021-09-28 34/week @ 2021-10-05 46/week @ 2021-10-12 25/week @ 2021-10-19

81 downloads per month
Used in 3 crates

MIT license

275KB
5.5K SLoC

Russell Lab - Matrix-vector laboratory including linear algebra tools

This crate is part of Russell - Rust Scientific Library

This repository is a "rust laboratory" for vectors and matrices.

Documentation:

Installation

Install some libraries:

sudo apt-get install \
    liblapacke-dev \
    libopenblas-dev

Add this to your Cargo.toml (choose the right version):

[dependencies]
russell_lab = "*"

Number of threads

By default OpenBLAS will use all available threads, including Hyper-Threads that make the performance worse. Thus, it is best to set the following environment variable:

export OPENBLAS_NUM_THREADS=<real-core-count>

Furthermore, if working on a multi-threaded application, it is recommended to set:

export OPENBLAS_NUM_THREADS=1

Examples

Compute the pseudo-inverse matrix

use russell_lab::*;

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

    // 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(())
}

Compute eigenvalues

use russell_lab::*;
use russell_chk::*;

fn main() -> Result<(), &'static str> {
    // set matrix
    let data = [
        [2.0, 0.0, 0.0],
        [0.0, 3.0, 4.0],
        [0.0, 4.0, 9.0],
    ];
    let mut a = Matrix::from(&data);

    // allocate output arrays
    let m = a.nrow();
    let mut l_real = vec![0.0; m];
    let mut l_imag = vec![0.0; m];
    let mut v_real = Matrix::new(m, m);
    let mut v_imag = Matrix::new(m, m);

    // perform the eigen-decomposition
    eigen_decomp(
        &mut l_real,
        &mut l_imag,
        &mut v_real,
        &mut v_imag,
        &mut a,
    )?;

    // check results
    let l_real_correct = "[11.0, 1.0, 2.0]";
    let l_imag_correct = "[0.0, 0.0, 0.0]";
    let v_real_correct = "┌                      ┐\n\
                          │  0.000  0.000  1.000 │\n\
                          │  0.447  0.894  0.000 │\n\
                          │  0.894 -0.447  0.000 │\n\
                          └                      ┘";
    let v_imag_correct = "┌       ┐\n\
                          │ 0 0 0 │\n\
                          │ 0 0 0 │\n\
                          │ 0 0 0 │\n\
                          └       ┘";
    assert_eq!(format!("{:?}", l_real), l_real_correct);
    assert_eq!(format!("{:?}", l_imag), l_imag_correct);
    assert_eq!(format!("{:.3}", v_real), v_real_correct);
    assert_eq!(format!("{}", v_imag), v_imag_correct);

    // check eigen-decomposition (similarity transformation) of a
    // symmetric matrix with real-only eigenvalues and eigenvectors
    let a_copy = Matrix::from(&data);
    let lam = Matrix::diagonal(&l_real);
    let mut a_v = Matrix::new(m, m);
    let mut v_l = Matrix::new(m, m);
    let mut err = Matrix::filled(m, m, f64::MAX);
    mat_mat_mul(&mut a_v, 1.0, &a_copy, &v_real)?;
    mat_mat_mul(&mut v_l, 1.0, &v_real, &lam)?;
    add_matrices(&mut err, 1.0, &a_v, -1.0, &v_l)?;
    assert_approx_eq!(err.norm(NormMat::Max), 0.0, 1e-15);
    Ok(())
}

Dependencies

~18MB
~194K SLoC