{'fraction_non_alphanumeric': 0.05750056793800641, 'fraction_numerical': 0.018451674786076684, 'mean_word_length': 4.337195204095379, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The role of optical loss and tight bounds in polariton-mediated near-field heat transfer

We introduce an analytical framework for near-field radiative heat transfer in bulk plasmonic and polar media. Considering material dispersion, we derive a closed form expression for the radiative thermal conductance, which disentangles the role of optical loss from other material dispersion characteristics, such as the spectral width of the Reststrahlen band in polar dielectrics, as well as from the temperature. We provide a universal condition for maximizing heat transfer that defines the optimal interplay between a material's optical loss and polariton resonance frequency, based on which we introduce tight bounds to near-field heat transfer. With this formalism, one can quantitatively evaluate all polaritonic materials in terms of their performance as near-field thermal emitters.

Accessing the thermal near-field, where the separation distance between objects is smaller than the relevant thermal wavelength, typically on the order of few microns, yields ultra-high radiative heat transfer (on the order of tens of W/cm 2 [1,5,[25][26][27][28]). In this range, radiative heat is optimally transferred by evanescent modes, such as surface plasmon polaritons (SPPs) and surface phonon polaritons (SPhPs) [27,29], that tunnel across a nanometric vacuum gap. SPPs are resonantly excited on a plasmonic metal surface [30], e.g. Au or Ag, at frequencies ranging from near-infrared (IR) to ultraviolet [31]. On the other hand, SPhPs are surface resonant modes occurring on polar dielectrics and semiconductors, e.g. hexagonal BN, SiC or GaAs, at frequencies ranging from mid-IR to < 10 THz [29].
Near-field radiative heat transfer (NFRHT) is rigorously described within the framework of fluctuational electrodynamics [1]. Computing the total NFRHT between bodies requires carrying out a spectral as well as a spatial integration of the exchanged thermal radiation over the entire frequency spectrum and for all relevant wavenumbers, respectively. Although various sophisticated concepts have already been considered [32,33], this integration remains challenging in practice even in the simplest case of NFRHT between semi-infinite planar layers. Hence numerical approaches are usually employed [34,35]. Recent works have analytically estimated upper bounds to the thermal emission spectrum per frequency (termed Φ(ω) henceforth) [14,[36][37][38][39], as well as upper bounds to the total radiative thermal conductance [36,38]. In [40,41], analytical solutions to the total integrated heat transfer for polaritonic media were reported.
Nonetheless, critical material parameters such as the polariton resonance frequency, the optical loss, as well as other dispersion characteristics, have remained mathematically intertwined. This prohibits a deep physical understanding of the role of each in NFRHT. To the authors knowledge, in no previous works have these parameters been disentangled to offer new insights into this otherwise well-studied problem in physics. Building a general framework for the classification of polaritonic materials in terms of NFRHT performance has remained a challenge.
In this Letter, we provide an analytical framework for NFRHT mediated by surface polaritons in planar bulk systems that describe the explicit dependence of NFRHT on material loss, quality factor, and temperature. Our results offer insights into polariton-mediated NFRHT and allow for a quantitative classification of all polaritonic media in terms of NFRHT performance.
We describe NFRHT at a mean temperature, T , by evaluating the radiative thermal conductance per unit area [42]: where θ(ω, T ) = ω/ [exp ( ω/k B T ) − 1] is mean energy per photon and Φ(ω) is the thermal emission spectrum, expressed in m −2 . We consider a vacuum gap of size d that separates two planar semi-infinite bodies exchanging heat. From fluctuational-electrodynamics [1,2,43,44], Φ(ω) is given by: where ξ p,s (ω, β) is the probability for a photon at frequency ω and in-plane wavenumber β to tunnel accross the gap. The subscripts p and s denote polarization, corresponding to TM (transverse magnetic) and TE (transverse electric), respectively. For gap sizes smaller than the thermal wavelength, λ T = b Wien /T , where b Wien = 2989 µm K [45], thermally excited SPPs and SPhPs dominate NFRHT in plasmonic and polar materials, respectively. Since these can only be excited in p-polarization [30], the emission spectrum, Φ, can be approximated as the contribution from ppolarization alone. At sufficiently small vacuum gaps, the dispersion of surface polaritons approaches the quasistatic limit, for which β k 0 [14], where k 0 = ω/c is the free-space wavenumber. In this limit, it can be shown (Supplementary Material) that the maximum in-plane wavenumber, β max , that satisfies the perfect photon tunneling condition, i.e. ξ p = 1, occurs near the surface polariton resonance frequency, Ω, at which Re is the Fresnel coefficient in the quasistatic regime, or: For low-loss materials, i.e. for Im{χ(Ω)} 1, Eq. (3) reduces to the more common expression Re {χ(Ω)} = −2 [30].
To obtain the emission spectrum, we carry out the integration of Eq. (2), over all available wavenumbers, β. Upon assuming β k 0 , this integration yields: where Li 2 is the dilogarithm or Spence's function [46]. This expression agrees with Rousseau et al. [40]. Its derivation, along with the more general expression for heat exchange between dissimilar materials, can be found in the Supplementary Material, where we also showcase the validity of Eq. (4). Eq. (4) directly computes the thermal emission spectrum provided the Fresnel coefficients, and is valid for any material's frequency dispersion. In the low-loss limit, the emission spectrum is maximum at ω = Ω, where Eq. (4) reduces to the result by Miller et al. in [47], as shown in the Supplementary Material, where we discuss an approach to quantitatively distinguish the low-loss from the high-loss regime in NFRHT (Eqs. (S27)-(S28)). The dielectric function of plasmonic and polar materials can be described by Drude and Lorentz oscillators, respectively: where ε ∞ is the high-frequency relative permittivity, and γ is the optical loss [48,49]. For plasmonic metals, ω p is the plasma frequency, near which the SPP mode occurs. For polar materials, ω TO and ω LO are the transverse and longitudinal optical phonon frequencies, respectively [29].
The spectral range [ω TO , ω LO ] defines the Reststrahlen band [50], within which the SPhP mode occurs. The resonance frequency, Ω, is found by solving Eq. (3). Since Ω does not vary significantly as γ increases (Supplementary Material), we evaluate it for γ → 0 as: Henceforth, we assume that Ω plasm and Ω polar are γindependent.
Radiative thermal conductance, h -In evaluating NFRHT performance, it is useful to introduce a material quality factor for plasmonic and polar media [29,51,52]: To obtain the radiative thermal conductance (Eq. (1)), we use contour integration in the complex frequency plane of Φ (section Radiative thermal conductance h in the Supplementary Material), similar to [40]. Considering that the Planck distribution varies slowly with respect to Φ(ω), we obtain: The functions Π and Ψ are given by: These functions are both bounded by unity, hence h max in Eq. (8) defines the maximum thermal conductance that a polaritonic material can reach in a planar configuration, and is given by: As can be seen, h max is temperature-and lossindependent. Furthermore, it expresses the well-known ∝ d −2 dependence of NFRHT from the vacuum gap size [34,36,40,41]. The parameter B is termed material residue henceforth, and it is defined as B = Q/Im{r p (Ω)}, by setting Q → ∞, which reduces to: for Drude and Lorentz materials, respectively. Eq. (8) is the key contribution of this paper. Unlike expressions presented in previous works [40,41], Eq. (8) distinctly separates the role of the optical loss, described by the quality factor Q, in NFRHT, from other dispersion parameters that are captured by B, and temperature. The decoupling of temperature, material quality factor, and material residue, in Eq. Eq. (8) is a very good approximation of the exact result obtained via fluctuational electrodynamics, for Q considerably larger than unity, which is satisfied by all relevant materials for NFRHT [53]. Further, Eq. (8) is exact for ε ∞ = 1 in either polar or plasmonic media. We stress that Eq. (11) represents a tight bound to NFRHT that accounts for material dispersion. This is to be contrasted to previous results that derived upper bounds to h with the idealized assumption of a dispersion-less perfect blackbody in the near-field (i.e. ξ = 1) [36], thus yielding a thermal conductance that is orders of magnitude larger than our result in Eq. (8) (Fig. 3).
With Eqs. (8), (10), one can identify the optimal material characteristics, independent of temperature, that maximize NFRHT. In particular, Ψ describes how NFRHT changes with optical loss. Seeking for the maximum of Ψ (Eq. (10)), we obtain: Hence, NFRHT is maximized when the material quality factor is 4.5 times the material residue function, given in Eq. (12). The work in [36] yielded a universal optimal quality factor, namely Q * opt = 2.72, which, however, is independent of B, thus suggesting that all materials that have the same Q should perform identically in terms of NFRHT. In contrast, Eq. (13) demonstrates that other dispersion characteristics, beyond the quality factor, are critical in evaluating NFRHT response. From Eq. (12a), the material residue for plasmonic materials depends only on ε ∞ . Typically, ε ∞ 10, hence B plasm remains well below 10. Thus, from Eq. (13), Q opt for plasmonic materials is relatively low, namely Q opt 50. Hence, plasmonic materials with good NFRHT performance have high-loss (γ) and modest Q, and NFRHT is enhanced due to the broadband nature of the plasmonic resonance. In contrast, the material residue for polar materials, B polar (Eq. (12b)) is inversely proportional to the spectral width of the Reststrahlen band, (ω LO − ω TO ). The Reststrahlen band of most polar materials is narrow, hence B polar > B plasm , therefore Q opt for polar media is higher than for plasmonic ones. In contrast to plasmonic media, in polar ones, it is the narrowband nature of SPhPs that enhances NFRHT (See Supplementary Material, Eq. (S5)).
In Fig. 1, Eq. (13) is shown with the solid line. We also evaluate the NFRHT performance of several relevant polaritonic emitters considered in literature. These include polar materials such as Silicon Carbide (SiC), hexagonal Boron Nitride (hBN), and doped semiconductors, e.g. Gallium Arsenide (GaAs), Indium Arsenide (InAs) [29,55,56], as well as plasmonic materials such as standard noble metals, e.g. Gold (Au), Silver (Ag), and heavily doped oxides, e.g. IZO and GZO [29,31,57]. The distance between each point in Fig. 1 and the solid curve representing Eq. (13) expresses how far from the ideal material performance each material falls. Interestingly, an ultra-high Q does not necessarily yield optimal NFRHT. By contrast, it is the interplay between Q and B that is critical, making, for instance, GaAs, AZO and GaN near-optimal materials for NFRHT as compared to Ag or 3C-SiC, even though the latter exhibit ultra-high quality factors. This demonstrates the importance of the material residue, B, in evaluating NFRHT performance.
The parameter h max in Eq. (11) is the maximum radiative thermal conductance achievable for each material, if one adjusted their quality factor such that Ψ Q B → 1, and in the limit of infinite temperature, for which Π Ω T → 1. In Fig. 2 (a), we calculate h max for the materials considered in Fig. 1. To quantify the degree to which the loss of each material deviates from the optimal value as defined in Eq. (13), we plot these points against the ratio Q/B, which is inversely proportional to γ, the optical loss. Fig. 2   Au as compared to that of Ag (see Fig. 1). In other words, small deviations of ε ∞ from unity in plasmonic metals (Eq. (12a)), and, similarly, sub-optimal Reststrahlen band spectral widths with respect to Ω in polar materials (Eq. (12b)), can considerably affect the optimal point of NFRHT.
The dependence of NFRHT from the optical loss is captured explicitly in Ψ (Eq. (10)), and is shown graphically in Fig. 2 (b). As described in Eq. (13), Ψ is maximum at Q opt = 4.5 B, depicted with the vertical dashed line. The horizontal distance between this line and each point in Fig. 2 (a) indicates how close each material is to the ideal optical loss, for its particular resonance frequency, Ω. For example, despite the comparable material residue, B, of Ag and Au, the loss (γ) of Au yields a value of Ψ that is much closer to unity as compared to Ag, hence Au presents overall better NFRHT performance, which is consistent with Fig. 1.
In Fig. 2 (b), we also append points that correspond to exact calculations with fluctuational electrodynamics, for few commonly used materials in NFRHT. These calculations are performed in the limit of infinite temperature, for the sake of a meaningful comparison with our formalism in Eq. (8). These exact results, represented as markers, are in very good agreement with our theory (solid line) for all considered materials. Small discrepan-cies occur in the range of relatively low-Q, for example in the case of IZO [31], since our formalism assumes resonant material response, hence its accuracy improves as the material quality factor increases (see Supplementary Material for details).
The temperature dependence of NFRHT is described via Π Ω T in Eq. (9), which is the only temperaturedependent term in Eq. (8), and agrees with previous analytical results [36,40]. In contrast to Wein's displacement law in the far-field, where h scales as T 3 , in the near-field, Π scales as ∼ T −2 . On the other hand, from Eq. (11), h max scales with the resonance frequency, therefore materials supporting polaritons at high frequencies (high-Ω) will, in principle, reach higher NFRHT rates. However, for this to occur in practice, they ought to operate at dramatically higher temperatures. Specifically, since Π decreases exponentially with the ratio Ω/T , to avoid a dramatic damping in h, a higher resonance frequency should be compensated by a higher operating temperature, as expected. This can be seen in Fig. 3, showing the total radiative thermal conductance, h, computed via our analytical result (Eq. (8)), as well as the exact result (fluctuational electrodynamics), where the wavenumber and frequency integrations are carried out numerically. We consider a set of plasmonic materials, i.e. IZO and Ag, and a set of polar ones, i.e. 4H-SiC and AlN. It is clear that plasmonic materials reach higher NFRHT than their polar counterparts, however this occurs at very high temperatures. This is expected since the resonance frequency, Ω, of plasmonic media is significantly higher than that of polar ones. Using Eq. 9, one can estimate the optimal temperature of operation for each material (Fig. S5 in the Supplementary Material).
In Fig. 3, we also append the exact results with fluctuational electrodynamics (dotted). These are in excellent agreement with our analytical formalism, except for small deviations that occur only for materials with relatively low Q. This is expected, since a low-Q suggests a spectrally broadband response, whereas our formalism applies to polaritonic resonances (Supplemental Material, Fig. S4 (a)). To conclude, the vast majority of polar and plasmonic materials, one can compute exactly their NFRHT properties with Eqs. (8)(9)(10)(11).
As a reference, in Fig. 3, we also show the fundamental limit to the radiative thermal conductance h, i.e. h U max = k 2 B T 3 d 2 , as derived by Ben-Abdallah et al. [36] (dashed line), and the upper bound derived by Venkataram et al. [38] for one of the considered materials, viz. 4H-SiC, denoted with h opt (4H-SiC) (orange curve). By comparing with our exact results, both h U max and h opt (4H-SiC) represent loose bounds to NFRHT. This is to be contrasted with the expression in Eq. (11), which is the limit to which h actually saturates at high temperatures for op-  timal loss, i.e. Q = Q opt , for every material (see right y-axis in Fig. 3).
Conclusion -We present a simple analytical framework that describes NFRHT in polaritonic bulk systems. We derive a universal closed form expression (Eq. (8)) for the thermal conductance that is valid for any plasmonic or polar material. This expression clarifies what the role of optical loss (γ) and material quality factor (Q ∝ γ −1 ) are in NFRHT, as well as their interplay with other material dispersion characteristics. We show that the quality factor of a material's polariton resonance alone is not sufficient to accurately describe NFRHT. In contrast, we introduced the material residue parameter, B that completes the analytical framework for the classification of all plasmonic and polar materials for NFRHT.
We derive a material-dependent optimal condition that maximizes NFRHT, namely Q = 4.5 B, where the quality factor of the polaritonic resonance is inversely proportional to the optical loss, and the material residue is loss-independent and encompasses critical properties in polaritonic materials, i.e. the resonance frequency and the spectral width of the Reststrahlen band.
In previous works, upper bounds to the spectral emissivity [38,39,58] and loose upper bounds to the total near-field thermal conductance have been determined [14,36]. In contrast, here, we provide a tight bound to the thermal conductance, h max . Other than the wellknown dependence from the gap-size d −2 , h max also rigorously demonstrates the role of other material dispersion characteristics in NFRHT.
can be written as  We derive the analytical expression of the emission spectrum and the heat-transfer conductance for two semi-infinite layers made of non-magnetic, isotropic, homogeneous polaritonic media. We clarify the approximations leading to the compact heat-transfer conductance expression, given in Eq. (8) of the main manuscript, and we estimate a material-specific optimal operating temperature in the near-field, which differs significantly from Wein's displacement law, relevant in the far-field regime. We also compare our results with existing literature [1][2][3].
For gap sizes smaller than the thermal wavelength, , evanescent waves (β > k 0 ) dominate the heat flux. In this range, the corresponding transmission probability, per polarization, can be expressed as [5][6][7][8]: where k 0 = ω/c 0 is the free-space wavenumber, c 0 is the speed of light in vacuum, η 0 = β 2 − k 2 0 is the out-ofplane wavenumber in vacuum, and r p,s are the Fresnel coefficients at the vacuum-material interface for p-and s-waves, respectively [9]. In plasmonic and polar homogeneous non-magnetic isotropic media, thermally excited SPPs and SPhPs, respectively, dominate NFRHT. These can only be excited in p-polarization [10], hence, the emission spectrum, Φ, can be approximated as the contribution from p-polarization alone, without loss of generality. Furthermore, at sufficiently small vacuum gaps, the dispersion of surface polaritons approaches the quasistatic limit, for which β k 0 [11]. In this limit, η 0 ≈ β, and one can approximate the Fresnel coefficient with The transmission probability ξ p for p-polarization (Eq. S1) in the electrostatic limit β k 0 can be written as: Perfect photon tunneling occurs at ξ p (ω, β) = 1. From Eq. (S3), this occurs at an in-plane wavenumber of: Eq.
(S4) is valid for frequencies ω such that |1 + 2/χ(ω)| < 1, and defines a curve in the (ω, β) parameter space, near which the NFRHT is maximal [1]. The maximum β res that satisfies Eq. (S4) occurs near the surface polariton resonance frequency, Ω, such that The logarithmic dependence of β max from the imaginary part of the Fresnel coefficient showcases the role of material loss in NFRHT [1]. The emission spectrum Φ is therefore given by where we have assumed k 0 d 1 and we have made the substitution βd → x. By plugging the expression of ξ p given in Eq. (S3) in Eq. (S6), and by making the substitution e −2x → y, we obtain: We now solve the last two integrals in Eq. (S7). We note that the function log y 1−y is analytical inside the contour L in the complex plane depicted in Fig. S1. Thus, by applying Cauchy integral's theorem [12], hence the first integral is solved as: where Li 2 is the dilogarithm or Spence's function, defined as [13] Li The second integral is simply . Therefore, by plugging this result and the one of Eq. (S8), in Eq. (S7), we can express the thermal emission spectrum, Φ, as: where we have used the dilogarithm identity Li 2 [z] + Li 2 [1 − z] = π 2 6 − log z log (1 − z), ∀z ∈ C \ {1} [14]. Eq. (S10) coincides with Eq. (7) in the main manuscript, and applies to any material dispersion.
Eq. S12 agrees exactly with the result by Miller et al.
in [3] (Eq. (10)), derived for planar configurations. By contrast, in the high-loss limit, Eq. (S14) may have no real solutions, Ω, and the maximum of Φ needs to be calculated by maximizing the right hand side in Eq. (S10). The range of frequencies for which Eq. (S12) is valid and the threshold between low-loss and high-loss regimes is discussed in the following section. We now showcase the validity of Eq. (S10) via comparison with fluctuational electrodynamics.
The dielectric function of plasmonic and polar materials can be adequately described by the Drude and Lorentz oscillator models, respectively: Here, ε ∞ is the high-frequency relative permittivity, and γ is the damping rate [15,16] in both plasmonic and polar media. For plasmonic metals, ω p is the plasma frequency, near which the SPP mode occurs. For polar materials, ω TO and ω LO are the transverse and longitudinal optical phonon frequencies, respectively [17]. The spectral range [ω TO , ω LO ] is known as the Reststrahlen band [18], and the SPhP mode occurs within this range. Drude metals and polar dielectrics support surface polaritons at the resonance frequencyΩ such that Re r p (Ω) = 0 or Re 2 χ(Ω) = −1. (S14) As an example, in Fig. S2, we plot Φ(ω) normalized with d −2 for silicon carbide (SiC), a widely used polar dielectric in the NFRHT literature [20][21][22]. We consider a representative Lorentz model for its permittivity, with ε ∞ = 6.7, ω TO = 1.49 × 10 14 rad/s, ω LO = 1.83 × 10 14 rad/s, γ = 8.97 × 10 11 rad/s [23]. The exact emission spectrum for p-polarization (Φ p ), obtained via fluctuational electrodynamics (FE), viz. via numerical integration of Eq. (S6) with no approximations, is shown with the red curve. Its s-polarization counterpart as well as their sum, Φ = Φ p + Φ s , are shown  S2. Thermal emission spectrum, Φ, normalized by d −2 as a function of frequency, in the neighborhood of the SPhP resonance Ω for two bulk planar layers made of SiC exchanging thermal radiation in the near-field. Φp, Φs, and Φ, shown with the red, blue, and black curve, respectively, correspond to the p-polarization, s-polarization, and total spectrum as computed via fluctuational electrodynamics, whereas the analytical prediction in Eq. (S10) as well as in [2] is shown with the green dashed line. We also show the spectral upper bound Φopt from [19] (orange curve) and the emission spectrum at the resonance frequency from [3] (purple triangle), given in Eq. (S12).
with the blue and black curve, respectively. As anticipated, the p-polarization component dominates the emission spectrum in almost the entire frequency range near Ω, and thus coincides with the total Φ. Importantly, the green dashed curve shows Φ p obtained with the analytical solution in Eq. (S10) and in [2]. This curve overlaps nearly perfectly with FE for frequencies nearresonance. We also display an upper bound to the spectrum of NFRHT, Φ opt , as derived by Venkataram et al. [19], obtained through singular value decomposition of the Maxwell Green's tensor. This result accurately estimates the response of SiC only on resonance. Similarly, we show with a triangle-shaped marker the upper bound to NFRHT on resonance, derived by Miller et al. [3], and given in Eq. (S12).

RADIATIVE THERMAL CONDUCTANCE h
The radiative thermal conductance for two closely spaced semi-infinite layers is defined as: where where θ(ω, T ) = ω/ [exp ( ω/k B T ) − 1] is the mean energy per photon and Φ(ω) is the emission spectrum, whose closed form expression has been derived in the previous section and is given in Eq. (S10). We now particularize this derivation to plasmonic and polar media, whose dispersion relations are given in Eqs. (S13) and (S13b), respectively. By applying Eq. (S14), for γ → 0, the surface plasmon polariton (SPP) and surface phonon polariton (SPhP) resonance frequency can be expressed as: We assume that the function Π ω T (Eq. (9) of the main text) is slowly varying with respect to the emission spectrum Φ(ω), which peaks at the polariton resonance frequencyΩ. Therefore, we make the first step toward the analytical integration of Eq. (S15) by sampling Π at the frequencyΩ, i.e.: We now carry out the frequency integration of the emission spectrum in Eq. (S18). By using its expression in Eq. (S10), we have (S19) Since r p (ω), defined in Eq. (S2), inherits the hermiticity (or PT-symmetry) from the permittivity function ε(ω), i.e. r * p (ω) = −r p (−ω) ( * is the complex-conjugate operator), the integrand in Eq. (S19) is an even function of ω. Therefore, we can extend the integration to include the negative frequency axis, namely where the complex valued function f (z) is We now tackle the complex integration of f (ω) in Eq. (S20) by means of contour integration in the complex plane, following a similar strategy to the one empolyed in [2]. Specifically, we intend to perform the integration on the closed contour L shown in Fig. S3, composed by the real axis and a semicircular contour of positive radius R lying the lower half-plane, in the limit R → ∞. Since the integrand is vanishing on the semicircular contour in the limit lim R→∞ , from Jordan's lemma [24], the integral on this contour is also zero. Hence, we can rewrite the integrated emission spectrum as:

Re
The function f (z) is the analytical extension of the realvalued function f (ω) to the complex plane. It must be noted that Im{r p }(z), Re {r p } (z) are no longer constrained to be real valued functions, and therefore the notation Im{·}, Re {·} no longer refers to the real and imaginary part operators. Nevertheless, we keep using the same notation in the following calculations for the sake of simplicity, while accounting that Im{r p }(z), Re {r p } (z) are functions derived for real variable, and extended to the complex plane. For instance, the step Im{r p (ω)} → Im{r p }(z) for plasmonic dispersions is the following: We now carry the complex integration using the Cauchy's residue theorem [24]. The first step is identifying the poles of the integrand f (z). By inspecting the expression of f (z) in Eq. (S21), it is clear that the poles are exactly the frequenciesΩ solving the resonance condition in Eq. (S14). By solving Eq. (S14), using the expression of the plasmonic and polar permittivity given in Eq. (S13), we have: where Ω is the polariton resonance frequency in the absence of optical losses, given in Eq. (S17), and Q is the quality factor of the polaritonic material resonance, evaluated at frequency ω = Ω, namely [25]: The parameter F for plasmonic and polar cases can be written as: .
The parameter B, which in the main manuscript we have labeled as material residue function, is independent of the material losses (independent of γ), and is calculated as B = Q Im{r p (Ω)} , in the limit Q → ∞. Therefore, for plasmonic and polar dispersions, respectively, B can be written as: Finally, the parameters Q th and Q 2 have the following expressions: where F for plasmonic and polar dispersions is given in Eq. (S26). From Eq. (S24), in the case ε ∞ = 1, it can be shown thatΩ 1,2 are real numbers only if Q < Q 2 and Q > Q th . It can be proven that Q 2 < 1, and hence we can focus only on the cases Q < Q th and Q > Q th . Therefore, Q th represents a threshold for the material quality factor below which the poles of f (z) become complex and move away from the real axis, as shown in Figs. S3b-c.
We can now apply the residue theorem to carry out the complex integration in Eq. (S22). It is important to notice that, since for ε ∞ = 1 and Q > Q th the poles of f (z) are on the integration contour, we have to add a 1 2 factor to the standard residue theorem formula, for which the poles are in the interior of the integration contour L [24]. On the other hand, the frequencyΩ in the case ε ∞ = 1 does not depend on the Q factor, and the poles are always real.
Case ε ∞ = 1. Via algebraic manipulation, it can be shown that the function f (z) for both plasmonic and polar dispersions in the case ε ∞ = 1 has the following expression: Therefore, f (z) has two first-order real poles, i.e. where Res(f, w) is the residue of f at w. Here, we have used the parity of f , and the minus sign comes from having chosen a clockwise (negative) orientation of the contour L, shown in Fig. S3a. By taking the limit lim z→Ω (z − Ω)f (z), we calculate the residue Res(f, +Ω), which has the following expression: By plugging this result in Eq. (S31), we can finally write the expression for the radiative thermal conductance as: where Π Ω T is defined in Eq. (S16), and Eq. (S33) coincides with Eq. (8) in the main manuscript.
It must be noted that all the redundant scaling factors, e.g. B at the denominator of h max and Ψ, have been introduced such that Ψ would be bounded above by 1.
The function Ψ reaches its maximum at: Eq.
Case ε ∞ = 1. The function f (z) for both plasmonic and polar dispersions in the case ε ∞ = 1 has the following expression: where F is given in Eq. (S26). For ε ∞ = 1, we have different resuts according to the position of Q with respect to Q th . From Eq. (S24), it is clear that if Q > Q th , then the function f (z) has four first-order real poles {±Ω 1 , ±Ω 2 }, as shown in Fig. S3b; if Q < Q th then the function f (z) has four first-order complex poles {±Ω 1 , ±Ω 2 }, with Ω 2 = −Ω * 1 as shown in Fig. S3c. Thus, if Q < Q th , the only poles contributing to the integral are the two in the interior of the contour L, viz. {Ω 1 , Ω 2 }. Therefore, in both cases Q > Q th and Q < Q th we can solve Eq. (S22) by applying the residue theorem as follows: By making the limit lim z→Ωj (z − Ω j )f (z), we calculate the residues Res(f, Ω 1 ), Res(f, Ω 2 ) which have the following expression: where F is given in (S26), and the complex function Im{r p }(z) is the analytical extension of the imaginary part of the Fresnel coefficient in the complex plane (e.g., see Eq. (S23) for the plasmonic dispersion). By inserting the residues in Eq. (S38) into Eq. (S37), and in turn plugging this into Eq. (S22), one can finally get the expression for the heat conductance.

APPROXIMATIONS FOR h AND THE POLARITON RESONANCE FREQUENCYΩ
We now simplify the exact expressions for h, derived in the previous section, for the three scenarios, viz. {ε ∞ = 1, ∀Q}, {ε ∞ = 1, Q < Q th } and {ε ∞ = 1, Q > Q th }, aiming at providing a single expression valid in any regime of Q and ε ∞ .
We make the following approximations: (i ) we neglect the contribution from the second pole Ω 2 ; (ii ) we assume Q Q th . Under these assumptions, we can approximate Ω 1 ≈ Ω, with Ω given in Eq. (S17). It can be shown that the resulting residue Res(f, Ω) has the same form as the case ε ∞ = 1 in Eq. (S32), and the heat transfer conductance expression is the same as in Eq. (S33) and Eq. (13) in the main manuscript.
Even if this expression is derived in the high-Q limit, it represents a good approximation also in the low-Q case, as shown in Fig. S4 (a) for a case of study. Specifically, we show the integral in Eq. (S22) for a polar dispersion with Lorentz parameters ε ∞ = 4, ω TO = 1.49 × 10 14 rad/s, ω LO = 1.83×10 14 rad/s, Ω = 1.77×10 14 rad/s as a function of the quality factor Q = Ω/γ, and there is good agreement with the exact solution over all the considered Q range.
It must be noted that the polariton resonance fre-quencyΩ 1 given in Eq. (S24) for both plasmonic and polar dispersions is weakly dependent from the optical loss, i.e. γ or the quality factor, Q. In Fig. S4 (b) We show this by monitoring the peak position of the emission spectrum, Φ(ω), for the same polar dispersion used in panel (a), for decreasing values of the quality factor, Q, or equivalently for increasing value of γ, starting from Q 1 ≈ 200 and arriving to Q = Q 1 /50 ≈ 4. It is clear that, even in the lowest-Q case, the emission spectrum peaks very closely to the resonance frequency Ω calculated in the limit γ → 0 or Q → ∞, given in Eq. (S17). (a) Integral in Eq. (S22), equivalent to 8π 2 d 2 ∞ 0 Φ(ω)dω, for a polar dispersion with Lorentz parameters ε∞ = 4, ωTO = 1.49 × 10 14 rad/s, ωLO = 1.83 × 10 14 rad/s, and resonance frequency Ω = 1.77 × 10 14 rad/s. We compare the exact solution (black curve) calculated using Eq. (S37), with the approximated one (red dashed line) calculated using Eq. (S31), and used as final result in the main manuscript. We also show with a green line the threshold Q th = 10.85 between the low-and high-loss regimes, calculated using Eq. (S28). For the same dispersion parameters, in (b) we plot the emission spectrum normalized by d −2 , for increasing values of Q = Ω/γ, namely Q = Q1/50, Q1/45, Q1/40, . . . , Q1, being Q1 ≈ 200. We also mark the peak of each emission spectrum with a red marker: in every case, the resonance occurs at frequencies very close to Ω, calculated for Q → ∞.
Thus, assuming Ω as the polariton resonance frequency in any material loss condition, as considered in the main manuscript, represents a good approximation.

TEMPERATURE DEPENDENCE IN THE NEAR-FIELD VS WIEN'S LAW
We now discuss the temperature dependence of NFRHT in bulk systems. In the inset of Fig. S5, we display Π, given in Eq. (S16), which decreases exponentially as a function of the ratio Ω/T . Hence, to avoid a dramatic damping of h, a higher resonance frequency  [26], while the Drude parameter for the other plasmonic materials and the Lorentz parameters for polar materials are taken from [17] (tables 1-2). The superscript o and e stand for the ordinary and extraordinary principal axes of the corresponding material, respectively. In the inset, we plot Π (in Eq. (S16)), which expresses the normalized thermal conductance in the optimal loss condition Q = Qopt = 4.5 B (see Eq. (S35)). We identify the optimal operating temperature for a material by setting Π = 0.9 (green marker).
should be compensated by a higher operating temperature, as expected. This is well-understood in the far-field regime with Wien's displacement law that estimates the optimal resonance frequency of a thermal emitter at a given temperature for maximizing the power emitted in the farfield. One can similarly estimate the optimal temperature, T opt , of a polaritonic thermal emitter in a planar near-field configuration, by maximizing Π Ω T . Since Π reaches its maximum Π = 1 in the limit of infinite temperature, we compute the optimal temperature of operation as a function of resonance frequency, Ω, by setting the term Π Ω T to 0.9, as shown with the green lines in the inset of Fig. S5. The resulting optimal temperature is expressed as: where b NF = 12729µm K ≈ 4.4 × b Wien (in which the subscript stands for Near-Field) and λ = 2πc 0 /Ω is the resonance wavelength. This dependence of T opt from Ω is shown with the solid line in Fig. S5. As a reference, we also display with the dashed line Wein's displacement law, relevant in the far-field. One can therefore see that in the near-field, considerably higher temperatures are required to reach optimal material performance, as com-pared to the far-field. This is expected, since far-field thermal emission is generally more broadband with respect to the near-field thermal emission, especially for polaritonic materials [20,27]. Hence, the maximal spectral overlap between the mean energy per photon, θ(ω, T ) in Eq. (S15), and a blackbody spectrum is achieved at much lower temperatures as compared to the overlap between a narrowband near-field resonance and θ(ω, T ), since θ broadens as T increases. We display the optimal temperature, calculated using Eq. (S39), as a function of the resonance wavelength 2πc o /Ω, for polar materials such as Silicon Carbide (SiC), hexagonal Boron Nitride (hBN), and doped semiconductors, e.g. Gallium Arsenide (GaAs), Indium Arsenide (InAs) [17,28,29], as well as plasmonic materials such as standard noble metals, e.g. Gold (Au), Silver (Ag), and heavily doped oxides, e.g. IZO and GZO [17,30,31], as in Figs. 1-3 of the main text. It can be seen that most polar media, with resonance frequencies mainly in the mid-IR, will achieve optimal performance at temperatures that are up to two orders of magnitude lower than their plasmonic counterparts, since plasmonic resonances occur mainly in the near-IR, visible and UV regimes.

COMPARISON WITH LITERATURE [1]
We now compare the near-field radiative thermal conductance derived in this work with the expression for polar dielectrics provided by Ben-Abdallah et al. in [1]. In [1], the authors considered a polar material dispersion with high-frequency dielectric permittivity ε ∞ = 1 (see Eq. (S13)). In their derivation (Eq. (14)), the radiative thermal conductance that they obtained, h , can be written as: where Π Ω T is the same as in Eq. (S33), given in Eq. (S16), and the functions h max and Ψ are given by: and Ψ (Q) = log Q 0.37 Q .
In both our result (Eq. (S33)) and the result from [1] (Eq. (S40)), since (Ψ, Ψ ) and Π are functions bounded above by 1, h max and h max represent the maximum heattransfer conductance achievable in this configuration. As we shown in the previous sections, the material residue, B, is greater than unity, i.e. B > 1, therefore we can compare h max and h max as follows: [1] given in Eq. (S42), respectively, as a function of the polariton resonance quality factor Q. We considered a Lorentz dispersion with parameters ε∞ = 1, ωTO = 1.49 × 10 14 rad/s, ωLO = 1.83 × 10 14 rad/s, for which B = 4.93, calculated using Eq. (S27b). The different optimal quality factors at which the curves peak are also marked with a dashed line.
Therefore our analytical estimation for the maximum radiative thermal conductance is at least 4.4 times smaller than the one predicted in [1] under the assumption of blackbody-like thermal emission in the near-field (ξ = 1). We now compare the functions (Ψ, Ψ ), taking into account the dependence of heat transfer from optical loss. The function Ψ depends only on the quality factor Q of the polariton resonance, and neglects the dependence from the other features of a material's dispersion, such as the size of the Reststrahlen band (ω LO − ω TO ) and its position. In our work, these are included in the material residue, B. According to the definition of Ψ (Q), in Eq. (S42), this function is maximized at the optimal quality factor Q * opt = e, being e ≈ 2.72 the Napier's constant, valid for any polar dielectric with a Lorentz dispersion relation (Eq. (S13)), with ε ∞ = 1. Conversely, from our derivation, the function Ψ Q B depends on the dispersion's characteristics through the factor B, and the optimal quality factor Q is given by Q opt = 4.5 B.