Article pubs.acs.org/JPCA

Distributed Polarizability Models for Imidazolium-Based Ionic Liquids Claude Millot,*,† Alain Chaumont,‡ Etienne Engler,‡ and Georges Wipff‡ †

Université de Lorraine, CNRS, SRSMC, UMR 7565, Equipes TMS/ReSolve, Faculté des Sciences et Technologies, Boulevard des Aiguillettes, BP 70239, Vandoeuvre-lès-Nancy F-54506, France ‡ Université de Strasbourg, CNRS, Institut de Chimie, Laboratoire MSM, UMR 7177, 1 rue Blaise Pascal, Strasbourg F-67000, France S Supporting Information *

ABSTRACT: Quantum chemical calculations are used to derive distributed polarizability models sufficiently accurate and compact to be used in classical molecular dynamics simulations of imidazolium-based room temperature ionic liquids. Two distributed polarizability models are fitted to reproduce the induction energy of three imidazolium cations (1,3dimethyl-, 1-ethyl-3-methyl-, and 1-butyl-3-methylimidazolium) and four anions (tetrafluoroborate, hexafluorophosphate, nitrate, and thiocyanate) polarized by a point charge located successively on a grid of surrounding points. The first model includes charge-flow polarizabilities between firstneighbor atoms and isotropic dipolar polarizability on all atoms (except H), while the second model includes anisotropic dipolar polarizabilities on all atoms (except H). For the imidazolium cations, particular attention is given to the transferability of the distributed polarizability sets. The molecular polarizability and its anisotropy rebuilt by the distributed models are found to be in good agreement with the exact ab initio values for the three cations and 23 additional conformers of 1-ethyl-3-methyl-, 1-butyl-3-methyl-, 1-pentyl-3-methyl-, and 1-hexyl3-methylimidazolium cations. the dipolar polarizability29,30 or to the dipolar plus quadrupolar contributions.31 However, in polyatomic ions, the molecular polarizability is generally distributed over selected sites (atoms, bonds, lone pairs, and functional groups). The resulting “distributed polarizabilities” (noted “DPs” for short) thus depend on the choice of sites, of the multipolar rank (monopole, dipole, etc.), as well as on the partitioning protocol. With regard to the partitioning procedure, several empirical approaches have been proposed to assign isotropic or anisotropic dipolar polarizabilities on atoms, bonds, or functional groups,32−44 sometimes adding intramolecular charge-flow.41 In the simple Applequist−Thole’s model, atomic isotropic dipolar polarizabilities are used, and the molecular polarizability anisotropy is rebuilt via intramolecular interactions between polarizable sites.37 Gu and Yan have used that approach to derive atomic polarizabilities for 216 cations and 80 anions of RTILs.26 The fluctuating charge model, describing polarization effects by intramolecular charge transfer, is based on the electronegativity equalization principle.45 This model originally parametrized from experimental data can also be parametrized from quantum chemical calculations. Other approaches to dissect molecular polarizability are based on a direct derivation of distributed polarizabilities from

1. INTRODUCTION In the last 15 years, the modeling of room temperature ionic liquids (RTILs) has rapidly developed,1−3 generally based on molecular dynamics (MD) simulations using two-body additive force fields for convenience.4−8 However, electron reorganization effects like polarization and charge transfer can be important to fine-tune the calculated physicochemical properties of ionic liquids. For instance, the explicit inclusion of polarizability in the force field of RTILs has been found to reduce viscosity and to increase ionic diffusion, improving the agreement with experimental results.9−11 Likewise, charge reduction of ions from ±1 e to smaller values, as suggested by ab initio static or dynamic [ab initio molecular dynamics (AIMD)] studies, improves the description of water dragging by hydrophobic RTILs at aqueous interfaces.12,13 This strategy of charge reduction has also been used in simulation studies of RTILs14−20 and compared to the use of a polarizable model.21 Both approaches can lead to good results for collective properties, but the reduced charge approach is not so efficient for describing interactions at short range.21 As regards to the inclusion of polarization effects in MD simulations of RTILs or molten salts, there are so far several reports, generally using atomic dipolar polarizabilities, corresponding to the standard protocol22 or to the Thole’s model.9,11,23−27 In other MD studies, polarization has been mimicked by the charge on springs model (Drude oscillator or shell model).10,21,28 For monatomic ions, the “molecular polarizability” corresponds to © 2014 American Chemical Society

Received: June 4, 2014 Revised: August 8, 2014 Published: August 18, 2014 8842

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

where Q̂ sl11,m1 is the operator of an electric multipole moment of a site s1 located at the position s1, |0> is the molecular electronic ground-state with a total energy E0 and |n > ≠ |0> is an excited state with energy En. When several sites are chosen, one obtains charge−charge (l1 = l2 = 0), charge−dipole (l1 = 0, l2 = 1), dipole−dipole (l1 = l2 = 1), ..., local (s1 = s2), and nonlocal (s1 ≠ s2) polarizabilities. Instead of complex multipole moments Q̂ sljj,mj (mj = −lj,...,lj) in eq 1, their real counterpart defined with respect to a frame linked to the site sj are used in practice:

quantum chemical (QM) calculations. Some use the sum over states formula with localized molecular orbitals to define local polarizabilities.46−48 Another method, also based on localized molecular orbitals, uses a polarization propagator approach.49 The Stone general formalism describes the reorganization of charge distribution at selected polarizable sites induced by an external potential and its successive derivatives.50,51 By including local and nonlocal polarizabilities, this model affords a detailed description of the polarization effects, including intramolecular charge flow, generally absent in most other schemes. Calculations of atomic dipolar polarizabilities52,53 and DPs defined by Stone54 have been implemented in the Bader’s theory of atoms in molecules. Recently, density fitting schemes have been proposed to get DPs.55−57 The Hirshfeld scheme of partitioning electron density58 has also been used.59−61 Generally, most of the DP models yield a large number of local and nonlocal multipole−multipole (charge−charge, charge−dipole, dipole−dipole, etc.) polarizability parameters that may not be easily transferable and are furthermore CPUtime-consuming in MD simulations. Aternatively, polarizability models can be obtained from routine QM calculations via fitting of induction energies or of perturbed potential maps62−66 or by fitting of distributed multipoles of a polarized molecule.67 Perturbed potential maps have been used to parametrize polarizable force fields based on the fluctuating point charge model alone68 or with added isotropic atomic dipolar polarizabilities69,70 or using the classical Drude oscillator model.71 In this paper, we decided to fit distributed polarizabilities on induction energies and to optimize two DP models compact and simple enough to be used in MD simulations of RTILs, thereby focusing on atomic dipolar polarizability and charge flow components. The first model, noted CFDiso, includes charge-flow polarizabilities between first-neighbor atoms and isotropic dipolar polarizability on all atoms (except H), while the second one, noted Daniso, includes anisotropic dipolar polarizabilities on all atoms (except H). They are developed on typical RTIL ions, namely 1-alkyl-3-alkylimidazolium cations and four anions (tetrafluoroborate BF−4 , hexafluorophosphate PF−6 , nitrate NO−3 , and thiocyanate SCN−) . For the cations, 1,3-dimethylimidazolium (DMIM+), 1-ethyl-3-methylimidazolium (EMIM+), and 1butyl-3-methylimidazolium (BMIM+) have been chosen as test cases to establish transferable models. For each ion, induction energies are calculated by QM for different positions of the surrounding polarizing +1 e charge, and the relevant parameters are fitted for the two models. We also test the transferability of these parameters in the alkyl-imidazolium cation series.

