1. Introduction

A stack of graphene sheets, separated by subwavelength dielectric spacers, can be regarded as a composite material with uniaxial electric properties under certain conditions. Graphene layers strongly affect the complex effective permittivity of the composite for electric field components polarized parallel to the graphene sheets. Uniaxial anisotropic materials in general offer a variety of interesting electromagnetic properties. In particular, we investigate here a subcategory denoted as hyperbolic metamaterials (HMs), named after the hyperbolic iso-frequency wavevector dispersion that arise due to the negative permittivity experienced by the electric field component along either the axis of anisotropy or a direction orthogonal to the axis of anisotropy [1, 2].

Strong interest in HMs is based on their specific property that enables propagation in a very wide spatial spectrum, that would be otherwise evanescent in free space, which is in principle unbounded for the ideal case of purely hyperbolic iso-frequency wavevector dispersion. In case of realistic hyperbolic-like dispersion, the spatial spectrum allowed for propagation can still be extremely wide, as shown in [3–5]. This property is shared with uniaxial anisotropic materials having elliptic iso-frequency dispersion diagram, with a very large axial ratio. Though in practical cases purely hyperbolic dispersion cannot be obtained, effective medium models based on hyperbolic dispersion proves to be a very useful tool for understanding the physics behind the interesting electromagnetic properties of these metamaterials with extremely subwavelength features.

Recently, metal-dielectric multilayers were proposed as candidates to realize HMs at optical frequencies [3, 6, 7]. The wide spatial spectrum of propagation supported by these HMs can lead to novel phenomena as increasing the power emitted by imposed dipoles [3,8] or scattered by nanoparticles [3,9] on HM surfaces, and this power is mostly directed into HMs. This exotic property enables features like focusing with very subwavelength resolution [10,11], controlling absorption [12], enhancement of spontaneous emission [13], increasing the decay rate of emitters [6], designing quantum and thermal emitters [14]. HMs can also host backwards waves and thus they are utilized for achieving negative refraction as in [15]. As stated in [4], HMs are considered to be promising materials for advancement in the fields of tunable surface plasmon polaritons, super Planckian thermal emission [5], radiative decay engineering [16], and nano-imaging.

HMs attract attention also because they are easy to fabricate using metal-dielectric multilayers or metallic nanowires embedded in dielectric substrates. Also, doped semiconductors [17, 18] and conductive oxides used for generating surface plasmons can be used for HM designs in near- and mid-IR frequency bands [15]. In multilayer structures, metal is used as a negative permittivity layer-spaced by dielectric layers, overall creating a negative permittivity effect for the electric field components tangential to the layers. This effect does not rely on any resonant behavior and thus is a very wide-frequency-band.

In this paper we investigate a multilayer HM design based on stacking dielectric layers and graphene sheets. A graphene sheet (a one-atom-thick carbon layer) is able to support surface plasmon modes at terahertz frequencies [19–23]. A recent work in [24] a discussed HM based on a semi-infinite stack of graphene-dielectric multilayers, studied at a temperature of 4 K, thus assuming graphene losses negligible, leading to purely real permittivity and wavenumbers. Our analysis instead accounts for losses at a room temperature of 300 K and their effect on (i) the effective permittivity that assumes complex values, (ii) the hyperbolic-like dispersion (it is not exactly hyperbolic), hence (iii) wavenumbers that may assume complex values, and (iv) on the radiated power by a dipole near the HM surface, where losses play an important role. Moreover, we investigate practical cases with a finite number of graphene-dielectric layers, and quantitatively show tunability aspects of graphene-based HMs using electrostatic biasing. We also show a detailed study about the dependence of power spectrum emitted by a dipole source in the proximity of the graphene-based HM on the number of layers, as well as on frequency. Alternative to the HM implemented as a graphene-dielectric multilayer studied here, in [25] graphene stripes analogous to a metallic wire medium, are utilized for realizing hyperbolic dispersion in cylindrical wavenumber coordinates, with the aim of designing a hyperlens. The peculiar electronic properties of graphene [26, 27], have been investigated for different electromagnetic applications such as lensing [28], transformation optics [29], nanomechanical resonators [30], and solar cells [31]. Moreover, based upon periodic patterning of graphene, bi-periodic and/or multilayered graphene structures were extensively studied for enhanced transmission [32], optical absorption [33], and tunable metasurfaces in [34], as well as isolators and polarizers in the microwave regime [35–37]. It was shown that crystalline Graphite (the 3D parent of graphene) possess indefinite permittivity at UV frequencies [38].

A graphene sheet has properties at THz frequencies similar to those of a thin metal film at optical frequencies, as shown in our work. In principle, any inductive infinitesimally-thin layer can be used to realize HMs, however, in the terahertz regime, designing inductive layers is difficult due to metallic losses and spatial and frequency dispersion introduced by periodically patterned conductive layers. For this aim, metallic meshes are mainly effective only in lower microwave regime [39–41]. On the other hand, the use of highly dispersive metals is practical the optical frequencies below the plasma frequency. Here we show that stacking graphene sheets can be utilized for designing tunable HMs in a wide frequency spectrum ranging from millimeter-waves up to tens of terahertz frequencies, encompassing the whole far-infrared band. In Section 2, we develop an effective medium approximation (EMA) to facilitate the characterization of graphene-based multilayer structure and its use in possible devices and we assess its limits. Moreover, in Section 2 we describe the basic properties of graphene and the associated hyperbolic wavenumber dispersion of the multilayer structure composed of graphene-dielectric layers. In Section 3 we investigate plane wave transmission and reflection for a thin slab made of several graphene sheets, as well as their tunability features, and we also show how EMA is able to describe these properties. In Section 4 we investigate the radiation of a dipole at the interface of a finite thickness HM and we show that the HM is able to enhance the total power radiated by several orders of magnitude, reporting an enhancement of 5 × 10² at 2 THz. We also show that most of the power is directed into the HM, offering a viable route for wide band and wide incidence-angle super absorption interfaces at far-infrared frequencies, as previously discussed in [3, 9, 13, 42] for optical frequencies. We show that this large enhancement of power emission is associated to the wide spatial spectrum being able to propagate inside the HM, that would be otherwise evanescent in free space.

2. Tunable hyperbolic metamaterial made of graphene-dielectric multilayers

A graphene monolayer is electrically characterized by its surface conductivity σ(ω, μ_c), where μ_c is the chemical potential related to the electrostatic biasing, which quantifies the electronic transport properties [43]. The frequency dependent conductivity follows the interband (bound-electrons) and the intraband (free-electrons) sum rules [44, 45]. Spatial dispersion of graphene has negligible effects, since the graphene lattice constant a ≈0.246 nm is extremely subwavelength at THz frequencies [34]. Therefore, a graphene layer is modeled by the local isotropic surface conductivity σ(ω, μ_c) = σ′ + jσ″ (assuming time-harmonic fields of the form e^jωt) that is calculated by the Kubo formula

Since a graphene sheet is infinitesimally thin with respect to the dielectric thickness, one may assume that ε_z = ε_d, because z–directed electric field would not excite any current in the graphene sheet. The transverse effective relative permittivity ε_t is determined as follows, considering a unit cell made of a dielectric layer with dielectric constant ε_d between z = 0 and z = d, and a graphene sheet at z = 0. Within this unit cell one can write ∇ × H = jωε₀ε_dE + J = jωε₀ε̱_d · E, where the current density J [A/m²] in the graphene sheet is reduced to the surface current along the sheet J = δ(z)σE_t, and E_t is the transverse component of the electric field. Therefore, one has ${\underset{\_}{\mathit{\epsilon }}}_d={\epsilon }_d\underset{\_}{\mathbf{I}}-j\frac{\sigma }{\omega {\epsilon }_0}\delta \left(z\right)\left(\widehat{\mathbf{x}}\widehat{\mathbf{x}}+\widehat{\mathbf{y}}\widehat{\mathbf{y}}\right)$, where I̱ is the unit tensor, that when averaged over a period along z leads to the effective relative “transverse” permittivity ε_t,

