Linear Algebra

We will use two key crates:

  • nalgebra, a general-purpose linear algebra library with transformations and statically-sized or dynamically-sized matrices. However it supports only vectors (1d) and matrices (2d) and not higher-dimensional tensors.
  • ndarray is less featureful than nalgebra↗ but supports arbitrarily dimensioned arrays.

Add Matrices

ndarray ndarray~crates.io ndarray~github ndarray~lib.rs cat~data-structures cat~science

Creates two 2-D matrices with ndarray::arr2↗ and sums them element-wise.

Note that the sum is computed as let sum = &a + &b. The & operator is used to avoid consuming a and b, making them available later for display. A new array is created containing their sum.

use ndarray::arr2;

/// This example demonstrates how to add two matrices using the `ndarray` crate.
fn main() {
    // Define the first matrix.
    let a = arr2(&[[1, 2, 3], [4, 5, 6]]);

    // Define the second matrix.
    let b = arr2(&[[6, 5, 4], [3, 2, 1]]);

    let sum = &a + &b;

    println!("{a}");
    println!("+");
    println!("{b}");
    println!("=");
    println!("{sum}");
}

Multiply Matrices

ndarray ndarray~crates.io ndarray~github ndarray~lib.rs cat~data-structures cat~science

Creates two matrices with ndarray::arr2↗ and performs matrix multiplication on them with ndarray::ArrayBase::dot↗.

use ndarray::arr2;

/// This example demonstrates matrix multiplication using the `ndarray` crate.
fn main() {
    // Define the first matrix 'a' as a 2x3 array.
    let a = arr2(&[[1, 2, 3], [4, 5, 6]]);

    // Define the second matrix 'b' as a 3x2 array.
    // Note that the number of columns in 'a' must match the number of rows in
    // 'b' for matrix multiplication.
    let b = arr2(&[[6, 3], [5, 2], [4, 1]]);

    println!("{}", a.dot(&b));
}

Multiply a Scalar with a Vector and a Matrix

ndarray ndarray~crates.io ndarray~github ndarray~lib.rs cat~data-structures cat~science

Creates a 1-D array (vector) with ndarray::arr1↗ and a 2-D array (matrix) with ndarray::arr2↗.

First, a scalar is multiplied by the vector to get another vector. Then, the matrix is multiplied by the new vector with ndarray::Array2::dot↗ (Matrix multiplication is performed using ndarray::Array2::dot↗, while the * operator performs element-wise multiplication.)

In ndarray↗, 1-D arrays can be interpreted as either row or column vectors depending on context. If representing the orientation of a vector is important, a 2-D array with one row or one column must be used instead. In this example, the vector is a 1-D array on the right-hand side, so ndarray::Array2::dot↗ handles it as a column vector.

//! This example demonstrates how to multiply a scalar by a vector and then
//! multiply a matrix by the resulting vector.

use ndarray::Array1;
use ndarray::arr1;
use ndarray::arr2;

fn main() {
    let scalar = 4;

    let vector = arr1(&[1, 2, 3]);

    let matrix = arr2(&[[4, 5, 6], [7, 8, 9]]);

    let new_vector: Array1<_> = scalar * vector;
    println!("{new_vector}");

    let new_matrix = matrix.dot(&new_vector);
    println!("{new_matrix}");
}

Compare Vectors

ndarray ndarray~crates.io ndarray~github ndarray~lib.rs cat~data-structures cat~science

The ndarray↗ crate supports a number of ways to create arrays -- this recipe creates ndarray::Array↗ from std::vec::Vec using std::convert::From↗. Then, it sums the arrays element-wise.

This recipe contains an example of comparing two floating-point vectors element-wise. Floating-point numbers are often stored inexactly, making exact comparisons difficult. However, the approx::assert_abs_diff_eq↗ macro from the approx↗ crate allows for convenient element-wise comparisons. To use the approx↗ crate with ndarray↗, the approx↗ feature must be added to the ndarray↗ dependency in Cargo.toml. For example, ndarray = { version = "0.15.6", features = [ "approx" ] }.

This recipe also contains additional ownership examples. Here, let z = a + b consumes a and b, updates a with the result, then moves ownership to z. Alternatively, let w = &c + &d creates a new vector without consuming c or d, allowing their modification later. See Binary Operators With Two Arrays↗ for additional detail.

//! This example demonstrates how to compare vectors using the `approx` crate.

use approx::assert_abs_diff_eq;
use ndarray::Array;

fn main() {
    // Create vectors.
    let a = Array::from(vec![1., 2., 3., 4., 5.]);
    let b = Array::from(vec![5., 4., 3., 2., 1.]);
    let mut c = Array::from(vec![1., 2., 3., 4., 5.]);
    let mut d = Array::from(vec![5., 4., 3., 2., 1.]);

    let z = a + b;
    let w = &c + &d;

    // Check that the sum of the vectors is correct.
    assert_abs_diff_eq!(z, Array::from(vec![6., 6., 6., 6., 6.]));

    println!("c = {c}");

    // Modify the vectors.
    c[0] = 10.;
    d[1] = 10.;

    // Assert approximate equality (using the absolute difference).
    assert_abs_diff_eq!(w, Array::from(vec![6., 6., 6., 6., 6.]));
}

