Temporal interface in dispersive hyperbolic media

1 Introduction

Electromagnetic wave interaction with time-varying media has recently gained considerable attention and growing interest [1], [2], [3]. This topic offers the notion of four-dimensional (4D) metamaterials [4], wherein material parameters such as relative permittivity can undergo rapid changes in time, either independently or in conjunction with spatial variations. Interest in spatiotemporal modulation, mostly as applied to circuits, has roots dating back to the 1950s [5], [6], [7], [8], [9], yet it has recently garnered significant attention in various research groups worldwide. This newfound interest is mainly due to its promising potential for unconventional wave manipulation and a diverse range of potential applications [1], [2], [4], [10], [11], [12], [13], [14], [15], [16], [17], [18], [19], [20], [21] that include artificial magnetic field for photons [22], optically induced negative refraction [23], frequency conversion [24], amplification [25], [26], Doppler shift [27], [28], Fresnel drag [29], camouflage [30], [31], and nonreciprocity [32], [33], [34], [35], [36], [37] to name a few.

A simple, while fundamental scenario one may devise is a temporal interface. In analogy with a spatial interface between two different semi-infinite media, a temporal interface is defined when a spatially unbounded medium in which a wave propagates is abruptly transformed in time into another medium with different material parameters, e.g., when the relative permittivity of a medium is sharply altered in time [7]. Unlike conventional spatial interfaces, temporal interfaces exhibit three distinctive properties: immutability of momentum (wave vector) accompanied by a frequency shift, lack of electromagnetic energy conservation, and generation of a backward wave, which, due to causality, propagates in the medium after the temporal interface [7], [38], [39]. These features have been experimentally validated, for example, in plasma physics [40], [41], with water waves [42], and in transmission lines that operate at megahertz frequencies [43]. Also, according to these fundamental characteristics and considering complex electromagnetic systems, such as anisotropic and bianisotropic media [44] and metasurfaces [45], a myriad of possibilities and opportunities for manipulation of classical and quantum fields have been uncovered. These advancements encompass: the creation of “wiggler mode” [46], [47], [48], temporal aiming [12], direction-dependent wave manipulation [49], inverse prism [50], antireflection temporal coatings [51], [52], [53], polarization engineering [49], [54], [55] and polarization-dependent analog computing [56], wave freezing and thawing [57], [58], the transformation of surface waves into free-space radiation [57], [59], photon-pair generation [60], [61], angular-dependent inhibition of photon production [62], photon-pair destruction, and vacuum state generation [61].

Despite the growing body of research in this area, the sudden creation of hyperbolic media, resulting in a temporal interface, remains largely unexplored. Here, we extend the notion of temporal interface to dispersive hyperbolic media. Specifically, we explore how a monochromatic electromagnetic uniform plane wave undergoes changes when a host medium, assumed to be a simple isotropic dispersionless dielectric (e.g., free space), is rapidly transformed into a dispersive hyperbolic medium, formed by a stack of many bilayers made of metal and dielectric layers. The frequency dispersion of the hyperbolic medium is taken into account by considering the Drude dispersion for those metal layers. In particular, we reveal that a temporal interface in the presence of such anisotropy and frequency dispersion causes the splitting of the initial frequency of the wave into three pairs, which propagate mainly along the optical axes of the crystal, exhibiting canalization.

The paper is organized as follows. In Section 2, we discuss the mechanism of the corresponding temporal interface. In Sections 3 and 4, we explain the frequency conversion and the evolution of electromagnetic fields as the result of the temporal interface. In Section 5, we demonstrate numerical simulation results, and, finally, in Section 6, we conclude the work.

2 Description of the problem

To start, let us consider our hyperbolic medium as an infinitely extended collection of identical bilayers, each formed by a dielectric layer of thickness d_d and relative permittivity ϵ_d and a metallic layer of thickness d _m and relative permittivity ϵ _m . We assumed all these bilayers are parallel with the xz − plane of a Cartesian coordinate system, with its y axis being normal to these bilayers. According to the effective medium theory, the elements of the relative permittivity tensor of such a medium can be written [63] as

with f = d _m /(d_d + d _m ). Let us assume the relative permittivity of the dielectric layers not frequency dispersive, and the relative permittivity of metallic layers lossless Drude-dispersive,

where ω_p is the plasma frequency and ϵ_∞ is the relative permittivity at infinite frequency. Plugging Eqs. (2) into (1a) and (1b) and assuming ϵ_d = ϵ_∞, we can get frequency-dependent expressions for ϵ _xx , ϵ _yy , and ϵ _zz ,

where

Dispersion in Eq. (3a) is of the Drude type, whereas dispersion in Eq. (3b) is of the Lorentzian type, with resonance frequency ω_0,yy (Eq. (4b)). Interestingly, and as expected, one can engineer effective media parameters by properly selecting the relative thicknesses of the layers.

In the equations above, as mentioned, it is assumed that ϵ_∞ = ϵ_d. This choice ensures that, when the plasma frequency is zero, we have an isotropic nondispersive effective medium ϵ_d (see Eqs. (3) and (4) when ω_p = 0). Clearly, as shown, introducing a nonzero plasma frequency (ω_p ≠ 0) imparts dispersive and anisotropic properties to the effective medium. Thus, a sudden change in plasma frequency presents a unique opportunity to establish a temporal interface between two markedly distinct media. This is illustrated by Figure 1.

Figure 1:

Temporal interface between a simple isotopic medium and a rapidly generated dispersive anisotropic medium: different colors of the waves indicate different frequencies after temporal jump. ⟨S₁⟩, ⟨S₂₁⟩, ⟨S₂₂⟩, ⟨S₂₃⟩ denote the time-averaged Poynting vectors of the corresponding FW waves. We do also have time-averaged Poynting vector for the BW waves (not shown here for the sake of brevity).

3 Frequency conversion due to temporal interface

For simplicity, let us initially consider a monochromatic uniform plane wave with angular frequency ω₁ propagating inside a medium described by Eq. (2) but initially parameterized with a zero plasma frequency, which makes the medium a simple isotropic and dispersionless dielectric. At t = t₀, we abruptly ionize the metallic layers by increasing ω_p from zero to, say, 2ω₁, making the medium both anisotropic and frequency dispersive. Frequency 2ω₁ is chosen to show a large enough contrast of the effects, and, in principle, it can be chosen arbitrarily. The dispersion relation for stationary (i.e., time-invariant) anisotropic media with the electric field components in the x − y plane and the magnetic field parallel with the z axis reads

where the wave vector components k _x and k _y should be conserved quantities across the temporal boundary at t = t₀. In the absence of material dispersion, this conservation of wave vectors allows us to find, from Eq. (5), the new converted frequency, ω₂, after t₀ [12]. In the case of frequency dispersive isotropic media, one obtains dispersion relations from a transcendental equation, as discussed in [64]. In the present work, however, we have a combination of anisotropy and frequency dispersion, requiring the medium to be characterized by frequency-dependent ϵ _xx and ϵ _yy , leading to another transcendental equation, which reduces to

Due to anisotropy, the converted frequencies after t₀, ω₂, depend (through k _x and k _y ) on the initial incidence angle of propagation ψ. Equation (6) can, therefore, be recast as the following 6th-order characteristic equation

In our lossless scenario, this expression can be made 3rd order by considering ω 2 2 as the new variable, with each of the three solution pairs ±ω₂ indicating a forward- and a backward-propagating wave, which hereafter will be denoted as FW and BW waves, respectively.

It is interesting to note that a temporal interface in an unbounded dispersionless media (either isotropic or anisotropic) only shifts the incident frequency, resulting in a single pair of FW and BW waves with a single converted frequency (see, e.g., [7], [12]). It has also been shown that in lossless isotropic medium with Lorentz dispersion, one obtains two positive solution pairs [64], providing two pairs of FW and BW waves with two converted frequencies. In our specific scenario here, which involves anisotropy and frequency dispersion, we observe the emergence of three solution pairs. Figure 2(a) shows the polar plots of converted frequencies after the temporal jump as a function of the direction of propagation angle of the initial wave ψ. Propagation along the axes is equivalent to propagation in a medium with the corresponding Lorentz or Drude dispersion. Specifically, propagation along x axis (ψ = 0 and ψ = π) or y axis (ψ = π/2 and ψ = 3π/2) is equivalent to propagation in a medium with permittivity ϵ _yy (see Eq. (3b)) or ϵ _xx (see Eq. (3a)), respectively. Figure 2(a) shows three frequencies for all ψ except of ψ = mπ/2 with m is arbitrary integer number, where the amplitude of one of the three pairs of FW and BW waves is zero. More will be said below.

Figure 2:

Converted frequencies and relative permittivity values after temporal interface: (a) polar plots of the normalized converted frequencies versus incidence angle, after sudden temporal transition from isotropic medium to anisotropic hyperbolic medium, assuming d_d = d _m , ω_p = 2ω₁, and ϵ_∞ = 1. (b) and (c) Effective relative permittivity in x and y directions for different converted frequencies in (a): (b) hyperbolic medium ω_2,1, (c) hyperbolic medium ω_2,2, (d) elliptic medium at ω_2,3. Note that the circles represent vanishing coupling to the corresponding mode.

Figure 2(b)–(d) shows the xx and yy elements of relative permittivity tensors of the medium after the temporal interface, as evaluated for each of the three converted frequencies. (The ϵ _zz is the same as ϵ _xx .) It is evident that the medium retains its anisotropic nature; however, for each of the frequencies, the nature of anisotropy is different. For ω_2,1 and ω_2,2, the components of the permittivity tensor exhibit opposite signs, indicating a hyperbolic nature of the medium. In contrast, at frequency ω_2,3, all components of the permittivity tensor are positive (and less than unity), implying that the anisotropy at this frequency is of elliptic type. Additionally, at frequencies ω_2,1 and ω_2,2, permittivity tensors have different negative components, i.e., ϵ_xx,1 and ϵ_yy,2 are negative for ω_2,1 and ω_2,2, respectively, indicating that the isofrequency curves of the material for these two frequencies are rotated by π/2.

3.1 Band topology

From the previous, one can see that ω_2,2(ψ = π/2) and ω_2,3(ψ = 0) result from pure Drude and Lorentz dispersion, respectively, and thus ω_2,1(ψ = 0) is the second solution pair expected from a Lorentzian response (as shown in [64], ω_2,3 > ω₁ can be connected to the resonant frequency ω_0,yy > ω₁, whereas ω_2,1 < ω₁ approaches the solution in a nondispersive time interface [7]). On the other hand, ω_2,1(ψ = π/2) = 0 gives us the missing “wiggler” mode from the Drude plasma. For arbitrary ψ, these three pure frequencies are hybridized into three different solutions, considering that sum of them squared is a conserved quantity independent of ψ, reduced when ϵ_∞ = 1 to:

(8) ω 2,1 2 + ω 2,2 2 + ω 2,3 2 = ω 1 2 + ( 1 + f ) ω p 2 .

The dispersion diagram of the new homogenized material in Figure 3(a) reveals three distinct frequency bands. (i) A hyperbolic lower band in the range [0, min(ω_0,yy, ω_p,eff)]. (ii) A mid band in the range [ω_p,eff, ω_p], which is elliptic (hyperbolic) below (above) resonance ω_0,yy. (This requires that ω_0,yy > ω_p,eff, only satisfied if d _m < d _d ; otherwise, if ω_0,yy < ω_p,eff, a bandgap opens between the first and the second bands. In Figure 3(a), d _m = d _d , so this band is only hyperbolic and there is no gap underneath. Coincidentally, moreover, ω_0,yy = ω_p,eff equate the surface plasmon frequency ω p / 2 of the actual metallic layers.) (iii) An elliptic upper band above ω_p. In this figure, the cylinder’s surface encompasses all the possible angular directions of the initial wave in vacuum, so its intersection with the three bands provides the necessary momentum matching condition (such three intersections (blue, red, and yellow lines) are precisely the curves in Figure 2(a)). Moreover, the intersection of the dispersion surface with the horizontal planes ω = ω_2,l (chosen here for ψ = π/4) yields the three corresponding isofrequency contours (black curves).

Figure 3:

Dispersion topology of the homogenized hyperbolic material. (a) Dispersion surface of the three bands when d _m = d _d : from momentum conservation, their intersection with ω₁ light in vacuum (cylindrical surface) reveals the excited states after the temporal interface; corresponding isofrequency plots are indicated in black for ψ = π/4. (b) Continuum of isofrequency contours (d _m = d _d /2) that correspond with the blueshifting evolution of the excited states ω_2,l as the plasma frequency goes from 0 to 2ω₁. Each ω_2,l renders a different surface, with its own colormap. The black isofrequency contours indicate a transition from hyperbolicity to ellipticity within the second band (ω_2,2). The zenithal view of the inset shows how the forbidden directions of the first and second bands are perpendicular. (c) Coalescence of excited states from bands 2 and 3 when ψ = π/2 (left), and from bands 1 and 2 when ψ = 0 (right). This happens when f = 1 − κ² (e.g., d _m = 3d _d for ω_p = 2ω₁) and f = 1 + 2 κ 2 − 1 + 4 κ 4 2 , respectively, with κ = ω₁/ω_p. (d) Time evolution of the FW amplitudes when a perfectly abrupt step-function temporal interface is replaced by a deeply subcycle sigmoid function, with ψ = 1°: during this transition, the matched wave in band 2 (red) first feels a hyperbolic medium and grazes the corresponding forbidden direction, so it is barely coupled to. As ω_p keeps increasing, the medium becomes elliptic in the new excited portion of band 2 (out-of-plane dielectric tensor element shown in the inset), and the wave amplitude |A₂| begins to increase before tunneling out to band 1 (this tunneling effect is explained in Section 5).

