Wavefunction-less, density matrix-based approach to computational quantum chemistry is briefly discussed. Implementation of second-order M oller-Plesset Perturbation Method energy and dipole moment calculations within the new paradigm is presented. Efficiency and reliability of the method is analyzed

Kowalczyk, Krzysztof

Makowski, Marcin

Mazur, Grzegorz

Sumera, Jakub

English

Computing and Informatics (E-Journal - Institute of Informatics, SAS, Bratislava)

Computing and Informatics, Vol. 31, 2012, 665–673EFFICIENT ENERGY AND ELECTROSTATICPROPERTIES CALCULATIONS AT THE MP2 THEORYLEVEL: A CASE STUDY OF DENSITY MATRIX-BASEDCOMPUTATIONAL QUANTUM CHEMISTRYGrzegorz MazurDepartment of Computational Methods in ChemistryJagiellonian Universityul. R. Ingardena 3, 30-060 Kraków, Polande-mail: mazur@chemia.uj.edu.plMarcin MakowskiDepartment of Theoretical ChemistryJagiellonian Universityul. R. Ingardena 3, 30-060 Kraków, Polande-mail: makowskm@chemia.uj.edu.plJakub Sumera, Krzysztof KowalczykFaculty of ChemistryJagiellonian Universityul. R. Ingardena 3, 30-060 Kraków, PolandCommunicated by Jacek KitowskiAbstract. Wavefunction-less, density matrix-based approach to computationalquantum chemistry is briefly discussed. Implementation of second-order Møller-Plesset Perturbation Method energy and dipole moment calculations within thenew paradigm is presented. Efficiency and reliability of the method is analyzed.Keywords: Computational chemistry, MP2, Laplace transform, linear scaling666 G. Mazur, M. Makowski, J. Sumera, K. Kowalczyk1 INTRODUCTIONSignificantly growing demand for high quality computational chemistry results fornano-materials and bio-chemical systems can be recently observed. The standard,wavefunction based computational quantum chemistry methods often cannot ful-fill the expectations. The main reason is inherent high computational complexityof these methods directly stemming from the bad locality of molecular orbitals,which are the main building block of any wavefunction-based quantum-chemicalformalism. This gives rise to the paradigm shift towards alternative, purely den-sity matrix-based computational quantum chemistry methods which promise linearcomputational complexity. Supported by high-performance hardware and efficientalgorithms, the new approach may provide viable toolkit for treating nano- andbiosystems.The aim of the paper is twofold. We highlight the crucial algorithmic aspects ofdensity-matrix based approach, and study its efficiency and reliability. We choosethe implementation of the second order Møller-Plesset perturbation method (MP2)as an illustrative example.2 LAPLACE-TRANSFORM MP2To keep the paper self-contained, we briefly introduce the fundamental elements ofthe formalism we implemented and analyzed.The second order energy correction in the Møller-Plesset perturbation method(MP2) takes the form [1, 2]E2 = −occ∑ijvirt∑ab(ia|jb)[2(ia|jb)− (ib|ja)]ǫa + ǫb − ǫi − ǫj(1)where i, j denote occupied orbitals, a, b denote virtual orbitals, ǫ stands for orbitalenergy and (ia|jb) represents a two-electron integral in chemists’ notation.From the computational complexity point of view, the crucial step is the trans-formation of the two-electron integrals from the atomic orbital (AO) to the molecularorbital (MO) basis set(ia|jb) =N∑µνλσCµiCνaCλjCσb(µν|λσ) (2)where Greek letters denote AOs, C stands for the molecular orbital coefficient matrixand N represents the dimension of the Hilbert subspace spanned by the atomicorbitals. The complexity stems from the fact that the C matrix is dense. Hence,even though its AO basis counterpart shows sparse structure with only O(N 2) non-negligible elements, the integral tensor in the MO basis is dense. This results in thetransformation having time complexity of O(N 5) and memory complexity of O(N 4),which basically precludes calculations for large systems.Efficient Properties Calculations at the MP2 Theory Level 667The breakthrough in reducing computational complexity of MP2 calculationswas made by Häser [3]. Noticing that the formula for energy correction is a LaplacetransformL[f ](x) =∫∞0f(t)e−txdt (3)of a constant functionL[1](x) =∫∞0e−txdt =1x(4)he proposed an alternative formulation of the energy correctionE2 = −∫ occ∑ijvirt∑ab(ia|jb)[2(ia|jb)− (ib|ja)]e−t(ǫa+ǫb−ǫi−ǫj)dt. (5)Introducing energy-weighted molecular orbitals|i′〉 = e12tǫi|i〉 |a′〉 = e−12tǫa|a〉 (6)the energy correction can be expressed asE2 = −∫ occ∑ijvirt∑ab(i′a′|j ′b′)[2(i′a′|j ′b′)− (i′b′|j ′a′)]dt. (7)The main feature of the Laplace-transform-based MP2 energy expression is thatit is invariant with respect to unitary transformations [3]. This allows for choosinga basis in which the complexity of the calculations is optimal. And since the com-plexity is directly related to the locality of the orbitals, the most natural choice isthe atomic orbital basis set.Introducing orbital-energy-weighted pseudo-density matricesP µ′µ =occ∑iCµ′iCµieεit (8)P ν′ν =virt∑aCν′aCνae−εat (9)where the C matrix represents the canonical molecular orbitals, and defining trans-formed AO orbitals as|µ〉 =N∑νP µν |ν〉 |µ〉 =N∑νP µν |ν〉 (10)yields the workhorse of the Laplace-transform Atomic Orbital MP2 (LT-AO MP2)668 G. Mazur, M. Makowski, J. Sumera, K. Kowalczykmethod [3]E2 = −∫∞0e2(t)dt = −∫∞0N∑µνλσ(µν|λσ)[2(µν|λσ)− (µσ|λν)]dt. (11)The formula can be efficiently implemented owing to the sparsity of the pseudo-density matrices, strictly following the sparsity pattern of the density matrix. Na-mely, the number of non-zero elements of a density matrix in an insulator is scalinglinearly with the size of the system. Moreover, the non-zero elements may occur onlyin a relatively narrow diagonal band. The sparsity, complemented with standard pre-screening techniques results in a significant reduction in number of the transformedtwo-electron integrals which have to be calculated and stored in memory. This,in turn, allows for a relatively straightforward implementation of calculations hav-ing both time and memory complexity of O(N 2), where N stands for the numberof orbitals. Furthermore, stricter prescreening should lead to the linearly scalingmethod. Low complexity should be achievable despite the fact that the integrationin Equation (11) has to be performed numerically. (Analytical integration wouldbring us back to the canonical MP2 method.)While the LT-AO MP2 formulation by Häser allows to efficiently calculate theenergy correction, it has some drawbacks. The main one is that the definition ofthe pseudo-density matrices involves the canonical MO coefficients and orbital ener-gies. This prevents complementing the LT-AO MP2 method with linearly-scaling,diagonalization-free Hartree-Fock calculations. To make things worse, implemen-tation of even basic energy-derivative based molecular properties computations in-volves dense matrix of MO coefficients derivatives, significantly degrading the effi-ciency.The issue has been recently solved by Súrjan, who proved [4] that the LT-AOMP2 pseudo-density matrices may be expressed purely in terms of the density andFock matrices, without any reference to the MO coefficients or orbital energies. Theresulting expressions [4]P (t) = etPFP P (t) = e−tQFQ (12)where P and Q denote the electron and hole density matrices, are the basis ofdensity matrix-based LT-AO MP2 formalism.This formulation, while keeping all the advantages of the original LT-AO MP2method, allows for efficient calculations of the energy-derivatives based molecularproperties like dipole moment or polarizability.3 IMPLEMENTATIONCalculations of the LT-AO MP2 energy and dipole moment have been implementedin the Niedoida computational chemistry package [5].Efficient Properties Calculations at the MP2 Theory Level 669The main aspect of the algorithm is an efficient representation of the transformedtwo-electron integral tensor. In our implementation, a four-level tree is used. Theintegral values are stored in contiguous memory regions attached to the leaves ofthe tree. Such a structure allows for optimal memory consumption, and at the sametime maintains efficiency of the algebraic operations performed on the integrals.The process of tree construction involves multi-level prescreening procedure, whichallows for skipping branches containing only negligibly small integrals.To improve efficiency, the computational procedure was parallelized. We adopt-ed Single Program Multiple Data (SPMD) model implemented by means of theMessage Passing Interface (MPI). The scheduling of subtasks was implemented ac-cording to [6].4 RESULTSIn order to verify applicability of the LT-AO MP2 method, we performed a seriesof calculations for selected model systems, for which canonical MP2 results areavailable. We assessed both the accuracy and computational complexity of themethod.The numerical integration was performed in a transformed coordinate system.This way it was possible to make the integration range finite and at the same timeto reduce the variation of the integrand, allowing for good accuracy even with verylow number of quadrature points. Specifically, all LT-AO MP2 calculations reportedhere were performed using the coordinate transformation defined by [7]t(r) =∑k akrk(1− r)m. (13)The optimal parameter values m = 2, k = 3 . . .6, ak = 2.64,−6.09, 6.07,−2.7 werefound by minimizing the integration error of exponential function e−αx with varyingvalue of the α factor. As suggested in [7], we used the Euler-Maclaurin quadrature.All LT-AO MP2 calculations reported here employed 10-point grid.The integral prescreening was performed according to the Häser protocol [3].The threshold used for skipping negligible individual integral contributions was setto 10−8 a.u. The accuracy was calculated with respect to the reference results fromcanonical MP2 method as implemented in the GAMESS package [8].The first system under study is a linear chain of water molecules. While notvery realistic, such a model is well suited (and often used) for testing methodsapplicable to large systems. The other testing set consists of selected nucleobasesand their pairs, which allows for direct evaluation of applicability of LT-AO MP2 tocalculations for molecules of importance in biochemistry.The results collected in Figure 1, Tables 1 and 2 show that even with modestintegration grid size the error introduced by the LT-AO MP2 approximation is wellwithin the bounds of so called chemical accuracy.670 G. Mazur, M. Makowski, J. Sumera, K. Kowalczyk 1e-06 1e-05 0.0001 1  2  3  4  5  6  7  8  9  10Absolute error (au)System sizedipole momentenergyFig. 1. Accuracy of LT-AO MP2 energy and dipole moment calculations for linear chainsof water molecules in the 6-31G** basis set. System size denotes the number ofmolecules in chain.Nucleobase Error (10−3 a.u.)Adenine (A) 1.8Cytosine (C) 0.3Guanine (G) 1.7Thymine (T) 1.0A–T 7.0C–G 2.0Table 1. Absolute error of LT-AO MP2/6-31G** energy calculations with respect to thecanonical MP2 resultsNucleobaseError (10−4 a.u.)x y zAdenine (A) 1.5 0.2 0.5Cytosine (C) 4.0 0.0 1.6Guanine (G) 4.7 4.8 4.1Thymine (T) 3.3 2.8 0.9A–T 2.0 1.8 4.0C–G 5.7 0.4 0.7Table 2. Absolute error of LT-AO MP2/6-31G** dipole moment calculations with respectto the canonical MP2 resultsEfficient Properties Calculations at the MP2 Theory Level 671Basis setComplexityLT-AO MP2 Canonical MP2STO-3G 2.81 5.546-31G** 2.60 5.44Table 3. Overall computational complexity of energy calculations for linear chains of watermoleculesBasis setComplexityLT-AO MP2 Canonical MP2STO-3G 2.82 4.733-21G 2.49 4.806-31G** 2.65 –Table 4. Overall computational complexity of dipole moment calculations for linear chainsof water molecules. (For canonical MP2 calculations in the 6-31G** basis not enoughdata-points were collected to reliably fit the complexity.)Efficiency data in Tables 3 and 4 show that the computational complexity of theLT-AO MP2 method is more than two orders lower with respect to the canonicalcounterpart. For actual computational time, the crossover in the case of linear chainsof water molecules was observed at the system size of roughly 20 molecules.5 CONCLUSIONSWe have shown that the density matrix-based quantum chemistry stands up to theexpectations with respect to improved computational complexity, at the same timemaintaining adequate accuracy.While the new paradigm introduces some algorithmic issues which require fur-ther studies, it is already possible to obtain quadratic scaling with the system sizeat the MP2 level, and linear scaling seems achievable.REFERENCES[1] Møller, C.—Plesset, M. S.: Note on an Approximation Treatment for Many-Electron Systems. Phys Rev, Vol. 46, 1934, p. 618.[2] Helgaker, T.—Jorgensen, P.—Olsen, J.: Molecular Electronic-Structure The-ory. Wiley 2000.[3] Häser, M.: Moller-Plesset (MP2) Perturbation Theory for Large Molecules. TheorChim Acta, Vol. 87, 1993, p. 147.[4] Surjan, P.R.: The MP2 Energy As a Functional of the Hartree-Fock Density Ma-trix. Chem Phys Lett, Vol. 406, 2005, p. 318.[5] Mazur, G.—Makowski, M.: Development and Optimization of ComputationalChemistry Algorithms. Comput Informat, Vol. 28, 2009, p. 115.672 G. Mazur, M. Makowski, J. Sumera, K. Kowalczyk[6] Mazur, G.—Makowski, M.—Brela, M.: Effective Resource Allocation in Pa-rallel Quantum-Chemical Calculations. Comput Informat, Vol. 30, 2011, p. 761.[7] Ayala, P.Y.—Scuseria, G. E.: Linear Scaling Second-order Moller-Plesset The-ory in the Atomic Orbital Basis for Large Molecular Systems. Vol. 110, 1999, No. 8,p. 3660.[8] Schmidt, M.W.—Baldridge, K.K.—Boatz, J.A.—Elbert, S. T.—Gor-don, M. S.—Jensen, J.H.—Koseki, S.—Matsunaga, N.—Nguyen, K.A.—Su, S.—Windus, T. L.—Dupuis, M.—Montgomery, J.A.: General Atomic andMolecular Electronic Structure System. J Comput Chem, Vol. 14, 1993, p. 1347.Grzegorz Mazur is an Assistant Professor at the Departmentof Computational Methods in Chemistry of Jagiellonian Univer-sity in Cracow (Poland). He received his Ph.D. in chemistry in2001 from the same university. His current research interests in-clude theoretical description of the excited states properties andoptimization of quantum chemistry algorithms.MarcinMakowski is an Assistant Professor at the Departmentof Theoretical Chemistry of Jagiellonian University in Cracow(Poland). He received his M. Sc. degree in 2001 and his Ph.D. in2004, both in chemistry, from the same university, and his B. Sc.degree in computer science from the University of Mining andMetallurgy in Cracow in 2004. His current research interests in-clude theoretical molecular spectroscopy, linear scaling methodsin quantum chemistry and optimization of quantum chemistryalgorithms. He participates in several research projects relatedwith the development of theoretical chemistry formalisms andnumerical methodologies that allow to calculate electronic structure of large molecularsystems efficiently.Jakub Sumera has recently received his M. Sc. in chemistry.His field of interest includes theoretical and computational che-mistry.Efficient Properties Calculations at the MP2 Theory Level 673Krzysztof Kowalzyk has recently received his M. Sc. in che-mistry. His field of interest includes theoretical and computa-tional chemistry.

Efficient Energy and Electrostatic Properties Calculations at the MP2 Theory Level: A Case Study of Density Matrix-Based Computational Quantum Chemistry

http://www.cai.sk/ojs/index.php/cai/article/download/1013/453

Efficient Energy and Electrostatic Properties Calculations at the MP2 Theory Level: A Case Study of Density Matrix-Based Computational Quantum Chemistry

Abstract

Similar works

Full text

Available Versions

Computing and Informatics (E-Journal - Institute of Informatics, SAS, Bratislava)