Abstract

In part I of this series of papers, we introduced an efficient pixel-by-pixel optimization technique, utilizing the Dyson’s equation in a Green’s function formalism. This technique directly incorporates material and minimum feature size constraints in practical photonic devices. The implementation in part I, however, requires full storage of the Green’s function, which is costly. In this paper, we further develop this optimization technique by showing that only the diagonal part of the Green’s function needs to be stored, and by showing that the optimization process can be implemented using standard commercially available electromagnetic solvers. The development here should enable a wider use of this optimization technique.

© 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. INTRODUCTION

For on-chip photonic device optimization, there are a few important constraints that need to be satisfied. First, these devices are typically constructed out of a discrete set of materials with a corresponding discrete set of permittivities, and the final optimized design must use only this set of permittivities. In many cases, this set in fact consists of two materials, and the resulting design must be a binary design [1]. Second, the design must be amenable to nanofabrication technology; in particular, the design needs to satisfy the minimum feature size requirement.

While continuous optimization methods, in particular those based on the adjoint variable method [28], have been widely used for the optimization of photonic devices, it remains a considerable challenge to incorporate the constraints mentioned above, in a systematic fashion, while maintaining computational efficiency and high performance of the final design [9]. Discrete optimization methods, e.g., a pixel-by-pixel optimization strategy, can directly incorporate these constraints from the outset, and in fact can be implemented such that the search occurs only in the parameter space that satisfies these constraints; however, the computational cost of pixel-by-pixel optimization has been high [10,11]. In particular, in traditional implementation of pixel-by-pixel optimization, where one updates a single pixel in each iteration, one needs to simulate the performance of all single-pixel variations in order to accurately decide the optimal pixel to update in each iteration, leading to substantial computational costs.

In part I of this series of papers [12], we have introduced a method, based on Dyson’s equation [13] in a Green’s function formalism, to overcome the challenge associated with discrete pixel-by-pixel optimization. We show that starting from the Green’s function of a device configuration, the change in the figure of merit (FOM) for all single-pixel variations can be efficiently computed. This then allows us to update the single pixel that leads to greatest improvement of the FOM to arrive at a new device configuration. Moreover, the Green’s function of the new device configuration can also be efficiently determined. Thus, the process can continue iteratively until an optimal configuration is reached. In part I [12], we provide an implementation of this method using the discrete dipole approximation (DDA) [1315] as the simulation tool.

There are several drawbacks in the implementation of part I [12]. First of all, the implementation requires storage of the Green’s function. Thus the memory requirement for storage could be substantial: the Green’s function in real space is fully dense. Thus, the memory requirement scales as N2, where N is the number of pixels in the design region. This can become prohibitive if the design region is large. Second, many standard electromagnetic (EM) solvers that are commercially available are based on methods such as finite-difference time-domain (FDTD) [16,17], or iterative finite-element methods [18,19]. These methods efficiently compute the field distribution for a given source configuration in a structure, but usually do not compute the Green’s function. To use these methods to compute the Green’s function directly would require the calculation for a large number of source configurations, which can become very expensive.

In this paper, we show that the basic idea for efficient pixel-by-pixel optimization, as outlined in part I, can alternatively be implemented with the use of standard commercially available electromagnetic solvers, with a memory requirement that scales linearly with N. For this purpose, we note that to determine the position of the pixel, the modification of which will lead to the best improvement in the FOM, one requires only the diagonal part of the Green’s function, the storage of which scales linearly with N. The diagonal part of the Green’s function, in turn, can be updated with only a few full simulations of the structures. This idea can be generalized for an efficient block-by-block optimization, and therefore enables the direct incorporation of fabrication constraints.

The rest of the paper is organized as follows: in Section 2, we discuss the algorithms. In Section 3, we illustrate the efficiency of this algorithm with a few optimization examples. We conclude in Section 4.

2. MATHEMATICAL BACKGROUND

A. Single Pixel Modification

The same as in part I [12], we consider an electric field E(r) generated by an excitation for a structure as described by a given permittivity ε(r). For illustration purposes, we suppose we would like to optimize a simple figure of merit that consists of the intensity of the electric field at a target point r0:

FOM=|E(r0)|2.
The optimization then typically requires one to evaluate the following change in the FOM between two successive iterations:
ΔFOM=|E(new)(r0)|2|E(old)(r0)|2,
where E(old)(r) and E(new)(r) are the “old” and “new” electric fields, corresponding to the permittivity distribution ε(old) and ε(new) before and after the iteration, respectively.

Define the change in permittivity distribution Δεε(new)ε(old) between two iterations. Using the Lippmann–Schwinger equation [15], we can determine the new electric field in the whole space:

E(new)(r)=E(old)(r)+Vk02G(old)(r,r)Δε(r)E(new)(r)dr,
where k0 is the free-space wavevector, and G(old) is the electric Green’s functions of the system as described by ε(old). The integration is over only the volume V, where Δε0, i.e., only in the volume where there is variation of the permittivity between two successive iterations.

In the case where we consider a single-pixel modification in each iteration, Δε is non-zero only in a small volume V centered at r, and the change in FOM consists of

ΔFOM=|E(old)(r0)+k02VG(old)(r0,r)Δε(r)E(new)(r)|2|E(old)(r0)|2.
Applying Lorentz reciprocity,
G(old)(r0,r)=G(old)(r,r0)T,
and introducing the adjoint field Eadj(old)(r), which is excited by a dipole of the amplitude E(old)(r0)* located at r0, for the permittivity distribution ε(old), i.e.,
Eadj(old)(r)=G(old)(r,r0)E(old)(r0)*,
the change in FOM can then be expressed, neglecting the second-order term, as [12]
ΔFOM=2k02VRe(Eadj(old)(r)·Δε(r)E(new)(r)).
The approximation leading to Eq. (7) is justified as long as the size of the volume V is much smaller than the wavelength.

As we show in part I [12], the new electric field at r after a pixel has been updated in r can be efficiently computed, provided that the old Green’s function at this location is known, according to the Lippmann–Schwinger equation [12,13,15] applied in r:

E(new)(r)=E(old)(r)+k02VG(old)(r,r)Δε(r)E(new)(r).
Therefore, using (7) and (8), the best pixel to update can be rigorously and efficiently determined even for large and discrete permittivity changes, as required in binary systems [12]. We refer to this pixel as the “best-improvement pixel” below.

In the DDA-type approach described in part I [12], for the “new” structure where the dielectric constant is modified at the best-improvement pixel, the location of which is denoted as rb, the field was updated everywhere using the Lippmann–Schwinger equation:

E(new)(r)=E(old)(r)+k02VG(old)(r,rb)Δε(rb)E(new)(rb).
The full Green’s function was also updated, for use in the next iteration, by Dyson’s equation [12]:
G(new)(r,r)=G(old)(r,r)+k02VG(old)(r,rb)Δε(rb)G(new)(rb,r).
We solved Eq. (10) first for G(new)(rb,r):
G(new)(rb,r)=G(old)(rb,r)+k02VG(old)(rb,rb)Δε(r)G(new)(rb,r).
Once G(new)(rb,r) was determined, the Green’s function of the new structure G(new)(r,r) for an arbitrary pair of coordinates r and r was then determined using Eq. (10).

In Ref. [12], when implementing the algorithms above, we stored the full Green’s functions G(r,r) in each iteration process. When using a numerical solver based on the discrete dipole approximation, this approach is natural. On the other hand, for most other algorithms, this is costly, since G(r,r)is fully dense. Thus, it is useful to develop a scheme in which one needs not store the full Green’s function. From Eqs. (7) and (8), it is clear that only the diagonal part of the Green’s function G(old)(r,r), where r can vary over the entire design region, is needed in the iteration process in order to determine the best improvement pixel. The required storage for G(old)(r,r) scales as N, where N is the number of pixels in the design region. Also, instead of using Eq. (9), which involves the full Green’s function G(r,r) to update the field, the field for the new structure can alternatively be determined by using a standard EM solver, without the need of knowing the Green’s function of the old structure. Moreover, within each iteration, the diagonal part of the Green’s function can be updated for the new structure. Indeed, from Eq. (10), we have

G(new)(r,r)=G(old)(r,r)+k02VG(old)(r,rb)Δε(rb)G(new)(rb,r).
To determine G(new)(rb,r), we utilize the reciprocity theorem:
G(new)(rb,r)=G(new)(r,rb)T.
G(old)(r,rb) and G(new)(r,rb) can in principle be determined using the standard EM solver, by placing appropriate sources at rb, respectively, in the old and new structures [15].

In practice, the number of simulations required to determine G(old)(r,rb) and G(new)(r,rb) can be significantly reduced, by taking advantage of the simulations of the fields that have been already performed.

To start, we show a more efficient way to determine G(old)(r,rb).

In the case of transverse electric (TE) polarization in two dimensions, G(old)(r,rb) is scalar, and therefore Eq. (9) is sufficient to determine it. The required information is: E(old)(r) before the best-improvement pixel in rb is updated, and E(new)(r) and E(new)(rb) after the best-improvement pixel is updated. All of this information has already been calculated. Thus, G(old)(r,rb) can be straightforwardly determined by solving a simple linear equation:

E(new)(r)=E(old)(r)+k02VG(old)(r,rb)Δε(rb)E(new)(rb).
There is no longer a need for a separate simulation to determine G(old)(r,rb).

In the case of the transverse magnetic (TM) polarization in two dimensions, the Green’s function is a 2×2 tensor. Utilizing both E(new)(r) as well as Eadj(new)(r), which will be needed in order to compute the gradient of the FOM in Eq. (7), we have

