Limits to the Optical Response of Graphene and 2D Materials

TL;DR

This paper establishes fundamental upper limits on the optical responses of 2D materials like graphene, guiding the design of advanced nanophotonic devices by identifying optimal materials and configurations.

Contribution

It derives general bounds on light--matter interactions for 2D materials based solely on optical conductivity, independent of shape or structure.

Findings

Highly doped graphene is optimal in the infrared.

Single-atomic-layer silver is optimal in the visible.

Bounds can be approached in some geometries, indicating potential for design improvements.

Abstract

2D materials provide a platform for strong light--matter interactions, creating wide-ranging design opportunities via new-material discoveries and new methods for geometrical structuring. We derive general upper bounds to the strength of such light--matter interactions, given only the optical conductivity of the material, including spatial nonlocality, and otherwise independent of shape and configuration. Our material figure of merit shows that highly doped graphene is an optimal material at infrared frequencies, whereas single-atomic-layer silver is optimal in the visible. For quantities ranging from absorption and scattering to near-field spontaneous-emission enhancements and radiative heat transfer, we consider canonical geometrical structures and show that in certain cases the bounds can be approached, while in others there may be significant opportunity for design improvement. The bounds can encourage systematic improvements in the design of ultrathin broadband absorbers, 2D antennas, and near-field energy harvesters.

Equations

$$P_{\mathrm{abs}}$$

$$P_{\mathrm{ext}}$$

$$f(\mathbf{E})$$

$$P_{\mathrm{abs}}(\mathbf{E})\le P_{\mathrm{ext}}(\mathbf{E}).$$

$$L\mathbf{K}=\mathbf{E},$$

$$P_{\alpha }\le \frac{1}{2}{\beta }_{\alpha }{\int }_A{\mathbf{E}}_{\mathrm{inc}}^{\ast }\cdot {\left(\mathrm{Re}L\right)}^{-1}{\mathbf{E}}_{\mathrm{inc}}\mathrm{d}A,$$

