The random phase approximation applied to ice

Machine Learning Density Functionals from the Random-Phase Approximation

Kohn-Sham density functional theory (DFT) is the standard method for first-principles calculations in computational chemistry and materials science. More accurate theories such as the random-phase approximation (RPA) are limited in application due to their large computational cost. Here, we use machine learning to map the RPA to a pure Kohn-Sham density functional. The machine learned RPA model (ML-RPA) is a nonlocal extension of the standard gradient approximation. The density descriptors used as ingredients for the enhancement factor are nonlocal counterparts of the local density and its gradient. Rather than fitting only RPA exchange-correlation energies, we also include derivative information in the form of RPA optimized effective potentials. We train a single ML-RPA functional for diamond, its surfaces, and liquid water. The accuracy of ML-RPA for the formation energies of 28 diamond surfaces reaches that of state-of-the-art van der Waals functionals. For liquid water, however, ML-RPA cannot yet improve upon the standard gradient approximation. Overall, our work demonstrates how machine learning can extend the applicability of the RPA to larger system sizes, time scales, and chemical spaces.

Identification of synthesisable crystalline phases of water – a prototype for the challenges of computational materials design

We discuss the identification of experimentally realisable crystalline phases of water to outline and contextualise some of the diverse building blocks of a computational materials design process.

Overcoming the difficulties of predicting conformational polymorph energetics in molecular crystals via correlated wavefunction methods

Molecular crystal structure prediction is increasingly being applied to study the solid form landscapes of larger, more flexible pharmaceutical molecules. Despite many successes in crystal structure prediction, van der Waals-inclusive density functional theory (DFT) methods exhibit serious failures predicting the polymorph stabilities for a number of systems exhibiting conformational polymorphism, where changes in intramolecular conformation lead to different intermolecular crystal packings. Here, the stabilities of the conformational polymorphs of o-acetamidobenzamide, ROY, and oxalyl dihydrazide are examined in detail. DFT functionals that have previously been very successful in crystal structure prediction perform poorly in all three systems, due primarily to the poor intramolecular conformational energies, but also due to the intermolecular description in oxalyl dihydrazide. In all three cases, a fragment-based dispersion-corrected second-order Møller-Plesset perturbation theory (MP2D) treatment of the crystals overcomes these difficulties and predicts conformational polymorph stabilities in good agreement with experiment. These results highlight the need for methods which go beyond current-generation DFT functionals to make crystal polymorph stability predictions truly reliable.

Extension and evaluation of the D4 London-dispersion model for periodic systems

We present an extension of the DFT-D4 model [J. Chem. Phys., 2019, 150, 154122] for periodic systems. The main new ingredients are additional reference polarizabilities for highly-coordinated group 1-5 elements derived from pseudo-periodic electrostatically-embedded cluster calculations. To illustrate the performance of the updated method, several test cases are considered, for which we compare D4 to its predecessor D3(BJ), as well as to a comprehensive set of other dispersion-corrected methods. The largest improvements are observed for solid-state polarizabilities of 16 inorganic salts, where the D4 model achieves an unprecedented accuracy, surpassing its predecessor as well as other, computationally much more demanding approaches. For cell volumes and lattice energies of two sets of chemically diverse molecular crystals, the accuracy gain is less pronounced compared to the already excellently performing D3(BJ) method. For the challenging adsorption energies of small organic molecules on metallic as well as on ionic surfaces, DFT-D4 provides values in good agreement with experimental and/or high-level references. These results suggest the application of the proposed D4 model as a physically improved yet computationally efficient dispersion correction for standard DFT calculations as well as low-cost approaches like semi-empirical or even force-field models.

Assessment of Density Functional Theory in Predicting Interaction Energies between Water and Polycyclic Aromatic Hydrocarbons: from Water on Benzene to Water on Graphene

