Anomalous Reflection Phase of Graphene Plasmons and its Influence on Resonators

The phase picked up by a graphene plasmon upon scattering by an abrupt edge is commonly assumed to be $-\pi$. Here, it is demonstrated that for high plasmon momenta this reflection phase is $\approx -3\pi/4$, virtually independent on either chemical potential, wavelength or dielectric substrate. This non-trivial phase arises from a complex excitation of highly evanescent modes close to the edge, which are required to satisfy the continuity of electric and magnetic fields. A similar result for the reflection phase is expected for other two-dimensional systems supporting highly confined plasmons (very thin metal films, topological insulators, transition polaritonic layers, etc.). The knowledge of the reflection phase, combined with the phase picked up by the plasmon upon propagation, allows the estimation of resonator properties from the dispersion relation of plasmons in the infinite monolayer.

The phase picked up by a graphene plasmon upon scattering by an abrupt edge is commonly assumed to be −π. Here, it is demonstrated that for high plasmon momenta this reflection phase is ≈ −3π/4, virtually independent on either chemical potential, wavelength or dielectric substrate. This non-trivial phase arises from a complex excitation of highly evanescent modes close to the edge, which are required to satisfy the continuity of electric and magnetic fields. A similar result for the reflection phase is expected for other two-dimensional systems supporting highly confined plasmons (very thin metal films, topological insulators, transition polaritonic layers, etc.). The knowledge of the reflection phase, combined with the phase picked up by the plasmon upon propagation, allows the estimation of resonator properties from the dispersion relation of plasmons in the infinite monolayer.
Introduction.-Graphene can pave new ways for the development of nanoscale photonic and optoelectronic devices 1 . Graphene plasmons 2-10 (GPs) are especially interesting due to their ultra-strong confinement, that may lead to a strong enhancement of light-matter interaction [11][12][13][14][15] . Due to current limitations of the mobility of charge carriers in graphene samples, plasmon propagation lengths are not large compared to the free space wavelength λ. However, they can still be large compared to the graphene plasmon wavelength λ p which can be as small as λ/100. For this reason interesting resonant effects occur in graphene structures of deep subwavelength scales. For instance, resonant enhancement of absorption due to excitation of GPs has recently been proposed [16][17][18] and demonstrated 19,20 .
In particular, resonant plasmonic effects have been measured in arrays of ribbons in an ample range of ribbons widths 18,[20][21][22][23][24] . Often, as evident in these works, the analysis of the resonant absorption peaks has been based on the comparison with the dispersion relation of two-dimensional (2D) GPs via the relation k p = π/W , where W is the "active" part of the graphene ribbon (i.e. the region having a sufficiently high concentration of charge carriers). This simple relation between the modes inside the ribbon and 2D GPs corresponds to a Fabry-Perot cavity with ideally-reflecting walls and is justified by the high value of the GP reflection coefficient at a graphene edge 25,26 . Such an interpretation is however only valid when the phase shift of the GP reflected from the ribbon termination is −π.
In this paper we show an effect that has been neglected up to now: the influence of the phase picked up by a GP upon reflection at a graphene edge. We will demonstrate that this is a very important contribution, which affects both the resonator frequency and the nature of the electromagnetic modes.
1: (Color online) (a) Schematic of the reflection of a graphene plasmon by an edge. A graphene plasmon propagating from left to right (shown by a red arrow) acquires a phase ΦR and then propagates back (which is shown by a blue arrow). (b) Computed absolute value of the magnetic field |H| created by the interference between the incident graphene plasmon and that reflected from the edge. In the calculation graphene lies on a SiO2 substrate, and the considered free-space wavelength is 12µm. The parameters used for computing the graphene conductivity are: temperature 300 K, Fermi level 0.35 eV.
Phase picked up by reflection at an edge.-It is now well known that when a GP reaches an edge, it is reflected with virtually 100% amplitude 25,26 . This effect can be attributed to the large density of states on the GP channel (much larger than the density of photonic radiating modes).
In order to investigate the phase of the reflection coefficient, let us start by analyzing an idealized system: an homogeneously-doped lossless free-standing graphene plasmon resonator (ribbon). As will be shown later, despite its simplicity this model contains very useful information. We consider the situation, schematically represented in Fig. 1(a), where a GP is launched (from x < 0) at normal incidence towards a graphene edge, placed at x = 0, and compute 27 the full electromagnetic field. In these calculations, and throughout the manuscript, graphene is modeled by its 2D conductivity σ, obtained from the random phase approximation 4,5,28 . Fig. 1(b) renders the absolute value of the total magnetic field. The single-period interference pattern with regions of zero field in the plane confirm that the scattered field is essentially a back-reflected GP.
The different contributions to the magnetic field are represented in Fig 2(a), clearly showing that, beyond a transition region at very small distances to the edge, the incident and reflected wave have the same period and virtually the same amplitude, but present a phase shift.
To characterize this phase shift we define the auxiliary function Φ R (x) through the reflected field H r : Defined in this way, Φ R (x) is just a reparametrization of the phase of the reflected field but, at sufficiently long distances from the edge, it converges to the phase pick up by the GP upon the reflection at the edge, Φ R ≡ Φ R (x → −∞) . Fig 2(b) shows the computed Φ R (x) for several wavelengths, both for free standing graphene and for graphene on a SIO 2 substrate. Two regions are apparent in this figure. First, at small distances to the edge (smaller than 50 − 100 nm), Φ R (x) varies spatially, due to the excitation of deeply evanescent modes generated by the termination. These modes are reminiscent of the "quasi-cylindrical" waves originated by subwavelength emitters in both metals and graphene 29,30 and, in this problem, are responsible for a change in phase of the order of 15%. Second, at larger distances, Φ R (x) converges to a GP reflection phase Φ R ≈ −0.75π. Notably, this value is virtually independent of the intrinsic properties of graphene (Fermi level, temperature, etc) and dielectric environment.
It is possible to obtain a good estimation of this non-trivial reflection phase, and further insight into its origin, using a simple analytical model. Following an approach that was used to study the reflection of the surface plasmons in metal slabs 31 , we neglect the complex behavior of the EM field close to the edge, and approximate the field at the graphene side by the sum of the incident and reflected GPs. The field is expanded in vacuum (x > 0) in its eigenstates, thus having the form: x > 0 (2) where r E and h(k z ) are expansion coefficients, ± stay for z > 0 and z < 0 respectively (the magnetic field for a GP is antisymmetric with respect to the graphene position). The wave vector components are k zp = −g/α (with g = ω/c and α = 2πσ/c) and k p = g 2 − k 2 pz , k x = k x (k z ) = g 2 − k 2 z , with Im(k p ), Im(k x ) ≥ 0 to satisfy the radiation conditions. The electric field components can be found from Eqs. (2) using Maxwell's equations. The unknown coefficients r and h(k z ) are obtained by imposing continuity of the magnetic fields and the z-component of the electric field along the plane x = 0 (see the Supplementary Information for the details). The final result for the reflection coefficient is where This integral can be computed numerically (the results can be found in the Suppl. Mat.), but an excellent analytical approximation can be found for the case with large values of k p (which is the case of interest, due to strong confinement of the GP). Noticing that the main contribution to the integral originates from k z ∼ |k pz | 1, we can approximate k x ik z , and obtain a = −2i/π, yielding a reflection phase Φ R = arctan(−4π/(4 + π 2 )) ≈ −0.64π. Thus, this analytical result coincides with the phase extracted from full-wave calculations with 13% error, which is consistent with the neglected contributions from the highly confined modes in the graphene region that lead to the formation of the quasi cylindrical wave close to the edge.
In order to get a better description of experimentally observed samples, let us extend the previous discussion to situations where the carrier density, and thus the conductivity σ = σ(x), changes spatially close to the edge. In this case the reflection coefficient corresponds to the reflection from the whole inhomogeneous region, not only from the abrupt termination x = 0. However, taking into account that GPs can follow adiabatically variations in conductivity in length scales even smaller than the GP wavelength 25 , we can hypothesize that the GP scattering can then be split into 3 "events": (i) GP propagating through the inhomogeneity to the edge; (ii) reflection of the GP from the edge with the reflection coefficient r E ; (iii) backward propagation of the reflected GP r through the inhomogeneous region. In other words, the full reflection coefficient is simply the product of the reflection coefficient for the edge r E and phase due to the propagation of the GP along the inhomogeneity region (see details in Suppl. Mat.): with k p0 = k p (−∞). In order to check the validity of this approximation, we have carried out simulations for the two different conductivity profiles, rendered in the inset to Fig 3(a). These profiles present a ∼ 1/ √ x-like dependency for F (x) (see details in the Suppl. Mat.) since similar behaviour of the Fermi level was found in electrostatically doped graphene ribbons 32 . Fig 3(a) shows that the argument of the reflection coefficient is very different for both profiles (rendered by the continuous curves). However, when corrected by the optical phase picked up by propagation in the non-uniform region this argument converges to approximately −0.7π. This phase is very close to the previously found value for the edge of a homogeneous sheet, thus confirming that, to a good approximation, the GP follows adiabatically the conductivity variation up to the edge.
GP Modes in ribbons.
-Several experiments have analyzed the modes in arrays of graphene ribbons, illuminated at normal incidence with an electric field pointing perpendicularly to the ribbon axis, in terms of a Fabry-Perot model, and arg(r) − ∆ϕ ("corrected" phase, discontinuous curves) are rendered. The parameters used for graphene are the same as in Fig. 2 where excited plasmons propagate back and forth between the ribbon edges. In a single homogeneous ribbon, with width W , the modes should satisfy All analysis have, to the best of our knowledge, consider that the field vanishes at the edges, thus assuming Φ R = −π. In this case, the fundamental mode occurs for n = 0, leading to the commonly used expression k p = π/W . However, the existence of the non-trivial reflection phase analyzed in this manuscript modifies this picture. For instance, notice that the fundamental mode occurring for n = 0 satisfies k p = −Φ R /W . In order to illustrate this point, we have carried out simulations (using modal expansion 16 ) of the transmission spectra T of normal-incidence electromagnetic waves, with electric field polarized normal to the ribbons, impinging into arrays of ribbons placed on an insulating substrate. Resonances in the relative transmission, i.e. with respect to the graphene-free substrate case, correspond to enhanced absorption due to the excitation of graphene plasmons. In order to compare with experimental data 20 , we have assumed an array period L = 2W , a temperature of 300 K, a typical carrier relaxation time of 0.1 ps for CVD graphene, and SiO 2 substrate whose dielectric function is taken from Ref. 20 . Fig. 4 (a) shows the intensity plot of the relative transmission computed for ribbon arrays with different widths. The white discontinuous lines provide the expected position of the resonances, based on Eq. (6), where k p is given by the dispersion relation for graphene plasmons on a SiO 2 substrate. The result of the mapping for the first three even modes n = 0, 2, 4 is shown by the white discontinuous curves (modes char-  Fig. 4 (a) might be attributed to the growing influence of non-homogeneously doped edges for very thin ribbons. In Fig. 4b we render the computed results for a freestanding graphene sheet. In this plot we show the full calculation for a ribbon array, together with prediction from the traditionally used expression k = π/W (labeled "previous, n = 0") and that with the anomalous reflection phase using Eq. (6). Also in this case, the traditional model captures the trends but fails to reproduce accurately the spectral position for maximum transmission. In contrast, our simple Fabry-Perot model with non-trivial reflection phase Φ R provides a good approximation to the computed resonant spectra for ribbons, specially for small momentum (larger ribbons widths). Notice that, for the considered ratio between period and ribbon width, the ribbonribbon interaction affects mainly the lowest frequency resonance (n = 0) 16,23 and provides a slight red shift (around a few percent) of the resonance peak compared to the value predicted by Eq. (6). Calculations show that already for slightly larger inter-ribbon distance (L = 3W ) Eq. (6) provides a virtually exact value for the resonance spectra (see Suppl. Mat.).
Conclusions.-We have shown that when plasmon in graphene is reflected from an abrupt termination, the phase of the reflected plasmon is largely independent on dielectric environment, wavelength and local doping provided that the GP is strongly confined. The reflection phase has a nontrivial value of Φ R ≈ −3π/4. We reached our results using three different methods (finite elements calculation, quasi-analytical model and a fully analytical approximation). Taking this phase into account is essential to design and analyze optical properties of graphene resonators (ribbons, discs, etc.). By appropriately remapping of the modes inside a single graphene ribbon, we have successfully fit the experimentally observed resonances in periodic arrays of graphene ribbons. Our results can be also further extended to other two-dimensional systems supporting plasmons (like very thin metal films, topological insulators, transition polaritonic layers, etc.), where phase of the reflection coefficient is normally, and perhaps incorrectly, assumed to be −π.
Acknowledgements.-We express our sincere gratitude to F. Guinea for stimulating and fruitful discussions. We acknowledge support by the Spanish Ministerio de Economia y Competitividad within project MAT2011-28581-C02. LMM acknowledges support by the EC under Graphene Flagship (contract no. CNECT-ICT-604391).

Appendix A: Homogeneously doped graphene
We will solve the problem of the graphene plasmon (GP) reflection from the edge by assuming that a single mode (GP) from the region x < 0 couples to the "electromagnetic continuum" in the region x > 0. For simplicity and symmetry we assume that graphene is free-standing. Let us represent the magnetic field of the incident GP in the region x < 0 as follows (we use the notation H ≡ H y ) where + and − correspond to the regions z > 0 and z < 0 respectively and with g = ω/c being the free-space wavevector and α = 2πσ/c being the normalized conductivity of graphene. Using the Maxwell equations, z-component of the electric field reads We will seek for the reflected GP in the following form: where we have introduce the reflection coefficient r. Then the total fields in the region x < 0 read Now let us represent the field in the region x > 0 in the form of the plane waves expansion where The boundary conditions along the line x = 0 read Using (A6) and (A5) we have where ψ(z) is the GP profile function At x = 0 the z-component of the electric field must be continuous, so we project the second line of (A9) by e −ikz z 2π . However, the magnetic field does not need to be continuous, due to the existence of surface currents, so we project the first line of (A9) by ψ(z). After integration, we obtain: Substituting h(k z ) from the first line into the second one, we find the linear equation for r: The integrals in z, z are trivially performed Then (A11) becomes From here we find the reflection coefficient explicitly In the last expression the dimensionless components of the wavevectors read The results of the numeric calculations for the reflection coefficient phase as a function of the GP momentum are shown in Fig.2 (blue curve). The integral in (A15) has been calculated by using Cauchy theorem to avoid the branch point singularity q x = 0 (q z = 1).
The integral (A15) can be strongly simplified if we take into account that the major contribution comes from large wavevectors, i.e from those with q z 1. In this case q in the denominator of the integrand can be replaced by q iq z . Then the integration is analytic and the result reads simply so that the result is independent upon the GP wavevector and thus properties of graphene. According to Eqs. (A17),(A14) the reflection phase reads arg(r) = arctan(−4π/(4 + π 2 )) ≈ −0.64π. (A18) This result is shown in Fig. 2 (red line) and is in an excellent agreement with the numerical simulations.

(B6)
Once the total field is known (for instance from the full wave simulations), both the local reflection coefficient r(x) and reflection coefficient for the edge r E can be extracted from Eqs. (B4)-(B6).

The shape of the conductivity profile
For the calculations of the reflection phase in case of the inhomogeneous profile of the conductivity, we will use the inverse square root type dependency in the vicinity of the edge ∝ 1/ √ x. In order to avoid the singularity at x = 0, we displace the singular point to positive values of x, so that the profile has the following form where σ 0 is the conductivity far away from the perturbation (or "background" conductivity), σ 0 = σ(−∞) under the assumption that F (−∞) = 1. The function F (x) reads with f (x) = 1 + a |x − δ| , where δ and a are the parameters controlling the maximal value of the conductivity at the position x = 0 (at the edge) and the growth, while x 0 presents the distance from the edge where the conductivity saturates to the constant value. In the calculations for the paper we take x 0 = −1 µm and δ = 10nm.
Notice that we do not follow any special model for the shape of the conductivity, but take the above profile for the proof of principle.
Appendix C: The ribbon-ribbon interaction Fig. 4 shows the effect of the distance between the ribbons in the array. One can clearly see that for a larger separation (period) L the lowest-frequency resonance (n = 0) in the relative transmission is much closer to the curve predicted by a model (the white line).