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 *μ*m^{2}. 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, *I _{b}*, is particularly sensitive to this angle.

^{7}A voltage,

*V*, applied between top and bottom graphene layers [Fig. 1(b)] induces a charge density,

_{b}*ρ*

_{B,T}, in each layer and causes

*I*to flow through the hBN barrier. The graphene layers, with in-plane sheet resistance,

_{b}*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 electrostatics

^{4}are governed by the equation $eVb=\mu B\u2212\mu T+\varphi b$, where $\varphi b=eFbd$ and

*F*is the electric field in the barrier,

_{b}*e*is the magnitude of the electronic charge, and

*μ*

_{B,T}are the two Fermi levels [see Fig. 1(b)].

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 *I _{b}*

as a function of time, *t*, and *V _{b}*. The summation is over all initial and final states, with wavevectors,

**k**

_{B,T}, measured relative to the position of the nearest Dirac point in the bottom layer, $K\xb1=(\xb14\pi /3a0,0)$, where ± distinguishes the two non-equivalent Dirac points in the Brillouin zone and

*a*

_{0}= 2.46 Å is the graphene lattice constant. The Fermi function in each electrode is $fB,T(EB,T)=[1+e(EB,T\u2212\mu B,T)/kT]\u22121$, where $EB,T=sB,T\u210fvFkB,T$ is the electron energy,

*v*= 10

_{F}^{6}ms

^{−1}is the carrier speed, and

*s*

_{B,T}= ±1 labels electrons in the conduction (+) and valence (–) bands, at temperature

*T*= 300 K. Tunneling between equivalent valleys gives the same contribution to

*I*, so we consider transitions between

_{b}**K**

^{+}points only. In Eq. (1) the matrix element,

*M*, is

where Ξ = *ξe*^{−κd}, *ξ* is a normalisation constant determined by comparison with recent measurements^{7} of *I _{b}*,

*γ*(

*θ*) is the spatial overlap integral of the cell-periodic part of the wavefunction, $g(\phi B,\phi T)$ describes electron chirality,

*V*is the elastic scattering potential, and

_{S}**q**=

**k**

_{B}–

**k**

_{T}(see below). The decay constant of the wavefunction in the barrier is $\kappa =2m\Delta b/\u210f$, where the barrier height, Δ

_{b}= 1.5 eV, and the effective electron mass in the barrier

*m*= 0.5

*m*.

_{e}^{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, $\Delta K=(R(\theta )\u2212I)K+$, where *R*(*θ*) is the 2D rotation matrix. When $\theta <2\xb0,\u2009|q|\u2248|\Delta K|=\Delta K$, and electrons tunnel with conservation of in-plane momentum. However, tunneling electrons can scatter elastically from impurities and defects, broadening the features in *I _{b}*(

*V*).

_{b}^{18,19}Therefore, we use a scattering potential $VS(q)=V0/(q2+qc2)$, with amplitude

*V*

_{0}= 10 meV and length scale 1/

*q*= 15 nm, which gives the best fit to the measured

_{c}*I*(

_{b}*V*).

_{b}^{7}The misorientation also reduces the spatial overlap integral,

*γ*(

*θ*). The chiral wavefunctions give rise to the term $g(\phi B,\phi T)=1+sBei\phi B+sTe\u2212i\phi T+sBsTei(\phi B\u2212\phi T)$, where $\phi =tan\u22121(ky/kx)$ is the wavevector orientation.

Fig. 2 shows the equilibrium (static) *I _{b}*(

*V*) curve (blue), where

_{b}*V*≈

_{b}*V*and

*I*=

_{b}*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

*I*(

_{b}*V*) curve reproduces the measured line-shape, position of the resonant peak and current amplitude.

_{b}^{6,7}The peak occurs when many electrons can tunnel with momentum conservation, i.e., $q\u2212\Delta K\u223c0$, corresponding to a resonant increase in the matrix element

*M*, i.e., when $\varphi b=\u210fvF\Delta Kfor\theta $ close to 1°. Temperature has negligible effect on the

*I*(

*V*) curve when

_{b}*V*>

_{b}*kT*/

*e*∼ 30 mV.

^{6,7}

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 electrodes^{4} 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 solving^{20} self-consistently the current-continuity equations: *dρ*_{B,T}/*dt* = ±(*I _{b}* –

*I*)/

*A*, where the + (–) sign is for the bottom (top) graphene layers, see Fig. 1(b),

*ρ*

_{B,T}are related by Poisson's equation:

*ϵF*=

_{b}*ρ*–

_{B}*ρ*= −(

_{BD}*ρ*–

_{T}*ρ*), in which

_{TD}*ϵ*=

*ϵ*

_{0}

*ϵ*and

_{r}*ϵ*= 3.9 (Refs. 4 and 21) is the permittivity of the barrier, and

_{r}*ρ*(

_{BD}*ρ*) are the doping densities in each layer. The voltages across the inductor and resistor,

_{TD}*V*and

_{L}*V*, are given by $dI/dt=VL/L,VR=IR$, and

_{R}*V*=

*V*+

_{R}*V*+

_{b}*V*.

_{L}Following initial transient behavior, *I*(*t*) either decays to a constant value or oscillates with a frequency, *f*, and time-averaged current, $\u3008I(t)\u3009t$. 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 $\u3008I(t)\u3009t$ versus *V* (green) and *I _{b}*(

*V*) (blue curve) for an undoped device, with

_{b}*θ*= 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), $\Delta VL=VLmax\u2212VLmin$ (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 $\u3008I(t)\u3009t$ versus

*V*curve (green) diverges from the static current,

*I*(

_{b}*V*), (blue) in the NDC region due to asymmetric rectification of

_{b}*I*(

*t*) in the strongly nonlinear NDC region of

*I*(

_{b}*V*). When the device is biased in regions of positive differential conductance, i.e.,

_{b}*V*< 0.55 V or

*V*> 0.8 V, oscillations are suppressed and $\u3008I(t)\u3009t$ converges to

*I*(

_{b}*V*).

_{b}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 analysis^{16} provides insight into how *L*, *R*, and the form of *I _{b}*(

*V*) affect the circuit response and gives an approximate frequency

_{b}where *R _{N}* is the maximum negative differential resistance of the equilibrium

*I*(

*V*) curve (

*I*(

_{b}*V*)), the circuit factor $QN=RNC/L$, and $f0=1/2\pi LC$. Here,

_{b}*R*is large, therefore

_{N}*f*

^{s}^{ }≈

*f*

_{0}. 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

Consequently, *R*, and the form of *I _{b}*(

*V*) are also important for optimising the HF performance.

_{b}We now consider the self-consistent simulation of the charge dynamics. Fig. 3(a) shows the *f _{max}*(

*R*) curve calculated for the diode parameters, which compare well to recent measurements,

^{7}used to produce the

*I*(

_{b}*V*) curves in Fig. 2. We determine

_{b}*f*(

_{max}*R*) by finding the smallest

*L*value for self-sustained oscillations. The solid part of the curve in Fig. 3(a) shows

*f*over the range of

_{max}*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,23}

*f*= 1.8 GHz.

_{max}Fig. 3(a), inset, reveals the power law *f _{max}* ∝

*R*

^{−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*,

*R*, and

_{N}*C*.

^{16}For this case

which compares well with the full signal analysis.

To increase *f _{max}*, we can also modify

*I*(

_{b}*V*). Reducing the number of layers,

_{b}*N*, in the hBN tunnel barrier increases

_{L}*I*(∼20× for each layer removed

_{b}^{24}), thus reducing

*R*and increasing

_{N}*f*, see Eq. (5). Fig. 3(b) shows

_{max}*f*(

_{max}*R*) calculated for a device with

*N*= 4 (blue), 3 (green), and 2 (red). Reducing

_{L}*d*produces a large gain in

*f*for all

_{max}*R*. For example,

*f*for a device with

_{max}*N*= 2 is at least an order of magnitude higher than when

_{L}*N*= 4 (e.g., for

_{L}*R*= 50 Ω,

*f*= 26 GHz when

_{max}*N*= 2, compared to

_{L}*f*= 1.8 GHz when

_{max}*N*= 4).

_{L}The *I _{b}*(

*V*) characteristics can also be modified by doping the graphene chemically

_{b}^{25,26}or, equivalently, by applying a gate voltage,

*V*, to shift the current peak and, thereby, change

_{g}*R*and the PVR.

_{N}^{6,7}In Fig. 4(a), we show

*I*(

_{b}*V*) curves calculated when

_{b}*N*= 2 for an undoped (red curve) and an asymmetrically doped device with

_{L}*ρ*/

_{BD}*e*= 10

^{13 }cm

^{−2}and

*ρ*/

_{TD}*e*= 0 (green curve). When

*ρ*> 0, the resonant peak occurs at higher

_{BD}*V*than when

_{b}*ρ*= 0, and the current peak magnitude is higher, raising the PVR from 1.5 to 3.5.

_{BD}The shoulder of the green curve in Fig. 4(a), (arrowed) when *ρ _{BD}*/

*e*= 10

^{13 }cm

^{−2}, arises from the low density of states around the Dirac point. This gives rise to an additional quantum capacitance,

^{6,27}

*C*, whose effect is prominent when the chemical potential aligns with the Dirac point. The total capacitance is given by $C\u22121=CG\u22121+CQ\u22121$, where

_{Q}*C*=

_{G}*ϵ*

_{0}

*ϵ*/

_{r}A*d*is the geometric capacitance. When

*μ*

_{B,T}passes through the Dirac point,

*C*→ 0 and, hence,

_{Q}*C*→ 0, suggesting that the RC time constant could be reduced. In practice,

*C*is small for only a small fraction of the oscillation period and so its effect on the fundamental frequency of

_{Q}*I*(

*t*) is negligible.

Fig. 4(b) shows *f _{max}*(

*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

*f*increases with

_{max}*ρ*/

_{BD}*e*when

*R*= 50 Ω;

*f*increases by a factor of 1.3 when

_{max}*ρ*/

_{BD}*e*is increased to 10

^{13 }cm

^{−2}(and

*f*= 32 GHz) from

_{max}*ρ*/

_{BD}*e*= 0 (

*f*= 26 GHz).

_{max}To quantify the possible benefits of lattice alignment, Fig. 5(a) shows the effect of changing *θ* on *I _{b}*(

*V*). As

_{b}*θ*increases, the position of the current peak shifts to higher

*V*. The peak current amplitude,

_{b}*I*, decreases as

_{peak}*θ*increases due to increasing misorientation of the spatial parts of the wavefunction, see Fig. 5(b), so that

*I*could be ∼10× larger for an aligned device. However, for undoped samples, the PVR increases with increasing

_{peak}*θ*, 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*= 10

^{13 }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,

*R*( $\u221d(fmaxs)\u2009\u22122$; Eq. (5)) decreases with decreasing

_{N}*θ*, Fig. 5(c) inset, and with increasing

*ρ*, meaning that oscillation frequencies are highest for

_{BD}*θ*= 0° and when

*ρ*= 10

_{BD}/e^{13}cm

^{−2}. Fig. 5(d) shows that perfect alignment could increase

*f*by a factor of ∼2, i.e., for

_{max}*R*= 50 Ω,

*f*= 65 GHz when

_{max}*θ*= 0° compared to 32 GHz when

*θ*= 0.9°. The numerical results diverge from the small signal analysis power law of

*f*∝

_{max}*R*

^{−0.5}as

*R*becomes small, see black curve of Fig. 5(d), and it becomes necessary to vary

_{N}*V*to induce oscillations.

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 *I _{b}*, 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 RTDs

^{1}have similar

*I*and

_{b}*V*values as the GRTD reported here. We therefore expect that the GRTD will produce similar EM emission power (∼10

_{b}*μ*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.