The interactions of water with polycyclic aromatic hydrocarbons, from benzene to graphene, are investigated using various exchange-correlation functionals selected across the hierarchy of density functional theory (DFT) approximations. The accuracy of the different functionals is assessed through comparisons with random phase approximation (RPA) and coupled-cluster with single, double, and perturbative triple excitations [CCSD(T)] calculations. Diffusion Monte Carlo (DMC) data reported in the literature are also used for comparison. Relatively large variations are found in interaction energies predicted by different DFT models, with GGA functionals underestimating the interaction strength for configurations with the water oxygen pointing toward the aromatic molecules. The meta-GGA B97M-rV and range-separated hybrid, meta-GGA ωB97M-V functionals provide nearly quantitative agreement with CCSD(T) values for the water-benzene, water-coronene, and water-circumcoronene dimers, while RPA and DMC predict interaction energies that differ by up to ∼1 kcal/mol and ∼0.4 kcal/mol from the corresponding CCSD(T) values, respectively. Similar trends among GGA, meta-GGA, and hybrid functionals are observed for larger polycyclic aromatic hydrocarbons. By performing absolutely localized molecular orbital energy decomposition analyses (ALMO-EDA), it is found that, independently of the number of carbon atoms and exchange-correlation functional, the dominant contributions to the interaction energies between water and polycyclic aromatic hydrocarbon molecules are the electrostatic and dispersion terms while polarization and charge transfer effects are negligibly small. Calculations carried out with GGA and meta-GGA functionals indicate that, as the number of carbon atoms increases, the interaction energies slowly converge to the corresponding values obtained for an infinite graphene sheet.

Ab initio thermodynamics of liquid and solid water

Thermodynamic properties of liquid water as well as hexagonal (Ih) and cubic (Ic) ice are predicted based on density functional theory at the hybrid-functional level, rigorously taking into account quantum nuclear motion, anharmonic fluctuations, and proton disorder. This is made possible by combining advanced free-energy methods and state-of-the-art machine-learning techniques. The ab initio description leads to structural properties in excellent agreement with experiments and reliable estimates of the melting points of light and heavy water. We observe that nuclear-quantum effects contribute a crucial [Formula: see text] to the stability of ice Ih, making it more stable than ice Ic. Our computational approach is general and transferable, providing a comprehensive framework for quantitative predictions of ab initio thermodynamic properties using machine-learning potentials as an intermediate step.

Predicting vapor liquid equilibria using density functional theory: A case study of argon

Predicting vapor liquid equilibria (VLE) of molecules governed by weak van der Waals (vdW) interactions using the first principles approach is a significant challenge. Due to the poor scaling of the post Hartree-Fock wave function theory with system size/basis functions, the Kohn-Sham density functional theory (DFT) is preferred for systems with a large number of molecules. However, traditional DFT cannot adequately account for medium to long range correlations which are necessary for modeling vdW interactions. Recent developments in DFT such as dispersion corrected models and nonlocal van der Waals functionals have attempted to address this weakness with a varying degree of success. In this work, we predict the VLE of argon and assess the performance of several density functionals and the second order Møller-Plesset perturbation theory (MP2) by determining critical and structural properties via first principles Monte Carlo simulations. PBE-D3, BLYP-D3, and rVV10 functionals were used to compute vapor liquid coexistence curves, while PBE0-D3, M06-2X-D3, and MP2 were used for computing liquid density at a single state point. The performance of the PBE-D3 functional for VLE is superior to other functionals (BLYP-D3 and rVV10). At T = 85 K and P = 1 bar, MP2 performs well for the density and structural features of the first solvation shell in the liquid phase.

Assessing Density Functionals Using Many Body Theory for Hybrid Perovskites

Which density functional is the "best" for structure simulations of a particular material? A concise, first principles, approach to answer this question is presented. The random phase approximation (RPA)-an accurate many body theory-is used to evaluate various density functionals. To demonstrate and verify the method, we apply it to the hybrid perovskite MAPbI_{3}, a promising new solar cell material. The evaluation is done by first creating finite temperature ensembles for small supercells using RPA molecular dynamics, and then evaluating the variance between the RPA and various approximate density functionals for these ensembles. We find that, contrary to recent suggestions, van der Waals functionals do not improve the description of the material, whereas hybrid functionals and the strongly constrained appropriately normed (SCAN) density functional yield very good agreement with the RPA. Finally, our study shows that in the room temperature tetragonal phase of MAPbI_{3}, the molecules are preferentially parallel to the shorter lattice vectors but reorientation on ps time scales is still possible.

Raman and IR Spectra of Ice Ih and Ice XI with an Assessment of DFT Methods

IR and Raman spectroscopic technology can be directly used to identify the occurrence of ferroelectric ice XI in laboratory or extraterrestrial settings. The performance of 16 different DFT methods applied on the ice Ih, VIII, IX, and XI crystal phases is evaluated. Based on a selected DFT computational scheme, the IR and Raman spectra of ice Ih and XI are derived and compared. When the spectra, both IR and Raman, of ice Ih and ice XI are compared, the librational vibrations are found to be the most affected by the proton ordering. The spectroscopic fingerprint of ice XI can be used to distinguish ferroelectric ice XI from ice Ih in the universe. Furthermore, the existence of only one kind of H-bond in ice Ih is demonstrated from the overlapping subspectra for different types of H-bonded pair configurations in 16 isomers of ice Ih, which provides an illustration to the historic debate on whether one or two kinds of H-bonds existed in ice.

