Open Access
Add to BibSonomy
Abstract
An efficient gradient-based topology optimal design approach using beam propagation method (BPM) for a 3-D semi-vectorial design problem of optical waveguide devices is proposed. A semi-vectorial finite-difference BPM (SVFD-BPM) based on an alternating direction implicit method (ADIM) is employed for the wave propagation analysis in our design approach. In comparison to conventional 3-D topology optimal design approaches, it is expected that our approach can reduce its computational cost and design an optical waveguide which is too long to analyze using a finite element method or finite-difference time/frequency-domain method. In this paper, we design several waveguide devices (S-bend, power splitter, and mode converter), and it is shown that our topology optimal design approach surely works.
© 2017 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
1. Introduction
A topology optimal design method has great potentials to achieve further miniaturization of optical waveguide devices and to improve device performance. The topology optimization has much design freedom than a sizing optimization or a shape optimization, and it can optimize the refractive index distribution itself. Thus, it is expected that optical devices with higher performance can be obtained by using the topology optimization. The applications of the topology optimal design to optical devices have been reported and the effectiveness has been demonstrated in previous studies [1–11].
Topology optimal design methods reported so far employ a finite element method (FEM) [1–8], a finite-difference time-domain (FDTD) method [9,10], or a finite-difference frequency-domain (FDFD) method [11] for wave propagation analysis. A beam propagation method (BPM) [12–15] is widely used as an efficient numerical analysis method of optical waveguide devices. An envelope of the unknown electromagnetic field varies slowly for the propagation direction in waveguides without backward reflection. Thus the BPM realizes the cost effective numerical analysis because larger step size is allowed in the propagation direction. It is expected that the devices whose device length is too long to analyze characteristics with FEM, FDTD method, or FDFD method, can be designed by employing the BPM. We have already reported the topology optimal design approach employing the BPM in our previous study [16]. However, this approach is applied only for 2-D design problem. 3-D design is necessary to consider vertical confinement and radiations to a substrate and a cover.
In this paper, our design approach using the BPM is extended for the case of a semi-vectorial 3-D design problem. We employ a semi-vectorial finite-difference BPM (SVFD-BPM) based on an alternating direction implicit method (ADIM) [14, 15] for efficient calculation. Sensitivity analysis is necessary in gradient-based topology optimization, then we employ an adjoint variable method (AVM) as the sensitivity analysis method. Also, we employ a density method which is a way to represent refractive index distribution in the design region. In the density method, the design parameter called a normalized density parameter is allocated to each of discretized domains, and the refractive index in the domain is determined depending on the value of the normalized density parameter. The AVM in the case of using the 2-D scalar or 3-D full-vectorial BPM for a sizing optimization has already been reported in [17,18]. We extend the AVM reported in [17,18] for the case of using the density method and the SVFD-BPM based on the ADIM.
In section 2, we begin with a 3-D scalar Helmholtz equation, then we review the basic equations for SVFD-BPM based on the ADIM. After that, we describe the sensitivity analysis method in the case of using the SVFD-BPM and the density method. In section 3, we design several weakly guiding waveguide devices to verify the validity of the sensitivity analysis and our design approach. Section 4 is the conclusion.
2. Topology optimal design
2.1. Formulation of the SVFD-BPM based on ADIM
We assume core and cladding materials are linear, lossless, isotropic, and frequency-independent. From Maxwell’s equations in frequency domain, the Helmholtz equation of light in the waveguide shown in Fig. 1 can be obtained,
where k0 is a free space wavenumber, n(x, y, z) is a distribution of refractive index, Φ(x, y, z) is a scalar electric or magnetic field. Ex or Hy is used for x-polarization and Ey or Hx is for y-polarization. The definitions of the operators, Dxx and Dyy, are shown in Table 1.
Table 1. The definitions of Dxx and Dyy.
In the BPM, Φ is expressed in the form [12],
where n0 is a reference refractive index. By substituting Eq. (2) into Eq. (1) and employing Fresnel approximation (∂2ϕ/∂z2 = 0), the basic equation of the BPM is obtained, where By discretizing z-direction using Crank-Nicolson scheme and employing ADIM, the update equations of the BPM can be expressed in following form [15]: where m is a step number in z-direction, Δz is the step size, and ϕm denotes ϕ(x, y, mΔz). The matrix-vector form of Eqs. (5) and (6) is obtained by descritizing the transverse direction as follows: where [Γi](i = 1, 2, 3, 4) can be expressed by square matrices. In this paper, we employ Stern’s method [19] to discretize the transverse directions.2.2. Sensitivity analysis
In the density method, the refractive index distribution in the design region is determined by the distribution of a normalized density parameter ρ(x, y, z). We assume uniform core height in the design region, thus in this case, n and ρ do not depend on y in the design region. The refractive index distribution in the design region is represented by using the density method as follows:
where k is a step number in x-direction, n1 and n2 are core and cladding refractive indices respectively, nk,m and ρk,m denote n(kΔx, mΔz) and ρ(kΔx, mΔz) respectively (Δx is a step size in x-direction), and H(ρ) is a modified Heaviside function defined in [16]. The gray region, which has an intermediate refractive index between n1 and n2, can be controlled by a penalty parameter. The refractive index distribution is binarized by making the penalty parameter infinite. In the design examples shown in this paper, the penalty parameter is taken to be a lower value in the initial phase of the optimization, then the value of the parameter is increased as the optimization process is iterated. The penalty parameter is eventually large enough, thus the effect of the binarization is of no significant.The characteristics of optical waveguide devices can be usually expressed by S-parameters. A normalized output power into port-n is represented by following overlap integral in semivectorial problems,
where (ϕ exp(−jk0n0z), ψ exp(−jk0n0z)) = (Ei, Hj), (Hi, Ej) [(i, j) = (x, y), (y, x)], Nz is a step number at an output port in z-direction,* denotes complex conjugate, and ψn is an eigenmode field in port-n.This overlap integral can be expressed in following matrix-vector form by using rectangle rule.
where † denotes Hermitian transpose, Δy is a step size in y-direction. From Eqs. (7) and (8), Eq. (11) can be rewritten as where T denotes transpose, and {λn,m+1}T is defined by By differentiating Eq. (12) with respect to ρk,m, we get the following relation: When m ≠ Nz − 1, the vectors {λn,m} and {λn,m+1/2} can be calculated as follows: If m = Nz − 1, Eq. (15) is replaced with The derivative of the matrices [Γi] (i = 1, 2, 3, 4) can be calculated analytically because the refractive index is represented by Eq. (9).In weakly guiding waveguides, the S-parameters can be expressed approximately by only one component.
where and ϕn is an eigenmode field in port-n. In this case, ψn is replaced with pϕn in Eq. (11).Since the devices treated in this paper are weakly guiding waveguides, we employ Eq. (18) to evaluate the S-parameters.
3. Design examples
We show the validity of our design approach and the sensitivity analysis by designing several weakly guiding waveguide devices.
In the design examples, the density parameters are updated every iteration by a steepest descent method.
3.1. S-bend waveguide
The design model of an S-bend waveguide is shown in Fig. 2. The structural parameters shown in Fig. 2 are the following: n1 = 1.45, n2 = n3 = 1.445, w = 8 μm, h = 4 μm, l = 100 μm, d1 = d2 = 30 μm, s = 15 μm, Lx = 80 μm, and Lz = 900 μm. In order to maximize the normalized output power in port 2 when the fundamental TE-like or TM-like wave with the wavelength of 1.55 μm is launched into port 1, the objective function is set as follows:
The transverse and longitudinal step sizes are taken to be Δx = Δy = 0.2 μm and Δz = 1 μm, respectively. n0 is taken to be the effective index of the fundamental TE-like or TM-like mode in port 1.First, we show that our design approach surely works regardless of which field component is used as ϕ, and the BPM is available in a very low-contrast-index waveguide even though a waveguide structure is complicated to some extent. In this optimal design, a uniform medium (ρ = 0.3 in the design region) is given as an initial structure. In the topology optimization using the density method, fine structures typified by a checker-board pattern often emerge. Such fine structures can be removed by a structural smoothing filter (the detail is described in our previous paper [16]). In this design example, however, the structural smoothing filter is not used in order to show the results without an influence of the smoothing filter.
The normalized output power as a function of the iteration number is shown in Fig. 3. According to Fig. 3, we can see that the normalized output power is improved every iteration, and the power converges in several iterations in the case that ϕ = Ex, Hy, Ey, or Hx. Moreover, an inset in Fig. 3 shows the output power converges quickly. The convergence speed may be strongly affected by the property of an objective function. The composition of the objective function shown in this paper is relatively simple, thus the normalized output power seems to converge relatively quickly. For instance, in the design of an optical triplexer, the convergence speed may be lower because composition of the objective function is more complicated [8]. The objective function is expressed as sum of objectives at three different wavelengths when the triplexer is optimized. Since the structural variation to improve the device performance at a certain wavelength may degrade the performance at another wavelength, the objective function may have complicated dependency on the design variables. The binarized optimized structures and the propagation waves are shown in Fig. 4. The term "binarized" in this paper means that ρ is rounded off to the closest whole number, then the refractive index in the design region is determined to be n1 (core) or n2 (cladding). The normalized output power in the binarized optimized structure in the case that ϕ = Ex, Hy (TE-like wave), Ey, or Hx (TM-like wave) are 0.929, 0.930, 0.924, 0.929, respectively. It is noted that the binarized optimized structures shown in Fig. 4 are almost similar, and the difference of the output properties is not significant even though mosaic-like structures emerge slightly. The normalized power |Sn1|2 is calculated using Eq. (18) in this paper. If the approximation in the formulation of the BPM is not valid for a specified waveguide, the output power is largely different depending on which field component is used as an unknown variable (ϕ) even though the polarization is the same [20]. Therefore, the numerical results indicate that the BPM is available when refractive index difference is very low (about 0.34%) even though waveguide structure is complicated to some extent.
Fig. 3 Normalized output power as a function of the iteration number in the design of an S-bend waveguide in the case that ϕ = Ex, Hy, Ey, or Hx.
Fig. 4 Optimized structure and propagation wave: ϕ is (a),(b) Ex, (c),(d) Hy, (e),(f) Ey, and (g),(h) Hx. The refractive index distribution in the design region is binarized in these structures. Insets (i)–(iii) are magnified views of the optimized structure.
After this design example, the smoothing filter is employed every iteration in order to avoid a very complex profile and obtain a more practical structure. The smoothing filter itself has been widely used in the topology optimization. Topology optimization using the smoothing filter may lead to a local optimum structure that is simple and has good property to some extent. Although the smoothing filter certainly hinders optimization, considering practical manufacturing, it is necessary to employ this filter.
Next, we will show the dependence of optimization results on an initial structure. A uniform medium (ρ = 0.3 in the design region) or a tilted waveguide (ρ = 0.8 in core region, 0.3 in cladding region) is given as an initial structure.
Figure 5 shows that the normalized output power as a function of the iteration number in the case that a uniform medium or a tilted waveguide is given as an initial structure. In Fig. 5, although there are ripples due to the smoothing filter, the output power is improved almost every iteration. The binarized optimized structures are shown in Fig. 6. Although the profiles of the optimized structures are not identical, the operation principles are probably the same. According to Figs. 6(a) and (c), we can see that offsets emerge in the vicinity of the input and output ports in order to reduce mode mismatching. Moreover, in the design region, the waveguide widths are wider than those of the input and output ports to enhance the confinement of light. The normalized output power in the binarized optimized structure shown in Figs. 6(a) and (c) are 0.947 and 0.941, respectively. Certain optimized structures based on the same principle of the operation are automatically obtained whether a special structure is given as an initial structure, or not. In the binarized initial tilted waveguide, the output power is 0.769. Thus the performance of the optimized S-bend waveguides is superior to the binarized initial structures.
Fig. 5 Normalized output power as a function of the iteration number in the design of an S-bend waveguide in the case that a uniform medium or a tilted waveguide is given as an initial structure.
Fig. 6 Optimized structure and propagation wave (|Hy|): The initial structure is (a),(b) a uniform medium and (c),(d) a tilted waveguide. The refractive index distribution in the design region is binarized in these structures. Insets (i)–(iii) are magnified views of the optimized structure.
3.2. Y-branch type power splitter
The design model of a Y-branch waveguide is shown in Fig. 7. The structural parameters shown in Fig. 7 are as follows: n1 = 1.45, n2 = n3 = 1.445, w = 8 μm, h = 4 μm, l = 50 μm, d1 = 40 μm, d2 = 20 μm, s = 15 μm, Lx = 80 μm, and Lz = 900 μm. We assume that the fundamental TE-like wave with the wavelength of 1.55 μm is launched into port 1. The objective function to realize the a: b power splitter is as follows:
In this design, a = b = 0.5 and a symmetric condition is imposed to prevent imbalance of power. The dashed line shown in Fig. 7 represents the symmetric axis.The transverse and longitudinal step sizes are taken to be Δx = Δy = 0.2 μm and Δz = 1 μm, respectively. n0 is taken to be the effective index of the fundamental TE-like mode in port 1. A uniform medium (ρ = 0.3 in the design region) or a straight Y-branch waveguide (ρ = 0.8 in core region, 0.3 in cladding region) is given as an initial structure, then we will show that the dependence of optimization results on an initial structure.
Figure 8 shows that the normalized output power as a function of the iteration number in the case that a uniform medium or a straight Y-branch waveguide is given as an initial structure. Only the normalized output power in port 2 is shown in Fig. 8 because there is no imbalance of the output power. According to Fig. 8, we can see that the normalized power is improved almost every iteration to realize the 1:1 power splitter in both cases. The binarized optimized structures and the propagation wave are shown in Fig. 9. Although the optimized structures shown in Figs. 9(a) and (c) differ, the offset waveguide junctions appear in the vicinity of the output ports and Y-branch structures which can reduce radiation loss into forward direction emerge. A shape of Y-branch like the structures shown in Figs. 9(a) and (c) has been presented, e.g. in [21]. Therefore, the operation principle of both structures seems to be the same. The normalized output power in the binarized structure shown in Figs. 9(a) and (c) are 0.471 and 0.478, respectively. The output power in the binarized initial straight Y-branch waveguide is 0.413, thus our design approach provides the superior structure compared to the initial structures.
Fig. 8 Normalized output power as a function of the iteration number in the design of a 1:1 power splitter in the case that a uniform medium or a straight Y-branch waveguide is given as an initial structure.
Fig. 9 Optimized structure and propagation wave (|Hy|): The initial structure is (a),(b) a uniform medium and (c),(d) a straight Y-branch waveguide. The refractive index distribution in the design region is binarized in these structures.
3.3. Mode converter
Figure 10 shows the design model of a mode converter. This optimal design aims to obtain a converter which converts the 1st higher order TE-like (TE1) wave into the fundamental TE-like (TE0) wave. The structural parameters shown in Fig. 10 are as follows: n1 = 1.45, n2 = n3 = 1.445, w1 = 18 μm, w2 = 6 μm, h = 4 μm, l = 50 μm, d1 = d2 = 40 μm, s = 15 μm, Lx = 80 μm, and Lz = 900 μm. Port 2 is a single mode waveguide at 1.55 μm wavelength. The objective function is set as follows:
where |S21(TE1→TE0)|2 is a normalized output power of the TE0 wave in port 2 when the TE1 wave is launched into port 1.The transverse step sizes are taken to be Δx = Δy = 0.2 μm. The longitudinal step size is taken to be Δz = 2 μm to make calculation efficient. n0 is taken to be the effective index of TE1 mode in port 1. Although the fundamental and the higher order mode are guided in this device, since the relative index difference is very small (Δ≈ 0.3%) and it is expected that the propagation angle may be small enough, the Fresnel approximation is employed also in this design example.A uniform medium (ρ = 0.3 in the design region) or a linear taper (ρ = 0.8 in core region, 0.3 in cladding region) is given as an initial structure.
Figure 11 shows that the output power of the TE0 wave increases almost every iteration. The binarized optimized structure and the propagation wave are shown in Fig. 12. The normalized output power of the TE0 wave is 0.911 in Fig. 12(a), and 0.800 in Fig. 12(b). An optimized structure depends on an initial structure in this design problem. The dependence on an initial structure is mainly determined by the property of the objective function. However, it is difficult to predict the property and a better initial structure. Nevertheless, in weakly guiding waveguides, the dependence is not usually significant if the composition of the objective function is simple, and single-mode and single-wavelength transmission is taken into account in the objective function. Actually, the dependence is not significant in the design of an S-bend and a 1:1 power splitter. Although the better structure is obtained in this design example when a uniform medium is given, it may be better to give the certain initial structure that is based on a known operation principle if the dependence on initial structure is strong.
Fig. 11 Normalized output power of the TE0 wave as a function of the iteration number in the mode converter.
Fig. 12 Optimized structure and propagation wave (|Hy|) in the case that (a),(c),(e) a uniform medium or (b),(d),(f) a linear taper is given as an initial structure. (c),(d) TE1-to-TE0 conversion, and (e),(f) TE0-to-TE1 conversion. In the propagation analysis, Δz is taken to be 1 μm to analyze the device more accurately. The refractive index distribution in the design region is binarized.
This device can be used for a TE0-to-TE1 converter by exchanging input and output ports, however, in this conversion, crosstalk may occur in port 1 because port 1 is not a single mode waveguide. In the superior one, when the TE0 wave is launched into port 2, the normalized output power of the TE1 wave and the crosstalk are 0.912 and −21.32 dB, respectively. Although the insertion loss is −0.4 dB, a low crosstalk TE0-to-TE1 converter is realized.
We investigate the influence of discretization on optimization results because there is possibility that the optimized structure depends on the step size. Although we tried designing with several different step sizes, almost same structure and performance are obtained in this problem, probably because of weakly guiding waveguide.
Finally, we compare the mode converter designed by using our approach with a conventional one. In previous studies on a mode order converter, an asymmetric branching waveguide has been used [22]. Using the sizing optimization based on the steepest descent method, we design a TE1-to-TE0 converter based on an asymmetric branching waveguide. Figure 13 is the design model for a mode order converter. We optimize three parameters: w11, w12, and d12. The other parameters, an objective function, and analysis condition of the BPM are the same as the former design. Figure 14 is the optimized structure and the propagation wave. The optimized parameters are as follows: w11 = 15 μm, w12 = 8.6 μm, d12 = 12 μm. The normalized output power of the TE0 mode in the optimized structure is 0.846. Thus, it is shown that a mode converter that has higher performance can be obtained by using our design approach than that based on an asymmetric branching waveguide.
Fig. 14 Optimized structure of a TE0-to-TE1 converter based on an asymmetric branching waveguide and propagation wave.
4. Conclusion
An efficient topology optimal design approach using the SVFD-BPM and the AVM has been proposed for the first time. The sensitivity analysis method reported in [16] has been extended for the case that the SVFD-BPM based on the ADIM is employed. S-bend waveguides, Y-branch type power splitters, and TE1-to-TE0 mode converters are designed by using our design approach, then it is shown that the properties of the devices are improved every iteration in the optimization process. Therefore, the validity of our design approach is confirmed in the case of weakly guiding waveguides.
To apply our design approach to strongly guiding waveguides, it may be necessary to develop a constraint condition which prohibits appearance of the structure which is difficult to analyze by using the BPM. The applicability of our design approach to strongly guiding waveguides will be discussed in future works.
Funding
Japan Society for the Promotion of Science (JSPS) KAKENHI Grant Number 15K06009.
References and links
1. J. S. Jensen and O. Sigmund, “Systematic design of photonic crystal structures using topology optimization: Low-loss waveguide bends,” Appl. Phys. Lett. 84(12), 2022–2024 (2004). [CrossRef]
2. P. I. Borel, A. Harpøth, L. H. Frandsen, M. Kristensen, P. Shi, J. S. Jensen, and O. Sigmund, “Topology optimization and fabrication of photonic crystal structures,” Opt. Express 12(9), 1996–2001 (2004). [CrossRef] [PubMed]
3. J. S. Jensen, O. Sigmund, L. H. Frandsen, P. I. Borel, A. Harpøth, and M. Kristensen, “Topology design and fabrication of an efficient double 90 photonic crystal waveguide bend,” IEEE Photon. Technol. Lett. 17(6), 1202–1204 (2005). [CrossRef]
4. P.I. Borel, L.H. Frandsen, A. Harpøth, M. Kristensen, J. S. Jensen, and O. Sigmund, “Topology optimised broadband photonic crystal Y-splitter,” Electron. Lett. 41(2), 69–71 (2005). [CrossRef]
5. Y. Tsuji, K. Hirayama, T. Nomura, K. Sato, and S. Nishiwaki, “Design of optical circuit devices based on topology optimization,” IEEE Photon. Technol. Lett. 18(7), 850–852 (2006). [CrossRef]
6. Y. Tsuji and K. Hirayama, “Design of optical circuit devices using topology optimization method with function-expansion-based refractive index distribution,” IEEE Photon. Technol. Lett. 20(12), 982–984 (2008). [CrossRef]
7. T. Yasui, Y. Tsuji, J. Sugisaka, and K. Hirayama, “Design of three-dimensional optical circuit devices by using topology optimization method with function-expansion-based refractive index distribution,” J. Lightw. Technol. 31(23), 3765–3770 (2013). [CrossRef]
8. Z. Zhang, Y. Tsuji, T. Yasui, and K. Hirayama, “Design of ultra-compact triplexer with function-expansion based topology optimization,” Opt. Express 23(4), 3937–3950 (2015). [CrossRef] [PubMed]
9. R. Matzen, J. S. Jensen, and O. Sigmund, “Systematic design of slow-light photonic waveguides,” J. Opt. Soc. Amer. B 28(10), 2374–2382 (2011). [CrossRef]
10. L. H. Frandsen, Y. Elesin, L. F. Frellsen, M. Mitrovic, Y. Ding, O. Sigmund, and K. Yvind, “Topology optimized mode conversion in a photonic crystal waveguide fabricated in silicon-on-insulator material,” Opt. Express 22(7), 8525–8532 (2014). [CrossRef] [PubMed]
11. J. Lu and J. Vuckkovic, “Nanophotonic computational design,” Opt. Express 21(11), 13351–13367 (2013). [CrossRef] [PubMed]
12. M. D. Feit and J. A. Fleck, “Light propagation in graded-index optical fibers,” Appl. Opt. 17(24), 3990–3998 (1978). [CrossRef] [PubMed]
13. Y. Chung and N. Dagli, “An assessment of finite difference beam propagation method,” IEEE J. Quantum. Electron. 26(8), 1335–1339 (1990). [CrossRef]
14. J. Yamauchi, T. Ando, and H. Nakano, “Beam-propagation analysis of optical fibers by alternating direction implicit method,” Electron. Lett. 27(18), 1663–1665 (1991). [CrossRef]
15. P. L. Liu and B. J. Li, “Semivectorial beam-propagation method for analyzing polarized modes of rib waveguides,” IEEE J. Quantum. Electron. 28(4), 778–782 (1992). [CrossRef]
16. A. Iguchi, Y. Tsuji, T. Yasui, and K. Hirayama, “Topology optimization of optical waveguide devices based on beam propagation method with sensitivity analysis,” J. Ligthw. Technol. 34(18), 4214–4220 (2016). [CrossRef]
17. M. A. Swillam, M. H. Bakr, and X. Li, “Efficient adjoint sensitivity analysis exploiting the FD-BPM,” J. Ligthw. Technol. 25(7), 1861–1869 (2007). [CrossRef]
18. M. A. Swillam, M. H. Bakr, and X. Li, “Full vectorial 3-D sensitivity analysis and design optimization using BPM,” J. Lightw. Technol. 26(5), 528–536 (2008). [CrossRef]
19. M. S. Stern, “Semivectorial polarised H field solutions for dielectric waveguides with arbitrary index profiles,” IEE Proc. J 135(5), 333–338 (1988).
20. J. Yamauchi, Y. Nito, and H. Nakano, “A modified semivectorial BPM retaining the effects of the longitudinal field component and its application to the design of a spot-size converter,” J. Lightw. Technol. 27(13), 2470–2476 (2009). [CrossRef]
21. Z. Weissman, E. Marom, and A. Hardy, “Very low-loss Y-junction power divider,” Opt. Lett. 14(5), 293–295 (1989). [CrossRef] [PubMed]
22. W. Chen, P. Wang, and J. Yang, “Optical mode interleave based on the asymmetric multimode Y junction,” Photon. Technol. Lett. 26(20), 2043–2046 (2014). [CrossRef]
