Intrinsic Plasmon-Phonon Interactions in Highly Doped Graphene: A Near-Field Imaging Study.

: As a two-dimensional semimetal, graphene o ﬀ ers clear advantages for plasmonic applications over conventional metals, such as stronger optical ﬁ eld con ﬁ ne-ment, in situ tunability, and relatively low intrinsic losses. However, the operational frequencies at which plasmons can be excited in graphene are limited by the Fermi energy E F , which in practice can be controlled electrostatically only up to a few tenths of an electronvolt. Higher Fermi energies open the door to novel plasmonic devices with unprecedented capabilities, particularly at mid-infrared and shorter-wave infrared frequencies. In addition, this grants us a better understanding of the interaction physics of intrinsic graphene phonons with graphene plasmons. Here, we present FeCl 3 -intercalated graphene as a new plasmonic material with high stability under environmental conditions and carrier concentrations corresponding to E F > 1 eV. Near- ﬁ eld imaging of this highly doped form of graphene allows us to characterize plasmons, including their corresponding lifetimes, over a broad frequency range. For bilayer graphene, in contrast to the monolayer system, a phonon-induced dipole moment results in increased plasmon damping around the intrinsic phonon frequency. Strong coupling between intrinsic graphene phonons and plasmons is found, supported by ab initio calculations of the coupling strength, which are in good agreement with the experimental data.

D oped graphene is of fundamental importance for a plethora of applications, including graphene plasmonic device technology development. 1−7 Control over the doping level enables surface plasmons (SPs), whose wavelength scales with the fourth root of the Fermi energy E F . In addition, SPs interactions with the zone-center, in-plane optical phonons in graphene are intriguing, as they are predicted to exhibit unconventional mixing of plasmon and phonon polarizations due to the two-dimensional (2D) character of this material. 5 Although only Raman-active in single-layer graphene, the optical phonon becomes infrared (IR) active in bilayer graphene, that is, it acquires a net dipole moment due to broken inversion symmetry and variation of the electron− phonon coupling under an electric field. 9 This makes it amenable to many other advanced applications, such as phonon lasers 9 as well as nonlinear optical, quantum communication, and slow-light devices. 10 Recently, calculations by Low et al. 11 predicted novel plasmonic behavior in bilayer graphene at higher carrier concentration levels and in particular giant plasmonic enhance-ment of IR phonon absorption. Moreover, taking advantage of plasmon−phonon interactions in bilayer nanoribbons, Yan et al. 12 experimentally demonstrated Fano-like phonon-induced transparency. However, many questions regarding the nature of interactions between surface plasmons and intrinsic phonons in graphene remain unexplored, such as the strength of the plasmon−phonon coupling and the effects of dipolar moment contributions to plasmon decay. For instance, optical phonons in monolayer graphene are regarded by many as a source of plasmon damping, 13 despite the lack of a net dipole moment of those phonons, which one would expect to result in reduced plasmon−phonon interaction.
In this Letter, we study for the first time propagating plasmons in highly doped intercalated graphene. In particular, we use near-field plasmon imaging to study the effects of electron−phonon interactions on the plasmon dispersion and damping. We show that the phonon-induced dipole in bilayer graphene is distinctly different from that in monolayer graphene and leads to hybridization of the plasmon mode at the phonon frequency in the former, as evidenced by anticrossing behavior in the dispersion of the material. In monolayer graphene, inversion symmetry leads to a cancellation of the induced-charge contributions coming from each of the two atoms in the unit cell, thus resulting in a vanishing dipole moment within the unit cell for zero parallel (in-plane) wavevector. In contrast, the lack of inversion symmetry in bilayer graphene leads to a net phonon-induced dipole per unit cell, therefore enabling hybridization between plasmonic and vibrational modes. This is depicted in the diagram in Figure 1a, where the atomic structure of monolayer (top) and bilayer graphene (bottom) is presented, together with the phononinduced charge densities for each of the represented normal modes. In addition, we compare our experimental results to theoretical calculations from which we obtain the value of the phonon-induced dipole moment per unit cell by means of quantum chemistry methods using the Gaussian package. We show that experiment agrees well with theoretical results.
In order to investigate the effects of electron−phonon interactions, we performed plasmon nanoimaging measurements in graphene at unprecented high intrinsic carrier concentrations (corresponding to E F > 1 eV), achieved by intercalating graphene with FeCl 3 . 14 This new material has been recently used in optoelectronic applications 15 and it showed an extraordinary linear dynamic range, 16 together with an unexpected resilience to environmental conditions. 17 To our knowledge, this is the first near-field observation of propagating plasmons with energies exceeding 0.2 eV.
We studied samples with two different types of intercalated graphene structures: (1) one in which FeCl 3 intercalates in between each graphene layer of a flake, which was previously a bilayer, to form two separated single layers and (2) one in which a bilayer flake remains free of FeCl 3 in between its layers after fabrication. A diagram of the atomic structure of these two cases can be seen in Figure 1b. The FeCl 3 -intercalation resulted in carrier concentrations n ≥ 4.8 × 10 13 and 5.4 × 10 13 cm −2 for the two-monolayer and bilayer graphene samples, respectively. This corresponds to a Fermi energy E F of 1.21 eV (∼0.8 eV for each of the two monolayers) on the two-monolayer sample, and 1.4 eV for the bilayer graphene. The high carrier densities result in an upshift of the G-bands in the Raman spectra of the samples with respect to pristine graphene, as well as a variation in the shape of the 2D bands, as can be seen in Figures S1 and S2 of the Supporting Information, thus allowing for quantitative extraction of carrier concentrations from the careful analysis of the Raman features. 18,19 The carrier concentrations obtained using this method are also in excellent agreement with the plasmon dispersion measured on each corresponding sample.
Scattering-scanning near-field optical microscopy (s-SNOM) 20,21 was used to image the near-field properties of highly doped graphene with ∼20 nm resolution and at a broad frequency range. Figure 1c presents a diagram of the

Nano Letters
Letter DOI: 10.1021/acs.nanolett.7b01603 Nano Lett. 2017, 17, 5908−5913 experimental scheme where a tunable IR/MIR laser is focused onto the AFM tip of our s-SNOM setup in order to study the near-field on a FeCl 3 -intercalated graphene sample. The graphene lies on a SiO 2 /Si substrate. The confined electromagnetic field at the AFM tip provides the momentum needed to excite propagating plasmons in the graphene. These plasmons can be reflected at the edge of the flake 22,23 and subsequently interact with the tip. This results in a scattered optical field from the tip, which is detected in the far-field. As such, and as previously reported, 20,21,24 these type of measurements allow for the detection and real-space imaging of graphene plasmons from which the SPs wavelength (i.e., momentum ∼1/wavelength) can be measured, as shown in Figure 1c. All the measurements were performed at different instances over a period of time exceeding 18 months and under ambient conditions with no observable signs of sample degradation. The plasmon dispersion is extracted from the abovementioned near-field images by probing the plasmon wavelength (obtained by fitting to an oscillatory response) for a wide range of excitation frequencies. Figure 2 shows the dispersion of the measured plasmons for the two-monolayer (blue dots) and bilayer system (red dots). The error bars are smaller than the size of the symbols. The dashed magenta curve represents the plasmon dispersion resulting from the simple expression for the in-plane wavevector k ∥ , parallel to the graphene layer ω πσ where σ is the 2D conductivity, ω is the frequency, and ϵ s is the substrate permittivity. This corresponds to an uncoupled (freestanding) graphene monolayer. Because of strong coupling between the graphene plasmons and the SiO 2 substrate phonons, the plasmon dispersion is modified, which is captured by including the frequency dependence of the substrate permittivity ϵ s (ω). A comparison between the model and the experimental data shows clear evidence of splitting of the plasmon dispersion due to the hybridization of the plasmon mode with SiO 2 substrate surface phonons, 13,25,26 similar to previous reports on SPs in nanoribbons 27 and nanodots 8 on a SiO 2 substrate. In general, the data and model are in reasonable agreement when σ is substituted by a Drude conductivity with Fermi energy E F = 1.21 eV. This confirms that the Drude weight of the plasmon is very high, thanks to the efficient intercalation process. In this study, plasmons were probed up to energies of ∼0.28 eV (see Supporting Information, Figure S5). We remark that the interband region and the intraband region (indicated by the dashed green line) are "far" from the plasmon dispersion indicating the potential for this highly doped graphene to carry plasmons up to energies of 1 eV. Further studies are required to demonstate this, as it is very challenging to couple light with high energy graphene plasmons because the plasmon wavelength is very small (below 20 nm). While the dispersion of the bilayer graphene is expected to be shifted due to the higher effective Fermi energy in comparison to the two-monolayer system, it is also important to note the modal anticrossing feature present in the bilayer graphene at the graphene phonon frequency of ∼1585 cm −1 . This feature is not observed on the two-monolayer sample. This can be seen more clearly in Figure 3, which displays s-SNOM 2D images of the two-monolayer (top) and bilayer (bottom) samples. The scanned areas shown for each sample are approximately the same. The images correspond to incident frequencies of 1590, 1586, and 1569 cm −1 for the twomonolayer sample and frequencies of 1600, 1585, and 1570 cm −1 for the bilayer sample (additional s-SNOM images for a wider frequency range are shown in the Supporting Information, Figure S3). For the two-monolayer sample, the plasmon interference fringes show a very similar pattern for the three frequencies and exhibit a gradual decrease in spacing for increasing frequencies, as expected. The bilayer sample behaves distinctly different: as the incident frequency approaches the phonon frequency, the visibility of the fringes 24 decreases and eventually almost no interference fringes are observed at ∼1585 cm −1 . As the incident frequency further increases, away from the graphene phonon frequency, the fringes begin to gradually reappear and their amplitude increases, showing prominent, muliple fringes in each image on both sides of the phonon frecuency (see also Figure S3).
The left panels in Figure 3 correspond to the dispersion data within the dashed-rectangle region in Figure 2 for the twomonolayer (top) and bilayer (bottom) samples. The data highlights clear diferences between the dispersive behavior of the two samples. Here, an anticrossing behavior at the graphene intrinsic phonon frequency is clearly observed in the dispersion of the bilayer sample (bottom) whereas this feature is absent on the two-monolayer sample (top). The frequency interval between measurements around the graphene phonon frequency was not larger than 4 cm −1 , which is smaller than the energy width of the anticrossing (Δω ∼ 37 cm −1 ). This leads to a quantitative value for the electron−phonon coupling strength Δω/ω 0 ∼ 2.3% in units of uncoupled frequency ω 0 . 5 We attribute the splitting in the plasmon dispersion of the bilayer graphene to the breaking of inversion symmetry in this system, causing a finite dipole moment and rendering the graphene phonon mode IR-active. The background color plots in Figure 3 (left panels) correspond to the calculated imaginary part of the Fresnel reflection coefficient, which is a suitable magnitude to reveal the dispersion and strength of optical modes in graphene. 6 This is in turn obtained from the 2D conductivity, where we include contributions from intrinsic phonons in the latter, as obtained from ab initio calculations. 28 More precisely, an additional phonon term is incorporated, which is proportional to the squared dipole associated with the carbon phonon mode per unit cell. Our calculations produce a net dipole in the bilayer structure, which is in contrast to a vanishing dipole in the monolayer for vertical optical transitions. This readily translates into the presence of a phonon-induced splitting of the plasmon mode in the carbon bilayer, which is absent in the twomonolayer sample. Good agreement between the experimental data and the calculated splitting is observed.
We now turn our attention to the impact of the plasmon− phonon coupling on plasmon lifetime. Given a plasmon wavevector with real and imaginary parts k 1 and k 2 , respectively, depending on the different cases, the damping γ = k 2 /k 1 of graphene SPs is typically determined by possible decay channels such as impurity and defect scattering, electron− hole pair creation, as well as electron−phonon scattering. 5,13 Although radiative and Landau damping can also play a role in some particular cases, their influence is expected to be negligible in our experiments. 13 Here, we obtained the experimental value of γ for each incident frequency by fitting the line profile of s-SNOM 2D maps to an exponential function using the least-squares fit method (see Supporting Information, Figure S4).
A comparison between the frequency dependence of the modal lifetime of the two samples is shown in Figure 4. The data corresponding to the bilayer sample show a sharp lifetime decrease at the phonon frequency. This is in stark contrast to the plasmon lifetime of the two monolayer sample at the same frequency. The lifetimes of both samples experience a gradual decrease above the phonon frequency, possibly due to the opening of an additional decay channel beyond the phonon frequency, as reported by Yan et al. 13 However, at these frequencies, the lifetime of the two-monolayer sample is larger than the corresponding lifetime of the bilayer. It is important to note that the measured lifetime corresponds to a hybrid mode, which is plasmonic and phononic in character. 23,29,30 Additionally, however, our definition of the lifetime does not take into account variations of the group velocity of the mode. Thus, although the plotted values in Figure 4 underestimate phononic contributions to the lifetime, our approach allows for a quantitative analysis of the lifetime near the graphene phonon resonance, particularly around the phonon frequency. The results also show that electronic damping of the measured mode is affected by electron−phonon interactions, particularly at the phonon frequency and also beyond. The data for the In conclusion, we have observed propagating plasmons in ultrahighly doped graphene with Fermi energy exceeding 1.2 eV and a splitting in the plasmon dispersion due to the interactions of plasmons with the finite dipolar moment of intrinsic optical phonons in bilayer graphene. This effect is not observed in two-monolayer graphene because in such case there is no effective phonon dipole moment at the zone center due to inversion symmetry. We further present ab initio calculations in excellent quantitative agreement with the measured dispersion relations, confirming the existence of a net dipole associated with phonons in bilayer graphene but not in monolayer graphene. Our results correlate with lifetime measurements and show that electron−phonon interactions play a substantial role in plasmonic decay at and beyond the graphene intrinsic phonon frequency. These results are important for the development of future graphene plasmonic devices, particularly as the excitation of graphene plasmons at energies of 0.2 eV, the optical phonon energy, and above may be amenable to novel applications such as photodetectors and sensors.    Section S5. Dispersion of high-energy surface plasmons in graphene Figure 1: Raman spectroscopy and stacking order of intercalated graphene. a, Representative spectra of pristine (black) and FeCl 3 -intercalated (blue) three-layers graphene. The E 2g phonon of graphene is responsible for the resonant G 0 peak at 1580 cm -1 , while the A 1g mode is responsible for the double-resonant 2D peaks around 2650 cm -1 (inset). Upon doping with FeCl 3 , the G 0 peak upshifts and splits into two peaks G 1 and G 2 , each corresponding to a different stage of intercalation (see text), while the reduced coupling within the layers transforms the 2D band from a multi-Lorentzian to a single-Lorentzian (inset) peak. b, d, Optical images with annotated contrast and structure of the studied three-monolayer and bilayer flakes, respectively. c, e, Raman maps showing the normalized heights of the G-bands for the two flakes. The Silicon peak at 520 cm -1 is taken as reference. Yellow-dashed lines mark the region where plasmonic response is observed. Scale bars are 4µm in both panels.

Supplementary Information for Intrinsic
Using a combination of Optical microscopy, Raman spectroscopy and scattering-Scanning Near-field Optical Microscopy (s-SNOM) measurements, we were able to determine the stacking order of FeCl 3 -Few-Layer Graphene (FeCl 3 -FLG) and the total charge concentration of the stacked layers. The Raman spectrum of graphene intercalated with FeCl 3 shows two prominent features: an upshift of the G-bands, with respect to pristine graphene, and a change in the shape of the 2D band, from the convolution of multiple Lorentzian peaks to a single-Lorentzian fit, 1 as shown in Figure S1a. Charge transfer from the FeCl 3 molecules induces hole doping in graphene, which results in the stiffening of the E 2g mode 2,3 responsible for the resonant G-band at ∼ 1580 cm -1 in the Raman spectrum. Two intercalation stages are possible: 1) when a graphene layer is in contact with FeCl 3 from one side (known as stage-2 intercalation) the charge doping shifts the G-band from ∼ 1580 cm -1 to ∼ 1610 cm -1 (G 1 peak in Figure S1a). 2) Conversely, when a graphene layer is sandwiched between two FeCl 3 layers (known as stage-1 intercalation) the G band shifts up to ∼ 1625 cm -1 (G 2 peak in Figure S1a). 1 The increased distance between graphene layers reduces their coupling, causing the double-resonant 2D band, originating from the A 1g mode, to change from a multi-peak band, characteristic for two and three layer graphene, 4 to a single-Lorentzian peak, typical of stacked non-interacting graphene monolayers. 1 The stacking order of the flakes presented in the main text can be determined using optical microscopy and Raman spectroscopy. The pristine flakes are a three-and two-layer graphene ( Figure   S1b,d), as measured by optical contrast on Si/SiO2 (in Figure S1d a green filter was used to enhance the contrast of the bilayer, see ref. 1 Supplementary). After intercalation, the Raman maps ( Figure   S1c-e) of the two flakes show a prominent G 1 peak and no G 2 or pristine G 0 peaks, as shown also in Figure S2a. Hence, in this case, the intercalation is of stage-2 and all graphene layers must be in contact with FeCl 3 , resulting in the following stacking structures: the three-monolayers sample is intercalated so that one bilayer stacking is maintained (thus we refer to this sample as bilayer, Figure S1b); in the two-monolayers sample, one layer of FeCl 3 intercalates between the graphene, leaving two monolayers (thus we refer to this sample as two-monolayers, Figure S1c). Further confirmation of the bilayer structure being present in the first flake comes from the shape of the 2D band: the observation of two adjacent graphene layers shows up as a two Lorentzian peaks, instead of a single one, a situation not observed in the two-monolayer case, Figure S2a.
It has been shown by Bointon et al. 5 that the shift of the G-band can be related to the charge concentration in the graphene layers, using the model developed in ref. 3 , with an accuracy of 10% if compared to magneto-transport measurements. 1,5 Following the same procedure, we extract the charge concentration from the Raman maps shown in Figure S2b: statistics of the position of the G 1 peak show that, across the flakes, we have hole concentrations per layer of ρ ′ BL = 4.8x10 13 cm -2 and ρ ′ ML = 5.4x10 13 cm -2 for the bilayer (BL) and the two-monolayers (ML) samples, respectively.  and increases as the incident frequency moves away from it. We note particularly interesting features of the near-field images observed, specifically near the phonon frequency, where, even though plasmon waves propagating away from the edges gradually disappear, a pattern of circular "mounds" gradually appears throughout the scanned area. Although a direct comparison between all images reveals the presence of these "mounds" at most probed frequencies, albeit in some with less prominence, these features are most prominent exactly at the phonon frequency. At that frecuency, the first plasmon fringe at the edge of the graphene decreases in amplitude significantly, yet remains slightly visible. In addition, it is important to note that although the first plasmon fringe at the edge of the graphene decreases in intensity to a minimum near the phonon frequency, it is also always present on all measurements, in contrast to measurements performed at substrate phonon resonances, where the first fringe disappears completely (not shown), pressumably due to significantly larger plasmon damping at the latter.
Section S3: Profile line fitting  The open red circles and the blue curve correspond to data and our least-square fit, respectively.
The equation used for the fit describes the motion of a surface plasmon launched by the AFM tip of the s-SNOM, which propagates along the graphene surface and is reflected at the edge of the flake, similarly as described in the literature. 7,8 This fitting function can be written as, 7 where k p is the plasmon wavevector and γ is the inverse decay rate of the surface plasmon. The real coefficients A and B are fitting constants. The agreement between the data and the fit is excellent. In addition to the fitting, the empirical plasmon wavelength can also be obtained directly by measuring the fringe-to-fringe spacing (i.e., λ = 2x(fringe-to-fringe spacing)). 6-9 Figure 5: High-energy surface plasmons in graphene. Surface plasmon dispersion in graphene at high energies. The red circles correspond to experimental data obtained from the two-monolayer graphene sample. The inset shows a 2D s-SNOM image obtained under a 2240 cm -1 incident wavelength.

