Energy streamlines in near-field radiative heat transfer between hyperbolic metamaterials T. J. Bright, X. L. Liu, and Z. M. Zhang* George W. Woodruff School of Mechanical Engineering, Georgia Institute of Technology, Atlanta, Georgia 30332, USA * [email protected]

Abstract: Metallodielectric photonic crystals having hyperbolic dispersions are called indefinite materials because of their ability to guide modes with extremely large lateral wavevectors. While this is useful for enhancing nearfield radiative heat transfer, it could also give rise to large lateral displacements of the energy pathways. The energy streamlines can be used to depict the flow of electromagnetic energy through a structure when wave propagation does not follow ray optics. We obtain the energy streamlines through two semi-infinite uniaxial anisotropic effective medium structures, separated by a small vacuum gap, using the Green functions and fluctuation-dissipation theorem. The lateral shifts are determined from the streamlines within two penetration depths. For hyperbolic modes, the predicted lateral shift can be several thousand times of the vacuum gap width. ©2014 Optical Society of America OCIS codes: (050.5298) Photonic crystals; (160.1190) Anisotropic optical materials; (160.3918) Metamaterials; (260.2160) Energy transfer.

References and links 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14.

S. A. Ramakrishna, J. B. Pendry, M. C. K. Wiltshire, and W. J. Stewart, “Imaging the near field,” J. Mod. Opt. 50(9), 1419–1430 (2003). R. Wangberg, J. Elser, E. E. Narimanov, and V. A. Podolskiy, “Nonmagnetic nanocomposites for optical and infrared negative-refractive-index media,” J. Opt. Soc. Am. B 23(3), 498–505 (2006). J. Schilling, “Uniaxial metallo-dielectric metamaterials with scalar positive permeability,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 74(4), 046618 (2006). H. Shin and S. Fan, “All-angle negative refraction and evanescent wave amplification using onedimensional metallodielectric photonic crystals,” Appl. Phys. Lett. 89(15), 151102 (2006). D. de Ceglia, M. A. Vincenti, M. G. Cappeddu, M. Centini, N. Akozbek, A. D’Orazio, J. W. Haus, M. J. Bloemer, and M. Scalora, “Tailoring metallodielectric structures for superresolution and superguiding applications in the visible and near-IR ranges,” Phys. Rev. A 77(3), 033848 (2008). A. J. Hoffman, L. Alekseyev, S. S. Howard, K. J. Franz, D. Wasserman, V. A. Podolskiy, E. E. Narimanov, D. L. Sivco, and C. Gmachl, “Negative refraction in semiconductor metamaterials,” Nat. Mater. 6(12), 946–950 (2007). A. Fang, T. Koschny, and C. M. Soukoulis, “Optical anisotropic metamaterials: Negative refraction and focusing,” Phys. Rev. B 79(24), 245127 (2009). H. N. S. Krishnamoorthy, Z. Jacob, E. Narimanov, I. Kretzschmar, and V. M. Menon, “Topological transitions in metamaterials,” Science 336(6078), 205–209 (2012). J. Yao, Z. Liu, Y. Liu, Y. Wang, C. Sun, G. Bartal, A. M. Stacy, and X. Zhang, “Optical negative refraction in bulk metamaterials of nanowires,” Science 321(5891), 930 (2008). M. A. Noginov, Y. A. Barnakov, G. Zhu, T. Tumkur, H. Li, and E. E. Narimanov, “Bulk photonic metamaterial with hyperbolic dispersion,” Appl. Phys. Lett. 94(15), 151105 (2009). A. Orlov, I. Iorsh, P. Belov, and Y. Kivshar, “Complex band structure of nanostructured metal-dielectric metamaterials,” Opt. Express 21(2), 1593–1598 (2013). A. Narayanaswamy and G. Chen, “Thermal emission control with one-dimensional metallodielectric photonic crystals,” Phys. Rev. B 70(12), 125101 (2004). S.-A. Biehs, P. Ben-Abdallah, F. S. S. Rosa, K. Joulain, and J.-J. Greffet, “Nanoscale heat flux between nanoporous materials,” Opt. Express 19(S5 Suppl 5), A1088–A1103 (2011). S.-A. Biehs, M. Tschikin, and P. Ben-Abdallah, “Hyperbolic metamaterials as an analog of a blackbody in the near field,” Phys. Rev. Lett. 109(10), 104301 (2012).

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1112

