Research Update: Ca doping effect on the Li-ion conductivity in NASICON-type solid electrolyte LiZr₂(PO₄)₃: A first-principles molecular dynamics study

Much scientific attention has been devoted to developing all-solid-state Li-ion batteries (LIBs) as alternatives to batteries that contain liquid electrolytes. To practically institute alternative batteries, solid-state electrolytes are required to improve their safety (e.g., against liquid leakage, flammability, explosion, and volatilization).^1–4 Major advances in solid-state LIBs require the discovery and characterization of new solid electrolyte materials.

A promising oxide-based solid electrolyte is Na super ionic conductor (NASICON)-type materials that are similar to perovskite-type and garnet-type ionic conductors. Some of the Li-ion-conducting NASICON-type undoped LiM₂(PO₄)₃ and heterovalent cation-doped Li_1+yM_2−xM′_x(PO₄)₃ (e.g., M = Ti⁴⁺, Ge⁴⁺, Zr⁴⁺, Hf⁴⁺; M′ = another metal; x and y = arbitrary values) solid electrolytes have been investigated.^5–21 Rhombohedral-phase LiZr₂(PO₄)₃ (LZP), shown in Fig. 1(a), is attractive because it forms an insulating interface consisting of Li₃P and Li₈ZrO₆ when in contact with Li metal, with no redox reactions up to 5.5 V.¹⁴ Li_1+2xCa_xZr_2−x(PO₄)₃ (LCZP) has exhibited a better Li-ion conductivity than pristine LZP.^15,16 The highest bulk and total ionic conductivities were 1.2 × 10⁻⁴ and 4.9 × 10⁻⁵ S/cm, respectively, at room temperature for Li_1.2Ca_0.1Zr_1.9(PO₄)₃. However, the origin of the enhanced Li-ion conductivity has not been completely clarified. It is necessary to gain fundamental knowledge of the structural properties and ionic conductivity of LCZP at an atomic level to develop solid-state electrolyte materials.

In this paper, we investigate the dependence of the Li-ion conductivity in LCZP on Ca substitution using first-principles molecular dynamics (FPMD) simulations and elucidate the Li-ion conduction mechanism in the system by analyzing the Li-ion trajectories. We show the methodology of first-principles calculations based on density functional theory (DFT). We describe the discussions of the resulting Li-ion conduction properties and analyze the Li-ion migration behavior in LCZP systems. We focus on the Li-ion behaviors at both the M1 and M2 sites [see Figs. 1(b) and 1(c)] which correspond to 6b and 18e Wyckoff positions in R$3¯$c space-group crystal structures, respectively. The target materials in this paper exhibit Li-excess compositions by replacing Zr⁴⁺ with 2 Li⁺ + Ca²⁺, and there is a possibility of Li-ion occupation at M1 and M2 sites.

We prepared five samples with large supercell structures of random Ca atoms and excess Li ions for each LCZP composition [Li_16+2nCa_nZr_32−nP₄₈O₁₉₂ with 1 ≤ n ≤ 6 corresponding to Li_1+2xCa_xZr_2−x(PO₄)₃ with 0.063 ≤ x ≤ 0.375] by referring to a previous study.²² We optimized their structures and evaluated their total energies using the Vienna ab initio simulation package (VASP)^23–26 using the projector augmented-wave method^27,28 and a plane-wave basis set. A generalized gradient approximation-type exchange-correlation functional developed by Perdew, Burke, and Ernzernhof and modified for solid materials (PBEsol)²⁹ was used. The cutoff energy for the plane-wave basis was set to 500 eV. A 2 × 2 × 2 Monkhorst-Pack k-point grid was employed for the optimization calculations. The energetically most stable structure in the five samples for each composition was adopted as the initial structure for our FPMD simulations. Next, we performed DFT-based FPMD simulations using VASP to estimate the Li-ion conductivity of the LCZP systems. The cutoff energy was set to 350 eV, and a 1 × 1 × 1 k-point grid (only Γ point) was employed to reduce the computational cost. The time step was set to 1 fs, and our FPMD simulations were carried out in the canonical (NVT) ensemble using a Nosé thermostat³⁰ at 873–1473 K for 50 ps. From our FPMD simulations, the mean square displacements (MSD) of all elements were calculated, and then the Li-ion conductivity at room temperature (298 K) and activation energy for Li-ion diffusion were determined from the Arrhenius plot of the diffusion coefficients.

Figure 2(a) shows the MSD plots at 1173 K for x = 0.188. The MSD profile of the Li ions shows a linear increase with respect to the MD step, indicating a good Li-ion diffusion. By contrast, the MSD profiles of the other ions (Ca, Zr, P, and O) remained constant; the ions remained near their original positions during the MD run, and the MSD values represent the magnitude of thermal vibration. Similar results were also obtained for other situations (different compositions and/or different temperatures for MD). Figure 2(b) shows the Arrhenius plot of the Li-ion diffusion coefficients for x = 0.188. From a slope and an intercept of the linear fitting line, we estimated the room temperature Li-ion conductivity and activation energy for each composition. The dependences of the Li-ion conductivity and the activation energy on the Ca content (x) are shown in Fig. 3. The Li-ion conductivity results [Fig. 3(a)] are substantially different from the experimental results because our computational results include only the bulk contribution, while the experimental results show the total ionic conductivities. Nevertheless, the Li-ion conductivity is generally maximized when x ≈ 0.100, which is the same as the experimental results. The computational results agree with the experimental results for the activation energy, especially for large amounts of Ca substitution. The difference in the activation energy (∼0.14–0.20 eV) between pristine LZP and LCZP [Fig. 3(b)] is almost the same as the formation energy of Frenkel-like defects in pristine LZP (∼0.13 eV).²² Since intrinsic defect (interstitial Li⁺) formation is required for stoichiometric LZP for Li-ion migration in the lattice, the calculated activation energy of LCZP contains the defect formation energy and the migration energy of Li ions.

We estimated a Li-ion probability density during 50 ps (50 000 steps) of FPMD simulations to clarify the Li-ion migration behavior (Fig. 4). Figure 4(a) shows that the three-dimensional pathways of Li-ions are clearly visible in the LCZP crystal structure and in the result of our previous work for pristine LZP. The probability density near the Ca atoms is still visible as a higher isosurface level is applied [Fig. 4(b)]. In detail, the Li-ion probability densities at the M1 and M2 sites near Ca [yellow and magenta points in Fig. 4(b)] are higher than those at the other sites. Li localization is not observed near the M2 sites in Ca-undoped LZP (not shown); therefore, Li ions in the M2 sites are stabilized by Ca doping. Therefore, the Li trapping effect by Ca ions is indicated, and we infer that the effect is an interaction between positively charged interstitial Li ions (⁠$Lii⋅$⁠) and negatively charged impurity defects of Ca ions (⁠$CaZr′′$⁠).

Furthermore, we estimated the RDFs to understand the Ca doping effect on the Li-ion behaviors. Figure 5 displays the RDFs of Li–Li, Ca–Li, and Zr–Li combinations for each composition. Figure 5(a) clearly shows that the RDF of Li–Li at around 3.5 Å increases as the amount of excess Li ions from Ca doping increases. This interaction agrees with the reported distance between the 36f sites in the same cavity (M1 sites).²² These Li ions trigger a so-called pushing-out mechanism of Li-ion conduction, which has been discussed in detail in our previous study.²² In the conduction mechanism, strong Coulombic repulsion between these Li ions pushes out one of the Li ions to neighboring M1 sites and causes another pushing-out event at the neighboring site. Focusing on the Ca–Li RDFs shown in Fig. 5(b), we find that there is little dependence of the RDF values on the Ca content. Moreover, the Li trapping effect by Ca atoms is confirmed because the Ca–Li RDFs [the peak value shown in Fig. 5(b) is about 0.030] are almost twice as large as the Zr–Li RDFs [the peak value shown in Fig. 5(c) is about 0.014] at around 3.3 Å; it shows higher occupancies of Li ions at both M1 and M2 sites near Ca than those at the other sites. The Ca content dependence of the RDF values can be clarified by examining the Zr–Li RDFs [Fig. 5(c)]. The Zr–Li RDF at around 4.0 Å gradually increases because of the increase in Li ions by Ca doping. The peak shift of Zr–Li RDF from 3.2 (x = 0.000) to 3.5 Å (x = 0.375) is visible, while that of Ca–Li RDF is smaller. It can be implied that the shift of the RDF peak is ascribed to the local Ca–Li interaction, which is stronger than the Zr–Li interaction, and not to the expansion of the lattice volume by Ca doping.

To investigate the Li-ion behaviors in FPMD simulations and clarify the Li-ion conduction mechanism in LCZP, we also computed SDF values at 64 sites (16 M1 sites and 48 M2 sites) using Eq. (2). The details of the SDF analysis are shown in Fig. 6. Figure 1(c) shows two MO₆ octahedra (M = Ca or Zr) near the M1 and M2 sites. Here, we discuss the different Li-ion behaviors of the surrounding environment with CaO₆ and that with ZrO₆. For M1 sites with two neighboring Zr ions, several Li-ion migration (hopping) events can occur, as indicated by the different color lines in the left panels of Fig. 6(a). In this respect, more frequent hopping events occur for M1 sites sandwiched between Ca and Zr atoms. In addition, simultaneous double Li occupation is visible at the similar M1 sites sandwiched between Ca and Zr atoms, while such double Li occupation at M1 sites sandwiched between two Zr atoms is hardly observed. For a clearer understanding, the number of Li ions at the same M1 sites is shown in the right panels of Fig. 6(a), and the average Li-ion occupancies are evaluated for each site using FPMD simulations. Accordingly, the average Li-ion occupancies at M1 sites near Ca ions are larger than those at the other M1 sites. A more pronounced tendency is observed for the SDF analysis at M2 sites, as presented in the left panels of Fig. 6(b), in which each Li-ion behavior can be verified by different color lines as well as the left panels of Fig. 6(a). In the figure, Li ions tend to stay at the M2 sites near Ca ions for more than several picoseconds, whereas Li ions pass through the other M2 sites that are surrounded only by Zr⁴⁺ ions. Therefore, quasi-stable sites are formed at the M2 sites near the Ca ions (other SDF results are shown in Figs. S5 and S6 of the supplementary material). This is quantitatively supported by the results for the average Li-ion occupancies at the M2 sites. The dependence on the average Li-ion occupancies at the M1 and M2 sites with/without neighboring Ca ions on the Ca content (x) is shown in Fig. 7. We found that the occupancies at sites near Ca atoms are larger than those at the other sites for both M1 and M2 sites, indicating the Li trapping effect of Ca ions. In addition, occupancies at M1 sites gradually decrease, whereas occupancies at M2 sites gradually increase as the Ca content increases, regardless of the neighboring Ca ions. Thus, the site potentials of M1 and M2 get close to each other by Ca doping. Moreover, we discussed the relationship between the bottleneck size of the Li-ion pathways and Li-ion conductivity, which in general show strong relationship, as discussed in previous papers.^11,34 However, we confirmed the less significant effect of the bottleneck size of the Li-ion pathways to the Li-ion conduction in LCZP (see Figs. S7 and S8 of the supplementary material). Therefore, the Li trapping effect of Ca ions strongly leads to the worsened Li-ion conduction in LCZP with the larger Ca content (x ≥ 0.188 in our theoretical simulations).

We estimated the Li-ion conductivity at room temperature and the activation energy of Li-ion diffusion for Li-excess Ca-doped LZP (LCZP) solid electrolytes, Li_1+2xCa_xZr_2−x(PO₄)₃ (0.063 ≤ x ≤ 0.375), using DFT-based FPMD simulations. The calculated ionic conductivities of Li in LCZP were much higher than those in pristine LZP because interstitial Li ions were created by Ca doping. When x ≤ 0.125 in Li_1+2xCa_xZr_2−x(PO₄)₃, the calculated ionic conductivity of Li increased because of the formation of interstitial Li ions at the 36f sites. The activation energy was reduced to ∼0.14–0.20 eV, which agrees with intrinsic defect formation energy (∼0.13 eV).²² However, further Ca doping (x ≥ 0.188) gradually decreased the ionic conductivity despite the increase in the interstitial Li concentration and good agreement with the experimental results from earlier studies.^15,16 The trajectory of the Li ions, RDF, and SDF analyses were adopted to clarify the effect of Ca doping on the Li-ion migration behavior. We confirmed that Li-ion migration was induced by a pushing-out mechanism triggered by simultaneous double Li occupancy at the M1 sites (interstitial Li ions at the 36f sites), as reported in our previous paper.²² In addition, we confirmed a strong Li trapping effect of Ca ions because of the negative charge of impurity defects of $CaZr′′$⁠, which also reduced the site potential of Li at the M2 sites. Therefore, increasing the Ca doping enhanced the Li trapping effect at the M1 and M2 sites with neighboring Ca ions, which decreased the Li-ion conductivity.

See supplementary material for detailed information on theoretical simulations and results.

This work was supported by the “Materials research by Information Integration” Initiative (MI²I) project of the Support Program for Starting Up Innovation Hub from the Japan Science and Technology Agency (JST) and by the “Elements Strategy Initiative: To Form Core Research Centers” of the Ministry of Education, Culture, Sports, Science and Technology (MEXT) of Japan (since 2012). We also thank the Information Technology Center of Nagoya University for providing computing resources (CX400, FX100). The crystal structures were drawn using the Visualization for Electronic and Structural Analysis (VESTA) program.³⁵