## Abstract

We find the conditions for the existence of trapped modes in planar periodic particle arrays. Confined excitations of TE and TM symmetry are observed in symmetric environments, originating in lattice resonances that are signalled by the onset of new diffraction beams. This mechanism of mode formation is shown to be inhibited by the presence of a dielectric interface in an asymmetric configuration. Modes can still exist above a threshold finite distance from the interface. Both rigorous numerical simulation and analytical modeling are used to elucidate the origin and systematics of this unexpected difference in the behavior of trapped modes in self-standing and supported particle arrays.

©2009 Optical Society of America

## 1. Introduction

The visible and near-infrared optical properties of individual nanoparticles have recently attracted considerable interest, partly triggered by the increasing degree of control over size, composition, and morphology for structures ranging from a few nanometers to several microns, and produced by colloidal chemistry [1] and epitaxy [2] (bottom-up approach) or via lithography [3] (top-down). Metal nanoparticles are specially attractive because they are capable of sustaining plasmon resonances that are accompanied by large enhancement of the near field, which finds applications to nonlinear optics [4] and molecule sensing [5].

Nanoparticle assemblies offer additional ways of manipulating optical fields that have been exploited to realize interesting phenomena such as light waveguiding [6, 7] and narrowing of plasmon bands [8, 9, 10], relying on the existence of confined modes in periodic particle arrays. The latter have been shown to produce full reflection even at extreme dilute concentrations, provided the wavelength is close to the period [11, 12], as recently demonstrated by experiment [10]. Understanding and mastering trapped modes in particle arrays thus gives access to new ways of controlling the flow of light in confined systems and might find applications to micro lasers, sensing, and optical signal processing.

In this work, we investigate the conditions for the existence of trapped modes in planar periodic arrays of small particles. Modes of both TE and TM symmetry are observed in particle arrays embedded in a homogeneous environment. These modes are degenerate for isotropic particles, but they have different dispersion relation for anisotropic ones, which we illustrate by considering arrays of ellipsoids. The origin of the modes lies in the long-range dipole-dipole interaction that leads to lattice divergences signalled by the onset of new diffraction orders. The main finding of this work is that these divergences disappear when the particles are sitting in an asymmetric environment, for example, near the interface between two dielectric materials, so that the long-range interaction is suppressed. The exception to this rule is the perfect-conductor surface, for which TM excitations are still present in an array supported on it, although TE modes are not allowed in this case. We derive these results from both rigorous electromagnetic calculations and analytical modeling.

## 2. Analytical model

In making emphasis on understanding trapped modes in self-standing and supported particle arrays, we introduce an analytical model to describe the response of these systems. The modes of a planar periodic array formed by small particles can be obtained from a direct extension of the methods described elsewhere for dealing with planar disks [12]. Electric and magnetic dipoles (**p** and **m**, respectively) can be induced in response to an incident light plane wave. They admit the closed-form expression

where **E**
^{ext} and **H**
^{ext} are the external electric and magnetic fields; the matrix

$\alpha =\left[\begin{array}{cc}\multicolumn{1}{c}{{\alpha}^{E}}& \multicolumn{1}{c}{0}\\ \multicolumn{1}{c}{0}& \multicolumn{1}{c}{{\alpha}^{M}}\end{array}\right]$

is formed by the 3×3 electric and magnetic polarizability tensors, *α ^{E}* and

*α*, respectively; and the coefficients of

^{M}$G=\left[\begin{array}{cc}\multicolumn{1}{c}{{G}^{\mathrm{EE}}}& \multicolumn{1}{c}{{G}^{\mathrm{EM}}}\\ \multicolumn{1}{c}{{G}^{\mathrm{ME}}}& \multicolumn{1}{c}{{G}^{\mathrm{MM}}}\end{array}\right]$

are lattice sums that describe the interaction among particles. In a homogeneous environment, *G ^{EE}*=

*G*is a symmetric matrix, whereas

^{MM}*G*=-

^{EM}*G*is antisymmetric. More precisely,

^{ME}${G}_{\mathrm{ij}}^{\mathrm{EE}}=\sum _{\mathbf{R}\ne 0}{e}^{-i{\mathbf{k}}_{\parallel}\xb7\mathbf{R}}\left({k}^{2}{\delta}_{\mathrm{ij}}+{\partial}_{i}{\partial}_{j}\right)\genfrac{}{}{0.1ex}{}{{e}^{i\mathrm{kR}}}{R},\text{}$ ${G}_{\mathrm{ij}}^{\mathrm{EM}}=-ik\sum _{\mathbf{R}\ne 0}{e}^{-i{\mathbf{k}}_{\parallel}\xb7\mathbf{R}}\sum _{l}{\epsilon}_{\mathrm{ijl}}{\partial}_{l}\genfrac{}{}{0.1ex}{}{{e}^{i\mathrm{kR}}}{R},$,

where the sums run over lattice sites **R**, the indices *i, j, l* denote Cartesian directions, *k* is the light wavevector, we have used the permutation symbol *εijl*, and *δ _{ij}* is 1 if

*i*=

*j*and 0 otherwise. The parallel wavevector of the incident plane wave

**k**

_{‖}introduces an implicit exp(i

**k**

_{‖}·

**R**) dependence of

**p**and m on

**R**.

## 3. Self-standing particle arrays

For simplicity, we first focus on non-absorbing, perfectly-conducting particles possessing axial symmetry around the lattice-plane normal, so that their polarizability tensors are diagonal and have components *α*
_{‖} and *α*
_{⊥} along directions parallel and perpendicular to that plane. Then, the magnetostatic response is given by [13] *α*
^{M}
_{0⊥}=-*α*
^{E}
_{0‖}/2 and *α*
^{M}
_{0‖}=-*α*
^{E}
_{0⊥}/2 in terms of the electrostatic polarizability, which is in turn described by Gans’ theory [14] (the 0 subscript refers to the long-wavelength limit). The quantities *α*
^{E}
_{0⊥} and *α*
^{E}
_{0‖} are shown in the inset of Fig. 1 as functions of the particle aspect ratio (AR). Besides, direct application of the optical theorem leads to Im{-1/*α ^{E}*}=Im{-1/

*α*}≥2

^{M}*k*

^{3}/3, with the equal sign applying only to non-absorbing particles [15] (e.g., our perfectly-conducting ellipsoids). This allows incorporating radiative corrections in the single-particle response by approximating [16]

*α*=1/(

*α*

^{-1}

_{0}-2i

*k*

^{3}/3).

The surface-trapped modes of the array are signalled by the zeros of Δ=|*α*
^{-1}-*G*|, which correspond to solutions of Eq. (1) with finite polarization and no external field. These modes must lie in the evanescent region defined by *k*
_{‖}>*k*. Actually, they can have infinite lifetime only for **k**
_{‖} and *k* lying below the diffraction threshold [see region I in Fig. 2(b)]: the imaginary part of D cancels out in that region for non-absorbing particles because Im{*G*}=-2*k*
^{3}/3 exactly compensates Im{*α*
^{-1}} [12]. However, the polarizability scales with the cube of the particle diameter for fixed AR, and therefore, *α*
^{-1} diverges for small sizes, so that Δ=0 requires large values of *G* to match this divergence. In fact, some of the lattice sums are singular at the onset of new diffraction orders [17, 12], and more precisely, they diverge as

for *k*
_{‖}≥*k*. The rest of non-symmetry-related components of *G* are nonsingular. Here, *A* is the area of the lattice unit cell, while $\kappa =\sqrt{{k}_{\parallel}^{2}-{k}^{2}}$ characterizes the exponential decay of the modes with the distance *b* from the plane of the particles [*E*~exp(-*κb*)]. Solving Eq. (1) and using the approximate expressions of Eq. (2), we can obtain the specular reflection coefficients of the array by adding the contribution of all induced dipoles, as explained in Ref. [12]. More precisely, we can approximate the specular reflection coefficients of the array for TE and TM

polarization as

$${r}_{\mathrm{TM}}^{a}=\genfrac{}{}{0.1ex}{}{{\gamma}_{\parallel}^{M}}{{\Sigma}^{-1}-\left({\gamma}_{\perp}^{E}+{\gamma}_{\parallel}^{M}\right)}\approx \genfrac{}{}{0.1ex}{}{-{\alpha}_{0\perp}^{E}}{2{\Sigma}^{-1}-{\alpha}_{0\perp}^{E}},$$

where we have defined *γ*=1/Re{*α*
^{-1}}≈*α*
_{0}. The poles of *r ^{a}* (or equivalently, the solutions of Δ=0) give rise, after some algebra, to the mode dispersion relations contained in the first row of Table 1, subject to the condition stated in the caption of that table. The dispersion relations of both TE and TM modes have the common form

where *S*=*πD*
^{2}/4 is the projected area of the particles and D is the diameter of the ellipsoids. The dimensionless constant Γ only depends on the particle AR, as shown in Fig. 1. These modes are degenerate for spheres (AR=1), but otherwise, they exhibit very different behavior. Incidentally, the magnetic polarizability vanishes for small dielectric particles, but their modes behave similar to Fig. 1, and they are still degenerate for spheres.

In the flat-disk limit (AR=0), the TM modes disappear (Γ=0, so they are completely delocalized away from the array), whereas TE modes converge to a finite value of Γ=16/9*π*
^{3}. Interestingly, in virtue of Babinet’s principle, these are also the TM modes of the complementary thin-film drilled with a hole array [11, 12], which were experimentally resolved for millimeter waves by Ulrich and Tacke [19]. For large AR, the TM modes become strongly confined as a consequence of the divergence in *α*
^{E}
_{0⊥} when the ellipsoids evolve towards needles of constant *S*.

## 4. Supported particle arrays

Trapped modes change dramatically when the particles are sitting near a dielectric interface. For a perfectly-conducting surface, normal-electric and parallel-magnetic dipoles are doubled by their images, while parallel-electric and normal-magnetic dipoles are suppressed. Consequently, the modes obey the self-standing-array dispersion relations with the transformations *γ ^{E}*

_{⊥}→2

*γ*

^{E}_{⊥},

*γ*

^{M}_{‖}→2

*γ*

^{M}_{‖},

*γ*

^{E}_{‖}→0, and

*γ*

^{M}_{⊥}→0. As a result, TE modes disappear and TM modes are more tightly bound to the array (see the corresponding dispersion relation in the second row of Table 1).

It is important to stress that these modes are driven by lattice singularities, which ultimately produce the required divergences in G close to the onset of new diffraction orders [for example, for *k*
_{‖}≥*k* in the expression (2)]. However, for an interface separating real materials, the Fresnel coefficients satisfy *r ^{s}*

_{TM}=

*r*

^{s}_{TE}=-1 at grazing incidence, so that light emerging from each particle along a grazing direction is exactly canceled by the reflected wave, thus suppressing the interaction between distant particles along the surface, and therefore averting the divergences of

*G*and the possibility of having trapped modes. In contrast,

*r*

^{s}_{TM}=-

*r*

^{s}_{TE}=1 in a perfect-conductor surface, so that TM grazing waves are reinforced upon reflection, and consequently, TM modes can still exist, as noted above.

The asymmetric configuration is discussed in more detail in Fig. 2, which shows the disappearance of TM modes when an array of dielectric particles is placed close to an interface separating two dielectrics 1 and 2. Similar behavior is observed for TE modes (not shown). We present converged, rigorous numerical solutions of Maxwell’s equations obtained by using the Korringa-Kohn-Rostoka layer method [20, 21]. For comparison, the reflection coefficient of the array far from the interface is represented in Fig. 2(b) as a function of *k*
_{‖} and *k*
_{1} (here, *k _{j}* refers to the light momentum in medium

*j*). TM modes appear as a maximum in the reflectivity immediately below the light line (in black), as shown in the zoomed inset. The maximum is actually a divergence in the diffraction-free evanescent region I, in which the mode has infinite lifetime. This occurs when the particles are surrounded by a medium of lower permittivity, so that

*γ*>0 (and |

^{E}*γ*|≪

^{M}*γ*), and therefore, the conditions stated in the caption of Table 1 are satisfied. In contrast, there are no modes when the permittivity of the particles is lower, in which case their polarizability becomes negative (

^{E}*γ*<0).

^{E}The divergence in *r ^{a}*

_{TM}for the self-standing array eventually disappears when the distance t from the interface is sufficiently small. This is clearly illustrated by the reflection-coefficient curves represented in Fig. 2(c) as a function of

*k*

_{‖}within region I for

*ε*=6

_{p}*ε*

_{1},

*ε*

_{2}=

*ε*

_{1}/2, and

*k*

_{1}

*a*=

*π*/2, where

*a*is the lattice constant. The mode of the self-standing array is clearly visible as a divergence in

*r*

^{a}_{TM}(

*t*→∞ curve). For a finite separation between the array and the surface (

*t*=4

*a*), the maximum in

*r*

^{a}_{TM}is shifted closer to the light line (

*k*

_{‖}→

*k*). Finally, at sufficiently small separation (

*t*=2

*a*), the reflectivity does not have a sharp maximum, so that the mode no longer exists.

This behavior is summarized in Fig. 2(d), which shows the maximum reflection coefficient *r ^{a}*

_{max}(left scale) and the mode parallel-wavevector

*k*

^{a}_{‖}(right scale, solid curve) as a function of

*t*. The figure corroborates that the modes disappear at small distances, but surprisingly, this effect suddenly takes place below a certain threshold when approaching the interface from medium 1 (notice the sudden drop in

*r*

^{a}_{max}), thus resembling a first-order phase transition, with the distance

*t*acting as a phase parameter. For the particles inside medium 2 (

*t*<0), with

*ε*

_{2}<

*ε*

_{1}, the condition

*k*

^{a}_{‖}>

*k*

_{2}prevents propagation in the

*ε*

_{2}medium away from the array, while

*k*‖<

^{a}*k*

_{1}allows mode leakage into the medium 1 light cone, and therefore,

*r*

^{a}_{max}does not diverge, but simply increases with increasing -

*t*below a certain threshold.

We can seek further insight into this effect by examining the conditions for the existence of modes confined in the array-interface cavity, as given by the Fabry-Perot formula *F*(*k*
_{‖})=|*r ^{a}r^{s}* exp(-2κ

*t*)-1|=0. Using Eq. (3) for

*r*

^{a}_{TM}, together with the

*r*

^{s}_{TM}Fresnel coefficient for the planar interface reflection [20, 21], we obtain a minimum of

*F*(

*k*

_{‖}) as shown by the dashed curve in Fig. 2(d). Both themode threshold and the actual value of

*k*

^{a}_{‖}predicted by this analytical model are in excellent agreement with the full electromagnetic calculations (solid curve). It should be stressed that the minimum of

*F*is 0 for

*t*>0 (medium 1), indicating that the mode has infinite lifetime, in contrast to a finite value of

*F*for

*t*<0 (medium 2).

Noticing that *r ^{s}*≈-1 near the light line, one can work out from

*F*the following expression for the minimum distance that allows the existence of a fully trapped TM mode:

*t*

_{min}=

*A*/(4

*πk*

^{2}

*ε*

_{1}

*γ*

^{E}_{‖}), subject to the condition

*ε*

_{1}>

*ε*

_{2}, and neglecting the magnetic polarization. This expression becomes

*t*

_{min}=

*A*/(4

*πk*

^{2}

*ε*

_{1}

*γ*

^{E}_{⊥}) for TE modes.

From a different perspective, Fig. 2(e) analyzes the robustness of the modes in the presence of an interface separating twomedia near their index-matching point (*ε*
_{1}=*ε*
_{2}). Themode suddenly appears for *ε*
_{2}≥0.88*ε*
_{1}, and it stays bound up to a value *ε*
_{2}≈1.02*ε*
_{1}, for which *k ^{a}*

_{‖}matches the light momentum into the

*ε*

_{2}medium, so that the mode becomes leaky, losses strength, and its wavevector slowly goes back to the light cone of medium 1. It should be noted that the disappearance of the mode for

*ε*

_{2}<0.88

*ε*

_{1}is not due to leakage from the medium into which the particles are located (

*ε*

_{1}) into the surrounded medium (

*ε*

_{2}) of lower permittivity, but rather this disappearance is due to the fact that the mode is not well defined when it does not fit in the gap between the particles and the interface. These results are also confirmed by an analytical Fabry-Perot model similar to the one described above.

All the results discussed up to now are obtained assuming non-absorbing materials. For example, the configuration of Fig. 2 for *t*>0 could correspond to silicon particles embedded in a semi-infinite silica medium exposed to air in the non-absorbing near-infrared region of the spectrum. We have briefly analyzed the effect of absorption in the particles by adding an imaginary part (0.25*ε*
_{1}i and 0.25*ε*
_{1}i) to their dielectric function in Fig. 2(c) (dotted and dashed curves, respectively). This addition of absorption prevents the divergences observed in the solid curves of the figure for *t*=4*a* and *t*→∞, which now become resonances associated to a finite peak in the reflection coefficient. Even with such large imaginary parts in the dielectric function, the resonances are still well defined. Furthermore, the absence of a resonance for *t*=2*a* is maintained with absorption. Therefore, we conclude that the addition of a moderate amount absorption does not change our conclusions regarding the presence or absence of lattice resonances depending on the proximity of a substrate, although these resonances involve modes of finite lifetime when absorption is switched on, in contrast to the infinite lifetime predicted for non-absorbing systems.

## 5. Systems formed by real metals

An interesting scenario arises when the substrate has intrinsic modes, for instance surface-plasmon polaritons (SPPs). Figure 3 shows the reflection coefficient of an array formed by dielectric spheres embedded in a lower-index medium as a function of their distance from the interface with a silver surface [18]. It is known that lattice modes and SPPs can coexist at large separations between the array and the planar surface [12]. However, at shorter distances, SPPs hybridize with lattice modes to produce a complex evolution in the dispersion relation as a function of the distance from the planar array to the surface of a plasmon-supporting metal. The mode of the array in the absence of a substrate (*k ^{a}*

_{‖}) and the SPP (

*k*

^{SP}_{‖}) are indicated by vertical lines in Fig. 3. The twomodes interact with each other as the array approaches the metal surface, and eventually they converge to a single excitation for a distance of approximately twice the period, below which both modes seem to disappear: the intrinsic surface excitations of the substrate (e.g., surface plasmons) are modified by the presence of the array, and eventually disappear when the interaction is sufficiently strong.

Incidentally, the SPP feature cannot be observed at large *t* in Fig. 3 because we are representing the reflection coefficient for evanescent waves relative to the plane of the particle array; for large *t*, the incoming field intensity is decreasing evanescently at the planar metal surface.

It is interesting to see that a maximum of the reflection coefficient is achieved near the point where SPPs and lattice modes merge each other (*t* ~2*µ*m), that is, where the resulting mode becomes rather flat in *k*
_{‖}. This points at the presence of an intense, nearly dispersionless resonance (vanishing group velocity). It should be interesting to study whether this type of behavior is also obtained in non-absorbing systems, such as guided modes in a planar-dielectric waveguide rather than a planar-metal surface as a source of additional surface excitations.

## 6. Conclusion

In summary, we have assessed the conditions for the existence of confined collective excitations in particle arrays, either sitting in a symmetric environment or supported by a substrate. When placed close to an interface, the long-range nature of the inter-particle interaction is suppressed, thus preventing the existence of lattice resonances that are otherwise present in the symmetric configuration. Consequently, surface modes are only possible in (i) a homogeneous environment, (ii) near a perfect-conductor surface, or (iii) at finite but sufficiently large distance from a dielectric interface. For small particles, the existence of the modes is conditioned to the sign of their polarizability, which needs to be positive. Likewise, corrugated surfaces can be analyzed in a first instance as consisting of a particle array (the corrugations) near a surface (the remaining planar substrate). In fact, particle-array modes bear resemblance to trapped modes in gratings [22] and photonic crystal slabs [23, 24], although their nature differs substantially from the modes considered here, since the interaction is strongly mediated by the slab itself for narrow holes or it relies on strong inter-hole interaction for higher filling fractions of the apertures, thus giving rise to mixing of TE and TM modes (typical of 2D photonic crystals), in contrast to the modes of small-particle arrays. Larger particles should probably produce these type of photonic-crystal-specific effects. For example, surface modes in drilled perfect conductors can be described this way [12]. Additionally, we have shown that the interaction between collective particle excitations and SPPs can lead to mode annihilation at short array-metal distances. Work is currently in progress to apply these methods to other types of corrugated surfaces. We hope that our results can be used as a navigation chart into the modes of particle arrays for on-demand optical surface-mode design.

## Acknowledgments

This work has been supported by the Spanish MICINN (MAT2007-66050 and Consolider NanoLight.es) and the EU (NMP4-2006-016881-SPANS and NMP4-SL-2008-213669-ENSEMBLE).

## References and links

**1. **L. M. Liz-Marzán, “Tailoring surface plasmon through the morphology and assembly of metal nanoparticles,” Langmuir **22**, 32–41 (2006).
[CrossRef]

**2. **C. Didiot, S. Pons, B. Kierren, Y. Fagot-Revurat, and D. Malterre, “Nanopatterning the electronic properties of gold surfaces with self-organized superlattices of metallic nanostructures,” Nat. Nanotech. **2**, 617–621 (2007).
[CrossRef]

**3. **P. Ghenuche, S. Cherukulappurath, T. H. Taminiau, N. F. van Hulst, and R. Quidant, “Spectroscopic mode mapping of resonant plasmon nanoantennas,” Phys. Rev. Lett. **101**, 116,805 (2008).
[CrossRef]

**4. **M. Danckwerts and L. Novotny, “Optical frequency mixing at coupled gold nanoparticles,” Phys. Rev. Lett. **98**, 026,104 (2007).
[CrossRef]

**5. **L. Rodríguez-Lorenzo, R. A. Álvarez-Puebla, I. Pastoriza-Santos, S. Mazzucco, O. Stéphan, M. Kociak, L. M. Liz-Marzán, and F. J. García de Abajo, “Zeptomol Detection through controlled ultrasensitive surface-enhanced Raman scattering,” J. Am. Chem. Soc. **131**, 4616–4618 (2009).
[CrossRef] [PubMed]

**6. **M. Quinten, A. Leitner, J. R. Krenn, and F. R. Aussenegg, “Electromagnetic energy transport via linear chains of silver nanoparticles,” Opt. Lett. **23**, 1331–1333 (1998).
[CrossRef]

**7. **S. A. Maier, P. G. Kik, H. A. Atwater, S. Meltzer, E. Harel, B. E. Koel, and A. A. G. Requicha, “Local detection of electromagnetic energy transport below the diffraction limit in metal nanoparticle plasmon waveguides,” Nat. Mater. **2**, 229–232 (2003).
[CrossRef] [PubMed]

**8. **S. Zou, N. Janel, and G. C. Schatz, “Silver nanoparticle array structures that produce remarkably narrow plasmon lineshapes,” J. Chem. Phys. **120**, 10,871–10,875 (2004).
[CrossRef]

**9. **E. M. Hicks, S. Zou, G. C. Schatz, K. G. Spears, R. P. Van Duyne, L. Gunnarsson, T. Rindzevicius, B. Kasemo, and M. Käll, “Controlling plasmon line shapes through diffractive coupling in linear arrays of cylindrical nanoparticles fabricated by electron beam lithography,” Nano Lett. **5**, 1065–1070 (2005).
[CrossRef] [PubMed]

**10. **B. Auguie and W. L. Barnes, “Collective resonances in gold nanoparticle Arrays,” Phys. Rev. Lett. **101**, 143,902 (2008).
[CrossRef]

**11. **F. J. García de Abajo, R. Gómez-Medina, and J. J. Sáenz, “Full transmission through perfect-conductor subwavelength hole arrays,” Phys. Rev. E **72**, 016,608 (2005).
[CrossRef]

**12. **F. J. García de Abajo, “Light scattering by particle and hole arrays,” Rev. Mod. Phys. **79**, 1267–1290 (2007).
[CrossRef]

**13. **S. A. Podosenov, A. A. Sokolov, and S. V. Al’betkov, “Method for determining the electric and magnetic polarizability of arbitrarily shaped conducting bodies,” IEEE Trans. Electr. Compatibility **39**, 1–10 (1997).
[CrossRef]

**14. **R. Gans, “The shape of ultra microscopic gold particles,” Ann. Phys. (Leipzig) **37**, 881–900 (1912).

**15. **H. C. van de Hulst, *Light Scattering by Small Particles* (Dover, New York, 1981).

**16. **B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. **333**, 848–872 (1988).
[CrossRef]

**17. **J. W. S. Rayleigh, “Note on the remarkable case of diffraction spectra described by Prof. Wood,” Philos. Mag. **14**, 60–65 (1907).

**18. **We are concerned with the coefficient multiplying the zero-order reflected evanescent wave, normalized to the amplitude of an incident evanescent wave at the plane of the array.

**19. **R. Ulrich and M. Tacke, “Submillimeter waveguiding on periodic metal structure,” Appl. Phys. Lett. **22**, 251–253 (1973).
[CrossRef]

**20. **N. Stefanou, V. Yannopapas, and A. Modinos, “Heterostructures of photonic crystals: Frequency bands and transmission coefficients,” Comput. Phys. Commun. **113**, 49–77 (1998).
[CrossRef]

**21. **N. Stefanou, V. Yannopapas, and A. Modinos, “MULTEM 2: A new version of the program for transmission and band-structure calculations of photonic crystals,” Comput. Phys. Commun. **132**, 189–196 (2000).
[CrossRef]

**22. **R. Zengerle, “Light propagation in singly and doubly periodic planar waveguides,” J. Mod. Opt. **34**, 1589–1617 (1987).
[CrossRef]

**23. **S. G. Johnson, S. Fan, P. R. Villeneuve, and J. D. Joannopoulos, “Guided modes in photonic crystal slabs,” Phys. Rev. B **60**, 5751–5758 (1999).
[CrossRef]

**24. **S. G. Tikhodeev, A. L. Yablonskii, E. A. Muljarov, N. A. Gippius, and T. Ishihara, “Quasiguided modes and optical properties of photonic crystal slabs,” Phys. Rev. B **66**, 045,102 (2002).
[CrossRef]