Computational Studies of Solubilities of LiO₂ and Li₂O₂ in Aprotic Solvents

All density functional theory (DFT) and MP2 calculations on solvated Li₂O₂ and LiO₂ clusters were done with Gaussian 09.²⁵ Structures were optimized at the B3LYP/6-31G(2df,p) level and frequencies were scaled by 0.9854 as in G4(MP2) theory. Single-point calculations with PBE/6-31G(2df,p) were performed at the B3LYP geometries for comparison. The PCM,²⁶ CPCM,^27,28 and SMD²⁹ continuum solvation models implemented in Gaussian09 were used to calculate solvation energies implicitly. Single-point MP2 calculations were done in addition to the DFT calculations. The MP2 calculations were done with frozen core according to the G4 convention³⁰ using the G3MP2LargeXP basis set as in G4(MP2) theory.³⁰

AIMD calculations were done using DFT methods with plane wave basis sets as implemented in the VASP code.^31–34 All the calculations were spin-polarized and carried out using the gradient corrected exchange-correlation functional of Perdew, Burke and Ernzerhof (PBE)³⁵ under the projector augmented wave (PAW)^36,37 method, with plane wave basis sets up to a kinetic energy cutoff of 300 eV to accommodate the large number of atoms (up to 324 in the simulation). A comparison of total radial distribution function vs radial distance calculated using 300 eV and 400 eV shows that the lower cutoff energy is sufficient for the system (Figure S4). The PAW method was used to represent the interaction between the core and valence electrons, and the Kohn-Sham valence states (i.e. 1s for H, 2s for Li, 2s 2p for C and O) are expanded in plane wave basis sets. The convergence criterion of the total energy was set to be within 1 × 10⁻⁴ eV within a single K-point grid (i.e. Γ-point). For the simulations, supercells of 12.0 Å × 12.0 Å × 12.0 Å containing 10 DME molecules and 15.0 Å × 15.0 Å × 15.0 Å containing 20 dimethyoxyethane (DME) molecules are adopted throughout the work. The DME densities corresponding to such a model are ∼0.88 and 0.90 g/cm³, respectively, close to the reported DME solvent density, i.e. 0.87 g/cm³. To investigate the thermodynamic stability of the system at room temperature, the 20 DME solvent molecules were initially randomly arranged inside the simulation box and were then thermally equilibrated at T = 300 K for about 2 ps using AIMD simulations based on an NVT-ensemble with a time step of 1 fs. For the production run, the sampling is obtained from a 4 ps AIMD simulation run.

Solubility is determined by the intrinsic stability of the species in its solid form and solvation by the particular solution. The solubility of substance A (S_A) can be calculated using the following equation:³⁸

where the dissolution free energy ΔG*_{dissolution, A} is the summation of sublimation free energy ΔG*_{sub, A} and solvation free energy ΔG*_{solv, A}. This is a thermodynamic cycle in which molecules of the crystal are first sublimed into the gas phase and then dissolved in the solvent. The ΔG*_{sub, A} is the free energy change of one unit of A dissociating from a bulk entity and it is typically endothermic. The ΔG*_{solv, A} is the solvation free energy of one unit of A and can be calculated using either an implicit continuum model, or an explicit solvent model where the solvent-solute interaction is calculated using ab initio methods. The computational details of these two terms are discussed further below. R and T are the gas constant and the temperature, respectively. S_A and V_A are the solubility (mol/L) and the molar volume (L/mol) of A in its solid structure, respectively. The molar volume V_A of LiO₂ in this work was obtained from a theoretical crystal structure³⁹ since there is no available experimental data in literature due to the instability of the structure under normal condition, and that of Li₂O₂ from the density and molecular weight.

Solvation energy

In this work, we use both implicit and explicit models to calculate solvation free energies of LiO₂ and Li₂O₂ in DME. With implicit solvation models, the free energy of the molecule in solution can be calculated directly from a dielectric continuum model. The solvation energy is calculated by subtracting the gas phase electronic energy from the energy in solution. In the continuum models, because DME is not a standard solvent for the method, we specified the solvent as dichloroethane which has the same shape as DME,⁴⁰ and the dielectric constant is the literature value of 7.2.

With the explicit models including 1∼4 DME solvent molecules, several different cluster structures constructed based on chemical intuition of the molecular interactions were tried and the lowest energy ones were used. It has been shown before that the complete of the first solvation shell is satisfactory for calculating solvation free energies of lithium-oxygen species.⁴¹ Since 4 DME molecules complete the first solvation, we did not calculate structures with more than 4 solvent molecules. These structures are shown in Figure 1. The first solvation shell of the Li₂O₂ from a snapshot from AIMD is also shown in Figure 1. Two different thermodynamic cycles, the monomer cycle and the cluster cycle,^41,42 were used to derive the solvation free energy. The two models use different reference energies. In the monomer cycle the solute is treated as binding to n solvent molecules while in the cluster cycle the solute binds to a cluster of n solvent molecules. Correction factors are added so that all reactants and products are in the same standard states, namely 1 M ideal gas and 1 M solution. For the monomer cycle, the solvation free energy is calculated as:

where ΔG^o_{bind, g, A} is the gas phase binding free energy between the solute and n DME molecules. Solvation free energies of the solute A with n surrounding explicit solvent DME molecules ΔG*_{solv, [A(DME)n]} and the solvation energy of DME molecule ΔG*_{solv, DME} were calculated with PCM by subtracting the gas phase electronic energy from the corresponding SCRF energy.C_DME is the concentration of 1 mol of DME and it is 9.6 mol/L. ΔG^{o − >}* is the free energy change to convert one mol of ideal gas from 1 atmosphere (24.46 L/mol) to 1M (1mol/L), which is 1.9 kcal/mol at 298 K.⁴² For all continuum models used (PCM, CPCM and SMD), this ΔG^{o − >}* correction is added.

Zoom In Zoom Out Reset image size

For the cluster cycle, the solvation free energy is:

where ΔG*_{solv, (DME)n} is the solvation energy of a cluster of DME surrounding the solute calculated using a continuum model and C_DME/n is the concentration of the solvent cluster.

We also calculated the solvation energy using the AIMD method. The computed solvation energy is approximated as the formation energy ΔE*_{solv, A}(T) from direct AIMD sampling at a given temperature (T) defined as the following:

where E^tot_{sol(DME, A)}(avg, T) is the statistical average of the total energy of a solution from AIMD sampling of the solute A with 10 or 20 surrounding DME solvent molecules thermally equilibrated at T = 298 K. Similarly, the E^tot_sol(DME)(avg, T) and E^tot_A(avg, T) are the total energies of a solution from AIMD sampling of the 10 or 20 DME solvent molecules and one A solute molecule thermally equilibrated at T = 298 K, respectively. More statistical information on the AIMD simulations can be found in the Supplemental Materials. In this work, we compared all other methods against the AIMD results as the model captures more complete solvation shell information and thus is the most physical and likely the most accurate. It has also been shown in literature that the solvation energies calculated using PBE AIMD are in good agreement with the experimental values.⁴³

Sublimation energy

The standard sublimation free energy ΔG^o_{sub, A} must be corrected for isothermal expansion of an ideal gas:⁴²

where P_o is the standard condition pressure. Molar volumes V_A for LiO₂ and Li₂O₂ are 0.0173 and 0.0199 L/mol, respectively.

The ΔG^o_{sub, A} of bulk is the summation of ΔH^o_{sub, A} and ΔS^o_{sub, A}. The ΔH^o_{sub, A} can be calculated according to:

where U_latt is the lattice energy calculated from crystal structures defined as U_latt = (E^tot_{crystal, A} − nE_{mol, A}^tot)/n with E^tot_{crystal, A} being the total energy of the A crystal that consists of n basic A molecular units in a crystal cell, and E^tot_{mol, A} being the corresponding total energy of a basic A molecular unit from a DFT calculation. Based on this relation, we obtained −1.85 and −3.26 eV for U_latt for LiO₂ and Li₂O₂, respectively. The 2RT term arises by including lattice vibrational energy and energy of the vapor.³⁸ The ΔS^o_{sub, A} is approximated by the difference between the rotation and translation contributions to the entropy of the gas phase at 298 K and the intermolecular vibrational contribution to the entropy of the crystal at 298 K. Details of these calculations were described in a previous publication.³⁹ With these entropy terms, ΔG^o_{sub, A} of LiO₂ and Li₂O₂, are calculated to be 1.47 and 2.83 eV, respectively.

Calculated free energies of solvation for the various methods, including the AIMD results, are reported in Table I with the optimized structures of Li⁺-DME_n cluster binding shown in Figure 1. The B3LYP and PBE energies for the cluster cycles (Table I, Rows 5–7) and the continuum models (Row 8–10) are very comparable (∼1 kcal/mol). With a monomer cycles (Table I, Row 1–4), the differences between ΔG*_{solv, A} with different numbers of DME molecules are noticeable (6–12 kcal/mol for B3LYP). When clusters of DME are used instead (cluster cycle), ΔG*_{solv, A} (Table I, Row 5–7) calculated using different cluster sizes are more similar. However, ΔG*_{solv, A} calculated using monomer cycles (Table I, Row 1–4) are in better agreement with the AIMD results (Row 11–12) with the 4DME monomer cycle (Row 4) giving the best agreement. For cluster and monomer cycles (Table I, row 1–7), MP2 solvation energies for LiO₂ are slightly more negative than the B3LYP results except for the 1DME monomer model (Row 1), while for Li₂O₂ the MP2 energies are less negative and the deviation is larger.

For the implicit solvation models PCM and CPCM, the results with both DFT and MP2 methods (Table I, Row 8 and 9) are similar to the AIMD results (Row 11 and 12), while for the SMD model (Row 10) results deviate more. This is because SMD is not parameterized for Li. The AIMD simulations using a larger box of 20 DME (Table I, Row 12) gives a slightly more negative ΔG*_{solv, A} for Li₂O₂ than the 10 DME (Row 11) box (−40.3 vs −38.4 kcal/mol). Li₂O₂ has a higher (more negative) ΔG*_{solv, A} than LiO₂ due to its stronger binding to DME. Considering the computational cost, the PCM and CPCM methods are preferable to cluster or monomer cycles for solvation energy calculation. In addition, it is not always clear if the lowest energy structure has been found.

Solubilities were calculated and are reported in Table II by combining theoretical ΔG*_{sub, A} and ΔG*_{solv, A} as outlined in the methods section. The AIMD solubility (Table II, Row 11 and 12) of LiO₂ is 1.8 × 10⁻² mol/L and that of Li₂O₂ is 9.1 × 10⁻¹⁹ mol/L (2.0 × 10⁻¹⁷ mol/L with the larger 20 DME box), indicating that bulk LiO₂ is more soluble than Li₂O₂ in DME. The PCM and CPCM values as well as the monomer cycle in Table II based on 4MDE are in agreement with these values. We note that because the solubility is an exponential function of solvation energies, small errors in solvation energy calculation will lead to a significant difference in the calculated absolute solubility. Our calculations showing that Li₂O₂ is very insoluble in DME and other solvents is consistent with some studies^19,20 and in disagreement with others^21–24 which are many orders of magnitudes larger (details on these experimental measurements were discussed earlier). For example, the solubility of 10⁻¹⁷ M reported by Lodge et al.¹⁹ and the insolubility of Li₂O₂ in ethanol reported in the Handbook of Inorganic Compouds²⁰ are consistent with our calculated results. Our calculations indicate that although the solvation energy of Li₂O₂ is quite large, it is the stability of Li₂O₂ (large lattice energy) that makes it insoluble. Although the solvation energies of the differing levels of theory differ, they all are in agreement on this qualitative result with values of 10⁻²¹ to 10⁻¹⁷ mol/L. The reports of large solubility of Li₂O₂ in the literature^21–24 could be for several reasons. One possibility is that since Li₂O₂ is known to react with H₂O, trace amounts of water contamination in solution might result in erroneous solubility measurement and in one case the solubility closely matched the ppm of water reported in the solvent.²² Another possibility is that solubility measurements may depend on particle size or surface effects that increase solubility.

As mentioned earlier, the experimental solubility of LiO₂ is not available. However, a recent investigation that uses LiO₂ solubility of 1 mM to model particle growth successfully simulated the particle size evolution observed experimentally.¹⁸ The solubility used for modeling is close to our AIMD result of 1.8 × 10⁻² mol/L, supporting our calculation result.

Solvation free energies of Li_xO₂ in some other common organic solvents with dielectrics ranging from 2 to 109 were also calculated using the CPCM model and are reported in Table III. Also included in Table III are the computed solubilities of LiO₂ and Li₂O₂ in these solvents. The solubilities in various other solvents are found to vary greatly (up to 15 orders of magnitude) depending on the dielectric constant. The solubility of Li₂O₂ ranges from 10⁻¹⁶ to 10⁻³² and that of LiO₂ ranges from 10⁻¹ to 10⁻¹³.

The solvation energies and solubilities of LiO₂ and Li₂O₂ in aprotic solvents were calculated using density functional theory and various solvent models for dimethoxy ethane (DME). The following conclusions are drawn:

• (1)  
  The AIMD simulations give free energies of solvation of −29.4 and −40.3 kcal/mol for LiO₂ and Li₂O₂, respectively. Solvation energies calculated using the monomer cycles with 4 DME solvent molecules and continuum models PCM and CPCM were found to be comparable with the AIMD results. The PCM and CPCM methods are recommended to be used for solvation energy calculation due to their relatively lower computational cost.
• (2)  
  The solvation energies were used with calculated free energies of sublimation to determine solubilities. The calculated solubilities of LiO₂ (1.8 × 10⁻² mol/L) and Li₂O₂ (2.0 × 10⁻¹⁷ mol/L) in DME indicates that bulk LiO₂ is slightly soluble and Li₂O₂ is very insoluble.
• (3)  
  The computed solubilities of Li_xO₂ (x = 1 or 2) in various organic solvents vary significantly depending on the dielectric constant of the solvent. For example, the computed solubility of LiO₂ in 1,4-dioxane (ɛ = 2.2), DME (ɛ = 7.2) and DMSO (ɛ = 46.5) is 10⁻¹³, 10⁻⁴ and 10⁻² mol/L, respectively. The solubilities of Li₂O₂ in the three solvents above are 10⁻³², 10⁻²⁰ and 10⁻¹⁷ mol/L, respectively.
• (4)  
  Experimental solubilities for Li₂O₂ range from 10⁻¹⁷ mol/L or "insoluble" to 10⁻³ mol/L. Our calculations agree with the former. The results giving larger solubilities may be due to impurities or surface effects. The results for LiO₂ agree with a mesoscale model for particle growth that successfully simulated particle size evolution observed experimentally.

The calculated LiO₂ and Li₂O₂ solubilities provide important information for fundamental studies of discharge and charge chemistries in Li-O₂ batteries. Since solubility is an important parameter that determines nucleation and growth of the discharge product,¹⁸ tuning the solubility of the LiO₂ and Li₂O₂ using different electrolytes by design can help achieving desired product morphology.

This work was supported by the Joint Center for Energy Storage Research (JCESR), an Energy Innovation Hub funded by the U. S. Department of Energy, Office of Science, Basic Energy Sciences. We gratefully acknowledge the computing resources provided on Blues and Fusion, a high-performance computing cluster operated by the Laboratory Computing Resource Center at Argonne National Laboratory. This research also used resources of the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility supported by the Office of Science of the U. S. Department of Energy.