M-Labs
/
nalgebra
3
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Added support for dynamic matrices

This commit is contained in:
Fredrik Jansson 2020-04-12 11:59:06 +02:00
parent 0a3ee99cdb
commit e914afe2af
2 changed files with 42 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

23
src/linalg/exp.rs
View File

@ -3,7 +3,8 @@
use crate::{ use crate::{
base::{ base::{
allocator::Allocator, allocator::Allocator,
dimension::{DimMin, DimMinimum, DimName}, dimension::{Dim, DimMin, DimMinimum, U1},
storage::Storage,
DefaultAllocator, DefaultAllocator,
}, },
convert, try_convert, ComplexField, MatrixN, RealField, convert, try_convert, ComplexField, MatrixN, RealField,
@ -13,7 +14,7 @@ use crate::{
struct ExpmPadeHelper<N, D> struct ExpmPadeHelper<N, D>
where where
N: RealField, N: RealField,
D: DimName + DimMin<D>, D: DimMin<D>,
DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>>, DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>>,
{ {
use_exact_norm: bool, use_exact_norm: bool,
@ -40,13 +41,14 @@ where
impl<N, D> ExpmPadeHelper<N, D> impl<N, D> ExpmPadeHelper<N, D>
where where
N: RealField, N: RealField,
D: DimName + DimMin<D>, D: DimMin<D>,
DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>>, DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>>,
{ {
fn new(a: MatrixN<N, D>, use_exact_norm: bool) -> Self { fn new(a: MatrixN<N, D>, use_exact_norm: bool) -> Self {
let (nrows, ncols) = a.data.shape();
ExpmPadeHelper { ExpmPadeHelper {
use_exact_norm, use_exact_norm,
ident: MatrixN::<N, D>::identity(), ident: MatrixN::<N, D>::identity_generic(nrows, ncols),
a, a,
a2: None, a2: None,
a4: None, a4: None,
@ -321,10 +323,11 @@ fn factorial(n: u128) -> u128 {
fn onenorm_matrix_power_nonm<N, D>(a: &MatrixN<N, D>, p: u64) -> N fn onenorm_matrix_power_nonm<N, D>(a: &MatrixN<N, D>, p: u64) -> N
where where
N: RealField, N: RealField,
D: DimName, D: Dim,
DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>, DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>,
{ {
let mut v = crate::VectorN::<N, D>::repeat(convert(1.0)); let nrows = a.data.shape().0;
let mut v = crate::VectorN::<N, D>::repeat_generic(nrows, U1, convert(1.0));
let m = a.transpose(); let m = a.transpose();
for _ in 0..p { for _ in 0..p {
@ -337,7 +340,7 @@ where
fn ell<N, D>(a: &MatrixN<N, D>, m: u64) -> u64 fn ell<N, D>(a: &MatrixN<N, D>, m: u64) -> u64
where where
N: RealField, N: RealField,
D: DimName, D: Dim,
DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>, DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>,
{ {
// 2m choose m = (2m)!/(m! * (2m-m)!) // 2m choose m = (2m)!/(m! * (2m-m)!)
@ -367,7 +370,7 @@ where
fn solve_p_q<N, D>(u: MatrixN<N, D>, v: MatrixN<N, D>) -> MatrixN<N, D> fn solve_p_q<N, D>(u: MatrixN<N, D>, v: MatrixN<N, D>) -> MatrixN<N, D>
where where
N: ComplexField, N: ComplexField,
D: DimMin<D, Output = D> + DimName, D: DimMin<D, Output = D>,
DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>>, DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>>,
{ {
let p = &u + &v; let p = &u + &v;
@ -379,7 +382,7 @@ where
fn one_norm<N, D>(m: &MatrixN<N, D>) -> N fn one_norm<N, D>(m: &MatrixN<N, D>) -> N
where where
N: RealField, N: RealField,
D: DimName, D: Dim,
DefaultAllocator: Allocator<N, D, D>, DefaultAllocator: Allocator<N, D, D>,
{ {
let mut col_sums = vec![N::zero(); m.ncols()]; let mut col_sums = vec![N::zero(); m.ncols()];
@ -396,7 +399,7 @@ where
impl<N: RealField, D> MatrixN<N, D> impl<N: RealField, D> MatrixN<N, D>
where where
D: DimMin<D, Output = D> + DimName, D: DimMin<D, Output = D>,
DefaultAllocator: DefaultAllocator:
Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>> + Allocator<N, D>, Allocator<N, D, D> + Allocator<(usize, usize), DimMinimum<D, D>> + Allocator<N, D>,
{ {

30
tests/linalg/exp.rs
View File

@ -2,7 +2,7 @@
mod tests { mod tests {
//https://github.com/scipy/scipy/blob/c1372d8aa90a73d8a52f135529293ff4edb98fc8/scipy/sparse/linalg/tests/test_matfuncs.py //https://github.com/scipy/scipy/blob/c1372d8aa90a73d8a52f135529293ff4edb98fc8/scipy/sparse/linalg/tests/test_matfuncs.py
#[test] #[test]
fn exp() { fn exp_static() {
use nalgebra::{Matrix1, Matrix2, Matrix3}; use nalgebra::{Matrix1, Matrix2, Matrix3};
{ {
@ -98,4 +98,32 @@ mod tests {
assert!(relative_eq!(f, m.exp(), epsilon = 1.0e-7)); assert!(relative_eq!(f, m.exp(), epsilon = 1.0e-7));
} }
} }
#[test]
fn exp_dynamic() {
use nalgebra::DMatrix;
let m = DMatrix::from_row_slice(3, 3, &[1.0, 3.0, 0.0, 0.0, 1.0, 5.0, 0.0, 0.0, 2.0]);
let e1 = 1.0_f64.exp();
let e2 = 2.0_f64.exp();
let f = DMatrix::from_row_slice(
3,
3,
&[
e1,
3.0 * e1,
15.0 * (e2 - 2.0 * e1),
0.0,
e1,
5.0 * (e2 - e1),
0.0,
0.0,
e2,
],
);
assert!(relative_eq!(f, m.exp(), epsilon = 1.0e-7));
}
} }