Dyakonov Surface Waves: Anisotropy-Enabling Confinement on the Edge

Abstract

The title “Dyakonov surface waves: anisotropy enabling confinement on the edge” plainly sets the scope for this chapter. The focus here is on the formation of bounded waves at the interface of two distinct media, at least one of them exhibiting optical anisotropy, which are coined as Dyakonov surface waves (DSWs) in recognition to the physicist who reported their existence for the first time. First, the general aspects of the topic are discussed. It also treats the characterization of bounded waves in isotropic-uniaxial multilayered structures, allowing not only the derivation of the dispersion equation of DSWs but also that of surface plasmons polaritons (SPPs), for instance. Furthermore, the interaction of such surfaces waves, with the possibility of including guided waves in a given planar layer and external sources mimicking experimental setups, can be accounted for by using the transfer matrix formalism introduced here. Finally, special attention is devoted to hyperbolic media with indefinite anisotropy-enabling hybridized scenarios integrating the prototypical DSWs and SPPs.

1. Introduction

The planar interface of two dissimilar materials plays a relevant role in many optical phenomena. In recent years, particularly, evanescent waves have been used in newly developing technologies such as near-field spectroscopy. The electromagnetic surface wave, which is intimately tied to the interface, travels in a direction parallel to the interface but, on either side of the interface, its amplitude is imperceptible after a certain distance from the interface. The notion of an electromagnetic surface wave made a significant appearance in 1907 when Zen-neck [1] authored a theoretical paper exploring the possibility of a wave guided by the interface of the atmosphere and either Earth or a large body of water. His focus was on radio waves, a region of the electromagnetic spectrum far from the optical regime in which we are particularly interested, but the principles involved are the same, owing to the scale invariance of the Maxwell postulates.

Yet, nearly a century later, a unique type of wave, the surface plasmon polariton (SPP) wave, that dominates the nanotechnology scene, at least at optical frequencies, resulted in wonderful developments with the creation of extremely sensitive bio/chemical sensors, and improvements in this mature technology continue to this day [2]. Even in this highly developed application, the two partnering materials which meet at the interface may be simple: one is a typical metal, a plasmonic material at optical frequencies, and the other is a homogeneous, isotropic, dielectric material. While the interface of a plasmonic material and a polarizable material supports SPPs, a variety of other types of surface waves can be supported by the interface of two polarizable materials. Since polarizable materials such as dielectric materials are less dissipative, in general, than plasmonic materials such as metals, the advantage of these materials for long-range propagation of surface waves is apparent.

The interface of two homogeneous dielectric materials of which at least one is anisotropic may support surface-wave propagation of another type, even though the real parts of all components of the permittivity dyadics of both materials are positive. Interest in surface waves guided by the interface of two dielectric materials began to take after Dyakonov in 1988 [3] explored the propagation of a surface wave guided by the interface of a uniaxial dielectric material and an isotropic dielectric material. The Dyakonov surface wave (DSWs) is the focus of this chapter.

In this chapter, we perform a thorough analysis of DSWs taking place in semi-infinite anisotropic media. Basic concepts related to the propagation of electromagnetic waves in homogeneous media will be introduced, including isotropic and anisotropic materials. Birefringent metal-dielectric (MD) lattices will be also considered as a contribution of meta-materials in the development of DSWs [4]. Special emphasis will be put when the effective-medium approaches induce unsatisfactory results, which occur in most experimental configurations. Practical cases will be analyzed including dissipative effects due to Ohmic losses of the metal.

2. Wave propagation in bulk media and interfaces

In this section we introduce the basic concepts related to the propagation of electromagnetic waves in homogeneous media, including isotropic and anisotropic materials. We describe in detail complex multilayered structures. For that purpose, we introduce a transfer matrix formulation that applies to isotropic and uniaxial media simultaneously. Finally, we discuss the conditions that give rise to surface waves at the interface of two isotropic media; the case of dealing with anisotropic media is considered in Section 5. Moreover, we obtain the dispersion equation for SPPs, which appears at the interface between a dielectric and a metal.

2.1. Wave propagation in isotropic media

In this section, we consider the propagation of electromagnetic waves in linear, homogeneous, and isotropic dielectrics. Under these conditions, the relative permittivity ϵ relating E and D is a scalar constant. Considering that the medium is free of electric charges and currents, and taking into account the medium equation $\mathbf{B}={\mu }_0\mathbf{H}$, Maxwell’s equations can be written as

$$\mathbf{\nabla }\mathbf{\times }\mathbf{H}=\mathrm{ϵ}{\mathrm{ϵ}}_0\frac{\partial \mathbf{E}}{\mathrm{\partial t}},$$E1a

$$\mathbf{\nabla }\mathbf{\times }\mathbf{E}=\mathbf{-}{\mathrm{\mu }}_0\frac{\partial \mathbf{H}}{\mathrm{\partial t}}\mathbf{,}$$E1b

$$\mathbf{\nabla }·\mathbf{E}=0,$$E1c

$$\mathbf{\nabla }·\mathbf{H}=0.$$E1d

Now, each of the scalar components of $\mathbf{E}$ and $\mathbf{H}$ satisfies the wave equation ${\nabla }^{2}u-c^{-2}{\partial }_t^{2}u=0$, where $u$ represents any one of the six scalar components of the electromagnetic field and $c$ is the speed of the waves in the medium.

When the electromagnetic wave is plane and monochromatic, all components of the electric and magnetic fields are harmonic functions in time and space at the time frequency $\omega$ and spatial frequency $\mathbf{k}=\left(k_x,k_y,k_z\right),$ respectively. Particularly they may be set as

$$\mathbf{E}\left(\mathbf{r},\mathrm{t}\right)={\mathbf{E}}_0\text{exp}\left(\mathbf{ik}·\mathbf{r}\right)\text{exp}\left(-\mathrm{i\omega t}\right),$$E2a

$$\mathbf{H}\left(\mathbf{r},\mathrm{t}\right)={\mathbf{H}}_0\text{exp}\left(\mathbf{ik}·\mathbf{r}\right)\text{exp}\left(-\mathrm{i\omega t}\right),$$E2b

where ${\mathbf{E}}_0\left(\mathbf{r}\right)$ and ${\mathbf{H}}_0\left(\mathbf{r}\right)$ are the complex amplitudes of the electric and magnetic fields. If we substitute the vector wave fields of Eqs. (2) in Maxwell’s equations (1), with the help of the transformations $\nabla \to i\mathbf{k}$ and ${\partial }_t\to -\mathit{i\omega }$, we might attain a simplified expression of $\mathbf{k}\times \left(\mathbf{k}\times {\mathbf{E}}_0\right)$, enabling to obtain the following wave equation: $\mathbf{M}·{\mathbf{E}}_0=0$, where $\mathbf{M}\equiv \mathbf{k}\mathbf{\otimes }\mathbf{k}-k^2\mathbf{I}+k_0^2ϵ\mathbf{I}$. Here $\mathbf{I}$ is the 3 × 3 identity matrix, and $k$ is the modulus of the wavevector $\mathbf{k}$ and $k_0=\omega /c_0$. In order to obtain the dispersion equation, we look for nontrivial solutions of the electric field ${\mathbf{E}}_0$ by imposing that $\text{det}\left(M\right)=0$. Its solution leads to $k^2={\omega }^2/c^2$. The electric field amplitude ${\mathbf{E}}_0$ can be written as a linear combination of the following vectors

$${\widehat{e}}_1=\left(0,k_z,-k_y\right),$$E3a

$${\widehat{e}}_2=\left(k_y^2+k_z^2,-k_y{k}_x,-k_z{k}_x\right).$$E3b

We point out that ${\widehat{e}}_1$ is associated with ${\mathrm{TE}}^{\mathrm{x}}$ modes and that ${\widehat{e}}_1.{\widehat{e}}_2=0$. Although not demonstrated here, ${\widehat{e}}_2$ is related to ${\mathrm{TM}}^{\mathrm{x}}$-polarized plane waves. Note also that Eqs. (1c) and (1d) lead to the following orthogonality relations-hips: $\mathbf{k}·{\widehat{e}}_1=\mathbf{k}·{\widehat{e}}_2=0$. As a result, the vectors $\left\{{\widehat{e}}_1,{\widehat{e}}_2,\mathbf{k}\right\}$ form an orthogonal trihedron.

2.2. Wave propagation in uniaxial media

Uniaxial crystals are media with certain symmetries that make them have two equal principal refractive indices: $n_x=n_y\equiv n_o$ (ordinary index) and $n_z\equiv n_e$ (extraordinary index). The crystal is to be a positive uniaxial if $n_e>n_o$ and negative uniaxial if $n_e<n_o$. The $z$ axis of a uniaxial crystal is called the optic axis. We call it the ordinary polarization direction if the wave has an eigenindex of refraction $n_o$ and extraordinary, if the wave has an eigenindex of refraction $n_e$. In Table 1 we show the values of $n_o$ and $n_e$ for some natural birefringent materials (uniaxial crystals).

Birefringent material	$n_o$	$n_e$	${\Delta n}=n_e-n_o$
Crystal quartz	1.547	1.556	0.009
MgF	1.3786	1.3904	0.0118
YVO4	1.9929	2.2154	0.2225
Rutile (TiO2)	2.65	2.95	0.3
E7 liquid crystal	1.520	1.725	0.205
Calomel (Hg2Cl2)	1.96	2.62	0.68

Table 1.

Birefringence $\mathit{\Delta n}$ of some natural materials [5].

Let us now consider the propagation of wave planes in uniaxial media. To obtain the eigenvalues associated with plane-wave propagation, we proceed in a similar way as in Section 2.1 for isotropic media, but taking into account that now the relative permittivity $\mathbf{ϵ}={ϵ}_x\left(\mathbf{x}\mathbf{\otimes }\mathbf{x}+\mathbf{y}\mathbf{\otimes }\mathbf{y}\right)+{ϵ}_z\left(\mathbf{z}\mathbf{\otimes }\mathbf{z}\right)$ is a tensor. Therefore, we solve $\text{det}\left(\mathbf{M}\right)=0$, where the matrix $\mathbf{M}\equiv \mathbf{k}\mathbf{\otimes }\mathbf{k}-k^2\mathbf{I}+k_0^2\mathbf{ϵ}$. Then, after solving the abovementioned determinant, we obtain two solutions. The first of them is $k_x^2+k_y^2+k_z^2={ϵ}_x{k}_0^2$ which corresponds to the dispersion equation of the ordinary plane waves. The electric field for this kind of plane waves is proportional to the vector ${\widehat{e}}_o=\left(-k_y,k_x,0\right)$. As a consequence, ordinary plane waves are ${\mathrm{TE}}^{\mathrm{z}}$-polarized waves. The second solution gives us the dispersion equation of the extraordinary plane waves:

