Search

We present an effective computational protocol (cLR2) to describe both solvatochromism and fluorosolvatochromism. This protocol, which couples the polarizable continuum model to time-dependent density functional theory, simultaneously accounts for both linear-response and state-specific solvation effects. A series of test cases, including solvatochromic and fluorosolvatochromic compounds and excited-state intramolecular proton transfers, are used to highlight that cLR2 is especially beneficial for modeling bright excitations possessing a significant charge-transfer character, as well as cases in which an accurate balance between states of various polarities should be restored.

The capability to determine core-electron binding energies (CEBEs) is vital in the analysis of X-ray photoelectron spectroscopy, and the continued development of light sources has made inner shell spectroscopy of heavier elements increasingly accessible. Density functional theory is widely used to determine CEBEs of lighter elements (boron-fluorine). It is shown that good performance of exchange-correlation functionals for these elements does not necessarily translate to the calculation of CEBEs for the heavier elements from the next row of the periodic table, and in general, larger errors are observed. Two strategies are explored that improve the accuracy of the calculated CEBEs. The first is to apply element and functional dependent energy corrections, and the second is a reparametrization of a short-range corrected functional. This functional is able to reproduce experimental phosphorus and sulfur K-edge CEBEs with an average error of 0.15 eV demonstrating the importance of reducing the self-interaction error associated with the core electrons and represents progress toward a density functional theory calculation that performs equally well for ionization at the K-edge of all elements.

Dinitrogen pentoxide (N2O5) is an important intermediate in the atmospheric chemistry of nitrogen oxides. Although there has been much research, the processes that govern the physical interactions between N2O5 and water are still not fully understood at a molecular level. Gaining a quantitative insight from computer simulations requires going beyond the accuracy of classical force fields while accessing length scales and time scales that are out of reach for high-level quantum-chemical approaches. To this end, we present the development of MB-nrg many-body potential energy functions for nonreactive simulations of N2O5 in water. This MB-nrg model is based on electronic structure calculations at the coupled cluster level of theory and is compatible with the successful MB-pol model for water. It provides a physically correct description of long-range many-body interactions in combination with an explicit representation of up to three-body short-range interactions in terms of multidimensional permutationally invariant polynomials. In order to further investigate the importance of the underlying interactions in the model, a TTM-nrg model was also devised. TTM-nrg is a more simplistic representation that contains only two-body short-range interactions represented through Born-Mayer functions. In this work, an active learning approach was employed to efficiently build representative training sets of monomer, dimer, and trimer structures, and benchmarks are presented to determine the accuracy of our new models in comparison to a range of density functional theory methods. By assessing the binding curves, distortion energies of N2O5, and interaction energies in clusters of N2O5 and water, we evaluate the importance of two-body and three-body short-range potentials. The results demonstrate that our MB-nrg model has high accuracy with respect to the coupled cluster reference, outperforms current density functional theory models, and thus enables highly accurate simulations of N2O5 in aqueous environments.

