Nonlocal optical response is one of the emerging effects on the nanoscale for particles made of metals or doped semiconductors. Here we classify and compare both scalar and tensorial nonlocal response models. In the latter case the nonlocality can stem from either the longitudinal response, the transverse response, or both. In phenomenological scalar models the nonlocal response is described as a smearing out of the commonly assumed infinitely localized response, as characterized by a distribution with a finite width. Here we calculate explicitly whether and how tensorial models, such as the hydrodynamic Drude model and generalized nonlocal optical response theory, follow this phenomenological description. We find considerable differences, for example that nonlocal response functions, in contrast to simple distributions, assume negative and complex values. Moreover, nonlocal response regularizes some but not all diverging optical near fields. We identify the scalar model that comes closest to the hydrodynamic model. Interestingly, for the hydrodynamic Drude model we find that actually only one third (1/3) of the free-electron response is smeared out nonlocally. In that sense, nonlocal response is stronger for transverse and scalar nonlocal response models, where the smeared-out fractions are 2/3 and 3/3, respectively. The latter two models seem to predict novel plasmonic resonances also below the plasma frequency, in contrast to the hydrodynamic model that predicts standing pressure waves only above the plasma frequency.
© 2015 Optical Society of America
In plasmonics research, conducting materials are described on many levels, from classically all the way to atomistically [1, 2]. Here we focus on effects of nonlocal response (also known as ‘spatial dispersion’) that for metals emerge on the few-nanometer scale, which is the intermediate regime where deviations from classical electrodynamics occur while atomistic descriptions are not needed yet. Nonlocal response means that the polarization field at a certain point depends not just on the electric field in that point, but rather on the electric field in its neighborhood. In Fourier space it is characterized by a wavevector-dependent dielectric response.
In microscopic descriptions of light-matter interactions, nonlocal optical response may ultimately occur due to the quantum mechanical position uncertainty of the particles or material excitations that the light interacts with. The optical properties of conductors is typically dominated by the free-carrier response, which is a collective effect (except for semiconductors at very low doping ) that allows semi-classical descriptions of nonlocal response, both for metals [4–6], for doped semiconductors [3, 7–9], and recently also for graphene [10, 11]. Indeed various nonlocal-response models for conductors have been proposed, predicting different phenomena. This article is not an attempt to review all of these models but rather to classify some of them, after introducing a few developments in this active subfield of plasmonics.
Probably the best known semiclassical model to describe nonlocal collective response is the hydrodynamic Drude model (HDM) that goes back to Bloch [4,12] who included the dynamics of the electron gas in the Thomas-Fermi approximation. Although the HDM has been around for a while, it is relevant for current research both because new experiments become possible and because new effects are predicted with the HDM. Three nonclassical phenomena that are well known to be explained by the HDM are (i) size-dependent blueshifts of the dipolar resonance of noble metal nanoparticles as their sizes get smaller than typically 10 nanometers, see for example Ref. , (ii) standing bulk-plasmon wave resonances above the plasma frequency, first observed as anomalous absorption in Ag thin films , and (iii) somewhat reduced plasmonic field enhancements [15, 16].
A natural extension of the HDM to include nonlocality due to surface diffusion is the Generalized Nonlocal Optical Response (GNOR) model [6, 17]. Recently it was recognized that the HDM predicts that surface-plasmon resonances of higher multipolar order are blueshifted more [18, 19], and in the GNOR model also broadened more . The latter nonlocal broadening was used very recently to explain experimental electron-loss spectra of embedded silver spheres, in particular the disappearance of higher-order plasmon resonances for spheres of radii smaller than 4 nanometers .
Very actively studied nowadays are the optical properties of plasmonic dimer systems with (sub-)nanometer-size gaps. Again an important role is played by semi-classical models. Even though they cannot explain all phenomena, to get insight it is important to find out which phenomena can and which cannot be explained by these relatively simple models. In Ref.  the HDM was used to show that large dimers with small gaps respond non-classically even though classical electrodynamics suffices to describe the response of the individual monomers. An intermediate nonlocal regime of gap sizes was anticipated to be well described by semiclassical models, with gaps too small for classical electrodynamics to be accurate, yet too large to allow short-circuiting due to electronic spill-out at both interfaces . In some experiments this intermediate nonlocal regime was not identified , while in others it was . The observed broadening of optical spectra for smaller-gap dimers  is typically interpreted as due to quantum tunneling , so it is rather interesting that the semiclassical GNOR theory systematically neglects tunneling but still predicts an analogous spectral broadening of smaller-gap dimers as due to diffusive nonlocal response [6, 17, 24]. What is more, without further ingredients the GNOR theory also explains the observed spectral broadening for smaller monomers, which of course cannot be explained as due to quantum tunneling. Here we have the interesting situation that a semiclassical nonlocal-response theory suggests to reinterpret microscopic quantum calculations, because besides quantum tunneling there must be at least one other broadening mechanism at work in small-gap dimers. This is confirmed in Ref.  where it is shown that besides damping due to electron-hole pair creation in forward-scattering processes across the gap (i.e. tunneling) also backward scattering processes at metal-air interfaces are important, the latter processes both for monomers and dimers.
Another very recent generalization of the standard hydrodynamic Drude model is the self-consistent hydrodynamic model (with electron gas dynamics beyond the Thomas-Fermi approximation) that can describe electronic spill-out semiclassically . It predicts nonlocal blueshifts for noble metals while for simple metals redshifts are predicted for smaller nanoparticles , both in agreement with experiments. By contrast, the usual (hard-wall) hydrodynamic model systematically excludes spill-out at the interfaces and consequently always predicts nonlocal blueshifts. Furthermore, the self-consistent hydrodynamic model is the first semiclassical theory to describe so-called multipole surface plasmons, also known as Bennett resonances, and which have been experimentally observed, as reviewed in .
In this article we will classify several types of nonlocal-response models, both scalar and dyadic, as summarized in Table 1. We will derive their nonlocal response functions in real space, and compare the results with phenomenological scalar models for nonlocal response [32, 33]. Aim is to obtain a sharper intuition about the models considered and the phenomena that they describe, especially the analogies and differences between scalar and dyadic models.
2. Preliminaries: nonlocal response in real space and in wavevector space
2.1. From real space to wavevector space: Fourier relations
In general, nonlocal response means that the response of a medium at a certain point r depends on the exciting field in a neighborhood of that point, rather than on the exciting field at r only. Here for simplicity we only consider dielectric (nonmagnetic) response of a medium to light, and assume the response to be linear. Then the dielectric function (or rather: the 3 × 3 dielectric tensor ε) depends on two positions rather than one, so that the displacement field at a certain point depends on the electric field in a neighborhood of that point:(1). The nonlocal response function of infinite systems so obtained is also useful for the analysis of finite systems, as an approximation for example , or as input for a phenomenological specular reflection model for conduction electrons at the metal surface .
2.2. Transverse and longitudinal response
The dielectric function is written as a tensor, and indeed in a few important models that we consider below, the response function is tensorial (or dyadic, meaning 3 × 3 tensorial) rather than scalar. For isotropic bulk media there are only two fundamentally different directions in wavevector space, namely the longitudinal direction parallel to a given wavevector, and the transverse directions orthogonal to it. More specifically, given a wavevector k of length k = |k| in Fourier space, k̂ ≡ k/k is the corresponding longitudinal unit vector, where ‘longitudinal’ means that it points in the direction of k. In addition, we can find two mutually orthogonal transverse unit vectors e1,2(k) that by definition are also orthogonal to k. These three unit vectors together thus form an orthonormal basis for k-space. An arbitrary vector v has a projection of length v ·k̂ in the longitudinal direction k̂. The longitudinal projection of the vector is thus given by (v ·k̂)k̂, which is the same as the inner product of v with the dyadic product k̂ ⊗ k̂. Sometimes the symbol ⊗ that signifies dyadic product is omitted, but here we keep it. The 3 × 3 unit tensor I in k-space can then be expanded as
2.3. Real-space nonlocal response: from k-space back to real space
For a given dielectric tensor in k-space, we can find the corresponding dielectric function in real space by inverse Fourier transformation:Appendix of Ref. ): Appendix A and Ref. , but here we prefer instead to take the derivatives also out of the k1-integral and with εL(−k1, ω) = εL(k1, ω) obtain
We are interested in fully explicit forms of the real-space response functions εL,T(r, ω). The procedure to arrive at that in subsequent sections will be to assume specific model functions for εL,T(k1, ω), to perform the k1-integrals in Eqs. (9a) and (9b), and then to evaluate the double spatial derivatives. Here we make the general point that the result of doing the k1-integral will be a radial function, i.e. depending on the length r rather than on the full position difference vector r. That makes doing the spatial derivatives easier, since they follow immediately from the dyadic identityAppendix A. There is a sum rule for the real-space response functions that we will derive, namely that their integrals over all of space should give Eq. (6) and in the second one we identified δ(k1) in between the square brackets. The sum rule (11) will be a convenient consistency test of our results. On the right-hand side we see the k-space response in the limit k → 0. This limit typically coincides with the local-response limit, and then the sum rule describes that the integral over the nonlocal response function equals that of the corresponding local-response model. In this sense, nonlocal-response models ‘smear out’ the response of the corresponding local model. There are infinitely many different functions that have the same k → 0 limit, and equally many ways to smear out the delta-function response of local models. Below we focus on specific response functions that in most cases are physically motivated.
2.4. Relation between dyadic nonlocal response function and dyadic Green function
The dyadic Green function can conveniently be determined using the longitudinal and transverse projectors of Eq. (3). The retarded (dyadic) Green function for an infinite medium with nonlocal response is defined by the equationEq. (5) for ε, the Green function G(k, ω) for isotropic translationally invariant media can be written as the sum of a transverse and a longitudinal part, where GL(k, ω) equals GL(k, ω)Lk and GT(k, ω) equals GL,T(k, ω)Tk, and where the GL,T(k, ω) are scalar functions depending only on the magnitude of the wave vector. Inserting this expansion into Eq. (13) and using the properties of the projectors Lk and Tk of Eq. (4), we obtain two scalar identities for the transverse and the longitudinal parts of the resulting equation separately, whereby the scalar Green functions can be expressed in terms of the dielectric response functions as
The corresponding real-space representations of the Green functions (14) are found by inverse Fourier transformation, where we can use the isotropy of the medium, analogous to Eqs. (9a) and (9b) for the longitudinal and transverse response functions. The result is therefore
3. k-independent response [εL(ω) and εT(ω)]
3.1. Local response [εL(ω) = εT(ω)]
When the scalar functions εL,T(k, ω) actually have no k-dependence, it may be expected that the corresponding response in real space will be local. In general this is incorrect, as shown below, but let us first consider the simplest case for which it is true: if the longitudinal and transverse responses are equal and independent of k, εL(k, ω) = εT(k, ω) = ε(ω), then ε(k1, ω) = ε(ω)I. In this simple situation we do not need to consider longitudinal and transverse response separately since from Eqs. (2) and (6) it follows immediately that
Local response models are employed in most of plasmonics research, and different metals are described by different permittivities ε(ω). No need to elaborate much on this here. For later reference we mention the widely used Drude-like models, with scalar local response functions of the form
3.2. Nonlocal response model with local trace [εL(ω) and εT(ω)]
A slightly more complex situation occurs if the scalar transverse and longitudinal response functions are different, but still are independent of the wavevector, that is if we assumeEq. (2). The real-space response for this model will nevertheless be nonlocal, because the basis of Eq. (2) itself does depend on the wavevector, on its direction to be precise. We now calculate the spatial dispersion explicitly, starting with εL(r, ω): Appendix and in Refs. [35,36]. Hence the longitudinal response does not solely give rise to a delta-function, but also to a “dipolar” nonlocal term scaling as 1/r3. Similarly, for the transverse part we find 35,36]. Again, besides the delta-function term that describes the local part of the response, there is a dipolar nonlocal term scaling with distance as 1/r3. From the above two expressions it is clear that the adjectives “longitudinal” and “transverse” refer to being parallel or perpendicular to the wavevector k in Fourier space, rather than parallel or perpendicular to the position vector r in real space. The total response of this model is the sum of its longitudinal and transverse parts, and is given by Eq. (21) reduces as it should to the scalar local response found earlier in Eq. (16). On the other hand, when εT(ω) and εL(ω) differ, then the dyadic tensor ε(r, ω) is a non-scalar dyadic tensor. What does this mean in a bulk system? It means that the response in the direction defined by r differs from the response in directions perpendicular to r. In particular, along r and perpendicular to it we find the diagonal tensor components Eq. (21). The result is [2εT(ω) + εL(ω)]δ3(r)/3, which is a local quantity. In other words, the nonlocality of the dyadic response of this model does not show up in its trace. Related to this vanishing angle average, the nonlocal term also disappears when taking the volume average ∫ dr of the dielectric response in Eq. (21), because the 4π solid-angle integral over (I − 3r̂ ⊗ r̂) vanishes. In general one could hope that dyadic nonlocal response functions allow an approximately correct scalar description by taking the angle-averaged response, but the present model shows a qualitative change where taking the average turns a nonlocal model into a local one.
To summarize, in a model in which neither the k-longitudinal nor the k-transverse dielectric functions depend on the wavevector, we find that in real space there is nevertheless a concomitant spatial dispersion that scales as 1/r3, unless the two dielectric functions εL(ω) and εT(ω) are identical. This result also puts the commonly employed local response models in another perspective: in local-response models, neither the longitudinal nor the transverse responses are local; rather, the overall local response is the result of an exact cancellation between the longitudinal and the transverse nonlocal response functions.
3.3. Examples of nonlocal response models with local trace?
We have not come across situations in which different εL(ω) and εT(ω) are employed, so in that sense the nonlocal response model with local trace can be considered a toy model that at least gave some insight. It is a rather special model in the sense that the nonlocality does not involve a new length scale. Below we will see instead examples where either εL or εT has an explicit k-dependence, and where in the small-k limit the two become equal. It is actually an interesting open question whether nonlocal response models with local trace may describe a physical state of affairs in some limit of other dyadic nonlocal response models. Toy models or not, below we will make use of the results obtained for this model in order to understand other nonlocal-response models that do have a clear physical motivation.
4. Longitudinal nonlocal response [εL(k, ω) and εT(ω)]
Although we have just seen that even without explicit k-dependence in the longitudinal and transverse scalar functions the real-space response is typically still nonlocal, we now look at models where there is a k-dependence in the scalar longitudinal response response function. Let’s call such models longitudinal nonlocal response models.
4.1. Hydrodynamic Drude model
The hydrodynamic Drude model (HDM) in the Thomas-Fermi approximation in general describes a response of metals to light that is nonlinear [6, 12], but here we will only consider the linear-response limit. The linearized dynamics of the HDM can be described by two coupled equations of motion, in which the electric field is driven by the free-carrier current density and vice versa,17]. Mathematically, by taking the limit |β| → 0 the equation (23b) reduces to Ohm’s law J = σE. By inserting this into the first equation, we would obtain the Maxwell wave equation in a medium with local response as described by εD(ω) of Eq. (17). From the non-classical term proportional to β2∇ ⊗ ∇· J in Eq. (23b) one can identify the nonlocal response in the HDM as due to non-classical pressure waves in the electron density; the observation of corresponding standing pressure waves in finite structures as in  is an important example of the predictive power of the HDM .
The above coupled equations indeed also describe finite systems when combined with the right boundary conditions. For infinite systems on the other hand, we can find the longitudinal and transverse response functions in k-space, where the above coupled equations turn intoEq. (4) we can decompose E(k, ω) uniquely into the sum of EL = Lk ·E and ET = Tk ·E. This is essentially the Helmholtz decomposition carried out in k-space. Projection properties such as Lk · EL = EL and Tk · EL = 0 etcetera follow immediately. An analogous decomposition exists for J(k, ω). By multiplying both equations in (24) from the left with Tk, one obtains two coupled equations involving only the transverse fields ET and JT, the solution of which gives the dispersion relation k2 = εT(k, ω) (ω/c)2. Multiplying Eq. (24) from the left instead with Lk gives the longitudinal dispersion relation εL(k, ω) = 0. This way one obtains these response functions in the HDM as Eq. (25) one can appreciate that the response in the HDM is both a generalization of the Drude-like local response of Eq. (17) and a special case of the general nonlocal response model (5). We note in passing that the self-consistent hydrodynamic model  that can describe electronic spill-out also belongs to the class of longitudinal nonlocal response models.
4.2. Transverse response in real space in hydrodynamic model
From Eq. (25) it can be seen that the transverse dielectric function has no dependence on k, exactly as in the result Eq. (20) for the “nonlocal model with local trace” considered in Sec. 3.2. In fact, the integral is identical, and so is the result, namely that in the hydrodynamic Drude model the spatial dependence of the transverse response function is given byEq. (17). Thus in the HDM we find terms in the spatial dependence of ε(r, ω) that fall off as 1/r3, originating from the transverse hydrodynamic response. Similar rather singular terms appeared in the total response of “nonlocal model with local trace”, and in the local-response limit they cancelled exactly against similar terms originating from the longitudinal response. It is thus interesting to see whether such cancellation also exists in the HDM, so let us now also consider the real-space version of the longitudinal response in the HDM.
4.3. Longitudinal response in real space in hydrodynamic model
In order to obtain εL(r, ω) for the HDM, we insert into the integral Eq. (9a) the specific hydrodynamic form of εL(k, ω) as given in Eq. (25b), which has a bound-carrier part and a free-carrier part. For the first part we need the integralEq. (17).
4.4. Total response in real space in hydrodynamic modelEq. (17). All effects of spatial dispersion are therefore described by the other term in Eq. (30). Sometimes nonlocal response is said to “smear out” the local delta-function response over a finite volume. We can now determine whether that is indeed what happens for the HDM, the archetypical model with nonlocal response. “Smearing out” would only occur if the second term of Eq. (30) besides smooth functions of r contains a delta-function against which some of the local free-carrier response of the first term cancels. Shortly we can and will work out the spatial derivatives of the second term using Eq. (10), but that identity is valid everywhere except in r = 0, and thus is blind for any delta-function contribution. But if for the moment it is only the delta-function contribution that we are interested in, then we can zoom in on ∇ ⊗ ∇[eiqr/(4πr)] so close to the origin that |q|r ≪ 1 and find that the sought delta-function contribution indeed exists and equals that of the longitudinal delta function of Eq. (19). Thus we find that the second term in Eq. (30) was hiding a delta-function term −χD(ω)δ(r)I/3, which partly cancels the local response described by the first term of Eq. (30). It is interesting that the HDM thereby apparently “smears out” only one third of the local free-electron response over a finite spatial region, while leaving two thirds of the free-carrier local response and (of course) the entire local core response unaffected. We arrive at one of our main results, the explicit spatial response of the HDM as Eq. (31). Now the solid-angle integral equals I/3, so that the integral over terms proportional to (I − 3r̂ ⊗ r̂) vanishes. The last term in Eq. (31), the one proportional to r̂ ⊗ r̂, does contribute to the volume integral, and it does so in a neat way since it ensures that indeed ∫ drε(r, ω) equals εD(ω)I, the same outcome as for the local-response Drude-like model. This agrees with the sum rule (11) and it also implies that the nonlocal term in Eq. (30) had a vanishing volume integral. But most importantly, it entails that indeed in the HDM one third (no more and no less) of the free-carrier response is smeared out nonlocally. This interesting result is a direct consequence of the dyadic nature of the HDM. In k-space the response is only k-dependent in the longitudinal direction, not in the two other directions, which makes the HDM only one third non-local. This is a rather general argument to explain the factor 1/3 and we did not use any specifics of the hydrodynamic model other than that it belongs to the longitudinal nonlocal class. We therefore anticipate that the argument holds for the entire class, which for example includes models with a fourth-order polynomial of k in the denominator of εL(k, ω) instead of the parabolic k-dependence for the HDM in Eq. (25b).
Intuitively one might have expected that nonlocal response gives rise to a smooth real-space response function, but for the HDM the expression Eq. (31) shows several diverging near-field terms. In the near field (i.e. for qr ≪ 1), the nonlocal terms are proportional to the longitudinal delta function, thus scaling as r−3, both in directions parallel and perpendicular to r. By contrast, in the far field (qr ≫ 1) the dyadic response along r scales with distance as eiqr/r, while perpendicular to r the response falls off faster, namely as eiqr/r2. So in that sense, the k-dependence of the k-longitudinal response εL(k, ω) gives rise to a nonlocal response in real space that is also predominantly r-longitudinal, stronger along r than perpendicular to it.
The response function ε(r, ω) in Eq. (31) is an integration kernel in Maxwell’s equations, in particular in the constitutive relation Eq. (1). Because of the diverging terms in the kernel, evaluating the integral may be numerically much more costly than in case of a Gaussian kernel, say. Instead of solving nonlocal-response problems as an integro-differential equation, it is therefore often simpler to solve coupled differential equations such as in Eq. (23). But one could try to approximate the kernel to make the integro-differential approach more manageable. A first approximation could be to neglect the dipole term out of the kernel (31)), since its volume integral vanishes. A further approximation would be to map the dyadic HDM response onto the scalar response function that captures it best. We propose to use the angle-averaged response, obtained by taking one third of the trace of Eq. (31) and multiplying by I. The result is32] would be interesting for further study. For more scalar models see Sec. 6.1 and for a warning about using them we refer to the Discussion in Sec. 7.
5. Transverse nonlocal response [εL(ω) and εT(k, ω)]
While in the previous section we considered models with k-dependent longitudinal and k-independent transverse scalar response functions, here we consider just the opposite: models defined by a εL(ω) in combination with a function εT(k, ω) that has an explicit dependence on the wavevector. We discuss some physical predictions of such models in Sec. 7 below, but first let us specify one such a model and study its response in real space. A main question is whether the real-space response function is qualitatively different from the one of the HDM.
We study the model from the recent work of Ref. , which in our notation becomesEq. (19). The real-space transverse response requires a bit more work, and after adding up the longitudinal and transverse contributions one finds the total response function Eq. (34) indeed gives the local-response result I[εcore(ω) + χfree(ω)], in accordance with the sum rule of Eq. (11).
6. Both longitudinal and transverse nonlocal response [εL(k, ω) and εT(k, ω)]
6.1. Scalar models [εL(k, ω) = εT(k, ω)]
Scalar models for nonlocal response have also been proposed in the literature [7, 37], and most phenomenological models are also scalar [32, 33]. In the light of the previous models, the specialty about scalar models is that their longitudinal and scalar response functions are identical,7] assumes ε(k, ω) to have the functional form of εT(k, ω) in Eq. (33b). Here we focus instead on the interesting scalar model of Ref.  where ε(k, ω) is assumed to be of the same hydrodynamic form as εL(k, ω) in Eq. (25b). The crucial difference with the real HDM is that in the scalar model also an explicit k-dependence of the transverse response is assumed. The model can be derived from two coupled equations analogous to the two equations for the HDM, the only difference being that in the second equation, that is in Eq. (23b), the ∇ ⊗ ∇ is replaced by a Laplacian ∇2I: Eq. (25b). This leads to the scalar real-space nonlocal response function Eq. (37) indeed equals εD(ω)I, in accordance with the general sum rule (11). The scalar model (37) that we just obtained can be contrasted with the traced-out HDM model of Eq. (32). Only the latter has the correct amount of smearing out of the local response, whereas Eq. (37) overestimates the amount of nonlocality by a factor of three, at least when assuming that the HDM is accurate.
6.2. Non-scalar models with L- and T-nonlocal response [εL(k, ω) ≠ εT(k, ω)]
In the general case that both the transverse and the longitudinal scalar response functions εL,T(k, ω) are explicitly k-dependent but different, the total real-space response functions can be calculated by adding Eqs. (9a) and (9b). In this case all of the local free-space response will be smeared out, one third as determined by the typical length scales in εL(k, ω), and two thirds on the length scales of εT(k, ω). When assuming a hydrodynamic form of the response functions, this would correspond to different β parameters for the longitudinal and for the transverse scalar responses. Then again the total integral of ε(r, ω) will give the local Drude εD(ω)I, because the sum rule (11) holds for longitudinal and transverse parts separately.
The microscopic RPA calculations for a free-electron gas by Lindhard  also lead to response functions εL,T(k, ω) that are both nonlocal (and different). The main difference with the hydrodynamic Drude model is obviously that the latter has a local transverse response. But also the wavevector dependence of the longitudinal response is different, with more pronounced differences for larger k. In the long-wavelength approximation on the other hand, Lindhard’s εL,T(k, ω) both tend to the same classical Drude form, so we do not obtain a nonlocal model with local trace (see Sec. 3.2) by taking the k → 0 limit. For applications of the Lindhard response functions also for finite systems, see Ref. .
We have compared scalar, longitudinal, and transverse models for nonlocal response and already for translationally invariant systems found qualitatively different predictions, for example whether all local free-carrier response gets nonlocally distributed or not. If we suppose that longitudinal nonlocal models describe metal nanospheres accurately and transverse nonlocal models apply to low-doped semiconductors, can one then still call scalar nonlocal response models phenomenological? It all depends how well they describe physical phenomena.
For finite structures the HDM predicts novel resonances (standing pressure waves of the electron density) due to nonlocal response only above the plasma frequency, and indeed these have been observed experimentally . It was found for metals that scalar nonlocal models with hard-wall boundary conditions exhibit novel nonlocal standing-wave resonances also below the plasma frequency, for example in Ref.  (and follow-up papers) where a scalar hydrodynamic-like response function was assumed that in real space becomes Eq. (37). Likewise in Ref.  with assumed Gaussian scalar nonlocal response functions, near the dipole resonance new standing-wave patterns were predicted on a scale determined by the strength of the nonlocality. Since experimental observations of such resonances in metallic nanoparticles to date have not been reported, the longitudinal hydrodynamic model seems to be the better simple model of nonlocal response in metals. Further discussion and graphical comparisons of the resonances in scalar models as compared to the HDM are given in Ref. . It should be mentioned that non-classical plasmonic resonances of another kind have been observed below the plasma frequency for some metals, namely the so-called multipole surface plasmons (reviewed in ), but they are related to electronic spill-out at interfaces which classical electrodynamics and the standard HDM both neglect; the first semiclassical model to describe these non-classical resonances is the self-consistent hydrodynamic model .
In Eq. (32) we constructed the scalar nonlocal response function that comes closest to the standard HDM. It is an interesting question for future studies how well it captures the predictions of the HDM. However, one concern should be that ‘scalarizing’ the HDM response may again lead to unphysical predictions of novel nonlocal standing-wave resonances in nanostructures. If so, then one should perhaps give up attempts to model nonlocal response phenomena in metals by scalar functions. This would still not rule out transverse nonlocal response in metals, but only the assumption that the transverse response function equals the longitudinal one.
The for metals not-so-physical novel nonlocal resonances below the plasma frequency calculated in Refs. [32, 37] are related to the fact that also the transverse response was assumed to be nonlocal . Let us assume that Ref.  is right in modeling the nonlocal response of ultra-low doped semiconductor particles with a transverse nonlocal model. (But here the issue is not settled either, Ref.  should be compared with the scalar nonlocal model for doped semiconductors in Ref.  and with the longitudinal nonlocal model in Ref. ). Then the very interesting question arises whether this transverse model gives rise to new resonances below the plasma frequency, just like in the scalar models of Refs. [32, 37]. If so, then it is to be expected that new nonlocal resonances will show up below the plasma frequency when miniaturizing doped semiconductor nanoparticles. That would constitute a qualitative difference in the nonlocal response of doped semiconductors as compared to metals. Indeed, Ref.  features spectral asymmetries and oscillatory behavior on the high-energy side that could very well be due to such novel nonlocal resonances below the plasma frequency.
We considered different classes of scalar and dyadic nonlocal response models, and in particular their real-space behavior. We were inspired by the intuitive scalar models of Refs. [32, 33] that stressed the smearing out of the local response onto a finite region of space. That is what we also found here, both for scalar and for dyadic models. The simple sum rule Eq. (11) describes exactly this: the spatial integral of the nonlocal response functions equals that of the corresponding local-response model. So whatever delta-function response gets ‘missing’ due to nonlocality is balanced by the sum of the response elsewhere. This equality holds for longitudinal and transverse response functions separately.
What the sum rule does not describe is how much of the local response is smeared out in this way. The phenomenological models Refs. [32, 33] assumed implicitly that all free-carrier response is redistributed. That is indeed what we find for scalar models, but for longitudinal nonlocal models including the hydrodynamic Drude model we found to our initial surprise that only one third of the free-electron response is redistributed nonlocally, while for the class of transverse nonlocal models we find two thirds. In this sense scalar models are more nonlocal than transverse ones, which in turn are more nonlocal than longitudinal ones. We explained these properties of real-space response functions back in k-space, where longitudinal nonlocal response modified only one out of three directions (i.e. only along k), and transverse nonlocal response changed only the two other directions (only normal to k). These general arguments apply to the whole classes of models.
Another clear difference with phenomenological models is that the nonlocal response functions ε(r1 − r2) that we found for some important dyadic models are interesting and rather non-Gaussian: no simple non-negative distributions but instead complex-valued functions that diverge as r−3 for r → 0. Even for the scalar version of the hydrodynamic model a 1/r divergence remains. Although these divergencies are integrable, as guaranteed by the sum rule Eq. (11), it does explain that solving these nonlocal response problems as integro-differential equations can be slow, and slower for dyadic models such as the HDM than for scalar ones.
We finally suggested that novel resonances below the plasma frequency may arise in transverse nonlocal models. If so, and if the transverse nonlocal model of Ref.  indeed applies to doped semiconductors, then effects of novel collective nonlocal resonances may be observed in semiconductor nanoparticles that have no counterpart in metal nanoparticles.
A. Derivation of Eq. (10) and connection with spherical harmonics
By the chain rule, one obtains the identity for radial functions ∇f(r) = f′(r)∇r, and by repeating this chain rule once more one obtains ∇ ⊗ [∇f(r)] = ∇ ⊗ [f′(r)∇r] = f″(r)∇r ⊗ ∇r + f′(r)∇ ⊗ ∇r. By writing out in Cartesian coordinates it is easy to show that ∇r = r/r = r̂, and after only a little more work that ∇ ⊗ ∇r = (I − r̂ ⊗ r̂)/r. From the previous identities it then follows thatEq. (10) of the main text. Applying this formula to f(r) = −1/(4πr) gives Eq. (38) that was derived for r ≠ 0, thereby missing that a delta-function term must be added, as shown in a few lines in Ref. . To make the connection with spherical harmonics, notice that Eq. (38) can be rewritten as 34],
Stimulating discussions in the Structural Electromagnetic Materials group at DTU Fotonik are gratefully acknowledged, in particular with Asger Mortensen, Thomas Christensen, Søren Raza, Nicolas Stenger, Giuseppe Toscano, and Wei Yan. Funding is acknowledged from the Danish Council for Independent Research, Grant No. 1323-00087. The Center for Nanostructured Graphene is sponsored by the Danish National Research Foundation, Project DNRF58.
References and links
1. M. Barbry, P. Koval, F. Marchesin, R. Esteban, A. G. Borisov, J. Aizpurua, and D. Sánchez-Portal, “Atomistic near-field nanoplasmonics: reaching atomic-scale resolution in nanooptics,” Nano Lett. 15(5), 3410–3419 (2015). [CrossRef] [PubMed]
2. X. Chen, J. E. Moore, M. Zekarias, and L. Jensen, “Atomistic electrodynamics simulations of bare and ligand-coated nanoparticles in the quantum size regime,” Nature Commun. 6, 8921 (2015). [CrossRef]
3. H. Zhang, V. Kulkarni, E. Prodan, P. Nordlander, and A. O. Govorov, “Theory of quantum plasmon resonances in doped semiconductor nanocrystals,” J. Phys. Chem. C 118, 16035–16042 (2014). [CrossRef]
4. E. Forstmann and R. R. Gerhardts, Metal Optics Near the Plasma Frequency (Springer, 1986). [CrossRef]
5. O. Keller, “Screened electromagnetic propagators in nonlocal metal optics,” Phys. Rev. B 34(6), 3883–3899 (1986). [CrossRef]
6. S. Raza, S. I. Bozhevolnyi, M. Wubs, and N. A. Mortensen, “Nonlocal optical response in metallic nanostructures,” J. Phys.: Condens. Matter 27, 183204 (2015).
7. R. Ruppin, Optical properties of a spatially dispersive cylinder,” J. Opt. Soc. Am. B 6(8), 1559–1563 (1989). [CrossRef]
8. R. Carmina Monreal, T. J. Antosiewicz, and S. P. Apell, “Diffuse surface scattering in the plasmonic resonances of ultralow electron density nanospheres,” J. Phys. Chem. Lett. 6, 1847–1853 (2015). [CrossRef] [PubMed]
9. A. Moradi, “Infrared absorption spectra of a spatially dispersive polar semiconductor nanowire,” Solid State Commun. 212, 10–13 (2015). [CrossRef]
10. T. Christensen, W. Wang, A.-P. Jauho, M. Wubs, and N. A. Mortensen, “Classical and quantum plasmonics in graphene nanodisks: role of edge states,” Phys. Rev. B 90(24), 241414(R) (2014). [CrossRef]
11. D. Correas-Serrano, J. S. Gomez-Diaz, M. Tymchenko, and A. Alù, “Nonlocal response of hyperbolic metasurfaces,” Opt. Express 23(23), 29434 (2015). [CrossRef]
12. F. Bloch, “Bremsvermögen von Atomen mit mehreren Elektronen,” Z. Physik A 81, 363–376 (1933). [CrossRef]
13. 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(3), 161–166 (2013). [CrossRef]
14. I. Lindau and P. O. Nilsson, “Experimental verification of optically excited longitudinal plasmons,” Phys. Scr. 3, 87–92 (1971). [CrossRef]
15. A. I. Fernández-Domínguez, A. Wiener, F. J. García-Vidal, S. A. Maier, and J. B. Pendry, “Transformation-optics description of nonlocal effects in plasmonic nanostructures,” Phys. Rev. Lett. 108(10), 106802 (2012). [CrossRef] [PubMed]
16. G. Toscano, S. Raza, A.-P. Jauho, N. A. Mortensen, and M. Wubs, “Modified field enhancement and extinction in plasmonic nanowire dimers due to nonlocal response,” Opt. Express 20(4), 4176–4188 (2012). [CrossRef] [PubMed]
17. N. A. Mortensen, S. Raza, M. Wubs, T. Søndergaard, and S. I. Bozhevolnyi, “A generalized non-local optical response theory for plasmonic nanostructures,” Nature Commun. 5, 3809 (2014). [CrossRef]
18. W. Yan, N. A. Mortensen, and M. Wubs, “Green’s function surface-integral method for nonlocal response of plasmonic nanowires in arbitrary dielectric environments,” Phys. Rev. B 88(15), 155414 (2013). [CrossRef]
19. T. Christensen, W. Yan, S. Raza, A.-P. Jauho, N. A. Mortensen, and M. Wubs, “Nonlocal response of metallic nanospheres probed by light, electrons, and atoms,” ACS Nano 8(2), 1745–1758 (2014). [CrossRef] [PubMed]
20. S. Raza, S. Kadkhodazadeh, T. Christensen, M. Di Vece, M. Wubs, N. A. Mortensen, and N. Stenger, “Multipole plasmons and their disappearance in few-nanometer silver nanoparticles,” Nature Commun. 6, 8788 (2015). [CrossRef]
22. H. Jung, H. Cha, D. Lee, and S. Yoon, “Bridging the nanogap with light: continuous tuning of plasmon coupling between gold nanoparticles,” ACS Nano, Article ASAP (2015). [CrossRef]
23. R. Esteban, A. G. Borisov, P. Nordlander, and J. Aizpurua, “Bridging quantum and classical plasmonics with a quantum-corrected model,” Nature Commun. 3, 825 (2012). [CrossRef]
26. 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,” Nature Commun. 6, 7132 (2015). [CrossRef]
27. W. Yan, M. Wubs, and N. A. Mortensen, “Hyperbolic metamaterials: nonlocal response regularizes broadband supersingularity,” Phys. Rev. B 86(20), 205429 (2012). [CrossRef]
29. A. Moradi, Maxwell-Garnett effective medium theory: quantum nonlocal effects,” Phys. of Plasmas 22(4), 042105 (2015). [CrossRef]
30. A. Yanai, N. A. Mortensen, and U. Levy, “Absorption and eigenmode calculation for one-dimensional periodic metallic structures using the hydrodynamic approximation,” Phys. Rev. B 88(20), 205120 (2013). [CrossRef]
31. F. Intravaia and K. Busch, “Fluorescence in nonlocal dissipative periodic structures,” Phys. Rev. A 91(5), 053836 (2015). [CrossRef]
33. N. A. Mortensen, “Nonlocal formalism for nanoplasmonics: phenomenological and semi-classical considerations,” Phot. Nanostr. 11(4), 303 (2013). [CrossRef]
34. P. de Vries, D. V. van Coevorden, and A. Lagendijk, “Point scatterers for classical waves,” Rev. Mod. Phys. 70(2), 447 (1998). [CrossRef]
35. D. P. Craig and T. Thirunamachandran, Molecular Quantum Electrodynamics. An Introduction to Radiation Molecule Interactions (Academic, 1984); see Sec. 3.1: Transverse and longitudinal δ-dyadics.
36. C. Cohen-Tannoudji, J. Dupont-Roc, and G. Grynberg, Photons and Atoms. Introduction to Quantum Electrodynamics (Wiley, 1989); see pp. 36–42.
38. J. Lindhard, “On the properties of a gas of charged particles,” Kgl. Dansk Videnskab. Selskab, Mat. Fys. Medd. 28(8), 1–57 (1954).
39. 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(12), 121412(R) (2011). [CrossRef]