testing complex eigenvalues

This commit is contained in:
Marc Haubenstock 2022-10-23 20:36:06 +02:00
parent 4a1fd605e4
commit aa89cda0cd
5 changed files with 49 additions and 2 deletions

View File

@ -22,7 +22,7 @@ proptest-support = [ "nalgebra/proptest-support" ]
arbitrary = [ "nalgebra/arbitrary" ]
# For BLAS/LAPACK
default = ["netlib"]
default = ["openblas"]
openblas = ["lapack-src/openblas"]
netlib = ["lapack-src/netlib"]
accelerate = ["lapack-src/accelerate"]
@ -36,6 +36,7 @@ simba = "0.7"
serde = { version = "1.0", features = [ "derive" ], optional = true }
lapack = { version = "0.19", default-features = false }
lapack-src = { version = "0.8", default-features = false }
openblas-src = {version = "*", features = ["static"]}
# clippy = "*"
[dev-dependencies]
@ -44,3 +45,4 @@ proptest = { version = "1", default-features = false, features = ["std"] }
quickcheck = "1"
approx = "0.5"
rand = "0.8"

View File

@ -0,0 +1,24 @@
extern crate nalgebra as na;
extern crate nalgebra_lapack;
#[macro_use]
extern crate approx; // for assert_relative_eq
use na::Matrix3;
use nalgebra_lapack::Eigen;
use num_complex::Complex;
fn main() {
let m = Matrix3::<f64>::new(4.0/5.0, -3.0/5.0, 0.0, 3.0/5.0, 4.0/5.0, 0.0, 1.0, 2.0, 2.0);
let eigen = Eigen::new(m,true,true).expect("Eigen Creation Failed!");
let (some_eigenvalues, some_left_vec, some_right_vec) = eigen.get_complex_elements();
let eigenvalues = some_eigenvalues.expect("Eigenvalues Failed");
let left_eigenvectors = some_left_vec.expect("Left Eigenvectors Failed");
let eigenvectors = some_right_vec.expect("Right Eigenvectors Failed");
assert_relative_eq!(eigenvalues[0].re, Complex::<f64>::new(4.0/5.0,3.0/5.0).re);
assert_relative_eq!(eigenvalues[0].im, Complex::<f64>::new(4.0/5.0,3.0/5.0).im);
assert_relative_eq!(eigenvalues[1].re, Complex::<f64>::new(4.0/5.0,-3.0/5.0).re);
assert_relative_eq!(eigenvalues[1].im, Complex::<f64>::new(4.0/5.0,-3.0/5.0).im);
}

View File

@ -171,7 +171,8 @@ where
/// Returns a tuple of vectors. The elements of the tuple are the complex eigenvalues, complex left eigenvectors and complex right eigenvectors respectively.
/// The elements appear as conjugate pairs within each vector, with the positive of the pair always being first.
pub fn get_complex_elements(&self) -> (Option<Vec<Complex<T>>>, Option<Vec<OVector<Complex<T>, D>>>, Option<Vec<OVector<Complex<T>, D>>>) where DefaultAllocator: Allocator<Complex<T>, D> {
match !self.eigenvalues_are_real() {
match self.eigenvalues_are_real() {
true => (None, None, None),
false => {
let (number_of_elements, _) = self.eigenvalues_re.shape_generic();

View File

@ -0,0 +1,19 @@
use std::cmp;
use na::{Matrix3};
use nalgebra_lapack::Eigen;
use crate::proptest::*;
use proptest::{prop_assert, proptest};
proptest! {
//#[test]
// fn complex_eigen() {
// let n = cmp::max(1, cmp::min(n, 10));
// let m = DMatrix::<f64>::new_random(n, n);
// let eig = SymmetricEigen::new(m.clone());
// let recomp = eig.recompose();
// prop_assert!(relative_eq!(m.lower_triangle(), recomp.lower_triangle(), epsilon = 1.0e-5))
// }
}

View File

@ -7,3 +7,4 @@ mod real_eigensystem;
mod schur;
mod svd;
mod symmetric_eigen;
mod complex_eigen;