Hermitian Eigenproblems Research Articles

Refined isogeometric analysis for generalized Hermitian eigenproblems

We use refined isogeometric analysis (rIGA) to solve generalized Hermitian eigenproblems (Ku=λMu). rIGA conserves the desirable properties of maximum-continuity isogeometric analysis (IGA) while it reduces the solution cost by adding zero-continuity basis functions, which decrease the matrix connectivity. As a result, rIGA enriches the approximation space and reduces the interconnection between degrees of freedom. We compare computational costs of rIGA versus those of IGA when employing a Lanczos eigensolver with a shift-and-invert spectral transformation. When all eigenpairs within a given interval [λs,λe] are of interest, we select several shifts σk∈[λs,λe] using a spectrum slicing technique. For each shift σk, the factorization cost of the spectral transformation matrix K−σkM controls the total computational cost of the eigensolution. Several multiplications of the operator matrix (K−σkM)−1M by vectors follow this factorization. Let p be the polynomial degree of the basis functions and assume that IGA has maximum continuity of p−1. When using rIGA, we introduce C0 separators at certain element interfaces to minimize the factorization cost. For this setup, our theoretical estimates predict computational savings to compute a fixed number of eigenpairs of up to O(p2) in the asymptotic regime, that is, large problem sizes. Yet, our numerical tests show that for moderate-size eigenproblems, the total observed computational cost reduction is O(p). In addition, rIGA improves the accuracy of every eigenpair of the first N0 eigenvalues and eigenfunctions, where N0 is the total number of modes of the original maximum-continuity IGA discretization.

GPU-acceleration of the ELPA2 distributed eigensolver for dense symmetric and hermitian eigenproblems

The solution of eigenproblems is often a key computational bottleneck that limits the tractable system size of numerical algorithms, among them electronic structure theory in chemistry and in condensed matter physics. Large eigenproblems can easily exceed the capacity of a single compute node, thus must be solved on distributed-memory parallel computers. We here present GPU-oriented optimizations of the ELPA two-stage tridiagonalization eigensolver (ELPA2). On top of cuBLAS-based GPU offloading, we add a CUDA kernel to speed up the back-transformation of eigenvectors, which can be the computationally most expensive part of the two-stage tridiagonalization algorithm. We benchmark the performance of this GPU-accelerated eigensolver on two hybrid CPU–GPU architectures, namely a compute cluster based on Intel Xeon Gold CPUs and NVIDIA Volta GPUs, and the Summit supercomputer based on IBM POWER9 CPUs and NVIDIA Volta GPUs. Consistent with previous benchmarks on CPU-only architectures, the GPU-accelerated two-stage solver exhibits a parallel performance superior to the one-stage counterpart. Finally, we demonstrate the performance of the GPU-accelerated eigensolver developed in this work for routine semi-local KS-DFT calculations comprising thousands of atoms.

Chiral two-dimensional periodic blocky materials with elastic interfaces: Auxetic and acoustic properties

Two novel chiral block lattice topologies are here conceived having interesting auxetic and acoustic behavior. The architectured chiral material is made up of a periodic repetition of square or hexagonal rigid and heavy blocks connected by linear elastic interfaces, whose chirality results from an equal rotation of the blocks with respect to the line connecting their centroids. The governing equation of the Lagrangian model is derived and a hermitian eigenproblem is formulated to obtain the frequency band structure. An equivalent micropolar continuum is analytically derived through a standard continualization approach in agreement with the procedure proposed by Bacigalupo and Gambarotta (2017) from which an approximation of the frequency spectrum is obtained. Moreover, the overall elastic moduli of the equivalent Cauchy continuum are obtained in closed form via a proper condensation procedure. The parametric analysis involving the overall elastic moduli of the Cauchy equivalent continuum model and the frequency band structure is carried out to catch the influence of the chirality angle and of the ratio between the tangential and normal stiffness of the interface. Finally, it is shown how chirality and interface stiffness may affect strong auxeticity and how the equivalent micropolar model provides dispersion curves in excellent agreement with the current ones for a wide range of the wave vector magnitude.

