## Abstract

The pseudospectral method, proposed in our previous work, has not yet been constructed for optical waveguides with leaky modes or anisotropic materials. Our present study focuses on antiresonant reflecting optical waveguides (ARROWS) made by anisotropic materials. In contrast to the fields in the outermost subdomain expanded by Laguerre-Gaussian functions for guided mode problems, the fields in the high-index outermost subdomain are expanded by the Chebyshev polynomials with Mur’s absorbing boundary condition (ABC). Accordingly, the traveling waves can leak freely out of the computational window, and the desirable properties of the pseudospectral scheme, i.e., provision of fast and accurate solutions, can be preserved. A number of numerical examples tested by the present approach are shown to be in good agreement with exact data and published results achieved by other numerical methods.

© 2006 Optical Society of America

## 1. Introduction

Building optical integrated circuits for practical use in an optical communication system requires a variety of photonic components such as optical interconnectors, modulators, switches, filters, splitters, polarizers, and couplers. The diverse optical waveguides by using distinct materials and refractive index profiles, according to specifications, produce functional optical integrated circuits. To precisely compute propagation characteristics (i.e., dispersion, attenuation, mode distribution, and birefringence) of general optical waveguides with both guided and leaky modes, we may precisely tailor the expected optical devices. A low-loss waveguide antiresonant reflecting optical waveguide (ARROW), which utilizes two thin interference cladding films made on semiconductor substrates [1], has been proposed to construct a wide range of functional devices [2–5]. Unlike total internal reflection mechanisms for guided modes, ARROW structure supports leaky waves through high reflectivity corresponding to the antiresonant condition of a Fabry-Perot resonator. This kind of waveguide has several advantages [6,7] such as low loss, relatively large spot size for effectively coupling into optical fibers, high fabrication tolerance for thicknesses of interference claddings, high loss discrimination between the fundamental and higher-order modes, high polarization sensitivity, and polarization insensitivity [8].

In the viewpoint of computation, ARROW-type devices may be regarded as multilayer waveguide structures, and thus various analytical or numerical approaches including transverse resonance method (TRM) [7,9] and transfer matrix method (TMM) [10–12] have been developed to find the dispersion equations. To calculate propagation characteristics, except some traditional zero-search algorithms like Newton’s and scan methods, which cannot find all of the guided and leaky modes simultaneously, a more efficient and reliable technique is proposed for extracting the poles of dispersion equations in complex plane by the argument principle method (APM) [13,14]. However, the TRM and TMM are only accurate enough to be used with planar waveguides that have step-index profiles rather than graded-index profiles, because inhomogeneous materials are often approximated by a stack of dielectric layers with constant refractive indices. Moreover, the two schemes mentioned above are also unsuitable for the 2-D cross-section structures like rib ARROW providing also lateral confinement on mode profiles. Accordingly, it is necessary and desirable to pay close attention to the rigorous numerical solution methods. The schemes used most often such as the finite difference modal analysis [15], finite element modal analysis [16,17], and imaginary distance beam propagation method [18], have been built up to tackle various leaky waveguides.

Recently, an efficient and accurate numerical method, based on a number of orthogonal basis functions termed as the pseudospectral method [19], has been proposed to solve dielectric optical waveguides supporting only the guided modes [20–22]. The purpose of this paper is to extend the ability of the pseudospectral scheme to the ARROW structures with leaky modes. At the present time, to modify the consideration representing outermost subdomains by the Laguerre Gaussian (LG) functions for guided mode problems with the exponential decay feature in our previous work [20–22], the outermost subdomains occupied by oscillatory leaky waves are represented by the Chebyshev polynomials with Mur’s absorbing boundary condition (ABC) [23]. Accordingly, the oscillatory behaviors of leaky waves extending to infinity can be well represented, and the superior performances of pseudospectral scheme for accuracy and efficiency may also be preserved. A number of isotropic and anisotropic ARROW structures are analyzed to verify the performances of the proposed solution method. Other issues of the paper are organized as follows. In Section 2, we formulate the wave equations for an anisotropic medium. The description of the pseudospectral scheme with Mur’s ABC is stated in Section 3. Section 4 simulates several numerical examples to validate the accuracy and efficiency of the present approach by comparing exact data. Finally, Section 5 concludes this work.

## 2. Formulations of planar waveguides

Considering an anisotropic medium if the principal dielectric axes of the crystal are found to be parallel to the waveguide coordinate system, the dielectric tensor *ε̿* is expressed as follows:

where *ε*
_{0} is the dielectric constant in free space, and *n*_{xx}
, *n*_{yy}
, and *n*_{zz}
are the refractive indices along the *x*, *y*, and *z* axes of the waveguide coordinate system, respectively. In this paper, we aim at analyzing leaky modes in planar ARROW of arbitrary index profiles with a diagonal dielectric tensor shown in expression (1). Assuming that *ε̿* is only a function of the independent *x* (*i.e.*, planar waveguides with finite and infinite extension in *x* and *y* directions, respectively), for the source-free, time-harmonic of exp(*iωt*) form, and nonmagnetic media, the modal equations in frequency domain derived from Maxwell’s two curl equations are given as follows:

for transverse-electric (TE) and that

for transverse-magnetic (TM) polarizations. Here, *k*
_{0}=2*π*/*λ*
_{0}, *λ*
_{0} denotes the wavelength in free space, *n*_{eff}
denotes the complex effective refractive index, and *E*_{y}
(*x*) and *H*_{y}
(*x*) represent, respectively, the components of the electric and magnetic fields in the *y* -direction.

## 3. Numerical scheme

#### 3.1 Pseudospectral scheme and interfacial patching conditions

The framework of the pseudospectral scheme is to divide computational domain into a number of subdomains determined by discontinuous material interfaces. In each subdomain, the optical field *φ*(*x*) is expanded by a set of proper interpolation functions *C*_{k}
(*x*) and the unknown grid point values *φ*_{k}
as follows:

where

Here, *ψ*
_{n+1}(*x*)denotes an arbitrary basis function of order *n*+1 and is determined according to the feature of that subdomain, the prime denotes the first derivative of *ψ*
_{n+1}(*x*) with respect to *x*, and *x*_{k}
denotes the corresponding collocation point fulfilling the condition of *C*_{l}
(*x*_{k}
)=*δ*_{lk}
, where *δ*_{lk}
denotes the Kronecker delta. The explicit forms of *C*_{k}
(*x*) for distinct basis functions are represented below. If we take the Chebyshev polynomials as basis functions, we have [19]

where *T*_{n}
(*x*) denotes the Chebyshev polynomial of order *n*, *c*
_{0}=*c*_{n}
=2, and *c*_{k}
=1 (1≤*k*≤*n*-1). For the LG functions, namely *ψ*
_{n+1}(*αx*)=exp(-*αx*/2)(*αx*)*L*_{n}
(*αx*) where*L*_{n}
(*αx*)denotes the Laguerre polynomial of order *n*, the explicit form of*C*_{k}
(*x*)is given as follows [19]:

Here, in LG functions, exists an important parameter *α* called as a scaling factor that can influence the accuracy for a given number of terms of basis functions. The definite procedure for determining *α* has been proposed in detail in our previous work [20, 21]. From Eqs. (2) and (3), we have an eigenvalue problem in the form

where **I** denotes a unit matrix, **Φ** denotes the vector of unknown optical field, and **A** denotes the transverse differential operator. For obtaining the entries of the transverse operator to the *e* th subdomain, we substitute Eq. (4) with a suitable basis function into Eqs. (2) and (3) and obtain

for TE and that

for TM polarizations. In Eqs. (9) and (10), ${C}_{k}^{\left(n\right)}$(*x*_{l}
) means the *n* th derivative result of *C*_{k}
(*x*) with respect to *x* at collocation point *x*_{l}
. The global transverse operator **A** can be obtained by assembling the entries of the transverse operators to all subdomains. However, in matrix **A**, the rows located in the dielectric interfaces *x*_{r}
({*r*=1, 2,…*m*}, where *m* denotes the number of discontinuous nodes) shared by two adjacent subdomains are imposed to be replaced by the interfacial boundary conditions

