## Abstract

We analyzed surface-wave propagation that takes place at the boundary between a semi-infinite dielectric and a multilayered metamaterial, the latter with indefinite permittivity and cut normally to the layers. Known hyperbolization of the dispersion curve is discussed within distinct spectral regimes, including the role of the surrounding material. Hybridization of surface waves enable tighter confinement near the interface in comparison with pure-TM surface-plasmon polaritons. We demonstrate that the effective-medium approach deviates severely in practical implementations. By using the finite-element method, we predict the existence of long-range oblique surface waves.

© 2013 Optical Society of America

## 1. Introduction

Artificial nanostructured materials can support electromagnetic modes that do not propagate in conventional media making them attractive for photonic devices with capabilities from nanoscale waveguiding to invisibility [1, 2]. The availability of metamaterials may lead to enhanced electromagnetic properties such as chirality, absorption, and anisotropy [3–5]. Engineered spatial dispersion is established as an essential route for diffraction management and subwavelength imaging [6, 7]. In particular, a giant birefringence also creates proper conditions for excitation of nonresonant hybrid surface waves with potential application in nanosensing [8, 9].

In its pioneering paper, Dyakonov theoretically demonstrated the existence of nondissipative surface waves at the boundary of a dielectric material and a transparent uniaxial medium [10]. However, the first experimental observation of Dyakonov surface waves came into sight more than 20 years later, most specially caused by a weak coupling with external sources [11, 12]. Indeed Dyakonov-like surface waves (DSWs) also emerge in the case that a biaxial crystal [13, 14] or a structurally chiral material [15, 16] takes the place of the uniaxial medium. The case of metal-dielectric (MD) multilayered media is specially 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 [17–19]. We also include the relative ease of nanofabrication, bulk three-dimensional response, and broadband nonresonant tunability.

For near-infrared and visible wavelengths, nanolayered MD compounds behave like plasmonic crystals enabling a simplified description of the medium by using the long-wavelength approximation, which involves an homogenization of the structured metamaterial [20–22]. Under certain conditions, the second-rank tensor denoting permittivity in the medium include elements of opposite sign, leading to extremely anisotropic metamaterials [23, 24]. This class of nanostructured media with hyperbolic dispersion are promising metamaterials with a plethora of practical applications from biosensing to fluorescence engineering [25]. In this context, Jacob *et al.* showed for the first time the existence of DSWs when considering anisotropic media with indefinite permittivity [8]. This study was devoted mainly to surface waves enabling subdiffraction imaging in magnifying superlenses [26], where hyperbolic DSWs exist at the interface of a metal and an all-dielectric birefringent metamaterial. When considering hyperbolic media, however, the authors provided only an elusive analysis of DSWs.

In this paper we retake the task and we perform a thorough analysis of DSWs taking place in semi-infinite MD lattices showing hyperbolic dispersion. In the first part of our study, our approach puts emphasis on the effective-medium approximation (EMA). Under these conditions, different regimes are found and they are thoroughly analyzed. These regimes include DSWs with non-hyperbolic dispersion. Validation of our results is carried out when put into practice using numerical simulations based on the finite-element method (FEM). The major points of our interest are nonlocal effects caused by the finite size of the layers and dissipative effects due to ohmic losses in the metals. Finally, the main conclusions are outlined.

## 2. The hyperbolic regimes of the plasmonic crystal

The system under analysis is depicted in Fig. 1. An isotropic material of dielectric constant *ε* is set in the semi-space *x* > 0. Filling the complementary space, *x* < 0, we consider a periodic bilayered structure made of two materials alternatively stacked along the *z* axis. Specifically, a transparent material of dielectric constant *ε _{d}* and slab width

*w*is followed by a metallic layer, the latter characterized by the permittivity

_{d}*ε*and the width

_{m}*w*. For simplicity, we assume that dielectric materials are nondispersive; indeed we set

_{m}*ε*= 1 and

*ε*= 2.25 in our numerical simulations. Furthermore, if a Drude metal is included, its permittivity may be written as

_{d}*ω/ω*.

_{p}Form anisotropy of this type of plasmonic devices may be modelled in a simple way by employing average estimates of the dyadic permittivity *ε̿*[27]. Under the conditions in which the EMA can be used (see details in Appendix A), the plasmonic lattice behaves like a uniaxial crystal whose optical axis is normal to the layers (the *z* axis in Fig. 1). The relative permittivity is set as *ε̿* = *ε*_{⊥}(**xx**+**yy**)+*ε*_{||}**zz**. We point out that the engineered anisotropy of the 1D lattice is modulated by the filling factor of the metal, *f*, but also by its strong dispersive character [28].

In Fig. 2(a) we represent permittivities *ε*_{||} and *ε*_{⊥} of the plasmonic crystal shown in Fig. 1 for a wide range of frequencies. In practical terms, the filling factor of the metal governs the dissipative effects in the metamaterial, thus low values of *f* are of great convenience. Indeed *f* = 1/4 in our numerical simulations. For low frequencies, Ω ≪ 1, we come near the following expressions: *ε*_{⊥} ≈ *fε _{m}* < 0 and 0 <

*ε*

_{||}≈

*ε*/(1−

_{d}*f*). Therefore, propagating TE

*modes (*

^{z}*E*= 0) cannot exist in the bulk crystal since it behaves like a metal in these circumstances. On the other hand, TM

_{z}*waves (*

^{z}*B*= 0) propagate following a spatial dispersion curve,

_{z}*k*=

_{p}*ω*/

_{p}*c*. Note that Eq. (2) denotes an hyperboloid of one sheet (see Fig. 2(b) at Ω = 0.20). Furthermore, the hyperbolic dispersion exists up to a frequency for which

*ε*

_{⊥}= 0; in our numerical example it occurs at Ω

_{1}= 0.359. For slightly higher frequencies, both

*ε*

_{||}and

*ε*

_{⊥}are positive and Eq. (2) becomes an ellipsoid of revolution; Figure 2(b) at Ω = 0.45 is associated with such a case. Since its minor semi-axis is $\mathrm{\Omega}\sqrt{{\epsilon}_{\perp}}$, the periodic multilayer simulates an anisotropic medium with positive birefringence. Raising the frequency even more,

*ε*

_{||}diverges at leading to the so-called canalization regime [3]. For the plasmonic lattice in Fig. 1 it happens at Ω

_{2}= 0.756. In general, Ω

_{1}< Ω

_{2}provided that

*f*< 1/2. Beyond Ω

_{2}, Eq. (2) turns to an hyperboloidal shape. In the range Ω

_{2}< Ω < 1, however, the dispersion curve has two sheets. Figure 2(b) illustrate this case at Ω = 0.80. Note that the upper limit of this hyperbolic band is determined by the condition

*ε*

_{||}= 0, or equivalently

*ε*= 0, occurring at the plasma frequency. To conclude,

_{m}*ε*

_{⊥}< (1 −

*f*)

*ε*in this spectral range.

_{d}## 3. Dyakonov-like surface waves

From an analytical point of view, our surface waves are localized at *x* = 0, their amplitudes decay as |*x*| → ∞, and ultimately they must satisfy Maxwell’s equations. For that purpose we follow Dyakonov by considering a modal treatment of our problem [10], also shown in Appendix B. This is a simplified procedure that is based on the characterization of the plasmonic lattice as a uniaxial crystal, enabling to establish the diffraction equation that gives the in-plane wave vector **k*** _{D}* = [0,

*k*,

_{y}*k*] of the DSW.

_{z}In the isotropic medium, we consider a superposition of TE* ^{x}* (

*E*= 0) and TM

_{x}*(*

^{x}*B*= 0) space-harmonic waves whose wave vectors have the same real components

_{x}*k*and

_{y}*k*in the plane

_{z}*x*= 0. These fields are evanescent in the isotropic medium, in direct proportion to exp(−

*κx*), being the attenuation constant

*k*. On the other side of the boundary, the ordinary wave (

_{p}*o*-wave) and extraordinary wave (

*e*-wave) existing in the effective-uniaxial medium also decay exponentially with rates given by and

**k**

*.*

_{D}In the special case of the surface wave propagation perpendicular to the optical axis (*k _{z}* = 0), Eq. (8) reveals the following solution:

*ε*

_{⊥}< 0 and

*ε*< |

*ε*

_{⊥}|, this equation has the well-known solution: which resembles the dispersion of conventional surface plasmon polaritons [29]. Indeed here we have purely TM

*polarized waves, as expected. It is worth noting that no solutions in the form of surface waves can be found from the Eq. (8) in the case of propagation parallel to the optical axis (*

^{x}*k*= 0) for hyperbolic metamaterials:

_{y}*ε*

_{⊥}

*ε*

_{||}< 0. That means that there is a threshold value of

*k*for the existence of surface waves. However, for the frequencies and the filling factors when both

_{y}*ε*

_{⊥}< 0;

*ε*

_{||}< 0, the solutions of Eq. (8) appear in the form of Bloch surface waves [22], i.e. for

*k*= 0.

_{y}## 4. DSWs in hyperbolic media

The analysis of DSWs takes into consideration effective-anisotropic metamaterials with indefinite permittivity, that is, provided that *ε*_{⊥}*ε*_{||} < 0. In this case, DSWs may be found in different regimes, which depend not only on the elements of *ε̿* characterizing the metamaterial, but also that of the surrounding isotropic material *ε*. Next we describe distinct configurations governing DSWs, first subject to a low value of the refractive index
$n=\sqrt{\epsilon}$, namely *ε* < *ε*_{||} (*ε* < *ε*_{⊥}) for low frequencies (in the neighborhood of the plasma frequency), and latter focused on a high index of refraction.

#### 4.1. Low index of refraction n

First we analyze the case when *ε* < *ε*_{||} in the range Ω < Ω_{1}, where in addition *ε*_{⊥} < 0 provided that *f* < 1/2. Note that this is satisfied if *ε* < *ε _{d}*/(1 −

*f*). In the effective-uniaxial medium,

*o*-waves are purely evanescent, and it is easy to see that

*κ*<

*κ*and also

_{o}*κ*<

_{e}*κ*. Under these circumstances, all brackets in Eq. (8) are positive provided that

_{o}*ε*

_{⊥}<

*ε*, in this case we cannot find a stationary solution of Maxwell’s equations satisfying Dyakonov’s equation (8). This happens within the spectral band Ω

_{0}< Ω < Ω

_{1}, where Note that Ω

_{0}= 0.292 in our numerical simulation. For instance, in the limiting case

*ε*= −

*ε*

_{⊥}, the unique solution of Eq. (8) is found for

*k*→ ∞ and

_{y}*k*= 0, as can be deduced straightforwardly from Eq. (10).

_{z}In Fig. 3(a) and (b) we illustrate the dispersion equation of DSWs for two different frequencies in the range 0 < Ω < Ω_{0}. In these cases, DSW dispersion curve approaches a hyperbola. Contrarily to what is shown in Fig. 3(b), we find a bandgap around *k _{z}* = 0 in (a). In general terms it occurs if Ω < 0.271, whose limiting frequency is determined by the condition

*k*= 0 are additionally constrained to the condition ${k}_{y}\ge \mathrm{\Omega}\sqrt{{\epsilon}_{\left|\right|}}$ [see also Eq. (10)], which is a necessary condition for

_{z}*κ*to exhibit real and positive values. Finally, a case similar to that shown in Fig. 3(b) was first reported in Fig. 5(b) from Ref. [8], though some discrepancies are evident.

_{e}In order to determine the asymptotes of the hyperbolic-like DSW dispersion curve, we consider the quasi-static regime (Ω → 0) since |**k*** _{D}*| =

*k*≫ Ω. Under this approximation,

_{D}*κ*=

*k*,

_{D}*κ*=

_{o}*k*, and

_{D}*κ*= Θ

_{e}*k*, where

_{D}*k*=

_{y}*k*cos

_{D}*θ*, and

*k*=

_{z}*k*sin

_{D}*θ*. Note that 0 ≤ Θ ≤ 1. By inserting all these approximations into Eq. (8), and performing the limit

*k*→ ∞, we attain the equation

_{D}*ε*+

*ε*

_{⊥}Θ = 0 straighforwardly. The latter equation indicates that hyperbolic-like solutions of Dyakonov’s equation may be found provided that

*ε*

_{⊥}< 0 and additionally

*ε*< −

*ε*

_{⊥}, occurring in the range Ω < Ω

_{0}. In this case, the asymptotes follow the equation

*k*=

_{z}*k*tan

_{y}*θ*, where

_{D}*e*-waves dispersion curve, in the

*k*-plane, have slopes given by tan2

_{y}k_{z}*θ*= −

_{e}*ε*

_{⊥}/

*ε*

_{||}. As a consequence

*θ*<

_{D}*θ*, as illustrated in Fig. 3(b), and in the limit Ω → 0 (

_{e}*ε*

_{⊥}→ −∞) we obtain

*θ*.

_{D}→θ_{e}Moving into the high-frequency band Ω_{2} < Ω < 1, we now find that *ε*_{||} < 0 < *ε*_{⊥}. The plot shown in Fig. 3(c) corresponds to this case. In a similar way found in Fig. 3(a) and (b), note the relevant proximity of DSW dispersion curve to *κ _{e}* = 0. Opposedly it crosses the

*e*-wave hyperbolic curve at two different points, where solutions of Dyakonov’s equation begin and end respectively. In comparison, the angular range of DSWs turns to be significantly low. Apparently the

*z*-component of

**k**

*tends to approach $\mathrm{\Omega}\sqrt{{\epsilon}_{\perp}}$ caused by the simultaneous dominance of*

_{D}*o*- and

*e*-waves. In general, an slight increase of the refractive index in the isotropic medium pushes the wave vector

**k**

*to higher values, leading to an enormous shortening in the dispersion curve of the surface waves. As a consequence, high-*

_{D}*n*materials give rise to adverse conditions for the excitation of DSWs in the neighborhood of the plasma frequency.

Figure 4 shows the magnetic field for the points *A*, *B*, and *C*, all highlighted in Fig. 3. Also we represent the *z*-component of the field **B** that is associated with the point *SP* appearing in Fig. 3(b), and that corresponds to a surface plasmon (*B _{x}* = 0). The wave field is tightly confined near the surface

*x*= 0, in a few units of 1/

*k*, for the cases

_{p}*A*and

*B*. Such a wave localization is even stronger than the confinement of the surface plasmon appearing at Ω = 0.28 (at

*k*= 0.625). This is caused by the large in-plane wavenumber of the DSW, being

_{y}*k*= 1.44 and 1.03 for the points

_{D}*A*and

*B*, respectively. Exceptionally, the lowest confinement is produced at Ω = 0.85 when making the choice

*C*, in spite of considering a DSW with large wave vector

**k**

*= [0, 0.2, 1.07]. In this case, the interplay of slowly-decaying*

_{D}*o*- and

*e*-waves counts against localization of the surface wave.

#### 4.2. High index of refraction n

Next we analyze the case when *ε*_{||} < *ε* in the spectral domain Ω < Ω_{1}. Therefore, the curve *κ* = 0 characterizing the isotropic medium crosses the TM* ^{z}* dispersion curve

*κ*= 0 of the uniaxial metamaterial. Note that such a curve crossing is mandatory by considering materials with all-positive permittivities, alike the pioneer paper by Dyakonov [10, 13]. Our current case discloses some similarities with DSWs analyzed in Sec. 4.1. For instance, no solutions of Dyakonov’s equation (8) are found if additionally −

_{e}*ε*

_{⊥}<

*ε*, occurring at Ω

_{0}< Ω < Ω

_{1}. Nevertheless, there are some distinct features which are worthy to mention.

In our numerical simulations we used a dielectric material with *ε* = 10, leading to Ω_{0} = 0.145. In Fig. 5(a) we illustrate the dispersion equation of DSWs at Ω = 0.1. Here DSW curve also approaches a hyperbola. In our instance, however, on-axis bandgaps are not found even at lower frequencies, and Eq. (8) provides solutions for every real value of *k _{z}*. Figs. 5(b) and 5(c) show the profile of the magnetic field along the

*x*axis for two different points (

*D*and

*SP*) of the dispersion curve. Once again, hybrid surface waves (case

*D*) exhibit a tighter confinement near the boundary

*x*= 0 than that offered by the solution of Eq. (8) and that is attributed to surface plasmons with pure TM

*polarization (case*

^{x}*SP*).

Finally by assuming *ε*_{⊥} < *ε* at a given frequency of the spectral window Ω_{2} < Ω < 1, we have not found solutions of Eq. (8). As discussed in Sec. 4.1, high values of *ε* goes in prejudice of the appearance of hybrid surface waves.

## 5. Validity of the effective-medium approximation

In the previous sections we have utilized EMA to represent MD multilayered metamaterial as a uniaxial plasmonic crystal. However, it was shown recently that EMA does not describe properly neither nonlocalities (even in the case of negligible losses) [31], nor it takes into account correctly metallic losses [32]. In fact, the MD layered nanostructure demonstrates strong optical nonlocality due to excitation of surface plasmon polaritons, depending on the thicknesses of the layers. Dispersion and diffraction properties of the periodic nanolayered metamaterial can be dramatically affected by strong coupling among surface plasmon polaritons excited at metal-dielectric interfaces. In many cases, EMA produce results that are so different from the exact transfer-matrix method [27] that cannot be taken into account as small corrections [30]. In particular, the existence of additional extraordinary wave has been revealed [31] which enable double refraction of TM* ^{z}*-polarized beam into negatively and positively refracted parts, within the favourable frequency range. Here, however, we deal with hyperbolic metamaterials that require frequency ranges out of the double refraction range. Confining ourselves to the wave fields of TM

*polarization which approach*

^{z}*e*-wave in the regime of validity of EMA, the exact map of wave vectors

**k**

*characterizing Bloch waves is calculated via:*

_{D}*φ*=

_{q}*k*, and ${k}_{x}^{2}+{k}_{y}^{2}+{k}_{qz}^{2}={\epsilon}_{q}{k}_{0}^{2}$ is the dispersion equation for bulk waves inside the dielectric (

_{qz}w_{q}*q*=

*d*) and Drude metal (

*q*=

*m*).

Equation (16) is graphically represented in Fig. 6 for the MD crystal displayed in Fig. 1 at *x* < 0, considering different widths *w _{m}* of the metal but maintaining its filling factor,

*f*= 0.25; ohmic losses in the metal are neglected once again. We use the same configurations appearing in Figs. 3(a)–3(c) and Fig. 5(a). In the numerical simulations, we set

