Inversely design a silicon-on-insulator grating coupler using the built-in parametrized geometry optimization module (PO) for FDTD. Compared to other optimization methods such as particle swarm optimization (PSO), this optimization algorithm enables obtaining the best solution in just a few iterations. In this example, the etch depth will be fixed at 80% while pitch and duty cycle will be optimized for 20 gratings. The optimal design will then be exported into a GDS file for further simulation and/or fabrication.

Related: Grating Coupler

Minimum product version: 2019b r2

Updating the model with your parameters

Appendix: Parametric Optimization (Adjoint method)

Appendix: S-parameter-extraction

Understand the simulation workflow and key results

Lumerical's inverse design capability provides unparalleled optimization performance by combining the power of gradient based optimization routine with efficiencies found in fundamental properties of Maxwell equations.

This example will demonstrate how to use the inverse design method to generate a TE silicon on insulator (SOI) grating coupler design with maximized coupling efficiency. Moreover, we'll demonstrate how to modify the example with your parameters so you can reuse this approach for your own design.

This example draws extensively from the LumOpt framework:

•Git

The goal for this initial step is to define the basic parameters of the FDTD simulation that will be used for the adjoint optimization. Here, a 2D base simulation file for the grating coupler that includes monitors for field gradient and figure of merit (FOM) calculations, source and mesh are defined to run the optimization routine.

The second step is to define a polygon geometry that will represent the body of the grating coupler (with 80% etch depth) to be optimized. The number of gratings is fixed at 20 and pitch and duty cycles are optimized for maximum transmission/coupling. This is defined as a Python function that accepts a set of optimization parameters and generates corresponding polygon geometry appended to the base FDTD simulation geometry defined in step 1. Initial parameter values and their limits will also get defined in this step.

Run the adjoint optimization Python script with the base simulation file and parametrized polygon object as an input. This will be a 2D simulation for fast delivery of a near optimal design which will be used for the 3D simulations and s-parameter extraction in the next step. To find the best set of parameters that define the optimal grating coupler shape, the optimization routine will use the specified field monitors in the base simulation to calculate:

•the field gradients: the gradients of the field due to perturbation of permittivity resulting from slight changes in geometry of the optimizable region

•the gradient of the figure of merit (FOM): the gradient of the mode overlaps to the fundamental TE mode of the waveguide as a result of the shape changes

After optimization, the optimized grating coupler component shape with curved gratings will be exported into GDS II format and into 3D_grating_coupler.fsp simulation file which can be used for further simulation and/or fabrication (mask design). We used circular gratings that cover the source area to make sure that all the fiber light is collected. We also used a short linear taper to connect to 500nm wide waveguide. The geometry of the gratings and taper can be further optimized by the designer.

The grating coupler geometry will be automatically updated after finishing the 2D optimization step. Run the 3D_grating_coupler.fsp simulation file to compare the results with 2D and extract the s-parameter for circuit simulations.

Instructions for running the model and discussion of key results

Step 1: Define base simulation parameters

This step is not necessary if you want to run the example as is, however, it is required for modification of base simulation file based on your needs.

1.Open the base simulation file (grating_base.fsp) in FDTD

2.Modify the desired base simulation parameters

3.Save the file.

The resulting base simulation which includes the FDTD simulation region, source, monitors, mesh override for optimizable geometry, and angled single mode fiber is shown below.

This step is not necessary if you want to run the example as is, however it is required for modification of the optimizable geometry based on your needs.

1.Open the adjoint optimization python script file (grating_opt.py) in FDTD script editor

2.Define the initial values of the input parameters of the optimization (the nodes of the optimizable polygon) such as etch depth based on fabrication goals

3.Define the bounds of the optimization parameters within the array bounds

4.Make any other necessary modifications to the optimization parameters. See “Important model settings” for more information about these parameters

5.Save the script

1. Before running the optimization, make sure the current working directory for FDTD is set to the location where the example files are located

2. Run the optimization script file from the script editor in FDTD

The geometry of the grating coupler at the start of the optimization with default parameters is shown in the figure below:

A command prompt (in windows operating system) or terminal (in macOS or Linux) window will open showing the progress of the optimization. A screenshot of the command prompt window at the end of the optimization is shown below which contains information about the optimization progress, FOM and parameter values:

During the optimization, a window will report the optimization progress providing various plots for evolution of the optimization geometry, figure of merit, parameters and other related information. The figure below shows the optimization report window at the end of the optimization process.

From the plot, it can be seen that the optimization started with an initial design providing an FOM of ~0.20 and ended up with a geometry with a FOM of about 0.45, resulting in more than 100% improvement through the optimization process in just a few iterations which highlights the advantage of this method compared to other optimization algorithms such as particle swarm which generally need more iterations to achieve an optimum solution. Also, note that the FOM will have little improvement after about 15 iterations which suggests higher number of iterations probably won't be necessary. While the number of iterations is set to 200, simulation automatically stops after 50 iteration (2 hours of simulations) as it reaches to optimum value.The field gradients (due to the changes in geometry) and forward direction field distributions within the device’s structure can also be seen in the plot. In addition, the changes in the optimizable parameter values become less noticeable after a few iterations as the FOM gets closer to the optimum value. The changes in FOM gradient (based on changes in the parameter values) also decrease on average as seen in gradient evolution plot as a function of the iterations.

At the end of the optimization, the final parameter values will be saved in a text file (2D_parameters.txt) in the same folder as the example files which will be used in the 3D optimization for the next step. In addition, copies of the simulation files for each iteration, a text file containing the history of parameters and FOM values for each iteration and still images of the optimization report window for each iteration will be saved in a folder (with a name starting with “opts”) in the same folder as example files for future reference. Each optimization run creates a new folder with the same name but a numbered suffix to avoid overwriting previous runs data.

Here is the final geometry optimized in 2D:

The final geometry will be created in 3D with curved gratings. The optimization Python script file will also export the optimal design to a GDS file at the end of optimization. The GDS file will be saved in the same folder as the example files. The example files come with encrypted script file with .lsfx extension which are needed for the GDS export operation. It should be saved in the same location as the Python script file.

1.Go to the “Optimizations and Sweeps” window and run the sweep_source.

2.Once the sweep is finished, visualize the transmission to obtain the optimal fiber position.

3.Update the fiber and port1 x-position based on the results from step 2.

4.Once again, go to the “Optimizations and Sweeps” window, select the sweep item named “S-parameter sweep” and run the sweep. The sweep will launch two simulations that should complete in about an hour.

5.Once the sweep has completed, the scattering parameters of the device become available. To see them, left click on the S-parameters sweep item and choose “Visualize” followed by “S parameters” on the context menu. Select scalar operation “Abs^2” to see the power s-parameters

6.To export the S-parameter results, left click on the S-parameters sweep item and choose “Export to INTERCONNECT” on the context menu. In the subsequent “Export to INTERCONNECT” window, enter a name and location for the data file and click on “Save”.

The obtained s-parameter spectra indicate a peak power coupling efficiency of about 40% is illustrated below:

For further information about s-parameter extraction, see the appendix section of this example.

Description of important objects and settings used in this model

Polarization: The chosen refractive index values are representative of an SOI chip manufacturing process. Because of the high refractive index contrast between silicon and silicon oxide, there is a large difference between effective indices of the two fundamental modes – TE and TM – of the integrated waveguide. For this reason, SOI grating couplers are strongly polarization selective. The presented design excites a TE mode since this is the most common choice, however, it is also possible to set up an optimization targeting the TM mode.

Tilting angle: The coupling efficiency depends heavily on how the fiber meets the top silicon oxide cladding layer. In this example, the end of the fiber is assumed to be polished at a small angle so that the fiber tilts as it is mounted on the top cladding. Such tilting prevents reflections into the fiber.

Etch depth: The coupling efficiency is highly sensitive to the grating’s pitch, duty cycle and etch depth. For simplicity, a fixed etch depth is employed here, however, it can also be varied if the available fabrication process provides that degree of freedom.

Substrate: If a silicon substrate is present in the fabricated device, it should be included in the simulations. The substrate will have a noticeable impact on how light is coupled and cannot be omitted as is often done in other device designs.

Structure groups: Structure groups are used here to update the geometric primitives of the fiber. Any changes to the fiber parameters must be applied using the interface of the structure group. The advantage of using structure groups is that they can apply a single parameter change to multiple primitives, however, they can also override manual changes on individual geometries. The fiber structure group is setup to mimic a fiber polished at an angle. This is done by applying a slight rotation to the objects in the group and the use of “mesh order” setting for the objects so that the fiber only extends down to some point above the grating coupler not all the way through the fiber as might be depicted from the layout view. This can be verified by visualizing the index profile of the structure as shown below:

Optimization figure of merit (FOM): Since the purpose of the design is to have the best possible coupling at a desired wavelength, the optimization figure of merit is chosen to be the average transmission through the grating coupler from 1530 – 1570 nm, and the optimization algorithm will try to maximize this value.

Definition of optimizable geometry: The optimizable geometry is defined as a polygon which has fixed values for the y coordinates of its boundary points and the points’ x coordinates can be modified to obtain the optimal geometry.

Optimization fields monitor: The optimization fields DFT monitor (opt_field) is used to collect the field data within the optimizable geometry which is then used to calculate field gradients used in optimization algorithm. As such, the location of this monitor is of great importance and should cover the entire optimizable geometry.

Figure of merit fields monitor: Since this monitor is used to calculate the figure of merit required for optimization (mode overlap to the fundamental TE mode of the output waveguide), it should be located within the output waveguide of the grating coupler with appropriate dimensions.

Height of the geometry objects: The height (depth) of the geometry objects is chosen based on fabrication/foundry process. This is especially important for the 3D simulation to ensure correct geometry setup. This is done by adjusting the “wg_height” and “etch depth” in the script as the input argument when calling the grate function to form the optimizable geometry. These value should also be passed to the GDS export portion of the script as shown below (highlighted in red):

gds_export_script = "\

gds_filename = 'Grating_Coupler.gds';\

top_cell = 'model';\

layer_def = [1,0.132e-6,0.176e-6 ; 2,0.11e-6, 0.22e-6];\

n_circle = 64;\

n_ring = 64;\

n_custom = 64;\

n_wg = 64;\

round_to_nm = 1;\

grid = 1e-9;\

Lumerical_GDS_auto_export;\

"

...

fdtd.eval(gds_export_script)

fdtd.close()

Optimization parameters bounds and initial values: The range over which the optimization algorithm is allowed to vary the parameters and also the initial values for those parameters as a starting point are defined by initial_params and bounds arrays as input arguments when calling the function FunctionDefinedPolygonI.

Optimization spectral range: The spectral range of the optimization can be specified as the array “wavelengths” as well as the number of frequency points considered in optimization. Note that choosing a large number of frequency points will make the optimization much slower and can cause issues when transferring large amount of data between FDTD and python environment. It is recommended to keep this number as low as possible especially for 3D simulations. For a single frequency optimization, simply choose the same start and stop values with points set to 1.

Maximum number of iterations: While the algorithm has the capability to stop the optimization once the gradient of the figure of merit falls below a certain threshold, the “max_iter” variable used for defining the optimization algorithm can be used to limit the number of iterations the algorithm is able to perform.

Target FOM: The target values for FOM at different wavelengths can be specified with target_T_fwd input argument when defining the FOM. This should be an array with the same length as wavelength data and contain target FOM values (forward transmission) for each wavelength which can be less than 1. The default is 1 for all wavelengths. This gives users the option to aim for a desired value of transmission within the design spectral range.

Instructions for updating the model based on your device parameters

Geometry: If you need to have your own geometry defined for the grating coupler, including SOI device layer thickness (grating height) and etch depth, corresponding changes should be made in base simulation file and/or optimization setup Python script. This might need changes in objects span and location in the base script, the polygon’s “etch_depth” parameter in python script which its value is also passed to the GDS export portion of the script. The source, simulation region, mesh override, and field gradient and FOM monitors dimensions should also be adjusted accordingly to make sure they cover the entire structure properly.

Materials: The permittivity of the materials included in the simulation (the material making the optimizable geometry (core) and the material surrounding it (cladding)) for the desired simulation wavelength should be passed to the optimizer when defining the geometry by calling the function FunctionDefinedPolygon (here they are presented as refractive index (n) squared). The refractive index (not permittivity) of the waveguides (which is the same as the optimizable geometry) should also be defined in the base simulation setup script.

Grating geometry: We used circular shape gratings for 3D simulations and GDS export. One can use more sophisticated grating shapes similar to focusing grating couplers.

Information and tips for users that want to further customize the model

S-parameter: The S-parameter extraction step generates scattering parameter results only for the TE polarization. To add the TM polarization, simply select both the TE and TM modes on the integrated waveguide port and re-run the S-parameter sweep. This change will add one additional simulation and one additional result to the S-parameter sweep.

Taper optimization: The 3-D coupler model uses a linear taper section to connect to the integrated waveguide at the start of the grating. The coupling efficiency can also be improved by optimizing the taper shape (see SOI taper design).

High efficiency grating couplers: Couplers with an efficiency higher than 90% over a large bandwidth have been designed with FDTD using more sophisticated gratings and a mixed 2-D/3-D optimization strategy (see references).

Additional documentation, examples and training material

•Inverse design of grating coupler: Neil V. Sapra et.al, “Inverse design and demonstration of broadband grating couplers” IEEE Journal of Selected Topics in Quantum Electronics ( Volume: 25 , Issue: 3 , May-June 2019 )

•Basic grating coupler design: D. Taillaert, F. Van Laere, M. Ayre, W. Bogaerts, D. Van Thourhout, P. Bienstman and R. Baets, “Grating Couplers for Coupling between Optical Fibers and Nanophotonic Waveguides,” Japanese Journal of Applied Physics, vol. 45, no. 8a, pp. 6071-6077, 2006.

•Advanced optimization: R. Marchetti, C. Lacava, A. Khokhar, X. Chen, I. Cristiani, D. J. Richardson, G. T. Reed, P. Petropoulos and P. Minzioni, “High-efficiency grating-couplers: demonstration of a new design strategy,” Scientific Reports, Article number: 16670, 2017. ; T. Watanabe, M. Ayata, U. Koch, Y. Fedoryshyn and J. Leuthold, “Perpendicular Grating Coupler Based on a Blazed Antiback-Reflection Structure,” Journal of Lightwave Technology, vol. 35, no. 21, pp. 4663- 4669, 2017.

•MATLAB driven grating coupler optimization

•Related Lumerical University courses

Additional background information and theory

In a typical device optimization, the goal is to minimize or maximize a figure or merit \( F(p) \) that depends on a set of \( N \) design parameters \( p=(p_1,...,p_n) \). The figure of merit can be any quantity used to benchmark the device’s performance and the design parameters can be any quantity that alters the device’s response. To evaluate \( F(p) \) for a given set of values of \( p \), a physics based modeling tool - such as FDTD - is used to map values of \( p \) to values of \( F \). The model usually contains a description of the device’s geometry, materials and all relevant environmental conditions. In a parameterized shape optimization, the goal is to allow some parts of the device’s geometry to vary to minimize or maximize the figure of merit. For this purpose, a shape function \( f(p) \) is typically defined to specify part of the device’s geometry; an example is shown in the following figure:

The shape function is typically user defined and must be fed to the modelling tool to construct the device’s geometry. Once the figure of merit \( F(f(p)) \) can be evaluated for a given shape function, a minimizer can make multiple evaluations of the figure of merit to find the optimal shape. There is a very large body of minimization tools and techniques that can be used to try to find the optimal shape parameters. When selecting a minimization technique, it is important to consider that running a device simulation to evaluate the figure of merit can be a slow operation. So, the best suited minimization techniques are those that limit the number of evaluations of the figure of merit. With a good starting point, gradient based methods can quickly find the optimal parameters with a relatively small number of figure of merit evaluations provided that it is also possible to evaluate the gradient of the figure of merit:

$$∇_pF=(\frac{\partial F}{\partial p_1},⋯,\frac{\partial F}{\partial p_n})$$

To evaluate the gradient without having to perform one simulation for each partial derivative, a technique called the adjoint method can be employed. This technique exploits the properties of linear partial differential equations to evaluate the gradient using only two simulations independently of the number of optimization parameters. The Python based shape optimization module in FDTD (lumopt) allows users to define a shape function and run the corresponding device simulation to evaluate the figure of merit. The module also implements the adjoint method to provide the gradient. In this way, any of the gradient based minimizers in SciPy’s optimize module can be used to perform a fully automated optimization.

Additional background information and theory

The S-parameter matrix formalism is a common approach to build compact models of photonic devices to be used in circuit-level simulations. Assuming the response of the device to optical signals is linear, it can be modeled by a network (black box) with multiple network ports, where each of them receives an incoming signal and scatters or reflects an outgoing signal.

We assume the device has N physical ports (which typically are input/output waveguide channels) and each of these supports a certain number of modes of interest. Each element of the S-parameter matrix is the ratio between the complex envelopes of the optical modes in two physical ports. For example, to model the behavior of the fundamental TE and TM modes at the input/output waveguide channels of a 1x2 MMI we need a 6x6 S-parameter matrix, given that there are two modes at each of the three physical ports. Each column of the S-parameter matrix corresponds to a fixed input physical port and mode and each row corresponds to a fixed output physical port and mode.

The approach to extract the S-parameter matrix of a photonic device depends on the solver being used. In both FDTD and EME solvers, port objects can be setup at each of the physical ports. The main differences between the solvers for the S-parameter extraction are:

•FDTD injects one mode in one input port object per simulation, while EME can solve for all the elements of the S-matrix in one simulation;

•FDTD is broadband so many frequency points can be calculated in one simulation, while results from EME are for a single frequency per simulation.

Therefore, in both solvers it is necessary to run sweeps to obtain the full S-parameter matrix at multiple frequency values, either by sweeping over the input physical ports and modes or over the frequency points.

Use the examples in the Application Gallery to determine the appropriate solver for your component. This section describes the recommended approach for S-parameter extraction and export to INTERCONNECT once you have determined the appropriate solver.

In this section we describe how to:

1.For FDTD, calculate the S-parameter matrix and save it to a file with the appropriate format to be imported in an Optical N-port S-parameter primitive element in INTERCONNECT.

2.For EME, calculate the S-parameter matrix and save it to a file with the appropriate format to be imported in an Optical N-port S-parameter primitive element in INTERCONNECT.

3.Optionally, more accurately account for group delay of the device. This is important in elements that will be part of phase-sensitive circuits (for example, cavities, resonating structures and interferometric devices). We will use the group delay specification in the S-parameter file, required for the group delay option in the digital filter of the Optical N-port S-parameter element in INTERCONNECT.

FDTD solver:

1.Check that the port objects in the simulation have been set up correctly to collect all the data required for the S-parameter matrix. In particular, make sure that all the desired modes are selected at the ports; the “user select” option in the Modal Properties allows you to pick multiple modes (for example, fundamental TE and TM modes) when necessary for the full S-parameter matrix.

2.Create a S-parameter matrix sweep item in the “Optimizations and Sweeps” window with the following settings:

▪S-Matrix Setup tab:

•The modes selected at each port are displayed in this tab. Make sure the list is consistent with the physical ports and modes you want to include in the S-parameter matrix (see step 1).

•The option “Excite all ports” is enabled by default, which means that a simulation will be run for every port and mode in the list to calculate each column of the S-parameter matrix. To reduce the number of simulations you can map between columns taking advantage of symmetries in the structure. If that is the case, disable the “Excite all ports” option and click on the “Auto symmetry” button to let the solver determine the appropriate mapping based on the port locations. In some cases adjustments to the mapping are required so it is always recommended to review the settings after applying auto symmetry.

▪INTERCONNECT Export Setup tab:

•Select the “Custom define” option to define the mode ID and port location properties.

•The file export automatically assigns a mode ID to each mode. It is often necessary to manually set the mode ID so they are consistent with the orthogonal identifier convention to be used in INTERCONNECT. This is important, for example ,when connecting the Optical N-port S-parameter element to an Optical Network Analyzer (ONA) in INTERCONNECT. The setting for the excitation mode in the ONA is the “orthogonal identifier” property in the “Waveguide” section. The typical convention is to use “orthogonal identifier” = 1 for the fundamental TE mode and “orthogonal identifier” = 2 for the fundamental TM mode.

3.After running the S-parameter sweep, the “Visualize” and “Export to INTERCONNECT” options become available in the context menu of the S-parameter sweep item. The “Export to INTERCONNECT” option allows you to generate a text file with the S-parameters. This file has the format required by the Optical N-port S-parameter primitive in INTERCONNECT without group delay specification, which is appropriate when the phase or group delay of the device are not critical for the compact model.

4.If phase and group delay are critical, the group delay can be estimated by calculating the slope of the phase as a function of frequency at the center frequency/wavelength. The script add_group_delay_FDTD_Spar.lsf takes as input the S-parameter text file created in step 3 and modifies it to include the group delay. For more information see the section S-parameter extraction for group delay option in INTERCONNECT below.

EME solver:

1.Create a parameter sweep in the “Optimizations and Sweeps” window with “model::EME::wavelength” as a parameter and “::model::EME::user s matrix” as a result; this sweep will run multiple simulations with different wavelengths. An alternative is to use the “Wavelength sweep” feature in the EME Analysis window, which will use only one simulation to quickly estimate the S-parameters over a given wavelength range; however, this feature should not be used to obtain final results as the accuracy of the calculation depends on the type of simulated structure.

2.Examples that use EME include a script that exports S-parameters in the format required by the Optical N-port S-parameter primitive. In the script you need to provide the INTERCONNECT port definitions (names and locations) and the S-parameter file name. The following actions are performed by the script:

▪Run the parameter sweep for wavelength, collect the results and write them in the text file with the appropriate format.

▪In situations where the phase and group delay are critical, the script will also run the simulation and the EME propagation for the center wavelength with the option “calculate group delays” enabled. The group delay is added to the text file and the S-parameter phase is corrected as explained below in the section S-parameter extraction for group delay option in INTERCONNECT.

▪The output of the script is the S-parameter text file in the correct format ready to be imported in INTERCONNECT.

Some photonic circuits rely on controlling the phase of the optical signal traveling through the different components; therefore, the phase of the S-parameters can be very important. There are a couple of challenges to capture the phase correctly in the Optical N-port S-parameter primitive in INTERCONNECT:

•Spectral sampling of the phase in the S-parameter data: Typically, the variation of the phase as a function of frequency is mostly linear with a slope that is proportional to the group delay (the delay of the optical signal between two ports). Therefore, in long devices the phase can change very quickly with frequency, and a fine sampling in frequency would be necessary to avoid phase changes greater than \( 2 \pi \) so that the phase has the correct slope after unwrapping it. This is important for both frequency- and time-domain simulations in INTERCONNECT.

•Capturing the correct group delay and frequency dependence of the S-parameters in INTERCONNECT time-domain simulations: For time-domain simulations the most common type of filter used in the Optical N-port S-parameter primitive is the finite impulse response (FIR) one, which has multiple options for the estimation of taps. The choice depends on what is the most significant aspect of the element response to be captured in the circuit simulations. In elements where the phase is important, typically it is necessary to capture both the group delay and the frequency dependence of the S-parameters correctly.

The “group delay” option for the “number of taps estimation” in the digital filter settings of the Optical N-port S-parameter primitive provides solutions to both challenges:

•The group delay is captured correctly without requiring a fine spectral sampling of the S-parameters by using group delay estimates provided by the user for each S-parameter matrix element at the center frequency/wavelength.

•In most cases, the FIR filter in time-domain simulations can capture the group delay and frequency dependence of the element response more accurately with the "group delay" option, compared to other options for tap estimation. With the "group delay" option, the number of taps available to the filter depends on the sample rate of the INTERCONNECT simulation and the group delay estimate provided by the user: for a given sample rate, there are more taps available for larger group delays. Since most photonic devices are long compared to the wavelength of the signal, typically the number of available taps is not a concern. Only in cases where the device is relatively short, it might be necessary to increase the sample rate.

Compared to the regular file format for S-parameter files in INTERCONNECT, the group delay option requires a slightly different format, which includes the group delay estimate in the header of each S-parameter data block; furthermore, the phase provided in the file needs to be adjusted to remove the linear contribution associated with the group delay estimate. The previous section describes how this type of file can be generated when using FDTD or EME for the S-parameter extraction. For more information on the group delay option see the discussion here.

•KB page: S-parameter simulator (SPS)

•KB page: S-parameter matrix sweep feature in KB

•KB page: S-parameter file formats