The formula for ε_t in Eq. (3) could be obtained alternatively by following the method used for homogenization involving thin metal-dielectric layers [49]. Accordingly, a graphene sheet may be treated as a layer with extremely subwavelength, but finite, thickness with bulk properties. Exploiting the continuity of the electric field, along x and y, at the boundaries between graphene and dielectric layers, and averaging the transverse component of the effective displacement current over a period d also leads to Eq. (3) (here the effective displacement current account for the displacement current in the dielectric and the conduction current in the graphene sheets). In Eq. (3), we have highlighted that the graphene conductivity is strongly dependent on the frequency and chemical potential. It is important to note that if a graphene sheet has a sufficiently large inductive susceptance, i.e., if σ″ < −ωε₀ε_dd, then the effective relative permittivity term ε_t = ε′_t − jε″_t has a negative real part, i.e., ε′_t < 0. Under this condition, and recalling that ε_z > 0, extraordinary waves, with TM^z polarization (magnetic field transverse to z), are allowed to propagate inside the HM, with wavevectors exhibiting iso-frequency hyperbolic dispersion as explained in [3]; whereas ordinary waves with TE^z polarization (electric field transverse to z) are mainly evanescent. This allows for the propagation of TM^z waves with transverse wavenumber $k_t>\omega \sqrt{{\mu }_0{\epsilon }_0{\epsilon }_d}$, that would be otherwise evanescent in free space.

 

Fig. 1 Composite multilayer material made by stacking graphene sheets and dielectric layers. Under certain conditions it exhibits hyperbolic-like iso-frequency wavevector dispersion as depicted in the inset, where v_g indicates the direction of the group velocity.

Download Full Size | PDF

To better illustrate the possible homogenized parameters that can be obtained, we report in Fig. 2 the real and imaginary parts of ε_t = ε′_t − jε″_t using the EMA formula Eq. (3), assuming that graphene sheets are spaced by silica layers with permittivity ε_d = 2.2 and thickness d= 0.1 μm. Note that ε′_t approaches ε_d at higher frequencies due to saturation of the graphene conductivity to its universal value πe²/(2h) [44] which gives a negligible contribution compared to ε_d. Moreover, the imaginary part of the effective permittivity term (ε″_t) is remarkably small in a certain frequency band, showing that waves in this composite material can propagate with minimal attenuation. Indeed, it is known that a graphene layer may support weakly attenuated plasmonic modes in at terahertz frequencies [45]. We note that ε_t is very sensitive to the chemical potential μ_c, and we show that the zero-crossing frequency of ε′_t, occurring when σ″ = −ωε₀ε_dd, occurs at 6.6, 11, and 24.6 THz, for different values of μ_c = 0, 0.1, and 0.4 eV, respectively. The tunability of the proposed HM structure is then quantified, as seen in Figs. 3(a) and 3(b), where the chemical potential and stacking density (thickness of dielectric layer d) are varied assuming a frequency of 12 THz. When the dielectric thickness d is increased, one should note that ε′_t increases toward 0. Note also that as the chemical potential (tuned by electrostatic biasing) increases ε′_t takes smaller values. Hence, at 12 THz, ε′_t is positive (≈2.2) for zero chemical potential, and at higher chemical potential values, ε′_t becomes negative: for example at μ_c =0.5 eV, with d= 0.1 μm, one has ε′_t ≈−9.5. In summary, a composite material made by layered graphene sheets possesses desirable performance in terms of losses, inductive response, and tunability from millimeter-waves up to mid-infrared frequencies, hard to find in any other known material. This makes it a good candidate for realizing HM designs in the THz range.

 

Fig. 2 Effective medium complex relative permittivity term ε_t = ε′_t − jε″_t for biased and unbiased graphene multilayer configuration.

Download Full Size | PDF

 

Fig. 3 Relative effective medium complex permittivity term ε_t = ε′_t− jε″_t of graphene HM versus the graphene sheets’ chemical potential for various spacer thicknesses at 12 THz.

Download Full Size | PDF

When considering plane waves propagating in a homogeneous uniaxial anisotropic medium it is useful to decompose them into the modal polarizations TE^z and TM^z. A description using EMA, when valid, gives a neat physical insight into wave propagation in this structure. Wavevector dispersion relations for ordinary (TE^z) and extraordinary (TM^z) waves in a uniaxial anisotropic materials are written as

With the aim of assessing the validity of EMA in Eq. (4), we calculate the iso-frequency wavevector dispersion with the more accurate Bloch theory for the periodic structures. This is done by treating each graphene layer as a complex lumped admittance Y_s = σ (where the subscript “s” denotes surface) as a shunt load in a transverse equivalent network (TEN, see Chapters 2 and 3 in [1]). This leads to the dispersion relation as (see Appendix A for more details)

