# From Graphite to Graphene: The Confinement of the Fermi Surface to the Line KH

Authors: Helmut Bross
Journal: ISRN Condensed Matter Physics
 (2013)

---

## Abstract

The graphite structure is self-consistently calculated by use of the all electron Modified Augmented Plane Wave (MAPW) scheme with lattice constantsc considerably enlarged above the experimental value of graphite. Overall, the band structures of the series are found to be quite similar: the energy levels of the highly symmetric states K and H almost coincide, essentially fixing the Fermi level of the semimetallic solid. The dispersion along lines parallel to the atomic planes, already small in graphite at the experimental value of c, continues to flatten with increasing value of c. The structure with an interlayer distance enlarged by the factor 3 over the experimental value provides a good approximation of the behaviour of a monoatomic sheet. In this context, the unusual behaviour of graphene appears in a new light.


---

## Body

## 1.   Introduction

Because of its layered structure with relatively large separation between the layers, graphite has been modelled as a two-dimensional solid in the first theoretical investigations based on various methods [1]. In the following decade, the full three-dimensional crystal structure has been investigated by all modern band structure schemes: pseudo-potential [2–5], full-potential linearised augmented plane wave [6, 7], full-potential linear muffin-tin-orbital [8], full-potential linear combination of Gaussian-type orbitals fitting-function technique (LCGTO-FF) [9] which all, due to inherent systematic restrictions, produce different degrees of accuracy. To understand the extraordinary behaviour of graphene we use the reverse strategy and study how the electronic structure of the graphite structure changes with increasing interlayer distance. We find that the semimetallic behaviour is drastically reduced as the Fermi surface shrinks to the surrounding of the KH line accompanied with a loss of the dispersion in the relevant bands. Only highly accurate schemes can properly describe this phenomenon. By a proper choice of the intrinsic parameters and with some specific care the MAPW method [10, 11] satisfies this criterion. The present investigations raise doubts on previous investigations of graphite and graphene based on the tight binding scheme with the nearest and second nearest neighbours couplings [12–14].The outline of the paper is as follows: in Section2 some issues of the crystallography of graphite and graphene are described. Then we sketch some features and the proper choice of intrinsic parameters of MAPW as well as further details of the calculation. Section 3 reports the results of self-consistent calculations done with increasing interlayer distance. Especially, we discuss some common features in the band structure. Then we describe in Section 4 the strategies to deal with the semimetallic behaviour. The construction of the Fermi surface (FS) clearly illustrates the confinement in k→ space, which has far reaching consequences. Finally, Section 5 gives a short outlook on other electronic features of graphite and graphene. Unfortunately, a detailed comparison with experiments is not possible since their results are mostly discussed in context of the Slonczewski-Weiss-McClure [13, 14] model.


## 2. Electronic Structure

### 2.1. Crystallography of Graphite and Graphene

The Bravais lattice of graphite is primitive hexagonal with the basisa→1=(0,-1,0)a, a→2=(3/2,1/2,0)a, a→3=(0,0,1)c [15]. Its crystalline structure consists of equidistant planes of carbon atoms, each forming a hexagonal net, stacked in a manner such that half of the carbon atoms (type A) are located directly above each other in adjacent planes, while the other half (type B) are located above the centres of the hexagons in the adjacent plane [1]. There are four atoms per unit cell. The space group of graphite is P63/mmc=D6h4, No. 194 [15]. The corresponding point group has 24 elements and is isomorphic to D6h. Half of the operations of the point group have nonprimitive translations v→=a→3/2 associated with them, carrying from an atom on one plane to an atom of the same type on the adjacent plane.The numerical work is considerably reduced if we choose the origin to be located halfway between adjacent layers on ac-axis through atoms of type A. It is an inversion centre of the point group with the advantage that the Fourier transforms of the charge density, the crystal potential, and other quantities invariant under the operations of the space group are real-valued. In this choice of the elementary cell the atoms of the basis are located at
(1)Aatoms:R→1=a→34,Batoms:R→2=(a→1+2a→2)3+a→34
and two further atoms in inversion positions.In graphene we have only the basis vectorsa→1 and a→2. The pattern of carbon atoms looks the same as that of one plane in graphite. But there exist no crystallographic reasons to distinguish closely neighboured atoms on a hexagon ring. Again its 2-dimensional space group is composed of the point group D6h where the nonprimitive translation v→=(a→1+2a→2)/3 is associated with half of its operations.


### 2.2. Some Details of the Band Structure Calculations

