Symmetry Breaking in Few Layer Graphene Films

Recently, it was demonstrated that the quasiparticle dynamics, the layer-dependent charge and potential, and the c-axis screening coefficient could be extracted from measurements of the spectral function of few layer graphene films grown epitaxially on SiC using angle-resolved photoemission spectroscopy (ARPES). In this article we review these findings, and present detailed methodology for extracting such parameters from ARPES. We also present detailed arguments against the possibility of an energy gap at the Dirac crossing ED.


Isolation of graphene
Recently, much attention has been given to the electronic and other properties of graphene. Following the isolation and dramatic transport measurements of individual graphene flakes by exfoliation [1,2,3], there has been an explosion of theoretical and experimental interest in graphene. Among the interesting properties found are the massless, relativistic nature of graphene's carriers, and the impact of Berry's topological phase factor on the transport properties of single and bilayer graphene. Especially interesting from a technological point of view is the extremely long lifetime of carriers, due to weak backscattering arising from their chiral nature [4,5]. This chiral nature derives from special symmetry properties of the graphene lattice.
Exploitation of these effects for electronic devices requires the precise and scalable control of graphene nanostructures, which cannot as yet be achieved with exfoliated flakes. Therefore, much attention has been given to the epitaxial growth of graphene on various substrates. Forbeaux et al. were the first to demonstrate that high-quality epitaxy of single and few-layer graphene (FLG) could be achieved on the silicon-rich SiC(0001) surface [6]. Transport measurements and demonstration of the feasibility of patterned graphene devices were demonstrated by Berger et al. [7,8]. Fig. 1 shows the atomic arrangement in monolayer and bilayer of graphene. The unit cell consists of two equivalent C atoms, labelled A and B with bond length 1.42Å. Jones proved that for a closed-packed hexagonal lattice, the energy gap along the zone boundary disappears where bands from adjacent unit cells cross [9]. This is illustrated in Fig. 2, which shows the tight-binding (TB) band structure E(k) of graphene, evaluated to third nearest neighbor using the parameters of Reich [10]. (Here we restrict consideration to the π and the π * states, which are derived from the p z orbitals of the carbon atoms [11]). Quantitative fits of the TB model to experimentally determined bands were presented by Bostwick et al. (Ref. [12]). These states meet at the so-called Dirac crossing point at energy E D in agreement with Jones' theorem. For neutral (undoped) graphene, the Fermi energy (the energy of the least-bounds states) Many of the interesting properties of graphene revolve around the fact that the band crossing at E D is strictly gapless. This means that at zero doping and zero temperature, graphene is a gapless semiconductor or a zero-overlap semimetal. Upon doping the graphene by either deposition of foreign atoms [13,14], molecules [15] or in a gated geometry [1,2,3], the carrier density can be easily manipulated. With this control, we can systematically study the many-body interactions in graphene as a function of doping using angle-resolved photoemission spectroscopy.

Angle-resolved Photoemission Spectroscopy
The Fermi surface is defined as a constant energy surface E(k)| E=E F , and determines all the transport properties of conducting materials. While transport measurements on doped graphene can determine the relevant properties such as group velocity and lifetime Figure 2. Theoretical tight-binding band structure for graphene, based on third NN parameters due to Reich [10].
of carriers on the Fermi surface, angle-resolved photoemission spectroscopy (ARPES) is a useful complementary tool. It can determine the electronic band structure, so it is capable of measuring not only the group velocity and Fermi surface, but also the constant-energy surfaces for all occupied states and the full occupied bandstructure E(k). Furthermore, the technique also accesses important information about manybody effects [16]. When there is sufficient energy and momentum resolution, the experimentally determined spectral width of the Fermi contours can be taken to be the inverse of the mean free path, and the measurement of E(k) is taken as a measure of the many-body spectral function A(k, ω).
This spectral function is in turn related to an electronic self-energy Σ(k, ω) as follows (see Ref. [17] and therein): where ω is the measured binding energy and ω b (k) is another energy defined below (whereh = 1). We make the approximation that Σ(k, ω) is k-independent. In this form, we see that A(k, ω), when evaluated at constant energy ω, is a Lorentzian function whose width is given by ImΣ representing the inverse lifetime (proportional to the inverse mean free path). Eq. 1 is valid when the scattering rate of the charge carriers (expressed in energy units) is not too large compared to their energy. In this situation, we refer to the charge carriers as quasiparticles (QPs). In our measurements, the QPs are holes which have been injected as part of the photoemission process. Their binding energy ω (hereh = 1) is taken as a negative number, and we speak of increasing energy as an increase in this negative value.
One can draw an analogy between QPs propagating in a scattering medium and light traveling in a lossy optical medium. Such a medium is characterized by a complex dielectric function, and the effects on the light propagation are not only through its absorption but also a dispersion. To satisfy causality, the real and imaginary parts of the dielectric function are related by a Hilbert transformation. Similarly, the propagation of QPs in a scattering medium leads not only to inelastic scattering (whose lifetime is encoded in ImΣ) but also renormalization of the carrier's energy, encoded in the real part of Σ(k, ω). These real and imaginary parts of Σ(k, ω) are also related by a Hilbert transform, and the function ReΣ is defined as the difference between the measured carrier energy ω and the "bare" band energy ω b (k) (that is, in the absence of scattering interactions), as indicated in Eq. (1). Following this formalism, ARPES can determine the energy-dependent lifetime due to scattering from other excitations in the system.
For a valid spectral function analysis, the ARPES spectra must be acquired with sufficient resolution and the samples must be of high quality so that defect scattering is negligible. They must also be well-characterized in thickness to ensure that the pure graphene signal is accessed.
The first ARPES measurements on FLG on SiC were from thick films [18,19,20] aimed toward studying the properties of graphite. Later, Rollings et al. [21] measured the Fermi surface and other constant energy surfaces around E D for a film with thicknesses around 2-3 layers, determined by core level shifts of C 1s electrons. Systematic core level and valence band offset studies were carried out around the same time by Seyller et al. [22]. Because of the contribution of carbon from the SiC substrate to the core level signal, such measurements give a rough measure of the film thickness but cannot give a precise thickness measurement.
As shown below, the ARPES measurements themselves can give not only a precise thickness determination, but also determine other crucial quantities. The initial formation of the graphene valence band from the silicon-rich SiC surface through to the first monolayer graphene was by Emtsev et al. [23]. Valence band measurements to discriminate film thicknesses greater than 1 monolayer were first shown by Ohta et al. for bilayer [13] and later systematically for monolayer-quadlayer graphene films [24].
These studies also demonstrated the crucial role of substrate preparation for good quality valence band measurements. The first detailed spectral function by ARPES from graphene were published by Bostwick et al. [14] and could show a rich spectrum dominated by electron-electron, electron-phonon, and electron-plasmon scattering.

Experimental
Here we briefly review the growth method of graphene on SiC in our work [13,14,24]. Films were grown on n-type (∼ 10 18 cm −3 N atoms) SiC(0001) wafers which were precleaned by annealling in 1 bar of Hydrogen gas at 1550 • C for around 30 minutes. The role of this cleaning step is essential, as by etching it removes the polishing scratches while maintaining bulk SiC stoichiometry. As-cleaned substrates were found to be atomically flat with wide terraces (Ohta et al., unpublished). Formation of graphene layers by heating to around 1200 • in ultrahigh vacuum was monitored with low energy electron diffraction (LEED) following Forbeaux [6] and ARPES as described below. The base pressure of our system was 1-2 × 10 −11 T, and graphene growth was always performed at pressures better than 1 × 10 −10 T. All measurements were obtained at phonon energy hν =94 eV unless otherwise noted.

Monolayer Graphene: a gap at E D due to symmetry breaking?
As Forbeaux et al. showed, FLG formation is accompanied by a 6 √ 3 × 6 √ 3 reconstruction at the graphite-SiC interface, which was initially attributed to the interference between the graphene and SiC lattice constants [6]. We now know from photoemission [23], theoretical calculations [25,26], and scanning tunnelling microscopy (STM) [27] that the 6 √ 3 × 6 √ 3 represents a non-interacting "0 th " graphene layer whose electronic structure resembles graphene only insofar as it has an intact σ-like bands (derived from sp 2 -hybridized in-plane bonds) but lacking the π bands characteristic of the out-of-plane p z states of graphene. The presence of such a 0 th layer is important because it acts like a dead layer, saturating or interacting with the underlying SiC bonds while forming a template for a subsequent first graphene overlayer. From symmetry considerations, it is known that the π bands from the latter and the σ bands of the former cannot interact. Therefore, the first graphene layer's chemical interaction with the substrate is very weak, and therefore we expect the π bands of graphene on SiC are to a very good approximation the same as those of freestanding doped graphene. In the following, we do not count this dead 0 th layer among the active graphene layers in our FLG film.
It is well-known that the Hamiltonian of one layer of graphene near the K point of the Brillouin zone can be approximated [4,5,28,29] by where the wavefunctions Ψ = (ψ A , ψ B ) are written in terms of p z orbitals centered on the A and B atoms in the graphene basis set. The parameter ∆ represents a possible asymmetry between the A and B sites. For ordinary graphene, ∆ = 0 since the atoms are indistinguishable. The off-diagonal terms represent the hopping between A and B sublattices, and v is the band velocity around E D . The Hamiltonian in Eq. 2 leads to a conical bandstructure E(k)= vk when ∆ = 0. Here k is the momentum relative to one of the K points at the corner of the graphene Brillouin zone (see Fig. 2). Experimental Fermi surfaces and underlying bandstructures for clean and alkali-dosed graphene are shown in Fig. 3(a-b), adapted from Ref. [14]. We can immediately recognize the expected nearly linear dispersions as well as the Dirac crossing points (middle panels) in the bands at the Dirac energy E D . We also see that there is a non-trivial change in intensity when traversing around the Fermi contour. This will be discussed in detail below, but for now we regard it as a photoemission cross section effect. Because of this effect, when we sample the bandstructure in the y-direction (relative to Fig. 3), we see only one of the two expected bands; the other is extinguished (right panels).
We also observe that even the clean, as-grown graphene films have a Fermi level E F significantly above (by around 0.45 eV) the Dirac energy E D . This in-built doping was first reported by Rollings et al. [21] and can be attributed to the greater electron affinity of graphene compared to the substrate. Our experiments have shown that this intrinsic n-doping is independent of whether the substrate dopants have been frozen out (at T ∼ 40K). Since its discovery by ARPES, this intrinsic n-doping has also been predicted theoretically [25,26]. An important feature of the one monolayer data is the appearance of kinks in the energy band structure below E F [14]. These kinks occur at two energy scales. First we see a slight kink at ∼ 200 meV below E F . This kink is hardly visible on the large energy scale plotted in Fig. 3, but it is accompanied by pronounced sharpening between 200 meV and E F that is readily observed. This kink is similar to ones which have been observed at similar energy scale in graphite [30,31] and bilayer graphene [13] that have been attributed to electron-phonon scattering. We can understand the presence of the kink within the spectral function formalism in Eq. 1, noting that there is an observable increase in linewidth of the band at binding energies greater than 200 meV, signifying a decrease in the lifetime of the states as electrons absorb or emit phonons. We will discuss this formalism further below but for now it is sufficient to identify this feature with phonons for two reasons: first because of the energy scale, which corresponds to the in-plane LO and TO phonons, and second, because the kink feature's energy scale remains constant with doping, as expected for the small doping levels considered here.
There is a second anomaly in the dispersion around the Dirac crossing point in Fig. 3. In the middle panels, where both bands have equal intensity, the region of the crossing of the bands seems spread out in energy. In the right panels, where one band is extinguished, it appears that this spread is associated with a second kink feature which is at the Dirac crossing point of the bands. Unlike the phonon kink, this anomaly moves to higher binding energy with doping, and must therefore be somehow associated with the Dirac energy E D . Similar to the phonon kink, it is stronger at higher doping, and it is associated with a change in linewidth-the bands are locally broadened around E D .
What causes this second feature? In Bostwick et al., it was proposed to be a kink due to electron-plasmon interaction [14] but it has been recently proposed that the observed spreading of the bands around E D is associated with substratedependent energy gap at E D [32,33]. Such a gap would be interesting because it suggests an electronic or chemical control of the electronic character (2D semimetal vs. semiconductor) and is proposed on the basis of possible symmetry breaking. First, we discuss this idea and then present the evidence against it followed by evidence in favor of the electron-plasmon scattering model.
Within the simple Hamiltonian (Eq. 2) a gap of magnitude ∆ appears at the Dirac crossing energy E D when the parameter ∆ = 0. A physical interpretation of this gap is the symmetry breaking of the A and B atoms. This occurs for, e.g. replacement of C atoms with B and N in hexagonal boron nitride. It also occurs in a scenario where the bonding of A and B atoms to the 0 th layer is asymmetric as recently proposed [32,33].
We present arguments against this scenario in graphene on SiC as follows. (1) The interaction between the 1st and 0th layer is very weak. This was established by ARPES [23], theory [25,26], and STM [27]. A possible argument against the weak interlayer attraction is the appearance in monolayer graphene films of replica π bands with 6 √ 3 × 6 √ 3 symmetry, ascribed by us as due to diffraction of the outgoing photoelectrons [14], similar to other nearly incommensurate systems [34]. These satellite bands lead to replicas of the constant energy contours, illustrated for the Dirac crossing energy in Fig. 4. With a linear grey scale in (a) the replica bands are hardly noticable but with a highly non-linear grey scale (b), they can be emphasized.
It is tempting to ascribe the replica bands to a possible 6 √ 3 × 6 √ 3 superlattice potential felt by the first graphene layer. If this were true, additional energy minigaps would appear where the replica and main bands cross [35];however, no such gaps have been observed [13]. Furthermore, the replica bands, very weak at low temperature (about a factor 40 reduced intensity compared to the primary band) do not appear at room temperature [23] and at 100K are dramatically broadened [J. L. McChesney, unpublished]. This observation violates the hypothesis that the first graphene layer has 6 √ 3 × 6 √ 3 symmetry potential which would demand the linewidths of the replica and main bands to be identical by symmetry. Instead, we can easily understand the broadening of the replica bands as due to a Debye-Waller factor, confirming the origin of these replica bands as due to final-state diffraction.
(2) The doping dependence shows a clear increase in the spread of the states at the Dirac crossing. If this spread were due to a gap from coupling to the substrate, the coupling strength should be independent of the doping density (or become smaller due to enhanced screening). (3) We observed that the bands above and below E D are misaligned [14], so that the projections of the π states below E D do not pass through the π * states above E D . This is illustrated by the dashed lines in Fig. 5(a), which reproduces the clean graphene bandstructure. This misalignment does not occur in the energy gap scenario, but comes naturally when many-body interactions are present.
(4) The density of states (DOS) does not show a gap at E D . This is illustrated in Fig. 5(b) for the momentum-integrated DOS. In a gap scenario one expects a decreased DOS but we see a peak (expected for crossed bands).
(5) The energy distribution curve (EDC) at the Dirac crossing shows only a single peak, not a split peak as expected in a gap scenario (see Fig. 5(c)).
(6) The intensity distribution along the Fermi surface provides a stringent test for A-B atom symmetry breaking. It is observed that one side of the Fermi contours is very weak or absent. In the strictly symmetric case ∆ = 0, the intensity on one side of the Fermi contour is strictly zero. Rather than a simple vanishing photoemission matrix element, the cancellation results from the interference between emission from A and B sites, as shown by Shirley [36]. This cancellation, like the Dirac nature of the quasiparticles, and the lack of backscattering, follow from the strict A-B atom symmetry. If we break the A-B atom symmetry, we not only open a gap at E D (thereby destroying the massless character), but also destroy the phase cancellation affecting the Fermi surface intensity.
These effects are illustrated in Fig. 6. In (a) we show as polar maps the measured angular distribution of the band intensity taken about the K point for monolayer and bilayer graphene (closed and open circles, resp.). These data were obtained by fitting the momentum distribution curves taken along radial cuts for an energy window ∼ 75 meV below E F . The monolayer and bilayer Fermi surfaces are practically identical, but as indicated in the figure, the bilayer signal is not completely extinguished in any direction. In contrast, for the monolayer, the intensity is completely extinguished in one direction, apart from a very weak minority contribution from bilayer regions. This residual bilayer signal is easy to subtract since it is well separated from the monolayer bands below E D [13]. After subtraction, we determined that the monolayer band intensity is zero within a very low noise floor (about 0.15% as indicated by the central yellow circle).
Shirley derived a simple formula for the symmetric case ∆ = 0 for monolayer graphene; we extended this model to the case of finite ∆ and show in Fig. 6(a) the expected angular distributions for a ∆ = 0.0, 0.1, 0.2 eV (leading to energy gaps at E D of the same values). This plot shows that we are fairly sensitive to the possible symmetry breaking (and this sensitivity can be enhanced simply by acquiring the bands with better statistics). Fig. 6(b) shows a plot of the intensity reduction as a function of ∆, which can be compared to our noise floor (< .015%). From this comparison, we can conservatively estimate the maximum gap at E D to be under 60 meV. Since the apparent kink at E D (with a resulting spreading of the states there) is much wider in energy than this, we can rule out the symmetry breaking as being the dominant factor to explain the anomalous dispersion at E D .
As an aside, the reason the bilayer is not completely extinguished is that A-B atom symmetry breaking is indeed violated for the bilayer. That is because only one atom (A, say) in the outer layer is directly over an atom in the inner graphene layer (see Fig.  1). This symmetry breaking also explains the well-known symmetry breaking in STM images of bilayer and thicker films [27,36]. (A complete model of the bilayer angular intensity profile is outside the scope of this paper and will be presented elsewhere.) (7) It is worth pointing out the very high momentum resolution and accuracy of positioning of the sample that is required to obtain spectra precisely at E D . In Fig.  5(a), we see that the entire span of the Fermi bands is only about 0.1Å −1 . Only a small misalignment on the order of 0.05 • could result in an apparent gap in the bands.

The case for self-consistency.
Having ruled out the gap scenario, we can now consider many-body interactions to explain the kinked dispersions around E D . The first issue is whether a self-consistent model is possible even in principle. We will first establish that the kinks and the linewidth variations are consistent with each other. As discussed above and in the literature [17,37], we analyze the spectral function data in terms of real and imaginary parts of the self-energy function Σ(k, ω). Fig. 7(b) shows an experimentally acquired spectral function A(k, ω) for relatively highly doped graphene (n = 5.6 × 10 13 cm −2 ). The dressed band position ω =ω b (k)+ ReΣ(k, ω) is determined by fitting momentum distribution curves (MDCs, that is, individual constant-energy slices) to Lorentzian functions. The positions are plotted in Fig. 7(a) (black line) and the Lorentzian width as a function of ω is plotted in Fig. 8(a). In order to converge to a self-consistent interpretation, it is necessary to iteratively apply the following procedure. We take a second order polynomial as a trial bare band ω b (k). Given this ω b (k), we can easily scale the MDC widths (units ofÅ) into the function ImΣ(k, ω) (units of eV), shown in 8(b). This function is smoothed and then Hilbert transformed into a trial ReΣ(k, ω) function. We can also extract the ReΣ(k, ω) function by subtracting the trial bare band from the fitted band position. These two ReΣ(k, ω) functions (Fig. 7(c)) are compared and the trial bare band adjusted until the model ReΣ(k, ω) and ImΣ(k, ω) are in good agreement with the experimentally extracted curves as plotted in Fig. 7(b-c).
As a final check of self-consistency, we can generate a trial spectral function A(k, ω) from only the fitted MDC widths and the mathematically transformed ReΣ(k, ω), shown in Fig. 7(b). It is in excellent agreement with the experimental function in Fig. 7(a). Having demonstrated this self-consistency, we can say with high degree of confidence that the kink anomalies must be attributed to many-body interactions, and not any details of the single-particle bandstructure. That is, we can safely rule out not only the . Comparison of calculated and measured MDC widths. (a) measured MDC widths (dots) for the highest-doping sample (n = 5.6 × 10 13 cm −2 ) are compared to the total scattering rate contribution from Bostwick et al. [14](solid). (b) the calculated contributions to the scattering rate due to electron-hole pair generation, electron-phonon scattering, and electron-plasmon scattering [14]. (c-e) experimental MDC widths for n =1.2, 3.0, and n = 5.6×10 13 cm −2 in comparison to the calculations of Hwang et al. [38]. Adapted from Ref. [14].
superlattice gap scenario outlined above, but also strain, defects and other initial-state effects.

