Realizing negative refraction of highly squeezed polaritons, especially that supported by two-dimensional (2D) materials [1–6], such as graphene plasmon polaritons, is an important step toward the active manipulation of light at the extreme nanoscale [7–9] and can promise many photonic and optoelectronic applications [10–17]. In 2017, the phenomenon of all-angle negative refraction between highly squeezed isotropic graphene plasmons and hexagonal boron nitride’s (BN) phonon polaritons, with their in-plane polaritonic wavelengths squeezed by a factor over 100, is theoretically shown possible in the graphene-BN heterostructures [18]. By following [18], henceforth we use the parameter of squeezing factor to define the ratio between the photon wavelength in free space and the in-plane polaritonic wavelength. The extreme spatial confinement of these highly squeezed isotropic polaritons, however, also limits our ability to tailor their dispersion relations [18, 19]; in other words, their effective negative refractive index exists only within a certain frequency range. Consequently, the phenomenon of all-angle negative refraction of these highly squeezed isotropic polaritons is restricted to only a narrow frequency range (near 23 THz) with a bandwidth of less than 1 THz [18]. The narrow bandwidth unavoidably hinders the demonstration of this exotic phenomenon in experiments and thus the potential applications. It remains a challenge to realize the all-angle negative refraction of highly squeezed polaritons within a broad frequency range.

Here we predict the phenomenon of broadband all-angle negative refraction of highly squeezed hyperbolic polaritons, with their squeezing factors over 100, in the infrared regime. The highly squeezed hyperbolic polaritons can be supported by nanostructured 2D materials such as graphene nanoribbon arrays [20, 21] or BN nanoribbon arrays [22, 23], by 2D materials supported by nanostructured substrates [13, 24], or by naturally existing hyperbolic 2D materials such as black phosphorous [2, 25]. In addition, we note that the hyperbolic polaritons can also be supported by metal based metasurfaces at visible [26, 27] and even microwave [28] regimes. The peculiar isofrequency contour of hyperbolic polaritons [2, 20, 29] gives us an extra freedom to tailor their in-plane propagation direction and thus the flexibility to realize the all-angle negative refraction of highly squeezed polaritons in a wide frequency range. As an example, for the all-angle negative refraction of hyperbolic polaritons shown below, the working bandwidth can vary from several tens of THz to over a hundred of THz by simply tuning the chemical potential of 2D materials. Our work thus indicates that 2D materials are a versatile platform for the design of advanced nanometasurfaces and nanoimaging elements, which are of fundamental importance for the active control of light at the nanoscale.

We note that there are other types of negative refraction studied in the platform of 2D materials, including the negative refraction of electrons [30, 31], negative refraction of light propagating through a monolayer graphene [32], and plasmonic [33, 34] and optical [35, 36] negative refraction in 3D bulk materials (i.e., 2D material-based periodic structures). However, these electromagnetic refractions [32–36] occur out of the plane of 2D materials and not in plane, thus without taking advantage of the high spatial confinement of polaritons that is critical in this work. In addition, the all-angle negative refraction of graphene plasmons was studied in the hybrid graphene-photonic crystal structures [37], but the working bandwidth is still far less than 1 THz for a given chemical potential of graphene, along with the squeezing factor less than 20.

To highlight the underlying physics, we begin with the dispersion of hyperbolic polaritons supported by a uniaxial metasurface. The metasurface, such as that in the left region of Figure 1(a), can be modelled by an effective anisotropic surface conductivity of {\stackrel{̿}{\sigma }}_l=[{\sigma }_{xx,l},{\sigma }_{yy,l}]. Contrary to the isotropic metasurface (i.e., {\sigma }_{xx,l}={\sigma }_{yy,l}) which supports the propagation of either transverse-magnetic (TM) [38, 39] or transverse-electric (TE) polaritons [40–43], the uniaxial metasurface supports the hybrid TM-TE polaritons [20, 24, 44–46]. From the electromagnetic theory [47], the dispersion relation for hybrid polaritons [44–46] is derived as

In the above, {\epsilon }_{r\mathrm{1,2}} is the relative permittivity of region above/below the metasurface; k_{z\mathrm{1,2}}=\sqrt{\left({\omega }^2/c^2\right){\epsilon }_{r\mathrm{1,2}}-k_x^2-k_y^2} and {\stackrel{-}{k}}_{\left|\right|}=\widehat{x}{k}_x+\widehat{y}{k}_y are the out-of-plane and in-plane wavevectors of the hybrid polariton, respectively; \varphi is the angle between {\stackrel{-}{k}}_{\left|\right|} and \widehat{x}; \omega is the angular frequency; {\epsilon }_0, {\mu }_0, and c are the permittivity, permeability, and light speed in free space, respectively. It is worth noting that the left side of (1) is only relevant to the TM field components of hybrid polaritons, while the denominator of the right side of (1) only corresponds to the TE field components. To some extent, the numerator of the right side of (1) denotes the coupling strength between TM and TE field components.

Importantly, for the highly squeezed hybrid polaritons (i.e., \left|{\stackrel{-}{k}}_{\left|\right|}\right|\gg \omega /c), the right side of (1) is shown to be much smaller than 1, due to the difference of the spatial confinement between pure TE and TM polaritons; see Supplementary Materials. This way, the dispersion relation of hybrid polaritons in (1) can be compactly reduced to [20, 24]

In other words, the highly squeezed hybrid TM-TE polaritons are dominant by the TM field components. The isofrequency contour governed by (2) is hyperbolic, when .

Figure 1(b) shows the hyperbolic isofrequency contour of highly squeezed polaritons supported by a graphene metasurface at 15 THz. As a conceptual demonstration, here we use the nanostructured graphene, i.e., graphene nanoribbon arrays, to create the hyperbolic metasurface; see the schematic structure in Figure 1(a). According to [20], when the pitch L of nanoribbons is much smaller than the polaritonic wavelength {\lambda }_{\mathrm{p}\mathrm{o}\mathrm{l}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{t}\mathrm{o}\mathrm{n}}, i.e., L\ll {\lambda }_{\mathrm{p}\mathrm{o}\mathrm{l}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{t}\mathrm{o}\mathrm{n}}, the effective medium theory can be applied to describe the graphene metasurface. Then the effective anisotropic surface conductivity of graphene metasurface can be described by {\sigma }_{xx,l}=L{\sigma }_s{\sigma }_C/\left(W{\sigma }_C+\left(L-W\right){\sigma }_s\right) and {\sigma }_{yy,l}={\sigma }_s\left(W/L\right), where W is the width of nanoribbon, {\sigma }_C=-i\left(\omega {\epsilon }_{0}L/\pi \right)\ln [\mathrm{c}\mathrm{s}\mathrm{c}(\pi (L-W)/2L)] is an equivalent conductivity associated with the near-field coupling between adjacent nanoribbons, and {\sigma }_s is the surface conductivity of monolayer graphene modelled by the Kubo formula [20, 48]. Here the nanostructured graphene has a chemical potential of 0.1 eV, a conservative electron mobility of 10000 cm²V⁻¹s⁻¹ [38, 39], a pitch of L=\mathrm{30} nm, and a width of W=\mathrm{20} nm. The region below/above the graphene metasurface is the dielectric substrate (e.g., SiO₂) and air, respectively.

For the emergence of refraction phenomenon, the graphene metasurface in the right region in Figure 1(a) shall be different from the left region. One simple way is to rotate the graphene metasurface in the right region with a certain angle with respect to the left region, such as 90⁰ shown in Figure 1(a). This way, the graphene metasurface in the right region can be characterized by a surface conductivity of {\stackrel{̿}{\sigma }}_r=[{\sigma }_{xx,r},{\sigma }_{yy,r}], where {\sigma }_{xx,r}={\sigma }_{yy,l} and {\sigma }_{yy,r}={\sigma }_{xx,l}. By applying the conservation law for wavevectors parallel to the boundary between left and right regions (i.e., parallel to \widehat{y}), we can determine the eigenmodes for the incident polariton in the left region and the transmitted polariton in the right region; see Figure 1(b). Through judicious design (e.g., by rotating the right region 90° with respect to the left region), we can create the case that the directions of the component of group velocity parallel to the interface (i.e., the y-component) are opposite for the eigenmodes of incident and transmitted polaritons, as shown in Figure 1(b). Since the group velocity (instead of the phase velocity) determines the direction of power flow, the directions of the component of power flow parallel to the interface are also opposite for the incident and transmitted polaritons. This enables the negative refraction of hyperbolic polaritons at the interface between left and right regions. The above underlying mechanisms for negative refraction of polaritons are the same as [18, 49]. In addition, we note that the hyperbolic polaritons in the left region can only propagate within a certain range of directions. Figure 1(b) indicates that, for these hyperbolic polaritons (i.e., incident from the left region with arbitrary angles), negative refraction can always happen at the boundary. Here we denote this phenomenon as the all-angle negative refraction [18, 49, 50] of hyperbolic polaritons (i.e., graphene plasmons).

Figure 1(c) numerically demonstrates the all-angle negative refraction of hyperbolic polaritons at 15 THz, by using the finite-element method (COMSOL Multiphysics). Hyperbolic polaritons are excited by a dipole source in the left region and propagate directionally towards the boundary. At the boundary, the polaritonic beams are negatively refracted. Moreover, these beams converge to form an image in the right region, which computationally validates the all-angle negative refraction. We note that the full width at half maximum (FWHM) of the image is only 0.035 μm [Fig. S2], which is less than 1/100 of the working wavelength (i.e., 20 μm) in free space and can enable deep-subwavelength imaging. Here the squeezing factor \left(\mathrm{R}\mathrm{e}(|{\stackrel{-}{k}}_{\left|\right|}|)\right)/\left(\omega /c\right)>\mathrm{100} [see Figure 1(b)] indicates that when compared with the wavelength in free space, the polaritonic wavelength is squeezed at least by a factor over 100. It is worth emphasizing that the all-angle negative refraction of highly squeezed hyperbolic polaritons with directional propagation, as an important advantage over that of isotropic polaritons which propagate omnidirectionally [18, 46, 47], might facilitate the design of novel compact guidance. We note that the imaging mechanism here is not a perfect image recovery. This is because the reflection at the boundary is unavoidable due to the impedance mismatch (caused by the modal mismatch) between the two hyperbolic metasurfaces, and the propagation loss of hyperbolic polaritons will degrade the imaging quality; see Supplementary Materials. To optimize the quality of formed image in the right region, one shall consider both the reflection and the propagation loss of polaritons [18, 46, 47].

Figure 2 shows the effective anisotropic conductivity of graphene metasurface, which can help to infer the working frequency range of the all-angle negative refraction of hyperbolic polaritons. As shown in Figure 2(a), the frequency range that has and thus supports the hyperbolic polaritons spans from 0 to 47 THz, when the chemical potential of graphene is 0.1 eV. It shall be noted that the phenomenon of all-angle negative refraction of hyperbolic polaritons can happen at arbitrary frequency within this frequency range, as long as the effective medium theory for graphene metasurface is valid (i.e., when L\ll {\lambda }_{\mathrm{p}\mathrm{o}\mathrm{l}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{t}\mathrm{o}\mathrm{n}}). To guarantee the validity of the effective medium theory for metasurfaces based on 2D materials, a small value of pitch L, which although might increase the complexity in structural fabrication [22, 23], can be adopted. This way, the working frequency range for the all-angle negative refraction of highly squeezed polaritons revealed here is no longer limited by the frequency range supporting negative-index polaritons (such as [18]). Consequently, the working frequency range can be broadband and actively controllable via tuning the chemical potential of 2D materials. In addition, the relative bandwidth for negative refraction of hyperbolic polaritons can be up to 2, if we define the relative bandwidth as \left(f_{\mathrm{m}\mathrm{a}\mathrm{x}}-f_{\mathrm{m}\mathrm{i}\mathrm{n}}\right)/\left((f_{\mathrm{m}\mathrm{a}\mathrm{x}}+f_{\mathrm{m}\mathrm{i}\mathrm{n}})/2\right), where f_{\mathrm{m}\mathrm{a}\mathrm{x}} and f_{\mathrm{m}\mathrm{i}\mathrm{n}} (f_{\mathrm{m}\mathrm{i}\mathrm{n}}\to 0) are the maximum and minimum frequencies supporting the negative refraction of hyperbolic polaritons.

Figure 2(b) shows the working bandwidth of all-angle negative refraction in graphene metasurfaces. The bandwidth of graphene metasurface having changes from 40 THz to over 100 THz, by increasing the chemical potential from 0.1 eV to 0.5 eV. Therefore, in principle, the bandwidth of all-angle negative refraction of hyperbolic polaritons can vary from several tens of THz to even over a hundred THz, by simply increasing the chemical potential of graphene. Interestingly, Figure 2(b) shows that the loss can increase the bandwidth when {\mu }_c is smaller than 0.18 eV; see analysis in Supplementary Materials.

To numerically validate the all-angle negative refraction of hyperbolic polaritons in a broad bandwidth, Figure 3 demonstrates this phenomenon at other frequencies, i.e., at 10 THz in Figure 3(a) and at 20 THz in Figure 3(b). The squeezing factors for hyperbolic polaritons at these two frequencies are both over 100; see Supplementary Materials. Therefore, by using hyperbolic metasurfaces based on 2D materials, we can extend the working bandwidth of all-angle negative refraction of highly squeezed polaritons to at least several tens of THz, which is favored for practical applications.

To further extend the bandwidth of all-angle negative refraction of highly squeezed polaritons, one may adopt the naturally anisotropic 2D materials to support tunable hyperbolic polaritons, such as those described in [2] including black phosphorous. From [2], these 2D materials can be directly characterized by an anisotropic surface conductivity, i.e.,where n is the concentration of electrons, m_j is the electron’s effective mass along the j direction, \tau is the relaxation time, {\omega }_j is the frequency of the onset of interband transitions for the j component of conductivity, s_j accounts for the strength of interband component, and \Theta \left(\omega -{\omega }_j\right) is a step function. Equation (3) circumvents the requirement (i.e., L\ll {\lambda }_{\mathrm{p}\mathrm{o}\mathrm{l}\mathrm{a}\mathrm{r}\mathrm{i}\mathrm{t}\mathrm{o}\mathrm{n}} discussed in the above) for the validity of effective medium theory for metasurfaces based on nanostructured 2D materials. For these anisotropic 2D materials, the value of n can be flexibly tunable via electrostatic gating, just like graphene, and {\omega }_x can be different from {\omega }_y. These, along with the broad class of anisotropic materials, give us the flexibility to realize the all-angle negative refraction of highly squeezed polaritons in a broad frequency range.

As a concrete example, Figure 4 shows the all-angle negative refraction of hyperbolic polaritons supported by the naturally anisotropic 2D materials. By following [2], here we set n={3\times 10}^{13}\hspace{0.17em}\hspace{0.17em}{cm}^{-2}, m_x={0.2m}_0, m_y=m_0, {\omega }_x=1\hspace{0.17em}\hspace{0.17em}\mathrm{e}\mathrm{V}, {\omega }_y=0.35\hspace{0.17em}\hspace{0.17em}\mathrm{e}\mathrm{V}, \tau =0.4\hspace{0.17em}\hspace{0.17em}\mathrm{p}\mathrm{s}, s_x=1.7s_0 and s_y=3.7s_0, where s_0=e^2/4\hslash and m_0 is the free-electron mass. The frequency range, which has and thus supports the hyperbolic polaritons, spans from 24 to 143 THz in Figure 4(a), with a bandwidth of 119 THz. This indicates that the negative refraction of hyperbolic polaritons can happen within the above frequency range; see the negative refraction of hyperbolic polaritons at 60 THz in Figure 4(b), for example.

In conclusion, we have revealed a viable way to realize the all-angle negative refraction of highly squeezed polariton in a broadband infrared regime, by utilizing hyperbolic metasurfaces based on 2D materials or naturally anisotropic 2D materials. Due to the combined advantages of highly directional propagation, active tunability, low loss, and ultrahigh confinement provided by hyperbolic polaritons in 2D materials, the broad class of 2D materials can provide a versatile platform for the manipulation of light-matter interaction at the extreme nanoscale and for the design of highly compact nanodevices and circuits.

The finite element simulation is implemented via the frequency domain simulation in the commercial software of COMSOL Multiphysics. To enable the high calculation accuracy, the 2D material is modeled as a surface, where an anisotropic surface conductivity is used to fulfill the conditions for discontinuities in the electromagnetic fields. The meshing resolution in the plane of 2D material is 5 nm. For Figures 1(c) and 3, a z-polarized dipole source and for Figure 4 a y-polarized electric dipole is placed in the left region at 5 nm above the 2D material; the phenomenon of all-angle negative refraction is pronounced, independent of the vertical position of the source. The fields in Figures 1(c), 3, and 4 are obtained at the plane with 2 nm above the plane of 2D material. For the clarity of conceptual demonstration, the value of \mathrm{R}\mathrm{e}({\sigma }_{yy,l}) is artificially set to be equal to \mathrm{R}\mathrm{e}({\sigma }_{xx,l}) in Figures 1(c), 3, and 4. For the cases considering the realistic material loss, the phenomenon of all-angle negative refraction is shown in Fig. S3. The phenomenon of all-angle negative refraction is also independent of the relative permittivity of the dielectric substrate {\epsilon }_{r2}, as shown in Fig. S6; in the main text, we set {\epsilon }_{r2}=3.6.