Modal analysis for nanoplasmonics with nonlocal material properties

TL;DR

This paper develops contour integral methods to compute eigenfrequencies and resonant states in nanoplasmonic systems with nonlocal material properties, enabling detailed modal analysis of nanowires.

Contribution

It introduces a novel approach for modal analysis in nonlocal plasmonic systems using contour integral techniques, advancing the understanding of nanoplasmonic resonances.

Findings

Computed all complex eigenfrequencies for a nanowire across a broad frequency range.

Identified resonant states relevant for the nanowire's extinction cross section.

Demonstrated the effectiveness of the method for nonlocal hydrodynamic Drude models.

Abstract

Plasmonic devices with feature sizes of a few nanometers exhibit effects which can be described by the nonlocal hydrodynamic Drude model. We demonstrate how to exploit contour integral methods for computing eigenfrequencies and resonant states of such systems. We propose an approach for deriving the modal expansion of relevant physical observables. We use the methods to perform a modal analysis for a metal nanowire. All complex eigenfrequencies in a large frequency range and the corresponding resonant states are computed. We identify those resonant states which are relevant for the extinction cross section of the nanowire.

Equations

$$\nabla \times$$

$$=i\omega {\mathbf{J}}_{\mathrm{hd}}(\mathbf{r},\omega )+i\omega \mathbf{J}(\mathbf{r},\omega ),$$

$${\beta }^2\nabla$$

$$=i\omega {\omega }_{\mathrm{p}}^2{ϵ}_0\mathbf{E}(\mathbf{r},\omega ),$$

$$\mathbf{E}(\mathbf{r},{\omega }_0)=\frac{1}{2\pi i}\underset{C_0}{\oint }\frac{\mathbf{E}(\mathbf{r},\omega )}{\omega -{\omega }_0}\text{ d}\omega ,$$

$$\mathbf{E}(\mathbf{r},{\omega }_0)=$$

$$+\frac{1}{2\pi i}\underset{C_{\text{nr}}}{\oint }\frac{\mathbf{E}(\mathbf{r},\omega )}{\omega -{\omega }_0}\text{ d}\omega ,$$

$${\tilde{\mathbf{E}}}_k(\mathbf{r},{\omega }_0)=-\frac{1}{2\pi i}\underset{{\tilde{C}}_k}{\oint }\frac{\mathbf{E}(\mathbf{r},\omega )}{\omega -{\omega }_0}\text{ d}\omega ,$$

$${\mathbf{E}}_{\text{nr}}(\mathbf{r},{\omega }_0)=\frac{1}{2\pi i}\underset{C_{\text{nr}}}{\oint }\frac{\mathbf{E}(\mathbf{r},\omega )}{\omega -{\omega }_0}\text{ d}\omega ,$$

$$T(\omega )v=f(\omega ),$$

$$v^{\circ }(\omega )=T^{\circ }(\omega {)}^{-1}{f}^{\circ }(\omega ).$$

$$\sigma ({\omega }_0)=\frac{1}{P_{\mathrm{pw}}}$$

$$+{\int }_{{\Omega }_{\mathrm{nw}}}\frac{1}{2}\mathrm{Re}\left({\mathbf{E}}^{\ast }(\mathbf{r},{\omega }_0)\cdot {\mathbf{J}}_{\mathrm{hd}}(\mathbf{r},{\omega }_0)\right)\mathrm{d}V],$$

Full text

I Introduction

Nanofabrication technologies allow for a rapid progress in engineering nano-optical devices Lindquist et al. (2012). Plasmonic resonances are the center of attention for many topical applications exploring new regimes of physics. Examples comprise the demonstration of plasmonic lasers Oulton et al. (2009), tailoring light emission of nanoantennas Curto et al. (2010); Giannini et al. (2011), probing single molecules and nanoparticles by Raman scattering Nie and Emory (1997), plasmonic photochemistry Zhang et al. (2018), and quantum emitters interacting with metal nanoresonators Chikkaraddy et al. (2016).

