Low-cost evaluation of the exchange Fock matrix from Cholesky and density fitting representations of the electron repulsion integrals

Psi4 1.1: An Open-Source Electronic Structure Program Emphasizing Automation, Advanced Libraries, and Interoperability

Psi4 is an ab initio electronic structure program providing methods such as Hartree-Fock, density functional theory, configuration interaction, and coupled-cluster theory. The 1.1 release represents a major update meant to automate complex tasks, such as geometry optimization using complete-basis-set extrapolation or focal-point methods. Conversion of the top-level code to a Python module means that Psi4 can now be used in complex workflows alongside other Python tools. Several new features have been added with the aid of libraries providing easy access to techniques such as density fitting, Cholesky decomposition, and Laplace denominators. The build system has been completely rewritten to simplify interoperability with independent, reusable software components for quantum chemistry. Finally, a wide range of new theoretical methods and analyses have been added to the code base, including functional-group and open-shell symmetry adapted perturbation theory, density-fitted coupled cluster with frozen natural orbitals, orbital-optimized perturbation and coupled-cluster methods (e.g., OO-MP2 and OO-LCCD), density-fitted multiconfigurational self-consistent field, density cumulant functional theory, algebraic-diagrammatic construction excited states, improvements to the geometry optimizer, and the "X2C" approach to relativistic corrections, among many other improvements.

Analytic energy gradients for the coupled-cluster singles and doubles method with the density-fitting approximation

An efficient implementation is presented for analytic gradients of the coupled-cluster singles and doubles (CCSD) method with the density-fitting approximation, denoted DF-CCSD. Frozen core terms are also included. When applied to a set of alkanes, the DF-CCSD analytic gradients are significantly accelerated compared to conventional CCSD for larger molecules. The efficiency of our DF-CCSD algorithm arises from the acceleration of several different terms, which are designated as the "gradient terms": computation of particle density matrices (PDMs), generalized Fock-matrix (GFM), solution of the Z-vector equation, formation of the relaxed PDMs and GFM, back-transformation of PDMs and GFM to the atomic orbital (AO) basis, and evaluation of gradients in the AO basis. For the largest member of the alkane set (C10H22), the computational times for the gradient terms (with the cc-pVTZ basis set) are 2582.6 (CCSD) and 310.7 (DF-CCSD) min, respectively, a speed up of more than 8-folds. For gradient related terms, the DF approach avoids the usage of four-index electron repulsion integrals. Based on our previous study [U. Bozkaya, J. Chem. Phys. 141, 124108 (2014)], our formalism completely avoids construction or storage of the 4-index two-particle density matrix (TPDM), using instead 2- and 3-index TPDMs. The DF approach introduces negligible errors for equilibrium bond lengths and harmonic vibrational frequencies.

An integral-factorized implementation of the driven similarity renormalization group second-order multireference perturbation theory

We report an efficient implementation of a second-order multireference perturbation theory based on the driven similarity renormalization group (DSRG-MRPT2) [C. Li and F. A. Evangelista, J. Chem. Theory Comput. 11, 2097 (2015)]. Our implementation employs factorized two-electron integrals to avoid storage of large four-index intermediates. It also exploits the block structure of the reference density matrices to reduce the computational cost to that of second-order Møller-Plesset perturbation theory. Our new DSRG-MRPT2 implementation is benchmarked on ten naphthyne isomers using basis sets up to quintuple-ζ quality. We find that the singlet-triplet splittings (ΔST) of the naphthyne isomers strongly depend on the equilibrium structures. For a consistent set of geometries, the ΔST values predicted by the DSRG-MRPT2 are in good agreements with those computed by the reduced multireference coupled cluster theory with singles, doubles, and perturbative triples.

Multireference Perturbation Theory with Cholesky Decomposition for the Density Matrix Renormalization Group

We present a second-order N-electron valence state perturbation theory (NEVPT2) based on a density matrix renormalization group (DMRG) reference wave function that exploits a Cholesky decomposition of the two-electron repulsion integrals (CD-DMRG-NEVPT2). With a parameter-free multireference perturbation theory approach at hand, the latter allows us to efficiently describe static and dynamic correlation in large molecular systems. We demonstrate the applicability of CD-DMRG-NEVPT2 for spin-state energetics of spin-crossover complexes involving calculations with more than 1000 atomic basis functions. We first assess, in a study of a heme model, the accuracy of the strongly and partially contracted variant of CD-DMRG-NEVPT2 before embarking on resolving a controversy about the spin ground state of a cobalt tropocoronand complex.

Truncated Transition Densities for Analysis of (Nonlinear) Optical Properties of carbo-Chromophores

The optical properties of several quadrupolar carbo-benzene derivatives are investigated at various levels of calculation (TDDFT and CASPT2) and analyzed using a new theoretical tool here disclosed: The "visualization" of the transition dipole moment from the transition density truncated to the main monoexcitations involved in the electronic transition (TTD). The experimental or calculated one-photon UV-visible absorption spectra of the carbo-benzene derivatives fit with the Gouterman model originally proposed for porphyrins, where the first four excited states involve linear combinations of monoexcitations of the same four frontier molecular orbitals. The relative intensities of the absorption bands are analyzed from the transition dipole moments calculated from the TTDs and an analogy between porphyrins and carbo-benzenes is argued. The two-photon absorption (TPA) cross section related to the third-order nonlinear optical response is calculated for each two-photon-allowed excited state |f⟩ from the contribution of all possible intermediate excited states |i⟩ using the "sum-over-state" (SOS) scheme. The quadrupolar carbo-benzene derivatives fit into the three-level model, as their TPA cross section exhibits a dominant contribution of one of the intermediate excited states. The origin of TPA efficiency (enhancement) upon carbo-merisation of the C-C link to the para-substituents is discussed from the excitation energies of the intermediate and final excited states and from the two corresponding transition dipole moments (μ0i and μif). The latter may be calculated from the TTDs.

Molecular and Electronic Structure of Re2Br4(PMe3)4

The dinuclear rhenium(II) complex Re2Br4(PMe3)4 was prepared from the reduction of [Re2Br8](2-) with (n-Bu4N)BH4 in the presence of PMe3 in propanol. The complex was characterized by single-crystal X-ray diffraction (SCXRD) and UV-visible spectroscopy. It crystallizes in the monoclinic C2/c space group and is isostructural with its molybdenum and technetium analogues. The Re-Re distance (2.2521(3) Å) is slightly longer than the one in Re2Cl4(PMe3)4 (2.247(1) Å). The molecular and electronic structure of Re2X4(PMe3)4 (X = Cl, Br) were studied by multiconfigurational quantum chemical methods. The computed ground-state geometry is in excellent agreement with the experimental structure determined by SCXRD. The calculated total bond order (2.75) is consistent with the presence of an electron-rich triple bond and is similar to the one found for Re2Cl4(PMe3)4. The electronic absorption spectrum of Re2Br4(PMe3)4 was recorded in benzene and shows a series of low-intensity bands in the range 10 000-26 000 cm(-1). The absorption bands were assigned based on calculations of the excitation energies with the multireference wave functions followed by second-order perturbation theory using the CASSCF/CASPT2 method. Calculations predict that the lowest energy band corresponds to the δ* → σ* transition, while the next higher energy bands were attributed to the δ* → π*, δ → σ*, and δ → π* transitions.

Octafluorodirhenate(III) Revisited: Solid-State Preparation, Characterization, and Multiconfigurational Quantum Chemical Calculations

A simple method for the high-yield preparation of (NH4)2[Re2F8]·2H2O has been developed that involves the reaction of (n-Bu4N)2[Re2Cl8] with molten ammonium bifluoride (NH4HF2). Using this method, the new salt [NH4]2[Re2F8]·2H2O was prepared in ∼90% yield. The product was characterized in solution by ultraviolet-visible light (UV-vis) and (19)F nuclear magnetic resonance ((19)F NMR) spectroscopies and in the solid-state by elemental analysis, powder X-ray diffraction (XRD), and infrared (IR) spectroscopy. Multiconfigurational CASSCF/CASPT2 quantum chemical calculations were performed to investigate the molecular and electronic structure, as well as the electronic absorption spectrum of the [Re2F8](2-) anion. The metal-metal bonding in the Re2(6+) unit was quantified in terms of effective bond order (EBO) and compared to that of its [Re2Cl8](2-) and [Re2Br8](2-) analogues.