15. S.-A. Biehs, M. Tschikin, R. Messina, and P. Ben-Abdallah, “Super-Planckian near-field thermal emission with phonon-polaritonic hyperbolic metamaterials,” Appl. Phys. Lett. 102(13), 131106 (2013). 16. Y. Guo, C. L. Cortes, S. Molesky, and Z. Jacob, “Broadband super-Planckian thermal emission from hyperbolic metamaterials,” Appl. Phys. Lett. 101(13), 131106 (2012). 17. Y. Guo and Z. Jacob, “Thermal hyperbolic metamaterials,” Opt. Express 21(12), 15014–15019 (2013). 18. X. L. Liu, R. Z. Zhang, and Z. M. Zhang, “Near-field thermal radiation between hyperbolic metamaterials: graphite and carbon nanotubes,” Appl. Phys. Lett. 103(21), 213102 (2013). 19. X. L. Liu, R. Z. Zhang, and Z. M. Zhang, “Near-field radiative heat transfer with doped-silicon nanostructured metamaterials,” Int. J. Heat Mass Transfer 73, 389–398 (2014). 20. 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-4), 59–112 (2005). 21. Z. M. Zhang, Nano/Microscale Heat Transfer (McGraw-Hill, New York, 2007). 22. K. Park and Z. M. Zhang, “Fundamentals and applications of near-field radiative energy transfer,” Frontiers Heat Mass Transfer 4(1), 013001 (2013). 23. S. Shen, “Experimental studies of radiative heat transfer between bodies at small separations,” Annu. Rev. Heat Transfer 16(1), 327–343 (2013). 24. Z. M. Zhang and B. J. Lee, “Lateral shift in photon tunneling studied by the energy streamline method,” Opt. Express 14(21), 9963–9970 (2006). 25. B. J. Lee, K. Park, and Z. M. Zhang, “Energy pathways in nanoscale thermal radiation,” Appl. Phys. Lett. 91(15), 153101 (2007). 26. B. J. Lee and Z. M. Zhang, “Lateral shifts in near-field thermal radiation with surface phonon polaritons,” Nanoscale Microscale Thermophys. Eng. 12(3), 238–250 (2008). 27. S. Basu, L. P. Wang, and Z. M. Zhang, “Direct calculation of energy streamlines in near-field thermal radiation,” J. Quant. Spectrosc. Radiat. Transf. 112(7), 1149–1155 (2011). 28. S. Basu and Z. M. Zhang, “Ultrasmall penetration depth in nanoscale thermal radiation,” Appl. Phys. Lett. 95(13), 133104 (2009). 29. S. Basu and M. Francoeur, “Penetration depth in near-field radiative heat transfer between metamaterials,” Appl. Phys. Lett. 99(14), 143107 (2011). 30. X. L. Liu and Z. M. Zhang, “Metal-free low-loss negative refraction in the mid-infrared region,” Appl. Phys. Lett. 103(10), 103101 (2013). 31. M. Tschikin, S.-A. Biehs, R. Messina, and P. Ben-Abdallah, “On the limits of the effective description of hyperbolic materials in the presence of surface waves,” J. Opt. 15(10), 105101 (2013). 32. X. L. Liu, T. J. Bright, and Z. M. Zhang, “Application conditions of effective medium theory in near-field radiative heat transfer between multilayered metamaterials,” J. Heat Transfer. in press. 33. S. Barkeshli, “On the electromagnetic dyadic Green’s functions for planar multi-layered anistropic uniaxial material media,” Int. J. Infrared Millim. Waves 13(4), 507–527 (1992). 34. A. Eroglu, Y. H. Lee, and J. K. Lee, “Dyadic Green’s functions for multi-layered uniaxially anisotropic media with arbitrarily oriented optic axes,” IET Microwaves Antennas Propag. 5(15), 1779–1788 (2011). 35. S. M. Rytov, Y. A. Kravtsov, and V. I. Tatarskii, Principles of Statistical Radiophysics 3: Elements of Random Fields (Springer, Berlin, 1989). 36. M. Francoeur, M. Pinar Mengüç, and R. Vaillon, “Solution of near-field thermal radiation in one-dimensional layered media using dyadic Green’s functions and the scattering matrix method,” J. Quant. Spectrosc. Radiat. Transf. 110(18), 2002–2018 (2009). 37. J. E. Sipe, “New Green-function formalism for surface optics,” J. Opt. Soc. Am. B 4(4), 481–489 (1987). 38. J. Lee and J. Kong, “Dyadic Green’s functions for layered anisotropic medium,” Electromagnetics 3(2), 111–130 (1983). 39. X. L. Liu, L. P. Wang, and Z. M. Zhang, “Wideband tunable omnidirectional infrared absorbers based on dopedsilicon nanowire arrays,” J. Heat Transfer 135(6), 061602 (2013). 40. A. Narayanaswamy and G. Chen, “Direct computation of thermal emission from nanostructures,” Annu. Rev. Heat Transfer 14(14), 169–195 (2005). 41. K. Park, S. Basu, W. P. King, and Z. M. Zhang, “Performance analysis of near-field thermophotovoltaic devices considering absorption distribution,” J. Quant. Spectrosc. Radiat. Transf. 109(2), 305–316 (2008). 42. L. P. Wang, S. Basu, and Z. M. Zhang, “Direct and indirect methods for calculating thermal emission from layered structures with nonuniform temperatures,” J. Heat Transfer 133(7), 072701 (2011). 43. D. M. F. Chapman, “Using streamlines to visualize acoustic energy flow across boundaries,” J. Acoust. Soc. Am. 124(1), 48–56 (2008). 44. O. A. Godin, “Wave refraction at an interface: Snell’s law versus Chapman’s law,” J. Acoust. Soc. Am. 125(4), EL117–EL122 (2009). 45. H. F. Schouten, T. D. Visser, and D. Lenstra, “Optical vortices near sub-wavelength structures,” J. Opt. B Quantum Semiclassical Opt. 6(5), S404–S409 (2004). 46. M. V. Bashevoy, V. A. Fedotov, and N. I. Zheludev, “Optical whirlpool on an absorbing metallic nanoparticle,” Opt. Express 13(21), 8372–8379 (2005). 47. S. Basu, B. J. Lee, and Z. M. Zhang, “Infrared radiative properties of heavily doped silicon at room temperature,” J. Heat Transfer 132(2), 023301 (2010).

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1113

