Accurate Potential Energy Surfaces Using Atom-Centered Potentials and Minimal High-Level Data.

We demonstrate that a Δ-density functional theory (Δ-DFT) approach based on atom-centered potentials (ACPs) represents a computationally inexpensive and accurate method for representing potential energy surfaces (PESs) for the HONO and HFCO molecules and vibrational frequencies derived therefrom. Using as few as 100 CCSD(T)-F12a reference energies, ACPs developed for use with B3LYP/def2-TZVPP are shown to produce PESs for HONO and HFCO with mean absolute errors of 27.7 and 5.8 cm-1, respectively. Application of the multiconfigurational time-dependent Hartree (MCTDH) method with ACP-corrected B3LYP/def2-TZVPP PESs produces vibrational frequencies for cis- and trans-HONO with mean absolute percent errors (MAPEs) of 0.8 and 1.1, compared to 0.8 obtained for the two isomers with CCSD(T)-F12a/cc-pVTZ-F12/MCTDH. For HFCO, the vibrational frequencies obtained using the present (Δ-DFT)/MCTDH approach give a MAPE of 0.1, which is the error obtained with CCSD(T)-F12a/cc-pVTZ-F12/MCTDH. The ACP approach is therefore successful in representing a PES calculated at a high level of theory (CCSD(T)-F12a) and a promising method for the development of a general protocol for the representation of accurate molecular PESs and the calculation of molecular properties from them.

Minimal Basis Set Hartree-Fock Corrected with Atom-Centered Potentials for Molecular Crystal Modeling and Crystal Structure Prediction.

Crystal structure prediction (CSP), determining the experimentally observable structure of a molecular crystal from the molecular diagram, is an important challenge with technologically relevant applications in materials manufacturing and drug design. For the purpose of screening the randomly generated candidate crystal structures, CSP protocols require energy ranking methods that are fast and can accurately capture the small energy differences between molecular crystals. In addition, a good ranking method should also produce accurate equilibrium geometries, both intramolecular and intermolecular. In this article, we explore the combination of minimal-basis-set Hartree-Fock (HF) with atom-centered potentials (ACPs) as a method for modeling the structure and energetics of molecular crystals. The ACPs are developed for the H, C, N, and O atoms and fitted to a set of reference data at the B86bPBE-XDM level in order to mitigate basis-set incompleteness and missing correlation. In particular, ACPs are developed in combination with two methods: HF-D3/MINIs and HF-3c. The application of ACPs greatly improves the performance of HF-D3/MINIs for lattice energies, crystal energy differences, energy-volume and energy-strain relations, and crystal geometries. In the case of HF-3c, the improvement in the crystal energy differences is much smaller than in HF-D3/MINIs, but lattice energies and particularly crystal geometries are considerably better when ACPs are used. The resulting methods may be useful for CSP but also for quick calculation of molecular crystal lattice energies and geometries.

Finding critical points and reconstruction of electron densities on grids.

The quantum theory of atoms in molecules (QTAIM), developed by Bader and co-workers, is one of the most popular ways of extracting chemical insight from the results of quantum mechanical calculations. One of the basic tasks in QTAIM is to locate the critical points of the electron density and calculate various quantities (density, Laplacian, etc.) on them since these have been found to correlate with molecular properties of interest. If the electron density is given analytically, this process is relatively straightforward. However, locating the critical points is more challenging if the density is known only on a three-dimensional uniform grid. A density grid is common in periodic solids because it is the natural expression for the electron density in plane-wave calculations. In this article, we explore the reconstruction of the electron density from a grid and its use in critical point localization. The proposed reconstruction method employs polyharmonic spline interpolation combined with a smoothing function based on the promolecular density. The critical point search based on this reconstruction is accurate, trivially parallelizable, works for periodic and non-periodic systems, does not present directional lattice bias when the grid is non-orthogonal, and locates all critical points of the underlying electron density in all tests studied. The proposed method also provides an accurate reconstruction of the electron density over the space spanned by the grid, which may be useful in other contexts besides critical point localization.

A density-functional benchmark of vibrational free-energy corrections for molecular crystal polymorphism.

Many crystal structure prediction protocols only concern themselves with the electronic energy of molecular crystals. However, vibrational contributions to the free energy (Fvib) can be significant in determining accurate stability rankings for crystal candidates. While force-field studies have been conducted to gauge the magnitude of these free-energy corrections, highly accurate results from quantum mechanical methods, such as density-functional theory (DFT), are desirable. Here, we introduce the PV17 set of 17 polymorphic pairs of organic molecular crystals, for which plane wave DFT is used to calculate the vibrational free energies and free-energy differences (ΔFvib) between each pair. Our DFT results confirm that the vibrational free-energy corrections are small, having a mean value of 1.0 kJ/mol and a maximum value of 2.3 kJ/mol for the PV17 set. Furthermore, we assess the accuracy of a series of lower-cost DFT, semi-empirical, and force-field models for computing ΔFvib that have been proposed in the literature. It is found that calculating Fvib using the Γ-point frequencies does not provide ΔFvib values of sufficiently high quality. In addition, ΔFvib values calculated using various approximate methods have mean absolute errors relative to our converged DFT results of equivalent or larger magnitude than the vibrational free-energy corrections themselves. Thus, we conclude that, in a crystal structure prediction protocol, it is preferable to forego the inclusion of vibrational free-energy corrections than to estimate them with any of the approximate methods considered here.

Experimental and Theoretical Study of SbPO4 under Compression.

SbPO4 is a complex monoclinic layered material characterized by a strong activity of the nonbonding lone electron pair (LEP) of Sb. The strong cation LEP leads to the formation of layers piled up along the a axis and linked by weak Sb-O electrostatic interactions. In fact, Sb has 4-fold coordination with O similarly to what occurs with the P-O coordination, despite the large difference in ionic radii and electronegativity between both elements. Here we report a joint experimental and theoretical study of the structural and vibrational properties of SbPO4 at high pressure. We show that SbPO4 is not only one of the most compressible phosphates but also one of the most compressible compounds of the ABO4 family. Moreover, it has a considerable anisotropic compression behavior, with the largest compression occurring along a direction close to the a axis and governed by the compression of the LEP and the weak interlayer Sb-O bonds. The strong compression along the a axis leads to a subtle modification of the monoclinic crystal structure above 3 GPa, leading from a 2D to a 3D material. Moreover, the onset of a reversible pressure-induced phase transition is observed above 9 GPa, which is completed above 20 GPa. We propose that the high-pressure phase is a triclinic distortion of the original monoclinic phase. The understanding of the compression mechanism of SbPO4 can aid to improve the ion intercalation and catalytic properties of this layered compound.

Interplay between local structure, vibrational and electronic properties on CuO under pressure.

The electronic and local structural properties of CuO under pressure have been investigated by means of X-ray absorption spectroscopy (XAS) at Cu K edge and ab initio calculations, up to 17 GPa. The crystal structure of CuO consists of Cu motifs within CuO4 square planar units and two elongated apical Cu-O bonds. The CuO4 square planar units are stable in the studied pressure range, with Cu-O distances that are approximately constant up to 5 GPa, and then decrease slightly up to 17 GPa. In contrast, the elongated Cu-O apical distances decrease continuously with pressure in the studied range. An anomalous increase of the mean square relative displacement (EXAFS Debye-Waller, σ2) of the elongated Cu-O path is observed from 5 GPa up to 13 GPa, when a drastic reduction takes place in σ2. This is interpreted in terms of local dynamic disorder along the apical Cu-O path. At higher pressures (P > 13 GPa), the local structure of Cu2+ changes from a 4-fold square planar to a 4+2 Jahn-Teller distorted octahedral ion. We interpret these results in terms of the tendency of the Cu2+ ion to form favorable interactions with the apical O atoms. Also, the decrease in Cu-O apical distance caused by compression softens the normal mode associated with the out-of-plane Cu movement. CuO is predicted to have an anomalous rise in permittivity with pressure as well as modest piezoelectricity in the 5-13 GPa pressure range. In addition, the near edge features in our XAS experiment show a discontinuity and a change of tendency at 5 GPa. For P < 5 GPa the evolution of the edge shoulder is ascribed to purely electronic effects which also affect the charge transfer integral. This is linked to a charge migration from the Cu to O, but also to an increase of the energy band gap, which show a change of tendency occurring also at 5 GPa.