An adequate description of material dispersion plays an important role for the investigation of light-matter interaction in plasmonic structures Murray and Barnes (2007). In many cases, the material dispersion can be described by the Drude-Lorentz model or by a rational function fit to measured material data Sehmi et al. (2017); Garcia-Vergara et al. (2017). Such models are based on spatially local interactions between the light and the free electron gas of the plasmonic scatterers Johnson and Christy (1972). When the scatterers are at the size of a few nanometers, nonlocal material models are required Raza et al. (2015). These models lead to additional resonances of the electromagnetic field with sub-nm wavelengths. Recently, surface plasmon resonance blueshifts have been observed in metal nanoparticles Scholl et al. (2012); Raza et al. (2013) which could be explained Christensen et al. (2014) using the nonlocal hydrodynamic Drude model (HDM) Boardman (1982). This model assumes that the motion of the electron gas behaves as a hydrodynamic flow and allows for the investigation of nonlocal physical effects Ruppin (2001); Palomba et al. (2008); García de Abajo (2008); Intravaia and Busch (2015); Toscano et al. (2015); Schnitzer et al. (2016); Moeferdt et al. (2018).

For the study of physical phenomena in nanoplasmonic systems, a deeper understanding of the effects based on the HDM is required. A modal description is the most instructive approach Kristensen and Hughes (2014); Lalanne et al. (2018). In the case of local material models, numerically computed resonant states of plasmonic systems have been successfully used to derive modal expansions Sauvan et al. (2013); Kamandar Dezfouli et al. (2017a); Zschiedrich et al. (2018); Yan et al. (2018). However, in the case of the HDM, a coupled system of equations has to be solved McMahon et al. (2009); Hiremath et al. (2012); Toscano et al. (2012); Kamandar Dezfouli et al. (2017b). To the best of our knowledge, for this system, the computations of eigenfrequencies in a large frequency range with corresponding resonant states and modal expansions have not yet been reported.

In this work, we investigate plasmonic resonances based on the HDM. We present a contour-integral-based framework for a modal analysis. Typical physical observables are sesquilinear forms which involve a complex conjugation of the solution fields. We propose a general approach for the computation of modal sesquilinear quantities. The framework is applied to calculate the eigenfrequencies and corresponding resonant states of a metal nanowire. Furthermore, the modal extinction cross section of the nanowire illuminated by plane waves is computed. This allows to classify the resonant states of the nanowire into states which couple to the light sources and into states which have no contribution to the extinction cross section.

This work is structured as follows. Section II introduces a coupled system of equations describing the HDM and summarizes numerical methods for modal analysis. In the subsection Modal expansion of sesquilinear quantities, we extend the framework of the Riesz projection expansion (RPE) Zschiedrich et al. (2018) in order to obtain modal expansions of physical observables, such as the extinction cross section. Section III applies these methods for an investigation of the resonances of a metal nanowire. Section IV concludes the study.

II Plasmonic resonances based on the hydrodynamic Drude model

The HDM is based on the interaction of a nonlocal polarization current and its resulting electric field. In the frequency domain and for nonmagnetic materials, this is described by the coupled system of equations,

[TABLE]

for the electric field $\mathbf{E}(\mathbf{r},\omega)$ and the nonlocal hydrodynamic current density $\mathbf{J}_{\mathrm{hd}}(\mathbf{r},\omega)$, where $\mathbf{J}(\mathbf{r},\omega)$ is a given impressed current density, $\omega$ is the frequency, $\epsilon_{\mathrm{loc}}(\mathbf{r},\omega)$ is the permittivity resulting from the local material response, $\epsilon_{0}$ is the vacuum permittivity, and $\mu_{0}$ is the vacuum permeability. The damping constant $\gamma$ and the plasma frequency $\omega_{\mathrm{p}}$ correspond to the local Drude model $\epsilon_{\mathrm{d}}(\omega)=\epsilon_{0}(\epsilon_{\infty}-\omega_{\mathrm{p}}^{2}/(\omega^{2}+i\gamma\omega))$, where $\epsilon_{\infty}$ is the relative permittivity at infinity. The factor $\beta=\sqrt{3/5}\,v_{\mathrm{F}}$ relates to the Fermi velocity $v_{\mathrm{F}}$ Boardman (1982).

The nonlocal material response is caused by $\mathbf{J}_{\mathrm{hd}}(\mathbf{r},\omega)$, which affects the permittivity function for the free electron gas. If $\beta\rightarrow 0$, then the coupled system simplifies to Maxwell’s equations for the local Drude model. For an illustration of the effect of the HDM, a nanowire excited by a plane wave is sketched in Fig. 1. While, for the local Drude model, the electric field intensity inside of the nanowire is nearly constant, the electric field pattern is radially oscillating considering the HDM [see Figs. 1(a) and 1(b), respectively]. The reader is referred to Hiremath et al. (2012); Toscano et al. (2012) for a detailed derivation of Eqs. (1) and (2) including the applied assumptions and approximations.