for TE and that

for TM polarizations, where the ${x}_{r}^{+}$ and ${x}_{r}^{-}$ are, respectively, referred to the locations at the infinitesimal right and left of the interface *x*_{r}
.

#### 3.2 Determining basis functions and boundary conditions

The best choice is that the interior and outermost subdomains are expanded, respectively, by the Chebyshev polynomials and LG functions, if only guided modes exist [20–22]. However, the outermost subdomains represented by the LG functions are inappropriate while tackling the waveguide structures with leaky modes. In this work, we mainly propose an efficient solution method for handling leaky waveguides. In general, the optical fields in the outermost subdomains exhibit two conditions: one is with evanescent wave, and the other is with traveling wave. Certainly, the outermost subdomains with the evanescent decay characteristic can be efficiently represented by the LG functions as well as those considered in Refs. [20–22]. As for that with traveling field, since the oscillatory behavior resulting from leaky modes penetrating to substrate, which have a larger refractive index than that in the core layer, this kind of subdomain is expanded by the Chebyshev polynomials at collocation points but has to be fulfilled by the Mur’s ABC [23] at the boundary point (*i. e.*, outmost collocation point), which allows the fields to radiate freely out of the computational window. The perfectly matched layer (PML) condition may be a more efficient scheme than Mur’s ABC, but the implementation of PML to the present method is very complicated. In this paper, we mainly aim to validate the accuracy and efficiency of the pseudospectral approach to leaky mode problems. As a result, we prefer choosing a simpler boundary condition like Mur’s ABC to accomplish the work.

In this work, we adopt the second-order formula of Mur’s schemes in frequency domain, and it can be given as follows:

where
$i=\sqrt{-1}$
, *φ* denotes the unknown field and is expanded by the Chebyshev polynomials, *r* denotes the outward direction normal to the boundary, and *k*_{r}
denotes the wavenumber in *r*-direction. In the paper, *r* and *k*_{r}
can be replaced by *x* and *k*_{x}
in the present coordinate system, respectively. Further, Eq. (13) can be expanded as follows:

where *x*_{bp}
denotes the boundary point in the outermost subdomain. Note that the chosen computational interval of the outermost subdomain expanded by the Chebyshev polynomials slightly affects the convergence. Observing Eq. (14), we find that the plane wave solution of *φ*(*x*)=*φ*(0)exp(-*ik*_{x}*x*)is fulfilled [17]. After that, we substitute the plane wave solution into Eqs. (2) and (3) for the outermost subdomain with homogeneous refractive index; we have

for TE, and

for TM polarizations. Through Eqs. (2) and (3), the first term in Eq. (14) can be replaced by the terms of -${k}_{x}^{2}$
*φ* for both TE and TM polarizations. As a result, Eq. (14) can be further modified as follows:

for TE and TM polarizations. For ensuring only the outgoing plane wave, we must restrict that the real part of *k*_{x}
is positive. In numerical implementation, in order to keep the unit matrix **I** in Eq. (8), the term (${k}_{0}^{2}$
${n}_{\mathit{\text{eff}}}^{2}$)*φ* is added simultaneously in two sides of Eq. (17). Therefore a matrix eigenvalue problem is obtained by evading a generalized matrix eigenvalue problem. Since the introduction of the Mur’s ABC, the eigenvalue *n*_{eff}
resides within the global matrix and is required to be solved by iteration scheme. The computational procedure is that we first assign the initial guess of ${n}_{\mathit{\text{eff}}}^{0}$ to **A**, and then a new value of ${n}_{\mathit{\text{eff}}}^{1}$ is found by solving a general linear matrix eigenvalue. Subsequently, ${n}_{\mathit{\text{eff}}}^{1}$ is used as a new value to calculate the next one ${n}_{\mathit{\text{eff}}}^{2}$. We do not accept the convergent results until the differences of real and imaginary part of *n*_{eff}
between the adjacent iterations are simultaneously less than the order of 10^{-6}.

