Open Access
Add to BibSonomyAbstract
The reflection coefficient is one important parameter of the perfectly matched layer (PML). Here we investigate its effect on the modal analysis of leaky waveguide modes by examining three different leaky waveguide structures, i.e., the holey fiber, the air-core terahertz pipe waveguide, and the gain-guided and index-antiguided slab waveguide. Numerical results reveal that the typical values 10−8 ~10−12 are inadequate for obtaining the imaginary part of the complex propagation constant, and the suggested reflection coefficient would be much smaller, for example, 10−50 or 10−100. With such a small coefficient, both the computational window size and the PML thickness can be significantly reduced without loss of stability. Moreover, in some cases, the modal field profiles can only be accurately obtained with such a small coefficient.
©2011 Optical Society of America
1. Introduction
Modal analysis of optical waveguides plays an important role in the design of optical systems and photonic integrated circuits. Usually, the waveguide structure is simulated within a finite space, and an artificial boundary condition is imposed to truncate the computational window. The easily used zero boundary condition or the transparent boundary condition [1] might be satisfactory for the analysis of guided modes. While for the leaky-mode analysis, owing to the outgoing characteristics of the modal fields, a more powerful boundary condition is required. Since the effectiveness of the perfectly matched layer (PML) [2–4] has been confirmed to provide reflectionless absorption of the incident fields for the finite-difference time-domain (FDTD) method [5], in recent years, the PML is popularly employed as the boundary condition for the eigenmode analysis of the leaky waveguide modes [6–9].
To solve the waveguide modes with PML, in the PML region, the derivative with respect to α (α=x or y) appearing in Maxwell’s equations can be transformed by the complex coordinate stretching as [3,4]:
where the coordinate stretching factor s is given by [7,8]in which m is usually taken as 2, λ is the free-space wavelength, n is the refractive index of the adjacent medium in the computational window, d is the thickness of the PML, ρ is the distance from the beginning of the PML, and R is the theoretical reflection coefficient at the interface between the computational window and the PML for the normal incident wave [4,10]. The reflection coefficient R is an important parameter of the PML, especially when the leaky modes are to be treated. Typically, R is explicitly assigned and ranges from 10−8 to 10−12 [7–9]. However, its effect on the modal analysis of the leaky waveguide modes, including calculating the complex propagation constants β and solving the modal field distributions, is seldom addressed in the literature. Recently, with the modified PML [11], the effect of R on calculating β for guided modes was illustrated using a pseudospectral mode solver [12], but the effect for leaky modes was not discussed in detail.In this paper, we investigate the effect of PML reflection coefficient on leaky-mode analysis by using the finite-difference frequency-domain (FDFD) method [8]. Three different leaky structures are examined, which are a photonic crystal fiber (holey fiber) with a ring of six air holes [13], an air-core terahertz (THz) pipe waveguide [14], and a gain-guided and index-antiguided (GG + IAG) slab waveguide [15]. Numerical results reveal that the typical values 10−8 ~10−12 of R are only adequate for obtaining the real part of β; but for the imaginary part, to which the leakage loss is related, a larger computational window size or a thicker PML thickness is required to obtain the converged magnitude. The ideal R value for calculating the imaginary part of β is found to be much smaller than the typical ones and would be as small as, for example, 10−50 or 10−100. With such a small R value, both the PML thickness and the computational window size can be significantly reduced to minimum degrees to obtain the converged magnitude without loss of stability. Moreover, it is shown that in some cases the modal field profiles can only be accurately obtained with such a small R value
2. The case of holey fiber
We first investigate the holey fiber with a ring of six air holes symmetrically arranged around the core. The structure is shown in Fig. 1 , with the hole diameter D=5 μm, the hole pitch Λ=6.75 μm, and the refractive indices ns=1.45 and na=1. The computational window size is denoted as W and the PML thickness is denoted as d. For the FDFD simulation, the grid sizes Δx=Δy=0.05 μm and λ=1.45 μm are assumed.
We calculate the complex effective index neff, defined as neff=β/k 0, where k 0 is the free-space wavenumber. Calculated results of Re(neff) and Im(neff) of the fundamental mode are shown in Figs. 2(a) and 2(b), respectively, as a function of W with a fixed d=1 μm. Obviously, results using R=10−8 are the least stable because they vary and oscillate with W most seriously, and the oscillation amplitude becomes smaller when W increases. However, results of Re(neff) are consistent over the whole simulation range of W with up to the first eight digits (1.4453975), while those of Im(neff) differ even in the first digit. Hence the typical value R=10−8 might be acceptable to calculate Re(neff), but is inadequate for Im(neff). The latter is so highly affected by the computational window size that it is difficult to decide what size of W should be adopted to obtain the reasonable converged magnitude. (Of course a very large computational window size far above those shown in Fig. 2 might lead to the converged result, but it is unpractical if the limited computing resource is taken into account.) Results using R=10−12 show similar characteristics but vary with a less significant oscillation. In general, as the value of R decreases, the variation and oscillation become gradually diminished. In the case of R=10−50, it is clear both Re(neff) and Im(neff) are extremely consistent over the whole W range. Note that in Fig. 2, W starts from 20 μm, which is very close to the required minimum computation window size 18.5μm in this case (2Λ + D, the distance between the outmost boundaries of the two air holes located in the central horizontal line). In other words, if a sufficiently small value of R like 10−50 is used for simulation, only a small W near the minimum size requirement would be enough to obtain the converged magnitude of neff, which significantly saves the computation cost in terms of memory and time.
Fig. 2 Calculated effective index neff of the holey fiber as a function of W with d=1 μm. (a) Real part. (b) Imaginary part.
Now we fix W=24 μm and then change d to calculate neff. Results are shown in Figs. 3(a) and 3(b) for Re(neff) and Im(neff), respectively. Again, the results using R=10−8 exhibit the worst performance, and the variation and oscillation get diminished as R decreases. In the case of R=10−50, both Re(neff) and Im(neff) are consistent over the simulated range of d. Note that for R=10−50, a very thin PML thickness d=0.2 μm, which is equivalent to only 4 grids in PML, is sufficient to obtain the converged magnitude. Therefore, the computation cost can be saved when such a small R value is used for simulation.
Fig. 3 Calculated effective index neff of the holey fiber as a function of d with W=24 μm. (a) Real part. (b) Imaginary part.
Figure 4 shows Im(neff) as a function of R, where W and d are fixed to 24 μm and 1 μm, respectively. It is observed that Im(neff) will converge to a stable magnitude when R is smaller than some threshold value, and in this case of holey fiber it is about 10−50. The calculated result of neff with R=10−50 is 1.445397590 − j3.23×10−8 under the adopted grid size of 0.05 μm (it is 1.445395151 − j3.1965×10−8 when the grid size decreases to 0.0125 μm), which is close to 1.445395345 − j3.15×10−8 obtained by the multipole method [13] and 1.445395232 − j3.1945×10−8 by the integral equation method [16]. Note that the degree of accuracy depends on the inherent property of the numerical scheme used for simulation, and is out of the scope of this paper. (Interested readers might refer to [17], which compares several mode solvers by using a photonic wire structure as a common example. As a reference, our result to that example is 2.412290 − j2.9202×10−8 with the parameters R=10−50 and the grid size of 0.005 μm.) What we aim to emphasize here is the effect of R, as discussed above and in the following, which is illustrated in this work using the finite difference method. This effect would be general to other numerical techniques as long as the PML is employed for leaky-mode analysis.
3. The case of air-core pipe waveguide
The second leaky structure is the air-core pipe waveguide recently proposed for THz waveguiding [14]. It consists of a large air core and a thin low-index dielectric cladding, as shown in Fig. 5 , and its guiding mechanism is similar to that of the antiresonant reflecting optical waveguide (ARROW) [18,19].
By assuming the core diameter D=9 mm, the cladding thickness t=1 mm, the refractive indices n 1=1 and n 2=1.4, Δx=Δy=0.025 mm, and the operating frequency being 840 GHz, we calculate neff of the fundamental mode by changing W and d, respectively. Since Re(neff) is usually much easier to converge, only the results of Im(neff) are shown in Figs. 6(a) and 6(b), respectively. Calculated results exhibit the same characteristics as that of the previously discussed holey fiber. For larger R such as 10−10, significant variation and oscillation occur, and the oscillation becomes smaller as W or d increases; therefore, a larger W or a thicker d is required to obtain the converged Im(neff). While for R=10−100, calculated Im(neff)’s are stable over the simulation range, so that both W or d can be significantly reduced to obtain the converged magnitude without loss of stability. From Fig. 7 , which shows Im(neff) as a function of R, it is clear Im(neff) converges if R is smaller than the threshold value of 10−100.
Fig. 6 (a) Im(neff) of the pipe waveguide as a function of W with d=0.4 mm. (b) Im(neff) of the pipe waveguide as a function of d with W=14 mm.
4. The case of gain-guided and index-antiguided slab waveguide
Previous two examples illustrate the effect of R on calculating the complex effective indices (or propagation constants) of the leaky waveguide modes. Now the third example is to demonstrate the effect on solving the leaky field distributions. Being with a high potential for large mode area laser applications, the GG + IAG waveguide [15,20] has a large gain in the core and the core has a lower refractive index than that of the cladding. Here we consider the GG + IAG slab waveguide with a core width of 100 μm, a cladding index of 1.5, and a complex core index of 1.4998 + jg/2k 0, where g is the small signal power gain coefficient. By assuming W=800 μm, d=8 μm, Δx=0.5 μm, and λ=1.55 μm, we first calculate the Im(neff) of the fundamental transverse-electric (TE) mode as a function of R for g=0 (IAG case) and g=1 cm−1 (GG + IAG case), and the results are shown in Figs. 8(a) and 8(b), respectively. Note that Im(neff) for the IAG case is negative, indicating the IAG waveguide mode is leaky (or radiative); but that for the GG + IAG case is positive, which is due to the gain in the core region. From Fig. 8, it is interesting to find that Im(neff) converges as R gets small, and the threshold value is about 10−200 in either case.
Fig. 8 Im(neff) of the gain-guided and index-antiguided slab waveguide as a function of R. (a) g=0 (IAG case). (b) g=1 cm−1 (GG + IAG case).
We then examine the modal field profiles. Calculated modal field distributions for g=0 and g=1 cm−1 are shown in Figs. 9(a) and 9(b), respectively. Clearly, for both cases, using the typical R value 10−10 is unable to obtain the correct modal profiles. R must be sufficiently small so that the stable and converged modal fields could be solved, even though there is gain in the waveguide core. Here, R should be at most 10−200.
Fig. 9 Fundamental TE modal field profiles of the gain-guided and index-antiguided slab waveguide. (a) g=0 (IAG case). (b) g=1 cm−1 (GG + IAG case).
5. Conclusion
In conclusion, we have investigated the effect of the reflection coefficient of PML on the modal analysis of leaky waveguide modes. By examining three different leaky structures, it reveals that the typical values of the reflection coefficient 10−8 ~10−12 are inadequate for leaky-mode analysis, and the satisfactory coefficient should be much smaller. Using such a small coefficient, both the computational window size and the PML thickness can be significantly reduced to minimum degrees to obtain the converged complex propagation constant. Moreover, in some cases, the modal field profiles can only be accurately obtained with such a small coefficient. Note that the threshold value of the reflection coefficient is case dependent and could be as small as 10−100 or 10−200; thus we suggest assigning the coefficient as small as possible under the constraint of the computer software. The software we used in this study is the MATLAB, and the ARPACK program which is accessible via the function eigs() is employed to solve the complex eigenvalue problems. The smallest R value acceptable in this study is about 10−320, and applying such a small R value to computation does not increase additional CPU time and memory space. Finally, although our simulation is based on the finite difference scheme, it is believed that the conclusion would be also applicable to other numerical techniques such as the finite element method, the pseudospectral method, etc.
Acknowledgement
This work was supported in part by the National Science Council of the Republic of China under grant NSC97-2221-E-002-043-MY2, in part by the Excellent Research Projects of National Taiwan University under grant 98R0062-07, and in part by the Ministry of Education of the Republic of China under “The Aim of Top University Plan” grant.
References and links
1. G. R. Hadley, “Transparent boundary condition for beam propagation,” Opt. Lett. 16(9), 624–626 (1991). [CrossRef] [PubMed]
2. J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J. Comput. Phys. 114(2), 185–200 (1994). [CrossRef]
3. W. C. Chew and W. H. Weedon, “A 3-D perfectly matched medium from modified Maxwell’s equation with stretched coordinates,” Microw. Opt. Technol. Lett. 7(13), 599–604 (1994). [CrossRef]
4. C. M. Rappaport, ““Perfectly matched absorbing boundary conditions based on anisotropic lossy mapping of space,”IEEE Microw. Guid. Wave Lett. 5(3), 90–92 (1995). [CrossRef]
5. A. Taflove, and S. C. Hagness, Computational Electrodynamics: The Finite-Difference Time-Domain Method, 3rd ed. (Artech House, 2005).
6. W. P. Huang, C. L. Xu, W. Lui, and K. Yokoyama, “The perfectly matched layer boundary condition for modal analysis of optical waveguides: leaky mode calculations,” IEEE Photon. Technol. Lett. 8(5), 652–654 (1996). [CrossRef]
7. Y. Tsuji and M. Koshiba, “Guided-mode and leaky-mode analysis by imaginary distance beam propagation method based on finite element scheme,” J. Lightwave Technol. 18(4), 618–623 (2000). [CrossRef]
8. C.-P. Yu and H.-C. Chang, “Yee-mesh-based finite difference eigenmode solver with PML absorbing boundary conditions for optical waveguides and photonic crystal fibers,” Opt. Express 12(25), 6165–6177 (2004), http://www.opticsinfobase.org/oe/abstract.cfm?URI=OPEX-12-25-6165. [CrossRef] [PubMed]
9. W.-P. Huang and J. Mu, “Complex coupled-mode theory for optical waveguides,” Opt. Express 17(21), 19134–19152 (2009), http://www.opticsinfobase.org/oe/abstract.cfm?URI=OPEX-17-21-19134. [CrossRef]
10. P. L. Ho and Y. Y. Lu, “A mode-preserving perfectly matched layer for optical waveguides,” IEEE Photon. Technol. Lett. 15(9), 1234–1236 (2003). [CrossRef]
11. B. Chen, D. G. Fang, and B. H. Zhou, “Modified Berenger PML absorbing boundary condition for FD-TD meshes,” IEEE Microw. Guid. Wave Lett. 5(11), 399–401 (1995). [CrossRef]
12. P.-J. Chiang and Y.-C. Chiang, “Pseudospectral frequency-domain formulae based on modified perfectly matched layers for calculating both guided and leaky modes,” IEEE Photon. Technol. Lett. 22(12), 908–910 (2010). [CrossRef]
13. T. P. White, B. T. Kuhlmey, R. C. McPhedran, D. Maystre, G. Renversez, C. M. de Sterke, and L. C. Botten, “Multipole method for microstructured optical fibers. I. Formulation,” J. Opt. Soc. Am. B 19(10), 2322–2330 (2002). [CrossRef]
14. C.-H. Lai, Y.-C. Hsueh, H.-W. Chen, Y.-J. Huang, H.-C. Chang, and C.-K. Sun, “Low-index terahertz pipe waveguides,” Opt. Lett. 34(21), 3457–3459 (2009). [CrossRef] [PubMed]
15. T.-H. Her, X. Ao, and L. W. Casperson, “Gain saturation in gain-guided slab waveguides with large-index antiguiding,” Opt. Lett. 34(16), 2411–2413 (2009). [CrossRef] [PubMed]
16. H. Cheng, W. Y. Crutchfield, M. Doery, and L. Greengard, “Fast, accurate integral equation methods for the analysis of photonic crystal fibers I: Theory,” Opt. Express 12(16), 3791–3805 (2004), http://www.opticsinfobase.org/oe/abstract.cfm?URI=OPEX-12-16-3791. [CrossRef] [PubMed]
17. P. Bienstman, S. Selleri, L. Rosa, H. P. Uranus, W. C. L. Hopman, R. Costa, A. Melloni, L. C. Andreani, J. P. Hugonin, P. Lalanne, D. Pinto, S. S. A. Obayya, M. Dems, and K. Panajotov, “Modelling leaky photonic wires: a mode solver comparison,” Opt. Quantum Electron. 38(9−11), 731–759 (2006). [CrossRef]
18. M. A. Duguay, Y. Kokubun, T. L. Koch, and L. Pfeiffer, “Antiresonant reflecting optical waveguides in SiO2–Si multilayer structures,” Appl. Phys. Lett. 49(1), 13–15 (1986). [CrossRef]
19. C.-H. Lai, B. You, J.-Y. Lu, T.-A. Liu, J.-L. Peng, C.-K. Sun, and H.-C. Chang, “Modal characteristics of antiresonant reflecting pipe waveguides for terahertz waveguiding,” Opt. Express 18(1), 309–322 (2010), http://www.opticsinfobase.org/oe/abstract.cfm?URI=OPEX-18-1-309. [CrossRef] [PubMed]
20. A. E. Siegman, “Propagating modes in gain-guided optical fibers,” J. Opt. Soc. Am. A 20(8), 1617–1628 (2003). [CrossRef]
