Plasmon modulation in three-dimensional periodic structure of graphene ribbons

Recently, plasmon dynamics in graphene has become an emerging subfield, benefitting from the remarkable electrical and optical properties [1–3], such as the relatively low loss compared to customary matals [4, 5] and the tunable ability utilizing electrical gating or chemical doping [6]. For instance, plasmon-polaritons (PP) in graphene were experimentally and theoretically [7–17] studied extensively and were also applied in optical modulators [18], photodectectors [19], biochemical sensors, and so on.

So far the most issues on the plasmon dynamics in graphene mainly focused on planar structures, for instance, on the planar graphene superlattices [14], on the periodic grid of graphene ribbons [3, 6, 20]. To improve the efficiency of photo-detectors or biochemical sensors, or to improve the ability of optical modulation, we performed a preliminary study of the plasmon dynamics based on graphene in a three-dimensional structure at the first time.

Here we study the PP spectra in a three-dimensional periodic structure of graphene ribbons embedded in a medium. For a fixed carrier concentration, the influence on the PP spectra from ribbon width, interface distance and period in the interface were given. Furthermore, we compared our results with those in a planar periodic structure of graphene ribbons and sharp differences were also shown. We believe that the structure may have applications in graphene-based photodetector and sensor.

For practical photodetectors, we always face a challenge, mismatch. For a certain wave vector, the frequency of PP, at least for graphene system, and that of electromagnetic(EM) radiation in the medium can hardly be equal. In fact, frequency of PP is generally lower than that of EM radiation. We argue that in the discussed structure there always exists a match point where the frequencies of plasmon-polarition and electromagnetic radiation are the same.

Consider a three-dimensional periodic structure of graphene ribbons embedded in a medium with relative dielectric constant, , as shown in figure 1. For simplification, we assume that each graphene ribbon, parallel to the planar defined by z = 0, has infinite length in the y direction and width d in the x direction. The periods in x direction and interface distance along z direction are a and b, respectively.

Zoom In Zoom Out Reset image size

We here suppose that the plasmon propagating in xz plane and only consider the transverse magnetic (TM) mode, that is, the magnetic field has only y component and correspondingly the electric field has only x and z components. We also require the time dependence of the fields in the form e^{−i ω t} and the absence of free volume currents and charges for simplification. The Maxwell equations read as

$$\begin{eqnarray}\begin{array}{rcl}{\partial }_{z}{E}_{x}-{\partial }_{x}{E}_{z} & = & i\omega {B}_{y}\\ {\partial }_{z}{B}_{y} & = & \displaystyle \frac{i\omega \epsilon }{{c}^{2}}{E}_{x}\\ {\partial }_{x}{B}_{y} & = & -\displaystyle \frac{i\omega \epsilon }{{c}^{2}}{E}_{z}\end{array}\end{eqnarray} \tag{ 1 }$$

Notice that all the field components are invariant along y direction, that is, ∂_y E_x  = ∂_y E_z  = ∂_y B_y  = 0.

From above equations one can obtain the Helmholtz equations for each field components. For B_y , for example, the Helmoholtz equation is,

$$\begin{eqnarray}&&{\partial }_{{xx}}^{2}{B}_{y}+{\partial }_{{zz}}^{2}{B}_{y}+{k}^{2}{B}_{y}=0,\end{eqnarray} \tag{ 2 }$$

where ${k}^{2}=\tfrac{{\omega }^{2}\epsilon }{{c}^{2}}$ is the wave number of electromagnetic (EM) radiation with angle frequencyω in the medium. Equations for E_x and E_z are similar.

Since the structure is periodic in z direction and there is no volume current and charge, the boundary conditions are

$$\begin{eqnarray}\begin{array}{rcl}{E}_{x}(z=0) & = & {E}_{x}(z=b)\\ {B}_{y}(z=0)-{B}_{y}(z=b) & = & -{\mu }_{0}\sigma {E}_{x}(z=0),\end{array}\end{eqnarray} \tag{ 3 }$$