Section S4: Calculation of the Fresnel reflection coefficient
The Fresnel reflection coefficient is calculated as defined in the literature 10 in terms of the conductivity of free-standing graphene. Here, we include interaction effects between the graphene and the substrate surface phonons by incorporating these into the relationship, 8 where ω and k p represent the angular frequency and electron parallel wavevector, respectively, and κ = (ϵ s + 1)/2. The expresion for the graphene conductivity σ is described in detail in ref. 11 It is important to note here that substrate contributions to the dielectric function ϵ s can be parametrized as, with the longitudinal and transverse phonon frequenciesν LO andν T O , respectively, and damping rate Γ for the various phonon modes of the SiO 2 substrate. [12][13][14][15][16][17] Given the amorphous nature of the SiO 2 substrate, the values in the literature for the frequency of these modes can be found to slightly vary. 13,[16][17][18] In our case, the values obtained were ϵ ∞ = 2.52 withν T O = 820 cm -1 ,ν LO = 550 cm -1 and Γ = 10 cm -1 for substrate surface phonon 1 andν T O = 1110 cm -1 ,ν LO = 1290 cm -1 and Γ = 50 cm -1 for substrate surface phonon 2. Note that another substrate phonon mode, typically found at frequencies approximately below 460 cm -1 , was not included in our calculations, as this frequency lies far away from those studied in our experiments. 16 Section S5: Dispersion of high-energy surface plasmons in graphene Figure S5 presents experimental data for the dispersion of the two-monolayer graphene sample at higher energies. The data was obtained using the s-SNOM imaging technique, as discussed previously for measurements at lower energies. The graph shows results for incident wavelengths varying from 2105 to 2250 cm -1 . The inset displays a 2D s-SNOM image taken at 2240 cm -1 incident wavelength. Using a similar analysis as that described in the main text for the two-monolayer sample, the results in Fig. S5 follow the expected dispersion for higher energy plasmons. However, it is important to note that in addition to the significant decrease in plasmon wavelength, the fringe visibility and number of visible fringes also decrease, suggesting a significant increase in plasmonic losses at these frequencies. Such a plasmon damping increase can be the results of additional decay channels well above the optical phonon frequency, as described by the literature. 19,20 Further investigations will be necessary order to get a better understanding of this phenomena, however, a