Toward a Quantitative Description of Crystal ... - ACS Publications

L. Holland , Arnold L. Rheingold , Curtis E. Moore , and Joseph M. O'Connor .... Alistair R. Lennie , W. D. Samuel Motherwell , Simon Parsons , El...
0 downloads 0 Views 310KB Size
Toward a Quantitative Description of Crystal Packing in Terms of Molecular Pairs: Application to the Hexamorphic Crystal System, 5-Methyl-2-[(2-nitrophenyl)amino]-3-thiophenecarbonitrile†

CRYSTAL GROWTH & DESIGN 2005 VOL. 5, NO. 6 2180-2189

J. D. Dunitz*,# and A. Gavezzotti*,§ Chemistry Department OCL, ETH-Ho¨ nggerberg HCI H333, Zurich, Switzerland, and Dipartimento di Chimica Strutturale e Stereochimica Inorganica, University of Milano, Via Venezian 21, 20133 Milano, Italy Received March 18, 2005;

Revised Manuscript Received June 28, 2005

ABSTRACT: Elements of a theory of crystal packing are presented in the form of a systematic analysis of crystal packing molecular pairs, i.e., neighboring molecule-molecule pairs in the crystal, rather than in terms of selected intermolecular atom-atom contacts. Intermolecular energies, based on the molecular electron density distribution and partitioned over Coulombic, polarization, dispersion, and repulsion contributions, are calculated for such pairs by the Pixel-semiclassical density sums (SCDS) method, recently updated for a better treatment of the dependence between electron density overlap and repulsion energy. The advantages of the pairs treatment are illustrated by a study of the six known polymorphs of 5-methyl-2-[(2-nitrophenyl)amino]-3-thiophenecarbonitrile, in which diffuse interactions between aromatic rings over a wide range of stacking modes are shown to contribute more to the cohesive energy than do lateral interactions between polar moieties, often described as weak hydrogen bonds. Introduction As the title of this journal might encourage us to believe, the analysis, prediction, and control of crystal architecture is a promising and even productive area of contemporary science. As far as organic compounds are concerned, the prediction of crystal structure, given the molecular structure,1,2 has been the subject of a series of three blind tests, organized by the Cambridge Crystallographic Data Centre,3-5 in which many leading research groups have participated. However, only modest success has been achieved so far.6,7 The main problem seems to be not so much a matter of generating stable crystal structures but rather that of selecting one or more possible structures from very many almost equienergetic candidates. At a more descriptive level, following a long and fruitful chemical approach, the analysis of crystal structures is usually carried out in terms of intermolecular atom-atom contacts, interpreted as intermolecular bonds or, more circumspectly, as bonding interactions. However, the problem is complicated, because packing forces are weaksfar weaker than standard chemical bondssand scarcely selective, so that any proposed partitioning of the lattice energy into discrete intermolecular atom-atom interactions must be controversial, to say the least.8 The apparent analogy between such weak intermolecular atom-atom bonds in crystals and the many times stronger chemical bonds that hold molecules together may be appealing, but it can be misleading. Rather than focusing on interactions between atoms in neighboring molecules, we utilize a procedure based on the evaluation of the interaction energy between †

For Michael McBride, scientist, scholar, teacher, in long friendship. * To whom correspondence should be addressed. E-mail: dunitz@ org.chem.ethz.ch (J.D.D.); [email protected] (A.G.). # ETH, Zurich. § University of Milano.

pairs of entire molecules. Quantum mechanical calculations of intermolecular energies, including a full account of electron correlation, are still problematic, partly because of matters of principle, and also because of extreme computational demands. In our simplified approach, intermolecular energies are calculated using a method that relies on integration over the molecular charge densitysthe Pixel-semiclassical density sums (SCDS) method,9,10 successfully applied to molecular crystals.11,12 To be sure, every quantum mechanical approximation creates its own energy partitioning, and the Pixel method yields a partitioning into Coulombic, polarization-dispersion and repulsion contributions. This partitioning, even if to a certain extent arbitrary and in large part (but not entirely) parameter dependent, as discussed below, is a key element of the theoretical construction, being readily applicable to organic crystals made of molecules of substantial size and having the additional advantage that individual energy contributions are easily interpretable in terms of familiar physicochemical concepts that are useful in the comparison of different packing modes for a given compound. Pixel total cohesive energies for many molecular dimers have been shown to be reliable at least at the same level as MP2-type quantum mechanical calculations.13 As a specific example, we take up the analysis of a hexamorphic crystal system, that of 5-methyl-2-[(2nitrophenyl)amino]-3-thiophenecarbonitrile, 1. The crystallography of 1 has been extensively studied.14 The compound has been crystallized as six solvent-free polymorphs, designated Y, YN, ON, OP, R, and ORP from the initials of the crystal color (Y, yellow; O, orange; R, red) and the crystal shape (N, needles; P, plates). All six polymorphs can coexist at room temperature, and all have one molecule in the asymmetric unit, with different conformations produced by rotation around

10.1021/cg050098z CCC: $30.25 © 2005 American Chemical Society Published on Web 09/01/2005

Quantitative Description of Crystal Packing

the two single bonds joining the central NH group to the aromatic rings, leading to the different colors. We thus have six different crystal packing arrangements for the same molecule, with relative stabilities estimated from solid-state conversions and calorimetric data.14 This experimental crystal structure landscape thus contains a wide selection of possible approach modes accessible to the constituent molecules. We calculate the lattice energies of the six polymorphs and also the binding energies of pairs of neighboring molecules in these structures and ask if any conclusions can be drawn about what kinds of interactions are most effective in making molecules stick together. Scheme 1a

a

N, green, O, red; S, yellow.

Computational Methods The Pixel-SCDS formulation9,10 requires first a calculation of the valence molecular electron density, here carried out for each of the six polymorphs of 1, based on the molecular geometries as found in the crystal structures, retrieved from the Cambridge Structural Database (refcodes are QAXMEH for ON, QAXMEH01 for Y, QAXMEH02 for R, QAXMEH03 for OP, QAXMEH04 for YN, and QAXMEH05 for ORP). Hydrogen atom positions were adjusted to standard geometries (C-H, 1.08 Å; N-H, 1.00 Å) since calculations and comparisons based on X-ray determined H-atom positions are notoriously unreliable. The electron density calculations were carried out at the MP2/6-31G** level with the Gaussian15 package; in this calculation, generous box limits and a step of 0.08 Å were used. The density was then condensed into 4 × 4 × 4 (condensation level n ) 4) super-pixels. Those with a charge smaller than 10-6 electrons were discarded, and the remaining pixel contents were renormalized to balance the sum of the nuclear charges. After this screenout, the number of remaining pixels per molecule of 1 was about 20 000, thus representing a delocalized description of the electron distribution. For the formation of dimers or molecular clusters in the various crystal structures, the nuclear positions and the electron densities are repeated in space according to the desired dimer geometry or to crystal symmetry operations, generating a collection of juxtaposed, undeformed, free-molecule electron densities. The Coulombic interaction energy (ECOUL) between any two molecules is then calculated by direct summation over pixel-pixel, pixel-nucleus, and nucleus-nucleus Coulombic contributions. The other terms are estimated along similar lines to those described previously11,12 but with a few minor adjustments suggested by further experience.13 As

Crystal Growth & Design, Vol. 5, No. 6, 2005 2181 Table 1. Atomic Polarizabilities, Atomic Radii, and Pauling Electronegativities atom

polarizability Å3 a

atomic radius Åb

electronegativity

Cl C (unsaturated) C (methyl) N O F H

2.30 1.35 1.05 0.95 0.75 0.50 0.39

1.76 1.77 1.77 1.64 1.58 1.46 1.10

3.0 2.5 2.5 3.0 3.5 4.0 2.1

a A selection with minor adjustments from values given in Miller, K. J. J. Am. Chem. Soc. 1990, 112, 8533. b Rowland, R. S.; Taylor, R. J. Phys. Chem. 1996, 100, 7384.

before, the polarization energy (EPOL) is estimated by calculating the electric field at each pixel and applying the linear polarization assumption. For this purpose, a local polarizability needs to be allotted to each pixel. Of course, there is no rigorous way of doing this. The procedure adopted is to assign each pixel to a particular atom and take the pixel polarizability to be proportional to the associated atom polarizability, scaled as the ratio of pixel charge to nuclear charge. This assumption may seem paradoxical in that it assigns lower polarizabilities to outer electrons rather than to inner electrons. However, although charge elements away from the nucleus are less restrained by the positive charge, they generate smaller dipoles for a given field because of their smaller charge. In any case, while the polarization distribution has admittedly9,10 a questionable physical basis, our procedure does at least ensure that the integrated volume polarizability of each atom is correctly reproduced. In principle, the assignment of pixels to atoms could be done by using atomic basins as defined in the atomsin-molecules description,16 but this would require a topological analysis of the electron density distribution, a considerable effort for a minor gain since it has been checked that the results are scarcely sensitive to the details of pixel allotment to atoms. We use a simpler procedure for assigning electron charge pixels to particular atoms, which differs slightly from the one described earlier.11 Let p be the number of atoms for which the nucleus-pixel distance is smaller than the atomic radius (these atomic radii were as previously described11). If p ) 1, the charge pixel is within one atomic sphere only, and it is assigned to that atom. If p > 1, the pixel is assigned to the atom from which the distance is the smallest fraction of the atomic radius. If p ) 0, the pixel is assigned to the atom to whose atomic surface it is closest. Table 1 shows the latest list of atomic polarizabilities, which have been slightly readjusted for a more systematic assignment. The charge distribution is not readjusted after polarization, so static polarization energies only are calculated. Note that the polarization energy is a many-body contribution, while the other terms are two-body ones. Dispersion energies (EDISP) are evaluated by a Londontype formula, again as a summation over pixel pairs. The repulsion energy is estimated from the electron density overlap. The overlap integral between the electron densities of any two molecules, SAB, is calculated by numerical integration. A new modification has been introduced to increase the flexibility of the approximation used to calculate repulsion energies from

2182

Crystal Growth & Design, Vol. 5, No. 6, 2005

overlap, for a better description of the wide variety of chemical environments found in organic crystals, ranging from pure hydrophobic interactions to hydrogen bonds. The total overlap is subdivided into contributions from the various atomic species, this having been made possible by the assignment of each electron density pixel to one of the atoms in the molecule as described above. For species i and j, the repulsion energy is evaluated as EREP,ij ) (K1 - K2 ∆pij)SAB, where ∆pij is the difference in Pauling electronegativity between atomic species i and j. The total molecule-molecule repulsion is then a sum over all these contributions from pairs of atomic species. The numerical values of the two adjustable parameters K1 and K2 have been optimized at 4800 and 1200, respectively (for energies in kJ mol-1 with pixel densities in electrons Å-3). The rationale for this approximation depends on the principle that any reorganization of the straightforward overlap density of the separate molecular electron density distributions is always stabilizing. The greater the electronegativity difference between the overlapping atoms, the greater should be the reorganization of the overlap electron density and also the greater the concomitant stabilization, thus lowering the Pauli repulsion energy contribution associated with the simple overlap model and resulting in a smaller repulsion energy. The total repulsion energy is a sum of molecule-molecule terms, although in some extreme cases there might be simultaneous overlap of the densities of more than two molecules. Undoubtedly, this revised formulation of the Pixel method includes more elaborate partitionings that depend on the nature of the atoms involved. The benefits appear to be worthwhile, however, judging from a comparison between calculated lattice energies and sublimation enthalpies for 91 organic crystals, and between Pixel and quantum mechanical results for 55 molecular dimers representing the most common recognition modes in organic crystals.13 On a standard personal computer, a calculation of the interaction energy for a molecular dimer takes about a minute, and a lattice energy calculation takes about half an hour. Parallelization of the computer code, presently under way,17 suggests the possibility of a 100-fold speed increase. The Significance of Electron Energy Partitioning. Any partitioning of the intermolecular energy is to some extent arbitrary; each method defines its own partitioned energies. In terms of basic physics, any chemical bonding, either intra- or intermolecular, is the result of electrostatic Feynman forces at atomic nuclei, these forces being zero when the whole nuclear and electronic structure has reached its equilibrium configuration. Thus, in principle, if the correct electron density were available, a purely Coulombic calculation would suffice; in practice, such a calculation is not feasible for an extended system because an extremely accurate wave function and extensive optimization would be needed. As a compromise, the Pixel approach uses an accurate electron density for the individual molecules and then approximates the interaction energy due to reorganization of the electron density in the extended system by introducing separate polarization, dispersion, charge transfer, and repulsion terms. The Pixel Cou-

Dunitz and Gavezzotti

lombic terms are parameter-free, and the results of our integration procedure are quite similar to those obtained by much more sophisticated methods.18 Since the electron distribution is delocalized over thousands of pixels, Pixel Coulombic energies include penetration effects, a crucial ingredient of molecular interaction at close contact, which, properly considered, may overturn many popular ideas about the details of intermolecular stabilization at a short distance. For the rest, Pixel calculations give a clear and realistic representation of the relative importance of the partitioned terms over different chemical categories: for example, dispersion energies are large in contacts among hydrocarbon molecules, especially polarizable, aromatic ones, and are less important in chemical environments where contacts between atomic species of different electronegativity are expected (e.g., anhydrides, quinones, hydrogen-bonding systems). Polarization energies are important in hydrogen-bonded systems, where they compensate, so to speak, for neglect of the electron density reorganization (partial covalent character of strong hydrogen bonds). Polarization energies are always stabilizing, even when Coulombic energies are destabilizing. The repulsion term is large whenever electron clouds overlap significantly, as a partial surrogate for the missing effects of Pauli exclusion, with ensuing Coulombic repulsion due to exposure of nuclei. Crystal Packing Molecular Pair Analysis. In a crystal, some chemical unit repeats itself with translational periodic symmetry. Here, for simplicity, we consider only the case in which the repeating unit consists of one entire molecule, but extension to mixed crystals, ionic crystals, solvates, and crystals with more than one molecule in the asymmetric unit is straightforward. The entire crystal is represented in terms of the positions of the atomic nuclei and the position and charge of the electron charge pixels of one molecule, the reference molecule, together with the cell parameters and crystal symmetry operations relating this molecule to all the others in the crystal. In an extension of a previous description,19 we here define a molecular pair as a pair of molecules in the crystal, the reference molecule plus another, usually a neighbor, characterized by the symmetry operator that relates the two molecules, the distance between their centers of mass, and the Coulombic, polarization, dispersion, repulsion, and total interaction energies (ECOUL, EPOL, EDISP, EREP, ETOT) of the pair. Some symmetry operators (center of symmetry, dyad axis, mirror plane) generate only a dimer, while others (those with translation components, such as screw axes and glide planes) generate molecular ribbons. In the latter case, the reference molecule and only one of its two equivalent nearest neighbors along the ribbon are considered; the total contribution to the energy of the crystal is then twice the cohesion energy of such a pair (assuming no many-body effects along the ribbon). The polarization energy in a crystal is a many-body term, and cannot, in principle, be partitioned over molecular pairs. Nevertheless, in this approach, we calculate the polarization energy of a single pair of molecules as if it were an isolated gas-phase dimer. As the distance between paired molecules increases, the pair energy tends to decrease. For our purpose, it is sufficient

Quantitative Description of Crystal Packing

