pub struct Matrix<T> { /* private fields */ }
Implementations§
source§impl<T: Clone + Number> Matrix<T>
impl<T: Clone + Number> Matrix<T>
sourcepub fn set_col(&mut self, col: usize, vec: Vector<T>)
pub fn set_col(&mut self, col: usize, vec: Vector<T>)
Set a column of the matrix using a vector
sourcepub fn delete_row(&mut self, row: usize)
pub fn delete_row(&mut self, row: usize)
Delete a row from the matrix
sourcepub fn multiply(&self, vec: Vector<T>) -> Vector<T>
pub fn multiply(&self, vec: Vector<T>) -> Vector<T>
Multiply the matrix by a (column) vector and return a vector
sourcepub fn resize(&mut self, n_rows: usize, n_cols: usize)
pub fn resize(&mut self, n_rows: usize, n_cols: usize)
Resize the matrix (empty entries are appended if necessary)
sourcepub fn transpose_in_place(&mut self)
pub fn transpose_in_place(&mut self)
Transpose the matrix in place
sourcepub fn swap_elem(
&mut self,
row_1: usize,
col_1: usize,
row_2: usize,
col_2: usize
)
pub fn swap_elem( &mut self, row_1: usize, col_1: usize, row_2: usize, col_2: usize )
Swap two elements of the matrix
sourcepub fn fill_diag(&mut self, elem: T)
pub fn fill_diag(&mut self, elem: T)
Fill the leading diagonal of the matrix with specified elements
sourcepub fn fill_band(&mut self, offset: isize, elem: T)
pub fn fill_band(&mut self, offset: isize, elem: T)
Fill a diagonal band of the matrix with specified elements ( offset above main diagonal +, below main diagonal - )
sourcepub fn fill_tridiag(&mut self, lower: T, diag: T, upper: T)
pub fn fill_tridiag(&mut self, lower: T, diag: T, upper: T)
Fill the main three diagonals of the matrix with specified elements
source§impl<T: Clone + Copy + Number + Signed + PartialOrd> Matrix<T>
impl<T: Clone + Copy + Number + Signed + PartialOrd> Matrix<T>
sourcepub fn solve_basic(&mut self, b: Vector<T>) -> Vector<T>
pub fn solve_basic(&mut self, b: Vector<T>) -> Vector<T>
Solve the system of equations Ax=b where b is a specified vector using Gaussian elimination (the matrix A is modified in the process)
Examples found in repository?
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
fn main() {
println!( "------------------ Sparse linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let n: usize = 4;
let mut a = Mat64::new( n, n, 0.0 );
let mut triplets = Vec::new();
for i in 0..n {
a[i][i] = 2.0;
triplets.push( ( i, i, 2.0 ) );
if i > 0 {
a[ i - 1 ][ i ] = 1.0;
triplets.push( ( i - 1, i, 1.0 ) );
}
if i < n - 1 {
a[ i + 1 ][ i ] = 1.0;
triplets.push( ( i + 1, i, 1.0 ) );
}
}
println!( " * A ={}", a );
let b = Vec64::random( n );
println!( " * b^T ={}", b );
println!( " * as both dense and sparse linear systems.");
let dense_x = a.clone().solve_basic( b.clone() );
println!( " * The dense system gives the solution vector");
println!( " * x^T ={}", dense_x );
let mut sparse = Sparse::from_triplets( triplets );
let sparse_x = sparse.solve( b ).unwrap();
println!( " * The sparse system gives the solution vector");
println!( " * x^T ={}", sparse_x );
let diff = dense_x - sparse_x;
println!( " * The maximum difference between the two is {:+.2e}", diff.norm_inf() );
println!( "-------------------- FINISHED ---------------------" );
}
More examples
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
fn main() {
println!( "------------------ Linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let mut a = Mat64::new( 3, 3, 0.0 );
a[0][0] = 1.0; a[0][1] = 1.0; a[0][2] = 1.0;
a[1][0] = 0.0; a[1][1] = 2.0; a[1][2] = 5.0;
a[2][0] = 2.0; a[2][1] = 5.0; a[2][2] = - 1.0;
println!( " * A ={}", a );
let b = Vec64::create( vec![ 6.0, -4.0, 27.0 ] );
println!( " * b^T ={}", b );
let x = a.clone().solve_basic( b );
println!( " * gives the solution vector");
println!( " * x^T ={}", x );
println!( "---------------------------------------------------" );
println!( " * Lets solve the complex linear system Cy = d where" );
let mut c = Matrix::<Cmplx>::new( 2, 2, Cmplx::new( 1.0, 1.0 ) );
c[0][1] = Cmplx::new( -1.0, 0.0 );
c[1][0] = Cmplx::new( 1.0, -1.0 );
println!( " * C ={}", c );
let mut d = Vector::<Cmplx>::new( 2, Cmplx::new( 0.0, 1.0 ) );
d[1] = Cmplx::new( 1.0, 0.0 );
println!( " * d^T ={}", d );
println!( " * Using LU decomposition we find that ");
let y = c.clone().solve_lu( d );
println!( " * y^T ={}", y );
println!( "---------------------------------------------------" );
println!( " * We may also find the determinant of a matrix" );
println!( " * |A| = {}", a.clone().determinant() );
println!( " * |C| = {}", c.clone().determinant() );
println!( " * or the inverse of a matrix" );
let inverse = a.clone().inverse();
println!( " * A^-1 =\n{}", inverse );
println!( " * C^-1 =\n{}", c.clone().inverse() );
println!( " * We can check this by multiplication" );
println!( " * I = A * A^-1 =\n{}", a * inverse );
println!( " * which is sufficiently accurate for our purposes.");
println!( "-------------------- FINISHED ---------------------" );
}
sourcepub fn lu_decomp_in_place(&mut self) -> (usize, Matrix<T>)
pub fn lu_decomp_in_place(&mut self) -> (usize, Matrix<T>)
Replace the matrix with its LU decomposition and return the number of pivots and a permutation matrix
sourcepub fn solve_lu(&mut self, b: Vector<T>) -> Vector<T>
pub fn solve_lu(&mut self, b: Vector<T>) -> Vector<T>
Solve the system of equations Ax=b where b is a specified vector using LU decomposition (the matrix A is modified in the process)
Examples found in repository?
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
fn main() {
println!( "------------------ Linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let mut a = Mat64::new( 3, 3, 0.0 );
a[0][0] = 1.0; a[0][1] = 1.0; a[0][2] = 1.0;
a[1][0] = 0.0; a[1][1] = 2.0; a[1][2] = 5.0;
a[2][0] = 2.0; a[2][1] = 5.0; a[2][2] = - 1.0;
println!( " * A ={}", a );
let b = Vec64::create( vec![ 6.0, -4.0, 27.0 ] );
println!( " * b^T ={}", b );
let x = a.clone().solve_basic( b );
println!( " * gives the solution vector");
println!( " * x^T ={}", x );
println!( "---------------------------------------------------" );
println!( " * Lets solve the complex linear system Cy = d where" );
let mut c = Matrix::<Cmplx>::new( 2, 2, Cmplx::new( 1.0, 1.0 ) );
c[0][1] = Cmplx::new( -1.0, 0.0 );
c[1][0] = Cmplx::new( 1.0, -1.0 );
println!( " * C ={}", c );
let mut d = Vector::<Cmplx>::new( 2, Cmplx::new( 0.0, 1.0 ) );
d[1] = Cmplx::new( 1.0, 0.0 );
println!( " * d^T ={}", d );
println!( " * Using LU decomposition we find that ");
let y = c.clone().solve_lu( d );
println!( " * y^T ={}", y );
println!( "---------------------------------------------------" );
println!( " * We may also find the determinant of a matrix" );
println!( " * |A| = {}", a.clone().determinant() );
println!( " * |C| = {}", c.clone().determinant() );
println!( " * or the inverse of a matrix" );
let inverse = a.clone().inverse();
println!( " * A^-1 =\n{}", inverse );
println!( " * C^-1 =\n{}", c.clone().inverse() );
println!( " * We can check this by multiplication" );
println!( " * I = A * A^-1 =\n{}", a * inverse );
println!( " * which is sufficiently accurate for our purposes.");
println!( "-------------------- FINISHED ---------------------" );
}
sourcepub fn determinant(&self) -> T
pub fn determinant(&self) -> T
Calculate the determinant of the matrix ( via LU decomposition )
Examples found in repository?
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
fn main() {
println!( "------------------ Linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let mut a = Mat64::new( 3, 3, 0.0 );
a[0][0] = 1.0; a[0][1] = 1.0; a[0][2] = 1.0;
a[1][0] = 0.0; a[1][1] = 2.0; a[1][2] = 5.0;
a[2][0] = 2.0; a[2][1] = 5.0; a[2][2] = - 1.0;
println!( " * A ={}", a );
let b = Vec64::create( vec![ 6.0, -4.0, 27.0 ] );
println!( " * b^T ={}", b );
let x = a.clone().solve_basic( b );
println!( " * gives the solution vector");
println!( " * x^T ={}", x );
println!( "---------------------------------------------------" );
println!( " * Lets solve the complex linear system Cy = d where" );
let mut c = Matrix::<Cmplx>::new( 2, 2, Cmplx::new( 1.0, 1.0 ) );
c[0][1] = Cmplx::new( -1.0, 0.0 );
c[1][0] = Cmplx::new( 1.0, -1.0 );
println!( " * C ={}", c );
let mut d = Vector::<Cmplx>::new( 2, Cmplx::new( 0.0, 1.0 ) );
d[1] = Cmplx::new( 1.0, 0.0 );
println!( " * d^T ={}", d );
println!( " * Using LU decomposition we find that ");
let y = c.clone().solve_lu( d );
println!( " * y^T ={}", y );
println!( "---------------------------------------------------" );
println!( " * We may also find the determinant of a matrix" );
println!( " * |A| = {}", a.clone().determinant() );
println!( " * |C| = {}", c.clone().determinant() );
println!( " * or the inverse of a matrix" );
let inverse = a.clone().inverse();
println!( " * A^-1 =\n{}", inverse );
println!( " * C^-1 =\n{}", c.clone().inverse() );
println!( " * We can check this by multiplication" );
println!( " * I = A * A^-1 =\n{}", a * inverse );
println!( " * which is sufficiently accurate for our purposes.");
println!( "-------------------- FINISHED ---------------------" );
}
sourcepub fn inverse(&self) -> Matrix<T>
pub fn inverse(&self) -> Matrix<T>
Return the inverse of the matrix ( via LU decomposition )
Examples found in repository?
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
fn main() {
println!( "------------------ Linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let mut a = Mat64::new( 3, 3, 0.0 );
a[0][0] = 1.0; a[0][1] = 1.0; a[0][2] = 1.0;
a[1][0] = 0.0; a[1][1] = 2.0; a[1][2] = 5.0;
a[2][0] = 2.0; a[2][1] = 5.0; a[2][2] = - 1.0;
println!( " * A ={}", a );
let b = Vec64::create( vec![ 6.0, -4.0, 27.0 ] );
println!( " * b^T ={}", b );
let x = a.clone().solve_basic( b );
println!( " * gives the solution vector");
println!( " * x^T ={}", x );
println!( "---------------------------------------------------" );
println!( " * Lets solve the complex linear system Cy = d where" );
let mut c = Matrix::<Cmplx>::new( 2, 2, Cmplx::new( 1.0, 1.0 ) );
c[0][1] = Cmplx::new( -1.0, 0.0 );
c[1][0] = Cmplx::new( 1.0, -1.0 );
println!( " * C ={}", c );
let mut d = Vector::<Cmplx>::new( 2, Cmplx::new( 0.0, 1.0 ) );
d[1] = Cmplx::new( 1.0, 0.0 );
println!( " * d^T ={}", d );
println!( " * Using LU decomposition we find that ");
let y = c.clone().solve_lu( d );
println!( " * y^T ={}", y );
println!( "---------------------------------------------------" );
println!( " * We may also find the determinant of a matrix" );
println!( " * |A| = {}", a.clone().determinant() );
println!( " * |C| = {}", c.clone().determinant() );
println!( " * or the inverse of a matrix" );
let inverse = a.clone().inverse();
println!( " * A^-1 =\n{}", inverse );
println!( " * C^-1 =\n{}", c.clone().inverse() );
println!( " * We can check this by multiplication" );
println!( " * I = A * A^-1 =\n{}", a * inverse );
println!( " * which is sufficiently accurate for our purposes.");
println!( "-------------------- FINISHED ---------------------" );
}
source§impl Matrix<f64>
impl Matrix<f64>
source§impl<T: Clone + Number> Matrix<T>
impl<T: Clone + Number> Matrix<T>
sourcepub fn new(rows: usize, cols: usize, elem: T) -> Self
pub fn new(rows: usize, cols: usize, elem: T) -> Self
Create a new matrix of specified size
Examples found in repository?
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
fn main() {
println!( "------------------ Sparse linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let n: usize = 4;
let mut a = Mat64::new( n, n, 0.0 );
let mut triplets = Vec::new();
for i in 0..n {
a[i][i] = 2.0;
triplets.push( ( i, i, 2.0 ) );
if i > 0 {
a[ i - 1 ][ i ] = 1.0;
triplets.push( ( i - 1, i, 1.0 ) );
}
if i < n - 1 {
a[ i + 1 ][ i ] = 1.0;
triplets.push( ( i + 1, i, 1.0 ) );
}
}
println!( " * A ={}", a );
let b = Vec64::random( n );
println!( " * b^T ={}", b );
println!( " * as both dense and sparse linear systems.");
let dense_x = a.clone().solve_basic( b.clone() );
println!( " * The dense system gives the solution vector");
println!( " * x^T ={}", dense_x );
let mut sparse = Sparse::from_triplets( triplets );
let sparse_x = sparse.solve( b ).unwrap();
println!( " * The sparse system gives the solution vector");
println!( " * x^T ={}", sparse_x );
let diff = dense_x - sparse_x;
println!( " * The maximum difference between the two is {:+.2e}", diff.norm_inf() );
println!( "-------------------- FINISHED ---------------------" );
}
More examples
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
fn main() {
println!( "------------------ Linear system ------------------" );
println!( " * Solving the linear system Ax = b, for x, where" );
let mut a = Mat64::new( 3, 3, 0.0 );
a[0][0] = 1.0; a[0][1] = 1.0; a[0][2] = 1.0;
a[1][0] = 0.0; a[1][1] = 2.0; a[1][2] = 5.0;
a[2][0] = 2.0; a[2][1] = 5.0; a[2][2] = - 1.0;
println!( " * A ={}", a );
let b = Vec64::create( vec![ 6.0, -4.0, 27.0 ] );
println!( " * b^T ={}", b );
let x = a.clone().solve_basic( b );
println!( " * gives the solution vector");
println!( " * x^T ={}", x );
println!( "---------------------------------------------------" );
println!( " * Lets solve the complex linear system Cy = d where" );
let mut c = Matrix::<Cmplx>::new( 2, 2, Cmplx::new( 1.0, 1.0 ) );
c[0][1] = Cmplx::new( -1.0, 0.0 );
c[1][0] = Cmplx::new( 1.0, -1.0 );
println!( " * C ={}", c );
let mut d = Vector::<Cmplx>::new( 2, Cmplx::new( 0.0, 1.0 ) );
d[1] = Cmplx::new( 1.0, 0.0 );
println!( " * d^T ={}", d );
println!( " * Using LU decomposition we find that ");
let y = c.clone().solve_lu( d );
println!( " * y^T ={}", y );
println!( "---------------------------------------------------" );
println!( " * We may also find the determinant of a matrix" );
println!( " * |A| = {}", a.clone().determinant() );
println!( " * |C| = {}", c.clone().determinant() );
println!( " * or the inverse of a matrix" );
let inverse = a.clone().inverse();
println!( " * A^-1 =\n{}", inverse );
println!( " * C^-1 =\n{}", c.clone().inverse() );
println!( " * We can check this by multiplication" );
println!( " * I = A * A^-1 =\n{}", a * inverse );
println!( " * which is sufficiently accurate for our purposes.");
println!( "-------------------- FINISHED ---------------------" );
}
Trait Implementations§
source§impl<T: Clone + Number> AddAssign<T> for Matrix<T>
impl<T: Clone + Number> AddAssign<T> for Matrix<T>
source§fn add_assign(&mut self, rhs: T)
fn add_assign(&mut self, rhs: T)
Add the same value to every element in a mutable matrix
source§impl<T: Clone + Number> AddAssign for Matrix<T>
impl<T: Clone + Number> AddAssign for Matrix<T>
source§fn add_assign(&mut self, rhs: Self)
fn add_assign(&mut self, rhs: Self)
Add a matrix to a mutable matrix and assign the result ( += )
source§impl<T: Clone + Number> DivAssign<T> for Matrix<T>
impl<T: Clone + Number> DivAssign<T> for Matrix<T>
source§fn div_assign(&mut self, rhs: T)
fn div_assign(&mut self, rhs: T)
Divide a mutable matrix by a scalar (matrix /= scalar)
source§impl<T: Clone + Number> MulAssign<T> for Matrix<T>
impl<T: Clone + Number> MulAssign<T> for Matrix<T>
source§fn mul_assign(&mut self, rhs: T)
fn mul_assign(&mut self, rhs: T)
Multiply a mutable matrix by a scalar (matrix *= scalar)