High-efficiency, large-area, topology-optimized
metasurfaces

Thaibao Phan, David Sell, Evan Wang, Sage Doshay, Kofi Sényo Edee, Jianji

Yang, Jonathan Fan

► To cite this version:

Thaibao Phan, David Sell, Evan Wang, Sage Doshay, Kofi Sényo Edee, et al. High-efficiency, large-area, topology-optimized metasurfaces. Light: Science and Applications, 2019, 8 (1), 10.1038/s41377-019-0159-5. [https://doi.org/10.1038/s41377-019-0159-5](https://doi.org/10.1038/s41377-019-0159-5)

HAL Id: hal-03049167

https://uca.hal.science/hal-03049167v1

Submitted on 9 Dec 2024

**HAL** is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# High-efficiency, large-area, topology-optimized metasurfaces

Thaibao Phan<sup>1</sup>, David Sell<sup>2</sup>, Evan W. Wang<sup>1</sup>, Sage Doshay<sup>2</sup>, Kofi Edee<sup>3,4</sup>, Jianji Yang<sup>1</sup> and Jonathan A. Fan<sup>1</sup>

## Abstract

Metasurfaces are ultrathin optical elements that are highly promising for constructing lightweight and compact optical systems. For their practical implementation, it is imperative to maximize the metasurface efficiency. Topology optimization provides a pathway for pushing the limits of metasurface efficiency; however, topology optimization methods have been limited to the design of microscale devices due to the extensive computational resources that are required. We introduce a new strategy for optimizing large-area metasurfaces in a computationally efficient manner. By stitching together individually optimized sections of the metasurface, we can reduce the computational complexity of the optimization from high-polynomial to linear. As a proof of concept, we design and experimentally demonstrate large-area, high-numerical-aperture silicon metasurface lenses with focusing efficiencies exceeding 90%. These concepts can be generalized to the design of multifunctional, broadband diffractive optical devices and will enable the implementation of large-area, high-performance metasurfaces in practical optical systems.

## Introduction

Metasurfaces are optical devices that utilize subwavelength-scale structuring to shape and manipulate electromagnetic waves<fnref n="1" />. They are powerful complements to bulk refractive and scalar diffractive optics and have a wide range of potential applications in imaging<fnref n="2" />, lithography³, sensing⁴, and computing platforms⁵. To date, a broad range of metasurface design concepts have emerged. The most widely used methods, which we will term “conventional methods,” sample the desired phase profile using discrete phase-shifter elements to form a nanoscale phased array (Fig. 1a). These methods utilize a library of simple, physically intuitive building blocks, including anisotropic waveguides⁶, Mie resonators⁷, plasmonic resonators⁸, and dielectric transmit arrays⁹, and can quickly produce macroscale device designs. However, these approaches lack the necessary degrees of freedom for realizing high-efficiency in devices that are

designed for large-angle deflections, multiple functions, and broadband responses¹⁰, thereby preventing meta-surfaces from being practically applied in many contexts.

Topology optimization is an alternative design method that can overcome the efficiency limitations of conventional approaches¹¹. Inverse design methods, including objective-first and adjoint-based topology optimization¹², have produced photonic crystals¹³,¹⁴, optical demultiplexers¹⁵, and spectral splitters¹⁶ with outstanding performance resulting from the devices’ unusual geometries and nonintuitive optical dynamics. More recently, topology optimization has led to high-performance meta-surfaces that have a broad range of capabilities, such as high-efficiency light deflection with nearly arbitrary input and output angles¹⁷, spectral sorting of plane waves to distinct diffraction orders¹⁸, and light focusing with field-curvature correction¹⁹. While devices that are designed using topology optimization support enhanced efficiencies, they require extensive computational resources that dramatically increase with the size of the device. Hence, current topology-optimized devices are either microscale in dimension or limited to periodic structures with microscale unit cells.

<figure><img src="image_2.png" /></figure>

<footnote>Correspondence: Jonathan A. Fan (jonfan@stanford.edu)</footnote>

<footnote n="1">¹Department of Electrical Engineering, Stanford University, Stanford, CA 94305, USA</footnote>

<footnote n="2">²Department of Applied Physics, Stanford University, Stanford, CA 94305, USA
Full list of author information is available at the end of the article.</footnote>
<figure><img src="image_3.png" /></figure>

We report a conceptually new approach to metasurface design that extends the high performance of topology-optimized devices to macroscopic areas in a computationally efficient manner. The concept is schematically illustrated in Fig. 1b and is comprised of three parts. First, we discretize the desired phase profile into a series of wavelength-scale, linear sections. Next, we use topology optimization to design metasurface elements that optimally scatter light as required by each section of the phase profile. Finally, we stitch all the elements together to form a full metasurface. Our approach produces devices that operate more efficiently than conventional designs by accounting for and optimizing the near-field optical coupling between neighboring nanostructures. Our design method is also substantially more computationally efficient than existing topology optimization methods.

# Results

## Assessment of computational efficiency

To benchmark the improvements in computational efficiency that are afforded by our approach, we perform adjoint-based topology optimization on metagratings that are made of silicon ridges. These periodic metasurfaces are designed to deflect light of a specified wavelength $λ$ to the +1 diffraction order. We perform electromagnetic simulations for topology optimization using rigorous coupled wave analysis (RCWA), which is also known as the Fourier modal method$^{20–23}$, on a personal computer and plot the time that is required for optimizing a full device as a function of its width $L$ in Fig. 1c. The simulation time scales approximately as $O(L^{2.4})$, which corresponds to the general scaling trend for electromagnetic solvers that utilize standard matrix multiplication and inversion algorithms$^{24}$. These trends demonstrate that prohibitively large computational resources are required

for directly optimizing devices that are many times larger than the wavelength.

If we instead divide the metagating into sections of width $d$, we have $L/d$ sections to optimize. The total computation time now scales as $O(d^{2.4} \cdot L/d) = O(d^{1.4} \cdot L)$, which is a linear function of the total size. This expression also shows that the computation time decreases with decreasing section size. In practice, we found that using 3λ-wide sections minimizes the total computation time, as the benefits from using smaller sections are outweighed by the computational overhead in the electromagnetic solver. We plot the computation time for a metasurface that has been optimized using 3λ-wide sections in Fig. 1c. The observed reduction in the computation time enables us to optimize devices with dimensions that are many times larger than the wavelength using realistic computational resources. For example, via our method, we can produce millimeter-scale topology-optimized metasurfaces that operate at near-infrared wavelengths in less than 1 day using a personal computer. In contrast, optimizing the full device at once would take nearly 1 year and require intractably large amounts of memory.

We realize additional savings in computation time via our method when using multiple computing cores because individual segments can be optimized on different cores. Our design concept allows the segments to be treated independently, thereby enabling the optimizations to be parallelized without issues concerning race conditions and synchronization$^{25}$. With $N$ computing nodes, metasurfaces that are subdivided into $N$ sections can be optimized in the same amount of time as it takes to optimize a single segment.

These reductions in computational complexity also apply to fully three-dimensional topology-optimized metasurfaces, which exhibit even more severe scaling
<figure><img src="image_4.png" /><figcaption>Fig. 2 Impact of phase profile linearization on the cylindrical lens performance. a Cylindrical lenses of focal length 36λ and NA of 0.7 are constructed using linear sections of length d for various values of d. b Line scans of the field intensity at the focal planes of the lenses. Lenses that are linearized with sections that are smaller than 4λ have peak intensities that are within 1% of that of the ideal lens</figcaption></figure>

trends: the time it takes to optimize a metasurface of size L × L all at once scales approximately as O(L$^{4.8}$). This trend indicates the necessity of linear sectioning for these more intricate design problems, which would reduce the runtime to a more reasonable O(L$^{2}$).

## Sectioning of a curvilinear phase profile

The process of approximating a curvilinear phase profile with a series of linear sections introduces wavefront error, which is denoted as $\epsilon$ (Fig. 1b, inset). However, this error has a negligible impact on the overall metasurface performance if the sections are sufficiently small. To analyze the effect of linearizing a general curvilinear phase profile $\phi(x)$, we locally describe each section of $\phi(x)$ at location $x_0$ using a 2nd-order Taylor series expansion as follows:

$$ \phi(x) \approx \phi(x_0) + \phi'(x_0)(x - x_0) + \frac{1}{2}\phi''(x_0)(x - x_0)^2 \quad (1) $$

This section can be approximated as a line of slope $\phi'(x_0)$ with a phase offset of $\phi(x_0) + \Delta\phi$, which incurs an error of $\epsilon(x) = \phi''(x_0)(x - x_0)^2/2 - \Delta\phi$. Given a section of length d at position $x_0$, the root-mean-square (RMS) wavefront error is minimized when $\Delta\phi = \phi''(x_0)d^2/24$ and is expressed as follows:

$$ \epsilon_{rms}(\Delta x) = \frac{1}{12\sqrt{5}}\phi''(x_0)d^2 \quad (2) $$

We use this result to analyze the impact of linearization on a focusing cylindrical lens, which enables us to benchmark the device performance using well-established metrics in lens design. To quantify the performance of lenses that are constructed with linear sections, we use the Strehl ratio$^{26}$, which is a metric that compares the diffraction efficiency of our lenses to that of an ideal lens. For a lens that focuses normally incident light, the ideal phase

profile is $\phi(x) = (2\pi/\lambda)(f - \sqrt{f^2 + x^2})^{27}$, where f is the focal length and $\lambda$ is the wavelength. A lens with a Strehl ratio of 0.98, which corresponds to an RMS wavefront error of $\lambda/50$, can be realized if we use linear sections that are no larger than the following:

$$ d < 0.73\sqrt{f\lambda} \quad (3) $$

This equation provides a practical and quantitative guide for linearizing phase profiles in a manner that minimizes the phase error. As an example, consider a cylindrical lens with a focal length of $f=36\lambda$ and a numerical aperture (NA) of 0.7 (Fig. 2a). Using Eq. 3, we expect that linearizing the phase profile with segments that are smaller than 4.4λ will have a negligible impact on the performance. To verify, we simulate lenses that are linearized with various segment lengths and calculate the field intensities at the focal plane. These field intensity profiles are plotted in Fig. 2b, according to which the lenses that are linearized with section lengths of 2λ and 4λ are nearly indistinguishable from the ideal lens.

Our approach to sectioning readily extends to three-dimensional phase profiles, which can be approximated as series of planar tiles. This approximation is discussed in detail in the Supplementary Section. For a hyperboloid that is discretized into tiles with dimensions $d \times d$, an RMS wavefront error of $\lambda/50$ can be realized if:

$$ d < 0.61\sqrt{f\lambda} \quad (4) $$

## Topology optimization for finite-sized, isolated devices

To design metasurface elements that have the desired linear phase profiles, we utilize adjoint-based topology optimization$^{28}$. Adjoint-based optimization is an iterative algorithm that modifies the device's dielectric constant distribution, namely, $\epsilon(x)$, to maximize a figure of merit (FoM). Our objective is to optimize a device that scatters
normally incident electromagnetic waves in a desired direction with electric field amplitude $E_{tgt}$ and phase $\phi_{tgt}$. To compute the FoM, we run a forward simulation in which waves that are incident onto the metasurface element scatter in the desired direction with field amplitude $E_{fwd}$ and phase $\phi_{fwd}$. Near-to-far-field transformations from the forward simulations are used to evaluate $E_{fwd}$ and $\phi_{fwd}^{29}$. The FoM describes the difference between the current and desired responses and has the following form:

$$ FoM = -A_1 \left[ |E_{tgt}|^2 - |E_{fwd}|^2 \right]^2 - A_2 \left[ \arg\left(e^{i(\phi_{tgt} - \phi_{fwd})}\right) \right]^2 \quad (5) $$

The terms $A_1$ and $A_2$ are weights that balance how strongly the FoM is biased toward optimizing the amplitude and phase, respectively. To determine how $\varepsilon(x)$ should be modified to improve the FoM each iteration, we perform a pair of forward and adjoint simulations and record the electric fields in the device for each excitation condition. These fields are used to calculate $\delta FoM$, which is the gradient of the FoM with respect to the dielectric constant at each position $x$:

$$ \delta FoM = 2A_1 (|E_{tgt}|^2 - |E_{fwd}|^2) \mathcal{R}e\{E_{fwd} \cdot \delta E^*\} - 2A_2 (\phi_{tgt} - \phi_{fwd}) \frac{1}{|E_{fwd}|^2} \mathcal{I}m\{E_{fwd} \cdot \delta E^*\} \quad (6) $$

where $\delta E$ is a function of the adjoint field and represents the variations of the field in the target direction in response to variations of the refractive index within the device. A more detailed discussion of the adjoint optimization method that is applied here is provided in the Supplementary Section.

## Aperiodic Fourier modal method

To apply these concepts to the design of isolated, finite-sized device elements, we have developed an aperiodic Fourier modal method (AFMM), which is a hybrid method that combines a solver for periodic systems with perfectly matched layers (PMLs). The key challenge of implementing PMLs involves describing both the periodic incident plane wave (the input field) and the aperiodic scattered field (the output field) of the isolated device within the same formalism. To address this challenge, we introduce a hybrid method that combines a Fourier basis, Maxwell's equations in complex coordinates, and the Stratton-Chu integral formalism$^{30–34}$.

A metasurface is typically composed of a single layer of patterned material. The patterned material can be expressed as a distribution of the relative permeability, namely, $\varepsilon(x, y)$, and the permeability, namely, $\mu(x, y)$, on the $xy$-plane. Along the thickness of the device in the $z$-direction, the device cross-section is constant. In this case, it can be shown from Maxwell's equations that the transverse electric fields satisfy the following eigenvalue

equation:

$$ -\gamma^2 \begin{bmatrix} E_x \\ E_y \end{bmatrix} = \mathcal{L}_{EH} \mathcal{L}_{HE} \begin{bmatrix} E_x \\ E_y \end{bmatrix} \quad (7) $$

where $\mathcal{L}_{EH}$ and $\mathcal{L}_{HE}$ are differential operators that are based on $\varepsilon$ and $\mu$ and the electric field solution can be expressed as $\Phi(x, y)e^{-ikyz}$. The Fourier modal method can easily be used to expand the input fields, which are assumed to be periodic, into a Fourier basis:

$$ \Phi_{in}(x, y, z) = \sum_{p} A_{p} e^{-ik_{p}y} \sum_{nm} \delta_{nm,p} e^{-ik_{n}x} e^{-ik_{p}y} \quad (8) $$

Next, the PML boundary conditions are described similarly. The PMLs can be introduced via a change of coordinates $(x, y, z) \rightarrow (\tilde{x}, \tilde{y}, \tilde{z})^{30,31}$:

$$ \left\{ \begin{array}{l} \tilde{x}(x) = (\chi_x - i\eta_x)(x - x_0) + x_0 \\ \tilde{y}(y) = (\chi_y - i\eta_y)(y - y_0) + y_0 \end{array} \right\} \quad (9) $$

Here, the parameter $\chi$ controls the scaling of the PML layers, while the parameter $\eta$ controls the PML attenuation. This transformation is also useful for expressing the output scattered fields, which are computed as eigenfunctions of Maxwell's equations in complex coordinates:

$$ \Phi_s(x, y, z) = \sum_p B_p e^{-ik_p y} \sum_{nm} \Phi_{nm,p} e^{-ik_n x} e^{-ik_p y} \quad (10) $$

Finally, the Stratton-Chu integral equation computes the radiated field in all of space$^{34}$. In this study, we design devices that are comprised of nanoridges, which can be described in this analysis by a one-dimensional Fourier basis.

## Section optimization demonstration

To demonstrate adjoint-based optimization in the design of metasurface elements, we design a 2.5$\lambda$-wide element that scatters incident TM-polarized light at a 20° angle with a phase response of $\pi/2$. The dielectric distribution and the amplitude and phase scattering profiles of the element at successive iterations in the optimization process are shown in Fig. 3a. The optimization begins with a random dielectric continuum with values between air and silicon. After a few iterations, the continuum begins to strongly scatter light at the desired angle. The final metasurface section is a binary structure of silicon in air that possesses a peak scattering amplitude and phase response matching with the targeted values. The full width at half maximum of the scattering peak is consistent with that expected from light that is diffracting from a 2.5$\lambda$-wide aperture; hence, the element is performing directional scattering near its physical limits.
<figure><img src="image_5.png" /></figure>

An intensity plot of the scattered fields from the fully optimized section (Fig. 3b) shows strong near-field coupling between neighboring nanostructures; hence, optimal near-field coupling is responsible for mediating strong scattering in the desired direction. As a method of gradient descent, adjoint-based topology optimization is a local optimizer and is sensitive to the initial dielectric distribution¹⁰. To obtain high-performance elements for a desired scattering angle and phase target, we perform ten optimizations with various initial dielectric distributions and select the best result.

## Optimized metalens simulation results

Combining everything into a proof-of-concept demonstration, we stitch together optimized metasurface elements to construct cylindrical metalenses. We design these metalenses to focus TM-polarized light at a wavelength of 640 nm. To enable device operation at visible wavelengths, we use 250-nm-thick crystalline silicon, which has relatively low absorption compared to polycrystalline and amorphous silicon, but much higher index contrast than materials such as titanium dioxide³⁵,³⁶.

First, we design and simulate 64-µm-wide metalenses with NAs that range from 0.2 to 0.9. We divide the metalenses into sections that are 2 µm wide, which is below the phase error limit of Eq. 3 and near the optimal size for efficient computation (Fig. 1c). Further reductions in section size lead to degradation of the device performance. The reason can be traced to our design of each metasurface element, which is optimized in isolation with PML boundary conditions. When the elements are stitched together to produce a device, the optical fields that are guided by a single element have evanescent tails and can couple to a neighboring element in a parasitic manner. Smaller section sizes require more elements to be

stitched together to produce a desired metasurface, thereby resulting in more boundaries and more parasitic coupling. Below a section size of 2λ, the device performance begins to degrade and below a section size of 1λ, the aperiodic boundary conditions in our optimizer are no longer valid.

There are a few approaches for addressing the issue of stitching error. One is to keep the section size relatively large compared to the wavelength. Another is to perform boundary optimization on the stitched regions to eliminate the stitching error. A third approach, which we use here, is to separate silicon structures from other sections by a gap of at least 0.2λ, thereby reducing the near-field coupling between sections. To ensure a reduction in the stitching error with this scheme, we simulate stitched sections to check for spurious diffraction and redesign the sections in the event of excess error.

The efficiencies of our simulated metalenses are summarized in Fig. 4a. The absolute efficiency is defined as the amount of power that is contained in the principal lobe of the focus compared to that of an ideal lens with 100% transmission. The relative efficiency, or focusing efficiency, compares the power in the principal lobe to that of an ideal lens that transmits the same amount of power as the device³⁷. This efficiency corresponds to the efficiency of the diffraction process, as it removes the effects of absorption from the material and reflection at the metalens interface.

The efficiency plots show that the relative efficiencies are consistently high, namely, above 93%, with minimal drop-off in performance as the NA increases. This trend is unlike that of conventional metalenses, where the efficiency decreases with increasing NA because conventional designs cannot efficiently deflect light at large angles²⁷,³⁸. The absolute efficiencies of the metalenses all
<figure><img src="image_6.png" /><figcaption>**Fig. 4** Simulation results of optimized metalenses. **a** Relative and absolute efficiencies for metalenses that were designed with various numerical apertures. **b** A full-field simulation of a metalens with an NA of 0.9</figcaption></figure>

exceed 75%, with approximately 10% of the light reflected from the metalens and 10% absorbed by the silicon. Reflection losses can be reduced via the use of more intricate three-dimensional silicon nanostructures, while absorption losses can be minimized by designing silicon-based devices for longer wavelength operation²⁸. A simulated field profile of a metalens with an NA of 0.9 is shown in Fig. 4b, which demonstrates that the lens focuses strongly with minimal spurious diffraction.

## Performance of fabricated metalenses

We design, fabricate, and characterize 200-µm-wide metalenses with NAs of 0.2, 0.5, and 0.8. To prepare crystalline silicon thin films on glass, we use hydrogen silesquioxane to bond silicon-on-insulator wafers onto Pyrex wafers under high temperature and pressure³⁶,³⁹. After removing the silicon handle wafer and the buried oxide layer, we pattern and etch the devices via standard electron beam lithography and dry etching techniques. We characterize the metalenses by collimating polarized, monochromatic light from a tunable white-light laser onto the devices and imaging the light at the focal plane with a ×100 objective (NA = 0.9) and a CCD sensor.

Scanning electron microscope images of the center of a representative device are shown in Fig. 5a and show silicon nanostructures that exhibit smooth and vertical sidewalls. The metalenses all have relative efficiencies that exceed 89% and absolute efficiencies that exceed 67% for all NAs, which are within 10% of the simulated values (Fig. 5b). All the metalenses exhibit diffraction-limited performance, as shown by the theoretical and experimental intensity plots in Fig. 5c–e. The device with an NA of 0.8 can focus light at a wavelength of 640 nm to a spot with a beam waist of 340 nm. The central lobes of the foci are all much stronger than the side lobes; hence, the focusing efficiency is high.

The metalenses maintain reasonably high efficiencies for wavelengths that range from 580 to 700 nm, as shown

in Fig. 5f–h. As these lenses are not designed to be achromatic, the focal length shifts with the wavelength. The simulated shifts are shown in the Supplementary section in Fig. S3. Future work will focus on generalizing our design approach to include achromatic functionality, which can be addressed by modifying the optimizer's figure of merit to include multiple wavelengths¹⁸. The figure of merit can be specified so that each metasurface section deflects all wavelengths in the same direction and realizes the correct dispersion for ensuring constructive interference at the focus⁴⁰,⁴¹.

## Discussion

In summary, we present a computationally efficient method for designing large-area, topology-optimized metasurfaces. Our approach of approximating a desired phase profile with a series of linear segments renders the design problem computationally tractable while minimally impacting the device performance. Wavelength-scale scattering elements that are designed via topology optimization can strongly scatter light at a desired angle and phase, thereby serving as high-performance building blocks for larger metasurfaces.

The metalenses that are demonstrated here have limited functionality and other heuristic design methods are available that can produce high-performance, high-numerical-aperture metalenses⁴². What differentiates our method is that our design platform can be generalized systematically to multifunctional, broadband, and even multilayer metasurfaces with large areas and high efficiencies. To practically extend topology optimization in these directions, continued improvements in electromagnetic solver speed will be required, particularly for fully three-dimensional devices. For multilayer devices, new innovations will be needed to address the stitching error, which can compound due to interactions between neighboring elements. We envision that these design methods will enable the realization of compact and
<figure><img src="image_7.png" /></figure>

monolithic electromagnetic devices that exhibit high efficiencies and practical functionality across the electromagnetic spectrum.

## Materials and methods

### Sample fabrication

Crystalline silicon on glass wafers was prepared via an approach that is similar to that done by Sell, et al.³⁶. A Pyrex and silicon-on-insulator (SOI) wafer were bonded under high temperature and pressure using hydrogen silesquioxane. The silicon handle wafer was etched away using SF₆ plasma and the buried oxide layer was etched

away using HF solution. The wafer was diced into pieces for later use. The metasurfaces were patterned onto a piece in the AR-P-6200 resist via electron beam lithography. Aluminum oxide was evaporated onto a piece and removed using liftoff to form a hard mask. The silicon was etched using a Cl₂ and HBr plasma. The aluminum oxide was removed using a solution of HCl at 50 °C. Additional details can be found in the supplementary information.

### Measurement procedure

A collimated beam of light from a tunable white light laser is filtered with a longpass filter and polarized with a
linear polarizer. Then, the beam passes through a weak cylindrical lens to focus it onto the metasurface aperture without significantly changing the angle of incidence. The metasurface focuses the light. Next, the focal plane is imaged onto a CMOS detector using a 0.9 NA, 100x objective and a tube lens.

To evaluate the lens efficiencies, a line scan of the focal plane image is obtained. The lens relative efficiency is calculated by integrating over a small aperture around the main focal lobe and comparing it to what would be expected for an ideal lens of the same numerical aperture. The lens absolute efficiency is calculated by multiplying the relative efficiency by the measured transmission. A schematic diagram of the optical setup and additional characterization details can be found in the supplementary information.

## Simulations

We used Reticolo, a rigorous coupled-wave analysis solver, as a basis to perform optimizations²¹. The final metasurface lens designs were simulated using Lumerical FDTD. The efficiencies were calculated by determining the optical power that is contained in the main lobe of the focus and comparing it to that expected from an ideal lens.

## Acknowledgements

This work was supported by the U.S. Air Force under Award Number FA9550-18-1-0070, the Office of Naval Research under Award Number N00014-16-1-2630, and the David and Lucile Packard Foundation. DS was supported by the National Science Foundation (NSF) through an NSF Graduate Research Fellowship. SD was supported by the Department of Defense (DoD) through an NDSEG Research Fellowship. This work was performed in part at the Stanford Nanofabrication Facility (SNF) and the Stanford Nano Shared Facilities (SNSF), which are supported by the National Science Foundation as part of the National Nanotechnology Coordinated Infrastructure under award ECCS-1542152.

## Author details

¹Department of Electrical Engineering, Stanford University, Stanford, CA 94305, USA. ²Department of Applied Physics, Stanford University, Stanford, CA 94305, USA. ³Université Clermont Auvergne, Institut Pascal, BP 10448, F-63000 Clermont-Ferrand, France. ⁴CNRS, UMR 6602, Institut Pascal, F-63177 Aubière, France

## Author contributions

T.P., D.S., K.E. and J.Y. formulated the design algorithms. T.P., D.S., E.W. and S.D. fabricated the devices. D.S. performed the measurements and analyzed the data. T.P. and J.F. wrote the manuscript with input from all authors.

## Conflict of interest

The authors declare that they have no conflict of interest.

**Supplementary information** is available for this paper at https://doi.org/10.1038/s41377-019-0159-5.

Received: 9 December 2018 Revised: 30 April 2019 Accepted: 5 May 2019
Published online: 29 May 2019

## References

1. Yu, N. F. & Capasso, F. Flat optics with designer metasurfaces. *Nat. Mater.* **13**, 139–150 (2014).

2. Pahlevaninezhad, H. et al. Nano-optic endoscope for high-resolution optical coherence tomography in vivo. *Nat. Photonics* **12**, 540–547 (2018).

3. Chen, K. et al. Moiré nanosphere lithography. *ACS Nano* **9**, 6031–6040 (2015).

4. Tittl, A. et al. Imaging-based molecular barcoding with pixelated dielectric metasurfaces. *Science* **360**, 1105–1109 (2018).

5. Silva, A. et al. Performing mathematical operations with metamaterials. *Science* **343**, 160–163 (2014).

6. Mueller, J. P. B. et al. Metasurface polarization optics: independent phase control of arbitrary orthogonal states of polarization. *Phys. Rev. Lett.* **118**, 113901 (2017).

7. Iyer, P. P., Pendharkar, M. & Schuller, J. A. Electrically reconfigurable metasurfaces using heterojunction resonators. *Adv. Opt. Mater.* **4**, 1582–1588 (2016).

8. Ding, F. et al. Broadband high-efficiency half-wave plate: a supercell-based plasmonic metasurface approach. *ACS Nano* **9**, 4111–4119 (2015).

9. Arbabi, A. et al. Dielectric metasurfaces for complete control of phase and polarization with subwavelength spatial resolution and high transmission. *Nat. Nanotechnol.* **10**, 937–943 (2015).

10. Yang, J. & Fan, J. A. Topology-optimized metasurfaces: impact of initial geometric layout. *Opt. Lett.* **42**, 3161–3164 (2017).

11. Lalau-Keraly, C. M. et al. Adjoint shape optimization applied to electromagnetic design. *Opt. Express* **21**, 21693–21701 (2013).

12. Molesky, S. et al. Inverse design in nanophotonics. *Nat. Photonics* **12**, 659–670 (2018).

13. Frandsen, L. H. et al. Broadband photonic crystal waveguide 60° bend obtained utilizing topology optimization. *Opt. Express* **12**, 5916–5921 (2004).

14. Borel, P. I. et al. Topology optimization and fabrication of photonic crystal structures. *Opt. Express* **12**, 1996–2001 (2004).

15. Piggot, A. Y. et al. Inverse design and demonstration of a compact and broadband on-chip wavelength demultiplexer. *Nat. Photonics* **9**, 374–377 (2015).

16. Xiao, T. P. et al. Diffractive spectral-splitting optical element designed by adjoint-based electromagnetic optimization and fabricated by femtosecond 3D direct laser writing. *ACS Photonics* **3**, 886–894 (2016).

17. Sell, D. et al. Ultra-high-efficiency anomalous refraction with dielectric metasurfaces. *ACS Photonics* **5**, 2402–2407 (2018).

18. Sell, D. et al. Periodic dielectric metasurfaces with high-efficiency, multi-wavelength functionalities. *Adv. Opt. Mater.* **5**, 1700645 (2017).

19. Lin, Z. et al. Topology-optimized multilayered metaoptics. *Phys. Rev. Appl.* **9**, 044030 (2018).

20. Lalanne, P. & Morris, G. M. Highly improved convergence of the coupled-wave method for TM polarization. *J. Opt. Soc. Am. A* **13**, 779–784 (1996).

21. Hugonin, J. P., & Lalanne, P. Reticolo Software for Grating Analysis. (Institut d'Optique, 2005).

22. Knop, K. Rigorous diffraction theory for transmission phase gratings with deep rectangular grooves. *J. Opt. Soc. Am.* **68**, 1206–1210 (1978).

23. Li, L. F. Use of Fourier series in the analysis of discontinuous periodic structures. *J. Opt. Soc. Am.* **A** **13**, 1870–1876 (1996).

24. Coppersmith, D. & Winograd, S. Matrix multiplication via arithmetic progressions. *J. Symb. Comput.* **9**, 251–280 (1990).

25. Patterson, D. A. & Hennessy, J. L. Computer Organization and Design. 5th edn. (Morgan Kaufmann, Waltham: MK, 2014).

26. Strehl, K. Aplanatische und fehlerhafte abbildung im femrohr. *Z. für Instrum.* **15**, 362–370 (1895).

27. Arbabi, A. et al. Subwavelength-thick lenses with high numerical apertures and large efficiency based on high-contrast transmitarrays. *Nat. Commun.* **6**, 7069 (2015).

28. Sell, D. et al. Large-angle, multifunctional metagratings based on freeform multimode geometries. *Nano Lett.* **17**, 3752–3757 (2017).

29. Yang, J. J., Hugonin, J. P. & Lalanne, P. Near-to-far field transformations for radiative and guided waves. *ACS Photonics* **3**, 395–402 (2016).

30. Chew, W. C. & Weedon, W. H. A 3D perfectly matched medium from modified Maxwell's equations with stretched coordinates. *Microw. Opt. Technol. Lett.* **7**, 599–604 (1994).

31. Chew, W. C., Jin, J. M. & Michielssen, E. Complex coordinate stretching as a generalized absorbing boundary condition. *Microw. Opt. Technol. Lett.* **15**, 363–369 (1997).
32. Teixeira, F. L. & Chew, W. C. Differential forms, metrics, and the reflectionless absorption of electromagnetic waves. *J. Electromagn. Waves Appl.* **13**, 665–686 (1999).

33. Edee, K., Granet, G. & Plumey, J. P. Complex coordinate implementation in the curvilinear coordinate method: application to plane-wave diffraction by nonperiodic rough surfaces. *J. Opt. Soc. Am. A* **24**, 1097–1102 (2007).

34. Stratton, J. A. & Chu, L. J. Diffraction theory of electromagnetic waves. *Phys. Rev.* **56**, 99–107 (1939).

35. Yang, J. J. & Fan, J. A. Analysis of material selection on dielectric metasurface performance. *Opt. Express* **25**, 23899–23909 (2017).

36. Sell, D. et al. Visible light metasurfaces based on single-crystal silicon. *ACS Photonics* **3**, 1919–1925 (2016).

37. Palmer, C. & Loewen, E. Diffraction Grating Handbook. 6th edn. (New York: Newport Corporation, Rochester, 2005).

38. Yang, J. J., Sell, D. & Fan, J. A. Freeform metagratings based on complex light scattering dynamics for extreme, high efficiency beam steering. *Ann. der Phys.* **530**, 1700302 (2018).

39. Chung, J. W. et al. Seamless on-wafer integration of Si(100) MOSFETs and GaN HEMTs. *IEEE Electron Device Lett.* **30**, 1015–1017 (2009).

40. Shrestha, S. et al. Broadband achromatic dielectric metalenses. *Light: Sci. Appl.* **7**, 85 (2018).

41. Chen, W. T. et al. A broadband achromatic metalens for focusing and imaging in the visible. *Nat. Nanotechnol.* **13**, 220–226 (2018).

42. Paniagua-Dominguez, R. et al. A metalens with a near-unity numerical aperture. *Nano Lett.* **18**, 2124–2132 (2018).