Open Access
Add to BibSonomy
Abstract
Adjoint optimization is an effective method in the inverse design of nanophotonic devices. In order to ensure the manufacturability, one would like to have control over the minimal feature sizes. Here we propose utilizing a level-set method based on b-spline surfaces in order to control the feature sizes. This approach is first used to design a wavelength demultiplexer. It is also used to implement a nanophotonic structure for artificial neural computing. In both cases, we show that the minimal feature sizes can be easily parameterized and controlled.
© 2020 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
1. Introduction
Adjoint method has become a very useful tool in the inverse design of photonics [1–16] and among other optimization methods such as genetic algorithms [17–19], swarm optimization [20], or nonlinear optimization [21], adjoint method is particularly efficient as it computes the gradient of a cost function w.r.t. a large number of design parameters.
The inverse design process mostly involves the optimization of the dielectric constant distribution $\varepsilon _{(r)}$. In a gradient-based approach, $\varepsilon$ is treated as a continuous variable. However, in practice, the permittivity can only take discrete values determined by the types of materials used. To bridge this gap, one can apply the method of level-set [3] or density-based [11] topology optimization. However, even with level-set, or density-based topology optimization to address the issue of discrete material parameters, we still end up with structural features that might be too small for today’s nanofabrication. Additional constraints need to be applied to control the size of the smallest features, such as adding new terms to the cost function, using filters to remove small features throughout optimization, or fixing the size of the pixels that make up the device, above critical values [5,11,15,21–24]. In this work, we use b-splines to achieve effective control of feature sizes. This approach has been a highly successful and proven approach for size control in structural shape optimization of mechanical structures [24,25]. And splines in general have previously been used to define curved shapes in photonic structures [26]; however, to the best of our knowledge, our work is the first to adapt b-spline surfaces in the context of inverse designing photonic devices.
There are many advantages to using b-splines for topology optimization as we present it here. B-splines can make anisotropic designs in the $x$ and $y$ axes, and varying spatial resolution (different feature sizes) throughout the structure possible, by merely varying some simple parameters that makeup the splines across the simulation domain. Another approach that b-splines can be used to design devices while taking fabrication constraints into account are digitized devices. This design approach makes up the device using elements of fixed size that satisfy minimal feature-size fabrication constraints. For example [27] proposes a pseudo adjoint optimization method to design such devices, or [21] uses fixed square pixel size elements along with a type of genetic optimization algorithm to optimize a nanophotonic device. Similar structures can be easily achieved using b-splines while still enjoying all the benefits of a straightforward adjoint method implementation. Moreover, b-splines are easily implementable and memory efficient since, as is explained in the following sections, a $D$ dimensional b-spline is implemented as the tensor-product of $D$ univariate(1D) b-splines. Furthermore, defining the level-set function that is being optimized using a b-spline, provides the surface with control parameters that directly affect the feature sizes of the final device. This allows the optimization to be done in a single stage. Of course, we can also adopt a multi-stage optimization scheme and fine tune our design by adding extra considerations to our cost function and using the conventional implementation of level-set method in the second stage, similar to [15]. Finally, we should point out that a similar control over the level-set function can be achieved by using a different set of basis functions such as radial basis functions [5]. However, b-spline functions are numerically more stable and easier to implement.
2. B-splines
Splines are functions generally used for interpolating data and designing curves and surfaces in computer aided designs and computer graphics. These functions are generated by a linear combination of piece-wise polynomial functions. Their strength is in the fact that they are capable of accurate interpolation of data with low degree polynomial components, that would otherwise require a high degree non-local polynomial. The degree of the spline is the degree of its highest order polynomial.
One way of representing a spline is by defining a set of basis functions, called basis-spline functions or b-spline functions. In order to define any one of these basis functions, the domain over which the spline is going to be defined, has to be segmented with a set of knots $U = [u_0,\ldots, u_m]$. In this work we only work with knots that appear only once in this set, and are equidistant from one another. Another property that is required to define a b-spline function is its degree (k). With both knots and degree defined, basis-spline functions can be expressed using the Cox-de Boor recursion formula [28,29].
Fig. 1. A domain of 1 to 31 has been segmented by 4 knots in the set $U=[1,11,21,31]$ (the knots have been have shown with red circles on the x-axis). (a)-(c) show the b-splines of degree 0 each of which spans 1 interval. (d) and (e) show the b-splines of degree 1 that span 2 intervals, and (f) shows the b-spline of degree 2 that spans 3 intervals.
With the basis spline functions, a curve can be generated using a linear combination of these basis functions. In Eq. (2), the $P(i)$’s are the set of coefficients that determine the contribution of each b-spline component to the curve.
We now discuss the control of feature sizes. The degree of the basis function(k) plays a critical role in controlling the minimal feature size. This can be clearly seen in Fig. 1: as the degree of the basis function is increased, the minimum area that it affects also increases and so the minimum size of the features increases. Another way to control the minimal feature size is to tune the spatial density of the knots. A similar effect happens when the knots are distributed more sparsely over the domain. This effect can be seen in Fig. 2: the denser the knots are packed together, the smaller the area each of these functions can affect becomes. Consequently, more nuances can be generated in a curve. This is the approach we have taken in this work to control the critical dimensions of our designs.Fig. 2. The domain from 1 to 31 has been segmented by (a)8, (b)7, and (c)6 knots to create cubic b-splines(b-splines of degree $k=3$). The location of each of the knots has been shown with a red circle on the x-axis.
3. Formulation
B-splines can be used for a density based approach or a level-set one. Next we use 2 dimensional (2D) photonic design to illustrate this method in the context of the level-set method. In 2D, we use basis spline functions to construct a surface $\phi _{(x,y)}$ using the following equation:
In order to use the level-set method for topology optimization, we first need to define a surface $\phi _{(x,y)}$ over the area that the optimized device will occupy (The value of each point on this surface can vary across a limited range around zero). Then the medium can be shaped according to Eq. (1) based on this surface.
Fig. 3. (a) The structures resulting from the knots of a cubic b-spline being distributed at every $2^{nd}$ and $8^{th}$ (b) $2^{nd}$ and $2^{nd}$ (c) $8^{th}$ and $2^{nd}$ grid points on the vertical and horizontal axes respectively.
We also mentioned that the spatial resolution of the device can be tuned using the b-spline. In the previous section it was explained that the knot distribution and the degree of the splines can be used for this purpose. Since in the next section of the paper we will focus on the implementation of the former, here we merely show an example of the latter. Figure 4 demonstrates how using splines of different degrees changes the feature sizes of the nanophotonic device.
Fig. 4. nanophotonic structures resulting from b-splines of degree (a) 0, (b) 2, (c) 4. The degree of splines on both axes are set equal to each other in all these cases.
The next case that we will briefly discuss here is the method that can provide us with digitized structures. In order to design a device made up of square pixels, we can simply set the degree of the splines to zero. This allows each b-spline coefficient to control one of the basis functions and consequently one of the pixels. The size of the pixels can then be tuned by changing the knot distribution on the simulation grid. Figure 5 depicts examples of such structures made up of square pixels with three different sizes.
Fig. 5. nanophotonic structures resulting from b-splines of degree (a) 0, (b) 2, (c) 4. The degree of splines in $x$ and $y$ axes are set equal to each other in all these cases.
Now, optimizing the medium with a level-set method consists of evolving the boundaries between the two materials (all the points where $\phi (r)=0$) in a manner that improves the performance of the device. With a b-spline surface for this task, we have direct control on this surface, and consequently on the boundary between the constituent materials. Having defined the level-set surface as a b-spline surface, we can now directly optimize a cost function J with respect to the control parameters of the b-spline surface.
Once the gradient of the cost function is calculated, gradient descent can be used to update the b-spline coefficients. However, in order to ensure the stability of the process, in each update iteration the coefficients are normalized by the highest value among them [31]. The mentioned steps can be summarized as shown in Eq. (6), where $P^t$ is the matrix of all the b-spline coefficients at iteration $t$, $||P^t||_\infty$ is the maximum absolute value in this matrix, and $\gamma$ is a constant (called the learning rate) with which gradient descent updates the coefficients.
4. Application
In order to demonstrate the application of this approach, we utilized it to design two nanophotonic devices: a wavelength demultiplexer [4] and a nanophotonic structure for artificial neural computing [32]. Of course, such binary devices have been designed and fabricated in many other works before; here we only wish to use these designs to demonstrate the utility of our approach.
The wavelength demultiplexer is a three port device that guides the light coming from the input port to one of the two output ports based on its wavelength. In our design, this device is of size $5\,\mu m\times 2.5\,\mu m$ designed to demultiplex $1500\,nm$ and $1550\,nm$ wavelengths. We optimize three versions of this device with the intention of showcasing how the critical dimensions of the device can be tuned by the basis spline functions.
The b-spline surface is defined on the same grid as the one used for the electromagnetic simulation and over the area to be optimized. We can use either the distribution of the knots on the grid, or the degree of the b-splines to control the critical dimension. In this work we opt for the former. We define our basis functions as cubic b-splines, and use three different knot distributions to realize three version of the devices with different feature sizes. The first design has the knots placed at every second point on the simulation grid, the second design at every fifth point, and the third design at every ninth point. This effectively provides us with 4559, 629, and 152 b-spline coefficients respectively. Once the b-spline surface is set up, we use Eq. (4) to define the medium, where $\varepsilon _1$ is set to the permittivity of $SiO_2$ and $\varepsilon _2$ is set to that of $Si$, and the structures are simulated using Finite-Difference Frequency-Domain(FDFD), with a spatial resolution of $25\,nm$s as we show in Code 1 (Ref. [33]). Finally, a cost function is necessary to optimize the structure. The cost function we define here tries to maximize the the transmission of the time-averaged energy flux to the output port corresponding to the input wavelength.
The devices were optimized with respect to the cost function mentioned above using Eq. (4) and Eq. (5). The results for the optimization of these three instances have been depicted in Fig. 6. As it was mentioned earlier, since placing the knots more sparsely results in an increase in the area that each basis-spline function affects, the critical dimension of the device increases. This can clearly be seen in Fig. 6. As we move from Fig. 6(a) to Fig. 6(c) the feature sizes distinctly expand. However, even for the case of the sparsely distributed knots, the occasional nuances in the b-spline surface can cause the appearance of small features in the device. This problem is mitigated by applying erosion and dilation processes (morphological operations generally used in image processing) on the device in each update iteration to eliminate these occasional smaller features. Ref [3] uses these operations to eliminate the narrow bridges that appear between small features in their devices throughout the design process. We should point out that this control over the feature sizes of the device comes at a price: there is a trade off between the performance of the device and the critical dimensions of the structure. As we use larger critical dimensions, the degree of freedom, i.e. the number of b-Spline functions used, decreases, and the optimization is further constrained, consequently the performance generally degrades. This is also evident in Fig. 6, where the transmission drops moving from the first device to the last one.
Fig. 6. The wavelength demultiplexers acting on $1500\,nm$ and $1550\,nm$ wavelengths. These devices are designed with P matrices (matrices of the b-spline coefficients) with sizes (a) $97\times 47$, (b) $37\times 17$, and (c) $19\times 8$. As shown here, the higher the number of coefficients(corresponding to packing the b-spline knots closer together), the smaller the details of the designed device, and the higher the confinement of the field to the device. (d) shows the transmission for each of these three instances. Transmission for $1500\,nm$ is shown with continuous lines, while the dashed lines represent the transmission for $1550\,nm$. The degrading effect of reducing the number of b-spline coefficients on transmission can be clearly seen here.
Nanophotonic artificial neural computing This device can accomplish the same task of neural inference as in [32] . The objective is to recognize the value of a handwritten digit 0-9. The input light from the handwritten digit is focused by a nanophotonic structure to different locations according to the value of the digit. The medium then has to map all the different writing styles of the same digit to the specific location corresponding to that digit; furthermore, the medium must perform well on handwritten digits that it has not yet seen.
The setup for this task is as follows: the image of the handwritten digit is first vectorized and mapped to an array of light sources placed on the left side of the nanophotonic structure. Then a set of ten receivers corresponding to each of the ten classes 0-9 are placed on the right side of the medium. Inside the device, the input light scatters and reflects off of the numerous material interfaces, and is then focused on a spatial location on the output side. On the output side, the receivers measure the light intensity at their respective locations and the class with the maximum measured value is selected as the correct class for the input image.
The device is made possible by a much more extensive optimization process than that for the wavelength demultiplexer. This optimization process is equivalent to the training process in digital neural networks. We define a cost function and use adjoint method to calculate the gradient. At the time of training, we measure the intensity of light at the location of all the receivers. Then we use a cross-entropy cost function between the normalized output array and a $1\times 10$ one-hot vector(a vector that is all zeros except at the index of the correct class), to calculate the cost value for that specific instance. Details of the working principle can be found in Ref [32]. In this work we implement the level-set function based on a b-spline surface to control the minimal critical dimension. In our previous work [32], the level-set method was implemented based on a surface set as a signed distance function. Here, the initial medium is created by generating a random set of values for the b-spline coefficients, and then the medium is created based on the generated surface with the surface evolving afterwards by updating the b-spline coefficients. On the other hand, in Ref [32] the medium was generated by randomly distributing inclusions made up of a second material inside the host medium, and then a surface $\phi _{(r)}$ is generated based on this initial medium (similar to the previous section, the surface relates to the medium with Eq. (4)). The structure then evolves by updating this signed-distance surface at each training iteration using a modified version of the following equation.
In this equation $v(x, y)$ is the velocity with which each point on the zero crossing curves of the level-set function moves normal to each of those curves. This velocity is set equal to the gradient acquired with the adjoint method. This implementation of the level-set method for the image classification task turned out to be quite dependant on the initial distribution of the inclusions, and therefore we had to initialize the medium with many small inclusions to get a good performance from the device. Doing so resulted in small feature sizes in the device which make the fabrication of the device rather difficult. However, with the b-spline approach to level-set function this problem does not happen as the optimization process is less dependant on the initial shape of the medium.We optimize two versions of this device, one with the signed distance function approach, and the other with the b-spline approach. Both devices were optimized for the wavelength $1\,\mu m$ with the size $42\,\mu m\times 30\,\mu m$ and made up of $SiO_2$ as the host material and $air$ as the inclusions. For the signed-distance implementation 2000 inclusions of size $5\,\mu m\times 5\,\mu m$ are spread throughout the device for initialization, and for the b-spline implementation the knots are distributed at every tenth point on the simulation grid in both axial dimensions. It should be pointed out that, similar to the previous section, both devices are simulated using FDFD, with an spatial resolution of $100\,nm$s. The accuracy on the test set for the two devices turned out to be $76.8\%$ and $77.2\%$. The performance of the b-spline approach is slightly better in this case which is not of great consequence as the gap between the two values can be reduced with better initialization of the signed-distance approach. However, the important point that is apparent in Fig. 7 is the difference between the shape of the two structures. As it can be seen, the first device shown in Fig. 7(a) has much smaller features whereas the second device depicted in Fig. 7(b) has a much larger critical dimension. This makes the second device much easier to fabricate, and thus the b-spline approach much more desirable. Finally, Fig. 7(c) and Fig. 7(d) show the device in action. As it can be seen, light from different class of images is focused at distinct locations corresponding to that class. Moreover, different image samples of the same digit produce different field distributions; however, these different field distributions result in the same output class.
Fig. 7. Nanophotonic structures can perform artificial computing. The structures are optimized by adjoint method. The two structures can realize similar performance levels but the b-spine method (b) produces structures that are much easier to fabricate. (a) The structure acquired by following the evolution of a signed distance level-set surface. (b) A medium achieved by following the b-spline surface scheme with larger minimal feature sizes than (a). (c) and (d) show the operation of the device shown in (b) for different input samples. (c) Two image samples of the same digit 2 generate different field distributions but they are both recognized as the same digit. (d) The same phenomenon can be seen for two image samples of handwritten digit 8, where different field distributions result in recognition of the same digit.
5. Conclusion
In this work, we adopt a proven technique to size control, that has been a useful tool in the field of mechanics, to the inverse design of photonic devices. It uses b-spline basis to describe the distribution of the dielectric constant. The minimal feature size is controlled by the knots’ distribution and degre. This method is very easy to implement and interpret. We demonstrated the proposed method for a wavelength demultiplexer, and a nanophotonic structure for artificial neural computing. The first class of nanophotonic devices we implemented was used to demonstrate how the trade-off between their complexity of and their performance is handled. For the second class, two variations of the level-set method were used to achieve the same goal. The comparison between the two cases shows how having a well-defined set of control parameters can help us in designing a nanophotonic device in one optimization stage, without the need to have a good initial guess of the variables we are trying to optimize.
Funding
National Science Foundation (1561917); Defense Advanced Research Projects Agency (HR00111820046).
Acknowledgments
The authors thank W. Shin for his help on improving the computational speed of the implementation of this method. The views, opinions and/or findings expressed are those of the author and should not be interpreted as representing the official views or policies of the Department of Defense or the U.S. Government.
Disclosures
The authors declare no conflicts of interest.
References
1. C. M. Lalau-Keraly, S. Bhargava, O. D. Miller, and E. Yablonovitch, “Adjoint shape optimization applied to electromagnetic design,” Opt. Express 21(18), 21693–21701 (2013). [CrossRef]
2. F. Callewaert, V. Velev, P. Kumar, A. Sahakian, and K. Aydin, “Inverse-designed broadband all-dielectric electromagnetic metadevices,” Sci. Rep. 8(1), 1358 (2018). [CrossRef]
3. A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7(1), 1786 (2017). [CrossRef]
4. A. Y. Piggott, J. Lu, K. G. Lagoudakis, J. Petykiewicz, T. M. Babinec, and J. Vučković, “Inverse design and demonstration of a compact and broadband on-chip wavelength demultiplexer,” Nat. Photonics 9(6), 374–377 (2015). [CrossRef]
5. W. Frei, D. Tortorelli, and H. Johnson, “Geometry projection method for optimizing photonic nanostructures,” Opt. Lett. 32(1), 77–79 (2007). [CrossRef]
6. L. Su, R. Trivedi, N. V. Sapra, A. Y. Piggott, D. Vercruysse, and J. Vučković, “Fully-automated optimization of grating couplers,” Opt. Express 26(4), 4023–4034 (2018). [CrossRef]
7. R. Pestourie, C. Pérez-Arancibia, Z. Lin, W. Shin, F. Capasso, and S. G. Johnson, “Inverse design of large-area metasurfaces,” Opt. Express 26(26), 33732–33747 (2018). [CrossRef]
8. L. F. Frellsen, Y. Ding, O. Sigmund, and L. H. Frandsen, “Topology optimized mode multiplexing in silicon-on-insulator photonic wire waveguides,” Opt. Express 24(15), 16866–16873 (2016). [CrossRef]
9. J. S. Jensen and O. Sigmund, “Topology optimization for nano-photonics,” Laser Photonics Rev. 5(2), 308–321 (2011). [CrossRef]
10. A. Michaels and E. Yablonovitch, “Inverse design of near unity efficiency perfectly vertical grating couplers,” Opt. Express 26(4), 4766–4779 (2018). [CrossRef]
11. T. W. Hughes, M. Minkov, I. A. Williamson, and S. Fan, “Adjoint method and inverse design for nonlinear nanophotonic devices,” ACS Photonics 5(12), 4781–4787 (2018). [CrossRef]
12. L. H. Frandsen and O. Sigmund, “Inverse design engineering of all-silicon polarization beam splitters,” in Photonic and Phononic Properties of Engineered Nanostructures VI, vol. 9756 (International Society for Optics and Photonics, 2016), p. 97560Y.
13. W. R. Frei, H. Johnson, and K. D. Choquette, “Optimization of a single defect photonic crystal laser cavity,” J. Appl. Phys. 103(3), 033102 (2008). [CrossRef]
14. S. Molesky, Z. Lin, A. Y. Piggott, W. Jin, J. Vucković, and A. W. Rodriguez, “Inverse design in nanophotonics,” Nat. Photonics 12(11), 659–670 (2018). [CrossRef]
15. D. Vercruysse, N. V. Sapra, L. Su, R. Trivedi, and J. Vučković, “Analytical level set fabrication constraints for inverse design,” Sci. Rep. 9(1), 8999 (2019). [CrossRef]
16. O. D. Miller and E. Yablonovitch, Inverse Optical Design (Springer, Berlin Heidelberg, 2015), pp. 729–732.
17. A. Hakansson, J. Sánchez-Dehesa, and L. Sanchis, “Inverse design of photonic crystal devices,” IEEE J. Select. Areas Commun. 23(7), 1365–1371 (2005). [CrossRef]
18. Z. Yu, H. Cui, and X. Sun, “Genetic-algorithm-optimized wideband on-chip polarization rotator with an ultrasmall footprint,” Opt. Lett. 42(16), 3093–3096 (2017). [CrossRef]
19. Z. Yu, H. Cui, and X. Sun, “Genetically optimized on-chip wideband ultracompact reflectors and fabry–perot cavities,” Photonics Res. 5(6), B15–B19 (2017). [CrossRef]
20. J. C. Mak, C. Sideris, J. Jeong, A. Hajimiri, and J. K. Poon, “Binary particle swarm optimized 2× 2 power splitters in a standard foundry silicon photonic platform,” Opt. Lett. 41(16), 3868–3871 (2016). [CrossRef]
21. B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm 2 footprint,” Nat. Photonics 9(6), 378–382 (2015). [CrossRef]
22. M. Zhou, B. S. Lazarov, F. Wang, and O. Sigmund, “Minimum length scale in topology optimization by geometric constraints,” Comput. Methods Appl. Mech. Eng. 293, 266–282 (2015). [CrossRef]
23. L. Su, A. Y. Piggott, N. V. Sapra, J. Petykiewicz, and J. Vuckovic, “Inverse design and demonstration of a compact on-chip narrowband three-channel wavelength demultiplexer,” ACS Photonics 5(2), 301–305 (2018). [CrossRef]
24. X. Qian, “Topology optimization in b-spline space,” Comput. Methods Appl. Mech. Eng. 265, 15–35 (2013). [CrossRef]
25. M. Wang and X. Qian, “Efficient filtering in topology optimization via b-splines,” J. Mech. Des. 137(3), 031402 (2015). [CrossRef]
26. W. Bogaerts and S. K. Selvaraja, “Compact single-mode silicon hybrid rib/strip waveguide with adiabatic bends,” IEEE Photonics J. 3(3), 422–432 (2011). [CrossRef]
27. X. Ren, W. Chang, L. Lu, M. Yan, D. Liu, and M. Zhang, “Digitized adjoint method for inverse design of digital nanophotonic devices,” arXiv preprint arXiv:1902.00654 (2019).
28. C. De Boor, “On calculating with b-splines,” J. Approx. theory 6(1), 50–62 (1972). [CrossRef]
29. M. G. Cox, “The numerical evaluation of b-splines,” IMA J. Appl. Math. 10(2), 134–149 (1972). [CrossRef]
30. D. J. Kroon, “Isocontour Webpage,” (2011). https://www.mathworks.com/matlabcentral/fileexchange/30525-isocontour.
31. O. Bernard, D. Friboulet, P. Thévenaz, and M. Unser, “Variational b-spline level-set: a linear filtering approach for fast deformable model evolution,” IEEE Trans. on Image Process. 18(6), 1179–1191 (2009). [CrossRef]
32. E. Khoram, A. Chen, D. Liu, L. Ying, Q. Wang, M. Yuan, and Z. Yu, “Nanophotonic media for artificial neural inference,” Photonics Res. 7(8), 823–827 (2019). [CrossRef]
33. E. Khoram, “Controlling the minimal feature sizes in adjoint optimization of nanophotonic devices using b-spline surfaces(code),” (2019). https://github.com/erfan6511/b-spline.
