# Eigenpulses of dispersive time-varying media

S. A. R. Horsley<EMAIL_ADDRESS>School of Physics and Astronomy,
University of Exeter, Stocker Road, Exeter, UK, EX4 4QL  E. Galiffi
<EMAIL_ADDRESS>Photonics Initiative, Advanced Science Research Center,
City University of New York, 10027 New York, USA  Y.-T. Wang  Photonics
Initiative, Advanced Science Research Center, City University of New York,
10027 New York, USA

###### Abstract

We develop a compact theory that can be applied to a variety of time-varying
dispersive materials. The continuous wave reflection and transmission
coefficients are replaced with equivalent operator expressions. In addition to
comparing this approach to existing numerical and analytical techniques, we
find that the eigenfunctions of these operators represent pulses that do not
change their spectra after interaction with the time-varying, dispersive
material. In addition, the poles of these operators represent the non–time
harmonic bound states of the system.

Typical electromagnetic parameters (e.g. refractive index, or impedance) are
constant in time, deriving from both the dispersive response and _spatial_
arrangement of the constituent atoms, or meta–atoms. Designing such composite
materials (metamaterials) with specified behaviour has been the subject of
intense research for decades, with recent developments including one–way
propagation [1], PT–symmetric lasing [2], and a variety of methods for
invisibility cloaking [3]. Yet many fundamental limitations can be overcome if
the material properties are additionally structured in _time_ as well as space
[4]. For example, a lossless time-varying grating can amplify waves [5],
something impossible with a static grating. The Kramers–Kronig relations [6]
also restrict the available palette of static materials (limiting e.g. the
thickness to bandwidth ratio of an absorber [7]), but do not apply to
time–varying media (TVM) [8]. In addition, TVM exhibit fundamentally new
phenomena such as ‘time-refraction’ [9, 10], ‘time-reflection’ [11, 12, 13,
14], and ‘temporal aiming’ [15], where wave energy can be re–directed in the
absence of spatial inhomogeneity.

While many of the properties of TVM are yet to be experimentally explored, one
promising platform is the conducting compound Indium Tin Oxide (ITO) [16].
Close to the epsilon–near–zero (ENZ) frequency of ITO, it exhibits a large
non-linear susceptibility [17], leading to an effective permittivity that is
switchable on a sub–picosecond timescale [18]. This rapid switchability has
been used to demonstrate time-refraction [9, 10] and time-diffraction [19] as
well as to implement both time-varying metasurfaces [20], and spectrum-
modifying mirrors [21].

As most materials are static, the theoretical tools for treating TVM are not
well developed. For instance, typical TVM are highly dispersive in the
frequency regimes where their nonlinearities are strongest, which makes the
susceptibility a two–time function. This is difficult to incorporate into
existing analytical results or numerical schemes. Several works have found
methods to circumvent these difficulties. For instance, Zhou et al. [9]
introduced an effective non-dispersive time-varying response that was
optimized to fit the experimental data. Many current predictions alternatively
integrate Maxwell’s equations in full over time using e.g. COMSOL Multiphysics
[22], and such an approach can be found in e.g. [15]. However this becomes
computationally intensive for the dispersive structured materials discussed in
e.g. [20, 21].

In this Letter we develop a compact semi–analytical theoretical approach that
can be applied to a variety of time varying dispersive materials and is
particularly well suited to describe the thin film ITO experiments discussed
above. We find that the continuous wave reflection and transmission
coefficients can be replaced with equivalent operator expressions that are
simple to evaluate numerically and act on the spectrum of the incident wave.
Although distinct from this work, aspects of our operator-based approach
(where e.g. the wave–vector is treated as an operator) appear in the recent
extension of Mie theory given by Ptitcyn et al. [23]. We find that the
eigenfunctions of these operators represent spectra of incident waves that are
e.g. unchanged on reflection from a TVM. Furthermore when the eigenvalue of
e.g. the reflection operator is zero or infinite, we have respectively the
time-varying analogues of a reflectionless medium or a bound state, which we
term the ‘eigenpulses’ of the system. In addition, we compare this approach to
existing numerical and semi-analytical techniques, evaluating the superior
performance of this method in terms of both accuracy and efficiency.