$$\frac{k_x^2+k_y^2}{{ϵ}_z}+\frac{k_z^2}{{ϵ}_x}=k_0^2.$$E4

In this case, the electric field is proportional to the vector ${\widehat{e}}_e=\left(k_x{k}_z,k_y{k}_z,k_z^2-{ϵ}_x{k}_0^2\right)$.

2.3. Matrix formulation for multilayered media

In this section we look at the case of multilayered media composed of different nonmagnetic materials which are separated by planar parallel interfaces, displaced on $x=x_i$ as shown in Figure 1. In particular, we deal with either uniaxial or isotropic materials, including metals. We make a detailed description of the electromagnetic fields inside a given medium “i”, which lies in $x_{i-1}<x<x_i$. Our objective is the analysis of the appropriate conditions for the propagation of surface waves at the abovementioned interfaces.

Figure 1.

Schematic arrangement of the multilayered media. The amplitudes $A_i$ and $A{\prime }_i$ ($B_i$ and $B{\prime }_i$) correspond to waves propagating along the positive (negative) $x$-axis. These amplitudes characterize a given state of polarization: For uniaxial media we may deal with ordinary (for example, $A_{\mathit{oi}}$) and extraordinary ($A_{\mathit{ei}}$) waves. For isotropic materials we have TE-polarized ($A_{\mathit{TEi}}$) and TM-polarized ($A_{\mathit{TMi}}$) waves.

2.3.1. Electromagnetic fields in a uniaxial elementary layer

Let us first consider a multilayered media of uniaxial materials. For simplicity, let us take into account only relative permittivities of the form ${\mathbf{ϵ}}_i={ϵ}_{\mathit{xi}}\left(\mathbf{x}\mathbf{\otimes }\mathbf{x}+\mathbf{y}\mathbf{\otimes }\mathbf{y}\right)+{ϵ}_{\mathit{zi}}\left(\mathbf{z}\mathbf{\otimes }\mathbf{z}\right)$ for a given medium “$i$.” As we demonstrated in Section 2.2, the dispersion equation for extraordinary waves propagating in bulk uniaxial media is given by Eq. (4). Due to the boundary conditions at the interfaces between different media, the components of the wave vector $k_y$ and $k_z$ are conserved, but not its projection upon the $x$-axis. More specifically, if we rename $k_{\mathit{xi}}\equiv k_{\mathit{oi}}$ for ordinary waves, we can write

$$k_{\mathit{oi}}=\sqrt{{ϵ}_{\mathit{xi}}{k}_0^2-\left(k_y^2+k_z^2\right)},$$E5

and for extraordinary waves (here $k_{\mathit{xi}}\equiv k_{\mathit{ei}}$)

$$k_{\mathit{ei}}=\sqrt{{ϵ}_{\mathit{zi}}{k}_0^2-\left(k_y^2+\frac{{ϵ}_{\mathit{zi}}{k}_z^2}{{ϵ}_{\mathit{xi}}}\right)},$$E6

taken from Eq. (4).

The total electric field of the elementary layer “$i$” is set as

$${\mathbf{E}}_{\mathit{tot}}^{\left(i\right)}={\mathbf{E}}^{\left(i\right)}\left(x\right)\text{exp}\left(i{k}_{y}y+i{k}_{z}z-\mathit{i\omega t}\right).$$E7

The part of the electric field that varies along with the spatial coordinate $x$ can be written as

$$\begin{array}{ll}{\mathbf{E}}^{\left(i\right)}\left(x\right)=& A_{\mathit{oi}}{\widehat{a}}_{\mathit{oi}}\text{exp}\left[i{k}_{\mathit{oi}}\left(x-x_i\right)\right]+B_{\mathit{oi}}{\widehat{b}}_{\mathit{oi}}\text{exp}\left[-i{k}_{\mathit{oi}}\left(x-x_i\right)\right]\hfill \\ & +A_{\mathit{ei}}{\widehat{a}}_{\mathit{ei}}\text{exp}\left[i{k}_{\mathit{ei}}\left(x-x_i\right)\right]+B_{\mathit{ei}}{\widehat{b}}_{\mathit{ei}}\text{exp}\left[-i{k}_{\mathit{ei}}\left(x-x_i\right)\right],\hfill \end{array}$$E8

where the amplitude $A_{\mathit{oi}}$ (and $A_{\mathit{ei}}$) corresponds to propagating ordinary (and extraordinary) waves, and $B_{\mathit{oi}}$ (and $B_{\mathit{ei}}$) is related with counter-propagating ordinary (and extraordinary) waves. Note that all these amplitudes have zero dephase at $x=x_i$. Finally, the vectors ${\widehat{a}}_{\mathit{oi}}$ and ${\widehat{a}}_{\mathit{ei}}$ are rewritten as:

$${\widehat{a}}_{\mathit{oi}}=\left(-k_y,k_{\mathit{oi}},0\right),$$E9a

$${\widehat{a}}_{\mathit{ei}}=\left(k_{\mathit{ei}}{k}_z,k_y{k}_z,k_z^2-{ϵ}_{\mathit{xi}}{k}_0^2\right).$$E9b

In the case of counter-propagating waves, we take into consideration that $k_{\mathit{xi}}=-k_{\mathit{oi}}$ for ordinary waves and $k_{\mathit{xi}}=-k_{\mathit{ei}}$ for extraordinary waves. This fact leads us to introduce the vector fields

$${\widehat{b}}_{\mathit{oi}}=\left(-k_y,-k_{\mathit{oi}},0\right),$$E10a

$${\widehat{b}}_{\mathit{ei}}=\left(-k_{\mathit{ei}}{k}_z,k_y{k}_z,k_z^2-{ϵ}_{\mathit{xi}}{k}_0^2\right).$$E10b

The vector fields ${\widehat{a}}_{\mathit{oi}}$ and ${\widehat{b}}_{\mathit{oi}}$ can be given in units of $k_0$, whereas ${\widehat{a}}_{\mathit{ei}}$ and ${\widehat{b}}_{\mathit{ei}}$ can be expressed in units of $k_0^2$.

For convenience, the field function ${\mathbf{E}}^{\left(i\right)}$ can be set in terms of the wave amplitudes $A{\prime }_{\mathit{oi}}$, $A{\prime }_{\mathit{ei}}$ and $B{\prime }_{\mathit{oi}}$ and $B{\prime }_{\mathit{ei}}$ with zero dephase at $x=x_{i-1}$, namely:

$$\begin{array}{ll}{\mathbf{E}}^{\left(i\right)}\left(x\right)& =A{\prime }_{\mathit{oi}}{\widehat{a}}_{\mathit{oi}}\text{exp}\left[i{k}_{\mathit{oi}}\left(x-x_{i-1}\right)\right]+B{\prime }_{\mathit{oi}}{\widehat{b}}_{\mathit{oi}+1}\text{exp}\left[-i{k}_{\mathit{oi}}\left(x-x_{i-1}\right)\right]\hfill \\ & +A{\prime }_{\mathit{ei}}{\widehat{a}}_{\mathit{ei}}\text{exp}\left[i{k}_{\mathit{ei}}\left(x-x_{i-1}\right)\right]+B{\prime }_{\mathit{ei}}{\widehat{b}}_{\mathit{ei}}\text{exp}\left[-i{k}_{\mathit{ei}}\left(x-x_{i-1}\right)\right].\hfill \end{array}$$E11

The complete set of amplitudes is $A_{\mathit{qi}}$, $A{\prime }_{\mathit{qi}}$ and $B_{\mathit{qi}}$ and $B{\prime }_{\mathit{qi}}$ where $q=\left\{o,e\right\}$ satisfies the following relationships:

$$A{\prime }_{\mathit{qi}}=A_{\mathit{qi}}\text{exp}\left[-i{k}_{\mathit{qi}}\left(x_i-x_{i-1}\right)\right],$$E12a

$$B{\prime }_{\mathit{qi}}=B_{\mathit{qi}}\text{exp}\left[i{k}_{\mathit{qi}}\left(x_i-x_{i-1}\right)\right].$$E12b

For completeness we calculate the magnetic field in a given elementary layer “$i$.” By using the Maxwell’s equation ${\mathbf{H}}_{\mathit{tot}}^{\left(i\right)}={\left(i\omega {\mu }_0\right)}^{-1}\mathbf{\nabla }\mathbf{\times }{\mathbf{E}}_{\mathit{tot}}^{\left(i\right)}$, and considering that the magnetic field may be written as ${\mathbf{H}}_{\mathit{tot}}^{\left(i\right)}={\mathbf{H}}^{\left(i\right)}\left(x\right)\text{exp}\left(i{k}_{y}y+i{k}_{z}z-\mathit{i\omega t}\right)$, we finally obtain the following expression for the variation of the field along the $x$ direction, namely:

$$\begin{array}{ll}\omega {\mu }_0{\mathbf{H}}^{\left(i\right)}\left(x\right)& =A_{\mathit{oi}}{\widehat{c}}_{\mathit{oi}}\text{exp}\left[i{k}_{\mathit{oi}}\left(x-x_i\right)\right]+B_{\mathit{oi}}{\widehat{d}}_{\mathit{oi}}\text{exp}\left[-i{k}_{\mathit{oi}}\left(x-x_i\right)\right]\hfill \\ & +A_{\mathit{ei}}{\widehat{c}}_{\mathit{ei}}\text{exp}\left[i{k}_{\mathit{ei}}\left(x-x_i\right)\right]+B_{\mathit{ei}}{\widehat{d}}_{\mathit{ei}}\text{exp}\left[-i{k}_{\mathit{ei}}\left(x-x_i\right)\right].\hfill \end{array}$$E13

In the previous equations, the new vector fields are set as

$${\widehat{c}}_{\mathit{oi}}=-\left(k_{\mathit{oi}}{k}_z,k_y{k}_z,k_z^2-{ϵ}_{\mathit{xi}}{k}_0^2\right),$$E14a

$${\widehat{c}}_{\mathit{ei}}={ϵ}_{\mathit{xi}}{k}_0^2\left(-k_y,k_{\mathit{ei}},0\right),$$E14b

$${\widehat{d}}_{\mathit{oi}}=-\left(-k_{\mathit{oi}}{k}_z,k_y{k}_z,k_z^2-{ϵ}_{\mathit{xi}}{k}_0^2\right),$$E14c

$${\widehat{d}}_{\mathit{ei}}=-{ϵ}_{\mathit{xi}}{k}_0^2\left(k_y,k_{\mathit{ei}},0\right).$$E14d

We point out that the vector fields ${\widehat{c}}_{\mathit{oi}}$ and ${\widehat{d}}_{\mathit{oi}}$ can be given in units of $k_0^2$, whereas ${\widehat{c}}_{\mathit{ei}}$ and ${\widehat{d}}_{\mathit{ei}}$ can be expressed in units of $k_0^3$.

2.3.2. Boundary conditions for anisotropic layered media

Once we have a complete description of the wave fields in every elementary layer of our metamaterial, we have to impose some boundary conditions at the interfaces $x=x_i$. The components of the electric field ${\mathbf{E}}_{\mathit{tot}}^{\left(i\right)}$ and the magnetic field ${\mathbf{H}}_{\mathit{tot}}^{\left(i\right)}$ lying on the planes $x=x_i$ must be continuous. Particularly these boundary conditions may be expressed in a matrix form as

$${\mathbf{D}}_i{\mathbf{v}}_i={\mathbf{D}}_{i+1}\mathbf{v}{\prime }_{i+1}.$$E15

We defined the following matrix, explicitly given as

$${\mathbf{D}}_i=\left[\begin{array}{cccc}{k}_{\mathit{oi}}& -k_{\mathit{oi}}& k_y{k}_z& k_y{k}_z\\ 0& 0& k_z^2-{ϵ}_{\mathit{xi}}{k}_0^2& k_z^2-{ϵ}_{\mathit{xi}}{k}_0^2\\ -k_y{k}_z& -k_y{k}_z& {ϵ}_{\mathit{xi}}{k}_0^2{k}_{\mathit{ei}}& -{ϵ}_{\mathit{xi}}{k}_0^2{k}_{\mathit{ei}}\\ -k_z^2+{ϵ}_{\mathit{xi}}{k}_0^2& -k_z^2+{ϵ}_{\mathit{xi}}{k}_0^2& 0& 0\end{array}\right].$$E16

On the other hand, we introduced the amplitude column vectors

$${\mathbf{v}}_i=\left[\begin{array}{c}{A}_{\mathit{oi}}\\ B_{\mathit{oi}}\\ A_{\mathit{ei}}\\ B_{\mathit{ei}}\end{array}\right],\mathbf{v}{\prime }_i=\left[\begin{array}{c}A{\prime }_{\mathit{oi}}\\ B{\prime }_{\mathit{oi}}\\ A{\prime }_{\mathit{ei}}\\ B{\prime }_{\mathit{ei}}\end{array}\right].$$E17

The matrix formulation can also be used to relate the amplitude vector ${\mathbf{v}}_i$ with zero-phase shift at $x=x_i$ with the amplitude vector $\mathbf{v}{\prime }_i$ exhibiting zero dephase at $x=x_{i-1}$, previously established in Eq. (12a), (12b). For that purpose, we introduce the propagation matrix ${\mathbf{P}}_i$, which takes into account the amplitude dephasing between the boundaries of each layer. Explicitly we may write

$$\mathbf{v}{\prime }_i={\mathbf{P}}_i·{\mathbf{v}}_i,$$E18

being

$${\mathbf{P}}_i=\left[\begin{array}{cccc}{e}^{-i{k}_{\mathit{oi}}{w}_i}& 0& 0& 0\\ 0& e^{i{k}_{\mathit{oi}}{w}_i}& 0& 0\\ 0& 0& e^{-i{k}_{\mathit{ei}}{w}_i}& 0\\ 0& 0& 0& e^{i{k}_{\mathit{ei}}{w}_i}\end{array}\right],$$E19

where $w_i=x_i-x_{i-1}$ denotes the width of the elementary layer “$i.$” We note that equivalent matrix formulations for anisotropic multilayered media can be found elsewhere [6, 7].

2.3.3. Electromagnetic fields in layered isotropic media

At this point, once we have described the electromagnetic fields in uniaxial media, let us study a multilayered media composed of isotropic materials. Considering an isotropic medium of relative permittivity ${ϵ}_i$, once again, the projection of the wave vector along the positive $x$ direction is set as

$$k_{\mathit{TEi}}=k_{\mathit{TMi}}=\sqrt{{ϵ}_i{k}_0^2-\left(k_y^2+k_z^2\right)}.$$E20

Formally, $k_{\mathit{TMi}}$ applies to TM-polarized waves and $k_{\mathit{TEi}}$ corresponds to TE-polarized waves. The total electric field of the elementary layer “$i$” can be set, again, as given in Eq. (7). The part of the electric field that varies along with the spatial coordinate $x$ can be written using the amplitude $A_{\mathit{TEi}}$ (and $A_{\mathit{TMi}}$), which corresponds to propagating TE (and TM) waves and $B_{\mathit{TEi}}$ (and $B_{\mathit{TMi}}$) that is related with counter-propagating TE (and TM) waves. This finally reads as

$$\begin{array}{ll}{\mathbf{E}}^{\left(i\right)}\left(x\right)& =A_{\mathit{TEi}} {\widehat{a}}_{\mathit{TEi}}\text{exp}\left[i{k}_{\mathit{TEi}}\left(x-x_i\right)\right]+B_{\mathit{TEi}} {\widehat{b}}_{\mathit{TEi}}\text{exp}\left[-i{k}_{\mathit{TEi}}\left(x-x_i\right)\right]\hfill \\ & +A_{\mathit{TMi}} {\widehat{a}}_{\mathit{TMi}}\text{exp}\left[i{k}_{\mathit{TMi}}\left(x-x_i\right)\right]+B_{\mathit{TMi}} {\widehat{b}}_{\mathit{TMi}}\text{exp}\left[-i{k}_{\mathit{TMi}}\left(x-x_i\right)\right].\hfill \end{array}$$E21

Note that all these amplitudes have zero dephase at $x=x_i$. Finally, the vectors ${\widehat{a}}_{\mathit{TEi}}$ and ${\widehat{a}}_{\mathit{TMi}}$ are given in Eq. (3a), (3b) by ${\widehat{e}}_1$ and ${\widehat{e}}_2$, respectively, which we rewrite as

$${\widehat{a}}_{\mathit{TEi}}={\widehat{b}}_{\mathit{TEi}}=\left(0,k_z,-k_y\right),$$E22a

$${\widehat{a}}_{\mathit{TMi}}=\left(k_y^2+k_z^2,-k_y{k}_{\mathit{TMi}},-k_z{k}_{\mathit{TMi}}\right),$$E22b

$${\widehat{b}}_{\mathit{TMi}}=\left(k_y^2+k_z^2,k_y{k}_{\mathit{TMi}},k_z{k}_{\mathit{TMi}}\right).$$E22c

In Eq. (22a) and (22c) we have included the field vectors ${\widehat{b}}_{\mathit{TEi}}$ and ${\widehat{b}}_{\mathit{TMi}}$ which are associated with counter-propagating waves. Again, we point out that the vector fields ${\widehat{a}}_{\mathit{TEi}}$ and ${\widehat{b}}_{\mathit{TEi}}$ can be given in units of $k_0$, and the vectors ${\widehat{a}}_{\mathit{TMi}}$ and ${\widehat{b}}_{\mathit{TMi}}$ can be expressed in units of $k_0^2$. Similarly as performed in the previous section, the field function ${\mathbf{E}}^{\left(i\right)}$ can also be set in terms of the wave amplitudes $A{\prime }_{\mathit{TEi}}$, $A{\prime }_{\mathit{TEi}}$, $B{\prime }_{\mathit{TMi}},$ and $B{\prime }_{\mathit{TMi}}$ with zero dephase at $x=x_{i-1}$. As in the previous section, the amplitudes $A_{\mathit{qi}}$, $A{\prime }_{\mathit{qi}}$, $B_{\mathit{qi}},$ and $B{\prime }_{\mathit{qi}}$, where $q=\left\{\mathrm{TE},\mathrm{TM}\right\}$, are related by $A{\prime }_{\mathit{qi}}=A_{\mathit{qi}}\text{exp}\left(-i{k}_{\mathit{qi}}{w}_i\right)$ and $B{\prime }_{\mathit{qi}}=B_{\mathit{qi}}\text{exp}\left(i{k}_{\mathit{qi}}{w}_i\right)$. Following the same procedure seen earlier, we calculate the magnetic field in every elementary layer “$i$.” We finally obtain the following expression for the variation of the field along the $x$ direction, namely

$$\begin{array}{ll}\omega {\mu }_0{\mathbf{H}}^{\left(i\right)}\left(x\right)& =A_{\mathit{TEi}} {\widehat{c}}_{\mathit{TEi}}\text{exp}\left[i{k}_{\mathit{TEi}}\left(x-x_i\right)\right]+B_{\mathit{TEi}} {\widehat{d}}_{\mathit{TEi}}\text{exp}\left[-i{k}_{\mathit{TEi}}\left(x-x_i\right)\right]\hfill \\ & +A_{\mathit{TMi}} {\widehat{c}}_{\mathit{TMi}}\text{exp}\left[i{k}_{\mathit{TMi}}\left(x-x_i\right)\right]+B_{\mathit{TMi}} {\widehat{d}}_{\mathit{TMi}}\text{exp}\left[-i{k}_{\mathit{TMi}}\left(x-x_i\right)\right].\hfill \end{array}$$E23

In the previous equation we introduced the vector fields

$${\widehat{a}}_{\mathit{TEi}}={\widehat{b}}_{\mathit{TEi}}=\left(0,k_z,-k_y\right),$$E24a

$${\widehat{a}}_{\mathit{TMi}}=-{\widehat{c}}_{\mathit{TEi}}=-\left(-k_y^2-k_z^2,k_y{k}_{\mathit{TEi}},k_z{k}_{\mathit{TEi}}\right),$$E24b

