{'fraction_non_alphanumeric': 0.04160562554937006, 'fraction_numerical': 0.014259205000488329, 'mean_word_length': 4.422418358340688, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Quantization of graphene plasmons

Beatriz A. Ferreira, B. Amorim, A. J. Chaves, and N. M. R. Peres1,2∗ Department of Physics, Center of Physics, and QuantaLab, University of Minho, Campus of Gualtar, 4710-057, Braga, Portugal International Iberian Nanotechnology Laboratory (INL), Av. Mestre José Veiga, 4715-330 Braga, Portugal CeFEMA, Instituto Superior Técnico, Universidade de Lisboa, Av. Rovisco Pais, 1049-001 Lisboa, Portugal and Department of Physics, Instituto Tecnológico de Aeronáutica, DCTA, 12228-900 São José dos Campos, Brazil


I. INTRODUCTION
In many cases, light-matter interaction can be understood in a semi-classical picture, where matter is quantized and the electromagnetic field (EM) is treated classically. This semi-classical approach holds when the number of photons in the field is large or the light source is coherent. On the other hand, in order to understand the properties of a small number of photons the quantization of the EM field is required. Typical phenomena where the quantization of the EM field is necessary involve entanglement, squeezed light, cavity electrodynamics, interaction of photons with nano-mechanical resonators, and near-field radiative effects [1].
In near-field radiative effects, plasmonics emerges as a promising candidate to observe quantum effects of the electromagnetic radiation, an example being the Hong-Ou-Mandel interference of plasmons [2]. Many other quantum effects in plasmonics exist, such as the survival of entanglement, particle-wave duality, quantum size effects due to reduced dimensions of metallic nanostructures, quantum tunneling of plasmons (which are simultaneously light and matter), and coupling of quantum emitters to surface plasmons [3][4][5][6][7][8][9][10][11].
The exploration of quantum effects in plasmonics in unusual spectral ranges, such as the THz and the mid-IR, has been deterrent by the poor plasmonic response of noble metals in these regions of the electromagnetic spectrum. However, with the emergence of graphene plasmonics [12,13] the possibility of exploring quantum Figure 1. Schematic representation of the three systems considered in this paper: a single graphene layer (left), and graphene double layer (center), and a graphene sheet near a perfect metal. The quantities n refer to the dielectric permittivity of the medium n and σg refers to the optical conductivity of graphene. effects in these yet unexplored spectral regions is a possible prospect. Despite the fact that, at the time of writing, quantum effects involving graphene plasmons remain illusive, the fact that graphene plasmons are characterized by low losses [14][15][16] boosts the above hope. In addition, graphene plasmons screened by a nearby metal (also called screened or acoustic plasmons) can be confined down to the atomic limit [17], which certainly opens the prospects of finding rich grounds for quantum plasmonics and nonlocal effects [18,19]. Indeed, the idea of developing quantum optics with plasmons has already a long history [20] and quantization of localized plasmons in metallic nanoparticles was recently performed [21,22].
The development of quantum theory of the electromagnetic field in the presence of dieletric media has a long history and several approaches have been developed [23][24][25][26][27][28][29][30][31][32][33][34][35]. These methods are typically based either on the quantization of the macroscopic classic energy [27], on the classical dyadic Green's function for the electric field [28] or on the diagonalization of Hopfield-type Hamiltonians [26]. The quantization of evanescent EM waves [36,37] and of the EM field in the vicinity of a metal [38] have also been considered in the past.
In this paper, we perform the quantization of graphene plasmons, obtaining the plasmonic electromagnetic field mode-functions and, importantly, their normalization, when losses are neglected. These mode-functions are then used to study the interaction of graphene plasmons with nearby quantum emitters and determine, in a very intuitive way, using Fermi's golden rule the spontaneous decay rate of the emitter due to plasmon emission. Thee quantization of graphene plasmons is performed in two ways: (i) A macroscopic approach, which starts from the classical time-averaged energy of the electromagnetic field in a dielectric [24,27,39,40]. This method allows for the inclusion of dispersion in the optical response of graphene. (ii) A hydrodynamic approach, where graphene charge carriers are described in terms of an electronic fluid, which couples to the eletromagnetic field [41,42]. The hydrodynamic approach allows not only for the inclusion of dispersion, but also for the inclusion of non-local effects in the optical response of graphene.
The paper is organized as follows: in Section II, we present the general macroscopic approach for the quantization of the eletromagnetic field in dispersive, lossless media and a normalization condition for the modefunctions is determined. In Section III, we present the quantum hydrodynamic model or graphene. We will see that when non-local effects are neglected, the result of the macroscopic approach is recovered. In Section IV, the plasmon dispersion relations, mode-functions and modefunction normalization for a single graphene layer and for a graphene-metal structures are obtained. In Section V, we use the quantized plasmon fields to compute the decay rate of a quantum emitter due to the spontaneous emission of plasmons. The plasmon emission rate is compared with the total decay rate of the emitter, which is obtained from the complete dyadic Green's function for the electric field. The role of non-local response of graphene is analyzed. Finally, we conclude with Section VI, commenting the obtained results and discussing future research paths. A set of appendices detailing the calculations is also presented.

II. MACROSCOPIC QUANTIZATION OF THE PLASMON ELECTROMAGNETIC FIELD
In this section, we will describe how the plasmon fields can be quantized using the macroscopic classical energy of the electromagnetic field in a dispersive, lossless dielectric medium, a method first used in Refs. [24,27,40]. For electric and magnetic fields with a harmonic time dependence, (1) close to a central frequency ω, the time-averaged classical electromagnetic energy in the presence of a dispersive, lossless dieletric is given by [39,43] where¯ r (r, ω) is the relative dielectric tensor of the medium, 0 and µ 0 are, respectively, vacuum permittivity and permeability. The idea of this method is to take the above equation as the quantum mechanical energy of a EM field eigenmode with frequency ω. We will work in the Weyl gauge, in which the scalar potential is set to zero φ = 0, such that the electric and magnetic fields are obtained only from the vector potential A as The vector potential is then expanded in modes as where α λ are amplitudes for the mode λ, with modefunction A λ (r) and corresponding frequency ω λ . The mode-functions and frequencies are determined by solving the non-linear eigenvalue problem with c vacuum's speed of light, which is just Ampère's law in the dielectric medium for the mode-function A λ (r). Next, we assume that the total time-averaged energy for the vector potential Eq. (6) is given by The quantization of the theory is achieved by promoting the amplitudes α λ to quantum mechanical operators whereâ λ (â † λ ) are bosonic annihilation (creation) operators, which obey the usual equal-time commutation relations and L λ is a mode-length which is determined by demanding that the quantum mechanical Hamiltonian which is obtained from Eq. (8) by performing the replacement of Eqs. (9) and (10), coincides with the Hamiltonian for a collection of quantum harmonic oscillatorŝ Imposing this condition, we have that the mode-length is given by L λ = U (ω λ )/ 2 0 ω 2 λ . Using the mode-function equation (7) to write The mode-length can be written as Notice that in the absence of dispersion, the second term vanishes and L λ reduces to the usual norm in the presence of a position dependent dielectric tensor¯ r (r, ω λ ).
Although this phenomenological method appears to be unjustified, it has actually been shown to be correct for the case in when the dielectric is modeled by a Lorentz oscillator [40]. We will see in the next section, that Eq. 15 remains valid within a quantum hydrodynamic model of graphene, even when non-local effects are included in the optical response. As a matter of fact Eq. (15) for the mode-length remains valid for any linear optical medium (including effects of dispersion, non-locality, inhomogeneity and anisotropy) as long as losses can be neglected [44].

III. PLASMON QUANTIZATION WITHIN A QUANTUM HYDRODYNAMIC MODEL
In this section, we will perform the quantization of graphene surface plasmons employing a hydrodynamic model. The hydrodynamic model treats both the electron gas and the electromagnetic field using a classical picture and provides a simple and elegant way of including non-local effects [45]. Non-local effects are taken into account by including a pressure term in the Boltzmann equation, that arises due to Pauli's exclusion principle. A detailed derivation of the hydrodynamic model for graphene can be found in [41,42]. To illustrate the method, we choose the simple case of a single graphene sheet located at the plane z = 0, embedded by a static dielectric medium with relative dielectric constant¯ d (r).

A. Classical hydrodynamic Lagrangean and Hamiltonian
The classical Lagrangean density for the hydrodynamic model of graphene can be written as: where L EM is the Lagrangean density of the electromagnetic field and L 2D hyd describes the electronic fluid of graphene and its coupling to the electromagnetic field. Using once again the Weyl gauge, L EM is given by where¯ d (r) is allowed to be position dependent, but is frequency independent. Within the hydrodynamic model, the electronic fluid of the graphene layer is described by the fluctuation, n, of the density around the equilibrium density, n 0 , and the displacement vector υ = (υ x , υ y , 0), which should not be confused with the velocity field. In the Weyl gauge, L 2D fluid is written as [41,42] where the δ-function δ(z) restricts the electronic fluid to the z = 0 2D plane, m is the Drude mass and β appears from the relation between the degeneracy pressure and the electronic density and depends on the band dispersion for the carriers (see Ref. [41]). In the approximation of the linear dispersion for graphene electrons, the hydrodynamic parameters are given by [41]: where k F is the Fermi wavenumber and v F the Fermi velocity of graphene. Equation (18) is the 2D equivalent Lagrangian for the hydrodynamic model presented in [46].
Using the Euler-Lagrange equations for Eq. (18) with respect to A, we obtain which is nothing more than Ampère's law in the presence of a surface current given by Using the Euler-Lagrange equations for Eq. (18) with respect to the fluid variables n and υ, we obtain the continuity and Newton's second law with a diffusion term, which read respectively from which we recognize the fluid electronic surface density Equations (19)- (22) correspond to the linearized hydrodynamic model for graphene [41] (see also [47]). Notice that Eq. (21) has no dynamics. Therefore, we can use it to eliminate the field n, thus obtaining a new Lagrangean density with which can be obtained by combining Eqs. (21) and (22). From the Lagrangean density Eq. (24), we define the canonical momenta conjugate to A and υ, respectively, as In terms of the variables A, Π, υ and π, the classical Hamiltonian obtained from Eq. (24) is given by

B. Canonical quantization of hydrodynamic model
In order to quantize the classical Hamiltonian Eq. (29), we start by introducing the eigenmodes of the coupled electronic fluid + electromagnetic field. Assuming, we have in-plane translational invariance, we write the vector potential and fluid displacement variables as where S is the area of the graphene layer, α q,λ (t) = α q,λ e −iω q,λ t are mode amplitudes with, ω q,λ the mode frequency, and A q,λ (z) and υ q,λ are mode-functions, which, following from Eqs. (19) and (26), obey the equations where we defined the differential operator D q = iq +ẑ∂ z . From Eq. (33), we can write whereσ hyd g (q, ω) is the conductivity within the hydrodynamic model, which we separate into transverse and longitudinal components as respectively given by where we identify D = e 2 n 0 /m as the Drude weight, which for graphene is given by D = e 2 v F k F /π. Notice that in the the limit q → 0,σ hyd g (q, ω) reduces to the Drude model. Inserting Eq. (34) into Eq. (32), we obtain with the dieletric function, including screening by graphene electrons, being given bȳ in agreement with Eq. (7). Inserting the expansions Eq. (30) and (31) into Eq. (29), and using the orthogonality properties of the mode-functions (A q,λ (z), υ q,λ ) it is possible to write the Hamiltonian for the hydrodynamic model as (see Appendix (B)) with the mode-length defined here as Promoting the mode amplitudes to quantum mechanical operators as whereâ † q,λ (â q,λ ) are creation (annihilation) operators for the plasmon-polaritons, satisfying the usual equal-time commutation relations Eq. (11), the quantum Hamiltonian for the hydrodynamic model becomeŝ We will now see that Eq. (41) can be recast in the same form as Eq. (41). In order to do so, we use Eq. (34) which allows to write Eq. (41) as It is easy to see that the above equation can also be written as with¯ r (q, z, ω) given by Eq. (39).

IV. DISPERSION RELATIONS AND MODE-FUNCTIONS OF GRAPHENE PLASMON IN TWO DIFFERENT STRUCTURES
We now wish to determine the dispersion relation and mode-functions of graphene plasmons in two different structures (see Fig. 1): a single graphene layer and a graphene sheet in the vicinity of a perfect metal. As in the previous section, we make use of the in-plane translational invariance of the structures being considered. Therefore, the non-linear eigenvalue problem for the mode-functions, Eq. 7, can be written as wherē (48) where¯ d (z) is the dielectric constant of the medium, which we assume to be isotropic and a piecewise homogeneous function of z, and labels the graphene layers which are located at the planes z = z , with conductivityσ g (q, ω). We model the conductivity of each graphene layer with Eq. (35), which when including losses becomes where D = 4E F σ 0 / (π ) is the Drude weight, with E F graphene's Fermi energy, σ 0 = πe 2 / (2h) and γ is a relaxation rate. When determining mode-functions we will set γ = 0, but we allow for γ = 0, for the situation when ohmic losses are included in Section (V). The presence of the δ-functions in Eq. (48), implies that boundary conditions at each interface located at z = z : where E q,λ (z) = iω q, A q,λ (z) and B q,λ (z) = D q × A q,λ (z) are the electric and magnetic fields corresponding to mode A q,λ (z), and J s,q,λ (z ) =σ g (q, ω)·E q,λ (z ) is the surface current in the graphene layer . In addition to the boundary conditions Eqs. (51) and (52), to determine of the plasmon modes one must impose that the fields decay for z → ±∞. Having determined the plasmon mode-function, A q,sp (z), and dispersion, ω q,sp , the mode-length can be obtained from Eqs. (46) and (48) as We will now analyse the different structures in detail.

A. Single layer graphene
We first discuss the simplest case of a single graphene sheet (see left panel of Fig. 1). The problem of finding the spectrum of the surface plasmons in a graphene sheet was first considered in [48] and a detailed derivation can be found in Refs. [49,50]. We assume that the single layer of graphene is located at z = 0, with a encapsulating dielectric medium n = 2 for z > 0 and a medium n = 1 for z < 0, such that In order to determine the plasmon mode, we look for ppolarized solutions of the electric field (the electric field lies in the plane of incidence) in the form of evanescent waves for z → ±∞. In each piecewise homogeneous region we have that D q · E q,λ (z) = 0. The mode-function must then take the form where with c n = c/ √ n the speed of light in medium n, and we introduced the vectors Imposing the boundary conditions Eqs. (51) and (52), we obtain the following implicit relation for the surface plasmon dispersion: In general, Eq. (58) has no analytical solution, except in the simple case where 1 = 2 , in which case its solution reduces to solving a quadratic equation. The corresponding mode-function can be written as and the mode-length is obtained according to Eq. (53) as where the last term is due to the dispersion in the graphene layer. We point out, that within the hydrodynamic model used for conductivity of graphene, this contribution is only non-zero if non-local effects are also included, that is if β = 0.

B. Graphene-metal structure
We now move to the more complex case of a graphene sheet near a perfect metal (see right panel of Fig. 1).
We assume that the metal interface is located at z = −d and the graphene layer is located at z = 0. The dieletric constant is given by Noticing that the plasmon field should decay for z → ∞, the mode-function should have the form (62) Notice that the presence of a perfect metal at z = −d implies that the tangential component of the electric field should vanish. Imposing the previous condition and the boundary conditions Eqs. (51) and (52) at z = 0, we obtain a homogeneous system of equations where . The dispersion relation of the screened plasmons is obtained by looking for the zero the determinant of the previous matrix, which leads to the condition for the dispersion relation The boundary conditions imply that the mode-function is given by , z > 0 i q q sinh (κ 1 (z + d)) + + q κ1,qẑ cosh (κ 1 (z + d)) , 0 > z > −d The mode-length can be determined from Eq. 53, and we obtain where, as in Eq. (60), the last term is due to the dispersion of graphene. We note that the dispersion relation for the screened plasmons, Eq. (64), is the same one that would be obtained for the acoustic plasmons in a symmetric graphene double layer structure (center panel of Fig. 1). In this structure, we have two graphene layers located at z = 0 and z = −d gdl . The dieletric constant of the encapsulating medium is given by Since the plasmonic modes should decay for z → ±∞, the plasmon mode-function must have the form (67) Imposing the boundary conditions Eqs. (51) and (52) at z = 0 and z = d, we obtain the following homogeneous system of equations . The dispersion relation is obtained from zeroing the determinant of the previous matrix. Since the system is composed of two graphene sheets the double layer structure has two dispersion branches, a low energy one -the acoustic mode -and a high energy one -the optical mode. In the optical mode, the charge oscillations in the two graphene sheets are in-phase; whereas in the acoustic mode, the charge oscillations are out-of-phase. In the particular case where, we have a symmetry structure, with 2 = 3 and the conductivities of the top and bottom graphene layers are the same σ gtop,L = σ gbot,L = σ g,L , the zeroing of the determinant factorizes into two independent expressions for the optical mode, and for the acoustic one. Notice that the relation for the acoustic mode dispersion Eq. (70) coincides with the equation for the screened plasmon Eq. (64) provided d dlg = 2d. This fact can be understood in terms of image charges as depicted in Fig. 2.
In the bottom panel of Fig. 3 we depict the loss function for the graphene-metal system. Clearly, only one branch is seen, which coincides with the acoustic branch of the double layer graphene upon considering d equal to half that of the double layer system, as explained in Fig. 2.
An alternative way to obtain the plasmon dispersion relation is to look for poles (or resonances in the presence of losses) in the so called loss function [50] , which is defined as where r p (q, ω) is the reflection coefficient of the structure in consideration for the p-polarization, and q and ω are Figure 2. Comparison of the double layer graphene system and the graphene-metal case. Due to the image charges in the metal-graphene structure, the spectrum of the double layer graphene is equivalent to that of the graphene-metal system if we take into account that graphene is at at distance d gdl /2 from the metal, where d gdl is the interlayer distance in the double layer case. Therefore, the graphene-metal distance has to be d/2 to obtain the same spectrum as for the graphene double layer . For a full numerical equivalence of the two spectra it is also necessary to have 1 = 3 in the double layer geometry.
the in-plane wavevector and frequency of the impingent radiation, respectively (see Appendix A). For a symmetric graphene double layer ( 1 = 3 and σ gtop,L = σ gbot,L ) and neglecting losses γ → 0, this coefficient has poles at the solutions of Eqs. (70) and (69), as can be seen comparing those equations with Eq. (A5). The loss function for the double layer graphene is depicted in the left panel of Fig. 3, as function of ω and of a dimensionless parameter s = qc/ω which defines the dispersion relation of the single graphene layer, clad by two different dielectrics of dielectric functions 1 and 2 , in the electrostatic limit by where E F is the Fermi energy of graphene and α is the fine structure constant of vacuum. In the top panel of Fig. 3 two branches are clearly seen: a high energy onethe optical branch -and the acoustic branch at lower energies. At high energies and high s the two branches merge and converge to the single layer branch. The reader may wonder why the lower branch starts at finite momentum. This happens due to the definition of the s parameter, which involves both the frequency and the real wavenumber q. This choice allows to clearly separate the two branches in the (ω, s) plane.

V. APPLICATION: QUANTUM EMISSION CLOSE TO GRAPHENE STRUCTURES
We will now apply the quantization of the plasmon modes in grahene structures to the problem of spontaneous emission by a quantum emitter which is located above the structure. We model the quantum emitter as a two-level system embedded in medium 2 at position r 0 = (0, 0, z 0 ). The quantum emitter couples to the plasmonic field via dipolar coupling:Ĥ sp-d = −d ·Ê sp (r 0 ), whereĉ † g/e (ĉ g/e ) is the creation (annihilation) operator for the ground/excited state of the two-level system and d ge = g d e is the dipole matrix element.
The transition rate of an emitter due to emission of surface plasmons in graphene is given by Fermi's golden rule [51,52]: where ω 0 is the transition frequency, |g; n q,sp + 1 represents a final state with one more surface plasmon and the emitter in the ground state and |e; n q,sp represents an initial state with n q,sp plasmons in graphene and the emitter in the excited state. The transition matrix element for spontaneous plasmon emission, when there are no surface plasmons in the initial state,is given by With this result the transition rate reads Using the mode-function we can write ψ is the angle the dipole makes with the axis perpendicular to graphene (z-axis), and φ is the azimuthal angle. The prefactor N q is define as N q = 1 for the singlelayer graphene case and as N q = sinh 2 (κ 1,q d) for the graphene+metal structure. Using the in-plane isotropy of the system, the momentum integration in Eq. (77) can be trivially performed, yielding: where q 0 is the momentum of a surface plasmon, with frequency ω 0 , i.e., ω 0 = ω q0,sp . So far the expression for Γ sp is general. The differences arise from the particular forms of the dispersion ω q,sp , the mode length L q,sp and the prefactor N q .

A. Decay rate for local conductivity
We will now study the plasmon emission rate, when non-local effects are neglected, β = 0 in Eq. (50).
We first focus on the case of a quantum emitter close to a single graphene layer and the same dielectric above and bellow the graphene layer, 1 = 2 = . Using the analytic solution of Eq. (58) for this case, we can write Eq. (79) as In Figs. 4 and5, we plot the ratio Γ gsl sp /Γ 0 , where Γ 0 = d 2 ge ω 3 0 / 3π 0 c 3 the total decay rate of an emitter in vacuum, for different graphene-quantum emitter distances z 0 and dipole orientations. For comparison, we will display the ratio Γ full /Γ 0 , where Γ full is the total decay rate (including both plasmon and photonic losses) of the quantum emitter. The total decay rate can be computed from the knowledge of the reflection coefficients, which are incorporated into the dyadic EM Green's function [50,51,53] (see Ref. [54] for a detailed study of the properties of an emitter near graphene using Dyadic Green's functions). For a dipole in medium n = 2, the total decay rate is given by where r s/p are reflection coefficients for the s/ppolarization and k n = ω/c n and k z,n are defined in Appendix A. A number of details are worth mentioning. There are two distinct behaviors of the decay rate Γ full . At low frequency the curves shoot up due to Ohmic losses, which are not included in Eq. (80). At intermediate frequencies, the curves develop a clear resonance due to the excitation of surface plasmons in graphene. Also the maximum of the resonances blue-shifts with the decrease of the distance of the dipole to the graphene sheet. This behavior is easily explained remembering that the dispersion of the surface plasmons in single layer graphene is proportional to √ q. Since the distance z 0 introduces a momentum scale q ∼ 1/z 0 , smaller z 0 values correspond to higher q-values and higher energies of the resonant maximum. Equation (80) for the plasmon emission rate produces exactly the same resonance (same magnitude and same maximum of the resonance position) as Γ full , indicating that in this region, the decay rate of the quantum emitter is dominated by plasmon emission. This is further shown in Fig. 5, where losses were arbitrarily reduced in the evaluation of Γ full . We can also avoid the superposition of the Ohmic and surface plasmon contributions choosing either a larger Fermi energy or a larger distance from graphene to the metal. In Fig. 5 the ratio Γ/Γ 0 is  smaller than the ones in Fig. 4 due to the larger distance of the dipole to the graphene sheet.
An analytic expression for the plasmon emission rate can be also obtained for the case of graphene-metal structure assuming that κ 1 d 1. In this limit, Eq. 64 can be approximately solved, yielding ω q,sp = 4αdE F c/ q. Plugging this result in Eq. 79, we obtain This plasmon emission rate is shown in Fig. 6 as a function of ω and the graphene-metal separation d. As in the case of a single graphene layer, for each d there is a well defined peak as function of ω.

