Source-linked AI summary
Inverse design of photonic crystals through automatic differentiation
Momchil Minkov, Ian A. D. Williamson, Lucio C. Andreani, Dario Gerace, Beicheng Lou, Alex Y. Song, Tyler W. Hughes, Shanhui Fan
TL;DR
Mode-expansion methods widely used for photonic crystals lacked a straightforward adjoint formulation for gradient-based optimization. This paper uses automatic differentiation to implement plane-wave and guided-mode expansion, then demonstrates optimization of waveguide dispersion and an ultra-small lithium-niobate cavity. The approach extends gradient-based inverse design to these simulations while supporting flexible objectives and parameterizations.
Problem
Mode-expansion methods used for periodic optical structures had not yet been incorporated into gradient-based optimization because their adjoint variable method could not be defined explicitly as in standard field solvers.
Method
The paper implements differentiable plane-wave and guided-mode expansion methods using automatic differentiation, which generates gradients through arbitrary computational dependencies.
Results
The implementation computes gradients efficiently with respect to all input parameters and optimizes photonic-crystal waveguide dispersion and an ultra-small lithium-niobate cavity’s quality factor.
Takeaways & Limitations
Automatic differentiation extends gradient-based inverse design to widely used photonic-crystal mode-expansion simulations and enables flexible structure parameterizations and objective functions.
Abstract
from arXiv · showhide
Gradient-based inverse design in photonics has already achieved remarkable results in designing small-footprint, high-performance optical devices. The adjoint variable method, which allows for the efficient computation of gradients, has played a major role in this success. However, gradient-based optimization has not yet been applied to the mode-expansion methods that are the most common approach to studying periodic optical structures like photonic crystals. This is because, in such simulations, the adjoint variable method cannot be defined as explicitly as in standard finite-difference or finite-element time- or frequency-domain methods. Here, we overcome this through the use of automatic differentiation, which is a generalization of the adjoint variable method to arbitrary computational graphs. We implement the plane-wave expansion and the guided-mode expansion methods using an automatic differentiation library, and show that the gradient of any simulation output can be computed efficiently and in parallel with respect to all input parameters. We then use this implementation to optimize the dispersion of a photonic crystal waveguide, and the quality factor of an ultra-small cavity in a lithium niobate slab. This extends photonic inverse design to a whole new class of simulations, and more broadly highlights the importance that automatic differentiation could play in the future for tracking and optimizing complicated physical models.
I. INTRODUCTION
Photonic inverse design uses gradients to optimize complex optical devices, but mode-expansion simulations make explicit adjoint implementations difficult. The paper applies automatic differentiation to these computational graphs and demonstrates gradient-based optimization with differentiable plane-wave and guided-mode expansion methods.
- Photonic inverse design exploits wavelength-scale device degrees of freedom to automate compact, efficient optical-device design.
- Gradient-based optimization computes objective gradients with respect to all tunable parameters and updates device parameters along the gradient direction.The adjoint variable method efficiently computes these gradients for many parameters when simulations solve a linear system.
- Eigenvalue-based photonic-crystal and mode-expansion simulations create complex parameter dependencies that make explicit backward-simulation implementations difficult or infeasible.The system matrix itself can depend non-trivially on structural parameters, while objectives may depend on eigenvalues and eigenmodes.
- Automatic differentiation applies the adjoint variable method to arbitrary computational graphs, generating backward computations from a defined forward computation.Gradient-aware elementary functions can be flexibly composed to provide end-to-end gradient support for complex computations.
- The paper implements differentiable 2D plane-wave expansion and guided-mode expansion methods, computing gradients of simulation outputs with respect to all structural input parameters.The implementation supports objectives based on dispersion, field profiles, or loss rates without explicitly defining adjoint fields.
- Two demonstrations optimize photonic-crystal waveguide dispersion and the quality factor of an ultra-small lithium-niobate-slab cavity.The waveguide targets dispersion curves, while the cavity optimization addresses a small-volume resonator.
B. Mode expansion methods
Mode-expansion methods represent Maxwell eigenstates in the eigenbasis of a simpler operator, producing a matrix eigenproblem for expansion coefficients. Their accuracy depends on basis completeness and truncation, with guided-mode expansion remaining approximate but computationally faster.
- The formulation assumes linear, isotropic, lossless, nondispersive, and non-magnetic media, with periodic or spatially decaying fields for the inner product.For open systems, correctly defining the inner product is beyond the scope of the work.
- Mode-expansion methods express eigenstates of an operator in the basis of eigenstates from a different, simpler operator.
- The basis modes solve Maxwell’s equations for a typically simple starting structure, such as free space, and form an orthonormal set for the Hermitian problem.
- The expansion coefficients satisfy a matrix eigenproblem with elements Hνµ = (Hν, ˆΘHµ).
- Mode-expansion accuracy requires a complete basis and all modes, but practical calculations truncate the basis to a finite subset.
- Guided-mode expansion restricts the basis to fully guided slab modes, making it approximate while retaining substantially faster computation than first-principles simulations.The plane-wave basis is complete apart from finite truncation, whereas the guided-mode basis is incomplete.
III. PLANE-WAVE EXPANSION
The paper reviews the plane-wave expansion method as a standard approach for simulating periodic structures and as a foundation for understanding guided-mode expansion and its automatic-differentiation implementation.
- The plane-wave expansion method is a well-known approach for simulating periodic structures in photonics and quantum mechanics.The review supports understanding both guided-mode expansion and the automatic-differentiation implementation of the two methods.
- The paper reviews plane-wave expansion fundamentals to clarify its relationship to guided-mode expansion and differentiable implementation.
A. Method description
The 2D plane-wave expansion models periodic photonic crystals by expanding their modes in reciprocal-lattice plane waves and solving the resulting eigenvalue problem. For 2D structures, the formulation uses Bloch periodicity, TE/TM polarization separation, Fourier-permittivity matrix elements, and a truncated reciprocal-space basis.
- 2D plane-wave expansion: The 2D PWE formulation assumes periodicity in the xy-plane and translational invariance in z, with modes represented using Bloch momentum k.The periodic part of each mode has the 2D lattice periodicity; the discussion focuses on kz = 0, although nonzero kz is straightforward to generalize.
- 2D plane-wave expansion: Reflection symmetry separates the modes into transverse electric and transverse magnetic polarizations, even for arbitrary in-plane permittivity.The plane-wave basis is obtained from the free-space system with relative permittivity equal to one.
- Reciprocal-space representation: The plane-wave expansion uses reciprocal lattice vectors, while Fourier components of the inverse permittivity provide the matrix elements coupling basis states.Only reciprocal-lattice-compatible plane-wave combinations contribute because non-commensurate Fourier components vanish.
- Numerical implementation: The reciprocal-space basis is truncated at |G| ≤ Gmax to make the eigenproblem numerical.For discontinuous permittivity, convergence improves when the permittivity Fourier transform is computed first and the inverse-permittivity matrix is obtained by matrix inversion.
- Numerical implementation: The computational graph for the full PWE consists of three main operations required to compute eigenmodes of an arbitrary 2D photonic crystal.The figure presents this graph alongside the reciprocal-lattice construction and the primitive cell.
B. Automatic differentiation
Automatic differentiation extends adjoint-style gradient propagation to the computational graph of mode-expansion simulations. The implementation differentiates eigensolves, matrix inverses, Fourier transforms, and shape-based parameterizations, while eigenvector gradients remain approximate when only a subset of modes is computed.
- Differentiating the computational graph: Reverse-mode AD propagates derivatives backward through the Hermitian eigensolve, which maps a Hermitian matrix to eigenvalues and associated eigenvectors.The eigensolve required extending Autograd because its existing eigh operation supported only real symmetric matrices.
- Differentiating the computational graph: The eigensolve derivative assumes nondegenerate eigenvalues, although degenerate cases can be handled with additional perturbation-theory considerations.Eigenvalues may be labeled, for example, by sorting them according to magnitude.
- Eigensolve derivatives: Eigenvalue derivatives depend only on their corresponding eigenvectors, whereas eigenvector derivatives depend on all other eigenvectors.This distinction determines whether a truncated eigensolve yields an exact or approximate backward pass.
- Eigensolve derivatives: When an objective depends on eigenvectors, computing only a subset of modes gives exact eigenvalue gradients but only approximate eigenvector gradients.The approximation becomes exact when all eigenvectors are computed, and modes with widely separated eigenvalues contribute less through the eigenvalue-denominator dependence.
- Other graph operations: The matrix-inverse derivative uses dA^-1/dp = -A^-1(dA/dp)A^-1, while Fourier-transform derivatives can be propagated through a discrete Fourier transform.The implementation instead parameterizes permittivity with circles and polygons, whose Fourier transforms are differentiable through algebraic operations and an analytic Bessel-function derivative.
C. Waveguide optimization
The PWE and AD framework is applied to engineer photonic-crystal-waveguide dispersion toward progressively richer target curves. The first two targets are matched closely, while the highest-spatial-frequency target is harder to reach with the chosen parameterization and supercell.
- Target formulation: The waveguide is a missing row of holes in a hexagonal dielectric-rod lattice, and optimization targets three prescribed dispersion curves.The curves are generic, while dispersion engineering is motivated by nonlinear phase matching and frequency-comb generation.
- Target formulation: The objective is the mean-square error between the optimized reduced-frequency band and each target, allowing an arbitrary frequency offset.The optimized parameters are positions and radii of rods surrounding the waveguide, with 45 free parameters under mirror symmetry.
- Optimization results: The first optimization converges quickly to the target dispersion with a final MSE below 10^-6, while gradients for all 45 parameters take approximately one forward-simulation time.LBFGS is used to optimize the waveguide structure.
- Optimization results: The second, more complicated dispersion pattern requires more optimization iterations but still matches its target very well.The corresponding panels show the optimized band, structure and field intensity, and objective evolution.
- Optimization results: The third target is not matched perfectly: initialization from the second optimized structure improves the result, but higher-k components remain difficult to reproduce.The authors attribute this to smaller-scale real-space permittivity fluctuations and suggest arbitrary rod shapes or a larger supercell as possible improvements.
IV. GUIDED-MODE EXPANSION
Guided-mode expansion efficiently simulates many three-dimensional layered photonic-crystal slabs using a basis that analytically captures z-dependence. The paper extends the AD framework to GME and applies it to quality-factor optimization in an ultra-small lithium-niobate cavity.
- Method and application: GME simulates three-dimensional layered periodic photonic-crystal slabs with computational complexity similar in many cases to 2D PWE.Its basis modes capture the z-dependence analytically, while in-plane periodic permittivity enters through Fourier components.
- Method and application: GME is approximate because its basis modes are incomplete, but it agrees well with first-principles methods for modes well confined to the photonic-crystal region.The paper discusses AD implementation for GME and demonstrates quality-factor optimization for a lithium-niobate slab cavity with ultra-small volume.
A. Method description
The guided-mode expansion models layered periodic photonic structures using guided modes of an effective multilayer basis, then computes photonic bands and radiative losses.
- Layered structure: The method represents a layered photonic crystal with semi-infinite claddings and N finite layers sharing a common in-plane periodicity.Each layer has an in-plane permittivity distribution and thickness d_j.
- Basis construction: The expansion basis is derived from a homogeneous multilayer whose layer permittivities are elementary-cell averages of the original structure.Fourier components of the actual layer permittivities enter the subsequent matrix construction.
- Guided modes: Guided modes are found by enforcing evanescence in both claddings and solving D22(ω)=0 below the higher-index cladding light line.The guided-mode frequency and associated coefficients are then computed numerically.
- Eigenmode expansion: The photonic-crystal eigenmode is expanded in fully guided modes indexed by reciprocal vector, band, and polarization, with TE and TM components generally mixed.The resulting matrix is normalized and diagonalized to obtain Bloch-band frequencies and eigenvectors.
- Radiative loss: Radiative coupling is added through perturbative overlaps with outgoing modes, yielding the imaginary frequency component and quality factor Q = R(ω_k)/(2I(ω_k)).The homogeneous basis contains both discrete guided bands and radiative continua above the relevant light lines.
B. Automatic differentiation
Automatic differentiation extends the guided-mode expansion computational graph through its implicit guided-mode solve and lengthy algebraic operations, enabling efficient gradients of the full simulation.
- Computational graph: The guided-mode expansion is differentiated as a computational graph, with Fourier transforms and matrix inversions performed separately for every layer.The nontrivial added operation is differentiation through the nonlinear guided-mode solve.
- Implicit differentiation: Guided-mode frequencies depend implicitly on average permittivities and layer thicknesses, so reverse-mode differentiation propagates through the scalar equation D22(ω)=0.Autograd computes the required derivatives of D22 with respect to frequency, permittivity, and thickness.
- Gradient accumulation: Automatic differentiation handles the many algebraic operations needed for diagonalization and radiative overlaps without manually deriving their long derivative expressions.Reverse-mode differentiation also provides a speed advantage when many inputs contribute to a scalar objective.
C. Cavity optimization
The authors optimize a lithium-niobate L4/3 photonic-crystal cavity by moving hole positions, increasing its quality factor while reducing mode volume relative to the starting and prior designs.
- Motivation: Lithium niobate cavity design is challenging because its lower refractive index produces weaker index contrast than silicon near λ = 1.55µm.A prior design had theoretical Q = 1.5×106 and mode volume V = 2.43(λ/n)^3.
- Optimization setup: The optimized cavity is based on an L4/3 defect in a triangular lattice, with 70 hole-position parameters and quality factor Q as the objective.The underlying structure uses a = 620nm, slab thickness d = 270nm, and hole radius r = 145nm.
- Gradient validation: The automatic-differentiation gradients match finite-difference values, while requiring one reverse-mode computation rather than 70 extra simulations.The reverse-mode computation after the initial quality-factor evaluation took approximately the same time as the forward simulation.
- Optimization result: After about 100 LBFGS epochs, the optimized cavity reaches a GME-computed quality factor of Q = 9.7 × 106.The optimization is described as close to converged after 100 epochs.
- Validation and compactness: The final cavity has mode volume V = 0.49(λ/n)^3 and finite-difference validation gives Q = 2.4 × 106.The lower finite-difference value is attributed to the approximate nature of the guided-mode expansion.
V. DISCUSSION AND CONCLUSION
The paper extends gradient-based inverse design to photonic-crystal mode-expansion methods using differentiable implementations. It demonstrates practical optimization benefits and identifies applications enabled by the guided-mode expansion.
- V. DISCUSSION AND CONCLUSION: The guided-mode expansion improved the quality factor of a small-volume photonic-crystal cavity by more than two orders of magnitude.The method is approximate but has agreed well with first-principles simulations while remaining computationally faster.
- V. DISCUSSION AND CONCLUSION: Guided-mode expansion is useful for quasi-guided modes above the light line, which are difficult to isolate with first-principles finite-difference methods.The method is also described as extremely fast: a converged 16x10-cell simulation took about ten minutes on a personal computer.
- V. DISCUSSION AND CONCLUSION: The supplementary discussion identifies bound states in the continuum and fully 3D photonic-crystal slab waveguide dispersion as potential applications.The authors also note that optimization can include radiation losses and is not restricted to modes below the light line.
- V. DISCUSSION AND CONCLUSION: The work implements 2D plane-wave expansion and guided-mode expansion with reverse-mode automatic differentiation for efficient gradient computation.The implementation supports flexible structure parametrizations and objective functions, and the software was made publicly available.
SUPPLEMENTARY INFORMATION
The supplementary information specifies normalization procedures for guided and radiative modes in the guided-mode expansion. Radiative normalization uses outgoing cladding components and handles semi-radiative modes according to which cladding supports propagation.
- S2. RADIATIVE MODE NORMALIZATION: Guided-mode normalization is imposed by rescaling the modal coefficients by the square root of the normalization terms for TE and TM polarizations.The supplied text states this after the TE and TM normalization expressions.
- S2. RADIATIVE MODE NORMALIZATION: Radiative-mode normalization uses a box extending from −L to L, with the normalization integral modified accordingly.As L increases, contributions from inner layers become negligible and the cladding layers dominate.
- S2. RADIATIVE MODE NORMALIZATION: Only outgoing radiative components contribute to normalization because they determine the energy flux far from the waveguide.The outgoing-wave boundary condition is imposed by adding a small negative imaginary part to the squared frequency, which is then set to zero.
- S2. RADIATIVE MODE NORMALIZATION: If the 1/L factor is omitted, intermediate normalization and matrix-element quantities acquire compensating L-dependence while the final intrinsic losses remain unchanged.The supplied text states that expressions, matrix elements, and the density of states change with L but the final losses do not.
- S2. RADIATIVE MODE NORMALIZATION: For semi-radiative modes, only the coefficient associated with the propagating cladding component appears in the normalization.This occurs when unequal cladding dielectric constants place the mode between the two cladding light lines.
S3. MATRIX ELEMENTS
The matrix-element derivation transforms the scalar product into an evaluable form using the paper’s Fourier and modal field expressions. The resulting expressions apply to both guided-guided and guided-radiative overlaps.
- S3. MATRIX ELEMENTS: The scalar product for the matrix elements is transformed and evaluated using Eq. (12) together with the modal field expressions in Eqs. (36) and (39).These expressions are used for the various polarization combinations.
- S3. MATRIX ELEMENTS: The resulting matrix-element equations apply to overlaps between two guided modes and between a guided mode and a radiative mode.The supplied text explicitly states this scope for Eqs. (S8-S10).