Embedded and DFT Calculations on the Crystal Structures of Small

Article pubs.acs.org/crystal

Embedded and DFT Calculations on the Crystal Structures of Small Alkanes, Notably Propane A. Daniel Boese*,‡,† and Joachim Sauer‡ ‡

Humboldt Universität zu Berlin, Institut für Chemie, Unter den Linden 6, D-10099, Berlin, Germany Karl-Franzens Universität Graz, Institut für Chemie, Heinrichstrasse 28/IV, A-8010, Graz, Austria

†

Downloaded via UNIV OF SUSSEX on August 8, 2018 at 16:35:07 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: A hybrid method is introduced and applied that combines second-order Møller−Plesset perturbation theory (MP2) for cluster models with density functional theory for periodic models to obtain structures and energetics for molecular crystals of the alkane species, notably propane. Calculated lattice energies elucidate the fact that propane has the lowest melting temperature of all organic compounds. The hybrid method provides an alternative to currently existing methods, which use either periodic MP2 for the determination of periodic crystals or periodic DFT+D. We also compare the hybrid method to a number of modern density functionals which are commonly used for the computation of molecular crystals. Depending on the functional used, we get a large range of lattice energies. From our QM:QM results, we can deduce that the achievements of DFT+D for computing molecular crystals largely relies on error cancellation effects between the functional and its dispersion coefficients. be embedded into an empirical force field. Here, even MP2 or CCSD(T) has been embedded successfully, and almost any desired accuracy can be achieved. However, the embedded ab initio method has to go to large distances, and in some cases, where many-body interactions become important, these have to be included, making the method computationally demanding. In this work, we suggest an embedding in the spirit of the former QM:MM methods (although we use a subtractive rather than additive scheme), but into a QM:QM manner, using DFT +D as the embedding method. This way, we treat the manybody interactions by the already rather accurate DFT+D, combining it with more exact post-Hartree−Fock methods. By doing so, we can achieve an increased accuracy for these calculations. Currently, there are some, but not many benchmark sets for molecular crystals: Among others, the X23 benchmark of Reilly and Tkachenko,36 which is an extension of the C21 benchmark of Otero de la Roza and Johnson37 as well as the benchmark set of Maschio, Civalleri, Ugliengo, and Gavezotti38 are the best known. In this study, although benchmarking some density functionals as well as the newly proposed method, we just restrict ourselves to the crystal structures of the first alkanes: methane, ethane, propane, and n-butane. Propane has the lowest melting temperature of all organic compounds.39−41 When analyzing the packing of the shortchain alkanes it becomes obvious that this effect can be

1. INTRODUCTION In recent years, advances in method developments have been leading to increasingly accurate predictions of structures of molecular crystals. This has been highlighted by the impressive results which were made in the prediction of crystal structures for organic molecules and, recently, also organic salts.1,2 Much of this progress can be attributed to the increased use of density functional theory (DFT) including dispersion coefficients. Nowadays, a plethora of methods treating the nonincluded dispersion in density functionals exists.3−17 Many of these methods have thus proven extremely successful. However, to predict polymorphism (i.e., the possibility that two crystal forms exist for the same molecule, possibly even at the same pressure and temperature), we probably need to attain an accuracy which may be beyond that of DFT including such dispersion coefficients (DFT+D).18,19 Noteworthy, especially for hydrogen bonded systems, DFT does not perform as well as may be necessary to describe molecular crystals with the required accuracy. Here, the errors of DFT+D are at the minimum 3 kJ/mol per hydrogen bond for the best functionals,20,21 which implies that, for example, for molecular crystals with many hydrogen bonds present in the unit cell, these errors add up.22 These problems have been observed, for example, for ice phases.23 In order to go beyond the accuracy of DFT+D, there are several methods at hand: Local MP2 methods for periodic systems, which, however, do not yet give the possibility to optimize periodical structures.24−27 Another, also very promising approach, is the use of embedding techniques like QM:MM for solids,28−35 in which any quantum mechanical method can © 2017 American Chemical Society

Received: November 15, 2016 Revised: March 15, 2017 Published: March 15, 2017 1636

DOI: 10.1021/acs.cgd.6b01654 Cryst. Growth Des. 2017, 17, 1636−1646

Crystal Growth & Design

Article

In all cases of D3, Becke-Johnson damping has been employed.52,53 For methane, 2 × 2 × 2 k-points have been used; for ethane, 5 × 3 × 5 (since the MBD dispersion is computed in reciprocal space, for PBEMBD and PBE0-MBD, we employed a 8 × 5 × 8 k-point grid); for propane, 6 × 2 × 3 (MBD: 10 × 4 × 6); and for n-butane, 4 × 3 × 3 (MBD: 6 × 6 × 5). In order to compute sublimation enthalpies, the single molecules have been optimized in a cell of 17 × 17 × 17 Å3. As a limitation of VASP, all frequencies are calculated at the Γ-point only. Hence, we used the phonopy program suite54 in order to test the effect of more k-points on the zero-point, enthalpy, and, especially, entropy contributions utilizing PBE-D3. In all cases, cell volumes, densities, packing coefficients, and lattice energies have been corrected for thermal effects. Cell volumes, densities, and packing coefficients have been decreased by experimental estimates, and lattice energies by computed frequency contributions. 2.2. Embedding Calculations. For the embedding calculations, we need to treat both the cluster calculations as well as the periodic calculations at an equal footing. The periodic values are again obtained using the VASP suite of programs using a 1200 eV cutoff and hard potentials. Here, we utilized the PBE+D2 functional for embedding. The cluster calculations are performed with TURBOMOLE,55 using both RI-MP2 with the aug-cc-pVTZ56,57 and an extrapolated aug-ccpV{D,T}Z basis set. For the latter calculations, the Hartree−Fock energy was extrapolated with an l−5 formula,58 whereas the correlation energy was extrapolated with l−3,59 with l the cardinal angular momentum quantum number of the basis set (in this case basis set of a double-ζ quality corresponds to l equal two and the basis set of tripleζ quality to three). The interactions of the molecular fragments are calculated at the highest level of theory possible and are embedded into a larger “sea” treated at DFT+D as a more approximate, low level method,60,61 which have been linked by a modified version of the QMPOT program:62

understood by simply representing the alkane molecules via their projection of their outer, molecular shape as polygons arranged in layers: Propane can be represented by a pentagon, whereas methane, ethane, and n-butane would be triangles, rectangles, and inclined hexagons, respectively.40 The latter three molecules acquire a denser packing than the (C2symmetric) pentagons. Thus, the packing coefficient would first increase from methane to the heavier ethane and next for propane, drops dramatically even below that of methane. Then it follows the well-known even/odd melting temperature alternation for all heavier n-alkanes.39,42 Recently, the low melting temperature of propane was attributed to a vast number of effects: “The unfavorable shape, low electrostatic potential on molecular surface, and electrostatic mismatch between loosely packed molecules, the overall cohesion forces are weaker, and thermal vibrations break them more easily”.41 From a theoretical point of view, we have clearly two effects responsible for the lattice energy of a molecular crystal: one is density, and the other one is the interaction energy of a given single molecule with all its surrounding congeners. In addition, entropy effects can play a large role for molecular crystals. The paper thus comprises the following three objectives: (1) to introduce a new QM:QM method to compute molecular crystals, using post-Hartree−Fock methods embedded into DFT+D; (2) to give a benchmark of DFT and the QM:QM method for the small alkanes, which could be used to extend the currently existing benchmark sets for molecular crystals; (3) to analyze the interactions and discern as to why propane has the lowest melting point of all organic compounds.

2. COMPUTATIONAL AND THEORETICAL DETAILS 2.1. Periodic DFT Calculations. All periodic calculations have been performed with VASP, employing a large 1200 eV cutoff and hard potentials.43 For comparing different treatments of dispersion, we employed the PBE44 and hybrid PBE045 functional with the four D2,3, D3,5, TS6, and MBD10,11 dispersion corrections. The PBE functional is the main computational tool used in solid state physics and has been constructed using perturbations of the uniform electron gas (UEG). PBE0 mixes in 25% of Hartree−Fock exchange. Concerning the different dispersion contributions, D2 uses C6 coefficients and damping functions for a pairwise term in order to include the missing dispersion into DFT, whereas D3 also uses C8 coefficients. TS is a density based method using two-body dispersion contributions, which can be extended to MBD which incorporates many-body effects. For B3LYP,46,47 we employed the available D2, D3, and D2*4 dispersion corrections together with BLYP-D3.48,49 BLYP and B3LYP are the main functional tools used when computing gas phase molecules; for the BLYP functional, Becke 88 exchange is combined with the LYP correlation functional which is based on the correlation hole obtained by the second-order Hartree−Fock density matrix; B3LYP additionally adds Hartree−Fock exchange and correlation from the local spin density approximation. Furthermore, we tested HSE-D2 (which is similar to PBE0, with the difference that the HF exchange tends to zero for large interelectron distances; thus we also used the D2 dispersion coefficients of PBE0),50 RPBE-D3 (which is a variation of the PBE functional),51 optB88vdW,13,14 and the vdW-DF212 functional/dispersion combinations. The latter two functionals combine different exchange (a changed B88 and RPBE) functionals with a long-range two-body correlation functional (vdW) which includes dispersion directly.

N = Monomers

∑

Large E HL:LL = E LL −

Small Small [E LL (n) − E HL (n)]

n N = Monomers

∑

=

Small (E HL (n)) + ΔE LR

n

(1)

The first calculation ELarge implies a periodic computation of the LL solid (implying periodic boundary conditions and LL for low-level). This contrasts ESmall HL , the high-level molecular calculation on the embedded, small cluster (superscript C). The low-cluster calculation ESmall of the cluster represents a calculation of the embedded cluster LL with the periodic code used for the host calculation within a large box. The long-range correction in eq 1 emerges entirely from the periodic calculations employed. The same previously mentioned, subtractive scheme can be applied to gradients, too. The dominating part stems from the fragment dimers (or, particularly, dimers of molecules) N = Dimers

∑

Large E HL:LL = E LL −

Small Small [E LL (i) − E HL (i)]

i N = Monomer ∈ Dimers

−

∑

Small Small [E LL (j) − E HL (j)]

j N = Monomer ∉ Dimers

+

∑ k

Small Small [E LL (k) − E HL (k)]

(2)

Another notation of eq 2 is 1637


Crystal Growth & Design N = Dimers Large E HL:LL = E LL −

∑ i

Article

⎤ ⎡ (E Small (i) − E Small (i)) LL HL Dimer ⎥ ⎢ ⎥ ⎢ (E Small (i) E Small (i)) − HL Monomer1 ⎥ ⎢ − LL ⎥ ⎢ Small Small ⎣− (E LL (i) − E HL (i))Monomer2 ⎦

structure of the deuterated methane in the P4m2 space group is available.65 This is probably because the nondeuterated methane crystal structure is disordered. The crystal structure of methane has been previously calculated with PBE+D by one of us in the space group Fm3m.66 Ethane (P21/n) has been determined more accurately.67,68 Data on propane (P21/n) and n-butane (P21/c) are more recent.39 The melting temperatures of the n-alkanes, sublimation enthalpies (probably accurate to one kJ/mol69),70−72 cell volumes, and the number of molecules in the unit cell have been summarized in ref 40 and are listed here in Table 1 together with the cell volumes. Since the cell volumes of methane, ethane, and n-butane have been measured at higher temperatures (with temperatures around 80−90 K) using liquid nitrogen, only propane has been measured using liquid helium (at a temperature around 20−30 K); we corrected the cell volume for the former three species. The expansion of the cell volume from 20 to 30 to 80 K has been measured for propane,39 and we assume a similar expansion for the other four alkanes computed. Furthermore, we expect the expansion from 0 to 30 K to be negligible. The volume was thus scaled by a factor of 0.928 resulting at column 5 of Table 1, and we expect the final cell volumes to be accurate by ±3%. Concerning the experimental sublimation energies ΔHsublimation of column 2, these are easier to compare to theory when converted into lattice energies ΔH lattice. As the sublimation enthalpy is defined by (s denotes the solid, g the gas phase)

N = Monomers

+

∑ j

Small Small [E LL (j) − E HL (j)]

(3)

Using eq 2 rather than eq 3, several monomer calculations are skipped, as only the monomers which were counted more than once (or not counted at all) have to be added to the equation. Whenever the fragments approach each other, they are included as dimers, for which a certain threshold distance is chosen by the user. As we are already treating the periodic system by a rather accurate method, we use a cutoff distance of 5 Å. Moreover, the fragments are allowed to interact periodically and exceed the boundary of the periodic box. As an example, Figure 1 illustrates this situation for the 1dimensional case of the water dimer. Here, the first water fragment is

g−s g−s g−s g−s ΔHsublimation = ΔE lattice + ΔEvib + ΔErot + ΔEtrans − pV g−s g−s = ΔE lattice + ΔEvib + 4RT

Figure 1. Interaction of the one-dimensional water chain within a periodic box of 6.0 Å. The dimer fragments are marked in light red/ green.

(4)

we have to take into account temperature contributions to the enthalpy, zero-point energies ΔEg−s vib , the missing rotational and translational degrees of freedom as well as pV = RT. As temperatures, we used the ones reported for the sublimation enthalpies: 60 K for methane (as a value between 50 and 90 K), 90 K for ethane and propane, and 110 K for n-butane. With these corrections, we arrive at the lattice energies reported in column 3 by the calculations performed in Table 6 (see below). 3.2. Density Functional Theory. As a second step, we consider the performance of various DFT+D methods to compute sublimation enthalpies, cell volumes, densities ρ 1.66 × Z × M molecule (which are computed by the formula ρ = ) and V

interacting twice with the second fragment: once within the box boundaries and the second exceeding the box boundaries. This way, we are able to effectively employ a local post-Hartree−Fock method embedded into DFT+D.

3. RESULTS AND DISCUSSION 3.1. Experimental Reference Data. The methane crystal structure has been determined by X-ray crystallography (Figure 2).63,64 In the Cambridge Structural Data Base (CSD), only the

cell

Figure 2. Crystal structures of methane (space group P4m2), ethane (P21/n), propane (P21/n), and n-butane (P21/c) (from left to right). In the unit cell of methane, 32 methane molecules (Z = 32) are present, for ethane Z = 2, for propane Z = 4, and for n-butane Z = 2. Hydrogen atoms are depicted in white, carbon atoms in brown. 1638



Article

Table 1. Experimentally Determined Cell Volumes and Melting Temperatures of the Alkanes (n = 1−4) melting temp. [K]

sublimation enthalpy [kJ/mol]

derived lattice energy [kJ/mol]

cell volume [Å3]

cell volume V at 0 K [Å3]

molecules in cell (Z)

90.7 89.9 85.5 134.9

9.7 (53−91 K) 20.5 (90 K) 28.5 (86 K) 35.9 (107 K)

10.3 22.1 28.3 38.2