*k*= 0 in Eq. (16) accordingly. We observe that the EMA is extremely accurate for

_{x}*w*= 0.1 (in units of 1/

_{m}*k*) in the region of interest. However, deviations among the contours are evident for wider metallic layers. Apparently, Eq. (16) is in good agreement with the EMA in the vicinity of

_{p}*k*= 0, but rising

_{z}*w*makes the dispersion curve to sheer in direction to the

_{m}*z*axis. As a general rule, given a value of the on-axis frequency

*k*, Eq. (16) yields lower values of

_{z}*k*. As shown in the next section, this is also the cause of a spectral shift of

_{y}**k**

*along the same direction.*

_{D}The Bloch waves involved in the formation of the DSW and caused by nonlocality in the structured medium have an impact on the wave vector **k*** _{D}* at the interface

*x*= 0, as seen above, but also put additional conditions on this boundary. Since the introduction of an effective permittivity requires some kind of field averaging normally to the metal-dielectric layers, the excitation of evanescent fields in the isotropic medium would be fundamentally governed by the value of

**k**

*which determines the attenuation constant*

_{D}*κ*given in Eq. (5). However, spatial dispersion also leads to strong field oscillations across the system [30–32]. This means that the conventional boundary conditions imposed by the equation det(

*M̿*) = 0 are not valid anymore. Such a strong variation of the field is set on the scale of a single layer. Consequently, evanescent fields with spatial frequencies much higher than

*k*will participate vigorously in the isotropic medium extremely near the boundary. As we will see in the FEM simulations appearing in the next section, predominance of these high-frequency components of the field lies crucially by the edge of the metallic layers adjoining the isotropic medium.

_{D}## 6. Analysis of a practical case

Dissipation in metallic elements is a relevant issue of plasmonic devices. Taking into account ohmic losses, permittivities in Eqs. (1)–(18) become complex valued. Consequently Dyakonov’s equation (8) provides complex values of the wave vector **k*** _{D}*. This procedure has been discussed recently by Vuković et al in [18]. Nevertheless, practical implementations are out of the long-wavelength limit, thus we also direct our efforts toward nonlocal effects. Numerical techniques to solve Maxwell’s equations seem to be convenient tools in order to provide a conclusive characterization of DSWs at the boundary of realistic metallodielectric lattices.

In order to tackle this problem, we evaluate numerically the value of *k _{y}* for a given Bloch wavenumber

*k*. Since the imaginary part of

_{z}*ε*is not neglected anymore,

_{m}*k*becomes complex. This means that the DSW cannot propagate indefinitely; Im(

_{y}*k*) denotes the attenuation factor of the surface wave along the metallic-film edges. In our numerical simulations we consider a dissipative DSW propagating on the side of a Ag-PMMA lattice at a wavelength of

_{x}*λ*= 560 nm (normalized frequency Ω = 0.28), where the surrounding isotropic medium is air. Accordingly

*ε*= 1,

*ε*= 2.25, and

_{d}*ε*= −11.7 + 0.83

_{m}*i*(being Ω

*= 12.0 rad/fs) [33]. Bearing in mind a practical setting in the plasmonic lattice with current nanotechnology, we apply*

_{p}*w*= 9 nm and

_{m}*f*= 0.25.

Our computational approach lies on the FEM by means of COMSOL Multiphysics software. Thus given *k _{z}* = 0.25, which is associated with the point

*B*in Fig. 3(b), we finally estimate the complex propagation constant:

*k*= 0.70 + 0.06

_{y}*i*. We point out that Eq. (8) predicts a value

*k*= 0.85 + 0.24

_{y}*i*; in addition we obtain

*k*= 1.00 by neglecting losses in Dyakonov’s equation, as shown in Fig. 3(b). Our numerical experiment proves a “red shift” in the propagation constant caused by nonlocal effects. Furthermore, Im(

_{y}*k*) decreases sharply (roughly by a ratio 1/4) in comparison with EMA estimates. This major result enable DSWs propagating along distances significantly longer than those predicted by the long-wavelength approach.

_{x}Figure 7 shows the magnetic field **B** of the DSW in the *xz*-plane. The calculated pattern in one cell reveals the effects of retardation clearly. Along the *z* axis, an abrupt variation of **B** is evident inside the nanostructured material, in contrast with assumptions involving the effective-medium approach. The wave field cannot penetrate in the metal completely, and it is confined not only on the silver-air interface but also in the Ag-PMMA boundaries near *x* = 0. Indeed, from FEM simulations, the ratio of max|*B _{x}*| over max|

