Function eigen_symmetric_tridiagonal 

pub fn eigen_symmetric_tridiagonal(
    d_in: Vec<f64>,
    e_in: Vec<f64>,
) -> (Vec<f64>, Vec<f64>)

Symmetric Tridiagonal Eigenvalue Solver

This module solves the eigenvalue problem $Ax = \lambda x$ for a real symmetric tridiagonal matrix $A$.

Mathematical Formulation

A symmetric tridiagonal matrix $A$ has the form:

$$ A = \begin{pmatrix} d_1 & e_1 & 0 & \cdots & 0 \ e_1 & d_2 & e_2 & \cdots & 0 \ 0 & e_2 & d_3 & \ddots & \vdots \ \vdots & \vdots & \ddots & \ddots & e_{n-1} \ 0 & 0 & \cdots & e_{n-1} & d_n \end{pmatrix} $$

where $d_i$ are the diagonal elements and $e_i$ are the off-diagonal elements.

Algorithm: Implicit QL/QR

The solver uses the Implicit QL (or QR) algorithm, which is the standard efficient method for this problem. The algorithm works by iteratively applying orthogonal similarity transformations to diagonalize the matrix.

1. Shift: At each step, a shift $\sigma$ is chosen (typically the Wilkinson shift) to accelerate convergence.
2. Transformation: An orthogonal matrix $Q$ is constructed such that $A’ = Q^T A Q$ is still tridiagonal but with smaller off-diagonal elements.
3. Convergence: The process repeats until the off-diagonal elements are negligible, at which point the diagonal elements are the eigenvalues.

The eigenvectors are accumulated by multiplying the transformation matrices: $Z = Q_1 Q_2 \cdots Q_k$.

Implementation Details

This implementation wraps the LAPACK routine DSTEV.

• Input: Diagonal vector $d$ (length $N$) and off-diagonal vector $e$ (length $N-1$).
• Output: Eigenvalues in ascending order and the corresponding eigenvectors.
• Complexity: $O(N^2)$ for eigenvalues and eigenvectors.

The eigenvectors are returned in a row-major format Vec<Vec<f64>>, where z[i][j] corresponds to the $i$-th component of the $j$-th eigenvector (or row $i$, column $j$ of the eigenvector matrix $Z$ where $A = Z \Lambda Z^T$).

Arguments

• d_in - The diagonal elements of the matrix.
• e_in - The off-diagonal elements of the matrix. e[i] is the element at (i, i+1).

Returns

A tuple containing:

• A vector of eigenvalues.
• A flattened vector of eigenvectors (column-major). The element at row i and column j is at index i + j * n.

Examples

use solver::linalg::eigen::eigen_symmetric_tridiagonal;

// Matrix:
// [ 2 -1 ]
// [-1  2 ]
// Eigenvalues should be 1 and 3.

let d = vec![2.0, 2.0];
let e = vec![-1.0]; // Off-diagonal

let (evals, evecs) = eigen_symmetric_tridiagonal(d, e);

assert!((evals[0] - 1.0).abs() < 1e-6);
assert!((evals[1] - 3.0).abs() < 1e-6);