An Analysis of Solid-State Electrodeposition-Induced Metal Plastic Flow and Predictions of Stress States in Solid Ionic Conductor Defects

The use of inorganic solid-ionic conductors with a metal electrode, has been proposed as a way to increase energy density, decrease capacity loss and prevent failure from metal propagation. Current observations of Li-metal electrodes causing cell shorting in solid-state systems have been identi ﬁ ed as main obstacles limiting the development of this technology. However, many aspects of the involved phenomenon have not been fully addressed theoretically. In this work, we derive a mathematical model of electrodeposition-induced plastic ﬂ ow in metal/inorganic solid-conductor systems. We use a semi-analytical solution to derive pressure increase expressions at metal protrusions and assess the possibility of fracture. The results give ﬂ ow solutions analogous to laminar channel ﬂ ow. The solutions also show how taking into account a boundary traction potential from built up pressure, leads to ionic redistribution and effectively screens isolated ﬂ aws, making local current focusing an incomplete explanation for observed electrolyte fracture. We show that the boundary traction potential sets a maximum value for the pressure increase that can occur from deposition at an isolated ﬂ aw. We derive conditions under which fracture can occur, and quantify the role of ionic conductivity and electrolyte fracture toughness in extending safe operating regimes of solid-state electrolytes with metal electrodes.

The use of solid inorganic ion conductors paired with a reactive metal electrode, such as Li or Na, is being increasingly studied as a promising path to increase energy density of current Li-ion battery technology. Substantial progress has been done in finding chemically stable and highly conductive materials. 1,2 However successfully cycling against a Li/Na metal electrode still remains an elusive challenge. Understanding electrodeposition from a solid electrolyte to a reactive metal is a vital step toward resolving the underlying problems observed from cycling reactive metal solid state batteries (MSSB).
Current research has suggested that conductor-metal interface contact, parasitic interface reactions, deposited metal nucleation and metal propagation through the solid conductor are the main challenges limiting current development in MSSBs. 3 These unresolved challenges stem from a limited understanding of solid state electrodeposition, where electrochemical and mechanical behavior can be strongly coupled.
Recent experiments have shown remarkable behavior of deposited Li propagating through solid ionic conductors and causing cell shorting. [4][5][6][7][8] Ex situ SEM and EDS imaging have shown apparent Li deposits at grain boundaries of LLZO after single charge cell shorting. 4 Further studies using a gold tip electrode on various ionic conductors (LPS, LLZO, LS-PS, accounting for crystalline, polycrystalline and amorphous phases), where used to propose the role of critical size defects in the ionic conductor in Li propagation. 5 The most recent extension of this work shows that current inhomogeneities caused at edges of electrodes can serve as Li penetration initiation sites in single crystals of LLZO. 6 Based on experimental observations, a few mechanisms of Li nucleation and propagation have been proposed. Based on X-ray tomography the degree of pore connectivity in ionic conductors has been correlated with the critical current, suggesting that Li propagation tends to follow pore networks. 7 Fracture of the conductor from Li propagation has also been proposed as a propagation mechanism. The mechanism is described by increased current densities at sharp defects and edges leading to high deposition rates which cause large enough stress concentrations to drive fracture. 5,6 Based on DFT computed surface energy and densities of states, Li reduction at grain interfaces without direct contact with an electrode has also been proposed. 9 This mechanism allows lithium deposits in the bulk to eventually connect to form a short circuit without directly driving fracture of the conductor.
Theoretical and mathematical modeling of Lithium deposition and growth has been mostly limited to liquid and polymer electrolyte systems. [10][11][12][13][14][15] Results from liquid and polymer electrolyte modeling suggest that a stiff enough separator will block lithium propagation and reduce the risk of cell failure. These results however, cannot be directly translated to ionic conductor systems since the mode of ion transport and deposition are fundamentally different. The popular yet misunderstood criterion regarding electrolyte and metal shear moduli, derived from modeling of surface stability in polymer electrolyte systems, 14,15 does not apply to Lithium propagation in ceramic ionic conductors due to the very different interface morphology, the requirement to consider electrode plasticity, and the fundamentally different mode of ionic conduction. Recently, extensions of electrochemical and mechanical effects in metal deposition from ceramic solid electrolytes have been developed, [16][17][18] however a specific analysis of current focusing and potential mechanical failure in metal filled conductor flaws has not been directly addressed.
Decades old work on failure in liquid Na β-Alumina attributing the resulting high pressures to internal viscous flow of the liquid Na is remarkably relevant for the current problem in MSSBs. [19][20][21][22][23][24] However, justification for treating the molten Na flow as fully developed channel flow was never carefully considered, which led to incorrect use of the governing equations. 22 This type of analysis needs to be strongly reconsidered and appropriately developed before using it to rationalize or predict failure in solid electrolyte with solid metal electrode systems.
In this work, we present a continuum scale model for electrodeposition induced plastic flow in reactive metal electrode and solid ionic conductor systems, similar to what was done in work on molten Na and β-alumina. Subsequently we derive semi-analytic results to predict stress and hydrostatic pressure concentrations to assess the possibility of ionic conductor fracture. We show that when the effects of a boundary traction on deposition are neglected, and no other ionic redistribution methods are considered, the metal flow is exactly analogous to laminar and fully developed channel flow; and the pressure increase is spatially linear as given by Armstrong. 24 However, using this analysis we find that due to the high pressure predicted there is no size threshold for a flaw to become critical, rather, any flaw is a critical flaw. In contrast, we find that in the limit of fast electrode kinetics, the effect of a boundary traction on the metal ion potential, limits the maximum value of the pressure increase at an isolated flaw. The value of this saturation pressure allows us to determine thresholds for flaw size, and operating conditions for which flaws can become critical and fracture is initiated. We show that electrolyte fracture toughness and ionic conductivity are key materials properties that control the extent of sub-critical regions.

Theory of Electrodeposition-Induced Plastic Flow
In the present theoretical framework, we consider electrodeposition in pre-existing flaws in solid ionic conductors, and the mechanical response in the surrounding area of the bulk conductor and the deposited metal. Figure 1 shows a rendition of the theoretical domain: an idealized planar crack filled with the electrode metal. It has been previously proposed that metal deposition will initially fill voids and cracks based on a lower reduction overpotential at these sites. 5 Equivalently, we assume the flaw is filled by the electrode metal, and proceed to derive the mechanical response from further metal deposition.
The current density in the neighborhood of a flaw (see inset of Fig. 1) is driven by the local ion flux, which is locally homogeneous within this region but away from the metal/conductor boundary. The locally homogeneous current density near the flaw may be different than the driving current density of the cell due to a lack of perfect contact between the conductor and the electrode or contact loss upon cycling. 25 Local electrodeposition is initially enhanced at high curvature defects, such as crack tips. The local enhancement of electrodeposition is due to ion flux, driven by the electric field which concentrates near sharp conductor tips. 26 The higher deposition rates lead to stress building up at these sites triggering some stress release mechanisms. Possible stress release mechanisms include metal plastic flow, conductor fracture and ion redistribution. These stress release mechanisms do not necessarily occur independently but can occur simultaneously, in sequence or in any combination thereof. However, we contend that some degree of plastic flow must occur to explain observed critical currents, 4,5 below which no failure is observed. Barring some other rate dependent mechanism overlooked in this work, without metal flow, MSSB failure would be a function only of total charge passed, contrary to rate dependent observations. Furthermore, any spatial inhomogeneities in deposition can lead to high enough deviatoric strain where plastic deformation is inevitable.
We show that for isolated flaws such as the geometry shown in Fig. 1, flux concentration at the tip monotonically decreases along the sides of the protrusion and initially leads to plastic deformation and outward metal flow. The flow is driven by a pressure gradient along the length of the protrusion. The pressure increase couples to the ion flux in the conductor as a diffusion potential boundary condition, which in turn causes ion flux redistribution in the conductor. We find that the levels of ion redistribution predicted can homogenize ion flux and effectively screen sharp flaws in the ionic conductor. The extent to which ion redistribution homogenizes the diffusion potential in the ionic conductor is formalized by introducing an electro-plastic screening length that depends on the ionic conductivity of the electrolyte and the plastic behavior of the metal.
For screened defects, the pressure increase is predicted to be equal to the saturation pressure, which is determined by the electrolyte conductivity and driving current density. Furthermore, a locally homogeneous diffusion potential causes negligible deviatoric stress, such that plastic deformation within the flaw becomes irrelevant and the possibility of fracture can be determined from the built up saturation pressure alone.
Governing equations.-To quantify the phenomenon previously described, we use continuum ionic diffusive transport, linearelasticicty, viscous plastic flow, and linear-elastic fracture mechanics. We first list the general governing equations and describe suitable approximations. Subsequently we elaborate on the boundary conditions where key coupling between the governing equations occur.
Ionic diffusion in ceramic single-ion conductors under stress can be described using a generalized diffusion potential, Where μ is the ion chemical potential, z + is the valence of the mobile ions. ϕ is the electric potential. β and Ω v are the relative volume change of the host conductor lattice during ion diffusion and the molar volume of the conductor respectively. 27 e s and  e are the elastic stress and strain in the conductor.
The ionic flux in the conductors is proportional to the gradient of the generalized diffusion potential, M corresponds to the mobility of the ionic species. c is the ion concentration field within the conductor. T is the absolute temperature and R is the universal gas constant. shows where the considered domain fits within the total metal electrode and conductor system. The current density within the domain (which due to imperfect interface contact may not be the same as the total cell current density) drives electrodeposition. In addition the the operating stress due to the cell' s stack pressure, local stress is generated by electrodeposition, which triggers metal flow, ionic redistribution and possible fracture as local stress release mechanisms. Accordingly, key materials properties involved are the metal's pseudo-viscosity η, the conductor's ionic conductivity σ + , and the conductor's fracture toughness K IC .
At steady state the governing equation for ionic transport takes the form of a conservation of ions as follows, For single ion conductors with negligible electronic conductivity, the ionic concentration gradient must be zero everywhere within the conductor since the divergence of charge must be zero in order to satisfy local electroneutrality. Furthermore, the gradient of the elastic energy term, which accounts for stress effects on migration energy, is assumed to be negligible compared to the potential gradient term near high curvature flaws.
To an appropriate first order approximation, the mechanical response in the ionic conductor can be described with linear elasticity. For an isotropic homogeneous media the relationship between stress e s and  e is, Where E is the elastic modulus, ν is Poisson's ratio, and  is used for the identity tensor. Similarly, the mechanical governing equation for the conductor is given by a momentum conservation equation, Plastic flow of relatively soft metal electrodes such as Li and Na is treated using viscous plastic flow with negligible elastic strain (i.e. incompressible). The corresponding constitutive stress p s relation for the metal electrode is, Where P is the hydrostatic pressure.  p  is the plastic strain rate, and η is the metal pseudo-viscosity.
The pseudo-viscosity η is obtained following the Perzyna formulation of strain rate for associated flow, 28 using Von-Mises type flow. The pseudo-viscosity can be expressed in terms of the yield stress σ Y of the deforming metal, a plasticity relaxation rate γ, and the effective plastic strain rate  , which is compactly expressed using the entrywise product,    2 : The hardening exponent n quantifies the degree of work-hardening during deformation. For most metals 0.1 < n < 0.5.
For a perfectly plastic material with no hardening g  ¥, the pseudo-viscosity simplifies to, The plastic flow governing equation is also a steady state momentum conservation equation, Finally, fracture occurs when the internal loading reaches the critical stress that can be supported by the solid medium, in this case the solid electrolyte bulk. The critical stress is a material property that is usually measured experimentally in terms of a fracture toughness K IC . For linear elastic fracture, where plasticity contributes minimally to increased toughness (i.e. via crack tip blunting), the material toughness measured experimentally contains effects of fracture mode and crack geometry (stress concentration), in addition to surface and elastic energy terms. 29 For the case of solid-state electrolyte failure, the fracture mode from an internally pressurized crack is mode I fracture. Fracture occurs when the stress intensity factor (SIF) at a crack is equal to or greater than the material mode I fracture toughness. Possible fracture can be evaluated by computing the SIF for internally pressurized cracks using the pressure values obtained from solving the electrodeposition/plastic flow model. An expression for the SIF for internally pressurized cracks is, 30 Here L is the total length of the crack. l represents the spatial dimension along the length of the crack. P is the pressure inside the crack. And f(·) is a numerically fitted weight function, derived and reported by Hartranft. 30 For standard operation in confined cells, the effective SIF must include the corresponding cell confinement stress σ α . This is accounted for by using the principle of superposition, Mode I fracture toughness values have been previously used to investigate failure in garnet 31,32 and sulfide 33 systems. However direct application of experimentally measured fracture toughness values may not necessarily reflect accurate failure criteria, since material toughness for the case of an alkali metal being deposited into a crack involves the surface energy between the deposited metal and the electrolyte-in contrast to crack growth in the bare electrolyte where surface energy of the electrolyte only plays a role. As a point of comparison, differences have been reported in measured fracture toughness for β-Alumina in contact with Na under varying current densities compared to ex situ measurements. 23 No studies have been reported that try to quantify these effects for current MSSB systems.
We can in principle include these effects into a SIF based fracture criteria by using a lower effective SIF K IC,eff compared to the ex situ measured SIF K IC,m to determine fracture.

[ ]
Coupling boundary conditions.-Although explicit coupling between the ion transport 3, conductor 5 and metal 9 mechanical governing equations occurs in Eq. 3, where the ion flux expression includes the stress effects on the migration energy, in what follows we focus on the coupling between governing Eqs. 3 and 9 through the boundary conditions at the electrolyte/metal interface.
The boundary terms at the metal/conductor interface are shown graphically in Fig. 2. The boundary values of all terms will be labeled with a zero subscript for vectors and scalars and a zero superscript for tensors. The interface boundary condition for the diffusion potential includes an interface chemical potential term- which accounts for chemical interface reaction kinetics-plus a term associated with the work required to deposit ions against a normal surface traction. The diffusion potential boundary condition at the metal/electrolyte interface is, Where μ 0 is the aforementioned metal ion chemical potential at the interface which accounts for interface reaction kinetics. For fast electrode kinetics, this is constant and is set to zero. V m is the molar volume of the deposited metal, p 0 ( ) s is the metal boundary stress, and n is the surface normal vector. When μ 0 is constant, Eq. 13 accounts for the effect of mechanical stress on the ion diffusion potential.
Additionally an interface momentum balance couples Eqs. 9 and 5, Finally, the interface boundary strain rate condition depends on the deposition rate, which is proportional to the gradient of the boundary ion flux J 0 , The boundary conditions above lead to a coupled system of governing equations that are summarized for brevity in Table I. Additionally, a summary of the symbols used in the system of governing equations and the consecutive analysis is given in Table II.

Semi-Analytical Solution
We obtain two semi-analytical solutions to the set of governing equations presented. A set of simplifying assumptions with varying complexity regarding domain geometry and materials properties are made to obtain tractable solutions. Although the actual phenomena may exhibit more complex behavior, the results obtained provide qualitatively accurate results, and quantitative bounds.
In the most simple model we solve, boundary traction effects from mechanical stresses on the ion deposition are neglected. In other words the second term in the diffusion potential boundary condition 13 is ignored. Neglecting the surface traction term, makes the ion transport problem become a standard problem for an electrostatic field around a grounded conductor, and hence we refer to it as the purely electrostatic problem. We additionally solve a model which includes the effects of the surface traction term in Eq. 13, and which thereby allows us to quantify the effects of mechanical stresses on the ionic transport and deposition. We will present the model assumptions and the main results from solving these models. Further details on the solutions and all involved equations are given in the Appendix.
Model and assumptions.-The solutions are derived for a planar crack of length L and tip radius r 0 with circular tip as shown in Fig. 3. The ion concentration is constant and homogeneous within the conductor in order to satisfy local electroneutrality-hence does not play a role in the ion flux Eq. 2. We assume all materials properties to be uniform and isotropic. Metal is electrodeposited at steady state galvanostatic conditions. The metal chemical potential at the interface is set to zero under the assumption of fast and spatially homogeneous interface kinetics and negligible side reactions. Elastic terms in the ion flux Eq. 2 are neglected by assuming a negligible effect of elastic stress on ion migration energy compared with the electric potential effects. The current density in the local environment surrounding the flaw but sufficiently far from the interface and other inhomogeneities is assumed to be homogeneous with a value i 0 . The value of i 0 contains effects of imperfect contact and possible contact loss. 34 The local current density i 0 can be related to the cell current density by considering the effective interface contact area. 25 The ion transport continuity Eq. 3 simplifies to a Laplace equation of the generalized diffusion potential, Furthermore, by assuming the metal pseudo-viscosity η to be nearly constant for the plastic strain rates involved, the plastic flow Eq. 6 can be expressed as the Stokes equation from low Reynolds number hydrodynamics, Table I. Governing equations and interface boundary conditions for electrodeposition-induced plastic flow.

