working on issue 1106
This commit is contained in:
parent
3d52327f82
commit
d61f9e29ba
|
@ -7,9 +7,8 @@ use num_complex::Complex;
|
||||||
use simba::scalar::RealField;
|
use simba::scalar::RealField;
|
||||||
|
|
||||||
use crate::ComplexHelper;
|
use crate::ComplexHelper;
|
||||||
use na::allocator::Allocator;
|
use na::dimension::{Const, Dim, DimName};
|
||||||
use na::dimension::{Const, Dim};
|
use na::{DefaultAllocator, Matrix, OMatrix, OVector, Scalar, allocator::Allocator};
|
||||||
use na::{DefaultAllocator, Matrix, OMatrix, OVector, Scalar};
|
|
||||||
|
|
||||||
use lapack;
|
use lapack;
|
||||||
|
|
||||||
|
@ -148,7 +147,7 @@ where
|
||||||
eigenvalues_re: wr,
|
eigenvalues_re: wr,
|
||||||
eigenvalues_im: wi,
|
eigenvalues_im: wi,
|
||||||
left_eigenvectors: vl,
|
left_eigenvectors: vl,
|
||||||
eigenvectors: vr,
|
eigenvectors: vr
|
||||||
})
|
})
|
||||||
}
|
}
|
||||||
|
|
||||||
|
@ -168,6 +167,64 @@ where
|
||||||
|
|
||||||
det
|
det
|
||||||
}
|
}
|
||||||
|
|
||||||
|
/// 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() {
|
||||||
|
true => (None, None, None),
|
||||||
|
false => {
|
||||||
|
let number_of_elements = self.eigenvalues_re.nrows();
|
||||||
|
let number_of_complex_entries = self.eigenvalues_im.iter().fold(0, |acc, e| if !e.is_zero() {acc + 1} else {acc});
|
||||||
|
let mut eigenvalues = Vec::<Complex<T>>::with_capacity(2*number_of_complex_entries);
|
||||||
|
let mut eigenvectors = match self.eigenvectors.is_some() {
|
||||||
|
true => Some(Vec::<OVector<Complex<T>, D>>::with_capacity(2*number_of_complex_entries)),
|
||||||
|
false => None
|
||||||
|
};
|
||||||
|
let mut left_eigenvectors = match self.left_eigenvectors.is_some() {
|
||||||
|
true => Some(Vec::<OVector<Complex<T>, D>>::with_capacity(2*number_of_complex_entries)),
|
||||||
|
false => None
|
||||||
|
};
|
||||||
|
|
||||||
|
let eigenvectors_raw = self.eigenvectors;
|
||||||
|
let left_eigenvectors_raw = self.left_eigenvectors;
|
||||||
|
|
||||||
|
for mut i in 0..number_of_elements {
|
||||||
|
if self.eigenvalues_im[i] != T::zero() {
|
||||||
|
//Complex conjugate pairs of eigenvalues appear consecutively with the eigenvalue having the positive imaginary part first.
|
||||||
|
eigenvalues.push(Complex::<T>::new(self.eigenvalues_re[i].clone(), self.eigenvalues_im[i].clone()));
|
||||||
|
eigenvalues.push(Complex::<T>::new(self.eigenvalues_re[i].clone(), -self.eigenvalues_im[i].clone()));
|
||||||
|
|
||||||
|
if eigenvectors.is_some() {
|
||||||
|
let mut r1_vec = OVector::<Complex<T>, D>::zeros(number_of_elements);
|
||||||
|
let mut r1_vec_conj = OVector::<Complex<T>, D>::zeros(number_of_elements);
|
||||||
|
|
||||||
|
for j in 0..number_of_elements {
|
||||||
|
r1_vec[j] = Complex::<T>::new(self.eigenvectors.unwrap()[(i,j)].clone(),self.eigenvectors.unwrap()[(i,j+1)].clone());
|
||||||
|
r1_vec_conj[j] = Complex::<T>::new(self.eigenvectors.unwrap()[(i,j)].clone(),-self.eigenvectors.unwrap()[(i,j+1)].clone());
|
||||||
|
}
|
||||||
|
|
||||||
|
eigenvectors.unwrap().push(r1_vec);
|
||||||
|
eigenvectors.unwrap().push(r1_vec_conj);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
if left_eigenvectors.is_some() {
|
||||||
|
//TODO: Do the same for left
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
i += 1;
|
||||||
|
}
|
||||||
|
|
||||||
|
}
|
||||||
|
(Some(eigenvalues), left_eigenvectors, eigenvectors)
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
}
|
||||||
|
|
||||||
}
|
}
|
||||||
|
|
||||||
/*
|
/*
|
||||||
|
|
Loading…
Reference in New Issue