Scaling Atomic Partial Charges of Carbonate ... - ACS Publications

Oct 21, 2016 - ... Applied Science and Technology, Politecnico di Torino, Turin 10129, Italy ... Tulane University, New Orleans, Louisiana 70118, Unit...
0 downloads 0 Views 1MB Size
Article pubs.acs.org/JCTC

Scaling Atomic Partial Charges of Carbonate Solvents for Lithium Ion Solvation and Diffusion Mangesh I. Chaudhari,*,† Jijeesh R. Nair,‡ Lawrence R. Pratt,¶ Fernando A. Soto,§ Perla B. Balbuena,§ and Susan B. Rempe*,† †

Center for Biological and Engineering Sciences, Sandia National Laboratories, Albuquerque, New Mexico 87185, United States Department of Applied Science and Technology, Politecnico di Torino, Turin 10129, Italy ¶ Department of Chemical and Biomolecular Engineering, Tulane University, New Orleans, Louisiana 70118, United States § Department of Chemical Engineering, Texas A&M University, College Station, Texas 77843, United States ‡

S Supporting Information *

ABSTRACT: Lithium-ion solvation and diffusion properties in ethylene carbonate (EC) and propylene carbonate (PC) were studied by molecular simulation, experiments, and electronic structure calculations. Studies carried out in water provide a reference for interpretation. Classical molecular dynamics simulation results are compared to ab initio molecular dynamics to assess nonpolarizable force field parameters for solvation structure of the carbonate solvents. Quasi-chemical theory (QCT) was adapted to take advantage of fourfold occupancy of the near-neighbor solvation structure observed in simulations and used to calculate solvation free energies. The computed free energy for transfer of Li+ to PC from water, based on electronic structure calculations with cluster-QCT, agrees with the experimental value. The simulation-based direct-QCT results with scaled partial charges agree with the electronic structure-based QCT values. The computed Li+/PF6− transference numbers of 0.35/0.65 (EC) and 0.31/0.69 (PC) agree well with NMR experimental values of 0.31/0.69 (EC) and 0.34/0.66 (PC) and similar values obtained here with impedance spectroscopy. These combined results demonstrate that solvent partial charges can be scaled in systems dominated by strong electrostatic interactions to achieve trends in ion solvation and transport properties that are comparable to ab initio and experimental results. Thus, the results support the use of scaled partial charges in simple, nonpolarizable force fields in future studies of these electrolyte solutions.



INTRODUCTION

The performance of lithium-ion batteries (LiBs) depends on efficient movement of Li+ within and between heterogeneous materials. Performance also depends on the interaction of Li+ with nonaqueous electrolytic solvents. The development of new materials to solvate and transport Li+ efficiently is a prime research focus in the battery field.1 Molecular simulations can guide the design of new materials2 and also facilitate the study of molecular-level phenomena that are difficult to probe experimentally. Those phenomena include ion solvation and transport in complex battery materials. We focus here on the solvation and transport of Li+ in carbonate solvents, specifically (Figure 1) ethylene carbonate (EC) and propylene carbonate (PC). Primary challenges for molecular simulations include the description of interionic forces important to ion solvation and the study of long time-scale phenomena like diffusion that require statistical sampling. Polarizable force fields have been used in preliminary studies of Li+ solvation structure and transport properties.3 Although the interionic forces may be © 2016 American Chemical Society

Figure 1. EC (left) and PC (right) structures and atom names.

reasonable, those force fields are specialized, not available in high-performance molecular simulation codes, and elude long time-scale investigations. Received: August 20, 2016 Published: October 21, 2016 5709

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation Table 1. Partial Charges on EC and PC Atomsa

Beyond force-field models with atomic polarizability, ab initio molecular dynamics (AIMD) methods provide higher resolution results at a higher computational cost and thus a more limited extent. Molecular simulation data are available on Li+ solvation structure in water4−6 and carbonate solvents7−10 using AIMD techniques. Although AIMD results provide meaningful structural data, the calculations are sufficiently costly that reliable evaluation of transport properties in these carbonate solvents is inaccessible. A compromise between highly accurate interionic forces and accessibility to long time-scale studies exists in the form of force field molecular dynamics (FFMD) simulations. Current nonpolarizable force fields, however, only partially reproduce experimental structural data.11 Ions and solvent molecules have strong electrostatic interactions that can affect the structure of solvent molecules around an ion. By understanding this behavior, we may be able to exploit it to tune the structure and transport properties of Li+ ion in organic carbonate electrolytes described by nonpolarizable force field parameters. Below we follow a simple plan to identify a force field for long time-scale FFMD calculations. We scale a standard set of partial charges for EC and PC and then compare the results obtained against an array of test quantities. The test quantities are obtained specifically from AIMD structural characteristics, from experiments for transport and other properties and from electronic structure assessments of ion solvation free energies. Scaling of standard partial charges has been followed before in studies of ionic liquids.13−16 That idea recognizes that electrostatic interactions dominate intermolecular forces for these systems. The array of test quantities obtained here is unusually diverse. For example, we measure Li+ transport properties with impedance spectroscopy to provide experimental results for low ion concentrations. Important and different test quantities are the ion solvation free energies obtained from electronic structure calculations on the basis of quasi-chemical theory (QCT). The QCT results obtained are remarkably accurate in comparison with experimental values for water and PC. Those QCT/electronic structure results give an independent and distinct comparison, particularly in the absence of experimental data on EC, for identification of the force field we seek. We also computed the density, molecular dipole, and dielectric constants of the carbonate solvents described by force fields to compare with experimental values. Additionally, as noted above, we obtained AIMD structural results and analyzed charge transfer for these systems to expand the comparison.

partial charges EC

PC

atom

100%

90%

80%

C3/C2 H Oe Oc Co C3 H Oe1 Co Oc Oe2 C2 C1

0.033 0.1041 −0.4684 −0.6452 1.0996 −0.0040 0.1165 −0.4509 1.0489 −0.6378 −0.4120 0.0832 −0.3264

0.0297 0.09369 −0.42156 −0.58068 0.98964 −0.0036 0.10485 −0.40581 0.94401 −0.57402 −0.3708 0.07488 −0.29376

0.0264 0.08328 −0.37472 −0.51616 0.87968 −0.0032 0.0932 −0.36072 0.83912 −0.51024 −0. 3296 0.06656 −0.26112

