Quantum oscillations in 2D electron gases with spin-orbit and Zeeman interactions

Shubnikov-de Haas (SdH) oscillations have served as a paradigmatic experimental probe and tool for extracting key semiconductor parameters such as carrier density, effective mass, Zeeman splitting with g-factor $g^*$, quantum scattering times and spin-orbit (SO) coupling parameters. Here, we derive for the first time an analytical formulation for the SdH oscillations in 2D electron gases (2DEGs) with simultaneous Rashba, Dresselhaus, and Zeeman interactions. Our analytical and numerical calculations allow us to extract both Rashba and Dresselhaus SO coupling parameters, carrier density, quantum lifetimes, and also to understand the role of higher harmonics in the SdH oscillations. More importantly, we derive a simple condition for the vanishing of SO induced SdH beatings for all harmonics in 2DEGs: $\alpha/\beta= [(1-\tilde \Delta)/(1+\tilde \Delta)]^{1/2}$, where $\tilde \Delta$ is a material parameter given by the ratio of the Zeeman and Landau level splitting. We also predict beatings in the higher harmonics of the SdH oscillations and elucidate the inequivalence of the SdH response of Rashba-dominated ($\alpha>\beta$) vs Dresselhaus-dominated ($\alpha<\beta$) 2DEGs in semiconductors with substantial $g^*$. We find excellent agreement with recent available experimental data of Dettwiler ${\it et\thinspace al.}$ Phys. Rev. X $\textbf{7}$, 031010 (2017), and Beukman ${\it et\thinspace al.}$, Phys. Rev. B $\textbf{96}$, 241401 (2017).


I. INTRODUCTION
The spin-orbit (SO) interaction couples the orbital and spin degrees of freedom, not only forms the basis for a range of spin related effects such as the spin Hall effect 1-4 and the persistent spin helix [5][6][7] , but also underlies the physical mechanisms of new phases of matter, e.g., topological insulators, quantum spin Hall materials [8][9][10] , and Majorana [11][12][13] , Dirac and Weyl fermions 14 .Accordingly, advancing techniques and methods to measure and extract SO couplings from experimental data are crucial for the development of these fields.
Shubnikov-de Haas (SdH) oscillations 16,17 are among the best techniques to probe simultaneously spin-and charge-related quantities associated to electrons in semiconductors, including effective masses, gyromagnetic ratios, quantum scattering times, densities and SO couplings.Most recently, they have been crucial to the study and understanding of new materials, as for example, 2D-materials, transition metal dichalcogenides, van der Waals heterostructures [18][19][20][21][22][23][24][25] , and also materials hosting new phases of matter e.g., topological insulators 26 , unconventional superconductivity 27 and correlated insulator behavior 28 .It has also been used to establish the presence of nodal-lines 29 , Berry's phase 30,31 , and differ-ent topology of Fermi surfaces 32 .SdH oscillations are magneto-oscillations in the resistivity and originate from the sequential crossings of the discrete Landau Levels (LLs) through the Fermi energy.Without SO coupling and in the low-field regime, the period of the SdH oscillations can be related to the density of the electron gas 33 .In the presence of SO interaction, on the other hand, the energy spectrum changes dramatically thus leading to additional frequencies in the magnetoresistivity and hence beatings, Figs.1(a).This was first theoretically described semiclassically by Das et.al., 34 .In the socalled Onsager's picture, different sub-bands possess different Sommerfeld quantized orbits (playing the role of the LLs), which cross the Fermi energy with different frequencies in B −1 .The spin-split bands give rise to two distinct oscillating frequencies in the magnetotransport.The standard experiment relies on Fourier analyzing the measured SdH oscillations.An experimental method introduced in Refs.35-37 has often been used to estimate the strength of the Rashba coupling via the splitting of the Fourier frequency peaks.However, these methods have been criticized for not accounting for the Zeeman splitting (through the g-factor g * ) nor for the additional Dresselhaus SO coupling 15 .
There have been some attempts to analyze the SdH os- cillations taking into account both α, β and g * .However, these mostly involved qualitative comparison with the energy spectrum of pure Rashba and pure Dresselhaus 38,39 .In Ref. 40, fully numerical calculations of magnetooscillations were performed but for relatively high magnetic fields and low electron densities, far away from the regime of recent experimental works 41 .Moreover, it was realized that in the absence of the Zeeman interaction, important features are absent.More specifically, without accounting for the spin mixing generated by the magnetic field (via the Zeeman interaction), predictions become imprecise 42 , and even fail to describe phenomena such as magnetic inter-subband scattering 43 and magnetic breakdown 44 .In general, full quantum mechanical numerics are generally done in order to check agreement with experiments, which are neither very practical nor elucidate much of the physics happening in those systems 41,45 .Finally, all the previous works have neglected the influence of higher harmonics, recently seen experimentally 46 .
Here, we present a detailed investigation of SdH oscillations in the presence of SO couplings of both Rashba α and Dresselhaus β types and Zeeman interaction with g-factor g * .Our main result is the derivation, for the first time in the literature, of a simple analytical expression for the SdH oscillations in the presence of simultaneous arbitrary couplings α and β in addition to g * .We note that earlier analytical descriptions of SdH magnetoresistivity oscillations considered particular cases, namely, when either only one of the parameters α, β or g * was nonzero or any two of these parameters were nonzero, with the exceptions (α = β = 0, g * = 0) and (α = β, g * = 0).Interestingly, our analytical formula generalizes previous results 44 for g * = 0 and predicts a new condition for the vanishing of the SdH magneto-oscillation beatings in all harmonics [e.g., Figs.1(c)] in Rashba-Dresselhaus coupled 2DEGs with substantial Zeeman splittings, namely, where ∆ is a material parameter given by the ratio between the Zeeman splitting and the Landau-level spacing.
As we discuss later on, Eq. ( 1) is not associated with a conserved quantity in our system; this contrasts with the persistent-spin-helix condition α = β, which predicts spin conservation along particular axes [5][6][7] .We stress that this case with α = β and generic g * = 0 leads to beating in the frequency spectrum of our system, Figs.1(d), as opposed to our new condition in Eq. 1.As we discuss below, our numerical and analytical approaches show excellent agreement with available data from Refs.41 and 46.
Our approach combines a semiclassical formulation for the resistivity of 2DEGs with a trace formula for the density of states (DOS) in a quantizing magnetic field.The trace formula expresses the DOS using the usual Poisson summation formula 47 .This formulation brings out the oscillatory part of the DOS, thus allowing us to clearly identity the higher harmonics of the SdH oscillations.It enables us to conveniently separate the frequency scales into "fast" and "slow" oscillations thus allowing for a clearer interpretation of the underlying physical phenomena, e.g., the slow beating SdH oscillations due to the SO coupling.
Our main results for the oscillatory part of magnetoresistivity δρ xx (1/B) and its frequency spectra I(f ) [panel insets] are show in Fig. 1.For pure Rashba [α = 0, β = 0, Fig. 1a)] and pure Dresselhaus [α = 0, β = 0, Fig. 1b)], but non-zero Zeeman term (g * = 0), the frequency spectra, as usual, show two main peaks, which correspond to the first two Fourier components of δρ xx (1/B).These two cases, however, exhibit a marked contrast: while the pure Rashba shows a peak splitting at the fundamental frequency, the pure Dresselhaus exhibits a peak splitting in the second harmonic.As we explain in detail in Sec.V D, this contrasting behavior arises from the interplay between the Zeeman and SO interactions, which makes the SdH magneto-responses with nonzero gfactors g * inequivalent for Rashba-dominated (α > β) vs. Dresselhaus-dominated (α < β) 2DEGs.For g * = 0, the pure Rashba and pure Dresselhaus cases give identical results.
Figure 1(c) illustrates our prediction in Eq. (1) thus showing no peak splitting in the frequency spectraat any harmonic -when this condition is satisfied.To emphasize this condition emulates a situation with no SO coupling (i.e., no beating), we plot in Fig. 1(c) the α = β = 0 (with g * = 0) case [dashed curve in 1(c)], which shows complete overlap with the case satisfying Eq. (1).In contrast and for completeness, Fig. 1d) shows the α = β = 0 case with g * = 0, which exhibits peak splitting in the second harmonic.
We have applied our analytical description to lowdensity GaAs-based quantum wells for which there are experimental data 46 showing several harmonics in the SdH magneto-oscillations. Figure 2 shows the excellent agreement obtained, thus illustrating that our semiclassical formulas can satisfactorily capture the higher harmonics of the SdH oscillations.Moreover, we have applied our analytical approach to low-density InSb-based 2DEGs 15,39 where, unlike GaAs-based 2DEGs, a strong SO coupling manifests itself as beatings in the measured SdH oscillations, and find good agreement.We have also implemented a detailed numerical calculation for highdensity InAs-based 2DEGs for which an analytical description is not adequate.Here again we find very good agreement with available data 41 and are able to extract SO coupling parameters.
Next (Sec II), we present a description of the Hamiltonian of our system.In Sec.III we discuss how to obtain the "F -function", the central quantity in our formulation, from the Landau-quantized energy spectrum of our system and its connection with the density of states (DOS).The formalism for obtaining the Shubnikov-de Haas oscillations in terms of the Poisson summation formula and the F-function is described in Sec.IV.Finally, in Sec.V we present and analyze different particular cases of SdH oscillations and, more important, derive the new condition in Eq. (1) for the complete absence of beatings (all harmonics) in the SdH oscillations, for 2DEGs with non-zero Rashba, Dresselhaus, and Zeeman couplings.The appendices present relevant details of our theoretical formulation.

II. 2DEG HAMILTONIAN
Our starting point is the Hamiltonian for a 2DEG confined in a quantum well (xy plane) grown along the [001] crystallographic direction, taken as z axis.In the presence of a perpendicular external magnetic field B = (0, 0, B) and both Rashba 48 and Dresselhaus 49 spin orbit interactions, the Hamiltonian reads where g * is the g-factor, m * is effective mass, Π = p−qA is the canonical momentum, q is the electric charge, µ B is the Bohr magneton, the reduced Planck's constant and σ x , σ y , σ z denote the usual Pauli matrices.The parameters α and β denote the linear-in-k Rashba and Dresselhaus SO couplings, respectively.The β coupling includes a density dependent correction arising from the cubic Dresselhaus term.As we discuss later on [Sec.VI A], our numerical results will account for the full cubic Dresselhaus term.
Let us introduce the annihilation and creation operators associated to the Landau level |n , with the magnetic length and the center of the Landau orbit denoted by c = |qB| and y 0 = ekx |qB| , respectively.In this work, we have q = −e, where e > 0 is the absolute value of the elementary electronic charge, and we choose B > 0, yielding ζ = 1.Using Eqs. ( 3) and (4), our Hamiltonian [Eq.( 2)] becomes where the cyclotron frequency is ω c = eB/m * , ∆ = g * µ B B, which inherits its sign from g * , and σ ± = σ x ± iσ y , with σ x and σ y denoting Pauli matrices.We now perform the canonical transformation H = UHU † with U = e −i π 4 ( σz 2 +a † a) , which yields H where we have introduced the real valued, dimensionless Analytical solutions for the above Hamiltonian [Eq.(9)] can be found for the cases with either pure Rashba or pure Dresselhaus 48,50 .The specific cases of α = ±β turn out to be of great physical interest, where persistent spin helix (PSH) 5,6,46 and persistent Skyrmion lattice (PSL) 7 were predicted.Interestingly, the case with α = ±β maps to the Rabi model in quantum optics and was recently solved exactly 51 .The exact solution relies on obtaining zeros of a transcendental function.Moreover, previous studies of the Rabi model have important implications for our system.For instance, we have shown that the Rabi parity symmetry 51,52 remains valid in our problem for arbitrary α and β (See Appendix B).This enables us to separate the Hilbert space in two subspaces with different parities, which can be individually analyzed and compared.As for general couplings α and β, similar systems have been studied before in the framework of Landau levels, using either variational (Hartree-Fock) methods 53 , second order perturbation 54,55 or obtaining the spectrum in terms of solutions of transcendental equations 56 .A perturbation scheme based on 4th order Schrieffer-Wolff transformation has also been used to find approximate analytical solutions 57 .However, we are unaware of any exact analytical solution for general Rashba, Dresselhaus and Zeeman coupling.

III. F -FUNCTION AND ITS CONNECTION WITH THE ENERGY SPECTRUM AND DOS
For our 2DEG in the presence of perpendicular magnetic field, the low magnetic field regime corresponds to having a very large number of Landau levels below the Fermi energy ε F (taken as constant and equal to its zero-field value), i.e., many occupied states.The system is thus assumed to be far away from the integer quantum Hall regime where few Landau levels are occupied and the effects of electron-electron interaction become important 33 .Let us denote the eigenenergies of our dimensionless Hamiltonian Eq. ( 9) by ε n,s , where n ∈ N 0 represents the LL number and s = ± represents a pseudospin associated to the presence of two spin-split bands (due to the Zeemann and SO interactions).With this notation, the density of states (DOS) reads where D = A/2π 2 c is the LL degeneracy and A the 2DEG area.This LL degeneracy is the same for all 2DEGs studied here in the presence or absence of Zeemann and SO interactions.
As we show in the next section, the magnetotransport properties of the system can be determined by the Landau levels sequentially crossing ε F .The rate at which these crossings happen will determine a periodic behavior of the magnetotransport properties of the system as the magnetic field is varied.In order to describe this periodicity, we introduce the F -function 33 (see Appendix A for details), which is defined by the relation The F -function gives the Landau level index n of the state that has energy ε and pseudo-spin s at magnetic field B. 47,58 It is important to notice that the equation for n [Eq.(11)] can also return non-integer values for n.
In such cases the F -function provides a measure of how close a Landau level n is to the energy ε, for a given pseudo-spin s and magnetic field B.
Since one can relate transport phenomena with the density of states, we rewrite the DOS of our system in a way that highlights its oscillatory behavior dependence on both α and β.First we introduce the F s function into Eq.( 10) we obtain where D 0 = m * 2π 2 is the (constant) density of states per spin for the 2DEG at zero magnetic field (see Appendix A for details).As we are going to see later, F + contains the fast oscillations with respect to 1/B, which is is proportional to the electron density n 2D .On the other hand, F − contains the slow oscillations that are determined by the spin-orbit coupling terms, α and β.Moreover, the fast oscillations arising from the terms with l > 1 correspond to the higher harmonics, and have be seen in experiments 46 .

