In principle any decomposition Eq. (\[decompos\]) of $H$ into real skew-symmetric parts will do. Largely guided by previous work [@DeRaedt87; @DeRaedt94; @Michielsen98; @DeRaedt00], we have adopted a decomposition that is efficient, flexible, sufficiently accurate and easy to program.

One dimension {#subsec41}
-------------

We consider a 1D system along the $x$-direction. Accordingly, Maxwell’s equations contain no partial derivatives with respect to $y$ or $z$ and ${\varepsilon}$ and $\mu$ do not depend on $y$ or $z$. Under these conditions, the TDME reduce to two independent sets of first-order differential equations [@BornWolf]. The solutions to these sets are known as the transverse electric (TE) mode and the transverse magnetic (TM) mode [@BornWolf]. As the equations of the TE- and TM-mode only differ by a sign we can restrict our considerations to the TM-mode and obtain the result for the TE-mode by reversing the time.

From Eq. (\[eqn:mtxeqn\]) it follows that the magnetic field $H_y(x,t)=X_y(x,t)/\sqrt{\mu(x)}$ and the electric field $E_z(x,t)=Y_z(x,t)/\sqrt{{\varepsilon}(x)}$ of the TM-mode are solutions of $$\begin{aligned}
{\frac{\partial}{\partial t}}X_y(x,t) & = &
\frac{1}{\sqrt{\mu(x)}}{\frac{\partial}{\partial x}}\left(\frac{Y_z(x,t)}{\sqrt{{\varepsilon}(x)}}\right),
\label{eqn:TM1Dx} \\
{\frac{\partial}{\partial t}}Y_z(x,t) & = &
\frac{1}{\sqrt{{\varepsilon}(x)}}{\frac{\partial}{\partial x}}\left(\frac{X_y(x,t)}{\sqrt{\mu(x)}}\right).
\label{eqn:TM1Dy}\end{aligned}$$ Note that in 1D, the divergence of $H_y(x,t)$ and $E_z(x,t)$ is zero. Hence Eqs.(\[max3\]) and (\[max4\]) are automatically satisfied.

Using 