Article pubs.acs.org/JCTC
Cite This: J. Chem. Theory Comput. XXXX, XXX, XXX-XXX
Modeling Halogen Bonds in Ionic Liquids: A Force Field for Imidazolium and Halo-Imidazolium Derivatives Carlos E. S. Bernardes* and José N. Canongia Lopes Centro de Química Estrutural, Instituto Superior Técnico, 1049-001 Lisboa, Portugal S Supporting Information *
ABSTRACT: In this work, a force field for molecular dynamics and Monte Carlo simulations of ionic liquids containing imidazolium and halo-imidazolium derivatives is presented. This force field is an extension of the well-known CL&P and OPLS-AA models and was validated by comparing predicted crystalline structures for 22 ionic liquid compounds with the corresponding data deposited at the Cambridge Structural Database. The obtained results indicate that the proposed force field extension allows the reproduction of the crystal data with an absolute average deviation lower than 2.4%. Finally, it was also established that the halogen atoms covalently bound to the studied imidazolium cations are positively charged and do not exhibit a so-called σ-hole feature. For this reason, the formation of halogen bonds in the proposed force field appears naturally from the parametrized atomic pointcharge distribution, without the necessity of any extra interaction sites.
■
INTRODUCTION The term “halogen bond” usually refers to a directional noncovalent interaction between a covalently bonded halogen atom, X in a R−X molecule, and a Lewis base species, B. Although not as well-known as their hydrogen-bond counterparts, these R−X···B interactions have been known for a long time,1−3 particularly in the context of donor−acceptor interactions in the solid phase. More recently halogen bonding has been recognized as a major player in different chemical and biological processes.4,5 Halogen bonds are somewhat counterintuitive in the sense that one does not expect an electron-donor species (the Lewis base B) to interact with a covalently bound halogen atom with a complete valence shell and an overall negative character. However, R−X···B interactions are easy to understand if one takes into account the depletion of electronic charge along the extension of the R−X bond (the so-called σ-hole) caused by the geometry and interactions of the molecular orbitals involved in the R−X covalent bond.6 Such depletion results in a narrow region of positive electrostatic potential at the outer side of the halogen atom opposite to its covalent bond, which can interact with the electrons from the B species and trigger the formation of the halogen bond. Some molecular ions that constitute ionic liquids can contain covalently bound halogen atoms such as chlorine, bromine, or iodine. For instance, halogen substitutions on the imidazolium cationone of the most common cations in ionic liquid formulationsyield highly dense and hydrophobic ionic liquids with very interesting potential applications.7,8 Furthermore, halogenation is also a route to analyze in a systematic way the © XXXX American Chemical Society
nature and strength of cation−anion interactions and intermolecular forces in different ionic liquids, with the underlying motivation of assisting on the rational design of new ionic liquids with specific properties.9,10 Nowadays, it is accepted that the properties of ionic liquids emerge from the interplay and balance between electrostatic and van der Waals interactions. Hydrogen bonding can also play a significant role in such context.11 It is time to add halogen bonding to this list. The present work deals with the extension of the well-known CL&P force field12used in the systematic modeling of ionic liquid homologous seriesto incorporate halogenated (or otherwise substituted) imidazolium cations. Those cations will interact with different anions (Lewis bases) and have the potential to yield halogen bonds. Of particular interest will be ionic liquid formulations where covalently bound halogen atoms will form halogen bonds with halide anions.
■
FORCE FIELD
The force field presented in this work is an extension to the CL&P parametrization,12 to include halogen atoms directly bonded to the imidazolium ring and/or different types of hydrocarbon substitutions (Figure 1). The potential function is defined in the general form: Received: June 20, 2017 Published: November 1, 2017 A
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
Table 1. Intermolecular Force Field Parameterization Established According to Equation 1a atoms
Imidazolium Derivatives 0.15 NA NAp 0.08 NAM 0.04 NA1 0.07 NA2 −0.04 NA3 −0.07 CR −0.11 CRM 0.19 CRCl −0.03 CRBr −0.09 CRI −0.15 CW −0.13 CWM 0.06 CWCl −0.01 CWBr −0.06 CWI −0.10 CW1 −0.25 CWM1 0.15 CWCl1b −0.13 CWBr1b −0.18 CWI1b −0.22 C1 −0.17 C1p 0.13 C1i 0.15 C2 0.01 CE −0.05 CS −0.12 CT −0.18 CCM −0.26 CM1 −0.035 CM −0.23 HM 0.115 H1 0.13 H1p 0.04 H1i 0.07 HC 0.06 HCW 0.21 HCR 0.21 ClCR 0.13 ClCW 0.09 BrCR 0.19 BrCW 0.14 ICR 0.25 ICW 0.18 Anions Cl− −1 Br− −1 I− −1
Figure 1. Nomenclature for the interaction sites of the imidazolium derivatives investigated in this work. The highlighted blue areas indicate that, after this point, the OPLS-AA force field should be used.13,14 In the case of the halogen substitution, X corresponds to Cl, Br, and I. bonds
U=
∑ ij
kr , ij 2
dihedrals
+
angles
(rij − ro , ij)2 +
ijk 4
∑ ∑ ijkl
∑
n=1
Vn , ijkl 2
kθ , ijk 2
(θijk − θo , ijk)2
[1 + (− 1)n cos(nφijkl)]
⎧ ⎡⎛ ⎞12 ⎛ ⎞6 ⎤ qiqj ⎫ σij ⎥ ⎪ ⎢ σij ⎪ ⎬ + ∑ ∑ ⎨4εij⎢⎜⎜ ⎟⎟ − ⎜⎜ ⎟⎟ ⎥ + r ⎝ rij ⎠ ⎦ 4πεorij ⎪ i j>i ⎪ ⎩ ⎣⎝ ij ⎠ ⎭
q, e
(1)
In the previous equation, rij and ro,ij are the distance and equilibrium distance between atoms i and j, respectively; θijk and θo,ijk are the angle and equilibrium angle between atoms i− j−k, respectively; kr,ij and kθ,ijk are the bond and angle harmonic force constants, respectively; φijkl is the dihedral angle i−j−k−l and Vn,ijkl the corresponding Fourier constants; σij and εij are the 12-6 Lennard-Jones (LJ) parameters; εo is the vacuum permittivity; and qi and qj and the atomic point charges (APCs) of i and j. In eq 1 the three initial summations describe the intramolecular interactions (covalent bonds, valence angles, and torsion dihedrals), while the last two summations represent the van der Waals (VDW) and Coulombic interactions, used to describe the intermolecular forces. Within the same molecule, the electrostatic and VDW interactions between atoms separated by one or two bonds are ignored, and those separated by three bonds scaled by a 0.5 factor. As in the case of the OPLS-AA force field,13 the LJ parameters between different types of atoms are computed using standard combination rules based on geometrical means. Intermolecular Potential. The atomic point charges and 12-6 Lennard-Jones (LJ) parameters proposed in this work are given in Table 1. The LJ parametrization used for the cations was retrieved from previously published CL&P data12,14,15 and from the OPLS-AA force fields13 that are frequently used for the simulation of ionic liquids and organic compounds, respectively. For the halogen atoms attached to the imidazolium ring, it was assumed that the parametrization reported for halobenzene molecules16 could be used. Regarding the anions, the parametrization in CL&P was followed,12 but
ε, kJ·mol−1
σ, Å
0.71128 0.71128 0.71128 0.71128 0.71128 0.71128 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.29288 0.27614 0.27614 0.27614 0.27614 0.27614 0.27614 0.27614 0.27614 0.31798 0.31798 0.12552 0.12552 0.12552 0.12552 0.12552 0.12552 0.12552 1.25600 1.25600 1.96780 1.96780 2.51040 2.51040
3.25 3.25 3.25 3.25 3.25 3.25 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.55 3.50 3.50 3.50 3.50 3.50 3.50 3.50 3.50 3.55 3.55 2.42 2.50 2.50 2.50 2.50 2.42 2.42 3.40 3.40 3.47 3.47 3.75 3.75
0.72000 0.86000 1.05000
3.71 3.97 4.31
a
The data in bold were obtained in this work. The values in normal font were taken from the OPLS-AA or CL&P force fields.13,14,16,19. b CW connected to a halogen atom, in an α position to a substituted CW carbon atom.
for the case of the halogen anions, a revision of the previously recommended data was performed. This procedure included fitting the Born−Huggins−Mayer (BHM) potential for chloride, bromide, and iodide, reported by Tosi and Fumi,17,18 to the 12-6 LJ relation in eq 1 (see Supporting B
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
Table 2. Intramolecular Force Field Parametrization Established in This Work, According to the Nomenclature in Figure 1 and Equation 1a bond C1p−CM1 CM−CMc CM−HMc NA−CWb NA−CRb CW−CWb CR/w−Clb CR/W−Brb CR/W−Ib CR/W−Hb CR/W−CCMb NA−C*b,d C*−H*d C*−C*d angles NA−CR−Cl/Br/Ib NA−CW−Cl/Br/Ib NA−CR−HCRb NA−CW−HCWb CR−NA−CWb NA−CR−NAb C*−NA−CWb,d C*−NA−CRb,d NA−C*−H*b,d NA−C*−C*b,d NA−C1p−CM1b C1p−CM1−CM C1p−CM1−HM CM−CM−HMc CM1−C1p−H1p HM−CM−HM H*−C*−H*d H*−C*−C*d C*−C*−C*d CW/R−CCM−H1b CW/R−CCM−C*b,d NA−CW−CCMb CW−CW−CCMb CW−CW−Cl/Br/Ib CW−CW−HCWb
ro, Å 1.510 1.340 1.080 1.378 1.315 1.341 1.725 1.870 2.080 1.080 1.510 1.466 1.090 1.529 θ0, deg 125.1 122.0 125.1 122.0 108.0 109.8 125.6 126.4 110.7 112.7 109.5 124.0 119.1 120.0 109.5 117.0 107.8 110.7 112.7 109.5 114.0 122.0 130.9 130.9 130.9
kr, kJ·mol−1·Å−2
NA−CW−CW NA−CR−CCMb dihedral b
2654 4597 2847 3574 3992 4351 2512 2512 2303 2847 2654 2820 2847 2243 kθ, kJ·mol−1·rad−2
107.1 125.1 V1, kJ·mol−1
CW−NAp−C1p−CM1b −7.1535 CR−NAp−C1p−CM1b −5.2691 NAp−C1p−CM1−CM 1.2991 H1p−C1p−CM1−CM 0.0 NAp−C1p−CM1−HMb 0.6484 H1p−C1p−CM1−HM 0.0 CR−NA−C1i−H1i 2.7941 X−CM−CM−Xc,e 0.0 X−NA−CW−Xb,e 0.0 X−NA−CR−Xb,e 0.0 X−CW−CW−Xd,e 0.0 CR−NA−C*−H*b,d 0.0 CW−NA−C*−H*b,d 0.0 CR−NA−C*−C*b,d −5.2690 CW−NA−C*−C*b,d −7.1540 NA−C*−C*−H*b,d 0.0 NA−C*−C*−C*b,d −7.4800 NA−CW/R−CCM−C*b,d 0.0 NA−CW/R−CCM−H*b,d 0.0 CW−CW−CCM−H1b 0.0 CW−CW−CCM−C*b,d 0.0 CW/R−C*−C*−H*b,d 0.0 CW/R−C*−C*−C*b,d 5.4428 C*−C*−C*−C*d 7.2800 C*−C*−C*−H*d 0.0 H*−C*−C*−H*d 0.0 improper V1, kJ·mol−1
628 628 293 293 585 585 585 585 314 418 335 586 293 293 293 293 276 314 488 293 528 585 585 628 293
kθ, kJ·mol−1·rad−2
θ0, deg
angles
X−CM−X−Xc,e X−NA−X−Xb,e X−CW/R−X−Xb,e
0.0 0.0 0.0
V2, kJ·mol−1
585 585 V3, kJ·mol−1
6.1064 0.0 −3.6108 0.0 1.6600 0.0 5.9164 58.6152 12.5500 19.4600 44.9800 0.0 0.0 0.0 6.1060 0.0 3.1640 0.0 0.0 0.0 0.0 0.0 −0.2093 −0.6569 0.0 0.0 V2, kJ·mol−1
0.7939 0.0 −5.9253 −1.5575 −0.1942 1.3314 0.0 0.0 0.0 0.0 0.0 0.0 0.5190 0.0 0.7940 0.0 −1.2030 0.0 0.0 0.0 0.0 1.9343 0.8374 1.1673 1.5313 1.3305 V3, kJ·mol−1
125.604 8.370 9.200
0.0 0.0 0.0
a
The data in bold was established in this work. bNA, CR, and CW, valid for any group attached to these atoms. cCM, valid for CM and CM1. dC* represents a generic aliphatic carbon, C1, C1i, C1p, CCM, CE, C2, CS, or CT. H* represents H1, HC, H1p, or H1i (to be used unless specifically included in the list). eX denotes any atom.
respectively.32 A small survey to test the effect of these atomic radii in the computed APCs was performed by increasing and decreasing the iodine value by 5%. The obtained results indicated a variation ≤0.01 e in the charges of the carbon and iodine atoms at the CR position of the imidazolium ring. This value is well within the APC variations normally accepted during the development of a generic force field, thus supporting the use of the radii values previously recommended. More details regarding the APC calculations are given as Supporting Information. From the analyses of the obtained data, it was found that, when an isopropyl group replaces an alkyl chain connected to a NA atom, the APC of the imidazolium ring remains unchanged. However, the presence of double bonds in the alkyl chain leads to a modification of the NA charge from 0.15 e to 0.08 e (Table 1). A trend similar to that of the latter case is found when the hydrogen atoms of the imidazolium ring are replaced by halogen elements, X. In this case, the ab initio results indicated that only the APC, q, of the carbon atom
Information). For bromide, it was found that the previously recommended parameters19 accurately fitted the BHM function, and for this reason, they were kept in this work. The atomic point charges were established considering those previously reported for ionic liquids in the CL&P force field,12,14,15 the general parametrization for organic compounds, OPLS-AA force field,13 and quantum chemical calculations performed with Gaussian 09.20 The APC calculations were obtained by the ChelpG methodology,21 applied to MP2 results.22−25 The basis set aug-cc-pVDZ26−29 was used for all molecules, except if they contained iodine. For the latter case, the DGDZVP30,31 full-electron basis set was selected for all atoms. The APC’s calculations were performed over fully optimized geometries at the corresponding level of theory. The implementation of the ChelpG methodology was based on the default atomic radius for carbon, hydrogen, nitrogen, and chlorine, in Gaussian 09. In the case of bromine and iodine atoms, the values of 2.095 and 2.250 Å were selected, C
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
corresponded to the energy minimum at each evaluated dihedral, three independent PES were made, where the initial conformation of the substituent was modified. The starting molecular conformations investigated in this study are given as Supporting Information. (ii) A similar procedure was then performed using the molecular dynamics force field with the target dihedral Fourier constants set to zero. For this purpose, a series of MD runs performed using DL_POLY40 and the dihedral constrained at each angle were used to anneal the molecule from a temperature of 1000 to 10 K. The annealing process consisted of 40 ps simulation runs with time step of 1 fs, performed at 1000 K, 500 K, 100 K, and 10 K. The final conformation energy was then assumed to correspond to the energy of the molecule with the annealed conformation, after a 40 ps “zero” temperature molecular dynamics run, performed as defined in DL_POLY (i.e., a low temperature simulation where the molecule moves according to the forces and torques, but their velocity cannot exceed that corresponding to a temperature of 10 K). Five independent runs were performed to obtain the energy minimum at each dihedral angle. (iii) The energy difference between the DFT and MD results was fitted to the Fourier series in eq 1. (iv) To check the validity of the obtained coefficients, procedure ii was repeated using the newly fitted data, and the obtained energies were compared with the DFT calculations. This result is presented in Figure 2.
attached to the halogen changes, and it can be assumed that qCW/R − qHCW/R = qCWX/RX − qXCW/R (see examples in Figure S2 of the Supporting Information). An important consequence of the presence of halogens in organic compounds, as discussed above, is the possibility for the formation of halogen bonds.33 These structures result from the appearance of an electrophilic region (σ-hole) along the R−X axis (where R can be any atom attached to X), at the negatively charged halogen element. This allows the halogen atoms to form intermolecular bonds with nucleophilic atoms or group atoms (e.g., anions).34 To include this feature in classical MD simulations, Jorgensen and Schyman35 proposed the inclusion of this effect into the OPLS-AA force field, by adding a positive charge at the position of the σ-hole. In the present work, as will be demonstrated below, it was found that this modification is not necessary for the imidazolium ions. In fact, as can be observed in Table 1, because the halogen atoms are bonded to the positively charged cation rings, the APCs of these atoms vary from +0.09 e to +0.25 e. Thus, unlike in neutral molecules, where the APCs of halogen elements are negative and, for this reason, mainly nucleophilic, the fact that in imidazolium derivatives their charges are positive allows a single charge to be used to describe their electrostatic field. The substitution of the hydrogen atoms bonded to CR and CW by alkyl chains has, however, a larger impact in the APCs distribution of the imidazolium ring. In this case, in addition to a change in the charges of the CR and CW atoms, the APC of the atoms at the α position is also affected. This leads to the existence of five possible APC distributions, according to the number of hydrogen substitutions by alkyl groups (Figure 1). In a previous work, a force field for imidazolium cations with a methyl group at the CR position was proposed.36 In the sequence of the current work, a typo in the APC distribution previously reported was identified for the first time. In Chart 1 of the previous work, an HC atom (with charge 0.06 e) is attached to the CCM atom (the carbon atoms of the methyl group; see Figure 1). By using this APC distribution, a net charge of +0.79 for the cation is obtained, instead of the expected +1. To obtain the correct value, instead of an HC, an H1 atom should be used. Thus, the typo in the previous communication36 was caused by an error in the preparation of the labels in Chart 1. It should also be highlighted that the APC distribution in the alkyl chain is independent of the carbon atom position to which they are attached (CR or CW). Finally, the parametrization proposed in this work is completely compatible with the OPLS-AA force field.13 This is depicted in Figure 1, by the blue regions that demark the molecular segment where this force field can be used. Intramolecular Potential. The intramolecular force field for the imidazolium derivatives is given in Table 2. Apart from the isopropyl and propene derivatives, all data included in the table was compiled from the existent parametrization in the OPLS-AA13,16 and CL&P14,36 force fields. However, the available dihedral potentials were insufficient to describe the rotation along the isopropyl NA−C1i bond and propene C1p− CM bond. These torsion potentials were, therefore, computed from potential energy surface (PES) scan results, obtained from density functional theory (DFT) calculations, as follows: (i) Initially a PES for the target dihedral was computed at the B3LYP/6-31G(d,p) level of theory37−39 using Gaussian 09,20 determining the configuration energy of the molecule with the dihedral angle incremented of 20°. To attest that the obtained electronic energies found for the propene substituents
Figure 2. Comparison between the potential energy surface scan results obtained from DFT (blue dots) and MD (red line) calculations using the final force field, for (a) isopropyl and (b) propene substituents. D
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
Figure 3. Compounds investigated for evaluation of the atomistic force field parametrization presented in this work. The code in parentheses refers the Cambridge Structural Database reference code.43 The asterisk after the reference code indicates the compounds that can form halogen bonds. In the images, the following color code was used: carbon, gray; hydrogen, white; nitrogen, blue; chlorine, green; bromine, brown; iodine, purple; fluorine, light green; sulfur, yellow.
C1p−CM1−HM as in the OPLS-AA force field.13 The obtained parameters are given in bold in Table 2, and the comparison of the DFT data with those computed using the complete force field are given in Figure 2b. For the bonds containing halogen atoms, the parametrization was adapted from the data available in the OPLS-AA force field for halobenzene compounds.16 The harmonic force constants of the angle functions were also retrieved from this force field, but the equilibrium angles changed to match those observed in the imidazolium ring.
For an isopropene substituent attached to a NA atom of an imidazolium ring, four types of dihedrals including the NA−C1i bond can be defined. While verifying the rotation barrier along this bond, it was found that the parameters previously used for CR/W−NA−C1−H1 and CR/W−NA−C1−CE/S14 poorly describe the energy variation of this substituent (see comparison between the MD results before the reparameterization and the DFT data, given as Supporting Information). To correct this problem, the dihedral CR −N A−C1i−H 1i, that was previously set with all Fourier coefficients as zero, was reparametrized in this work. For this purpose, the energy variation was fitted to the Fourier series in eq 1, leading to the parameters given in bold in Table 2. The comparison of the DFT data with those computed using the complete force field is shown in Figure 2a. This figure indicates that the refitted parameters can now accurately describe the DFT data. A similar verification was performed for the propene substituent bonded to the imidazolium NA atom. In this case, a correction for the dihedral torsions around the C1p−CM1 bonds was found necessary, relative to the previously published values. For this purpose, the torsion angles of the dihedrals NA−C1p−CM1−CM and NA−C1p−CM1−HM were reparameterized, keeping the Fourier constants of the angles H1p−C1p−CM1−CM and H1p−
■
FORCE FIELD VALIDATION The validation of a force field can be achieved comparing experimental data (e.g., liquid density, cohesive energies, and single crystal structures) with predictions made with the proposed parametrization.14,41,42 To the best of our knowledge, no experimental data for liquid densities and cohesive energies are available for the compounds evaluated in this work. But, a survey at the Cambridge Structural Database (CCD 5.37)43 revealed 22 single crystal structures of ionic liquids containing the functional groups investigated in this work. The validation was, therefore, made through the evaluation of the proposed force field to reproduce these single crystal structures. By using E
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
Table 3. Comparison between the Experimental and Computed Parameters Obtained for the Single Crystal Structures Available in the Cambridge Structural Database43 CSD ref. code
compound
space group
T, K
a, Å
b, Å
c, Å
α, deg
β, deg
γ, deg
ρ, g·cm−3
ref
8.284 8.501 7.847 7.580 11.537 11.407 16.162 15.830 8.450 8.544 9.674 9.756 9.206 8.927 7.134 6.695 7.367 7.430 6.350 6.487 23.614 23.585 7.543 7.671 7.429 7.626 12.706 12.439 6.729 6.585 15.044 15.571 7.672 7.390 7.587 7.464 6.250 5.978 8.364 8.431 11.188 10.957 9.560 9.342
16.116 16.367 14.858 15.251 8.630 9.241 7.080 7.120 12.057 11.742 13.830 13.903 9.614 9.018 7.465 7.800 10.002 10.079 7.638 7.500 23.614 23.352 15.549 15.448 9.206 9.452 15.103 15.051 13.147 13.291 10.159 10.158 8.018 8.322 10.736 10.770 9.775 10.097 13.660 13.252 12.405 12.050 14.566 14.672
26.336 25.444 10.918 11.181 15.377 14.781 9.592 9.750 15.132 15.206 13.769 13.167 10.692 10.982 8.789 8.861 16.077 16.275 8.359 8.414 18.837 19.524 12.755 12.945 20.878 19.878 7.846 8.390 16.047 16.258 17.765 17.977 12.327 12.208 11.185 11.294 8.362 7.930 15.139 15.132 15.472 16.117 11.330 11.283
90.00 90.00 90.00 90.00 90.00 89.95 90.00 90.02 90.00 90.00 90.00 90.00 101.78 99.40 90.00 90.00 90.00 89.87 90.00 89.93 90.00 89.69 90.00 90.07 90.00 90.00 90.00 90.00 90.00 90.00 90.00 85.62 93.36 88.66 90.00 90.00 90.00 90.00 90.00 90.00 90.00 90.00 90.00 90.00
90.00 90.00 98.60 97.82 107.14 107.27 90.00 90.00 101.21 101.65 106.98 108.63 105.43 106.31 92.13 86.48 90.00 90.00 90.93 86.54 90.00 90.09 90.00 92.54 90.00 90.00 90.00 90.00 90.00 90.00 96.37 95.74 105.94 107.38 90.00 90.00 108.91 107.45 103.21 102.73 90.00 90.00 97.77 96.30
90.00 90.00 90.00 90.00 90.00 89.96 90.00 90.00 90.00 90.00 90.00 90.00 111.35 106.58 90.00 90.00 90.00 90.02 90.00 90.03 90.00 90.47 90.00 89.93 90.00 90.01 90.00 90.00 90.00 90.00 90.00 90.50 116.97 116.56 90.00 90.00 90.00 90.00 90.00 90.00 90.00 90.00 90.00 90.00
2.430 2.413 2.237 2.199 2.270 2.330 2.196 2.193 1.959 1.983 2.070 2.155 2.241 2.287 2.080 2.106 2.113 2.054 1.651 1.638 1.585 1.548 1.712 1.671 1.375 1.371 1.915 1.836 1.300 1.297 1.237 1.183 2.124 2.122 1.495 1.501 1.719 1.819 1.654 1.689 1.646 1.661 1.429 1.453
45 MD 46 MD 47 MD 48 MD 49 MD 45 MD 49 MD 50 MD 48 MD 51 MD 52 MD 52 MD
Z′/Z
AGOCIK
[C2HC1II im][tf2N]
Pbca
173
1/8
VUSNIH
[C1HC3BrBr im][I]
P21/n
173
1/4
PUZZIU
[C1HC3II im][PF6]
P21/c
173
1/4
YOXWAK
[C1HC3BrBr im][Br]
Pnma
173
0.5/4
VURKUP
[C1HC4BrBr im][tf]
P21/c
173
1/4
AGOCOQ
[C1HC2BrBr im][tf2N]
P21/c
173
1/4
VURLAW
[C1HC4II im][tf]
P1̅
173
1/2
WOMYED
[C1HC1ClCl im][I]
P21
100
1/2
YOXVUD
[C1HC4BrBr im][Br]
P212121
173
1/4
HIPYUB
[C1HC1ClCl im][Cl]
P21
133
1/2
CIHJEK
[C1Br(C3H5)BrBr im][BF4]
I41/a
173
2/32
CIHJIO
[C2Br(C3H5)BrBr im][PF6]
P21/n
173
1/4
AWENOF
[iC3 Br iC3C1C1 im][Br]
P212121
173
1/4
AWEPEX
[iC3 I iC3C1C1 im][I]
Pnma
293
0.5/4
FALCIF
[iC3 Br iC3C1C1 im][NO3]
Cmcm
223
0.25/4
XOMMER
[iC3 Cl iC3C1C1 im][Cl]
C2/c
173
1/8
HARXAZ
[C2 I C2C1C1 im][I]
P1̅
150
1/2
JIRCAP
[C1H iC3HH im][Br]
P212121
110
1/4
KEXBEV
[C1H(C3H5)HH im][I]
P21
100
1/2
EJIVAV
[C1C1C1C1C1 im][tf2N]
P21/c
173
1/4
SOTWED
[C1C1C1C1C1 im][I]
Pbca
293
1/8
IGAZAR
[iC3C2iC3C1C1 im][I]
Cc
173
1/4
this procedure, it is intended to evaluate if the interaction forces and atomic sizes rising from the force field can reproduce and maintain the spatial arrangements of the molecules in crystal structures, thus indicating the proposed parametrization is able to realistically model the ions in this class of ionic liquids. All simulations were performed using DL_POLY 4.08,40 with input files prepared with DLPGEN 2.0.2.41 A cutoff distance of 1.5 nm was used for the Coulomb and van der Waals interactions, and Ewald summation corrections were applied to account for the electrostatic interactions beyond this limit. All solid simulations were performed at the temperatures reported in the single crystal structure datasheet and at a pressure of 0.1 MPa. The temperature and pressure were controlled using a Nosé−Hoover thermostat and barostat, with relaxation time constants of 1 and 4 ps, respectively. An anisotropic
53
MD 53 MD 54 MD 55 MD 56 MD 57 MD 58 MD 59 MD 60 MD 61 MD
isothermal−isobaric ensemble (N−σ−T) was considered in all solid simulations. Equilibration stages of 1 ns were used before the production stage of 4.0 ns. A time step of 2 fs was applied in all runs. Because the initial configurations were close to equilibrium, it was previously demonstrated that these simulation times are sufficient to give reliable results.41 All initial conformations were prepared by staking several crystal unit cells of each phase, to obtain approximately squared boxes with at least 4 nm sides. The supercell dimensions used to prepare the simulation boxes are given in Table S1 of Supporting Information. Figure 3 gives the image of an ion pair retrieved for each compound found in CSD using the program Mercury 3.8.44 Also included in this figure is the corresponding reference code in the database and the compound acronym attributed in this F
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation work. The nomenclature used in the imidazolium cations was as follows: (i) the name contains 6 entries, five substituents attached to the imidazolium ring plus the termination “im”, that identifies the cation; (ii) the first substituent is bonded to one of the NA atoms, followed by that attached to CR and so on; (iii) alkyl chains are defined by the number of carbon atoms (e.g., methyl = C1; ethyl = C2); (iv) the isopropyl and propene substituents are specified as iC3 and C3H5, respectively; and (v) the hydrogen and halogen atoms are indicated by their corresponding atomic symbol. The obtained MD results for the 22 crystal structures are compared in Table 3 with the corresponding experimental data. Also given in this table are the experimental details of the single crystal structures (temperature, space group, and number of molecules in the asymmetric unit, Z′, and unit cell, Z). The deviations between theoretical and experimental data are given in Table S2 of Supporting Information. The analysis of the results reveals a good agreement between the data: (i) on average, a deviation of −0.1% is found between the MD results obtained using the proposed parametrization and the experimental data; (ii) from the 176 evaluated parameters in Table 1, only 8 exhibit a deviation larger than 4% to the experimental and, if this criterion is reduced to deviations larger than 2%, only 25 values are outside this threshold; and (iii) no clear relation is observed between the data that exhibit larger deviations, which could suggest systematic problems with the parametrization. It should also be pointed out that the absolute average (2.4%) and maximum (7.1%) differences to the experimental values found in this work are similar to those previously observed while using the OPLS-AA force field, to reproduce crystal structural data (average and maximum deviations of 2.3 and 7%, respectively).41 This result suggests, therefore, a good reliability of the proposed force field. As mentioned above, the inclusion of halogen atoms in the imidazolium ring enables the possibility for the formation of halogen bonds between the halogen atom and, e.g., the anions. In the test set used in this work, 17 of the selected structures can form this sort of interaction (please check Figure 3, for the compounds marked with an asterisk). The good agreement between the theoretical and the experimental data of the single crystal structures of these materials suggests, therefore, that the present parametrization can reproduce this sort of interaction, without the need for the inclusion in the model of a “σ-hole”, as recently proposed.35 In fact, the analysis of the electrostatic potential mapped around the imidazolium derivatives, as illustrated in Figure 4 for [iC3 Br iC3C1C1 im]+, reveals that, in the case of these compounds, no σ-hole is present in the region of the bromine atom. Thus, in the case of halogenated imidazolium cations, no conventional halogen bonds are expected to be observed, i.e., interactions between electronegative atoms and the positive σ-hole present in the halogen atom of a molecules.34 As mentioned above, experimental data, like liquid densities, for the compounds modeled in this work are scarce or nonexistent. Thus, the proposed force field can be used to make some predictions about the liquid structure in these substances. Particularly interesting is to evaluate if the structural organization found in the crystal structure is preserved in the liquid state, namely, their ability to form halogen interactions. This study was performed for the ionic liquids [C1HC1ClCl im][I] and [iC3 Br iC3C1C1 im][Br]. For this purpose, 300 ion pairs were randomly distributed inside a simulation box using Packmol,62 and simulations were made using GROMACS
Figure 4. Electrostatic potential mapped onto an electron density isosurface, obtained at the MP2/aug-cc-pVDZ level of theory for [iC3 Br iC3C1C1 im]+. In the images, the following color code was used: carbon, gray; hydrogen, white; nitrogen, blue; bromine, brown.
5.1.4.63 All input files for these programs were prepared using DLPGEN.41 A cutoff of 1.6 nm was used for all types of interactions. Beyond this boundary, the particle-mesh Ewald technique was used to correct the electrostatic interactions. The temperature was maintained at 300 K using a v-rescale thermostat (relaxation time constant of 1 ps), and the pressure was controlled at 0.1 MPa, using a Parrinello−Rahman barostat (relaxation time of 5 ps; compressibility of 4.5 × 10−5). The simulation boxes were initially equilibrated by performing several 5 ns simulation runs, until a constant density had been achieved. Finally, a production run was performed, collecting the trajectory of the molecules each 100 ps for 40 ns. Figure 5 gives the comparison of the radial distribution function (RDF) between the imidazolium bromine atom and the bromide anion of [iC3 Br iC3C1C1 im][Br], in the crystal and liquid phases. Due to the similarity between the results
Figure 5. Comparison between the radial distribution functions obtained between the imidazolium bromine atom and the bromide anion of the [iC3 Br iC3C1C1 im][Br] ionic liquid, computed from MD simulation results of the crystalline material at 173 K (CSD reference code AWENOF) and the liquid phase at 300 K. G
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
(2) Ouvrard, C.; Le Questel, J. Y.; Berthelot, M.; Laurence, C. Halogen-Bond Geometry: a Crystallographic Database Investigation of Dihalogen Complexes. Acta Crystallogr., Sect. B: Struct. Sci. 2003, 59, 512−526. (3) Legon, A. C. The Halogen Bond: an Interim Perspective. Phys. Chem. Chem. Phys. 2010, 12, 7736−7747. (4) Metrangolo, P.; Neukirch, H.; Pilati, T.; Resnati, G. Halogen Bonding Based Recognition Processes: A World Parallel to Hydrogen Bonding. Acc. Chem. Res. 2005, 38, 386−395. (5) Parisini, E.; Metrangolo, P.; Pilati, T.; Resnati, G.; Terraneo, G. Halogen Bonding in Halocarbon-Protein Complexes: a Structural Survey. Chem. Soc. Rev. 2011, 40, 2267−2278. (6) Politzer, P.; Murray, J. S.; Clark, T. Halogen Bonding: an Electrostatically-Driven Highly Directional Noncovalent Interaction. Phys. Chem. Chem. Phys. 2010, 12, 7748−7757. (7) Ye, C. F.; Shreeve, J. M. Syntheses of Very Dense Halogenated Liquids. J. Org. Chem. 2004, 69, 6511−6513. (8) Tan, K. L.; Li, C. X.; Meng, H.; Wang, Z. H. Improvement of Hydrophobicity of Ionic Liquids by Partial Chlorination and Fluorination of the Cation. Chin. J. Chem. 2009, 27, 174−178. (9) Fumino, K.; Reichert, E.; Wittler, K.; Hempelmann, R.; Ludwig, R. Low-Frequency Vibrational Modes of Protic Molten Salts and Ionic Liquids: Detecting and Quantifying Hydrogen Bonds. Angew. Chem., Int. Ed. 2012, 51, 6236−6240. (10) Li, H. Y.; Lu, Y. X.; Wu, W. H.; Liu, Y. T.; Peng, C. J.; Liu, H. L.; Zhu, W. L. Noncovalent Interactions in Halogenated Ionic Liquids: Theoretical Study and Crystallographic Implications. Phys. Chem. Chem. Phys. 2013, 15, 4405−4414. (11) Skarmoutsos, I.; Dellis, D.; Matthews, R. P.; Welton, T.; Hunt, P. A. Hydrogen Bonding in 1-Butyl- and 1-Ethyl-3-methylimidazolium Chloride Ionic Liquids. J. Phys. Chem. B 2012, 116, 4921−4933. (12) Canongia Lopes, J. N.; Padua, A. A. H. CL&P: A Generic and Systematic Force Field for Ionic Liquids Modeling. Theor. Chem. Acc. 2012, 131, 1129. (13) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (14) Canongia Lopes, J. N.; Deschamps, J.; Padua, A. A. H. Modeling Ionic Liquids Using a Systematic All-Atom Force Field. J. Phys. Chem. B 2004, 108, 2038−2047. (15) Canongia Lopes, J. N.; Padua, A. A. H.; Shimizu, K. Molecular Force Field for Ionic Liquids IV: Trialkylimidazolium and Alkoxycarbonyl-imidazolium Cations; Alkylsulfonate and Alkylsulfate Anions. J. Phys. Chem. B 2008, 112, 5039−5046. (16) Jorgensen, W. L.; Ulmschneider, J. P.; Tirado-Rives, J. Free Energies of Hydration From a Generalized Born Model and an ALLAtom Force Field. J. Phys. Chem. B 2004, 108, 16264−16270. (17) Fumi, F. G.; Tosi, M. P. Ionic Sizes and Born Repulsive Parameters in NaCl-Type Alkali Halides - I: the Huggins-Mayer and Pauling Forms. J. Phys. Chem. Solids 1964, 25, 31−43. (18) Tosi, M. P.; Fumi, F. G. Ionic Sizes and Born Repulsive Parameters in NaCl-Type Alkali Halides - II: the Generalized HugginsMayer Form. J. Phys. Chem. Solids 1964, 25, 45−52. (19) Canongia Lopes, J. N.; Padua, A. A. H. Molecular Force Field for Ionic Liquids III: Imidazolium, Pyridinium, and Phosphonium Cations; Chloride, Bromide, and Dicyanamide Anions. J. Phys. Chem. B 2006, 110, 19586−19592. (20) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.;
obtained for this ionic liquid and those of [C1HC1ClCl im][I], the data for the latter substance are given as Supporting Information. The RDF obtained for the crystalline material (Figure 5) reflects the regularity of the arrangement of the atoms (and molecules) in the solid state, with the formation of direct interactions between the halogen atoms of the cation and the anion Br− (peak centered at approximately 0.35 nm). The comparison of this data with the RDF obtained in the liquid phase indicates that much of the local organization found in the solid is preserved in the fluid phase, including the presence of strong and specific Br···Br interactions.
■
CONCLUSIONS In the case of ionic liquids, halogen bond interactions can be modeled using a “natural” extension of the currently used CL&P force field, without the need to parametrize ad-hoc σhole interaction sites. The systematic extension of the CL&P force field to halogenated and otherwise substituted imidazolium cations was successfully validated against 22 different crystalline structures deposited at the CSD. Molecular Dynamics simulations in the liquid phase show that many of the structural features present in the solid are still present in the liquid, including the existence of halogen−halide interactions.
■
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.7b00645. Tables S1 and S2 give additional data regarding the MD simulations. Figure S1 contains the scheme of the cations used in the ab initio calculation to obtain APCs. Figure S2 shows the effect of the APC distribution and of the substitution of the hydrogens in the imidazolium ring by halogen atoms. Figure S3 with the starting molecular conformations used in the DFT potential energy surface scans. Figure S4 shows the comparison between MD PES scan results before dihedral parametrization and the DFT data. Figure S5 shows the fitting results of the Tosi and Fumi functions for halogen anions to a 12-6 LennardJones potential. Figure S6 contains the radial distribution functions obtained between the imidazolium chlorine atoms and the iodide anion of the [C1HC1ClCl im][I] ionic liquid (PDF)
■
AUTHOR INFORMATION
Corresponding Author
*(C.E.S.B.) E-mail:
[email protected]. ORCID
Carlos E. S. Bernardes: 0000-0003-1490-9728 Funding
This research was supported by Fundaçaõ para a Ciência e Tecnologia (FCT) through project UID/QUI/00100/2013. C. E. S. Bernardes acknowledge postdoctoral grant SFRH/BPD/ 101505/2014. Notes
The authors declare no competing financial interest.
■
REFERENCES
(1) Bent, H. A. Structural Chemistry of Donor-Acceptor Interactions. Chem. Rev. 1968, 68, 587−648. H
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation
(41) Bernardes, C. E. S.; Joseph, A. Evaluation of the OPLS-AA Force Field for the Study of Structural and Energetic Aspects of Molecular Organic Crystals. J. Phys. Chem. A 2015, 119, 3023−3034. (42) Bernardes, C. E. S.; Canongia Lopes, J. N.; da Piedade, M. E. M. All-Atom Force Field for Molecular Dynamics Simulations on Organotransition Metal Solids and Liquids. Application to M(CO)n (M = Cr, Fe, Ni, Mo, Ru, or W) Compounds. J. Phys. Chem. A 2013, 117, 11107−11113. (43) Allen, F. H. The Cambridge Structural Database: a Quarter of a Million Crystal Structures and Rising. Acta Crystallogr., Sect. B: Struct. Sci. 2002, 58, 380−388. (44) Macrae, C. F.; Edgington, P. R.; McCabe, P.; Pidcock, E.; Shields, G. P.; Taylor, R.; Towler, M.; van de Streek, J. Mercury: Visualization and Analysis of Crystal Structures. J. Appl. Crystallogr. 2006, 39, 453−457. (45) Mukai, T.; Nishikawa, K. 4,5-Dihaloimidazolium-Based Ionic Liquids: Effects of Halogen-Bonding on Crystal Structures and Ionic Conductivity. RSC Adv. 2013, 3, 19952−19955. (46) Mukai, T.; Nishikawa, K. Zigzag Sheet Crystal Packing in a Halogen-bonding Imidazolium Salt: 1-Butyl-4,5-dibromo-3-methylimidazolium Iodide. X-Ray Struct. Anal. Online 2010, 26, 31−32. (47) Mukai, T.; Nishikawa, K. Halogen Bonding and Hydrogen Bonding in 4,5-Diiodo-3-methyl-1-propylimidazolium Hexafluorophosphate. X-Ray Struct. Anal. Online 2010, 26, 39−40. (48) Mukai, T.; Nishikawa, K. Halogen-Bonded and HydrogenBonded Network Structures in Crystals of 1-Propyl- and 1-Butyl-4,5dibromo-3-methylimidazolium Bromides. Chem. Lett. 2009, 38, 402− 403. (49) Mukai, T.; Nishikawa, K. Syntheses and Crystal Structures of Two Ionic Liquids With Halogen-Bonding Groups: 4,5-dibromo- and 4,5-diiodo-1-butyl-3-methylimidazolium Trifluoromethanesulfonates. Solid State Sci. 2010, 12, 783−788. (50) Hindi, K. M.; Siciliano, T. J.; Durmus, S.; Panzner, M. J.; Medvetz, D. A.; Reddy, D. V.; Hogue, L. A.; Hovis, C. E.; Hilliard, J. K.; Mallet, R. J.; Tessier, C. A.; Cannon, C. L.; Youngs, W. J. Synthesis, Stability, and Antimicrobial Studies of Electronically Tuned Silver Acetate N-heterocyclic Carbenes. J. Med. Chem. 2008, 51, 1577−1583. (51) Scheele, U. J.; Dechert, S.; Meyer, F. A Versatile Access to Pyridazines With Tethered Imidazolium Groups - New Precursors for Mono- and Binucleating NHC/Pyridazine Hybrid Ligands. Tetrahedron Lett. 2007, 48, 8366−8370. (52) Froschauer, C.; Kahlenberg, V.; Laus, G.; Schottenberger, H. Halogen Interactions in 2,4,5-Tribromoimidazolium Salts. Crystals 2012, 2, 1017−1033. (53) Kuhn, N.; Abu-Rayyan, A.; Eichele, K.; Schwarz, S.; Steimann, M. Weak Interionic Interactions in 2-Bromoimidazolium Derivatives. Inorg. Chim. Acta 2004, 357, 1799−1804. (54) Kuhn, N.; Richter, M.; Steimann, M.; Strobele, M.; Sweidan, K. Hydrogen Bonding in Imidazolium Nitrates. Z. Anorg. Allg. Chem. 2004, 630, 2054−2058. (55) Kuhn, N.; Abu-Rayyan, A.; Gohner, M.; Steimann, M. Synthesis and Crystal Structure of the Dichlorine Adduct of 2,3-dihydro-1,3diisopropyl-4,5-dimethylimidazol-2-ylidene. Z. Anorg. Allg. Chem. 2002, 628, 1721−1723. (56) Kuhn, N.; Kratz, T.; Henkel, G. A Stable Carbene Iodine Adduct: Secondary Bonding in 1,3-Diethyl-2-Iodo-4,5-Dimethylimidazolium Iodide. J. Chem. Soc., Chem. Commun. 1993, 1778−1779. (57) Golovanov, D. G.; Lyssenko, K. A.; Vygodskii, Y. S.; Lozinskaya, E. I.; Shaplov, A. S.; Antipin, M. Y. Crystal Structure of 1,3dialkyldiazolium Bromides. Russ. Chem. Bull. 2006, 55, 1989−1999. (58) Fei, Z. F.; Kuang, D. B.; Zhao, D. B.; Klein, C.; Ang, W. H.; Zakeeruddin, S. M.; Gratzel, M.; Dyson, P. J. A Supercooled Imidazolium Iodide Ionic Liquid as a Low-Viscosity Electrolyte for Dye-Sensitized Solar Cells. Inorg. Chem. 2006, 45, 10407−10409. (59) Roth, C.; Peppel, T.; Fumino, K.; Kockerling, M.; Ludwig, R. The Importance of Hydrogen Bonds for the Structure of Ionic Liquids: Single-Crystal X-ray Diffraction and Transmission and Attenuated Total Reflection Spectroscopy in the Terahertz Region. Angew. Chem., Int. Ed. 2010, 49, 10221−10224.
Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, O.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, revision D.01; Gaussian, Inc.: Wallingford, CT, 2009. (21) Breneman, C. M.; Wiberg, K. B. Determining Atom-Centered Monopoles from Molecular Electrostatic Potentials. The Need for High Sampling Density in Formamide Conformational Analysis. J. Comput. Chem. 1990, 11, 361−373. (22) Frisch, M. J.; Head-Gordon, M.; Pople, J. A. A Direct MP2 Gradient Method. Chem. Phys. Lett. 1990, 166, 275−280. (23) Frisch, M. J.; Head-Gordon, M.; Pople, J. A. Semi-Direct Algorithms for the MP2 Energy and Gradient. Chem. Phys. Lett. 1990, 166, 281−289. (24) Head-Gordon, M.; Head-Gordon, T. Analytic MP2 Frequencies without Fifth-Order Storage. Theory and Application to Bifurcated Hydrogen-Bonds in the Water Hexamer. Chem. Phys. Lett. 1994, 220, 122−128. (25) Head-Gordon, M.; Pople, J. A.; Frisch, M. J. MP2 Energy Evaluation by Direct Methods. Chem. Phys. Lett. 1988, 153, 503−506. (26) Kendall, R. A.; Dunning, T. H.; Harrison, R. J. ElectronAffinities of the 1st-Row Atoms Revisited - Systematic Basis-Sets and Wave-Functions. J. Chem. Phys. 1992, 96, 6796−6806. (27) Wilson, A. K.; Woon, D. E.; Peterson, K. A.; Dunning, T. H. Gaussian Basis Sets for Use in Correlated Molecular Calculations. IX. The Atoms Gallium Through Krypton. J. Chem. Phys. 1999, 110, 7667−7676. (28) Woon, D. E.; Dunning, T. H. Gaussian-Basis Sets for Use in Correlated Molecular Calculations. III. The Atoms Aluminum through Argon. J. Chem. Phys. 1993, 98, 1358−1371. (29) Woon, D. E.; Dunning, T. H. Gaussian-Basis Sets for Use in Correlated Molecular Calculations. IV. Calculation of Static Electrical Response Properties. J. Chem. Phys. 1994, 100, 2975−2988. (30) Sosa, C.; Andzelm, J.; Elkin, B. C.; Wimmer, E.; Dobbs, K. D.; Dixon, D. A. A Local Density Functional-Study of the Structure and Vibrational Frequencies of Molecular Transition-Metal Compounds. J. Phys. Chem. 1992, 96, 6630−6636. (31) Godbout, N.; Salahub, D. R.; Andzelm, J.; Wimmer, E. Optimization of Gaussian-Type Basis-Sets for Local Spin-Density Functional Calculations. Part I. Boron through Neon, Optimization Technique and Validation. Can. J. Chem. 1992, 70, 560−571. (32) Jaguar 7.7- Quick Start Guide; Schrödinger Press: New York, 2010; p 236. (33) Cavallo, G.; Metrangolo, P.; Milani, R.; Pilati, T.; Priimagi, A.; Resnati, G.; Terraneo, G. The Halogen Bond. Chem. Rev. 2016, 116, 2478−2601. (34) Cavallo, G.; Terraneo, G.; Monfredini, A.; Saccone, M.; Priimagi, A.; Pilati, T.; Resnati, G.; Metrangolo, P.; Bruce, D. W. Superfluorinated Ionic Liquid Crystals Based on Supramolecular, Halogen-Bonded Anions. Angew. Chem., Int. Ed. 2016, 55, 6300−6304. (35) Jorgensen, W. L.; Schyman, P. Treatment of Halogen Bonding in the OPLS-AA Force Field; Application to Potent Anti-HIV Agents. J. Chem. Theory Comput. 2012, 8, 3895−3901. (36) Canongia Lopes, J. N.; Padua, A. A. H.; Shimizu, K. Molecular Force Field for Ionic Liquids IV: Trialkylimidazolium and Alkoxycarbonyl-Imidazolium Cations; Alkylsulfonate and Alkylsulfate Anions. J. Phys. Chem. B 2008, 112, 5039−5046. (37) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (38) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron-Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (39) Frisch, M. J.; Pople, J. A.; Binkley, J. S. Self-Consistent Molecular-Orbital Methods 25. Supplementary Functions for Gaussian-Basis Sets. J. Chem. Phys. 1984, 80, 3265−3269. (40) Todorov, I. T.; Smith, W.; Trachenko, K.; Dove, M. T. DL_POLY_3: New Dimensions in Molecular Dynamics Simulations via Massive Parallelism. J. Mater. Chem. 2006, 16, 1911−1918. I
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
Article
Journal of Chemical Theory and Computation (60) Kuhn, N.; Henkel, G.; Kreutzberg, J. Synthesis and Reactions of 1,2,4,5-Tetramethylimidazole: the Crystal-Structure of Pentamethylimidazolium Iodide. Z. Naturforsch., B: J. Chem. Sci. 1991, 46, 1706− 1712. (61) Kuhn, N.; Gohner, M.; Steimann, M. Imidazole Derivatives. 47. On the alkylation of 2,3-dihydro-1,3-diisopropyl-4,5-dimethylimidazol2-ylidene. Z. Naturforsch. B 2002, 57, 631−636. (62) Martínez, L.; Andrade, R.; Birgin, E. G.; Martínez, J. M. PACKMOL: A Package for Building Initial Configurations for Molecular Dynamics Simulations. J. Comput. Chem. 2009, 30, 2157− 2164. (63) Abraham, M. J.; Murtola, T.; Schulz, R.; Páll, S.; Smith, J. C.; Hess, B.; Lindahl, E. GROMACS: High Performance Molecular Simulations Through Multi-Level Parallelism from Laptops to Supercomputers. SoftwareX 2015, 1−2, 19−25.
J
DOI: 10.1021/acs.jctc.7b00645 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX
