The Role of Local Inhomogeneities on Dendrite Growth in LLZO-Based Solid Electrolytes

The entire lower length scale model framework has been developed based on Density Functional Theory (DFT), Ab-Initio Molecular Dynamics (AIMD), Monte Carlo (MC) and classical Molecular Dynamics (MD) techniques. Appropriate methodologies have been evoked for predicting different parameters of interest. For example, a combination of DFT and AIMD has been used for determining the LLZO surface facet with lowest energy with respect to lithium metal.⁴¹ Classical MD has been used to predict the Young's modulus at the GB regions of the electrolyte,⁹ which is an almost impossible task with the AIMD techniques due to the computational limitations associated with it. A description of the various atomistic computational techniques mentioned above is provided below.

First principles calculations have been performed using the Vienna Ab-initio Simulation Package (VASP)⁴² code with plane wave basis sets and projector-augmented wave (PAW) pseudo-potentials.⁴³ The exchange-correlation functional has been treated within the generalized gradient approximation (GGA) of Perdew-Burke-Ernzerhof (PBE).⁴⁴ The cutoff energy for the plane-wave basis has been set at 600 eV, which was tested and applied for all super-cells. The convergence criterion of the total energy has been set to be within 1 × 10⁻⁵ eV within the K-point integration, and all the geometries were optimized until the residual forces becomes less than 1 × 10⁻² eV Å⁻¹. For the interface structure with extended super-cells in each phase, a 5 × 1 × 1 Monkhorst-Pack k-points mesh has been used in the self-consistent calculations, whereas 3 × 3 × 3 Monkhorst-Pack k-points mesh has been applied to the bulk and surface calculations.

To investigate the interaction of a Li metal anode and LLZO surface, ab initio molecular dynamics (AIMD) simulations have been performed using the VASP code.⁴² The exchange-correlation functional is treated within the generalized gradient approximation (GGA) of Perdew-Burke-Ernzerhof (PBE)⁴³ with Γ-point Brillouin zone sampling and a 600 eV plane-wave energy cutoff. A time step of 1 femtosecond is used for the integration of the equations of motion and a 10⁻⁶ eV energy convergence criterion. An NVT-ensemble simulation at the temperature T = 300 K is employed. A supercell containing 6 layers of Li (100) (78 atoms) and LLZO (100) with Li termination (192 atoms) is developed. The fully optimized structures at T = 0 K have been used as initial structures to simulate AIMD calculations. The final Li/LLZO interface used for the AIMD structure is shown in Fig. 1b. The stability of the interface between Li electrode and the lowest energy LLZO surface ((100) Li terminated) has been investigated with ab initio molecular dynamics (AIMD) simulations. All production runs of AIMD trajectories have been obtained after 20 picosecond of thermal equilibration. It is worth noting that the (100) surface of LLZO has already been pointed out as the lowest energy surface in earlier studies.⁴⁵

Zoom In Zoom Out Reset image size

To establish the atomic structure/stoichiometry of LLZO GBs, Monte Carlo and Molecular dynamics simulations have been conducted using LAMMPS with soft BV potentials.^12,46–48 The GB structure has been determined through a series of calculations including zero Kelvin geometry optimization, NPT MD heating to 1200 K, NPT equilibration at target temperature, NVT equilibration at target temperature for 3 ns, and NVT Monte Carlo (MC) simulations. These calculations have provided the structure of low angle grain boundaries $\left[\displaystyle {\sum }_{3}\left.\left(112\right)\right)/\left(1-10\right)\right],$ where the tilt angle is 70.5°. It reveals that the equilibrium composition near the $\displaystyle {\sum }_{3}\left(112\right)/\left(1-10\right)$ GB in LLZO is enriched with Li. The structure of the grain boundaries has been used to calculate various properties, such as, Li ion conductivity, and mechanical properties that are needed for the continuum level modeling of dendrite growth. Approximately 1.9 M atoms have been used for mechanical property calculations.

All the computational analysis in the atomic scale were conducted in length and time scales in the range of nanometers and nanoseconds.^9,12 Continuum level analysis usually involve length and time scales in the range of hundreds of microns and thousands of seconds, for example growth of dendrite and complete short circuit of a LLZO pellet.^7,30 In order to bridge these atomic scale estimations with the predictions conducted at the continuum level, it is necessary to successfully capture the physicochemical phenomena that occurs at the intermediate length and time scales.⁴⁹ In the present context, ion transport in the GI/GB region of the LLZO SSE, lithium deposition at the electrode/LLZO interface, mechanical stress generation within Li and LLZO, have been modeled in detail. Several prior research activities characterized the computational modeling of physical phenomena at GI/GB domains as mesoscale models.^50,51 Hence, the intermediate length (sub-micron) and time (seconds) scale models developed here to understand the mechanical stress generation and charge transport processes within Li and GI/GB of LLZO, can be characterized as mesoscale models.⁴⁹

The mesoscale analysis has been conducted in a two-step fashion, where the stress induced electrochemical potential at the GB and the bulk has been evaluated using a Localized Model, and the propensity of current focusing and dendrite growth in the electrolyte have been estimated using the Generalized Model.²¹ Suppression of Li protrusions due to the plastic deformation and creep flow of metal deposits have been estimated using a modified version of the Localized Model, which were later incorporated into the Generalized Model as an effective parameter. Details of this two-step computational methodology for evaluating the possibility of current focusing and fracture initiation at the GB domain has already been described in an earlier manuscript published by the authors.²¹ For the sake of completion, the most relevant governing equations and boundary conditions will be briefly discussed below.

As is evident from the name itself, the Localized Model focuses at a very small domain at the triple junction point of lithium electrode, GB and bulk grain interior, where the domain width is approximately 50 nm.²¹ Based on the evolution of mechanical stress at the Li/LLZO interface, the magnitude of stress induced electrochemical potential is estimated, which impacts the overall reaction current density.^6,22 Due to differences in the elastic properties at the GB and bulk of LLZO, different magnitudes of stress induced electrochemical potential evolves there, which results in variation in reaction current density at the GB and bulk domains.²¹ The magnitude of stress induced electrochemical potential has been transferred from the Localized Model to the Generalized Model, where the extent of current focusing along with possibility of fracture initiation and propagation, have been investigated.

During deposition, the total amount of Li that should deposit at a particular point can be estimated from the local current density using the Faraday's Law.^52,53 Due to plastic deformation of Li, all the metal that deposits at the Li/LLZO interface, or within a cracked region of the LLZO electrolyte, does not cause deformation of the solid electrolyte. Actual deformation of LLZO is obtained from a modified version of the Localized Model, where the fractional vertical deformation of LLZO has been estimated during Li deposition at the electrode/electrolyte interface. Similarly, fractional horizontal deformation of LLZO has been estimated for metal deposits at open cracked regions. These fractional horizontal and vertical deformations of LLZO have been used in the Generalized Model while analyzing the initiation and propagation of lithium dendrites.

