## Abstract

Plasmonic response of the metallic structure characterized by sub-nanometer dielectric gaps can be strongly affected by nonlocal or quantum effects. In this paper, we investigate these effects in spherical Na and Au nanomatryoshka structures with sub-nanometer core-shell separation. We use the state-of-the-art quantum hydrodynamic theory (QHT) to study both near-field and far-field optical properties of these systems: results are compared with the classical local response approximation (LRA), Thomas–Fermi hydrodynamic theory (TF–HT), and the reference time-dependent density functional theory (TD–DFT). We find that the results obtained using the QHT method are in a very good agreement with TD–DFT calculations, whereas other LRA and TF–HT significantly overestimate the field-enhancements. Thus, the QHT approach efficiently and accurately describes microscopic details of multiscale plasmonic systems whose sizes are computationally out-of-reach for a TD–DFT approach; here, we report results for Na and Au nanomatryoshka with a diameter of 60 nm.

© 2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

## 1. Introduction

Plasmonic systems have been a topic of substantial research due to their promising applications such as optical sensing [1], nanolasing [2], surface-enhanced Raman scattering [3] and enhancement of nonlinear optical processes [4], to name a few. Metallic structures with sub-nanometer inter-distances make an important class of plasmonic systems presenting unique properties [5–8]. At such distances in fact conventional classical theories, such as Drude-type models, fail to correctly describe electron dynamics inside the metallic system. This deviation is due to microscopic details that are neglected in the traditional approach but that give rise to a nonlocal response or quantum effects, which become observable in the macroscopic properties of the system [9–12].

A correct theoretical description of the optical response of these systems can be obtained by first-principle quantum mechanical calculations, such as time-dependent density functional theory (TD–DFT) [13]. However, TD–DFT methods are applicable only to very small systems, a few nanometers in size. Most of the plasmonic components involved in experiments are tens to hundreds of nanometer in size and the TD–DFT is unable to handle such scales as their computational cost increases cubically with the number of electrons.

A computationally less demanding alternative to the first-principle approach is the hydrodynamic theory (HT) [14–17], in which the quantum electron dynamics is solved via macroscopic observable quantities, such as the charge density *n*(*r*, *t*) and the current density **J**(*r*, *t*). In the HT quantum properties are taken into account via complex expressions of the electron internal energy or pressure. In the local response approximation (LRA) in fact, induced charges can be squeezed into an infinitesimal thin layer at the surface of the metallic particle. The role of the electron pressure is then to smoothen the singularity so that the induced charges spread out from the metal surface.

The most common and simple approximation in the HT is obtained using the Thomas–Fermi (TF) expression for the electron pressure [18–22], which is usually combined with the hard-wall boundary condition (electrons cannot escape the metal boundary) and a uniform equilibrium electron density. Although the TF–HT gives a good qualitative description of gap-dependent resonance shifts in plasmonic systems [8, 9], it neglects important quantum effects such as electron spill-out and quantum tunneling. These effects can be incorporated in HT by taking into account the spatial dependence of the electron density. However, in order to avoid unphysical plasmonic modes [23,24], the TF kinetic energy needs to be enriched with a new ∇*n*-dependent term: the von Weizsäcker kinetic energy. This method is referred in the literature as the quantum hydrodynamic theory (QHT). The QHT has been applied to study optical properties of spherical and cylindrical monomer and dimer structures and it has been seen to efficiently describe the nonlocal or quantum effects [23,25–29].

Concentric core-shell nanostructures—also known as nanomatryoshka—offer an unusual large tunability due to strong coupling between core and shell elements. These structures could be fabricated in very small sizes with subnanometer core-shell spacing and might be useful in surface enhanced spectroscopy and biomedical applications [30,31]. These structures with dielectric or metallic core have been extensively studied in classical framework [32–36]. A full quantum mechanical treatment for spherical and cylindrical nanomatryoshka (NM) structures has been presented in [37,38]. They showed that at core-shell spacing below 5 Å, the charge transfer between the core and shell takes place and the the tunnelling effect substantially influences the absorption cross-section and field enhancement. Reference [39] presents an *ab initio* study of bimetallic core-shell nanostructures emphasizing the quantum effects on the optical response of these nanoparticles. However, the plasmonic systems considered in these works are much smaller than those of practical interest.

