Spin-orbit coupling in a graphene bilayer and in graphite

The intrinsic spin-orbit interactions in bilayer graphene and in graphite are studied, using a tight binding model, and an intraatomic LS coupling. The spin-orbit interactions in bilayer graphene and graphite are larger, by about one order of magnitude, than the interactions in single layer graphene, due to the mixing of pi and sigma bands by interlayer hopping. Their value is in the range 0.1 - 1K. The spin-orbit coupling opens a gap in bilayer graphene, and it also gives rise to two edge modes. The spin-orbit couplings are largest, 1-4K, in orthorhombic graphite, which does not have a center of inversion.


INTRODUCTION
The isolation and control of the number of carriers in single and few layer graphene flakes [1,2] has lead to a large research activity exploring all aspects of these materials [3]. Among others, the application of graphene to spintronic devices [4][5][6][7][8][9][10] and to spin qubits [11][12][13] is being intensively studied. The understanding of these devices requires a knowledge of the electronic spin-orbit interaction. In principle, this interaction turns single layer graphene into a topological insulator [14], which shows a bulk gap and edge states at all boundaries. The magnitude of the spin-orbit coupling in single layer graphene has been studied [15][16][17][18]. The calculated couplings are small, typically below 0.1K. The observed spin relaxation [8,19] suggests the existence of stronger mechanisms which lead to the precession of the electron spins, like impurities or lattice deformations [20][21][22].
Bilayer graphene is interesting because, among other properties, a gap can be induced by electrostatic means, leading to new ways for the confinement of electrons [23]. The spin-orbit interactions which exist in single layer graphene modulate the gap of a graphene bilayer [24]. The unit cell of bilayer graphene contains four carbon atoms, and there are more possible spin-orbit couplings than in single layer graphene.
We analyze in the following the intrinsic and extrinsic spin-orbit couplings in bilayer graphene, using a tight binding model, and describing the relativistic effects responsible for the spin-orbit interaction by a L S intraatomic coupling. We use the similarities between the electronic bands of a graphene bilayer and the bands of three dimensional graphite with Bernal stacking to generalize the results to the latter.

THE MODEL
. We describe the electronic bands of a graphene bilayer using a tight binding model, with four orbitals, the 2s and the three 2p orbitals, per carbon atom. We consider hoppings between nearest neighbors in the same plane, and nearest neighbors and next nearest neighbors between adjacent layers, see [25]. The couplings between each pair of atoms is parametrized by four hoppings, V ss , V sp , V ppπ and V ppσ . The model includes also two intraatomic levels, s and p , and the intraatomic spin-orbit coupling The parameters used to describe the π bands of graphite [26,27], γ 0 , γ 1 , γ 2 , γ 3 , γ 4 , γ 5 and ∆, can be derived from this set of parameters. We neglect the difference between different hoppings between atoms which are next nearest neighbors in adjacent layers, which are responsible for the difference between the parameters γ 3 and γ 4 . We also set the difference in onsite energies between the two inequivalent atoms, ∆ to zero. The parameters γ 2 and γ 5 are related to hoppings between next nearest neighbor layers, and they do not play a role in the description of the bilayer. The total number of parameters is 15, although, without loss of generality, we set p = 0. We do not consider hoppings and spin orbit interactions which include d levels, although they can contribute to the total magnitude of the spin-orbit couplings [18,28]. The effects mediated by d orbitals do not change the order of magnitude of the couplings in single layer graphene, and their contribution to interlayer effects should be small. The main contribution to the effective spin-orbit at the Fermi level due to the interlayer coupling is due to the hoppings between p orbitals in next nearest neighbor atoms in different layers. This interaction gives rise to the parameters γ 3 and γ 4 in the parametrization of the bands in graphite. For simplicity, we will neglect couplings between s and p orbitals in neighboring layers. The non zero hoppings used in this work are listed in Table I.
The hamiltonian can be written as a 32 × 32 matrix for each lattice wavevector. We define an effective hamiltonian acting on the π, or p z , orbitals, by projecting out the rest of the orbitals:  [29,30], and the spinorbit coupling from [31]. Superindices 0,1, and 2 correspond to atoms in the same layer, nearest neighbors in different layers, and next nearest neighbors in different layers. We isolate the effect of the spin-orbit coupling by defining: Note that H so π depends on the energy, ω.
We analyze H so π at the K and K points. The two matrices have a total of 16 entries, which can be labeled by specifying the sublattice, layer, spin, and valley. We define operators which modify each of these degrees of freedom using the Pauli matricesσ,μ,ŝ, andτ . The unit cell is described in Fig. 1.
The hamiltonian has inversion and time reversal symmetry, and it is also invariant under rotations by 120 • . These symmetries are defined by the operators: where K is complex conjugation. The possible spin dependent terms which respect these symmetries were listed in [32], in connection with the equivalent problem of three dimensional Bernal graphite (see below). In the notation described above, they can be written as The first term describes the intrinsic spin-orbit coupling in single layer graphene. The other three, which involve the matrices µ i , are specific to bilayer graphene. The term proportional to λ 3 can be viewed as a Rashba coupling with opposite signs in the two layers.
The energy dependence of the four couplings in eq. 5 is shown in Fig. 2. The values of the couplings scale linearly with ∆ so . This dependence can be understood by treating the next nearest neighbor interlayer coupling and the intratomic spin-orbit coupling as a perturbation. The spin-orbit coupling splits the spin up and spin down states of the σ bands in the two layers. The interlayer couplings couple the π band in one layer to the σ band in the other layer. Their value is of order γ 3 . The π states are shifted by: where 0 σ is an average value of a level in the σ band. The model gives for the only intrinsic spin-orbit coupling in single layer graphene the value This coupling depends quadratically on ∆ so , δ π± ∼ ±∆ 2 so / 0 σ [15]. The band dispersion of bilayer graphene at low energies, in the absence of spin-orbit couplings is given by four Dirac cones, because of trigonal warping effects associated with γ 3 [23]. Hence, we must to consider the couplings for wavevectors k slightly away from the K and K points. We have checked that the dependence of the couplings λ i on momentum, in the range where trigonal warping is relevant, is comparable to the changes with energy shown in Fig. 2.
A gap, E g , between the two layers breaks inversion symmetry, and can lead to new couplings. The calculations show no new coupling greater than 10 −6 meV for gaps in the range −0.1eV ≤ E g ≤ 0.1eV. The dependence of the couplings on the value of the gap is shown in Fig. 3. This calculation considers only the effect in the shift of the electrostatic potential between the two layers. The existence also of an electric field will mix the p z and s orbitals within each atom, leading to a Rashba term similar to the one induced in single layer graphene [15,16].
The effect of λ 1 is to open a gap of opposite sign in the two valleys, for each value of s z . The system will become a topological insulator [14,33]. The number of edge states is two, that is, even. The spin Hall conductivity The terms with λ 3 and λ 4 describe spin flip hoppings which involve a site coupled to the other layer by the parameter γ 1 . The amplitude of the wavefunctions at these sites is suppressed at low energies [23]. The shifts induced by λ 3 and λ 4 in the low energy electronic levels will be of order λ 2 3 /γ 1 , λ 2 4 /γ 1 .

Bulk graphite
. The hamiltonian of bulk graphite with Bernal stacking can be reduced to a set of bilayer hamiltonians with interlayer hoppings which depend on the momentum along the direction perpendicular to the layers, k z . We neglect in the following the (small) hoppings which describe hoppings between next nearest neighbor layers, γ 2 and γ 5 , and the energy shift ∆ between atoms in different sublattices. At the K and K points of the three dimensional Brillouin Zone (2k z c = 0, where c is the interlayer distance) the hamiltonian is that of a single bilayer where the value of all interlayer hoppings is doubled. At the H and H points, where 2k z c = π, the hamiltonian reduces to two decoupled layers, and in the intermediate cases the interlayer couplings are multiplied by |2 cos(k z c)|. Carrying out the calculations described in the previous section, k z dependent effective couplings, λ i (k z ), can be defined. These couplings are shown in Fig. 4. The results for bilayer graphene correspond to k z c = 2π/3, 4π/3. The layers are decoupled for k z c = π. In this case, the only coupling is λ 1 , which gives the coupling for a single layer, given in eq. 7.
The significant dispersion as function of momentum parallel to the layers shown in Fig. 4 implies the existence of spin dependent hoppings between layers in different unit cells. This is consistent with the analysis which showed that the spin-orbit coupling in a bilayer has a contribution from interlayer hopping, see eq. 6.
The spin-orbit couplings can be larger in bulk graphite than in a graphene bilayer. The bands in Bernal graphite do not have electron-hole symmetry. The shift in the Fermi energy with respect to the Dirac energy is about E F ≈ 20meV λ 1 , λ 3 [34]. Hence, the spin-orbit coupling is not strong enough to open a gap throughout the entire Fermi surface, and graphite will not become an insulator.
A similar analysis applies to orthorhombic graphite, which is characterized by the stacking sequence ABCABC · · · [35]. The electronic structure of this allotrope at low energies differs markedly from Bernal graphite [36,37], and it can be a model for stacking defects [36][37][38]. If hoppings beyond nearest neighbor layers are neglected, the hamiltonian can be reduced to an effective one layer hamiltonian where all sites are equivalent. The effective hamiltonian which describes the K and K valleys contains eight entries, which can be described using the matrices σ i , s i , and τ i . Orthorhombic graphene is not invariant under inversion, and a Rashba like spin-orbit coupling is allowed. The spin-orbit coupling takes the form: As in the case of Bernal stacking, the couplings have a significant dependence on the momentum perpendicular to the layers, k z , and interlayer hopping terms are induced. For ω = 0, k = 0 and k z = 0, we find: In orthorhombic graphite the Fermi level is away from the K and K points, in the vicinity of a circle defined by | k| = γ 1 /v F [36,37]. The variation of the couplings as function of wavevector is shown in Fig. 5.

CONCLUSIONS
We have studied the intrinsic spin-orbit interactions in a graphene bilayer and in graphite. We assume that the origin of the couplings is the intraatomic L S interaction, and we use a tight binding model which includes the 2s and 2p atomic orbitals.
The intrinsic spin-orbit couplings in a graphene bilayer and in graphite are about one order of magnitude larger than in single layer graphene, due to mixing between the π and σ bands by interlayer hoppings. Still, these couplings are typically of order 0.01 − 0.1meV, that is, 0.1 − 1K.
Bilayer graphene becomes an insulator with an even number of edge states. These states can be mixed by perturbations which do not break time reversal symmetry. These perturbations can only arise from local impurities with strong spin-orbit coupling, as a spin flip process and intervalley scattering are required.
The interplay of spin-orbit coupling and interlayer hopping leads to spin dependent hopping terms. The spinorbit interactions are largest in orthorhombic graphite, which does not have inversion symmetry.