The structure, equation of state, IR, Raman, and inelastic neutron scattering (INS) spectra of high-pressure, proton-ordered phase VIII of ice are calculated by the second-order many-body perturbation and coupled-cluster singles and doubles methods. Nearly all the observed features of the pressure-dependence of the structures and spectra are reproduced computationally up to 60 GPa insofar as the anharmonic effects can be neglected. The calculations display no sign of the hypothetical isostructural transition in 2–3 GPa to phase VIII′, the existence of which has been a matter of controversy for over a decade, while they do not contradict the interpretation of the spectral anomaly at 10–14 GPa as a precursor of the VIII-X phase transition. The calculated INS spectra correct a systematic error in the peak positions of the observed spectra

Molecular crystals are chemists’ solids in the sense that their structures and properties can be understood in terms of those of the constituent molecules merely perturbed by a crystalline environment. They form a large and important class of solids including ices of atmospheric species, drugs, explosives, and even some organic optoelectronic materials and supramolecular assemblies. Recently, surprisingly simple yet extremely efficient, versatile, easily implemented, and systematically accurate electronic structure methods for molecular crystals have been developed. The methods, collectively referred to as the embedded-fragment scheme, divide a crystal into monomers and overlapping dimers and apply modern molecular electronic structure methods and software to these fragments of the crystal that are embedded in a self-consistently determined crystalline electrostatic field. They enable facile applications of accurate but otherwise prohibitively expensive ab initio molecular orbital theories such as Møller–Plesset perturbation and coupled-cluster theories to a broad range of properties of solids such as internal energies, enthalpies, structures, equation of state, phonon dispersion curves and density of states, infrared and Raman spectra (including band intensities and sometimes anharmonic effects), inelastic neutron scattering spectra, heat capacities, Gibbs energies, and phase diagrams, while accounting for many-body electrostatic (namely, induction or polarization) effects as well as two-body exchange and dispersion interactions from first principles. They can fundamentally alter the role of computing in the studies of molecular crystals in the same way ab initio molecular orbital theories have transformed research practices in gas-phase physical chemistry and synthetic chemistry in the last half century.

In this Account, after a brief summary of formalisms and algorithms, we discuss applications of these methods performed in our group as compelling illustrations of their unprecedented power in addressing some of the outstanding problems of solid-state chemistry, high-pressure chemistry, or geochemistry. They are the structure and spectra of ice Ih, in particular, the origin of two peaks in the hydrogen-bond-stretching region of its inelastic neutron scattering spectra, a solid–solid phase transition from CO2-I to elusive, metastable CO2-III, pressure tuning of Fermi resonance in solid CO2, and the structure and spectra of solid formic acid, all at the level of second-order Møller–Plesset perturbation theory or higher.

Despite its terrestrial abundance and astrochemical significance, many aspects of the phase diagram of solid carbon dioxide remain uncertain or unknown. The observed transition pressures from cubic to orthorhombic phase range widely from 2.5 GPa at 80 K to above 18 GPa at room temperature. The vibrational Raman bands that appear at higher pressure and serve as a decisive proof of the existence of the orthorhombic phase have never been assigned. Here we introduce a general ab initio computational method that can predict the Gibbs free energies and thus phase diagrams of molecular crystals. Using this with second-order Møller—Plesset perturbation theory, we obtain the transition pressure of 13 GPa at 0 K with small temperature dependence, which is in line with many experiments. We also computationally reproduce the vibrational Raman bands and explain the pressure dependence of the structure parameters and Raman band positions of both phases quantitatively.

The symmetric-stretching fundamental (ν1) and the bending first overtone (2ν2) of CO2, which are accidentally degenerate with the same symmetry, undergo a Fermi resonance and give rise to two Raman bands with a frequency difference of 107 cm−1 and an intensity ratio of 2.1. Both the frequency difference and intensity ratio can be varied by pressure applied to CO2 in condensed phases, which has been utilized as a spectroscopic geobarometer for minerals with CO2 inclusion. This study calculates the pressure dependence of the Fermi dyad frequency difference and intensity ratio by combining the embedded-fragment second-order Møller–Plesset perturbation calculations of harmonic frequencies of solid CO2 under pressure and the coupled-cluster singles and doubles with noniterative triples and vibrational configuration-interaction calculations of anharmonic frequencies of molecular CO2. It reproduces frequency difference quantitatively and intensity ratio qualitatively up to 10 GPa. The analysis of the results is shown to render strong support for one particular order of unperturbed frequencies, ν1 > 2ν2, in both the gas and solid phases, which has been a matter of controversy for decades.

Ice Ih is arguably the most important molecular crystal in nature, yet our understanding of its structural and dynamical properties is still far from complete. We present embedded-fragment calculations of the structures and vibrational spectra of the three-dimensional, proton-disordered phase of ice Ih performed at the level of second-order many-body perturbation theory with a basis-set superposition error correction. Our calculations address previous controversies such as the one related to the O–H bond length as well as the existence of two types of hydrogen bonds with strengths differing by a factor of two. For the latter, our calculations suggest that the observed spectral features arise from the directionality or the anisotropy of collective hydrogen-bond stretching vibrations rather than the previously suggested vastly different force constants. We also report a capability to efficiently compute infrared and Raman intensities of a periodic solid. Our approach reproduces the infrared and Raman spectra, the variation of inelastic neutron scattering spectra with deuterium concentration, and the anomaly of heat capacities at low temperatures for ice Ih.

