|
NumericAL
A header-only C++20 lexer built on REAL
|
The eigenvalues and eigenvectors of a general square matrix. More...
#include <algorithm>#include <array>#include <cmath>#include <complex>#include <cstddef>#include <exception>#include <limits>#include <type_traits>#include <utility>#include <vector>#include "decomposition.hpp"Go to the source code of this file.
Classes | |
| struct | numerical::detail::balancing |
How a matrix was balanced: swaps, in the order they were applied, then the scaling of the active block low … high − 1 by scale. More... | |
| struct | numerical::detail::real_operand |
| A real vector that rows are summed against by dot_product. More... | |
| struct | numerical::detail::complex_operand |
A complex vector that rows are summed against, kept as (re, −im) and (im, re) pairs: a complex row, read as its parts in order ([complex.numbers]), against each is the real or the imaginary part of the sum, one dot_product each. More... | |
| struct | numerical::detail::hessenberg_panel< T > |
One panel of hessenberg_blocked, Q_panel = I − V T Vᴴ. More... | |
| struct | numerical::detail::standard_block |
| The rotation standardize_block found, and the eigenvalues of the block it standardised. More... | |
| struct | numerical::detail::schur_target< T > |
Where a QR sweep applies its reflections: the rows and columns it updates, and Z. More... | |
| struct | numerical::detail::one_based< T > |
A row-major matrix read through indices from 1, (1, 1) its first entry. More... | |
| struct | numerical::detail::plane_rotation< T > |
A plane rotation [c s; −s̄ c] and what it leaves of its first argument. More... | |
| struct | numerical::detail::deflation_count |
| What early_deflation found: the eigenvalues left as shifts, and those deflated. More... | |
| struct | numerical::detail::scaled_quotient< S > |
A quotient of divided and the factor its dividend was scaled by: a type of its own, where a std::pair of them made GCC on arm64 print, at every caller, a note on its ABI change of GCC 10.1. More... | |
| struct | numerical::detail::real_components |
| A real eigenvector during its back-substitution, the sums of rows against it by dot_product. More... | |
| struct | numerical::detail::split_components |
| A real matrix's complex eigenvector during its back-substitution, its real and imaginary parts apart: a row of the real Schur form against it is two dot_product. More... | |
| struct | numerical::detail::interleaved_components |
A complex matrix's eigenvector during its back-substitution, kept as (re, −im) and (im, re) pairs: a complex row, read as its parts in order ([complex.numbers]), against each is the real or the imaginary part of their product, one dot_product each. More... | |
| struct | numerical::detail::real_eigenvectors |
A real matrix's eigenvectors before their back-transformation: Xᵀ's real and imaginary parts. More... | |
| class | numerical::eigen_decomposition< T > |
The eigenvalues, and optionally the right eigenvectors, of a general square matrix, as numpy.linalg.eig. More... | |
Namespaces | |
| namespace | numerical |
| The NumericAL public API. | |
| namespace | numerical::detail |
Typedefs | |
| using | numerical::detail::lapack_index = std::ptrdiff_t |
| A signed index from 1, as LAPACK's: the multishift sweep starts its bulges in phantom columns above the matrix, before its first row. | |
Functions | |
| double | numerical::detail::modulus_1 (double x) |
|x| of a real x: what modulus_1 of a complex one stands for. | |
| double | numerical::detail::modulus_1 (const std::complex< double > &x) |
|re x| + |im x|, the cheap modulus LAPACK's complex routines judge with. | |
| template<typename A , typename B > | |
| auto | numerical::detail::times (const A &a, const B &b) |
a · b. | |
| std::complex< double > | numerical::detail::times (const std::complex< double > &a, const std::complex< double > &b) |
a · b of two complex numbers by the schoolbook formula, without C Annex G's recovery of infinities, which a large translation unit leaves out of line (a call per product: complex eig of order 64 took 2.7 times as long in the Python binding as alone). The entries here are finite. | |
| double | numerical::detail::largest_part (double x) |
|x| of a real x. | |
| double | numerical::detail::largest_part (const std::complex< double > &x) |
max(|re x|, |im x|): within √2 of |x|, without its square root. | |
| template<typename T > | |
| void | numerical::detail::exchange (buffer< T > &a, std::size_t n, std::size_t i, std::size_t j) |
Exchanges rows i and j, then columns i and j, of the row-major n × n a: a similarity. | |
| template<typename T > | |
| balancing | numerical::detail::balance (buffer< T > &a, std::size_t n) |
Balances the row-major n × n a in place, as LAPACK's ‘dgebal('B’)`: rows with no off-diagonal entry in the active block go to its bottom and columns with none to its top, each isolating an eigenvalue on the diagonal; then each row and column of the remaining block is scaled by a power of two until their 2-norms are within a factor of 2, which leaves every entry exact. | |
| template<typename T > | |
| buffer< T > | numerical::detail::hessenberg_unblocked (buffer< T > &a, std::size_t n, std::size_t low, std::size_t high, bool form) |
Reduces the active block low … high − 1 of the balanced row-major n × n a to upper Hessenberg form by Householder reflections one column at a time, as LAPACK's dgehd2: A = Q H Qᴴ, the entries below the subdiagonal set to zero. | |
| bool | numerical::detail::already_reduced (double, double rest) |
Whether (α, x), of tail norm rest, is already (β, 0) with β real: for a real one, its tail zero. | |
| bool | numerical::detail::already_reduced (const std::complex< double > &alpha, double rest) |
The same for a complex one: its tail zero and α real. | |
| template<typename T > | |
| T | numerical::detail::reflector_tau (T &alpha, T *x, std::size_t count) |
The Householder reflection H = I − τ v vᴴ, v = (1, v₁…), with Hᴴ (α, x) = (β, 0) and β real, as LAPACK's zlarfg (dlarfg for a real T): alpha becomes β and the count entries at x become v₁…. | |
| real_operand | numerical::detail::operand_of (const double *v, std::size_t stride, std::size_t count) |
The count entries stride apart at v as a real_operand. | |
| complex_operand | numerical::detail::operand_of (const std::complex< double > *v, std::size_t stride, std::size_t count) |
| The same as a complex_operand. | |
| template<typename T > | |
| buffer< T > | numerical::detail::hessenberg_blocked (buffer< T > &a, std::size_t n, std::size_t low, std::size_t high, bool form) |
hessenberg_unblocked by panels of hessenberg_block columns, as LAPACK's dgehrd and dlahr2: a panel's reflectors are formed one column at a time, each column first updated by those before it through Y = A V T; the rest of the matrix then takes the panel from both sides through the matrix product kernel, A ← (I − V Tᴴ Vᴴ)(A − Y Vᴴ). Its arguments and answer are hessenberg_unblocked's. | |
| template<typename T > | |
| buffer< T > | numerical::detail::hessenberg (buffer< T > &a, std::size_t n, std::size_t low, std::size_t high, bool form) |
| hessenberg_blocked from an active order of hessenberg_crossover, else hessenberg_unblocked. | |
| double | numerical::detail::small_reflector (double(&v)[3], std::size_t count) |
The Householder reflection I − τ v vᵀ, v = (1, x…), that takes (α, x…) to (β, 0…), as LAPACK's dlarfg for two or three entries: on return v[0] = β and v[1…] holds x's part of v. | |
| std::complex< double > | numerical::detail::complex_reflector (std::complex< double >(&v)[2]) |
The Householder reflection I − τ v vᴴ, v = (1, x), whose adjoint takes (α, x) to (β, 0), β real, as LAPACK's zlarfg for two entries: on return v[0] = β and v[1] holds x's part of v. | |
| standard_block | numerical::detail::standardize_block (double &a, double &b, double &c, double &d) |
The Schur factorisation of the real 2 × 2 block [a b; c d], as LAPACK's dlanv2: [a b; c d] ← [cs −sn; sn cs]ᵀ [a b; c d] [cs −sn; sn cs], upper triangular where its eigenvalues are real, else with equal diagonal entries and b c < 0. | |
| template<typename T > | |
| void | numerical::detail::rotate (T *x, T *y, std::size_t count, std::size_t stride, double cs, double sn) |
(x, y) ← (cs x + sn y, cs y − sn x) over count pairs stride apart: a plane rotation. | |
| template<typename T > | |
| void | numerical::detail::reflect_three_rows (T *r1, T *r2, T *r3, std::size_t from, std::size_t to, const T &tau, const T &v2, const T &v3, bool adjoint) |
A reflection I − τ v vᴴ, v = (1, v₂, v₃), on the entries from … to − 1 of three rows: its adjoint from the left with adjoint (rows of H), else from the right (rows of Zᵀ, which are columns of Z). A reflection of two passes r3 = r2 and v3 = 0. One loop over the entries, which vectorises. | |
| template<typename T > | |
| void | numerical::detail::reflect_two_rows (T *r1, T *r2, std::size_t from, std::size_t to, const T &tau, const T &v2, bool adjoint) |
reflect_three_rows for a reflection of two entries, v = (1, v₂), on two rows: its own loop, which vectorises where the third row aliasing the second would not. | |
| template<typename T > | |
| void | numerical::detail::reflect_three_columns (T *first, std::size_t rows, std::size_t ld, const T &tau, const T &v2, const T &v3) |
A reflection I − τ v vᴴ, v = (1, v₂, v₃), from the right on rows rows (ld apart) of the three columns at first: x ← x − τ (x₁ + v₂ x₂ + v₃ x₃) (1, v̄₂, v̄₃) on each. The rows are short and apart: one call for them all. | |
| template<typename T > | |
| void | numerical::detail::reflect_two_columns (T *first, std::size_t rows, std::size_t ld, const T &tau, const T &v2) |
The same for a reflection of two entries, v = (1, v₂), on two columns. | |
| void | numerical::detail::real_schur (const schur_target< double > &t, buffer< double > &re, buffer< double > &im, std::size_t steps) |
The real Schur form of the upper Hessenberg t's active block, as LAPACK's dlahqr: Francis's implicit double-shift QR, with the deflation criterion of Ahues and Tisseur and exceptional shifts after 10 and 20 iterations without deflation. The eigenvalues are left in re and im, a complex pair's in two consecutive places, the positive imaginary part first. | |
| void | numerical::detail::complex_schur (const schur_target< std::complex< double > > &t, buffer< std::complex< double > > &values, std::size_t steps) |
The complex Schur form of the upper Hessenberg t's active block, as LAPACK's zlahqr: the shifted QR algorithm with Wilkinson's shift, the deflation criterion of Ahues and Tisseur, and exceptional shifts after 10 and 20 iterations without deflation. The eigenvalues are left in values. | |
| std::size_t | numerical::detail::count_of (lapack_index count) |
| A count as an unsigned size (never negative where it is taken). | |
| constexpr lapack_index | numerical::detail::multishift_from (double) |
| The order from which a real Schur form is taken by multishift_schur, below which real_schur is faster (eight arm64 cores: order 128 1.85 against 2.60 ms; 224, 8.9 against 8.6; with the vectors, 192, 8.1 against 7.1). | |
| constexpr lapack_index | numerical::detail::multishift_from (const std::complex< double > &) |
| The same for a complex one (order 76 1.8 against 2.1 ms, 128 7.5 against 6.5; with the vectors, 100, 6.4 against 5.6). | |
| plane_rotation< double > | numerical::detail::givens (double f, double g) |
The rotation that takes (f, g) to (r, 0), as LAPACK's dlartg: r of f's sign. | |
| plane_rotation< std::complex< double > > | numerical::detail::givens (const std::complex< double > &f, const std::complex< double > &g) |
The rotation that takes (f, g) to (r, 0), as LAPACK's zlartg: c real, r of f's phase. | |
| template<typename T > | |
| void | numerical::detail::rotate_pair (T *x, T *y, std::size_t count, std::size_t stride, double c, const T &s) |
(x, y) ← (c x + s y, c y − s̄ x) over count pairs stride apart. | |
| template<typename T > | |
| void | numerical::detail::exchange_diagonal (one_based< T > t, one_based< T > qt, lapack_index n, lapack_index k) |
Exchanges the diagonal entries k and k + 1 of the upper triangular (or quasi-triangular, both 1 × 1 there) t by a plane rotation, as LAPACK's ztrexc step and dlaexc for two 1 × 1 blocks; the rows of qt (Qᵀ) turn with it. | |
| void | numerical::detail::reflect_rows_3 (const std::array< double, 3 > &u, double tau, double *t, std::size_t ld, std::size_t count) |
H = I − τ u uᵀ of three entries applied from the left to the 3 × count rows at t (row stride ld). | |
| void | numerical::detail::reflect_columns_3 (const std::array< double, 3 > &u, double tau, double *t, std::size_t ld, std::size_t count) |
The same from the right, to the count × 3 columns at t. | |
| void | numerical::detail::small_sylvester (const one_based< double > &tl, const one_based< double > &tr, const one_based< double > &b, lapack_index n1, lapack_index n2, double(&x)[4]) |
The solution X (n1 × n2, row-major in x) of TL X − X TR = B, n1 and n2 at most 2, by Gaussian elimination with complete pivoting on its Kronecker form, a pivot below ε times the largest coefficient taken as that, as LAPACK's dlasy2 (the blocks' entries near 1: no scaling against overflow). | |
| void | numerical::detail::standardize_at (one_based< double > t, one_based< double > qt, lapack_index n, lapack_index j) |
The standard form of the 2 × 2 block at j of the quasi-triangular t, its rotation carried to the rest of T and to the rows of qt, as dlaexc ends. | |
| bool | numerical::detail::exchange_blocks (one_based< double > t, one_based< double > qt, lapack_index n, lapack_index j1, lapack_index n1, lapack_index n2) |
Exchanges the adjacent diagonal blocks of the real quasi-triangular t at j1, of orders n1 then n2, by an orthogonal similarity carried to the rows of qt, as LAPACK's dlaexc. | |
| lapack_index | numerical::detail::block_order (one_based< double > t, lapack_index n, lapack_index j) |
The order of the diagonal block of the quasi-triangular t that starts at j: 2 where a subdiagonal entry below it is nonzero. | |
| bool | numerical::detail::reorder (one_based< double > t, one_based< double > qt, lapack_index n, lapack_index ifst, lapack_index &ilst) |
Moves the diagonal block of the real quasi-triangular t at ifst to ilst by exchanges of adjacent blocks, the rows of qt (Qᵀ) turning with it, as LAPACK's dtrexc. ilst becomes where the block ended: short of its goal where an exchange was refused (then. | |
| bool | numerical::detail::reorder (one_based< std::complex< double > > t, one_based< std::complex< double > > qt, lapack_index n, lapack_index ifst, lapack_index &ilst) |
Moves the diagonal entry of the complex upper triangular t at ifst to ilst, as LAPACK's ztrexc. | |
| void | numerical::detail::bulge_column (one_based< double > h, lapack_index count, const std::complex< double > &s1, const std::complex< double > &s2, double *v) |
A multiple of the first column of (H − s₁ I)(H − s₂ I) for the leading 2 × 2 or 3 × 3 (count) of the real Hessenberg h, the shifts a conjugate pair or two reals, as LAPACK's dlaqr1. | |
| void | numerical::detail::bulge_column (one_based< std::complex< double > > h, lapack_index count, const std::complex< double > &s1, const std::complex< double > &s2, std::complex< double > *v) |
The same for a complex h and any two shifts, as LAPACK's zlaqr1. | |
| void | numerical::detail::pair_shifts (std::complex< double > *s, lapack_index count, double) |
Arranges a real matrix's shifts so that each pair is a conjugate pair or two reals, as dlaqr5 begins (conjugate shifts already adjacent). | |
| void | numerical::detail::pair_shifts (std::complex< double > *, lapack_index, const std::complex< double > &) |
| A complex matrix's shifts pair as they come. | |
| double | numerical::detail::as_scalar (const std::complex< double > &s, double) |
The shift s as the matrix's scalar: a real matrix's are real where they are used one by one. | |
| std::complex< double > | numerical::detail::as_scalar (const std::complex< double > &s, const std::complex< double > &) |
| A complex matrix's shift as it is. | |
| template<typename T > | |
| void | numerical::detail::vigilant_deflation (one_based< T > h, lapack_index k, lapack_index ktop, lapack_index kbot, double smallest) |
H(k + 1, k) set to zero where the deflation criterion of Ahues and Tisseur holds, as the sweep of dlaqr5 tests after each step of a bulge (vigilant deflation). | |
| template<typename T > | |
| void | numerical::detail::reflect_from_left (T *x, std::size_t stride, const T &tau, const T &v2, const T &v3) |
The adjoint I − τ̄ v vᴴ of a reflection of three entries from the left on one column (x three entries stride apart). | |
| template<typename T > | |
| void | numerical::detail::chase_bulges (const schur_target< T > &target, lapack_index ktop, lapack_index kbot, std::complex< double > *s, lapack_index count, bool accumulate) |
One small-bulge multishift QR sweep over the active block ktop … kbot of target's Hessenberg matrix with the count shifts at s, as LAPACK's dlaqr5 and zlaqr5: chains of 3 × 3 bulges, two shifts each, chased down the diagonal by slabs of 4 · count / 2 columns; with accumulate, each slab's reflections are accumulated in U and applied to the rest of the matrix and to Z through the matrix product kernel, else applied to them directly (LAPACK's kacc22, which accumulates from 14 shifts). A subdiagonal entry that becomes negligible on the way is set to zero (vigilant deflation). | |
| template<typename T > | |
| void | numerical::detail::multishift_schur (const schur_target< T > &target, std::complex< double > *w, bool recursive, lapack_index steps) |
The Schur form of target's active block by the multishift QR algorithm with aggressive early deflation, as LAPACK's dlaqr0 and zlaqr0 (dlaqr4 without recursive): each iteration deflates what it can from a window at the bottom of the active block (early_deflation) and, unless that deflated enough, sweeps the block with the window's other eigenvalues as shifts (chase_bulges). The eigenvalues are left in w (one per row; a real matrix's complex pairs adjacent, the positive imaginary part first). | |
| void | numerical::detail::plain_schur (const schur_target< double > &t, std::complex< double > *w) |
The Schur form of t's active block by real_schur, its eigenvalues to w (one per row of t). | |
| void | numerical::detail::plain_schur (const schur_target< std::complex< double > > &t, std::complex< double > *w) |
| The same by complex_schur. | |
| lapack_index | numerical::detail::undeflated (one_based< double > t, one_based< double > vt, lapack_index jw, double s, double smallest, std::complex< double > *w) |
The eigenvalues of the real Schur form t (jw × jw) of a deflation window that do not deflate, moved to its top, as dlaqr3: from the bottom, a block whose spike entries s v₁ⱼ are negligible deflates, any other is moved up; the undeflated blocks are then sorted by decreasing modulus (which helps graded matrices), and every eigenvalue of the window written to w. | |
| lapack_index | numerical::detail::undeflated (one_based< std::complex< double > > t, one_based< std::complex< double > > vt, lapack_index jw, const std::complex< double > &s, double smallest, std::complex< double > *w) |
The same for a complex Schur form, as zlaqr3: 1 × 1 blocks only, sorted by selection. | |
| template<typename T > | |
| deflation_count | numerical::detail::early_deflation (const schur_target< T > &target, lapack_index ktop, lapack_index kbot, lapack_index nw, std::complex< double > *w, bool recursive) |
Aggressive early deflation, as LAPACK's dlaqr3 and zlaqr3: the trailing nw × nw window of the active block ktop … kbot is brought to Schur form; its eigenvalues whose coupling to the rest (the spike, s times the first row of the window's Schur vectors) is negligible deflate; the window is reduced back to Hessenberg form above them, and its transformation carried to the rest of H and to Z through the matrix product kernel. Every eigenvalue of the window is written to w: the deflated ones final, the others the next shifts. recursive: the window itself by multishift_schur past its crossover. | |
| lapack_index | numerical::detail::recommended_shifts (lapack_index nh) |
The shifts LAPACK's iparmq recommends for an active block of order nh. | |
| lapack_index | numerical::detail::recommended_window (lapack_index nh, lapack_index ns) |
The deflation window LAPACK's iparmq recommends for an active block of order nh and ns shifts. | |
| void | numerical::detail::exceptional_shifts (one_based< double > h, std::complex< double > *w, lapack_index ks, lapack_index kbot, lapack_index ktop) |
Exceptional shifts for a real matrix, as dlaqr0's after kexsh iterations without deflation: from the subdiagonal entries, as conjugate pairs of the 2 × 2 blocks dlahqr's exceptional shifts build. | |
| void | numerical::detail::exceptional_shifts (one_based< std::complex< double > > h, std::complex< double > *w, lapack_index ks, lapack_index kbot, lapack_index) |
Exceptional shifts for a complex matrix, as zlaqr0's. | |
| void | numerical::detail::pair_for_sweep (std::complex< double > *w, lapack_index ks, lapack_index kbot, double) |
A real matrix's shifts brought into pairs, each a conjugate pair or two reals, as dlaqr0 does before its sweep (conjugate shifts already adjacent). | |
| void | numerical::detail::pair_for_sweep (std::complex< double > *, lapack_index, lapack_index, const std::complex< double > &) |
| A complex matrix's shifts need no pairing. | |
| bool | numerical::detail::two_real_shifts (const std::complex< double > *w, lapack_index ks, lapack_index kbot, double) |
| Whether a real matrix's sweep is down to two shifts, both real: then the one nearer the corner serves twice. | |
| bool | numerical::detail::two_real_shifts (const std::complex< double > *, lapack_index ks, lapack_index kbot, const std::complex< double > &) |
| The same for a complex matrix: two shifts of any kind. | |
| void | numerical::detail::schur_form (const schur_target< double > &t, buffer< double > &re, buffer< double > &im, std::size_t steps) |
The real Schur form of t's active block: real_schur below multishift_from, else multishift_schur, as LAPACK's dhseqr chooses. steps bounds real_schur's iterations. | |
| void | numerical::detail::schur_form (const schur_target< std::complex< double > > &t, buffer< std::complex< double > > &values, std::size_t steps) |
| The complex Schur form likewise: complex_schur or multishift_schur. | |
| template<typename S > | |
| double | numerical::detail::solve_block (const S(&m)[4], S(&r)[2], double floor) |
Solves [a b; c d] x = r in place of r by Gaussian elimination with complete pivoting, a pivot below floor taken as floor, as LAPACK's dlaln2; r is scaled first, and. | |
| template<typename S , typename P > | |
| scaled_quotient< S > | numerical::detail::divided (S r, P pivot, double floor) |
r / pivot, a pivot below floor taken as floor, r scaled first where the quotient would pass back_substitution_limit. | |
| void | numerical::detail::normalize_vector (std::complex< double > *v, std::size_t n, bool signed_real) |
Scales v to unit Euclidean norm and turns it so that its first component of largest modulus is real: positive for a complex matrix (zgeev), of its real part's sign for a real one's pair (dgeev). | |
| template<typename T , typename X > | |
| void | numerical::detail::solve_upward (const buffer< T > &t, std::size_t n, std::size_t k, std::size_t top, typename X::scalar lambda, double floor, X &x) |
Solves (T − λ I) x = 0 upward from the component k, rows k − 1 … 0 of the row-major quasi-triangular t (its 2 × 2 blocks where a subdiagonal entry is nonzero), x [k … top] given: as LAPACK's dtrevc and ztrevc, a pivot below floor taken as floor, the components scaled down where they would grow past back_substitution_limit. The sums of rows against x are X 's, vectorised. | |
| buffer< std::complex< double > > | numerical::detail::complex_vectors (const buffer< std::complex< double > > &t, std::size_t n) |
The eigenvectors of the complex upper triangular t (row-major n × n) as the rows of Xᵀ: row k is x with (T − t_kk I) x = 0, x_k = 1, zero past k, as LAPACK's ztrevc. | |
| real_eigenvectors | numerical::detail::real_vectors (const buffer< double > &t, const buffer< double > &re, const buffer< double > &im, std::size_t n) |
The eigenvectors of the real quasi-triangular t (row-major n × n, its 2 × 2 blocks standard) as the rows of Xᵀ, as LAPACK's dtrevc: real for a real eigenvalue, solved in real arithmetic; for a pair re ± i im (im [k] > 0 the first) the first's vector, the second's its conjugate. | |
| template<typename T > | |
| buffer< T > | numerical::detail::back_transformed (const buffer< T > &xt, const buffer< T > &zt, std::size_t n) |
Vᵀ = Xᵀ Zᵀ, row-major n × n, through the matrix product kernel. | |
| template<typename T > | |
| vector< std::complex< double > > | numerical::eigenvalues (const matrix< T > &a) |
The eigenvalues alone of the square a, as numpy.linalg.eigvals. | |
Variables | |
| constexpr double | numerical::detail::schur_ulp {std::numeric_limits<double>::epsilon()} |
| The relative spacing LAPACK's ‘dlamch('P’)` names: the roundoff a QR sweep judges a subdiagonal by. | |
| constexpr std::size_t | numerical::detail::schur_steps_per_order {30} |
The iterations one eigenvalue (or pair) may take per unit of the active order, as dlahqr's 30. | |
| constexpr double | numerical::detail::balance_limit {0x1p256} |
| The largest power of two a balancing scales a row or a column by, either way: the matrix enters at a largest magnitude near 1, so no entry overflows, and the scaling stops where LAPACK's would be one factor among dozens. | |
| constexpr double | numerical::detail::back_substitution_limit {0x1p400} |
| The modulus an eigenvector's components are kept under during the back-substitution: the sums stay finite for any order a matrix can have, and a component that would grow past it scales the vector down. | |
| constexpr std::size_t | numerical::detail::hessenberg_block {32} |
| The columns of a panel of hessenberg_blocked. | |
| constexpr std::size_t | numerical::detail::hessenberg_crossover {128} |
| The active order from which hessenberg reduces by panels. | |
| constexpr lapack_index | numerical::detail::multishift_crossover {75} |
LAPACK's nmin: the active order below which dlaqr0's iterations skip a sweep after enough deflation. | |
| constexpr lapack_index | numerical::detail::window_growth {5} |
The iterations without deflation after which the window grows (kexnw). | |
| constexpr lapack_index | numerical::detail::exceptional_after {6} |
The iterations without deflation after which the shifts are exceptional (kexsh). | |
| constexpr lapack_index | numerical::detail::nibble {14} |
The share of the window, in percent, whose deflation makes the sweep unnecessary (nibble). | |
The eigenvalues and eigenvectors of a general square matrix.
numerical::eigen_decomposition follows LAPACK's dgeev and zgeev: the matrix is balanced (permuted to isolate the eigenvalues that need no iteration, then scaled by powers of two so that its rows and columns have comparable norms), reduced to upper Hessenberg form by Householder reflections, and brought to Schur form by the shifted QR algorithm: Francis's implicit double shift for a real matrix (dlahqr), whose Schur form is quasi-triangular with a 2 × 2 block for each pair of complex conjugate eigenvalues, and a single Wilkinson shift for a complex one (zlahqr). The eigenvectors are those of the Schur form, by back-substitution (dtrevc, ztrevc), carried back through the reductions; each has unit Euclidean norm and its component of largest modulus real, as LAPACK's.
The eigenvalues come in the order the Schur form holds them, which is LAPACK's: no order is promised. A real matrix's complex eigenvalues come in conjugate pairs, the one of positive imaginary part first, and their vectors are conjugate to the bit.
Definition in file eigen.hpp.