Physical scattering solutions $\mathbf{E}(\mathbf{r},\omega_{0})$ and $\mathbf{J}_{\mathrm{hd}}(\mathbf{r},\omega_{0})$ of the coupled system can be obtained for real frequencies $\omega_{0}\in\mathbb{R}$. The eigenfrequencies are defined as the complex resonance poles $\tilde{\omega}_{k}\in\mathbb{C}$ of the analytical continuation of $\mathbf{E}(\mathbf{r},\omega_{0})$ and $\mathbf{J}_{\mathrm{hd}}(\mathbf{r},\omega_{0})$ into the complex plane yielding $\mathbf{E}(\mathbf{r},\omega)$ and $\mathbf{J}_{\mathrm{hd}}(\mathbf{r},\omega)$, where $\omega\in\mathbb{C}$ Zschiedrich et al. (2018). The resonant states, also called eigenmodes, of the coupled system correspond to these eigenfrequencies.

A Numerical methods for modal analysis

The contour integral method Beyn’s Algorithm Beyn (2012) is applied to numerically solve the nonlinear eigenproblem Güttel and Tisseur (2017) corresponding to the coupled system given by Eqs. (1) and (2). Contour integral methods for such problems require the definition of an integration path in the complex frequency plane which encloses the eigenfrequencies corresponding to the eigenmodes of interest. The numerical integration along this contour projects vector fields onto the space spanned by these eigenmodes. In this way, an approximate eigenspace is constructed. Then, e.g., the methods proposed in Asakura et al. (2009); Beyn (2012) apply a singular-value decomposition (SVD) to this approximate eigenspace and solve a linear eigenproblem of small dimension. The approach presented in Gavin et al. (2018) applies the Rayleigh-Ritz method to the approximate eigenspace and solves a nonlinear eigenproblem of small dimension. The common property of these methods is that they essentially require the solution of scattering problems for the integration points on the chosen contour. This is in contrast to standard approaches for solving nonlinear eigenproblems, such as the Arnoldi method, which are based on linearization of the nonprojected problems using auxiliary fields Saad (2011); Brûlé et al. (2016).

For the modal expansion of scattering problems, an unconjugated scalar product can be used Yan et al. (2018). In this context, it is an open problem how to deal with the expansion of nonholomorphic quantities, e.g., the extinction cross section. The contour-integral-based RPE Zschiedrich et al. (2018) allows one to perform a modal expansion without a scalar product. A solution $\mathbf{E}(\mathbf{r},\omega_{0})$ to the coupled system given by Eqs. (1) and (2) can be expanded into a weighted sum of eigenmodes yielding the coupling of the modes to specific sources $\mathbf{J}(\mathbf{r},\omega_{0})$ with $\omega_{0}\in\mathbb{R}$. Cauchy’s integral formula,

[TABLE]

is exploited, where $\mathbf{E}(\mathbf{r},\omega)$, $\omega\in\mathbb{C}$, is the analytical continuation of $\mathbf{E}(\mathbf{r},\omega_{0})$ into the complex plane and $C_{0}$ is a closed integration path around $\omega_{0}$ so that $\mathbf{E}(\mathbf{r},\omega)$ is holomorphic inside of $C_{0}$. Deforming the integration path and applying Cauchy’s residue theorem yield

[TABLE]

where $\tilde{C}_{1},\dots,\tilde{C}_{K}$ are contours around the eigenfrequencies $\tilde{\omega}_{1},\dots,\tilde{\omega}_{K}$ and $C_{\text{nr}}$ is a contour including $\omega_{0}$, the eigenfrequencies $\tilde{\omega}_{1},\dots,\tilde{\omega}_{K}$, and no additional eigenfrequencies. The Riesz projections,

[TABLE]

corresponding to $\tilde{\omega}_{k}$ describe the coupling of the eigenmodes to the considered source field. The field,

[TABLE]

contains nonresonant components as well as components corresponding to eigenfrequencies outside of the contour $C_{\text{nr}}$. For this modal expansion approach, instead of projecting random vectors as for Beyn’s Algorithm, the numerical integration is performed by solving the coupled system using physical source fields at the integration points.

Equations (1) and (2) are spatially discretized with the finite element method (FEM) Monk (2003); Weiser (2016). The FEM solver JCMsuite is used to solve scattering problems. Perfectly matched layers (PMLs) are applied to realize outgoing radiation conditions Berenger (1994). High order polynomial ansatz functions and mesh refinements are used to reach a sufficient numerical accuracy Lalanne et al. (2019). We write

[TABLE]

for the coupled system given by Eqs. (1) and (2), where $T({\omega})\in\mathbb{C}^{n\times n}$ is the system matrix resulting from the FEM discretization and $v\in\mathbb{C}^{n}$ is the vector corresponding to $\mathbf{E}(\mathbf{r},\omega)$ and $\mathbf{J}_{\mathrm{hd}}(\mathbf{r},\omega)$. The dimension $n$ results from the spatial mesh and from the degrees of the polynomial ansatz functions of the FEM discretization. The right-hand side $f(\omega)$ corresponds to the impressed current density $\mathbf{J}(\mathbf{r},\omega)$ and incoming source fields. In this notation, $T({\tilde{\omega}_{k}})\tilde{v}_{k}=0$ holds for an eigenfrequency $\tilde{\omega}_{k}$ and an eigenmode $\tilde{v}_{k}$. Solving $T(\omega)v=f(\omega)$ with $f(\omega)\neq 0$ corresponds to solving a scattering problem.

B Modal expansion of sesquilinear quantities

Typical physical quantities are quadratic forms associated with a sesquilinear map $q(v,v^{*})$ for solution fields $v$ and their complex conjugates $v^{*}$. Examples include the electromagnetic absorption and the electromagnetic energy flux. For two reasons, the construction of a meaningful modal expansion of sesquilinear forms $q(v,v^{*})$ is not straightforward. First, the missing orthogonality $q(\tilde{v}_{k},\tilde{v}^{*}_{l})\neq 0$ yields cross terms in the expansion. Secondly, the conjugation $v^{*}(\omega_{0})$ renders $q(v(\omega_{0}),v^{*}(\omega_{0}))$ nonholomorphic and the evaluation of this expression for complex eigenfrequencies $\tilde{\omega}_{k}$ is problematic.

To derive a modal expansion of sesquilinear quantities with well-defined expansion coefficients, we extend the framework of the RPE. The method is based on an analytical continuation of the sesquilinear form $q(v(\omega_{0}),v^{*}(\omega_{0}))$ from the real axis $\omega_{0}\in\mathbb{R}$ into the complex plane $\omega\in\mathbb{C}$. We remark that $v^{*}(\omega_{0})$ is the solution to ${T}^{*}({\omega_{0}})v^{*}(\omega_{0})=f^{*}(\omega_{0})$. The system matrix ${T}^{*}({\omega_{0}})$ and the right-hand side $f^{*}({\omega_{0}})$ have analytical continuations, which we denote by ${T}^{\circ}({\omega})$ and $f^{\circ}({\omega})$. Consequently, the analytical continuation of $v^{*}({\omega_{0}})$ reads as

[TABLE]

Finally, this gives the analytical continuation $q(v(\omega),v^{\circ}(\omega))$ into the complex plane and the modal expansion can be computed.

Note that if a solution of the coupled system given by Eqs. (1) and (2) has a pole in $\omega=\tilde{\omega}_{k}$, then its complex conjugate has a pole in $\omega={\tilde{\omega}^{*}_{k}}$. Thus, $q(v(\omega),v^{\circ}(\omega))$ has poles in $\tilde{\omega}_{k}$ and also in ${\tilde{\omega}^{*}_{k}}$. This has to be taken into account for the RPE. The calculation of a modal quantity corresponding to a specific $\tilde{\omega}_{k}$ involves the summation of the Riesz projections for $\tilde{\omega}_{k}$ and $\tilde{\omega}^{*}_{k}$.

As the derivation of $v^{\circ}(\omega)$ is given formally, we remark, for a better physical understanding, that the complex conjugation of the system matrix and the right-hand side corresponds to solving the coupled system for $\omega=-\omega_{0}$ with sign-inverted radiation conditions.

III Resonances of a nanowire

We consider a specific setup, a cylindrical metal nanowire which has also been investigated in the literature, to study HDM-based effects theoretically Ruppin (2001). For typical nanoplasmonic applications, a quantity of interest is the extinction cross section. In the following, we first compute eigenfrequencies and eigenmodes of the nanowire. Based on this, we then investigate the extinction cross section in a modal sense, i.e., it is shown which of the eigenmodes scatter and absorb an incoming source field and which of the modes do not interact with the light source. When the nonlocal HDM is replaced by a local Drude model, only a single resonance is observed in the extinction cross section Ruppin (2001); Schnitzer et al. (2016); Hiremath et al. (2012).

The investigated sodium nanowire of radius $R=2~{}\mathrm{nm}$, infinitely extended in the $z$ direction [see Fig. 1(a)], is described by $\epsilon_{\infty}=1$, $\omega_{\mathrm{p}}=8.65\times 10^{15}~{}\mathrm{s}^{-1}$, $\gamma=0.01\,\omega_{\mathrm{p}}$, and $\epsilon_{\mathrm{loc}}=\epsilon_{0}\epsilon_{\infty}$. The Fermi velocity is given by $v_{\mathrm{F}}=1.07\times 10^{6}~{}\mathrm{m}\mathrm{s}^{-1}$. The nanowire is surrounded by free space with refractive index equal to one. The source field is a $y$ polarized plane wave with unit amplitude propagating in the $x$ direction. For the FEM discretization, a mesh containing about $2000$ triangles with edge lengths from about $0.05~{}\mathrm{nm}$ to $1~{}\mathrm{nm}$ is applied. The polynomial degree of the finite elements is set to $p=3$.

The frequency range $0.4\,\omega_{\mathrm{p}}<\omega_{0}<1.4\,\omega_{\mathrm{p}}$ is selected for the modal analysis. To compute eigenmodes $\tilde{v}_{k}$ using Beyn’s Algorithm, an integration contour around this range is defined. The parameters for the algorithm are $N=160$ integration points, $l=200$ random vectors, and, for the rank drop detection within the SVD, a tolerance of $\mathrm{tol}_{\mathrm{rank}}=10^{-8}$ is chosen. The SVD and the solution of the resulting small linear eigenproblem are performed within Matlab. We obtain $118$ eigenfrequencies inside the integration contour. The imaginary parts of these eigenfrequencies are $\mathrm{Im}(\tilde{\omega}_{k})=-0.0050\,\omega_{\mathrm{p}}$, except for $\tilde{\omega}_{1}=(0.7313-0.0054i)\omega_{\mathrm{p}}$. We note that the eigenmodes corresponding to eigenfrequencies with $\mathrm{Im}(\tilde{\omega}_{k})=-0.0050\,\omega_{\mathrm{p}}$ are localized in the nanowire material, which is modeled with a constant damping $\gamma$. Other loss channels are not significant for these modes. This results in the very similar imaginary parts of the eigenfrequencies. To numerically assess the quality of the approximations of the eigenfrequencies and eigenmodes, we compute the residuals $\mathrm{res}(\tilde{\omega}_{k})=||T(\tilde{\omega}_{k})\tilde{v}_{k}||_{2}/||T(\tilde{\omega}_{k})||_{F}$, where $||\tilde{v}_{k}||_{2}=1$. The residuals for eigenfrequencies within the integration contour are smaller than $6\times 10^{-15}$. The residuals for computed eigenfrequencies outside the integration contour increase with the distance to the contour. The integration points, all computed eigenfrequencies, and the residuals are shown in Fig. 2(a). Plots of the electric field intensity of an exemplary selection of eigenmodes corresponding to eigenfrequencies in frequency ranges below and beyond the plasma frequency are shown in Figs. 2(b) and 2(c), respectively. Note that these eigenfrequencies are semi-simple with an algebraic and geometric multiplicity of two. The chosen indicies of the eigenfrequencies and eigenmodes are increasing with increasing real parts of the eigenfrequencies and are intended to guide the reader through the figures.

Based on the computed spectrum, we investigate the extinction cross section,

[TABLE]

where the first term is the power flux across the boundary of the entire computational domain, denoted by $\delta\Omega$, and the second term is the energy loss in the domain where the nanowire exists, denoted by $\Omega_{\mathrm{nw}}$ Hiremath et al. (2012). The incoming plane wave with real frequencies $\omega_{0}$ is normalized so that the power flux through the geometrical cross section of the nanowire is $P_{\mathrm{pw}}=4\times 10^{-9}~{}\mathrm{W}$. To quantify the coupling of the light source to specific eigenmodes, the RPE is applied. This requires the holomporphic evaluation of sesquilinear quantities from Eq. (3) and yields the modal extinction cross section $\tilde{\sigma}_{k}(\omega_{0})$ corresponding to an eigenfrequency $\tilde{\omega}_{k}$. The direct solution of the coupled system given by Eqs. (1) and (2) yields the total extinction cross section $\sigma_{\mathrm{tot}}(\omega_{0})$.