The electronic structure was self-consistently determined by use of the modified augmented plane-wave (MAPW) formalism [10]. Common to most modern schemes, the constituents within the atomic polyhedron are surrounded by nonoverlapping spheres (APW spheres). Outside these spheres the Bloch functions 〈r→∣nk→〉 are approximated by superposition of plane waves, whereas within each APW sphere the plane waves are suitably augmented to properly describe the atomic character of the wave function. Specifically, in the MAPW scheme, the angular decomposition of 〈r→∣nk→〉 is a sum of products of spherical harmonics Ylm and a set of radial functions Rsl where s is associated with various orbital energies. The Rsl are generated in a spherical averaged potential, mostly by requiring that their logarithmic derivative dRsl/dr is either 1 or −1. To avoid any truncation error, the augmentation of the plane waves is only performed for the leading angular momenta l≤lmax, whereas for l>lmax the spherical Bessel functions jl are kept. In contrast to all other versions of the APW scheme [6], the full-wave function and its derivative are made exactly continuous on the surface of the APW spheres by use of additional constraints in the extremal problem, with the consequence that the corresponding eigenvalue problem is a generalised one, but is real symmetric.Test investigations in the full atomic configuration 1s2, 2s2, 2p2 of C had the results that the two lowest bands are well separated from the other bands by more than 18 Ry and do not significantly influence the valence electrons. Thus the 1s electrons are considered to constitute the core which, however, is not kept fixed in the self-consistent cycles. Thanks to the strong localisation of the core wave functions the valence states are almost fully orthogonal on the core states. Due to the strong s-p character of graphite we choose lmax=1 and a set of radial functions Rsl consisting of 4 elements for each l=0,1 with alternating sign of the logarithmic derivative.In the present study which focuses on the influence of increasing interlayer distance the choice of the plane waves used to describe a Bloch statek→  is rather subtle. We have found that the inequality
(2)(k→∥+K→∥)2+(k⊥+K⊥)2≤QMXSQ,
which limits for a given value of QMXSQ the kinetic energy of the plane waves and guarantees that the correct symmetry behaviour in k→ space is most appropriate. ∥ denotes the components of k→ and K→ parallel to the graphite sheets and ⊥ their perpendicular components. This requirement has the effect that the wave functions in the interplanar space are properly described as the number of reciprocal vectors K→ increasing due to their perpendicular components when the interlayer distance grows. For any value of the wave vector k→ the actual value of QMXSQ fixes the rank of the MAPW eigenvalue problem and thus controls the accuracy of the Bloch eigenvalues Enk→. Cursory investigations have been done with the value QMXSQ=20.16 a.u. corresponding to a maximal kinetic energy of the same value in Rydberg.The density of the electronsρ(r→) obtained from the occupied states is, as a consequence of the MAPW ansatz, a symmetric combination of plane waves everywhere, superimposed by an angular-dependent contributions that are treated in the accuracy wanted by use of a suitably chosen fine mesh of r→ points [16].The Brillouin zone (BZ) integrations over all occupied states yielding the electron densityρ(r→) and the total energy Etot are approximated by a sum over k→ points properly chosen in the irreducible wedge with the concept of Monkhorst and Pack [17] applied to the hcp structure. They are part of a hexagonal Bravais lattice spanned by the basis vectors (π/MG1a,π/MG2a,π/MG3c) which avoids locating the points H and K in the centre of the parallelepipeds. Here, the MGi,i=1,2,3 are specially chosen integers. Thus the influence of possible energy degeneracies is kept small. Within the self-consistent cycles the same value of QMXSQ is used which has been chosen to keep the fluctuation of the rank of the eigenvalue problems small.The characteristic ground-state properties are calculated within the local-density approximation (LDA) by use of the exchange-correlation functionals [18, 19]. As we found no significant difference we restrict ourselves in the following to results obtained with the functional [19].


## 2.1. Crystallography of Graphite and Graphene

The Bravais lattice of graphite is primitive hexagonal with the basisa→1=(0,-1,0)a, a→2=(3/2,1/2,0)a, a→3=(0,0,1)c [15]. Its crystalline structure consists of equidistant planes of carbon atoms, each forming a hexagonal net, stacked in a manner such that half of the carbon atoms (type A) are located directly above each other in adjacent planes, while the other half (type B) are located above the centres of the hexagons in the adjacent plane [1]. There are four atoms per unit cell. The space group of graphite is P63/mmc=D6h4, No. 194 [15]. The corresponding point group has 24 elements and is isomorphic to D6h. Half of the operations of the point group have nonprimitive translations v→=a→3/2 associated with them, carrying from an atom on one plane to an atom of the same type on the adjacent plane.The numerical work is considerably reduced if we choose the origin to be located halfway between adjacent layers on ac-axis through atoms of type A. It is an inversion centre of the point group with the advantage that the Fourier transforms of the charge density, the crystal potential, and other quantities invariant under the operations of the space group are real-valued. In this choice of the elementary cell the atoms of the basis are located at
(1)Aatoms:R→1=a→34,Batoms:R→2=(a→1+2a→2)3+a→34
and two further atoms in inversion positions.In graphene we have only the basis vectorsa→1 and a→2. The pattern of carbon atoms looks the same as that of one plane in graphite. But there exist no crystallographic reasons to distinguish closely neighboured atoms on a hexagon ring. Again its 2-dimensional space group is composed of the point group D6h where the nonprimitive translation v→=(a→1+2a→2)/3 is associated with half of its operations.


## 2.2. Some Details of the Band Structure Calculations