ChASE

Solving dense Hermitian eigenproblems arranged in a sequence with direct solvers fails to take advantage of those spectral properties that are pertinent to the entire sequence and not just to the single problem. When such features take the form of correlations between the eigenvectors of consecutive problems, as is the case in many real-world applications, the potential benefit of exploiting them can be substantial. We present the Chebyshev Accelerated Subspace iteration Eigensolver (ChASE), a modern algorithm and library based on subspace iteration with polynomial acceleration. Novel to ChASE is the computation of the spectral estimates that enter in the filter and an optimization of the polynomial degree that further reduces the necessary floating-point operations. ChASE is written in C++ using the modern software engineering concepts that favor a simple integration in application codes and a straightforward portability over heterogeneous platforms. When solving sequences of Hermitian eigenproblems for a portion of their extremal spectrum, ChASE greatly benefits from the sequence’s spectral properties and outperforms direct solvers in many scenarios. The library ships with two distinct parallelization schemes, supports execution over distributed GPUs, and is easily extensible to other parallel computing architectures.

Non-linear Least-Squares Optimization of Rational Filters for the Solution of Interior Hermitian Eigenvalue Problems

Rational filter functions can be used to improve convergence of contour-based eigensolvers, a popular family of algorithms for the solution of the interior eigenvalue problem. We present a framework for the optimization of rational filters based on a non-convex weighted Least-Squares scheme. When used in combination with a contour based eigensolvers library, our filters out-perform existing ones on a large and representative set of benchmark problems. This work provides a detailed description of: (1) a set up of the optimization process that exploits symmetries of the filter function for Hermitian eigenproblems, (2) a formulation of the gradient descent and Levenberg-Marquardt algorithms that exploits the symmetries, (3) a method to select the starting position for the optimization algorithms that reliably produces effective filters, (4) a constrained optimization scheme that produces filter functions with specific properties that may be beneficial to the performance of the eigensolver that employs them.

On local quadratic convergence of inexact simplified Jacobi–Davidson method for interior eigenpairs of Hermitian eigenproblems

For the Hermitian eigenproblems, under proper assumption on an initial approximation to the desired eigenvector, we prove local quadratic convergence of the inexact simplified Jacobi–Davidson method when the involved relaxed correction equation is solved by a standard Krylov subspace iteration, which particularly leads to local cubic convergence when the relaxed correction equation is solved to a prescribed precision proportional to the norm of the current residual. These results are valid for the interior as well as the extreme eigenpairs of the Hermitian eigenproblem and, hence, generalize the results by Bai and Miao (2017) from the extreme eigenpairs to the interior ones.

On local quadratic convergence of inexact simplified Jacobi–Davidson method

For the Hermitian eigenproblems, we prove local quadratic convergence of the inexact simplified Jacobi–Davidson method when the involved relaxed correction equation is solved by a standard Krylov subspace iteration. This method then shows local cubic convergence rate when the relaxed correction equation is solved to a prescribed precision proportional to the norm of the current residual. As a by-product, we obtain local cubic convergence of the simplified Jacobi–Davidson method. These results significantly improve the existing ones that show only local linear convergence for the inexact simplified Jacobi–Davidson method, which lead to local quadratic convergence for the simplified Jacobi–Davidson method when the tolerance of the inexact solve is particularly set to be zero. Numerical experiments confirm these theoretical results.

Preconditioned eigensolvers for large-scale nonlinear Hermitian eigenproblems with variational characterizations. I. Extreme eigenvalues

Efficient computation of extreme eigenvalues of large-scale linear Hermitian eigenproblems can be achieved by preconditioned conjugate gradient (PCG) methods. In this paper, we study PCG methods for computing extreme eigenvalues of nonlinear Hermitian eigenproblems of the form T ( λ ) v = 0 T(\lambda )v=0 that admit a nonlinear variational principle. We investigate some theoretical properties of a basic CG method, including its global and asymptotic convergence. We propose several variants of single-vector and block PCG methods with deflation for computing multiple eigenvalues, and compare them in arithmetic and memory cost. Variable indefinite preconditioning is shown to be effective to accelerate convergence when some desired eigenvalues are not close to the lowest or highest eigenvalue. The efficiency of variants of PCG is illustrated by numerical experiments. Overall, the locally optimal block preconditioned conjugate gradient (LOBPCG) is the most efficient method, as in the linear setting.

Preconditioned Eigensolvers for Large-Scale Nonlinear Hermitian Eigenproblems with Variational Characterizations. II. Interior Eigenvalues

We consider the solution of large-scale nonlinear algebraic Hermitian eigenproblems of the form $T(\lambda)v=0$ that admit a variational characterization of eigenvalues. These problems arise in a variety of applications and are generalizations of linear Hermitian eigenproblems $Av\!=\!\lambda Bv$. In this paper, we propose a preconditioned locally minimal residual (PLMR) method for efficiently computing interior eigenvalues of problems of this type. We discuss the development of search subspaces, preconditioning, and eigenpair extraction procedures based on the refined Rayleigh--Ritz projection. Extension to the block methods is presented, and a moving-window-style soft deflation is described. Numerical experiments demonstrate that PLMR methods provide a rapid and robust convergence toward interior eigenvalues. The approach is also shown to be efficient and reliable for computing a large number of extreme eigenvalues, dramatically outperforming standard preconditioned conjugate gradient methods.

Improved Accuracy and Parallelism for MRRR-Based Eigensolvers---A Mixed Precision Approach

The real symmetric tridiagonal eigenproblem is of outstanding importance in numerical computations; it arises frequently as part of eigensolvers for standard and generalized dense Hermitian eigenproblems that are based on a reduction to real tridiagonal form. For its solution, the algorithm of multiple relatively robust representations (MRRR) is among the fastest methods. Although fast, the solvers based on MRRR do not deliver the same accuracy as competing methods like Divide & Conquer or the QR algorithm. In this paper, we demonstrate that the use of mixed precisions leads to improved accuracy of MRRR-based eigensolvers with limited or no performance penalty. As a result, we obtain eigensolvers that are not only as accurate as or more accurate than the best available methods, but also---under most circumstances---faster and more scalable than the competition.

High-Performance Solvers for Dense Hermitian Eigenproblems

We introduce a new collection of solvers---subsequently called EleMRRR---for large-scale dense Hermitian eigenproblems. EleMRRR solves various types of problems: generalized, standard, and tridiagonal eigenproblems. Among these, the last is of particular importance as it is a solver in its own right, as well as the computational kernel for the first two; we present a fast and scalable tridiagonal solver based on the algorithm of multiple relatively robust representations---referred to as PMRRR. Like the other EleMRRR solvers, PMRRR is part of the freely available Elemental library, and is designed to fully support both message-passing and multithreading parallelism. As a result, the solvers are equally effective in message-passing environments with and without multithreading. We conducted a thorough performance study of EleMRRR and ScaLAPACK's solvers on two supercomputers. Such a study, performed with up to $8{,}192$ cores, provides precise guidelines for assembling the fastest solver within the ScaLAPACK framework; it also indicates that EleMRRR outperforms even the fastest solvers built from ScaLAPACK's components.

Divide and Conquer on Hybrid GPU-Accelerated Multicore Systems

With the raw computing power of graphics processing units (GPUs) being more widely available in commodity multicore systems, there is an imminent need to harness their power for important numerical libraries such as LAPACK. In this paper, we consider the solution of dense symmetric and Hermitian eigenproblems by the LAPACK divide and conquer algorithm on such modern heterogeneous systems. We focus on how to make the best use of the individual strengths of the massively parallel manycore GPUs and multicore CPUs. The resulting algorithm overcomes performance bottlenecks faced by current implementations that are optimized for a homogeneous multicore. On a dual socket quad-core Intel Xeon 2.33 GHz with an NVIDIA GTX 280 GPU, we typically obtain up to about a tenfold improvement in performance for the complete dense problem. The techniques described here thus represent an example of how to develop numerical software to efficiently use heterogeneous architectures. As heterogeneity becomes more common in the architecture design, the significance of and need for this work are expected to grow.