$${\widehat{b}}_{\mathit{TMi}}=-{\widehat{d}}_{\mathit{TEi}}=\left(k_y^2+k_z^2,k_y{k}_{\mathit{TEi}},k_z{k}_{\mathit{TEi}}\right),$$E24c

$${\widehat{c}}_{\mathit{TMi}}={\widehat{d}}_{\mathit{TMi}}={ϵ}_i{k}_0^2\left(0,k_z,-k_y\right).$$E24d

Note that ${\widehat{c}}_{\mathit{TMi}}={ϵ}_i{k}_0^2{\widehat{a}}_{\mathit{TEi}}$. We conclude that the amplitudes $A_{\mathit{TEi}}$ and $B{\prime }_{\mathit{TEi}}$ are associated with TE-polarized waves along the $x$-axis, that is, ${\mathrm{TE}}^x$ waves, and $A_{\mathit{TMi}}$ and $B{\prime }_{\mathit{TMi}}$ are field amplitudes of ${\mathrm{TM}}^x$ waves.

2.3.4. Application of the boundary conditions

At a given interface $x=x_i$, the electromagnetic fields should accomplish the continuity boundary conditions. Note that we can write the four equations derived from the boundary conditions in the following matrix form:

$${\mathbf{D}}_i{\mathbf{v}}_i={\mathbf{D}}_{i+1}\mathbf{v}{\prime }_{i+1}.$$E25

In the previous matrix equation, we introduced the element

$${\mathbf{D}}_i=\left[\begin{array}{cccc}{k}_z& k_z& -k_y{k}_{\mathit{TMi}}& k_y{k}_{\mathit{TMi}}\\ -k_y& -k_y& -k_z{k}_{\mathit{TMi}}& k_z{k}_{\mathit{TMi}}\\ k_y{k}_{\mathit{TEi}}& -k_y{k}_{\mathit{TEi}}& {ϵ}_i{k}_0^2{k}_z& {ϵ}_i{k}_0^2{k}_z\\ k_z{k}_{\mathit{TEi}}& -k_z{k}_{\mathit{TEi}}& -{ϵ}_i{k}_0^2{k}_y& -{ϵ}_i{k}_0^2{k}_y\end{array}\right].$$E26

Finally, the amplitude vectors now are represented as

$${\mathbf{v}}_i=\left[\begin{array}{c}{A}_{\mathit{TEi}}\\ B_{\mathit{TEi}}\\ A_{\mathit{TMi}}\\ B_{\mathit{TMi}}\end{array}\right],\mathbf{v}{\prime }_i=\left[\begin{array}{c}A{\prime }_{\mathit{TEi}}\\ B{\prime }_{\mathit{TEi}}\\ A{\prime }_{\mathit{TMi}}\\ B{\prime }_{\mathit{TMi}}\end{array}\right].$$E27

The matrix formulation can also be used to relate the amplitude vector ${\mathbf{v}}_i$ with zero-phase shift at $x=x_i$ with the amplitude vector $\mathbf{v}{\prime }_i$ exhibiting zero dephase at $x=x_{i-1}$. For that purpose, we introduce the propagation matrix for the ${\mathrm{TE}}^x$ and ${\mathrm{TM}}^x$ modes

$${\mathbf{P}}_i=\left[\begin{array}{cccc}{e}^{-i{k}_{\mathit{TEi}}{w}_i}& 0& 0& 0\\ 0& e^{i{k}_{\mathit{TEi}}{w}_i}& 0& 0\\ 0& 0& e^{-i{k}_{\mathit{TMi}}{w}_i}& 0\\ 0& 0& 0& e^{i{k}_{\mathit{TMi}}{w}_i}\end{array}\right],$$E28

which takes into account the amplitude phase shift due the finite width $w_i=x_i-x_{i-1}$ of each layer. Explicitly we will write $\mathbf{v}{\prime }_i={\mathbf{P}}_i·{\mathbf{v}}_i$, which is formally the same as Eq. (18), previously derived for uniaxial media.

3. Surface modes in isotropic media

The main purpose of this chapter is the analysis of Dyakonov surface waves, which originally was formulated for an isotropic medium and a uniaxial crystal. However, this analysis is developed in Section 5. Here we introduce the most well-known surface waves arisen at the interface between isotropic media of different dielectric constants. In addition, these surface waves will play a relevant role when dealing with metal-dielectric multilayered structures.

The so-called surface plasmon polaritons are waves that propagate along the surface of a conductor, usually a metal [2, 8, 9]. These are essentially light waves that are trapped on the surface, evanescently confined in the perpendicular direction and caused by their interaction with the free electrons of the conductor, the latter oscillating in resonance with the electromagnetic field. To describe these wave fields, we use the matrix formalism applied in the vicinity of a single interface between two isotropic media with different dielectric permittivities.

Let us consider the propagation of bound waves on the interface between two semi-infinite media, which are denoted as medium 1 and medium 2 with dielectric permittivities ${ϵ}_1$ i ${ϵ}_2,$ respectively. This interface is located at $x_1=0$. For medium 1, the electric and magnetic fields varying along the $x$-axis are given by Eqs. (22a), (22b), (22c) and (24a), (24b), (24c), (24d) respectively.

As we are only interested in bound states; the elements of the remaining field vectors read as

$${\mathbf{v}}_1=\left[\begin{array}{c}0\\ B_{\mathit{TE}1}\\ 0\\ B_{\mathit{TM}1}\end{array}\right]\text{and}\mathbf{v}{\prime }_2=\left[\begin{array}{c}A{\prime }_{\mathit{TE}2}\\ 0\\ A{\prime }_{\mathit{TM}2}\\ 0\end{array}\right].$$E29

The values of the amplitudes $A_{\mathit{TE}1}$, $A_{\mathit{TM}1}$, $B{\prime }_{\mathit{TE}2}$, and $B{\prime }_{\mathit{TM}2}$ are identically zero in case of lack of interaction with external sources, as we assume here. Since we are dealing with bound states, the wavenumbers $k_{\mathit{TE}1}$, $k_{\mathit{TM}1}$, $k_{\mathit{TE}2,}$, and $k_{\mathit{TM}2}$ are purely imaginary. The application of the boundary conditions, ${\mathbf{D}}_1{\mathbf{v}}_1={\mathbf{D}}_2\mathbf{v}{\prime }_2$, gives us the following two equations:

$$0=\frac{{ϵ}_2{k}_{\mathit{TM}1}+{ϵ}_1{k}_{\mathit{TM}2}}{2{ϵ}_1{k}_{\mathit{TM}1}}A{\prime }_{\mathit{TM}2},$$E30a

$$B_{\mathit{TM}1}=\frac{{ϵ}_2{k}_{\mathit{TM}1}-{ϵ}_1{k}_{\mathit{TM}2}}{2{ϵ}_1{k}_{\mathit{TM}1}}A{\prime }_{\mathit{TM}2}.$$E30b

To accomplish Eq. (30a), the following equation must be satisfied, $k_{\mathit{TM}2}/k_{\mathit{TM}1}=-{ϵ}_2/{ϵ}_1$. This is the dispersion equation of the ${\mathrm{TM}}^x$-polarized surface modes. As $k_{\mathit{TM}1}$ and $k_{\mathit{TM}2}$ have a vanishing real part and a positive imaginary part, this requires that ${ϵ}_2/{ϵ}_1<0$, that is, one of the relative permittivities must be negative. Once we have obtained the relationship between the relative permittivities and the purely imaginary wavenumbers, we may obtain a new ligature involving the field amplitudes, $B_{\mathit{TM}1}=\left({ϵ}_2/{ϵ}_1\right)A{\prime }_{\mathit{TM}2}$, which have been derived from Eq. (30b). We may rewrite the dispersion equation, making use of the definition of $k_{\mathit{TMi}}$ given in Eq. (20), resulting:

$$k_{\mathit{SPP}}=k_0\sqrt{\frac{{ϵ}_1{ϵ}_2}{{ϵ}_1+{ϵ}_2}},$$E31

where $k_{\mathit{SPP}}=\sqrt{k_y^2+k_z^2}$ is the SPP wavenumber. In Figure 2 we represent Eq. (31) giving the dispersion relation for TM-bounded modes at frequencies lower than the surface plasmon frequency:

$${\omega }_{\mathit{SPP}}=\frac{{\omega }_p}{\sqrt{1+{ϵ}_1}}.$$E32

Figure 2.

SPP dispersion relation at the interface between a lossless Drude metal (${ϵ}_2\left(\omega \right)=1-{\omega }_p^2/{\omega }^2$) and a dielectric: air in blue (${ϵ}_1=1$) and GaAs in red (${ϵ}_1=11.55$). Note that $k_p={\omega }_p/c$.

At $\omega \to 0$ the SPP wavenumber tends to zero; however, when $\omega \to {\omega }_{\mathit{SPP}}$ we find that $k_{\mathit{SPP}}\to \infty$. In addition, radiative modes may arise at higher frequencies, typically $\omega \ge {\omega }_p$; in such cases, the wave field is not confined near the interface and it will lose its energy by radiation.

4. Multilayered plasmonic lattices

Wave propagation in periodic media can be treated as the motion of electrons in crystalline solids. In fact, formulation of the Kronig-Penney model used in the elementary band theory of solids is mathematically identical to that of the electromagnetic radiation in periodic layered media. Thus, some of the physical concepts used in material physics such as Bloch waves, Brillouin zones, and forbidden bands can also be used here. A periodic layered medium is equivalent to a one-dimensional lattice that is invariant under lattice translation.

Here we will treat the propagation of electromagnetic radiation in a simple periodic layered medium that consists of alternating layers of transparent nonmagnetic materials with different electric permittivities. The layers are set in a way that the $x$-axis points along the perpendicular direction of the layers. The permittivity profile is given by ${ϵ}_1$, for $x_0<x<x_1$, and ${ϵ}_2$, for $x_1<x<x_2$. In addition, the relative permittivity satisfies the condition of periodicity, $ϵ\left(x\right)=ϵ\left(x+\Lambda \right)$, where $w_1=x_1-x_0$ $\left(w_2=x_2-x_1\right)$ is the thickness of the layers of permittivity ${ϵ}_1$ $\left({ϵ}_2\right)$ and $\Lambda =w_1+w_2$ represents the period of the structure.