Calculate Vector Norms

ndarray ndarray~crates.io ndarray~github ndarray~lib.rs cat~data-structures cat~science

This recipe demonstrates use of the ndarray::Array1↗ type, ndarray::Array1↗ type, ndarray::ArrayBase::fold method, and ndarray::ArrayBase::dot↗ method in computing the l1↗ and l2↗ norms of a given vector.

The l2_norm↗ function is the simpler of the two, as it computes the square root of the dot product of a vector with itself. The l1_norm↗ function is computed by a ndarray::ArrayBase::fold↗ operation that sums the absolute values of the elements. (This could also be performed with x.mapv(f64::abs).scalar_sum(), but that would allocate a new array for the result of the mapv.)

Note that both l1_norm↗ and l2_norm↗ take the ndarray::ArrayView1↗ type. This recipe considers vector norms, so the norm functions only need to accept one-dimensional views, hence ndarray::ArrayView1↗. While the functions could take a parameter of type &Array1<f64> instead, that would require the caller to have a reference to an owned array, which is more restrictive than just having access to a view (since a view can be created from any array or view, not just an owned array).

ndarray::Array↗ and ndarray::ArrayView↗ are both type aliases for ndarray::ArrayBase↗. So, the most general argument type for the caller would be &ArrayBase<S, Ix1> where S: Data, because then the caller could use &array or &view instead of x.view(). If the function is part of a public API, that may be a better choice for the benefit of users. For internal functions, the more concise ArrayView1<f64> may be preferable.

use ndarray::Array1;
use ndarray::ArrayView1;
use ndarray::array;

/// Calculates the L1 norm (Manhattan norm) of a vector.
///
/// The L1 norm is the sum of the absolute values of the vector's elements.
fn l1_norm(x: ArrayView1<f64>) -> f64 {
    x.fold(0., |acc, elem| acc + elem.abs())
}

/// Calculates the L2 norm (Euclidean norm) of a vector.
///
/// The L2 norm is the square root of the sum of the squares of the vector's
/// elements.
fn l2_norm(x: ArrayView1<f64>) -> f64 {
    x.dot(&x).sqrt()
}

/// Normalizes a vector to have a unit L2 norm.
fn normalize(mut x: Array1<f64>) -> Array1<f64> {
    let norm = l2_norm(x.view());
    x.mapv_inplace(|e| e / norm);
    x
}

fn main() {
    let x = array![1., 2., 3., 4., 5.];
    println!("||x||_2 = {}", l2_norm(x.view()));
    println!("||x||_1 = {}", l1_norm(x.view()));
    println!("Normalizing x yields {:?}", normalize(x));
}

Invert a Matrix

nalgebra nalgebra~crates.io nalgebra~github nalgebra~lib.rs cat~mathematics cat~no-std cat~science cat~wasm

Creates a 3x3 matrix with nalgebra::Matrix3↗ and inverts it, if possible.

use nalgebra::Matrix3;

/// This example demonstrates how to invert a 3x3 matrix using the `nalgebra`
/// crate.
fn main() {
    // Create a 3x3 matrix.
    let m1 = Matrix3::new(2.0, 1.0, 1.0, 3.0, 2.0, 1.0, 2.0, 1.0, 2.0);
    // Print the matrix.
    println!("m1 = {m1}");
    // Try to invert the matrix.
    match m1.try_inverse() {
        Some(inv) => {
            println!("The inverse of m1 is: {inv}");
        }
        None => {
            println!("m1 is not invertible!");
        }
    }
}

(De)serialize a Matrix

ndarray ndarray~crates.io ndarray~github ndarray~lib.rs cat~data-structures cat~science

Serialize and deserialize a matrix to and from JSON. Serialization is taken care of by serde_json::to_string↗ and serde_json::to_string↗ performs deserialization.

Note that serialization followed by deserialization gives back the original matrix.

use nalgebra::DMatrix;

/// This example demonstrates how to serialize and deserialize a `DMatrix` using
/// `serde_json`.
///
/// It creates a 50x100 matrix, serializes it to a JSON string, deserializes it
/// back, and then verifies that the deserialized matrix is equal to the
/// original matrix.
fn main() -> Result<(), std::io::Error> {
    let row_slice: Vec<i32> = (1..5001).collect();
    let matrix = DMatrix::from_row_slice(50, 100, &row_slice);
    println!("{matrix}");

    // Serialize the matrix.
    let serialized_matrix = serde_json::to_string(&matrix)?;

    // Deserialize the matrix.
    let deserialized_matrix: DMatrix<i32> =
        serde_json::from_str(&serialized_matrix)?;

    // Verify that `deserialized_matrix` is equal to `matrix`.
    assert!(deserialized_matrix == matrix);

    Ok(())
}