*B*| yields 0.80, considerably higher that its value on the basis of the EMA (equal to 0.37). This proves a field enhancement on the walls of the metallic films and inside the dielectric nanolayers, minimizing dissipative effects in the lossy metamaterial. Finally, the distributed field along the

_{z}*x*axis is analogous in all cases, also by comparing with EMA-based results.

## 7. Conclusions

We showed that excitation of DSWs at the boundary of an isotropic dielectric and a hyperbolic metamaterial enable distinct regimes of propagation. At low frequencies, DSWs exhibit the well-known hyperbolic-like dispersion [8]. Exclusively, this occurs provided that 0 < *ε* < −*ε*_{⊥}. Importantly its vertex is placed on the *k _{y}* axis, coinciding with a non-hybrid TM

*surface wave. If ${\epsilon}_{\left|\right|}^{-1}-{\epsilon}_{\perp}^{-1}<{\epsilon}^{-1}$ is satisfied, a bandgap around*

^{x}*k*= 0 arises in the DSW dispersion curve. However, these bandgaps disabling negative refraction would disappear by using an isotropic material with a sufficiently-high index of refraction,

_{z}*n*.

Conversely, severe restrictions are found at high frequencies where *ε*_{||} < 0 < *ε*_{⊥} (two-sheet hyperboloidal dispersion). DSWs dispersion curve crosses the *e*-wave hyperbolic curve at two different points, where solutions of Dyakonov’s equation begin and end. In comparison with the previous case, the spectral range of DSWs turns to be significantly low. Furthermore, high-*n* materials give rise to adverse conditions for the excitation of DSWs.

We point out that the wave field of DSWs keeps tightly confined near the surface *x* = 0. Such a wave localization is even stronger than the confinement of the (TM* ^{x}*) surface plasmon appearing at

*k*= 0.

_{z}In the last part of this paper we analyzed retardation effects and dissipation effects caused by the finite width and ohmic losses of the metallic nanolayers, respectively. As a practical implementation we considered a dissipative DSW propagating at the boundary of an Ag-PMMA lattice and air. Our FEM-based numerical experiment proves a “red shift” in the propagation constant caused by nonlocal effects. Furthermore, the field is enhanced on the walls of the metallic films and inside the dielectric nanolayers, minimising dissipative effects in the lossy metamaterial. As a consequence, DSWs propagate along distances significantly longer than those predicted by the long-wavelength approach.

Finally, let us point out that the properties of the resulting bound states change rapidly with the refractive index of the surrounding medium, suggesting potential applications in chemical and biological sensing and nanoimaging.

## Appendix A: Effective-medium approximation

The EMA involves a material homogenization, a procedure that requires the metallic elements to had 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, *L _{D}*[22, 30, 34, 35]. In the Drude approach,

*L*may be roughly estimated by the inverse of

_{D}*k*, that is deeply subwavelength in practice. In reference to this point, recent development of microfabrication technology makes it possible to create such subwavelength structures [36, 37]. Under these conditions, the plasmonic lattice behaves like a uniaxial crystal characterized by a relative permittivity tensor

_{p}*ε̿*=

*ε*

_{⊥}(

**xx**+

**yy**) +

*ε*

_{||}

**zz**, where

## Appendix B: Derivation of the Dyakonov’s equation

This appendix shows the procedure to derive the dispersion equation for the wave vector **k*** _{D}* of the DSW. We start with the analytical expressions of the electric and magnetic fields involved in the boundary problem. In the isotropic medium, the bimodal electric field may be written in the complete form as

**E**(

**r**

*,t*) =

**E**(

*x*)

*f*(

*y,z,t*), where

*x*> 0. Here time- and space-coordinates are normalized to the inverse of

*ω*and

_{p}*k*, respectively. In addition

_{p}*A*and

_{TE}*A*are complex-valued amplitudes, and the vectors Note that the polarization vector

_{TM}**a**

_{1}lies on the plane

*x*= 0 where the DSW propagates, but it is perpendicular to its direction of propagation,

**k**

*. For the magnetic field we apply the Faraday’s law of induction,*

_{D}*i*Ω

*c*

**B**= ∇ ×

**E**. In this case

**B**(

**r**

*,t*) =

**B**(

*x*)

*f*(

*y,z,t*), where

*A*and

_{TE}*A*represent amplitudes of the TE

_{TM}*and TM*

^{x}*modes, respectively.*

^{x}On the other side of the boundary, the dependence along the *x* direction of the multimodal electric field in *x* < 0 may be written as

*A*and

_{o}*A*stand for the amplitudes of the

_{e}*o*-wave and the

*e*-wave, respectively. Also we define the vectors

**b**

*are perpendicular to the optical axis. The part of the magnetic field which provides its variation normally to the isotropic-uniaxial interface is set now as*

_{o,e}*o*-wave corresponds to a TE

