Ultrasmall metal-insulator-metal nanoresonators: impact of slow-wave effects on the quality factor

We study the quality factor variation of three-dimensional Metal-Insulator-Metal nanoresonators when their volume is shrunk from the diffraction limit (λ/2n)3 down to a deep subwavelength scale (λ/50)3. In addition to rigorous fully-vectorial calculations, we provide a semi-analytical expression of the quality factor Q obtained with a Fabry-Perot model. The latter quantitatively predicts the absorption and radiation losses of the nanoresonator and provides an in-depth understanding of the mode lifetime that cannot be obtained with brute-force computations. In particular, it highlights the impact of slow-wave effects on the Q-factor as the size of the resonator is decreased. The Fabry-Perot model also evidences that, unexpectedly, wave retardation effects are present in metallic nanoparticles, even for deep subwavelength dimensions in the quasi-static regime.

In view of potential applications of MIM resonators in a variety of areas ranging from metamaterials to photovoltaics, it is important to obtain a deep understanding of their optical properties, in particular when their size is scaled down to deep subwavelength dimensions. Indeed, ultra-small resonators are of major importance both for metamaterials, where metaatoms should be much smaller than the wavelength, and for photodetector applications, where several resonators can be associated within a subwavelength cell to engineer multi-resonant structures [24]. The tunability of the resonance wavelength of MIM resonators with the resonator dimensions being well-known [11,12,[25][26][27], we will mainly focus on the analysis of the quality factor (Q-factor) variation as the resonator size is scaled down.
We study the resonance with the lowest energy supported by a single 3D MIM resonator, see Fig. 1. This fundamental resonance is crucial in the context of ultra-small resonators because it has no cut-off and can be scaled down to deep subwavelength dimensions in the quasi-static regime. Moreoever, this resonance presents a "magnetic" character and can be used to engineer negative-index metamaterials [14,15,28,29]. We especially quantify the variation of the Q-factor when the resonator volume is shrunk from the diffraction limit V = (λ /2n) 3 down to a deep subwavelength scale V = (λ /50) 3 . As shown by 3D fully-vectorial calculations, the main trend is a significant increase of the Q-factor by one order of magnitude when the size is reduced. To explain this increase, we use an approximate Fabry-Perot model, which allows us to derive semi-analytical expressions of the quality factor, of the absorption and of the radiation losses. The model thus provides a real understanding of the physics governing the mode lifetime at deep subwavelength scales. In particular, it evidences the crucial role played by the slowness of the SPP modes that are bouncing back and forth inside the resonator, a physical effect that is completely hidden in brute-force computations. Over the past ten years, Fabry-Perot models have been successfully applied to a variety of micro-and nanoresonators, including photonic-crystal microcavities and micropilars [30,31], and plasmonic nanowire resonators with a longitudinal length of a few wavelengths [32][33][34][35][36]. Hereafter we evidence that such a useful picture remains valid and helpful even for analyzing nanoresonators that are operating far below the diffraction limit in the quasi-static regime. We note that a Fabry-Perot model has already been used for studying two-dimensional (2D) MIM resonators [26]. However, the authors in [26] overlooked the importance of slow-wave effects and consequently derived a Q-factor expression that is not consistent with the asymptotic value derived in the quasi-static limit for metallic nanoparticles of arbitrary shape [37]. We correct this discrepancy by properly taking into account slow-wave effects in the Fabry-Perot model. Indeed, as we show in Sections 3 and 4, a correct Fabry-Perot expression for the Q-factor is highly accurate and fully consistent with other asymptotic expressions derived in the quasistatic limit [37]. In Section 3 we highlight three different regimes for the Q-factor variation with the resonator size. As the volume shrinks from the diffraction limit (λ /2n) 3 to (λ /20) 3 , the Q-factor first increases because of a reduction of the radiation losses. Then, from (λ /20) 3 to (λ /35) 3 , the radiation losses reduction is balanced by an increase of the absorption and one would intuitively expect a drop of the Q-factor as predicted in [26]. However, the Q-factor keeps on increasing because of the slowdown of the SPP mode that bounces back and forth inside the resonator. This slow-wave effect, which is similar to the one classically encountered in photonic crystal microcavities [30,31], results in an increase of the cavity-mode lifetime. Finally, for ultrasmall resonators, the Fabry-Perot model correctly predicts the quality factor saturation toward the asymptotic value obtained in the quasi-static limit [37].

Magnetic resonance of a single 3D MIM resonator
The 3D MIM resonator under study consists of a dielectric rectangular nanoparticle (width w, length L and thickness t d ) sandwiched between a metal substrate and a metal layer with the same width and length and with a thickness denoted by t m , see Fig. 1(a). Calculations have been performed for silver, whose permittivity has been taken from the data tabulated in [38]. In order to maximize the confinement and the Q-factor, the dielectric material has been chosen to be a semiconductor (such as GaAs) with a high refractive index (n = 3.5). Since it has only a weak impact on the resonator properties as long as it remains larger than the skin depth, the thickness t m of the top metallic layer is fixed in the following and we take t m = 25 nm. All rigorous fully-vectorial calculations are performed with a 3D frequency-domain fully-vectorial modal method known as the aperiodic Fourier Modal Method (a-FMM) [39]. Stretching of the numerical space through coordinate transforms is additionally incorporated for improving the computational accuracy, see the method MM3 in the benchmark article [40] for more details.
We study the resonance with the lowest energy supported by the MIM resonator -the fundamental mode. Figure 1(b) presents the spectrum of the electric field enhancement inside the MIM resonator (dimensions are provided in the figure caption) when the latter is illuminated from air by a plane wave at normal incidence and polarized along the z-direction. Figure 1(c) shows the distribution of the magnetic field |H x | 2 and of the electric current density J(r) induced in the nanoresonator by the incident plane wave at resonance, λ = 947 nm. From the curl Maxwell equations, the induced current is proportional to the total electric field E(r), with ω the frequency, ε 0 the vacuum permittivity, ε(r) the relative permittivity distribution that defines the resonator and ε ref (r) the relative permittivity of a reference background. We have chosen a uniform background ε ref (r) = 1, which implies that the field radiated by the current density J(r) placed in air includes both the reflection of the incident plane wave by the metallic plane and the scattering by the resonator. As can be seen from Fig. 1(c), at the resonance wavelength, the electric current density J(r) induced by the incident plane wave is essentially localized around the resonator and forms a loop in the (y, z) plane, Re(J x ) Re(J y ) and Re(J x ) Re(J z ). Moreover, the resonant magnetic field is perpendicular to the loop. This highlights the magnetic character of the fundamental MIM resonance, as previously discussed in the literature [14,15,28,29]. In the following, we thoroughly study the quality factor Q of this resonance. Note that, throughout the report, the resonance wavelength λ 0 = 950 nm is maintained constant. In order to fulfill this condition as the two transverse resonator dimensions (the width w and the dielectric thickness t d ) are varied, one simply adjusts the length L of the MIM particle. The quality factor can be directly derived from the full-width-at-half-maximum of a Lorentzian fit of the intensity enhancement spectrum shown in Fig. 1(b). However, in order to get more physical insight, we analyze the MIM resonator as a Fabry-Perot cavity.

Fabry-Perot model of the resonance
Within the Fabry-Perot picture, the resonance is described as a standing wave pattern along the z-direction created by the bouncing between the nanoparticle facets at z = −L/2 and z = L/2 of a single mode, the fundamental plasmonic mode of the MIM waveguide, whose cross-section in the (x, y) plane is represented in Fig. 2. The single mode approximation amounts to assume that all higher-order modes of the MIM waveguide play a negligible role to build the resonance. The validity of the assumption is discussed in more details in the following; as will be shown, it depends mostly on the transversal size of the MIM waveguide. The key parameters of the Fabry-Perot model are the propagation constant β = k 0 n eff + i α 2 of the fundamental plasmonic mode supported by the MIM waveguide and its complex reflection coefficient r = √ R exp(iφ r ) at the air/MIM interface.
We first calculate and discuss these two important physical quantities in Section 3.1. Then we provide in Section 3.2 analytical expressions for the resonance length L and the quality factor Q. These formulae are validated against fully-vectorial calculations and the gained physical understanding for the Q-factor variation with the resonator volume is discussed in Section 3.3.

Fabry-Perot parameters: fundamental mode of the MIM waveguide and its reflectivity
The three main field components of the fundamental mode of the MIM waveguide, |E y | 2 , |E z | 2 and |H x | 2 , are shown in Figs. 2(a)-2(c) for w = 40 nm and t d = 20 nm. For w → ∞ (planar MIM stack) and t m → ∞, this mode results from the coupling with a symmetric magnetic field H x of two metal/insulator SPPs [41]. The H x -symmetric mode of a planar MIM stack has no cut-off as the insulator thickness decreases, unlike the antisymmetric mode that is cut-off at roughly λ /(2n). Therefore, for the thicknesses of interest, t d < 100 nm, the planar MIM stack is monomode. Reducing the top metal thickness to t m = 25 nm does not change the nature of the fundamental mode, except that the symmetry is slightly broken, see Figs. 2(a)-2(c). As the dielectric thickness t d decreases, the symmetric mode of the MIM stack interacts more strongly with the metal and presents a rapid increase of both its effective index n eff = Re(β )/k 0 and its attenuation α = 2Im(β ), as shown by the thin-solid curves in Figs. 2(d)-2(e). The fundamental mode of the MIM waveguide with a finite width follows a similar trend. As can be seen from Figs. 2(d)-2(e), the impact of the width w on the propagation constant is much weaker than the impact of the thickness t d ; the propagation constant of the 2D waveguide is almost equal to that of the planar stack, especially for small thickness. Let us emphasize that the mode group index n g = c/v g , with v g the group velocity, also increases as the dielectric thickness is reduced, see Fig. 2(f). As will be shown hereafter, the slowdown of the plasmonic mode supported by thin MIM waveguides directly impacts the resonator lifetime. Finally, it is noteworthy that for finite widths higher-order modes with the same profile along y as the fundamental plasmonic mode may exist. They are determined by a total-internal-reflection condition at the air/dielectric interfaces and possess an increasing number of nodes along x, see the insets in Fig. 3(c).
The large increase of the effective index is a key ingredient to build ultra-small resonators since the resonator length L corresponding to the fundamental resonance is of the order of λ 0 /(2n eff ). However, the price to pay for a reduction of the resonator volume is an enhanced dissipation in the metal. In order to accurately quantify the impact on the quality factor of the absorption increase, we need to thoroughly evaluate the variation of the radiation losses with the size. In the Fabry-Perot picture, the normalized power dissipated by radiation per roundtrip is simply equal to 1 − R, with R the reflectivity of the fundamental mode at the air/MIM interfaces, z = −L/2 and z = L/2. We have performed fully-vectorial calculations of the modal reflection coefficient r = √ R exp(iφ r ) as a function of w and t d with the a-FMM [39,40]. The results are presented in Fig. 3(a). The reflectivity strongly increases as the waveguide cross-section is reduced because of the increased mode confinement. For w = t d = 100 nm, the facets only weakly reflect the fundamental plasmonic mode, R = 0.7, whereas for w = 40 nm and t d = 10 nm R is as large as 0.999. Figure 3(b) shows the dependence of the reflectivity with the thickness t d for w = 40, 100, 350 nm and ∞. The four curves show the same main trend, namely a significant enhancement of R as t d is reduced. The width dependence is highlighted in Fig. 3(c) for t d = 20, 35 and 50 nm. It is important to note the non-monotonous variation of the reflectivity. The narrow dips correspond to a strong backscattering into higher-order MIM waveguide modes that are cutoff at the dip wavelength. The insets in Fig. 3(c) show the profile of the first two x-symmetric higher-order modes [42]. When a higher-order mode passes its cut-off and becomes propagative, the total density of states is modified by the appearance of a new mode (often with a small group velocity) and this causes a sharp dip in the modal reflectivity of the fundamental mode.
The dip is especially pronounced for the appearance of the first higher-order mode. Finally, it is noteworthy that even extremely tiny waveguides with transverse dimensions as small as w < 100 nm and t d < 10 nm are not monomode. This unexpected observation can be intuitively understood from the large value of the effective index n 1D eff of planar MIM stacks (w infinite) with a small dielectric thickness t d ; the cut-off width of higher-order modes is indeed linked to the value of λ /(2n 1D eff ), which becomes smaller and smaller as t d is decreased.

Fabry-Perot equations: phase-matching condition and quality factor
Since we have fixed the resonance wavelength to λ 0 = 950 nm, the length L of the resonator is easily obtained from the phase-matching condition with m = 1, 2, . . ., n eff the effective index of the waveguide mode given in Fig. 2(d) and φ r the phase of the modal reflection coefficient r. We are interested hereafter in the fundamental resonance m = 1. The reflection phase varies between φ r = π/8 for the largest waveguide cross-section (w = 350 nm and t d = 100 nm) and φ r = −π/20 for the smallest cross-section (w = 40 nm and t d = 5 nm) and therefore, with a good approximation, one may consider that Eq. (2) leads to L = λ 0 /(2n eff ), corresponding to the classical half-wavelength condition for the resonator length [43]. Figure 4(a) shows the resonator length obtained from Eq. (2) as a function of the dielectric thickness t d for w = 40, 100 and 350 nm. Because of the strong increase of n eff , L rapidly decreases as t d decreases. On the other hand, the resonator width w is found to weakly impact the length. Within the Fabry-Perot picture, and under the assumption of a narrow resonance, the quality factor Q can be derived analytically as [31] with k 0 = 2π/λ 0 , n g = n eff − λ ∂ n eff ∂ λ the group index of the MIM waveguide mode, R eff = R exp(−αL) the effective reflectivity that includes the absorption loss over one-half roundtrip, and L eff = L + 2L p the effective resonator length that includes the penetration length L p = − λ 2 4πn g ∂ φ r ∂ λ , which is due to the dispersive nature of the reflection. The penetration length represents less than 10% of the effective length L eff and it can be safely neglected for all values of w and t d considered in the present work [44].
The predictions of Eq. (3) are shown in Fig. 4(b) for w = 40, 100 and 350 nm together with fully-vectorial calculations of the Q-factor, extracted from Lorentzian fits of the intensity enhancement spectra and shown with various markers. The analytical formula of Eq. (3) accurately predicts the Q-factor increase. When the resonator volume is shrunk from roughly the diffraction limit V = (λ /2n) 3 (w = 120 nm, t d = 100 nm and L = 100 nm) down to a deep subwavelength scale V = (λ /50) 3 (w = 40 nm, t d = 5 nm and L = 30 nm), the Q-factor is significantly enhanced by one order of magnitude.
Strictly speaking, the Fabry-Perot model should be applied only when all higher-order modes are evanescent. For instance, for w = 100 nm, the MIM is monomode only for t d > 15 nm, and for w = 350 nm, the MIM is multimode whatever the value of t d . Hopefully, the presence of higher-order modes does not alter the model predictions, except close to their cut-off where a drop of Q is inaccurately predicted in Fig. 4(b) because of the rapid change of the facet reflectivity. Note that the deviation between the model predictions and rigorous calculations vanishes as the waveguide width increases.

Analysis of the Q-factor increase
Because of its analytical treatment, the model allows us to provide a comprehensive analysis of the Q-factor limitations. As shown by Eq. (3), the cavity-mode lifetime is impacted both by absorption (α in R eff ) and radiation (1 − R). In order to separate the contributions of these two loss channels and to clarify the impact of the SPP slowdown, Fig. 4(c) shows the Q-factor predicted by Eq. (3) for w = 40 nm (thick solid blue curve) together with the Q-factor values that would have been achieved for α = 0, R = 1 and n g = 5. Taking α = 0 in Eq. (3) amounts to remove the absorption; the black dashed-dotted curve thus gives an approximation of the radiation-limited Q. Then, considering the particle facets as perfect reflectors, R = 1 in Eq. (3), amounts to neglect the radiation losses and gives an approximation of the absorption-limited Q.
Consistently with the fact that the absorption and the scattering of small particles respectively scale as the particle volume and as the particle volume squared [45], we find that the Q-factor of large MIM resonators with V = (λ /2n) 3 is mostly limited by radiation, Q ≈ Q α=0 , whereas the Q-factor of small resonators in the quasi-static limit with V = (λ /50) 3 is completely limited by absorption, Q ≈ Q R=1 . Radiation and absorption losses are balanced for t d = 35 nm. Finally, the impact on the resonance lifetime of the slowdown of the MIM waveguide mode is unraveled by considering Eq. (3) with a constant group index n g = 5, a value corresponding to thick MIM waveguides with t d = 100 nm (red dashed curve). Three different regimes of Q-factor variation can be considered. They are marked by vertical dashed lines in Fig. 4(c). The first regime concerns large resonator thicknesses, t d > 35 nm. The absorption αL and the group index n g are roughly constant and the increase of the quality factor is purely due to a reduction of the radiation losses (increase of the reflectivity R). In the second regime, 10 < t d < 35 nm, the Q-factor increases (solid thick blue curve) due to the slowdown of the SPP mode bouncing inside the resonator. Indeed, it is noticeable that a resonator with the actual absorption and radiation but with a constant group index (red dashed curve) would see its Q-factor decreasing because the radiation-loss reduction is balanced by an increase of the absorption. Finally, for extremely thin dielectric layers, t d < 10 nm, we enter the quasi-static regime where the resonator performance is purely limited by absorption. The validity of the Fabry-Perot model in this last regime is discussed in the next Section.

Fabry-Perot model in the quasi-static limit
In the quasi-static limit, i.e., for dimensions much smaller than the wavelength, the quality factor Q s of any metallo-dielectric resonant nanoparticle composed of a single lossy metal is completely determined from the relative permittivity of the metal ε m = ε m + iε m [37], In this limit, the specific structure of the particle is not important and the Q-factor depends neither on the geometric shape nor on the dielectric media, provided that the latter are dispersionless and lossless. This universal result has been derived in [37] by assuming a nanoparticle built in a noble metal (|ε m | ε m ) that supports a purely electrostatic resonance (∇ × E = 0) whose damping is solely due to absorption (no radiation losses). This last assumption can be justified by the fact that the scattering of small particles scales as the particle volume squared whereas the absorption scales only as the particle volume [45]. The quasi-static regime is usually associated to the absence of wave retardation effects because k 0 L 1 for small size-to-wavelength ratios. On the other hand, the Fabry-Perot model describes a resonance as the result of a round-trip phase accumulation of 2mπ, m = 1, 2, . . . [43]. Thus, at first sight, the Fabry-Perot picture seems inappropriate to describe a resonance in the quasi-static limit. This intuitive conclusion is in contradiction with Fig. 4(b) that evidences the accuracy of the Fabry-Perot predictions even for ultrasmall resonators with V = (λ /50) 3 . In order to clarify whether a nanoparticle can be described by a non-vanishing wave retardation in the quasi-static limit, we now analytically study the asymptotic behavior of the Fabry-Perot model as resonator sizes tend toward zero. We emphasize that retardation effects are present in very small objects and that MIM nanoparticles, even in the quasi-static limit, are equivalent to half-wavelength antennas. This equivalence has been recently discussed for the resonance wavelength [46], and hereafter, we additionally evidence that it applies also quantitatively to the quality factor of the resonance. In particular, we show that the closed-form expression of the Q-factor given by Eq. (3) reduces to Eq. (4) in the quasi-static limit.
We consider a Fabry-Perot resonance with a Q-factor given by Eq. (3). As in [37], we assume that the resonance damping in the quasi-static limit is purely due to absorption. This amounts to consider R = 1 in Eq. (3), an assumption legitimated for the geometry under study by the numerical results of Fig. 4(c). Moreover, the resonator length is small and the absorption term in the Q-factor can be approximated by exp(−αL) ∼ 1 − αL. With these two approximations, Eq. (3) leads to To obtain this result, we have considered that L p L and L eff ≈ L [44]. Equation (5) emphasizes that, in the absence of radiation losses, the Q-factor is purely driven by the guided SPP mode bouncing back and forth inside the particle. The important parameter is the group-index-toattenuation ratio, and not only the attenuation as concluded in [26]. Next we examine the asymptotic value of the propagation constant of MIM waveguides as their transverse dimensions become much smaller than the wavelength. Actually, the limit as t d /λ → 0 of the normalized propagation constant n 1D eff of planar MIM stacks is known analytically [26] with ε d the relative permittivity of the dielectric layer. We note from Fig. 2 that the normalized propagation constant of the fundamental mode of 2D MIM waveguides with finite widths w becomes independent of the width as t d /λ → 0 and asymptotically tends toward n 1D eff . This is understood by realizing that for any fixed value of the width w, we can have w λ /(2n 1D eff ) as t d is reduced. Thus Eq. (6) can be safely used to derive the limit of the n g /α ratio for any 2D MIM waveguide. The expressions of n g = n 1D eff − λ ∂ n 1D eff ∂ λ and α = 2k 0 Im(n 1D eff ) are easily derived from Eq. (6), and by inserting them into Eq. (5), the asymptotic value Q s FP of the Fabry-Perot Q-factor in the quasi-static limit is found to be given by which is exactly the Q s value of Eq. (4) derived in [37]. Note that to derive Eq. (7), we assumed that |ε m | ε m . The prediction of Eq. (7) is represented by the horizontal arrow in Fig. 4(b). It is noteworthy that the 3D MIM nanoparticle under study enters the quasi-static regime only for very small dielectric thicknesses, t d < 10 nm. The derivation of Eq. (7) leads to two important conclusions. First, it evidences that the Fabry-Perot resonator model can be appropriately used to describe a resonance in the quasistatic limit. In other words, the quasi-static regime does include wave retardation effects linked to the propagation inside the particle of plasmonic modes whose wavevector diverges as the length shrinks (actually the product k 0 n eff L = π is fixed by the Fabry-Perot phase-matching condition). The second important conclusion concerns the properties of plasmonic waveguides in the quasi-static limit. Indeed, the Q-factor expressions of Eqs. (4) and (5) should be identical as the transverse cross-section shrinks, whatever the shape of the waveguide cross-section. We can therefore deduce that, in the limit of small transverse dimensions compared to the wavelength, the propagation constant of any plasmonic waveguide should satisfy regardless of the geometric shape of the waveguide cross-section or of the dielectric media composing it, provided that they are dispersionless and lossless. Note that Eq. (8) holds in a regime where the scattering of the plasmonic mode at an interface can be neglected. The waveguide dimensions corresponding to this regime can be different from one waveguide to the other [47].

Conclusion
In summary, we have studied the quality factor of the fundamental resonance supported by MIM nanoresonators when their volume is shrunk from the diffraction limit V = (λ /2n) 3 down to a deep subwavelength scale V = (λ /50) 3 . The ten-fold increase of the Q-factor given by rigorous fully-vectorial calculations has been accurately predicted with a Fabry-Perot model, which provides analytical expressions of the quality factor, of the absorption and of the radiation losses. These expressions allow for a comprehensive analysis of the increase of the mode lifetime as the dimensions of the resonator are shrunk. The main mechanisms are a reduction of the radiation losses and a slowdown of the plasmonic mode bouncing back and forth inside the resonator.
The slow-wave effect is responsible for the Q-factor increase when the radiation-loss decrease is fully balanced by an increase of the absorption. A second important result of this work is to show that the Fabry-Perot model remains quantitatively valid down to very small dimensions far below the diffraction limit. This evidences that the localized plasmon resonances of MIM nanoparticles in the quasi-static limit can be quantitatively analyzed with the same wave retardation effects as the delocalized resonances of nanowires with a length of a few wavelengths. The nanoparticle is indeed an half-wavelength antenna but with a vanishing effective wavelength λ 0 /(2n eff ) [33].
We hope that the Fabry-Perot analytical treatment of the magnetic resonance supported by MIM nanoparticles will be helpful for further "bottom-up" approaches for the design of optical metamaterials starting from the optical properties of single meta-atoms.