Origin of Nanoscale Friction Contrast between Supported Graphene, MoS2, and a Graphene/MoS2 Heterostructure.

Ultralow friction can be achieved with 2D materials, particularly graphene and MoS2. The nanotribological properties of these different 2D materials have been measured in previous atomic force microscope (AFM) experiments sequentially, precluding immediate and direct comparison of their frictional behavior. Here, friction is characterized at the nanoscale using AFM experiments with the same tip sliding over graphene, MoS2, and a graphene/MoS2 heterostructure in a single measurement, repeated hundreds of times, and also measured with a slowly varying normal force. The same material systems are simulated using molecular dynamics (MD) and analyzed using density functional theory (DFT) calculations. In both experiments and MD simulations, graphene consistently exhibits lower friction than the MoS2 monolayer and the heterostructure. In some cases, friction on the heterostructure is lower than that on the MoS2 monolayer. Quasi-static MD simulations and DFT calculations show that the origin of the friction contrast is the difference in energy barriers for a tip sliding across each of the three surfaces.

Oscillatory motion in layered materials: graphene, boron nitride, and molybdenum disulfide.

Offset-driven self-retraction and oscillatory motion of bilayer graphene has been observed experimentally and is potentially relevant for nanoscale technological applications. In a previous article, we showed that friction between laterally offset graphene layers is controlled by roughness and proposed a simple reduced-order model based on density-functional theory (DFT) and molecular dynamics (MD) data, with which predictions on the experimental size-scale could be made. In this article, we extend our study to other layered materials, with emphasis on boron nitride (BN) and molybdenum disulfide (MoS2). Using MD and DFT simulations of these systems and a generalized version of the reduced-order model, we predict that BN will exhibit behavior similar to graphene (heavily-damped oscillation with a decay rate that increases with roughness) and that MoS2 shows no oscillatory behavior even in the absence of roughness. This is attributed to the higher energy barrier for sliding in MoS2 as well as the surface structure. Our generalized reduced-order model provides a guide to predicting and tuning experimental oscillation behavior using a few parameters that can be derived from simulation data.

The role of roughness-induced damping in the oscillatory motion of bilayer graphene.

A multi-scale theoretical model is presented that is the first to offer quantitative agreement with experimental measurements of self-retraction and oscillation of bilayer graphene. The model integrates density-functional theory calculations of the energetics driving flake retraction and molecular-dynamics simulations capturing the dynamic response of laterally-offset rough surfaces. We demonstrate that nanoscale roughness explains self-retraction motion and propose a recipe for tuning that motion by controlling friction.

Density-functional description of electrides.

Electrides are a unique class of ionic solids in which the anions are stoichiometrically replaced by electrons localised within the crystal voids. There are only nine electrides with known crystal structures and their study represents a challenge for theory. A systematic investigation of their electronic structure is conducted using semilocal density-functional theory (DFT) in this article. The band structure is calculated for each electride and a range of techniques including plots of the valence densities, procrystal densities, non-covalent interaction isosurfaces, and Bader's quantum-chemical topology are applied. All of these methods provide consistent results, confirming the presence of localised interstitial electrons and demonstrating that these crystals display a characteristic electronic structure.

Metallophilic interactions from dispersion-corrected density-functional theory.

In this article, we present the first comprehensive study of metallophilic (aurophilic) interactions using dispersion-corrected density-functional theory. Dispersion interactions (an essential component of metallophilicity) are treated using the exchange-hole dipole moment (XDM) model. By comparing against coupled-cluster benchmark calculations on simple dimers, we show that LC-ωPBE-XDM is a viable functional to study interactions between closed-shell transition metals and that it performs uniformly better than second-order Møller-Plesset theory, the basic computational technique used in previous works. We apply LC-ωPBE-XDM to address several open questions regarding metallophilicity, such as the interplay between dispersion and relativistic effects, the interaction strength along group 11, the additivity of homo- and hetero-metallophilic effects, the stability of [E(AuPH3)4](+) cations (E = N, P, As, Sb), and the role of metallophilic effects in crystal packing. We find that relativistic effects explain the prevalence of aurophilicity not by stabilizing metal-metal contacts, but by preventing gold from forming ionic structures involving bridge anions (which are otherwise common for Ag and Cu) as a result of the increased electron affinity of the metal. Dispersion effects are less important than previously assumed and their stabilization contribution is relatively independent of the metal.

