OSA's Digital Library

Optics Express

Optics Express

  • Editor: C. Martijn de Sterke
  • Vol. 15, Iss. 15 — Jul. 23, 2007
  • pp: 9681–9691

Optics InfoBase > Optics Express > Volume 15 > Issue 15 > Page 9681

« Show journal navigation

The complex Bloch bands of a 2D plasmonic crystal displaying isotropic negative refraction

Marcelo Davanco, Yaroslav Urzhumov, and Gennady Shvets  »View Author Affiliations


Optics Express, Vol. 15, Issue 15, pp. 9681-9691 (2007)
http://dx.doi.org/10.1364/OE.15.009681


View Full Text Article

Acrobat PDF (769 KB)


Advanced Search

Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

The propagation characteristics of a subwavelength plasmonic crystal are studied based on its complex Bloch band structure. Photonic crystal bands are generated with an alternative 2D Finite Element Method formulation in which the Bloch wave problem is reduced to a quadratic eigenvalue system for the Bloch wavevector amplitude k. This method constitutes an efficient and convenient alternative to nonlinear search methods normally employed in the calculation of photonic bands when dispersive materials are involved. The method yields complex wavevector Bloch modes that determine the wave-scattering characteristics of finite crystals. This is evidenced in a comparison between the band structure of the square-lattice plasmonic crystal and scattering transfer-functions from a corresponding finite crystal slab. We report on a wave interference effect that leads to transmission resonances similar to Fano resonances, as well as on the isotropy of the crystal’s negative index band. Our results indicate that effective propagation constants obtained from scattering simulations may not always be directly related to individual crystal Bloch bands.

© 2007 Optical Society of America

1. Introduction

The study of the optical properties of Photonic Crystals (PCs) has relied on the generation of Bloch-mode photonic band structures. Information such as the existence of forbidden band gaps, phase and group velocity, group velocity dispersion, photonic densities of states, propagation anisotropy, etc. may be conveniently obtained by inspection or manipulation of the computed propagation bands [1

J. D. Joannopoulos, R. D. Meade, and J. N. Winn, Photonic Crystals: Molding the Flow of Light . (Princteon University Press, 1995).

, 2

K. Sakoda, Optical Properties of Photonic Crystal, ser. Optical Sciences. (New York: Springer, 2001).

]. A commonly adopted approach for the calculation of PC band structures using frequency-domain methods consists in obtaining the allowed propagation frequencies ω for specific Bloch k-vectors as eigenvalues of an electromagnetic wave equation defined in a unit-cell of the periodic medium. If the PC incorporates only dispersionless materials (i.e., the dielectric constant ε(ω) constant for all frequencies), the problem is reduced to a linear generalized eigenvalue system of the form Akxk =ω 2 Bkxk [1

J. D. Joannopoulos, R. D. Meade, and J. N. Winn, Photonic Crystals: Molding the Flow of Light . (Princteon University Press, 1995).

–4

B. P. Hiett, B. D. H., S. J. Cox, J. M. Generowicz, M. Molinari, and K. S. Thomas, “Aplication of finite element methods to photonic crystal modelling,” IEE Proc - Sci. Meas. Technol. 149, 293–296 (2002). [CrossRef]

].

An important class of metamaterials composed of a periodic arrangement of metallic inclusions in a dielectric host, or Subwavelength Plasmonic Crystals (SPCs) [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

–7

C. Enkrich, M. Wegener, S. Linden, S. Burger, L. Zschiedrich, F. Schmidt, J. F. Zhou, T. Koschny, and C. M. Soukoulis, “Magnetic metamaterials at telecommunication and visible frequencies,” Phys. Rev. Lett. 95, 203901 (2005), doi:10.1103/PhysRevLett.95.203901. [CrossRef] [PubMed]

], are particularly useful in the UV and optical frequency ranges, where the permittivity of the constitutive metal inclusions displays a plasma-like dispersion. In [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

], such structures were shown to display a negative refractive index propagation band in a subwavelength regime (~λ/10), and simulations in fact predicted the possibility of subwavelength-resolution focusing. In this case, the full variation of the metallic permittivities must be taken into account, which results in a nonliear equation of the form Ak (ω)xk =ω 2 Bk (ω)xk , or, more generally, Tk (ω)xk =0. Solutions to this type of problem are generally obtained through iterative algorithms -normally extensions of algorithms for linear problems [8

A. Ruhe, “Algorithms for the nonlinear eigenvalue problem,” SIAM J. Numer. Anal. 10, 674–689 (1973). [CrossRef]

], or developments from Newton’s Method [9

A. Spence and C. Poulton, “Photonic band structure calculations using nonlinear eigenvalue techniques,” J. Comput. Phys. 204, 65–81 (2005). [CrossRef]

]- in which successive approximations of ω are computed within a chosen convergence domain. Successful and efficient convergence usually requires a good initial guess; in addition matrix T(ω) and, depending on the method, its derivatives, must be recalculated at several frequencies within the search domain, causing the solution process to be time-consuming and computationally intensive. Hence, for this reason, band structure diagrams are not abundant in the metamaterials and plasmonic crystals literature.

Here, we revisit the Finite Element Method (FEM) formulation for the calculation of 2D photonic crystal bands developed by Hiett et al. [4

B. P. Hiett, B. D. H., S. J. Cox, J. M. Generowicz, M. Molinari, and K. S. Thomas, “Aplication of finite element methods to photonic crystal modelling,” IEE Proc - Sci. Meas. Technol. 149, 293–296 (2002). [CrossRef]

]. This method is shown to yield a quadratic eigenvalue equation in the Bloch k-vector magnitude; frequency in this case is a parameter, so that material dispersion is readily taken into account. While the quadratic eigenvalue equation is nonlinear, it is a more tractable problem than the general nonlinear case above, and can yield solutions more efficiently. In addition, this formulation inherently yields bands of purely imaginary- and complex-wavevector Bloch modes, which may be particularly hard to obtain with nonlinear search routines. Imaginary k modes are generally available at the photonic bandgap frequencies and play an important role in representing the evanescent field inside a finite or semi-infinite PC slab under external excitation. This has been suggested in [10

E. Istrate, A. A. Green, and E. H. Sargent, “Behavior of light at photonic crystal interfaces,” Phys. Rev. B 71, 195122 (2005). [CrossRef]

], in which an equivalent formulation to that presented here, based on the plane-wave expansion method [3

S. G. Johnson and J. D. Joannopoulos, “Block-iterative frequency-domain methods for maxwell’s equations in a planewave basis,” Opt. Express 8, 173–190 (2001), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-8-3-173. [CrossRef] [PubMed]

], was used to calculate imaginary k bandgap modes in PCs composed of nondispersive, purely dielectric materials. On the other hand, modes with complex, rather than purely imaginary, k-vectors may also exist within, and even outside, bandgap regions that have considerable influence in the transmission and reflection transfer functions of finite plasmonic crystals.

In Sec. 2, we derive the FEM formulation for the plasmonic crystal eigenvalue problem following similar steps as [4

B. P. Hiett, B. D. H., S. J. Cox, J. M. Generowicz, M. Molinari, and K. S. Thomas, “Aplication of finite element methods to photonic crystal modelling,” IEE Proc - Sci. Meas. Technol. 149, 293–296 (2002). [CrossRef]

]. This formulation is implemented in a simple way with the COMSOL Finite Elements package [11

[Online]. Available: http://www.comsol.com

] and used to produce the photonic bands presented in Sec. 3, for a square-lattice SPC displaying a negative refraction band [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

,6

G. Shvets and Y. Urzhumov, “Electric and magnetic properties of sub-wavelength plasmonic crystals,” J. Opt. A: Pure Appl. Opt. 7, S23–S31 (2005). [CrossRef]

]. In these references, the negative refraction band was calculated only along one of the high-symmetry crystal directions; in the present article, the photonic bands are obtained over the entire irreducible first Brillouin zone, such that the isotropy of the negative refraction band near k=0 is revealed. A comparison between scattering data from a finite SPC and its band diagram is also performed, in which the role of the multiple existing Bloch modes (including those with imaginary and complex wavevectors) is determined. In particular, we report on the existence of resonances in transmission through a finite crystal slab that stem from the coexistence of two SPC bands in the same frequency range, being similar in origin and character to Fano resonances [12

U. Fano, “Effects of configuration interaction on intensities and phase shifts,” Phys. Rev. 124, 1866–1878 (1961). [CrossRef]

]. In Sec. 4 we present conclusions.

2. Finite-element method

We consider a periodic 2D medium with electromagnetic waves propagating in the xy plane. Solutions toMaxwell’s equations are classified as either TE (H=Hz ẑ, where H is the magnetic field) or TM (E=Ez ẑ, where E is the electric field), with the respective z-components obeying the wave equation

∇·p ∇ϕ+ ω2 c2qϕ=0.
(1)

For TE waves, ϕ=Hz , p=1/ε, q=1 and ε is the electric permittivity; for TM waves, ϕ=Hz , p=1, q=ε. In Eq.(1), ω/c is the vacuum angular wavenumber.

The medium is periodic, and hence ε(r+T)=ε(r), where T is a lattice translation vector. From Bloch’s theorem we consider solutions ϕ(r)=u(r)exp(-i k·r), with u(r+T)=u(r) and k in the first Brillouin zone. Along with Eq.(1), this leads to:

∇· (p∇u) −i∇· (pku) −ik·p∇ u− k2pu= − ω2 c2qu
(2)

Following the Galerkin procedure to reduce Eq.(1) to a weighted-residuals expression [13

J. Jin, The Finite Element Method in Electromagnetics , (2nd ed. Wiley, 2002).

], we multiply that equation on both sides by a weight function w(r), and note that

∇·(p∇u) = ∇· [p (∇u)w] −p∇ u·∇w
(3)

and

∇· (pku)w=∇·[pkuw] −k · (pu ∇w)
(4)

to arrive at

∇· [p (∇u)w] −p∇ u·∇ w−i ∇· (pkuw)+

ik· {pu∇ w− [p (∇u)w]} − k2 puw=

− ω2 c2quw
(5)

This expression is integrated over a closed domain Γ with boundary δ Γ that covers one unit cell of the periodic medium (the integration domain may actually extend over an integer number of unit cells). Making use of the Divergence theorem, the following weak-form expression is obtained:

∫Γ p∇u· ∇wdΓ+ k2 ∫ΓpuwdΓ+

−ik { ∫Γpu k̂ ·∇wd Γ− ∫Γp k̂ · (∇u)wdΓ}=

= ω2 c2 ∫Γ quwdΓ+ ∮ δΓpw ( ∇u−iku) · n̂d δΓ
(6)

In Eq.(6), n̂ is the outward normal unit vector to the boundary δ Γ. Since u is periodic, the line integral vanishes. The resulting integrodifferential equation may be transformed into matrix format by following the usual FEM discretization procedure [4

B. P. Hiett, B. D. H., S. J. Cox, J. M. Generowicz, M. Molinari, and K. S. Thomas, “Aplication of finite element methods to photonic crystal modelling,” IEE Proc - Sci. Meas. Technol. 149, 293–296 (2002). [CrossRef]

, 13

J. Jin, The Finite Element Method in Electromagnetics , (2nd ed. Wiley, 2002).

]: the domain Γ is divided into several triangular subdomains (elements) in which locally supported expansion functions are defined; u is expanded in terms of such functions within each element; w is taken to be each one of the local expansion functions inside each element; and the material parameters p and q are allowed to be constant inside each element. Then the following quadratic matrix eigenvalue equation in k results:

[ A− ω2 c2D]u= ik ( C−B)u+ k2Cu.
(7)

Here, u is a vector containing the (complex) coefficients of the expansion of u, and matrices A through D may be individually related to each integral in Eq.(6) by inspection. The results in Sec. 3 are all calculated using second-order Lagrange elements. Explicit expressions for the matrices in this case can be found in [13

J. Jin, The Finite Element Method in Electromagnetics , (2nd ed. Wiley, 2002).

]. The eigenvalue equation may be solved at a fixed frequency ω (and thus fixed ε(ω)) and k-vector direction. The eigenvalue itself is the k-vector amplitude, k. Due to the periodicity of the system, Re{k} is expected to be periodic in k-space. Note that, for real ε, all matrices are Hermitian, so that the eigenvalues k are either purely real or are complex conjugate pairs (k,k̄) [14

F. Tisseur and K. Meerbergen, “The quadratic eigenvalue problem,” SIAM Rev. 43, 235–286 (2001). [CrossRef]

]. The most common way of solving the Quadratic Eigenvalue Problem is by linearization, which results in a (linear) system twice the original size [10

E. Istrate, A. A. Green, and E. H. Sargent, “Behavior of light at photonic crystal interfaces,” Phys. Rev. B 71, 195122 (2005). [CrossRef]

,14

F. Tisseur and K. Meerbergen, “The quadratic eigenvalue problem,” SIAM Rev. 43, 235–286 (2001). [CrossRef]

]. Other possibilities are to use nonlinear iterative methods such as Nonlinear Inverse Iteration or Rayleigh Coefficient Iteration [8

A. Ruhe, “Algorithms for the nonlinear eigenvalue problem,” SIAM J. Numer. Anal. 10, 674–689 (1973). [CrossRef]

].

3. Rectangular array plasmonic crystal

Here, the TE polarization (magnetic field perpendicular to the crystal plane) band structure of the plasmonic crystal studied in [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

, 6

G. Shvets and Y. Urzhumov, “Electric and magnetic properties of sub-wavelength plasmonic crystals,” J. Opt. A: Pure Appl. Opt. 7, S23–S31 (2005). [CrossRef]

], consisting of a 2D rectangular array of infinite metallic cylinders in air, is analyzed. The metal permittivity was considered to follow the Drude model where ε=1-ω 2 p /ω(ω-iωc ), with ωp the plasma frequency and ωc the damping frequency related to the mean electron collision rate. The square unit cell has a lattice parameter a=c/ωp and the cylinder radii are R=0.45a. All of our results below are given in terms of the normalized frequency Ω=ω/ωp .

For comparison, the transmission (t) and reflection (r) coefficients for plane-waves incident on a finite (i.e., with N of unit cells) crystal slab were obtained using the same technique as described in [15

G. Shvets and Y. Urzhumov, “Negative index meta-materials based on two-dimensional metallic structures,” J. Opt. A: Pure Appl. Opt. 8, S122–S130 (2006). [CrossRef]

]. The incident plane wavew is in the Γ-X direction of the square-lattice Brillouin zone. Furthermone, an effective index n eff is calculated using expressions from the scattering-parameter technique of Smith et al. [15

G. Shvets and Y. Urzhumov, “Negative index meta-materials based on two-dimensional metallic structures,” J. Opt. A: Pure Appl. Opt. 8, S122–S130 (2006). [CrossRef]

, 16

D. R. Smith, D. C. Vier, T. Koschny, and C. M. Soukoulis, “Electromagnetic parameter retrieval from inhomogeneous metamaterials,” Phys. Rev. E 71, 036617 (2005). [CrossRef]

], i.e.:

neff= lnY (ω) iωΔ⁄c
(8)

with

X≡ ( 1− r2+ t2)⁄2t
(9)

and

Y≡ X± 1− X2
(10)

These expressions model a 1D Fabry-Pérot cavity for waves propagating with a single wavevector at each frequency. This in principle precludes their direct application in situations when multiple waves associated with different k-vectors exist simultaneously - for instance in systems displaying Fano-type resonances as shown below.

Figure 1 shows the calculated band diagram in the Γ-X direction of the square lattice irreducible Brillouin zone, with ωc =0, in terms of both kr and ki (k=kr +iki ). Modes for which ki =0 are marked with crosses, while complex-k modes are marked with circles and triangles. Duplicates of such eigensolutions are also produced by the formulation, with k′ r =kr ±2m·π/a (kr in the first B.Z., m=1,2, …) and identical ki . Such solutions only arise because they satisfy the established periodic boundary conditions at the edges of the unit-cell; they do not participate in the Bloch-mode expansion of the total field in the crystal. Notice as well that since the calculated bands are for a crystal with zero damping losses (wc =0), the existence of complex eigenvalues is completely unrelated to power dissipation.

Same as in the case of completely dielectric photonic crystals, real-k modes extend over the entire crystal, with a phase evolution given by the k(Ω) dispersion and with no amplitude variation apart from that of the periodic Bloch envelope. In Fig. 1, purely real eigenvalues are found in the ranges Ω≈0.2-0.36 and Ω≈0.54-0.6.

It has been suggested that imaginary k modes play an important role in representing the decaying field inside a finite or semi-infinite PC upon wave incidence at bandgap frequencies [10

E. Istrate, A. A. Green, and E. H. Sargent, “Behavior of light at photonic crystal interfaces,” Phys. Rev. B 71, 195122 (2005). [CrossRef]

]. The inference that the imaginary part of the Bloch wavevector relates to the field decay rate in this situation is supported by our analysis. In Fig. 1, non-real eigenvalues may be either purely imaginary (kr =0), or complex. Modes with purely imaginary k (found in the intervals Ω=0.43-0.52 and Ω=0.6-0.61) have no phase variation across a unit-cell, and are non-propagating. Complex k modes, on the other hand, present a phase variation concomitant with an exponential amplitude decay (or growth), which translate into interference effects in scattering from a finite crystal. In Fig. 1, complex k modes with |kr | <π/a are found between Ω=0.36 and 0.43, and between Ω=0.52 and 0.54. Bands with ki marked with triangles in Fig. 1(c) have real parts at the Brillouin zone boundaries, kr =±π/a, such that the phase varies in full cycles across a unit cell. These solutions will henceforth be referred to as zone-boundary modes.

Figure 2(a) shows the transmission and reflection coefficients for a crystal with N=4 and ωc /ωp =0.001. In Figs. 2(b) and (c), the band structure with kr and ki on the horizontal axis are presented once again. The superimposed thick curves are the effective k-vectors k eff=Ω·n eff, with effective index n eff obtained from from Eq. (8). A clear correspondence exists between the high-reflectivity band in the frequency range Ω=0.36-0.54 and the exclusive existence of non-real-k modes in Figs. 2(b) and 2(c). The correspondence between the the purely real band at Ω<≈ 0.36 and k eff is also apparent - in this case, disregarding high-decay-rate kr =±π modes, only single propagating modes are available at each frequency, which is an assumption of the n eff expression, as pointed out above. In the range Ω≈0.43-0.56, k eff coincides almost exactly with the lowest-decay-rate purely imaginary Bloch band. We infer that the waves travelling through the slab suffer an exponential attenuation given by the lower-amplitude imaginary Bloch k. In addition, waves related to the larger ki modes either are not as strongly excited, or simply decay too quickly to be meaningful in the determination of k eff.

Fig. 1. (a)Square-lattice plasmonic crystal band structure as function of kr . (b) Band structure as function of ki . Crosses indicate bands with ki =0; circles, ki ≠0; and triangles, kr =±π/a, ki ≠0
Fig. 2. (a)Amplitude-squared transmission and reflection coefficients, (b) Band structure as function of kr and (c) Band structure as function of ki for a.

The k eff resonances observed in the interval Ω ≈ 0.53-0.6, are associated with Fabry-Pérot resonances of waves related to the purely real band within this range. It is apparent, however, that the complete transmission transfer function also involves the coexisting purely imaginary band. Evidence for this is that the imaginary part of k eff follows the ki band (see Fig. 2(c)), rather than tending to zero. Furthermore, the transmission transfer functions for the N=4 structure in Fig. 3(b) can be roughly approximated by a function of the form t(Ω)=α(Ω)+i·β (Ω)·t FP(Ω), with t FP(Ω) the Fabry-Pérot cavity transmission coefficient. Examples of this may be seen in Fig. 3(a), in which the transmission coefficient amplitude |t| for 4- and 3-unit cell structures are compared to those obtained from the expression above, under the following conditions: the Fabry-Pérot cavity single-pass phase varied as kr ·N · a, with kr from the purely real band (N is the number of unit cells); the cavity reflection coefficient was taken as r=0.92; β=0.5; the term α(Ω) was made to vary as exp(-ki ·N·a), where ki is from the purely imaginary band; the total transmission was normalized to its maximum value. General features such as sharp peaks (related to resonances of t FP) and alternating maxima and minima in-between sharp resonances (related to minima of t FP with phases of respectively odd and even multiples of π/2), are similar in both numerical and analytical curves, and for both cavity lengths. Thus, the terms α(Ω) and β (Ω) may be related to the amplitudes of two waves at the output side of the slab: one which is transmitted without a phase change, corresponding to the purely imaginary k-vector Bloch mode, and a second that propagates dispersively given by the real-k band, leading to the Fabry-Pérot resonances. In this case, both purely real and imaginary bands have relevant contributions to the field within the finite slab and thus a strong influence on k eff. Conversely, n eff obtained from Eq.(8) does not allow one to discern the propagation characteristics of the different SPC modes. Notice as well that from Fig. 3(b), transmission resonances are quite sharp with null damping losses (ωc =0), with peaks approaching unity. Introduction of even small damping (ωc =0.001ωp ) causes strong degradation of the resonance quality factor. The asymmetric shape of the transmission transfer function around resonant peaks is similar to that of Fano absorption lineshapes, which in quantum systems arise from the interference of a discrete state with a continuum [12

U. Fano, “Effects of configuration interaction on intensities and phase shifts,” Phys. Rev. 124, 1866–1878 (1961). [CrossRef]

]. The transmission resonances presented here have a similar origin to these, namely the interference between waves forming the Fabry-Pérot resonant peaks and the continuum formed by the evanescent waves.

Fig. 3. (a) Comparison between numerically calculated transmission coefficient amplitudes for 3- and 4-unit-cell crystal slabs (thin lines) and the analytical expression |t(Ω)|= A ·| exp(-ki (Ω)·N·a)+i·β· t FP(ki (Ω)·N·a)| (thick lines), with ki (Ω) and kr (Ω) taken from the purely real and imaginary Bloch bands respectively. (b)Amplitude and phase of the transmission transfer function for the 4-unit-cell crystal with ωc /ωp =0.001 (thick line) and ωc =0 (thin line).

Finally, it is important to point out that the upper real-k bands discussed above are negative refraction bands, with group velocity vg of opposite sign to the phase velocity ω/k; however, as shown in [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

,6

G. Shvets and Y. Urzhumov, “Electric and magnetic properties of sub-wavelength plasmonic crystals,” J. Opt. A: Pure Appl. Opt. 7, S23–S31 (2005). [CrossRef]

], the smallest number of rods that can be meaningfully analyzed in the search for optical magnetism is four, since each rod must be surrounded by other rods to create a magnetically active resonance.

In the range Ω≈0.36-0.43, the transmitted power is extremely low, yet the transfer function presents a series of anti-resonances corresponding to the observed peaks in Im{k eff} with concomitant step-transitions in Re{k eff}. At the same time, disregarding the resonances, both real and imaginary parts of the effective wavevector closely follow the corresponding photonic bands. Hence, this behavior has an origin similar to the resonances discussed above, given the existence of two evanescent Bloch bands, as seen in Figs. 1(a) and 1(b).