Multiconfigurational Second-Order Perturbation Theory with Frozen Natural Orbitals Extended to the Treatment of Photochemical Problems

A new flavor of the frozen natural orbital complete active space second-order perturbation theory method (FNO-CASPT2, Aquilante et al., J. Chem. Phys. 131, 034113) is proposed herein. In this new implementation, the virtual space in Cholesky decomposition-based CASPT2 computations (CD-CASPT2) is truncated by excluding those orbitals that contribute the least toward preserving a predefined value of the trace of an approximate density matrix, as that represents a measure of the amount of dynamic correlation retained in the model. In this way, the amount of correlation included is practically constant at all nuclear arrangements, thus allowing for the computation of smooth electronic states surfaces and energy gradients-essential requirements for theoretical studies in photochemistry. The method has been benchmarked for a series of relevant biochromophores for which large speed-ups have been recorded while retaining the accuracy achieved in the corresponding CD-CASPT2 calculations. Both vertical excitation energies and gradient calculations have been carried out to establish general guidelines as to how much correlation needs to be retained in the calculation for the results to be consistent with the CD-CASPT2 findings. Our results feature errors within a tenth of an eV for the most difficult cases and have been validated to be used for gradient computations where an up to 3-fold speed-up is observed depending on the size of the system and the basis set employed.

Detection and Electronic Structure of Naked Actinide Complexes: Rhombic-Ring (AnN)2 Molecules Stabilized by Delocalized π-Bonding

The major products of the reaction of laser ablated and excited U atoms and N2 are the linear N≡U≡N dinitride molecule, isoelectronic with the uranyl dication, and the diatomic nitride U≡N. These molecules form novel cyclic dimers, (UN)2 and (NUN)2, with complex electronic structures, in matrix isolation experiments, which increase on UV photolysis. In addition, (NUN)2 increases at the expense of (UN)2 upon warming the codeposited matrix samples into the 20-40 K range as attested by additional nitrogen and argon matrix infrared spectra recorded after cooling the samples back to 4 or 7 K. These molecules are identified through matrix infrared spectra with nitrogen isotopic substitution and by comparing the observed matrix frequencies with those from electronic structure calculations. The dimerization is strong (theory predicts the dimer to be on the order of 100 kcal/mol more stable than the monomers), since the ground state involves 12 bonding electrons, 8 in the σ-system, and 4 in the delocalized π-system. This delocalized π bonding is present in the U, Th, La, and Hf analogues further demonstrating the interesting interplay between the 5f and 6d orbitals in actinide chemistry. The (UN)2(+) cation is also observed in solid argon, and calculations indicate that the bonding in the ring is preserved. On the other hand, the NUN dimer is of lower C2h symmetry, and the initial NUN molecules are recognizable in this more weakly bonded (ΔE = -64 kcal/mol) structure. The NThN molecules bind more strongly in the (NThN)2 dimer than the NUN molecules in (NUN)2 since NUN itself is more stable than NThN.

TD-DFT study of the light-induced spin crossover of Fe(iii) complexes

Two light-induced spin-crossover Fe(iii) compounds have been studied with time-dependent density functional theory (TD-DFT) to investigate the deactivation mechanism and the role of the ligand-field states as intermediates in this process.

Orbital-optimized linearized coupled-cluster doubles with density-fitting and Cholesky decomposition approximations: an efficient implementation

An efficient implementation of the orbital-optimized linearized coupled-cluster double method with the density-fitting (DF-OLCCD) and Cholesky decomposition (CD-OLCCD) approximations is presented.

Simulations of iron K pre-edge X-ray absorption spectra using the restricted active space method

The intensities and relative energies of metal K pre-edge features are sensitive to both geometric and electronic structures.

Cholesky Decomposition-Based Multiconfiguration Second-Order Perturbation Theory (CD-CASPT2): Application to the Spin-State Energetics of Co(III)(diiminato)(NPh)

The electronic structure and low-lying electronic states of a Co(III)(diiminato)(NPh) complex have been studied using multiconfigurational wave function theory (CASSCF/CASPT2). The results have been compared to those obtained with density functional theory. The best agreement with ab initio results is obtained with a modified B3LYP functional containing a reduced amount (15%) of Hartree-Fock exchange. A relativistic basis set with 869 functions has been employed in the most extensive ab initio calculations, where a Cholesky decomposition technique was used to overcome problems arising from the large size of the two-electron integral matrix. It is shown that this approximation reproduces results obtained with the full integral set to a high accuracy, thus opening the possibility to use this approach to perform multiconfigurational wave-function-based quantum chemistry on much larger systems relative to what has been possible until now.

Molcas 8: New capabilities for multiconfigurational quantum chemical calculations across the periodic table

In this report, we summarize and describe the recent unique updates and additions to the Molcas quantum chemistry program suite as contained in release version 8. These updates include natural and spin orbitals for studies of magnetic properties, local and linear scaling methods for the Douglas-Kroll-Hess transformation, the generalized active space concept in MCSCF methods, a combination of multiconfigurational wave functions with density functional theory in the MC-PDFT method, additional methods for computation of magnetic properties, methods for diabatization, analytical gradients of state average complete active space SCF in association with density fitting, methods for constrained fragment optimization, large-scale parallel multireference configuration interaction including analytic gradients via the interface to the Columbus package, and approximations of the CASPT2 method to be used for computations of large systems. In addition, the report includes the description of a computational machinery for nonlinear optical spectroscopy through an interface to the QM/MM package Cobramm. Further, a module to run molecular dynamics simulations is added, two surface hopping algorithms are included to enable nonadiabatic calculations, and the DQ method for diabatization is added. Finally, we report on the subject of improvements with respects to alternative file options and parallelization.

Pushing the Limits of Delta Bonding in Metal-Chromium Complexes with Redox Changes and Metal Swapping

Into the metalloligand Cr[N(o-(NCH2P((i)Pr)2)C6H4)3] (1, CrL) was inserted a second chromium atom to generate the dichromium complex Cr2L (2), which is a homobimetallic analogue of the known MCrL complexes, where M is manganese (3) or iron (4). The cationic and anionic counterparts, [MCrL](+) and [MCrL](-), respectively, were targeted, and each MCr pair was isolated in at least one other redox state. The solid-state structures of the [MCrL](+,0,-) redox members are essentially the same, with ultrashort metal-metal bonds between 1.96 and 1.74 Å. The formal shortness ratios (r) of these interactions are between 0.84 and 0.74 and are interpreted as triple to quintuple metal-metal bonds with the aid of theory. The trio of (d-d)(10) species [Cr2L](-) (2(red)), MnCrL (3), and [FeCrL](+) (4(ox)) are S = 0 diamagnets. On the basis of M-Cr bond distances and theoretical calculations, the strength of the metal-metal bond across the (d-d)(10) series increases in the order Fe < Mn < Cr. The methylene protons in the ligand are shifted downfield in the (1)H NMR spectra, and the diamagnetic anisotropy of the metal-metal bond was calculated as -3500 × 10(-36), -3900 × 10(-36), and -5800 × 10(-36) m(3) molecule(-1) for 2(red), 3, and 4(ox) respectively. The magnitude of diamagnetic anisotropy is, thus, affected more by bond polarity than by bond order. A comparative vis-NIR study of quintuply bonded 2(red) and 3 revealed a large red shift in the δ(4) → δ(3)δ* transition energy upon swapping from the (Cr2)(2+) to the (MnCr)(3+) core. Complex 2(red) was further investigated by resonance Raman spectroscopy, and a band at 434 cm(-1) was assigned as the Cr-Cr bond vibration. Finally, 4(ox) exhibited a Mössbauer doublet with an isomer shift of 0.18 mm/s that suggests a primarily Fe-based oxidation to Fe(I).

