Interfacial studies on the effects of patterned anodes for guided lithium deposition in lithium metal batteries

Long-range electric vehicles and long-duration renewable energy storage are critical to moving toward a more energy-efficient future. However, due to the limited energy density of current battery technologies, such as lithium-ion (Li-ion) batteries, higher energy density battery technologies are needed to meet the expected future energy storage demands. Lithium (Li) metal batteries have emerged as a potential alternative due to their high-density anode material with a theoretical capacity (∼3860 mA h g⁻¹), which is about ten times larger than that of a graphite anode (∼372 mA h g⁻¹) used in Li-ion batteries.¹ While Li metal batteries show promise for higher energy density storage, the growth of Li dendrites during the plating and stripping processes decreases Coulombic efficiency and creates safety concerns.^1–3 These dendrites can break off, forming dead Li, and potentially cause a short circuit in the battery and thermal runaway, which can lead to catastrophic failure.^4,5 Therefore, it is critical to suppress dendrite growth.

Different methods have been proposed to aid in dendrite suppression to varying degrees of success. Proposed methods include electrolyte addatives,⁶ 3D current collectors,⁷ porous separators,⁸ solid composite electrolytes,⁹ and artificial solid electrolyte interphase (SEI) layers.¹⁰ Mechanical modification to the anode surface has also been researched as a way of limiting growth. It has been shown that the Li metal possesses a natural roughness that has a direct impact on the dendrite growth and morphology.¹¹ Becking et al.¹² decreased this natural roughness by roll pressing the anode and showed a smooth surface that helped to decrease heterogeneous deposition, leading to a smoother surface after cycling.

In comparison to mechanically smoothing the anode surface, several groups have taken an alternate route and created a patterned anode surface to help guide Li deposition. Modifying the anode surface via nanochannels, for example, allows the manipulation of lithium ion (Li⁺) flux above the anode. Liu et al.¹³ proposed modifying the anode surface via a thin polymer coating that contained vertically aligned nanochannels. The channels cause the Li⁺ to move normal to the anode surface instead of parallel, helping to control the flux near the Li deposits that sit inside the channels. The size of these nanochannels also has an impact on the morphology of dendrites. Sharon et al.⁵ found, in a recent study, that at smaller channel widths, the dendrites were denser and less branched. Conversely, as the confinement width increased, the dendrites were longer and more branched.

While some researchers have used the addition of a nonreactive layer to the anode surface, others have attempted to control growth through modifications to the anode itself, such as creating patterns in the Li metal surface. Guided Li deposition via mechanical modification on the anode surface has experimentally shown that the shape and size of the patterned imprints have an impact on how Li deposits on the surface. In the work of Kim et al.,¹⁴ two different sized square hole patterns were evaluated, showing that for larger hole patterns (∼40 μm wide), Li deposited initially at the bottom surface and top corners of the hole. The cavities were then fully filled with precipitates, followed by additional deposition on the top surface, leading to a flat anode. Although the surface area of the smaller holes is higher than that of the larger holes, the increased resistivity difference between the top and bottom of the well leads to preferential deposition on the bottom surface. In addition, Li et al.¹⁵ used circular patterns and showed that Li will deposit in areas where the curvature is higher such as the top corners of the wells, agreeing with the work of Kim et al.¹⁴ This phenomenon occurs because at a point of high convexity, such as the corner of the hole, the charge density is significantly higher than elsewhere along the anode. This increase in the charge density causes increased Li⁺ concentrations at the corners, leading to more Li depositing at these locations.^15,16

These experimental interfacial studies indicate that controlling the interface morphology can have a significant impact on dendrite growth in Li metal batteries by influencing the local conditions at the interface. It is difficult to experimentally determine variations in concentrations, current density, or other properties along the interface under in situ or operando conditions. Since the dendrite growth and morphology are largely influenced by local current densities, overpotentials, and Li⁺ concentrations,^3,17,18 resolving these complex interfacial phenomena is critical to overcoming dendrite growth. Computational modeling at the appropriate scales can give insights into the chemical and physical phenomena contributing to dendrite growth and can isolate critical physics to understand their influence on dendrite growth and battery performance. Recent computational efforts ranging from atomistic¹⁹ to continuum scale²⁰ have provided insights into the effects that various physical parameters have on the electrodeposition morphology. These studies include how overpotential,²¹ thermal effects,²² and mechanical effects, such as an applied external pressure,²³ impact the morphological evolution of dendrite growth.

This paper presents a computational model of the electrode–electrolyte interface [the extended Butler–Volmer-smoothed particle hydrodynamics (eBV-SPH) model], which implements an extended Butler–Volmer (eBV) equation^4,17,24,25 and uses the smoothed particle hydrodynamics (SPH) computational method^27–28 to simulate the dynamic moving interface, capturing the dendrite nucleation and morphology along the electrode surface. The Butler–Volmer (BV) equation describes the rate of reaction at the electrode–electrolyte interface. By utilizing the expanded version of this equation, the rate at which Li is locally deposited on the anode surface can be carefully tracked along with other important phenomena, i.e., current density, overpotentials, and Li⁺ concentration. In this work, the eBV-SPH model is used to investigate the effects of overpotential on the dendrite morphology and how different sized micropatterns in the Li metal anode influence Li deposition at the interface for the design of better batteries.

The SPH computational method^27–28 is used to simulate the electrode–electrolyte interface. SPH is a Lagrangian, particle-based computational fluid dynamics method that is well suited to moving and deforming boundaries, such as dendrite growth, due to its particle nature. SPH is able to capture moving boundaries inherently without the need for front tracking methods and can capture the complex geometry of dendrites.¹⁷ The model used in this work builds off of previous modeling of Li metal batteries²⁵ and includes Butler–Volmer kinetics.^30–31 

To understand the effects of the surface morphology on dendrite growth, the simulation domain includes the anode surface, an electrolyte layer near the interface where diffusion occurs, and dendritic structures; see Fig. 1. The concentration of Li⁺ and anions is assumed to be initially constant throughout the domain and constant in the bulk of the electrolyte, a distance L from the anode surface (outside the diffusion layer),

As a result of the applied potential during charging, the anions migrate away from the anode, while the cations (Li⁺) migrate toward the anode. Li⁺ builds up at the anode surface and reacts with electrons, depositing as the Li metal. This causes a Li⁺ concentration gradient and the subsequent diffusion flux. The general reaction between the Li⁺ and electrons is⁴ 

The eBV-SPH model includes the transport of Li⁺ through the electrolyte and the plating of Li on the surface, which leads to dendrite structures. The eBV-SPH model extends previous work on dendrite growth and dissolution by including the effects of the potential and current on transport and Li deposition.^17,25,32 This is accomplished by including an extended Butler–Volmer (eBV) equation that relates the local concentrations to the local current density (i),²⁹ 

where α_c and α_a are the cathodic transfer and anodic transfer coefficients, F is Faraday’s constant, R is the gas constant, T is the temperature, η is the overpotential, $Cc(rs⃗,t)$ and $Ca(rs⃗,t)$ are the time dependent concentrations of the cations and anions at the surface of the anode, C_c,0 and C_a,0 are the constant concentrations in the bulk electrolyte, i_o is the exchange current density, which is a function of the bulk concentration, and r_s is a point on the electrode–electrolyte interface, Γ.

The exchange current density, i₀, is dependent on the bulk concentration, given by³⁰ 

where k⁰ is the reaction rate constant.

The overpotential, η in Eq. (1), is

where ϕ_APP is the applied potential and ϕ_eq is the equilibrium potential.

The conventional BV equation assumes that the concentration of the reactants and products near the electrode is constant. In this model, the local concentrations are governed by the mass transport equations. The eBV equation considers changes in those concentrations and is related to the reaction kinetics at the electrode–electrolyte interface via Faraday’s law,^33,34

where k_J is the reaction rate flux and $irs⃗,t$ is the local current density given by Eq. (1).

The Nernst–Planck equations solve for mass transport of ions in the electrolyte,

where the subscript i denotes the ion type, cation (c) or anion (a). The first term on the right-hand side of the equation is the diffusion term, where D_i is the coefficient of diffusivity. The second term tracks migration of these ions where μ_i is the migration mobility relating how fast an ion moves with respect to the electric field, and ϕ is the local potential in the electrolyte. The electrolyte fluid domain is denoted as Ω_F.

