OSA's Digital Library

Optics Express

Optics Express

  • Editor: C. Martijn de Sterke
  • Vol. 18, Iss. 5 — Mar. 1, 2010
  • pp: 4088–4102
« Show journal navigation

Arbitrary-order interface conditions for slab structures and their applications in waveguide analysis

Yih-Peng Chiou and Cheng-Han Du  »View Author Affiliations


Optics Express, Vol. 18, Issue 5, pp. 4088-4102 (2010)
http://dx.doi.org/10.1364/OE.18.004088


View Full Text Article

Acrobat PDF (222 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

We extend the continuity relations of field derivatives across an abrupt interface to arbitrary orders for transverse electric and magnetic waves in slab structures. Higher-order finite-difference formulation is then obtained by combining the systematically-obtained interface conditions with Taylor series expansion. Generalized Douglas scheme is also adopted to further enhance the convergence of truncation errors by two orders. We apply the derived finite-difference formulation, up to nine-points in this paper, to solve the guided modes in simple a slab waveguide and multiple quantum well waveguides. The results shows the truncation error is much higher, up to tenth order, as expected. Using those higher-order schemes, accurate results are obtained with much fewer sampled points, and hence with tremendously less computation time and memory.

© 2010 Optical Society of America

1. Introduction

Semivectorial [1

1. M. S. Stern, “Semivectorial polarised finite difference method for optical waveguides with arbitrary index profiles,” Inst. Elect. Eng. Proc. J. 135, 56–63 (1988).

] and full-vectorial [2

2. C. L. Xu, W. P. Huang, M. S. Stern, and S. K. Chaudhuri, “Full-vectorial mode calculations by finite difference method,” Proc. Inst. Elect. Eng. J. 141, 281–286 (1994). [CrossRef]

] formulations for graded-index structure were proposed and applied in commercial tools. Such techniques were further extended to beam propagation analysis [3–5

3. M. S. Stern, “Rayleigh quotient solution of semivectorial field problems for optical waveguides with arbitrary index profiles,” Proc. Inst. Elect. Eng. J. 138, 185–190 (1991).

]. Although graded-index approximation is applicable to structures with minimal index difference, it has difficulty in accurately modeling field behaviors near abrupt interfaces. Moreover, recent advancement of nanoprocess makes step-index structures practical solutions for next-generation optical design, indicating that a modified modeling formulation dedicated to step-index problem is required. Vassallo [6

6. C. Vassallo, “Improvement of finite difference methods for step-index optical waveguides,” Proc. Inst. Elect. Eng. J. 139, 137–142 (1992).

] proposed formulation for one-dimensional (1D) and two-dimensional (2D) guiding structure in scalar and semivectorial cases, which takes boundary continuity relations of field and derivatives into consideration. Chiou et al. have extended this technique to full-vectorial waveguide mode solver and photonic crystal analysis [7–9

7. Y.-C. Chiang, Y.-P. Chiou, and H.-C. Chang, “Improved full-vectorial finite-difference mode solver for optical waveguides with step-index profiles,” J. Lightwave Technol. 20, 1609–1618 (2002). [CrossRef]

]. Other modeling algorithm considering abrupt interface conditions are developed as well.

One critical property of finite-difference method is its convergence of truncation error. The convergence of commonly used three-point formulation is normally no more than second order (O(h 2)), since only two adjacent sampled points are referenced to approximate second-order derivative in Helmholtz equation. Higher-order convergence is desired in certain situations to achieve higher accuracy and less computation. One effective method to increase the order is using generalized Douglas (GD) scheme. Hadley [10

10. G. R. Hadley, “Low-truncation-error finite difference equations for photonics simulation I: beam propagation,” J. Lightwave Technol. 16, 134–141 (1998). [CrossRef]

] uses generalized Douglas scheme to reduce error of beam propagation method on step-index problems, while modification for graded-index problem has also been reported [11

11. J. Yamauchi, J. Shibayama, and H. Nakano, “Application of the generalized Douglas scheme to optical waveguides analysis,” Optic. Quantum Eletron. 31, 675–687 (1999). [CrossRef]

]. Chiou et al. [12

12. Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, “Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices,” J. Lightwave Technol. 18, 243–251 (2000). [CrossRef]

] has developed three-point formulation considering interface conditions with generalized Douglas scheme and successfully achieved fourth-order convergence (O(h 4)) for mode solving of step-index structures. Such formulation is adopted in waveguide discontinuities [13

13. H. Zhang, Q. Guo, and W. Huang, “Analysis of waveguide discontinuities by a fourth-order finite-difference reflective scheme,” J. Lightwave Technol. 25, 556–561 (2007). [CrossRef]

] and extended three-dimensional semivectorial formulation, by including two-transverse directions, to solve wide-angle beam-propagation method [14

14. J. Shibayama, T. Takahashi, J. Yamauchi, and H. Nakano, “A three-dimensional multistep horizontally wide-angle beam-propagation method based on the generalized Douglas scheme,” Photon. Tech. Lett. 18, 2535-2537 (2006). [CrossRef]

] and surface plasmon metallic waveguides [15

15. B. Hu, P. Sewell, J. G. Wykes, A. Vukovic, and T. M. Benson, “Fourth-order accurate subsampling for finite-difference analysis of surface plasmon metallic waveguides,” Microwave Opt. Technol. Lett. 50, 995–1000 (2008). [CrossRef]

].

Convergence more than O(h 4) may be desired when we want to achieve even higher accuracy or even less computation. Hadley proposed a finite-difference formulation magically free of truncation error [16

16. G. R. Hadley, “Low-truncation-error finite difference equations for photonics simulation II: vertical-cavity surface-emitting lasers,” J. Lightwave Technol. 16, 142–151 (1998). [CrossRef]

]. It is an exact solution for eigenmode in one-dimensional problem by placing the sampled points at interfaces. It appears about 40% faster than matrix method. However, it can not be directly applied to beam propagation analysis, since the eigen value n̄2 from ∂2/∂z 2 in Helmholtz equation is inside the trigonometric functions. The derivation of finite-difference formulation is back to the commonly used Taylor series expansion. To achieve convergence of higher orders, we may include more sampled points to evaluate higher-order derivatives. However, the formulation is usually derived in homogeneous region, and it does not consider the discontinuity across abrupt interfaces. Therefore, when there are abrupt interfaces, the convergence is still of low order. It might not be problematic when the number of interfaces is small, since it is still of high order in the rest homogeneous regions. When the number of interfaces are large, like multiple quantum well structures, the convergence will be lowered, especially in structures of higher index difference. Therefore, a generalized higher-order finite-difference formulation considering the interface conditions is required in such cases to improve the numerical performance.

2. Derivation of Higher-Order Interface Conditions

For monochromatic electromagnetic wave propagating of angular frequency ω in source-free, isotropic, and piecewisely homogeneous, Maxwell equations in derivatives form are expressed as

×E=B
(1)
×H=D
(2)
·D=0
(3)
·B=0.
(4)

The Helmholtz equations for electrical and magnetic fields are satisfied in a piecewisely homogeneous region as

2E+k02n2E=0
(5)

and

2H+k02n2H=0,
(6)

E=E¯exp(jβz)
(7)

and

H=H¯exp(jβz)
(8)

respectively. Helmholtz equations then become

2E¯x2=ω2μεE¯=β2E¯
(9)

and

2H¯x2=ω2μεH¯=β2H¯.
(10)

The interface conditions are the continuity of tangential components of electrical and magnetic field intensities, Ey, Ez, Hy, and Hz, and the continuity of normal components of electrical displacement Dx = εEx and magnetic flux density Bx = μHx as in Fig. 1. From the continuity of the field components, we will derive the continuity of their derivatives.

Fig. 1. Continuity of fields on interface and sample points.

2.1. TE case

In transverse electric (TE) modes, Ey, Hx, and Hz are not zero. From the continuity of the tangential electrical and magnetic fields, we have

EyR=EyL
(11)

and

HzR=HzL
(12)

where subscript L and R represent infinitesimally close to the interface on the left and right, respectively, as shown in Fig. 1. From (1), we have

Eyx=jωμHz.
(13)

Since Hz is continuous across the interface, we have continuity relation for the first-order derivative

EyRx=μRμLEyLx.
(14)

The continuity of field and its first derivative in (11) and (14) are commonly used. To achieve higher-order convergence in the finite-difference method, the continuity conditions of higher-order derivatives should be derived. From one-dimensional Helmholtz equation

2Eyx2+ω2μεEy=β2Ey
(15)

we have

Ey=1β2[2x2+ω2με]Ey.
(16)

From (11), we have the continuity relation for the second-order derivative

[2x2+ω2μLεL]EyR=[2x2+ω2μRεR]EyL
(17)

or

2EyRx2=[2x2+ω2μLεL]EyLω2μRεREyR=[2x2+(ω2μLεLω2μRεR)]EyL.
(18)

Substituting (18) and (16) into (11) and (14) recursively, we have generalized continuity relations for even- and odd-order derivatives expressed as

2nx2nEyR=[2x2+(ω2μLεLω2μRεR)]nEyL.
(19)

and

2n+1x2n+1EyR=μRμLx[2x2+(ω2μLεLω2μRεR)]nEyL.
(20)

respectively. We prove (19), and (20) can be proved similarly. First, the statement is true n = 0 from (11). Second, we assume the statement is true for n = k, that is,

2kx2kEyR=[2x2+(ω2μLεLω2μRεR)]kEyL.
(21)

Substituting (16) into both sides of (21) we have

2kx2k{1β2[2x2+ωR2μRεR]EyR}=1β2{2k+2EyRx2k+2+ω2μRεR2kEyRx2k}
=[2x2+(ω2μLεLω2μRεR)]k{1β2[2x2+ω2μLεL]EyL}.
(22)

Substituting (21) into (22), we have

2k+2EyRx2k+2=[2x2+(ω2μLεLω2μRεR)]k{[2x2+ω2μLεL]EyL}
ω2μRεR{[2x2+(ω2μLεLω2μRεR)]kEyL}
=[2x2+(ω2μLεLω2μRεR)]k+1EyL.
(23)

The statement is also true for n = k + 1. Therefore, by the principle of mathematical induction, the statement holds for all non-negative integers.

Equations (19) and (20) up to 2N-th order can be written in matrix form as

[ψRψRψRψR(3)ψR(4)ψR(2N1)ψR(2N)]=[1000000a0000b010000ab0a00b202b0100C0N1abN10C1N1abN200C0NbN0C1NbN10C2NbN2CNN][ψLψLψLψL(3)ψL(4)ψR(2N1)ψL(2N)]
(24)

or

ΨR=MR:LΨL,
(25)

where a = μRL and bk = ω2LεL − μRεR). We will use this equation to drive higher-order finite-difference formulation.