Can the zero-point vibrational energies (ZPVE) of molecular clusters and crystals be evaluated as sums of ZPVE of constituent molecular fragments embedded in the cluster or crystal electrostatic environment? What is the appropriate unit of fragmentation: monomers or overlapping dimers? Can the contributions of acoustic phonons, which are fundamentally delocalized, be recuperated at satisfactory accuracy? These questions are answered by this study applying embedded monomer- and dimer-fragmentation methods to the harmonic ZPVE of hydrogen fluoride clusters, hydrogen fluoride crystal, and water clusters. Our findings are as follows: (1) ZPVE are reproduced accurately by both fragmentation schemes within a few percents of exact values or a few tenths of 1 kcal mol−1 per molecule even for crystalline hydrogen fluoride, which has acoustic phonons. (2) Both the monomer- and dimer-based fragmentation are nearly equally accurate and useful for the absolute values of ZPVE, but the latter is more reliable than the former in reproducing the relative ZPVE of cluster isomers of the same size. (3) The embedding field is essential as it renders nonzero frequencies to the translational and rotational motions of monomers and dimers, accounting for the pseudo-translational and librational motions of the entire clusters or crystals. (4) Some of these low-frequency modes of fragments are calculated to have imaginary frequencies because the fragments are not at their equilibrium geometries, causing ZPVE to be complex. The imaginary part of ZPVE, which is nonphysical and is guaranteed to vanish in the exact limit of the many-body expansion, is nonetheless a useful estimate of errors in the real part.

A linear-scaling, embedded-fragment, second-order many-body perturbation (MP2) method with basis sets up to aug-cc-pVTZ is applied to the antiparallel structure of solid hydrogen fluoride and deuterium fluoride under 0–20 GPa of ambient pressure. The optimized structures, including the lattice parameters and molar volume, and phonon dispersion as well as phonon density of states (DOS), are determined as a function of pressure. The basis-set superposition errors are removed by the counterpoise correction. The structural parameters at 0 GPa calculated by MP2 agree accurately with the observed, making the predicted values at higher pressures a useful pilot for future experiments. The corresponding values obtained by the Hartree–Fock method have large, systematic errors. The MP2/aug-cc-pVDZ frequencies of the infrared- and Raman-active vibrations of the three-dimensional solids are in good agreement with the observed and also justify previous vibrational analyses based on one-dimensional chain models; the non-coincidence of the infrared and Raman mode pairs can be explained as factor-group (Davydov) splitting. The exceptions are one pair of modes in the librational region, for which band assignments based on a one-dimensional chain model need to be revised, as well as the five pseudo-translational modes that exist only in a three-dimensional treatment. The observed pressure dependence of Raman bands in the stretching region, which red-shift with pressure, is accounted for by theory only qualitatively, while that in the pseudo-translational region is reproduced with quantitative accuracy. The present calculation proves to be limited in explaining the complex pressure dependence of the librational modes. The hydrogen-amplitude-weighted phonon DOS at 0 GPa is much less structured than the DOS obtained from one-dimensional models and may be more realistic in view of the also broad, structureless observed inelastic neutron scattering spectra. All major observed peaks can be straightforwardly assigned to the calculated peaks in the DOS. With increasing pressure, MP2 predicts further broadening of bands and breach of the demarcation between the pseudo-translational and librational bands.

A linear-scaling, local-basis, electron-correlation method based on a truncated many-body expansion of energies has been applied to crystalline hydrogen fluoride in three dimensions. The energies, equilibrium atomic positions, lattice constants, and dipole moments of the two structures (polar and nonpolar) have been determined, taking account of one- and two-body Coulomb (electrostatic), exchange, and correlation interactions exactly and three-body and higher-order Coulomb interactions approximately within certain truncation radii. The longer-range two-body Coulomb interactions are also included to an infinite distance by computing the Madelung constant. The second-order Møller−Plesset perturbation method has been used in conjunction with the aug-cc-pVDZ and aug-cc-pVTZ basis sets for correlation. Counterpoise corrections of the basis-set superposition errors have also been made. Predicted relative energies show that the nonpolar arrangement is considerably more stable than the polar one, establishing the precise three-dimensional structure of this crystal and finally resolving the controversy. The computed lattice constants of the nonpolar configuration agree with the observed to within 0.3 Å.

A linear-scaling electron-correlation method based on a truncated many-body expansion of the energies of molecular crystals has been applied to solid hydrogen fluoride. The energies, structures, harmonic, and anharmonic frequencies of the infrared- and/or Raman-active vibrations, phonon dispersions, and inelastic neutron scattering (INS) of the solid have been simulated employing an infinite, periodic, one-dimensional zigzag hydrogen-bonded chain model. The Hartree–Fock, second-order Møller–Plesset (MP2), coupled-cluster singles and doubles (CCSD), and CCSD with a noniterative triples correction [CCSD(T)] methods have been combined with the aug-cc-pVDZ and aug-cc-pVTZ basis sets and, in some instances, the counterpoise corrections of the basis-set superposition errors. The computed structural parameters agree with the observed within 0.1–0.2 Å and a few degrees, and the anharmonic frequencies obtained by vibrational MP2 allowing two-phonon couplings reproduce the observed frequencies semiquantitatively if the potential energy surface is obtained by a correlated theory. They support the revised infrared and Raman band assignments of librational modes made by Hirata and Iwata (J Phys Chem A 1998, 102, 8426) and provide more detailed assignments of the observed INS features.