Derivation of general analytic gradient expressions for density-fitted post-Hartree-Fock methods: an efficient implementation for the density-fitted second-order Møller-Plesset perturbation theory

General analytic gradient expressions (with the frozen-core approximation) are presented for density-fitted post-HF methods. An efficient implementation of frozen-core analytic gradients for the second-order Møller-Plesset perturbation theory (MP2) with the density-fitting (DF) approximation (applying to both reference and correlation energies), which is denoted as DF-MP2, is reported. The DF-MP2 method is applied to a set of alkanes, conjugated dienes, and noncovalent interaction complexes to compare the computational cost of single point analytic gradients with MP2 with the resolution of the identity approach (RI-MP2) [F. Weigend and M. Häser, Theor. Chem. Acc. 97, 331 (1997); R. A. Distasio, R. P. Steele, Y. M. Rhee, Y. Shao, and M. Head-Gordon, J. Comput. Chem. 28, 839 (2007)]. In the RI-MP2 method, the DF approach is used only for the correlation energy. Our results demonstrate that the DF-MP2 method substantially accelerate the RI-MP2 method for analytic gradient computations due to the reduced input/output (I/O) time. Because in the DF-MP2 method the DF approach is used for both reference and correlation energies, the storage of 4-index electron repulsion integrals (ERIs) are avoided, 3-index ERI tensors are employed instead. Further, as in case of integrals, our gradient equation is completely avoid construction or storage of the 4-index two-particle density matrix (TPDM), instead we use 2- and 3-index TPDMs. Hence, the I/O bottleneck of a gradient computation is significantly overcome. Therefore, the cost of the generalized-Fock matrix (GFM), TPDM, solution of Z-vector equations, the back transformation of TPDM, and integral derivatives are substantially reduced when the DF approach is used for the entire energy expression. Further application results show that the DF approach introduce negligible errors for closed-shell reaction energies and equilibrium bond lengths.

Electronic Spectroscopy of [FePAH](+) Complexes in the Region of the Diffuse Interstellar Bands: Multireference Wave Function Studies on [FeC6H6](+)

The low-energy states and electronic spectrum in the near-infrared-visible region of [FeC6H6](+) are studied by theoretical approaches. An exhaustive exploration of the potential energy surface of [FeC6H6](+) is performed using the density functional theory method. The ground state is found to be a (4)A1 state. The structures of the lowest energy states ((4)A2 and (4)A1) are used to perform multireference wave function calculations by means of the multistate complete active space with perturbation at the second order method. Contrary to the density functional theory results ((4)A1 ground state), multireference perturbative calculations show that the (4)A2 state is the ground state. The vertical electronic spectrum is computed and compared with the astronomical diffuse interstellar bands, a set of near-infrared-visible bands detected on the extinction curve in our and other galaxies. Many transitions are found in this domain, corresponding to d → d, d → 4s, or d → π* excitations, but few are allowed and, if they are, their oscillation strengths are small. Even though some band positions could match some of the observed bands, the relative intensities do not fit, making the contribution of the [Fe-C6H6](+) complexes to the diffuse interstellar bands questionable. This work, however, lays the foundation for the studies of polycyclic aromatic hydrocarbons (PAHs) complexed to Fe cations that are more likely to possess d → π* and π → π* transitions in the diffuse interstellar bands domain. PAH ligands indeed possess a larger number of π and π* orbitals, respectively, higher and lower in energy than those of C6H6, which are expected to lead to lower energy d → π* and π → π* transitions in [FePAH](+) than in [FeC6H6](+) complexes.

Charge-constrained auxiliary-density-matrix methods for the Hartree-Fock exchange contribution