D.S.P. and D.J.F. thank the EPSRC for funding via Prosperity Partnership EP/S035990/1. D.S.P. and D.J.F. thank the ARCHIE-WeSt High-Performance Computing Centre (www.archie-west.ac.uk) for computational resources. The UCL authors thank Prof. Keith Refson for guidance with the phonon calculations, which used the ARCHER U.K. National Supercomputing Service (http://www.archer.ac.uk) as part of the U.K. HEC Materials Chemistry Consortium, which is funded by the EPSRC (EP/L000202, EP/R029431). R.G. was funded by MagnaPharm, a project funded by the European Union’s Horizon 2020 Research and Innovation programme under grant agreement number 736899. We demonstrate that physics-based calculations of intrinsic aqueous solubility can rival cheminformatics-based machine learning predictions. A proof-of-concept was developed for a physics-based approach via a sublimation thermodynamic cycle, building upon previous work that relied upon several thermodynamic approximations, notably the 2RT approximation, and limited conformational sampling. Here, we apply improvements to our sublimation free-energy model with the use of crystal phonon mode calculations to capture the contributions of the vibrational modes of the crystal. Including these improvements with lattice energies computed using the model-potential-based Ψmol method leads to accurate estimates of sublimation free energy. Combining these with hydration free energies obtained from either molecular dynamics free-energy perturbation simulations or density functional theory calculations, solubilities comparable to both experiment and informatics predictions are obtained. The application to coronene, succinic acid, and the pharmaceutical desloratadine shows how the methods must be adapted for the adoption of different conformations in different phases. The approach has the flexibility to extend to applications that cannot be covered by informatics methods. Publisher PDF

International audience; Tungsten trioxide (WO3)-derived nanostructures have emerged recently as feasible semiconductors for photocatalytic purposes due to their visible-light harvesting that overcomes the drawbacks presented by TiO2-derived nanoparticles (NPs). However, applications are still limited by the lack of fundamental knowledge at the nanoscale due to poor understanding of the physical processes that affect their photoactivity. To fill this gap, we report here a detailed computational study using a combined density functional theory (DFT)-GW scheme to investigate the electronic structure of realistic WO3 NPs containing up to 1680 atoms. Different phases and morphologies are considered to provide reliable structure–property relationships. Upon proper benchmark of our DFT-GW methodology, we use this highly accurate approach to establish relevant rules for the design of photoactive WO3 nanostructures by pointing out the most stable morphologies at the nanoscale and the appropriate size regime at which the photoactive efficiency is enhanced

A comprehensive overview of the performance of local hybrid functionals for molecular properties like excited states, ionization potentials within the GW framework, polarizabilities, magnetizabilities, NMR chemical shifts, and NMR spin-spin coupling constants is presented. We apply the generalization of the kinetic energy, τ, with the paramagnetic current density to all magnetic properties and the excitation energies from time-dependent density functional theory. This restores gauge invariance for these properties. Different ansatze for local mixing functions such as the iso-orbital indicator, the correlation length, the Gorling-Levy second-order limit, and the spin polarization are compared. For the latter, we propose a modified version of the corresponding hyper-generalized gradient approximation functional of Perdew, Staroverov, Tao, and Scuseria (PSTS) [Phys. Rev. A2008, 78, 052513] to allow for a numerically stable evaluation of the exchange-correlation kernel and hyperkernel. The PSTS functional leads to a very consistent improvement compared to the related TPSSh functional. It is further shown that the "best" choice of the local mixing function depends on the studied property and molecular class. While functionals based on the iso-orbital indicator lead to rather accurate excitation energies and ionization energies, the results are less impressive for NMR properties, for which a considerable dependence on the considered molecular test set and nuclei is observed. Johnson's local hybrid functional based on the correlation length yields remarkable results for NMR shifts of compounds featuring heavy elements and also for the excitation energies of organic compounds.

We present a time-dependent density functional theory (TDDFT) approach to compute the light-matter couplings between two different manifolds of excited states relative to a common ground state in the context of 4d transition metal systems. These quantities are the necessary ingredients to solve the Kramers-Heisenberg (KH) equation for resonant inelastic X-ray scattering (RIXS) and several other types of two-photon spectroscopies. The procedure is based on the pseudo-wavefunction approach, where the solutions of a TDDFT calculation can be used to construct excited-state wavefunctions, and on the restricted energy window approach, where a manifold of excited states can be rigorously defined based on the energies of the occupied molecular orbitals involved in the excitation process. Thus, the present approach bypasses the need to solve the costly TDDFT quadratic-response equations. We illustrate the applicability of the method to 4d transition metal molecular complexes by calculating the 2p4d RIXS maps of three representative ruthenium complexes and comparing them to experimental results. The method can capture all the experimental features in all three complexes to allow the assignment of the experimental peaks, with relative energies correct to within ∼0.6 eV at the cost of two independent TDDFT calculations.

We present a practical approach to treat static and dynamical correlation accurately in large multiconfigurational systems. The static correlation is taken into account by using the spin-flip approach, which is well-known for capturing static correlation accurately at low-computational expense. Unlike previous approaches to add dynamical correlation to spin-flip models which use perturbation theory or coupled-cluster theory, we explore the ability to use the on-top pair-density functional theory approaches recently developed by Gagliardi and co-workers (J. Comput. Theor. Chem., 2014, 10, 3669). External relaxations are performed in the spin-flip calculations through a restricted active space framework for which a truncation scheme for the orbitals used in the external excitation is presented. The performance of the approach is demonstrated by computing energy gaps between ground and excited states for diradicals, triradicals, and linear polyacene chains ranging from naphthalene to dodecacene. Accurate results are obtained using the new approach for these challenging open-shell molecular systems.

This paper presents a new theory called multiconfiguration density coherence functional theory (MC-DCFT). This theory provides a new route to define density functionals for multiconfiguration wave functions, in particular by using the one-particle density matrix in the coordinate representation. The theory is illustrated by calculating the dissociation curve of four heteronuclear and homonuclear diatomic molecules, namely, H2, F2, N2, and HF, using density coherence functionals converted from PBE, BLYP, and PW91. By introducing two parameters in the converted density functionals, we are able to calculate bond dissociation energies of comparable accuracy as those calculated by multiconfiguration pair-density functional theory (MC-PDFT) and complete active space second-order perturbation theory (CASPT2). This demonstrates that it would be possible to build a successful multiconfiguration density functional theory based on density coherence.

A new algorithm called Fast and Flexible CrystAl Structure Predictor (FFCASP) was developed to predict the structure of covalent and molecular crystals. FFCASP is massively parallel and able to handle more than 200 atoms in the unit cell (in other terms, it allows global optimization around 100 individual parameters). It uses a global optimizer specialized for Crystal Structure Prediction (CSP) which combines particle swarm and simulated annealing optimizers. Three different molecular crystals, including diverse intermolecular interactions, namely, cytosine, coumarin, and pyrazinamide, have been selected to evaluate the performance of FFCASP. While cytosine polymorphs have been searched by employing two different force fields (a DFT-SAPT based intermolecular potential and generalized amber force field (GAFF)) up to Z = 16, only GAFF has been used both in coumarin and pyrazinamide polymorph searches up to Z = 4. For these three molecular crystals, FFCASP generated more than 20 000 crystal structures, and the unique ones have been further treated by DFT-D3. A combination of data mining and a machine learning approach was introduced to determine the unique structures and their distribution into different clusters, which ultimately gives an opportunity to retrieve the common features and relations between the resulting structures. There are two known experimental crystal structures of cytosine, and both were successfully located with FFCASP. Two of the reported crystal structures of coumarin have been reproduced. Similarly, in pyrazinamide, three known experimental structures have been rediscovered. In addition to finding the experimentally known structures, FFCASP also located other low-energy structures for each considered molecular crystals. These successes of FFCASP offer the possibility to discover the polymorphic nature of other important molecular crystals (e.g., drugs) as well.

Collagen proteins are spread in almost every vertebrate's tissue with mechanical function. The defining feature of this fundamental family of proteins is its well-known collagen triple-helical domain. This helical domain can have different geometries, varying in helical elongation and interstrands contact, as a function of the amino acidic composition. The helical geometrical features play an important role in the interaction of the collagen protein with cell receptors, but for the vast majority of collagen compositions, these geometrical features are unknown. Quantum mechanical (QM) simulations based on the density functional theory (DFT) provide a robust approach to characterize the scenario on the collagen composition-structure relationships. In this work, we analyze the role of the adopted computational method in predicting the collagen structure for two purposes. First, we look for a cost-effective computational approach to apply to a large-scale composition-structure analysis. Second, we attempt to assess the robustness of the predictions by varying the QM methods. Therefore, we have run geometry optimization on periodic models of the collagen protein using a variety of approaches based on the most commonly used DFT functionals (PBE, HSE06, and B3LYP) with and without dispersion correction (D3ABC). We have coupled these methods with several different basis sets, looking for the highest accuracy/cost ratio. Furthermore, we have studied the performance of the composite HF-3c method and the semiempirical GFN1-xTB method. Our results identify a computational recipe that is potentially capable of predicting collagen structural features in line with DFT simulations, with orders of magnitude reduced computational cost, encouraging further investigations on the topic.

Ab initio molecular dynamics (MD) with hybrid density functionals and a plane wave basis is computationally expensive due to the high computational cost of exact exchange energy evaluation. Recently, we proposed a strategy to combine adaptively compressed exchange (ACE) operator formulation and a multiple time step integration scheme to reduce the computational cost significantly [J. Chem. Phys.2019,151, 151102 ]. However, it was found that the construction of the ACE operator, which has to be done at least once in every MD time step, is computationally expensive. In the present work, systematic improvements are introduced to further speed up by employing localized orbitals for the construction of the ACE operator. By this, we could achieve a computational speedup of an order of magnitude for a periodic system containing 32 water molecules. Benchmark calculations were carried out to show the accuracy and efficiency of the method in predicting the structural and dynamical properties of bulk water. To demonstrate the applicability, computationally intensive free-energy computations at the level of hybrid density functional theory were performed to investigate (a) methyl formate hydrolysis reaction in neutral aqueous media and (b) proton-transfer reaction within the active-site residues of the class C β-lactamase enzyme.

For the large and chemically diverse GMTKN55 benchmark suite, we have studied the performance of density-corrected density functional theory (HF-DFT), compared to self-consistent DFT, for several pure and hybrid GGA and meta-GGA exchange-correlation (XC) functionals (PBE, BLYP, TPSS, SCAN) as a function of the percentage of HF exchange in the hybrid. The D4 empirical dispersion correction has been added throughout. For subsets dominated by dynamical correlation -- particularly noncovalent interaction subsets -- HF-DFT is highly beneficial, particularly at low HF exchange percentages. For subsets with significant static correlation (i.e., where a Hartree-Fock determinant is not a good zero-order wavefunction), HF-DFT may do more harm than good. While the self-consistent series show optima at or near 37.5% (i.e., 3/8) for all four XC functionals -- consistent with Grimme's proposal of the PBE38 functional -- HF-BnLYP-D4, HF-PBEn-D4, and HF-TPSSn-D4 all exhibit minima nearer 25% (i.e., 1/4). Intriguingly, for HF-SCANn-D4, the minimum is near 10%, but the weighted mean absolute error (WTMAD2) for GMTKN55 is only barely lower than that of HF-SCAN-D4 (i.e., where the post-HF step is a pure meta-GGA). The latter becomes an attractive option, only slightly more costly than pure Hartree-Fock, and devoid of adjustable parameters other than the three in the dispersion correction. Moreover, its WTMAD2 is only surpassed by the highly empirical M06-2X and by the combinatorically optimized empirical range-separated hybrids wB97X-V and wB97M-V., Comment: Journal of Chemical Theory and Computation (2021)

The new R739×5 data set from the Non-Covalent Interactions Atlas series (www.nciatlas.org) focuses on repulsive contacts in molecular complexes, covering organic molecules, sulfur, phosphorus, halogens, and noble gases. Information on the repulsive parts of the potential energy surface is crucial for the development of robust empirically parametrized computational methods. We use the new data set of highly accurate CCSD(T)/CBS interaction energies to test selected density functional theory (DFT) and semiempirical quantum-mechanical methods. The double-hybrid functionals were the best performing, with the revDSD-PBEP86-D3 being the most accurate DFT method, followed by the range-separated ωB97X functionals. Out of semiempirical methods, GFN2-xTB yielded the best results. On the example of the PM6 method, we analyze the source of error and its relation to the difficulties in the description of conformational energies, and we also devise an immediately applicable correction that fixes the most serious uncorrected issues previously encountered in practical calculations.

Cyclic voltammograms (CVs) are a central experimental tool for assessing the structure and activity of electrochemical interfaces. Based on a mean-field ansatz for the interface energetics under applied potential conditions, we here derive an ab initio thermodynamics approach to efficiently simulate thermodynamic CVs. All unknown parameters are determined from density functional theory (DFT) calculations coupled to an implicit solvent model. For the showcased CVs of Ag(111) electrodes in halide-anion-containing solutions, these simulations demonstrate the relevance of double-layer contributions to explain experimentally observed differences in peak shapes over the halide series. Only the appropriate account of interfacial charging allows us to capture the differences in equilibrium coverage and total electronic surface charge that cause the varying peak shapes. As a case in point, this analysis demonstrates that prominent features in CVs do not only derive from changes in adsorbate structure or coverage but can also be related to variations of the electrosorption valency. Such double-layer effects are proportional to adsorbate-induced changes in the work function and/or interfacial capacitance. They are thus especially pronounced for electronegative halides and other adsorbates that affect these interface properties. In addition, the analysis allows us to draw conclusions on how the possible inaccuracy of implicit solvation models can indirectly affect the accuracy of other predicted quantities such as CVs.

Cerium (Ce) plays important roles in catalysis. Its position in the sixth period of the periodic table leads to spin-orbit coupling (SOC) and other open-shell effects that make the quantum mechanical calculation of cerium compounds challenging. In this work, we investigated the low-lying spin states of Ce+ and the bond energy of CeH+, both by multiconfigurational methods, in particular, SA-CASSCF, MC-PDFT, CASPT2, XMS-PDFT, and XMS-CASPT2, and by single-configurational methods, namely, Hartree-Fock theory and unrestricted Kohn-Sham density functional theory with 34 choices of the exchange-correlation functional. We found that only CASPT2, XMS-CASPT2, and SA-CASSCF (among the five multiconfigurational methods) and GAM, HCTH, SOGGA11, and OreLYP (among the 35 single-configuration methods) successfully predict that the SOC-free ground spin state of Ce+ is a doublet state, and CASPT2 and GAM give the most accurate multireference and single-reference calculations, respectively, of the excitation energy of the first SOC-free excited state for Ce+. We calculated that the ground doublet state of Ce+ is an intra-atomic hyper-open-shell state. We calculated the spin-orbit energy (ESO) of Ce+ by the five multiconfigurational methods and found that ESO calculated by CASPT2 is the closest to the experimental value. Taking advantage of the availability of an experimental D0 for CeH+ as a way to provide a unique test of theory, we showed that all the multiconfigurational methods overestimate D0 by at least 246 meV (5.7 kcal/mol), and only three functionals, namely, SOGGA, MN15, and GAM, have an error of D0 that is less than 200 meV (5 kcal/mol).

The calculation of harmonic vibrational frequencies (HVF) to interpret infrared (IR) spectra and to convert molecular energies to free energies is one of the essential steps in computational chemistry. A prerequisite for accurate thermostatistics so far was to optimize the molecular input structures in order to avoid imaginary frequencies, which inevitably leads to changes in the geometry if different theoretical levels are applied for geometry optimization and frequency calculations. In this work, we propose a new method termed single-point Hessian (SPH) for the computation of HVF and thermodynamic contributions to the free energy within the modified rigid-rotor-harmonic-oscillator approximation for general nonequilibrium molecular geometries. The key ingredient is the application of a biasing potential given as Gaussian functions expressed with the root-mean-square-deviation (RMSD) in Cartesian space in order to retain the initial geometry. The theory derived herein is generally applicable to quantum mechanical (QM), semiempirical QM, and force-field (FF) methods. Besides a detailed description of the underlying theory including the important back-correction of the biased HVF, the SPH approach is tested for reaction paths, molecular dynamics snapshots of crambin, and supramolecular association free energies in comparison to high-level density functional theory (DFT) values. Furthermore, the effect on IR spectra is investigated for organic dimers and transition-metal complexes revealing improved spectra at low theoretical levels. On average, DFT reference free energies are better reproduced by the newly developed SPH scheme than by conventional calculations on freely optimized geometries or without any relaxation.

Protein-embedded chromophores are responsible for light harvesting, excitation energy transfer, and charge separation in photosynthesis. A critical part of the photosynthetic apparatus are reaction centers (RCs), which comprise groups of (bacterio)chlorophyll and (bacterio)pheophytin molecules that transform the excitation energy derived from light absorption into charge separation. The lowest excitation energies of individual pigments (site energies) are key for understanding photosynthetic systems, and form a prime target for quantum chemistry. A major theoretical challenge is to accurately describe the electrochromic (Stark) shifts in site energies produced by the inhomogeneous electric field of the protein matrix. Here, we present large-scale quantum mechanics/molecular mechanics calculations of electrochromic shifts for the RC chromophores of photosystem II (PSII) using various quantum chemical methods evaluated against the domain-based local pair natural orbital (DLPNO) implementation of the similarity-transformed equation of motion coupled cluster theory with single and double excitations (STEOM-CCSD). We show that certain range-separated density functionals (ωΒ97, ωΒ97X-V, ωΒ2PLYP, and LC-BLYP) correctly reproduce RC site energy shifts with time-dependent density functional theory (TD-DFT). The popular CAM-B3LYP functional underestimates the shifts and is not recommended. Global hybrid functionals are too insensitive to the environment and should be avoided, while nonhybrid functionals are strictly nonapplicable. Among the applicable approximate coupled cluster methods, the canonical versions of CC2 and ADC(2) were found to deviate significantly from the reference results both for the description of the lowest excited state and for the electrochromic shifts. By contrast, their spin-component-scaled (SCS) and particularly the scale-opposite-spin (SOS) variants compare well with the reference DLPNO-STEOM-CCSD and the best range-separated DFT methods. The emergence of RC excitation asymmetry is discussed in terms of intrinsic and protein electrostatic potentials. In addition, we evaluate a minimal structural scaffold of PSII, the D1-D2-CytB559 RC complex often employed in experimental studies, and show that it would have the same site energy distribution of RC chromophores as the full PSII supercomplex, but only under the unlikely conditions that the core protein organization and cofactor arrangement remain identical to those of the intact enzyme.

Funding Information: This work was supported by the Academy of Finland (project numbers 308647, 314298, and 316601) and through their Flagship program: Finnish Center for Artificial Intelligence FCAI. We thank CSC, the Finnish IT Center for Science, and Aalto Science IT for computational resources. This work is supported by COST (European Cooperation in Science and Technology) Action 18234. L.F. thanks Guoxu Zhang, Marc Dvorak, Jingrui Li, and Annika Stuke for the help with FHI-Aims. He also acknowledges financial support from the Chinese Scholarship Council (grant no. [2017]3109). 1 Publisher Copyright: © 2020 American Chemical Society. All rights reserved. Copyright: Copyright 2021 Elsevier B.V., All rights reserved. Finding low-energy molecular conformers is challenging due to the high dimensionality of the search space and the computational cost of accurate quantum chemical methods for determining conformer structures and energies. Here, we combine active-learning Bayesian optimization (BO) algorithms with quantum chemistry methods to address this challenge. Using cysteine as an example, we show that our procedure is both efficient and accurate. After only 1000 single-point calculations and approximately 80 structure relaxations, which is less than 10% computational cost of the current fastest method, we have found the low-energy conformers in good agreement with experimental measurements and reference calculations. To test the transferability of our method, we also repeated the conformer search of serine, tryptophan, and aspartic acid. The results agree well with previous conformer search studies.

Photoredox catalysts (PCs) have contributed to the advancement of organic chemistry by accelerating conventional reactions and enabling new pathways through the use of reactive electrons in excited states. With a number of successful applications, chemists continue to seek new promising organic PCs to achieve their objectives. Instead of labor-intensive manual experimentation, quantum chemical simulations could explore the enormous chemical space more efficiently. The reliability and accuracy of quantum chemical simulations have become important factors for material screening. We designed a theoretical protocol capable of predicting redox properties in excited states with high accuracy for a selected model system of dihydroquinoxalino[2,3-b]quinoxaline derivatives. Herein, three factors were considered as critical to achieving reliable predictions with accurate physics: the solvent medium effect on excited-state geometries, an adequate amount of Hartree-Fock exchange (HFX), and the consideration of double-excitation character in excited states. We determined that it is necessary to incorporate solvent medium during geometry optimizations to obtain planar excited-state structures that are consistent with the experimentally observed modest Stokes shift. While density functionals belonging to the generalized gradient approximation family perform well for the prediction of photoelectrochemical properties, an incorrect description of exciton boundedness (spontaneous dissociation of excitons or extremely weak boundedness) on small organic molecules was predicted. The inclusion of an adequate amount of Hartree-Fock exchange was suggested as one approach to obtain bound excitons, which is physically reasonable. The last consideration is the double-excitation character in S1 states. As revealed by the second-order algebraic diagrammatic construction theory, non-negligible double excitations exist in S1 states in our model systems. Time-dependent density functional theory (TDDFT) is blind to doubly excited states, and this motivated us to use spin-flip DFT (SF-DFT). We established a theoretical protocol that could provide highly accurate estimations of photophysical properties and ground-/excited-state redox properties, focusing on the three factors mentioned above. Geometry optimization with DFT and TDDFT employing the B3LYP functional (20% HFX) in solution and energy refinement by SF-DFT reproduced the experimental redox properties in the excited and ground states remarkably well with mean signed deviations (MSDs) of 0.01 and -0.15 V, respectively. This theoretical protocol is expected to contribute to the understanding of exciton behavior in organic PCs and to the efficient design of new promising PC candidates.

The many-body expansion (MBE) of energies of molecular clusters or solids offers a way to detect and analyze errors of theoretical methods that could go unnoticed if only the total energy of the system was considered. In this regard, the interaction between the methane molecule and its enclosing dodecahedral water cage, CH$_4$(H$_2$O)$_{20}$, is a stringent test for approximate methods, including density-functional theory (DFT) approximations. Hybrid and semilocal DFT approximations behave erratically for this system, with three- and four-body nonadditive terms having neither the correct sign nor magnitude. Here we analyze to what extent these qualitative errors in different MBE contributions are conveyed to post-Kohn-Sham random-phase approximation (RPA). The results reveal a correlation between the quality of the DFT input states and the RPA results. Moreover, the renormalized singles energy (RSE) corrections play a crucial role in all orders of MBE. For dimers, RSE corrects the RPA underbinding for every tested Kohn-Sham model: generalized-gradient approximation (GGA), meta-GGA, (meta-)GGA hybrids, as well as the optimized effective potential at the correlated level. Remarkably, the inclusion of singles in RPA can also correct the wrong signs of three- and four-body nonadditive energies as well as mitigate the excessive higher-order contributions to the MBE. The RPA errors are dominated by the contributions of compact clusters. As a workable method for large systems, we propose to replace those compact contributions with CCSD(T) energies and to sum up the remaining many-body contributions up to infinity with supermolecular or periodic RPA. As a demonstration of this approach, we show that for RPA(PBE0)+RSE it suffices to apply CCSD(T) to dimers and 30 compact, hydrogen-bonded trimers to get the methane-water cage interaction energy to within 1.6% of the reference value., Comment: Data available at https://doi.org/10.5281/zenodo.4429677

While CCSD(T) with spin-restricted Hartree-Fock (RHF) orbitals has long been lauded for its ability to accurately describe closed-shell interactions, the performance of CCSD(T) on open-shell species is much more erratic, especially when using a spin-unrestricted HF (UHF) reference. Previous studies have shown improved treatment of open-shell systems when a non-HF set of molecular orbitals, like Brueckner or Kohn-Sham density functional theory (DFT) orbitals, is used as a reference. Inspired by the success of regularized orbital-optimized second-order Møller-Plesset perturbation theory (κ-OOMP2) orbitals as reference orbitals for MP3, we investigate the use of κ-OOMP2 orbitals and various DFT orbitals as reference orbitals for CCSD(T) calculations of the corrected ground-state harmonic vibrational frequencies of a set of 36 closed-shell (29 neutrals, 6 cations, 1 anion) and 59 open-shell diatomic species (38 neutrals, 15 cations, 6 anions). The aug-cc-pwCVTZ basis set is used for all calculations. The use of κ-OOMP2 orbitals in this context alleviates difficult cases observed for both UHF orbitals and OOMP2 orbitals. Removing two multireference systems and 12 systems with ambiguous experimental data leaves a pruned data set. Overall performance on the pruned data set highlights CCSD(T) with a B97 orbital reference (CCSD(T):B97), CCSD(T) with a κ-OOMP2 orbital reference (CCSD(T):κ-OOMP2), and CCSD(T) with a B97M-rV orbital reference (CCSD(T):B97M-rV) with RMSDs of 8.48 cm-1, and 8.50 cm-1, and 8.75 cm-1 respectively, outperforming CCSD(T):UHF by nearly a factor of 5. Moreover, the performance on the closed- and open-shell subsets shows these methods are able to treat open-shell and closed-shell systems with comparable accuracy and robustness. CCSD(T) with RHF orbitals is seen to improve upon UHF for the closed-shell species, while spatial symmetry breaking in a number of restricted open-shell HF (ROHF) references leads CCSD(T) with ROHF reference orbitals to exhibit the poorest statistical performance of all methods surveyed for open-shell species. The use of κ-OOMP2 orbitals has also proven useful in diagnosing multireference character that can hinder the reliability of CCSD(T).

A proper treatment of electron correlation effects is indispensable for accurate simulation of compounds. Various post-Hartree-Fock methods have been adopted to calculate correlation energies of chemical systems, but time complexity usually prevents their usage in a large scale. Here, we propose a density functional approximation, based on machine learning using neural networks, which can be readily employed to produce results comparable to second-order Møller-Plesset perturbation (MP2) ones for organic compounds with reduced computational cost. Various systems have been tested and the transferability across basis sets, structures, and nuclear configurations has been evaluated. Only a small number of molecules at the equilibrium structure has been needed for the training, and generally less than 5% relative error has been achieved for structures outside the training domain and systems containing about 140 atoms. In addition, this approach has been applied to make predictions on correlation energies of nuclear configurations extracted from density functional theory-based molecular dynamics trajectories with only one or two structures as training data.

A simple and robust range-separated (RS) double-hybrid (DH) time-dependent density functional approach is presented for the accurate calculation of excitation energies of molecules within the Tamm-Dancoff approximation. The scheme can be considered as an excited-state extension of the ansatz proposed by Toulouse and co-workers [J. Chem. Phys. 2018, 148, 164105], which is based on the two-parameter decomposition of the Coulomb potential, for which both the exchange and correlation contributions are range-separated. A flexible and efficient implementation of the new scheme is also presented, which facilitates its extension to any combination of exchange and correlation functionals. The performance of the new approximation is tested for singlet excitations on several benchmark compilations and thoroughly compared to that of representative DH, RS hybrid, and RS DH functionals. The one-electron basis set dependence and computation times are also assessed. Our results show that the new approach improves on standard DHs in most cases, and it can provide a more robust and accurate alternative. In addition, on average, it noticeably surpasses the existing RS hybrid and RS DH functionals.

The premise for Kohn-Sham density functional theory is that the exact electron density can be generated by a set of orbitals in a single Slater determinant. While this is ensured in a complete basis set, it has been shown that it cannot hold in small basis sets. The present work probes how accurately a reference electron density of the full-CI type can be reproduced by a set of orbitals in a single Slater determinant, as a function of the basis set used for the fitting electron density. The key finding is that the fitting error may be significant for basis sets of double- or triple-ζ quality. It is also shown that it is important that the fitting basis set includes the same basis functions as used for generating the reference electron density. The main limitation in a given basis set is the lack of higher order polarization functions. The error for practical purposes becomes insignificant for basis sets of quadruple-ζ or better quality, and this should be the choice when assessing the accuracy of exchange-correlation functionals by comparing electron densities to accurate reference results generated by wave function methods. The methodology in the present work can be used to transform an electron density from a multideterminant wave function into a set of orbitals in a single Slater determinant, and this may be useful for developing and testing new exchange-correlation functionals.

Connections between the Görling-Levy (GL) perturbation theory and the parameters of double-hybrid (DH) density functional are established via adiabatic connection formalism. Moreover, we present a more general DH density functional theory, where the higher-order perturbation terms beyond the second-order GL2 one, such as GL3 and GL4, also contribute. It is shown that a class of DH functionals including previously proposed ones can be formed using the present construction. Based on the proposed formalism, we assess the performance of higher-order DH and long-range corrected DH formed on the Perdew-Burke-Ernzerhof (PBE) semilocal functional and second-order GL2 correlation energy. The underlying construction of DH functionals based on the generalized many-body perturbation approaches is physically appealing in terms of the development of the non-local forms using more accurate and sophisticated semilocal functionals.

Interpreting response properties such as the polarizability, optical rotation (OR), or hyperpolarizabilities is a complex task for which a uniform strategy would be desirable. We propose a response analysis procedure called the RespA approach with two interrelated schemes to describe molecular optical response properties in terms of natural response orbitals (NROs) and chemical fragment response for convenient elucidation of structure-(optical)property relationships. These quantities can be easily extracted from the frequency-dependent perturbed one-electron transition/current density matrix obtained from any quantum mechanical response function calculation. NROs provide the most compact representation of the virtual excitations occurring in the (hyper)scattering process. It is decomposed in hole and electron NRO pairs providing a simple exciton picture. For a chemist, it is natural to interpret a property by decomposing it into functional groups or fragment contributions. In this spirit, the response is partitioned into on-site and between-site fragment responses, allowing a property mapping into real space. The new RespA procedure was implemented and tested at the simplified time-dependent density functional theory (sTD-DFT) level enabling calculations for large systems. The RespA strategy is a method-independent route for the understanding of a wide variety of response properties. We showcase how the chemically intuitive RespA approach extracts easily structure-property relationships for the particularly difficult case of OR. As examples, we demonstrate how to enhance the OR of [5]helicene and norbornenone, provide an analysis of the change of the OR observed for camphor and fenchone, and finally investigate the case of a (P,P)-bis-helicenic 2,2'-bipyridine chiroptical switch.

Inelastic neutron scattering (INS) is a widely used technique to study atomic and molecular vibrations. With the increasing complexity of materials and thus the INS spectra, being able to simulate the spectra from various atomistic models becomes an essential step and also a major bottleneck for INS data analysis. The conventional approach using density functional theory and lattice dynamics often falls short when the materials of interest are complex (e.g., defective, disordered, heterogeneous, amorphous, large-scale), for which molecular dynamics driven by an interatomic force field is a more common approach. In this paper, we demonstrate a method to directly convert molecular dynamics trajectories into simulated INS spectra, including not only fundamental but also higher order excitations. The results are compared with data collected on various representative samples from different neutron spectrometers. This development will open great opportunities by providing the key tool to perform in-depth analysis of INS data and to validate and optimize computer models.

This paper is a study of density coherences in multiconfiguration self-consistent field theory and Kohn-Sham (KS) density functional theory. We visualize and compare the nondiagonal elements of the first-order reduced density matrix in the electronic coordinate representation. This is the electronic density coherence, and it also appears in Hartree-Fock (HF) theory in the integrand of the exchange integral. The density coherence is calculated as a function of the internuclear distance for three diatomic molecules (H2, F2, and HF) using both restricted and unrestricted KS and HF theory, as well as the complete active space self-consistent field method. We identify a group of closely associated peaks corresponding to the coherence of electrons on opposite sides of the center of the molecule, and we call this the exchange massif. We find that Slater-determinant methods with a higher percentage of HF exchange tend to underestimate the density coherence at the exchange massif. We explain the trends in terms of a multireference diagnostic and the number of unpaired electrons.

The development of variational density functional theory approaches to excited electronic states is impeded by limitations of the commonly used self-consistent field (SCF) procedure. A method based on a direct optimization approach as well as the maximum overlap method is presented and the performance compared with previously proposed SCF strategies. Excited-state solutions correspond to saddle points of the energy as a function of the electronic degrees of freedom. The approach presented here makes use of a preconditioner determined with the help of the maximum overlap method to guide the convergence on a target nth-order saddle point. The method is found to be more robust and to converge faster than previously proposed SCF approaches for a set of 89 excited states of molecules. A limited-memory formulation of the symmetric rank-one method for updating the inverse Hessian is found to give the best performance. A conical intersection for the carbon monoxide molecule is calculated without resorting to fractional occupation numbers. Calculations on excited states of the hydrogen atom and a doubly excited state of the dihydrogen molecule using a self-interaction corrected functional are presented. For these systems, the self-interaction correction is found to improve the accuracy of density functional calculations of excited states.

The binding energies for cation-π complexation are underestimated by traditional fixed-charge force fields owing to their lack of explicit treatment of ion-induced dipole interactions. To address this deficiency, an explicit treatment of cation-π interactions has been introduced into the OPLS-AA force field. Following prior work with atomic cations, it is found that cation-π interactions can be handled efficiently by augmenting the usual 12-6 Lennard-Jones potentials with 1/r4 terms. Results are provided for prototypical complexes as well as protein-ligand systems of relevance for drug design. Alkali cation, ammonium, guanidinium, and tetramethylammonium were chosen for the representative cations, while benzene and six heteroaromatic molecules were used as the π systems. The required nonbonded parameters were fit to reproduce structure and interaction energies for gas-phase complexes from density functional theory (DFT) calculations at the ωB97X-D/6-311++G(d,p) level. The impact of the solvent was then examined by computing potentials of mean force (pmfs) in both aqueous and tetrahydrofuran (THF) solutions using the free-energy perturbation (FEP) theory. Further testing was carried out for two cases of strong and one case of weak cation-π interactions between druglike molecules and their protein hosts, namely, the JH2 domain of JAK2 kinase and macrophage migration inhibitory factor. FEP results reveal greater binding by 1.5-4.4 kcal/mol from the addition of the explicit cation-π contributions. Thus, in the absence of such treatment of cation-π interactions, errors for computed binding or inhibition constants of 101-103 are expected.

Understanding charge transport through molecular wires is important for nanoscale electronics and biochemistry. Our goal is to establish a simple first-principles protocol for predicting the charge transport mechanism in such wires, in particular the crossover from coherent tunneling for short wires to incoherent hopping for longer wires. This protocol is based on a combination of density functional theory with a polarizable continuum model introduced by Kaupp et al. for mixed-valence molecules, which we had previously found to work well for length-dependent charge delocalization in such systems. We combine this protocol with a new charge delocalization measure tailored for molecular wires, and we show that it can predict the tunneling-to-hopping transition length with a maximum error of one subunit in five sets of molecular wires studied experimentally in molecular junctions at room temperature. This suggests that the protocol is also well suited for estimating the extent of hopping sites as relevant, for example, for the intermediate tunneling-hopping regime in DNA.

In this work, we explore the accuracy of post-Hartree–Fock (HF) methods and double-hybrid density functional theory (DFT) for the computation of solid-state NMR chemical shifts. We apply an embedded cluster approach and investigate the convergence with cluster size and embedding for a series of inorganic solids with long-range electrostatic interactions. In a systematic study, we discuss the cluster design, the embedding procedure, and basis set convergence using gauge-including atomic orbital (GIAO) NMR calculations at the DFT and MP2 levels of theory. We demonstrate that the accuracy obtained for the prediction of NMR chemical shifts, which can be achieved for molecular systems, can be carried over to solid systems. An appropriate embedded cluster approach allows one to apply methods beyond standard DFT even for systems for which long-range electrostatic effects are important. We find that an embedded cluster should include at least one sphere of explicit neighbors around the nuclei of interest, given that a sufficiently large point charge and boundary effective potential embedding is applied. Using the pcSseg-3 basis set and GIAOs for the computation of nuclear shielding constants, accuracies of 1.6 ppm for 7Li, 1.5 ppm for 23Na, and 5.1 ppm for 39K as well as 9.3 ppm for 19F, 6.5 ppm for 35Cl, 7.4 ppm for 79Br, and 7.5 ppm for 25Mg as well as 3.8 ppm for 67Zn can be achieved with MP2. Comparing various DFT functionals with HF and MP2, we report the superior quality of results for methods that include post-HF correlation like MP2 and double-hybrid DFT.

We re-evaluate two key steps in the mechanism of CO2 reduction on copper at a higher level of theory capable of correcting inherent errors in density functional theory (DFT) approximations, namely,...

Computational simulation of nonadiabatic molecular dynamics is an indispensable tool for understanding complex photoinduced processes such as internal conversion, energy transfer, charge separation, and spatial localization of excitons, to name a few. We report an implementation of the fewest-switches surface-hopping algorithm in the NWChem computational chemistry program. The surface-hopping method is combined with linear-response time-dependent density functional theory calculations of adiabatic excited-state potential energy surfaces. To treat quantum transitions between arbitrary electronic Born-Oppenheimer states, we have implemented both numerical and analytical differentiation schemes for derivative nonadiabatic couplings. A numerical approach for the time-derivative nonadiabatic couplings together with an analytical method for calculating nonadiabatic coupling vectors is an efficient combination for surface-hopping approaches. Additionally, electronic decoherence schemes and a state reassigned unavoided crossings algorithm are implemented to improve the accuracy of the simulated dynamics and to handle trivial unavoided crossings. We apply our code to study the ultrafast decay of photoexcited benzene, including a detailed analysis of the potential energy surface, population decay timescales, and vibrational coordinates coupled to the excitation dynamics. We also study the photoinduced dynamics in trans-distyrylbenzene. This study provides a baseline for future implementations of higher-level frameworks for simulating nonadiabatic molecular dynamics in NWChem.

Proton transfer at metal oxide/water interfaces plays an important role in electrochemistry, geochemistry, and environmental science. The key thermodynamic quantity to characterize this process is the surface acidity constant. An ab initio method that combines density functional theory-based molecular dynamics (DFTMD) and free energy perturbation theory has been established for computing surface acidity constants. However, it involves a reversible proton insertion procedure in which frequent proton hopping, e.g., for strong bases and some oxide surfaces (e.g., SnO2), can cause instability issues in electronic structure calculation. In the original implementation, harmonic restraining potentials are imposed on all O-H bonds (denoted by the VrH scheme) to prevent proton hopping and thus may not be applicable for systems involving spontaneous proton hopping. In this work, we introduce an improved restraining scheme with a repulsive potential Vrep to compute the surface acidities of systems in which proton hopping is spontaneous and fast. In this Vrep scheme, a Buckingham-type repulsive potential Vrep is applied between the deprotonation site and all other protons in DFTMD simulations. We first verify the Vrep scheme by calculating the pKa values of H2O and aqueous HS- solution (i.e., strong conjugate bases) and then apply it to the SnO2(110)/H2O interface. It is found that the Vrep scheme leads to a prediction of the point of zero charge (PZC) of 4.6, which agrees well with experiment. The intrinsic individual pKa values of the terminal five-coordinated Sn site (Sn5cOH2) and bridge oxygen site (Sn2ObrH+) are 4.4 and 4.7, respectively, both being almost the same as the PZC. The similarity of the two pKa values indicates that dissociation of terminal water has almost zero free energy at this proton hopping interface (i.e., partial water dissociation), as expected from the acid-base equilibrium on SnO2.

By using a combination of classical Hamiltonian replica exchange with high-level quantum mechanical calculations on more than one hundred drug-like molecules, we explored here the energy cost associated with binding of drug-like molecules to target macromolecules. We found that, in general, the drug-like molecules present bound to proteins in the Protein Data Bank (PDB) can access easily the bioactive conformation and in fact for 73% of the studied molecules the "bioactive" conformation is within 3kBT from the most-stable conformation in solution as determined by DFT/SCRF calculations. Cases with large differences between the most-stable and the bioactive conformations appear in ligands recognized by ionic contacts, or very large structures establishing many favorable interactions with the protein. There are also a few cases where we observed a non-negligible uncertainty related to the experimental structure deposited in PDB. Remarkably, the rough automatic force field used here provides reasonable estimates of the conformational ensemble of drugs in solution. The outlined protocol can be used to better estimate the cost of adopting the bioactive conformation.

The study of photoinduced dynamics in chemical systems necessitates accurate and computationally efficient electronic structure methods, especially as the systems of interest grow larger. The linea...

Quadrupolar NMR relaxation rates were computed for aqueous 133Cs+, 131Xe, and 127I– via Kohn–Sham (KS) density functional theory-based ab initio molecular dynamics and KS calculations of the electr...

We systematically investigate an impact of the size and content of a quantum (QM) region, treated at the density functional theory level, in embedding calculations on one- (OPA) and two-photon absorption (TPA) spectra of the following fluorescent proteins (FPs) models: Aequorea victoria green FP (avGFP) with neutral (avGFP-n) and anionic (avGFP-a) chromophore as well as Citrine FP. We find that amino acid (a.a.) residues as well as water molecules hydrogen-bonded (h-bonded) to the chromophore usually boost both OPA and TPA processes intensity. The presence of hydrophobic a.a. residues in the quantum region also non-negligibly affects both absorption spectra but decreases absorption intensity. We conclude that to reach a quantitative description of OPA and TPA spectra in multiscale modeling of FPs, the quantum region should consist of a chromophore and most of a.a. residues and water molecules in a radius of 0.30-0.35 nm (ca. 200-230 atoms) when the remaining part of the system is approximated by the electrostatic point-charges. The optimal size of the QM region can be reduced to 80-100 atoms by utilizing a more advanced polarizable embedding model. We also find components of the QM region that are specific to a FP under study. We propose that the F165 a.a. residue is important in tuning the TPA spectrum of avGFP-n but not other investigated FPs. In the case of Citrine, Y203 and M69 a.a. residues must definitely be part of the QM subsystem. Furthermore, we find that long-range electrostatic interactions between the QM region and the rest of the protein cannot be neglected even for the most extensive QM regions (ca. 350 atoms).

Orbital-free approaches might offer a way to boost the applicability of density functional theory by orders of magnitude in system size. An important ingredient for this endeavor is the kinetic energy density functional. Snyder et al. [Phys. Rev. Lett.2012, 108, 25300223004593] presented a machine learning approximation for this functional achieving chemical accuracy on a one-dimensional model system. However, a poor performance with respect to the functional derivative, a crucial element in iterative energy minimization procedures, enforced the application of a computationally expensive projection method. In this work we circumvent this issue by including the functional derivative into the training of various machine learning models. Besides kernel ridge regression, the original method of choice, we also test the performance of convolutional neural network techniques borrowed from the field of image recognition.

The fundamental ideas for a non-local density functional theory -- capable of reliably capturing van der Waals interaction -- were already conceived in the 1990's. In 2004, a seminal paper introduced the first practical non-local exchange-correlation functional called vdW-DF, which has become widely successful and laid the foundation for much further research. However, since then, the functional form of vdW-DF has remained unchanged. Several successful modifications paired the original functional with different (local) exchange functionals to improve performance and the successor vdW-DF2 also updated one internal parameter. Bringing together different insights from almost two decades of development and testing, we present the next-generation non-local correlation functional called vdW-DF3, in which we change the functional form while staying true to the original design philosophy. Although many popular functionals show good performance around the binding separation of van der Waals complexes, they often result in significant errors at larger separations. With vdW-DF3, we address this problem by taking advantage of a recently uncovered and largely unconstrained degree of freedom within the vdW-DF framework that can be constrained through empirical input, making our functional semi-empirical. For two different parameterizations, we benchmark vdW-DF3 against a large set of well-studied test cases and compare our results with the most popular functionals, finding good performance in general for a wide array of systems and a significant improvement in accuracy at larger separations. Finally, we discuss the achievable performance within the current vdW-DF framework, the flexibility in functional design offered by vdW-DF3, as well as possible future directions for non-local van der Waals density functional theory.

The projector augmented wave (PAW) method is a powerful numerical algorithm that serves as a backend, enabling efficient density functional theory (DFT) calculations through the smoothing of valence electronic descriptions. Although it is mainly used in conjunction with plane-wave basis for solid-state systems, its generality permits the combination with other types of basis functions. In the previous study, we proposed a scheme to incorporate the PAW method into the conventional quantum chemical DFT implementation based on Gauss-type function (GTF) basis (Xiong et al.

Efficient and high-capacity storage materials are indispensable for a hydrogen-based economy. In silico tools can accelerate the process of discovery of new adsorbent materials with optimal hydrogen adsorption enthalpies. Density functional theory is well-poised to become a very useful tool for enabling high-throughput screening of potential materials. In this work, we have identified density functional approximations that provide good performance for hydrogen binding applications following a two-pronged approach. First, we have compiled a data set (H2Bind275) that comprehensively represents the hydrogen binding problem capturing the chemical and mechanistic diversity in the binding sites encountered in hydrogen storage materials. We have also computed reference interaction energies for this data set using coupled-cluster theory. Second, we have assessed the performance of 55 density functional approximations for predicting H2 interaction energies and have identified two hybrid density functionals (ωB97X-V and ωB97M-V), two double hybrid density functionals (DSD-PBEPBE-D3(BJ) and PBE0-DH), and one semilocal density functional (B97M-V) as the best performing ones. We have recommended the addition of empirical dispersion corrections to systematically underbinding density functionals such as revPBE, BLYP, and B3LYP for improvements in performance at negligible additional cost. We have also recommended the usage of the def2-TZVPP basis set as it represents a good compromise between accuracy and cost, limiting the finite basis set errors to less than 1 kJ/mol.

Hematite is a possible photoanode for photoelectrochemical cells (PECs) that is widely studied using density functional theory (DFT). In this paper we perform more accurate calculations of the absorption spectrum of hematite using the one-shot Green's function (G0W0) and Bethe-Salpeter equation (BSE) methods, which take excited states into account and compare the spectrum to experimental data. We found a match between our calculations and the observed absorption spectra in peak locations. Furthermore, there is anisotropy of the absorption spectra that is concurrent with the crystal structure. We also calculated the absorption spectrum of hematite intermediates during catalysis of the oxygen evolution reaction to better understand which intermediate is dominant during the reaction and the contribution of excited states to catalysis. The *O intermediate was found to be the most optically and chemically dominant species during catalysis.

We present a graph theoretic approach to adaptively compute contributions from many-body approximations in an efficient manner and perform accurate hybrid density functional theory (DFT) electronic structure calculations for condensed-phase systems. The salient features of the approach are ONIOM-like. (a) The full-system calculation is performed at a lower level of theory (pure DFT) by enforcing periodic boundary conditions. (b) This treatment is then improved through a correction term that captures many-body interactions up to any given order within a higher (in this case, hybrid DFT) level of theory. (c) In the spirit of ONIOM, contributions from the many-body approximations that arise from the higher level of theory [i.e., from (b) above] are included through extrapolation from the lower level calculation. The approach is implemented in a general, system-independent, fashion using the graph-theoretic functionalities available within Python. For example, the individual one-body components within the unit cell are designated as "nodes" within a graph. The interactions between these nodes are captured through edges, faces, tetrahedrons, and so forth, thus building a many-body interaction hierarchy. Electronic energy extrapolation contributions from all of these geometric entities are included within the above-mentioned ONIOM paradigm. The implementation of the method simultaneously uses multiple electronic structure packages. Here, for example, we present results which use both the Gaussian suite of electronic structure programs and the Quantum ESPRESSO program within a single calculation. Thus, the method integrates both plane-wave basis functions and atom-centered basis functions within a single structure calculation. The method is demonstrated for a range of condensed-phase problems for computing (i) hybrid DFT energies for condensed-phase systems at pure DFT cost and (ii) large, triple-zeta, multiply polarized, and diffuse atom-centered basis-set energies at computational costs commensurate with much smaller sets of basis functions. The methods are demonstrated through calculations performed on (a) homogeneous water surfaces as well as heterogeneous surfaces that contain organic solutes studied using two-dimensional periodic boundary conditions and (b) bulk simulations of water enforced through three-dimensional periodic boundary conditions. A range of structures are considered, and in all cases, the results are in good agreement with those obtained using large atom-centered basis and hybrid DFT calculations on the full periodic systems at significantly lower cost.

The absorption spectra of the octahedral [UX6]n- (X = Cl, Br; n = 1, 2) complexes in the near-infrared (NIR) and UV-vis spectral regions were studied theoretically, using a relativistic restricted active space second-order perturbation theory (RASPT2) wavefunction framework, with the spin-orbit (SO) coupling treated by state interaction, in conjunction with Kohn-Sham density functional theory calculations for determining the vibrational normal modes. The electric-dipole-allowed ligand-to-metal charge-transfer (LMCT) and 5f-to-6d transitions, and the electric-dipole-forbidden 5f-to-5f ligand field (LF) transitions, are thereby obtained within the same theoretical framework. For the 5f-to-5f LF transitions, the observed absorption intensity is mostly due to vibronic coupling with low-energy electric-dipole-allowed transitions, but in some cases, the magnetic dipole intensity of the purely electronic transition has comparable intensity to the vibronic transitions. Experimental LF spectra of 5f2 open-shell systems have been reported decades back, but ab initio calculations of their vibronic intensity have not yet been reported in the literature. Although the LF spectra for the 5f2 systems can be assigned in detail, based on the calculations, the spectra are very complex and the underlying electronic states are strongly multiconfigurational. Therefore, the usefulness of the LF spectra beyond serving as a "fingerprint" of the LF and the metal oxidation state appears to be limited.

We present the details of a GPU capable exchange correlation (XC) scheme integrated into the open source QUantum Interaction Computational Kernel (QUICK) program. Our implementation features an octree based numerical grid point partitioning scheme, GPU enabled grid pruning and basis/primitive function prescreening and fully GPU capable XC energy and gradient algorithms. Benchmarking against the CPU version demonstrated that the GPU implementation is capable of delivering an impres- sive performance while retaining excellent accuracy. For small to medium size protein/organic molecular systems, the realized speedups in double precision XC energy and gradient computation on a NVIDIA V100 GPU were 60 to 80-fold and 140 to 780- fold respectively as compared to the serial CPU implementation. The acceleration gained in density functional theory calculations from a single V100 GPU significantly exceeds that of a modern CPU with 40 cores running in parallel.

The relationships between Kohn-Sham (KS) and generalized KS (GKS) density functional orbital energies and fundamental gaps or optical gaps raise many interesting questions including the physical meanings of KS and GKS orbital energies when computed with currently available approximate density functionals (ADFs). In this work, by examining three diverse databases with various ADFs, we examine such relations from the point of view of the exciton shift of quasiparticle theory. We start by calculating a large number of excitation energies by time-dependent density functional theory (TDDFT) with a large number of ADFs. To relate the exciton shift implicit in TDDFT to the exciton shift that is explicit in Green's function theory, we define the exciton shift in TDDFT as the difference of the response shift and the quasiparticle shift. We found a strong correlation between the response shift and the amount of Hartree-Fock exchange included in the density functional, with the response shift varying between -1 and 5 eV. This range is an order of magnitude larger than the mean errors of the TDDFT excitation energies. This result suggests that, with currently available functionals, the KS or GKS orbital energies should be treated as intermediate mathematical variables in the calculation of excitation energies rather than as the energies of independent-particle reference states for quasiparticle theory.

Computation of the ionization energies and of the respective Dyson orbitals based on the use of the extended Koopmans theorem (EKT) is implemented in connection with an ensemble density functional theory (eDFT) method, the state-interaction state-averaged spin-restricted ensemble-referenced Kohn–Sham (SI-SA-REKS or SSR) method. The new methodology enables fast computation of the ionization energies and evaluation of the respective Dyson orbitals, the square norms of which are related with the ionization probabilities, in the ground and excited electronic states of molecules. As the application of EKT recycles the intermediate quantities from the SSR analytical energy gradient, evaluation of the ionization energies and probabilities can be carried out on-the-fly during the nonadiabatic molecular dynamics simulations. This opens up a perspective for fast theoretical simulation of the time-resolved photoelectron spectroscopy observations. In the present work, the new methodology is tested in the computation of the ionization energies and Dyson orbitals of several molecules in the ground and excited electronic states, including strongly correlated species, such as the ozone molecule, dissociating chemical bonds, and conical intersections.