The electronic structure was self-consistently determined by use of the modified augmented plane-wave (MAPW) formalism [10]. Common to most modern schemes, the constituents within the atomic polyhedron are surrounded by nonoverlapping spheres (APW spheres). Outside these spheres the Bloch functions 〈r→∣nk→〉 are approximated by superposition of plane waves, whereas within each APW sphere the plane waves are suitably augmented to properly describe the atomic character of the wave function. Specifically, in the MAPW scheme, the angular decomposition of 〈r→∣nk→〉 is a sum of products of spherical harmonics Ylm and a set of radial functions Rsl where s is associated with various orbital energies. The Rsl are generated in a spherical averaged potential, mostly by requiring that their logarithmic derivative dRsl/dr is either 1 or −1. To avoid any truncation error, the augmentation of the plane waves is only performed for the leading angular momenta l≤lmax, whereas for l>lmax the spherical Bessel functions jl are kept. In contrast to all other versions of the APW scheme [6], the full-wave function and its derivative are made exactly continuous on the surface of the APW spheres by use of additional constraints in the extremal problem, with the consequence that the corresponding eigenvalue problem is a generalised one, but is real symmetric.Test investigations in the full atomic configuration 1s2, 2s2, 2p2 of C had the results that the two lowest bands are well separated from the other bands by more than 18 Ry and do not significantly influence the valence electrons. Thus the 1s electrons are considered to constitute the core which, however, is not kept fixed in the self-consistent cycles. Thanks to the strong localisation of the core wave functions the valence states are almost fully orthogonal on the core states. Due to the strong s-p character of graphite we choose lmax=1 and a set of radial functions Rsl consisting of 4 elements for each l=0,1 with alternating sign of the logarithmic derivative.In the present study which focuses on the influence of increasing interlayer distance the choice of the plane waves used to describe a Bloch statek→  is rather subtle. We have found that the inequality
(2)(k→∥+K→∥)2+(k⊥+K⊥)2≤QMXSQ,
which limits for a given value of QMXSQ the kinetic energy of the plane waves and guarantees that the correct symmetry behaviour in k→ space is most appropriate. ∥ denotes the components of k→ and K→ parallel to the graphite sheets and ⊥ their perpendicular components. This requirement has the effect that the wave functions in the interplanar space are properly described as the number of reciprocal vectors K→ increasing due to their perpendicular components when the interlayer distance grows. For any value of the wave vector k→ the actual value of QMXSQ fixes the rank of the MAPW eigenvalue problem and thus controls the accuracy of the Bloch eigenvalues Enk→. Cursory investigations have been done with the value QMXSQ=20.16 a.u. corresponding to a maximal kinetic energy of the same value in Rydberg.The density of the electronsρ(r→) obtained from the occupied states is, as a consequence of the MAPW ansatz, a symmetric combination of plane waves everywhere, superimposed by an angular-dependent contributions that are treated in the accuracy wanted by use of a suitably chosen fine mesh of r→ points [16].The Brillouin zone (BZ) integrations over all occupied states yielding the electron densityρ(r→) and the total energy Etot are approximated by a sum over k→ points properly chosen in the irreducible wedge with the concept of Monkhorst and Pack [17] applied to the hcp structure. They are part of a hexagonal Bravais lattice spanned by the basis vectors (π/MG1a,π/MG2a,π/MG3c) which avoids locating the points H and K in the centre of the parallelepipeds. Here, the MGi,i=1,2,3 are specially chosen integers. Thus the influence of possible energy degeneracies is kept small. Within the self-consistent cycles the same value of QMXSQ is used which has been chosen to keep the fluctuation of the rank of the eigenvalue problems small.The characteristic ground-state properties are calculated within the local-density approximation (LDA) by use of the exchange-correlation functionals [18, 19]. As we found no significant difference we restrict ourselves in the following to results obtained with the functional [19].


## 3.  Results: Some Common Features of the Band Structure

Self-consistent calculations were performed for the lattice constanta=2.456Å=4.641167 Bohr [1] and the lattice constant c varies from 6.696Å=12.653605 Bohr [1] up to 36.0 Bohr.Figure1 shows the band structure along different special directions of the Brillouin zone characteristic for the experimental value of c. Up to c=36 Bohr the band structures look quite similar to previous investigations apart from a rigid shift caused by the Ewald potential and a global compression by 5.9% against the LCGTO-FF calculations [9] in the energy scale.Figure 1
Band structure along different lines in the Brillouin zone for the valuec=12.6536 Bohr. The notation of Bradley-Cracknell [15] is used. The broken line shows the Fermi level. In the outermost panel the twofold degenerated bands are marked red. All energies are in Ry.(1)
The low lying bands are well separated from the excited states. Both groups are only connected at the special points H and K which have almost coinciding energies. Table1 lists the values of these energies for different values of the lattice parameter c. With increasing value of c the energies EK and EH next to the Fermi energy get closer without any order in the sign. The low lying bands completely accommodate the 2s and 2p electrons of graphite for all values of c considered. The well-separated excited states are not occupied. The Fermi level is close to the eigenvalues of the H and K states quoted in Table 1. Its relative location to these values determines the semimetallic behaviour.
(2)
Parallel to the carbon planes, for example, along the linesΔ, U, and P, the dispersion of the bands weakens with increasing lattice constant c and is almost negligible at c≥24.0 Bohr. To make this more clear, the four bands closest to the Fermi level ϵF are displayed in Figure 2 for different values of c on an enlarged scale. Both outermost curves mark nondegenerate states, whereas the curve inbetween mark’s two-fold degenerate states. In previous investigations these bands were denoted by E1,E2, and E3, respectively. The twofold degeneracy is lifted as we move away from the zone edge; one of these bands will host the Fermi level. All energies are given as offsets to the energy values of the twofold states at the point K. Their scales have been chosen in such a way that their shapes become similar. It is worth mentioning that the scales vary by more than two order of magnitude when the interlayer distance is doubled. For the values of c considered the twofold generated band is cut by the Fermi level and thus separates occupied and non occupied states. Further details are discussed in the next section.
(3)
The metallicity becomes distinctly apparent in the density-of-states-curve, displayed in Figure4, by a fine dip near the energies of the twofold generated band (red curve) which sharpens with increasing value of the lattice constant c. Within the accuracy of the plotting our band structure fairly well agrees with that obtained with a special version of the Gaussian-type orbitals technique (LCGTO-FF) [9]. Especially, both calculations show a splitting at the bottom of the occupied bands at the point Γ in contrast to the FLMTO investigations [8].Table 1
The energy eigenvalues of the point K and the eigenvalues of the twofold degenerated state K in units of Rydberg with increasing interlayer distance. For each value ofc, by proper choice of QMXSQ, the set of the plane waves is generated which minimises the fluctuations of the number of plane waves NPW along the line KH.
c (Bohr)
NPW
E
K
NPW
E
H
-
E
K
12.6536
339
1.295390
336
0.002543
14.0
378
1.165958
378
0.000442
16.0
441
1.004455
438
0.000085
18.0
483
0.872945
483
0.000449
20.0
537
0.762965
540
-
0.000219
22.0
615
0.669224
612
0.000198
24.0
669
0.589686
666
-
0.000064
26.0
708
0.521203
708
0.000010
28.0
771
0.461227
768
0.000020
30.0
813
0.408885
816
0.000233
32.0
885
0.362247
882
0.000189
34.0
921
0.321406
918
0.000167
36.0
981
0.284721
978
0.000060Band structure along the line KH for the experimental lattice constantc=12.6536 Bohr, and the enlarged values c=16.0 Bohr, c=20.0 Bohr and c=24.0 Bohr. The middle, almost flat, curve shows the twofold degenerate band, whereas the other mark the nondegenerate bands. Horizontal axis in units of 2π/c. The notation of the relevant bands is in accord with previous work. Note the close agreement in the case of c=12.6536 Bohr with the result displayed in Figure 7(a) [20] or  Figure 2 [5] obtained by LDA and GW calculations.
(a)
(b)


## 4. Strategies to Deal with the Semimetallic Behaviour

With increasing value of the lattice constantc, the surface dividing occupied and nonoccupied states (Fermi surface) almost completely contracts to the line KH. In the energy region of interest the H point hosts two twofold degenerate eigenstates with energies differing by 0.034 mRy, whereas in the K point four eigenvalues cover the range of 0.160 mRy with a twofold degenerate state inbetween. Both values refer to c=24.0 Bohr with the tendency to shrink in the limit c→∞.Because of the confinement ink→ space all schemes to locate the Fermi surface, originally designed for ordinary metals, fail even if the number of sampling points of the wave vector is increased. On the other hand, it allows to use the k→·p→ approximation, which has already been applied within the initial investigations of graphite [12, 13], based on the four almost degenerate eigenstates on the KH line. For fixed value of k⊥, 0≤k⊥≤2π/c, we get an eigenvalue problem of rank 4
(3)ΣjHijxj=εxi,
where the leading contribution to the matrix elements Hij originates from the components of the momentum operator p→. Apart from k⊥ the elements Hij depend on the two-dimensional vector κ→ perpendicular to the vector characterising the KH line, for example, (2π/a3,2π/3a). As already mentioned, it is useful to define the energy eigenvalues ε relative to the energy of the degenerate eigenstate on the line κ→=0.Nontrivial solutions of the equation are defined by the characteristic equation(4)ε4-ε3A1+ε2A2-ε13A3+A4=0,
which defines four eigenvalues εi. The scalar magnitudes Ai are approximated by polynomials in the two-dimensional vector κ→. Symmetry operations of the small group with respect to the line KH (see Bradley-Cracknell [15, page 373]) leave their values invariant. Thus the k→·p→ approximation only allows even powers of the norm of the 2d vector κ→ with the important consequence that the Fermi surface remains unchanged under any rotation around the line KH. The evaluation of the matrix elements Hij in the MAPW scheme is straightforward but quite tedious. We can avoid it by the fact that the characteristic equation (4) may be rewritten by use of the eigenvalues εi(5)(ε-ε1)(ε-ε2)(ε-ε3)(ε-ε4)=0,
which conversely allows to express the scalars Ai by the four roots εi. Well known are the examples A1=ε1+ε2+ε3+ε4 and A4=ε1ε2ε3ε4. Analogous relations for A2 and A3 are easily derived.Thus we can proceed in the following way.(1)
For a fixed value ofk⊥ and a set of values |κ→| determine the roots εi by solving the MAPW problem and evaluate the scalars Ai by use of (4).(2)
In accord with the rotation invariance use the power expansions:(6)A1=α1+β1|κ→|2,A2≈α2+β2|κ→|2,A3=β3|κ→|4,A4≈β4|κ→|4.
for any value of  k⊥ and determine the coefficients αi and βi by a least squares fit. Since the energy eigenvalues have been chosen relative to the energy value of the twofold degenerate state on the line KH, the power expansion of A3 and A4 starts with |κ→|4. The omitted terms in A2 and A4 are by a factor |κ→|2 smaller than the included terms and are negligible for the modest values of c.With this background, the strategy to find ink→ space the points of equal energy is established. The characteristic equation (4) combined with the power expansions equations (6) is an implicit equation which determines the radius |κ→| of the circle connecting points of equal energy. As a consequence of the nonzero values β3 and β4 it is a quadratic equation defining either two different sheets or no sheet at all. Depending on the sign of the derivative d|κ→|2/dε, which can also be obtained from the characteristic equation, these circles run on a electron or hole sheet. Finally, the number of states Z+(ε) per unit volume with energies lower than ε is obtained by integrating over the area of the circles with positive sign of d|κ→|2/dε:
(7)Z+(ε)=±2(2π)3∫02π/c‍π|κ→|2dk⊥.This 1d integral is mostly restricted to subdomains where the integrand is defined. The analogous expressionZ-(ε) counts the hole states with energies greater than ε.The semimetallic behaviour of graphite requires that at least one electron and one hole-like sheet exist. In the case that the energyEK is larger than EH, the hole-like surface originates from the point K, whereas the electron-like surface develops from the point H. The exact compensation of the electron- and hole-like states defines the correct Fermi energy by the implicit equation
(8)Z+(εF)+Z-(εF)=0,
which is suitably solved near an estimated value ε0(9)εF=ε0-Z+(ϵ0)+Z-(ϵ0)N+(ϵ0)+N-(ϵ0).HereN±=dZ±/dε is nothing else than the density of electron and hole-like states per unit volume. This strategy produces the value the Fermi energy with the same accuracy as the one-particle energies Ek→.At this point the extraordinary small energy scale used in Figure2 needs a special consideration since we cannot expect that the MAPW scheme guarantees this high accuracy even when more than 600 plane waves are chosen according to (2). Thus it may happen that with slightly different values k⊥, distinctly different sets of plane waves are used with the consequence that the corresponding radii |κ→(ε,k⊥)| show distinctly visible jumps on the relevant length scale in k→ space. Such a jump may be seen in the case of c=20.0 Bohr near k⊥=0.4π/c. It allows to estimate the accuracy of the Bloch energies to be less than 0.25 mRy. Consequently the bands in the case c=24.0 Bohr show distinct fluctuations. To avoid them the bands displayed have been found by use of the same set of plane waves along the whole line.The result of this refined analysis for different values of the interlayer distancec is compiled in Table 2. It is remarkable that the number of electron or hole-like states Z± is decay by three orders of magnitude when the interlayer distance is doubled. The density of states at the Fermi level N± shows a weaker decay, but the decrease of the metallicity is still pronounced.Table 2
Fermi energy, volume of elementary cellVc, number of the states with energies below Z- and above the Fermi energy Z+,  and density of electron and hole states per unit volume at the Fermi level N±(ϵF). All magnitudes are expressed in corresponding atomic units.
c
V
c
ϵ
F
Z
±
(
ϵ
F
)
N
+
(
ϵ
F
)
N
-
(
ϵ
F
)
12.6536
236.0474
1.295935
±
0.6401
1
0
-
5
0.1464
1
0
-
1
0.1464
1
0
-
1
16.0
298.4730
1.00513582
±
0.1269
1
0
-
6
0.5041
1
0
-
2
0.3372
1
0
-
2
20.0
373.0912
0.76237693
±
0.1295
1
0
-
7
0.8873
1
0
-
3
0.8841
1
0
-
3
24.0
447.7094
0.59982141
±
0.1564
1
0
-
7
0.2324
1
0
-
3
0.1224
1
0
-
3Only forc=12.6536 Bohr a comparison with measurements is possible. Even our sophisticated treatment of the metallicity overestimates the number of conduction electrons and holes per C atom, which is defined by Vc|Z±(ϵF)|/4, by a factor of 4. McClure [13] found that his calculated carrier densities 2.4×10-5  per atom for electrons and 1.8×10-5  per atom for holes are in rough agreement with estimates derived from galvanomagnetic data.Measurements of the heat capacity [22, 23] give the value N(ϵF)=0.145 states per Ry and cell. This discrepancy can certainly be ascribed to many body effects beyond the LDA approximation which are completely neglected in the present work. According to [5], the electron-electron correlation produces the main corrections. It is suspected that in graphite the phonon enhancement is negligible.


