Insights Gained from Refined Force-Field for Pure and Aqueous

Article Cite This: J. Phys. Chem. B XXXX, XXX, XXX−XXX

pubs.acs.org/JPCB

Insights Gained from Refined Force-Field for Pure and Aqueous Ethylene Glycol through Molecular Dynamics Simulations Supreet Kaur, Shobhna, and Hemant K. Kashyap* Department of Chemistry, Indian Institute of Technology Delhi, Hauz Khas, New Delhi 110016, India

Downloaded via BUFFALO STATE on July 25, 2019 at 05:14:58 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: Recent experimental and first-principles simulation studies have shown that in liquid phase of ethylene glycol (EG), an equilibrium between both more prevalent gauche conformers and less probable trans conformers of EG molecule exists. Gleaning into the complexities faced during classical molecular dynamics (MD) simulations of condensed phase EG due to its conformational richness and considering the aforesaid observations, here we propose a refined force-field for EG molecule for atomistic MD studies. By employing the refined parameters, we have thoroughly investigated the structure and dynamics of pure EG liquid and its aqueous mixtures and compared the results with the available experimental data. The proposed force-field justifies the important role played by intra- and intermolecular hydrogen bonding rendered by EG molecules. The simulated X-ray scattering structure function for pure EG liquid at 298 K is found to be in excellent agreement with experimental X-ray scattering structure function which precisely confirms the ability of the proposed force-field to mimic the structure of liquid phase EG. Additionally, the accuracy of the refined force-field for the microscopic dynamics and self-diffusion coefficient of pure as well as aqueous EG were also assessed here. Temperature dependence of hydrogen bonding interactions and their dynamics reveals that with increasing temperature the intermolecular hydrogen bonds in pure ethylene glycol becomes weaker and consequently render faster dynamics. In aqueous mixture, intermolecular hydrogen bonding interaction between EG molecules tends to decrease with decrease in ethylene glycol mole fraction due to invasion of water. conformation.10,11 Similar inference has been given for the crystalline state, but ethylene glycol exhibits distinct behavior in liquid state. This nonidentical behavior of ethylene glycol in liquid phase has been a matter of debate in the recent past, and the origin is found to be the existing competition between intra- and intermolecular hydrogen bonding interactions rendered by ethylene glycol. Taking the conformational aspect into account, various experimental and theoretical studies have observed the preference of gauche isomer of liquid ethylene glycol with a scanty probability of trans isomer as well.12,13 The greater preference of gauche isomer in pure liquid as well as in the presence of various weakly polar solvents is also seen in early NMR and Raman study.14,15 The presence of trans isomer with a probability of ∼20% due to competing intermolecular hydrogen bonding and intramolecular hyperconjugation has been reported through recent ab initio molecular dynamics (MD) and nuclear magnetic resonance (NMR) studies.12,13 However, all these investigations emphasize establishing conformational equilibrium present in the liquid phase ethylene glycol irrespective of the overall molecular arrange-

1. INTRODUCTION Ethylene glycol (or 1,2-ethanediol) molecule serves as a prototype for study of various polyols of this family having multiple hydroxyl groups. It has a large number of applications, such as cryogenic liquid and heat transfer agent, in synthesis of various polymers and pharmaceuticals.1 In polymeric form, polyethylene glycol (PEG) serves as a macromolecular crowder.2 A plethora of investigations have been performed to address this phenomenon and also to study the effect of PEG on hydration structure and dynamics of proteins.3,4 In its aqueous form EG is also found to maintain the biological activity of DNA and thus exhibits the tendency for storing DNA for longer time periods.5−7 The application domain of ethylene glycol is not only restricted to its pure form, but also as a major component in deep eutectic solvents (DESs) that have appeared as designer solvents.8,9 Although ethylene glycol is a small molecule as compared to other polar molecules, it is found to be conformationally very rich, and this makes the thorough understanding of its liquid, vapor, and crystalline phases more complex. In the past, multiple experimental and theoretical techniques were employed to study the gas phase as well as condensed phase structure of EG, and it has been found that the structure of ethylene glycol differs in both these phases. Extensive investigation of vapor phase ethylene glycol reports the preference of this molecule to adapt a gauche © XXXX American Chemical Society

Received: April 27, 2019 Revised: June 20, 2019 Published: July 9, 2019 A

DOI: 10.1021/acs.jpcb.9b03950 J. Phys. Chem. B XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry B ment manifested in its bulk phase. Through the DFT calculations and X-ray and neutron scattering studies, Bakó and co-workers elucidated the microscopic structural features of ethylene glycol and depicted the preference for gauche conformers due to presence of intramolecular hydrogen bonding interactions in the molecule.16 Capturing all these macroscopic and microscopic characteristics of liquid ethylene glycol through molecular simulations has been a challenge for quite some time. Distinct classical models have been reported in the past that have predicted different behavior of condensed phase ethylene glycol in terms of microscopic structure and dynamics as well as conformational isomerism.7,17−21 Recently, by employing a modified OPLS-AA force-field, Kaiser et al. investigated structural and dynamic properties of pure ethylene glycol and its aqueous mixture.19,22 Saiz and co-workers compared different variants of other force-fields for liquid ethylene glycol, addressing the ability of different force-fields to capture different features of ethylene glycol in structure and dynamics.17 Recently Zhang et al. have reported the AMBER force-field parameters for the whole composition range of aqueous mixture of ethylene glycol, where intermolecular interactions were taken into consideration.7 Due to lack of reliable force-field parameters for ethylene glycol-based DESs, the computational investigations for these systems are scarce. However, it is very crucial to understand the microscopic structure and physicochemical properties of these novel systems due to the increasing demand of alternative solvents. Recently, the compatibility of CHARMM-based force-field parameters for the modeling of DESs have been shown.23,24 As pointed out earlier, ethylene glycol, being a notable hydrogen bond donor molecule, has gained popularity as a major component in DESs. However, the lack of reliable classical force-field for pure ethylene glycol is causing obstacles in obtaining a fundamental understanding of various interactions present in DESs. Despite several attempts to get an accurate model of liquid ethylene glycol for DESs modeling, one does not find any reliable force-field for the same that can be used for a quantitative comparison with experimental properties, such as X-ray scattering structure functions of ethylene glycol-based DESs.25,26 Since molecular simulation acts as a tool to gain a microscopic level understanding and to acquire physical insights with higher precision, it prompted us to propose refined force-field (FF) for ethylene glycol molecule by keeping both X-ray scattering structure and microscopic dynamics of liquid ethylene glycol at base. In the present study, we have justified the refined FF by comparing our results with the experimental data, wherever available. We have also assessed the accuracy of the proposed force-field for aqueous mixtures of ethylene glycol covering the whole composition range.