Besides the consideration of choosing basis functions in the present work, the Chebyshev polynomials with Mur’s ABC can also be used to represent the evanescent field having decaying characteristics at the boundary as well as that used in the finite element scheme with Sommerfeld’s ABC [17]. However, in our numerical experiments concerning computational efficiency (the required number of basis functions), the preferable choice is still the LG functions. Moreover, no supplementary boundary conditions are demanded for evanescent fields [20–22], and thus the numerical implementation may be significantly simplified.

## 4. Numerical results and discussion

In order to verify the performances of the present scheme as applied to various ARROW structures, we first analyze in detail the propagation characteristics of an isotropic example, which has been studied by Baba and Kokubun [6] using the system interference matrix and by Chen *et al.* [14] using TMM with APM. Except for the anisotropy, the second example is the same as the first example. Finally, we study the coupling length and field distributions of even and odd modes for an ARROW-based directional coupler device. The accuracy of these examples calculated by the proposed approach is discussed and compared with other published results.

#### 4.1 Isotropic ARROW structure

The schematic diagrams of the geometry and refractive index profile of the ARROW structure are depicted in Fig. 1. The relevant parameters used are those with an operating wavelength of λ=0.6328µm, the refractive indices are n_{a}=1for air and n_{s}=3.85 for substrate, and the thicknesses (refractive indices) for the first cladding, second cladding, and core layers are, respectively, d_{1}=0.142λ (n_{1}=2.3), d_{2}=3.15λ (n_{2}=1.46), and d_{c}=6.3λ (n_{c}=1.46) under the antiresonant condition. The computational domain is divided into five subdomains. The calculated results of the complex effective indices for the first *i* time iterations ${\mathrm{n}}_{\text{eff}}^{\mathrm{i}}$ are shown in Table 1, which also shows that ${\mathrm{n}}_{\text{eff}}^{0}$=1.459 denotes the chosen initial guess of effective index. In the example, the number of terms of LG basis functions for the outermost subdomain occupied by air is N_{a}=10 while the computational interval estimated according to our criterion [21] for spreading of guided mode is *S*_{a}
=1*µm*, the number of Chebyshev polynomials for the first cladding, second cladding, and core layers are N_{1}=N_{2}=N_{c}=20, and the Chebyshev polynomials for substrate layer is N_{s}=40 while the computational interval of the outermost subdomain chosen is *S*_{sub}
=1*µm*. More terms of basis functions in the substrate layer are necessary due to the sharply oscillatory behavior of leaky waves. In addition, considering the effect caused by the computational interval of the outermost subdomain with higher refractive index, the differences of the convergent results between *S*_{sub}
=1*µm* and *S*_{sub}
=2*µm* for all of the modes are negligible while using N_{s}=50 for *S*_{sub}
=2*µm*. We conclude that the major factor in obtaining accurate results is use of the pseudospectral scheme. Table 1 also shows that with only three iterations, our method achieves the convergent values. In addition, the exact four-digit values for n_{eff} can be obtained even when only one iteration is executed. In our numerical experiments, the present approach obtains the same convergent rate for using different ${\mathrm{n}}_{\text{eff}}^{0}$, as long as the chosen value of ${\mathrm{n}}_{\text{eff}}^{0}$ is slightly smaller than n_{c}. The solutions of this work, obtained by the third iteration, and the results calculated by Chen *et al.* [14] are shown in Table 2. We can see that the two schemes are in excellent agreement.

The initial guess of effective index is ${\mathrm{n}}_{\text{eff}}^{0}$=1.459, and the number of terms of basis function in each subdomain is N_{a}=10, N_{1}=N_{2}=N_{c}=20, and N_{s}=40.