A larger computational domain has been adopted in the Generalized Model, where the domain width is around 500 nm, and contains multiple grain/GB microstructural features.²¹ The total number of grain/GB/Li triple points depends on the grain size adopted in the analysis. Along with the mechanical stress induced electrochemical potentials, the relative conductivity and Li concentration of the bulk, with respect to the GB, impacts the presence or absence of current focusing at the GBs.^9,12,21 To simplify the computational analysis, no lithium electrode has been analyzed in the Generalized Model. In the presence of current focusing, non-uniform deposition of lithium should occur, and that can lead to application of point loads within LLZO electrolyte. It has also been demonstrated by the authors in their earlier research work, that current focusing at the GB domain can eventually lead to initiation of crack fronts at the same place.²¹ Propensity of these crack fronts to propagate and form spanning cracks through which lithium dendrites can easily grow, will also be predicted by the Generalized Model.

To capture the ion transport within lithium electrode and LLZO electrolyte, the charge balance equation, as demonstrated below, has been solved in both the Localized as well as the Generalized Models.^21,52,53

$$\begin{eqnarray}&&\mathop{{\rm{\nabla }}}\limits^{ \rightharpoonup }\cdot \left(\kappa \mathop{{\rm{\nabla }}}\limits^{ \rightharpoonup }\phi \right)=0\end{eqnarray} \tag{ 1 }$$

where, $\mathop{{\rm{\nabla }}}\limits^{ \rightharpoonup }$ is the gradient operator, $\kappa$ indicates conductivity, and $\phi$ denotes potential. The magnitude of conductivity is substantially different in lithium metal, bulk LLZO and GB domains, which has been adopted from the lower length scale simulations. Assuming that ceramic LLZO behaves as a single ion conductor, no equation has been solved to capture the evolution of concentration gradient within the solid electrolyte.^21,52 Reaction current at the Li/LLZO interface is given by the modified version of the Butler-Volmer equation, as shown below^6,22,54 (see Appendix A and Fig. A·1 for a physical understanding):

$$\begin{eqnarray}&&\begin{array}{l}{i}_{BV}={i}_{0,ref}\cdot \exp \left(\displaystyle \frac{\left(1-{\alpha }_{a}\right){\rm{\Delta }}{\mu }_{{e}^{-}}}{RT}\right)\cdot \left[\exp \left(\displaystyle \frac{{\alpha }_{a}F{\eta }_{s}}{RT}\right)\right.\\ \left.\,\,-\,\exp \left(-\displaystyle \frac{{\alpha }_{c}F{\eta }_{s}}{RT}\right)\right]\end{array}\end{eqnarray} \tag{ 2 }$$

where, ${i}_{0,ref}$ is the reference exchange current density, ${\alpha }_{a}$ and ${\alpha }_{c}$ are the anodic and cathodic transfer coefficients, respectively, ${\rm{\Delta }}{\mu }_{{e}^{-}}$ indicates the mechanical stress induced electrochemical potential, and $F,$ $R,$ and $T$ are Faraday's constant, the universal gas constant and the absolute temperature. The surface overpotential, ${\eta }_{s},$ for this lithium-electrolyte interface is given by ${\eta }_{s}={\phi }_{s}-{\phi }_{e}+\left({\rm{\Delta }}{\mu }_{{e}^{-}}/F\right)-U,$ where ${\phi }_{s}$ and ${\phi }_{e}$ are the potentials in the electrode and electrolyte phases, and $U$ indicates the open circuit potential, which has been assumed to be zero for the reduction reaction of lithium metal. The reference exchange current density, ${i}_{0,ref},$ is given by the following expression^22,54:

$$\begin{eqnarray}&&{i}_{0,ref}=F\cdot {k}_{ref}\cdot {\left({c}_{e}\right)}^{{\alpha }_{a}}\end{eqnarray} \tag{ 3 }$$

where, ${k}_{ref}$ indicates the reference rate of reaction, and ${c}_{e}$ denotes the concentration of lithium salt within the electrolyte adjacent to the electrode. As mentioned earlier, in solid ceramic electrolytes lithium ions are part of the crystal lattice structure, and hence, the concept of salt concentration does not make any sense.^12,55 Even though concentration gradients do not evolve within the LLZO solid electrolytes, the concentration of Li ions at the Li/LLZO interface affects the exchange current density through its impact on the jumping frequency. There may be a variation in the concentration of lithium ions within the bulk LLZO and GB domains.^9,12 This difference in ion concentration can affect the reference exchange current density, which will again have an impact on the current focusing observed within LLZO electrolytes.^22,54 The appropriate equilibrium concentration of lithium ions at the bulk and GB of LLZO has been estimated using the lower length scale simulations. The stress induced electrochemical potential can be expressed in terms of the mechanical stress $\left(\mathop{\sigma }\limits^{=}\right)$ and the partial molar volume of lithium metal $\left({\bar{V}}_{Li}\right)$ as,⁶

$$\begin{eqnarray}&&{\rm{\Delta }}{\mu }_{{e}^{-}}=-\displaystyle \frac{{\bar{V}}_{Li}}{2}\left\{\hat{n}\cdot \left[\hat{n}\cdot \left({\mathop{s}\limits^{=}}^{Li}-{\mathop{s}\limits^{=}}^{SSE}\right)\right]\right\}+\displaystyle \frac{{\bar{V}}_{Li}}{2}\left\{{p}^{Li}+{p}^{SSE}\right\}\end{eqnarray} \tag{ 4 }$$

where, $\hat{n}$ is the unit normal vector at the lithium/electrolyte interface pointing into the solid state electrolyte layer, ${\bar{V}}_{Li}$ is the molar volume of metallic lithium, $\mathop{s}\limits^{=}$ is the deviatoric stress tensor and $p$ is the hydrostatic component of the stress tensor. The $Li$ and $SSE$ superscripts indicate lithium electrode and solid-state electrolyte, respectively. Please note that the compressive hydrostatic stress has been assumed to be positive here.⁶ Since the interface between Li-metal and the SSE is taken to be flat, there is no curvature-induced component. Also, in SSEs, the lithium cations are a part of the crystal lattice structure, so their partial molar volume within the SSE is taken to be zero.^18,22 It should be pointed out that the present formulation of stress induced electrochemical potential, as reported in Eq. 4, is only applicable to ceramic based inorganic single ion conducting solid state electrolytes, where only lithium ions carry the current within the bulk of the electrolyte. Since, the anions are part of the lattice structure, they do not participate in the charge carrying process. The transference number of cations in 1.0, and the transference number for anions is 0.0, within these LLZO SSEs. Following Eq. (20) in the article Monroe and Newman (JES, 2004, A880),²² negligible mobility of anions results in zero partial molar volume of cations within LLZO. In polymer electrolytes, the anions actively participate in the ion transport process. Hence, the transference number of anions in polymers is greater than zero, and results in positive magnitudes of lithium partial molar volume.²² It is evident from Eqs. 2–4 that the reaction current density at the Li/LLZO interface is affected by both lithium ion concentration and stress induced electrochemical potential.