where μ₀ is the magnetic permeability and σ the conductivity. One can write them in the form of B_y ,

$$\begin{eqnarray}\begin{array}{rcl} & & {\partial }_{z}{B}_{y}(z=0)={\partial }_{z}{B}_{y}(z=b),\\ & & {B}_{y}(z=0)-{B}_{y}(z=b)=\displaystyle \frac{i\sigma (x,\omega )}{\omega \epsilon {\epsilon }_{0}}{\partial }_{z}{B}_{y}(z=0),\end{array}\end{eqnarray} \tag{ 4 }$$

where ₀ is the vacuum electric permittivity, and σ(x, ω) is the conductivity in the interface.

Now we assume the wave vector of plasmon is (q_x , 0, q_z ). From the Fourier-Floquet decomposition, B_y can be written as

$$\begin{eqnarray}&&{B}_{y}={e}^{i({q}_{x}x+{q}_{z}z)}G(x,z),\end{eqnarray} \tag{ 5 }$$

where G(x, z) is a quasiperiodic function. In other words, G(x + a, z) = G(x, z) and G(x, z + b) = G(x, z) almost everywhere except at the interfaces, where G(x, 0) and G(x, b) are not equal but determined by boundary conditions. With the preparations, G(x, z) has the form $G(x,z)={\sum }_{m=-\infty }^{\infty }{X}_{m}{Z}_{m}(z)$, where ${X}_{m}={e}^{i\tfrac{2\pi m}{a}x}$.

Instituting above arguments into equation (2) one found that

$$\begin{eqnarray}&&\left[-{\left({q}_{x}+\displaystyle \frac{2\pi m}{a}\right)}^{2}+{k}^{2}+{\left({{iq}}_{z}+{\partial }_{z}\right)}^{2}\right]{Z}_{m}(z)=0.\end{eqnarray} \tag{ 6 }$$

To study the boundary condition, we first deal with the conductivity. In the first period, 0 < x < a, the conductivity can be written as σ(x, ω) = σ(ω)f(x)θ(d − x), where θ(x) denotes the Heaviside step function, σ(ω) is determined by different processes, such as intraband and interband contributions [6, 21], and f(x) is determined by charge contribution in the ribbon. Due to the periodic arrangement in x direction, we have

$$\begin{eqnarray}&&\sigma (x,\omega )=\displaystyle \sum _{m=-\infty }^{\infty }{\sigma }_{m}{X}_{m},\end{eqnarray} \tag{ 7 }$$

where ${\sigma }_{m}=\tfrac{\sigma (\omega )}{a}{\int }_{0}^{d}f(x){X}_{m}^{* }{dx}$. A simply choice f(x) = 1, leaded to ${\sigma }_{m}=\tfrac{i\sigma (\omega )}{2m\pi }(1-{e}^{2{im}\pi d/a})$. In the [22], authors used a spatial modulation of conductivity. But their coefficients are twice as the ones obtained here except that the two coefficients at the case of m = 0are the same, σ(ω)d/a. However, such choice does not satisfy the edge condition of each ribbon and furthermore, it takes the coupling between different m'th mode, which will complicate our analysis and make the calculation unstable. The validity of the edge condition [3, 23] assumes that $f(x)=\chi \sqrt{{d}^{2}/4-{\left(x-d/2\right)}^{2}}$ where χ was determined by a lengthy function and also should be adjusted to preserve the stability in calculations. Here, to make calculations stable and simple, and to compare our results with those of monolayer graphene grid, we simply choose ${\sigma }_{m}=\sigma (\omega )\tfrac{d}{a}{\delta }_{m0}$. In the long wavelength approximation, we do not concern the detailed distribution of electromagnetic field and the choice is suitable.

Now boundary conditions change to

$$\begin{eqnarray}&&({{iq}}_{z}+{\partial }_{z}){Z}_{m}(0)={e}^{{{iq}}_{z}b}({{iq}}_{z}+{\partial }_{z}){Z}_{m}(b),\end{eqnarray} \tag{ 8 }$$

