Optical properties of two-dimensional honeycomb crystals graphene, silicene, germanene, and tinene from first principles

We compute the optical conductivity of 2D honeycomb crystals beyond the usual Dirac-cone approximation. The calculations are mainly based on the independent-quasiparticle approximation of the complex dielectric function for optical interband transitions. The full band structures are taken into account. In the case of silicene, the influence of excitonic effects is also studied. Special care is taken to derive converged spectra with respect to the number of k points in the Brillouin zone and the number of bands. In this way both the real and imaginary parts of the optical conductivity are correctly described for small and large frequencies. The results are applied to predict the optical properties reflection, transmission and absorption in a wide range of photon energies. They are discussed in the light of the available experimental data.


Introduction
Graphene, the flat two-dimensional (2D) honeycomb lattice of carbon atoms with sp 2 bonding, has attracted much attention due to its unique properties ( [1] and references therein). Because of its linear conical bands near the six corners of the hexagonal Brillouin zone (BZ), a constant dynamical conductivity and optical absorbance πα s , with α = 1 137.076 s the Sommerfeld fine structure constant, are observed experimentally in the infrared spectral region [1,2]. Besides graphene, the 2D honeycomb structure of silicon, silicene, has recently been discovered [3]. The stability of freestanding silicene and properties similar to those of graphene, despite the tendency for sp 3 bonding, have been confirmed theoretically [4,5]. The same holds for other group-IV materials within the 2D pristine geometry, germanene and tinene, or stanene [6]. The name of the pristine layer based on one monolayer of tin atoms is not unique. One may use the English name 'tin' or the Latin name 'stannum' and combine it with the suffix used in graphene [6,7]. Unique optical properties in the long-wavelength limit have been confirmed theoretically for graphene [8,9] but also predicted for silicene, germanene and tinene [10,11] if spin-orbit coupling (SOC) does not play an important role [6].
Chemical treatments by oxygen, hydrogen and halogen atoms may drastically change the electronic and accompanying properties of the 2D crystals with zero gap and Dirac-cone-like band structures. For instance, alternately hydrogenated and hence fully sp 3 -bonded graphene, so-called graphane, opens the fundamental gap at least up to 5.4 eV [12,13]. At the corresponding absorption edge strongly bound excitons occur. Strong excitonic effects near the quasiparticle absorption edge in the yellow or blue spectral region have also been predicted for hydrogenated silicene and germanene, silicane and germanane [13,14]. The symmetry break to open a substantial fundamental gap of a 2D honeycomb crystal also happens if the pristine geometry of graphene is conserved but each second C atom is replaced by a Si one. It results in 2D silicon carbide, which however possesses a large direct gap at the BZ corners due the charge transfer between Si and C atoms [15]. A further interesting group of 2D materials with honeycomb symmetry is isolated layers of van der Waals-bonded three-dimensional (3D) crystals, e.g. transition-metal dichalcogenides. A prototypical example is a monolayer molybdenum disulfide MoS 2 . It has also a fundamental direct gap at the BZ corners, as a result of inversion symmetry breaking in its honeycomb lattice structure [16,17]. Strong excitonic effects influence the optical properties below the fundamental gap region [17].
The treatment of the optical properties of isolated 2D honeycomb crystals but also in combination with other materials [1,2,18] requires particular attention, since for a macroscopic measurement technique such as the optical reflection, transmission or absorption spectroscopy the layer thickness d has to be studied in the limit → d 0. As a consequence, a macroscopic dielectric function of such a layered material cannot be defined, e.g. in a spectral ellipsometry technique. This is only possible by introducing a finite thickness such as d = 0.3 or 0.7 nm for graphene [19]. Another experimental approach is the assumption of a negligible reflection of a one-atomic-layer thick system such as graphene. Then, an absorbance of such a layer can be defined, which takes the universal value of 2.3 % in the spectral range from near-infrared to visible [1,2,8,10,11]. In order to avoid the thickness problem, frequently the optical conductivity σ ω ( ) 2D of a 2D system such as graphene is theoretically investigated [1,8,20,21] and also used to characterize experimental spectra [1,2,22]. The optical conductivity can be used to characterize the optical properties of any system containing 2D sheets [23]. The description of the optical conductivity σ ω ( ) 2D of electrons and holes in Dirac-cone-like band structures is relatively easy to derive [1,8,20]. However, the majority of first-principles calculations use a supercell or slab method to model the properties of single 2D sheets [9][10][11][12][13][14][15]. They allow the computation of the macroscopic complex, frequency-dependent dielectric function ϵ ω ( ) of the 3D periodic arrangement of almost one-atom-layer thick 2D crystals. The question arises how the calculated frequency dependences can be related to the optical properties of isolated sheet crystals or combinations of such 2D sheets with other materials.
In this paper, we show how such calculations of 3D arrangements can be related to the 2D optical conductivity of a single sheet crystal and how the result can be used to predict frequency-dependent macroscopic optical properties of systems containing such 2D crystals. As prototypical examples the undoped honeycomb crystals graphene, silicene, germanene, and tinene are studied.

Ab initio methodology
We employ total-energy and electronic-structure calculation methods which are based on the density functional theory (DFT) as implemented in the Vienna ab initio simulation package (VASP) [24]. We compute the atomic geometries of the four freestanding group-IV honeycomb crystals, in particular the lattice constant a and the buckling parameter Δ, within the generalized-gradient approximation of Perdew, Burke, and Ernzerhof for exchange and correlation (XC) [25]. The spin-orbit interaction (SOI) is generally taken into account in order to describe correctly the electronic states of heavy elements such as Ge and Sn. The geometry results have been published elsewhere [6].
Neglecting lattice vibrations the frequency dependence of the optical properties is mainly determined by the electronic structure. This is calculated including partially quasiparticle effects by describing XC within the screened hybrid density functional HSE06 according to Heyd, Scuseria, and Ernzerhof [26][27][28]. It yields to a reasonable zero-order approach to the quasiparticle band structures [29]. Pseudopotentials for C1s 2 , Si1s s p 10 cores and all-electron-like wave functions are generated within the projector-augmented wave (PAW) method [30]. The wave functions between the cores are expanded in plane waves up to a kinetic energy cutoff of 400 eV (graphene), 245 eV (silicene), 175 eV (germanene), and 103 eV (tinene) depending on the core radius.
The 2D crystals are simulated within a supercell approach with L = 20 Å of vacuum in the direction perpendicular to the 2D crystal surface in order to avoid artificial interactions between the periodic images of the 2D sheet crystals. The linear optical properties of this artificial 3D material are described by a frequency-dependent, complex dielectric function calculated in the framework of the independent-quasiparticle approximation [31]. Most convenient is the description of the optical transition matrix elements within the longitudinal gauge [32]. In this way, the influence of SOI on the optical oscillator strength and the inclusion of possible spinflip transitions are easily described by the corresponding dipole matrix elements with Pauli spinors.
The computation of the dielectric function ϵ ω ( ) contains sums over the Bloch states ν | 〉 k and hence the BZ of the artificial 3D crystal. Because of the large lattice constant L the BZ extent in the normal direction is extremely small and only one layer of k-points is sufficient. The k-point sampling is only necessary for the 2D hexagonal BZ of the corresponding honeycomb crystal. Near the six corner points K and ′ K of the BZ Dirac cones with linear bands and vanishing gap appear. Only the SOI opens a small gap of 0 (C), 2 (Si), 33 (Ge) or 101 (Sn) meV [6]. Therefore, for small frequencies the correct description of ϵ ω ( ) requires a fine sampling of the BZ around the K and ′ K points because of the ω −1 divergence of ϵ ω ( ) without SOI and hence the significant ω variation for small frequencies in the presence of SOC. In order to make the computations numerically tractable we apply a hybrid k-mesh strategy to sample the 2D BZ sufficiently. Around the K and ′ K points a much finer k-point mesh is used. In general, a 128 × 128 Γ-centered grid is used for the optical calculations, which however can be reduced to a 32 × 32 mesh to compute the screened exchange interaction in the HSE06 XC functional [28]. Two additional refined grids are applied, each including 869 symmetry-reduced k points within a circle of the radius π × a 0.05 2 and 0.005 × π a 2 , respectively, around one single K (or ′ K ) point. In contrast to the hydrogenated group-IV honeycomb crystals [13] the excitonic effects on the optical absorption of graphene have been shown to be small [9,14]. We verify this conclusion for silicene within the time-dependent (TD) DFT scheme for the inclusion of electron-hole interactions. We perform the calculations in the framework of the adiabatic extension of the local density approximation for the XC functional [33]. Since, as expected, this approximation fails in infinite systems, we also use the so-called long-range kernel α − q 2 in momentum space, which has the right asymptotic behavior of the exact TDDFT kernel [34,35]. Its validity to describe excitonic effects has been demonstrated in optical spectra of bulk materials, e.g. Si [34], but also of 2D systems, such as the Si(111)2 × 1 surface [36]. Even for such low-dimensional systems it has been demonstrated that the empirical kernel appropriately describes excitonic effects in optical absorption spectra.