s s s 1 Q̂ l ,jm = ( −1)mj (Q̂ l ,jm c + iQ̂ l ,jm s) j j j j j j 2

Within this convention, one obtains DPs αl1,κ1,l2,κ2(s1,s2) where κj stands for mjc or mjs, related to the local axes at site sj. The molecular polarizability αL,κ,L′,κ′(x) defined with respect to the position x is obtained by translating the DPs involving site positions sj according to the formula: αL , κ , L

=

∑ n≠0

s

2

s

2

2

1

1 1

1 1 2

2

s1, s2 l1, κ1 l 2 , κ2

QM < QM − ,QM − qV 0,QM ind, i = , i 0 i

(4)

where ,QM is the energy of the total system (ion plus point i QM charge q), , 0 is the energy of the unpolarized ion alone, and VQM 0,i is the molecular electrostatic potential at point i created by the unperturbed molecular charge distribution. With a chosen polarizability model {αl1,κ1,l2,κ2(s1,s2)}, the polarizing charge q creates induced moments ΔQsl11,κ1 at each polarizable site, s1, in the molecule:

s

2

∑ ∑ ∑ WL ,κ ,l ,κ (x − s1)αl ,κ , l ,κ (s1, s2)

(x − s 2 ) (3) ′ , κ ′ , l2 , κ2 where Wa,b,c,d(r) is a translation function defined by Stone.72 For all calculations presented here, the molecular dipolar polarizabilities are calculated at the origin of the laboratory (x, y, z) frame (see Figure 1). The atomic coordinates are given in the Supporting Information. The distributed polarizability parameters α’s are obtained from a fitting procedure aiming to reproduce induction energies of the RTIL ion polarized by a point charge q of +1e (1.602177 × 10−19 C) placed at different positions on a grid around it. The grids take into account the symmetry of the ion with a spacing in the x, y, and z directions of 2 Å. For the nitrate and the thiocyanate ions, a tighter spacing of 1.5 and 1.0 Å, respectively, has been chosen to get enough independent points. Then the grid is symmetrized to cover the whole space around the ion. Selected points sit in spherical shells around each atom with radii ranging from 3 to 5 times its van der Waals radius,73 excluding points located at less than 3 times the van der Waals radius of any atom. This mimimum distance is sufficiently large to avoid penetration and hyperpolarizability effects.64 Finally, the number Np of selected points ranges between 248 (for BF−4 ) and 768 (for BMIM+). For every point i, the induction energy, < QM ind, i , of a given ion polarized by the point charge q is calculated by QM as follows:

+ 1

(x) = ′,κ′

WL

αl1, m1, l2 , m2(s1, s 2) s

(2)

s s s 1 Q̂ l ,j−m = (Q̂ l ,jm c − iQ̂ l ,jm s) j j j j j j 2

2. DERIVATION OF DISTRIBUTED POLARIZABILITIES FROM QUANTUM CHEMICAL CALCULATIONS 2.1. Definition of Distributed Polarizabilities and Fitting Method. We use the definitions and spherical tensor notation introduced by Stone.50 The distributed polarizabilities α’s are defined by

1

mj > 0

1

ΔQ ls1, κ = −

En − E0

1 1

∑ s2 , l 2 , κ2

(1) 8843

αl1, κ1, l2 , κ2(s1, s 2)Tls22,,κi2 ,0,0q

(5)

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article model < ind, i =

1 2

ΔQ ls1, κ Tls11,,κi1,0,0q

∑

1 1

(6)

s1, l1, κ1

From eqs 5 and 6, one obtains: 1 2 model < ind, i = − q 2

αl1, κ1, l2 , κ2(s1, s 2)Tls22,,κi2 ,0,0Tls11,,κi1,0,0

∑ ∑ s1, l1, κ1 s2 , l 2 , κ2

(7)

Equations 5 to 7 assume that nonlinear effects are negligible. The double sum in eq 7 contains Nc terms corresponding to the number of nonzero unknown polarizability components (α’s). To obtain these, two algorithms have been proposed. The first one is the standard root-mean-square deviation (RMSD) algorithm.64 It finds the optimum set of Nc α’s which minimizes the RMSD between the fitted and the exact QM induction energies, over all Np points of a grid. The second one, used here, is based on a statistical analysis of α’s (SADP).65 A subset of Nc points i is first selected at random among the Np points of the total grid (with Nc ≪ Np). The corresponding set of Nc polarizability components that exactly reproduces the induction energies induced by the q charge at the i positions is obtained by solving the following linear system of Nc equations: 1 2 < QM ind, i = − q 2

∑ ∑

αl1, κ1, l2 , κ2(s1, s 2)Tls22,,κi2 ,0,0Tls11,,κi1,0,0

s2 , l 2 , κ2 s1, l1, κ1

i = 1, ..., Nc

(8)

This procedure is repeated a large number of times Ne (here, Ne = 107), and a distribution function is built for every parameter of the model, whose optimal value is taken as the most probable value of a Cauchy distribution fitting the distribution function, as done in ref 65. When an intramolecular charge flow is included in the model, as the total charge of the molecule or ion must be conserved, the charge-multipole polarizabilities of each polarizable site s1 must satisfy the constraints:

∑ αl ,κ ,0,0(s1, s2) = 0 1 1

s2

(9)

where the summation over s2 includes all sites exchanging some charge with s1, plus s1 itself. Such constraints are introduced to reduce the number of unknown α parameters in the linear system of equations (eq 8). 2.2. Selected CFDiso and Daniso Models. As our aim is to get real polarizability models simple enough to be tractable in MD simulations of RTILs, we have selected two previously considered models, 64,65 noted here CFD iso and Daniso, respectively. In both models, the polarizability in eq 7 is distributed over atomic sites (s1 = A and s2 = B). The CFDiso model contains charge-flow (CF) polarizabilities α0,0,0,0(A,B), noted αAB 00 , between first neighbors A and B, plus local isotropic dipolar (Diso) polarizabilities α1,0,1,0(A,A) = α1,1c,1,1c(A,A) = α1,1s,1,1s(A,A), noted αAdip, on nonhydrogen A atoms. The indices 1c and 1s in spherical tensor notation correspond to the Cartesian x and y components, respectively, and κi = 0 corresponds to the z component. The Daniso model includes local anisotropic dipolar polarizabilities α1,κ1,1,κ2(A,A), noted αAdip,ab, on nonhydrogen A atoms, where a and b correspond to two of the (x, y, z) axes (in the laboratory frame) or two of the (x′, y′, z′) axes (in the local frame at site A, diagonalizing the local dipolar polarizabilty tensor). See Figure 1 for the studied ions. Note that a combination of fluctuating charges and atomic

