Ultrahigh mobility and Rashba spin splitting in Sb-substituted bismuth telluride and bismuth selenide

Abstract

Topological insulators (TIs) such as Sb-doped Bi₂Te₃ and Bi₂Se₃ exhibit promising phenomena for advanced spintronics. While previous studies explored isolated doping levels; a systematic understanding of how Sb concentration influences topological behavior, Rashba-type spin splitting, and surface state formation is lacking. Here, we use density functional theory to investigate the structural, electronic, topological and transport properties of (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃ thin films across 0 ≤ x ≤ 1. We identify pronounced Rashba spin splitting in Bi₂Te₃ at x = 0.5, 0.6, and 0.9 with in-plane helical spin textures. We identified the orbital origins of topological surface states and demonstrate that band inversion persists across the Sb doping range. At x = 0.2, 0.4, and 0.8, calculated surface electron mobilities are consistent with experiments and increase an order of magnitude over Sb₂Te₃, with minimal impact on bulk mobilities. These insights advance our understanding of TIs for spintronic and quantum device applications.

---

Introduction

Topological insulators (TIs) are quantum materials characterized by an insulating bulk and topologically protected metallic surface states exhibiting spin-momentum locking.^1–3 The nontrivial topology arises from spin–orbit coupling (SOC)-driven band inversion, leading to Dirac-like surface states with significant potential for spintronic and quantum applications.^2,4,5

Beyond the fundamental topological surface states, another SOC-driven phenomenon plays a crucial role in TI functionality:^6,7 Rashba-type spin splitting, which emerges under broken inversion symmetry.⁸ Unlike bulk band inversion, Rashba splitting is typically induced by surface reconstruction⁹ or interface effects,^10–12 resulting in momentum-dependent spin textures that can be systematically controlled through chemical doping.^13–15

While the effects of individual dopants on topological insulators have been widely studied,^16–19 the interplay between doping-induced modifications and Rashba-type effects on surface state properties remains largely unexplored. Previous studies have primarily focused on isolated doping scenarios, neglecting a critical aspect. When dopant atoms exhibit correlations or mutual interactions, the physical behavior fundamentally changes, and Rashba effects become significantly more pronounced.

This correlated doping regime offers control over TI properties, potentially enabling fine-tuning of both topological protection and spin–orbit phenomena in ways not achievable through conventional random doping approaches.^14,15,20 This research gap is particularly crucial, as optimal TI functionality likely depends on a delicate balance between topological robustness and tunable spin–orbit effect, where dopant correlations and enhanced Rashba contributions could serve as a critical pathway for achieving advanced quantum control. Such a balance can only be realized through systematic compositional engineering that incorporates dopant interactions and their amplified impact on surface state manipulation.

Among known TIs, Bi₂Te₃ ²¹ and Bi₂Se₃ have emerged as benchmark systems due to their well-defined topological properties,^22–24 layered crystal structures, and strong SOC. These materials have been widely studied both theoretically and experimentally, serving as model platforms for exploring the manifestation of topological surface states (TPSS) and their manipulation via external fields and chemical modifications.^4,25,26 A common strategy to modulate their electronic structure is isovalent substitution, particularly by introducing antimony (Sb) in place of bismuth (Bi). Sb shares the same valence electron configuration with Bi, making it a chemically compatible dopant that allows for electronic tuning without causing severe distortions to the crystal lattice.^27,28

Sb substitution in Bi₂Te₃ and Bi₂Se₃ has been extensively explored through experimental and density functional theory (DFT) studies. For 20% Sb-doped Bi₂Se₃ thin films, Sb shifts the Dirac cone below the Fermi level and enhances surface–bulk hybridization.²⁹ In bulk (Bi_1−xSb_x)₂Te₃ (0 ≤ x ≤ 1), DFT revealed that Sb tuning modifies the Dirac cone, spin–orbit coupling strength, and defect formation energies.³⁰ Experimental and DFT studies on ultrathin Bi_1.5Sb_0.5Te_1.7Se_1.3 showed a thickness-dependent hybridization gap and topological phase transitions.³¹ A theoretical DFT-based study on 3D Bi₂Se₃-type topological insulators (Bi₂Te₃, Bi₂Se₃, Sb₂Te₃) developed an effective four-band k·p Hamiltonian to describe bulk electronic structures and surface Dirac states, providing a foundation for modeling topological features in these materials.³² Another DFT investigation on Bi₂Se₃ nanosheets with Sb₂Se₃ and SbBiSe₃ surface capping layers analyzed the role of interface asymmetry and biaxial strain (±5%), showing that Sb-based capping opens a band gap by breaking inversion symmetry, while strain can restore or tune the topological surface states and spin-momentum locking.³³

While Sb doping^27,28,34 has been extensively used to shift the Fermi level, the full extent of its impact on the electronic topology, including the evolution of band inversion, the emergence or suppression of Rashba-type splitting, and the modification of SOC strength, is not yet fully clarified. In particular, previous studies have often focused on isolated doping levels or specific material systems, lacking a comprehensive and systematic investigation across a continuous doping range. Additionally, the surface versus bulk character of emerging states under doping has not been thoroughly analyzed, leaving questions about the robustness of topological protection under high substitution levels.

Recent work on (Bi_1−xSb_x)₂Te₃ demonstrated that Sb concentration can be used to suppress bulk conduction and reveal surface-dominated transport at x ≈ 0.5, enabling effective gating and ambipolar field effects.³⁴ At the same time, giant and tunable Rashba-type spin splitting has been observed in polar Bi-based systems such as BiTeI and surface-doped Bi₂Se₃, where strong SOC and inversion asymmetry induce spin-resolved electronic states suitable for spintronic applications.^35,36 Despite these advances, the systematic emergence and evolution of Rashba-type features across the full Sb doping range remain unresolved. Moreover, the competition and possible coexistence between Rashba-split states and topological surface states in doped Bi₂Te₃ and Bi₂Se₃ are not yet fully understood.

In this study, we address these critical gaps by performing a systematic DFT analysis on two-dimensional (2D) Sb-doped Bi₂Te₃ and Bi₂Se₃ across a wide doping range, from pristine compound to the fully Sb-substituted case. By examining the evolution of the electronic band structure, orbital-projected contributions, and symmetry-breaking effects, we aim to clarify how Sb doping modifies key topological and spin-related features. Special attention is devoted to identifying the existence of band inversion features, assessing the onset and strength of Rashba-type spin splitting, and examining how these topological features respond to varying doping levels. This work provides a detailed theoretical framework to understand the doping-dependent transformation of topological insulators, offering insights that are essential for their integration into next-generation spintronic and quantum technologies. We calculate the electron effective masses, bulk and surface electronic mobilities as a function of Sb concentration. Sb-alloyed Bi₂Te₃ derivatives exhibit strong spin–orbit coupling, giving rise to Rashba-type spin splitting and high spin-charge conversion. To describe this situation, semiclassical Boltzmann transport calculations combined with DFT-derived parameters provide an efficient route to link SOC. In this work, we focus on the (Bi_0.4Sb_0.6)₂Te₃ composition, analyzing its Rashba-split bands and transport properties using the Rashba model.³⁷ In particular, existing Rashba-type splitting and spin-momentum locking in these systems highlights their potential utility in spin field-effect transistors, spin–orbit torque-based memory, and topological quantum computing architectures.

Computational methodology

All calculations were carried out using density functional theory (DFT) as implemented in the Vienna ab initio simulation package (VASP).^38,39 The exchange–correlation interactions were treated within the generalized gradient approximation (GGA) using the Perdew–Burke–Ernzerhof (PBE) functional.⁴⁰ The projector augmented-wave (PAW) method was employed to construct pseudopotentials, ensuring an accurate description of the electronic interactions of individual atoms.⁴¹ To account for relativistic effects in the electronic structure calculations, SOC was incorporated alongside GGA. A plane-wave basis set with a kinetic energy cutoff of 500 eV was used. The Gaussian smearing method was applied to total-energy calculations with a smearing width of 0.01 eV. A vacuum spacing of approximately 20 Å was introduced, isolating the slab consisting of five quintuple layers (QL) of Bi₂Te₃ and Bi₂Se₃. The equilibrium interlayer distances between QLs were found to be approximately 2.62 Å for Bi₂Te₃ and 2.70 Å for Bi₂Se₃. To systematically investigate the effects of Sb doping, Bi atoms were substituted with Sb at varying concentrations while preserving the TI character of the structures. All configurations were generated using consistent substitution patterns to maintain structural uniformity across the doping series. The Sb concentration was varied in increments of 0.1, resulting in eleven distinct compositions ranging from x = 0 to x = 1.0 for both (Bi_1−xSb_x)₂Te₃ (BST) and (Bi_1−xSb_x)₂Se₃ (BSS) systems. During ionic relaxations, the Brillouin zone (BZ)⁴² was sampled using a 12 × 12 × 2 Gamma-centered k-point mesh for structural optimization and a 9 × 9 × 1 Monkhorst–Pack k-point mesh, which was refined to 18 × 18 × 1 for more precise electronic density of states calculations. The self-consistent DFT iterations were considered converged when the energy difference between consecutive steps fell below 10⁻⁵ eV. Additionally, the residual stress in the unit cell was minimized to less than 1 kbar along all three directions. Additionally, ab initio molecular dynamics simulations (AIMD) were performed with a 3 × 3 × 1 supercell (225 atoms) at 300 K of 5 ps to examine thermal stability of (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃ as a function of Sb concentration. Because of the huge amount of time consumed in running an AIMD calculation, a k-point mesh of 2 × 2 × 1 was used for the BZ integration. The transport properties of two-dimensional electron gas (2DEG) were analyzed using the Boltzmann Transport Equation (BTE) under the Relaxation-Time Approximation (RTA), including the Gaussian impurity scattering model.

To cross-check Rashba-related features and corresponding transport response we can use the Hamiltonian and DFT-derived parameters as an input to the Boltzmann Transport Equation (BTE), validating that the observed spin splitting behavior originate from structural inversion asymmetry. The velocity and spin expectation values were computed over the Fermi surface to calculate the corresponding transport tensors, including Edelstein and conductivity tensors.⁴³ A numerical solver was implemented using single-particle Hamiltonian which includes Rashba spin–orbit coupling and warping effect, expressed as:

where m* is the effective mass, λ_w is the hexagonal warping term (Even though the material is anisotropic, the warping term was set to zero in our model) that arises from the C_3v crystal symmetry, and α_R is the Rashba parameter which originates from the structural inversion asymmetry. The eigenvalues and eigenvectors of the Hamiltonian were calculated numerically to obtain energy and velocity components. The Fermi surface was detected using an adaptive refinement algorithm that increases the sampling density around the target energy level.⁴⁴ This method enables more precise energy-dependent line or surface integrals in reciprocal space. To cross-check the Rashba-type splitting observed in DFT band structures of (Bi_1−xSb_x)₂Te₃ for doping levels with x = 0.6, the parameters α_R and m* were calculated using DFT data near the Γ-point (k = 0). These parameters were used as inputs to the model to calculate band structures, electron mobility, and transport properties.

Results and discussion

Threshold thickness for TI properties in pristine Bi₂Te₃ and Bi₂Se₃

We perform the electronic band structures of Bi₂Te₃ with orbital projections onto the p orbitals of constituent atoms in 1. In these plots, blue dots represent the contribution from Bi (p) orbitals, while red dots indicate Te (p) orbitals. In Fig. 1, the band structure of Bi₂Te₃ reveals an energy gap of 408 meV at the Γ-point for 1QL, with no visible band inversion. At 2QL, the band gap drops to 114 meV, reflecting the strong sensitivity of the band structure to quantum confinement. By 3QL, a topological phase transition becomes evident through the closure of the hybridization gap and the onset of band inversion, establishing 3QL as the critical thickness for the emergence of topological surface states in Bi₂Te₃. The electronic band structure of Bi₂Se₃ exhibits a clear thickness-dependent evolution from 1 to 6QLs, as shown in Fig. 1 in the SI. At 1QL, a large trivial band gap of 818 meV is present at the Γ-point, with no indication of band inversion. As the thickness increases, the gap gradually narrows, decreasing to 139 meV at 2QL, 44 meV at 3QL, and 26 meV at 4QL. This critical thickness (4QL) marks the onset of topological insulating behavior. Beyond this point, at 5QL and 6QL, the band gaps reduce further to 12 meV and 5 meV, respectively, and clear signatures of band inversion emerge, signifying the transition to a topologically non-trivial phase. These findings emphasize the crucial role of thickness in governing the topological phase transitions in Bi₂Te₃ and Bi₂Se₃. As the number of QLs increases, the hybridization-induced gap steadily decreases, eventually allowing the intrinsic band inversion driven by SOC to prevail. In this study, a thickness of 5QL was chosen to ensure the emergence of topologically non-trivial states while also providing a practical framework for exploring a wide range of Sb doping concentrations.

Fig. 1 Thickness-dependent electronic band structures of Bi₂Te₃ from 1 to 6QLs, and projected onto Bi (blue) and Te (red) p-orbitals. The color intensity indicates the orbital contribution to each band. Qualitatively similar behavior has been observed for the thickness-dependent Bi₂Se₃ films, whose electronic band structures and band inversion have been provided in the Fig. S1 in the SI.

The emergence of TI properties in Bi₂Te₃ and Bi₂Se₃ is fundamentally governed by two key mechanisms: band inversion and band gap closing, both of which occur as a function of thickness. In their bulk form, these materials exhibit strong SOC, which induces an inverted band structure and supports the formation of gapless Dirac surface states. However, when these compounds are confined to ultrathin films, the interaction between the top and bottom surfaces becomes significant. This inter-surface hybridization opens a gap at the Dirac point, effectively suppressing the TI phase and rendering the system a trivial insulator.⁴⁵ This hybridization gap dominates the electronic structure in films with very small thicknesses, where the surface wavefunctions strongly overlap. With increasing film thickness, the interaction between opposite surfaces weakens, allowing the emergence of band inversion. Once the hybridization gap weakens sufficiently, the system undergoes a topological phase transition, restoring the gapless Dirac surface states that are the hallmark of TIs. The critical thickness at which this transition occurs differs between Bi₂Te₃ and Bi₂Se₃, owing to intrinsic differences in their electronic band structures and the strength of their SOC.^46–48

The first-principles calculations⁴⁶ confirm that the decay length of surface states is approximately 1 nm, resulting in strong hybridization in 1QL films. This hybridization effect becomes progressively weaker beyond 2QL, leading to reduced gap opening. Experimentally, this manifests as the emergence of a single Dirac cone surface state at 2QL.⁴⁶

The structural properties of (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃

The V–VI compound family (including Bi₂Te₃, Sb₂Te₃, Bi₂Se₃, and Sb₂Se₃) share a common rhombohedral primitive structure (space group R3̄m) and a hexagonal conventional cell, with the notable exception of Sb₂Se₃, which naturally crystallizes in an orthorhombic phase (Pnma).²¹ While the first three compounds are well-established TIs, Sb₂Se₃ is a trivial insulator in its equilibrium phase.⁴⁵ However, theoretical and experimental studies have demonstrated that Sb₂Se₃ can undergo a topological phase transition under specific conditions such as via external strain⁴⁹ or by forming a heterointerface with a topological phase like Bi₂Se₃ ⁵⁰ where interfacial hybridization induces nontrivial surface states.

Fig. 2 illustrates the side-view atomic structures of the 5QL slabs used in this study for both (Bi_1−xSb_x)₂Te₃ (BST) and (Bi_1−xSb_x)₂Se₃ (BSS) systems. As shown in Fig. 2(a), each QL consists of five atomic planes arranged in a Te(Se)–Bi–Te(Se)–Bi–Te(Se) sequence, with an inter-QL spacing.

Fig. 2 Side-view atomic structures of 5QL Sb-doped (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃ systems, modeled with identical atomic configurations at each composition. Panels (a)–(k) show the structural evolution with increasing Sb concentration from x = 0 to x = 1.0 in increments of 0.1. Panel (l) displays the corresponding two-dimensional Brillouin zone and the high-symmetry path used in the electronic structure calculations.

To systematically investigate the effect of Sb concentration on structural and electronic properties, we substituted Bi atoms with Sb in 0.1 increments, generating eleven configurations from x = 0 to x = 1.0. The goal was to identify the doping level at which TI behavior emerges, as shown in Fig. 2(a)–(k). For each doping level, the configuration with the lowest formation energy was selected to ensure structural stability and physical relevance in subsequent analyses (SI Table S1). All structures were modeled within the same rhombohedral crystal framework, using consistent substitution patterns to ensure structural uniformity across the series. This controlled setup allows for a clear comparison of doping effects without interference from phase transitions or distortions. Electronic band structure and spin texture were calculated along the high-symmetry path K̅–[capital Gamma, Greek, macron]–M̅ in the two-dimensional surface Brillouin zone, as illustrated in Fig. 2(l).

To evaluate structural stability, we computed both the formation energy (E_for) and cohesive energy (E_coh) per atom using the following expressions:

where E_BST(BSS) denotes the total energy of the BST or BSS structure, E_{(Bi,Sb,Te,Se)} represent the energies of isolated single atoms Bi, Sb, Te, and Se, and n_{(Bi,Sb,Te,Se)} are the numbers of each respective atom type in the structure. The chemical potentials, μ_{(Bi,Sb,Te,Se)}, are determined from their 3D bulk phases, where for each element, the atomic chemical potential μ_atom is estimated by the equation μ_atom = E_bulk/N, with E_bulk being the total energy of the bulk structure and N the number of atoms in the bulk unit cell.

Fig. 3 shows the systematic evolution of the in-plane lattice parameters and the out-of-plane lattice constant, represented by the 5QL slab thickness, as a function of Sb concentration for both BST and BSS systems. In BST, the in-plane and out-of-plane lattice constants a and c decrease from 4.302 and 49.658 Å at x = 0.1 to 4.247 and 48.828 Å at x = 0.9. A similar pattern emerges in BSS, with a and c shrinking from 4.107 and 45.792 Å to 4.025 and 45.476 Å. This behavior arises from the smaller atomic radius of Sb compared to Bi, leading to a lattice shrinkage upon substitution. The nearly linear variation of lattice parameters with composition in both systems is consistent with Vegard's law. Additionally, our ab initio molecular dynamics (AIMD) simulation results also demonstrate that these materials remain thermally stable at room temperature as shown in Fig S2 within the SI.

Fig. 3 Vegard's law: variation of in-plane (a) and out-of-plane (c) lattice constants of 5QL (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃ as a function of Sb concentration x. Solid markers represent DFT-calculated values, while dashed lines indicate linear trends based on Vegard's law.

Tables 1 and 2 show that interlayer chalcogen–chalcogen (Te–Te or Se–Se) distances do not follow a monotonic trend. In BST, Te–Te separations fluctuate between 2.854 and 2.925 Å, while in BSS, Se–Se distances vary from 2.669 to 2.760 Å. These deviations likely result from local structural relaxations or distortions induced by Sb substitution. Importantly, both E_for and E_coh decrease progressively with increasing Sb concentration, suggesting a decline in thermodynamic stability. Although the extent of destabilization is comparable between systems, the lower absolute cohesive energy in BST is consistent with the softer bonding nature of telluride systems. These results indicate that while Sb substitution is energetically feasible, higher doping levels reduce lattice cohesion and compromise bonding integrity.

Table 1 Calculated formation and cohesive energies per atom (in eV), average interlayer distance (d_int-ave in Å), Γ-point bandgap (E_Γ-gap in eV), and Rashba parameter (α_R in eV Å) for 5QL (Sb_xBi_1−x)₂Te₃ at varying Sb content x

Compound	Efor (eV)	Ecoh (eV)	dint-ave (Å)	EΓ-gap (eV)	αR (eV Å)
Bi2Te3	−0.343	3.121	2.870	0.009	—
(Bi0.9Sb0.1)2Te3	−0.329	3.110	2.880	0.030	—
(Bi0.8Sb0.2)2Te3	−0.314	3.101	2.877	0.069	—
(Bi0.7Sb0.3)2Te3	−0.299	3.091	2.889	0.031	—
(Bi0.6Sb0.4)2Te3	−0.285	3.081	2.854	0.009	—
(Bi0.5Sb0.5)2Te3	−0.270	3.071	2.892	0.082	0.328
(Bi0.4Sb0.6)2Te3	−0.257	3.061	2.920	0.131	2.098
(Bi0.3Sb0.7)2Te3	−0.243	3.052	2.912	0.031	—
(Bi0.2Sb0.8)2Te3	−0.229	3.043	2.925	0.024	—
(Bi0.1Sb0.9)2Te3	−0.215	3.034	2.905	0.047	0.719
Sb2Te3	−0.202	3.025	2.925	0.004	—

---

Table 2 Calculated formation and cohesive energies per atom (in eV), average interlayer distance (d_int-ave in Å), Γ-point bandgap (E_Γ-gap in eV), and Rashba parameter (α_R in eV Å) for 5QL (Sb_xBi_1−x)₂Se₃ at varying Sb content x

Compound	Efor (eV)	Ecoh (eV)	dint-ave (Å)	EΓ-gap (eV)	αR (eV Å)
Bi2Se3	−0.513	3.372	2.682	0.010	—
(Bi0.9Sb0.1)2Se3	−0.495	3.358	2.694	0.020	—
(Bi0.8Sb0.2)2Se3	−0.476	3.343	2.718	0.064	0.304
(Bi0.7Sb0.3)2Se3	−0.458	3.330	2.711	0.034	—
(Bi0.6Sb0.4)2Se3	−0.440	3.317	2.702	0.023	—
(Bi0.5Sb0.5)2Se3	−0.422	3.303	2.736	0.035	—
(Bi0.4Sb0.6)2Se3	−0.405	3.290	2.738	0.079	—
(Bi0.3Sb0.7)2Se3	−0.388	3.278	2.757	0.067	—
(Bi0.2Sb0.8)2Se3	−0.371	3.265	2.731	0.080	—
(Bi0.1Sb0.9)2Se3	−0.354	3.253	2.760	0.079	—
Sb2Se3	−0.337	3.241	2.771	0.078	—

---

Electronic band structure and Rashba-type spin splitting

To elucidate the impact of Sb doping on the electronic properties of (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃, we computed their electronic band structures with SOC for Sb concentrations ranging from x = 0 to x = 1 in increments of 0.1. The results are presented in Fig. 4, with panels (a) and (b) corresponding to BST and BSS, respectively. All calculations were performed on 5QL slabs to ensure the emergence of topologically protected surface states and to facilitate the systematic arrangement of doping concentrations.

Fig. 4 Electronic band structures of 5QL (a) (Bi_1−xSb_x)₂Te₃ and (b) (Bi_1−xSb_x)₂Se₃ for varying Sb concentrations (x = 0 to 1), showing the evolution of surface states and band dispersion with doping.

In the BST system [Fig. 4(a)], we observe a progressive evolution of surface states with increasing Sb concentration. Notably, a Dirac gap opens at the Γ-point, which coexists with Rashba-type splitting of surface bands around this region. The Rashba splitting exhibits a non-monotonic trend with increasing Sb concentration. The most pronounced splitting occurs at x = 0.6, which also corresponds to the maximum Dirac gap opening at the Γ-point (see Table 1). The finite bandgap opening can thus be attributed to residual coupling between top and bottom surface states in thin films, where Sb-induced local structural asymmetry modifies the decay length and overlap of surface Bloch states, leading to inter-surface hybridization and the formation of a Dirac gap at Γ even without breaking time-reversal symmetry.³³ Additional peaks in Rashba-type separation are identified at x = 0.5 and x = 0.9, suggesting that structural inversion asymmetry (SIA) is modulated in a non-linear manner by Sb doping. This behavior likely originates from complex lattice relaxations and local potential gradients induced by Sb substitution.

In contrast, the BSS system [Fig. 4(b)] exhibits a notably different trend. Despite following identical doping configurations, significant band splitting is observed only at x = 0.2. For all other doping levels, the Dirac-like surface states remain essentially degenerate near the Γ-point. This disparity between BST and BSS highlights the material-specific sensitivity of Rashba splitting to doping, which can be attributed to differences in lattice stiffness, intrinsic SOC strength, and atomic hybridization. Specifically, the heavier Te atom in BST contributes to stronger SOC, amplifying spin-splitting effects.

To provide direct evidence of Rashba-type interaction, we examine the spin-resolved band structure and spin textures for the x = 0.6 case in BST, where the band splitting is maximal (Fig. 5). Fig. 5(a) shows the spin-split band structure near the Γ-point, with clear formation of inner and outer bands, resembling canonical Rashba behavior. Fig. 5(b) depict the in-plane isoenergetic spin textures for the outer and inner bands at a constant energy. The outer band exhibits counterclockwise (CCW) spin rotation, while the inner band shows clockwise (CW) rotation, confirming helical spin-momentum locking typical of Rashba systems.^35,51

Fig. 5 Spin-resolved electronic structure of 5QL (Bi_0.4Sb_0.6)₂Te₃. (a) Spin-split band structure near the Γ-point showing inner and outer Rashba-type bands. (b) Isoenergy spin textures at E ≈0.3 eV for both bands. (c) Schematic representation of the Rashba split bands. (d) Rashba band dispersion obtained from Rashba model, confirming the consistency between the analytical Rashba model and DFT-derived spin-split bands.

To further support this, Fig. 5(e) illustrates the spin textures at E = 0.3 eV, where both the inner and outer bands remain visible, reinforcing the non-degenerate, spin-polarized nature of the states. These results are well described by the Rashba Hamiltonian, which captures the spin splitting induced by structural inversion asymmetry and spin–orbit interaction:^32,52

where α_R is the Rashba parameter, quantifying the strength of the linear-in-momentum spin–orbit coupling. Here, σ_x and σ_y are the Pauli matrices, and k = (k_x, k_y) is the in-plane wave vector.

Diagonalizing this Hamiltonian yields the energy dispersion relation:

which describes two spin-split parabolic bands. The splitting in momentum space is given by:where m* is the effective mass, and the energy separation between the two branches is:

ΔE = 2α_Rk (7)

As illustrated in Fig. 4(a) for x = 0.6, the Rashba effect results in a characteristic band structure featuring a momentum offset (k₀ = Δk/2) and a Rashba energy E_R, defined as the energy difference between the band minimum and the degeneracy point at k = 0. The Rashba parameter α_R can be conveniently expressed in terms of E_R and k₀ as:

This pronounced spin splitting, particularly evident in BST at intermediate Sb concentrations, arises from the strong spin–orbit coupling and structural inversion asymmetry, manifesting as Rashba-type band splitting with helical spin textures, consistent with the spin-split surface states observed in systems such as Au(111) due to broken inversion symmetry and SOC.⁵³

A quantitative analysis of the extracted α_R values (provided in Tables 1 and 2) reveals a peak at x = 0.6 in BST, where α_R = 2.098 eV Å, indicating a remarkably strong linear Rashba interaction. In contrast, the BSS system shows a significantly weaker response, with a maximum value of only α_R = 0.304 eV Å at x = 0.2. This disparity highlights the stronger structural inversion asymmetry (SIA) and spin–orbit coupling in Te-based systems compared to their Se-based counterparts under similar doping configurations. The pronounced Rashba spin splitting observed in BST, especially at intermediate Sb doping, parallels the large and tunable Rashba effect reported in Bi₂Se₃-based 2DEGs,³⁶ highlighting how intrinsic SIA and strong SOC can induce spin-polarized, non-degenerate states even without external gating, offering a robust platform for spintronic applications. The calculated dipole moment (μ_z) and the LDOS distribution (SI Table S2 and Fig. S3) both show polarity reversal and charge asymmetry at x = 0.6, confirming strong SIA, which in turn leads to a pronounced Rashba splitting.

The observed Rashba parameter in BST not only exceeds that of BSS, but also surpasses many well-known Rashba-active materials. For instance, quantum-confined near-surface states in K-deposited Bi₂Se₃ yield α_R = 0.79 ± 0.03 eV Å (RB1) and α_R = 0.35 eV Å (RB2),⁵⁴ while magnetotransport studies on epitaxial Bi₂Te₃ films report α_R = 1.58 eV Å.¹³ Additionally, metallic systems such as Bi(111),⁵⁵ Ir(111),⁵⁶ and PtPb₄ ⁵⁷ exhibit lower values of approximately 0.55, 0.35, and 1.0 eV Å, respectively. Even the Bi/Ag(111) interface, engineered to maximize Rashba splitting, reaches approximately 3.0 eV Å under ideal conditions.⁵⁸ The fact that BST approaches these engineered values through doping alone demonstrates the effectiveness of Sb substitution in enhancing SIA. The pronounced Rashba splitting in BST, particularly at intermediate Sb concentrations, reflects robust spin–momentum locking and strong SOC, hallmark features of topological materials. These characteristics position Sb-doped Bi₂Te₃ as a promising platform for spintronic technologies, including spin field-effect transistors (SpinFETs), spin–orbit torque (SOT)-based memory devices, and future topological quantum computing architectures.

Deformation potential theory and carrier mobility in Sb-doped Bi₂Te₃

As originally proposed by Bardeen and Shockley,^59,60 the coupling between charge carriers and acoustic phonons can be effectively described using deformation potential theory, which has since become a widely adopted method for evaluating carrier mobility in two-dimensional (2D) materials. In this work, we employed the effective mass approximation to estimate the acoustic-phonon-limited carrier mobility, using the following expression:where e is the elementary charge, ħ is the reduced Planck constant, k_B is the Boltzmann constant, T is the absolute temperature (298 K), C is the in-plane elastic modulus, m* is the carrier effective mass, m_d is the density-of-states effective mass and E_d is the deformation potential constant.

Due to the intrinsic anisotropy of topological insulator (TI) materials, calculations were performed independently along the in-plane x- and y-directions. To enable directional strain application, the original rhombohedral crystal structures of (Bi_1−xSb_x)₂Te₃ with x = 0, 0.2, 0.6, 0.8, and 1 were transformed into orthorhombic supercells (see Fig. S4(a) in the SI). After structural relaxation, a series of uniaxial strains ranging from −6% to +6% (in 2% increments) were applied along each direction. The total energy corresponding to each strain state was computed using first-principles calculations, and the resulting energy–strain relationship was fitted using a quadratic function.

The elastic modulus C was obtained from the second derivative of the total energy with respect to strain as shown in Fig. S4(b–f) in the SI:

where S₀ is the equilibrium surface area of the relaxed system.

The deformation potential constant E₁ was extracted by tracking the shift in the energy of the conduction band minimum (CBM) of both surface and bulk states as a function of strain, shown in Fig. S5 and S6 in the SI. The vacuum level was used as a reference to determine the absolute positions of the band edges under strain. The energy shifts were linearly fitted, and the slope of the fit yielded the deformation potential:

The carrier effective mass m* was calculated by fitting the curvature of the electronic bands near the CBM or VBM using the parabolic approximation:

High-resolution band structure data along the x direction were analyzed using VASPkit. When necessary, the density-of-states effective mass was calculated using the geometric mean: .

In this study, the electron mobilities and related transport parameters of Sb-doped Bi₂Te₃ were analyzed across the doping range x = 0.0 to x = 1.0 at room temperature (298 K), with a particular focus on surface and bulk contributions. Hole mobilities were not evaluated separately for surface and bulk states due to the difficulty in distinguishing the corresponding bands. As shown in Table 3, surface mobility μ_Surface dominates over bulk mobility μ_Bulk across all compositions, confirming the expected behavior in topological insulator systems where robust surface states play a key role in conduction. A notable enhancement in surface mobility is observed around x = 0.4 and x = 0.8 reaching 1.04 × 10⁵ and 9.76 × 10⁴ cm² V⁻¹ s⁻¹, respectively, whereas the bulk mobilities remain significantly lower (<1.1 × 10⁴ cm² V⁻¹ s⁻¹). This behavior correlates with the minimized effective mass and moderate deformation potential values in this doping range, indicating reduced carrier-phonon scattering and enhanced carrier transport along the surface channels.

Table 3 Surface and bulk electron deformation potential constants (E₁), in-plane elastic constants (C_2D), effective masses (m*), density-of-states masses (m_d), and electron mobilities (μ) along the x and y directions at 298 K for various compositions of (Bi_1−xSb_x)₂Te₃

Compound	E^Surf1 (eV)	E^Bulk1 (eV)	C2D (N m^−1)	mSurf^a (m0)	mBulk^a (m0)	md,Surf (m0)	md,Bulk (m0)	μSurf (×10^4 cm^2 V^−1 s^−1)	μBulk	Method
a The film was grown on a buffer layer of In2Se3/(Bi0.5In0.5)2Se3 on an Al2O3 substrate.
Bi2Te3	8.18	7.35	253.03	0.05	0.31	0.04	0.24	3.38	0.13	In this study
(Bi0.8Sb0.2)2Te3	4.25	3.61	260.81	0.06	0.21	0.05	0.20	9.03	1.03	In this study
(Bi0.6Sb0.4)2Te3	5.94	5.17	266.48	0.04	0.30	0.04	0.33	10.43	0.21	In this study
(Bi0.2Sb0.8)2Te3	5.98	5.39	280.29	0.04	0.42	0.04	0.45	9.76	0.11	In this study
Sb2Te3	6.74	6.09	292.46	0.13	0.27	0.07	0.22	1.43	0.27	In this study
Bi2Te3	—	—	—	—	—	—	—	0.60	0.09	ME^61
Bi2Te3	—	—	—	—	—	—	—	—	0.05	MBE^62
Bi2Se3 (uncapped)^a	—	—	—	—	—	—	—	1.60	—	MBE^63
Bi2Se3 (Se-capped)^a	—	—	—	—	—	—	—	0.36	—	MBE^63
Bi2Se3 (MoO3-capped)^a	—	—	—	—	—	—	—	0.33	—	MBE^63
Graphene	—	—	—	—	—	—	—	—	1.00	CVD^64
Graphene	—	—	—	—	—	—	—	—	35.00	CVD^65

---

Our theoretical predictions, based on DFT calculations at 0 K for ideal 2D (Bi_1−xSb_x)₂Te₃ systems, demonstrate surface mobilities that are significantly higher than typical experimental values. For instance, we find that the calculated surface mobility peaks at x = 0.4, exceeding 10⁵ cm² V⁻¹ s⁻¹ which are consistent with the trends reported in experimental studies. For example, exfoliated Bi₂Te₃ ultrathin films have shown mobilities of approximately 6030 cm² V⁻¹ s⁻¹ at 2 K,⁶¹ considerably exceeding those obtained from MBE-grown (521 cm² V⁻¹ s⁻¹)⁶² and chemically synthesized samples. Although our theoretical results inherently overestimate absolute mobility values⁶⁶ due to the absence of phonon scattering, defects, and substrate effects; they capture the non-monotonic doping dependence, with maximum mobility occurring near x = 0.4–0.6. This trend is in good agreement with previous experimental observations of enhanced Rashba spin splitting and surface-dominated conduction in this doping regime.^13,67–69 The observed dependence on Sb content further reflects the underlying changes in electronic band structure, effective mass, and electron–phonon coupling, all of which are sensitively captured in first-principles deformation potential calculations. Thus, while quantitative differences remain between theory and experiment, the agreement in doping-dependent mobility behavior strongly supports the predictive power of DFT-based models for guiding the optimization of topological materials. When compared to high-mobility 2D materials such as graphene, where mobilities range from 10 000 cm² V⁻¹ s⁻¹ (ref. 64) in ambient conditions to 350 000 cm² V⁻¹ s⁻¹ for CVD-grown ultraclean samples,⁶⁵ Sb-doped Bi₂Te₃ films show competitive surface transport properties while maintaining topologically protected states. These results demonstrate the tunability of surface transport in (Bi_1−xSb_x)₂Te₃ through controlled doping and provide a pathway for engineering high-performance topological devices, particularly in regimes where surface conduction can be maximized and bulk contributions can be minimized.

For doping levels x = 0.6, the DFT-derived Rashba parameters correspond to a significant spin splitting. The BTE results reveals that the largest Edelstein response, which describes the efficiency and the direction of charge-to-spin conversion, and charge conductivity (shown in SI) occur near at x = 0.6, consistent with the DFT-predicted maximum Rashba-splitting at similar doping levels. Edelstein response quantifies how efficiently an applied charge current generates a net spin density, providing a direct measure of charge-to-spin conversion in Rashba systems, which is essential for spintronic applications. This confirms that moderate Sb substitution optimizes the balance between spin–orbit interaction strength and band inversion. Furthermore, the similarity between the analytical Rashba model and the DFT band structure confirms that the simplified approach captures the essential physics of spin–orbit-driven band structures while enabling much faster computation transport tensors compared to full DFT calculations. The observed trends can be traced to the strength of Rashba coupling and the nature of impurity scattering. Increasing α_R enhances Rashba splitting, which enhances the in-plane Edelstein response. The numerical results demonstrate that spin–charge conversion peaks near x = 0.6, instead of DFT-derived parameters within the semiclassical BTE framework. While the model neglects electron–phonon coupling and interfacial roughness effects, it provides a reliable first-order description of spin–orbit-driven transport at moderate impurity concentrations.

Effect of Sb doping on band inversion and surface-state localization in BST and BSS systems

To assess whether the observed band modifications influence the topological character, Fig. 6 presents the p-orbital projections for Bi, Sb, and Te in BST and BSS at Sb doping levels x = 0.1, 0.5, and 0.9. Band inversion persists across all doping levels [see Fig. S7 within the SI], although its strength gradually diminishes with increasing Sb concentration, particularly near the Γ-point. This weakening suggests a progressive evolution toward a trivial insulating phase, yet without a complete suppression of the TI characteristics. A similar trend is observed in BSS, where band inversion is retained but noticeably weakened upon Sb substitution. These results indicate a continuous modulation of topological order rather than an abrupt topological phase transition in both systems. At the same time, we computed the time-reversal Z₂ index of the 5-QL slabs using the Wilson-loop (hybrid Wannier charge-center) approach implemented in Z2Pack,⁷⁰ interfaced with VASP including spin–orbit coupling and spinor wavefunctions. The calculations were performed by fixing k_z = 0 and evaluating closed loops along k_x while sweeping k_y from 0 to 1/2. We first established the topological character of the 5-QL monolayers. The Z₂ index distinguishes a topological insulator (ν = 1) from a trivial insulator (ν = 0) under time-reversal symmetry. In thin slabs, a two-dimensional Z₂ invariant is well defined when the in-plane Brillouin zone remains gapped. Our results indicate that representative (Bi_1−xSb_x)₂Te₃ with higher Sb concentrations (x = 0.5, 0.6, 0.7, 0.8, 0.9) exhibits a Z₂ invariant equal to 1, confirming that they are topological insulator character.

Fig. 6 Orbital-projected electronic band structures of 5QL (a) (Bi_1−xSb_x)₂Te₃ and (b) (Bi_1−xSb_x)₂Se₃ for x = 0.1, 0.5, and 0.9. The contributions from Bi, Sb, and chalcogen (Te or Se) p-orbitals are shown in blue, green, and red, respectively.

To determine whether the weakening of band inversion manifests in the topological properties, the contributions of surface and bulk atoms to the formation of surface states (SS) were examined. The crystal structures were spatially partitioned: the top and bottom QLs were designated as the surface region, while the central three QLs represented the bulk. This classification enables a focused evaluation of the Dirac cone's spatial origin, based on the assumption that surface states predominantly arise from the outermost atomic layers.

Fig. 7(a) displays the projected p-orbital band structures for the surface and bulk regions of BST at x = 0.1, 0.5, and0.9. The corresponding results for other doping levels in both BST and BSS are provided in Fig. S8 and S9 within the SI for BST and BSS, respectively. As shown, the surface atoms clearly dominate the formation of surface states at all concentrations. This finding is further corroborated by Fig. 7(b), which illustrates the projections onto the bulk atoms, revealing minimal or negligible contribution to the SS.

Fig. 7 Orbital-weighted band structures of 5QL (Bi_1−xSb_x)₂Te₃ for x = 0.1, 0.5, and 0.9, showing (a) surface contributions from the top and bottom quintuple layers (QLs) and (b) bulk contributions from the three middle QLs. Red dots indicate states with dominant surface or bulk character, as visualized in the structural schematics on the right.

This distinction underscores a defining characteristic of topological insulators: while the bulk remains insulating, the metallic surface states are confined to the system's boundaries. The consistent absence of significant bulk contribution across all doping levels affirms that the Dirac-like surface states originate from topological protection enforced by crystal symmetry and strong spin–orbit coupling.

To better understand the transport behavior of the studied materials, carrier mobilities were calculated separately for surface and bulk states based deformation potential theory. The results reveal a clear and systematic distinction between surface and bulk mobilities across all compositions, with surface mobilities consistently exceeding their bulk counterparts by a significant margin. This indicates that charge transport is predominantly governed by surface states rather than bulk conduction. The enhanced mobility at the surface suggests more favorable conditions for carrier movement, likely due to reduced scattering and a distinct band structure. In contrast, the bulk region exhibits much lower mobility values, reflecting limited carrier transport efficiency. These findings consistently highlight the dominant role of surface channels in the overall electronic transport behavior.

Conclusions

In this study, we performed a systematic first-principles investigation of the structural, electronic, and topological properties of Sb-doped Bi₂Te₃ and Bi₂Se₃ thin films. By exploring a comprehensive range of doping concentrations within a 5QL framework, we elucidated how Sb substitution modulates lattice parameters, band inversion, Rashba-type spin splitting, and the spatial origin of surface states. Our results reveal that Sb doping introduces non-monotonic changes in structural asymmetry and spin–orbit interaction strength, particularly in Bi₂Te₃, where a pronounced Rashba effect emerges at intermediate doping levels (notably at x = 0.6) with a peak Rashba parameter of α_R = 2.098 eV Å.

Furthermore, we demonstrated that the topological surface states (TSS) persist across the full doping range, albeit with a gradual weakening of band inversion, indicating a continuous rather than abrupt topological phase evolution. Spin texture analysis confirmed helical spin-momentum locking and clear signatures of inner and outer spin-split bands, which are typical of Rashba systems, substantiating the coexistence of Rashba-type interaction and topological protection in these materials.

This study offers new perspectives on Sb-induced modifications in well-known topological insulators but also provide a foundation for controlling key properties such as Rashba spin splitting and the stability of topological surface states across varying compositions. Our carrier mobility analysis further confirms that surface transport channels dominate over bulk conduction at all doping levels, reinforcing the critical role of surface states in practical charge transport. Together, these findings position (Bi_1−xSb_x)₂Te₃ and (Bi_1−xSb_x)₂Se₃ as highly tunable platforms, paving the way for advanced applications in spin-orbitronics and topological quantum computing.

Author contributions

M. C. O. conceived and supervised the project. R. K. and B. A. H. performed the calculations. R. K. and B. A. H. and K. A. and A. I. S. and A. J. performed the theoretical analysis. All authors contributed to the writing of the manuscript.

Conflicts of interest

The authors declare no competing interests.

Data availability

Our article does not require legal or ethical confidentiality requirements.

Data for this article, including Fig. 2 are available at VESTA at https://jpminerals.org/vesta/en/download.html.

Data for this article, including Fig. 1 and 4–7 are available at VASP (Vienna Ab initio Simulation Package) at https://www.vasp.at/.

The illustration shown in Fig. 5(d) was drawn using Adobe Illustrator at https://www.adobe.com/products/illustrator.html.

The videos used in the supplementary information (SI) were recorded with ASE (Atomic Simulation Environment) at https://ase-lib.org/about.html.

The data supporting this article have been included as part of the SI. Supplementary Information provides details on stable atomic configurations, additional electronic band structures, Rashba-related transport (Edelstein tensor and conductivity), dipole moment and LDOS analyses, deformation-potential calculations, AIMD thermal stability tests, and surface-bulk orbital projections supporting the main text. Supplementary information is available. See DOI: https://doi.org/10.1039/d5nr03682h.

Acknowledgements

This work was performed using the KUACC HPC Cluster and the High Performance Grid Computing Center (TR-Grid e-Infrastructure) at TÜBİTAK ULAKBİM. This study was supported by the European Research Council (ERC) under the Starting Grant SKYNOLIMIT (Grant No. 948063), the ERC Proof of Concept Grant SuperPHOTON (Grant No. 101100718), the European Association of National Metrology Institutes (EURAMET) under contract number 23FUN07 and project acronym QuAHMET, and the Scientific and Technological Research Institution of Türkiye (TÜBİTAK) within TÜBİTAK–Chinese Academy of Sciences (CAS) bilateral collaboration program under contract number 122N485 and project acronym 2DThermoelectrics. This article is based upon work from COST Action POLYTOPO CA23134, supported by COST (European Cooperation in Science and Technology).

References

1. M. Z. Hasan and C. L. Kane, Rev. Mod. Phys., 2010, 82, 3045–3067.
2. Z.-K. Lin, Q. Wang, Y. Liu, H. Xue, B. Zhang, Y. Chong and J.-H. Jiang, Nat. Rev. Phys., 2023, 5, 483–495.
3. X.-L. Qi and S.-C. Zhang, Rev. Mod. Phys., 2011, 83, 1057–1110.
4. P. Singha, S. Das, N. Rana, S. Mukherjee, S. Chatterjee, S. Bandyopadhyay and A. Banerjee, J. Appl. Phys., 2024, 135, 025001.
5. O. V. Yazyev, J. E. Moore and S. G. Louie, Phys. Rev. Lett., 2010, 105, 266806.
6. A. Tarasov, D. Estyunin, A. Rybkin, A. Frolov, A. Sergeev, A. Eryzhenkov, V. Anferova, T. Estyunina, D. Glazkova and K. Kokh, et al., Phys. Rev. B, 2025, 111, 165115.
7. Y. Xia, D. Qian, D. Hsieh, L. Wray, A. Pal, H. Lin, A. Bansil, D. H. Y. S. Grauer, Y. S. Hor, R. J. Cava and M. Z. Hasan, Nat. Phys., 2009, 5, 398–402.
8. Y. A. Bychkov and É. I. Rashba, JETP Lett., 1984, 39, 78.
9. G. Bihlmayer, O. Rader and R. Winkler, New J. Phys., 2015, 17, 050202.
10. Q. Peng, Y. Lei, X. Deng, J. Deng, G. Wu, J. Li, C. He and J. Zhong, Phys. E: Low-Dimens. Syst. Nanostructures, 2022, 135(114944).
11. G. Lazrak, B. Göbel, A. Barthélémy, I. Mertig, A. Johansson and M. Bibes, Phys. Rev. Res., 2024, 6, 023074.
12. W. Xue, J. Li, X. Peng, C. He, T. Ouyang, X. Qi, C. Zhang, C. Luo, J. Deng and Q. Peng, et al., J. Electron. Mater., 2022, 51, 5142–5149.
13. D. Holgado, K. Bolaños, S. De Castro, H. Monteiro, F. Pena, A. Okazaki, C. Fornari, P. Rappl, E. Abramof and D. Soares, et al., Appl. Phys. Lett., 2020, 117, 102108.
14. C. Mondal, C. Barman, A. Alam and B. Pathak, APS March Meeting Abstracts, 2021, pp. E45-014.
15. H. M. Benia, C. Lin, K. Kern and C. R. Ast, Phys. Rev. Lett., 2011, 107, 177602.
16. W. Liu, L. He, Y. Xu, K. Murata, M. C. Onbasli, M. Lang, N. J. Maltby, S. Li, X. Wang and C. A. Ross, et al., Nano Lett., 2015, 15, 764–769.
17. Y. Nakazawa, T. Akiho, K. Kanisawa, H. Irie, N. Kumada and K. Muraki, Appl. Phys. Lett., 2024, 125, 083101.
18. J.-M. Zhang, W. Ming, Z. Huang, G.-B. Liu, X. Kou, Y. Fan, K. L. Wang and Y. Yao, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 88, 235131.
19. B. Y. Yavorsky, N. Hinsche, I. Mertig and P. Zahn, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 84, 165208.
20. M. Michiardi, F. Boschini, H.-H. Kung, M. Na, S. Dufresne, A. Currie, G. Levy, S. Zhdanovich, A. Mills and D. Jones, et al., Nat. Commun., 2022, 13, 3096.
21. Y. Chen, J. G. Analytis, J.-H. Chu, Z. Liu, S.-K. Mo, X.-L. Qi, H. Zhang, D. Lu, X. Dai and Z. Fang, et al., Science, 2009, 325, 178–181.
22. H. Zhang, C.-X. Liu, X.-L. Qi, X. Dai, Z. Fang and S.-C. Zhang, Nat. Phys., 2009, 5, 438–442.
23. Y. Xia, D. Qian, D. Hsieh, L. Wray, A. Pal, H. Lin, A. Bansil, D. Grauer, Y. S. Hor and R. J. Cava, et al., Nat. Phys., 2009, 5, 398–402.
24. T. Kato, H. Kotaka and F. Ishii, Mol. Simul., 2015, 41, 892–895.
25. H. Osterhage, J. Gooth, B. Hamdou, P. Gwozdz, R. Zierold and K. Nielsch, Appl. Phys. Lett., 2014, 105, 123117.
26. I. Aguilera, C. Friedrich, G. Bihlmayer and S. Blügel, Phys. Rev. B: Condens. Matter Mater. Phys., 2013, 88, 045206.
27. G. Wang, X.-G. Zhu, Y.-Y. Sun, Y.-Y. Li, T. Zhang, J. Wen, X. Chen, K. He, L.-L. Wang and X.-C. Ma, et al., Adv. Mater., 2011, 23, 2929.
28. X. He, H. Li, L. Chen and K. Wu, Sci. Rep., 2015, 5, 8830.
29. P. Das, T. Whitcher, M. Yang, X. Chi, Y. Feng, W. Lin, J. Chen, I. Vobornik, J. Fujii and K. Kokh, et al., Phys. Rev. B, 2019, 100, 115109.
30. C. Niu, Y. Dai, Y. Zhu, Y. Ma, L. Yu, S. Han and B. Huang, Sci. Rep., 2012, 2, 976.
31. S. K. Chong, L. Liu, K. Watanabe, T. Taniguchi, T. D. Sparks, F. Liu and V. V. Deshpande, Nat. Commun., 2022, 13, 6386.
32. C.-X. Liu, X.-L. Qi, H. Zhang, X. Dai, Z. Fang and S.-C. Zhang, Phys. Rev. B: Condens. Matter Mater. Phys., 2010, 82, 045122.
33. N. Bahadursha, B. Sadhukhan, T. Nag, S. Bhattacharya and S. Kanungo, Mater. Sci. Semicond. Process., 2025, 186, 109042.
34. D. Kong, Y. Chen, J. J. Cha, Q. Zhang, J. G. Analytis, K. Lai, Z. Liu, S. S. Hong, K. J. Koski and S.-K. Mo, et al., Nat. Nanotechnol., 2011, 6, 705–709.
35. K. Ishizaka and M. Bahramy, et al., Nat. Mater., 2011, 10, 521–526.
36. P. King, R. C. Hatch, M. Bianchi, R. Ovsyannikov, C. Lupulescu, G. Landolt, B. Slomski, J. Dil, D. Guan and J. Mi, et al., Phys. Rev. Lett., 2011, 107, 096802.
37. A. Johansson, J. Phys.: Condens. Matter, 2024, 36, 423002.
38. G. Kresse and J. Hafner, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 558.
39. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169.
40. J. P. Perdew, K. Burke and M. Ernzerhof, Phys. Rev. Lett., 1996, 77, 3865.
41. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953.
42. H. J. Monkhorst and J. D. Pack, Phys. Rev. B, 1976, 13, 5188.
43. N. Bovenzi and M. Diez, Phys. Rev. B, 2017, 95, 205430.
44. A. Johansson, J. Henk and I. Mertig, Phys. Rev. B, 2016, 93, 195440.
45. Z. Zhu, Y. Cheng and U. Schwingenschlögl, Phys. Rev. B: Condens. Matter Mater. Phys., 2012, 85, 235401.
46. Y.-Y. Li, G. Wang, X.-G. Zhu, M.-H. Liu, C. Ye, X. Chen, Y.-Y. Wang, K. He, L.-L. Wang and X.-C. Ma, et al., Adv. Mater., 2010, 22, 4002–4007.
47. J. M. Crowley, J. Tahir-Kheli and W. A. Goddard III, J. Phys. Chem. Lett., 2015, 6, 3792–3796.
48. S. Murakami, Phys. E, 2011, 43, 748–754.
49. W. Liu, X. Peng, C. Tang, L. Sun, K. Zhang and J. Zhong, Phys. Rev. B: Condens. Matter Mater. Phys., 2011, 84, 245105.
50. Q. Zhang, Z. Zhang, Z. Zhu, U. Schwingenschlogl and Y. Cui, ACS Nano, 2012, 6, 2345–2352.
51. V. Kumar and N. Ganguli, Phys. Rev. B, 2022, 106, 125127.
52. D. C. Vaz, F. Trier, A. Dyrdał, A. Johansson, K. Garcia, A. Barthelemy, I. Mertig, J. Barnaś, A. Fert and M. Bibes, Phys. Rev. Mater., 2020, 4, 071001.
53. S. LaShell, B. McDougall and E. Jensen, Phys. Rev. Lett., 1996, 77, 3419.
54. Z.-H. Zhu, G. Levy, B. Ludbrook, C. Veenstra, J. Rosen, R. Comin, D. Wong, P. Dosanjh, A. Ubaldini and P. Syers, et al., Phys. Rev. Lett., 2011, 107, 186405.
55. Y. M. Koroteev, G. Bihlmayer, J. Gayone, E. V. Chulkov, S. Blügel, P. M. Echenique and P. Hofmann, Phys. Rev. Lett., 2004, 93, 046403.
56. A. Varykhalov, D. Marchenko, M. Scholz, E. Rienks, T. Kim, G. Bihlmayer, J. Sánchez-Barriga and O. Rader, Phys. Rev. Lett., 2012, 108, 066804.
57. K. Lee, D. Mou, N. H. Jo, Y. Wu, B. Schrunk, J. M. Wilde, A. Kreyssig, A. Estry, S. L. Bud'ko and M. C. Nguyen, et al., Phys. Rev. B, 2021, 103, 085125.
58. C. R. Ast, J. Henk, A. Ernst, L. Moreschini, M. C. Falub, D. Pacilé, P. Bruno, K. Kern and M. Grioni, Phys. Rev. Lett., 2007, 98, 186807.
59. J. Bardeen and W. Shockley, Phys. Rev., 1950, 80, 72.
60. J. Xiao, M. Long, X. Zhang, J. Ouyang, H. Xu and Y. Gao, Sci. Rep., 2015, 5, 9961.
61. D. Mo, W. Wang and Q. Cai, Nanoscale Res. Lett., 2016, 11, 354.
62. H.-T. He, G. Wang, T. Zhang, I.-K. Sou, G. K. L. Wong, J.-N. Wang, H.-Z. Lu, S.-Q. Shen and F.-C. Zhang, Phys. Rev. Lett., 2011, 106, 166805.
63. N. Koirala, M. Brahlek, M. Salehi, L. Wu, J. Dai, J. Waugh, T. Nummy, M.-G. Han, J. Moon and Y. Zhu, et al., Nano Lett., 2015, 15, 8245–8249.
64. A. Tyagi, V. Mišeikis, L. Martini, S. Forti, N. Mishra, Z. M. Gebeyehu, M. A. Giambra, J. Zribi, M. Frégnaux and D. Aureau, et al., Nanoscale, 2022, 14, 2167–2176.
65. L. Banszerus, M. Schmitz, S. Engels, J. Dauber, M. Oellers, F. Haupt, K. Watanabe, T. Taniguchi, B. Beschoten and C. Stampfer, Sci. Adv., 2015, 1, e1500222.
66. Z. Li, J. Wang and Z. Liu, J. Chem. Phys., 2014, 141, 144107.
67. M. T. Pettes, J. Kim, W. Wu, K. C. Bustillo and L. Shi, APL Mater., 2016, 4, 104810.
68. R. Sun, S. Yang, X. Yang, E. Vetter, D. Sun, N. Li, L. Su, Y. Li, Y. Li and Z.-z. Gong, et al., Nano Lett., 2019, 19, 4420–4426.
69. Y. Wang, B. Liu, Y.-X. Huang, S. V. Mambakkam, Y. Wang, S. A. Yang, X.-L. Sheng, S. A. Law and J. Q. Xiao, Phys. Rev. B, 2022, 106, L241401.
70. D. Gresch, G. Autes, O. V. Yazyev, M. Troyer, D. Vanderbilt, B. A. Bernevig and A. A. Soluyanov, Phys. Rev. B, 2017, 95, 075146.

---

Footnote

† These authors contributed equally to this work as co-first authors.

---