MR3-SMP: A symmetric tridiagonal eigensolver for multi-core architectures

The computation of eigenvalues and eigenvectors of symmetric tridiagonal matrices arises frequently in applications; often as one of the steps in the solution of Hermitian and symmetric eigenproblems. While several accurate and efficient methods for the tridiagonal eigenproblem exist, their corresponding implementations usually target uni-processors or large distributed memory systems. Our new eigensolver MR^3-SMP is instead specifically designed for multi-core and many-core general purpose processors, which today have effectively replaced uni-processors. We show that in most cases MR^3-SMP is faster and achieves better speedups than state-of-the-art eigensolvers for uni-processors and distributed-memory systems.

A quadrature-based eigensolver with a Krylov subspace method for shifted linear systems for Hermitian eigenproblems in lattice QCD

We consider a quadrature-based eigensolver to find eigenpairs of Hermitian matrices arising in lattice quantum chromodynamics. To reduce the computational cost for finding eigenpairs of such Hermitian matrices, we propose a new technique for solving shifted linear systems with complex shifts by means of the shifted CG method. Furthermore using integration paths along horizontal lines corresponding to the real axis of the complex plane, the number of iterations for the shifted CG method is also reduced. Some numerical experiments illustrate the accuracy and efficiency of the proposed method by comparison with a conventional method.

Two-sided Grassmann–Rayleigh quotient iteration

The two-sided Rayleigh quotient iteration proposed by Ostrowski computes a pair of corresponding left–right eigenvectors of a matrix C. We propose a Grassmannian version of this iteration, i.e., its iterates are pairs of p-dimensional subspaces instead of one-dimensional subspaces in the classical case. The new iteration generically converges locally cubically to the pairs of left–right p-dimensional invariant subspaces of C. Moreover, Grassmannian versions of the Rayleigh quotient iteration are given for the generalized Hermitian eigenproblem, the Hamiltonian eigenproblem and the skew-Hamiltonian eigenproblem.

Controlling Inner Iterations in the Jacobi–Davidson Method

The Jacobi-Davidson method is an eigenvalue solver which uses an inner-outer scheme. In the outer iteration one tries to approximate an eigenpair while in the inner iteration a linear system has to be solved, often iteratively, with the ultimate goal to make progress for the outer loop. In this paper we prove a relation between the residual norm of the inner linear system and the residual norm of the eigenvalue problem. We show that the latter may be estimated inexpensively during the inner iterations. On this basis, we propose a stopping strategy for the inner iterations to maximally exploit the strengths of the method. These results extend previous results obtained for the special case of Hermitian eigenproblems with the conjugate gradient or the symmetric QMR method as inner solver. The present analysis applies to both standard and generalized eigenproblems, does not require symmetry, and is compatible with most iterative methods for the inner systems. It can also be extended to other types of inner-outer eigenvalue solvers, such as inexact inverse iteration or inexact Rayleigh quotient iteration. The effectiveness of our approach is illustrated by a few numerical experiments, including the comparison of a standard Jacobi-Davidson code with the same code enhanced by our stopping strategy.

Detecting and Solving Hyperbolic Quadratic Eigenvalue Problems

Hyperbolic quadratic matrix polynomials $Q(\lambda) = \lambda^2 A + \lambda B + C$ are an important class of Hermitian matrix polynomials with real eigenvalues, among which the overdamped quadratics are those with nonpositive eigenvalues. Neither the definition of overdamped nor any of the standard characterizations provides an efficient way to test if a given Q has this property. We show that a quadratically convergent matrix iteration based on cyclic reduction, previously studied by Guo and Lancaster, provides necessary and sufficient conditions for Q to be overdamped. For weakly overdamped Q the iteration is shown to be generically linearly convergent with constant at worst $1/2$, which implies that the convergence of the iteration is reasonably fast in almost all cases of practical interest. We show that the matrix iteration can be implemented in such a way that when overdamping is detected a scalar $\mu<0$ is provided that lies in the gap between the n largest and n smallest eigenvalues of the $n \times n$ quadratic eigenvalue problem (QEP) $Q(\lambda)x = 0$. Once such a $\mu$ is known, the QEP can be solved by linearizing to a definite pencil that can be reduced, using already available Cholesky factorizations, to a standard Hermitian eigenproblem. By incorporating an initial preprocessing stage that shifts a hyperbolic Q so that it is overdamped, we obtain an efficient algorithm that identifies and solves a hyperbolic or overdamped QEP maintaining symmetry throughout and guaranteeing real computed eigenvalues.

