1. Introduction

Photonic topology optimization (TO) or “inverse design” is a powerful device-design methodology in which performance (e.g. transmitted power) can be optimized over thousands or millions of degrees of freedom characterizing every “pixel” of the fabricated device, allowing completely unexpected designs with unprecedented performance to emerge [1–11]. The vast design freedom of TO must be constrained, however, by the limitations of practical manufacturing processes, and in previous TO work this has taken the form of various techniques to constrain the minimum lengthscales and curvature present in the optimized design [12–15]. In order to be useful within the commercial semiconductor foundry ecosystem, however, additional design constraints must be developed [16,17]: foundry lithography requires a minimum area and a minimum enclosed area, in addition to minimum curvature and lengthscales (linewidths and linespacings). The theoretical challenge is to formulate these constraints in a manner that is mathematically compatible with practical TO algorithms.

In this paper, we present a unified framework for density-based topology optimization which produces integrated photonic devices that obey complicated commercial foundry design rule checks (DRC). To begin with, we consolidate previous techniques to impose minimum linewidth, linespacing, and curvature constraints [13,18,19] (Sec. 3.). In addition, we present a new TO-compatible formulation of minimum-area and minimum enclosed-area constraints (Sec. 4.). Our methodology relies only on classical density-based projection methods [20] and does not require any additional reparameterization [12,21–24] (e.g. level sets or splines). Furthermore, we combine our manufacturing constraints with robust optimization [14,25,26], incorporating under/over-etch uncertainty into the TO design process (Sec. 5.), whereas previous work typically implemented robustness without manufacturing constraints [25] or vice-versa. We demonstrate our approach on three broadband silicon-photonic test problems (Sec. 6.)—dielectric waveguide bends, T-splitters, and reflectors—all optimized for broadband operation using time-domain topology optimization [27,28] developed using the free/open-source software package Meep [29].

Density-based TO methods, so named because the design parameters directly correspond to the local material properties, typically require a sequence of standard image-processing steps, such as filtering and thresholding, to map the design parameters to the physical structure, allowing the optimization algorithm to act on a differentiable function (Sec. 2.) [30]. The sensitivity of this resulting structure is efficiently calculated using an adjoint-variable method requiring just two electromagnetic simulations, regardless of the number of design variables [31]. In the time domain, a short pulse can excite a broadband response that is then used to calculate the respective sensitivities at multiple frequency points simultaneously [27,28]. Since the optimizer is gradient-based, manufacturing constraints need to be formulated in the form $g(\mbox {design}) \le 0$ where $g$ is some function that is (mostly) differentiable. Our approach follows this general framework exactly, allowing it to be directly incorporated into pre-existing TO schemes without any need to reparameterize [32] during optimization, effectively automating the design process for fabrication.

Ensuring that a final design is “DRC clean” is standard practice in the foundry ecosystem, where thousands of checks are performed to verify that the geometries can be reliably fabricated and pose no risk to surrounding designs or the manufacturing tooling [16]. While circuit-level DRC violations can be algorithmically fixed using design post-processing [33], devices designed with TO must integrate these design rules into the optimization process itself. Important “1D” constraints such as the minimum linewidth, linespacing, and curvature are well understood within the context of TO [12,13]. More complicated “2D” constraints, like minimum areas and enclosed areas, require new developments, especially since these constraints are often more restrictive than their 1D counterparts and currently are only enforceable through static B-spline parameterizations [23], limiting the problem to shape-optimization methods that require a pre-determined material topology (connectivity, number of holes, etcetera). Figure 1 illustrates each of these fundamental design rules on three example problems produced using TO. Indeed many other rules, including the maximum area, run length, and prescribed density [16] have also been overlooked in TO. Our experiments demonstrate that satisfying the five aforementioned fundamental constraints, however, often translates implicitly to satisfying these “higher order” constraints as well.

 

Fig. 1. Five fundamental design rule constraints for semiconductor foundries. 1D constraints consist of minimum linewidth, minimum linespacing and minimum curvature (a). 2D constraints include minimum area and minimum enclosed area (b). Example DRC violations in a broadband mirror, broadband bend, and broadband T-splitter are highlighted in respective colors (c).

Download Full Size | PDF

Even after DRC violations are cleared, integrated photonic devices may exhibit sensitivity to even minor process variations. Robust optimization generally means optimizing the worst case of an objective over some set of uncertainties [34]. Correspondingly, in Sec. 5 we optimize the worst of the design, an under-etched design, and an over-etched design [34]. In contrast to previous robust TO [14,25,26,34], our approach decouples the DRC constraints from over/under-etch uncertainties, meaning that the topologies between the eroded, dilated, and blueprint design fields need no longer remain consistent. Similarly, the erosions and dilations can now be applied by a simple formula, independent of the dynamic range of the underlying design field. Our methodology, validated by dozens of numerical experiments on nine unique foundry rulebooks, provides a simple and cohesive framework for practical photonic TO targeting high-yield foundry applications.

2. Photonic topology optimization overview

In this section, we review photonics TO, following Ref. [20], in the context of semiconductor foundry design and discuss various projection methods. We first formulate broadband design as a minimax optimization problem over $N$ distinct objective functions subject to Maxwell’s equations spanning $M$ frequency points:

In order to parameterize the permittivity $\boldsymbol {\epsilon }_r(\boldsymbol {\rho })$, we use a simple density-based interpolation scheme where the original (latent) design parameters $\boldsymbol {\rho }$, are first filtered using a linear density filter,

The second filter kernel is a conic filter, as the strength of the filter linearly decays toward the edge of the filter support and is defined by

Once the design variables are filtered using the appropriate kernel, either Eq. (4) or Eq. (5), the resultant field is projected onto a binary value using a differentiable and nonlinear function. In this work, we use

 

Fig. 2. Design variables, $\boldsymbol {\rho }$, evolution as they are convolved with various filter kernels, $w$, forming a filtered field, $\tilde {\boldsymbol {\rho }}$, and subsequently a projected field, $\bar {\boldsymbol {\rho }}$. The filter kernels from left to right are a uniform cylinder, a uniform ellipse rotated by 45$^{\circ }$, a uniform square rotated by 45$^{\circ }$, a conic filter, and a conic filter with twice the radius. The uniform filter kernels (marked in pink) are used for morphological transforms (Sec. 3.2) that enable systematic robust optimization (Sec. 5.). The conic filters (marked in yellow) are used for the geometric constraints (Sec. 3.1) that enforce minimum linewidth and linespacing constraints.

Download Full Size | PDF

The final permittivity is then interpolated from the design parameters using

Using this methodology, we can systematically erode and dilate the design field using two different approaches, as illustrated in Fig. 3. The first approach adjusts the value of $\eta$ used in the projection step. Dilation is produced by $\eta <0.5$ and erosion is produced by $\eta >0.5$. While often used in practice [14], the resulting erosions and dilations are directly dependent on the rate of change of the filtered design field $\widetilde {\rho }$. In other words, using the same value for $\eta$ on two different filtered profiles does not guarantee the same amount of erosion/dilation. This is inconvenient when trying to robustly design for a specific amount of over- and under-etching, as is further discussed in Sec. 5.

 

Fig. 3. Comparison of dilated (left) and eroded (right) design fields using the standard projection method (a) and harmonic filters (b). The harmonic filters can explicitly perform uniform erosions and dilations on the original design field, but are limited by the design field resolution and require an extra nonlinear filter step. In contrast, the standard projection method has more flexibility at the expense of predictability. The erosions and dilations are uniform for each contour and depend on the underlying design variable distribution.

Download Full Size | PDF

To overcome this issue, we propose a new approach which assumes $\eta =0.5$ for all elements of $\boldsymbol {\rho }$, but involves a nonlinear filter step after the projection. The harmonic erosion filter [37] is defined by

The filter $w$ is typically a uniform kernel and its support determines the amount of erosion/dilation that occurs in a particular direction. $\alpha$ is a regularization parameter to ensure differentiability; the exact dilation/erosion operation is obtained in the $\alpha \to 0$ limit, and we find that $\alpha = 10^{-3}$ is sufficient for usage in TO. For example, a filter support defined by a circle with a radius of 20 nm will systematically erode/dilate the design by 20 nm uniformly in all directions. The advantage of this approach is that it performs the same morphological transformation regardless of the rate of change of the underlying filtered field, as the nonlinear filters are applied after the projection step (and the resulting structures are already quasi-binary).

3. Minimum linewidth, spacing, and curvature constraints

Here, we consolidate existing methodologies for lengthscale constraints from several papers [13,18,26,36] to ensure that our designs obey minimum-linewidth, minimum-linespacing, and minimum-curvature requirements.

3.1 Geometric constraints

Reference [13] introduced two geometric constraints that enforce a minimum lengthscale on both the solid and void regions of a topology. These particular constraints correspond exactly to the minimum-linewidth and minimum-linespacing constraints of a foundry’s design rulebook.

The constraints work by identifying inflection regions that dictate whether or not the design topology remains unchanged after a strategically chosen set of erosions and dilations. If the topology is consistent (i.e. no additional “holes” or “islands”), then the lengthscale is satisfied. Enforcing an optimization constraint that implements this process requires an indicator function capable of identifying the inflection regions and a known relationship between the filter kernel and the necessary erosions/dilations.

The minimum-linewidth constraint, $g_{LW} \le 0$, is described by the function

To determine $\eta _e$ and $\eta _d$ from the foundry minimum linewidth ($l_w$) and minimum spacing ($l_s$), the following relations hold [26]:

For the damping coefficient $c$, we use $c=r^{4}$ where $r$ is the design grid resolution (i.e. $\frac {1}{\Delta x}$). As we discuss in Sec. 6., we impose the constraint gradually during optimization by employing a constraint $g_k \le G_k$, where $G_k$ is reduced in a systematic way as the optimization progresses. The choice of tolerances $G_k$ allows some flexibility in the tradeoff between the geometric constraints and the area constraints.

Finally, we note that both the linewidth and linespacing constraints also place implicit constraints on the corresponding minimum curvature for both solid and void regions. Assuming a circular filter is used and that the constraint is satisfied, the resulting topology cannot contain circular elements whose diameters are smaller than the corresponding lengthscale (e.g. linewidth for solid regions and linespacing for void regions). Consequently, the corresponding minimum radius of curvature $\kappa _{w,s}$ is defined by

In some cases it may be desirable to place a smaller radius of curvature constraint than is dictated by the minimum linewidth or spacing (e.g. for platforms that are able to resolve sharper angles). In such cases, a morphological transform may be a better approach.

3.2 Morphological-transform constraints

As an alternative to the geometric constraints described above, it has been proposed that one can enforce minimum linewidth, linespacing, and curvature constraints using differentiable morphological transforms [18]. In practice, we found these morphological constraints to be inconvenient for generating binary designs as described below. However, we review them here because we found them to be useful, instead, in implementing robust optimization as described in Sec. 5.

First, from the erosion and dilation operators of Eqs. (8–9), one defines an open operator ($\mathcal {O}_{\mathcal {N}}$)

While elegant, this approach still requires a suitable constraint function that is differentiable. One obvious implementation might be an $L^{2}$ norm

4. Minimum area and enclosed area constraints

We now define new minimum-area and minimum-enclosed-area constraints along with their respective gradients. The basic procedure for enforcing each design rule follows the methodology of the geometric constraints rather closely. First, indicator functions identify regions that violate the minimum-area constraint (“islands”) and regions that violate the minimum-enclosed-area constraint (“holes”). Next, each constraint will encourage the optimizer to either eliminate the violating regions via a localized erosion (for islands) or a localized dilation (for holes). Figure 4 illustrates the two constraint functions on a sample design with their corresponding gradients.

 

Fig. 4. Minimum area and minimum enclosed area constraint functions. (a) The nominal design under consideration is first parsed using an image segmentation algorithm to identify both the (b) islands and (e) holes that correspond to regions violating either constraint (highlighted in red). These contour masks serve as indicator functions for the optimization constraint. Gradients for the (c) minimum area constraint and the (f) minimum enclosed area constraints are computed by backpropagating through the rest of the mapping functions using the chain rule. The minimum area constraint seeks to eliminate violating regions, which corresponds to a negative gradient (red). The minimum hole constraint seeks to fill holes that are too small, which corresponds to a positive gradient (blue). All images were generated using resolutions that match the corresponding simulation resolution, which corresponds to small, pixelated artifacts in the design geometry, field patterns, and constraint gradients.

Download Full Size | PDF

We define the new minimum-area constraint function, $g_A$, as