The momentum balance equation, as shown below, has been solved at the electrolyte as well as within the lithium electrode domain to determine the mechanical equilibrium within the solids.^56,57

$$\begin{eqnarray}&&\mathop{{\rm{\nabla }}}\limits^{ \rightharpoonup }\cdot \mathop{\sigma }\limits^{=}=0\end{eqnarray} \tag{ 5 }$$

The magnitude of mechanical stress $\left(\mathop{\sigma }\limits^{=}\right)$ is directly proportional to the strain $\left(\overline{\overline{{\epsilon }}}\right),$ where the elastoplastic modulus acts as proportionality constant^.56,58 The elastic moduli of lithium, bulk LLZO, and GB domain vary substantially, and their magnitudes have been estimated using the lower length calculations. The elastoplastic properties of Li metal, such as, yield strength, hardening modulus, have been obtained from existing literature^{.5,23,27,32,59} In the Localized Model, since estimation of the magnitude of stress is important, the momentum balance equation is solved using the finite element method (FEM).⁵⁷ Also, determination of the fractional horizontal and vertical deformations of LLZO electrolyte has been conducted in modified Localized Model using FEM. Whereas, in the Generalized Model, capturing the evolution of mechanical degradation in the form of crack formation and propagation within the SSE is more important, which encouraged the usage of a lattice spring methodology (LSM) for solving the equilibrium equation.⁶⁰ Please refer to earlier publications by the authors for more detailed description of the computational techniques used for estimating the mechanical equilibrium within the solid electrode and electrolyte.²¹

In the Generalized Model, propensity of rupture within the SSE was estimated based on the evolution of strain energy $\left({\rm{\Psi }}\right)$ within each of the lattice spring elements⁶¹:

$$\begin{eqnarray}&&{\Psi }=\displaystyle \frac{1}{2}\cdot \mathop{f}\limits^{ \rightharpoonup }\cdot \mathop{{\rm{\Delta }}u}\limits^{ \rightharpoonup }\end{eqnarray} \tag{ 6 }$$

where $\mathop{f}\limits^{ \rightharpoonup }$ is the local force vector and $\mathop{{\rm{\Delta }}u}\limits^{ \rightharpoonup }$ indicates the local displacement vector. If the local force is tensile and the magnitude of strain energy exceeds the fracture threshold of that lattice element $\left({{\rm{\Psi }}}_{t}\right),$ it is considered broken and irreversibly removed from the network.⁶¹ For each lattice element, the fracture threshold value is separately sampled from a uniform distribution with mean ${{\Psi }}_{t},$ adopted from atomistic calculations, and ranging between $0.9\cdot \overline{{{\Psi }}_{t}}$ and 1.1$\cdot \overline{{{\Psi }}_{t}}.$^60–62 If the lattice element resides within the bulk grains, the fracture threshold obtained from DFT calculations have been used $\left({\left.\overline{{{\Psi }}_{t}}\right|}_{bulk}=\overline{{{\Psi }}_{t}}\right).$ Whereas, for lattice elements residing within the GB domain, fracture threshold energy have been assumed to be half of that observed in the bulk domains $\left({\left.\overline{{{\Psi }}_{t}}\right|}_{GB}=0.5\times \overline{{{\Psi }}_{t}}\right).$⁶³

As an initial part of this study we investigated the seven possible surfaces of LLZO to determine which surface is most stable. Since the LLZO surface can have different terminations, we also investigated four different terminations, Li, La, Zr and O. The calculations have been carried out with the PBE functional with a plane wave basis. A summary of the results is given in Fig. 1a for 28 possible surfaces. The lowest energy surface is found to be the LLZO (100) surface with Li termination, which has a surface energy of 0.8 J m⁻² (see Fig. 1a). Earlier studies have also revealed that the (100) surface of LLZO is the lowest energy surface.⁴⁵ In addition, another four Li terminated surfaces [(011), (101), (110), (111)] have low energies. The zirconium (Zr) terminated surfaces are generally quite high in energy. Based on the LLZO surface studies, an interfacial supercell between Li anode and LLZO has been constructed (as shown in Fig. 1b). Please note that, to be consistent with experimental measurements of partial occupancy on the Li sublattice, we generated LLZO structures where Li ions have been quasi-randomly distributed on the 24d and 96h sites, resulting in occupancies of 0.542 and 0.448, respectively.⁶⁵

To establish the atomic structure/stoichiometry of the LLZO grain boundaries (GB), Monte Carlo and molecular dynamics simulations have been conducted using LAMMPS with soft BV potentials.^12,46–48 Several experimentalists reported existence of Li metal at the GB region, and attribute that to the propagation of Li dendrites.⁷ There may be a possibility of segregation of Li ions at the GB of as synthesized LLZO pellets,^12,66 which can help in the deposition of Li ions at the LLZO GBs due to electrochemical reduction of Li.¹⁷ The structure of the grain boundaries has been used to calculate various properties such as Li ion conductivity and mechanical properties that are needed for the continuum level modeling of dendrite growth. Based on our MD and MC simulations, the lithium ion concentrations within the bulk and GB domains of LLZO are $0.0108\,{\rm{atoms}}\,{{\rm{A}}}^{-3}$ and $0.0195\,{\rm{atoms}}\,{{\rm{A}}}^{-3},$ respectively. These numbers translate into lithium ion concentrations of ${c}_{L{i}^{+},Bulk}=18,012\,{\rm{mol}}\,{{\rm{m}}}^{-3}$ and ${c}_{L{i}^{+},GB}=32,422\,{\rm{mol}}\,{{\rm{m}}}^{-3}$ within the bulk and GB regions, respectively.

Around 1.9 million atoms have been simulated in the MD technique for mechanical property calculations. Fig. 2a shows the stress-strain curves as obtained by applying uniaxial load in the computational domain. Different numbers of atoms have been considered in the MD calculations, all of which leads to almost similar stress-strain patterns. At higher magnitudes of strain (strain > 0.1), the stress values drop due to brittle fracture of the LLZO samples.⁶⁷ It is possible to estimate the Young's modulus in the bulk and GB of the LLZO solid electrolyte from the slope of the linear portion of the stress-strain curves (as shown in Fig. 2b). Figures 2b (i) and 2b (ii) depicts the stress-strain data in bulk grains obtained from two different samples and the corresponding linear fits for estimating the magnitude of Young's modulus. Similarly, Figs. 2b (iii) and 2b (iv) shows the stress-strain data obtained within the GB region along with a linear fit. By measuring the slope of these linear fits, the exact magnitude of Young's modulus has been estimated to be 141 GPa within the bulk and 120 GPa in the GB regions, which is very similar to that observed experimentally.⁶⁸ It is also possible to estimate the fracture threshold energy at the bulk of LLZO by calculating the area under the stress-strain curve shown in Fig. 2a. Based on the present atomistic calculations, the fracture threshold energy of bulk LLZO is around 5.3 J m⁻², which is approximately half of that observed in experiments.⁶⁷ Since, the formation of a grain boundary requires some excess energy, the total energy required to create fracture within a GB region is less than that needed to form a crack in the grain-interior. Hence, the fracture energy of GB can be assumed to be smaller than the bulk.⁶³ All these mechanical parameters for LLZO have been used in the mesoscale model for analyzing nucleation and propagation of lithium dendrites.

