Recursive Factorization of the Inverse Overlap Matrix in Linear-Scaling Quantum Molecular Dynamics Simulations

Abstract

We present a reduced complexity algorithm to compute the inverse overlap factors required to solve the generalized eigenvalue problem in a quantum-based molecular dynamics (MD) simulation. Our method is based on the recursive, iterative refinement of an initial guess of Z (inverse square root of the overlap matrix S). The initial guess of Z is obtained beforehand by using either an approximate divide-and-conquer technique or dynamical methods, propagated within an extended Lagrangian dynamics from previous MD time steps. With this formulation, we achieve long-term stability and energy conservation even under the incomplete, approximate, iterative refinement of Z. Linear-scaling performance is obtained using numerically thresholded sparse matrix algebra based on the ELLPACK-R sparse matrix data format, which also enables efficient shared-memory parallelization. As we show in this article using self-consistent density-functional-based tight-binding MD, our approach is faster than conventional methods based on the diagonalization of overlap matrix S for systems as small as a few hundred atoms, substantially accelerating quantum-based simulations even for molecular structures of intermediate size. For a 4158-atom water-solvated polyalanine system, we find an average speedup factor of 122 for the computation of Z in each MD step.

Canonical-ensemble extended Lagrangian Born–Oppenheimer molecular dynamics for the linear scaling density functional theory

“…In the extended Lagrangian BOMD method, on the other hand, the conserved quantity remained stable with no systematic energy drift ( figure 1(b)), in agreement with previous work [12][13][14][15][16][17][18][19][20]31]. The number of iterations in the DMM procedure was moderate (the maximum number of iterations for ε = 10 −5 was 7 cycles per time step).…”

Localized inverse factorization

“…L 2 gm (Z; E) = L gm (Z; E), and thus has bounded powers. Let (38) h be the mapping associated with (31). Let Z * SZ = I.…”

“…However, the multifrontal approach could also be combined with the AINV algorithm for direct computation of the inverse Cholesky factor. A starting guess for iterative refinement is in [31] built up from the inverse factors of overlapping principal submatrices. The convergence has not been proven but this approach could possibly lead to a lower initial factorization error and reduced number of iterative refinement iterations compared to our binary partition at the expense of more computations to construct the starting guess.…”

“…Hermitian. Thus,(31) is stable if non-Hermitian perturbations of δ are disallowed, which can simply be achieved by storing only the upper triangle of δ.The stability properties of the matrix iterations above are demonstrated inFigure 1. Again, we follow[13] and apply the iterations to the Wilson matrix…”

Three‐body expansion of the fragment molecular orbital method combined with density‐functional tight‐binding

“…The diagonalization step is the bottleneck, resulting in the cubic scaling on the system size. DFTB has been combined with several low scaling methods, hybrid approaches, and FMO, eliminating this bottleneck.…”