In this model, the eBV equation is used to model the reduction reaction at the electrode–electrolyte interface,

where the subscript c denotes the cations.

Combining Eqs. (1)–(4) gives the reaction term for Eq. (6),

There is a zero-flux boundary condition for anions at the electrode–electrolyte interface,

The SPH implementation of the governing equations includes the discretization scheme of Tartakovsky et al.³² for Fickian diffusion and precipitation, the discretization scheme of Cannon et al.²⁵ for the migration term of the Nernst–Planck equation, and the continuum surface reaction (CSR) method³⁵ is used to implement the heterogeneous reaction boundary condition. The final forms of the discretized equations are summarized below, and further details on their implementations and verification can be found in the references.

The Nernst–Planck equation is discretized as^25,27,32,35

where m, ρ, D, and C are the mass, density, diffusion coefficient, and concentration, respectively. The subscript i represents the species (anions or cations), superscripts a and b represent individual SPH fluid particles, and k represents the solid particles. μ is the mobility, ϕ is the potential, and $rab⃗$ is the distance between the particles a and b. $n⃗$ is the normal vector, β is the characteristic function used to identify the anode interface,³⁵ k_J is the concentration dependent BV equation, and W is the SPH weighting function.

In our system, the concentration of the Li⁺ is the oxidized species, and the mass of the solid is the reduced species such that the discretized form of Eq. (7) is

where $Mk(rs⃗,t)$ and $Cca(rs⃗,t)$ are the time-dependent mass and Li⁺ concentration at the surface of the anode, and M₀, C_a,0, and C_c,0 are the constant mass, anion, and cation concentrations, respectively, in the bulk electrolyte and anode. The subscripts a and c represent the species (anion or cation), and the superscripts a and k represent individual SPH fluid and solid particles, respectively.

Furthermore, the precipitation is modeled via the discretization scheme of Tartakovsky et al.³² and is implemented into the eBV-SPH model as

where the superscripts k and a represent solid and liquid SPH particles, respectively. During a charging cycle, the rate of loss of the solute in the liquid phase must be balanced by the rate of gain of mass of the solid (M^k) due to precipitation. As the local Li⁺ concentration at the interface decreases, the mass of the solid Li particles will increase. Once the mass of a Li particle is twice the initial mass, the nearest fluid particle will “precipitate” and become a new solid particle. Therefore, the change in volume of solid particles will exactly match the volume change in liquid particles. This method has been implemented and verified by Tartakovsky et al.³² and in several previous studies from our group.^4,18,25

The governing equations [Eqs. (9)–(11)] are implemented in an SPH module in the open source Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) code base.²⁷ For all simulations used in the verification of the eBV-SPH model, the simulation domain is 5 × 5 µm² and is made up of 360 000 discrete SPH particles with a density of 30 particles μm⁻². It contains the anode–electrolyte interface and a finite diffusion boundary layer, as shown in Fig. 1.

The implementation of the eBV equation is an extension of a previous SPH model that has been verified and presented by Cannon et al.²⁵ However, the use of the eBV equation to model the reduction reaction at the electrode–electrolyte interface has not been verified. To accomplish this, the current density and applied potential are compared to the standard Butler–Volmer equation and the Tafel equation.

First, we compare the predictions of the eBV-SPH model to the BV equations using a simplified setup that includes diffusion, migration, and a reactive surface but does not model plating/stripping. From our eBV-SPH model, we obtain a local mass flux (N_i) at the surface via Eq. (10),

where m_Li is the atomic mass of Li and k_J is the eBV equation [Eq. (10)]. The local mass flux is then used to calculate the overall current density as³¹ 

where F is the Faraday constant, z_i is the charge number, and the subscript i represents the species. This value is then compared to the predictions of the standard BV equation,

To obtain the data used in the verification, the applied electric field was varied over multiple simulations to vary the surface overpotential. The mass flux was then calculated following the application of the charging voltage. The mass flux was obtained early in the simulation when the system was still in a reaction-limited regime, allowing analysis of how the reaction equation was behaving in the eBV-SPH model. The data were used to calculate the predicted current density through Eq. (13). As shown in Fig. 2(a), the current density predicted by the eBV-SPH model compares well with the current density of the standard BV equations with an average error of 6%.

Based on the data from Fig. 2(a), an exponential relationship was observed as the applied potential is increased agreeing with the expected result based on the BV equation. As the applied surface overpotential becomes more positive, the second term in the BV equation becomes negligible, which leads to the exponential relationship shown in Fig. 2(a). The exchange current density of the standard BV equation was calculated from the Tafel plot of the data. The log of the current density was plotted against the overpotential in Fig. 2(b), resulting in a Tafel plot. The exchange current density can be solved via the Tafel equation,²⁹ 

where A is the Tafel slope given by a linear fit of the data in Fig. 2(b) and η is the overpotential. The intercept of the fit gives the exchange current density.²⁹ This method for obtaining the exchange current density was used because the eBV-SPH model is not a constant current system. This means that as the charging voltage is applied, there is an initial spike in flux near the surface that decreases as the simulation proceeds, leading to the same trend in the simulated current density. Therefore, measuring the exchange current density from the Tafel plot allows for a better comparison of the simulation data to the trends from the BV equation.

In addition to verifying the implementation of the governing equations, the accuracy of the eBV-SPH model to represent a physical system can be assessed by a comparison to the experimental data. Many researchers have studied and visualized the dendrite growth and morphology under varying conditions.^5,36,37 The challenge in experimental validation of the eBV-SPH model is the lack of quantitative data at the interfacial level. As such, the best data available allow for a qualitative comparison of the dendrite growth and morphology. In particular, two experimental studies are compared to the eBV-SPH model, one considering the effect of the current density on the dendrite morphology³⁶ and another looking at the effect of nanochannels.⁵ 

Bai et al.³⁶ studied the effects of current density on the dendrite morphology. They found that at higher current densities, the system enters a transport-limited regime earlier than at low current densities. In this transport-limited regime, the Li⁺ concentration decreases significantly at the base of the dendrites, and the dendrites grow from the tips, leading to a longer and more branched morphology [Fig. 3(a)]. At lower current densities, the system will remain in a reaction-limited regime for long periods of time, and thus, the Li⁺ concentration at the base of the dendrites is not depleted as quickly, allowing for the growth throughout the entire dendrite as compared to the tips of the dendrites. Thus, denser and less branched dendrites were observed.

To study the effects of current density, four different simulations with varying applied electric fields, in turn, adjusting the applied potential, were run using the eBV-SPH model. The applied electric fields ranged from 16 to 50 m V μm⁻¹ resulted in current densities ranging from 0.26 to 6.94 mA cm⁻², respectively. Additional simulation parameters are in Table I. This range was chosen because it encompasses the experimental range used by Bai et al.³⁶ (1.09 to 6.78 mA cm⁻²). The current densities were calculated from the simulations as discussed previously in the verification of the eBV-SPH model [Eqs. (12) and (13)]. While the current densities were consistent between the experiments and the eBV-SPH model, the exact experimental conditions of Bai et al.³⁶ were not simulated, namely, the time scales are different, and this leads to different lengths of dendrites but the same morphology trends. Figure 4 shows the dendrite growth predicted by the eBV-SPH model under different electric fields. The solid particles, including the anode and dendrite growth, are colored gray to distinguish them from the electrolyte so that the morphology of the dendrites can be visualized. From the simulation results, as the current density and electric field increase, the dendrites get longer with more of a branching morphology making them less dense [Fig. 4(d)]. The simulations using a low electric field and current density [Figs. 4(a) and 4(b)] showed a more even distribution of growth across the anode surface with denser dendrites.

Looking at Fig. 4(a), the current density is the smallest, and the concentration of Li⁺ near the interface is higher than in the simulations with a high current density [Fig. 4(d)]. Under these conditions, the system remains in a reaction-limited regime, leading to dendrite growth that is dense and more uniform. As the current density increases, the system transitions into a transport-limited regime early in the simulation, and the Li⁺ are depleted at the base of the dendrites. This causes growth to occur at the tips of the dendrites, resulting in long and thin dendrite structures [Fig. 4(d)], which agrees with the morphology of the experimental results [Fig. 3(b), 6.78 mA cm⁻²].

In addition to the changes in the current densities, the spatial distribution of the local potential varies throughout the electrolyte based on the applied electric field. Near the anode, the local potential remains constant and then increases linearly through the electrolyte. At low applied electric fields [Fig. 5(a)], the local potential is small. Under these conditions, the dendrites grow in a uniform fashion and will remain in the portion of the electrolyte where there is a low potential. However, when a large electric field is applied, the dendrites grow quickly with a more branched morphology [Fig. 5(b)]. Due to the linear increase in the potential from the electrode to the electrolyte region, as the dendrites increase in size, they grow toward the higher local potential region. The tips will be surrounded by a larger potential, thus leading to more Li⁺ moving toward these regions and growth from the tips of the dendrites, resulting in the branched morphology exhibited at higher applied electric fields [Fig. 5(b)].

Further experimental validation of the ability of the eBV-SPH model to predict the dendrite morphology compares the simulation to experimental studies on nanochannels at the anode interface. Sharon et al.⁵ studied the effect that polymeric nanostructures have on the morphology of silver dendrites. Their results show that as the nanochannels decrease in size, the dendrites grow in a linear fashion normal to the anode. It was also observed that as the nanochannels increase in size, the morphology of the growth becomes more branched (Fig. 6).

To compare the behavior of the eBV-SPH model to the experimental results, three cases were run with different width nanochannels. As in the experiments, the channel walls were nonreactive and the channel widths were 0.07, 0.25, and 0.5 µm. The dendrite growth for all three cases is presented in Fig. 7. At very small channel widths [Fig. 7(a)], the Li⁺ is guided to move normal to the anode surface, causing linear growth. When the channel size is large, as shown in Fig. 7(c), the Li⁺ is free to deposit in all directions, resulting in a branched dendrite morphology. The predictions of the eBV-SPH model compare well with the experimental images of Sharon et al.⁵ 

The eBV-SPH model was used to investigate the effect that modifications to the anode surface have on the dendrite growth and morphology. Various groups have studied different types of surface modifications including mechanical modifications to the anode surface.^12,14 The theory behind mechanically modifying the anode surface is that there is a natural roughness on the anode that results in dendritic growth. Despite the attempts of some researchers,¹² this inherent roughness cannot be completely eliminated, but only minimized. Instead of trying to smooth the surface, we can control the surface geometry and use these modifications to guide the deposition of Li.¹⁴ 

In the work of Kim et al.,¹⁴ the anode surface was roll pressed. Then, using a square pillar polymer mold, the surface was modified to have square cavities (Fig. 8). The experiments consisted of three different anodes, a control, a 6 µm square pattern, and a 40 µm square pattern. The control shows normal dendritic growth, while the 6 and 40 µm square patterns are intended to guide the Li deposition. The smaller holes show growth on the corners and tops of the anode, whereas the larger cavities show growth on the bottom surface and then on the top surface of the anode (Fig. 8).

To model a patterned anode surface, various sized square cavities were tested using the eBV-SPH model. Potential deformations to the SEI layer were neglected, and a uniform reaction rate across the entire anode surface was used. The depth of the square imprints remained constant at 3 µm, and the widths varied from 1.5 to 9 µm, as shown in Fig. 9.¹⁴ A pattern size of 6 µm was chosen as the starting point and then adjusted to observe how adjusting the width would change dendrite growth. The reaction rate across the anode surface was set to 1.0 × 10⁻³ μm/s. In addition, the overpotential will remain constant for all simulations, and thus, the widths of the channels will not have a direct impact on the eBV reaction equation [Eq. (10)], allowing the guiding capabilities of the channels to be analyzed.

Looking at Fig. 9, when the reaction rate is uniform across the entire anode surface, there is a significant growth on the top corners of the cavities as expected from the experimental study (Fig. 8). In Fig. 9(a), where the width of the square cavities is the smallest, the Li⁺ concentration is quickly depleted, limiting the possible dendrite growth on the bottom surface. As the width increases, Li⁺ more easily move toward the bottom surface, leading to a slower depletion of Li⁺. This allows for the guided deposition of Li on the bottom surface instead of the top.