a

The default charges (100%) were taken from Soeten et al.12 and subsequently reduced to 90 and 80%. Refer to Figure 1 for atom identifiers.

chosen to be consistent with experimental solutions studied here. All FFMD simulations used the GROMACS molecular simulation package.19 The simulation procedure consisted of minimization of the system energy followed by 1 ns of dynamics to equilibrate solvent density and 50 ns of molecular dynamics production run. Prior to the simulations, solvent molecules were randomly placed in a periodic simulation cell. A Parrinello−Rahman20 barostat maintained atmospheric pressure conditions (P = 1 atm) in a constant particle number, pressure, and temperature (NPT) ensemble. A Nosé−Hoover21,22 thermostat maintained T = 300 K for simulations involving water and PC solvents. A higher temperature of T = 320 K was used for EC solvent to avoid freezing of EC molecules. In separate simulations, a single Li+···PF6− ion pair was solvated in 249 EC or PC molecules to calculate diffusion constants. After a 50 ns simulation, an additional 1 ns simulation was run to record high-resolution trajectory data (200 000 frames/ns). Properties were computed in postprocessing steps. Ion solvation structure was determined in the standard way by computing the radial distribution of water or carbonyl oxygen atoms around Li+. Solvent density, molecular dipole moment, and dielectric constant were calculated using utilities in the GROMACS software package. The approaches taken to compute diffusion constants, transference numbers, and solvation free energies are described below. Ab Initio Molecular Dynamics. Ab initio molecular dynamics simulations were carried out to investigate solution structures of a single Li+ ion in solvent. The experimental densities of the pure solvents were matched in a cubic cell of edge length between 12 and 15 Å with 30 carbonate or 64 water molecules. Specifically, 31 PC molecules and a single Li+ occupied a cubic box of ∼15 Å per side, and 32 EC occupied a box of ∼12 Å per side. The cubic box with water solvent also measured ∼12 Å per side. The number of solvent molecules was chosen to keep the number of atoms around ∼400 and thus the expensive AIMD simulations were computationally accessible. VASP23 software was used for these simulations. The core− valence electronic interactions were described using the improved projector augmented-wave (PAW) method.24 The



SIMULATION METHODS Molecular Dynamics. Standard force field molecular dynamics (FFMD) simulations were used to assess ion solvation structure and free energy and solvent properties. For those studies, we used the nonpolarizable all-atom optimized potentials for liquid simulations (OPLS-AA)17 to represent EC, PC, the Li+ cation, and the hexafluorophosphate anion PF6−. The partial charges on EC and PC molecules were those of Soeten et al.12 initially, which were scaled by 90% and then 80% to study the effect of partial charges on ion solvation structure and transport properties. The partial charges are tabulated in Table 1. Three systems were studied: a single Li+ ion solvated by 249 EC solvent molecules, 249 PC carbonate solvent molecules, or 999 SPC/E18 water molecules. Solution concentrations were 0.06 M in water and EC solvent and 0.05 M in PC solvent, 5710

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

Figure 2. Radial distribution function of water oxygens around Li+ using (left) FFMD and (right) AIMD simulations. The position of the first peak is closer to Li+ in the case of classical FFMD simulations (1.78 Å) compared with that of AIMD (1.95 Å), indicating a tighter first hydration shell. The first four near-neighbor contributions are shaded in gray, showing that the first solvation shell in both cases consists of the first four water neighbors. The FFMD results also resemble AIMD in that the next two waters occupy the second solvation shell (shaded red). The running coordination r numbers (⟨n(r)⟩ = 4πρ ∫ g(x)x2dx) also confirm four near neighbors in the first solvation shell from both simulations (red dashed line).



0

PBE generalized gradient approximation25 was used to describe the exchange-correlation energy in density functional theory. The valence electronic orbitals were expanded in plane waves with a high kinetic energy cutoff of 400 eV for accuracy. The temperature was maintained using the Nosé−Hoover thermostat21,22 at 300 K for simulations of PC and water solvents and 320 K for the EC solvent. Simulations of 100 ps were performed in a constant volume (NVT) ensemble with a time step of 1 fs. The last 20 ps of data were used for calculating radial distribution functions, g(r). To assess solvent polarizability, 50 configurations were sampled from the last 10 ps of the AIMD simulations, and charge transfer between solvent molecules and Li+ was evaluated by a Bader analysis.26 Electronic Structure Calculations. Clusters of Li+ ion with four neighboring solvent molecules were extracted from the AIMD trajectories as input for electronic structure calculations of ion solvation free energy. The Gaussian 0927 suite of programs was employed for these calculations. Using the PBE1PBE (PBE0) hybrid exchange-correlation density functional and the 6-31+G(d,p) basis set, the clusters were optimized first in isolation. PBE is a parameter-free generalized gradient approximation functional. The PBE0 version additionally contains a predefined amount of exact exchange.28 Recent work identified the PBE functional as reliable for predicting ion−water clustering free energies in vacuum.29 With a larger 6-311++G(3df,2pd) basis set, the cluster structures were reoptimized, and frequency calculations30 implemented to obtain accurate thermal free energy data. All vibrational frequencies were positive, confirming that the structures represented energy minima (Figure S1). The environment was modeled using experimental bulk dielectric constant values (Table S1) for the respective solvents in a polarizable continuum solvation model (PCM).31 Previous work on barium (Ba2+) ion hydration32 describes further details on the electronic structure simulation methods.

EXPERIMENTAL METHOD The Li+ ion transference number (t+) was measured using combined alternating current (AC) impedance and direct current (DC) polarization measurements of a Li/electrolyte− glasswool/Li cell.33,34 The test cells (EL−cell−Std) were purchased from EL-CELL GmbH (Germany), and the separators used were made of glass fiber with average thickness of 300 μm. The area of the cell was 2.54 cm2. While the lithium/lithium symmetric cells were being prepared, two separators soaked in 0.1 M solution of LiPF6 in EC or PC were used. For the EC-based symmetric cell, the cell was kept at 50 °C overnight to achieve a stable interface between the electrolyte and the lithium electrode. The same procedure was followed for PC, which was kept at 25 °C for the same amount of time. At those respective temperatures, a DC potential (ΔV = 10 mV) was applied until a steady current was obtained (generally 2−3 h). The electrochemical impedance spectroscopy (EIS) as well as the chronoamperometry measurements of the cells were performed using a PARSTAT-2273 potentiostat/galvanostat/ FRA (frequency response analyzer) instrument purchased from Princeton Applied Research (United States). The initial (I0) and steady (Is) currents that flow through the cell were measured. Simultaneously, impedance spectra were recorded (100 kHz and 0.1 Hz) with an oscillating potential of 10 mV before and after DC polarization. Subsequently, the initial (R0) and final (Rs) bulk resistances of the electrolyte and the initial (RC0) and final (RCs) charge transfer resistances of the interfacial layers of the Li metal electrode/electrolyte were derived. Using these measured values, t+ was calculated t+ =

IsR s(ΔV − I0 RC0) I0R 0(ΔV − Is RCs)

(1)

It is well-known that this method gives slightly higher values for transference numbers compared to NMR results.35,36



STRUCTURAL ANALYSES Water is the reference solvent for testing our methods and comparing structural results. Toward that end, we first analyze 5711

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

Figure 3. Radial distribution (g(r)) of carbonyl oxygens (Oc) in (top) EC and (bottom) PC around Li+ using (left) FFMD and (right) AIMD simulations. In the FFMD case, the partial charges on EC and PC molecules were reduced from 100 to 90% and subsequently to 80%. The first solvation shell becomes broader as the charges decrease, which allows additional solvent molecules to interact with the ion. Corresponding running coordination numbers (⟨n(r)⟩) are also plotted to show how the number of solvent molecules changes with a change in partial charges. As many as 2 additional solvent molecules are allowed into the first solvation shell when the charges are reduced to 80%. In all cases, contributions from the first four solvent molecules are sufficient to define the first peak solvation structure around Li+. AIMD results show that Li+ coordinates with only four solvent molecules.

our AIMD results for Li+ hydration in water and compare them with recent structural data37 to assess the accuracy of the AIMD predictions. A favorable comparison leads to new AIMD studies of Li+ solvation by EC and PC. Because our AIMD structural results for Li+ in PC agree well with recent X-ray absorption spectroscopy38 and neutron diffraction39 results, as described below, we use these validated AIMD results to compare solvent structure predicted by the nonpolarizable force field with scaled charges. Recent experimental neutron scattering studies of lithium chloride (LiCl) in water reported that Li+ ion is surrounded by 4−5 solvent molecules with an average occupancy in the first hydration shell of 4.85 (water) in solutions of 1 M concentrations.37 Our current AIMD results (Figure 2) describe 1 M solution concentrations and concur with the experimental characterization while also emphasizing a strong fourfold coordination structure for the first hydration shell. A slightly lower hydration number of 4.0 from AIMD simulations may be attributed to short simulation times (10−10 s) relative to the time of water exchange between first and second hydration shells (10−11 s). Nevertheless, the position of the first maximum

(1.95 Å) in the radial distribution of oxygen atoms around Li+ from AIMD matches the neutron scattering result (1.97 Å).37 FFMD studies generally result in tighter structuring of the first hydration shell around ions compared with that in the AIMD results.40,41 In the case studied here of Li+ in water (Figure 2), we also observe that the position of the first maximum in g(r) occurs closer (1.8 Å) to the ion in FFMD than AIMD or experimental results (1.95 to 1.97 Å). The stronger association of the first hydration shell in FFMD is also apparent from the higher peak height compared with that in the AIMD results (Figure 2). Nevertheless, both FFMD and AIMD simulations predict the same fourfold water occupancy in the first hydration shell with the next two waters in a nearby second hydration shell. Given the favorable comparison between Li+ hydration structure by AIMD and experimentally, we next computed Li+ solvation structure by the two carbonate solvents using AIMD simulations. The ion solvation structure by both EC and PC resembles the hydration structure (Figure 3). In particular, four carbonate solvents occupy the first solvation shell defined by a maximum in g(r) at 2.0 Å, similar to earlier AIMD studies.10 In comparison, studies of PC solvation structure by neutron 5712

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

Figure 4. Diffusion constants were calculated from FFMD using MSD of Li+ and PF6− ions in EC (left) and PC (right). MSD was calculated for 25 ps. The average diffusion constant was calculated using the portion where the derivative of MSD becomes constant over time (Figure S2). Experimental values are taken from Hayamizu et al.44 As evident from these plots, the diffusion constants match experimental values best for 80% of partial charges on EC and 90% on PC.

diffraction39 reported a maximum in g(r) of 2.04 Å. Additionally, both neutron diffraction39 and X-ray spectroscopy studies38 reported an average solvation structure of 4.5 PC molecules surrounding Li+. One difference between hydration and solvation by carbonates appears in the second solvation shell. Due to the larger size of EC and PC solvents relative to water, the second solvation shell is pushed out to larger distances from the ion. Similar to the hydration case, simulation by FFMD results in a closer (1.78 Å) and tighter structuring of the first solvation shell compared with AIMD (Figure 3). Yet, four carbonate molecules still occupy the first solvation shell, as in the case of water. Different Lennard−Jones force field parameters for the Li+ ion in the FFMD model produce a first peak at 2 Å, in agreement with experiment and our AIMD results, but also allow two additional carbonate molecules to interact with the ion.11 To mitigate the problem of solvation structure and peak positions in FFMD, we modified the partial charges on EC and PC molecules while keeping the partial charge on Li+ fixed. In separate calculations, Bader charge analysis26 of AIMD configurations suggests that solvent molecules donate 0.1 electrons to the Li+ ion. This observation emphasizes the natural approach to scale partial charges of solvents to incorporate polarization effects indirectly. The FFMD radial distributions of solvent carbonyl oxygens around Li+ in EC and PC change with alterations to solvent partial charges (Figure 3). As the partial charges decrease, the first solvation shell becomes less structured and shifts to larger separations, increasing agreement with the AIMD results. Less structuring of the solvation shell also allows additional solvent molecules to interact with the Li+ ion, lessening agreement with the AIMD results in this respect. Specifically, the first peak position of g(r) expands from 1.78 Å for the default 100% charges to 1.86 and 1.88 Å at 80% charge for EC and PC solvents, respectively. Running coordination numbers (Figure 3) show that the first four solvent molecules maintain strong interactions with the ion even when additional solvent molecules interact with the ion as the partial charges are lowered. Also, the carbonyl oxygen atoms orient toward Li+ ion to form a strong solvation structure

