In contrast to conventional lenses, which have a constant refractive index, gradient-index (GRIN) lenses have a refractive index that changes throughout their volume^[1]. The distribution of these GRIN lenses provides unique optical properties, such as reduced aberrations, light weight, and compactness. Examples of these GRIN lenses include the Luneburg lens^[2] and the Gutman lens^[3]. Each of these lenses features a radially varying refractive index profile and is designed to be aberration-free. The Luneburg lens focuses planar waves to a point on its opposite surface. In contrast, the Gutman lens can dynamically adjust its focal distance from the center to the surface, allowing it to focus planar waves within itself. Thus, while the Luneburg lens consistently focuses waves to a specific point on its surface, the Gutman lens offers variable focal lengths, making it a generalized form of the Luneburg lens. Therefore, the Luneburg lens can be regarded as a special case of the Gutman lens when the focal distance is at the surface of the lens.

In 2012, Dyachenko et. al. introduced a discrete design of a two-dimensional photonics Luneburg lens engineered using rods of different radii for producing an effective GRIN distribution. In that work, a new dodecagonal quasicrystal array was introduced, and its efficiency was compared against a square array^[4]. Two-dimensional GRIN Luneburg lenses have garnered significant importance in the field of plasmonics^[5–7]. In this case, wave propagation occurs through surface plasmon polaritons (SPPs), which are collective oscillations of electrons on the surface of a metal^[8]. As a result, these lenses have been designed as metasurfaces, which allow for the efficient control and manipulation of the two-dimensional nature of SPP waves^[9].

An effective GRIN distribution in metasurfaces can be engineered through the creation of an array of nanoscale-sized holes, each with a distinct size, patterned within a thin dielectric film. The proportion of hole sizes to dielectric film surroundings leads to specific GRIN distributions for a given wavelength^[10,11]. The geometrical array of holes often involves organizing the nanoscale-sized holes in a rectangular pattern with a designated periodicity^[6,7]. This can allow for the design of different kinds of GRIN lenses like the well-studied Luneburg lenses and the modified Eaton lenses^[5]. We note that the Luneburg lens is a particular case more generally known as a Gutman GRIN lens. Finding the geometry for this lens in a metasurface is a challenging problem.

In this paper, we present for the first time, to the best of our knowledge, a plasmonic metasurface Gutman lens (PMGL) using a carefully designed arrangement of nanometer-sized holes in thin dielectric films. We present a detailed comparison of two distinct designs for a PMGL optimized for a wavelength of 950 nm. The first design is a “12-fold symmetry design,” constructed based on a dodecagonal quasicrystal array, while the second utilizes a periodic rectangular array of air holes that will be henceforth called the “square lattice design.” To evaluate the performance of each design, we have created separate metasurfaces with varying focal lengths, thoroughly examining their focal properties. A comprehensive comparative analysis between the two designs has been conducted. Despite the theoretical nature of this work, all calculations are meticulously executed with a focus on practical experimental details that should be easily reproduced empirically.

The design of the PMGL incorporates three distinct media: gold (Au), characterized by a refractive index of $n_{\mathrm{Au}}=0.20152+6.0961i$^[12]; polymethyl methacrylate (PMMA), with a refractive index of $n_2=1.4846$^[13]; and air, with a refractive index of $n_1=1$, at a wavelength of $\lambda =950\mathrm{nm}$. In this configuration, two effective refractive indices can be identified for the SPPs. The first effective refractive index, $n_{\mathrm{eff},1}$, corresponds to the region devoid of PMMA, encompassing solely the Au-air interface, referred to as Region I. The second effective refractive index, $n_{\mathrm{eff},2}$, characterizes the region comprising air, PMMA, and Au, labeled as Region II. This design is depicted in Fig. 1(a). Utilizing the established dispersion relations of the SPPs^[8], we can derive the effective refractive index for Region I as follows: $$n_{\mathrm{eff},1}=\sqrt{\frac{n_1^2{n}_{\mathrm{Au}}^2}{n_1^2+n_{\mathrm{Au}}^2}},$$(1)obtaining $n_{\mathrm{eff},1}=1.013$.

On the other hand, the effective refractive index for Region II can be calculated using $$n_{\mathrm{eff},2}=\frac{{\beta }_2}{k_0},$$(2)where $k_0=\frac{2\pi }{\lambda }$ is the free space propagation constant and ${\beta }_2$ is the corresponding one to the SPPs^[5]. To determine the propagation constant ${\beta }_2$, it is necessary to solve the following transcendental equation: $$\tanh (k_2{\epsilon }_{2}d)=-\left(\frac{k_2{k}_1+k_{\mathrm{Au}}{k}_2}{k_2^2+k_{\mathrm{Au}}{k}_1}\right),$$(3)where $k_i=\frac{\sqrt{{\beta }_2^2-{\epsilon }_i{\omega }^2/c^2}}{{\epsilon }_i}$ is the wavenumber, $\epsilon$ is the permittivity, $\omega$ is the frequency of the SPPs, $c$ is the speed of light in vacuum, and $d$ is the PMMA thickness. In our design, $d=390\mathrm{nm}$ and $i$ represents each medium. Applying these equations, we determine that $n_{\mathrm{eff},2}=1.495$.

In order to create a metasurface with controlled refractive index variation through an array of nanoscale-sized holes, the Maxwell–Garnett effective medium theory has been employed^[14]. This theory relies on the principles of composite materials, enabling manipulation of the effective refractive index by adjusting the volume fraction of each component. In this scenario, introducing nanoscale air holes onto the PMMA slab alters the value of $n_{\mathrm{eff},2}$ within a specified region. The resultant effective refractive index variation is determined by $$n_{\mathrm{eff}}^2=n_{\mathrm{eff},2}^2(1-f)+n_{\mathrm{eff},1}^{2}f,$$(4)where $f$ is the volume fraction of the air holes. For a 12-fold symmetry with a fixed hole radius $r$, the value of $f=\frac{N\pi r^2}{S}$, where $S$ represents the area of the metasurface and $N$ denotes the number of the holes present. For a rectangular pattern with a specific periodicity, i.e., for a square lattice, the value of $f$ is calculated with $f=\frac{\pi r^2}{l^2}$, where $l$ is the lattice period.

If the desired value of $n_{\mathrm{eff}}$ is known, it is possible to determine the corresponding radius of the hole in a plasmonic metasurface using Eq. (4). For the case of a 12-fold symmetry, $$r=\sqrt{\frac{S}{N\pi }\left(\frac{n_{\mathrm{eff}}^2-n_{\mathrm{eff},2}^2}{1-n_{\mathrm{eff},2}^2}\right)}.$$(5)

For the case of a square lattice, the radius is given by $$r=l\sqrt{\frac{n_{\mathrm{eff}}^2-n_{\mathrm{eff},2}^2}{\pi (1-n_{\mathrm{eff},2}^2)}}.$$(6)

These equations enable the calculation of the area proportion of the metasurface occupied by the holes, providing valuable insights into the spatial distribution and arrangement of the structures.

Experimentally, a GRIN metasurface with approximated radial symmetry is constructed using a multilayer system consisting of glass, gold, PMMA, and air. The fabrication of the PMGL involves precise steps. First, a 70-nm-thick gold film is deposited onto a glass substrate using a physical vapor deposition technique such as electron-beam evaporation, sputtering, or thermal evaporation. Next, a 390-nm-thick PMMA layer is uniformly deposited onto the gold film through spin coating, as shown in Fig. 1(a). Once the desired multilayer structure is obtained, nanometer-sized cylindrical holes and the delimiting lens area are etched into the PMMA film using electron-beam lithography (EBL), which typically has a resolution limit of approximately 50 nm, as shown in Fig. 1(b).

The effective refractive index of the PMGL is sensitive to both material and geometric parameters, including the refractive indices of gold, PMMA, and air, as well as the nanoscale hole configuration in the dielectric film^[7]. Errors in the effective refractive index can degrade the lens performance by affecting focal properties, causing discrepancies in the focal length and position, and introducing aberrations and irregularities in field intensity^[5].

A proposed experimental process for SPP excitation and propagation along plasmonic metasurfaces involves using a dielectric grating set on top of the gold thin film^[6]. The grating’s period matches the SPP wavelength along the gold/air interface, efficiently coupling light into SPPs. Positioned a few microns apart from the PMGL, this grating excites and propagates SPPs into the lens. Techniques such as leakage-radiation microscopy can be used to monitor and characterize the SPP propagation and interaction with the PMGL^[5].

The theory developed here provides a way to calculate the radii at each array location for any known GRIN distribution. We are interested in the Gutman lens whose refractive index distribution is given by^[3]$$n_{\mathrm{eff}}^2=\frac{R_L^2+F^2-{\rho }^2}{F^2},$$(7)where ${\rho }^2=x^2+z^2$ is the radial distance in Cartesian coordinates, $R_L$ is the radius of the lens, and $F$ is the focal length. Notice that if $F\le R_L$, the Gutman lens functions as an aberration-free lens^[3], which is our case of interest in this work, but if $F>R_L$, the lens introduces spherical aberration^[15]. A special case occurs when $F$ is equal to $R_L$, transforming the Gutman lens into a Luneburg lens with a GRIN distribution of $n=\sqrt{2-(\rho /R_L)}$^[2].

Equation (7) indicates that the refractive index decreases from its maximum value at the center, $n^2(\rho =0)=\frac{R_L^2+F^2}{F^2}$, to a minimum value at the lens surface, $n^2(\rho =R_L)=1$. The minimum value in the GRIN distribution is always 1, regardless of the value of $F$. However, the maximum value of the GRIN occurs at the center of the lens. The variation of this maximum value for different values of $F$ is shown in Fig. 2(a). Notice that the maximum value increases as the value of $F$ decreases.

The value of the effective refractive index in the fabrication of metasurfaces is subject to two limiting factors: the effective refractive index of Region II ($n_{\mathrm{eff},2}=1.495$ for $\lambda =950\mathrm{nm}$) and the resolution limit of the EBL ($\sim 50\mathrm{nm}$). As a result of these limitations, we can only fabricate PMGLs with $F\ge 0.90801R_L$ for the 12-fold symmetry case and $F\ge 0.90959R_L$ for the square lattice configuration. The provided values of $F$ ensure that there are no cylinders with radii smaller than the resolution limit of the EBL, as shown in Fig. 2(b).

Figure 3 shows the dodecagonal and square metasurface arrays with $F\ge 0.90801R_L$ and $F\ge 0.90959R_L$, respectively. For our study, we take $R_L=5\mathrm{\mu m}$. It is evident that for both arrangements, starting from $\rho =4.35\mathrm{\mu m}$ (depicted by the red circle), there is a noticeable overlap among the holes. To address this concern, we remove these overlapping holes, resulting in a refined lens with a maximum radius of 4.35 µm. Nevertheless, the value of $R_L$ will remain constant at 5 µm. We will observe that the induced aberrations are negligible for this small distance.

We have created, for both designs, separate metasurfaces for different values F, namely $0.90801R_L$ (only for the 12-fold symmetry design), $0.90959R_L$ (only for the square lattice design), $0.93R_L$, $0.95R_L$, $0.97R_L$, and $R_L$ (for both designs), to evaluate their performance.

The simulations were carried out using the commercial platform COMSOL Multiphysics 5.1, employing the frequency-domain solver. The simulation area consisted of a symmetrical box with a size of $21.05\lambda$ for each side. The plasmonic metasurface lenses were placed at the center of the simulation domain. An incident plane wave, with an out-of-plane magnetic field, was propagated from the left boundary to the right, as depicted in Fig. 4. To avoid artificial reflections, scattering boundary conditions were applied on all sides of the simulation box. A triangular mesh was implemented with a minimum triangle size of $0.0016\lambda$ and a maximum size of $0.052\lambda$.

The first and second columns of Fig. 4 exhibit the field intensity of SPPs through metasurfaces with a 12-fold symmetry design. The third and fourth columns display the analogous evolution through metasurfaces with a square lattice design. The field intensity along the $z$-axis is illustrated in Figs. 5(a) and 5(b) for metasurfaces with a 12-fold symmetry design and a square lattice design, respectively.

We have found by a detailed comparison of the field intensities that the two designs show the following: 1.The square lattice design shows a neater propagation than that generated from the 12-fold symmetry design, as can be seen in Fig. 4, where a large amount of scattering is observed on the 12-fold symmetry, which disturbs the propagation of SPPs. This result seems to contrast with the results by Dyachenko et. al. for a photonic Luneburg lens^[4]. They investigated the photonic lens in a normalized frequency range of 0.34–0.38 ($l/\lambda$). In our plasmonic scenario, with a hole separation value $l$ of 300 nm and a $\lambda$ of 950 nm, the resulting normalized frequency is approximately 0.32. This value falls outside the operational range specified for the 12-fold symmetry design. Also, they investigated a Luneburg lens, and in our case, it is a Gutman lens. However, our results are consistent with frequencies below 0.33.2.Figure 5(b) shows how the square lattice reproduces the focusing properties of the Gutman lens, i.e., for each selected value of $F$, a different focal point is generated in an orderly manner. In contrast, this variability is not observed for the 12-fold symmetry design, as depicted in Fig. 5(a).3.Figure 5(c) compares the field intensities of Gutman lenses with $F=R_L$: a continuous lens (the black curve), one with a square lattice design (the magenta curve), and another with a 12-fold symmetry design (the blue curve). The field intensity along the $z$-axis of a metasurface with a square lattice more closely resembles that of a continuous Luneburg lens, while it differs significantly for the 12-fold symmetry design. Note that field intensities for both designs are calculated on the $x–z$ plane, as shown in Figs. 5(e) and 5(f). The same holds true for the location of the point with the highest peak intensity, considered as the focal point. It can be observed, in Fig. 5(c), that the focal point of a metasurface with a square lattice is $0.9382R_L$ (the gray dot), which closely approximates the focal point of a continuous Luneburg lens at $R_L$ (the green dot), while it differs significantly for the 12-fold symmetry design (the yellow dot), localized at $0.5545R_L$.4.The focal points generated by metasurfaces with a square lattice for each value of $F$ are located at offsets from the specified $F$ values: $F_{0.90959R_L}=0.8813R_L$, $F_{0.93R_L}=0.8885R_L$, $F_{95R_L}=0.9011R_L$, $F_{0.97R_L}=0.9219R_L$, and $F_{R_L}=0.9382R_L$. These offsets are due to the discretization of the lens and the transition from a radius of 5 to 4.35 µm due to hole interceptions. For example, from Fig. 5(d), the field intensity of a continuous Gutman lens with $F=R_L$ cut at $\rho =4.35\mathrm{\mu m}$ generates a focal point at $0.9434R_L$. In comparison, the focal point generated by a Gutman lens with a square lattice design, where $F=R_L$, is located at $0.9382R_L$. The difference between these focal points is $0.0052R_L$, equivalent to 25.9 nm. This difference is negligible compared to the wavelength.5.The full width at half-maximum (FWHM) was measured at the focal plane for each Gutman lens with the square lattice design. The resulting FWHM values were $0.682\lambda$, $0.711\lambda$, $0.713\lambda$, $0.709\lambda$, and $0.725\lambda$ for F values of $0.90959R_L$, $0.93R_L$, $0.95R_L$, $0.97R_L$, and $R_L$, respectively. To determine the quality of the focus, a propagation was generated in a homogeneous lens, which is without GRIN [see Fig. 1(a)], yielding an FWHM value of $1.83\lambda$ [see Figs. 5(g) and 5(h)]. In other non-plasmonic works, it has been demonstrated that effective medium theory can be utilized to achieve super-resolution focusing or imaging^[16,17]. We can observe that our plasmonic work aligns with these studies, as we achieve super-resolution focusing using effective medium theory. The resulting FWHM values from the 12-fold symmetry design are presented in Table 1. While in some instances, the FWHM is smaller in the 12-fold symmetry design compared to that in the square lattice design, this is an isolated focal point that does not meet the characteristics of a Gutman lens. Additionally, in Table 1, the focal points are provided to facilitate a clearer comparison of the field intensities.

In conclusion, we have designed two plasmonic metasurface Gutman lenses using the square lattice and the 12-fold symmetry designs. Upon comparison, evidence shows that the square lattice design exhibits smoother propagation and more closely resembles the continuous Gutman lens properties than the 12-fold symmetry design. We demonstrate that the irregularity in the propagation of the 12-fold symmetry design is influenced by the hole positions rather than their quantity because both designs have a comparable number of holes (567 for the 12-fold symmetry design and 665 for the square lattice design). These findings emphasize the significance of element arrangement and pose critical considerations in the design of metasurfaces for optical applications.