*mode and the*

^{z}*e*-wave is associated with a TM

*mode.*

^{z}Next we apply the standard electromagnetic-field boundary conditions at *x* = 0, that is, continuity of the *y*- and *z*-components of the fields **E** and **B** at the planar interface. This problem may be set in matrix form as

**A**= [

*A*,

_{TE}*A*,

_{TM}*A*,

_{o}*A*] includes the amplitudes of the four modes that the DSW integrates. In addition, the 4 × 4 matrix

_{e}*M̿*vanishes. In this way, Dyakonov derived equation (8), which provides the spectral map of wave vectors

**k**

*.*

_{D}## Acknowledgments

This research was funded by the Spanish Ministry of Economy and Competitiveness under the project TEC2011-29120-C05-01, by the Serbian Ministry of Education and Science under the project III 45016, and by the Qatar National Research Fund under the project NPRP 09-462-1-074. CJZR gratefully acknowledges a financial support from the Generalitat Valenciana (grant BEST/2012/060).

## References and links

**1. **W. Cai, U. K. Chettiar, A. V. Kildishev, and V. M. Shalaev, “Optical cloaking with metamaterials,” Nat. Photon. **1**, 224–227 (2007) [CrossRef] .

**2. **S. Lal, S. Link, and N. J. Halas, “Nano-optics from sensing to waveguiding,” Nat. Photon. **1**, 641–648 (2007) [CrossRef] .

**3. **P. A. Belov and Y. Hao, “Subwavelength imaging at optical frequencies using a transmission device formed by a periodic layered metal-dielectric structure operating in the canalization regime,” Phys. Rev. B **73**, 113110 (2006) [CrossRef] .

**4. **E. Plum, V. A. Fedotov, A. S. Schwanecke, N. I. Zheludev, and Y. Chen, “Giant optical gyrotropy due to electromagnetic coupling,” Appl. Phys. Lett. **90**, 223113 (2007) [CrossRef] .

**5. **J. Hao, L. Zhou, and M. Qiu, “Nearly total absorption of light and heat generation by plasmonic metamaterials,” Phys. Rev. B **83**, 165107 (2011) [CrossRef] .

**6. **M. Conforti, M. Guasoni, and C. D. Angelis, “Subwavelength diffraction management,” Opt. Lett. **33**, 2662–2664 (2008) [PubMed] .

**7. **C. J. Zapata-Rodríguez, D. Pastor, M. T. Caballero, and J. J. Miret, “Diffraction-managed superlensing using plasmonic lattices,” Opt. Commun. **285**, 3358–3362 (2012) [CrossRef] .

**8. **Z. Jacob and E. E. Narimanov, “Optical hyperspace for plasmons: Dyakonov states in metamaterials,” Appl. Phys. Lett. **93**, 221109 (2008) [CrossRef] .

**9. **C. J. Zapata-Rodríguez, J. J. Miret, J. A. Sorni, and S. Vuković, “Propagation of dyakonon wave-packets at the boundary of metallodielectric lattices,” IEEE J. Sel. Top. Quant. Electron. **19**, 4601408 (2013) [CrossRef] .

**10. **M. I. D’yakonov, “New type of electromagnetic wave propagating at an interface,” Sov. Phys. JETP **67**, 714–716 (1988).

**11. **O. Takayama, L.-C. Crasovan, S. K. Johansen, D. Mihalache, D. Artigas, and L. Torner, “Dyakonov surface waves: A review,” Electromagnetics **28**, 126–145 (2008) [CrossRef] .

**12. **O. Takayama, L. Crasovan, D. Artigas, and L. Torner, “Observation of Dyakonov surface waves,” Phys. Rev. Lett. **102**, 043903 (2009) [CrossRef] [PubMed] .

**13. **D. B. Walker, E. N. Glytsis, and T. K. Gaylord, “Surface mode at isotropic-uniaxial and isotropic-biaxial interfaces,” J. Opt. Soc. Am. A **15**, 248–260 (1998) [CrossRef] .

**14. **M. Liscidini and J. E. Sipe, “Quasiguided surface plasmon excitations in anisotropic materials,” Phys. Rev. B **81**, 115335 (2010) [CrossRef] .

**15. **J. Gao, A. Lakhtakia, J. A. Polo Jr., and M. Lei, “Dyakonov-Tamm wave guided by a twist defect in a structurally chiral material,” J. Opt. Soc. Am. A **26**, 1615–1621 (2009) [CrossRef] .

**16. **J. Gao, A. Lakhtakia, and M. Lei, “Dyakonov-Tamm waves guided by the interface between two structurally chiral materials that differ only in handedness,” Phys. Rev. A **81**, 013801 (2010) [CrossRef] .