We next analyze the effects or small damping losses in the calculated band structure. In Fig. 4(a) and 4(b), the real and imaginary k bands for a lossless and a lossy structure with ωc /ωp =0.001 are compared. A noticeable difference between the two cases is the disappearance of the the small bandgap centered at Ω ≈ 0.36. Whereas in the lossless case the purely real kr -bands reach the Brillouin zone boundary at gap frequencies, in the lossy structure the corresponding bands do not reach kr =π/a, turning around at some k max <π/a where the group velocity defined as vg =dω/dkr becomes infinite. Notice that in the lossless case the group velocity as defined above is infinite throughout the photonic bandgap, since kr =π/a everywhere in this range. It is apparent that the band turn-around observed in the lossy crystal case is caused by a less-efficient Bragg wave-interference in the presence of damping losses. This effect was also observed in [17

K. C. Huang, E. Lidorikis, X. Jiang, J. D. Joannopoulos, K. A. Nelson, P. Bienstman, and S. Fan, “Nature of lossy bloch states in polaritonic photonic crystals,” Phys. Rev. B 69, 195111 (2004). [CrossRef]

], for photonic crystals containing polaritonic inclusions. It must be pointed out that within the bandgap frequency region, the group velocity as defined above does not correspond to the velocity with which an incident wavepacket would traverse a finite length of the SPC, so that vg >c does not imply superluminal propagation [18

H. G. Winful, “The meaning of group delay in barrier tunnelling: a re-examination of superluminal group velocities,” New J. Phys., Phys. 8, 101 (2006). [CrossRef]

].

Fig. 4. (a) Band structure for (a)kr and (b)ki . Dots correspond to ωc =0, crosses to ωc /ωp =0.001.

As mentioned above, in the lossless case, the matrices involved are Hermitian, such that eigenvalues are either purely real or come in complex conjugate pairs. The latter case is observed in Ω=0.36-0.43, where two degenerate modes at each kr exist. The degeneracy is evident from the opening of bifurcations in ki at Ω≈ 0.36 and splitting of ki bands in Ω=0.36-0.43 due to small damping losses. The singled-out band in Fig. 4(a) is a perturbed version of a zone-boundary band that exists at Ω⪷0.365 in the lossless case. Above Ω⪷0.365, the latter ceases to be a zone-boundary band, and becomes a degenerate complex-k band. This transition is evident in the departure of the perturbed band from kr =π/a and its merging with a second complex-k band around Ω=0.36.

Finally, to show that the method may be conveniently applied to the generation of bands in arbitrary directions in the k-vector space, Figs. 5 and 6 display, respectively, the low-frequency and two high-frequency purely real-k bands. These diagrams were obtained simply by changing the direction of the unit k-vector in Eq. (6).

Fig. 5. Low-frequency, purely real bands of the square-lattice plasmonic crystal over the first Brillouin zone.

Figure 5 shows that the lowest purely real band is highly isotropic at lower frequencies, however becomes anisotropic as the photonic bandgap is approached. Figure 6 shows slices of the upper real bands, which in reality form a closed surface. The negative group velocity band Ω=0.57-0.6 is highly isotropic near its top, thus making the square array of circular rods a perfectly isotropic negative-index metamaterial for the range of frequencies Ω=0.58-0.6.

Fig. 6. High-frequency, purely real bands of the square-lattice plasmonic crystal over the first Brillouin zone.

On the other hand, the longitudinal-wave, or Bulk Plasmon, band near Ω=0.0.61 predicted in [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

] displays very strong anisotropy near the Γ-point. This anisotropy is apparent in Fig. 7, where the Bulk Plasmon band was plotted for the Γ-M (0°) and Γ-X (45°) directions.

Fig. 7. Bulk Plasmon band along Γ-M (0°) and Γ-X (45°) directions.

It is important to notice that, in [5

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

] and [6

G. Shvets and Y. Urzhumov, “Electric and magnetic properties of sub-wavelength plasmonic crystals,” J. Opt. A: Pure Appl. Opt. 7, S23–S31 (2005). [CrossRef]

], the origin of the negative refraction band above was related to the hybridization of quasi-static multipole resonances of the plasmonic rods in close proximity to each other. Over the negative-index band, the lattice constant is ~λ/10, where λ is the vacuum wavelength. Thus we believe that negative refraction occurs in the longwavelength rather than in the Bragg diffraction regime.

4. Conclusions

The Finite-Element formulation presented in Section 2 is convenient for the calculation o photonic band structures of periodic metamaterials composed of dispersive materials. The formulation yields a quadratic eigenvalue equation system that may be easily implemented with commercial, generic Finite-Element packages such as COMSOLMultiphysics, and that may be more efficiently solved than the general nonlinear eigenvalue problem that is normally solved in the metamaterials community. In addition, calculation of bands in arbitrary directions in k-space simply involves determination of the components of the unit k-vector, which is convenient for the analysis of medium anisotropy.

An important advantage associated with solving the quadratic eigenvalue problem is that it yields Bloch eigenmodes with not only real but also complex wavevectors, which, as suggested in [10

E. Istrate, A. A. Green, and E. H. Sargent, “Behavior of light at photonic crystal interfaces,” Phys. Rev. B 71, 195122 (2005). [CrossRef]

], may be associated with many effects observed in the scattering of waves from finite crystals, especially at, but not limited to, bandgap frequency ranges. A clear indication of this in fact is the existence of the Fano-type resonances in transmission as shown in Section 3.

It is clear from the analysis in Section 3 that the effective index n eff determined from Eq.(8) via scattering from finite crystals will not completely coincide with effective indices of particular Bloch modes when multiple such modes exist that influence the scattering process at similar levels. At the same time, just the Bloch band structure by itself does not give an immediately clear picture of wave scattering and therefore should be used with some attention in predicting the scattering response of finite structures. At any rate, both methods may be used in conjunction to form a more complete understanding of the propagation characteristics of a metamaterial.

Finally, it must be pointed out that alternative formulations to the present finite element-based one presented in this article may also lead to quadratic eigenvalue equations, see for instance the plane-wave expansion formulation used in [10

E. Istrate, A. A. Green, and E. H. Sargent, “Behavior of light at photonic crystal interfaces,” Phys. Rev. B 71, 195122 (2005). [CrossRef]

].

Acknowledgements

The authors thank Dr. Stephen Forrest for support, encouragement and useful comments. This work was supported by the Air Force Office of Scientific Research under contract FA 9550-06-01-0279 through the Multidisciplinary University Research Initiative Program, as well as ARO MURI W911NF-04-01-0203, and DARPA contract HR0011-05-C-0068.

References and links

1.

J. D. Joannopoulos, R. D. Meade, and J. N. Winn, Photonic Crystals: Molding the Flow of Light . (Princteon University Press, 1995).

2.

K. Sakoda, Optical Properties of Photonic Crystal, ser. Optical Sciences. (New York: Springer, 2001).

3.

S. G. Johnson and J. D. Joannopoulos, “Block-iterative frequency-domain methods for maxwell’s equations in a planewave basis,” Opt. Express 8, 173–190 (2001), http://www.opticsexpress.org/abstract.cfm?URI=OPEX-8-3-173. [CrossRef] [PubMed]

4.

B. P. Hiett, B. D. H., S. J. Cox, J. M. Generowicz, M. Molinari, and K. S. Thomas, “Aplication of finite element methods to photonic crystal modelling,” IEE Proc - Sci. Meas. Technol. 149, 293–296 (2002). [CrossRef]

5.

G. Shvets and Y. A. Urzhumov, “Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances,” Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]

6.

G. Shvets and Y. Urzhumov, “Electric and magnetic properties of sub-wavelength plasmonic crystals,” J. Opt. A: Pure Appl. Opt. 7, S23–S31 (2005). [CrossRef]

7.

C. Enkrich, M. Wegener, S. Linden, S. Burger, L. Zschiedrich, F. Schmidt, J. F. Zhou, T. Koschny, and C. M. Soukoulis, “Magnetic metamaterials at telecommunication and visible frequencies,” Phys. Rev. Lett. 95, 203901 (2005), doi:10.1103/PhysRevLett.95.203901. [CrossRef] [PubMed]

8.

A. Ruhe, “Algorithms for the nonlinear eigenvalue problem,” SIAM J. Numer. Anal. 10, 674–689 (1973). [CrossRef]

9.

A. Spence and C. Poulton, “Photonic band structure calculations using nonlinear eigenvalue techniques,” J. Comput. Phys. 204, 65–81 (2005). [CrossRef]

10.

E. Istrate, A. A. Green, and E. H. Sargent, “Behavior of light at photonic crystal interfaces,” Phys. Rev. B 71, 195122 (2005). [CrossRef]

11.

[Online]. Available: http://www.comsol.com

12.

U. Fano, “Effects of configuration interaction on intensities and phase shifts,” Phys. Rev. 124, 1866–1878 (1961). [CrossRef]

13.

J. Jin, The Finite Element Method in Electromagnetics , (2nd ed. Wiley, 2002).

14.

F. Tisseur and K. Meerbergen, “The quadratic eigenvalue problem,” SIAM Rev. 43, 235–286 (2001). [CrossRef]

15.

G. Shvets and Y. Urzhumov, “Negative index meta-materials based on two-dimensional metallic structures,” J. Opt. A: Pure Appl. Opt. 8, S122–S130 (2006). [CrossRef]

16.

D. R. Smith, D. C. Vier, T. Koschny, and C. M. Soukoulis, “Electromagnetic parameter retrieval from inhomogeneous metamaterials,” Phys. Rev. E 71, 036617 (2005). [CrossRef]

17.

K. C. Huang, E. Lidorikis, X. Jiang, J. D. Joannopoulos, K. A. Nelson, P. Bienstman, and S. Fan, “Nature of lossy bloch states in polaritonic photonic crystals,” Phys. Rev. B 69, 195111 (2004). [CrossRef]

18.

H. G. Winful, “The meaning of group delay in barrier tunnelling: a re-examination of superluminal group velocities,” New J. Phys., Phys. 8, 101 (2006). [CrossRef]

OCIS Codes
(000.4430) General : Numerical approximation and analysis
(160.4670) Materials : Optical materials
(260.2030) Physical optics : Dispersion
(260.2110) Physical optics : Electromagnetic optics

ToC Category:
Physical Optics

History
Original Manuscript: June 1, 2007
Revised Manuscript: July 15, 2007
Manuscript Accepted: July 16, 2007
Published: July 19, 2007

Citation
Marcelo Davanco, Yaroslav Urzhumov, and Gennady Shvets, "The complex Bloch bands of a 2D plasmonic crystal displaying isotropic negative refraction," Opt. Express 15, 9681-9691 (2007)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-15-15-9681


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. J. D. Joannopoulos, R. D. Meade, and J. N. Winn, Photonic Crystals: Molding the Flow of Light. (Princteon University Press, 1995).
  2. K. Sakoda, Optical Properties of Photonic Crystal, ser. Optical Sciences. (New York: Springer, 2001).
  3. S. G. Johnson and J. D. Joannopoulos, "Block-iterative frequency-domain methods for maxwell’s equations in a planewave basis," Opt. Express 8, 173-190 (2001). [CrossRef] [PubMed]
  4. B. P. Hiett, B. D. H., S. J. Cox, J. M. Generowicz, M. Molinari, and K. S. Thomas, "Aplication of finite element methods to photonic crystal modelling," IEE Proc - Sci. Meas. Technol. 149, 293-296 (2002). [CrossRef]
  5. G. Shvets and Y. A. Urzhumov, "Engineering the electromagnetic properties of periodic nanostructures using electrostatic resonances," Phys. Rev. Lett. 93, 243902 (2004). [CrossRef]
  6. G. Shvets and Y. Urzhumov, "Electric and magnetic properties of sub-wavelength plasmonic crystals," J. Opt. A: Pure Appl. Opt. 7, S23-S31 (2005). [CrossRef]
  7. C. Enkrich, M. Wegener, S. Linden, S. Burger, L. Zschiedrich, F. Schmidt, J. F. Zhou, T. Koschny, and C. M. Soukoulis, "Magnetic metamaterials at telecommunication and visible frequencies," Phys. Rev. Lett. 95, 203901 (2005). [CrossRef] [PubMed]
  8. A. Ruhe, "Algorithms for the nonlinear eigenvalue problem," SIAM J. Numer. Anal. 10, 674-689 (1973). [CrossRef]
  9. A. Spence and C. Poulton, "Photonic band structure calculations using nonlinear eigenvalue techniques," J. Comput. Phys. 204, 65-81 (2005). [CrossRef]
  10. E. Istrate, A. A. Green, and E. H. Sargent, "Behavior of light at photonic crystal interfaces," Phys. Rev. B 71, 195122 (2005). [CrossRef]
  11. [Online]. Available: http://www.comsol.com
  12. U. Fano, "Effects of configuration interaction on intensities and phase shifts," Phys. Rev. 124, 1866-1878 (1961). [CrossRef]
  13. J. Jin, The Finite Element Method in Electromagnetics, (2nd ed. Wiley, 2002).
  14. F. Tisseur and K. Meerbergen, "The quadratic eigenvalue problem," SIAM Rev. 43, 235-286 (2001). [CrossRef]
  15. G. Shvets and Y. Urzhumov, "Negative index meta-materials based on two-dimensional metallic structures," J. Opt. A: Pure Appl. Opt. 8, S122-S130 (2006). [CrossRef]
  16. D. R. Smith, D. C. Vier, T. Koschny, and C. M. Soukoulis, "Electromagnetic parameter retrieval from inhomogeneous metamaterials," Phys. Rev. E 71, 036617 (2005). [CrossRef]
  17. K. C. Huang, E. Lidorikis, X. Jiang, J. D. Joannopoulos, K. A. Nelson, P. Bienstman, and S. Fan, "Nature of lossy bloch states in polaritonic photonic crystals," Phys. Rev. B 69, 195111 (2004). [CrossRef]
  18. H. G. Winful, "The meaning of group delay in barrier tunnelling: a re-examination of superluminal group velocities," New J. Phys., Phys. 8, 101 (2006). [CrossRef]

Cited By

 Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.


« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited