1 Introduction

Electronic structure and transport in graphene: quasi-relativistic Dirac – Hartry – Fock self-consistent field approximation

H. V. Grushevskaya and G. G. Krylov
E-mail: grushevskaja@bsu.by
Physics Department, Belarusian State University, 4 Nezalezhnasti Ave., 220030 Minsk, BELARUS

## 1 Introduction

Graphene and graphene-like materials are considered today as a prominent candidates to be used in new devices with functionality based on quantum effects and(or) spin-dependent phenomena in low-dimensional systems. Technically, the main obstacle to such devices implementation is the lack of methods that provide minimizatioon of distortion of these material unique properties in bulk nanoheterostructures. Theoretical approaches and computer simulation play an important role in systematic search of this kind of nanostructures for nanoelectronic applications.

However, significant part of all theoretical consideration of graphene-like materials as well as modern model representations of charge transport in these systems on pseudo-Dirac massless fermion model, originally based on tight binding approximation and applied to the description of graphite [1, 2, 3, 4], which is a bulk material.

According to this approach [2], -electrons in graphene are massless fermion type quasiparticle excitations moving with the Fermi velocity. The approach has been seriously developed and successfully applied to a number of experimental situation, some related review papers on the topics and further references are in [17, 18]. There are few known key points where (at our knowledge) one could expect the necessity of some generalized consideration. The first one is the cyclotron mass dependence upon the carriers concentration. Due to weakness of the signal, modern experimental techniques can register cyclotron mass of charge carriers which is just a little smaller than 0.02 of a free electron mass [5, 6, 7, 8]. Assessments are absent whether this mechanism of conductivity prevails in the region of very small values of charge carriers concentration.

The another point is the experimentally observable carrier asymmetry in graphene. According to modern theoretical concepts the bands for pseudo-relativistic electrons and holes in graphene must be symmetrical. With this in mind in the paper [9] in the generalized gradient approximation there were modelled the hexagonal Si and Ge, with the same structure as in graphene. But the electrons and holes bands near the Dirac points in the Brillouin zone turned out to be strongly asymmetric ones for both cases. Firstly, a Dirac cone deformation takes place far away from a circular shape, as for the second, the Dirac velocities for the valent and conduction bands are different [9]. Therefore, one can assume the existence of some asymmetry in the behavior of pseudo-relativistic electrons and holes of the graphene as well. Since the value of the asymmetry seems to be very small, for its experimental observation one should use a highly sensitive method, such e.g., as based on the measurement of noises [10]. For graphene, such a method is based on large amplitudes of non-universal fluctuations of charge carriers current in the form of nonmonotonic noise in the crossover region of the scattering [11]. At high charge densities, the contribution to the resistance of clean graphene basically gives the scattering of charge carriers on long-range impurities at ordinary (symplectic) diffusion. Regime of pseudo-diffusion with a charge carriers scattering on short-range impurities is realized in the vicinity of the Dirac points of the Brillouin zone. In the papers [12, 13] measurements of quantum interference noise in a crossover between a pseudodiffusive and symplectic regime and magnetoresistance measurements in graphene p-n junctions have been performed, which established asymmetric behaviour of pseudo-relativistic electrons and holes based on asymmetric form of non-monotonic dependence of noise and magnetoresistance. And the third known key point needed theoretical explanation is the replicas existence. A weakly interacting epitaxial graphene on the surface of Ir(111) has a non-distorted hexagonal symmetry due to the weakness of the interaction with the substrate in temperature range up to a room temperature. Therefore in ARPES (angle-resolved photo-electron spectroscopy) spectra, the perturbation of the band structure is manifested in the form of replica of the inverted Dirac cone and mini-gaps in places of quasi-crossings of replicas and the Dirac cone [14]. Authors of paper [14] propose replicas existence explanation such that replicas are produced only in the areas of convergence of C and Ir atoms, that explains the weak intensity of photoelectron emission replicas and brightness of the main cone. However, besides different intensity, the asymmetry of photoelectron emission spectra also manifests itself in the fact that maxima from the replicas are below then zero maximum () of Dirac cone in the ARPES spectrum at the same angle of incidence of photons and with are equal to maxima of ARPES spectrum of replicas oppositely arranged on the hexagon. And conversely, if the corners of oppositely disposed replicas are in the neighborhood , the Dirac cone in the ARPES spectrum is located lower. The above described is possible if axes of the Dirac cone and its replicas are not parallel. It means that the top of the replica does not correspond to the corners of the hexagonal mini Brillouin zone, centered at the Dirac cone corner for the epitaxial graphene on the surface of Ir(111). It has been also demonstrated that epitaxial graphene on SiC(0001) holds a hexagonal mini Brillouin zone near the Dirac points [15, 16].

And at last, a bit more philosophical but also important comment. The majority of modern software for ab initio band structure simulations uses models being some variant of the Dirac equation or at least take into account the known relativistic corrections to the Schrödinger equation when attacking the problems. The quasi-Dirac massless fermion approach based purely on tight binding non-relativistic Hamiltonian seems to be oversimplified and hardly extendable.

With the goal to investigate the balance of exchange and correlation interactions the ab initio band structure simulations of loosely-packed solids(it was used the developed generalization of the LMTO method) have been performed and demonstrated that the strong exchange leads to appearance of an energy gap in the spectrum whereas strong correlation interactions leads to tightening of this gap [20]. The spin-unpolarized ab initio simulations of partial electron densities of two-dimensional graphite have shown that the material is a semiconductor. Interlayer correlations tighten the energy gap that results in semi-metal behaviour of three-dimensional graphite [20]. This means that in the absence of correlation holes, the correlation interaction in a monoatomic carbon layer (monolayer) is weak in comparison with the exchange. This theoretical prediction for spin-nonpolarized graphene were confirmed experimentally in [15, 16], where it was demonstrated a bandgap in bilayer graphene on SiC(0001) and its diminishing up to vanishing in multilayer graphene.

By the way, in a monolayer graphene on SiC(0001) one observes the dispersion of the Dirac cone apexes [15, 16].

Experimentally manufactured quasi-two-dimensional systems, such as graphene, carbon armchair nanotubes and ribbons as well as some types of zigzag carbon nanotubes manifest metallic properties (see, for example, [17, 19]). In this regard, there are discrepancies between theoretical predictions and experimental data.

Therefore, based on the results of ab initio spin-unpolarized simulations of two-dimensional and three-dimensional graphite [20], we can make the following assumption. Carbon low-dimensional systems having the properties similar to graphite-like materials should possess spin-polarized electronic states with the correlation holes (a magnetic ordering). Enhancing of correlation interaction due to correlation-hole contribution leads to tightening of the energy gap and, as a consequence, the emergence of semi-metal conductivity.

The approach we use has been developed earlier in [22] and applied for graphene-like material in [23].

The goal of this chapter is to represent a Dirac – Hartree – Fock self-consistent field quasirelativistic approximation for quasi-two-dimensional systems and to describe the origin of asymmetry of electron – correlation hole carriers in graphene-like materials.

## 2 Graphene model bands with correlation holes

Single atomic layer of carbon atoms (two-dimensional graphite) is called monolayer graphene. Its hexagonal structure can be represented by two triangular sublattices A, B [21]. The primitive unit cell of the graphene contains two carbon atoms C and C belonging to the sublattices A and B respectively.

The carbon atom has four valent electrons s, p, p, p. Electrons s, p, p are hybridized in the plane of the monolayer, p-electron orbitals form a half filled band of -electron orbitals on a hexagonal lattice.

Let electrons with spin ”down”  (”up”) are placed on the sublattice , and the electrons with spin ”up”  (”down”) – on the sublattice , as shown in Fig. 1.

With such a symmetry of the problem, all the relevant bands of sublattices are half-filled and are formed due to correlation holes. In the representation of secondary quantization and Hartree – Fock self-consistent field approximation, when not accounting for the electron density fluctuation correlation, the hole energy is simply added to the electron energy [30]:

 [H(→r)+^Vsc(kr)−^Σx(kr)]ψm(kr)=(ϵm(0)−n∑j=1^ϵ†Pj)ψm(kr), (1) ^ϵ†=ϵ(k)^I (2)

because the sum of projection operators in in parentheses equals to the identity operator : .

We denote spinor wave functions of the valent electrons of graphene as and . From Fig. 1 it follows that the spinor quantum fields and are transformed into each other under the mirror reflection . Therefore, a quasi-particle excitation in the proposed model of graphene is a pair of an electron and a correlation hole. As an electron-hole pairs at the same time represent themselves their proper antiparticle, the wave functions belong to the space of Majorana bispinors , and upper and lower spin components , are transformed via different representations of the Lorenz group

 Ψ′=(ψ′σ˙ψ′−σ)=(eκ2→σ⋅→nψσeκ2(−→σ)⋅→n˙ψ−σ). (3)

It means that behaves as a component , and – as a component of bispinor (3). Using the expression (3) and properties of these operators:

 ˆχ†−σA(→r)|0,+σ⟩\lx@stackreldef==0, ˆχ†+σB(→r)|0,−σ⟩\lx@stackreldef==0

