Graphene plasmons in the presence of adatoms

We theoretically investigate graphene plasmons in the presence of a low density of adatoms on the graphene surface. The adatoms can significantly modify the conductivity and plasmonic properties of graphene and may produce a level splitting with the plasmon mode, resulting in two plasmon branches. The high energy branch exhibits large losses and the low energy branch exhibits low losses. Our model may also be considered as a simple model for molecules on graphene and we show that graphene plasmons are sensitive to such changes in the environment. Our microscopic treatment of plasmons and adatoms shows the sensitivity of plasmons and highlights the potential of graphene plasmons for sensing purposes.

However, the quality of graphene still limits many proposed applications [21], and high quality graphene devices are labor-intensive to fabricate. Even the cleanest graphene samples exhibit some momentum relaxation [22], and thus a theoretical analysis of various loss mechanisms is of much interest [23][24][25][26][27][28]. The high frequency relaxation is of fundamental importance since the plasmons of interest are in this regime. Graphene conductivity in a relaxation time approximation [29,30] was investigated by Rana [31] and by Jablan et al [5] who found substantial plasmon losses for realistic relaxation times. This was also found in the experiments by Chen et al [9] and by Fei et al [10].
Graphene has previously been considered for sensing purposes, see for instance [32,33], and the large surface-to-volume ratio is one of the main advantages of graphene in this regard. Chemical sensing has been explored in the mid-infrared part of the spectrum where plasmons have been exploited to detect changes in refractive indices [34,35] and vibrational states in biomolecules [18,36]. These applications show much promise for plasmonic-based sensing in the future. However, also electronic transitions in atoms and molecules can couple to graphene plasmons but have so far not been analyzed in detail. Electronic properties of molecules and atoms that adsorb on the graphene surface have been studied extensively using various computational methods [37][38][39][40][41][42] and different substances have different coupling strengths and energy of the electronic levels. This variability between substances makes it possible to consider graphene plasmonic-based sensors which have the ability to selectively detect various compounds. We use these previous results to investigate adatom effects on the plasmonic properties of graphene. Our approach is complementary to previous works on plasmonic-based sensing in graphene as we study the plasmon response to microscopic degrees of freedom, rather than to changes in the dielectric environment [35] or to vibrational modes [18,36,43].
In this article we develop and analyze a model for uncorrelated adatoms, coupled to the graphene surface by tunneling (see figure 1). Since the adatoms are not the only imperfections in graphene we also include an electron relaxation time which we include in a number conserving manner following the Mermin prescription [5,29]. This relaxation time describes, phenomenologically, all sources of damping except the adatoms we are investigating. We explore the effects of adatoms on the single-particle properties of graphene as well as on the conductivity. We focus on the graphene surface plasmon mode and investigate its dispersion and the related losses. We find that plasmons close to resonance with the transition from the adatom energy level to above the Fermi energy become lossy. Furthermore, depending on the density of adatoms, their presence can split the plasmon mode into two separate branches, one low energy branch which experiences low losses and one high energy branch experiencing high degree of losses. We discuss how this can be used for ultra-sensitive sensing under the right conditions.
The article is organized as follows: in section 2 we treat the graphene plasmon dispersion and the graphene loss function. In section 3 we develop a manybody description of the system and derive an expression for the nonlocal longitudinal conductivity q, s w ( ). Finally, in section 4 we analyze the effects of the adatoms on plasmons, in particular on the plasmon dispersion relation and damping. In appendix A we present a derivation of two central equations in the article, equations (3) and (4). Appendix B gives details of the microscopic model and in appendix C a simplified expression for the susceptibility tensor is derived, from which the conductivity is obtained. Throughout the article we use c 1  = = .

Graphene plasmons
Longitudinal plasmons confined at a conducting interface between two dielectrics with relative dielectric constants 1 e and 2 e satisfy the dispersion relation [5,6] q q q , , This is a good approximation as long as the wavelength of the mode is much smaller than the free space wavelength, which is the regime we are investigating. In general equation (2) is a complex equation and for any given ω it can be solved by allowing complex wave vectors, q q q i 1 2 = + . Physically this means that the corresponding oscillation of the mode is damped. For weak damping, we can expand equation (2) in small q q 2 1 , and separating the real and imaginary parts of the Figure 1. Schematic view of the studied system. An infinite sheet of graphene is decorated by a dilute random distribution of adatoms, coupled to graphene by tunneling. The adatoms, here represented by the red spheres, are described by three phenomenological parameters: the energy level, 0  , measured with respect to the Dirac point, an intrinsic lifetime 1 d -, and a matrix element t describing the electron hopping between graphene and an individual adatom. The dielectric material on the two sides of the graphene have the relative permittivities 1 e and 2 e , respectively. The back gate V bg makes it possible to tune the Fermi energy in the graphene. Inset: zoom on the adatom with the hopping element represented with a red vertical bond.
) . The plasmon losses are given by equation (4) and in the second equality we have used equation (2) to express q, 2 s w ( )on the plasmon mode. Equations (3) and (4) show explicitly that in the low-loss limit the graphene plasmon dispersion is determined by q, ). The smaller the effective velocity, the higher the dissipative loss.
An alternative way to describe plasmons is to analyze the imaginary part of the current-current correlation function evaluated in the random-phase approximation. This is the spectral function of current fluctuations and describes where in q w space it is possible to deposit energy. For this reason it is sometimes called the loss function and it is defined as [30,44] S q v q q , 1 Im 1 , = e e e w + ( ) . The loss function has peaks where equation (2) is satisfied and these peaks are interpreted as signatures of the plasmons [30].
The plasmon dispersion, i.e., solutions to equations (3) and (4) as well as the loss function are discussed in section 4. First, the conductivity of graphene with adatoms needs to be calculated.