to consider only neighbor pairs, i.e., those that constitute the coordination sphere of the reference molecule. Molecular pairs found in a given crystal structure may or may not correspond to energy minima for the free dimers because the molecular arrangement in the crystal results from the cooperative action of forces exercised by all the molecules in the coordination sphere. Of course, energies of molecular pairs also can be calculated by other methods, for example, by atomatom or point-charge force fields or by quantum mechanical methods. A few comparisons between Pixel results and those obtained by other methods will be presented. In the molecular pair type of analysis, it is convenient to arrange the various molecular pairs in a crystal structure in descending order of energy. As more and more molecular pairs are considered, the sum of the pair energies, with allowance for the multiplicity of the symmetry operators, approaches the lattice energy of the crystal but usually only the first three or four types of molecular pairs provide significant stabilizing energies, the rest merging into a continuum of small stabilizing and also sometimes slightly destabilizing contributions. The top contributors may be regarded as the fundamental structural units for the analysis of a crystalline aggregate, in terms of their cohesive energies. In the molecular pair type of analysis, all intermolecular interactions in the crystal are analyzed systematically and objectively, without any a priori assumptions about their importance. In the energy partitioning scheme, polar (Coulombic plus polarization) and nonpolar (dispersion) terms can be related to molecular electron distributions and polaritiessthat is, with chemistrysin what promises to reach a better understanding of crystal structure and hence improve the possibilities of crystal structure prediction and design. Of course, a molecular packing pair analysis is not limited to actual crystal structures. It is easily extendable to virtual crystal structures. With present computational techniques, it is relatively straightforward to generate a few hundred possible crystal structures in the most common space groups for any given compound.3-5 Such a collection of possible crystal structures can be expected to contain many if not all of the favorable molecular pairs, the low-energy structural association modes, available for a given molecular constitution, even though the task becomes more difficult with increasing degrees of torsional freedom. The estimated lattice energies of the structures in the sample need not be highly accurate and can therefore be calculated by any empirical, inexpensive force field. What one learns from a molecular pair analysis of a collection of real and virtual structures is an inventory of which geometrical motifs are unlikely (as they do not appear in any of the stable crystal structures), which motifs are common to many structures (and hence cannot be discriminating), and which motifs are obligatory (since no low-energy, close-packed crystal structure without them appears in the collection). Such a panorama of virtual crystal structures tells us much more than can be learned from any single structure and may warn us against assuming that an association mode found in some particular crystal structure, especially

Crystal Growth & Design, Vol. 5, No. 6, 2005 2183

an experimentally observed one, represents a specially favorable association of molecules. Results Illustrations of the crystal structures of the six polymorphs of 1 are shown in Figure 1. It is not obvious from the figure whether or how the structures are related, nor is it evident that they are roughly equienergetic. There might appear to be a vague similarity between the packings in the ON and YN polymorphs, but the molecular conformations are clearly different, the nitro and nitrile substituents being respectively syn and anti. There is some similarity between the OP and ORP polymorphs, but, again, the molecular conformations are not the same. Indeed, the figure tells us very little about the relationships among the structures. Table 2 shows the calculated lattice energies for the six polymorphs according to the Pixel method and also from the parametric UNI atom-atom force field.20 For all six polymorphs, by far the largest contribution to crystal stabilization comes from dispersion. Polymorph Y is calculated to have the lowest lattice energy, in agreement with the thermodynamic measurements of Yu et al.14 This preference is not altered when conformational energy differences are taken into account (Table 3). The energy differences are relatively small and could easily change by a few kJ mol-1, depending on minor changes in the positioning of the hydrogen atoms or minor differences in bonding geometry due to experimental inaccuracies. As an alternative, optimization of all molecular geometries within the same theoretical framework could be carried out, but this would undoubtedly mask at least a good part of the conformational changes occurring upon crystal packing. In any case, there does not seem to be any obvious correlation between molecular energy and the angle between the two rings or the intramolecular hydrogen-bonding distance. For the same reasons, lattice energy differences could change slightly with minor changes in intermolecular bond geometry, especially in the position of the hydrogenbonding hydrogen atom. These considerations indicate that caution is called for in taking the small energy differences among polymorphs as a strict ranking of structural stability. Indeed, this counsel applies to practically all computational studies of crystal polymorphism based on molecular simulation of crystal energies. In this particular case, we do not believe that a detailed discussion of relative stabilities in comparison with experiments is warranted, given the small values of the implied differences and the lack of temperature-dependent terms in the simulation methods. Table 2 shows that although dispersion provides the largest contributions to lattice energies, the extra stability of polymorph Y arises mainly from its favorable Coulombic energy. Yu et al. ascribed major importance to the structure-determining role of molecular dipole interactions in the crystal structures of 1. Our analysis does not support this interpretation. In fact, the Y polymorph, constructed from a nearly perpendicular, low-dipole conformer, has the largest Coulombic energy contribution to its lattice energy. We emphasize once again that at intermolecular distances comparable to or smaller than the dimensions of the molecules concerned, energies estimated from the point-dipole ap-

2184

Crystal Growth & Design, Vol. 5, No. 6, 2005

Dunitz and Gavezzotti

Figure 1. Projections of the six crystal structures down selected axes (left to right, top to bottom: Y, ON, YN, OP, R, ORP). There are no very striking similarities among them. Table 2. Partitioned Packing Energies and Unit Cell Data14 of the Six Polymorphs of 1a polymorph

Ecoul

Epol

Edisp

Erep

Etot

Eatom-atom20

density9 g cc-1

Y ON YN OP R ORP

-52 -46 -49 -47 -45 -44

-21 -19 -17 -18 -15 -14

-138 -136 -136 -132 -130 -124

86 84 87 84 77 69

-125 -117 -113 -113 -113 -112

-160 -154 -144 -151 -150 -149

1.447 1.428 1.431 1.435 1.438 1.429

Y ON YN OP R ORP a

a

b

c

8.500 3.945 4.592 7.976 7.492 13.177

16.413 18.685 11.249 13.319 7.790 8.021

8.537 16.395 12.315 11.676 11.911 22.801

R

β

71.1

91.8 93.8 89.9 104.7 77.8

75.4

γ

88.2 63.6

space group P21/n P21/c P1 h P21/n P1 h Pbca

Energies in kJ mol-1, cell parameters in Å and (°).

proximation are meaningless. Our subsequent analysis (see below) shows that the higher stability of the Y polymorph is not to be associated with the presence of any particular structural feature, not even the hydrogen bond, but rather with its overall and not easily dissectable Coulombic energy. From the most stable to the least stable polymorph, the Pixel energy difference is only 13 kJ mol-1 (roughly the same from the UNI calculation). In the total attractive energy (sum of Coulombic, polarization, and dispersion), the difference is 29 kJ mol-1, but this is countered by a reduction of 17 kJ mol-1 in the repulsion

Table 3. Relative Energies of the Molecular Conformations and of the Crystals of the Six Polymorphs of 1 polymorph ∆E (lattice)a ∆E (conf)b ∆E (total) R (NH‚‚‚ON)c φd Y R ON OP YN ORP

0 +12 +8 +12 +12 +13

0 -5 0 +4 +6 +9

0 +7 +8 +16 +18 +22

1.90 1.87 1.80 1.98 1.88 1.81

75 22 53 46 76 39

a Relative lattice energies, see Table 2. b Relative conformational energies, MP2/6-31G** calculation based on the conformation present in the crystal. c Intramolecular hydrogen bond distance. d Absolute values (°), standardized < 90°, of the angle between the mean molecular planes of the nitrophenyl and thiophenecarbonitrile rings. Energies in kJ mol-1, distances in Å.

energy, which is of course smaller for the less tightly packed polymorphs. The lattice energy calculations confirm that molecule 1 finds many sites of almost equal stability in its crystal energy hypersurface thanks to its conformational flexibility and its various chemical functionalities. The energy differences among the polymorphs come from complex factors, delocalized over the whole molecular electron distribution. We now proceed to the analysis of the crystal packing of the six polymorphs in terms of molecular packing pairs (Table 4). In assigning a structural designation to these pairs, it is convenient to refer to relevant molecular regions, such as the nitrophenyl (NP) and thiophenecarbonitrile (TC) groups, hydrogen-bonding sites, and so on. This helps to describe and classify them.

Quantitative Description of Crystal Packing

Crystal Growth & Design, Vol. 5, No. 6, 2005 2185

Table 4. Interaction Energies and Interaction Types for the Most Important Molecular Pairs in the Six Polymorphs of 1a determinantb

Ecoul

Epol

Edisp

Erep

Etot

R (coc) and symmetryc

approximate interaction typed

YA YN A OP A ORP A OP B ON A RA RB RC YB ORP B RD YN B YN C YN D ORP C RE OP C YN E YC YD YE

-14 -10 -31 -10 -12 -4 -14 -9 -9 -11 -3 -15 -14 -22 -18 -22 -21 -6 -14 -1 -16 -3

-6 -6 -11 -4 -4 -8 -3 -5 -4 -4 -4 -5 -4 -7 -8 -6 -6 -3 -3 -3 -8 -2

-58 -61 -27 -43 -41 -76 -34 -44 -44 -33 -41 -22 -25 -20 -27 -11 -9 -29 -8 -28 -14 -16

33 34 28 17 19 50 14 25 26 17 21 18 17 25 31 14 12 17 4 13 22 7

-44 -44 -42 -40 -38 -38 -37 -33 -32 -32 -26 -25 -26 -25 -22 -25 -25 -20 -22 -19 -17 -13

5.83 I 4.59 T 6.80 I 6.42 I 6.26 G 3.95 T 6.42 I 7.03 I 6.61 I 7.25 G 5.21 G 6.20 I 7.95 I 6.50 I 5.70 I 9.75 I 10.24 I 7.92 G 9.37 I 6.63 I 7.10 G 8.50 T

NP‚‚‚NP stack TC‚‚‚TC and NP‚‚‚NP stack double CN‚‚‚H-N hydrogen bond TC‚‚‚TC stack NP‚‚‚TC stack TC‚‚‚TC and NP‚‚‚NP stack NP‚‚‚NP stack TC‚‚‚TC stack NP‚‚‚NP stack NP‚‚‚TC stack TC‚‚‚TC perpendicular, S‚‚‚TC nitro-nitro NP‚‚‚NP offset stack N-H‚‚‚ON hydrogen bond TC‚‚‚NP perpendicular, S‚‚‚NP CN‚‚‚H-C ring, N‚‚‚H 2.55 CN‚‚‚H-C ring, N‚‚‚H 2.55 NP‚‚‚TC offset stack CN‚‚‚H-C ring, N‚‚‚H NP...TC perpendicular, C-H‚‚‚π CN‚‚‚H-N single hydrogen bond TC‚‚‚NP perpendicular, S‚‚‚NP

a Energies in kJ mol-1, distances in Å. b For each molecular pair: polymorph name and label, A, B, C, ... in descending order of energy contribution in each polymorph. c Distance between centers of coordinates (the equivalent of center of mass but with unit masses) of the two molecules in the pair; the symmetry operators are translation (T), inversion (I), screw (S), or glide (G). d NP: nitrophenyl aromatic ring; TC: thiophenecarbonitrile aromatic ring. Stack: parallel rings on top of one another, the two polar substituents (nitro or carbonitrile group) pointing in opposite directions; offset stack: stack as above but with a horizontal offset; perpendicular: nearly perpendicular rings; S-NP or S‚‚‚TC: one sulfur atom pointing approximately toward the center of the second ring..