$$\begin{eqnarray}&&{Z}_{m}(0)-{e}^{{{iq}}_{z}b}{Z}_{m}(b)=\displaystyle \sum _{l}\displaystyle \frac{i{\sigma }_{m-l}}{\omega \epsilon {\epsilon }_{0}}({{iq}}_{z}+{\partial }_{z}){Z}_{l}(0).\end{eqnarray} \tag{ 9 }$$

From equation (6), in the first period, one has,

$$\begin{eqnarray}&&{Z}_{m}(z)={c}_{1m}{e}^{-({p}_{m}+{{iq}}_{z})z}+{c}_{2m}{e}^{({p}_{m}-{{iq}}_{z})z},\end{eqnarray} \tag{ 10 }$$

where ${p}_{m}=\sqrt{{\left({q}_{x}+\tfrac{2\pi m}{a}\right)}^{2}-{k}^{2}}$ and c_1m and c_2m are constants to be determined. From equation (8), we have ${c}_{2m}=-{c}_{1m}{e}^{-{p}_{m}b}$. Therefore, we wrote Z_m (z) in the first period as

$$\begin{eqnarray}&&{Z}_{m}(z)={c}_{1m}{e}^{-{{iq}}_{z}z}({e}^{-{p}_{m}z}-{e}^{{p}_{m}(z-b)}).\end{eqnarray} \tag{ 11 }$$

Above expression guarantees the evanescence of B_y . In fact, we have Z_m (b/2) ≡ 0 and Z_m (0) ≡ − Z_m (b) for each m.

Then, we have from equation (9),

$$\begin{eqnarray}&&\displaystyle \sum _{l=-\infty }^{\infty }{A}_{{ml}}{c}_{1l}=0,\end{eqnarray} \tag{ 12 }$$

where ${A}_{{ml}}={\delta }_{{ml}}+\tfrac{i{\sigma }_{m-l}{p}_{l}}{2\omega \epsilon {\epsilon }_{0}}\tfrac{1+{e}^{-{{bp}}_{l}}}{1-{e}^{-{{bp}}_{m}}}$. The requirement that there are nonzero solutions to c_1m means that the  ∞ × ∞  matrix A satisfies $\det A=0$. Therefore, $\det A=0$ determines the dispersion relation of PP in the periodic structure of graphene ribbons.

Since matrix A is independent on q_z , the plasmon dispersion only depends on q_x . In other words, plasmon frequencies at different q_z are the same, as long as q_x is at the same. We thus can generate not only PP propagating in x direction, but also degenerate one propagating in xz plane. The interesting result is favourable to manipulate plasmon in the structure.

It seems that the phase factors ${e}^{{{iq}}_{z}z}$ in equation (5) and ${e}^{-{{iq}}_{z}z}$ in equation (11) should be counteracted by each other. Then, there is a puzzle that what role q_z plays. In fact, the factor ${e}^{{{iq}}_{z}z}$ in equation (5) plays its role in the whole z field while ${e}^{-{{iq}}_{z}z}$ in equation (11) only in the region 0 < z < b, the first period of z. In other word, only in each one period of z, these two factors counteract each other. The phase difference between points (x,z) and (x, z + nb), for instance, is nq_z b, where n is an integer. Therefore, q_z supplies a mechanism of plasmon waveform modulation. Since q_z does not enter in the dispersion relation of plasmon, we simplify assume q_z = 0 in the following.

Here we consider that the environment temperature T is zero. As for the graphene conductivity σ(ω), it is determined by doping concentrate or Fermi energy (${E}_{F}={\hslash }{v}_{F}\sqrt{\pi {n}_{e}}$, where v_F  ∼ 1 × 10⁶ m s⁻¹ is the Fermi velocity of carriers and n_e is electronic density). If a ∼ 10 μm, E_F  ∼ 0.1 eV and we focus on the THz spectral range (ℏ ω ∼ 1 meV), the dominant contribution to conductivity is Drude term [6] and the conductivity reads as

$$\begin{eqnarray}&&\sigma (\omega )=i{\sigma }_{0}\displaystyle \frac{4{E}_{F}}{\pi {\hslash }\omega },\end{eqnarray} \tag{ 13 }$$