Eigenvector approximation leading to exponential speedup of quantum eigenvalue calculation.

We present an efficient method for preparing the initial state required by the eigenvalue approximation quantum algorithm of Abrams and Lloyd. Our method can be applied when solving continuous Hermitian eigenproblems, e.g., the Schrödinger equation, on a discrete grid. We start with a classically obtained eigenvector for a problem discretized on a coarse grid, and we efficiently construct, quantum mechanically, an approximation of the same eigenvector on a fine grid. We use this approximation as the initial state for the eigenvalue estimation algorithm, and show the relationship between its success probability and the size of the coarse grid.

Scaled total least squares fundamentals

The standard approaches to solving overdetermined linear systems $Bx \approx c$ construct minimal corrections to the vector c and/or the matrix B such that the corrected system is compatible. In ordinary least squares (LS) the correction is restricted to c, while in data least squares (DLS) it is restricted to B. In scaled total least squares (STLS) [22], corrections to both c and B are allowed, and their relative sizes depend on a real positive parameter $\gamma$ . STLS unifies several formulations since it becomes total least squares (TLS) when $\gamma=1$ , and in the limit corresponds to LS when $\gamma\rightarrow 0$ , and DLS when $\gamma\rightarrow \infty$ . This paper analyzes a particularly useful formulation of the STLS problem. The analysis is based on a new assumption that guarantees existence and uniqueness of meaningful STLS solutions for all parameters $\gamma >0$ . It makes the whole STLS theory consistent. Our theory reveals the necessary and sufficient condition for preserving the smallest singular value of a matrix while appending (or deleting) a column. This condition represents a basic matrix theory result for updating the singular value decomposition, as well as the rank-one modification of the Hermitian eigenproblem. The paper allows complex data, and the equivalences in the limit of STLS with DLS and LS are proven for such data. It is shown how any linear system $Bx \approx c$ can be reduced to a minimally dimensioned core system satisfying our assumption. Consequently, our theory and algorithms can be applied to fully general systems. The basics of practical algorithms for both the STLS and DLS problems are indicated for either dense or large sparse systems. Our assumption and its consequences are compared with earlier approaches.

The Padé-Rayleigh-Ritz method for solving large hermitian eigenproblems

We make use of the Pade approximants and the Krylov sequencex, Ax,,...,A m−1 x in the projection methods to compute a few Ritz values of a large hermitian matrixA of ordern. This process consists in approaching the poles ofR x(λ)=((I−λA)−1 x,x), the mean value of the resolvant ofA, by those of [m−1/m]Rx(λ), where [m−1/m]Rx(λ) is the Pade approximant of orderm of the functionR x(λ). This is equivalent to approaching some eigenvalues ofA by the roots of the polynomial of degreem of the denominator of [m−1/m]Rx(λ). This projection method, called the Pade-Rayleigh-Ritz (PRR) method, provides a simple way to determine the minimum polynomial ofx in the Krylov subspace methods for the symmetrical case. The numerical stability of the PRR method can be ensured if the projection subspacem is “sufficiently” small. The mainly expensive portion of this method is its projection phase, which is composed of the matrix-vector multiplications and, consequently, is well suited for parallel computing. This is also true when the matrices are sparse, as recently demonstrated, especially on massively parallel machines. This paper points out a relationship between the PRR and Lanczos methods and presents a theoretical comparison between them with regard to stability and parallelism. We then try to justify the use of this method under some assumptions.