Hyperbolic metasurfaces exhibit unique dispersion and polarization properties, making them a promising platform for a plethora of photonic applications. At the same time, the ability to engineer the hyperbolicity via the predefined spectral positions of the metasurface resonances remains a notable challenge. Here, we analyze the dependencies of the spectral positions of the resonances corresponding to the limits of the hyperbolic regime for the metasurfaces based on square arrays of the rectangular nanopatches. We show that the spectral difference between the resonances increases linearly with stretching of the nanopatch, but this dependence becomes quadratic when the length of the stretched nanopatch exceeds 85% of the lattice constant, indicating the regime of extreme anisotropy. Finally, we demonstrate the characteristic feature of the engineered resonances by showing the canalization (divergenceless propagation) of the surface plasmon-polariton along the anisotropic nanopatch-based metasurface in the vicinity of the resonance. The results obtained may be used for the engineering of the anisotropic nanoparticle-based metasurfaces for a plethora of photonic applications.

Two-dimensional (2D) materials offer unique opportunities for photodetection, light emission, energy harvesting, and enhanced light–matter interactions.1,2 Even more interest brings the artificially engineered 2D micro- and nanostructures with on-demand properties, such as surface conductivity or surface polarizability, paving the way towards a plethora of specific applications and devices including lensing, holography, imaging, polarimetry, biosensing, etc.3–5 The rapidly developing use of 2D nanostructures poses new challenges for their proper engineering.

The particular interest arises for hyperbolic metasurfaces (HMSs), which are described by the 2D anisotropic tensor of the negative determinant within the local effective medium approximation.6 In other words, the components of the diagonalized 2D tensor have different signs leading to the topological transition from elliptical to hyperbolic shapes of the isofrequency contours.7–9 The 2D tensor of the surface conductivity is usually studied for the plasmonic metasurfaces,10,11
σ ^ = ( σ x 0 0 σ y ) .
(1)
In this case, the hyperbolic regime of an anisotropic metasurface is characterized by the negative product of the imaginary parts of the surface conductivity tensor main-diagonal components,
Im ( σ x ) Im ( σ y ) < 0.
(2)
The HMSs have shown a number of applications in both far- and near-field, including negative refraction, planar hyperlensing, enhanced emission, sensing, etc.12,13 In the near-infrared and visible spectrum, HMSs are represented by metallic gratings14 and lattices of nanoparticles.15–17 

While the HMS have the clear advantages for many applications, the selection of their specific designs requires the extensive numerical simulation, which can be time-consuming and require significant effort. Some specific applications requiring the proper engineering of HMS include, for instance, the high-directional excitation of surface waves,18,19 twisted hyperbolic optics,20–22 spin–orbit photonics,23–25 and mode engineering for high-Q resonators26,27 and sensors.28,29 For all these applications, it is highly important to set two well-separated resonances of the surface conductivity, namely, to adjust their wavelength and amplitude independently. These resonances may be characterized by the peaks in the reflectance spectra17,30 or the plasmon canalization regime, which is the collinear divergenceless propagation of the surface plasmon-polariton along the specific direction12,31 (Fig. 1). The plasmon canalization has been previously observed in hyperbolic metasurfaces,10,11,32 natural anisotropic materials,33–35 and twisted hyperbolic layers.20–22 

FIG. 1.

The artistic view of the surface plasmon canalization (shown by electric field distribution) along the hyperbolic metasurface based on gold nanopatches (marked by yellow bricks) in the vicinity of the resonant wavelength excited by the point vertical dipole (marked by the green arrow).

FIG. 1.

The artistic view of the surface plasmon canalization (shown by electric field distribution) along the hyperbolic metasurface based on gold nanopatches (marked by yellow bricks) in the vicinity of the resonant wavelength excited by the point vertical dipole (marked by the green arrow).

Close modal

In this work, we analyze the spectral positions of the resonances and the bandwidth of the hyperbolic regime for dozens of plasmonic metasurface designs based on gold nanopatches. Starting with nanopatches of a square cross section, we gradually stretch them into gratings, introducing anisotropy. Thus, a nanopatch-based design allows one to study the full dynamics of the transition between isotropic and extreme anisotropic cases for the HMS resonances. We show that the resonant wavelength for the case of the electric field oriented along the stretching direction and the spectral bandwidth between two resonances are initially linearly proportional to the stretching factor, but at a certain degree of the stretching (anisotropy), their dependence becomes quadratic. Then, we demonstrate the plasmon canalization in the vicinity of one of the resonances highlighting the relevance of the metasurface engineering for the in-plane optical signal transferring. Our results bring the significant contribution to the engineering of the resonances for the anisotropic 2D materials and nanostructures.

We consider the metasurface based on the periodic array of gold nanopatches with a height of h = 20 nm and a rectangular cross section ( a x × a y ), packed in a square array with a lattice constant of p = 300 nm at the interface of a quartz substrate ( n s = 1.45). The dispersion of gold is taken from the experimental data for the thin films.36 Initially, we consider a plasmonic metasurface based on the nanopatches with a square cross section with a side of a = 150 nm [the unit cell is shown in Fig. 2(a 1)]. Then, we stretch the nanopatches along the x-direction transforming the square cross section into the rectangular one [Figs. 2(a 2)2(a 4)]. The total area of the cross section is saved during the stretching; thus, the long and short sides of a rectangle correspond to a x = a η and a y = a / η, respectively, where η is the stretching factor characterizing the degree of the anisotropy for a metasurface. For further analysis, we introduce the parameter of a filling factor along the direction of the stretching defined as f x = a x / p, i.e., the long side of a nanopatch normalized per period. The unit cells, shown in Figs. 2(a 1)2(a 4), correspond to the stretching factors η = 1 , 1.4 , 1.8 , 2 and filling factors along the x-axis f x = 0.5 , 0.7 , 0.9 , 1, respectively. The first ( a x = a y = a) and last ( a x = p) unit cells, shown in Figs. 2(a 1) and 2(a 4), correspond to the isotropic (within effective local medium approach) and extreme anisotropic cases, respectively. Namely, the nanopatch-based metasurface in extreme anisotropic case represents a grating with a period p and a width a y.

FIG. 2.

(a 1)–(a 4) The metasurfaces’ unit cells of a period p = 300 nm consisting of the gold nanopatches with (a 1) a square cross section of a side a = 150 nm (isotropic case), (a 2) and (a 3) a rectangular cross section of the long and short sides a x = a η and a y = a / η (anisotropic case), and (a 4) grating of a width a y = a / η (extreme anisotropic case) characterized by the stretching factors (a 1) η = 1, (a 2) η = 1.4, (a 3) η = 1.8, and (a 4) η = 2. (b 1)–(b 4) The reflectance (solid lines) and transmittance (dashed lines) spectra of the metasurfaces shown in (a 1)–(a 4), respectively. The purple and blue lines in (b 2)–(b 4) correspond to the electric field orientation along x and y axes, respectively. (c 1)–(c 4) The dependencies of the real (dashed lines) and imaginary (solid lines) parts of the effective surface conductivity on the wavelength retrieved via the fitting of Eq. (4) using the data from (b 1)–(b 4), respectively. The blue, red, and green regions correspond to (i) Im ( σ x ) > 0 , Im ( σ y ) > 0, (ii) Im ( σ x ) < 0 , Im ( σ y ) < 0, and (iii) the hyperbolic regime Im ( σ x ) Im ( σ y ) < 0, respectively. The resonant wavelengths for the isotropic and anisotropic cases are shown by green ( λ 0), purple ( λ x), and blue ( λ y) dashed lines, respectively.

FIG. 2.

(a 1)–(a 4) The metasurfaces’ unit cells of a period p = 300 nm consisting of the gold nanopatches with (a 1) a square cross section of a side a = 150 nm (isotropic case), (a 2) and (a 3) a rectangular cross section of the long and short sides a x = a η and a y = a / η (anisotropic case), and (a 4) grating of a width a y = a / η (extreme anisotropic case) characterized by the stretching factors (a 1) η = 1, (a 2) η = 1.4, (a 3) η = 1.8, and (a 4) η = 2. (b 1)–(b 4) The reflectance (solid lines) and transmittance (dashed lines) spectra of the metasurfaces shown in (a 1)–(a 4), respectively. The purple and blue lines in (b 2)–(b 4) correspond to the electric field orientation along x and y axes, respectively. (c 1)–(c 4) The dependencies of the real (dashed lines) and imaginary (solid lines) parts of the effective surface conductivity on the wavelength retrieved via the fitting of Eq. (4) using the data from (b 1)–(b 4), respectively. The blue, red, and green regions correspond to (i) Im ( σ x ) > 0 , Im ( σ y ) > 0, (ii) Im ( σ x ) < 0 , Im ( σ y ) < 0, and (iii) the hyperbolic regime Im ( σ x ) Im ( σ y ) < 0, respectively. The resonant wavelengths for the isotropic and anisotropic cases are shown by green ( λ 0), purple ( λ x), and blue ( λ y) dashed lines, respectively.

Close modal
The HMSs under consideration are described via the effective surface conductivity tensor (1). We assume that the dispersion of the surface conductivity obeys the Drude–Lorentz model as follows:30,
σ x , y = i σ 0 x , y β x , y + i γ x , y ,
(3)
where σ 0 x , y is the normalized dimensionless amplitude of surface conductivity, β x , y = 1 Ω x , y 2; Ω x , y and γ x , y are the dimensionless resonant frequency and the bandwidth of the resonance, respectively, normalized per the angular frequency ω = 2 π c / λ, c is the speed of light, and λ is the operational wavelength. Deriving the Fresnel equations with a 2D conducting layer described via the surface conductivity (3), one can express the reflectance through the parameters of a surface conductivity as follows:17,30
R i = ( σ 0 i ) 2 + ( β i 2 + γ i 2 ) ζ 2 + 2 γ i σ 0 i ζ ( σ 0 i ) 2 + ( β i 2 + γ i 2 ) ζ + 2 + 2 γ i σ 0 i ζ + ,
(4)
where the index i means x and y components and ζ + , = n s ± 1. One can notice that the condition σ 0 i = 0 in Eq. (4) leads to the conventional formula for the Fresnel equations used for normal incidence.37 The reflectance and transmittance spectra of a normally incident plane wave with electric field oriented along x- and y-axis were calculated by the finite-element method with COMSOL Multiphysics for a single unit cell with imposed periodic boundary conditions [Figs. 2(b 1)2(b 4)]. Then, we retrieve the corresponding components of the surface conductivity tensor σ x and σ y by fitting the calculated reflectance spectra using the least-squares method according to Eq. (4); see Figs. 2(c 1)2(c 4). The effect of the finite thickness of a metasurface can be compensated by properly correcting the reflection phase.30 
The nanopatch shown in Fig. 2(a 1) exhibits the same properties along the x- and y-direction. The surface conductivity of the corresponding metasurface is isotropic in terms of the effective local medium approximation and experiences the resonance at the wavelength λ 0 = 790 nm. The imaginary part of the surface conductivity is positive and negative at the wavelengths smaller ( λ < λ 0) and larger ( λ > λ 0) than the resonant one [Fig. 2(c 1)]. After the nanoparticle is stretched, the reflectance and transmittance become different along the mutually orthogonal directions leading to the emergence of two new resonances at the wavelengths λ x and λ y [Fig. 2(b 2)]. When the electric field is directed along the x and y axes, the resonant wavelength redshifts and blueshifts, and the amplitude of the resonance is larger and smaller, respectively, than those for the isotropic case [Figs. 2(c 2) and 2(c 3)]. As a consequence, the hyperbolic regime characterized by the condition (2) emerges between the resonances in the wavelength range λ y < λ < λ x. We introduce the spectral difference between the resonances (SDBR) defined as the difference between the resonant wavelengths of the surface conductivity tensor components,
Δ λ = λ x λ y .
(5)
The SDBRs are Δ λ = 315 nm and Δ λ = 745 nm for the cases of f x = 0.7 and f x = 0.9, respectively [Figs. 2(c 2) and 2(c 3)]. For the extreme stretching, which means the transformation into the grating, the resonance along x-axis infinitely redshifts, signifying the hyperbolic regime for λ > λ y [Fig. 2(c 4)].

The corresponding electric field distributions at the resonances for all the cases from Figs. 2(a 1)2(a 4) are shown in Fig. 3. Figures 3(a 1) and 3(b 1) correspond to the isotropic case (square nanopatch) at the wavelength λ 0 = 790 nm and, as a result, are the same with respect to the rotation by 90 °. Figures 3(a 2) and 3(a 3) and 3(b 2)3(b 4) correspond to the wavelengths λ x (1000, and 1345 nm) and λ y (685, 600 and 580 nm), respectively. The field distribution for the non-resonant case, when the electric field is directed along the grating at λ = 1500 nm, is shown in Fig. 3(a 4).

FIG. 3.

The spatial distribution of the electric field amplitude within the unit cell at the distance of 20 nm from a top of nanopatches along (a 1)–(a 4) and across (b 1)–(b 4) the stretching direction ( x-axis) for the cases, shown in Figs. 2(a 1)2(a 4), respectively, at the corresponding resonant wavelengths: (a 1) and (b 1) 790 nm, (a 2) 1000 nm, (a 3) 1345 nm, (a 4) 1500 nm, (b 2) 685 nm, (b 3) 600 nm, and (b 4) 580 nm. The electric field amplitude changes from 0 (blue color) to the local maximum value E 0 (red color). The inset in (a 4) is shown for the same case, but for the color bar limits from 0 to E 0 / 3.

FIG. 3.

The spatial distribution of the electric field amplitude within the unit cell at the distance of 20 nm from a top of nanopatches along (a 1)–(a 4) and across (b 1)–(b 4) the stretching direction ( x-axis) for the cases, shown in Figs. 2(a 1)2(a 4), respectively, at the corresponding resonant wavelengths: (a 1) and (b 1) 790 nm, (a 2) 1000 nm, (a 3) 1345 nm, (a 4) 1500 nm, (b 2) 685 nm, (b 3) 600 nm, and (b 4) 580 nm. The electric field amplitude changes from 0 (blue color) to the local maximum value E 0 (red color). The inset in (a 4) is shown for the same case, but for the color bar limits from 0 to E 0 / 3.

Close modal

We analyze the spectral positions of the resonances and the SDBR during the stretching of the nanopatches up to the extreme near-grating case characterized by f x = 0.99. First, we set the value of the period p = 300 nm and the stretching factor η to 1 (isotropic case), 1.2, 1.5, and 2 and change the side of a square nanopatch a from 100 to 297 nm, 125 to 247.5 nm, from 100 to 198 nm, and from 75 to 148.5 nm, respectively, so that f x changes from 0.5 to 0.99 for all three cases [Figs. 4(a) and 4(d)]. Second, we set the value of the period p = 300 nm and the side of a square nanopatch a to 100, 150, and 180 nm, that correspond to the resonant wavelength λ 0 = 705 , 790, and 835 nm, and change the stretching factor η from 1.5 to 2.97, from 1 to 1.98, and from 1 to 1.65, respectively [Figs. 4(b) and 4(e)]. For the last case, a = 180 nm, p = 300 nm, the filling factor along the stretching direction f x changes from 0.6 to 0.99. Third, we fix the long side of a rectangular nanopatch as a x = 225 nm via three cases a = 180 nm and η = 1.25, a = 150 nm and η = 1.5, a = 125 nm and η = 1.8, and change the period p from 227.3 to 450 nm, so that f x changes from 0.5 to 0.99 for all three cases [Figs. 4(c) and 4(f)]. The latter three cases are compared with the isotropic case of a = a x = 150 nm with the period varying from 151.5 to 300 nm [gray line in Fig. 4(c)]. The corresponding values of the resonant wavelengths λ x, λ y and the SDBR Δ λ are shown in Figs. 4(a)4(c) and 4(d)4(f), respectively.

FIG. 4.

The dependencies of (a)–(c) the resonant wavelengths λ 0 (gray lines), λ x (upper solid lines), and λ y (lower dashed lines), and (d)–(f) the spectral difference between the resonances Δ λ on the filling factor along the stretching direction f x in the range from f x = 0.5 to f x = 0.99 (extreme anisotropic case) for the fixed (a) and (d) period ( p = 300 nm) and stretching factors [ η = 2 (green lines), η = 1.5 (red lines), and η = 1.2 (blue lines)], (b) and (e) period ( p = 300 nm) and sides of a square nanopatch [ a = 100 (green lines), a = 150 (red lines), and a = 180 nm (blue lines)], and (c) and (f) stretching factors and sides of a square nanopatch [ a = 125 nm, η = 1.8 (green lines), a = 150 nm, η = 1.5 (red lines), and a = 180 nm, η = 1.25 (blue lines)]. The blue regions correspond to the values of f x > 0.85. Here, the dots correspond to the numerical calculations, while the lines to their fitting with linear and quadratic (for λ 0 , λ x and Δ λ) curves for f x change from 0.5 to 0.85 and from 0.85 to 0.99, respectively.

FIG. 4.

The dependencies of (a)–(c) the resonant wavelengths λ 0 (gray lines), λ x (upper solid lines), and λ y (lower dashed lines), and (d)–(f) the spectral difference between the resonances Δ λ on the filling factor along the stretching direction f x in the range from f x = 0.5 to f x = 0.99 (extreme anisotropic case) for the fixed (a) and (d) period ( p = 300 nm) and stretching factors [ η = 2 (green lines), η = 1.5 (red lines), and η = 1.2 (blue lines)], (b) and (e) period ( p = 300 nm) and sides of a square nanopatch [ a = 100 (green lines), a = 150 (red lines), and a = 180 nm (blue lines)], and (c) and (f) stretching factors and sides of a square nanopatch [ a = 125 nm, η = 1.8 (green lines), a = 150 nm, η = 1.5 (red lines), and a = 180 nm, η = 1.25 (blue lines)]. The blue regions correspond to the values of f x > 0.85. Here, the dots correspond to the numerical calculations, while the lines to their fitting with linear and quadratic (for λ 0 , λ x and Δ λ) curves for f x change from 0.5 to 0.85 and from 0.85 to 0.99, respectively.

Close modal

The calculations with square nanopatches (isotropic case) are shown by the gray lines in Figs. 4(a) and 4(c). One can notice that resonant wavelength λ 0 increases linearly up to f x 0.85 and then quadratically for f x > 0.85 for both types of the induced anisotropy. The quadratic dependence of λ x ( f x ) in the range 0.85 < f x < 0.99 is associated with the strong mutual interaction between the nanoparticles leading to the substantial modification of the eigenmodes dispersion.30,38 The wavelength of the smaller resonance λ y slightly linearly redshifts [Fig. 4(a)] and blueshifts [Figs. 4(b) and 4(c)] for the cases of the fixed p and η, p and a, and a x, respectively, with the increase of f x. The wavelength of the larger resonance λ x increases linearly up to f x 0.85 and then quadratically for f x > 0.85 for all three cases [Figs. 4(a)4(c)] as well as the resonant wavelength λ 0 in the isotropic case. It leads to the significant change of the SDBR from 145, 270, and 415 nm at f x = 0.5 to 300, 520, and 720 nm at f x = 0.85, and then to 965, 1155, and 1620 nm at f x = 0.99 for the case of the fixed period p = 300 nm and the stretching factors η = 1.2 , 1.5, and 2, respectively [Fig. 4(d)]. For the case of the fixed period p = 300 nm and the side of a square nanopatch a = 150 nm, the SDBR changes from 0 nm at f x = 0.5 ( η = 1, isotropic case) to 610 nm at f x = 0.85 ( η = 1.7, anisotropic case) and then to 1620 nm at f x = 0.99 ( η = 1.98, extreme anisotropic case) as it is shown by the red line in Fig. 4(e). While for the case of the fixed long sides of a rectangular nanopatch a x = 225 nm and a y = 100 nm, the SDBR changes from 360 nm at f x = 0.5 ( p = 450 nm) to 480 nm at f x = 0.85 ( p = 264.7 nm) and then to 1420 nm at f x = 0.99 ( p = 227.3 nm) as it is shown by the red line in Fig. 4(f). The complete picture can be represented by the color map in Fig. 5 demonstrating the dependence of the SDBR on the stretching factor η and the side of a square nanopatch a at the period of p = 300 nm. The strong hyperbolicity is observed within the region 0.85 < f x < 0.99 and corresponds to the extremely anisotropic case. The transition from linear to quadratic dependence of the resonant wavelength and SDBR on the anisotropy degree is counter-intuitive and important to engineer properly the ultra-anisotropic systems and the plasmon canalization phenomena, where both resonances must be tuned simultaneously in extreme anisotropic case.

FIG. 5.

The dependence of the spectral difference between the resonances Δ λ shown by color on the side of the square nanopatch a and the stretching factor η. The black dashed lines mark the range of f x between 0.85 and 0.99. Here, the period of metasurface is equal to 300 nm.

FIG. 5.

The dependence of the spectral difference between the resonances Δ λ shown by color on the side of the square nanopatch a and the stretching factor η. The black dashed lines mark the range of f x between 0.85 and 0.99. Here, the period of metasurface is equal to 300 nm.

Close modal

Finally, we demonstrate the typical near-resonant spatial field distribution of the surface plasmon-polariton excited by the vertical ( z-oriented) dipole and propagating in the canalization regime along the plasmonic nanopatch-based metasurface (Fig. 6). We take the metasurface consisting of 15 × 15 nanopatches taken from Fig. 2(a 2), i.e., p = 300 nm, a = 150 nm, η = 1.4. The canalization regime is observed when one of the surface conductivity components is large, while the orthogonal one is near-zero. Therefore, we draw the spectral dependence of Im ( σ x ) / Im ( σ y ) characterizing the efficiency of plasmon canalization in Fig. 6(a) by the blue line. One can notice the small peak in the vicinity of λ y = 700 nm and the pronounced peak in the vicinity of λ x = 1000 nm. For the small peak, both components of the surface conductivity are near-zero meaning inefficient propagation with a divergence, while for the large peak, σ y is near-zero and σ x is resonant leading to the divergenceless propagation of the surface wave. Nevertheless, the most efficient plasmon canalization takes place for the wavelength slightly shifted from the resonance that is explained by the high absorption losses in the vicinity of the resonance characterized by the large real part of the surface conductivity component along the canalization direction, as it is shown by the red dashed line in Fig. 6(a). As a consequence, we choose the wavelength λ = 1100 nm for the optimum plasmon canalization for this specific metasurface design.

FIG. 6.

(a) The dependence of the real part of σ x and the ratio of the imaginary parts [ Im ( σ x ) / Im ( σ y )] on the wavelength. (b) and (d) The spatial in-plane distribution of the (b) amplitude and (d) real part of the normal component of the magnetic field for the surface plasmon-polariton excited by the vertical magnetic dipole located in the center of the structure (marked by the green point) at the wavelength λ = 1100 nm marked by the orange star and dashed line in (a). The metasurface consists of the nanopatches shown in Fig. 2(a 2) with a period of p = 300 nm, and the total size of a metasurface is 15 × 15 periods ( 4.5 μ m × 4.5 μ m). (c) The Fourier spectrum of the magnetic field spatial distribution shown in (b) and (d) within the first Brillouin zone. The maximum values correspond to the isofrequency contour at the wavelength λ = 1100 nm. The yellow circle corresponds to the light cone in the air. The blue lines correspond to the isofrequency contour, calculated analytically using the surface conductivity from Fig. 2(c 2).

FIG. 6.

(a) The dependence of the real part of σ x and the ratio of the imaginary parts [ Im ( σ x ) / Im ( σ y )] on the wavelength. (b) and (d) The spatial in-plane distribution of the (b) amplitude and (d) real part of the normal component of the magnetic field for the surface plasmon-polariton excited by the vertical magnetic dipole located in the center of the structure (marked by the green point) at the wavelength λ = 1100 nm marked by the orange star and dashed line in (a). The metasurface consists of the nanopatches shown in Fig. 2(a 2) with a period of p = 300 nm, and the total size of a metasurface is 15 × 15 periods ( 4.5 μ m × 4.5 μ m). (c) The Fourier spectrum of the magnetic field spatial distribution shown in (b) and (d) within the first Brillouin zone. The maximum values correspond to the isofrequency contour at the wavelength λ = 1100 nm. The yellow circle corresponds to the light cone in the air. The blue lines correspond to the isofrequency contour, calculated analytically using the surface conductivity from Fig. 2(c 2).

Close modal

The spatial distributions of the normal component of the magnetic field H z are shown in Figs. 6(b) and 6(d) explicitly demonstrating the collinear divergenceless propagation of surface plasmon-polariton with a small decay along the x-axis. The canalization regime is characterized by the flat isofrequency contour, which also indicates the topological transition between the elliptical and hyperbolic shapes of the isofrequency contours. In order, to extract the isofrequency contour, we apply the two-dimensional Fourier transform to the calculated distribution of the H z-field [Fig. 6(c)] following the standard procedure.32,39 On the other hand, we calculated the corresponding isofrequency contour analytically [blue dashed line in Fig. 6(c)] using the dispersion equation for hyperbolic plasmon-polaritons at metasurfaces10,11 with the retrieved surface conductivity tensor, shown in Fig. 2(c 2). We claim a good coincidence between the isofrequency contours calculated analytically with the effective surface conductivity and extracted from the full-wave numerical simulations for the full-scale metasurface consisting of 15 × 15 nanopatches.

The experimental implementation of the surface waves excitation in the hyperbolic and canalization regimes includes the loop and probe near-field technique in the microwave range,32 the prism coupling,15,30 the defect-induced structures,14,40 the metallic41 and dielectric nanoantennas,42,43 the metal tip of the microscope,33,41,44 and the edge-oriented launching44–46 in the visible and infrared ranges.

We have investigated a number of plasmonic metasurface designs based on gold nanopatches by varying the side of the nanoparticle, its stretching (anisotropy) degree and the period of a metasurface from the isotropic (square nanopatches) to extreme anisotropic (gratings) cases. We define the quadratic dependence of the spectral bandwidth for one of the resonances on the anisotropy degree when the electric field is oriented along the stretching direction. It leads to the broadband spectral bandwidth between the resonances characterizing the hyperbolic regime of a metasurface. As a possible application verifying the eligibility of the discussed results, we demonstrated the plasmon canalization in the vicinity of the engineered resonance of the plasmonic metasurface.

The results obtained may be extended to different meta-atom geometries, lattice types, and materials (e.g., silver, aluminum). We believe that these results may be used by optical researchers and engineers for the efficient designing, engineering, and optimization of HMS plasmon resonances for a wide range of applications.

The numerical calculations have been conducted with the support of the U.S. Civilian Research and Development Foundation (CRDF Global) (No. G-202401-71609). A.H. thanks the IEEE Microwave Theory and Technology Society (IEEE MTT-S) for the support within the Undergraduate/Pre-Graduate Scholarship. S.P. and O.D. acknowledge the EURIZON Fellowship Programme “Remote Research Grants EU#3003-EURIZON” (Agreement No. 871072). O.Y. acknowledges the support of the Alexander von Humboldt Foundation within the framework of the Humboldt Research Fellowship for postdoctoral researchers.

The authors have no conflicts to disclose.

A. Hrinchenko: Data curation (lead); Formal analysis (lead); Investigation (lead); Software (lead); Validation (lead); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal). S. Polevoy: Formal analysis (equal); Software (equal); Visualization (equal); Writing – review & editing (equal). O. Demianyk: Software (equal); Visualization (equal); Writing – review & editing (equal). O. Yermakov: Conceptualization (equal); Formal analysis (equal); Methodology (lead); Supervision (lead); Visualization (equal); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

1.
A. K.
Geim
and
I. V.
Grigorieva
, “
Van der Waals heterostructures
,”
Nature
499
,
419
425
(
2013
).
2.
F.
Xia
,
H.
Wang
,
D.
Xiao
,
M.
Dubey
, and
A.
Ramasubramaniam
, “
Two-dimensional material nanophotonics
,”
Nat. Photonics
8
,
899
907
(
2014
).
3.
H.-T.
Chen
,
A. J.
Taylor
, and
N.
Yu
, “
A review of metasurfaces: Physics and applications
,”
Rep. Prog. Phys.
79
,
076401
(
2016
).
4.
H.-H.
Hsiao
,
C. H.
Chu
, and
D. P.
Tsai
, “
Fundamentals and applications of metasurfaces
,”
Small Methods
1
,
1600064
(
2017
).
5.
C.-W.
Qiu
,
T.
Zhang
,
G.
Hu
, and
Y.
Kivshar
, “
Quo vadis, metasurfaces?
,”
Nano Lett.
21
,
5461
5474
(
2021
).
6.
Z.
Guo
,
H.
Jiang
, and
H.
Chen
, “
Hyperbolic metamaterials: From dispersion manipulation to applications
,”
J. Appl. Phys.
127
,
071101
(
2020
).
7.
H. N.
Krishnamoorthy
,
Z.
Jacob
,
E.
Narimanov
,
I.
Kretzschmar
, and
V. M.
Menon
, “
Topological transitions in metamaterials
,”
Science
336
,
205
209
(
2012
).
8.
A.
Poddubny
,
I.
Iorsh
,
P.
Belov
, and
Y.
Kivshar
, “
Hyperbolic metamaterials
,”
Nat. Photonics
7
,
948
957
(
2013
).
9.
L.
Ferrari
,
C.
Wu
,
D.
Lepage
,
X.
Zhang
, and
Z.
Liu
, “
Hyperbolic metamaterials and their applications
,”
Prog. Quantum Electron.
40
,
1
40
(
2015
).
10.
O.
Yermakov
,
A.
Ovcharenko
,
M.
Song
,
A.
Bogdanov
,
I.
Iorsh
, and
Y. S.
Kivshar
, “
Hybrid waves localized at hyperbolic metasurfaces
,”
Phys. Rev. B
91
,
235423
(
2015
).
11.
J. S.
Gomez-Diaz
,
M.
Tymchenko
, and
A.
Alu
, “
Hyperbolic plasmons and topological transitions over uniaxial metasurfaces
,”
Phys. Rev. Lett.
114
,
233901
(
2015
).
12.
J.
Gomez-Diaz
and
A.
Alu
, “
Flatland optics with hyperbolic metasurfaces
,”
ACS Photonics
3
,
2211
2224
(
2016
).
13.
O.
Takayama
and
A. V.
Lavrinenko
, “
Optics with hyperbolic materials
,”
J. Opt. Soc. Am. B
36
,
F38
F48
(
2019
).
14.
A. A.
High
,
R. C.
Devlin
,
A.
Dibos
,
M.
Polking
,
D. S.
Wild
,
J.
Perczel
,
N. P.
De Leon
,
M. D.
Lukin
, and
H.
Park
, “
Visible-frequency hyperbolic metasurface
,”
Nature
522
,
192
196
(
2015
).
15.
A.
Samusev
,
I.
Mukhin
,
R.
Malureanu
,
O.
Takayama
,
D. V.
Permyakov
,
I. S.
Sinev
,
D.
Baranov
,
O.
Yermakov
,
I. V.
Iorsh
,
A. A.
Bogdanov
, and
A. V.
Lavrinenko
, “
Polarization-resolved characterization of plasmon waves supported by an anisotropic metasurface
,”
Opt. Express
25
,
32631
32639
(
2017
).
16.
V. G.
Kravets
,
A. V.
Kabashin
,
W. L.
Barnes
, and
A. N.
Grigorenko
, “
Plasmonic surface lattice resonances: A review of properties and applications
,”
Chem. Rev.
118
,
5912
5951
(
2018
).
17.
A.
Hrinchenko
and
O.
Yermakov
, “
Designing optical hyperbolic metasurfaces based on gold nanodisks
,”
J. Phys. D: Appl. Phys.
56
,
465105
(
2023
).
18.
S. A. H.
Gangaraj
,
G. W.
Hanson
,
M. G.
Silveirinha
,
K.
Shastri
,
M.
Antezza
, and
F.
Monticone
, “
Unidirectional and diffractionless surface plasmon polaritons on three-dimensional nonreciprocal plasmonic platforms
,”
Phys. Rev. B
99
,
245414
(
2019
).
19.
A.
Nemilentsau
,
T.
Stauber
,
G.
Gómez-Santos
,
M.
Luskin
, and
T.
Low
, “
Switchable and unidirectional plasmonic beacons in hyperbolic two-dimensional materials
,”
Phys. Rev. B
99
,
201405
(
2019
).
20.
G.
Hu
,
A.
Krasnok
,
Y.
Mazor
,
C.-W.
Qiu
, and
A.
Alù
, “
Moiré hyperbolic metasurfaces
,”
Nano Lett.
20
,
3217
3224
(
2020
).
21.
G.
Hu
,
Q.
Ou
,
G.
Si
,
Y.
Wu
,
J.
Wu
,
Z.
Dai
,
A.
Krasnok
,
Y.
Mazor
,
Q.
Zhang
,
Q.
Bao
,
C.-W.
Qiu
, and
A.
Alu
, “
Topological polaritons and photonic magic angles in twisted α-MoO 3 bilayers
,”
Nature
582
,
209
213
(
2020
).
22.
A.
Girich
,
L.
Ivzhenko
,
A.
Hrinchenko
,
S.
Tarapov
, and
O.
Yermakov
, “
Manipulation over surface waves in bilayer hyperbolic metasurfaces: Topological transition and multidirectional canalization
,”
IEEE Microw. Wireless Technol. Lett.
33
,
367
370
(
2023
).
23.
P. V.
Kapitanova
,
P.
Ginzburg
,
F. J.
Rodríguez-Fortuño
,
D. S.
Filonov
,
P. M.
Voroshilov
,
P. A.
Belov
,
A. N.
Poddubny
,
Y. S.
Kivshar
,
G. A.
Wurtz
, and
A. V.
Zayats
, “
Photonic spin Hall effect in hyperbolic metamaterials for polarization-controlled routing of subwavelength modes
,”
Nat. Commun.
5
,
3226
(
2014
).
24.
O. Y.
Yermakov
,
A. I.
Ovcharenko
,
A. A.
Bogdanov
,
I. V.
Iorsh
,
K. Y.
Bliokh
, and
Y. S.
Kivshar
, “
Spin control of light with hyperbolic metasurfaces
,”
Phys. Rev. B
94
,
075446
(
2016
).
25.
Y.
Mazor
and
A.
Alù
, “
Routing optical spin and pseudospin with metasurfaces
,”
Phys. Rev. Appl.
14
,
014029
(
2020
).
26.
Y.
Liang
,
K.
Koshelev
,
F.
Zhang
,
H.
Lin
,
S.
Lin
,
J.
Wu
,
B.
Jia
, and
Y.
Kivshar
, “
Bound states in the continuum in anisotropic plasmonic metasurfaces
,”
Nano Lett.
20
,
6351
6356
(
2020
).
27.
Y.
Liang
,
H.
Lin
,
S.
Lin
,
J.
Wu
,
W.
Li
,
F.
Meng
,
Y.
Yang
,
X.
Huang
,
B.
Jia
, and
Y.
Kivshar
, “
Hybrid anisotropic plasmonic metasurfaces with multiple resonances of focused light beams
,”
Nano Lett.
21
,
8917
8923
(
2021
).
28.
O.
Kotov
and
Y. E.
Lozovik
, “
Enhanced optical activity in hyperbolic metasurfaces
,”
Phys. Rev. B
96
,
235403
(
2017
).
29.
G.
Palermo
,
K. V.
Sreekanth
,
N.
Maccaferri
,
G. E.
Lio
,
G.
Nicoletta
,
F.
De Angelis
,
M.
Hinczewski
, and
G.
Strangi
, “
Hyperbolic dispersion metasurfaces for molecular biosensing
,”
Nanophotonics
10
,
295
314
(
2020
).
30.
O. Y.
Yermakov
,
D. V.
Permyakov
,
F. V.
Porubaev
,
P. A.
Dmitriev
,
A. K.
Samusev
,
I. V.
Iorsh
,
R.
Malureanu
,
A. V.
Lavrinenko
, and
A. A.
Bogdanov
, “
Effective surface conductivity of optical hyperbolic metasurfaces: From far-field characterization to surface wave analysis
,”
Sci. Rep.
8
,
14135
(
2018
).
31.
D.
Correas-Serrano
,
A.
Alù
, and
J. S.
Gomez-Diaz
, “
Plasmon canalization and tunneling over anisotropic metasurfaces
,”
Phys. Rev. B
96
,
075436
(
2017
).
32.
O.
Yermakov
,
V.
Lenets
,
A.
Sayanskiy
,
J.
Baena
,
E.
Martini
,
S.
Glybovski
, and
S.
Maci
, “
Surface waves on self-complementary metasurfaces: All-frequency hyperbolicity, extreme canalization, and TE-TM polarization degeneracy
,”
Phys. Rev. X
11
,
031038
(
2021
).
33.
P.
Li
,
G.
Hu
,
I.
Dolado
,
M.
Tymchenko
,
C.-W.
Qiu
,
F. J.
Alfaro-Mozaz
,
F.
Casanova
,
L. E.
Hueso
,
S.
Liu
,
J. H.
Edgar
,
S.
Velez
,
A.
Alu
, and
R.
Hillenbrand
, “
Collective near-field coupling and nonlocal phenomena in infrared-phononic metasurfaces for nano-light canalization
,”
Nat. Commun.
11
,
3663
(
2020
).
34.
P.-H.
Chang
,
C.
Lin
, and
A. S.
Helmy
, “
Field canalization using anisotropic 2D plasmonics
,”
npj 2D Mater. Appl.
6
,
5
(
2022
).
35.
Q.
Xing
,
J.
Zhang
,
Y.
Fang
,
C.
Song
,
T.
Zhao
,
Y.
Mou
,
C.
Wang
,
J.
Ma
,
Y.
Xie
,
S.
Huang
,
M.
Mu
,
Y.
Lei
,
W.
Shi
,
F.
Huang
, and
H.
Yan
, “
Tunable anisotropic van der Waals films of 2M-WS2 for plasmon canalization
,”
Nat. Commun.
15
,
2623
(
2024
).
36.
D. I.
Yakubovsky
,
A. V.
Arsenin
,
Y. V.
Stebunov
,
D. Y.
Fedyanin
, and
V. S.
Volkov
, “
Optical constants and structural properties of thin gold films
,”
Opt. Express
25
,
25574
25587
(
2017
).
37.
M.
Born
and
E.
Wolf
,
Principles of Optics: Electromagnetic Theory of Propagation, Interference and Diffraction of Light
(
Elsevier
,
2013
).
38.
M. A.
Gorlach
and
P. A.
Belov
, “
Effect of spatial dispersion on the topological transition in metamaterials
,”
Phys. Rev. B
90
,
115136
(
2014
).
39.
Y.
Yang
,
L.
Jing
,
L.
Shen
,
Z.
Wang
,
B.
Zheng
,
H.
Wang
,
E.
Li
,
N.-H.
Shen
,
T.
Koschny
,
C. M.
Soukoulis
, and
H.
Chen
, “
Hyperbolic spoof plasmonic metasurfaces
,”
npg Asia Mater.
9
,
e428
(
2017
).
40.
J.
Lin
,
J. B.
Mueller
,
Q.
Wang
,
G.
Yuan
,
N.
Antoniou
,
X.-C.
Yuan
, and
F.
Capasso
, “
Polarization-controlled tunable directional coupling of surface plasmon polaritons
,”
Science
340
,
331
334
(
2013
).
41.
P.
Li
,
I.
Dolado
,
F. J.
Alfaro-Mozaz
,
F.
Casanova
,
L. E.
Hueso
,
S.
Liu
,
J. H.
Edgar
,
A. Y.
Nikitin
,
S.
Vélez
, and
R.
Hillenbrand
, “
Infrared hyperbolic metasurface based on nanostructured van der Waals materials
,”
Science
359
,
892
896
(
2018
).
42.
I. S.
Sinev
,
A. A.
Bogdanov
,
F. E.
Komissarenko
,
K. S.
Frizyuk
,
M. I.
Petrov
,
I. S.
Mukhin
,
S. V.
Makarov
,
A. K.
Samusev
,
A. V.
Lavrinenko
, and
I. V.
Iorsh
, “
Chirality driven by magnetic dipole response for demultiplexing of surface waves
,”
Laser Photonics Rev.
11
,
1700168
(
2017
).
43.
F.
Walla
,
F.
Bürkle
,
I.
Sinev
,
M. M.
Wiecha
,
N.
Mecklenbeck
,
K.
Ladutenko
,
R.
Malureanu
,
F.
Komissarenko
,
A.
Lavrinenko
,
A.
Bogdanov
,
A.
Soltani
, and
H. G.
Roskos
, “
Near-field observation of guided-mode resonances on a metasurface via dielectric nanosphere excitation
,”
ACS Photonics
5
,
4238
4243
(
2018
).
44.
S.
Dai
,
Q.
Ma
,
Y.
Yang
,
J.
Rosenfeld
,
M. D.
Goldflam
,
A.
McLeod
,
Z.
Sun
,
T. I.
Andersen
,
Z.
Fei
,
M.
Liu
, and
Y.
Shao
, “
Efficiency of launching highly confined polaritons by infrared light incident on a hyperbolic material
,”
Nano Lett.
17
,
5285
5290
(
2017
).
45.
W.
Ma
,
P.
Alonso-González
,
S.
Li
,
A. Y.
Nikitin
,
J.
Yuan
,
J.
Martín-Sánchez
,
J.
Taboada-Gutiérrez
,
I.
Amenabar
,
P.
Li
,
S.
Vélez
, and
C.
Tollan
, “
In-plane anisotropic and ultra-low-loss polaritons in a natural van der Waals crystal
,”
Nature
562
,
557
562
(
2018
).
46.
Z.
Dai
,
G.
Hu
,
G.
Si
,
Q.
Ou
,
Q.
Zhang
,
S.
Balendhran
,
F.
Rahman
,
B. Y.
Zhang
,
J. Z.
Ou
,
G.
Li
,
A.
Alù
,
C.-W.
Qiu
, and
Q.
Bao
, “
Edge-oriented and steerable hyperbolic polaritons in anisotropic van der Waals nanocavities
,”
Nat. Commun.
11
,
6086
(
2020
).