where σ₀ = e²/(4ℏ).

To make conclusion more general, we introduce dimensionless quantities, b₀ = b/a, q₀ = q_x a, d₀ = d/a, ω₀ = ω a/c, ${k}_{0}={ka}={\omega }_{0}\sqrt{\epsilon }$, ${p}_{0l}={p}_{l}a=\sqrt{{\left({q}_{0}+2\pi l\right)}^{2}-{\omega }_{0}^{2}\epsilon }$ where c is the vacuum light speed. Then, we have

$$\begin{eqnarray}&&i\displaystyle \frac{\sigma (\omega ){p}_{l}}{2\omega {\epsilon }_{0}\epsilon }=-{g}_{0}\displaystyle \frac{\sqrt{N}}{{\omega }_{0}^{2}\epsilon }{p}_{0l}\end{eqnarray} \tag{ 14 }$$

where we have used ${g}_{0}\equiv \tfrac{{e}^{2}{v}_{F}}{2\sqrt{\pi }{\hslash }{\epsilon }_{0}{c}^{2}}=8.63\times {10}^{-5}$ and N = n_e a². We notice that the validity of Drude expression, ℏ ω/E_F  ≪ 1, means ${\omega }_{0}\ll \tfrac{\sqrt{\pi }{v}_{F}}{c}\sqrt{N}\simeq 5.91\times {10}^{-3}\sqrt{N}$.

At last, with the notion ${f}_{{ml}}=\tfrac{1+{e}^{-{p}_{0l}{b}_{0}}}{1-{e}^{-{p}_{0m}{b}_{0}}}$, we have,

$$\begin{eqnarray}&&{A}_{{ml}}={\omega }_{0}^{2}\epsilon {\delta }_{{ml}}-{g}_{0}\sqrt{N}{p}_{0l}{d}_{0}{f}_{{ml}}{\delta }_{{ml}}\end{eqnarray} \tag{ 15 }$$

which determines PP dispersion by $\det A=0$. Notice that if d₀ = 1 and p_0m b₀ ≫ 1 (leading to f_ml  → 1), one shall obtain the famous result of surface PP where there is only one interface [3, 6] from m = l = 0. We thus nominate the fraction f_ml as the correction factor due to the periodic arrangement along z direction. If d₀ ≠ 1, frequency of PP with only one interface, ${\omega }_{0}^{{\prime} }$, reads as

$$\begin{eqnarray}&&{\omega }_{0}^{{\prime} 2}{\epsilon }_{r}=\displaystyle \frac{{g}_{0}^{2}{d}_{0}^{2}N}{2}\left(\sqrt{1+\displaystyle \frac{4{q}_{0}^{2}}{{g}_{0}^{2}{d}_{0}^{2}N}}-1\right).\end{eqnarray} \tag{ 16 }$$

Our choice of σ_m decouples PP mode with different m and makes the mode with m = 0 the ground one in the first Brillouin zone, ∣q₀∣ < π. Such choice therefore simplify our discussions. In the following we focus on the correction to the ground states due to f_ml .

For the structure the carrier concentration is always fixed and not tunable after the device was fabricated. However, the theoretical study on the effect of different concentration or equivalent N is meaningful. Results of N effect were listed in figure 2. In the figure, ${\omega }^{{\prime} }$ was given by equation (16). We found that, as expected, PP frequency gets larger as N and/or q₀ increase . However, there is a coupling between different interface. The smaller the q₀, the stronger the coupling. Therefore, ω₀ becomes more smoother than ${\omega }_{0}^{{\prime} }$ and PP spectrum is lifted. The lifting makes that ω₀ does not tend to zero at q₀ → 0 although ${\omega }_{0}^{{\prime} }$ does. p₀₀ is imagine number and Z₀ behaves as combination of cosine functions, ${e}^{-i\sqrt{{\omega }_{0}^{2}{\epsilon }_{r}-{q}_{0}^{2}}z}-{e}^{i\sqrt{{\omega }_{0}^{2}{\epsilon }_{r}-{q}_{0}^{2}}(z-b)}$. In fact, ${\omega }_{0}^{{\prime} }$ behaves as convex function but ω₀ as concave one.

