Hybrid Extended Lagrangian, Post-Hartree Fock Born Oppenheimer Ab . - IU
Articlepubs.acs.org/JCTCHybrid Extended Lagrangian, Post-Hartree Fock Born Oppenheimer ab Initio Molecular Dynamics Using Fragment-BasedElectronic StructureJunjie Li, Cody Haycraft, and Srinivasan S. Iyengar*Department of Chemistry and Department of Physics, Indiana University, 800 E. Kirkwood Ave., Bloomington, Indiana 47405,United StatesS Supporting Information*ABSTRACT: We present a hybrid ab initio molecular dynamics scheme that includes bothDFT and Hartree Fock-based extended Lagrangian and converged post-Hartree FockBorn Oppenheimer components, combined within the framework of a molecularfragmentation-based electronic structure. The specific fragmentation algorithm used here isderived from ONIOM but includes multiple, overlapping “model” systems. The interactionbetween the various overlapping model systems is approximated by invoking the principle ofinclusion exclusion at the chosen higher level of theory and within a “real” calculationperformed at the chosen lower level of theory. Furthermore, here, the lower level electronicstructure of the full system is propagated through an extended Lagrangian formalism, whereasthe fragments, treated using post-Hartree Fock electronic structure theories, are computedusing the normal converged Born Oppenheimer treatment. This conservative dynamicalapproach largely reduces the computational cost to approximate on-the-fly dynamics usingpost-Hartree Fock electronic structure techniques and can be very beneficial for large systemswhere SCF convergence may be challenging and time consuming. Benchmarks are providedfor medium-sized protonated water clusters, H9O 4 and H13O 6 , and polypeptide fragments, including a proline tripeptidefragment, and alanine decamer. Structural features are in excellent agreement between the hybrid approach using an MP2:B3LYPfragment-based electronic structure and BOMD using MP2 for the full system. Vibrational properties derived from dynamicalcorrelation functions do show a small redshift for the extended Lagrangian treatments, especially at higher frequencies. Strategiesare discussed to improve this redshift. The computational methodology works in parallel using both MPI and OpenMP andshows good scaling with the processor number. The timing benchmarks are provided for the alanine decamer. A powerful featureof the computational implementation is the fact that it is completely decoupled from the electronic structure package beingemployed and thus allows for an integrated approach that may include several different packages. These computational aspectswill be further probed in future publications.I. INTRODUCTIONSince the seminal work of Karplus1 and Leforestier2 and theextended Lagrangian generalization by Car and Parrinello,3 abinitio molecular dynamics (AIMD) has become a central toolfor studying reactive events and vibrational properties beyondthe harmonic approximation and for biochemical modelingwhen combined with hybrid techniques. Here, electronicstructure and gradients are simultaneously computed in-stepwith nuclear motion. Due to this, AIMD is in general muchmore resource intensive as compared to both classical, forcefield-based, molecular dynamics methods as well as commonlyused electronic structure applications such as geometryoptimization and frequency calculations. This limitation hasthus far restricted AIMD mostly to DFT-based application formedium-sized systems. Improvements in electronic structuremethodology, certainly, have critical impact on AIMD, andone of the ways to reduce the cost of the on-the-fly electronicstructure calculations is through the utilization of the recentlypopularized fragment-based electronic structure schemes,4 25where a molecular system is divided into multiple regions and 2016 American Chemical Societythe total energy of the system is assembled from the energy ofindividual pieces. Detailed discussion of these methods can befound in several recently published review articles.26 30To explore the application of fragment-based electronicstructure in AIMD and potentially ab initio quantum nucleardynamics,31 we recently32 (a) considered a generalization tothe ONIOM5,13,33 41 energy extrapolation scheme using theset-theoretic principle of inclusion and exclusion42 to treatmultiple overlapping model systems32 (abbreviated in thispublication as the PIE-ONIOM energy expression) and (b)applied this method to obtain accurate potential energysurfaces and to compute ab initio dynamical trajectories withcorrelation functions that agree well with higher levels oftheory such as MP2 and coupled cluster. Specifically, meanabsolute errors in potential energy surfaces are (a) less than0.05 kcal/mol when the PIE-ONIOM results that combineB3LYP and MP2 are compared with MP2 calculations on theReceived: January 1, 2016Published: May 10, 20162493DOI: 10.1021/acs.jctc.6b00001J. Chem. Theory Comput. 2016, 12, 2493 2508
ArticleJournal of Chemical Theory and ComputationFigure 1. Figure 1 and Figure 2 briefly summarize the developments from ref 32. Here, the vibrational ground state (ρ(R)) weighted meanabsolute error in the potential surface is provided for H9O 4 . The error is computed between the fragment-based calculations and full CCSD(T)calculations as dR{ρ(R) ECCSD(T) ECCSD(T):B3LYP }. The statistics are obtained from 693 different geometries chosen along normal modecoordinates where the horizontal axis in the figure represents mode frequencies. Thus, the figure depicts the error in the critical regions of thesurface sampled by the ground state wavepacket and is clearly very small. When the full surface is included (that is, without vibrational groundstate weighting), the error rises to the order of 0.1 kcal/mol. All calculations use the 6-31 G(d,p) basis function.Figure 2. Solvated Zundel cation, H13O 6 : The vibrational density of states (VDOS) are computed from the fragment-based AIMD trajectories(frag-BOMD MP2:B3LYP), and these agree quite well with full MP2-based AIMD. All calculations use the 6-31 G(d,p) Gaussian basis set.and low levels. Thus, the paper is organized as follows:Section II introduces the new formalism with facilitatingtechnical discussions provided in Appendix A. The conservative properties of the approach are then tested in SectionIII. The approach is used to compute structural properties ofprotonated water clusters in Section IV, dynamical trajectories,and correlation functions for protonated water clusters andpolypeptide fragments in Section V. In all cases, thetrajectories are compared with MP2-based AIMD trajectories,and the results are in very good agreement. We alsodemonstrate numerical benchmarks in Section VI that showa large reduction in computational effort. Furthermore, thealgorithm is written in parallel using MPI (also demonstratedin Section VI) and through a C package that invokesexternal electronic structure software. Concluding remarks aremade in Section VII.full system and (b) less than 0.07 kcal/mol when a similarhybrid calculation is compared with CCSD(T) on the fullsystem. This data was obtained in ref 32 by considering 693different geometries of protonated water clusters sampledalong normal mode coordinates, and the associated computational effort is a very small fraction of the computational costincurred in the full calculations. The AIMD trajectoriescomputed using this approach reproduce vibrational density ofstates from regular MP2 calculations as well, and these resultsare summarized in Figures 1 and 2. As seen in Figure 1(b),and noted in the caption, the potential energy surface errorsare indeed extremely small and vibrational density of statespeaks agree quite well when the hybrid approach is comparedto higher level theory as seen from Figure 2(b).Despite the preliminary success of our approach in ref 32,several challenges remain: (1) During dynamical simulations,the molecular structural variation can lead to changes in thefragmentation topologies that must be handled carefully tomaintain conservative dynamics. (2) The computationalexpense may be dominated by the lower level calculationon the full system since our approach is based on ONIOM.This second hurdle is also seen in other hybrid fragmentationschemes, such as the Generalized Energy Based Fragmentation (GEBF),11 Molecular Tailoring Approach (MTA),14 theeXtended ONIOM method,15 and Molecule-in-Molecule(MIM).16 (3) Analytical gradient for the embedded chargescan be expensive.43,44 In the current publication, we propose asolution to the second problem listed above by constructing ahybrid method that combines Extended Lagrangian treatmentof the full system at the lower level of theory with BornOppenheimer dynamics for the fragments computed at highII. HYBRID EXTENDED LAGRANGIANBORN OPPENHEIMER DYNAMICS USINGFRAGMENT-BASED ELECTRONIC STRUCTUREFROM PIE-ONIOM FOR POST-HARTREE FOCKTREATMENTClassical ground state ab initio molecular dynamics primarilyhas two flavors: one is the traditional Born Oppenheimermolecular dynamics (BOMD), where the electronic energy isconverged, and the other is the extended Lagrangianmolecular dynamics (ELMD)3,45 57 formalism, where theelectronic degrees of freedom, either single particle densitymatrix or molecular orbital coefficients, are assigned fictitiousinertia and propagated simultaneously with the classical nucleithrough a simple adjustment of time scales.58 The first ELMD2494DOI: 10.1021/acs.jctc.6b00001J. Chem. Theory Comput. 2016, 12, 2493 2508
ArticleJournal of Chemical Theory and Computationmethod, the Car Parrinello molecular dynamics(CPMD)3,45,47 approach was influenced by earlier workfrom Anderson59 and Parrinello and Rahman45,60 and wasdeveloped in 1985. In CPMD, the wavefunction is expandedin a plane-wave basis, and pseudopotentials are used todescribe the core electrons. In 2001, Atom-centered DensityMatrix Propagation (ADMP),48 54 an ELMD formalism thatemploys the single particle electronic density matrix written inan all-electron atom-centered Gaussian basis was developed.The ADMP method has been demonstrated throughapplications to vibrational spectroscopy61 66 and biologicalion-channel studies54 and the study of atmospheric reactions,where pump probe techniques have been incorporated intothe dynamics.67 69 Critically, the ADMP formalism has hadan important role in the prediction of the amphiphilic natureof the hydrated proton70,71 that was later confirmed throughseveral experiments.70 76 The ADMP method has manyattractive features. Systems can be simulated by accuratelytreating all electrons or by using pseudopotentials. Throughthe use of smaller values for the tensorial fictitious mass,relatively large time steps can be employed, and lighter atomssuch as hydrogen nuclei are routinely used. Atom-centeredfunctions can be used with the appropriate physical boundaryconditions for molecules, polymers, surfaces, and solids,without the need to treat replicated images to impose 3dperiodicity. Consequently, charged systems and QM/MMmodels of biological systems54 can be treated as easily asneutral molecules. Deviations from the Born Oppenheimersurface and the mixing of fictitious and real kinetic energiescan also be rigorously controlled on-the-fly in ADMP. ADMPtrajectories of the order of picoseconds show stable dynamics,and the adiabaticity can be controlled effectively in thesesystems without thermostats. The important conceptual andcomputational differences between ADMP and other Gaussianbasis set-based implementations77 80 of the Car Parrinelloformalism have been discussed in detail in ref 50.The elimination of SCF in extended Lagrangian MD is acomputational advantage, especially for large molecules andsystems where SCF convergence may be challenging.However, thus far, a critical restriction for ELMD (Car Parrinello as well as ADMP) has been that it can only beapplied to single particle formalisms such as DFT andHartree Fock. Here, we combine ELMD and post-Hartree Fock-based BOMD through PIE-ONIOM and allow extendedLagrangian (Car Parrinello and ADMP) techniques tobecome applicable to higher level electronic structuremethods, such as MP2 and potentially CCSD, and alsoimprove the computational efficiency of fragment-basedAIMD. However, before we embark into a discussion ofthis dynamical formalism, we must first summarize the PIEONIOM energy extrapolation equations. The overall PIEONIOM molecular energy expression32 is given bywhen two levels of theory are used. (This approach is easilygeneralized to more levels,32 but this generalization is notconsidered here.) In eq 1, :(i j) represents theextrapolation term for the derivative fragment (i j), formedfrom the intersection of the ith and jth primitive fragments.Every fragment is treated at two levels of theory as written ineq 2, and the entire system (0th fragment) is only consideredat the lowest level of theory (with energy Elow(0)). Thisenergy extrapolation scheme is pictorially demonstrated inFigure 3. A powerful feature of the algorithm is the fact thatFigure 3. Illustration for the energy expression in eq 1. Theintersections between primitive fragments are subtracted out toremove overcounting according to the principle of inclusion andexclusion.the identity of the overlapping model systems or fragmentsare efficiently evaluated either (i) through an automated bitmanipulation algorithm presented in ref 32 that recomputesfragments at each instant in time during dynamics or (ii)through user specification. In the bit-manipulation algorithm,fragments are represented as bits in an integer, andoverlapping regions are computed using binary ANDoperations. The use of bitwise arithmetic improves efficiency,but a large integer library may be required, depending on thenumber of fragments involved and the platform being used forthe calculation. One such library is available at http://www.boost.org/.Clearly, for a large enough system the computationalexpenses may be dominated by the lower level calculation,Elow(0). Thus, in the scheme introduced here, it is electronicdegrees of freedom of the entire system, Elow(0) that we willpropagate using an extended Lagrangian technique. The restof the fragments, where more accurate calculations maybeneeded, will be treated using converged Born Oppenheimertreatments with post-HF electronic structure. While wechoose to demonstrate this approach using the ADMPLagrangian48,49 as the choice for ELMD, extensions to employCPMD (that is, using molecular orbital coefficients instead ofthe single particle density matrices with plane wave basis setsinstead of the atom-centered Gaussian basis functionsemployed here) can be done in a similar fashion. TheLagrangian of the entire system is111/41/4 2Wlow, 0 μ low,0)]3 Tr[VTMV] Tr[(μ low,022 EPIE‐ONIOM(R, Plow, 0) Tr[Λlow, 0(P2low, 0 Plow, 0)]nEPIE‐ONIOM E low (0) :(i) i 1 :(i j)1 i j n:(i j k) ···1 i j k n ( 1)n 1 :(1 ··· n)(1)(3)where the extrapolation term is defined as:(i) E high(i) E low (i)In this equation, R, V, and M are the nuclear geometry,velocity, and mass, respectively, and the single particle density(2)2495DOI: 10.1021/acs.jctc.6b00001J. Chem. Theory Comput. 2016, 12, 2493 2508
ArticleJournal of Chemical Theory and Computationmatrix Plow,0 that determines Elow(0) has fictitious inertiatensor μlow,0, and velocity Wlow,0. The quantity Λlow,0 is aLagrangian multiplier that maintains the N-representability ofPlow,0. The energy functional EPIE‑ONIOM(R,Plow,0) is identicalto that in eq 1 but differs with the additional dependence ofElow(0) Elow(0)(R,Plow,0) since the lower level full systemelectronic structure is now to be propagated rather thanconverged. This leads to the equations of motion for nucleiand density matrixM E 2R 2 tPIE‐ONIOMcalculations, which appears to be a problem for Born Oppenheimer-type dynamics calculations in ref 82.Proceeding further with the discussion of the equations ofmotion, since Plow,0 is propagated according to eq 5, thegradient term(4)Plow, 0and1/2μ low,0 2Plow, 0 t 2 1/2μ low,0 EPIE‐ONIOM(R, Plow, 0) Plow, 0 Λlow, 0Plow, 0R Plow, 0Λlow, 0 Λlow, 0(5)where EPIE‐ONIOM(R, Plow, 0) Rn i 1 :(i) R Plow, 0 1 i j n E low (0) RPlow, 0 :(i j) ··· R(6) E low (0) Plow, 0(7)and EPIE‐ONIOM(R, Plow, 0) Plow, 0 RRin eq 6 contains forces from thePlow, 0unconverged density matrix, Plow,0, in addition to the standardBorn Oppenheimer gradients, and these aspects have beenwell described in refs 48, 49, and 51. We summarize thesegradients in Appendix A since the major portion of theseforces have already been discussed in other publications notedabove. It is only the interpretation of these gradients withrespect to the fragment-based technology that is new here.The equations of motion are integrated using the velocityVerlet algorithm,83 while the Lagrangian multiplier matricesare determined through an iterative procedure.48,49 In thesections below, we probe the accuracy and efficiency of thishybrid scheme by studying protonated water clusters andpolypeptides fragments. In all cases, consistency of thedynamics is checked by computing the conserved quantities,and these conserved quantities (the conjugate Hamiltonian forthe Lagrangian in eq 3) are outlined in Appendix B.However, before we proceed further, we must note anddistinguish our approach from two other recent developmentsin Born Oppenheimer and extended Lagrangian dynamics. Inref 84, the author has carefully analyzed the time scale ofoscillations of the different constituents of energy in a postHartree Fock (MP2) gradient expression. The author infersthrough many computational studies that the post-Hartree Fock portion of the gradient is in fact a slowly varyingquantity as a function of nuclear motion, which allows theauthors to construct a multiple time step approach thatappears to be quite powerful. This extension is certainly anoption available to us here and will be considered in futurepublications. In ref 82, the authors address the critical driftproblem in standard Born Oppenheimer dynamics bypropagating a surrogate electronic density through anextended Lagrangian formalism, where such a density isharmonically bound to a converged density. This latterapproach has extended Lagrangian and Born Oppenheimercomponents, much like the current formalism, but differs inthat the ideas in ref 82 are thus far implemented for DFT andDFTB, whereas the approach here is specifically gearedtoward computing post-Hartree Fock dynamics. Thus, todifferentiate between the two approaches, from hereon, themethod developed in this paper has both extended Lagrangianand post-Hartree Fock Born Oppenheimer dynamics components funneling from the ADMP Lagrangian in Eq. (3), andis hence referred to as ADMP-pHF, and the PIE-ONIOMenergy functional, eq 1, is implied. When the extendedLagrangian treatment is not invoked and Born Oppenheimertreatment is used for all components of the energy expressionin eq 1, the associated method is referred to here as fragBOMD.(R, Plow, 0) R E low (0) RIn eq 6, the constraint, “ Plow,0”, is only included on the firstterm on the right side since that is the only term that is a :(i j)function of Plow,0. The terms :(i) ,, etc., involve R Rstandard electronic structure gradients that are already presentin most electronic structure packages. This is because, asnoted in eq 2, these terms involve gradients from higher levelelectronic structure methods (MP2 for the calculationspresented in here), and when DFT is used for the lowerlevel in eq 2, the SCF is converged to a high degree ofaccuracy. In this publication, all SCF calculations involved in :(i) :(i j),, are well converged requiring the RMS change in R Rdensity matrix to be less than 10 8 and the maximum changein the density matrix to be less that 10 6. These are thedefault threshold conditions used in the version of theGaussian quantum chemistry package81 used here to computethe electronic energy and nuclear gradients for the fragments.The energy change during SCF is not used to testconvergence, but a 10 N RMS change in density matrixtypically corresponds to a 10 2N change in energy in atomicunits. While in principle this translates to an energyconvergence threshold of 10 16 Hartree, in practice we findthat SCF convergence is always better than 10 12 Hartree inmost calculations presented here, with very few at 10 10. (Amore detailed analysis on this aspect is presented in SectionIII.) As a result of this, there is practically no drift in theseIII. ACCURACY AND EFFICIENCY OF THE HYBRIDEXTENDED LAGRANGIAN, POST-HARTREE FOCKBORN OPPENHEIMER DYNAMICS (ADMP-PHF)TREATMENTWe benchmark the dynamics formalism by studying (smalland medium-sized) protonated water clusters and twodifferent polypeptide fragments. We choose two challenging2496DOI: 10.1021/acs.jctc.6b00001J. Chem. Theory Comput. 2016, 12, 2493 2508
ArticleJournal of Chemical Theory and ComputationFigure 4. Protonated water clusters shown in Figures 1(a) and 2(b) are fragmented as depicted here, and the electronic energy is assembledusing eq 1. The proline tripeptide (c) is similarly fragmented using overlapping dipeptides (d). Dangling bonds are saturated with link atoms(hydrogen) as in standard ONIOM.protonated water clusters, the Eigen cation (H9O 4 ) and thesolvated Zundel cation (H13O 6 ), as our benchmark systemswith their structures shown in Figure 4(a) and (b). Thesetype of systems have significance in biological,85 92 atmospheric,69,93 95 and condensed phase chemistry and have beenstudied extensively by both experiments96 104 and bytheories.61,62,64,70,71,97,105 120 Here, accurate treatment ofhydrogen bonding is crucial to model structural anddynamical behavior. In addition, polypeptide fragments haverecently been studied in some detail using mass spectrometrywith a goal to elucidate conformational transformations andkinetics.121,122 For example, in refs 121 and 123, the authorsobserve a step-by-step conformational transition from an allcis polyproline helix chain to an all-trans polyproline helixchain. While our primary benchmarks here that test theaccuracy of ADMP-pHF and frag-BOMD through comparisonwith full system MP2 calculations are performed for theprotonated water clusters stated above, we also presentpreliminary benchmarks on an all-trans proline tripeptidefragment (Figure 4(c)) to show the consistency of ourapproach and its potential for future applications to the studyof biomolecules. In addition, we also present an analysis ofthe computational gain from the fragment-based ADMP-pHFand frag-BOMD formalisms by computing trajectories for anall-trans polyalanine fragment.At this point, it is critical to note that for the peptidesystems, chemical bonds intersect the boundary between twofragments, and hence, link atoms are added to saturate thedangling valencies as in standard ONIOM.31,54,124 Thepositions of link atoms are uniquely determined based onthe connectivity of the system. Specifically, the placement oflink atoms is constrained by the positions of the two atomson either ends of the bond intersecting the fragmentboundary according to rlink rbond g(rsub rbond), whererbond is the position of the bonded atom (B) inside thefragment in question to which the link is bound to substitutethe atom (S) belonging to the surrounding region and locatedat rsub. The quantity g is a scale factor normally set to R(B L)/R(B S), where R(B L) is a typical bond length betweenthe bonded atom within the fragment and the link atom(hydrogen in all cases in this work), and R(B S) is a typicalbond length for the pair of atoms describing the bond beingdissected to create the fragment. Careful tests have shownthat the ONIOM results are quite insensitive to the value ofg.37 It must be stressed that since the link atom positions arefully defined by the coordinates of the full system, as instandard ONIOM, the gradients with respect to link atomcoordinates can be transferred to coordinates of the fullsystem for this fragment-based methodology.54All calculations are performed using a C package, that is,MPI parallelized across nodes and OpenMP parallelizedwithin each node. For the current publication, the packageinvokes a development version of the Gaussian quantumchemistry program81 to compute the electronic energy andnuclear gradients for the fragments and the full system.However, our software package can be easily generalized toinclude other electronic structure software or combinations ofthese. These computational aspects and associated generalizations will be discussed in future publications.The dynamical simulations are performed five differentways for each system with all calculations using the 631 g(d,p) atom-centered Gaussian basis set. The five differentsets of trajectories are (1) BOMD using pure B3LYP(abbreviated as BOMD B3LYP), (2) ADMP using pureB3LYP (abbreviated as ADMP B3LYP), (3) BOMD usingpure MP2 (abbreviated as BOMD MP2), and (4) BOMDusing two-layer PIE-ONIOM fragment-based electronicstructure. MP2 is used for the high level electronic structure,while the low level is calculated at B3LYP. This calculation isabbreviated as frag-BOMD MP2:B3LYP and was fullybenchmarked in ref 32. Also see Figures 1 and 2 where wepresent the accuracy of this calculation. Finally, (5) the two2497DOI: 10.1021/acs.jctc.6b00001J. Chem. Theory Comput. 2016, 12, 2493 2508
ArticleJournal of Chemical Theory and ComputationTable 1. Energy Conservation Properties for Dynamical Simulations (NVE Simulations)systemlevel of theoryatimebH9O 4BOMD B3LYPADMP B3LYPBOMD MP2frag-BOMD MP2:B3LYPADMP-pHF MP2:B3LYPBOMD B3LYPADMP B3LYPBOMD MP2frag-BOMD MP2:B3LYPADMP-pHF MP2:B3LYPADMP B3LYPBOMD B3LYPfrag-BOMD MP2:B3LYPADMP-pHF MP2:B3LYP20.0 ps20.0 ps20.0 ps5.2 ps8.6 ps20.0 ps13.4 ps5.2 ps6.0 ps5.0 ps4.33 ps2.18 ps2.97 ps2.90 psH13O 6proline trimerfave temp 2151.22156.44237.96238.09246.30246.07 .522.5Δ/ c (in kcal/mol)driftd (in kcal/mol)ΔESCFe (in 74 10 116.384.083.342.75 10 1210 1210 1210 111.32 10 116.05 10 122.93 10 128.41 10 122.96 10 121.64 10 11aWater clusters: 6-31 g(d,p) basis is used for B3LYP and MP2. Proline trimer: 6-31 g(d) basis is used for B3LYP and MP2. bTime step 0.2 fs for allsimulations. cRoot mean square deviation of the total energy in kcal/mol. The total energy corresponds to the Hamiltonian in eq B4. dDrift in totalenergy is computed as a difference between the average total energies for the first and last 100 fs of dynamics (in kcal/mol). In addition, in Figure SI1 of the Supporting Information, we provide the evolution of total energy. This is provided in eV to facilitate comparison with results from otherwork such as ref 82. eWhile SCF convergence thresholds are placed on the density matrix (see text for details) rather than energy, in this column wenote the average absolute change in energy (in Hartrees) during the last interation before convergence is achieved for all steps in the trajectory. fCterminus is capped with a carboxylic acid group (see Figure 4(c)).Figure 5. Comparison of forces around an optimized geometry (a) and during dynamics (b) for H9O 4 . See text for details.in density matrix be less than 10 8 and the maximum changein the density matrix be less that 10 6. The correspondingSCF energy convergence is very tight and is noted in Table 1.In addition to gauging conserved quantities, we alsocompare forces in Figure 5. Specifically, in Figure 5(a) wecompute the forces for 693 different geometries chosen alongnormal mode coordinates where the horizontal axis in thefigure represents mode frequencies. The forces are computedat the level of MP2 and also using the MP2:B3LYP fragment FMP2:B3LYP FMP2 2method. The differenceis presented as alayer PIE-ONIOM fragment-based ADMP-pHF schemeproposed in this publication. The methods MP2 and B3LYPare used for high and low level of theory, respectively, and thecalculation is abbreviated as ADMP-pHF MP2:B3LYP.The stability and consistency of the dynamics trajectoriescan be gauged from the total energy conservation data forthese classical trajectories. The total energy corresponds tothe conjugate Hamiltonian in eq B4 and contains the nuclearkinetic energy, the fictitious density matrix kinetic energy,EPIE‑ONIOM, and any energy penalty that is added as a result ofresidual loss of N-representability. As shown in Table 1, alltrajectories have the total energy well conserved to within 2 3 orders of magnitude less than 1 kcal/mol. The ADMP-pHFscheme produces energy conservation comparable to theBOMD simulations. We also inspect the drift in the conjugateHamiltonian in eq B4, and as noted in Table 1, this quantityis also extremely small and well within the acceptable range.82In addition, in Table 1, we also present the change in SCFenergy in the last iteration before convergence is achieved.This quantity in some sense represents a convergencethreshold of the electronic energy. As noted previously,convergence is enforced on the density matrix and not on theSCF energy. In this case, it is required that the RMS change3NAtomsfunction of harmonic frequency. In addition, for the analysisin Figure 5(b), we have retrieved geometries at fixed intervalsfrom the extended Lagrangian Born Oppenheimer trajectoryfor H9O 4 and computed forces using MP2. The vertical axisin Figure 5(b) is the mean absolute deviation between theseforces. As shown, the forces in both figures are well withinthat expected for chemical accuracy.IV. STRUCTURAL FEATURES FROM DYNAMICSSIMULATIONSNext, we examine the structural properties of thesetrajectories by calculating radial distribution functions for2498DOI: 10.1021/acs.jctc.6b00001J. Chem. Theory Comput. 2016, 12, 2493 2508
ArticleJournal of Chemical Theory and ComputationFigure 6. Radial distribution function for O O an
Classical ground state ab initio molecular dynamics primarily has two flavors: one is the traditional Born Oppenheimer molecular dynamics (BOMD), where the electronic energy is converged, and the other is the extended Lagrangian molecular dynamics (ELMD)3,45 57 formalism, where the electronic degrees of freedom, either single particle density
hartree partners purchase order terms and conditions: 12-02-20 1 hartree partners purchase order terms and conditions . the terms and conditions contained herein shall apply only in the event that contractor and hartree partners, lp or any hartree affiliate (collectively referred to herein as “hartree”) are not parties to any other agreement in
2.1 Lagrangian Mechanics 2.1.1 The Lagrangian The core of Lagrangian Mechanics is the Lagrangian, a function of positions xAand velocities x_A of all the particles, which summarizes the dynamics of a system. Any function which generates the correct equations of motion can be taken
Hartree Partners Hartree Partners, LP is a global merchant commodities firm specialising in energy and its associated industries. Founded originally as Hess Energy Trading Company LLC (HETCO), Hartree’s business spans 10 offices around the world. Objective Hartree Partners occupied circa 8,000 sq ft offices in a building
Finite temperature Hartree-Fock Mermin (1963) examined finite temperature Hartree-Fock, showing in the grand ensemble that there is an ensemble of Slater determinants that minimizes the Hartree-Fock free energy and provides an upper bound as at zero temperature. with respective orbitals and Single particle equations FHF T THF SHF U .
The Hartree and Hartree-Fock eigenproblems provide quantum mechanical models of atoms and molécules which are more tractable than the full Schrodmger équations for these Systems. They have been extensively used for computational modelmg since their introduction by Hartree, Fock and Slater [10], [7] and [21] in the early days of quantum theory.
Hartree and Hartree-Fock Equations Starting with the many-body Schr odinger equation i@ t H n; for a system of n identicalbosonsorfermionsand restricting it to theHartreeandHartree-Fock states n 1 and n 1 i; we obtain theHartreeand theHartree-Fockequations. There is a considerable literature on I the derivation of the Hartree and Hartree .
2. Extended Lagrangian: Nose-Hoover thermostat [1, 7]. Equations of motion derived from the Lagrangian of the system conserve the total energy of the system. One can write an extended Lagrangian, by adding fictitious degrees of freedom, such that the overall total energy is conserved but the atomic subsystem can span ensembles .
Paediatric Anatomy Paediatric ENT Conditions Paediatric Hearing Tests and Screening. 1 Basic Sciences HEAD AND NECK ANATOMY 3 SECTION 1 ESSENTIAL REVISION NOTES medial pterygoid plate lateral pterygoid plate styloid process mastoid process foramen ovale foramen spinosum jugular foramen stylomastoid foramen foramen magnum carotid canal hypoglossal canal Fig. 1 The cranial fossa and nerves .
