Fullwave Maxwell inverse design of axisymmetric, tunable, and multi-scale multi-wavelength metalenses
Rasmus E. Christiansen, Zin Lin, Charles Roques Carmes, Yannick Salamin, Steven E. Kooi, John D. Joannopoulos, Marin Soljačić, Steven G. Johnson
FFullwave Maxwell inverse design of axisymmetric, tunable, and multi-scale multi-wavelengthmetalenses
Rasmus E. Christiansen, , , , ∗ Zin Lin, , Charles Roques-Carmes, , Yannick Salamin, ,Steven E. Kooi, , John D. Joannopoulos, , , Marin Soljaˇci´c, , , and Steven G. Johnson , Department of Mechanical Engineering, Technical University of Denmark,Nils Koppels All´e, Building 404, 2800 Kongens Lyngby, Denmark NanoPhoton—Center for Nanophotonics, Technical University of Denmark,Ørsteds Plads 345A, DK-2800 Kgs. Lyngby, Denmark. Department of Mathematics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA Research Laboratory of Electronics, Massachusetts Institute of Technology, 50 Vassar St., Cambridge, MA Department of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA and Institute for Soldier Nanotechnologies, 500 Technology Square, Cambridge, MA
We demonstrate new axisymmetric inverse-design techniques that can solve problems radically differentfrom traditional lenses, including reconfigurable lenses (that shift a multi-frequency focal spot in responseto refractive-index changes) and ultra-broadband multi-wavelength lenses ( λ = 1 µ m and µ m). We alsopresent experimental validation for an axisymmetric inverse-designed monochrome lens in the near-infraredfabricated via two-photon polymerization. Axisymmetry allows fullwave Maxwell solvers to be scaled up tostructures hundreds or even thousands of wavelengths in diameter before requiring domain-decomposition ap-proximations, while multilayer topology optimization with ∼ degrees of freedom can tackle challengingdesign problems even when restricted to axisymmetric structures. I. INTRODUCTION
In this paper, we demonstrate that axisymmetric metalenses can be designed with fullwave Maxwell simulations (as opposed tothe scalar-diffraction [1] or domain-decomposition approximations [2] used in prior metasurface designs), for > wavelengths( λ ) in diameters, combined with multilayer variable-height topology optimization (TO) with > ∼ degrees of freedom ( ∼ per λ per layer) as shown in Fig. 1. The capability and flexibility of our design approach is demonstrated by solving two challengingnew design problems with 10-layer metalenses. First (in Sec. III), we design a multi-scale metalens that simultaneously focuses λ = 1 µ m and λ = 10 µ m at the same diffraction-limited focal point (numerical aperature NA = 0 . , Strehl ratios . and . , efficiencies % and %). Second (in Sec. IV), we design an active metalens that shifts its achromatic multi-wavelength(three λ s over a 6% bandwidth) focal spot from NA = 0 . to . as the index of the material (GSS4T1 [3]) is changed from n = 3 . to n = 4 . (thermally or optically), in contrast to previous work that showed only monochromatic reconfigurability [4].As a proof of concept, we also show (in Sec. V) an experimental realization of single-layer axisymmetric TO-designed metalensfor λ = 1550 nm, fabricated by two-photon polymerization 3D-printing (Nanoscribe Professional GT), demonstrating thatsuch variable-height surfaces are manufacturable. As discussed in Sec. VI, our approach could easily be scaled to much largerdiameters and number of layers, and the vast number of design degrees of freedom coupled with the lack of approximationsmakes it a uniquely attractive method for the most difficult metasurface inverse designs.Flat-optics metalenses have received widespread attention due to their potential for achieving multiple functionalities withinan ultra-compact form factor [5–9]. Prior work on metalens design has largely focused on exploiting local resonant conditions [5]under locally periodic approximation (LPA), using rather small unit cells ( < ∼ λ ) [2, 10]. Recently, it has been shown that theunit-cell approach with LPA can lead to fundamental limitations on metalens performance [11, 12], whereas some of theselimitations may be mitigated by using overlapping boundaries, perfectly matched layers or larger domains [13, 14] some maynot. Meanwhile, axisymmetric multi-level diffractive lenses (MDL) have been proposed as an alternative to metalenses forachieving enhanced functionalities [15, 16]; however, MDL designs utilize scalar diffraction theory subject to locally uniformapproximation, neglecting multiple scatterings or resonant phenomena, and, thus, have limited design complexity and physicalbehavior [17]. In contrast to previous works, our approach considers axisymmetric multilayered freeform metalenses which canbe modelled by rigorous fullwave Maxwell equations without any uncontrolled approximations.The prospect of fabricating single-layer metasurfaces with traditional single-step lithography processes is very promising forlarge-scale high-throughput integration [6, 7]. However, single-layer metasurfaces have been limited in their functionality tonarrow angular and spectral bandwidths of operation, with some progress towards chromatic [9] and geometrical aberrationcorrection [18]. Achieving truly multifunctional metasurfaces requires more advanced designs, such as closely-packed multi-layer structures [10, 11, 13, 19]. Recently, there has been a surge in interest in fabricating multilayer metasurfaces [20, 21], ∗ [email protected] a r X i v : . [ phy s i c s . op ti c s ] J u l bolstered by advances in inverse-designed nanophotonics. However, these designs can only be fabricated with more advancedfabrication techniques, such as multi-step lithography, or multi-photon lithography. Multi-photon lithography/polymerization isa technique that enables the fabrication of sub-micron (down to ∼
150 nm) arbitrary three-dimensional structures. Two-photonpolymerization enabled the demonstration of three-dimensional chiral/helical structures [22, 23] and was more recently appliedto the fabrication of supercritical lenses [24, 25] and to the demonstration of full three dimensional control of optical fields usinga metasurface [26]. Nonetheless, the possibility of 3D printing inverse-designed metasurfaces with two-photon lithography pro-cesses remains largely unexplored. In this work, we realize a proof-of-concept experiment with an inverse-designed, freeform,single-layer metalens working at λ = 1550 nm fabricated via two-photon polymerization.The use of TO as a tool for inverse design in nanophotonics has increased steadily over the last two decades [27, 28] with arecent application to metalens design [19]. Our proposed multilayer metalens design framework utilizes density-based topologyoptimization [29] as the inverse design tool. Rather than allowing fully free-form 3D designs, not amenable to nano-scalefabrication, we propose using TO to control the radial height-variation of the N -layers constituting the lens. In addition, wepropose using a filtering technique [30, 31] combined with a threshold operation to limit the gradient of the height variations,making it possible to ensure that they comply with fabrication constraints.FIG. 1: A) Sketch of the multilayer design domain [left] and model domain [right]. B) Illustration of an active metalens,operating at λ = 10 µ m at n = 3 . [left] and n = 4 . [right], showing the max-normalized transmitted | E | -field (thermal) inthe (x,y)-plane through the center of the lens. C)
3D rendering of the active metalens geometry.
II. MODEL AND DESIGN PROBLEM
The design problem is modelled in an axisymmetric domain, Ω , sketched in Fig. 1A(right). The model domain contains adesignable region, Ω D (gray), where the metalens is placed on a solid material surface (dark gray). The sketch also indicates theplane, Γ FF (magenta line), used for computing the far-field transformation in Eq. 1, the focal plane(s), Γ FP (blue and red lines),and the focal spot(s), r FP (black dots), of the lens. Finally, the model domain is truncated using a perfectly matching layer in Ω PML [32, 33]. The lens-design itself consists of N layers of material (Fig. 1A(left)), each with a variable height controlled bythe design field, ¯ ξ L ( r ) . Each designable layer is separated from the next by a layer of air (light gray) and a layer of solid material(dark gray) of fixed thicknesses.The physics is modelled in Ω using the Maxwell Equations [34] assuming time-harmonic behavior. Doing so, we capture thefull wave-behaviour of the electromagnetic field without simplifying assumptions, thus enabling the exploitation of the full wave-behavior to design metalenses exerting precise control of the electromagnetic field. To make the model problem numericallytractable for large design domains, we assume an axisymmetric geometry. This enables the reduction of the full 3D MaxwellEquations to their 2D axisymmetric counterparts [35], therefore significantly reducing the computational effort required to solvethe model problem at the cost of a geometric restriction on the design.The model problem is solved using the scattered-field formulation, E = E b + E s , where the background field, E b , is takento be a planewave propagating along the z-axis (decomposed into two counter-rotating circularly polarized waves). The modelproblem is discretized using the Finite Element Method (FEM) [35] and solved using COMSOL Multiphysics v5.5 [36].A far-field transformation [37] may be used to compute the electric field at any point in space given knowledge of the fieldin the plane Γ FF above the lens, see Fig. 1. The use of this transformation removes the need for simulating the spatial domainbetween the lens and focal point, hereby significantly reducing computational cost. The far field transformation may be writtenas, E FF ( r ) = (cid:90) Γ FF G E ( r , r (cid:48) ) K ( r (cid:48) ) + G H ( r , r (cid:48) ) J ( r (cid:48) ) d r (cid:48) . (1)Here E FF ( r ) denotes the electric far field at the point r , G E ( r , r (cid:48) ) and G H ( r , r (cid:48) ) denotes the electric and magnetic field Green’sfunctions, respectively. Finally K ( r (cid:48) ) and J ( r (cid:48) ) denote the equivalent magnetic and electric surface currents computed from theelectric and magnetic near field obtained by solving the model problem.The figure of merit (FOM) used in the design process is the electric field-intensity at the focal point, r FP . The design problemis formulated as the following continuous constrained optimization problem, max ξ ( r ) ∈ [0 , Φ( ξ ) = | E FF ( r FP , ξ ) | , (2) s . t . A L ≤ ξ L ( r ) ≤ B L , L ∈ { , , ..., N } , N ∈ N (3)Here ξ L ( r ) denotes a radially-varying design field, which controls the thickness of the L ’th layer of the metalens. The electricfield at the focal point, E FF ( r FP , ξ ) , is computed using the solution to the physical model problem and Eq. 1.We propose using a standard PDE-filter [31] to limit the layer-thickness gradient, by applying it to ξ L ( r ) through the choiceof filter radius, r f . After filtering we propose using ξ L ( r ) to control the layer height through the smoothed threshold operation[38] as, ¯ ξ L = 1 − tanh( β · ξ L ) + tanh( β · ( z L − ξ L ))tanh( β · ξ L ) + tanh( β · ( B L − ξ L )) , β ∈ [1 , ∞ [ , ξ L ∈ [ A L , B L ] , z L ∈ [ A L , B L ] . (4)Here z L denotes the spatial position inside each designable layer. The value z L = A L corresponds to the bottom of thedesignable region in the L ’th layer, and the value z L = B L corresponds to the top of the designable region in the L ’th layer. Thethreshold sharpness is controlled by β . In the limit of β → ∞ the field ¯ ξ L takes the value 0 when z L > ξ L and the value 1 when z L < ξ L . A continuation approach may be used to gradually increase β during the inverse design process to enforce a 0/1 finaldesign.The field ¯ ξ L is used to interpolate the relative permittivity, ε r ( r ) , in space between the background material and the materialconstituting the metalens using a linear scheme, ε r ( r ) = ε r,bg + ¯ ξ L ( ε r,lens − ε r,bg ) . (5)Here ε r,bg (resp. ε r,lens ) denotes the relative permittivity of the background (resp. lens).The design problem, Eqs. (2-3), is solved using the Method of Moving Asymptotes [39], for which the sensitivites of theFOM are computed using adjoint sensitivity analysis [40]. Details regarding the modelling and optimization process as well asthe parameter choices for each example are found in Appendix A and Appendix B, respectively.The final designs are all evaluated numerically using a high resolution model by exciting the lens using a linearly polarizedplanewave decomposed into two counter-rotating circularly polarized waves introduced in the model using a first order scatteringboundary condition. An example of a reconfigurable metalens operating at λ = 10 µ m for two different refractive indices isshown in Fig. 1 B . III. MULTI-SCALE MULTI-WAVELENGTH MULTILAYER METALENS
As the first example of our framework we tailor a 10-layer silicon ( n = 3 . ) in air metalens to focus λ = 1 µ m light(Fig. 2 A) and λ = 10 µ m light (Fig. 2 B ) simultaneously at the same focal spot (NA = 0 . ). The lens is 100 µ m in diameterand has a thickness of 10 µ m. The inversely-designed lens is presented in Fig. 2 E with the insert showing an example of thelayer-height variations.FIG. 2: A-B)
Max-normalized | E | -field (thermal) and focal plane (green line) with design overlay (black) in the (x,z)-planethrough the center of the lens for A) λ = 1 µ m and B) λ = 10000 nm planewave excitation. C-D)
Powerflow in the z-directionthrough the Focal plane normalized to the maximum of the Airy disc for C) λ = 1 µ m and D) λ = 10 µ m planewave excitation. E)
3D rendering of the metalens design.From Fig. 2 A -2 B it is clear that the lens exhibits the desired numerical aperture (green line). The focusing capability of thelens is found to reach the diffraction-limit for both wavelengths, when measured in terms of the Full Width at Half Maximum(FWHM) of the main lobe in the focal plane (Fig. 2 C-D ). The Strehl ratio (SR) at the two targeted wavelengths, λ = 1 µ mand λ = 10 µ m, is computed to SR = ≈ . and SR = ≈ . , respectively. The SR is computed based on the power flowthrough the focal plane (blue lines) and the corresponding Airy discs (dashed red lines), shown in Fig. 2 C -2 D . The absolutepower transmission from the substrate of silicon through the lens is computed to T A ,λ ≈ and T A ,λ ≈ , relative tothe incident power in the silicon substrate within the lens diameter. Appendix C includes an additional design example targetingNA = 0 . rather than NA= . while keeping all other parameters fixed, demonstrating the methods versatility. For that secondexample we also achieve diffraction-limited focusing and attain Strehl ratios of SR ≈ . and SR ≈ . for λ = 1 µ m and λ = 10 µ m, respectively. IV. TUNABLE MULTI-WAVELENGTH MULTILAYER METALENS
As a second example of our framework, we design of a 10-layer tunable three-wavelength metalens (see Fig. 1 C ) capable ofshifting the numerical aperture of the lens from NA = 0 . (see Fig. 3[Left Column]) to NA = 0 . (see Fig. 3[Right Column]) bychanging the refractive index of the active material (GST41T1 [3]) from n = 3 . to n = 4 . .The lens is 625 µ m in diameter and has a thickness of 25 µ m. The lens is designed to operate in the mid-infrared region atwavelengths, λ = 9 . µ m (Fig. 3[Row 1]), λ = 10 µ m (Fig. 3[Row 2]) and λ = 10 . µ m (Fig. 3[Row 3]). From Fig. 3 it isFIG. 3: | E | -field normalized to the largest value across the six cases (thermal) at the [Rows] three targeted wavelengths for the[Columns] two targeted values of the refractive index with the focal plane (green line) and design (black) overlaid.observed that the lens exhibits the desired numerical aperture at all three wavelengths for both values of the refractive index. TheStrehl ratio, absolute power transmission and FWHM of the main lobe at the focal point for the three targeted wavelengths andtwo refractive indices are presented in Tab. I. In brief, a Strehl ratio of approximately . is achieved across all six cases withthe spatial focusing being at most from the diffraction limit. Finally, a T A of ≈ . for n = 3 . and of ≈ . for n = 4 . isachieved. V. EXPERIMENTAL VALIDATION OF SINGLE-LAYER VARIABLE-HEIGHT METALENS
Finally, as a proof of concept, we demonstrate experimentally that the proposed method can be used to design variable-heightmetasurfaces for given fabrication specifications (details about the fabrication and experiment are given in Appendix D andAppendix E). Figure 4 G shows a 3D rendering of the designed single-layer varying-height metalens. The metalens is fabricatedvia 3D two-photon polymerization in IP-Dip, a low-refractive-index polymer [41] that can be printed in voxel sizes with in-planefeature sizes ∼ nm and fixed voxel aspect ratio of ∼ to 3. This example is not aimed at designing the largest area lens λ µ m 10.0 µ m 10.3 µ m n = 3 . , Strehl ratio [ · ] ≈ . ≈ . ≈ . n = 4 . , Strehl ratio [ · ] ≈ . ≈ . ≈ . n = 3 . , FWHM main lobe (cid:2) λ (cid:3) ≈ . ≈ . ≈ . n = 4 . , FWHM main lobe (cid:2) λ (cid:3) ≈ . ≈ . ≈ . n = 3 . , T A (cid:104) P lens P inc (cid:105) ≈ . ≈ . ≈ . n = 4 . , T A (cid:104) P lens P inc (cid:105) ≈ . ≈ . ≈ . TABLE I: Strehl ratio, FWHM of main lobe in the focal plane and the absolute power transmission relative to the incident powerin the Si substrate for the lens in Fig. 3.possible nor at achieving the highest possible numerical performance, but at designing a metalens that complies with fabricationconstraints. In this respect, the design is restricted to a diameter of 200 µ m with a 300 nm radial pixel size and a varying heightwith a maximum height of 900 nm, restricted to height-variations in 100 nm increments. The height of the individual radial pixelis allowed to vary independently of its neighbors (i.e. no filtering is applied to ξ L ).The lens is designed to focus λ = 1550 nm light at normal incidence with a numerical aperture of . . The numericallycomputed electric-field intensity at 1550 nm for planewave illumination of the lens at normal incidence is shown in Fig. 4 A ,clearly showing that the targeted numerical aperture (green line) is achieved. Numerically the lens achieves near diffraction-limited focusing in terms of the FWHM of main lobe of the power flow in the z-direction through the focal plane. A FWHMof ≈ nm is computed numerically, corresponding to ≈ . above the diffraction limit (Using the theoretical limit λ ≈ nm).The absolute power transmission from the substrate of IP-Dip through the lens is computed at T A ≈ , relative to theincident power in the IP-Dip substrate within the lens diameter. A Strehl ratio of SR ≈ . is computed by numerical integrationof the power flow over the focal plane. This SR value reveals that a significant fraction of the power is not flowing through thefocal point. From a design point of view, the Strehl ratio is easy to improve using our framework by increasing the designfreedom, either by changing the metasurface material; by decreasing the radial pixel size; by increasing the number of heightincrements; by increasing the total height of the lens and/or by introducing multiple-layers in the lens. All of these weredemonstrated in the two previous examples.Experimentally the Strehl ratio is estimated to be ≈ . by integrating the power flow over an µ m × . µ m region centeredat the focal spot. Computing the SR numerically using the same integration area we obtain SR ≈ . showing that a majorityof the power transmitted through the lens is not focused at the focal spot but flows through the focal plane outside this area. Thediscrepancy between the experimentally measured and numerically computed SR suggests that the experiment overestimatesthe SR, due to the camera’s limited field of view. This is supported by the relatively low measured absolute focusing efficiencyof ≈ . The measured focal spot (Fig. 4(D-F)) exhibits FWHMs of . ± . µm (resp. . ± . µm ) along thehorizontal (resp. vertical) direction, corresponding to ± (resp. ± ) above the diffraction limit. These experimentalresults validate the feasability of freeform axisymmetric metasurfaces experimentally. While this proof-of-concept experimentwas limited to a single-layer metasurface, the radially-varying height of the structure can, to the authors knowledge, only beimplemented with fabrication techniques such as 2.5D lithography or multi-photon polymerization. This is a first step towardsrealizing the full potential of the freeform axisymmetric inverse design technique presented in this work. VI. CONCLUSION
In this paper, we demonstrated that fullwave Maxwell Equation based inverse design of axisymmetric structures can tacklechallenging new design problems involving radically different wavelengths or active materials. We believe that the proposeddesign framework opens the way to many new applications whose functionality goes far beyond traditional lenses, such asend-to-end design [42], hyperspectral imaging [43], depth sensing [44] and nonlinear imaging [45].Computationally, there are many ways to scale our algorithm to much larger designs. The simplest improvement would be toincorporate near-to-farfield transformations [37] to omit the homogeneous region above the lens from the computation, whichwould allow us to increase the radial size by a factor of ∼ . Approximate domain-decomposition could be used to partitiona larger lens into overlapping subdomains solved in parallel (but optimized together) [13]. To increase design freedom, theaxisymmetry could be relaxed to various forms of N -fold or other rotational symmetries. One could also explore fully free-formtopology optimization for 3D-printed structures with manufacturability constraints [46–49].Experimentally, we have shown a proof-of-concept fabricated structure using a two-photon 3D-lithography process. TheFIG. 4: 3D-printed single-layer circular symmetric metasurface. A) Max-normalized electric field intensity | E | (thermal) andfocal plane (green line) with design overlay (black). B) Scanning Electron Micrograph of the full lens. Scale bar = 30 µ m. C) Scanning Electron Micrograph of a smaller area, showing the height variation along the radial direction. Scale bar = 4 µ m. D) Horizontal cut of the focal spot, showing a Gaussian fit to the spot and the corresponding Airy disk (which defines the StrehlRatio). E) Vertical cut of the focal spot. The inset shows the focal spot recorded by the imaging setup (measured on the NIRimaging camera). F) Focal spot measured at various positions along the optic axis. G)
3D rendering of metalens design.inverse-designed metasurface achieved focusing at the telecommunication wavelength of 1550 nm, close to the diffraction limit,with a numerical aperture of 0.4. In the future, we will develop multilayer fabrication of these structures, in order to realizethe full potential of the design technique developed in this work. A key challenge is to realize mechanically stable multilayeredstructures from which unpolymerized resist can be extracted. Application-specific two-photon polymerization setups [24] canachieve more height levels and some control over the voxel aspect ratio. For devices operating at shorter wavelengths, thusrequiring proportionately smaller feature sizes, the design process would shift to multilayer structures with piecewise-constantcross-section [10, 50]. Conversely, at longer wavelengths such as for microwave wavefront shaping, multilayer structures couldbe straightforwardly fabricated, for instance, by stacking multiple stacks of 3D-printed resins or drilled materials [51].
FUNDING
This work was supported in part by Villum Fonden through the NATEC (NAnophotonics for TErabit Communications)Centre (grant no. 8692); the Danish National Research Foundation through NanoPhoton Center for Nanophotonics (grantno. DNRF147); the U. S. Army Research Office through the Institute for Soldier Nanotechnologies (award no. W911NF-18-2-0048); and the MIT-IBM Watson AI Laboratory (challenge no. 2415).
DISCLOSURES
The authors declare no conflicts of interest. [1] M. Meem, S. Banerji, C. Pies, T. Oberbiermann, A. Majumder, B. Sensale-Rodriguez, and R. Menon, “Large-area, high-numerical-aperture multi-level diffractive lens via inverse design,”
Optica , vol. 7, no. 3, pp. 252–253, 2020.[2] R. Pestourie, C. P´erez-Arancibia, Z. Lin, W. Shin, F. Capasso, and S. G. Johnson, “Inverse design of large-area metasurfaces,”
OpticsExpress , vol. 26, no. 26, pp. 33732–33747, 2018.[3] Y. Zhang, J. B. Chou, J. Li, H. Li, Q. Du, A. Yadav, S. Zhou, M. Y. Shalaginov, Z. Fang, H. Zhong, C. Roberts, P. Robinson, B. Bohlin,C. R ˜Aos, H. Lin, M. Kang, T. Gu, J. Warner, V. Liberman, K. Richardson, and J. Hu, “Broadband transparent optical phase changematerials for high-performance nonvolatile photonics,”
Nature Communications , vol. 10, p. 4279, 2019.[4] M. Y. Shalaginov, S. An, Y. Zhang, F. Yang, P. Su, V. Liberman, J. B. Chou, C. M. Roberts, M. Kang, C. Rios, et al. , “Reconfigurableall-dielectric metalens with diffraction limited performance,” arXiv preprint arXiv:1911.12970 , 2019.[5] N. Yu and F. Capasso, “Flat optics with designer metasurfaces,”
Nature Materials , vol. 13, no. 2, p. 139, 2014.[6] M. Khorasaninejad, W. T. Chen, R. C. Devlin, J. Oh, A. Y. Zhu, and F. Capasso, “Metalenses at visible wavelengths: Diffraction-limitedfocusing and subwavelength resolution imaging,”
Science , vol. 352, no. 6290, pp. 1190–1194, 2016.[7] M. Khorasaninejad and F. Capasso, “Metalenses: Versatile multifunctional photonic components,”
Science , vol. 358, no. 6367,p. eaam8100, 2017.[8] F. Aieta, M. A. Kats, P. Genevet, and F. Capasso, “Multiwavelength achromatic metasurfaces by dispersive phase compensation,”
Science ,vol. 347, no. 6228, pp. 1342–1345, 2015.[9] W. T. Chen, A. Y. Zhu, V. Sanjeev, M. Khorasaninejad, Z. Shi, E. Lee, and F. Capasso, “A broadband achromatic metalens for focusingand imaging in the visible,”
Nature nanotechnology , vol. 13, no. 3, pp. 220–226, 2018.[10] Z. Lin, V. Liu, R. Pestourie, and S. G. Johnson, “Topology optimization of freeform large-area metasurfaces,”
Optics Express , vol. 27,no. 11, pp. 15765–15775, 2019.[11] H. Chung and O. D. Miller, “High-NA achromatic metalenses by inverse design,”
Optics Express , vol. 28, no. 5, pp. 6945–6965, 2020.[12] F. Presutti and F. Monticone, “Focusing on bandwidth: achromatic metalens limits,”
Optica , vol. 7, no. 6, pp. 624–631, 2020.[13] Z. Lin and S. G. Johnson, “Overlapping domains for topology optimization of large-area metasurfaces,”
Optics Express , vol. 27, no. 22,pp. 32445–32453, 2019.[14] T. Phan, D. Sell, E. W. Wang, S. Doshay, K. Edee, J. Yang, and J. A. Fan, “High-efficiency, large-area, topology-optimized metasurfaces,”
Light: Science & Applications , vol. 8, no. 1, pp. 1–9, 2019.[15] S. Banerji, M. Meem, B. Sensale-Rodriguez, and R. Menon, “Metalenses or diffractive lenses for imaging?,” in
Imaging Systems andApplications , pp. ITu4B–3, Optical Society of America, 2019.[16] M. Meem, S. Banerji, A. Majumder, J. C. Garcia, P. W. Hon, B. Sensale-Rodriguez, and R. Menon, “Imaging from the visible to thelongwave infrared wavelengths via an inverse-designed flat lens,” arXiv preprint arXiv:2001.03684 , 2020.[17] J. Engelberg and U. Levy, “The advantages of metalenses over diffractive lenses,”
Nature Communications , vol. 11, no. 1, pp. 1–4, 2020.[18] B. Groever, W. T. Chen, and F. Capasso, “Meta-lens doublet in the visible region,”
Nano letters , vol. 17, no. 8, pp. 4902–4907, 2017.[19] Z. Lin, B. Groever, F. Capasso, A. W. Rodriguez, and M. Lonˇcar, “Topology-optimized multilayered metaoptics,”
Physical ReviewApplied , vol. 9, no. 4, p. 044030, 2018.[20] M. Mansouree, H. Kwon, E. Arbabi, A. McClung, A. Faraon, and A. Arbabi, “Multifunctional 2.5 d metastructures enabled by adjointoptimization,”
Optica , vol. 7, no. 1, pp. 77–84, 2020.[21] Y. Zhou, I. I. Kravchenko, H. Wang, J. R. Nolen, G. Gu, and J. Valentine, “Multilayer noninteracting dielectric metasurfaces for multi-wavelength metaoptics,”
Nano letters , vol. 18, no. 12, pp. 7529–7537, 2018.[22] J. K. Gansel, M. Thiel, M. S. Rill, M. Decker, K. Bade, V. Saile, G. von Freymann, S. Linden, and M. Wegener, “Gold helix photonicmetamaterial as broadband circular polarizer,”
Science , vol. 325, no. 5947, pp. 1513–1515, 2009.[23] S. Wu, S. Xu, T. L. Zinenko, V. V. Yachin, S. L. Prosvirnin, and V. R. Tuz, “3d-printed chiral metasurface as a dichroic dual-bandpolarization converter,”
Optics letters , vol. 44, no. 4, pp. 1056–1059, 2019.[24] W. Fang, J. Lei, P. Zhang, F. Qin, M. Jiang, X. Zhu, D. Hu, Y. Cao, and X. Li, “Multilevel phase supercritical lens fabricated by synergisticoptical lithography,”
Nanophotonics , vol. 1, no. ahead-of-print, 2020.[25] J. Coompson, M. Liang, C. Auginash, A. Gin, M. Yang, Z. Qu, I. B. Djordjevic, and H. Xin, “3d-printed phase controlled focusingmetalens at 1550 nm wavelength,” in , pp. 1685–1686, IEEE, 2018.[26] J. W. E. S. J. R. H. Alan Zhan, Ricky Gibson and A. Majumdar, “Controlling three-dimensional optical fields via inverse mie scattering,”
Science Advances , vol. 5, no. 10, pp. 1–5, 2019.[27] J. S. Jensen and O. Sigmund, “Topology optimization for nano-photonics,”
Laser & Photonics Reviews , vol. 5, pp. 308–321, 2011.[28] S. Molesky, Z. Lin, A. Y. Piggott, W. Jin, J. Vuckovic, and A. W. Rodriguez, “Inverse design in nanophotonics,”
Nature Photonics ,vol. 12, pp. 659–670, 2018.[29] M. P. Bends ˜A¸e and O. Sigmund,
Topology Optimization . Springer, 2003.[30] B. B, “Filters in topology optimization,”
International Journal for Numerical Methods in Engineering , vol. 50, pp. 2143–2158, 2001.[31] B. S. Lazarov and O. Sigmund, “Filters in topology optimization based on helmholtz-type differential equations,”
International Journalfor Numerical Methods in Engineering , vol. 86, pp. 765–781, 2011. [32] J.-P. Berenger, “A perfectly matched layer for the absorption of electromagnetic waves,”
Journal of Computational Physics , vol. 114,pp. 185–200, 1994.[33] W. C. Chew, J. M. Jin, and E. Michielssen, “Complex coordinate stretching as a generalized absorping boundary condition,”
Microwaveand Optical Technology Letters , vol. 15(6), p. 363, 1997.[34] D. J. Griffiths,
Introduction to Electrodymanics - Fourth Edition . Pearson Education Limited, 2014.[35] J.-M. Jin,
The Finite Element Method in Electromagnetics - Third Edition . Wiley-IEEE, 2014.[36] “Comsol multiphysics R (cid:13) Computational Electrodynamics: the Finite-Difference Time-Domain Method . Artech house, 2005.[38] F. Wang, B. S. Lazarov, and O. Sigmund, “On projection methods, convergence and robust formulations in topology optimization,”
Structural Multidiciplinary Optimization , vol. 43, pp. 767–784, 2011.[39] K. Svanberg, “A class of globally convergent optimization methods based on conservative convex separable approximations,”
SiamJournal on Optimization , vol. 12, pp. 555–573, 2002.[40] D. A. Tortorelli and P. Michaleris, “Design sensitivity analysis: Overview and review,”
Inverse Problems in Engineering , vol. 1, pp. 71–105, 1994.[41] T. Gissibl, S. Wagner, J. Sykora, M. Schmid, and H. Giessen, “Refractive index measurements of photo-resists for three-dimensionaldirect laser writing,”
Optical Materials Express , vol. 7, no. 7, pp. 2293–2298, 2017.[42] Z. Lin, C. Roques-Carmes, R. Pestourie, M. Soljaˇci´c, A. Majumdar, and S. G. Johnson, “End-to-end inverse design for inverse scatteringvia freeform metastructures,” arXiv preprint arXiv:2006.09145 , 2020.[43] K. Monakhova, K. Yanny, N. Aggarwal, and L. Waller, “Spectral DiffuserCam: Lensless snapshot hyperspectral imaging with a spectralfilter array,” arXiv preprint arXiv:2006.08565 , 2020.[44] Q. Guo, Z. Shi, Y.-W. Huang, E. Alexander, C.-W. Qiu, F. Capasso, and T. Zickler, “Compact single-shot metalens depth sensors inspiredby eyes of jumping spiders,”
Proceedings of the National Academy of Sciences , vol. 116, no. 46, pp. 22959–22965, 2019.[45] C. Schlickriede, S. S. Kruk, L. Wang, B. Sain, Y. Kivshar, and T. Zentgraf, “Nonlinear imaging with all-dielectric metasurfaces,”
NanoLetters , 2020.[46] R. E. Christiansen, B. S. Lazarov, J. S. Jensen, and O. Sigmund, “Creating geometrically robust designs for highly sensitive problemsusing topology optimization - acoustic cavity design,”
Structural and Multidisciplinary Optimization , vol. 52, pp. 737–754, 2015.[47] M. Zhou, B. S. Lazarov, F. Wang, and O. Sigmund, “Minimum length scale in topology optimization by geometric constraints,”
ComputerMethods in Applied Mechanics and Engineering , vol. 293, pp. 266–282, 2015.[48] F. Wang, R. E. Christiansen, Y. Yu, J. Mørk, and O. Sigmund, “Maximizing the quality factor to mode volume ratio for ultra-smallphotonic crystal cavities,”
Applied Physics Letters , vol. 113, p. 241101, 2018.[49] Q. Li, W. Chen, S. Liu, and L. Tong, “Structural topology optimization considering connectivity constraint,”
Structural and Multidisci-plinary Optimization , vol. 54, pp. 971–984, May 2016.[50] A. Y. Piggott, J. Petykiewicz, L. Su, and J. Vuˇckovi´c, “Fabrication-constrained nanophotonic inverse design,”
Scientific reports , vol. 7,no. 1, pp. 1–7, 2017.[51] P. Camayd-Mu˜noz, C. Ballew, G. Roberts, and A. Faraon, “Multifunctional volumetric meta-optics for color and polarization imagesensors,”
Optica , vol. 7, no. 4, pp. 280–283, 2020.
APPENDICESAPPENDIX A. OPTIMIZATION AND NUMERICAL MODELLING
The physics is modelled in COMSOL Multiphysics [36] and the optimization problem is solved using the Globally Conver-gent Method of Moving Asymptotes (GCMMA) [39].In the design process Ω D and the solid material regions in Ω are discretized using a structured quadrilateral mesh, while thesurrounding air regions are discretized using an unstructured triangular mesh, both of which uses ≥ elements per λ/n . Thefinite element method with a linear Lagrangian basis is used to discretize the physics [35].The following stopping criterion is used to terminate the iterative solution of the optimization problem: if i ≥ i min thenif | Φ i − Φ i − n | / | Φ i | ≤ . ∀ n { , , ..., } then Terminate optimization. end ifend if
Here i denotes the current optimization iteration, i min = 70 denotes the minimum number of design iterations taken. Φ i denotesthe objective function value at the i ’th iteration and n ∈ N + .0 APPENDIX B. STUDY PARAMETERS
The parameters used in setting up the models and associated optimization problems for the three examples follow here.
B. 1. Multi-scale multi-wavelength multilayer metalens
For the problem treated in Sec. III the following parameter values are used:The axisymmetric model domain Ω has a width of 57 µ m in the r-direction and a height of 82 µ m in the z-direction. Ω issurrounded on three of four sides by a perfectly matched layer with a depth of 1500 nm (Fig. 1). The metalens design domain Ω D is taken to have a radius of 50 µ m and a height of 10 µ m and is separated into ten layers of equal height. Each layer hasa total height of 1 µ m with the designable region having a height of 600 nm and the fixed air and silicon regions each havingheights of 200 nm. It is placed on a slab of material of 2 µ m thickness placed at the bottom edge of the model domain.The radial design pixel size is restricted to a minimum of 200 nm and the height-variation is restricted to 25 nm increments.The two wavelengths of the incident field are taken to be λ = 1 µ m and λ = 10 µ m. The lens is taken to be made of silicon inan air background. The refractive index of air are taken to be n air = 1 . . The refractive index of silicon is taken to be n si = 3 . at both operating wavelengths. The speed of light is taken to be c = 3 · m/s. The numerical aperture is taken to be NA = 0 . .The initial guess for the design field is ξ L, initial ( r ) = 0 . ∀ r ∈ Ω D for all 10 layers. A filter radius of r f = 400 nm is used tolimit the gradient of the heigh variation in each layer to avoid rapid pixel-by-pixel oscillations in the design. The value of thethresholding sharpness parameters is β = 40 . B. 2. Tunable multi-wavelength multilayer metalens
For the problem treated in Sec. IV the following parameter values are used:The axisymmetric model domain Ω has a width of 342.5 µ m in the r-direction and a height of 380 µ m in the z-direction. Ω issurrounded on three of four sides by a perfectly matched layer with a depth of 15 µ m (Fig. 1). The metalens design domain Ω D is taken to have a radius of 312.5 µ m and a height of 25 µ m and is separated into ten layers of equal height with a 2000 nmdesignable region and 250 nm fixed air region and 250 nm fixed solid region. It is placed on a slab of material of 5 µ m thicknessplaced at the bottom edge of the model domain.The radial design pixel size is restricted to a minimum of 600 nm and the height-variation is restricted to 100 nm increments.The three wavelengths of the incident field are taken to be λ = 9 . µ m, λ = 10 µ m and λ = 10 . µ m. The lens is taken to bemade of GST41T1 in an air background. The refractive index of air are taken to be n air = 1 . . The refractive index of the activematerial is taken to be n GST , = 3 . in the first configuration and n GST , = 4 . in the second at all operating wavelengths. Thespeed of light is taken to be c = 3 · m/s. The numerical aperture of the lens is taken to be NA = 0 . in the first configurationand NA = 0 . in the second.The initial guess for the design field is ξ L, initial ( r ) = 0 . ∀ r ∈ Ω D for all 10 layers. A filter radius of r f = 3 µ m is used to limitthe gradient of the height-variation in each layer (see the insert in Fig. 1 C ). The value of the thresholding sharpness parametersis β = 40 . B. 3. Single-layer variable-height metalens
For the problem treated in Sec. V the following parameter values are used:The axisymmetric model domain Ω has a width of 106 µ m in the r-direction and a height of 301.8 µ m in the z-direction. Ω issurrounded on three of four sides by a perfectly matched layer with a depth of 3 µ m (Fig. 1). The metalens design domain Ω D is taken to have a radius of 100 µ m and a height of 900 nm and comprises a single layer constituting the designable region. Thedesign domain is placed on a slab of material of 500 nm thickness placed at the bottom edge of the model domain.1The design is discretized into 300 nm radial increments and 100 nm height increments.The wavelength of the incident field is taken to be λ = 1550 nm. The lens is taken to be made of IP-Dip in an air background.The refractive index of air are taken to be n air = 1 . . The refractive index of IP-Dip is taken to be n si = 1 . at both operatingwavelengths. The speed of light is taken to be c = 3 · m/s. The numerical aperture is taken to be NA = 0 . .The initial guess for the design field is ξ L, initial ( r ) = 0 . ∀ r ∈ Ω D . No smoothing filter is applied. The value of the thresholdingsharpness parameters is β = 40 . APPENDIX C. SECOND EXAMPLE OF A MULTI-SCALE MULTI-WAVELENGTH MULTILAYER METALENS DESIGN
We tailor a 10-layer silicon ( n = 3 . ) in air metalens to focus λ = 1 µ m light (Fig. 2 A) and λ = 10 µ m light (Fig. 2 B )simultaneously at the same focal spot (NA = 0 . ). The lens has identical dimensions and design resolution as the lens in Sec. III.The final lens design is presented in Fig. 5 E with the insert showing an example of the layer-height variations.Figures 2 A -2 B show that the lens exhibits the desired numerical aperture at both wavelengths (green line). Further, thefocusing capability of the lens is diffraction-limited for both wavelengths. The Strehl ratio (SR) at the two targeted wavelengths, λ = 1 µ m and λ = 10 µ m, is computed to SR = ≈ . and SR = ≈ . , respectively, from the data in Fig. 2 C -2 D . APPENDIX D. FABRICATION
The metalens was fabricated using a commercial two-photon polymerization system (Nanoscribe Photonic Professional GT)on a -micron-thick fused silica substrate, where the structures are written in circles with height increments of nm . Forthis purpose, piezo actuators move the sample in the out-of-plane direction after fabricating each layer. Geometrical parametersand dose (scanning speed and laser power) are optimized with a dose test on this specific machine. In the in-plane direction thelaser beam is guided by galvanometric mirrors parallel to the substrate. After printing, the structures are put in a developer bath(PGMEA 5 min) and dried in IPA with a critical point dryer Auto Samdri 815 Series A. APPENDIX E. EXPERIMENT
For the proof-of-concept experimental results presented in Fig. 4, we used the imaging setup shown in Fig. 6(a). A AndoAQ4321D Tunable Laser Source produces a fiber-coupled output at 1550 nm. The fiber output is collimated with a set of lenses.In the measuring configuration Fig. 6, the collimated beam is focused by the metasurface, and the focal spot is imaged by anobjective - tube lens - IR imaging camera system. For this measurement, we used a 100X Mitutoyo Plan Apo NIR HR InfinityCorrected Objective, a ThorLabs f = 200 mm tube lens, and a EC MicronViewer 7290A. The imaging setup was first calibratedusing the configuration shown in Fig. 6(b), where the equivalent pixel size on the detector is evaluated by imaging a USAF1951target. To evaluate the efficiency of the metasurface, we measured the equivalent power going through a µ m diameterpinhole with the configuration shown in Fig. 6(c).To estimate the metasurface efficiency, we use the intensity-voltage relation of the NIR camera provided by the vendor.It has the form I = KV /gs , where I is the incident optical power on a pixel, V s the generated voltage at that pixel, and g the characteristic nonlinear slope of the intensity-voltage relation, which is given to be g ∼ . . We first calibrate theproportionality constant K by measuring the signal produced by the camera of a known beam power. This allows us to translatethe measured voltage on a pixel to an incident power (in W). We also measure the incident intensity on the metasurface area withthe experimental configuration shown in Fig. 6(c). The efficiency is then calculated asEff = KL (cid:80) i ∈ pixels in focal spot V /gi P ref , where L is the estimated optical loss through the objective and tube lens, which is 0.55 (objective) × A-B)
Max-normalized | E | -field (thermal) and focal plane (green line) with design overlay (black) in the (x,z)-planethrough the center of the lens for A) λ = 1 µ m and B) λ = 10000 nm planewave excitation. C-D)
Powerflow in the z-directionthrough the Focal plane normalized to the maximum of the Airy disc for C) λ = 1 µ m and D) λ = 10 µ m planewave excitation. E)
3D rendering of the metalens design.3100x, NA 0.7objectivemetasurface tube lens IR imaging camera1550 nm collimated source 100x, NA 0.7objectiveUSAF1951 Target tube lens IR imaging camera1550 nm collimated source 200 μm pinhole1550 nm collimated source Power meter(a)(b)(c)