According to the Floquet theorem, solutions of the wave equation for a periodic medium may be set in the form ${\mathbf{E}}^{\left(i+2\right)}\left(x+\Lambda \right)={\mathbf{E}}^{\left(i\right)}\left(x\right)·\text{exp}\left(\mathit{iK\Lambda }\right)$. The constant $K$ is known as the Bloch wavenumber. The problem is thus that of determining $K$ and ${\mathbf{E}}^{\left(i\right)}\left(x\right)$. Finally, the equation

$$\text{cos}\left(K_q\Lambda \right)=\frac{1}{2}{\mathrm{tr}}_q.$$E33

represents the dispersion relation for TE and TM modes, written in a compact way. Representing $K_{\mathit{TE}}$ ($K_{\mathit{TM}}$), the Bloch wavenumber $K$ associated with the mode TE (TM), and writing

$${\mathrm{tr}}_{\mathit{TE}}=2\text{cos}\left(w_1{k}_{\mathrm{TE}1}\right)\text{cos}\left(w_2{k}_{\mathrm{TE}2}\right)-\frac{\left(k_{\mathrm{TE}1}^2+k_{\mathrm{TE}2}^2\right)}{k_{\mathrm{TE}1}{k}_{\mathrm{TE}2}}\text{sin}\left(w_1{k}_{\mathrm{TE}1}\right)\text{sin}\left(w_2{k}_{\mathrm{TE}2}\right),$$E34a

$${\mathrm{tr}}_{\mathit{TM}}=2\text{cos}\left(w_1{k}_{\mathrm{TM}1}\right)\text{cos}\left(w_2{k}_{\mathit{TM}2}\right)-\frac{\left({ϵ}_2^2{k}_{\mathit{TM}1}^2+{ϵ}_1^2{k}_{\mathit{TM}2}^2\right)}{{ϵ}_1{ϵ}_2{k}_{\mathit{TM}1}{k}_{\mathit{TM}2}}\text{sin}\left(w_1{k}_{\mathit{TM}1}\right)\text{sin}\left(w_2{k}_{\mathit{TM}2}\right).$$E34b

we find the dispersion equation of a binary periodic medium for each polarization.

Neglecting losses in the materials, regimes where $\mid t{r}_q\mid <2$ correspond to real $K_q$ and thus to propagating Bloch waves, when $\mid t{r}_q\mid >2$; however, $K_q=\mathit{m\pi }/\Lambda +i{K}_{\mathit{qi}}$, where $m$ is an integer and $K_{\mathit{qi}}$ is the imaginary part of $K_q$, which gives an evanescent behavior to the Bloch wave. These are the so-called forbidden bands of the periodic medium. The band edges are set for $\mid t{r}_q\mid =2$.

In Figure 3, the transverse wavenumber reads as

$$k_t=\sqrt{k_y^2+k_z^2}.$$E35

Figure 3.

Exact dispersion curves derived from Eq. (33) for a lossless Drude metal/GaAs composite medium with $\lambda =1.55\mu$ m, ${\omega }_p=13.94{\mathrm{fs}}^{-1}$ and ${ϵ}_1=11.55$, for (a) TM modes and (b) TE modes. GaAs layer is always 12 times thicker than the metallic layer.

For ultra-thin metallic layers, for instance, $w_m=3$ nm, the curves resembles ellipses and circumferences for TM and TE modes, respectively. In the particular case of TM modes, a secondary curve surges, in relation with the excitation of SPPs. For increasing values of the metal width, some deviations are evident. Reaching a given value of $w_m$, a large band gap centered at $k_t=0$ surges. Particularly for TM modes, at $K_{\mathit{TM}}=0$, the two curves collapse. This fact may be understood as that the symmetric and antisymmetric surface modes in the metallic layer approach for SPP in a single metal-dielectric interface, giving $k_t\to k_{\mathit{SPP}}$ [see Eq. (31)]. For illustration, the wavenumber of the SPP propagating on the interface of our Drude metal/GaAs materials yields $k_{\mathit{SPP}}=0.746·2\pi /\Lambda$, assuming that $\Lambda =325\mathrm{nm}$ (associated with the period of a multilayer with metal thickness $w_m=25$ nm).

4.1. Effective medium approach

For near-infrared and visible wavelengths, nanolayered metal-dielectric compounds enable a simplified description of the medium by using the long-wavelength approximation, which involves a homogenization of the structured metamaterial [10, 11]. The effective medium approach (EMA), as Rytov exposed in his seminal paper [12], involves representing MD multilayered metamaterial as an uniaxial plasmonic crystal, whose optical axis is normal to the layers (in our case, the x-axis is the optical axis), a procedure that requires the metallic elements to have a size of a few nanometers. This is caused by the fact that transparency of noble metals is restrained to a propagation distance not surpassing the metal skin depth. In Ref. to this point, recent development of nanofabrication technology makes it possible to create such subwavelength structures. Under this condition, the plasmonic lattice behaves as a uniaxial crystal characterized by a relative permittivity tensor $\mathbf{ϵ}={ϵ}_{\parallel }\left(\mathbf{x}\mathbf{\otimes }\mathbf{x}\right)+{ϵ}_{\perp }\left(\mathbf{y}\mathbf{\otimes }\mathbf{y}+\mathbf{z}\mathbf{\otimes }\mathbf{z}\right)$, where

$${ϵ}_{\parallel }=\frac{{ϵ}_1{ϵ}_2}{\left(1-f\right){ϵ}_2+f{ϵ}_1},$$E36

gives the permittivity along the optical axis, and ${ϵ}_{\perp }=\left(1-f\right){ϵ}_1+f{ϵ}_2$ corresponds to the permittivity in the transversal direction. In the previous equations, $f=w_2/\left(w_1+w_2\right)$, denotes the filling factor of medium “2” providing the metal rate in a unit cell. The dispersion equations given by the EMA are $k_x^2+k_t^2={ϵ}_{\perp }{k}_0^2$, corresponding to TE ($o$-) waves, where $k_x$ represents the Bloch wavenumber $K_{\mathit{TE}}$, and for TM ($e$-) waves, we have

$$\frac{k_x^2}{{ϵ}_{\perp }}+\frac{k_t^2}{{ϵ}_{\parallel }}=k_0^2,$$E37

being now $k_x$, the Bloch wavenumber $K_{\mathit{TM}}$.

The validity of the EMA is related on the assumption that the period $\Lambda$ is much shorter than the wavelength, that is, $\Lambda \ll {\lambda }_0$. Apparently, Eq. (33) is in good agreement with the EMA in the vicinity of $k_x=0$ for ${\mathrm{TM}}^x$ waves only. In contrast, propagation along the $x$-axis, where $k_y=k_z=0$, results in large discrepancies. Even small-filling factors of the metallic composite lead to enormous birefringences. Such metamaterials enlarge the birefringence of the effective-uniaxial crystal at least in one order of magnitude in comparison with values shown in Table 1. However, the size of birefringence displayed by extraordinary waves is reduced if $w_2$ increases. On the other hand, the isotropy of the isofrequency curve is practically conserved for ordinary waves.

4.2. Hyperbolic media

As we have seen in Section 3.3, nanolayered metal-dielectric compounds behave like plasmonic crystals enabling a simplified description of the medium by using the long-wavelength approximation [10, 11, 12]. Under certain conditions, the permittivity of the medium set in the form of a second-rank tensor includes elements of opposite signs, leading to a metamaterial of extreme anisotropy [13, 14]. This class of nanostructured media with hyperbolic dispersion is promising metamaterials with a plethora of practical applications from biosensing to fluorescence engineering [15].

Type I hyperbolic media refers to a special kind of uniaxially anisotropic media, that can be described by a permittivity tensor where element ${ϵ}_{\parallel }$ is negative and ${ϵ}_{\perp }$ is positive. In this case, Eq. (37) leads to a two-sheet hyperboloid. Type II hyperbolic media lead to positive ${ϵ}_{\parallel }$ and negative ${ϵ}_{\perp }$, and Eq. (37) gives us a one-sheet hyperboloid [16]. The fulfillment of hyperbolic dispersion allows wave propagation over a wide spatial spectrum that would be evanescent in an ordinary isotropic dielectric. At the optical range, hyperbolic media can be manufactured with metal-dielectric multilayers or metallic nanowires. Multilayered hyperbolic metamaterials at an optical range take advantage of the wide frequency band in which metals exhibit negative permittivity and support plasmonic modes.

The system under analysis is a periodic binary medium, where we take as medium 1 a transparent dielectric medium that is ideally nondispersive. In our three numerical simulations we take a lossless Drude metal where its permittivity is ${ϵ}_2=1-{\Omega }^{-2}$ and dielectric media with permittivities: (a) ${ϵ}_1=1$, (b) ${ϵ}_1=2.25$ and (c) ${ϵ}_1=11.55$. Note that frequencies can be expressed in units of the plasma frequency, $\Omega =\omega /{\omega }_p$.

In Figure 4(a), we represent the permittivities ${ϵ}_{\parallel }$ and ${ϵ}_{\perp }$ of our plasmonic crystals for a wide range of frequencies. Note that the metal-filling factor takes control on the dissipative effects in the metamaterial; accordingly low values of $f$ are of great convenience. We set $f=1/4$ in our numerical simulations. For low frequencies, $\Omega \ll 1$, the following approximations can be used: ${ϵ}_{\perp }\approx f{ϵ}_2<0$ and $0<{ϵ}_{\parallel }\approx {ϵ}_1/\left(1-f\right)$. Therefore, propagating ${\mathrm{TE}}^x$ modes ($E_x=0$) cannot exist in the bulk crystal since it behaves like a metal in these circumstances. On the other hand, ${\mathrm{TM}}^x$ waves propagate following the spatial dispersion curve of Eq. (37). This is a characteristic of Type II hyperbolic media. As mentioned earlier, Eq. (37) denotes a hyperboloid of one sheet (see Figure 4(b) for $\Omega =0.20$).

Figure 4.

(a) Variation of relative permittivities ${ϵ}_{\parallel }$ (blue solid line) and ${ϵ}_{\perp }$ (magenta solid line) as a function of normalized frequency $\Omega$, for the plasmonic crystal including a lossless Drude metal and fused silica (${ϵ}_1=2.25$) as dielectric material. ${\Omega }_1$ and ${\Omega }_2$ yield 0.359 and 0.759, respectively. (b) and (c) plot Eq. (37) in the $k_x{k}_y$ plane for extraordinary waves (${\mathrm{TM}}^x$ modes) for a plasmonic effective crystal, including a dielectric of permittivity ${ϵ}_1=2.25$, in the range $\Omega <1$. Solid line corresponds to $k_z=0$ and shaded regions are associated with harmonic waves with $k_z>0$ (nonevanescent fields).

Furthermore, the hyperbolic dispersion exists up to a frequency

$${\Omega }_1=\frac{1}{\sqrt{1+{ϵ}_1\left(1-f\right)/f}},$$E38

for which ${ϵ}_{\perp }=0$. For slightly higher frequencies, both ${ϵ}_{\parallel }$ and ${ϵ}_{\perp }$ are positive and Eq. (37) becomes an ellipsoid of revolution. Since its minor semi-axis is $\Omega \sqrt{{ϵ}_{\perp }}$, the periodic multilayer simulates a uniaxial medium with positive birefringence. Raising the frequency even more, ${ϵ}_{\parallel }$ diverges at

$${\Omega }_2=\frac{1}{\sqrt{1+{ϵ}_{1}f/\left(1-f\right)}},$$E39

leading to the so-called canalization regime. In general, ${\Omega }_1<{\Omega }_2$ provided that $f<1/2$. Beyond ${\Omega }_2$, Eq. (37) turns to a hyperboloidal shape. In the range ${\Omega }_2<\Omega <1;$ however, the dispersion curve has two sheets (Type I hyperbolic medium). Figure 4(c) illustrates this case. Note that the upper limit of this hyperbolic band is determined by the condition ${ϵ}_{\parallel }=0$ or in an equivalent way, ${ϵ}_2=0$, occurring at the plasma frequency.

5. Dyakonov surface waves

Dyakonov surface waves (DSWs) are another kind of surface waves, supported at the interface between an optically isotropic medium and a uniaxial-birefringent material. In the original work by Dyakonov (English version was reported in 1988 [3]), the optical axis of the uniaxial medium was assumed in-plane with respect to the interface. This is the case we deal with here.

The importance of DSWs for integrated optical applications, such as sensing and nano-waveguiding, was appreciated in a series of papers [17, 18]. Indeed Dyakonov-like surface waves also emerge in the case that a biaxial crystal [19] or a structurally chiral material [20] takes the place of the uniaxial medium. The case of metal-dielectric (MD) multilayers as structurally anisotropic media is especially convenient since small-filling fractions of the metallic inclusions enable metamaterials with an enormous birefringence, thus enhancing density of DSWs and relaxing their prominent directivity [21, 22, 23].

5.1. Dispersion equation of DSWs

The system under study is the plotted in Figure 5, where we have two semi-infinite media, one of them is isotropic and the second one is an MD lattice. In our case, the indices “1” and “2” make reference to the plasmonic lattice and the isotropic medium, respectively. We have previously reported a comprehensive analysis of this case in [4]. As we have seen earlier, the plasmonic lattice can be taken as an effective uniaxial crystal. In this case, the permittivity along its optical axis, ${ϵ}_{z1}={ϵ}_{\parallel ,}$ is given by Eq. (36); also, the permittivity in the transverse direction ${ϵ}_{x1}={ϵ}_{y1}={ϵ}_{\perp }$ may be appropriately averaged. From hereon, the permittivity ${ϵ}_2$ of the isotropic medium in $x>0$ will be denoted by $ϵ$. Note that our analysis serves for natural birefringent materials characterized by permittivities ${ϵ}_{\parallel }$ and ${ϵ}_{\perp }$.

Figure 5.

Schematic setup under study, consisting in a semi-infinite dielectric-metal superlattice $\left(x<0\right)$ and an isotropic substrate $\left(x>0\right)$.

Since we treat the plasmonic lattice as a uniaxial crystal, we may establish analytically the diffraction equation that gives the 2D wave vector ${\mathbf{k}}_D=\left[0,k_y,k_z\right]$ in $x=0$. For that purpose, we follow Dyakonov [3] by considering hybrid-polarized surface modes. In the isotropic medium we consider ${\mathrm{TE}}^x$ ($E_x=0$) and ${\mathrm{TM}}^x$ ($H_x=0$) waves whose wave vectors have the same real components $k_y$ and $k_z$ in the plane $x=0$. Therefore the electric field in both media may be set as

$${\mathbf{E}}_{\mathit{tot}}=\mathbf{E}\left(x\right)\text{exp}\left(i{k}_{y}y+i{k}_{z}z-\mathit{i\omega t}\right).$$E40

Moreover, these fields are evanescent in the isotropic medium and in the superlattice. In the anisotropic medium ( $x<0$) the evanescent electric amplitude can be written as

$$\mathbf{E}\left(x\right)=B_{o1} {\widehat{b}}_{o1}\text{exp}\left(-i{k}_{o1}x\right)+B_{e1} {\widehat{b}}_{e1}\text{exp}\left(-i{k}_{e1}x\right),$$E41

where the ordinary and extraordinary waves in the effective uniaxial medium decay exponentially with rates given by ${\kappa }_o=-i{k}_{o1}$ and ${\kappa }_e=-i{k}_{e1}$, respectively. Taking the formulation given in Section 2.3.1, the amplitudes $A_{o1}$ and $A_{e1}$ are identically zero. In the isotropic medium ($x>0$) the amplitude of the electric field is

$$\mathbf{E}\left(x\right)=A{\prime }_{\mathit{TE}2} {\widehat{a}}_{\mathit{TE}2}\text{exp}\left(i{k}_{\mathit{TE}2}x\right)+A{\prime }_{\mathit{TM}2} {\widehat{a}}_{\mathit{TM}2}\text{exp}\left(i{k}_{\mathit{TM}2}x\right),$$E42

where the evanescent decay for TE and TM modes is $\kappa =-i{k}_{\mathit{TE}2}=-i{k}_{\mathit{TM}2}$. Now the amplitudes $B{\prime }_{\mathit{TE}2}$ and $B{\prime }_{\mathit{TM}2}$ are zero.

Once we have the amplitudes in both sides of the interface, we apply the boundary conditions at $x=0$, ${\mathbf{D}}_1·{\mathbf{v}}_1={\mathbf{D}}_2·\mathbf{v}{\prime }_2$, where ${\mathbf{D}}_1$ is provided by Eq. (16), ${\mathbf{D}}_2$ from Eq. (26), ${\mathbf{v}}_1$ from Eq. (17), and $\mathbf{v}{\prime }_2$ from Eq. (27). This equation reduces to:

$$\left[\begin{array}{c}0\\ B_{o1}\\ 0\\ B_{e1}\end{array}\right]={\mathbf{M}}_h\left[\begin{array}{c}A{\prime }_{\mathit{TE}2}\\ 0\\ A{\prime }_{\mathit{TM}2}\\ 0\end{array}\right],$$E43

where the transmission matrix ${\mathbf{M}}_h={\mathbf{D}}_1^{-1}·{\mathbf{D}}_2$ establishes a relationship between the amplitudes of hybrid polarization modes. Using the elements $M_{\mathit{ij}}$ of the matrix ${\mathbf{M}}_h$, and defining ${\mathbf{M}}_i$ and ${\mathbf{M}}_a$ as

$${\mathbf{M}}_i=\left[\begin{array}{cc}{M}_{11}& M_{13}\\ M_{31}& M_{33}\end{array}\right]=\left[\begin{array}{cc}-\frac{i{k}_z\left({\kappa }_o+\kappa \right)}{2\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}& \frac{k_0^2{k}_y\left(ϵ{\kappa }_o+\kappa {ϵ}_{\perp }\right)}{2{\kappa }_o\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}\\ -\frac{\left({\kappa }_e+\kappa \right)k_y}{2{\kappa }_e\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}& -\frac{i{k}_z\left(\kappa {\kappa }_e{ϵ}_{\perp }+ϵ{\kappa }_o^2\right)}{2{\kappa }_e{ϵ}_{\perp }\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}\end{array}\right],$$E44a

$${\mathbf{M}}_a=\left[\begin{array}{cc}{M}_{21}& M_{23}\\ M_{41}& M_{43}\end{array}\right]=\left[\begin{array}{cc}\frac{i{k}_z\left({\kappa }_o-\kappa \right)}{2\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}& \frac{k_0^2{k}_y\left(ϵ{\kappa }_o-\kappa {ϵ}_{\perp }\right)}{2{\kappa }_o\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}\\ \frac{\left(\kappa -{\kappa }_e\right)k_y}{2{\kappa }_e\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}& -\frac{i{k}_z\left(\kappa {\kappa }_e{ϵ}_{\perp }-ϵ{\kappa }_o^2\right)}{2{\kappa }_e{ϵ}_{\perp }\left(k_z^2-k_0^2{ϵ}_{\perp }\right)}\end{array}\right].$$E44b

Eq. (43) can be rewritten as a set of two independent matrix equations, namely

$$\left[\begin{array}{c}0\\ 0\end{array}\right]={\mathbf{M}}_i·\left[\begin{array}{c}{A}_{\mathit{TE}2}^{\prime }\\ A_{\mathit{TM}2}^{\prime }\end{array}\right],$$E45a

$$\left[\begin{array}{c}{B}_{o1}\\ B_{e1}\end{array}\right]={\mathbf{M}}_a·\left[\begin{array}{c}{A}_{\mathit{TE}2}^{\prime }\\ A_{\mathit{TM}2}^{\prime }\end{array}\right].$$E45b

Note that ${\mathbf{M}}_i$ governs the amplitudes $A{\prime }_{\mathit{TE}2}$ and $A{\prime }_{\mathit{TM}2}$ of the isotropic medium, and ${\mathbf{M}}_a$ (also ${\mathbf{M}}_h$) may be used to determine the amplitudes $B_{o1}$ and $B_{e1}$ of the anisotropic medium.

Dyakonov equation is obtained by means of letting the determinant of ${\mathbf{M}}_i$ equal to zero, giving

$$k_0^2{k}_y^2{ϵ}_{\perp }\left(\kappa +{\kappa }_e\right)\left(ϵ{\kappa }_o+{ϵ}_{\perp }\kappa \right)={\kappa }_o{k}_z^2\left(\kappa +{\kappa }_o\right)\left(ϵ{\kappa }_o^2+{ϵ}_{\perp }\kappa {\kappa }_e\right),$$E46