Optical conductivity of 2D sheets
A 3D arrangement of the honeycomb sheets represents an artificial hexagonal crystal with two independent components of the dielectric tensor. Since we are interested in the optical properties of an isolated sheet with a vanishing (macroscopic) thickness → d 0, we study only the in-plane component ϵ ω ( ) of the dielectric tensor, i.e., a light polarization perpendicular to the sheet normal, the hexagonal axis. For a restriction to normal incidence of the polarized light optical properties are dominated by the in-plane optical conductivity σ ω ( ). However, also for arbitrary angles of incidence only the in-plane field-induced current modifies the boundary condition for the transverse components of the magnetic field vectors.
According to the Maxwell equation [36] it holds for a system without background dielectric constant. The quantity ε 0 is the permittivity of the vacuum. The 3D conductivity can be related to the 2D one σ ω ( ) 2D of an isolated sheet by changing from a 3D to a 2D sample. For the superlattice arrangement the relation is mediated by the lattice constant L as

2D
Within the independent-quasiparticle approach, for which the dielectric function ϵ ω ( ) can be directly related to the eigenstates ν | 〉 k and eigenvalues ε ν k ( ) via the Ehrenreich-Cohen formula [11], and without lifetime broadening effects, one finds for the real part of the 2D conductivity (using international units) with A as the sheet area and k ∈ of the 2D BZ. The interband electric-dipole matrix element is defined as [32]  for the four different group-IV honeycomb crystals [6,11] (and references therein).
The Kramers-Kronig relations for the dielectric function [37] are used to derive the imaginary part of the 2D optical conductivity from expression (3). By means of equations (1) and (2) with  as the Cauchy principal value. Expressions (1), (2), (3), and (5) indicate that there are two main ways to compute the real and imaginary part of the 2D optical conductivity. The first way is related to the direct computations of ϵ ω ( ) for the 3D sheet arrangement via the Ehrenreich-Cohen formula [11] and an application of equation (1) to obtain the conductivity σ ω ( ). The result for graphene is illustrated in figure 1 and compared with curves derived within a second method.
The direct calculation of ϵ ω ( ) with a homogeneous 128 × 128 k-point mesh using the VASP implementation yields a reasonable description of the dielectric function above photon energies of 1 eV. Below this energy the calculations are not converged because of too few k points employed to sample the K and ′ K regions of the 2D BZ. The second way begins with a direct calculation of σ ω Re ( ) via the Kramers-Kronig relation (5). Thereby, special care is taken to treat the singularity of the integrand and the reasonable quality of σ ω Re ( ) varying with the frequency as a consequence of the employed hybrid k mesh. First, we remove the pole replacing (5) by [38] This is possible since the integral over the principal value of the denominator vanishes. Second, since σ ω Re ( ) 2D has been computed only for a finite number of not equally distributed frequencies the integrand in (6) is interpolated using cubic splines [39]. This procedure allows for an easy performance of the integral in (6). The resulting conductivity σ ω ( ) is represented by smooth functions for both the imaginary and the real part and has the correct behavior for vanishing frequencies.
The influence of excitonic effects on the optical absorption is studied in figure 2 without taking quasiparticle shifts into account. In the framework of the non-local, energy-independent TDDFT kernel α ∼ − q 2 [34], we observe an almost vanishing influence of the direct screened electron-hole attraction. We have used α = 0.2 since this value works well for silicon surfaces. Test calculations with other values of α did not change the conclusions. Optical local-field effects which are described within the TDDFT by a Hartree response are more important. Anyway, since the peak positions are practically not affected and only a minor redistribution of spectral strength is visible, we neglect the excitonic effects in the following studies and take into account only quasiparticle shifts in the framework of the hybrid HSE approach.
The described way to perform converged computations has been applied also to the other three honeycomb crystals silicene, germanene, and tinene. The results for the real and imaginary part of the optical conductivity are displayed for the four group-IV 2D sheet crystals in figure 3.  Results of an independent-particle approximation using DFT-LDA (black line) are compared with those of the TDDFT including a nonlocal, energy-independent XC kernel −α q 2 with α = 0.2 (red line) [34]. The influence of optical local-field effects is illustrated by the green line. 0 has been discussed elsewhere [6,11]. Apart from variations due to the SOI-induced small gap the spectra tend to 1, while the curves for σ ω Im ( ) 2D vanish in the long-wavelength limit. For arbitrary frequencies, even in a wide energy range up to 20 eV (C) or 10 eV (Si, Ge, Sn) the most important intensity variations with the photon energy can be simply described by two (C, Si) or three (Ge, Sn) damped harmonic oscillators which may be traced back to some van Hove singularities in the 2D interband structure [6,11]. The low-energy oscillators have been discussed in the literature. For instance, the first peak in σ ω Re ( ) 2D of graphene near ω = 5.0 eV is a consequence of 2D saddle points in the π π * − interband structure located at the six (or better three counting their weight) M points at the BZ boundary. This peak has been also clearly identified in optical measurements with a position of the peak maximum at ω = 4.62 eV [40]. It should have a Fano lineshape which we, however, cannot confirm since no electron-hole coupling effects are taken into account. The first peaks near 2 eV in silicene, germanene and tinene have almost the same reason. The higher peaks near 14/16 eV (C), 5 eV (Si), 4.0/4.7 eV (Ge), and 3.1/4.0 eV (Sn) are related to σ σ → * transitions mainly at the Γ point of the 2D BZ. These peaks can basically also be related to van Hove singularities in the interband structure (not shown). For instance, for graphene the flat, i.e., dispersionless contributions of conduction and valence bands appear along the high-symmetry lines ΓK , ΓM but also along the zone boundary MK . The corresponding interband transitions are basically responsible for the shoulder and the peak visible in the high-energy spectrum of graphene.
Interestingly, the two main spectral features in the real part of the optical conductivity of pristine graphene agree well with the π and π σ + plasmon peaks found in experimental and theoretical energy loss spectra at loss energies of 4.7 eV and 14.6 eV [41,42]. These findings are roughly in agreement with the spectral variation of the in-plane optical conductivity in figure 3(a).
σ ω = Im ( ) 0 D 2 near ω = 5 eV and 16 eV is a consequence of ϵ ω ≈ Re ( ) 1 for the superlattice arrangement of graphene sheets. As a consequence it nearly holds ϵ ω ϵ ω ≈ − Im Im ( ) (1/ ( )). That means that the spectral behavior of the optical absorption and the energy loss function are rather similar.

Reflection, transmission and absorption
We study a single graphene, silicene, germanene, and tinene layer surrounded by two media = j 1, 2 which are characterized by dielectric functions ϵ ω ( ) j or by the corresponding complex index of refraction j j j j with its real part, the index of refraction ω n ( ) j , and negative imaginary part, the extinction coefficient κ ω ( ) j . The optical response of the group-IV layer is described by a normalized conductivity σ ω ε σ ω μ ϵ σ ω . For vanishing imaginary parts κ ω ≡ ( ) 0 j the situation is illustrated in figure 4. The layer normal is parallel to the z-axis. Optical isotropy is assumed in the xy-plane. In the three-layer system we distinguish between s-and p-polarization. We follow the description of the non-normal light propagation in [23,43]. With the convention Both quantities can be complex. According to [23,43] it holds for reflectance R, transmittance T and absorbance A σ σ =˜−˜−+˜+R n n n n , Thereby, the formulas (10) and (11) have been generalized for the case that the second medium 2 is absorbing. In the special case of two dielectrics with dielectric constants ϵ j expressions (10), (11) and (12) agree with the predictions in [23]. Many experimental situations can be studied by variation of the materials in the two halfspaces, the angle of incidence, and the light polarization. We give explicit results in figure 5 only for one of the simplest cases. We investigate the optical properties of a freestanding group-IV honeycomb crystal in air for normal incidence. It holds˜=˜= n n 1 1 2 and θ θ = = 0 1 2 . The expressions (10) These expressions show how the frequency-dependent conductivity influences the optical properties. According to the monatomic character of the sheet crystals the transmission is only weakly influenced by an isolated layer. This can be easily underlined in the long-wavelength limit ω → 0. Neglecting SOI the optical conductivity is known to be σ πα = (0) s [44]. It results

Summary and conclusions
In summary, ab initio DFT calculations with a hybrid exchange-correlation function and the supercell method are applied to calculate the 2D optical conductivity, real and imaginary parts, of an isolated single-atom-thick layer consisting of a group-IV honeycomb crystal. Special care is taken to compute converged spectra, in particular in the long-wavelength limit, with respect to the sampling of the Brillouin zone. Test calculations have been performed to study the possible influence of excitonic effects on the optical conductivity. Independent of the material its lineshape can be modeled within a picture of damped oscillators. Important oscillators can be related to interband transitions, for instance, between π and π * bands at the M points and σ σ → * transitions near the Γ point. Their peak positions are not too far from those of the π and π σ + plasmons in graphene. The influence of an infinitely thin 2D crystal on the optical reflectance, transmittance, and absorbance is modeled by the normalized 2D conductivity σ ω ε c ( ) 2D 0 . For freestanding honeycomb crystals the reflectance is vanishingly small outside the resonance frequencies. As a consequence ω σ ω ε = A Re c ( ) ( ) 2D 0 holds. Despite the small layer thickness the absorption of such a 2D crystal is not vanishing. Rather, near the optical resonances the absorbance ω A ( ) approaches large values, one order of magnitude larger than the approximate zero-frequency value πα = A (0) s . We suggest future experiments on how single or multiple layers of group-IV honeycomb crystals are suitable to tailor the optical properties of a given substrate in a wide energy range.