Zoom In Zoom Out Reset image size

MD simulations have been conducted to estimate the diffusion coefficient of Li within the bulk and GB of the LLZO lattice structure. The results are shown in Fig. 3a, where the magnitudes of Li diffusivity in the bulk (black squares) and GB (red squares) of LLZO are shown with respect to temperature. The corresponding black and red solid lines indicate linear fits to the obtained data, the slope of which can be used to estimate the activation energy associated with Li transport in LLZO. It is evident that with increasing temperature, the transport of Li becomes easier within the solid electrolytes. It can also be concluded that the diffusivity of Li in bulk LLZO is larger in magnitude than the diffusivity in GB. Also, the estimated diffusion coefficients correlate very well with the calculation conducted by Yu and Siegel, which are denoted by the blue (bulk LLZO) and cyan (GB) circles.¹² The conductivity of Li ions within LLZO has been estimated from the following equation: $\kappa =\left({F}^{2}{D}_{+}{c}_{+}\right)/\left(RT\right),$ where, $\kappa$ is the ionic conductivity of LLZO, $F$ indicates Faraday constant, ${D}_{+}$ is the Li diffusivity in LLZO, ${c}_{+}$ denotes Li concentration, $R$ is the universal gas constant, and $T$ is temperature.⁵² The temperature dependent conductivity of Li ions within the bulk and GB of LLZO is shown in Fig. 3b by the black and red squares, respectively. It is obvious that with decreasing temperature, the conductivity of Li ions decrease. The blue circles indicate experimentally observed conductivities of un-doped LLZO at room temperature conditions.⁶⁹ An extremely good correlation between the computational predictions and experimental observations serve as a good validation for the developed atomistic model. The exact values of LLZO conductivity at room temperature, shown in Fig. 3b, are around three orders of magnitude smaller than that of doped LLZO reported by experimentalists.^13,30,39,40 Since, no dopants are used in the atomistic calculations to stabilize the cubic phase at room temperature conditions, the predicted conductivity corresponds to that of tetragonal LLZO, and correlates well with the experimentally observed conductivities of tetragonal LLZO pellets.⁶⁹

Zoom In Zoom Out Reset image size

The magnitude of exchange current density is a very important parameter required for successfully capturing the current distribution at the Li solid electrolyte interface.^52,53 In the present research, the magnitude of exchange current density $\left({i}_{0}\right)$ at the bulk LLZO has been estimated from the charge transfer resistance $\left({R}_{ct}\right)$ data obtained from experimental literature,^39,40,53,70 according to the equation, ${i}_{0}=\left(RT/F{R}_{ct}\right).$⁵³ The exchange current density at the GB/Li interface has been obtained by scaling the bulk/LLZO exchange current density by a fraction obtained from the relative Li concentration at bulk and GB of LLZO $\left({i}_{0,GB}={i}_{0,Bulk}\cdot {\left({c}_{L{i}^{+},GB}/{c}_{L{i}^{+},Bulk}\right)}^{0.5}\right).$^22,54 Detailed estimation of the exchange current density at Li/LLZO interface is in Table II.

Table II.  Various experimental studies from where the magnitude of exchange current density has been estimated according to the relation, ${i}_{0}=\left(RT/F{R}_{ct}\right).$

Charge transfer resistance $\left({R}_{ct}\right)$ $\left[{\rm{\Omega }}.{{\rm{cm}}}^{2}\right]$	Exchange current density $\left({i}_{0}\right)$ $\left[{\rm{A}}/{{\rm{m}}}^{2}\right]$	References	Average exchange current density $\left({i}_{0}\right)$ $\left[{\rm{A}}/{{\rm{m}}}^{2}\right]$
400	0.64	40	1.52
512	0.5	70
75	3.44	39

Next, using the parameters obtained from the atomistic calculations, current focusing, Li nucleation and dendrite growth through the GB region of LLZO solid electrolytes were investigated using the mesoscale model. The charge balance (shown in Eq. 1) and mechanical equilibrium (given by Eq. 5) equations are solved as part of the mesoscale analysis. A highly nonlinear Butler-Volmer equation (given by Eq. 2) has been solved at the electrode/electrolyte interface to estimate the distribution of reaction current in bulk grain and GB regions.^22,54 The stress field at the Li/LLZO interface also affects the reaction current density through the mechanical stress induced electrochemical potential $\left({\rm{\Delta }}{\mu }_{{e}^{-}}\right)$ term (see Eqs. 2 and 4).²¹ The entire mesoscale analysis has been conducted in two different length scales: i) In the lower length scale, a very small portion around the GB region has been analyzed (domain size 50 nm—100 nm), which is characterized as the Localized Model, and ii) In the higher length scale, a larger domain size (∼500 nm) has been simulated with various grain/GB microstructure with grain size ranging from 25 nm–200 nm, which is characterized as the Generalized Model. The Localized Model provides information regarding the current focusing and fractional distortion occurring at GB region, whereas, the Generalized Model informs about the time required for nucleation and propagation of Li dendrites.

The localized version of the mesoscale model provides information regarding the distribution of current at the Li/LLZO interface by taking into account the effect of conductivity, Li concentration and mechanical stress induced change in electrochemical potential. Figure 4a shows the schematic diagram of a LLZO microstructure (adopted from Cheng et al., ACS Appl. Mater. Interfaces (2015) 7 2073–2081),⁸ with a Li metal electrode. The zoomed in view of the vicinity of a GB region is shown in Fig. 4b, which has been adopted to run the Localized Model simulations. Conductivity of ions in a particular region depends on the concentration of ions as well as the energy barrier associated with ion transport within the domain.⁷¹ It is interesting to note that even though the concentration of Li ions within the GB region is larger than the bulk, the activation energy for Li transport in GB is substantially larger than the activation energy within the bulk domain.¹² The combined effect leads to lower conductivity of Li ions within the GB domain as compared to the bulk (see Fig. 3). Due to lower conductivity of Li in the GB region,^9,12,13 the electrical driving force tends to deposit more Li in the bulk grain interiors. However, higher Li concentration at the GB increases the propensity of Li deposition at the GB region because of its impact on the exchange current density (as pointed out in Eq. 3).^12,66 Also, lower magnitudes of elastic modulus associated with GB leads to the development of smaller stress, and effectively produce smaller magnitudes of ${\rm{\Delta }}{\mu }_{{e}^{-}}$ in the GB than in bulk grain interiors $\left({\rm{\Delta }}{\mu }_{{e}^{-},GB}\lt {\rm{\Delta }}{\mu }_{{e}^{-},Bulk}\right).$^9,21 This leads to enhanced current focusing at the GB region solely due to the inhomogeneity of stress field between bulk and GB of LLZO. Figure 4c demonstrates the current distribution at the Li/LLZO interface under an applied current density of 1 A m⁻². The overall current distribution is shown by the black line, which has contributions from reduced GB conductivity, enhanced Li concentration at GB, and inhomogeneity in stress distribution at the Li/LLZO interface. The GB region has been assumed to be 10 nm thick, and resides between 20 nm and 30 nm in Fig. 4c. The current focusing observed in the GB region can be broadly divided into two components:

• i)  
  Current focusing induced by higher Li concentration at GB $\left({c}_{L{i}^{+},GB}\gt {c}_{L{i}^{+},Bulk}\right),$ which has been denoted by the blue dashed line in Fig. 4c.
• ii)  
  Current focusing induced by lower Young's modulus at GB $\left({Y}_{GB}\lt {Y}_{Bulk}\right),$ which has been denoted by the red dash-dot line in Fig. 4c.

Zoom In Zoom Out Reset image size

It is worthwhile to point out that the increase in Li ion concentration at the GB region is definitely observed far away from the Li/LLZO interface. But, near the Li/LLZO interface, there exists an uncertainty regarding the increment in Li ion concentration at the GB domain, because the presence of Li metal can impact the distribution of ions. Overall, the GB region experiences higher current focusing, which is captured by the lithium concentration (denoted as ${c}_{e}$ in Eq. 3) and stress induced electrochemical potential term (denoted as ${\rm{\Delta }}{\mu }_{{e}^{-}}$ in Eq. 2). This variation in ${\rm{\Delta }}{\mu }_{{e}^{-}}$ and ${c}_{L{i}^{+}}$ between the grain interior and GB of LLZO is transferred to the Generalized Model, where propensity for dendrite nucleation and growth is investigated.

In the present analysis, it has been assumed that Li dendrites can only nucleate by creating a crack within LLZO at the Li/electrolyte interface.^9,21 Unlike several other mechanisms, electronic conductivity through LLZO is not needed for growing Li dendrites in this particular process.^16,17,20 Since the current focusing at the GB is higher than the bulk, it is assumed that initial nucleation of Li dendrite will happen by forming a crack at the GB domain.²¹ Also, lower fracture energy of the GB region helps to initiate crack fronts.⁶³ A schematic diagram of the Li/electrolyte system used in this modified version of the Localized Model is shown in Fig. 5a. The strain energy required to rupture the LLZO-GB comes from the non-uniform mechanical deformation of the Li/LLZO interface during Li deposition. Due to the possibility of extreme plastic deformation of Li metal, it has been argued that under smaller operating currents (∼1 A m⁻²), Li should flow freely (yield strength ∼0.3 MPa) under the application of a small amount of stress, and should never accumulate substantially to cause fracture within the LLZO GBs.²⁷ Due to its viscoplastic nature, the yield strength of Li increases with increased rate of deposition, and is supposed to demonstrate larger yield strengths (∼100 MPa) while depositing under higher current densities (∼100 A m⁻²).^25,27 Also, the size of the propagating Li dendrites is extremely small (∼10 nm), and because of the Hall-Petch effect, can demonstrate higher yield strengths (∼100 MPa).^31,32 Hence, it is obvious that depending on the rate of deposition and size of the dendritic protrusion, the yield strength of Li can vary from 0.3 MPa to 100 MPa.

Zoom In Zoom Out Reset image size

Irrespective of the parameters, it is well understood that not all the Li that deposits at the interface result in deformation of the LLZO electrolyte. As shown in Fig. 5b, due to the current focusing at GB, according to Faraday's law the amount of Li deposited can be denoted as ${h}_{dep}.$ However, during mechanical equilibration of Li and LLZO interface, and plastic deformation of Li metal, the equilibrated height of the deposit will be smaller (as shown in Figs. 5c and 5d), and can be denoted as ${h}_{fin}$ (please note that ${h}_{fin}\lt {h}_{dep}$). The difference between ${h}_{fin}$ and ${h}_{dep}$ depends on the yield strength of Li (see Figs. 5c and 5d). Let's define the ratio between ${h}_{fin}$ and ${h}_{dep}$ as the "fractional vertical distortion" of the LLZO solid electrolyte. It is also expected that both the yield strength of Li, as well as the elastic modulus of the GB region $\left({Y}_{GB}\right),$ will impact the magnitude of fractional vertical distortion.^9,21 Accordingly, Fig. 5e demonstrates the variation of this fraction of vertical displacement experienced by LLZO with respect to Li yield strength, which has been plotted for three different magnitudes of GB elastic modulus $\left({Y}_{GB}/{Y}_{Bulk}\right).$ It is evident that, with increasing yield strength of Li, the depositing metal becomes stronger, and can deform the electrolyte to a larger extent, effectively increasing the magnitude of fractional vertical distortion. Similarly, decreasing the elastic modulus of the GB region makes the electrolyte softer so that it can deform easily, which effectively helps to increase the magnitude of fractional vertical distortion. It is worth mentioning that the Young's modulus of Li is smaller than the Young's modulus of LLZO by an order of magnitude. Hence, even at very high yield strength of Li, even though there is no plastic deformation, the fractional vertical distortion is much smaller than unity (around 0.1 in Fig. 5e) due to the smaller magnitude of elastic modulus. Please note that, estimation of the fractional vertical distortion has been conducted prior to any fracture initiation within LLZO.

Incorporating this fractional vertical distortion into the Generalized Model, it is possible to estimate whether the depositing Li will be completely squeezed out of the GB domain, or will it accumulate and lead to fracture of the LLZO electrolyte. If the depositing Li keeps accumulating at the GB, the local strain energy density can exceed its fracture threshold, and the LLZO electrolyte can fail. Figure 6 shows the evolution of normalized strain energy density within a LLZO microstructure with grain size 197 nm and applied current density of 1 A m⁻². In Fig. 6a, the yield strength of Li is 0.7 MPa, whereas, in Fig. 6b, the yield strength of Li is 100 MPa. It is evident that even for low magnitudes of Li yield strength (∼0.7 MPa), it is possible to accumulate substantial amount of Li at the GB domain, such that the local strain energy density exceeds the fracture energy threshold, and crack fronts start to evolve within LLZO. However, as shown in Fig. 6a, for low yield strength Li, around 0.7 MPa, the time taken to initiate a crack is large, approximately 4.8 h. Whereas, Fig. 6b shows that when the Li yield strength is much larger, around 100 MPa, lithium deposition for a smaller time, only 1.6 h, can lead to fracturing of LLZO. It is also worth noting that under lower yield strengths, the excess strain energy is very much concentrated (see Fig. 6a), whereas, for high yield strength Li, excess strain energies are more distributed in the LLZO SSE (see Fig. 6b). This clearly indicates that for higher yield strengths, Li dendrites can be more branched, and under lower yield limit, the dendrites will be more unidirectional.

Zoom In Zoom Out Reset image size

Next, using the Generalized Model, the dependence of the current focusing at GBs and lithium dendrite initiation time on the LLZO grain size is investigated. The current focusing at the GB region (induced by variation in Li concentration and stress induced electrochemical potential) and magnitudes of fractional horizontal distortion, has been extracted from the Localized Model simulations, and incorporated into the Generalized Model. The domain size used in the Generalized Model is around 500 nm, and the applied current density has been maintained at 1 A m⁻². The GB elastic modulus is around 85% of the bulk $\left(\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85\right).$ As the LLZO grain size decreases from 200 nm to 20 nm, the variation in current focusing at the GB region is shown in Fig. 7a. For grain sizes larger than 100 nm, changing the grain size has no impact on current focusing. Such sizes have been characterized as the "bulk zone." For grain sizes smaller than 100 nm, decreasing grain size helps to minimize the extent of current focusing observed at the GB region. These sizes can be denoted as the "GB affected zone." For grains with size smaller than 100 nm, the GBs are located close enough to each other, and one GB can see its neighbor from an electrochemical perspective. As the grain size decreases, the number of grain boundaries increases within a unit length of the Li/LLZO interface (2D scenario). With more GBs, the total current at the Li/LLZO interface has multiple places to get distributed, and current at each of the GBs are effectively minimized. This minimization in current focusing with decreasing grain size has already been discussed in an earlier manuscript by the authors.²¹ For grains with size larger than 100 nm, the GBs are located so far away from each other, that from an electrochemical standpoint they cannot see each other, and changing the grain size does not alter the current focusing by any significant amount.

Zoom In Zoom Out Reset image size

Assuming the magnitude of Li yield strength to be around 0.7 MPa, the time required to initiate fracture within LLZO microstructures with various grain sizes have been investigated, which is shown in Fig. 7b. It has been assumed that initiation of Li dendrites is only possible through formation of fracture at the LLZO grain boundaries. Hence, fracture initiation time can be considered as a good approximation for the initiation of dendrites. In the "bulk zone," as the grain size decreases from 200 nm to 100 nm, the time required for nucleating Li dendrites at the Li/LLZO interface decreases from 4 h to 3.7 h. However, in the "GB affected zone" as the grain size is decreased even more, to 25 nm or so, time for dendrite nucleation increases to 4.2 h. The corresponding thickness of deposited Li is denoted by the axis on the right side. There exists an intermediate grain size (∼100 nm) where fracture initiation occurs the earliest. For grains with size smaller and larger than this optimum magnitude, it takes longer time to initiate fracture within the LLZO electrolyte.

Grains with size smaller than 100 nm falls under the "GB affected zone," where decreasing grain size helps to minimize the current focusing at the GB regions (as shown in Fig. 7a). Nonuniform deposition of Li at the GB region provides the driving force for crack initiation within the LLZO electrolytes. The amount of Li deposited at the GB region dictates how easily cracks can form. In the "GB affected zone," decreasing grain size helps to minimize the amount of Li deposited at the GB region. Hence, decreasing the grain size should increase the time required to nucleate dendrites within the LLZO SSE. Next, the reason behind increase in dendrite nucleation time with increasing grain size, for grains larger than 100 nm, will be discussed.

Two different LLZO microstructures have been demonstrated in Figs. 7c and 7d with average grain sizes of 197 nm and 120 nm, respectively. Both of them reside within the "bulk zone," because their grain size is larger than 100 nm. Enhanced deposition of Li at the GB region leads to application of mechanical force at the GBs. This eventually causes rupture in LLZO SSEs and dendrites can nucleate there. For LLZO grains with size larger than 100 nm, the extent of current focusing at the GB is independent of the grain size, which is evident from Fig. 7a. Let's assume that Li deposition induced point load acting at the GB of Fig. 7c is given by ${{\rm{\Sigma }}}_{1}$ (grain size ∼197 nm), and the same for Fig. 7d is given by ${{\rm{\Sigma }}}_{2}$ (grain size ∼120 nm). Since the current focusing for grains with size larger than 100 nm is independent of the grain size, one can conclude that the Li deposition induced point load that acts at the GB region is also independent of the grain size, or ${{\rm{\Sigma }}}_{1}\approx {{\rm{\Sigma }}}_{2}\approx {\rm{\Sigma }}.$ The total load acting at the Li/LLZO interface is given by the product of the number of grain boundaries and the load at each of the GB region $\left({\rm{\Sigma }}\right).$ In the "bulk zone," where grain size is larger than 100 nm, LLZO microstructures with smaller grain size contain more GB regions. For example, in Fig. 7c the average grain size is 197 nm, and there exist only two GBs at the Li/LLZO interface. Whereas, in Fig. 7d, with grain size 120 nm, there are three GBs at the Li/LLZO interface. Hence, the total load experienced by the LLZO microstructure with grain size of 120 nm is around $3{\rm{\Sigma }},$ whereas, the total load experienced by the LLZO microstructure with grain size of 197 nm is around $2{\rm{\Sigma }}.$ As a result, the total load experienced by the LLZO microstructure with grain size 197 nm is less than that experienced by the LLZO microstructure with grain size 120 nm. Time required to nucleate cracks should be proportional to the total amount of load acting at the Li/LLZO interface. Hence, the total time required to initiate cracks in LLZO with larger grains is higher than that required for LLZO microstructures with smaller grains.

Similar results have been reported by Sharafi et al., where they observed an increase in critical current density by increasing the grain size.³⁰ Even though the experimental study is not same as that being reported here, delayed dendrite nucleation with increasing LLZO grain size, is consistent in both the studies. It is worthwhile to note that, the overall improvement in dendrite nucleation time obtained by increasing the grain size is not significant. It is evident from Fig. 7b, as the LLZO grain size increases from 100 nm to 200 nm, the time required for initiation of Li dendrites increases from 3.7 h to 4 h, a mere improvement of 20 min only.

Variation in dendrite initiation time with increasing Li yield strength is reported in Fig. 8a for $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85$ and average LLZO grain size of 197 nm. Since, Li metal with higher yield strength behaves as a stronger material, it can cause a larger deformation within LLZO, and that can lead to early fracture within the solid electrolyte. Hence, dendrite nucleation time decreases substantially as the lithium yield strength is increased from 0.4 MPa to 100 MPa. The right side of the y-axis indicates that for high yield strength of Li (∼100 MPa), dendrites should start to grow after deposition of only 200 nm thick Li film. On the other hand, for Li yield strength of 0.4 MPa, around 2.3 μm thick Li film can be deposited prior to nucleation of dendrites. It is also interesting to note that even for very low Li yield strengths (∼0.4 MPa), dendrite nucleation can occur in LLZO under applied current densities of 1A m⁻². To understand the impact of the GB elastic modulus on the dendrite initiation time, it has been plotted in Fig. 8b with respect to the Li yield strength in a log-log scale. The GB elastic modulus has been changed between $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85$ to $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.1.$⁹ Grain size of the LLZO microstructure has been kept constant at 197 nm. As the magnitude of ${Y}_{GB}$ decreases, the LLZO microstructure becomes softer, and it becomes easier to deform and rupture under Li deposition induced load. Hence, the time required for dendrite nucleation with ${Y}_{GB}\sim 0.1{Y}_{Bulk}$ is much smaller than that needed for nucleating dendrites with ${Y}_{GB}\sim 0.85{Y}_{Bulk}.$ By comparing the dendrite nucleation time in Figs. 7b and 8b, it is evident that the impact of lithium yield strength and GB elastic modulus is more significant than grain size in determining the structural integrity of LLZO electrolytes. Also note that, in order to keep the analysis simple and tractable using the present computational technique, Li deposition in a particular direction has been assumed in all these simulations. If deposition in the reverse direction is also taken into account, it should lead to more complex phenomena, such as, nonuniform dissolution, formation of dead Li, etc, which is out of the scope of the developed computational methodology.

Zoom In Zoom Out Reset image size

After nucleation of Li at the GB region, dendrites grow to short the cell. It is well accepted that growth of Li dendrites occurs mostly through the GB regions.^7,8,30 Growth of Li dendrite through the GB region should involve deposition of Li at the fractured GB, and making the crack front propagate, until the growing dendrite shorts the cell. Under elastic deformation of Li, it has been demonstrated by researchers that Li dendrite within a cracked region prefers to propagate by increasing the crack length, rather than getting extruded from the crack opening.²⁹ However, it has been argued by other researchers that under externally applied loads due to plastic deformation of Li, metal deposited at the GB cracks should flow out without causing any accumulation or further rupture of the solid electrolyte.²⁷ In order to reconcile the extent of deformation experienced by LLZO due to Li deposition in a cracked region, fractional horizontal distortion of the electrolyte has been analyzed for different yield strengths of Li metal.

Figure 9a schematically demonstrates an updated version of the Localized Model where the bulk and GB region of an electrolyte is shown by blue and black domains, respectively. In the same figure, a thin crack has been assumed in the middle of the GB region, which is denoted by white color. The total amount of Li deposit at the cracked region, as predicted by the Faraday's Law, is shown in Fig. 9b (Li denoted by orange color). Prior to mechanical stress equilibrium, the thickness of the Li deposit is denoted as ${\delta }_{dep}$ (see Fig. 9b). After mechanical equilibrium, because of the very high modulus of LLZO and plastic deformation of Li, the deposited metal experiences a decrease in thickness, and some amount of Li is squeezed out of the crack region. This is schematically shown in Figs. 9c and 9d, for high and low yield strength of Li metal, respectively. The final thickness of the Li deposit is denoted as ${\delta }_{fin},$ where ${\delta }_{fin}\lt {\delta }_{dep}.$ The ratio between ${\delta }_{fin}$ and ${\delta }_{dep}$ is denoted as the fractional horizontal distortion, and this parameter indicates the extent of horizontal deformation experienced by LLZO during Li deposition at a cracked region. Figure 9e demonstrates how the fractional horizontal distortion changes with Li yield strength. Li metal with higher yield limit behaves as a stronger material and can deform the LLZO electrolyte to a larger extent, which results in higher magnitudes of fractional horizontal distortion. As the elastic modulus of the LLZO GB domain decreases, deformation of the solid electrolyte becomes easier, resulting in larger magnitudes of fractional horizontal distortion (also demonstrated in Fig. 9e). The magenta line with triangular symbols in Fig. 9e is for $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 1,$ which indicates the Li deposition induced fractional horizontal distortion experienced by crack openings that exist in the bulk grain interior. It is worth mentioning that for extremely high yield strength of Li, the metal should deform elastically, and fractional horizontal distortion should saturate to a certain magnitude. This trend is supposed to be similar to that observed in Fig. 5e, observed in very high yield strength Li. However, the magnitude of fractional vertical distortion in Fig. 5e is much smaller than the magnitude of fractional horizontal distortion observed in Fig. 9e. This difference can be attributed to the type of deformation reported in Fig. 5e, where the LLZO did not have any cracks, whereas, in Fig. 9e, the GB of LLZO already contains a crack. Presence of a pre-existing crack front makes the LLZO GB behave like a much softer material. Please note that, during this analysis, no propagation of the existing crack front at LLZO GB has been taken into consideration.

Zoom In Zoom Out Reset image size

Transferring the magnitude of fractional horizontal distortion from the Localized Model into the Generalized one, crack propagation, and subsequently Li dendrite growth, through LLZO electrolyte has been analyzed. Figures 10a and 10b demonstrate the crack path, and hence Li dendrite, through the LLZO electrolyte for Li yield strength of 0.7 MPa and 100 MPa, respectively. The grain size of the LLZO microstructure is 88 nm, GB elastic modulus is 85% of bulk $\left(\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85\right),$ and dendrite propagation has been demonstrated under applied current density of 1 A m⁻². It is evident that the majority of the cracks reside within the GB domain, which can be attributed to the smaller fracture energy density of GB considered in the present study. However, growth of Li dendrite through GB is consistent with other experimental observations, where the cross section of shorted LLZO reveals propagation of Li dendrites only within the GB region.^7–9,17 It is also clear from Figures 10a and 10b that Li nucleates at almost all the GBs, but only a few of them propagate due to the statistical nature of the fracture energy threshold. Under both low and high yield strengths of Li, the dendrites take almost the same pathway, which can be attributed to the distribution of grains and their orientation within the LLZO microstructure. However, at high yield strengths of Li, as shown in Fig. 10b, the unsuccessful dendrites (see the GBs located around 150 nm, 250 nm and 350 nm along the x-direction) within the GB region tend to grow larger. Also, note that under low yield strength of Li, shown in Fig. 10a, a particular crack front, around 400 nm along the x-axis, tries to penetrate through the bulk grains, but fails to propagate further due to larger load requirements for growing cracks within the grains. Such attempts of propagation through the bulk grain interior have not been observed at higher Li yield strengths (see Fig. 10b). Figure 10c shows the fraction of broken elements that reside within the GB (denoted by red circles) and bulk (depicted by black squares) of LLZO electrolyte. Dendrite propagation through various LLZO microstructures with grain size ranging between 25 nm to 200 nm has been analyzed here. It is evident that, irrespective of grain size majority of the Li dendrite reside within the GB region because of its lower fracture energy density. For grain sizes larger than 100 nm, enhanced cracking within the bulk grains have been observed, which can be attributed to GBs residing extremely far away from each other and minor cracking within the bulk grains. These dendrites through the grain-interior eventually fail to propagate, as shown in Fig. 10a around 400 nm along x-axis.

Zoom In Zoom Out Reset image size