Three new variants of the auxiliary-density-matrix method (ADMM) of Guidon, Hutter, and VandeVondele [J. Chem. Theory Comput. 6, 2348 (2010)] are presented with the common feature that they have a simplified constraint compared with the full orthonormality requirement of the earlier ADMM1 method. All ADMM variants are tested for accuracy and performance in all-electron B3LYP calculations with several commonly used basis sets. The effect of the choice of the exchange functional for the ADMM exchange-correction term is also investigated.

Efficient calculation of integrals in mixed ramp-Gaussian basis sets

Algorithms for the efficient calculation of two-electron integrals in the newly developed mixed ramp-Gaussian basis sets are presented, alongside a Fortran90 implementation of these algorithms, RampItUp. These new basis sets have significant potential to (1) give some speed-up (estimated at up to 20% for large molecules in fully optimised code) to general-purpose Hartree-Fock (HF) and density functional theory quantum chemistry calculations, replacing all-Gaussian basis sets, and (2) give very large speed-ups for calculations of core-dependent properties, such as electron density at the nucleus, NMR parameters, relativistic corrections, and total energies, replacing the current use of Slater basis functions or very large specialised all-Gaussian basis sets for these purposes. This initial implementation already demonstrates roughly 10% speed-ups in HF/R-31G calculations compared to HF/6-31G calculations for large linear molecules, demonstrating the promise of this methodology, particularly for the second application. As well as the reduction in the total primitive number in R-31G compared to 6-31G, this timing advantage can be attributed to the significant reduction in the number of mathematically complex intermediate integrals after modelling each ramp-Gaussian basis-function-pair as a sum of ramps on a single atomic centre.

Efficient implementation of the pair atomic resolution of the identity approximation for exact exchange for hybrid and range- separated density functionals

An efficient new molecular orbital (MO) basis algorithm is reported implementing the pair atomic resolution of the identity approximation (PARI) to evaluate the exact exchange contribution (K) to self-consistent field methods, such as hybrid and range-separated hybrid density functionals. The PARI approximation, in which atomic orbital (AO) basis function pairs are expanded using auxiliary basis functions centered only on their two respective atoms, was recently investigated by Merlot et al. [J. Comput. Chem. 2013, 34, 1486]. Our algorithm is significantly faster than quartic scaling RI-K, with an asymptotic exchange speedup for hybrid functionals of (1 + X/N), where N and X are the AO and auxiliary basis dimensions. The asymptotic speedup is 2 + 2X/N for range separated hybrids such as CAM-B3LYP, ωB97X-D, and ωB97X-V which include short- and long-range exact exchange. The observed speedup for exchange in ωB97X-V for a C68 graphene fragment in the cc-pVTZ basis is 3.4 relative to RI-K. Like conventional RI-K, our method greatly outperforms conventional integral evaluation in large basis sets; a speedup of 19 is obtained in the cc-pVQZ basis on a C54 graphene fragment. Negligible loss of accuracy relative to exact integral evaluation is demonstrated on databases of bonded and nonbonded interactions. We also demonstrate both analytically and numerically that the PARI-K approximation is variationally stable.

Communication: resolving the three-body contribution to the lattice energy of crystalline benzene: benchmark results from coupled-cluster theory

Coupled-cluster theory including single, double, and perturbative triple excitations [CCSD(T)] has been applied to trimers that appear in crystalline benzene in order to resolve discrepancies in the literature about the magnitude of non-additive three-body contributions to the lattice energy. The present results indicate a non-additive three-body contribution of 0.89 kcal mol(-1), or 7.2% of the revised lattice energy of -12.3 kcal mol(-1). For the trimers for which we were able to compute CCSD(T) energies, we obtain a sizeable difference of 0.63 kcal mol(-1) between the CCSD(T) and MP2 three-body contributions to the lattice energy, confirming that three-body dispersion dominates over three-body induction. Taking this difference as an estimate of three-body dispersion for the closer trimers, and adding an Axilrod-Teller-Muto estimate of 0.13 kcal mol(-1) for long-range contributions yields an overall value of 0.76 kcal mol(-1) for three-body dispersion, a significantly smaller value than in several recent studies.