Predicting the relative solubilities of racemic and enantiopure crystals by density-functional theory.

Isolation of chiral molecules as pure enantiomers remains a fundamental challenge in chemical research. Enantioselective enrichment through preferential crystallization is an efficient method to achieve enantiopure compounds, but its applicability depends on the relative stability of the enantiopure and racemic crystal forms. Using a simple thermodynamic model and first-principles density-functional calculations, it is possible to predict the difference in solubility between the enantiopure and racemic solid phases. This approach uses dispersion-corrected density functionals and is capable of accurately predicting the solution-phase entantiomeric excess to within about 10 % of experimental measurements on average. The accuracy of the exchange-hole dipole moment (XDM) model of dispersion enables the viability of the proposed method.

Polymorphism and Phase Stability of Hydrated Magnesium Carbonate Nesquehonite MgCO3·3H2O: Negative Axial Compressibility and Thermal Expansion in a Cementitious Material.

The P-T phase diagram of the hydrated magnesium carbonate nesquehonite (MgCO3·3H2O) has not been reported in the literature. In this paper, we present a joint experimental and computational study of the phase stability and structural behavior of this cementitious material at high-pressure and high-temperature conditions using in situ single-crystal and synchrotron powder X-ray diffraction measurements in resistive-heated diamond anvil cells plus density functional theory calculations. Our results show that nesquehonite undergoes two pressure-induced phase transitions at 2.4 (HP1) and 4.0 GPa (HP2) at ambient temperature. We have found negative axial compressibility and thermal expansivity values, likely related to the directionality of the hydrogen bonds. The equations of state of the different phases have been determined. All the room-temperature compression effects were reversible. Heating experiments at 0.7 GPa show a first temperature-induced decomposition at 115 °C, probably into magnesite and a MgCO3·4H2O phase.

Performance of conventional and dispersion-corrected density-functional theory methods for hydrogen bonding interaction energies.

The approximate CCSD(T)/CBS binding energies for the set of 23 hydrogen-bonded dimers (HB23) of the S66 set reported by Rezác et al. (J. Chem. Theory Comput. 2011, 7, 2427-2438) were expected to be under-estimated based on the known under-binding tendency of the counterpoise correction combined with small basis sets. In this work, we present binding energies for the HB23 set of dimers obtained using a composite approach recently described by Mackie and DiLabio (J. Chem. Phys. 2011, 135, 134318) that averages the counterpoise- and non-counterpoise-corrected energies, while utilizing standard approaches to obtain CCSD(T)/CBS-type energies. The binding energies for the HB23 set are revised upward by an average of 0.12 kcal mol(-1) and as much as 0.35 kcal mol(-1). We use these improved benchmark-level binding energies to evaluate the ability of pure, hybrid, long-range-corrected, and dispersion-corrected density-functional theory (DFT) methods to accurately predict the binding energies of hydrogen-bonded dimers. We find that, in general, the inclusion of dispersion into the DFT approach is required in order to obtain reasonable results for the HB23 set. We find that the dispersion-corrected DFT methods we tested produce results of variable quality, as measured by mean absolute deviation relative to the revised reference values we computed: B97D, 0.57 kcal mol(-1); B3LYP-D3, 0.44 kcal mol(-1); ωB97XD, 0.25 kcal mol(-1); LC-ωPBE-D3, 0.24 kcal mol(-1); M06-2X, 0.21 kcal mol(-1); B3LYP-DCP, 0.23 kcal mol(-1); B971-XDM, 0.18 kcal mol(-1).

Efficient basis sets for non-covalent interactions in XDM-corrected density-functional theory.

In the development and application of dispersion-corrected density-functional theory, the effects of basis set incompleteness have been largely mitigated through the use of very large, nearly-complete basis sets. However, the use of such large basis sets makes application of these methods inefficient for large systems. In this work, we examine a series of basis sets, including Pople-style, correlation-consistent, and polarization-consistent bases, for their ability to efficiently and accurately predict non-covalent interactions when used in conjunction with the exchange-hole dipole moment (XDM) dispersion model. We find that the polarization-consistent 2 (pc-2) basis sets, and two modifications thereof with some diffuse functions removed, give performance of comparable quality to that obtained with aug-cc-pVTZ basis sets, while being roughly 12 to 23 times faster computationally. The behavior is explained, in part, by the role of diffuse functions in recovering small density changes in the intermolecular region. The general performance of the modified basis sets is tested by application of XDM to standard intermolecular benchmark sets at, and away from, equilibrium.

XDM-corrected hybrid DFT with numerical atomic orbitals predicts molecular crystal lattice energies with unprecedented accuracy.

Molecular crystals are important for many applications, including energetic materials, organic semiconductors, and the development and commercialization of pharmaceuticals. The exchange-hole dipole moment (XDM) dispersion model has shown good performance in the calculation of relative and absolute lattice energies of molecular crystals, although it has traditionally been applied in combination with plane-wave/pseudopotential approaches. This has limited XDM to use with semilocal functional approximations, which suffer from delocalization error and poor quality conformational energies, and to systems with a few hundreds of atoms at most due to unfavorable scaling. In this work, we combine XDM with numerical atomic orbitals, which enable the efficient use of XDM-corrected hybrid functionals for molecular crystals. We test the new XDM-corrected functionals for their ability to predict the lattice energies of molecular crystals for the X23 set and 13 ice phases, the latter being a particularly stringent test. A composite approach using a XDM-corrected, 25% hybrid functional based on B86bPBE achieves a mean absolute error of 0.48 kcal mol-1 per molecule for the X23 set and 0.19 kcal mol-1 for the total lattice energies of the ice phases, compared to recent diffusion Monte-Carlo data. These results make the new XDM-corrected hybrids not only far more computationally efficient than previous XDM implementations, but also the most accurate density-functional methods for molecular crystal lattice energies to date.

Structural Behavior of Minrecordite Carbonate Mineral upon Compression: Effect of Mg → Zn Chemical Substitution in Dolomite-Type Compounds.

We report the structural behavior and compressibility of minrecordite, a naturally occurring Zn-rich dolomite mineral, determined using diamond-anvil cell synchrotron X-ray diffraction. Our data show that this rhombohedral CaZn0.52Mg0.48(CO3)2 carbonate exhibits a highly anisotropic behavior, the c axis being 3.3 times more compressible than the a axis. The axial compressibilities and the equation of state are governed by the compression of the [CaO6] and [ZnO6] octahedra, which are the cations in larger proportion in each layer. We observe the existence of a dense polymorph above 13.4(3) GPa using Ne as a pressure-transmitting medium, but the onset pressure of the phase transition decreases with the appearance of deviatoric stresses in nonhydrostatic conditions. Our results suggest that the phase transition observed in minrecordite is strain-induced and that the high-pressure polymorph is intimately related to the CaCO3-II-type structure. A comparison with other dolomite minerals indicates that the transition pressure decreases when the ratio Zn/Mg in the crystal lattice of pure dolomite is larger than 1. Density functional theory (DFT) calculations predict that a distorted CaCO3-II-type structure is energetically more stable than dolomite-type CaZn(CO3)2 above 10 GPa. However, according to our calculations, the most stable structure above this pressure is a dolomite-V-type phase, a polymorph not observed experimentally.

Revealing non-covalent interactions in solids: NCI plots revisited.

In this article, the NCI method [Johnson et al., J. Am. Chem. Soc., 2010, 132, 6498] for plotting and analysing non-covalent interactions (NCI) is extended to periodic (solid-state) electron densities and implemented in the critic program. The new code uses self-consistent electron densities from a variety of electronic structure methods (pseudopotentials/plane-wave, FP-LAPW, local orbitals, etc.), and it can also build the promolecular density from the crystal geometry alone. As an example of the new code, intermolecular interactions in several molecular crystals are presented and analyzed. The connection with QTAIM studies is established and a reinterpretation of the NCI domains is given regarding the current knowledge of the field. The connection between NCI domains and intermolecular vibrations is made apparent, as well as the ability of the method to reveal the locality of bonding.

Fast and Accurate Quantum Mechanical Modeling of Large Molecular Systems Using Small Basis Set Hartree-Fock Methods Corrected with Atom-Centered Potentials.

There has been significant interest in developing fast and accurate quantum mechanical methods for modeling large molecular systems. In this work, by utilizing a machine learning regression technique, we have developed new low-cost quantum mechanical approaches to model large molecular systems. The developed approaches rely on using one-electron Gaussian-type functions called atom-centered potentials (ACPs) to correct for the basis set incompleteness and the lack of correlation effects in the underlying minimal or small basis set Hartree-Fock (HF) methods. In particular, ACPs are proposed for ten elements common in organic and bioorganic chemistry (H, B, C, N, O, F, Si, P, S, and Cl) and four different base methods: two minimal basis sets (MINIs and MINIX) plus a double-ζ basis set (6-31G*) in combination with dispersion-corrected HF (HF-D3/MINIs, HF-D3/MINIX, HF-D3/6-31G*) and the HF-3c method. The new ACPs are trained on a very large set (73 832 data points) of noncovalent properties (interaction and conformational energies) and validated additionally on a set of 32 048 data points. All reference data are of complete basis set coupled-cluster quality, mostly CCSD(T)/CBS. The proposed ACP-corrected methods are shown to give errors in the tenths of a kcal/mol range for noncovalent interaction energies and up to 2 kcal/mol for molecular conformational energies. More importantly, the average errors are similar in the training and validation sets, confirming the robustness and applicability of these methods outside the boundaries of the training set. In addition, the performance of the new ACP-corrected methods is similar to complete basis set density functional theory (DFT) but at a cost that is orders of magnitude lower, and the proposed ACPs can be used in any computational chemistry program that supports effective-core potentials without modification. It is also shown that ACPs improve the description of covalent and noncovalent bond geometries of the underlying methods and that the improvement brought about by the application of the ACPs is directly related to the number of atoms to which they are applied, allowing the treatment of systems containing some atoms for which ACPs are not available. Overall, the ACP-corrected methods proposed in this work constitute an alternative accurate, economical, and reliable quantum mechanical approach to describe the geometries, interaction energies, and conformational energies of systems with hundreds to thousands of atoms.

Small-Basis Set Density-Functional Theory Methods Corrected with Atom-Centered Potentials.

Density functional theory (DFT) is currently the most popular method for modeling noncovalent interactions and thermochemistry. The accurate calculation of noncovalent interaction energies, reaction energies, and barrier heights requires choosing an appropriate functional and, typically, a relatively large basis set. Deficiencies of the density-functional approximation and the use of a limited basis set are the leading sources of error in the calculation of noncovalent and thermochemical properties in molecular systems. In this article, we present three new DFT methods based on the BLYP, M06-2X, and CAM-B3LYP functionals in combination with the 6-31G* basis set and corrected with atom-centered potentials (ACPs). ACPs are one-electron potentials that have the same form as effective-core potentials, except they do not replace any electrons. The ACPs developed in this work are used to generate energy corrections to the underlying DFT/basis-set method such that the errors in predicted chemical properties are minimized while maintaining the low computational cost of the parent methods. ACPs were developed for the elements H, B, C, N, O, F, Si, P, S, and Cl. The ACP parameters were determined using an extensive training set of 118655 data points, mostly of complete basis set coupled-cluster level quality. The target molecular properties for the ACP-corrected methods include noncovalent interaction energies, molecular conformational energies, reaction energies, barrier heights, and bond separation energies. The ACPs were tested first on the training set and then on a validation set of 42567 additional data points. We show that the ACP-corrected methods can predict the target molecular properties with accuracy close to complete basis set wavefunction theory methods, but at a computational cost of double-ζ DFT methods. This makes the new BLYP/6-31G*-ACP, M06-2X/6-31G*-ACP, and CAM-B3LYP/6-31G*-ACP methods uniquely suited to the calculation of noncovalent, thermochemical, and kinetic properties in large molecular systems.