To further demonstrate the accuracy of our method, we calculate the dispersion and radiation loss characteristics of the ARROW structure versus the thickness of each layer while the thicknesses of other layers are at corresponding antiresonant conditions. First, Figs. 2(a) and 2(b) show the dispersion and loss characteristics of the first six ARROW modes of TE polarization, labeled TE_{1}-TE_{6}, versus the thickness of the first cladding (d_{1}/λ), respectively, while other parameters are fixed. Here, referring to Ref. [6], we call TE_{0} as the first cladding mode. In Fig. 2(a), the ranges of d_{1}/λ with flat portions of effective index for each mode correspond to the antiresonances of Fabry-Perot cavities formed in the two interference claddings. The same portions represent the low-loss portions while referring to Fig. 2(b). It is clear that the low-loss portion of TE1 is broad so as to allow large fabrication tolerance. In contrast to the flat portions, the transitional portions in Fig. 2(a) illustrate the resonances subject to high-loss portions in Fig. 2(b). For comparing the order of losses between different polarized waves, Fig. 2(b) also illustrates the first two ARROW modes of TM polarization. It is clear that the losses of TM modes are larger than TE ones at a corresponding mode number. As a result, the higher-order TE and TM modes can be easily filtered out to obtain single-mode (TE_{1}) propagation. Besides, in Fig. 2(b), the loss characteristics display periodic variety as a function of (d_{1}/λ). Likewise, the propagation characteristics versus the thickness of the second cladding (d_{2}/λ) and core (d_{c}/λ) layers are shown in Figs. 3 and 4, respectively. In Figs. 3(b) and 4(b), at corresponding first antiresonant conditions, high loss discrimination as that observed in Fig. 2(b) is also shown.

Furthermore, Figs. 5(a)–5(f) show the real part of relative field profiles of various TE (solid curve) and TM (dotted curve) polarized waves simultaneously at the first antiresonant condition (*i. e.*, d_{1}=0.142λ, d_{2}=3.15λ, and d_{c}=6.3λ). Here, “relative” means that the field profiles are normalized by their maximum values [6]. In this paper, Figs. 5(a)–5(f) are sequentially labeled as TE_{0} (TM_{0}) to TE_{6} (TM_{6}). We can observe that the fields in Fig. 5(a) are confined to the vicinity of the first cladding layer, and the confinement of the TE mode is tighter than that of the TM mode. Figure 5(b) illustrates that the lowest ARROW mode (TE_{1}) is most localized in the core layer, which is of practical interest, and that it is similar to the fundamental mode of conventional waveguides. The oscillatory characteristics in the substrate layer are also clearly shown, and the larger loss of TM mode is observed, as well as the values listed in Table 2. From Figs. 5(b)–5(f), excluding Figs. 5(b) and 5(e), the fields penetrating the interference claddings are fairly large.

#### 4.2 Anisotropic ARROW structure

Following the same structure analyzed in Subsection 4.1, but with anisotropy in inner layers sandwiched between air and substrate layers, the anisotropic material with a diagonal form of dielectric tensor is expressed in Eq. (1). The refractive indices are n_{zzc}=1.46 for the core, n_{zz1}=2.3for the first cladding, and n_{zz2}=1.46 for the second cladding. Other values follow the relations of n_{xxj}=n_{yyj}=1.03n_{zzj}, where (j=c, 1,2). Here, the initial guess of effective index of ${\mathrm{n}}_{\text{eff}}^{0}$=1.459 is chosen, and three iterations are executed to achieve convergent solutions. The results obtained by our method, using 140 unknowns (N_{a}=10,N_{1}=N_{2}=N_{c}=30,N_{s}=40), are shown in Table 3 with the six-order accuracy of the finite element method [17] and TMM using APM [14]. The number of unknowns in Ref. [17] is nearly 1300 because of the effective mesh size of 0.005µm. This shows that the computational efficiency of the present scheme is superior to that of the finite element method with six-order accuracy [17]. In addition, the results of complex effective indices obtained by our work show excellent agreement with exact solutions [14]. However, the present scheme obtains a full matrix. The computational time is 0.18 sec (each iteration) executed on Pentium IV PC with a CPU clock rate of 3.0 GHz in a double precision. We obtain that the effective indices and losses as a function of thickness of each inner layer are similar to its isotropic counterpart. As for that the principal dielectric axes of the crystal are not parallel to the waveguide coordinate system, the more complicated hybrid modes and the propagation characteristics depend strongly on the included angle between the optics axis and the waveguide coordinate system.

#### 4.3 ARROW-based directional coupler

The directional coupler, which consists of two parallel waveguides, is an interesting photonic structure to build various optical communication devices such as a modulator, filter, power divider, and polarizer. The coupling mechanism of a conventional directional coupler is well studied by the coupled mode theory [24], which considers merely the weakly coupled approximation. In addition, the coupling strength decreases as an exponential function of waveguide separation. Because of the widespread applications of ARROW-based structures, a novel dual ARROW-based directional coupler utilizing antiresonant reflection as the guiding mechanism has been proposed to investigate the propagation characteristics [25]. In comparison with a conventional directional coupler that utilizes total internal reflection, the ARROW structures have strong coupling, which is due to their having leaky waves rather than the evanescent waves found in conventional two parallel waveguides; furthermore, the variation of the coupling length versus waveguide separation displays periodicity [25].

In Fig. 6, the configuration of the dual ARROW structure grown on a Si substrate with a refractive index n_{s}=3.5 consists of two ARROW structures separated by a separation cladding layer with a refractive index n_{sep}=1.46 and thickness d_{sep}=2µm [26]. At the operating wavelength λ=0.6328µm, the two waveguide cores are with n_{c1}=n_{c2}=1.46 and d_{c1}=d_{c2}=4µm. The core 1 is sandwiched by the two thin cladding layers with n_{h1}=n_{h2}=2.3 and d_{h1}=d_{h2}=0.089µm, and the same circumstance is also encountered by core 2, which is sandwiched by two cladding layers with n_{h3}=n_{h4}=2.3and d_{h3}=d_{h4}=0.089µm. For the upper cladding layer, which is grown over the thin cladding layer of d_{h1} and in contact with air having a n_{a}=1, the parameters are n_{11}=1.46 and d_{11}=2µm. As for the lower cladding layer grown over the Si substrate, it is with n_{12}=1.46 and d_{12}=2µm. The purpose of the different interference claddings is to accomplish the antiresonant conditions. According to the conclusion found by Chen and Huang [26], the coupling efficiency can be flexibly controlled through the adjustment of outermost cladding thickness d_{11}. If we consider the conditions of d_{11}=d_{12}=2µm and d_{11}=0,d_{h1}=0.03µm, the maximum coupling efficiency of C_{o}=99.87% (while the coupling length is L_{c}=λ/(N_{s}-N_{a})=59mm, where N_{s} and N_{a} denote lowest order symmetric (even) and asymmetric (odd) modes, respectively) and decoupled phenomenon can be achieved [26], respectively.