Microscopic model for graphene conductivity
The system considered consists of a pristine infinite graphene sheet [45,46] where momentum relaxation is added to model losses in graphene [5,29]. We introduce the convention that clean graphene means graphene with the finite momentum relaxation. To the clean graphene, we add a dilute distribution of identical adatoms that are modeled as non-magnetic Fano-Anderson localized states [30] coupled to graphene by tunneling as sketched in figure 1. When the coverage is dilute, correlations between adatoms are unimportant, and each adatom can be considered independently. Here, dilute means n 1 imp  , where n imp is the fraction of adatoms per lattice site. The total system can be described by a tight binding Hamiltonian that includes the tight binding Hamiltonian of graphene, the Hamiltonian of the adatoms, and the hopping between graphene and adatoms. Effective hopping Hamiltonians for different adatoms on graphene are obtained in [39][40][41][42], by use of DFT modeling of the composed system. In the following we assume a spin degenerate system, where the spin degree of freedom is included as a spin degeneracy factor g s = 2. An individual adatom is situated on the graphene atom at site x connected to a single graphene lattice site by a hopping parameter t. The adatom has a single energy level at 0  , measured relative to the charge neutrality point of clean graphene, and has an intrinsic lifetime 1 d -.
These are phenomenological parameters that are inputs of the model and in this section we explore the general features of the model for various adatom parameters. In section 4 we obtain parameters for hydrogen from [40] to examine a simple adatom. The Hamiltonian of the system with a single adatom is where dˆis the annihilation operator for the electron on the adatom and i ŷ is the annihilation operator of graphene electrons on site i. The sum in equation (6) is over all nearest neighbor sites of the graphene lattice and t 0 is the corresponding hopping parameter. The Hamiltonian is quadratic in the operators and can be diagonalized also for the case of a dilute density of adatoms, see appendix B for details.
In the regime of dilute density of adatoms and in the Dirac approximation, the Green's function of the Hamiltonian is given by [47] G g a p p , , 1 , a r c t a n , 7 where the renormalization factors Z p l  ( ) are the residues of the poles, i.e. Z G p p Res , figure 2. Due to the hybridization between the graphene bands and adatoms there is a level repulsion around 0  , close to which the hybridization is strong and gives rise to a finite lifetime of the order 1 dto the energy bands. Far from the level splitting, the energy states approach their uncoupled behaviors. The right panel of figure 2 shows the total density of states of the coupled system which exhibits a significant deviation close to the level splitting compared to the pristine graphene case. Specifically, there is a significant increase in the density of states close to 0  due to the coupling to the adatoms. The change of the bands compared to pristine graphene opens up new possible electronic transitions that alter the conductive and plasmonic properties of graphene, as discussed in the following.
The conductivity of the system can be computed using the Green's function in equation (7a). This is achieved by calculating the current response to an electric field ) (t is time in this paragraph, not to be confused with the coupling above). We restrict our analysis to the response to a longitudinal electric field which in the temporal gauge, t We set the electric field, E, and the momentum, q, to be parallel to the x-axis. According to minimal substitution [30], the perturbation given by the [45,46,48]. The diamagnetic current is zero in the Dirac approximation [48]. The conductivity of the system can be obtained from the current-current response function (longitudinal susceptibility) which relates the average value of the current to the vector potential to linear order j A q q , , [30,49]. From this expression and the relationship between the vector potential and the electric field, the conductivity can be seen to be q q , , The current-current response can be expressed in terms of the Green's function as where g g s v is the total spin-valley degeneracy factors, see appendix C for details. In pristine graphene, i0 , the energy integral in equation (8) can be performed analytically to arrive at the expression considered in [50][51][52][53]. In the temporal gauge the zero frequency component of the response is unphysical and has to be removed to avoid having a response to a static vector potential [53]. In the following we present results for electron doped graphene and when adatoms are present we take 0  to be positive (with respect to the Dirac point). However, our model can also be applied for hole doped graphene and/or negative 0  . Before analyzing the conductivity in the presence of the adatoms it is useful to examine the conductivity for pristine graphene. In the pristine case there exists a Pauli blocked triangle, inside which plasmon losses vanish at zero temperature due to the real part of the conductivity being identically zero [3,6,50,51]. This triangle in q, w ( )-space is shown in the inset of figures 3 and 4. For non-zero temperatures or when momentum relaxation is included, e.g. through a finite relaxation time, the triangle is no longer completely lossless but for moderate temperatures and relaxation times it is still the region where plasmons with low losses are expected to exist [5]. In this article we take the momentum relaxation time ( 1 G -) to be 1 ps, as reported in [54]. This relaxation time accounts for all intrinsic relaxation channels of the graphene and we include this in a number conserving way following the Mermin prescription [5,29]. As already introduced, we refer to graphene with the finite relaxation time as clean graphene to distinguish from graphene together with the adatoms (and also a relaxation time). The calculated conductivity in the presence of adatoms is shown in figure 3 for various values of the hopping parameter t and an impurity energy fixed at E 0.8 F 0  = , i.e., the energy level of the adatoms is close the Fermi energy. The presence of adatoms with energies close to the Fermi energy has an effect on the conductivity for frequencies close to the transition frequency E F 0  -| | between the adatom energy level and the Fermi energy. In particular there are peaks that appear inside the originally lossless triangle, which will give rise to larger plasmons losses, see equation (4). The imaginary part of the conductivity is also changed which will lead to changes in the plasmon dispersion as can be seen from equation (3). Figure 4 shows the conductivity for an increasing adatom energy level detuning from the Fermi energy and a fixed density of adatoms. As the energy level moves further from the Fermi energy the conductivity becomes more and more like the conductivity of clean graphene and the effect of the adatoms becomes negligible. From this we conclude that for the adatoms to have a large effect, the energy level of the adatoms needs to be close to the Fermi energy.
The new features that are present in the conductivity, as shown in figures 3 and 4, arise from the modification of the graphene bands caused by the presence of adatoms which is visible in figure 2. This modification of the conduction band around 0  changes the possible electronic transitions and in particular the allowed intraband transitions within the conduction band. These new transitions start playing a role around energies |, which is where there is enough energy for electrons to transition from the modified part of the band to unoccupied parts of the band above E F . This is the frequency around which the changes in conductivity start occurring. =˜is small for the densities we consider. Hydrogen is chosen since it is a simple atom and can serve as a typical atom adsorbed on the graphene. We emphasize that our model can be applied for other types of atoms and even simple molecules as well.

Graphene plasmons in the presence of a dilute density of adatoms
As stated previously, the loss function exhibits peaks where the dispersion equation has solutions. Figure (3), for the same density of adatoms as the loss function and the solid red line is the obtained loss, q 2 , from equation (4). The level splitting is accompanied by large plasmons losses and an emergence of two separate plasmon branches. The low energy plasmon branch exhibits low losses and the high energy branch has a large amount of accompanying loss. The larger loss in the upper branch can be understood by considering the Fermi golden rule, a new loss channel is opened for the plasmons in this branch. The loss channel is the excitation of a single electron around the adatom energy (of which there are many, see the DoS in figure 2) to above the Fermi energy where there are unoccupied electron states. The plasmons in the low energy branch do not have enough energy to lose energy through this channel. On resonance with this transition there is a very pronounced plasmons loss which separates the plasmon branches.
The right panel of figure 6 shows the evolution of a loss function cut at q k 0.013 F = as the density of adatoms is varied. This particular cut is chosen to show the loss function evolution for this particular adatom species (hydrogen) as clearly as possible. For small values of adatom density, there is only one plasmon peak visible, but as the density increases, the splitting into two branches is visible in the two emergent peaks. The separation between the two peaks grows as the adatom density is increased even further. To explore the adatom effect on plasmons in a large frequency range, figure 7 shows the plasmon propagation length along the dispersion until it crosses into the single-particle continuum. The plasmon propagation length (defined as the distance covered until the intensity of the plasmon falls by e 1 - [6]) in units of the plasmon wavelength, L p p l , can be obtained from the ratio q q 1 2 in equation (4) as L q q 4 Both panels of figure 7 show the plasmon propagation length as a function of frequency for different adatom densities on the graphene. The left panel is calculated for E 0.2 eV F = (k F  1/(3.3 nm)) and the right panel is calculated for E 0.4 eV F = (k F  1/ (1.65 nm)). The propagation lengths for clean graphene in both cases are shown in figure 7 as the magenta lines. The reason for the different propagation lengths in the two panels is that the relevant parameter for damping in the clean case is E F G . The plasmons are significantly affected by the presence of adatoms, in particular the damping is increased for energies above the transition frequency to excite electrons from the adatom energy level to above the Fermi energy. In the left panel of figure 7, this energy is roughly E 0.2 F , which corresponds to 0.16 eV 0  = (hydrogen), and E 0.2 eV F = . By changing the Fermi energy, the energy needed to make a transition to an unoccupied state changes. Therefore the energy at which the propagation length shows a step is different in the left and right panels of figure 7. For large enough densities, this step is where the plasmon is split into two branches.  The left and right panel of figure 7 are obtained for identical parameters except for the Fermi energy. The left panel clearly shows a larger effect for the same density of adatoms compared to the right panel. The conclusion is that for sensing purposes the Fermi energy should be tuned close to the adatom energy level for the sensitivity to be large.

Discussion
In this article we have focused on the graphene plasmon properties and how they are influenced by small densities of adatoms on the graphene surface. The induced plasmon losses and the level splitting that we find can be probed by for example light scattering in a grating environment [11], or on patterned graphene microribbons [1]. For a properly dimensioned grating or microribbon array, it would be possible to perform a laser frequency sweep and measure optical signatures of the presence of the adatoms as seen in figure 6. For the doping levels considered in this article, the typical dimension of the grating periodicity or the microribbon arrays needed are on the order of a few hundred nanometers. An alternative route to investigate the plasmons is by nanotip experiments such as in [9,10], where the number of plasmon oscillations are measured. Obtaining such data for different frequencies could reveal the presence of small amounts of adatoms and their energy levels.
Our analysis is restricted to the presence of a single kind of adatom, i.e., all the adatoms on the surface are characterized by the same 0  , t, and δ. In our model, the values of 0  , t, and the adatom density determine the characteristic level splitting that separates the two plasmon branches, see figures 6 and 7. Also, the sensitivity of the plasmons to the adatoms is found to be large when the adatom energy is close to the Fermi energy. Hence, by measuring the plasmonic properties and taking advantage of the tunability of E F offered by graphene, it is possible to determine 0  and thus discriminate between different adatoms. We stress that different adatoms exhibit different values of 0  as indicated by DFT calculations, see [37][38][39][40][41][42]55], and that our model is general enough to handle various adatoms and simple molecules. Our model thus enables selective sensing of various adatoms and molecules by probing the plasmonic properties of graphene. It should be noted that the adatom densities involved, 30 300 m 2 m --, is enough to increase plasmon losses and create a level splitting, making sensing of minute amounts of substances possible using graphene plasmons coupled to electron energy levels of adatoms. This is an increase of 2-3 orders of magnitude in sensitivity compared to experimental results for biomolecules obtained in [36], where plasmon coupling to vibrational modes of molecules was utilized for sensing. In our model, the adatom densities needed to produce a measurable plasmon response depends on the adatom coupling strength, t, which is considered as an input in our model. Thus, the sensitivity of the proposed sensing scheme is different for different adatoms. Adatoms that couple strongly to graphene will give rise to larger plasmon response for a given density than weakly coupled adatoms.
In this article we have taken a view towards sensing of the adatoms on the graphene surface. However, the adatoms may also be considered as imperfections on the graphene that impedes electron propagation by allowing the electrons to tunnel onto the adatom. For the purpose of plasmonics, long propagation lengths are often sought after and such damping is unwanted. We find that even small amounts of adatoms may have a significant effect on the plasmon damping. This is potentially one of several mechanisms that induces the large plasmon damping found in experiments [9, 10]. mrespectively. The sharp drop around E 1.3 F is due to the plasmon dispersion crossing into the particle-hole continuum where the plasmon becomes heavily damped.

