OSA's Digital Library

Optics Express

Optics Express

  • Editor: J. H. Eberly
  • Vol. 8, Iss. 3 — Jan. 29, 2001
  • pp: 167–172
« Show journal navigation

Photonic band structure computations

Daniel Hermann, Meikel Frank, Kurt Busch, and Peter Wölfle  »View Author Affiliations


Optics Express, Vol. 8, Issue 3, pp. 167-172 (2001)
http://dx.doi.org/10.1364/OE.8.000167


View Full Text Article

Acrobat PDF (187 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

We introduce a novel algorithm for band structure computations based on multigrid methods. In addition, we demonstrate how the results of these band structure calculations may be used to compute group velocities and effective photon masses. The results are of direct relevance to studies of pulse propagation in such materials.

© Optical Society of America

1 Introduction

Over the last years, tremendous progress in microfabricating two- and three-dimensional compounds with a spatially periodic dielectric function has lead to a dramatically increased interest in studying these so-called photonic crystals [1

1. C.M. Soukoulis (Ed.), Photonic Band Gap Materials, NATO ASI Series E315, Kluwer Academic Publishers1996 [CrossRef]

]. In such structures, Bragg-scattering of electromagnetic waves leads to the formation of a photonic band-structure, the engineering of which offers enormous potential for applications both in basic science as well as in telecommunication. Due to the advanced state in various planar lithography techniques, the fabrication of two-dimensional photonic crystals has matured to a point which suggests that the first practical applications are imminent [2

2. A. Birner et al., “Macroporous Silicon: A Two-Dimensional Photonic Bandgap Material Suitable for the Near-Infrared Spectral Range,” phys. stat. sol (a) 165, 111–117 (1998) [CrossRef]

]. The theoretical description of photonic crystals is generally based on methods of electronic bandstructure theory and is playing a key role in characterizing existing as well as designing novel structures. For instance, the seminal paper of Ho et al. [3

3. K.M. Ho, C.T. Chan, and C.M. Soukoulis, “Existence of a photonic band gap in periodic structures,” Phys. Rev. Lett. 65, 3152–3155 (1990) [CrossRef] [PubMed]

] described the first structure that exhibits a complete photonic band gap (PBG), i.e., a range of frequencies for which propagation is forbidden irrespective of the direction of propagation. This work suggested that the formation of a complete PBGis favored by structures exhibiting diamond or hexagonal symmetry and has subsequently lead to the fabrication of the so-called “Yablonovite” and Lincoln-log structures.

In this paper we describe in detail the calculation of photonic band and associated electromagnetic mode structures of two-dimensional photonic crystals. In Section 2, we introduce a novel algorithm for efficient computation of the electromagnetic mode structures. In addition, we outline the calculation of group velocities and effective photon masses in section 3.

2 Two-dimensional photonic crystals

We consider a two-dimensional photonic crystal consisting of a periodic lattice of infinitely extended parallel rods with diameter d. The rods are made from an isotropic dielectric material (dielectric constant a ) and are embedded in a matrix of dielectric material characterized by an isotropic dielectric constant b . The direction of the rods defines the z-axis. For electromagnetic waves propagating in the xy-plane, the two transverse polarizations decouple leading to two separate scalar problems. For the TM-polarization (E-field parallel to the rods, i.e., E⃗ (r⃗)≡(0, 0, E(r⃗)), we obtain from Maxwell’s equations

1(r)(x2+y2)E(r)+ω2c2E(r)=0,
(1)

while for the TE-polarization (H-field parallel to the rods, i.e., H⃗ (r⃗)≡(0, 0,H(r⃗)), we have

x+(1(r)xH(r))+y(1(r)yH(r))+ω2c2H(r)=0.
(2)

Here, x/∂x,y/y , and r⃗ (x, y) is a two-dimensional vector. The lattice defined through the centers of cylinders is given by the set R={R⃗=n 1 a1+n 2 a⃗2,n 1, n 2N} of two-dimenisonal lattice vectors R⃗ that are generated by the primitive translations a1 and a2. The corresponding reciprocal (dual) lattice is defined through G={G⃗:G⃗·R⃗=2πn, nN,R∊R}. The dielectric function (r⃗)=b -(a -b )∑ Rθ(d/2-|r⃗-R⃗|) is a lattice periodic function which contains all the information of the photonic crystal. Here, θ(r) denotes the Heaviside step function. Eq. (1) and (2) represent differential equations with periodic coefficients. Therefore, their solutions obey the Bloch-Floquet theorem

Ek(r+ai)=eikaiEk(r)
(3)
Hk(r+ai)=eikaiHk(r),
(4)

where i=1, 2. The wave vector k⃗∊1.BZ is a vector of the first Brillouin zone (BZ) known as the crystal momentum.

A straightforward way to solving Eqs. (1) and (2) is to expand all the periodic functions into Fourier series over the reciprocal lattice G, thereby transforming the differential equations into a infinite matrix eigenvalue problem, which is then truncated and solved numerically. Details of this plane wave method (PWM) for isotropic systems can, for instance, be found in [4

4. K. Busch and S. John, “Photonic band gap formation in certain self-organizing systems,” Phys. Rev. E 58, 3896–3908 (1998) [CrossRef]

] and for anisotropic systems in [5

5. K. Busch and S. John. “Liquid crystal photonic band gap materials: The tunable electromagnetic vacuum,” Phys. Rev. Lett. 83, 967–970 (1999) [CrossRef]

].

3 Efficient computation of the mode structure

While the PWM provides a straightforward approach to computing the band structure of photonic crystals, it also exhibits a number of shortcomings. The truncation of Fourier series in the presence of discontinuous changes in the dielectric function may lead to serious convergence problems. In particular, we find that an accurate computation of the electromagnetic mode structure through PWM generally requires more than 2000 plane waves. This leads to very large matrices for which not only eigenvalues but also eigenvectors have to be computed. Finally, the Fourier series for the Bloch functions have to be evaluated. This clearly is a formidable task. In what follows, we describe the details of a new approach to efficiently compute the electromagnetic mode structure of two-dimensional photonic crystals.

The pairs of Eqs. (1) and (3) as well as Eqs. (2) and (4) comprise elliptic partial differential eigenproblems with Bloch-type boundary conditions at the boundaries of the unit cell. Upon a real-space discretization of the differential operators, both problems for any wave vector k⃗ may be cast into the matrix form of a linear system of equations

kUi=ΛiUi,
(5)

where Λ i =ωi2 /c 2 and Ui represents the discretized fields E k(r⃗) and H k(r⃗) in case of TM- and TE-polarization, respectively. The precise form of the operator matrix L k and the vectors Ui depend on the type of polarization (TE or TM) and the type of discretization (finite difference or finite element method). Mathematically, the problem is then to find approximations to the first few eigenfrequencies (bands) Λ1≤Λ2≤Λ3≤… and associated Bloch functions U 1; U 2; U 3…. Due to the local nature of the underlying differential operators, the operator matrix L k is sparse, which suggests that Eq. (5) should be solved iteratively by Rayleigh quotient iteration or Lanczos. Such approaches of treating the discrete problem as a purely algebraic one can result in loss of valuable information, especially concerning the smoothness of the Bloch functions. In general, the Bloch functions corresponding to the desired first bands are very smooth, so that they are well approximated on coarser grids. Certain multigrid methods take full advantage of this smoothness and are therefore very effective for solving such problems [6

6. P. Wesseling, An Introduction to Multigrid Methods, John Wiley & Sons (1992)

]. In addition, multigrid methods can be applied to nonlinear problems such as Eq. (5), where both eigenfrequency and Bloch functions are unknown. To do so, Eq. (5) has to be amended by a ortho-normalization condition [7

7. A. Brandt, S. McCormick, and J. Ruge, “Multigrid methods for differential eigenproblems,” SIAM J. Sci. Stat. Comput. 4, 244–260 (1983). [CrossRef]

]

1Vwscd2rEnk*(r)(r)Emk(r)=δnm
(6)
1Vwscd2rHnk*(r)Hmk(r)=δnm,
(7)

where V denotes the volume of the Wigner-Seitz cell (WSC) and δnm is the Kronecker symbol for the band indices n and m.

The efficiency of multigrid methods results from the fact that, although iteration (relax-ation) of Eq. (5) using a simplified version of the matrix operator L k (Jacobi- or Gauss-Seidel iteration [6

6. P. Wesseling, An Introduction to Multigrid Methods, John Wiley & Sons (1992)

]) is usually slow to converge, it is quick to reduce high-frequency error components [6

6. P. Wesseling, An Introduction to Multigrid Methods, John Wiley & Sons (1992)

]. This allows the problem to be transfered to a coarser grid where the error can be resolved with much less work. Not only is the relaxation cheaper per sweep on coarser grids, but the solution process is also much more effective. The coarse grid equation can be solved by relaxation and appeal to still coarser grids. The coarsest grid used is chosen so that solution of the problem there is inexpensive compared to the work performed on the fine grid. The number of relaxation sweeps needed to smooth the error on each grid is generally small.

We discretize Eqs. (1) and (2) using a finite volume discretization scheme on a uniform mesh and employ a simple white-black ordered Gauss-Seidel relaxation method [6

6. P. Wesseling, An Introduction to Multigrid Methods, John Wiley & Sons (1992)

]. The multigrid (MG) iteration follows a V-cycle as illustrated in Fig. 1. The relaxation steps on all but the coarsest level are performed with fixed eigenfrequencies Λ i i=1,2,…, q. The eigenfrequencies are updated after every relaxation step on the base level using the Rayleigh quotient

Λi=k1Ui1,Ui1Ui1,Ui1,
(8)
Fig. 1. Graphical illustration of the V-cycle used in the multigrid iteration of Eq. (5) depicting four levels of grids. A movie (881.mov(0.8M)) shows the multigrid iteration using two levels for the third band at the x-point in the TM-polarized case (see Fig. 2).

where the superscript 1 indicates the representation of L k and Ui on the lowest level and 〈.,.〉 stands for the discrete version of scalar products defined in Eq. (6). This V-cycle scheme is repeated until stable values for eigenfrequencies Λ i and Bloch-functions Ui emerge (generally two V-cycles are sufficient) [7

7. A. Brandt, S. McCormick, and J. Ruge, “Multigrid methods for differential eigenproblems,” SIAM J. Sci. Stat. Comput. 4, 244–260 (1983). [CrossRef]

].

As an example, we consider a square lattice of dielectric rods (a =13) in air (b =1). Then, the primitive translations are given by a1=(a, 0) and a2=(0, a), where a is the lattice constant. The radius r of the rods is r/a=0.45. The photonic band structure of this structure is shown in Fig. 2. In Tab. 1, we compare the band structure data

Fig. 2. Photonic band structure for a square lattice of dielectric rods (a =13, r/a=0.45) in air (b =1) for TM-polarization (blue) and TE-polarization (green).

from the MG-method for the first four bands of TM-polarization using different mesh sizes with corresponding data using PWM with different numbers N of plane waves. We observe, that for large numbers of plane waves, the PWM aproaches the values of the MG-method to within 1% but requires about 100% more CPU-time. In addition, the result from the MG-computations suggest that a 128×128 mesh is sufficient to obtain converged results. It should be noted that in this example PWM was used to compute the band structure only. If an evaluation of the electromagnetic mode structure were

Table 1. Comparison of numerical data obtained from PWM (upper part) and the multigrid method (lower part) using different numbers N of plane waves and different mesh-sizes (128×128 and 256×256), respectively. The values for the defect give an estimate of the absolute error of the Bloch-functions [7]. The computations have been performed at the X-point for TM-polarization using a 466 MHz Celeron processor.

table-icon
View This Table

required, the computational effort would increase substantially leading to much longer CPU-times as quoted above. In contrast, the MG-method obtains the Bloch functions without any additional numerical effort.

4 Computation of group velocities and effective photon masses

In order to understand pulse propagation in photonic crystals, additional calculations have to be performed. In this section we demonstrate how to obtain group velocities and effective masses characterizing the speed of pulse propagation and the broadening of pulses. This approach closely follows the so-called kp-perturbation theory (kp-PT) of electronic band structure theory and has been applied to one-dimensional systems [8

8. C. Martijn de Sterke and J.E. Sipe, “Envelope-function approach for the electrodynamics of nonlinear periodic media,” Phys. Rev. A 38, 5149–5165 (1988) [CrossRef]

]. For simplicity, we will only discuss the TM-polarization.

We begin by writing the Bloch-Floquet theorem, Eq. (3), as

Ek(r)=eikruk(r),
(9)

where u k(r⃗) is a lattice periodic function. The equation of motion for E k(r⃗), Eq. (1) may now be transformed into an equation of motion for u k(r⃗)

(Δ+2i·kk2)uk(r)+ωk2c2(r)uk(r)))=0.
(10)

Introducing the operator Ĥ (k⃗)=Δ+2i∇·k⃗+k2, where Δ=x2 +y2, allows to cast the equation of motion, Eq. (10), for u k⃗+q(r⃗) (|q⃗≪|π/a) into the form

Ĥ(k)uk+q(r)+(2q·Ωq2)uk+q(r)+ωk+q2c2(r)uk+q(r)=0.
(11)

Here, we have introduced Ω⃗=-i(∇+ik⃗). Eq. (11) together with the fact that |q⃗|≪π/a suggests that we treat the second and third part on its l.h.s as a perturbation on Ĥ(k⃗) such that the eigenfrequency ω k⃗+q⃗ may be regarded as a perturbed eigenvalue and compared to a Taylor-expansion of ω k⃗+q⃗ around k

ωk+q=ωk+(1storderkpPT)+(2ndorderkpPT)+
(12)
=ωk+q·υk+q·Mk·q+.
(13)
Fig. 3. Group velocity of the first three bands (first band: black; second band: red; third band: green) for a square lattice of dielectric rods (a =13, r/a=0.45) in air (b =1) for TM-polarized light. The corresponding band structure is displayed in Fig. 2.

The computation of group velocities υk=∂ kω k and the inverse effective photon mass tensors (group velocity dispersion tensors) Mkij=ki kj ω k, i=1, 2 are thus reduced to second order perturbation theory [8

8. C. Martijn de Sterke and J.E. Sipe, “Envelope-function approach for the electrodynamics of nonlinear periodic media,” Phys. Rev. A 38, 5149–5165 (1988) [CrossRef]

].

In Fig. 3 we show the group velocity for the first three TM-polarized bands along the high-symmetry lines of the 1.BZ. As expected from the corresponding band structure in Fig. 2, the higher bands exhibit rather low values of the group velocity. In particular, the values of group velocity for band 3 along the Γ-X direction are about an order of magnitude smaller than the group velocity in band 1, making band 3 an excellent candidate for realizing nonlinear optical properties in this material.

5 Discussion

In summary, we have introduced a novel scheme for computing photonic band and associated electromagnetic mode structures based on multigrid methods. This approach is much more efficient than the plane wave method especially when the electromagnetic mode structure is required as, for instance, is the case when calculating group velocities and effective photon masses in photonic crystal. The MG-approach may be further refined using finite elements instead of finite differences as well as adaptive instead of uniform grids. An extension of this approach to three-dimensional photonic band structures is straightforward. In addition, we have shown how group velocities and effective photon masses may be calulated directly from photonic band structure computations through the photonic analogue of the well-known kp-perturbation theory. These results are of direct relevance to studies of pulse propagation in photonic crystals.

6 Acknowledgments

K.B. acknowledges the financial support of the Deutsche Forschungsgemeinschaft under Bu 1107/2-1 (Emmy-Noether program).

References and links

1.

C.M. Soukoulis (Ed.), Photonic Band Gap Materials, NATO ASI Series E315, Kluwer Academic Publishers1996 [CrossRef]

2.

A. Birner et al., “Macroporous Silicon: A Two-Dimensional Photonic Bandgap Material Suitable for the Near-Infrared Spectral Range,” phys. stat. sol (a) 165, 111–117 (1998) [CrossRef]

3.

K.M. Ho, C.T. Chan, and C.M. Soukoulis, “Existence of a photonic band gap in periodic structures,” Phys. Rev. Lett. 65, 3152–3155 (1990) [CrossRef] [PubMed]

4.

K. Busch and S. John, “Photonic band gap formation in certain self-organizing systems,” Phys. Rev. E 58, 3896–3908 (1998) [CrossRef]

5.

K. Busch and S. John. “Liquid crystal photonic band gap materials: The tunable electromagnetic vacuum,” Phys. Rev. Lett. 83, 967–970 (1999) [CrossRef]

6.

P. Wesseling, An Introduction to Multigrid Methods, John Wiley & Sons (1992)

7.

A. Brandt, S. McCormick, and J. Ruge, “Multigrid methods for differential eigenproblems,” SIAM J. Sci. Stat. Comput. 4, 244–260 (1983). [CrossRef]

8.

C. Martijn de Sterke and J.E. Sipe, “Envelope-function approach for the electrodynamics of nonlinear periodic media,” Phys. Rev. A 38, 5149–5165 (1988) [CrossRef]

OCIS Codes
(260.2110) Physical optics : Electromagnetic optics
(350.3950) Other areas of optics : Micro-optics

ToC Category:
Focus Issue: Photonic bandgap calculations

History
Original Manuscript: November 13, 2000
Published: January 29, 2001

Citation
Daniel Hermann, Meikel Frank, Kurt Busch, and Peter Wolfle, "Photonic band structure computations," Opt. Express 8, 167-172 (2001)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-8-3-167


Sort:  Journal  |  Reset  

References

  1. C.M. Soukoulis (Ed.), Photonic Band Gap Materials, NATO ASI Series E 315, (Kluwer Academic Publishers, 1996). [CrossRef]
  2. A. Birner et al., "Macroporous Silicon: A Two-Dimensional Photonic Bandgap Material Suitable for the Near-Infrared Spectral Range," phys. stat. sol (a) 165, 111-117 (1998) [CrossRef]
  3. K.M. Ho, C.T. Chan, and C.M. Soukoulis, "Existence of a photonic band gap in periodic structures," Phys. Rev. Lett. 65, 3152-3155 (1990). [CrossRef] [PubMed]
  4. K. Busch and S. John, "Photonic band gap formation in certain self-organizing systems," Phys. Rev. E 58, 3896-3908 (1998). [CrossRef]
  5. K. Busch and S. John. "Liquid crystal photonic band gap materials: The tunable electromagnetic vacuum," Phys. Rev. Lett. 83, 967-970 (1999). [CrossRef]
  6. P. Wesseling, An Introduction to Multigrid Methods, (John Wiley & Sons, 1992).
  7. A. Brandt, S. McCormick, and J. Ruge, "Multigrid methods for differential eigenproblems," SIAM J. Sci. Stat. Comput. 4, 244-260 (1983). [CrossRef]
  8. C. Martijn de Sterke and J.E. Sipe, "Envelope-function approach for the electrodynamics of non-linear periodic media," Phys. Rev. A 38, 5149-5165 (1988) [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.

Figures

Fig. 1. Fig. 2. Fig. 3.
 

Multimedia

Multimedia FilesRecommended Software
» Media 1: MOV (813 KB)      QuickTime

« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited