Jacobi eigenvalue decomposition for small symmetric matrices. More...
#include "config.h"#include "utilities/expected.h"#include "matrix_traits.h"#include "fused/fused_matrix.h"#include "fused/fused_vector.h"#include "math/math_utils.h"
Include dependency graph for eigen.h:
Go to the source code of this file.
Classes | |
| struct | matrix_algorithms::EigenResult< T, N > |
| Result of eigenvalue decomposition. More... | |
Namespaces | |
| namespace | matrix_algorithms |
Functions | |
| template<typename T , my_size_t N> | |
| Expected< EigenResult< T, N >, MatrixStatus > | matrix_algorithms::eigen_jacobi (const FusedMatrix< T, N, N > &A, my_size_t max_iters=100, T tol=T(PRECISION_TOLERANCE)) |
| Compute eigenvalues and eigenvectors of a symmetric matrix via Jacobi. | |
Detailed Description
Jacobi eigenvalue decomposition for small symmetric matrices.
Computes all eigenvalues and eigenvectors of a symmetric matrix A by iteratively applying Givens rotations to zero the largest off-diagonal element until convergence. The result is A = V·D·Vᵀ where:
- D is diagonal (eigenvalues)
- V is orthogonal (eigenvectors as columns)
ALGORITHM
Repeat until off-diagonal norm < tolerance:
- Find the largest off-diagonal element |A(p,q)|
- Compute Givens rotation angle θ to zero A(p,q): If A(p,p) == A(q,q): θ = π/4 Else: τ = (A(p,p) − A(q,q)) / (2·A(p,q)) Compute t = sign(τ) / (|τ| + √(1+τ²)) (smaller root for stability) c = 1/√(1+t²), s = t·c
- Apply rotation: A' = Jᵀ·A·J (only affects rows/cols p and q)
- Accumulate: V = V·J
Complexity: O(N²) per sweep, typically 5–10 sweeps for convergence. Total: O(5N²) to O(10N²) for small matrices.
NOTES
- Only for symmetric matrices. Returns NotSymmetric otherwise.
- Convergence is guaranteed for symmetric matrices (classical Jacobi).
- Best suited for small matrices (N ≤ ~20). For larger matrices, tridiagonalization + QR iteration (not implemented) is preferred.
- Eigenvalues are returned unsorted. Caller can sort if needed.
FAILURE MODES
- MatrixStatus::NotSymmetric — input fails isSymmetric() check
- MatrixStatus::NotConverged — max iterations exceeded