For a static material, the electric current $\bm{j}$ is linked to the electric
field $\bm{E}$ through the conductivity $\sigma(t-t^{\prime})$ that represents
the movement of charge in response to the past behaviour of the electric field
and depends only on the time _difference_ $t-t^{\prime}$. When the material is
explicitly time-dependent, due to e.g. a pump pulse [17] (at optical
frequencies) or electronic modulation (at radio frequencies), the conductivity
can be replaced with a two–time function such that,

$\bm{j}(t)=\int_{-\infty}^{\infty}{\rm
d}t^{\prime}\sigma(t,t-t^{\prime})\,\bm{E}(t^{\prime}).$ (1)

As in the static case, causality requires $\sigma(t,t-t^{\prime})=0$ when
$t^{\prime}>t$. We understand $\sigma(t,t-t^{\prime})$ as the response at time
$t$ to the electric field at the earlier time $t^{\prime}$, where this
response is modulated as a function of the observation time $t$. We could also
develop the same formalism by taking the first argument of $\sigma$ as
$t^{\prime}$ instead of $t$. As the two times are related by
$t^{\prime}=t-(t-t^{\prime})$, our results can be applied to either form of
Eq. (1), with only minor modifications. We could also equally develop the
formalism in terms of the permittivity and/or permeability instead of the
conductivity.

Performing a Fourier transform of (1), $\tilde{\bm{j}}(\omega)=\int\,{\rm
d}t\,\bm{j}(t)\exp({\rm i}\omega t)$, the frequency dependent current can be
written as

$\displaystyle\tilde{\bm{j}}(\omega)$
$\displaystyle=\int_{-\infty}^{\infty}\frac{{\rm
d}\omega^{\prime}}{2\pi}\tilde{\bm{E}}(\omega^{\prime})\int_{\infty}^{\infty}{\rm
d}t\,\sigma(t,\omega^{\prime}){\rm e}^{{\rm i}(\omega-\omega^{\prime})t}$ (2)
$\displaystyle=\int_{-\infty}^{\infty}{\rm
d}\omega^{\prime}\tilde{\bm{E}}(\omega^{\prime})\hat{\sigma}(-{\rm
i}\partial_{\omega},\omega^{\prime})\delta(\omega-\omega^{\prime})$
$\displaystyle=\hat{\sigma}(-{\rm
i}\partial_{\omega},\omega)\tilde{\bm{E}}(\omega)$ (3)

where $\hat{\sigma}$ is the operator obtained by replacing the first argument
$t^{\prime}$ with the operator $-{\rm i}\partial_{\omega}$. In the final line
of (3), all derivatives $\partial_{\omega}$ within $\hat{\sigma}$ must be
ordered such that they appear to the _left_ of all the frequency dependence of
$\hat{\sigma}$ 111The opposite (normal) operator ordering (where the
$\partial_{\omega}$ sit on the _right_ of the factors of frequency in
$\hat{\sigma}$) is required if we specify the conductivity in Eq. (1) as
$\sigma(t^{\prime},t-t^{\prime})$.. This prescription is reminiscent of the
anti–normal operator ordering adopted in quantum mechanics [25]. To use Eq.
(3) we write the operator as e.g. $\hat{\sigma}=\sum_{n}a_{n}(-{\rm
i}\partial_{\omega})b_{n}(\omega)$. The derivative $\partial_{\omega}$ is
numerically constructed as an $N\times N$ matrix acting on $N$ frequency
points, via the finite difference approximation or a Fourier transform. The
operators $a_{n}$ are then evaluated as matrix valued functions, and
$b_{n}(\omega)$ is a diagonal matrix. This idea of using an operator valued
function is similar to the exponential function of the Hamiltonian operator
used as the time evolution operator in quantum mechanics [26].

In this work we assume the magnetic permeability is unity, and use the Drude
model with a time-varying plasma frequency $\omega_{p}$ (see Supplementary
Information)

$\hat{\sigma}(-{\rm i}\partial_{\omega},\omega)=\omega_{p}^{2}(-{\rm
i}\partial_{\omega})\frac{{\rm i}\epsilon_{0}}{\omega+{\rm i}\gamma},$ (4)

where we have imposed the aforementioned anti–normal ordering, and $1/\gamma$
is the collision time. Note that throughout this work we use the symbol
$\omega_{0}=\omega_{p}(-\infty)$, i.e. the plasma frequency before the time
variation.

The simplest application of Eq. (3) is where the medium is homogeneous,
propagation is along $x$, and the field is polarized such that
$\tilde{\bm{H}}=H\bm{e}_{z}$. In this case Maxwell’s equations become
$\partial_{x}^{2}H+\hat{K}_{p}^{2}H=0$, the solutions to which are

