The study of near field thermal radiation is gaining renewed interest thanks in part to their great potential in energy harvesting applications. It is well known that plasmonic or polaritonic materials exhibit strongly enhanced fields near the surface, but it is not trivial to quantitatively predict their impact on thermal radiation intensity in the near field. In this paper, we present a case study for a metamaterial that supports a surface plasmon mode in the terahertz region and consequently exhibits strongly enhanced near field thermal radiation at the plasmon resonance frequency. We implemented a finite-difference time-domain method that thermally excites the metamaterial with randomly fluctuating dipoles according to the fluctuation-dissipation theorem. The calculated thermal radiation from the metamaterial was then compared with the case of optical excitation by the plane wave incident on the metamaterial surface. The optical excitation couples only to the mode that satisfies the momentum matching condition while thermal excitation is not bound by it. As a result, the near field thermal radiation exhibits substantial differences compared to the optically excited surface plasmon modes. Under thermal excitation, the near field intensity at 1 µm away from metal surface of the metamaterial reaches a maximum enhancement of 43 fold over the far field at the frequency of the Brillouin zone boundary mode while the near field intensity under optical excitation reaches a maximum enhancement of 24 fold at the frequency of the Brillouin zone center mode. In addition, the peak near field intensity under thermal excitation shows a 4-fold enhancement over blackbody radiation with linear polarization radiation in the far field. The ability to precisely predict the local field intensity under thermal excitation is critical to the development of advanced energy devices that take advantage of this near field enhancement and could lead to the development of new generation of novel energy technology.
© 2017 Optical Society of America
Thermal radiation has been an intriguing topic over the past century. The study of blackbody radiation stimulated the development of quantum theory, which opened a new era of human history. The advance of nanotechnology in the past few decades has enabled researchers to study thermal radiation in the near field which exhibits many new phenomena, including coherent thermal radiation , surface wave excitation by a thermal source , and near-field radiation exceeding the blackbody radiation limit [3,4]. This new research frontier has inspired a variety of applications, including thermal energy conversion , radiative cooling , nano-fabrication [7,8], and near-field imaging .
The thermal radiation research in nanostructures calls for a rigorous and efficient way of numerical modeling. In some studies , optical excitation by the far field plane wave is used as the excitation source. For studying far field thermal radiation, it is sufficient and computationally efficient to calculate absorption under far field illumination and invoke the Kirchhoff’s law of thermal radiation. But the far field optical excitation method has the fundamental limitation that it can only excite the modes that satisfy the momentum matching condition. Modes that are not excited by optical excitation but that can be excited thermally could lead to significant deviations of the near field thermal radiation from the results under optical excitation. Even for the modes that are excited optically, the coupling efficiency may vary from geometry to geometry and is generally different from the case of thermal excitation. To accurately simulate the thermal radiation, therefore, one must simulate the radiation resulting from thermal fluctuations using the fluctuation-dissipation theorem . This approach has been implemented in both time domain [11,12] and frequency domain  to directly calculate thermal radiation from complex structures but these techniques have not yet been fully utilized to investigate the near field thermal radiation in artificially structured materials.
In this paper, we present the results of numerical study of the near field thermal radiation from a simple metamaterial thermal source using both far field plane wave optical excitation and thermal excitation. The metamaterial is designed to support surface plasmon polariton (SPP) modes at 5 THz at the Brillouin zone center and 4.6 THz at the Brillouin zone boundary, respectively. In the optical excitation model, only the Brillouin zone center mode can be excited, and exhibits a near field intensity enhancement peak at 5 THz with an enhancement factor of 24. However, the actual thermal radiation extracted from the thermal excitation model reveals a peak near field intensity enhancement at 4.6 THz with an enhancement factor of 43. From the thermal excitation simulations, we also calculate the energy density radiating from this metamaterial thermal source at a temperature of 600 K. The near field intensity was found to exceed the same temperature blackbody radiation intensity by approximately a factor of 4 in the frequency band of the surface plasmon mode. Our study confirms that a metamaterial structure supporting surface plasmon mode can exhibit near field thermal radiation intensity exceeding the blackbody radiation. It also shows the differences between optical and thermal excitation modeling and thereby the need for implementing direct thermal modeling for quantitative study of near field thermal radiation phenomena and their device applications.
Results and discussion
To design a 1D periodic metamaterial thermal source, we first use the transverse magnetic (TM) polarized far field plane wave as the excitation source, and tune the reflectance spectrum so that the dip occurs at 5 THz. We choose this frequency band with an eye on eventual application in the rectenna technology . A rectenna directly rectifies the ac electric field of an incoming electromagnetic wave and produces dc electricity. One of the main challenges for the energy harvesting application is that the rectenna efficiency tends to be low at high operating frequencies while the thermal radiation intensity is low at low frequencies. A metamaterial thermal emitter offers a solution by significantly increasing the thermal radiation intensity in the near field at low frequencies. In this context, 5 THz is considered a good operating frequency. The simulation is conducted by the finite-difference time-domain (FDTD) method using the commercial software Lumerical. The metamaterial structure schematically shown in Fig. 1(a) consists of 1D periodic trenches in copper. The trenches are filled with SU-8 and there is also a thin layer of SU-8 on the metal surface. In the FDTD simulation of far field optical excitation, we set up a two-dimensional (2D) computational cell containing one unit cell along the direction of periodicity and terminate with periodic boundary condition (PBC). Above the metamaterial structure is an air region, which is terminated with the perfect matched layer (PML) boundary condition that absorbs all out-going waves exiting the computational domain. As for the material properties, we use the Drude dielectric function for copper , and a constant refractive index of 1.7 for SU-8 . In the final design, the period is 24 μm with a trench depth of 7.3 μm and width of 2 μm. The SU-8 overlayer is 1 μm thick. In the future integration with a photovoltaic (PV) converter, the PV device will be placed on the SU-8 overlayer and therefore we evaluate the field in the following discussion is at the SU-8 surface or 1 μm above the metal surface.
The reflectance spectrum of the metamaterial under TM polarized plane wave excitation is shown in Fig. 1(c). It has a strong resonance at 5 THz with reflectance suppressed down to 20%. This indicates good coupling between the plane wave excitation source and the surface plasmon mode. To further confirm the excitation of the surface plasmon mode, we calculate the photonic band structure of the metamaterial also using the FDTD method. In this calculation, we place multiple random electric dipoles near the metamaterial surface, and monitor the modes that persist long after the dipoles are turned off. The result is shown in Fig. 1(d). The wave number ky is for the direction along the periodicity (y direction as shown in Fig. 1(a)). The maximum value of ky is limited to the first Brillouin zone edge, , where = 24 μm is the periodicity. The straight line, which bounces back at the Brillouin zone boundary, is the light line in free space. The flat line that crosses the y-axis at 5 THz is the surface plasmon mode and the y-intercept corresponds to the resonance dip observed in the reflectance spectrum under normal incidence or ky = 0. At the Brillouin zone edge, the surface plasmon mode is located at 4.6 THz. The curve overlaps with the light line at low frequencies, and bends to the right side at higher frequencies until it reaches the resonance frequency. This surface plasmon polariton branch below the light line has the ability to produce strongly localized optical modes at desired frequencies and is sometimes called spoof surface plasmon . Near the resonance frequency, the slope of the dispersion curve decreases close to zero, indicative of a large local density of states and local electromagnetic field. However, due to the momentum mismatch, these modes cannot be excited by coupling light from free space. Thus, no resonance dip is observed around the spoof surface plasmon frequency of 4.6 THz in the reflectance spectrum in Fig. 1(c).
To investigate the near field thermal radiation in a rigorous and quantitative manner, we implemented a direct thermal modeling method based on fluctuation-dissipation theorem following the approach proposed by Luo et al [11,18,19]. Briefly, we introduced a random thermal fluctuation term into the copper polarization response P to local electric field E20], the time correlation of in frequency domain is given as follows:11, 19]. has the uniform distribution with a range, where N is the FDTD time steps, and is the mesh volume element used in the simulation. The final simulated thermal radiation is calibrated by requiring the far field frequency spectrum to be consistent with the Kirchhoff’s law. To implement this scheme on the commercial FDTD software Lumerical, we worked with the Lumerical engineers to develop a new material model using their “flexible material plugin framework”. The new material model simulates thermal radiation by adding random perturbation to the polarization at each FDTD step. Because the random excitation term should be implemented in 3-dimensional (3D) space, we used 3D FDTD modeling even though the structure is invariant in z direction. In the z direction, we extended the length to 2 μm, and set boundary condition to continuity. The cross-section on the x-y plane is shown in Fig. 1(b). To faithfully represent the random thermal fluctuations, we set up five simulations with different Bloch phase shifts along the y direction. The phase shift was varied at uniform interval from 0 to, where Λ = 24 μm is the period of the metamaterial structure. For each Bloch phase shift, we averaged the results of seven simulations with different random number seeds to reduce the stochastic noise inevitable in this method . The final result is an average of all simulations with different Bloch phase shifts and random number seeds used to generate the term. Taking into consideration the computational time and the FDTD simulation frequency resolution, we chose the FDTD time long enough to reach the temporal steady state and to have 0.2 THz resolution in frequency. The simulation is first tested for thermal emission from a flat aluminum surface for which an analytical solution exists. The simulation showed an excellent agreement with the analytical result. The details are presented in the Appendix.
For quantitatively comparisons between the optical and thermal excitation of surface plasmon mode, we first present the optical excitation results. In the following discussion, the effect of surface plasmon will be assessed by comparing the optical energy density, ε|E|2, and will be loosely termed as intensity enhancement, as it is proportional to the actual intensity in the far field, although one cannot rigorously define intensity for evanescent modes. Figure 2(a) shows the profile of intensity enhancement over the excitation source at 5 THz, the frequency of the surface plasmon at ky = 0. The intensity is largest near the trench and decays away from the metal surface and also from the trench center. At the largest intensity point, the intensity enhancement factor over the incident intensity reaches about 300. Inside the metal, the field is almost zero because copper has a very large negative permittivity and thus behaves like a perfect electric conductor in this frequency range. Since the local field intensity varies strongly along the y-direction, we averaged the field over one unit cell length for comparison’s sake. The averaged intensity enhancement at 1 μm away from the copper surface is shown in Fig. 2(c) from frequency 3 THz to 7 THz. It has a clear peak at 5 THz with an enhancement factor of 24 over the far field intensity. The peak position aligns well with the frequency of the reflectance dip as well as the zone center mode (ky = 0) in the band structure. This is the well-known excitation of surface plasmon mode at ky = 0 via grating coupling.
Similarly to the optical excitation case, we calculated the enhancement of near field optical intensity under thermal excitation. Once again, the near and far field intensities were calculated by averaging the intensity over the unit cell along y direction as shown in Fig. 1(b). The averaged near field intensity was evaluated at 1 μm above the metamaterial surface and the averaged far field intensity at 100 μm away. We confirmed that all near field modes decay sufficiently at 100 μm distance away from the metamaterial. Figure 2(b) shows the field profile of the thermally excited surface plasmon mode at 5.1 THz. Since no visual difference of intensity is observed in the x dimension from 30 μm and farther, the figure shows only the x range from −10 μm to 30 μm. Limited by our FDTD frequency resolution, 5.1 THz was the closest frequency to 5 THz, the zone center mode frequency at ky = 0. It is apparent that the mode shape is the same as that shown in Fig. 2(a) under far field TM polarized plane wave excitation. In both cases, most of the energy is confined near the trench and decays away from it. In direct thermal excitation, the plasmon mode is excited thermally, resulting in a different near field enhancement factor. Figure 2(b) shows strong field with randomness in intensity and spatial distribution in the copper region, which arise from the randomly fluctuating dipoles. In contrast, under far field plane wave excitation, no field resides inside the metal as shown in Fig. 2(a). Also, since the thermal source can couple into all available modes, the intensity enhancement profile in Fig. 2(b) must contain all surface plasmon modes with different ky values and polarizations.
Figure 2(d) shows the averaged near field intensity enhancement over the far field intensity. At the frequency of the zone center mode of 5 THz, the enhancement is around 6 fold. It is much weaker than that in the optical exaction case shown in Fig. 2(c). This means that the far field TM polarized plane wave excitation model overestimates the available thermal energy at the 5 THz zone center mode frequency. In contrast, the largest enhancement happens at the zone boundary mode frequency of 4.5 THz with an enhancement factor of 43. This peak enhancement at 4.6 THz does not show up in the optical excitation model because the far field plane wave excitation source cannot couple to the zone boundary mode due to momentum mismatch. In contrast, a thermal source can excite surface plasmon modes with various ky values and produce large near field thermal radiation over the whole surface plasmon frequency band. As shown in the band structure in Fig. 1(d), the first branch of the surface plasmon dispersion curve lies below the light line with wave number larger than light line. Therefore, the excited modes are all surface waves. The thermal energy coupled into these waves is trapped near the metamaterial surface and eventually gets absorbed back by copper. The mode profile of the zone boundary mode is shown in Fig. 3(b). The mode shape is very similar to the 5 THz zone center mode in Fig. 2(b), but the near field intensity enhancement over far field is much stronger.
All the previous discussions are focused on the relative thermal radiation intensity enhancement of the near field over the far field. It is also important to investigate the actual radiative energy that is produced by the metamaterial thermal source, and compare it with the same temperature blackbody. The blackbody radiation is the upper limit of thermal radiation to the far field from any object. It is ultimately limited by the free space density of state of propagating waves which convey thermal energy into the far field. In the near field, thanks to the additional surface modes, the density of states can exceed that of the free space and consequently the radiative energy exceeds that of blackbody radiation [4,21]. It is stressed that this phenomenon happens only in the near field region and in the far field region the thermal radiation still obeys the Kirchhoff’s law. To make a specific comparison with blackbody, we chose a temperature of 600 K. The thermal excitation we implement here does not explicitly include temperature information in it. We therefore use Kirchhoff’s law of thermal radiation to specify the temperature. Specifically, using the emissivity calculated by the far field optical excitation simulation, we calculate the far field radiation intensity spectrum for the temperature 600 K. By equating this value with the far field intensity extracted from direct thermal excitation modeling, we obtain the frequency dependent scaling factors for this specified temperature over all frequencies of interest. We then use the same scaling factors to calibrate the near field. Once again, we use ε|E|2 as relative intensity to directly compare the optical field in the near and far field regions. In Fig. 3(a), the black dashed line is the blackbody radiation intensity at 600 K. The red line is the far field radiation intensity from the metamaterial thermal source, which is the product of emissivity and the blackbody radiation intensity. Since the emissivity cannot exceed unity at any frequency, the far field radiation intensity is always weaker than blackbody radiation over the whole frequency band. It has a peak at around 5 THz due to the large emissivity arising from the out-coupling of the zone center mode. The blue line is the averaged relative intensity at 1 μm away from the metamaterial surface. It exceeds the blackbody radiation over the frequency band from around 4.2 THz to 5.2 THz due to the excitation of surface plasmon modes. The peak intensity is located near the zone boundary mode frequency of 4.6 THz with approximately 4 times higher intensity than the blackbody radiation. Judging from the intensity enhancement profile in Fig. 3(b), we expect even stronger radiation at positions closer to the copper surface.
We also studied the polarization of the thermal emission from the metamaterial source. Here, we consider only the periodic boundary condition along y direction for the 0th diffraction order, which corresponds to the Brillouin zone center mode (ky = 0) in the band structure. As shown in Fig. 4(b), the far field radiation is strongly polarized in the y direction, especially near the surface plasmon resonance frequency of 5 THz. The intensity in the y polarization direction is over 1000 times stronger than the intensity in the x polarization direction. It indicates that the thermal radiation at this frequency is dominated by the excitation and subsequent out coupling of the surface plasmon mode. This is in good agreement with its reciprocal optical excitation, where a TM polarized source has a very strong coupling with the Brillouin zone center mode and consequently strong absorption at the surface plasmon frequency while a TE polarized source has very poor absorption. The polarization of the thermal radiation in the near field as shown in Fig. 4(a) is very different due to the presence of bound evanescent waves near the metamaterial surface. The ratio of intensity in y polarization over that in x polarization also reaches the peak at around 5 THz due to the Brillouin zone center mode. However, the value of the intensity ratio evaluated at 1 μm above the metal surface is much smaller, ranging between 1.5 and 2.5. This is in excellent agreement with the optical excitation simulations which show the polarization intensity ratio of 2.5 at 5 THz and confirms that the thermally excited plasmon modes retain the same polarization characteristics as the optically excited modes.
In this paper, we present a numerical study on the near and far field thermal radiation by a metamaterial, which can generate near field thermal radiation about 4 times larger than the same temperature blackbody, and emit linearly polarized thermal radiation in the far field. The large thermal radiation originates from the surface plasmon mode supported by the metamaterial. In designing the metamaterial, we used far field TM polarized plane wave excitation to target the mode at 5 THz frequency. The photonic band structure calculation reveals the two branches of surface plasmon modes with the Brillouin zone boundary mode at 4.6 THz and the zone center mode at 5 THz. The surface plasmon mode in the lower branch lying beneath the light line does not show up in the far field plane wave excitation because of the momentum mismatch. However, in thermal excitation, the energy comes from the randomly fluctuating dipoles inside the material, which can excite all available modes. To quantify the differences between these two types of excitation sources, we implemented a FDTD method containing randomly fluctuating dipole sources, and compared the resulting optical intensity profiles between the thermal and optical excitation. We discovered that the thermal excitation can excite all surface plasmon modes irrespective of the momentum matching conditions and produces the maximum field at the frequencies of surface plasmon modes not accessible with optical excitation. The near field enhancement factors are also different from those under optical excitations. This study shows the need for thermal excitation modeling for studies of thermally excited surface plasmons. The implementation of thermal excitation based on the fluctuation-dissipation theorem allows us to quantitatively determine the thermally excited radiative energy from complex artificial geometry in both near and far field regions, making it a useful tool for designing novel thermal sources for energy harvesting devices.
The analytical solution of radiation field from thermal source is limited to only a few simple geometries with high symmetry. To validate our simulation with analytical solution, here we simulate the thermal radiation from a semi-infinite hot aluminum plate using exactly the same simulation method we presented in the manuscript. We will extract the enhancement factor (EF) of near field thermal energy density to the far field, and compare it with analytical solution.
At the equilibrium temperature T, the density of electromagnetic energy is expressed as :
The simulation of this hot aluminum plate is straightforward. We use the same simulation method as described in the original manuscript. In the FDTD computational domain, we include a 50 nm thick (z direction) aluminum plate with size of 500 nm by 500 nm (x-y plane). On top of the aluminum plate (along z direction), we leave 150 nm thick vacuum before the perfect matching layers (PML) boundary. The back of aluminum plate is also terminated with a PML boundary. All other directions are implemented with periodic boundary condition to represent the semi-infinite aluminum plate. For the dielectric function of aluminum, we use the Drude model: rad s−1, rad s−1 . From the simulation, we extract electromagnetic energy density along the z direction at the wavelength of 500 nm or 600 THz. The energy density decays exponentially away from the surface and remains nearly constant after around 100 nm from the aluminum surface. We consider it as the far field energy density, and use it to calculate the enhancement factor of energy density in the near field. The results are shown in Fig. 5. Excellent agreement between the simulation and analytical results are observed, validating our simulation method.
Redwave Energy Inc.; Advanced Research Projects Agency – Energy (ARPA-E) (DE-AR0000676); National Research Foundation (NRF), Ministry of Information and Communication Technologies and Future Planning of Korea (NRF-2016K1A1A2912758).
References and links
2. M. Kreiter, J. Oster, R. Sambles, S. Herminghaus, S. Mittler-Neher, and W. Knoll, “Thermally induced emission of light from a metallic diffraction grating, mediated by surface plasmons,” Opt. Commun. 168(1–4), 117–122 (1999). [CrossRef]
3. O. Ilic, M. Jablan, J. D. Joannopoulos, I. Celanovic, and M. Soljacić, “Overcoming the black body limit in plasmonic and graphene near-field thermophotovoltaic systems,” Opt. Express 20(10), A366–A384 (2012). [CrossRef] [PubMed]
4. L. Hu, A. Narayanaswamy, X. Chen, and G. Chen, “Near-field thermal radiation between two closely spaced glass plates exceeding Planck’s blackbody radiation law,” Appl. Phys. Lett. 92(13), 133106 (2008). [CrossRef]
5. M. Laroche, R. Carminati, and J. J. Greffet, “Near-field thermophotovoltaic energy conversion,” J. Appl. Phys. 100(6), 063704 (2006). [CrossRef]
9. Y. De Wilde, F. Formanek, R. Carminati, B. Gralak, P.-A. Lemoine, K. Joulain, J.-P. Mulet, Y. Chen, and J.-J. Greffet, “Thermal radiation scanning tunnelling microscopy,” Nature 444(7120), 740–743 (2006). [CrossRef] [PubMed]
10. L. Novotny and B. Hecht, Principles of Nano-Optics (Cambridge University Press, 2012).
12. D. L. Chan, M. Soljacić, and J. D. Joannopoulos, “Direct calculation of thermal emission for three-dimensionally periodic photonic crystal slabs,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 74(3), 036615 (2006). [CrossRef] [PubMed]
13. A. W. Rodriguez, M. T. H. Reid, and S. G. Johnson, “Fluctuating-surface-current formulation of radiative heat transfer: theory and applications,” Phys. Rev. B 88(5), 054305 (2013). [CrossRef]
14. G. Moddel and S. Grover, Rectenna Solar Cells (Springer, 2013).
15. E. J. Zeman and G. C. Schatz, “An accurate electromagnetic theory study of surface enhancement factors for Ag, Au, Cu, Li, Na, Al, Ga, In, Zn, and Cd,” J. Phys. Chem. 91(3), 634–643 (1987). [CrossRef]
16. S. Arscott, F. Garet, P. Mounaix, L. Duvillaret, J.-L. Coutaz, and D. Lippens, “Terahertz time-domain spectroscopy of films fabricated from SU-8,” Electron. Lett. 35(3), 243–244 (1999). [CrossRef]
18. D. L. Chan, M. Soljacić, and J. D. Joannopoulos, “Thermal emission and design in one-dimensional periodic metallic photonic crystal slabs,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 74(1), 016609 (2006). [CrossRef] [PubMed]
19. C. R. Otey, L. Zhu, S. Sandhu, and S. Fan, “Fluctuational electrodynamics calculations of near-field heat transfer in non-planar geometrics: a brief overview,” J. Quant. Spectrosc. Radiat. Transf. 132(C), 3–11 (2014). [CrossRef]
20. S. M. Rytov, Y. A. Kravsov, and V. I. Tatarskii, Principles of Statistical Radiophysics: Elements of Random Fields (Springer Verlag, 1989).
21. A. Narayanaswamy, S. Shen, L. Hu, X. Chen, and G. Chen, “Breakdown of the Planck blackbody radiation law at nanoscale gaps,” Appl. Phys., A Mater. Sci. Process. 96(2), 357–362 (2009). [CrossRef]
22. K. Joulain, R. Carminati, J.-P. Mulet, and J.-J. Greffet, “Definition and measurement of the local density of electromagnetic states close to an interface,” Phys. Rev. B 68(24), 245405 (2003). [CrossRef]
23. M. Blaber, M. D. Arnold, and M. J. Ford, “Search for the ideal plasmonic nanoshell: the effects of surface scattering and alternatives to gold and silver,” J. Phys. Chem. C 113(8), 3041–3045 (2009). [CrossRef]