# Tunable enhancement of the cylindrical Luneburg lens focusing ability with the aid of a conformal graphene strip

Iryna O. Mikhailikova, Denys M. Natarov, Sergii V. Dukhopelnykov, Ronan Sauleau, Mario Lucido, Oleksandr I. Nosych

► To cite this version:

Iryna O. Mikhailikova, Denys M. Natarov, Sergii V. Dukhopelnykov, Ronan Sauleau, Mario Lucido, et al.. Tunable enhancement of the cylindrical Luneburg lens focusing ability with the aid of a conformal graphene strip. Optics Express, 2024, 32 (23), pp.41726-41740. 10.1364/OE.529269. hal-04835640

HAL Id: hal-04835640

https://hal.science/hal-04835640v1

Submitted on 26 May 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Tunable enhancement of the cylindrical
Luneburg lens focusing ability with the aid of a
conformal graphene strip

IRYNA O. MIKHAILIKOVA,<fnref n="1" />,<fnref n="2" />,* ID DENYS M. NATAROV,<fnref n="2" />
SERGII V. DUKHOPELNYKOV,<fnref n="1" />,<fnref n="2" />,<fnref n="3" /> ID RONAN SAULEAU,<fnref n="3" />
MARIO LUCIDO,<fnref n="4" /> ID AND OLEKSANDR I. NOSYCH<fnref n="2" /> ID

**Abstract:** We consider the plane wave focusing characteristics of the layered cylindrical Luneburg lens equipped with a conformal strip of graphene, in the H-polarization case. The angular width and location of the strip is arbitrary, and its surface impedance is characterized with the aid of the quantum-physics Kubo formalism. We use a mathematically accurate full-wave analytical regularization technique, which is based on the explicit inversion of the problem static part and yields a Fredholm second-kind matrix equation. This guarantees the convergence of the resulting meshless numerical algorithm. We compute the focusing ability of a microsize lens as a function of the frequency in the wide range up to 60 THz. This analysis shows that a graphene strip, placed into the focal area of the Luneburg lens, enhances its focusing ability at the resonance frequency of the strip plasmon mode proportionally to the quality factor. This frequency is defined by the strip width and is tunable with the aid of graphene's chemical potential.

© 2024 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

## 1. Introduction

The Luneburg lens (LL) is a spherical or circular-cylindrical dielectric lens of the radius $R_M$ with the relative dielectric permittivity depending on the radial coordinate $r$ as $\epsilon(r) = 2 - (r/R_M)^2$ so that it takes the value of 1 at the lens surface. In the geometrical-optics (GO) approximation, such a lens has the focal point just at its outer surface [1-3]. In the realistic scenario, LL is finite in size and discrete [2-7], for instance, concentrically layered as in Fig. 1(a), and the outer lens layer has non-vanishing contrast with the free space. Under the plane-wave illumination, such a lens has not a focal point but finite-size focal area, which shrinks to the half-wavelength limit if the lens size gets larger.

The lens focusing ability can be estimated as the maximum field magnitude within the focal area (normalized by the plane-wave amplitude); it grows proportionally to $kR_M$, where $k$ is the wavenumber. Note that in real-life systems the focusing ability is spoiled each time when the frequency coincides with a high-Q eigenfrequency of a whispering-gallery mode of the outer lens layer [4]. To shift these frequencies off the operational range, the lens should have the outer layer permittivity as close to 1 as possible.

Without a size increment, a raise of the focusing ability of a quasi-optical antenna is possible with the aid of the resonance sub-wavelength elements placed into the focal area. Recently, this opportunity has been pointed out to in [8] (see Figs. 4 and 5 with brief discussion after them) for the simplest dielectric lens antenna shaped as a uniform circular dielectric rod decorated with a

<footnote n="1">¹Department of Mathematics, V. N. Karazin Kharkiv National University, Kharkiv, Ukraine</footnote>

<footnote n="2">²Laboratory of Micro and Nano Optics, Institute of Radio-Physics and Electronics NASU, Kharkiv 61085, Ukraine</footnote>

<footnote n="3">³IETR - UMR 6164, Université de Rennes, Rennes, France</footnote>

<footnote n="4">⁴Department of Electrical and Information Engineering “Maurizio Scarano”, University of Cassino and Southern Lazio, 03043 Cassino, FR, Italy</footnote>

<footnote>*irinamihaylikova@gmail.com</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Cross sections of the bare discrete *M*-layer concentric Luneburg lens (a) and the same lens decorated with a conformal strip of graphene (b), illuminated by the plane wave. The strip is shown as a black arc centered at $\phi = 0$.</figcaption></figure>

conformal graphene strip. The achieved enhancement is explained by the interplay of the GO focusing, a.k.a. “electromagnetic jet effect” [9,10], and the plasmon-mode resonances of the graphene strip, provided that the latter is placed into the focal area.

Indeed, in the case of the H (TE) polarization, a strip of graphene is able to support the transversal plasmon modes with low natural frequencies and moderate Q-factors. Both the frequencies and the associated Q-factors can be controlled by the graphene chemical potential and electron relaxation time, as follows from the quantum-physics Kubo theory [11] of the graphene conductivity. To help the readers, we remind the corresponding expressions in Appendix A.

Note that the plasmon modes of dielectric scatterers with graphene elements are the object of active research today [12–21]. These studies show that such composite scatterers are complicated open resonators, supporting both the dielectric-part modes and the plasmon modes of graphene elements.

It should be noted, however, that the focusing ability of a layered LL and opportunity of its improvement with the aid of the resonance graphene elements has not been studied so far. Instead, LL’s ability to collimate the main beam of the radiation pattern, when fed by a localized source placed on the lens surface, has been the main object of research in [4,5]; therefore, LL-based antennas are famous, first of all, for their high directivity.

Still, another important application of LL is its use as enhancer of the backward (monostatic radar) scattering cross section. In such applications, discrete LL is equipped with conformal metal (in the modeling, perfectly electrically conducting (PEC)) “cap” [22–24]. In contrast to more common metal right-corner reflector, such a radar target shows high reflection in a wider range of the incidence angles.

Replacement of PEC reflector with graphene one, which supports plasmon modes, still leads to the enhanced backscattering [25]. Besides, as graphene is lossy, this configuration shows enhanced absorption, especially at the plasmon resonance frequency [26].

The main challenge in the accurate full-wave modeling of such a composite scatterer as LL with conformal strip of graphene is that the lens should have quasi-optical size in terms of the free-space wavelength while the graphene strip width, to resonate on the lowest plasmon mode, should be deeply sub-wavelength. This requires, for trusted analysis, more efficient, i.e., economic and accurate, numerical techniques than conventional numerical approximations and commercial codes. Guided by these considerations, we consider the scattering from a circular layered dielectric rod, on the outer surface of which a conformal strip of graphene is placed, using the code based on the Method of Analytical Regularization (MAR) [27–29], able to overcome the above indicated difficulties and provide guaranteed convergence.
This paper builds on the research presented in the preceding conference papers [25,26], however, greatly extended and deepened. Unlike these works, where the analysis was concentrated on the far-field scattering and absorption characteristics, here we study the focusing, i.e., the near-field behavior of layered cylindrical LL, equipped with conformal graphene strip. As a figure-of-merit, we consider the total field magnitude at the GO focal point, i.e., at the outer boundary of the lens, however, from the inside of the lens.

On the possibility of layered LL fabrication, it should be stated that today such lenses are market available as the components of millimeter-wave antennas. As for the THz and infrared waves, where the graphene properties and the plasmon-mode resonances are the most promising, the existing technologies allow controlled manufacturing of graphene micro and nanotubes with dielectric filling [30]; continuing technological progress inspires optimism.

## 2. Formulation and method

Suppose that a plane H-polarized time-harmonic ($e^{-i\omega t}$) wave is incident onto a concentrically layered dielectric rod decorated with a conformal longitudinal graphene strip, normally to its axis – see Fig. 1(b). Here, $R_p$ and $\varepsilon_p$ are the outer radius and the relative dielectric permittivity of the p-th layer with $p = 1, \dots, M$, $2\delta$ is the angular width and $L = 2\delta R_M$ is the physical width of the graphene strip cross-section; then $2\theta = 2\pi - 2\delta$ is the angular width of the slot. Cylindrical coordinates $(r, \phi, z)$ are introduced co-axially with the rod. Then, the H-polarized electromagnetic field is fully characterized by its longitudinal component, $H_z$. The black arrow in Fig. 1 shows the direction of the plane wave propagation, and $\phi_0$ is the angle of incidence counted from the x-axis.

Then, the following 2-D wave-scattering problem is formulated: find the function $H_z(r, \phi)$ that satisfies

(i) the Helmholtz equation with the piecewise wavenumber $k_p$ in the domains bounded by the circles of radii $R_p$ and $R_{p+1}$, $p = 1, \dots, M$, i.e.,

$$ \Delta H_z(r, \phi) + (k_p)^2 H_z(r, \phi) = 0, \quad (1) $$

(ii) the dual conditions at $r = R_M$, namely, on the graphene arc $L: \{r = R_M, |\phi| \le \delta\}$,

$$ E_{\phi}^{M} + E_{\phi}^{M+1} = 2ZZ_{0}(H_{z}^{M} - H_{z}^{M+1}), \quad E_{\phi}^{M} = E_{\phi}^{M+1} \quad (2) $$

and on the slot arc $S: \{r = R_M, |\phi| \le \theta\}$,

$$ H_{z}^{M} = H_{z}^{M+1}, \quad E_{\phi}^{M} = E_{\phi}^{M+1}, \quad (3) $$

where $M + 1$ denotes the domain outside of the lens,

(iii) the tangential components continuity conditions at $r = R_p, p = 1, 2, \dots, M - 1$,

$$ H_{z}^{p} = H_{z}^{p+1}, \quad E_{\phi}^{p} = E_{\phi}^{p+1}, \quad (4) $$

(iv) the Sommerfeld radiation condition at infinity,

(v) the local power finiteness condition.

Here, $k_p = (\omega/c)\sqrt{\varepsilon_p} = k\sqrt{\varepsilon_p} = k\alpha_p$, where $c$ is the light velocity, $Z$ is the surface impedance of graphene (i.e. the frequency dependent complex resistivity), normalized by the free-space impedance $Z_0$ – see Appendix A, and $E_{\phi}^p = (Z_0/ik\varepsilon_p)\partial H_z^p/\partial r$. Note that $|Z| \gg 1$ at all frequencies above several GHz [11].

As known, the conditions (i) – (v) guarantee the uniqueness of the boundary-value problem solution for all real-valued $k$.
To mimic the Luneburg lens, we use the same approach as in [4], i.e. assume the layer thicknesses being identical and their relative dielectric permittivities taking the following values:

$$ \epsilon_p = 2 - (p - 1/2)^2 / M^2, \qquad (5) $$

Further, both the incident plane wave, $H_z^{inc} = e^{ik(x \cos \phi_0 + y \sin \phi_0)}$, and the total field can be presented as the angular Fourier series,

- in the domain outside of the lens,

$$ H_z^{M+1} = H_z^{inc} + \sum_{n=-\infty}^{+\infty} d_n \frac{H_n(k_{M+1}r)}{H'_n(k_{M+1}R_M)} e^{in\phi}, \quad r \ge R_M \qquad (6) $$

- in the central domain of the lens,

$$ H_z^1 = \sum_{n=-\infty}^{+\infty} c_n \frac{\alpha_1 J_n(k_1 r)}{J'_n(k_1 R_1)} e^{in\phi}, \quad r<R_1 \qquad (7) $$

- in the p-th layer domain, $R_{p-1} < r < R_p$, $p = 2, \dots, M$,

$$ H_z^p = \sum_{n=-\infty}^{+\infty} \left[ a_n^p \frac{\alpha_p J_n(k_p r)}{J'_n(k_p R_p)} + b_n^p \frac{\alpha_p H_n(k_p r)}{H'_n(k_p R_p)} \right] e^{in\phi}, \qquad (8) $$

Here, $J_n(\cdot)$ and $H_n(\cdot)$ are the Bessel and Hankel (first kind) functions of the order $n$ and $c_n, d_n, a_n^p, b_n^p$ are unknown coefficients, which should be found. Note that the functions (6)–(8) satisfy the Helmholtz equation and the radiation condition. The condition of the power finiteness in each domain leads to the following inequalities for the sets of unknowns:

$$ \sum_{n=-\infty}^{+\infty} \{ |a_n|^2, |c_n|^2 \} n^{-2} < \infty, \quad \sum_{n=-\infty}^{+\infty} \{ |b_n|^2, |d_n|^2 \} n^2 < \infty \qquad (9) $$

Substituting (7) and (8) into the continuity conditions (4) with $p=1$, we obtain the following equations:

$$ -b_n^2 L_n^2 = c_n G_n^1, \quad a_n^2 K_n^2 = c_n \Gamma_n^1. \qquad (10) $$

Considering the p-th boundary between the rings, we obtain the recurrent formulas for unknown coefficients,

$$ a_n^p G_n^p + b_n^p S_n^p = -b_n^{p+1} L_n^{p+1}, \qquad (11) $$

$$ a_n^p \Gamma_n^p + b_n^p D_n^p = a_n^{p+1} K_n^{p+1}. \qquad (12) $$

Here, the following notations are used:

$$ L_n^p = 2i[\pi k R_{p-1} J_n'(k \alpha_p R_{p-1}) H_n'(k \alpha_p R_p)]^{-1}, \qquad (13) $$

$$ K_n^p = 2i[\pi k R_{p-1} J_n'(k \alpha_p R_p) H_n'(k \alpha_p R_{p-1})]^{-1}, \qquad (14) $$

$$ \Gamma_n^p = \left[ \alpha_p \frac{J_n(k \alpha_p R_p)}{J'_n(k \alpha_p R_p)} - \alpha_{p+1} \frac{H_n(k \alpha_{p+1} R_p)}{H'_n(k \alpha_{p+1} R_p)} \right], \qquad (15) $$

$$ S_n^p = \left[ \alpha_p \frac{H_n(k \alpha_p R_p)}{H'_n(k \alpha_p R_p)} - \alpha_{p+1} \frac{J_n(k \alpha_{p+1} R_p)}{J'_n(k \alpha_{p+1} R_p)} \right], \qquad (16) $$

$$ D_n^p = \left[ \alpha_p \frac{H_n(k \alpha_p R_p)}{H'_n(k \alpha_p R_p)} - \alpha_{p+1} \frac{H_n(k \alpha_{p+1} R_p)}{H'_n(k \alpha_{p+1} R_p)} \right], \qquad (17) $$
$$G_n^p = \left[ \alpha_p \frac{J_n(k\alpha_p R_p)}{J'_n(k\alpha_p R_p)} - \alpha_{p+1} \frac{J_n(k\alpha_{p+1} R_p)}{J'_n(k\alpha_{p+1} R_p)} \right]. \quad (18)$$

In similar way, denoting the Fourier coefficients of the incident wave and its derivative as

$$f_n = i^n e^{-in\phi_0} J_n(kR_M), \quad f'_n = i^n e^{-in\phi_0} k J'_n(kR_M). \quad (19)$$

and considering the continuity conditions of the tangential fields on the outer lens boundary, we obtain that

$$d_n + f'_n = a_n^M + b_n^M. \quad (20)$$

Finally, using the dual boundary conditions (2) and (3), Eq. (20), and introducing new notations and unknowns,

$$A_n = a_n^M \left[ \alpha_M \frac{J_n(k\alpha_M R_M)}{J'_n(k\alpha_M R_M)} - \frac{H_n(kR_M)}{H'_n(kR_M)} \right] - f_n + \\ b_n^M \left[ \alpha_M \frac{H_n(k\alpha_M R_M)}{H'_n(k\alpha_M R_M)} - \frac{H_n(kR_M)}{H'_n(kR_M)} \right]_n + f_n \frac{H_n(kR_M)}{H'_n(kR_M)}, \quad (21)$$

$$F_n = \left[ f'_n \frac{H_n(kR_M)}{H'_n(kR_M)} - f_n \right] W_n, \quad (22)$$

we obtain a dual series equation for the unknowns $A_n$:

$$\left\{ \begin{array}{ll} \displaystyle \sum_{n=-\infty}^{+\infty} A_n e^{in\phi} = 0, & |\phi| \le \theta \\ \\ \displaystyle \sum_{n=-\infty}^{+\infty} (W_n - iZ)A_n e^{in\phi} = \sum_{n=-\infty}^{+\infty} F_n e^{in\phi}, & \theta \le |\phi| \le \pi \end{array} \right. \quad (23)$$

The coefficients $W_n$ are found explicitly from (23) and Eqs. (9)–(12) with the aid of recurrent formulas – see [25]. They have the following asymptotic representation at $|n| \to \infty$:

$$W_n = B_1 |n| + O(|n|^{-1}), \quad B_1 = [kR_M(\varepsilon_M + \varepsilon_{M+1})]^{-1}. \quad (24)$$

This behavior shows that the static part of Eq. (23) forms the Riemann-Hilbert Problem on the arc of the unit circle on the complex plane. Using the exact analytical solution of RHP in the form given in [31,32], we perform the analytical regularization of the dual series Eq. (23) and recast it to the following infinite matrix equation:

$$A_m = \sum_{n=-\infty}^{+\infty} [(|n| - B_1^{-1}W_n + iB_1^{-1}Z)A_n + B_1^{-1}F_n]T_{mn}(\theta), \quad (25)$$

where $m = 0, \pm 1, \pm 2, \dots$, and $T_{mn}(\theta)$ are the combinations of the Legendre polynomials of the argument $-\cos \theta$, see [32].

As can be shown, (25) is fully equivalent to the original boundary-value problem (i)-(v) of Section 2 and its solution satisfies (9). Note that the estimation (24) for $W_n$, together with (22) and the large-index behavior of $T_{mn}(\theta)$ that can be found in [32], prove that (25) is the Fredholm second kind matrix equation in the space of numerical sequences $l_2$. Therefore, its exact solution exists in this space for all real values of the wavenumber $k$ (as the purely real eigenvalues of $k$ are prohibited by the Poynting theorem). Moreover, the Fredholm nature of (25) guarantees that the solutions of the truncated equation tend to the exact solution of the infinite equation if the order of truncation is taken larger, $N \to \infty$, by the norm in the space $l_2$.
### 3. Numerical results and discussion

Before making systematic computations of the lens characteristics, we perform the check of the convergence of our code. Although it is guaranteed as explained above, the rate of convergence and the error dependence on the parameters are of interest. Here, we introduce the computational error as a function of the truncation number, N,

$$e_{ACS}(N) = |S_{ACS}(N) - S_{ACS}(N_{\max})|/S_{ACS}(N_{\max}), \qquad (26)$$

where $S_{ACS}$ is the absorption cross-section. This latter quantity corresponds to the power lost in the graphene strip of our scatterer and can be found, using the Optical Theorem, as

$$S_{ACS}(N) = -(4/k)\mathrm{Re}D_{H}(\phi_{0}, N) - S_{TSCS}(N), \qquad (27)$$

where the far-field angular scattering pattern is given by

$$D_{H}(\phi, N) = \sum_{n=-N}^{+N} d_{n}(-i)^{n}(H'_{n})^{-1} e^{in\phi}, \qquad (28)$$

and the total scattering cross-section is computed as

$$S_{TSCS}(N) = (4/k) \sum_{n=-N}^{+N} |d_n|^2 |H'_n|^{-2}, \qquad (29)$$

with coefficients $d_n$ calculated from the solutions to (25) by the expression

$$d_n = \left[ A_n + f_n - f_n' \frac{H_n(kR_M)}{H'_n(kR_M)} \right] W_n - f_n'. \qquad (30)$$

Our calculations, presented in Fig. 2, show that for any lens and strip parameters, the error decreases with N in exponential manner, with the number of layers having no significant effect on the error; in contrast, smaller values of $\delta$ or $\theta$ entail larger N to provide the same accuracy.

<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Computational error defined by (26) versus the matrix truncation number, for $N_{max} = 100$.</figcaption></figure>

It can be estimated that, to obtain D correct digits (assuming that $D \ge 2$) in the surface current (i.e. have the error at $10^{-D}$), this order should be taken as $N = |g_{plasm}|R_M + \pi/2\delta + 30D$, where $g_{plasm}$ ($|g_{plasm}| \gg k$) is the propagation constant of the plasmon guided wave supported by the
sheet of graphene [13]. This places N in the hundreds for a 20-lambda in diameter lens in the whole THz range even if D is wanted at the machine precision level. As visible in and discussed below Fig. 4 of [33] (note that considered there Nystrom-type code is numerically equivalent to MAR, as shown in [34]) commercial codes generally agree with the MAR solutions at low frequencies, however, start to deviate from them if the frequency exceeds a few THz, and finally hit a “numerical wall” because the number of unknowns reaches millions and the computation time becomes unacceptable. These and other merits of MAR-based numerical solutions are highlighted in [27–29].

Also note that, according to our results, the Optical Theorem is satisfied with machine precision (i.e., at the $10^{-16}$ level) if the field expansion coefficients are found from (25).

As known, in the GO approximation a continuous-epsilon LL has the focal point, with infinite field magnitude, at the lens outer surface $r = R_M$. That point, however, turns into a finite-area focal domain with a finite-value field maximum in the full-wave analysis of discrete LLs. When estimating the lens focusing ability, one has to take into account both the focal domain size and the field magnitude maximum. However, a reasonable compromise can be achieved if the field magnitude in the GO focus, i.e., at $r = R_M$, is taken as a simplified figure-of-merit.

In Appendix B, we present the color maps of the magnetic field magnitude along the x-axis as a function of the lens electric size, $kR_M$ (Fig. 8(a)-(c)) for discrete LLs with $M = 2, 5$ and $7$, without the graphene strip. The maps show that the focal domain is indeed close to the lens outer boundary and shrinks at higher frequencies, provided that the number of layers is $M \ge 5$.

Additionally, comparison of panels (a), (b), and (c) shows that the field absolute value maximum in the focal domain increases with larger numbers of the lens layers, $M$. Indeed, we see that the field magnitude at $r = R_M$ (i.e. the lens focusing ability) shows a linear growth with this parameter (Fig. 8(d)). This growth is spoiled, however, by the whispering-gallery mode effect, which appears above 25 THz for $M = 2$ but is moved away by taking $M$ larger and, hence, providing a better matching of the LL outer layer with the free space.

In the case of LL equipped with a graphene strip, the magnetic field is not continuous across the strip – see (2). Therefore, we define the lens-with-strip focusing ability, FA, as the field magnitude from inside, i.e., at $r = R_M - 0, \phi = \phi_0$,

$$FA = \left| \sum_{n=-N}^{+N} \left[ a_n^M \frac{\alpha_M J_n(k_M R_M)}{J'_n(k_M R_M)} + b_n^M \frac{\alpha_M H_n(k_M R_M)}{H'_n(k_M R_M)} \right] e^{in\phi_0} \right| \quad (31)$$

To study the effect of the graphene strip on the focal domain of discrete LL, we plot the magnetic field magnitude along the x-axis in the case of the normal incidence of the plane wave, $\phi_0 = 0$, see Fig. 3. Here, the graphene strip is located right in the focal domain. We present the results of the calculations for two different angular widths of the graphene strip, $2\delta = 5^\circ$ on the left panels and $2\delta = 1^\circ$ on the right panels, and take the number of layers as $M = 2, 5$ and $7$.

The focusing property of the lens is not yet visible on the left panels, where the strip is $5^\circ$ wide and the electric size of the lens is small, $0 \le kR_M \le 7.33$. Instead, the plasmon resonance $P_1$ of the graphene strip is visible at 6.13 THz. At this frequency, the field absolute value displays maximum at $r = R_M - 0$ with a standing wave inside the lens and a shadow zone behind the lens. Note that if $M = 7$, then $\epsilon_M = 1.1378$ and analytical estimations (37) and (38) in Appendix A yield the mode frequency as 6.52 THz and the Q-factor as 26.7, i.e. show good agreement with full-wave computations.

In contrast, the focusing is very well seen on the right panels, where the strip is $1^\circ$ wide and $12.56 \le kR_M \le 16.76$, as a bright area along $r = R_M$. The plasmon resonance $P_1$ of the $1^\circ$ wide strip is now around 16 THz and has larger Q-factor of 59.7, see (37), (38); the standing wave inside the lens and the shadow outside are even better visible than for a wider strip.

In Fig. 4, we plot the dependences of FA (31) versus the frequency for two values of the graphene strip angular width, $\delta = 5^\circ$ (a) and $1^\circ$ (b). Off the plasmon resonance, these curves
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** The near field magnitude along the x-axis for two values of graphene strip's angular width, $2\delta = 5^\circ$ (left column) and $1^\circ$ (right column). The incident plane wave propagates along the LL symmetry line. The number of the layers is $M = 2$ for the top raw, 5 for the middle raw, and 7 for the bottom raw.</figcaption></figure>

show approximately linear growth with the frequency, as expected for a quasi-optical antenna thanks to the GO focusing mechanism.

However, this dependence overlaps with a remarkable peak at the first plasmon mode $P_1$ frequency. The enhancement of FA is proportional to the Q-factor of that mode that depends on the graphene parameters as discussed in Appendix A.

Note that on the curves for the wider strip, of $2\delta = 5^\circ$, the next resonance, on the plasmon mode $P_3$, has the shape of a drop in FA.

The in-resonance near magnetic field patterns are shown in Figs. 4(c) and (e), respectively, while their zooms in the focal domain are given in Figs. 4(d) and (f). Their inspection tells that, at the resonance frequencies, the “focal domain” is dominated by the plasmon-mode field, i.e. is completely different from the GO case.
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** The focusing ability as a function of the frequency for two cases of the graphene strip angular width, $2\delta=5^\circ$ (a) and $1^\circ$ (b), under symmetric plane-wave excitation, $\phi_0 = 0$, and the amplitude near field patterns and their zooms in the $P_1$ (c), (d) and $P_3$ (e), (f) resonances for $2\delta=5^\circ$. The matrix truncation order is $N=160$.</figcaption></figure>

The near field in Fig. 4(f) allows to explain the drop in FA. Indeed, the central among the three hot field spots now has the opposite phase, in comparison to $P_1$, so that a field minimum appears at the strip center from inside the lens (the maximum appears outside). The resonances on the higher-order plasmon modes, for both strip widths, are not observable.

Besides, in both cases the curves for $M=2$ show periodic ripples at the frequencies above 25 THz, explained by the whispering-garden mode effect – see Appendix B.
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** The color map (a) and the plots (b) of the focusing ability versus the frequency and the chemical potential for the strip width $2\delta = 1^\circ$ under symmetric plane-wave excitation, $\phi_0 = 0$. White dashed line corresponds to Eq. (37).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Fig. 6.** The near field magnitude along the x-axis for the graphene strip width, $2\delta = 5^\circ$ (top row) and $1^\circ$ (bottom row), versus the plane wave incidence angle, $\phi_0$. The number of layers is $M = 7$.</figcaption></figure>

Thus, both the frequency and the height of the FA peak associated to the plasmon-mode resonance can be electrically tuned with the aid of the graphene chemical potential, i.e. the DC bias – this tunability is demonstrated in Fig. 5. As known, today the realistic values of the chemical potential of the best CVD samples of graphene are limited to 1 eV. This means that the plasmon-mode frequency can be tuned within the factor of 3 and even larger.
<figure><img src="image_8.png" /><figcaption>**Fig. 7.** Focusing ability as a function of the plane wave incidence angle for two cases of graphene strip angular width, $2\delta=5°$ (a) and $1°$ (b). The frequency, in each case, corresponds to the strip plasmon-mode $P_1$ resonance frequency.</figcaption></figure>

To highlight the importance of precise positioning of the graphene strip, we have built the maps of the near-field cuts along the plane-wave propagation direction (left column) and along the x-axis (right column), for different values of the incidence angle, $\phi_0$ (Fig. 6).

Here, $M = 7$ and the frequencies are 6.33 THz and 15.65 THz, which correspond to the $P_1$ resonance of the $5°$ and $1°$ graphene strip, i.e., the maximum FA values in Figs. 4(a) and (b), respectively.

These color maps show the influence of the incident wave angle of arrival, counted from the line of symmetry, which passes through the strip's middle point. If the angle of incidence shifts behind the strip edge, this leads to significant reduction in the near-field amplitude. Thus, precise placement of the strip is important for the highest near-field enhancement.

This observation is supported by the plots of FA (31) as a function of the incidence angle of the plane wave for the lens with $M = 7$ layers, at the same plasmon-mode resonance $P_1$ frequencies, as before – see Fig. 7.

**4. Conclusions**

We have shown that a narrow strip of graphene, placed on the surface of the microsized layered cylindrical Luneburg lens, is able to provide a dramatic, although narrow-band, enhancement of the lens focusing ability in the THz range. Moreover, this enhancement is electrostatically tunable in wide range, which, at the current level of technology, corresponds to the factor of 3. To study accurately the focusing by such a complicated scatterer, with the diameter up to 20 free-space wavelengths however decorated with a deeply sub-wavelength strip, we use the advanced RHP technique that belongs to the family of the MAR techniques [27–29]. MAR-RHP inverts analytically the static part of the H-polarized plane-wave scattering problem and yields the algorithm that has the mathematically guaranteed convergence and easily controlled accuracy. This accuracy is equally high off and at the plasmon-mode resonances of the graphene strip that manifest themselves most brightly in the THz range, where the surface impedance of graphene is reasonably high and dominated by the inductive reactance. As the reported here effect should be present in the other lens configurations, we believe that it can be exploited for boosting the sensitivity of various THz receiver systems.
## Appendix A

The Kubo formalism is widely recognized today as a frequency-dependent model of the quantum electron mobility in graphene [11]. Here, infinite graphene sheet thickness is considered zero and its surface conductivity, $\sigma(\omega, \mu_c, \tau, T)$, consists of two contributions, intraband and interband conductivities. The first has the analytic expression known as the Drude term [11],

$$ \sigma_{\text{intra}} = \frac{iq_e^2 k_B T}{\pi \hbar^2 (\omega + i\tau^{-1})} \left\{ \frac{\mu_c}{k_B T} + 2 \ln \left[ 1 + \exp \left( -\frac{\mu_c}{k_B T} \right) \right] \right\}, \quad (32) $$

where $\mu_c$ is the chemical potential, $\tau$ is the electron relaxation time. $q_e$ is the charge of the electron, $T$ is the temperature, $k_B$ is the Boltzmann constant, and $\hbar$ is the reduced Planck constant. If $\mu_c >> k_B T$ (note that $k_B T = 0.026$ eV at $T = 300$ K), then the $\ln(\cdot)$ term in (32) can be neglected.

The second contribution is expressed as an integral of the known functions (see [11]). If $\mu_c >> k_B T$, that integral can be reduced to a simple expression,

$$ \sigma_{\text{inter}} = \frac{iq_e^2}{4\pi\hbar} \ln \frac{2\mu_c - (\omega + i\tau^{-1})\hbar}{2\mu_c + (\omega + i\tau^{-1})\hbar} \quad (33) $$

Then, the normalized surface impedance of graphene, involved into (2), is

$$ Z(\omega) = Z_0^{-1} (\sigma_{\text{intra}} + \sigma_{\text{inter}})^{-1} \quad (34) $$

where $Z_0 = \sqrt{\mu_0/\epsilon_0}$ is the impedance of the free space.

Note that the hexagonal fine structure of graphene has a cell size less than 10 nm. Therefore, infinite-sheet expressions for graphene's conductivity are applicable (have good accuracy), if the size of the patterned graphene is around 30 nm or larger. Besides, graphene's anisotropy can be safely neglected at all frequencies below the X-ray range. The relative contribution of two terms into (34) depends on the frequency and chemical potential. According to [11], the interband conductivity, (33), is smaller in the absolute value than the intraband one, (32), in a wide range from the statics to a certain high frequency, which scales with the chemical potential, due to the dominance of the factor $\mu_c(k_B T)^{-1}$ in the Drude term. For instance, if $\tau = 0.5$ ps, $T = 300$ K and $\mu_c = 0.5$ eV, then $|\sigma_{\text{inter}}| \le 0.1|\sigma_{\text{intra}}|$ at the frequencies below 80 THz. Neglecting the interband conductivity, we obtain

$$ Z(\omega) \approx (Z_0 \sigma_{\text{intra}})^{-1} = -i(\omega + i\tau^{-1})\Omega^{-1}, \quad (35) $$

where the quantity

$$ \Omega = \frac{q_e^2 k_B T Z_0}{\pi \hbar^2} \left\{ \frac{\mu_c}{k_B T} + 2 \ln \left[ 1 + \exp \left( -\frac{\mu_c}{k_B T} \right) \right] \right\} \quad (36) $$

does not depend on the frequency. Still, at the near infrared and visible light frequencies, the description of graphene should take into account both types of conductivity.

Due to the fact that $\text{Im}Z<0$ (i.e. the reactance is inductive), infinite flat sheet of graphene can support the propagation of the plasmon wave [11]. On the patterned graphene samples, this wave bounces between the edges and creates standing waves, i.e., natural modes. According to derivations presented in [13,18], the real and imaginary parts of the complex wavenumbers of the plasmon modes of the flat strip of graphene, lying on the boundary between two material half-spaces with the relative dielectric permittivities $\epsilon_I$ and $\epsilon_{II}$, are given by, respectively,

$$ k'_m \approx \left[ \frac{\pi(m - 0.25)\Omega}{(\epsilon_I + \epsilon_{II})cL} \right]^{1/2}, \quad k''_m \approx -\frac{1}{2c\tau}. \quad (37) $$

where $m = 1, 2, \dots$ is the mode index, $L$ is the strip width, and $\Omega$ is given in (36).
Then, using the same line of reasoning as in [12], the quality factors of the plasmon modes of graphene strip as a surface-wave Fabry-Perot open resonator are obtained as

$$Q_m \approx \tau \left[ \frac{\pi(m - 0.25)c\Omega}{(\varepsilon_I + \varepsilon_{II})L} \right]^{1/2} . \quad (38)$$

Note that (37) and (38) are derived neglecting the interband contribution to the graphene conductivity. Therefore, these equations are not valid at the very high frequencies (lying in the visible range [18]) where the interband conductivity cannot be neglected. Besides, in (37) the radiation losses of the plasmon modes are neglected in favor of the losses in graphene.

## Appendix B

In this appendix, we summarize the results of the study of the focal domain and the focusing ability of the “bare” discrete M-layer cylindrical LL (i.e. without the graphene strip), cross-section of which is shown in Fig. 1(a).

In Fig. 8, we present the color maps of the field magnitude along the x-axis, as explained in the text, for M = 2, 5 and 7 – see panels (a), (b) and (c), respectively. On panel (d), we show the frequency dependences of FA with these layer numbers

<figure><img src="image_9.png" /><figcaption>Fig. 8. Near field magnitude along the x-axis of the “bare” layered LL for M = 2 (a), 5 (b) and 7 (c) and the focusing ability as a function of the frequency (d).</figcaption></figure>

Note that if M = 2 then the focusing ability shows periodic ripples at higher frequencies, which become more and more pronounced and sharper. This is the effect of the whispering-gALLERY modes, studied in [4] – they are excited at the lower frequencies for the lens with larger optical
contrast between the outer layer and the host medium. Here, we recall that $\varepsilon_2 = 1.4375$ for $M=2$, $\varepsilon_5 = 1.19$ for $M=5$, and $\varepsilon_7 = 1.1378$ for $M=7$. The ripples can be eliminated in continuous-epsilon LL, which has $\varepsilon(R_M) = 1$. Such a lens can be studied using the approach developed in [35], however, its fabrication is far beyond the existing technologies.

**Funding.** National Research Foundation of Ukraine (2020-02-0150).

**Acknowledgments.** IOM acknowledges support of the IEEE Antennas and Propagation Society via E. F. Knott Pre-Doctoral Research Grant Award, the IEEE Photonics Society via Graduate Student Scholarship, and the N. I. Akhiyezer Foundation via Student Scholarship. SVD acknowledges support of the University of Rennes, "Fondation Université de Rennes" and "Rennes Métropole" via the PAUSE-Ukraine Program of the Ministry of Higher Education, Research and Innovation, France. OIN acknowledges support of the Institute of International Education's Scholar Rescue Fund. SVD and OIN are grateful to the Institute of Electronics and Numerical Technologies of the University of Rennes for the hospitality.

The authors also thank two anonymous reviewers for their valuable comments and encouragement.

**Disclosures.** The authors declare no conflicts of interest.

**Data availability.** The authors have presented the equations, which can be directly programmed by the reader, delivering all the results reported in the paper.

**References**

1. R. K. Luneburg, *The Mathematical Theory of Optics* (Brown University Press, 1941).

2. J. A. Lock, "Scattering of an electromagnetic plane wave by a Luneburg lens III. Finely stratified sphere model," *J. Opt. Soc. Am. A* **25**(12), 2991-3000 (2008).

3. A. Parfitt, J. Graeme, J. Kot, et al., "A case for the Luneburg lens as the antenna elements for the square-kilometre array radio telescope," *Radio Sci. Bulletin* **293**, 32-37 (2000).

4. A. V. Boriskin and A. I. Nosich, "Whispering-gallery and Luneburg lens effects in a beam-fed circularly-layered dielectric cylinder," *IEEE Trans. Antennas Propag.* **50**(9), 1245-1249 (2002).

5. Z. Sipus, D. Bojanjac, and T. Komljenovic, "Electromagnetic modeling of spherically stratified lenses illuminated by arbitrary sources," *IEEE Trans. Antennas Propag.* **63**(4), 1837-1843 (2015).

6. Y.-Y. Zhao, Y. L. Zhang, M. L. Zheng, et al., "Three-dimensional Luneburg lens at optical frequencies," *Laser Photonics Rev.* **10**(4), 665-672 (2016).

7. V. Kaschten, D. Tihon, D. Lederer, et al., "Mismatched Luneburg lenses: numerical analysis and design rules," *IEEE Trans. Antennas Propag.* **72**(5), 4487-4500 (2024).

8. S. V. Dukhopelnykov, R. Sauleau, M. Garcia-Vigueras, et al., "Combined plasmon-resonance and photonic-jet effect in the THz wave scattering by dielectric rod decorated with graphene strip," *J. Appl. Phys.* (Melville, NY, U. S.) **126**(2), 023104 (2019).

9. A. Heifetz, S. C. Kong, A. V. Sahakian, et al., "Photonic nanojets," *J. Comput. Theor. Nanosci.* **6**(9), 1979-1992 (2009).

10. S. V. Dukhopelnykov, M. Lucido, R. Sauleau, et al., "Circular dielectric rod with conformal strip of graphene as tunable terahertz antenna: interplay of inverse electromagnetic jet, whispering gallery and plasmon effects," *IEEE J. Select. Topics Quantum Electron.* **27**(1), 1-8 (2021).

11. G. W. Hanson, "Dyadic Green's functions and guided surface waves for a surface conductivity model of graphene," *J. Appl. Phys.* **103**(6), 064302 (2008).

12. Z. Ullah, G. Witjaksono, I. Nawi, et al., "A review on the development of tunable graphene nanoantennas for terahertz optoelectronic and plasmonic applications," *Sensors* **20**(5), 1401 (2020).

13. S. V. Dukhopelnykov, R. Sauleau, and A.I. Nosich, "Integral equation analysis of terahertz backscattering from circular dielectric rod with partial graphene cover," *IEEE J. Quantum Electron.* **56**(6), 1-8 (2020).

14. M. Cuevas, "Graphene coated subwavelength wires: a theoretical investigation of emission and radiation properties," *J. Quant. Spectrosc. Radiat. Transfer* **200**, 190-197 (2017).

15. D. O. Herasymova, S. V. Dukhopelnykov, and A. I. Nosich, "Infrared diffraction radiation from twin circular dielectric rods covered with graphene: plasmon resonances and beam position sensing," *J. Opt. Soc. Am. B* **38**(9), C183-C190 (2021).

16. F. O. Yektushenko, S. V. Dukhopelnykov, Y. G. Rapoport, et al., "Spoiling of tunability of on-substrate graphene strip grating due to lattice-mode-induced transparency," *RSC Adv.* **12**(8), 4589-4594 (2022).

17. M. Lucido, M. V. Balaban, and A. I. Nosich, "Terahertz range plasmon and whispering gallery mode resonances in the plane wave scattering from thin microsize dielectric disk with graphene covers," *Proc. R. Soc. A.* **478**(2262), 20220126 (2022).

18. D. O. Herasymova, S. V. Dukhopelnykov, D. M. Natarov, et al., "Threshold conditions for transversal modes of tunable plasmonic nanolasers shaped as single and twin graphene-covered circular quantum wires," *IOP Nanotechnology* **33**(49), 495001 (2022).

19. M. Yasir, F. Peinetti, and P. Savi, "Enhanced graphene based electronically tunable phase shifter," *Micromachines* **14**(10), 1877 (2023).
20. Y. Jeyar, M. Antezza, and B. Guizal, "Electromagnetic scattering by a partially graphene-coated dielectric cylinder: efficient computation and multiple plasmonic resonances," *Phys. Rev. E* **107**(2), 025306 (2023).

21. F. O. Yektushenko, S. V. Dukhopelnykov, Y. G. Rapoport, et al., "Tunability of non-plasmon resonances in E-polarized terahertz wave scattering from microsize graphene strip-on-substrate grating," *Opt. Mater. Express* **13**(8), 2274-2287 (2023).

22. H. Sakurai, T. Hashidate, M. Ohki, et al., "Electromagnetic scattering by the Luneberg lens with reflecting cap," *IEEE Trans. Electromagn. Compat.* **40**(2), 94-96 (1998).

23. N. Nikolic, J. S. Kot, and S. Vinogradov, "Scattering by a Luneburg lens partially covered by a metallic cap," *J. Electromagn. Waves Applic.* **21**(4), 549-563 (2007).

24. P. Kadera, A. Jiménez-Sáez, and J. Lacik, "Sub-THz Luneburg lens enabled wide-angle frequency-coded identification tag for passive indoor self-localization," *Int. J. Microw. Wireless Technol.* **15**(1), 59-73 (2023).

25. I. O. Mikhailikova and S. V. Dukhopelnykov, "MAR algorithm for cylindrical Luneburg lens equipped with conformal graphene strip," *Proc. European Conf. Antennas Propagat.* (EuCAP 2023), Florence, 2023, pp. 1-4.

26. I. O. Mikhailikova and S. V. Dukhopelnykov, "Cylindrical Luneburg lens equipped with conformal graphene strip as efficient THz absorber," *Proc. European Microwave Conf.* (EuMC-2023), Berlin, 2023 pp. 484-487.

27. A. I. Nosich, "Method of Analytical Regularization in wave-scattering and eigenvalue problems: foundations and review of solutions," *IEEE Antennas Propag. Mag.* **41**(3), 34-49 (1999).

28. A. I. Nosich, "Method of analytical regularization in computational photonics," *Radio Sci.* **51**(8), 1421-1430 (2016).

29. M. Lucido, K. Kobayashi, F. Medina, et al., "Guest Editorial: Method of analytical regularisation for new frontiers of applied electromagnetics," *IET Microwaves, Antennas Propagat.* **15**(10), 1127-1132 (2021).

30. C. Dai, K. Agarwal, H. A. Bechtel, et al., "Hybridized radial and edge coupled 3D plasmon modes in self-assembled graphene nanocylinders," *Small* **17**(14), 2100079 (2021).

31. A. I. Nosich, Y. Okuno, and T. Shiraishi, "Scattering and absorption of E and H-polarized plane waves by a circularly curved resistive strip," *Radio Sci.* **31**(6), 1733-1742 (1996).

32. D. Colak, A. I. Nosich, and A. Altintas, "Radar cross-section study of cylindrical cavity-backed apertures with inner or outer material coating: the case of H-polarization," *IEEE Trans. Antennas Propagat.* **43**(5), 440-447 (1995).

33. O. V. Shapoval, J. S. Gomez-Diaz, J. Perruisseau-Carrier, et al., "Integral equation analysis of plane wave scattering by coplanar graphene-strip gratings in the THz range," *IEEE Trans. THz Sci. Technol.* **3**(5), 666-674 (2013).

34. S. V. Dukhopelnykov, "Comparison of two convergent numerical methods for solving the problem of wave-scattering by a dielectric rod with a conformal strip of graphene," *Proc. Int. Conf. Microwaves, Communications, Antennas, Electronic Syst.* (COMCAS-2019), Tel Aviv, 2020, pp. 1-4.

35. J. L. Tsalamengas, "Oblique scattering from radially inhomogeneous dielectric cylinders: An exact Volterra integral equation formulation," *J. Quant. Spectrosc. Radiat. Transfer* **213**, 62-73 (2018).