The d-DFT energy minimizations were carried out with the computer program GRACE
, which uses the computer program VASP
for single-point pure DFT calculations. GRACE
implements an efficient minimization algorithm to reduce the number of expensive single-point DFT calculations, and GRACE
augments the pure DFT energies with a dispersion correction from hybridization-dependent isotropic atom–atom potentials. The details are given in Neumann & Perrin (2005
); we mention here only that we use the Perdew–Wang-91 functional and a plane-wave energy cut-off of 520 eV. The dispersion-correction parameters for iodine were kindly provided by Dr J. Kendrick of the Institute for Pharmaceutical Innovation in Bradford; the dispersion-correction parameters for boron and bromine came from in-house parameterizations.1
All dispersion-correction parameters were parameterized against low-temperature (2–130 K) crystal structures and the d-DFT method was intended to reproduce unit-cell parameters at essentially 0 K. No dispersion-correction parameters were available for charged atoms: the parameters of the corresponding neutral species were used. The convergence criteria for the minimization were < 0.003 Å for the maximum Cartesian displacement (including H atoms), < 2.93 kJ mol−1
for the maximum force and < 0.00104 kJ mol−1
per atom for the energy difference between the last two minimization steps.
The energy optimizations were divided into two steps: first an energy optimization with the unit cell fixed, followed by a second step with the unit cell free, starting from the energy-minimized crystal structure from the first step
This two-step procedure has a computational advantage. From a numerical perspective, the energy of certain strong interactions such as chemical bonds is very sensitive to the atomic positions and small experimental uncertainties can result in large initial forces. At the beginning of the minimization procedure, when the optimization algorithm has no or only approximate information about the anisotropy of the curvature of the potential energy hypersurface, such forces can result in a large step in the wrong direction, ultimately leading to the structure getting trapped in a less favourable side minimum. The robustness of the minimization procedure is improved if ‘hard’ degrees of freedom, in practice intramolecular degrees of freedom, are minimized first. With respect to separating hard and soft degrees of freedom, the above scheme is not perfect since in the first minimization the soft molecular translations and rotations are adjusted together with the intramolecular degrees of freedom. In fact, in order to avoid getting trapped in a side minimum, for one crystal structure it turned out to be necessary to apply a three-step optimization procedure, with the unit-cell parameters, the molecular positions and the molecular orientations being held fixed for the first minimization
The three-step procedure requires more CPU time than the two-step procedure, and the three-step procedure should only be used if there are reasons to suspect that the crystal structure may have ended up in a side minimum.
Since pure DFT optimizations, without dispersion correction, are common in the crystallographic literature, almost invariably with the experimental unit cell kept fixed during the optimization, the calculations were repeated with pure DFT with the experimental unit cell kept fixed for comparison. The pure DFT calculations with fixed unit cell were carried out merely to reassure other authors that such calculations are indeed meaningful, and these calculation will only be mentioned briefly as part of the discussion.
Unless otherwise indicated, the experimental space-group symmetry was used, which imposes certain constraints on unit-cell parameters, atomic positions and Z.
For validation, two test sets will be used:
- (i) a test set of crystal structures that can be assumed to be correct;
- (ii) a test set of crystal structures that are known to be ambiguous or wrong.
For a test set of correct crystal structures, all 249 organic crystal structures from the August 2008 issue of Acta Cryst. Section E were downloaded, with permission. Acta Cryst. Section E is an open access journal, making the test set publicly available to all. Two crystal structures contained silicon and six contained phosphorus, two elements for which the dispersion correction has not yet been parameterized. These eight crystal structures had to be omitted from the test set, leaving 241 crystal structures. These 241 crystal structures cover a wide spectrum of molecular crystal structures including sugars, a high-energy material, drug molecules, chiral molecules, disordered structures, hydrates, solvates, salts and a range of space groups, functional groups and elements (C, H, B, Br, Cl, F, I, N, O and S). Three crystal structures are polymorphs of earlier determinations, but the test set contains no pairs of polymorphs. There were 16 disordered crystal structures which had to be adjusted manually before minimization. These disordered structures were not included in the validation set and will be discussed separately, leaving 225 crystal structures for the validation set.
For the test set of incorrect crystal structures, we took four structures that were known to be wrong. Two were from the literature (examples 1 and 6) and two turned up among the 225 structures in the Acta Cryst.
Section E test set (examples 8 and 9). Four more crystal structures were added as examples where structure solution from powder diffraction data had yielded ambiguous results (examples 3, 4, 5 and 7). These structures require individual discussion and they are described below in §4
Each crystal structure was energy-optimized in two ways: with the experimental unit-cell parameters kept fixed and with the unit cell allowed to vary. This provides us with a set of 225 times three crystal structures: the experimental crystal structure plus the two optimized structures. By comparing any two out of those three crystal structures and calculating the volume difference, the energy difference, the r.m.s. or the maximum Cartesian displacement with or without H atoms etc., a large number of possible quality measures can be calculated. Moreover, two quality measures can be plotted against each other to generate two-dimensional scatterplots, quadratically increasing the number of plots. Several one-dimensional quality measures were explored in some detail, but one turned out to be the most relevant one for the purpose of discriminating between correct and incorrect crystal structures: the r.m.s. Cartesian displacement between the experimental crystal structure and the fully optimized crystal structure (including unit cell), excluding H atoms.
‘Cartesian displacement’ is not uniquely defined when the unit cells of the two crystal structures to be compared are different, as is the case when we compare the experimental crystal structure to the d-DFT optimized structure with the unit cell allowed to vary. In this work the Cartesian displacement for an atom in two crystal structures (1) and (2) is
are the fractional coordinates of the atoms in crystal structure i
, and G
is the transformation matrix from fractional to Cartesian coordinates for crystal structure i
. This definition of Cartesian displacement has the advantages that it is symmetric with respect to the two structures to be compared, that it varies smoothly upon smooth distortions of either or both of the two structures to be compared, and that there is no need for a user-defined parameter such as the number of molecules used for the comparison.