B. Decay rate with non-local effects
We now focus on the role of non-local effects in the graphene conductivity play in the emitter decay rate. We will restrict ourselves to the case of a quantum emitter close to a single layer of graphene. In Fig. 7, we compare the transition rate of a quantum emitter calculated taking into account non-local effects β = 0 in Eq. (37), with the local case β = 0. The left panel shows the result for E F = 0.6 eV, z 0 = 7 nm, while the right panel shows the result forE F = 0.6 eV, z 0 = 2 nm. We can see that the nonlocal effects play an important role at smaller emitter-graphene distances. As the distance between quantum emitter and the graphene layer, z 0 , decreases, non-locality leads to an increase and blueshift of the resonance in the transition rate as a function of the emitter frequency. We have also verified (not shown) that the non-local effects are also more prominent in the case of lower electronic densities. Figure 7. Comparison between the plasmon emission rate of a quantum emitter close to a single layer of graphene using the local (β = 0) and non-local (hydr.) (β = 0) models for the graphene conductivity (see Eq. (37)). The used EF = 0.6 eV, 1 = 2 = 3.9, ψ = 0. Left panel : z0 = 7 nm. Right panel: z0 = 2 nm. The nonlocal effects causes an increase and blueshift of the resonance in the transition rate as the distance z0 between the quantum emitter and the graphene decreases.

