Trait peroxide::structure::matrix::LinearAlgebra [−][src]
pub trait LinearAlgebra {}Show methods
fn back_subs(&self, b: &Vec<f64>) -> Vec<f64>; fn forward_subs(&self, b: &Vec<f64>) -> Vec<f64>; fn lu(&self) -> PQLU; fn waz(&self, d_form: Form) -> Option<WAZD>; fn qr(&self) -> QR; fn svd(&self) -> SVD; fn rref(&self) -> Matrix; fn det(&self) -> f64; fn block(&self) -> (Matrix, Matrix, Matrix, Matrix); fn inv(&self) -> Matrix; fn pseudo_inv(&self) -> Matrix; fn solve(&self, b: &Vec<f64>, sk: SolveKind) -> Vec<f64>; fn solve_mat(&self, m: &Matrix, sk: SolveKind) -> Matrix;
Expand description
Linear algebra trait
Required methods
Implementors
impl LinearAlgebra for SPMatrix
[src]
impl LinearAlgebra for SPMatrix
[src]Linear algebra for sparse matrix
Caution : In every ops in this trait, there is converting process to dense matrix
fn back_subs(&self, _b: &Vec<f64>) -> Vec<f64>
[src]
fn forward_subs(&self, _b: &Vec<f64>) -> Vec<f64>
[src]
fn lu(&self) -> PQLU
[src]
fn waz(&self, _d_form: Form) -> Option<WAZD>
[src]
fn qr(&self) -> QR
[src]
fn det(&self) -> f64
[src]
fn block(&self) -> (Matrix, Matrix, Matrix, Matrix)
[src]
fn inv(&self) -> Matrix
[src]
fn pseudo_inv(&self) -> Matrix
[src]
fn rref(&self) -> Matrix
[src]
fn solve(&self, _b: &Vec<f64>, _sk: SolveKind) -> Vec<f64>
[src]
fn solve_mat(&self, _m: &Matrix, _sk: SolveKind) -> Matrix
[src]
fn svd(&self) -> SVD
[src]
impl LinearAlgebra for Matrix
[src]
impl LinearAlgebra for Matrix
[src]fn lu(&self) -> PQLU
[src]
fn lu(&self) -> PQLU
[src]LU Decomposition Implements (Complete Pivot)
Description
It use complete pivoting LU decomposition. You can get two permutations, and LU matrices.
Caution
It returns Option<PQLU>
- You should unwrap to obtain real value.
PQLU
has four field - p
, q
, l
, u
.
p
, q
are permutations.
l
, u
are matrices.
Examples
#[macro_use] extern crate peroxide; use peroxide::fuga::*; fn main() { let a = matrix(vec![1,2,3,4], 2, 2, Row); let pqlu = a.lu(); let (p,q,l,u) = (pqlu.p, pqlu.q, pqlu.l, pqlu.u); assert_eq!(p, vec![1]); // swap 0 & 1 (Row) assert_eq!(q, vec![1]); // swap 0 & 1 (Col) assert_eq!(l, matrix(c!(1,0,0.5,1),2,2,Row)); assert_eq!(u, matrix(c!(4,3,0,-0.5),2,2,Row)); }
fn qr(&self) -> QR
[src]
fn qr(&self) -> QR
[src]QR Decomposition
Translation of RosettaCode#Python
Example
extern crate peroxide; use peroxide::fuga::*; fn main() { let a = ml_matrix("12 -51 4;6 167 -68; -4 24 -41"); let qr = a.qr(); let r = ml_matrix("-14 -21 14; 0 -175 70; 0 0 -35"); #[cfg(feature="O3")] { assert_eq!(r, qr.r); } qr.r.print(); }
fn svd(&self) -> SVD
[src]
fn svd(&self) -> SVD
[src]Singular Value Decomposition
Examples
extern crate peroxide; use peroxide::fuga::*; fn main() { let a = ml_matrix("3 2 2;2 3 -2"); #[cfg(feature="O3")] { let svd = a.svd(); assert!(eq_vec(&vec![5f64, 3f64], &svd.s, 1e-7)); } a.print(); }
fn rref(&self) -> Matrix
[src]
fn rref(&self) -> Matrix
[src]Reduced Row Echelon Form
Implementation of RosettaCode
fn det(&self) -> f64
[src]
fn det(&self) -> f64
[src]Determinant
Examples
#[macro_use] extern crate peroxide; use peroxide::fuga::*; fn main() { let a = matrix!(1;4;1, 2, 2, Row); assert_eq!(a.det(), -2f64); }
fn block(&self) -> (Self, Self, Self, Self)
[src]
fn block(&self) -> (Self, Self, Self, Self)
[src]Block Partition
Examples
#[macro_use] extern crate peroxide; use peroxide::fuga::*; fn main() { let a = matrix!(1;16;1, 4, 4, Row); let (m1, m2, m3, m4) = a.block(); assert_eq!(m1, matrix(c!(1,2,5,6), 2, 2, Row)); assert_eq!(m2, matrix(c!(3,4,7,8), 2, 2, Row)); assert_eq!(m3, matrix(c!(9,10,13,14), 2, 2, Row)); assert_eq!(m4, matrix(c!(11,12,15,16), 2, 2, Row)); let b = matrix!(1;16;1, 4, 4, Col); let (m1, m2, m3, m4) = b.block(); assert_eq!(m1, matrix(c!(1,2,5,6), 2, 2, Col)); assert_eq!(m3, matrix(c!(3,4,7,8), 2, 2, Col)); assert_eq!(m2, matrix(c!(9,10,13,14), 2, 2, Col)); assert_eq!(m4, matrix(c!(11,12,15,16), 2, 2, Col)); }
fn inv(&self) -> Self
[src]
fn inv(&self) -> Self
[src]Inverse of Matrix
Caution
inv
function returns Option<Matrix>
Thus, you should use pattern matching or unwrap
to obtain inverse.
Examples
#[macro_use] extern crate peroxide; use peroxide::fuga::*; fn main() { // Non-singular let a = matrix!(1;4;1, 2, 2, Row); assert_eq!(a.inv(), matrix(c!(-2,1,1.5,-0.5),2,2,Row)); // Singular //let b = matrix!(1;9;1, 3, 3, Row); //let c = b.inv(); // Can compile but.. //let d = b.det(); //assert_eq!(d, 0f64); }
fn pseudo_inv(&self) -> Self
[src]
fn pseudo_inv(&self) -> Self
[src]Moore-Penrose Pseudo inverse
Description
$X^\dagger = (X^T X)^{-1} X^T$
Examples
#[macro_use] extern crate peroxide; use peroxide::fuga::*; fn main() { let a = matrix!(1;4;1, 2, 2, Row); let inv_a = a.inv(); let pse_a = a.pseudo_inv(); assert_eq!(inv_a, pse_a); // Nearly equal }
fn solve(&self, b: &Vec<f64>, sk: SolveKind) -> Vec<f64>
[src]
fn solve(&self, b: &Vec<f64>, sk: SolveKind) -> Vec<f64>
[src]Solve with Vector
Solve options
- LU: Gaussian elimination with Complete pivoting LU (GECP)
- WAZ: Solve with WAZ decomposition
Reference
- Biswa Nath Datta, Numerical Linear Algebra and Applications, Second Edition
- Ke Chen, Matrix Preconditioning Techniques and Applications, Cambridge Monographs on Applied and Computational Mathematics