48. E. D. Palik, Handbook of Optical Constants of Solids (Academic Press, San Diego, 1998). 49. S. Lang, M. Tschikin, S.-A. Biehs, A. Yu. Petrov, and M. Eich, “Large penetration depth of near-field heat flux in hyperbolic media,” Appl. Phys. Lett. 104(12), 121903 (2014).

1. Introduction Hyperbolic metamaterials are artificial nanostructures that support radiative modes with extremely large in-plane (lateral) wavevectors [1–10]. Their anisotropic properties create a hyperbolic dispersion that deviates from the usual elliptical dispersion of uniaxial medium. This unique property leads to potential applications in optical imaging systems, waveguides, near-field filters, and thermal radiation sources, which can be engineered by modifying the structural dimensions [5–12]. The energy exchange in the near field between two semi-infinite anisotropic media separated by a vacuum gap has also been investigated [13–19]. Both surface modes and hyperbolic modes can enhance radiative transfer through tunneling of photons with large lateral wavevectors [18–22]. Significantly enhanced near-field radiative transfer at nanometer distances has been demonstrated by recent experiments conducted worldwide as reviewed in [22,23]. The potentially large penetration depth and lateral shift of photons associated with these structures can lead to both lateral and vertical characteristic dimensions being quite large. While the lateral shift of energy streamlines and the penetration depth have been investigated for near-field radiative transfer between isotropic materials [24– 29], the characteristic lateral dimensions at which three-dimensional effects become important for hyperbolic modes in these structures are unknown. The present study focuses on the pertinent characteristic dimensions of such structures with a nanometer vacuum spacing between two metallodielectric photonic crystals, by analyzing the energy streamlines using the method of Green’s functions and the fluctuationdissipation theorem. The energy streamlines trace the direction of the Poynting vector inside the structure to depict the propagation of evanescent waves that cannot be obtained simply by ray optics (i.e., geometric optics). The characteristic dimensions are then determined based on the overall lateral displacement of the electromagnetic energy within a few penetration depths. 2. Near-field radiative transfer in uniaxial media Typically, uniaxial metamaterials are constructed from alternating metallic and dielectric layers or from metallic nanorods arrays surrounded by a dielectric matrix [5–11]. The term metallic here refers to a material whose dielectric function has a negative real part in a certain frequency region. For a multilayered structure, if the period is small enough that the electric displacement does not vary significantly over several layers (e.g., the period is subwavelength), the structures can be treated as an effective homogeneous medium with anisotropic properties [2, 3, 6, 10, 30]. This criterion may not be sufficient for near-field radiative transfer if the presence of surface polaritons between the metallic layers becomes important [31]. In the near field, due to the large transverse wavevector, the requirement is that the period be much smaller than the vacuum gap distance, which is usually much smaller than the characteristic wavelength. Detailed discussions and quantitative criteria are given in a separate publication by the present authors [32]. Consider two of the structured metamaterials separated by a vacuum gap of distance d as shown in Fig. 1. The dielectric tensor in Cartesian coordinates (x,y,z) for a uniaxial medium with a vertical optical axis is given by [33, 34]

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1114

Fig. 1. Illustration of near-field radiative transfer between two hyperbolic metamaterials whose optical axis is normal to the surface. A vacuum gap d separate the two semi-infinite media. Cylindrical coordinates are used in which ρ specifies the component parallel to the plane. The wavevector k is decomposed into a parallel component β and a perpendicular component γ. The source (medium 1) is at 300 K and the receiver (medium 2) is at 0 K.

ε t   ε = ε t I − ( ε t − ε z ) zz =  0  0

0

εt 0

0 0  ε z 

(1)

where Ι is the unit dyad. In most of the following discussion, cylindrical coordinates are adopted for convenience, ρ represents the radial coordinate in the x-y plane. The dielectric tensor remains diagonal with two transverse components being equal to εt because of the rotational symmetry. For plane waves in a uniaxial anisotropic medium, the dispersion relations are [3, 9]

γ o2 + β 2 = ε t k02 , for ordinary waves

(2)

γ e2 β 2 + = k02 , for extraordinary waves εt ε z

(3)

Here, k0 = ω / c (with c being the speed of light) is the wavevector in vacuum, and β and γ are respectively the in-plane and vertical components of the wavevector in the medium. Note that Eq. (2) is for an ordinary wave (s polarization) for which the electric field lies in the plane of the stratified layers; while Eq. (3) is for extraordinary waves (p polarization) for which the magnetic field lies in the x-y plane. The subscripts o and e denote properties related to ordinary or extraordinary waves, respectively. The phase-matching boundary condition requires the parallel (transverse) wavevector β to be the same in all layers [21]. For ordinary waves, the uniaxial medium behaves the same as an isotropic material with a dielectric function ( ε t ). Care must be taken in dealing with extraordinary waves. When the spacing (gap width) d between the two media is comparable to the characteristic wavelength of thermal radiation, near-field effects become important. Coupled evanescent waves can enhance the heat flux beyond the far-field limit. Thermal radiation in both the near and far fields as described by the fluctuation-dissipation theorem is due to random charge fluctuations (or currents) inside of a medium at non-zero temperatures [35]. The fluctuational currents in the source region at location r ′ generate an electric field at a location r . A timevarying magnetic field is also generated according to Maxwell’s equations. In the frequency domain, the fields due to a current source can be evaluated by integrations with the assistance

of the Green functions G and Γ as follows [20, 35]:

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1115

E ( r , ω ) = iωμ0  G ( r , r ′, ω ) J r ( r ′, ω ) d 3r ′

(4)

H ( r , ω ) =  Γ ( r , r ′, ω ) J r ( r ′, ω ) d 3r ′

(5)

V′

and V′

where Ε and H are the electric and magnetic field vectors, ω is the angular frequency, V ′ is the volume containing the source, and μ0 is the vacuum permeability since only nonmagnetic media are considered here. The cross-spectral density of the random currents, J r , for an anisotropic media are described by [20, 35]:

Θ (ω , T ) ε 0ε ik′′ (ω ) ω

J r,i ( r, ω ) J r,* k ( r ′, ω ) =

π

δ ( r − r′)

(6)

where Θ is the mean energy of a Planck oscillator, ε 0 is the vacuum permittivity, ε ik′′ is the imaginary part of the dielectric tensor component, and δ is the Dirac delta function. The Green function for the magnetic field Γ is related to that for the electric field G by Γ = ∇ × G for a nonmagnetic medium. The Green function G can be expressed by the Weyl

transform with the Fourier component g according to [13, 36] G ( r , r ′, ω ) =

1 2π



0 g ( β , z, z ′, ω ) e

i β ( ρ − ρ ′)

βdβ

(7)

where the integral is over the wavevector plane considering axial symmetry. Maxwell’s equations can be solved by assuming forward and backward plane waves in each layer and using the Fresnel coefficients to obtain g [37]. The Fourier component for Γ is h , and their relationship is similar to Eq. (7) between G and g . Furthermore, h = ik × g . While the Green functions for a current source embedded in a uniaxial anisotropic media with arbitrary number of layers were obtained previously [33, 34, 38], they have not been applied to nearfield radiation between anisotropic media. For a three-layer structure descripted in Fig. 1, g in each medium is expressed in the following: ± g1 ( β , z , z ′, ω ) =

(

)

i z − z′ iγ z ′− z ± iγ ˆˆ e o,1 ( ) + Rs e o,1 ( ) oo 2γ o,1

+ iF ( β , ω )

g 0 ( β , z, z ′, ω ) =



Ts e o,1 2γ o,1 t02,s i

z′

(

(e

+ iF ( β , ω ) Z10 ( β , ω ) g 2 ( β , z , z ′, ω ) =

i 2γ o,1

Ts e

iγ o ,1 z ′ iγ o ,2 ( z − d )

e

z − z′ ± iγ eˆ1± eˆ1± e e,1 ( )

iγ 0 ( z − d )

Tp e

iγ e,1 z ′

t02,p

+

iγ z′− z Rp e e,1 ( ) eˆ1− eˆ1+

(8)

)

)

ˆˆ + r02,s e−iγ 0 z oo

(

eiγ 0 ( z − d ) eˆ 0+

+ r02,p e−iγ 0 z eˆ 0−

)

(9)

eˆ1+

ˆ ˆ + iF ( β , ω ) Z12 ( β , ω ) Tp eiγ e ,1 z eiγ e ,2 ( z − d ) eˆ +2 eˆ1+ (10) oo ′

where

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1116

F (β , ω ) =

2 k02 ( ε t ,1 + ε z ,1 ) − ke,1

2γ e,1k02ε z ,1

and Z10 ( β , ω ) = k0

2 γ e,1 γ2 β2 β2 + 2 ; Z12 ( β , ω ) = e,2 + 2 2 2 ε t ,1 ε z ,1 ε t ,2 ε z ,2

2 γ e,1 β2 + 2 2 ε t ,1 ε z ,1

Here, T and R are the transmission and reflection coefficients from medium 1 to medium 2 considering the three-layer structure of Fig. 1 [21], rik and tik are Fresnel’s reflection and transmission coefficients at the interface between media i and k [39], and subscript s or p denotes s or p polarization, respectively. When z is inside the source region (medium 1), there exists a primary component that describes waves from the source z ′ that reach location z directly, without being reflected or refracted by any interface. This primary wave travels forward (+) when z > z ′ and backward (–) when z < z ′ . Equation (8) contains superscript ± that should be chosen such that the integration of Eq. (7) uses the positive sign from −∞ to z and the negative sign from z to 0 to properly account for the forward and backward primary waves, respectively. Since the temperature of the receiver is assumed to be absolute zero, there are no primary waves in medium 2. The extraordinary vector (p polarization) for medium l, where l is 1 or 2, is given in cylindrical coordinates as eˆ l± =

γ e,l ε z ,l ρˆ + βε t ,l zˆ

(γ e,l ε z ,l ) + ( βε t ,l ) 2

(11)

2

For the vacuum region, the direction of the electric field for p polarization is eˆ 0± =

γ 0ρˆ + β zˆ k0

(12)

The ordinary wavevector in cylindrical coordinates is simply oˆ = −θˆ . The calculation of the spectral Poynting vector for isotropic media is described in [36, 40– 42]. For a prescribed frequency and lateral wavevector, the Poynting vector in the vertical and radial direction can be obtained from S z ( z, ω, β ) = S ρ ( z, ω, β ) =

k02 Θ (ω, T ) β 2π

3

k02 Θ (ω, T ) β 2π

3

(

)

(

)

0 * * * Re i  ε t′′,1 g11h21 dz ′ (13) + ε z′′,1 g13h23 − ε t′′,1 g 22 h12  −∞  0 * * * Re i  ε t′′,1 g 22 h32 dz ′ (14) − ε t′′,1 g31h21 − ε z′′,1 g33h23  −∞ 

where * represents the complex conjugate. The tensors g and h are for the region (0, 1 or 2) where z is located. Subscripts for g and h indicate the tensor component. It can be shown that the terms involving g22 in Eqs. (13) and (14) are for s polarization, while the rest are for p polarization. For s polarization, the near-field enhancement is rather limited, because it cannot support hyperbolic modes with high β values. All the calculations below are for p polarization only, while the s polarization terms are included for completeness. The energy streamlines are determined by tracing a line tangent to the Poynting vectors inside of the medium [24–27], viz.

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1117

m=



(15)

Sz

where m = tan(θs ) , with θs being the polar angle of the Poynting vector. The energy streamlines are calculated by finding the line tangent at discrete points in the ρ-z plane and then tracing along the points to generate a curve. The concept of energy streamlines has been applied to acoustic waves and optical diffraction [43–46]. A structure must have dimensions which are much larger than the characteristic lateral displacement in order for it to be treated as one-dimensional. If the lateral dimensions are too small, near-field radiative heat transfer may be limited since the presence of larger lateral wavevectors may be limited or altered by the dimensions of the structure. The method described above can be readily extended to multilayered media with each medium being either isotropic or uniaxial. In this case, the Green function in each sub-layer must be determined and Eqs. (13) and (14) need to be integrated within each single layer and then summed up for all the source layers. It should be noted that S z is independent of z in vacuum and decays exponentially in the receiver. In the emitter, it increases as z approaches the interface. The behavior of the zcomponent Poynting vector allows the definition of a near-field penetration depth [28, 29]. The heat flux between the two uniaxial media can be obtained by integrating S z in vacuum over the wavevector plane and over all frequencies, giving ′′ = qnet

1 4π



2

0



Θ (ω , T1 ) − Θ (ω , T2 )  d ω  β 0

 ξ j (ω , β ) d β

(16)

j = s, p

Here, ξ j indicates the power transmission factor for either s or p polarization. For two media with the same dielectric tensor, the transmission factor is given as [13, 16, 18]

(

)

2 2  2iγ 0d 2 2 − 1 r 1 − r10, , j je 10,  ξ j (ω , β ) =  4  Im( r )  2 e 2iγ 0d 1 − r 2 e 2iγ 0d 2 , 10, j  10, j  



β < k0   β > k0  

(17)

where r10, j is the reflection coefficient for j polarization at the interface between the uniaxial media and vacuum, using the anisotropic Fresnel coefficients [39]. When β > k0 the waves become evanescent in the vacuum region. Biehs et al. [13] showed that the Green’s function method is equivalent to the approach based on scattering theory. However, the Fourier components of Green’s function in all media are needed in order to obtain the streamlines. 3 Energy streamlines in hyperbolic metamaterials

Consider a periodic multilayered structure consisting of alternating metallic and dielectric layers. Here, the metallic layer is either doped Si (D-Si) or SiC, whose dielectric function becomes negative at certain frequencies due to free carriers or lattice vibrations (or phonons). Ge is chosen as the dielectric layers for both structures. Assuming that the layer thicknesses are sufficiently thin for the effective medium theory to be valid, the two components of the dielectric tensor in Eq. (1) can be expressed as follows [3, 11]

ε t = f ε m + (1 − f ) ε d

#209344 - $15.00 USD (C) 2014 OSA

(18)

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1118

εz =

ε dε m

f ε d + (1 − f ) ε m

(19)

where f is the volume filling fraction of the metal component, and subscript m or d refers to bulk metal or dielectric, respectively. In the present study, the dielectric function of n-doped Si with a doping concentration of 1020 cm−1 is taken from Basu et al. [47]. The optical constants of SiC are taken from Palik’s handbook [48]. The dielectric function of Ge is taken to be a constant, ε d = 16 , in the region of interest [48]. Note that room-temperature properties are used for these materials, even though the temperature of the receiver is set to 0 K. In all the calculations, medium 1 and medium 2 as shown in Fig. 1 are assumed to have the same dielectric properties. The real part of the dielectric function for the two types of multilayered structures is shown in Fig. 2. The filling ratio for D-Si is 0.5 and that for SiC is 0.3. In both cases, two hyperbolic bands are shown in the shaded regions: I for ε t′ > 0 > ε z′ and II for ε z′ > 0 > ε t′ . Negative refraction or backwards bending of the Poynting vector can only exist with type I hyperbolic dispersion, but both hyperbolic modes can support propagating waves with large β that may enhance near-field radiative transfer [14, 18, 19]. The hyperbolic bands in the D-Si are due to free electron absorption, and thus the type II band appears fairly broad toward long wavelengths. The hyperbolic modes of SiC are much narrower, as they are caused by the infrared phonon mode associated with the Restrahlen band of SiC [21, 48]. The predicted dielectric functions are used for future calculations with fixed filling ratio for the two multilayered structures, D-Si/Ge and SiC/Ge.

Fig. 2. Real part of the dielectric function predicted by EMT: (a) doped Si and Ge with f = 0.5; (b) SiC and Ge with f = 0.3. The multilayered structure is illustrated by the inset for each case. These filling ratios are used in the property calculations throughout this paper.

The power transmission factor for extraordinary waves between two identical semi-infinite metamaterials is shown in Figs. 3(a) and 3(b), which corresponds to a D-Si/Ge and SiC/Ge multilayered structure, respectively. The contour plots show ξp versus frequency and in-plane wavevector β, divided by k0, for a vacuum gap distance of 10 nm. Figures 3(c) and 3(d) show the spectral heat flux corresponding to the two structures. As seen from Fig. 3(a), ξp is broad in both type I and type II hyperbolic bands, although the values in the type II region are not as high as those in the type I region. Furthermore, due to loss, the distinction between the regions is not so clear. When ε z′ > 0 > ε t′ , evanescent waves exist in the metamaterials when

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1119

β < βcr = k0 ε z as can be seen from Eq. (3). Thus, the transmission factor is rather small for the type II band below the critical wavevector β cr [18]. Interestingly, the spectral heat flux for the D-Si/Ge structure shows a broad continuous peak, which can be seen in Fig. 3(c), presumably due to loss. The type I and type II hyperbolic bands can also be seen for SiC/Ge structures as shown in Fig. 3(c), although the magnitude of ξp is much higher and extends to greater β /k0 values. The spectral heat flux for D-Si/Ge multilayers is several orders higher than that for blackbodies in the far field shown in Fig. 3(c) as the dashed line. For the SiC/Ge metamaterial, the distinction between the hyperbolic bands turned out to be clearer and two narrowband peaks can be seen corresponding to type I and type II bands in the spectral heat flux curve shown in Fig. 3(d). The frequencies of the two peaks correspond well to the transverse and longitudinal phonon modes in SiC [21]. Compared with the total heat flux of 460 W/m2 between blackbodies, the near-field heat flux at d = 10 nm is 640 and 1270 times higher with D-Si/Ge and SiC/Ge multilayered structures, respectively.

Fig. 3. Power transmission factor and spectral heat flux with a gap spacing d = 10 nm. Contour plots of the transmission factor for p-polarization between two semi-infinite hyperbolic metamaterials made of (a) D-Si/Ge and (b) SiC/Ge multilayered structures; (c) Spectral heat flux between D-Si/Ge multilayered structures plotted together with that between blackbodies in the far field; (d) Spectral heat flux between SiC/Ge multilayered structures. The hyperbolic bands for the multilayered structures are highlighted. Note that the unit of Sω is heat flux [W/m2] per frequency [rad/s] interval.

To study the lateral displacement, it is noted that the energy streamlines depend on the particular frequency and given parallel wavevector. For the SiC/Ge structure, the two peak frequencies at 1.82 ×1014 and 1.54 ×1014 rad/s can be readily picked. A frequency in between at 1.59 ×1014 rad/s is also chosen since it represents elliptical dispersion where the spectral heat flux is much lower. For the D-Si/Ge structure, five representative frequencies are selected in the region of interest from 0.5 ×1014 to 3 ×1014 rad/s (wavelength λ from 37.7 to

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1120

6.3 μm). Two of the chosen frequencies are 0.5 ×1014 and 1.0 ×1014 rad/s, which fall inside the type II hyperbolic region. The third at 1.8 ×1014 rad/s is at the edge of the type II band with ε t′ close to zero. The fourth frequency of choice at 2.7 ×1014 rad/s falls in the type I hyperbolic band. The fifth frequency location is at 2.9 × 1014 rad/s, where ε z′ is near zero. These frequencies were chosen to demonstrate different and unique energy propagation characteristics. The selection of suitable and representative parallel wavevector is illustrated with Fig. 4 for the D-Si/Ge and SiC/Ge structures, respectively. Figures 4(a) and 4(b) plot the z-component of the Poynting vector calculated from Eq. (13) at each selected frequency for the two structures, respectively, as functions of the parallel wavevector. The magnitude corresponds well to the transmission factor shown in Figs. 3(a) and 3(b) for the given frequency. In the type II hyperbolic band, the Poynting vector is small but increases quickly and reaches a peak at β = βcr, where kz = 0. For the elliptical mode in SiC/Ge structure at 1.59 ×1014 rad/s, as shown in Fig. 4(b), the Poynting vector is relatively large and increases

Fig. 4. Poynting vector and cumulative heat flux: (a) z-component of the Poynting vector versus lateral wavevector at select frequencies for (a) D-Si/Ge and (b) SiC/Ge multilayer metamaterials; Integration of the total heat flux over the wavevector space for (c) D-Si/Ge and (d) SiC/Ge. The horizontal line in (c) and (d) shows where 50% of the heat flux falls above and below the median wavevector for each frequency. Note that line styles in (c) and (d) correspond to the line styles and frequencies in (a) and (b), respectively. The unit of Sz is the heat flux [W/m2] per unit frequency [rad/s], wavevector [rad/m], and azimuthal angle [rad].

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1121

with β for small β values. However, drops quickly at β = k0 ε z′ ≈ 6k0 , beyond which only evanescent modes exist in the metamaterial. The enhancement at β > k0 is still due to photon tunneling because the waves are evanescent in vacuum. Figures 4(c) and 4(d) show the cumulative spectral heat flux for the selected frequencies for both structures. The 50% lines intercept with each curve at the median wavevector β median , which can be used as a characteristic value to plot the energy streamlines for the determination of the lateral displacement. Note that for the elliptical mode, as shown in Fig. 4(d), β * ≈ β median / k0 = 4.4 is much smaller than other modes since most of the energy is transferred at β < 6k0 . The penetration depth is an important characteristic of the vertical dimension for the medium to be considered opaque [21]. Furthermore, to evaluate the lateral shift, one needs to choose a certain depth in the medium away from the vacuum gap. For near-field radiation, the power penetration depth can be calculated according to δ = 1 ( 2γ e′′ ) , where γ e′′ is the imaginary part of the z-component wavevector for the extraordinary wave, since the contribution of the ordinary wave is negligibly small [27, 28]. The penetration depth at each of the chosen frequencies is shown in Figs. 5(a) and 5(b) for the two structures, respectively.

Fig. 5. Penetration depth δ (ω , β ) and spectral penetration depth δω (ω ) for d = 10 nm: (a,b)

δ / d for selected frequencies as a function of β/k0 for D-Si/Ge and SiC/Ge, respectively; (c,d) δω / d as a function of frequency for D-Si/Ge and SiC/Ge structures, respectively. Note that the hyperbolic bands are highlighted in (c,d).

#209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1122

In the type I hyperbolic band, δ decreases with β. A sharp peak in some of the curves corresponds to the critical wavevector in the type II hyperbolic band, where γ e′′ has a minimum. At 1.8 × 1014 rad/s, where ε t′ is near zero for the D-Si/Ge structure, the penetration depth is nearly a constant at small β , but decreases at large values of β . At 1.59 × 1014 rad/s for the SiC/Ge structure, the penetration depth is very high for small β values due to the existence of propagating waves with low loss. For large values of β, the lateral distance that the electromagnetic energy might travel in the medium may be very large because the slab acts as a waveguide. However, due to loss and the relation between γ e and β through Eq. (3), a large lateral wavevector results in a smaller penetration depth. At β >> k0, the penetration depth for an isotropic medium is given by δ = 1 ( 2 β ) [27, 28]. For a uniaxial medium, it can be shown that the asymptote of the penetration depth is

δ (ω ,β ) =

± Re

(

1

εt ε z

)

1 2β

Note that δ should always be positive and the negative sign should be used if Re

(20)

(

εt ε z

)
d 

#209344 - $15.00 USD (C) 2014 OSA

(22)

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1123

such that inside of the media an increment of one is a distance of one penetration depth. The penetration depth δ(ω,β) can be found from Fig. 5. The lateral shift can be estimated according to the energy streamlines within two penetration depths since only approximately 14% of the energy flux remains beyond 2δ. For all the streamlines, the origin (z = 0 and ρ = 0) is defined at the interface between the emitter (source) and the vacuum region. At a given frequency and lateral wavevector, the lateral displacement is defined as Δ (ω , β ) = ρ ( z 2 ) − ρ ( z1 )

(23)

where ρ ( z ) is the lateral location of the energy streamline at z. Note that z = 0 corresponding

to the interface between the emitter and the vacuum gap, and ρ ( 0 ) is defined to be at the

origin. For the lateral shift in the vacuum gap, z1 = 0 and z2 = d ; for the lateral shift in the emitter, z1 = 0 and z2 = 2δ (two penetration depths inside of the emitter); for the receiver,

z1 = d and z2 = 2δ + d . Negative penetration depth and lateral shift are allowed. Figure 6(a) shows the streamlines for three frequencies for D-Si/Ge structure in the emitter, while the inset shows the streamlines in the vacuum gap and the receiver. Very large lateral displacements are found inside the emitter especially at the two lowest frequencies in the hyperbolic type II bands, where the lateral displacement of the streamline at two penetration depths is Δ / d = 3730 at ω = 0.5 × 1014 rad/s and 1600 at ω = 1.0 × 1014 rad/s. At ω = 1.8 × 1014 rad/s in the source region, the lateral shift is near 200d. As shown in Fig. 6(b), the type I hyperbolic mode (2.7 × 1014 rad/s) undergoes a negative refraction in terms of the Poynting vectors and Δ ≈ 23d inside the emitter. At ω = 2.9 × 1014 rad/s where ε z′ = 0 , Δ is relatively small (about 3.7d). Usually, the lateral shift is the largest inside the emitter, while the shift inside the vacuum gap and receiver could also be significant in some cases. Similar results are shown in Fig. 6(c) for SiC/Ge structure, where the lateral shift for the type II mode becomes extremely large, close to 4900d at ω = 1.54 × 1014 rad/s. For propagating waves at ω = 1.59 × 1014 rad/s, Δ / d = 433 in the emitter. Negative refraction can be seen from Fig. 6(d) for the type I mode at ω = 1.82 × 1014 rad/s with a moderate Δ around 95d. Considering the large lateral shift for hyperbolic type II modes, the lateral dimension of the multilayer structures must be sufficiently large in order for the near-field radiative transfer to be treated as one-dimensional using the multilayer Green’s functions. Note that except ω = 1.59 × 1014 rad/s for the SiC/Ge structure, the selected frequencies all have a relatively high z-component of the Poynting vector, as shown in Figs. 3(c) and 3(d). Therefore, the large lateral displacements at these frequencies are relevant to the requirement of the lateral dimension of the structures in order for them to be considered one-dimensional. Furthermore, the penetration depth will scale with the vacuum gap distance as well as the lateral displacement because the transmission coefficient varies with βd in the near field. A weighted average of the lateral displacement can be defined using the energy streamline method, considering the transmission coefficient, according to ∞

Δω

 (ω ) = 0

βξ p Δ (ω , β )) d β ∞

0

βξ p d β

(24)

The lateral displacements across the vacuum gap and inside the emitter are plotted in Fig. 7 as a function of frequency. The lateral displacement varies significantly at different frequencies. For D-Si/Ge, the lateral displacement in the emitter is negative for 2.3 × 1014 rad/s < ω < 2.8 × 1014 rad/s, although the magnitude is relatively small. This corresponds to negative refraction in the type I hyperbolic band. For the SiC/Ge structure at certain frequencies, the magnitude of the lateral displacement can be very large even in the type I hyperbolic region. For example, |Δω|/d exceeds 4000 at 1.7 × 1014 rad/s; see Fig. 7(b). #209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1124

Fig. 6. Energy streamlines at select frequencies based on the median wavevector for (a,b) DSi/Ge and (c,d) SiC/Ge multilayer structures at d = 10 nm. Here, β * = β median / k0 that corresponds to 50% of the cumulative heat flux for each frequency as shown in Figs. 4(c) and 4(d).

Inside of the emitter or receiver, the lateral displacement is related to the eccentricity of the hyperbolic dispersion, since the Poynting vector is perpendicular to the hyperbola. In the case of negligible loss, the eccentricity for hyperbolic dispersions are given below

ε z′ for type I hyperbolic dispersion ε t′

(25)

ε t′ for type II hyperbolic dispersion ε z′

(26)

ec = 1 +

and ec = 1 +

In the type II region, if the eccentricity is near one, then S z >> S ρ

and the line tangent

m defined in Eq. (15) becomes smaller and the lateral displacement is also expected to be smaller. This agrees with the trend for the D-Si/Ge at frequencies up to 1.8 × 1014 rad/s, where #209344 - $15.00 USD (C) 2014 OSA

Received 1 Apr 2014; revised 14 May 2014; accepted 15 May 2014; published 2 Jun 2014 30 June 2014 | Vol. 22, No. S4 | DOI:10.1364/OE.22.0A1112 | OPTICS EXPRESS A1125

ε t′ ≈ 0 . Although it is expected that the penetration depth will have a maximum at this frequency according to Eq. (20), this did not happen due to loss. The lateral displacement in SiC/Ge between 1.50 × 1014 rad/s and 1.55 × 1014 rad/s in the narrow type II hyperbolic band also agrees with that predicted by Eq. (26).

Fig. 7. Average lateral displacements in vacuum and emitter for (a) D-Si/Ge structure and (b) SiC/Ge structure, with a gap spacing d = 10 nm. The highlighted regions correspond to the hyperbolic band (type I or type II).

The trend of the lateral displacement for type I hyperbolic dispersion is more complicated. According to Eq. (25) and the shape of the hyperbola, see [3] for example, m is large when ec → 1 or ε z′ ε t′

Energy streamlines in near-field radiative heat transfer between hyperbolic metamaterials.

Metallodielectric photonic crystals having hyperbolic dispersions are called indefinite materials because of their ability to guide modes with extreme...
2MB Sizes 0 Downloads 3 Views