Negative differential resistance and characteristic nonlinear electromagnetic response of a Topological Insulator

Materials exhibiting negative differential resistance have important applications in technologies involving microwave generation, which range from motion sensing to radio astronomy. Despite their usefulness, there has been few physical mechanisms giving rise to materials with such properties, i.e. GaAs employed in the Gunn diode. In this work, we show that negative differential resistance also generically arise in Dirac ring systems, an example of which has been experimentally observed in the surface states of Topological Insulators. This novel realization of negative differential resistance is based on a completely different physical mechanism from that of the Gunn effect, relying on the characteristic non-monotonicity of the response curve that remains robust in the presence of nonzero temperature, chemical potential, mass gap and impurity scattering. As such, it opens up new possibilities for engineering applications, such as frequency upconversion devices which are highly sought for terahertz signal generation. Our results may be tested with thin films of Bi2Se3 Topological Insulators, and are expected to hold qualitatively even in the absence of a strictly linear Dirac dispersion, as will be the case in more generic samples of Bi2Se3 and other materials with topologically nontrivial Fermi sea regions.

The enhanced nonlinearity of the response of TI surface states can be understood as follows. A TI heterostructure has two conducting surfaces, the top and bottom surfaces, while a Graphene sheet only has one. Due to the existence of the substrate in the TI heterostructure, structural inversion symmetry has to be broken, leading to a breaking of the degeneracy of the two TI surface states which opens up a Rashba-type spin splitting 26 . This results in an unique Dirac ring bandstructure which exhibits a much stronger nonlinear electromagnetic response than a Dirac cone alone, thereby opening up a venue for interesting physics as well as potential applications.
In this work, we model a TI heterostructure as a Dirac ring system, and analytically and numerically study its semiclassical nonlinear response. We first consider the case of an ideal Dirac ring, i.e. at zero mass gap, temperature and impurities. Next we consider deviations from these ideal conditions, and crucially show that the characteristic features of the response curve remain robust. We further discuss how this semiclassical analysis can be generalized to a more general setting with scattering and/or Rashba-like dispersion, and its implications for the output signal. Finally, we discuss some experimental proposals and engineering applications.

Results
We first introduce some basic theory on the semi-classical electromagnetic response of a generic Hamiltonian. With that, we present our main results on the characteristic nonlinear response curve of Dirac ring systems. Such systems have been detected in the surface states of thin films of Bi 2 Se 3 TIs via ARPES experiments 26 , and we will return to discussing the experimental signatures of our results after developing its theory.
Theory of semiclassical response. Consider a generic system described by a Hamiltonian ( ) H p under the influence of a driving field ( ) E t . At the semi-classical level, the field shifts the crystal momenta p of the partially occupied bands, leading to an induced current is the canonical velocity and ε ( ) p is the eigenenergy for a particular band. In other words,  J is the expectation value of the current  e v over states weighted by the time-dependent occupation function ( , ) g p t . As will be shown rigorously in Sect. 6, the time-dependent occupation function takes the form of the Fermi-Dirac distribution ( , ) = ( − ( )) g p t F p p t 0 0 , but with momentum shifted by an effective driving impulse . We shall elaborate on exactly how E eff depends on E and the scattering time τ in Sect. 6. For now, we shall proceed by treating p 0 as an external influence, and note that = E E eff in the limit of zero scattering. This approach considers only intra-band transport processes, and is valid when p 0 originates from an oscillatory electric field with  µ Ω , ,  T m max { }, i.e. with frequencies under 20 THz in typical applications. Henceforth, we shall work in units where  π = 2 1 and k B = 1 for notational simplicity. If the bands of the hamiltonian are isotropic with eigenenergies ε ε ( ) = ( ) p p where = p p , F 0 depends only on p via ε(p) and we can further express  J as is the contribution from the p-momentum shell. This decomposition is particularly useful since ( )  j p can often be expressed in closed-form, although  J usually cannot be. We now specialize to the Hamiltonians with a Dirac ring, which is the focus of this work. The energy dispersion of a Dirac ring system takes the form (Fig. 1) . Note that it reduces to the Dirac cone in Graphene when M = 0, which was systematically studied in ref. 19. In TI heterostructure realization reported and analyzed in refs 26-29, m is the mass gap induced by interlayer coupling and M measures the extent of inversion symmetry breaking.
Nonlinear response of Ideal Dirac rings. We shall first study the proptotypical case of ideal Dirac rings, where the gap m and temperature T are both zero. In the limit of small chemical potential µ  M, which can be tuned by varying the gate voltage 29 , the filled states on the valence band form a very thin ring bounded by inner and outer Fermi momenta with radial thickness 2μ/v F . Being such a thin ring, its total current is thus well-approximated by ( / )  j M v F in Eq. 3, which is simple enough to visualize schematically and exactly evaluate analytically (Eq. 15).
The canonical velocity is simply given by which is a vector of constant magnitude v F . It is always pointing in the radial direction, and is positive(negative) outside(inside) the Dirac ring. Such a scenario occurs when the filled states forms a non simply-connected 'ring-like' shape instead of a 'blob-like' shape. In this case, the lower dimensionality of the ring may allow for some novel kind of 'destructive interference' to occur between the contributions on both sides of the ring, and hence lead to negative differential resistance. Before proceeding with the calculations, let us attempt to understand that more intuitively.
Without a driving field, = p 0 0 and the ring of filled states lie exactly above the ring of Dirac nodes. On either side of it, the vector field ( )  v p points in equal and opposite directions, thereby resulting in a zero net current in Eq. 33. and 3. Upon a small impulse = | | p p 0 0 from the driving field, the ring of filled states will be slightly displaced in reciprocal (momentum) space, leading to an imbalance between the contributions of ( ) p from inside and outside the ring. As illustrated in Fig. 2,  J is relatively large for small p 0 because its arises from  v contributions that point in the same horizontal direction. As | | p 0 increases to more than half the radius of the ring, contributions from inside the ring oppose those from outside the ring, thereby leading to a decrease in  J . Finally, for p 0 larger than the radius of the ring, the filled states disentangle from the Dirac ring completely, and the  v contributions point in same horizontal direction again, adding up more strongly than before. As p 0 continue to increase,  v will eventually become parallel, leading to a maximal current  J proportional to the size of the ring. The above arguments suggest a response curve =  J J that rises sharply to a moderately large value when p 0 is very small, decreases when < < M v p M 2 F 0 , and increases to an even larger value for larger p 0 . We identify the decreasing region as the region of negative differential resistance ∝ < 0 . This agrees exactly with analytic expression derived in section and plotted in Fig. 3: Left) The bandstructure of a Dirac ring system ±ε(p) given by Eq. 4, with asymptotic canonical velocity (slope) given by v k F . Right) 3D view of the same band structure, with the axis of rotation the energy axis and the plane of rotation the p x , p y plane. The occupied (red) band is separated from the partially occupied (blue) band by a mass gap 2m. They can be obtained by gapping out two intersecting Dirac cone with energy difference of 2M. A ring of filled states (red) at chemical potential μ > m lies in the ring of minima on the valence (blue) band.
Scientific RepoRts | 5:18008 | DOI: 10.1038/srep18008 π α π ϕ α π ϕ α π ϕ α π ϕ α π α µ Distortion of a sinusoidal signal. Below, we show how an ideal Dirac ring system distorts periodic signals of different amplitudes. Due to the segment of negative differential resistance in the response curves e v is relatively large as it is contributed by  v whose horizontal components add constructively (light blue). Middle) As p 0 increases to the range , the ring of filled states disentangles with the ring of Dirac nodes (grey). The contributions from  v (light blue) add even more strongly, resulting in an even greater J. Very importantly, the behavior described here still holds for a ring thickened by (not too large) chemical potential or temperature, and even remain qualitatively true in the presence nonuniformities in  v , i.e. from nonzero mass. We now clearly see the asymptotic properties of the varies curves. As μ decreases, the response curve becomes 'sharper' due to the decreasingly thickness of the ring of filled states, and eventually converges nonuniformly to G(Q, Q) (from Eq. 16) as µ0.
Scientific RepoRts | 5:18008 | DOI: 10.1038/srep18008 in Fig. 3, additional kinks and lobes are introduced in the output signal  J . These lead to larger high-frequency components than what can be obtained with Graphene 19 , whose distorted output signal is a square-wave which also appears in the ∞  A 0 limit of the Dirac ring (Fig. 4). The extent of the lobes in the output can be quantified by the Fourier coefficients of the output current . When an energy scale is introduced by either nonzero gap or temperature, f 2n−1 decays exponentially. Even then, the frequency multiplication factor still outperforms that of Graphene for the the first few harmonics.
Response of imperfect Dirac ring materials. So far, we have considered Dirac ring systems with vanishing mass gap m, such as those realized in Bi 2 Se 3 TI thin films heterostructures thicker than 6 quintuple layers (QLs) 26 . To demonstrate the robustness of the regime of negative differential resistance, we shall now analyze Dirac ring systems which possess a gap and correspondingly a departure from perfect linear dispersion. They arise in sufficiently thin TI films where hybridization of the surface states on either side of the TI opens up a gap due to wavefunction overlap 27,28 . This had been predicted 27,28 to occur and was indeed observed 26 in Bi 2 Se 3 heterostructures thinner than 6 QLs.
Real materials under experimental conditions furthermore experience non-negligible temperature effects and impurity scattering, both of which can undermine the preceding Dirac ring interference analysis. But very importantly, we shall show that the qualitative features of the response curve, particularly  (Fig. 3). The negative differential resistance part of the response curve produces the interesting lobes in the output current. They will be smoothed out at larger μ, nonzero gap or nonzero temperature. In the → ∞ A 0 limit, the Dirac ring becomes irrelevant and the response curve reduce to that of Graphene (∝ sgnp 0 ) in the same limit.
the region of negative differential resistance, remain robust. As long as the occupied states still occupy a ring-like region in reciprocal space, we indeed observe: 3. An increasing, even larger J for larger p 0 .
To understand why, it is useful to examine Fig. 2 again. The non-monotonicity of the response (Center) is a consequence of the destructive interference of the contributions of  v from inside and outside the ring. This is a generic feature for a ring of energy minima, since the ε = ∇ ( ) v p p , which is always of opposite sides of the ring. In particular, note that it does not depend on the dispersion being linear, although the destructive interference will be less pronounced when the ring is thickened by large μ or m, or fuzzied by nonzero temperature T.
Below, we shall substantiate the above schematic arguments with detailed analyses and numerical results.
Effect of nonzero band gap m. In a gapped Dirac ring system, Eq. 4 gives the canonical velocity in a gapped Dirac ring system as which vanishes linearly near the Dirac ring. As such, we expect a more 'rounded' response curve, as shown in Fig. 6 (Left). Note that Eq. 8 also describes the small p 0 response in a Rashba system. The latter, however, has a asymptotically quadratic dispersion which leads to a (rather uninteresting) linear response for large > / p m v F 0 , which also suppresses the region of negative differential resistance. The results for nonzero mass gap m are shown in Fig. 6 (Left). The response curves exhibit the same qualitative shape, but with a broader linear response regime for small v p M F 0 , which is studied in more detail in Sect. . Examples of such cases include the bulk states of HgTe 34 and GaAs 35 quantum wells with inversion symmetry breaking.

Effects of nonzero temperature T.
In the presence of nonzero temperature, the Dirac ring is smeared out over a radius of ∆ ≈ / p T v F . Furthermore, hole carrier from the ε − ( ) p band also participate in current temperature T, which has the maximum possible frequency multiplication factor without negative differential resistance. It is exceeded by that of the gapless, zero-temperature Dirac ring at several different harmonics. This is expected from the appearance of additional, higher-frequency lobes in Fig. 4, where the sharper lobes in the A 0 = 2 curve now manifest themselves as a suppression of the first higher harmonic f 3 relative to the A 0 = 1 case. At nonzero T or m, the higher harmonics are exponentially suppressed. However, the n = 2 or n = 3 multiplication factor can still exceed the maximal possible without negative differential resistance (black line), as well as that of Graphene at T = 300 K, μ = 0 under the same input signal as given in the discussion section.
transport. The smearing of the Dirac ring results in qualitatively similar modifications to the response curve as increasing μ or m. However, the participation of the hole carriers undermines the destructive interference of  v described in Fig. 2, and can destroy the characteristic non-monotonicity and thus negative differential resistance of the response curve at sufficiently high T. Nonetheless, we numerically find (Fig. 6 (Right)) that the ideal response curve retains its qualitative shape till room temperature (300K) for ≈ . Here, we shall model the effect of scattering by a classical relaxation time τ, and show that its effects can be incorporated into our calculations by replacing the input driving field E with a 'renormalized' effective driving field E eff .
The electronic states are distributed according to a time-dependent occupation function g (p(t), t) that obeys the Boltzmann equation where = ( ( )) = ( ( )) = ( + ) p is the local equilibrium electron distribution, F 0 being the Fermi-Dirac distribution which can deviate from g = g(p(t), t), the true, non-equilibrium distribution. Being the equilibrium distribution, the functional form of ( ) f p remains unchanged, and all fluctuations in time occurs only through the time-varying momentum of an electron ( ) p t . By contrast, g can deviate from this equilibrium through its explicit time dependence, although it will relax towards f with a characteristic time τ. We will assume spatial homogeneity throughout.
Eq. 9 has an explicit solution though elementary calculus: where the dependence on ( ) p t has been made implicit. Intuitively, dg dt is a Laplace transform of df dt , with contributions from earlier times exponentially suppressed due to scattering. As derived in Sect. , Eq. 10 has an explicit solution: ( ) There is essentially no difference in the curves at T = 1 K and T = 77 K (liquid nitrogen temperature).
Scientific RepoRts | 5:18008 | DOI: 10.1038/srep18008 is the effective electric field that contains the effect of scattering. Essentially, E eff is contributed by past increments of the electric field E that are exponentially suppressed with a characteristic time τ. As τ increases, the increments have more time to constructively interfere. Physically, | | E eff increases as τ increases because a longer ballistic motion contributes greater to the momentum shift of each particle.
The manifest causality structure results from the time-reversal asymmetry of the system cf. Eq. 9.
For a simple sinusoidal driving field ( ) = Ω E t E t sin 0 , Eq. 12 (or Eq. 34) can be solved for the exact effective field . The latter is trivial for a single frequency mode, but will lead to interference effects when there is a mixture of modes. That will further studied in Sect. , in the limits of weak and strong scattering. Note that an expression very similar to Eq. 14 also appeared in ref. 19 in the context of radiation damping. But here τ has a more generic interpretation, encompassing all generic damping mechanisms. The DC limit can be recovered by taking Ω − π  t 1 2 , so that ( ) ≈ E t E 0 . The effective impulse is equal to eE t 0 in the ballistic regime, but saturates at τ eE 0 when t is increased to τ and effect of scattering is felt.

Discussion
Through an approach based on the Boltzmann equation, we have analytically and numerically obtained the novel nonlinear response curve of a Dirac ring system representing Topological Insulator thin films. Due to the special topology of a ring-shaped region of filled states, the response curve is intrinisically non-monotonic, leading to a regime with negative differential resistance which heightens frequency upconversion of a periodic signal. This property remains robust in real, experimentally fabricated samples in the terahertz gap, with nonzero mass gap, chemical potential, temperature and impurity scattering rate.
One test of our results will be to reproduce the output signals in Fig. 4. To produce the A 0 = 1 curve, for instance, we can input either an electrical or optical sinusoidal signal 37  and frequency Ω = 100 GHz. The output current should closely agree with that in Fig. 4 at the temperature of liquid nitrogen T = 77 K, and deviate slightly from it at room temperature T = 300 K (cf. Fig. 6(Right)). All other output signal shapes in Fig. 4 may be reproduced by varying the input signal amplitude. Since the effective impulse from the field increases linearly with the scattering time τ, TIs with much larger τ will exhibit the aforesaid negative differential resistance and nonlinear response behavior at much weaker physical electric fields. This will be very attractive for low power applications, and is well on track to becoming a reality with the rapid development of TI material growth techniques as well as the discovery of novel TI materials.

Methods
The current response of a gapless Dirac ring. Here we derive the analytic expression for the response of an (ideal) gapless Dirac ring, for which the µ0 result in Eq. 6 is a special case.
Substituting the canonical velocity which is just Eq. 6. For larger μ, numerical integration yields the plots in Fig. 3 which has the following properties: • The current is notably proportional to v F . This property is unique to the quasi-1D shape of the ring. While the velocity operator  v is proportional to v F , the area of the ring is independent of it. This is because its radius is proportional to it, while its thickness is inversely proportional to it.
• The response at any finite μ does not converge uniformly to Eq. 6: At small perturbations p 0 , the ring feels opposite velocity fields at both sides of the Fermi surface (FS), and there is a resultant linear regime. From Eq. 6, one can show that to a very high degree of accuracy, as shown in Fig. 3. Next comes a nonlinear regime where the response is negative approximately in the range μ < 0.6v F . It can also be shown that The effective driving field due to scattering. Here we show that the non-equilibrium state distribution g is of the form Care has to be taken in handling these PDEs: while g depends on t both explicitly and implicitly through ( ) p t , the dependence of f on t is only implicit through p and p d , respectively before and after integrating out the damping effect. To motivate the solution to Eq. 10, we first solve it for a periodic driving field. For each fourier component Ω , we have = − ( ) = − Ω Ω eE t eE e dp dt i t , so that