An ab-initio coupled mode theory for near field radiative thermal transfer

We investigate the thermal transfer between finite-thickness planar slabs which support surface phonon polariton modes (SPhPs). The thickness-dependent dispersion of SPhPs in such layered materials provides a unique opportunity to manipulate and enhance the near field thermal transfer. The key accomplishment of this paper is the development of an ab-initio coupled mode theory that accurately describes all of its thermal transfer properties. We illustrate how the coupled mode parameters can be obtained in a direct fashion from the dispersion relation of the relevant modes of the system. This is illustrated for the specific case of a semi-infinite SiC substrate placed in close proximity to a thin slab of SiC. This is a system that exhibits rich physics in terms of its thermal transfer properties, despite the seemingly simple geometry. This includes a universal scaling behavior of the thermal conductance with the slab thickness and spacing. The work highlights and further increases the value of coupled mode theories in rapidly calculating and intuitively understanding near-field transfer. © 2014 Optical Society of America OCIS codes: (030.5620) Radiative transfer; (230.4555) Coupled resonators; (240.5420) Polaritons. References and links 1. E. Rousseau, A. Siria, G. Jourdan, S. Volz, F. Comin, J. Chevrier, and J.-J. Greffet, “Radiative heat transfer at the nanoscale,” Nature Photonics 3(9), 514–517 (2009). 2. S. Shen, A. Narayanaswamy, and G. Chen, “Surface phonon polaritons mediated energy transfer between nanoscale gaps,” Nano Lett. 9(8), 2909–13 (2009). 3. S. Basu, Y.-B. Chen, and Z. M. Zhang, “Microscale radiation in thermophotovoltaic devices-a review,” International J. Energy Res. 31(6), 689–716 (2007). 4. S. Y. Lin, J. Moreno, and J. G. Fleming, “Three-dimensional photonic-crystal emitter for thermal photovoltaic power generation,” Appl. Phys. Lett. 83(2), 380 (2003). 5. J. A. Schuller, T. Taubner, and M. L. Brongersma, “Optical antenna thermal emitters,” Nature Photonics 3(11), 658–661 (2009). 6. N. Dahan, A. Niv, G. Biener, Y. Gorodetski, V. Kleiner, and E. Hasman, “Enhanced coherency of thermal emission: beyond the limitation imposed by delocalized surface waves,” Phys. Rev. B 76(4), 045427 (2007). 7. D. Polder and M. Van Hove, “Theory of radiative heat transfer between closely spaced bodies,” Phys. Rev. B 4(10), 3303–3314 (1971). 8. A. V. Shchegrov, K. Joulain, R. Carminati, and J.-J. Greffet, “Near-field spectral effects due to electromagnetic surface excitations,” Phys. Rev. Lett. 85(7), 1548–1551 (2000). 9. J.-P. Mulet, K. Joulain, R. Carminati, and J.-J. Greffet, “Enhanced radiative heat transfer at nanometric distances,” Microscale Thermophysical Engineering 6(3), 209–222 (2002). #222793 $15.00 USD Received 11 Sep 2014; revised 27 Oct 2014; accepted 14 Nov 2014; published 24 Nov 2014 (C) 2014 OSA 1 December 2014 | Vol. 22, No. 24 | DOI:10.1364/OE.22.030032 | OPTICS EXPRESS 30032 10. M. Planck, The Theory of Heat Radiation (Blakiston’s Son & Co, 1914). 11. N. Shitrit, I. Yulevich, E. Maguid, D. Ozeri, D. Veksler, V. Kleiner, and E. Hasman, “Spin-optical metamaterial route to spin-controlled photonics,” Science 340(6133), 724–726 (2013). 12. N. Dahan, Y. Gorodetski, K. Frischwasser, V. Kleiner, and E. Hasman, “Geometric doppler effect: spin-split dispersion of thermal radiation,” Phys. Rev. Lett. 105(13), 136402 (2010). 13. M. Krüger, G. Bimonte, T. Emig, and M. Kardar, “Trace formulas for nonequilibrium casimir interactions, heat radiation, and heat transfer for arbitrary objects,” Phys. Rev. B 86(11), 115423 (2012). 14. P. Ben-Abdallah, K. Joulain, J. Drevillon, and G. Domingues, “Near-field heat transfer mediated by surface wave hybridization between two films,” J. Appl. Phys. 106(4), 044306 (2009). 15. M. Francoeur, M. P. Mengüç, and R. Vaillon, “Spectral tuning of near-field radiative heat flux between two thin silicon carbide films,” J. Phys.D: Appl. Phys. 43(7), 075501 (2010). 16. M. Francoeur, M. P. Mengüç, and R. Vaillon, “Coexistence of multiple regimes for near-field thermal radiation between two layers supporting surface phonon polaritons in the infrared,” Phys. Rev. B 84(7), 075436 (2011). 17. M. Krüger, T. Emig, and M. Kardar, “Nonequilibrium electromagnetic fluctuations: heat transfer and interactions,” Phys. Rev. Lett. 106(21), 210404 (2011). 18. C. Otey and S. Fan, “Numerically exact calculation of electromagnetic heat transfer between a dielectric sphere and plate,” Phys. Rev. B 84(24), 245431 (2011). 19. A. P. McCauley, M. T. H. Reid, M. Krüger, and S. G. Johnson, “Modeling near-field radiative heat transfer from sharp objects using a general three-dimensional numerical scattering technique,” Phys. Rev. B 85(16), 165104 (2012). 20. M. T. H. Reid, A. W. Rodriguez, and S. G. Johnson, “Fluctuation-induced phenomena in nanoscale systems: harnessing the power of noise,” Proceedings of the IEEE 101(2), 531–545 (2013). 21. K. Sasihithlu and A. Narayanaswamy, “Proximity effects in radiative heat transfer,” Phys. Rev. B 83(16), 161406 (2011). 22. C. R. Otey, L. Zhu, S. Sandhu, and S. Fan, “Fluctuational electrodynamics calculations of near-field heat transfer in non-planar geometries: a brief overview,” J. Quantitative Spectroscopy and Radiative Transfer 132, 3–11 (2014). 23. H. Haus and W. Huang, “Coupled-mode theory,” Proceedings of the IEEE 79(10), 1505–1518 (1991). 24. C. R. Otey, W. T. Lau, and S. Fan, “Thermal rectification through vacuum,” Phys. Rev. Lett. 104(15), 154301 (2010). 25. J. C. Swihart, “Field solution for a thin-film superconducting strip transmission line,” J. Appl. Phys. 32(3), 461 (1961). 26. E. Economou, “Surface plasmons in thin films,” Phys. Rev. 182(2), 539–554 (1969). 27. W. Eckhardt, “First and second fluctuation-dissipation-theorem in electromagnetic fluctuation theory,” Opt. Commun. 41(5), 305–309 (1982). 28. K. Joulain, J.-P. Mulet, F. Marquier, R. Carminati, and J.-J. Greffet, “Surface electromagnetic waves thermally excited: radiative heat transfer, coherence properties and casimir forces revisited in the near field,” Surf. Sci. Rep. 57(3), 59–112 (2005). 29. S. Biehs, E. Rousseau, and J. Greffet, “Mesoscopic description of radiative heat transfer at the nanoscale.” Phys. Rev. Lett. 105(23), 234301 (2010). 30. P. Ben-Abdallah and K. Joulain, “Fundamental limits for noncontact transfers between two bodies,” Phys. Rev. B 82(12), 121419 (2010). 31. W. Spitzer, D. Kleinman, and D. Walsh, “Infrared properties of hexagonal silicon carbide,” Physical Review 113(1), 127–132 (1959). 32. W. Spitzer, D. Kleinman, and C. Frosch, “Infrared properties of cubic silicon carbide films,” Phys. Rev. 113(1), 133–136 (1959). 33. R. D. Kekatpure, A. C. Hryciw, E. S. Barnard, and M. L. Brongersma, “Solving dielectric and plasmonic waveguide dispersion relations on a pocket calculator.” Opt. Express 17(26), 24112–29 (2009).


Introduction
Obtaining new methods to control the thermal transfer is of great concern for a variety of important applications, including magnetic heat assisted lithography and recording [1,2] to energy conversion [3,4] and thermal emitters [5,6]. It was demonstrated that near-field effects [7][8][9] can be harnessed to achieve dramatic increases in thermal transfer compared over conventional farfield radiative transfer following Planck's law [10]. These enhancements are facilitated through opening new channels for transfer in closely-spaced bodies that support evanescent optical modes [7]. Fundamentally new levels of control can be achieved on thermal transfer by tailoring the dispersive properties of such systems [11,12]. Even though, near-field thermal transfer d t SiC SiC Vacuum Fig. 1. (a) Schematic of the Swihart geometry which is composed of a thick slab of SiC at a temperature T 1 separated by a vacuum gap of width d from a slab of SiC of thickness t and at a temperature T 2 . The profiles of the three supported SPhP modes are shown in this structure in the limit of negligible coupling (b) Schematic of the mode profiles at a frequency right below the SPhP resonance frequency for three special cases of (i) t = 0.5d (ii) t = d (iii) t = 2d. Curves show the normalized absolute value of the magnetic field.
for planar structures have been calculated for different thicknesses and the law of its variation has been considered Rigorous theory exists to quantify thermal transfer between systems of arbitrary shape and number of components [13]. However, exact numerical solutions for arbitrary structures requires tremendous amount of computational power. As a result, near-field thermal transfer calculations have been limited to relatively basic geometries such as planar-to-planar [7,[14][15][16] as well as planar structures to a sphere [17][18][19], a cylinder [19], and even a cone [19]. Reference [20] reviews some of the recent results in different considered structures. Approximate theories have been developed to make thermal transfer calculations more tractable and the proximity method is perhaps the most popular [1,21,22]. Unfortunately, this method cannot be applied to sub wavelength structures and fails in the far-field regime. Another approximation proposed for calculating the near field thermal transfer is coupled mode theory [23]. Whereas this method provides some valuable intuition on the relevant modes that control the transfer, the coupled mode parameters were obtained from best fits to a full fluctuation dynamics calculation [24]. Here, we illustrate how these parameters can be obtained in a direct fashion from the dispersion relation of the relevant modes of the system. This significantly increases the value of coupled mode theories in rapidly calculating and intuitively understanding near-field transfer.
In this letter, we consider the specific case of a SiC slab that is closely-spaced to an underlying substrate slab (see Fig. 1(a)). This system has been studied in the context of superconductivity and is termed the Swihart geometry [25]. SiC is a polar material and planar structures support surface waves termed surface phonon polaritons (SPhPs). These modes provide the dominant channels for near-field thermal transfer [8].

Exact theory
Near-field transfer is critically-dependent on the nature of the relevant modes and their dispersive properties, which are well-established for the Swihart geometry [26]. Figure 1(b) illustrates the mode profiles at a frequency just below the SPhP frequency ω SPhP and for three different thicknesses. These mode profiles results from a hybridization of the individual interface modes shown in Fig. 1(a). The nature of the hybridization evolves as the ratio of t/d is increased. Whereas for small t/d the mode profile resembles that of a slab-mode, for large t/d the mode profile closely mimics a gap mode. When t is approximately equal to d, there is strong hybridization of these modes. This hybridization can lead to an increased thermal transfer, which will be shown. The presented coupled mode theory will provide insights into this fact.
Using the green function method, combined by the fluctuation dissipation theorem [27,28] (Appendix A), thermal transfer between the two bodies can be written as: where in that V is given by: in which R p,s and T p,s are reflection and transmission coefficients corresponding to p and s polarized waves and Θ (ω, T ) =hω/ (exp (hω/k B T ) − 1) is the mean energy of a harmonic oscillator with angular frequency of ω at temperature T . Moreover, β corresponds to the transverse component of the wave vector. Finally, k z = εω 2 /c 2 − β 2 and k zv = ω 2 /c 2 − β 2 are z components of the wave vector in SiC and vacuum, respectively.
These relations can be simplified in the near field regime. In this regime the contribution by the s polarized waves becomes negligible compared with the p polarized ones. This comes back to the fact that only p polarized waves support SPhP resonances. In this regime, if we define the ratio ζ as ζ = (1 − ε) / (1 + ε), the thermal transfer can be written as: Where S (ω, β ) quantifies the contribution to the transfer at different ω and β values and is given by: This can be viewed as a generalization of the earlier results obtained for the case of two semi-infinite structures [29,30]. Using the above expression for the S (ω, β ), it is easy to show that the thermal transfer can be written as a function f of the ratio t/d: This generalizes the previously reported results for the near-field thermal transfer between semi-infinite slabs of SiC which is proportional to d −2 [29,30]. The asymptotic behavior of the near-field thermal transfer for planar structures as a function of distance for different regimes of thickness to distance ratio has been considered recently [16]. This expression, however, shows that the only dependence on the thickness for the Swihart geometry comes from the t/d ratio. Equation (5) shows that the thermal transfer obeys a universal scaling with t, for a given vacuum gap. Moreover, it achieves its maximum value for a thickness which is proportional to the gap. This comes from the fact that the modal hybridization and the associated near field transfer is dependent on the thickness through the ratio between thickness and vacuum gap. Note should be added here that these results hold for the case of near-field regime, where the contribution of the s polarized waves becomes negligible.
In Fig. 2(a), we plot the normalized thermal transfer rate as a function of the ratio between slab thickness and vacuum gap for four different sizes of vacuum gaps. For these calculations, the permittivity of SiC is taken to be in the form of with ε ∞ = 6.7, ε s = 10, δ /ω 0 = 0.006 and ω 0 /(2π) = 2.38 × 10 13 sec −1 (12.6µm) [31,32]. In addition, the temperature that is assumed in the numerical calculations is T = 315K.
As we expect per Eq. (5), all the curves follow the same universal trend in this plot. The maximum value is achieved for a slab thickness nearly equal to gap width. This leads to about 10% enhancement of the thermal transfer relative to the case of semi-infinite sized slabs. It is interesting to note that the thermal transfer decreases with increasing d for d > t, despite the fact that the slab volume (i.e. the amount of material) increases. Figure 2(b) shows the expected d −2 scaling of the maximum achievable thermal transfer rate with gap size on a log-log scale. Figure 3 shows the variation of S (ω, β ) as a function of ω, β for three different thicknesses of t = 0.5d, t = d, and t = 2d. In each panel, the calculated dispersion relationships for the relevant surface modes are shown with dashed white lines. These correspond to the real ω vs. β dispersion profiles of the coupled surface modes (CSMs). It is clear that the main contribution to the thermal transfer is coming from those CSMs. For the case of Swihart geometry, there is an extra channel for thermal transfer at exactly ω = ω SPhP , contrary to the case of thermal transfer between semi-infinite substrates [29]. In continuation, by using the coupled mode theory, we show how the individual surface resonances at ω = ω SPhP , branch into three coupled surface modes.
The spectral contributions to the near-field thermal transfer are maximized near the SPhP resonance frequency. At this frequency, if we write the epsilon as ε = −1 + jε , the coefficient ζ becomes ζ (ω SPhP ) = −1 − 2 j/ε . As it is the case for SiC, we assume that |ε | 2. As a consequence, S (ω SPhP , β ) can be written as:

Coupled mode theory
Next, we develop a coupled mode theory [23] to understand the unique thermal transfer properties of the Swihart geometry. This structure supports three basic surface modes (see Fig. 1(a)). The coupling between these three modes affords thermal transport between two bodies in the near-field regime. First of all, using the above mentioned form for the permittivity, it is easy to show that each of the isolated planar objects in the limit of large β resonates at frequency of ω SPhP = ω 0 (1 + ε s ) 1 2 (1 + ε ∞ ) − 1 2 and exhibits a decay rate of γ = δ /2. For the gap and slab-like modes, the resonance frequencies are given by ω SPhP ± k. Here, the parameter k quantifies the coupling between two adjacent surface modes and can be calculated as . Also, it is easy to show that δ can be written in terms of the imaginary part of epsilon at ω SPhP as δ = −ω 0 Γε . Details of these calculations are done in the appendix B.
For the Swihart geometry, we have three interfaces leading to three surface modes. From the relevant dispersion relations [26], we find that in the limit of large β , their Eigen frequencies are given by ω = ω SPhP and ω = ω SPhP ± Γω 0 e −2β d + e −2βt − e −2β d e −2βt 1 2 . Each of them has a corresponding decay rate of γ = δ /2. Three modes can be seen to be formed as a result of coupling of three original isolated modes each of them supports a SPhP resonance at frequency of ω SPhP and the decay rate of γ = δ /2. In the appendix C, it is shown that by taking the amplitudes of the original modes to be a 1 , a 2 , and a 3 , the coupling between these modes can be expressed as: Where A = a 1 a 2 a 3 is the matrix of amplitudes and N = n 1 n 2 n 3 is the  matrix of noise fluctuations. Here, κ is the Hermitian coupling matrix and is given as: Using the above coupling matrix, the Eigen modes of the system are consistent with the previously calculated ones based on dispersion relations. As such, these relations extend the coupled mode theory obtained for two semi-infinite SiC slabs [24] to the case of Swihart geometry. Whereas the required coupling coefficients and decay rates in the coupled mode theory were previously obtained by parameter fitting, here they are derived from the mode analysis and dispersion relations of the CSMs. Note that the thermal transfer between two bodies is carried out by the coupling of the surface modes located at two sides of the vacuum gap (resonators 1 and 2). Moreover, there is a coupling between two sides of the thin slab of SiC (resonators 2 and 3). However, no direct coupling exists between resonators 1 and 3. Since the first resonator is kept at temperature T 1 and the second and third resonators are kept at temperature T 2 , then [23]: Following the coupled mode theory, the power transfer between two modes can be written as S (t) = 2ℑ [κ 12 a * 1 (t) a 2 (t) ]. By using the above-mentioned expressions and taking the Fourier transform of the coupled differential equations, we can then calculate the thermal transfer between two bodies. After summing over all transfer channels at different β s and ωs, we obtain the total thermal transfer between the two bodies. Details of those calculations are described in the appendix D. With the total thermal transfer expressed by Eq. (3), the final result for the S (ω, β ) using couple mode theory can be expressed as: At ω SPhP this gives exactly the same result as the one obtained by the exact method in the limit of |ε | 2. Figure 3 shows calculations of Eq. (10) as an inset for the three different t/d values considered. The contribution to the thermal transfer coming from the β channels at frequency of ω obtained with coupled mode theory closely resembles the exact analysis. Quantitatively, the coupled mode theory becomes increasingly accurate with increasing β . In fact, the magnitudes of the spectral peaks in S (ω, β ) exhibit perfect quantitative correspondence, even in the small β regime. We observe small shifts in the peak locations due to second order effects that cannot be handled by the first-order coupled mode theory. However, these have a negligible effect on the magnitude of the thermal transfer as its dominant contributions come from the large β regime.
To highlight the difference between the two methods of thermal transfer calculation, we plot the thermal transfer variation as a function of ω in Fig. 4. These graphs are for the structure with thickness t = d at three different values of β = 0.1/d, β = 2/d, and β = 5/d. Even for small values of β , magnitudes of the peaks of S (ω, β ) calculated by the coupled mode theory, are pretty close to the ones obtained by the exact analysis. Moreover, they essentially give the same result for the thermal transfer contribution at the frequency of ω = ω SPhP .

Conclusion
In conclusion, we have calculated the thermal transfer in a Swihart's geometry composed of a SiC slab of infinite thickness and a thin SiC slab. We have shown that a maximum in the nearfield thermal transfer is achieved for a slab thickness that is nearly equal to gap width. Moreover, we have applied coupled mode theory to accurately quantify the thermal conductance and to obtain intuitive insights into modes that govern the thermal transport. The presented approach of calculating near-field thermal transfer can also be generalized to other systems in which a limited number of near-field optical modes dominate the transfer.

Appendix
In the following appendices, we start by presenting the exact analysis in more details. The reflection and transmission coefficients used in order to calculate the exact thermal transfer are expressed here. Also, a rigorous proof is presented explaining the thermal transfer variation as a function of SiC slab thickness and vacuum gap. In order to explain the coupled mode theory, we need to investigate the dispersion relations of different planar structures made of SiC and vacuum. The next section is devoted for this investigation in basic planar geometries. We start by the simple IM structure, composed of SiC neighbored by the vacuum. Next the dispersion relation of the MIM structure, composed of a vacuum slab sandwiched by SiC layers, is considered. In continuation, a SiC slab as an IMI structure is studied. By these analysis, we are able to introduce the effect of coupling between their two supported IM modes. In the next section, we investigate the dispersion relations of the Swihart geometry. We also obtain the coupling matrix describing interactions between three IM modes supported in this geometry. Finally, in the last section by using the obtained coupling matrix, we calculate the thermal transfer rate based on the coupled mode theory.
Appendix A: Details of exact theory Swihart geometry as described in the main manuscript is composed of a thick SiC substrate separated by a vacuum gap from a thin SiC slab. The thermal transfer rate in this geometry can be calculated by using the green functions and fluctuation dissipation theorem, as: where V 12 = V 1 −V 2 with the following expressions for V 1 and V 2 : As mentioned in the main manuscript, R p,s and T p,s are reflection and transmission coefficients corresponding to p and s polarized waves and Θ (ω, T ) =hω/ (exp (hω/k B T ) − 1) is the mean energy of a harmonic oscillator with angular frequency of ω at temperature T . Moreover, β corresponds to the transverse component of the wave vector. In addition, k z = εω 2 /c 2 − β 2 and k zv = ω 2 /c 2 − β 2 are z components of the wave vector in SiC and vacuum, respectively. If we further simplify this expression we arrive at the Eq. (1), expressed in the main manuscript.
For its derivation, we should first calculate the total thermal radiation from SiC substrate to outside. We should then adjust it by subtracting from it the part that corresponds to the thermal transfer rate from SiC substrate to the vacuum at far distances which is not absorbed by the thin layer of SiC. These two contributions are represented by the factors V 1 and V 2 , respectively. There exists a rigorous theory expressing the thermal transfer in terms of general Green functions. For obtaining the above expressions, relevant green functions are calculated using the reflection and transmission coefficients. We consider the unity amplitude waves inside the SiC substrate right at its interface with the vacuum gap. Reflection coefficients correspond to the amplitude of these waves after reflection from the interface with the vacuum gap. Whereas, transmission coefficients correspond to the amplitude of these waves after traversing the thin SiC slab and the vacuum gap. These coefficients for the p polarized waves are expressed as: in which α, ξ are defined as α = 1 + εk zv /k z and ξ = 1 − εk zv /k z . Moreover, the same relationship holds for the R s and T s for s polarized waves with the exception of modified definitions for α = 1 + k zv /k z and ξ = 1 − k zv /k z . These are in fact the modified Fresnel coefficients for the Swihart geometry.
In the near-field regime, the contribution of the p polarized waves become dominant in the thermal transfer. In this regime, as it was explained in the main manuscript, the thermal transfer rate can be expressed as Eq. (3) with the S (ω, β ) expressed as Eq. (4). If we define the ratio, r = t/d and b = β d, the special form of S (ω, β ) allows us to express it as S (ω, β ) = f 1 (ω, r, b). Using this form, the thermal transfer rate can be written as: in which f 2 (r) is a function describing the variation of thermal transfer with the slab thickness to vacuum gap ratio. This form of thermal transfer explicitly describes why the maximum achievable thermal transfer at each vacuum gap varies with distance as d −2 . Moreover, it describes why thermal transfer follows a universal curve as a function of thickness to vacuum gap ratio at different vacuum gaps.
Using the same kind of approximation, we arrive at similar results as the ones for MIM structure in the limit of large β . The resonance frequencies are given again by ω SPhP ±k with the same value of k defined before. However, in this case the + sign corresponds to the asymmetric mode while the − sign corresponds to the symmetric mode.
By using this expression for the t (ω, β ), S (ω, β ) is calculated to be as Eq. (10). This completes the calculation of thermal transfer based on coupled mode theory.