one gets the following expression for the bispinor wave function of an electron in graphene:

 (4)

where is a vacuum vector which consists of uncorrelated vacuum states with spin “down”  and “up” : .

The density matrix is expressed through the components of bispinor (4) as

 ^ρABrr′=⎛⎜ ⎜⎝ˆχ†−σA(→r)ˆχσA(→r)ˆχ†−σA(→r)ˆχ−σB(→r′)ˆχ†σB(→r′)ˆχσA(→r)ˆχ†σB(→r′)ˆχ−σB(→r′)⎞⎟ ⎟⎠. (5)

## 3 Equation for the density matrix

In description we will consider only valent electrons. We denote by the number of atoms in two sublattices. For valent electrons,the Dirac hamiltonian has the following form:

 HD=∑L=A,BN/2∑i=14∑v=1⎧⎪⎨⎪⎩c→α⋅→pivL+βmec2−N∑k=1Ze2|→rivL−→Rk|+∑L
 (7)

Here is a set of Dirac matrices, is the set of Pauli matrices, indices and enumerate s-, p-, p and p electron orbitals, indices and enumerate sublattices and atoms within them respectively, is the electron radius-vector, is the radius-vector of -th carbon atom without valent electrons (atomic core), is the charge of the atomic core, is the electron charge, is the free electron mass, is the speed of light.

The operator (matrix) of the electron density in the mean field approximation when neglecting correlation interactions between electrons, satisfies the equation [29, 30]

 (ı∂∂t−(^h+Σx+Vsc))∑n^ρnn′;rr′=(−ϵn(0))NvNδrr′δ(t−t′), (8)

where is the kinetic energy operator for a single-particle state, is the self-consistent potential, is the exchange interaction, () is the eigenvalue of a non-excited single-particle state (energy of an electron orbital for an isolated atom), , is the number of valent electrons. At the equation (8) can be considered as the equation for the Green function of the quasiparticle excitations [29, 30]:

 (ı∂∂t−(^h+Σx+Vsc))∑n′^ρnn′;rr′=δ(→r−→r′)δ(t−t′). (9)

Further, the ”electron” will be uses in the sense of a quasiparticle.

A Dirac – Hartree – Fock Hamiltonian for quasiparticle excitations in graphene can be obtained by the procedure of secondary quantization of the Dirac Hamiltonian (LABEL:Dirac-Hamiltonian):

 1NvN(ı∂∂t−^hDFH)∑n^ρABnn′;rr′=(−ϵn(0))δrr′δ(t−t′), (10) HDFH=^hD+Σxrel+Vscrel, (11)

where are relativistic analogs of operators .

The equation (10) can be rewritten for the quasiparticle field as

 [E(p)−(^hD+Σxrel+Vscrel)]χ(→r,t)=0. (12)

Here is the energy of the quasiparticle excitation.

Now, one can write the relativistic equation (12) in an explicit form:

 ⎛⎜ ⎜⎝mec2−∑Nk=1Ze2|→r−→Rk|c→σ⋅→pc→σ⋅→p−mec2−∑Nk=1Ze2|→r−→Rk|⎞⎟ ⎟⎠⎛⎜ ⎜⎝ˆχ†−σA(→r)ˆχ†σB(→r)⎞⎟ ⎟⎠|0,−σ⟩|0,σ⟩ +1NvNVNvN∑i,i′=1∫∫d→rid→ri′ ×⎛⎜ ⎜⎝⟨0,−σi|ˆχ†−σAi(→ri)Vri,rˆχσAi(→ri)|0,−σi⟩⟨0,−σi|ˆχ†−σAi(→ri)Vri,rˆχ−σBi′(→ri)|0,−σi′⟩⟨0,σi′|ˆχ†σBi′(→ri′)Vri′,rˆχσAi(→ri′)|0,σi⟩⟨0,σi′|ˆχ†σBi′(→ri′)Vri′,rˆχ−σBi′(→ri′)|0,σi′⟩⎞⎟ ⎟⎠ ×⎛⎜⎝ˆχ†−σi′A(→r)ˆχ†σiB(→r)⎞⎟⎠|0,−σ⟩|0,σ⟩=E(p)I⎛⎜ ⎜⎝ˆχ†−σA(→r)ˆχ†σB(→r)⎞⎟ ⎟⎠|0,−σ⟩|0,σ⟩, (13) ˆχ†−σA(→r)=1NvNNvN∑i′=1ˆχ†−σi′A(→r),ˆχ†σB(→r)=1NvNNvN∑i=1ˆχ†σiB(→r) (14)