{E(new)(r)=E(old)(r)+k02VG(old)(r,rb)Δε(rb)E(new)(rb),Eadj(new)(r)=Eadj(old)(r)+k02VG(old)(r,rb)Δε(rb)Eadj(new)(rb).
Thus, in this case, G(old)(r,rb) can be determined by solving a 4×4 linear system. (Notice that all E fields are two vectors in this case.)

Finally, for 3D problems, where electric fields are three dimensional and hence G(old)(r,rb) is a 3×3 tensor, one has to run a third simulation (that we will denote with a subscript Green), in addition to forward and adjoint simulations. This simulation can be done with any excitation source as long as such an excitation source is linearly independent of the sources used in the simulations for E(new)(r) and Eadj(new)(r). We can use a dipole source located in rb, and with the dipole moment orthogonal to both E(old)(rb) and Eadj(old)(rb), which provides a linear system with good conditioning [15]. The Lippmann–Schwinger equations for the three new fields then read

{E(new)(r)=E(old)(r)+k02VG(old)(r,rb)Δε(rb)E(new)(rb),Eadj(new)(r)=Eadj(old)(r)+k02VG(old)(r,rb)Δε(rb)Eadj(new)(rb),EGreen(new)(r)=EGreen(old)(r)+k02VG(old)(r,rb)Δε(rb)EGreen(new)(rb),
which defines a 9×9 linear system from which G(old)(r,rb) can be determined.

We now switch to the determination of G(new)(r,rb). We note that Eqs. (14)–(16) are all in the form of the Lippmann–Schwinger equation, which has the following properties between the “old” and the “new” fields: one can simply switch the subscript “new”, and “old,” and replace Δε(rb) with Δε(rb), to arrive at a new set of equations that involve G(new)(r,rb). As a result, Eq. (14) transforms to

E(new)(r)=E(old)(r)+k02VG(new)(r,rb)Δε(rb)E(old)(rb).
Equations (15) and (16) transform to
{E(new)(r)=E(old)(r)+k02VG(new)(r,rb)Δε(rb)E(old)(rb),Eadj(new)(r)=Eadj(old)(r)+k02VG(new)(r,rb)Δε(rb)Eadj(old)(rb),
{E(new)(r)=E(old)(r)+k02VG(new)(r,rb)Δε(rb)E(old)(rb),Eadj(new)(r)=Eadj(old)(r)+k02VG(new)(r,rb)Δε(rb)Eadj(old)(rb),EGreen(new)(r)=EGreen(old)(r)+k02VG(new)(r,rb)Δε(rb)EGreen(old)(rb).
G(new)(rb,r) can then be determined by solving the small linear systems in Eqs. (17)–(19).

We therefore arrive at the following algorithm: at the start of each iteration, we assume that we know the field distribution E(old)(r), the adjoint field distribution Eadj(old)(r), as well as the diagonal part of the Green’s function G(old)(r,r), where r takes the value over the design region. We use Eqs. (7) and (8) to determine the best-improvement pixel at which the modification of the structure leads to the best improvement of the FOM. We then update the structure by modifying the best-improvement pixel, arriving at the new structure. Afterwards we determine the field distributions E(new)(r) and Eadj(new)(r) of the new structure using a standard electromagnetic solver, and determine the diagonal part of the Green’s function G(new)(r,r) using Eqs. (12)–(19), which may involve using the standard solver for one additional field excited by a source placed at the best-improvement pixel in the 3D case. This completes a single iteration. This process iterates until the FOM no longer improves, at which point the algorithm converges to an optimal solution. In this algorithm, the storage scales linearly as the number of pixels of the system. The main computational costs for each iteration are two or three full simulations of the fields in each iteration with a standard electromagnetic solver.

In the algorithm above, we assume that at the initial step of the iteration, we have the information about the diagonal part of the Green’s function G(r,r) over the entire design region. In some cases, e.g., if the initial structure consists of only vacuum, G(r,r) is known analytically. On the other hand, there are common situations where the Green’s function for the initial structure is not known analytically. Instead of numerically determining G(r,r) for the initial structure, which can be costly with many standard EM solvers, we can adopt the following procedure where we build up the knowledge about G(r,r) as the iteration progresses. Starting from an initial structure labeled “0,” we can then randomly flip a pixel at r1 to arrive at the structure labeled “1.” Using Eqs. (17)–(19), part of the Green’s function for structure 1, G(1)(r1,r1) can be determined using the information of the fields for structures 0 and 1, without the knowledge of G(0)(r1,r1). At this point, we approximate G(1)(r,r)G(1)(r1,r1), and proceed to determine the best-improvement pixel r2 and generate the structure labeled “2.” We note that G(2)(r2,r2) can be determined rigorously in the same process as described above using Eqs. (17)–(19), whereas G(2)(r1,r1) can be determined accurately from G(1)(r1,r1) using Eq. (12). We then construct an approximation of G(2)(r,r) from G(2)(r1,r1) and G(2)(r2,r2). We repeat this process as outlined above in the iterations. As the iterations continue, the approximation to G(r,r) becomes increasingly accurate, since G(r,r) at the best-improvement pixels at all previous iterations are accurately determined. We also note that the accurate information of G(r,r) is not crucial for the beginning part of the iteration process, but becomes increasingly important as the iteration progresses. Thus, the algorithm above allows us to progressively build upon the information about G(r,r) without compromising the optimization process.

Here, for illustration purposes, we focus on a simple FOM as described in Eq. (1). The generalization to other types of FOMs has been described in Ref. [12], which involves modification of the adjoint simulation. The procedure for updating the diagonal part of the Green’s function, as described above, is unchanged.

B. Single-Block Modification

As already stated in Ref. [12], instead of updating a single pixel in each iteration, one can generalize the algorithm above to describe the case of single-block modification, where each block contains multiple pixels, and in the optimization process, a single block is modified in every step. Such single-block modification can be used to handle technological constraints, such as minimum feature size requirement, when the minimum feature size exceeds the size of one pixel (we recall that the size of the pixel needs to be sufficiently small as compared to the wavelength so that the discretization of the field is accurate). It is also useful in a 3D situation for the design of photonic circuits in a guiding layer. In this case, a block consists of multiple pixels within, and stacked in the direction perpendicular to, the guiding layer.

In this case, the derivation largely parallels the case of single-pixel modification. Equation (7) becomes

ΔFOM=2k02blockdrRe(Eadj(old)(r)·Δε(r)E(new)(r)),
where the integration occurs over the block under consideration. And Eq. (8) becomes
E(new)(r)=E(old)(r)+k02blockdr1G(old)(r,r1)Δε(r1)E(new)(r1).
Thus, to compute the change in the FOM over single-block modification, one needs to have the information about the block-diagonal part of the Green’s function G(r2,r1), where both r1 and r2 belong to the same block. With such information, to evaluate Eq. (20) involves solving Eq. (21), which is a small linear system with the dimension linearly proportional to the number of pixels (Npixel) inside a block. With this information, we then modify the block that leads to the best improvement of FOM, to arrive at the new structure; below, we refer to this block as the “best-improvement block.”

The electric field of the new structure can be determined using a standard electromagnetic solver. The remaining task to complete the iteration is to update the block-diagonal part of the Green’s function G(r2,r1). This is done by noting that updating a single block is equivalent to updating sequentially the pixels inside the block. Thus, we can update G(r2,r1) in a pixel-by-pixel fashion. Denoting the pixel being updated as rb, the block-diagonal part of the Green’s function is updated using

G(new)(r2,r1)=G(old)(r2,r1)+k02VG(old)(r2,rb)Δε(rb)G(new)(rb,r1).
G(old)(r2,rb) and G(new)(rb,r1) are determined with a similar procedure shown by Eqs. (13)–(19).

3. DESIGN EXAMPLES

We now provide examples of implementation for the algorithm as described above. We use a commercial code RSOFT-Fullwave [20], which is based on the FDTD method. The FDTD method is one of the most popular EM methods: it is versatile, able to deal with arbitrary geometries and many different kinds of materials, and it is also very memory efficient in 3D, as its memory requirement scales linearly with the system size. We implement the optimization process as a code in MATLAB. At each iteration, the optimization code communicates with RSOFT-Fullwave to get the updated electromagnetic fields, and provides the new permittivity distribution to RSOFT-Fullwave. Communications between the two codes are done by exchanging text data files [20].

A. 2D TE Example in Homogeneous Background

We start by reproducing the results of part I of the paper. We optimize in 2D a device that maximizes light intensity at a target point excited from the source at a distance of one wavelength away [Fig. 1(a)]. The pixel size is 50nm×50nm. We consider the TE polarization where the electric field is perpendicular to the 2D plane.

 

Fig. 1. 2D TE example in homogeneous background. (a) Schematic of the setup of the optimization, which aims to maximize light intensity at a target point, separated from an emitting source by a distance of one wavelength (λ=1.55μm). We consider a binary system with two permittivities of 2.25 and 6.25, respectively. The design space is a 2λ-sized square, divided in 50-nm-sized square pixels, with the source at the center. (b) Optimized structure obtained in part I [12] with the DDA approach. (c) FOM as a function of iteration number. (d) Optimized structure. (c)–(d) are obtained using the approach described in this paper. (e) Electric field of the optimized structure, as determined using FDTD simulation.

Download Full Size | PPT Slide | PDF

In the optimization process, we start from a silica homogeneous background, for which the Green’s function can be analytically obtained [15]. At each iteration, we follow the algorithms as discussed in Section 2.A, with the required fields determined from FDTD simulations.

We reproduce in Fig. 1(b) the optimized structure as obtained using the DDA approach and presented in part I. In Fig. 1(c), we show the FOM evolution as a function of iteration number, and in Figs. 1(d) and 1(e), we show the optimized structure and the field distribution of the optimized structure. The FOM of the optimized structure, as well as the number of iterations required to reach convergence, is very similar to those from the DDA approach in Ref. [12]. The final optimized structure [Fig. 1(d)] is also very similar to the structure obtained using the DDA approach [Fig. 1(b)]. However, the memory required by the optimization code and the RSOFT code [90*90 50-nm-sized square pixels including the perfectly matched layers (PML)s] together is only 257 kB, compared to 236 MB required for the DDA approach in Ref. [12], which represents a significant decrease in memory usage. In the total time spent, the optimization code accounts for only 0.38% and is thus negligible. Almost all simulation times are spent on FDTD simulations.

B. 2D TE Example in Complex Background

We then reproduce the results obtained in part I for the 2D TE device design problem with complex background [12]. The device operates as a two-wavelength demultiplexer, which separates two incident waves with wavelengths λ1=1.31μm and λ2=1.55μm from the fundamental mode of a 1-μm-wide input multimode waveguide, to two 0.3-μm-wide single-mode output guides [Fig. 2(a)]. The FOM is the average of power coupling efficiencies to the two output waveguides [12].

 

Fig. 2. 2D TE example in complex background. (a) Schematic of the setup of the optimization: light incident as the fundamental mode of a 1-μm-wide multimode input waveguide is routed to two 0.3-μm-wide single-mode output waveguides in a wavelength-dependent fashion. Light at wavelength λ1=1.31μm is routed to the left, and light at the wavelength λ2=1.55μm is routed to the right. (b) Optimized structure obtained in part I [12] with the DDA approach. (c) FOM as a function of iteration number, starting from a uniform silica background in the design region. (d) Optimized structure. (c)–(d) are obtained using the approach described in this paper. (e) Corresponding electric fields at wavelengths λ1=1.31μm (left panel) and λ2=1.55μm (right panel), as obtained using FDTD simulations of the optimized structure.

Download Full Size | PPT Slide | PDF

In this structure, the initial structure consists of the input and output waveguides, with the design region being uniform. The Green’s function of this structure is not known analytically. In contrast to the DDA-based method in part I, for which the Green’s function needed to be initialized, requiring time-costly pre-calculations or matrix inversion [12], we use the method discussed towards the end of Section 2.A, where we approximate the diagonal part of the Green’s function and gradually improve the approximation as the iteration progresses. The final structure produced in this process [Fig. 2(d)] does operate as intended [Fig. 2(e)], and is very similar to the structure produced using the DDA approach [Fig. 2(b)], with a comparable FOM of 78% [Fig. 2(c)]. The small discrepancies between the structures in Figs. 2(b) and 2(d) are expected, since the Green’s function used in the two methods are not identical. The results here indicate that the approximate Green’s function used is sufficient to reach high-performing structures. The memory required by the optimization code and the RSOFT code (110*110 50-nm-sized square pixels including PMLs) together is only 316 kB, compared to 41 MB required for the DDA approach.

C. 2D TM Example in Complex Background

Here we perform the same 2D problem as example 3.B, but for TM-polarized light where the electric field is an in-plane (Fig. 3). In this case, in each iteration, we use both forward and adjoint simulations to determine the diagonal part of the Green’s function at the best-improvement pixel, as described by Eq. (18) with r=rb. Then, the update of the diagonal part of the Green’s functions for all former best-improvement pixels positions is done with Eq. (12), using Eqs. (15) and (18). An approximate form of the diagonal part of the Green’s function is built in the process as described in Section 2.A. We confirm, here again, this process produces an optimal structure with a high FOM. The memory required by the optimization code and the RSOFT code (110*110 50-nm-sized square pixels including PMLs) together is only 392 kB, compared to 164 MB required for the DDA approach.

 

Fig. 3. 2D TM example in complex background. Same as in Fig. 2, but for TM polarization. (a) FOM as a function of iteration number, starting from a uniform silica background in the design region. (b) Optimized structure. (c) Corresponding magnetic fields at wavelengths λ1=1.31μm (left panel) and λ2=1.55μm (right panel), as obtained using FDTD simulations of the optimized structure.

Download Full Size | PPT Slide | PDF

D. 3D Example for Optimization of a Slab Structure

We then perform a similar problem as example 3.A, but in a 3D configuration [Fig. 4(a)]. We consider a silicon [ε=12.25) slab with 0.3 μm thickness, surrounded by silica (ε=2.25)]. The aim of the optimization is to maximize light intensity emitted from a y-polarized source dipole, toward a target point at half-wavelength distance away. Both the source and target are located in the middle of the slab along the y direction [Fig. 4(a)]. In the optimization process, we consider structures where various patterns are introduced into the slab. Since these patterns are generally fabricated by lithography and etching, we assume that the patterning results in structures that are uniform in y within the slab. Therefore, we use the algorithms as discussed in Section 2.B, where each block contains a line of pixels aligned along the y direction.

 

Fig. 4. 3D example for optimization of a slab structure. (a) Schematic of the setup of the optimization, which aims to maximize light intensity at a target point, separated from a y-polarized emitting source dipole by a distance of a half-wavelength (λ=1.55μm), with the source at the center. Source and target point are located in the middle of the slab along the y direction. We consider a binary system with two permittivities of 2.25 (silica) and 12.25 (silicon), respectively. The design space is a 1.55μm×0.3μm×1.55μm slab, divided in 50-nm-sized cubic pixels. (b) FOM as a function of iteration number, starting from a uniform silica background. (c) Optimized structure, using the approach described in this paper. (d) Ey electric field of the optimized structure in the middle of the slab, as determined using FDTD simulation.

Download Full Size | PPT Slide | PDF

For this optimization, we choose the initial structure consisting of uniform silica regions. As the initial background is homogenous, made of silica, one can initialize analytically the Green’s function. However, this initialization should be consistent with the way the Green’s function is updated with 3D-FDTD. In 3D-FDTD, as other 3D finite-difference methods, the electric field components are calculated at edges of the Yee cells [16,17]. Therefore, one has to analytically initialize the Green’s function correspondingly, in particular, the diagonal part, in a way different from the standard procedure defined in literature where the spatial coordinates of the Green’s function is assumed to be at the center of the cell [15,21].

Alternatively, it is also possible to initialize numerically the Green’s function for a homogeneous background. Here, we use Eq. (16) and we refer to the initial uniform silicon structure as the “old” structure, and the structure where we add a first pixel at rb of a given block as the “new” structure. Equation (16) can then be used to deduce G(old)(r2,rb) for all r2 within the same block as rb. Since the initial structure is uniform, from this we can then deduce G(old)(r2,r1) for r1 and r2 belonging to the same block, which concludes the initialization process for the block-diagonal part of the Green’s function. This process is directly compatible with 3D-FDTD simulations and results in a rigorous initialization of the block-diagonal part of the Green’s function. After the initialization process, the iteration process can then proceed according to the description in Section 2.

For this optimization problem, the FOM as a function of iteration number is shown in Fig. 4(b). As an important characteristic of this class of methods we propose here, the FOM improves for each interaction until convergence. The cross section of the resulting structure is shown in Fig. 4(c). The structure can be fabricated by etching such a cross-sectional pattern into the slab, as intended for this optimization tool. The optimal structure shows strong field concentration at the target point [Fig. 4(d)].

The reduction in memory use obtained with our approach, as compared to the DDA stand-alone method of part I [12], is here again large. The memory required by the optimization code and the RSOFT code (50*38*50 50-nm-sized cubic pixels including PMLs) together is only 9.5 MB, whereas the DDA approach would require 4.8 GB, for final results with similar performance. The additional time for the optimization algorithm accounts for 35.6% of the simulation time, due mainly to the additional Green’s simulation (33%).

E. 3D Example of a Slab Structure in Complex Background with Minimum Feature Size Constraints

We finally end with a 3D slab structure example in a complex background, with minimum feature size constraints. The design domain is 1 μm wide and 0.3 μm thick, made of silicon and silica, optimized to convert TE00 mode into TE01 mode, for identical 1-μm-wide and 0.3-μm-thick input and output Si waveguides [Fig. 5(a)]. We consider three minimum feature sizes: 50 nm, 100 nm, and 200 nm, leading, respectively, to vertical blocks of six pixels, 24 pixels, and 96 pixels. The Green’s function is not analytical for this problem in a complex background, and is thus progressively determined according to the procedure described earlier. As shown in Fig. 5, as the minimum feature size increases, the device performance naturally decreases, but the method finds the optimal device, within a given minimum feature size constraint. The memory required by the optimization code and the RSOFT code (50*38*50 50-nm-sized cubic pixels including PMLs) together is 6.6 MB, 12.9 MB, and 37.6 MB, respectively, for the 50 nm, 100 nm, and 200 nm minimum feature size cases, whereas it would be 829 MB for the DDA approach.

 

Fig. 5. 3D example of a slab structure in complex background with minimum feature size technological constraints. (a) Schematic of the setup of the optimization, which aims to maximize conversion of TE00 mode into TE01 mode, for 1-μm-wide 0.3-μm-thick input and output waveguides (λ=1.55μm). We consider a binary system with two permittivities of 2.25 (silica) and 12.25 (silicon), respectively. The design space is a 1μm×0.3μm×1μm slab, divided in 50-nm-sized cubic pixels. (b) Optimized structures, using the approach described in this paper, for minimum feature sizes of (i) 50 nm, (ii) 100 nm, and (iii) 200 nm. (c) Ex electric field of the optimized structures in the middle of the slab, as determined using FDTD simulations, with indicated conversion efficiencies.

Download Full Size | PPT Slide | PDF

4. CONCLUSION

In part I [12] of the paper, we introduced an efficient pixel-by-pixel optimization technique, based on Dyson’s equation in a Green’s function formalism. This technique naturally incorporates material and minimum feature size constraints in practical devices. In this paper, we have taken a step further, by showing that this technique can be implemented with standard commercially available electromagnetic solvers, with a memory requirement that scales linearly with the number of pixels in the design region. Our work should lead to further development of discrete optimization techniques for design and implementation of nanophotonic devices.

Funding

Air Force Office of Scientific Research (AFOSR) (FA9550-17-1-0002); FP7 People: Marie-Curie Actions (PEOPLE) (600382).

REFERENCES

1. D. J. Lockwood and L. Pavesi, “Silicon fundamentals for photonics applications,” in Silicon Photonics (Springer, 2004), pp. 1–50.

2. N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002). [CrossRef]  

3. G. Veronis, R. W. Dutton, and S. Fan, “Method for sensitivity analysis of photonic crystal devices,” Opt. Lett. 29, 2288–2290(2004). [CrossRef]  

4. C. M. Lalau-Keraly, S. Bhargava, O. D. Miller, and E. Yablonovitch, “Adjoint shape optimization applied to electromagnetic design,” Opt. Express 21, 21693–21701 (2013). [CrossRef]  

5. O. D. Miller, “Photonic design: from fundamental solar cell physics to computational inverse design,” arXiv:1308.0212 (2013).

6. J. Lu and J. Vučković, “Objective-first design of high-efficiency, small-footprint couplers between arbitrary nanophotonic waveguide modes,” Opt. Express 20, 7221–7236 (2012). [CrossRef]  

7. 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, 374–377 (2015). [CrossRef]  

8. V. Liu, Y. Jiao, D. A. Miller, and S. Fan, “Design methodology for compact photonic-crystal-based wavelength division multiplexers,” Opt. Lett. 36, 591–593 (2011). [CrossRef]  

9. A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017). [CrossRef]  

10. M. A. Seldowitz, J. P. Allebach, and D. W. Sweeney, “Synthesis of digital holograms by direct binary search,” Appl. Opt. 26, 2788–2798 (1987). [CrossRef]  

11. B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015). [CrossRef]  

12. S. Boutami and S. Fan, “Efficient pixel-by-pixel optimization of photonic devices utilizing the Dyson’s equation in a Green function formalism: Part I: implementation with the method of discrete dipole approximation,” J. Opt. Soc. Am. B 36, 2378–2386 (2019). [CrossRef]  

13. O. J. Martin, A. Dereux, and C. Girard, “Iterative scheme for computing exactly the total field propagating in dielectric structures of arbitrary shape,” J. Opt. Soc. Am. A 11, 1073–1080 (1994). [CrossRef]  

14. B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. 333, 848–872 (1988). [CrossRef]  

15. O. J. Martin and N. B. Piller, “Electromagnetic scattering in polarizable backgrounds,” Phys. Rev. E 58, 3909 (1998). [CrossRef]  

16. K. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966). [CrossRef]  

17. A. Taflove, “Application of the finite-difference time-domain method to sinusoidal steady-state electromagnetic-penetration problems,” IEEE Trans. Electromagn. Compat. EMC-22, 191–202 (1980). [CrossRef]  

18. J. M. Jin, The Finite Element Method in Electromagnetics (Wiley, 2015).

19. W. Shin and S. Fan, “Choice of the perfectly matched layer boundary condition for frequency-domain Maxwell’s equations solvers,” J. Comput. Phys. 231, 3406–3431 (2012). [CrossRef]  

20. Fullwave® and Synopsis®, “Rsoft Design,”.

21. A. D. Yaghjian, “Electric dyadic Green’s functions in the source region,” Proc. IEEE 68, 248–263 (1980). [CrossRef]  

References

  • View by:
  • |
  • |
  • |

  1. D. J. Lockwood and L. Pavesi, “Silicon fundamentals for photonics applications,” in Silicon Photonics (Springer, 2004), pp. 1–50.
  2. N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
    [Crossref]
  3. G. Veronis, R. W. Dutton, and S. Fan, “Method for sensitivity analysis of photonic crystal devices,” Opt. Lett. 29, 2288–2290(2004).
    [Crossref]
  4. C. M. Lalau-Keraly, S. Bhargava, O. D. Miller, and E. Yablonovitch, “Adjoint shape optimization applied to electromagnetic design,” Opt. Express 21, 21693–21701 (2013).
    [Crossref]
  5. O. D. Miller, “Photonic design: from fundamental solar cell physics to computational inverse design,” arXiv:1308.0212 (2013).
  6. J. Lu and J. Vučković, “Objective-first design of high-efficiency, small-footprint couplers between arbitrary nanophotonic waveguide modes,” Opt. Express 20, 7221–7236 (2012).
    [Crossref]
  7. 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, 374–377 (2015).
    [Crossref]
  8. V. Liu, Y. Jiao, D. A. Miller, and S. Fan, “Design methodology for compact photonic-crystal-based wavelength division multiplexers,” Opt. Lett. 36, 591–593 (2011).
    [Crossref]
  9. A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
    [Crossref]
  10. M. A. Seldowitz, J. P. Allebach, and D. W. Sweeney, “Synthesis of digital holograms by direct binary search,” Appl. Opt. 26, 2788–2798 (1987).
    [Crossref]
  11. B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
    [Crossref]
  12. S. Boutami and S. Fan, “Efficient pixel-by-pixel optimization of photonic devices utilizing the Dyson’s equation in a Green function formalism: Part I: implementation with the method of discrete dipole approximation,” J. Opt. Soc. Am. B 36, 2378–2386 (2019).
    [Crossref]
  13. O. J. Martin, A. Dereux, and C. Girard, “Iterative scheme for computing exactly the total field propagating in dielectric structures of arbitrary shape,” J. Opt. Soc. Am. A 11, 1073–1080 (1994).
    [Crossref]
  14. B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. 333, 848–872 (1988).
    [Crossref]
  15. O. J. Martin and N. B. Piller, “Electromagnetic scattering in polarizable backgrounds,” Phys. Rev. E 58, 3909 (1998).
    [Crossref]
  16. K. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966).
    [Crossref]
  17. A. Taflove, “Application of the finite-difference time-domain method to sinusoidal steady-state electromagnetic-penetration problems,” IEEE Trans. Electromagn. Compat. EMC-22, 191–202 (1980).
    [Crossref]
  18. J. M. Jin, The Finite Element Method in Electromagnetics (Wiley, 2015).
  19. W. Shin and S. Fan, “Choice of the perfectly matched layer boundary condition for frequency-domain Maxwell’s equations solvers,” J. Comput. Phys. 231, 3406–3431 (2012).
    [Crossref]
  20. Fullwave® and Synopsis®, “Rsoft Design,”.
  21. A. D. Yaghjian, “Electric dyadic Green’s functions in the source region,” Proc. IEEE 68, 248–263 (1980).
    [Crossref]