## 5. Further Results

### 5.1. The Fermi Surface

According to the preceding section the Fermi surfaces are rotationally symmetric around the line KH even in the case of the experimental value ofc. This is in contradiction to a series of previous investigations [1, 2, 5]. With our choice of the parameter QMXSQ the electron-like surfaces are found to be located at the point K, whereas the compensating hole-like states are located around the point H in the case c=12.6536, 16.0 and 24.0 a.u. This is in accordance with previous investigations [24–26]. In the case c=20.0 a.u. the arrangement is reversed.Figure3 shows sections of a plane through the line KH with surfaces of equal energies magnified in radial direction that yield similarly shaped curves. The black lines describe electron-like states, and hole-like states are marked red. The Fermi surfaces around the plane KΓM have approximately the shape of extremely thin spindles. As a consequence of the existence of a mirror plane, the operation σ of the point group D6h, those around the plane HLA,  are double spindles with minimal calibres near the mirror plane. Values of the extremal cross sections are listed in Table 3. Again they are decayed by 3 orders of magnitude with the doubling of the interlayer distance. They clearly demonstrate how the Fermi surfaces shrink to the line HK with increasing value of c, with the consequence of extremely small de-Haas-van-Alphen frequencies. The cyclotron masses mc as well as the inverse transport masses mtrns-1 also listed in the table show a distinctly weaker decrease. The latter is defined by the sum in k→ space
(10)mtrns-1=13NSpur∑nk→‍∂ϵnk→∂kμ∂ϵnk→∂kν∂f∂ϵ(ϵnk→),
where N is the number of electron per cell. A highly simplified treatment of the dc conductivity based on the assumption of a relaxation time τ describing the scattering of electrons by lattice defects or phonons gives for an element of the dc tensor
(11)σμν=e2∑nk→‍τnk→∂ϵnk→∂kμ∂ϵnk→∂kν∂f∂ϵ(ϵnk→),
which constitutes, as consequence of (10), nothing else than a generalised Drude formula. Thus the transport mass allows to estimate the magnitude of the dc conductivity. It also determines the plasma frequency ωpl2=4πe2Nmtrns-1/ℏ.Table 3
Extremal areas, cyclotron masses, transport masses of electron (e) and hole  (h) Fermi surfaces. The areas arein units of (2π/a)2, and the masses are multiples of the free electron mass.
c
A
ext
,
e
m
c
,
e
mas
s
trns
,
e
A
ext
,
h
m
c
,
h
mas
s
trns
,
h
12.6536
0.7742
1
0
-
4
0.1195
1
0
-
1
0.7069
1
0
-
2
0.4467
1
0
-
3
-
0.1023
1
0
-
1
-
0.8727
1
0
-
2
16.0
0.1293
1
0
-
5
0.3965
1
0
-
3
0.3557
1
0
-
2
0.8211
1
0
-
4
-
0.2316
1
0
-
1
-
0.7092
1
0
-
2
20.0
0.1132
1
0
-
6
0.5056
1
0
-
3
0.3442
1
0
-
3
0.1771
1
0
-
6
-
0.9011
1
0
-
3
-
0.3442
1
0
-
3
24.0
0.2490
1
0
-
7
0.1491
1
0
-
3
0.1818
1
0
-
3
0.9565
1
0
-
8
-
0.7589
1
0
-
4
-
0.2519
1
0
-
3Intersections of planes through line KH with the surface of equal energies. Black line: electron states, red line: hole states. Dimensions in units of2π/a, 2π/c; respectively; however, the radii are enlarged by a factor of 10000. Rotation around the line KH generates the corresponding Fermi surface.
(a)
(b)Figure 4
The density of states (DOS) for graphite. The full line shows the calculated DOS, while the broken line is obtained by convoluting with a Gaussian ofFWHM=0.05 Ry−2. The dotted line shows the experimental X-ray photo emission spectra (XPS) [21].Measurements of the oscillatory magneto resistance (Shubnikov-de Haas Effect) of a high-quality graphite single crystal yield effective mass values in the basal plane of 0.039 free electron masses for electrons and 0.057 free electron masses for holes [27] which are significantly larger than theoretical values quoted in Table 3. The de-Haas-van-Alphen periods are 1.55 and 2.1×10-5 G-1, respectively, electrons at K and holes at H [28].


### 5.2. Density of States (DOS)

The electronic density of states as function of the energyN(E) has been calculated using the Gilat-Raubenheimer scheme [29] using 150k→ points in the Brillouin zone. For details see [30]. The typical step function in Figure 4 is due to the lack of dispersion along the line M–L. The pronounced peak at E≈0.85 Ry can be associated with the highest occupied state at the point L. In accordance with the previous discussion we find a deep notch at the Fermi level. In a crude manner we took into account the correlation enhancement by folding N(E) with a Gaussian of width 0.1 Ry. This broadening results in a curve is compared favourably with the spectra obtained by X-ray photo emission spectroscopy [21].


### 5.3. Charge Density

Figure5 shows the charge-density contours of the valence electrons in the fundamental triangle consisting of the two basis atoms in the lower corners and the centre of the hexagonal ring in the upper corner. The charge-density ρval assumes pronounced values near the two basis atoms; that is, ρval=4.9568 and 4.3869 electrons Bohr−3 at atoms of type A and type B, respectively, decreases rapidly and reaches minimal values near the line through the centre. Outside the circle which is the cut of the APW sphere with the plane: of the hexagonal ring, solely the plane wave part of the wave functions contributes to the charge density. Between the two basis atoms we find a secondary maximum with the value 0.28 electrons Bohr−3 which explains the covalent bonds between the different atoms of the hexagonal ring. As consequence of the symmetry D6h the contours around the centre of the hexagonal ring are almost spherically symmetric. These results qualitatively agree with previous investigations [8]. The charge density in the plane spanned by the vectors (a→1-a→2)/3 and a→3/2 is displayed in Figure 6. The quasi two-dimensional behaviour of graphite is obvious from the charge-density distribution in the plane containing the hexagonal ring, since most of the charge density is found in the planes, with very little overlap between the different planes. Again the secondary maximum of the density distribution on the line connecting the atoms of type A and type B is quite pronounced. By comparing Figure 7 which also shows the charge-density contours in the fundamental quadrangle for the enlarged c value 36.0 Bohr with Figure 6 we see that the overlap and interaction between the different planes decrease substantially with growing values of c. Thus the results derived for the biggest value of c well describe the charge-density distribution of graphene.Figure 5
Charge-density contours for hexagonal graphite at the experimental value of the lattice constants,a=4.6412 Bohr, c=12.6536 Bohr. The vector (a→1+2a→2)/3 points from the A-type atom, right-hand side corner to the B-type atom, left-hand side corner. The vector (a→1-a→2)/3 points from the A-type atom to the centre of the hexagonal ring. The charge-density contours ρval=0.025, 0.050, 0.075, 0.100, 0.125, 0.25, and 0.28 electrons Bohr−3 are indicated by the colours red, bright red, yellow, green, bright blue, blue, and black, respectively.Figure 6
Charge-density contours for hexagonal graphite at the experimental value of the lattice constants in the fundamental quadrangle spanned by the vectors(a→1+2a→2)/3 pointing in the horizontal direction and the vector a→3/2 pointing in the vertical direction. In both directions the same scale, in Bohr, is used. The median of the quadrangle is the cut with the hexagonal ring. The charge-density contours ρval=0.005, 0,025, 0.05, 0.1, 0.2, 0.3, 0.4, and 0.8 electrons Bohr−3 are indicated by the colours black, red, light-red, yellow, green, bright-blue, black, respectively. Note that the vertical midline is not a cut with the reflection plane.Figure 7
Charge-density contours for hexagonal graphite at the lattice constantsa=4.6411 Bohr and c=36.0 Bohr in the fundamental quadrangle. To show the weak dependence on c in the limit of large values the horizontal scale has been reduced by the factor 12.6526/36.0. For further details see Figure 6.


## 5.1. The Fermi Surface

According to the preceding section the Fermi surfaces are rotationally symmetric around the line KH even in the case of the experimental value ofc. This is in contradiction to a series of previous investigations [1, 2, 5]. With our choice of the parameter QMXSQ the electron-like surfaces are found to be located at the point K, whereas the compensating hole-like states are located around the point H in the case c=12.6536, 16.0 and 24.0 a.u. This is in accordance with previous investigations [24–26]. In the case c=20.0 a.u. the arrangement is reversed.Figure3 shows sections of a plane through the line KH with surfaces of equal energies magnified in radial direction that yield similarly shaped curves. The black lines describe electron-like states, and hole-like states are marked red. The Fermi surfaces around the plane KΓM have approximately the shape of extremely thin spindles. As a consequence of the existence of a mirror plane, the operation σ of the point group D6h, those around the plane HLA,  are double spindles with minimal calibres near the mirror plane. Values of the extremal cross sections are listed in Table 3. Again they are decayed by 3 orders of magnitude with the doubling of the interlayer distance. They clearly demonstrate how the Fermi surfaces shrink to the line HK with increasing value of c, with the consequence of extremely small de-Haas-van-Alphen frequencies. The cyclotron masses mc as well as the inverse transport masses mtrns-1 also listed in the table show a distinctly weaker decrease. The latter is defined by the sum in k→ space
(10)mtrns-1=13NSpur∑nk→‍∂ϵnk→∂kμ∂ϵnk→∂kν∂f∂ϵ(ϵnk→),
where N is the number of electron per cell. A highly simplified treatment of the dc conductivity based on the assumption of a relaxation time τ describing the scattering of electrons by lattice defects or phonons gives for an element of the dc tensor
(11)σμν=e2∑nk→‍τnk→∂ϵnk→∂kμ∂ϵnk→∂kν∂f∂ϵ(ϵnk→),
which constitutes, as consequence of (10), nothing else than a generalised Drude formula. Thus the transport mass allows to estimate the magnitude of the dc conductivity. It also determines the plasma frequency ωpl2=4πe2Nmtrns-1/ℏ.Table 3
Extremal areas, cyclotron masses, transport masses of electron (e) and hole  (h) Fermi surfaces. The areas arein units of (2π/a)2, and the masses are multiples of the free electron mass.
c
A
ext
,
e
m
c
,
e
mas
s
trns
,
e
A
ext
,
h
m
c
,
h
mas
s
trns
,
h
12.6536
0.7742
1
0
-
4
0.1195
1
0
-
1
0.7069
1
0
-
2
0.4467
1
0
-
3
-
0.1023
1
0
-
1
-
0.8727
1
0
-
2
16.0
0.1293
1
0
-
5
0.3965
1
0
-
3
0.3557
1
0
-
2
0.8211
1
0
-
4
-
0.2316
1
0
-
1
-
0.7092
1
0
-
2
20.0
0.1132
1
0
-
6
0.5056
1
0
-
3
0.3442
1
0
-
3
0.1771
1
0
-
6
-
0.9011
1
0
-
3
-
0.3442
1
0
-
3
24.0
0.2490
1
0
-
7
0.1491
1
0
-
3
0.1818
1
0
-
3
0.9565
1
0
-
8
-
0.7589
1
0
-
4
-
0.2519
1
0
-
3Intersections of planes through line KH with the surface of equal energies. Black line: electron states, red line: hole states. Dimensions in units of2π/a, 2π/c; respectively; however, the radii are enlarged by a factor of 10000. Rotation around the line KH generates the corresponding Fermi surface.
(a)
(b)Figure 4
The density of states (DOS) for graphite. The full line shows the calculated DOS, while the broken line is obtained by convoluting with a Gaussian ofFWHM=0.05 Ry−2. The dotted line shows the experimental X-ray photo emission spectra (XPS) [21].Measurements of the oscillatory magneto resistance (Shubnikov-de Haas Effect) of a high-quality graphite single crystal yield effective mass values in the basal plane of 0.039 free electron masses for electrons and 0.057 free electron masses for holes [27] which are significantly larger than theoretical values quoted in Table 3. The de-Haas-van-Alphen periods are 1.55 and 2.1×10-5 G-1, respectively, electrons at K and holes at H [28].


