OSA's Digital Library

Optics Express

Optics Express

  • Editor: J. H. Eberly
  • Vol. 9, Iss. 9 — Oct. 22, 2001
  • pp: 461–475
« Show journal navigation

Propagation in planar waveguides and the effects of wall roughness

J. Merle Elson  »View Author Affiliations


Optics Express, Vol. 9, Issue 9, pp. 461-475 (2001)
http://dx.doi.org/10.1364/OE.9.000461


View Full Text Article

Acrobat PDF (551 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

We consider the propagation of guided waves in a planar waveguide that has smooth walls except for a finite length segment that has random roughness. Maxwell’s equations are solved in the frequency domain for both TE and TM polarization in 2-D by using modal expansion methods. Obtaining numerical solutions is facilitated by using perfectly matched boundary layers and the R-matrix propagator.Varying lengths of roughness segments are considered and numerical results are obtained for guided wave propagation losses due to roughness induced scattering. The roughness on each waveguide boundary is numerically generated from an assumed Gaussian power spectrum. The guided waves are excited by a Gaussian beam incident on the waveguide aperture. Considerable numerical effort is given to determine the stability of the algorithm.

© Optical Society of America

1 Introduction

The work presented here has been motivated by the development of waveguide components for use in fiber optic gyros. The fiber optic gyro requires that two counter rotating guided waves be the same single polarization and mode. These “clean” guided modes are created by a specially designed waveguide that supports and allows only a single mode and polarization to emerge from the output side. The output intensity must be sufficient and this is where loss mechanisms are important. Since the waveguide has wall roughness, radiative scattering losses can occur and the question is how much roughness can be tolerated. Another possibility is that roughness scattering is secondary to scattering due to bulk material inhomogeneity of the waveguide core. In this work scattering from material inhomogeneity is not considered, rather we calculate radiative losses per unit length of wall roughness.

The method used in this work is flexible so that many planar waveguide problems can be considered with either metallic or dielectric walls. This includes step discontinuities, periodic or random roughness, multiple waveguide configurations, multiplexing, beam combining, etc. Using this method, it is straightforward to consider the case of additional scattering due to waveguide core inhomogeneity and this work is in progress.

We numerically solve Maxwell’s equations in the frequency domain for propagation in planar waveguides. Waveguide modes are initiated by a Gaussian beam incident on the waveguide aperture. This is illustrated by the schematic in Fig. 1. Part of the Gaussian incident energy that is not reflected is coupled into one or more guided wave modes that propagate down a smooth portion of the waveguide that extends from z=L t to z=L r. At the end of the smooth segment begins a segment of roughness that has length extending from z=L t to z=0. The length of the smooth portion is chosen so that at z=L r, the only remaining electromagnetic energy propagating downward is a guided wave mode. While transiting the roughness region, some of the guided wave energy is lost to scattering and after exiting the roughness region at z=0, the remaining guided wave mode energy continues to propagate unimpeded in the semi-infinite smooth substrate region (z<0) of the waveguide. The 1-D roughness is numerically generated with a Gaussian correlation function.

Many authors have investigated waveguide scattering losses [8

8. D. Marcuse, “Mode conversion caused by surface imperfections in a dielectric slab waveguide,” Bell Sys. Tech. J. 48, 3187–3215 (1969).

]–[15

15. J. Rodríguez, R. D. Crespo, S. Fernández, J. Pandavenes, J. Olivares, S. Carrasco, I. Ibáñez, and J. Virgós, “Radiation losses on discontinuities in integrated optical waveguides,” Opt. Engr. 381896–1906 (1999). [CrossRef]

]. Early theoretical work by Marcuse [8

8. D. Marcuse, “Mode conversion caused by surface imperfections in a dielectric slab waveguide,” Bell Sys. Tech. J. 48, 3187–3215 (1969).

, 9

9. D. Marcuse, “Radiation losses of dielectric waveguides in terms of the power spectrum of the wall distortion function,” Bell Sys. Tech. J. 48, 3233–3242 (1969).

] considered the coupling between guided and radiative wave modes. Based on these works, Payne and Lacey [10

10. F. P. Payne and J. P. R. Lacey, “A theoretical analysis of scattering loss from planar optical waveguides,” Opt. and Quantum Elec. 26977–986 (1994). [CrossRef]

, 11

11. J. P. R. Lacey and F. P. Payne, “Radiation loss from planar waveguides with random wall imperfections,” IEE Proc. J. 137282–288 (1990).

] also provided a theoretical analysis of scattering losses from planar optical waveguides. Their analytical results are dependent on numerous waveguide characteristics including the power spectrum of wall roughness that is available to allow coupling into the radiative scattering field. Some experimental results are given in [12

12. K. K. Lee, D. R. Lim, H. Luan, A. Agarwal, J. Foresi, and L. Kimerling, “Effect of size and roughness on light transmission on a Si/SiO2 waveguide: Experiment and model,” Appl. Phys. Lett. 771617–1619 (2000). [CrossRef]

]–[14

14. F. Ladouceur, J. D. Love, and T. J. Senden, “Measurement of surface roughness in buried channel waveguides,” Electron. Lett. 281321–1322 (1992). [CrossRef]

] and these data are also compared with theory similar to that given in [10

10. F. P. Payne and J. P. R. Lacey, “A theoretical analysis of scattering loss from planar optical waveguides,” Opt. and Quantum Elec. 26977–986 (1994). [CrossRef]

, 11

11. J. P. R. Lacey and F. P. Payne, “Radiation loss from planar waveguides with random wall imperfections,” IEE Proc. J. 137282–288 (1990).

].

Fig. 1. Schematic of waveguide with wall roughness. The waveguide channel, having nominal width ω and permittivity 2, is bounded by media with permittivity 1. In the z-direction, the waveguide channel consists of a smooth segment (L rzL t), a roughness segment (0≤zL r), and a semi-infinite substrate region z≤0. The x-dimension of the computational region is bounded by xL x/2 and the shaded regions at the extreme left and right denote the PML absorption layers.

2 Theory

2.1 Maxwell’s Equations

We begin with Maxwell’s equations for the electric E⃗ and magnetic B⃗ fields in the frequency domain.

×E(x,y,z)=i(ω/c)μ(x,y,z)H(x,y,z)
(1)
×H(x,y,z)=i(ω/c)(x,y,z)E(x,y,z)
(2)

Reducing Eqs. (1) and (2) to two dimensions and incorporating the PML features [6

6. J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J. Comp. Phys. 114185–200 (1995). [CrossRef]

, 7

7. J. P. Berenger, “Three-dimensional perfectly matched layer for the absorption of electromagnetic waves,” J. Comp. Phys. 127363–379 (1996). [CrossRef]

] yields the following coupled first order equations that can be written in real space as

i(ω/c)Ex(x,z)=αz(x,z)Hy(x,z)z
(3)
i(ω/c)εx(x,z)Ez(x,z)=Hy(x,z)x
(4)
i(ω/c)Hy(x,z)=βx(x,z)Ez(x,z)xβz(x,z)Ex(x,z)z
(5)
i(ω/c)Hx(x,z)=βz(x,z)Ey(x,z)z
(6)
i(ω/c)μx(x,z)Hz(x,z)=Ey(x,z)x
(7)
i(ω/c)Ey(x,z)=αx(x,z)Hz(x,z)xαz(x,z)Hx(x,z)z
(8)

For reasons that will become clear, Eqs. (3)(8) have been specifically written in the form as shown and this will be discussed in more detail below. We have defined the quantities

μj(x,z)=1+4πiσj*(x,z)/ω;εj(x,z)=(x,z)μj(x,z)
(9)
βj(x,z)=1/μj(x,z);αj(x,z)=1/εj(x,z)
(10)

where j=x or z. These position-dependent functions describe the details of the waveguide structure by their dependence on the scalar permittivity (x, z), the magnetic conductivity σj(x, z)* and the electric conductivity σj(x, z). The PML conductivity terms, σj(x, z) and σj(x, z)*, control absorption for propagation of the electric and magnetic fields in the j-direction, respectively. The electric conductivity σj(x, z) is not explicit in Eqs. (9) and (10) because the relationship σj(x, z)=∊(x, z)σj* (x, z) is required for impedance matching. Ideally, impedance matching eliminates reflection when the PML absorbing layers are encountered. Outside of a PML absorption region, the σj=σj*=0, and then Eqs. (3)(8) revert to the usual form of Maxwell’s equations.

Equations Eqs. (3)(5) can be obtained from equations Eqs. (6)(8) (and vice versa) by making the replacements

HE;EH
(11)
(μx,μz)(εx,εz);(εx,εz)(μx,μz)
(12)

For this reason, we will limit the remaining theory discussion to TM polarization only. The case for TE polarization is obtained with Eqs. (11) and (12).

We now convert Eqs. (3)(5) to a k-space representation and make the transition to a finite sized problem. When these conversions are properly done, equations Eqs. (3)(5) assume a form that has been shown to have excellent numerical convergence properties [16

16. P. Lalanne and G. Morris, “Highly improved convergence of the coupled-wave method for TM polarization,” J. Opt. Soc. Am. 13779–784 (1996). [CrossRef]

, 17

17. L. Li, “Use of Fourier series in the analysis of discontinuous periodic structures,” J. Opt. Soc. Am. 131870–1876 (1996). [CrossRef]

]. This is in contrast to directly applying the same numerical approach as described here to equations Eqs. (3)(8) where convergence and numerical stability are very poor. Before making the conversion to k-space, we make two approximations to Eqs. (3)(5). First, we truncate the x-dimension of the solution space to a finite length L x and discretize the x-coordinate. Second, we assume that the material parameters in Eq. (9) and (10) are z-invariant.

We discretize the x-coordinate into N points as x m=mΔx where integer m=-N/2 → N/2 and Δx=L x/N. Likewise, we discretize the x-component of the wavevector k n=nΔk where n=-N/2 → N/2 and Δk=2π/L x. Now that we have discretized and truncated the problem, we use matrix-vector notation to represent the Fourier transform and inverse Fourier transform: F and F -1. Such transform relations are written as

f(x,z)=F(x,k)f(k,z);f(k,z)=F1(x,k)f(x,z)
(13)

where the F is a N×N matrix and f is a N-element column vector. Note that these are one dimensional transforms applied to the k x-x transform pair where k xk. Applying these transform operators to Eqs. (3)(5) yields

i(ω/c)Ex(k,z)=αz(k,k)Hy(k,z)z
(14)
(ω/c)εx(k,k)Ez(k,z)=kHy(k,z)
(15)
i(ω/c)Hy(k,z)=iβx(k,k)kEz(k,z)+βz(k,k)Ex(k,z)z
(16)

where k is interpreted as a diagonal matrix with elements k n and the H y, E x, and E z are N-element column vectors with each element corresponding to k n. We have defined the N×N square matrices

μj(k,k)=F1(x,k)μj(x)F(x,k);εj(k,k)=F1(x,k)εj(x)F(x,k)
(17)
βj(k,k)=F1(x,k)βj(x)F(x,k);αj(k,k)=F1(x,k)αj(x)F(x,k)
(18)

where the quantities µ j(x), ε j(x), β j(x), and αj(x) are interpreted as diagonal matrices with elements µ j(x n), etc.

Using Eqs. (14)(16), we obtain the second-order differential equation

2Hy(k,z)z2=αz1βz1[βxkεx1kI(ω/c)2]Hy(k,z)=MHy(k,z)
(19)

where I is the identity matrix. Since there are N discrete k values, Eq. (19) represents N coupled differential equations. When M is independent of z (this is the case when the material parameters are z-invariant), solution of this equation system is straightforward by diagonalization as S -1 MS=ξ 2=ξ·ξ. The columns of the N×N matrix S are the eigenvectors of M and the diagonal matrix ξ 2 contains the eigenvalues of M. The solution to Eq. (19) is

Hy(k,z)=S(eξzC++eξzC)
(20)

and the electric field follows from Eq. (14) as

Ex(k,z)=i(c/ω)αzSξ(eξzC+eξzC)
(21)

The C ± are scalar constants and the e ±ξz represent diagonal matrices where the n-th diagonal element is an exponential term as exp(±ξ n z). The ξ n is the n-th diagonal root-eigenvalue element of ξ. In Eqs. (20) and (21) we now have solutions to Maxwell’s equations for an inhomogeneous solution space that has an x-dimension of length L x and has material parameters that are z-invariant.

2.2 R-matrix Propagator

In order to consider solution spaces that are not z-invariant, we assume that we can build up such structures from numerous z-invariant layers. Thus, for all the z-invariant layers in the computational domain, we need to relate the solutions for each layer as given in Eqs. (20) and (21). A numerically stable relationship is given by the R-matrix algorithm which is but one of a class of several types of transfer matrices. First, for solutions within a single z-invariant layer bounded from zz+δ, a relationship is assumed to be given by

(Ex(k,z)Ex(k,z+δ))=(r11(δ)r12(δ)r21(δ)r22(δ))(Hy(k,z)Hy(k,z+δ))
(22)

where δ is the thickness of the layer. Inserting Eqs. (20) and (21) into Eq. (22) yields the N×N r-matrices as

r11(δ)=r22(δ)=(ic/ω)αzSξ(eξδ+eξδ)(eξδeξδ)1S1
(23)
r12(δ)=r21(δ)=(2ic/ω)αzSξ(eξδeξδ)1S1
(24)

These r-matrices relate the electric and magnetic fields at the boundaries across a single z-invariant layer. Generalizing to the case of two or more contiguous layers bounded from zz+z t, we use the R-matrices which have a similar form as

(Ex(k,z)Ex(k,z+zt))=(R11(zt)R12(zt)R21(zt)R22(zt))(Hy(k,z)Hy(k,z+zt))
(25)

The parameter z t is the cumulative thickness of all the z-invariant layers as described in Eq. (22). The N×N R-matrices are calculated recursively by

R11(zt)=R11(ztδ)+R12(ztδ)[r11(δ)R22(ztδ)]1R21(ztδ)
(26)
R12(zt)=−R12(ztδ)[r11(δ)R22(ztδ)]1r12(δ)
(27)
R21(zt)=r21(δ)[r11(δ)R22(ztδ)]1R21(ztδ)
(28)
R22(zt)=r22(δ)r21(δ)[r11(δ)R22(ztδ)]1r12(δ)
(29)

2.3 Superstrate and Substrate

Referring to Fig. 1, we want to relate solutions within 0≤zL t to solutions in the semi-infinite superstrate zL t and substrate z≤0 and we can use Eq. (25) to provide such a relation. Setting z=0 and z t=L r in Eq. (25), we can relate to the superstrate and substrate fields by invoking the continuity of the tangential components. If we have incident and reflected fields in the superstrate and a transmitted field in the substrate, then these continuity conditions are

Ex(k,0)=Ext(k,0);Hy(k,0)=Hyt(k,0)
(30)
Ex(k,Lt)=Exr(k,Lt)+Exi(k,Lt);Hy(k,Lt)=Hyr(k,Lt)+Hyi(k,Lt)
(31)

where the superscripts i, r, and t denote incident, reflected, and transmitted, respectively.

In Fig. 1 we see that the superstrate and substrate regions are z-invariant. In addition, the superstrate region zL t is homogeneous. In the superstrate region, we let a Gaussian beam be incident on the waveguide aperture and our goal will be to calculate the resulting fields at the z=0 and z=L t boundaries. This will allow us to further calculate the guided wave transmitted field after it has emerged from the roughness region and propagates into z<0. Although we do not do so here, we can also calculate the reflected field for z>L t.

2.3.1 Transmitted field z≤0

In this section, we discuss the Ext(k, z) and Hyt(k, z) transmitted fields. Solutions of the form given by Eqs. (20) and (21) are still valid for the transmitted region by setting C-=0 and this yields

Hyt(k,z)=SeξzC+Hyt(k,z)=SeξzS1Hyt(k,0)
(32)
Ext(k,z)=i(c/ω)αzSξeξzC+Ext(k,z)=i(c/ω)αzSξeξzS1Hyt(k,0)
(33)

Since the region z≤0 has only downward propagating waves, the eigenvalues associated with Eqs. (32) and (33) must have ℜξ≥0 and ℑξ≤0. The latter equations in Eqs. (32) and (33) give the z-dependence after the guided wave exits the roughness region. These equations can be examined numerically after we have obtained the solution for Hyt(k, 0). Equations Eqs. (32) and (33) may be combined to give the relationship between the electric and magnetic transmitted fields in the inhomogeneous region z≤0 which is

Ext(k,z)=i(c/ω)αzSξS1Hyt(k,z)=THyt(k,z)
(34)

and this will be used later.

2.3.2 Incident and Reflected Fields: zLt

In the homogeneous superstrate region, the Gaussian incident field is assumed to be a superposition of plane waves [18

18. A. A. Maradudin, T. Michel, A. R. McGurn, and E. Mendez, “Enhanced backscattering of light from a random grating,” Annals of Physics 203225–307 (1990). [CrossRef]

] with half-width determined by ωσ/c and the angle of incidence centered about θ i. We write the x and y-components of the incident field as

(Exi(x,z)Hyi(x,z))=ωσ2cππ/2π/2dθ(cosθ1)
×exp[(ωσ2c)2(θθi)2]exp[i(ω/c)(xsinθzcosθ)]
(35)

where the superstrate medium is assumed to be vacuum. As in Eqs. (30) and (31), we want the Fourier transforms of these fields and Fourier inversion of these components yields the column vectors Exi(k, z) and Hyi(k, z) where the n-th element is given by

Exi(kn,z)=σπexp[(ωσ2c)2(θnθi)2]exp[i(ω/c)zcosθn]
(36)
Hyi(kn,z)=Exi(kn,z)/cosθn
(37)

when k n<(ω/c) and zero when k n>(ω/c). The wavevector k n=2πn/L x and sin θ n=k n/(θ/c).

The reflected fields are written as a superposition of plane waves as

(Exr(x,z)Hyr(x,z))=12πdk(Exr(k)Hyr(k))exp[i(kx+qz)]
=12πdk(Exr(k,z)Hyr(k,z))exp(ikx)
(38)

where the z-component of the wavevector is q=(ω/c)2k2. In the latter equation, the z-dependence has been absorbed into the Fourier coefficient. It is straightforward to show that Exr(k, z)=(qc/ω)Hyr(k, z). With this and consistent with the notation of Eq. (31), we relate the column vectors associated with the reflected field as

Exr(k,z)=Z(k)Hyr(k,z)
(39)

where the n-th element of the diagonal matrix Z(k) is

Z(kn)=(ω/c)2kn2(ω/c)
(40)

2.3.3 Surface Fields: z=0 and z=Lt

We now have the relationships to calculate the fields at the upper and lower boundaries of the computational domain as shown in Fig. 1. Using Eqs. (30), (31), (34), (39), and (25) yields

(TR11R12R21ZR22)(Hyt(k,0)Hyr(k,Lt))=(R12Hyi(k,Lt)R22Hyi(k,Lt)Exi(k,Lt))
(41)

and this linear equation system may be solved for the surface fields Hyt(k, 0) and Hyr(k,L t). Along with the eigenvectors S and eigenvalues ξ for the substrate region, the first of these two solutions may be then used in Eqs. (32) and (33) or Eq. (34) to obtain the transmitted fields for z<0.

3 Numerical Results

The theory discussed in Section 2i s applied to the problem of guided wave propagation in rough waveguides. Numerical results of transmission losses due to the wall roughness are given in Figs. 35. In Figs. 68 we present some results showing numerical convergence characteristics of the algorithm used in this work. In Fig. 9, we illustrate how effective the PML layers can be in preventing reflection from the hard limits of the computational domain. Finally, in Fig. 10, the transmission curves are fit to analytical model to estimate power loss per unit length of wall roughness.

3.1 General Comments

There are some parameters common to the numerical data and reference to Fig. 1 and Fig. 2 may be helpful in parts of the following discussion. All realizations for the waveguide wall roughness were numerically generated for a Gaussian autocorrelation function with a correlation length of 1.0λ and various values of rms roughness: 0.05λ, 0.1λ and 0.2λ. Four realizations have been used and numerical results for each realization are displayed separately rather than an attempt at ensemble averaging. Except in Fig. 10, the initial smooth portion of the waveguide has length L t-L r=3000λ and nominal waveguide channel width ω=1λ. This length L t-L r is sufficient to allow the remaining transmitted energy that reaches the z=L r plane to only be one or more guided wave modes. With this, only a guided wave mode is incident on the roughness region. The length L r of the roughness region is a parameter of study in this work that varies over many values. In the numerical algorithm, the length L r=N r δ is the cumulative thickness of N r z-invariant layers each having thickness δ. This is illustrated in Fig. 2. Within each layer of thickness δ, the waveguide channel width is dependent on the roughness at each waveguide boundary. We have let δ=0.1λ for all calculations presented here. The waveguide permittivity values are ∊1=(2.25, 0.) and ∊2=(2.50, 0.). In most cases, we assume that the Gaussian beam is normally incident (θ i=0) and the ratio σ/λ=5 (see Eq. (35)). In Figs. 35, the extent of the computational region in the x-dimension is L x=16.3λ and this length is divided into N=299 segments of length Δx=L x/N=0.055λ.

Fig. 2. Schematic showing modeling of waveguide roughness region. The waveguide channel has nominal width ω. For this example, the roughness region is divided into N r=5 sublayers of thickness δ and the channel width within each sublayer is varied in accordance with the roughness. By doing this, each sublayer is z-invariant and solutions of the form given in Eqs. (34) and (35) apply throughout the waveguide. In the numerical results of this work, the number of sublayers ranged from 0 to 1600, δ=0.1λ, and 0λ≤L r≤160λ. The permittivity of the light-shaded and dark-shaded regions are ∊1=(2.25, 0.) and ∊2=(2.50, 0.), respectively.

In the theory described in Section 2, some of the PML absorption parameters have a z-subscript that identifies them as pertaining to absorption for propagation in the z-direction. However, in the problem described here, absorption for propagation in the z-direction is not applicable and this means that we set σz*=0 everywhere. However, in order to prevent reflection and simulate infinite extent in the x-direction, absorption of fields having propagation components in the x-direction is applicable. As shown in Fig. 1, there are two PML absorption regions when |x| in the vicinity of L x/2eac h region has a total thickness denoted by δ PML. If η is the magnitude of the penetration depth into the PML region, where η ranges over 0 → δ PML, we assume that the PML absorption increases quadratically with penetration depth as

4πσx*(η)/ω=APML(η/δPML)2
(42)

which will be used in Eq. (9). The PML regions are divided into N PML layers each with thickness Δx and the total thickness may then be written as δ PML=N PMLΔx. In all cases except the data shown in Fig. 9, we have set N PML=24 and the value of A PML=8. Since the x coordinate is discretized, we calculate the N PML values of Eq. (42) with η set to the middle value of the appropriate layer in a PML region.

3.2 Transmission Loss

The normalized transmitted power is calculated as the ratio of P t/P i where

Pi=n=N/2N/2Exi(kn,Lt)[Hyi(kn,Lt)]*;Pt=n=N/2N/2Ext(kn,0)[Hyt(kn,0)]*
(43)

Shown in Figs. 35 is transmitted power lost as L r increases. These three plots consider both TE and TM polarization and each plot includes four roughness realizations. The difference between the three plots is the rms roughness used to generate the realizations where the rms roughness values are 0.05λ, 0.10λ, and 0.20λ, respectively. The transmitted power is only plotted for every 10δ=1λ increase in L r. Thus, each curve in Figs. 35 contain 160 data points. However, since the numerical algorithm actually calculates in increments of δ=0.1λ, the total number of sublayers over the roughness region ranges from 0 to 1600 over the range of L r. To include the smooth region L rzL t, we use four sublayers with δ=750λ for a total thickness L t-L r=3000λ. As would be expected, it is seen in Figs. 35 that the rate of scattered power lost as L r increases becomes greater as the rms roughness increases. These curves are summarized later in Section 4.

Fig. 3. Normalized transmitted power versus length of roughness calculated at z=-10000λ. There are four groups of curves where each group has four curves. In each group the four curves, labeled 1–4, correspond to a different roughness realization. The upper and lower two sets of curves are for TE and TM polarization, respectively. The upper two sets of data have been displaced upward 0.4 units for clarity. The horizontal curves are for the normalized reflected power. The downward sloping curves are the normalized transmitted power. The number of discretization points N=299 and the nominal waveguide channel width ω=1λ. All realizations are generated with rms roughness 0.05λ and correlation length 1λ.
Fig. 4. Normalized transmitted power versus length of roughness calculated at z=-10000λ. These data are analogous to those in Fig. 3 except that all realizations are generated with rms roughness 0.10λ.
Fig. 5. Normalized transmitted power versus length of roughness calculated at z=-10000λ. These data are analogous to those in Fig. 3 except that all realizations are generated with rms roughness 0.20λ.

3.3 Convergence

We now look at the stability of the numerical solution as certain parameters are varied. In Fig. 6, the number of digitization points parameter N is varied with all other parameters held constant. Figure 6 is similar to Fig. 3 through Fig. 5 where the residual transmitted power versus L r is shown, but only for realization 1. The main purpose of Fig. 6 is to show convergence as the number of digitization points N increases for constant L x=16.3λ. The data in Figs. 35 were generated for N=299 whereas in Fig. 6, the seven TM curves cover N=199 → 449 and the five TE cover N=199 → 399. These data indicate numerical convergence as N increases, or equivalently, as Δx decreases.

In Fig. 7 convergence is now considered for the case where the increment Δx is held fixed. We choose L x and N such that L x/Nx≈0.055λ. Values of L r range from 14.3λ to 24.3λ and it is seen that as L r increases, the transmission at z=-10000λ is noticeably greater for the larger values of N and L x. The short answer for this discrepancy is that for larger L x, the plane wave solutions are absorbed by the PML layers at a lower rate as z → -∞. As an example, when L x=16.3λ one of the 299 calculated complex eigenvalues is ξ=(ω/c)(0.9340×10-4-1.5i). When L x=2 4.3λ there are 441 calculated eigenvalues and the corresponding eigenvalue is now ξ=(ω/c)(0.1217×10-4,-1.5i). In the transmission region, the solutions are proportional to exp(ξz) (see Eqs. (32) and (33)) and it is clear that these particular eigenvalues pertain to the downward propagating plane wave solution in the =2.25 bounding medium of the waveguide structure. The absorptive real part of ξ is generated by the presence of the PML layers. Note that at z=-10000λ and the two ξ values above, the products ℜ(ξ)z are -0.9340(ω/c) and -0.1217(ω/c). For these values it is clear that exp(ξz) is still contributing a non-negligible solution value to the transmitted field. In hindsight, this discrepancy could have been avoided simply by choosing z<<-10000λ such that ℜ(ξ)z<<0 and all plane wave solutions would have been sufficiently decayed to be a negligible contribution to the transmitted field. We conclude that the discrepancy noted in Fig. 7 is not a convergence problem.

Fig. 6. Convergence of normalized transmitted power for various values of discretization points N. For the same roughness realization 1 as used in Figs. 35, the transmitted power is calculated at z=-10000λ versus length of roughness. The upper set of curves, which have for clarity been displaced upward by 0.2, are for TM polarization and the values of N vary from 199 to 449. The lower set of curves are for TE polarization and the values of N vary from 199 to 399. For the corresponding polarization, the curves are essentially converged except for the N=199 curves. The rms roughness is 0.2λ.
Fig. 7. Convergence of normalized TM transmitted power for various values of discretization points N and length L x. The ratio of L x/Nx is adjusted so that Δx≈0.055λ is nearly constant. For the same roughness realization 1 that was used in Figs. 35, the transmitted power is calculated at z=-10000λ versus length of roughness and the rms roughness is 0.2λ. Note that as L x increases, there is an increase in transmitted power due to a decrease in the absorption rate for the plane wave solutions. This is not a convergence related issue.

Shown in Fig. 8 are convergence characteristics as the number of PML layers N PML and the absorption coefficient A PML are varied. For curves labeled N PML=1 through 24, we set A PML=8 and vary N PML as indicated. For the curve labeled N PML=0, we have for numerical reasons actually set N PML=1 with A PML=0.1. This latter case does not completely eliminate PML absorption, but the effect is relatively minimal. We see that the N PML=0 curve shows little absorption at z=-10000λ and examination of the plane wave eigenvalues as was done in Fig. 7 shows a very low rate of absorption coefficient. As the PML absorption is increased, the curves quickly converge to a common result.

Fig. 8. Normalized TM transmitted power for various values of PML absorption layers N PML. For 1≤N PML≤24, the value A PML=8. For numerical reasons, the results shown by the curve labeled N PML=0 was actually calculated with N PML=1 and A PML=0.1 (see Eq. (42)). Thus, the N PML=0 curve only approximates the absence of PML layers. The incident beam is at normal incidence, L x=16.3λ, N=299, and rms roughness is 0.2λ. These data are for the same roughness realization 1 that was used in Figs. 35. The group of horizontal lines are the normalized reflected power. It is clear that if there is sufficient PML absorption, then the solutions converge to a common result.

In Fig. 10, the TM curves in Figs. 35 are averaged and modeled by fitting a curve of the form y(x)=m 010m1x where y is the normalized transmission and x is L r/λ. The m 0 and m 1 are fitting parameters. In Figs. 35, the fitting parameters are m 0=0.293, 0.286, and 0.277, and m 1=-0.00045, -0.00168, and -0.00542, respectively. These analytical curves are plotted in Fig. 10 and the loss can be expressed in dB as -10Log10[y(x)/y(0)] which yields 0.0045L r/λ, 0.0168L r/λ, and 0.0542L r/λ, respectively.

4 Concluding Remarks

We have presented a frequency domain calculation method and numerical results that predict scattering losses in planar waveguides due to wall roughness. The wall roughness is numerically generated assuming a Gaussian autocorrelation function with root mean square values of 0.05λ, 0.10λ, and 0.20λ and a correlation length of 1.0λ. We have also considered several aspects regarding convergence of the analytical/numerical method where various numerical parameters are varied.

Fig. 9. TM E-field intensity profile at various depths from the waveguide opening. There is no roughness region in this example and the 8 indicated z values are relative to the waveguide input aperture plane. The width L x=18.3λ and the angle of incidence is θ i=30°. For clarity, all curves, except the z=-205λ curve, have been displaced by multiples of 0.01 units. The z=0λ curve is the electric field intensity across the waveguide input aperture plane and the remaining z-value curves show the progression of the intensity curves at indicated depths into the waveguide. It is clear that the refracted part of the incident Gaussian profile is ultimately absorbed by the right-hand set of PML layers without reflection. In addition, the final z=-205λ curve is beginning to take shape as the residual guided wave since much of the transmitted energy has disappeared. The two vertical lines at x=±0.5λ indicate the boundaries of the waveguide channel. The two vertical lines at x=±7.65λ indicate the start of the PML absorbing layers.

The numerical trend of the transmission energy loss results obtained here are intuitive and entirely predictable: the rougher the random surface and the longer its length, the greater the guided wave propagation losses. However, the actual numerical amount of loss per unit length of roughness is not predictable and this can be an important question that depends on many parameters. For example, suppose the intrinsic material absorption is negligible and the measured scattering losses are unacceptably large. It may not be obvious as to what the dominant source or sources of scattering losses are. Should the focus be on reducing the wall roughness or decreasing the inhomogeneity of the waveguide core? These kinds of questions can apply to many electro-optic devices and theoretical/analytical tools can, among other things, help interpret experimental results and aid in development/research decisions. Regarding scattering due to waveguide core inhomogeneity, the work presented here can, with some minor modifications, be used to predict losses resulting from this source.

The objective of this paper is part methodological as well. The method is sufficiently general that many other interesting problems can be examined, such as photonic crystal structure, for example. The main criteria is that the structure of interest be described by spatially varying linear media. With the use of PML absorbing layers, the computational domain is finite and not dependent on periodic boundary conditions. Reflection from the “hard limits” at the truncation boundary of the computational domain is prevented by the PML algorithm which is commonly used in finite difference time domain algorithms. Of course, the method is easily adaptable to be used for infinitely periodic structure as well. The R-matrix algorithm allows very large computational domains to be considered with numerical stability.

Fig. 10. Analytical curves derived from the data shown in Figs. 35. The four curves in each of these figures are each compiled into one analytical curve as a power of 10. This yields an estimate of the dB loss as a function of the length of roughness.

References and links

1.

J. M. Elson and P. Tran, “R-matrix propagator with perfectly matched layers for the study of integrated optical components,” J. Opt. Soc. Am. A 162983–2989 (1999). [CrossRef]

2.

J. M. Elson and P. Tran, “Coupled-mode calculation with the R-matrix propagator for the dispersion of surface waves on a truncated photonic crystal,” Phys. Rev. B 541711–1715 (1996). [CrossRef]

3.

J. M. Elson and P. Tran, “Band structure and transmission of photonic media: a real-space finite-difference calculation with the R-matrix propagator,” NATO ASI Series E: Applied Sciences on Photonic Band Gap MaterialsVol. 315 341–354 Crete, Greece June 15–29 (1995).

4.

L. Li, “Multilayer modal method for diffraction gratings of arbitrary profile, depth, and permittivity,” J. Opt. Soc. Am. A 112816–2828 (1993). [CrossRef]

5.

L. Li, “Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings,” J. Opt. Soc. Am. 131024–1035 (1996). [CrossRef]

6.

J. P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,” J. Comp. Phys. 114185–200 (1995). [CrossRef]

7.

J. P. Berenger, “Three-dimensional perfectly matched layer for the absorption of electromagnetic waves,” J. Comp. Phys. 127363–379 (1996). [CrossRef]

8.

D. Marcuse, “Mode conversion caused by surface imperfections in a dielectric slab waveguide,” Bell Sys. Tech. J. 48, 3187–3215 (1969).

9.

D. Marcuse, “Radiation losses of dielectric waveguides in terms of the power spectrum of the wall distortion function,” Bell Sys. Tech. J. 48, 3233–3242 (1969).

10.

F. P. Payne and J. P. R. Lacey, “A theoretical analysis of scattering loss from planar optical waveguides,” Opt. and Quantum Elec. 26977–986 (1994). [CrossRef]

11.

J. P. R. Lacey and F. P. Payne, “Radiation loss from planar waveguides with random wall imperfections,” IEE Proc. J. 137282–288 (1990).

12.

K. K. Lee, D. R. Lim, H. Luan, A. Agarwal, J. Foresi, and L. Kimerling, “Effect of size and roughness on light transmission on a Si/SiO2 waveguide: Experiment and model,” Appl. Phys. Lett. 771617–1619 (2000). [CrossRef]

13.

F. Ladouceur, J. D. Love, and T. J. Senden, “Effect of side wall roughness in buried channel waveguides,” IEE Proc.-Optoelectron. 141242–248 (1994). [CrossRef]

14.

F. Ladouceur, J. D. Love, and T. J. Senden, “Measurement of surface roughness in buried channel waveguides,” Electron. Lett. 281321–1322 (1992). [CrossRef]

15.

J. Rodríguez, R. D. Crespo, S. Fernández, J. Pandavenes, J. Olivares, S. Carrasco, I. Ibáñez, and J. Virgós, “Radiation losses on discontinuities in integrated optical waveguides,” Opt. Engr. 381896–1906 (1999). [CrossRef]

16.

P. Lalanne and G. Morris, “Highly improved convergence of the coupled-wave method for TM polarization,” J. Opt. Soc. Am. 13779–784 (1996). [CrossRef]

17.

L. Li, “Use of Fourier series in the analysis of discontinuous periodic structures,” J. Opt. Soc. Am. 131870–1876 (1996). [CrossRef]

18.

A. A. Maradudin, T. Michel, A. R. McGurn, and E. Mendez, “Enhanced backscattering of light from a random grating,” Annals of Physics 203225–307 (1990). [CrossRef]

OCIS Codes
(230.7390) Optical devices : Waveguides, planar
(290.5880) Scattering : Scattering, rough surfaces

ToC Category:
Research Papers

History
Original Manuscript: September 28, 2001
Published: October 22, 2001

Citation
J. Merle Elson, "Propagation in planar waveguides and the effects of wall roughness," Opt. Express 9, 461-475 (2001)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-9-9-461


Sort:  Journal  |  Reset  

References

  1. J. M. Elson and P. Tran, "R-matrix propagator with perfectly matched layers for the study of integrated optical components," J. Opt. Soc. Am. A 16, 2983-2989 (1999). [CrossRef]
  2. J. M. Elson and P. Tran, "Coupled-mode calculation with the R-matrix propagator for the dispersion of surface waves on a truncated photonic crystal," Phys. Rev. B 54, 1711-1715 (1996). [CrossRef]
  3. J. M. Elson and P. Tran, "Band structure and transmission of photonic media: a real-space finite-difference calculation with the R-matrix propagator," NATO ASI Series E: Applied Sciences on Photonic Band Gap Materials Vol. 315 341-354 Crete, Greece June 15-29 (1995).
  4. L. Li, "Multilayer modal method for diffraction gratings of arbitrary profile, depth, and permittivity," J. Opt. Soc. Am. A 11, 2816-2828 (1993). [CrossRef]
  5. L. Li, "Formulation and comparison of two recursive matrix algorithms for modeling layered diffraction gratings," J. Opt. Soc. Am. 13, 1024-1035 (1996). [CrossRef]
  6. J. P. Berenger, "A perfectly matched layer for the absorption of electromagnetic waves," J. Comp. Phys. 114, 185-200 (1995). [CrossRef]
  7. J. P. Berenger, "Three-dimensional perfectly matched layer for the absorption of electromagnetic waves," J. Comp. Phys. 127, 363-379 (1996). [CrossRef]
  8. D. Marcuse, "Mode conversion caused by surface imperfections in a dielectric slab waveguide," Bell Sys. Tech. J. 48, 3187-3215 (1969).
  9. D. Marcuse, "Radiation losses of dielectric waveguides in terms of the power spectrum of the wall distortion function," Bell Sys. Tech. J. 48, 3233-3242 (1969).
  10. F. P. Payne and J. P. R. Lacey, "A theoretical analysis of scattering loss from planar optical waveguides," Opt. and Quantum Elec. 26, 977-986 (1994). [CrossRef]
  11. J. P. R. Lacey and F. P. Payne, "Radiation loss from planar waveguides with random wall imperfections," IEEE Proc. J. 137, 282-288 (1990).
  12. K. K. Lee, D. R. Lim, H. Luan, A. Agarwal, J. Foresi and L. Kimerling, "Effect of size and roughness on light transmission on a Si/SiO2 waveguide: Experiment and model," Appl. Phys. Lett. 77, 1617-1619 (2000). [CrossRef]
  13. F. Ladouceur, J. D. Love and T. J. Senden, "Effect of side wall roughness in buried channel waveguides," IEEE Proc.-Optoelectron. 141, 242-248 (1994). [CrossRef]
  14. F. Ladouceur, J. D. Love and T. J. Senden, "Measurement of surface roughness in buried channel waveguides," Electron. Lett. 28, 1321-1322 (1992). [CrossRef]
  15. J. Rodrguez, R. D. Crespo, S. Fernandez, J. Pandavenes, J. Olivares, S. Carrasco, I. Ibanez, J. Virgos, "Radiation losses on discontinuities in integrated optical waveguides," Opt. Engr. 38, 1896-1906 (1999). [CrossRef]
  16. P. Lalanne and G. Morris, "Highly improved convergence of the coupled-wave method for TM polarization," J. Opt. Soc. Am. 13, 779-784 (1996). [CrossRef]
  17. L. Li, "Use of Fourier series in the analysis of discontinuous periodic structures," J. Opt. Soc. Am. 13, 1870-1876 (1996). [CrossRef]
  18. A. A. Maradudin, T. Michel, A. R. McGurn, and E. Mendez, "Enhanced backscattering of light from a random grating," Annals of Physics 203, 225-307 (1990). [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

OSA is a member of CrossRef.

CrossCheck Deposited