A sheet of graphene under magnetic bias attains anisotropic surface conductivity, opening the door for realizing compact devices such as Faraday rotators, isolators and circulators. In this paper, an accurate and analytical method is proposed for a periodic array of graphene ribbons under magnetic bias. The method is based on integral equations governing the induced surface currents on the coplanar array of graphene ribbons. For subwavelength size ribbons subjected to an incident plane wave, the current distribution is derived leading to analytical expressions for the reflection/transmission coefficients. The results obtained are in excellent agreement with full-wave simulations and predict resonant spectral effects that cannot be accounted for by existing semi-analytical methods. Finally, we extract an analytical, closed form solution for the Faraday rotation of magnetically-biased graphene ribbons. In contrast to previous studies, this paper presents a fast, precise and reliable technique for analyzing magnetically-biased array of graphene ribbons, which are one of the most popular graphene-based structures.
© 2019 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
Graphene is a two-dimensional (2D) nanomaterial that is composed of a single layer of carbon atoms organized into a hexagonal lattice . In the past decade, graphene-based devices have become a prime focus of research due to their unique electrical, optical, thermal, and mechanical properties ,,. The electrical conductivity of graphene can change by an external DC voltage or chemical doping. This tunability provides additional flexibility in the design of diverse functionalities. Therefore, numerous graphene-based applications have been proposed from microwave to terahertz regions [5–14]. Under an external magnetic bias, the surface conductivity of graphene becomes a tensor, and it will also have non-reciprocal properties . These extraordinary properties make graphene a promising material for many novel applications, such as [16, 17] , nonreciprocal couplers , optical modulators , and Faraday rotators , in the terahertz regime.
Due to strong interaction of graphene with electromagnetic (EM) fields, patterned magnetically-biased graphene has been assessed up to now for several applications [21–23]. For example, in  it has been shown that the Faraday rotation of graphene micro ribbons can be controlled dynamically and that a giant Faraday rotation can be obtained, the frequency of which corresponds to the width of ribbons. But a closed-form relationship between the width of the graphene ribbons and the resonance frequency was not presented.
Consequently, the analytical and numerical study of continuous and patterned anisotropic graphene is currently of high interest. A continuous graphene sheet can be easily analyzed by applying appropriate boundary conditions. Sounas and Caloz studied the transmission and reflection properties of magnetically-biased graphene sheet using the general anisotropic conductivity tensor . To the best of the authors’ knowledge, no rigorous and analytical solution has yet been presented on patterned anisotropic graphene under magnetic bias. Although certain numerical methods, including the Fourier modal-based , the finite difference time-domain method , and discontinuous Galerkin time domain (DGTD)  have been used to numerically solve this problem, they suffer from a slow convergence rate and more importantly, are not analytic results. Furthermore, some efforts have been made to decrease the computation time [27, 28]; however, they still require higher acceleration. In addition, a simple method for the analysis of magnetic-biased graphene ribbons based on the effective medium approach was proposed in ; however, as it will be shown, this method does not have accuracy and is only valid for the extreme subwavelength regime. Hence a rigorous and analytical analysis of magnetically-biased Graphene-based meta-surfaces is of paramount importance.
In this paper, a rigorous analysis is performed through solving the integral equations governing the induced surface currents on the graphene ribbons to accurately predict the electromagnetic performance of an array of magnetically-biased graphene ribbons. The paper is organized as follows: In Section (2) we shall first study the scattering of a TE/TM electromagnetic wave by a single magnetically-biased graphene ribbon. Then an analytical expression is extracted for the surface current density induced within a periodic array of magnetically-biased ribbons in the subwavelength regime. In Section (3), the reflection/transmission coefficient of the structure is calculated and the proposed method is validated against full-wave simulations. We show how plasmonic resonances across the graphene ribbons enter the reflection and transmission coefficients of such a structure. The limitations of our proposed method is discussed. Moreover, an analytical expression is presented for Faraday rotation of magnetically-biased graphene ribbons, as one of the important phenomena in optics and photonics.
2. Surface currents on a perodic array of anisotropic graphene ribbons
We shall consider a periodic array of 1D graphene ribbons under magnetic bias as shown in Fig. 1. Each ribbon has a width of w. The periodicity of the array is D in the x-direction. The graphene ribbons are placed on the x − y plane (infinite along y) and a static magnetic field B0 is applied perpendicular to the ribbons. A normal EM plane wave is incident on the arrangement. Due to the incident wave, a surface current density with components Jx and Jy will be induced onthe ribbons. Since there is no variation in y-direction, these induced currents are functions of x only (i.e. ).
2.1. Electromagnetic parameters of magnetically-biased graphene
As noted earlier, graphene is a 2-D nanomaterial consisting of a single layer of carbon atoms organized within a hexagonal lattice. Assume that the graphene sheet is on the x − y plane and that a magnetic field is also established perpendicular to the graphene sheet (). The surface conductivity of this graphene sheet is represented by a tensor () that is obtained from microscopic, semi-classical, and quantum mechanical considerations as (1) can be simplified by a Drude-like model as  24]. Carrier mobility is assumed to be within the range of to , which is practical and achievable for graphene with different substrates at room temperature [31, 32]. It should also be noted that the time harmonic of the form is assumed throughout this paper.
2.2. Scattering from a single ribbon
Let us first consider a single anisotropic graphene ribbon in free space. The integral equation governing the induced surface currents on the single ribbon can be expressed as 
Assuming the strip to be narrow (), we next adopt the quasi-static approximation . As a result, the second term on the right hand side of Eq. (3) becomes negligible compared to the first term. Moreover, the Green’s function Eq. (4) can be approximated by
It must be noted that σ0 is the graphene conductivity in the absence of the external magnetic field.
Assume, for the moment, that is known. Equation (6) can then be viewed as can integral equation for the function . This integral equation has a complicated, but exact solution. However, in the limit whereEq. (6) becomes negligible and the solution is trivially given by
To see why Eq. (8) is satisfied, note thatEq. (8) is satisfied.
Let us now return to Eq. (3) and combine it with Eq. (3) to arrive at an equation for . However, one may simplify the resulting equation by neglecting the second term on the right hand side of Eq. (3). This term is small compared to the first term for a narrow ribbon. The integral equation governing Jx thus becomes
2.3. Scattering from an infinite array of ribbons
Now consider a periodic array of coplanar graphene ribbons with period D as shown in Fig. 1. The ribbon array is subject to a plane electromagnetic wave of arbitrary polarization that is normally incident on the array plane. Equation (3) that locally relates the ribbon surface currents to the total electric field is still applicable. However, on each ribbon, must now account for the total electric field that is produced by all ribbons. Fortunately, in case of normal incidence on an infinite number of array elements, the distributions of surface currents on all ribbons are identical. As a result, on each individual ribbon, Eq. (3) and Eq. (3) remain provided that the Green’s function Eq. (4) is replaced by
As in the case of a single ribbon, we first treat the equation for Jy. However, the interaction between ribbons, manifested in the periodic Green’s function Eq. (12), leads to an important modification of the approximation Eq. (5). For sub-wavelength arrays where , a lengthy calculation presented in the appendix shows that
The interaction between sub-wavelength ribbons leads to a Green’s function whose leading term is of the order . Unlike the case of a single ribbon, the second term in the integral equation Eq. (6) cannot be neglected anymore. By replacing G0 with Eq. (13) in Eq. (6), and keeping the leading term alone, one ends up with the equation
Note that, since the external field is now a normally incident plane wave, , do not depend on x. This simple equation is solved by
This equation is identical to the one describing an array of graphene ribbons in the absence of a magnetic field , except for F replacing the external field . Using perturbation theory, the solution of Eq. (18) can be written as
Note that in Eq. (22) equals the admittance of a capacitance in series with a conductance .
The constant F self-consistently depends on Jx through Ix [see Eqs. (16) and (19)]. To calculate Ix, an integration is carried out over Eq. (21), and Eq. (19) is used. After basic calculations one obtains
Substitution in Eq. (19) then yields
In all the aforementioned equations, the graphene ribbons are assumed to be suspended in free space. In a different medium with a relative dielectric constant of εr, all the expressions obtains remain valid provided that ε0 and k0 are replaced by and , respectively.
These calculations can be readily generalized to the cases of oblique incidence. This generalization can be obtained by the correction of . Consider an array of graphene ribbons when illuminated by an oblique incidence angle (θi). Equation (13) can be re-expressed as
2.4. Scattering parameters
The distribution of surface current on graphene ribbons is next used to calculate the reflection and transmission coefficients of the structure. To that end, the electric field in the region is expressed as the Rayleigh expansion
Similarly, the electric field in the region is written as
In these relationships
In the sub-wavelength regime where only the zero-order mode (n = 0) is propagating. All other modes are evanescent. For sub-wavelength ribbon arrays the reflection and transmission coefficient are, therefore, defined using the zero-order mode amplitudes alone, according to the relationships
In order to compute the reflection and transmission matrices, the usual boundary conditions on electric and magnetic fields is applied at z = 0, and an integration is carried out over the resulting equations on a single period, i.e., from x = 0 to . It then follows that
The next step is to express in terms of the incident field. To that end, we carry out an integration over Eq. (15) from to . This results inEq. (7). Combination of Eqs. (26),(28), (42), and (39)–(37) then yields
For simplification, in all the equations of this subsection, we assume that graphene is suspended in free space. These equations can be generalized to the case where the ribbon is sandwiched between two homogenous media with permittivities ε1 and ε2, respectively as shown in Fig. 1. However, this method has complicated calculations. In this case, we employ another simple method to solve the integral equations of Eq. (14) and Eq. (18). In this method we omit the term but account for its effect by including the reflected fields in the external fields [34, 36]. and now not only include the true external field but also the field reflected by the ribbon array. For example, for an incident TM wave with unity amplitude, the total amplitude of the electric field must be replaced by and Ryx in the x and y directions, respectively.
3. Result, discussion and application
3.1. Numerical results
In this section we validate and verify the accuracy of the proposed method through some examples. Consider a periodic array of graphene ribbons in free space with for the first example. For a normally incident plane wave, the magnitude (absolute value) of reflection coefficients Eqs. (43), (44) and (45) is plotted as function of frequency in Fig. 2. For the sake of comparison the results obtained using effective medium theory  (dotted blue), and finite integration technique (FIT) (dashed black) are also plotted. The FIT results are obtained using CST Microwave Studio 2018. Excellent agreement can be observed between the FIT results and our analytical method. By contrast, the effective medium approach  can not reproduce the electromagnetic response of magnetically-biased graphene ribbons.
The resonance features observed in Fig. 2 can be explained as follows. Equations (21) and (22) imply that when becomes very small for some n, becomes largeand a resonance occurs. Near the corresponding resonance frequency, the distribution of is predominantly given by . Apart from constant terms, the distribution of is identical to as can be seen from Eq. (15). Figure 3 shows the magnitude of across a graphene ribbon at the first (n = 1) and second (n = 3) resonances, that occur at 9.78- and 19.13 THz, respectively. It must be noted that is an odd function of x with respect to x = 0 (center of a ribbon) for even values of n. As a result Sn = 0 for even values of n [see Eq. (23)] so that only odd values of n enter the summation Eq. (21). The resonances observed are thus due to odd modes . Near these resonance frequencies, peaks in reflection coefficients are observed. Note, however, that Ryy does not show any resonance features in the absence of an external magnetic bias, i.e., when . This can be seen from Eq. (45).
As a second example, consider a periodic array of graphene ribbons that are illuminated by an oblique incident wave. The periodic graphene ribbon structure has , , , . The result of the magnitude and phase and of Ryx vs frequency is presented in Fig. 4. For brevity, only Ryx has been plotted here. A very good agreement is seen between the analytical method and the FIT. It can be observed from Fig. 4 that oblique incidence does not change the resonance frequency, and this point is clear from Eqs. (21),(22),(32) and (33). Oblique incidence does not change the eigenfunction and eigenvalue of the magnetically-biased graphene ribbons, and as previously mentioned each resonance frequency corresponds to an eigenfunction and eigenvalue. It can be further seen from the figure that if the angle of incidence increases, the magnitude of Ryx can increase Eq. (44).
With regards to the limitations of the presented theory, it should be noted that our analysis assumed operation in the sub-wavelength regime where . As shown in Fig. 2, when the frequency increases, the accuracy of the proposed method decreases. For the position of the first resonance, we have a relative error of while the relative error for the fourth resonance is . Moreover, solution Eq. (21) is based on perturbation theory in which the interaction between neighboring ribbons is assumed not to affect the eigenvalue Eq. (24) . Increasing the fill factor (), will lead to a smaller gap and a stronger interaction between neighboring ribbons. One would, therefore, expect the accuracy of our method to decrease. The frequency of the first resonance is plotted as function of fill factor in Fig. 5 for . It is observed that by increase in the fill factor, the relative error in predicting the resonance frequency is slightly increased. Finally, derivation of Eq. (9) was based on Eq. (8). According to Eqs. (2) and (7), when Ef is large, Eq. (8) cannot be satisfied. If the Fermi level energy increases, the accuracy of the proposed method can slightly decrease. The effect of this approximation, as well as perturbation theory, is small compared to the sub-wavelength approximation. The latter, however, does not pose a strict limitation, since anisotropic graphene ribbons have almost always been utilized in the sub-wavelength regime [21–23, 29]. Consequently, the proposed analysis is much more accurate and affordable than its counterparts already reported in the literature [23, 28, 29].
3.2. Design of a Faraday rotator
Faraday rotation is one of the most attractive phenomena in optics, with many applications in optical diodes, sensing, magnetic microscopy, optical communications, data storage, optical isolators, and phase modulators [20, 23]. For the array of graphene ribbons under magnetic bias, Faraday rotation angle can be directly computed as follows :Eqs. (43) and (44). Large Faraday rotation angles are observed at resonance. To achieve a giant Faraday rotation at 10 THz, the structure parameters are designed as . The results are shown in Fig. 6. A good agreement between the results of our method and numerical FIT is observed. Also, we compare our results with the semi-analytical results . It is clear that our method yields much more accurate results for the Faraday rotation angle of the magnetically-biased graphene ribbons. This is because the analysis in  is based on the low-relaxation-time assumption. Besides, it assumes that period D is much smaller than the free-space wavelength. Neither condition is satisfied here. Finally, we investigate the effect of magnetic bias on performance of the Faraday rotator in Fig. 6(b). It can be seen that increasing the magnetic bias leads to larger Faraday rotation angles. Also, giant Faraday rotation occurs at higher frequencies.
Based on integral equations governing the field and the surface current density on magnetically-biased graphene ribbons, this paper proposed a analytical, fast, and accurate method for the analysis of such structures. Unlike previous studies, the proposed method does not suffer from extremely slow convergence rates of numerical methods, and more importantly is able to predict spectral resonance features, not possible with previous analytic solutions. The expressions were derived under a quasi-static approximation and extended to periodic arrays of graphene ribbons using perturbation theory. By Rayleigh expansion, we determined the reflection/transmission coefficient using a rigorous method. The proposed method was validated against commercial software equipped with finite integration technique (FIT), showing excellent agreement between the two solutions. In addition, a Faraday rotator was designed and investigated as an important application. The present study suggests a new method for the analysis of magnetically-biased graphene ribbons, which is more precise, reliable, and affordable than the methods already reported in the literature. The proposed method might also be useful for the analysis and development of different devices such as polarizer, isolators, and so on.
Consider the Green function Eq. (12). The Hankel function is written as37] in the complex θ plane. Since , one has:
Equation (48) can be rewritten as
Since is small (and so is ) the integrand in the first integral can be expanded in powers of ,. The result can be easily shown to be of the order . To compute the second integral it is rewritten asEq. (49) leads to Eq. (13). Take note that Eq. (20).
1. A. K. Geim and K. S. Novoselov, “The rise of graphene,” in Nanoscience and Technology: A Collection of Reviews from Nature Journals, (World Scientific, 2010), pp. 11–19.
2. K. S. Novoselov, A. K. Geim, S. V. Morozov, D. Jiang, Y. Zhang, S. V. Dubonos, I. V. Grigorieva, and A. A. Firsov, “Electric field effect in atomically thin carbon films,” science 306, 666–669 (2004). [CrossRef] [PubMed]
5. M. Rahmanzadeh, H. Rajabalipanah, and A. Abdolali, “Multilayer graphene-based metasurfaces: robust design method for extremely broadband, wide-angle, and polarization-insensitive terahertz absorbers,” Appl. optics 57, 959–968 (2018). [CrossRef]
6. A. Momeni, K. Rouhi, H. Rajabalipanah, and A. Abdolali, “An information theory-inspired strategy for design of re-programmable encrypted graphene-based coding metasurfaces at terahertz frequencies,” Sci. reports 8, 6200(2018). [CrossRef]
7. P.-Y. Chen, J. Soric, Y. R. Padooru, H. M. Bernety, A. B. Yakovlev, and A. Alù, “Nanostructured graphene metasurface for tunable terahertz cloaking,” New J. Phys. 15, 123029 (2013). [CrossRef]
8. M. Rahmanzadeh, H. Rajabalipanah, and A. Abdolali, “Analytical investigation of ultrabroadband plasma–graphene radar absorbing structures,” IEEE Transactions on Plasma Sci. 45, 945–954 (2017). [CrossRef]
10. S. AbdollahRamezani, K. Arik, A. Khavasi, and Z. Kavehvash, “Analog computing using graphene-based metalines,” Opt. letters 40, 5239–5242 (2015). [CrossRef]
11. R. Filter, M. Farhat, M. Steglich, R. Alaee, C. Rockstuhl, and F. Lederer, “Tunable graphene antennas for selective enhancement of thz-emission,” Opt. express 21, 3737–3745 (2013). [CrossRef] [PubMed]
12. K. Rouhi, H. Rajabalipanah, and A. Abdolali, “Multi-bit graphene-based bias-encoded metasurfaces for real-time terahertz wavefront shaping: From controllable orbital angular momentum generation toward arbitrary beam tailoring,” Carbon 149, 125–138 (2019). [CrossRef]
13. B. Sensale-Rodriguez, S. Rafique, R. Yan, M. Zhu, V. Protasenko, D. Jena, L. Liu, and H. G. Xing, “Terahertz imaging employing graphene modulator arrays,” Opt. express 21, 2324–2330 (2013). [CrossRef] [PubMed]
15. G. W. Hanson, “Dyadic green’s functions for an anisotropic, non-local model of biased graphene,” IEEE Transactions on Antennas Propag. 56, 747–757 (2008). [CrossRef]
16. M. Tamagnone, C. Moldovan, J.-M. Poumirol, A. B. Kuzmenko, A. M. Ionescu, J. R. Mosig, and J. Perruisseau-Carrier, “Near optimal graphene terahertz non-reciprocal isolator,” Nat. communications 7, 11216 (2016). [CrossRef]
17. X. Lin, Z. Wang, F. Gao, B. Zhang, and H. Chen, “Atomically thin nonreciprocal optical isolation,” Sci. reports 4, 4190 (2014). [CrossRef]
19. M. Tamagnone, A. Fallahi, J. R. Mosig, and J. Perruisseau-Carrier, “Fundamental limits and near-optimal design of graphene modulators and non-reciprocal devices,” Nat. photonics 8, 556 (2014). [CrossRef]
20. I. Crassee, J. Levallois, A. L. Walter, M. Ostler, A. Bostwick, E. Rotenberg, T. Seyller, D. Van Der Marel, and A. B. Kuzmenko, “Giant faraday rotation in single-and multilayer graphene,” Nat. Phys. 7, 48 (2011). [CrossRef]
21. J.-M. Poumirol, P. Q. Liu, T. M. Slipchenko, A. Y. Nikitin, L. Martin-Moreno, J. Faist, and A. B. Kuzmenko, “Electrically controlled terahertz magneto-optical phenomena in continuous and patterned graphene,” Nat. communications 8, 14626 (2017). [CrossRef]
22. A. Fallahi and J. Perruisseau-Carrier, “Manipulation of giant faraday rotation in graphene metasurfaces,” Appl. Phys. Lett. 101, 231605 (2012). [CrossRef]
24. D. L. Sounas and C. Caloz, “Electromagnetic nonreciprocity and gyrotropy of graphene,” Appl. Phys. Lett. 98, 021911 (2011). [CrossRef]
25. X.-H. Wang, W.-Y. Yin, and Z. Chen, “Matrix exponential fdtd modeling of magnetized graphene sheet,” IEEE Antennas Wirel. Propag. Lett. 12, 1129–1132 (2013). [CrossRef]
26. P. Li and L. J. Jiang, “Modeling of magnetized graphene from microwave to thz range by dgtd with a scalar rbc and an ade,” IEEE Transactions on Antennas Propag. 63, 4458–4467 (2015). [CrossRef]
27. P. K. Khoozani, M. Maddahali, M. Shahabadi, and A. Bakhtafrouz, “Analysis of magnetically biased graphene-based periodic structures using a transmission-line formulation,” JOSA B 33, 2566–2576 (2016). [CrossRef]
28. M. Feizi, V. Nayyeri, and O. M. Ramahi, “Modeling magnetized graphene in the finite-difference time-domain method using an anisotropic surface boundary condition,” IEEE Transactions on Antennas Propag. 66, 233–241 (2018). [CrossRef]
29. J. S. Gomez-Diaz and A. Alù, “Magnetically-biased graphene-based hyperbolic metasurfaces,” in 2016 IEEE International Symposium on Antennas and Propagation (APSURSI), (IEEE, 2016), pp. 359–360.
30. A. Ferreira, J. Viana-Gomes, Y. V. Bludov, V. Pereira, N. Peres, and A. C. Neto, “Faraday effect in graphene enclosed in an optical cavity and the equation of motion method for the study of magneto-optical transport in solids,” Phys. Rev. B 84, 235410 (2011). [CrossRef]
31. L. Ju, B. Geng, J. Horng, C. Girit, M. Martin, Z. Hao, H. A. Bechtel, X. Liang, A. Zettl, Y. R. Shen, et al., “Graphene plasmonics for tunable terahertz metamaterials,” Nat. Nanotechnology 6, 630 (2011). [CrossRef]
32. C. R. Dean, A. F. Young, I. Meric, C. Lee, L. Wang, S. Sorgenfrei, K. Watanabe, T. Taniguchi, P. Kim, K. L. Shepard, et al., “Boron nitride substrates for high-quality graphene electronics,” Nat. Nanotechnology 5, 722 (2010). [CrossRef]
33. C.-T. Tai, Dyadic Green Functions in Electromagnetic Theory (Institute of Electrical & Electronics Engineers (IEEE), 1994).
34. A. Khavasi and B. Rejaei, “Analytical modeling of graphene ribbons as optical circuit elements,” IEEE J. Of Quantum Electron. 50, 397–403 (2014). [CrossRef]
35. M. Rahmanzadeh, A. Abdolali, A. Khavasi, and H. Rajabalipanah, “Adopting image theorem for rigorous analysis of a perfect electric conductor–backed array of graphene ribbons,” JOSA B 35, 1836–1844 (2018). [CrossRef]
36. S. Barzegar-Parizi, B. Rejaei, and A. Khavasi, “Analytical circuit model for periodic arrays of graphene disks,” IEEE J. Quantum Electron. 51, 1–7 (2015). [CrossRef]
37. I. S. Gradshteyn and I. M. Ryzhik, Table of Integrals, Series, and Products (Academic press, 2014).