## 1 Introduction
Moreover, the optical response of metals clearly depends on the frequency, and therefore, we have to take into account the non-locality in time and space by generalizing the linear relationships to:  
\[\mathbf{D}(\mathbf{r},t)=\epsilon_{0}\int dt^{\prime}d\mathbf{r}^ {\prime}\varepsilon(\mathbf{r}-\mathbf{r}^{\prime},t-t^{\prime})\mathbf{E}( \mathbf{r}^{\prime},t^{\prime})\] \[\mathbf{J}(\mathbf{r},t)=\int dt^{\prime}d\mathbf{r}^{\prime} \sigma(\mathbf{r}-\mathbf{r}^{\prime},t-t^{\prime})\mathbf{E}(\mathbf{r}^{ \prime},t^{\prime}) \tag{1}\]  
It should be noticed that we have implicitly assumed that all length scales are significantly larger than the lattice spacing of the material, i.e., the impulse response functions do not depend on absolute spatial and temporal coordinates, but only on their differences. A fundamental relationship between the dielectric function and the conductivity is given by:  
\[\varepsilon(\mathbf{k},\omega)=1+\frac{i\sigma(\mathbf{k},\omega)}{\epsilon_{ 0}\omega} \tag{2}\]  
The general form of the dielectric response \(\varepsilon(\mathbf{k},\omega)\), in the interaction between the light and metals, can be simplified to the limit of a spatially local response through \(\varepsilon(\mathbf{k},\omega)=\varepsilon(\omega)\). The simplification is valid as long as the wavelength in the material is significantly longer than all characteristic dimensions such as the size of the unit cell or the mean free path of the electrons. In general, \(\varepsilon(\mathbf{k},\omega)\) is a complex valued function, and the imaginary part of the dielectric function determines the amount of absorption inside the medium. Most importantly, including quantum interlayer contributions leads to increasing the imaginary part of the dielectric function, and it turns out that the effects of quantum mechanics are very vital in systems in which interlayer transitions play an important role. It is basically known that the traveling-wave solutions of [PERSON]'s equations in the absence of external stimuli is given by:  
\[k^{2}\mathbf{E}-\mathbf{k}(\mathbf{k}\cdot\mathbf{E})=\varepsilon(\mathbf{k}, \omega)\frac{\omega^{2}}{c^{2}}\mathbf{E} \tag{3}\]  
Two cases, depending on the polarization direction of the electric field vector, need to be distinguished. For transverse waves, where \(\mathbf{k}\cdot\mathbf{E}=\mathbf{0}\), yielding the generic dispersion relation and more intriguingly for longitudinal waves, where \(\mathbf{k}(\mathbf{k}\cdot\mathbf{E})-k^{2}\mathbf{E}=0\), it requires the particular condition where \(\varepsilon(\mathbf{k},\omega)=0\), signifying that longitudinal collective oscillations can only occur at frequencies corresponding to zeros of the dielectric function. Plasmon modes in doped graphene, which are obtained from the condition in which \(\varepsilon(\mathbf{k},\omega)=0\), show many special properties, and some of them are listed here. * Their momentum is larger than the light momentum with the same energy [4];
* They can be actively tuned through chemical doping or electrical gating in real time (tuning charge carriers) [5,6];
* They illustrate higher levels of confinement (\(\lambda_{SPP}/\lambda_{light}\simeq 0.025\) in the normal case) [4];
* They have a longer lifetime and propagating lengths (\(\tau\simeq 500\) fs) [7];
* They occur in the terahertz and mid-infrared modes, which are absent in normal metals [4];
* They can be coupled with quasiparticles (for instance, generating plasmarons) [8];
* They are used in the quasi-zero dimension as emitters [9,10];
* They show very particular properties in hybrid structures (combining graphene with other 2D crystalline materials) [11,12];
* At long wavelength limits, they behave like \(\sqrt{n^{1/2}q}\), which is proportional to the charge career as \(n^{1/4}\)[4]. Those properties of plasmon modes in doped graphene have been measured by using the near-field optical microscopic technique with nanometer resolution.