In the simulations of Fig. 9, it is assumed that there are no heterogeneities in the SEI layer along the anode surface, and a uniform reaction rate was used along the entire anode surface. However, when mechanical modifications are made, it can create deformations in the native SEI layer, such as cracks, resulting in the exposed Li metal and an increase in the local reaction rate. To model the effect that the imprints have on the SEI layer, different reaction rates were assigned across the anode. The varying reaction rates have a direct impact on the eBV equation [Eq. (10)], thus increasing the reaction flux. On the bottom surface of the cavities, the reaction rate was significantly larger than on the tops and side walls, representing the exposed Li metal. Following the initial nucleation of Li on the surface, the new solid particles adopt a constant reaction rate regardless of their location. This new reaction rate will be the higher reaction rate in the ratio to resemble the exposed Li. This assumption is based on the work of Cohen et al.,³⁸ showing that as Li deposits, the SEI layer is unable to accommodate the volume changes that occur, causing it to crack. Preferential deposition occurs at these cracks due to the exposed Li, thus leading to dendrite formation. The width of the holes and the ratio of the reaction rates were varied, as shown in Table II, to understand the interplay between changes to the surface geometry and surface chemistry.

For the cases with a smaller ratio of high to low reaction rate (A–F), growth was seen on the sides and top of the cavities. This growth can be seen throughout the entire cavity (Fig. 10, cases A–F). The growth seen in these cases compared well with the experimental results by Kim et al.,¹⁴ showing dendrite growth that is larger at points of high convexity and dendrite growth along the tops, sides, and bottom of the holes (Fig. 8). Like the uniform reaction rate cases (Fig. 9), the imprints with a wider width allow for more Li⁺ to move into the holes and deposit on the bottom surface.

In both the low and high reaction ratio cases, similar dendrite growth was seen at the bottom of the cavities. For the low reaction ratio cases (A–F), there is a substantial growth on the bottom surface of the wide width channel (Fig. 10, case A) and a very little growth on the bottom of the narrow channels (Fig. 10, case F). However, for the high reaction ratio cases, similar growth is seen on the bottom surface for all cases (Fig. 10, cases G–I). This occurs because for the low reaction ratio cases, the dendrite growth at the top corners and sides of the patterns quickly hinders the migration of Li⁺ into the cavities, and thus, the concentration is depleted, and the dendrites at these regions stop growing. At a high reaction ratio (Fig. 10, cases G–I), Li⁺ has more time to move into the cavities before the dendrites at the corners block their path, leading to a higher concentration of Li⁺ inside the cavities and subsequently more dendrite growth on the bottom surface. In addition, limited growth on the side walls is observed (Fig. 10, cases G–I). This observation does not align with the expected growth that is seen experimentally, indicating that the experimental reaction ratio does not vary by as much as the high reaction ratio used in these simulations (Fig. 10, cases G–I).

Along with the reaction rate ratio, the size of the square patterns was adjusted. For both cases of high and low ratios, it was observed that as the pattern widths decreased the dendrite growth at the top corners and sides of the holes caused the migration of Li⁺ into the cavities to be hindered, causing a decrease in the concentration within the cavities. Eventually, the Li⁺ concentration within the hole will be depleted, and the dendrites on the bottom surface will stop growing. This results in preferential deposition on the top of the anode and has a negative impact on the lifetime of the battery. For the wider widths, depletion is seen less severely, and the Li⁺ concentration toward the bottom of the holes remains higher, allowing for dendrites to continue to grow at the bottom of the cavities. In Fig. 7, the comparison between different sized nanochannels that are nonreactive were studied. Here, it was observed that when the channels consisted of a nonreactive material, the smaller width channel suppressed the dendrite growth. However, in cases where mechanical modifications to the anode were made and the sides and tops of the patterns were reactive, the narrow channel widths did not suppress the dendrite growth (Fig. 10). Specifically, preferential deposition on the top surfaces was seen. Thus, it can be concluded that when the modifications are made to the anode itself and the channel walls are reactive, large widths are more successful in the guided Li deposition, but when the channels consist of a nonreactive material, smaller channel sizes are better.

Based on the computational modeling, in conjunction with the experimental work by Kim et al.¹⁴ and Sharon et al.⁵ when the channel walls are reactive, a pattern width of greater than 9 µm has the largest potential to suppress dendrite growth. The larger the width is, the more guided the deposition is, as seen from the experimental results when the width is increased to 40 µm (Fig. 8). In addition, when the channel walls are nonreactive, channel widths smaller than 0.10 µm are most effective in guiding Li⁺, resulting in more uniform growth.

In this paper, a novel SPH model that details the microscopic and mesoscopic phenomena at the electrode–electrolyte interface was presented with the goal of helping to design better interfaces. An extended Butler–Volmer equation is implemented to describe the rate of reaction at this interface expanding on previous modeling of dendrite growth in Li batteries.^{4,17,24,25,32} The eBV-SPH model was verified by a comparison to the standard Butler–Volmer equation. In addition, the eBV-SPH model was compared against the experimental results, looking at the effect that the current density has on the morphology of Li dendrites³⁶ and comparing how different sized nanochannels guide the deposition of Li.⁵ In both cases, the eBV-SPH model compared well to the experimental results, supporting its capabilities to accurately simulate the dendritic growth and morphology.

The focus of the eBV-SPH model is on the electrode–electrolyte interface and the physics at this location and allows tracking of important phenomena at the interface including the overpotential, which is a proposed factor in increased nucleation and growth of dendrites. In this work, it is shown that as the electric field is increased, dendrite growth correspondingly escalates. At higher electric fields, and, in turn, higher overpotentials, more dendritic structures and growth were observed. While at low electric fields, more uniform deposition and less growth were observed.

The Li metal anodes have an inherent surface roughness that can promote nucleation. Patterning the anode surface holds potential for better control over where the nucleation is taking place. Using various sized square patterns, it was shown that the patterns can direct the location of dendrite growth and that the larger hole widths promote dendrite growth on the bottom surfaces of the cavities. While significant growth occurred at points of high convexity, such as the corners of the imprints, due to a larger charge density at these locations than at other locations along the anode,^15,16 larger widths allowed for Li⁺ to migrate into the hole and continue depositing in these regions. Growth on the top corners of the narrow channels led to the depletion of cation concentration and additional growth on the top surface, which is detrimental to the health of the battery.

Finally, the effect that deformations in the SEI layer have on guided Li deposition was studied in order to understand the interplay of the anode geometry and surface chemistries. Two different reaction rate ratios were used. One ratio implementing a higher reaction rate along the bottom surface of the anodes resulted in dendrites initially growing at the bottom of the hole in the anode before growing on the sides and corners of the imprints. However, the results from this study suggested that the experimental ratio was smaller than the one being used. Using the smaller ratio, it was observed that dendrites were formed on all walls of the cavities. At the corners, the dendrites grew larger than elsewhere and blocked the transport of Li⁺ to the bottom surface in the narrow channel cases. As the channel width increased, we observe a higher concentration of Li⁺ at the bottom of the holes, allowing for more of a guided deposition. From this study, it can be concluded that when mechanical modifications are made to the anode surface and the cavity walls are reactive, wider geometries are more successful at guiding Li deposition. However, when nonreactive channels are used, the narrower the width, the more it suppresses the dendrite growth. Narrower channels will not work under reactive conditions because the dendrite growth along the walls cuts off access to the channel, thus, leading to more growth along the top of the anode. Under nonreactive conditions, this is not an issue, thus, narrower channels help to direct Li⁺ transport, promoting a more even deposition at the base of the channels.

Financial support for this research was provided by the National Science Foundation through Award Nos. 1911698, 1727316, and 2034154, the Boston University Undergraduate Research Opportunity Program, and the Hariri Institute for Computing and Computational Science and Engineering at Boston University through their Focused Research Program.

The authors have no conflicts to disclose.

All data were obtained via the open source LAMMPS code base modified for our system. The data that support the findings of this study are available from the corresponding author upon reasonable request.