which provides a spectral map of allowed values $\left(k_y,k_z\right)$. After fairly tedious algebraic transformations we can reduce Eq. (46) to a more convenient form [3].

$$\left(\kappa +{\kappa }_e\right)\left(\kappa +{\kappa }_o\right)\left(ϵ{\kappa }_o+{ϵ}_{\perp }{\kappa }_e\right)=\left({ϵ}_{\parallel }-ϵ\right)\left(ϵ-{ϵ}_{\perp }\right)k_0^2{\kappa }_o.$$E47

Assuming that ${ϵ}_{\parallel }$, ${ϵ}_{\perp }$ and all decay rates are positive, the additional restriction ${ϵ}_{\perp }<ϵ<{ϵ}_{\parallel }$ can be deduced for the existence of surface waves. As a consequence, positive birefringence is mandatory to ensure a stationary solution of Maxwell’s equations. Therefore, layered superlattices supporting Dyakonov-like surface waves cannot be formed by all dielectric materials [21].

5.2. DSWs in nano-engineered materials

To illustrate the difference between using conventional birefringent materials and plasmonic crystals, we solve Eq. (47) for liquid crystal E7 with ${ϵ}_{\parallel }=2.98$ i ${ϵ}_{\perp }=2.31$ at a wavelength of ${\lambda }_0=1.55\mu$ m and N-BAK1 substrate of dielectric constant $ϵ=2.42$. In this case, DSWs propagate in a narrow angular region $\mathit{\Delta \theta }={\theta }_{\mathit{max}}-{\theta }_{\mathit{min}}$, where $\theta$ stands for the angle between the in-plane vector $\left(k_y,k_z\right)$ and the optical axis. More specifically, the angular range yields $\mathit{\Delta \theta }={0.92}^{\circ }$ around the mean angle $\overline{\theta }={26.6}^{\circ }$. It appears that the resulting angular range $\mathit{\Delta \theta }$ is short. But this range would become much smaller when using other optical crystals like quartz, exhibiting a common birefringence. In order to gain in angular extent $\mathit{\Delta \theta },$ we consider a GaAs-Ag crystal (${ϵ}_1=12.5$ and ${ϵ}_2=-103.3$, where we neglect losses) that leads to values of ${ϵ}_{\parallel }=14.08$, derived from Eq. (36), and ${ϵ}_{\perp }=0.92$ with a metal-filling factor $f=0.10$. Form birefringence now yields $\mathit{\Delta n}=2.79$. Moreover, solutions to Eq. (47) can be found in the region of angles comprised between ${\theta }_{\mathit{min}}={39.0}^{\circ }$ and ${\theta }_{\mathit{max}}={71.3}^{\circ }$. It is important to point out that the total angular range of existence of DSWs, $\mathit{\Delta \theta }={32.3}^{\circ }$, grows by more than an order of magnitude. Figure 6 shows the dispersion curve for DSWs; one can observe that ${\theta }_{\mathit{min}}$ is attained under the condition $\kappa =0$ (red solid line), where ${\mathrm{TE}}^x$ and ${\mathrm{TM}}^x$ waves are uniform in the substrate $x>0$. Looking at the other side of the dispersion curve, ${\theta }_{\mathit{max}}$ is established by $k_e=0$, shown as a black solid line, for which the extraordinary wave will not decay spatially at $x\to -\infty$.

Figure 6.

Dispersion Eq. (47) for DSWs (dotted-dashed line) propagating on the interface of a semi-infinite dielectric-metal lattice with metal filling factor f = 0.1. The solid elliptical line and the black dashed line are associated with homogeneous extraordinary waves $\left({\kappa }_e=0\right)$ and homogeneous ordinary waves $\left({\kappa }_o=0\right)$, respectively. The isofrequency curve $\left(\kappa =0\right)$ of isotropic N-BAK1 are represented in by the red solid line, which applies for ${\mathrm{TE}}^x$ and ${\mathrm{TM}}^x$ waves. Adapted with permission from [4] of copyright ©2013 IEEE photonics society.

Consequently, the solution for Eq. (47) can be traced near the curves $\kappa =0$ and ${\kappa }_e=0$; thus, DSWs are always found close to the crosspoint $P_0\left(k_{y0},k_{z0}\right)$ of both curves.

5.2.1. Nonlocal effects

As we discussed in Section 4.1, the EMA is limited to metallic slabs’ width $w_m\ll {\lambda }_0$. However, this condition must be taken into account with care, since the skin depth of noble metals is extremely short, $\delta \approx c/{\omega }_p$. For instance, we estimate $\delta =24$ nm in the case of silver. If the metal thickness is comparable to its skin depth, the EMA will substantially deviate from exact calculations. Note that experimental studies from multilayer optics rarely incorporate metallic slabs with a thickness below 10 nm.

We emphasize that moderate changes in the birefringence of the plasmonic crystal will substantially affect the existence of DSWs. More specifically, an enlargement of ${ϵ}_{\perp }$ driven by increasing $w_m$, provided $f$ is fixed (see Section 4.1), will lead to a significant modification of the DSW dispersion curves. Ultimately, this phenomenon is clearly attributed to nonlocal effects in the effective-medium response of nanolayered metamaterials [24], which is associated with a strong variation of the fields on the scale of a single layer.

We conclude that, in order to excite DSWs, one may counterbalance the decrease of birefringence in the plasmonic lattice by means of a dielectric substrate of higher index of refraction. To illustrate this matter, the dispersion Eq. (47) for DSWs is represented in Figure 7, in addition to using the values of ${ϵ}_{\parallel }$ and ${ϵ}_{\perp }$ from nonlocal estimators [4]. When $w_m$ grows but $f$ is kept fixed, the dispersion curve of the Dyakonov surface waves tends to approach the optic axis.

Figure 7.

Solutions to Dyakonov equation, drawn in dotted-dashed lines, for a MD lattice with the same filling factor $f=0.10$ but different $w_m$, using estimates from the nonlocal birefringence approach described in [4]. Note that the red solid line designates the isofrequency curve of isotropic substrate (a) N-BAK1 and (b) P-SF68. Adapted with permission from [4] of copyright ©2013 IEEE photonics society.

For an N-BAK1 substrate, as shown in Figure 7(a), ${\theta }_{\mathit{max}}={68.2}^{\circ }$ and ${58.7}^{\circ }$ for $w_m=3$ nm and 6 nm, respectively. Also, ${\theta }_{\mathit{min}}={37.6}^{\circ }$ and ${32.1}^{\circ }$ for these two cases. As a consequence the angular range $\mathit{\Delta \theta }$ shrinks when $w_m$ increases. In the limit ${ϵ}_{\perp }\to ϵ$, which occurs for $w_m=10.3$ nm using a N-BAK1 substrate, DSWs are not supported at the interface of the MD lattice and the isotropic dielectric. A substrate with greater relative permittivity $ϵ$ would be necessary. For example, if we use a substrate with greater relative permittivity $ϵ$ as P-SF68 [see Figure 7(b)], DSWs exist for $w_m=12$ nm with an angular range $\mathit{\Delta \theta }={12.5}^{\circ }$.

5.2.2. Dissipative effects

Up to now, we have avoided another important aspect of plasmonic devices namely dissipation in metallic elements. In this regard, effective permittivities are fundamentally complex, and consequently the Dyakonov Eq. (47) is expected to give complex values of $\left(k_y,k_z\right)$. This procedure has been discussed by Sorni et al. [25] recently. In order to tackle this problem, we evaluate numerically the value of the Bloch wavenumber $k_z$ for a given real value $k_y$. The spatial frequency $k_z$ becomes complex since $\text{Im}\left[{ϵ}_m\right]=8.1$. As a consequence, the surface wave cannot propagate indefinitely, undergoing an energy attenuation given by $l={\left(2\text{Im}\left[k_z\right]\right)}^{-1}$. Furthermore, we naturally assume that the real part of the parameters $\kappa ,{\kappa }_{o,}$ and ${\kappa }_e$ are all positive. These positive values correlate with a decay at $\mid x\mid \to \infty$ and thus with a confinement of the wave near $x=0$.

Figure 8(a) depicts the dispersion curve corresponding to dissipative DSWs, for the case of a plasmonic MD lattice with $f=0.10$ and $w_m=12$ nm. We used a commercial software (COMSOL Multiphysics) based on the finite-element method (FEM) in order to perform our numerical simulations. We cannot observe surface waves by setting an N-BAK1 substrate with $n=1.56$, suggesting that this is a retardation effect. More specifically, Figure 8(a) shows the isofrequency curves when $n=1.95$, corresponding to P-SF68.

Figure 8.

(a) Isofrequency curve that corresponds to hybrid surface waves existing at the boundary between a semi-infinite P-SF68 substrate and a lossy MD superlattice of $f=0.10$ and $w_m=12$ nm. (b) Ratio of Im($k_z$) over re($k_z$) representing dissipation effects in the propagation of DSWs. (c) Three contour plots of the magnetic field $\mid H_x\mid$ computed using the finite-element method. The superlattice is set on the left, for which only one period is represented. Capital letters a, B, and C, designate the transverse spatial frequencies $k_y=0.8k_0$, $1.2k_0$, and $1.6k_0$, respectively. Adapted with permission from [4] of copyright ©2013 IEEE photonics society.

We observe that the dispersion curve for dissipative DSWs is flatter and larger than the curve obtained by neglecting losses. Specifically ${\theta }_{\mathit{max}}={49.9}^{\circ }$ and ${\theta }_{\mathit{min}}={23.7}^{\circ }$, giving an angular range $\mathit{\Delta \theta }={26.2}^{\circ }$. Figure 8(b) shows Im($k_z$)/Re($k_z$) in the range of existence of the surface waves. In these two figures, capital letters A, B, and C designate the transverse spatial frequencies $k_y=0.8k_0$, $1.2k_0$, and $1.6k_0$, respectively. Figure 8(c) shows the magnetic field $\mid H_x\mid$ for the three different cases denoted by capital letters A, B, and C. Note that in the case of paraxial surface waves, for which $k_y$ reaches a minimum value (case A), one achieves $\text{Im}\left(k_z\right)\ll \text{Re}\left(k_z\right)$ as depicted in Figure 8(b). This is induced by an enormous shift undergone by the field maximum in the direction to the isotropic medium, as shown in Figure 8(c), where dissipation effects are barely disadvantageous on surface-wave propagation. Moreover, this would be consistent with a condition $\text{Re}\left(\kappa \right)\ll \text{Re}\left({\kappa }_e\right)$. On the other hand, for nonparaxial waves, having the largest values of $k_y$, the fields show slow energy decay inside the plasmonic superlattice. In case C, the magnetic field $\mid H_x\mid$ is localized around the metallic layer and takes significant values far from the boundary of the substrate. As a consequence, losses in the metal translate into a significant rise in the values of Im($k_z$).