where , is the identity matrix, index () enumerates all valent electrons of graphene. From eqs. (13) and (14) the expressions follow for relativistic self-consistent Coulomb potential

 Vscrel⎛⎜ ⎜⎝ˆχ†−σA(→r)ˆχ†σB(→r)⎞⎟ ⎟⎠|0,−σ⟩|0,σ⟩=((Vscrel)AA00(Vscrel)BB)⎛⎜ ⎜⎝ˆχ†−σA(→r)ˆχ†σB(→r)⎞⎟ ⎟⎠|0,−σ⟩|0,σ⟩, (15) (Vscrel)AA=NvN∑i=1∫d→ri⟨0,−σi|ˆχ†−σAi(→ri)V(→ri−→r)ˆχσAi(→ri)|0,−σi⟩, (16) (Vscrel)BB=NvN∑i′=1∫d→ri′⟨0,σi′|ˆχ†σBi′(→ri′)V(→ri′−→r)ˆχ−σBi′(→ri′)|0,σi′⟩; (17)

and exchange interaction term [22]

 Σxrel⎛⎜⎝ˆχ†−σA(→r)ˆχ†σB(→r)⎞⎟⎠|0,−σ⟩|0,σ⟩=(0(Σxrel)AB(Σxrel)BA0) ×⎛⎜⎝ˆχ†−σA(→r)ˆχ†σB(→r)⎞⎟⎠|0,−σ⟩|0,σ⟩, (18) (Σxrel)ABˆχ†σB(→r)|0,σ⟩ =NvN∑i=1∫d→riˆχ†σiB(→r)|0,σ⟩⟨0,−σi|ˆχ†−σAi(→ri)V(→ri−→r)ˆχ−σB(→ri)|0,−σi′⟩, (19) (Σxrel)BAˆχ†−σA(→r)|0,−σ⟩ =NvN∑i′=1∫d→ri′ˆχ†−σAi′(→r)|0,−σ⟩⟨0,σi′|ˆχ†σBi′(→ri′)V(→ri′−→r)ˆχσA(→ri′)|0,σi⟩. (20)

Substitution of the expressions (15) and (18) into eq.  (13) gives

 (21)

Let us perform a variable change and write down the system (21) in components

 [−N∑k=1Ze2|→r−→Rk|+(Vscrel)AA−E(p)]ˆχ†−σA(→r)|0,−σ⟩ +[c→σ⋅→p+(Σxrel)AB]ˆχ†σB(→r)|0,σ⟩=0, (22) [c→σ⋅→p+(Σxrel)BA]ˆχ†−σA(→r)|0,−σ⟩ +[−2mec2−N∑k=1Ze2|→r−→Rk|+(Vscrel)BB−E(p)]ˆχ†σB(→r)|0,σ⟩=0. (23)

From the last equation of the system (2223) we find the equation for the component

 ˆχ†σB(→r)|0,σ⟩=12mec2⎧⎪ ⎪⎨⎪ ⎪⎩1+⎡⎢ ⎢⎣∑Nk=1Ze2|→r−→Rk|−(Vscrel)BB+E(p)2mec2⎤⎥ ⎥⎦⎫⎪ ⎪⎬⎪ ⎪⎭−1 ×[c→σ⋅→p+(Σxrel)BA]ˆχ†−σA(→r)|0,−σ⟩. (24)

## 4 Quasirelativistic corrections

In quasirelativistic limit it is possible to neglect lower components of bispinor on respect to upper ones, as the components of the bispinor have an order of . So, it is sufficient to find upper components to describe the behavior of the system.

With this in mind we eliminate the small lower components in the equation (22), expressing small components through large ones with the help of (24):

 [−N∑k=1Ze2|→r−→Rk|+(Vscrel)AA−E(p)]ˆχ†−σA(→r)|0,−σ⟩+12mec2[c→σ⋅→p+(Σxrel)AB] ×⎧⎪ ⎪⎨⎪ ⎪⎩1−⎡⎢ ⎢⎣−∑Nk=1Ze2|→r−→Rk|+(Vscrel)BB−E(p)2mec2⎤⎥ ⎥⎦⎫⎪ ⎪⎬⎪ ⎪⎭−1[c→σ⋅→p+(Σxrel)BA]ˆχ†−σA(→r)|0,−σ⟩=0.

Expanding the factor in curly brackets in a power series on a small parameter

we obtain the quasirelativistic Dirac – Hartree – Fock approximation for graphene:

 [−N∑k=1Ze2|→r−→Rk|+(Vscrel)AA−E(p)]ˆχ†−σA(→r)|0,−σ⟩+12mec2[c→σ⋅→p+(Σxrel)AB] ×[c→σ⋅→p+(Σxrel)BA]ˆχ†−σA(→r)|0,−σ⟩ +12mec2[c→σ⋅→p+(Σxrel)AB]⎡⎢ ⎢⎣−∑Nk=1Ze2|→r−→Rk|+(Vscrel)BB−E(p)2mec2⎤⎥ ⎥⎦ ×[c→σ⋅→p+(Σxrel)BA]ˆχ†−σA(→r)|0,−σ⟩=0. (26)

Let us find the non-relativistic limit. With this goal in eq.  (26) we write down

 (Σxrel)ABˆχ†−σA(→r)|0,−σ⟩→(Σxrel)ABˆχ†σB(→r)|0,σ⟩ (27)

and leave only first order terms on :

 [−N∑j=1Ze2|→r−→Rj|+(Vscrel)AA−E(p)]ˆχ†−σA(→r)|0,−σ⟩ +12mec2[c→σ⋅→p c→σ⋅→p+c→σ⋅→p(Σxrel)BA+(Σxrel)AB(Σxrel)BA]ˆχ†−σA(→r)|0,−σ⟩ +12mec2(Σxrel)ABc→σ⋅→p ˆχ†σB(→r)|0,σ⟩=0. (28)

After some elementary algebra, we transform the equation (28) to the form

 [→p22me−N∑k=1Ze2|→r−→Rk|+(Vscrel)AA−E(p)]ˆχ†−σA(→r)|0,−σ⟩ +12[(Σxrel)BAˆχ†−σA(→r)|0,−σ⟩+(Σxrel)AB ˆχ†σB(→r)|0,σ⟩] +12mec2(Σxrel)AB(Σxrel)BAˆχ†−σA(→r)|0,−σ⟩=0. (29)
 (30)

Since at replacements and , first two terms do not change the form of equation, then they give the non-relativistic contributions. Quadratic summand

 (Σxrel)AB(Σxrel)BA (31)

is a quasirelativistic correction, because its form is sensitive to the above mentioned change.

Since in non-relativistic limit the quasirelativistic quadratic correction (31) should be omitted, the substitution of the expressions (15) and (18) into eq. (30) leads to non-relativistic equation

 [→p22me−N∑k=1Ze2|→r−→Rk|+NvN∑i=1∫d→ri⟨0,−σi|ˆχ†−σAi(→ri)VVri,rˆχσAi(→ri)|0,−σi⟩−E(p)] ×ˆχ†−σA(→r)|0,−σ⟩+12[NvN∑i′=1∫d→ri′ˆχ†−σAi′(→r)|0,−σ⟩⟨0,σi′|ˆχ†σBi′(→ri′)Vri′,rˆχσA(→ri′)|0,σi⟩ (32)

Presenting as a difference of -th energy eigenvalue for one-electron non-excited state and the energy eigenvalue for the hole : and taking into account the chain of equalities

 ≡NvN∑i′=1∫d→ri′ˆχ†−σAi′(→r)|0,−σ⟩⟨0,σi′|ˆχ†σBi′(→ri′)V(→ri′−→r)ˆχ−σB(→ri′)|0,−σi′⟩ =NvN∑i′=1∫d→ri′ˆχ†σBi′(→r)|0,−σ⟩⟨0,−σi′|ˆχ†−σAi′(→ri′)V(→ri′−→r)ˆχ−σB(→ri′)|0,−σi′⟩, (33)

one can rewrite the equation (32) as

 [→p22me−N∑k=1Ze2|→r−→Rk| +NvN∑i=1∫d→ri⟨0,−σi|ˆχ†−σAi(→ri)V(→ri−→r)ˆχσAi(→ri)|0,−σi⟩−(ϵ(0)m−ϵ(p))] ×ˆχ†−σA(→r)|0,−σ⟩+NvN∑i=1∫d→ri (34)

As mentioned above, in non-relativistic limit the indices can be omitted and eq. (34) can be written in a final form

 ϵ(0)mˆψ†σm(→rm)|0,σm⟩−⟨0,−σi|ˆϵ† ˆI|0,−σi⟩ˆψ†σm(→rm)|0,σm⟩ =^h(→rm)ˆψσm(→rm)|0,σm⟩ −n∑i=1∫d→riˆψ†σi(→rm)|0,σm⟩V(→ri−→rm)⟨0,−σi|ˆψ†−σi(→ri)ˆψσm(→ri)|0,−σi⟩ +n∑i=1∫d→riˆψ†σm(→rm)|0,σm⟩V(→ri−→rm)⟨0,