**17. **O. Takayama, D. Artigas, and L. Torner, “Practical dyakonons,” Opt. Lett. **37**, 4311–4313 (2012) [CrossRef] [PubMed] .

**18. **S. M. Vuković, J. J. Miret, C. J. Zapata-Rodríguez, and Z. Jaks̆ić, “Oblique surface waves at an interface of metal-dielectric superlattice and isotropic dielectric,” Phys. Scripta **T149**, 014041 (2012) [CrossRef] .

**19. **J. J. Miret, C. J. Zapata-Rodríguez, Z. Jaks̆ić, S. M. Vuković, and M. R. Belić, “Substantial enlargement of angular existence range for Dyakonov-like surface waves at semi-infinite metal-dielectric superlattice,” J. Nanophoton. **6**, 063525 (2012) [CrossRef] .

**20. **S. M. Rytov, “Electromagnetic properties of layered media,” Sov. Phys. JETP **2**, 466 (1956).

**21. **A. Yariv and P. Yeh, “Electromagnetic propagation in periodic stratified media. II. Birefringence, phase matching, and x-ray lasers,” J. Opt. Soc. Am. **67**, 438–448 (1977) [CrossRef] .

**22. **S. M. Vukovic, I. V. Shadrivov, and Y. S. Kivshar, “Surface Bloch waves in metamaterial and metal-dielectric superlattices,” Appl. Phys. Lett **95**, 041902 (2009) [CrossRef] .

**23. **D. R. Smith and D. Schurig, “Electromagnetic wave propagation in media with indefinite permittivity and permeability tensors,” Phys. Rev. Lett. **90**, 077405 (2003) [CrossRef] [PubMed] .

**24. **I. I. Smolyaninov, E. Hwang, and E. Narimanov, “Hyperbolic metamaterial interfaces: Hawking radiation from Rindler horizons and spacetime signature transitions,” Phys. Rev. B **85**, 235122 (2012) [CrossRef] .

**25. **Y. Guo, W. Newman, C. L. Cortes, and Z. Jacob, “Applications of hyperbolic metamaterial substrates,” Advances in OptoElectronics **2012**, ID 452502 (2012) [CrossRef] .

**26. **I. I. Smolyaninov, Y.-J. Hung, and C. C. Davis, “Magnifying superlens in the visible frequency range,” Science **315**, 1699–1701 (2007) [CrossRef] [PubMed] .

**27. **P. Yeh, *Optical Waves in Layered Media* (Wiley, 1988).

**28. **B. Wood, J. B. Pendry, and D. P. Tsai, “Directed subwavelength imaging using a layered metal-dielectric system,” Phys. Rev. B **74**, 115116 (2006) [CrossRef] .

**29. **S. A. Maier, *Plasmonics: Fundamentals and Applications* (Springer, 2007).

**30. **J. Elser, V. A. Podolskiy, I. Salakhutdinov, and I. Avrutsky, “Nonlocal effects in effective-medium response of nanolayered metamaterials,” Appl. Phys. Lett. **90**, 191109 (2007) [CrossRef] .

**31. **A. A. Orlov, P. M. Voroshilov, P. A. Belov, and Y. S. Kivshar, “Engineered optical nonlocality in nanostructured metamaterials,” Phys. Rev. B **84**, 045424 (2011) [CrossRef] .

**32. **A. Orlov, I. Iorsh, P. Belov, and Y. Kivshar, “Complex band structure of nanostructured metal-dielectric metamaterials,” Opt. Express **21**, 1593–1598 (2013) [CrossRef] [PubMed] .

**33. **E. D. Palik and G. Ghosh, *The Electronic Handbook of Optical Constants of Solids* (Academic, 1999).

**34. **E. Popov and S. Enoch, “Mystery of the double limit in homogenization of finitely or perfectly conducting periodic structures,” Opt. Lett. **32**, 3441–3443 (2007) [CrossRef] [PubMed] .

**35. **A. V. Chebykin, A. A. Orlov, A. V. Vozianova, S. I. Maslovski, Y. S. Kivshar, and P. A. Belov, “Nonlocal effective medium model for multilayered metal-dielectric metamaterials,” Phys. Rev. B **84**, 115438 (2011) [CrossRef] .

**36. **P. Chaturvedi, W. Wu, V. J. Logeeswaran, Z. Yu, M. S. Islam, S. Y. Wang, R. S. Williams, and N. X. Fang, “A smooth optical superlens,” Appl. Phys. Lett. **96**, 043102 (2010) [CrossRef] .

**37. **H. N. S. Krishnamoorthy, Z. Jacob, E. Narimanov, I. Kretzschmar, and V. M. Menon, “Topological Transitions in Metamaterials,” Science **336**, 205–209 (2012) [CrossRef] [PubMed] .