array_matrix/matrix/
eig.rs1use num_complex::{Complex};
2use num_traits::{Float};
3
4use crate::{MMul, Diag, QRHouseholder, Matrix};
5
6const ITERATIONS: usize = 1000;
7
8pub trait Eig: Matrix
9{
10 type Output;
11
12 fn eig(&self) -> Self::Output;
26}
27
28
29impl<F: Float, const L: usize, const H: usize> Eig for [[Complex<F>; L]; H]
30where
31 [[Complex<F>; L]; H]: QRHouseholder,
32 <[[Complex<F>; L]; H] as QRHouseholder>::OutputR: Diag
33 + MMul<<[[Complex<F>; L]; H] as QRHouseholder>::OutputQ, Output = [[Complex<F>; L]; H]>
34{
35 type Output = <<[[Complex<F>; L]; H] as QRHouseholder>::OutputR as Diag>::Output;
36 fn eig(&self) -> Self::Output
37 {
38 let mut a = self.clone();
39 for _ in 0..ITERATIONS
40 {
41 let (q, r) = a.qr_householder();
42 a = r.mul(q)
43 }
44 let r = a.qr_householder().1;
45 r.diag()
46 }
47}
48
49impl<const L: usize, const H: usize> Eig for [[f32; L]; H]
50where
51 Self: Matrix,
52 [[Complex<f32>; L]; H]: Eig
53{
54 type Output = <[[Complex<f32>; L]; H] as Eig>::Output;
55 fn eig(&self) -> Self::Output
56 {
57 self.map(|ar| ar.map(|arc| Complex::from(arc))).eig()
58 }
59}
60
61impl<const L: usize, const H: usize> Eig for [[f64; L]; H]
62where
63 Self: Matrix,
64 [[Complex<f64>; L]; H]: Eig
65{
66 type Output = <[[Complex<f64>; L]; H] as Eig>::Output;
67 fn eig(&self) -> Self::Output
68 {
69 self.map(|ar| ar.map(|arc| Complex::from(arc))).eig()
70 }
71}