An efficient approach for the calculation of Franck-Condon integrals of large molecules

A general and efficient approach for the calculation of Franck-Condon integrals (FCIs) of large molecules is presented. In a first step, by exploiting the diagonally dominant and sparse structure of the Duschinsky matrix, a model system is constructed for which the Duschinsky matrix takes a block-diagonal form. For each of these blocks separately, the FCIs are calculated discarding all below a certain threshold. From those integrals retained the FCIs of the model system are obtained by simple multiplication. These serve as an estimate for the FCIs of the exact system which are calculated for those integrals which lie above a certain threshold. By systematically decreasing the threshold, the simulation can be reliably converged to the exact result with an arbitrary accuracy. Using this scheme, a considerable reduction of the number of FCIs which have to be calculated is achieved which leads to an improved scaling behavior of the computational effort with system size. The approach has been tested thoroughly for a set of molecules including difficult cases. For the larger systems a speedup of up to three orders of magnitude compared to an exact calculation is observed while the errors can be kept negligible. With this approach accurate calculations of FCIs are feasible also for large molecules encountered in "real-life" chemistry, especially biochemistry and material science.     

Systematic sparse matrix error control for linear scaling electronic structure calculations

Efficient truncation criteria used in multiatom blocked sparse matrix operations for ab initio calculations are proposed. As system size increases, so does the need to stay on top of errors and still achieve high performance. A variant of a blocked sparse matrix algebra to achieve strict error control with good performance is proposed. The presented idea is that the condition to drop a certain submatrix should depend not only on the magnitude of that particular submatrix, but also on which other submatrices that are dropped. The decision to remove a certain submatrix is based on the contribution the removal would cause to the error in the chosen norm. We study the effect of an accumulated truncation error in iterative algorithms like trace correcting density matrix purification. One way to reduce the initial exponential growth of this error is presented. The presented error control for a sparse blocked matrix toolbox allows for achieving optimal performance by performing only necessary operations needed to maintain the requested level of accuracy.     

Hartree-Fock calculations with linearly scaling memory usage

We present an implementation of a set of algorithms for performing Hartree-Fock calculations with resource requirements in terms of both time and memory directly proportional to the system size. In particular, a way of directly computing the Hartree-Fock exchange matrix in sparse form is described which gives only small addressing overhead. Linear scaling in both time and memory is demonstrated in benchmark calculations for system sizes up to 11 650 atoms and 67 204 Gaussian basis functions on a single computer with 32 Gbytes of memory. The sparsity of overlap, Fock, and density matrices as well as band gaps are also shown for a wide range of system sizes, for both linear and three-dimensional systems.     

Representations of the symmetric group generated by projected spin functions: A graphical approach

The representation matrices generated by the projected spin functions have some very interesting properties. All the matrix elements are integers and they are quite sparse. A very efficient algorithm is presented for the calculation of these representation matrices based on a graphical approach and a new indexing scheme for representation of primitive spin functions is introduced. Test calculations show that the method is very fast and suited for calculations on vector computers. © 1996 John Wiley & Sons, Inc.     

Gaussian basis sets for accurate calculations on molecular systems in gas and condensed phases

We present a library of Gaussian basis sets that has been specifically optimized to perform accurate molecular calculations based on density functional theory. It targets a wide range of chemical environments, including the gas phase, interfaces, and the condensed phase. These generally contracted basis sets, which include diffuse primitives, are obtained minimizing a linear combination of the total energy and the condition number of the overlap matrix for a set of molecules with respect to the exponents and contraction coefficients of the full basis. Typically, for a given accuracy in the total energy, significantly fewer basis functions are needed in this scheme than in the usual split valence scheme, leading to a speedup for systems where the computational cost is dominated by diagonalization. More importantly, binding energies of hydrogen bonded complexes are of similar quality as the ones obtained with augmented basis sets, i.e., have a small (down to 0.2 kcal/mol) basis set superposition error, and the monomers have dipoles within 0.1 D of the basis set limit. However, contrary to typical augmented basis sets, there are no near linear dependencies in the basis, so that the overlap matrix is always well conditioned, also, in the condensed phase. The basis can therefore be used in first principles molecular dynamics simulations and is well suited for linear scaling calculations.     

How to choose one-dimensional basis functions so that a very efficient multidimensional basis may be extracted from a direct product of the one-dimensional functions: energy levels of coupled systems with as many as 16 coordinates

