Graphene-Dielectric Composite Metamaterials: Evolution from Elliptic to Hyperbolic Wavevector Dispersion and The Transverse Epsilon-Near-Zero Condition
aa r X i v : . [ c ond - m a t . m e s - h a ll ] M a y Graphene-Dielectric Composite Metamaterials: Evolution from Elliptic toHyperbolic Wavevector Dispersion and The Transverse Epsilon-Near-ZeroCondition
Mohamed A. K. Othman, Caner Guclu, and Filippo Capolino ∗ Department of Electrical Engineering and Computer Science,University of California, Irvine, CA, , USA
We investigated a multilayer graphene-dielectric composite material, comprising graphene sheets separatedby subwavelength-thick dielectric spacer, and found it to exhibit hyperbolic isofrequency wavevector dispersionat far- and mid-infrared frequencies allowing propagation of waves that would be otherwise evanescent in adielectric. Electrostatic biasing was considered for tunable and controllable transition from hyperbolic to ellipticdispersion. We explored the validity and limitation of the effective medium approximation (EMA) for modelingwave propagation and cutoff of the propagating spatial spectrum due to the Brillouin zone edge. We found thatEMA is capable of predicting the transition of the isofrequency dispersion diagram under certain conditions.The graphene-based composite material allows propagation of backward waves under the hyperbolic dispersionregime and of forward waves under the elliptic regime. Transition from hyperbolic to elliptic dispersion regimesis governed by the transverse epsilon-near-zero (TENZ) condition, which implies a flatter and wider propagatingspectrum with higher attenuation, when compared to the hyperbolic regime. We also investigate the tunabletransparency of the multilayer at that condition in contrast to other materials exhibiting ENZ phenomena.
I. INTRODUCTION
Hyperbolic metamaterial (HM) refers to a subcategory of uniaxially anisotropic metamaterial, that canbe modeled by a diagonal permittivity tensor (in Cartesian coordinates) comprising entries with bothpositive and negative real parts. The realization of hyperbolic dispersion allows wave propagation overa wide spatial spectrum (infinite for an ideal HM), that would be evanescent in a common isotropicdielectric [1]. HMs are realized at optical frequencies using metal-dielectric multilayers [2–4], or metallicnanowires [5], and at terahertz and infrared frequencies using semiconductor-dielectric multilayers [6, 7]or carbon nanotubes [8]. In multilayer HMs, the emergence of hyperbolic dispersion does not rely on anyresonant feature, thus it poses a potential for broadband enhancement of the local density of states (LDOS)[9], subwavelength imaging [10, 11], and lensing [12]. Spontaneous emission rate of an emitter, as wellas the radiative decay of dye molecules, is proportional to the LDOS [11], hence it can be substantiallyenhanced in the proximity of a hyperbolic metamaterial [13, 14]. It was demonstrated in [2] that thepower scattered by a passive nanosphere located in the proximity of a metal-dielectric HM is enhancedby orders of magnitude, while the HM absorbs most of the scattered power, opening a new frontier insuper absorbers designs based on near-fields transformation from evanescent to propagating regimes. Awide band absorption was devised in [15] using tilted carbon nanotubes.Multilayer HMs at optical frequencies take advantage of the wide frequency band in which metalsexhibit negative permittivity and support plasmonic modes [2, 3]. At infrared frequencies, graphene asa tunable inductive layer constitutes a potential building block for multilayer HM realizations. Further-more tunability of HMs can be achieved using static fields to bias graphene [16, 17]. It is a remarkablematerial with a wide operational frequency band starting from microwave regime [18], through terahertzfrequencies [19], and optical frequencies [20]. Graphene was utilized in design of metasurfaces in manydifferent applications, such as polarizers and absorbers [21, 22], and cloaking devices [23]. ∗ Electronic address: [email protected]
In this paper we investigate a graphene-dielectric multilayer material that shows promising propertiesas tunable HM at far- and mid-infrared frequencies, that was predicted to provide a large enhancementin the Purcell factor [24, 25]. In that recent work, the enhancement of emitted power by electrically-small emitter near the interface of graphene-based HM as well as the near-field absorption propertieswere developed using effective medium approximation (EMA) and transfer matrix methods, where thelimitations and validity of EMA were established [25]. Here we show how the wavevector dispersiondiagram can be controlled and even transformed between hyperbolic and elliptic curves at mid- and far-infrared regime. Moreover, we demonstrate the design guidelines of the graphene-based HM in terms ofthe physical parameters for the purpose of engineering the evolution from hyperbolic to elliptic dispersioncondition. . In the last part of the paper we explore the transverse epsilon near zero (TENZ) condition, itsrelation to the dispersion diagram and the enhanced transparency of a thin film made of TENZ graphene-dielectric layers for TM waves with a wide range of incidence angle. The fabrication of the metamaterialcomprising as few as ten graphene-dielectric layers, which were shown to have characteristics that resem-bles those of a semi-infinite stack [25], could be realized utilizing commercially available, high qualitychemical-vapor-deposition-(CVD)-grown graphene monolayer on a transition metallic (Ni or Cu) foil[26, 27]. from which graphene can be transferred onto a SiO /Si substrate using an intermediate hostsuch as a thermoplastic polymethyl-methacrylate (PMMA) for enhancing the transfer process efficiency[28]. This process is followed by depositing a thin film of SiO or SiC on the graphene flake using CVD.However, it was shown that a graphene monolayer on SiO can become highly disordered and increasesscattering losses [29]. The transfer of few-layer graphene (FLG) [26] on other compatible materials suchas Boron-Nitride (h-BN) might be of interest toward realizing the metamaterial, since h-BN shares thesame hexagonal structure with graphene [30]. II. EFFECTIVE MEDIUM ANALYSIS OF GRAPHENE-DIELECTRIC MULTILAYERS
Graphene is a one-atom-thick layer of hexagonal arrangement of carbon atoms with a lattice constant of . nm, hence spatial dispersion effects introduced by graphene periodicity can be in general neglectedat terahertz frequencies. Although the existence of extremely slow surface modes can trigger spatialdispersion effects [18, 31], those modes are essentially highly evanescent due to the periodicity of themultilayer structure studied here, as it will be shown in Sec. 3. Graphene is electrically modeled bythe local isotropic sheet conductivity σ = σ ′ + jσ ′′ (assuming time-harmonic variation of e jωt ), whichaccounts for both interband and intraband contributions to the total electronic transport [32, 33]. The sheetconductivity σ is computed by the Kubo formula [34], which yields a function of frequency, chemicalpotential µ c , phenomenological scattering rate Γ , and temperature T . Here we assume for graphene Γ = 0 . meV (using the same notation as in [34]), which corresponds to a mean electron scattering timeof about ps, at room temperature T = 300 K. Graphene supports relatively low loss TM plasmonicmodes [16] (dictated by the negative imaginary part of the surface conductivity σ ′′ < ). As such, σ ′′ ,modeling the reactive response of graphene, plays a fundamental role in the manifestation of hyperbolicdispersion in multilayer graphene-dielectric materials, as described in the following. We aim at analyzingan infinite periodic multilayer structure depicted in Fig. 1 whose unit cell is composed of a graphene sheetand a dielectric layer of subwavelength thickness d and relative permittivity ǫ d . A physical understandingof wave propagation in such multilayers with subwavelength period can be established by using theeffective medium approximation (EMA) approach, which is a quasi-static or local approximation formetamaterials, often adopted for metal-dielectric multilayers [2, 3, 31]. According to EMA, the periodicmultilayer is regarded as an anisotropic homogeneous medium with effective relative permittivity tensor ǫǫǫ eff = ǫ t (ˆ x ˆ x + ˆ y ˆ y )+ ǫ z ˆ z ˆ z , where the relative effective transverse permittivity ǫ t is found by averaging thetransverse effective displacement current over the associated electric field in a unit cell (here, the effectivedisplacement current is defined as a quantity that includes both displacement current in the dielectric slaband conduction current in the infinitesimally-thin graphene sheet). Then the relative effective permittivityparameter for transversely polarized field is Figure 1: Graphene-dielectric multilayer HM topology, modeled by a periodically-loaded transmission line. The unitcell is indicated on the right and the graphene sheet is represented as a shunt admittance, and we denote the referenceplane for evaluating the Bloch impedance. At far- and mid-infrared frequencies, TM z waves exhibit hyperbolicisofrequency wavevector dispersion. ǫ t = ǫ ′ t − jǫ ′′ t = ǫ d − j σωǫ d . (1)Since an individual graphene sheet is infinitesimally-thin, the conduction current is always along thesheet, hence the permittivity experienced by z − directed electric field is not affected by graphene, leadingto ǫ z = ǫ d . The relation in Eq. (1) implies that when the graphene sheet is adequately inductive, inparticular when σ ′′ < − ωǫ ǫ d d , we obtain ǫ ′ t < and in turn the isofrequency wavevector dispersion ishyperbolic [2], as demonstrated next. Let us consider plane waves propagating inside the metamaterialwith the spatial dependence e − j k . r where k = k x ˆ x + k y ˆ y + k z ˆ z is the wavevector. A plane wave analysisis particularly useful in understanding the multilayer’s response to sources because the radiation of adipole inside or close to the metamaterial can be represented as a spatial spectral sum of plane waves.Due to the symmetry of the multilayer metamaterial with respect to the z axis, we will use k t = q k x + k y for denoting the transverse wavenumber component and in the following k t is taken real representing thespatial spectrum of TE z (electric field transverse to z ) and TM z (magnetic field transverse to z ) waves.The z -directed wavenumber k z = β z − jα z can assume complex values modeling propagation andattenuation, accounting also for natural losses in the material constituents. Accordingly, the wavevectordispersion of TE z and TM z waves inside the effective medium is given as k z = ǫ t k − k t , TE z (2) k z = ǫ t k − ǫ t ǫ d k t , TM z (3)where k = ω √ µ ǫ is the wavenumber in free space. When the losses are neglected (i.e., if σ ′ → ) onewould obtain purely real ǫ t , hence k z (obtained via Eq.(2) and Eq. (3)) assumes either purely real values,denoting the propagating spectrum, or purely imaginary values, denoting the evanescent spectrum. Inthis lossless case, hyperbolic dispersion occurs when ǫ t < , and the HM uniaxial medium allows forpropagation (i.e., k z is purely real) of extraordinary waves (TM z ) with a large transverse wavenumber k t > √ ǫ d k ; these waves with k t > √ ǫ d k would be otherwise evanescent (i.e., k z is purely imaginary)either in a isotropic dielectric with permittivity ǫ d , or in a generic uniaxial anisotropic media with ǫ t > . This unusual phenomenon implies that high k t spectrum emanating from sources, which would beevanescent in free space, can be converted to propagating waves at HM interfaces. Ordinary waves (TE z )are, however, evanescent for any k t when ǫ t < . On the other hand, when ǫ t > we have real k z onlyfor limited spectrum of TM z waves with k t < √ ǫ d k , which leads to the elliptic isofrequency wavevectordispersion. Therefore the transition between hyperbolic to elliptic regimes is associated to the condition ǫ t = 0 .Instead, for realistic lossy cases, k z is complex and the wavevector isofrequency dispersion becomeselliptic-like and hyperbolic-like (for ǫ ′ t > and ǫ ′ t < , respectively), as shown in the examples in nextsection. However, the interpretations regarding propagation of power are still valid provided that lossesare relatively small, and we will show that moderate propagation losses is a major advantage of graphene-based HMs at far- and mid-infrared frequencies. When applying EMA, the dispersion relation β z − k t ishyperbolic-like for k t > √ ǫ d k when ǫ ′ t < , and it converges to the asymptote | β z | ≈ | ǫ ′ t k t /ǫ d | for largespatial wavenumber k t , i.e., the β z − k t dispersion becomes linear, with a slope of | σ ′′ / ( ωǫ ǫ d d ) | .To validate our EMA hypothesis, we obtain a more accurate representation of the wavevector dispersionrelation by employing Bloch theory [35] for a periodically loaded transmission line whose unit cell isillustrated in Fig. 1. When each graphene sheet is modeled with a complex admittance Y s = σ = σ ′ + jσ ′′ , the dispersion relation for TM z or TE z waves in the periodic structure is cast in the form cos k z d = cos κ d d + j Y s Z d sin κ d d, (4)where κ d = p ǫ d k − k t is the z -directed wavenumber of a wave inside the dielectric spacer, Z TM d = κ d / ( ωǫ ǫ d ) and Z TE d = ωµ /κ d are the characteristic wave impedances for TM z and TE z waves, respec-tively. This relation in (4) is yet accurate for arbitrary d and k t , i.e., accounts for transverse wavenumberdispersion. For the spectrum in which the dielectric layer’s thickness is much smaller than the Blochwavelength and the wavelength inside the dielectric itself ( | k z d | ≪ , | κ d d | ≪ ), we can apply thefollowing small argument approximations cos x ≈ − x and sin x ≈ x , the dispersion relation in Eq.(4) simplifies to the one obtained via EMA in Eq. (2) and Eq. (3) using the same definitions for ǫ t and ǫ z [25]. As we will discuss thoroughly in Sec. 3, Bloch theory proves that the propagating spectrum ofTM z waves is limited due to the periodicity, manifested by the Brillouin zone edge at which β z = ± π/d ,and therefore the propagating spectrum in realistic HMs has an upper bound even in lossless cases. Nev-ertheless, the Brillouin zone edge (i.e., β z = ± π/d ) is reached in general at higher values of k t , providedthat the period d is extremely subwavelength .In the following we report some aspects that demonstrate the merits of graphene-based HM: Grapheneconductivity σ = σ ′ + jσ ′′ is tunable with chemical potential variation via electrostatic biasing, hence ǫ ′ t is also tunable through negative or positive values, at a fixed frequency. This implies a possible transitionbetween hyperbolic to elliptic wavevector dispersion. The realization of HMs using graphene is alsoprone to graphene’s frequency response. For instance, graphene sheets are mainly capacitive in mid-and near-infrared frequencies, because intraband contributions in graphene are dominant, and the TM z surface modes on a single graphene sheet become on the improper Riemann sheet [16]. On the otherhand, at very low frequencies (GHz regime), the interband conductivity dominates leading to high losses.Hence a proper frequency range for realizing hyperbolic dispersion extends from far-infrared up to lowmid-infrared frequencies. Furthermore, the dielectric thickness also plays role on the frequency rangeof HM design. As the dielectric thickness is increased, the frequency range of negative ǫ ′ t shifts tolower frequencies which are undesirable due to significant losses in graphene. Moreover, thicker spacersrequire a larger biasing electrostatic potential between layers to achieve a moderate chemical potentiallevel in graphene sheets. On the other hand, when considering smaller periods (in the range of severalnanometers), it is expected that graphene sheets are no longer electronically isolated for such quantum-scale interspacing, and a tight binding model for graphene layers must be taken into account in orderto evaluate the conductivity of graphene sheets [36, 37]. Therefore, for very small thicknesses, bothEMA relation, reported in (3), and transfer matrix analysis must be modified to account for quantumtunneling between graphene sheets. In the next section we will explore and provide illustrative examplesfor graphene-based HM designs in terms of frequency response, losses and tunability and we will assessthe validity of the EMA in predicting hyperbolic or elliptic dispersion regimes. III. HYPERBOLIC AND ELLIPTIC WAVEVECTOR DISPERSION
Let us consider a multilayer stack depicted in Fig. 1, that comprises graphene sheets and dielectriclayers with ǫ d = 2 . and thickness d. In our illustrations we only adopt positive values for graphenechemical potential owing to the assumed reciprocity in the multilayers, and consider a typical range for µ c up to 0.5 eV in individual graphene sheets as suggested in [23]. We plot in Fig. 2 the relative transverse permittivity ǫ t = ǫ ′ t − jǫ ′′ t versus frequency, for various chemical potential levels ( µ c = 0 , 0.25, and 0.5eV) and dielectric thickness ( d = ǫ ′ t ,where σ ′′ = − ωǫ ǫ d d , is primarily defined by the period d and it can be tuned via the chemical potential;in turn the frequency of transition between the hyperbolic and the elliptic dispersion regimes can becontrolled. Assuming d = 100 nm (solid lines) in Fig. 2(a) we show that the frequency at which ǫ ′ t = 0 shifts from . THz to . THz by increasing the chemical potential from 0 eV to 0.5 eV. For d = ǫ ′ t = 0 is observed by varying µ c . Moreover when µ c = 0 ,we see that ǫ ′ t = 0 occurs at 8.7 THz for d =
50 nm, a higher frequency than the d = 100 nm case whosezero-crossing frequency is around 6.6 THz. Graphene sheets become capacitive at higher frequencies( σ ′′ = 0 denotes the transition from inductive to capacitive, for instance, σ ′′ = 0 at ≃ THz when µ c = 0 eV), however its contribution to ǫ ′ t becomes negligible because of both ω in the denominatorof (1) and graphene conductivity saturates to πe / (2 h ) ≈ µ S with a very small imaginary part, andhence ǫ ′ t approaches ǫ d .We show a relative variation in ǫ ′′ t when µ c is increased, indicating a possible way to tune losses. Notethat when the frequency dependent transverse permittivity ǫ ′ t turns positive and becomes close to unity,satisfying σ ′′ ≈ ωǫ d (1 − ǫ d ) , for instance at 15.6 THz when µ c = d =
100 nm, a finitegraphene-dielectric multilayer becomes almost transparent to TE z and TM z plane waves in free spacewith k t ≪ k , and all waves would travel with k z ≈ k , as seen from (3) when ǫ ′ t ≈ .In order to address some design considerations and tuning opportunities of graphene-based HM, weshow in Fig. 3(a,b), the real and imaginary parts of ǫ t as a colormap versus µ c and d . We also indicate the ǫ ′ t = 0 contour denoting the transition between hyperbolic and elliptic dispersion regimes. The selectionof d determines the range of chemical potential levels in which hyperbolic/elliptic dispersion occurs. Forinstance, when d = 0 . µ m, a tuning range for hyperbolic dispersion starts at µ c = 0 . eV, while for d = 0 . µ m it begins at µ c = 0 . eV; this illustrates the need for thinner dielectric spacers due to (cid:1) (cid:1) (cid:2)(cid:3)(cid:4)(cid:3)(cid:5)(cid:6) (cid:2)(cid:3)(cid:4)(cid:7)(cid:8)(cid:9)(cid:3)(cid:5)(cid:6) (cid:2)(cid:3)(cid:4)(cid:7)(cid:9)(cid:3)(cid:5)(cid:6) (cid:1) (cid:1) (cid:1) (cid:1) (cid:2) (cid:3) (cid:21)(cid:4) (cid:4) (cid:21)(cid:4) (cid:21) (cid:26)(cid:8)(cid:4)(cid:26)(cid:21)(cid:9)(cid:26)(cid:21)(cid:4)(cid:26) (cid:9)(cid:4)(cid:9) (cid:10)(cid:11)(cid:5)(cid:12)(cid:13)(cid:5)(cid:14)(cid:1)(cid:15)(cid:3)(cid:16)(cid:17)(cid:18)(cid:19)(cid:20) (cid:21)(cid:4) (cid:4) (cid:21)(cid:4) (cid:21) (cid:8)(cid:21)(cid:7)(cid:9)(cid:21)(cid:4)(cid:7)(cid:9)(cid:4) (cid:10)(cid:11)(cid:5)(cid:12)(cid:13)(cid:5)(cid:14)(cid:1)(cid:15)(cid:3)(cid:16)(cid:17)(cid:18)(cid:19)(cid:20) (cid:2) (cid:4) (cid:5) (cid:2) (cid:4) (cid:5)(cid:5) (cid:22)(cid:23)(cid:24) (cid:22)(cid:25)(cid:24) Figure 2: Real and imaginary parts of the effective relative transverse permittivity ǫ t = ǫ ′ t − jǫ ′′ t for graphene-basedmultilayer HM for two possible designs with d = 100 nm (solid lines) and d = 50 nm (dashed lines). (cid:1) (cid:1) (cid:2) (cid:2)(cid:3) (cid:2) (cid:3) (cid:1) (cid:4) (cid:5) (cid:3) (cid:6)(cid:7)(cid:8) (cid:6)(cid:7)(cid:9) (cid:6)(cid:7)(cid:10) (cid:6)(cid:7)(cid:11) (cid:12)(cid:6)(cid:6)(cid:7)(cid:12)(cid:6)(cid:7)(cid:8)(cid:6)(cid:7)(cid:13)(cid:6)(cid:7)(cid:9)(cid:6)(cid:7)(cid:14) (cid:1) (cid:1) (cid:2) (cid:2)(cid:3) (cid:2) (cid:3) (cid:1) (cid:4) (cid:5) (cid:3) (cid:6)(cid:7)(cid:8) (cid:6)(cid:7)(cid:9) (cid:6)(cid:7)(cid:10) (cid:6)(cid:7)(cid:11) (cid:12)(cid:6)(cid:6)(cid:7)(cid:12)(cid:6)(cid:7)(cid:8)(cid:6)(cid:7)(cid:13)(cid:6)(cid:7)(cid:9)(cid:6)(cid:7)(cid:14)(cid:15)(cid:8)(cid:6)(cid:15)(cid:12)(cid:14)(cid:15)(cid:12)(cid:6)(cid:15)(cid:14)(cid:6) (cid:6)(cid:6)(cid:7)(cid:8)(cid:6)(cid:7)(cid:9)(cid:6)(cid:7)(cid:10) (cid:4) (cid:5) (cid:6) (cid:16)(cid:17)(cid:6) (cid:4) (cid:5) (cid:6) (cid:4) (cid:5) (cid:6) (cid:6) (cid:6)(cid:7)(cid:18)(cid:14) (cid:8)(cid:7)(cid:8) (cid:4) (cid:5) (cid:6) (cid:16)(cid:17)(cid:6)(cid:19)(cid:20)(cid:21) (cid:19)(cid:22)(cid:21) (cid:1) Figure 3: Contour plot exploring the tuning capabilities of ǫ t for graphene-based HM via chemical potential µ c anddielectric thickness d at 10 THz. the limitations on the chemical potential levels’ adjustability, up to 0.5 eV in this paper. On the otherhand, the choice of a thinner dielectric spacer, i.e., smaller d , effectively induces higher ǫ ′′ t , so the lossesembodied in ǫ ′′ t are larger at the same frequency and bias. For example, when d = 0 . µ m, ǫ ′′ t ≃ . but when d = 0 . µ m we notice that ǫ ′′ t ≃ . , with larger negative ǫ ′ t in the former case than in thelatter. Nonetheless, a thin dielectric spacer allows feasible biasing by standard values of static potential[21]. This demonstrates a basic trade-off in graphene-dielectric HM design, between the tuning ranges,losses, and effective negative values of ǫ ′ t , and leads to a broad interpretation of the respective wavevectordispersion, as described next.The TM z wavevector dispersion diagrams according to EMA Eq. (3) and Bloch theory for the multilay-ered medium Eq. (4) are shown in Fig. 4. Here we report one of the two solutions of Eq. (3) and Eq. (4)for k z = β z − jα z that corresponds to a wave whose Poynting vector is directed towards the + z direction,noting that the other root − k z is also a solution of (3) and (4), not reported for symmetry reasons. Accord-ingly, the attenuation constant α z has positive sign, associated to the field decay (due to possible losses)along the + z direction. On the other hand, for the hyperbolic regime one observes β z < indicating back-ward wave propagation because it satisfies the backward wave condition β z α z < explained in [38], for k t > √ ǫ d k . In general, for the elliptic case, when k t < √ ǫ d k the valid k z = β z − jα z solution withpositive α z is the one with β z > , indicating that waves under the elliptic dispersion regime are forwardwaves because they satisfy the condition β z α z > . In Fig. 4(c,d) we show the dispersion diagrams in amuch wider spatial spectrum than in Fig.4(a,b) for the same cases. In the reported cases, all with d = 100 nm, β z curves in Fig. 4 keep either an overall hyperbolic or elliptic shape due to limited losses. When µ c = ǫ ′ t > ) the medium exhibits elliptic dispersion, moreover β z is nonzerofor k t > √ ǫ d k where α z exhibits a dramatic increase, i.e., waves become mostly evanescent. On theother hand, when µ c = ǫ ′ t < leading to hyperbolic dispersion. We emphasizethat EMA is fully capable of predicting the hyperbolic and elliptic wavevector dispersion regimes in thespatial spectrum reported in Fig. 4(a,b) in perfect agreement with the Bloch wavenumber. In a muchwider range of the spatial spectrum k t as in Fig. 4(c,d) the EMA-based normalized wavenumber β z /k starts to deviate from Bloch theory. Bloch theory predicts the band edge where β z approaches − π/d and α z exhibits a dramatic increase, denoting a bandgap. However EMA assumes infinite growth of β z /k following the asymptotic linearized β z − k t relation, given by β z ≈ − ǫ ′ t k t /ǫ d when k t ≫ k . For highernegative values of ǫ ′ t , (corresponding to higher µ c ), the Brillouin zone band edge is met at smaller k t dueto steeper β z − k t curves, as seen from Fig. 4(c) ( ǫ ′ t ≃ − and ǫ ′ t ≃ − for µ c = 0 . and 0.5 eV). Al-though the effective permittivity parameters are important for fast characterization of graphene-dielectriccomposites and providing physical interpretation of the evolution from elliptic to hyperbolic dispersion,they do not account for transverse wavenumber dispersion [31, 39]. Accordingly, EMA predicts an in-definite propagating spatial spectrum in HMs (that is indeed limited by Brillouin zone edge according to Figure 4: Wavevector dispersion diagram of (a) β z and (b) α z versus k t (both normalized by k ) at 10 THz and d = 100 nm. In (c) and (d) a wider spatial spectrum of the wavevector dispersion is provided in order to identify k t values where β z approaches the Brillouin zone edge ( β z = − π/d ) denoted by a horizontal dotted line in (c). Thishappens when k t ≈ k and k t ≈ k for µ c = 0 . eV and µ c = 0 . eV, respectively. Calculations are basedon both EMA (dash-dotted lines) and Bloch theory (solid lines).Figure 5: Real and imaginary parts of the Bloch (solid lines) and effective (dashed lines) impedance of graphene-dielectric multilayers with d = 100 nm when µ c = 0 . eV for (a) k t = 0 and (b) k t = 5 k . Bloch model), and consequently overestimate the LDOS and the near-field power absorption in HMs asalready discussed in [2–4, 25].We provide in Fig. 5 both the Bloch impedance of graphene-dielectric multilayers at the referenceplane shown in Fig. 1, with d = 100 nm. In addition, we report the effective wave impedance of themetamaterial obtained via EMA for TM z plane wave, Z eff = k z / ( ωǫ ǫ t ) where k z is evaluated using (cid:20) (cid:21) (cid:2) (cid:22)(cid:23) (cid:2) (cid:24) (cid:25) (cid:26) (cid:17) (cid:27) (cid:28) (cid:29) (cid:19) (cid:15) (cid:15)(cid:30)(cid:13) (cid:15)(cid:30)(cid:31) (cid:15)(cid:30)(cid:32) (cid:15)(cid:30)(cid:33) (cid:15)(cid:30)(cid:34)(cid:13)(cid:15)(cid:31)(cid:15)(cid:32)(cid:15)(cid:33)(cid:15)(cid:34)(cid:15) (cid:34)(cid:15)(cid:15)(cid:13)(cid:15)(cid:15)(cid:15)(cid:13)(cid:34)(cid:15)(cid:15) (cid:5) (cid:6) (cid:17)(cid:2)(cid:18)(cid:19) (cid:20) (cid:21) (cid:2) (cid:22)(cid:23) (cid:2) (cid:24) (cid:25) (cid:26) (cid:17) (cid:27) (cid:28) (cid:29) (cid:19) (cid:15) (cid:15)(cid:30)(cid:13) (cid:15)(cid:30)(cid:31) (cid:15)(cid:30)(cid:32) (cid:15)(cid:30)(cid:33) (cid:15)(cid:30)(cid:34)(cid:13)(cid:15)(cid:31)(cid:15)(cid:32)(cid:15)(cid:33)(cid:15)(cid:34)(cid:15) (cid:34)(cid:15)(cid:15)(cid:13)(cid:15)(cid:15)(cid:15)(cid:13)(cid:34)(cid:15)(cid:15) (cid:1)(cid:2)(cid:3) (cid:1) (cid:4) (cid:5) (cid:6)(cid:7)(cid:3) (cid:1) (cid:4) (cid:5) (cid:8)(cid:9)(cid:10) (cid:8)(cid:11)(cid:10) (cid:12)(cid:13) (cid:13) (cid:2) (cid:3) (cid:4) (cid:14)(cid:12)(cid:15) (cid:2) (cid:3) (cid:4) (cid:16)(cid:12)(cid:15) (cid:2) (cid:3) (cid:4) (cid:16)(cid:12)(cid:15) (cid:2) (cid:3) (cid:4) (cid:14)(cid:12)(cid:15) (cid:5) (cid:6) (cid:17)(cid:2)(cid:18)(cid:19) Figure 6: (a) Real and (b) imaginary parts of the Bloch impedance for d = 100 nm and k t = 0 . Eq. (3), see [40]. The two impedances are close to each other for k t = 0 case (Fig. 5(a)) whereas for k t = 5 k the effective impedance shows noticeable difference for both real and imaginary parts fromthe Bloch calculations. Nonetheless, the effective impedance provides a good prediction regarding thetransition frequency between propagating and evanescent spectra. Moreover, we notice that the real partof the impedance is negligible at low frequencies in Fig. 5(a), whereas it peaks at the frequency where ǫ ′ t = 0 . From Fig. 6(a) one can see that after ǫ ′ t turns positive, the impedance becomes dominantly real,with relatively small reactive part, owing to the presence of a mainly propagating plane wave in ellipticdispersion regime for k t = 0 . On the contrary for k t = 5 k case, at lower frequencies , wave propagatesin the hyperbolic dispersion regime while having ǫ ′ t < , and the impedance real part is relatively large,as depicted in Fig. 5(b), whereas the impedance becomes almost purely reactive after ǫ ′ t turns positive,denoting a mainly evanescent wave. At higher frequencies, the impedance for k t = 0 case becomesmatched to free space at ≈
37 THz at which ǫ ′ t ≈ as shown in Fig. 6(a). At much higher frequencyranges, the impedance approaches the impedance in isotropic lossless dielectric where ǫ t ≈ ǫ d in both Fig.5(a) and (b). For clarification, we report the Bloch impedance as a color plot showing the dependanceon frequency and chemical potential in Fig. 6, where the impedance peaking is observed as a clearmanifestation of the TENZ condition, as it will be demonstrated also in Sec. 4. Based on the conclusionsin [25], in order to guarantee the validity of EMA for each spectral component of propagating planewaves with k t < k , the dielectric thickness should be electrically-small, i.e., d < . λ for accuraterepresentation of the impedance and wavevector using the homogenized model derived above.We report in Fig. 7 the frequency dependance of the quantity | β z /α z | where α z and β z are calculatedby Bloch theory, for graphene-dielectric multilayers with d = 100 nm. The ratio | β z /α z | constitutes afigure of merit for understanding if a wave is mainly propagating or attenuating. The horizontal whitedash-dotted line marks the transition frequency from hyperbolic to elliptic dispersion (the latter occurringalways above the transition frequency) and the transition happens when the real part ǫ ′ t crosses zero andturns positive causing the elliptic regime. For k t < √ ǫ d k , β z is relatively very small compared to α z , which implies mainly evanescent spectrum (purely evanescent in absence of losses), for hyperbolicdispersion frequencies ω < − σ ′′ / ( ǫ ǫ d d ) . However, for k t > √ ǫ d k , wavevector dispersion has ahyperbolic-like shape, with attenuation α z moderately low (and slightly increasing as seen in Fig. 7) dueto the losses in graphene, and therefore | β z /α z | exhibits an overall increase, where it reaches a maximumvalue ≃ as in µ c = 0 .5 eV yielding a wide propagating spectrum √ ǫ d < k t /k . at 10 − k t , the propagation constant β z tends to − π/d while α z experiencesan abrupt increase, as shown in Fig. 4(d), denoting the beginning of a strong evanescent spectrum. Inthe elliptic dispersion regime, occurring at higher frequencies such that ω > − σ ′′ / ( ǫ ǫ d d ) , the trend for β z and α z is reversed. Elliptic dispersion arises at 6.6 THz for µ c = 0 eV, as depicted in Fig. 7, andthe propagating spectrum with k t < √ ǫ d k is allowed in the composite multilayer. For higher chemicalpotentials, as for example µ c = 0 . eV, hyperbolic wavevector dispersion is supported for frequencies (cid:1) (cid:2) (cid:1)(cid:2)(cid:3)(cid:2)(cid:4)(cid:5)(cid:6)(cid:7)(cid:8) (cid:6)(cid:9)(cid:8) (cid:10) (cid:11) (cid:4) (cid:12)(cid:13) (cid:4) (cid:14) (cid:15) (cid:16) (cid:17) (cid:18) (cid:19) (cid:20) (cid:21) (cid:22)(cid:3)(cid:23)(cid:3)(cid:24)(cid:3)(cid:25)(cid:3)(cid:26)(cid:3) (cid:10) (cid:11) (cid:4) (cid:12)(cid:13) (cid:4) (cid:14) (cid:15) (cid:16) (cid:17) (cid:18) (cid:19) (cid:20) (cid:21) (cid:22)(cid:3)(cid:23)(cid:3)(cid:24)(cid:3)(cid:25)(cid:3)(cid:26)(cid:3)(cid:26)(cid:22)(cid:3)(cid:22)(cid:26)(cid:23)(cid:3) (cid:26)(cid:3)(cid:22)(cid:3)(cid:3)(cid:22)(cid:26)(cid:3) (cid:22)(cid:3) (cid:3) (cid:4) (cid:5)(cid:3) (cid:6) (cid:22)(cid:3) (cid:3) (cid:22) (cid:22)(cid:3) (cid:23) (cid:22)(cid:3) (cid:3) (cid:4) (cid:5)(cid:3) (cid:6) (cid:22)(cid:3) (cid:3) (cid:22) (cid:22)(cid:3) (cid:23) (cid:27)(cid:28)(cid:28)(cid:29)(cid:30)(cid:31)(cid:29)(cid:15) (cid:19)(cid:16)(cid:30)(cid:4)(cid:11)(cid:9)(cid:32)(cid:28)(cid:29)(cid:15) (cid:19)(cid:16)(cid:30)(cid:4)(cid:11)(cid:9)(cid:32)(cid:28)(cid:29)(cid:15) (cid:27)(cid:28)(cid:28)(cid:29)(cid:30)(cid:31)(cid:29)(cid:15) (cid:3) (cid:1) (cid:2) (cid:1)(cid:2)(cid:3)(cid:33)(cid:26)(cid:2)(cid:4)(cid:5) Figure 7: The figure of merit | β z /α z | versus frequency and spatial wavenumber k t , for both hyperbolic and ellipticregimes. Two chemical potential levels are considered: (a) µ c = µ c = 0 . eV. (cid:1) (cid:2) (cid:3)(cid:1) (cid:4) (cid:5) (cid:1) (cid:2) (cid:3) (cid:4) (cid:5)(cid:6)(cid:7)(cid:8)(cid:8)(cid:6)(cid:7)(cid:9) (cid:5) (cid:1) (cid:2) (cid:3) (cid:4) (cid:5)(cid:6)(cid:7)(cid:8)(cid:8)(cid:6)(cid:7)(cid:9)(cid:7)(cid:5)(cid:8)(cid:5)(cid:5)(cid:8)(cid:7)(cid:5) (cid:7)(cid:5)(cid:8)(cid:5)(cid:5)(cid:8)(cid:7)(cid:5) (cid:8)(cid:5) (cid:5) (cid:8)(cid:5) (cid:8) (cid:8)(cid:5) (cid:9) (cid:8)(cid:5) (cid:5) (cid:8)(cid:5) (cid:8) (cid:8)(cid:5) (cid:9) (cid:1) (cid:2) (cid:3)(cid:1) (cid:4) (cid:10)(cid:11)(cid:12)(cid:13)(cid:14)(cid:15)(cid:16)(cid:17)(cid:18)(cid:19) (cid:20)(cid:17)(cid:17)(cid:18)(cid:12)(cid:21)(cid:18)(cid:19)(cid:20)(cid:17)(cid:17)(cid:18)(cid:12)(cid:21)(cid:18)(cid:19) (cid:5) (cid:5)(cid:5)(cid:6)(cid:8) (cid:5)(cid:6)(cid:8) (cid:22)(cid:23)(cid:24) (cid:22)(cid:15)(cid:24) (cid:6) (cid:7) (cid:25)(cid:26)(cid:5)(cid:26)(cid:13)(cid:27) (cid:6) (cid:7) (cid:25)(cid:26)(cid:5)(cid:6)(cid:7)(cid:26)(cid:13)(cid:27) Figure 8: The figure of merit | β z /α z | versus dielectric thickness d and spatial wavenumber k t , at 10 THz, for bothhyperbolic and elliptic regimes. Two chemical potential levels are considered: (a) µ c = µ c = 0 . eV. up to . THz, and the dispersion becomes elliptic thereafter. Notice that at frequencies less than 1THz, waves poorly propagate due to higher losses in graphene sheets, i.e., wave propagation has a lowfigure of merit. On the other hand, elliptic dispersion regime, occurring for frequencies greater than 30THz, has small attenuation constant for k t < √ ǫ d k due to relatively low loss in graphene, and thus ahigh figure of merit | β z /α z | > . Note that the lowest operational frequency for hyperbolic dispersionregime with high | β z /α z | is limited by graphene losses, whereas the highest frequency is tunable by thechemical potential.We now examine the how the figure of merit | β z /α z | varies versus the transverse wavenumber k t ,assuming different design values for the dielectric spacing d . In Fig. 8(a) we observe | β z /α z | at 10 THzvarying d , for µ c = 0 eV, where only elliptic dispersion regime is observed for any thickness d considered.However, hyperbolic dispersion is supported when appropriate chemical potential is achieved, as shownin Fig. 8(b) for µ c = 0 . eV. In this latter case, when d = 1 µ m, TM z waves are mainly evanescentfor large transverse wavenumber k t > √ ǫ d k , irrespective of the chemical potential levels reported here.Consequently, a typical dielectric thickness in the range of 50 −
100 nm is deemed appropriate to utilizein graphene-dielectric multilayers for tunable HM designs.0
IV. TRANSVERSE ǫ -NEAR-ZERO CONDITION Finally, we describe an interesting frequency region at which ǫ ′ t changes sign and it assumes values veryclose to zero. We denote this regime as transverse epsilon near zero (TENZ), which is manifested underthe condition σ ′′ ≈ − ωǫ ǫ d d, i.e., when graphene sheet’s inductive susceptance compensates for thesmall capacitive susceptance of each dielectric layer. We show in Fig. 9(a) and (b), the level of biasingpotential ( µ c ) required to achieve the TENZ condition at a given frequency, and the corresponding ǫ ′′ t ,respectively. We note that the required bias voltage for TENZ at a certain frequency decreases for thinnerunit cells, i.e., smaller d , however losses become larger due to increased graphene sheet density, especiallyat low frequencies. For example when d = 50 nm, we require µ c to be tuned to 0.1 eV in order to achievethe TENZ condition at 15 THz, and we have ǫ ′′ t ≈ . , whereas if the metamaterial is designed with d = 200 nm, the amount of bias required to realize TENZ condition at the same frequency is about 0.2eV and the losses are lower ǫ ′′ t ≈ . . In view of such observations one can easily identify the tuningranges and show that for smaller unit cell thickness the tuning range is larger but one must tolerate thelosses in such design.When considering wave propagation at that particular condition, and if losses are to be neglected with-out compromising the generality of the conclusions, the quasi-static approximation derived from EMAEq. (3) reveals a β z − k t dispersion relation with very small slope, i.e., the dispersion curve is almostflat. However, at higher k t the EMA approximations become inaccurate, and β z grows until it reachesthe Brillouin zone edge − π/d . The accurate wavevector dispersion of TM z waves according to Blochtheory, using Eq. (4) and Z TM d = κ d / ( ωǫ ǫ d ) , is given by cos k z d = cos κ d d + j ( σ ′ + jσ ′′ )2 κ d ωǫ ǫ d sin κ d d. (5)The condition ǫ ′ t ≈ is satisfied when ωǫ ǫ d d ≈ − σ ′′ , and it leads to cos k z d ≈ cos κ d d + κ d d κ d d + j (cid:12)(cid:12)(cid:12)(cid:12) σ ′ σ ′′ (cid:12)(cid:12)(cid:12)(cid:12) κ d d sin κ d d. (6)This latter dispersion equation is further simplified under the small argument approximation, | κ d d | ≪ as cos k z d ≈ j ( κ d d ) (cid:12)(cid:12)(cid:12)(cid:12) σ ′ σ ′′ (cid:12)(cid:12)(cid:12)(cid:12) + O ( | κ d d | ) . (7) (cid:2) (cid:10)(cid:2) (cid:11)(cid:2) (cid:12)(cid:2) (cid:13)(cid:2)(cid:2)(cid:2)(cid:14)(cid:10)(cid:2)(cid:14)(cid:11)(cid:2)(cid:14)(cid:12)(cid:2)(cid:14)(cid:13)(cid:2)(cid:14)(cid:15) (cid:1) (cid:2) (cid:3) (cid:4) (cid:1) (cid:2) (cid:3)(cid:4)(cid:5)(cid:6)(cid:7) (cid:5) (cid:6) (cid:3) (cid:8) (cid:9) (cid:7) (cid:2) (cid:10)(cid:2) (cid:11)(cid:2) (cid:12)(cid:2) (cid:13)(cid:2)(cid:2)(cid:14)(cid:15)(cid:2)(cid:14)(cid:13)(cid:2)(cid:14)(cid:12)(cid:2)(cid:14)(cid:11)(cid:2)(cid:14)(cid:10)(cid:2) (cid:2) (cid:4) (cid:3)(cid:3) (cid:7) (cid:16)(cid:1)(cid:16)(cid:15)(cid:2)(cid:16)(cid:17)(cid:18) (cid:7) (cid:16)(cid:1)(cid:16)(cid:10)(cid:2)(cid:2)(cid:16)(cid:17)(cid:18) (cid:7) (cid:16)(cid:1)(cid:16)(cid:11)(cid:2)(cid:2)(cid:16)(cid:17)(cid:18) (cid:7) (cid:16)(cid:1)(cid:16)(cid:15)(cid:2)(cid:2)(cid:16)(cid:17)(cid:18) (cid:1) (cid:2) (cid:3) (cid:4) (cid:1) (cid:2) (cid:3)(cid:4)(cid:5)(cid:6)(cid:7) (cid:19)(cid:20)(cid:21) (cid:19)(cid:22)(cid:21) Figure 9: The zero-crossing frequency of ǫ ′ t evaluated according to EMA f ǫ ′ t =0 = σ ′′ / (2 πǫ ǫ d d ) varying thechemical potential, for various thicknesses d . (b) Imaginary part of the transverse permittivity ǫ ′′ t evaluated at f ǫ ′ t =0 . (cid:20)(cid:18) (cid:18) (cid:20)(cid:18) (cid:20) (cid:20)(cid:18) (cid:25) (cid:18)(cid:18)(cid:23)(cid:26)(cid:20)(cid:20)(cid:23)(cid:26)(cid:25) (cid:1) (cid:2) (cid:3) (cid:4) (cid:5) (cid:20)(cid:18) (cid:18) (cid:20)(cid:18) (cid:20) (cid:20)(cid:18) (cid:25) (cid:30)(cid:26)(cid:30)(cid:31)(cid:30)(cid:32)(cid:30)(cid:25)(cid:30)(cid:20)(cid:18) (cid:4) (cid:6) (cid:3)(cid:4) (cid:5) (cid:7) (cid:2) (cid:3) (cid:4) (cid:5) (cid:1)(cid:2)(cid:3)(cid:4)(cid:5)(cid:6)(cid:7)(cid:8)(cid:9)(cid:10)(cid:11)(cid:8)(cid:8)(cid:9)(cid:3)(cid:12)(cid:9)(cid:10) (cid:13)(cid:11)(cid:14)(cid:15)(cid:16) (cid:8)(cid:9) (cid:12) (cid:17)(cid:18) (cid:13)(cid:11)(cid:14)(cid:15)(cid:16) (cid:8)(cid:9) (cid:12) (cid:19)(cid:18) (cid:20)(cid:18)(cid:21)(cid:13)(cid:1)(cid:22) (cid:20)(cid:20)(cid:23)(cid:24)(cid:21)(cid:13)(cid:1)(cid:22) (cid:20)(cid:25)(cid:21)(cid:13)(cid:1)(cid:22) (cid:20)(cid:26)(cid:21)(cid:13)(cid:1)(cid:22) (cid:27)(cid:28)(cid:29) (cid:27)(cid:28)(cid:29) (cid:4) (cid:6) (cid:3)(cid:4) (cid:5) Figure 10: Isofrequency wavevector dispersion in the TENZ, hyperbolic, and elliptic regimes, showing both (a) β z and (b) α z calculated by Bloch theory at four different frequencies (10, 11.9, 12, 15 THz), when µ c = 0 . eV. The imaginary term in Eq. (7) is negligible since (cid:12)(cid:12)(cid:12) ( κ d d ) σ ′ / σ ′′ (cid:12)(cid:12)(cid:12) ≪ for graphene-dielectric multi-layer with a subwavelength period, and therefore one simply obtains k z ≈ , far enough from the Bril-louin zone edge. Therefore, the TENZ condition ǫ ′ t ≈ , implies a flat isofrequency dispersion diagramwith small k z over a wide range of k t . We report in Fig. 10(a,b) the isofrequency wavevector dispersionat four different frequencies, at which we show hyperbolic dispersion (10 THz with ǫ t ≃ − . − j . ),elliptic dispersion (15 THz with ǫ t ≃ . − j . ), and the TENZ transitional state (at 11.9 THz and12 THz, with ǫ t ≃ − . − j . and ǫ t ≃ . − j . , respectively), where both β z and α z forall cases are normalized by k . In Fig. 10(a) one can observe that the slope of the β z − k t dispersion isreduced when | ǫ ′ t | is much smaller than unity, as also predicted analytically in Eq. (7), still preservinglimited values of the attenuation constant α z . Note that the elliptic regime (at 15 THz) also shows avery low slope of the β z − k t dispersion, however the attenuation constant α z is large, because wavesare mainly evanescent for large k t . Fig. 10(a) shows that the TENZ regimes are responsible for almostflat propagation constant ( | β z /k | < ) up to k t ≃ k , with a moderately low attenuation constant α z . However, for larger k t , we observe that β z experiences a sharp increase towards the Brillouin zoneedge, together with an increase of the attenuation constant α z . In Fig. 10(b) we observe that the atten-uation constant exhibits significant difference for HM and TENZ regimes that requires some importantconsideration. Although the two TENZ cases have smaller ǫ ′′ t than the hyperbolic one (at 10 THz), theyexperience a higher attenuation than HM case for k t > √ ǫ d k , whereas the opposite relation is validfor k t < √ ǫ d k . Therefore we can observe the two trends: on one hand TENZ allows flatter β z − k t relation and a wider k t spectrum than a fully hyperbolic regime, on the other hand the hyperbolic regimeexhibits smaller attenuation constant α z than the TENZ cases. Note also that the TENZ is a transitionalstate toward elliptic dispersion, at which the attenuation α z becomes even higher for k t > √ ǫ d k , andforward waves ( β z α z > ) can propagate for k t < √ ǫ d k with low attenuation constant.It has been shown in [41, 42] that isotropic epsilon-near-zero (IENZ) material inside a waveguidesupporting TE modes is able to tunnel electromagnetic waves. Here we elaborate on TENZ materialsat far- and mid-infrared frequencies designed using graphene-dielectric multilayers and explore theircapabilities of tunneling electromagnetic waves [43]. Consider an electrically-thin slab of thickness h made by either a TENZ ( ǫ t ≈ , ǫ z = 0 ) or an IENZ ( ǫ t = ǫ z = ǫ r ≈ ) material in free space. Under TE z wave incidence, TENZ and IENZ slabs provide an identical response and the reflection from such slabscan be set arbitrarily small by decreasing their thickness, as reported in [44]. However, for TM z obliqueplane waves impinging on a lossless IENZ semi-infinite material, total reflection occurs for angles greater2than the critical angle k ct /k = sin θ ci = √ ǫ r ≈ . For an electrically-thin IENZ slab, transmission ofTM z plane wave takes place for small angles of incidence ( < θ i < θ ci , where θ ci is considerably small)due to evanescent waves exhibiting frustrated multiple reflections at the slab interfaces. By including theeffect of losses in IENZ slabs, absorption and local electric field enhancement were reported for specificincident angles θ i > θ ci in [45]. Instead, we provide here the TM z reflection and transmission coefficients( R TENZTM and T TENZTM ) for a thin TENZ slab R TENZTM = ζ Z + ζ , T TENZTM = 2 Z Z + ζ , (8)where ζ = jh ( k − k t /ǫ z ) ωǫ , Z = p k − k t ωǫ . (9)Therefore upon having a thin slab of TENZ material, ζ can be made small enough (due to the existenceof finite, non vanishing ǫ z ) in order to observe complete transmission for oblique TM z waves with awide range of incidence angles. This is in contrast to what happens for the IENZ case with ǫ z assumingnear-zero values; which implies that transmission only occurs around k t ≈ . We show in Fig. 11 thereflection and transmission at 37 THz, by a TENZ material with ǫ t = − . and ǫ z = 2 . , and by anIENZ material with ǫ r = − . , assuming in both cases negligible losses. It is clear that the IENZmaterial exhibits a very narrow transmission around θ i ≈ ◦ only due to evanescent waves (permittivityhas a negative value) tunneling through the subwavelength slab [43, 45], and the transmission windowdramatically diminishes as ǫ r approaches zero or h increases, in accordance with the trend observed in[44]. On the contrary, the TENZ slab exhibits large and stable transmission over a wide range of incidenceangles, inherently complying with the flat wavevector dispersion relation in Eq. (7). Also, one shouldpoint out that the TM z transmission in TENZ materials occurs up to much larger incidence angles thanTE z transmission, which is identical to an IENZ slab’s TE z transmission discussed in [44]. In principlethe different properties illustrated in the preceding simple example reveal the advantage of TENZ materialover conventional IENZ material in enhancing transmission under oblique TM z plane wave incidence.For a more practical comparison, we report in Fig. 12 the transmission and reflection for two possibleTENZ and IENZ materials at mid-infrared. We consider a TENZ made of graphene-dielectric multilayerbiased with µ c = 0 . eV, accounting for losses, and having total thickness of h = N d where d = 50 (cid:1)(cid:2) (cid:1)(cid:2)(cid:3)(cid:4)(cid:5)(cid:2)(cid:6)(cid:7)(cid:2) (cid:1)(cid:2) (cid:1)(cid:2)(cid:5)(cid:2)(cid:6)(cid:7)(cid:2) (cid:1)(cid:2) (cid:1)(cid:2)(cid:8)(cid:2)(cid:6)(cid:7)(cid:2) (cid:15)(cid:16)(cid:3) (cid:15)(cid:17)(cid:3) (cid:15)(cid:18)(cid:3) (cid:3) (cid:18)(cid:3) (cid:17)(cid:3) (cid:16)(cid:3)(cid:3)(cid:3)(cid:4)(cid:8)(cid:3)(cid:4)(cid:19)(cid:3)(cid:4)(cid:17)(cid:3)(cid:4)(cid:20)(cid:5) (cid:3) (cid:4) (cid:9)(cid:10)(cid:11)(cid:12)(cid:4)(cid:13) (cid:14) (cid:5) (cid:6) (cid:15)(cid:16)(cid:3) (cid:15)(cid:17)(cid:3) (cid:15)(cid:18)(cid:3) (cid:3) (cid:18)(cid:3) (cid:17)(cid:3) (cid:16)(cid:3)(cid:3)(cid:3)(cid:4)(cid:8)(cid:3)(cid:4)(cid:19)(cid:3)(cid:4)(cid:17)(cid:3)(cid:4)(cid:20)(cid:5) (cid:3) (cid:4) (cid:9)(cid:10)(cid:11)(cid:12)(cid:4)(cid:13) (cid:14) (cid:7) (cid:6) (cid:21)(cid:22)(cid:23) (cid:21)(cid:24)(cid:23) Figure 11: Different characteristics of TM z plane wave (a) reflection and (b) transmission from a thin slab made by aTENZ material (solid lines) and IENZ material (dashed lines) at 37 THz. Material losses in this example are assumednegligible. The TENZ material exhibits much wider and flatter parameters varying angle of incidence than the IENZmaterial. (cid:1)(cid:2) (cid:1)(cid:2)(cid:3)(cid:4)(cid:5)(cid:2)(cid:6)(cid:7)(cid:2) (cid:1)(cid:2) (cid:1)(cid:2)(cid:5)(cid:2)(cid:6)(cid:7)(cid:2) (cid:1)(cid:2) (cid:1)(cid:2)(cid:8)(cid:2)(cid:6)(cid:7)(cid:2) (cid:13)(cid:14)(cid:3) (cid:13)(cid:15)(cid:3) (cid:13)(cid:16)(cid:3) (cid:3) (cid:16)(cid:3) (cid:15)(cid:3) (cid:14)(cid:3)(cid:3)(cid:3)(cid:4)(cid:8)(cid:3)(cid:4)(cid:17)(cid:3)(cid:4)(cid:15)(cid:3)(cid:4)(cid:18)(cid:5) (cid:19) (cid:3) (cid:4) (cid:5) (cid:6) (cid:20)(cid:21)(cid:22)(cid:23)(cid:4)(cid:24)(cid:13)(cid:14)(cid:3) (cid:13)(cid:15)(cid:3) (cid:13)(cid:16)(cid:3) (cid:3) (cid:16)(cid:3) (cid:15)(cid:3) (cid:14)(cid:3)(cid:3)(cid:3)(cid:4)(cid:8)(cid:3)(cid:4)(cid:17)(cid:3)(cid:4)(cid:15)(cid:3)(cid:4)(cid:18)(cid:5) (cid:19) (cid:7) (cid:4) (cid:5) (cid:6) (cid:20)(cid:21)(cid:22)(cid:23)(cid:4)(cid:24) (cid:9)(cid:10)(cid:11) (cid:9)(cid:12)(cid:11) Figure 12: TM z plane wave (a) reflection and (b) transmission from a slab made by graphene-dielectric layers with d = 50 nm and h = Nd (solid lines, using transfer matrix analysis) and an isotropic InAsSb slab of thickness h (dashed lines) at 37 THz. nm, at 37 THz. Under these conditions EMA estimates ǫ t ≈ − . − j . as seen from Fig. 2.The IENZ material is assumed to be a heavily n-doped InAsSb semiconductor [46], which is engineeredvia doping to exhibit low loss IENZ in this frequency range, i.e., ǫ InAsSb ≈ − . − j . at ≈ ǫ r is increased, angular transmission is slightly broadened, especially as h increases.This indicates an advantage of using the graphene-based TENZ materials in tuning and enhancing TM z plane wave transmission for wide angles of incidence. On the other hand, losses in natural materials orengineered metamaterials that exhibit IENZ behavior degrades the performance considerably, and mayrequire integration of gain materials as in [45]. V. CONCLUSION
We have reported a HM implementation at far- and mid-infrared frequencies that comprises graphene-dielectric layers, and showed that EMA describes the hyperbolic wavevector dispersion as well as thetransition to elliptic regime for specific conditions. Hyperbolic dispersion have manifested mainly atfar-infrared frequencies, where we have investigated the propagating spectrum properties and discussedthe effect of losses. We also showed that hyperbolic and elliptic dispersion regimes are associated tobackward and forward wave propagation, respectively. We have explored the tuning opportunities anddesign considerations of the structure, as well as the translation from hyperbolic to elliptic wavevectordispersion, and demonstrated a transitional state, TENZ, at which the wavevector dispersion diagrambecomes very flat. Furthermore, we have demonstrated that a thin slab made by a TENZ material becomestransparent to both TE z and TM z plane wave, with the interesting characteristic that the transmissionand reflection of TM z waves are stable with respect to the incident angle, in contrast to what happensin conventional IENZ materials. This property can be utilized in designing ultra-thin films for tunableinfrared applications. [1] D. R. Smith and D. Schurig, “Electromagnetic wave propagation in media with indefinite permittivity andpermeability tensors,” Phys. Rev. Lett. (7), 077405 (2003). [2] C. Guclu, S. Campione, and F. Capolino, “Hyperbolic metamaterial as super absorber for scattered fields gen-erated at its surface,” Phys. Rev. B (7), 205130 (2012).[3] O. Kidwai, S. V. Zhukovsky, and J. E. Sipe, “Effective-medium approach to planar multilayer hyperbolic meta-materials: Strengths and limitations,” Phys. Rev. A (5), 053842 (2012).[4] O. Kidwai, S. V. Zhukovsky, and J. E. Sipe, “Dipole radiation near hyperbolic metamaterials: applicability ofeffective-medium approximation,” Opt. Lett. (5), 2530–2532 (2011).[5] E. Narimanov, M. A. Noginov, H. Li, and Y. Barnakov, “Darker than black: Radiation-absorbing metamate-rial,” in The Conference on Lasers and Electro-Optics(CLEO)/The International Quantum Electronics Confer-ence(QELS) (OSA, 2010), p. QPDA6.[6] C. Rizza, A. Ciattoni, E. Spinozzi, and L. Columbo, “Terahertz active spatial filtering through optically tunablehyperbolic metamaterials,”
Opt. Lett. (16), 3345–3347 (2012).[7] G. Naik and A. Boltasseva, “Semiconductors for plasmonics and metamaterials,” physica status solidi RPL ,295–297 (2010).[8] I. Nefedov and S. Tretyakov, “Ultrabroadband electromagnetically indefinite medium formed by aligned carbonnanotubes,” Phys. Rev. B (11), 113410 (2011).[9] Z. Jacob, J. Kim, G. Naik, A. Boltasseva, E. Narimanov, and V. Shalaev, “Engineering photonic density of statesusing metamaterials,” Appl. Phys. B (1), 215–218 (2010).[10] K. J. Webb and M. Yang, “Subwavelength imaging with a multilayer silver film structure,”
Opt. Lett. (14),2130–2132 (2006).[11] C. L. Cortes, W. Newman, S. Molesky, and Z. Jacob, “Quantum nanophotonics using hyperbolic metamaterials,” J. Opt. (6), 063001 (2012).[12] J. Bénédicto, E. Centeno, and A. Moreau, “Lens equation for flat lenses made with hyperbolic metamaterials,” Opt. Lett. (22), 4786–4788 (2012).[13] J. Kim, V. Drachev, Z. Jacob, G. Naik, A. Boltasseva, E. Narimanov, and V. Shalaev, “Improving the radiativedecay rate for dye molecules with hyperbolic metamaterials,” Opt. Express (7), 8100–8116 (2012).[14] H. N. S. Krishnamoorthy, Z. Jacob, E. Narimanov, I. Kretzschmar, and V. M. Menon, “Topological transitionsin metamaterials,” Science (6078), 205–209 (2012).[15] S. M. Hashemi and I. S. Nefedov, “Wideband perfect absorption in arrays of tilted carbon nanotubes,”
Phys.Rev. B (19), 195411 (2012).[16] G. W. Hanson, A. B. Yakovlev, and A. Mafi, “Excitation of discrete and continuous spectrum for a surfaceconductivity model of graphene,” J. Appl. Phys. (11), 114305 (2011).[17] A. Vakil and N. Engheta, “Transformation optics using graphene,”
Science (6035), 1291–1294 (2011).[18] G. Lovat, “Equivalent circuit for electromagnetic interaction and transmission through graphene sheets,”
IEEETrans. Electromag. Compat. (1), 101 –109 (2012).[19] M. Tamagnone, J. Gomez-Diaz, J. Mosig, and J. Perruisseau-Carrier, “Analysis and design of terahertz antennasbased on plasmonic resonant graphene sheets,” J. Appl. Phys. (11), 114915 (2012).[20] S. Thongrattanasiri, F. H. L. Koppens, and F. J. Garcia de Abajo, “Complete optical absorption in periodicallypatterned graphene,”
Phys. Rev. Lett. (4), 047401 (2012).[21] A. Fallahi and J. Perruisseau-Carrier, “Design of tunable biperiodic graphene metasurfaces,”
Phys. Rev. B (19), 195408 (2012).[22] A. Fallahi and J. Perruisseau-Carrier, “Manipulation of giant faraday rotation in graphene metasurfaces,” Appl.Phys. Lett. (23), 231605–231605 (2012).[23] P.-Y. Chen and A. Alu, “Atomically thin surface cloak using graphene monolayers,”
ACS Nano (7), 5855–5863(2011).[24] I. V. Iorsh, I. S. Mukhin, I. V. Shadrivov, P. A. Belov, and Y. S. Kivshar, “Hyperbolic metamaterials based onmultilayer graphene structures,” Phys. Rev. B (7), 075416 (2013).[25] M. A. K. Othman, C. Guclu, and F. Capolino, “Graphene-based tunable hyperbolic metamaterials and enhancednear-field absorption,” Opt. Express (6), 7614–7632 (2013).[26] A. Reina, X. Jia, J. Ho, D. Nezich, H. Son, V. Bulovic, M. S. Dresselhaus, and J. Kong, “Large area, few-layergraphene films on arbitrary substrates by chemical vapor deposition,” Nano lett. (1), 30–35 (2008).[27] X. Li, W. Cai, J. An, S. Kim, J. Nah, D. Yang, R. Piner, A. Velamakanni, I. Jung, E. Tutuc et al. , “Large-areasynthesis of high-quality and uniform graphene films on copper foils,” Science (5932), 1312–1314 (2009).[28] X. Li, Y. Zhu, W. Cai, M. Borysiak, B. Han, D. Chen, R. D. Piner, L. Colombo, and R. S. Ruoff, “Transferof large-area graphene films for high-performance transparent conductive electrodes,”
Nano Lett. (12), 4359–4363 (2009). PMID: 19845330.[29] C. Dean, A. Young, I. Meric, C. Lee, L. Wang, S. Sorgenfrei, K. Watanabe, T. Taniguchi, P. Kim, K. Shepard et al. , “Boron nitride substrates for high-quality graphene electronics,” Nature nanotechn. (10), 722–726(2010).[30] J. Wang, Y. Xu, H. Chen, and B. Zhang, “Ultraviolet dielectric hyperlens with layered graphene and boronnitride,” J. Mater. Chem. (31), 15863–15868 (2012).[31] J. J. Miret, C. J. Zapata-Rodriguez, Z. Jaksic, S. Vukovic, and M. R. Belic, “Substantial enlargement of angularexistence range for dyakonov-like surface waves at semi-infinite metal-dielectric superlattice,” J. Nanophoton. (1), 063525–063525 (2012).[32] K. Novoselov, A. Geim, S. Morozov, D. Jiang, Y. Zhang, S. Dubonos, I. Grigorieva, and A. Firsov, “Electricfield effect in atomically thin carbon films,” Science (5696), 666–669 (2004).[33] V. P. Gusynin, S. G. Sharapov, and J. P. Carbotte, “Sum rules for the optical and hall conductivity in graphene,”
Phys. Rev. B (16), 165407 (2007).[34] G. W. Hanson, “Dyadic green’s functions and guided surface waves for a surface conductivity model ofgraphene,” J. Appl. Phys. (6), 064302 (2008).[35] D. Pozar,
Microwave engineering (John Wiley & Sons, 2009).[36] L. A. Falkovsky and S. S. Pershoguba, “Optical far-infrared properties of a graphene monolayer and multilayer,”
Phys. Rev. B (15), 153410 (2007).[37] L. Hao and L. Sheng, “Optical conductivity of multilayer graphene,” Solid State Communications (43),1962–1966 (2009).[38] S. Campione, S. Steshenko, M. Albani, and F. Capolino, “Complex modes and effective refractive index in 3dperiodic arrays of plasmonic nanospheres,”
Opt. Express (27), 26027–26043 (2011).[39] A. Chebykin, A. Orlov, C. Simovski, Y. Kivshar, and P. Belov, “Nonlocal effective parameters of multilayeredmetal-dielectric metamaterials,” Phys. Rev. B (11), 115420 (2012).[40] L. Felsen and N. Marcuvitz, Radiation and Scattering of Waves , pp.745-747, IEEE Press, NJ (1994).[41] M. Silveirinha and N. Engheta, “Tunneling of electromagnetic energy through subwavelength channels andbends using ε -near-zero materials,” Phy. rev. lett. (15), 157403 (2006).[42] B. Edwards, A. Alu, M. E. Young, M. Silveirinha, and N. Engheta, “Experimental verification of epsilon-near-zero metamaterial coupling and energy squeezing using a microwave waveguide,” Phys. rev. lett. (3), 33903(2008).[43] J. Luo, Y. Xu, H. Chen, B. Hou, W. Lu, and Y. Lai, “Oblique total transmissions through epsilon-near-zerometamaterials with hyperbolic dispersions,”
EPL (Europhys. Lett.) (4), 44001 (2013).[44] A. Alu, M. G. Silveirinha, A. Salandrino, and N. Engheta, “Epsilon-near-zero metamaterials and electromag-netic sources: Tailoring the radiation phase pattern,”
Phys. Rev. B (15), 155410 (2007).[45] S. Campione, D. de Ceglia, M. A. Vincenti, M. Scalora, and F. Capolino, “Electric field enhancement in epsilon-near-zero slabs under tm-polarized oblique incidence,” Phys. Rev. B (3), 035120 (2013).[46] D. C. Adams, S. Inampudi, T. Ribaudo, D. Slocum, S. Vangala, N. A. Kuhta, W. D. Goodhue, V. A. Podolskiy,and D. Wasserman, “Funneling light through a subwavelength aperture with epsilon-near-zero materials,” Phys.Rev. Lett.107