Diagonalization of complex symmetric matrices: Generalized Householder reflections, iterative deflation and implicit shifts

Abstract

We describe a matrix diagonalization algorithm for complex symmetric (not Hermitian) matrices, A̲=A̲T, which is based on a two-step algorithm involving generalized Householder reflections based on the indefinite inner product 〈u̲,v̲〉∗=∑iuivi. This inner product is linear in both arguments and avoids complex conjugation. The complex symmetric input matrix is transformed to tridiagonal form using generalized Householder transformations (first step). An iterative, generalized QL decomposition of the tridiagonal matrix employing an implicit shift converges toward diagonal form (second step). The QL algorithm employs iterative deflation techniques when a machine-precision zero is encountered “prematurely” on the super-/sub-diagonal. The algorithm allows for a reliable and computationally efficient computation of resonance and antiresonance energies which emerge from complex-scaled Hamiltonians, and for the numerical determination of the real energy eigenvalues of pseudo-Hermitian and PT-symmetric Hamilton matrices. Numerical reference values are provided. Program summaryProgram Title: HTDQLSProgram Files doi:http://dx.doi.org/10.17632/x24wjxtrsg.1Licensing provisions: GPLv3Programming language: Fortran 90 using fixed form notationNature of problem: Calculating the eigenvalues and optionally the eigenvectors of complex symmetric (non-Hermitian), densely populated matrices.Solution method: The complex symmetric (not Hermitian) input matrix is diagonalized in two steps. First step: The matrix is tridiagonalized via a series of (n−2) generalized Householder reflections, where n is the rank of the input matrix. Second step: The tridiagonal matrix is diagonalized via a generalization of the “chasing the bulge” technique, which is an iterative process utilizing an implicitly shifted initial rotation followed by (n−2) Givens rotations. This technique is an implementation of QL factorization, and converges roughly as [(λi−σi)∕(λi+1−σi)]j where λi is the eigenvalue located in the (i,i) position of the final diagonal matrix and the eigenvalues are ordered (|λ1|<|λ2|<…<|λn|), and j is the iteration. The “educated guess” σi for the eigenvalue λi is obtained from the analytic determination of the eigenvalues of (k×k)-submatrices of A̲, in the vicinity of the ith element, where k=0,1,2,3 (here, k=0 means that the implicit shift vanishes, σi=0). The routine optionally calculates the rotation matrix Z̲, such that Z̲−1A̲Z̲=D̲ where A̲ is the input matrix and D̲ is the diagonal matrix containing the eigenvalues. The ith column of Z̲ then is the eigenvector of A̲ corresponding to the eigenvalue found at the element D̲(i,i), in the ith position on the diagonal of the matrix D̲.Unusual features:For simplicity, the “wrapper” program which contains an example application and the HTDQLS routine are distributed in the same file.