We assess the potential of two-terminal graphene-hexagonal boron nitride-graphene resonant tunneling diodes as high-frequency oscillators, using self-consistent quantum transport and electrostatic simulations to determine the time-dependent response of the diodes in a resonant circuit. We quantify how the frequency and power of the current oscillations depend on the diode and circuit parameters including the doping of the graphene electrodes, device geometry, alignment of the graphene lattices, and the circuit impedances. Our results indicate that current oscillations with frequencies of up to several hundred GHz should be achievable.

Resonant tunneling diodes (RTDs) operating at 1.4 THz and 10 μW output power have been demonstrated recently.1–3 An addition to the family of RTDs is the graphene tunnel transistor,4–15 in which negative differential conductance (NDC), with a room temperature peak-to-valley ratio (PVR) of 2:1, arises from constraints imposed by energy and momentum conservation of Dirac Fermions, which tunnel through a hexagonal boron nitride (hBN) barrier.6,7

Here, we analyse how the device and circuit parameters can be tuned to increase the operating frequency of graphene RTDs (GRTDs). Our model device, shown schematically in Fig. 1(a), comprises two graphene layers separated by a hBN tunnel barrier of thickness, d. The bottom (B) and top (T) graphene electrodes are arranged in an overlapping cross formation, resulting in a tunneling area, A = 1 μm2. We consider the general case when the two graphene crystalline lattices are slightly misorientated by a twist angle, θ, see Fig. 1(a). The tunnel current, Ib, is particularly sensitive to this angle.7 A voltage, Vb, applied between top and bottom graphene layers [Fig. 1(b)] induces a charge density, ρB,T, in each layer and causes Ib to flow through the hBN barrier. The graphene layers, with in-plane sheet resistance, R, carry current, I, (black arrows) from two pairs of Ohmic contacts [orange in Fig. 1(a)] to the tunneling region, i.e., I/2, flows to/from each contact. The electrostatics4 are governed by the equation eVb=μBμT+ϕb, where ϕb=eFbd and Fb is the electric field in the barrier, e is the magnitude of the electronic charge, and μB,T are the two Fermi levels [see Fig. 1(b)].

FIG. 1.

(a) Schematic diagram: device comprises bottom (red) and top (blue) graphene lattices, misaligned by an angle θ and separated by a hBN tunnel barrier (dark green). The current, I, passes through the tunnel barrier between the graphene electrode layers to/from Ohmic contacts (orange). The diode is mounted on a hBN layer (light green) and an insulating substrate (purple). (b) Schematic diagram of the resonant circuit incorporating the diode (in box) showing the voltage applied, V, circuit inductance, L, and resistance, R. The band diagram of the tunnel region is shown in box.

FIG. 1.

(a) Schematic diagram: device comprises bottom (red) and top (blue) graphene lattices, misaligned by an angle θ and separated by a hBN tunnel barrier (dark green). The current, I, passes through the tunnel barrier between the graphene electrode layers to/from Ohmic contacts (orange). The diode is mounted on a hBN layer (light green) and an insulating substrate (purple). (b) Schematic diagram of the resonant circuit incorporating the diode (in box) showing the voltage applied, V, circuit inductance, L, and resistance, R. The band diagram of the tunnel region is shown in box.

Close modal

A device with NDC provides instability that can generate self-sustained current oscillations when placed in an RLC circuit.16,17 To investigate the frequency response of the GRTD, we solve the time-dependent current continuity and Poisson equations self-consistently, using the Bardeen transfer Hamiltonian method to calculate Ib

Ib=8πekB,kT|M|2[fB(EB)fT(ET)]δ(EBETϕb),
(1)

as a function of time, t, and Vb. The summation is over all initial and final states, with wavevectors, kB,T, measured relative to the position of the nearest Dirac point in the bottom layer, K±=(±4π/3a0,0), where ± distinguishes the two non-equivalent Dirac points in the Brillouin zone and a0 = 2.46 Å is the graphene lattice constant. The Fermi function in each electrode is fB,T(EB,T)=[1+e(EB,TμB,T)/kT]1, where EB,T=sB,TvFkB,T is the electron energy, vF = 106 ms−1 is the carrier speed, and sB,T = ±1 labels electrons in the conduction (+) and valence (–) bands, at temperature T = 300 K. Tunneling between equivalent valleys gives the same contribution to Ib, so we consider transitions between K+ points only. In Eq. (1) the matrix element, M, is

M=Ξγ(θ)g(φB,φT)VS(qΔK),
(2)

where Ξ = ξeκd, ξ is a normalisation constant determined by comparison with recent measurements7 of Ib, γ(θ) is the spatial overlap integral of the cell-periodic part of the wavefunction, g(φB,φT) describes electron chirality, VS is the elastic scattering potential, and q = kBkT (see below). The decay constant of the wavefunction in the barrier is κ=2mΔb/, where the barrier height, Δb = 1.5 eV, and the effective electron mass in the barrier m = 0.5 me.4 

In recently studied GRTDs,7 the crystal lattices of the graphene layers are only misorientated by an angle θ ≈ 1°. Nevertheless, this gives rise to a significant misalignment of the Dirac cones of the two layers, ΔK=(R(θ)I)K+, where R(θ) is the 2D rotation matrix. When θ<2°,|q||ΔK|=ΔK, and electrons tunnel with conservation of in-plane momentum. However, tunneling electrons can scatter elastically from impurities and defects, broadening the features in Ib(Vb).18,19 Therefore, we use a scattering potential VS(q)=V0/(q2+qc2), with amplitude V0 = 10 meV and length scale 1/qc = 15 nm, which gives the best fit to the measured Ib(Vb).7 The misorientation also reduces the spatial overlap integral, γ(θ). The chiral wavefunctions give rise to the term g(φB,φT)=1+sBeiφB+sTeiφT+sBsTei(φBφT), where φ=tan1(ky/kx) is the wavevector orientation.

Fig. 2 shows the equilibrium (static) Ib(Vb) curve (blue), where Vb ≈ V and Ib = I, calculated for an undoped device with θ = 0.9° and d = 1.3 nm (4 layers of hBN), similar to that studied in Ref. 7. The calculated Ib(Vb) curve reproduces the measured line-shape, position of the resonant peak and current amplitude.6,7 The peak occurs when many electrons can tunnel with momentum conservation, i.e., qΔK0, corresponding to a resonant increase in the matrix element M, i.e., when ϕb=vFΔKforθ close to 1°. Temperature has negligible effect on the I(Vb) curve when Vb > kT/e ∼ 30 mV.6,7

FIG. 2.

Equilibrium and non-equilibrium current-voltage curves calculated for θ = 0.9°, L = 140 nH, and R = 50 Ω. Blue: equilibrium current-voltage characteristic Ib(Vb). Note, in equilibrium, Vb ≈ V and Ib = I. Green: time-averaged current I(t)t vs V. Red: peak-to-peak voltage amplitude (right scale) of the stable current oscillations. Inset: I(t) plot showing stable oscillations with f = 4.2 GHz.

FIG. 2.

Equilibrium and non-equilibrium current-voltage curves calculated for θ = 0.9°, L = 140 nH, and R = 50 Ω. Blue: equilibrium current-voltage characteristic Ib(Vb). Note, in equilibrium, Vb ≈ V and Ib = I. Green: time-averaged current I(t)t vs V. Red: peak-to-peak voltage amplitude (right scale) of the stable current oscillations. Inset: I(t) plot showing stable oscillations with f = 4.2 GHz.

Close modal

We now consider the non-equilibrium charge dynamics when the device is in a series circuit with inductance, L, and resistance, R, see Fig. 1(b); the diode has an in-built capacitance, C. The device has no in-built inductance and therefore, to oscillate, requires a series inductance. Recently, self-excited plasma oscillations have also been shown to cause instabilities and oscillations in GRTDs.15 The primary contribution to R arises from the graphene electrodes4 and depends on their charge densities, ρB,T. This dependence does not significantly affect the high-frequency (HF) response: for most of the oscillation period, changes in ρB,T do not greatly affect R. Therefore, we take R to be independent of t. However, R can be changed by altering the device geometry, e.g., by reducing the length of the electrodes, and we consider this effect on the performance of the GRTD. We also consider how L affects the frequency, which could be controlled by careful design of the microwave circuit, e.g., by using a resonant cavity or integrated patch antennas.2 

We determine the current, I(t), in the contacts and external circuit by solving20 self-consistently the current-continuity equations: B,T/dt = ±(IbI)/A, where the + (–) sign is for the bottom (top) graphene layers, see Fig. 1(b), ρB,T are related by Poisson's equation: ϵFb = ρBρBD = −(ρTρTD), in which ϵ = ϵ0ϵr and ϵr = 3.9 (Refs. 4 and 21) is the permittivity of the barrier, and ρBD (ρTD) are the doping densities in each layer. The voltages across the inductor and resistor, VL and VR, are given by dI/dt=VL/L,VR=IR, and V = VR + Vb + VL.

Following initial transient behavior, I(t) either decays to a constant value or oscillates with a frequency, f, and time-averaged current, I(t)t. Fig. 2, inset, shows a typical I(t) curve, for V = 0.48 V, exhibiting oscillations with f = 4.2 GHz. In Fig. 2, we show I(t)t versus V (green) and Ib(Vb) (blue curve) for an undoped device, with θ = 0.9°, placed in a resonant circuit with R = 50 Ω and L = 140 nH. The plot reveals that when V is tuned in the NDC region (0.55 V < V < 0.8 V), ΔVL=VLmaxVLmin (red curve) becomes non-zero indicating that self-sustained oscillations are induced. Here, VLmax/min is the maximum/minimum voltage dropped across the inductor during an oscillation period. Also, the I(t)t versus V curve (green) diverges from the static current, Ib(Vb), (blue) in the NDC region due to asymmetric rectification of I(t) in the strongly nonlinear NDC region of Ib(Vb). When the device is biased in regions of positive differential conductance, i.e., V < 0.55 V or V > 0.8 V, oscillations are suppressed and I(t)t converges to Ib(Vb).

This behavior is similar to that recently measured in a GRTD, where oscillations with f ∼ 2 MHz were reported.7 That device had high circuit capacitance due to large-area contact pads and coupling to the doped Si substrate (gate). This effect can be modelled by placing a capacitor in parallel with the GRTD. Including this large capacitance (65 pF) limits the maximum observed f value.7 When parasitic circuit capacitances are minimised, using the geometry shown in Fig. 1(a), the only significant contribution to the total capacitance is from the overlap area of the graphene electrodes, as described by the charge-continuity equation. This enables us to investigate the potential of GRTDs optimised for HF applications.

A small signal analysis16 provides insight into how L, R, and the form of Ib(Vb) affect the circuit response and gives an approximate frequency

fs=f0(1R/RN)QN2(1QN2R/RN)2/4,
(3)

where RN is the maximum negative differential resistance of the equilibrium I(V) curve (Ib(Vb)), the circuit factor QN=RNC/L, and f0=1/2πLC. Here, RN is large, therefore f s≈ f0. For a given C (that depends on A and d), f can be increased by reducing L. The decay parameter of the small signal analysis reveals that the circuit will oscillate only if

(RN/RQN2)>0.
(4)

Consequently, R, and the form of Ib(Vb) are also important for optimising the HF performance.

We now consider the self-consistent simulation of the charge dynamics. Fig. 3(a) shows the fmax(R) curve calculated for the diode parameters, which compare well to recent measurements,7 used to produce the Ib(Vb) curves in Fig. 2. We determine fmax(R) by finding the smallest L value for self-sustained oscillations. The solid part of the curve in Fig. 3(a) shows fmax over the range of R values that can be achieved with small modifications to the design of existing devices, e.g., by reducing the length of the graphene between the tunnel area and the Ohmic contacts, or by doping the electrodes. The dashed part is calculated for R values that may be possible in future configurations. The curve reveals that for a recently attained R = 50 Ω,22,23fmax = 1.8 GHz.

FIG. 3.

(a) fmax(R) calculated when NL = 4. Inset: Log-log plot. (b) fmax vs R when NL = 2 (red), 3 (green), and 4 (blue). Inset: fmax vs NL calculated when R = 50 Ω. Curves are solid for R values presently obtainable in GRTDs and dashed for R values that could be achieved by future device designs. All curves are for undoped devices.

FIG. 3.

(a) fmax(R) calculated when NL = 4. Inset: Log-log plot. (b) fmax vs R when NL = 2 (red), 3 (green), and 4 (blue). Inset: fmax vs NL calculated when R = 50 Ω. Curves are solid for R values presently obtainable in GRTDs and dashed for R values that could be achieved by future device designs. All curves are for undoped devices.

Close modal

Fig. 3(a), inset, reveals the power law fmaxR−0.505, which can be derived by setting Eq. (4) equal to zero and rearranging to find the smallest L value for a given R, RN, and C.16 For this case

fmaxs=(2πCRRN)1R0.5,
(5)

which compares well with the full signal analysis.

To increase fmax, we can also modify Ib(Vb). Reducing the number of layers, NL, in the hBN tunnel barrier increases Ib (∼20× for each layer removed24), thus reducing RN and increasing fmax, see Eq. (5). Fig. 3(b) shows fmax(R) calculated for a device with NL = 4 (blue), 3 (green), and 2 (red). Reducing d produces a large gain in fmax for all R. For example, fmax for a device with NL = 2 is at least an order of magnitude higher than when NL = 4 (e.g., for R = 50 Ω, fmax = 26 GHz when NL = 2, compared to fmax = 1.8 GHz when NL = 4).

The Ib(Vb) characteristics can also be modified by doping the graphene chemically25,26 or, equivalently, by applying a gate voltage, Vg, to shift the current peak and, thereby, change RN and the PVR.6,7 In Fig. 4(a), we show Ib(Vb) curves calculated when NL = 2 for an undoped (red curve) and an asymmetrically doped device with ρBD/e = 1013 cm−2 and ρTD/e = 0 (green curve). When ρBD > 0, the resonant peak occurs at higher Vb than when ρBD = 0, and the current peak magnitude is higher, raising the PVR from 1.5 to 3.5.

FIG. 4.

(a) Ib(Vb) characteristics calculated for a doped (green, ρBD/e = 1013 cm−2) and undoped (red) device, with NL = 2. The arrow shows the shoulder that arises due to the quantum capacitance effect. (b) fmax vs R curves for the devices in (a). Inset: fmax vs (ρBD/e) calculated when R = 50 Ω, with ρTD/e = 0. Curves in (b) are solid for R values presently obtainable in GRTDs and dashed for R values that could be achieved by future device designs.

FIG. 4.

(a) Ib(Vb) characteristics calculated for a doped (green, ρBD/e = 1013 cm−2) and undoped (red) device, with NL = 2. The arrow shows the shoulder that arises due to the quantum capacitance effect. (b) fmax vs R curves for the devices in (a). Inset: fmax vs (ρBD/e) calculated when R = 50 Ω, with ρTD/e = 0. Curves in (b) are solid for R values presently obtainable in GRTDs and dashed for R values that could be achieved by future device designs.

Close modal

The shoulder of the green curve in Fig. 4(a), (arrowed) when ρBD/e = 1013 cm−2, arises from the low density of states around the Dirac point. This gives rise to an additional quantum capacitance,6,27CQ, whose effect is prominent when the chemical potential aligns with the Dirac point. The total capacitance is given by C1=CG1+CQ1, where CG = ϵ0ϵrA/d is the geometric capacitance. When μB,T passes through the Dirac point, CQ → 0 and, hence, C → 0, suggesting that the RC time constant could be reduced. In practice, CQ is small for only a small fraction of the oscillation period and so its effect on the fundamental frequency of I(t) is negligible.

Fig. 4(b) shows fmax(R) curves calculated for undoped (red) and doped (green) devices and reveals that the doped device is faster for all R. Fig. 4(b) inset shows that fmax increases with ρBD/e when R = 50 Ω; fmax increases by a factor of 1.3 when ρBD/e is increased to 1013 cm −2 (and fmax = 32 GHz) from ρBD/e = 0 (fmax = 26 GHz).

To quantify the possible benefits of lattice alignment, Fig. 5(a) shows the effect of changing θ on Ib(Vb). As θ increases, the position of the current peak shifts to higher Vb. The peak current amplitude, Ipeak, decreases as θ increases due to increasing misorientation of the spatial parts of the wavefunction, see Fig. 5(b), so that Ipeak could be ∼10× larger for an aligned device. However, for undoped samples, the PVR increases with increasing θ, see inset in Fig. 5(b), converging to a value of 3.4 as θ approaches 2°: at higher θ, more states are available to tunnel resonantly at the current peak.10 For the doped samples (ρBD/e = 1013 cm−2), the valley current is small for all θ, thus the PVR is consistently large, see Fig. 5(c). Consequently, the increase in current magnitude, which results from alignment, leads to higher f values without the power reduction associated with undoped samples. We find that, generally, RN ((fmaxs)2; Eq. (5)) decreases with decreasing θ, Fig. 5(c) inset, and with increasing ρBD, meaning that oscillation frequencies are highest for θ = 0° and when ρBD/e = 1013 cm−2. Fig. 5(d) shows that perfect alignment could increase fmax by a factor of ∼2, i.e., for R = 50 Ω, fmax = 65 GHz when θ = 0° compared to 32 GHz when θ = 0.9°. The numerical results diverge from the small signal analysis power law of fmaxR−0.5 as RN becomes small, see black curve of Fig. 5(d), and it becomes necessary to vary V to induce oscillations.

FIG. 5.

(a) Ib(Vb) curves calculated for samples with misalignment angles θ = 0° (black), 0.5° (blue), 0.9° (green), and 2° (magenta), taking ρBD/e = 0 cm−2 and NL = 2. (b) Current amplitude at the peak vs misalignment angle, θ. Inset: PVR of Ib(Vb) vs θ. (c) Ib(Vb) calculated when θ = 0° and ρBD/e = 1013 cm−2. Inset: RN(θ) for undoped (upper) and ρBD/e = 1013 cm−2 (lower) diodes. (d) fmax vs R curves calculated for an aligned sample (black) and misaligned sample with θ = 0.9° (green), when ρBD/e = 1013 cm−2. Curves are shown solid over the range of R presently obtainable in GRTDs and dashed for R values that could be achieved by future device designs. For all curves, ρTD/e = 0 cm−2.

FIG. 5.

(a) Ib(Vb) curves calculated for samples with misalignment angles θ = 0° (black), 0.5° (blue), 0.9° (green), and 2° (magenta), taking ρBD/e = 0 cm−2 and NL = 2. (b) Current amplitude at the peak vs misalignment angle, θ. Inset: PVR of Ib(Vb) vs θ. (c) Ib(Vb) calculated when θ = 0° and ρBD/e = 1013 cm−2. Inset: RN(θ) for undoped (upper) and ρBD/e = 1013 cm−2 (lower) diodes. (d) fmax vs R curves calculated for an aligned sample (black) and misaligned sample with θ = 0.9° (green), when ρBD/e = 1013 cm−2. Curves are shown solid over the range of R presently obtainable in GRTDs and dashed for R values that could be achieved by future device designs. For all curves, ρTD/e = 0 cm−2.

Close modal

In conclusion, we have investigated the performance of GRTDs as the active element in RLC oscillators. These devices could oscillate at mid-GHz frequencies, by careful design of the RLC circuit. We have also quantified the effect of changing the parameters of the GRTD. Reducing the barrier width (a modest change to the structure of existing devices) increases Ib, and thus raises the oscillation frequency by an order of magnitude. Adjusting the doping of the electrodes can enhance f. Finally, we have considered the effect of misalignment of the graphene electrodes: in devices with aligned lattices, frequencies approaching 1 THz may be attainable. Double barrier (GaIn)As/AlAs RTDs1 have similar Ib and Vb values as the GRTD reported here. We therefore expect that the GRTD will produce similar EM emission power (∼10 μW). Our results illustrate the potential of graphene tunnel structures in HF graphene electronics.

This work was supported by the EU Graphene Flagship Programme. K.S.N. and M.T.G. acknowledge the support of the Royal Society and of The Leverhulme Trust, respectively.

1.
S.
Suzuki
,
M.
Asada
,
A.
Teranishi
,
H.
Asugiyama
, and
H.
Yokoyama
,
Appl. Phys. Lett.
97
,
242102
(
2010
).
2.
Y.
Koyama
,
R.
Sekiguchi
, and
T.
Ouchi
,
Appl. Phys. Express
6
,
064102
(
2013
).
3.
M.
Feiginov
,
H.
Kanaya
,
S.
Suzuki
, and
M.
Asada
,
Appl. Phys. Lett.
104
,
243509
(
2014
).
4.
L.
Britnell
,
R. V.
Gorbachev
,
R.
Jalil
,
B. D.
Belle
,
F.
Schedin
,
A.
Mishchenko
,
T.
Georgiou
,
M. I.
Katsnelson
,
L.
Eaves
,
S. V.
Morozov
,
N. M. R.
Peres
,
J.
Leist
,
A. K.
Geim
,
K. S.
Novoselov
, and
L. A.
Ponomarenko
,
Science
335
,
947
(
2012
).
5.
T.
Georgiou
,
R.
Jalil
,
B. D.
Belle
,
L.
Britnell
,
R. V.
Gorbachev
,
S. V.
Morozov
,
Y. J.
Kim
,
A.
Gholinia
,
S. J.
Haigh
,
O.
Makarovsky
,
L.
Eaves
,
L. A.
Ponomarenko
,
A. K.
Geim
,
K. S.
Novoselov
, and
A.
Mishchenko
,
Nat. Nanotechnol.
8
,
100
103
(
2013
).
6.
L.
Britnell
,
R. V.
Gorbachev
,
A. K.
Geim
,
L. A.
Ponomarenko
,
A.
Mishchenko
,
M. T.
Greenaway
,
T. M.
Fromhold
,
K. S.
Novoselov
, and
L.
Eaves
,
Nat. Commun.
4
,
1794
(
2013
).
7.
A.
Mishchenko
,
J. S.
Tu
,
Y.
Cao
,
R. V.
Gorbachev
,
J. R.
Wallbank
,
M. T.
Greenaway
,
V. E.
Morozov
,
S. V.
Morozov
,
M. J.
Zhu
,
S. L.
Wong
,
F.
Withers
,
C. R.
Woods
,
Y.-J.
Kim
,
K.
Watanabe
,
T.
Taniguchi
,
E. E.
Vdovin
,
O.
Makarovsky
,
T. M.
Fromhold
,
V. I.
Fal'ko
,
A. K.
Geim
,
L.
Eaves
, and
K. S.
Novoselov
,
Nat. Nanotechnol.
9
,
808
813
(
2014
).
8.
B.
Fallahazad
,
K.
Lee
,
S.
Kang
,
J.
Xue
,
S.
Larentis
,
C.
Corbet
,
K.
Kim
,
H. C. P.
Movva
,
T.
Taniguchi
,
K.
Watanabe
,
L. F.
Register
,
S. K.
Banerjee
, and
E.
Tutuc
,
Nano Lett.
15
,
428
(
2015
).
9.
S.
Kang
,
B.
Fallahazad
,
L.
Kayoung
,
H.
Movva
,
K.
Kyounghwan
,
C. M.
Corbet
,
T.
Taniguchi
,
K.
Watanabe
,
L.
Colombo
,
L. F.
Register
,
E.
Tutuc
, and
S. K.
Banerjee
,
IEEE Electron Device Lett.
36
(
4
),
405
407
(
2015
).
10.
R. M.
Feenstra
,
D.
Jena
, and
G.
Gu
,
J. Appl. Phys.
111
,
043711
(
2012
).
11.
P.
Zhao
,
R. M.
Feenstra
,
G.
Gu
, and
D.
Jena
,
IEEE Trans. Electron Devices
60
,
951
957
(
2013
).
12.
L.
Brey
,
Phys. Rev. Appl.
2
,
014003
(
2014
).
13.
F. T.
Vasko
,
Phys. Rev. B
87
,
075424
(
2013
).
14.
V.
Ryzhii
,
A. A.
Dubinov
,
V. Y.
Aleshkin
,
M.
Ryzhii
, and
T.
Otsuji
,
Appl. Phys. Lett.
103
,
163507
(
2013
).
15.
V.
Ryzhii
,
A.
Satou
,
T.
Otsuji
,
M.
Ryzhii
,
V.
Mitin
, and
M. S.
Shur
,
J. Phys. D: Appl. Phys.
46
,
315107
(
2013
).
16.
M. E.
Hines
,
Bell Syst. Tech. J.
39
,
477
(
1960
).
17.
H.
Mizuta
and
T.
Tanoue
,
The Physics and Applications of Resonant Tunnelling Diodes
(
Cambridge University Press
,
1995
).
18.
S.
Adam
,
E. H.
Hwang
,
E.
Rossi
, and
S.
Das Sarma
,
Solid State Commun.
149
,
1072
(
2009
).
19.
Q.
Li
,
E. H.
Hwang
,
E.
Rossi
, and
S.
Das Sarma
,
Phys. Rev. Lett.
107
,
156601
(
2011
).
20.
M. T.
Greenaway
,
A. G.
Balanov
,
E.
Schöll
, and
T. M.
Fromhold
,
Phys. Rev. B
80
,
205318
(
2009
).
21.
K.
Kim
,
A.
Hsu
,
X.
Jia
,
S.
Kim
,
Y.
Shi
,
M.
Dresselhaus
,
T.
Palacios
, and
J.
Kong
,
ACS Nano
6
(
10
),
8583
8590
(
2012
).
22.
A. V.
Kretinin
,
Y.
Cao
,
J. S.
Tu
,
G. L.
Yu
,
R.
Jalil
,
K. S.
Novoselov
,
S. J.
Haigh
,
A.
Gholinia
,
A.
Mishchenko
,
M.
Lozada
,
T.
Georgiou
,
C. R.
Woods
,
F.
Withers
,
P.
Blake
,
G.
Eda
,
A.
Wirsig
,
C.
Hucho
,
K.
Watanabe
,
T.
Taniguchi
,
A. K.
Geim
, and
R. V.
Gorbachev
,
Nano Lett.
14
,
3270
3276
(
2014
).
23.
L.
Wang
,
I.
Meric
,
P. Y.
Huang
,
Q.
Gao
,
Y.
Gao
,
H.
Tran
,
T.
Taniguchi
,
K.
Watanabe
,
L. M.
Campos
,
D. A.
Muller
,
J.
Guo
,
P.
Kim
,
J.
Hone
,
K. L.
Shepard
, and
C. R.
Dean
,
Science
342
,
614
617
(
2013
).
24.
L.
Britnell
,
R. V.
Gorbachev
,
R.
Jalil
,
B. D.
Belle
,
F.
Schedin
,
M. I.
Katsnelson
,
L.
Eaves
,
S. V.
Morozov
,
A. S.
Mayorov
,
N. M. R.
Peres
,
A. H.
Castro Neto
,
J.
Leist
,
A. K.
Geim
,
L. A.
Ponomarenko
, and
K. S.
Novoselov
,
Nano Lett.
12
(
3
),
1707
1710
(
2012
).
25.
A.
Das
,
S.
Pisana
,
B.
Chakraborty
,
S.
Piscanec
,
S. K.
Saha
,
U. V.
Waghmare
,
K. S.
Novoselov
,
H. R.
Krishnamurthy
,
A. K.
Geim
,
A. C.
Ferrari
, and
A . K.
Sood
,
Nat. Nanotechnol.
3
,
210
215
(
2008
).
26.
H.
Liu
,
Y.
Liu
, and
D.
Zhu
,
J. Mater. Chem.
21
,
3335
3345
(
2011
).
27.
S.
Luryi
,
Appl. Phys. Lett.
52
(
6
),
501
503
(
1988
).