To validate the characteristics studied in Ref. [26] and demonstrate the ability of our scheme, by the present scheme, the dual ARROW structure is considered as an 11-layer waveguide structure, and the computational window is divided into 11 subdomains. The numbers of terms of LG functions used for air is N_{a}=10, of Chebyshev polynomials for substrate is N_{s}=40, and of others are N_{i}=20 (where *i* represents other layers). Compared with the coupling length of L_{c}=59mm, obtained in Ref. [26], the calculated coupling length is L_{c}=58.60mm (where N_{s}=1.45785857 and N_{a}=1.45785317). The coupling length of TM mode obtained by our scheme is L_{c}=9.51mm (where N_{s}=1.45787059 and N_{a}=1.45783732), which is far smaller than that of the TE mode. This is because the coupling resonances of leaky waves of TM modes are larger than that of TE modes. In our calculation, the imaginary parts (radiation loss) of complex effective indices of TM modes are two orders larger than that of TE modes (namely, TE:≈10^{-6} dB/cm, TM:≈10^{-4} dB/cm). From the dispersion characteristics of TE modes illustrated in Fig. 7 it can be clearly seen that the maximum coupling is located at d_{11}=2µm, making the dual ARROW acts as a symmetric coupler, and the decoupling portions are close to d_{11}=0µm and d_{11}=4µm. Additionally, the field profiles of the even and odd TE modes for three cases of d_{11}=2µm (maximum coupling), d_{11}=0.5µm (half coupling), and d_{11}=0 (decoupling, while d_{h1}=0.03*µ*m is used [26]) are clearly shown by the decreased order of coupling in Figs. 8(a)–8(c), respectively. Finally, in Fig. 9, the coupling characteristics of the TE and TM modes as a function of the waveguide separation are also illustrated. Obviously, the coupling lengths of the TE and TM modes appear to be a periodic function, as predicted in Ref. [25]. Consequently, this shows that we may apply a dual ARROW-based coupler to easily accomplish a variety of applications in a remote coupler [25].

## 5. Conclusion

We successfully demonstrate the ability of the present method to deal with the leaky modes supported by the isotropic and anisotropic ARROW structures. The main consideration is to represent the outermost subdomains with strongly oscillatory characteristics by the efficient Chebyshev polynomials with Mur’s ABC allowing outgoing waves to freely penetrate into substrate, evading most reflected waves. As a result, the oscillatory characteristics are accurately represented. The results calculated by the present scheme show excellent agreement with exact solutions and the computational efficiency is far higher than the finite element method. The application of our scheme to a dual ARROW structure is considered as an ARROW-based directional coupler; the accurate values of the coupling length are obtained, and periodicity of the coupling length versus waveguide separation is also validated. The extension of the present method to three-dimensional ARROW structures of practical interests will be reported elsewhere.

## Acknowledgments

This work was supported by the National Science Council, Taiwan, Republic of China under contract No. NSC 95-2221-E-275-005.

## References and links

**1. **M. A. Duguay, Y. Kokubun, T. L. Koch, and L. Pfeiffer, “Antiresonant reflecting optical waveguides in SiO_{2}-Si multilayer structures,” Appl. Phys. Lett. **49**, 13–15 (1986). [CrossRef]

**2. **T. Baba, Y. Kokubun, T. Sakaki, and K. Iga, “Loss reduction of an ARROW waveguide in shorter wavelength and its stack configuration,” J. Lightwave Technol. **6**, 1440–1445 (1988). [CrossRef]

**3. **M. Mann, U. Trutschel, C. Wachter, L. Leine, and F. Lederer, “Directional coupler based on antiresonant reflecting optical waveguide,” Opt. Lett. **16**, 805–807 (1991). [CrossRef]

**4. **Z. M. Mao and W. P. Huang, “An ARROW optical wavelength filter: design and analysis,” J. Lightwave Technol. **11**, 1183–1188 (1992). [CrossRef]

**5. **F. Prieto, A. Llobera, D. Jimenez, C. Domenguez, A. Calle, and L. M. Lechuga, “Design and analysis of silicon antiresonant reflecting optical waveguides for evanescent field sensor,” J. Lightwave Technol. **18**, 966–972 (2000). [CrossRef]

**6. **T. Baba and Y. Kokubun, “Dispersion and radiation loss characteristics of antiresonant reflecting optical waveguides-numerical results and analytical expressions,” IEEE J. Quantum Electron. **28**, 1689–1700 (1992). [CrossRef]

**7. **W. P. Huang, R. M. Shubair, A. Nathan, and Y. L. Chow, “The modal characteristics of ARROW structures,” J. Lightwave Technol. **10**, 1015–1022 (1992). [CrossRef]

**8. **T. Baba and Y. Kokubun, “New polarization-insensitive antiresonant reflecting optical waveguide (ARROW-B),” IEEE Photon. Technol. Lett. **1**, 232–234 (1989). [CrossRef]