2019 (1)

2017 (1)

A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
[Crossref]

2015 (2)

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, 374–377 (2015).
[Crossref]

B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
[Crossref]

2013 (1)

2012 (2)

J. Lu and J. Vučković, “Objective-first design of high-efficiency, small-footprint couplers between arbitrary nanophotonic waveguide modes,” Opt. Express 20, 7221–7236 (2012).
[Crossref]

W. Shin and S. Fan, “Choice of the perfectly matched layer boundary condition for frequency-domain Maxwell’s equations solvers,” J. Comput. Phys. 231, 3406–3431 (2012).
[Crossref]

2011 (1)

2004 (1)

2002 (1)

N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
[Crossref]

1998 (1)

O. J. Martin and N. B. Piller, “Electromagnetic scattering in polarizable backgrounds,” Phys. Rev. E 58, 3909 (1998).
[Crossref]

1994 (1)

1988 (1)

B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. 333, 848–872 (1988).
[Crossref]

1987 (1)

1980 (2)

A. D. Yaghjian, “Electric dyadic Green’s functions in the source region,” Proc. IEEE 68, 248–263 (1980).
[Crossref]

A. Taflove, “Application of the finite-difference time-domain method to sinusoidal steady-state electromagnetic-penetration problems,” IEEE Trans. Electromagn. Compat. EMC-22, 191–202 (1980).
[Crossref]

1966 (1)

K. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966).
[Crossref]

Allebach, J. P.

Babinec, T. M.

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, 374–377 (2015).
[Crossref]

Bakr, M. H.

N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
[Crossref]

Bandler, J. W.

N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
[Crossref]

Bhargava, S.

Boutami, S.

Dereux, A.

Draine, B. T.

B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. 333, 848–872 (1988).
[Crossref]

Dutton, R. W.

Fan, S.

Georgieva, N. K.

N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
[Crossref]

Girard, C.

Glavic, S.

N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
[Crossref]

Jiao, Y.

Jin, J. M.

J. M. Jin, The Finite Element Method in Electromagnetics (Wiley, 2015).

Lagoudakis, K. G.

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, 374–377 (2015).
[Crossref]

Lalau-Keraly, C. M.

Liu, V.

Lockwood, D. J.

D. J. Lockwood and L. Pavesi, “Silicon fundamentals for photonics applications,” in Silicon Photonics (Springer, 2004), pp. 1–50.

Lu, J.

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, 374–377 (2015).
[Crossref]

J. Lu and J. Vučković, “Objective-first design of high-efficiency, small-footprint couplers between arbitrary nanophotonic waveguide modes,” Opt. Express 20, 7221–7236 (2012).
[Crossref]

Martin, O. J.

Menon, R.

B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
[Crossref]

Miller, D. A.

Miller, O. D.

C. M. Lalau-Keraly, S. Bhargava, O. D. Miller, and E. Yablonovitch, “Adjoint shape optimization applied to electromagnetic design,” Opt. Express 21, 21693–21701 (2013).
[Crossref]

O. D. Miller, “Photonic design: from fundamental solar cell physics to computational inverse design,” arXiv:1308.0212 (2013).

Pavesi, L.

D. J. Lockwood and L. Pavesi, “Silicon fundamentals for photonics applications,” in Silicon Photonics (Springer, 2004), pp. 1–50.

Petykiewicz, J.

A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
[Crossref]

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, 374–377 (2015).
[Crossref]

Piggott, A. Y.

A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
[Crossref]

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, 374–377 (2015).
[Crossref]

Piller, N. B.

O. J. Martin and N. B. Piller, “Electromagnetic scattering in polarizable backgrounds,” Phys. Rev. E 58, 3909 (1998).
[Crossref]

Polson, R.

B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
[Crossref]

Seldowitz, M. A.

Shen, B.

B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
[Crossref]

Shin, W.

W. Shin and S. Fan, “Choice of the perfectly matched layer boundary condition for frequency-domain Maxwell’s equations solvers,” J. Comput. Phys. 231, 3406–3431 (2012).
[Crossref]

Su, L.

A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
[Crossref]

Sweeney, D. W.

Taflove, A.

A. Taflove, “Application of the finite-difference time-domain method to sinusoidal steady-state electromagnetic-penetration problems,” IEEE Trans. Electromagn. Compat. EMC-22, 191–202 (1980).
[Crossref]

Veronis, G.

Vuckovic, J.

A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
[Crossref]

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, 374–377 (2015).
[Crossref]

J. Lu and J. Vučković, “Objective-first design of high-efficiency, small-footprint couplers between arbitrary nanophotonic waveguide modes,” Opt. Express 20, 7221–7236 (2012).
[Crossref]

Wang, P.

B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
[Crossref]

Yablonovitch, E.

Yaghjian, A. D.

A. D. Yaghjian, “Electric dyadic Green’s functions in the source region,” Proc. IEEE 68, 248–263 (1980).
[Crossref]

Yee, K.

K. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966).
[Crossref]

Appl. Opt. (1)

Astrophys. J. (1)

B. T. Draine, “The discrete-dipole approximation and its application to interstellar graphite grains,” Astrophys. J. 333, 848–872 (1988).
[Crossref]

IEEE Trans. Antennas Propag. (1)

K. Yee, “Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media,” IEEE Trans. Antennas Propag. 14, 302–307 (1966).
[Crossref]

IEEE Trans. Electromagn. Compat. (1)

A. Taflove, “Application of the finite-difference time-domain method to sinusoidal steady-state electromagnetic-penetration problems,” IEEE Trans. Electromagn. Compat. EMC-22, 191–202 (1980).
[Crossref]

IEEE Trans. Microwave Theory Tech. (1)

N. K. Georgieva, S. Glavic, M. H. Bakr, and J. W. Bandler, “Feasible adjoint sensitivity technique for EM design optimization,” IEEE Trans. Microwave Theory Tech. 50, 2751–2758 (2002).
[Crossref]

J. Comput. Phys. (1)

W. Shin and S. Fan, “Choice of the perfectly matched layer boundary condition for frequency-domain Maxwell’s equations solvers,” J. Comput. Phys. 231, 3406–3431 (2012).
[Crossref]

J. Opt. Soc. Am. A (1)

J. Opt. Soc. Am. B (1)

Nat. Photonics (2)

B. Shen, P. Wang, R. Polson, and R. Menon, “An integrated-nanophotonics polarization beamsplitter with 2.4× 2.4 μm2 footprint,” Nat. Photonics 9, 378–382 (2015).
[Crossref]

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, 374–377 (2015).
[Crossref]

Opt. Express (2)

Opt. Lett. (2)

Phys. Rev. E (1)

O. J. Martin and N. B. Piller, “Electromagnetic scattering in polarizable backgrounds,” Phys. Rev. E 58, 3909 (1998).
[Crossref]

Proc. IEEE (1)

A. D. Yaghjian, “Electric dyadic Green’s functions in the source region,” Proc. IEEE 68, 248–263 (1980).
[Crossref]

Sci. Rep. (1)

A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vučković, “Fabrication-constrained nanophotonic inverse design,” Sci. Rep. 7, 1786 (2017).
[Crossref]

Other (4)

D. J. Lockwood and L. Pavesi, “Silicon fundamentals for photonics applications,” in Silicon Photonics (Springer, 2004), pp. 1–50.

O. D. Miller, “Photonic design: from fundamental solar cell physics to computational inverse design,” arXiv:1308.0212 (2013).

J. M. Jin, The Finite Element Method in Electromagnetics (Wiley, 2015).

Fullwave® and Synopsis®, “Rsoft Design,”.

Cited By

OSA participates in Crossref's Cited-By Linking service. Citing articles from OSA journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (5)

Fig. 1.
Fig. 1. 2D TE example in homogeneous background. (a) Schematic of the setup of the optimization, which aims to maximize light intensity at a target point, separated from an emitting source by a distance of one wavelength ( λ = 1.55 μm ). We consider a binary system with two permittivities of 2.25 and 6.25, respectively. The design space is a 2 λ -sized square, divided in 50-nm-sized square pixels, with the source at the center. (b) Optimized structure obtained in part I [12] with the DDA approach. (c) FOM as a function of iteration number. (d) Optimized structure. (c)–(d) are obtained using the approach described in this paper. (e) Electric field of the optimized structure, as determined using FDTD simulation.
Fig. 2.
Fig. 2. 2D TE example in complex background. (a) Schematic of the setup of the optimization: light incident as the fundamental mode of a 1-μm-wide multimode input waveguide is routed to two 0.3-μm-wide single-mode output waveguides in a wavelength-dependent fashion. Light at wavelength λ 1 = 1.31 μm is routed to the left, and light at the wavelength λ 2 = 1.55 μm is routed to the right. (b) Optimized structure obtained in part I [12] with the DDA approach. (c) FOM as a function of iteration number, starting from a uniform silica background in the design region. (d) Optimized structure. (c)–(d) are obtained using the approach described in this paper. (e) Corresponding electric fields at wavelengths λ 1 = 1.31 μm (left panel) and λ 2 = 1.55 μm (right panel), as obtained using FDTD simulations of the optimized structure.
Fig. 3.
Fig. 3. 2D TM example in complex background. Same as in Fig. 2, but for TM polarization. (a) FOM as a function of iteration number, starting from a uniform silica background in the design region. (b) Optimized structure. (c) Corresponding magnetic fields at wavelengths λ 1 = 1.31 μm (left panel) and λ 2 = 1.55 μm (right panel), as obtained using FDTD simulations of the optimized structure.
Fig. 4.
Fig. 4. 3D example for optimization of a slab structure. (a) Schematic of the setup of the optimization, which aims to maximize light intensity at a target point, separated from a y -polarized emitting source dipole by a distance of a half-wavelength ( λ = 1.55 μm ), with the source at the center. Source and target point are located in the middle of the slab along the y direction. We consider a binary system with two permittivities of 2.25 (silica) and 12.25 (silicon), respectively. The design space is a 1.55 μm × 0.3 μm × 1.55 μm slab, divided in 50-nm-sized cubic pixels. (b) FOM as a function of iteration number, starting from a uniform silica background. (c) Optimized structure, using the approach described in this paper. (d) Ey electric field of the optimized structure in the middle of the slab, as determined using FDTD simulation.
Fig. 5.
Fig. 5. 3D example of a slab structure in complex background with minimum feature size technological constraints. (a) Schematic of the setup of the optimization, which aims to maximize conversion of TE 00 mode into TE 01 mode, for 1-μm-wide 0.3-μm-thick input and output waveguides ( λ = 1.55 μm ). We consider a binary system with two permittivities of 2.25 (silica) and 12.25 (silicon), respectively. The design space is a 1 μm × 0.3 μm × 1 μm slab, divided in 50-nm-sized cubic pixels. (b) Optimized structures, using the approach described in this paper, for minimum feature sizes of (i) 50 nm, (ii) 100 nm, and (iii) 200 nm. (c) Ex electric field of the optimized structures in the middle of the slab, as determined using FDTD simulations, with indicated conversion efficiencies.

Equations (22)

Equations on this page are rendered with MathJax. Learn more.

FOM = | E ( r 0 ) | 2 .
Δ FOM = | E ( new ) ( r 0 ) | 2 | E ( old ) ( r 0 ) | 2 ,
E ( new ) ( r ) = E ( old ) ( r ) + V k 0 2 G ( old ) ( r , r ) Δ ε ( r ) E ( new ) ( r ) d r ,
Δ FOM = | E ( old ) ( r 0 ) + k 0 2 V G ( old ) ( r 0 , r ) Δ ε ( r ) E ( new ) ( r ) | 2 | E ( old ) ( r 0 ) | 2 .
G ( old ) ( r 0 , r ) = G ( old ) ( r , r 0 ) T ,
E adj ( old ) ( r ) = G ( old ) ( r , r 0 ) E ( old ) ( r 0 ) * ,
Δ FOM = 2 k 0 2 V Re ( E adj ( old ) ( r ) · Δ ε ( r ) E ( new ) ( r ) ) .
E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( old ) ( r , r ) Δ ε ( r ) E ( new ) ( r ) .
E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E ( new ) ( r b ) .
G ( new ) ( r , r ) = G ( old ) ( r , r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) G ( new ) ( r b , r ) .
G ( new ) ( r b , r ) = G ( old ) ( r b , r ) + k 0 2 V G ( old ) ( r b , r b ) Δ ε ( r ) G ( new ) ( r b , r ) .
G ( new ) ( r , r ) = G ( old ) ( r , r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) G ( new ) ( r b , r ) .
G ( new ) ( r b , r ) = G ( new ) ( r , r b ) T .
E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E ( new ) ( r b ) .
{ E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E ( new ) ( r b ) , E adj ( new ) ( r ) = E adj ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E adj ( new ) ( r b ) .
{ E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E ( new ) ( r b ) , E adj ( new ) ( r ) = E adj ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E adj ( new ) ( r b ) , E Green ( new ) ( r ) = E Green ( old ) ( r ) + k 0 2 V G ( old ) ( r , r b ) Δ ε ( r b ) E Green ( new ) ( r b ) ,
E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( new ) ( r , r b ) Δ ε ( r b ) E ( old ) ( r b ) .
{ E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( new ) ( r , r b ) Δ ε ( r b ) E ( old ) ( r b ) , E adj ( new ) ( r ) = E adj ( old ) ( r ) + k 0 2 V G ( new ) ( r , r b ) Δ ε ( r b ) E adj ( old ) ( r b ) ,
{ E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 V G ( new ) ( r , r b ) Δ ε ( r b ) E ( old ) ( r b ) , E adj ( new ) ( r ) = E adj ( old ) ( r ) + k 0 2 V G ( new ) ( r , r b ) Δ ε ( r b ) E adj ( old ) ( r b ) , E Green ( new ) ( r ) = E Green ( old ) ( r ) + k 0 2 V G ( new ) ( r , r b ) Δ ε ( r b ) E Green ( old ) ( r b ) .
Δ FOM = 2 k 0 2 block d r Re ( E adj ( old ) ( r ) · Δ ε ( r ) E ( new ) ( r ) ) ,
E ( new ) ( r ) = E ( old ) ( r ) + k 0 2 block d r 1 G ( old ) ( r , r 1 ) Δ ε ( r 1 ) E ( new ) ( r 1 ) .
G ( new ) ( r 2 , r 1 ) = G ( old ) ( r 2 , r 1 ) + k 0 2 V G ( old ) ( r 2 , r b ) Δ ε ( r b ) G ( new ) ( r b , r 1 ) .

Metrics