VI. CONCLUSIONS
In this paper, we have performed the quantization of graphene plasmons, in the absence of losses, and applied the field quantization to the interaction of an emitter with doped graphene. The quantization was performed using both a macroscopic energy approach and a quantum hydrodynamic model, which allows for the inclusion of non-local effects in the EM response of graphene. The quantization approaches used allow for the determination of the plasmon EM field mode-functions and, importantly, their normalization, which becomes non-trivial when dispersion is included.
When comparing the total decay rate of a quantum emitter (as obtained using the full EM dyadic Green's function) with the decay rate due to quantum emission, it was shown that plasmon emission completely dominates the decay rate, for typical emitter-graphene separations and emitter frequencies. It was shown that nonlocal effects in the graphene response, become increasingly important for smaller graphene-emitter separations and smaller Fermi energies.
The advantage of the quantization method developed in this work lies in its simplicity. The mode-functions are obtained from the solution of the Maxwell equations for the vector potential, and the normalization of the modefunctions can be expressed as a simply integral, which only involves the dielectric function of the medium. For situations when only a few modes contribute significantly to the physics, as in the case of quantum emission dominated by plamons, mode-functions allows to a much simpler and physically transparent description than the full EM dyadic Green's function. Since, the determination of the mode-functions only involves the solution of the classical Maxwell equations, the quantization of the electromagnetic field of more complexes structures can be easily achieved. We also note that the procedure gives both the quantized form of the electric and magnetic fields.
Therefore, we could also study the enhancement of the spin relaxation. The only change would be the interaction Hamiltonian, which in this case would be of the form H I = −µ · B, where µ is the magnetic dipole moment of the emitter and B the quantized magnetic field of the surface plasmon.
In possession of a quantized theory for graphene plasmons, we have set the stage for the future discussion of other quantum effects involving these collective excitation made simultaneously of light and matter.