Figure 1. Anions and cations atom numbering and definition of laboratory (x, y, z) frame and local (x′, y′, z′) frames. For the central atom of a given anion, (x′, y′, z′) coincides with (x, y, z). Note that for the O(NO3) oxygens, the y′ axis is parallel to the z′ = z axis of the N atom. For each atom of the imidazolium ring, the y′ axis is perpendicular to the ring plane. For a (X−)CH3 group, the z′ axis sits along the X−C bond and a local C3v symmetry is assumed. For CH2 groups, the z′ axis bissects the H−C−H angle, and the x′ axis sits in the N−C(H2)−C or the C−C(H2)−C plane.

where, Tsl22,κ,i 2,0,0 are elements of the electrostatic tensor. The latter, multiplied by q, afford the potential (when l2 = 0), the field (with inversed sign) (when l2 = 1), or higher derivatives (when l2 > 1) created at site s2 by the charge q located at point i.74,75 The induction energy of the molecule is then obtained by 8844

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

and the anisotropic component, αF2 dip,xz. The resulting polarizabilities parallel and perpendicular to the B−F bond, obtained by diagonalization of the polarizability tensor, are tabulated in Table 1. For NO−3 , only three parameters are fitted, namely the

dipolar polarizabilities has already been used in MD simulations of liquid water.70,76 For that liquid, comparison of different models has emphasized the importance of extending polarizability terms beyond single-site dipolar contributions.76 Regarding the Daniso approach, it has been used, for instance in the molecular mechanics MM3 force field, where the polarizability is empirically partitioned in bond contributions77 or in the QM-based SIBFA force field where the polarizability is distributed over bond barycenters and lone pairs.47,78 For both models, the molecular dipolar polarizability tensor will be noted αa,b with a,b = x, y, or z. The corresponding isotropic polarizability α̅ is defined by α̅ = (αx,x + αy,y + αz,z)/3, while the anisotropy can be defined by the anisotropy factor κaniso: κaniso =

(αx

2

′,x′

2

Table 1. Fitted Distributed Polarizabilities (in a.u.) for the CFDiso and Daniso Models of BF−4 , PF−6 , NO−3 , and SCN− Anionsa CFDiso BF−4

PF−6

2

− α̅ ) + (αy , y − α̅ ) + (αz , z − α̅ ) ′ ′ ′ ′ 6α̅ 2

NO−3

(10)

In eq 10, the molecular frame (x′, y′, z′) is the frame where the molecular dipolar polarizability tensor is diagonal. Geometry optimization of the anions and of the DMIM+, EMIM+, and BMIM+ cations and all QM calculations have been performed at the MP2/6-311+G(2df) level, using the Gaussian03 package.79 For EMIM+ and BMIM+, the conformers of Cs symmetry with all-trans alkyl chains (Figure 1), although not the most stables ones, have been chosen because they have large polarizability anisotropy. The root-mean-square deviation (RMSD = [(1/ 2 1/2 Np)∑i =Np 1(Umodel − UQM ind,i ind,i) ] ) and the percentage rootmean-square deviation (RMSD % = 100[(1/Np)∑iN=p 1(Umodel ind,i − 2 QM 2 1/2 UQM ind,i) /(Uind,i) ] ) between the reference induction energies at all Np points of the grids and those obtained from the derived distributed polarizabilities are used to test the quality of the proposed models. To test the transferability of the DP model in the alkylimidazolium series, the molecular polarizability has been calculated by QM for the optimized elongated conformation of 1-pentyl-3-methylimidazolium (PMIM+) and 1-hexyl-3-methylimidazolium (HMIM+) cations and for two additional conformers of EMIM+, seven new conformers of BMIM+ and seven of PMIM+ and HMIM+, obtained by turning the dihedral angles of the alkyl chains. 2.3. Fitting Protocols of the Selected Distributed Polarizability Models. As observed previously, fitting of CFDiso and Daniso models to reproduce induction energies sometimes leads to counterintuitive parameters.64,65 For example, atomic dipolar polarizability components αAdip and charge-flow polarizabilities αAA 00 , which are naturally positive can be found negative for interior atoms, leading to additional problems of transferability. To circumvent this feature and reduce the number of independent parameters, we decided to impose some proportionality constraints between selected atomic dipolar polarizabities, using as much as possible ratio’s from the literature,80 as specified below. For the anions, we imposed α̅ Bdip/α̅Fdip = 5.4, α̅Pdip/α̅ Fdip = 6.5, N α̅dip/α̅ Odip = 1.4, α̅ Sdip/α̅Cdip = 1.65, and α̅Ndip/α̅ Cdip = 0.62. Thus, with the CFDiso model, one charge-flow parameter α00 AA is fitted for the central atom of BF−4 , PF−6 , and NO−3 , and for the S and C atoms of SCN−, in addition to one isotropic dipolar polarizability for the F (BF−4 ), F(PF−6 ), O(NO−3 ), and C(SCN−) atoms. With the Daniso model, two parameters are fitted for BF−4 , namely the average dipolar polarizability of the F2 atom, α̅F2 dip

SCN−

αBB 00 αBdip αF,∥ dip αF,⊥ dip αPP 00 αPdip αF,∥ dip αF,⊥ dip αNN 00 N αdip,x′x′ N αdip,y′y′ αNdip,z′z′ αOdip,z′z′ O αdip,x′x′ O αdip,y′y′ αSS 00 αCC 00 αS,∥ dip αS,⊥ dip αC,∥ dip αC,⊥ dip αN,∥ dip αN,⊥ dip

2.95 8.10 1.50 1.50 5.5 6.5 1.0 1.0 6.80 5.53 5.53 5.53 3.95 3.95 3.95 1.65 6.10 15.923 15.923 9.65 9.65 5.983 5.983

Daniso 12.42 3.90 1.50 15.058 4.15 1.40 11.34 11.34 5.67 11.15 5.05 4.05

37.785 15.675 22.90 9.50 14.198 5.89

a (∥): parallel to the bond and (⊥) perpendicular to the bond; (x′, y′, z′): local site frame (see Figure 1).

three polarizability components of O atom in its local frame. The N parameters are obtained taking αNdip,x′x′ = αNdip,y′y′ and using a scaling factor of 1.4 with respect to the O atom for inplane and out-of-plane components. For the C atom of SCN−, two polarizability parameters are fitted: αCdip,z′z′ and αCdip,x′x′ = αCdip,y′y′, respectively, parallel and perpendicular to the symmetry axis z′ = z. The corresponding parameters of S and N atoms are obtained using the previously mentioned scaling factor. For the imidazolium cations, dipolar polarizabilities on N and C atoms of different hybridization states have been 2 C,sp2 C,sp3 C,sp2 distinguished by imposing α̅N,sp dip /α̅ dip = 0.62 and α̅ dip /α̅ dip 35 = 0.82. For the CFDiso model, the charge-flow polarizabilities αAA 00 (A = C, N, and H) have been fitted, assuming the n H atoms of a CHn group to be equivalent and that the CH2 groups at C8 and C9 positions of BMIM+ are also equivalent (see Figure 1). The resulting CFDiso models contain 8, 15, and 18 fitted parameters for the DMIM+, EMIM+, and BMIM+ ions, respectively. Such an approach gave reasonable parameters for each ion, but these turned out to be not sufficiently transferable from a cation to the other. Thus, a new set of parameters have been fitted by constraining some of them to keep the same values for the DMIM+, EMIM+, and BMIM+ ions. First, the charge-flow polarizability components αAA 00 of the H atoms at C1, C3, and C4 positions and at each CH3 group have been fixed to 0.7, 0.75, 0.75, and 0.95, respectively. These values have been chosen from the results of test fits. For the Daniso model, some constraints to the atomic or group components located 8845

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

on the N and C atoms have been applied. For the CH3 carbons, the local C3v symmetry imposes αCdip,x′x′ = αCdip,y′y′ and we also imposed αCdip,z′z′ = 1.05αCdip,x′x′ (this relation stems from test calculations on DMIM+). For CH2 carbons, we assumed αCdip,x′x′ C = 1.65αdip,z′z′ , also inspired by the observed change of polarizability by adding CH2 groups in the studied cation series. For the aromatic atoms, we used αNdip,z′z′ = αNdip,x′x′ for the N atoms and αCdip,z′z′ = 0.5αCdip,x′x′ for the CH carbons. The latter relationship is consistent with a previous study on benzene.64 Test calculations with the above-mentioned constraints led, however, to a negative value of αC1 dip,y′y′. We therefore further C1 assumed that αC1 dip,y′y′ = 0.3αdip,x′x′. These two components were in fact found to be proportional by a factor close to 0.3 in a series of fits on BMIM+, where C1, C3, and C4 were taken equivalent. With this set of constraints, the Daniso model finally required the fitting of 4 independent parameters for DMIM+ and 6 parameters for the EMIM+ and BMIM+ ions. For both CFDiso and Daniso models, the parameters have been averaged over the three cations.

Table 3. Fitted Distributed Polarizabilities (in a.u.) for the CFDiso and Daniso Models of 1-Alkyl,3-alkylimidazolium cationsa

3. RESULTS AND DISCUSSION The distributed polarizability parameters obtained with the CFDiso and Daniso models are given in Table 1 for anions, while Table 2 reports comparisons with the reference QM results

a

CFDiso αC1C1 00 αN2N2 00 αC3C3 00 αC4C4 00 αN5N5 00 ) αC6C6(CH 00 3 C7C7(CH) α00 2 αHH(C1) 00 αHH(C3,C4) 00 ) αHH(CH 00 3 HH(CH) α00 2

CFDiso

PF−6

NO−3

SCN−

α̅ RMSD RMSD% α̅ RMSD RMSD% α∥ α⊥ α̅ RMSD RMSD% α∥ α⊥ α̅ RMSD RMSD%

21.04 0.031 5.0 29.67 0.015 4.2 36.62 17.38 30.21 0.042 3.4 70.22 31.56 44.44 0.031 9.0

Daniso 21.62 0.017 2.3 28.96 0.016 4.9 35.64 17.82 29.70 0.047 4.1 74.88 31.07 45.67 0.031 8.4

CFDiso

Daniso

5.5 5.5 5.5 5.5 5.5 5.5 3.41 3.41 3.41 4.51 4.51 4.51 4.51 4.51 4.51

12.4 3.72 6.2 15.9 5.0 7.95 7.394 3.1 7.394 13.3 13.3 13.965 16.1 8.8 9.758

(x′, y′, z′): local site frame (see Figure 1).

Table 4. Molecular Dipolar Polarizabilities (in a.u.) of 1Alkyl-3-Alkylimidazolium DMIM+, EMIM+, and BMIM+ Cations Rebuilt from CFDiso and Daniso Models and Obtained Directly by QM MP2/6-311+G(2df) Calculationsa

Table 2. Molecular Dipolar Polarizabilities (in a.u.) of BF−4 , PF−6 , NO−3 , and SCN− Anions Rebuilt from CFDiso and Daniso Models and Obtained Directly by QM MP2/6-311+G(2df) Calculationsa BF−4

3.1 3.8 2.9 2.9 3.8 4.7 5.7 0.7 0.75 0.95 0.85

αC1 dip,x′x′ αC1 dip,y′y′ αC1 dip,z′z′ αC3,C4 dip,x′x′ αC3,C4 dip,y′y′ αC3,C4 dip,z′z′ αN2,N5 dip,x′x′ αN2,N5 dip,y′y′ αN2,N5 dip,z′z′ 3) αC(CH dip,x′x′ C(CH3) αdip,y′y′ 3) αC(CH dip,z′z′ 2) αC(CH dip,x′x′ C(CH2) αdip,y′y′ 2) αC(CH dip,z′z′

DMIM+

QM 20.75

29.03

EMIM+

35.85 17.32 29.67 BMIM+

73.96 31.98 45.97

a α̅ : isotropic polarizability; α∥: parallel to the bond or to the molecular plane; α⊥: perpendicular to the bond or to the molecular plane (see Figure 1). RMSD (in 10−3 a.u.) and RMSD % are defined in the text.

αx,x αy,y αz,z αx,z α̅ RMSD RMSD % αx,x αy,y αz,z αx,z α̅ RMSD RMSD % αx,x αy,y αz,z αx,z α̅ RMSD RMSD %

CFDiso

Daniso

QM

84.14 43.12 65.94 0.000 64.40 0.020 2.7 100.20 52.35 78.22 −1.92 76.92 0.020 2.4 132.31 70.67 102.59 −2.22 101.86 0.027 3.2

79.99 46.52 70.41 0.000 65.64 0.024 3.6 95.753 55.320 80.509 −1.56 77.19 0.025 3.0 127.49 72.92 100.48 −3.95 100.30 0.032 3.3

81.73 43.74 66.82 0.000 64.10

97.32 52.78 76.49 −1.21 75.53

131.44 71.28 96.51 −3.06 99.74

α̅ : isotropic polarizability; (x, y, z): global laboratory frame (see Figure 1). RMSD (in 10−3 a.u.) and RMSD % are defined in the text. a

regarding the calculation of the induction energies and of the rebuilt molecular dipolar polarizabilities from the DP models. Results for cations are reported similarly in Tables 3 and 4. For the CFDiso model, Tables 1 and 3 contain only the fitted charge-flow polarizabilities. Those resulting from the charge conservation constraint are obtained from eq 9. The full sets of charge-flow polarizabilities of the DMIM+, EMIM+, and BMIM+ cations are listed in the Supporting Information. With regard to the comparison of the induction energies of the polarized ions obtained from the DP models and from the QM calculations, it can be observed that both models perform similarly and are quite good in this respect, especially for the

cations for which the RMSD amounts to about 3%. For the anions, the RMSD is around 3−4% for BF−4 , PF−6 , and NO−3 and 8−9% for SCN−. Another quality test of the models concerns the rebuilding of the components of the molecular polarizability tensor αa,b (a,b = x, y, or z). Globally, both models rebuild a molecular average dipolar polarizability α̅ in good agreement with the target ab initio value. The difference is typically an increase of about 1− 3% with respect to the ab initio value, except for the SCN− anion, for which a decrease by 3% is obtained for the CFDiso model. These models reproduce also quite well the polar8846

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

Figure 2. Diagonal molecular polarizability components (αx′,x′, αy′,y′, and αz′,z′) and anisotropy factor κaniso for 26 conformers of DMIM+, EMIM+, BMIM+, PMIM+, and HMIM+ ions: comparison of the CFDiso to the ab initio values.

Figure 3. Diagonal molecular polarizability components(αx′,x′, αy′,y′, αz′,z′) and anisotropy factor κaniso for 26 conformers of DMIM+, EMIM+, BMIM+, PMIM+, and HMIM+ ions: comparison of the Daniso to the ab initio values.

EMIM+, BMIM+, PMIM+, and HMIM+. The CFDiso model (Figure 2) reproduces very well αx′,x′ but slightly overestimates the αy′,y′ and αz′,z′ components; the κaniso factor is thus slightly underestimated on the average. The Daniso model (Figure 3) leads to less scattered polarizability components with respect to

izability anisotropy. Figures 2 and 3 present the three components (αx′x′, αy′y′, and αz′z′) of the diagonalized molecular polarizability tensor and the κaniso anisotropy factor for CFDiso and Daniso models correlated to the reference ab initio QM values for the three cations and the additional 23 conformers of 8847

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

components of RTILs have been fitted from ab initio induction energy calculations. These models describe the polarization effects by either charge flow between neighbor atoms plus isotropic dipolar polarizabilities of non-hydrogen atoms (CFDiso model) or by anisotropic dipolar polarizabilities of nonhydrogen atoms (Daniso model). For the cations, particular attention has been given to obtain parameters transferable in the 1-alkyl-3-alkylimidazolium series. Both models account for the anisotropy of the molecular dipolar polarizability in a natural way. Both also quite well reproduce induction energies of polarized RTIL ions as well as the main characteristics of the polarizability tensor. With regard to the computational demands, both CFDiso and Daniso models are similar, as seen by MD tests performed on the [DMIM][BF4] ionic liquid (343 ion pairs, assuming rigid ions): adding the charge flow polarizabilities (i.e., from the Daniso to the CFDiso model) increases the cost by only 6%. Further comparisons between the two models and parameter development should emerge from molecular dynamics simulations of ionic liquids based on different ions combinations.

the ab initio values. The trend is to slightly underestimate αx′,x′ and to slightly overestimate αy′,y′ and αz′,z′, thus the κaniso anisotropy factor is systematically somewhat underestimated. It is interesting to compare our models for cations to Thole’s models fitted by Gu and Yan to reproduce molecular polarizabilities obtained from QM calculations at the MP2/ aug-cc-pVDZ level.26 For our studied ions (DMIM+, EMIM+, BMIM+, and HMIM+), the average ionic polarizabilities α̅ are lower than theirs by 2−4% with a maximum difference of 5.3% for HMIM+ with the Daniso model. Nitrate ions have been simulated by molecular dynamics with polarizable models.9,25,81−87 Salvador et al. studied an aqueous nitrate ion using an anisotropic polarizability model with an α̅ of 26.3 a.u.82 Yan et al. used the Thole’s model for nitrate in ionic liquids, with α̅ equal to 23,9 32.8,84 or 30.4 a.u.25 Note that the latter value is quite close to ours (29.8 for CFDiso and 30.3 a.u. for Daniso). The BF−4 anion has also been simulated in RTILs with polarization using α̅ values of 24.088 and 18.4 a.u.,89 bracketing our values (21.2 and 21.6 a.u., respectively, for the CFDiso and Daniso models). For the SCN− anion, Petersen et al. reported gas phase polarizabilites of αz,z = 78.3 a.u. and αx,x = αy,y = 39.1 a.u. and used isotropic atomic polarizabilities of 13.5 a.u. on S, C, and N atoms in a MD simulation of sodium thiocyanate in liquid water.90 We find an anisotropy factor κaniso in the range of 0.29 (CFDiso)−0.32 (Daniso), a little bit larger than the value determined by Petersen et al. in the gas phase (0.25). Our average polarizability α̅ for SCN− around 45−46 a.u. is intermediate between the values used by Petersen et al. for gas phase (52.2 a.u.) and in MD simulations of aqueous solutions (40.5 a.u.). Regarding alkylimidazolium-based RTILs, Bica et al. reported ionic polarizabilities obtained from refractive index and density measurements.91 These authors also proposed a set of atomic polarizabilities, reproducing experimental ionic polarizabilities. Their resulting α̅ for the anions are 18.9 (BF−4 ), 28.2 (PF−6 ), 20.6 (NO−3 ), and 32.8 (SCN−) a.u., and those for the cations are 74.1 (DMIM+), 86.1 (EMIM+), 110.3 (BMIM+), 122.3 (PMIM+) and 134.4 (HMIM+) a.u. When compared to these values, ours obtained in the gas phase (Tables 1 and 2 and in the Supporting Information) are larger by a factor of 1.05−1.45 for the anions and lower by a factor of 0.85−0.9 for the cations. Note that in condensed phases, the influence of the environment on ion polarizability can be significant, as suggested by ab initio simulations.92−95 For instance, in water, a reduction by 30−50% of the polarizability has been estimated for monatomic anions, whereas for monatomic cations the reduction is much smaller (at most 15% in the case of Cs+).94 Another study reported, however, a smaller reduction of α̅ for hydrated chloride ion with respect to the gas phase value.96 What happens with RTILs in the condensed phase deserves further study.

■

ASSOCIATED CONTENT

S Supporting Information *

Atomic Cartesian coordinates, MP2/6-311+G(2df) dipolar polarizabilies, CFDiso and Daniso ionic dipolar polarizabilities of the four anions and of the 26 conformers of imidazolium cations, a full list of nonzero charge-flow polarizabilities of DMIM+, EMIM+, and BMIM+ ions are listed, and a full ref 79 is reported. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Tel: ++33 (0)3 83 68 43 84. Fax: ++33 (0)3 83 68 43 71. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We thank the University of Lorraine, the University of Strasbourg, and the CNRS for their support. C.M. thanks Marius Retegan for providing some polarizability calculations of imidazolium cations and Jean-Christophe Soetens for performing a MD simulation test on liquid [DMIM][BF4].

■

REFERENCES

(1) Maginn, E. J. Atomistic Simulation of the Thermodynamic and Transport Properties of Ionic Liquids. Acc. Chem. Res. 2007, 40, 1200− 1207. (2) Hunt, P. A.; Maginn, E. J.; Lynden-Bell, R. M.; Del Pópolo, M. G. In Ionic Liquids in Synthesis; Wasserscheid, P., Welton, T., Eds.; WileyVCH Verlags GmBH & Co KGaA: Weinheim, 2008; pp 206−249. (3) Dommert, F.; Wendler, K.; Berger, R.; Delle Site, L.; Holm, C. Force Fields for Studying the Structure and Dynamics of Ionic Liquids: A Critical Review of Recent Developments. ChemPhysChem 2012, 13, 1625−1637. (4) Hanke, C. G.; Price, S. L.; Lynden-Bell, R. M. Intermolecular Potentials for Simulations of Liquid Imidazolium Salts. Mol. Phys. 2001, 99, 801−809. (5) de Andrade, J.; Böes, E. S.; Stassen, H. Computational Study of Room Temperature Molten Salts Composed of 1-Alkyl-3-methylimidazolium Cations: Force Field Proposal and Validation. J. Phys. Chem. B 2002, 106, 13344−13351.

4. CONCLUSION Accounting for polarization in force fields used to simulate condensed phases or biological systems by classical MD is still a developing field.97−99 The use of relatively simple models is required for an efficient implementation of polarizable models in simulations, as implemented so far in major force field groups like AMBER, CHARMM, OPLS, GROMOS, and SIBFA to cite a few. An excellent review can be found in ref 99. In the present work, two models, CFDiso and Daniso, of distributed polarizabilities of ions (BF−4 , PF−6 , NO−3 , SCN− anions and 1-alkyl-3-alkylimidazolium cations) that are typical 8848

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

Point Induced Dipoles at Various Densities. J. Chem. Phys. 2013, 138, 204119. (28) Schröder, C.; Steinhauser, O. Simulating Polarizable Molecular Ionic Liquids with Drude Oscillators. J. Chem. Phys. 2010, 133, 154511. (29) Madden, P. A.; Wilson, M. ‘Covalent Effects’ in Ionic’ Systems. Chem. Soc. Rev. 1996, 25, 339−350. (30) Salanne, M.; Madden, P. A. Polarization Effects in Ionic Solids and Melts. Mol. Phys. 2011, 109, 2299−2315. (31) Wilson, M.; Madden, P. A.; Costa Cabral, B. J. Quadrupole Polarization in Simulations of Ionic Systems: Application to AgCl. J. Phys. Chem. 1996, 100, 1227−1237. (32) Le Fèvre, R. J. W. Molecular Refractivity and Polarizability. Adv. Phys. Org. Chem. 1965, 3, 1−90. (33) Applequist, J.; Carl, J. R.; Fung, K.-K. An Atom Dipole Interaction Model for Molecular Polarizability. Application to Polyatomic Molecules and Determination of Atomic Polarizabilities. J. Am. Chem. Soc. 1972, 94, 2952−2960. (34) Applequist, J. An Atom Dipole Interaction Model for Molecular Optical Properties. Acc. Chem. Res. 1977, 10, 79−85. (35) Miller, K. J.; Savchik, J. A. A New Empirical Method to Calculate Average Molecular Polarizabilities. J. Am. Chem. Soc. 1979, 101, 7206−7213. (36) Birge, R. R. Calculation of Molecular Polarizabilities Using an Anisotropic Atom Point Dipole Interaction Model which Includes the Effect of Electron Repulsion. J. Chem. Phys. 1980, 72, 5312−5319. (37) Thole, B. T. Molecular Polarizabilities Calculated with a Modified Dipole Interaction. Chem. Phys. 1981, 59, 341−350. (38) Miller, K. J. Additivity Methods in Molecular Polarizability. J. Am. Chem. Soc. 1990, 112, 8533−8542. (39) No, K. T.; Cho, K. H.; Jhon, M. S.; Scheraga, H. A. An Empirical Method to Calculate Average Molecular Polarizabilities from the Dependence of Effective Atomic Polarizabilities on Net Atomic Charge. J. Am. Chem. Soc. 1993, 115, 2005−2014. (40) Voisin, C.; Cartier, A. Determination of Distributed Polarizabilities to Be Used for Peptide Modeling. THEOCHEM 1993, 286, 35−45. (41) Applequist, J. Atom Charge Transfer in Molecular Polarizabilities. Application of the Olson-Sundberg Model to Aliphatic and Aromatic Hydrocarbons. J. Phys. Chem. 1993, 97, 6016−6023. (42) Stout, J. M.; Dykstra, C. E. Static Dipole Polarizabilities of Organic Molecules. Ab Initio Calculations and a Predictive Model. J. Am. Chem. Soc. 1995, 117, 5127−5132. (43) Bode, K. A.; Applequist, J. A New Optimization of Atom Polarizabilities in Halomethanes, Aldehydes, Ketones and Amides by Way of the Atom Dipole Interaction Model. J. Phys. Chem. 1996, 100, 17820−17824. (44) van Duijnen, P. T.; Swart, M. Molecular and Atomic Polarizabilities: Thole’s Model Revisited. J. Phys. Chem. A 1998, 102, 2399−2407. (45) Rappé, A. K.; Goddard, W. A., III Charge Equilibration for Molecular Dynamics Simulations. J. Phys. Chem. 1991, 95, 3358−3363. (46) Karlström, G. Local Polarizabilities in Molecules, Based on Ab Initio Hartree-Fock Calculations. Theor. Chim. Acta 1982, 60, 535− 541. (47) Garmer, D. R.; Stevens, W. J. Transferability of Molecular Distributed Polarizabilities from a Simple Localized Orbital Based Method. J. Phys. Chem. 1989, 93, 8263−8270. (48) Åstrand, P.-O.; Karlström, G. Local Polarizability Calculations with Localized Orbitals in the Uncoupled Hartree Fock Approximation. Mol. Phys. 1992, 77, 143−155. (49) Giribet, C. G.; Demarco, M. D.; Ruiz de Azúa, M. C.; Contreras, R. H. Ab Initio ‘CLOPPA’ Decomposition of the Static Molecular Polarizability Tensor. Mol. Phys. 1997, 91, 105−111. (50) Stone, A. J. Distributed Polarizabilities. Mol. Phys. 1985, 56, 1065−1082. (51) Le Sueur, C. R.; Stone, A. J. Practical Schemes for Distributed Polarizabilities. Mol. Phys. 1993, 78, 1267−1291.

(6) Liu, Z.; Huang, S.; Wang, W. A Refined Force Field for Molecular Simulation of Imidazolium-Based Ionic Liquids. J. Phys. Chem. B 2004, 108, 12978−12989. (7) Liu, Z.; Wu, X.; Wang, W. A Novel United-Atom Force Field for Imidazolium-Based Ionic Liquids. Phys. Chem. Chem. Phys. 2006, 8, 1096−1104. (8) Canongia Lopes, J. N.; Pádua, 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. (9) Yan, T.; Burnham, C. J.; Del Pópolo, M. G.; Voth, G. A. Molecular Dynamics Simulation of Ionic Liquids: The Effect of Electronic Polarizability. J. Phys. Chem. B 2004, 108, 11877−11881. (10) Lynden-Bell, R. M.; Youngs, T. G. A. Simulations of Imidazolium Ionic Liquids: When Does the Cation Charge Distribution Matter ? J. Phys.: Condens. Matter 2009, 21, 424120. (11) Yan, T.; Wang, Y.; Knox, C. On the Dynamics of Ionic Liquids: Comparison Between Electronically Polarizable and Nonpolarizable Models II. J. Phys. Chem. B 2010, 114, 6886−6904. (12) Chaumont, A.; Schurhammer, R.; Wipff, G. Aqueous Interfaces with Hydrophobic Room-Temperature Ionic Liquids: A Molecular Dynamics Study. J. Phys. Chem. B 2005, 109, 18964−18973. (13) Bühl, M.; Chaumont, A.; Schurhammer, R.; Wipff, G. Ab initio Molecular Dynamics of Liquid 1,3-Dimethylimidazolium Chloride. J. Phys. Chem. B 2005, 109, 18591−18599. (14) Chevrot, G.; Schurhammer, R.; Wipff, G. Molecular Dynamics Simulations of the Aqueous Interface with the [BMI][PF6] Ionic Liquid: Comparison of Different Solvent Models. Phys. Chem. Chem. Phys. 2006, 8, 4166−4174. (15) Bhargava, B. L.; Balasubramanian, S. Refined Potential Model for Atomistic Simulations of Ionic Liquid [BMIM][PF6]. J. Chem. Phys. 2007, 127, 114510. (16) Zhao, W.; Eslami, H.; Cavalcanti, W. L.; Müller-Plathe, F. A Refined All-Atom Model for the Ionic Liquid 1-n-Butyl-3-methylimidazolium Bis(trifluoromethylsulfonyl)imide [bmim][Tf2N]. Z. Phys. Chem. 2007, 221, 1647−1662. (17) Youngs, T. G. A.; Hardacre, C. Application of Static ChargeTransfer within an Ionic-Liquid Force Field and Its Effect on Structure and Dynamics. ChemPhysChem 2008, 9, 1548−1558. (18) Chaban, V. V.; Voroshylova, I. V.; Kalugin, O. N. A New Force Field Model for the Simulation of Transport Properties of Imidazolium-Based Ionic Liquids. Phys. Chem. Chem. Phys. 2011, 13, 7910−7920. (19) Chaban, V. Polarizability versus Mobility: Atomistic Force Field for Ionic Liquids. Phys. Chem. Chem. Phys. 2011, 13, 16055−16062. (20) Zhang, Y.; Maginn, E. J. A Simple AIMD Approach to Derive Atomic Charges for Condensed Phase Simulation of Ionic Liquids. J. Phys. Chem. B 2012, 116, 10036−10048. (21) Schröder, C. Comparing Reduced Partial Charge Models with Polarizable Simulations of Ionic Liquids. Phys. Chem. Chem. Phys. 2012, 14, 3089−3102. (22) Chang, T.-M.; Dang, L. X. Computational Studies of Structures and Dynamics of 1,3-Dimethylimidazolium Salt Liquids and Their Interfaces Using Polarizable Potential Models. J. Phys. Chem. A 2009, 113, 2127−2135. (23) Borodin, O. Polarizable Force Field Development and Molecular Dynamics Simulations of Ionic Liquids. J. Phys. Chem. B 2009, 113, 11463−11478. (24) Bedrov, D.; Borodin, O.; Li, Z.; Smith, G. D. Influence of Polarization on Structural, Thermodynamic, and Dynamic Properties of Ionic Liquids Obtained from Molecular Dynamics Simulations. J. Phys. Chem. B 2010, 114, 4984−4997. (25) Yan, T.; Wang, Y.; Knox, C. On the Structure of Ionic Liquids: Comparison Between Electronically Polarizable and Nonpolarizable Models I. J. Phys. Chem. B 2010, 114, 6905−6921. (26) Gu, Y.; Yan, T. Thole Model for Ionic Liquid Polarizability. J. Phys. Chem. A 2013, 117, 219−227. (27) Taylor, T.; Schmollngruber, M.; Schröder, C.; Steinhauser, O. The Effect of Thole Functions on the Simulation of Ionic Liquids with 8849

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

(52) Laidig, K. E.; Bader, R. F. W. Properties of Atoms in Molecules: Atomic Polarizability. J. Chem. Phys. 1990, 93, 7213−7224. (53) Bader, R. F. W.; Keith, T. A.; Gough, K. M.; Laidig, K. E. Properties of Atoms in Molecules: Additivity and Transferability of Group Polarizabilities. Mol. Phys. 1992, 75, 1167−1189. (54) Á ngyán, J. G.; Jansen, G.; Loos, M.; Hättig, C.; Heß, B. A. Distributed Polarizabilities Using the Topological Theory of Atoms in Molecules. Chem. Phys. Lett. 1994, 219, 267−273. (55) Misquitta, A. J.; Stone, A. J. Distributed Polarizabilities Obtained Using a Constrained Density-Fitting Algorithm. J. Chem. Phys. 2006, 124, 024111. (56) Rob, F.; Szalewicz, K. Asymptotic Dispersion Energies from Distributed Polarizabilities. Chem. Phys. Lett. 2013, 572, 146−149. (57) Rob, F.; Szalewicz, K. Distributed Molecular Polarizabilities and Asymptotic Intermolecular Interaction Energies. Mol. Phys. 2013, 111, 1430−1455. (58) Hirshfeld, F. L. Bonded-Atom Fragments for Describing Molecular Charge Densities. Theor. Chim. Acta 1977, 44, 129−138. (59) Yang, M.; Senet, P.; Van Alsenoy, C. DFT Study of Polarizabilities and Dipole Moments of Water Clusters. Int. J. Quantum Chem. 2005, 101, 535−542. (60) Krishtal, A.; Senet, P.; Yang, M.; Van Alsenoy, C. A Hirshfeld Partitioning of Polarizabilities of Water Clusters. J. Chem. Phys. 2006, 125, 034312. (61) Marenich, A. V.; Cramer, C. J.; Truhlar, D. G. Reduced and Quenched Polarizabilities of Interior Atoms in Molecules. Chem. Sci. 2013, 4, 2349−2356. (62) Nakagawa, S.; Kosugi, N. Polarized One-Electron Potentials Fitted by Multicenter Polarizabilities and Hyperpolarizabilities. Ab initio SCF-CI Calculation of water. Chem. Phys. Lett. 1993, 210, 180− 186. (63) Alkorta, I.; Bachs, M.; Perez, J. J. The Induced Polarization of the Water Molecule. Chem. Phys. Lett. 1994, 224, 160−165. (64) Celebi, N.; Á ngyán, J. G.; Dehez, F.; Millot, C.; Chipot, C. Distributed Polarizabilities Derived from Induction Energies: A Finite Perturbation Approach. J. Chem. Phys. 2000, 112, 2709−2717. (65) Dehez, F.; Soetens, J.-C.; Chipot, C.; Á ngyán, J. G.; Millot, C. Determination of Distributed Polarizabilities from a Statistical Analysis of Induction Energies. J. Phys. Chem. A 2000, 104, 1293−1303. (66) Williams, G. J.; Stone, A. J. Distributed Dispersion: A New Approach. J. Chem. Phys. 2003, 119, 4620−4628. (67) Lillestolen, T. C.; Wheatley, R. J. First-Principles Calculation of Local Atomic Polarizabilities. J. Phys. Chem. A 2007, 111, 11141− 11146. (68) Banks, J. L.; Kaminski, G. A.; Zhou, R.; Mainz, D. T.; Berne, B. J.; Friesner, R. A. Parametrizing a Polarizable Force Field from Ab Initio Data. I. The Fluctuating Point Charge Model. J. Chem. Phys. 1999, 110, 741−754. (69) Stern, H. A.; Kaminski, G. A.; Banks, J. L.; Zhou, R.; Berne, B. J.; Friesner, R. A. Fluctuating Charge, Polarizable Dipole, and Combined Models: Parametrization from Ab Initio Quantum Chemistry. J. Phys. Chem. B 1999, 103, 4730−4737. (70) Stern, H. A.; Rittner, F.; Berne, B. J.; Friesner, R. A. Combined Fluctuating Charge and Polarizable Dipole Models: Application to a Five-Site Water Potential Function. J. Chem. Phys. 2001, 115, 2237− 2251. (71) Anisimov, V. M.; Lamoureux, G.; Vorobyov, I. V.; Huang, N.; Roux, B.; MacKerell, A. D., Jr. Determination of Electrostatic Parameters for a Polarizable Force Field Based on the Classical Drude Oscillator. J. Chem. Theory Comput. 2005, 1, 153−168. (72) Stone, A. J.; Tong, C.-S. Local and Non-Local Dispersion Models. Chem. Phys. 1989, 137, 121−135. (73) Pauling, L. The Nature of the Chemical Bond, 3rd ed.; Cornell University Press: Ithaca, NY, 1960; p 260. (74) Stone, A. J. In Theoretical Models of Chemical Bonding; Maksić, H., Ed.; Springer−Verlag: Berlin, 1991; Vol. 4; pp 103−31. (75) Chipot, C.; Á ngyán, J. G.; Millot, C. Statistical Analysis of Distributed Multipoles Derived from the Molecular Electrostatic Potential. Mol. Phys. 1998, 94, 881−895.

(76) Soetens, J.-C.; Millot, C. Effect of Distributing Multipoles and Polarizabilities on Molecular Dynamics Simulations of Water. Chem. Phys. Lett. 1995, 235, 22−30. (77) Ma, B.; Lii, J.-H.; Allinger, N. L. Molecular Polarizabilities and Induced Dipole Moments in Molecular Mechanics. J. Comput. Chem. 2000, 21, 813−825. (78) Gresh, N.; Cisneros, G. A.; Darden, T. A.; Piquemal, J.-P. Anisotropic, Polarizable Molecular Mechanics Studies of Inter- and Intramolecular Interactions and Ligand-Macromolecule Complexes. A Bottom-Up Strategy. J. Chem. Theory Comput. 2007, 3, 1960−1986. (79) 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.; et al. GAUSSIAN 03, revision C.02; Gaussian Inc.: Wallingford, CT, 2004. (80) Nagle, J. K. Atomic Polarizability and Electronegativity. J. Am. Chem. Soc. 1990, 112, 4741−4747. (81) Ribeiro, M. C. C.; Almeida, L. C. J. Validating a Polarizable Model for the Glass-Forming Liquid Ca0.4K0.6(NO3)1.4 by Ab Initio Calculations. J. Chem. Phys. 2000, 113, 4722−4731. (82) Salvador, P.; Curtis, J. E.; Tobias, D. J.; Jungwirth, P. Polarizability of the Nitrate Anion and Its Solvation at the Air/ Water Interface. Phys. Chem. Chem. Phys. 2003, 5, 3752−3757. (83) Urahata, S. M.; Ribeiro, M. C. C. Molecular Dynamics Simulation of Molten LiNO3 with Flexible and Polarizable Anions. Phys. Chem. Chem. Phys. 2003, 5, 2619−2624. (84) Yan, T.; Li, S.; Jiang, W.; Gao, X.; Xiang, B.; Voth, G. A. Structure of the Liquid-Vacuum Interface of Room-Temperature Ionic Liquids: A Molecular Dynamics Study. J. Phys. Chem. B 2006, 110, 1800−1806. (85) Thomas, J. L.; Roeslová, M.; Dang, L. X.; Tobias, D. J. Molecular Dynamics Simulations of the Solution-Air Interface of Aqueous Sodium Nitrate. J. Phys. Chem. A 2007, 111, 3091−3098. (86) Borodin, O.; Smith, G. D.; Sewell, T. D.; Bedrov, D. Polarizable and Nonpolarizable Force Fields for Alkyl Nitrates. J. Phys. Chem. B 2008, 112, 734−742. (87) Ribeiro, M. C. C. Polarization Effects in Molecular Dynamics Simulations of Glass-Formers Ca(NO3)2.nH2O, n=4,6 and 8. J. Chem. Phys. 2010, 132, 134512. (88) Borodin, O.; Smith, G. D.; Douglas, R. Force Field Development and MD Simulations of Poly(ethylene oxide)/LiBF4 Polymer Electrolytes. J. Phys. Chem. B 2003, 107, 6824−6837. (89) Wick, C. D.; Chang, T.-M.; Dang, L. X. Molecular Mechanism of CO2 and SO2 Molecules Binding to the Air/Liquid Interface of 1Butyl-3-Methylimidazolium Tetrafluoroborate Ionic Liquid: A Molecular Dynamics Study with Polarizable Potential Models. J. Phys. Chem. B 2010, 114, 14965−14971. (90) Petersen, P. B.; Saykally, R. J.; Mucha, M.; Jungwirth, P. Enhanced Concentration of Polarizable Anions at the Liquid Water Surface: SHG Spectroscopy and MD Simulations of Sodium Thiocyanide. J. Phys. Chem. B 2005, 109, 10915−10921. (91) Bica, K.; Deetlefs, M.; Schröder, C.; Seddon, K. R. Polarizabilities of Alkylimidazolium Ionic Liquids. Phys. Chem. Chem. Phys. 2013, 15, 2703−2711. (92) Jë mmer, P.; Fowler, P. W.; Wilson, M.; Madden, P. A. Environmental Effects on Anion Polarizability: Variation with Lattice Parameter and Coordination Number. J. Phys. Chem. A 1998, 102, 8377−8385. (93) Salanne, M.; Vuilleumier, R.; Madden, P. A.; Simon, C.; Turq, P.; Guillot, B. Polarizabilities of Individual Molecules and Ions in Liquids from First Principles. J. Phys.: Condens. Matter 2008, 20, 494207. (94) Molina, J. J.; Lectez, S.; Tazi, S.; Salanne, M.; Dufrêche, J.-F.; Roques, J.; Simoni, E.; Madden, P. A.; Turq, P. Ions in Solutions: Determining Their Polarizabilities from First-Priciples. J. Chem. Phys. 2011, 134, 014511. (95) Bauer, B. A.; Lucas, T. R.; Krishtal, A.; Van Alsenoy, C.; Patel, S. Variation of Ion Polarizability from Vacuum to Hydration: Insights from Hirshfeld Partitioning. J. Phys. Chem. A 2010, 114, 8984−8992. 8850

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851

The Journal of Physical Chemistry A

Article

(96) Masia, M. Estimating Chloride Polarizability in a Water Solution. J. Phys. Chem. A 2013, 117, 3221−3226. (97) Yu, H.; van Gunsteren, W. F. Accounting for Polarization in Molecular Simulation. Comput. Phys. Commun. 2005, 172, 69−85. (98) Warshel, A.; Kato, M.; Pisliakov, A. V. Polarizable Force Fields: History, Test Cases and Prospects. J. Chem. Theory Comput. 2007, 3, 2034−2045. (99) Cieplak, P.; Dupradeau, F.-Y.; Duan, Y.; Wang, J. Polarization Effects in Molecular Mechanical Force Fields. J. Phys.: Condens. Matter 2009, 21, 333102.

8851

dx.doi.org/10.1021/jp505539y | J. Phys. Chem. A 2014, 118, 8842−8851