In this paper, we apply the QHT to investigate microscopic and macroscopic optical response of so-called nanomatryoshkas (NMs) which are composed of a solid metallic spherical core surrounded by a concentric metallic spherical shell. A systematic analysis within local, TF–HT and QHT approximations has been presented and the impact of core-shell spacing on the optical properties of Au and Na NMs has been studied. For Na NMs, we performed TD–DFT calculations which show that the QHT approach is very accurate both in reproducing the optical spectrum as well as the induced density. It is worth noticing however that the QHT does not fully incorporate the wave nature of an electron and, therefore, it does not capture the Friedel (quantum) oscillations. Similar findings are found for Au NM nanostructures, where TD–DFT results can be found in literature [37]. Finally, NMs with sizes much larger than those which can be handled by the TD–DFT have also been considered both for Au and Na metals. We take full advantage of symmetry of the structure and use a 2.5D technique for numerical implementation, which reduces the computational burden to a great extent [40].

## 2. Theory and computational details

#### 2.1. System details

We consider a spherical structure composed of a metallic solid core encapsulated by a thin metallic shell, as shown in Fig. 1. For simplicity we assume that the dielectric medium between the core and shell is vacuum. The dimensions of the NM are represented by (*R*_{1}, *R*_{2}, *R*_{3}) where *R*_{1} represents the radius of the core while *R*_{2} and *R*_{3} indicate respectively the inner and outer radii of the encapsulating shell. We model the NM using Jellium approximation [41,42] where the electrons in the metal are confined by uniformly-distributed positive background charges with density ${n}^{+}={\left({r}_{s}^{3}4\pi /3\right)}^{-1}$ inside and zero outside the metal boundary. The Wigner-Seitz radius *r _{s}* = 3.02 a.u. for Au and

*r*= 4 a.u. for Na metal.

_{s}For comparison, we perform calculation within the LRA. In this case we use a Drude-type dielectric function for the metal description:

where*ε*

_{∞}represents the background core permittivity,

*ω*is the spatially constant plasma frequency and

_{p}*γ*is the electron collision rate. For Au we choose

*ε*

_{∞}= 8,

*ω*= 9.07 eV and

_{p}*γ*= 0.135 eV, as given in [37] and for Na we use

*ε*

_{∞}= 1,

*ω*= 6.04 eV and

_{p}*γ*= 0.16 eV, taken from [43].

#### 2.2. Quantum hydrodynamic theory

Hydrodynamic equation of motion of an electronic system when coupled to Maxwell equations leads to the following system in frequency domain [28]:

**E**is the electric field,

**P**is the polarization and

*μ*

_{0}and

*c*represent the magnetic permeability and speed of light in free-space respectively;

*n*

_{0}(

**r**) is the spatially dependent equilibrium charge density,

*γ*is the damping rate,

*e*and

*m*are the electron charge and mass respectively. In the QHT the internal energy of the electronic system is: where

*T*

_{TF}and

*T*

_{vW}are the Thomas–Fermi and von Weizsäcker kinetic energy functional and

*E*

_{xc}indicates the exchange and correlation (XC) energy.

*σ*

^{(kxc)}in Eq. (3) is the viscoelastic kinetic-exchange-correlation tensor which takes into account the nonlocal broadening of the plasmon energies. Explicit expressions for the terms in Eq. (3) can be found in [27, 28]. The TF–HT approximation can be recovered if the XC and von Weizsäcker energy functional are neglected in Eq. (4) and the equilibrium electron density is assumed as constant, i.e.,

*n*

_{0}(

**r**) =

*n*

_{0}.

We solve Eqs. (2) and (3) under plane-wave excitation for **E** and **P** fields, which are further used to compute optical properties of the plasmonic systems. In particular, Eq. (3) is only imposed in a region that contains the NM. On the boundary of this region, far enough form the ion edge, we impose **P** = **0**; here the equilibrium density and hence **P** is practically zero (although not exactly since the electron density decays exponentially). Commercial software Comsol Multiphysics [44], which allows a flexible implementation of these equations, has been used. We exploit the symmetry of the geometry and use the 2.5D simulation technique [40], which significantly reduces the computational efforts in terms of memory and processing time. 2.5D method requires all fields to be written in terms of an azimuthal mode number *m*, such that for a vector field **v** can be expressed as: $\mathbf{v}(\rho ,\varphi ,z)={\sum}_{m\in \mathbb{Z}}{\mathbf{v}}^{(m)}(\rho ,z)\text{exp}[-\mathit{im}\varphi ]$, where *m* ∈ ℤ. The advantage of this method is that an initially three-dimensional problem is reduced to a few (2*m*_{max} + 1) two-dimensional problems. For sub-wavelength structures *m*_{max} < 3 is usually enough to accurately describe the problem.

The equilibrium electron density is computed self-consistently from the following nonlinear differential equation [25,28],

*ε*

_{0}is the permittivity of the free-space and

*n*

^{+}represent the uniform background positive charge. The ground-state density

*n*

_{0}(

**r**), which include the spatial dependence of the electron density, has a strong impact on the QHT optical response. To calculate the unperturbed electron density

*n*

_{0}(

**r**) of the NM self-consistently, we consider the following expression for the initial guess (as used in [27,28] for spheres):

*r*is the distance from the center of the NM and

*κ*represents an asymptotic decay of the electron density. The solution is then found solving Eq. (5) iteratively.

#### 2.3. DFT and TD–DFT calculations

DFT data are obtained by solving self-consistently the Kohn-Sham (KS) equations and computing the TD–DFT optical spectra of jellium NM systems with the Wigner-Seitz radius *r _{s}*=4. Calculations have been performed using an in-house developed code using a uniform numerical grid. The implemented TD–DFT algorithm follows the literature [41,45,46]. We used an empirical broadening of 40 meV. In this work we only consider jellium NM systems with filled shells. A NM system described by an arbitrary (

*R*

_{1},

*R*

_{2},

*R*

_{3}) may not have filled shells nor an integer number of electrons. Thus, for a given NM system described by (

*R*

_{1},

*R*

_{2}=

*R*

_{1}+ Δ,

*R*

_{3}) we first considered only the core of the jellium nanomatryoshka system, i.e. a jellium nanosphere with

*N*electron and we searched for the filled-shell jellium nanosphere with a radius

_{c}*R′*

_{1}closest to the target one

*R′*

_{1}≈

*R*

_{1}. Then we added the shell with

*N*electrons considering the exact gap region (i.e.

_{s}*R′*

_{2}=

*R*

_{1}+ Δ) and an external radius giving a magic-number total number (

*N*=

_{t}*N*+

_{c}*N*) of electron and an external radius

_{s}*R′*

_{3}closest to the target one

*R′*

_{3}≈

*R*

_{3}).

## 3. Numerical results and discussion

#### 3.1. Optical response of scaled NM structures

In this section, we consider five NM structures (see Fig. 1) with dimensions SF ×(*R*_{1}, *R*_{2}, *R*_{3}) with *R*_{1} = 8.5 Å, *R*_{2} = 9.5 Å and *R*_{3} = 15.9 Å where SF is the scaling factor taking an integer value from 5 to 1 corresponding to a core-shell separation from 5 to 1 Å. According to classical electrodynamics, scaling of a plasmonic system does not change its optical properties: only the amplitude of the spectra is affected while the plasmon resonance remains exactly at the same position [47]. Since the particle sizes under study are much smaller than the operating wavelength, therefore, the retardation effects can be entirely neglected [37]. The NM structure supports three hybridized plasmon modes arising from the interaction of the core mode and the two bonding and antibonding modes of the shell. [37,38]. In the classical framework the Au (Na) NMs support a lower energy mode (LEM), localized in the gap, approximately at 1.77 (1.38) eV and a higher energy mode (HEM) at 2.86 (3.46) eV in the absorption spectra, whereas the third mode is not visible in the spectrum due to low dipole moment.

Absorption efficiency (i.e. the photoabsorption cross-section divided by the geometrical cross-section, $\pi {R}_{3}^{2}$) for the five NMs computed within the LRA, TF–HT and QHT methods for Au (left panels) and Na (right panels) are shown in Fig. 2. It is evident from the LRA calculations that the scaling of the structure affects only the relative amplitude of the absorption spectrum for both Au and Na. The positions of lower and the higher energy resonant modes are unchanged after such scaling. The TF–HT predicts a blue-shift with respect to local calculations for both lower and higher energy modes with a larger shift for smaller SF (or gaps). This shift of the HEM peak towards the higher energies has already been found in dimer structures with sub-nanometers gaps [48].

The QHT, on the other hand, predicts spectra whose shift direction (with respect to LRA) depends on the material. For Au NMs, as in the TF–HT case, HEM peaks undergo a blue-shift (however smaller) as the SF decreases from 5 to 1. Conversely, for Na NMs, the QHT predicts a red-shift with respect to LRA. When the gap between the core and shell decreases the LEM starts vanishing due to electron spill-out, and for the gap size below 3 Å for Au and below 4Å for Na (i.e. with SF smaller than 3 and 4 respectively) the LEM completely disappears in the QHT spectrum, as shown in the Fig. 2(e)–2(f).

It is clear that for sodium, quantum mechanical effects are more pronounced compared to gold. This is due to two reasons: i) Na has a lower work function so that electrons are more relaxed at the surface and hence the spill-out region has a slightly larger extension; ii) Na does not have *d*-band electrons, which in Au systems counteract the effects of free electrons [49].

The results above are in a very good agreement with the previously published TD–DFT calculations presented in [37] for Au nanomatryoshkas. A more detailed comparison for the Na systems is reported in the next section.

#### 3.2. Comparison with TD–DFT results

TD–DFT calculations have been performed for filled-shell Na nanomatryoshka systems, corresponding to SF=1–4. The properties of considered systems are reported in Tab. 1. Note that the *N _{t}* increase very fast (cubically) with SF: due to some technical limitation of the current implementation (no fractional occupation numbers implemented) we couldn’t compute systems with SF larger than 4. We found that the filled-shell NM are very close to the target systems, with deviation in the radius less than 0.3Å.

In Fig. 3 we report (from top to bottom) the ground-state density *n*_{0}, the photoabsorption spectrum, the induced charge density *n*_{1} at the resonance peak, for SF=4 (left column) and SF=1 (right column).

As previously discussed for jellium spheres [27], the QHT density reproduces the DFT density quite accurately and the agreement increases for large systems, see Figs. 3(a) and 3(b). Clearly, (quantum) oscillations are neglected, due to the approximate KE functional: however, their amplitude decreases with the system size.

Concerning the main plasmon peak (HEM) in photoabsorption spectrum (see Figs. 3c and 3d) QHT is in very good agreement with the TD–DFT one, whereas LRA and TF–HT are significantly blue-shifted. We recall that in TD–DFT the line-shape originates from the convolution of many electronic transitions. This is quite evident for SF=1. For large systems, all electronic transitions collapse into a single plasmon peak. On the other hand, in the QHT approach, only one collective (plasmon) excitation is computed, which is broadened [28] with the viscoelastic term from Eq. (3). Concerning the LEM peak, both QHT and TD–DFT do not present any significant peak, in contrast to LRA and TF–HT, see insets in Figs. 3(c) and 3(d).

A more detailed comparison between QHT and TD–DFT is illustrated in Figs. 3(e)–3(f) (induced polarization density). Interestingly, the QHT approach reproduces very accurately the main oscillations at the external surface of the system as well as inside the gap. On the other hand, not all (quantum) oscillations can be reproduced by the QHT, which results averaged out inside the metal [see e.g. Figs. 3(e) and 3(f)].

#### 3.3. Field enhancements

Having validated the accuracy of the QHT approach, in this section we will discuss in more details the field enhancements, i.e. |**E*** _{tot}*|/|

**E**

*| = |*

_{inc}**E**

