ENZ materials and Anisotropy: Enhancing Nonlinear Optical Interactions at the Nanoscale in Metal/Conducting-Oxide Multilayer Stacks

Epsilon-near-zero materials are exceptional candidates for studying electrodynamics and nonlinear optical processes at the nanoscale. We demonstrate that by alternating a metal and a highly doped conducting-oxide, the epsilon-near-zero regime may be accessed resulting in an anisotropic, composite nanostructure that significantly enhances nonlinear interactions. Using two independent and different computation techniques we show that the structure can enhance the local field intensity by nearly two orders of magnitudes, in large part due to the onset of the effective anisotropy. The investigation of the multilayer nanostructure using a microscopic, hydrodynamic approach also sheds light on the roles of two competing contributions that are for the most part overlooked, but that can significantly modify linear and nonlinear responses of the structure: nonlocal effects, which blueshift the resulting resonance, and the hot electron nonlinearity, which redshifts the plasma frequency as the effective mass of free electrons increases as a function of incident power density. Finally, we show that, even in absence of second order bulk nonlinearity, second order nonlinear processes are also significantly enhanced by the layered structure.


Introduction
In the last few decades the scientific community has focused significant effort on the design and fabrication of artificial nanostructures, with the intent of circumventing limitations imposed by natural materials, like absorption, and provide novel functionalities by harnessing processes that may occur only at the nanoscale, while simultaneously reducing the size of optical devices. Among the plethora of proposed novel materials and nanostructures, those operating near their epsilonnear-zero (ENZ) region have been shown to possess a set of peculiar properties [1][2][3]. These materials may be exploited to manipulate radiation direction and polarization [4][5][6][7][8], but may also serve as platforms to study nonlinear optical interactions thanks to their ability to enhance local fields under specific circumstances, i.e., TM-polarized field incident at oblique incidence [3].
Among the explored nonlinear processes [9][10][11][12][13][14][15][16] that have been studied both theoretically and experimentally [17][18][19][20][21][22][23][24][25][26], second and third harmonic generation stand out.. The ENZ condition is known to manifest itself in any natural material. For example, GaAs, Si, and GaP have their zero crossing for the real part of the permittivity in the UV, metals like Au, Ag and Cu in the visible, and highly-doped oxides like ITO and AZO in the infrared [27]. However, it is possible to design artificial nanostructures and tailor the frequency response at will [28] or to compensate losses [18,29,30], which are thought to represent the main limitation for field enhancement in ENZ materials [31]. Still, loss compensation techniques based on the inclusion of active materials might be challenging for fabrication and not practical for certain applications.
An alternative route to achieving high local field enhancement without resorting to loss compensation techniques is to exploit the effective anisotropy of artificial nanostructures via the realization of longitudinal epsilon-near-zero, or LENZ, materials [32]. Indeed, it has been shown that LENZ possess exceptional abilities that span from perfect light bending [33,34] and angular filtering and polarization control [35], to coherent perfect absorption [36] and control of leaky wave radiation [37]. Moreover, since the effective anisotropy tends to enhance the local field intensity and improve tolerances in terms of acceptance angles and a weaker sensitivity to material thickness [32], LENZ turn out to be good candidates for the investigation of nonlinear processes.
Recently, LENZ have been studied in the context of second order nonlinear optical processes arising from symmetry breaking at each interface of a layered structure. It was shown that they can boost the enhancement of second harmonic generation (SHG) beyond their isotropic counterparts [38]. Here we first investigate third harmonic generation (THG) with the intent to further examine the role of anisotropy on nonlinear processes in ENZ media. We simulate an effective anisotropic medium and demonstrate that, regardless of the thickness of the medium, THG can be enhanced dramatically by properly increasing the material's degree of anisotropy. We then analyze a realistic scenario, where a metal/conducting-oxide multilayer stack composed of Au and Dy:CdO (dysprosium-doped Cadmium Oxide) is used to test the impact of anisotropy on the same third order process. We compare the response of the multilayer with and without the effects of nonlocality, and introduce the contribution of hot electrons, two processes that can be triggered in both metals and conducting oxides, and compete to either blueshift or redshift the resonant linear and nonlinear responses of the structure. Finally, we demonstrate that the high degree of anisotropy achieved in a multilayer environment also enhances second order nonlinear processes even in absence of bulk nonlinearities.

Improving nonlinear phenomena through anisotropy
Anisotropic ENZ or LENZ have been shown to increase SHG conversion efficiencies in the absence of dipole-allowed quadratic nonlinearity [38]. This improvement is possible thanks to the high local fields that can be achieved in anisotropic ENZ [32] even without the aid of geometric resonances. Although it is clear from previous investigations that anisotropy facilitates nonlinear optical processes in the presence of the ENZ condition, to date it is still not clear if nanostructures with a low degree of anisotropy can perform as well as structures with a high degree of anisotropy, and whether or not the improvement associated with the anisotropy also applies to bulk nonlinearities.
In order to shed light on these two aspects we first investigate THG in a homogenous slab of material having the dielectric permittivity of Dy:CdO [39]. Then we artificially modify the transverse components of the dielectric permittivity (εz = εDy:CdO, εx = εy = εDy:CdO +Δε) to simulate a Dy:CdO-like material and introduce an arbitrary degree of anisotropy. For simplicity, in this section we assume the material is isotropic at the third harmonic frequency. The geometry of the structure under consideration is depicted in Fig.1: a slab of thickness d is illuminated by a TMpolarized field incident at a variable angle ϑi. compatible with those measured for other doped oxides, such as ITO and AZO [14,17]. No dipoleallowed second-order nonlinearity is considered. The TM-polarized field is assumed to be a superposition of two monochromatic signals at the FF and TH harmonic frequencies. The TH signal is initially null. Therefore, the time-independent, nonlinear complex current densities at FF and TH may be written as [40]: In order to understand the impact of anisotropy we begin our investigation by calculating the field intensity enhancement, FIE =  As one departs from the isotropic condition (Δε = 0), FIE increases for a certain range of angles of incidence, regardless of slab thickness. We note that the sign of the anisotropy also appears to impact FIE, which attains slightly larger values when   is positive. The plots also confirm that the tolerance to the angle of incidence improves as the thickness of the slab is increased, as indicated in Ref. [32]. Although FIE improvement is generally a good guideline when attempting to boost nonlinear optical processes, the maps of total (forward plus backward) THG for thicker The reason behind the discrepancy between the FIE and THG trends in thicker slabs resides in the definition of FIE, which considers the average value of the longitudinal electric field inside the slab. In fact, the maximum value of the longitudinal electric field component Ez (not shown), which occurs at the first interface between vacuum and the slab, is largest for large negative anisotropy ( 0    ), and produces a stronger TH signal. Put another way, while a large degree of anisotropy is desirable in order to improve FIE, the sign of the anisotropy may also be relevant to determine the efficiency of the nonlinear optical process, especially in thicker slabs.

Third-Harmonic Generation from Au/Dy:CdO multilayer stack
A careful look at Fig. 2 reveals that relatively high |Δε| are generally needed to significantly improve the efficiency of nonlinear processes compared to the isotropic case (Δε = 0). For example, a 100nm thick slab shows a maximum THG efficiency occurring near 65, when Δε -100. At first sight, the problem of assembling a practical device may appear to be daunting. The requirement is to achieve high anisotropy while preserving a zero-crossing in the real part of the effective permittivity in the longitudinal direction. However, if we aim to realize an anisotropic material with 0    , one may resort to the inclusion of metals because they naturally display large, negative permittivities in practically all spectral ranges of interest. For this reason we direct our attention to structures similar to that depicted in Fig.3: a five periods metal-dielectric stack composed of Au [27] and Dy:CdO [39], each 10 nm thick (a = b = 10 nm), for a total thickness d = 100 nm, comparable in size to the slab of Fig. 2 (c) and (e). The multilayer is situated on top of a SiO2 semi-infinite substrate, and it is illuminated by a TM-polarized field incident at a variable angle ϑi.
The effective permittivities of the stack calculated at λFF = 2147 nm (the crossing point of Dy:CdO) are εz,FF = 2.110 -7 + 0.34i, εx,FF = εy,FF = -101.34 + 11.34i at the FF, and εz,TH = 12.46 + 0.51i, εx,TH = εy,TH = -6.52 + 0.92i at the TH. While in the previous section the eventual anisotropy at the TH frequency was ignored, a realistic nanostructure will indeed display anisotropic behavior also at the TH. To boot, in our preliminary analysis we also assumed a purely real Δε. The introduction of metals in the multilayer stack adds a realistic degree of damping in the system, as revealed by the imaginary part of the effective permittivity in the transverse direction.
Metals like gold and silver typically require one free and two bound electron contributions in order to accurately reproduce the local dielectric constant over a broad range that includes UV wavelengths, where interband transitions significantly contribute to the dielectric response.
Conducting oxides, on the other hand, typically require one Drude and one Lorentz electron species. Therefore, gold data were modeled by using a combination of a Drude (  . Both approaches yield qualitatively similar results. We stress that, while gold has been reported to have a nonlinear susceptibility of the order of 10 -16 m 2 /V 2 in the visible range, its magnitude generally depends on pulse duration and may decrease substantially away from the plasma frequency of the material [14], a fact that can easily be verified using the Drude-Lorentz oscillator approach. In addition, the fields are strongly localized inside the CdO layers, a circumstance that greatly diminishes the nonlinear response of gold. Pump irradiance is assumed once again to be 1 GW/cm 2 . FIE and THG maps are shown in Fig. 4 as functions of incident wavelength and angle of incidence. In the multilayer FIE = propagation z E inside all layers of the stack, keeping in mind, however, that field localization is localized mostly inside the CdO layers. Both FIE and THG peak in the vicinity of the zerocrossing for the effective permittivity in the propagation direction. The peaks are also relatively wide as a function of the angle of incidence, in contrast with the response of isotropic ENZ materials [38]. We note that the TH conversion efficiency shown in Fig. 4  and may only be used as a guideline to design multilayers. Care should be exercised when predicting the nonlinear response from the structure using effective medium approaches. In the multilayer scenario, three important aspects need to be taken into account for both linear and nonlinear simulations: (i) the thickness of each layer is 10 nm, so that nonlocal contributions should be taken into account in both linear and nonlinear regimes; (ii) it is known that hot carrier contributions may be triggered in both Au and Dy:CdO, significantly altering the linear and nonlinear responses [42]; (iii) a second-order response that arises from symmetry breaking at the interfaces in not negligible in this kind of structure [38], may lead to cascaded THG [26], and must be included in the calculations [42,44,45].
Therefore, in addition to the nonlinear contributions that arise from bound electrons, which peak in a different wavelength range [42] and are taken into account directly into Maxwell's equations, the nonlinearities associated with equations (1) and (2) are now related directly to hot electrons, and are treated in the Drude portion of the simulation. Adding the usual nonlocal terms to equations (1) and (2), with the understanding that we are now dealing with free electrons, we obtain the following novel equations of motion for the nonlinear, free currents: F  TH  TH  TH  FF  FF  FF  TH  ,TH  TH  TH  TH 3 / 5 ( ) .
p  (3) hot  measures the change of the free electron's effective mass as the medium is irradiated. The temperature of the free electron gas generally depends on the rate of absorption, and thus the fluence of the incident pulse. In turn, if temperature excursions are limited to a few thousand degrees Kelvin, the effective mass changes linearly with temperature, so that to first order the nonlinearity may be assumed to be a function of peak power density, as outlined in Eqs. (3)(4).
Therefore, (3) hot  will depend on equilibrium plasma frequency, electron temperature, frequencydependent conductivity, pulse duration, and will be of order 10 -18 m 2 /V 2 at the TH and 10 -19 m 2 /V 2 at FF.
The second order response is accounted for by evaluating Coulomb, magnetic, and convective contributions on the surface and volume occupied by the free electrons. In the pulsed regime, these nonlinear sources are described in [44,45]. The second harmonic electromagnetic problem is also solved as outlined in Refs. [13,46] using the finite element method. We thus analyze the impact of equations (3) and (4) and discuss the second order nonlinear effects separately.
The introduction of nonlocal and hot electron contributions alters significantly both the linear and nonlinear responses of the nanostructure. We begin our linear analysis by monitoring linear and nonlinear absorption of the multilayer stack. Results are shown in Fig. 5(a) and (b), respectively. so as not to trigger nonlinear effects. Remarkably, the nonlocal terms, whose strength is predetermined by the Fermi energy and effective free electron mass, blueshift the resonance by nearly 200 nm. The effect of pulse duration may also be ascertained from the figure: the 120 fs pulse couples more strongly with the structures, leading to larger local fields, stronger light-matter coupling, and more pump absorption. In Fig.5(b) we show the effect of increasing input power density to 1 GW/cm 2 , and simultaneously triggering hot electron contributions. We monitor pump absorption within the multilayer in the following scenarios: local, with (red line, square markers) and without hot electron contributions (black line, triangle markers); and nonlocal, with the hot electron nonlinearity included (blue line, circle markers) [see Fig. 5(b)]. The inclusion of the hot electron contribution in Fig.5(b) mitigates the blueshift displayed in the linear regime ( Fig.5(a), broadening, redshifting, and smearing the resonance as the plasma frequency decreases (red line, square markers). The significance and complexity of this dynamical, time-dependent interplay cannot be overstated, and it impacts the fundamental and its harmonics. All absorption curves in Fig. 5 were calculated at ϑi=65°.  Fig. 6 with those in Fig.4 (a), two macroscopic differences may be discerned: (i) the nonlocality blue-shifts the structure' spectral response [47,48]; (ii) FIE values drop for all wavelengths and angles of incidence. Moreover, while hot carrier contributions do not seem to alter the maximum FIE that may be achieved in the multilayer, significant distortion of both spectral and angular features are noticeable when hot electrons contributions are introduced [ Fig. 6(b)]. The degree of distortion is proportional to pump irradiance value (not shown) [42].  As mentioned earlier, and as previously reported [38], the second order nonlinear response in these structures is significant. Additional effects occur at the metal/oxide interfaces due to the large field confinement and discontinuities in the free electron density [44,49]. In Fig.8 we limit ourselves to depicting the predicted, reflected SHG conversion efficiency spectrum for ϑi = 65° incident angle, and different pulse durations (50 fs pulses -blue, solid line, 100 fs pulses -red, dashed line, 150 fs pulses, green, dotted line.) SHG undergoes shifts similar to THG as pulse duration is increased, but the maximum conversion efficiency is not altered significantly because SHG is more sensitive to field discontinuities than field localization. Finally, we compare the SH signal generated by the anisotropic, metal/conducting-oxide multilayer stack with the SH signal generated by a single, 20 nm-thick Dy:CdO layer set in a Kretschmann configuration [17], assuming a MgO prism and a 150 fs pulse incident at ϑi = 36° (black line, diamond markers). The figure suggests that at least one order of magnitude and 100nm separate the efficiency of the SH signals generated under the two different circumstances.

Conclusions
In conclusion, we have attempted to clarify the role anisotropy plays in third order nonlinear processes. We demonstrated that a high degree of anisotropy can enhance nonlinear processes thanks to the high local field values that can be achieved. We also found that the sign of the anisotropy plays a different role in FIE and TH conversion efficiency when thicker slabs are investigated, suggesting that it may possible to exploit metals to achieve a large degree of anisotropy. Finally, we demonstrated that while an effective medium approach may be useful as a general guideline to design appropriate nanostructures, care should be exercised when nonhomogeneous field distributions are achieved, as in the proposed multilayer. The introduction of nonlocal effects and hot electron contributions in fact reveals how these two effects compete and interact, leading not only to significant spectral and angular shifts, but also to the possible suppression/enhancement of the generated signal that depends on pump excitation conditions.