Partially ordered state of ice XV

Most ice polymorphs have order-disorder "pairs" in terms of hydrogen positions, which contributes to the rich variety of ice polymorphs; in fact, three recently discovered polymorphs- ices XIII, XIV, and XV-are ordered counter forms to already identified disordered phases. Despite the considerable effort to understand order-disorder transition in ice crystals, there is an inconsistency among the various experiments and calculations for ice XV, the ordered counter form of ice VI, i.e., neutron diffraction observations suggest antiferroelectrically ordered structures, which disagree with dielectric measurement and theoretical studies, implying ferroelectrically ordered structures. Here we investigate in-situ neutron diffraction measurements and density functional theory calculations to revisit the structure and stability of ice XV. We find that none of the completely ordered configurations are particular favored; instead, partially ordered states are established as a mixture of ordered domains in disordered ice VI. This scenario in which several kinds of ordered configuration coexist dispels the contradictions in previous studies. It means that the order-disorder pairs in ice polymorphs are not one-to-one correspondent pairs but rather have one-to-n correspondence, where there are n possible configurations at finite temperature.

Modeling Polymorphic Molecular Crystals with Electronic Structure Theory

Interest in molecular crystals has grown thanks to their relevance to pharmaceuticals, organic semiconductor materials, foods, and many other applications. Electronic structure methods have become an increasingly important tool for modeling molecular crystals and polymorphism. This article reviews electronic structure techniques used to model molecular crystals, including periodic density functional theory, periodic second-order Møller-Plesset perturbation theory, fragment-based electronic structure methods, and diffusion Monte Carlo. It also discusses the use of these models for predicting a variety of crystal properties that are relevant to the study of polymorphism, including lattice energies, structures, crystal structure prediction, polymorphism, phase diagrams, vibrational spectroscopies, and nuclear magnetic resonance spectroscopy. Finally, tools for analyzing crystal structures and intermolecular interactions are briefly discussed.

Benchmarking DFT and semiempirical methods on structures and lattice energies for ten ice polymorphs

Water in different phases under various external conditions is very important in bio-chemical systems and for material science at surfaces. Density functional theory methods and approximations thereof have to be tested system specifically to benchmark their accuracy regarding computed structures and interaction energies. In this study, we present and test a set of ten ice polymorphs in comparison to experimental data with mass densities ranging from 0.9 to 1.5 g/cm(3) and including explicit corrections for zero-point vibrational and thermal effects. London dispersion inclusive density functionals at the generalized gradient approximation (GGA), meta-GGA, and hybrid level as well as alternative low-cost molecular orbital methods are considered. The widely used functional of Perdew, Burke and Ernzerhof (PBE) systematically overbinds and overall provides inconsistent results. All other tested methods yield reasonable to very good accuracy. BLYP-D3(atm) gives excellent results with mean absolute errors for the lattice energy below 1 kcal/mol (7% relative deviation). The corresponding optimized structures are very accurate with mean absolute relative deviations (MARDs) from the reference unit cell volume below 1%. The impact of Axilrod-Teller-Muto (atm) type three-body dispersion and of non-local Fock exchange is small but on average their inclusion improves the results. While the density functional tight-binding model DFTB3-D3 performs well for low density phases, it does not yield good high density structures. As low-cost alternative for structure related problems, we recommend the recently introduced minimal basis Hartree-Fock method HF-3c with a MARD of about 3%.

Probing the structural and dynamical properties of liquid water with models including non-local electron correlation

Water is a ubiquitous liquid that displays a wide range of anomalous properties and has a delicate structure that challenges experiment and simulation alike. The various intermolecular interactions that play an important role, such as repulsion, polarization, hydrogen bonding, and van der Waals interactions, are often difficult to reproduce faithfully in atomistic models. Here, electronic structure theories including all these interactions at equal footing, which requires the inclusion of non-local electron correlation, are used to describe structure and dynamics of bulk liquid water. Isobaric-isothermal (NpT) ensemble simulations based on the Random Phase Approximation (RPA) yield excellent density (0.994 g/ml) and fair radial distribution functions, while various other density functional approximations produce scattered results (0.8-1.2 g/ml). Molecular dynamics simulation in the microcanonical (NVE) ensemble based on Møller-Plesset perturbation theory (MP2) yields dynamical properties in the condensed phase, namely, the infrared spectrum and diffusion constant. At the MP2 and RPA levels of theory, ice is correctly predicted to float on water, resolving one of the anomalies as resulting from a delicate balance between van der Waals and hydrogen bonding interactions. For several properties, obtaining quantitative agreement with experiment requires correction for nuclear quantum effects (NQEs), highlighting their importance, for structure, dynamics, and electronic properties. A computed NQE shift of 0.6 eV for the band gap and absorption spectrum illustrates the latter. Giving access to both structure and dynamics of condensed phase systems, non-local electron correlation will increasingly be used to study systems where weak interactions are of paramount importance.

