Resonant third-harmonic generation driven by out-of-equilibrium electron dynamics in sodium-based near-zero index thin films

1 Introduction

Nonlinear (NL) radiation-matter interactions give rise to a wide range of diverse physical phenomena, e.g., frequency conversion [1], all-optical self-effects [2], generation of non-classical radiation [3], and many others. In particular, harmonic generation in photonic materials provides fundamental insights into quantum mechanical processes, and further offers a promising platform to devise compact ultraviolet (UV) radiation sources [4], [5]. Artificial photonic materials with low dielectric permittivity, known as epsilon-near-zero (ENZ) metamaterials [6]–[8], enhance the inherently weak NL interaction producing extreme NL dynamics [9], second and third-harmonic generation (THG) [10], and provide active control of tunneling [11], and optical switching [12]. Furthermore, they naturally bypass phase-matching requirements owing to the large effective wavelength ensuing at ENZ conditions, thus leading to high THG efficiency [13]. Moreover, ENZ materials can also naturally exist in the form of plasmas, transparent conductors, and metals near their bulk plasma frequency [14]. ENZ media are currently adopted for vortex generation applications [15] and for sensing, guiding, trapping, and emission of visible/infrared (IR) radiation [16]. The ENZ extraordinary boosting of third-order nonlinearity in ENZ thin films [17] ensues from the combined effect of (i) field enhancement for transverse magnetic (TM) excitation, (ii) increased averaged field intensity due to amplified effective wavelength, and (iii) slowdown of light propagation enabling nonlinearity accumulation over time [18]–[22]. Furthermore, amplified Kerr nonlinearity in ENZ conditions enables light-induced “metal–dielectric” transition [23], [24] producing self-organization of frozen light into still solitary spatial doughnuts [25]. Currently adopted ENZ media for enhanced NL optics mainly focus in the near-IR part of the spectrum, where doped semiconductors, oxides, and nitrides display ENZ features [14]. Moreover, poor metals exhibit marked ENZ response in the UV [14] and are promising for manipulation and guidance of UV radiation [26]. Sodium (Na) is a particularly promising UV-ENZ material displaying “metal–dielectric” crossover wavelength at λ _ENZ ≃ 230 nm and an imaginary susceptibility of two orders of magnitude smaller than typical IR-ENZ media like indium–tin–oxide (ITO) [27], thus undergoing mitigated absorption. As a drawback, because Na interacts strongly with air and water, its practical implementation in NL optics applications is challenging. Recently, however, thermo-assisted spin-coating [28] and phase-shift photo-lithography [29] have been adopted to fabricate high-quality Na films, enabling surface plasmon polaritons [30] and thermosensitive plasmonic color [31].

Here, we investigate the potential of Na-based thin films for THG applications at the nanoscale. We model NL electron dynamics in Na by the Fokker–Planck–Landau (FPL) theoretical framework [32], a well-established classical approach in out-of-equilibrium (OOE) statistical mechanics to model plasmas. We calculate analytically the Landau collisional integral in the weak coupling limit [33], obtaining a set of hydrodynamical equations (HDEs) accounting for NL electron dynamics in Na. We emphasize that our HDEs account for damping and heating saturation occurring when OOE instantaneous electron velocities become higher than the thermal velocity owing to collision quenching [33]. By solving HDEs perturbatively, we derive an analytical expression for the dependence over the pump field angular frequency ω = 2πc/λ of the collision-driven THG NL-susceptibility χ coll ( 3 ω ) ( ω ) , where c is the speed of light in vacuum and λ is the pump wavelength. Finally, by solving the scattering problem of the pump field in the undepleted pump approximation, we calculate the THG NL-polarization within the Na-based thin film and model its radiation of λ _THG = λ/3 forward and backward waves, accounting for multiple reflections of the THG signal. We observe that, when the THG wavelength λ _THG matches the angle-dependent near-zero index (NZI) dispersion curve of the Na-based thin film θ _NZI(λ _THG), the THG process becomes resonant for TM pump excitation due to the enhanced radiation amplitude of microscopic dipoles oscillating at λ/3 in an ENZ environment. Moreover, for both TM and transverse-electric (TE) pump excitation, we observe THG enhancement for λ _THG < λ _ENZ due to Fabry–Perot resonances of the Na-based thin film, which becomes transparent in such a spectral range. We illustrate the dependence of the THG process over the Na-based thin film thickness d and over the pump impinging angle, wavelength and intensity. Such a systematic analysis unveils that THG is maximised at a peculiar thin film thickness d ̄ ≃ 40 nm where evanescent back-reflection enhances the average field intensity in the thin film.

2 OOE electron dynamics in Na