Instead of the material dispersion in its final state, Figure 3(b) superimposes the isofrequency contours that correspond to each of the ω_2,l eigenstates, as they increase with increasing plasma frequency (from 0 to its final value after the transition) for a fixed propagating angle (ψ = π/4). Each of the three surface contours ω_2,l is plotted with a different colormap and shows how the second band makes a transition from hyperbolic to elliptic in its excited frequency ω_2,2 when d _m < d _d (we here choose d _m = d _d /2). In the limit of ω_p → ∞, ω 2,1 → [ 0 , 1 − f ω 1 ] depending on ψ, ω 2,2 → f ω p , and ω_2,3 → ω_p, which satisfies Eq. (8).

Now that the interplay between ω₁/ω_p and f has been visualized, further insight is gained by noticing that ω_2,2(ψ = 0) = ω_p,eff, i.e., at ψ = 0 this mode of zero field amplitude becomes the effective bulk plasmon resonance ϵ_xx,2 = 0, which is thus not coupled onto. Likewise, ω_2,3(ψ = π/2), of zero field amplitude also, becomes ω p,eff 2 + ω 0 , y y 2 (which is simply ω_p when ϵ_∞ = 1) and, therefore, ϵ_yy,3 = 0. However, the described behavior holds as long as f < 1 − κ² (below resonance, i.e., ω₁ < ω_0,yy) and f > 1 + 2 κ 2 − 1 + 4 κ 4 2 , κ being ω₁/ω_p. Right at the upper threshold for f, ω_2,2 and ω_2,3 become degenerate and equal to ω_p when ψ = π/2 (Figure 3(c), left panel). Above this threshold, the incident wave does couple to ω_2,3 when ψ = π/2, but instead there is no coupling to ω_2,2 when ψ is either 0 or ψ = π/2. In the limit f → 1, coupling to ω_2,2 vanishes for all directions: the Lorentz mode ω_2,3 plays the role of an effective isotropic plasma, since ω_0,yy → 0, and ω_2,1 → 0 describes the associated DC mode for all ψ. Similarly, at the lower threshold, ω_2,1 and ω_2,2 coalesce into ω_p,eff (Figure 3(c), right panel). Below this threshold, when ψ = 0, there is no coupling to ω_2,1 but there is to ω_2,2. Taking the limit f → 0 (no temporal boundary at all), there is only FW coupling to ω_2,2 → ω₁, regardless of ψ, given that

(9) ω 2,1 → 0 , ω 2,2 → ω 1 , ω 2,3 → 2 ω 1 .

In short, the incoming wave feels vacuum at ω_2,2 in the new medium parameterized by ω_p,eff → 0 and ω_0,yy → 2ω₁, again assuming ϵ_∞ = 1. This vacuum limit is different if, e.g., one recovers d_d = d _m and makes ω_p → 0, in which case

(10) ω 2,1 → 0 , ω 2,2 → 0 , ω 2,3 → ω 1 .

In such scenario, the incoming wave couples to the new medium – where the resonance frequency ω_0,yy now also vanishes – through ω_2,3. These limiting cases boil down to the vanishing width of the first and second bands as f → 0 and f → 1, respectively.

From Figure 3(a) and (b), propagation of the ω_2,2 eigenmode is forbidden along the x direction in its hyperbolic phase. In order to see this elliptic-to-hyperbolic transition, we choose d _m = d _d /2 and ψ = 1° and solve the temporal transition in the eigenstate basis versus time, assuming a subcycle sigmoid function of transition rate R for the parameterization of the time-varying ω_p. As seen in Figure 3(d), only after the transition of ϵ_yy,2 from positive to negative, at t ≈ − 2/R (see inset), does the FW mode amplitude |A₂| (red dashed) begin to increase (a detailed derivation of the complex field amplitudes is found in the next section).

4 Field amplitudes after temporal interface

To get a better insight of the wave phenomena after the abrupt emergence of anisotropy, one needs to determine amplitudes of the fields after the temporal jump. The fields after an abrupt temporal interface consist of the three pairs of FW and BW waves. To determine their six amplitudes, we need to have six independent temporal boundary conditions. From our homogenized model of anisotropy in Eq. (2), the second-order nature of the time-varying Drude and Lorentzian response in the form of two separate second-order differential equations for the x and y components of the polarization density vector P and electric field vector E can be written as

where ϵ_∞ is assumed to be unity. Therefore, temporal continuity of P _x and P _y and their time derivatives at the temporal interface is implied. Together with temporal continuity of electric and magnetic flux densities, D = ϵ₀E + P and B = μ₀H, we have a set of six temporal boundary conditions. Worth mentioning that continuity of D _x and D _y yields the same equations (as an aside, it is worth noting that the temporal continuity of P and D automatically guarantees the temporal continuity of E). Let us assume that before the temporal jump at t = t₀, we have

where subscript “1” indicates the quantities before t₀. For the sake of simplicity, let us assume t₀ = 0, then quantities after t₀ are denoted with subscript “2” and they can be written as

Using equations in (12) and (13) and the six temporal boundary conditions, one can write the following expressions:

Here, each line of the matrix from top to bottom corresponds to continuity of B _z , D _x (or D _y , same equation), E _x (or P _x ), E _y (or P _y ), d P x d t and d P y d t , respectively. Solving this system of equations provides field amplitudes for FW and BW waves. Figure 4(a) and (b) shows linear plots of field amplitudes as a function of the direction of propagation angle ψ. Figure 4(c) and (d) shows the amplitude and direction of the time-averaged Poynting vector of FW and BW for each converted frequency. Directions of arrows indicate the direction of the Poynting vector. Colored contours and length of the arrows indicate the amplitude of the Poynting vector of the corresponding frequency branch.

Figure 4:

Field amplitudes and time-averaged Poynting vectors after the temporal interface: field amplitudes for (a) FW waves, A _l , and (b) BW waves, B _l ; (c) time-averaged Poynting vector of FW for each converted frequency as well as the initial Poynting vector divided by 2; (d) time-averaged Poynting vector of BW for each converted frequency as well as the initial Poynting vector divided by 10. (Since the Poynting vectors for FW and BW waves are small compared to S₁, in order to show their details, we scale down S₁ in the figure.)

From Figures 2 and 4, one can make several observations: First – for the limiting cases ψ = 0 and ψ = π/2 (and obviously also for ψ = π and ψ = 3π/2), only two sets of wave have nonzero amplitudes, since A₂ = B₂ = 0 at ψ = 0 (and ψ = π) and A₃ = B₃ = 0 at ψ = π/2 (and ψ = 3π/2). This results in zero time-averaged Poynting vector S_2,2 for FW and BW waves for ψ = 0 (and ψ = π) and S_2,3 for ψ = π/2 (and ψ = 3π/2). For all other angles of propagation, we always get three pairs of FW and BW waves. Second – the first set for ω_2,1 is always nonzero, i.e., A₁ ≠ 0 and B₁ ≠ 0 for all ψ. Third – for the first converted frequency branch, while ϵ_xx,1 becomes negative (and infinitely large for ψ = π/2 (and ψ = 3π/2)), the corresponding wave always exists and its converted frequency ω_2,1 approaches zero for ψ = π/2 (and ψ = 3π/2), which essentially means DC field. However, we note that while the field amplitudes A₁ and B₁ at angle ψ = π/2 (and ψ = 3π/2) are not zero, the x component of the electric field is zero since the expression ω_2,1ϵ_xx,1 in the denominator approaches infinitely large value as ω_2,1 approaching zero at ψ = π/2 (and ψ = 3π/2) (see Eqs. (13a) and (13b)). Effectively for ψ = π/2 (and ψ = 3π/2), we have a Drude medium, where the first solution is a DC (effectively “frozen”) magnetic field (while the xx element of relative permittivity for this set is infinitely negative, thus behaving as a perfect electric conductor causing the x component of the electric field to be zero) and the second solution is a propagating wave. Fourth, it is worth noting that for waves propagating along ψ = π/2 (and ψ = 3π/2), the medium after the temporal jump behaves as a Drude medium (with A₁ ≠ 0, B₁ ≠ 0, A₂ ≠ 0 and B₂ ≠ 0), but for waves propagating along ψ = 0 (and ψ = π), the medium behaves as a Lorentzian medium (with A₁ ≠ 0, B₁ ≠ 0, A₃ ≠ 0 and B₃ ≠ 0). For all other angles, after the temporal jump, the medium behaves as a medium with a dispersion not resembling solely Drude or Lorentzian type, with all three pairs of amplitudes A₁ ≠ 0, B₁ ≠ 0, A₂ ≠ 0, B₃ ≠ 0, A₃ ≠ 0, and B₃ ≠ 0. Finally, the Poynting vector in Figure 4(c) and (d) demonstrates another interesting insight on the wave properties. One notices that for all three converted frequencies, the corresponding energy flows mainly along the optical axes, which can be attributed to a certain level of canalization of energy in this system. We point out that along the directions ψ = π/2 and ψ = 3π/2, the Poynting vector for the converted frequency ω_2,1 (which is near zero) is identically zero (since the x component of the electric field is zero for this set). Therefore, the blue-color plots of the Poynting vector in Figure 4(c) and (d) have a zero value along those directions, but attain nonzero values when the direction of propagation deviates from ψ = π/2 and ψ = 3π/2 with its flow being primarily along the x axis.

5 Numerical simulation

One may argue that, rigorously speaking, the continuity of the homogenized dP/dt adopted above might still be, at least in principle, open to debate for the following reasons: (i) it is true that, within each metal layer with the Drude dispersion, a temporal discontinuity in the plasma frequency leads to temporally continuous nonhomogenized dP/dt when the Drude current response follows; (ii) it is also clear that Eq. (11) correctly models the homogenized anisotropic response once ω_p (and thus ω_p,eff and ω_0,yy) is time-invariant. But the assumption that these same equations are still indeed valid across the time interface, when ω_p(t) varies in time, should be tested. To the best of our knowledge, a theory of time-varying homogenization in the presence of dispersion has not been done yet. If Eq. (11) were not to be exact (in the following we prove they are), different temporal boundary conditions might be applicable for the homogenized quantities. In light of this, in order to validate our analytic derivations, we numerically solve homogenized model described by Eq. (11) and the actual nonhomogenized problem with deeply subwavelength Drude layers. We present the results only for the latter case, since the results for the former case of homogenized model are identical with the analytical results, if the layers are thin enough (∼λ/10³ in our case). Regardless of the model, the in-plane momentum k _x allows for a dimensionality reduction, so only y needs to be parameterized. We can thus write

with phasor H ̃ z ( y , t ) = e j ω 1 t e − j k y y before the temporal transition, in an unbounded vacuum. We can thus focus on the time evolution of the k-phasors, which abides by the curl equations

and the nonhomogenized time-varying current response

where ω p 2 ( y , t ) is nonzero only inside the Drude layers, depending on y. Using finite differences in both y and t and a Yee grid, Eqs. (16) and (17) can be solved in leapfrog fashion [65]. Figure 5 shows the temporal evolution of the fields at x = y = 0 when ψ = π/4, with perfect overlap of the analytic and numerical curves, demonstrating that our assumption of temporal continuity used above for the homogenized quantities is indeed valid. An exceedingly small unit cell of a = λ/10³ (y step of a/200) is chosen in order to keep the spectral content within the long-wavelength portion of the first band in the dispersion diagram, such that the homogenization behind our analytic expressions is applicable. (Anyhow, regardless of a, we should point out that out-of-plane momentum k _y is partially conserved too, in that H ̃ z ( y + a , t ) = e − j k y a H ̃ z ( y , t ) at all times. That is, the crystal momentum – i.e., the portion of k _y within the first Brillouin zone – is conserved, and any temporal variation of this crystal, as long as it is structured in y with a periodicity, will only contribute to the total out-of-plane momentum with integer multiples of the reciprocal lattice vector 2π/a. The solution could thus be expanded as a Bloch wave, or else, as in our case, one can simply simulate one unit cell in the finite-differences solver by imposing the Bloch boundary condition.) Accordingly, the numerical curves for E _y and P_x/y, discontinuous across the successive air/metal interfaces, follow from y-averaging over a few unit cells. Figure 5(d) illustrates these discontinuities in E _y (t = 6T) near y = 0.

Figure 5:

Comparison of analytical and numerical results. (a–c) Show magnetic field (H _z ), x and y components of electric field (E _x and E _y ), and x and y components for polarization density (P _x and P _y ) obtained analytically and numerically. Insets in (a–c) highlight continuity of the plotted quantities across the temporal interface. (d) Electric field E _y along y axis when x = 0 ( e − j k x x dependence, as seen in Eq. (15)), which exhibits discontinuities, and its homogenized version ⟨E _y ⟩ (resulting from y-averaging over a distance of λ/100 or, equivalently, about 10 unit cells) agreeing with effective medium theory.

Finally, the effects of the finite switching time and loss have been studied and reported in Figure 6. Panels (c) and (d) illustrate how, as we reduce the rate of change R of the plasma frequency, which follows an arbitrarily fast (though continuous) sigmoid – particularly, the results in panel (c) are already very close to the previous analytic results for a perfectly abrupt time boundary – both temporal reflection (BW terms) and frequency splitting tend to decrease, as expected from the adiabatic limit picture. In contrast, evidently, the allowed frequencies across the transition (panel (a)), and therefore their associated permittivity tensors (panel (b)) are only function of ω_p and not of the rate R. Panel (e) depicts the temporal decay of the waves, following the imaginary part of the eigenfrequencies, when the Drude model of the metallic layers in Eq. (2) is extended with a damping term γ.

Figure 6:

Impact of finite switching times and loss. (a) Evolution of the three eigenfrequencies as ω_p rises from 0 to 2ω₁ versus time (normalized to the inverse of the sigmoid’s transition slope R), assuming ψ = 89° and d _m = d _d . When ω_p satisfies f = 1 − ( ω 1 / ω p ) 2 , ω_2,2 and ω_2,3 become nearly degenerate (see black circle in the figure) and strongly interact: given that the initial mode ω_2,3 in its final state is short-circuited by the metal–air interfaces, complete interband transition takes place from the third to the second level, as shown in panels (c) and (d). (b) Dielectric tensor elements of modes 2 and 3: at the degeneracy point, both out-of-plane permittivities ϵ_yy,2 and ϵ_yy,3 (dashed lines) become epsilon-near-zero (ENZ). (c) and (d) Temporal change of the FW (solid lines) and BW (dashed lines) field amplitudes for two different switching times: quasi-instantaneous (R = 100/T, panel (c)) and slow (R = 1/T, panel (d)). As R decreases, so does temporal reflection, with backward propagation tending to disappear. Still, the longer time of interaction between the two degenerate modes 2 and 3 favors frequency mixing down to mode 1 (DC), for which FW or BW waves lose their propagation character: we rather have a stationary pattern (dictated by the original wavenumber) that slowly oscillates in time. (e) Effect of adding some loss to the Drude model of the metallic layers: the evolution of the magnetic field phasor (x = y = 0) for several values of the nonradiative electron collision rate γ shows, over the first periods, very similar results to the ideal lossless material discussed so far. The decaying character of the (now complex) frequencies ω_2,l is nonetheless very clear after a few periods in the lossiest case, for which the actual metallic layer has a dielectric function of roughly −0.02i at the ENZ condition.

6 Conclusions

In summary, we have extended the notion of temporal interfaces to hyperbolic frequency dispersive media. Particularly, our analysis focused on a temporal interface between vacuum and a layered structure that exhibits hyperbolic properties. The interplay between anisotropy and frequency dispersion in this system was analyzed, which, together with a temporal interface, resulted in the splitting of the original wave into three different frequency pairs. We discussed the topological features in the dispersion of the effective time-varying anisotropic medium and its final properties at the new (converted) frequencies and noticed that it leads to energy canalization, i.e., propagation of power primarily along optical axes of the structure. Finally, full-wave simulations were conducted to corroborate our theoretical findings.

On the experimental side, several implementation scenarios can be envisioned. One possibility involves a periodic arrangement of dielectric and semiconductor layers, optically pumped to induce a temporal interface. The dielectric and semiconductor materials must be carefully selected so that the optical pumping leaves one medium unaffected while reconfiguring the other. As an alternative to the optical approach, one could consider a low-frequency, two-dimensional geometry – such as a metasurface. For instance, a hyperbolic metasurface operating at radio frequencies could be realized using lumped elements and switches, providing a feasible platform for proof-of-concept experiments.