Gain-assisted harmonic generation in near-zero permittivity metamaterials made of plasmonic nanoshells

We investigate enhanced harmonic generation processes in gain-assisted, near-zero permittivity metamaterials composed of spherical plasmonic nanoshells. We report the presence of narrow-band features in transmission, reflection and absorption induced by the presence of an active material inside the core of the nanoshells. The damping-compensation mechanism used to achieve the near-zero effective permittivity condition also induces a significant increase in field localization and strength and, consequently, enhancement of linear absorption. When only metal nonlinearities are considered, second- and third-harmonic generation efficiencies obtained by probing the structure in the vicinity of the near-zero permittivity condition approach values as high as 10- 7 for irradiance value as low as 10 MW cm- 2 . These results clearly demonstrate that a relatively straightforward path now exists for the development of exotic and extreme nonlinear optical phenomena in the kW cm−2 range.


Introduction
Finite thickness metamaterials composed of sub-wavelength constitutive elements have stimulated researchers to engineer intriguing homogenized material properties, such as negative refractive index, artificial permeability and near-zero permittivity (NZP), which in turn have found applications in many practical scenarios. For example, artificial magnetism has been achieved via the split ring resonator (SRR) structure [1], initially introduced at microwave frequencies and then extrapolated to infrared frequencies by scaling the dimensions of the SRR [2,3]. Plasmonic nanoshells arranged in periodic lattices have been investigated for their unique optical properties, such as strongly directional band gaps [4], negative refractive index [5], NZP [6][7][8] and various applications, such as surface enhanced Raman scattering, surface enhanced infrared absorption [9] and negative refraction [10]. Moreover, artificial composite materials exhibiting NZP capabilities [11] have attracted the attention of the scientific community in view of their potential applications, including tunneling of electromagnetic energy [12], boosting of optical nonlinearities [13], low-threshold nonlinear effects [14,15] and producing narrow directive beams [16,17].
The main advantage of exploiting the NZP condition to access nonlinear processes is offered by the huge field enhancement that occurs at the interface of an NZP medium. For a lossy NZP medium large field enhancement values due to the continuity requirement of the longitudinal component of the displacement field at the interface are achieved for a transverse magnetic (TM)-polarized incident wave and specific angles of incidence [12,14]. Recently, it was demonstrated that efficient second-harmonic (SH) and third-harmonic (TH) generation may be achieved in the proximity of the zero crossing point of the real part of the dielectric permittivity without resorting to any resonant mechanism [15]. Although Vincenti et al [15] showed several solutions to achieve singularity-driven nonlinear processes in natural materials, losses represent a significant impediment to the triggering of low-threshold, nonlinear phenomena. An efficient way to abate damping in composite, artificial nanostructures consists in the inclusion of active media, such as dyes or quantum dots [6,8,18]. For example, in [6] a relative permittivity value of Im(ε) ≈ 10 −4 is estimated in the vicinity of the zero crossing point (Re(ε) = 0) for an array of gold nanoshells, where silica-like cores including 10 mM of Rhodamine 800 fluorescent dyes pumped at λ ≈ 680 nm provided damping compensation. 3 Full-wave linear simulations of the metamaterial nanostructure proposed in [6] are shown in section 2 as a function of frequency and incident angle of the incoming plane wave. Linear and nonlinear responses arising from nanoshells similar to those described in [6] that do not resort to damping compensation have been reported in [19]. The simultaneous availability of the NZP condition and damping compensation, however, complicates both linear and nonlinear responses of the nanostructure. In this paper, we aim to exploit the enhanced narrow-band processes arising from the metamaterial damping-compensation mechanism to boost the already promising theoretical results reported in [15]. We investigate the nonlinear contributions arising from the metal (section 3) and from dielectric media separately by either filling with or hosting the nanoparticles in a nonlinear medium (section 4). Although metals do not possess intrinsic, dipolar, quadratic nonlinear terms, they display an effective second-order nonlinearity that arises from a combination of symmetry breaking at the surface, magnetic response (Lorentz force), inner-core electrons, convective nonlinear sources and electron gas pressure [20]. We also consider third-order nonlinearity of the gold shells that together with effective second-order nonlinear sources contribute significantly to the generated signals [21,22]. A detailed study of each nonlinear contribution, without making any a priori assumption about the relative weights of surface and volume sources, allows us to predict SH and TH conversion efficiencies of the order of 10 −7 when probe irradiance is as low as 10 MW cm −2 . The inclusion of a dielectric nonlinear medium characterized by relatively low second-and third-order susceptibilities shows how to further improve conversion efficiencies, paving the way for as yet unexplored lowthreshold nonlinear phenomena. Finally, in section 5 we analyze the critical role of gain-assisted harmonic generation by investigating the effect of different concentrations of the active material. This study serves also as a guide for future experimental works on the subject, where either dye concentration or pump intensity variations can modify significantly the nonlinear response of the metamaterial.

Linear response of metamaterial slabs composed of plasmonic nanoshells
We start our analysis by considering a composite material having finite thickness h along the z direction, made of four layers of arrayed nanoshells with period a = 100 nm, embedded in a homogeneous host medium with relative permittivity ε h = 2.25 (see figure 1). The composite material extends to infinity in the x and y directions, also with period a. To mitigate metamaterial losses, the core of the nanoshells is filled with glass doped with Rhodamine 800 fluorescent molecules, homogeneously modeled with relative permittivity according to the relation ε c = 2.25 + ε c ( pump , N 0 ), with ε c ( pump , N 0 ) as in [6], assuming that the concentration of dye molecules is N 0 = 10 mM (other concentration values will be considered in section 5), pumping rate pump = 6.5 × 10 9 s −1 and core radius r c = 30 nm. Each shell is made of gold, with relative permittivity ε s modeled through the Drude model as shown in [6], and external radius r s = 35 nm. The relationship between ε c ( pump , N 0 ) and the pump light intensity I pump is not straightforward when the gain medium is either enclosed by a nanoshell or is close to other nano scatterers (as in the proposed case), whereas it is rather simple for an open gain medium. As a first approximation, one may consider the relation in open gain media where I pump = pump h f 30 /σ abs , where σ abs is the dye's absorption cross section, h is Planck's constant and f 30 is the pump frequency [23]. For Rhodamine 800 dye molecules, σ abs = 3.14 × 10 −16 cm 2 and f 30 = 441 THz [24,25], which results in I pump = 6.5 MW cm −2 . Further improvements to the model may be based on a more extensive analysis of the electrodynamic system at the pump Figure 1. Sketch of the composite material with finite thickness h made of four layers of nanoshells with distance a embedded in a homogeneous medium with relative permittivity ε h . The composite material extends to infinity in the x and y directions, with period a. The nanoshell core is made of glass doped with Rhodamine 800 fluorescent molecules, with homogenized relative permittivity ε c and radius r c . The shell is made of gold, with relative permittivity ε s and external radius r s . FF stands for fundamental frequency and SH and TH stand for second and third harmonic, respectively. frequency and on determining local pump fields and local absorption. We note that the analysis carried out here is performed in the steady state. A more extensive analysis that includes pulsed dynamics and rate equations [23] will be implemented in a future work.
We assume a pump/probe optical setup: the pump, used to activate the dye molecules and set the metamaterial in gain condition, is tuned at λ = 680 nm (corresponding to 441 THz) and leads to the linear response shown in figure 2; we will define the details of the probe signal for the harmonic processes in sections 3 and 4. More details of the design of the metamaterial under investigation may be found in [6]. The main feature of the proposed structure is the achievement of efficient damping compensation in the vicinity of the zero crossing point of the real part of the metamaterial effective dielectric constant. These two conditions, in fact, are not coincident in frequency: an effective Im (ε) ≈ 10 −4 is obtained for λ ≈ 712 nm, while Re (ε) ≈ 0 at λ ≈ 710.6 nm. The effective permittivity has been calculated by means of the Nicolson-Ross-Weir method applied to a metamaterial slab of finite thickness under normal incidence conditions [26,27].
The metamaterial slab is illuminated with a TM-polarized plane wave, incident at an angle θ with respect to the z-axis, with electric field and wavevector lying in the x-z plane (figure 1). Transmission, reflection and absorption as functions of wavelength and incident angle are shown in figure 2, calculated using two different methods but yielding similar results: a single dipole approximation (SDA; figures 2(a), (c) and (e)) and a full-wave simulation based on a finiteelement method in the frequency domain (COMSOL Multiphysics; figures 2(b), (d) and (f)). The SDA approach [28][29][30] assumes that each nanoshell behaves as an electric dipole, with induced electric dipole moment p = α ee E loc , where α ee is the electric polarizability modeled via the Mie theory [28] and E loc is the local field acting on it. The SDA is a good approximation when the plasmonic nanoshells are used close to their fundamental plasmonic frequency, when particle dimensions are much smaller than the incident wavelength, and when the period a 3r s . In general, more accurate results may be obtained by including multipole field contributions [28,[31][32][33].
As shown in figure 1, we model N = 4 layers of nanoshells embedded in a homogeneous environment. Let us consider a reference nanoshell located at the nth layer r n = x 0x + y 0ŷ + z 0,nẑ (we can assume that x 0 = 0, y 0 = 0 and z 0,n = na − a/2, with n = 1, . . . , 4). The local electric field acting on it is given by E loc (r n ) = E inc (r n ) + 4 m=1 G layer n,m · p m , where E inc (r) = E 0 e ik·r is the TM incident electric field with the wavevector k lying in the x-z plane, and p m , m = 1, . . . , N , are the electric dipole moments at positions r m . Furthermore, the term G layer n,m ≡ G ∞ (r n , r m ), for n = m, represents the periodic dyadic Green's function (PDGF) for a two-dimensional periodic array of dipoles [34,35], whereas G layer n,n ≡ G ∞ (r n , r n ) represents its regularized counterpart (without the contribution from the dipole at r n ) [34,35]. Both G ∞ and G ∞ are evaluated through the Ewald method [36,37] that provides: (i) rapid converging summations (only a handful of summation terms is needed to achieve convergence) and (ii) analytic continuation to the complex wavenumber plane. The Ewald method for the twodimensional periodic array of scalar Green's functions is detailed in [38][39][40]. Following the formulation in sections III.C and III.D of [35], we then calculate the induced dipole moments p n , n = 1, . . . , 4, by solving the four combined linear equations 4 m=1 A n,m · p m = α ee E inc (r n ), where A n,m = −α ee G layer n,m , for n = m, and A n,n = I − α ee G layer n,n . Once the solution p n has been determined, reflection and transmission coefficients are simply evaluated using again the PDGFs G ∞ (r, r n ). A slight blue shift in the results retrieved via SDA with respect to the full-wave COMSOL result (see figure 2) is observed. Also, transmission, reflection and absorption values are slightly different in the SDA case because we are not taking into account the multipolar response. However, similar spectral and angular features are clearly observed in both cases in figure 2. The presence of a damping-compensation mechanism, or more generally NZP condition in the presence of low losses [41], allows the occurrence of narrow-band features and boosts electric field values. We stress that we do not treat the metamaterial as a homogenized slab with effective parameters, since for any incident angle and wavelength the field localization inside the metamaterial slab varies substantially and thus leads to different nonlinear responses.
In order to investigate the more general behavior of the structure in the proximity of the Re (ε) ≈ 0 condition, we calculated the induced dipole moment p = p xx + p zẑ in the x and z directions at the section z = 3a/2, which corresponds to the second bottom layer of nanoshells (figure 3). Inspection of figure 3 reveals that in the vicinity of the NZP condition, the magnitude of the longitudinal component of the induced dipole moment p z (figure 3(a)) is much larger than the transverse component p x ( figure 3(b)), which remains small. By looking at the phase of the induced dipole moment in the same spectral/angular region, one may also note that the Re (ε) ≈ 0 condition causes a π -shift in the phase of the longitudinal component, while no phase shift is observed for the transverse component. These conditions may be important in phase-sensitive processes such as harmonic generation and optical bistability.
A key aspect for the purposes of efficient harmonic generation and low-threshold optical bistability processes is electric field enhancement inside the structure. We used a full-wave approach to calculate the electric field intensity enhancement, defined as EF = max(|E| 2 /|E inc | 2 ), where |E| is the total electric field magnitude inside the metamaterial slab and |E inc | is the electric field magnitude of the incident plane wave. Our calculations predict the following electric field enhancements: EF ≈ 1000 at λ FF 1 = 709.5 nm and θ = 13 • ; EF ≈ 1900 at λ FF 2 = 710.6 nm and θ = 7 • ; and EF ≈ 3000 at λ FF 3 = 711.35 nm and θ = 4 • , with angular trends that match well the shape of the absorption profile in the structure (see figure 2(f)). Although high field enhancement values are promising for efficient nonlinear processes, one should keep in mind that nonlocal effects may intervene for small metal thicknesses, as may be the case for a 5 nm thick gold shell. However, in what follows, we will assume that nonlocal effects are either negligible or may be compensated for with higher concentrations of gain material, and we postpone this discussion to a future effort.

Nonlinear response of metamaterial slabs arising from metal shells
We now consider a pump-probe scenario, where we assume that a continuous wave (CW) pump signal keeps the system in the damping-compensation regime described in the previous section. We then assume a time harmonic probe signal with irradiance 10 MW cm −2 and calculate its second and third generated harmonics. This relatively small irradiance guarantees that the probe will remain undepleted so that we may neglect any self-phase modulation and down-conversion processes in our calculations. Here we also neglect any self-phase modulation of the pump signal and cross-phase modulation between pump and probe signals. Moreover, we expect that harmonic generated signals triggered by the pump field will not significantly alter the harmonic generation conversion efficiencies relative to the probe frequency for two reasons: (i) both the pump and the probe are CW signals so that their harmonics will fall at different frequencies and (ii) the pump signal is tuned in a spectral region where the metamaterial absorbs and no substantial field enhancement is expected. For these reasons, in what follows, we will neglect the pump contribution to harmonic generation associated with the probe signal.

Second-harmonic generation from gold nanoshells
The nonlinear contributions to second-harmonic generation (SHG) arising from the metal come mostly from free electrons [20,42,43]. The hydrodynamic model is used to describe freeelectron motion, and we do not make any assumptions on the relative weight of surface and volume contributions. A nonlinear equation of motion for the current density of free electrons 8 may be written in the time domain as follows (a tilde (∼) is used to denote quantities in the time domain) [20,42,43]: whereJ f is the free-electron current density, γ f = 7.28 × 10 13 s −1 is the damping coefficient associated with gold losses, e is the electron charge, m * = 9.109 × 10 −31 kg is the effective electron mass, µ 0 is the vacuum magnetic permeability, n 0 = 4.963 × 10 22 cm −3 is the conduction electron density andH andẼ are the local magnetic and electric fields, respectively.
We transform equation (1) into the frequency domain by expressing all the variables as superposition of two components, one oscillating at the fundamental wavelength and the other oscillating at the SH frequency. We may write the following frequency-domain expression for the current density at the SH frequency: The first term on the right-hand side of equation (2) is the linear response at the SH angular frequency ω SH due to both free and bound electrons, and depends on the experimental permittivity of gold at the SH frequency ε s (ω SH ). The second (J vol ) and third (J surf ) terms are bulk and surface nonlinear contributions due to free electrons. These source terms may be derived within the context of the hydrodynamic model in equation (1), following the method illustrated in [44,45] for vacuum-metal interfaces, which is adapted to take into account absorption losses in the metal and the effects of a generic dielectric with permittivity ε B in contact with metal as detailed in [46], n · J surf = − iω SH ε 0 χ (2) ⊥ E 2 FF,⊥ . In expressions (3) and (4) ω FF is the fundamental angular frequency, ε s,f (ω FF ) = −(e 2 n 0 /m * )/(ω 2 FF + iω FF γ f ) is the free-electron permittivity, α = ω FF /(ω FF + iγ f ), β = ω SH /(ω SH + iγ f ) and E FF is the electric field at the fundamental frequency;t andn define a local, boundary coordinate system and represent unit vectors pointing in the directions tangential and outward normal to the metallic surfaces, respectively; E FF, and E FF,⊥ are the magnitudes of the tangential and normal electric field components at the fundamental frequency in this local coordinate system. The surface, quadratic nonlinear susceptibilities are then given by where ε B assumes the value c or h depending on whether the inner or outer interface of the shell is considered for χ (2) and χ (2) ⊥ calculations. In figure 4 we show the forward, backward and total (i.e. forward plus backward) SH conversion efficiencies calculated as η SH = P 2ω /P ω , where P ω is the input probe power and P 2ω is the radiated power (in the forward, backward or forward plus backward directions) at the SH [47]. We analyzed three different wavelengths where either absorption or transmission peaked (λ FF 1 = 709.5 nm-figure 2(f); λ FF 2 = 710.6 nm-figure 2(b); λ FF 3 = 711.35 nm-figure 2(f)). The results in figures 4(b)-(d) reveal a dramatic increase of the generated signal by at least five orders of magnitude when compared with a flat gold film ( figure 4(a)) having the same thickness as the shell (5 nm). Moreover, we note that forward and backward harmonic efficiencies are strongly modulated as the incident angle varies. In particular, while one may certainly surmise a relationship between the trends of linear absorption and total SH generation, it is not straightforward to establish a one-to-one relation between the maxima of the forward and backward harmonic processes and absorption or transmission profile. Field localization and phase accumulation at the probe (the probe is tuned in a region where phase accumulation is close to zero) and harmonic frequencies assume completely different character for each peak of transmission/absorption for the probe, so that forward and backward harmonic generation arising from both surface and volume contributions should be examined on a case-by-case basis and cannot be merely inferred from the linear response of the structure. A rough estimation of the effective second-order nonlinear susceptibility that comes from symmetry breaking at the surface, magnetic dipoles, inner-core electrons and convective nonlinear sources may be obtained by replacing these terms with a χ (2) eff ∼ 5000 pm V −1 . This value gives the reader an idea of the kind of field enhancement that is achieved inside each gold nanoshell.

Third-harmonic generation from gold nanoshells
Like all other metals, gold possesses a relatively large third-order nonlinearity. However, different material configurations, growth techniques and pulse duration have led to a disparate set of values for the third-order susceptibility, ranging from 10 -14 to 10 -19 m 2 V -2 [48][49][50]. Since we are not specifically interested in quantifying or optimizing conversion efficiency values and instead wish to explore the nonlinear optical properties triggered by the NZP and dampingcompensation regimes, in what follows we assume that χ (3) x x x x = χ (3) yyyy = χ (3) zzzz = 10 −17 m 2 V −2 (all other tensor components are assumed to be zero). As we did for SHG, here too we calculate forward, backward and total TH signals at λ FF 1 = 709.5 nm ( figure 5(b)), λ FF 2 = 710.6 nm ( figure 5(c)), and λ FF 3 = 711.35 nm ( figure 5(d)). Conversion efficiency is calculated as η TH = P 3ω /P ω , where P 3ω and P ω are the radiated power at the TH frequency and the input probe power [47], respectively. A cursory examination of figure 5 reveals an increase of the generated TH signal from the nanocomposite metamaterial by about three orders of magnitude when compared with a flat gold film of the same thickness as the shell (5 nm) in figure 5(a). As mentioned above, care should be exercised when relating the features of linear and nonlinear processes. The highly selective angular and spectral features found in both linear and harmonic responses may indeed be ascribed to different mechanisms taking place inside the structure. An investigation of the field distribution inside the structure reveals that a number of resonant modes (three in this case) and an impedance matching condition may be identified in this structure, all leading to different degrees of field localization across the slab and consequently different harmonic responses. Moreover, different phase accumulations in the structure due to the NZP regime where the probe is tuned, together with the difference in parity and phases of the excited resonant modes, may in turn either maximize or abate field overlap, leading to complicated harmonic angular responses and differences between forward and backward generated signals. Moreover, a different angular dependence of TH signals (figure 5) when compared to SH (figure 4) is also expected based on the different trends of the two nonlinear processes for the 5 nm thick flat gold film (compare figure 4(a) with figure 5(a)) and on the fact that they arise either from surface and volume terms or solely from bulk contributions.

Nonlinear response of metamaterial slabs arising from bulk nonlinearities
One may now ask the following questions: (i) how does the nonlinear response change if dielectric bulk nonlinearities are also considered in the calculations? (ii) What is the relative importance of dielectric versus metal nonlinearities, and how much does each contribute to the generated signals? In this section, we will attempt to answer these questions by simulating the effect of bulk nonlinearities associated with the dielectric materials in the system, i.e. host medium and shells' core, and by comparing these results with the radiated SH and TH arising only from the gold shells.
To account for nonlinear effects in the host and core dielectrics, we express the leading contributions of the nonlinear polarization densities in the j direction at second (P SH, j ) and third (P TH, j ) harmonic frequencies as [47] P SH, j = ε 0 3 l,m=1 where j, l, m, n are the Cartesian coordinates, 0 is the vacuum electric permittivity and χ (2) jlm and χ (3) jlmn are the instantaneous second-and third-order susceptibility tensors' components, respectively, and ω TH is the TH angular frequency. In what follows, we will consider χ (3) x values compatible with commercial glasses [47], where the subscript 'h' and 'c' stand for the host and core, respectively, and we assume a relatively low second-order susceptibility value, i.e. χ (2) x x x = χ (2) yyy = χ (2) zzz = χ (2) h = χ (2) c = 10 pm V −1 (all other tensor components are assumed to be zero).

Effects of second-and third-order susceptibilities in the host medium
We begin our study by first introducing second-and third-order susceptibilities in the host medium (indicated by h in figure 1). The nonlinear susceptibilities are included only in the fourperiod-thick region that actually contains the nanoshells and are set to zero outside. In order to identify the relative nonlinear contribution, we neglect all nonlinear terms in the metal described in the previous section. We consider the same probe irradiance (10 MW cm −2 ) and incident wavelength, λ FF 2 = 710.6 nm, where Re( ) = 0. Our nonlinear calculations (figure 6) reveal that both SH and TH have different angular trends when compared with figure 4(c). Moreover, we note that the SH signal radiated by the metal (figure 4(c)) exceeds the SH generated from the host medium (shown in figure 6(a)) by one order of magnitude, while the third-harmonic generation (THG) that arises from the host medium is negligible compared to the TH signal that emanates from the metal (figure 5(c)).

Effects of second-and third-order susceptibilities in the nanoshells' core
We repeat the calculations shown in section 4.1 by considering non-zero second-and thirdorder susceptibilities in the cores of the nanoshells: χ (2) x x x = χ (2) yyy = χ (2) zzz = χ (2) c = 10 pm V −1 and χ (3) x the probe irradiance is 10 MW cm −2 and wavelength is λ FF 2 = 710.6 nm. Just as before, in this scenario SH and TH signals manifest yet a different trend when compared with the same signals arising from the metal and from bulk nonlinearities of the host medium external to the nanoshell. In this case, SH generation becomes less efficient when only the core of the nanoshells is assumed to be nonlinear (see figure 7, in contrast to the case analyzed in section 4.1 where the host medium is nonlinear). On the other hand, the core-nonlinearity boosts THG by approximately one order of magnitude compared to the case of section 4.1 (see figure 7 for comparison). These contrasting trends in SHG and THG are most likely due to the substantial difference in field localization. For example, at λ FF 2 = 710.6 nm the electric field is well localized inside the gold shell, which makes the SH signal significantly larger when the metal nonlinearity is retained. The same field localization arguments cannot be applied directly to the improvements noted for the TH signal because the third-order nonlinearity of the metal is at least four orders of magnitude larger compared to the host or core dielectric's third-order nonlinearity, a discrepancy that could account for most of the improvements that we calculate.

The role of gain
In order to stress the key role of damping-compensation processes in shaping the linear properties to enhance electric field localization to lower nonlinear thresholds, we also calculate transmission, reflection and absorption values for the same structure sketched in figure 1 by assuming different gain conditions. We first compare the absorption maps obtained with full-wave simulations for nanoshells that contain a silica core doped with Rhodamine 800 with decreasing concentrations: 10 mM (figure 8(a)), 7 mM (figure 8(b)), 3 mM (figure 8(c)) and 0 mM ( figure 8(d)). Under the same probing conditions the absorption profiles change dramatically while the Re( ) = 0 condition red-shifts from λ ≈ 710.6 nm at 10 mM to λ ≈ 714.3 nm at 7 mM and from λ ≈ 715.9 nm at 3 mM to λ ≈ 717.7 nm for 0 mM. We stress that the Re( ) = 0 conditions for different concentrations are determined by means of the Nicolson-Ross-Weir method applied to a metamaterial slab of finite thickness under normal incidence conditions [26,27]. By comparing the results in figure 8 one may clearly infer that the majority of the features that characterize absorption (and implicitly also reflection and transmission-not shown here) arise from the ability to compensate for damping in the metamaterial. Although a predominant, wide-band absorption feature is present even without damping compensation ( figure 8(d)), this spectral/angular feature narrows and bends when dye concentration is increased and damping is compensated for. On the other hand, while the effective damping of the structure is being compensated for (Im( ) goes from ≈0.4 for 0 mM to ≈10 -4 for 10 mM, at their respective zero-crossing point for Re( )), absorption is still dramatically high because field localization is boosted. This conclusion is straightforward if one considers that the absorption is proportional to the product Im( eff )|E| 2 . However, high concentrations of dye molecules may impact overall damping compensation due to the presence of fluorescence quenching and other nonradiative phenomena [51]. The quenching effect leads to a reduction of gain in the system, although an analytical treatment of the molecule-nanoparticle interaction is probably needed to estimate its real impact. Alternative designs that, for example, use quantum dots [8] may be employed here.
An immediate consequence of the differences in field localization among the four structures with variable damping compensation is a variation in the strength of harmonic signal. In order to understand how critical is the introduction of proper gain (or alternatively, proper pumping condition) to the nonlinear response, we calculate forward and backward TH signals (arising from the metal shells only) from the four structures with absorption profiles in figure 8 where Re( ) = 0: λ ≈ 710.6 nm for 10 mM (figure 9(a)), λ ≈ 714.3 nm for 7 mM (figure 9(b)), λ ≈ 715.9 nm for 3 mM (figure 9(c)) and λ ≈ 717.7 nm for 0 mM ( figure 9(d)). The detrimental effect induced by the reduction or complete elimination of the active material from the cores of the nanoshells is quite evident from an inspection of figure 9. A dramatic drop in TH efficiencies by about three orders of magnitude is observed when gain is eliminated. These results should also be compared with THG from a 5 nm thick gold film, which exhibits very similar angular behavior and conversion efficiency ( figure 5(a)). The results shown in figure 9 demonstrate quite clearly how the introduction of gain in the metamaterial not only gives a different angular trend to the generated signals but also improves significantly conversion efficiencies between two and three orders of magnitude. Similar results are also expected for SH generation.

Conclusions
We have studied SH and TH generation from metamaterials composed of plasmonic nanoshells working in a damping-compensation regime, to stress the ability of engineered materials exhibiting NZP properties to boost nonlinear processes. The introduction of a gain material into the cores of the nanoshells not only results in peculiar, narrow transmission, reflection and absorption features but also shapes dramatically the nonlinear response increasing both SH and TH efficiencies, which may be evaluated on a case-by-case basis for each resonant feature of the metamaterial. These results demonstrate how efficient, singularity-driven, low-threshold SH and TH generation processes may originate from the strong enhancement of the electric field achievable in metamaterials with NZP properties. Moreover, the investigation of both linear and nonlinear regimes performed by reducing the concentration of gain material or pump intensity proves the extreme sensibility of this regime to pump conditions and serves as a guide for future experimental works within the framework of damping compensation.