We model OOE electron dynamics in Na by classical kinetic theory, where the system is regarded as a dilute electron gas with number density n(r, t) immersed in an idealized homogeneous background of positive uniform density n ₀, physically accounting for lattice ions of mass M staying at equilibrium temperature T ₀. The time-dependent electron distribution function f(r, w, t) is defined from the infinitesimal probability d⁶ P = (1/N)f(r, w, t)d³ rd³ w that an electron in Na with effective mass m ≃ 1.184 × 10⁻³⁰ [34] at time t is placed within a volume element d³ r around r and has a velocity within a velocity-space element d³ w around w, where N is the total electron number within the system volume V. Owing to the large inertia of the lattice ions, we neglect their motion upon external electromagnetic (EM) excitation, so that their distribution function f 0 ( r , w ) = n 0 ( M / 2 π k b T 0 ) 3 / 2 ⁡ exp − M w 2 / 2 k b T 0 stays at equilibrium, where k _b is the Boltzmann constant. The electron distribution normalization ∫f(r, w, t)d³ rd³ w = N accounts for the system charge neutrality dictated by e∫d³ rd³ w(f ₀ − f) = 0, where −e is the electron charge. Within the kinetic theory framework, the temporal evolution of f(r, w, t) is governed by the Boltzmann equation

where F _w (r, t) = −e E(r, t) − e w × B(r, t) is the external EM force exerted by the electric E(r, t) and magnetic B(r, t) pump fields, and ( ∂ t f ) coll = ( ∂ t f ) coll el − el + ( ∂ t f ) coll el − ph is the total collision rate resulting from electron–electron (el–el) ( ∂ t f ) coll el − el and electron–phonon (el–ph) ( ∂ t f ) coll el − ph collisions. While in principle the numerical solution of the Boltzmann equation provides the system’s EM response, it is computationally demanding and lacks insights into NL dynamics. To address electron collisions, the relaxation time approximation (RTA) simplifies the collision integrals but overlooks the NL dependence of ( ∂ t f ) coll over f(r, w, t). In order to model collision-driven NL electron dynamics, we evaluate ( ∂ t f ) coll in the weak coupling assumption where only grazing-angle collisions are considered and improbable frontal collisions are disregarded [33], obtaining

where the C _ee,ep parameters account for el–el and el–ph collisions, D ( β ) = ∇ w ∇ w β and D 0 = D ( β 0 ) are diffusion tensors, α(r, w, t) = ∫d³ w ₁ f(r, w ₁, t)/|w − w ₁| and β(r, w, t) = ∫d³ w ₁ f(r, w ₁, t)|w − w ₁| are the Rosenbluth potentials [32], α ₀ = α(f ₀), and β ₀ = β(f ₀). In this assumption, we solve Eq. (1) by the method of moments and we truncate the ensuing hierarchy of equations at the second moment, obtaining the approximate solution

where n(r, t) = ∫f(r, w, t)d³ w, n(r, t)v(r, t) = ∫w f(r, w, t)d³ w is the current density, and (3/2)n(r, t)k _b T _e (r, t) = (m/2)∫|w − v|² f(r, w, t)d³ w is the OOE energy density. In turn, the OOE temperature of the electron gas T _e (r, t) is defined starting from the OOE energy density, and the moments n(r, t) (zero-order), v(r, t) (first-order), and T _e (r, t) (second-order) satisfy the hierarchy of HDEs

where F _eff = −e E − e v × B is the external effective force, γ ≃ 24.6 ps⁻¹ (depending only over the el–ph collision parameter C _ep) is the linear depolarization rate of Na [27], v T ( T e ) = 2 k b ( T 0 / M + T e / m ) is the electron thermal velocity, v _T0 = v _T(T ₀), γ _th = 2mγ/(m + M) is the hot-electron linear relaxation rate, and

Note in the expressions above that F γ ( v , T e ) and Q(v, T _e ) accounting for NL current damping and electron heating have been calculated by the analytical integration of Rosebluth potentials and of ( ∂ t f ) coll in the weak coupling assumption. Note also that, for weak excitation T _e − T ₀ ≪ T ₀, |v/v _T0|≪ 1, one recovers the RTA limit F γ ≃ 1 . Moreover, in such a limit Q ≃ − γ _th(T _e − T ₀), and thus γ _th represents the hot electron relaxation rate to equilibrium via el–ph collisions.

3 NL response of Na thin films

We consider quasi-monochromatic pump ultrashort pulses with duration of the order of 500 fs, carrier angular frequency ω and associated electric field E 0 ( r , t ) = R e A 0 ( t ) e i k 0 ⋅ r − i ω t impinging over a Na-based thin film with thickness d immersed in air, see Figure 1a. We emphasize that, in view of the quasi-monochromatic regime of such ultrafast excitation, we adopt the slowly-varying envelope approximation (SVEA), where the pump vectorial envelope A ₀(t) is temporally modulated over a duration much longer than the single-cycle period 2π/ω. In turn, in the calculations below, the field envelopes are treated as time-independent quantities and the reported field intensities refer to peak intensities of the considered ultrashort pulses. Without any loss of generality, we assume that k 0 = k x e ̂ x + ω 2 / c 2 − k x 2 e ̂ z lies in the x − z plane, where e ̂ x , y , z indicate the x, y, z unit vectors, k _x = (ω/c)sinθ indicates the in-plane wave-vector component, and θ is the pump incidence angle, see Figure 1a. Owing to multiple reflections by the thin film interfaces, the pump electric field within the ENZ slab E in ( r , t ) = R e A in ( r ) e − i ω t is composed of both forward (FW) and backward (BW) waves, both incorporated within the vector field profile A _in(r). Neglecting relativistic magnetic effects, we consider the weak EM excitation limit v ≃ v ₁ + v ₃, where |A _in(r)| ≃ o(ɛ), |v ₁| ≃ o(ɛ), ΔT = T _e − T ₀ ≃ o(ɛ ²), |v ₃| ≃ o(ɛ ³), and ɛ is a small dummy parameter. Moreover, we focus on the local response where n(r, t) = n ₀, spatial derivatives in Eqs. (4a)–(4c) are neglected, and the electron temperature and mean velocity depend parametrically over the pump electric field at position r. In turn, developing a multiple scale expansion of HDEs, at first-order ≃o(ɛ) the electron mean velocity v ₁(r, t) remains unaffected by heating and satisfies the Drude model ∂ _t v ₁ = −γ v ₁ − (e/m)E _in, providing the solution v 1 ( r , t ) = R e ( − e / m ) ( γ − i ω ) − 1 A in ( r ) e − i ω t . At second order ≃o(ɛ ²), the OOE electron temperature variation ΔT(r, t) is governed by ∂ t Δ T = ( m / 3 k b ) ( 2 γ − γ th ) v 1 2 − γ th Δ T , whose stationary solution arising from the balance between heating (first term) and cooling (second term) provides

where C H = e 2 ( 2 γ − γ th ) / 6 k b m . Note that the OOE electron temperature oscillates with angular frequency 2ω around the intensity-dependent temperature shift produced in such a dynamical equilibrium between collision-induced heating upon illumination and relaxation via el–ph collisions. Finally, at third-order ≃o(ɛ ³), the electron mean velocity correction v ₃(r, t) is governed by

providing the solution v 3 = R e s 3 ( 3 ω ) A in 2 A in e − 3 i ω t + ( s 3 , + ( ω ) | A in | 2 A in + s 3 , − ( ω ) A in 2 A in * ) e − i ω t , where

Figure 1:

Na EM response. (a) Schematic of the considered Na-based thin film with thickness d, illuminated by a pump field with intensity I ₀ and vacuum wavelength λ, producing forward and backward THG at λ _THG = λ/3. (b) Dependence of real (red curve) and imaginary (dashed blue curve) parts of Na relative dielectric permittivity ϵ(λ). The ENZ wavelength of Na λ _ENZ where Re[ϵ(λ _ENZ)] = 0 is indicated by the black solid circle in the figure inset. (c) Pump vacuum wavelength-dependence of the real (red curve) and imaginary (dashed blue curve) parts of the THG susceptibility χ coll ( 3 ω ) ( λ ) . (d) Wavelength dependence of the THG NL susceptibility modulus | χ coll ( 3 ω ) ( λ ) | for several distinct electron thermalization rates γ _th = 1/300, 1/100, 1/50 fs⁻¹, indicated by full, dashed, and dotted curves, respectively.

Note that, owing to OOE electron dynamics, the electron mean velocity NL correction oscillates at ω (Kerr nonlinearity) and 3ω (THG). Kerr nonlinearity of the pump is highly amplified at ENZ conditions [17]. In our THG calculations below, we focus on pump excitation around 600 nm, where the real part of the relative dielectric permittivity of sodium is largely negative Re ϵ ≃ −10, see Figure 1b, and Kerr nonlinearity does not play a crucial role. In turn, in our theoretical framework we neglect Kerr nonlinearity of the pump, obtaining accurate results for λ > 330 nm (where Re ϵ < − 1, see Figure 1b), while at shorter pump wavelengths Kerr-induced spectral shifts of the pump absorbance produce NL modulations of THG efficiency that are not accounted in our theoretical framework. In such assumptions, the macroscopic polarization field P(r, t) can be calculated from the current density J(r, t) = −en ₀ v(r, t) = ∂ _t P(r, t), providing P ( r , t ) = ϵ 0 R e ( ϵ − 1 ) e − i ω t + χ coll ( 3 ω ) A in 2 e − 3 i ω t A in , where ϵ ( ω ) = 1 − ω p 2 / ω ( ω + i γ ) is the Drude relative permittivity, ω p = n 0 e 2 / ϵ 0 m ≃ 8.2 × 1 0 15 rad/s is the sodium plasma frequency [27], ϵ ₀ is the dielectric permittivity of vacuum, and χ coll ( 3 ω ) ( ω ) = m ω p 2 / 3 i e ω s 3 ( 3 ω ) ( ω ) is the collision-induced THG NL susceptibility. Note that the ion mass M enters only in the relaxation rate γ _th. In Figure 1b and c we depict the wavelength dependence of (b) real (red curve) and imaginary (dashed blue curve) parts of the relative dielectric permittivity of Na [ϵ(λ)] and (c) real (red curve) and imaginary (dashed blue curve) parts of its collision-induced THG NL susceptibility χ coll ( 3 ω ) ( λ ) . Note that Na displays marked ENZ behavior at λ _ENZ = 2πc/ω _p ≃ 230 nm, where Im[ϵ(λ _ENZ)] ≃ 10⁻³, see Figure 1b. In turn, conversely to ITO where ENZ behavior is attained in the IR and | ϵ ( λ ENZ ( I T O ) ) | ≃ 0.5 [17], for Na such a quantity becomes much smaller |ϵ(λ _ENZ)| ≃ 10⁻³, showing great theoretical potential for ENZ functionality, see below. Note that the collision-induced THG NL susceptibility is highly dispersive, ranging from | χ coll ( 3 ω ) | ≃ 1 0 − 24 m²/V² at λ = 200 nm to | χ coll ( 3 ω ) | ≃ 1 0 − 21 m²/V² at λ = 800 nm. We emphasize that band anharmonicity and interband transitions (negligible for Na, for which linear optical response is well described by the Drude model) provide other contributions to the third-order NL susceptibility that we do not account here because the scope of the present work focuses on collision-induced NL dynamics. Furthermore, in our calculations we assume γ _th = (300 fs)⁻¹, implying that hot electron relaxation in Na occurs over the same timescale of gold [35], [36]. More precise evaluation of γ _th requires pump-probe reflectivity measurements on Na by experimental schemes similar to the ones adopted for gold [35]. Because electron thermalization arises from electron-phonon scattering, it is highly dependent on sample pristinity. In turn, in practical experimental realizations, the precise electron thermalization rate will depend over the particular sample to be adopted. It is worth emphasizing that electron thermalization rate influences the third-order NL response of the system, see Eqs. (8), while it does not affect its linear response. As a consequence, the resonant THG mechanism theoretically described below in principle depends over the particular electron thermalization rate of the sample. However, we find that THG NL susceptibility is only weakly modulated by the electron thermalization rate, changing only by a factor 1.5 for (50 fs)⁻¹ < γ _th < (300 fs)⁻¹, see Figure 1d, where we depict the wavelength dependence of the THG NL susceptibility modulus | χ coll ( λ ) | for different electron thermalization rates γ _th. Thus, the precise γ _th value does not lead to major quantitative THG NL susceptibility deviations.

4 THG in Na thin films

In order to account for both TM and TE pump excitation of the Na-based thin film, we set A 0 = A 0 n ̂ 0 , where n ̂ 0 = α I ⁡ cos ⁡ θ e ̂ x + β I e ̂ y − α I ⁡ sin ⁡ θ e ̂ z is the impinging field unit vector, A ₀ is the electric field amplitude of the impinging wave, β _I = 1 − α _I and α _I = 0, 1 are dimension-less coefficients enabling selective TE (α _I = 0, β _I = 1) or TM (α _I = 1, β _I = 0) excitation. In turn, in the undepleted pump approximation, the pump field over all space is given by E p = R e A p ( r ) e − i ω t , where A p ( r ) = A 0 e i k 0 ⋅ r + A r e i k r ⋅ r Θ ( − z ) + A in ( r ) Θ in ( z ) + A t e i k 0 ⋅ r Θ ( z − d ) , A r = A r n ̂ r , n ̂ r = α I ⁡ cos ⁡ θ e ̂ x + β I e ̂ y + α I ⁡ sin ⁡ θ e ̂ z , A _r is the reflected field amplitude, k r = k x e ̂ x − ω 2 / c 2 − k x 2 e ̂ z is the reflected wave-vector, A in = ∑ s = ± 1 A s e i k s ⋅ r is the electric field vector amplitude within the thin film, k s = k x e ̂ x + s ( ω / c ) 2 ϵ − k x 2 e ̂ z are the FW (s = 1) and BW (s = −1) wave-vectors of the EM field with vector amplitudes A s = α I A s , x e ̂ x + β I A s , y e ̂ y − s α I ( sin ⁡ θ / n eff ) A s , x e ̂ z , n eff ( ω ) = ϵ − sin 2 ⁡ θ , A _s,{x,y} are the transverse field components in Na, A t = A t n ̂ 0 , A _t is the transmitted field amplitude, Θ_in(z) = Θ(z) − Θ(z − d), and Θ(z) is the Heaviside step function. The pump induction magnetic field over all space is calculated from the Faraday law in local form, providing B p = R e ( 1 / i ω ) ∇ × A p ( r ) e − i ω t . Thus, the unknown transmitted, reflected, FW and BW pump field amplitudes are obtained from boundary conditions (BCs) for the continuity of the magnetic field, the tangential electric field components, and the normal component of the displacement vector at the interfaces z = 0, d. Such BCs provide two uncoupled 4 × 4 inhomogeneous systems of equations for TM and TE amplitudes, explicitly given by

where V TM = [ A r A + , x A − , x A t ] T , W _TM = −A ₀[cosθ 1 0 0]^T, V TE = [ A r A + , y A − , y A t ] T , W _TE = A ₀[−1 cosθ 0 0]^T, η ± ( ω , θ ) = e ± i ω n eff d / c , η ₀(ω, θ) = e^iωcosθd/c and

We emphasize that, for the pump field, all the elements of the TM M ̂ TM ( ω , θ ) and TE M ̂ TE ( ω , θ ) scattering matrices are evaluated at the pump angular frequency ω and incident angle θ, indicated in their superscripts. Also, note that the formalism above refers to single TE (α _I = 0) or TM (α _I = 1) excitation (β _I = 1 − α _I in both cases), but is not suited to treat mixed polarization states, requiring a more general treatment accounting for the polarization dynamics of reflected and transmitted waves. Inverting the system above, one gets all the pump field amplitudes for TE and TM excitation. Figure 2 illustrates the dependence of the pump (a, c) absorbance A = 1 − R − T and (b, d) reflectance R = 1 − |A _r |²/|A ₀|² and transmittance T = 1 − |A _t |²/|A ₀|² over the pump wavelength λ and incidence angle θ for (a, b) TE (α _I = 0, β _I = 1) and (c, d) TM (α _I = 1, β _I = 0) excitation of a Na-based thin film of thickness d = 300 nm. Note that absorbance, reflectance, and transmittance heavily depend over both the pump wavelength and incidence angle when λ < λ _ENZ owing to NZI resonance occurring when 0 < Re[ϵ(λ)] < 1 and R e n eff 2 = 0 , which provides the dispersion relation θ NZI ( λ ) = a s i n R e [ ϵ ( λ ) ] , indicated by the white dashed curves in Figure 2a and c. Note that the maximum absorbance is attained for TM excitation exactly at the NZI resonance, see Figure 2c. Moreover, due to the “metal-dielectric” transition occurring at λ _ENZ, for λ < λ _ENZ Fabry–Perot absorption resonances are excited by both TE and TM impinging polarization, see Figure 2a and c, for θ < θ _NZI(λ). The pump field THG NL polarization within the Na-based thin film is given by P nl ( r , t ) = R e p nl ( z ) e 3 i k x x − 3 i ω t , where

and ϕ _ω (z) = ωn _eff(ω)z/c. Such NL polarization generates 3ω radiation due to the associated polarization (i) current density J _nl = ∂ _t P _nl, and (ii) volume ρ _nl = −∇⋅P _nl and (iii) surface σ nl , ± = ± P nl ( z = d [ 1 ± 1 ] / 2 ) ⋅ e ̂ z charge densities appearing as inhomogeneous terms for the THG signal in (i, ii) Maxwell’s equations within the thin film and (iii) BCs at the two interfaces z = 0, d. We take the Ansatz E THG ( r , t ) = R e e THG ( z ) e 3 i k x x − 3 i ω t , where e THG = a r e 3 i k r , z z n ̂ r Θ ( − z ) + A THG in ( z ) Θ in ( z ) + a t e 3 i k 0 , z z n ̂ 0 Θ ( z ) . Hence, the THG induction magnetic field is obtained from the Faraday law in local form, providing B THG ( r , t ) = R e ( 1 / 3 i ω ) ∇ × e THG ( z ) e 3 i k x x e − 3 i ω t . Inserting this Ansatz in Maxwell’s equations accounting for J _nl and ρ _nl, we obtain an inhomogeneous ordinary differential equation for the THG field within the Na-based thin film, explicitly given by

where

Ω s ( ω , θ ) = − 9 ω 2 χ coll ( 3 ω ) A s 2 A s / c 2 , n ̂ nl = [ 3 ⁡ sin ⁡ θ e ̂ x + s n eff ( ω ) e ̂ z ] , Γ s = A s 2 A − s + 2 ( A s ⋅ A − s ) A s , and ϒ _s = Γ _s − (A _s ⋅A _−s )A _s . Thus, we obtain the solution A THG in = ∑ s = ± 1 Λ ̃ s e i s ϕ ω + Ω ̃ s e 3 i s ϕ ω + a s e 3 i s ϕ 3 ω , where Λ ̃ s = Λ s c 2 / ω 2 9 n eff 2 ( 3 ω ) − n eff 2 ( ω ) and Ω ̃ s = Ω s c 2 / 9 ω 2 n eff 2 ( 3 ω ) − n eff 2 ( ω ) , a s = α I a s , x e ̂ x + β I a s , y e ̂ y − s α I [ sin ⁡ θ / n eff ( 3 ω ) ] a s , x e ̂ z , and a _s,{x,y} are homogeneous amplitudes yet to be determined.

Figure 2:

Linear response of Na-based thin films. Dependence over the pump wavelength λ and incidence angle θ of the (a, c) absorbance A ^(TE,TM) and (b, d) reflectance R ^(TE,TM) (full curves) and transmittance T ^(TE,TM) (dahsed curves) for (a, b) TE and (c, d) TM excitation of a Na-based thin film with thickness d = 300 nm. The dashed white curves in (a, c) indicate the NZI dispersion curve θ _NZI(λ).

In order to obtain the unknown THG field amplitudes, we impose BCs at 3ω (accounting also for the polarization surface charge densities produced by p _nl), explicitly given by

where p z 0 = ( 2 z 0 − d ) / d , z ₀ = 0, d and d _THG = ϵ ₀{1 + [ϵ(3ω) − 1]Θ_in(z)}e _THG. Such BCs provide two uncoupled 4 × 4 inhomogeneous systems of equations for TM and TE amplitudes of the THG signal, explicitly given by

where

K = E, M is an index labelling TE/TM excitation, W ̃ TE = W ̃ y ( 1 ) − W ̃ y ( 2 ) W ̃ y ( 3 ) W ̃ y ( 4 ) T , W ̃ TM = W ̃ x ( 1 ) ( W ̃ x ( 2 ) − sin ̃ θ W ̃ z ( 1 ) ) W ̃ x ( 3 ) ( W ̃ x ( 4 ) + sin ⁡ θ W ̃ z ( 3 ) ) T , W ̃ ( 1 ) = ∑ s = ± 1 [ Ω ̃ s + Λ ̃ s ] , W ̃ ( 2 ) = n eff ( ω ) ∑ s = ± 1 s [ Ω ̃ s + ( 1 / 3 ) Λ ̃ s ] , W ̃ ( 3 ) = − ∑ s = ± 1 Ω ̃ s e 3 i s ϕ ω ( d ) + Λ ̃ s e i s ϕ ω ( d ) , and W ̃ ( 4 ) = n eff ( ω ) ∑ s = ± 1 s Ω ̃ s e 3 i s ϕ ω ( d ) + ( 1 / 3 ) Λ ̃ s e i s ϕ ω ( d ) .

Note that, analogously to the pump, see Eq. (9), the THG field amplitudes satisfy a similar inhomogeneous system, see Eq. (14). However, we emphasize that in Eq. (14) all scattering matrices [ M ̂ TM,TE ( 3 ω , θ ) ] elements are calculated at angular frequency 3ω, as labelled in the superscripts, and account for the THG radiation produced by P _nl within the film. Inverting Eq. (14) we obtain all the THG field amplitudes for either TE or TM polarization of the impinging pump with intensity I ₀ = (1/2)ϵ ₀ c|A ₀|². In Figure 3 we illustrate the dependence of the FW intensity of the THG signal I THG T = ( 1 / 2 ) ϵ 0 c | a t | 2 over the pump wavelength λ, incidence angle θ and intensity I ₀. We observe that, in spite of the reduced THG NL susceptibility χ coll ( 3 ω ) by ≃3 orders of magnitude at λ ≃ 200 nm with respect to λ ≃ 800 nm, see Figure 1c, the THG process is amplified by ≃10² at λ ≃ 200 nm owing to ENZ field enhancement of the pump for TM polarization, see Figure 3a, implying ≃10⁵ THG amplification by the NZI resonance of a Na-based thin film of thickness d = 300 nm. Because the FW THG intensity I THG T ∝ I 0 2 scales with the square of the impinging pump intensity, see Figure 3d and h, this implies an effective averaged field enhancement > 10 by a thin film of thickness d = 300 nm at λ _ENZ. Moreover, for TM polarization we observe THG amplification also for λ < 3λ _ENZ due to NZI resonance of the THG signal field itself, see the dashed white curve in Figure 3a indicating the NZI dispersion relation θ _NZI(λ/3). We emphasize that such signal amplification does not arise from the non-resonant pump field enhancement, but rather from the boost of THG radiation by microscopic dipoles immersed in an ENZ environment. This is mathematically accounted by the Gauss law in local form for the THG signal ∇ ⋅ e _THG = −∇ ⋅ p _nl/ϵ ₀ ϵ(3ω), accounting for the boost of THG polarization volume charge density by ≃3 orders of magnitude at the ENZ resonance, where |ϵ(3ω)| ≃ 10⁻³, see Figure 1b. In addition to NZI resonance, several Fabry–Perot resonances of the THG signal for λ < 3λ _ENZ produce broadband THG amplification ×50 for both TE and TM polarization of the impinging pump, see Figure 3b, c, f, and g. Such a broadband functionality is enabled by the “metal–dielectric” transition [23], [24] of the ENZ material (for THG wavelength λ/3 < λ _ENZ) that, thanks to the reduced absorption of Na, can support Fabry–Perot resonances with quality factor almost comparable to transparent photonic materials. Indeed, we observe a strong dependence of the THG conversion efficiency over the Na-based film thickness d. We investigate systematically the FW THG intensity I THG T produced by the thin film, see Figure 4, where we illustrate its dependence over the thin film thickness d, the pump wavelength λ and the incidence angle θ for both TM and TE polarization of the pump field. Note that enhanced THG at the NZI dispersion relation θ _NZI(λ/3) is heavily modulated by the thin film thickness d for TM polarization of the pump, see Figure 4a and b. Indeed, at λ ≃ 200 nm for every fixed incidence angle θ, we observe an optimal film thickness d ̄ ≃ 20 − 50 nm that decreases when the incidence angle increases. This is due to destroying interference of THG by FW and BW waves, which gets minimised when d is reduced. Moreover, for both TM, see Figure 4a–d, and TE, see Figure 4e–h, pump polarization and λ < 3λ _ENZ, we observe that Fabry–Perot resonances at λ/3 produce optimal THG conversion efficiency at several thicknesses d ̄ n that depend over the incidence angle θ owing to the modified optical path of THG waves. Moreover, we observe optimal THG amplification also for θ > θ _NZI(λ/3) where the Na-based thin film optically behaves as a metal, see the white dashed vertical lines in Figure 4c and g indicating the NZI critical angle θ _NZI(λ/3). Such behavior is not produced by THG resonance, but rather by maximisation of the pump field intensity within the thin film. Indeed, when λ > λ NZI = θ NZI − 1 ( θ ) , see Figures 2 and 3a and e, the Na-based thin film is highly reflective for the pump field. In turn, in such conditions, reducing the film thickness does not reduce the pump intensity within Na, but rather amplifies it due to back-reflection from the z = d interface. In order to account for this effect, we calculate the averaged NL current density per unit length J nl ( Λ ) ( x , t ) = R e j nl e 3 i k x x − 3 i ω t = ∫ 0 d R e − 3 i ω p nl ( z ) e 3 i k x x − 3 i ω t d z , where j _nl is the effective vectorial amplitude of J nl ( Λ ) . In Figure 5a–c, we depict the dependence of |j _nl| over the film thickness d and the pump incidence angle θ and wavelength λ. Note that the d dependence of THG amplification for λ ≃ 800 nm, see Figure 4a, b, e, and f, coincides with the |j _nl| enhancement dependence observed in Figure 5a–c. In turn, THG amplification for λ > 3λ _ENZ at small film thicknesses arises from the enhancement of the effective NL current density per unit length J nl ( Λ ) . The physical origin of such unusual behavior stems from the high reflectivity of Na in such spectral range, producing evanescent pump field within the slab, see Figure 5d, where we illustrate the spatial dependence of the pump field within the film. Note that, when the thickness d decreases, the pump field within the film increases due to enhanced back-reflection of evanescent waves and in turn |j _nl| is increased up to a critical thickness d ̄ ≃ 40 nm, where the field becomes uniform and |j _nl| ≃ |A _in|d starts decreasing due to the reduced film thickness, vanishing for d = 0. Overall, for both TE and TM pump polarization, we observe two distinct incidence angles where resonant THG amplifies the FW signal by ≃10⁴, and a similar behavior is observed also for the BW THG signal.

Figure 3:

THG by Na-based thin films. (a, e) Dependence of the FW intensity of the THG signal I THG T over the pump wavelength λ and incidence angle θ for fixed pump intensity I ₀ = 100 GW/cm². (b, c, f, g) Dependence of I THG T over the pump incidence angle θ for (b, f) fixed pump wavelength λ = 630 nm and several intensities I ₀ = 50, 75, 100 GW/cm², and (c, g) fixed pump intensity I ₀ = 100 GW/cm² and several wavelengths λ = 615, 630, 650 nm [highlighted by the gray shaded area in (a, e)]. (d, h) Dependence of I THG T over the pump intensity I ₀ for fixed incidence angle θ = 11° and several wavelengths λ = 615, 630, 650 nm [highlighted by the gray shaded area in (a, e)]. All plots refer to either (a–d) TM or (e–h) TE polarization of the pump field for fixed thin film thickness d = 300 nm.

Figure 4:

THG manipulation by Na-based thin films. (a, b, e, f) Dependence of the FW intensity of the THG signal I THG T over the thin film thickness d and the pump wavelength λ for incidence angle (a, e) θ = 20° and (b, f) θ = 60°. (c, d, g, h) Dependence of I THG T over the thin film thickness d and the pump incidence angle θ for wavelength (c, g) λ = 630 nm and (d, h) λ = 500 nm. All plots refer to either (a–d) TM or (e–h) TE polarization of the pump field. The white dashed vertical lines in (c, d, g, h) indicate the NZI critical angle θ _NZI(λ/3). All plots are obtained for fixed pump intensity I ₀ = 100 GW/cm².

Figure 5:

Surface-like NL interaction. (a–c) Dependence of the effective NL current per unit length modulus |j _nl| over the thin film thickness d and the pump incidence angle θ and wavelength λ. (a, b) Dependence of |j _nl| over θ and d for fixed pump wavelength λ = 800 nm and (a) TE and (b) TM polarization. (c) Dependence of |j _nl| over λ and d at normal incidence. (d) Spatial dependence of the pump electric field profile modulus within the Na-based thin film rescaled to the impinging amplitude |A _in(z)|/|A ₀| at transverse electromagnetic (TEM) normal incidence for fixed wavelength λ = 800 nm and several thicknesses d = 10, 20, 50, 100, 200 nm. All plots are obtained for fixed pump intensity I ₀ = 100 GW/cm².

Such a surface-like NL interaction arises from thermal nonlinearity, producing THG by 2ω oscillations of the OOE electron temperature. In order to illustrate the ultrafast response of such an NL process, we consider an impinging pump optical envelope A 0 ( t ) = 2 I 0 / ϵ 0 c e − 2 ⁡ ln ⁡ 2 t 2 / τ 2 e ̂ x with peak intensity I ₀ and temporal full width at half maximum (FWHM) τ at the z = 0 position producing within the Na-based thin film both FW and BW waves accounted by the solutions of Eq. (9) in the SVEA. We solve numerically Eq. (4) under such a driving field at every z position by an adaptive step fourth-order Runge–Kutta algorithm, enabling us to calculate the spatio-temporal dependence of the OOE electron temperature, illustrated in Figure 6a for I ₀ = 100 GW/cm², τ = 0.1 ps and carrier wavelength λ = 690 nm. Note that, in such excitation conditions, the OOE electron temperature increases only within the initial propagation distance ≃50 nm, confirming the surface-like NL interaction. In Figure 6b and c we illustrate the temporal evolution of the OOE electron temperature T _e for (b) several distinct pulse duration FWHM τ and fixed peak intensity I ₀ = 100 GW/cm², and (c) several distinct peak intensities and fixed pulse duration FWHM τ = 0.1 ps. We observe that the ultrafast electron heating produces OOE electron temperature increase of only ≃40 K owing to the slow damping rate γ ≃ 24.6 ps⁻¹, producing reduced EM radiation absorption and heating.

Figure 6:

Thermal dynamics in Na-based thin films. (a) Dependence of the electron temperature T _e over the longitudinal position within the thin film z and the excitation time t for fixed pump peak intensity I ₀ = 100 GW/cm² and pulse duration τ = 0.1 ps. (b) Temporal evolution of T _e for fixed pump peak intensity I ₀ = 100 GW/cm², longitudinal position z = 0 and several pulse durations τ = 0.1, 0.5, 1 ps. (c) Temporal evolution of T _e for several pump peaks intensities I ₀ = 0.1, 0.5, 1 TW/cm² and fixed pulse duration τ = 0.1 ps and longitudinal position z = 0. All plots are obtained at TEM normal incidence and for fixed pump carrier wavelength λ = 690 nm.

5 Conclusions

In summary, collision-induced third-order nonlinearity produces efficient THG in Na-based thin films thanks to highly dispersive NL THG susceptibility χ coll ( 3 ω ) ranging from | χ coll ( 3 ω ) | ≃ 1 0 − 24 m²/V² at λ ≃ 200 nm to | χ coll ( 3 ω ) | ≃ 1 0 − 21 m²/V² at λ ≃ 800 nm, where λ is the optical pump angular frequency. Moreover, Na displays great theoretical potential for highly efficient NZI features at λ/3 < λ _ENZ ≃ 230 nm thanks to mitigated absorption producing Im[ϵ(λ _ENZ)] ≃ 10⁻³, about 500 times smaller than ITO. We find that, thanks to such extreme NZI behavior and the excitation of diverse Fabry–Perot THG resonances, it is possible to attain resonant THG over a broad spectral range by tuning the pump incidence angle θ for λ < 3λ _ENZ. NZI resonances produce a giant THG enhancement by ≃10⁴ at peculiar thicknesses d ̄ ≃ 40 nm dictated by optimal pump field penetration within the highly reflective Na-based thin film. For pump intensities I ₀ ≃ 100 GW/cm² at λ ≃ 800 nm, we observe FW THG intensities I THG T ≃ 1 kW/cm². Our results indicate that Na-based thin films are promising for the development of integrated UV sources and innovative spectroscopy schemes at the nanoscale.