Modified rigorous coupled-wave analysis for grating-based plasmonic structures with a delta-thin conductive channel: far- and near-field study

Yurii M. Lyaschuk; Serhii M. Kukhtaruk; Serhii M. Kukhtaruk; Vytautas Janonis; Vadym V. Korotyeyev

doi:10.1364/JOSAA.410857

1. INTRODUCTION

Nowadays, structures with spatially periodical lateral structurization/metasurfaces are a focus of studies as the key elements of the many opto- and optoelectronics devices with broad application areas, including spectroscopy, imaging, holography, optical lithography, biochemical sensing, and military application (for more information, see recent review in Ref. [1]). For example, diffractive gratings are widely utilized in different spectral ranges from microwaves to deep ultraviolet as dispersive optical components [2], antenna elements [3], polarizers [4], etc.

Recently, great attention has been paid to exploitation of the subwavelength metasurfaced structures in terahertz (THz) and far-infrared spectral ranges. Particularly, the semiconductor structures with surface-relief gratings [5,6], quantum well (QW) heterostructures [7,8] or graphene-based structures [9] incorporated with metallic gratings are widely discussed as potential efficient emitters and detectors of the THz radiation [10]. In such structures, the gratings play a role of the coupler, facilitating the resonant interaction between incident electromagnetic ($\textit{em}$) waves and charge density waves such as surface plasmon-polaritons or two-dimensional (2D) plasmons. Moreover, detailed investigations of the resonant properties of grating-based plasmonic structures can serve as additional tool for basic characterization of the 2D electron gas (2DEG) in QWs [11,12], graphene [13–17], and other novel 2D layered materials [18].

Mentioned research is faced with the problem of the rigorous electrodynamic simulations of optical characteristics of grating-based structures containing very thin or even atomically thin conductive layers. In the simulations, these layers should be treated as delta-thin with (2D) parameters such as the sheet conductance of the electron channel. For such structures, the technique of integral equations (IEs) is often applied for the solutions of the Maxwell’s equations. This technique uses Green function formalism and is based on a reduction of the Maxwell’s system of equations to the linear IEs. The latter can be solved, for example, using Galerkin schemes with guaranteed convergence. The IE technique has been exploited for the different problems including investigations of 2D plasmon instabilities under the grating [8,19,20], detection of THz radiation [21–23], and interaction of THz radiation with conductive-strip gratings [24,25]. In spite of apparent advantages with respect to fast convergence, this method is typically formulated for modeling structures with simple geometries, and grating is treated as delta-thin.

Another powerful method for solution of grating-related electrodynamic problems is known as rigorous coupled-wave analysis (RCWA) [26–28]. This method belongs to the matrix-type methods and operates with systems of algebraic equations, which are formulated for coefficients of the Fourier expansion of the actual components of the $\textit{em}$ field. The RCWA can be applied for arbitrary complexity of the grating-based structures with any geometry of gratings. Typically, the realization of the RCWA relates to the structures where each layer is described by the bulk parameters. The RCWA method for structures with delta-thin conductive channels requires essential modifications, which were recently discussed in Refs. [29,30] for plasmonic structures with graphene. The aim of this paper is to present such modifications in more detail, focusing on the study of convergence of the proposed method on the example of calculations of optical characteristics of a metallic grating-based resonant plasmonic structure with a QW in THz frequency range. The effect of the grating depth on the plasmon resonance, analysis of the near-field pattern, and comparison to other methods are discussed. A total absorption of the plasmonic structures separating contributions of the grating-gated 2DEG and metallic grating couplers is also simulated using spatial distributions of the $\textit{em}$ field in the near-field region.

Mathematical formalism of the modified RCWA method is presented in Section 2. RCWA will be formulated for the planar diffraction problem (plane of incidence is perpendicular to grating strips) including the cases of TM and TE polarizations. The investigations of the convergence of the proposed method and comparison with the IE technique will be illustrated in Section 3 on the example of the calculations of transmission, reflection, and absorption spectra of a QW plasmonic structure with a deeply subwavelength metallic grating. The effect of the finite thickness of the metallic grating will be studied in detail. In Section 4, we will perform the analysis of the near-field patterns under conditions of the plasmon resonances. The main results will be summarized in Section 5.

2. MATHEMATICAL FORMALISM

Let us assume that multilayered structure with the grating is illuminated by plane $\textit{em}$ wave of TM polarization with frequency, $\omega$, and angle of incidence, $\theta$ (see Fig. 1). The grating of period ${a_g}$ is formed by the infinitely long in $y$-direction rectangular bars with width, ${w_g}$, and height, ${h_g}$. The structure consists of $N$ layers with thicknesses ${d_j} = {z_j} - {z_{j - 1}}$, $j = 1\ldots N$. Each $j$-layer, including the grating region, is described by the own dielectric permittivity, ${\epsilon _{\omega ,j}}(x)$. The delta-thin conductive channel is placed between $j - 1$ and $j$ layers and is described by high-frequency conduction current, ${\vec J^{2D}}(x,t)\delta (z - {z_j})$. The whole structure is set between two non-absorbing half-spaces with ${\epsilon _0}$ (at $z \lt 0$) and ${\epsilon _{N + 1}}$ (at $z \gt D$, where $D$ is the total thickness of the structure).

Fig. 1. Schematic sketch of the geometry of the multilayered plasmonic structure with 2DEG.

Download Full Size | PDF

Assuming that all components of the $\textit{em}$ field oscillate in time as $\exp (- i\omega t)$, the Maxwell’s equations written for amplitudes take the form

(1)$$\begin{split}{\rm rot}{\vec H_\omega} &= - i{k_0}{\epsilon _\omega}(x,z){\vec E_\omega} + \frac{{4\pi}}{c}\vec J_\omega ^{2D}(x)\delta (z - {z_j}), \\ {\rm rot}{\vec E_\omega}& = i{k_0}{\vec H_\omega},\end{split}$$

where ${k_0} = \omega /c$, ${\epsilon _\omega}(x,z) = {\epsilon _{\omega ,j}}(x)$ at $z \in [{z_{j - 1}},{z_j}]$, and each non-zero component of the vectors $\vec E$ and $\vec H$ is the function of $x$ and $z$ coordinates. For simplicity, we consider non-magnetic structure with magnetic permittivity equal to 1. For the case of planar diffraction and TM polarization, the non-zero components are ${E_x}, {E_z}$, and ${H_y}$. Then, Eq. (1) can be rewritten as

(2)$$\begin{split}\frac{{\partial\! {H_{\omega ,y}}}}{{\partial z}}& = i{k_0}{\epsilon _\omega}(x,z){E_{\omega ,x}} - \frac{{4\pi}}{c}J_{\omega ,x}^{2D}(x)\delta (z - {z_j}), \\ \frac{{\partial\! {H_{\omega ,y}}}}{{\partial x}}& = - i{k_0}{\epsilon _\omega}(x,z){E_{\omega ,z}}, \\ \frac{{\partial\! {E_{\omega ,x}}}}{{\partial z}} &- \frac{{\partial\! {E_{\omega ,z}}}}{{\partial x}} = i{k_0}{H_{\omega ,y}}.\end{split}$$

According the Floquet theorem, we can search for the solutions in the form of Fourier expansion,

(3)$$\!\!\!\left({\begin{array}{*{20}{l}}{{H_{\omega ,y}}(x,z)}\\{{E_{\omega ,\{x,z\}}}(x,z)}\\[3pt]{J_{\omega ,x}^{2D}(x)}\end{array}}\! \right) = \sum\limits_{m = - M}^M \left({\begin{array}{*{20}{l}}{{H_{\omega ,m,y}}(z)}\\{{E_{\omega ,m,\{x,z\}}}(z)}\\[3pt]{J_{\omega ,m,x}^{2D}}\end{array}}\! \right)\exp (i{\beta _m}x),\!$$

where ${\beta _m} = {k_0}\sqrt {{\epsilon _0}} \sin \theta + {q_m}$ with ${q_m} = 2\pi m/{a_g}$. The truncation rank $M$ (actual number of Fourier harmonics) is selected in such way to provide the convergence of the solution with a given accuracy. Using Eq. (3), Eq. (2) written for each $j $-layer in Fourier representation reads as

(4)$$\begin{split}\frac{{\partial {\textbf{H}_{y,j}}}}{{\partial z}} &= i{k_0}\hat {\cal E}_{{\rm inv},j}^{- 1}{\textbf{E}_{x,j}}, \\ \hat \beta {\textbf{H}_{y,j}}& = - {k_0}{\hat {\cal E}_j}{\textbf{E}_{z,j}}, \\ \frac{{\partial {\textbf{E}_{x,j}}}}{{\partial z}}& = i{k_0}{\textbf{H}_{y,j}} + i\hat \beta {\textbf{E}_{z,j}}.\end{split}$$

Here, we introduced the matrix notations, where Fourier vectors ${\textbf{H}_{y,j}}$, ${\textbf{E}_{\{x,z\} ,j}}$ contain $2M + 1$ corresponding Fourier components, ${H_{\omega ,m,y}}$, ${E_{\omega ,m,\{x,z\}}}$: $\hat \beta$ is the diagonal matrix formed by the elements ${\beta _m}{\delta _{m,m^\prime}}$ (here $\delta$ denotes the Kronecker delta symbol). Elements of the matrix ${\hat {\cal E}_j}$ and ${\hat {\cal E}_{{\rm inv},j}}$ are expressed through the spatial profile of the dielectric permittivity ${\epsilon _{\omega ,j}}(x)$,

(5)$$\begin{split}{[{\hat {\cal E}_j}]_{m,m^\prime}}& = \int_0^1 {\epsilon _{\omega ,j}}(\bar x)\exp (- 2\pi i[m - m^\prime]\bar x){\rm d}\bar x, \\ {[{\hat {\cal E}_{{\rm inv}j}}]_{m,m^\prime}}&= \int_0^1 \epsilon _{\omega ,j}^{- 1}(\bar x)\exp (- 2\pi i[m - m^\prime]\bar x){\rm d}\bar x,\end{split}$$