**9. **W. Jiang, J. Chrostowski, and M. Fontaine, “Analysis of ARROW waveguides,” Opt. Commun. **72**, 180–186 (1989). [CrossRef]

**10. **J. Chilwell and I. Hodgkinson, “Thin-film field-transfer matrix theory for planar multilayer waveguides and reflection from prism-loaded waveguides,” J. Opt. Soc. Am. A **1**, 742–753 (1984). [CrossRef]

**11. **J. Kubica, D. Uttamchandani, and B. Culshaw, “Modal propagation within ARROW waveguides,” Opt. Commun. **78**, 133–136 (1990). [CrossRef]

**12. **J. M. Kubica, “Numerical analysis of InP/InGaAsP ARROW waveguides using transfer matrix approach,” J. Lightwave Technol. **10**, 767–771 (1992). [CrossRef]

**13. **E. Anemogiannis and E. N. Glytsis, “Multilayer waveguides: efficient numerical analysis of general structures,” J. Lightwave Technol. **10**, 1344–1351 (1992). [CrossRef]

**14. **C. K. Chen, P. Berini, D. Feng, S. Tanev, and V. P. Tzolov, “Efficient and accurate numerical analysis of multilayer planar waveguides in lossy anisotropic media,” Opt. Express **7**, 260–272 (2000). [CrossRef]

**15. **W. P. Huang, C. L. Xu, W. Lui, and K. Yokoyama, “The perfectly matched layer boundary conditions for modal analysis of optical waveguides: leaky mode calculations,” IEEE Photon. Technol. Lett. **8**, 652–654 (1996). [CrossRef]

**16. **J. C. Grant, J. C. Beal, and N. J. P. Frenette, “Finite element analysis of the ARROW leaky optical waveguide,” IEEE J. Quantum Electron. **30**, 1250–1253 (1994). [CrossRef]

**17. **H. P. Uranus, H. J. W. M. Hoekstra, and E. V. Groesen, “Simple high-order Galerkin finite scheme for the investigation of both guided and leaky modes in anisotropic planar waveguides,” Opt. Quantum Electron. **36**, 239–257 (2004). [CrossRef]

**18. **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**, 618–623 (2000). [CrossRef]

**19. **J. P. Boyd, “Chebyshev and Fourier Spectral methods,” in *Lecture Notes in Engineering*, 2nd ed. (Springer Verlag, 2001).

**20. **C. C. Huang, C. C. Huang, and J. Y. Yang, “An efficient method for computing optical waveguides with discontinuous refractive index profiles using spectral collocation method with domain decomposition,” J. Lightwave Technol. **21**, 2284–2296 (2003). [CrossRef]

**21. **C. C. Huang, C. C. Huang, and J. Y. Yang, “A full-vectorial pseudospectral modal analysis of dielectric optical waveguides with stepped refractive index profiles,” IEEE J. Sel. Top. Quantum Electron. **11**, 457–465 (2005). [CrossRef]

**22. **C. C. Huang and C. C. Huang, “An efficient and accurate semivectorial spectral collocation method for analyzing polarized modes of rib waveguides,” J. Lightwave Technol. **23**, 2309–2317 (2005). [CrossRef]

**23. **G. Mur, “Absorbing boundary conditions for the finite difference approximation of the time-domain electromagnetic field equations,” IEEE Trans. Electromagn. Compat. **23**, 377–382 (1981). [CrossRef]

**24. **A. Hardy and W. Streifer, “Coupled mode theory of parallel waveguides,” J. Lightwave Technol. **3**, 1135–1146 (1985). [CrossRef]

**25. **M. Mann, U. Trutschel, C. Wachter, L. Leine, and F. Lederer, “Directional coupler based on an antiresonant reflecting optical waveguide,” Opt. Lett. **16**, 805–807 (1991). [CrossRef]

**26. **Y. H. Chen and Y. T. Huang, “Coupling-efficiency analysis and control of dual antiresonant reflecting optical waveguides,” J. Lightwave Technol. **14**, 1507–1513 (1996). [CrossRef]