Deep learning assisted modeling of Li-ion transport in SEI: a graph neural network based study

Understanding Li-ion transport through the solid electrolyte interphase (SEI) is essential for improving the stability of lithium-metal batteries, as nonuniform ionic transport through the inorganic SEI can promote spatially localized Li deposition and dendrite formation. In this work, we develop a deep-learning-assisted framework to model Li-ion transport across the inorganic SEI by combining density functional theory (nudged elastic band) calculations with graph neural network learning. A systematic diffusion dataset was first generated for eight major inorganic SEI components, namely, LiF, LiCl, LiBr, LiI, Li 2 O, Li 2 S, Li 3 N, and Li 2 CO 3, covering both bulk (grain) diffusion and grain-boundary diffusion over representative low-energy surfaces and interfaces. This dataset includes homogeneous and heterogeneous interfaces, providing a unified design space for Li-ion migration in SEI environments. A path-aware graph variational autoencoder (GVAE) was then used to learn latent representations of NEB trajectories, and the learned embeddings were incorporated one into a predictive GNN model for minimum-energy-path and migration-barrier estimation. The combined GVAE-GNN framework achieved strong predictive performance for both grain and grain-boundary diffusion, with test-set R 2 values of 0.93 and 0.94, respectively. Feature-importance and latent-space analyses further show that migration behavior is governed not only by composition, but also by local structural factors such as exposed surface, saddle-point character, and grain-boundary energetics. The results reveal clear transport trends across SEI chemistries, with halide-rich systems generally associated with narrower low-barrier distributions, while Li 3 N, Li 2 CO 3, and structurally mismatched heterogeneous interfaces exhibit broader and higher-barrier landscapes. This study give a physics-informed machine-learning framework for mapping and predicting Li-ion migration in complex SEI structures and provides insight into how SEI chemistry and microstructure jointly control interfacial ion transport.

Lithium-metal batteries (LMBs) are widely regarded as one of the most promising next-generation energy-storage technologies because they can deliver significantly higher energy densities than conventional lithium-ion batteries. Lithium metal possesses both the highest theoretical specific capacity (3860 mAh g − 1) and the lowest electrochemical potential (− 3.04 V vs. SHE) among known anode materials. Despite these advantages, the practical deployment of LMBs is severely limited by the instability of the lithium/electrolyte interface. Lithium as an electrode lies at the most reducing end of the electrochemical series and is therefore highly reactive toward electrolyte components. As a result, spontaneous electrolyte decomposition occurs during cycling, leading to the formation of a Solid Electrolyte Interphase (SEI), a chemically and structurally complex passive layer on the electrode surface (Peled, 1979; Aurbach et al., 1999). Ideally, the SEI should electronically insulate the anode while allowing fast Li-ion transport, thereby suppressing continued electrolyte decomposition. n practice, however, the SEI formed on Li metal is often heterogeneous and mechanically fragile, and it distributes Li ions nonuniformly across the electrode surface, undermining electrochemical performance. This results in uneven Li stripping and plating, low Coulombic efficiency, continuous interfacial side reactions, and eventual capacity decay. More importantly, spatial inhomogeneity in the SEI promotes nonuniform current distribution and localized Li deposition, which can trigger uneven growth of Li-metal, interfacial failure, and leads to cell short circuit. Although lithium electrode is an ideal anode from a theoretical standpoint, the unstable nature of the SEI remains one of the principal barriers to the practical realization of LMBs. Therefore, understanding and controlling Li-ion transport through the SEI—especially through its inorganic sublayer—is critical for enabling stable LMBs.

SEI is a thin, complex, multilayered structure, typically several hundred angstroms thick, that commonly forms on the electrode surface during the initial charging and discharging cycles. It originates from electrochemical reduction and decomposition of the electrolyte at the electrode surface. After its initial formation, the SEI continues its evolution chemically and structurally during battery operation mainly due to the active Li-ion transport, electron tunneling, and ongoing interfacial reactions induced by electrochemical cycling. Owing to this dynamic nature, the SEI is generally described as a heterogeneous, mosaic-like bi-layered interphase (Wang et al., 2018).

Despite the highly disordered nature of practical SEI layers, crystalline inorganic phases and partially ordered nano-crystals are frequently observed experimentally and are commonly used as atomistic model systems in first-principles investigations of SEI transport behavior (Das Goswami et al., 2023). In particular, inorganic SEI components such as LiF, Li 2 O, Li 2 CO 3, Li 2 S, and Li 3 N often exhibit locally ordered domains that can strongly influence Li-ion transport pathways and interfacial stability.

In LMBs, SEI is generally modeled as a bi-layer mosaic architecture consisting of an inorganic-rich inner layer adjacent to the electrode and a more porous, heterogeneous outer layer in contact with the electrolyte (Shi et al., 2012). The inner layer is typically enriched in inorganic species such as lithium halides, oxides, sulfides, and other salt-derived compounds, whereas the outer layer contains a larger proportion of organic decomposition products. In lithium-metal systems, this inner layer is often particularly LiF-rich, with other commonly reported components including LI 2 CO 3, Li 2 O, and related anion-derived inorganic salts such as Li 3 N and Li 2 S. The SEI formed on graphite anodes in conventional lithium-ion batteries has also been reported to contain similar inorganic species, including LiF, Li 2 CO 3, Li 2 O, and other lithium salt–derived species, within its inorganic sublayer (Kulathuvayal and Su, 2023; Beheshti et al., 2022; Adenusi et al., 2023).

Morphologically, the inorganic SEI sublayer can be described as a densely packed grains composed of the aforementioned inorganic components, with numerous grain boundaries between them. An illustration of SEI is given in Figure 1. Earlier microscopic characterization studies have revealed this grain-boundary-rich architecture of the SEI (Golozar et al., 2019). This model of SEI is further supported by previous computational and experimental investigations, including the ‘two-layer–two-mechanism’ diffusion framework proposed by Shi et al. (2012) and related studies by Lu et al. (2014), both of which corroborate the mosaic, multilayered, grain boundary rich nature of the SEI.

The failure of the SEI is generally classified into three major categories: thermal, chemical, and mechanical degradation. Thermal and chemical degradation are the dominant failure modes for the SEI formed on graphite anodes in conventional lithium-ion batteries, whereas mechanical failure is more prevalent in the SEI formed on silicon and lithium-metal electrodes (Ding et al., 2021). Since the present study focuses on lithium-metal systems, particular attention is given to the mechanical failure of the SEI associated with lithium dendritic growth. In general, dendrite formation proceeds through three stages: nucleation, growth, and penetration. Among these, the nucleation stage is especially critical, as achieving spatially uniform nucleation is essential to prevent dendrites from becoming formed (Ding et al., 2021). One important factor that can trigger dendrite nucleation is the accumulation of low-energy Li-ion diffusion channels at specific regions of the electrode surface (Liu et al., 2026; Sun et al., 2013). This phenomenon with the repeated cycling amplifies deposition inhomogeneities, leading to interfacial instability, and ultimately the growth of filamentary lithium dendrites. Therefore, uniform Li-ion diffusion through the grains and grain boundaries of the inorganic SEI is essential for stable lithium deposition and long-term battery operation.

Earlier computational studies on dendrite formation and SEI behavior have largely focused on mapping Li-ion transport in a limited set of representative SEI components. Most of these investigations rely on density functional theory (DFT), particularly nudged elastic band (NEB) calculations, to determine Li-ion migration pathways and activation barriers. A notable study by Ramasubramanian et al. (2019) examined LiF, Li 2 O, and Li 2 CO 3 and demonstrated Li-ion transport through the grain boundaries of these phases. Their results showed that grain boundaries create multiple open diffusion channels that facilitate faster Li-ion transport, with the highest diffusion rate observed for the heterogeneous LiF/ Li 2 O grain boundary compared with the corresponding homogeneous LiF/LiF and Li 2 O / Li 2 O boundaries. However this work include a few selected diffusion pathways and limited surface terminations for a given grain boundary. Similarly, Zheng et al. (2021) investigated several ionic migration mechanisms in bulk LiF, Li 2 CO 3, LiOH, and Li 2 O. While this work attributed the relatively fast Li-ion transport in LiF compared to other components, the broader transport picture within the SEI remains unaddressed, as the SEI is inherently an amalgam of crystalline grains and grain boundaries. Soto et al. (2018) reported, in their NEB study of Li-ion diffusion mechanisms in the grains of LiF and Li 2 CO 3, that among the mechanisms considered—including knock-off, concerted exchange, vacancy diffusion, and direct hopping—vacancy diffusion and knock-off are the most energetically favorable in these two components. Existing studies of Li-ion transport in SEI materials largely focus on one or a few components particularly about LiF, Li 2 O and Li 2 CO 3, with relatively limited attention to the halide-rich SEI chemistry observed in modern electrolytes (Shi et al., 2012; Kulathuvayal and Su, 2023; Iddir and Curtiss, 2010; Chen et al., 2026). Also the studies incorporating grain boundary rich SEI model concentrate only on a limited number of migration pathways within selected grains or grain boundaries. This demands the need of an collective SEI modeling that include grain and grain boundary diffusion in SEI through all reported SEI components and their possible grain boundary configurations.

As a first phase of this study, we consider the major inorganic SEI components LiF, LiCl, LiBr, LiI, Li 2 O, Li 2 S, Li 3 N, and Li 2 CO 3. Li-ion diffusion pathways in these materials were mapped using NEB calculations. As grain formation in the SEI is strongly influenced by the exposed surface and its crystallographic orientation, grain-diffusion calculations were performed over a broad range of low to medium energy surfaces identified on the basis of surface energy. For grain-boundary diffusion, both homogeneous grain boundaries and heterogeneous interfaces were included using the listed SEI components. These grain-boundary pairs were constructed on the basis of grain-boundary energy and Coincidence Site Lattice (CSL) theory. For both grains and grain boundaries, all possible symmetry-unique Li-ion diffusion pathways were identified and considered for NEB calculations, with more details provided in the Methodology section. This procedure resulted in a large dataset describing Li-ion migration across the SEI, including diffusion coordinates, crystallographic orientations, migration energetics, surface energies, and grain-boundary energies. To the best of our knowledge, this is the first study to systematically map Li-ion migration across such a broad set of major inorganic SEI phases, their low-medium energy surface orientations, and representative grain-boundary configurations within a unified SEI design space.

In the second stage of this work, the Li-ion diffusion information generated from the DFT–NEB calculations is used to develop a graph neural network (GNN)-assisted predictive model for estimating Li-ion migration pathways in a given grain or grain-boundary environment of SEI components. This step is motivated by the high computational cost of DFT-NEB calculations, which makes exhaustive screening of the large number of possible SEI chemistries, surface orientations, and grain-boundary configurations impractical. For instance, in LiF, three stable surfaces are reported in literature with surface energy ranging from 20–50 meV/Å 2 (Ramasubramanian et al., 2019; Jand and Kaghazchi, 2014). This alone can results in three distinct grains of LiF with respective growth surfaces. Considering the symmetrically distinct diffusion sites in this grain which further demands respective number of times NEB calculations required to completely map the potential energy terrain of Li-ion intra-grain diffusion. As the slab thickness used in this work was intentionally kept small to enable high-throughput DFT-NEB calculations across a large number of SEI components and interfaces, the corresponding intra-grain diffusion should be interpreted mainly as surface or near-surface diffusion rather than fully bulk-like diffusion in large supercells. Therefore, a balance between computational cost and dataset diversity was adopted to generate diffusion data for multiple SEI chemistry, surfaces, and interface configurations. Although thicker slabs generally provide more bulk-like environments and quantitatively better-converged migration barriers, the thinner slabs used in this work were intentionally adopted to enable high-throughput generation of a diverse diffusion dataset spanning numerous SEI components, surfaces, and interface configurations.

Given the SEI components and their diverse local atomic distributions, crystallographic orientation, and interfacial environment involved in the diffusion dataset, crystal graph data is extensively used throughout the modeling in this work. A GNNs based modeling approach is therefore particularly well suited to this modeling, as they can naturally encode atomic connectivity, local coordination, and geometric features that govern migration energetics, unlike traditional machine-learning models that rely primarily on fixed feature vectors and often struggle to capture such structural complexity. With recent advances in machine learning and its applications to materials science, state-of-the-art increasingly focused on the use of graph neural networks (GNNs) for the prediction of structural and material properties, the conditional generation of novel materials and molecules, and the design of new materials tailored to desired target properties (Abramson et al., 2024; Jiao et al., 2023; Watson et al., 2023; Schneuing et al., 2024; Banik et al., 2023; Dwivedi et al., 2023).

When its comes to SEI and machine learning perspectives, a few works can be exemplified. One earlier work by Soleymanibrojeni et al. (2024) used a kinetic Monte Carlo method to model SEI formation as a function of reaction barriers. It built an initial database using a design-of-experiments approach to train a Gaussian process classification model, which was later improved through active learning. However, a key limitation of this approach is that the migration barriers are not learned or predicted directly from atomistic structure; instead, they are treated as predefined coarse-grained input parameters whose ranges are approximated from first-principles calculations and literature data. In a similar work Li et al. (2025) developed active-learning ML interatomic potentials trained on DFT data from amorphous SEI structures to predict Li transport and migration barriers in Li 2 CO 3 and Li 2 EDC which covers only a limited subset of SEI chemistries, and remains computationally expensive for long-timescale simulations when using the more accurate neural-network potential. In general, there have been several attempts to use machine-learning-assisted simulations and predictions for battery interfaces, including SEI formation, lithium dendrite growth, interfacial ion transport, and atomistic modeling with machine-learned interatomic potentials (Sun et al., 2025). However, these studies still rely largely on limited, chemistry-specific datasets—often focused on only one or two major SEI components—and their accuracy and transferability for migration-barrier prediction remain constrained by the limited availability of high-quality DFT- or simulation-derived training data. Based on these observations, we note that a GNN-based modeling approach has not yet been applied specifically to understanding Li-ion diffusion in the SEI, nor have DFT–NEB datasets spanning a wide range of SEI components been systematically integrated with machine learning to establish design rules for SEI architecture. In this scenario, the present study, which provides a large NEB dataset encompassing both grain and grain-boundary diffusion together with a GNN-based modeling approach, establishes a rational design space for SEI structures by the prediction of Li-ion migration behavior across a broad range of SEI components.

The remainder of this article is organized as follows. Section presents the methodology, beginning with data generation using DFT-NEB, with separate discussions of grain and grain-boundary model construction. It then describes the graph construction scheme, design of the graph variational autoencoder and graph-level predictive model, and training and validation schemes. section presents the results and discussion, including the general trends in the DFT-NEB migration-barrier dataset, followed by the GNN-based predictive modeling results for grain diffusion and grain-boundary diffusion. Finally, section concludes the article.

This study consists of two phases: (i) generation of Li-ion DFT-NEB migration-barrier data for selected SEI components, namely, LiF, LiCl, LiBr, LiI, Li 2 O, Li 2 S, Li 3 N, and Li 2 CO 3 and (ii) development of GNN-based predictive models trained using DFT-NEB data. The methodology as a whole is illustrated as a flowchart as shown in Figure 2.

The structures of the SEI components were obtained from the Materials Project repository using the MPRester API, with the lowest-energy-above-hull entries selected for each component (Jain et al., 2013). These bulk structures were subsequently geometrically optimized through DFT calculations. DFT calculations were performed using the Vienna Ab-initio Simulation Package (VASP) (Kresse and Furthmüller, 1996). The projector augmented-wave method was used to describe core–valence interactions, together with the Perdew–Burke–Ernzerhof (PBE) exchange–correlation functional within the generalized gradient approximation (GGA) (Blöchl, 1994; Perdew et al., 2008; Perdew et al., 1996). A plane-wave kinetic-energy cutoff of 500 eV was applied throughout. Geometrical optimizations were continued until the forces on all atoms were below 0.1 eV/Å. Brillouin-zone integrations employed Monkhorst–Pack k -point meshes of 7 × 7 × 7 for bulk and 3 × 3 × 1 for slab initial optimization calculations. The structural and energy information after the optimization step is given in Supplementary Table S1.

Diffusion data generation

This study considered Li-ion vacancy diffusion to estimate migration barriers, following prior reports that vacancy-mediated diffusion is the most favorable mechanism in SEI components compared with concerted exchange and knock-off mechanisms (Shi et al., 2012; Soto et al., 2018; Yildirim et al., 2015). The MinimumDistanceNN() class in Pymatgen was used to identify candidate diffusion sites and pathways in the grain structures while excluding symmetrically equivalent sites. Initial and final configurations for each diffusion event were then constructed and geometrically relaxed before NEB calculations. Intermediate images representing the Li-ion migration coordinates were then generated using the VTST scripts maintained by Henkelman Research Group (SCRIPTS, 2025). For each diffusion process, five intermediate images were introduced between the initial and final states to reduce the overall computational cost of the DFT calculations, with the images obtained through linear interpolation along the migration path. All images were then optimized simultaneously to determine the minimum energy path (MEP). Migration barriers were calculated using the NEB method as implemented in VASP, with the climbing-image NEB (CI-NEB) scheme employed to locate the saddle point accurately (Henkelman et al., 2000). The migration barrier was defined as the maximum energy along the MEP relative to that of the initial state.

Grain diffusion

The geometry-optimized bulk structures were used to construct slab models for migration-barrier calculations in grains. For each slab, a vacuum layer of 10 Å was introduced along the z-direction, while periodicity was preserved along the x and y directions. Slab thicknesses of 5–7.5 Å were chosen to balance physical accuracy with the computational cost associated with the multiple DFT calculations required for each NEB simulation. This vacuum-slab configuration was employed in the DFT-NEB calculations to reduce spurious interactions between adjacent images. Grain-diffusion slabs were subsequently generated using Pymatgen (Ong et al., 2013) by selecting exposed surfaces with surface energies below 250 meV/Å 2, as lower-energy surfaces are generally more stable and more likely to occur naturally in grain structures (Ramasubramanian et al., 2019). The surface energy, (γ surf), was calculated using Equation 1: where E slab is the total energy of the relaxed slab, E bulk is the bulk energy per unit structure, N is the number of units, and S is the surface area of the slab.

The resulting structures of grain, together with their corresponding surface terminations and surface energies, are summarized in Table 1. Diffusion sites were then identified from these grain models using nearest-neighbor analysis to determine plausible Li migration pathways. When multiple diffusion sites were symmetrically equivalent, only one representative site was retained to avoid redundancy in the dataset. On average, 10–15 NEB calculations were performed for each grain structure, enabling systematic mapping of the diffusion-energy landscape.

Grain boundary diffusion

Homogeneous grain boundaries: Grain-boundary migration-barrier data were generated using homophase grain-boundary models constructed on the basis of Coincidence Site Lattice (CSL) theory. In this approach, the degree of lattice coincidence between two adjoining grains is quantified by the CSL parameter, defined as the ratio of the area of the coincidence site lattice to that of the fundamental lattice. To sample representative low- Σ boundaries, CSL values of 3, 5, 7, 9, 11, and 13 were considered. The corresponding rotation axes, rotation angles, grain-boundary planes, number of grain-boundary planes, and interfacial spacings used in this work are summarized in Table 2.

The grain-boundary structures were generated using the open-source Python package aimsgb, which enables construction of periodic grain-boundary models for both cubic and non-cubic crystals. For each SEI component, grain-boundary supercells were built from the parent crystal by joining grains A and B according to the prescribed CSL relationship. Default aimsgb settings were used for model generation, except 10 Å vacuum thickness, 1.75 Å interfacial separation beyond default construction, primitive-cell reduction enabled, and an atomic-plane tolerance factor of 0.25 Å. A few atomic layers were removed along the ± z directions to reduce the system size to approximately 100–150 atoms, thereby keeping the models computationally less demanding. The resulting grain-boundary slabs were constructed with thicknesses in the range of 5–10 Å. A few sample grain boundary structures constructed in such a way is illustrated in Figure 3.

To exclude energetically unfavorable grain boundaries, the grain-boundary energy of each model was evaluated using Equation 2: where γ GB is the grain-boundary energy, E AB GB is the total energy of the relaxed grain-boundary structure, E A and E B are the bulk energies of structures A and B, N A and N B are the numbers of unit cells of A and B, respectively, and S is the interfacial area. Grain boundaries with energies greater than 1,000 mJ/m 2 were excluded from the dataset, since such high-energy interfaces are generally much less stable, less likely to form or persist under realistic conditions, and may therefore introduce physically less relevant diffusion pathways into the analysis.

To focus the diffusion calculations on the grain-boundary region, selective dynamics were employed to constrain atoms far from the interface, while atoms in and near the boundary region were allowed to relax. Candidate Li migration pathways were then identified by selecting diffusion sites across the grain boundary on the basis of local coordination analysis using Pymatgen (Ong et al., 2013). In particular, the destination site for Li diffusion was chosen as a neighboring Li-vacancy site in the vicinity of the interface. Symmetrically equivalent initial and final diffusion sites were excluded to avoid redundant NEB calculations. The resulting grain-boundary diffusion pathways were subsequently used in climbing image-NEB calculations to determine the corresponding migration barriers.

Heterogeneous interfaces: In addition to homogeneous grain boundaries, heterogeneous interface models were also incorporated into the training dataset to account for interfacial diffusion across distinct SEI components. As classical CSL theory is not directly applicable to hetero-interfaces, these structures were constructed using a coherent lattice-matching approach based on the Zur–McGill algorithm (Zur and McGill, 1984), which identifies commensurate in-plane supercells between two dissimilar crystals under prescribed mismatch and angular tolerances. We used CoherentInterfaceBuilder() class from Pymatgen to build these systems. To maintain computational feasibility, only interface models with small supercell areas, average lattice mismatch below 5%, and atom counts compatible with DFT calculations were retained. After lattice matching, the film was rotated and strained to match the substrate, and the interface was assembled with an appropriate interfacial gap (≈ 1.5 Å) and vacuum spacing (10 Å). Only a limited number of such heterogeneous interfaces were included, as many candidate interfaces were found to be structurally less stable and associated with relatively high grain-boundary energies. A few representative heterogenous interfaces are given in Supplementary Figure S1 As in the homogeneous grain-boundary models, diffusion pathways were selected across the interface by excluding symmetrically equivalent sites and choosing nearby Li-vacancy destinations for NEB calculations.

Training data

The training dataset was generated from two categories of diffusion data for SEI components: (i) grain diffusion and (ii) grain-boundary diffusion. For both categories, the target property is the Li-ion migration barrier, determined from the MEP of the diffusion process connecting the initial and final Li diffusion coordinates. The set of energies corresponding MEP of images in the NEB calculation, hereafter referred to as the profile energy along the MEP, was also included as a target, given that the migration barrier corresponds to the maximum energy value within this set.

The overview of grain diffusion training data is given in Figure 4, where the first subplot depicts the perpendicular nudged elastic force along the MEP, showing that the force distributions are generally centered near zero for most SEI components, which suggests reasonably well-converged NEB pathways. However, Li 3 N and Li 2 CO 3 exhibit noticeably broader spreads and larger extrema, indicating greater variability in local migration environments. The second subplot presents the profile energies along the diffusion path, which remain relatively narrow and close to zero for most components, whereas Li 3 N shows a substantially wider distribution with pronounced negative excursions, showing a more diverse set of intermediate-state energetics. The third subplot summarizes the migration barriers, one of the target data, revealing that LiCl, LiBr, and LiI are associated with comparatively lower and more compact barrier distributions, while Li 2 O, Li 2 S, and especially Li 3 N display broader distributions and higher upper ranges, indicative of more high-energy diffusion pathways. LiF and Li 2 CO 3 also show moderate-to-high barriers with appreciable spread. And, the last subplot shows the number of NEB runs used for model training and testing. Total number of NEB data generated for grain diffusion is 280 where LiF contributes the largest number of samples, while the remaining components maintain reasonably balanced train–test splits. Overall, the dataset captures both comparatively simple and highly variable grain-diffusion behaviors across the SEI components, providing a diverse basis for model training.

The grain-boundary diffusion dataset consists of 198 homogeneous grain boundary and 104 heterogeneous interfaces NEB calculations. Figure 5a summarizes these data on a per-component basis. For homogeneous grain boundaries, each case is assigned to the corresponding SEI component since both grains belong to the same material. For heterogeneous interfaces, component-wise involvement counting was used. For example, if LiF forms grain boundaries with both Li 2 S and Li 2 O, its involvement is counted twice, once for the LiF/ Li 2 S boundary and once for the LiF/ Li 2 O boundary. Hence, the bars corresponding to the heterogeneous interfaces in Figure 5a represent the frequency with which a given SEI component participates in heterogeneous interfaces. For homogeneous grain-boundary cases, Figure 5b further resolves the data according to the CSL value of the interface, classified by the corresponding Σ values. This subplot therefore illustrates how the homogeneous grain-boundary dataset is being distributed across different boundary types for each SEI component.

The heterogeneous grain-boundary dataset is summarized in Figure 5c as a pairwise component matrix, where each entry denotes the number of heterogeneous grain-boundary configurations formed between two SEI components. The planes involved such grain boundaries are given in Supplementary Table S2.

Together, the two types of grain boundaries contribute a total of 2,114 data points as each NEB calculation consists of seven intermediate structures. Figure 5d shows the train–test split of this dataset for each SEI component. The labels above the bars indicate the corresponding numbers of training and test samples, which highlights the component-wise fraction of data used for training and testing. Figure 6 depicts the distribution of (a) NEB force along band, (b) profile energy, (c) migration barrier and (d) grain boundary energy in homogeneous grain boundaries.

Featurization: The final dataset contains both the scalar migration barrier and the profile energies along the diffusion pathway for each NEB calculation. To describe grain and grain-boundary diffusion environments, two categories of descriptors were considered: (i) composition-based chemical features and (ii) structure-based features.

For composition-based descriptors, we employed Magpie-preset featurization implemented in the Matminer library (Ward et al., 2018). Specifically, the compositional feature set was generated using a MultipleFeaturizer combining Stoichiometry(), ElementProperty(), ValenceOrbital(), and IonProperty(). These descriptors capture complementary chemical information, including elemental fractions, statistical elemental-property attributes, average valence-orbital characteristics, and ionic-property-related features. The initial feature filtering was performed using Pearson’s correlation coefficient; specifically, when two features exhibited a correlation greater than 90%, the feature with the larger mean absolute deviation was retained. For predictive modeling in materials science, the use of fundamental physicochemical properties of a system as input features has proven promising in previous studies (Alghadeer et al., 2024; Kulathuvayal et al., 2024; Chan et al., 2022; Jung et al., 2024). The detailed information of the composition features used in this study are tabulated in Supplementary Table S3. The space-group information included in the feature space was determined from the local crystal structures associated with the NEB images and diffusion environments rather than from composition-based estimations. In addition to these composition-derived descriptors, several structural descriptors were incorporated in train data. These include the exposed surface of the grain slab model (or the grain-boundary surface for grain-boundary systems), surface energy (or grain boundary energy for grain-boundary slabs) and a binary indicator specifying whether a given NEB image corresponds to the saddle-point configuration. Importantly, slab models associated with the NEB pathway, including the initial, intermediate, and final images, were represented as graphs to preserve the local atomic connectivity and geometric environment during diffusion. More information about graph construction is given in section. Together, these composition and structure-based descriptors provide a unified representation of the chemical identity and atomistic diffusion environment of each SEI system.

It should be noted that some composition-derived descriptors, such as melting temperature and space-group number, are not interpreted here as direct causal variables governing Li-ion transport. Rather, these descriptors act as proxy representations of broader structural and bonding characteristics of the SEI components. In contrast, local structural descriptors such as exposed surface orientation, grain-boundary energetics, and saddle-point environment are more directly connected to the underlying diffusion physics because they influence Li-site connectivity, local coordination environment, bottleneck geometry, and interfacial structural disorder along the migration pathway.

Feature-dependence analysis was also performed to identify the most important features among the original descriptors that are directly correlated with the MEP. For this, permutation importance was used. In this approach, the trained model kept fixed, and one feature is shuffled across the test samples while the prediction error for the MEP is recalculated. A large increase in error indicates that the shuffled feature is important for predicting the target property, MEP. This method is particularly suitable here because it measures the model’s actual predictive reliance on each feature and remains valid even for nonlinear models.

Graph neural network

Crystal graph construction: The encoding process begins with the construction of a periodic crystal graph from a POSCAR and CONTCAR structures (VASP input/output files) corresponding to NEB images. Each structure is represented as a graph, G = (V, E), where atoms correspond to nodes, V = { v i } i = 1 N, and edges, E ⊂ V × V, represent interatomic interactions within a cutoff radius. The featurizer extracts fractional atomic coordinates, r i fractional ∈ R 3, atom types, Z i, lattice parameters (lengths a and angles α), and periodic boundary information. The resulting PyTorch Geometric object contains node features, edge indices, and periodic image offsets, thereby enabling graph construction under periodic boundary conditions. This study uses ‘Directional Message Passing Neural Network’ (DimNet++) to generate graph embedding of grain and grain boundary structures (Gasteiger et al., 2020a; Gasteiger et al., 2020b). According to review articles published in recent years on GNN applications in materials science and drug design, DimeNet has been widely used in these domains (Al-Maeeni et al., 2024; Reiser et al., 2022; Zhang et al., 2025). Additional information about its architecture and message-passing scheme, including the fractional-to-Cartesian coordinate transformation, periodic distance computation, triplet construction for angular interactions, radial and spherical basis expansions, edge embedding, edge-to-node projection, angular interaction blocks, accumulated node representations, and graph-level readout, are briefly provided in the Supplementary Material.

Predictive model: A multi layered GNN model is used to predict diffusion energy profiles and migration barriers directly from initial NEB structures. Figure 7a depicts the architecture of the model. Training data, the NEB trajectory consists of seven structures, including the two end point geometries, along the reaction coordinate, and each structures, as mentioned earlier, is represented as a periodic crystal graph. Individual structures are encoded using a GNN encoder using DimeNet++, producing embeddings g i ∈ R 256. Additionally, a pretrained graph variational autoencoder (GVAE) was used to incorporate SEI component distribution information into each NEB data point, thereby enabling each data point to retain awareness of its position within the SEI component space. A summery of GVAE used is provided in section.

The role of the GVAE pretraining in the present framework is primarily to learn embedded latent representations of the crystal graph structures and their configurational relationships (from the latent space) within the diverse SEI diffusion dataset. Although the GVAE is pretraining includes NEB trajectory information along with structure information, the resulting latent embeddings are intended to capture generalized structural characteristics of diffusion environments rather than directly encode explicit migration barriers or NEB energy profiles. During inference, the predictive GNN model operates on graph-based structural representations together with the pretrained latent embeddings, without directly requiring NEB energies, force profiles, or converged minimum-energy pathway information as input features.

To incorporate diffusion pathway information, each structural embedding is augmented with projected scalar physical descriptors (f i) namely, exposed surface information, surface or grain boundary energy, saddle-point indicator, force along the band, DFT energy per atom and compositional features and a learnable positional embedding (p i) specifying the position of intermediate structure along the diffusion path. Here, NEB associated features are included mainly in GVAE encoder model, whereas GNN-predictive model only interoperates structural and chemical descriptive. Thus the resulting per-image representation (e i) is given as Equation 3:

which is then grouped into a sequence E = [ e 0, …, e 6 ] and processed using a multi-layer Transformer encoder. Self-attention across images enables the model to capture energy curvature, saddle localization, and collective structural evolution along the diffusion pathway.

From the contextualized embeddings, the model predicts per-image profile energies, y ̂ i, which are shifted to enforce the VTST (Henkelman et al., 2000) convention.

The migration barrier is then derived directly from the predicted profile. Depending on configuration, the barrier is computed as a differentiable softmax pooling as in Equation 4:

where parameter τ controls smoothness. The model is trained using a multi-task objective combining profile energy regression, barrier regression, and a consistency regularizer enforcing agreement between the predicted barrier and the maximum predicted profile energy. This physics-informed formulation ensures that the learned representation captures both local structural manifold and global diffusion energetics within a unified pathway-level model.

Graph variational auto encoder

A path-aware GVAE design logic is used in this study to learn structural latent representations of NEB data points. Unlike conventional crystal VAEs that treat configurations independently, VAE used here encodes all seven structures of an NEB trajectory jointly as a single data point, allowing the model to learn collective geometric information along the reaction coordinate and positional awareness in SEI component space. This design enables the learning of diffusion processes in individual SEI components while preserving an intuitive understanding of the SEI as a whole.

The GVAE architecture illustrated in Figure 7b, is adapted from the Crystal Diffusion Variational Autoencoder developed by Xie et al. (2021). In contrast to a conventional reconstruction-based VAE, this framework employs a denoiser-based decoder, which is trained to recover the original input from slightly perturbed representations rather than simply reconstructing an unaltered ground-truth structure. For continuous attributes, such as atomic coordinates, a perturbation is introduced through Gaussian noise. For discrete attributes, such as atom types, direct addition of continuous noise is not meaningful; therefore, the training procedure instead perturbs the data by interpolating between the true and predicted probability distributions. Such a denoising objective improves the strength of the learned latent space and captures chemically and structurally relevant variations. Detailed information about training scheme using GVAE is given in supplementary material.

Training and validation

Grain diffusion: The GVAE model developed in this work serves as a supportive module for the predictive GNN model. Its main function is to learn a latent representation for each NEB calculation in a higher-dimensional space. This latent representation encodes the relative position of each NEB calculation, and therefore each MEP, within the broader configurational space defined by all SEI components and their corresponding diffusion pathways.

The GVAE was tuned by varying the training batch size, latent-space dimension, and feature/target normalization strategy as shown in Supplementary Figures S1a–c. Based on this analysis, the best overall performance was obtained with a batch size of 16, a latent dimension of 128, and a learning rate of 0.004. Since the GVAE used here follows a semi-latent-diffusion VAE scheme, the effect of latent diffusion strength was also examined by varying the perturbation applied to atomic coordinates and lattice parameters; these results are also summarized in Supplementary Figure S4.

The performance of the fine-tuned GVAE was further analyzed by varying the KL-divergence weight parameter, β, from 0.1 to 1.0, as shown in Figure 8a. The parameter β controls the extent of latent-space regularization in the GVAE model (Burgess et al., 2018). In the present case, this regularization determines how the NEB data corresponding to different SEI components are distributed in the learned latent space. The corresponding KL-divergence convergence behavior is shown in Figure 8b. Among the tested values, β = 0.6 provides the best balance between regularization and reconstruction quality. A visualization of this latent space in 3-D with t-Distributed Stochastic Neighbor Embedding (t-SNE) based dimensionality reduction technique is shown in Figure 11a. A clear clustering of SEI components are visible in this basis with LiF forms comparatively bigger hull among the components. Furthermore, as shown in Figure 8c, this GVAE gives the most accurate reconstruction of lattice parameters, atomic coordinates, and atom types. Representative grain structures reconstructed by the GVAE are shown in Supplementary Figure S5, demonstrating that the model reproduces the lattice parameters, atomic coordinates, and atom types with good accuracy.

After training the GVAE with β = 0.6, the learned latent embeddings of the NEB structures were incorporated into the predictive GNN model. To evaluate the contribution of these embeddings to the predictive performance, the GNN was tested under three different input schemes: (1) using the original descriptors (graph representations of the NEB images and the composition-based chemical descriptors) together with the GVAE-derived latent embeddings, (2) using only the original descriptors, and (3) using the GVAE-derived latent embeddings. In the first setting, original descriptors were concatenated with the latent embeddings learned by the GVAE. In the second setting, each NEB calculation was represented only by its original descriptors, without latent-space information. In the third setting, the GNN was trained using only the GVAE-derived latent representation. The performance of these three models was compared using the test R 2 score and mean absolute error (MAE) for MEP prediction, as shown in Supplementary Figures S3b–d.

Among the three training schemes, the GNN model that incorporated both the original descriptors and the GVAE-derived latent embeddings achieved the best predictive performance, with a test-set R 2 score of 0.93 and an MSE of 1.02 × 1 0 − 2, as shown in Figure 8d. This result indicates that the latent representation learned by the GVAE contains additional physically meaningful information from SEI component space about the configurational relationships among diffusion pathways that is not fully captured by the original descriptors alone. Therefore, this model was selected for further benchmarking.

To assess its predictive reliability across the SEI components, representative test cases were selected from all components considered in this study with their low-index exposed surfaces. The corresponding minimum energy paths predicted by the GNN are compared with those obtained directly from NEB calculations in Figures 9a–h. Overall, the predicted MEPs reproduce the main features such as saddle points of the NEB calculated profiles with good agreement across the different grain systems, demonstrating that the combined GVAE-GNN framework is able to capture the Li-ion diffusion pathways in grain diffusion in selected SEI components.

To validate the model a separate sample set of grain models constructed for all SEI components described in this study. This set of data never been used in any part of training or testing. The validations results are provided in Supplementary Figure S9.

Feature-importance analysis was further carried out using the trained GNN model, and the corresponding importance scores are presented in Figure 10. Among the composition-derived descriptors shown in Figure 10a, melting temperature, space-group number, and ground-state volume emerge as statistically important features. However, these quantities should primarily be interpreted as proxy descriptors reflecting broader structural and bonding characteristics of the SEI components rather than direct causal variables governing Li-ion transport. For example, higher melting temperatures generally correlate with stronger ionic bonding and reduced lattice flexibility, which can suppress local structural relaxation during Li migration and thereby increase activation barriers. Similarly, the space-group symmetry influences the connectivity and anisotropy of Li-ion hopping pathways through its control over neighboring Li-site arrangement and diffusion-channel topology.

More physically direct insight is obtained from the structural descriptors shown in Figure 10b. The exposed surface orientation strongly affects the local Li coordination environment, surface relaxation behavior, and migration bottleneck geometry, all of which influence the smoothness of the migration-energy landscape. In particular, the relatively higher importance of the (210) and (111) surfaces suggests that local crystallographic topology plays a significant role in determining accessible low-energy diffusion channels. Likewise, the strong contribution of saddle-point and surface-energy descriptors indicates that Li-ion migration in SEI grains is governed not only by composition, but also by the local energetic roughness and structural continuity of the diffusion pathway. These trends are further consistent with the relatively smooth and narrow migration-barrier distributions observed for lithium halides compared with Li 3 N and Li 2 CO 3.

Latent-space analysis further indicates that the GVAE learns chemically distinct regions within the generative manifold. As shown by the PCA and t-SNE based projections in Supplementary Figures S6a,b, the latent representations exhibit clear component-dependent clustering for the major SEI phases, demonstrating that the model captures meaningful differences in their underlying structural and Li-ion diffusion characteristics. To visualize cluster separation, the original 256-dimensional latent vectors were projected onto a three-dimensional Cartesian basis using t-SNE, as shown in Figure 11a. This reduced representation clearly highlights both the separation and the partial overlap among the various SEI components, reflecting similarities and differences in their learned latent descriptors. Sampling points for migration barrier generated from each component-specific hull through the barrier head produced characteristic distributions for each SEI component. These data distributions are given as violin plots in Figures 11b,c which further show that halides such as LiCl and LiI are associated with relatively narrow, low-barrier distributions, whereas Li 3 N and Li 2 CO 3 exhibit substantially broader and higher-barrier distributions, with Li 2 O and Li 2 S lying in between. Overall, these results suggest that the GVAE latent space is not only structured by composition, but also encodes physically relevant variations in ionic migration behavior across different SEI components.

Grain boundary diffusion: Following the procedure adopted for grain diffusion, the GVAE model for grain boundary diffusion was optimized through hyperparameter tuning for latent space dimension, batch quantity, learning rate, and weight decay, as shown in Supplementary Figures S7a–d. Among the tested configurations, a batch size of 8, a latent dimension of 128, and a learning rate of 0.001 yielded the best overall performance of GVAE model. Evolution of total test loss of the tuned model for KL-weight, β = 0.2, 0.4, 0.6, 0.8 and 1.0, shown in Figure 12a, indicates stable convergence. The corresponding KL-divergence also shown in Figure 12b show that β = 0.4 provides the best balance between latent-space regularization and reconstruction fidelity. Then the reconstruction quality of GVAE model was assessed using the R 2 scores of reconstructed lattice parameters, atomic positions, and atomic types, as presented in Figure 12c. Overall, the model reconstructs the grain boundary slabs with high fidelity, achieving an average R 2 score of approximately 0.85. A representative set of samples from distinct grain boundaries are given in Supplementary Figures S3a–h for a visual comparison of reconstruction efficiency.

Using the optimized GVAE latent representation, the predictive GNN model was subsequently trained for grain boundary diffusion. As in the grain diffusion case, three training schemes with distinct feature sets were considered. In addition to the NEB descriptors used in the case of grain diffusion, the feature space for grain boundary diffusion was expanded to include grain boundary type (homogeneous or heterogeneous), grain boundary energy, and the exposed surface of the grain boundary. These additional descriptors were introduced to better capture the grain boundary configurations in such systems.

The GNN model for grain-boundary diffusion showed the best performance when the composition-derived features were combined with the latent embeddings learned by the GVAE model. A comparison of the three training settings is presented in Supplementary Figures S9a–c. The best-performing model Figure 12d achieved an R 2 score of 0.94 with an MAE of 1.21 × 1 0 − 2. This improvement indicates that the GVAE-derived latent embeddings transfer useful information about the distribution of different grain-boundary environments in the latent space, enabling the GNN model to predict the MEP more accurately compared to other two training approaches.

Using the best-performing GNN model, feature-importance analysis was further carried out, and the results are shown in Figure 13. Similar to grain diffusion, composition-derived descriptors such as average melting temperature and space-group number appear among the statistically important features. However, these descriptors should again be interpreted mainly as indirect representations of the underlying structural and bonding characteristics of the SEI phases rather than direct mechanistic descriptors of diffusion itself.

More physically meaningful trends emerge from the structural descriptors shown in Figure 13b. The dominant importance of saddle-point environment, grain-boundary energy, CSL Σ, and exposed grain-boundary surfaces indicates that Li-ion transport across interfaces is governed primarily by local interfacial structure and energetic disorder. Grain-boundary energy is particularly important because it reflects the degree of lattice mismatch and structural distortion at the interface. Higher-energy interfaces generally exhibit more disrupted Li-site connectivity and a more heterogeneous local energy landscape, which increases migration barriers.

Similarly, the influence of CSL Σ suggests that crystallographic coherence between neighboring grains plays a key role in maintaining continuous low-energy diffusion pathways. Lower- Σ grain boundaries tend to preserve better atomic registry across the interface, whereas higher- Σ boundaries introduce more disordered migration environments. The importance of exposed surface orientation further indicates that the local coordination geometry at the interface strongly influences Li-ion mobility. These observations collectively suggest that interfacial transport in SEI is controlled not only by chemistry, but also by the structural compatibility and energetic stability of neighboring grains forming the interface.

Using the optimized GVAE model, the latent-space distribution was projected into three dimensions using t-SNE, as shown in Figure 14a. The visualization reveals tightly packed clusters corresponding to various grain boundaries, with only negligible overlap, indicating that the model effectively captures the diversity of grain-boundary environments. To further analyze migration-barrier trends, each cluster was enclosed within a hull and sampled within the latent space. The resulting migration-barrier distributions generated through the GVAE decoder are shown in Figure 14b. These distributions indicate that heterogeneous interfaces generally exhibit higher migration barriers than homogeneous grain boundaries. Among the homogeneous grain boundaries considered, LiF–LiF and Li 2 O - Li 2 O show the lowest barrier ranges, approximately 0.58–1.20 eV and 0.20–0.86 eV, respectively. In contrast, among the heterogeneous interfaces, Li 2 O - Li 2 CO 3 gives a comparatively higher barrier range of about 1.05–1.65 eV. Furthermore, grain boundaries involving lithium halides tend to produce migration barriers in the range of 0.5–1.5 eV. These results suggest that chemical and structural mismatch at heterogeneous interfaces creates a more resistive and energetically heterogeneous landscape for Li-ion transport than that found in homogeneous grain-boundary systems. Such mismatch disrupts the continuity of low-energy Li diffusion channels and introduces localized variations in coordination environment and interfacial energetics.

In this work, we developed a unified DFT-NEB and GNN-based framework to describe Li -ion transport across both grains and grain boundaries formed in SEI. Widely reported components in SEI, namely, LiF, LiCl, LiBr, LiI, Li 2 O, Li 2 S, Li 3 N, and Li 2 CO 3 are considered in this study. The generated diffusion dataset spans 280 grain-diffusion NEB calculations and grain-boundary dataset comprising 198 homogeneous and 104 heterogeneous interfaces, corresponding to a total of 302 diffusion data points. Each NEB diffusion data point consists of seven intermediate structures. The predictive framework achieved high accuracy for both transport regimes, with test R 2 = 0.93 for grain diffusion and R 2 = 0.94 for grain-boundary diffusion, showing that the combination of graph-based structural encoding and GVAE latent representations can reliably capture the minimum-energy-path landscape across a chemically and structurally diverse SEI design space.

Beyond predictive accuracy, the present results provide quantitative design insights into how SEI chemistry and microstructure jointly govern Li-ion migration. For grain diffusion, the latent-space sampling results show that lithium halides generally define the lowest-barrier and most spatially uniform transport landscape. In particular, LiCl and LiI exhibit mean migration barriers of only ≈ 0.12 and ≈ 0.11 eV, respectively, with narrow standard deviations of 0.08 and 0.10 eV. Li 2 O also remains comparatively favorable, with μ ≈ 0.14 eV and σ ≈ 0.10 eV. Further, Li 3 N and Li 2 CO 3 show much broader and higher barrier distributions, with Li 3 N reaching μ ≈ 0.93 eV and σ ≈ 0.69 eV, and Li 2 CO 3 yielding μ ≈ 0.50 eV and σ ≈ 0.40 eV. These differences indicate that grain diffusion in the SEI is not controlled solely by composition, but also by the roughness of the underlying migration-energy landscape: halides and Li 2 O provide both lower average barriers and better transport uniformity, whereas Li 3 N and Li 2 CO 3 introduce stronger local variability in Li -ion mobility. This is consistent with the observed feature-importance trends, where the dominant composition-based descriptors are melting point temperature, space-group number, and ground-state volume, while the saddle-point indicator and exposed surface orientation are the most influential structural factors. In the grain diffusion model, the saddle-point flag alone carries an importance score of about 0.03, far exceeding other structural descriptors, and among the exposed surfaces the (210) and (111) planes are more influential than (001), which contributes almost negligibly.

The grain-boundary results further show that interfacial transport is systematically more restrictive than intra-grain transport, and that the transport penalty depends strongly on interface type and boundary energetics. Even the most favorable sampled homogeneous grain boundaries, such as Li 2 O - Li 2 O and LiF - LiF, have mean migration barriers of 0.62 and 0.71 eV, which are already higher than the best grain-diffusion barriers. More resistive interfaces, including Li 2 CO 3 - Li 2 CO 3, Li 2 O - Li 2 CO 3, Li 2 O - Li 2 S, Li 2 S - Li 2 S, Li 3 N - Li 3 N, and LiCl - LiCl, cluster around 1.08–1.28 eV, with the largest spreads observed for Li 3 N - Li 3 N (σ ≈ 0.33 eV) and Li 2 S - Li 2 S (σ ≈ 0.25 eV). These results indicate that grain boundaries can act as major bottlenecks for ionic transport even when the parent phases themselves are not intrinsically poor conductors. The feature-importance analysis reinforces this interpretation: in the grain-boundary model, the saddle-point indicator has the largest structural importance score, followed by grain-boundary energy, grain-boundary type, DFT energy per atom, and CSL Σ. Thus, unlike grain diffusion, grain-boundary transport depends explicitly on interfacial stability and crystallographic matching, not merely on bulk material identity. This distinction highlights that favorable SEI transport cannot be predicted solely from bulk chemical composition, but must also account for local interfacial topology and structural coherence between neighboring phases.

Another important outcome of this study is that lower-energy and lower-complexity interfaces tend to support more favorable Li-ion transport. The results indicate that lower CSL Σ boundaries are generally associated with lower migration barriers, whereas higher Σ interfaces tend to produce a more rugged barrier landscape. Similarly, lower-index grain-boundary planes tend to be more favorable, especially in the lithium halides. Among the important grain-boundary surfaces identified, (1 ̄ 11), (110), (1 1 ̄ 2), (0 3 ̄ 2), and (001) are also consistent with surfaces reported in experimentally observed stable SEI configurations in earlier studies. In addition, the present results show that heterogeneous interfaces are generally more resistive than homogeneous ones. This trend is especially significant for halide-carbonate combinations, for which the barrier distributions can exceed 1.50 eV, indicating that chemical mismatch across the interface can severely disrupt low-energy Li-ion pathways.

Despite the SEI design scope of the present framework, a few limitations of this model should be noted. First, the diffusion analysis is restricted to vacancy-mediated Li-ion transport, while other mechanisms such as interstitial and knock-off diffusion were not considered. Second, the study focuses only on crystalline inorganic SEI components and therefore does not capture the role of organic species, mixed inorganic-organic interfaces that are often present in SEI structures particularly found close to electrolyte. In addition, the calculations and model predictions are based on static DFT-NEB energetics and do not explicitly include temperature-dependent effects, dynamic structural disorder, or long-time interfacial evolution during electrochemical cycling. These factors should be incorporated in future work to extend the present structure-transport framework toward a more complete description of realistic SEI behavior.

Taken together, these results suggest a physical picture for SEI transport design. Fast and spatially uniform Li-ion conduction is most likely to arise in SEI morphologies enriched in halide and oxide grains with relatively low-barrier surfaces, while high fractions of Li 3 N, Li 2 CO 3, and energetically unfavorable heterogeneous interfaces are expected to increase transport heterogeneity. Since localized low-energy pathways can concentrate ionic flux and thereby promote nonuniform lithium deposition, the microstructural distribution of grains and grain boundaries is expected to be just as important as the SEI chemistry itself for suppressing dendrite nucleation and growth. In this sense, the present framework does not merely provide a predictive model for migration barriers; it establishes a quantitative structure-transport map for the inorganic SEI, showing that SEI performance should be engineered through control of phase selection, exposed surface orientation, grain-boundary character, and interface stability.

A limitation of the present work is that the SEI was modeled using crystalline slabs and idealized interfaces, whereas experimentally formed SEI layers are often amorphous, nanocrystalline, and dynamically evolving (Wang et al., 2018). This approximation was adopted to enable systematic high-throughput DFT-NEB calculations across a broad range of inorganic SEI components and interfacial environments within feasible computational cost. Although such ordered models cannot fully capture the structural complexity of real SEI layers, they nevertheless provide physically meaningful insight into how local coordination environment, surface orientation, crystallographic matching, and interfacial energetics influence Li-ion migration behavior. Future work extending the present framework, toward amorphous and dynamically evolving SEI structures using larger-scale simulations and machine-learning interatomic potentials would further improve the realism of the model.