$H(x,\omega)=\exp(\pm{\rm i}\hat{K}_{p}x)\,H_{0}(0,\omega)$ (5)

where $k_{0}=\omega/c$ and $\hat{K}_{p}$ is the matrix square root of
$\hat{K}^{2}_{p}={\rm i}\eta_{0}\hat{\sigma}k_{0}+k_{0}^{2}$ 222It is
important to determine which sign of the matrix square root $\pm\hat{K}_{p}$
describes right–going waves. We take right–going waves to be those where the
eigenvalues have a positive imaginary part. and $H(0,\omega)$ is the Fourier
amplitude of the wave at $x=0$, which we are free to choose. The quantity
$\eta_{0}=\sqrt{\mu_{0}/\epsilon_{0}}$ is the impedance of free space.
Although the solutions (5) have the appearance of plane waves, the operator
$\exp(\pm{\rm i}\hat{K}_{p}x)$ modifies the spectral content of the wave as
the observation point $x$ is changed, describing the reshaping of the pulse
during propagation. Eq. (5) shows that those spectra $H(0,\omega)$ that are
eigenfunctions of $\hat{K}_{p}$ with eigenvalue $\lambda$ have a plane wave
spatial dependence $\exp(\pm{\rm i}\lambda x)$ and retain the same spectrum
during propagation, despite the time variation of the material parameters.
Note also that Eq. (5) is similar to the aforementioned time evolution
operator in quantum mechanics, where a state $|\psi\rangle$ evolves in time as
$|\psi(t)\rangle=\exp(-{\rm i}\hat{H}t/\hbar)|\psi(0)\rangle$, where $\hat{H}$
is the Hamiltonian operator. Fig. 1 illustrates the basic idea of our
formalism, and gives the form of the operator $\hat{K}_{p}$ for a typical time
variation of the material parameters (4) (See SI for further details).

Figure 1: Scattering from a time varying dispersive medium. (a) Our theory
treats this as a generalization of a time independent problem, with operators
$\hat{r}$ and $\hat{t}$ replacing the usual reflection and transmission
coefficients. (b) Here we assume the Drude model permittivity (4) with an
asymmetric time variation of the plasma frequency, shown here in white. (c)
The operator $\hat{K}_{p}$, determines the spatial evolution of the spectral
content of the wave via Eq. (5). For the time dependence shown in panel (b),
the operator causes a spectral broadening (the smearing around
$\omega=\omega^{\prime}$, with this largest close to the plasma frequency) and
reversal of propagation direction (the line around $\omega=-\omega^{\prime}$).

## Fresnel coefficients for a dispersive, time–varying interface:

Consider a pulse incident from vacuum onto a dispersive time varying
half–space ($x>0$). Using the operator formalism described above, we calculate
the reflection $\hat{r}$ and transmission $\hat{t}$ operators for an incident
pulse. Just as for static materials, this calculation must be done separately
for transverse electric (TE) and transverse magnetic (TM) polarizations.

Assuming incidence in the $x$–$y$ plane with in–plane wave–vector
$\bm{k}_{\parallel}$, TE polarized waves have an electric field
$\tilde{\bm{E}}=E\bm{e}_{z}$ obeying the operator Helmholtz equation
$\partial_{x}^{2}E+\hat{K}^{2}_{s}E=0$ where $\hat{K}_{s}^{2}={\rm
i}\eta_{0}k_{0}\hat{\sigma}+k_{0}^{2}-k_{\parallel}^{2}$. Inside the TVM
($x>0$) the solution is given by Eq. (5), $E(x>0,\omega)=\exp({\rm
i}\hat{K}_{s}x)C_{s}(\omega)$. Meanwhile, on the entrance side the field is a
sum of plane waves for each frequency $E(x<0,\omega)=A_{s}(\omega)\exp({\rm
i}k_{x}x)+B_{s}(\omega)\exp(-{\rm i}k_{x}x)$ where
$k_{x}=[k_{0}^{2}-k_{\parallel}^{2}]^{1/2}$. The spatial boundary conditions
are the same as for static media, with both electric $E$ and magnetic
$\eta_{0}H_{y}={\rm i}k_{0}^{-1}\partial_{x}E$ fields continuous across $x=0$.
Substituting the forms of the fields in the respective regions leads to the
following reflection and transmission operators

$\displaystyle\hat{r}_{s}$
$\displaystyle=\left(1-\hat{Z}_{s}\right)\left(1+\hat{Z}_{s}\right)^{-1}$
$\displaystyle\hat{t}_{s}$ $\displaystyle=2\left(1+\hat{Z}_{s}\right)^{-1}$
(6)

where $\hat{Z}_{s}=k_{x}^{-1}\hat{K}_{s}$, $B_{s}=\hat{r}_{s}A_{s}$ and
$C_{s}=\hat{t}_{s}A_{s}$. Eqns. (6) are the TE Fresnel coefficients [28], with
an operator replacing the usual expression for the wave–vector in the
material.

The derivation is slightly different for TM polarization, revealing the
importance of operator ordering in these calculations. Taking the magnetic
field as $\bm{H}=H\bm{e}_{z}$, it obeys the operator Helmholtz equation
$\partial_{x}^{2}H+\hat{K}_{p}^{2}H=0$. The square of the wave–vector is as
above $\hat{K}_{p}^{2}={\rm
i}\eta_{0}\hat{\sigma}k_{0}+k_{0}^{2}-k_{\parallel}^{2}$, differing from the
TE expression due to the non–commuting nature of $\hat{\sigma}$ and $\omega$.
Applying the continuity of the magnetic $H$, and electric $E=-{\rm
i}[k_{0}+{\rm i}\eta_{0}\hat{\sigma}]^{-1}\,\,\partial_{x}\eta_{0}H$ fields at
the $x=0$ interface leads to the reflection and transmission operators

$\displaystyle\hat{r}_{p}$
$\displaystyle=\left(1-\hat{Z}_{p}\right)\left(1+\hat{Z}_{p}\right)^{-1}$
$\displaystyle\hat{t}_{p}$ $\displaystyle=2\left(1+\hat{Z}_{p}\right)^{-1}$
(7)

where $\hat{Z}_{p}=k_{x}^{-1}[1+{\rm
i}\eta_{0}\hat{\sigma}k_{0}^{-1}]^{-1}\hat{K}_{p}$. These are again operator
analogues of the TM Fresnel coefficients, although in this case the operator
ordering would not be obvious without applying the boundary conditions.
Importantly, at normal incidence the two wave–vector operators differ by a
similarity transformation: $\hat{K}_{p}=k_{0}^{-1}\hat{K}_{s}k_{0}$. The two
impedance operators are then simply related by
$\hat{Z}_{p}=\hat{K}_{s}^{-1}k_{0}=\hat{Z}_{s}^{-1}$, making the reflection
operators (6–7) differ by a minus sign $\hat{r}_{s}=-\hat{r}_{p}$ as expected
for the two polarizations at normal incidence 333We define the TM reflection
coefficients in terms of the _magnetic_ field amplitude, and the TE reflection
coefficient in terms of the _electric_ field amplitude, hence their difference
by a sign at normal incidence..

As discussed above, it is again interesting to examine the eigenvalues and
eigenvectors of the reflection and transmission operators. These reveal that
there are pulse spectra (‘eigenpulses’) that retain an identical spectrum
after interaction with the TVM (bar an overall multiplicative factor).
Alternatively, taking a singular value decomposition of the reflection and
transmission operators, we can find pulse spectra that are scaled by a set
level (the singular value), but have a different spectral content after
interaction with the material. We can see from expressions (6–7) that, in the
case of a single interface, the eigenpulses are the eigenfunctions of the
impedance operators $\hat{Z}_{s,p}$ and thus _both_ transmitted and reflected
spectra are unchanged. Fig. 2 shows a comparison between the reflection of a
Gaussian pulse and an eigenpulse from a TVM (plots show the incident field
just before the interface). While the Gaussian pulse is significantly
broadened and reshaped by the interaction with the TVM, the eigenpulse reduces
in frequency in tandem with the plasma frequency, retaining an identical
spectrum upon reflection. In this case (modulus of eigenvalue $|r|=1$), the
eigenpulse is also entirely reflected by the medium, as if it were a mirror.

Figure 2: Reflection of an eigenpulse from a time varying dispersive
half–space. (a) We compare the reflection of two different incident pulses; a
Gaussian pulse (upper curve), and an eigenpulse (lower curve) computed from
the reflection operator (7), with eigenvalue $|r|=1.00$ (the zero level is
displaced to aid visualization). (b) Time variation of
$N(t)=\omega_{p}^{2}(t)/\omega_{0}^{2}$ in Eq. (4). (c–d) Magnitude of
incident and reflected Fourier spectra computed via a numerical integration of
Maxwell’s equations (See Supplementary Information), for an incident (c)
eigenpulse, and (d) Gaussian pulse.