Conclusions
We have investigated how graphene plasmons are affected by adatoms by comparing plasmons in realistic quality graphene with and without adatoms. We found that adatoms with energy levels close to the Fermi energy induces a strong level splitting between the bare plasmon mode and the adatom energy level. This level splitting is accompanied by large plasmon losses and depending on the adatom density may separate the plasmon mode into two separate branches, one low energy branch and one high energy branch. The low energy branch is virtually unaffected by the presence of the adatoms, whereas the high energy branch experiences larger losses. This is due to a new plasmon decay channel opening up, namely the excitation of an adatom electron to an unoccupied state above the Fermi energy.
Furthermore, we studied the sensitivity of the plasmon losses to the presence of adatoms. As a typical atom, we considered hydrogen and we found that a density of 300 adatoms per m 2 m is enough to give rise to a significant level splitting, and already 30 adatoms per m 2 m is enough to damp the upper branch. These effects could be measured in various light scattering experiments using dielectric gratings as well as using nanotips, making it possible to envision ultra-sensitive devices that measure the plasmon dispersion and losses to infer the presence of adatoms and molecules on the graphene surface.
Our results highlight the sensitivity of graphene plasmons to microscopic degrees of freedom and the possibility to use this effect in applications. Microscopic models for coupling various degrees of freedom to the plasmons is a very rich field and has the potential to further increase the already large sensing potential of graphene plasmons.

Acknowledgments
Financial support from the Knut and Alice Wallenberg foundation (KAW), the Swedish Research Council (VR), and the Swedish Foundation for Strategic Support (SSF) is gratefully acknowledged. We wish to thank Philippe Tassin, Daniel J Persson and Samuel L Avila for interesting discussions. G V and T W have contributed equally to this project.

Appendix A. Derivation of equations (3) and (4) in the main article
This appendix details the derivation of equations (3) and (4) starting from equation (2). Equation (2) defines the plasmon dispersion relation in the non-retarded limit. In the presence of losses, the dispersion relation is obtained by allowing complex wave vectors, q q q i The real part and zeroth order in q q

Appendix B. Fano-Anderson model in graphene
This appendix gives details on the microscopic model used in the main text. In particular the full Green's function in equation (7a) is derived.
We first consider pristine graphene, described in [45,46,48], coupled to adatoms by tunneling. The adatoms are modeled here as Fano-Anderson localized states [30] as described in the main text. To solve the system in the case of many adatoms, the following approximations are used: (i) all the impurities are identical, i.e., all of them are characterized by the same parameters , 0  d and t, (ii) the adatoms are uncorrelated and far apart so an average on position can describe the system. These assumptions are also used in the T-matrix formalism for weakly interacting electron systems in the presence of low densities of impurities [30].
In the main text the Hamiltonian of graphene and a single adatom is presented in equation (6). The Hamiltonian for many adatoms is [45,48]  where p is the momentum, t t e l x p p , , i l l = l˜· , l =  is the graphene band index and x l is the position of the lth impurity. Hence, for many impurities with vanishing hopping between them, the total Hamiltonian is obtained by adding single impurity contributions which is performed by the sum on l.
We introduce the notation P p, l = ( )and E E . Now, substitute these solutions into equations (B.4) and (B.5) we obtain o proceed further we take advantage of the assumption of identical impurities and we may thus perform the sums on l. We assume that the positions of the impurities are independent and randomly distributed, so the average on the impurity position is t t t l P l P N P P l l , ,

|˜| ||
. N imp is the number of impurities while n imp is the number of impurities per lattice site. This is standard procedure in T-matrix formalism for low density impurities [30]. ), in the temporal gauge since the adatoms carry no in-plane current so it is possible to avoid including terms related to the adatoms.
Here we compute the average value of the longitudinal current j ev q,

( )
. For E 0 = and in equilibrium, the values of currents and density deviation are vanishing. In Keldysh formalism, the current to linear order in the perturbation is written j r G q q , T , . As shown in [57] , , , , , C . 4 Here P c and T c are given by the poles of the Green's function and Fermi function respectively. We underline that to separate q, )in these two terms, the only assumption is that the self-energy is a smooth function of ò. Hence, it can be extended to include also other contributions to the electronic self-energy such as the phonon and electron-electron interaction. We underline that the susceptibility contributions in equations (C.6) and (C.7) are for one valley and one spin and thus needs to be multiplied by g g s v , as is done in the main text, to obtain the total result for graphene.
The momentum integrals in j j