Fundamental limits to multi-functional and tunable nanophotonic response

Tunable and multi-functional nanophotonic devices are used for applications from beam steering to sensing. Yet little is understood about fundamental limits to their functionality. The difficulty lies with the fact that it is a single structure that must exhibit optimal response over multiple scenarios. In this article, we present a general theoretical framework for understanding and computing fundamental limits to multi-functional nanophotonic response. Building from rapid recent advances in bounds to light-matter interactions, we show that after rewriting the design problems in terms of polarization fields, the introduction of suitable cross-correlation constraints imposes the crucial"single-structure"criteria. We demonstrate the utility of this approach for two applications: reflectivity contrast for optical sensing, and maximum efficiency for optical beam switching. Our approach generalizes to any active or multi-functional design in linear optics.


I. INTRODUCTION
Nanophotonic devices that offer multiple functionalities, from liquid-crystal devices for beam steering [1][2][3][4] to polychromatic metasurface lenses [5][6][7][8][9][10][11][12], have tremendous design complexity due to the need for a single structure to simultaneously optimize multiple objectives.In this Letter, we develop a general framework for identifying fundamental limits to achievable response in multifunctional devices.Despite significant recent interest in identifying nanophotonic bounds [13][14][15][16][17][18][19][20][21][22][23][24][25][26][27][28][29], such works have almost exclusively applied only to devices with only a single function.Here, building from recent discoveries that design problems can be transformed to quadratically constrained quadratic programs [26,30], we introduce a simple mechanism for constructing "cross-correlation" constraints that encapsulate simultaneous requirements on a single device.Such constraints can be utilized for maximum functionality across multiple frequencies, incident fields, and constituent-material properties, as well as active modulation.We identify two prototypical examples where such bounds illuminate the limits to what is possible: multi-frequency reflection control in optical filters, for applications such as optical sensing, and optimal beam switching via liquid-crystal-based nanophotonic metagratings.Our framework promises to identify the limits to complex optical functionality for wideranging applications.
Fundamental limits to optical response have dictated technological selection for decades, ranging from the Wheeler-Chu limits to broadband antenna design for mobile devices [31][32][33][34] to photovoltaic energy-conversion efficiency approaching the canonical Shockley-Queisser limits [35][36][37][38].Similarly, imaging techniques are mea-sured against the Abbe diffraction limit, while allangle, broadband absorbers are designed against the Yablonovitch limit [39].Yet none of these approaches account for the full complexity of wave-scattering physics allowed by Maxwell's equations, instead utilizing smallsize (Wheeler-Chu), free-space propagation (Abbe), or ray-optical (Yablonovitch) regimes in which the physics is dramatically simpler.
Recent progress has suggested the possibility for identifying bounds in more complex, "full-wave" regimes.After developing the idea of "communication channels," which limit information capacities transmitted between known scattering bodies [40][41][42], D. A. B. Miller utilized single-frequency scattering sum rules to identify bounds on corralling and separating pulses, and more complex functionality [43,44].Those bounds, however, take specialized consideration of currents non-orthogonal to "pass-through" and "single-pass" waves in a way that appears difficult to extend beyond one-dimensional structures.All-frequency sum rules utilize complex-frequency contour integrals to connect all-frequency response to relatively simple constants (related to either electrostatic response or electron densities) [19,[45][46][47][48][49][50][51][52], but can only be utilized for very special response functions such as extinction, and they cannot meaningfully bound response over smaller frequency ranges.
Analytical bounds to single-frequency response have been developed using a variety of constraints that essentially distill to energy-conservation constraints over the relevant scattering bodies [13-18, 20, 22, 25, 29, 53-57].One can accommodate multiple functionalities in the sense of, for example, maximizing scattering while constraining absorption [17,29], but multiple scenarios cannot be accounted for.It is possible to combine the contour-integral, broad-bandwidth approach with the energy-conservation-based, single-frequency approach to identify bounds over any frequency range [19], but such bounds are feasible only for response functions with very special optical-theorem-like structure, such as local den-sity of states [19] and extinction [26].Tighter singlefrequency bounds can be achieved with a recently developed computational approach in which many "local" conservation laws are identified that must necessarily be satisfied by any design [26,30].The mathematical framework of this approach appears to generalize to physical design problems across many domains; one successful application has been to quantum optimal control [58].A unique instance of multiple functionalities, that of multiple incident fields (or initial wavefunctions in the quantum case) was proposed in Refs.[58,59], while an alternative formulation of multi-scenario design was proposed in Ref. [60].However, the approach of Ref. [60] to multiple scenarios does not actually incorporate any constraints between scenarios, and hence the bounds for many scenarios simplify to a weighted average of the single-scenario bounds.
In this article, we combine the recently developed integral-equation framework for optical-response bounds [26,30] with a new approach to capturing constraints that must be satisfied in multi-functional devices.The key is to identify generalized conservation laws that relate the polarization fields induced in one scenario to the incident and scattered fields in another (Sec.II).These conservation laws are necessary and sufficient conditions for any design solution, and therefore one can replace the Maxwell-equation constraints, which are nonlinear (and non-polynomial) in the structural design variables, with a set of quadratic constraints in the polarization fields.From this quadratically constrained quadratic program (QCQP), one can "lift" the degrees of freedom to a higher-dimensional space, and relax the constraints, such that the problem becomes a semidefinite program (SDP) whose global optimum can be computed efficiently [61].Equation (5) represents the culmination of the transformations to a QCQP; the solutions of the corresponding SDP represent fundamental limits to what is possible is multi-functional systems.We apply this approach to two classes of applications in Sec.III, where we show that the limits illuminate fundamental constraints on traits such as size, sensitivity, optical contrast, and switching efficiency in multi-functional devices.

II. A FRAMEWORK FOR MULTI-FUNCTIONAL BOUNDS
We start with a conventional statement of the multifunctional Maxwell design problem, of maximizing some objective subject to the Maxwell-equation constraints.We show that the problem can be transferred to one with a different mathematical structure, a QCQP, without any approximation or relaxation, through identification of power-conservation and cross-correlation constraints embedded in the structure of Maxwell's equations.Given the QCQP formulation, one can "relax" the problem through now-standard techniques (e.g.semidefinite relaxations), which leads to problems that are ef-ficiently solvable and whose solutions represent global bounds on what is possible.
In a Maxwell design problem, one maximizes (or minimizes) some objective subject to the constraint of Maxwell's equations over all scenarios of interest.For notational simplicity, we will assume here nonmagnetic, isotropic materials in free space.(The generalization to anisotropic materials in arbitrary backgrounds is straightforward).For N scenarios of interest (incident fields, frequencies, etc.), if we denote the field of scenario s as E s (x), and take the designable geometry as a spatially dependent electric susceptibility χ(x), then one can mathematically write the design problem as max.
for all s from 1 through N , where we have taken e −iωst time-dependencies and units ε 0 = µ 0 = 1.The key mathematical difficulty in solving such design problems is that the Maxwell equations are nonconvex in the geometrical degrees of freedom χ s (x); in fact, the design problem is in the class of NP hard problems [60].We will show, however, that under suitable transformations, bounds can be identified via efficient algorithms.
Our first step is to transform the differential form of Eq. ( 1) to an equivalent integral-equation form, which is then amenable to the transformation to a QCQP.The volume-integral formulation arises by separating the total fields E into its incident and scattered components, E inc and E scat , respectively, and then re-writing the total and scattered fields instead in terms of the polarization currents that generate them.First, the definition of susceptibility implies that E(x) = χ −1 (x)P(x), for polarization field P and at any point x in the material (such that χ(x) = 0).Second, the scattered field is the field radiated by the polarization currents as though they were in free space, which is given by the convolution of the freespace (or background) dyadic Green's function with the polarization field: E scat (x) = ´V G 0 (x, x )P(x ), where V is the volume of the designable region.Hence one can replace the differential form of Maxwell's equation in Eq. ( 1) with the equivalent integral equations [62], again over all scenarios s.Again, however, the expressions are nonconvex and hard to work with due to the presence of both variables χ(x) and P(x) in the expressions.Also, Eq. ( 2) is only valid at points within the material for which χ(x) = 0, which is an unknown domain before the design problem has been solved.We can simplify the integral constraint into a simpler quadratic form that extends to all design space V through the following idea.In the designable region, the desired material is either present or absent.If the material is absent, Eq. ( 2) simplifies to P s (x) = 0, as there is no polarization field.If the material is present, then at that point we know χ s (x) equals the constant value χ s , without position argument, and Eq. ( 2) simplifies to a linear equation in P s .Hence at any point x in the design region, either P s (x) = 0 (in the absence of material), or χ −1 s P s (x) − ´V G 0,s (x, x )P s (x ) dx − E inc,s (x) = 0 (in the presence of material).We can combine these possibilities into a single, quadratic constraint, by multiplying them together (and conjugating the first equation for mathematical convenience): If we enforce the constraint of Eq. ( 3) at all points x in the designable region, it is equivalent to enforcing the integral (or differential) Maxwell constraints.One can show that the real and imaginary parts of Eq. ( 3) represent conservation of real and reactive power flow through an infinitesimal bounding surface centered at x.These constraints (or equivalent versions of them [22]) were used in Ref. [26] to identify bounds in single-function photonic structures.Their key drawback, however, is that they do not enforce any restrictions between different scenarios.For example, if one wants a single structure to exhibit optimal response at two frequencies, there is no restriction in Eq. ( 3) that the polarization fields in the two scenarios must arise from the same structure.Hence, such bounds (including those discussed in Ref. [60]) yield multi-functional design bounds that are often trivial in some way.

A. Cross-correlation constraints
Here, we show that one can naturally enforce the desired constraints through the creation of new, "crosscorrelation" constraints.The key idea is to consider the two terms multiplied together in Eq. ( 3), and to take similar products of the terms across all combinations of different scenarios.Mathematically, this means adding constraints of the form for all combinations of scenarios, i = j.Why does this enforce the desired constraints?Consider a given i and j.Equation (4) requires that either the material in scenario i is absent (hence P si = 0), or the material in scenario j is present.This achieves half of the necessary requirement: if the material in scenario i is present, and hence the first condition fails, then the material in scenario j must also be present.The reverse condition, that the absence of

With cross-correlation constraints
One structure for all scenarios FIG.1: Previous approaches to identifying fundamental limits could not account for the constraint that a single device must account for response across all scenarios with a multi-functional objective (left).We introduce cross-correlation constraints in Eq. ( 4) that correctly account for a single device (right), yielding meaningful bounds on what is possible.one material must imply the absence of the other, arises from the equivalent set of constraints in which i and j are reversed.Figure 1 illustrates the difference between imposing the cross-correlation constraints and not: without the cross-correlation constraints, different optimal structures would implicitly be considered for each scenario.Using the cross-correlation constraints, one can enforce a single structure that must be utilized for each scenario.

B. Optimization problem
The transformations discussed above lead us to a new formulation of the Maxwell design problem, in terms of the polarization fields for all scenarios, P s (x): max.
Equation ( 4) satisfied for all i, j and all x in the design region. ( Any typical objective f of interest in linear optics (scattered/absorbed power, spontaneous-emission enhancement, mode overlap, etc.) will be a quadratic (or linear) function of the polarization fields P s , in which case Eq. ( 5) is a quadratic objective subject to quadratic constraints, and is thereby a QCQP.Formulating a problem as a QCQP does not necessarily make it easy to solve; generically, QCQPs are NP hard to solve [63].However, they arise in wide-ranging problems across engineering and physics [63][64][65][66][67][68], and numerous tools have been developed to aid in their solution.In particular, semidefinite relaxations provide a mechanism for relaxing QCQPs to semidefinite programs (SDPs), which can be solved efficiently using interior-point methods [69].
The general prescription to relax a QCQP is to "lift" the variable to a higher-dimensional space in which the quadratic constraints are transformed to linear constraints [63].In our case, let us assume that we stack the space-, polarization-, and scenario-dependent polarization-field values into a single vector p.Then our constraints are of the form p † Ap + y † p = 0.The difficult part is the nonconvex quadratic term p † Ap.That term can be re-written as Tr App † ; if we define a new matrix variable X = pp † , then the constraint now has the linear form Tr (AX) + y † x.The variable X is not free to take any form, however, as it must be a rank-one positive semidefinite matrix.The rank-one constraint is the only nonconvex constraint, and it can be dropped to transform the problem to a convex SDP.This is the only necessary relaxation, and it is the point at which any looseness (e.g.non-enforcement of the single-structure constraint) may arise.Finally, the corresponding SDP can be solved using standard packages [70,71], and the resulting solution is a fundamental limit on the multifunctional design problem of interest.

III. APPLICATIONS
To illustrate the generality and utility of the bound framework established in Sec.II, we demonstrate its ability to predict fundamental limits for two prototype problems: the design of maximum-reflectivity-contrast optical filters, as may be useful for sensing applications, and the design of a liquid-crystal-based beam-switching device.In each case we show that previous bounds trivialize in some way, whereas our bounds can be quite close to (and always above) the performance of real theoretical designs.These applications demonstrate that our bounds identify fundamental limits to what is possible, and that they tease out the key scaling laws of multi-functional systems.

A. Maximum reflectivity contrast
A key problem in sensing applications [72][73][74][75] is to detect changes, however minuscule, in the incident frequency of light.To do so, the scatterer must be designed to maximize the variation in its response with respect to changes in incident frequency.While there are several candidates for such optical response, we focus here on reflectance, which naturally arises in many applications [74,75], and on multilayer films, i.e., optical filters.The problem is depicted schematically in Fig. 2, with any combination of angles and frequencies for the incident, reflected, and transmitted waves, while the primary objective might be the normal-incident contrast in reflectivity between two nearby frequencies.We take the multilayer film to comprise alternating layers of Al 2 O 3 (as is commonly used for multilayer thin films [76,77]) and vacuum.Reflectance is a quadratic function of the polarization fields.Specifically, the reflection coefficient is the amplitude of the back-scattered wave, which can be evaluated through a convolution of the free-space Green's function with the polarization field, and the reflectance R(ω) is the square of this quantity.Our objective is the contrast in reflectivity at any two frequencies of interest, FIG.2: The first application we consider is that of maximum reflectivity contrast between two nearby frequencies, ω 1 and ω 2 , on a multilayer structure.Our bounds capture tradeoffs between the key parameters such as designable-domain size, refractive index, frequency separation, and reflectivity contrast.
which we can write as R(ω 1 ) − R(ω 2 ).Our set of constraints is those of Eq. ( 5), including the cross-correlation constraints between the polarization fields at the two frequencies.
The solid black lines in Fig. 3(a) show the global bounds for reflectivity contrast as a function of the maximum thickness L of the designable region (normalized to wave number k = ω/c).Also included are the bounds without the cross-correlation constraints, as well as actual reflectivity contrasts achieved by designs identified through a local-optimization routine.For each thickness of the designable region, we used gradient descent, a standard optimization technique [61], to optimize over many initial conditions (random designs) and select the bestperforming design.In each case, one can see that the results of the local optimizations are able to quite closely approach, though not surpass, the global bounds.This suggests that the bounds have very little "slack," i.e., are nearly the best possible bounds.The bounds without the cross-correlation constraints, shown in grey, are not as close to the local-optimization results.In fact, the bounds without the cross-correlation constraints are trivial in the sense that they always identify zero reflectance as possible for R(ω 2 ) (since an all-vacuum structure would be reflectionless), and there are no constraints requiring the same structure to also provide large reflectance at a nearby frequency.One can see that for two relatively close frequencies, it is vital to incorporate the cross-correlation constraints to achieve meaningful bounds.
The inability of previous bound approaches to identify key features of multi-functional bounds is further illustrated in Fig. 3(b), which sweeps across a variety of relative frequency values (holding ω 1 fixed).In each case, the solid line depicts the global bound, which tracks quite closely with the best designs from the local-optimization computations.By contrast, the dashed lines show the bounds without the cross-correlation constraints, which trivially do not change as a function of relative frequency because R(ω 2 ) is always bounded below by 0. Using the cross-correlation constraints is crucial to capturing the difficulty of achieving high contrast for small changes in frequency.
Finally, we can use this approach to identify a bound on the minimum thickness of any design that could possibly achieve a desired reflectivity contrast ∆R as a function of the relative bandwidth of the two frequencies of interest, as shown in Fig. 3(c).As expected, the minimum thickness increases both as a function of the desired reflectivity difference as well as as a function of decreasing relative bandwidth.But the scaling lines of Fig. 3(c) that indicate precise tradeoffs between sensitivity and size would not be possible with any other approach.

B. Beam switching via liquid crystals
Another prototypical example of a multi-functional optical device is one in which a voltage is applied to change the refractive index and modulate the optical response.Here, we consider an example of liquid-crystalbased beam switching, in which the target functionality is to direct light in one direction for one voltage state, and another direction for another voltage state.The design of such unit cells could be useful not only for twostate switching in periodic structures, but also for manystate switching in a larger composition of varying unit cells [78,79].
We consider periodic grating structures with the active designable medium comprising a liquid crystal (LC).(More realistic models with a LC filling a patterned material such as silicon can be similarly modeled [80].)For a period Λ and grating thickness L (as shown in the top of Fig. 4), we consider any possible pattern of the LC material.We choose E7 [81] as the LC material, which has refractive indices of about 1.7 and 1.5 in the voltage-on and off state respectively (for purposes of demonstration, we allow for a small Im n of 0.1 for both cases).Given a monochromatic field at normal incidence, we choose our objective to be the power switching efficiency, defined to be the sum of the power in the target directions for the voltage-on and voltage-off.We choose diffraction orders m and −m, with angles θ m = sin −1 (2πm/kΛ) and θ −m = −θ m with respect to the normal, to define the target directions for the on and off states, respectively.The power diffracted into the two states, P m and P −m , are quadratic functions of the polarization field and naturally fit our bound framework.
Figure 4(a) shows bounds and designs for the liquidcrystal beam-switching problem with the unit cell taken to be about 5 free-space wavelengths wide (kΛ = 10π+δ), and take m = 1.(The constant δ = 0.1 is added to avoid the singularities and infinite-Q resonances possible via bound-state-in-continuum modes [82], which cause numerical instabilities.)Hence the angular deflection is small (±11.5 • ), and one might expect a perfect switching efficiency to be relatively easy to achieve.Yet the solid black line shows the bound computed using the crosscorrelation constraints, which even for significant thickness (kL = 8) shows the possibility to reach only about 70% of the maximum possible switching efficiency (which is 2).We also consider three possible design strategies: a fill-fraction optimization (red), in which the optimal fill fraction of a simple grating structure (fill fraction here specifies the fraction of the material occupying the unit cell) is computed, a "const-z" optimization (blue), in which the all air holes must have equal depth to be compatible with lithography, and a freeform optimization, in which the permittivity is allowed to take either material value at any point in the domain (teal).The freeform optimizations have the most degrees of freedom, and come closest to the bounds, while the fill-fraction optimizations are feasible at small thicknesses but deteriorate in quality at larger thicknesses.The freeform and const-z approaches both show similar trendlines to the computed bounds.By contrast, the bounds without cross-correlation constraints quickly approach 2, which is a trivial value that would suggest 100% possible efficiency in both states.Figure 4(b) isolates a single thickness kL = 3 for a larger unit-cell period kΛ = 20π and sweeps over the target diffraction order m, with angular deflections increasing to 64 • for m = 9.Perhaps surprisingly, the bound suggests that at this large unit-cell period, the switching efficiency can increase as the angular of deflection increases.This is borne out by the constz optimizations, which show a similar trend (though the noisiness of local optimizations makes the trend less clear than the bounds do).

IV. EXTENSIONS
The essence of our approach is as follows: starting from a linear partial differential equation (PDE) that governs the response of field variables under multiple sce-narios, one can transform the PDE to an integral equation in terms of polarization-field variables.Next, one can transform the linear but domain-dependent integral equations to quadratic, domain-independent equations.Crucially, those quadratic constraints largely comprise cross-correlation equations that are necessary to describe multi-scenario problems.Upon reaching this QCQP, standard methods can be used to compute bounds via semidefinite programming.
No part of our approach requires the linear PDE to be Maxwell's equations.And, indeed, this approach can be applied to applications such as quantum optimal control [58,[83][84][85][86][87], and potentially to applications in acoustics, elasticity, and more [88].Of course, within nanophotonics, the choice of an electric susceptibility and electricfield variable as written in Eq. ( 1) is not unique, and this approach can be applied for arbitrary materials and for electric and/or magnetic fields.
Looking forward, a key opportunity is in the development of faster computational algorithms for solving the semidefinite programs that arise in this approach.Interior-point algorithms are highly effective but scale at least with the cube of the number of degrees of freedom of the design [89].First-order methods such as ADMM [90] have been shown to achieve near-linear scaling for certain instances of SDPs [91], and would come with an additional benefit: integral operators have a special hi-erarchical rank structure [92], oft-utilized for fast multipole methods [93,94], for example, that could provide significant further speed improvements of matrix multiplications with the integral operators.
In addition to improving SDP computation times, another avenue is to reduce the total number of constraints that are necessary.One approach, as developed in Ref. [26], is to iteratively choose only the constraints that are most violated by a given solution.Another approach, unique to the multi-scenario case, is to choose a minimal set of pairs of indices for which to impose the cross-correlation constraints: instead of enforcing crosscorrelation across every pair, for example, enforcing it be-tween pairs 1 and 2, then 2 and 3, and so on, until N and N − 1, would actually enforce exactly the same "singlestructure" conditions but with ∼ N constraints instead of ∼ N 2 .The QCQP would be unchanged, though it is unclear whether the resulting SDP formulations would be.This is connected with a well-known opportunity in SDPs, in which at times it can be beneficial to impose "redundant" constraints in a quadratic formulation but lead to better SDP bounds [95].
Connecting nanophotonics design problems to QCQPs in optimization theory offers a new lens with which to view them, and a new approach for understanding what is possible.

FIG. 3 :
FIG. 3: (a) Comparison of bounds on reflectivity contrast, both with and without cross-correlation constraints, to designs achieved via gradient-based local optimization.The relative frequency ω 2 /ω 1 is set to 0.9, and kL is the normalized thickness of the designable region k = 2π/λ.The inset shows a locally-optimized design at kL = 9.3, with reflectance of 0.79 at ω 1 and 0.005 at ω 2 .(b) Similar to (a), but as a function of relative frequency ω 2 /ω 1 , where kL is fixed at 3 and 5. (c) The minimum thickness required for different values of desired reflectivity contrast, as predicted by our bounds with cross-correlation constraints.Since reflectance ranges from 0 to 1, the reflectivity contrast ranges from -1 to 1.

FIG. 4 :
FIG.4:(Top) Liquid-crystal-based beam beam-switching device.(a,b) Bounds on switching efficiency, both with and without cross-correlation constraints, compared to 3 different designs: simple grating structure with optimal fraction of material (fill-fraction optimization), lithography-friendly optimization with equal depth of air holes (const-z), and arbitrary permittivity at any point in the designable region (freeform).The switching efficiency is plotted against (a) thickness of the designable region for unit cell period kΛ = 10π and diffraction order m = 1 and (b) diffraction order for unit cell period kΛ = 20π and and thickness kL = 3.The switching efficiency is the sum of the power into diffraction orders m and −m as illustrated at the top, where the power diffracted into each state is normalized to 1 (so that the maximum switching efficiency is 2).The freeform designs are not shown in (b) because they are identical to const-z optimization in terms of switching efficiency at the (subwavelength) thickness of kL = 3.