Valley permitted Klein tunneling and magnetoresistance in ferromagnetic monolayer MoS2

Valley permitted Klein tunneling and magnetoresistance in ferromagnetic monolayer MoS2

Superlattices and Microstructures 86 (2015) 243–249 Contents lists available at ScienceDirect Superlattices and Microstructures journal homepage: ww...

2MB Sizes 0 Downloads 70 Views

Superlattices and Microstructures 86 (2015) 243–249

Contents lists available at ScienceDirect

Superlattices and Microstructures journal homepage: www.elsevier.com/locate/superlattices

Valley permitted Klein tunneling and magnetoresistance in ferromagnetic monolayer MoS2 M. Khezerlou, H. Goudarzi ⇑ Department of Physics, Faculty of Science, Urmia University, P.O. Box: 165, Urmia, Iran

a r t i c l e

i n f o

Article history: Received 26 June 2015 Received in revised form 1 July 2015 Accepted 24 July 2015 Available online 26 July 2015 Keywords: Molybdenum disulfide Klein tunneling Spin–orbit coupling Exchange field

a b s t r a c t Recently, effect of proximity-induced ferromagnetism or superconductivity in two dimensional monolayer molybdenum disulfide on the energy dispersion of relativistic charge carriers has attracted much effort. Inequivalent nondegenerate K and K 0 valleys resulted from strong spin–orbit coupling in monolayer MoS2 give rise to obtain novel behaviors in some physical phenomena, such as Andreev reflection or spin–valley transport properties, comparing with corresponding structure, graphene. Considerable direct band-gap in monolayer MoS2 makes it potentially useful for electronics applications. Motivated by these, we study, in particular, the influence of various combinations of valley and spin indices on charge transmission in the MoS2 -based ferromagnetic/insulator/ferromagnetic junction. It is shown that the tunneling process is suppressed with arbitrary range of chemical potential for all incident electrons. We determine the exact magnitude of Fermi energy by chemical potential. Topological and mass-related terms ða; bÞ are taken into account in the Dirac-like Hamiltonian, and conductance resulted from the transmission probability is investigated in terms of physical parameters of structure. We, further, study the dependence of the magnetoresistance on exchange field for various cases of valley-coupling. Magnetoresistance of system is changed from a positive value to negative for K valley, only. The Klein tunneling is resulted in different exhibitations for different allowed spin and valley indices. As expected, the exact oscillations of transmission coefficient is found in terms of barrier parameter. Ó 2015 Elsevier Ltd. All rights reserved.

1. Introduction From the point of view of Klein paradox [1,2], recently experimentally discovered two-dimensional exotic condensed matter structures, such as graphene [3] and monolayer molybdenum disulfide (ML-MDS) [4–6] may gain more attention, since charge carriers in these structures can exhibit either electron-like or hole-like quasiparticles belonging to two inequivalent valleys (K and K 0 Dirac points) located at the corners of the hexagonal Brillouin zone. In the last few years, experimentally detecting relativistic effects, such as Klein tunneling became feasible [7–9]. For nanoelectronics applications, existence of a band-gap in low-energy band structure of graphene or ML-MDS can play a crucial role, and, of course, strongly can affect Klein tunneling [10]. In this regards, MoS2 as a member of transition-metal dichalcogenides can be shaped into monolayers, which exhibits distinct physical properties from its bulk counterpart [4,11]. It is experimentally shown that the ML-MDS has a direct band-gap in the visible frequency range ( 0:95 eV), and relativistic charge carriers have a Fermi velocity ⇑ Corresponding author. E-mail addresses: [email protected] (M. Khezerlou), [email protected], [email protected] (H. Goudarzi). http://dx.doi.org/10.1016/j.spmi.2015.07.054 0749-6036/Ó 2015 Elsevier Ltd. All rights reserved.

244

M. Khezerlou, H. Goudarzi / Superlattices and Microstructures 86 (2015) 243–249

v F  0:53  106 ms1

and a room-temperature mobility over 200 cm2 =ðV sÞ [6]. Analogous to graphene, there is a valley index s ¼ 1, which is expected to be robust against scattering by smooth deformations and long wavelength phonons due to a large valley splitting. Moreover, regarding spin, valley and electronics transport applications, monolayer MoS2 is distinguished from graphene by two important features: (i) breaking inversion symmetry, which, in the presence of an in-plan electric field, leads to the valley Hall effect, where quasiparticles in different valleys are transported in opposite transverse edges [12], (ii) strong spin– orbit coupling (SOC) resulted from the heavy metal atoms molybdenum disulfide [13]. The second above characteristics may actually provide a realizable way to investigate spin-related applications, so that, ML-MDS in contact with a ferromagnetic section can be considered as a good system to demonstrate the role of strong SOC, which may result in a major distinction in tunneling conductance between spin-up and spin-down carriers giving rise to magnetoresistance. In Refs. [14–19], it is shown that a ferromagnetic property can be achieved in monolayer MoS2 . Also, existence of defects can induce the ferromagnetic behavior to ML-MDS [20,21]. Recently, the physics of spin and valley coupling [22], ferromagnetic/superconductor/fer romagnetic [23], and normal metal/superconducting [24] in the ML-MDS structures have been studied. In this paper, we study effect of an exchange field on transmission coefficient (caused by Klein tunneling) and resulted magnetoresistance in a ferromagnetic/insulator/ferromagnetic ðF=I=FÞ monolayer MoS2 junction, in which we are allowed to choose two ferromagnetic sections with parallel or antiparallel spin polarizations. Furthermore, one can control the spin and valley indices to obtain the permitted magnetoresistance from equivalent or inequivalent valleys related to the different ferromagnetic sections separated by a thin barrier insulator. Xiao et al. [22] have stated that, at the valence band due to the existence of valley-contrasting spin splitting (0.1–0.5 eV) caused by inversion symmetry breaking spin and valley relaxation are suppressed, and flip of each index alone is forbidden. Consequently, it seems that, in contrast to the similar structure of graphene, the large SOC and valley degeneracy breaking in ML-MDS gives rise to appear a distinct magnetoresistance corresponding to spin-up and spin-down carriers. However, tunneling relative to the mentioned above forbidden valley exchange can actually be checked (note that, in our system the different valleys are taken from two distinct sections, which one can consider it as an inter-valley transport). The effective Hamiltonian of ML-MDS is taken with topological term (b) [25], and its effect is investigated on the charge carriers behavior. This paper is organized as follows. Section 2 is devoted to present the proposed model and formalism to obtain the exact form of transmission coefficient for both parallel and antiparallel exchange fields using the normalized Dirac spinors of ferromagnetic and insulator sections. The numerical results of tunneling conductance and corresponding magnetoresistance of structure are presented, and their main characteristics are discussed in Sections 3 and 4. Finally, we close with a brief summary in Section 5. 2. Model of structure and formalism

A typical F=I=F structure on top of a monolayer MoS2 sheet is sketched in Fig. 1. The ferromagnetic regions are extended from x ¼ 1 to x ¼ 0 and from x ¼ L to x ¼ þ1 for all y. A gate voltage covers the central insulator region ð0 < x < LÞ, that the electrostatic potential is taken V 0 in this region, and zero otherwise. The low-energy band structure of ML-MDS can be described by the modified Dirac Hamiltonian. This Hamiltonian in addition to the first order term of momentum for 2D massive fermions, contains the quadratic terms originated from the difference between electron and hole masses [25]. The strong spin–orbit coupling leads to distinct spin splitting at the valence band for different valleys. In the presence of an exchange field h, Dirac-like Hamiltonian is given by:

0

1    h2 k2  D þ 2m0F a2 þ 2b  sh v F h skx  iky A: H¼@    h2 k2  v F h skx þ iky D þ 2kss þ 2m0F a2  2b  sh

ð1Þ

The spin-up and spin-down carriers are labeled by s ¼ 1, and valley index s ¼ 1 denotes the K and K 0 valleys. The bare electron mass is m0 ¼ 0:05  1010 ðeV s2 =m2 Þ, and topological band parameters are evaluated by b ¼ 2:21; a ¼ 0:43. D is the direct gap mentioned in the previous section, k  0:04 eV denotes the spin–orbit coupling. From the Hamiltonian (1), the energy dispersion (relative to the Fermi energy EF tuned by chemical potential) can be obtained as below:

vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi ! u 2 2  2 u  2 k2F a  k b h h 2 2 F   EF ¼ kss þ þ v 2F h kF :  sh  t D  kss þ 2m0 2 2m0 2

ð2Þ

M o S2 Fig. 1. A schematic plot of ferromagnetic/insulator/ferromagnetic molybdenum disulfide junction for which directions of the polarization of two valleys in ferromagnetic regions can be parallel or antiparallel.

M. Khezerlou, H. Goudarzi / Superlattices and Microstructures 86 (2015) 243–249

245

Fig. 2. The excitation energy of inequivalent valleys in the ferromagnetic molybdenum disulfide segment. Red curves indicate spin-up subbands while the green ones stands for spin-down subbands. Due to the exchange field, the splitting of up and down spin subbands is different in two distinct valleys. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

The energy dispersion at the presence of an exchange splitting of Zeeman type is shown for inequivalent Dirac points in Fig. 2. It is notable that in presence of exchange field, spin-splitting of valence band strongly depends on the valley index, so that, as we observe in Fig. 2 this splitting can considerably differ from one valley (s ¼ þ1) to the other (s ¼ 1) relative to the exchange field. The corresponding Fermi wavevector of ferromagnetic MoS2 can be acquired from energy eigenvalue:

0

x1 a2 þ x2 2b þ m0 v 2F 

kF ¼ @

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi112 2 2 ðx1 2b þ x2 a2Þ þ ðm0 v 2F Þ þ 2m0 v 2F ðx1 a2 þ x2 2bÞ A;   h2 a2  b 2 2 2m0 2

ð3Þ

where we define x1 ¼   EF  kss þ sh; x2 ¼ kss þ D. In the insulator region the eigenstates can be found by substituting h ¼ 0. We can introduce a dimensionless barrier strength Z ¼ kI L. This parameter is defined as characteristics of thin barrier when the electrostatic potential and width of barrier are taken to be V 0 ! 1; L ! 0, respectively. The charge and current density of quasiparticles may be defined by nonrelativistic and relativistic terms based on the Lorentz covariant continuity equation. Using the modified Hamiltonian (1), continuity equation results in:

@  y  @ ww þ @t @x



vF



 wy ðsrx Þw þ

 

h @ @ ¼ 0: wy ða þ brz Þ w  wy ða þ brz Þw @x @x 4m0 i

ð4Þ

The wavefunctions of the reflected and transmitted quasiparticles in the left and right ferromagnetic sections can be described by:

    1 1 1 r iðskFx xþky yÞ p ffiffiffiffiffiffi e eiðskFx xþky yÞ ; wFL ¼ pffiffiffiffiffiffi þ N seish AF N seish AF   1 t eiðskFx xþky yÞ ; wFR ¼ pffiffiffiffiffiffi N seisht AF  h2 k2

where we have defined AF ¼  hv F kF =ð  EF þ D  2kss  2m0F

a 2

ð5Þ

  2b þ shÞ. The scattering coefficients r and t denote the normal

reflection and transmission amplitude, respectively. Importantly, the spin and, in particular, valley indices can be in the same or opposite sign in two ferromagnetic regions. Accordingly, we could allowed to access the possibility of intervalley Klein tunneling. The transmitted electron angle is determined from conservation of y-component of Fermi wavevector. The normalization factor N ensure that the particle current density of states is same, and by using Eq. (4) is given by:

N ¼ AF cosðshÞ þ

shkF  4m0 v F

 ða þ bÞ þ A2F ða  bÞ cosðshÞ:

ð6Þ

Taking into account the amplitudes of incoming and outgoing electrons (a and b) in the barrier region, and regarding the effect of gate potential, the wavefunction reads:

    1 1 iðskIx xþky yÞ e eiðskIx xþky yÞ : wI ¼ a þ b s eihI s eihI

ð7Þ

To evaluate the transmission probability, we match wavefunctions at the boundaries wFL jx¼0 ¼ wI jx¼0 ; wI jx¼L ¼ wFR jx¼L , and hence the spin and valley-dependent transmission amplitude is found as following:

246

M. Khezerlou, H. Goudarzi / Superlattices and Microstructures 86 (2015) 243–249 1

0.5 (a)

0.8

0.3

0.6

0.2

0.4

0.1

0.2

θ/π

0.4

0 0

1

2 Z/π

3

1

0.5

0.4

0.8

0.4

0.3

0.6

0.3

0.2

0.4

0.2

0.2

0.1

0.5

4

0

0.9 (c)

(b)

0.8 0.7 0.6

θ/π

θ/π

0.5

0.1

0.4 0.3 0.2 0.1

0

0

1

2 Z/π

3

0

0

4

0.5

0

1

2 Z/π

3

4

0.5 (d)

1 (e)

0.8

0.4

0.7

0.4

0.8

0.3

0.6

0.2

0.4

0.1

0.2

0.6 0.3

0.4

0.2

θ/π

θ/π

0.5

0.3 0.2

0.1

0.1 0

0

1

2 Z/π

3

0

4

0.6

0.2

0.4

0.1