Zoom In Zoom Out Reset image size

To study the effect of the coupling between different interfaces we also listed our results of PP dispersion at different b₀ in figure 3. As expected, the smaller the b₀ is, the more obvious the lifting becomes. On the other hand, if b₀ ≫ 1, the lifting at larger q₀ is neglected and ${\omega }_{0}\to {\omega }_{0}^{{\prime} }$, as shown by the case of b₀ = 4 in the figure. However, for small q₀ (or larger wavelength), the lifting can not be neglected and we have still no ω₀ → 0 when q₀ → 0.

Zoom In Zoom Out Reset image size

To study the effect of width of ribbon, d₀, we listed PP dispersion at different d₀ in figure 4. We found that for a fixing q₀, ω₀ becomes smaller with decreasing d₀. The phenomenon can be understood from equation (15). The equation told us that a smaller d₀ is equivalent to a smaller N, which leading to a smaller ω₀.

Zoom In Zoom Out Reset image size

In the section we discussed the PP spectrum in the 3D periodic graphene structure at different ribbon width, interface distance, period in the interface and carrier concentrate. The influences of different parameter were given.

Because of the correction factors, there is always a lifting in the PP spectrum. The smaller q₀ becomes, the more obvious the lifting gets. When q₀ → 0, we have ${\omega }_{0}\to \tfrac{{k}_{0}}{\sqrt{\epsilon }}\ne 0$, where k₀ is determined by ${k}_{0}\tan \tfrac{{k}_{0}{b}_{0}}{2}={g}_{0}\sqrt{N}{d}_{0}$ and becomes larger with increasing N, d₀ and decreasing b₀, as shown in figures 2, 3 and 4. where EM radiation spectra in the medium were also shown. It was found that due to the lifting, there exists an intersection between PP dispersion and EM radiation dispersion, where the two dispersions are the same and a match between PP and EM radiation occurs. To estimate the matching point, we notice that the match always occurs at long wave region, where the frequency of PP hardly changes, therefore, match occurs roughly at the wave vector k₀. If $\tfrac{2{g}_{0}\sqrt{N}{d}_{0}}{{b}_{0}}$ is not large, ${k}_{0}\simeq \sqrt{\tfrac{2{g}_{0}\sqrt{N}{d}_{0}}{{b}_{0}}}$. Therefore, one can easily excite PP in the structure by EM radiation at the match point.

In conclusion, to improve the efficiency of photodetectors, we first designed a three-dimensional periodic structure of graphene ribbons embedded in a medium and investigated the dispersion of TM plasmon-polaritons in the structure. It was shown that for a plasmon-polariton propagating in xz-plane, the wave vector of which is, for instance, (q_x , q_z ), the component of wave vector perpendicular to graphene ribbons, q_z , is irrelevant to the PP spectrum but supplies a wave function modulation.

We further considered the spectra of plasmon-polaritons in the structure and compared results with those in monolayer graphene interface embedded in the medium. There are sharp differences between these results. Firstly, wave functions are not exponential decaying but combinations of cosine functions due to the periodic boundary conditions. Secondly, the spectra are lifted and behave as concave functions but not convex ones with respect to q_x . The larger the carrier concentration, the more obvious the lifting. Thirdly, although for the case q_x b ≫ 1, where b is the interface distance, the two results approximately equal, for the case of q_x  → 0, the frequencies of plasmon-polaritons in our structure do not tend to zero. The matching condition, where the frequencies of PP and EM radiation are equal, occurs at the wave vector k₀. The match makes it easier to excite plasmon polariton with electromagnetic radiations in the three-dimensional period structure than that in a two-dimensional period structure. Thus, the three-dimensional structure supplies a high efficiency of photoelectric detection.

The National Natural Science Foundation of China (21773181, 21 573 172) has financially supported this work.