To determine which regions (if any) violate the area and enclosed area rules, we use the marching-squares algorithm [38], a versatile image-segmentation technique often used to extract contours from images. The algorithm relies on a user-defined threshold point (between 0 and 1) such that “clusters” of pixels above the threshold point form part of the contour. In many ways, the procedure follows the same logic as the filter-threshold mapping used to evolve the design’s parameter field during optimization. For simplicity, we used an out-of-the-box implementation in Python [39] and a threshold parameter of 0.6.

Once the contours are identified, we calculate the area of each contour using a discrete summation of all the density values inside the contour, identified using morphological dilations. We used another out-of-the-box implementation in Python [40] to quickly and accurately compute each region’s area. If a particular contour has an area smaller than the prescribed design rule, the filled contour region is dilated by 1 pixel (e.g. using Eq. (9)) and subsequently added to the indicator function. As explained later, the minimal dilation conditions the gradient computation.

The constraints are then defined by

The constraint gradients w.r.t. the filtered design parameters are calculated as

To calculate the gradients w.r.t. the actual design parameters, we use the chain rule and backpropagate [30] through the thresholding and filtering steps

At first glance, it may seem that the indicator functions $I_A$ and $I_{EA}$ pose a severe challenge, because they are complicated to compute and are not always differentiable. However, it turns out that this difficulty is solved by the geometric constraints from Sec. 3., which ensure that the “islands” are well separated. In consequence, we can choose the boundaries of the indicator function’s support so that its boundaries occur only where $\bar {\rho }$ (for $I_A$) or $1-\bar {\rho }$ (for $I_{EA}$) is exponentially small, and hence changes in the indicator function make a negligible contribution to the constraint gradients. By slightly dilating each violating contour within the indicator function, we avoid any singularities associated with the edges of each contour. We can therefore safely omit this term from our analysis and still produce sufficiently accurate gradients. The gradients used for all numerical experiments in Sec. 6 are calculated using

Occasionally, the minimum-lengthscale constraints and the minimum-area constraints compete, effectively stalling the optimization. We overcome this process by weighting the two sets of constraints within the minimax optimization problem, as further explained in Sec. 6. This new formulation allows the optimizer to move freely between feasible regions without the usual issues that accompany stiff optimization (where the corresponding Hessian is ill-conditioned).

5. Robust optimization

In this section, we present a comprehensive approach to robust optimization that includes geometric constraints, area constraints, and an explicit under/over-etch tolerance. Commercial foundry lithography introduces yield variability, both random and systematic, throughout each step of the fabrication process. Current density-based robust optimization methodologies account for over/under-etch variability by running a worst-case optimization across three separate design fields [14,25,34]. The eroded and dilated design fields (which are geometric representations of over and under etching, respectively) are typically determined by choosing threshold parameters ($\eta$) that implicitly enforce quasi-lengthscales [14]. This approach can be problematic because all three design fields (eroded, dilated, and the nominal “blueprint”) must have the same topology [34], a strong requirement that may not hold in TO with small features. While some double-filtering techniques have shown potential for alleviating this restriction [15], predicting the actual amount of over- and under-etching relative to the desired lengthscale constraints is also difficult as it strongly depends on the dynamic range of the latent design parameters ($\rho$). In short, a reliable and deterministic method that decouples the lengthscale constraints from the prescribed etch tolerance has not yet been demonstrated.

We propose a straightforward approach that obeys arbitrary lengthscale and area design constraints, without any restrictions on the amount of over/under etching or on the topology itself. First, we place geometric constraints and minimum-area constraints on just the blueprint design. To ensure a predictable over/under-etch perturbation, we use the morphological transforms described in Sec. 2., such that the final erosions and dilations do not require any hyperparameter tuning (e.g. $\eta$) as they no longer depend on the dynamic range of the latent design field (illustrated in Fig. 3). Furthermore, this technique allows the topology to change between design fields without compromising the minimum length scales or requiring intermediary filter steps. The design rule constraints and robustness to manufacturing uncertainty are now completely decoupled, allowing the user to freely activate one or the other (although typically both will be applied).

Figure 5 compares a non-robust splitter to a robust splitter designed to be tolerant to roughly $\pm 20$ nm of etch variation. Both splitters were optimized with the geometric constraints and area constraints described above. The robust variant requires two extra parallel Maxwell solves per device. As expected, the nominal design between the robust and non-robust variant show topological differences. The broadband performance (indicated in errorbars) for the robust splitter demonstrates greater tolerance than the non-robust variant even beyond the expected $\pm 20$ nm variations.

 

Fig. 5. Comparison of a robust T-splitter to a nominal T-splitter. (a) The nominal, non-robust design. (b) The broadband performance as a function of etching variation for both the robust and nominal design. The (c) eroded, (d) blueprint, (e) and dilated designs for the robust optimization. The final pixelated device geometries were smoothed (but not otherwise altered) to better illustrate the erosions and dilations

Download Full Size | PDF

The devices were designed on a grid with 17 nm resolution. The resulting $\pm 20$ nm target variation, along with the additional $\pm 30$ nm and $\pm 40$ nm perturbations, were approximated using the nonlinear harmonic filters described in Sec. 3.2. Like any spatial parameter, these values are only approximate once the grid discretization is included, but our 17 nm resolution is sufficient to clearly resolve distinct $\pm 20$ nm, $\pm 30$ nm, and $\pm 40$ nm perturbations (filter radii of 1.2, 1.8, and 2.4 pixels), as shown in Fig. 5. Moreover, the point of performing robust optimization is that the results become less sensitive to small details such as the discretization resolution or the exact perturbation amount.

Note that the nominal design here actually has a certain degree of robustness simply by virtue of being broadband, which eliminates the extreme sensitivities that can arise in single-frequency optimization due to resonant and interference effects [41]. Therefore, it is not surprising that in Fig. 5 the nominal design still performs some splitting even for over/under-etching, and is qualitatively not too unlike the robust design. However, because the robust design is designed specifically to tolerate over/under-etching, its performance is quantitatively much better (its error bars are up to $24 \times$ smaller for variations $\le 20$ nm) for that specific imperfection.

While this approach produces robust devices that simultaneously obey length scales, further work is needed to extend this methodology to random perturbations on the design field. Indeed existing methods that randomly vary the thresholding parameter $\eta$ could be used, but they suffer from the same challenges discussed earlier.

6. Numerical examples

In this section, we present numerical results for three broadband silicon photonic examples to demonstrate the effectiveness and flexibility of our proposed methodology. Among the many photonic TO examples already discussed in the literature (1D grating couplers, demultiplexers, directional couplers, etc.) we chose (2D) devices that typically produce small islands and holes. Specifically, we demonstrate the design of several different waveguide reflectors, bends, and T-splitters because they each route light within a highly confined space (3 $\mu$m $\times$ 3 $\mu$m). Each example uses the geometric and area constraints described earlier for nine different foundry specifications (for twenty-seven unique examples in total). While these specific design rules are synthetic, they are comparable to existing commercial-foundry DRC constraints.

Each example problem has an objective function depending on a waveguide modal overlap, which is proportional to an entry of the scattering matrix ($S$ matrix):

The fabrication constraints must be enforced “gradually” in the optimization problem, so that at early stages of the optimization the geometry can change freely (since changes in topology, such as the appearance of a new hole, temporarily violate the constraints [13]. For a fabrication constraint $g_k \le 0$, this is accomplished by changing the constraint to $g_k \le G_k$ for some constants $G_k \ge 0$ that are decreased as optimization proceeds:

In many photonics-optimization problems, such as the examples here, the objective can be formulated so that $t \to \max _n{f(\omega _n)} \to 0$ as the optimization proceeds, and this can be used to simplify the choice of constraint bounds $G_k$. For example, if we are minimizing reflection or $1-\textrm {transmission}$, then a high-performance structure will an achieve near-zero optimum. In such cases, we can simply set $G_k = a_k t$ for a suitable constant $a_k$, so that $G_k \to 0$ as the optimization proceeds. Here, $a_k=10^{-5}$ for $k \in \{LS,LW,A,EA\}$ worked for all the examples presented below. In the future, we are hopeful that such hyperparameters as $G_k$ and the projection scale factor $\beta$ can be automatically inferred from the convergence rate of the optimization algorithm. For now, however, they offer a strategic “knob” to turn when the optimization stalls or underperforms (further discussed in Sec. 6.3).

For the forward and adjoint solves, we used Meep, a free/open-source FDTD software package [29]. Each simulation ran with a resolution of 30 pixels/$\mu$m ($\Delta x=33$ nm) and a design-region ($\rho$) resolution of 60 pixels/$\mu$m ($\Delta x=17$ nm). Solid regions represent silicon with a dispersionless refractive index of $n=3.4$. Void regions represent SiO$_2$ with a dispersionless refractive index of $n=1.44$. We note that our method generalizes to materials with dispersion, and that dispersionless approximations were used here for simplicity. All examples were optimized using 70 iterations at $\beta =8,16,$ and $32$ for a total of 210 iterations. The foundry constraints $g_k \le G_k$ were not included until $\beta =32$. We used a free/open-source implementation [42] of the globally convergent method of moving asymptotes (MMA) [35] optimization algorithm, restarting the algorithm each time $\beta$ changed values. We ran each simulation on four cores of an Intel Xeon Gold 6226 2.7 GHz CPU using resources provided by the Partnership for an Advanced Computing Environment (PACE) at the Georgia Institute of Technology. Each example finished after about 6–8 hours of runtime.

6.1 Mirror

In this example, we designed a silicon photonic broadband mirror such that the forward-propagating mode is reflected, transferring energy to the backward-propagating mode of the same waveguide. Figure 6 illustrates the design evolution of the broadband mirror, along with the steady state field pattern of the final device at the center of the band and the broadband reflection response. It is clear that a topology with several small holes and islands is favored by the optimization algorithm until the area and geometric constraints are activated during the last thresholding epoch (iteration 140, $\beta =32$). At this point, the performance drops significantly and is gradually recovered while the foundry constraints are simultaneously satisfied. The final device exhibits no small holes or islands and obeys the prescribed linewidth, linespacing, and curvature rules, while still performing well across the band with $<1$% variation in reflectivity.

 

Fig. 6. Evolution and performance of the broadband mirror designed to obey a minimum linewidth of 90 nm, a minimum spacing of 90 nm, a minimum area of 0.08 $\mu$m$^{2}$, and a minimum enclosed area of 0.2 $\mu$m$^{2}$. (a) Projected design field ($\bar {\rho }$) evolution for iterations 1, 31, 81, 146, and 210. The final iteration (210) also displays the steady-state field pattern of the device at $\lambda$=1.55 $\mu$m. The two iterations enclosed in red indicate when the constraint was active at $\beta =32$. (b) The evolution of the device broadband reflection as a function of iteration. The mean performance across the band is depicted by a solid line, while the bounds are illustrated by a blue shaded region. All constraints were enforced starting with iteration 141. (c) The broadband performance of the device from 1.50 $\mu$m to 1.60 $\mu$m.

Download Full Size | PDF

The objective function used to design the mirror is

Figure 7 illustrates the final designs: nine different broadband mirrors, each of which obey a different combination of fabrication constraints. The devices are presented such that those with the least-restrictive constraints are on the upper-left side of the figure and the constraints become more restrictive toward the bottom-right. The least-restrictive rulebook, for example, imposes a minimum linewidth/linespacing of 60 nm and a minimum area/enclosed area of 0.005 $\mu$m$^{2}$, and produced a device shown in the upper-left corner. In contrast, the design shown in the lower-right corner obeys a minimum linewidth/linespacing of 120 nm and a minimum area/enclosed area of 0.2 $\mu$m$^{2}$. The devices located between these two boundary cases are intermediate versions of these design rules, including cases where linewidth $\ne$ linespacing and/or area $\ne$ enclosed area.

 

Fig. 7. Comparison of various broadband mirror designs optimized to perform under different combinations of fabrication constraints. (top) Each of the minimum linewidth, minimum spacing, minimum area, and minimum enclosed area constraints drawn to scale. The minimum area constraints, for example, are demonstrated by a circle that shares the same area as the constraint. The minimum linewidth and spacing constraints are demonstrated by circles whose diameters align with the corresponding constraint. By extension, the minimum radius of curvature also aligns with this depiction. (bottom) Each broadband mirror design is drawn alongside the constraints used to produce it. The average reflection of each device from $\lambda =$1.5 $\mu$m to $\lambda =$1.6 $\mu$m, along with the corresponding minimum and maximum is depicted above each design.

Download Full Size | PDF

All but two of the resulting devices reflect $>90$% of the light with very little variation across the band (the most restrictive rulebook produced a device with just $\pm 0.3\%$ of variation across the band). Because these algorithms produce only a local optimum, not a global optimum, the devices that have less restrictive constraints (upper-left) do not always outperform the devices designed with tighter restrictions (lower-right). For example, the device with a linewidth/linespacing constraint of 90 nm and a minimum area/enclosed area constraint of $0.005\mu$m$^{2}$ only achieve an average of 68.8% transmission with a rather high variability of 18.7% across the band. A global optimizer subject to this particular rulebook could have also produced the same design as the device on its right, as it shares the same linewidth/linespacing constraint (and thus the same filter radius) but with slightly stricter area constraints (which still satisfy a less-strict constraint by definition). However, as seen here, a local optimizer may happen to find a poorer local optimum in a case with looser constraints, and it is well known that some experimentation (e.g. different starting points or different constraint weights) is often needed to find a satisfactory local optimum in high-dimensional nonconvex optimization (see Sec. 6.3). Nevertheless, as is often observed in topology optimization [20], in all cases the algorithm converged to a local optimum with good performance.

6.2 Bend

In the next example, we design a silicon photonic broadband $90^{\circ }$ bend in a single-mode waveguide. Figure 8 illustrates the design evolution of the broadband bend, along with the steady-state field pattern of the final device at the center of the band and the broadband transmission response. As expected, the design converges to a binary layout that satisfies the prescribed DRC rules. Unlike the mirror examples earlier, the optimization favored larger, elongated structures/gaps that satisfied the area rules even before the constraints were applied. As such, the drop in performance at the start of the final epoch is only about 10%.

 

Fig. 8. Evolution and performance of the broadband bend designed to obey a minimum linewidth of 90 nm, a minimum spacing of 90 nm, a minimum area of 0.08 $\mu$m$^{2}$, and a minimum enclosed area of 0.2 $\mu$m$^{2}$. (a) Projected design field ($\bar {\rho }$) evolution for iterations 1, 31, 81, 146, and 210. The final iteration (210) also displays the steady-state field pattern of the device at $\lambda =$1.55 $\mu$m. The two iterations enclosed in red indicate when the constraint was active at $\beta =32$. (b) The evolution of the device broadband reflection as a function of iteration. The mean performance across the band is depicted by a solid line, while the bounds are illustrated by a blue shaded region. All constraints were enforced starting with iteration 141. (c) The broadband performance of the device is from $\lambda =$1.50 $\mu$m to $\lambda =$1.60 $\mu$m.

Download Full Size | PDF

The objective function used to design the bend is

Figure 9 shows the final designs for nine different broadband bends, each of which obeys a different combination of fabrication constraints (just as for the mirror designs above). Again, the devices perform well across the various rulebooks, with the average transmission always $>89\%$. We note that the most constrained device (lower-right corner) exhibits a transmission of 98.2% with practically no variation across the band. Interestingly, the area constraints removed the small Bragg mirror-like features on the boundaries as soon as the constraints were activated. The resulting device naturally took a rather different optimization evolution than the other devices constrained by different rulebooks.

 

Fig. 9. Comparison of various broadband bend designs optimized to perform under different combinations of fabrication constraints. (top) Each of the minimum linewidth, minimum spacing, minimum area, and minimum enclosed area constraints drawn to scale. The minimum area constraints, for example, are demonstrated by a circle that shares the same area as the constraint. The minimum linewidth and spacing constraints are demonstrated by circles whose diameters align with the corresponding constraint. By extension, the minimum radius of curvature also aligns with this depiction. (bottom) Each broadband mirror design is drawn alongside the constraints used to produce it. The average transmission of each device from $\lambda =$1.5 $\mu$m to $\lambda =$1.6 $\mu$m, along with the corresponding minimum and maximum is depicted above each design.

Download Full Size | PDF

6.3 T-splitter

In the final example, we design a symmetric, silicon-photonic broadband T-splitter. Unlike conventional splitters that route the outputs along the same direction as they entered the splitting structure, the outputs of this design are each routed at a $90^{\circ }$ angle (forming a "T"). Consequently, we intuitively expect structures containing small islands and holes in order to reflect the light within the device and route it as intended. Figure 10 illustrates the design evolution of the broadband splitter, along with the steady-state field pattern of the final device at the center of the band and the broadband transmission response. Similar to the mirror example above, we first see several small islands forming in the initial two epochs ($\beta =8,16$). Once the constraints are applied, however, these violating regions are modified such that they no longer pose any DRC issues. Specifically, the geometric constraints “bridged” the regions that were within the rulebook’s defined linespacing, while the area constraints eliminated the remaining well-separated islands. The enclosed-area constraints ensured that no new small holes would form as a result of this bridging process. The resulting device exhibits $~0.2\%$ of variation across the band with a mean transmission of 49%.

 

Fig. 10. Evolution and performance of the T-splitter designed to obey a minimum linewidth of 90 nm, minimum spacing of 90 nm, a minimum area of 0.08 $\mu$m$^{2}$, and minimum enclosed area of 0.08 $\mu$m$^{2}$. (a) Projected design field ($\bar {\rho }$) evolution for iterations 1, 30, 80, 145, and 210. The two iterations enclosed in red indicate when the constraint was active at $\beta =32$. (b) The evolution of the device splitting ratio as a function of iteration. All constraints were enforced starting with iteration 140. (c) The steady state field pattern of the final device at $\lambda$=1.55 $\mu$m. (d) The broadband performance of the device from $\lambda$=1.50 $\mu$m to $\lambda$=1.60 $\mu$m.

Download Full Size | PDF

The objective function used to design the T-splitter is

We imposed a geometric symmetry constraint, forcing the final splitting ratio to also be symmetric.

Figure 11 illustrates the final designs for nine different broadband splitters, each of which obeys a different combination of fabrication constraints (just as for the mirror and bend designs above). In this case, we see that the least-restricted device (upper-left corner) happens to exhibit rather poor performance across the band ($27.2\% \pm 5.7\%$), which turns out to simply be an unlucky local optimum. We reran the exact same optimization problem for this rulebook, but modified the constraint weight factors to $a_k=10^{-3}$, effectively relaxing the strength of each constraint (without modifying the constraint function itself). The resulting device still obeyed the prescribed design rule constraints, but exhibited far superior performance across the band ($49.5\% \pm 0.05\%$). In general, we found that if the constraint weights were too strict, the optimizer could become stuck in a local minimum. The constraint weights ($a_k$) offer additional tunability to combat slow convergence.

 

Fig. 11. Comparison of various broadband T-splitter designs optimized to perform under different combinations of fabrication constraints. (top) Each of the minimum linewidth, minimum spacing, minimum area, and minimum enclosed area constraints drawn to scale. The minimum area constraints, for example, are demonstrated by a circle that shares the same area as the constraint. The minimum linewidth and spacing constraints are demonstrated by circles whose diameters align with the corresponding constraint. By extension, the minimum radius of curvature also aligns with this depiction. (bottom) Each broadband mirror design is drawn alongside the constraints used to produce it. The average transmission of each device from $\lambda =$1.5 $\mu$m to $\lambda =$1.6 $\mu$m, along with the corresponding minimum and maximum is depicted above each design.

Download Full Size | PDF

7. Conclusion

We presented a comprehensive photonic TO design methodology that includes minimum linewidth, linespacing, curvature, area, and enclosed-area “DRC” constraints. The linewidth, linespacing, and curvature constraints are adapted from previous work and combined with novel area constraints. All the constraints rely on density-based projection steps without any additional re-parameterizations. By additionally incorporating differentiable morphological transforms, we provide a method to perform systematic robust optimization while simultaneously imposing the fabrication constraints. We illustrated these techniques using three distinct test devices constrained by nine different design rulebooks.

There remain many opportunities for future work on these topics. As we noted in Sec. 4., it is desirable to have a more flexible area constraint that allows the topology-violating region to either expand or contract, instead of only contracting as in this work. Similarly, a more rigorous analysis that relates the hyperparameters of the geometric constraints ($G_k$) is needed, such that the need to tune each optimization run could be reduced. One possible path forward involves deriving an improved constraint function based on the morphological transforms described in Sec. 3.2. While our new robust-optimization paradigm decouples systematic under/over-etch variations from DRC constraints, future efforts should seek to handle a larger class of manufacturing uncertainties, such as over-etch in some regions and under-etch in others, and surface roughness, within a similar morphological-transform framework. Similarly, our methodology does not explicitly incorporate calibration steps that are typically performed by the foundry, such as optical proximity correction (OPC). While recent work has demonstrated a variant of TO that uses filters to mitigate the need for OPC [43,44], most semiconductor foundry pipelines (and their respective DRC rulebooks) are not set up for this workflow.