0.2

0

0

1

2 Z/π

3

4

0

3

4

0

0.25

0.4

0.2

0.3

0.15

0.2

0.1

0.1

0.05

θ/π

0.3 θ/π

0.8

2 Z/π

(g)

(f) 0.4

1

0.5

1

0.5

0

0

0

1

2 Z/π

3

4

Fig. 3. (a–g) Plots of transmission probability for distinct ferromagnetic sections, as a function of incident angle and barrier strength when h ¼ 0:5k and EF ¼ 1:05 eV. There are 7 permitted cases that are relevant to transition from (a) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ, (b) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ, (c) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ, (d) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ, (e) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ, (f) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ and (g) ðs ¼ 1; s ¼ 1Þ to ðs ¼ 1; s ¼ 1Þ.

ts;s ¼

NL 4sL sI AFL cosðsL hÞ NR M 2 M 3 eiZ  M1 M 4 eiZ

where

M 1 ¼ sL AFL eisL h  sI ; is 3 h t

M 3 ¼ sR AFR e

þ sI ;

M2 ¼ sL AFL eisL h þ sI M4 ¼ sR AFR eis3 ht  sI :

This parameter can be calculated for different valley and spin indexes in two ferromagnetic regions.

ð8Þ

M. Khezerlou, H. Goudarzi / Superlattices and Microstructures 86 (2015) 243–249

2

247

α=0.43, β=2.21 α=0.43, β=0 α=0, β=2.21

G/G

0

1.5

1

0.5

0

0

0.5

1

1.5 Z/π

2

2.5

3

Fig. 4. The charge conductance versus barrier strength is depicted when h ¼ 0:5k and EF ¼ 1:05 eV. The role of the topological terms on the conductance has been demonstrated.

65 60

MR (%)

55 50 45 40 35

0

0.5

1

1.5 Z/π

2

2.5

3

Fig. 5. The magnetoresistance as a function of barrier strength is plotted. We have set h ¼ 0:5k and EF ¼ 1:05 eV. The parallel and antiparallel conductance includes 4 and 3 permitted cases, respectively.

3. Valley-dependent Klein tunneling In order to consider the transport of charge as a function of thin barrier with a gate voltage V 0 between two ferromagnetic regions separated by a insulator, we propose a F=I=F junction on top of a monolayer molybdenum disulfide. The inequivalent valence spin subbands in two valleys can be considered as an interesting platform to study the effect of spin flips with valley exchanging. To clarify this issue, the dispersion relation in momentum space near the K and K 0 valleys of the p-doped ferromagnetic region is shown in Fig. 2 for different spin subbands. It is shown that the exchange field h results in a distinct valley-related energy splitting in the valence band, and, consequently, effective energy gap varies notably from one valley (s ¼ 1) to the other (s ¼ 1). Therefore, we start by defining critical chemical potential in valence band edges for each set of spin and valleys, since the transmission of quasiparticles may be strongly suppressed apart from these values. The exchange field is taken in units of the spin–orbit coupling h ¼ 0:5k for the numerical results. Thus, for the possible spin and valley parameters ðs ¼ 1; s ¼ 1Þ; ðs ¼ 1; s ¼ 1Þ; ðs ¼ 1; s ¼ 1Þ and ðs ¼ 1; s ¼ 1Þ the critical values of chemical potential in order to determine the Fermi level EFc ¼ D þ 2kss  sh are given by 0.89, 1.01, 1.05 and 0.85, respectively. Because the transport of charge can include various combinations of valley and spin indices, the magnitude of Fermi energy in this structure needs to necessarily be determined by the incident electron from K 0 valley with spin-up. Otherwise, if the Fermi energy is taken to be EF ¼ 1:01 eV, the incoming electrons from K 0 valley with spin up are filtered. By regarding four degree of freedom for each electron in the two distinct ferromagnetic sections, there are 16 cases to calculate the transmission amplitude. From these cases, only seven cases are allowed to have a permissible tunneling and contribute to the charge current. This important fact can be understood from the band structure in the ferromagnetic regions. Accordingly, valley and spin can not be simultaneously flipped from the left ferromagnetic region to the right. Also, one of the 9 forbidden cases is relevant to transition from indices s ¼ 1; s ¼ 1 to the s ¼ 1; s ¼ 1, because of the Fermi level intersects the valence band edge in the case s ¼ 1; s ¼ 1. Our numerical results indicate that we have no valley exchange tunneling in the system.

248

M. Khezerlou, H. Goudarzi / Superlattices and Microstructures 86 (2015) 243–249

The behavior of transmission probability for seven possible cases in terms of barrier strength and electron incident angle is shown in Fig. 3. Similar to other relativistic structures such as graphene, transmission probability oscillates with barrier strength, and period of these oscillations is a value of p, independently on spin and valley indices. Remarkably, Klein effect occurs in Z ¼ np (n is integer) only for parallel alignment of magnetizations in both K and K 0 valley as shown in Fig. 3(a), (b), (e) and (f). The transmission probability for antiparallel case, for similar conditions, is smaller than one, and depends on valley index (see Fig. 3(c), (d) and (g)). These results can be described by two important properties of ML-MDS as spin–orbit coupling and valley degeneracy breaking. Also, the incident angle dependence of the transmission probability is different. 2 The minimum value of ts;s occurs when the incident electron from K 0 valley with spin-up tunnels from the potential barrier as spin-down (Fig. 3(g)), while the magnitude of this probability is noticeable in K 0 valley for smaller incident angles (see, Fig. 3(c)). 4. Equivalent valleys magnetoresistance The dimensionless normalized valley and spin resolved conductance is defined as:

Gs;s ¼ G0ðs;sÞ

Z p=2 2 t s;s cos hdh;

ð9Þ

0

where we define G0ðs;sÞ ¼ e2 =hNðEF Þ with the density of state being NðEF Þ ¼ kF W=p, in which W being the width of the ML-MDS sheet. Regarding the tunneling behavior of structure, the charge conductance for all incoming electrons versus barrier strength is demonstrated in Fig. 4. As expected, the period of p is observed for conductance, and it has a constant maximum in Z ¼ np. Moreover, effect of the mass asymmetry term and topological characteristic of effective ML-MOS Hamiltonian on charge conductance can be considered. We find that presence of the a term in the Hamiltonian has no effect on the conductance, while the conductance amplitude decreases by decreasing the b term. Using conductance for parallel and antiparallel spin polarizations, one can calculate the magnetoresistance (MR) of the system:

100

(a)

90 80

τ=−1

MR (%)

70 60 50 40

τ=+1

30 20 10 0

0

100

0.5

1

1.5 Z/π

2

2.5

3

2.5

3

(b)

80

τ=−1

MR (%)

60 40

τ=+1

20 0 −20

0

0.5

1

1.5 Z/π

2

Fig. 6. The magnetoresistance for different valleys versus barrier strength is separately depicted. We have set EF ¼ 1:05 eV and (a) h ¼ 0:5k (b) h ¼ k.

M. Khezerlou, H. Goudarzi / Superlattices and Microstructures 86 (2015) 243–249

  G s0 ;s0 : MRs;s % ¼ 100  1  Gs;s

249

ð10Þ

Difference between magnitudes of conductance in four possible parallel cases and three possible antiparallel once leads to an important physical quantity of magnetoresistance. The MR effect in Dirac materials-based junctions has been investigated in many previous experimental and theoretical studies [26–30]. Due to the relatively weakness of spin–orbit coupling in carbon-based nanostructures, the spin coherence times is long. So, magnitude of magnetoresistance in MoS2 -based junction is noticeable in presence of strong spin–orbit coupling. In Fig. 5, behavior of magnetoresistance of system is presented as a function of barrier strength. It can be clearly seen that MR exhibits the oscillation feature, and maximum magnitude occurs at Z ¼ np. This feature corresponds to the same period of conductance of parallel and antiparallel cases. In addition, to observe the effect of valley degeneracy breaking, we calculate MR for different valleys, separately, as shown in Fig. 6(a) and (b). It is remarkable that the MR of junction for K 0 valley has high values comparing with that for K valley. This is reasonable due to increasing the gap relative to spin–orbit coupling in the valence band. However, by increasing the exchange field upto h ¼ k, the MR can appear with a negative value in K valley (Fig. 6(b)). Varying exchange field affects both magnitude and sign of the MR, but changes the shape of MR curves. 5. Conclusion In summary, the valley and spin transport characteristics of a monolayer MoS2 F=I=F junction have been studied. Using the modified Dirac Hamiltonian, that in addition to first order term of momentum contains the quadratic terms corresponding to the topological and difference mass between electron and hole terms, and matching the wavefunctions in boundaries of junction, we have obtained analytical expression for the transmission probability at the interface. Relative to the effective spin–orbit coupling in valence band, the allowed values of spin and valley indices have been fixed in the two ferromagnetic sections. The transmission probability of an incident electron from left ferromagnetic region with spin number s from valley s extremely depends on the Fermi level of valence band. We have predicted the requirements to occur Klein tunneling. Also, we have proposed seven possibilities to have tunneling when the chemical potential of first section is in valence band edge of K 0 valley with spin-up. In particular, simultaneously spin and valley flip is energetically forbidden in two ferromagnetic regions. We have observed that the amplitude of conductance oscillations varies with topological term b of MoS2 Hamiltonian. The necessity of equivalent valleys to have tunneling has been found to affect significantly the magnetoresistance of junction comparing with similar graphene-based junctions. References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30]