IV. SDH OSCILLATIONS IN THE MAGNETORESISITIVITY
As already mentioned, the oscillations in the magnetoresistivity as a function of the magnetic field are called SdH oscillations 33 .They appear as a consequence of the sequential depopulation of the LLs when the magnetic field is increased.For low magnetic fields where multiple LL are occupied, i.e., far from the integer quantum Hall regime 33 , a semi-classical description of the magnetooscillations can be used.
In experiments, the measurement of the SdH oscillations is accessed via the longitudinal differential resistivity.In general, the resistivity tensor is defined as the inverse matrix of the conductivity tensor, The relevant magnetoresistivity component for us is where where f 0 (ε) is the Fermi-Dirac distribution.Using a semi-classical approach, we account for the magnetic field dependence of the conductivity via the electron scattering time τ (ε, B), which is proportional to the DOS D(ε, B) via Fermi's golden rule.Accordingly, up to linear order on the deviation of the DOS, we obtain with D 0 (ε) = D(ε, B = 0) and τ 0 (ε) = τ (ε, B = 0).Using the Drude semi-classical equations for the frequencyindependent current 33 , the normalized longitudinal resistivity reads For the DOS in the presence of Landau level broadening due to scattering processes, the relation in Eq. ( 10) is replaced by where L Γ (x) describes the broadening function, e.g., Lorentzian or Gaussian, and Γ is parameter defining the broadening of the levels (see Appendix A for details).
After applying the Poisson summation formula, we obtain a result that resembles Eq. ( 14), apart from the appearance of the the cosine Fourier transform of L Γ (x), denoted with LΓ (x), The so-called Dingle factor LΓ (x) 33 sets the limit of validity of the semi-classical approximation, i.e., that the oscillatory part of the resistivity should be much smaller than the constant term.It also gives the regime where it is valid to consider only the lowest harmonic.
Higher harmonics have been observed in magnetoresistivity measurements 46 in GaAs-based 2DEGs.The F − function can be related to the envelope of the SdH oscillations.The general form of the temperature-dependent normalized resistivitity reads × cos(2πlF − ) cos(2πlF + ).
Even though we only consider the zero-temperature limit in the present work, for completeness, below we present the temperature-dependence of δρ xx (B, T ) valid in the relevant parameter range considered in this work and for all the systems studied here.As show in Appendix G, we find where the temperature-dependent coefficient accounts for the temperature dependence of the SdH oscillations.In the limit of vanishing α and β this reduces to the result in Ref. 59, and in the case of both vanishing β and broadening (Γ = 0) gives the result in Ref. 50  [Eq.(9.28)].Here we assume that ε is close to the zero magnetic field Fermi energy ε F = 2 k 2 F /2m.A widely used method to extract spin-orbit couplings and electronic densities is to analyze the oscillations by calculating the quantity which defines the power spectrum of the normalized magneto-resistivity with a trivial background value ρ xx (B = B 1 ) removed.Note that B 1 should be small enough such that the semiclassical regime of a constant ρ xx (B → 0) is reached.
In Fig. 2 the power spectrum is shown for data from Fig. S11a in Ref. 46, where magnetoresistivity SdH oscillations were measured in a GaAs 2DEG over a magnetic field interval [0.20, 1.5] T. The power spectrum shows a SdH peak at f ≈ 10.5 T (the fundamental frequency), and higher harmonics are clearly visible at 21.0 T and 31.5 T, corresponding to the first and second harmonic, respectively.The experimental data was fitted with Eq. (30) with one fit parameter: τ q .The resulting fit matches very well the harmonics of the SdH signal.To acccount for the small background shift in the experimental data as seen in the inset a more elaborate modeling of the data FIG.2. The power spectrum I(f ) for δρxx measurements on a GaAs 2DEG in Ref. 46 obtained using Eq. ( 26).The calculated results used Eq. ( 30) with one fitting parameter: τq.
The inset shows the magneto-resistivity data and the corresponding calculated δρxx.
would be required.The fitting was done using six harmonics, and resulted in τ q = 0.8 ps, using standard GaAs parameters m = 0.067m 0 and g * = −0.44.Note that we have used Eq.( 30), which does not include SO coupling, for our fitting procedure here.This is justifiable because GaAs-based 2DEGs have relatively small SO couplings, not accessible via SdH measurements.Weak antilocalization measurements can access the SO parameter in these systems 41 .However, GaAs-based 2DEGs have relatively high mobilities thus making it possible to see many harmonics.

V. RESULTS AND DISCUSSIONS
In this section we present the energy spectrum, Ffunction and magnetoresistivity SdH oscillations for different parameter regimes of our Hamiltonian, Eq. ( 9).Additionally, we discuss in detail the interpretation of the SdH oscillations within the trace formula description (e.g., contribution of higher harmonics) and show how to extract relevant spin-orbit couplings from it.The results are presented in order of simplicity, i.e., from the simplest to the more complex case.

A. Landau levels with only Zeeman interaction
In the presence of Zeeman and no Rashba and Dresselhaus SO couplings, i.e., α = β = 0, the eigenenergies of our Hamiltonian [Eq.( 9)] are given by with n ∈ N 0 and s = 1 (s = −1) representing the pure spin state |↑ (|↓ ).In Fig. 3 we plot the four energy levels corresponding to n = 4, 5 and s = ±1, along with ε F / ω c , using the following InSb QW parameters from Refs. 15 and 39: m * = 0.019m o , g * = −34 and electron density n 2D = 3.3 × 10 −3 nm −2 .For these parameters, the ordering of the energies obeys 3 shows how successive levels cross the Fermi energy as a function of the magnetic field.This, in turn, will reflect on the oscillations of the resistivity once for ε F ≈ ε n,s , an increase on the conductivity will happen due to the resonance condition between the corresponding LL and the Fermi energy.
From the energy expressions above [Eq.( 27)], we can obtain the F -functions through Eq. ( 11), namely, yielding the fast and slow components [Eq.( 13)] At ε = ε F these can be expressed (to a very good approximation) as m0 , where we assume that The corresponding resistivity can now be determined through Eq. ( 22) and reads where f SdH = hn 2D 2e and we have assumed a Lorentzian form for the L Γ broadening.For small magnetic fields, both effective mass and g-factor nominal values do not depend on the magnetic field 60 .As a result, the 1/Bdependence of the resistivity in a 2DEG with only Zeeman coupling, displays oscillations with frequencies multiple of f SdH , and absence of beating.This can be seen from Fig. 4, where we plot δρ xx (B) vs 1/B for the harmonics l = 1, 2, 3 and clearly see oscillations with the respective frequencies f SdH ,2f SdH , and 3f SdH .The solid (dotted) curves correspond to g * = −34 and m * = 0.019m o (g * = 0 and m * = 0.019m o ) 15,39 .Note that the higher harmonics have smaller resistivity amplitudes.This occurs due to the Dingle factor ∝ e −l/B , which suppresses the higher harmonic components.We should stress that the effects of the Zeeman coupling within the plot of δρ xx (B) are not immediately obvious.For instance, it can be seen that for g * = 0 and g * = 0, the corresponding δρ l=1 xx (B) (blue curves depicting the first harmonic) only differ from themselves by the amplitude of the oscillation.For ∆ = −0.323,cos(2π ∆/2) is smaller than one, thus yielding a reduction of the total amplitude for g * = 0 as compared to g = 0.As a consequence, the presence of Zeeman is not readily evident from the oscillations of δρ l=1 xx (B).Conversely, the Zeeman is only manifested within the resistivity when one considers many harmonics, as we discuss below.
The definition of DOS in Eq. ( 21) gives broadened Landau levels separated by ω c , which are in turn spin split by the Zeeman term ∆ [See Eq. ( 27) and Fig. 3].This spin splitting can only be seen in the resistivity [Eq.(30)] when the contributions from the first and second harmonics, cos(2πf SdH /B − π) cos(2π ∆/2) and cos(4πf SdH /B − 2π) cos(4π ∆/2), respectively, have opposite signs.For the parameters of Fig. 3 ∆ = −0.323 the Zeeman term significantly affects the maximum of the resistivity.This can be seen in Fig. 4, where the resistivities associated to harmonics l = 1 and l = 2 (blue and cyan solid curves, respectively), interfere in a destructive way, producing the double-peak feature in the total resistivity (purple solid lines), characteristic of the incipient spin splitting in such data.We emphasize, however, that this feature can be absent depending on the broadening of the energy levels (due to the overlap of the spin-split levels).This is the reason why the double-peak feature is not seen on the other maximum peaks.
Although the g * -factor term does not depend explicitly on magnetic field, it can manifest itself in the magnetooscillations.More specifically, Zeeman-only effects can have a pronounced effect on the magneto-oscillations, controlling the amplitude and sign of how subsequent harmonics are added, either constructively or destructively, before being damped by the quantum life time.Furthermore, it is important to say that the Zeeman can give rise to interesting features and affect drastically the understanding of the magneto-oscillations.For instance, if one could engineer a material 61 such that ∆ = g * 2 m * mo = 0.5 + m with m ∈ Z, then the main weight of the resistivity would be due to the second harmonic with SdH frequently 2f SdH as cos(lπ ∆) = 0 for l = 1.

B. Landau Levels with Zeeman and Rashba interactions
We now analyze the case where we have the presence of both Zeeman and Rashba terms, i.e., ∆ = 0, α = 0 and no Dresselhaus coupling β = 0 in Eq. ( 9).In the spin basis {|↑ , |↓ }, the corresponding Hamiltonian assumes the following matrix form Interestingly, the operator As a consequence, a linear combination of |n, ↑ and |n + 1, ↓ is also an eigenstate of our Hamiltonian Eq. ( 31).This motivates us to rewrite the total Hamiltonian as a direct sum of 2 × 2 block Hamiltonians in the basis {|n, ↑ , |n + 1, ↓ }( H|n,↑ ;|n+1,↓ ), in addition to the non-degenerate decoupled Hamiltonian ( H|0,↓ ), namely with H|0,↓ = ε 0,↓ and The diagonalization of the Hamiltonian Eq. ( 33) yields energies with s = ± and n ∈ N 0 , which already incorporates the energy of the decoupled state . These LLs are plotted in Fig. 5 as a function of 1/B for parameters α = 10 meV nm, m * = 0.019m o and g * = −34 15,39 .Due to the spin-orbit coupling, the energy levels ωc are no longer equidistant, and their separation changes as function of 1/B.On this scale, the energy dispersion appears linear in 1/B.In fact, for ∆ < 0 ( ∆ > 0) the spin-splitting is enhanced (suppresses) relative to the case with α = 0 (See Fig. 3).This can be seen through the expansion of the term (1− ∆) 2 within the square root of Eq. ( 34), yielding −2 ∆, which enhances the Zeeman splitting in the presence of Rashba SO coupling. 39cordingly, for this case we obtain Differently from the results in the previous section, here both F ± functions depend on the magnetic field.As a consequence, we will have more complex oscillations in ρ xx (B) as compared to the case without Rashba coupling [Fig.4].In Fig. 6, we plot the total differential magnetoresistivity δρ xx (B), and the independent contributions from harmonics l = 1, 2 and l = 3.Here we use α = 10 meV nm, m * = 0.019m o , g * = −34 and n 2D = 3.3 × 10 −3 nm −215,39 .Similarly to the case with α = 0 (dashed line in Fig. 6), here we also see oscillations for the l = 1, 2, 3 harmonics with respective frequencies f SdH , 2f SdH and 3f SdH .However, for l = 1 we observe beating, which can be expected as both F − (ε, B) and F + (ε, B) frequencies now depend on 1/B.More specifi- cally, this beating appears here because in the magnetic range considered we have 2πlF − (B) = π 2 , which leads to a node in δρ xx as δρ xx ∝ cos [2πlF − (B)].Note that this only occurs for l = 1, since for higher harmonics this condition is not satisfied.Due to the larger amplitude of the harmonic l = 1, this beating is also seen in the total magneto-resistivity.

C. Landau Levels with Zeeman and Dresselhaus interaction
In the case of Zeeman with pure Dresselhaus, i.e., ∆ = 0, α = 0 and β = 0, the Hamiltonian Eq. ( 9) in the spin basis is given by Differently from the case of pure Rashba, here the operator As a consequence, a linear combination of |n, ↓ and |n + 1, ↑ is also an eigenstate of our Hamiltonian.Therefore, differently from the previous case here the Hamiltonian reads,

H|n,↓
The diagonalization of the Hamiltonian Eq. ( 39) yields energies with s = ± and n ∈ N 0 , which already incorporates the energy of the decoupled state |0, ↑ , ε 0, Here, it is important to notice the opposite sign of s with respect to Eq. ( 33).This happens because the pure Dresselhaus Hamiltonian Eq. ( 39) has opposite basis ordering of the spin states as compared to the pure Rashba Hamiltonian Eq. ( 33).Accordingly, the F ± functions change slightly and read 1.0 1.5 2.0 2.5 3.0 0.5 Due to the: i) similarity of the Dresselhaus expression Eqs. ( 40), ( 41) and ( 42) to the ones arising from the pure Rashba case, Eqs. ( 34), ( 35) and (36); ii) cosine dependence of the F ± functions within the resistivity Eq. ( 22); all the results and equations in the last section also holds here by making α B → β B , ∆ → − ∆ and s → −s.This can also be seen on the level of the Hamiltonian in Eq. (2) where applying the unitary transformation W = e i π 2 σx e i π 4 σz results in which is the expected result.This mapping from (α, ∆) to (β, − ∆) has visible consequences on the energy levels.In Fig. 7 we plot the corresponding LLs [Eq.( 40)] as a function of 1/B for parameters β = 10 meV nm, m * = 0.019m o and g * = −34 15,39 .Due to the spin-orbit coupling, the energy levels ωc are no longer equidistant, and their separation changes as function of 1/B.However, differently from the pure Rashba case, now the Dresselhaus competes with the Zeeman coupling, even leading to LL-dependent crossings.This can be seen through the expansion of (1 + ∆) 2 within the square root [Eq.(40)], which will give rise to 2 ∆ < 0, thus suppressing the spin splitting in the presence of Dresselhaus SO coupling.
In Fig. 8 we plot the total differential magnetoresistivity δρ xx (B), and the individual contributions from the harmonics l = 1, 2 and l = 3.We use β = 10 meV nm, m * = 0.019m o , g * = −34 and n 2D = 3.3 × 10 −3 nm −215,39 .First, similarly to the previous cases, here we can also clearly see oscillations with frequencies f SdH , 2f SdH , 3f SdH .Differently from the pre-vious case with α = 10 meV nm and β = 0, now we see no beating for the l = 1 harmonic but find beating for l = 2.This happens as 2πlF − (B) = π 2 -the condition to observe beating -is only satisfied for l = 2.Even though the beating appears within the second harmonic, it is not manifested in the total differential magnetoresistivity δρ xx (B) for our choice of parameters.This is due to the smaller oscillation amplitude of l = 2 with respect to l = 1.

D. Beatings in the SdH oscillations with nonzero
Zeeman and in the presence of either Rashba or Dresselhaus: a unified description In this section we will discuss more thoroughly the conditions for the appearance of beatings.The two functions F + and F − , Eq. ( 13), determine the fast and slow component, respectively, of the SdH oscillations.To highlight this point and its connection to the power spectrum in Eq. ( 26), we start by rewriting Eqs. ( 35)-( 36), and Eqs.(42)-( 42) in a unified way where we have introduced the magneto-oscillation frequencies where the R (D) index refers to either pure Rasbha (Dresselhaus) case, with k R = mα f R(D) /B 1, the beating frequency takes the stan- , in which case ∆ becomes irrelevant for the magnitude of the beating frequency 36 .
The frequency f SdH R(D) [Eq.(46)] is the main SdH frequency of the magneto-resistance oscillations, usually extracted from experiments to infer the 2D electronic density n 2D .On the other hand, the frequency f R(D) [Eq.( 47)] is the one allowing for possible beatings in the magneto-oscillation.As previously discussed in the last two sections, the presence of beating happens when 2πlF − (B) = π 2 is satisfied, which depends on the value of both f R(D) and ∆.
The presence or absence of beatings can also be visualized through the power spectrum defined by Eq. (26).From interference of waves, we know that the presence of beatings correspond to sum of cosines waves with slightly different frequencies.Accordingly, the power spectrum for this case would show two peaks located at slightly different frequencies.In Fig. 9 we plot I(f ) for m * = 0.019m o and n 2D = 3.3 × 10 −3 nm −2 , using different spin-orbit parameters and g-factor values.For all different sets of parameters, we always have the presence of two main peaks located at both 1/B ≈ 6.8 T −1 and 1/B ≈ 13.6 T −1 .These correspond to the main SdH frequencies for the first and second harmonics, and 2f SdH R(D) , respectively.In the absence of both Rashba, Dresselhaus and g-factor (dashed yellow curve), we observe no beating in the δρ xx (Fig. 4).
On the other hand, for the case of pure Rashba α = 10 meV nm with g = −34 (solid red curve), the presence of the beating in Fig. 6 is made clear by the splitting of the peak of the power spectrum around f = f SdH R in Fig. 9. Interestingly, for α = 10 meV nm with g = 0 (dashed red curve), the splitting of the peak is not seen anymore, thus highlighting the important role of the Zeeman on the visualization of beatings.For the pure Dresselhaus case with β = 10 meV nm and g = −34 (solid blue line), we do not see a peak splitting at the f = f SdH D but rather at f = 2f SdH D , which is consistent with the presence of the beating seen on the second harmonic in Fig. 8. Similarly to the pure Rashba case, for β = 10 meV nm with g = 0 (dashed blue line), the splitting of the peak is not seen anymore, corroborating again the role of the Zeeman term on the presence of beatings.
The apparent "asymmetry" in having peak-splitting for Rashba spin-orbit coupling but not for Dresselhaus (even when they have same SO strength) can be understood from the behavior of the lF − -function vs 1/B, shown in Fig. 10.As already discussed previously in Secs.V B and V C, the condition for beating happens when cos(2lπF − ) = 0 or equivalently, lF − = ±1/4 (±1/4 plotted as gray lines).In the case of Rashba (green lines) one has (1 − ∆) > 1, and the condition for a beating node, cos(2lπF − ) = 0, is reached in the interval of 1/B for l = 1 (solid purple) (gray circles).In the Dresselhaus case, (1 + ∆) < 1, such that lF − for l = 1 only crosses −1/4 for large values of 1/B, where the amplitude of the SdH has already been suppressed.Conversely, lF − crosses 1/4 for l = 2 at smaller values of 1/B, thus guaranteeing the presence of a beating within the magnetic field range, as shown in Fig. 8.