After substituting the characteristic impedance by its corresponding value for TE^z and TM^z waves, Eq. (7) leads to

 

Fig. 4 Iso-frequency wavevector dispersion (k_z = β_z − jα_z) versus k_t computed by both Bloch theory (dashed-dotted lines) and EMA (solid lines), for different chemical potential levels at 2 THz (a,b), and at 12 THz (c,d).

Download Full Size | PDF

It is clear from the analytical analysis above and from Fig. 4 that EMA well describes the wavevector dispersion when |κ_dd| ≪ 1 and |k_zd| ≪ 1, which can be verified when the transverse wavenumber k_t is not too large. For larger and larger values of k_t, the two assumptions would not be valid anymore. In Fig. 4, we also report the evolution of the dispersion curves by varying the chemical potential μ_c at 2 THz and 12 THz. For example at 2 THz, when μ_c = 0.4 eV, the spatial spectrum is bounded by k_t ≈ 36k₀ after which α_z increases dramatically. However, by increasing the chemical potential, ε′_t assumes larger negative values and β_z − k_t dispersion in Fig. 4 evolves to a flatter curve, thus the Brillouin zone edge is reached at larger values of k_t. Note that at larger spatial spectrum, the attenuation constant α_z increases due to finite losses, as shown in Fig. 4(b). By tuning the chemical potential, the dispersion characteristics can be controlled, for example, at 12 THz the dispersion for unbiased graphene is elliptic as well as when μ_c = 0.1 eV. This behavior appears since ε′_t exhibits zero crossing and becomes positive at 6.6 THz and 11 THz, when μ_c = 0 and μ_c = 0.1 eV, respectively. However, when the chemical potential is increased to 0.4 eV, hyperbolic dispersion arises at 12 THz, as shown in Figs. 4(c) and 4(d) where the inset of Fig. 4(c) shows the elliptic behavior for $k_t<\sqrt{{\epsilon }_d}{k}_0$. At high frequencies where the conductivity saturates to its universal value, the TM^z plasmonic modes are extremely confined to graphene layers (σ″ becomes very small) and higher frequencies, once σ″ > 0, graphene layers are incapable of supporting those modes [21]. Hence, the spectrum $k_t>\sqrt{{\epsilon }_d}{k}_0$ becomes mainly evanescent at frequencies with ε′_t > 0. In other words, after ε′_t exhibits a zero-crossing and becomes positive, the wavevector dispersion becomes elliptic.

3. Plane wave reflection and transmission by a thin film made of graphene-dielectric multilayers

A finite thickness graphene-dielectric multilayer film is considered comprising N graphene sheets stacked with silica SiO₂ dielectric spacer, such that a graphene sheet is at the topmost layer. The thickness of each SiO₂ spacer is 0.1 μm, and the total multilayered film thickness is D = Nd. For simplicity, all graphene layers are biased equally using a constant electrostatic potential [34]. For practical consideration, suppose that a thin film of silica (in the order of 100 nm) is deposited on an epitaxially-grown graphene monolayer repeatedly until creating an N layer stack; though larger thicknesses could be considered, it is rather simple to achieve the biasing range (μ_c up to 0.5 eV) using relatively lower electrostatic potential for smaller thicknesses [34, 53].

We investigate reflection and transmission under normal plane wave incidence, and at 30° oblique incidence for both TE^z and TM^z polarizations (here k_t = k₀ sin(θⁱ), where θⁱ is the incidence angle). Reflection and transmission coefficients are reported using the transfer matrix method (solid lines), and using EMA (circles) as well, for the two cases with N = 10 (D = 1 μm) and N = 20 (D = 2 μm), when μ_c= 0 eV and 0.4 eV in Fig. 5 and Fig. 6, respectively. As discussed in Section 2, at lower frequencies, TE^z waves are evanescent for any k_t when ε′_t < 0, while for TM^z waves the iso-frequency wavenumber dispersion is hyperbolic (when ε′_t < 0), consequently, a plane wave impinging on the structure with $k_t<\sqrt{{\epsilon }_d}{k}_0$ is very weakly transferred, specifically by evanescent coupling. This property is demonstrated in Fig. 5 at frequencies lower than 6 THz, and in Fig. 6 at frequencies lower than 24.6 THz. However, after ε′_t exhibits a zero-crossing and becomes positive, the plane wave is able to propagate, and the hyperbolic wavevector dispersion turns into elliptic (thus waves with $k_t<\sqrt{{\epsilon }_d}{k}_0$ can propagate). The transmission peak for normal incidence occurs when the effective ε′_t is near unity (matched to the free space, where the losses are negligible) for μ_c= 0 eV at ≈ 8 THz (reported in Fig. 5) and for μ_c= 0.4 eV at ≈ 33 THz (reported in Fig. 6); this is in accordance with the effective ε′_t plotted in Fig. 2(a). It is clear that changing the chemical potential of the graphene layers offers great tunability and possibility to control the transmission peak and spectrum. EMA is a good tool to describe plane wave interaction with a graphene-dielectric multilayer thin film, for small dielectric thickness d. In order to explore the validity of EMA for thicker dielectric spacers, we report in Fig. 7 the reflection and transmission coefficients for 10 layers of graphene-dielectric layers with varying thickness d, at 10 THz, assuming μ_c = 0.4 eV. It is shown that EMA yields a noticeable deviation from transfer matrix calculations when d > 0.2λ₀. Note that in Fig. 7 the transition from hyperbolic to elliptic dispersion occurs at d = 0.02λ₀, which implies ε′_t ≈ 0.

 

Fig. 5 Reflection and transmission versus frequency for a finite thickness graphene-silica multilayered HM, at normal, and oblique incidence for both TE^z and TM^z polarizations, calculated by transfer matrix method (solid lines) and EMA (circles) when graphene layers are unbiased, i.e., μ_c = 0 eV.

Download Full Size | PDF

 

Fig. 6 Reflection and transmission versus frequency, for the same set of parameters as Fig. 5, except that now graphene layers are biased with μ_c = 0.4 eV.

Download Full Size | PDF

 

Fig. 7 Reflection and transmission of a 30° TM^z wave from a 10 layer graphene-dielectric stack at 10 THz with variable spacer d, based on transfer matrix (solid lines) and EMA (circles).

Download Full Size | PDF

In particular, these results show two remarkable facts: (i) EMA agrees well with transfer matrix calculations for a wide range of frequencies and dielectric thicknesses, (ii) transmission and reflection by the graphene-based multilayered structure can be effectively tuned by electrostatic biasing. It is evident that graphene layers despite controlling transmission with such small thicknesses, at the same time can be designed to be almost transparent to plane wave excitation [32].

In this Section we have investigated reflection and transmission, and how this is predicted by EMA, for an incident plane wave, however a source or scatterer near the HM interface is able to generate a very wide spatial spectrum of plane waves, including the spectrum with k_t > k₀, which would be evanescent in free space. In the next Section we show how this wide spectrum is able to propagate inside the HM, similarly to what was done in [3] for a HM at optical frequencies made of dielectric and metallic layers.

4. Enhancement of emitted power by an impressed dipole at the surface of a graphene-based HM film

We investigate the power emitted by a transverse dipole located at a distance h above the graphene-silica multilayered HM as depicted in Fig. 8(a), over a a silicon substrate (sufficiently thick to be assumed infinitely long, with relative permittivity ε_Si=11.7). We assume here a unit cell of the HM consisting of a 0.1-μm thick silica layer stacked with a sheet of graphene sheet on top. We calculate the power emitted by the transverse dipole located at z = 0 as in Fig. 8(a), by using the spatial spectral formalism of TE^z and TM^z waves as outlined in [1]. The total power P_tot = P_up + P_down emitted by the transverse dipole illustrated in Fig. 8 is decomposed into the power terms directed toward the +z and −z directions (P_up and P_down, respectively) that are found by the spectral integrals

 

Fig. 8 (a) Dipole near-field emission over a finite thickness multilayer graphene HM over of a substrate, and (b) the its transverse equivalent network (TEN) for every spectral wave (both TE^z and TM^z).

Download Full Size | PDF

 

Fig. 9 (a) Ratio between power emitted in the lower space with the one in the upper space, P_down/P_up, and (b) the ratio P_tot/P_{free space} related to the transverse dipole located near the interface of free space and graphene-based HM made by N graphene layers on top of Si substrate. Calculations done via multilayer transfer matrix method (lines) and via EMA (markers) when chemical potential is μ_c = 0 eV.

Download Full Size | PDF

 

Fig. 10 (a) Ratio P_down/P_up and (b) ratio P_tot/P_{free space} for the same set of parameters in Fig. 9, but when chemical potential is μ_c = 0.4 eV.

Download Full Size | PDF

The ratio P_tot/P_{free space} also represents the increase of the local density of states (LDOS) with respect to LDOS at a point in free space [5, 50, 54], and this is also referred to as Purcell effect [6, 24]. We model the HM with thickness Nd via both the more accurate transfer matrix method (denoted by lines in the figure) and EMA (denoted by circles), and provide the results in Fig. 9(a) and (b) for unbiased graphene (μ_c= 0 eV), and in Fig. 10 for biased graphene (μ_c= 0.4 eV). In Fig. 9 one can observe that, at the lowest frequency 0.1 THz, there is a clear trend showing that when the number of layers (N) increases the ratio P_down/P_up also increases from ≈ 10⁶ up to ≈ 2.7 × 10⁸, when using calculations based on the transfer matrix method. Moreover in the lower frequency range, EMA overestimates P_down/P_up by almost one order of magnitude, as also discussed in [3] for a different HM configuration; however, as the frequency increases EMA and the transfer matrix method agree well. In Fig. 9(b), we observe the same disagreement of the transfer matrix calculations and EMA at lower frequencies, and it is clearly seen that the normalized emitted power (P_tot/P_{free space}) is ≈ 4 × 10³ at the lowest frequency 0.1 THz and drops linearly as the frequency increases. P_down/P_up exhibits a very sharp drop for N = 1 after around 1 THz, for N = 3 after about 2 THz, for N = 10 after 4 THz, and for N → ∞ after 6 THz (note that the hyperbolic to elliptic dispersion curve transformation occurs at 6.6 THz when μ_c = 0 eV obtained via EMA, see Fig. 2, in very good agreement with the N → ∞ case). When the chemical potential is increased to μ_c = 0.4 eV, the “transverse” permittivity ε′_t decreases to further negative values and the frequency of hyperbolic to elliptic dispersion curve transformation shifts from 6.6 THz (μ_c = 0 eV) to 24.6 THz (μ_c = 0.4 eV). In Fig. 10(a), at lower frequencies one can observe that P_down/P_up is increased by one order of magnitude for all cases whereas P_tot/P_{free space} decreases by almost one order of magnitude when compared to the case with μ_c = 0 eV. Moreover the frequency where P_tot/P_{free space} exhibits a sharp decrease shifts to a higher frequency when the chemical potential is increased to 0.4 eV in agreement with the change in the frequency of hyperbolic to elliptic dispersion curve transformation, as illustrated from Fig. 4(c). As seen in Figs. 10(a) and 10(b) both the enhancement of emitted power and the ratio of power directed to the –z direction are much larger than the Si-substrate case for a wide frequency band (1–6 THz) in the presence of graphene-dielectric. The interesting features in Figs. 9–10 are related to the power spectrum in Eq. (15), which is described in the following.

We report the emitted power spectrum for TM^z waves (solid lines), $p^{\text{TM}}\left(k_t\right)=p_{\text{up}}^{\text{TM}}\left(k_t\right)+p_{\text{down}}^{\text{TM}}\left(k_t\right)$, versus normalized transverse wavenumber k_td/π at 0.1 in Fig. 11(a) and 3 THz in Fig. 11(b), varying the number of graphene-dielectric layers assuming μ_c = 0 eV, and for a better visualization we provide Figs. 11(a) and 11(b) in both logarithmic and linear scales for the horizontal axis, in the left and right panels, respectively. For comparison we also show the power spectrum $p^{\text{TE}}\left(k_t\right)=p_{\text{up}}^{\text{TE}}\left(k_t\right)+p_{\text{down}}^{\text{TE}}\left(k_t\right)$ for N =1 and N → ∞ (dashed lines). At these two frequencies the composite multilayer exhibits hyperbolic dispersion for TM^z waves and propagation inside the HM occurs for $k_t>\sqrt{{\epsilon }_d}{k}_0$. We observe in Fig. 11(b) that in the high k_t spectrum, there are a larger number of spectral peaks when N increases; eventually yielding a continuous distribution of large spectral intensities when N → ∞. This explain the advantage of having a large number of layers. Moreover when N → ∞, one can notice that the power spectrum starts to rise strongly after $k_t=\sqrt{{\epsilon }_d}{k}_0$, in agreement with the propagating spectrum’s lower limit in the hyperbolic dispersion diagrams in Fig. 4. We would like to emphasize that plots in linear scale in Fig. 11 clearly show that the propagating power spectrum in the large k_t region is very wide and therefore it strongly contributes to the spectral integral in Eq. (14). Note that the upper limit of the power spectrum, cut-off at k_t,max, can be determined by evaluating the evanescent decay in free space between the dipole at z = 0 and the surface of the composite material at z = −h, given by $\text{exp}\left(-\sqrt{k_t^2-k_0^2}h\right)$. For example, by setting $\text{exp}\left(-\sqrt{k_{t,\text{max}}^2-k_0^2}h\right)=\xi$, where ξ is a predetermined small number, we can consider the power spectrum negligible when k_t > k_t,max. It is important to note that for h ≪ λ₀ this upper wavenumber limit k_t,max is independent of the operating frequency when k_t,max ≫ k₀, because in this case $\text{exp}\left(-\sqrt{k_{t,\text{max}}^2-k_0^2}h\right)\approx \text{exp}\left(-k_{t,\text{max}}h\right)$. These considerations explain why all spectral curves decay with very similar profile for very large k_t and therefore it is imposed mainly by the dipole distance h, for both frequencies examined in Fig. 11(a) and (b), i.e., at 0.1 and 3 THz.

 

Fig. 11 Emitted power spectrum $p^{\text{TM}}\left(k_t\right)=p_{\text{up}}^{\text{TM}}\left(k_t\right)+p_{\text{down}}^{\text{TM}}\left(k_t\right)$ in Eq. (15), solid lines, versus normalized traverse wavenumber k_td/π at (a) 0.1 THz, and (b) 3 THz, for different number of graphene-dielectric layers:N = 1, 10, and N → ∞. For comparison we also show the power spectrum $p^{\text{TE}}\left(k_t\right)=p_{\text{up}}^{\text{TE}}\left(k_t\right)+p_{\text{down}}^{\text{TE}}\left(k_t\right)$ for N =1 and N → ∞ (dashed lines). The points A, B, and C denote the spectrum points k₀d/π, $\sqrt{{\epsilon }_d}{k}_{0}d/\pi$, and $\sqrt{{\epsilon }_{Si}}{k}_{0}d/\pi$, respectively. Left panel plots have a horizontal logarithmic scale whereas right panel plots have a horizontal linear scale.

Download Full Size | PDF

At low frequency, in Fig. 11(a) all curves with different N, tend to exhibit the same behavior at the large k_t, in particular when (k_td/π) > 10⁻². The reason of this low frequency property, that does not occur at higher frequency in Fig. 11(b), is explained as follows. Propagation inside such multilayer stack consists in strong evanescent coupling between adjacent graphene sheets [23, 32], that is approximately proportional to $\text{exp}\left(-\sqrt{k_t^2-{\epsilon }_d{k}_0^2}d\right)$. Considering now the large k_t region where the power spectrum p^TM(k_t) is intense in Fig. 11(a), the exponential interlayer decay becomes stronger, at fixed k_t, when the frequency decreases (i.e., when ${\epsilon }_d{k}_0^2$ decreases). Thus, at low frequencies, less power is coupled to lower graphene layers, whereas most of the power is coupled to losses in the first graphene sheet closest to the dipole, implying that the number of layers becomes less effective on the power spectrum for large k_t and hence on the total emitted power integral in Eq. (14). For example, the case with 0.1 THz in Fig. 11(a), the total TM^z power emitted by the dipole is dominated by the wide power spectrum region (5 × 10⁻³π/d) < k_t < k_t,max, which is weakly dependent on the number of layers. Under this low frequency condition, we observe that the total TM^z emitted power becomes proportional to ω², independently on the number of layers, in agreement with the findings in [21] for a single graphene layer. A similar trend occurs for the power emitted as TE^z waves at these lower frequencies, though it is several orders of magnitude weaker than TM^z cases for large k_t (see the dashed curves in Fig. 11). Note that, instead, the free space emitted power by a dipole |p_t|²ω⁴η₀/ (12πc²) is proportional to ω⁴. Comparing the low frequency trends of the power emitted in free space with the one in presence of the HM one can explain the strong increase of the power ratio (P_tot/P_{free space}) ∝ ω⁻² in Fig. 10(b) as frequency decreases.