In addition, the scattering operators can also exhibit poles. For example in
Eqns. (6–7) these poles occur where ${\rm det}(1+\hat{Z}_{s,p})=0$. The
vectors in the null–space of $(1+\hat{Z}_{s,p})$ then represent non–time
harmonic modes that are—in this case—confined to the interface of the
material. In the Supplementary Information we find the surface plasmon–like
eigenpulses that are confined to the interface of a TVM,

## Time varying layer:

We can straightforwardly extend this approach to any multilayer and any simple
geometry (e.g. a spherical, cylindrical, or ellipsoidal object) that admits an
analytic solution to Maxwell’s equations in the static case. Broadly speaking,
the results for the scattering operators will have an identical form but with
an operator replacing the material parameters. To illustrate this in a
non–trivial case we calculate the reflection and transmission operators for a
slab of thickness $d$, which is relevant to the experiments reported in [9,
18, 10, 21, 19].

Figure 3: Comparison between the proposed theory and COMSOL multiphysics [22].
The figures demonstrate normalized Fourier transmitted spectra for a time-
modulated Drude slab excited by a 45 degree incident probe wave. (See
Supplementary Information for more details)

Assuming TM polarization, the magnetic field within the layer is taken to be
of the form

$H(0<x<d,\omega)={\rm e}^{{\rm i}\hat{K}_{p}x}C_{p}(\omega)+{\rm e}^{-{\rm
i}\hat{K}_{p}x}D_{p}(\omega)$ (8)

with the field in the external regions equal to $H(x<0)=\exp({\rm
i}k_{x}x)A_{p}(\omega)+\exp(-{\rm i}k_{x}x)B_{p}(\omega)$ and
$H(x>d)=\exp({\rm i}k_{x}(x-d))F_{p}(\omega)$. Imposing the same boundary
conditions described above we obtain the reflection and transmission
coefficients for the slab

$\hat{r}_{\rm slab}=\left[\hat{A}_{+}{\rm e}^{{\rm
i}\hat{K}_{p}d}\hat{A}_{-}-\hat{A}_{-}{\rm e}^{-{\rm
i}\hat{K}_{p}d}\hat{A}_{+}\right]\\\ \times\left[\hat{A}_{-}{\rm e}^{{\rm
i}\hat{K}_{p}d}\hat{A}_{-}-\hat{A}_{+}{\rm e}^{-{\rm
i}\hat{K}_{p}d}\hat{A}_{+}\right]^{-1}$ (9)

and

$\hat{t}_{\rm slab}=4\hat{Z}_{p}\left[\hat{A}_{+}{\rm e}^{-{\rm
i}\hat{K}_{p}d}\hat{A}_{+}-\hat{A}_{-}{\rm e}^{{\rm
i}\hat{K}_{p}d}\hat{A}_{-}\right]^{-1}$ (10)

where $\hat{A}_{\pm}=1\pm\hat{Z}_{p}$. Expressions (9–10) reduce to the
familiar reflection and transmission coefficients of a dielectric slab [28]
when the operators are replaced with their scalar counterparts. When $d=0$ the
reflection operator (9) is identically zero, and the transmission operator
(10) becomes the identity, as they should.

In Fig. 3 we give a comparison between results obtained using COMSOL
Multiphysics (see SI), and calculations made using the reflection and
transmission operators (9–10). We plot the normalized transmitted spectra as a
function of pulse delay time $\Delta t$. As shown in the lower panel of this
figure, there is excellent agreement between the finite element calculation
and our operator approach. Additional comparisons to an adiabatic multiple-
timescale approach used to model past experiments [18, 21] are also available
in the SI. Importantly, these tests demonstrate the advantage of this method
for the efficient modelling of structures that feature extremely subwavelength
layers, circumventing the need for expensive numerical calculations.

## Summary and Conclusions:

We have developed a compact theoretical approach for treating the problem of
scattering from dispersive TVM. We have shown that our analytic expressions
match full wave numerical simulations well. Although the theory is formally
similar to the case of static materials, the TVM parameters are given in terms
of operators that depend on both the frequency and frequency derivatives,
which must be carefully ordered. The advantage of our theory is that it is
semi–analytical, allowing us to give explicit operator expressions for
scattering coefficients from the TVM, and thus determine conditions for e.g.
incoming modes that are bound, not reflected, or completely reflected by the
material. We have numerically constructed these operators and found the
‘eigenpulses’ of a time–modulated Drude half–space, numerically verifying that
there are input pulse spectra that e.g. reflect as if the TVM was a static
mirror. This approach may be readily extended to other areas of wave physics
such as pressure acoustics [30] and elasticity and may be of interest to those
working on TVM as well as multiple scattering, where our reflectionless
eigenpulses are analogous to the concept of open scattering channels in
disordered media (see e.g. [31]).

###### Acknowledgements.

## Acknowledgements:

SARH acknowledges the Royal Society and TATA for financial support through
grant URF\R\211033, and thanks Riccardo Sapienza, Euan Hendry, James Capers,
and Dean Patient for useful conversations. EG acknowledges funding from the
Simons Foundation through a Junior Fellowship of the Simons Society of Fellows
(855344,EG).

## References

  * Segev and Bandres [2021] M. Segev and M. A. Bandres, Topological photonics: Where do we go from here?, Nanophotonics 10, 425 (2021).
  * El-Ganainy  _et al._ [2018] R. El-Ganainy, K. G. Makris, M. Khajavikhan, Z. H. Musslimani, S. Rotter, and D. N. Christodoulides, Non-hermitian physics and pt symmetry, Nature Phys. 14, 11 (2018).
  * Lee  _et al._ [202] K.-T. Lee, C. Ji, H. Iizuka, and D. Banerjee, Optical cloaking and invisibility: From fiction toward a technological reality, Journal of Applied Physics 129, 231101 (202).
  * Galiffi  _et al._ [2022a] E. Galiffi, R. Tirole, S. Yin, H. Li, S. Vezzoli, P. A. Huidobro, M. G. Silveirinha, R. Sapienza, A. Alú, and J. B. Pendry, Photonics of time-varying media, Adv. Photon. 4, 014002 (2022a).
  * Galiffi  _et al._ [2022b] E. Galiffi, P. Huidobro, and J. Pendry, An archimedes’ screw for light, Nat. Commun. 13, 2523 (2022b).
  * Landau and Lifshitz [2005] L. D. Landau and E. M. Lifshitz,  _Statistical Physics: Part 1_ (Butterworth-Heinemann, 2005).
  * Rozanov [2000] K. N. Rozanov, Ultimate thickness to bandwidth ratio of radar, IEEE Trans. Ant. Prop. 48, 1230 (2000).
  * Solis and Engheta [2021] D. M. Solis and N. Engheta, Functional analysis of the polarization response in linear time-varying media: A generalization of the kramers-kronig relations, Phys. Rev. B 103, 144303 (2021).
  * Zhou  _et al._ [2020] Y. Zhou, M. Z. Alam, M. Karimi, J. Upham, O. Reshef, C. Liu, A. E. Willner, and R. W. Boyd, Broadband frequency translation through time refraction in an epsilon-near-zero material, Nat. Commun. 11, 2180 (2020).
  * Bohn  _et al._ [2021a] J. Bohn, T. S. Luk, S. A. R. Horsley, and E. Hendry, Spatiotemporal refraction of light in an epsilon-near-zero indium tin oxide layer: frequency shifting effects arising from interfaces, Optica 8, 1532 (2021a).
  * Morgenthaler [1958] F. R. Morgenthaler, Velocity modulation of electromagnetic waves, IRE Transactions on microwave theory and techniques 6, 167 (1958).
  * Mendonça and Shukla [2002] J. T. Mendonça and P. K. Shukla, Time refraction and time reflection: Two basic concepts, Physica Scripta 65, 160 (2002).
  * Bacot  _et al._ [2016] V. Bacot, M. Labousse, A. Eddi, M. Fink, and E. Fort, Time reversal and holography with spacetime transformations, Nature Physics 12, 972 (2016).
  * Moussa  _et al._ [2022] H. Moussa, G. Xu, S. Yin, E. Galiffi, Y. Radi, and A. Alù, Observation of temporal reflections and broadband frequency translations at photonic time-interfaces, arXiv preprint arXiv:2208.07236 (2022).
  * Pacheco-Peña and Engheta [2020] V. Pacheco-Peña and N. Engheta, Temporal aiming, Light: Science & Applications 9, 129 (2020).
  * Jaffray  _et al._ [2022] W. Jaffray, S. Saha, V. M. Shalaev, A. Boltasseva, and M. Ferrera, Transparent conducting oxides: from all-dielectric plasmonics to a new paradigm in integrated photonics, Adv. Opt. Photon. 14, 148 (2022).
  * Alam  _et al._ [2016] M. Z. Alam, I. De Leon, and R. W. Boyd, Large optical nonlinearity of indium tin oxide in its epsilon-near-zero region, Science 352, 795 (2016).
  * Bohn  _et al._ [2021b] T. S. Bohn, J. ans Luk, C. Tollerton, S. W. Hutchings, I. Brener, S. A. R. Horsley, W. L. Barnes, and E. Hendry, All-optical switching of an epsilon-near-zero plasmon resonance in indium tin oxide, Nat. Commun. 12, 1017 (2021b).
  * Tirole  _et al._ [2022a] R. Tirole, S. Vezzoli, E. Galiffi, I. Robertson, D. Maurice, B. Tilmann, S. A. Maier, J. B. Pendry, and R. Sapienza, Double-slit time diffraction at optical frequencies, arXiv preprint arXiv:2206.04362 (2022a).
  * Liu  _et al._ [2021] C. Liu, M. Z. Alam, K. Pang, K. Manukyan, O. Reshef, Y. Zhou, S. Choudhary, J. Patrow, A. Pennathurs, H. Song, Z. Zhao, R. Zhang, F. Alishahi, A. Fallahpour, Y. Cao, A. Almaiman, J. M. Dawlaty, M. Tur, R. W. Boyd, and A. E. Willner, Photon acceleration using a time-varying epsilon-near-zero metasurface, ACS Photonics 8, 716 (2021).
  * Tirole  _et al._ [2022b] R. Tirole, E. Galiffi, J. Dranczewski, T. Attavar, B. Tilmann, Y.-T. Wang, P. A. Huidobro, A. Alú, J. B. Pendry, S. A. Maier, S. Vezzoli, and R. Sapienza, Saturable time-varying mirror based on an enz material, arXiv:2202.05937 (2022b).
  * Multiphysics [1998] C. Multiphysics, Introduction to comsol multiphysics®, COMSOL Multiphysics, Burlington, MA, accessed Feb 9, 2018 (1998).
  * Ptitcyn  _et al._ [2021] G. Ptitcyn, A. G. Lamprianidis, T. Karamanos, V. S. Asadchy, R. Alaee, M. Müller, M. Albooyeh, M. S. Mirmoosa, S. Fan, S. A. Tretyakov, and C. Rockstuhl, Scattering of light by spheres made from a time-modulated and dispersive material, 2021 Fifteenth International Congress on Artificial Materials for Novel Wave Phenomena (Metamaterials) , 347 (2021).
  * Note [1] The opposite (normal) operator ordering (where the $\partial_{\omega}$ sit on the _right_ of the factors of frequency in $\hat{\sigma}$) is required if we specify the conductivity in Eq. (1) as $\sigma(t^{\prime},t-t^{\prime})$.
  * Shewell [1959] J. R. Shewell, On the formation of quantum-mechanical operators, Am. J. Phys. 27, 16 (1959).
  * Sakuri and Napolitano [2018] J. J. Sakuri and J. Napolitano,  _Modern Quantum Mechanics_ (Cambridge University Press, 2018).
  * Note [2] It is important to determine which sign of the matrix square root $\pm\hat{K}_{p}$ describes right–going waves. We take right–going waves to be those where the eigenvalues have a positive imaginary part.
  * Landau  _et al._ [2004] L. D. Landau, E. M. Lifshitz, and L. P. Pitaevskii,  _Electrodynamics of Continuous Media_ (Butterworth-Heinemann, 2004).
  * Note [3] We define the TM reflection coefficients in terms of the _magnetic_ field amplitude, and the TE reflection coefficient in terms of the _electric_ field amplitude, hence their difference by a sign at normal incidence.
  * Cho  _et al._ [2020] C. Cho, X. Wen, N. Park, and J. Li, Digitally virtualized atoms for acoustic metamaterials, Nat. Commun. 11, 251 (2020).
  * Choi  _et al._ [2011] W. Choi, A. P. Mosk, Q.-H. Park, and W. Choi, Transmission eigenchannels in a disordered medium, Phys. Rev. B 83, 134207 (2011).