2.2. TM case

Derivation in transverse magnetic modes is quite similar to that in TE case. Hy, Ex and Ez are nonzero in TM mode. Continuity relations of the tangential fields

HyR=HyL
(26)

and

EzR=EzL,
(27)

while (2) leads to

Hyx=jωεEz.
(28)

From (27) and (28), we have continuity relation of the first-order derivative

HyRx=εRεLHyLx
(29)

Similar to (19) and (20), generalized continuity relations are expressed as

2nx2nHyR=[2x2+(ω2μLεLω2μRεR)]nHyL.
(30)

and

2n+1x2n+1HyR=εRεLx[2x2+(ω2μLεLω2μRεR)]nHyL.
(31)

Similarly, a matrix form expression for (30) and (31) is obtained as that for (24) and (25), except that a = εRL in TM case.

3. Implementation in Mode Solvers

3.1. Normal Finite-Difference Scheme

Considering a single abrupt interface and equispaced sampled points shown in Fig. 1, we first relate ψ i+1 and its derivatives with ψi and its derivatives up to 2N-th order. Sampled points in homogeneous medium can be described by Taylor series expansion up to 2N-th order, and we have the relation between ψL and ψi as

[ψLψLψLψL(2N)]=[1pp22!p2N(2N)!01pp2N1(2N1)!001p2N2(2N2)!0001][ψiψiψiψi(2N)]+HOT
(32)