E. Landau Levels with simultaneous Zeeman, Rashba and Dresselhaus interactions: Analytical results
As mentioned earlier, to the best of our knowledge, there are no general exact analytical results for the energies and SdH oscillations corresponding to the case with simultaneous and arbitrary Zeeman, Rashba and Dresselhaus couplings.Therefore, in this section we will outline how to derive an effective approximate solution that can be used to shed light on magnetotransport results for materials, e.g.GaAs or InAs, in which all the three couplings are present.For convenience, we define the sum and difference of the spin-orbit couplings [see definitions of α B and β B following Eq.( 9)] which allows us to rewrite Eq. ( 9) as Note that both the pure Rashba and pure Dresselhaus cases are recovered from the equation above for γ = δ and γ = −δ, respectively.Next, we define the Hamiltonian for γ = δ and γ = −δ which describes the pure Rashba (+) and pure Dresselhaus (−) cases in the presence of the Zeeman coupling.As we already discussed in the previous sec-tions, by defining the operator N ± = a † a ± 1 2 σ z , we obtain [ H± , N ± ] = 0, so the eigenstates of H± are also eigenstates of N ± .The eigenstates of N + (N − ) are then constructed from the pair {|n, ↑ , |n + 1, ↓ } ({|n, ↓ , |n + 1, ↑ }).The above statement is true except for the decoupled eigenstates |0, ↑ (|0, ↓ ) with corresponding eigenenergy ω c (1 − ∆)/2 [ ω c (1 + ∆)/2].The diagonalization of each two state subspace results in with s = + (−) and n ∈ N 0 .Note that this form is valid for both pure Rashba (δ = γ > 0) and Dresselhaus (δ = −γ < 0), Eqs. ( 34) and ( 40), respectively, thus also including the corresponding decoupled state with the lowest eigenvalues of N ± .Note that to recover the pure Zeeman case with no Rashba and Dresselhaus, we should take δ → 0 with δ/|δ| → 1.
When both Rashba and Dresselhaus are present, we can use second order perturbation theory with respect to δ, γ 1 (See Appendix E), to obtain the approximate eigenvalues of the Hamiltonian in Eq. ( 50), namely where the quantities Λ and Ω are defined as where we have introduced ε R / ω c = α 2 B and ε D / ω c = β 2 B .Our goal now is to rewrite Eq. ( 53) in a form that recovers the already obtained exact results for pure Rashba and pure Dresselhaus cases.First, we write Λ = Λ |Λ| |Λ| since Λ changes sign depending on the relative strengths of α and β, similarly to the sign of δ that enters into Eq.(52).By adding and subtracting a term s 2 Λ |Λ| in Eq. ( 53) and after some straightforward algebra we obtain In the case of pure Rashba we have Λ = Ω = δ 2 1− ∆ > 0 while for pure Dresselhaus Λ = −Ω = − δ 2 1+ ∆ < 0; these neatly reduce to the exact results when using second order Taylor expansion of Eq. ( 52).Note that Eq. ( 56) also reproduces the exact result for when α = β and g * = 0 62 , represented here by Λ → 0 with Λ/|Λ| → 1, ∆ = 0, and Ω = 2ε D/R / ω c .The mathematical expression of Eqs. ( 34) and ( 40) motivate us to rewrite Eq. ( 56) as where we have used 1 + x 2 ≈ √ 1 + x 63 .It is important to note that although |Λ| 1, Λ enters the square root multiplied by n, the Landau level index.This means that for high enough n, the product |Λ|n is not necessarily a small quantity.Accordingly, although the equation above becomes exact for either pure Rashba or Dresselhaus case, for α, β = 0 Eq. ( 57) is only valid when |Λ|n 1, besides α B , β B , δ, γ 1 already assumed in Appendix E to obtain Eq. ( 53).
We reiterate that Eq. ( 57) satisfies the exact results for (i) the Zeeman-only case [Eq.( 27)], (ii) the pure Rashba plus nonzero g * [Eq.(34)] and (iii) the pure Dresselhaus plus nonzero g * [Eq.(40)].The case α = β with g * = 0, for which there is also an exact solution 62 , is satisfied to leading order using 1.That is, as mentioned in the previous paragraph, the approximate solution given by Eq. ( 56) reproduces the exact solution for α = β with g * = 0 62 .
As in the case of pure Zeeman, Rashba or Dresselhaus, we can now calculate the F -function from Eq. ( 57).The corresponding results are presented in Appendix F, and by neglecting SO contributions higher or equal than second order in the spin-orbit parameters Λ and Ω (or fourth order in γ and δ), we obtain It is easy to see that these equations recover all the previous results: pure Zeeman [Eq.( 29)], Zeeman with pure Rashba [Eqs.(35) and ( 36)], and Zeeman with pure Dresselhaus [Eqs.(41) and (42)].Additionally, in the case of Λ ≈ 0, F − ≈ − ∆/2 , which reduces to the pure Zeeman case.Accordingly, here F − becomes independent of B (for B 1 T), and therefore, we expect the absence of beatings in the magneto-resistivity, previously seen for both pure Rashba and pure Dresselhaus cases.
F. Generalized SdH magneto-resistivity for arbitrary α, β and g * : new prediction for the absence of beatings.
Using the Eqs.( 58) and ( 59) in Eq. ( 24), we can derive the magnetoresistivity δρ xx (B) for the case with arbitrary Rashba and Dresselhaus couplings and simultaneous nonzero Zeeman field, From Eq. ( 60), we can derive the condition for the absence of beatings for any l by finding the condition for the second cosine being independent of 1/B ; this implies |Λ| = 0, which leads to thus yielding Eq. ( 1) presented in the introduction.For ∆ 1, the above condition is reduced to α ≈ β, corresponding to the situation where the total SO k-dependent effective field becomes unidirectional [5][6][7] .
Note that the above condition does not correspond to any fundamental symmetry, since there is no new conserved quantity in our Hamiltonian with both non-zero Zeeman (g * = 0) and Rashba-Dresselhaus couplings.We reiterate that Eq. ( 61) is entirely distinct from the persistent-spin-helix condition α = β.As shown in Fig. 1(d), the case α = β and g * = 0 does not show peak splitting in the first harmonic but ehxibits beating (or peak splitting) in the second harmonic.Only when g * = 0 (no Zeeman) and α = β there are peak splittings absent altogether 44,62 .

G. Beatings for both α and β non-zero
In the previous sections, we studied the effect of the Zeeman interaction on the frequency splitting of the power spectrum peaks, which represents the beatings in the SdH oscillations.Here we study the interplay of both the Dresselhaus and Rashba interactions on the beatings of the SdH oscillations.
Similarly to what we did leading up to Eq. ( 47), we can obtain the effective beating frequency from the F −function in Eq. ( 59) which results in where the effective SO momentum is We start with the pure Rashba case plus Zeeman, α = 7.0 meV nm and g * = −34.The corresponding power spectrum yields the bottom curve in Fig. 11, similar to the one plotted in Fig. 9.This curve shows two main peaks representing the first two harmonics, and the presence of a split main peak.We assume a Lorentzian broadening τ −1 q = 1.75 meV.When the Dresselhaus coupling β increases, we see the splitting of the main peak reduces until it vanishes for β = 5.0 meV nm (the frequency splitting from Eq. ( 62) is indicated by the gray circles).The absence of beating is indeed expected as predicted by the  61).For larger β, we see that the splitting of the main peak remains neglible.However, in the second harmonic a clear splitting opens up.The condition for having no peak splitting at any harmonics is indeed the condition in Eq. ( 61), where the effects of the SO couplings basically disappear [there are still small SO terms ε R , ε D in Eq. ( 60)].The power spectrum using full numerical calculations are also shown [black dashed], and for this parameter regime the analytical and numerical results agree well.
A similar analysis can be done for the case of pure Dresselhaus with Zeeman, β = 5.0 meV nm and g * = −34, see Fig. 12.Here the splitting is not observed in the main peak, but rather in the second harmonic.As α is increased from 0.0 to 9.0 meV nm, the splitting in the second harmonic decreases, and vanishes at α = 7.0 meV nm, which again corresponds to the condition in Eq. (61).Despite the good accuracy of the approximate analytical solution for α 8.0 meV nm, it starts to deviate from the exact one (full numerics) for higher values of β.This happens because for these values, the combined effect Rashba and Dresselhaus is more pronounced, producing an anti-crossing between different energy levels (see blue curves Fig. 14, discussed further below).While the approximate energies obtained here are always monotonic with respect to 1/B, around the anti-crossing the numerical ones are not.Accordingly, our F -function calculation will not be able to fully describe the SdH oscillations and frequencies around the anti-crossing regions, specifically the approximate solution misses a central peak that starts developing, which will be discussed in the next section.In terms of the Ffunction, the occurance of level anticrossings corresponds to |F − | ≈ 1/2.Since the power spectrum is obtained by integrating δρ xx over a range of 1/B, there is no simple condition determining the validity of the approximate solution.However, looking at the λ B term in Eq. ( 59) the condition yields a useful estimate for the 1/B values where the Dingle factor has not suppressed δρ xx .Equation ( 64) generalizes a similar condition derived in Ref. 57.It is also interesting to note that the analytical result is more accurate for higher harmonics, as the Dingle-factor helps diminishing the amplitude of the anti-crossing at higherfields (see Fig. 14).

VI. LANDAU LEVELS WITH ZEEMAN, RASHBA AND DRESSELHAUS INTERACTIONS:
NUMERICAL RESULTS In the previous section, we have derived an approximate analytical result for the magnetoresistance oscillations in the presence of both Rashba, Dresselhaus and Zeeman interactions.The assumptions and approximations underlying the derivation involved the relatively small SO coupling and the low number of occupied Landau levels.These are satisfied in the low electron density InSb-based 2DEGs of Refs.39 and 64.For higher electron density systems (but still with just a singly-occupied subband at B = 0), such as the InAs/GaSb wells in Ref. 41, a numerical approach is needed.Below we out-line the numerical procedure.The numerical approach also allows us to account for the full form of cubic Dresselhaus term, see Sec.VI A.
For the case of either pure Rashba or Dresselhaus with Zeeman, the absence of anti-crossing in the LL spectrum allow us to obtain exact analytical results for the problem.As we explain below, this does not hold in the presence of both Rashba and Dresselhaus with the Hamiltonian (in the spin basis) Eq. ( 9) Therefore, here we calculate the magnetotransport numerically via the diagonalization of the Hamiltonian above.The F -function method used for the analytical cases can be extended to allow for numerical methods for calculating the energy spectrum, see App.D. As opposed to both the pure Rashba and pure Dresselhaus cases, N ± do not commute with the Hamiltonian above, and therefore, the diagonal basis cannot be described by any linear combination of the previous degenerate eigenstates of N ± .However, there is still a unitary operator, P = exp iπ N ± − 1 2 that commutes with this Hamiltonian, called the parity operator 51,52 , which is discussed in detail in App.B. The corresponding unitary transformation gives PaP † = −a, Pa † P † = −a † and Pσ ± P † = −σ ± , which clearly makes the Hamiltonian Eq. ( 9) invariant due to presence of only a † a, a † σ ± and aσ ± terms.The eigenvalues of P, ±1, help analyze the energy spectrum behavior.
In terms of the spectrum, in the presence of only Rashba SO coupling, we observe multiple crossing between the Rashba eigenstates {|n, − , |n, + } for different n ∈ N 0 , with energy given by Eq. ( 34), obtained through the diagonalization of the Rashba blocks (red boxes within the Hamiltonian matrix in Fig. 13).This is shown by the red solid lines in Fig. 14(a) for α = 7.5 meV nm.In the presence of Dresselhaus SO coupling, the states |n, − and |n + ∆n, + with ∆n ∈ N odd belong to the same parity subspace and adding a Dresselhaus contribution will yield anti-crossing, which open up gaps in the spectrum (blue curves).Conversely, the decoupling between the different parity sets, i.e., |n, − and |n+∆n, + with ∆n ∈ N even , implies multiple crossing between their corresponding energy states.These features are shown by the blue curve in Figs.14(a) and (c), where we have used β = 3.0 meV nm.Other parameters are m * = 0.04, g * = −12 and n 2D = 17.6 × 10 −3 nm −2 .These parameters are for InAs/GaSb-based (double) quantum wells 41 in the electron regime.This regime, as emphasized in Ref. 41, corresponds to the configuration in which the GaSb well is depleted and the system is effectively a single InAs-based asymmetric quantum well with electrons only.Furthermore, we also observe that the effect of the Dresselhaus term is to simply shift the crossing point to a different magnetic field and energy (the crossing-point energy remains constant to lowest order in β but does in general shift for higher values of β).
The contrasting behavior of crossings vs. anti-crossings has direct consequences on the F -function, which will be analyzed in the next paragraphs.First we consider the crossing between states |n, − and |n + ∆n, + , with even ∆n (corresponding to states belonging to different parity subspaces).The F -function are and where we have explicitly added their dependence on α and β.This results in an F -function difference [see Eq.
Note that since the SdH oscillation is dependent on F ± in the form of cos(2πF ± ), we can re-define F − to lie within an unit interval, e.g., Next, we look at the crossing between states belonging to the same parity subspace, i.e., |n, − and |n + ∆n, + for odd ∆n.We recall that this crossing only exists for the pure Rashba case, shown in both Figs.14(a) and (c).Here the relations in Eqs. ( 66) and ( 67) still hold, the only difference being the value of ∆n, which results in Adding a non-zero Dresselhaus contribution will couple these states and lead to an anti-crossing, shown in Figs.14(a) and (c).The anti-crossing result in non halfinteger values of F ± in Eqs. ( 66) and (67) and will lead to a rounding of the sawtooth pattern as seen in Fig. ( The conditions in Eqs. ( 68) and ( 69) lead to values of cos(2πF − ) = 1 [filled circle and rectangle in Fig. 14b)] and cos(2πF − ) = −1 [open cirlce circle in 14b)], respectively, in the case of either pure Rashba or pure Dresselhaus.However, when both Rashba and Dresselhaus are present only the former condition cos(2πF − ) = 1 holds (crossing of states with opposite parity) but the latter condition changes such that cos(2πF − ) > −1 due to anticrossings of states with same parity eigenvalue [open rectangle in Fig. 14b)].This, in turn, affects the shape of the magneto-oscillations leading to an asymmetry in the maximum and minimum values of cos(2πF − ).In Fig. 15 this asymmetry is visible in the magnetoosillations.Here we assume Gaussian broadening with B q = 0.50 T which forms an envelope (black dashed curve).The lowest curve is the pure Rashba (α, β) = (7.5, 0.0) meV nm and there all maximas intersect the envelope [black circles].
The curve for (α, β) = (7.5, 3.0) meV nm shows that only some maxima intersect the envelope, the other maximas correspond to cos(2πF − ) > −1 do not (black circle).This is a direct consequence of the anti-crossing in the spectrum in Fig. 14.The curves for (α, β) = (5.5, 3.0) meV nm and (4.5, 3.0) meV nm show how the anti-crossing becomes larger, eventually leading to an absence of beatings.This can also be seen in the frequency spectrum shown in Fig. 16, for the f ≈ f SdH peak.The lowest curve (blue) corresponds to (α, β) = (7.5, 3.0) meV nm where the spectrum shows well separated peaks.How- ever, as the strength of the Rashba coupling is decreased all the way down to α = 0.5 meV nm for a fixed value of β = 3.0 meV nm a central peak develops and for α between 4.5 and 1.5 meV nm, the two split peaks are barely visible.

A. Extracting α and β from SdH data
The magneto-oscillations can be thought of as a fingerprint of the sample parameters, including Fermi energy ε F , effective mass m * , g * , and α and β.To better capture the influence of the spin-orbit couplings for higher electron density, the full form of the Dresselhaus interaction will be used.For non-zero magnetic fields, this corresponds to having Dresselhaus SO term in Eq. ( 9) replaced with where β 1 = γ k 2 z , γ is material-dependent parameter describing the SO interaction due to bulk inversion asymmetry, and k 2 z is the expectation value of the zcomponent of the square of momentum operator (divided by ), see App.C for details of full Dresselhaus coupling.Note that β in Eq. ( 2) is assumed to include the first harmonic of the cubic Dresselhaus 7,49 , which makes it linearly dependent on the electron density.For instance, if the potential confining the 2DEG is assumed to be an infinite well of width d QW then k 2 z = π 2 /d 2 QW .To model the magnetoresistance data we start from Eq. ( 22), which features (i) a sum over higher harmonics , (ii) rapid oscillations coming from F + , and (iii) damping due to Landau level broadening LΓ .The analysis introduced in the previous section was based on the study of the properties of cos(2πF − ), which forms an envelope on top of the rapid oscillations.Note that in the case having both Rashba and Dresselhaus coupling the rapid oscillations are still dominated by the normal SdH oscillations, i.e.
where B q = √ 2π m * Γ e and Γ is a constant Landau level broadening.For curve C1 in Fig. 17 the fitting with linear Dresselhaus yields values α = 7.2 meV nm and β = 3.0 meV nm.On the other hand, for fitting to the full model we obtain α = 7.6 meV nm, and γ = 85 meV nm 3 .We see that both fits produce equally good curves fitting the experimental data points, with comparable values for the extracted Rashba SO coupling.This indicates that when the Rashba coupling dominates the cubic Dresselhaus term (a 3 -term in Eq. ( 70)), fitting the data with the addition of the cubic term does not strongly affect the result.The results for curve C5 in Fig. 18 behave similarly, i.e. we find fitted values of the Rashba coefficient, α = 6.7 meV nm for the linear Dresselhaus with However, the story is different for the curve C10 shown in Fig. 19.Here the value of Rashba and Dresselhaus coupling are closer, and then the details of the linear vs. cubic Dresselhaus become relevant.Indeed, the linear Dresselhaus model fitting yields α = 5.5 meV nm and β = 2.5 meV nm while the cubic fit gives α = 4.9 meV nm.More importantly the error in the linear fit is quite high, and the fit [black dashed curve] fails to describe the data points.However, the cubic model gives a good fit , with γ = 80 meV nm 3 .This clearly shows the importance of the cubic contributions in samples with high density, where the Rashba and Dresselhaus contributions are of comparable magnitudes.
The fit results in Fig. 17 z , specific values of α β 1 , and β are obtained from the fit.The fit results for α and β for each curve remain relatively insensitive to k 2 z -variations.Note that as k 2 z varies β 1 changes quite rapidly via the fitted value of γ.This is to be expected since lower values k 2 z , correspond to the electron leaking out the InAs quantum well γ into the GaSb, which has a higher bulk value of γ.For higher values of k 2 z the system becomes more strongly confined in the InAs quantum well and the value of γ should tend to the value corresponding to bulk InAs.
The fact that the values of α and β change only slightly as function of k 2 z , as can be seen in Fig. 20, has important consequences on the fitting proceedure.For this reason a fitting with γ and k 2 z both being independent fitting parameters can not be performed, since if β 1 = γ k 2 z is the dominant contribution to the Dresselhaus couplings then there are multiple (infinite) solutions to the equation γ k 2 z = const.and fitting the data with γ and k 2 z independent will not converge 41 .

VII. SUMMARY
We investigated the SdH magneto-oscillations in the resistivity ρ xx of 2DEGs in the presence of spin-orbit (Rashba-Dresselhaus) and Zeeman couplings.We used a semiclassical approach for the resistivity combined with a Poisson summation formula for the Landau-quantized DOS.Our approach allows for an intuitive separation of the slow and fast quantum oscillations in terms of "F-functions", central quantities in our description, essentially being the inverse functions of the spin-resolved Landau-level structure of the system.We study a variety of exact cases such as the pure Zeemann, pure Dresselhaus, and pure Rashba cases -all of which provide analytical expressions for the magnetoresistivity.
More importantly, from our unified and general formulation we also derive, for the first time, an analytical solution for the case with arbitrary Rashba and Dresselhaus couplings and simultaneous non-zero Zeeman coupling (g * = 0).Interestingly, this allows us to derive a unique new condition for the vanishing of the SO-induced beatings in the SdH signals: α/β = [(1 − ∆)/(1 + ∆)] 1/2 , where ∆ = g * m * /2m 0 (i.e., ratio (Zeeman energy)/ ω c ).This new condition does not correspond to any conserved quantity in our Hamiltonian, unlike the persistent-spinhelix condition α = β which is associated with the conservation of spin along some particular axes.We emphasize that our new condition precludes beatings in all harmonics of the quantum oscillations.
We have applied our analytical formulation to describe low-density data for SdH oscillations showing many harmonics in GaAs-based 2DEGs (see SM in Ref. 46) and found an excellent agreement, Fig. 2. We have also applied our theory to low-density InSb-based 2DEGs 15,39 .In addition, we have also developed a detailed numerical calculation for high-density InAs-based 2DEGs, in which an analytical description is not satisfactory.We find excellent agreement with available data for highdensity InAs-based 2DEGs 41,46 .We have also pointed out an inequivalence between the Rashba-dominated + Zeeman vs the Dresselhaus-dominate + Zeeman cases, with only the former showing beatings.This follows from a distinct interplay between the SO and Zeeman terms in these two regimes.We hope our detailed study and unified general formulation will stimulate furhter experimental investigations aiming at verifying our theoretical predictions.