similar to the hydration structure for all values of solvent charges. Solvent properties of density, dipole moment, and dielectric constant decrease with lowered partial charges in the FFMD model (Table S1). At 80% partial charge on EC and 90% on PC, the computed dielectric constants of 68.8 and 72.3 come closest to the experimental values of 90 and 64,42 respectively, and to previous FFMD results obtained with different partial charges.43 For those same partial charge values, the dipole moments for EC and PC are similar (5.1 and 5.7 D), and the density of EC (1.3 kg/m3) is slightly higher than that of the larger PC molecule (1.2 kg/m 3 ), in agreement with experimental trends (Table S1). These individual changes are subtle. Though partial charges affect the dielectric constants substantially, all dielectric constant values are sufficiently high that the net effect is again subtle. As we discuss further below, the cluster-QCT results for Li+ solvation free energy put these various properties together in a convincing way. In summary, the entire FFMD first peak structure for all solvent charge cases results from the first four solvent molecules, in agreement with the AIMD results. This is an important result that will be used during solvation free energy calculations with quasi-chemical theory (see Solvation Free Energies).



DIFFUSION CONSTANTS Ion transport through solvent can be characterized by the selfdiffusion constant D. Ionic movement in the solvent depends on electrostatic interactions, dielectric constants, solvent viscosity, and ion size. Classical FFMD simulations can be used to calculate D of Li+ and the counterion PF6−, but computational expense prohibits AIMD-based estimates. Fortunately, experimental results from NMR measurements are available for comparison with FFMD predictions. We used the Einstein relation ⎛d⎞ 2 ⎜ ⎟⟨Δr (t ) ⟩ ∼ 6D ⎝ dt ⎠ i

(2)

to calculate D. Here, ⟨Δri(t)2⟩ is the mean-squared displacement (MSD) at time t for each ion, Li+ or PF6−, considered 5713

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

Figure 5. Transference numbers for Li+ and PF6− in EC (left) and PC (right) were calculated using FFMD and eq 3. The NMR diffusion constants at 1 M salt concentration were used to calculate experimental transference numbers (NMR),44 and the combined AC impedance and DC polarization results apply to 0.1 M salt concentration (EIS). The experimental results show little dependence on salt concentration, and the calculated results agree with the experimental numbers.

An analogous equation applies to the anion PF−6 , and the cationic and anionic transference numbers sum to unity. The computed transference values for the Li+ ···PF6− ions are plotted in Figure 5 for all solvent partial charges along with the experimental values measured earlier by NMR and here by impedance spectroscopy. Our previous conclusions that 90% charge on PC and 80% charge on EC result in the best ion diffusion constants holds also for ionic transference numbers. The computed Li+ ion transport numbers of 0.35 (EC) and 0.31 (PC) agree well with the NMR experimental values of 0.31 (EC) and 0.34 (PC).44 Similar values result from the impedance spectra measured for a solution 10-fold more dilute and emphasize the weak dependence of transference number on solution concentration. The diffusion constants of both ions change with the solvent partial charge due to changes in electrostatic interactions between ion and solvent molecules. This effect does not translate to ionic transference numbers. Nevertheless, the transference numbers are useful for validating simulation results on ion transport properties due to their experimental relevance. In summary, changes to the partial charges on PC and EC solvents alter solvation structure and transport properties of Li+ and PF−6 ions. On the basis of our results for radial distribution functions (Figure 3), diffusion constants (Figure 4), transference numbers (Figure 5), and solvent dielectric constants (Table S1), we chose 90% scaling of PC partial charges and 80% scaling for EC for the following analysis of Li+ solvation free energy.

individually. Calculation of diffusion constants can be sensitive to the time a particle is tracked or the method chosen to calculate MSD.8 In this study, we ran an additional 1 ns simulation to collect high-resolution data after 50 ns of production run. The Li+ or PF6− ions were then tracked for 25 ps each during MSD data collection. D was evaluated as the slope of the MSD with respect to time after diffusive behavior was verified (Figure S2). Diffusion constants for Li+ and PF6− in EC and PC solvents with different partial charges along with the experimental NMR results are plotted in Figure 4. The scaling of the partial charges clearly affects those coefficients. Ions move faster or slower depending on the equilibrium and kinetic characteristics of the solvent. In these simulations, both ions are present and sometimes come into contact with each other (Figure S3). Simulations of a single ion with scaled solvent partial charges show a trend not observed in the simulations that involve ion pairs: D increases as solvent partial charges decrease (Table S2). This trend suggests that ion−ion interactions affect D, but only in the prefactor. Although all the computed D values show the correct order of magnitude of 10−10 m2/s, the diffusion constants averaged between Li+ and PF6− in 80% charged EC and 90% charged PC come closest to the averaged experimental values. In ethylene carbonate, the experimental prefactors of D are 1.82 (Li+) and 3.63 (PF6−). In propylene carbonate, the analogous experimental values are 0.96 (Li+) and 1.84 (PF6−). Classical simulations with 80% scaling of the partial charges on EC predict diffusion constant prefactors of 2.12 (Li+) and 3.85 (PF6−) in EC, compared with 1.16 (Li+) and 2.55 (PF6−) with 90% partial charges on PC. Self-diffusion constants are slightly higher in EC relative to those in PC due to slightly smaller solvent coordination numbers in EC (Figure 3).



SOLVATION FREE ENERGIES Li hydration free energy, μ(ex), was studied previously using a cluster-continuum method,45 quasi-chemical theory (QCT),4,46 and experiments.47,48 In contrast, little data is available for Li+ solvation free energy in carbonate solvents. Experimental data is available only for the transfer free energy, ΔΔG, of Li+ from water to PC.49 Here, we apply two complementary QCT approaches to study Li+ solvation free energies in the carbonate solvents. Quasi-Chemical Theory. QCT is based on the potential distribution theorem50 and can be used to calculate the excess + chemical potential, μ(ex) Li+ , for a solution species such as Li . An 51−54 original QCT approach was built up directly from electronic structure calculation on clusters; here, we call that approach the cluster-QCT method. Cluster-QCT predicts +



TRANSFERENCE NUMBERS The ionic transference number, or transport number, is an experimentally accessible quantity that is directly related to the amount of current carried by an ion in the battery. The transference number for Li+, tLi+, can be calculated from ratios of D according to t Li+ =

D Li+ D Li+ + DPF−6

(3) 5714

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

Figure 6. Ion solvation free energy evaluated with n = 4 near-neighbor solvent molecules in the direct-QCT approach (eq 6). Solvent molecules are illustrated in blue, near-neighbor solvent molecules are in green, ions are in red, and the inner solvation shell of radius λ is in the dashed circle.

contrast to that common practice,60−62 here, we chose n = 4 to represent the first solvation shells and λ as the first minima in the ion−solvent g(r), as depicted in Figure 6. Note also that μ(ex) Li+ defines intrinsic solvation free energies that do not include any potentials of the phase.5,63 The potential disadvantage of this nλ = 4 choice (eq 6) is that evaluation of the outer shell contribution ln⟨eβε | nλ= 4⟩ might be challenging because the 4 near neighbors have to be included, too (Figure 6). The outer-shell contribution is typically approximated by

hydration free energies for monovalent and divalent metal ions in good agreement with tabulated experimental values.4,32,41,54−59 As an example, the cluster-QCT approach applied here produces a Li+ hydration free energy of μ(ex) Li+ = −116.6 kcal/mol, in close agreement with the experimental value adjusted for the same standard state, μ(ex) Li+ = −117.3 kcal/ mol.47 Here, we use the cluster-QCT method as a benchmark for comparison of QCT results for μ(ex) Li+ obtained from MD simulation with simple force fields. In cluster-QCT, the excess chemical potential is represented by the sum of three terms,

⎛ ⟨δε 2|n = 4⟩ ⎞ λ ⎟ ln⟨e βε|nλ = 4⟩ ≈ β⟨ε|nλ = 4⟩ + β 2⎜ 2 ⎝ ⎠

(0) n μLi(ex) + kT ln pLi+ (n) + = − kT ln K n ρ sol (ex) (ex) + (μLi(sol) ) + − nμ sol

(4)

n

which assumes a Gaussian distribution of binding energies, ε. Here, interestingly, the required conditional distribution of binding energies, P(ε|nλ = 4), is indeed accurately Gaussian (Figure S4). The binding energy ε here is calculated using 6−12 Lennard−Jones potentials and a generalized reaction field for Coulomb interactions64 with a cutoff at 1.2 nm. The physical conclusion of this Gaussian observation is that Li(sol)4+ clusters are sufficiently strongly bound so that only the quadratic well-bottom is required to describe these energetics of near-neighbor interactions. As usual, interactions longer-ranged than near-neighbor are numerous and weak, and that is additionally consistent with a Gaussian approximation for those energetics. The picture that results is a classic confirmation of the physical intuition upon which QCT theory is based. The advantage of using observed Li(sol)4+ clusters is that the required nλ = 4 condition permits a large data set to be applied for the outer-shell and chemical contributions of eq 6. Yet, the packing contribution of eq 6 is still manageable. Indeed, the initial contribution of eq 6 is still a “packing” contribution in the sense that occupancies nλ ≠ 4 are still excluded. Free Energy Results. Figure 7 compares free energies of Li+ transfer ΔΔGLi+ to carbonate solvent from water for the two QCT implementations: accurate electronic structure calculations and classical FFMD simulation with simple force fields. The cluster-QCT result for the free energy of Li+ transfer to PC from water agrees with the value tabulated by Marcus.49 No experimental value is available for EC. The direct-QCT evaluations of ΔΔGLi+ for the FFMD simulations agree reasonably well with cluster-QCT electronic structure calculations when the partial charges of the force fields are scaled by 80% (EC) or 90% (PC). The direct-QCT MD calculations agree to within 2 kcal/mol on the 10 kcal/mol difference in transfer free energies between PC and EC predicted by the cluster-QCT calculations.

each with clear physical meaning.51−54 The Boltzmann factor, k, and temperature, T, set the thermal energy scale. In the first term, K(0) is the equilibrium ratio for the Li+-solvent (sol) n association reaction Li+ + nsol ⇌ Li(sol)+n

(5)

treated as in an ideal gas phase, hence the superscript (0). The indicated solvent density, ρsol, describes the availability of solvent molecules serving as ligands in this association. In the second term, pLi+(n) is the probability of observing n solvent molecules in the defined inner solvent shell of radius λ. This population term contributes zero (0) if only the single coordination number n is realized in the given inner shell. The probability pLi+(n) is readily evaluated from AIMD simulations. The third term describes solvation of the Li(sol)n+ complex by the solvation environment external to it. That combination, (ex) μ(ex) Li(sol)n+ − nμsol , makes a favorable contribution to the free energy. For analyzing the results of MD simulation with simple force fields, we also use a direct-QCT approach,54 (0) βμLi(ex) (nλ) + ln⟨e βε|nλ⟩ + ln p(nλ) + = − ln p

(6)

−1

which is tautologically related to eq 4. Here, β = kT and ε are the binding energies of Li+. The three terms in eq 6, from left to right, are packing, outer-shell, and chemical contributions to the total free energy. Physical discussions of these terms are available elsewhere.54 The advantage of this simulation-based QCT is that it permits calculation of physically meaningful contributions to solvation free energy, μ(ex), for a range of λ defining the inner solvation shell, and this permits further insight into these thermodynamic characteristics. The QCT approaches (eqs 4 and 6) are independent of n and λ parameters. The choice of nλ = 0 (eq 6) produces a theory strongly connected to a van der Waals picture. In 52

(7)

5715

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

encouraging support for use of scaled partial charges in simple, nonpolarizable force fields in future studies of these electrolyte systems.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.6b00824. Solvent densities, dipole moments, dielectric constants, and Li+ diffusion constants calculated from simulation data, experimental values, MSD plots and Gaussian behavior of the outer-shell free energy contribution for one case, optimized structure (Li+(PC)4) used for cluster-QCT, plotted distances between Li+ and PF6 ions for different partial charges to show dissociation of ions during MSD calculations (PDF)



Figure 7. Transfer free energies, ΔΔGLi+, comparing FFMD directQCT results (left) and cluster-QCT results (right) using the Gaussian 09 electronic structure software package. The cluster-QCT results for the free energy of Li+ transfer to PC from water agree with tabulated experimental values to within 1 kcal/mol.49 Our experience using cluster-QCT to predict ion hydration free energies32,51,59 suggests that the Gaussian 09 results have accuracy comparable to that of other ab initio predictions,5,65−67 and hence provide a useful benchmark for FFMD results.

AUTHOR INFORMATION

Corresponding Authors

*E-mail: [email protected]. *E-mail: [email protected]. Funding

This work was supported by the Assistant Secretary for Energy Efficiency and Renewable Energy, Office of Vehicle Technologies of the U.S. Department of Energy under Contract DEAC02-05CH11231, Subcontract 7060634 under the Advanced Batteries Materials Research (BMR) Program and Sandia’s LDRD program (MIC and SBR).

Positive transfer free energies favor lithium-ion solvation by water compared with either carbonate solvent. From the perspective of the cluster-QCT calculations, the replacement free energy reflecting the availability of the solvent molecules as ligands is the foremost factor leading to that result. When EC and PC transfer free energies are compared, again from the perspective of cluster-QCT, the solvation of the bare EC/PC molecules serving as ligands is decisive in arriving at a positive free energy of transfer from PC to EC (EC being slightly smaller).

Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS Sandia National Laboratories is a multimission laboratory managed and operated by Sandia Corporation, a wholly owned subsidiary of Lockheed Martin Corporation, for the U.S. Department of Energy’s National Nuclear Security Administration under Contract DE-AC04-94AL85000. This work was performed, in part, at the Center for Integrated Nanotechnologies (CINT), an Office of Science User Facility operated for the U.S. DOE’s Office of Science by Los Alamos National Laboratory (Contract DE-AC52-06NA25296) and SNL.



CONCLUSIONS Electrostatic forces dominate Li+ solvation in carbonate-based solvents and in water. We tested three sets of scaled partial charges for EC and PC solvents, comparing them with respect to Li+ structural and transport properties. As might be expected, both structural and diffusion properties depend significantly on partial charges. Solvent partial charges scaled by 80% for EC and 90% for PC gave results that compare reasonably with experimental data. These force field models with scaled partial charges were used for calculating free energies of transfer of Li+ to the carbonate solvents from water, applying direct-QCT with n = 4. This variant of direct-QCT with nonzero occupancy allows gathering of voluminous statistical data. The perfect Gaussian fit to the binding energy and cancellation of packing and chemical terms during comparison makes it suitable for this analysis. Lacking the corresponding experimental thermodynamic data, solvation free energies based on electronic structure calculations with cluster-QCT were used for comparison. The computed free energy of transfer of Li+ to PC from water based on the electronic structure calculations with cluster-QCT agrees with the tabulated experimental value to within 1 kcal/ mol.49 The observed reasonable comparison between the electronic structure-based QCT value and the simulation-based direct-QCT with scaled partial charges provides further



REFERENCES

(1) Xu, K. Electrolytes and Interphases in Li-Ion Batteries and Beyond. Chem. Rev. 2014, 114, 11503−11618. (2) Leung, K.; Chaudhari, M. I.; Rempe, S. B.; Fenton, K. R.; Pratt, H. D., III; Staiger, C. L.; Nagasubramanian, G. Density Functional Theory and Conductivity Studies of Boron-Based Anion Receptors. J. Electrochem. Soc. 2015, 162, A1927−A1934. (3) 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. (4) Rempe, S. B.; Pratt, L. R.; Hummer, G.; Kress, J. D.; Martin, R. L.; Redondo, A. The Hydration Number of Li+ in Liquid Water. J. Am. Chem. Soc. 2000, 122, 966−967. (5) Leung, K.; Rempe, S. B.; von Lilienfeld, O. A. Ab Initio Molecular Dynamics Calculations of Ion Hydration Free Energies. J. Chem. Phys. 2009, 130, 204507−204517. (6) Alam, T. M.; Hart, D.; Rempe, S. L. B. Computing the 7Li NMR Chemical Shielding of Hydrated Li+ using Cluster Calculations and 5716

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation Time-averaged Configurations from Ab Initio Molecular Dynamics Simulations. Phys. Chem. Chem. Phys. 2011, 13, 13629−13637. (7) Silva, L. B.; Freitas, L. C. G. Structural and Thermodynamic Properties of Liquid Ethylene Carbonate and Propylene Carbonate by Monte Carlo Simulations. J. Mol. Struct.: THEOCHEM 2007, 806, 23− 34. (8) Ong, M. T.; Verners, O.; Draeger, E. W.; van Duin, A. C. T.; Lordi, V.; Pask, J. E. Lithium Ion Solvation and Diffusion in Bulk Organic Electrolytes from First-Principles and Classical Reactive Molecular Dynamics. J. Phys. Chem. B 2015, 119, 1535−1545. (9) Blint, R. J. Binding of Ether and Carbonyl Oxygens to Lithium Ion. J. Electrochem. Soc. 1995, 142, 696−702. (10) Ganesh, P.; Jiang, D.; Kent, P. R. C. Accurate Static and Dynamic Properties of Liquid Electrolytes for Li-Ion Batteries from Ab Initio Molecular Dynamics. J. Phys. Chem. B 2011, 115, 3085−3090. (11) Arslanargin, A.; Powers, A.; Beck, T. L.; Rick, S. W. Models of Ion Solvation Thermodynamics in Ethylene Carbonate and Propylene Carbonate. J. Phys. Chem. B 2016, 120, 1497−1508. (12) Soetens, J.-C.; Millot, C.; Maigret, B. Molecular Dynamics Simulation of Li+BF4− in Ethylene Carbonate, Propylene Carbonate, and Dimethyl Carbonate Solvents. J. Phys. Chem. A 1998, 102, 1055− 1061. (13) Schröder, C. Comparing Reduced Partial Charge Models with Polarizable Simulations of Ionic Liquids. Phys. Chem. Chem. Phys. 2012, 14, 3089−3102. (14) Chen, M.; Pendrill, R.; Widmalm, G.; Brady, J. W.; Wohlert, J. Molecular Dynamics Simulations of the Ionic Liquid 1-n-Butyl-3Methylimidazolium Chloride and Its Binary Mixtures with Ethanol. J. Chem. Theory Comput. 2014, 10, 4465−4479. (15) Youngs, T. G. A.; Hardacre, C. Application of Static Charge Transfer within an Ionic-Liquid Force Field and Its Effect on Structure and Dynamics. ChemPhysChem 2008, 9, 1548−1558. (16) Beichel, W.; Trapp, N.; Hauf, C.; Kohler, O.; Eickerling, G.; Scherer, W.; Krossing, I. Charge-Scaling Effect in Ionic Liquids from the Charge-Density Analysis of N,N′-Dimethylimidazolium Methylsulfate. Angew. Chem., Int. Ed. 2014, 53, 3143−3146. (17) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic Liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (18) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The Missing Term in Effective Pair Potentials. J. Phys. Chem. 1987, 91, 6269−6271. (19) 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. (20) Parrinello, M.; Rahman, A. Polymorphic Transitions in Single Crystals: A New Molecular Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (21) Nosé, S. A Molecular Dynamics Method for Simulations in the Canonical Ensemble. Mol. Phys. 1984, 52, 255−268. (22) Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions. Phys. Rev. A 1985, 31, 1695−1697. (23) Kresse, G.; Hafner, J. Ab Initio Molecular Dynamics for Liquid Metals. Phys. Rev. B: Condens. Matter Mater. Phys. 1993, 47, 558−561. (24) Perdew, J. P.; Chevary, J. A.; Vosko, S. H.; Jackson, K. A.; Pederson, M. R.; Singh, D. J.; Fiolhais, C. Atoms, Molecules, Solids, and Surfaces - Applications of the Generalized Gradient Approximation for Exchange and Correlation. Phys. Rev. B: Condens. Matter Mater. Phys. 1992, 46, 6671−6687. (25) Blöchl, P. Projector Augmented-Wave Method. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953−17979. (26) Tang, W.; Sanville, E.; Henkelman, G. A Grid-Based Bader Analysis Algorithm without Lattice Bias. J. Phys.: Condens. Matter 2009, 21, 084204. (27) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima,

T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, Jr, J. A.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö ; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09, revision E.01; Gaussian Inc.: Wallingford, CT, 2009. (28) Adamo, C.; Barone, V. Toward Reliable Density Functional Methods without Adjustable Parameters: The PBE0Model. J. Chem. Phys. 1999, 110, 6158−6170. (29) Soniat, M.; Rogers, D. M.; Rempe, S. B. Dispersion- and Exchange-Corrected Density Functional Theory for Sodium Ion Hydration. J. Chem. Theory Comput. 2015, 11, 2958−2967. (30) Rempe, S. B.; Jónsson, H. A Computational Exercise Illustrating Molecular Vibrations and Normal Modes. Chem. Educ. 1998, 3, 1−17. (31) Tomasi, J.; Mennucci, B.; Cammi, R. Quantum Mechanical Continuum Solvation Models. Chem. Rev. 2005, 105, 2999−3093. (32) Chaudhari, M. I.; Soniat, M.; Rempe, S. B. Octa-Coordination and the Aqueous Ba2+ Ion. J. Phys. Chem. B 2015, 119, 8746−8753. (33) Watanabe, M.; Nagano, S.; Sanui, K.; Ogata, N. Estimation of Li+ Transport Number in Polymer Electrolytes by the Combination of Complex Impedance and Potentiostatic Polarization Measurements. Solid State Ionics 1988, 28−30, 911−917. (34) Evans, J.; Vincent, C. A.; Bruce, P. G. Electrochemical Measurement of Transference Numbers in Polymer Electrolytes. Polymer 1987, 28, 2324−2328. (35) Zhao, J.; Wang, L.; He, X.; Wan, C.; Jiang, C. Determination of Lithium-Ion Transference Numbers in LiPF6-PC Solutions Based on Electrochemical Polarization and NMR Measurements. J. Electrochem. Soc. 2008, 155, A292−A296. (36) Zugmann, S.; Fleischmann, M.; Amereller, M.; Gschwind, R. M.; Wiemhöfer, H. D.; Gores, H. J. Measurement of Transference Numbers for Lithium Ion Electrolytes via Four Different Methods, a Comparative Study. Electrochim. Acta 2011, 56, 3926−3933. (37) Mason, P. E.; Ansell, S.; Neilson, G. W.; Rempe, S. B. Neutron Scattering Studies of the Hydration Structure of Li+. J. Phys. Chem. B 2015, 119, 2003−2009. (38) Smith, J. W.; Lam, R. K.; Sheardy, A. T.; Shih, O.; Rizzuto, A. M.; Borodin, O.; Harris, S. J.; Prendergast, D.; Saykally, R. J. X-Ray Absorption Spectroscopy of LiBF4 in Propylene Carbonate: A Model Lithium Ion Battery Electrolyte. Phys. Chem. Chem. Phys. 2014, 16, 23568−23575. (39) Kameda, Y.; Umebayashi, Y.; Takeuchi, M.; Wahab, M. A.; Fukuda, S.; Ishiguro, S.-i.; Sasaki, M.; Amo, Y.; Usuki, T. Solvation Structure of Li+ in Concentrated LiPF6-propylene Carbonate Solutions. J. Phys. Chem. B 2007, 111, 6104−6109. (40) Leung, K.; Rempe, S. B. Ab Initio Rigid Water: Effect on Water Structure, Ion Hydration, and Thermodynamics. Phys. Chem. Chem. Phys. 2006, 8, 2153−2162. (41) Sabo, D.; Jiao, D.; Varma, S.; Pratt, L. R.; Rempe, S. B. Case Study of Rb+(aq), Quasi-Chem. Theory of Ion Hydration, and the No Split Occupancies Rule. Annu. Rep. Prog. Chem., Sect. C: Phys. Chem. 2013, 109, 266−278. (42) Payne, R.; Theodorou, I. E. Dielectric Properties and Relaxation in Ethylene Carbonate and Propylene Carbonate. J. Phys. Chem. 1972, 76, 2892−2900. (43) You, X.; Chaudhari, M. I.; Rempe, S. B.; Pratt, L. R. Dielectric Relaxation of Ethylene Carbonate and Propylene Carbonate from Molecular Dynamics Simulations. J. Phys. Chem. B 2016, 120, 1849− 1853. (44) Hayamizu, K. Temperature Dependence of Self-Diffusion Coefficients of Ions and Solvents in Ethylene Carbonate, Propylene Carbonate, and Diethyl Carbonate Single Solutions and Ethylene Carbonate + Diethyl Carbonate Binary Solutions of LiPF6 Studied by NMR. J. Chem. Eng. Data 2012, 57, 2012−2017. 5717

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718

Article

Journal of Chemical Theory and Computation

(66) Ȧ qvist, J. Ion-Water Interaction Potentials Derived from Free Energy Perturbation Simulations. J. Phys. Chem. 1990, 94, 8021−8024. (67) Bhatt, M. D.; Cho, M.; Cho, K. Interaction of Li+ Ions with Ethylene Carbonate (EC): Density Functional Theory Calculations. Appl. Surf. Sci. 2010, 257, 1463−1468.

(45) Carvalho, N. F.; Pliego, J. R. Cluster-Continuum Quasichemical Theory Calculation of the Lithium Ion Solvation in Water, Acetonitrile and Dimethyl Sulfoxide: An Absolute Single-Ion Solvation Free Energy Scale. Phys. Chem. Chem. Phys. 2015, 17, 26745−26755. (46) Asthagiri, D.; Pratt, L. R.; Ashbaugh, H. S. Absolute Hydration Free Energies of Ions, Ion−Water Clusters, and Quasichemical Theory. J. Chem. Phys. 2003, 119, 2702−2708. (47) Marcus, Y. Thermodynamics of Solvation of Ions. Part 5. Gibbs Free Energy of Hydration at 298.15 K. J. Chem. Soc., Faraday Trans. 1991, 87, 2995−2999. (48) Tissandier, M. D.; Cowen, K. A.; Feng, W. Y.; Gundlach, E.; Cohen, M. H.; Earhart, A. D.; Coe, J. V.; Tuttle, T. R. The Proton’s Absolute Aqueous Enthalpy and Gibbs Free Energy of Solvation from Cluster-Ion Solvation Data. J. Phys. Chem. A 1998, 102, 7787−7794. (49) Marcus, Y. Thermodynamic Functions of Transfer of Single Ions from Water to Nonaqueous and Mixed Solvents: − Part 1  Gibbs Free Energies of Transfer to Nonaqueous Solvents. Pure Appl. Chem. 1983, 55, 977−021. (50) Beck, T. L.; Paulaitis, M. E.; Pratt, L. R. The Potential Distribution Theorem and Models of Molecular Solutions; Cambridge University Press: Cambridge, MA, 2006. (51) Pratt, L. R.; Rempe, S. B. Simulation and Theory of Electrostatic Interactions in Solution. In AIP Conf. Proc.; Hummer, G., Pratt, L. R., Eds.; AIP Press: New York, NY, 1999; Vol. 492, pp 172−201. (52) Asthagiri, D.; Dixit, P. D.; Merchant, S.; Paulaitis, M. E.; Pratt, L. R.; Rempe, S. B.; Varma, S. Ion Selectivity from Local Configurations of Ligands in Solutions and Ion Channels. Chem. Phys. Lett. 2010, 485, 1−7. (53) Rogers, D. M.; Rempe, S. B. Probing the Thermodynamics of Competitive Ion Binding Using Minimum Energy Structures. J. Phys. Chem. B 2011, 115, 9116−9129. (54) Rogers, D. M.; Jiao, D.; Pratt, L. R.; Rempe, S. B. Annu. Rep. Comput. Chem. 2012, 8, 71−127. (55) Rempe, S. B.; Pratt, L. R. The Hydration Number of Na+ in Liquid Water. Fluid Phase Equilib. 2001, 183, 121−132. (56) Rempe, S. B.; Asthagiri, D.; Pratt, L. R. Inner Shell Definition and Absolute Hydration Free Energy of K+(aq) on the Basis of QuasiChemical Theory and Ab Initio Molecular Dynamics. Phys. Chem. Chem. Phys. 2004, 6, 1966−1969. (57) Varma, S.; Rempe, S. B. Structural Transitions in Ion Coordination Driven by Changes in Competition for Ligand Binding. J. Am. Chem. Soc. 2008, 130, 15405−15419. (58) Jiao, D.; Leung, K.; Rempe, S. B.; Nenoff, T. M. First Principles Calculations of Atomic Nickel Redox Potentials and Dimerization Free Energies: A Study of Metal Nanoparticle Growth. J. Chem. Theory Comput. 2011, 7, 485−495. (59) Asthagiri, D.; Pratt, L. R.; Paulaitis, M. E.; Rempe, S. B. Hydration Structure and Free Energy of Biomolecularly Specific Aqueous Dications, Including Zn2+ and First Transition Row Metals. J. Am. Chem. Soc. 2004, 126, 1285−1289. (60) Sabo, D.; Varma, S.; Martin, M. G.; Rempe, S. B. Studies of the Thermodynamic Properties of Hydrogen Gas in Bulk Water. J. Phys. Chem. B 2008, 112, 867−876. (61) Jiao, D.; Rempe, S. B. CO2 Solvation Free Energy using QuasiChemical Theory. J. Chem. Phys. 2011, 134, 224506. (62) Chaudhari, M. I.; Sabo, D.; Pratt, L. R.; Rempe, S. B. Hydration of Kr(aq) in Dilute and Concentrated Solutions. J. Phys. Chem. B 2015, 119, 9098−9102. (63) You, X.; Chaudhari, M. I.; Pratt, L. R. In Aqua Incognita: Why Ice Floats on Water and Galileo 400 Years on; Nostro, P. L., Ninham, B., Eds.; Connor Court Press: Ballarat, 2014. (64) Tironi, I. G.; Sperb, R.; Smith, P. E.; van Gunsteren, W. F. A Generalized Reaction Field Method for Molecular Dynamics Simulations. J. Chem. Phys. 1995, 102, 5451. (65) Yanase, S.; Oi, T. Solvation of Lithium Ion in Organic Electrolyte Solutions and Its Isotopie Reduced Partition Function Ratios Studied by Ab Initio Molecular Orbital Method. J. Nucl. Sci. Technol. 2002, 39, 1060−1064. 5718

DOI: 10.1021/acs.jctc.6b00824 J. Chem. Theory Comput. 2016, 12, 5709−5718