In this paper we propose a scheme for choosing basis functions for quantum dynamics calculations. Direct product bases are frequently used. The number of direct product functions required to converge a spectrum, compute a rate constant, etc., is so large that direct product calculations are impossible for molecules or reacting systems with more than four atoms. It is common to extract a smaller working basis from a huge direct product basis by removing some of the product functions. We advocate a build and prune strategy of this type. The one-dimensional (1D) functions from which we build the direct product basis are chosen to satisfy two conditions: (1) they nearly diagonalize the full Hamiltonian matrix; (2) they minimize off-diagonal matrix elements that couple basis functions with diagonal elements close to those of the energy levels we wish to compute. By imposing these conditions we increase the number of product functions that can be removed from the multidimensional basis without degrading the accuracy of computed energy levels. Two basic types of 1D basis functions are in common use: eigenfunctions of 1D Hamiltonians and discrete variable representation (DVR) functions. Both have advantages and disadvantages. The 1D functions we propose are intermediate between the 1D eigenfunction functions and the DVR functions. If the coupling is very weak, they are very nearly 1D eigenfunction functions. As the strength of the coupling is increased they resemble more closely DVR functions. We assess the usefulness of our basis by applying it to model 6D, 8D, and 16D Hamiltonians with various coupling strengths. We find approximately linear scaling.     

A regionally contracted multireference configuration interaction method: general theory and results of an incremental version

This work proposes to take benefit of the localizability of both occupied and virtual inactive molecular orbitals (MOs) in the context of complete active space singles and doubles configuration interaction (CAS-SDCI). The doubly occupied MOs are partitioned into blocks, or regions, corresponding to a subset of adjacent bonds and lone pairs. The localized virtual MOs are attributed to these regions from a spatial criterion. Then a series of limited post-CAS-CI calculations is performed, using the same reference space, one for each block, and then one per pair of blocks. From these independent CI calculations contracted external functions are defined for each block or for each pair of blocks, and for each state. A general multistate formalism is proposed, the CI matrix being expressed in the space defined by the CAS and the contracted functions. Preliminary numerical studies, resting on the evaluation of single-block and two-block contributions to the dynamical correlation energy of each state, are presented. Provided that size-consistency corrections are taken into account the results of the procedure are shown to be in excellent agreement with those of the nonpartitioned post-CAS-CI. The computational benefits of this evidently parallelizable procedure are underlined.     

New algorithm for tensor contractions on multi‐core CPUs,GPUs, and accelerators enables CCSD and EOM‐CCSD calculations with over 1000 basis functions on a single compute node

A new hardware‐agnostic contraction algorithm for tensors of arbitrary symmetry and sparsity is presented. The algorithm is implemented as a stand‐alone open‐source code libxm . This code is also integrated with general tensor library libtensor and with the Q‐Chem quantum‐chemistry package. An overview of the algorithm, its implementation, and benchmarks are presented. Similarly to other tensor software, the algorithm exploits efficient matrix multiplication libraries and assumes that tensors are stored in a block‐tensor form. The distinguishing features of the algorithm are: (i) efficient repackaging of the individual blocks into large matrices and back, which affords efficient graphics processing unit (GPU)‐enabled calculations without modifications of higher‐level codes; (ii) fully asynchronous data transfer between disk storage and fast memory. The algorithm enables canonical all‐electron coupled‐cluster and equation‐of‐motion coupled‐cluster calculations with single and double substitutions (CCSD and EOM‐CCSD) with over 1000 basis functions on a single quad‐GPU machine. We show that the algorithm exhibits predicted theoretical scaling for canonical CCSD calculations, O (N ⁶), irrespective of the data size on disk. © 2017 Wiley Periodicals, Inc.     

Multidimensional time-dependent discrete variable representations in multiconfiguration Hartree calculations

In the multiconfiguration time-dependent Hartree (MCTDH) approach, the wave function is expanded in time-dependent basis functions, called single-particle functions, to increase the efficiency of the wave-packet propagation. The correlation discrete variable representation (CDVR) approach, which is based on a time-dependent discrete variable representation (DVR), can be employed to evaluate matrix elements of the potential energy. The efficiency of the MCTDH method can be further enhanced by using multidimensional single-particle functions. However, up to now the CDVR approach could not be used in MCTDH calculations employing multidimensional single-particle functions, since this would require a general multidimensional non-direct-product DVR scheme. Recently, Dawes and Carrington presented a practical scheme to implement general non-direct-product multidimensional DVRs [R. Dawes and T. Carrington, Jr., J. Chem. Phys. 121, 726 (2004)]. The present work utilizes their scheme in the MCTDH/CDVR approach. The accuracy is tested using the photodissociation of NOCl as example. The results show that the CDVR scheme based on multidimensional time-dependent DVRs allows for an accurate evaluation of the potential in MCTDH calculations with multidimensional single-particle functions.     

