Band structure of graphene, massless Dirac fermions as low-energy quasiparticles, Berry phase, and all that

Graphene as the first truly two-dimensional crystal

The surprising experimental discovery of a two-dimensional (2D) allotrope of carbon, termed graphene, has ushered unforeseen avenues to explore transport and interactions of low-dimensional electron system, build quantum-coherent carbon-based nanoelectronic devices, and probe high-energy physics of "charged neutrinos" in table-top experiments. Graphene represents one-atom–thick layer of carbon atoms tightly packed into a honeycomb crystal lattice. Its discovery is unexpected in the view of Mermin-Wagner theorem which prevents breaking of continuous symmetry in less than three-dimensions, which means that at first sight 2D crystals (of astronomical size though) would be unstable at any finite temperature.

Basic chemistry of graphene

FIG 1: Chemistry of graphene.

Carbon has four valence electrons, of which three are used for the $sp^{2}$ bonds in graphene, graphite, and carbon nanotubes. The origin of covalent bonding between two carbon atoms in the Fig. 1 is tunneling of electrons through the potential barrier between the atoms which removes the degeneracy of the original two-level system, making spatially symmetric wave function (the bonding one) to have lower energy eigenvalue of the corresponding Hamiltonian.

When many carbon atoms are brought into a crystal, each bonding or antibonding orbital will acquire dispersion to become an energy band. The analysis below is concerned with the band structure which originates from $p_{z}$ orbital hosting the fourth valence electron. The bands which correspond to the dispersion of bonding and antibonding molecular orbital (constructed from $p_{z}$ orbitals on two carbon atoms) are called $\pi$ and $\pi ^{*}$ bands.

Lattice structure of graphene

FIG 2: (a) Honeycomb lattice of graphene. (b) Primitive unit cell with primitive translation vectors $\mathbf {a}$ and $\mathbf {b}$ and the corresponding primitive vectors $\mathbf {G} _{1}$ and $\mathbf {G} _{2}$ of the reciprocal lattice. (c) Segment of graphene reciprocal lattice and the first Brillouine zone within it.

Real space honeycomb lattice: The honeycomb lattice of graphene show in Fig. 2(a) is not a Bravais lattice. Instead, it can be viewed as bipartite lattice composed of two interpenetrating triangular sublattices. Thus, such lattice with a basis contains two identical carbon atoms per primitive unit cell.

Reciprocal space honeycomb lattice: The corners of the Brillouin zone are at the points $\mathbf {K} _{i}$, which can be expressed as:

$\mathbf {K} _{1}=(\mathbf {G} _{2}-\mathbf {G} _{1})/3$,
$\mathbf {K} _{2}=(2\mathbf {G} _{2}+\mathbf {G} _{1})/3$,
$\mathbf {K} _{3}=(2\mathbf {G} _{1}+\mathbf {G} _{2})/3$,
etc.

Among six corners, only two are inequivalent. For example, notice that $\mathbf {K} _{3}=\mathbf {K} _{1}+\mathbf {G} _{1}$.

Symmetries: Because the lattice is two-dimensional, all translations commute with reflection in the plane of the lattice, so all electron (or vibrational) eigenstates can be chosen to be either even or odd under this reflection. For this reason, the single-particle electron states are rigorously separated into two classes, called $\sigma$ and $\pi$. The even $\sigma$ states are derived from carbon $s$ and $p_{x}$, $p_{y}$ orbitals (i.e., their hybridized $sp^{2}$ orbitals discussed above), while the odd $\pi$ states are derived from carbon $p_{z}$ orbitals. These latter are cylindrically symmetric in the $xy$ plane, lie near the Fermi level (half-filled) and are the electrically active states of interest in low energy experimental probing of graphene.

Tight-binding Hamiltonian of graphene

A useful picture of electron behavior can be derived by using the phenomenological nearest-neighbor tight-binding model to look at the $\pi$ electrons. The two sublattices are shown below in different colors, with the $\mathrm {A}$ sublattice at vectors $\mathbf {R} =n_{1}\mathbf {a} +n_{2}\mathbf {b}$, and the $\mathrm {B}$ sublattice at vectors $\mathbf {R} +\mathbf {\tau }$ with $\mathbf {\tau } =-(\mathbf {a} +\mathbf {b} )/3$.

FIG 3: Vectors specifying the position of carbon atoms on two triangular Bravais sublattices comprising the honeycomb lattice of graphene.

The tight-binding model is then defined as:

${\hat {H}}=-t\sum _{\mathbf {R} }\left(|\mathbf {R} \rangle \langle \mathbf {R} +\mathbf {\tau } |+|\mathbf {R} \rangle \langle \mathbf {R} +\mathbf {a} +\mathbf {\tau } |+|\mathbf {R} \rangle \langle \mathbf {R} +\mathbf {b} +\mathbf {\tau } |+\mathrm {H.c.} \right).$

Here $\langle \mathbf {r} |\mathbf {R} \rangle =\Psi _{p_{z}}(\mathbf {R} -\mathbf {r} )$ is a wave function of $p_{z}$ orbital on the sublattice $\mathrm {A}$ atom at site $\mathbf {R}$, while $\langle \mathbf {r} |\mathbf {R} +\mathbf {\tau } \rangle$ is similar state on the $\mathrm {B}$ sublattice atom, and $t$ is the "hopping integral" from a state to an adjacent similar state.

In the simplest model discussed here, the hopping matrix element couples states on the $\mathrm {A}$ sublattice only to states on the $\mathrm {B}$ sublattice, and vice versa. The hopping is chosen as $t\approx 2.7$ eV in order to match the first principles computed band structure around the corners of the Brillouin zone (to match the density functional theory computed band structure over the whole Brillouin zone, one has to use the third-nearest neighbor tight-binding model - see Lecture 5 for more details).

Diagonalization of the tight-binding model of graphene: LCAO method

The symmetry arguments above show why we were able to truncate the basis to a single $p_{z}$ orbital per atom, as long as we cay about the $\pi$ band near the Fermi level. Since there are two atoms in the unit cell, we can form two candidate Bloch orbitals:

$|\mathbf {k} A\rangle ={\frac {1}{\sqrt {N}}}\sum _{\mathbf {R} }e^{i\mathbf {k} \cdot \mathbf {R} }|\mathbf {R} \rangle ,$
$|\mathbf {k} B\rangle ={\frac {1}{\sqrt {N}}}\sum _{\mathbf {R} }e^{i\mathbf {k} \cdot \mathbf {R} }|\mathbf {R} +\mathbf {\tau } \rangle .$

The terms "Bloch" indicates that they are eingestates of the translation operator. Thus, for the basis functions from which we can construct approximate (in the sense that we assume that no other orbitals can mix into the energy eigenstates near the Fermi energy) Bloch eigenstates of the Hamiltonian:

$|\mathbf {k} \rangle =\alpha |\mathbf {k} A\rangle +\beta |\mathbf {k} B\rangle .$

where the expansion coefficients (normalized by $|\alpha |^{2}+|\beta |^{2}=1$ are eigenvectors of the $2\times 2$ Hamiltonian matrix:

${\hat {H}}(\mathbf {k} )={\begin{pmatrix}\langle \mathbf {k} A|{\hat {H}}|\mathbf {k} A\rangle &\langle \mathbf {k} A|{\hat {H}}|\mathbf {k} B\rangle \\\langle \mathbf {k} B|{\hat {H}}|\mathbf {k} A\rangle &\langle \mathbf {k} B|{\hat {H}}|\mathbf {k} B\rangle \\\end{pmatrix}}$

The diagonal matrix elements are zero:

$\langle \mathbf {k} A|{\hat {H}}|\mathbf {k} A\rangle =\langle \mathbf {k} B|{\hat {H}}|\mathbf {k} B\rangle =0,$

while the off-diagonal are given by:

$\langle \mathbf {k} A|{\hat {H}}|\mathbf {k} B\rangle ={\frac {1}{N}}\sum _{\mathbf {RR^{\prime }} }e^{i\mathbf {k} \cdot (\mathbf {R^{\prime }} -\mathbf {R} )}\langle \mathbf {R} |{\hat {H}}|\mathbf {R} +\mathbf {\tau } \rangle =-t\left(1+e^{i\mathbf {k} \cdot \mathbf {a} }+e^{i\mathbf {k} \cdot \mathbf {b} }\right)=e(\mathbf {k} ).$

The single particle Bloch eigenenergies are thus $\varepsilon (k_{x},k_{y})=\pm |e(\mathbf {k} )|$, where

$|e(\mathbf {k} )|/t={\sqrt {1+4\cos \left({\frac {{\sqrt {3}}k_{x}a}{2}}\right)\cos \left({\frac {k_{y}a}{2}}\right)+4\cos ^{2}\left({\frac {k_{y}a}{2}}\right)}}.$

This band structure can be plotted in a usual way (the so-called "spaghetti diagram") shown in Fig. 4(b), along lines defined by the high symmetry points on the edges of the "irreducible wedge" in Fig. 4(a). Energy at any other $\mathbf {k}$-vector outside of the wedge is obtained from $\varepsilon (U\mathbf {k} )=\varepsilon (\mathbf {k} )$, where $U$ is the point group symmetry operation of the crystal. Similarly, energy at any $\mathbf {k}$-vector outside of the First Brillouin zone is obtained from $\varepsilon (T\mathbf {k} )=\varepsilon (\mathbf {k} )$, where $T$ is the translation group symmetry operation of the crystal. In fact, since band dispersion $\varepsilon (k_{x},k_{y})$ for 2D crystal is just a function of two variables, we can plot its complete surface as shown in Fig. 4(c).

FIG 4: (a) Irreducible wedge within the first Brillouin zone of the honycomb lattice. (b) Conventionally plotted band structure of graphene along the lines connecting points of high symmetry. (c) Full band dispersion over the whole Brillouin zone for $\pi$ (lower surface or valence band) and $\pi ^{*}$ (upper surface or conduction band) bands of graphene.

Let us write $e(\mathbf {k} )=|e(\mathbf {k} )|e^{i\phi (\mathbf {k} )}$. The Schrödinger equation is then:

${\hat {H}}(\mathbf {k} )|\Psi _{\mathbf {k} }\rangle =\varepsilon (\mathbf {k} )|\Psi _{\mathbf {k} }\rangle$ 

with the Hamiltonian matrix being:

${\hat {H}}(\mathbf {k} )={\begin{pmatrix}0&|e(\mathbf {k} )|e^{i\phi (\mathbf {k} )}\\|e(\mathbf {k} )|e^{-i\phi (\mathbf {k} )}&0\\\end{pmatrix}}.$

The eigenvectors $|\Psi _{\mathbf {k} }\rangle$ are given by:

$|\mathbf {k} -\rangle ={\frac {1}{\sqrt {2}}}{\begin{pmatrix}e^{i\phi (\mathbf {k} )/2}\\-e^{-i\phi (\mathbf {k} )/2}\\\end{pmatrix}},|\mathbf {k} +\rangle ={\frac {1}{\sqrt {2}}}{\begin{pmatrix}e^{i\phi (\mathbf {k} )/2}\\e^{-i\phi (\mathbf {k} )/2}\\\end{pmatrix}},$

Note that the phase factors $1+e^{i\mathbf {k} \cdot \mathbf {a} }+e^{i\mathbf {k} \cdot \mathbf {a} }$ become the 1/3 roots of unity, $1+e^{i2pi/3}+e^{-i2pi/3}=0$, when $\mathbf {k}$ lies on corner points of the Brillouin zone, K and K' (the other four are equivalent to those as they can be reached from them by adding the reciprocal lattice vector). Therefore, the energy $e(\mathbf {k} )=0$ at the zone corners. Everywhere else in k-space, $e(\mathbf {k} )\neq 0$ and the splitting of two graphene bands [$\pi$ and $\pi ^{*}$ or conduction and valence, respectively, in Fig. 4(c)] lying symmetrically above and below the Fermi energy $\varepsilon _{F}=0$ is $2|e(\mathbf {k} )|$.

A perfect undoped graphene sheet has one electron per carbon in the $\pi$ level. Therefore, the Fermi level is between the two symmetrical bands, with zero excitation energy needed to excite an electron from just below the Fermi energy to just above at the $K$ points in the corners of the Brillouin zone. Note that the $(\pi \pi ^{*})$ degeneracy at isolated K points forming the Fermi "surface" is a general to the single particle description of graphene. It follows from symmetry, and it is not just accidental result of the simply nearest-neighbor tight-binding model. For example, in density functional theory calculations (see Lecture 5), there is no particle-hole symmetry, but the degeneracy at the K and K' points persists.

Massless Dirac fermions as low-energy quasiparticles and their Berry phase

The wavefunctios of low-energy quasiparticles (i.e., for energy close to the Fermi energy or charge neutrality point $\varepsilon _{F}=0$) have attracted a lot of interest. Let us consider a circular path in k-space around $K$ point. The energy is a linear function of $\delta \mathbf {k} =\mathbf {k} -\mathbf {K}$. If we move adiabatically in k-space around this point, the wavefunctions acquire a "Berry phase" $e^{i\pi }=-1$ when completing a circle. This can be seen by expanding the wavefunctions to first order in $\delta \mathbf {k}$:

$-{\frac {e(\mathbf {k} )}{t}}=1+e^{(\mathbf {K} +\delta \mathbf {k} )\cdot \mathbf {a} }+e^{(\mathbf {K} +\delta \mathbf {k} )\cdot \mathbf {b} }\approx 1+\left(-{\frac {1}{2}}-i{\frac {\sqrt {3}}{2}}\right)(1+i\delta \mathbf {k} \cdot \mathbf {a} )+\left(-{\frac {1}{2}}+i{\frac {\sqrt {3}}{2}}\right)(1+i\delta \mathbf {k} \cdot \mathbf {b} )\approx -{\frac {1}{2}}i\delta \mathbf {k} \cdot (\mathbf {a} +\mathbf {b} )+{\frac {\sqrt {3}}{2}}i\delta \mathbf {k} \cdot (\mathbf {a} -\mathbf {b} )={\frac {\sqrt {3}}{2}}|\delta \mathbf {k} |ae^{-i\pi /2}e^{-i\theta },$

where $\theta$ is the angle from the x-axis in the $k_{x},k_{y}$ plane. Here the formulas $\mathbf {a} +\mathbf {b} ={\sqrt {3}}a\mathbf {i}$ and $\mathbf {a} -\mathbf {b} =-a\mathbf {j}$ were employed. Thus, $|e(\mathbf {k} )|={\sqrt {3}}|\delta \mathbf {k} |at$, and $e(\mathbf {k} )$ has phase $\phi =-\theta -\pi /2$. Since the wavefunctions have phases $\pm \phi /2$, they change phase by $\pi$ when $\theta$ increase by $2\pi$, that is, when the $\mathbf {k}$-vector goes once around a loop surrounding a $K$ point.

The code word massless Dirac fermion is used to refer not only to linear $\varepsilon (\delta \mathbf {k} )=\pm v_{F}\hbar {\sqrt {\delta \mathbf {k} _{x}^{2}+\delta \mathbf {k} _{y}^{2}}}$ dispersion, but it also takes into account that the quasiparticle Hamiltonian around the Dirac points K and K' has a form of the Weyl Hamiltonian, i.e., the Dirac Hamiltonian of a relativistic particle with zero mass:

${\hat {H}}(\delta \mathbf {k} )=-{\frac {\sqrt {3}}{2}}at{\frac {1}{\hbar }}{\begin{pmatrix}0&\hbar (\delta k_{x}-i\delta k_{y})\\\hbar (\delta k_{x}+i\delta k_{y})&0\\\end{pmatrix}}=-v_{F}{\hat {\boldsymbol {\sigma }}}\cdot {\hat {\mathbf {p} }}$

where the Fermi velocity $v_{F}={\frac {\sqrt {3}}{2}}at\approx c/300$ replaces the velocity of light $c$ when the same Hamiltonian is used in quantum field theory of elementary particles (such as neutrinos).

Pseudospin, isospin and chirality of massless Dirac fermions

Although graphene linear spectrum is important, it is not the only essential feature. Above zero energy, the current carrying states in graphene are, as usual, electron-like and negatively charged. At negative energies, if the valence band is not full, unoccupied electronic states behave as positively charged quasiparticles (holes), which are often viewed as a condensed matter equivalent of positrons. Note, however, that electrons and holes in condensed matter physics are normally described by separate Schrödinger equations, which are not in any way connected (as a consequence of the so-called Seitz sum rule, the equations should also involve different effective masses).

In contrast, electron and hole states in graphene should be interconnected, exhibiting properties analogous to the charge-conjugation symmetry in quantum electrodynamics (QED). For the case of graphene, the latter symmetry is a consequence of the crystal symmetry, because graphene low-energy quasiparticles have to be described by two-component wave functions discussed above, which are needed to define the relative contributions of the $\mathrm {A}$ and $\mathrm {B}$ sublattices in the quasiparticles make-up. The two-component description for graphene is very similar to the spinor wave functions in QED, but the ‘spin’ index for graphene indicates the sublattice rather than the real spin of the electrons and is usually referred to as isospin.

FIG 6: Effective continuous (i.e., long wave length) Hamiltonian of graphene with both isospin and pseudospin degrees of freedom taken into account. Due to the isospin conservation, A-B symmetric potential cannot backward scatter chiral fermions.

This allows one to introduce chirality – formally a projection of isospin on the direction of motion – which is positive and negative for electrons and holes, respectively. The description of the electron spectrum of graphene in terms of Dirac massless fermions is a kind of continuum-medium description applicable for electron wavelengths much larger than interatomic distances. However, even at these length scales, there is still some retention of the structure of the elementary cell, that is, the existence of two sublattices. In terms of continuum field theory, this can be described only as an internal degree of freedom of the charge carriers, which is just the chirality.

The concepts of isospin and chirality, i.e., the fact that isospin direction is linked to the axis determined by the electron momentum, allows us to understand some transport experiments. For example, QED predicts the so-called Klein paradox - unimpeded penetration of relativistic particles through high and wide potential barriers - as one of the most exotic and counterintuitive phenomena. Although discussed in many contexts in particle, nuclear and astro-physics, direct observation of the Klein paradox using elementary particles have so far proved impossible. On the other hand, it has been observed in recent experiments on graphene.

The matching between directions of chirality illustrated in the Fig. 6 for quasiparticles inside and outside the barrier results in perfect tunneling. In the strictly one-dimensional case, such perfect transmission of Dirac fermions has been discussed in the context of electron transport in carbon nanotubes as well.

FIG 6: Tunnelling through a potential barrier in graphene: (a) Schematic diagrams of the spectrum of quasiparticles in single-layer graphene. The spectrum is linear at low Fermi energies (< 1 eV). The red and green curves emphasize the origin of the linear spectrum, which is the crossing between the energy bands associated with crystal sublattices $\mathrm {A}$ and $\mathrm {B}$. (b) Potential barrier of height $V_{0}$ and width $D$. The three diagrams in a schematically show the positions of the Fermi energy E across such a barrier. The Fermi level (dotted lines) lies in the conduction band outside the barrier and the valence band inside it. The blue filled areas indicate occupied states. The isospin denoted by vector $\sigma$ is parallel (antiparallel) to the direction of motion of electrons (holes), which also means that $\sigma$ keeps a fixed direction along the red and green branches of the electronic spectrum. Adapted from Nature Physics 2, 620 (2006).