Inverse design of 3D nanophotonic devices with structural integrity using auxiliary thermal solvers

Oliver Kuster; Yannick Augenstein; Roberto Narváez Hernández; Carsten Rockstuhl; Thomas Jebb Sturges

doi:10.1515/nanoph-2025-0067

Article Open Access

Inverse design of 3D nanophotonic devices with structural integrity using auxiliary thermal solvers

Oliver Kuster , Yannick Augenstein , Roberto Narváez Hernández , Carsten Rockstuhl and Thomas Jebb Sturges

Published/Copyright: April 15, 2025

Published by

Become an author with De Gruyter Brill

Submit Manuscript Author Information

From the journal Nanophotonics Volume 14 Issue 9

Abstract

3D additive manufacturing enables the fabrication of nanophotonic structures with subwavelength features that control light across macroscopic scales. Gradient-based optimization offers an efficient approach to design these complex and non-intuitive structures. However, expanding this methodology from 2D to 3D introduces complexities, such as the need for structural integrity and connectivity. This work introduces a multi-objective optimization method to address these challenges in 3D nanophotonic designs. Our method combines electromagnetic simulations with an auxiliary heat-diffusion solver to ensure continuous material and void connectivity. By modeling material regions as heat sources and boundaries as heat sinks, we optimize the structure to minimize the total temperature, thereby penalizing disconnected regions that cannot dissipate thermal loads. Alongside the optical response, this heat metric becomes part of our objective function. We demonstrate the utility of our algorithm by designing two 3D nanophotonic devices. The first is a focusing element. The second is a waveguide junction, which connects two incoming waveguides for two different wavelengths into two outgoing waveguides, which are rotated by 90° to the incoming waveguides. Our approach offers a design pipeline that generates digital blueprints for fabricable nanophotonic materials, paving the way for practical 3D nanoprinting applications.

Keywords: inverse design; 3D additive manufacturing; nanophotonics

1 Introduction

3D nanoprinting enables the fabrication of photonic devices with voxel sizes below the wavelength of light, achieving intricate structures that seamlessly integrate features from the nanoscale to millimeter or even centimeter scale [1], [2], [3], [4], [5]. Fabrication happens upon two-photon absorption in a photoresist within a tiny volume into which a writing laser is focused. Two-photon absorption triggers the polymerization of the photoresist, and a structure is formed. Micro-optical components are widely explored with such fabrication techniques [6], [7], [8], [9], [10], [11], [12], [13], [14], [15], [16], [17], [18], and complex photonic systems with a non-intuitive design have also been realized [19], [20], [21], [22], [23]. By raster scanning the laser across a volume, free-form designs in three dimensions can be realized, where every voxel constitutes a design parameter, i.e., it may or may not be polymerized. In principle, this high dimensionality promises structures with almost arbitrary functionalities due to the many degrees of freedom available.

However, the high dimensionality of the parameter space is simultaneously a blessing and a curse. All degrees of freedom must be efficiently adjusted while accommodating fabrication and application constraints in the design. Due to the complexity of such devices, efficient numerical methods are indispensable to solve the inverse problem [24], [25], [26], [27], [28], [29], [30], [31], [32]. For free-form 3D nanoprinting, global optimization techniques are no longer applicable, and only gradient-based inverse design methods can be used. Here, the adjoint formalism is particularly beneficial as it allows us to iteratively update a design with only two simulations of our system to show improved performance in relation to a chosen objective function. The adjoint formalism is the core of topology optimization, which is a highly efficient method in designing nanophotonic devices [33], [34], [35], [36], [37], [38], [39], [40], [41], [42], [43].

However, not just the optical functionality matters in the design. Besides the usual fabrication constraints, such as minimum feature size, a prime design challenge is ensuring the structural integrity of the 3D nanoprinted devices. We must ensure the final design is fully connected and does not collapse on itself. In other words, free-floating components are not an option for realistic devices. Simultaneously, we must ensure that the design does not include cavities. Cavities would trap the undeveloped photoresist, which is detrimental to the optical functionality.

One way to ensure the structural integrity of the material is by using physics-based constraints [44], [45], [46]. In this work, we present a method to enforce the structural integrity of 3D photonic structures by using an auxiliary heat-diffusion solver. Comparable approaches have already been discussed in the past to optimize 2D structures [47], [48]. While structural solvers exist, we require a framework which is differentiable and can easily be combined with an electromagnetic simulation. The heat solver provides a direct and computationally cheaper approach to enforcing global connectivity, which – at the nanoscale – is typically sufficient to also guarantee mechanical stability. In addition, the heat based solver allows us to easily enforce connectivity of the void, which is a requirement for the 3D additive manufacturing process of nanophotonic devices. Using structural solvers which rely on vectorial forces [44] would become difficult and restrict the design space too much.

To evaluate how well-connected a given structure is, we consider the written structure to be a heat source. By solving a heat-diffusion equation, we then study how well the structure can dissipate its heat to strategically placed heat sinks at the outer rims of the design region. As we permit diffusion only inside the written material, the material must be well-connected to reach a low final temperature. By calculating the total temperature of the device (upon integration of the temperature distribution across the design), we measure how well connected the structure is. This process is repeated for the void, which ensures that no cavities are formed within the design. To be clear, the thermal solver is purely fictitious and is used only within the optimization pipeline to ensure connectivity.

These two separate heat-diffusion problems are solved with a finite-element method. The electromagnetic simulation is solved with a finite-difference frequency-domain (FDFD) method [49]. For each of these three simulations, an adjoint problem is also solved to enable our multi-objective optimization workflow with gradient-based optimization. The overall objective function is constructed as the sum of three sub-objectives, each of which is normalized relative to a predefined threshold and processed through a softplus function. The normalization ensures balanced contributions from the electromagnetic, material, and void simulations [50] while the softplus function deactivates further thermal optimization once structural connectivity is established. The additional computational cost of the thermal simulations is comparatively low. Still, the outcome guarantees that the resulting devices possess structural integrity. Therefore, these designs provide a digital blueprint for direct fabrication. In addition, we find designs that are almost as optically performant as designs optimized only for their optical response. With this, we develop a critical component of a design pipeline for 3D nanoprinted photonic devices that will find applications in perceiving a future generation of novel structures that can serve societal needs (see Figure 1).

Figure 1:

Sketch of our design procedure. The initial design, expressed in terms of a material density ρ ∈ [0, 1], is mapped to three different material distributions. Firstly, it is mapped to the permittivity ɛ(ρ) of the dielectric material to be printed, from which we calculate the electromagnetic response. Additionally, a fictitious heat-diffusion through the system is studied in two dedicated auxiliary simulations by assuming the written material and the void as heat sources, respectively. The material density is mapped to a thermal conductivity κ(ρ) for the material, and κ(1 − ρ) for the void. From all three simulations, an objective function is evaluated. Each sub-objective is then renormalized, passed through a softplus-function, and summed to the final objective function. In combination, these sub-objectives aim to balance the goal of maximizing a given optical figure of merit while ensuring structural connectivity of the devices with no cavities. The softplus-function ensures that the fictitious thermal performance of the design stops being optimized once structural connectivity is achieved. The gradients of the objective function relative to the material density at every point in the design space are then used to update the design iteratively and to optimize the digital blueprint with structural integrity.

Besides this introduction, the paper is structured into three sections. In the following section, we outline our general design pipeline based on topology optimization. Afterward, in Section Results, we demonstrate the applicability of our approach in two carefully selected devices. Finally, we conclude on our work in Section Discussion.

2 Topology optimization

Topology optimization, as used in our design pipeline, is a gradient-based, local optimization method. Given a figure of merit (also called the objective function) L ( ρ ) , topology optimization can formally be written as

(1) min ρ L ( ρ ) s.t. c i ( ρ ) = 0 , s.t. c j ( ρ ) < 0 ,

with equality constraints c _i(ρ) and inequality constraints c _j(ρ), i , j ∈ N . By calculating ∇ ρ L ( ρ ) and using gradient-based optimizers, a local minimum for the figure of merit can be found for any initial condition chosen as the starting point for the iterative design process. In our case, ρ represents a material density that can be mapped to either the relative permittivity or the thermal conductivity of our designs. In each case, it represents the material distribution we aim to optimize.

A forward simulation is performed to evaluate the objective function. Then, using the adjoint method, the gradients of the objective function with respect to each degree of freedom can be calculated efficiently at the cost of only one additional simulation, called a backward simulation [33], [51]. Notably, the computational complexity of calculating the gradients using the adjoint method is independent of the number of degrees of freedom, enabling the free-form design of structures characterized by many millions of parameters.

2.1 Forward simulation

To perform the forward simulation, the material density must first be mapped to the corresponding material distribution. We will optimize the density in the design region D ⊆ S ⊆ R 3 , which is embedded into the simulation domain S . The continuous material density is represented as ρ≔ρ(x, y, z), with ρ ( x 0 , y 0 , z 0 ) ∈ [ 0,1 ] , ( x 0 , y 0 , z 0 ) ∈ D . Discretizing the spatial distribution of ρ divides the design region into voxels, each of which can then be independently parameterized, enabling free-form structural design [52], [53], [54]. To map ρ to the physical parameters required for our simulations, we use multiple, continuously differentiable transformations. The transformations map the material density to the electric permittivity or the thermal conductivity, respectively, while also incorporating fabrication limitations.

As such, we first apply a Gaussian filter to ρ to enforce a minimum feature size of 100 nm [14] for our structures.

(2) ρ ̃ = ρ ∗ w ( σ ) with w ( σ ) = e − x 2 + y 2 + z 2 2 ⋅ σ 2 .

The Gaussian filter is given by a discrete convolution of ρ with a Gaussian kernel w(σ), with a kernel size of l = 2 ⋅ 4 ⋅ σ + 1 . The minimum feature size itself can be approximated as 3 σ . Considering our minimum feature size of 100 nm and the resolution of our simulations, we use σ ≈ 2.3 px.

We are only interested in designs consisting of either “void” or “material” with no intermediate values. This corresponds to only two electric permittivities or thermal conductivities. The binarization of the material density is imposed by a soft threshold filter [24], [55]

(3) ρ ̃ ̂ = tanh ( β α ) + tanh ( β ( ρ ̃ − α ) ) tanh ( β α ) + tanh ( β ( 1 − α ) ) .

Here, β controls the steepness of the function and quantifies the degree of binarization, while α defines its center. We set α = 0.5 in our simulations and gradually increase β from 1 to 30 throughout the optimization. Initially, β is kept low to allow broad exploration of the parameter space. As the optimization progresses, we incrementally increase β to enforce binarization, ensuring convergence towards a fabricable device.

The last set of transformations interpolates the binarized density to the physical quantities using a linear mapping. For the electromagnetic simulation, this is given by

(4) ε r = ε min + ρ ̃ ̂ ( ε max − ε min ) ,

where ɛ _min and ɛ _max are the electric permittivity of the void and the material, respectively. For the temperature simulation, we use

(5) κ = κ min + ρ ̃ ̂ ( κ max − κ min ) ,

where κ _min and κ _max refer to the thermal conductivity of the void and the material, respectively [44]. Please note that the heat conductivities considered are purely fictitious and do not correspond to the actual material that will be printed. We choose to keep the values representing physical properties as constant and use simulation parameters as our hyperparameters for the optimization. The values are considered in the auxiliary heat solver only and are chosen to ensure connectivity among all domains. Once we have the material distributions, we can consider them in the respective simulation. From these simulations, we can judge the suitability of a given device for its anticipated purpose. Two types of simulations are performed: the electromagnetic simulation and the heat-diffusion simulations.

We simulate the electromagnetic response using an FDFD Maxwell solver [49] and evaluate the total electric field distribution E ⃗ ( x , y , z , ω ) at the frequency ω at every point in space for a given illumination. We take our design region D and embed it into the full simulation domain S . The total electric field distribution is then used to calculate the electromagnetic figure of merit L EM ≔ L EM ( ρ ̃ ̂ ) . For the purpose of readability, throughout this manuscript we omit the explicit dependency of our figures of merit on ρ ̃ ̂ . The electromagnetic figure of merit encodes the desired functionality of the device into a single scalar value that indicates how good the optical performance is. Therefore, by finding a local maximum of L EM , we get a density distribution of ρ ̃ ̂ which manipulates the electromagnetic wave in the desired way.

In addition to the electromagnetic simulation, two heat-diffusion simulations are performed, from which we judge the connectivity of the material and the voids. Specifically, we solve Poisson’s equation

(6) − κ ( ρ ̃ ̂ ) ∇ 2 u ( x , y , z ) = q ( ρ ̃ ̂ ) s.t. u ( x 0 , y 0 , z 0 ) = 0 for ( x 0 , y 0 , z 0 ) ∈ U ⊆ ∂ D , s.t. ∇ u ( x 0 , y 0 , z 0 ) = 0 for ( x 0 , y 0 , z 0 ) ∈ ∂ D \ U ,

where u(x, y, z) is the scalar temperature distribution of the system, κ ( ρ ̃ ̂ ) the thermal conductivity, and q ( ρ ̃ ̂ ) the heat generation rate, i.e. corresponding to the heat sources. A steady-state temperature distribution is obtained by solving equation (6) with an in-house finite element solver using an H8 element. For readability purposes, we implicitly assume the boundary conditions to hold true for our optimization problem and do not state them explicitly anymore. When optimizing for material connectivity, we designate material regions as heat sources, i.e. q ( ρ ̃ ̂ ) = ρ ̃ ̂ , and place heat sinks at points U along the boundary. We strategically position the heat sinks within regions where material is known to be desirable, like waveguide ports. We assume Neumann boundary conditions for ∂ D \ U . Specifying heat sinks will lead to Dirichlet boundary conditions for U. Note that the heat sinks do not need to be on the boundaries but can, in principle, be placed anywhere. We assume a small background heat conductivity everywhere in space to avoid the temperature diverging for non-connected structures. Such a setting allows us to reach an equilibrium with a high but finite temperature in spatial domains not connected to the sinks. The finite element method requires us to set a non-zero heat conductivity for κ _min. Our chosen value range for the thermal conductivity is κ min = 1 0 − 5 M W μ m ⋅ K and κ max = 1 M W μ m ⋅ K . Structural integrity is achieved by minimizing the total temperature in the system, encouraging the optimizer to connect all material regions to the heat sinks, thereby reducing the temperature within the structure.

The simulation that ensures connectivity of the voids uses q ( ρ ̃ ̂ ) = 1 − ρ ̃ ̂ as heat sources, enforcing connectivity of the void and thus avoiding the formation of cavities that trap unwritten photoresist in the final structure. The edges of the design domain that were not considered heat sinks in the material-connectivity simulation will now act as heat sinks in the void-connectivity simulation.

In contrast to the electromagnetic simulation, in the heat-diffusion simulations, we consider only the design region D and not the full simulation domain. Such a restriction saves computational time and memory, as only D is relevant for the structural integrity.

Our figure of merit for the thermal simulations is the integrated final temperature distribution within the material (or respectively void). We define these figures of merit as L material = L heat ( q = ρ ̃ ̂ ) and L void = L heat ( q = 1 − ρ ̃ ̂ ) respectively, where

(7) L heat = ∫ D u ( x , y , z ) d x d y d z .

This implies that we strive to minimize the final temperature within the design region in both simulations.

2.2 Optimization problem

Using the three forward simulations, we get three contributions to our figure of merit, which must be balanced so that the optimizer weights each one appropriately. Since all these simulations yield vastly different numerical values, we first renormalize each sub-objective to become comparable.

For L material , we define

(8) L m a t e r i a l n = L material − L material thresh L material thresh ,

to renormalize the contribution from the material heat-diffusion equation. We introduce a threshold value L material thresh where the heat-dissipation is good enough that the material has a high probability of being fully connected. Therefore, the sign of L material n provides a binary measure of whether the device is connected enough or not. This will allow us to formulate a multi-objective optimization scheme that deprioritizes those sub-objectives that have reached their threshold value. In a similar fashion, we renormalize L void to L void n using a threshold value L void thresh .

Since we aim to minimize both L material n and L void n , while the contribution from L EM has to be maximized, L EM is renormalized as

(9) L EM n = L EM thresh − L EM L EM max − L EM thresh .

We choose L EM thresh to be roughly the same as L EM of a design which was optimized only for its optical performance. Such a choice ensures that we can still find highly efficient designs, even with the more restricted parameter space, since the optimization will mostly focus on L EM once the threshold values for the material and the void have been reached. We choose L EM max to be larger than any L EM during the optimization, so we can ensure that L EM n stays larger than −1 for big values of L EM .

As discussed, we chose L material n and L void n such that they become negative once the threshold value is reached, which corresponds to fully connected material (or void). As a subsequent step, we apply a softplus function to the three sub-objectives that assigns a vanishing weight for negative values,

(10) softplus ( L n ) = ln ( 1 + e L n ) , if L n < 0.4 , L n , otherwise .

We also add a regularization term as a fourth sub-objective, which penalizes non-binarized voxels. This binarization penalty is given as

(11) L binary = ∫ D 4 ⋅ ρ ̃ ̂ ( x , y , z ) ⋅ ( 1 − ρ ̃ ̂ ( x , y , z ) ) d x d y d z ,

and renormalized to L binary n = L binary − L binary thresh L binary thresh using L binary thresh = 1 0 − 2.3 μ m 3 . We only turn on the binarization penalty once we are close to convergence. At this point, most of the voxels are binarized anyway. By turning on the binarization penalty, we can binarize the few remaining non-binarized voxels without majorly impacting the properties of the design.

To combine the contributions from each individual figure of merit into a single scalar, we use the l ² norm. Our full optimization problem is then given as

(12) min ρ L ( ρ ̃ ̂ ) = ∑ i | softplus ( L i n ( ρ ̃ ̂ ) ) | 2 s.t. ( ∇ × ∇ × − ω 2 ε ( ρ ̃ ̂ ) ) E ⃗ ( ω ) = − i ω J ⃗ ( ω ) , s.t. − κ ( ρ ̃ ̂ ) ∇ 2 u ( x , y , z ) = q ( ρ ̃ ̂ ) , s.t. − κ ( 1 − ρ ̃ ̂ ) ∇ 2 u ( x , y , z ) = q ( 1 − ρ ̃ ̂ ) , s.t. 0 ≤ ρ ̃ ̂ ( x , y , z ) ≤ 1 ,

with i ∈ {EM, material, void, binary}. We choose to formulate our figure of merit using a softplus function instead of simply using prefactors for two reasons. First, the softplus function has a built in cutoff. This means, that once a certain threshold value is crossed, the respective sub-objective will barely be considered in the full figure of merit. Since we are not interested in the thermal performance of our device, the difference between two values of L heat will not make a significant difference in the structural integrity once the connectivity threshold is crossed. This helps us to avoid getting stuck in a local minimum where the thermal performance is optimized but the electromagnetic performance is not. Second, since we are only interested in crossing certain threshold values for our subobjectives, our figure of merit is consequently less sensitive to the choice of our hyperparameters, avoiding the need for a proper hyperparameter optimization as we discuss later on.

The main disadvantage of this approach is that it requires some prior knowledge about the threshold values. Only with this prior knowledge can a proper renormalization be done. Still, the method is fairly robust. A rough estimate of the threshold values is usually sufficient, saving time compared to a dedicated hyperparameter optimization. The choice of the values for L heat thresh mostly depend on the surface the respective heat sinks occupy as well as the overall volume of the structure. If the heat sinks occupy more surface, the system has more efficient means to disperse the heat, making it less prone to the initial value of L heat n . A bigger volume for the structure means more heat sources. Nonetheless, we have found, that choosing L heat thresh , so that the initial value for L heat n is in (0, 1) works sufficiently well for most structures. A good initial value tends to be L heat n = 0.5 . To get these initial values for L heat n , one heat simulation for the material and one heat simulation for the void on the initial structure are required, which gives us L heat which we can then use to determine L heat thresh to give us our desired initial L heat n .

There are various options to choose L EM thresh . We want to choose L EM thresh so that L EM n can vary in a suitable range. If we already know the L EM max beforehand, for instance the normalized transmission in a waveguide, we can choose L EM thresh so, that L EM n = 0 for L EM = L EM thresh initially. If we have no prior knowledge about what constitutes a good value for L EM thresh , we usually have to do a full optimization without the heat solver to determine a good guess for L EM thresh . Then we use L EM max to tune the initial L EM n . Choosing an initial value of L EM n = 1 tends to be a good starting point.

All simulations are done using a single NVIDIA A100 GPU and run for a maximum of 120 optimization steps. The optimization is split into five stages with 20 optimization steps each, or until convergence (absolute change of less than 10⁻⁴ or relative change of less than 10⁻⁵ per iteration), during which we gradually increase β from 1 to 30 with every stage. Then, a final stage is done with β = 30, and we turn on the binarization penalty. We discretize our spatial extent using 40 px per micrometer. Our initial design is set to be ρ = 0.5 everywhere inside the design region. One full optimization takes roughly half a day to two days, depending on the simulation size and number of sources.

3 Results

We showcase the outlined design pipeline with two examples. First, we discuss a focusing element. Second, we discuss a specific waveguide coupler.

3.1 Focusing element

Our first example is a focusing element with an operational wavelength of 1.55 μm. A permittivity of ɛ _r = 2.25 characterizes the written polymer. Maximizing the normalized electromagnetic field at the focal point for a given illumination achieves the desired focusing behavior. Such a figure of merit can be expressed as

(13) L EM = | E x ( x f , y f , z f ) | 2 ∫ S | E x ( x , y , z ) | 2 d x d y d z ,

where E _x(x, y, z) is the electric field component in the x-direction and (x _f, y _f, z _f) is the focal point. The other electric field components are negligible because the incident field is an x-polarized plane wave propagating in the +z-direction. The volume within which we optimize the material distribution is a cuboid with side lengths of 4 μm × 4 μm × 2 μm, resulting in roughly 2 million parameters for our chosen resolution. The full simulation domain has a size of 5  μm × 5 μm × 6 μm with the design region placed at the center. Additionally, we use PMLs with width 0.5 μm to avoid any unphysical scattering effects. We also pad our design region with material in the x-y-direction, creating a “frame”. This padding acts as a base for fabrication and supports mechanical stability. So we want to connect our heat sinks for the material to said frame. The heat sinks for the void are placed on the top and on the bottom of the design in the z-direction. The focal point is chosen to be 1 μm behind the terminating interface on the optical axis.

To get insights into the optimization process, we perform two parameter sweeps. In both parameter sweeps we use L EM as a measure of the performance, as we are mostly interested in how a device (with structural integrity) optically performs in comparison to a structure optimized exclusively for its optical functionality.

First, we use a constant value for L material thresh , and optimize several devices for different values of L void thresh . We select our threshold value such that (given the initial distribution of ρ = 0.5 everywhere in the design region) the initial value of L material n is approximately 0.5. Second, we do the same by fixing sweeping through L material thresh whilst initializing L void n ≈ 0.5 in the same way. The results of the parameter sweeps can be seen in Figure 2. The figure shows the electromagnetic figure of merit ultimately obtained for each optimized device. As an orientation, we also plot the electromagnetic figure of merit for an optimized device that was not designed for structural integrity. This figure of merit (blue dashed line in Figure 2) serves as an indication of the largest possible value that is achievable when optimizing only the electromagnetic response.

$Figure 2: Different focusing devices designed with structural integrity. A parameter sweep is done by varying the threshold value for L material thresh ${\mathcal{L}}_{\text{material}}^{\text{thresh}}$ while L void threshold ${\mathcal{L}}_{\text{void}}^{\text{threshold}}$ is constant (red) and by varying the threshold value for L void thresh ${\mathcal{L}}_{\text{void}}^{\text{thresh}}$ while L material thresh ${\mathcal{L}}_{\text{material}}^{\text{thresh}}$ is constant (black). Their performance as focusing devices, evaluated using L EM ${\mathcal{L}}_{\text{EM}}$ is compared to a device designed by maximizing only L EM ${\mathcal{L}}_{\text{EM}}$ (blue dotted line). Some selected optimized devices are also shown, where the markers (colored square or triangle) indicate the corresponding threshold values. The displayed devices show the material structure and the field amplitude above a cut-off value. We also display the device optimized for its optical performance only. It can be seen that this device does not possess structural integrity (see Supplementary Material for a camera angle that better shows this). The best-performing design with structural integrity (red triangle) can be seen in Figure 3 in detail.$

Figure 2:

Different focusing devices designed with structural integrity. A parameter sweep is done by varying the threshold value for L material thresh while L void threshold is constant (red) and by varying the threshold value for L void thresh while L material thresh is constant (black). Their performance as focusing devices, evaluated using L EM is compared to a device designed by maximizing only L EM (blue dotted line). Some selected optimized devices are also shown, where the markers (colored square or triangle) indicate the corresponding threshold values. The displayed devices show the material structure and the field amplitude above a cut-off value. We also display the device optimized for its optical performance only. It can be seen that this device does not possess structural integrity (see Supplementary Material for a camera angle that better shows this). The best-performing design with structural integrity (red triangle) can be seen in Figure 3 in detail.

In addition, Figure 2 shows various focusing devices where different values of L heat thresh were chosen, as well as the focusing device that was optimized by maximizing only L EM . Notably, all designs using an auxiliary heat-dissipation solver possess structural integrity in the material and the void. In contrast, the device that considered only the electromagnetic response does not.

When looking into the results, we observe that if we choose the threshold value too small, the optimization process ends up in structures with ρ = 0 or ρ = 1 everywhere in the design region. This can be seen clearly in the exemplary devices with the red and black squares, where the optimized structure consists of only a void (red square) or only material (black square). This happens because normalizing to a tiny threshold value results in a thermal objective function that has a much higher value than the electromagnetic objective function by orders of magnitude. Since this happens in the linear regime of the softplus function, the gradient will also point towards a solution which optimizes the thermal objective over the electromagnetic objective. For example, when L material thresh ≪ 1 , we obtain L heat n ≫ 1 , while L EM n ≈ 1 . The optimizer strongly prioritizes minimizing L material n , effectively disregarding L EM n (and L void n ). As a result, the optimization converges to the trivial solution ρ = 0 everywhere (i.e., no material and therefore no heat sources). The opposite applies when L void thresh ≪ 1 . The optimization fills the entire design region with material (ρ = 1 everywhere), resulting in no heat generated by the void. Again, this is just a consequence of the fact that if there are no voids, no spatial domain acts as a heat source in the void-connectivity simulation, and thus the temperature stays low everywhere. But of course, these trivial solutions that only satisfy a single thermal sub-objective, whilst disregarding the EM functionality, are not what we want. Furthermore, the optimization gets stuck in the local minimum of the trivial solution if it converges towards it initially, even if better local minima are available (see Supplementary Material).

Increasing the threshold leads to non-trivial structures which can function as proper focusing devices. Interestingly, it seems to be possible to increase the threshold value for either L material thresh or L void thresh to high values, without losing the structural integrity of material and void. Examples of such optimized devices are also shown in Figure 2, and they are indicated with the red and black triangle, respectively. However, that observation is most likely problem-specific. In the case of choosing high L material thresh , it is even possible to find designs that are almost as good as a design that was optimized by maximizing only L EM .

In contrast, when considering only L EM in the optimization, the resulting structure is free-floating. We can achieve structural integrity only after adding the auxiliary thermal solvers. The best performing device with structural integrity can be seen in more detail in Figure 3. The electromagnetic object function reaches a value of L EM = 16.7 μ m − 3 . Compared to a value of L EM = 17.4 μ m − 3 for the design optimized only for its optical functionality, we can conclude that the enforcement of structural integrity led only to a minor degradation in optical performance. However, we have achieved a fully connected, self-sustaining structure free of voids.

Figure 3:

Best performing focusing device with enforced structural integrity, shown from different angles: Tilted (a), top view (b), and side view (c). The electric field distribution |E _x| is shown in (d).

3.2 Waveguide coupler

The second design is a waveguide coupler. We consider two incoming waveguides for two different wavelengths λ ₁ = 700 nm and λ ₂ = 600 nm. The waveguides are made from a polymer characterized by a permittivity of ɛ _r = 2.25. Each waveguide couples to an outgoing waveguide rotated by 90°. We place our heat sinks for the material at the waveguide ports, which is where we are certain that material is required for the optical performance. The heat sinks for the void are placed on the boundaries of the design domain except where the waveguides are. Our design region is defined by a cube with sidelength 2 μm, resulting in around half a million parameters for our chosen resolution. The full simulation domain has a size of 3 μm × 3 μm × 6 μm with the design region at its center. We use PMLs with width 0.5 μm. As a figure of merit, we choose to maximize the intensity of each field in their respective target waveguides

where E x 1 ( x , y , z ) and E x 2 ( x , y , z ) are the x-component of the electric field for wavelength λ ₁ and λ ₂. In each waveguide, an x-polarized TE₀-waveguide mode is propagating, which serves as the illumination. W 1 and W 2 are the cross-section of the outgoing waveguides where each wavelength should couple to. Each waveguide is 2 μm long and has a sidelength of 0.5 μm and is also displaced by 0.75 μm from the center in their respective direction. The fields are evaluated inside the outgoing waveguides.

While the optimization that exclusively optimizes the electromagnetic response already favors structures that mostly possess structural integrity, the final designs may still contain free-floating artifacts and cavities. As with the focusing device, we sweep through both L material thresh and L void thresh , while fixing the other. In contrast to the focusing device, we enforce a higher degree of connectivity by setting the threshold values such that the initial values of L material n and L void n are approximately 1 (twice as large as for the focusing device). Lower values still tend to produce free-floating artifacts, which we want to avoid. This is most likely due to the stepwise binarization, which can produce “islands” that are initially connected with non-binarized values, and their connections are slowly eroded.

The results of the parameter sweep can be seen in Figure 4. As with the focusing device, choosing the value for L material thresh too low results in structures where ρ = 0 almost everywhere in the design region. Therefore, as there is no structure which connects the input with the output waveguide, the electromagnetic objective function is rather low. There are non-trivial structures that appear for low threshold values of L material thresh during the parameter sweep. These non-trivial solutions appear due to the placement of the heat sinks, as both input and output waveguides are heat sinks, resulting in solutions which do not connect the waveguides but still minimize the figure of merit with respect to the heat and influence the light in a very non-efficient manner. Since waveguide-like structures tend to be the optimal solution to the problem, we will end up with fully connected structures for higher threshold values L material thresh . At the transition, the bumps begin to elongate towards one another, forming rudimentary guiding paths that outperform the trivial solution but still do not fully connect the input waveguides to the output waveguides. Only once the waveguide ports are connected we can actually guide the light from the input waveguides to the output waveguides in an efficient manner.

$Figure 4: Different waveguide couplers designed with structural integrity. A parameter sweep is done by varying the threshold value for L material thresh ${\mathcal{L}}_{\text{material}}^{\text{thresh}}$ while L void threshold ${\mathcal{L}}_{\text{void}}^{\text{threshold}}$ is constant (red) and by varying the threshold value for L void thresh ${\mathcal{L}}_{\text{void}}^{\text{thresh}}$ while L material thresh ${\mathcal{L}}_{\text{material}}^{\text{thresh}}$ is constant (black). Their performance as waveguide couplers evaluated using L EM ${\mathcal{L}}_{\text{EM}}$ is compared to a waveguide coupler designed by maximizing only L EM ${\mathcal{L}}_{\text{EM}}$ (blue dotted line). Some selected optimized devices are also shown, where the markers (colored square or triangle) indicate the corresponding threshold values. The displayed devices show the material structure and the field amplitude above a cut-off value. We also display the device optimized for its optical performance only. It can be seen that this device does not possess structural integrity (see Supplementary Material for a camera angle that better shows this). The best-performing design with structural integrity (red triangle) can be seen in Figure 5 in detail.$

Figure 4:

Different waveguide couplers designed with structural integrity. A parameter sweep is done by varying the threshold value for L material thresh while L void threshold is constant (red) and by varying the threshold value for L void thresh while L material thresh is constant (black). Their performance as waveguide couplers evaluated using L EM is compared to a waveguide coupler designed by maximizing only L EM (blue dotted line). Some selected optimized devices are also shown, where the markers (colored square or triangle) indicate the corresponding threshold values. The displayed devices show the material structure and the field amplitude above a cut-off value. We also display the device optimized for its optical performance only. It can be seen that this device does not possess structural integrity (see Supplementary Material for a camera angle that better shows this). The best-performing design with structural integrity (red triangle) can be seen in Figure 5 in detail.

The behavior is different when changing L void thresh . Unlike in the case of the focusing element, we do not encounter a solution with ρ = 1 everywhere in the design region. The heat sinks for the void make up the majority of the surface of the design region, allowing the heat to be dissipated more easily than for the material. This, in turn, means that the threshold values required for trivial solutions are magnitudes lower than the range we looked at. Using a high value for L void thresh slightly improves the performance while also using less material in the final design. The best-performing design can be seen in Figure 5. This design can reach L EM = 150 μ m − 1 in comparison to L EM = 152 μ m − 1 for a design which was optimized by maximizing only L EM . Therefore, by using an auxiliary heat-dissipation simulation, we can remove free-floating islands and cavities while reducing the required material, all with minimal impact on the optical performance.

$Figure 5: Best performing waveguide crossings with enforced structural integrity, shown from different angles: Tilted (a), top view (b), and side view (c). The electric field distributions | E x 1 | $\vert {E}_{x}^{1}\vert $ and | E x 2 | $\vert {E}_{x}^{2}\vert $ with the device are shown in (d).$

Figure 5:

Best performing waveguide crossings with enforced structural integrity, shown from different angles: Tilted (a), top view (b), and side view (c). The electric field distributions | E x 1 | and | E x 2 | with the device are shown in (d).

Interestingly, both with and without the auxiliary heat solver, we find solutions that are not just two waveguides connecting the waveguide ports. Instead, we get structures with additional features. These additional features help to guide the light through the curved waveguides by catching the field which is propagating outside of the waveguides [56]. This leads to an increased performance, as even the tail end of the electromagnetic waves can be captured. Once we enforce structural integrity using our heat solver, we favor solutions with less material. The resulting less bulky structures still retain some of the additional features which catch the evanescent field outside the central waveguides, which seemingly still help to improve the performance of our waveguides. This phenomenon is similar to additional features appearing in 2D waveguide bends [50], [57], where resonator-like structures appear to improve the performance of the waveguide bend.

4 Discussion

This work introduces a robust methodology to ensure the structural integrity of 3D-printed nanophotonic devices through the integration of an auxiliary heat-diffusion solver within a gradient-based topology optimization framework. Our approach not only prevents disconnected or structurally unstable designs but also preserves the desired optical functionalities with minimal performance trade-offs.

By leveraging the heat-diffusion solver as a fictitious physics-based soft-constraint, we effectively enforce connectivity in both material and void domains. The use of renormalised figures of merit and the softplus function ensures that thermal optimisation ceases once structural integrity is achieved, allowing the optimizer to focus on the optical performance. In combination with the iterative binarization strategy, our approach results in designs that are directly fabricable, offering a streamlined pipeline for 3D nanoprinting applications.

The proposed framework was demonstrated through the design of two devices – a focusing element and a waveguide crossing. These case studies validate the effectiveness of our method in achieving structurally robust and optically performant devices. In both cases, the devices with structural integrity performed only marginally worse than those optimized for their optical performance only.

This study lays the foundation for a broader application of physics-based constraints in the design of complex 3D nanostructures. Future work could explore extensions of this methodology to include additional auxiliary (or indeed real) physical constraints. Further investigation into adaptive weighting schemes and thresholding strategies could also enhance the robustness and versatility of the optimization process. By addressing the longstanding challenge of structural integrity in 3D nanophotonic design, this work represents a significant step forward in enabling practical, high-performance devices for real-world applications.

Corresponding author: Oliver Kuster, Institute of Theoretical Solid State Physics, Karlsruhe Institute of Technology (KIT), Kaiserstrasse 12, 76131 Karlsruhe, Germany, E-mail: oliver.kuster@kit.edu

Funding source: Alexander von Humboldt-Stiftung

Funding source: Deutsche Forschungsgemeinschaft

Award Identifier / Grant number: 390761711

Funding source: Carl-Zeiss-Stiftung

Acknowledgement

This work was performed on the HoreKa supercomputer funded by the Ministry of Science, Research and the Arts Baden-Württemberg and by the Federal Ministry of Education and Research.

Research funding: OK and CR acknowledge support from the German Research Foundation within the Excellence Cluster 3D Matter Made to Order (EXC 2082/1 under project number 390761711) and by the Carl Zeiss Foundation. TJS acknowledges funding from the Alexander von Humboldt Foundation.
Author contributions: All authors have accepted responsibility for the entire content of this manuscript and consented to its submission to the journal, reviewed all the results and approved the final version of the manuscript. CR, TJS and YA conceived the idea and surpervised the research. RNH has done initial investigations for 2D, OK wrote the code and designed the devices in 3D. OK wrote the manuscript. CR and TJS revised the paper.
Conflict of interest: Authors state no conflict of interest.
Data availability: The code to reproduce the datasets analyzed during the current study are available in the GitHub repository, https://github.com/OlloKuster/Structural_Integrity_3D/tree/main.

References

[1] J. Fischer and M. Wegener, “Three-dimensional optical laser lithography beyond the diffraction limit,” Laser Photonics Rev., vol. 7, no. 1, pp. 22–44, 2013, https://doi.org/10.1002/lpor.201100046.Search in Google Scholar

[2] C. Barner-Kowollik, et al.., “3D laser micro- and nanoprinting: challenges for chemistry,” Angew. Chem., Int. Ed., vol. 56, no. 50, pp. 15828–15845, 2017, https://doi.org/10.1002/anie.201704695.Search in Google Scholar PubMed

[3] V. Hahn, et al.., “Rapid assembly of small materials building blocks (voxels) into large functional 3d metamaterials,” Adv. Funct. Mater., vol. 30, no. 26, p. 1907795, 2020, https://doi.org/10.1002/adfm.201907795.Search in Google Scholar

[4] P. Kiefer, H. Vincent, S. Kalt, Q. Sun, Y. M. Eggeler, and M. Wegener, “A multi-photon (7×7)-focus 3D laser printer based on a 3D-printed diffractive optical element and a 3D-printed multi-lens array,” Light: Adva. Manuf., vol. 4, no. 1, pp. 28–41, 2024.10.37188/lam.2024.003Search in Google Scholar

[5] P. Somers, A. Münchinger, S. Maruo, C. Moser, X. Xu, and M. Wegener, “The physics of 3D printing with light,” Nat. Rev. Phys., vol. 6, no. 2, pp. 99–113, 2024.10.1038/s42254-023-00671-3Search in Google Scholar

[6] H. Gehring, et al.., “Low-loss fiber-to-chip couplers with ultrawide optical bandwidth,” APL Photonics, vol. 4, no. 1, p. 010801, 2019, https://doi.org/10.1063/1.5064401.Search in Google Scholar

[7] L. Siegle, D. Xie, C. A. Richards, P. V. Braun, and H. Giessen, “Diffractive microoptics in porous silicon oxide by grayscale lithography,” Opt. Express, vol. 32, no. 20, pp. 35678–35688, 2024, https://doi.org/10.1364/oe.538142.Search in Google Scholar

[8] P. Ruchka, et al.., “3D-printed micro-axicon enables extended depth-of-focus intravascular optical coherence tomography in vivo,” Advanced Photonics, vol. 7, no. 2, 2025.10.1117/1.AP.7.2.026003Search in Google Scholar

[9] M. D. Schmid, A. Toulouse, S. Thiele, M. Simon, A. M. Herkommer, and H. Giessen, “3D direct laser writing of highly absorptive photoresist for miniature optical apertures,” Adv. Funct. Mater., vol. 33, no. 39, p. 2211159, 2023, https://doi.org/10.1002/adfm.202211159.Search in Google Scholar

[10] S. Papamakarios, I Katsantonis, M Kafesaki, K. L. Tsakmakidis, and M Farsari, “Cross-shaped resonators for BICs in mid-IR fabricated using multiphoton lithography,” in 2024 Eighteenth International Congress on Artificial Materials for Novel Wave Phenomena (Metamaterials), Crete, Greece, IEEE, 2024, pp. 1–3.10.1109/Metamaterials62190.2024.10703218Search in Google Scholar

[11] O Tsilipakos, G Perrakis, M Farsari, and M Kafesaki, “Polymeric optical metasurfaces by two-photon lithography: practical designs for beam steering,” in 2024 Eighteenth International Congress on Artificial Materials for Novel Wave Phenomena (Metamaterials), Crete, Greece, IEEE, 2024, pp. 1–3.10.1109/Metamaterials62190.2024.10703264Search in Google Scholar

[12] S. Papamakarios, et al.., “Cactus-like metamaterial structures for electromagnetically induced transparency at thz frequencies,” ACS Photonics, vol. 12, no. 1, 2024. https://doi.org/10.1021/acsphotonics.4c01179.Search in Google Scholar PubMed PubMed Central

[13] D. Pereira, M. S. Ferreira, M. Zeisberger, and M. A. Schmidt, “Spatially controlled phase modulation - selective higher-order mode excitation in 3D nanoprinted on-Chip hollow-core waveguides,” ACS Photonics, vol. 11, no. 8, pp. 3178–3186, 2024, https://doi.org/10.1021/acsphotonics.4c00524.Search in Google Scholar

[14] J. Bauer, C. Crook, and T. Baldacchini, “A sinterless, low-temperature route to 3D print nanoscale optical-grade glass,” Science, vol. 380, no. 6648, pp. 960–966, 2023, https://doi.org/10.1126/science.abq3037.Search in Google Scholar PubMed

[15] Y. Dana, et al.., “Free-standing microscale photonic lantern spatial mode (De-) multiplexer fabricated using 3D nanoprinting,” Light: Sci. Appl., vol. 13, no. 1, p. 126, 2024, https://doi.org/10.1038/s41377-024-01466-6.Search in Google Scholar PubMed PubMed Central

[16] S.-F. Liu, Z.-W. Hou, L. Lin, Z. Li, and H.-B. Sun, “3D laser nanoprinting of functional materials,” Adv. Funct. Mater., vol. 33, no. 39, p. 2211280, 2023, https://doi.org/10.1002/adfm.202211280.Search in Google Scholar

[17] T. J. Sturges, et al.., “Inverse-designed 3D laser nanoprinted phase masks to extend the depth of field of imaging systems,” ACS Photonics, vol. 11, no. 9, pp. 3765–3773, 2024, https://doi.org/10.1021/acsphotonics.4c00953.Search in Google Scholar

[18] C. Su, et al.., “Ring-in-ring beam shaping based on a 3D nanoprinted microlens on fiber tip,” in 2024 14th International Symposium on Communication Systems, Networks and Digital Signal Processing (CSNDSP), Rome, Italy, IEEE/IET, 2024, pp. 351–353.10.1109/CSNDSP60683.2024.10636389Search in Google Scholar

[19] P. Mainik, C. A. Spiegel, and E. Blasco, “Recent advances in multi-photon 3D laser printing: active materials and applications,” Adv. Mater., vol. 36, no. 11, p. 2310100, 2024, https://doi.org/10.1002/adma.202310100.Search in Google Scholar PubMed

[20] J. Bürger, J. Kim, T. Weiss, S. A. Maier, and M. A. Schmidt, “On-chip twisted hollow-core light cages: enhancing planar photonics with 3D nanoprinting,” Arxiv Preprint, 2024.Search in Google Scholar

[21] H. Wang, et al.., “Two-photon polymerization lithography for imaging optics,” Int. J. Extreme Manuf., vol. 6, no. 4, p. 042002, 2024, https://doi.org/10.1088/2631-7990/ad35fe.Search in Google Scholar

[22] G. C. Marques, et al.., “Fully printed electrolyte-gated transistor formed in a 3D polymer reservoir with laser printed drain/source electrodes,” Adv. Mater. Technol., vol. 8, no. 22, p. 2300893, 2023, https://doi.org/10.1002/admt.202300893.Search in Google Scholar

[23] J. Weinacker, et al.., “Multi-photon 3d laser micro-printed plastic scintillators for applications in low-energy particle physics,” Adv. Funct. Mater., vol. 35, no. 2, p. 2413215, 2024. https://doi.org/10.1002/adfm.202413215.Search in Google Scholar

[24] R. E. Christiansen and O. Sigmund, “Inverse design in photonics by topology optimization: tutorial,” J. Opt. Soc. Am. B, vol. 38, no. 2, p. 496, 2021, https://doi.org/10.1364/josab.406048.Search in Google Scholar

[25] S. Molesky, Z. Lin, A. Y. Piggott, W. Jin, J. Vucković, and A. W. Rodriguez, “Inverse design in nanophotonics,” Nat. Photonics, vol. 12, no. 11, pp. 659–670, 2018, https://doi.org/10.1038/s41566-018-0246-9.Search in Google Scholar

[26] J. Noh, T. Badloe, C. Lee, J. Yun, S. So, and J. Rho, “Inverse design meets nanophotonics: from computational optimization to artificial neural network,” Intel. Nanotechnol., pp. 3–32, 2023, https://doi.org/10.1016/b978-0-323-85796-3.00001-9.Search in Google Scholar

[27] S. So, T. Badloe, J. Noh, J. Bravo-Abad, and J. Rho, “Deep learning enabled inverse design in nanophotonics,” Nanophotonics, vol. 9, no. 5, pp. 1041–1057, 2020, https://doi.org/10.1515/nanoph-2019-0474.Search in Google Scholar

[28] P. R Wiecha, A. Arbouet, C. Girard, and O. L Muskens, “Deep learning in nano-photonics: inverse design and beyond,” Photonics Res., vol. 9, no. 5, pp. B182–B200, 2021, https://doi.org/10.1364/prj.415960.Search in Google Scholar

[29] T. W. Radford, P. R. Wiecha, A. Politi, I. Zeimpekis, and O. L. Muskens, “Inverse design of unitary transmission matrices in silicon photonic coupled waveguide arrays using a neural adjoint model,” ACS Photonics, vol. 12, no. 3, 2025.10.1021/acsphotonics.4c02081Search in Google Scholar PubMed PubMed Central

[30] D. Peng, K. Sun, O. L. Muskens, C. H. de Groot, and R. Huang, “Inverse design of a vanadium dioxide based dynamic structural color via conditional generative adversarial networks,” Opt. Mater. Express, vol. 12, no. 10, pp. 3970–3981, 2022, https://doi.org/10.1364/ome.467967.Search in Google Scholar

[31] D. Lee, W. Chen, L. Wang, Y.-C. Chan, and W. Chen, “Data-driven design for metamaterials and multiscale systems: a review,” Adv. Mater., vol. 36, no. 8, p. 2305254, 2024, https://doi.org/10.1002/adma.202305254.Search in Google Scholar PubMed

[32] I. Tanriover, D. Lee, W. Chen, and K. Aydin, “Deep generative modeling and inverse design of manufacturable free-form dielectric metasurfaces,” ACS Photonics, vol. 10, no. 4, pp. 875–883, 2023.10.1021/acsphotonics.2c01006Search in Google Scholar

[33] J. S. Jensen and O. Sigmund, “Topology optimization for nano-photonics,” Laser Photonics Rev., vol. 5, no. 2, pp. 308–321, 2011, https://doi.org/10.1002/lpor.201000014.Search in Google Scholar

[34] M. P. Bendsøe and O. Sigmund, “Bibliographical notes,” in Topology Optimization: Theory, Methods, and Applications, M. P. Bendsøe, and O. Sigmund, Eds., Berlin, Heidelberg, Springer, 2004, pp. 305–318.10.1007/978-3-662-05086-6_6Search in Google Scholar

[35] N. V Sapra, et al.., “Inverse design and demonstration of broadband grating couplers,” IEEE J. Sel. Top. Quantum Electron., vol. 25, no. 3, pp. 1–7, 2019, https://doi.org/10.1109/jstqe.2019.2891402.Search in Google Scholar

[36] N. Murai, Y. Noguchi, K. Matsushima, and T. Yamada, “Multiscale topology optimization of electromagnetic metamaterials using a high-contrast homogenization method,” Compu. Methods Appl. Mech. Eng., vol. 403, p. 115728, 2023, https://doi.org/10.1016/j.cma.2022.115728.Search in Google Scholar

[37] 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, vol. 9, no. 6, pp. 374–377, 2015, https://doi.org/10.1038/nphoton.2015.69.Search in Google Scholar

[38] A. Nanda, M. Kues, and A. C. Lesina, “Exploring the fundamental limits of integrated beam splitters with arbitrary phase via topology optimization,” Opt. Lett., vol. 49, no. 5, pp. 1125–1128, 2024, https://doi.org/10.1364/ol.512100.Search in Google Scholar

[39] E. Hassan and A. C. Lesina, “Topology optimization of dispersive plasmonic nanostructures in the time-domain,” Opt. Express, vol. 30, no. 11, pp. 19557–19572, 2022, https://doi.org/10.1364/oe.458080.Search in Google Scholar PubMed

[40] J. Gedeon, I. Allayarov, A. C. Lesina, and E. Hassan, “Time-domain topology optimization of power dissipation in dispersive dielectric and plasmonic nanostructures,” IEEE Trans. Antennas Propag., vol. 10, no. 11, 2025. https://doi.org/10.1109/tap.2024.3517156.Search in Google Scholar

[41] V. Igoshin, A. Kokhanovskiy, and M. Petrov, “Inverse design of mie resonators with minimal backscattering,” Optics Letters, vol. 50, no. 5, 2025.10.1364/OL.552002Search in Google Scholar PubMed

[42] V. D. Igoshin, A. Y. Kokhanovskiy, and M. I. Petrov, “Inverse design of stacked cylinders scatterers by evolutionary optimization,” in 2024 Photonics & Electromagnetics Research Symposium (PIERS), Chengdu, China, IEEE, 2024, pp. 1–5.10.1109/PIERS62282.2024.10618079Search in Google Scholar

[43] W. Jin, S. Molesky, Z. Lin, K.-M. C Fu, and A. W Rodriguez, “Inverse design of compact multimode cavity couplers,” Optics express, vol. 26, no. 20, pp. 26713–26721, 2018, https://doi.org/10.1364/oe.26.026713.Search in Google Scholar

[44] Y. Augenstein and C. Rockstuhl, “Inverse design of nanophotonic devices with structural integrity,” ACS Photonics, vol. 7, no. 8, pp. 2190–2196, 2020, https://doi.org/10.1021/acsphotonics.0c00699.Search in Google Scholar

[45] B. M. de Aguirre Jokisch, R. E. Christiansen, and O. Sigmund, “Engineering optical forces through Maxwell stress tensor inverse design,” Journal of the Optical Society of America B, vol. 42, no. 4, 2025.10.1364/JOSAB.546272Search in Google Scholar

[46] B. M. de Aguirre Jokisch, R. E. Christiansen, and O. Sigmund, “Topology optimization framework for designing efficient thermo-optical phase shifters,” J. Opt. Soc. Am. B, vol. 41, no. 2, pp. A18–A31, 2024, https://doi.org/10.1364/josab.499979.Search in Google Scholar

[47] R. E Christiansen, “Inverse design of optical mode converters by topology optimization: tutorial,” J. Opt., vol. 25, no. 8, p. 083501, 2023, https://doi.org/10.1088/2040-8986/acdbdd.Search in Google Scholar

[48] Q. Li, W. Chen, S. Liu, and L. Tong, “Structural topology optimization considering connectivity constraint,” Struct. Multidiscip. Optim., vol. 54, no. 4, pp. 971–984, 2016, https://doi.org/10.1007/s00158-016-1459-5.Search in Google Scholar

[49] J. D. Fischbach, “Jan-david-fischbach/jaxwell,” 2024. Available at: https://github.com/jan-david-fischbach/jaxwell.Search in Google Scholar

[50] M. F. Schubert, A. K. C. Cheung, I. A. D. Williamson, A. Spyra, and D. H. Alexander, “Inverse design of photonic devices with strict foundry fabrication constraints,” ACS Photonics, vol. 9, no. 7, pp. 2327–2336, 2022, https://doi.org/10.1021/acsphotonics.2c00313.Search in Google Scholar

[51] C. M. Lalau-Keraly, S. Bhargava, O. D. Miller, and E. Yablonovitch, “Adjoint shape optimization applied to electromagnetic design,” Opt. Express, vol. 21, no. 18, p. 21693, 2013, https://doi.org/10.1364/oe.21.021693.Search in Google Scholar PubMed

[52] M. P. Bendsøe, “Optimal shape design as a material distribution problem,” Struct. optim., vol. 1, no. 4, pp. 193–202, 1989, https://doi.org/10.1007/bf01650949.Search in Google Scholar

[53] O. Sigmund and K. Maute, “Topology optimization approaches,” Struct. Multidiscip. Optim., vol. 48, no. 6, pp. 1031–1055, 2013, https://doi.org/10.1007/s00158-013-0978-6.Search in Google Scholar

[54] M. Zhou, B. S. Lazarov, F. Wang, and O. Sigmund, “Minimum length scale in topology optimization by geometric constraints,” Compu. Methods Appl. Mech. Eng., vol. 293, pp. 266–282, 2015, https://doi.org/10.1016/j.cma.2015.05.003.Search in Google Scholar

[55] F. Wang, B. Stefanov Lazarov, and O. Sigmund, “On projection methods, convergence and robust formulations in topology optimization,” Struct. Multidiscip. Optim., vol. 43, no. 6, pp. 767–784, 2011, https://doi.org/10.1007/s00158-010-0602-y.Search in Google Scholar

[56] M. Paszkiewicz-Idzik, M. Sukhova, W. Dörfler, and C. Rockstuhl, “Optimizing the trajectories of freeform waveguides using a modal theory,” in Computational Optics, vol. PC13023, D. G. Smith and A. Erdmann, Eds., San Diego, USA, SPIE/International Society for Optics and Photonics, 2024, p. PC1302305.10.1117/12.3014421Search in Google Scholar

[57] S. Irfan, J.-Y. Kim, and H. Kurt, “Ultra-compact and efficient photonic waveguide bends with different configurations designed by topology optimization,” Sci. Rep., vol. 14, no. 1, p. 6453, 2024, https://doi.org/10.1038/s41598-024-53881-9.Search in Google Scholar PubMed PubMed Central