The dynamics of dendrite propagation through the LLZO electrolyte have been analyzed in Fig. 11. The increase in height of Li dendrite within LLZO with respect to time is shown in Fig. 11a (indicated by the solid black line) during Li deposition at 1 A m⁻². A LLZO microstructure with grain size around 197 nm is used here, along with Li yield strength of 0.7 MPa and $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85.$ Nucleation of Li dendrites takes long time, around 4.5 h, but it propagates very fast, and can effectively short the cell within 60 min. The red squares denote instantaneous velocity of dendrite growth, which can be as large as 2 nm s⁻¹ at certain times. However, the average velocity of dendrite growth is much smaller. This discontinuous growth indicates that, propagation of Li dendrites through LLZO is not a gradual process, but occurs through bursts of fracture events within the brittle solid electrolyte. During this propagation of Li dendrite, variation in electrochemically active surface area and potential drop across the solid electrolyte has been shown in Fig. 11b by the magenta and blue lines, respectively. As the Li dendrite propagates, the electrochemically active surface area increases, which helps to reduce the activation overpotential, and the total potential drop across the LLZO electrolyte decreases. In Fig. 11b, the sudden decrease in electrochemically active surface area (magenta line) around 5 h can be attributed to isolation of electrolyte domain due to fracture though the GB. This kind of decrease in potential drop with propagating dendrite has also been observed by experimentalists during Li deposition with point electrodes.²⁹

Zoom In Zoom Out Reset image size

The average velocity of growing Li dendrites have been estimated and plotted with respect to grain size in Fig. 11c. The Young's modulus in the GB has been assumed to be 85% of the bulk, whereas Li yield strength is 0.7 MPa. The applied current density is maintained at 1 A m⁻². The average velocity of dendrite growth is defined as:

$$\begin{eqnarray}&&{v}_{Li}=\displaystyle \frac{{\rm{Final}}\,{\rm{height}}\,{\rm{of}}\,{\rm{Li}}\,{\rm{dendrite}}}{{\rm{Time}}\,{\rm{of}}\,{\rm{short}}\,{\rm{circuit}}-{\rm{Dendrite}}\,{\rm{initiation}}\,{\rm{itme}}}\end{eqnarray} \tag{ 7 }$$

This definition makes the average velocity independent of the effect of variation in dendrite initiation time (which is similar to fracture initiation time) with different grain sizes (as shown in Fig. 8a). It is evident from Fig. 11c that the average dendrite growth velocity decreases with decreasing grain size. Growth velocity observed in LLZO microstructures with 50 nm grains is almost half of that observed in LLZO with 200 nm grains. This decrease can be attributed to an enhanced amount of GBs and a more tortuous pathway taken by the propagating dendrites inside LLZO microstructures with smaller grain sizes. The errorbars are obtained by averaging over five samples for each of the data points. However, the large standard deviation observed at smaller grain sizes indicates that the decrease in dendrite growth velocity for smaller grains is not significant. Variation in dendrite growth velocity with Li yield strength is shown in Fig. 11d for $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85.$ Li metal with a higher yield limit behaves as a stronger material, and can result in much faster growth of the dendrite through LLZO. The error-bars at each value of Li yield strength have been obtained by averaging over all the grain sizes. As shown in Fig. 11c, changing the grain size can alter the dendrite growth velocity from 0.1 nms⁻¹ to 0.25 nms⁻¹, whereas, according to Fig. 11d, changing yield strength can vary the dendrite growth velocity from 0.1 nms⁻¹ to 0.4 nms⁻¹. Hence, Li yield strength has a slightly more prominent effect on determining the dendrite growth velocity than grain size.

It has been observed in various experiments that increasing the current density leads to early nucleation of dendrites and eventually a short circuit of the cell.^8,30 Since, GBs have lower ionic conductivity than the bulk, at higher current densities more Li should pass through the grain-interiors than the GB regions. This cannot explain the early appearance of a Li dendrite through the GB under higher current densities. In the present study, the authors argue that being a viscoplastic solid; Li demonstrates higher yield strength during operation at higher current densities. During Li deposition at 1 A m⁻² and 10 A m⁻², the Li yield strength can be assumed to be around 0.7 MPa and 1 MPa, respectively.²⁷ The Young's modulus of GB has been assumed to be 85% of the bulk domain. Taking into consideration these parameters, the predicted height of Li dendrite within an LLZO microstructure of grain size 197 nm is shown in Fig. 12a. Li deposition at two different current densities are considered here, 1 A m⁻² and 10 A m⁻², which are denoted by the solid and dash-dot lines, respectively. Irrespective of applied current density, dendrite growth is a much faster process than dendrite initiation. It is also evident from Fig. 12a that, given enough time for lithium deposition, it is possible to grow dendrites through LLZO even at low current densities of 1 A m⁻². Hence, there should not exist any critical current density for LLZO SSEs, and whether short circuit should occur or not, depends on a combination of applied current density and the total amount of Coulombs passed through the LLZO SSE.

Zoom In Zoom Out Reset image size

The average time for short circuit is shown in Fig. 12b by the black squares. The time for short circuit can be considered as the summation of dendrite initiation time and dendrite propagation time. In these simulations $\left({Y}_{GB}/{Y}_{Bulk}\right)\sim 0.85$ has been maintained. The errorbars are obtained by averaging over all the different grain sizes that ranges approximately between 25 nm to 200 nm. Due to viscoplasticity of Li, the yield strength of Li has been assumed to vary according to the applied current density.^25,27,28 It is evident that at higher current densities, it takes smaller time to nucleate and propagate Li dendrites, which is consistent with various experimental observations.^8,30,70 The solid blue line indicates a stability limit, below which safe deposition of Li is possible, whereas, operating conditions that reside above the blue line should lead to dendrite growth and short circuit of the cell. In Fig. 12b, the numbers in red parenthesis gives the total thickness of Li deposited prior to shorting. Note that, before initiation of dendrites, more lithium can be passed during operation at higher current densities. This trend is counterintuitive, and can be attributed to the lower conductivity of GB than bulk LLZO $\left({\kappa }_{GB}\lt {\kappa }_{Bulk}\right).$ During operation at higher current densities, the transport limitation become dominant, and the magnitude of current focusing at the GB region starts to diminish (similar trends were also reported by the authors in an earlier publication²¹). Hence, at higher current densities, total amount of Li deposition at GB with respect to the entire electrode is relatively less than that observed at lower current densities $\left({\left.\left({i}_{GB}/{I}_{app}\right)\right|}_{\begin{array}{c}high\\ current\end{array}}\lt {\left.\left({i}_{GB}/{I}_{app}\right)\right|}_{\begin{array}{c}low\\ current\end{array}}\right),$ which leads to more Li deposition at high current before short circuit. Along the same direction, at higher currents, the total amount of Li deposited is more than that at lower currents $\left({\left.{i}_{GB}\right|}_{\begin{array}{c}high\\ current\end{array}}\gt {\left.{i}_{GB}\right|}_{\begin{array}{c}low\\ current\end{array}}\right),$ which results in earlier short circuit of the cells under high current operation, as reported in Figs. 12a and 12b. This analysis clearly indicates that only the short circuit time is not a definitive descriptor of the interfacial stability between Li and LLZO. The amount of charge passed must also be reported to successfully determine the stability of a solid electrolyte against Li metal anode. These results also conclude that there does not exist any critical current density for LLZO electrolytes. Shorting of a cell due to lithium dendrite growth depends on a combination of applied current density and the Coulombs passed through the electrolyte.