First, we investigate the modal extinction cross section in a small frequency range including $\tilde{\omega}_{7},\dots,\tilde{\omega}_{12}$. Figure 3(a) shows $\tilde{\sigma}_{7}(\omega_{0}),\dots,\tilde{\sigma}_{12}(\omega_{0})$, and $\sigma_{\mathrm{tot}}(\omega_{0})$. The eigenmode $\tilde{v}_{10}$ has a significant contribution to $\sigma_{\mathrm{tot}}(\omega_{0})$. The contributions of the eigenmodes $\tilde{v}_{7},\tilde{v}_{8},\tilde{v}_{9},\tilde{v}_{11},$ and $\tilde{v}_{12}$ are negligible.

Secondly, in order to understand why a specific eigenmode couples to the incoming plane wave, a fast Fourier transform of the electric field intensities of the eigenmodes on a circle inside the nanowire is performed. This yields the number of intensity maxima of the eigenmodes along the boundary of the nanowire, which we denote by $n_{\varphi}(\tilde{\omega}_{k})$. In this way, it is possible to classify the eigenmodes. Figure 3(b) shows $n_{\varphi}(\tilde{\omega}_{k})$ for the frequency range $0.4\,\omega_{\mathrm{p}}<\mathrm{Re}(\tilde{\omega}_{k})<1.4\,\omega_{\mathrm{p}}$. The field intensities of the six eigenmodes with $n_{\varphi}(\tilde{\omega}_{k})=2$ are plotted in Fig. 3(c). It can be seen that these modes are dipolelike. Due to the relation of the radius of the nanowire and the wavelength of the plane wave, $R\ll\lambda_{0}$, the overlap integral of source field and eigenmode field has a significant contribution only for these modes.

Finally, the modal extinction cross sections $\tilde{\sigma}_{k}(\omega_{0})$ for the eigenfrequencies with $n_{\varphi}(\tilde{\omega}_{k})=2$ are computed. Figure 3(d) shows the sum of the modal extinction cross sections $\sum_{M}\tilde{\sigma}_{k}(\omega_{0}),M=\{k\,|\,n_{\varphi}(\tilde{\omega}_{k})=2\}$. For the investigated scattering of a plane wave, the agreement of the expansion with the total extinction cross section $\sigma_{\mathrm{tot}}(\omega_{0})$ demonstrates that the complex scattering behavior of the HDM-based nanowire is governed by a few eigenmodes only. Note that the total extinction cross section is in agreement with results from the literature Ruppin (2001); Hiremath et al. (2012).

For illumination with different types of source fields, e.g., dipole sources, also the remaining eigenmodes of the rich spectrum can be excited.

IV Conclusions

We investigated the light-matter interaction in nanoplasmonic systems described by the HDM. We presented a contour-integral-based framework for modal analysis, which enables the direct computation of the spectrum of nonlocal material systems. We introduced an approach for the modal expansion of sesquilinear quantities. This opens the possibility to investigate typical physical observables, e.g., the energy flux, the energy absorption, and overlap integrals for extraction efficiencies. Due to the generality of this approach, we expect that it will prove useful also in other fields of physics. Resonant states and the modal extinction cross section of a metal nanowire were calculated. While the spectrum of this system consists of many eigenfrequencies, only a few resonant states have a significant contribution to the extinction cross section. These resonant states were identified and used to expand the quantity of interest.

As demonstrated, nanoplasmonic systems on small length scales exhibit a large number of additional resonant states described by the HDM. A typical feature of these states is their high local field energy concentration. With precisely defined source fields, specific states can be excited. We expect that this will allow for additional degrees of freedom in tailoring light-matter interactions. A modal picture is a prerequisite for the understanding and for the design of corresponding nanoplasmonic devices.

Acknowledgements

We acknowledge Philipp-Immanuel Schneider and Fridtjof Betz for fruitful discussions. We acknowledge funding by the Deutsche Forschungsgemeinschaft under Germany’s Excellence Strategy - The Berlin Mathematics Research Center MATH+ (EXC-2046/1, Project No. 390685689, AA4-6). This work is partially funded through the project 17FUN01 (BeCOMe) within the Programme EMPIR. The EMPIR initiative is co-founded by the European Union’s Horizon 2020 research and innovation program and the EMPIR Participating Countries.