J. Phys. Chem. C 2010, 114, 20523–20530
20523
Molecular Crystals: A Test System for Weak Bonding† T. Todorova and B. Delley* Paul Scherrer Institute Switzerland, CH-5232 Villigen, Switzerland ReceiVed: May 31, 2010; ReVised Manuscript ReceiVed: August 20, 2010
Molecular crystals are classified by the leading order power law for the long-range interaction among the molecular constituents and additionally by hydrogen bonds. It is found that for molecular crystals formally dominated by ionic, electrostatic dipole-dipole or quadrupole-quadrupole, and even van der Waals interactions, reasonably accurate predictions of crystal structure can be made. The PBE density functional approximation leads to a systematic overestimation of lattice parameters of type length reaching 7% on average for the van der Waals class. This is consistent with the findings of earlier investigations. The PBEsol functional has a signed mean deviation of i
l
i,lj(λ)
- Di,lj)2W(Di,lj) + F0(λ)
j
1 + p1 + D4(p2 + p3D4) D4
D < Dcut
where Di,lj is the bond distance between atom i and l, j for the undistorted reference geometry and for the lattice deformation λ, which may involve any combination of ABCRβγ. The weighting function W was inspired from Halgren Lipscomb.22 The parameters p1p2p3 are taken to make W(D) vanish with continuous first and second derivative at the cutoff Dcut. This keeps the number of terms to be calculated for a solid at a finite number. The cutoff radius has to be taken large enough to include some intermolecular atom pair distances. As the function F only relates to internal coordinates, a small extra term is required to arrive at a well-defined minimization problem for all coordinates. This function F0 is conveniently taken as a quadratic with its minimum for all atoms at their original (Ri(λ) - Rλi )2, undistorted coordinates. F0(λ;{Ri(λ)}) ) Σcell0 i λ where Ri are the uniformly distorted atomic coordinates. The function F(λ;{Ri(λ)}) is minimized with respect to all atomic coordinates Ri(λ) to start optimizations at each deformation step. This should yield a fair initial guess of the slight deformations that are forced upon the molecular constituents by the lattice deformations. In simple cases, the initial guess from this minimization results in the uniformly distorted atomic coordinates. The Hessian for the atomic coordinates, as obtained by the previous undistorted lattice, serves as initialization for the optimization at the distorted lattice. This initialization, and its parameters have no effect on the final optimized structure, it only serves to accelerate calculations by an initial guess for the distorted crystal structures. III. Results and Discussion As suggested by the Introduction, the PBEsol functional was used to model all exchange and correlation effects that may be of importance for the intramolecular and intermolecular energy terms. The results of the present lattice optimizations are compared with experimentally determined25-30 lattice parameters. The experimental lattices are measured at low temperature. In the special case of cubic acetylene, the sample was held at 141 K to avoid the low temperature orthorhombic phase. No thermal corrections were applied to the experimental lattice parameters. The calculated lattice parameters are from zero temperature DFT, and no thermal expansion corrections were applied. Experimental lattice parameters and the input crystal structures for our calculations are presented as Supporting Information. For the present purpose, length type lattice parameters A, B, and C are compared and shown as percent deviations ∆A, ∆B, and ∆C in Table 2. Redundant deviations are not listed in the table nor used for statistics. A difficulty with monoclinic and triclinic lattices is to present the deviations of lattice angles in a fashion comparable to that for the lattice parameters with length dimension. Here it was opted to present the deviation of the cube root of the cell volume as additional parameter of length dimension under ∆L. In the case of triclinic
20526
J. Phys. Chem. C, Vol. 114, No. 48, 2010
Todorova and Delley
TABLE 2: Deviations (%) from Experimental Lattice Constants with PBEsol Functionala compound
crystal
∆A
∆B
∆C
∆L
U0
TTF-TCNQ histidine resolved cholamide dihydrate orotic acid
m o m a
Monopole -0.05 0.66 0.44 -0.42 -0.67 -0.54 -1.13 1.22 0.40 3.31 2.31 -0.40
acetanilide alanine aniline apha aspartic acid benzamide benzoic acid caprolactam cinnamide citric acid ethanol glucose glycine h2s ice imidazole mono paracet toluene urea
o o m m m m m m m m m o m o o m m m t
Dipole 8.14 -4.38 6.38 -5.37 -1.18 6.41 1.19 -0.71 -0.91 0.41 2.99 -0.35 2.68 0.39 3.17 0.59 -2.24 0.48 1.49 1.31 2.25 3.21 -2.13 -1.22 -0.75 -0.72 0.90 -6.82 -6.10 -5.04 -1.20 1.90 -0.20 0.62 1.42 2.84 -0.13
acetylene adipic acid anthracene benzene βsuccinic acid biphenyl indathrone indigo1 indigo2 naphtalene pigment red quinone sulfur tere acid1 tere acid2
c m m o m m m m m m a m o a a
Quadrupole -1.87 -0.88 7.84 9.83 2.63 0.22 2.76 2.04 0.87 4.40 0.32 0.68 1.41 -1.71 2.13 2.13 1.08 1.14 2.02 1.75 -1.46 -0.42 3.79 0.70 -0.27 3.42 2.95 5.41 4.01 0.75 1.47 2.16 -1.73 2.28 0.10 0.47 0.31 -0.51 0.46 0.88 -0.82 3.90 1.21
0.20 2.21 1.19 1.67 0.36 0.22 0.81 1.14 0.82 0.35 1.43 0.68 0.56 0.78 1.49 0.76 3.07 0.27 1.77 1.09 0.65 0.42 0.31 0.42 1.22 2.14 1.22
adamantane c60 cubane hexamine
t c r c
gOctupole, VdW 2.72 3.49 2.48 1.09 0.18
0.19 0.12 1.00 0.16 0.45
ne ar kr xe
c c c c
VdW 5.89 7.93 7.18 5.49
TABLE 3: Deviations (%) from Experimental Lattice with KS Functionala dH NH
0.12 1.87 8.10 163 4 0.12 3.40 170 9 1.01 2.16 152 5
1.87 0.66 186 1 0.06 1.28 168 3 0.08 1.69 0.45 4.88 1.75 0.45 0.63 0.01 1.70 160 4 -0.22 0.79 0.76 202 2 0.48 1.21 0.67 177 1 -2.01 0.69 0.64 184 1 4.73 0.50 0.80 202 2 -1.74 0.40 1.47 155 4 -2.35 1.18 0.53 188 1 -0.82 1.89 176 5 -0.08 -0.32 1.30 178 3 -5.99 0.34 -7.02 0.88 176 2 -0.87 1.97 0.78 183 1 -2.97 -0.18 1.00 180 2 3.46 2.15 0.20 -1.23 1.07 195 4
compound
crystal
TTF-TCNQ histidine resolved cholamide dihydrate orotic acid
m o m a
acetanilide alanine aniline apha aspartic acid benzamide benzoic acid caprolactam cinnamide citric acid ethanol glucose glycine h2s ice imidazole mono paracet toluene urea
o o m m m m m m m m m o m o o m m m t
acetylene adipic acid anthracene benzene βsuccinic acid biphenyl indathrone indigo1 indigo2 naphtalene pigment red quinone sulfur tere acid1 tere acid2
c m m o m m m m m m a m o a a
adamantane c60 cubane hexamine
gOctupole, VdW t 0.86 c 1.39 r -2.32 c -1.58
ne ar kr xe
c c c c
149 2
169 2
198 2 211 2 198 2
159 2 176 2
0.01 0.02 0.04 0.05
a
For monoclinic and lower symmetry lattices, the deviation for the cube root of the cell volume is an additional independent lattice parameter of dimension length. Also shown is the calculated cohesive energy (eV) of the molecular subunits. The shortest hydrogen bond length dH (pm) in the structure according to experiment and the number NH of hydrogen bonds are shown. Column 2 identifies crystal types: c, cubic; t, tetragonal; o, orthorhombic; m, monoclinic; a, triclinic.
(a: anorthic) lattices, represented by four cases in the present data set, deviations in six lattice parameters are possible, but only four independent deviations are evaluated statistically here. The TTF-TCNQ crystal shows real charge transfer. A necessary consequence of charge transfer in continuous density functional theory is that the LUMO orbital of the anion is moving up in energy approaching the HOMO of the cation. In the case of TTF-TCNQ, the molecular levels are broadened into overlapping bands and the compound becomes a narrow band metal. This property is qualitatively consistent with the well-
∆A
∆B
∆C
∆L
U0
Monopole -0.51 -0.55 -1.50 -1.07 2.93 -2.23 -0.87 0.27 8.75 -0.86 0.34 -2.07 -0.95 4.72 3.01 2.03 -4.48 -1.07 2.84 Dipole 0.91 0.05 -2.53 -0.28 -0.67 0.01 -1.39 -3.51 -2.01 -0.92 -1.91 -2.21 -0.93 -0.22 -5.30 -3.80 -1.82 -2.70 -1.65
-4.22 -3.26 0.35 -2.36 -1.87 -1.06 -1.36 -1.15 -0.95 -0.12 -0.60 -0.81 -2.75 -6.73 -5.70 -1.83 -2.41 -1.31
Quadrupole -4.40 -0.46 2.44 -1.36 -0.34 0.46 -0.95 -1.06 -3.74 -3.65 -1.53 -0.15 -0.46 -1.49 -1.14 -0.45 0.21 0.18 0.68 0.14 0.09 -5.58 -1.43 -0.03 0.31 -0.73 -0.98 -0.27 1.13
VdW -8.76 -1.60 -1.05 1.27
-1.55 0.72 -1.81 0.08 -0.42 -4.13 -1.72 -1.67 -0.43 -2.82 -3.12 -2.46 0.73 -5.90 -5.81 -0.43 -2.96 -0.55 -1.07
4.82 0.24 -0.19 -0.25 0.20 0.35 -1.43 -0.79 0.91 -0.61 -0.39 -0.97 -4.04 -1.55 -1.26
-1.30 -0.85 -0.79 -1.74 -1.55 -1.91 -1.12 -1.43 -2.00 -0.96 -1.65 -1.84 -1.44
1.05 1.64 0.76 1.11 2.19 1.13 1.02 1.00 1.24 2.10 0.78 2.56 1.67 0.50 1.05 1.01 1.47 0.49 1.36
0.34 -1.22 1.67 -0.83 0.78 0.48 -1.69 1.58 -1.66 0.77 -0.07 1.55 -1.00 1.42 -0.25 1.42 0.27 0.62 -0.21 1.75 -2.42 0.72 0.76 -1.84 1.69 -0.35 1.69 0.56 0.62 -2.16 0.52 0.82 0.05 0.10 0.13 0.14
a
Also shown is the calculated cohesive energy (eV) of the molecular subunits.
known experimental observation of near metallic conductivity of TTF-TCNQ. In the case of histidine, charge transfer remains small enough that a significant band gap of several electronvolts persists. This may seem in contradiction to an expected formal charge of -1 for the Cl atoms. However, histidine is not a strong electron donor. So its HOMO gets charge equilibrated after a fractional electron transfer to Cl. Nevertheless, the crystal cohesive energy is due in a large part to the monopole electrostatics. The other two compounds in the monopole class contain crystal water as the second molecular species. The monopole
Molecular Crystals
J. Phys. Chem. C, Vol. 114, No. 48, 2010 20527
Figure 2. Distribution of deviations of lattice parameters (A, B, C, V1/3) for the classes of molecular crystals: dash-dot, monopole; full line, dipole; dashed, quadrupole; dotted, VdW-octupole. Panels: (a) functional PBEsol; (b) KS; (c) PBE; (d) PBE with small atomic radial cutoff.
terms remain small, and the crystal cohesive energy is less pronouncedly larger than for the other classes of molecular crystals considered here. The electrostatic dipole-dipole class covers the largest fraction of molecular crystals containing a single species. Mostly, the deviations of the lattice constants remain in the low percentage range, with a tendency toward too large lattice parameters. Notable exceptions in the present collection are the perfectly ordered model-ice crystal and the H2S crystal with clearly too short lattice parameters. The largest deviation toward too long lattice parameters is found for toluene. It is significant that the binding energy for this molecule with rather low dipole moment is the lowest in the dipole-dipole class. The next class, where the leading long distance interactions are down to static quadrupole-quadrupole type show lower calculated cohesive energies and a more pronounced tendency toward too long lattice parameters. In the case of the terephthalic acid crystals 1 and 2, cohesion may still be enhanced by hydrogen bonds. These are not shown in the table as they exceed our somewhat arbitrary cutoff at 300 pm. The class of molecules with symmetry so high that the van der Waals interaction dominates at long distance contains a very limited number of molecules, which are of lesser importance for most applications. In this class, the deviations grow to larger values. The case of hexamine is notable. Because of its polar groups, this molecule still produces a strong electrostatic near field, which helps to significantly smaller deviations than the other members of this class. For completeness, also the rare gas crystals are shown. That the deviations are largest here is no surprise, as only the inaccurately modeled van der Waals interaction provides cohesion. It is perhaps a surprise that the deviations remain well below 10%. The table also shows the bare cohesive energy per molecule without vibrational and classical free energy contributions. In the case of multi component crystals, the energy is per molecule of the lead species. In histidine, one of the additional species in the crystal is the Cl ion. In principle, one should refer to its
atomic ground state as its dissociated reference energy. Atomic orbital occupation according to Hund’s rule leads to open shell atoms with nonspherical density, which may exhibit the lowest possible DFT atomic energy. For Cl, however, no lowering occurs for the functionals considered here, and the spherically averaged atom energy is recovered as reference. The basic KS functional performs perhaps unexpectedly well for lattice parameters, as shown in Table 3. In view of Figure 1, it is no surprise that KS predicts larger cohesive energies than PBEsol. For the weakest bound crystals, KS yields more than twice the cohesive energy predicted by PBEsol. For the crystals bound by stronger electrostatics, the fine points of electron correlation matter less and the two functionals agree to better than 10% for crystal of histidine and to about 25% for crystals like glucose. Figure 2 summarizes the statistics for the nonredundant deviations ∆A, ∆B, ∆C, and ∆L shown in the tables for each of the classes of molecular crystals. The rare gas crystals were omitted and only the 4 high symmetry molecules, with 6 independent deviations, retained for the VdW-octupole+ class. Obviously, the data set tested here is small. To obtain smooth looking distribution curves and to deemphasize the specifics our small statistical sample, a Gaussian broadening with rms width of 1% was applied for the figure. For the present data, the KS functional shows a quite consistent systematic underestimation of length lattice parameters by 1% and best consistency of the scatter around the signed mean error. Without the artificial broadening used for the figure smoothing, the scatter of KS deviations around the mean deviation remains very consistently below 2% for the four molecular compound classes considered here. PBEsol shows somewhat wider distributions of deviations for molecular compounds than KS, but up to the quadrupole class, the systematic mean error remains small and on the overestimation side. Figure 2 shows also the distribution of deviations for the PBE functional, which was shown to be successful for predicting molecular enthalpies.4 Similar to its tendency toward too long bond distances, significant overesti-
20528
J. Phys. Chem. C, Vol. 114, No. 48, 2010
mation of length lattice parameters appears on average. This is consistent with the reported behavior of pure DFT9,12 with the PBE functional. Also the probability distribution of deviations gets significantly wider than for PBEsol and KS. While the main interest of this article is on limitations due to particular density functional approximations, it was noted that significantly better lattice parameter predictions result from the PBE functional in conjunction with the short local orbital cutoffs: Obviously, there is partial error compensation between the too large lattice constants from the PBE functional and the spurious attraction from using short radial cutoffs for the local orbital basis sets. This may be pragmatically considered to be a quantum mechanical model (PBE-sc), which is primarily defined here by the density functional plus the specific basis set. With the small cutoff radius, the basis set errors grow to similar magnitude, but opposite sign as the errors from the PBE functional. Predictions of the PBE-sc model for lattice parameters that are similar to PBEsol results, with notably a smaller width of the distribution of deviations for the dipole class and the weaker bound classes. The predicted cohesive energies from this model as estimated from the U0 values in Table 4 are significantly larger than the ones from PBE (not shown) and fall intermediate between PBEsol and KS. The shortening of the basis function tails enhanced the basis set superposition errors. Even for the significantly longer tails used for the main part of this study, a small amount of basis set superposition errors remain. This was investigated by counterpoise calculations for the case of an acetylene molecule in its calculated crystal environment with its nearest 40 neighbor atoms. For PBEsol 28 meV BSSE was found, to be compared to uncorrected 200 meV for U0. For KS, the shorter intermolecular distance leads to a larger BSSE of 40 meV for a U0 of 340 meV. For the PBE-sc model, the BSSE energy term is an order of magnitude larger, a consequence of the exponential dependence on cutoff radius. The short cutoff of the atomic tails introduces a systematic average reduction of length lattice parameters for molecular crystals by 2-3%. This is beneficial for application of the PBE functional. For our standard cutoffs, the BSSE effect falls below about 0.25%, as verified with test calculations with large cutoff values. The DMol3 numerical basis sets have been known for a long time for low BSSE because of the accurate evanescent tails of the molecular orbitals. Other sources of error for low cohesive energies need to be carefully controlled. A technically difficult case is C60 with only 0.12 eV cohesion per molecular unit on a background of total energy of more than 60 000 eV and molecular dissociation energy still exceeding 476 eV. In the case of C60, the coordinates of the buckyballs in the experimental lattice were replaced by gas phase optimized buckyballs superposed with the original experimental coordinates. In this way, the lattice optimization starts to work on the ever so small distortions from the gas phase structure arising from lattice interactions. Without this preparation, there is a risk that the solid state optimization stops early while relaxing the experimental structure. In the extreme, no cohesion with respect to the gas phase dissociation products appears in the calculation. The performance seen for weakly bonded crystals when using the PBEsol, KS functionals or the PBE functional with short cutoff is very similar to the performance found for crystals with covalent, metallic, or ionic bonds. Milman31 reports a deviation range from -3.5% to +4.5% for a collection involving metals and other inorganic compounds. A previous study by one of the present authors32 reports deviations in the ranges of -0.4%
Todorova and Delley TABLE 4: Deviations (%) from Experimental Lattice with PBE Functional and Using a DND Basis Set with Short Cutoff compound
crystal
TTF-TCNQ histidine resolved cholamide dihydrate orotic acid
m o m a
Monopole 0.95 0.48 -0.14 7.16
acetanilide alanine aniline apha aspartic acid benzamide benzoic acid caprolactam cinnamide citric acid ethanol glucose glycine h2s ice imidazole mono paracet toluene urea
o o m m m m m m m m m o m o o m m m t
Dipole 0.23 1.62 -1.42 -0.05 0.13 2.73 0.87 4.74 -0.31 2.25 2.94 -0.84 0.99 3.68 -2.19 -0.42 -0.26 0.30 0.25
acetylene adipic acid anthracene benzene βsuccinic acid biphenyl indathrone indigo1 indigo2 naphtalene pigment red quinone sulfur tere acid1 tere acid2
c m m o m m m m m m a m o a a
adamantane c60 cubane hexamine
gOctupole, VdW t 3.27 c 2.98 r 3.88 c 1.32
ne ar kr xe
c c c c
∆A
∆B
∆C
∆L
U0
3.40 1.45 1.66 2.20 0.89 -0.76 8.02 1.98 0.11 0.68 3.84 4.80 -3.39 2.24 2.14 -0.61 0.72 6.87 1.74 0.96 0.35 1.73 -0.31 1.83 1.42 2.39 0.49 0.51 -1.76 -4.74 -0.33 -0.16 0.32
3.10 0.81 -1.89 2.42 1.84 0.86 0.76 -2.17 2.00 0.51 -0.66 0.68 1.12 -1.53 -3.22 0.90 0.55 1.20 0.46
0.79 1.39 0.98 1.29 1.11 0.84 1.06 1.51 1.66 0.96 0.05 0.48 1.03
0.86 1.43 0.67 0.89 1.82 0.95 0.80 0.85 1.06 1.53 0.63 1.97 1.45 0.34 0.85 0.89 1.16 0.43 1.17
Quadrupole -2.55 0.34 0.52 1.19 1.53 0.93 1.31 1.83 2.34 0.73 1.57 0.64 2.49 1.95 1.51 0.41 5.19 -4.76 4.95 1.62 1.24 -0.11 0.76 0.57 0.41 0.62 2.51 4.24 0.27 2.34 1.13 2.27 -0.59 -0.06 1.48 1.08 5.53 -1.75 0.41 2.79 1.06 2.07 3.22 0.88 2.18 0.37 5.03 0.96 0.61 2.01 1.32 1.35 0.76 0.38 0.84 0.58 3.30 3.43 4.21 0.35 0.05 0.37 1.59 0.86 1.29 0.57 4.23 -0.22 2.14 1.28
VdW -0.26 1.51 1.63 5.10
4.32
0.51 0.25 2.73 0.27 0.77 0.03 0.08 0.12 0.12
to +1.8% for metals, +0.2% to +2.9% for semiconductors, and -2.0% to +3.1% for salts. Molecular compounds are much softer than the conventionally bonded crystals. So it is quite an achievement, if the unavoidable approximations to the density functional do not result in much larger errors for the lattice parameters. Table 5 summarizes the present results in comparison with a previous study32 covering some elemental solids, zinc blende structure semiconductors and rock-salt structure ionic solids. The average deviation (mean signed error) for the lattice parameters of dimension length suggests that both KS and PBE overestimate bond lengths in strongly bound solids, while PBEsol is closer. For molecular solids PBE has a systematic bias toward too long weak bonds, while PBEsol
Molecular Crystals
J. Phys. Chem. C, Vol. 114, No. 48, 2010 20529
TABLE 5: Average Percent Deviations from Experimental Lattice Parametersa compound type 32
elemental semiconductor32 ionic32 molecular: monopole molecular: dipole molecular: quadrupole molecular: VdW-octupole
N
PBEsol
KS
PBE
PBE-sc
45 23 21 15 69 55 6
-0.63 (1.5, 1.6) 0.86 (1.0, 1.0) 0.06 (1.1, 0.8) 0.43 (1.2, 0.9) 0.17 (2.9, 2.1) 1.54 (2.1, 1.9) 1.83 (1.3, 1.8)
0.61 (1.5, 1.8) 1.75 (1.0, 1.8) 0.34 (1.0, 0.8) -0.70 (1.7, 1.5) -1.78 (1.6, 1.9) -0.72 (1.6, 1.2) -0.85 (1.6, 1.6)
0.62 (1.5, 1.8) 2.42 (0.9, 2.4) 2.28 (1.2, 2.5) 3.94 (2.0, 3.9) 3.94 (3.3, 4.3) 5.03 (2.9, 5.1) 6.81 (1.8, 6.8)
1.20 (1.7, 2.2) 2.34 (1.0, 2.3) 2.61 (1.3, 2.9) 1.43 (2.4, 2.0) 0.69 (1.7, 1.4) 1.51 (1.9, 1.9) 3.08 (1.0, 3.5)
a Deviations are evaluated for the subset of lattice parameters A,B,C,L which are not symmetry redundant. Values are MSE (rms, MUE) in (%) relative to experimental values. MSE average signed deviation, rms scatter of deviations around the signed average deviation, MUE mean unsigned error with respect to the experimental parameters. N: number of independent lattice parameters. Results with conservative cutoff for PBEsol KS PBE and under PBE-sc results with short cutoff as shown in Table 4.
shows this overestimation only for the weakest bound classes. The KS functional shows a modest systematic underestimation of length dimensions for molecular solids. This is reinforced by the effect of the residual BSSE in the present study. An underestimation is also suggested in Figure 1. The mean unsigned deviations with respect to the experimental parameters is an overall measure for the accuracy of the computational approximations combined with the functional. The width of the distribution, summarized by the rms scatter around the mean remains pleasantly low for KS and PBEsol, at about 2 times the width for the established, strongly bound, standard crystals. IV. Summary and Conclusion Molecular crystals were classified by the leading long-range power law for the interaction among the molecular constituents. Charge transfer compounds, and more generally compounds with more than one constituent species have monopole interaction terms as the leading term of intermolecular interactions. Most molecular compounds have static dipole interactions as the leading long-range interaction term. Higher symmetry compounds have static quadrupole interactions. Very high symmetry molecules have van der Waals interaction as the leading order interaction for large separation. Some density functional approximations, including PBEsol, PBE, and the Kohn-Sham local exchange model, find the binding energy D0 of a pure van der Waals test case within a factor of 2. These functionals have been used to calculate lattice parameters for a test set of molecular crystals. PBEsol and KS yield lattice parameters of dimension length for almost all compound classes considered here with an average systematic error of less than 1%. The rms scatter around the systematic mean error for molecular compounds remains less than 3% for PBEsol and less than 2% for KS. This is about twice as much scatter as for nonmolecular compounds like metals, covalent semiconductors and inorganic salts. The mean unsigned error with respect to experimental parameters is about 2. The PBE functional systematically overestimates lattice parameters on average by about 4% for the stronger bound molecular crystal classes, and even more for the weaker bound ones. Scatter is slightly larger than for PBEsol. The present findings for PBE confirm earlier reports. It was also found that a PBE-sc model involving short local orbital cutoff radii can lead to improved calculated lattice parameters as compared to an accurate PBE calculation. The already quite good predictions for molecular lattices with the PBEsol and the KS functional could certainly be improved by addition of semiempirical van der Waals corrections with optimization of the unavoidable empirical rolloff parameters. An interesting open question is, could the much improved PBEsol be the basis for an even more successful method with semiempirical corrections?
The present predictions for lattice parameters of molecular crystals by a universally applicable nonempirical density functional yield about 2 times larger deviations than a density functional augmented by van der Waals terms including optimized semiempirical roll-off parameters. While semiempirical fit parameters are not at the present focus of interest, it might be interesting to explore in future work if PBEsol or KS with tuned dispersion contribution can improve over the state reached by Neumann and followup work.10-15 It is interesting that the widely usable PBEsol functional, an approximation arising from restoring the density gradient expansion for exchange in solids and surfaces, performs more favorably than PBE for the weak interactions in molecular solids. This is an encouraging hint to expect that relatively simple functional approximations may exist, which could lead to accurate predictions for molecular crystals while retaining universal applicability. Acknowledgment. Support by Swiss NSF grant 200021109358 is gratefully acknowledged. Supporting Information Available: Text file of data for the molecular crystals studied here: name of compound; our input crystal structure data with Cartesian coordinates for the translation vectors and all atoms in the cell [Bohr], experimental lattice parameters [Å] and degrees; calculated lattice parameters for KS; percent deviations for length and (Vcalc/Vexp)1/3; singular value analysis of cell transformation exp f KS with singular values presented as percent deviation from unity; ∆ (ref 12). This material is available free of charge via the Internet at http:// pubs.acs.org. References and Notes (1) Byrn, S. R. Pfeiffer, R. R. Stowell, J. G. Solid State Chemistry of Drugs, 2nd ed.; SSCI Inc.: West Lafayette, IN, U.S., 1999. (2) Schmidt, G. M. J. Pure Appl. Chem. 1971, 27, 647. (3) Kohn, W.; Sham, L. J. Phys. ReV. A 1965, 34, 737. (4) Delley, B. J. Phys. Chem. A 2006, 110, 13632. (5) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. ReV. Lett. 1996, 77, 3865. (6) Perdew, J. P.; Ruzsinszky, A.; Csonka, G. I.; Vydrov, O. A.; Scuseria, G. E.; Constantin, L. A.; Zhou, X. L.; Burke, K. Phys. ReV. Lett. 2008, 100, 136406. (7) Zhang, Y.; Pan, W.; Yang, W. J. Chem. Phys. 1997, 107, 7921. (8) Elstner, M.; Hobza, P.; Frauenheim, T.; Suhai, S.; Kaxiras, E. J. Chem. Phys. 2001, 114, 5149. (9) Wu, Q.; Yang, W. J. Chem. Phys. 2002, 116, 515. (10) Grimme, S. J. Comput. Chem. 2004, 25, 1463. 2004. (11) Zhechkov, L.; Heine, T.; Patchokovskii, S.; Seifert, G.; Duarte, H. A. J. Chem. Theory Comput. 2005, 1, 841. (12) Neumann, M. A.; Perrin, M. A. J. Phys. Chem. B 2005, 109, 15531. (13) Asmadi, A.; Neumann, M. A.; Kendrick, J.; Girard, P.; Perring, M. A.; Leusen, F. J. J. J. Phys. Chem. B 2009, 113, 16303. (14) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. J. Chem. Phys. 2010, 132, 154104.
20530
J. Phys. Chem. C, Vol. 114, No. 48, 2010
(15) Sorescu, D. C.; Rice, B. M. J. Phys. Chem. B 2010, 114, 6734. (16) Rapcewicz, K.; Ashcroft, N. W. Phys. ReV. B 1991, 44, 4032. (17) Andersson, Y.; Langreth, D. C.; Lundqvist, B. I. Phys. ReV. Lett. 1996, 76, 102. (18) Dobson, J. F.; Dinte, B. P. Chem. Phys. Lett. 1996, 76, 1780. (19) Kohn, W.; Meir, Y.; Makarov, D. E. Phys. ReV. Lett. 1998, 80, 4153. (20) Delley, B. An All-Electron Numerical Method for Solving the Local Density Functional for Polyatomic Molecules. J. Chem. Phys. 1990, 92, 508. (21) Delley, B. From molecules to solids with the DMol3 approach. J. Chem. Phys. 2000, 113, 7756. (22) Halgren, T. A.; Lipscomb, W. N. Chem. Phys. Lett. 1977, 49, 225. (23) Andzelm, J.; King-Smith, R. D.; Fitzgerald, G. Chem. Phys. Lett. 2001, 325, 321. (24) Baker, J.; Kessi, A.; Delley, B. J. Chem. Phys. 1996, 105, 192.
Todorova and Delley (25) Accelrys Materials Studio database. (26) Cox, E. G.; Cruickshank, D. W. J.; Smith, J. A. S. Proc. R. Soc. 1958, 247, 1–21. (27) Anderson, M.; Bosio, L.; Bruneauspouille, J.; Fourme, R. J. Chim. Phys. 1977, 74, 68–73. (28) Jo¨nsson, P. G. Acta Crystallogr. 1976, B32, 232–235. (29) Van Nes, G. J. H.; van Bolhuis, F. Acta Crystallogr. 1979, B35, 2580–2593. (30) David, W. I. F.; Ibberson, R. M.; Matthewman, K. P.; Dennis, T. J. S.; Hare, J. P.; Kroto, H. W.; Taylor, R.; Walton, D. R. M. Nature 1991, 353, 147. (31) Milman, V.; Winkler, B.; White, J. A.; Pickard, C. J.; Payne, M. C.; Akhmatskaya, E. V.; Nobes, R. H. Int. J. Quantum Chem. 2000, 77, 895. (32) Delley, B. Phys. ReV. B 2002, 66, 155125.
JP1049759