V=

1 2

∑

kr(r − r0)2 +

bonds

1 2

kθ(θ − θ0)2 +

∑ angles

∑ ∑ kϕ, j[1 + cos(njϕ − δj)] dihedrals

1 + 2 +

1 2

j

∑

kω(ω − ω0)2

impropers

∑

k ub(rub − rub0)2 +

ÄÅ É ÅÅi y12 i y6ÑÑÑ ÅÅjj σij zz σ j ij z ÑÑ ∑ 4ϵÅÅÅÅjjjj zzzz − jjjjj zzzzz ÑÑÑÑ + ÅÅ rij r Ñ LJ k ij { ÑÑÑÖ ÅÅÇk { Urey−Bradley

∑ Coulomb

1 qiqj 4π ϵ0 rij

(1)

The final refined force-field parameters defined according to the atom types shown in Figure 1 are given in Tables 1, 2, 3,

Figure 1. Chemical structure along with labeled atom type notation of ethylene glycol molecule.

Table 1. Partial Charges (q in |e| Unit) and Nonbonded LJ Parameters (σ in nm and ϵ in kJ mol−1) for Ethylene Glycol Used in the Simulation

a

atom type

q

σ

ϵ

CGa CGb OGa,b HGa HGb HOa,b

0.05 0.05 −0.65 0.09 0.09 0.42

0.360 0.358 0.314 0.221 0.238 0.040

0.216 0.234 0.804 0.378 0.146 0.192

Parameters used in present work. bCGenFF.27,30−32

Table 2. Parameters Defining the Bond Potential Energy Functions for Ethylene Glycol Systema bond type

r0

kr

CG−CGb,c CG−OGb,c CG−HGb,c OG−HOb,c

0.1530 0.1420 0.1111 0.0960

186188.00 358150.40 258571.20 456056.00

The units of r0 and kr are in nm and kJ mol−1 nm−2, respectively. Kindly note that these parameters are taken in original form from CGenFF.27,30−32 bParameters used in present work. cCGenFF.27,30−32

a

and 4. Herein, the refined FF parameters are chosen to mimic the conformational isomerism in order to obtain an equilibrium between the gauche and trans isomers in the case of liquid ethylene glycol. The dihedral distribution about the CG−CG bond obtained from the original CHARMM general force-field (CGenFF) parameters gives maximum population of trans isomer in liquid phase (Figure S1 of the Supporting

2. FORCE-FIELD MODEL AND ITS REFINEMENT The functional form of the FF proposed here is based on CHARMM FF given as27 B


Article

The Journal of Physical Chemistry B

calculations for both OG−CG−CG−OG and CG−CG−OG− HO dihedral, we have modified the parameters for these two dihedrals in a way to achieve the previously mentioned population ratio of gauche and trans isomers of EG in its liquid form. Along with this, further refinement is aimed at reproducing the experimental properties such as density, isothermal compressibility, thermal expansion coefficient, and self-diffusion coefficient in order to fit the Lennard-Jones (LJ) parameters of the EG atoms. At this point we have used the experimental X-ray scattering structure function to fix the remaining inadequacy in the parameters in order to obtain the final set of refined force-field parameters. Such kind of uncommon refinement of LJ parameters proves to be beneficial in parametrization of complex condensed phase systems and has been performed in the past.28,29 The final form of refined FF parameters gives us probability of finding gauche and trans isomers of ethylene glycol molecules as ∼83% and ∼17%, respectively (see Figure 2a). The equilibrium snapshot of simulation box for EG liquid shown in Figure 2b clearly depicts the dominating gauche conformers in the liquid along with the less probable trans isomers. Note that in the present study the partial atomic charges are kept the same as that of CGenFF.27,30−32 In order to asses the validity of this model, the simulated densities of pure ethylene glycol at different temperatures were compared with the experimental data,33 and the trend is shown in Figure 3a. The proposed force-field model for liquid ethylene glycol shows a remarkable match between simulated and experimental densities of pure ethylene glycol at all the corresponding temperatures. The deviation in the two sets of densities is within 0.5−1.1% range. The tabulated density comparison can be seen in Table 5. Also we have shown the comparison of simulated densities of all the compositions of aqueous ethylene glycol mixtures with the experimental data in Figure 3b at 298 K, which again is found to be a good match.33 In addition to this, the validity of refined FF model has also been checked by comparing other physicochemical properties, such as isothermal compressibility (βT), thermal expansion coefficient (αP),18 and self-diffusion coefficient (D) of liquid

Table 3. Parameters Defining the Angle Potential Energy Functions for Ethylene Glycol Used in the Simulationa θ0

kθ

rub0

kub

110.10 110.10 108.89 109.00 106.00

633.4576 221.7520 384.0912 297.064 418.400

0.00 0.2179 0.00 0.1802 0.00

0.00 18853.10 0.00 4518.72 0.00

angle b,c

CG−CG−OG CG−CG−HGb,c OG−CG−HGb,c HG−CG−HGb,c CG−OG−HOb,c

a The units for θ0, kθ, rub0, and kub are deg, kJ mol−1 rad−2, nm, and kJ mol−1 nm−2, respectively. Kindly note that these parameters are taken in original form from CGenFF.27,30−32 bParameters used in present work. cCGenFF.27,30−32

Table 4. Parameters Defining the Dihedral Potential Energy Functions for Ethylene Glycol Used in the Simulationa proper dihedral

δj

kϕ,j

nj

OG−CG−CG−HGb,c OG−CG−CG−OGb

0.0 180 0.0 180 0.0 0.0 180 0.0 180 0.0 0.0 0.0 0.0

0.81588 9.00 −5.736 7.800 0.8368 0.92048 1.656 −3.139 4.30 1.00416 0.58576 4.72792 0.75312

3 3 2 1 3 3 3 2 1 3 2 1 3

OG−CG−CG−OGc HG−CG−CG−HGb,c CG−CG−OG−HOb

CG−CG−OG−HOc

HG−CG−OG−HOb,c

The units for δj and kϕ,j are deg and kJ mol−1, respectively. Parameters used in the present work. cCGenFF.27,30−32

a

b

Information), which was not commensurate with recent experimental and theoretical observations.12−14 As mentioned previously, the liquid ethylene glycol consists of ∼80% cis isomers due to dominating intramolecular interactions. Starting from the dihedral angle potential parameters given by Cordeiro et al.,18 which is based on quantum chemical

Figure 2. (a) Probability distribution of OG−CG−CG−OG dihedral angle for liquid ethylene glycol using the refined FF at 298 K. The trans conformation is defined as ϕ ranging between ±180° to ±145° of OG−CG−CG−OG dihedral angle. (b) Equilibrium simulation box snapshot for pure ethylene glycol liquid at 298 K. Here the blue isosurface represents the gauche conformers. The oxygen, carbon, and hydrogen atoms of trans conformers are shown in red, cyan, and white color, respectively. C


Article


Figure 3. Comparison of simulated and experimental densities, ρ, for (a) pure ethylene glycol at different temperatures and (b) ethylene glycol aqueous mixtures at 298 K.

0.6, 0.4, 0.2, and 0.0, where the total number of molecules were maintained to 5000. All these systems were equilibrated for 15 ns in NPT ensemble at 298 K temperature and 1 bar pressure. Nosé−Hoover38−40 thermostat and Parrinello−Rahman41 barostat were used to maintain the temperature and pressure of these systems. Apart from this, temperature-dependent study was carried out for pure ethylene glycol for three temperatures: 298, 318, and 338 K. Pure water simulation was performed by employing usual methodology for SPC/E water model present in the literature.36 All the force-field parameters used in this study for ethylene glycol are given in Tables 1−4. PACKMOL42 package was used to generate the initial configuration for each system. Periodic boundary conditions in all three directions and minimum image convention were applied. Electrostatic interactions were evaluated using particle-mesh-Ewald (PME)43,44 summation technique. The equations of motion were solved by using the leapfrog algorithm with 1 fs time step. The cutoff radius for the short-range interactions was set to 1.2 nm with a switching function used from 1.0 to 1.2 nm. For the calculation of offdiagonal Lennard-Jones (LJ) interactions, Lorentz−Berthelot mixing rule was applied for all the systems except pure water. All the structural properties were analyzed using the last 5 ns of the equilibration run saved at every 100 fs. To analyze the microscopic dynamics of all the systems, NVT simulations for all compositions were carried out for at least 20 ns after the NPT run. The NVT trajectories were saved at every 50 fs time step for the computation of properties. Nosé−Hoover38−40 thermostat was used to maintain the desired temperature. As already shown, the equilibrium values of the simulated bulk densities at three temperatures for the systems studied were found to be in excellent agreement with experimental data (see Figure 3). The maximum deviation between the simulated and experimental densities were about 1.1% for the pure ethylene glycol and 0.7% for its aqueous mixtures. The X-ray scattering static structure function, S(q), and its subcomponents were computed using the methodology proposed in the literature.45−50 We first compute the radial distribution function (RDF), gij(r), for the atoms of type i and j. The computed RDFs include both intra- and intermolecular terms. By using these gij(r) values, the X-ray scattering S(q) can be computed as

Table 5. Comparison of Simulated and Experimental Densities, ρ, in g cm−3 for Pure as Well as Aqueous Mixtures of Ethylene Glycola xEG

T (K)

simulated ρ

exptl ρ

1.0

298 318 338 298 298 298 298 298

1.122 1.104 1.085 1.113 1.101 1.082 1.051 0.998

1.110 1.095 1.080 1.106 1.098 1.085 1.059 0.997

0.8 0.6 0.4 0.2 0.0 a

The experimental data has been taken from ref 33.

ethylene glycol at room temperature with the available experimental data34,35 (see Table 6) as well as with the other Table 6. Simulated Bulk Density, Isothermal Compressibility, Thermal Expansion Coefficient, and SelfDiffusion Coefficient of Pure Ethylene Glycol at 298 K and Their Comparison with Corresponding Experimental Values property

simulated

exptl

ρ (g cm−3) βT (10−5 bar−1) αP (10−4 K−1) D (10−6 cm2 s−1)

1.122 2.6 8.4 1.94

1.110 3.3 6.7 0.90

models proposed in the past (see Table S1 in the Supporting Information). Good agreement of these chief properties with the literature values gave us confidence to use this model to explore the structure and dynamics of this molecule as well as its aqueous mixtures in condensed phase as depicted in the upcoming sections.

3. SIMULATION DETAILS Large-scale all-atom MD simulations were performed by using the GPU version of GROMACS-5.1.1 (single-precision) package.36,37 In order to equilibrate pure ethylene glycol and its aqueous mixtures, we took six mole fractions xEG1.0, 0.8, D



ÅÄÅ L /2 n n 4πr 2ÅÅÅgij(r ) − ρo ∑i = 1 ∑ j = 1 xixjfi (q)f j (q)∫ ÅÇ 0 S(q) = ÄÅ ÉÄ É ÅÅÅ∑n x f (q)ÑÑÑÑÅÅÅÅ∑n x f (q)ÑÑÑÑ ÅÅ i = 1 i i ÑÑÅÅ j = 1 j j ÑÑ ÅÅÇ Ö ÑÖÅÇ

Article

ÑÉÑ sin qr 1ÑÑÑ qr ω(r )dr ÑÖ

(2)

In eq 2, xi and f i(q) are mole fraction and X-ray atomic form factor of atom of type i, respectively.51 ρo = Natom/⟨V ⟩ is the total number density, and L and V are the box length and volume, respectively. ω(r) is a Lorch window function, ω(r) = sin(2πr/L)/(2πr/L),52,53 which could be used to reduce the effects of finite integration error. Hydrogen bond autocorrelation function (C(t)) was calculated using54,55 C(t ) = ⟨h(0)h(t )⟩/⟨h(0)2 ⟩

(3)

where h(t) is hydrogen bond population variable which corresponds to 1 if a particular hydrogen bond was hydrogen bonded at time t = 0 and is still hydrogen bonded at time t; otherwise, it corresponds to 0. Note that this kind of definition of hydrogen bond autocorrelation function accounts for hydrogen bonds that are broken for a short time and formed again.21 The geometric criteria followed for calculation of C(t) is ∠HDA ≤ 30° and r(HA) ≤ 0.25 nm, which corresponds to intermolecular interactions between the hydrogen-bonded pair in EG. In order to get an accurate estimate of average hydrogen bond lifetime, the calculated C(t) is fitted by using a combination of an exponential function and a stretched exponential function given as C(t ) = a1exp[( −t /τ1)] + a 2exp[( −t /τ2)α ]

Figure 4. Experimental and simulated X-ray scattering structure functions (S(q)) for pure ethylene glycol at 298 K. The experimental data has been taken from ref 16.

the experimental data given by Bakó et al.16 By closely examining the intermolecular as well as intramolecular regions, we observe that the simulated S(q) is able to nicely encapsulate all the features of experimental X-ray scattering structure function, and overall very good agreement of the S(q) values is displayed. Gleaning the intermolecular region, we observe two peaks in the simulated X-ray S(q) centered around q = 1.57 Å−1 (R = 4 Å) and 2.6 Å−1 (R = 2.4 Å). Here, we name the first peak as the principal peak. The major contributions to this peak are from CG−CG and CG−OG correlations. Minor contributions from OG−OG, OG−HO, OG−HG, CG−HO, and CG−HG correlations are also present. In contrast, the second peak at around q = 2.6 Å−1 originates majorly due to OG−OG neighbor correlations along with minor contribution from OG−HO and OG−HG correlations. We have also compared the simulated X-ray S(q) obtained using CGenFF with that obtained from our refined FF and the experimental X-ray scattering structure function in Figure S2 of the Supporting Information where we see more inconsistencies in low q as well as high q features of S(q) obtained using original CGenFF. In summary, the positions and intensities of both the low q peaks are consistent with the experimental S(q), giving us further confidence to use the refined FF parameters for modeling liquid ethylene glycol. 4.1.2. Real Space Correlation Functions and Hydrogen Bonding. In order to appreciate the role of hydrogen bonding in molecular arrangement, we have also analyzed the RDFs corresponding to OG−OG and HO−OG pairs at 298 K shown in Figure 5a−b. The positions of first maximum and minimum corresponding to these RDFs and respective coordination number are also provided in Table 7. From Figure 5a, we observe that the closest distance peak in the total (inter- and intramolecular) OG−OG RDF occurs at rmax = 0.278 nm which is in agreement with the distance obtained from the neutron and X-ray scattering experiments performed by Bakó et al. 16 Likewise, the simulated coordination number is also very close to what is observed in the experiment.16 In the total OG−OG RDF, we also observe a comparatively less intense peak at 0.37 nm, which is apparently due to intramolecular OG−OG interaction. A close inspection of these peaks reveals that the more profound peak at 0.278 nm corresponds to OG−OG inter- and intramolecular interaction in gauche conformers, and the less profound peak at

(4)

where α is stretching exponent and 0 ≤ α ≤ 1 and a1 + a2 = 1. τ1 and τ2 are two time constants. Analytical integration of eq 4 gives us average lifetime as ⟨τ ⟩ = a1τ1 +

a 2τ2 ij 1 yz Γjj zz α kα{

(5)

where Γ is the gamma function. Self-diffusion coefficients of the species for all the xEG were computed with the help of mean-square displacement (MSD) using the Einstein relation given as56 Di =

1 d 2 lim ⟨| ri (⃗ t ) − ri (0) ⃗ |⟩ 6 t →∞ dt

(6)

where ri (⃗ t ) is the position coordinate of the ith particle at time t. To confirm the diffusive regime in the MSD of the species for different systems, we have also computed the β(t) function which is given as β (t ) =

2 d log⟨| ri (⃗ t ) − ri (0) ⃗ |⟩ d log t

(7)

4. RESULTS AND DISCUSSION 4.1. Accuracy of Refined FF for Pure Ethylene Glycol Liquid. 4.1.1. X-ray Scattering Structure Function, S(q). In this section, we have investigated the structural morphology of liquid ethylene glycol at 298 K and discussed the observations with the corresponding experimental inferences. X-ray or neutron scattering studies have proven to be powerful techniques to explore the microscopic structure of liquids. In Figure 4 we have shown the comparison of simulated X-ray scattering structure function, S(q), calculated using eq 2, with E


Article


Figure 5. Radial distribution function, g(r), for (a) OG−OG and (b) HO−OG pair. In the HO−OG RDF shown here directly bonded HO and OG pairs are excluded. Radial angular distribution function, g(r, θ), for (c) ∠HO−OG−OG angle and rOG−OG distance and (d) ∠HO−OG−OG angle and rHO−OG distance for the hydrogen-bonded HO-OG pair in pure ethylene glycol molecule at 298 K.

in Figure 5b. While intramolecular hydrogen bonding is exhibited only by the gauche conformers, the intermolecular hydrogen bonding peak at 0.194 nm has contribution from both gauche and trans isomers. The peak corresponding to the second coordination sphere at around 0.364 nm in total HO− OG RDF also has contribution from intramolecular hydrogen−oxygen interaction exhibited by trans conformers. The intensity is reduced for intermolecular case; however, the distance broadly remains the same. To support these observations, again we have shown the RADF for the hydrogen-bonded pair (HO−OG) in Figure 5d. It is evident that while intermolecular H-bonding occurs at shorter HO− OG distance and smaller ∠HO−OG−OG angle, the intramolecular H-bonding occurs at longer HO−OG distance and larger ∠HO−OG−OG angle. The presence of two different types of hydrogen bonding is in corroboration with the previous observations. These observations clearly depict the important role played by intramolecular interactions in laying down the structure of liquid ethylene glycol. It is also observed that the trans conformers have more tendency to remain aggregated with other trans EG molecules. Thus, small clusters of trans conformers are formed which are also clearly visible in the equilibrium snapshot shown in Figure 2b. Since CG−CG and CG−OG correlations are found to contribute in the principal peak present in total S(q), we have

Table 7. Positions (in nm) of First Maximum and Minimum in Total and Only Intermolecular RDF Pairs in Pure Ethylene Glycol Liquid and Corresponding Coordination Number at 298 K RDF pair OG−OG HO−OG

total intermolecular total intermolecular

rmax

rmin

ncoord

0.278 0.288 0.200 0.194

0.350 0.366 0.280 0.262

3.12 2.64 1.50 0.82

0.37 nm corresponds to OG−OG intramolecular interaction in trans conformers of EG. The intermolecular OG−OG correlation for the first coordination sphere is observed at 0.288 nm. This slight shift in the peak position for inter- and intramolecular correlation is also clearly visible in the radial angular distribution function (RADF), g(r, θ), shown in Figure 5c. In the RADF plot, the probability of finding intermolecular OG−OG pair is maximum at distance 0.288 nm and ∠HO− OG−OG ≤ 30°. In OG−OG RDF, one can also observe two peaks separated by a distance of 0.1 nm in the second coordination sphere corresponding to only the intermolecular interactions, irrespective of the any specific isomer. The other significant feature of EG in liquid phase is the exhibition of both intra- as well as intermolecular hydrogen bonding interactions which are manifested in the RDF shown F


Article

The Journal of Physical Chemistry B shown the RDFs for these two pairs in Figure 6. Broadly, both these correlations appear at the same length scale at which the

Figure 8. Distribution of H-bond number (NH−bond) per ethylene glycol molecule for total (intra- and intermolecular) hydrogen-bonded HO−OG pair in pure ethylene glycol at 298, 318, and 338 K.

Figure 6. Intermolecular radial distribution functions, g(r), for CG− CG and CG−OG pair in pure ethylene glycol at 298 K.

close to those obtained in the pioneering work by Padró et al.21 and previous literature.12,19,21 As anticipated, temperature causes the decrease in the average number of H-bonds in the system due to overall depreciation in the hydrogen bonding interactions in the system. Significance of H-bonding in this system can also be manifested through hydrogen bond dynamics analysis. Figure 9 shows the hydrogen bond

principal peak is observed in total S(q). Apparently, the refined force-field model is reliable enough to be able to capture the microscopic structural arrangement of liquid ethylene glycol. 4.1.3. Temperature-Dependent Behavior of Structure and Dynamics of Pure Ethylene Glycol Liquid. Looking into the temperature-dependent behavior of pure ethylene glycol, we observe from Figure 7 that on increasing the temperature from

Figure 7. Dependence of total X-ray scattering structure function, S(q), on temperature. A consecutive offset of 0.1 on ordinate is maintained to clearly show the minor shift in the position of the peaks with temperature.

Figure 9. Temperature dependence of hydrogen bond autocorrelation function, C(t), for intermolecular hydrogen bonding interaction for HO−OG pair in liquid ethylene glycol. Note that the abscissa has been plotted on a logarithmic scale.

298 to 338 K, only slight shift in the principal peak of total Xray S(q) toward lower q region is observed; however, the overall structural morphology is not significantly affected. It has already been discussed that H-bonding plays a vital role in laying down the structure of liquid ethylene glycol, so here we have explored the temperature-dependent behavior of statics and dynamics of H-bonding using refined FF. Figure 8 depicts the distribution of H-bond number for the total (intra- and intermolecular) hydrogen bonding interactions for HO−OG per ethylene glycol molecule and its change with temperature. The corresponding average number of H-bonds per ethylene glycol molecule comes out to be ∼3.9 (at 298 K) which is also

autocorrelation function, C(t), as a function of time that has been calculated using eq 3 for the intermolecular HO−OG hydrogen-bonded pair present in ethylene glycol. Also we have looked into the change of this correlation with altering the temperature. We can observe that the rate of decay of H-bond autocorrelation function increases with increase in temperature as the intermolecular hydrogen bonding interactions tend to depreciate with temperature. The corresponding fitting parameters and average lifetime are provided in Table 8. The amplitudes of the time scales are independent of temperature. While the α does not change significantly, the shorter as well as the longer time scales are found to vary with temperature. One G


Article

The Journal of Physical Chemistry B Table 8. Fitting Parameters Used To Fit Hydrogen Bond Autocorrelation Function for Intermolecular Hydrogen Bond Interaction Present in Pure Ethylene Glycol System at All the Three Temperatures and Its Respective Average Lifetime, ⟨τ⟩ T (K)

a1

a2

τ1 (ps)

τ2 (ps)

α

⟨τ⟩ (ps)

298 318 338

0.399 0.399 0.393

0.601 0.601 0.606

18.75 10.33 6.72

2.43 1.27 0.77

0.266 0.268 0.268

32.14 16.37 10.15

can clearly see the decrease in the overall lifetime with increasing the system’s temperature. We infer that with increasing temperature hydrogen bonding becomes weaker and its dynamics becomes faster. Also we see that the dynamics is faster than that reported by Padró and co-workers in early past investigation of EG using different force field.21 However, it is important to note that this correlation only corresponds to the intermolecular interactions, and we can observe a similar trend with temperature in the overall H-bond C(t) wherein intramolecular hydrogen bonds are also included. Being an important dynamical property to be considered while validating a force-field, mean-square displacement (MSD) as a function of time for ethylene glycol at all three temperatures is been shown in Figure 10. As we can observe that the dynamics become faster with increasing the temperature, these inferences are consistent with the previously discussed results in this section.

Figure 11. Simulated total X-ray scattering structure functions (S(q)) for all aqueous mixtures of ethylene glycol at 298 K.

can observe that by addition of water not only the intensity of principal peak decreases but also its position is slightly shifted toward higher q values. As already discussed, the principal peak originates due to the CG−CG and CG−OG correlations mainly, and we observe a depreciation in these correlations for lower mole fraction of ethylene glycol. However, the second peak gets more profound and also shifts toward higher q region. Eventually, both the peaks resemble pure water characteristic peaks.57,58 The simulated intermolecular RDFs for OG−OG (Figure 12a) and HO−OG (Figure 12b) pairs of ethylene glycol for all xEG reveal that the intermolecular hydrogen bonding interaction between ethylene glycol molecules tends to decrease with decrease in ethylene glycol mole fraction due to invasion of water. Similar observations have also been reported through simulation studies carried out using different force-field parameters in the past.7,59 Thus, it is clear that water molecules significantly hinder the overall structural arrangement of EG molecules. For further support, equilibrium snapshots of the simulation boxes for the complete range of aqueous mixtures are shown in Figure 13. In order to compare the translational dynamics of the aqueous mixtures of ethylene glycol, we have shown the comparison of simulated self-diffusion coefficients (D) of EG molecule for the whole composition range with the experimentally reported data at 298 K in Figure 14. The simulated self-diffusion coefficients have been computed from the respective MSD, as shown in Figure 14a, using eq 6. Note that these diffusion coefficients have been computed using MSDs in their respective diffusive regimes, which corresponds to the respective β(t) = 1 (Figure 14a). The fair agreement of the simulated and experimental D values for all mole fractions further validates the use of refined force-field for pure EG as well as the aqueous mixtures of EG (see Table 9).60The increase in the diffusion coefficient of the ethylene glycol molecule on adding water is due to invasion of water molecules, which tend to decrease the intermolecular interactions.

Figure 10. Mean-square displacement (MSD) and its first derivative β(t) plot for pure ethylene glycol at all temperatures. Note that both abscissa and ordinate of the upper panel have been plotted in logarithmic scale.

4.2. Performance of Refined FF for Aqueous Ethylene Glycol. Although the validation of the proposed refined FF for pure EG liquid is requisite to use this FF for all-atom MD simulation of pure EG, the compatibility of this FF with its aqueous mixture also validates its reliability. As shown in Figure 3b, the simulated densities for the aqueous mixtures also show a good agreement with the experimental densities for whole composition range. Since the structure of these mixtures has not been explored completely in the past, here we have computed the total X-ray scattering S(q) along with RDFs for the whole composition range at 298 K. From Figure 11 one

5. CONCLUDING REMARKS In present study, we have carried out refinement in the atomistic CHARMM general force-field to well mimic the microscopic structure and dynamics of liquid ethylene glycol. All possible comparisons of results obtained from the present study with the available literature have been done to asses the accuracy of the proposed force-field. The refined force-field H


Article


Figure 12. Intermolecular radial distribution functions, g(r), for (a) OG−OG and (b) HO−OG pairs for all the aqueous mixtures of ethylene glycol at 298 K.

Figure 13. Equilibrium snapshots of the simulation boxes of pure and aqueous ethylene glycol for xEG = (a) 1.0, (b) 0.8, (c) 0.6, (d) 0.4, and (e) 0.2. The gauche and trans conformers of ethylene glycol are shown in blue and white isosurfaces, respectively. The oxygen and hydrogen atoms of water molecule are shown in red and white molecular representations, respectively.

Figure 14. Composition dependence of (a) MSD and its first derivative β(t) and (b) self-diffusion coefficient (D) of ethylene glycol for all aqueous mixtures of ethylene glycol at 298 K.

conformers of ethylene glycol in its liquid state is also ensured. The less probable trans conformers are observed to be segregated in pure ethylene glycol liquid. The refined forcefield also predicts the structure and dynamics of aqueous mixtures of ethylene glycol with reasonable agreement with the experiments wherever available. Moreover, the new force-field parameters also allowed the detailed investigation of other

parameters predict multiple structural and dynamical properties of pure ethylene glycol in very good agreement with corresponding experimental data. The simulated total X-ray scattering structure function for pure ethylene glycol liquid has been found to be in excellent agreement with its experimental counterpart. Overall, as observed in recent experimental and ab initio simulation studies, the desired ratio of trans and gauche I


Article


Association by Intermolecular Excluded Volume Interactions. J. Phys. Chem. B 2011, 115, 2683−2689. (5) Eliasson, R.; Hammarsten, E.; Lindahl, T. Use of Ethylene Glycol in Counter Flow Electrophoresis and in Cold Concentration of DNA and Protein Solutions. Biotechnol. Bioeng. 1962, 4, 53−56. (6) Lindahl, T. The World of DNA in Glycol Solution. Nat. Rev. Mol. Cell Biol. 2016, 17, 335. (7) Zhang, N.; Li, M.-R.; Zhang, F.-S. Structure and Dynamics Properties of Liquid Ethylene Glycol from Molecular Dynamics Simulations. Chem. Phys. Lett. 2019, 718, 12−21. (8) Abbott, A. P.; Capper, G.; McKenzie, K. J.; Ryder, K. S. ̀̀ Alloys from Deep Eutectic Solvents Electrodeposition of ZincâTin Based on Choline Chloride. J. Electroanal. Chem. 2007, 599, 288−294. (9) Hizaddin, H. F.; Sarwono, M.; Hashim, M. A.; Alnashef, I. M.; Hadj-Kali, M. K. Coupling the Capabilities Oxide of Different Complexing Agents into Deep Eutectic Solvents to Enhance the Separation of Aromatics from Aliphatics. J. Chem. Thermodyn. 2015, 84, 67−75. (10) Guvench, O.; MacKerell, A. D. Quantum Mechanical Analysis of 1,2-Ethanediol Conformational Energetics and Hydrogen Bonding. J. Phys. Chem. A 2006, 110, 9934−9939. (11) Cramer, C. J.; Truhlar, D. G. Quantum Chemical Conformational Analysis of 1,2-Ethanediol: Correlation and Solvation Effects on the Tendency To Form Internal Hydrogen Bonds in the Gas Phase and in Aqueous Solution. J. Am. Chem. Soc. 1994, 116, 3892−3900. (12) Jindal, A.; Vasudevan, S. Conformation of Ethylene Glycol in the Liquid State: Intra- versus Intermolecular Interactions. J. Phys. Chem. B 2017, 121, 5595−5600. (13) Ghanghas, R.; Jindal, A.; Vasudevan, S. Distinguishing Intraand Intermolecular Interactions in Liquid 1,2-Ethanediol by 1H NMR and Ab Initio Molecular Dynamics. J. Phys. Chem. B 2018, 122, 9757−9762. (14) Pachler, K.; Wessels, P. Rotational Isomerism: X. A Nuclear Magnetic Resonance Study of 2-fluoro-ethanol and Ethylene Glycol. J. Mol. Struct. 1970, 6, 471−478. (15) Murli, C.; Lu, N.; Dong, Z.; Song, Y. Hydrogen Bonds and Conformations in Ethylene Glycol under Pressure. J. Phys. Chem. B 2012, 116, 12574−12580. (16) Bakó, I.; Grósz, T.; Pálinkás, G.; Bellissent-Funel, M. C. Ethylene Glycol Dimers in the Liquid Phase: A Study by X-ray and Neutron Diffraction. J. Chem. Phys. 2003, 118, 3215−3221. (17) Saiz, L.; Padró, J. A.; Guàrdia, E. Structure of Liquid Ethylene Glycol: A Molecular Dynamics Simulation Study with Different Force Fields. J. Chem. Phys. 2001, 114, 3187−3199. (18) Szefczyk, B.; Cordeiro, M. N. D. S. Physical Properties at the Base for the Development of an All-Atom Force Field for Ethylene Glycol. J. Phys. Chem. B 2011, 115, 3013−3019. (19) Kaiser, A.; Ismailova, O.; Koskela, A.; Huber, S. E.; Ritter, M.; Cosenza, B.; Benger, W.; Nazmutdinov, R.; Probst, M. Ethylene Glycol Revisited: Molecular Dynamics Simulations and Visualization of the Liquid and its Hydrogen-bond Network. J. Mol. Liq. 2014, 189, 20−29. (20) de Oliveira, O. V.; Freitas, L. C. G. Molecular Dynamics Simulation of Liquid Ethylene Glycol and Its Aqueous Solution. J. Mol. Struct. THEOCHEM 2005, 728, 179−187. (21) Padró, J.; Saiz, L.; Guàrdia, E. Hydrogen Bonding in Liquid Alcohols: A Computer Simulation Study. J. Mol. Struct. 1997, 416, 243−248. (22) Kaiser, A.; Ritter, M.; Nazmutdinov, R.; Probst, M. Hydrogen Bonding and Dielectric Spectra of Ethylene Glycol−Water Mixtures from Molecular Dynamics Simulations. J. Phys. Chem. B 2016, 120, 10515−10523. (23) Kumari, P.; Shobhna; Kaur, S.; Kashyap, H. K. Influence of Hydration on the Structure of Reline Deep Eutectic Solvent: A Molecular Dynamics Study. ACS Omega 2018, 3, 15246−15255. (24) Kaur, S.; Sharma, S.; Kashyap, H. K. Bulk and Interfacial Structures of Reline Deep Eutectic Solvent: A Molecular Dynamics Study. J. Chem. Phys. 2017, 147, 194507.

Table 9. Comparison of Simulated and Experimental SelfDiffusion Coefficient (D) of Ethylene Glycol (in 10−6 cm2 s−1) for Pure as Well as Aqueous Mixtures of Ethylene Glycola simulated

experimental

xEG

D

xEG

D

1.0 0.8 0.6 − 0.4 0.2

1.94 2.41 3.04 − 4.19 6.57

1.0 0.8 0.7 0.5 0.4 0.2

0.9 1.16 1.47 2.07 2.62 4.87

a

The experimental data has been taken from refs 35 and 60.

important properties, such as hydrogen bonding interactions and their dynamics and lifetime. The success of present ethylene glycol force-field also encourages us to extend its use in modeling the ethylene glycol-based DESs and other biological systems in a more justified manner.

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.9b03950. Figures for comparison of dihedral angle distribution and X-ray scattering structure functions; table for comparison of physical properties with literature (PDF)

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Phone: +91-(0)1126591518. Fax: +91-(0)11-26581102. ORCID

Supreet Kaur: 0000-0002-3756-0007 Hemant K. Kashyap: 0000-0001-9124-2918 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We sincerely thank Dr. Imre Bakó for providing us the X-ray scattering experimental data. S.K. and Shobhna thank UGC, India for fellowship. The authors thank the IIT Delhi HPC facility for computational resources. This work is financially supported by the Department of Science and Technology (DST), India, through FIST grant awarded to Department of Chemistry, IIT Delhi.

■

REFERENCES

(1) Yue, H.; Zhao, Y.; Ma, X.; Gong, J. Ethylene Glycol: Properties, Synthesis, and Applications. Chem. Soc. Rev. 2012, 41, 4218−4244. (2) Senske, M.; Törk, L.; Born, B.; Havenith, M.; Herrmann, C.; Ebbinghaus, S. Protein Stabilization by Macromolecular Crowding through Enthalpy Rather than Entropy. J. Am. Chem. Soc. 2014, 136, 9036−9041. (3) Samanta, N.; Luong, T. Q.; Das Mahanta, D.; Mitra, R. K.; Havenith, M. Effect of Short Chain Poly(ethylene glycol)s on the Hydration Structure and Dynamics around Human Serum Albumin. Langmuir 2016, 32, 831−837. (4) Rosen, J.; Kim, Y. C.; Mittal, J. Modest Protein-Crowder Attractive Interactions can Counteract Enhancement of Protein J


Article


Scattering of Imidazolium-Based Room-Temperature Ionic Liquids? J. Phys. Chem. B 2010, 114, 16838−16846. (46) Kashyap, H. K.; Santos, C. S.; Annapureddy, H. V. R.; Murthy, N. S.; Margulis, C. J.; Castner, E. W., Jr. Temperature-dependent Structure of Ionic Liquids: X-ray Scattering and Simulations. Faraday Discuss. 2012, 154, 133−143. (47) Santos, C. S.; Annapureddy, H. V. R.; Murthy, N. S.; Kashyap, H. K.; Castner, E. W., Jr.; Margulis, C. J. Temperature-dependent Structure of Methyltributylammonium Bis(trifluoromethylsulfonyl)amide: X Ray Scattering and Simulations. J. Chem. Phys. 2011, 134, 064501. (48) Kashyap, H. K.; Santos, C. S.; Daly, R. P.; Hettige, J. J.; Murthy, N. S.; Shirota, H.; Castner, E. W.; Margulis, C. J. How Does the Ionic Liquid Organizational Landscape Change When Nonpolar Cationic Alkyl Groups are Replaced by Polar Isoelectronic Diethers? J. Phys. Chem. B 2013, 117, 1130−1135. (49) Kashyap, H. K.; Santos, C. S.; Murthy, N. S.; Hettige, J. J.; Kerr, K.; Ramati, S.; Gwon, J.; Gohdo, M.; Lall-Ramnarine, S. I.; Wishart, J. F.; et al. Structure of 1-alkyl-1-methylpyrrolidinium Bis(trifluoromethylsulfonyl)amide Ionic Liquids with Linear, Branched, and Cyclic Alkyl Groups. J. Phys. Chem. B 2013, 117, 15328−15337. (50) Gupta, A.; Sharma, S.; Kashyap, H. K. Composition Dependent Structural Organization in Trihexyl(tetradecyl)phosphonium Chloride Ionic Liquid-methanol Mixtures. J. Chem. Phys. 2015, 142, 134503. (51) International Tables for Crystallography; Prince, E., Ed.; International Union of Crystallography, 2006; Vol. C. (52) Lorch, E. Neutron Diffraction by Germania, Silica and Radiation-Damaged Silica Glasses. J. Phys. C: Solid State Phys. 1969, 2, 229−237. (53) Du, J.; Benmore, C. J.; Corrales, R.; Hart, R. T.; Weber, J. K. R. A Molecular Dynamics Simulation Interpretation of Neutron and Xray Diffraction Measurements on Single Phase Y2O3−Al2O3 Glasses. J. Phys.: Condens. Matter 2009, 21, 205102. (54) Rapaport, D. Hydrogen Bonds in Water. Mol. Phys. 1983, 50, 1151−1162. (55) Luzar, A.; Chandler, D. Hydrogen-bond Kinetics in Liquid Water. Nature 1996, 379, 55−57. (56) Hansen, J.; McDonald, I. Theory of Simple Liquids; Academic Press: Cambridge, MA, 1986. (57) Kaur, S.; Kashyap, H. K. Three-dimensional Morphology and X-ray Scattering Structure of Aqueous Tert-butanol Mixtures: A Molecular Dynamics Study. J. Chem. Sci. 2017, 129, 103−116. (58) Gupta, A.; Kaur, S.; Kashyap, H. K. How Water Permutes the Structural Organization and Microscopic Dynamics of Cholinium Glycinate Biocompatible Ionic Liquid. J. Phys. Chem. B 2019, 123, 2057−2069. (59) Gubskaya, A. V.; Kusalik, P. G. Molecular Dynamics Simulation Study of Ethylene Glycol, Ethylenediamine, and 2-Aminoethanol. 2. Structure in Aqueous Solutions. J. Phys. Chem. A 2004, 108, 7165− 7178. (60) Ambrosone, L.; D’Errico, G.; Sartorio, R.; Costantino, L. Dynamic Properties of Aqueous Solutions of Ethylene Glycol Oligomers As Measured by the Pulsed Gradient Spin-echo NMR Technique at 25◦. J. Chem. Soc., Faraday Trans. 1997, 93, 3961−3966.

(25) Perkins, S. L.; Painter, P.; Colina, C. M. Experimental and Computational Studies of Choline Chloride-Based Deep Eutectic Solvents. J. Chem. Eng. Data 2014, 59, 3652−3662. (26) Ferreira, E. S. C.; Voroshylova, I. V.; Pereira, C. M.; Cordeiro, M. N. D. S. Improved Force Field Model for the Deep Eutectic Solvent Ethaline: Reliable Physicochemical Properties. J. Phys. Chem. B 2016, 120, 10124−10137. (27) Vanommeslaeghe, K.; Hatcher, E.; Acharya, C.; Kundu, S.; Zhong, S.; Shim, J.; Darian, E.; Guvench, O.; Lopes, P.; Vorobyov, I.; Mackerell, A. D. CHARMM General Force Field: A Force Field for Drug-like Molecules Compatible With the CHARMM All-atom Additive Biological Force Fields. J. Comput. Chem. 2009, 31, 671− 690. (28) Köddermann, T.; Paschek, D.; Ludwig, R. Molecular Dynamic Simulations of Ionic Liquids: A Reliable Description of Structure, Thermodynamics and Dynamics. ChemPhysChem 2007, 8, 2464− 2470. (29) Borodin, O. Polarizable Force Field Development and Molecular Dynamics Simulations of Ionic Liquids. J. Phys. Chem. B 2009, 113, 11463−11478. (30) Vanommeslaeghe, K.; Raman, E. P.; MacKerell, A. D. Automation of the CHARMM General Force Field (CGenFF) II: Assignment of Bonded Parameters and Partial Atomic Charges. J. Chem. Inf. Model. 2012, 52, 3155−3168. (31) Vanommeslaeghe, K.; MacKerell, A. D. Automation of the CHARMM General Force Field (CGenFF) I: Bond Perception and Atom Typing. J. Chem. Inf. Model. 2012, 52, 3144−3154. (32) Yu, W.; He, X.; Vanommeslaeghe, K.; MacKerell, A. D., Jr. Extension of the Charmm General Force Field to Sulfonyl-containing Compounds and Its Utility in Biomolecular Simulations. J. Comput. Chem. 2012, 33, 2451−2468. (33) Sun, T.; Teja, A. S. Density, Viscosity, and Thermal Conductivity of Aqueous Ethylene, Diethylene, and Triethylene Glycol Mixtures between 290 and 450 K. J. Chem. Eng. Data 2003, 48, 198−202. (34) Chandrasekhar, N.; Krebs, P. The Spectra and the Relative Yield of Solvated Electrons Produced by Resonant Photodetachment of Iodide Anion in Ethylene Glycol in the Temperature Range 296?t? 453 K. J. Chem. Phys. 2000, 112, 5910−5914. (35) Jerie, K.; Baranowski, A.; Przybylski, J.; Glinski, J. Ultrasonic and Positron Annihilation Studies of Liquid Solutions of N-hexanol in 1,2-ethanediol. J. Mol. Liq. 2004, 111, 25−31. (36) Abraham, M. J.; Murtola, T.; Schulz, R.; Pá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. (37) van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: Fast, Flexible, and Free. J. Comput. Chem. 2005, 26, 1701−1718. (38) Nosé, S. A Molecular Dynamics Method for Simulations in the Canonical Ensemble. Mol. Phys. 1984, 52, 255−268. (39) Nosé, S. A Unified Formulation of the Constant Temperature Molecular Dynamics Methods. J. Chem. Phys. 1984, 81, 511−519. (40) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695. (41) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (42) 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. (43) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An N.log(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (44) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. A Smooth Particle Mesh Ewald Method. J. Chem. Phys. 1995, 103, 8577−8593. (45) Annapureddy, H. V. R.; Kashyap, H. K.; De Biase, P. M.; Margulis, C. J. What Is the Origin of the Prepeak in the X-Ray K


Insights Gained from Refined Force-Field for Pure and Aqueous

Recommend Documents