5.3. New families of DSWs in lossy media

In this section we carry out a thorough analysis of DSWs that takes place in lossy uniaxial metamaterials. Special emphasis is put when the effective-medium approach induces satisfactory results. The introduction of losses leads to a transformation of the isofrequency curves, which deviates from spheres and ellipsoids, as commonly considered by ordinary and extraordinary waves, respectively. As a consequence, one can find two different families of surface waves as reported by Sorni et al. [25]. One family of surface waves is directly related with the well-known solutions derived by Dyakonov [3]. Importantly, the existence of a new family of surface waves is revealed, closely connected to the presence of losses in the uniaxial effective crystal. We point out that the solutions to Dyakonov equation presented earlier are partial ones insofar as the $z$-component of the wavevector is kept real valued. Nevertheless, the whole set of solutions includes all possible wavevectors that feature a complex-valued $k_x$.

6. Dyakonov surface waves in hyperbolic media

In this section we perform a thorough analysis of DSWs taking place in semi-infinite MD lattices exhibiting hyperbolic dispersion. Part of this section was previously reported by Zapata-Rodríguez et al. [26]; we point out that recently further studies on DSW in hyperbolic metamaterials have been reported by other authors [27]. Our approach puts emphasis on the EMA. Under these conditions, different regimes can be found including DSWs with nonhyperbolic dispersion. The system under analysis is again as depicted in Figure 5. For simplicity, we assume that dielectric materials are nondispersive; indeed, we set $ϵ=1$ and ${ϵ}_d=2.25$ in our numerical simulations. Furthermore, Drude metals are included, and frequencies will be expressed in units of its plasma frequency, $\Omega =\omega /{\omega }_p$. Again, Dyakonov Eq. (47) provides the spectral map of wave vectors ${\mathbf{k}}_D=\left[0,k_y,k_z\right]$. Note that in this section, spatial frequencies will be expressed in units of $k_p$.

In the special case of the surface wave propagation perpendicular to the optical axis ($k_z=0$), Eq. (47) reveals the following solution: $ϵ{\kappa }_o+{ϵ}_{\perp }\kappa =0$. In the case: ${ϵ}_{\perp }<0$ and $ϵ<\mid {ϵ}_{\perp }\mid$, this equation has the well-known solution

$$k_y=\Omega \sqrt{\frac{ϵ{ϵ}_{\perp }}{ϵ+{ϵ}_{\perp }}},$$E48

which resembles the dispersion equation of conventional SPPs [see Eq. (31)]. Here we have purely ${\mathrm{TM}}^x$ polarized waves, as expected. Note that no solutions to Eq. (47) can be found, in the form of surface waves, considering wave propagation parallel to the optical axis ($k_y=0$) for hyperbolic metamaterials: ${ϵ}_{\perp }{ϵ}_{\parallel }<0$. That means that a threshold value of $k_y$ can be found for the existence of surface waves.

Next we describe a specific configuration governing DSWs, subject to a low value of the refractive index $n=\sqrt{ϵ}$, namely $ϵ<{ϵ}_{\parallel }$ ($ϵ<{ϵ}_{\perp }$) occurring at low and moderate frequencies; other cases are treated elsewhere [26]. In the effective-uniaxial medium, it is easy to realize that $\kappa <{\kappa }_o$ and also ${\kappa }_e<{\kappa }_o$. Under these circumstances, all brackets in Dyakonov Eq. (47) are positive provided $ϵ{\kappa }_o+{ϵ}_{\perp }{\kappa }_e>0$. This happens within the spectral band ${\Omega }_0<\Omega <{\Omega }_1$, where

$${\Omega }_0=\frac{1}{\sqrt{1+ϵ/f+{ϵ}_d\left(1-f\right)/f}}.$$E49

Note that ${\Omega }_0=0.292$ in our numerical simulation.

In Figure 9(a) and (b), we illustrate the dispersion equation of DSWs for two different frequencies within the spectral range $0<\Omega <{\Omega }_0$. In these cases, the dispersion curve approaches a hyperbola. We find a bandgap around $k_z=0$ in Figure 9(a), unlike what occurs in (b). Note that hybrid solutions near $k_z=0$ are additionally constrained to the condition $k_y\ge \Omega \sqrt{{ϵ}_{\parallel }}$ [see also Eq. (48)], which is a necessary condition for ${\kappa }_e$ to exhibit real and positive values. We consider the quasi-static regime ($\Omega \to 0$) where $\mid {\mathbf{k}}_D\mid =k_D\gg \Omega$ to determine the asymptotes of the hyperbolic-like DSW dispersion curve. Under this approximation, $\kappa =k_D$, ${\kappa }_o=k_D$ and ${\kappa }_e=\Theta k_D$, where

$$\Theta =\sqrt{{\text{cos}}^2\theta +\left({ϵ}_{\parallel }/{ϵ}_{\perp }\right){\text{sin}}^2\theta },$$E50

being $k_y=k_D\text{cos}\theta$ and $k_z=k_D\text{sin}\theta$. These asymptotes establish a canalization regime leading to a collective directional propagation of DSW beams [4, 28]. At this point it is necessary to remind that the asymptotes of the $e$-waves dispersion curve, in the $k_y{k}_z$ plane, have slopes satisfying the condition ${\theta }_D<{\theta }_e$, as illustrated in Figure 9(b).

Figure 9.

Solutions to Eq. (47), drawn in solid line, providing the spatial dispersion of DSWs which can exist in the arrangement of Figure 5, at different frequencies: (a) $\Omega =0.20$, (b) $\Omega =0.28$, and (c) $\Omega =0.85$. Here, the isotropic medium is air and the multi-layered metamaterial has a filling factor $f=0.25$ and. Also, we include equations $\kappa =0$ (dotted line) and ${\kappa }_e=0$ (dashed line). Adapted with permission from [26] of copyright ©2013 Optical Society of America.

In the high-frequency band ${\Omega }_2<\Omega <1$ we find that ${ϵ}_{\parallel }<0<{ϵ}_{\perp }$, as occurs in Figure 9(c). Note the relevant proximity of DSW dispersion curve to ${\kappa }_e=0$, the same way we also find in Figure 9(a) and (b). Conversely it crosses the $e$-wave hyperbolic curve at two different points, where solutions to the Dyakonov equation begin and end, respectively. It is clear that the angular range of DSWs now turns to be significantly low.

6.1. DSWs in band-gap hyperbolic media

In previous sections we demonstrated that the presence of metallic nano-elements leads to nonlocal effects and dissipation effects which reshape the propagation dynamics of the surface signal. Here, we briefly discuss the extraordinary favorable conditions which may appear in band-gap metal–insulator-layered media for the existence of DSWs. As reported thoroughly by Miret et al. [29], engineering secondary bands by tuning the plasmonic-crystal geometry may lead to a controlled optical anisotropy, which is markedly dissimilar to the prescribed hyperbolic regime that is derived by the EMA, however, assisting the presence of DSWs on the interface between such hyperbolic metamaterial and an insulator.

In particular, a surface wave propagating on an Ag-Ge grating was considered, where the environment medium that is set above the metallic grating is formed by ${\mathrm{SiO}}_2$. If the metal-filling factor was $f=0.25$, the effective permittivities of the anisotropic metamaterial would be estimated as ${ϵ}_{\perp }=-11.48+i2.05$ and ${ϵ}_{\Vert }=25.96+i0.14$ at a wavelength of ${\lambda }_0=1550$ nm. Disregarding losses, the DSW dispersion curve describes an incomplete hyperbolic curve, finding an endpoint under the condition ${\kappa }_e=0$, where the extraordinary wave breaks its confinement in the vicinities of the isotropic-uniaxial interface [26].

Considering now a realistic nanostructure consisting of Ag layers of $w_2=40$ nm interspersed between Ge layers of $w_1=120$ nm, thus maintaining a metal-filling factor of $f=0.25$ as analyzed earlier, a first TM band with hyperbolic-like characteristics dominates at high in-plane frequencies $k_t$. Additionally, a second band emerges for TM Bloch modes, which exhibits a moderate anisotropy, demonstrating near-elliptical dispersion curves (with positive effective permittivities) and positive birefringence. Furthermore, TE modal dispersion is roughly isotropic. Therefore, satisfactory conditions are found near the second TM band for the existence of Dyakonov-like surface waves. Finally, in order to numerically obtain the dispersion curves and wave fields associated with DSWs, one may follow the same computational procedure followed by Zapata-Rodríguez et al. and [4, 26].

7. Summary

In this chapter we provide several methods to analytically calculate and numerically simulate modal propagation of DSWs governed by material anisotropy. We focused on the spatial properties of DSWs at optical and telecom wavelengths, particularly using uniaxial metamaterials formed of dielectric and metallic nanolayers. We developed an electromagnetic matrix procedure enabling different aspects reviewed in this chapter, specially adapted to complex multilayered configurations. The EMA results are particularly appropriate for the characterization of the form birefringence of a multilayered nanostructure, though limitations driven by the layers width have been discussed. Through a rigorous full-wave analysis, we showed that hybrid-polarized surface waves may propagate obliquely at the boundary between a plasmonic bilayer superlattice and an isotropic loss-free material. We revealed that realistic widths of the slabs might lead to solutions which deviate significantly from the results derived directly from the EMA and Dyakonov analysis. Finally, we showed that excitation of DSWs at the boundary of an isotropic dielectric and a hyperbolic metamaterial enables a distinct regime of propagation. It is important to note that the properties of the resulting bound states change drastically with the index of refraction of the surrounding medium, suggesting potential applications in chemical and biological sensing and nanoimaging.

This chapter was supported by the Qatar National Research Fund (Grant No. NPRP 8-028-1-001) and the Spanish Ministry of Economy and Competitiveness (Grants No. TEC2014-53727-C2-1-R and TEC2017-86102-C2-1R).