Consistent structures and interactions by density functional theory with small atomic orbital basis sets

A density functional theory (DFT) based composite electronic structure approach is proposed to efficiently compute structures and interaction energies in large chemical systems. It is based on the well-known and numerically robust Perdew-Burke-Ernzerhoff (PBE) generalized-gradient-approximation in a modified global hybrid functional with a relatively large amount of non-local Fock-exchange. The orbitals are expanded in Ahlrichs-type valence-double zeta atomic orbital (AO) Gaussian basis sets, which are available for many elements. In order to correct for the basis set superposition error (BSSE) and to account for the important long-range London dispersion effects, our well-established atom-pairwise potentials are used. In the design of the new method, particular attention has been paid to an accurate description of structural parameters in various covalent and non-covalent bonding situations as well as in periodic systems. Together with the recently proposed three-fold corrected (3c) Hartree-Fock method, the new composite scheme (termed PBEh-3c) represents the next member in a hierarchy of "low-cost" electronic structure approaches. They are mainly free of BSSE and account for most interactions in a physically sound and asymptotically correct manner. PBEh-3c yields good results for thermochemical properties in the huge GMTKN30 energy database. Furthermore, the method shows excellent performance for non-covalent interaction energies in small and large complexes. For evaluating its performance on equilibrium structures, a new compilation of standard test sets is suggested. These consist of small (light) molecules, partially flexible, medium-sized organic molecules, molecules comprising heavy main group elements, larger systems with long bonds, 3d-transition metal systems, non-covalently bound complexes (S22 and S66×8 sets), and peptide conformations. For these sets, overall deviations from accurate reference data are smaller than for various other tested DFT methods and reach that of triple-zeta AO basis set second-order perturbation theory (MP2/TZ) level at a tiny fraction of computational effort. Periodic calculations conducted for molecular crystals to test structures (including cell volumes) and sublimation enthalpies indicate very good accuracy competitive to computationally more involved plane-wave based calculations. PBEh-3c can be applied routinely to several hundreds of atoms on a single processor and it is suggested as a robust "high-speed" computational tool in theoretical chemistry and physics.

Singles correlation energy contributions in solids

The random phase approximation to the correlation energy often yields highly accurate results for condensed matter systems. However, ways how to improve its accuracy are being sought and here we explore the relevance of singles contributions for prototypical solid state systems. We set out with a derivation of the random phase approximation using the adiabatic connection and fluctuation dissipation theorem, but contrary to the most commonly used derivation, the density is allowed to vary along the coupling constant integral. This yields results closely paralleling standard perturbation theory. We re-derive the standard singles of Görling-Levy perturbation theory [A. Görling and M. Levy, Phys. Rev. A 50, 196 (1994)], highlight the analogy of our expression to the renormalized singles introduced by Ren and coworkers [Phys. Rev. Lett. 106, 153003 (2011)], and introduce a new approximation for the singles using the density matrix in the random phase approximation. We discuss the physical relevance and importance of singles alongside illustrative examples of simple weakly bonded systems, including rare gas solids (Ne, Ar, Xe), ice, adsorption of water on NaCl, and solid benzene. The effect of singles on covalently and metallically bonded systems is also discussed.

Construction and application of a new dual-hybrid random phase approximation

The direct random phase approximation (dRPA) combined with Kohn-Sham reference orbitals is among the most promising tools in computational chemistry and applicable in many areas of chemistry and physics. The reason for this is that it scales as N(4) with the system size, which is a considerable advantage over the accurate ab initio wave function methods like standard coupled-cluster. dRPA also yields a considerably more accurate description of thermodynamic and electronic properties than standard density-functional theory methods. It is also able to describe strong static electron correlation effects even in large systems with a small or vanishing band gap missed by common single-reference methods. However, dRPA has several flaws due to its self-correlation error. In order to obtain accurate and precise reaction energies, barriers and noncovalent intra- and intermolecular interactions, we construct a new dual-hybrid dRPA (hybridization of exact and semilocal exchange in both the energy and the orbitals) and test the performance of this new functional on isogyric, isodesmic, hypohomodesmotic, homodesmotic, and hyperhomodesmotic reaction classes. We also use a test set of 14 Diels-Alder reactions, six atomization energies (AE6), 38 hydrocarbon atomization energies, and 100 reaction barrier heights (DBH24, HT-BH38, and NHT-BH38). For noncovalent complexes, we use the NCCE31 and S22 test sets. To test the intramolecular interactions, we use a set of alkane, cysteine, phenylalanine-glycine-glycine tripeptide, and monosaccharide conformers. We also discuss the delocalization and static correlation errors. We show that a universally accurate description of chemical properties can be provided by a large, 75% exact exchange mixing both in the calculation of the reference orbitals and the final energy.

The individual and collective effects of exact exchange and dispersion interactions on the ab initio structure of liquid water

In this work, we report the results of a series of density functional theory (DFT) based ab initio molecular dynamics (AIMD) simulations of ambient liquid water using a hierarchy of exchange-correlation (XC) functionals to investigate the individual and collective effects of exact exchange (Exx), via the PBE0 hybrid functional, non-local van der Waals/dispersion (vdW) interactions, via a fully self-consistent density-dependent dispersion correction, and an approximate treatment of nuclear quantum effects, via a 30 K increase in the simulation temperature, on the microscopic structure of liquid water. Based on these AIMD simulations, we found that the collective inclusion of Exx and vdW as resulting from a large-scale AIMD simulation of (H2O)128 significantly softens the structure of ambient liquid water and yields an oxygen-oxygen structure factor, SOO(Q), and corresponding oxygen-oxygen radial distribution function, gOO(r), that are now in quantitative agreement with the best available experimental data. This level of agreement between simulation and experiment demonstrated herein originates from an increase in the relative population of water molecules in the interstitial region between the first and second coordination shells, a collective reorganization in the liquid phase which is facilitated by a weakening of the hydrogen bond strength by the use of a hybrid XC functional, coupled with a relative stabilization of the resultant disordered liquid water configurations by the inclusion of non-local vdW/dispersion interactions. This increasingly more accurate description of the underlying hydrogen bond network in liquid water also yields higher-order correlation functions, such as the oxygen-oxygen-oxygen triplet angular distribution, POOO(θ), and therefore the degree of local tetrahedrality, as well as electrostatic properties, such as the effective molecular dipole moment, that are in much better agreement with experiment.

The accurate calculation of the band gap of liquid water by means of GW corrections applied to plane-wave density functional theory molecular dynamics simulations

Knowledge about the intrinsic electronic properties of water is imperative for understanding the behaviour of aqueous solutions that are used throughout biology, chemistry, physics, and industry.

Periodic MP2, RPA, and Boundary Condition Assessment of Hydrogen Ordering in Ice XV

Ice XV is the hydrogen-ordered form of the ice VI phase whose structure was predicted to be Cc and ferroelectric using periodic DFT approaches. However, neutron diffraction and Raman spectroscopy data show the structure to have P1̅ symmetry and to be antiferroelectric. Recent work1 using fragment-based MP2 and CCSD(T) approaches predicts the experimental structure as the ground state. We have reconsidered this problem using fully periodic MP2 and RPA approaches and find that the ferroelectric Cc structure is the lowest energy configuration. However, ubiquitously employed tinfoil boundary conditions stabilize polar structures. We suggest that ferroelectric Cc crystals can grow within a paraelectric ice VI matrix but may become unstable once a fraction of the matrix has become hydrogen-ordered. The reduction in dielectric constant causes P1̅ and other structures with small polarization to become favored, providing a possible resolution between observation and theoretical predictions.

Proton Ordering of Cubic Ice Ic: Spectroscopy and Computer Simulations

Several proton-disordered crystalline ice structures are known to proton order at sufficiently low temperatures, provided that the right preparation procedure is used. For cubic ice, ice Ic, however, no proton ordering has been observed so far. Here, we subject ice Ic to an experimental protocol similar to that used to proton order hexagonal ice. In situ FT-IR spectroscopy carried out during this procedure reveals that the librational band of the spectrum narrows and acquires a structure that is observed neither in proton-disordered ice Ic nor in ice XI, the proton-ordered variant of hexagonal ice. On the basis of vibrational spectra computed for ice Ic and four of its proton-ordered variants using classical molecular dynamics and ab initio simulations, we conclude that the features of our experimental spectra are due to partial proton ordering, providing the first evidence of proton ordering in cubic ice. We further find that the proton-ordered structure with the lowest energy is ferroelectric, while the structure with the second lowest energy is weakly ferroelectric. Both structures fit the experimental spectral similarly well such that no unique assignment of proton order is possible based on our results.