Transport properties of imidazolium based ionic liquid electrolytes from molecular dynamics simulations

Ionic liquids (ILs) are promising electrolytes for high-performance Li-ion batteries (LIBs), which can significantly improve the safety and energy storage capac-ity. Although extensive experimental and computational studies have reported, further exploration is needed to understand the properties of IL systems, their microscopic structures and dynamics, and the behavior of Li ions in ILs. We report here results of molecular dynamics simulations as a function of electric field for Li diffusion in two IL systems, [EMIM][TFSI] and [BMIM][TFSI] doped with various concentrations of LiTFSI. We find that the migration of each individual Li ion depends largely on its micro-environment, leading to differences by factors of up to 100 in the diffusivity. The structural and dynamical properties indicate that Li diffusion is affected significantly by the coordination and interaction with the oxygen species in the TFSI anions. Moreover, the IL cations also contribute to the Li diffusion mechanism by attenuating the Li–TFSI interaction.


INTRODUCTION
Rechargeable lithium-ion batteries (LIBs) have become one of the most important and widely-used electrical energy storage devices since their successful commercialization in the early 1990s. 1 Despite their excellent properties, there is a growing demand for a new generation of electrochemical storage systems with superior properties suitable for a number of applications including transportation. However, many challenges remain that Full article doi.org/10.1002/elsa.202100007 stability, non-volatility/flammability, high ionic conductivity, and a wide electrochemical stability range. 4,5 Bis(trifluoromethanesulfonyl)imide (TFSI) are the most widely used and well-studied IL anions. 6 The flexible S−N−S bonds in TFSI lowers both the melting point and viscosity while the electron-withdrawing CF 3 SO 2 groups provide a large electrochemical stability window. 7,8 ILs combining of TFSI with imidazolium-based cations, such as 1-butyl-3-methylimidazolium (BMIM) or 1-ethyl-3-methylimidazolium chloride (EMIM), have received increasing attention due to their attractive properties such as relatively low viscosity and high conductivity. 8 Since performance of electrochemical cells depends significantly on electrolyte, it is crucial to fully understand the physico-chemical processes in ILs under electrochemical conditions. Extensive experimental and computational studies have been conducted to elucidate the structural and ion transport characteristics of ILs, [9][10][11][12][13] because they are directly related to macroscopic properties such as viscosity, electrical conductivity, and chemical stability. It has been reported that ion self-diffusion coefficients in typical [BMIM] [TFSI] or [EMIM][TFSI] IL are of the order of 10 -11 m 2 /s. 14,15 This diffusion is affected by the addition of Li salts, such as LiTFSI, LiNO 3 , LiPF 6, and LiBF 4 , with a gradually decrease observed as the Li salt concentration increases. 16 This occurs because the small size of the Li ion compared to that of the IL cation induces stronger interaction with anions. However, the IL ion diffusivity increases again at very high Li concentrations. 10 Due to the complexity of ILs compared to common solvents, experimental approaches produce limited information on their structures and diffusion mechanisms, in particular at the atomistic level. In contrast, computational simulations are a powerful tool to investigate the properties of such complex systems to generate important information for developing new efficient electrolytes. 10,11,[17][18][19][20][21][22][23][24] Here, molecular dynamics (MD) provides a most powerful means to investigate structural and dynamical properties of ILs, including dependence of ion behavior on additive salt concentration and external electric field. 10,19 In this paper, we report the results of MD studies of two IL systems: [ 27,28 which were minimized and equilibrated by the following procedure. The equilibration step consisted of a. NVT (constant number of molecules, volume, and temperature) MD simulation at 10 K for 10 ps to generate initial velocities, b. heating from 10 K to 300 K for 100 ps. For activation energy calculations, we also heated to 375 K and 450 K over 100 ps, c. NPT (constant number of molecules, pressure, and temperature) MD simulation at the target temperature for 1 ns to optimize cell parameters and density of the system, and d. NVT MD simulation for 1 ns with the density obtained in step c to equilibrate the system.
After equilibration, 5 ns of NVT MD was performed on the IL model structures. For nonpolarizable force fields, simulations of ILs with the ±1e ion charge often fail to reproduce the transport properties because of neglecting electronic effects such as polarization and charge transfer. 29 Charge scaling is an efficient way to approximately include those effects without further computational cost. 30 In this study, we applied a total charge of +0.7e to each BMIM or EMIM cation, and -0.7e for the TFSI anion. Li charge was also scaled to +0.7e for the charge neutrality. The scaled point charges on each atom were assigned based on the QM-based electrostatic potential (ESP) calculation ( Figure S1), where [Li][TFSI] and [TFSI][BMIM/EMIM] pairs were used, and all geometries were fully optimized at the B3LYP/6-31G(d,p) level using the Gaussian 09 program. 31 The transport properties of ions in MD simulations are usually measured by the self-diffusion coefficient (D) from the mean square displacement (MSD): where r i (t) is the position of ion i (center of mass for cations/anions) at time t and t 0 is the reference time. We identified the Fickian regime for each MSD curve (Figure S2) to calculate the ion diffusion coefficients. The largest domain, τ, where the log−log slope is nearly unity, is selected as the fitting region. We also calculated the effect of electric fields by applying DC fields to the IL bulk structures from 0 to 0.1 V/Å along the x direction. Electric fields were applied as a [TFSI] as a function of 1000/T. Activation energies for ion diffusions were estimated from the slopes of the linear fits static force to each particle of F = qE, where q is the charge on the particle and E is the electric field strength (see more detailed information on this procedure in the LAMMPS manual 26 ). We used periodic boundary conditions, with the local E equal to the applied electric field strength.

Ion diffusion in ionic liquids with various Li concentrations
The densities (ρ) of the selected IL systems, [BMIM] [TFSI] and [EMIM][TFSI], were calculated as a function of the Li concentration (Table 1). They are higher than the experimental ones, by 2% and 1%, respectively, 14,32 and increase gradually with increasing the Li concentration. This is in good agreement with reported experimental and theoretical studies 10,33 and can be interpreted in terms of stronger interactions between Li + and [TFSI] -. Also, the shorter carbon chain of EMIM led to a density 4% higher than To investigate the transport property of each ion species, self-diffusion coefficients (D) were calculated from the MSD plots ( Figure S2). Figures 1A and 1B show the predicted diffusion coefficients for various Li concentrations at 300 K. We also calculated the mobilities (μ) as where q -charge, k B -Boltzmann's constant, and T -temperature, and listed the corresponding values in Table S1.
For the pure IL systems, the diffusion coefficients of the cation and anion are 1.74 and 1.03 × 10 -11 m 2 /s for [BMIM][TFSI], and 3.79 and 2.67 × 10 -11 m 2 /s for [EMIM][TFSI], respectively, which are in good agreement with experimental results, ranging from 2.4 to 4.8 × 10 -11 m 2 /s for cations and from 2.0 to 2.9 × 10 -11 m 2 /s for anions. 15,16,34 The diffusion coefficients of the IL cation [TFSI], leading to diffusion constants from 1.0 × 10 -11 m 2 /s (pure IL) to 2.4 × 10 -11 m 2 /s (23%Li) for BMIM and from 0.6 × 10 -11 m 2 /s (pure IL) to 2.1 × 10 -11 m 2 /s (23%Li) for TFSI, respectively. 16 For both systems, we find that Li diffusion increases slightly from 0 to 10%Li, the highest diffusion coefficient value, and then decreases with concentration at higher loading ( Figures 1A and 1B). This trend is consistent with reported experimental results, where a distinct Lidiffusion peak was observed at 10 wt% for [EMIM][FSI] 10 and a much less pronounced Li-diffusion peak could be distinguished at 5 wt% for [BMIM] [TFSI]. 16 Note that the trend changes at very high Li content, where the Li diffusivity increases again at Li concentrations over 30 wt%. 10 The relative self-diffusion coefficients are D cation > D anion > D Li , which is consistent with reported experimental and theoretical data. 10,15,16 The small Li ions interact strongly with anions, resulting in slower diffusivity.
We estimated the activation energies from the Arrhenius plots at three temperatures (300, 375, and 450 K) for each IL system doped with 10%Li ( Figures 1C and 1D

Ion diffusion in ionic liquids under electric field
In real batteries, an electrolyte is sandwiched between two electrodes, producing a voltage across a cell. Thus, it is important to understand how IL ions behave under external electric fields (E-field). To investigate the effect of the E-field on transport properties of ions, we performed MD simulations under various E-field strengths, from 0 to 0.1 V/Å, along the x direction. Note that the E-field in real batteries are generally <0.001 V/Å. However, modeling these materials under strong E-fields is of interest to provide a deeper insight into the ion behavior under extreme conditions. Figure 2 shows the predicted diffusion coefficients for each ion species as a function of E-field strength with the corresponding calculated mobilities listed in Table S3. We find that for E-fields less than 10 -2 V/Å, ion diffusion in both IL systems is barely affected by the magnitude of E, whereas ion diffusion increases drastically for E-fields stronger than 10 -2 V/Å. An analogous trend was observed for all Li concentrations studied here ( Figure S3). Similar results have been reported for theoretical calculations on other systems. 19,38 Note that using a standard thermostat may cause interference along the E-field direction. However, as it was shown by Petravic, 39 this effect becomes apparent only at extremely high E-fields (>0.5 V/Å) and imposing a thermostat only on motion perpendicular to the E-field direction provides controlling the temperature of a system without any interference with the development of current. To be confident that the interference is negligible in our simulations, we performed additional computations on the [BMIM][TFSI] system with 10%Li under the 0.01 and 0.1 V/Å E-field (the strongest in our simulations), in which the thermostat was applied only along Full article doi.org/10.1002/elsa.202100007

F I G U R E 3 Radial distribution function (RDF) of the Li-N TFSI pair for the (A, B) BMIM-TFSI and (C, D) EMIM-TFSI systems at various Li concentrations (A, C) without E-field and (B, D) under various E-field strengths (V/Å) for the 10%Li concentration
the y and z directions perpendicular to the E-field along the x direction. We obtained the essentially same diffusion coefficients as before using the standard thermostat ( Figure S4).
Our result demonstrates that the effect of E-field on ion diffusion can be neglected for real energy storage devices, although the complexity of heterogenous environments, in particular at interfaces, should be analyzed carefully.

Structural analysis of ionic liquids
To explore how interactions of Li with IL anions and cations are modified from adding Li salt or applying E-field, we calculated the radial distribution functions (RDFs) (Figures S5 and S6). Among the calculated RDFs, the Li-N TFSI (nitrogen atom in TFSI) RDF shows the most significant changes (Figure 3). The RDF exhibits three major peaks at around 3.1, 4.3, and 5.1 Å, corresponding to different configurations between Li and TFSI molecules. In particular, we find different interactions with the oxygen species in TFSI ( Figure S7). Typically, the middle peak at 4.3 Å indicates that the Li ion interacts with two oxygens of a TFSI molecule (bidentate coordination), whereas the other two peaks indicate that the Li ion interacts with one oxygen of a TFSI molecule (monodentate coordination). The intensity of the middle peak gradually decreases with increasing Li concentration or E-field strength, whereas the intensities of the other peaks at 3.1 Å and 5.1 Å increase slightly. This change implies that the number of bidentate configurations decreases as the Li concentration or E-field strength increases. In addition, this change in the Li-TFSI interactions can be related to a flexible anion conformation, that is, cis and trans isomers. The C-S-S-C dihedral angle distribution in TFSI shows that the ratio of cis-isomer increases at higher Li concentrations and for stronger E-fields ( Figure S7), which is consistent with our earlier conclusion that the anion conformation can be changed by both Li concentration 10 and E-field strength. 19 We observed different peak posi- Figure S8), indicating that the length of the IL cation side chain affects the TFSI conformation.

Diffusion of individual Li in ionic liquids
To understand the Li transport in ILs in more detail, the self-diffusion of each individual Li ion was analyzed. We showed that the transport properties obtained from our 5 ns MD simulations in terms of densities, activation energies, and self-diffusion coefficients are in good agreement with experimental data. 14 But we find that diffusion of individual Li ions varies by up to two orders of magnitude  (Figures 4C, 4D and S10). This indicates that the fast Li ions tend to interact with only one oxygen of a TFSI molecule, which is apparently more favorable for fast ion diffusion than the more rigid interaction with two oxygens of a TFSI molecule, which corresponds to the peak at 4.3 Å.
The 5 ns MSDs and trajectories of the most and least diffusive Li ions with the surrounding IL molecules are shown in Figure 5. The least diffusive Li ion occupies a position in a closed rigid cage formed by three TFSI molecules with bidentate coordination. In contrast, the most diffusive Li ion is in an open flexible cage of TFSI molecules, allowing it to move out of the cage at certain favorable conditions. This makes this cage more open due to interactions with TFSI anions and neighboring IL cations. During the 5 ns MD simulation, the least diffusive Li ion dangled within the closed cage of the TFSI molecules, moving slightly together with the cage (vehicular diffusion) ( Figure 5B), whereas the most diffusive Li ion undergoes two big jumps, corresponding to hopping from one cage to another. Figure 5C shows the trajectory of the most diffusive Li ion, which stretches over three different cages. It can be partitioned into three intervals: 0-1.8 ns, 1.8-2.9 ns, and 2.9-5 ns. Our results indicate that individual Li diffusion depends strongly on its microenvironment with two types of the diffusion mechanisms for Li ions in ILs: vehicular (slow diffusion) and hopping (fast diffusion). These mechanisms were first described in the study of solutes in ILs, 40,41 and then similar results have been reported in recent experimental and computational works. 20,42 We also observed very different diffusion behaviors of individual Li ions under E-field ( Figure 6). For the zero and 0.01 V/Å E-field cases, the least diffusive Li ion barely moved during the 5 ns MD simulations, while the most diffusive Li ions showed significant migration through hopping several times. It should be noted that the E-field was applied along the x-axis, whereas the most diffusive Li ion migrated roughly along the z-axis. This demonstrates that an E-field < 0.01 V/Å is too weak to affect the structural and dynamic properties of the systems. Under stronger Efields, 0.05 V/Å and 0.1 V/Å, the fast transport of the most diffusive Li ion occurs by frequent hopping along the Efield direction. Interestingly, even under these high E-field strengths, the least diffusive Li ion moves slowly, indicating that interactions between the Li and IL ions governed by the micro-environments significantly affect the Li diffusion.

CONCLUSION
We Our computational studies reveal that Li diffusion is affected significantly by interactions between the IL ions, in particular by the coordination and interaction of Li with the oxygen species of the IL anions. Moreover, the IL cations contribute to Li transport by interacting with the TFSI anions to create favorable conditions for Li hopping. The fast Li-ion diffusion occurs via the hopping mechanism, with Li-ions jumping from one cage to another, while the Li-ion diffusion via the vehicular mechanism is slow. The ion diffusion is affected only slightly under E-fields less than 0.05 V/Å.
The complexity of heterogenous environments and interfaces should be considered carefully to fully understand the effect of the E-field on the structure and properties of an electrochemical cell.
These insights into the Li diffusion process at the atomistic level provide essential information for developing advanced electrolyte materials for a new generation of high performance LIBs.

A C K N O W L E D G M E N T S
This work was initiated by a research grant from LG Chem. It was completed with support provided by InnoHK. We used the Extreme Science and Engineering Discovery Environment (XSEDE) system, which is