Spin-relaxation time in materials with broken inversion symmetry and large spin-orbit coupling

We study the spin-relaxation time in materials where a large spin-orbit coupling (SOC) is present which breaks the spatial inversion symmetry. Such a spin-orbit coupling is realized in zincblende structures and heterostructures with a transversal electric field and the spin relaxation is usually described by the so-called D’yakonov-Perel’ (DP) mechanism. We combine a Monte Carlo method and diagrammatic calculation based approaches in our study; the former tracks the time evolution of electron spins in a quasiparticle dynamics simulation in the presence of the built-in spin-orbit magnetic fields and the latter builds on the spin-diffusion propagator by Burkov and Balents. Remarkably, we find a parameter free quantitative agreement between the two approaches and it also returns the conventional result of the DP mechanism in the appropriate limit. We discuss the full phase space of spin relaxation as a function of SOC strength, its distribution, and the magnitude of the momentum relaxation rate. This allows us to identify two novel spin-relaxation regimes; where spin relaxation is strongly non-exponential and the spin relaxation equals the momentum relaxation. A compelling analogy between the spin-relaxation theory and the NMR motional narrowing is highlighted.

This Supplementary Material is organized as follows: we discuss the details of the O (L/EF) approximation, which we employed in the main text. We further discuss the relationship between S and the dynamic spin-susceptibility and we also give an analytic treatment of the spin-dynamics in the clean limit. We present the equivalence of the two-site motional narrowing NMR problem and the 2D Bychkov-Rashba induced spinrelaxation for spins aligned perpendicular to the plane. The spin-relaxation is discussed when both Dresselhaus and Bychkov-Rashba spin-orbit couplings are present. FIG. 1. Schematics of the O (L/EF) approximation for a one-dimensional quadratic band dispersion. For such a dispersion, the Bychkov-Rashba type SOC shifts the up/down dispersions to the right/left, respectively. At the Fermi surface intersection, the SOC induced splitting is |k| dependent. In contrast, our approximation is equivalent to substituting it by a hypothetical, Zeeman-like split band structure, where the SOC induced splitting is independent of |k|.
It was mentioned in the main text that both the Monte Carlo and the diagrammatic technique in Ref. 1 neglects the effect of the SOC on the Fermi surface. A proper calculation should consider that the bands are split due to the SOC and therefore the corresponding Ω(k) is |k| dependent due to this effect. In contrast, our approximation neglects this effect and Fig. 1 depicts this approximation for a one-dimensional band dispersion. The figure suggests that corrections to our approximation are on the order of O (L/E F ).

II. RELATIONSHIP BETWEEN S AND THE DYNAMIC SPIN-SUSCEPTIBILITY
We discussed in the main text that the time decay of a spin-polarized ensemble (described by s(t)) is calculated with a Monte Carlo approach for both the clean and dirty cases. We found that the real part of its Fourier transform, ReS(ω) can be conveniently displayed in order to demonstrate the spin-relaxation properties. Here, the i, i ′ components of S refer to the i ′ component of a spin ensemble which was originally polarized along i, with i, i ′ = x, y, or z.
We also showed that a relation between ReS(ω) and the dynamic spin-susceptibility, χ(ω) holds: We present the above relationship on a quantitative agreement which we numerically obtained by comparing our Monte Carlo results on S(ω) with analytic calculations of χ(ω) for a particular Hamiltonian which is available in the literature (Ref. 2). The result is shown in Fig. 2 where L R and L D are the strengths of the Bychkov-Rashba and the Dresselhaus type spin-orbit couplings, respectively. The dynamic spin-susceptibility was calculated in the absence of momentum relaxation and when L R ∼ L D ≪ E F . This limit matches the above discussed O (L/E F ) approximation. Eq. 15 in Ref. 2 gives the dynamic spin-susceptibility for the various components as: . demonstrates a good agreement between the two kinds of data which supports the statement in the main text concerning the connection between S and χ.

III. ANALYTIC TREATMENT OF THE SPIN DYNAMICS BY THE TIME EVOLUTION IN THE CLEAN LIMIT
We give an analytic description of the time evolution of the magnetization of an electron ensemble subjected to internal SOC fields to verify further the Monte Carlo method applied in the main text since this description is equivalent to the numerical method in the clean limit. We consider the SOC in the form H 0 = 1 2 Ω(k) · σ, where σ is a vector composed by the Pauli matrices, and Ω(k) is the k-dependent internal SOC field.
The time evolution of the state of an electron under the SOC is determined by the time evolution operator U (t) = exp(−iH 0 t/ ). Supposing that the electron is initially in spin-up state, i.e. its spin is polarized along the z direction, its state ket at time t is obtained as ] and E ± are the eigenkets and eigenenergies of the Hamiltonian H 0 . Finally, the time development of the z component of the electron spin is obtained as where Ω = Ω 2 x + Ω 2 y + Ω 2 z . The quantity S z (t), i.e. the z component of a spin ensemble, calculated by Monte Carlo method in the main text is obtained as within this approach, i.e. by integration over k on the Fermi surface. Arbitrary Ω(k) can be considered in the Hamiltonian H 0 such as the two-dimensional Bychkov-Rashba SOC or the three-dimensional Dresselhaus case discussed in the main text. We note that for a complicated distribution of the SOC fields the k integration in Eq. (7) might not be performed analytically. In such a case, the integration can be performed by choosing random k values on the Fermi surface. When this calculation is performed according to Eq. (6) for different model Hamiltonians such as those given in the main text, we obtained numerically identical results (data not shown) as for the Monte Carlo, and the analytic result for the dynamic spin susceptibility given in Eq. (3) is reproduced as well.
By taking the Bychkov-Rashba SOC, Ω(k) = L kF [k x , k y , 0], Ω(k) becomes k independent as Ω(k) = L/ = Ω, which means a single oscillating component in S z (t) as it is shown by dashed curve in Fig.2 of the main text, and two Dirac-delta peaks at ±Ω in the real part of the Fourier transform S(ω). In three-dimensional cases, there is always a t-independent non-zero term in S z (t) coming from the last term in Eq. (6), which explains the finite S z value in Fig.4 of the main text with Γ = 0. This non-zero and time independent term corresponds to a Dirac-delta function centered on ω = 0 in S(ω).
The origin of this effect is further supported by a geometric consideration which is depicted in Fig. 3. The presence of the Dirac-delta peak for S(ω = 0) is a generic feature and its absence for the two-dimensional electron gas and the Bychkov-Rashba SOC is an exeption. For the latter, when the spins are aligned perpendicular to the 2D plane, all SOC fields are in the plane, i.e. the precession of the spins around the built in Ω(k) results in a zero-averaged net magnetization. However, for a general distribution of the SOC fields and the corresponding Ω(k) vectors, the precession of the spins retains a finite positive s z component as Fig. 3. depicts. A straightforward geometric consideration shows that the s z , i.e. the Dirac delta function strength is given by Ω 2 z /Ω 2 (Ω z and Ω are the z component and t he magnitude of the Ω(k) vector, respectively) for a particular Ω(k) component. Similarly, we obtain that the amplitude of the oscillation goes as 1 − Ω 2 z /Ω 2 , in full agreement with Eq. (6). In Fig. 4., we present the comparison between the ReS(ω) and the histogram of the internal Larmor frequency distribution for the three-dimensional Dresselhaus model. The latter data integrated for the positive frequencies gives the average of the SOC fields which coincides with the frequency value where it is peaked and its weighted integral, i.e the integral of the histogram values multiplied by the related frequency gives unity since all weights are summed up in this way. FIG. 5. Simulated lineshapes for the two-site NMR motional narrowing problem. Note that two peaks are observed for smaller values of Γ whose linewidth increases with increasing Γ. In contrast, a single, motionally narrowed peak is observed for larger Γ values, whose linewidth decreases with increasing Γ.

IV. ANALOGY BETWEEN THE MOTIONAL NARROWING AND THE SPIN-RELAXATION FOR THE 2D BR MODELL
Abragam 3 considered the so-called two-site NMR motional narrowing problem: a nuclei is allowed to jump with the transition rate Γ c = 1/τ c between two sites with different local Larmor frequencies: ±Ω around a central Larmor frequency (defined as zero in this case). The resulting NMR lineshape is shown in Fig. 5 for a fixed Ω = 1 and different values of the jumping frequency, Γ c . The analogy between the spin-relaxation and the motional narrowing is clear: the ±Ω local Larmor frequencies correspond to the built-in Zeeman field distribution of the spin-relaxation problem and the jumping frequency (Γ c ) of the motional narrowing problem corresponds to the Γ momentum relaxation rate (besides a factor 2 which is discussed below). The analogy can be quantified for the simplest case as follows.
Our definition of Γ differs from that of Γ c in Abragam's work in a factor 2 as discussed herein, as Γ corresponds to the momentum relaxation rate and Γ c corresponds to the transition rate between the two states. The rate equations for the two state's populations: We can see that 1 τm = 2 τc , so Γ = 2Γ c . With this change, the real and imaginary values of the above defined ω 1,2 are shown in Fig. 6. The real part of the roots describe the position of the two peaks and the imaginary parts describe the linewidths in agreement with Fig. 6. Remarkably, this figure is identical to the spin-relaxation problem for the 2D Rashba model in the main paper with a straightforward identification of the correspondence of the two parameters.

V. SPIN-RELAXATION WHEN DRESSELHAUS AND BYCHKOV-RASHBA SPIN-ORBIT COUPLINGS ARE PRESENT
Both the Dresselhaus and Bychkov-Rashba spin-orbit couplings are present in III-V sermiconductor based quantum wells. The earlier is due to the bulk inversion symmetry breaking while the latter is due to a longitudinal electric field acting on the quantum well.
When the growth direction is along the [0, 0, 1], both the Dresselhaus and the Bychkov-Rashba Hamiltonian has the usual form which is mentioned in the main text and the resulting the Larmor-frequency distribution reads: When the growth direction is along the [1, 1, 1], we take the z-axis along this direction, therefore the Bychkov-Rashba Hamiltonian is retained and the Dresselhaus one is transformed. The resulting Larmor-frequency distribution reads: In our examples we use L D = L R = L. In Fig. 7., we show the components of ReS(ω) for the various SOC configurations and for different values of Γ. The first column is for the pure Dresselhaus SOC and is an identical result to that given in the main text. Γ values in the rows are 0, 0.02L, and 5L, respectively. The scales are the same for the first two rows of data but are different for the third one (with large Γ).
The important observations are: i) the overall SOC field becomes larger as expected with a singular (step-like) feature (indicated by an arrow in the figure), ii) the spin-relaxation becomes anisotropic when the Bychkov-Rashba term is also present and even an xy term is present for the [0, 0, 1] growth direction, which however vanishes in the D'yakonov-Perel' regime, iii) the TABLE I. The spin-relaxation broadening parameter, Γs in the D'yakonov-Perel' regime for the various SOC combinations in Fig. 7. for Γ = 5L. The definition of Γs is given in the text .

SOC model and direction
Γs in units of L anisotropy remains for the [1, 1, 1] growth direction but the xy term is zero, iv) a single Lorentzian component is observed in the Dyakonov-Perel regime but with varying widths which are summarized in Table I. The broadening parameter, Γ s is the HWFM of the Lorentzian curves in the figures.
These parameters indicate how strong the spin-relaxation in these models is. E.g. the 0.017 value is related to the corresponding spatial average of the S clean function.