## 5.2. Density of States (DOS)

The electronic density of states as function of the energyN(E) has been calculated using the Gilat-Raubenheimer scheme [29] using 150k→ points in the Brillouin zone. For details see [30]. The typical step function in Figure 4 is due to the lack of dispersion along the line M–L. The pronounced peak at E≈0.85 Ry can be associated with the highest occupied state at the point L. In accordance with the previous discussion we find a deep notch at the Fermi level. In a crude manner we took into account the correlation enhancement by folding N(E) with a Gaussian of width 0.1 Ry. This broadening results in a curve is compared favourably with the spectra obtained by X-ray photo emission spectroscopy [21].


## 5.3. Charge Density

Figure5 shows the charge-density contours of the valence electrons in the fundamental triangle consisting of the two basis atoms in the lower corners and the centre of the hexagonal ring in the upper corner. The charge-density ρval assumes pronounced values near the two basis atoms; that is, ρval=4.9568 and 4.3869 electrons Bohr−3 at atoms of type A and type B, respectively, decreases rapidly and reaches minimal values near the line through the centre. Outside the circle which is the cut of the APW sphere with the plane: of the hexagonal ring, solely the plane wave part of the wave functions contributes to the charge density. Between the two basis atoms we find a secondary maximum with the value 0.28 electrons Bohr−3 which explains the covalent bonds between the different atoms of the hexagonal ring. As consequence of the symmetry D6h the contours around the centre of the hexagonal ring are almost spherically symmetric. These results qualitatively agree with previous investigations [8]. The charge density in the plane spanned by the vectors (a→1-a→2)/3 and a→3/2 is displayed in Figure 6. The quasi two-dimensional behaviour of graphite is obvious from the charge-density distribution in the plane containing the hexagonal ring, since most of the charge density is found in the planes, with very little overlap between the different planes. Again the secondary maximum of the density distribution on the line connecting the atoms of type A and type B is quite pronounced. By comparing Figure 7 which also shows the charge-density contours in the fundamental quadrangle for the enlarged c value 36.0 Bohr with Figure 6 we see that the overlap and interaction between the different planes decrease substantially with growing values of c. Thus the results derived for the biggest value of c well describe the charge-density distribution of graphene.Figure 5
Charge-density contours for hexagonal graphite at the experimental value of the lattice constants,a=4.6412 Bohr, c=12.6536 Bohr. The vector (a→1+2a→2)/3 points from the A-type atom, right-hand side corner to the B-type atom, left-hand side corner. The vector (a→1-a→2)/3 points from the A-type atom to the centre of the hexagonal ring. The charge-density contours ρval=0.025, 0.050, 0.075, 0.100, 0.125, 0.25, and 0.28 electrons Bohr−3 are indicated by the colours red, bright red, yellow, green, bright blue, blue, and black, respectively.Figure 6
Charge-density contours for hexagonal graphite at the experimental value of the lattice constants in the fundamental quadrangle spanned by the vectors(a→1+2a→2)/3 pointing in the horizontal direction and the vector a→3/2 pointing in the vertical direction. In both directions the same scale, in Bohr, is used. The median of the quadrangle is the cut with the hexagonal ring. The charge-density contours ρval=0.005, 0,025, 0.05, 0.1, 0.2, 0.3, 0.4, and 0.8 electrons Bohr−3 are indicated by the colours black, red, light-red, yellow, green, bright-blue, black, respectively. Note that the vertical midline is not a cut with the reflection plane.Figure 7
Charge-density contours for hexagonal graphite at the lattice constantsa=4.6411 Bohr and c=36.0 Bohr in the fundamental quadrangle. To show the weak dependence on c in the limit of large values the horizontal scale has been reduced by the factor 12.6526/36.0. For further details see Figure 6.


## 6. Conclusions and Outlook

With the careful precautions described perviously it is shown that the Fermi surface (FS) of crystalline graphite, even in the case of the measured interlayer distance, is essentially located around the line KH in the reciprocal space. It evolves from a twofold degenerate eigenstate on this line. In the light of thek→·p→ approximation and the small extension of the FS in the plane perpendicular to the line KH it is evident that the angular decomposition of the Bloch functions defined on the Fermi surface consists mainly of a term depending on the z coordinate which is especially located at the atoms of B type. Since approximately only half of the charge density is found to be within the APW spheres we have great doubts whether investigations based on the tight binding scheme or on the concept of localised Wannier functions give reliable results. The expansion of the interlayer distance c favours the semimetallic property of graphite accompanied with a further shrinking of the Fermi surface, for example, for c=36.0a Bohr the FS dimensions are reduced by an order of magnitude. This insight will considerably simplify the evaluation of other ground state properties based on the band structure near the Fermi surface. For example in (11) we can assume the relaxation time τnk→ on the Fermi surface to depend only on k⊥. As the total charge density strongly is decayed outside the atomic layers the results derived for the biggest value of c provide a good description for the case of graphene which is planned to be investigated in the near future.



---

*Source: 327913-2013-05-07.xml*