Transient perturbative nonlinear responses of plasmonic materials
Mikko J. Huttunen, Jussi Kelavuori and Marco Ornigotti
Photonics Laboratory, Physics Unit, Tampere University, FI-33014 Tampere, Finland
(Dated: May 6, 2020)
Recent investigations on optical nonlinearities of plasmonic materials suggest their responses may
be even beyond the usual perturbative description. To better understand these surprisingly strong
responses, we develop here a simple but general approach to describe the nonlinear optical response
of plasmonic materials up to nth perturbation order. We apply the approach to understand spectral
broadening occurring in resonant metasurfaces and investigate the enhancement of high-harmonic
generation from multiply-resonant metasurfaces, predicting an over million-fold enhancement of
higher harmonics.
DOI:
I. INTRODUCTION
Photonic metamaterials, i.e., artificial materials that
enable sub-wavelength length scale control of light, play
a key role in modern nanophotonics [1, 2]. Besides the
study of the linear optical responses of metamaterials, a
growing interest has emerged in the recent years to under-
stand also their nonlinear responses [3–7], as they are es-
sential in many different applications including frequency
conversion, high-harmonic generation (HHG), ultrashort
pulse generation and frequency combs [8, 9]. Spec-
tral broadening and supercontinuum generation (SCG)
are particularly interesting nonlinear phenomena that
are utilized to make spectrally broad and coherent light
sources [10], which have found numerous applications for
example in gas sensing [11], generation of few-cycle pulses
[12, 13], optical metrology [14], spectroscopy [15, 16], and
optical coherence tomography [17]. Consequently, these
two nonlinear phenomena are of considerable scientific
and technological interest.
Spectral broadening and SCG occur ubiquitously in
solids, liquids and gases [18–20]. Unfortunately, the in-
trinsic material nonlinearities giving rise to these phe-
nomena are very weak, making it necessary to use strong
excitation pulses from amplified laser systems with peak
intensities of the order of 1014 W/cm2 to exploit those
nonlinear effects. Nonlinear optical fibers, for example,
provide a relatively flexible platform that can be also uti-
lized with moderate peak intensities [10]. However, their
big disadvantage is the requirement for long interaction
lengths, which makes them incompatible with the small
footprint typical of nanophotonic devices.
Integrated photonic supercontinuum sources, on the
other hand, are characterized by small footprint and com-
patibility with mass production, representing an ideal
platform for the realization of on-chip SCG [14, 21]. The
recent demonstrations of SCG from plasmonic nanopar-
ticles and metasurfaces, in particular, have raised the
interest to further investigate their properties, in partic-
ular the out-of-equilibrium dynamics of conduction elec-
trons occurring in such systems, as they are responsible
for both SCG and spectral broadening [22–24]. Many of
the experiments have been performed using only mod-
erate peak intensities (1010–1011 W/cm2), prompting to
investigate the possible origins of the responses, because
conventional theory predicts necessity to use orders of
magnitude higher intensities to realize efficient SCG from
such thin materials [24].
Moreover, because the complex dynamics of conduc-
tion electrons in plasmonic materials dictate their optical
properties [22], it seems natural to assume that the same
dynamics are pivotal also in SCG occurring in such sys-
tems. To this aim, the hydrodynamic model has recently
re-emerged as a powerful tool to understand the conduc-
tion electron dynamics [25–27], and has also been used
to propose that nonlocal and nonperturbative effects oc-
curring in plasmonic materials could play a major role in
the surprisingly strong SCG [24, 27]. However, it is not
yet clear whether other effects could also participate in
the process.
In this work, we investigate the transient nonlinear re-
sponses of plasmonic systems motivated by the fact that
modern experiments are often performed using ultrashort
laser pulses. We extend the simple anharmonic oscillator
model, commonly used to describe nonlinear responses,
by taking into account the transient response of the sys-
tem and by generalizing the treatment to include up to n
perturbation orders. We show that considerable spectral
broadening of ultrashort pulses may take place in res-
onant plasmonic systems, such as in metasurfaces. We
also predict that dramatic 6–7 orders of magnitude en-
hancement of higher-harmonic processes could occur in
plasmonic structures exhibiting multiple resonances.
Our work is organized as follows: in Sect. II we use
the Green function approach to extend the usual anhar-
monic oscillator model to include the transient dynamics
up to nth perturbation order. In Sect. III we discuss
how accounting for the transient response of the system
could result in a significant broadening of the incident
ultrafast pulse. We also discuss how the transient con-
tributions into the light–matter interaction could result
in a significant enhancement of the HHG process in mul-
tiresonant plasmonic structures. In Sect. IV, we draw the
conclusions.
ar
X
iv
:2
00
5.
02
08
0v
1
[p
hy
sic
s.o
pti
cs
]
5 M
ay
20
20
2II. THEORY
We start our analysis by considering the classical an-
harmonic oscillator model and seek to describe the time-
varying displacement of the conduction electron x˜(t). We
assume the electron to be forced into motion by an inci-
dent field
E˜(t) =
1√
2pi
∫
∆ω
E(ω′)e−iω
′t dω′ , (1)
where the Fourier integral extends over the incident pulse
frequencies ∆ω centered at the fundamental frequency
ω. We further assume E˜ to be strong enough, so that
the electron oscillation becomes noticeably anharmonic,
while we take E˜ to be weak enough to allow describ-
ing the anharmonic motion perturbatively. In this case,
the forced oscillation of the electron displacement can
be described using the classical anharmonic oscillator
model [28]
¨˜x+ 2γ ˙˜x+ ω20 x˜+ ax˜
2 = −eE˜(t)/m , (2)
where γ is the damping constant, ω0 is the resonance
frequency of the oscillator, e is the electric charge, m
is the effective mass of the electron and a is a parame-
ter describing the anharmonicity of the oscillator. The
tilde and over-dot notations denote time-varying quan-
tities and their time derivatives, respectively. Although
no analytical solutions of Eq. (2) are available, in the
case of weak nonlinearities (ax˜ ω20), an approximate
perturbative solution of form
x˜(t) = x˜1(t) + x˜2(t) + x˜3(t) + ...+ x˜n(t) , (3)
is commonly solved up to the first few orders in terms of
its steady-state solution [28]. Here, we extend the con-
ventional treatment to include the transient response of
the oscillator and by finding a general solution taking into
account perturbative corrections up to nth order. This
allows us to calculate the perturbative behavior of the
anharmonic oscillator for arbitrary input profiles E˜(t).
Substituting the ansatz above into Eq. (2) and equat-
ing the terms of the same perturbative order, we get the
following set of coupled differential equations, for the var-
ious perturbation orders:
¨˜x1 + 2γ ˙˜x1 + ω
2
0 x˜1 =− eE˜(t)/m , (4a)
¨˜x2 + 2γ ˙˜x2 + ω
2
0 x˜2 =− ax˜21 , (4b)
¨˜x3 + 2γ ˙˜x3 + ω
2
0 x˜3 =− 2ax˜1x˜2 , (4c)
...
¨˜xn + 2γ ˙˜xn + ω
2
0 x˜n =− a
∑
|α|=2
(
2
α
)
x˜α . (4d)
The right-hand side of the nth equation follows from the
multinomial theorem, and takes use of the multi-indices
α = (α1, α2, ..., αn) and x˜
α = x˜α11 x˜
α2
2 ...x˜
αn
n . Note, that
x˜3 can be seen to arise from cascaded nonlinear processes,
where light generated by the 2nd-order correction term
(x˜2) gives rise to higher-order terms via process of sum-
frequency generation. Similarly, the higher-order terms
can be seen to arise from cascaded processes of lower
order.
The complete linear response, i.e., the solution of
Eq. (4a), can be written as a convolution of the Green
function G˜ of the unperturbed oscillator and the incident
field E˜ [29]
x˜1(t) = − e
m
∫ ∞
−∞
G˜(t− t′)E˜(t′) dt′ = − e
m
G˜~ E˜ . (5)
In the last equality we have introduced ~ as a shorthand
to indicate the convolution operation. The time-domain
Green function of an under-damped oscillator is [29]
G˜(t) =
e−γt
ωtr
sin (ωtrt) Θ(t) , (6)
where ωtr =
√
ω20 − γ2 is the natural frequency of the
oscillator and Θ(t) is the step function imposed due to
causality of the system. Note, that Eq. (5) contains an
important result: while the motion of an electron de-
scribed by Eq. (5) will, at equilibrium, oscillate follow-
ing the driving field E˜(t), on a timescale of the order of
τ = γ−1, the electron oscillates at its natural frequency
ωtr, rather than following the impinging field. If this
transient timescale, which for bulk gold and silver, for
example, is around 10 fs and 30 fs, respectively [30], is
comparable with the characteristic timescale of the driv-
ing pulse, efficient coupling between the transient and
driven electron dynamics will take place.
The solution for the second- and third-order correction
terms x˜2 and x˜3 can then be found by solving Eqs. (4b)
and (4c) with the aid of the first-order solution. To do
that, it is convenient to work in the frequency domain,
where the convolution can be regarded as a simple mul-
tiplication, allowing us to rewrite Eq. (5) as
xˆ1(ω) = − e
m
Gˆ(ω)Eˆ(ω) , (7)
where the hat notation indicates that variable is de-
scribed in the frequency domain. The frequency-domain
Green function of the under-damped oscillator Gˆ(ω) is
explicitly given as
Gˆ(ω) =
A0
ωtr
(
1
γ + i(ω − ωtr) −
1
γ + i(ω + ωtr)
)
, (8)
where the scaling of the spectral amplitude A0 is dictated
by the investigated structure. Using Eq. (7) and the re-
sult above, the solutions of Eqs. (4b)–(4d) can be then
3written in the following compact form
xˆ2(ω) = −a Gˆ(ω) [xˆ1(ω)~ xˆ1(ω)] , (9a)
xˆ3(ω) = −2a Gˆ(ω) [xˆ1(ω)~ xˆ2(ω)] , (9b)
...
xˆn(ω) = −a Gˆ(ω)
∑
|α|=2
(
2
α
)
xˆα(ω) . (9c)
Here, the nth-order solution is concisely written using
multi-index notation xˆα = xˆα11 ~ xˆα22 ~ ... ~ xˆαnn . Note,
that because convolution in the frequency domain results
in a sum between the various frequencies involved in the
process, xˆ2(ω) contains terms that oscillate, in addition
to the expected frequency 2ω, also at frequencies 2ωtr,
and ω + ωtr, corresponding, respectively, to the second-
harmonic generation centered at the transient frequency,
and sum-frequency generation involving the driving fre-
quency and the transient frequency. Similarly, xˆ3(ω) will
contain terms oscillating, in addition to the expected 3ω,
also at 3ωtr, 2ω + ωtr, and 2ωtr + ω.
Equations (9a)–(9c) are the main result of our work,
and allow us a convenient way to find the general so-
lution for the driven anharmonic oscillator system gov-
erned by Eq. (2) up to the nth perturbative order. This
general solution includes the possible transient contribu-
tions to the overall response, which may be important
either when the driving fields are very short pulses (∼fs)
or when the system exhibits dynamics occurring at sim-
ilar/longer time scales.
Taking into account the transient dynamics of the con-
duction electrons in plasmonic systems gives rise to the
generation of new frequencies and modulation of existing
frequency components associated with the incident pulse.
Interestingly, we see that the characteristic transient fre-
quency of the system ωtr dictates the occurrence of the
nonlinear processes allowing means to engineer them. For
example, one could control the spectral broadening of an
incident pulse by designing a system where the center fre-
quency of the incident pulse and the transient frequency
ωtr coincide, where the strongest spectral broadening is
expected to occur near ωtr. This spectral broadening
mechanism is interesting because the frequency compo-
nents retain their phase coherence potentially being use-
ful also for ultrashort pulse generation [8, 31].
III. RESULTS AND DISCUSSION
To validate our findings, we first calculate the emis-
sion spectrum of a metasurface consisting of identical
plasmonic gold nanoparticles, by using the perturbative
solution to Eq. (3) as found above. The nanoparticles
are taken to exhibit a localized plasmon resonance, peak-
ing near 1.54 eV and associated with the relaxation time
τ = 10 fs similar to bulk gold [30]. The nonlinear coeffi-
cient a was taken to be 20 m−1s−2, resulting in calculated
SHG efficiencies that agree with earlier experiments on
gold nanoparticle arrays [32].
We start our analysis by taking an incident pulse with
a Gaussian temporal profile with a full width at half
maximum of 100 fs, centered at λ = 805 nm (1.54 eV),
and investigate the spectral broadening of the incident
field upon interaction. The peak intensity of the inci-
dent pulse is assumed to be 20 TW/cm2, which can be
readily achieved using amplified systems. A represen-
tative calculated emission spectrum obtained from the
metasurface defined above is shown in Fig. 1, where the
spectral broadening of the incident pump field (red solid
curve) is clearly visible. The calculations were repeated
while varying the number of included perturbative terms
n = 5, 10, and 15. Surprisingly high perturbative terms
(n > 10) were found to still visibly affect the predicted
spectral broadening. This demonstrates the power of the
introduced approach, which allows these contributions to
be included with ease. Overall, the result shown in Fig. 1
suggests that resonant metasurfaces could be used also
for spectral broadening of ultrashort pulses. The advan-
tages of metasurfaces compared to traditional materials
include compact size, wavelength tunability, and the pos-
sibility to engineer the wavefront of interacting free-space
beams [33, 34].
FIG. 1: Calculated emission spectrum showing spectral
broadening as a function of included perturbation terms
n. Input pump (red curve) is centered at 805 nm
(1.54 eV), and a single material resonance (τ = 10 fs) is
taken to coincide with the pump wavelength. Blue
curve acts as a guide to the eye and visualizes the
spectral profile of the Green function associated with
the metasurface.
As a second demonstration of the introduced approach,
we investigate HHG in plasmonic metasurfaces. It has
been recently realized that collective responses of peri-
odic nanostructures can support narrow resonances [35–
37]. The decay times of these narrow resonances can be
considerably longer than those of bulk metals. In addi-
tion, the possibility to fabricate metasurfaces support-
ing multiple resonances has also been recently demon-
strated [38, 39]. Here, we are interested to understand
how the emission spectrum and process of HHG could be
engineered by utilizing the above-mentioned plasmonic
4(d)
(b)(a)
2nd
3rd
4th
(c) 2nd
3rd
4th
5th
6th
FIG. 2: Calculated emission spectrum (for n = 10) using a pump beam centered at 1.54 eV. Blue curves act as
guides to the eye and visualize the Green functions associated with the metasurfaces. Emission spectra for
singly-resonant (a) and multiresonant (b) metasurfaces associated with resonances exhibiting relaxation times of
τ = 10 fs. HHG up to the fourth harmonic is visible. Emission spectra for singly-resonant (c) and multiresonant (d)
metasurfaces with τ = 100 fs. (d) HHG up to the 6th harmonic is clearly visible.
metasurfaces exhibiting multiple narrow resonances.
First, we consider the metasurface and incident pulse
as described above except for the assumed peak inten-
sity. Here, we use a lower peak intensity of 100 GW/cm2.
In the calculations, the first ten perturbative terms are
taken into account (n = 10). The calculated emission
spectrum is shown in Fig. 2(a) demonstrating HHG.
Then we study a multiresonant metasurface that ex-
hibits material resonances at frequencies of 1.54 eV,
3.08 eV, and 4.62 eV, which coincide with the funda-
mental, second- and third-harmonic peaks of the inci-
dent pulse. The relaxation times of these resonances are
taken to be the same as above (τ = 10 fs). As these
resonances occur at the harmonics of the incident pulse
frequency (1.54 eV), their presence is expected to en-
hance HHG. This is also seen in the calculated emission
spectrum shown in Fig. 2(b), where most notably the
field amplitude of the fourth-harmonic peak is increased
1560-fold compared to the singly-resonant metasurface
[Fig. 2(a)].
We then investigate singly- and multiply-resonant
metasurfaces that exhibit resonances with longer relax-
ation times (τ = 100 fs) than considered above. Such
metasurfaces could potentially be realized by utilizing
Fabry–Pe´rot resonances and collective responses occur-
ring in metasurfaces [38, 39]. The calculated emission
spectra for these singly- and multiply-resonant meta-
surfaces are shown in Figs. 2(c) and 2(d), respec-
tively. Comparing the calculated emission spectra of
singly-resonant metasurfaces [Figs. 2(a) and 2(c)], we see
that the narrower resonance results in a 7700-fold (40-
fold) enhancement of third-harmonic generation (fourth-
harmonic generation).
It is even more interesting to compare the emission
spectra of the singly-resonant metasurface against the
multiresonant case with narrow resonances [Figs. 2(a)
and 2(d)]. In this case, the field enhancements are calcu-
lated to be 3.94×106 (10.36×106) for the third-harmonic
(fourth-harmonic) peak. We also note for the case of the
multiresonant metasurface, even the fifth-harmonic and
sixth-harmonic peaks are now clearly visible. Together
these results prompt to investigate how such multireso-
nant metasurfaces could be designed and realized [38, 39].
Next, we discuss the origin of the higher-order pro-
cesses of spectral broadening and HHG, as it might seem
surprising that the anharmonic oscillator model [Eq. (2)],
exhibiting only a quadratic nonlinearity, predicts their
occurrence. However, their origin is understood by con-
sidering cascaded nonlinear processes, where the lower-
order nonlinear processes give rise to effective higher-
order nonlinear processes [40–44]. Here, at the begin-
ning of the light–matter interaction, the quadratic non-
linearity gives rise to second-order processes and asso-
ciated new frequency components [Eq. (9a)], that will
later on act as seeds in further quadratic nonlinear in-
teractions [Eq. (9b)]. These and subsequent cascaded
5second-order processes will then effectively manifest as
higher-order processes [Eq. (9c)]. Therefore, we see that
by using the introduced approach such cascaded higher-
order processes can be taken into account simply by in-
cluding higher-order perturbative terms. This is inter-
esting because our approach allows simple and intuitive
means to understand and even engineer the occurrence
of higher-order processes in metasurfaces.
Last, we consider the asymptotic behavior of the model
by taking an increasingly long relaxation time τ (i.e. van-
ishingly small damping term γ). In this case, the tran-
sient responses no longer quickly decay suggesting that
they remain non-negligible even under continuous-wave
illumination. Despite vanishingly small γ is strictly non-
physical, we note that recent experimental and numerical
investigations on collective responses of metal nanopar-
ticle arrays, known as surface lattice resonances, suggest
that systems associated with quite long relaxation times
(∼ps) can be fabricated [35, 36, 45].
IV. CONCLUSIONS
We have investigated the nonlinear optical responses
of plasmonic materials by extending the anharmonic os-
cillator model, commonly used to describe nonlinear re-
sponses, by taking into account also the transient re-
sponse of the system and by finding a general solution
including up to n perturbative orders. The approach al-
lows us to better understand the nonlinear responses of
these materials and provides means to engineer their re-
sponses for applications, such as for spectral broadening
or frequency conversion. Our results suggest, that con-
siderable spectral broadening of ultrashort pulses may
occur in metasurfaces due to the transient response of
the system. We also predict that higher-harmonic pro-
cesses could be enhanced over a million-fold by utilizing
multiply-resonant plasmonic structures.
ACKNOWLEDGEMENTS
We thank Ksenia Dolgaleva for her valuable feed-
back and comments. We acknowledge the support of
the Academy of Finland (Grant No. 308596) and the
Flagship of Photonics Research and Innovation (PREIN)
funded by the Academy of Finland (Grant No. 320165).
[1] A. Alu`, M. G. Silveirinha, A. Salandrino, and N. En-
gheta, Phys. Rev. B 75, 155410 (2007).
[2] C. M. Soukoulis and M. Wegener, Nat. Photonics 5, 523
(2011).
[3] M. Kauranen and A. V. Zayats, Nat. Photonics 6, 737
(2012).
[4] M. Lapine, I. V. Shadrivov, and Y. S. Kivshar, Rev.
Mod. Phys. 86, 1093 (2014).
[5] J. Butet, P. F. Brevet, and O. J. Martin, ACS Nano 9,
10545 (2015).
[6] G. Li, S. Zhang, and T. Zentgraf, Nat. Rev. Mater. 2, 1
(2017).
[7] E. Rahimi and R. Gordon, Adv. Opt. Mater. 6, 1 (2018).
[8] T. Brabec and F. Krausz, Rev. Mod. Phys. 72, 545
(2000).
[9] T. J. Kippenberg, R. Holzwarth, and S. A. Diddams,
Science 332, 555 (2011).
[10] J. M. Dudley, G. Genty, and S. Coen, Rev. Mod. Phys.
78, 1135 (2006).
[11] J. M. Langridge, T. Laurila, R. S. Watt, R. L. Jones,
C. F. Kaminski, and J. Hult, Opt. Express 16, 10178
(2008).
[12] J. M. Dudley and S. Coen, OSA Trends Opt. Photonics
Ser. 96 A, 1111 (2004).
[13] M. A. Foster, A. L. Gaeta, Q. Cao, and R. Trebino, Opt.
Express 13, 6848 (2005).
[14] A. R. Johnson, A. S. Mayer, A. Klenner, K. Luke, E. S.
Lamb, M. R. E. Lamont, C. Joshi, Y. Okawachi, F. W.
Wise, M. Lipson, U. Keller, and A. L. Gaeta, Opt. Lett.
40, 5117 (2015).
[15] P. Domachuk, N. A. Wolchover, M. Cronin-Golomb,
A. Wang, A. K. George, C. M. B. Cordeiro, J. C. Knight,
and F. G. Omenetto, Opt. Express 16, 7161 (2008).
[16] C. R. Petersen, U. Møller, I. Kubat, B. Zhou, S. Dupont,
J. Ramsay, T. Benson, S. Sujecki, N. Abdel-Moneim,
Z. Tang, D. Furniss, A. Seddon, and O. Bang, Nat. Pho-
tonics 8, 830 (2014).
[17] B. Schenkel, J. Biegert, and U. Keller, Opt. Lett. 28,
1987 (2003).
[18] R. R. Alfano and S. L. Shapiro, Phys. Rev. Lett. 24, 592
(1970).
[19] P. B. Corkum, C. Rolland, and T. Srinivasan-Rao, Phys.
Rev. Lett. 57, 2268 (1986).
[20] P.-a. Champert, V. Couderc, P. Leproux, S. Fe´vrier,
V. Tombelaine, L. Labonte´, P. Roy, and C. Froehly,
Opt. Express 12, 4366 (2004).
[21] D. Y. Oh, D. Sell, H. Lee, K. Y. Yang, S. A. Diddams,
and K. J. Vahala, Opt. Lett. 39, 1046 (2014).
[22] P. Muehlschegel, H. Eisler, O. J. F. Martin, B. Hecht,
P. Mu, and D. W. Pohl, Science 308, 1607 (2005).
[23] P. Biagioni, D. Brida, J. S. Huang, J. Kern, L. Duo`,
B. Hecht, M. Finazzi, and G. Cerullo, Nano Lett. 12,
2941 (2012).
[24] J. Chen, A. Krasavin, P. Ginzburg, A. V. Zayats, T. Pul-
lerits, and K. J. Karki, ACS Photonics 5, 1927 (2018).
[25] J. E. Sipe, V. C. Y. So, M. Fukui, and G. I. Stegeman,
Phys. Rev. A 21, 4389 (1980).
[26] M. Scalora, M. A. Vincenti, D. De Ceglia, V. Roppo,
M. Centini, N. Akozbek, and M. J. Bloemer, Phys. Rev.
A - At. Mol. Opt. Phys. 82, 1 (2010), arXiv:1006.0709.
6[27] A. V. Krasavin, P. Ginzburg, G. A. Wurtz, and A. V.
Zayats, Nat. Commun. 7, 1 (2016).
[28] R. W. Boyd, Nonlinear optics (Academic Press, San
Diego, 2003) p. 578.
[29] F. W. Byron and R. W. Fuller, Mathematics of Classical
and Quantum Physics (Addison-Wesley, Dover, 1992).
[30] P. B. Johnson and R. W. Christy, P. B. Johnson, and
R. W. Christy, Phys. Rev. B 6, 4370 (1972).
[31] G. Steinmeyer, D. H. Sutter, L. Gallmann, N. Matuschek,
and U. Keller, Science 286, 1507 (1999).
[32] R. Czaplicki, A. Kiviniemi, M. J. Huttunen, X. Zang,
T. Stolt, I. Vartiainen, J. Butet, M. Kuittinen, O. J. F.
Martin, and M. Kauranen, Nano Lett. 18, 7709 (2018).
[33] W. Ye, F. Zeuner, X. Li, B. Reineke, S. He, C. W. Qiu,
J. Liu, Y. Wang, S. Zhang, and T. Zentgraf, Nat. Com-
mun. 7, 1 (2016).
[34] F. Walter, G. Li, C. Meier, S. Zhang, and T. Zentgraf,
Nano Lett. 17, 3171 (2017).
[35] B. Auguie´ and W. L. Barnes, Phys. Rev. Lett. 101,
143902 (2008).
[36] M. J. Huttunen, K. Dolgaleva, P. To¨rma¨, and R. W.
Boyd, Opt. Express 24, 28279 (2016).
[37] V. G. Kravets, A. V. Kabashin, W. L. Barnes, and A. N.
Grigorenko, Chem. Rev. 118, 5912 (2018).
[38] M. J. Huttunen, O. Reshef, T. Stolt, K. Dolgaleva, R. W.
Boyd, and M. Kauranen, J. Opt. Soc. Am. B 36, E30
(2019).
[39] O. Reshef, M. J. Huttunen, G. Carlow, T. Sullivan, J.-M.
Me´nard, K. Dolgaleva, and R. W. Boyd, Nano Lett. 19,
6429 (2019).
[40] G. Assanto, G. I. Stegeman, E. VanStryland, and
M. Sheik-Bahae, IEEE J. Quantum Electron. 31, 673
(1995).
[41] X. Mu, X. Gu, M. V. Makarov, Y. J. Ding, J. Wang,
J. Wei, and Y. Liu, Opt. Lett. 25, 117 (2000).
[42] L. Misoguti, S. Backus, C. G. Durfee, R. Bartels, M. M.
Murnane, and H. C. Kapteyn, Phys. Rev. Lett. 87,
013601 (2001).
[43] K. Dolgaleva, H. Shin, and R. W. Boyd, Phys. Rev. Lett.
103, 113902 (2009).
[44] M. Celebrano, A. Locatelli, L. Ghirardini, G. Pellegrini,
P. Biagioni, A. Zilli, X. Wu, S. Grossmann, L. Carletti,
C. De Angelis, L. Duo`, B. Hecht, and M. Finazzi, Nano
Lett. 19, 7013 (2019).
[45] M. S. Bin-Alam, O. Reshef, Y. Mamchur, M. Z. Alam,
G. Carlow, J. Upham, B. T. Sullivan, J.-M. Me´nard,
M. J. Huttunen, R. W. Boyd, and K. Dolgaleva, arXiv
Prepr. (2020), arXiv:2004.05202.