1560 139 364.9 239.1

1448 129 364.9 221.9

32 2 4 2

methane ethane propane n-butane

Table 2. Performance of Several DFT Functionals for the Computation of the Marked Properties of the Methane Crystala methane

functional PBE

Dispersion

exp.

D2

D3

lattice energy cell volume density ρ packing coefficient

10.3 1448c 0.589 61.1e

14.2 1349 0.630 67.2

13.2 1558 0.546 58.1

Dispersion

exp. b


10.3 1448c 0.589 61.1e

TS 16.8 1509 0.563 60.0 B3LYP

BLYP

RPBE

optB88

DF2

MBD

D3

D3

vdW

vdW

13.1 1534 0.554 59.1

10.9 1492 0.570 60.5

17.4 1533 0.554 59.3

15.3 1523 0.558 59.4

13.2 1626 0.518 55.3 PBE0

D2d

D3

D3

TS

MBD

15.8 1202 0.707 74.4

11.7 1463 0.584 61.5

12.0 1528 0.558 58.9

14.9 1488 0.572 60.6

11.3 1480 0.576 60.8

The lattice energy is given in kJ/mol, the cell volume in Å3, the density in g cm−3 and the packing coefficient in %. bThe sublimation enthalpy has been corrected by the zero-point energies and thermal corrections, see Table 1. cCorrected to take thermal effects into account, see Table 1. d Calculation performed at the Γ-point. eCalculated at the estimated, corrected cell volume. a

Table 3. Performance of Several DFT Functionals for the Computation of the Marked Properties of the Ethane Crystala ethane

functional BLYP

RPBE

optB88

DF2

Dispersion

exp.

D2

D3

PBE TS

MBD

D3

D3

vdW

vdW


22.1b 129c 0.772 68.4e

28.7 113.8 0.875 79.2

24.2 132.0 0.755 68.7 B3LYP

34.1 127.7 0.780 70.7

25.7 129.4 0.770 69.9 HSE

23.8 125.8 0.792 71.8

23.2 138.0 0.723 65.7

33.5 128.3 0.776 70.6

28.6 132.6 0.751 68.0

Dispersion

exp.

D2d

D2*d

D3

D2d

D2d

D3

TS

MBD

32.1 105.5 0.944 83.9

14.3 121.20 0.822 74.3

26.9 118.5 0.841 76.0

21.2 114.2 0.872 78.4

20.0 114.2 0.872 78.4

23.1 127.5 0.781 70.4

31.5 127.0 0.784 70.6

23.7 125.3 0.797 71.8


b

22.1 129c 0.772 68.4e

PBE0

The lattice energy is given in kJ/mol, the cell volume in Å3, the density in g cm−3 and the packing coefficient in %. bThe sublimation enthalpy has been corrected by the zero-point energies and thermal corrections, see Table 1. cCorrected to take thermal effects into account, see Table 1. d6 × 6 × 6 k-points. eCalculated at the estimated, corrected cell volume. a

the Kitaigorodskii packing index (packing coefficient),73 in which the vdW-radii for carbon of 1.70 Å and for hydrogen of 1.20 Å are employed. The density ρ can be derived directly from the cell volume and the packing coefficient from the space filling of the molecule in the unit cell. The latter two quantities are measures often reported by crystallographers in order to assess the space used by the molecule in a molecular crystal. The most commonly used underlying functional when calculating molecular crystals is PBE. The older D2 parametrization of Grimme is well-known to overestimate the dispersion in molecular crystals, and has been reparametrized soon after D2 was published4,74 to give better descriptions of molecular crystals. The newer D3 parametrization yields values closer to experiment for methane, as shown in Table 2. Nevertheless, PBE-D3 now overestimates the cell volume. PBE-

TS, on the other hand, shows a somewhat surprising behavior, as both overestimate the cell volume (and consequently underestimate the density and packing coefficient) and the lattice energy. The newer MBD correction still overestimates both lattice energy and the cell volume; a method which is too strongly bound should underestimate the cell volume. When moving to other density functionals, we get a much diversified performance; the lattice energies vary from 10.9 (BLYP+D3) to 17.4 (optB88-vdW) kJ/mol, with all density functionals underestimating the experimental lattice energy. The cell volumes range from 1202 (B3LYP+D2) to 1626 Å3 (RPBE +D3). The same deviations similar to methane can be seen for ethane in Table 3: The different treatment of dispersion has a rather large impact for both functionals PBE and PBE0, 1639



Article

Table 4. Performance of Several DFT Functionals for the Computation of the Marked Properties of the Propane Crystala propane

functional exp.

Dispersion

PBE

BLYP

RPBE

optB88

DF2

D2

D3

TS

MBD

D3

D3

vdW

vdW


28.3b 364.9 0.801 63.3

35.1 324.4 0.901 76.0

29.5 375.2 0.779 66.1 B3LYP

43.1 351.3 0.832 70.5

32.4 359.3 0.813 69.0 HSE

30.5 348.6 0.838 71.0

28.3 390.5 0.749 63.7

42.7 356.0 0.821 69.9

35.7 369.2 0.791 67.2

Dispersion

exp.

D2c

D2*c

D3

D2c

D2c

D3

TS

MBD

41.1 302.5 0.966 80.2

38.3 343.1 0.851 72.0

30.9 344.5 0.850 71.5

27.4 326.1 0.896 75.2

25.2 324.3 0.901 75.5

28.6 356.9 0.820 66.1

39.1 349.1 0.839 70.3

29.7 351.5 0.840 70.5


b

28.3 364.9 0.801 63.3

PBE0

The lattice energy is given in kJ/mol, the cell volume in Å3, the density in g cm−3, and the packing coefficient in %. bThe sublimation enthalpy has been corrected by the zero-point energies and thermal corrections, see Table 1. c3 × 3 × 3 k-points. a

Table 5. Performance of Several DFT Functionals for the Computation of the Marked Properties of the n-Butane Crystala n-butane

Functional exp.

Dispersion lattice energy cell volume density ρ packing coefficient

38.2b 221.9c 0.870 64.1e

BLYP

RPBE

optB88

DF2

D2

D3

PBE TS

MBD

D3

D3

vdW

vdW

45.5 194.4 0.991 73.5

39.1 224.1 0.859 70.1 B3LYP

56.3 214.9 0.896 73.1

41.9 221.1 0.871 71.4 HSE

40.5 213.1 0.904 74.0

36.6 231.0 0.829 67.9

55.6 218.8 0.880 72.1

46.5 226.3 0.851 69.2

D3

D2

d

d

D3

TS

MBD

40.4 212.44 0.907 68.0

30.7 194.1 0.992 73.4

37.7 218.5 0.881 71.7

52.8 212.4 0.907 73.6

38.5 215.3 0.895 72.7

d

Dispersion

exp.

D2


38.2b 221.9c 0.870 64.1e

56.0 180.5 1.067 77.9

D2*

d

33.5 208.32 0.925 69.0

PBE0 D2

31.2 193.8 0.994 73.4

a The lattice energy is given in kJ/mol, the cell volume in Å3, the density in g cm−3 and the packing coefficient in %. bThe sublimation enthalpy has been corrected by the zero-point energies and thermal corrections, see Table 1. cCorrected to take thermal effects into account, see Table 1. d3 × 3 × 3 k-points. eCalculated at the estimated, corrected cell volume.

between 302.5 and 390.2 Å3 (propane) and 180.5 and 232.2 Å3 (n-butane) for the cell volume. Concerning the melting temperatures, the trend of these can only be deduced from the respective packing coefficients. This increases from 61.1 (methane) to 68.4 (ethane), and then goes down again to 63.3 (propane) and increases again to 64.1 (nbutane). The results for most functionals follow this trend. Nevertheless, this is in line with the packing coefficients following the known alternating trend in the melting temperatures of the n-alkanes.39,42 Summarizing the results of the DFT+D methods used, we can find the following ranking of the relative errors of the cell volume (see Supporting Information): PBE-MBD < optB88vdW ≈ PBE0-D3 < vdW-DF2 ≈ PBE-TS ≈ PBE0-MBD ≈ PBE0-TS ≈ PBE-D3 ≈ BLYP-D3 ≪ B3LYP-D3 ≪ RPBE-D3 ≪ PBE-D2 ≪ B3LYP-D2. The less approximate dispersion corrections vdW, MBD, and TS show a good performance, whereas the D3 methods (with the exception of PBE0-D3) are in general not as accurate. The use of the old D2 methods cannot be recommended anymore. In addition, DFT+D functionals based on PBE have smaller errors than B3LYP, BLYP, and RPBE, and hybrid functionals are not much more accurate than their GGA counterparts.

whereas the difference between the two functionals because of the exact exchange is rather small in comparison. Both PBE-D3 and PBE-MBD yield reasonable results for lattice energies and cell volumes. The older dispersion−parametrizations PBE-D2, PBE-TS, and PBE0-TS overestimate, whereas PBE0-D2 underestimates the sublimation enthalpy. Although PBE0-D2 underestimates the sublimation enthalpy (as well as HSE-D2), it also underestimates the cell volume. Whereas the differences between PBE0 and PBE are rather small, the differences between BLYP and B3LYP are large, and those functionals, as observed before,75,76 show larger deviations (in contrast to, for example, molecular calculations)in all cases the cell volume is underestimated. In general, we can note that the results can again vary between 105.5 and 138.0 Å3 for the cell volume and between 20.0 and 34.1 kJ/mol for the lattice energy. When scrutinizing the calculations of the propane and nbutane crystals in Tables 4 and 5, we get a similar picture with larger deviations: Again, PBE-D3, PBE-MBD, and their hybrid counterparts yield rather good lattice energies and cell volumes, whereas most of the other combinations show large differences between the different functional/dispersion combinations: They vary between 25.2 and 43.1 kJ/mol (propane) and 30.7 and 56.3 kJ/mol (n-butane) for the sublimation enthalpies and 1640



Article

Table 6. Zero Point Vibrational Energy (ZPVE), Thermal Effects on the Sublimation Enthalpy and the Entropies of the Alkanes in the Rigid-Rotor-Harmonic Oscillator (RRHO) Approximation per Moleculea molecule CH4

C2H6

C3H8

C4H10

method (frequency kpoints) PBE-D2 PBE-D3 PBE-D3 PBE-D2 PBE-D3 PBE-D3 PBE-D3 PBE-D2 PBE-D3 PBE-D3 PBE-D3 PBE-D2 PBE-D3 PBE-D3 PBE-D3

(2 × 2 × 2)

(3 × 2 × 3) (4 × 3 × 3)

(3 × 1 × 2) (4 × 2 × 2)

(2 × 2 × 2) (3 × 3 × 2)

ZPVE (kJ/mol) 3.5 2.6 2.7 5.1 2.0 3.4 3.3 4.5 2.3 2.2 2.0 5.0 3.7 3.8 4.2

ΔHT − ΔH0 (kJ/mol)

entropy ΔS (J/mol K)

entropy ΔS × T (kJ/mol)

ΔH − ΔS × T (kJ/mol)

ZPVE + ΔHT − ΔH0 − 4RT

−0.1

−23.67

−1.4

4.5 4.2

0.8

1.3 1.4

12.67 14.57

1.1 1.3

4.3 3.6 3.3

1.6 1.6

0.7 0.8

12.65 14.97

1.1 1.3

2.4 1.8 1.5

−0.1 −0.2

2.2 2.1

37.76 23.24

4.2 3.7

3.9 1.9 2.6

2.1 2.3

a

For methane as a hindered rotor, the approximation may be inaccurate. For the temperature in each of the calculations, we used 60 K for methane, 90 K for both ethane and propane, and 110 K for n-butane, close to the temperatures where the sublimation enthalpies of the alkanes were measured (see Table 1). All values were obtained from PBE+D2 with an energy-cutoff of 800 eV and from PBE+D3 with an energy-cutoff of 1200 eV.

For the lattice energies of the alkane crystals, the D3 methods show better performance: PBE0-MBD ≈ BLYP-D3 ≪ PBE0-D3 < PBE-D3 < RPBE-D3 ≈ B3LYP-D3 ≪ PBE-MBD ≪ PBE-D2 < vdW-DF2 ≪ PBE0-TS ≪ B3LYP-D2 ≪ optB88vdW ≈ PBE-TS. Here, PBE-TS, PBE0-TS, vdW-DF2, and optB88-vdW, which yielded good cell volumes, give much worse results than most of the other functionals. In contrast, B3LYP and BLYP seem to be a better method of choice. With these results in mind, PBE-D3 and PBE-MBD are most likely the best options when calculating dispersion-bound molecular crystals like the alkanes studied (for energies, PBE0-MBD and PBE0D3 seem to be slightly better), as we obtain low errors for both lattice energies and cell volumes with these methods. Finally, we need to scrutinize the zero-point energy, entropy, and thermal effects on the total energies which are less easy to compute. Unfortunately, these effects turn out to be nonnegligible, ranging from 2 to 7 kJ/mol when comparing to experimental sublimation enthalpies, as can be seen in Table 6. In general, we need to increase the k-point grid for small unit cells like the one of ethane, as the Γ-point VASP calculations for ethane differ by up to 1.4 kJ/mol from the values obtained using 3 × 2 × 3 k-points for the frequency calculation. When going further to 4 × 3 × 3 k-points, the difference is 1.3 kJ/mol, implying that this calculation is converged. We expect the calculations to be converged to 0.2 kJ/mol for the k-points, and the temperature effects to the enthalpies also to 0.1 kJ/mol. Effects of going beyond the rigid-rotor-harmonic-oscillator (RRHO) approximation are probably still larger than these contributions and a maximum of 10% of the zero-point vibrational energy with 0.5 kJ/mol. These errors are probably in the range of the uncertainties in the experimental enthalpies of sublimation. For the calculation of the lattice energy out of the sublimation enthalpy, we have to subtract 4RT (at 90 K) for the missing degrees of freedom of the gas-phase molecule. The lattice energy can be thus calculated from the (positive) sublimation energy by

exp exp calc ΔE lattice (0 K) = ΔHsublimation (T ) + ΔEvib − ΔH calc(T )

+ ΔH calc(0 K) − 4RT

(5)

with the vibrational contribution larger than zero (with g for the gas phase and s for the solid) g s ΔEvib = Evib − Evib

(6)

and the temperature contribution of the enthalpy ΔH calc = ΔH calc,s − ΔH calc,g

(7)

When adding all these terms together, we obtain the values of the last column in Table 6. These are the values used to calculate the lattice energies in Tables 2−5 together with the experimental sublimation energies from Table 1. Concerning the very low melting point of propane, we deduce from our data that this is caused by the molecular structure and the crystal packing. Whereas the lattice energies increase by a large amount when going from methane to ethane and from propane to n-butane, the lattice energy difference between ethane and propane is small (see Supporting Information). This is already visible when computing the dissociation energies of the methane, ethane, propane, and nbutane dimers.77 From an energetic point of view, this can be entirely attributed to the dispersion energy: When investigating functional and dispersion contribution separately, as it is done in the DFT+D scheme, the lattice energy is in all cases smaller than its dispersion correction contribution. Zero-point energies as well as temperature effects and electrostatics have, other than previously claimed,41 even the opposite effect and reduce the lattice energies. Propane is the only molecule in this series having a dipole moment. As the electrostatic contribution to the lattice energy is repulsive, the difference is larger when going from methane to ethane and from propane to n-butane than comparing ethane to propane, resulting in smaller differences for the dispersion and overall lattice energies. Keeping in mind especially the large deviations of DFT+D, we now turn to the method proposed in section 2.2 in order to see how such embedding methods perform. 1641



Article

3.3. QM:QM Method Proposed. As gradients for cell volumes have not been implemented yet, we start from the optimized PBE+D2 structure and alter the cell volumes obtained. When assuming that van der Waals interactions are the predominant interaction between the molecules, this is probably a good estimate of the cell dimensions. For the cubic cell of methane in Figure 3, this is of course not a problem.

Figure 4. Lattice energy of the ethane crystal calculated with MP2/ aug-cc-pVTZ and extrapolated aug-cc-pV{D,T}Z MP2 embedded into DFT+D2. A Morse potential curve has been fitted through the points obtained at the individual cell volumes.

Figure 3. Lattice energy of the methane crystal calculated with MP2/ aug-cc-pVTZ embedded into DFT+D2. A Morse potential curve has been fitted through the points obtained for the individual cell volumes.

Here, the minimum calculated by a Morse potential curve is at a lattice energy of 11.5 kJ/mol (experiment: 10.3 kJ/mol) and 1492 Å3, which is rather close to the experimentally determined value of 1444 Å3. Of course, we do not have the error cancellation effects which are in place in pure DFT+D, which does not necessarily mean that the values have to improve when going from PBE+D2 to embedded MP2/aug-ccpVTZ:PBE+D2. Also, when going from an aug-cc-pVTZ to an extrapolated aug-cc-pV{D,T}Z basis set for MP2, the dissociation energies do not necessarily improve because of error cancellation. When looking at the results of the ethane crystal in Figure 4, we already see some deviations from experiment (2.7 kJ/mol and 16 Å3, see Table 3); the lattice energies have been determined here not only with an aug-cc-pVTZ basis set, but also with a MP2:PBE+D2 geometry optimized with an extrapolated aug-cc-pV{D,T}Z basis set: The extrapolated lattice energies are 19.4 (aug-cc-pVTZ) and 15.4 (aug-ccpV{D,T}Z) kJ/mol, respectively. The 4.0 kJ/mol difference emerges partly from the basis set superposition error in the augcc-pVTZ, which we will discuss later. The cell volumes, however, are very little changed by the different basis sets, from 145.4 to 147.3 Å3. For propane in Figure 5, for which the crystal structure at 30 K is available and for which we thus expect very little thermal expansion, the results are closer to the experiment with 31.9 kJ/ mol for the lattice energy (experiment: 28.3 kJ/mol) and 358.3 Å3 for the cell volume (experiment: 364.9 Å3). Finally, for n-butane in Figure 6, the lattice energy is underestimated compared to the experimental value (37.7 compared to 38.2 kJ/mol), whereas the cell volume is about 20% too large (260.7 compared to 221.4 Å3). In order to estimate the accuracy of MP2/aug-cc-pVTZ for the dimers evaluated in the QM:QM method, we then compute the energies of the dimers in comparison to highly accurate

Figure 5. Lattice energy of the propane crystal calculated with MP2/ aug-cc-pVTZ embedded into DFT+D2. A Morse potential curve has been fitted through the points obtained at the individual cell volumes.

Figure 6. Lattice energy of the n-butane crystal calculated with MP2/ aug-cc-pVTZ embedded into DFT+D2. A Morse potential curve has been fitted through the points obtained at the individual cell volumes.

1642



Article

unfavorable molecular shape of propane as suggested with simplified models derived from the X-ray structure determinations. A more complete result would arise from the analysis of the intermolecular interactions like when using symmetryadapted perturbation theory (SAPT).86 This would have to be performed for the whole crystal lattices, but unfortunately the program codes for such calculations do not exist.

CCSD(T) calculations at the basis set limit using the structures optimized with the QM:QM method. For this purpose, we computed MP2 at the basis set limit using an aug-cc-pV{Q,5}Z basis set extrapolation (see Supporting Information). Then, the difference between CCSD(T) and MP2 is estimated with smaller basis sets, namely, those of double-ζ and triple-ζ quality. As both basis sets yield very similar values for the energy difference of CCSD(T) and MP2, we extrapolated these basis set differences and arrive at our best estimate. Then, we sum over all dimers computed by MP2, and compare this value to the MP2/aug-cc-pVTZ and extrapolated MP2/aug-ccpV{D,T}Z values which were used in the QM:QM calculations displayed in Figures 3−6. From these values, we can estimate the energy if we would calculate the energy of the molecular crystal with CCSD(T) rather than with MP2. The long-range correction from PBE+D2, coming mostly from the trimers of PBE+D2 and the long-range correction, reduces the lattice energy for both alkane crystals investigated. This implies that we have two effects counteracting each other: The first are the electrostatics of trimers, tetramers, and so on, which lowers the lattice energy and dominates especially for ethane. The second effect is the long-range dispersion beyond the used cutoff radius of 5 Å, which is raising the lattice energy: It gets larger for propane and is almost completely able to counter the repulsive many-body effect. The final best estimates of the lattice energies (12.9 and 23.4 kJ/mol) for ethane and propane, respectively are, however, deviating more from the experimental values of 22.1 and 28.3 kJ/mol than the MP2/aug-cc-pVTZ:PBE+D2 values (19.4 and 31.9 kJ/mol), mainly due to error cancellation effects. This implies that, even when improving the description of the dimers by itself, they do not necessarily yield better results when comparing to the experiment. Here, a more thorough analysis and larger clusters would have to be calculated, as we have done for the adsorption of CO, CH4, and C2H6 on the MgO (100) surface and CO on the NaCl(100) surface.78−81 Finally, we may come up with a more suitable method for embedding the dimers beyond PBE+D2 or a better embedding scheme using density embedding.82−85 This, however, also implies that for modern density functionals which yield rather good results a lot of error cancellation is in place between short-range and long-range dispersion and electrostatics or, in other words, there is a very fine balance between the density functional itself together with the extra dispersion term going with it. For example, the PBE functional includes some shortrange dispersion already, which is then adsorbed into the parametrization of the additional long-range dispersion. The BLYP functional does not include much of this short-range dispersion, having much larger long-range dispersion contributions (see Supporting Information). When replacing now in the short-range BLYP or PBE with an exact treatment of correlation and short-range dispersion, like it is done for example with CCSD(T), the error cancellation between functional and long-range dispersion is unbalanced, giving worse results despite a better short-range treatment. Referring to the recent discussion about the intermolecular interactions of solid propane responsible for the extremely low melting temperature,40,41 we compare the lattice energies of dimers to the neighboring molecules in the homologous series of n-alkanes. Accordingly, we consider ethane and n-butane plus methane for additional data. Analyzing the differences from data from the Supporting Information, we can state that the dominating lattice energy results from dispersion, i.e., the