It is also important to provide a physical insight into the effect of the distance h on the power emitted by the impressed transverse dipole. In Fig. 12, we provide the plots of P_tot/P_{free space} and P_down/P_up at 2 THz versus the dipole’s distance h for the semi-infinitely periodic (N → ∞) graphene-dielectric multilayered structure with μ_c= 0 eV and μ_c= 0.4 eV, obtained both via EMA (markers) and the transfer matrix method (lines). In Figs. 12(a) and 12(b) we report the power ratios P_down/P_up and P_tot/P_{free space} for both HMg and HMd configurations, denoting a HM with graphene (HMg) and dielectric (HMd) as topmost layer, respectively [3,49]. We notice that the responses of both HMg and HMd configurations are very similar while the HMg has a slightly larger P_down/P_up and P_tot/P_{free space} for smaller h, in agreement with the observations in [49]. For the smallest reported distance h = 0.2 μm, the emitted power (P_tot/P_{free space}) and P_down/P_up are maximum. However for small h, EMA overestimates the reported parameters by one to two orders of magnitude for small h, as demonstrated in [49], whereas for h > 1 μm both methods agree well at the given frequency. Using the transfer matrix method we find the maximum ratio P_down/P_up ≈ 2 × 10⁶ when μ_c=0 eV, and it decreases to P_down/P_up ≈ 5 × 10⁵ as a result of increasing $Y_{\text{down}}^{\text{TM}}$ when ε′_t possesses more negative values. A similar change is also observed such that the maximum ratio P_tot/P_{free space} becomes ≈ 2 × 10⁴ when μ_c= 0 eV, and it decreases to P_tot/P_{free space} ≈ 4 × 10² when μ_c= 0.4 eV. The total emitted power decreases as the distance h increases, due to the stronger evanescent decay of high k_t spectrum. However, since the distance h is subwavelength, still a lot of power is able to couple into the HM. The power ratio P_down/P_up also exhibits a decrease with increasing h showing the key role of the coupling of the evanescent spectrum in free space to the propagating spectrum in the HM. We finally note from the observations in Fig. 12 that the accuracy of EMA is influenced by changing the values of μ_c, implying the effect of ε_t on EMA’s validity.

 

Fig. 12 (a) Ratio P_down/P_up and (b) ratio P_tot/P_{free space} related to the transverse dipole located near the interface of free space and a semi-infinite graphene-based HM at 2 THz plotted versus dipole distance h, for different chemical potential μ_c values obtained via transfer matrix method (lines) and EMA (markers).

Download Full Size | PDF

In summary, we have shown that the power emitted by a dipolar source in the proximity of a graphene-based HM is strongly enhanced, and that it can be effectively tuned by electrostatically biasing the graphene sheets, which makes this HM a promising candidate for tunable applications in the far-infrared frequencies.

5. Conclusion

We have investigated a novel design of HM for far-infrared frequencies based on graphene layers. The multilayer structure has been analyzed using EMA which, based on a permittivity homogenization model, predicts the HM features at far-infrared frequencies. We have quantitatively shown the capability of tuning the composite material properties via chemical potential of graphene. We have investigated plane wave interaction with a thin film made by few graphene sheets, and showed how the transmission frequency spectrum can be tuned. We have assessed the validity of EMA for both plane wave incidence and near-field radiation from a dipole, and we have shown that under certain conditions EMA is in good agreement with the transfer matrix analysis. In the last part of the paper, we have shown that a very wide spatial spectrum emitted by an electric dipole is allowed to couple into the graphene-based HM, that would be otherwise evanescent in free space. This generates two interesting main features: (1) the power emitted by the dipole is strongly enhanced (up to several orders of magnitude) by the presence of the graphene-based HM, and (2) most of the power is directed into the HM, also for relative subwavelength HM thicknesses realized with only a few graphene sheets.These properties seem to enable the use of this tunable graphene-based HM to efficiently absorb mm-waves and terahertz frequencies, and give rise to other possible applications including super resolution lenses.