VIII. ACKNOWLEDGMENTS
down to the lowest transverse level results in where Π ± = Π x ± iΠ y , and Π 2 z = 2 k 2 z .Note that now the Dresselhaus spin-orbit coupling is parametrized by two parameters γ and k 2 z , while for the linear approximation, only the single parameter β = (−γ) k 2 z is required.Using the definition in Eqs. ( 3) and (4) the full Dresselhaus Hamiltonian becomes In the absence of spin-orbit interaction a † a can be replaced by its eigenvalue n, which in turn is related to the ratio of the Fermi energy and ω c (valid for ε In the presence of spin-orbit we can still formally rewrite Eq. (C2) as which reduces to the traditional definition of β for low density samples as considered in Sec.V.
The parity operator P introduced in App.B also commutes with the Hamiltonian in Eq. (C2), since the spinorbit terms involve odd powers of a, a † multiplied by either σ + or σ − , and the sign introduced the unitary transformation gets cancelled.For fixed parameter values, the eigenenergies of the Hamiltonian Eq. ( 9) take discrete values.They are obtained numerically by diagonalizing the Hamiltonian matrix using a large enough set of basis states.Finding the F -function as described in Eq. ( 11) is equivalent to a root finding problem for the function This requires the quantum number n to be a continuous variable.which leads to a minor modification of the Hamiltonian diagonlization procedure.The standard diagonalization proceedure is to construct a 2N L matrix from N L harmonic oscillator eigenstates, in addition to the spin degree of freedom.The Pauli matrices form 2×2 blocks that are coupled by the ladder operators a and a † , leading to block tri-diagonal matrix with 2 × 2 block matrices where l runs from 1 to N L (number of Landau levels used in the calculations).To obtain a continuous version of Eqs.(D2) and (D3) a variable x is added to the index l, resulting in The full block-tridiagonal matrix based on the submatrices in Eqs.(D4) and (D5) will then yield a spectrum ε n+x,s , for x ∈ [−1, 1].To further simply the calculations the basis states can be split into P = ±1 subspaces.Each P-subspace contains ordered states { 0+x , 1+x , . . .}.For each subspace, one chooses the two adjecent eigenenergies determined by the condition n+x < ε F ωc < n+x+1 .Subsequently the value of x is found by solving g s (n+x) = 0.
Appendix E: Perturbation theory and "Bogoliubov-de Gennes Hamiltonian" Here we solve the Hamiltonian Eq. ( 50) through a perturbative approach.As the Hamiltonian due to the spinorbit terms are generally much smaller than the Hamiltonian corresponding to free electron gas, we rewrite Eq. ( 50) as with corresponding unperturbed Hamiltonian and perturbation, H 0 and V, respectively.Using now the Schrieffer-Wolff transformation 66,67 , defined by e S , with the constraint V + [S, H 0 ] = 0, we obtain an effective Hamiltonian given by H ef f = H 0 + 1 2 [S, V] + O V 3 .For our system we find S = S γ + S δ , with Heff with The Hamiltonian Eq. (E3) can be rewritten in the Bogoliubov-de Gennes form as Heff which can be diagonalized by a 2 × 2 Bogoliubov-de Gennes transformation, and reads Heff with the diagonal operators ã and ã † .For most semiconductors, we have Ω, Λ, Γ 1.By neglecting the fourth order or higher spin-orbit terms, i.e., δ i γ j with i + j ≥ 4, we obtain Heff with energies For 1 − 2Λ > 0 we obtain which is Eq. ( 53) in the main text.
Appendix F: Approximations leading to Eqs. ( 58) and ( 59) Starting from Eq. ( 57) one can obtain the the Ffunction by inverting the energy levels to obtain l, for each value of s.The resulting equations are We will further simplify these equations by approximating Eqs.(F1) and (F2) up to second order in the spin-orbit parameters Λ and Ω (or fourth order in γ and δ).Accordingly, we rewrite these equations as where A = A 0 + A 1 + A 2 and B = B 1 , with ) Here, the nominal values of the subindices of A i and B j indicate their order on the spin-orbit terms Λ and Ω. Accordingly, we expand the square roots of Eqs.(F3) and (F4) and keep only terms up to second order in either Λ or Ω, yielding As a consequence, we can finally write which are Eqs.( 58) and (59), respectively.
Appendix G: Temperature dependence of the normalized differential resistivity In this section we derive the general temperature dependence of the normalized differential magnetoresistiv-ity in Eq. ( 24) for the systems studied in this work.(G2) being obviously temperature independent.When the temperature is finite but small, i.e., k B T µ ∼ ε F , we have a temperature dependent δρ xx (B).We now analyze the relevant case for low-density semiconductors, but with ε F ε R , ε D , and high number of populated Landau levels, i.e., ε F / ω c 1.With these conditions, all the different analyzed in this manuscript present F ± -functions constant or linearly dependent on the energy, so we write here, F ± ∝ ε + cte, see, for example Eqs. ( 29), ( 44), ( 45), (58), and (59).Using 2πlF + = 2Λ l + ε+φ l + , 2πlF − = 2Λ l − ε + φ l − , with φ l ± properly defined by comparison with these equations, and assuming an energy-independent Dingle factor (only true for Lorentzian broadening.), we need to calculate integrals of the following form,

FIG. 1 .
FIG. 1. Magnetoresistivity for a) pure Rashba α = 7.0 meV nm and b) pure Dresselhaus β = 7.0 meV nm with m * = 0.019mo and n2D = 3.3 × 10 11 cm −2 from Ref. 15.The curves in a) and b) are not the same due to the large gfactor g * = −34.The insets display the normalized FFT including the 2nd harmonic.The presence of beating nodes in δρxx are clearly visible in a) the fundamental and b) the 2nd harmonic, see Fig. 8.The condition for the absence of beatings (single peak for each harmonic) is α = 7.0 meV nm and β = 5.0 meV nm, shown in c), but not α = β = 7.0 meV nm, the persistant spin helix case, shown in d), clearly exhibiting a beating (here a splitting of 2nd harmonic peak).
13] at the crossing ε = ε c and B = B c Accordingly, integer values of F − are equivalent to F − = 0 and therefore, the vanishing of F − provides the field values where the crossing happens.The curves for F − are plotted in Fig. 14(b) for the same parameters as in Fig. 14(a).It presents a sawtooth pattern because values of |F − | > 1/2 are shifted back to the [−1/2, 1/2] interval.The role of the Dresselhaus coupling for these crossings is evident in Fig. 14(b), where the zeros of F − remain zeros for any value of β, but are simply shifted to new values of magnetic field, open circle moves to open rectangle Fig. 14(b).
Figures 17-19 show the experimental data from Ref. 41 for InAs/GaSb quantum wells in the electron regime) along with our theoretical fits [Eq.(72)].We focus on the experimental curves 1, 5 and 10, of Fig. S4 of Ref. 41 that we label as C1, C5 and C10 in Fig. 17-19.The data was fitted to δρ xx (B) in Eq. (72), where F − was calculated numerically.For the fitting we consider both the Dresselhaus coupling in Eq. (9) [black dashed lines], and also with the full Dresselhaus term in Eq. (70) [solid red lines].The black dots are reference points extracted from the data, which are used in the fitting of LΓ (B) cos(2πF − ).The best fittings were produced by assuming Gaussian broadening, namely.

- 19 .
Figures 17-19 show the experimental data from Ref. 41 for InAs/GaSb quantum wells in the electron regime) along with our theoretical fits [Eq.(72)].We focus on the experimental curves 1, 5 and 10, of Fig. S4 of Ref. 41 that we label as C1, C5 and C10 in Fig. 17-19.The data was fitted to δρ xx (B) in Eq. (72), where F − was calculated numerically.For the fitting we consider both the Dresselhaus coupling in Eq. (9) [black dashed lines], and also with the full Dresselhaus term in Eq. (70) [solid red lines].The black dots are reference points extracted from the data, which are used in the fitting of LΓ (B) cos(2πF − ).The best fittings were produced by assuming Gaussian broadening, namely.

FIG. 17 .
FIG. 17.The black dots are reference points for curve C1, solid black line.The black dashed curve is the linear Dresselhaus result and solid red curve full Dresselhaus result.Parameter values from fitting are shown in the inset.Other parameters 41 are m * = 0.019, g = −12 and n2D = 0.0176 nm −2 .

FIG. 19 .
FIG. 19.Similar to Fig. (17), but for curve C10, solid black line.The black dashed curve represents the linear Dresselhaus result, which fails to fit the data.However, the full cubic Dresselhaus term (solid red curve) results in a good fit.Extracted fitting parameters are shown in the inset.Other parameters 41 are m * = 0.019, g = −12 and n2D = 0.0176 nm −2 .
-19 were done for k 2 z = π 2 /d 2 QW where d QW = 12.5 nm 41 .To fully model the sample a self-consistent Poisson-Schrödinger calculation is required 7,46,65 , which is beyond the scope of this work.We can however use different values of k 2 z , which indirectly emulate self-consistent potential details, i.e. increasing the value of k 2 z suggests a stronger confinement in the InAs quantum well, and decreased value of k 2 z would correspond to wavefunctions being less localized in the InAs quantum well.In Fig. 20 the values of α, β 1 , and β are shown as a function of k 2 z from 0.75 π 2 from the three curves are indicated by different forms: C1: circle, C5: triangle, and C10: square.For each value of k 2
Appendix D: The numerical procedure for finding the F -function γ a † + a σ x + iδ a − a † σ y = H 0 / ω c + V/ ω c .

TABLE I .
Definitions of the Zeeman and SO-related quantities used is this work.