4. CONCLUSIONS Computing the crystal structures of small, dispersion-bound alkanes as model systems, we introduced an embedding method for the description of molecular crystals. The QM:QM embedding method uses high-level ab initio theory, and, unlike other contributions, embeds this into periodic DFT +D calculation. While it should be more accurate than DFT+D and MP2 embedded into force fields, the combination of MP2 with PBE+D2 yields less accurate results than the best functionals out of 16 functional and dispersion combinations tested. This suggests that the error cancellation which is in place for modern density functional and dispersion combinations is canceled, and other combinations than the one employed would yield better results. Despite these shortcomings, the new QM:QM method is inherently linear scaling when assuming a linear scaling of the underlying periodic DFT calculation and rather easy to understand and implement. It is also, in comparison to the aforementioned methods, systematically improvable by, for example, using the coupled-cluster method instead of Møller−Plesset perturbation theory, increasing the cutoff radius, computing trimers, and so forth. In addition, when finding a suitable functional instead of PBE +D2, more accurate results can be expected. For charged and hydrogen bonded species, for which modern DFT+D methods are not yielding as high-quality results as for dispersion-bound molecular crystals, methods such as the introduced one may become a viable alternative. Finally, focusing on the rather low melting point of propane, we can state that it can be entirely attributed to the dispersion effects and thus to the crystal packing caused by the molecular structure.

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.cgd.6b01654. All lattice energies, enthalpies, dispersion contributions, cell volumes, as well their differences to experiment together with the MP2 and CCSD(T) basis set extrapolations (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

A. Daniel Boese: 0000-0001-7388-778X Joachim Sauer: 0000-0001-6798-6212 Notes

The authors declare no competing financial interest.

■

REFERENCES

(1) Bardwell, D. A.; Adjiman, C. S.; Arnautova, Y. A.; Bartashevich, E.; Boerrigter, S. X. M.; Braun, D. E.; Cruz-Cabeza, A. J.; Day, G. M.; Della Valle, R. G.; Desiraju, G. R.; van Eijck, B. P.; Facelli, J. C.;

1643



Article

Ferraro, M. B.; Grillo, D.; Habgood, M.; Hofmann, D. W. M.; Hofmann, F.; Jose, K. V. J.; Karamertzanis, P. G.; Kazantsev, A. V.; Kendrick, J.; Kuleshova, L. N.; Leusen, F. J. J.; Maleev, A. V.; Misquitta, A. J.; Mohamed, S.; Needs, R. J.; Neumann, M. A.; Nikylov, D.; Orendt, A. M.; Pal, R.; Pantelides, C. C.; Pickard, C. J.; Price, L. S.; Price, S. L.; Scheraga, H. A.; van de Streek, J.; Thakur, T. S.; Tiwari, S.; Venuti, E.; Zhitkov, I. K. Towards crystal structure prediction of complex organic compounds–a report on the fifth blind test. Acta Crystallogr., Sect. B: Struct. Sci. 2011, 67, 535−551. (2) Reilly, A. M.; Cooper, R. I.; Adjiman, C. S.; Bhattacharya, S.; Boese, A. D.; Brandenburg, J. G.; Bygrave, P. J.; Bylsma, R.; Campbell, J. E.; Car, R.; Case, D. H.; Chadha, R.; Cole, J. C.; Cosburn, K.; Cuppen, H. M.; Curtis, F.; Day, G. M.; DiStasio, R. A., Jr; Dzyabchenko, A.; van Eijck, B. P.; Elking, D. M.; van den Ende, J. A.; Facelli, J. C.; Ferraro, M. B.; Fusti-Molnar, L.; Gatsiou, C.-A.; Gee, T. S.; de Gelder, R.; Ghiringhelli, L. M.; Goto, H.; Grimme, S.; Guo, R.; Hofmann, D. W. M.; Hoja, J.; Hylton, R. K.; Iuzzolino, L.; Jankiewicz, W.; de Jong, D. T.; Kendrick, J.; de Klerk, N. J. J.; Ko, H.Y.; Kuleshova, L. N.; Li, X.; Lohani, S.; Leusen, F. J. J.; Lund, A. M.; Lv, J.; Ma, Y.; Marom, N.; Masunov, A. E.; McCabe, P.; McMahon, D. P.; Meekes, H.; Metz, M. P.; Misquitta, A. J.; Mohamed, S.; Monserrat, B.; Needs, R. J.; Neumann, M. A.; Nyman, J.; Obata, S.; Oberhofer, H.; Oganov, A. R.; Orendt, A. M.; Pagola, G. I.; Pantelides, C. C.; Pickard, C. J.; Podeszwa, R.; Price, L. S.; Price, S. L.; Pulido, A.; Read, M. G.; Reuter, K.; Schneider, E.; Schober, C.; Shields, G. P.; Singh, P.; Sugden, I. J.; Szalewicz, K.; Taylor, C. R.; Tkatchenko, A.; Tuckerman, M. E.; Vacarro, F.; Vasileiadis, M.; Vazquez-Mayagoitia, A.; Vogt, L.; Wang, Y.; Watson, R. E.; de Wijs, G. A.; Yang, J.; Zhu, Q.; Groom, C. R. Report on the sixth blind test of organic crystal structure prediction methods. Acta Crystallogr., Sect. B: Struct. Sci., Cryst. Eng. Mater. 2015, 72, 439−459. (3) Grimme, S. Semiempirical GGA-type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787−1799. (4) Civalleri, B.; Zicovich-Wilson, C.; Valenzano, L.; Ugliengo, P. B3LYP augmented with an empirical dispersion term (B3LYP-D*) as applied to molecular crystals. CrystEngComm 2008, 10, 405−410. (5) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (6) Tkatchenko, A.; Scheffler, M. Accurate Molecular Van Der Waals Interactions from Ground-State Electron Density and Free-Atom Reference Data. Phys. Rev. Lett. 2009, 102, 073005. (7) Bucko, T.; Lebègue, S.; Hafner, J.; Á ngyán, J. G. Improved density dependent correction for the description of London dispersion forces. J. Chem. Theory Comput. 2013, 9, 4293−4299. (8) Bucko, T.; Lebègue, S.; Á ngyán, J. G.; Hafner, J. Extending the applicability of the Tkatchenko-Scheffler dispersion correction via iterative Hirshfeld partitioning. J. Chem. Phys. 2014, 141, 034114. (9) Bultinck, P.; Van Alsenoy, C.; Ayers, P. W.; Carbó Dorca, R. Critical analysis and extension of the Hirshfeld atoms in molecules. J. Chem. Phys. 2007, 126, 144111. (10) Tkatchenko, A.; Di Stasio, R. A.; Car, R.; Scheffler, M. Accurate and efficient method for many-body van der Waals interactions. Phys. Rev. Lett. 2012, 108, 236402. (11) Ambrosetti, A.; Reilly, A. M.; DiStasio, R. A.; Tkatchenko, A. Long-range correlation energy calculated from coupled atomic response functions. J. Chem. Phys. 2014, 140, 018A508. (12) Lee, K.; Murray, E. D.; Kong, L.; Lundqvist, B. I.; Langreth, D. C. Higher-Accuracy van der Waals Density Functional. Phys. Rev. B: Condens. Matter Mater. Phys. 2010, 82, 081101. (13) Klimes, J.; Bowler, D. R.; Michaelides, A. Chemical accuracy for the van der Waals density functional. J. Phys.: Condens. Matter 2010, 22, 022201. (14) Klimes, J.; Bowler, D. R.; Michaelides, A. Van der Waals density functionals applied to solids. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 83, 195131.

(15) Dion, M.; Rydberg, H.; Schröder, E.; Langreth, D. C.; Lundqvist, B. I. Van der Waals Density Functionals for General Geometries. Phys. Rev. Lett. 2004, 92, 246401. (16) Vydrov, O. A.; van Voorhis, T. Nonlocal van der Waals Density Functional Made Simple. Phys. Rev. A: At., Mol., Opt. Phys. 2010, 81, 062708. (17) Roman-Perez, G.; Soler, J. M. Efficient Implementation of a van der Waals Density Functional: Application to Double-Wall Carbon Nanotubes. Phys. Rev. Lett. 2009, 103, 096102. (18) Yang, J.; Hu, W.; Usvyat, D.; Matthews, D.; Schütz, M.; Chan, G. K. Ab initio determination of the crystalline benzene lattice energy to sub-kilojoule/mol accuracy. Science 2014, 345, 640−643. (19) Price, S. L. Lattice energy, nailed? Science 2014, 345, 619−620. (20) Boese, A. D. Density Functional Theory and Hydrogen Bonds: Are we there yet? ChemPhysChem 2015, 16, 978−985. (21) Boese, A. D. Basis Set Limit Coupled-Cluster studies of Hydrogen.Bonded Systems. Mol. Phys. 2015, 113, 1618−1629. (22) Friedrich, J.; Yu, H.; Leverentz, H. R.; Bai, P.; Siepmann, J. I.; Truhlar, D. G. Water 26-mers Drawn from Bulk Simulations: Benchmark Binding Energies for Unprecedentedly Large Water Clusters and Assessment of the Electrostatically Embedded ThreeBody and Pairwise Additive Approximations. J. Phys. Chem. Lett. 2014, 5, 666−670. (23) Brandenburg, J. G.; Maas, T.; Grimme, S. Benchmarking DFT and semiempirical methods on structures and lattice energies for ten ice polymorphs. J. Chem. Phys. 2015, 142, 124104. (24) Maschio, L.; Usvyat, D.; Schütz, M.; Civalleri, B. Periodic local Møller−Plesset second order perturbation theory method applied to molecular crystals: Study of solid NH3 and CO2 using extended basis sets. J. Chem. Phys. 2010, 132, 134706. (25) Del Ben, M.; Hutter, J.; VandeVondele, J. Second Order MøllerPlesset Perturbation Theory in the Condensed Phase: An Efficient and Massively Parallel Gaussian and Plane Waves Approach. J. Chem. Theory Comput. 2012, 8, 4177−4188. (26) Booth, G. H.; Grüneis, A.; Kresse, G.; Alavi, A. Towards an Exact Description of Electronic Wavefunctions in Real Solids. Nature 2012, 493, 365−370. (27) Suhai, S.; Ladik, J. Perturbation Theoretical Calculation of the Correlation-Energy in an Infinite Metallic Hydrogen Chain. J. Phys. C: Solid State Phys. 1982, 15, 4327−4337. (28) Wen, S.; Beran, G. J. O. Accurate Molecular Crystal Lattice Energies from a Fragment QM/MM Approach with On-the-Fly Ab Initio Force Field Parametrization. J. Chem. Theory Comput. 2011, 7, 3733−3742. (29) Nanda, K. D.; Beran, G. J. O. Prediction of organic molecular crystal geometries from MP2-level fragment quantum mechanical/ molecular mechanical calculations. J. Chem. Phys. 2012, 137, 174106. (30) Mörschel, P.; Schmidt, M. U. Prediction of molecular crystal structures by a crystallographic QM/MM model with full space-group symmetry. Acta Crystallogr., Sect. A: Found. Adv. 2015, 71, 26−35. (31) Beran, G. J. O.; Wen, S. H.; Nanda, K. D.; Huang, Y. H.; Heit, Y. Prediction and Calculation of Crystal Structures: Methods and Applications. Top. Curr. Chem. 2013, 345, 59−93. (32) Beran, G. J. O. A new era for ab initio molecular crystal lattice energy prediction. Angew. Chem., Int. Ed. 2014, 54, 396−398. (33) Heit, Y. N.; Nanda, K. D.; Beran, G. J. O. Predicting finitetemperature properties of crystalline carbon dioxide from first principles with quantitative accuracy. Chem. Sci. 2016, 7, 246−255. (34) Binnie, S. J.; Nolan, S. J.; Drummond, N. D.; Alfe, D.; Allan, N. L.; Manby, F. R.; Gillan, M. J. Bulk and surface energetics of crystalline lithium hydride: Benchmarks from quantum Monte Carlo and quantum chemistry. Phys. Rev. B: Condens. Matter Mater. Phys. 2010, 82, 165431. (35) Nolan, S. J.; Bygrave, P. J.; Allan, N. L.; Manby, F. R. Comparison of the incremental and hierarchical methods for crystalline neon. J. Phys.: Condens. Matter 2010, 22, 074201. (36) Reilly, A. M.; Tkatchenko, A. Understanding the role of vibrations, exact exchange, and many-body van der Waals interactions 1644



Article

in the cohesive properties of molecular crystals. J. Chem. Phys. 2013, 139, 024705. (37) Otero-de-la-Roza, A.; Johnson, E. R. A benchmark for noncovalent interactions in solids. J. Chem. Phys. 2012, 137, 054103. (38) Maschio, L.; Civalleri, B.; Ugliengo, P.; Gavezzotti, A. Intermolecular Interaction Energies in Molecular Crystals: Comparison and Agreement of Localized M?ller Plesset 2, DispersionCorrected Density Functional, and Classical Empirical Two-Body Calculations. J. Phys. Chem. A 2011, 115, 11179−11196. (39) Boese, R.; Weiss, H.-C.; Bläser, D. The melting point alternation in the short-chain n-alkanes: Single-crystal X-ray analyses of propane at 30 K and of n-butane to n-nonane at 90 K. Angew. Chem., Int. Ed. 1999, 38, 988−992. (40) Thalladi, V. R.; Boese, R. Why is the melting point of propane the lowest among n-alkanes? New J. Chem. 2000, 24, 579−581. (41) Podsiadlo, M.; Olejniczak, A.; Katrusiak, A. Why Propane? J. Phys. Chem. C 2013, 117, 4759−4763. (42) Baeyer, A. Ueber Regelmassigkeiten im Schmelzpunkt homologer Verbindungen. Ber. Dtsch. Chem. Ges. 1877, 10, 1286− 1288. (43) Kresse, G. et al. VASP 5.4.1; http://www.vasp.at, accessed 12.11.2016. (44) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (45) Adamo, C.; Barone, V. Toward reliable density functional methods without adjustable parameters: The PBE0 model. J. Chem. Phys. 1999, 110, 6158−6170. (46) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab-initio Calculation of Vibrational Absorption and CircularDichroism Spectra Using Density-Functional Force-Fields. J. Phys. Chem. 1994, 98, 11623−11627. (47) Becke, A. D. Density Functional Thermochemistry 0.3. The role of the exact exchange. J. Chem. Phys. 1993, 98, 5648−5652. (48) Becke, A. D. Density-Functional Exchange-Energy Approximation with correct Asymptotic-Behaviour. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (49) Lee, C.; Yang, W.; Parr, R. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron-Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (50) Krukau, A. V.; Vydrov, O. A.; Izmaylov, A. F.; Scuseria, G. E. Influence of the exchange screening parameter on the performance of screened hybrid functionals. J. Chem. Phys. 2006, 125, 224106. (51) Hammer, B.; Hansen, L. B.; Norskov, J. K. Improved adsorption energetics within density-functional theory using revised PerdewBurke-Ernzerhof functionals. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 7413−7421. (52) Johnson, E. R.; Becke, A. D. A post-Hartree-Fock model of intermolecular interactions: Inclusion of higher-order corrections. J. Chem. Phys. 2006, 124, 174104. (53) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem. 2011, 32, 1456−1465. (54) Togo, A.; Tanaka, I. First principles phonon calculations in materials science. Scr. Mater. 2015, 108, 1−5. (55) TURBOMOLE version 6.5, COSMOlogic; 2011 developed by R. Ahlrichs et al., see www.turbomole.de, accessed 12.11.2016. (56) Dunning, T. H. Gaussian-Basis Sets for Use in Correlated Molecular Calculations 0.1. the Atoms Boron through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (57) Kendall, R. A.; Dunning, T. H.; Harrison, R. J. Electron Affinities of the 1st-row Atoms Revisited - Systematic Basis-Sets and Wave-Functions. J. Chem. Phys. 1992, 96, 6796−6806. (58) Parthiban, S.; Martin, J. M. L. Assessment of W1 and W2 theories for the computation of electron affinities, ionization potentials, heats of formation, and proton affinities. J. Chem. Phys. 2001, 114, 6014−6029. (59) Halkier, A.; Helgaker, T.; Jorgensen, P.; Klopper, W.; Koch, H.; Olsen, J.; Wilson, A. K. Basis-set convergence in correlated calculations on Ne, N-2, and H2O. Chem. Phys. Lett. 1998, 286, 243−252.

(60) Tuma, C.; Sauer, J. A hybrid MP2/planewave-DFT scheme for large chemical systems: proton jumps in zeolites. Chem. Phys. Lett. 2004, 387, 388−394. (61) Tuma, C.; Sauer. Treating dispersion effects in extended systems by hybrid MP2: DFT calculations - protonation of isobutene in zeolite ferrierite. Phys. Chem. Chem. Phys. 2006, 8, 3955−3965. (62) Sierka, M.; Tuma, C.; Kerber, T. Humboldt Universität, Berlin, 2011. (63) Bolshutkin, D. N.; Gasan, V. M.; Prokhvatilov, A. I. J. Struct. Chem. 1971, 12, 734. (64) Hazen, R. M.; Mao, H. K.; Finger, L. W.; Bell, P. M. Structure and Compression of Crystalline Methane at High-Pressure and RoomTemperature. Appl. Phys. Lett. 1980, 37, 288. (65) Prokhvatilov, A. I.; Isakina, A. P An X-Ray-Powder diffraction study of Crystalline Alpha-Methane-D4. Acta Crystallogr., Sect. B: Struct. Crystallogr. Cryst. Chem. 1980, 36, 1576−1580. (66) Kerber, T.; Sierka, M.; Sauer, J. Application of semiempirical long-range dispersion corrections to periodic systems in density functional theory. J. Comput. Chem. 2008, 29, 2088−2097. (67) van Nes, G. J. H.; Vos, A. Single-crystal structures and ElectronDensity distributions of Ethane, Ethylene and Acetylene 0.1. SingleCrystal X-Ray Structure Determinations of 2 Modifications of Ethane. Acta Crystallogr., Sect. B: Struct. Crystallogr. Cryst. Chem. 1978, 34, 1947−1956. (68) Amoureux, J. P.; Foulon, F.; Muller, M.; Bee, M. Structures of Ethane, Hexamethylmethane and Hexamethyldisilane in their plastic phases. Acta Crystallogr., Sect. B: Struct. Sci. 1986, 42, 78−84. (69) Chickos, J. S. Enthalpies of Sublimation after a Century of Measurement: A View as Seen through the Eyes of a Collector. Netsu Sokutei 2003, 30, 116−124. (70) Stephenson, R. M.; Malanowski, S. in Handbook of Thermodynamics of Organic Compounds; Elsevier: New York, 1987. (71) Bondi, A. Heat of Siblimation of Molecular Crystals: A Catalog of Molecular Structure Increments. J. Chem. Eng. Data 1963, 8, 371− 381. (72) Geiseler, G.; Quitzsch, K.; Rauh, H.-J.; Schaffernicht, H.; Walther, H. Bildungsenthalpien und Mesomerienenergien von PiBindungssystemen 0.1. Bindungsenthalpien und Mesomerienenergien einiger mehrkerniger Aromaten und verschiedener Pseudoazulene. J. Ber. Phys. Chem. 1966, 70, 551. (73) Kitaigorodskii, A. I. in Organic Chemical Crystallography; Consultants Bureau: New York, 1961; p 107. (74) Neumann, M. A.; Perrin, M. A. Energy ranking of molecular crystals using density functional theory calculations and an empirical van der Waals correction. J. Phys. Chem. B 2005, 109, 15531−15541. (75) Boese, A. D.; Kirchner, M.; Echeverria, G. E.; Boese, R. Ethyl Acetate: X-Ray, Solvent and Computed Strutures. ChemPhysChem 2013, 14, 799−804. (76) Boese, A. D.; Boese, R. Tetrahydrothiophene and Tetrahydrofuran, Computational and X-ray Studies in the Crystalline Phase. Cryst. Growth Des. 2015, 15, 1073−1081. (77) Tsuzuki, S.; Honda, K.; Uchimaru, T.; Mikami, M. Estimated MP2 and CCSD(T) interaction energies of n-alkane dimers at the basis set limit: Comparison of the methods of Helgaker et al. and Feller. J. Chem. Phys. 2006, 124, 114304. (78) Boese, A. D.; Sauer, J. Accurate Adsorption Energies of Small Molecules on Oxide Surfaces: CO-MgO(001). Phys. Chem. Chem. Phys. 2013, 15, 16481−16493. (79) Tosoni, S.; Sauer, J. Accurate quantum chemical energies for the interaction of hydrocarbons with oxide surfaces: CH4/MgO(001). Phys. Chem. Chem. Phys. 2010, 12, 14330−14340. (80) Boese, A. D.; Sauer, J. Accurate Adsorption Energies for Small Molecules on Oxide Surfaces: CH4/MgO(001) and C2H6/ MgO(001). J. Comput. Chem. 2016, 37, 2374−2385. (81) Boese, A. D.; Saalfrank, P. CO Molecules on a NaCl(100) Surface: Structures, Energetics and Vibrational Davydov Splittings at Various Coverages. J. Phys. Chem. C 2016, 120, 12637−12653. 1645



Article

(82) Manby, F. R.; Stella, M.; Goodpaster, J. D.; Miller, T. F. A Simple, Exact Density-Functional-Theory Embedding Scheme. J. Chem. Theory Comput. 2012, 8, 2564−2568. (83) Jacob, C. R.; Neugebauer, J.; Visscher, L. A Flexible Implementation of Frozen-Density Embedding for Use in Multilevel Simulations. J. Comput. Chem. 2008, 29, 1011−1018. (84) Gomes, A. S. P.; Jacob, C. R.; Real, F.; Visscher, L.; Vallet, V. Towards systematically improvable models for actinides in condensed phase: the electronic spectrum of uranyl in Cs2UO2Cl4 as a test case Phys. Phys. Chem. Chem. Phys. 2013, 15, 15153−15162. (85) Heuser, J.; Höfener, S. Wave-function frozen-density embedding: Analytical nuclear ground-state gradients. J. Comput. Chem. 2016, 37, 1092−1101. (86) Jansen, G. Symmetry-adapted perturbation theory based on density functional theory for noncovalent interactions. WIREs Comput. Mol. Sci. 2014, 4, 127−144.

1646


Embedded and DFT Calculations on the Crystal Structures of Small

Recommend Documents