Contributions to scattering lifetime.
We now explain the physical origin of the measured ImΣ(k, ω) function in more quantitative fashion. For convenience we work with the fitted MDC widths, which are plotted in Fig. 9(a). The features to explain are, starting from zero energy, the monotonic increase in lifetime down to about −0.2 eV; the hump at around E D = −1.0 eV, and the remaining background rise. These were attributed [14] to electron-phonon (e-ph) coupling, electron-plasmon (e-pl ) coupling, and electron-hole (e-h) pair generation; computations of these contributions to the scattering rate are shown in Fig. 9(b). Schematic diagrams of these processes are shown in Fig. 10.
We can meaningfully consider only those excitations whose energy scale is greater than our energy resolution (∼ 25 meV). Considering the energy scale of the observed kink anomalies (≥ 200 meV) we can rule out any significant interactions between 25 to 200 meV, such as scattering from low-energy acoustical vibrational modes.
First, we qualitatively discuss the coupling to phonons at the 200 meV energy scale (a quantitative analysis has been presented by Bostwick et al. [12] for graphene, and for other surfaces in Refs. [39,40,41]). Since this energy is much larger than our temperature (k B T ∼ 2 meV, we can rule out phonon absorption and consider only decay of quasiparticles (QPs) by phonon emission (Fig. 10(a)). Such QP decays are forbidden for states at E F , but become available as the quasiparticle energy increases. Typically once the quasiparticle energy is greater than the phonon energy scale, the lifetime due to scattering is independent of QP energy. This change in QP lifetime is reflected in a monotonic increase in the imaginary part of the self energy ImΣ(k, ω). Because the real and imaginary parts of Σ(k, ω) are related by Hilbert transform, one expects to see a non-trivial change in the dispersion on the phonon energy scale, which is observed as a kink. Physically, we interpret the change of band slope between the kink and E F as a renormalization of the mass as the QPs become "dressed" with a virtual cloud of phonons. But we know that the QPs in graphene are effectively massless, so we speak in terms of a velocity renormalization (or equivalently a renormalization of the relativistic mass-equivalent energy).
The 200 meV kink is stronger for the K-covered graphene compared to the as-grown surface (see Fig. 3(a,b)) due to a phase-space argument. The density of electronic states in k-space is a constant, so that as the sample is doped, the bands span more electronic states near E F ; as these become available final states for phonon scattering, the decay probability increases. Left unexplained is the overall magnitude of the e-ph scattering rate at all dopings, which is about 6 times stronger than what is predicted by the normal deformation potential calculations [12,42,43].
The quantitative analysis of the phonon kink [12,14], which followed the standard formalism [44], is quite useful but does not perfectly describe the kink strength (it underestimates it slightly) and furthermore does not take into account the actual band structure of graphene: the actual phonon scattering rate should diminish near E D from the same phase-space argument just cited. A first-principles calculation of the phonon scattering rate should account for both of these effects.
In the case of the second kink near E D , the QP decay is through emission of plasmons ( Fig. 10(b)), which is favored over phonon scattering because of the kinematic constraints [14,45,38]. Whereas optical phonons are more or less delocalized in k-space with fixed energy scale, the plasmon spectrum in graphene has a fast energy dispersion in a narrow range of k. This follows from the dispersion relation for two-dimensional plasmons [46] in the long wavelength limit: where q is the plasmon momentum, m is the carrier mass, and ǫ ∼ 6 is the dielectric constant [12]. Although plasmons in principle exist in the domain 0 < q < ∞, in practice they propagate freely up to a critical momentum q < q c due to Landau damping (plasmon decay into electron-hole pairs) [47]. For graphene, the rest mass m is zero near E D but the relativistic mass equivalent to the kinetic energy, m r =E/v 2 (where v is the Fermi velocity), is on the order [2,3] of 0.1m e and can be used to set the plasmon energy scale ω pl . This means that more or less vertical interband decays by plasmon scattering are now the dominant factor determining the lifetime near E D .
Two other contributions to the scattering lifetime must be considered. First is defect or impurity scattering. Normally the defect scattering is taken to be a constant background to the imaginary self energy ImΣ(k, ω), which is a deconvolution of the residual momentum spread of the bands at E F and the instrumental function. In our case, the residual momentum spread is only about 0.005Å −1 , which is comparable to our instrumental resolution, so we can safely discard the defect scattering rate as negligible.
The remaining contribution to the scattering rate is the decay by e-h pair generation, which is the standard decay process in Fermi liquid theory. In this process (Fig. 10(c)), the decay of the quasiparticle is accompanied by an excitation of an electron above E F , creating a new hole below E F . For two-dimensional metals with circular Fermi surface, Hodges et al. proved the famous rule that e-h-pair scattering rate goes as ω 2 lnω [48]. This was determined by a phase space integration of all possible kinematically allowed e-h-scattering processes. For a 2D free electron gas this could be carried out analytically, but for graphene, we evaluated the appropriate integral numerically. This was done so that we could use the experimentally determined dispersion (although we assumed cylindrical symmetry and zero temperature to simplify the integration).
The most interesting finding is that for n-doped graphene, the e-h-scattering rate rises from E F down towards E D as it would be expected for any metal. Around E D , however, the scattering rate must necessarily drop, because in the vicinity of E D , the decays are mostly vertical transitions. Such a decay by e-h-pair generation is forbidden since we cannot find a momentum-conserving excitation near E F to satisfy the kinematic constraints. Only at energy scale around twice the Dirac energy do such excitations become available, and we see an associated rise in the scattering rate at high energy scales.
Considering the simplicity of the model, the total scattering rate function ( Fig.  9(b)) does a remarkably good job to model the data. Theoretical modelling of the e-pl and e-h scattering rates has also been performed by Hwang et al. [38]. Fig. 9(c-e)) shows a comparison of our measured MDC widths to their model for three different dopings. Although they overestimate the relative contribution of the e-pl to e-h processes, their calculation is in excellent qualitative agreement with the observed MDC widths. The main discrepancy is the failure of the model to account for the scattering rate from phonons, which was not included in their calculation.
The many-body effects we measure are present all the way down to zero doping and therefore may play a role in the transport of gated graphene devices. These are much more dilute carrier gases than we achieved by alkali metal doping. As the doping level decreases, the phonon and plasmon processes will overlap in in energy and therefore will not be separable. This is already seen in the lowest doping we probed (Fig. 9(c)). The plasmon and e-h pair scattering rates are reasonably separable at all dopings, but there is an energy overlap region just above E D where neither alone is a good description of the total electron-electron interaction. These observations imply that a proper description of the scattering rate will take into account much more complicated processes than in our simple treatment. In the language of Feynman diagrams, it means higher-order diagrams than are typically considered will be necessary to model the photoemission data. In addition, when E F is reduced to be comparable to the temperature, thermal excitation effects will increase in importance. This has already been discussed in relation to plasmons [45].

