Density Functional Theory screening of gas-treatment strategies for stabilization of high energy-density lithium metal anodes

To explore the potential of molecular gas treatment of freshly cut lithium foils in non-electrolyte based passivation of high energy-density Li anodes, density functional theory (DFT) has been used to study the decomposition of molecular gases on metallic lithium surfaces. By combining DFT geometry optimization and Molecular Dynamics, the effects of atmospheric (N2, O2, CO2) and hazardous (F2, SO2) gas decomposition on Li(bcc) (100), (110), and (111) surfaces on relative surface energies, work functions, and emerging electronic and elastic properties are investigated. The simulations suggest that exposure to different molecular gases can be used to induce and control reconstructions of the metal Li surface and substantial changes (up to over 1 eV) in the work function of the passivated system. Contrary to the other considered gases, which form metallic adlayers, SO2 treatment emerges as the most effective in creating an insulating passivation layer for dosages<= 1 mono-layer. The substantial Li->adsorbate charge transfer and adlayer relaxation produce marked elastic stiffening of the interface, with the smallest change shown by nitrogen-treated adlayers.


Introduction
The increasing demand for stable, high energy density rechargeable batteries for long-range electric vehicles motivates the growing interest in developing alternative chemistry and cell strategies to replace existing Li-ion insertion-based technologies [1][2]. Driven by the substantial theoretical increase in energy density, great efforts are currently being devoted to the development of Li-air and Li-sulfur batteries [3][4][5][6][7][8][9], which depend on the (to date hypothetical) availability of stable, highly reversible, lithium metal anodes, capable of delivering a nearly ten-fold increase in theoretical capacity (3,860 mAh g -1 ) over commercially used graphite anodes (372 mAh g -1 ) [2].
The highly electropositive nature of Li (-3.04 V redox-potential vs. SHE [10]) is responsible for its extreme (reducing) reactivity towards exposed molecular media. Similarly to graphite anodes, initial decomposition of the electrolyte and the ensuing formation of a protective solid-electrolyte interphase (SEI), which should be Li + -ion permeable yet mostly electronically insulating, is a beneficial process for the stabilization of Li anodes [1][2][3][4][5][6][7][8][9]. However, the repeated removal (stripping) and re-insertion (plating) of Li atoms beneath the SEI upon electrochemical cycling is known to cause cracks in the SEI and ensuing exposure of the electrolyte to metallic Li, leading to progressive electrolyte decomposition [11]. Another unresolved issue affecting the stabilization of Li metal anodes is the formation and growth, through cracks in the SEI, of highly reactive lithium metal protrusions (a.k.a. dendrites) during cycling. Lithium-dendrite growth eventually short-circuits the battery electrodes, which might cause an organic solvent electrolyte to ignite, leading to catastrophic failure of the battery [12][13][14][15][16][17][18][19][20][21][22][23][24]. Recent work indicates that besides cracking in the SEI during cycling, the presence of sub-surface impurities (nitrides and other compounds depending on the preparation/storage of the Li foil) at the Li anode can be critical for dendrite formation and growth [25].
Prevention of Li dendrites has so far focused on physical and chemical strategies to block their formation by controlling the SEI composition and morphology via use of pressure, SEI-stabilizing additives, ionic liquid-based electrolytes, as well as copolymer (solid) electrolytes and mixtures of liquid and polymeric electrolytes [14-15, 18-19, 23, 26-30]. Stabilization of Li anodes is made even more challenging by the simultaneous need for the electrolyte to decompose into electron-insulating SEI with sufficiently high diffusivity of Li + ions, allowing the flow of Li + ions to the cathode through a thermally and electrochemically stable electrolyte with a low boiling point. Although promising advances have been very recently shown to be possible via combination of carefully chosen liquid and polymer electrolytes with low reduction potential, high viscosity and large size anions [30], addition of halogenated salts (especially LiF) to the electrolyte [31], and hollow carbon nanosphere coating [32], stable cycling of Li anodes for several hundred cycles at room temperature at competitive (dis)charge rates (similar to those achieved with graphite) has not, to the best of our knowledge, been achieved yet.
The observed relationship between occurrence of cracks in the SEI and Li-dendrite formation suggests that creation of a tough (i.e. mechanical damage tolerant [33]) SEI should be beneficial in preventing dendrite formation. Recent research in damage tolerant natural and synthetic materials indicates that hierarchical multi-scale (nm to cm) structuring (extrinsic toughening [33]) of composite materials can be crucial for crack suppression [33][34]. However, the requirement of favorable Li + diffusivity through the SEI could be hardly meet by adoption of known extrinsic toughening strategies [33][34] leading to cm-thick SEI, which would exceed the thickness of commercial cells (both electrodes and electrolyte-soaked separator) by several orders of magnitude. These considerations make exploration of novel strategies towards formation of electrochemical and mechanical stable (nm-thick) SEI a necessity for viable stabilization of Li-anodes. To this end, the critical role of atomic relaxations for interface mechanical anomalies [35][36][37][38] and the expected substantial charge transfer involved in the SEI-formation call for atomistic insight into the structural and mechanical properties of Li anode SEI beyond available results from continuum models [24,29].
Apart from, to the best of our knowledge, one exception where gas (N2) pretreatment of metal Li was considered [39], the explored strategies for Li anode stabilization to date have invariably targeted formation of the SEI via decomposition of the cell electrolyte or components dissolved in it [14-15, 18-19, 23, 26-30]. Experimental work in the field has been complemented by a limited number of Density Functional Theory (DFT) studies of adatom energy and diffusion on vacuum-exposed [40] and implicitly solvated [41][42] Li surfaces, ionic liquid decomposition on defect-free Li(100) [43][44], force-field modelling of fractures in Li single crystal [45], and coarse-grained dynamic Monte Carlo studies of Li dendrites [46]. However, recent advances in the field indicate that major benefits can be achieved by pre-treating Li anodes before exposure to the cell electrolyte [39,47]. In addition, the recently established link between Li subsurface impurities and dendrite formation [25] suggests that controlled deposition of impurities in metal Li substrates could be effective in preventing dendrite formation and growth.
Qualitatively, the ideal SEI or, as we start to explore here, an alternative passivation layer created by pre-treatment of the metal Li anode, should fulfill the following conditions: i) it should be electronically insulating in order to prevent electron transfer from the Li anode to the electrolyte. ii) It should be thick enough to suppress electron tunneling from the (biased) electrode to the electrolyte.
To this end, we speculate that iii) the occurrence of a Li-SEI interface dipole opposing (zero-bias) electron-transfer from the passivated anode to the electrolyte may be beneficial. iv) The SEI should be tough [33] to adapt to the volume changes of the Li anode upon cycling (stripping during discharge and plating during charge) without cracking. Alternatively, v) a SEI capable of quickly self-healing [48][49][50] the cracks created during cycling may be also highly beneficial. vi) The SEI should allow good diffusivity of Li + ions. In this respect, nm-thick SEI (favoring Li-diffusion) may be preferable, provided the SEI is sufficiently thick to keep the Li-anode and the electrolyte electronically decoupled. Ideally, vii) it should be possible to tune a priori the Li + ion diffusivity of a given SEI to match the given cathode redox chemistry and (dis)charge rate, allowing for balanced battery assembly. viii) The SEI should be impermeable to (and insoluble in) the electrolyte solvent and other contaminants dissolved in it. Simultaneous fulfilment of all these conditions, and stability of the SEI over several (hundred to thousand) charge-discharge cycles is clearly a formidable challenge, which can be hardly met without a thorough understanding of the atomic-scale factors governing the SEI formation and evolution upon cycling.
The SEI formation via inherently out of equilibrium chemistry during the initial contact with the electrolyte and cycling of the Li anode, and the limited atomic scale control of the pristine Li surface present severe challenges to the characterization, thence understanding and eventual optimization, of the SEI formation in controlled and reproducible conditions. These considerations, encouraging results on the beneficial role of N2 treatment of Li metal anode [39], and the observed dependence of Li anode stability on the inert atmosphere (e.g., dehydrated air vs. argon) in which commercial Lifoils are made [47,51], make us wonder whether alternative gas-solid, equilibrium based, chemical strategies could be used to create a working (i.e. fulfilling conditions i-viii above) SEI or SEIprecursor layer on Li metal anodes, before contact with the liquid or polymeric electrolyte. To the best of this knowledge, this strategy has not been systematically studied, which motivates the present work.
Although the complexity and extension of the actual anode/SEI/electrolyte interfaces cycled at variable electrode bias is well beyond the current capabilities of standard Density Functional Theory (DFT) methods [52][53], previous success of DFT-based strategies in elucidating the reactivity of adsorbed molecules on Li metal surfaces [43][44]54] makes the approach convenient for exploration of some of the benefits which might be achieved by gas treatment of pristine, freshly cut Li anodes. Furthermore, the adopted computational approach allows preliminary assessment of the actual benefits of using hazardous gases (i.e. F2 and SO2, vide infra) without taking unnecessary experimental risks.
To explore the opportunities offered by molecular-gas-passivation of pristine Li metal surfaces, here we investigate 0 K and room temperature decomposition of different molecular gases on metal Li surfaces for coverages in the 0.25-1 mono-layer (ML) range. Using DFT geometry optimization and Molecular Dynamics we investigate the effects of atmospheric (N2, O2, CO2) and toxic (F2, SO2) gas decomposition on the relative energy of Li(bcc) (100), (110), and (111) surfaces, their reducing potential (approximated by the corresponding work function), and emerging electronic and elastic properties. 6 The presented results indicate that dosing of different gases can lead to passivation layers with profoundly different electronic and elastic properties. Depending on the dissociated gas, insulating or metallic adlayers, with surface elastic constants up to ten times stiffer or softer than the pristine Lisurfaces, can be obtained. We believe these results should be useful to inform future experimental efforts towards stabilization of high energy-density Li-metal anodes via gas pre-treatment of the Li electrodes.

Methods
DFT simulations were performed via the Projected Augmented Wave (PAW) method as implemented in the VASP program [55]. In all cases, the PBE exchange-correlation (XC) functional [56], a 400 eV plane wave energy cutoff, and 0.2 eV Gaussian smearing were used. At least 15 Å vacuum separation between periodic replicas of the slab models and dipole corrections were used for all the simulations. To prevent introduction of artificial dipoles perpendicular to the surfaces, molecules were adsorbed on both sides of the slabs. Geometry optimizations were performed without any atomic or symmetry constraints, with a force-convergence threshold of 0.05 eV Å -1 via the RMM-DIIS quasinewton algorithm [57]. The slab models of the Li(100), Li(110) and Li(111) surfaces were constructed using the DFT-optimized lattice constant for Li(bcc), 3.466 Å. In all cases, orthorhombic [Li(111)] grids of symmetry-irreducible k-points. The grids were chosen to maintain the same k-point spacing (≤ 0.003 Å -1 ), which was checked to yield energies converged to within less than 1 meV/atom for bulk Li(bcc) (13 symmetry irreducible k-point grid). This choice produced (100), (110) and (111) slabs with four top Li adsorption sites per exposed surface ( Figure 1).
Given their use for equilibration-only purposes, canonical (NVT) Molecular Dynamics (MD) simulations were run using the Verlet integration algorithm [58] and the Berendsen thermostat [59] as implemented in VASP. In all cases the time step was 1.5 fs. Geometry optimization and MD runs for all the molecularly decorated slabs were carried out allowing unconstrained spin-polarization in the system.
In analogy with previous studies of chemical bonding of organics at metallic surfaces [60][61], and as also discussed in [41], we found that use of van der Waals corrections (according to Grimme's parameterization [62]) negligibly affected the optimized geometry (< 0.01 Å changes in the adsorption length of the energetically favored systems).
Vibrational frequencies were calculated via symmetric finite displacements (± 0.05 Å) following further optimization (to within 0.01 eV Å -1 force-tolerance) of the selected systems with an increased plane wave energy cutoff of 600 eV. Elastic tensors accounted for ion-relaxation following the procedures described in [63] and [64], as implemented in VASP. Based on the numerical (non-zero) value of the elastic constants bound to be zero owing to the symmetry [65] (orthorhombic or tetragonal) of the cells, the error of the procedure is < 0.1 GPa for the bare Li slabs.
Bader charge analyses [66] were carried out on the basis of the total charge density i.e. accounting for both the electronic and ionic core charges.
Slab formation energies (Eform) were calculated as: where Nmol is the number of gas molecules initially present in the system and Emol is the energy of one molecule optimized in vacuo.
Work functions (W) were calculated from the difference between the vacuum-electrostatic plateau (Ev) and the computed Fermi energy (EF):

Choice of Li-surfaces and molecular gases
To explore the effects of molecular gas adsorption on Li(bcc) substrates, we considered three Li(bcc) surfaces with different surface-energy [67]. Specifically, we focused on the lowest surface-energy  [31], and speculation that dissociative adsorption of lone-pair rich systems, potentially leading to a lone-pair rich passivation layer, benefit Li + coordination and diffusivity across the layer.

Optimized molecular adlayers
For all the three considered crystallographic cuts [Li(100), Li(110), Li(111)] and gases, initial geometries were prepared for different coverage in the 0.25-1 ML range placing the undissociated adsorbate on different surface adsorption sites ( Figure 1) at distance of at least 1.8 Å from the topmost Li-atoms. For all considered gases, several different initial molecular orientations were explored, with more than 60 adsorption geometries being screened for each gas. Details on the initial geometry set up and energy screening after geometry optimization can be found in the supplementary material.
Tables S1-S6 in the supplementary material list all the considered initial geometries, together with the computed slab formation energy (Eform) after geometry optimization. In all cases we model strongly exothermic (Eform < 0) reaction between the molecular gases and the Li-surfaces, accompanied by significant rearrangement of the topmost Li-layers. Unsurprisingly, given their known large electronegativity and oxidizing chemistry [10], F2, O2 and SO2 yield the lowest Eform when reacted with Li-slabs. Reaction with N2 and CO2 is computed to be substantially less exothermic (less negative Eform). Figure 2 summarizes the computed lowest Eform for each molecular gas on Li(100), Li(110) and Li(111). The atomic structure of the lowest Eform systems for each gas is shown in Figure 3. The SI reports atomistic models of the lowest Eform system for all the considered Li surfaces. Alkali metal redox chemistry with aqueous and other reducible media is known to be extremely vigorous and potentially explosive depending on the mixing of the preliminary products [70], which is reflected in the computed very negative (< 10 eV or, equivalently, < 2.9 eV / adsorbate) Eform following dissociative adsorption of F2, O2, SO2. The substantial energy released upon dissociative adsorption of F2, O2, SO2 may be effective in promoting adsorbate induced Li-surface reconstruction of freshly cut Li-surfaces. Alternatively, co-dosing of small amount of F2, O2, SO2 during initial gas treatment of freshly cut Li anodes, and the ensuing energy release, could be used to activate and/or alter the surface chemistry with other, less reactive, molecular gases. Overall, these results clearly indicate that designing molecular-gas treatment of Li-slabs towards engineering of pristine passivation layers based on the experimental (or computed) EA of the molecular reactants could be highly misleading. Direct simulation of the reaction products turns out to be necessary for rational development of experimental gas-treatment strategies towards passivation of Li-substrates.
Inspection of the optimized geometry for the lowest Eform systems reveals dissociation for F2, O2, and SO2, molecular condensation for CO2, (formation of an adsorbed acetylenediolate, C2O2, species) and subsurface intercalation for N2, Nad, and SO2 (O-atoms). The supplementary material contains further analysis of the optimized geometry. It is worth to recall that these results have been obtained following geometry optimization, which indirectly points to the existence of barrier-less reaction and intercalation channels for the considered gases on Li surfaces (from the adopted initial geometries).
Consistent with recent DFT results on solvated alkali metal (Na) clusters [70], we find that atomic relaxation is not needed to trigger initial charge transfer at the immediate Li/adsorbate interface. This electron transfer strongly alters the potential energy surface governing the molecular and interface relaxation, leading to barrier-less reaction for all the considered adsorbates.
Whereas formation of oxide dissociation products following O2 adsorption is consistent with available XPS results for O2-treatment of metal Li films [73], the occurrence of an acetylenediolate C2O2 product (and oxide subsurface intercalation) from CO2 dissociation does not match experimental

XPS suggestions of oxalate (C2O4) intermediate formation on Li from the reaction of CO2 with metal
Li at 120 K [73][74][75][76][77][78]. While these deviations could be caused by biases in the simulations due to the limited size of the simulation cells and neglect of surface defects as well as temperature effects, we note that in [73] the Li-substrate was characterized (at 120-350 K) after substantially larger (30 Langmuir) molecular exposure then considered here, which may explain the observed differences.
Although the simulations suggest that formation of adlayers with isolated N-adatoms is energetically favored over N2 subsurface intercalation (Nad, Figs 2,3), it is interesting to note that, in spite of the substantial charge transfer (Fig. 4), intercalation of (markedly elongated to 1.34 Å, supplementary material) N2 molecules turns out to be favored over N2 dissociation (at 0 K) on defect-free substrate.
Based on the experimentally known occurrence of nitride (Li3N) contamination in N2-exposed Lifoils [25], we speculate that N2 dissociation may be triggered at surface defects (neglected in our models).  (111). The occurrence of metallic, therefore arguably conducting, Li-adsorbate reconstructions exposed to the medium suggests that larger molecular dosages (> 1 ML) are needed to grow thicker, expectedly insulating, passivation layers capable of electronically decoupling the anode and electrolyte. The SO2 case stands apart from the others since the adsorbate and topmost Li layers reveal a noticeably suppressed PDOS at EF. These results suggest that lowdosage SO2 treatment should be more effective than low-dosage F2, N2, O2 and CO2 exposure in creating ultra-thin insulating passivation layers, which may be beneficial for extremely fast Li + diffusion.
The characteristically strong reductive chemistry of metal Li (and other alkali metals) is intimately related to its high EF value, or equivalently, low work-function (W = 2.9 eV for polycrystalline metallic Li [79]) in comparison to more inert transition metals (> 4.5 eV [67,79]). Accordingly, increase of metal Li W by molecular passivation, resulting in an energetically more costly electron extraction, hence lower EF and expectedly lower reducing reactivity may be a rewarding strategy towards stabilization of Li-anodes. Figure  It is interesting to note that for higher Eform structures, which are therefore predicted to be less frequently observed, such as the lowest Eform CO2/(100) system (Figure 2), the calculated W increases by up to 2 eV with respect to the value for clean Li(100). This indicates that just by adsorption of ≤ 1 ML of different gases, and as a result of the different adlayer geometries and Li-adsorbate charge transfer, an engineered increase of metal Li W by more than 1.5 eV could be in principle possible.
Further work is in progress to investigate the evolution of the computed changes in metal Li W for larger dosages of molecular gases. These results will be reported elsewhere. shown as black continuous lines.

Optimized molecular adlayers following NVT MD equilibration
To investigate the occurrence of artefacts in the structural screening by geometry optimization of structures prepared starting from undissociated gas molecules, the lowest Eform system for each considered gas and Li-slab were subject to a short (> 1.5 ps) NVT MD equilibration at 300 K (supplementary material) followed by geometry optimization of the final MD snapshots. The optimized structures were then subject to structural and electronic characterization (Figure 2, 4-6).
We stress that rather than statistically robust insight into the real-time dynamics of the optimized adlayer (which would require larger simulation cells, longer MD trajectories, and rigorous canonical ensemble sampling [80]), the main target of this study was to use NVT MD equilibration, followed by structural relaxation, to identify lower Eform minima potentially missed by the initial screening.
This is a simple form of simulated annealing that permits the system to escape from high energy metastable configurations (local minima on the global potential energy surface).
The largest deviations in structure and Eform (> 6 eV) take place for F2/Li(100) and CO2/Li (100) (Figures 4 and 6). These changes, however, do not affect the conclusion that, apart from CO2, the computed Eform for adsorption of F2, O2, N2 and SO2 is not directly correlated with the charge transferred from the Li-slab.
Inspecting the optimized geometry for this second set of lowest Eform systems ( Figure 6) confirms dissociative adsorption for F2, O2, and SO2, as well as molecular condensation for CO2 with formation of an adsorbed acetylenediolate, C2O2 species and subsurface O-atom intercalation. Subsurface intercalation is consistently predicted (observed? predicted?) also for the energetically favored adlayer of N2 and Nad.
The computed work function (W) for these new lower Eform structures results is invariably smaller (up to more than 1 eV) values with respect to the pristine metal Li slabs ( Figure 5) for all considered gases, SO2 and CO2 included. These results confirm that that low-dosage molecular adsorption (0.25-1 ML), as considered here, is not effective in creating an interface dipole capable of lowering the slab EF, increasing its W, thence quenching the reducing activity of the metal Li substrates. To this end, larger dosage (> 1 ML) may be effective. Work in this respect is in progress and will be the subject of a forthcoming contribution.
Layer-resolved analysis of the PDOS for the lower Eform minima after NVT MD equilibration confirms that adsorption of SO2 on Li(111), leading to S-, O-and topmost Li-PDOS suppression at EF is still more effective then adsorption of F2, O2, N2 and CO2 in creation of a nearly insulating passivation layer capable of electronically decoupling the Li-subsurface from the exposed medium.
Overall, the modelled changes in structure, Eform and W for the optimized system before and after NVT MD equilibration suggest that extra care should be taken when modelling molecular adsorption on Li slabs in the absence of experimental structural input. For the considered systems, structural screening via geometry optimization alone is shown to be clearly not sufficient.

Vibrational and elastic properties of the passivated Li slabs
We now analyze the effects of the molecular dissociation on the elastic properties of Li slabs.
Although the limits of mechanical interface stabilities for Li-surface contacted to a polymer electrolyte have been previously studied via continuous elastic theory [29], the substantial rearrangements and charge transfer at the Li-adsorbate interface (Figures 3, 4, and 6), and ensuing likely change of the interface elasticity from the bulk counterparts (an aspect neglected in the continuous treatment [29]), calls for atomistic investigations of the elasticity of the passivation layers.
Surface and interface atomic relaxations are known to strongly affect the elastic properties of layered or nanostructured materials as well as grain-boundaries [35][36][37][38]. The observed correlation between Li dendrite growth and crack formation due to mechanical stress of the SEI [11][12][13][14][15][16][17][18][19][20][21][22][23][25][26][27][28][29][30] motivates our interest in the dependence of the Li-adlayer elastic properties on the composition of the dissociated gas, amount of interface re-organization, and Li-adsorbate charge transfer. We speculate that elastically compliant passivation layers with small elastic constants, leading to reduced propensity for irreversible plastic deformation and crack formation, should be beneficial for stabilization of Li-anodes. How to maximize the elastic compliance of Li-passivation layers has, to best of our knowledge, received no attention in the atomistic simulation literature, which prompts this section.

Harmonic Vibrational Frequencies
As a first approximation to the elastic properties of the passivated Li-slab, we computed the harmonic vibrational frequencies for the molecularly decorated slabs. To a first very qualitative and exploratory approximation (to be rigorously tested in the next section), one could suspect that the introduction of hard (high frequency) adlayer vibrations may lead to a stiffer (i.e. less elastically compliant) Liadsorbate interface.
For this analysis we consider the bare Li(100), (110) and (111) surfaces together with the lowest Eform system for each considered molecular gas. The computed (Γ-point) vibrational frequencies are shown (as wavenumbers) in Figure 7. The largest wavenumber vibrations for the bare Li slabs are below 340 cm -1 . In all cases, we compute a noticeable vibrational hardening, with at least 1.8-fold increase of the highest energy vibrational modes for F2/Li(100). While the computed vibrational hardening for O2, Nad and SO2 is somehow larger (× 1.9-2.5 increase), the acetylenediolate C-O and C-C stretchings for CO2/Li(110) (Figure 6) leads to a much larger (× 6.5) increase. Although the adopted approximated DFT functional (PBE [56]) is well known to generally yield underestimated vibrational frequencies [81], we expect the relative molecule-induced vibrational hardening to be qualitatively correct.
Using the vibrational hardness of the adlayer as an approximated measure of its elastic compliance, it is tempting to link the smallest computed vibrational hardening of F2/(100) with recent reports on the beneficial role that addition of LiF to the electrolyte, and ensuing formation of LiF-rich SEI, plays in the stabilization of Li-anodes [31]. Along the same line, and owing to the ×6.5 increase in the adlayer vibrations and expected reduction in the elastic compliance of the passivated Li-slab, the simulations suggest that initial low dosage CO2 treatment of Li substrates could be a very detrimental choice, which should be accordingly avoided. In the next section we test these results against explicit evaluation of the slab elastic constants.

Figure 7. Computed vibrational (Γ-point) wavenumbers (cm -1 ) for the bare Li slabs and the lowest
Eform systems (Figures 2 and 6) of each considered molecular gas after NVT MD equilibration.

Elastic constants
Investigation of the elastic properties of the passivated Li slab is extended by explicit evaluation of the surface elastic constants for the lowest Eform systems. We are particularly interested in the effects which adsorption of different molecular gases cause on the elastic compliance of the passivated slabs.
From a specialist perspective, and aiming at reducing the computational cost of further molecular screening, we are also interested in testing the (un)suitability of using the adlayer vibrational hardness as an approximated measure of its elastic compliance.
According to linear elasticity theory, the elastic constants of a system (collected in the stiffness tensor C) govern the proportionality between the stress (tensor, σ) generated in an isotropic material and the applied strain (tensor, ε) [ there exist five independent in-plane elastic constants Cijkl (Cxxxx, Cyyyy, Cxyxy, Cxxyy, and Cyyxx in extended Cartesian notation [37][38]65]). The number of independent in-plane elastic constants reduces to three (Cxxxx = Cyyyy, Cxxyy = Cyyxx, Cxyxy [37][38]65]) for the tetragonal Li(100) slab (4 mm plane point group [82]). Figure 8 displays the independent in-plane elastic constants for the bare Li slabs and the lowest Eform system for each considered gas.
In stark contrast to the results of the harmonic vibrations analysis, which suggests the smallest vibrational hardening for F2/Li(100) (Figure 7), the computed elastic constants reveal that F2 adsorption actually causes the largest elastic stiffening, with a nearly ten-fold increase of Cxxxx (= Cyyyy) and Cxyxy. Substantial and strongly anisotropic elastic stiffening is modelled also for SO2/Li(111), with increase in the elastic constants ranging from 20% (Cyyyy) to 65% (Cxxyy and Cyyxx). As shown in Figure 8, with the exception of F2/Li(100), the computed elastic stiffening (increase in Cijkl) and softening (decrease of Cijkl) do not correlate directly with the computed slab formation energy (Eform), which rules out also approximation to the adlayer elastic stiffening on the basis of computed formation energy or charge transfer (Figure 4).
Although it is tempting to link the smallest computed elastic stiffening of Nad/Li(110) with the measured increased cycling efficiency of Li metal anodes passivated with N2 gas-solid treatment [39], one must be cautious: The elastic constants of the composite Li-passivation layer will inevitably evolve with increased thickness (larger molecular gas dosage) of the adlayer. Accordingly, the results obtained for 0.25-1 ML coverage should not be taken as representative of a nm-thick SEI as present in reality [39]. Further work will focus on the study of the dependence of the passivation-layer elasticity on its thickness and structure.
To summarize this section, analysis of the elastic constants for the molecularly decorated slabs, strengthen previous conclusions on the complex role of atomic relaxations (and, based on this work, charge transfer) for the elastic properties of (Li-adsorbate) interfaces [35][36][37][38]. Perhaps unsurprisingly, the adlayer vibrational hardness, formation energy and Li→adsorbate charge transfer are found not to directly correlate with the adsorbate-induced elastic stiffening/softening of the slabs, suggesting that explicit evaluation of the elastic constants of a given passivation layer cannot be avoided to quantify its elastic compliance. viii) The introduction of high vibrational frequencies, strongly exothermic slab formation energies and large Li→adsorbate charge transfer was found not to directly correlate with the computed adlayer elastic stiffening.

S1.1. Surface Energies
For the clean Li-slabs, surface energies (Esurface) were calculated as: where Eslab is the energy of the optimized slab made up of NLi Li-atoms, ELi-bulk is the energy (per atom) of bulk Li(bcc) at the optimized lattice constant and A is the slab surface area. The factor 2 accounts for the occurrence of two relaxed surfaces in the slab.

S1.3. Initial Geometries
All investigated gases (N2, O2, CO2, F2, SO2) and single N-atoms (Nad) were placed vertically and horizontally on the different adsorptions sites on both sides of the Li(100), Li(110) and Li(111) slabs shown in Figure S3 for several coverages in the 0.25-1 ML range. Table S1-S6 reports details of the considered initial adsorption geometries. We recall that each slab models contained four topmost Liatoms per exposed surface (i.e. overall 8 topmost Li-atoms per slab). Where applicable, the molecules were also rotated around the z-axis (= slab normal), so that the inplane projection of the molecular axis was aligned with the different vectors shown in Figure S3. The initial closest distance between Li and adsorbate atoms was always 2.0 Å, except for O2 (1.8 Å).            Table S7.