An efficient solver for large structured eigenvalue problems in relativistic quantum chemistry

We report an efficient program for computing the eigenvalues and symmetry-adapted eigenvectors of very large quaternionic (or Hermitian skew-Hamiltonian) matrices, using which structure-preserving diagonalization of matrices of dimension N > 10000 is now routine on a single computer node. Such matrices appear frequently in relativistic quantum chemistry owing to the time-reversal symmetry. The implementation is based on a blocked version of the Paige-Van Loan algorithm [D. Kressner, BIT 43, 775 (2003)], which allows us to use the Level 3 BLAS subroutines for most of the computations. Taking advantage of the symmetry, the program is faster by up to a factor of two than state-of-the-art implementations of complex Hermitian diagonalization; diagonalizing a 12800 x 12800 matrix took 42.8 (9.5) and 85.6 (12.6) minutes with 1 CPU core (16 CPU cores) using our symmetry-adapted solver and Intel MKL's ZHEEV that is not structure-preserving, respectively. The source code is publicly available under the FreeBSD license.