Approaching the theoretical limit in periodic local MP2 calculations with atomic-orbital basis sets: the case of LiH

The atomic orbital basis set limit is approached in periodic correlated calculations for solid LiH. The valence correlation energy is evaluated at the level of the local periodic second order M?ller-Plesset perturbation theory (MP2), using basis sets of progressively increasing size, and also employing "bond"-centered basis functions in addition to the standard atom-centered ones. Extended basis sets, which contain linear dependencies, are processed only at the MP2 stage via a dual basis set scheme. The local approximation (domain) error has been consistently eliminated by expanding the orbital excitation domains. As a final result, it is demonstrated that the complete basis set limit can be reached for both HF and local MP2 periodic calculations, and a general scheme is outlined for the definition of high-quality atomic-orbital basis sets for solids.     

Direct inversion of the iterative subspace with contracted planewave basis functions

Ways to reduce the computational cost of periodic electronic structure calculations by using basis functions corresponding to linear combinations of planewaves have been examined recently. These contracted planewave (CPW) basis functions correspond to Fourier series representations of atom‐centered basis functions, and thus provide access to some beneficial properties of planewave (PW) and localized basis functions. This study reports the development and assessment of a direct inversion of the iterative subspace (DIIS) method that employs unique properties of CPW basis functions to efficiently converge electronic wavefunctions. This method relies on access to a PW‐based representation of the electronic structure to provide a means of efficiently evaluating matrix–vector products involving the application of the Fock matrix to the occupied molecular orbitals. These matrix–vector products are transformed into a form permitting the use of direct diagonalization techniques and DIIS methods typically employed with atom‐centered basis sets. The abilities of this method are assessed through periodic Hartree–Fock calculations of a range of molecules and solid‐state systems. The results show that the method reported in this study is approximately five times faster than CPW‐based calculations in which the entire Fock matrix is calculated. This method is also found to be weakly dependent upon the size of the basis set, thus permitting the use of larger CPW basis sets to increase variational flexibility with a minor impact on computational performance. © 2018 Wiley Periodicals, Inc.     

Bringing about matrix sparsity in linear-scaling electronic structure calculations

The performance of linear-scaling electronic structure calculations depends critically on matrix sparsity. This article gives an overview of different strategies for removal of small matrix elements, with emphasis on schemes that allow for rigorous control of errors. In particular, a novel scheme is proposed that has significantly smaller computational overhead compared with the Euclidean norm-based truncation scheme of Rubensson et al. (J Comput Chem 2009, 30, 974) while still achieving the desired asymptotic behavior required for linear scaling. Small matrix elements are removed while ensuring that the Euclidean norm of the error matrix stays below a desired value, so that the resulting error in the occupied subspace can be controlled. The efficiency of the new scheme is investigated in benchmark calculations for water clusters including up to 6523 water molecules. Furthermore, the foundation of matrix sparsity is investigated. This includes a study of the decay of matrix element magnitude with distance between basis function centers for different molecular systems and different methods. The studied methods include Hartree–Fock and density functional theory using both pure and hybrid functionals. The relation between band gap and decay properties of the density matrix is also discussed.     

Optimization of basis sets for isoelectronic series of closed-shell atoms in Hartree-Fock-Roothaan calculations

Optimization of the nonlinear parameters (orbital exponents) of basis functions in Hartree-Fock-Roothaan calculations may be canied out with a high degree of accuracy. A scheme using second-order methods is suggested for optimization of the orbital exponents of Slater type basis functions defining the ground state of closed-shell atoms. An exact fomula is derived for calculating the partial second derivatives of energy with respect to nonlinear parameters in temis of the density matrix. In bases of isoelectronic series, the orbital exponents are shown to be the linear functions of the charge of the ion nucleus. Optimization calculations are reported for He, Be, Ne, and Mg atoms and their isoelectronic series. Translated fromZhumal Struktumoi Khimii, Vol. 41, No. 2, pp. 217–228, March–April, 2000     

Using simultaneous diagonalization and trace minimization to make an efficient and simple multidimensional basis for solving the vibrational Schrodinger equation

In this paper we improve the product simultaneous diagonalization (SD) basis method we previously proposed [J. Chem. Phys. 122, 134101 (2005)] and applied to solve the Schrodinger equation for the motion of nuclei on a potential surface. The improved method is tested using coupled complicated Hamiltonians with as many as 16 coordinates for which we can easily find numerically exact solutions. In a basis of sorted products of one-dimensional (1D) SD functions the Hamiltonian matrix is nearly diagonal. The localization of the 1D SD functions for coordinate qc depends on a parameter we denote alphac. In this paper we present a trace minimization scheme for choosing alphac to nearly block diagonalize the Hamiltonian matrix. Near-block diagonality makes it possible to truncate the matrix without degrading the accuracy of the lowest energy levels. We show that in the sorted product SD basis perturbation theory works extremely well. The trace minimization scheme is general and easy to implement.     

Full optimization of Jastrow-Slater wave functions with application to the first-row atoms and homonuclear diatomic molecules

We pursue the development and application of the recently introduced linear optimization method for determining the optimal linear and nonlinear parameters of Jastrow-Slater wave functions in a variational Monte Carlo framework. In this approach, the optimal parameters are found iteratively by diagonalizing the Hamiltonian matrix in the space spanned by the wave function and its first-order derivatives, making use of a strong zero-variance principle. We extend the method to optimize the exponents of the basis functions, simultaneously with all the other parameters, namely, the Jastrow, configuration state function, and orbital parameters. We show that the linear optimization method can be thought of as a so-called augmented Hessian approach, which helps explain the robustness of the method and permits us to extend it to minimize a linear combination of the energy and the energy variance. We apply the linear optimization method to obtain the complete ground-state potential energy curve of the C(2) molecule up to the dissociation limit and discuss size consistency and broken spin-symmetry issues in quantum Monte Carlo calculations. We perform calculations for the first-row atoms and homonuclear diatomic molecules with fully optimized Jastrow-Slater wave functions, and we demonstrate that molecular well depths can be obtained with near chemical accuracy quite systematically at the diffusion Monte Carlo level for these systems.     

An efficient linear scaling procedure for constructing localized orbitals of large molecules based on the one-particle density matrix

We have developed a linear-scaling algorithm for obtaining the Boys localized molecular orbitals from the one-particle density matrix. The algorithm is made up of two steps: the Cholesky decomposition of the density matrix to obtain Cholesky molecular orbitals and the subsequent Boys localization process. Linear-scaling algorithms have been proposed to achieve linear-scaling calculations of these two steps, based on the sparse matrix technique and the locality of the Cholesky molecular orbitals. The present algorithm has been applied to compute the Boys localized orbitals in a number of systems including α-helix peptides, water clusters, and protein molecules. Illustrative calculations demonstrate that the computational time of obtaining Boys localized orbitals with the present algorithm is asymptotically linear with increasing the system size.     

Parallelization of the deMon2k code

The parallelization of the LCGTO-KS-DFT code deMon2k is presented. The parallelization of the three-center electron repulsion integrals, the numerical integration using a direct grid algorithm and the matrix multiplication and diagonalization are described. The efficiency of the parallelization is analyzed by selected benchmark calculations. It is shown that geometry optimizations of systems with more than 8,000 basis functions are feasible on cluster architectures.     

Direct energy functional minimization under orthogonality constraints

The direct energy functional minimization problem in electronic structure theory, where the single-particle orbitals are optimized under the constraint of orthogonality, is explored. We present an orbital transformation based on an efficient expansion of the inverse factorization of the overlap matrix that keeps orbitals orthonormal. The orbital transformation maps the orthogonality constrained energy functional to an approximate unconstrained functional, which is correct to some order in a neighborhood of an orthogonal but approximate solution. A conjugate gradient scheme can then be used to find the ground state orbitals from the minimization of a sequence of transformed unconstrained electronic energy functionals. The technique provides an efficient, robust, and numerically stable approach to direct total energy minimization in first principles electronic structure theory based on tight-binding, Hartree-Fock, or density functional theory. For sparse problems, where both the orbitals and the effective single-particle Hamiltonians have sparse matrix representations, the effort scales linearly with the number of basis functions N in each iteration. For problems where only the overlap and Hamiltonian matrices are sparse the computational cost scales as O(M2N), where M is the number of occupied orbitals. We report a single point density functional energy calculation of a DNA decamer hydrated with 4003 water molecules under periodic boundary conditions. The DNA fragment containing a cis-syn thymine dimer is composed of 634 atoms and the whole system contains a total of 12,661 atoms and 103,333 spherical Gaussian basis functions.