Region Equations
Conductor bulk Neglecting the elastic stress term in the diffusion potential, the elastic response in the conductor becomes fully decoupled from the ion transport in the conductor. Ion transport in the conductor and metal plastic flow remain coupled through the interface boundary condition, which setting the interface chemical potential μ 0 = 0 in Eq. 13 becomes, Such that the metal potential at the boundary is dictated by the normal stress that develops as Li ions are deposited and Li metal is pushed out of the flaw.
Using the definition of strain rate, the strain/flux boundary condition 15 can be re-expressed more simply in terms of the metal displacement rate, Hence, the ion transport 16 and metal flow 17 can be considered independently of the conductor elastic response. Subsequently the results can be used to analyze the conductor elastic response and/or potential fracture.
Electric field focusing effects.-To solve the model without boundary traction effects, the purely electrostatic problem, we obtain an expression for the field focusing at the defect tip. Analyzing this simpler problem allows a derivation of an upper bound on the deposition concentration and stress increase at sharp tips.
Similar to modeling work in field emission, 35,36 we quantify the degree of deposition inhomogeneity using a flux concentration factor (FCF) accounting for the amplification of the boundary ion flux at a sharp interfaces vs a flat interface, J 0 || || is the magnitude of the ion flux at the peak of a curved boundary, and J 0 || || is the magnitude of the ion flux at a reference flat surface, keeping all other boundary conditions equal.
We quantify the dependence of the FCF γ J on the crack geometry by expressing it in terms of a mean aspect ratio, defined as the product of the mean surface curvature at a peak and the length of the crack l k . For example, for the geometry shown in Fig. 3 with only one radius of curvature, the mean aspect ratio would simply be, Figure 4 shows this dependence obtained for a 3D flaw and a 2D planar flaw for values of the mean aspect ratio l k that can be expected for flaws in MSSB electrolytes. It is evident that the dimensionality of the flaw affects the FCF, with a 3D flaw having a higher FCF than a 2D flaw with the same mean aspect ratio. This observation allows us to estimate that real flaws with more complex geometries will exhibit FCFs somewhere between the 3D and 2D curves, as represented by the shaded region in Fig. 4. We see that field effects can increase local deposition by 3 to 4 orders of magnitude for the sharpest flaws expected (a flaw length of tens of micrometers with a tip radius of r 0 < 1 nanometers).
Metal flow at crack tips without a boundary traction potential.-We solve the metal flow Eq. 17, initially using the ion boundary flux for the purely electrostatic problem previously discussed. The solutions to the internal flow problem are obtained using a 2D stream function for a planar flaw (see Fig. 3). Further details of the solutions can be found in Appendix C.
The pressure field, ion flux lines and plastic displacement stream lines from the solution of the electrostatic problem are plotted in Figure 3. Simplified geometry used in semi-analytical derivations. Metal electrode/Inorganic Solid Conductor interfaces with metal filled planar flaw of length L and tip radius r 0 , held at stack pressure P 0 . Fig. 5. The corresponding metal displacement rate profile at a crosssection of the flaw tip is also shown. The ion flux lines plotted outside the flaw tip (white background) converge at the flaw tip boundary, exhibiting the concentration effects described in the previous section. We obtain a linear pressure dependence with the length along the flaw l. In addition we obtain a parabolic displacement rate profile, which along with the linear pressure drop are directly analogous to the pressure drop and velocity profile that occurs in internal channel flow at low Reynold's numbers. This suggests that using Poiseuille or channel flow based models such as those developed for molten Na electrodes, 19,20,22,24 may be applicable also for low yield metals, such as Li, at low deposition rates. However, we emphasize that due to the ion deposition boundary condition 15 the pressure always has a linear dependence at a flaw tip, even as the width of the corresponding channel shrinks to zero.
From our solution of the metal flow in the electrostatic problem, the viscous flow of metal, increases the hydro-static pressure increase along the length of a flaw as, Where α is a geometric factor of order unity dependent on the shape of the flaw. As could be expected, the pressure increase depends on the pseudo-viscosity η, and the metal mass deposited per unit area According to Eq. 22 high aspect or sharp flaws can-in theoryincrease pressure without bound for very sharp flaws. In this model, the pseudo-viscosity is a key material parameter controlling the pressure magnitudes that can occur-this is discussed in more detail in the discussion section. Finally, we estimate that the pressure increase strongly diverges as the flaw gets sharper ΔP ∼ (1/r 0 ) a where 2.5 ⩽ a ⩽ 3, as detailed in C Appendix C. Finally, we note that Eq. 22, giving the pressure increase, corresponds to the linear expression derived by Armstrong, 24 in contrast with the diverging expression incorrectly derived by Feldman, 22 for Na/β-Alumina cell degradation.
As previously mentioned, solving the purely electrostatic problem in the solid electrolyte provides an upper bound on the stress increase. Furthermore, the results from the purely electrostatic case and subsequent driven flow become a good approximation for moderate pressure increases from either low current densities, blunt flaws, soft metal or any combination thereof. In the next section we investigate in more detail the mitigating effect of mechanical stress induced potential changes which can redirect current.
Metal flow at crack tips with a boundary traction potential.-In the simplest model above, we obtained a solution from the electrostatic version of the model, which results in a pressure increase that is directly proportional to the product of the ion flux concentration factor times a geometrical concentration factor from the metal viscous deformation and flow. The pressure increase in that model diverges strongly as the tip radius gets smaller. It is clear that the pressure cannot grow unbound, and so it becomes imperative to consider mechanisms that accommodate and potentially mitigate the ion flux concentration and resulting pressure concentration factor. The boundary traction potential in Eq. 13, can give rise to ion redistribution, which serves as a stress relaxation mechanism effectively lowering the total pressure increase. At high enough pressures, the boundary traction potential term can lead to diffusion potential surface gradients that screen out the concentration effects at sharp flaws. We show that this screening effect is very significant for large pressures and that neglecting it leads to an overestimation of the flaw tip pressure.
We solve the ion transport and metal viscous flow in Eqs. 16 and 17, using the same techniques as used in the electrostatic problem, but now include the effects of boundary traction. Details of the solution are given in the Appendix. The solution to the pressure increase at a flaw tip given in Eq. 22 is given in terms of the FCF γ J , Figure 4. Flux concentration factor (FCF) arising from the electrostatic field with a constant potential boundary condition. The 2D curve is calculated using an elliptic cylindrical geometry and the 3D curve using a spheroidal geometry as detailed in the Appendix A. Shading corresponds to FCF range for more complex flaws. the tip radius r 0 , the metals pseudo-viscosity η and molar volume V m , and the operating current density i 0 as a function of the position along the flaw l, Where we identify an electro-plastic (EP) screening length ξ, defined in terms of only materials properties: the ionic conductivity of the electrolyte σ + , the molar volume of the metal V m , and the pseudo-viscosity of the metal η. The EP screening length, is a material system property that determines a length-scale at which all deposition concentration effects at defects of lower tip radius are effectively screened and deposition rate inhomogeneities become negligible. The stress at a flaw with a radius of curvature smaller than that determined by the EP screening length only increases proportional to the length of the flaw. For tip radii sufficiently smaller than the EP screening length, the pressure increase becomes, For this limiting case, the pressure increase can be recognized as simply caused by the Ohmic voltage drop in the electrolyte along the length of the flaw. Figure 6 shows the effects of the EP screening length on the diffusion potential in the solid electrolyte. Diffusion potential isocontours are plotted for different ratios of EP screening length to flaw tip radius. As the ratio increases, the isocountours become more and more uniform-horizontal equipotential lines. At values where ξ ≈ 2r 0 , the diffusion potential is essentially uniform, and the concentration effects from the flaw are fully screened. At this point, the pressure increase at the flaw is given only by that of the ohmic voltage drop in Eq. 25, which in contrast with Eq. 22 from the electrostatic solution, does not increase unbounded. In other words, the pressure does not increase with the sharpness of the flaw.
Using the pressure increase Eq. 23 and the SIF for internally pressurized cracks 10, we derive an expression for the SIF at flaw.
Similarly, for the saturated case, using Eq. 25, the value of the SIF for a screened flaw is, Where P 0 is the cell operation stack pressure. β 1 and β 2 are geometrical constants of order unity. For our model geometry, the values obtained at saturation are β 1 = 1.000 84 and β 2 = 1.090 87 from evaluating Eq. 10, using the parametrization given by Hartranft. 30

Discussion
We have shown that metal filled flaws in a solid ionic conductor can lead to concentration of deposition current density and as a consequence pressure increases along flaws. When no effect of stress on deposition is assumed, we obtain a pressure increase that diverges  as flaws become sharper. In theory this allows for pressure to grow until fracture is triggered. Versions of this model have been used to rationalize experimental results. 5,6 However, from a more sophisticated analysis in which the effects of deposition against an emergent boundary traction are included, it is clear that ion transport is redistributed and the pressure concentration at flaws can be partially or fully screened. From these results, unlike in previous levels of analysis, fracture is not guaranteed. The occurrence of fracture rests on the degree of screening and rates of ion redistribution possible. Before further discussing model results, we consider possible values of the metal pseudo-viscosity to better understand in which regime Li deposition from a solid electrolyte occurs. The previous derivations point at pseudo-viscosity as a critical property required to predict pressure increase values, the extent of deposition inhomogeneities, and subsequent failure criteria. We currently are unaware of any published studies quantifying the pseudoviscosity or any equivalent property of solid alkali metals. We estimate such values from recent studies of Li yield stress using Eq. 8. Figure 7 shows reported values of viscosity of molten Li at different temperatures. Figure 7 also shows an approximation of solid Li pseudo-viscosity based on direct application of Eq. 8 for multiple sets of values for Li yield stress, strain rate and temperature reported in each reference. [37][38][39][40] The first thing that stands out is the enormous range of magnitude for the approximated Li pseudoviscosities. Furthermore, it is important to note that although very different magnitudes of Li yield stress have been reported, 37-40 the highest value of the approximated pseudo-viscosity does not correspond to the highest yield stress value. This observation suggests that in addition to determining the yield strength of Li, the plastic behavior, included here as the pseudo-viscosity, at different strain rates is also crucial to understanding and predicting cycling behavior and failure in MSSBs.
Model results of the electrostatic problem-excluding the boundary traction term-predict that failure via fracture can occur for metals with high enough pseudo-viscosity-certainly for any of the values shown in Fig. 7 for Li metal. As was shown in Fig. 5, deposition onto a perfect conductor in a sharp flaw leads to a linear pressure distribution with a parabolic displacement rate profile along the flaw length. This is precisely analogous to Poiseuille or channel flow which occurs for fluid flows in pipes and channels where inertia terms are negligible compared to viscous terms. The channel flow framework has previously been used directly to study molten Na cells, 19,22,24 and alluded to in recent work for solid Li electrodes. 6 The pressure increase leads to values of the stress intensity factor (SIF) for internally pressurized cracks, 30 that easily meet fracture criteria based on reported conductor toughness values, and the estimated pseudo-viscosity values shown in Fig. 7. Figure 8 shows color maps of predicted values of tip pressure increase, and effective strain rates when the boundary traction term is neglected at a deposition current density of 0.1 mA cm −2 . The pressure has been normalized by the pseudo-viscosity. To determine the absolute value of the pressure increase, the magnitude shown needs to be scaled by the corresponding pseudo-viscosity. Critical regions-where flaws of those dimensions would lead to fracturefor given values of pseudo-viscosity are delimited from above by solid lines.
From the predicted critical dimensions shown in Fig. 8, it becomes evident that a pseudo-viscosity value above ∼100 Pa · s, leads to a non-negligible critical region where conductor fracture can occur for most flaw geometries expected in inorganic solid electrolytes. Furthermore any value of pseudo-viscosity in the range of those shown in Fig. 7 renders all defects critical, regardless of dimension.
The calculated strain rates for the simplified electrostatic model without any current redistribution are also shown in Fig. 8. Since the plastic flow rate is set by the deposition rate, the strain rate value is independent of the pseudo-viscosity. The values of strain rates that occur for the given current density are well above those reported in studies on Li yielding, with the exception those reported by Herbert. 37 Furthermore, it has recently been observed that at strain rates higher than 1 × 10 −2 s −1 Li hardening becomes important. 37,42 This suggests in depth studies on the plastic response and its possible strain rate dependence are critical for further understanding of Li deposition at flaws near the metal/electrolyte interface.
Our analysis shows that in a framework that does not incorporate the effects of the developed stress on ion diffusion potential, any sized flaw expected to occur in a ceramic ionic conductor will drive fracture even at low current densities (0.1 mA cm −2 ), and leads to extremely large pressure values in the flaw, which casts doubt on its predictive value. A purely electrostatic treatment of the conductor diffusion potential severely overestimates the current focusing on flaw tips. In a more realistic scenario, the pressure buildup in the flaw caused by the increased deposition will elevate the ion diffusion potential and redirect ion flux away from the tip. This physical situation is accounted for when including the boundary traction term in the diffusion potential boundary condition. In the analysis including the boundary traction term, the magnitude of the EP screening length is critically important to determine the effects of current focusing. This length scale depends both on mechanical properties of the metal and the conductivity of the electrolyte.
We have shown in Fig. 6 the effect of the EP screening length on the current concentration at a flaw tip. As evidenced by the diffusion potential contour lines, whenever the EP screening length exceeds the tip radius by more than 2 times, the geometrical focusing effects are fully screened out. Thus for any ionic conductor with a conductivity above 1 × 10 −2 mS cm −1 , and any metal with pseudo-viscosity close to the values shown in Fig. 7, geometrical current density concentration effects are predicted to be effectively screened for metal filled flaws that may occur within a solid electrolyte. It must be noted though, that predicted geometrical concentration effects are non-negligible for molten metal viscosity values, which means that failure of solid electrolytes can be very different when a molten metal is used.
As shown by Eqs. 23 and 24, screening of flaws improves as the conductivity of the conductor and the pseudo-viscosity of the metal increase. Hence, using a low conductivity value of σ + = 1 × 10 −3 mS cm −1 , and a lower bound for pseudo-viscosity of η = 1 × 10 −3 Pa · s of molten Li viscosity at 180°C, gives a lower estimate of the EP screening length of ξ ≈ 4 × 10 −4 nm. At this lower value we can expect strong effects of geometrical and ion flux concentration on the pressure increase at flaw tips. An upper bound estimate using σ + = 10 mS cm −1 and η = 1 × 10 7 Pa · s, gives ξ ≈ 130 μm. At the upper estimate for the screening length all flaws with radii below this size would be effectively screened, and pressure increase is bounded and will only increase along flaws according to Eq. 25. From these estimates it is again clear that obtaining a more precise value of the pseudo-viscosity is necessary to narrow the range of possible screening lengths. However, for solid metal electrodes, and current state of the art electrolytes (σ + ⩾ .1 mS cm −1 ), the screening lengths are expected to be of the order of tens of micrometers or higher. Figure 9 shows the maximum pressure increase which occurs at the tip of a 20 μm flaw for an electrolyte with conductivity σ + = 0.1 mS cm −1 . We note that for ionic conductivities of current state of the art solid electrolytes, the pressure increase quickly saturates to the value given by Eq. 25. We also restate that the saturation pressure is independent of how sharp the flaw is. The tip radius/curvature only affects how quickly the pressure saturates with respect to the EP screening length. Figure 9 also shows the ratio of the tip pressure increase with boundary traction correction (Eq. 23) over the pressure increase without this correction (Eq. 22). The abrupt drop in this ratio shows that the pressure increase values predicted with the pressure correction are of substantially smaller value and bounded by the saturation pressure 25 compared to the unbounded values predicted without considering the boundary traction correction.
Including the boundary traction term has a substantial effect on the magnitude of pressure difference obtained for flaw tips. In the fully screened case the pseudo-viscosity becomes irrelevant to the pressure increase since the pressure increase saturates at the maximum value shown. This occurs because with a uniform deposition flux, no plastic deformation occurs in the metal filling the flaw, since deposition will not generate shear strain according to Eq. 15. This behavior suggests that plastic deformation can only occur in the bulk of the electrode. Based on this behavior, the fundamentally different plastic response and yield stress of the metal in confined and reduced dimensions is not of much importance to this problem. Furthermore, as a result from pressure increase saturation, the SIF is also independent of the pseudo-viscosity as shown in Fig. 10. For the range of pseudo-viscosity values shown in Fig. 7, the SIF is constant. At saturation, the only mechanical variable that has an effect on the SIF is the applied stack pressure P 0 .
In contrast to the negligible effects of the pseudo-viscosity on SIF values, the ionic conductivity plays much more important role. Since the increase in pressure from deposition is set by the saturation pressure value. As given in Eq. 25, the ionic conductivity is the only materials property that directly dictates this maximum pressure. Therefore, higher ionic conductivity can effectively reduce the possible SIF values directly by lowering the saturation pressure value as shown in Fig. 10.
A final important feature to notice from Fig. 10 is that the applied stack pressure directly increases the SIF value and so can bring its magnitude above the solid electrolytes fracture toughness (0.23 MPa m 33 reported for LPS and 1 MPa m for LLZO 31,32 ). Even without any deposition, a high enough stack pressure can raise the SIF for a metal filled flaw to the critical value and initiate fracture from purely mechanical means as a result of the increase in hydrostatic pressure of the metal. Recent experimental work reports that shorting from purely mechanical means is observed in Li and argyrodite Li 6 PS 5 Cl cells at a stack pressure of 75 MPa. 8 Severe cracking was observed in these shorted samples, suggesting that fracture without any deposition can in fact occur. Finally, we note that the SIF values calculated do not include the external stress from cell constraints as given in Eq. 11, and therefore represent a worse case scenario calculation since an external compressive stress can act to reduce the total SIF.
Considering that the pressure increase and stress intensification caused by deposition from solid electrolytes occurs at their saturation values, we can estimate conditions under which a flaw becomes critical based on its length, the electrolyte conductivity, the local current density, and the stack pressure. The flaw becomes critical when the saturation SIF given in Eq. 27 reaches the electrolyte fracture toughness. Figure 11 shows sub-critical and critical regions for a 10 μm sized flaw for an electrolyte with MPa m 1 fracture toughness for different values of ionic conductivity. It is clear that increasing the ionic conductivity has a pronounced effect in increasing the size of sub-critical regions by directly affecting the slope of the critical region boundary. Similarly, Fig. 11 also shows sub-critical regions for different values of electrolyte fracture toughness for a given ionic conductivity of σ + = 0.5 mA cm −1 . Increasing the materials fracture toughness essential increases subcritical regions by setting the maximum stack pressure that a MSSB can operate at.
Ionic conductivity and fracture toughness are shown to be key electrolyte properties that can be optimized to increase the size of subcritical regions and therefore ensure stable deposition in MSSBs. The ionic conductivity has a stronger effect in increasing the maximum local current densities that can be accommodated. In contrast the fracture toughness plays a more important role in setting the maximum stack pressure that can be used during operation. In reality local current density and stack pressure are not decoupled, and their dependence and effects from ionic conductivity and fracture toughness involve other phenomena, such as interface reactions, and contact loss, that we have not considered.
In cases where the solid electrolyte is not stable against the metal electrode, interface reaction products may become relevant in determining the EP screening length. Interface reaction products still exhibit some level of ionic conductivity for most proposed and currently used electrolyte materials. 43,44 The extent to which screening changes will be dependent on the thickness of the interface layer. More precise quantification of the changes in screening and thickness thresholds require further development of our current model to account for this increased complexity. However we expect that for reaction products with conductivity above σ + = 1 × 10 −2 mS cm −1 the SIF due to deposition at flaws will still occur at or near the saturated value, albeit at a higher saturation value than that predicted for the electrolyte only.
Finally, we stress that plots such as those in Fig. 11 currently serve as materials and cell operation design guidance, but do not serve as predictors of the exact conditions for which fracture is initiated. In order to reach more precise predictive values, several other phenomena must be accounted for in the model. Notably, contact loss and material inhomogeneites (microstructure, interphase composition), which operate at a larger scales (compared to sharp flaws) but can lead to different values of the local current density near flaws.
For example, contact loss and void formation during striping, have been found to have an important effect on the critical current  value during plating below which metal propagation through the electrolyte is not observed. 34 Effects of contact loss would show up as an increase of the local current density i 0 near a flaw, which directly translates to a higher saturation pressure at the flaw tip based on Eq. 23. In addition, void formation and changes in interface contact during battery cycling are modified by the operating stack pressure. 34 This creates a dependence between the stack pressure and the value of the local current density near a flaw, which cannot be explicitly accounted for in the current model solutions, but should be considered when reading and using Fig. 11 as guidance.
Additionally, improved performance during cycling at higher temperatures 45 may be attributed to better interface contact since the metal electrode becomes more ductile. This ductility increase with temperature for Li has been recently shown to correspond to a decrease in Li creep stress. 42 These effects are also seen as a decrease in the pseudo-viscosity for a given strain rate, as observed in Fig. 7. However, the reduction in pseudo-viscosity does not have any appreciable effect on reducing the overall SIF, which can be seen in Fig. 10, where the SIF remains saturated well below the lowest pseudo-viscosity calculated in Fig. 7. In contrast, an increase in temperature also leads to improved ionic conductivity, which based on Fig. 10 will have a noticeable effect in reducing the SIF. Therefore, improvements in performance from elevated temperatures in addition to improved interface contact also involve lower SIF values at flaws from an increase in ionic conductivity, and therefore an overall lower chance of cracking in the electrolyte.
Furthermore, since fracture toughness directly affects the stable deposition regime in Fig. 11, tailoring the fracture toughness of the ceramic conductor can be attempted using insights from other toughening studies of ceramic materials. For instance, by adding particulates to the ceramic matrix, 46,47 or by tuning microstructure to exhibit elongated grains, 48 electrolyte fracture toughness could be improved. Tuning electrolyte particle size may also be an option, however some studies have shown no discernible effect on toughness. 46 Improving electrolyte densification and reducing porosity is expected to be beneficial since lower porosity has been shown improve fracture toughness, 49 in addition to reducing available paths for metal creep. Lastly, the precise value for materials fracture toughness must be considered under electrodeposition conditions. Previous studies on β-Alumina have shown that exposure to Na does indeed lower the fracture toughness. 23 Similar studies for Li systems would be highly informative.

Conclusion
We derive governing equations for the coupled transport, plastic deformation, and elastic response equations describing metal deposition from solid electrolytes. Using a set of simplifying assumptions and a simple domain involving a planar flaw filled with metal, we obtain a semi-analytic solution to the set of governing equations. Subsequently we calculate order of magnitude predictions on current focusing, strain-rates, and pressure buildup. An upper bound on the current density concentration factor of ∼10 3 at a sharp defect tip is obtained using a purely electrostatic version of the model (no stress boundary traction included in the ionic diffusion potential). Neglecting this boundary traction term, the model suggests that all defects become critical for even the lowest expected pseudo-viscosity values, even under low current densities of the order of 0.1 mA cm −2 .
When the boundary traction term in the diffusion potential is included, an electro-plastic screening length ξ emerges as a critical length scale. The EP screening length is set by mechanical properties of the metal and the ionic conductivity in the electrolyte. Our analysis predicts that for screening lengths exceeding 2 times the radius of curvature of a flaw, the current focusing effects at a defect are effectively screened out. The increase of the ion potential from the pressure buildup leads to redistribution of metal ions which effectively homogenizes the ion flux approaching the sharp tip. For a screened flaw, the pressure buildup at a tip is solely caused by the potential difference arising from an Ohmic drop in the electrolyte along the length of the flaw. The pressure increase is of the order of hundreds of MPa and the resulting stress intensity factors can be within the order of magnitude of experimental fracture toughness values measured for LLZO and LPS electrolytes.
We also find that negligible deformation occurs in the metal within a screened flaw, since only inhomogeneous deposition gives rise to shear strains. Thus the model suggests that any plastic deformation in the metal electrode should occur in the bulk of the electrode metal, and hence the mechanical properties of the bulk electrode are more important than deviations from these at small scales and confined dimensions.
Our analysis of the saturation pressure at a flaw tip highlights the importance of high ionic conductivity in the electrolyte, since high values of ionic conductivity increase the region of stable operation by accommodating much higher local current densities around sharp flaws. Hence a high ionic conductivity electrolyte is considerably more tolerant to flaws. The electrolyte fracture toughness on the other hand permits much higher stack pressures to be used during deposition. High stack pressures have been shown elsewhere to improve performance by reducing interface contact loss. 8,34 Therefore, continued efforts into developing faster superionic conductors as well as tailoring their fracture toughness is shown to be the two most crucial materials considerations to prevent MSSB failure from electrolyte fracture.
Finally, to further understand and more precisely determine relevant materials design and operating conditions to reduce or prevent MSSB degradation and failure from mechanical fracture, the plastic behavior of the bulk metal electrode, and the fracture toughness of the solid electrolyte need to be investigated under consistent conditions to those expected during cell operation. Specifically, understanding Li/Na metal plastic behavior beyond yielding becomes crucial to effectively determine a more precise value for the EP screening length. Similarly, the fracture toughness of the solid electrolyte in contact with the electrode, under operating overpotentials, and considering possible interface reactions needs to be further investigated to predict precise numerical values for failure criteria. Additionally, further complexities of solid state deposition, such as contact loss, deposition kinetics, interface reaction products, and specific values of materials properties based on microstructure should be investigated in order to successfully analyze and predict fracture in experimental samples. r r