template<typename RealScalar_>
class Eigen::DPR1EigenSolver< RealScalar_ >
Direct O(n^2) eigensolver for real symmetric diagonal-plus-rank-one matrices \( A = D + \rho\, z z^T \), via the secular equation.
This is the standalone version of the kernel at the heart of the divide-and-conquer symmetric eigensolvers (LAPACK's xLAED2/3/4): after deflation – entries with negligible \( |z_i| \) are eigenpairs of the diagonal already, and (nearly) equal diagonal entries are combined by Givens rotations whose dropped coupling is below a backward-stability threshold – the surviving eigenvalues are the roots of the secular equation
\[ f(\lambda) = 1 + \rho \sum_i \frac{z_i^2}{d_i - \lambda} = 0, \]
one in each interval between consecutive poles. Each root is bracketed and bisected in coordinates shifted to its nearest pole, so every distance \( \lambda - d_i \) is retained as an exact data difference plus a small offset instead of a cancellation-prone subtraction of close numbers. Eigenvectors are then built not from the original z but from the Gu-Eisenstat vector \( \hat z \) – the one for which the computed roots are exact secular eigenvalues – which is what makes the computed eigenvector matrix numerically orthogonal without any reorthogonalization.
The total cost is O(n^2): O(n log(1/eps)) per bisected root in the common case (the iteration cap is sized to the scalar's full exponent range, so even roots subnormally close to their pole resolve), O(n) per Gu-Eisenstat weight, and O(n) per eigenvector (the deflation rotations are replayed instead of accumulated into a dense matrix).
Both signs of \( \rho \) are supported (negative \( \rho \) is handled by negating the matrix), as are \( \rho = 0 \), zero z, repeated diagonal entries and any ordering of d. The problem is rescaled internally by the exact power of two that brings \( \max(\|D\|_\infty, |\rho| \|z\|^2) \) into [1/2, 1), which changes no bits of the result when \( \rho \|z\|^2 \) and the scaled poles are normal. InvalidInput is reported for non-finite input, for an update \( \rho \|z\|^2 \) that overflows, and for a computed eigenvalue that is not finite.
DPR1EigenSolver()=default
Matrix< double, Dynamic, Dynamic > MatrixXd
Matrix< double, Dynamic, 1 > VectorXd
- Template Parameters
-
| RealScalar_ | one of float, double, or long double. |
References:
- M. Gu and S. C. Eisenstat, "A stable and efficient algorithm for the
rank-one modification of the symmetric eigenproblem," SIAM J. Matrix Anal. Appl., 15(4):1266-1276, 1994.
- P. H. Sterbenz, "Floating-Point Computation", Prentice-Hall, 1974. Scaling by a power of two is exact, the property the problem scaling relies on.
- See also
- class SelfAdjointEigenSolver