respectively. Here $\bar x = x/{a_g}$. Note, that the emergence of the $\hat {\cal E}_{{\rm inv},j}^{- 1}$ matrix in the first equation of Eq. (4) relates to the Fourier factorization rule [31,32] of a product of two discontinues functions ${\epsilon _\omega}(x)$ and ${E_{\omega ,x}}(x)$ with concurrent jump in the region of the grating (product ${\epsilon _\omega}{E_{\omega ,x}}$ is the electrical induction, which is continuous in the $x $-direction). In the spatially uniform layers, matrices $\hat {\cal E}_{{\rm inv},j}^{- 1}$ and ${\hat {\cal E}_j}$ are diagonal and identical. The application of this rule considerably improves the convergence of the results compared with the old formulation of the RCWA method [26,28], where classical Laurent rule represents the product ${\epsilon _\omega}{E_{\omega ,x}}$ as a convolution-type sum in conventional form (${\hat {\cal E}_j}$ was used instead of $\hat {\cal E}_{{\rm inv},j}^{- 1}$).

Equation (4) can be rewritten in terms of the vector ${\textbf{H}_{y,j}}$,

(6)$$\begin{split}\frac{{{\partial ^2}{\textbf{H}_{y,j}}}}{{\partial {z^2}}} = k_0^2{\hat{\boldsymbol A}_j}{\textbf{H}_{y,j}},\quad {\rm where} \; {\hat{\boldsymbol A}_j} = \hat {\cal E}_{{\rm inv},j}^{- 1}\left[{\frac{{\hat \beta \hat {\cal E}_j^{- 1}\hat \beta}}{{k_0^2}} - \hat I} \right],\end{split}$$

and $\hat I$ is the identity matrix. Having ${\textbf{H}_{y,j}}$, we can find Fourier vectors of the electric field components,

(7)$${\textbf{E}_{x,j}} = - \frac{i}{{{k_0}}}{\hat {\cal E}_{{\rm inv},j}}\frac{{\partial {\textbf{H}_{y,j}}}}{{\partial z}},\quad {\textbf{E}_{z,j}} = - \frac{1}{{{k_0}}}\hat {\cal E}_j^{- 1}\hat \beta {\textbf{H}_{y,j}}.$$

Equation (6) should be solved consequently for each $j $-layer with appropriate boundary conditions on the $j $-interfaces:

(8)$$\begin{split}\!\!\!{\textbf{E}_{x,j}}({z_j}) = {\textbf{E}_{x,j + 1}}({z_j}), {\textbf{H}_{y,j}}({z_j}) - {\textbf{H}_{y,j + 1}}({z_j}) = \frac{{4\pi}}{c}\textbf{J}_{x,j}^{2D},\!\end{split}$$

where $\textbf{J}_{x,j}^{2D}$ is the Fourier vector formed by $J_{\omega ,m,x}^{2D}$ components. The first equation in Eq. (8) expresses the continuity of the tangential component of electric field, and the second equation describes the discontinuity of magnetic field component due to the presence of the conductive 2D layer [33,34]. In the frames of the linear response theory, $J_{\omega ,m,x}^{2D} = \sigma _{\omega ,m}^{2D}{E_{\omega ,m,x}}$, where high-frequency sheet conductivity, $\sigma _{\omega ,m}^{2D}$, takes into account both frequency and spatial dispersion of the 2DEG. Our approach is valid for an arbitrary form of the conductivity of 2DEG. Particular examples of $\sigma _{\omega ,m}^{2D}$ for the 2DEG with parabolic spectrum can be found in Refs. [8,19]. For the electrons with Dirac spectrum, as in graphene, $\sigma _{\omega ,m}^{2D}$ can be obtained using Kubo formalism (see Refs. [35–37]). In the case of doped graphene in steady-state applied electric field, $\sigma _{\omega ,m}^{2D}$ can be found in Ref. [38].

Thus, the second boundary condition in Eq. (8) can be rewritten in the form

(9)$${\textbf{H}_{y,j}}({z_j}) - \hat \Gamma _j^{2D}{\textbf{E}_{x,j}}({z_j}) = {\textbf{H}_{y,j + 1}}({z_j}),$$

where diagonal matrix $\hat \Gamma _j^{2D}$ is formed by the elements $4\pi /c \times \sigma _{\omega ,m}^{2D}{\delta _{m,m^\prime}}$.

Equation (6) composes the system of ordinary second-order differential equations with constant coefficients. The solution of this system can be expressed in terms of the eigenvalues, ${\lambda _{m,j}}$, and eigen vectors, ${\vec w_{m,j}}$, of the matrix $\hat{\boldsymbol A}$,

(10)$$\begin{split}{\textbf{H}_{y,j}}(z) &= \sum\limits_{\nu = 1}^{2M + 1} {\vec w_{\nu ,j}}\big[{C_{\nu ,j}^ + \exp (- {k_0}{{\bar \lambda}_{\nu ,j}}(z - {z_{j - 1}})) } \\ &\quad+ C_{\nu ,j}^ - \exp ({k_0}{\bar \lambda _{\nu ,j}}(z - {z_j})) \big],\end{split}$$

where ${\bar \lambda _{\nu ,j}}$ is the square root (with positive real part) of the eigenvalues ${\lambda _{\nu ,j}}$. Two terms in the square brackets describe two waves: transmitted (${+}$) and reflected (${-}$) into $j$ layer. The expression for Fourier vector, ${\textbf{E}_{x,j}}(z)$, can be obtained by means of the first equation of Eq. (7) and can be written as follows:

(11)$$\begin{split}{\textbf{E}_{x,j}}(z) &= \sum\limits_{\nu = 1}^{2M + 1} {\vec v_{\nu ,j}}\big[{C_{\nu ,j}^ + \exp (- {k_0}{{\bar \lambda}_{\nu ,j}}(z - {z_{j - 1}})) } \\[-2pt] &\quad- C_{\nu ,j}^ - \exp ({k_0}{\bar \lambda _{\nu ,j}}(z - {z_j}))\big],\end{split}$$

where vector ${\vec v_{\nu ,j}} = i{\hat {\cal E}_{{\rm inv},j}}{\bar \lambda _{\nu ,j}}{\vec w_{\nu ,j}}$.

Matching magnetic and electric fields on $j $-interface according to the boundary conditions in Eqs. (8) and (9), we come to the following recurrence relationship between constants of integration, $C_{\nu ,j}^ \pm$ and $C_{\nu ,j + 1}^ \pm$:

(12)$$\begin{split}&\sum\limits_{\nu = 1}^{2M + 1} \vec w_{\nu ,j}^{-}\exp (- {k_0}{\bar \lambda _{\nu ,j}}{d_j})C_{\nu ,j}^ + + \vec w_{\nu ,j}^{+}C_{\nu ,j}^ - \\ &\quad= \sum\limits_{\nu = 1}^{2M + 1} {\vec w_{\nu ,j + 1}}\big[{C_{\nu ,j + 1}^ + + C_{\nu ,j + 1}^ - \exp (- {k_0}{{\bar \lambda}_{\nu ,j + 1}}{d_{j + 1}})} \big]; \\& \sum\limits_{\nu = 1}^{2M + 1} {\vec v_{\nu ,j}}\big[{\exp (- {k_0}{{\bar \lambda}_{\nu ,j}}{d_j})C_{\nu ,j}^ + - C_{\nu ,j}^ -} \big] \\&\quad = \sum\limits_{\nu = 1}^{2M + 1} {\vec v_{\nu ,j + 1}}\big[{C_{\nu ,j + 1}^ + - C_{\nu ,j + 1}^ - \exp (- {k_0}{{\bar \lambda}_{\nu ,j + 1}}{d_{j + 1}})} \big],\end{split}$$

where $\vec w_{\nu ,j}^{{\pm}} = {\vec w_{\nu ,j}} \pm \hat \Gamma _j^{2D}{\vec v_{\nu ,j}}$ contains parameters of 2DEG entered into matrix $\hat \Gamma _j^{2D}$. Equation (12) can be written in the compact matrix form

(13)$$\begin{split}\left({\begin{array}{*{20}{l}}{\hat{\boldsymbol{W}}_j^ - {{\hat \Lambda}_j},}&{\hat{\boldsymbol{W}}_j^ +}\\[3pt]{{{\hat{\boldsymbol{V}}}_j}{{\hat \Lambda}_j},}&{- {{\hat{\boldsymbol{V}}}_j}}\end{array}} \right)\left({\begin{array}{*{20}{l}}{\vec C_j^ +}\\[3pt]{\vec C_j^ -}\end{array}} \right) = \left({\begin{array}{*{20}{l}}{{{\hat{\boldsymbol{W}}}_{j + 1}},}&{{{\hat{\boldsymbol{W}}}_{j + 1}}{{\hat \Lambda}_{j + 1}}}\\[3pt]{{{\hat{\boldsymbol{V}}}_{j + 1}},}&{- {{\hat{\boldsymbol{V}}}_{j + 1}}{{\hat \Lambda}_{j + 1}}}\end{array}} \right)\left({\begin{array}{*{20}{l}}{\vec C_{j + 1}^ +}\\[3pt]{\vec C_{j + 1}^ -}\end{array}} \right),\end{split}$$

where each vector, $\vec C_j^ \pm$, contains $2M + 1$ integration constants, and the matrices ${\hat{\boldsymbol{W}}_j}$ and ${\hat{\boldsymbol{V}}_j}$ are formed by elements of the corresponding eigen vectors, ${\vec w_{\nu ,j}}$ and ${\vec v_{\nu ,j}}$, defined for all eigenvalues, ${\lambda _{\nu ,j}}$. Matrices $\hat{\boldsymbol{W}}_j^{\,\pm} = {\hat{\boldsymbol{W}}_j} \pm \hat \Gamma _j^{2D}{\hat{\boldsymbol{V}}_j}$, and ${\hat \Lambda _j}$ is the diagonal matrix with the elements, $\exp (- {k_0}{\bar \lambda _{\nu ,j}}{d_j}){\delta _{\nu ,\nu ^\prime}}$.

Equation (13) allows us to couple amplitudes of the reflected wave (in region $z \lt {z_0}$) and transmitted wave (in the region $z \gt {z_N}$), and consequently to calculate transmission and reflection coefficients for different diffraction orders. Indeed, components of the Fourier vectors of $y $-magnetic and $x $-electric fields in the region $z \lt 0$ ($j = 0$) are

(14)$$\begin{split}{H_{\omega ,m,y}}& = {\delta _{m,0}}\exp (- {\bar \lambda _{m,0}}{k_0}z) + {r_{\omega ,m}}\exp ({\bar \lambda _{m,0}}{k_0}z), \\ {E_{\omega ,m,x}}& = \frac{{\cos \theta}}{{\sqrt {{\epsilon _0}}}}{\delta _{m,0}}\exp (- {\bar \lambda _{m,0}}{k_0}z) - \frac{{i{{\bar \lambda}_{m,0}}}}{{{\epsilon _0}}}{r_{\omega ,m}}\exp ({\bar \lambda _{m,0}}{k_0}z),\end{split}$$

and in the region $z \gt D$ ($j = N + 1$) are

(15)$$\begin{split}{H_{\omega ,m,y}}& = {t_{\omega ,m}}\exp (- {\bar \lambda _{m,N + 1}}{k_0}(z - {z_N})), \\ {E_{\omega ,m,x}} &= \frac{{i{{\bar \lambda}_{m,N + 1}}}}{{{\epsilon _{N + 1}}}}{t_{\omega ,m}}\exp (- {\bar \lambda _{m,N + 1}}{k_0}(z - {z_N})),\end{split}$$

where ${\bar \lambda _{m,\{0,N + 1\}}} = \sqrt {\beta _m^2 - {\epsilon _{\{0,N + 1\}}}k_0^2} /{k_0}$ if ${\beta _m} \gt \sqrt {{\epsilon _{\{0,N + 1\}}}} {k_0}$ and ${\bar \lambda _{m,\{0,N + 1\}}} = - i\sqrt {{\epsilon _{\{0,N + 1\}}}k_0^2 - \beta _m^2} /{k_0}$ if otherwise. The quantities ${r_{\omega ,m}}$ and ${t_{\omega ,m}}$ are the normalized magnetic field amplitudes of the $m$th backward-diffracted (reflected) and forward-diffracted (transmitted) waves, respectively. These amplitudes form the Fourier vectors $\textbf{R}$ and $\textbf{T}$.

Using Eqs. (14) and (15), boundary conditions in Eqs. (8) and (9), and Eqs. (13), we found that Fourier vectors $\textbf{R}$ and $\textbf{T}$ are coupled through following matrix equations:

(16)$$\begin{split}&\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I} - \frac{{\cos \theta}}{{\sqrt {{\epsilon _0}}}}\hat \Gamma _0^{2D}}\\[6pt]{\frac{{\cos \theta}}{{\sqrt {{\epsilon _0}}}}\hat{\boldsymbol I}}\end{array}} \right)\vec \delta + \left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I} - {{\hat{\boldsymbol Z}}_I}\hat \Gamma _0^{2D}}\\[6pt]{{{\hat{\boldsymbol Z}}_I}}\end{array}} \right) \\ &\textbf{R} =\prod\limits_{j = 1}^N \left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol{W}}}_j},}&{{{\hat{\boldsymbol{W}}}_j}{{\hat \Lambda}_j}}\\[6pt]{{{\hat{\boldsymbol{V}}}_j},}&{- {{\hat{\boldsymbol{V}}}_j}{{\hat \Lambda}_j}}\end{array}} \right){\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol{W}}_j^ - {{\hat \Lambda}_j},}&{\hat{\boldsymbol{W}}_j^ +}\\[6pt]{{{\hat{\boldsymbol{V}}}_j}{{\hat \Lambda}_j},}&{- {{\hat{\boldsymbol{V}}}_j}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{{{\hat{\boldsymbol Z}}_{\textit{II}}}}\end{array}} \right)\textbf{T},\end{split}$$

where matrices ${\hat{\boldsymbol Z}_I}$ and ${\hat{\boldsymbol Z}_{\textit{II}}}$ are diagonal with elements ${-}i{\bar \lambda _{m,0}}/{\epsilon _0}{\delta _{m,m^\prime}}$ and $i{\bar \lambda _{m,N + 1}}/{\epsilon _{N + 1}}{\delta _{m,m^\prime}}$, respectively, and $\vec \delta$ is the vector with elements ${\delta _{m,0}}$. Equation (16) is the master system of equations of the modified RCWA method providing the way for calculation of transmission and reflection coefficients for different diffraction orders. In contrast to the previous formulation of the RCWA method [27,28], account of 2DEG leads to the nontrivial modifications. Particularly, the master system Eq. (16) contains the matrices $\hat{\boldsymbol{W}}_j^{\,\pm}$ in the right-hand side and $\hat \Gamma _0^{2D}$ in the left-hand side. The latter term describes the possible existence of the 2DEG on the top of the grating.

It should be noted that usage of this system of equations, written in the present form, applying to the situation of deep surface grating and optically dense materials, can face the problem of the computational instability. This instability is associated with the procedure of the numerical inversion of the second matrix in the right-hand side of Eq. (16) when exponentially small terms, $\exp (- {k_0}{\bar \lambda _{\nu ,j}}{d_j})$ (standing in the ${\hat \Lambda _j}$), become smaller than machine precision. To avoid this obstacle, authors in Ref. [28] proposed to use the following decomposition of the badly inverted matrix:

(17)$$\begin{split}{\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol{W}}_j^ - {{\hat \Lambda}_j},}&{\hat{\boldsymbol{W}}_j^ +}\\[3pt]{{{\hat{\boldsymbol{V}}}_j}{{\hat \Lambda}_j},}&{- {{\hat{\boldsymbol{V}}}_j}}\end{array}} \right)^{- 1}} = {\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_j},}&{\hat{\boldsymbol 0}}\\{\hat{\boldsymbol 0},}&{\hat{\boldsymbol I}}\end{array}} \right)^{- 1}}{\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol{W}}_j^ - ,}&{\hat{\boldsymbol{W}}_j^ +}\\[3pt]{{{\hat{\boldsymbol{V}}}_j},}&{- {{\hat{\boldsymbol{V}}}_j}}\end{array}} \right)^{- 1}}.\end{split}$$

Note that second inversion matrix in the right-hand side of Eq. (17) is the regular and can be numerically inverted without any difficulties. Let the product of the last $N $-terms in Eq. (16) be

(18)$$\begin{split}\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol{W}}}_N},}&{{{\hat{\boldsymbol{W}}}_N}{{\hat \Lambda}_N}}\\{{{\hat{\boldsymbol{V}}}_N},}&{- {{\hat{\boldsymbol{V}}}_N}{{\hat \Lambda}_N}}\end{array}} \right){\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_N},}&{\hat{\boldsymbol 0}}\\{\hat{\boldsymbol 0},}&{\hat{\boldsymbol I}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol X}}_N}}\\{{{\hat{\boldsymbol Y}}_N}}\end{array}} \right)\textbf{T} \equiv \left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol f}}_N}}\\{{{\hat{\boldsymbol g}}_N}}\end{array}} \right)\textbf{T},\end{split}$$

where we introduced the following designation:

\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol X}}_N}}\\{{{\hat{\boldsymbol Y}}_N}}\end{array}} \right) \equiv {\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol{W}}_N^ - ,}&{\hat{\boldsymbol{W}}_N^ +}\\[3pt]{{{\hat{\boldsymbol{V}}}_N},}&{- {{\hat{\boldsymbol{V}}}_N}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{{{\hat{\boldsymbol Z}}_{\textit{II}}}}\end{array}} \right).

Making substitutions, $\textbf{T} = \hat{\boldsymbol X}_N^{- 1}{\hat \Lambda _N}{\textbf{T}_N}$, term

\begin{split}&{\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_N},}&{\hat{\boldsymbol 0}}\\{\hat{\boldsymbol 0},}&{\hat{\boldsymbol I}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol X}}_N}}\\{{{\hat{\boldsymbol Y}}_N}}\end{array}} \right)\textbf{T} = {\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_N},}&{\hat{\boldsymbol 0}}\\{\hat{\boldsymbol 0},}&{\hat{\boldsymbol I}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_N}}\\{{{\hat{\boldsymbol Y}}_N}\hat{\boldsymbol X}_N^{- 1}{{\hat \Lambda}_N}}\end{array}} \right) \\ &\quad\times {\textbf{T}_N} = {\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_N},}&{\hat{\boldsymbol 0}}\\{\hat{\boldsymbol 0},}&{\hat{\boldsymbol I}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{{{\hat \Lambda}_N},}&{\hat{\boldsymbol 0}}\\{\hat{\boldsymbol 0},}&{\hat{\boldsymbol I}}\end{array}} \right)\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{{{\hat{\boldsymbol Y}}_N}\hat{\boldsymbol X}_N^{- 1}{{\hat \Lambda}_N}}\end{array}} \right){\textbf{T}_N} \\ &\quad= \left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{{{\hat{\boldsymbol Y}}_N}\hat{\boldsymbol X}_N^{- 1}{{\hat \Lambda}_N}}\end{array}} \right){\textbf{T}_N},\end{split}

and we can obtain that

(19)$$\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol f}}_N}}\\{{{\hat{\boldsymbol g}}_N}}\end{array}} \right)\textbf{T} = \left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol{W}}}_N}\left[{\hat{\boldsymbol I} + {{\hat \Lambda}_N}{{\hat{\boldsymbol Y}}_N}\hat{\boldsymbol X}_N^{- 1}{{\hat \Lambda}_N}} \right]}\\[6pt]{{{\hat{\boldsymbol{V}}}_N}\left[{\hat{\boldsymbol I} - {{\hat \Lambda}_N}{{\hat{\boldsymbol Y}}_N}\hat{\boldsymbol X}_N^{- 1}{{\hat \Lambda}_N}} \right]}\end{array}} \right){\textbf{T}_N}.$$

Substituting Eq. (19) into Eq. (16) and sequentially performing above mentioned transformations for each $j$ th term in the product, we can rewrite master system Eq. (16) in the computationally stable form

(20)$$\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I} - \frac{{\cos \theta}}{{\sqrt {{\epsilon _0}}}}\hat \Gamma _0^{2D}}\\[6pt]{\frac{{\cos \theta}}{{\sqrt {{\epsilon _0}}}}\hat{\boldsymbol I}}\end{array}} \right)\vec \delta + \left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I} - \frac{{\cos \theta}}{{\sqrt {{\epsilon _0}}}}\hat \Gamma _0^{2D}}\\[6pt]{{{\hat{\boldsymbol Z}}_I}}\end{array}} \right)\,\textbf{R} = \left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol f}}_1}}\\{{{\hat{\boldsymbol g}}_1}}\end{array}} \right){\textbf{T}_1}.$$

The matrices ${\hat{\boldsymbol f}_1}$ and ${\hat{\boldsymbol g}_1}$ can be found from the following recurrence relationship:

(21)$$\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol f}}_{j - 1}}}\\{{{\hat{\boldsymbol g}}_{j - 1}}}\end{array}} \right) = \left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol{W}}}_{j - 1}}\left[{\hat{\boldsymbol I} + {{\hat \Lambda}_{j - 1}}{{\hat{\boldsymbol Y}}_{j - 1}}\hat{\boldsymbol X}_{j - 1}^{- 1}{{\hat \Lambda}_{j - 1}}} \right]}\\[6pt]{{{\hat{\boldsymbol{V}}}_{j - 1}}\left[{\hat{\boldsymbol I} - {{\hat \Lambda}_{j - 1}}{{\hat{\boldsymbol Y}}_{j - 1}}\hat{\boldsymbol X}_{j - 1}^{- 1}{{\hat \Lambda}_{j - 1}}} \right]}\end{array}} \right),$$

where

(22)$$\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol X}}_{j - 1}}}\\{{{\hat{\boldsymbol Y}}_{j - 1}}}\end{array}} \right) = {\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol{W}}_{j - 1}^ - ,}&{\hat{\boldsymbol{W}}_{j - 1}^ +}\\[3pt]{{{\hat{\boldsymbol{V}}}_{j - 1}},}&{- {{\hat{\boldsymbol{V}}}_{j - 1}}}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol f}}_j}}\\[3pt]{{{\hat{\boldsymbol g}}_j}}\end{array}} \right),$$

and $j$ is varied from $N$ to 2. The vector ${\textbf{T}_1}$ relates to Fourier vector of the transmission coefficient, $\textbf{T}$, as follows:

(23)$$\textbf{T} = \prod\limits_{j = N}^1 \hat{\boldsymbol X}_j^{- 1}{\hat \Lambda _j}{\textbf{T}_1}.$$

Components of the vectors $\textbf{T}$ and $\textbf{R}$ provide the transmission, ${T_m}$, reflection, ${R_m}$, and coefficients for any $m$th diffraction order as well as total absorption, $L$. Particularly,

(24)$$\begin{split}{T_m} &= \frac{{\sqrt {{\epsilon _0}} \sqrt {{\epsilon _{N + 1}}k_0^2 - \beta _m^2}}}{{{k_0}{\epsilon _{N + 1}}\cos \theta}}|\textbf{T}[m{]|^2},{R_m} = \frac{{\sqrt {{\epsilon _0}k_0^2 - \beta _m^2}}}{{{k_0}\sqrt {{\epsilon _0}} \cos \theta}}|\textbf{R}[m{]|^2}, \\ L &= 1 - \sum\limits_m {T_m} + {R_m}.\end{split}$$

Here, the summation is taken over all numbers of visible diffraction orders, i.e., over such values of $m$, which keep the positive expressions under the square roots in the nominators. In the case of subwavelength gratings, only zero diffraction order ($m = 0$) occurs, and we have that ${T_0} = |\textbf{T}{[0]|^2}$, ${R_0} = |\textbf{R}{[0]|^2}$ (if ${\epsilon _0} = {\epsilon _{N + 1}}$), and $L \equiv {L_0} = 1 - {T_0} - {R_0}$. Equations (20)–(23) together with Eq. (24) finalize the computationally stable realization of the modified RCWA method for characterization of the multilayered plasmonic structures with grating-gated delta-thin conductive channels in the case of TM polarization.

This method can be extended for the case of TE polarization of the incident radiation. Now, actual components of the $\textit{em}$ waves are ${E_y}$, ${H_x}$, and ${H_z}$. Matrix Eq. (6) is formulated for the Fourier vectors ${\textbf{E}_{y,j}}$, where

(25)$${\hat{\boldsymbol A}_j} = {\hat \beta ^2}/k_0^2 - {\hat {\cal E}_j}.$$

The Fourier vectors of the magnetic field components are given as follows:

{\textbf{H}_{x,j}} = \frac{i}{{{k_0}}}\frac{{\partial {\textbf{E}_{y,j}}}}{{\partial z}},\quad {\textbf{H}_{z,j}} = \frac{{\hat \beta}}{{{k_0}}}{\textbf{E}_{y,j}},

and boundary conditions read as

{\textbf{E}_{y,j}}({z_j}) = {\textbf{E}_{y,j + 1}}({z_j}),\; {\textbf{H}_{x,j}}({z_j}) + \hat \Gamma _j^{2D}{\textbf{E}_{y,j}}({z_j}) = {\textbf{H}_{x,j + 1}}({z_j}) .

Finally, the master system Eq. (16) takes the form

(26)$$\begin{split}\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{\hat \Gamma _0^{2D} - \sqrt {{\epsilon _0}} \cos \theta \hat{\boldsymbol I}}\end{array}} \right)\vec \delta + \left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{\hat \Gamma _0^{2D} + {{\hat{\boldsymbol Z}}_I}}\end{array}} \right)\\ \textbf{R} = \prod\limits_{j = 1}^N \left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol{W}}}_j},}&{{{\hat{\boldsymbol{W}}}_j}{{\hat \Lambda}_j}}\\[3pt]{{{\hat{\boldsymbol{V}}}_j},}&{{{\hat{\boldsymbol{V}}}_j}{{\hat \Lambda}_j}}\end{array}} \right){\left({\begin{array}{*{20}{c}}{{{\hat{\boldsymbol{W}}}_j}{{\hat \Lambda}_j},}&{{{\hat{\boldsymbol{W}}}_j}}\\[3pt]{\hat{\boldsymbol{V}}_j^ + {{\hat \Lambda}_j},}&{\hat{\boldsymbol{V}}_j^ -}\end{array}} \right)^{- 1}}\left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\[3pt]{{{\hat{\boldsymbol Z}}_{\textit{II}}}}\end{array}} \right)\textbf{T},\end{split}$$

where ${\hat{\boldsymbol{W}}_j}$ is formed by eigen vectors, ${\vec w_{\nu ,j}}$ is obtained for the each ${\lambda _{\nu ,j}}$ eigenvalue of the matrix Eq. (25). Matrix ${\hat{\boldsymbol{V}}_j}$ contains vectors, ${\vec v_{\nu ,j}} = - i{\bar \lambda _{\nu ,j}}{\vec w_{\nu ,j}}$ and $\hat{\boldsymbol{V}}_j^ \pm = {\hat{\boldsymbol{V}}_j} \pm \hat \Gamma _j^{2D}{\hat{\boldsymbol{W}}_j}$. The matrices ${\hat{\boldsymbol Z}_I}$ and ${\hat{\boldsymbol Z}_{\textit{II}}}$ are diagonal with elements $i{\bar \lambda _{m,0}}{\delta _{m,m^\prime}}$ and ${-}i{\bar \lambda _{m,N + 1}}{\delta _{m,m^\prime}}$, respectively. Applying the procedure (see above) for stable computation of the system of Eq. (26), we can find vectors $\textbf{T}$ and $\textbf{R}$, and calculate transmission, reflection coefficients for any $m$th diffraction order as well as total absorption for the case of TE-polarized incident radiation,

(27)$$\begin{split}{T_m} &= \frac{{\sqrt {{\epsilon _{N + 1}}k_0^2 - \beta _m^2}}}{{{k_0}\sqrt {{\epsilon _0}} \cos \theta}}|\textbf{T}[m{]|^2},\quad {R_m} = \frac{{\sqrt {{\epsilon _0}k_0^2 - \beta _m^2}}}{{{k_0}\sqrt {{\epsilon _0}} \cos \theta}}|\textbf{R}[m{]|^2}, \\ L &= 1 - \sum\limits_m {T_m} + {R_m}.\end{split}$$

The proposed modified RCWA method is applied for investigation of particular plasmonic structures with a grating-gated 2DEG channel. Such structures possess the resonant properties in the THz frequency range (wavelength of order of 100 µm) for TM-polarized incident radiation due to excitations of plasmons in conductive channel of 2DEG [7,11,19,39]. Below, we will study spectral characteristics of the AlGaN/GaN-based plasmonic structure with a deeply subwavelength (micron period) metallic grating, calculating transmission (${T_0}$), reflection (${R_0}$) and absorption (${L_0}$) coefficients, their convergence versus number of the Fourier harmonics, and the near-field mapping. Also, we will pay attention to the dependence of the plasmon resonances versus geometry of the grating.

3. FAR-FIELD CHARACTERISTICS AND THEIR CONVERGENCE

Here and below, we will study the case of a TM-polarized incident $\textit{em}$ wave with incidence angle $\theta = 0$. The structure under test is formed by $N = 3$ media, embedded into the air, including the rectangular metallic grating with ${\epsilon _{\omega ,1}}(x) = {\epsilon _{\omega ,M}}\Theta ({w_g} - x) + {\epsilon _0}\Theta (x - {w_g})$ (where $\Theta (x)$ stands the Heaviside step function, $x \in [0,{a_g}]$, ${\epsilon _{\omega ,M}} = 1 + 4\pi i{\sigma _M}/\omega$, and ${\sigma _M} = 4 \times {10^{17}}\,{{\rm s}^{- 1}}$ that corresponds to the gold), and the AlGaN barrier and GaN buffer layers with constant dielectric permittivities ${\epsilon _2} = 9.2$ and ${\epsilon _3} = 8.9$, respectively. The 2D conductive channel is formed in the plane $z = {z_2}$. The matrix ${[\Gamma _j^{2D}]_{m,m^\prime}} = 4\pi /c \times \sigma _{\omega ,m}^{2D}{\delta _{m,m^\prime}}{\delta _{j,2}}$, where for description of the high-frequency properties of 2DEG we used Drude–Lorentz model $\sigma _{\omega ,m}^{2D} = {e^2}{n_{2D}}{\tau _{2D}}/{m^*}(1 - i\omega {\tau _{2D}})$ with electron effective mass, ${m^*} = 0.22 \times {m_e}$, concentration of 2DEG, ${n_{2D}} = 6 \times {10^{12}}\,{{\rm cm}^{- 2}}$, and effective scattering time, ${\tau _{2D}} = 0.5 \;{\rm ps}$. Other geometrical parameters of the structure are listed in the caption for Fig. 2. The selected parameters are close to the parameters of the experimental structures recently studied in Ref. [11].

Fig. 2. Spectra of the (a) transmission, (b) reflection, and (c) absorption coefficients for zero diffraction order calculated at three values of the grating depth: ${h_g} \equiv {d_1} = 0.05$, 1, 5 µm. Grating period, ${a_g} = 1\;{\unicode{x00B5}{\rm m}}$; width of grating bars, ${w_g} = 0.5\;{\unicode{x00B5}{\rm m}}$. Thickness of AlGaN barrier, ${d_2} = 0.025\;{\unicode{x00B5}{\rm m}}$; GaN buffer layer, ${d_3} = 1\;{\unicode{x00B5}{\rm m}}$. Black dashed lines are the results of the IE method assuming delta-thin grating with 2D conductivity, $\sigma _g^{2D} = 2 \times {10^{12}}\,{\rm cm}/{\rm s}$. Dashed–dotted lines are the results for the structure without 2DEG. All spectra are obtained at $M = 100$. (d) and (e) Dependencies of the transmission coefficient, ${T_0}$, versus number of Fourier harmonics, $M$, for two selected frequencies. Insets: relative errors, ${\delta _M}$, are plotted in logarithmic scale as a function of $M$.

Download Full Size | PDF

The spectra of the far-field characteristics such as transmission, reflection, and absorption coefficients calculated for three depths of the metallic grating are illustrated in Fig. 2. As seen, all spectra possess a strong resonance at the frequency of 1.7 THz and much more weaker resonances at frequencies of 3.5 THz and 4.75 THz. The emergence of these resonances relates to the grating-assisted interaction of the incident $\textit{em}$ wave with the plasmons in the channel of the 2DEG. At the resonance, plasmon excitation of 2DEG with wave vectors determined by the grating period can effectively absorb energy of the incident $\textit{em}$ wave. Physics of 2D plasmons and their resonant interaction with $\textit{em}$ radiation are well-described in Refs. [7,19,39–41].

The considered plasmonic structure can provide considerable absorption of the THz radiation. For example, at 1.7 THz, the absorption coefficient ${L_0} \sim 40\%$. The RCWA calculations show that absorption coefficient of the plasmonic structure is almost independent of the depth of the metallic grating [see Fig. 2(c)]. Moreover, we show that results of the IE method [19,21,23], developed for the same plasmonic structure but with delta-thin grating (see black dashed lines), and the present RCWA for the grating with ${h_g} = 0.05\;{\unicode{x00B5}{\rm m}}$ almost coincide. An increase of the grating depth only leads to the additional dispersion of transmission and reflection coefficients in the higher frequency range. This dispersion is also observed in modeling structure without 2DEG (see dashed–dotted lines). Very weak dependence of the absorption of the plasmonic structure versus grating depth (even for deep grating with ${h_g} = 5\;{\unicode{x00B5}{\rm m}}$) indicates that subwavelength highly conductive grating plays the role of almost-lossless waveguide for incident TM-polarized wave of THz frequencies. For example, absorption of the deep grating in the structure without 2DEG does not exceed 2% in the considered spectral range. The calculations of the partial losses associated with grating and 2DEG gas in the grating-2DEG plasmonic structure can be performed using the pattern of the near-field, and one will be done below in the Section 4.

The convergence of the RCWA method versus number of Fourier harmonics, M, is illustrated in Figs. 2(d) and 2(e) on example of ${T_0}$ coefficient. The proposed method provides fast convergence, and results with reasonable accuracy can be already obtained using $M \approx 30 - 50$ for all considered cases in Figs. 2(d) and 2(e). In order to quantify the convergence of the RCWA method, we introduce relative error defined as follows: ${\delta _M} = |T_0^{(M)} - T_0^{(200)}|/T_0^{(200)}$ [see insets in Figs. 2(d) and 2(e)]. For example, for shallow grating (black circles), ${\delta _{30}} = 0.2\%$ for resonant frequency 1.7 THz, and ${\delta _{30}} = 0.01\%$ for the frequency 3.5 THz. Convergence becomes worse for the deep gratings: at resonant frequency 1.7 THz, ${\delta _{30}} = 0.5\%$ (for ${h_g} = 1\;{\unicode{x00B5}{\rm m}}$) and ${\delta _{30}} = 1.3\%$ (for ${h_g} = 5\;{\unicode{x00B5}{\rm m}}$); at non-resonant frequency 3.5 THz, ${\delta _{30}} = 0.1\%$ (for ${h_g} = 1\;{\unicode{x00B5}{\rm m}}$) and ${\delta _{30}} = 0.2\%$ (for ${h_g} = 5\;{\unicode{x00B5}{\rm m}}$). Thus, estimations show that convergence of the RCWA method exhibits dependence on grating depth and frequency of the incident radiation. The cases of the deep gratings and resonant frequencies of the plasmon excitation require account of the larger numbers of Fourier harmonics.

In the case of the plasmonic structure with narrow-slit grating, ${w_g} = 0.85\;{\unicode{x00B5}{\rm m}}$ and ${a_g} = 1\;{\unicode{x00B5}{\rm m}}$, our calculations predict much more pronounced features in the optical characteristics including intensity of the plasmon resonances versus grating depth (see Fig. 3). Moreover, narrow-slit grating provides more efficient coupling between incident radiation and plasmon excitations that leads to an emergence of well-pronounced multiple plasmon resonances, which are redshifted in comparison to the previous case. The redshift of the resonant frequency is the result of a larger contribution of the gated region of 2DEG where phase velocity of the plasmons is smaller than in the ungated region of 2DEG [42].

Fig. 3. The same as in Fig. 2 at ${w_g} = 0.85\;{\unicode{x00B5}{\rm m}}$.

Download Full Size | PDF

The first plasmon resonance occurs at frequency of 1.38 THz; at this, absorption of the THz waves reaches a value of ${\sim}50\%$. However, in this spectral range, the effect of the grating thickness is still a weak. Starting from the frequencies larger than 2 THz, spectral characteristics are essentially modified by grating thickness. As seen from Fig. 3(c), deep grating suppresses plasmonic mechanisms of the absorption of THz radiation. The absorption coefficient ${L_0}$ at resonant frequency of 3.78 THz is decreased from 28% for shallow grating (${h_g} = 0.05\;{\unicode{x00B5}{\rm m}}$) to 15% for the deepest grating (${h_g} = 5\;{\unicode{x00B5}{\rm m}}$). Apparently, this effect relates to an essential increase of the reflectivity of the plasmonic structures with thicker gratings as shown in Fig. 3(b). It means that for the deeper gratings, a smaller portion of the $\textit{em}$ energy is concentrated in 2DEG as it will be further illustrated in Section 4.

It should be noted that application of the RCWA methods for accurate calculations of far-field spectral characteristics of the plasmonic structure with narrower-slit grating requires a larger number of Fourier harmonics [see Figs. 3(d) and 3(e)]. Now, the relative errors ${\delta _{30}}$ for resonant frequency 1.38 THz are equal to 0.6%, 2.2%, and 2.5% for ${h_g} = 0.05\; \unicode{x00B5}{\rm m}$, 1 µm, and 5 µm, respectively. Relative errors less than 1% for deep gratings are achievable at $M \gt 60$. Similarly to the previous case, the convergence of the RCWA method is improved at higher frequencies. So, at frequency of 3.78 THz, accuracy of computation with relative errors ${\delta _M} \lt 1\%$ is achieved at $M \gt 30$. All spectra shown in Figs. 2 and 3 are obtained at $M = 100$.

Fig. 4. Spatial distributions of amplitudes of (a) ${H_y}$, (b) ${E_x}$, and (c) ${E_z}$ components of $\textit{em}$ field in the units of amplitude of incident wave, ${E_{\rm{ins}}}$, for the structure with ${a_g} = 1\;{\unicode{x00B5}{\rm m}}$, ${w_g} = 0.5\;{\unicode{x00B5}{\rm m}}$, ${h_g} = 0.05\;{\unicode{x00B5}{\rm m}}$ at frequency $\omega /2\pi = 1.7\,{\rm THz}$. Number of Fourier harmonics, $M = 450$.

Download Full Size | PDF

4. NEAR-FIELD STUDY

Together with calculations of optical characteristics relating to the far-field, the RCWA method allows us to study geometry of the near-field. Especially, we will pay attention to the spatial distributions of the ${E_x}(x,z)$ and ${E_z}(x,z)$ components of the $\textit{em}$ fields. Absolute values of these components determine the local absorption of the $\textit{em}$ wave and can be used for extraction of partial losses in metallic grating and 2DEG. Having RCWA data on Fourier vectors ${\textbf{T}_1}$ [taking from the solutions of master system Eq. (16)], we can find vectors of the integration constants $\vec C_j^ \pm$ in each $j$ layer of the structure,

(28)$$\left({\begin{array}{*{20}{c}}{\vec C_j^ +}\\[3pt]{\vec C_j^ -}\end{array}} \right) = \left({\begin{array}{*{20}{c}}{\hat{\boldsymbol I}}\\{{{\hat{\boldsymbol Y}}_j}\hat{\boldsymbol X}_j^{- 1}{{\hat \Lambda}_j}}\end{array}} \right){\textbf{T}_j},$$

where ${\textbf{T}_j}$ can be found recurrently, ${\textbf{T}_{j + 1}} = \hat{\boldsymbol X}_j^{- 1}{\hat \Lambda _j}{\textbf{T}_j}$. Substituting found constants into Eq. (10) (with known ${\hat{\boldsymbol{W}}_j}$ and ${\hat \Lambda _j}$ matrices), we can calculate Fourier vectors ${\textbf{H}_{y,j}}(z)$ and reconstruct a spatial distribution of the ${H_y}$ component in the each $\{x,z\}$ point inside plasmonic structures using Eq. (3). Equations (14) and (15) are used for reconstruction of the near-field distribution of ${H_y}$ and ${E_x}$ components outside the plasmonic structures.

The spatial distribution of $|{H_y}(x,z)|$ for the particular case of the plasmonic structure with shallow grating is shown in Fig. 4(a). This component is tangential to the grating sides and exhibits smooth behavior with partial penetration into the grating bar. Calculations give that skin depth, $\delta = c/\sqrt {2\pi {\sigma _M}\omega}$, of the gold at frequency of 1.7 THz is equal to 0.057 µm, which is comparable with height of the grating bar. The cold zone of ${H_y}$ component occupies the middle region of the bar near the bottom face. In the plane of 2DEG, $|{H_y}(x,z)|$ is a discontinuous quantity according boundary conditions in Eq. (8).

The electric components ${E_x}$ [Fig. 4(b)] and ${E_z}$ [Fig. 4(c)] show more interesting behavior with highly non-uniform distributions. The ${E_z}(x,z)$ can be directly obtained from the second relationship in Eq. (7) (using already found Fourier vectors, ${\textbf{H}_{y,j}}(z)$) and Eq. (3). For the correct reconstruction of the ${E_x}$ component, we follow the method discussed in Refs. [43,44]. In the region of the grating, $z \in [0,{h_g}]$, the $x$ component of the electric field is the normal to the grating bar’s sides, and one has a discontinuity. It is more effective to reconstruct a continuous quantity, the component of the displacement field, ${D_x}(x,z)$, which can be easily calculated from the derivative of the ${H_y}$ component with respect to the $z$ coordinate [see the first equation within Eq. (2)]. Then ${E_x}(x,y) = {D_x}(x,z)/{\epsilon _1}(x,z)$, where dielectric permittivity ${\epsilon _1}(x,z)$ is the known discontinuous function. Such method allows us to partially avoid an emergence of the unphysical spurious oscillations, known as Gibb’s phenomenon. Nevertheless, reconstruction of the near-field patterns requires account of much more Fourier harmonics than for calculations of the far-field characteristics. This circumstance was discussed in Ref. [44].

Fig. 5. Same as in Fig. 4 for ${h_g} = 1\;{\unicode{x00B5}{\rm m}}$.

Download Full Size | PDF

As seen, both ${E_x}$ and ${E_z}$ components demonstrate the field concentration effect. The energy of the $\textit{em}$ field is mainly concentrated near the ridges (hot zone I) of the metallic bars and in the region between grating bars and 2DEG (hot zone II). In the hot zones (I) and (II), both electric components are essentially enhanced. In the hot zone (II), the ${E_z}$ component predominantly dominates. The specific formation of the cold zone for the ${E_z}$ component at the vertical axis $x = {w_g}/2$ and $x = ({a_g} + {w_g})/2$ reflects the quadruple-related symmetry of the near-field (for details, see Ref. [39]). In the hot zones, amplitudes of the electric components can be in several tens of times larger than the amplitude of the incident wave. In spite of the magnetic component, the penetration of the electric components inside the metallic bar is strongly suppressed, which is the result of the edge effects. As seen, the ${E_x}$ component mainly penetrates to the grating’s bars from the upper and back faces, and ${E_z}$ from the side faces as the tangential ones for corresponding faces.

A similar geometry of the near-fields is realized for the case of the deep grating (see Fig. 5). The mappings of the ${E_x}$ and ${E_z}$ components show that the incident wave passes through the subwavelength grating in the form of TEM mode, i.e., in the grating slit, the wave has predominantly polarization along the $x$ direction with almost constant amplitude.

Additionally, we used COMSOL Multiphysics [45] to validate independently the obtained results by the finite element method. The Wave Optics module [45] is used to solve Maxwell’s equations for the system, which is shown in Fig. 1. The 2DEG was introduced as the surface current density at the interface of AlGaN and GaN. The uniform quadratic mesh with 5 nm size is used to resolve near-field components (see Fig. 6). COMSOL’s results of the electric component distributions for the case of the shallow grating with ${h_g} = 0.05\;{\unicode{x00B5}{\rm m}}$ are shown in Fig. 6. Excellent agreement between the modified RCWA and finite element methods is demonstrated in both far- and near-field studies.

Fig. 6. Same as in Fig. 4 obtained by usage of COMSOL Multiphysics.

Download Full Size | PDF

The spatial distributions of the $|{E_x}|(x,z)$ and $|{E_z}|(x,z)$ can be used for calculation of the partial losses in the grating, ${L_{\textit{gr}}}$ and 2DEG, ${L_{2DEG}}$,

(29)$${L_{\textit{gr}}} = \frac{{4\pi {\sigma _M}}}{c}\frac{{\int_0^{{w_g}} \int_0^{{h_g}} {\rm d}x{\rm d}z|{E_x}{|^2} + |{E_z}{|^2}}}{{\int_0^{{a_g}} {\rm d}x|{E_{\rm{ins}}}{|^2}}}$$

and

(30)$${L_{2DEG}} = \frac{{4\pi {\rm Re}[\sigma _\omega ^{2D}]}}{c}\frac{{\int_0^{{a_g}} {\rm d}x|{E_x}(x,{z_2}{{)|}^2}}}{{\int_0^{{a_g}} {\rm d}x|{E_{\textit{ins}}}{|^2}}},$$

where ${E_{\rm{ins}}}$ is the amplitude of the incident wave.

For the case in Fig. 4, we obtained that ${L_{\textit{gr}}} = 0.184\%$, which consists of $0.0741\%$ contribution of the ${E_x}$ component and 0.11% contribution of the ${E_z}$-component. Losses in 2DEG are considerably larger, ${L_{2DEG}} = 37.92\%$. Total losses from the near-field patterns are ${L_0} = {L_{2DEG}} + {L_{\textit{gr}}} = 38.11\%$. Calculations of the ${L_0}$ from far-field characteristics give the almost same value 38.21%. For the case of deep grating (see Fig. 5), we obtained the increase of absorption in the grating bars, ${L_{\textit{gr}}} = 0.56\%$ (with 0.072% and 0.49% contributions for ${E_x}$ and ${E_z}$ components, respectively), with almost same value of absorption in 2DEG ${L_{2DEG}} = 37.97\%$. The total losses ${L_0} = 38.53\%$ that almost coincide with the number 38.58% were obtained from far-field characteristics.

Calculations of the partial losses at the frequency of the first-order plasmon resonance indicate that the incident $\textit{em}$ wave is mainly absorbed by 2DEG and this absorption weakly depends on thickness of grating bars. This fact is illustrated by the spatial distribution of the amplitude of the $x $-component of the electric field, $|{E_x}(x,{z_2})|$, in the plane of 2DEG, calculated at the frequencies of the first [Fig. 7(a)] and second [Fig. 7(b)] plasmon resonances at three values of the grating depth. As seen, all three (grey, red, green) curves for lower frequency ($\omega /2\pi = 1.7\,{\rm THz}$) almost coincide, and all of them exhibit non-uniform, oscillating-like behavior in the gated region with almost flat distribution in the ungated region. The obtained distribution denotes that a larger part of $\textit{em}$ energy is absorbed in the gated region, i.e., under the metallic strip.

Fig. 7. Distribution of $|{E_x}(x,{z_2})|$ on one spatial period of the plasmonic structure with ${w_g} = 0.5\;{\unicode{x00B5}{\rm m}}$ and ${a_g} = 1\;{\unicode{x00B5}{\rm m}}$ at two resonant frequencies.

Download Full Size | PDF

For the higher frequency ($\omega /2\pi = 3.5\,{\rm THz}$), spatial distribution of the $|{E_x}(x,{z_2})|$ quantity [Fig. 7(b)] acquires a more complicated form with several spatial oscillations in the gated region. At this, the effect of the grating thickness becomes visible, i.e., the deep gratings start to screen the interaction of the $\textit{em}$ wave with 2DEG. The emergence of the several spatial oscillations in the distribution of $|{E_x}(x,{z_2})|$ leads to suppression of absorptivity of the plasmonic structures at higher order plasmon resonances. Also, according the Drude model, at higher frequencies, the response of electron gas on the $\textit{em}$ wave becomes weaker, which leads to a decrease of the prefactor standing in Eq. (30). This prefactor, $4\pi {\rm Re}[\sigma _\omega ^{2D}]/c$, for two considered frequencies $\omega /2\pi = 1.7$ and $\omega /2\pi = 3.5\,{\rm THz}$ is equal to 0.049% and 0.012%, respectively. Using the obtained distributions in Fig. 7(b), we found that for ${h_g} = 0.05, 1, 5\;{\unicode{x00B5}{\rm m}}$, ${L_{2DEG}} = 4.42, 4.3$ and 3.8%, and the corresponding values of total losses are calculated from far-field characteristics, ${L_0} = 4.8, 5.3, 6.3\%$. Note that, for the structure with the deepest grating, the absorptions in 2DEG and the grating bars become comparable.

The plasmonic structure with narrow-slit grating provides more efficient coupling between 2DEG and $\textit{em}$ radiation. The distributions of $|{E_x}(x,{z_2})|$ calculated for the structure with ${w_g} = 0.85\;{\unicode{x00B5}{\rm m}}$ at two resonant frequencies $\omega /2\pi = 1.38\,{\rm THz}$ (first-order plasmon resonance) and $\omega /2\pi = 3.78\,{\rm THz}$ (third-order plasmon resonance) are shown in Fig. 8. As seen, the geometry of the distributions obtained for the frequency of first-order plasmon resonance [Fig. 8(a)] is similar to the previous case depicted in Fig. 7(a). However, the wider gated region of 2DEG integrally provides larger contribution to the absorption of the $\textit{em}$ wave by 2DEG. The corresponding values of ${L_{2DEG}}$ are following: $47.5\%$ (for ${h_g} = 0.05\;{\unicode{x00B5}{\rm m}}$), $46.6\%$ (for ${h_g} = 1\;{\unicode{x00B5}{\rm m}}$), and $41.2\%$ (for ${h_g} = 5\;{\unicode{x00B5}{\rm m}}$). At this, ${L_0} = 48.2, 47.6, 43.9\%$, respectively. The distributions in Fig. 8(b) obtained at the frequency of third-order plasmon resonance demonstrate multiple spatial oscillations. The number of such oscillations is proportional to the order of plasmon resonances. Also, we see that the deepest grating with ${h_g} = 5\;{\unicode{x00B5}{\rm m}}$ essentially suppresses the plasmon absorption of the $\textit{em}$ wave. The corresponding values of ${L_{2DEG}}$ are the following: 27.2% (for ${h_g} = 0.05\;{\unicode{x00B5}{\rm m}}$), 24.1% (for ${h_g} = 1\;{\unicode{x00B5}{\rm m}}$), and 11.8% (for ${h_g} = 5\;{\unicode{x00B5}{\rm m}}$). At this, ${L_0} = 28.7\%$, 26.1%, and 15.1%, respectively.

Fig. 8. The same as in Fig. 7 for ${w_g} = 0.85\;{\unicode{x00B5}{\rm m}}$.

Download Full Size | PDF

5. SUMMARY

We have developed a computationally stable RCWA method for the solution of Maxwell’s equations in the case of multilayered plasmonic structures with a delta-thin grating-gated conductive channel. The method was formulated for a planar diffraction problem for both TM and TE polarization of an incident wave. The method was implemented for investigation of far- and near-field characteristics of the particular plasmonic structures based on AlGaN/GaN heterostructure with a deeply subwavelength metallic grating coupler.

The calculations of the far-field characteristics including transmission, reflection, and absorption coefficients for zero diffraction order were performed in THz frequency range where the considered structure has multiple resonances related to the excitations of 2D plasmons in the conductive channel of AlGaN/GaN heterostructure. The dependence of these characteristics versus grating parameters and their convergence versus the number of the Fourier harmonics were analyzed.

We found that spectra of transmission and reflection coefficients in the lower frequency range, 0.2 THz, have weak dependence on grating depth. Results for both shallow (${h_g}/{a_g} = 0.05 \;{\unicode{x00B5}{\rm m}}/1\;{\unicode{x00B5}{\rm m}}$) and deep (${h_g}/{a_g} = 1$, $5\, {\unicode{x00B5}{\rm m}}/1\;{\unicode{x00B5}{\rm m}}$) gold grating are almost identical and coincide with the results of the IE method where grating is treated as delta-thin. In a higher frequency range, 3.5 THz, increase of the grating depth suppresses transmission with the increasing of the reflection coefficients. At the same time, the absorption spectrum remain less sensitive to the grating depth. We showed that dispersion of far-field characteristics on ${h_g}$ becomes more pronounced for a narrower-slit grating with ${w_g}/{a_g} = 0.85\;{\unicode{x00B5}{\rm m}}/1\;{\unicode{x00B5}{\rm m}}$ than for a wide-slit grating with ${w_g}/{a_g} = 0.5\;{\unicode{x00B5}{\rm m}}/1\;{\unicode{x00B5}{\rm m}}$. We showed that convergence of the calculations depends on geometrical parameters of the grating and frequencies. Better convergence is achieved for shallow and wide-slit grating with relative errors of ${\sim}0.05 {-} 0.2\%$ (in dependence on frequency) with 30 Fourier harmonics. For deep and narrow-slit grating, relative errors of ${\sim}0.1..1\%$ are achieved at $M \sim 50$.

Procedure of the calculations of the near-field characteristics was discussed in detail. Analysis of spatial distribution of the amplitudes of the $\textit{em}$ wave’s components in the near-field zone reveals the main physical peculiarities of the interaction of the plasmonic structure with incident radiation. It was shown that subwavelength metallic grating plays the role of a perfect waveguide for the incident wave, concentrator of the $\textit{em}$ energy, and polarization rotator. In the region of the grating slit, the $\textit{em}$ wave has predominantly lateral polarization with amplitude close to the amplitude of the incident wave. The hot zone is formed in the region between grating bars and 2DEG where the $\textit{em}$ wave has predominantly vertical polarization with amplitudes that can in 100 times exceed the amplitude of the incident wave.

The pattern of the near-field also was used for the calculations of the partial losses related to the grating and 2DEG. It was shown that at the frequencies of the plasmon resonances, the structure can efficiently absorb THz radiation with absorption coefficient values in the order of 20–50% (independent of the grating filling factor and order of plasmon resonance). We found that the contribution of 2DEG to the total losses is dominant at low-frequency plasmon resonances with weak dependence on the grating depth. At high-frequency plasmon resonances, the effect of the grating depth becomes essential. The deep gratings can effectively screen interaction of the $\textit{em}$ waves with plasmon oscillation in 2DEG that leads to a decrease of the total absorption of THz radiation.

Also, it should be noted that the proposed modified RCWA has several advantages over conventional volumetric RCWA. First, our realization of the RCWA method allows us to avoid additional numerical manipulation with matrices that can reduce the computational time. For the considered structure, we have 25% in terms of computation time savings in comparison with conventional RCWA at the volumetric treatment of the conductive layer. This value can be increased in simulation of structures with a stack of 2D conductive layers. Second, we operate with one parameter, 2D concentration, ${n_{2D}}$, instead of two independent parameters of bulk concentration, ${n_{3D}}$, and thickness of the layer, $d$. It can be convenient for the metrology of the structures at the processing of the experimental data.

We suggest that the proposed modified RCWA algorithm can be effectively used for the modeling of the optical characteristics of various kinds of plasmonic structures with 2D conductive channels, including QW- or graphene-based structures, and the results of the paper provide deeper insight on physics of the interaction of THz radiation with grating-gated plasmonic structures.

Funding

Lietuvos Mokslo Taryba (DOTSUT-184); European Regional Development Fund (No. 01.2.2-LMT-K-718-03-0096; Bundesministerium für Bildung und Forschung (VIP+ “Nanomagnetron”).

Acknowledgment

The authors thank Prof. V. A. Kochelap (ISP NASU, Ukraine) and Dr. I. Kašalynas (FTMC, Lithuania) for fruitful discussions of the various aspects of this work. V. J. and V. V. K. acknowledge the support from the Research Council of Lithuania (Lietuvos mokslo taryba) through the “T-HP” Project (Grant No. DOTSUT-184) funded by the European Regional Development Fund according to the supported activity “Research Projects Implemented by World-class Researcher Groups” under the Measure No. 01.2.2-LMT-K-718-03-0096. S. M. K. would like to acknowledge the support of the Bundesministerium für Bildung und Forschung (VIP+ “Nanomagnetron”).

Disclosures

The authors declare no conflicts of interest.

REFERENCES AND NOTE

1. E. Popov, ed. Gratings: Theory and Numeric Applications, 1st ed. (Presses Universitaires de Provence, 2012).

2. W. Neumann, Fundamentals of Dispersive Optical Spectroscopy Systems (SPIE, 2014).

3. H. F. Hammad, Y. M. M. Antar, A. P. Freundorfer, and M. Sayer, “A new dielectric grating antenna at millimeter wave frequency,” IEEE Trans. Antennas Propag. 52, 36–44 (2004). [CrossRef]

4. S. Shena, Y. Yuana, Z. Ruana, and H. Tan, “Optimizing the design of an embedded grating polarizer for infrared polarization light field imaging,” Results Phys. 12, 21–31 (2019). [CrossRef]

5. G. A. Melentev, V. A. Shalygin, L. E. Vorobjev, V. Y. Panevin, D. A. Firsov, L. Riuttanen, S. Suihkonen, V. V. Korotyeyev, Y. M. Lyaschuk, V. A. Kochelap, and V. N. Poroshin, “Interaction of surface plasmon polaritons in heavily doped GaN microstructures with terahertz radiation,” J. Appl. Phys. 119, 093104 (2016). [CrossRef]

6. V. Janonis, S. Tumenas, P. Prystawko, J. Kacperski, and I. Kašalynas, “Investigation of n-type gallium nitride grating for applications in coherent thermal sources,” Appl. Phys. Lett. 116, 112103 (2020). [CrossRef]

7. V. V. Popov, D. V. Fateev, O. V. Polischuk, and M. S. Shur, “Enhanced electromagnetic coupling between terahertz radiation and plasmons in a grating-gate transistor structure on membrane substrate,” Opt. Express 18, 16771–16776 (2010). [CrossRef]

8. V. V. Korotyeyev, V. A. Kochelap, S. Danylyuk, and L. Varani, “Spatial dispersion of the high-frequency conductivity of two-dimensional electron gas subjected to a high electric field: collisionless case,” Appl. Phys. Lett. 113, 041102 (2018). [CrossRef]

9. V. Ryzhii, T. Otsuji, and M. Shur, “Graphene based plasma-wave devices for terahertz applications,” Appl. Phys. Lett. 116, 140501 (2020). [CrossRef]

10. T. Otsuji and M. Shur, “Terahertz plasmonics: good results and great expectations,” IEEE Microw. Mag. 15(7), 43–50 (2014). [CrossRef]

11. D. Pashnev, T. Kaplas, V. Korotyeyev, V. Janonis, A. Urbanowicz, J. Jorudas, and I. Kašalynas, “Terahertz time-domain spectroscopy of two-dimensional plasmons in AlGaN/GaN heterostructures,” Appl. Phys. Lett. 117, 051105 (2020). [CrossRef]

12. D. Pashnev, V. Korotyeyev, V. Janonis, J. Jorudas, T. Kaplas, A. Urbanowicz, and I. Kašalynas, “Experimental evidence of temperature dependent effective mass in AlGaN/GaN heterostructures observed via THz spectroscopy of 2D plasmons,” Appl. Phys. Lett. 117, 162101 (2020). [CrossRef]

13. B. Yan, J. Fang, S. Qin, Y. Liu, Y. Zhou, R. Li, and X.-A. Zhang, “Experimental study of plasmon in a grating coupled graphene device with a resonant cavity,” Appl. Phys. Lett. 107, 191905 (2015). [CrossRef]

14. M. M. Jadidi, A. B. Sushkov, R. L. Myers-Ward, A. K. Boyd, K. M. Daniels, D. K. Gaskill, M. S. Fuhrer, H. D. Drew, and T. E. Murphy, “Tunable terahertz hybrid metal–graphene plasmons,” Nano Lett. 15, 7099–7104 (2015). [CrossRef]

15. B. Zhao and Z. M. Zhang, “Strong plasmonic coupling between graphene ribbon array and metal gratings,” ACS Photon. 2, 1611–1618 (2015). [CrossRef]

16. H. Lu, J. Zhao, and M. Gu, “Nanowires-assisted excitation and propagation of mid-infrared surface plasmon polaritons in graphene,” J. Appl. Phys. 120, 163106 (2016). [CrossRef]

17. S. M. Kukhtaruk, V. V. Korotyeyev, V. A. Kochelap, and L. Varani, “Interaction of THz radiation with plasmonic grating structures based on graphene,” in International Conference on Mathematical Methods in Electromagnetic Theory, Lviv, Ukraine (2016), pp. 196–199.

18. T. Low, A. Chaves, J. D. Caldwell, A. Kumar, N. X. Fang, P. Avouris, T. F. Heinz, F. Guinea, L. Martin-Moreno, and F. Koppens, “Polaritons in layered two-dimensional materials,” Nat. Mater. 16, 182–194 (2017). [CrossRef]

19. S. A. Mikhailov, “Plasma instability and amplification of electromagnetic waves in low-dimensional electron systems,” Phys. Rev. B 58, 1517–1532 (1998). [CrossRef]

20. V. V. Korotyeyev and V. A. Kochelap, “Plasma wave oscillations in a nonequilibrium two-dimensional electron gas: electric field induced plasmon instability in the terahertz frequency range,” Phys. Rev. B 101, 235420 (2020). [CrossRef]

21. D. V. Fateev, V. V. Popov, and M. S. Shur, “Plasmon spectra transformation in grating-gate transistor structure with spatially modulated two-dimensional electron channel,” Semiconductor 44, 1406–1462 (2010). [CrossRef]

22. V. V. Popov, D. V. Fateev, T. Otsuji, Y. M. Meziani, D. Coquillat, and W. Knap, “Plasmonic terahertz detection by a double-grating-gate field-effect transistor structure with an asymmetric unit cell,” Appl. Phys. Lett. 99, 243504 (2011). [CrossRef]

23. Y. M. Lyaschuk and V. V. Korotyeyev, “Theory of detection of terahertz radiation in hybrid plasmonic structures with drifting electron gas,” Ukr. J. Phys. 62, 889–902 (2017). [CrossRef]

24. O. V. Shapoval, J. S. Gomez-Diaz, J. Perruisseau-Carrier, J. R. Mosig, and A. I. Nosich, “Integral equation analysis of plane wave scattering by coplanar graphene-strip gratings in the THz range,” IEEE Trans. Terahertz Sci. Technol. 3, 666–674 (2013). [CrossRef]

25. O. V. Shapovala and A. I. Nosich, “Finite gratings of many thin silver nano strips: optical resonances and role of periodicity,” AIP Adv. 3, 042120 (2013). [CrossRef]

26. M. G. Moharam and T. K. Gaylord, “Rigorous coupled-wave analysis of planar-grating diffraction,” J. Opt. Soc. Am. 71, 811–818 (1981). [CrossRef]

27. M. G. Moharam, E. B. Grann, D. A. Pommet, and T. K. Gaylord, “Formulation for stable and efficient implementation of the rigorous coupled-wave analysis of binary gratings,” J. Opt. Soc. Am. A 12, 1068–1076 (1995). [CrossRef]

28. M. G. Moharam, D. A. Pommet, E. B. Grann, and T. K. Gaylord, “Stable implementation of the rigorous coupled-wave analysis for surface-relief gratings: enhanced transmittance matrix approach,” J. Opt. Soc. Am. A 12, 1077–1086 (1995). [CrossRef]

29. S. Inampudi and H. Mosallaei, “Tunable wideband-directive thermal emission from SiC surface using bundled graphene sheets,” Phys. Rev. B 96, 125407 (2017). [CrossRef]

30. S. Inampudi, V. Toutam, and S. Tadigadapa, “Robust visibility of graphene monolayer on patterned plasmonic substrates,” Nanotechnology 30, 015202 (2019). [CrossRef]

31. L. Li, “Use of Fourier series in the analysis of discontinuous periodic structures,” J. Opt. Soc. Am. A 13, 1870–1876 (1996). [CrossRef]

32. E. Popov, M. Nevière, and N. Bonod, “Factorization of products of discontinuous functions applied to Fourier–Bessel basis,” J. Opt. Soc. Am. A 21, 46–52 (2004). [CrossRef]

33. E. Bleszynski, M. Bleszynski, and T. Jaroszewicz, “Surface-integral equations for electromagnetic scattering from impenetrable and penetrable sheets,” IEEE Antennas Propag. Mag. 35(6), 14–25 (1993). [CrossRef]

34. T. L. Zinenko and A. I. Nosich, “Plane wave scattering and absorption by resistive-strip and dielectric-strip periodic gratings,” IEEE Trans. Antennas Propag. 46, 1498–1505 (1998). [CrossRef]

35. V. P. Gusynin, S. G. Sharapov, and J. P. Carbotte, “Sum rules for the optical and Hall conductivity in graphene,” Phys. Rev. B 75, 165407 (2007). [CrossRef]

36. L. A. Falkovsky and S. S. Pershoguba, “Optical far-infrared properties of a graphene monolayer and multilayer,” Phys. Rev. B 76, 153410 (2007). [CrossRef]

37. M. V. Balaban, O. V. Shapoval, and A. I. Nosich, “THz wave scattering by a graphene strip and a disk in the free space: integral equation analysis and surface plasmon resonances,” J. Opt. 15, 114007 (2013). [CrossRef]

38. S. M. Kukhtaruk, V. A. Kochelap, V. N. Sokolov, and K. W. Kim, “Spatially dispersive dynamical response of hot carriers in doped graphene,” Physica E 79, 26–37 (2016). [CrossRef]

39. Y. M. Lyaschuk and V. V. Korotyeyev, “Interaction of a terahertz electromagnetic wave with the plasmonic system ‘grating–2D–gas.’ Analysis of features of the near field,” Ukr. J. Phys. 59, 495–504 (2014). [CrossRef]

40. A. V. Chaplik, “Absorption and emission of electromagnetic waves by two-dimensional plasmons,” Surf. Sci. Rep. 5, 289–335 (1985). [CrossRef]

41. M. Dyakonov and M. Shur, “Shallow water analogy for a ballistic field effect transistor: new mechanism of plasma wave generation by dc current,” Phys. Rev. Lett. 71, 2465–2467 (1993). [CrossRef]

42. Parametrical studies of the plasmon resonance at different configurations of the grating coupler are reported in Ref. [39] and in the recent experimental paper Ref. [11].

43. K.-H. Brenner, “Aspects for calculating local absorption with the rigorous coupled-wave method,” Opt. Express 18, 10369–10376 (2010). [CrossRef]

44. M. Weismann, D. F. G. Gallagher, and N. C. Panoiu, “Accurate near-field calculation in the rigorous coupled-wave analysis method,”J. Opt. 17, 125612 (2015). [CrossRef]

45. COMSOL AB, “COMSOL Multiphysics v. 5.5,” https://www.comsol.com.

Modified rigorous coupled-wave analysis for grating-based plasmonic structures with a delta-thin conductive channel: far- and near-field study

Abstract

1. INTRODUCTION

2. MATHEMATICAL FORMALISM

3. FAR-FIELD CHARACTERISTICS AND THEIR CONVERGENCE

4. NEAR-FIELD STUDY

5. SUMMARY

Funding

Acknowledgment

Disclosures

REFERENCES AND NOTE

Cited By

Figures (8)

Equations (34)

Journal of the Optical Society of America A