*|/*

_{tot}*E*

_{0}, which is together with the resonance shift, the most relevant feature for applications.

A comparison of the electric field distribution at LEM and HEM computed within LRA, TF–HT and QHT for Na NMs for SF=5 are shown in Fig. 4. It can be seen that the HEM produces a strong electric field inside the gap as well as outside the particle whereas at the LEM the induced field is confined only in the gap region.

The electric field enhancement inside the gap region, obtained by employing the LRA, TF–HT and QHT theories, is plotted as a function of distance in Fig. 5. Figures 5(a) and 5(b) report respectively results for Au and Na nanostructures for SF=5 at the resonance frequencies. The plot shows that LRA and TF–HT predict almost the same field enhancement with a linear profile inside the gap. On the other hand much lower field enhancement is predicted by the QHT approach and with a more parabolic profile, due to the spill-out effects.

Despite such more complicated profile, in the following we consider the field enhancement averaged in the gap, in order to have a global indicator to compare different methods and systems.

In Fig. 6 we report the electric field enhancement averaged inside the gap as a function of scaling factor (SF) at the lower and higher energy modes both for Au and Na. Note that for QHT spectra, LEM peaks result damped for small SFs. The field enhancement becomes then independent on the frequency within the associated region. The comparison shows that LRA and TF–HT predict similar field enhancements whereas much lower values are obtained from QHT due to quantum effects. More specifically, while LRA field enhancements are independent from the SF, TF–HT results show a weak size-dependence, converging towards the LRA results for larger SF. On the other had, QHT field enhancement increases quite rapidly with SF, for both modes.

It is worth noting that QHT field enhancement results are in good agreement with TD–DFT results, confirming the discussion of previous section. This is shown for Na NMs with SF=1–4 [see Fig. 6(b)]. Remark that TD–DFT LEM peaks cannot be identified for SF=3 or lower. This is true for QHT spectra as well, however in this case the field enhancement becomes independent on the frequency within the region associated to the LEM.

It is important to point out that the distribution of the electric-field enhancement (see Fig. 5) directly depends on spill-out effects. On the other hand, the absolute value of the electric-field enhancement (or its averaged values, reported in Fig. 6) depends on the absolute absorption efficiency at the resonance frequency, i.e. on the broadening of the plasmon resonance, which is treated differently by the different methods. In TD-DFT a myriad electronic transitions exist that translate into a broadening of the collective plasmon resonance. In QHT single-particle excitations cannot be exactly taken into account (since the basic assumption is that all the electrons lay in identical states) and the broadening is effectively described via a density-dependent viscous-like tensor [28]. In both cases, the broadening depends on the gap distance, contrarily to LRA and TF-HT in which the broadening only depends on the bulk electron collision rate.

Results in Fig. 6 are very important for applications. It is shown that TF–HT is even less accurate than LRA to described field enhancement in the gap. On the other hand, due to the importance of spill-out effects in these systems, QHT is the method of choice, being as accurate as TD–DFT.

#### 3.4. Optical response of NMs with D=60 nm

In this section, we consider Au and Na NM structures with a fixed particle size of diameter *D* = 60 nm. Systems of such size cannot be handled by TD–DFT methods. In particular, we consider a system with dimensions (*R*_{1}, *R*_{2}, *R*_{3}) where *R*_{1} = 17 nm is the fixed radius of the solid core, *R*_{2} = *R*_{1} + *g* nm is the inner radius and *R*_{3} = 30 nm is the outer radius of the shell. The gap between the core and shell is represented by *g* and it will vary from 4 nm to 0.1 nm.

The large size of the particle with sub-nanometer separation makes the problem multiscale and its numerical implementation to study both microscopic as well as macroscopic details becomes computationally virtually impossible for a TD–DFT approach but still possible, although every demanding, for QHT (and TF–HT for structures for which no analytical solution exists). As we mentioned before, the 2.5D formulation makes this task much easier. The absorption spectra for Au NMs calculated using the LRA, HT-TF and QHT methods are plotted in Fig. 7 (left panel). The change in core-shell spacing strongly influences the LEM whereas the plasmon peaks of HEM remain roughly at the same energy as the overall size of the particle remains the same (*D* = 60 nm). Results based on LRA and TF–HT show that as the core-shell distance gets smaller, the amplitude of the LEM decreases but there is always a distinct and well-defined peak for each gap size. However, in agreement with the previous study conducted in section Sec. 3.1 [Fig. 2 (c)], the QHT predicts that the amplitude of the LEM diminishes as the gap between the core and shell closes and it completely disappears for the gaps smaller than 0.3 nm (3Å) due to quantum effects.

Na NM absorption spectra show a third hybridized resonant mode, resulting from the combination of core resonant mode and antibonding resonant mode of the shell (see right panel of Fig. 7).We see that when the core-shell gap closes the spectra follow a similar trend as shown in Fig. 2 (f). That is, the LEM vanishes for the gaps smaller than 0.4 nm (4 Å). The QHT then efficiently describes near- and far-field optical properties of microscopic as well as macroscopic plasmonic systems.

Finally, in Fig. 8, we plot the average electric field in the gap for the LEM and HEM both for Au and Na. The insets present the results normalized with the average field in LRA for the gap sizes from 5 to 1Å. The plots show similar trend for average field as it was seen for the case of scaled Au and Na NMs, given in Fig. 6. It is interesting to observe that quantum effects, and particularly quantum tunneling, have the same impact on *macroscopic* systems as it was observed in small systems. This means that in sub-nm gap regimes, single-particle excitations might still drive the system optical response even when a system of millions of electrons is considered.

## 4. Conclusions

We have investigated near- and far-field optical properties of many Au and Na core-shell spherical nanomatryoshka structures, where the concentric core and shell are separated by a vacuum gap. The impact of particle size and core-shell spacing on the optical response was studied within the LRA, TF–HT, and QHT approximations and compared with TD–DFT results. We have found that LRA and TF–HT methods could be quite inaccurate to describe the microscopic behavior of these system, as essential quantum effects like electron spill-out and tunneling are neglected. The QHT, on the other hand, efficiently predicts their optical response and shows a very good agreement with full quantum TD–DFT calculations, not only for the absorption spectrum but also for the induced polarization density and local field enhancements. In fact, the field enhancement in the gap region is greatly reduced, particularly for smaller gaps, as compared to LRA and TF–HT results. Due to the importance of spill-out effects in these systems, QHT can be as accurate as TD–DFT.

Finally, we also considered large NM structures (*D* = 60 nm), with subnanometer gaps, whose sizes make them computationally prohibitive for TD–DFT. We found that when the gap size is fairly large, LRA, TF–HT and QHT theories approximate the same response because of absence of charge transfer between core and shell. However, as the core-shell separation is sufficiently decreased, one has to take into account quantum effects in order to correctly account for the correct response. The QHT can efficiently describe quantum spill-out effects as well as the plasmon resonances and near-field properties, thus providing an excellent tool to probe microscopic and macroscopic optical response of multiscale plasmonic systems.

## Funding

Air Force Office of Scientific Research (AFOSR) (FA9550-17-1-0177).

## References and links

**1. **J. J. Mock, R. T. Hill, Y.-J. Tsai, A. Chilkoti, and D. R. Smith, “Probing dynamically tunable localized surface plasmon resonances of film-coupled nanoparticles by evanescent wave excitation,” Nano Lett. **124**, 1757–1764 (2012). [CrossRef]

**2. **R. F. Oulton, V. J. Sorger, T. Zentgraf, R-M. Ma, C. Gladden, L. Dai, G. Bartal, and X. Zhang, “Plasmon lasers at deep subwavelength scale,” Nature **461**, 629–632 (2009). [CrossRef] [PubMed]

**3. **M. Moskovits, “Surface-enhanced spectroscopy,” Rev. Mod. Phys. **57**(3), 783–826 (1985). [CrossRef]

**4. **M. Scalora, M. A. Vincenti, D. de Ceglia, V. Roppo, M. Centini, N. Akozbek, and M. J. Bloemer, “Second- and third-harmonic generation in metal-based structures,” Phys. Rev. A **82**, 043828 (2010). [CrossRef]

**5. **E. Hao and G. C. Schatz, “Electromagnetic fields around silver nanoparticles and dimers,” J. Chem. Phys. **120**(1), 357–366 (2004). [CrossRef] [PubMed]

**6. **I. Romero, J. Aizpurua, G.W. Bryant, and F. J. García de Abajo, “Plasmons in nearly touching metallic nanoparticles: singular response in the limit of touching dimers,” Opt. Express **14**(21), 9988–9999 (2006). [CrossRef] [PubMed]

**7. **P. K. Jain and M. A. El-Sayed, “Plasmonic coupling in noble metal nanostructures,” Chem. Phys. Lett. **487**, 153–164 (2010). [CrossRef]

**8. **C. Ciracì, X. Chen, J. J. Mock, F. McGuire, X. Liu, S.-H. Oh, and D. R. Smith, “Film-coupled nanoparticles by atomic layer deposition: Comparison with organic spacing layers,” Appl. Phys. Lett. **104**, 023109 (2014). [CrossRef]

**9. **C. Ciracì, R. T. Hill, J. J. Mock, Y. A. Urzhumov, A. I. Fernandez-Dominguez, S. A. Maier, J. B. Pendry, A. Chilkoti, and D. R. Smith, “Probing the ultimate
limits of plasmonic enhancement,”
Science **337**,
1072–1074
(2012). [CrossRef] [PubMed]

**10. **K. J. Savage, M. M. Hawkeye, R. Esteban, and A. G. Borisov, “Revealing the quantum regime in tunnelling plasmonics,” Nature **491**, 574–577 (2012). [CrossRef] [PubMed]

**11. **J. A. Scholl, A. García-Etxarri, A. L. Koh, and J. A. Dionne, “Observation of quantum tunneling between two plasmonic nanoparticles,” Nano Lett. **13**(2), 564–569 (2013). [CrossRef]

**12. **S. Raza, N. Stenger, S. Kadkhodazadeh, S. V. Fischer, N. Kostesha, A-P. Jauho, A. Burrows, M. Wubs, and N. A. Mortensen, “Blueshift of the surface plasmon resonance in silver nanoparticles studied with EELS,” Nanophotonics **2**(2), 131–138 (2013). [CrossRef]

**13. **C. A. Ullrich, ed. *Time-dependent density functional theory: concepts and applications*, (Oxford University, 2011). [CrossRef]

**14. **A. Eguiluz, S. Ying, and J. Quinn, “Influence of the electron density profile on surface plasmons in a hydrodynamic model,” Phys. Rev. B **11**(6), 2118–2121 (1975). [CrossRef]

**15. **C. Schwartz and W. L. Schaich, “Hydrodynamic models of surface plasmons,” Phys. Rev. B **26**(12), 7008–7011 (1982). [CrossRef]

**16. **A. Boardman, *Electromagnetic surface modes hydrodynamic theory of plasmon-polaritons on plane surfaces*, (Wiley, 1982).

**17. **I. Tokatly and O. Pankratov, “Hydrodynamic theory of an electron gas,” Phys. Rev. B **60**(23), 15550–15553 (1999). [CrossRef]

**18. **A. Domps, P. G. Reinhard, and E. Suraud, “Time-dependent Thomas–Fermi approach for electron dynamics in metal clusters,” Phys. Rev. Lett. **80**(25), 5520–5523 (1988). [CrossRef]

**19. **S. Raza, G. Toscano, A-P. Jauho, M. Wubs, and N. A. Mortensen, “Unusual resonances in nanoplasmonic structures due to nonlocal response,” Phys. Rev. B **84**, 121412 (2011). [CrossRef]

**20. **G. Toscano, S. Raza, A-P. Jauho, N. A. Mortensen, and M. Wubs, “Modified field enhancement and extinction by plasmonic nanowire dimers due to nonlocal response,” Opt. Express **20**(4), 4176–4188 (2012). [CrossRef] [PubMed]

**21. **C. Ciracì, Y. A. Urzhumov, and D. R. Smith, “Effects of classical nonlocality on the optical response of three-dimensional plasmonic nanodimers,” J. Opt. Soc. Am. B **30**(10), 2731–2736 (2013). [CrossRef]

**22. **C. Ciracì, J. B. Pendry, and D. R. Smith, “Hydrodynamic model for plasmonics: a macroscopic approach to a microscopic problem,” Chem. Phys. Chem. **14**(6), 1109–1116 (2013). [CrossRef] [PubMed]

**23. **W. Yan, “Hydrodynamic theory for quantum plasmonics: linear-response dynamics of the inhomogeneous electron gas,” Phys. Rev. B **91**, 115416 (2015). [CrossRef]

**24. **E. Zaremba and H. C. Tso, “Thomas–Fermi–Dirac–von Weizsäcker hydrodynamics in parabolic wells,” Phys. Rev. B **49**(12), 8147–8162 (1994). [CrossRef]

**25. **G. Toscano, J. Straubel, A. Kwiatkowski, C. Rockstuhl, F. Evers, H. Xu, N. A. Mortensen, and M. Wubs, “Resonance shifts and spill-out effects in self-consistent hydrodynamic nanoplasmonics,” Nat. Commun. **6**, 7132 (2015). [CrossRef] [PubMed]

**26. **X. Li, H. Fang, X. Weng, L. Zhang, X. Dou, A. Yang, and X. Yuan, “Electronic spill-out induced spectral broadening in quantum hydrodynamic nanoplasmonics,” Opt. Express **23**(23), 29738–29745 (2015). [CrossRef] [PubMed]

**27. **C. Ciracì and F. D. Sala, “Quantum hydrodynamic theory for plasmonics: impact of the electron density tail,” Phys. Rev. B **93**(20), 205405 (2016). [CrossRef]

**28. **C. Ciracì, “Current-dependent potential for nonlocal absorption in quantum hydrodynamic theory,” Phys. Rev. B **95**(24), 245434 (2017). [CrossRef]

**29. **K. Ding and C. T. Chan, “Plasmonic modes of polygonal rods calculated using a quantum hydrodynamics method,” Phys. Rev. B **96**(12), 125134 (2017). [CrossRef]

**30. **D-K. Lim, K-S Jeon, J-H Hwang, H. Kim, Y. D. Suh, and J-M Nam, “Highly uniform and reproducible surface-enhanced Raman scattering from DNA-tailorable nanoparticles with 1-nm interior gap,” Nat. Nanotechnol. **6**, 452–460 (2011). [CrossRef] [PubMed]

**31. **L. Lin, M. Zapata, M. Xiong, Z. Liu, S. Wang, H. Xu, A. G. Borisov, H. Gu, P. Nordlander, J. Aizpurua, and J. Ye, “Nanooptics of Plasmonic Nanomatryoshkas: Shrinking the Size of a Core-Shell Junction to Subnanometer,” Nano Lett. **15**(10), 6419–6428 (2015). [CrossRef] [PubMed]

**32. **H. Xu, “Multilayered metal core-shell nanostructures for inducing a large and tunable optical field,” Phys. Rev. B **72**, 073405 (2005). [CrossRef]

**33. **Y. Hu, R. C. Flemming, and R. A. Drezek, “Optical properties of gold-silica-gold multilayer nanoshells,” Opt. Express **16**(24), 19579–19591 (2008). [CrossRef] [PubMed]

**34. **R. Bardhan, S. Mukherjee, N. A. Mirin, S. D. Levit, P. Nordlander, and N. J. Halas, “nanosphere-in-a-nanoshell: a simple nanomatryoshka,” J. Phys. Chem. C **114**(16), 7378–7383 (2010). [CrossRef]

**35. **J. Qian, W. Wang, Y. Li, J. Xu, and Q. Sun, “Optical extinction properties of perforated gold-silica-gold multilayer nanoshells,” J. Phys. Chem. C **116**(18), 10349–10355 (2012). [CrossRef]

**36. **C. Ayala-Orozco, J. G. Liu, M. W. Knight, Y. Wang, J. K. Day, P. Nordlander, and N. J. Halas, “Fluorescence enhancement of molecules inside a gold nanomatryoshka,” Nano Lett. **14**(5), 2926–2933 (2014). [CrossRef] [PubMed]

**37. **V. Kulkarni, E. Prodan, and P. Nordlander, “Quantum plasmonics: optical properties of a nanomatryoshka,” Nano Lett. **13**(12), 5873–5879 (2013). [CrossRef] [PubMed]

**38. **M. Zapata, A. S. C. Beltrań, A. G. Borisov, and J. Aizpurua, “Quantum effects in the optical response of extended plasmonic gaps: validation of the quantum corrected model in core-shell nanomatryoshkas,” Opt. Express **23**(6), 8134–8149 (2015). [CrossRef] [PubMed]

**39. **D-C. Marinica, J. Aizpurua, and A. G. Borisov, “Quantum effects in the plasmon response of bimetallic core-shell nanostructures,” Opt. Express **24**(21), 23941–23956 (2016). [CrossRef] [PubMed]

**40. **C. Ciracì, Y. A. Urzhumov, and D. R. Smith, “Far-field analysis of axially symmetric three-dimensional directional cloaks,” Opt. Express **21**(8), 9397–9406 (2013). [CrossRef] [PubMed]

**41. **W. Ekardt, “Size-dependent photoabsorption and photoemission of small metal particles,” Phys. Rev. B **31**(10), 6360–6370 (1985). [CrossRef]

**42. **M. Brack, “The physics of simple metal clusters: self-consistent jellium model and semiclassical approaches,” Rev. Mod. Phys. **65**(3), 677–732 (1993). [CrossRef]

**43. **S. Raza, S. I. Bozhevolnyi, M. Wubs, and N. A. Mortensen, “Nonlocal optical response in metallic nanostructures,” J. Phys.: Condens. Matter **27**, 183204 (2015).

**44. ** Comsol Multiphysics, http://www.comsol.com.

**45. **G. Bertsch, “An RPA program for jellium spheres,” Comput. Phys. Commun. **60**(2), 247–255 (1990). [CrossRef]

**46. **E. Prodan and P. Nordlander, “Electronic structure and polarizability of metallic nanoshells,” Chem. Phys. Lett. **352**(3), 140–146 (2002). [CrossRef]

**47. **B. K. P. Scaife, *Principles of dielectrics*, (Oxford Science, 1998).

**48. **J. Zuloaga, E. Prodan, and P. Nordlander, “Quantum description of the plasmon resonances of a nanoparticle dimer,” Nano Lett. **9**(2), 887–891 (2009). [CrossRef] [PubMed]

**49. **T. V. Teperik, P. Nordlander, J. Aizpurua, and A. G. Borisov, “Robust subnanometric plasmon ruler by rescaling of the nonlocal optical response,” Phys. Rev. Lett. **110**(26), 263901 (2013). [CrossRef] [PubMed]