$${\beta }_{\alpha }=\{\begin{array}{ll}1,& \alpha =\mathrm{absorption}\mathrm{or}\mathrm{extinction}\\ \frac{1}{4},& \alpha =\mathrm{scattering},\end{array}$$

$$\frac{{\sigma }_{\alpha }}{A}\le {\beta }_{\alpha }{Z}_0{\Vert {\mathbf{\sigma }}^{†}{\left(\mathrm{Re}\mathbf{\sigma }\right)}^{-1}\mathbf{\sigma }\Vert }_2$$

$$\frac{{\rho }_{\alpha }}{{\rho }_0}\le \frac{3{\beta }_{\alpha }}{8{\left(k_{0}d\right)}^4}{Z}_0{\Vert {\mathbf{\sigma }}^{†}{\left(\mathrm{Re}\mathbf{\sigma }\right)}^{-1}\mathbf{\sigma }\Vert }_2$$

$$\frac{\Phi }{{\Phi }_{\mathrm{BB}}}\le \frac{3}{2{\left(k_{0}d\right)}^4}{Z}_0^2{\Vert {\mathbf{\sigma }}^{†}{\left(\mathrm{Re}\mathbf{\sigma }\right)}^{-1}\mathbf{\sigma }\Vert }_2^2,$$

$$\displaystyle Z_{0}\left\|\bm{\sigma}^{\dagger}\left(\operatorname{Re}\bm{\sigma}\right)^{-1}\bm{\sigma}\right\|_{2}\approx Z_{0}\sigma\approx\begin{cases}4\alpha\big{(}\frac{E_{\rm F}}{\hbar\gamma}\big{)}&\omega\ll\gamma\\ \pi\alpha&\omega\gg 2E_{\rm F}/\hbar.\end{cases}$$

$$\underset{L}{\underbrace{\left[\frac{-i}{{ϵ}_0\omega {\omega }_{\mathrm{p}}^2}\left({\beta }^2+D\left(\gamma -i\omega \right)\right)\nabla \nabla \cdot +{\sigma }_{\mathrm{loc}}^{-1}\right]}}\mathbf{K}=\mathbf{E},$$

$$\frac{{\sigma }_{\mathrm{ext}}}{A}\le {\left[{\left(Z_0\frac{{\left|{\sigma }_{\mathrm{loc}}\right|}^2}{\mathrm{Re}{\sigma }_{\mathrm{loc}}}\right)}^{-1}+{\left(\frac{r^2}{{\ell }_D^2}\right)}^{-1}\right]}^{-1}.$$

Full text

2D Materials, Graphene, Upper bounds, Near-field optics, Nonlocality

2D materials Novoselov et al. (2005); Geim and Novoselov (2007) and emerging methods Pang et al. (2009); He et al. (2010); Thongrattanasiri et al. (2012); Zhan et al. (2012); Piper and Fan (2014); Cai et al. (2015) for patterning 2D layers and their surroundings are opening an expansive design space, exhibiting significantly different optical Koppens et al. (2011); Basov et al. (2016); Low et al. (2016) (and electronic) properties from their 3D counterparts. In this Letter, we identify energy constraints embedded within Maxwell’s equations that impose theoretical bounds on the largest optical response that can be generated in any 2D material, in the near or far field. The bounds account for material loss as encoded in the real part of a material’s conductivity—in the case of a spatially local conductivity tensor $\bm{\sigma}$, they are proportional to \big{\|}\bm{\sigma}^{\dagger}\left(\operatorname{Re}\bm{\sigma}\right)^{-1}\bm{\sigma}\big{\|}—and are otherwise independent of shape and configuration. We derive the bounds through convex constraints imposed by the optical theorem Newton (1976); Jackson (1999); Lytle et al. (2005) and its near-field analogue, leveraging a recent approach we developed for spatially local 3D materials Miller et al. (2016). In addition to accommodating nonlocal models, this work demonstrates starkly different near-field dependencies of 2D and 3D materials. For graphene, the 2D material of foremost interest to date, the bounds bifurcate into distinctive low- and high-energy regimes: the low-energy bounds are proportional to the Fermi level, whereas the high-energy bounds are proportional to the fine-structure constant, $\alpha$, for any geometrical configuration. We find that far-field bounds on the extinction cross-section can be approached by elliptical graphene disks, whereas the near-field bounds on the local density of states Novotny and Hecht (2012); Joulain et al. (2003); Martin and Piller (1998); D’Aguanno et al. (2004); Oskooi and Johnson (2013) and radiative heat transfer rate Polder and Van Hove (1971); Rytov et al. (1988); Pendry (1999); Mulet et al. (2002); Joulain et al. (2005); Volokitin and Persson (2007) cannot be approached in prototypical flat-sheet configurations. The bounds presented here provide a simple material figure of merit to evaluate the emerging zoo of 2D materials, and offer the prospect of greater optical response via computational design. The material figure of merit can guide ongoing efforts in 2D-material discovery, while the general bounds can shape and drive efforts towards new levels of performance and better optical components.

Plasmonics in 2D materials opens the possibility for stronger light–matter interactions, which may be useful for technological applications, including single-molecule imaging Nie and Emory (1997); Kneipp et al. (1997); van Zanten et al. (2009); Schermelleh et al. (2010) and photovoltaics Atwater and Polman (2010); Ilic et al. (2012a), as well as for basic-science discoveries, such as revealing forbidden transitions Rivera et al. (2016), and achieving unity optical absorption in graphene through optical impedance matching Thongrattanasiri et al. (2012); Jang et al. (2014); Zhu et al. (2016); Kim et al. (2017). Theoretical work towards understanding optical response in 2D materials has focused on analytical expressions using specific geometrical Koppens et al. (2011); Thongrattanasiri et al. (2012); de Abajo (2014) or metamaterial-based Tassin et al. (2012) models, but from a design perspective such assumptions are restrictive. Quasistatic sum rules can yield upper limits on the cross-section de Abajo and Manjavacas (2015); Miller et al. (2014), but have been restricted to far-field quantities and isotropic and spatially local materials. A well-known microwave-engineering bound, known as Rozanov’s theorem Rozanov (2000), offers a bandwidth limit as a function of material thickness, but its contour-integral approach requires perfectly conducting boundaries that are not applicable for 2D materials at optical frequencies. Here, we find constraints that do yield 2D-material optical-response bounds given only the material properties. We provide a general framework to derive limits to any optical-response quantity (including cross-sections, spontaneous-emission enhancements, and radiative-heat exchange), and we present computational results suggesting pathways to approach the new bounds. For a broad class of hydrodynamic nonlocal-conductivity models Ciracì et al. (2012); Mortensen et al. (2014), which capture several important nonclassical features at length scales approaching the quantum regime, we derive general bounds in terms of a constitutive-relation operator. We show that the nonlocal response is necessarily bounded above by the local-response bounds; further, by exploiting the quasistatic nature of interactions at nonlocal length scales, we show that the maximum response must be reduced in proportion to a ratio of the scatterer size to the effective “diffusion” length.

To derive general scattering bounds, consider a 2D scatterer embedded in a possibly heterogeneous background. Passivity, which implies the absence of gain and that polarization currents do no work Welters et al. (2014), requires that the powers absorbed ($P_{\rm abs}$) and scattered ($P_{\rm scat}$) by the target body are non-negative Miller et al. (2016). These almost tautological conditions in fact dictate bounds on the largest currents that can be excited at the surface of any 2D structure. The key is that their sum, extinction ($P_{\rm ext}=P_{\rm abs}+P_{\rm scat}$), is given by the imaginary part of a forward-scattering amplitude, which is a well-known consequence of the optical theorem Newton (1976); Jackson (1999); Lytle et al. (2005). For an arbitrarily shaped 2D scatterer with area $A$ that supports electric surface currents $\mathbf{K}$ (a magnetic-current generalization is given in the Supp. Info.), the absorbed and extinguished powers are given by Jackson (1999); Reid and Johnson (2015)

[TABLE]

where, in the latter expression, $\int_{A}\mathbf{E}^{*}_{\rm inc}\cdot\mathbf{K}\,{\rm d}A$ is a forward-scattering amplitude. A key feature of the optical theorem is that the extinction is the real part of an amplitude, which is linear in the induced currents. By contrast, absorption is a quadratic function of the currents/fields. Yet extinction must be greater than absorption (due to the $P_{\rm scat}\geq 0$ condition noted above), requiring the linear functional to be greater than the quadratic one, a condition that cannot be satisfied for large enough currents. The inequality $P_{\rm abs}\leq P_{\rm ext}$ thereby provides a convex constraint for any optical-response function. Any optical-response maximization can thus be formulated as an optimization problem subject to this convex passivity constraint Miller et al. (2016). For a generic figure of merit $f(\mathbf{E})$ of the fields (or, equivalently, currents), the design problem can be written

[TABLE]

Thanks to the convex nature of the constraint $P_{\rm abs}\leq P_{\rm ext}$ and the simple expressions for $P_{\rm abs}$ and $P_{\rm ext}$, Equation (2) can often be solved analytically—unlike the highly nonconvex Maxwell equations—thereby providing general upper-bound expressions without approximation.

To find bounds that solve Eq. (2), we must specify a relationship between the field $\mathbf{E}$ and the induced current $\mathbf{K}$. To maintain generality we assume only that they are related by a linear operator $\operatorname{\mathcal{L}}$,

[TABLE]

where in different size, material, and parameter regimes, $\operatorname{\mathcal{L}}$ may represent anything from a density-functional-theory operator Yamamoto et al. (2006) or a hydrodynamic model Mortensen et al. (2014); Raza et al. (2015), to a simple scalar conductivity. For a scalar conductivity $\sigma$, $\operatorname{\mathcal{L}}=1/\sigma$. Given this current–field relation, the quadratic dependence of absorption on induced current, per Eq. (1a), is made explicit: $P_{\rm abs}=(1/2)\operatorname{Re}\int_{A}\mathbf{K}^{*}\operatorname{\mathcal{L}}\mathbf{K}\,{\rm d}A$. If we choose the figure of merit to be the absorbed or scattered power, then straightforward variational calculus (see Supp. Info.) from Eq. (2) yields the bounds

[TABLE]

where $\alpha$ denotes absorption, scattering, or extinction. The variable $\beta$ takes the values

[TABLE]

which represent a power-balance asymmetry: absorption and extinction are maximized when $P_{\rm abs}=P_{\rm ext}$, whereas scattering is maximized when $P_{\rm scat}=P_{\rm abs}=P_{\rm ext}/2$, akin to conjugate-matching conditions in circuit theory Stutzman and Thiele (2012). Equation (4) sets a general bound, at any frequency, given only the incident field and the (material-driven) field–current relationship, dictated by the operator $\operatorname{\mathcal{L}}$. The bounds apply in the far field, where $\mathbf{E}_{\rm inc}$ might be a plane wave or Bessel beam, as well as the near field, where $\mathbf{E}_{\rm inc}$ might be the field emanating from dipolar sources. Further below, we show that $\left(\operatorname{Re}\operatorname{\mathcal{L}}\right)^{-1}$ can be considerably simplified in the case when $\operatorname{\mathcal{L}}$ is the differential operator arising in nonlocal hydrodynamic models. First, however, we simplify Eq. (4) for the important case of a spatially local conductivity.

A local conductivity $\bm{\sigma}$, relating currents at any point on the surface to fields at the same point, by $\mathbf{K}=\bm{\sigma}\mathbf{E}$, is the primary response model employed in the study of optical and plasmonic phenomena, in two as well as three dimensions. In 2D materials, it is common to have off-diagonal contributions to the conductivity (e.g. through magnetic-field biasing), and thus we allow $\bm{\sigma}$ to be a general $2\times 2$ matrix (implicitly restricting $\mathbf{E}$ to its two components locally tangential to the 2D surface). Given that $\operatorname{\mathcal{L}}=\bm{\sigma}^{-1}$, the term involving $\operatorname{\mathcal{L}}$ in the bound of Eq. (4) can be written: $\left(\operatorname{Re}\operatorname{\mathcal{L}}\right)^{-1}=\bm{\sigma}^{\dagger}\left(\operatorname{Re}\bm{\sigma}\right)^{-1}\bm{\sigma}$. In far-field scattering, the quantity of interest is typically not the total absorbed or scattered power, but rather the cross-section, defined as the ratio of the power to the average beam intensity. The scattering cross-section, for example, is given by $\sigma_{\rm scat}=P_{\rm scat}/I_{\rm inc}$, where $I_{\rm inc}=|\mathbf{E}_{\rm inc}|^{2}_{\rm avg}/2Z_{0}$. Then, the bound of Eq. (4) simplifies for the absorption, scattering, and extinction cross-sections to

[TABLE]

where $Z_{0}$ is the impedance of free space, $\beta_{\alpha}$ is defined above in Eq. (5), and $\|\cdot\|_{2}$ denotes the induced matrix 2-norm Trefethen and Bau (1997) (which is the largest singular value of the matrix). The power of Eq. (6) is its simplicity—the scattering efficiency of any 2D scatterer, whether it is a periodic array of circles Thongrattanasiri et al. (2012), a spherical coating Christensen et al. (2015), an isolated strip de Abajo (2014), or in any other configuration, has an upper bound defined solely by its material conductivity. We show below that simple ellipses can approach within $\approx 10\%$ of the bounds, and that structures with two additional degrees of freedom can approach within $<1\%$ of the bounds.

A key feature of the approach outlined here is that the optical response of a 2D material of interest can be cleanly delineated (without approximation) from the response of any “background” structures. Our formulation relies on the passivity constraints $P_{\rm scat},P_{\rm abs}>0$, and yet the choice of “incident” and “scattered” fields is arbitrary, as long as they sum to the total fields. As an example, there is significant interest in integrating 2D materials with photonic crystals Zhan et al. (2012); Piper and Fan (2014); we can define the incident field that controls the bounds in Eqs. (4,6) as the field in the presence of only the photonic crystal, and the scattered field as arising only from the addition of the 2D layer above it. The limits of Eqs. (4,6), as well as the limits derived below, then capture the maximum achievable enhancement due to the 2D material itself, subject to its inhomogeneous environment. Throughout this Letter we focus on free-standing graphene to understand its unique optical response, noting that generalization involving substrates and more complex surrounding structures can follow precisely this prescription.

Near-field optical response, in the presence of nearby emitters, is at least as important as far-field response. Here we find bounds to two important near-field quantities: (i) the local density of states (LDOS), which is a measure of the spontaneous-emission rate of a single excited dipole near the scatterer, and (ii) near-field radiative heat transfer, which is a measure of the radiation exchange between two bodies at different temperatures. The (electric) LDOS at a point $\mathbf{x}$ is proportional to the power radiated by an (orientation-averaged) electric dipole at that point, and is given by the expression $\rho=(1/\pi\omega)\operatorname{Im}\sum_{j}\mathbf{p}_{j}\cdot\mathbf{E}_{j}(\mathbf{x})$, where $\mathbf{E}_{j}$ is the electric field excited by the dipole with moment $\mathbf{p}_{j}$, and where the sum over $j=x,y,z$ accounts for orientation-averaging Novotny and Hecht (2012). The expression for $\rho$ shows that LDOS is dictated by a causal amplitude (not a squared amplitude), exhibiting similar mathematical structure to extinction. The source of the similarity is that both extinction and LDOS can be decomposed into radiative and nonradiative components, which for the LDOS we denote by $\rho_{\rm rad}$ and $\rho_{\rm nr}$, respectively. The nonradiative part of the LDOS is given by the absorption in the scattering body (which is often an antenna), and per Eq. (1a) is quadratic in the induced currents. Unlike far-field scattering, in the near field, the incident field increases rapidly at smaller distances $d$ ($|\mathbf{E}|\sim 1/d^{3}$). Thus, the same convex-optimization problem laid out in Eq. (2) leads to distance-dependent LDOS bounds via the replacements $P_{\rm ext}\rightarrow\rho$ and $P_{\rm abs}\rightarrow\rho_{\rm nr}$. For an arbitrarily shaped 2D surface separated from the emitter by some minimum distance $d$, the bounds are (Supp. Info.):

[TABLE]

where $\alpha$ in this context denotes the total, radiative, or nonradiative component of the LDOS, $k_{0}=\omega/c$, and $\rho_{0}$ is the free-space electric-dipole LDOS, $\rho_{0}=\omega^{2}/2\pi^{2}c^{3}$. Again $\beta_{\alpha}$ represents a power-balance (conjugate-matching) condition, and takes the value $1$ for nonradiative or total LDOS and $1/4$ for the radiative LDOS. Equation (7) includes the highest-order ($\sim 1/d^{3}$) term from the incident electric field; lower-order terms ($\sim 1/d^{2},1/d$) are generally negligible in the high-enhancement regimes of interest, as discussed quantitatively in Miller et al. (2016). The $3/8$ coefficient in Eq. (7) is for the common case in which the surface is separated from the emitter by a separating plane; if the scattering body surrounds the emitter across a solid angle $\Omega$, the bound in Eq. (7) is multiplied by $4\Omega$. Equation (7) provides a general answer to the question of how efficient and effective a 2D optical antenna can be.

Radiative heat transfer (RHT), in which a warm body transfers energy to a colder one via photon exchange, is also subject to optical-response bounds. It has long been known Polder and Van Hove (1971); Rytov et al. (1988); Pendry (1999) that near-field RHT can surpass the blackbody limit, as evanescent tunneling can outpace radiative exchange. Yet general limits to the process in conventional 3D materials had been unknown until our recent work Miller et al. (2015). The total RHT rate, $H$, is given by the net flux from one body at temperature $T_{1}$ to another at temperature $T_{2}$, typically expressed as (Joulain et al. (2005)) $H_{1\rightarrow 2}=\int_{0}^{\infty}\Phi(\omega)\left[\Theta(\omega,T_{1})-\Theta(\omega,T_{2})\right]\mathrm{d}\omega$, where $\Phi(\omega)$ is a temperature-independent energy flux and $\Theta$ is the Planck spectrum. The flux $\Phi$ is the power absorbed by the second body, having been emitted from the first, such that it is similar to the scattering problem bounded by Eq. (6). A key distinction is that the (incoherent) sources are in the interior of one of the scattering bodies, invalidating the conventional optical theorem. This difficulty can be circumvented by breaking the flux transfer into two scattering problems, connected by a generalized Kong (1972) reciprocity relation (the material conductivity does not need to be reciprocal), as outlined in Miller et al. (2015). The key distinction in the case of 2D materials is the dimensionality of the domain over which the field intensities are evaluated, which for bodies with identical conductivities $\bm{\sigma}$ leads to the bound

[TABLE]

where $d$ is the minimum separation distance between the arbitrarily shaped bodies, $\Phi_{\rm BB}=k_{0}^{2}A/4\pi^{2}$ is the blackbody limit (for infinite area $A$) Joulain et al. (2005), and the conductivity term is squared due to potential contributions from each body (see Supp. Info.). As for the LDOS bounds, Eq. (8) assumes a separating plane between the bodies; corrugated surfaces that are interlaced (but non-touching) have bounds of the same functional form but with different numerical prefactors. An interesting 2D-specific aspect of Eqs. (7,8) is that they exhibit identical $1/d^{4}$ distance dependencies, whereas for 3D bodies, RHT increases more slowly for smaller separations ($\sim 1/d^{2}$) than does the LDOS ($\sim 1/d^{3}$).

The fundamental limits of Eqs. (6–8) share a common dimensionless material “figure of merit” (FOM), Z_{0}\big{\|}\bm{\sigma}^{\dagger}\left(\operatorname{Re}\bm{\sigma}\right)^{-1}\bm{\sigma}\big{\|}_{2}. The FOM, which simplifies to $Z_{0}|\sigma|^{2}/\operatorname{Re}\sigma$ for a scalar conductivity, captures the intrinsic tradeoffs between high conductivity for large response and high losses that dissipate enhancement, and can be used to identify optimal materials. In Fig. 1 we plot the FOM across a range of frequencies, using experimentally measured or analytically modeled material data for common 2D materials of interest: graphene, for various Fermi levels Jablan et al. (2009), magnetic biasing Hanson (2008), and AA-type bilayer stacking Wang et al. (2016) (at $300\text{\,}\mathrm{K}$), hBN Brar et al. (2014), \ceMoS_2 Liu et al. (2014a), black phosphourous (BP) Low et al. (2016), \ceBi2Se3 (at THz frequencies Di Pietro et al. (2013)), and metals Ag, Al, and Au, all taken to have 2D conductivities dictated by a combination de Abajo and Manjavacas (2015) of their bulk properties and their interlayer atomic spacing. Strongly doped graphene ($E_{\rm F}=$0.6\text{,}\mathrm{eV}$$) offers the largest possible response across the infrared, whereas 2D Ag tends to be better in the visible. At THz frequencies, where graphene’s potential is well-understood Rana (2008); Ju et al. (2011); Low and Avouris (2014), the topological insulator \ceBi2Se3 shows promise for even larger response. More broadly, the simple material FOM, $|\sigma|^{2}/\operatorname{Re}\sigma$ or its anistropic generalization $\left\|\bm{\sigma}^{\dagger}\left(\operatorname{Re}\bm{\sigma}\right)^{-1}\bm{\sigma}\right\|$, offers a metric for evaluating emerging (e.g. silicene Vogt et al. (2012), phosphorene Xia et al. (2014); Liu et al. (2014b)) and yet-to-be-discovered 2D materials.

In the following we specialize our considerations to graphene, the standard-bearer for 2D materials, to examine the degree to which the bounds of Eqs. (6–8) can be attained in specific structures. We adopt the conventional local description, including intra- and interband dispersion. Appropriate modifications de Abajo (2014); Jablan et al. (2009) are included to account for a finite intrinsic damping rate, \gamma=1/\tau=\left({10}^{12}\text{,}\mathrm{eV}\mathrm{/}\mathrm{s}\right)/E_{\rm F}, which is taken as Fermi-level-dependent (corresponding to a Fermi-level-independent mobility), with a magnitude mirroring that adopted in de Abajo (2014). Figure 2 shows the cross-section bounds (dashed lines), per Eq. (6), alongside graphene disks (with $E_{\rm F}=$0.4\text{,}\mathrm{e}\mathrm{V}$$) that approach the bounds at frequencies across the infrared. For simplicity, we fix the aspect ratio of the disks at 2:1 and simply reduce their size to increase their resonant frequency; each disk achieves $\approx 85\%$ of its extinction cross-section bound. The disk diameters are kept greater than $10\text{\,}\mathrm{n}\mathrm{m}$ to ensure the validity of our local description. We employ a fast quasistatic solver Christensen (2017) to compute the response of the ellipses, which is verified with a free-software implementation Reid of the boundary element method (BEM) Harrington (1993) for the full electrodynamic problem with the surface conductivity incorporated as a modified boundary condition Reid and Johnson (2015). If edge scattering, or any other defect, were to increase the damping rate, such an increase could be seamlessly incorporated in the bounds of Eqs. (6–8) through direct modification of the conductivity. In the Supp. Info., we show that with two extra geometrical degrees of freedom (e.g., a “pinched ellipse”), one can reach $>99.6\%$ of the bound. The cross-section bounds can also be used as bounds on the fill fraction of graphene required for perfect absorption in a planar arrangement, and they suggest the potential for an order-of-magnitude reduction relative to the best known results Thongrattanasiri et al. (2012). Conversely, such room for improvement could be used to significantly increase the perfect-absorption bandwidth beyond the modern state-of-the-art.

The bounds simplify analytically at the low- and high-frequency extremes. In these regimes, graphene’s isotropic conductivity is real-valued and comprises simple material and fundamental constants, such that the material FOM is approximately

[TABLE]

The low-frequency proportionality to $E_{\rm F}/\hbar\gamma$ arises as a consequence of the intraband contributions to the conductivity, in contrast to the interband dominance at high frequencies. Interband contributions to the conductivity are often ignored at energies below the Fermi level, but even at those energies they are responsible for a sizable fraction of the loss rate, thus causing the quadratic roll-off (derived in Supp. Info.) of the maximum efficiency seen on the left-hand side of Fig. 2.

Famously, at high frequencies a uniform sheet of graphene has a scattering efficiency $\sigma/A\approx\pi\alpha$ (Refs. Nair et al. (2008); Kuzmenko et al. (2008); Mak et al. (2008)). Interestingly, Fig. 2 and Eq. (9) reveal that $\pi\alpha$ is the largest possible scattering efficiency, for any shape or configuration of graphene, at those frequencies. Per the incident-field discussion above, it is possible to increase the absolute absorption of a plane wave at those frequencies by structuring the background (e.g. with a photonic-crystal slab supporting the graphene), but the percentage of the background field intensity that can be absorbed by the graphene is necessarily $\leq\pi\alpha$, no matter how the graphene is structured. The right-hand side of Fig. 2 shows the bounds for each Fermi level converging to $\pi\alpha$, with the inset magnifying the high-energy region and showing that the response of a flat sheet indeed reaches the bound.

The near-field LDOS and RHT limits are more challenging to attain. We study the LDOS near a flat sheet of graphene, the most common 2D platform for spontaneous-emission enhancements to date Koppens et al. (2011); Gaudreau et al. (2013); Tielrooij et al. (2015), and show that there is a large performance gap between the flat-sheet response and the fundamental limits of Eq. (7). There are two key factors that control the near-field bounds (for both LDOS and RHT): the material FOM $|\sigma|^{2}/\operatorname{Re}\sigma$, and a “near-field enhancement factor” $1/d^{4}$, for emitter–sheet distance $d$. The $1/d^{4}$ near-field enhancement factor is particularly interesting, because it increases more rapidly than in 3D materials (for which the LDOS Miller et al. (2016) and RHT Miller et al. (2015) bounds scale as $1/d^{3}$ and $1/d^{2}$, resp.). In Fig. 3, we show the LDOS as a function of the emitter–graphene separation, for a fixed Fermi level $E_{\rm F}=$0.4\text{,}\mathrm{eV} and a range of frequencies (colored solid lines). The bounds for each frequency are shown in the colored dashed lines, and the ratio of the LDOS $\rho$ to the LDOS bound $\rho_{\rm bound}$ is shown in the inset. For low and moderate frequencies, there is an ideal distance at which the LDOS most closely approaches its frequency-dependent bound, whereas the high-frequency regime (e.g. $\hbar\omega=$0.7\text{\,}\mathrm{e}\mathrm{V}) is almost distance-insensitive due to high losses.

Figure 3 shows two asymptotic distance-scaling trends. First, at high frequencies and/or large separations ($50\text{\,}\mathrm{nm}$ to $1\text{\,}\mathrm{\SIUnitSymbolMicro m}$), the LDOS enhancement scales as $1/(k_{0}d)^{4}$. We show in the Supp. Info. that in this regime the LDOS further exhibits the material-enhancement factor $|\sigma|^{2}/\operatorname{Re}\sigma$, falling short of the bound only by a factor of 2. In this regime, the LDOS is dominated by a “lossy-background” contribution Gaudreau et al. (2013), which is insensitive to details of the plasmonic mode, and due instead predominantly to interband absorption in graphene (permitted even below $2E_{\rm F}$ for nonzero temperatures). Of more interest may be the opposite regime—higher frequencies at smaller separations—which are known Rodriguez-López et al. (2015) to have reduced distance dependencies. It is crucial to note that the bounds presented in this Letter are not scaling laws; instead, at each frequency and distance they represent independent response limits. We see in Fig. 3 that for each individual frequency, $\rho/\rho_{0}$ flattens towards a constant value at very small distances, because the corresponding plasmon surface-parallel wavenumber is smaller than $1/d$ and does not change; however, the envelope formed over many frequencies (for a given separation $d$) shows a $1/(k_{0}d)^{3}$ as higher-wavenumber plasmons are accessed at smaller distances. This suggests a simple potential approach to reach the bound: instead of finding a geometrical configuration that approaches the bound at all frequencies and separations, concentrate on finding a structure that reaches the bound at a single frequency and separation of interest. A “family” of structures that combine to approach the bounds over a large parameter regime may then naturally emerge.

Near-field RHT shows similar characteristics, in which the bounds may be approached with flat graphene sheets at specific energy, Fermi-level, and separation-distance parameter combinations. As a counterpart to the LDOS representation of Fig. 3, in Fig. 4 we fix the separation distance at $10\text{\,}\mathrm{nm}$ and plot the frequency-dependent RHT Ilic et al. (2012b) for three Fermi levels. The respective bounds, from Eq. (8), show the same “dip” as seen in the inset of Fig. 2(b), which occurs at the frequency where the imaginary part of the conductivity crosses zero. At these frequencies, the RHT between flat sheets can approach the bounds. However, at other frequencies, where the potential RHT is significantly larger, the flat sheets fall short by orders of magnitude, as depicted in Fig. 4 at $E_{\rm F}=$0.5\text{,}\mathrm{eV}$$. The flat-sheet case falls short due to near-field interference effects: as the sheets approach each other, the plasmonic modes at each interface interact with each other, creating a level-splitting effect that reduces their maximum transmission to only a narrow band of wavevectors Miller et al. (2016). By contrast, for two dipolar circles in a quasistatic approximation (Fig. 4 inset), the RHT between the two bodies can approach its respective bound. These examples suggest that patterned graphene sheets, designed to control and optimize their two-body interference patterns, represent a promising approach towards reaching the bounds and thereby unprecedented levels of radiative heat transfer. In the Supp. Info., we show that achieving RHT at the level of the bound, even over the narrow bandwidths associated with plasmonic resonances, would enable radiative transfer to be greater than conductive transfer through air at separations of almost $1\text{\,}\mu\mathrm{m}$, significantly larger than is currently possible Miller et al. (2016).

Having examined the response of graphene structures in the local-conductivity approximation, we now reconsider nonlocal conductivity models. For structures in the $210\text{\,}\mathrm{nm}$ size range, below the local-conductivity regime but large enough to not necessitate fully quantum-mechanical models, hydrodynamic conductivity equations Ciracì et al. (2012); Mortensen et al. (2014); Raza et al. (2015), or similar gradient-based models of nonlocality Hanson (2008); Fallahi et al. (2015), can provide an improved account of the optical response. In a hydrodynamic model, the currents behave akin to fluids with a diffusion constant $D$ and convection constant $\beta$ (both real-valued), with a current–field relation given by Mortensen et al. (2014)

[TABLE]

where $\sigma_{\rm loc}$, $\omega_{\rm p}$, and $\gamma$ are the local conductivity, plasma frequency, and damping rate of the 2D material, respectively. Per Eq. (4), the 2D-material response bounds depend only on the Hermitian part of the $\operatorname{\mathcal{L}}$ operator, denoted by an underbrace in Eq. (10). Before deriving bounds dependent on the hydrodynamic parameters, we note that the grad–div hydrodynamic term in Eq. (10) cannot increase the maximum optical response. The operator $-\nabla\nabla\cdot$ is a positive semidefinite Hermitian operator (for the usual $L^{2}$-space overlap-integral inner product), shown by integration by parts in conjunction with the no-spillout boundary condition. The Hermiticity of the grad–div operator means that the Hermitian part of $\operatorname{\mathcal{L}}$ is given by $\operatorname{Re}\operatorname{\mathcal{L}}=(\operatorname{\mathcal{L}}+\operatorname{\mathcal{L}}^{\dagger})/2=-\frac{D}{\omega_{\rm p}^{2}}\nabla\nabla\cdot+\operatorname{Re}\sigma_{\rm loc}^{-1}$. Because $-\nabla\nabla\cdot$ is a positive-semidefinite addition to the positive-semidefinite term $\operatorname{Re}\sigma_{\rm loc}^{-1}$, $\left\|\left(\operatorname{Re}\operatorname{\mathcal{L}}\right)^{-1}\right\|\leq\left(\operatorname{Re}\sigma_{\rm loc}^{-1}\right)^{-1}=|\sigma_{\rm loc}|^{2}/\operatorname{Re}\sigma_{\rm loc}$. Thus the nonlocal response is subject to the bound imposed by the underlying local conductivity, demonstrating that nonlocal effects of this type cannot surpass the local-conductivity response explored in depth above.

We can further show that hydrodynamic nonlocality necessarily reduces the maximum achievable optical response in a given 2D material, by exploiting the quasistatic nature of electromagnetic interactions at the length scales for which nonlocal effects manifest. The key insight required to derive bounds subject to the nonlocal current–field relation, Eq. (10), is that the absorbed power can be written as a quadratic form of both the currents $\mathbf{K}$ as well as $\nabla\cdot\mathbf{K}$ (proportional to the induced charge): $P_{\rm abs}=(1/2)\operatorname{Re}\int_{A}\mathbf{K}^{*}\cdot\mathbf{E}=1/2\int_{A}\left[a(\nabla\cdot\mathbf{K}^{*})(\nabla\cdot\mathbf{K})+b\mathbf{K}^{*}\cdot\mathbf{K}\right]$, where $a=D/\omega_{\rm p}^{2}$ and $b=\operatorname{Re}(\sigma_{\rm loc}^{-1})$. Similarly, the extinction can be written as a linear function of either $\mathbf{K}$ or $\nabla\cdot\mathbf{K}$ (exploiting the quasistatic nature of the fields), such that $P_{\rm abs}\leq P_{\rm ext}$ offers two convex constraints for the generalized nonlocal-conductivity problem. We defer to the Supp. Info. for a detailed derivation of general figures of merit under this constraint, and state a simplified version of the result for the extinction cross-section. The additional $\nabla\cdot\mathbf{K}$ constraint introduces a size dependence in the bound, in the form of a “radius” $r$ that is the smallest bounding sphere of the scatterer along the direction of the incident-field polarization. Defining a plasmonic “diffusion” length $\ell_{D}=\sqrt{cD/\omega_{\rm p}^{2}}$ (for speed of light $c$), the variational-calculus approach outlined above yields an analogue of Eq. (6) in the presence of a hydrodynamic nonlocality:

[TABLE]

Equation (11) has an appealing, intuitive interpretation: the cross-section of a scatterer is bounded above by a combination of the local-conductivity bound and a nonlocal contribution proportional to the square of the ratio of the size of the scatterer to the “diffusion” length. Thus as the size of the particle approaches $\ell_{D}$, and goes below it, there must be a significant reduction in the maximal attainable optical response. There is ambiguity as to what the exact value of $D$, or equivalently $\ell_{D}$, should be in 2D materials such as graphene; the bounds developed serve as an impetus for future measurement or simulation, to delineate the sizes at which the local/non-local transition occurs. Conversely, since the bound shows a dramatic reduction at sizes below $\ell_{D}$, Eq. (11) can serve as a means to extract this nonlocal property of any 2D material from experimental measurements.

General limits serve to contextualize a large design space, pointing towards phenomena and performance levels that may be possible, and clarifying basic limiting factors. Here we have presented a set of optical-response bounds for 2D materials, generalizing recent 3D-material bounds Miller et al. (2015, 2016) to incorporate both local and nonlocal models of 2D conductivities. We further studied the response of standard graphene structures—ellipses and sheets—relative to their respective bounds, showing that the far-field absorption efficiency bounds can be reliably approached within $10\%$, but that the near-field bounds are approached only in specific parameter regimes, suggesting the possibility for design to enable new levels of response. The figure of merit $\left\|\bm{\sigma}^{\dagger}\left(\operatorname{Re}\bm{\sigma}\right)^{-1}\bm{\sigma}\right\|$ can serve to evaluate new 2D materials as they are discovered, and their optical properties are measured. Our results point to a few directions where future work may further clarify the landscape for 2D-material optics. One topic of current interest is in patterned gain and loss Pick et al. (2017); Cerjan and Fan (2016) (esp. $\mathcal{PT}$-symmetry Guo et al. (2009); Regensburger et al. (2012); Cerjan et al. (2016)), which exhibit a variety of novel behaviors, from exceptional points to loss-induced transparency. Our bounds depend on passivity, which excludes gain materials, but in fact the bounds only require passivity on average, i.e., averaged over the structure. Thus Eqs. (4–8) should be extensible to patterned gain–loss structures. A second area for future work is in exploration of quantum models of the $\operatorname{\mathcal{L}}$ operator. We have shown here explicit bounds for the cases of local and hydrodynamic conductivities, but there is also significant interest in quantum descriptions of the response. Through, for example, density-functional theory Burke et al. (2005), analytical bounds in such cases may lead to a continuum of optical-response limits across classical, semi-classical, and quantum regimes.