Intermolecular Bonding Features in Solid Iodine - American Chemical

May 30, 2014 - Department of Chemistry, University of Turin, Via P. Giuria 7, 10125 ... (Institute of Molecular Science and Technology), University of...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/crystal

Intermolecular Bonding Features in Solid Iodine Published as part of the Crystal Growth & Design Mikhail Antipin Memorial virtual special issue Federica Bertolotti,*,† Anastasia V. Shishkina,*,‡ Alessandra Forni,§ Giuliana Gervasio,† Adam I. Stash,∥ and Vladimir G. Tsirelson‡ †

Department of Chemistry, University of Turin, Via P. Giuria 7, 10125 Turin, Italy D.I. Mendeleev University of Chemical Technology, Miusskaya Square 9, 125047 Moscow, Russia § CNR-ISTM (Institute of Molecular Science and Technology), University of Milan, Via Golgi 19, 20133 Milan, Italy ∥ L.Ya. Karpov Institute of Physical Chemistry, ul. Vorontsovo Pole 10, 103064 Moscow, Russia ‡

S Supporting Information *

ABSTRACT: A detailed description of the ability of halogen bonding to control recognition, self-organization, and self-assembly in I2 crystal, combining low-temperature X-ray diffraction experiments and theoretical DFT-D and MP2 studies of charge density, is reported. The bond critical point features were analyzed using the bonding ellipsoids, in order to make them more evident and easier to compare. We showed that one-electron potential, in contrast to Laplacian of electron density, allows the electron concentration and depletion regions in the valence shell of the iodine atoms to be revealed. Thus, it was demonstrated as an effective tool for understanding the molecular recognition processes in iodine crystal, describing the mutually complementary areas of concentration and depletion of electron density in adjacent molecules. This finding was also confirmed in terms of electrostatic potential, especially using the concept of σ-hole. The tiny features of the electrostatic component of halogen−halogen interactions were also visualized through the superposition of the gradient fields of electron density and electrostatic potential. The general picture provided significant arguments supporting the distinction between Type-I (van der Waals) and Type-II (Lewis molecular recognition mechanism) I··I interactions. The energies of these interactions, evaluated on the basis of empirical relationships with bond critical points parameters, have allowed estimating the lattice energy for crystalline I2, which has been found in reasonable agreement with the experimental sublimation energy.



transfer.10 Other studies11 highlighted that the electron density transfer, from the lone pair of the Lewis base to the σ*C−X antibonding orbital or to outer portions of the halogenated molecule, can be in some cases competitive with the electrostatic contribution. In the past decade, a lot of work has been done in order to analyze the features of ED in the XB intermolecular region in molecular complexes and crystals.12 The quantum theory of atoms in molecules and crystals (QTAIMC)13 plays a significant role in these studies because it offers a consistent way of reconstruction of the atomic interactions in manyelectron multinuclear systems, whose accurate wave function is computed or high resolution ED is measured. The QTAIMC features, such as the bounded atoms separated by the zero-flux surfaces of the gradient of ED, the critical points (CPs), and the lines of maximum density between nuclei (the bond paths, BPs),14 yield a crystal structure description at the level of the bonding details.15 In the majority of cases, BPs and associated bond critical points (BCPs) occur between pairs of atoms that

INTRODUCTION Intermolecular interactions are important in assembly of molecular and supramolecular systems.1 Among them halogen bonding (XB), a widely occurring type of noncovalent interaction, is now the focus of many studies.2 It is well established that a covalently bonded halogen atom is able to simultaneously interact with both the negative and positive sites of neighboring molecules due to its strong atomic electron density (ED) anisotropy.3 Therefore, the XB interactions play an important role in molecular recognition processes,4 supramolecular chemistry,5 and crystal engineering.6 The nature of XB has been discussed for several decades from different viewpoints.7 Conventional X-ray crystallography, Raman and microwave spectroscopy might provide a lot of empirical information about XB interactions that could be completed through quantum chemical calculations. However, owing to the very wide range (5−180 kJ/mol) of XB energies, where the weak Cl···Cl interactions in chlorocarbons and the very strong I−···I2 ones in I3− complexes are the extremes,8 there are contradictory data about the nature of these interatomic interactions. Some works show its primarily electrostatic character,9 with small contributions coming from second-order terms, such as polarization, dispersion, and charge © 2014 American Chemical Society

Received: April 14, 2014 Revised: May 24, 2014 Published: May 30, 2014 3587

dx.doi.org/10.1021/cg5005159 | Cryst. Growth Des. 2014, 14, 3587−3595

Crystal Growth & Design



could be typically associated through chemical bonds and molecular interactions.16 The same is valid for the van der Waals bonded molecules:17 in this case, the BPs reflect the directional nature of the dispersive interactions, which may be considered to play a decisive role in building molecular crystals.18 The molecular interactions in halogen crystals Cl2, Br2, and I2 are of special interest. Even if the structure of these crystals is well-known,19 the nature of molecular packing suggests deviations from common van der Waals type structures. These compounds crystallize in the orthorhombic space group Cmca, with X2 (X = Cl, Br, I) molecules lying in (100) planes. This feature leads to the formation of intermolecular contacts less than the sum of van der Waals radii among molecules belonging to the same crystallographic plane; however, these layered orthorhombic crystal structures cannot be anticipated to possess isotropic potentials with only quadrupole−quadrupole interactions,20 as in the structures of N2, NO, and CO.3 The theoretical QTAIMC analysis of solid Cl23 has shown that the intermolecular BPs, corresponding to the shorter contacts, link the regions of Cl atoms where the electron density is locally depleted, as revealed by the Laplacian of electron density,▽2ρ(r) > 0, with the atomic sites in adjacent molecules where it is locally concentrated (▽2ρ(r) < 0). Such complementary matching of holes (the electrophilic activity sites) and lumps (the nucleophilic sites) in adjacent molecules influences the observed three-dimensional architecture of many molecular crystals,3,21 as was predicted a long time ago by Kitaigorodsky et al.22 After this pioneer work of Tsirelson et. al,3 the interest of the scientific community moved toward halocarbons cocrystals and interhalogen compounds forgetting the pure halogen crystals. The main drawback of these systems is that the ED in the X···X interaction region is small and rather flat: it demands both highlevel ab initio quantum-chemical calculations and very accurate X-ray diffraction experiments. This work aims at the characterization of the XB interactions in crystalline I2 by low-temperature high-resolution X-ray diffraction data and DFT-D and MP2 calculations, using quantum-chemical periodic boundary condition codes. In principle, X-ray diffraction combined with modern data processing software is able nowadays to reconstruct the details of ED in crystals with an accuracy of 0.03−0.05 eÅ−3,23 if single crystals of excellent quality are available. However, the latter condition could be hardly satisfied in the case of I2, owing to unstable single crystal samples and related nontrivial extinction and absorption contributions, due to the presence of heavy atoms. Therefore, in the present study the experimental drawbacks, that we met dealing with this almost unsuitable system, will be as well described. Iodine is a fundamental element especially for its role in several biochemical processes and for some technological applications.24 Even if the structure of iodine at room temperature has already been determined in an earlier X-ray diffraction study by Harris et al.25a and refined later,25b,c the discovery of the mechanism through which XB occurs in this molecular crystal could be of some interest in the understanding of enzymatic reactions in human body or in general of its catalytic function and as redox mediator in many useful devices.

Article

MATERIALS AND METHODS

Experimental Section. Iodine single crystals have been obtained from commercial solid product by sublimation at room temperature. Seven experiments with different samples were performed at 110(1) K with a Gemini R Ultra diffractometer26 equipped with low temperature device (Oxford Cryojet with N2 stream), graphite-monochromatized λ(MoKα) = 0.71073 Å and Ruby CCD area detector. The goodness of each data set has been checked using Abrahams−Keve plots,27 implemented in XRDKplot program,28 and the most reliable one has been chosen (see Supporting Information, Figure S1). A total of 12 300 frame images in 36 scan sets with a scan angle Δω = 0.3° were measured; an integration time of 5.0 s was used for each scan set at the detector setting angle 2θ ≈ 40°, and an integration time of 35.0 and 50.0 s was used for scan sets at the detector setting angle 2θ ≈ 92°. The frame images were integrated using CrysAlis Pro software29 with a sufficiently large spot size (1.0), accounting for the Kα1−α2 splitting. The data were sorted and merged using SORTAV30 giving 1531 independent reflections with I > 2σ(I) (Rint = 0.0522) with a completeness of 100%. A spherical atom full-matrix least-squares refinement using SHELX-9731 was initially performed on F2 by using all data. Furthermore, the multipole hexadecapole-level Hansen− Coppens32a model, as implemented in the XD2006 program suite,33 was used, performing the required aspherical refinement. The atomic relativistic wave functions of Macchi & Coppens32b were used to describe both the core- and the valence electron shells. Atomic displacement parameters were described following the harmonic approximation. The Becke and Coppens isotropic secondary extinction correction34 was applied. Only 1281 reflections with I > 3σ(I) were included in the F2 refinement; the weighting scheme w = 1/[s2 + (ap)2], where s = σ(Fo2), p = f × Fo2 + (1 − f) × Fc2, f = 0.3333, a = −2.00, has been applied and figures of merit R = 0.0187, wR = 0.0444, S = 1.039 have been obtained at the convergence of the refinement. A noticeable improvement has been achieved including, in the leastsquares procedure, the third- and fourth-order Gram−Charlier expansion coefficients, to take into account the anharmonicity of atomic displacements parameters: R = 0.017, wR = 0.026, S = 1.040. The refined multipole parameters were used to calculate all the experimental-density-based features with WinXPRO program.35 More detailed experimental information is collected in Table 1. Computational Details. Kohn−Sham36 periodic-boundary method, using B3LYP hybrid functional augmented with the empirical Grimme dispersion correction37 (B3LYP-D, with coefficient C6 = 31.50 J·nm6/mol) was used in solid I2 calculations, as implemented in the CRYSTAL09 code.38 All-electron calculations have been performed with progressively extended basis sets. More accurate calculations have been carried out using tailored DZV quality basis set plus polarization functions.39 The truncation criteria for bielectronic integrals has been set to 7 7 7 7 15, and the shrinking factor of the commensurate reciprocal space grid was set to 4 corresponding to 21 independent kvectors in the irreducible Brillouin zone. The optimization of atomic positions has been carried out using analytical energy gradients; cell parameters have been fixed to the experimental values. This approximation provides reliable and consistent results in studying the intermolecular interactions.40 No imaginary frequencies have been found after computing the vibrational frequencies at Γ point.37a The XFAC module of the program suite was used in order to generate theoretical structure factors at the experimental resolution. The efficiency in description of long-range and dispersive interactions from solid state Kohn−Sham calculations was also compared with the results obtained applying the second-order local Møller−Plesset (LMP2) electron correlation correction to HF solution. While the former is easier and cost-effective because of the faster convergence of the results with respect to the basis set quality, the LMP2 is a full ab initio method, which provides robust and consistent results.41 This perturbative local approach, originally proposed by Pulay42 and further developed by Werner and co3588

dx.doi.org/10.1021/cg5005159 | Cryst. Growth Des. 2014, 14, 3587−3595

Crystal Growth & Design

Article

with A in Figure 1) and one at 3.9796(2) Å (B). Eight contacts linking the planes above and below are at slightly larger

Table 1. Crystallographic Data and Experimental Details compound formula

I2

Mr space group crystal system a/Å b/Å c/Å V/Å3 Z Dcalc/g cm−3 F(000) λ/Å μ/mm−1 crystal size/mm transmission factors (range) θ range/deg no. data used for merging no. unique reflections no. of data included in refinement (I > 3σI) hkl range

126.90 Cmca orthorhombic 7.1589(1) 4.6915(1) 9.8014(1) 329.190(9) 4 5.129 424.0 Mo−Kα/0.710 73 18.82 0.02 × 0.09 × 0.15 0.183−0.644 4.16−66.80 23 746 1531 1281 −16 < h < 18 −12 < k < 11 −25 < l < 22 0.0522, 0.0150 24 0.0177, 0.0172 0.0257 1.040 0.005003 0.826/−0.773 2σ(I) wR(F2) goodness of fit (S) refined extinction coefficient largest peaks in ED/e·Å−3 (all data)b max shift/esd in last cycle

Figure 1. Molecular structure of the iodine crystal. The thin lines show the intermolecular bond paths of different types A−E derived from the model ED refined from experimental data. The bonding ellipsoids for atomic interactions in solid I2 are depicted by green colors.

distances: 4.3499(1) Å (C), 4.2796(1) Å (D), and 4.4234(1) Å (E). I−I bond length from aspherical multipole-model refinement is 2.7179(2) Å (experiment) and 2.789 Å (theory). These values are in agreement with Bolhuis et. al (2.715 Å),25c and they are significantly longer than those computed for gas phase (2.662 Å),46 as a consequence of the presence of stronger XBs in solid and of the formation of a two-dimensional network within (100) planes. The geometrical features of XB in halocarbons (C−X···X−C, where X = Cl, Br, I) have been classified into two different types, based on the values of two angles, θ1(C−X···X) and θ2(X···X−C).47 In the case of iodine crystal, this classification can be extended from halocarbons to halogen molecules, as shown in Figure 2. Each iodine molecule forms only one

R = Σ(|Fo| − |Fc|)/Σ|Fo|; Rw = {Σ[w(Fo − Fc)2]/Σw(Fo)2}1/2 (where w is the XD2006 weighting scheme: w = 1/[s2 + (ap)2] (where s = σ(Fo2), p = f × Fo2 + (1 − f) × Fc2, f = 0.3333, a = −2.00)); Rint = Σ[n/ (n − 1)1/2 |Fo2 − Fo2 (mean)|/ΣFo2; Rσ = Σ[σ(Fo2)/Σ(Fo2)]; S = {Σ[w(Fo2 − Fc2)]/(nobs − npar)}1/2. bThe residual map ρo − ρmult is reported in Figure S2 (Supporting Information). a

workers,43 has been implemented recently in the CRYSCOR09 periodic code.44 Automatic full geometry optimizations at the correlated LMP2 level are not yet feasible within the present implementation of CRYSCOR09 code; therefore, single point calculations on Grimmecorrected optimized geometry have been performed. After that, a sizeextensive correlation correction to the HF result was performed in order to obtain the theoretical structure factor list, from which electron density distribution in crystals has been reconstructed. The hexadecapole-level Hansen−Coppens32a multipole model as implemented in the XD2006 program suite33 was used to perform the aspherical refinement on the theoretical structure factors, to be compared with the experimental ones. Unit weights were assigned to the theoretical structure factors during the refinement.45 The multipole parameters derived in this way were used to calculate the electron-density-based quantities through the WinXPRO program.35

Figure 2. Representation of I···I intermolecular interactions Type-I (e.g., I1···I3 in the picture) and Type-II (e.g., I1···I4).

interaction belonging to Type-II with the closer molecule in the same (100) plane, as a consequence of the presence of a 21 screw axis; it is also involved in other four symmetry independent Type-I I···I interactions, one in the same molecular plane and three among different planes, generated by the presence of inversion centers both on (100) planes and among them (see Supporting Information for further geometrical details, Table S2). Topological Analysis of Experimental and Theoretical Electron Densities. The bond critical points (BCPs) and the corresponding bond paths (BPs) in electron density were identified for all previously found geometrical contacts (Table



RESULTS AND DISCUSSION Crystal Structure. The iodine crystal is orthorhombic at 110 K with space group Cmca, and a = 7.1589(1), b = 4.6915(1), c = 9.8014(1) Å, V = 329.190(9) Å3, Z = 4. The molecules crystallize in layers parallel to (100), linked through other weaker intraplane I···I interactions. Each iodine atom has three close geometrical contacts with distances less than the sum of van der Waals radii (RvdW = 1.98 Å), laid in (100) planes: two of them are at distances of 3.5010(2) Å (labeled 3589

dx.doi.org/10.1021/cg5005159 | Cryst. Growth Des. 2014, 14, 3587−3595

Crystal Growth & Design

Article

Table 2. Theoretical B3LYP/DZVP/Model (First Line), B3LYP-D/DZVP/Model (Second Line), MP2/DZVP/Model (third line) and Experimental (Fourth Line) BCPs Features. The Electron Density Parameters Are Given in A.U.a contact I1−I2

I1···I3 A

I1···I4 B

I1···I5 C

I1···I6 D

I1···I7 E

R, Å

ρ(r)

▽2ρ(r)

λ1

λ2

λ3

2.798 2.789 2.789 2.710 3.468 3.484 3.484 3.501 3.940 3.915 3.915 3.980 4.357 4.343 4.343 4.280 4.280 4.280 4.280 4.350 4.398 4.410 4.410 4.423

0.054 0.055 0.055 0.050 0.018(0.019) 0.018 0.018 0.015 0.008(0.009) 0.009 0.009 0.009 0.004(0.005) 0.004 0.005 0.004 0.005(0.005) 0.005 0.005 0.005 0.004(0.004) 0.004 0.004 0.004

0.078 0.078 0.072 0.082 0.040(0.033) 0.039 0.037 0.036 0.020(0.020) 0.021 0.022 0.022 0.011(0.012) 0.011 0.011 0.011 0.012(0.013) 0.012 0.012 0.012 0.010(0.011) 0.010 0.010 0.009

−0.0346 −0.0348 −0.0357 −0.0393 −0.0101 −0.0098 −0.0096 −0.0074 −0.0040 −0.0043 0.0044 −0.0047 −0.0017 −0.0018 −0.0018 −0.0017 −0.0020 −0.0020 −0.0020 −0.0020 −0.0014 −0.0013 −0.0013 −0.0013

−0.0340 −0.0343 −0.0352 −0.0368 −0.0100 −0.0097 −0.0094 −0.0069 −0.0039 −0.0041 0.0043 −0.0044 −0.0009 −0.0009 −0.0009 −0.0009 −0.0016 −0.0016 −0.0016 −0.0018 −0.0013 −0.0012 −0.0013 −0.0012

0.1462 0.1467 0.1429 0.1568 0.0597 0.0582 0.0557 0.0500 0.0284 0.0298 0.0308 0.0311 0.0136 0.0140 0.0141 0.0131 0.0156 0.0157 0.0159 0.0158 0.0126 0.0123 0.0125 0.0121

a

For comparison, the values computed directly from the periodic wave-function are given in parentheses. The intermolecular contacts I···I are labeled by the same letters reported in Figure 1.

shif t bond, for which the charge is displaced from the internuclear space toward atomic cores.50 The intermolecular I···I bonds (noncovalent closed-shell interactions) are characterized by λ3 ≫ |λ1|,|λ2|, i.e., c |λ1|,|λ2|; therefore c < a,b. This result is opposite to that observed in typical covalent bonds, for which c > a, b13 and the corresponding ellipsoids are elongated parallel to the bonding interaction lines. In this respect, bond I−I resembles a specific kind of chemical bond called charge 3590

dx.doi.org/10.1021/cg5005159 | Cryst. Growth Des. 2014, 14, 3587−3595

Crystal Growth & Design

Article

Such a feature of the Laplacian prevents its application to study the bonding features in compounds containing heavy atoms. Thus, to characterize the bonding features in solid iodine from the ED viewpoint, other descriptors are needed. For this purpose the one-electron potential (OEP) has been computed:58

Table 3. B3LYP/DZVP (First Line), B3LYP-D/DZVP (Second Line) and MP2/DZVP Computed (Third Line) and Experimental (Fourth Line) Characteristics of the Intermolecular Interactions at the Bond Critical Points in the Electron Densitya type contact

G(rBCP)

−V(rBCP)

Eintb

Eintc

I1···I3 A

0.010 0.010 0.010 0.008 0.004 0.004 0.005 0.005 0.0022 0.0022 0.0022 0.002 0.0024 0.0024 0.0025 0.0024 0.0019 0.0019 0.0019 0.0018

0.010 0.010 0.010 0.008 0.004 0.004 0.004 0.004 0.0016 0.0016 0.0017 0.001 0.0018 0.0018 0.0018 0.0018 0.0014 0.0014 0.0014 0.0013

17.6 17.6 17.6 14.1 7.0 7.0 8.8 8.8 3.9 3.9 3.9 3.5 4.2 4.2 4.4 4.2 3.3 3.3 3.3 3.2

17.9 17.9 17.9 14.3 7.1 7.1 7.1 7.1 2.9 2.9 3.0 1.8 3.2 3.2 3.2 3.2 2.5 2.5 2.5 2.3

I1···I4 B

I1···I5 C

I1···I6 D

I1···I7 E

P(r) =

1 ∇2 ρ(r) 1 ∇ρ(r) ·∇ρ(r) − 4 ρ(r) 8 ρ2 (r)

The OEP defines exactly the potential governing the motion of a single electron in a many-electron system. Topological similarity was found59 between ▽2ρ(r) and OEP as reflected in the valence shell structures of light atoms, while OEP performs better for heavy atoms.60 OEP reveals minima and maxima associated with core atomic electronic shells and areas of electron concentration/depletion in the valence shells. They are indicated by negative and positive values of OEP, respectively. The negative areas of P(r) correspond to positive values of local one-electron kinetic energy meaning that the electron is classically allowed in these areas. On the contrary, positive values of P(r) reveal the potential barriers where electrons exhibit the quantum behavior. The first term in the expression for P(r) has the same sign as ▽2ρ(r), while the second term is always negative. Thus, classically allowed regions of P(r) correspond to regions of negative Laplacian of the electron density. The OEP maps for crystalline iodine computed from our MP2/DZVP calculation using the multipole model are shown in Figure 3. Unfortunately the same maps reconstructed from the experimental electron density yielded a picture with the electron lone-pair peaks inclined with respect to the I−I line (see Figure S4a, Supporting Information). Careful analysis of the data measurement, treatment, and refinement procedures did not allow location of the origin of this phenomenon to be determined. The attempts to apply the anisotropic extinction correction have also failed. Anyway, because of the good agreement of the experimental and theoretical BCPs features (Table 2), we could consider this anomalous behavior as an artifact, resulting from unrecovered local reciprocal-space distortion in the data. Therefore, from now on only maps obtained from theoretical calculations will be discussed. It is noteworthy, however, that a similar lone-pair behavior has already been observed in the experimental density study of a cocrystal of N-methylpyrazine iodide with I2.57a The shape of the OEP around I atom in solid iodine (Figure 3) shows the anisotropic distribution of the ED around the

a

ED parameters derived by multipole refinement are given in a.u. The interaction I···I energies, Eint, are in kJ/mol. The intermolecular contacts I···I are indicated by letters according to Figure 1. bEint = −γ· G(rBCP); γ = 0.67.54 cEint = β·V(rBCP); β = 0.68.54

shells.13 This property is used to locate the local electronic charge concentration and depletion regions in the outer electronic shells of bounded atoms, which are associated with the localized bonded and nonbonded Lewis electron pairs in the VSEPR model. 55 However, Kohout et al. 56 have demonstrated that beyond the third row elements, Laplacian of the ED is not able to resolve the outermost electronic shells from the atomic core, because the charge depletion of the d orbital from the penultimate shell overwhelms the charge concentration of the valence s and p orbitals. Figures S3b (experimental) and S6 (from calculation) of Supporting Information demonstrate that this situation occurs in the case of both a crystal (S3b and S6c) and a gas phase molecule (S6a) of iodine as well: Laplacian lacks regions with ▽2ρ(r) < 0.12c,57

Figure 3. One-electron potential distribution in iodine crystal (MP2/DZVP/model) in (a) the plane of contacts A and B; (b) the plane of the contacts C and D; (c) the plane of the contacts C and E. Line intervals are ±2,4,8 × 10n e·Å−5 (−2 ≤ n ≤ 2); solid lines correspond to negative oneelectron potential values. Bond and ring critical points are indicated by circles and triangles, respectively. The bond paths are marked by bold lines. 3591

dx.doi.org/10.1021/cg5005159 | Cryst. Growth Des. 2014, 14, 3587−3595

Crystal Growth & Design

Article

Figure 4. Theoretical (MP2/DZVP) electrostatic potential in the free iodine molecule (a) and molecule removed from the crystal (b) projected onto the 3D isosurfaces of the electron density (0.002 au) together with the intermolecular I···I bond paths. The σ-holes (blue regions) in A and A′ are associated with two different iodine atoms involved in XBs with the reference molecule.

Figure 5. Superposition of gradient fields in the theoretical (MP2/DZVP model calculation) electron density (blue) and the electrostatic potential (orange) in the planes containing intermolecular bond paths of type A and B (a), A, C, D (b) and B, C, and D (c). The bold lines correspond to the zero-flux atomic boundaries. The critical points (3, −1) in the electron density are in black.

Inner-Crystal Electrostatic Potential. Figure 4 shows the distribution of inner-crystal electrostatic potential (ESP) on an isosurface of electron density of 0.002 a.u. in a gas-phase I2 molecule (a) and in a molecule removed from the crystal (b). In the case of isolated molecule (Figure 4a), the negative “belt” of ESP is distributed in the equatorial plane and reflects ED concentration around the I atom. Simultaneously, an area of charge depletion behind the I atom on the extension of the intramolecular I−I bond forms the cone-like region of positive ESP. This feature, well-known as the σ-hole,2a,9c,65 has been found also in molecules with C−X (with X = Cl, Br, I) functional groups.12c,e,66 The axial symmetry of ESP distribution is present in the free molecule (Figure 4a); in contrast, the molecule removed from the crystal (Figure 4b) carries the f ingerprints of interactions with the nearest neighbors: axial symmetry of ESP is broken, and each σ-hole (region of maximum OEP) of I atoms is oriented toward the less positive ESP belt of the closer I atom in neighboring molecules (Type-II interaction, contact A and A′). It is worth noticing that a corresponding intermolecular I···I bond path runs through the σ-holes in the ESP. Thus, the amphoteric behavior of iodine is here explained, and additional information about the Lewis-type recognition mechanism of I2 molecules involved in intermolecular interactions within a crystal is provided. These Type-II interactions (XB) are driven, indeed, by the electrostatic factor, as discussed in ref 47a. To better understand the features of the electrostatic contributions to the atom−atom interactions in iodine crystal, it is possible to analyze the inner-crystal electrostatic field, E(r)

nucleus: the toroidal area of the negative OEP values localized in the equatorial plane can be associated with the nucleophilic part of the valence shell, while the electrophilic area is placed on the tip of σ-covalent bond. Figure 3a shows that the shortest I···I contacts labeled by A and A′ (I···I distance is 3.484 Å, the sum of van der Waals radii is 3.96 Å) are characterized by a Lewis-type mechanism,21c,61 with nucleophilic and electrophilic regions of the adjacent iodine atoms in the bc layer complementing each other along the bond paths. These interactions belong to Type-II and may be identified as XB.47a All the others interactions are arranged in order to increase the contact surface area of neighboring molecules, with internuclear distances exceeding the sum of the atomic van der Waals radii. It can be concluded that they are driven only by steric factors,62 and the molecular complementarity is realized in this case through a van der Waals mechanism, because the face-to-face molecular arrangement favors the dominant dispersion contribution to the interaction energy. It is reasonable to compare our OEP results with those obtained by the electron localization function (ELF).63 The ELF maps (Figure S4b from experimental data and Figure S8 from calculations), approximately computed according to ref 64, support the OEP conclusions: they distinctively show the closed-shell nature of intra- and intermolecular interactions in iodine crystal and point out the key-lock mechanism for interactions Type-II (contact A), against the van der Waals complementary interactions of Type-I (contacts B, C, D, and E). The latter do not show visible influence on the ELF distribution. 3592

dx.doi.org/10.1021/cg5005159 | Cryst. Growth Des. 2014, 14, 3587−3595

Crystal Growth & Design

Article

= −▽φ(r), where φ(r) is the electrostatic potential. Nuclei of neighboring atoms in any crystal are separated in the field E(r) by surfaces Pi, satisfying the zero flux condition E(r)·n(r) = 0, ∀r ∈ Pi(r), where n(r) is a unit vector normal to the surface Pi at r. Each surface Pi defines the φ-basin of the ith atom, inside of which the electrons are attracted to the nuclei i.67 The ▽φ(r) gradient field allows us to highlight the boundaries of the atomic-like φ-basins in solid iodine and to reveal the network of atomic electrostatic interactions in crystals. The superposition of the gradient fields of ▽φ(r) and ▽ρ(r) in solid iodine is shown in Figure 5. It allows revealing of the network of atomic electrostatic interactions in this crystal. Indeed, because the electron density within each φbasin is attracted to the corresponding nucleus, the overlap of the ρ-basins and φ-basins of adjacent atoms A and B leads to situation in which the part of the ED belonging to A atom falls into φ-basin of B atom. The combination of the ▽ρ(r) and ▽φ(r) gradient fields manifests the general physical mechanism of the electrostatic interactions in a molecule or a crystal. The significant mutual penetration of ρ- and φ-basins of neighboring molecules takes place along the bond paths only for Type-II interactions (contact A) supporting the considerations already made for ESP analysis: the Type-II interactions may be considered as driven by the electrostatic factor. For the other interactions (Type-I), the boundaries of φ-basins and ρbasins along the bond paths almost perfectly overlap. This behavior leads to the conclusion that molecular packing pattern in solid iodine is formed in such way to minimize repulsion along B, C, D, and E contacts by interfacing the neutral regions of their electrostatic potential surfaces. This general picture brings significant arguments in favor of the nature of the Type-I and Type-II halogen−halogen interactions as suggested in ref 47a. Estimation of the Molecular Interaction Energy in Solid Iodine. The total molecular interaction energy in homomolecular crystals is a well-defined property, which can be computed by quantum-chemical methods as Elatt = E(bulk)/Z − E(mol) + EBSSE.47d,61b,68 Here E(bulk) is the total energy of the unit cell, E(mol) is the total energy of the isolated gas-phase molecule calculated at the same level of theory, and EBSSE is the correction for basis set superposition error. B3LYP-D calculations for solid iodine led to Elatt = 67 kJ/mol (correction for zero-point energy has been applied). This value agrees well with the experimental value of sublimation energy at 0 K of crystalline I2 (69 kJ/mol),69 also proving the reasonable quality of the wave function used in this work. Unfortunately, the knowledge of Elatt does not allow the determination of the energy of each specific interatomic interaction. QTAIMC suggests a convenient approximate solution of this problem, because the energy of each specific interaction, corresponding to a BP, can be considered independently from the others, and Elatt may be approximated by a sum of the energies of the pairwise atomic X···Y interactions:70 ′ = E latt

densities at the X···Y bond critical points: Eint(X···Y) = β· V(rBCP) and Eint(X···Y) = −γ·G(rBCP). Considering the O···H hydrogen bond, it was shown that β = 0.571a and γ = 0.42971b (the atomic units are used). The coefficients β and γ are not universal, however, and they depend on the nature of X and Y atoms. The authors in ref 54 have considered a series of molecular complexes with XBs and have shown that β = 0.68 and γ = 0.67, for interactions involving Iodine atoms. The E′latt energy for solid iodine computed with these coefficients via the G(rBCP) and V(rBCP) values from Table 3 (B3LYP-D) are 70 and 63 kJ/mol, respectively. Both these values are close to the experimental absolute value of sublimation energy of the crystalline I2, 69 kJ/mol,69 with the value derived from the kinetic energy density looking preferable. Thus, the approach used in this work allowed evaluation of the lattice energy of a crystal using the quantum-topological features of intermolecular halogen interactions.



CONCLUSIONS The ability of halogen−halogen interactions to control recognition, self-organization, and self- assembly in I2 crystal has been extensively assessed through this work. The electron density features at the bond critical points derived using QTAIMC from both low-temperature X-ray diffraction experimental data and DFT-D and MP2 calculations enabled us to locate and even rank the intermolecular closed-shell interactions in solid iodine. To make the bonding features easier to compare, we represented the bond critical points properties by means of the bonding ellipsoids, whose sizes comparison allowed us to distinguish I···I interactions according to their strength. Using DFT-D and MP2 calculations, we were able to explain why previous QTAIM papers showed controversial results about Laplacian distribution around halogens atoms: the Laplacian is not able to resolve the outermost shell from the core, due to the influence of d orbitals. This effect is even increased by the relativistic expansion of the d orbitals and the contraction of the s shells going toward heavier atoms. We demonstrated that another descriptor, the one-electron potential, does not have this disadvantage, and its distribution showed the origin of iodine amphoteric behavior: it results from electron density concentration around I atom in the equatorial plane of the I2 molecule and electron density depletion on the extension of the I−I bond behind this atom. Thus, one-electron potential reveals the mutually complementary areas of concentration and depletions of electrons in adjacent I2 molecules, explaining the molecular recognition mechanism in the title crystal. This finding was also highlighted analyzing the electrostatic potential distribution, based on the concept of σhole. The tiny features of electrostatic component of halogen− halogen interactions in solid iodine have been also established through the superposition of the gradient fields of electrostatic potential, ▽φ(r), and of electron density, ▽ρ(r): this combination highlighted the features of the Type-II interactions (XB), corresponding to the shortest I···I contacts realized in the I2 crystal. The lattice energy of solid iodine was evaluated using the correlation between I···I interaction energy and the local potential/kinetic energy density at the corresponding bond critical points. The corresponding result, in good agreement with the computed total molecular interaction energy72 and with the experimental sublimation energy69 of iodine proved

∑ ∑ Eint ,Y ,X X Y