All top pairs but one (see Table 4) belong to the ring stacking type (Figures 2 and 3) and accordingly are stabilized by a major dispersion contribution. This contribution ranges between 60 and 30 kJ mol-1 per ring-ring contact, revealing a variety of offset arrangements, adopted no doubt in compromise with other energy requirements. The Coulombic-polarization contributions are smaller but relevant, especially because of the polarity of the ring substituents. Most of the stacked ring pairs adopt an antiparallel arrangement of the polar nitro or carbonitrile substituents across a center of symmetry. This kind of arrangement is, of course, not possible when double stacking of NP and TC rings occurs in molecular pairs formed by pure translation, and that is the reason the cohesive energy of double-ring stacking is about equal to, and not double that of, single-ring stacking (compare, for instance, the Y A and YN A pairs in Table 4). The one top pair without ring stacking (OP A, Figure 4a) shows an association mode based on an intermolecular double N-H‚‚‚NC hydrogen-bond cycle across a center of symmetry. The N‚‚‚H distance is 2.42 Å, which is on the long side for such an interaction. The stabilization here mainly arises from an interplay of Coulombic and polarization terms, and the hydrogen-bond interaction is probably supplemented by the Coulombic advantage of close antiparallel arrangement of the nitrile groups. In fact, the total stabilization energy, 42 kJ mol-1, is more than twice that of the pair containing a single N-H‚‚‚NC hydrogen bond (Y D, Figure 4b), 17 kJ mol-1. In both pairs, the hydrogen bond is bifurcated, the hydrogen atom being involved simultaneously in the intramolecular N-H‚‚‚ON and the intermolecular N-H...NC hydrogen bonds. Total lattice stabilization does not correlate with hydrogen bonding. In the stable Y polymorph, the CN‚‚‚H-N hydrogen-bonded chain (Y

Figure 2. Molecular pairs showing single-ring stackings with their Pixel energies (Table 4): top to bottom, Y (A) (44 kJ mol-1), ORP (A) (40 kJ mol-1), and OP (B) (38 kJ mol-1). Note the antiparallel arrangement of the polar substituents.

D) provides 17 × 2 ) 34 kJ mol-1 stabilization (see Table 4) and is by no means the energetically most important

2186

Crystal Growth & Design, Vol. 5, No. 6, 2005

Figure 3. Molecular pair showing double-ring stacking (Table 4): ON (A) (38 kJ mol-1). Thiophene and nitrophenyl rings stack parallel, related by translation vector of length 3.95 Å. The pair YN (A) (44 kJ mol-1) shows a rather similar arrangement.

Dunitz and Gavezzotti

Figure 5. Molecular pair YN (D), Table 4 (22 kJ mol-1), with the sulfur atom of one molecule pointing at the aromatic ring of its neighbor. The spheres of sulfur atoms have been artifically magnified.

Figure 6. Molecular pair YN (C), Table 4 (25 kJ mol-1), showing double N-H‚‚‚ON hydrogen bonding: the intermolecular H‚‚‚O distance is 2.39 Å, the intramolecular one 1.89 Å.

Figure 4. Molecular pairs (Table 4) showing double and single N-H‚‚‚NC hydrogen bonding: top, double, OP (A) (42 kJ mol-1) over a center of symmetry; the nitrile groups are in close antiparallel contact with N‚‚‚N 3.57 Å, C‚‚‚C 3.41 Å; bottom, single, Y (D) (17 kJ mol-1), forming a chain with glide-plane symmetry.

interaction in that crystal structure; the less stable OP polymorph forms a ring (OP A) with a double CN‚‚‚H-N hydrogen bond system providing a higher stabilization, 42 kJ mol-1. Crystal packing analyses based on the supposition that hydrogen bonds are the strongest intermolecular interactions can be misleading. At lower cohesive energies (15-25 kJ mol-1), a group of molecular pairs of a different nature appears, in which the dispersive contributions are smaller and the Coulombic-polarization contributions are larger. One such association mode (YN D, 22 kJ mol-1, Figure 5) has the sulfur atom of one molecule pointing at the

Figure 7. Molecular pair R (E), Table 4 (25 kJ mol-1), showing a cyclic CH‚‚‚NC ring structure.

aromatic ring of its neighbor. Another (YN C, 19 kJ mol-1, Figure 6), in the same polymorph, is based on a pair of bifurcated, intra- and intermolecular, classical N-H‚‚‚ON hydrogen bonds across a center of symmetry. Yet another (R E, 25 kJ mol-1, Figure 7) adopts a cyclic CH‚‚‚NC ring structure. Here one might well speak of a C-H‚‚‚N hydrogen bond, although part of the stabilization energy clearly arises from the antiparallel pair of nitrile groups, as already noted. Figure 8 shows the

Quantitative Description of Crystal Packing

Crystal Growth & Design, Vol. 5, No. 6, 2005 2187

Figure 8. Molecular pair Y (C), Table 4 (19 kJ mol-1), showing the almost perpendicular approach of a CsH bond of a nitrophenyl group of one molecule to the thiophene ring of its partner, the pattern being repeated across a center of symmetry. This could be classified as an example of two coupled CsH‚‚‚π interactions. Table 5. Interaction Energiesa for Molecular Pairs with Short Nonbonded Atom-Atom Distances in the Six Polymorphs of 1 pair

Ecoul Epol Edisp Erep Etot

OP D -9 ORP D -9 ON B -6 ON C -11 YF -2 RF -4 YG -4 RG 6 a

-3 -3 -1 -3 -3 -2 -1 0

-9 -12 -8 -6 -11 -6 -6 -1

7 10 3 8 8 8 5 0

-13 -14 -13 -11 -8 -4 -6 5

R (coc) 8.19 S 8.02 T 9.52 G 10.40 S 8.54 T 11.91 T 12.23 T 10.27 I

Figure 9. Molecular pairs showing intermolecular chains with short C-H‚‚‚O atom-atom contacts (Table 5): top, Y (F) (8 kJ mol-1), bottom, R (F) (4 kJ mol-1). The energies are hardly distinguishable from the noise level.

short atom-atom distances C-H‚‚‚N(C) 2.58 C-H‚‚‚ON 2.39 C-H‚‚‚N(C) 2.50 C-H‚‚‚ON 2.50 C-H‚‚‚ON 2.40 H3C‚‚‚ON 3.13 none, overall destabilizing

See footnotes to Table 4 for the meaning of the various labels.

Y C pair (19 kJ mol-1) with an almost perpendicular approach of a CsH bond of a nitrophenyl group of one molecule to the thiophene ring of its partner, the pattern being repeated across a center of symmetry. This could be classified as an example of two coupled CsH‚‚‚π interactions, but it may be merely a favorable compromise, given the existence of the much stronger Y A and Y B pairs in the same crystal structure. The R D pair (25 kJ mol-1), with contact between antiparallel nitro groups (O‚‚‚O 3.20 Å, N‚‚‚O 3.03 Å), does not seem to correspond to any of the usual stereotypes. From such examples, it should be clear that the morphological labeling of the molecular pairs in terms of specific interactions between particular groups of atoms is not always clear. The two partner molecules often appear to be no more than neighbors, without any special connotation in terms of chemical groups, so they would probably have been ignored in analyses based on a repertory of predetermined synthons. Still, even in such cases, the contribution of the molecular pair to crystal stabilization may be quite significant. This difficulty appears already for molecule 1, which would not be regarded as particularly complex in a synthetic, pharmacological, or technological context. Table 5 shows the energy breakdown for molecular pairs containing short intermolecular atom-atom distances, often associated with weak hydrogen bonding. Some are illustrated in Figures 9 and 10. A large part of the cohesive energies of such pairs arises from dispersion, which is difficult to characterize and most likely due to overall, nonlocalized interactions between the diffuse electron clouds of the paired

Figure 10. Molecular pairs showing intermolecular chains with short C-H‚‚‚NC atom-atom contacts (Table 5): top, OP (D) (13 kJ mol-1), bottom, ON (C) (11 kJ mol-1).

molecules. There is no apparent correlation between C-H‚‚‚O or C-H‚‚‚N distances and stabilization energy, and the energies of some of these neighbor pairs are hardly distinguishable from the noise of second-neighbor relationships. Among these is a molecular pair (Y G, 6 kJ mol-1, Table 5) with a short C‚‚‚O distance of 3.1 Å, which might have been expected to be associated with a strong Coulombic interaction. This is not the case, and indeed this pair is practically a mere spectator in the packing arena. This example confirms again, if any further conformation were needed, that short intermolecular atom-atom contact distances do not necessarily play a relevant role in crystal stabilization.

2188

Crystal Growth & Design, Vol. 5, No. 6, 2005

A final point concerns the comparison between the sum of Pixel Coulombic and polarization energies, representing the total electrostatic energies, and point charge Coulombic energies, EPC, here calculated as electrostatic potential derived charges by the “POP ) ESP” command in the GAUSSIAN package. Over the molecular pairs listed in Tables 4 and 5, there is a good linear relationship: ECOUL+POL ) 1.1EPC - 8.3. The Pixel Coulombic energies are exact, while polarization energies depend on one empirical damping parameter. In accordance with previous experience,8-13 the conclusion seems clear; Pixel Coulombic energies exceed point charge Coulombic energies because the Pixel calculation includes contributions from overlapping delocalized electron densities, sometimes called penetration energies (appropriate procedures having been devised to avoid singularities, as extensively discussed in previous work.9-13) For the same reason, in the complete crystal structures, point charge lattice energies are only 40% of the sum of Pixel Coulombic and polarization energies. For some low-energy molecular pairs, point charge energies may even be slightly destabilizing, while their Pixel counterparts are slightly stabilizing. Inferences based on point charge energies may thus be misleading. With the Pixel partitioning, Coulombic energies never exceed 50% of the total lattice energies of the polymorphs. Concluding Remarks New conceptual elements and computational methods, involving the delocalized molecular electron charge density, have been developed for a quantitative view of the crystal packing of organic molecules. It seems commonplace to state that to understand a complex system, it is necessary to understand the individual components and the way they interact with one another. How far is it useful to break the system down into smaller and smaller components? As new evidence keeps accumulating, 21-23 it becomes more and more evident that the cohesive energy in organic crystal structures is not amenable to simple reasoning based on intermolecular atom-atom bonding interactions. Rather, it is mainly the interactions between the molecules and not between the individual atoms that provide the best basis for analyzing and understanding the cohesive energy of organic crystals. Now that a reliable electron density for many organic molecules of interest can be calculated in a matter of minutes, this approach to the analysis of molecular association patterns is feasible and much to be preferred to the rough and ready arguments based on putative atom-atom bonding. So far our results concern the analysis of experimental crystal structures. Our methods inform us about the relative energetic relevance of various molecular arrangements. In particular, they can assess the energies of particular molecular association patterns, but they do not tell why a given molecule adopts certain association patterns in preference to others in its crystal structure(s). What is still missing is the consideration of the simultaneous demands of the various molecular assocations that come into play as translationally periodic symmetry is propagated into the long-range

Dunitz and Gavezzotti

order of the crystal. Although some progress has been made,24,25 the problem of crystal structure prediction from molecular structure remains to a large extent a challenge for the future. Acknowledgment. The figures were drawn using program SCHAKAL.26 Note Added in Proof After completion of this work, we became aware of the discovery of other polymorphs of ROY (see Chen, S.; Guzei, I. A.; Yu, L. J. Am. Chem. Soc. 2005, 127, 9881-9885, and references therein). These new structures tend to support the original conclusions by Yu et al.14 that hydrogen bonding is important in stabilizing the Y form. The calculated atom-atom lattice energy of the newly discovered YT04 form, -158 kJ mol-1, puts it at the second most stable place (see Table 2) in agreement with thermodynamic findings. Supporting Information Available: Table S1 with the atomic coordinates for the six polymorphic crystal structures and the symmetry operations corresponding to the structure determinants. This material is available free of charge via the Internet at http://pubs.acs.org.

References (1) Gavezzotti, A. J. Am. Chem. Soc. 1991, 113, 4622. (2) Gavezzotti, A. Acc. Chem. Res. 1994, 27, 309. (3) Lommerse, J. P. M.; Motherwell, W. D. S.; Ammon, H. L.; Dunitz, J. D.; Gavezzotti, A.; Hofmann, D. W. M.; Leusen, F. J. J.; Mooij, W. T. M.; Price, S. L.; Schweizer B.; Schmidt, M. U.; van Eijck, B. P.; Verwer, P.; Williams, D. E. Acta Crystallogr. 2000, B57, 697. (4) Motherwell, W. D. S.; Ammon, H. L.; Dunitz, J. D.; Dzyabchenko, A.; Erk, P.; Gavezzotti, A.; Hofmann, D. W. M.; Leusen, F. J. J.; Lommerse, J. P. M.; Mooij, W. T. M.; Price, S. L.; Scheraga, H.; Schweizer B.; Schmidt, M. U.; van Eijck, B. P.; Verwer, P.; Williams, D. E. Acta Crystallogr. 2002, B58, 647. (5) Day, G. M.; Motherwell, W. D. S.; Ammon, H.; Boerrigter, S. X. M.; Della Valle, R. G.; Venuti, E.; Dzyabchenko, A.; Dunitz, J. D.; Schweizer, B.; van Eijck, B. P.; Erk, P.; Facelli, J. C.; Bazterra, V. E.; Ferraro, M. B.; Hofmann, D. W. M.; Leusen, F. J. J.; Liang, C.; Pantelides, C. C.; Karamertzanis, P. G.; Price, S. L.; Lewis, T. C.; Nowell, H.; Torrisi, A.; Scheraga, H. A.; Arnautova, Y. A.; Schmidt, M. U.; Verwer, P. Acta Crystallogr., in press. (6) Dunitz, J. D. Chem. Commun. 2003, 545. (7) Dunitz, J. D.; Scheraga, H. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 14309-14311. (8) Dunitz, J. D.; Gavezzotti, A. Angew. Chem. Int. Ed. 2005, 44, 1766-1787. (9) Gavezzotti, A. J. Phys. Chem. 2002, B106, 4145. (10) Gavezzotti, A. J. Phys. Chem. 2003, B107, 2344. (11) Gavezzotti, A. CrystEngComm 2003, 5, 429. (12) Gavezzotti, A. CrystEngComm 2003, 5, 439. (13) Gavezzotti, A. J. Chem. Theor. Comput. 2005, in press. (14) Yu, L.; Stephenson, G. A.; Mitchell, C. A.; Bunnell, C. A.; Snorek, S. V.; Bowyer, J. J.; Borchardt, T. R.; Stowell, J. G.; Byrn, S. R. J. Am. Chem. Soc. 2000, 122, 585. (15) Gaussian 03, Revision A.l, Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A. Jr.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Kiene, M.; Li, X.; Knox, J. E.; Hratchian, H. P.; Cross, J. B.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.;

Quantitative Description of Crystal Packing Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malick, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Laham, M. A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian, Inc., Pittsburgh, PA, 2003. (16) Bader, R. F. W. Atoms in Molecules: A Quantum Theory, Oxford University Press: New York, 1990. (17) Novoa, J., University of Barcelona, 2005, personal communication. (18) See, for a discussion, Coppens, P.; Volkov, A. Acta Crystallogr. 2004, A60, 357.

Crystal Growth & Design, Vol. 5, No. 6, 2005 2189 (19) Gavezzotti, A.; Filippini, G. J. Am. Chem. Soc. 1995, 117, 12299. (20) Gavezzotti, A.; Filippini, G. J. Phys. Chem. 1994, 98, 4831. (21) Gavezzotti, A.; Dunitz, J. D. Helv. Chim. Acta 2002, 85, 3949 (22) Carlucci, L.; Gavezzotti, A. Chem. Eur. J., 2005, 11, 271. (23) Demartin, F.; Filippini, G.; Gavezzotti, A.; Rizzato, S. Acta Crystallogr. 2004, B60, 609. (24) Ouvrard, C.; Price, S. L. Cryst. Growth Des. 2004, 4, 1119. (25) Day, G. M.; Chisholm, J.; Shan, N.; Motherwell, W. D. S.; Jones, W. Cryst. Growth Des. 2004, 4, 1327. (26) Keller, E. SCHAKAL92, A program for the graphic representation of molecular and crystallographic models; University of Freiburg, 1993.

CG050098Z