# Using Graph Partitioning for Scalable Distributed Quantum Molecular Dynamics
## 1 Introduction  
The physical movements of multi-body systems on an atomistic level is at the core of molecular dynamics (MD) simulations. Those dynamics take place at the femtosecond (\(10^{-15}\) s) time scale and they are incorporated in a larger simulation which typically is of the order of pico- to nanoseconds (\(10^{-12}\) to \(10^{-9}\) s). A simple way to conduct MD simulations is to derive all forces from the potential energy surface for all interacting particles, and to compute molecular trajectories for the multi-particle system by solving [PERSON]'s equations numerically. In quantum-based molecular dynamics (QMD) simulations, the electronic structure is based on an underlying quantum mechanical description, from which interatomic forces are calculated. Several QMD methods are published in the literature for a variety of materials systems. So-called first principle methods can simulate a few hundred atoms over a picosecond range. Important approaches of this type include Hartree-Fock or density functional theory. Semiempirical methods such as self-consistent tight-binding techniques increase the applicability to systems of several thousand atoms. In contrast to regular first principle methods, approximate methods are often two to three orders of magnitude faster while still capturing the quantum mechanical behavior (for example, charge transfer, bond formation, excitations, and quantum size effects) of the system. One of the most efficient and widely used methods is density functional-based self-consistent tight-binding theory [1; 2; 3]. In this approach, the main computational effort stems from thediagonalization of a matrix, the so-called _Hamiltonian matrix_, which encodes the electronic energy of the system. The Hamiltonian matrix is needed to construct the _density matrix_ describing the electronic structure of the system. Before evaluating the forces at each time step of a QMD simulation, the self-consistent construction of the density matrix is carried out. Computing the density matrix requires matrix diagonalization with a computational cost of \(O(N^{3})\), where \(N\) is the dimension of the Hamiltonian matrix. This makes diagonalization only viable for small system sizes. For this reason, the last two decades have seen the development of several linear runtime (i.e., \(O(N)\)) algorithms. One such linear runtime approach is based on a recursive polynomial expansion of the density matrix [4]. A linear scaling with the system size for non-metallic systems is achieved by the sparse-matrix second-order spectral projection (SM-SP2) algorithm. On dense or sparse matrices, SM-SP2 competes with or outperforms regular diagonalization schemes with respect to both speed and accuracy [5]. SM-SP2 uses the expression  
\[D=\lim_{i\to\infty}f_{i}[f_{i-1}[\dots f_{0}[X_{0}]\dots]] \tag{1}\]  
to compute the density matrix \(D\) from the Hamiltonian \(H\), where \(f_{i}(X_{i})\) is a quadratic function (either \(X_{i}^{2}\) or \(2X_{i}-X_{i}^{2}\), depending on \(\mathit{Tr}(X_{i})\) or \(\mathit{Tr}(X_{i+1})\)), and the initial matrix \(X_{0}\) is a linearly modified version of \(H\). Usually, 20-30 iterations suffice to obtain a close approximation of \(D\). Additionally, thresholding is applied to further reduce the computational complexity, where small nonzero elements of the matrix (typically between \(10^{-5}\) to \(10^{-7}\)) are set to zero. Thresholding the matrix elements of the Hamiltonian is essential to arrive at a fast way to calculate the density matrix: without thresholding, the resulting matrix would quickly become dense due to fill-in. The cost of computing a matrix polynomial \(P\) (which is mainly due to the squaring of a matrix) dominates the computational cost of the SM-SP2 algorithm. Since we are interested in performing a large number of time steps (of the order of \(10^{4}\)-\(10^{6}\)) in a typical QMD simulation, we need to parallelize the evaluation of the matrix polynomials in order to keep the wall-clock time low.