Theoretical evidence of self-intercalated 2D materials for battery and electrocatalytic applications

Geometric structures and stabilities of the ic-2D materials

The ic-2D crystals consist of two layers of TMDs and one layer of transition metal intercalants, which are native. The MX₂ monolayer has a sandwiched structure, and the M layer is sandwiched between the two X layers, as shown in Figure 1A. The weak vdW bonding between the two TMD layers allows M atoms to insert and occupy the gap and locate at the octahedral vacant sites. Distinct stoichiometric phases are, therefore, produced depending on the different levels of intercalation. In this way, the chemical formula can be generally expressed as M_σM₂X₄ for ic-2D crystals, where σ represents the concentration of intercalated ions. In this work, we focus on seven different 1T phase 3d transition metal sulfides with σ = 33% in ic-2Ds: M₇S₁₂ (M = Ti, V, Cr, Mn, Fe, Co, and Ni). For M₇S₁₂, the vertical distance between the top and bottom sulfur layers is denoted as height h, while the bond lengths of M-S bonds are denoted as d₁-d₃ [Figure 1A], and the results for all the ic-2Ds are summarized in Supplementary Table 1. To screen out ic-2Ds with high stability for potential applications as Li-ion anodes and HER electrocatalysts, we utilize the screening procedure described in Figure 1B: first, thermodynamic, mechanical, and thermal stabilities of ic-2Ds are checked, followed by investigating their electrode and electrocatalytic performance.

Figure 1. (A) Representative geometric structure of the M₇S₁₂ ic-2Ds (σ = 33%) from the top and the side views. M atoms are denoted in blue, and S atoms are denoted in orange for the first layer and yellow for the second layer in each MS₂ layer for better visual distinction, respectively. (B) Procedure to screen potential ic-2D electrode materials and HER electrocatalysts.

Before studying the specific properties of a certain 2D material, stability is always the first thing we should consider. Herein, the thermodynamic, mechanical, and thermal stabilities of ic-2Ds are evaluated. First, negative cohesive energy (E_coh) is used as the condition to judge whether 2D materials are thermodynamically stable. From our results, all seven ic-2D systems can satisfy the thermodynamic criteria for stability, lying in the range between -4.08 and -6.45 eV/atom [Supplementary Table 2], indicating that they are feasible in experiments. Second, the mechanical stability is examined by calculating their 2D elastic constants by DFT. Using the Born-Huang criteria: C₁₁C₂₂ - C₁₂² > 0 and C₆₆ > 0 [Supplementary Table 3]^[34], we find that all the systems can satisfy the criteria and thus are mechanically stable. Finally, the thermal stability of the material is proved by performing a series of AIMD simulations at 300 K with a time step of 2 fs. After performing AIMD simulations for 10 ps, no significant structural collapse or deformation can be found, and the total energy and temperature of the system tend to oscillate around the equilibrium values, indicating thermal stability [Supplementary Figure 1]^[35-37].

The comprehensive assessment of stability discussed above demonstrates that all the seven systems are feasible for experimental synthesis (in fact, three of them have already been fabricated, which are V₇S₁₂, Cr₇S₁₂, and Fe₇S₁₂)^[22], and the optimized structures of these ic-2D systems are shown in Supplementary Figure 2. Due to the insertion of native transition metal atoms, ic-2Ds exhibit pseudo-2D structures, and some of them also exhibit anisotropic characters, such as V₇S₁₂, Cr₇S₁₂, and Mn₇S₁₂[Supplementary Table 1], which may originate from the lowering of symmetry due to their self-intercalation nature. Typically, the seven ic-2Ds have lattice constants (a, b), ranging from 5.5 to 6.0 Å, heights (h) between 7.2 and 8.7 Å, and metal-sulfur bond lengths (d) between 2.2 and 2.5 Å.

To further understand the origin of the high stability of ic-2Ds, differential charge density (DCD) and electron localization functions (ELF) are performed to investigate the bonding properties of ic-2Ds. DCD (∆ρ) is defined as $$\Delta \rho=\rho_{M_{7} s_{12}}-\rho_{M_{6} s_{12}}-\rho_{M}$$, and the DCD plots in Figure 2A and Supplementary Figure 3 indicate significant charge transfer between the M and the surrounding S atoms. Taking Co₇S₁₂ as an example [Figure 2A], the difference in local bonding arrangement can induce strong charge transfer from the cobalt atom to the nearby sulfur atom (for Co atoms in both the pristine CoS₂ layer and the intercalated layer). The ELF plots in Figure 2B further demonstrate that most electrons are located around the sulfur atoms, which reveals the ionic character of cobalt-sulfur bonds in Co₇S₁₂. The significant ionic bond contributes to its robust structural stability, and these phenomena can also be observed in other ic-2D systems [Supplementary Figures 3 and 4].

Figure 2. Bonding properties of Co₇S₁₂. (A) Differential charge density distributions of Co₇S₁₂, where the isosurface value is set as 0.012 e/Bohr³. Yellow and green colors represent charge accumulation and depletion, respectively. (B) ELF maps of Co₇S₁₂. The values of ELF > 0.5, ELF = 0.5, and ELF < 0.5 represent the covalent bond, metallic bond, and ionic bond, respectively. Co atoms and S atoms are denoted in purple and orange, respectively.

Electronic and magnetic properties of ic-2Ds

For the TMD materials, the electronic structure has been proven to be strongly dependent on the number of d-electrons and the coordination environment of the transition metal atoms^[38-40]. For example, 2H-MoS₂ exhibits semiconducting nature with a band gap of about 1.3 eV, and this can be attributed to the symmetry-induced splitting of the Mo 4d orbitals with the fully occupied Mo 4$$d_{z^{2}}$$ orbitals and unoccupied Mo 4$$d_{x y} / d_{x^{2}-y^{2}}$$^[41,42]. However, for the metallic feature 1T-MoS₂ monolayer, the Mo 4d orbitals split into the partially filled 4d_xy/d_yz/d_xz orbitals and non-occupied Mo 4$$d_{z^{2}}$$ and Mo 4$$d_{x^{2}-y^{2}}$$ levels. Thus, altering the coordination environment of the transition metal atoms may lead to the rearrangement of d electrons, giving rise to distinct electronic or magnetic properties.

Here, for the ic-2D systems, as revealed by the electronic density of states (DOS) and band structure calculations, all the ic-2D systems exhibit metallicity. Taking the case of Ti₇S₁₂ as an example, the metallic behavior is revealed by the Fermi level (E_F) crossing through the electronic bands, and this is mainly contributed by the Ti d orbitals [Figure 3A]. A similar phenomenon is found in the other six ic-2D systems [Supplementary Figure 5]. We also inspect the magnetic properties of ic-2Ds. Among the seven systems, non-zero magnetic moments can be found in five of them, i.e., V₇S₁₂, Cr₇S₁₂, Mn₇S₁₂, Fe₇S₁₂, and Co₇S₁₂, and they exhibit a ferromagnetic ground state (GS). Interestingly, the ferromagnetic property of Co₇S₁₂ is solely induced by the intercalated Co atom since the perfect monolayer 1T CoS₂ possesses a non-magnetic GS. In the pristine 1T CoS₂ system, the electronic states at E_F mainly consist of cobalt-d_z², cobalt-d_x²_-y², and cobalt-d_xy orbitals [Figure 3B]. For Co₇S₁₂, on the other hand, the intercalated cobalt atoms introduce additional bands (spin-split bands) across the E_F, and a magnetic GS, therefore, develops [Figure 3C and D]. Moreover, the introduced magnetic moments are mainly localized on the d orbitals of the intercalated cobalt atom, and the net magnetic moment is calculated as 1.17 μB/unit cell, which correlates with a large degree of charge transfer, nearly 0.6 |e| between the inserted Co and the surrounding S atoms as calculated by Bader charge analysis [Supplementary Figure 6]. On the other hand, V₇S₁₂, Cr₇S₁₂, Mn₇S₁₂, and Fe₇S₁₂ have ferromagnetic 1T-phase MS₂, and their ferromagnetism is contributed by both the intercalated ions and the pristine layers.

Figure 3. (A) Electronic band structures and DOS plots of Ti₇S₁₂. (B) The orbital-resolved band structure of pristine 1T-CoS₂. (C and D) The orbital-resolved and spin-resolved band structures of Co₇S₁₂ (left: spin-up; right: spin-down).

Performance as anode materials for LIBs

The predicted ic-2D crystals exhibit a morphology rather similar to TMDs, which endows them with similar characteristics such as high specific surface area and metallic features, encouraging us to investigate their capability as anodes for LIBs.

In order to study the properties of lithium adsorption on ic-2Ds, a 2 × 2 × 1 supercell is applied to investigate the possible lithium adsorption sites. Six high-symmetry adsorption sites on the surface, as well as the intercalated layer of ic-2D crystals, are considered [Supplementary Figure 7], and the adsorption property and favorable adsorption sites are investigated by comparing their adsorption energies E_ad quantitatively. The results show that a single Li atom prefers to be adsorbed on the S5 site for Ti₇S₁₂, V₇S₁₂, Mn₇S₁₂, and Ni₇S₁₂ and on the S1 site for Cr₇S₁₂, Fe₇S₁₂, and Co₇S₁₂ [Supplementary Table 4].

To explore the maximum theoretical capacity of ic-2D, layers of Li atoms are gradually added to the system. If we take Ti₇S₁₂ as an example, the first lithium layer prefers to be intercalated on the energetically most favorable S5 sites, followed by the redetermined S1 sites, and the third and fourth lithium layers are on the S6 and S4 sites, respectively [Figure 4A]. After layers of Li atoms are adsorbed on both sides of Ti₇S₁₂ crystals, the OCVs calculated with the number of Li layers ranging from one to seven are 3.8, 2.1, 1.7, 1.1, 0.9, 0.8, and 0.7 V, correspondingly [Supplementary Figure 8]. These positive OCV values suggest that seven lithium layers can be adsorbed on Ti₇S₁₂ ic-2D with high stability, and this corresponds to a high theoretical capacity value of 745.6 mA hg^-1 [Figure 4B]. As a comparison, the value for the commercialized graphite anode is only 372 mA hg^-1. In addition, due to the long distance between the third lithium layer and the ic-2D substrate, the chemical interactions between them are very weak, so the outer lithium layer is not further considered. A similar feature is observed in V₇S₁₂ [Supplementary Figure 9]. It should be noted that other systems cannot host up to seven layers of Li atoms owing to the collapse of the crystal structures. For instance, Co₇S₁₂ and Ni₇S₁₂ can only hold three and four layers, respectively.

Figure 4. (A) A schematic geometric structure of ic-2D with lithium atoms adsorbed at two layers on one side and two intercalated and (B) maximum theoretical capacity of ic-2Ds.

The metal ion diffusion is directly related to the cycle and rate capability of batteries. Herein, the climbing-image nudged elastic band method is used to study the Li diffusion process from one lithium adsorption site in the most energetically favorable configuration to the most neighboring stable site on both the surface and the intercalated layer [Figure 5]. For the surface diffusion process, we consider three different diffusion paths, which are denoted as path one, path two, and path three by red, black, and blue arrows, respectively [Figure 5A], which connect the two neighboring most preferred adsorption sites. For all the seven systems, the lowest diffusion barriers are calculated to be along path 1, and the values range from 0.17 to 0.69 eV; on the other hand, for the interlayer diffusion pattern, the calculated energy barriers are much higher, ranging from 0.58 to 1.58 eV, where Li diffusion on V₇S₁₂, Cr₇S₁₂, Mn₇S₁₂, Fe₇S₁₂, and Co₇S₁₂ is along path 5, diffusion on Ti₇S₁₂ is along path 6, and diffusion on Ni₇S₁₂ is along path 4 [Figure 5B]. This can be ascribed to the strong interlayer vdW interactions that impede the movement of Li ions along the interlayer. Therefore, Li ions prefer to diffuse along the surface instead of in the intercalation layer [Table 1].

Figure 5. Schematic illustration of the possible diffusion paths of Li on ic-2Ds. (A) Top views of the three possible surface diffusion paths and the corresponding barriers for ic-2Ds. (B) Top views of the three possible interlayer diffusion paths and the corresponding barriers for ic-2Ds.

To summarize, two ic-2Ds, i.e., Ti₇S₁₂ and V₇S₁₂, are theoretically predicted to be promising ic-2D anode materials for LIBs because of their relatively higher lithium capacity (745.6 mA hg^-1 for Ti₇S₁₂ and 723.9 mA hg^-1 for V₇S₁₂) and smaller Li diffusion barrier (0.17 eV surface barrier and 1.17 eV interlayer barrier for Ti₇S₁₂, 0.20 eV surface barrier and 0.92 eV interlayer barrier for V₇S₁₂). Besides, the metallicity of the two monolayers can be well maintained even after lithiation [Supplementary Figure 10], facilitating fast electron transfer.

Electrocatalytic performance of ic-2D crystals

The interlayer potential barriers of vdW gaps in 2D TMDs have been widely proven to cause poor interlayer charge transport^[43,44]. We believe that filling the vdW gap in TMDs with intercalants would be beneficial to address this bottleneck and thus enhance the HER performance of TMDs^[45]. To demonstrate the capability of ic-2Ds as electrocatalysts, we evaluated their catalytic performance toward HER at different θ levels.

Similar to the case of lithium adsorption, here, three adsorption sites (S1, S2, S3) are considered typical ic-2Ds [Figure 6A and B]. Since a 2 × 2 × 1 supercell of ic-2Ds contains 12 nonmetal atoms on the top layer of atoms, the hydrogen coverage (θ) is calculated to be 1/12 = 8.3%. At this low θ value, hydrogen prefers to be adsorbed on the S1 site for most ic-2Ds [Figure 6C], except for Ni₇S₁₂ [Supplementary Figure 11]. After hydrogen adsorption, the metallicity of ic-2Ds can also be maintained [Supplementary Figure 12], which is beneficial to electrocatalysis^[46,47].

Figure 6. Possible active sites of H adsorption on ic-2Ds from the (A) top and (B) side views. (C) Optimized structure for H on Co₇S₁₂ at 8.3% coverage. (D) Hydrogen adsorption Gibbs free energy on ic-2Ds.

Computational screening for the potential HER electrocatalyst is based on the descriptor: the hydrogen adsorption Gibbs free energy (ΔG_*H). According to the Sabatier principle^[48], an HER catalyst is promising when hydrogen is adsorbed on the catalyst surface neither too strongly nor too weakly, that is when ΔG_*H is close to zero^[49-51]. In this work, the energetically most favorable configurations of hydrogen adsorbed on ic-2Ds are adopted to calculate ΔG_*H [Supplementary Table 5]. Among them, Ti₇S₁₂, V₇S₁₂, Fe₇S₁₂, and Ni₇S₁₂ ic-2Ds bind H too weakly with large positive values of ΔG_*H [Figure 6D], while Mn₇S₁₂ binds H too strongly with the Gibbs free energy of -0.61 eV (which lead to difficult H desorption in the Heyrovsky/Tafel step). In comparison, Cr₇S₁₂ (0.09 eV) and Co₇S₁₂ (0.05 eV) bind H moderately with close-to-zero ΔG_*H values, indicating that Cr₇S₁₂ and Co₇S₁₂ are excellent HER electrocatalysts under this θ and the adsorption free energy is mainly determined by the nonmetal element sulfur that binds hydrogen.

Since the θ also affects ΔG_*H, here, we select Cr₇S₁₂ and Co₇S₁₂ to further study their catalytic performance at higher θ. According to our calculations, with θ increased from low (8.3%) to medium (25%) value, the ΔG_*H values of Cr₇S₁₂ and Co₇S₁₂ increase up to 0.45 eV and 0.13 eV, respectively. When θ is further increased to full coverage (100.00%), both Cr₇S₁₂ and Co₇S₁₂ exhibit no bond with H because of the large positive value of ΔG_*H [Supplementary Figure 13]. In general, Cr₇S₁₂ and Co₇S₁₂ can boost the HER at low θ but are not promising at high θ.