or

ΨL=ML:iΨi+HOT,
(33)

where HOT stands for truncated higher-order terms. Similarly, the relation between ψ i+1 and ψR as

or

[ψi+1ψi+1ψi+1ψi+1(2N)]=[1qq22!q2N(2N)!01qq2N1(2N1)!001q2N2(2N2)!0001][ψRψRψRψR(2N)]+HOT
(34)

or

Ψi+1=Mi+1:RΨR+HOT.
(35)

with (33), (33), and interface continuity conditions (25), the relation between Ψi+1 and Ψi is thus established as

Ψi+1=Mi+1:RΨR+HOT=Mi+1:RMR:LΨL+HOT=Mi+1:RMR:LML:iΨi+HOT
=Mi+1:iΨi+HOT.
(36)

Note that, there can be multiple interfaces between the sampled points ψ i+1 and ψi. For example, if there are T interfaces with left and right fields represented as ψLt and ψRt at tth interface, respectively, we can have the relation between Ψi+1 and Ψi as

Ψi+1=Mi+1:RTMRT:LTMLT:R(T1)ML(t+1):RtMRt:LtMR1:L1ML1:iΨi+HOT
=Mi+1:iΨi+HOT
(37)

where M Rt∙Lt’s are from matching the interface conditions, and M i+1∙RT , M L1:i, and M R(t+1)∙Lt’s are from Taylor series expansion. Applying the same procedure recursively, we can obtain the relation between Ψi+j and Ψi,

Ψi+j=Mi+j:i+j1Ψi+j1+HOT==Mi+j:i+j1Mi+1:iΨi+HOT
=Mi+j:iΨi+HOT,
(38)

where j = 1, 2, ,⋯ , N. Similarly, we can have relation for negative js.

We can express ψ i+j in terms of ψi and its derivatives from the first row in Ψi+j = M i+j:iΨi

ψi+j=[uj,0uj,1uj,2uj,2N]Ψi+O(h2N+1).
(39)

In the derivation, there is no approximation in matching the boundary conditions. All the errors are from the truncation of Taylor series expansion. Therefore, the truncation error is O(h 2N+1). After applying similar procedure on all the (2N + 1) sampled points, we have the relation between the field at sampled points and ψi and its derivatives up to 2Nth order as

[ψiNψiψi+N]=[uN,0uN,1uN,2NuN,0uN,1uN,2N][ψiψi(j)ψi(2N)]+O(h2N+1)
(40)

or

Fi=UiΨi+O(h2N+1).
(41)

If the higher-order terms are neglected, and the finite-difference approximation for the derivatives is expressed as Dxjψi, we have

Fi=UiDi,
(42)

where D i = [ψi, Dxψi, D 2 xψi, ⋯ , Dxjψi, ⋯ , D 2N x ψi]T. Finally, we obtain the finite-difference approximation of derivatives up 2N-th order considering the interface conditions at of ψi with (2N + 1) sampled points

Di=BiFi,
(43)

where B i = U -1 i. The third row in (43) is the approximation of second derivative expressed in terms of the sampled points, which is O(h 2N+1)/h 2 = O(h 2N-1) accurate. We use centered-difference scheme here, but it is not compulsory. We can also use forward-difference or backward-difference scheme. The advantage of using centered-difference scheme is that the truncation error can be enhanced by one order, from O(h 2N-1) to O(h 2N), when the medium is homogeneous and the grids are uniformly positioned [19

19. R. L. Burden and J. D. Faires, Numerical Analysis, 8th ed. (Brooks Cole, 2004).

].

To solve the eigen problem in mode solvers, we replace the second derivative ∂2 ψ/∂x 2 in (9) and (10) with D 2 x ψ. Then we have the discretized Helmholtz equation