O. Klein, Z. Phys. 53 (1929) 157. N. Dombey, A. Calogeracos, Phys. Rep. 315 (1999) 41. K.S. Novoselov, A.K. Geim, S.V. Morozov, D. Jiang, Y. Zhang, S.V. Dubonos, I.V. Grigorieva, A.A. Firsov, Science 306 (2004) 666. K.F. Mak, C. Lee, J. Hone, J. Shan, T.F. Heinz, Phys. Rev. Lett. 105 (2010) 136805; A. Splendiani, L. Sun, Y. Zhang, T. Li, J. Kim, C.Y. Chim, G. Galli, F. Wang, Nano Lett. 10 (2010) 1271. T. Korn, D. Stich, R. Schulz, D. Schuh, W. Wegscheider, C. Schller, Appl. Phys. Lett. 99 (2011) 102109. B. Radisavljevic, A. Radenovic, J. Brivio, V. Giacometti, A. Kis, Nat. Nanotech. 6 (2011) 147. M.I. Katsnelson, K.S. Novoselov, A.K. Geim, Nat. Phys. 2 (2006) 620. M.I. Katsnelson, K.S. Novoselov, Solid State Commun. 3 (2007) 143. J.P. Zhang, J.H. Yuan, Eur. Phys. J. B 85 (2012) 100. M.R. Setare, D. Jahani, J. Phys.: Condens. Matter 22 (2010) 245503. C. Lee, Q. Li, W. Kalb, X.-Z. Liu, H. Berger, R.W. Carpick, J. Hone, Science 328 (2010) 76. D. Xiao, W. Yao, Q. Niu, Phys. Rev. Lett. 99 (2007) 236809. Z.Y. Zhu, Y.C. Cheng, U. Schwingenschlogl, Phys. Rev. B 84 (2011) 153402. J. Zhang, J.M. Soon, K.P. Loh, J. Yin, J. Ding, M.B. Sullivian, P. Wu, Nano Lett. 7 (2007) 2370. Y. Li, Z. Zhou, S. Zhang, Z. Chen, J. Am. Chem. Soc. 130 (2009) 16739. S. Mathew et al, Appl. Phys. Lett. 101 (2012) 102103. Y. Ma, Y. Dai, M. Guo, C. Niu, Y. Zhu, B. Huang, ACS Nano 6 (2012) 1695. S. Tongay, S.S. Varnoosfaderani, B.R. Appleton, J. Wu, Appl. Phys. Lett. 101 (2012) 123105. R. Mishra, W. Zhou, S.J. Pennycook, S.T. Pantelides, J.-C. Idrobo, Phys. Rev. B 88 (2013) 144409. A. Vojvodic, B. Hinnemann, J.K. Norskov, Phys. Rev. B 80 (2009) 125416. C. Ataca, H. Sahin, E. Akturk, S. Ciraci, J. Phys. Chem. C 115 (2011) 3934. Di Xiao, G.-B. Liu, W. Feng, X. Xu, W. Yao, Phys. Rev. Lett. 108 (2012) 196802. L. Majidi, R. Asgari, Phys. Rev. B 90 (2014) 165440. L. Majidi, H. Rostami, R. Asgari, Phys. Rev. B 89 (2014) 045413. H. Rostami, A.G. Moghaddam, R. Asgari, Phys. Rev. B 88 (2013) 085440. A. Saffarzadeh, M. Ghorbani, Eur. Phys. J. B 67 (2009) 239. L. Brey, H.A. Fertig, Phys. Rev. B 76 (2007) 205435. F. Zhai, K. Chang, Phys. Rev. B 77 (2008) 113409. Y.-B. Zhou, B.-H. Han, Zh.-M. Liao, H.-Ch. Wu, D.-P. Yu, Appl. Phys. Lett. 98 (2011) 222502. H. Goudarzi, M. Khezerlou, A. Alav, Indian J. Phys. 89 (2015) 55.