Skip to content

Avoid unnecessary calculation of eigenvectors in EigVals::eigvals and relax the DataMut bound #286

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 3 commits into from
Jun 12, 2021
Merged

Avoid unnecessary calculation of eigenvectors in EigVals::eigvals and relax the DataMut bound #286

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
91b2be2
Include EigVals::eigvals in tests
jturner314 May 11, 2021
6eb1d85
Avoid calculating eigenvectors unless requested
jturner314 May 11, 2021
484cabe
Relax bound on EigVals from DataMut to Data
jturner314 May 12, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions ndarray-linalg/src/eig.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -66,13 +66,13 @@ pub trait EigVals {
impl<A, S> EigVals for ArrayBase<S, Ix2>
where
A: Scalar + Lapack,
S: DataMut<Elem = A>,
S: Data<Elem = A>,
{
type EigVal = Array1<A::Complex>;

fn eigvals(&self) -> Result<Self::EigVal> {
let mut a = self.to_owned();
let (s, _) = A::eig(true, a.square_layout()?, a.as_allocated_mut()?)?;
let (s, _) = A::eig(false, a.square_layout()?, a.as_allocated_mut()?)?;
Ok(ArrayBase::from(s))
}
}
55 changes: 37 additions & 18 deletions ndarray-linalg/tests/eig.rs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,8 +2,11 @@ use ndarray::*;
use ndarray_linalg::*;

// Test Av_i = e_i v_i for i = 0..n
fn test_eig<T: Scalar>(a: Array2<T>, eigs: Array1<T::Complex>, vecs: Array2<T::Complex>)
where
fn test_eig<T: Scalar>(
a: ArrayView2<'_, T>,
eigs: ArrayView1<'_, T::Complex>,
vecs: ArrayView2<'_, T::Complex>,
) where
T::Complex: Lapack,
{
println!("a\n{:+.4}", &a);
Expand Down Expand Up @@ -205,29 +208,37 @@ macro_rules! impl_test_real {
#[test]
fn [<$real _eigvals >]() {
let a = test_matrix_real::<$real>();
let (e, _vecs) = a.eig().unwrap();
assert_close_l2!(&e, &answer_eig_real::<$real>(), 1.0e-3);
let (e1, _vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
assert_close_l2!(&e1, &answer_eig_real::<$real>(), 1.0e-3);
assert_close_l2!(&e2, &answer_eig_real::<$real>(), 1.0e-3);
}

#[test]
fn [<$real _eigvals_t>]() {
let a = test_matrix_real_t::<$real>();
let (e, _vecs) = a.eig().unwrap();
assert_close_l2!(&e, &answer_eig_real::<$real>(), 1.0e-3);
let (e1, _vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
assert_close_l2!(&e1, &answer_eig_real::<$real>(), 1.0e-3);
assert_close_l2!(&e2, &answer_eig_real::<$real>(), 1.0e-3);
}

#[test]
fn [<$real _eig>]() {
let a = test_matrix_real::<$real>();
let (e, vecs) = a.eig().unwrap();
test_eig(a, e, vecs);
let (e1, vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
test_eig(a.view(), e1.view(), vecs.view());
test_eig(a.view(), e2.view(), vecs.view());
}

#[test]
fn [<$real _eig_t>]() {
let a = test_matrix_real_t::<$real>();
let (e, vecs) = a.eig().unwrap();
test_eig(a, e, vecs);
let (e1, vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
test_eig(a.view(), e1.view(), vecs.view());
test_eig(a.view(), e2.view(), vecs.view());
}

} // paste::item!
Expand All @@ -243,15 +254,19 @@ macro_rules! impl_test_complex {
#[test]
fn [<$complex _eigvals >]() {
let a = test_matrix_complex::<$complex>();
let (e, _vecs) = a.eig().unwrap();
assert_close_l2!(&e, &answer_eig_complex::<$complex>(), 1.0e-3);
let (e1, _vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
assert_close_l2!(&e1, &answer_eig_complex::<$complex>(), 1.0e-3);
assert_close_l2!(&e2, &answer_eig_complex::<$complex>(), 1.0e-3);
}

#[test]
fn [<$complex _eigvals_t>]() {
let a = test_matrix_complex_t::<$complex>();
let (e, _vecs) = a.eig().unwrap();
assert_close_l2!(&e, &answer_eig_complex::<$complex>(), 1.0e-3);
let (e1, _vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
assert_close_l2!(&e1, &answer_eig_complex::<$complex>(), 1.0e-3);
assert_close_l2!(&e2, &answer_eig_complex::<$complex>(), 1.0e-3);
}

#[test]
Expand All @@ -271,15 +286,19 @@ macro_rules! impl_test_complex {
#[test]
fn [<$complex _eig>]() {
let a = test_matrix_complex::<$complex>();
let (e, vecs) = a.eig().unwrap();
test_eig(a, e, vecs);
let (e1, vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
test_eig(a.view(), e1.view(), vecs.view());
test_eig(a.view(), e2.view(), vecs.view());
}

#[test]
fn [<$complex _eig_t>]() {
let a = test_matrix_complex_t::<$complex>();
let (e, vecs) = a.eig().unwrap();
test_eig(a, e, vecs);
let (e1, vecs) = a.eig().unwrap();
let e2 = a.eigvals().unwrap();
test_eig(a.view(), e1.view(), vecs.view());
test_eig(a.view(), e2.view(), vecs.view());
}
} // paste::item!
};
Expand Down