Dx2ψ+ω2μεψ=β2ψ.
(44)

As mentioned, the truncation error in D 2 x ψ is O(h 2N) if the medium is homogeneous between the (2N + 1) points and the grids are uniformly positioned [19

19. R. L. Burden and J. D. Faires, Numerical Analysis, 8th ed. (Brooks Cole, 2004).

]. If the grids are nonuniformly positioned, the truncation error is usually reduced to O(h 2N-1). If there is an abrupt interface or interfaces, the truncation error is usually O(h 2N-1). In waveguide problems, most grids are usually in homogeneous regions. Therefore, the eigen problem in (44) is O(h 2N) accurate in global sense.

3.2. Generalized Douglas Finite-Difference Scheme

Generalized Douglas scheme can also be implemented to magically enhance the accuracy by two orders. With two additional derivatives included in (41), we have

Fi=UiΨi+G1ψi(2N+1)+G2ψi(2N+2)+O(h2N+3),
(45)

where G 1 and G 2 are (2N + 1) × 1 vectors. Multiplying the inverse matrix B i = U -1 i as in (43), we have improved formulation

Di=BiFi=Ψi+Bi(G1ψi(2N+1))+G2ψi(2N+2)+O(h2N+1),
(46)

where the third row is expressed as

Dx2ψi=j=Nj=Nb3,jψi+jψi+v1ψi(2N+1)+v2ψi(2N+2)
=(1+v12N1x2N1+v22Nx2N)ψi,
(47)

where vj is the third element of the vector V = B i G j in (46). We further approximate 2N1x2N1 and 2Nx2N with D 2N-1 x and D 2N x in (43), and we have the generalized Douglas scheme to improve the accuracy of second-order derivative

Dx2ψi(1+v1Dx2N1+v2Dx2N)ψi,
(48)

which is sometimes written as

ψi=Dx21+v1Dx2N1+v2Dx2Nψi.
(49)

The accuracy is increased by two orders by including two additional terms.

Multiply the Helmholtz equation with (1 + v 1 D 2N-1 + v 2 D 2N x), we finally have improved discretized Helmholtz equation

Dx2ψ+ωμε(1+v1Dx2N1+v2Dx2N)ψ=β2(1+v1Dx2N1+v2Dx2N)ψ
(50)

based on generalized Douglas scheme, which considers the interface conditions. Note that matrix is still banded and the bandwidth is the same, (2N+1). The computation time and memory in the eigenvalue problem is the same, while the accuracy is enhanced by two orders. That is, the truncation error in is O(h 2N+2) if the medium is homogeneous between the (2N + 1) points and the grids are uniformly positioned. If the grids are nonuniformly positioned, the truncation error is usually reduced to O(h 2N+1). If there is an abrupt interface or interfaces, the truncation error is O(h 2N+1). In waveguide problems, most grids are usually in homogeneous region. Therefore, the eigen problem in (50) is usually O(h 2N+2) accurate in global sense.

4. Applications in Slab Waveguides

The convergence of propagation constants in slab waveguides is investigated in this section. First, a simple waveguide with a homogeneous core is investigated. Then multiple quantum waveguides (MQW) with more than than one-hundred interfaces in the guiding layer are investigated. The eigenvalues and vectors of the guided modes are solved. By slight modification, it can be directly applied to beam propagation analysis.

4.1. Slab Waveguide

Error analysis of slab waveguide simulation is demonstrated in this section. Refractive indices of core and cladding are ncore = √11.088 and nclad = √11.044, respectively. Waveguide width is 2 μm and wavelength is 1.55 μm. The computation domain is 40 μm in width and transparent boundary condition (TBC) is adopted on the outer boundary. The relative error of propagation constant β is defined as

Relativeerrorεr=βcalculatedβexactk0(ncorenclad),
(51)

where βcalculated is calculated from finite-difference method, and βexact is from the analytical solution of dispersion relation.

Figures 2 and 3 illustrate the relative errors in TE and TM modes, respectively, for different grid sizes ∆x. The thinnest black and blue lines are from three-point formulation based on graded-index (GI) formulation as in [1

1. M. S. Stern, “Semivectorial polarised finite difference method for optical waveguides with arbitrary index profiles,” Inst. Elect. Eng. Proc. J. 135, 56–63 (1988).

] with or without index averaging. Intermediate thick and thickest lines are from calculation considering the interface conditions without and with generalized Douglas (GD) scheme as in Section 3, respectively. The approximation of a step index interface with graded index is only O(h 0). The higher index difference the larger deviation of the approximation is. In this case, however, only two interfaces with small index difference exist, and most sampled points are in homogeneous region. Therefore, using conventional graded-index (GI) approach can still yield O(h 2) accuracy. The red lines are calculated from three-point formulation without and with GD scheme as proposed in [12

12. Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, “Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices,” J. Lightwave Technol. 18, 243–251 (2000). [CrossRef]

]. Similarly, the three-point formulation in (44) and (50) are usually O(h 1) and O(h 3) accurate, respectively. However, most sampled pointed are in homogeneous region, the accuracy are still O(h 2) and O(h 4) accurate, respectively. Furthermore, if the interfaces are at the center between two sampled points or at one sampled point, the three-point formulation in (44) and (50) can be O(h 2) and O(h 4) accurate, respectively.

The green, pink, and cyan lines are calculated from five-, seven-, and nine-point formulation, respectively. For (2N+1)-point calculation with boundary condition matching, the convergence is approximately O(h 2N) as the intermediate thick lines, and they can be further increased up to O(h 2N+2) by applying GD scheme as in the thickest lines. Overall, error analysis quite matches theoretical expectation. Including more points in evaluating the second derivative in Helmholtz equation, we can reach desired accuracy with much fewer grid points. For example, if we want to reach relative error εr = 10-8, the grid number in nine-point formulation with GD scheme is about 1/1000 of that in conventional three-point formulation. The required memory in the former is only about 3/1000 of that in the latter. The computation time in finding the eigenvalue in the former is almost negligible as compared with that in the latter.

Furthermore, the limitation of accuracy can also be put into lower value. Total numerical error is from truncation error and round-off error. Truncation error is due to the dropping of higher order terms (HOT) in the derivation. As we can see, the P-th order truncation error, O(hP), decreases polynomially with grid size h, Round-off error is due to finite digits of the floating point numbers in computers. It accumulates with the number of operations and increases with grid number polynomially. Therefore, the lower bound of of the accuracy is limited by the round-off error. As we can see in the figure, lower-order formulations require more operations and hence result in larger round-off error. For example, the lower bound of the relative error in three-point formulation is about 10-9 and that in nine-point formulation with GD scheme is about 10-14.

4.2. Multiple-Quantum-Well Waveguides

We extend to multiple quantum well (MQW) waveguides the finite-difference analysis based on higher-order approximation and GD scheme in this subsection. The guiding region consists of 56 barriers and 55 wells. The parameters are the same as those in [12

12. Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, “Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices,” J. Lightwave Technol. 18, 243–251 (2000). [CrossRef]

]. Widths of each barrier and well are 0.012 μm and 0.007 μm, respectively. Refractive indices of upper cover, wells, barriers, and lower cladding are n covr = 3.2224, n well = 3.3704 n well = 3.3704, n barr = 3.2874, and n clad = 3.2224, respectively. The relative error of propagation constant is defined as

Fig. 2. Relative error of TE mode propagation constant with respect to grid size in a simple slab waveguide.
Fig. 3. Relative error of TM mode propagation constant with respect to grid size in a simple slab waveguide.
Relativeerror=βcalculatedβexactk0(nwellnclad).
(52)

There are 111 thin layers in the guiding region. Ohke et al. proved that the MQW region could be approximately replaced with a homogeneous layer with an averaged or effective index [17

17. S. Ohke, T. Umeda, and Y. Cho, “Optical waveguides using GaAs-AlxGa1-xAs multiple quantum well,” Opt. Commun. 56, 235–239 (1985). [CrossRef]

, 18

18. S. Ohke, T. Umeda, and Y. Cho, “TM-mode propagation and form birefringence in a GaAs-AlGaAs multiple quantum well optical waveguide,” Opt. Commun. 70, 92–96 (1989). [CrossRef]

]. The averaged index is

navg2=iwini2iwi
(53)

for TE polarization and

1navg2=iwini−2iwi
(54)

for TM polarization, where wi is the thickness of i-th layer with refractive index ni.

Figures 4 and 5 show the relative errors of propagation constants of TE and TM fundamental modes for different grid sizes, respectively. The blue thinnest lines are from conventional method [1

1. M. S. Stern, “Semivectorial polarised finite difference method for optical waveguides with arbitrary index profiles,” Inst. Elect. Eng. Proc. J. 135, 56–63 (1988).

], and no index averaging is used in the guiding region. The refractive index at a sampled point in the layered region is either n well or n barr depending on the point position. The convergence is relatively slow and roughly O(∆x) with relatively large fluctuations. As mentioned, the graded index approximation is only O(h 0) accuracy. The results do not converge well since there are more than one-hundred layers.

The black thinnest lines are from conventional formulation [1

1. M. S. Stern, “Semivectorial polarised finite difference method for optical waveguides with arbitrary index profiles,” Inst. Elect. Eng. Proc. J. 135, 56–63 (1988).

], and index averaging in [17

17. S. Ohke, T. Umeda, and Y. Cho, “Optical waveguides using GaAs-AlxGa1-xAs multiple quantum well,” Opt. Commun. 56, 235–239 (1985). [CrossRef]

, 18

18. S. Ohke, T. Umeda, and Y. Cho, “TM-mode propagation and form birefringence in a GaAs-AlGaAs multiple quantum well optical waveguide,” Opt. Commun. 70, 92–96 (1989). [CrossRef]

] is used in guiding region. The convergence is about O(∆x 2) in TE polarization. It is worse in TM polarization and close to O(∆x). The electric field and its derivative in TE polariation are continuous across the interface, therefore, the index averaging gives a better approximation. However, the derivative of magnetic field in TM polarization is discontinuous across the interface, the index averaging does not working so well as that in TE polarization.

Intermediate thick and thickest lines are from calculation considering the interface conditions without and with generalized Douglas (GD) scheme as in Section 3, respectively. When the grid size is larger than the layer thickness, there are more than one interface between two sampled points. In such situation, the technique in (37) is adopted.

The red lines are calculated from three-point formulation without and with GD scheme as proposed in [12

12. Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, “Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices,” J. Lightwave Technol. 18, 243–251 (2000). [CrossRef]

]. Note that, the well width is 0.007 μm and the barrier width is 0.012 μm. Similarly, the three-point formulation in (44) and (50) are usually O(h 1) and O(h 3) accurate, respectively. As mentioned in [12

12. Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, “Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices,” J. Lightwave Technol. 18, 243–251 (2000). [CrossRef]

], if the interfaces are at the center between two sampled points or at one sampled point, the three-point formulation in (44) and (50) can be O(h 2) and O(h 4) accurate, respectively. Besides, the truncation error may be positive or negative depending on the interface position. Since there are many layers, the truncation error from different layers may be compensate each other to reach higher accuracy. Therefore, in large grid size, the accuracy fluctuates between O(h) and O(h 2) for the formulation without GD scheme, and it fluctuates between O(h 3) and O(h 4) for that with GD scheme. When the grid is smaller, more grids are in homogeneous region. The fluctuation becomes smaller with finer grid size and the accuracy is approximately O(h 2) and O(h 4) for the formulation without and with GD scheme, respectively.

Again, the green, pink, and cyan lines are calculated from five-, seven-, and nine-point formulation, respectively. Similar to three-point formulation, when the grid size is large, the truncation error for the formulation without GD scheme is between O(h 2N-1) and O(h 2N) as the intermediate thick lines, and that with GD scheme is between O(h 2N+1) and O(h 2N+2) as in the thickest lines. In seven-point formulation with GD scheme and nine-point formulation, when the round-off error dominates, the grid size is still larger than layer thickness of wells or barriers. The lines still fluctuate and no straight section exists.

Fig. 4. Relative error of TE mode propagation constant with respect to grid size in a MQW waveguide.
Fig. 5. Relative error of TM mode propagation constant with respect to grid size in a MQW waveguide.

We furthermore investigate a MQW waveguide with large index difference in guiding layers. Refractive indices are n well = 3.4, n barr = 1.6, and n cover = n clad = 1.5. Layer number and thickness of wells and barriers are the same. Relative error for TE and TM modes are illustrated in Figs. 6 and 7, respectively. The truncation error with (2N + 1)-point formulation is still between O(h 2N-1) and O(h 2N), and it is between O(h 2N+1) and O(h 2N+2) with generalized Douglas scheme. When the grid size is smaller, most points are in homogeneous region, the truncation error is almost O(h 2N) and O(h 2N+2) for formulation without and with GD scheme. The transverse variation is much larger due to larger index difference. The propagation constant diverges for grid size larger than 0.1 μm.

Fig. 6. Relative error of TE mode propagation constant with respect to grid size in a MQW waveguide with larger index difference.

4.3. Some Remarks

Although the derivation in this paper is for mode solving, the formulation can be directly applied to beam propagation method as commonly adopted three-point formulation. No extra limitation is imposed. By replacing β 2 with ∂2/∂z 2 and using ψ(x,z) = Ψexp(-jβz), we have

Dx2Ψ+ω2μεΨβ¯2Ψ+2Ψz2=j2β¯Ψz.
(55)

In paraxial approximation, ∂2Ψ/∂z 2 is dropped and Crank-Nicolson scheme is adopted. The resultant matrices are (2N + 1) banded instead of tridiagonal in three-point formulation. Note that wide angle techniques based on Padé approximants and multistep method [20

20. G. R. Hadley, “Multistep method for wide-angle beam propagation,” Opt. Lett. 17, 1743–1745 (1992). [CrossRef] [PubMed]

] can be implemented. In addition, absorption boundary conditions, reflection and transmission in waveguide discontinuities, and bidirectional beam propagation can all be implemented as usual three-point finite-difference beam propagation method. (Ref. [21

21. J. Yamauchi, Propagating Beam Analysis of Optical Waveguides, (Exeter, UK: Research Studies Press, 2003).

] gives details of propagation beam analysis.) Furthermore, it can be applied to semivectorial analysis of three-dimensional guiding structures as in [14

14. J. Shibayama, T. Takahashi, J. Yamauchi, and H. Nakano, “A three-dimensional multistep horizontally wide-angle beam-propagation method based on the generalized Douglas scheme,” Photon. Tech. Lett. 18, 2535-2537 (2006). [CrossRef]

] and [15

15. B. Hu, P. Sewell, J. G. Wykes, A. Vukovic, and T. M. Benson, “Fourth-order accurate subsampling for finite-difference analysis of surface plasmon metallic waveguides,” Microwave Opt. Technol. Lett. 50, 995–1000 (2008). [CrossRef]

], for both mode solver and beam propagation. With the help of higher-order formulation, the grid number is much less and computation time and memory are greatly reduced again in those numerical algorithms.

Fig. 7. Relative error of TM mode propagation constant with respect to grid size in a MQW waveguide with larger index difference.

5. Conclusion

