Silicon nitride-on-silicon bi-layer grating couplers were designed for the O-band using an optimization-based procedure that accounted for design rules and fabricated on a 200 mm wafer. The designs were sufficiently robust to fabrication variations to function well across the wafer. A peak fiber-to-chip coupling efficiency to standard single mode fiber of −2.2 dB and a 1-dB bandwidth of 72.9 nm was achieved in the representative device. Over several chips across the wafer, we measured a median peak coupling efficiency of −2.1 dB and median 1-dB bandwidth of 70.8 nm. The measurements had good correspondence with simulation.
© 2018 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
Silicon (Si) photonics can potentially realize low-cost, large-scale photonic integrated circuits (PICs) for applications such as optical switching [1–3], phased arrays [4, 5] and transceivers . However, waveguide loss and fiber-to-chip coupling loss limit the size Si photonic PICs. Multilayer silicon nitride (SiN)-on-Si foundry fabricated integrated photonic platforms address this problem by monolithically integrating lower optical loss, nonlinearity, and thermo-optic coefficient SiN waveguide layers atop the Si layer [7–10]. Such multi-layer platforms also support improved surface fiber-to-chip coupling through bi-layer grating couplers (GCs) with composite SiN and Si teeth. Bi-layer GCs allow for enhanced upwards directivity over single layer gratings by providing sufficient degrees of freedom to control phases and amplitudes of scattered waves, such that the upward emitted waves coupled into the fiber constructively interfere and other scattered waves and modes destructively interfere, as per the analysis in . A coupling efficiency or insertion loss (IL) of −1.3 dB with a 1-dB bandwidth (Δλ1dB) of 80 nm near center wavelength 1550 nm was demonstrated by the SiN-on-Si grating in , without the need for a back-reflector or other post-processing.
In this work, we expand upon  by demonstrating a systematic optimization design method for bi-layer GCs, using a new SiN-on-Si platform at CEA-LETI. The fabricated GCs have high performance consistent with simulation. Over the 9 chips measured from across the wafer, we measured a median IL of −2.1 dB and median Δλ1dB of 70.8 nm. To compare GCs with respect to both coupling efficiency and bandwidth at different center-wavelengths, we define a figure of merit efficiency-bandwidth product (EBWP) to be13–16]. The EBWP is lower than the Si-on-Si bi-layer GC reported in , but that GC requires additional post-processing, and wafer-scale measurements have not been reported.
We proceed in Section 2 by describing the optimization procedure to design the bi-layer GCs. Fabrication and measurements are described in Section 3. We compare our results with prior work in Section 4 and conclude after a discussion on avenues for improvement in Section 5.
2.1. Device parameterization and simulation setup
In our previous work , we heuristically found designs for uniform bi-layer GCs, which were subsequently apodized. Here, we present a systematic procedure for searching for uniform bi-layer GCs with the maximum coupling efficiency using a global optimization algorithm. The GCs were implemented in a new SiN-on-Si platform that is currently being developed at CEA-Leti. The layer thicknesses and material refractive indices were fixed and shown in Fig. 1(a). The GCs used the partially etched Si (150 nm thick) and the SiN (600 nm thick) layers, since GCs using the 300 nm thick Si layer were found to be less efficient.
For the optimization algorithm, the GCs were parameterized according to Fig. 1(b). All the parameters have units of µm, except angle parameters which have units of degrees. The GCs couple light from an angled polished SMF-28 fiber into a SiN waveguide on the chip. The SiN teeth have width, t, and gap, g; and the Si teeth have an offset ox and width tx. We include two transition teeth, which have a width of ta, a gap of ga, and offset of oa, at the start of the GC. These transition teeth have been found to improve coupling efficiency . The minimum feature sizes and clearances of this foundry process are accounted for by setting lower bounds on the parameters. The distance between the SiN waveguide and the center of the fiber is xs, and θ is the emission angle of the grating.
The optimization algorithm uses the value of the merit function computed using 2D Finite Difference Time Domain (FDTD) simulations of the 1D uniform GC in Fig. 1(b). For the FDTD calculations, the fiber input is approximated as a Gaussian source with a mode field diameter (MFD) of 9.2 um at 1310 nm, matching that of SMF-28. The merit function is coupling efficiency, η, at 1310 nm, which is extracted from the optical power coupled into the SiN waveguide mode. After the 1D uniform GC with the highest η is found, the grating teeth are curved to form a focusing GC [Fig. 1(c)] .
2.2. The DIRECT global optimization algorithm
For this work, we chose the DIviding RECTangles (DIRECT) optimization algorithm [18,19] for its reliability and ease of use. DIRECT is derivative-free (i.e., only evaluations of a merit function is needed). This is suitable here because FDTD does not directly provide the gradient of η. DIRECT takes advantage of a regularity condition of the merit-function (Lipschitz continuity) to prioritize search in the most promising regions of the domain, making it much more efficient than brute force parameter sweeping. This is done by recursively taking a rectangular domain, partitioning it, and estimating the optimal values of the sub-domains from the size of the domain and the maximum rate-of-change of the merit-function (Lipschitz constant). Sub-domains that are not promising are eliminated, and promising sub-domains are refined toward the optimum. For electromagnetic design problems, the merit function tends to be continuous functions of the geometric parameters (and hence Lipschitz continuous); thus, DIRECT is often compatible. This is certainly the case of finding the GCs with the maximum η.
A major advantage of DIRECT is that it works reliably without requiring extensive tuning of algorithm parameters to balance between global exploration and local search , unlike particle swarm  or genetic algorithms , which may get trapped in local optima without tuning of the algorithm parameters. However, a limitation of DIRECT is that its rate of convergence tends to slow down as the optimization progresses , making the search for larger domains disproportionately more time consuming than smaller domains. Reducing the size of the search domain greatly speeds up the convergence.
2.3. GC optimization
To speed up the convergence of the optimization for the GC, we first used a coarse sampling step that searched over a large feasible design domain to identify a smaller promising region, on which we applied DIRECT. Using this procedure, we were able to identify a promising design within 4630 simulations.
2.3.1. Coarse sampling
The sampling step was restricted to periodic designs parameterized as p(0) = (t, g, ox, tx, xs, θ), with the addition of 2 transition teeth below the SiN waveguide retaining the same teeth and gap widths as the Si teeth of the main grating to reduce the number of dimensions (i.e. ta = tx, ga = g + t − tx, and oa = ox + tx − t for a the gap from the first Si teeth in the grating and the first transition tooth is g + t − tx). We sampled over a wide yet physically realistic range for each design parameter. Due to factors such as substrate reflections and mode matching, it cannot be known a priori the emission angle, θ, that maximizes η. We let the coupling angle to vary in the range of θ ∈ (0°, 35°) corresponding to positive coupling angles with a safe margin less than the total internal reflection condition for the cladding-air interface (43.7°). The center of the source xs is bounded within (−2, 10) µm, since this is likely where there would be a non-trivial amount of emission from the grating. The g, t, tx, ta, and ga are set between minimum allowable features sizes (0.2 um for g and t, and 0.12 um for tx, ta, and ga) up to coarse features of around 1 um corresponding to an upper limit of roughly wavelength per period, beyond which the resultant grating period would not be efficient for the target wavelength. The offset ox is allowed to vary approximately over an entire period, with ox ∈ (−1, 1) µm. These ranges give rise to lower bounds l(0) = (0.2, 0.2, −1, 0.12, −2, 0) and upper bounds u(0) = (1, 1, 1, 1, 10, 35) for the sampling.
Instead of a grid based sweep, we used a low-discrepancy sequence (LDS) to select sampling points. LDS is a deterministic substitute for uniformly distributed random numbers with more even coverage that is useful for sampling high-dimensional functions and commonly applied to perform Monte-Carlo integration . LDS allows us to select the number of points that we want to sample without being constrained by a grid and increase resolution of the sampling as necessary. We sampled 4000 points from this parameter space using the Halton low discrepancy sequence, which struck a balance between adequate resolution and simulation time. The best 5 results with respect the coupling efficiency η, listed in Table 1, were then used to define a range of parameters to search for in the DIRECT optimization.
2.3.2. DIRECT refinement
After the sampling, the GC designs are refined using DIRECT. In this step, the transition teeth are allowed to vary. We parameterize the grating with p(1) = (t, g, ox, tx, oa, ta, ga, xs, θ), and initialize the bounds as l(1) = (0.5, 0.5, −1, 0.14, −1, 0.14, 0.14, 3, 10), u(1) = (1.2, 1.2, 1, 1, 1, 1, 1.5, 6, 32), using the minimum and maximum of Table 1 with some padding. We use the nlopt  implementation of DIRECT. Including the sampling, the entire process took around 16 hours to complete using a desktop with 8 core Intel i7 processor and 16 GB RAM. The improvement of coupling efficiency over 630 iterations is shown in Fig. 2(a). Although the best result was found after 335 iterations, the algorithm continues to search other partitions. In contrast to local optimization algorithms, where the figure of merit converges with increasing iterations, DIRECT is a global optimization algorithm which continues exploring other regions of the parameter space to test whether better solutions exist. Here, we terminated the algorithm when the figure-of-merit no longer improved.
We select the best design, p(1)∗ = (0.85, 0.617, −0.675, 0.304, −0.74, 0.283, 0.255, 5.83, 28), achieving coupling efficiency of 73%, to be our nominal 1D design. Figure 2(b) plots the transmission of the 1D design, and the corresponding focusing design simulated with 3D-FDTD. An IL of −1.5 dB at 1312 nm with Δλ1dB = 73 nm is expected for the focusing design, representing a EBWP of 3.9×10−2. Figure 2(c) shows that the GC has a high directionality with low back-reflection characteristic of bi-layer gratings, consistent with [11, 12]. Figure 2(d) shows the profile of the emitted power by the GC at the chip surface. The centering of this emission profile is consistent with where the fiber core is expected to be positioned (x = 5.83µm). Figures 3(a)–3(d) show the sensitivity of the nominal 1D design to variations in layer offset, duty cycle, and layer spacing and thickness. The peak coupling efficiency stays within 1 dB to layer offsets of ±60 nm, and layer spacing variations of ±60 nm. The worst case is due to a variation of −40 nm in the width of the Si teeth. The bandwidth is generally negligibly affected, but the center wavelength shifts by up to 50 nm for a layer spacing variation of +60 nm from the nominal.
3. Fabrication and measurement
Test structures consisting of pairs of GCs connected by SiN waveguides as in Fig. 4(a) were fabricated on a 200 mm silicon-on-insulator wafer using 193 nm deep ultraviolet (DUV) photolithography on the Si layer and 248 nm DUV photolithography for the SiN layer, as part of a trial for a SiN-on-Si platform that was under development at CEA-LETI. Fabrication targeted the layer thickness values in Fig. 1(a). While uniformity was good across the wafer for Si and SiN layers, the SiO2 interlayer spacing thickness σ had a large variability. Across the wafer, the layer thicknesses were characterized through ellipsometry: tSi had a range of 8 nm with median of 300 nm and a standard deviation of 2 nm, tetch had a range of 6 nm with median of 153 nm and a standard deviation of 1.4 nm, and tSiN had a range of 14 nm with median of 605 nm and a standard deviation of 4.4 nm. The intended values were 300 nm, 150 nm and 600 nm respectively. From examining Fig. 3(e)–3(f), we do not expect this range of variations to have a major impact on performance. In contrast, the interlayer spacing σ varied significantly. Over the wafer, σ had a range of 140 nm with median 229 nm and standard deviation 37.4 nm. The intended value was 200 nm. Based on the simulation in 3(d), we expect the variation of σ to be a dominant factor in the interdie variations across the measured samples.
We tested 9 samples from across the wafer [Fig. 4(b)], with the corresponding σ labeled. We used a fiber array polished at 29° with index fluid matching fluid (Norland IML150, refractive index n = 1.5) applied at the fiber-chip interface. The input laser polarization was set to transverse electric (TE) with a polarization controller to maximize the coupling efficiency. The transmission spectra are plotted in Fig. 4(c), and Table. 2 summarizes the performance of the measured samples. Over the 9 chips measured from across the wafer, we obtained a worst performance of EBWP = 1.5×10−2 (IL = −4.5 dB with Δλ1dB = 56.6 nm), a best performance of EBWP = 3.4×10−2 (IL = −2.2 dB with Δλ1dB = 72.9 nm), and a median EBWP of 3.3×10−2. The insertion loss and center wavelengths of the measurement data was plotted against interlayer spacing σ, and compared with simulated values in Fig. 4(d). The close correspondence of the trends confirm that the interdie variations are dominantly driven by the interlayer spacing variations in this set of samples.
The measurements showed an additional loss of approximately 0.7 dB compared to the simulation, which agreed well with simulations. The slight discrepancy may be caused by the lack of planarization of the oxide cladding, a slight mismatch in the refractive index between the index matching fluid (n = 1.5) and SiO2, variation or uncertainty in the refractive index, slopes in the sidewalls of the grating features, and the precise position and angle polish of the optical fiber.
For the chip with the grating having the best EBWP, which was also the one with σ = 206 nm closest to the nominal value of 200 nm, we studied the corners of the device for comparison with simulation. Figures 5(a) and 5(b) show the measured spectra for lithographically defined x and y direction offsets of ±30 and ±60 nm between the Si and SiN layer from the nominal design. A ±60 nm x-offset compromised coupling efficiency by at most 0.5 dB. The y-offset negligibly affected the performance. Figures 5(c) and 5(d) show the effect of measured spectra for lithographically defined teeth width (i.e., duty cycle) variation. GCs with −40 nm variation in SiN teeth width decreased the peak coupling efficiency by only 0.5 dB. The coupling efficiency is most sensitive to the width of the Si teeth, which are the smallest features in the device. A −40 nm deviation from the nominal width of 304 nm in the Si teeth led to a peak coupling efficiency of −3.1 dB but maintained a broad bandwidth of Δλ1dB = 69 nm. Overall, the measurements match the trends in the simulations and show reasonable robustness to variations in the interlayer offset and duty cycle.
Table 3 compares the present GC with past GCs in literature. This work achieves a comparable EBWP as the previous 1550 nm bi-layer GC, validating the effectiveness of the design procedure at a different λcenter while using a new platform. The polysilicon (p-Si) on c-Si GC in a CMOS platform has the highest EBWP to date . However, it requires post-processing steps, such as the removal of Si handle substrate. Compared to conventional single layer Si  and SiN GCs , Si GCs with back-reflectors , and subwavelength Si GCs , the EBWP of bi-layer GCs is substantially higher.
In the design procedure, the majority of the simulations were done for the sampling step. By using the grating equation to relate the emission angle to the grating average effective index and the grating period, the number of simulations may be reduced. Here, the sampling step was not time consuming and we let the emission angle to be a free parameter, so we did not use the grating equation as a constraint. Using the grating equation as a constraint may be helpful in designs that have a larger number of variables to optimize. For example, adding more parameters to the optimization, such as layer thicknesses, has the potential to yield even better solutions. On the other hand, it is not straightforward to apply the grating equation to the multilayer gratings. Adapting the grating equation to a bi-layer grating without an interlayer spacing is explored in , but the approach does not apply to our geometry which has a separation between independent waveguide layers. An incident waveguide mode from the SiN waveguide excites the symmetric and anti-symmetric TE modes in the bi-layer, so the effective indices used in the grating equation needs to account for two modes. In contrast, in , the bilayer region supports only one fundamental TE mode.
In this work, it was more effective to use a sampling step rather than the optimization algorithm over the feasible domain. To ascertain whether our design was a global optimum, we attempted to explore the entire design domain using DIRECT without limiting the number of iterations. The initial domain was l(1) = (0.2, 0.2, −1, 0.14, −1, 0.14, 0.14, −2, 0), u(1) = (1, 1, 1, 1, 1, 1, 1.5, 10, 35). DIRECT required around 6000 iterations to reach a similar solution as compared with the sample and refine strategy, and made negligible improvement afterwards up to at least 20,000 iterations. While our design focused on exploring all feasible coupling angles, if instead the design is aimed at a particular emission angle it may be possible to use the grating equation to generate a smaller range of bounds and skip the sampling step altogether. Newer variants of DIRECT  have better convergence behavior, and thus may be able to handle directly the entire feasible domain and bypass the sampling step.
Other methods of designing bi-layer gratings efficiently through custom solvers have also been reported, including a band-structure approach to reduce the simulation domain  and more recently adjoint methods [11, 26] which can converge more quickly by taking advantage of gradient information. Experimental demonstrations of fabricated devices are still pending for these methods. In our work, we opted to use a combination of available implementations of the optimization algorithms and a validated FDTD tool  to minimize debugging. In future work, alternative optimization and simulation methods can be also tried on the SiN-on-Si geometry.
With respect to the fabrication, the large variability in interlayer spacing was a major source of interdie variation. The platform used in this work was under development, and this was the first batch of chips that were made with the SiN and Si layers that allowed the testing of these grating design. The processes used for the planarization had not been fine tuned and resulted in the large variability in interlayer spacing; the uniformity of the interlayer spacing has since been significantly improved. In addition, a planarization step could also be added after the final SiO2 encapsulation, to reduce the losses from surface coupling.
Further apodization of the design to match the input mode is possible by increasing the number of freely varying teeth and successively using local optimization with the DIRECT designed GC as a starting point. A multi-objective design approach, beyond maximizing coupling efficiency, can be implemented by modifying the merit function to incorporate other metrics, such as the bandwidth. By doing so, we can generate designs which optimally trade off bandwidth and coupling efficiency.
In summary, we have demonstrated the effectiveness of the global optimization algorithm DIRECT in designing bi-layer GCs. To the best of our knowledge, this is the first time this algorithm has been applied to the design of photonic devices and bilayer GCs which have shown good performance across a 200 mm wafer. Starting from a large range of physically reasonable parameters, we used a sampling step to identify a smaller promising design space, which was then refined using DIRECT to reach a good design. This methodology was successfully applied to design a SiN-on-Si GC for the O-band. Over dies from across the wafer, the achieved EBWP ranged from 1.5×10−2 to 3.4×10−2, with a median of 3.3 × 10−2. Correspondingly, median insertion loss of −2.1 dB and median 1-dB bandwidth of 70.8 nm were measured. The main source of interdie variation was due to a large variability in the interlayer spacing, which has since been investigated and improved for future runs.
H2020 LEIT Information and Communication Technologies (688516).
References and links
1. E. Bernier, P. Dumais, D. J. Goodwill, H. Mehrvar, D. Celo, J. Jiang, C. Zhang, F. Zhao, X. Tu, C. Zhang, S. Yan, J. He, M. Li, W. Liu, Y. Wei, and D. Geng, “Large-scale silicon photonic switch,” in “Optical Fiber Communication Conference,” (Optical Society of America, 2018), pp. Th1J–1.
2. K. Suzuki, R. Konoike, J. Hasegawa, S. Suda, H. Matsuura, K. Ikeda, S. Namiki, and H. Kawashima, “Low insertion loss and power efficient 32 × 32 silicon photonics switch with extremely-high-Δ PLC connector,” in “Optical Fiber Communication Conference Postdeadline Papers,” (Optical Society of America, 2018), p. Th4B.5.
3. L. Qiao, W. Tang, and T. Chu, “32 × 32 silicon electro-optic switch with built-in monitors and balanced-status units,” Scientific Reports 7, 42306 (2017). [CrossRef]
4. J. Sun, E. Timurdogan, A. Yaacobi, Z. Su, E. S. Hosseini, D. B. Cole, and M. R. Watts, “Large-scale silicon photonic circuits for optical phased arrays,” IEEE Journal of Selected Topics in Quantum Electronics 20, 264–278 (2014). [CrossRef]
5. J. C. Hulme, J. K. Doylend, M. J. R. Heck, J. D. Peters, M. L. Davenport, J. T. Bovington, L. A. Coldren, and J. E. Bowers, “Fully integrated hybrid silicon two dimensional beam scanner,” Optics Express 23, 5861–5874 (2015). [CrossRef] [PubMed]
6. A. Ramaswamy, J. Roth, E. Norberg, R. S. Guzzon, J. Shin, J. Imamura, B. Koch, D. Sparacin, G. Fish, B. G. Lee, R. Rimolo-Donadio, C. Baks, A. Rylyakov, J. Proesel, M. Meghelli, and C. Schow, “A WDM 4×28 Gbps integrated silicon photonic transmitter driven by 32nm CMOS driver ICs,” in “Optical Fiber Communication Conference,” (Optical Society of America, 2015), pp. Th5B–5.
7. W. D. Sacher, Y. Huang, G.-Q. Lo, and J. K. S. Poon, “Multilayer silicon nitride-on-silicon integrated photonic platforms and devices,” Journal of Lightwave Technology 33, 901–910 (2015). [CrossRef]
8. W. D. Sacher, Z. Yong, J. C. Mikkelsen, A. Bois, Y. Yang, J. C. C. Mak, P. Dumais, D. Goodwill, C. Ma, J. Jeong, E. Bernier, and J. K. S. Poon, “Multilayer silicon nitride-on-silicon integrated photonic platform for 3D photonic circuits,” in “Conference on Lasers and Electro-Optics,” (Optical Society of America, 2016), p. JTh4C.3.
9. S. Malhouitre, B. Szelag, S. Brision, Q. Wilmart, D. Fowler, C. Dupré, and C. Kopp, “Heterogeneous and multi-level integration on mature 25Gb/s silicon photonic platform,” in “CPMT Symposium Japan (ICSJ), 2017 IEEE,” (IEEE, 2017), pp. 223–226.
10. C. Baudot, M. Douix, S. Guerber, S. Crémer, N. Vulliet, J. Planchot, R. Blanc, L. Babaud, C. Alonso-Ramos, D. Benedikovich, D. Pérez-Galacho, S. Messaoudène, S. Kerdiles, P. Acosta-Alba, C. Euvrard-Colnat, E. Cassan, D. Marris-Morini, L. Vivien, and F. Boeuf, “Developments in 300mm silicon photonics using traditional CMOS fabrication methods and materials,” in “Electron Devices Meeting (IEDM), 2017 IEEE International,” (IEEE, 2017), pp. 34.
12. W. D. Sacher, Y. Huang, L. Ding, B. J. F. Taylor, H. Jayatilleka, G.-Q. Lo, and J. K. S. Poon, “Wide bandwidth and high coupling efficiency Si3N4-on-SOI dual-level grating coupler,” Optics Express 22, 10938–10947 (2014). [CrossRef] [PubMed]
13. F. Van Laere, T. Claes, J. Schrauwen, S. Scheerlinck, W. Bogaerts, D. Taillaert, L. O’Faolain, D. Van Thourhout, and R. Baets, “Compact focusing grating couplers for silicon-on-insulator integrated circuits,” IEEE Photonics Technology Letters 19, 1919–1921 (2007). [CrossRef]
14. W. S. Zaoui, A. Kunze, W. Vogel, M. Berroth, J. Butschke, F. Letzkus, and J. Burghartz, “Bridging the gap between optical fibers and silicon photonic integrated circuits,” Optics Express 22, 1277–1286 (2014). [CrossRef] [PubMed]
15. C. R. Doerr, L. Chen, Y.-K. Chen, and L. L. Buhl, “Wide bandwidth silicon nitride grating coupler,” IEEE Photonics Technology Letters 22, 1461–1463 (2010). [CrossRef]
16. Q. Zhong, V. Veerasubramanian, Y. Wang, W. Shi, D. Patel, S. Ghosh, A. Samani, L. Chrostowski, R. Bojko, and D. V. Plant, “Focusing-curved subwavelength grating couplers for ultra-broadband silicon photonics optical interfaces,” Optics Express 22, 18224–18231 (2014). [CrossRef] [PubMed]
17. M. T. Wade, F. Pavanello, R. Kumar, C. M. Gentry, A. Atabaki, R. Ram, V. Stojanović, and M. A. Popović, “75% efficient wide bandwidth grating couplers in a 45 nm microelectronics CMOS process,” in “Optical Interconnects Conference (OI), 2015 IEEE,” (IEEE, 2015), pp. 46–47.
18. D. R. Jones, C. D. Perttunen, and B. E. Stuckman, “Lipschitzian optimization without the Lipschitz constant,” Journal of Optimization Theory and Applications 79, 157–181 (1993). [CrossRef]
19. J. M. Gablonsky and C. T. Kelley, “A locally-biased form of the DIRECT algorithm,” Journal of Global Optimization 21, 27–37 (2001). [CrossRef]
20. Q. Liu, G. Yang, Z. Zhang, and J. Zeng, “Improving the convergence rate of the direct global optimization algorithm,” Journal of Global Optimization 67, 851–872 (2017). [CrossRef]
21. Y. Shi and R. C. Eberhart, “Parameter selection in particle swarm optimization,” in “International Conference on Evolutionary Programming,” (Springer, 1998), pp. 591–600.
22. Á. E. Eiben, R. Hinterding, and Z. Michalewicz, “Parameter control in evolutionary algorithms,” IEEE Transactions on Evolutionary Computation 3, 124–141 (1999). [CrossRef]
23. W. J. Morokoff and R. E. Caflisch, “Quasi-Monte Carlo integration,” Journal of Computational Physics 122, 218–230 (1995). [CrossRef]
24. S. G. Johnson, “The NLopt nonlinear-optimization package,” http://ab-initio.mit.edu/nlopt. [Online].
25. J. Notaros and M. Popović, “Band-structure approach to synthesis of grating couplers with ultra-high coupling efficiency and directivity,” in “Optical Fiber Communication Conference,” (Optical Society of America, 2015), pp. Th3F–2.
27. Lumerical Inc., http://www.lumerical.com/tcad-products/fdtd/.