Out-of-plane symmetry breaking in multilayer graphene.
Multilayer graphene films grown on Silicon carbide have an obvious built-in symmetry breaking, because of the inequivalence of the two film interfaces (SiC and vacuum). Further symmetry breaking can be induced by either external field, or by growth of additional layers on top of the graphene films. Understanding these symmetry effects is important in order to exploit them for technological purposes. Extension of the Hamiltonian in Eq. 2 to multiple layers gives a simple framework to achieve this.
Extension to two layers is achieved by adding an additional hopping term between the B atoms of the first layer and the A atoms of the second layer [49,50,24]: Here α i acts with respect to the (A, B) sublattices of the i th layer, and β 0 is a 2 × 2 matrix where γ 1 is the hopping parameter between layers. The wave function now has four elements with basis set orbitals located at A 1 , B 1 , A 2 , B 2 atoms, where i is the layer number (1 or 2). There are two further generalizations of Eqn. 4. First, by adding more layers, and second by altering the stacking sequence. Adding a third layer, one couples the B atom of the second layer to the A atom of the third for conventional Bernal-type stacking (ABAB . . .) characteristic of bulk graphite. Repeating this sequence we come to the Bernal Hamiltonian for N layers, A useful generalization of Eq. 6 is where Now, if s = 0, then Eq. 7 is the Hamiltonian for Bernal stacking, while for s = 1, Eq. 7 is the Hamiltonian for Rhombohedral stacking (ABCABC . . .). We can further generalize this Hamiltonian to arbitrary stacking orders, by suitably choosing the different values of s in each block of the matrix.
In the above Hamiltonians, we have assigned to each layer its own on-site Coulomb energy E i . This allows for the possibility of a poorly screened electric field across the FLG film, which is reasonable in view of the predicted long screening lengths in this direction. It is straightforward to show that the Dirac crossing energy is E D =Tr H/2N where N is the number of layers. Fig. 11 shows the calculated bandstructures for one to four layer graphene. The calculations were for either Bernal (solid lines) or rhombohedral (dashed lines); the distinction is obviously meaningful only for films with N ≥ 3 layers. Far from E D , it turns out there is little distinction between rhombohedral and Bernal stacking. This is to our advantage, because as Fig. 11 shows, one can know the film thickness directly from band structure measurements by simply counting the number of π bands below E D . Near E D , the situation is quite different, since the two stacking types have quite different band dispersions. (Similar calculations have also been carried out with ab initio models [51,52,25,26,53]).
The detailed bandstructure around E D shows a strong sensitivity to the Coulomb energy terms E i that enter the Hamiltonian matrix (Eqs. 2, 7) [29]. This can be seen by comparing the upper and lower rows of Fig. 11 which were calculated for two cases. In the first case, the energies E i are all zero, and we find a gapless energy spectrum at E F =E D . For the lower row, we distributed a field change U=200 meV across the total film in uniform increments, which simulates FLG in a bias or inhomogeneosly doped geometry. This procedure opens gaps near E D ; for the special bilayer case N = 2, there is a complete gap at the Fermi level. This gap was proposed to be the basis of a new kind Figure 12. Gap Control of Bilayer Graphene. (a) An unbiased bilayer has a gapless spectrum, which we could observe for a doped sample which carefully balanced the field across the film. (b) For a bilayer with a field gradient, an energy gap is opened between π and π * states. of electronc switch, whereby lateral transport through the bilayer could be modulated by applying a modest field perpendicular to the film [13,54,55].
Systematic studies of the thickness and doping dependence by ARPES have been presented by Ohta et al. [13,24]. Fig. 12 shows the bilayer graphene bandstructure in two different field geometries, achieved by doping a bilayer graphene on SiC [13]. Simlar to monolayer graphene, the as-grown bilayer films have an intrinsic n-doping, which allows us to probe the states both above and below E D . Because the doped carriers are concentrated in the interface layer, the as-grown films have a field gradient across them and hence a gap at E D . Carefully balancing this charge imbalance allows us to close the gap ( Fig. 12(a)), while further doping of the surface layer allows us to create a net charge imbalance, thus reopening the gap ( Fig. 12(b)). Evidence for a similar gap opening was also presented for the surface layers of graphite when doped with Na [56].
Systematic thickness measurements at constant doping were presented by Ohta et al. [24]. For films of thickness N=1-4 layers, we found that the total charge density donated from the substrate was the same for all thicknesses. Similar to bilayer graphene, the charge was donated predominantly to the interface graphene layer. This is in accord with the metallic nature of the films, which can screen the interface layer from the rest of the film. The measured bandstructures for N = 1 − 4 layers are shown in Fig. 13.
These spectra are very rich in information: we could determine not only the number of layers (by counting the number of π states below E D straightforwardly) but also derive the stacking order. One can see easily for N=3 that there are two states (marked by circles) of equal intensity that can only be ascribed to equal populations of Bernal and Rhombohedral stacking. For quadlayer (and presumably thicker) careful analysis shows that Bernal stacking predominates. This can be taken as evidence for the role of second near-neighbor interactions to stabilize the Bernal stacking type in graphite.
The electronic information that could be derived from the data are equally rich: in analogy with the bilayer analysis, we could assign the different charge densities in Figure 13. Band structure of graphene films of thickness for (a-d) N = 1 − 4 layers, resp. Calculated bands for three configurations are shown: Bernal stackings ABAB and ABAC (blue and light blue, resp.) and Rhombohedral stackings (red).Adapted from Ref. [24]. each graphene layer, and determine the out-of-plane screening length. In the formation of the graphene films, about 1 × 10 13 carriers per cm 2 are donated to the film, with in general about 85% of the charge donated to the first layer, and most of the rest in the second layer [24].

Conclusions and outlook: Graphene, the simplest complex material.
In the last few years, there was an explosion of interest in graphene since isolation of high-quality samples was achieved and since its many novel properties were elucidated both experimentally and theoretically. Seldom does a new material come along that has such strong fundamental and practical interest. From an experimental point of view, graphene is highly attractive since unlike other low dimensional materials (such as high-mobility semiconductor two dimensional electron gases), graphene films are exposed to vacuum and can be directly probed by surface-sensitive techniques such as LEED, STM, and ARPES. ARPES has a special role to play because it is sensitive not only to the valence band energy structure but also its symmetry in k-space. Furthermore it can give direct information on the many-body interactions, such as mass renormalization. Through graphene's special sensitivity to symmetry, we could even derive much structural information such as stacking errors and electronic information such as charge density and screening length, which would be very hard to achieve with other probes.
In our opinion, graphene is unique in many ways. It is the first system to our knowledge to show electron-plasmon coupling in the ARPES signal, which suggests not only the exciting possibility of new coupling mechanisms, but also technological implications in the interaction with photons. Finally, it is a model system for correlation and many-body interactions which can supply stringent tests for theory.