We have proposed a novel technique to systematically derive the continuity of field derivatives at a slab interface up to arbitrary orders. The derived interface conditions can combine directly with Taylor series expansion to build finite-difference formulation for Helmholtz equation. We apply the proposed finite-difference formulation to the analysis of simple slab waveguides and multiple quantum waveguides. In simple slab waveguides, the truncation error with (2N + 1)-point formulation are O(h 2N), and it is enhanced to O(h 2N+2) with generalized Douglas scheme. In multiple quantum well waveguides, the truncation error with (2N + 1)-point formulation without and with generalized Douglas scheme are, at worst, O(h 2N-1) and O(h 2N+1), respectively. Adopting higher-order formulation, we can reach desired accuracy with much fewer grid points. Therefore, it can significantly reduce the computation time and memory. Furthermore, the ultimate accuracy is also lower due to smaller round-off error. Such higher-order formulation may find its applications in various optical and electromagnetic slab structures, especially when high accuracy is essential in their analysis. The formulation can be directly applied to beam propagation analysis. Three-dimensional semivectorial analysis is also possible by extending one-dimensional transverse direction to two-dimensional transverse directions.

Acknowledgement

This work was supported in part by the National Science Council of Taiwan under grant NSC98-2221-E-002-170-MY3 and in part by the Excellent Research Projects of National Taiwan University under grant 98R0062-07.

References and links

1.

M. S. Stern, “Semivectorial polarised finite difference method for optical waveguides with arbitrary index profiles,” Inst. Elect. Eng. Proc. J. 135, 56–63 (1988).

2.

C. L. Xu, W. P. Huang, M. S. Stern, and S. K. Chaudhuri, “Full-vectorial mode calculations by finite difference method,” Proc. Inst. Elect. Eng. J. 141, 281–286 (1994). [CrossRef]

3.

M. S. Stern, “Rayleigh quotient solution of semivectorial field problems for optical waveguides with arbitrary index profiles,” Proc. Inst. Elect. Eng. J. 138, 185–190 (1991).

4.

P.-L. Liu and B.-J. Li, “Semivectorial beam-propagation method for analyzing polarized modes of rib waveguides,” IEEE J. Quantum Electron. 28, 778–782 (1992). [CrossRef]

5.

W. Huang, C. Xu, S.-T. Chu, and S. K. Chaudhuri, “The finite-difference vector beam propagation method: analysis and assessment,” J. Lightwave Technol. 10, 295–305 (1992). [CrossRef]

6.

C. Vassallo, “Improvement of finite difference methods for step-index optical waveguides,” Proc. Inst. Elect. Eng. J. 139, 137–142 (1992).

7.

Y.-C. Chiang, Y.-P. Chiou, and H.-C. Chang, “Improved full-vectorial finite-difference mode solver for optical waveguides with step-index profiles,” J. Lightwave Technol. 20, 1609–1618 (2002). [CrossRef]

8.

Y.-P. Chiou, Y.-C. Chiang, C.-H. Lai, C.-H. Du, and H.-C. Chang, “Finite-difference modeling of dielectric waveguides with corners and slanted facets,” J. Lightwave Technol. 27, 2077–2086 (2009). [CrossRef]

9.

Y.-C. Chiang, Y.-P. Chiou, and H.-C. Chang, “Finite-difference frequency domain analysis of 2-D photonic cry-tals with curved dielectric interfaces,” J. Lightwave Technol. 26, 971–976 (2008). [CrossRef]

10.

G. R. Hadley, “Low-truncation-error finite difference equations for photonics simulation I: beam propagation,” J. Lightwave Technol. 16, 134–141 (1998). [CrossRef]

11.

J. Yamauchi, J. Shibayama, and H. Nakano, “Application of the generalized Douglas scheme to optical waveguides analysis,” Optic. Quantum Eletron. 31, 675–687 (1999). [CrossRef]

12.

Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, “Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices,” J. Lightwave Technol. 18, 243–251 (2000). [CrossRef]

13.

H. Zhang, Q. Guo, and W. Huang, “Analysis of waveguide discontinuities by a fourth-order finite-difference reflective scheme,” J. Lightwave Technol. 25, 556–561 (2007). [CrossRef]

14.

J. Shibayama, T. Takahashi, J. Yamauchi, and H. Nakano, “A three-dimensional multistep horizontally wide-angle beam-propagation method based on the generalized Douglas scheme,” Photon. Tech. Lett. 18, 2535-2537 (2006). [CrossRef]

15.

B. Hu, P. Sewell, J. G. Wykes, A. Vukovic, and T. M. Benson, “Fourth-order accurate subsampling for finite-difference analysis of surface plasmon metallic waveguides,” Microwave Opt. Technol. Lett. 50, 995–1000 (2008). [CrossRef]

16.

G. R. Hadley, “Low-truncation-error finite difference equations for photonics simulation II: vertical-cavity surface-emitting lasers,” J. Lightwave Technol. 16, 142–151 (1998). [CrossRef]

17.

S. Ohke, T. Umeda, and Y. Cho, “Optical waveguides using GaAs-AlxGa1-xAs multiple quantum well,” Opt. Commun. 56, 235–239 (1985). [CrossRef]

18.

S. Ohke, T. Umeda, and Y. Cho, “TM-mode propagation and form birefringence in a GaAs-AlGaAs multiple quantum well optical waveguide,” Opt. Commun. 70, 92–96 (1989). [CrossRef]

19.

R. L. Burden and J. D. Faires, Numerical Analysis, 8th ed. (Brooks Cole, 2004).

20.

G. R. Hadley, “Multistep method for wide-angle beam propagation,” Opt. Lett. 17, 1743–1745 (1992). [CrossRef] [PubMed]

21.

J. Yamauchi, Propagating Beam Analysis of Optical Waveguides, (Exeter, UK: Research Studies Press, 2003).

OCIS Codes
(000.4430) General : Numerical approximation and analysis
(040.4200) Detectors : Multiple quantum well
(130.2790) Integrated optics : Guided waves
(230.7400) Optical devices : Waveguides, slab

ToC Category:
Integrated Optics

History
Original Manuscript: January 11, 2010
Revised Manuscript: February 11, 2010
Manuscript Accepted: February 11, 2010
Published: February 16, 2010

Citation
Yih-Peng Chiou and Cheng-Han Du, "Arbitrary-order interface conditions for slab structures and their applications in waveguide analysis," Opt. Express 18, 4088-4102 (2010)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-18-5-4088


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. M. S. Stern, "Semivectorial polarised finite difference method for optical waveguides with arbitrary index profiles," Inst. Elect. Eng. Proc. J. 135, 56-63 (1988).
  2. C. L. Xu, W. P. Huang, M. S. Stern, and S. K. Chaudhuri, "Full-vectorial mode calculations by finite difference method," Proc. Inst. Elect. Eng. J. 141, 281-286 (1994). [CrossRef]
  3. M. S. Stern, "Rayleigh quotient solution of semivectorial field problems for optical waveguides with arbitrary index profiles," Proc. Inst. Elect. Eng. J. 138, 185-190 (1991).
  4. P.-L. Liu and B.-J. Li, "Semivectorial beam-propagation method for analyzing polarized modes of rib waveguides," IEEE J. Quantum Electron. 28, 778-782 (1992). [CrossRef]
  5. W. Huang, C. Xu, S.-T. Chu, and S. K. Chaudhuri, "The finite-difference vector beam propagation method: analysis and assessment," J. Lightwave Technol. 10, 295-305 (1992). [CrossRef]
  6. C. Vassallo, "Improvement of finite difference methods for step-index optical waveguides," Proc. Inst. Elect. Eng. J. 139, 137-142 (1992).
  7. Y.-C. Chiang, Y.-P. Chiou, and H.-C. Chang, "Improved full-vectorial finite-difference mode solver for optical waveguides with step-index profiles," J. Lightwave Technol. 20, 1609-1618 (2002). [CrossRef]
  8. Y.-P. Chiou, Y.-C. Chiang, C.-H. Lai, C.-H. Du, and H.-C. Chang, "Finite-difference modeling of dielectric waveguides with corners and slanted facets," J. Lightwave Technol. 27, 2077-2086 (2009). [CrossRef]
  9. Y.-C. Chiang, Y.-P. Chiou, and H.-C. Chang, "Finite-difference frequency domain analysis of 2-D photonic crytals with curved dielectric interfaces," J. Lightwave Technol. 26, 971-976 (2008). [CrossRef]
  10. G. R. Hadley, "Low-truncation-error finite difference equations for photonics simulation I: beam propagation," J. Lightwave Technol. 16, 134-141 (1998). [CrossRef]
  11. J. Yamauchi, J. Shibayama, and H. Nakano, "Application of the generalized Douglas scheme to optical waveguides analysis," Opt. Quantum Eletron. 31, 675-687 (1999). [CrossRef]
  12. Y.-P. Chiou, Y.-C. Chiang, and H.-C. Chang, "Improved three-point formulas considering the interface conditions in the finite-difference analysis of step-index optical devices," J. Lightwave Technol. 18, 243-251 (2000). [CrossRef]
  13. H. Zhang, Q. Guo, and W. Huang, "Analysis of waveguide discontinuities by a fourth-order finite-difference reflective scheme," J. Lightwave Technol. 25, 556-561 (2007). [CrossRef]
  14. J. Shibayama, T. Takahashi, J. Yamauchi, and H. Nakano, "A three-dimensional multistep horizontally wideangle beam-propagation method based on the generalized Douglas scheme," Photon. Tech. Lett. 18, 2535-2537 (2006). [CrossRef]
  15. B. Hu, P. Sewell, J. G. Wykes, A. Vukovic, and T. M. Benson, "Fourth-order accurate subsampling for finitedifference analysis of surface plasmon metallic waveguides," Microwave Opt. Technol. Lett. 50, 995-1000 (2008). [CrossRef]
  16. G. R. Hadley, "Low-truncation-error finite difference equations for photonics simulation II: vertical-cavity surface-emitting lasers," J. Lightwave Technol. 16, 142-151 (1998). [CrossRef]
  17. S. Ohke, T. Umeda, and Y. Cho, "Optical waveguides using GaAs-AlxGa1−xAs multiple quantum well," Opt. Commun. 56, 235-239 (1985). [CrossRef]
  18. S. Ohke, T. Umeda, and Y. Cho, "TM-mode propagation and form birefringence in a GaAs-AlGaAs multiple quantum well optical waveguide," Opt. Commun. 70, 92-96 (1989). [CrossRef]
  19. R. L. Burden and J. D. Faires, Numerical Analysis, 8th ed. (Brooks Cole, 2004).
  20. G. R. Hadley, "Multistep method for wide-angle beam propagation," Opt. Lett. 17, 1743-1745 (1992). [CrossRef] [PubMed]
  21. J. Yamauchi, Propagating Beam Analysis of Optical Waveguides, (Exeter, UK: Research Studies Press, 2003).

Cited By

Alert me when this paper is cited

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


« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited