CROSS-REFERENCE TO RELATED APPLICATION(S)

[001] This application claims priority, under 35 Ti.S.C. § 119(e), to U.S. Application No. 62/779,025, filed December 13, 2018, which is incorporated herein by reference in its entirety.

GOVERNMENT SUPPORT STATEMENT

[002] This invention was made with Government support under Grant No. DE-SC0001088 awarded by the Department of Energy (DOE), Grant No. FA9550-16-1-0382 awarded by the Air Force Office of Scientific Research (AFOSR), and Grant Nos. DMR-0819762, DMR- 1231319 and DMR-1809815 awarded by the National Science Foundation (NSF). The Government has certain rights in the invention.

BACKGROUND

[003] Due to the explosive growth of wireless technologies and portable devices, there is now increasing effort towards developing micro-scale devices that can harvest ambient energy into usable electrical energy. Because of fast-developing microwatt and nanowatt electronics and next-generation wireless networks, there is a large demand for energy harvesters that can harvest energy from electromagnetic fields in the microwave and terahertz (THz) frequency range. High-frequency rectifiers can also be used in sensor and detector technology for the infrared, far infrared and submillimeter bands, which has wide-ranging applications in medicine, biology, climatology, meteorology, telecommunication, astronomy, etc. However, there is a so-called terahertz gap (0.1 to 10 THz) between the operating frequencies of electrical diodes and photodiodes. At frequencies within this range, efficient detection technology remains to be developed.

[004] The physical process central to harvesting electromagnetic energy is rectification, which refers to the conversion from an oscillating electromagnetic field to a direct current (DC). Existing rectifiers operating at radio frequencies are mostly based on electrical circuits with diodes, where the built-in electric field in the semiconductor junction sets the direction of the DC current. But diodes face two fundamental limitations when used to rectify an alternating current (AC). First, rectification by a diode requires a threshold input voltage V _T = k _B T /e, where k _B is Boltzmann’s constant, T is temperature, and e is the charge of an electron, This voltage is known as the thermal voltage and is about 26 meV at room temperature. Second, the responsivity is limited by the transition time in diodes, which is typically on the order of nanoseconds and drops at high frequencies, limiting the operating frequency of diode-based rectifiers.

[005] Instead of using a semiconductor junction, rectification can be realized as the nonlinear electrical or optical response of a noncentrosymmetric crystal. In particular, the second-order nonlinearity c(w) is an intrinsic material property that characterizes the DC current generated by an external electric field oscillating at frequency w. Rectification in a single homogenous material is not limited by the thermal voltage threshold or the transition time that limits the responsivity of a semiconductor junction. Moreover, the nonlinear electrical and optical responses of metals and degenerate semiconductors are much faster than photothermal effects used in bolometers for thermal radiation detection. This means that a rectifier based on a noncentrosymmetric crystal could be used to sense high-frequency (e.g., THz) radiation at higher speeds than a bolometer.

[006] However, the second-order nonlinearity of most materials is small. While nonlinear optical properties of quantum materials, such as photocurrent and second-harmonic generation, have been studied extensively, much less is known about second-order responses at radio, microwave, or infrared frequencies. Recent works have predicted intra-band photocurrents and second-order nonlinear Hall effects due to“Berry curvature dipoles” in nonmagnetic materials at zero magnetic field. In particular, the nonlinear Hall effect is predicted to be prominent in materials with tilted Dirac or Weyl cones, which are sources of large Berry curvature dipoles. Very recently, this effect was observed in low-frequency (around 100 Hz) transport measurements on the two-dimensional (2D) transition metal dichalcogenide (TMD) 1T’- WTe ₂ . The second-order Hall conductivity of bilayer WTe ₂ is remarkably large, in agreement with its large Berry curvature dipole from the tilted Dirac dispersion. This and other types of nonlinear response from intra-band processes are also being explored in topological insulator surface states, Weyl semimetals, Rashba systems and heavy fermion materials.

SUMMARY

[007] Here, we present current rectifiers that use the nonlinear electric and terahertz response of noncentrosymmetric crystals. Unlike other current rectifiers, the inventive current rectifiers do not include diodes or other (semiconductor) junction structures. Instead, they operate by exploiting the second-order nonlinear response of a patch of noncentrosymmetric material. This patch of noncentrosymmetric material is coupled to an antenna, which receives an alternating current (AC) electric field. The antenna couples this AC electric field into the patch of noncentrosymmetric material, which generates a net DC current from the AC electric field via second-order skew scattering. A pair of electrodes electrically connected to the patch of noncentrosymmetric material conduct this DC current to a load.

[008] Without being bound by any particular theory, the large and highly tunable rectification effects in graphene multilayers, heterostructures, and 2H-TMD monolayers can be explained from the intra-band second-order response in time-reversal-invariant noncentrosymmetric materials using a semiclassical Boltzmann equation and taking account of electron Bloch wavefunctions via quantum-mechanical scattering rates. Skew scattering with nonmagnetic impurities contributes to rectification. This skew scattering arises from the chirality of Bloch wavefunctions in momentum space and does not require spin-orbit coupling. The contributions to rectification from skew scattering and from Berry curvature dipole scale differently with the impurity concentration, and skew scattering predominates at low impurity concentration. Moreover, skew scattering is allowed in all noncentrosymmetric crystals, whereas Berry curvature dipole requires more strict symmetry conditions.

[009] An example current rectifier may include a patch of a noncentrosymmetric material that is electrically coupled to an antenna and to a pair of electrodes. In operation, the antenna receives an alternating current (AC) electric field and couples it into the noncentrosymmetric material. The noncentrosymmetric material generates a direct current (DC) current via a second-order nonlinear response to the AC electric field. And the electrodes conduct the DC current out of the noncentrosymmetric material.

[0010] The AC electric field may oscillate at a frequency from about 0.1 THz to about 10 THz. The AC electric field may oscillate in a first direction in a plane of the noncentrosymmetric material and the DC current may flow from the noncentrosymmetric material in a second direction orthogonal to the first direction and in the plane of the noncentrosymmetric material. Alternatively, the AC electric field may oscillate and the DC current may flow in the same direction in a plane of the noncentrosymmetric material.

[0011] In either case, the noncentrosymmetric material may lack symmetry in a direction perpendicular to a plane of the noncentrosymmetric material. The noncentrosymmetric material may have a band gap of about 26 meV at room temperature (e.g., about 300 K). The noncentrosymmetric can generate the DC current from the AC electric field at room temperature.

[0012] The noncentrosymmetric material may include a graphene monolayer disposed on a dielectric layer that breaks an inversion symmetry of the graphene monolayer. Alternatively, the noncentrosymmetric material may include a graphene bilayer, in which case the current rectifier may further include first and second insulators and first and second gate electrodes. The first insulator layer is disposed on a first side of the graphene bilayer, and the second insulator layer is disposed on a second side of the graphene bilayer opposite the first side of the graphene bilayer. The first gate electrode is on a far side of the first insulator layer from the graphene bilayer and applies a first gate voltage to the graphene bilayer. And the second gate electrode is on a far side of the second insulator layer from the graphene bilayer and applies a second gate voltage to the graphene bilayer. The first and second gate voltages break an inversion symmetry of the graphene bilayer. The noncentrosymmetric material may also comprise a transition metal dichalcogenide.

[0013] All combinations of the foregoing concepts and additional concepts discussed in greater detail below (provided such concepts are not mutually inconsistent) are part of the inventive subject matter disclosed herein. In particular, all combinations of claimed subject matter appearing at the end of this disclosure are part of the inventive subject matter disclosed herein. The terminology used herein that also may appear in any disclosure incorporated by reference should be accorded a meaning most consistent with the particular concepts disclosed herein.

BRIEF DESCRIPTIONS OF THE DRAWINGS

[0014] The skilled artisan will understand that the drawings primarily are for illustrative purposes and are not intended to limit the scope of the inventive subject matter described herein. The drawings are not necessarily to scale; in some instances, various aspects of the inventive subject matter disclosed herein may be shown exaggerated or enlarged in the drawings to facilitate an understanding of different features. In the drawings, like reference characters generally refer to like features (e.g., functionally similar and/or structurally similar elements).

[0015] FIGS. 1A and IB show a high-frequency (e.g., THz) current rectifier with a noncentrosymmetric quantum material that converts an oscillating electric field into a transverse DC current. [0016] FIG. 2 shows a high-frequency (e.g., THz) current rectifier with a noncentrosymmetric quantum material that converts an oscillating electric field into a longitudinal DC current.

[0017] FIG. 3 A is a side view a noncentrosymmetric quantum material formed of a monolayer graphene/dielectric heterostructure where the monolayer graphene’s inversion symmetry is broken by a dielectric layer.

[0018] FIG. 3B is a top view of the monolayer graphene/di electric heterostructure of FIG. 3 A.

[0019] FIG. 4 is a side view of a noncentrosymmetric quantum material formed of a bilayer graphene stack whose inversion symmetry is broken by an applied bias electric field.

[0020] FIG. 5 A illustrates transverse ( X _yxx ) skew scattering on a honeycomb (e.g., graphene) lattice.

[0021] FIG. 5B illustrates the electric field and rectified current for transverse (x _yxx ) skew scattering on a honeycomb lattice.

[0022] FIG. 5C illustrates Fermi surface displacement at frequency w (green) for transverse ( X _y xx ) skew scattering on a honeycomb lattice.

[0023] FIG. 5D illustrates stationary Fermi surface displacement for transverse (X _yxx ) skew scattering on a honeycomb lattice.

[0024] FIG. 6 A illustrates longitudinal (x _yyy ) skew scattering on a honeycomb (e.g., graphene) lattice.

[0025] FIG. 6B illustrates the electric field and rectified current for longitudinal (x _yyy ) skew scattering on a honeycomb lattice.

[0026] FIG. 6C illustrates Fermi surface displacement at frequency w (green) for longitudinal ( X _yyy ) skew scattering on a honeycomb lattice.

[0027] FIG. 6D illustrates stationary Fermi surface displacement for longitudinal (x _yyy ) skew scattering on a honeycomb lattice.

[0028] FIG. 7A is a plot of the second-order conductivity c for a graphene monolayer heterostructure.

[0029] FIG. 7B is a plot of the reduced voltage responsivity h _n for a graphene monolayer heterostructure. [0030] FIG. 7C is a plot of the second-order conductivity c for a graphene bilayer heterostructure.

[0031] FIG. 7D is a plot of the reduced voltage responsivity for a graphene bilayer heterostructure.

[0032] FIG. 8A illustrates the lattice structure and crystal axes of monolayer graphene.

[0033] FIG. 8B illustrates the lattice structure and crystal axes of bilayer graphene.

[0034] FIG. 8C shows the Brillouin zone and energy contours for the monolayer graphene model of FIG. 8 A.

[0035] FIG. 9A is a plot of the Fermi energy e _F for monolayer graphene (upper trace) and bilayer graphene (lower trace).

[0036] FIG. 9B is a plot of the Fermi velocity v _F for monolayer graphene (flat trace) and bilayer graphene (sloped trace).

[0037] FIG. 9C is a plot of the linear DC conductivity s at zero temperature for monolayer graphene (sublinearly increasing trace) and bilayer graphene (linearly increasing trace).

[0038] FIG. 10A is a plot of the linear conductivity s of monolayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 (bottom trace) to 5 X 10 ¹² cm ^-2 (top trace) and a temperature of zero degrees Kelvin.

[0039] FIG. 10B is a plot of the second-order conductivity c of monolayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 to 5 X 10 ¹² cm ^-2 and a temperature of zero.

[0040] FIG. IOC is a plot of the reduced current responsivity 77 _/ of monolayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 to 5 X 10 ¹² cm ^-2 and a temperature of zero.

[0041] FIG. 10D is a plot of the reduced voltage responsivity h _n of monolayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 to 5 X 10 ¹² cm ^-2 and a temperature of zero.

[0042] FIG. 10E is a plot of the linear conductivity s of bilayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 (bottom trace) to 5 X 10 ¹² cm ^-2 (top trace) and a temperature of zero degrees Kelvin. [0043] FIG. 10F is a plot of the second-order conductivity c of bilayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 to 5 X 10 ¹² cm ^-2 and a temperature of zero.

[0044] FIG. 10G is a plot of the reduced current responsivity 77 _/ of bilayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 to 5 X 10 ¹² cm ^-2 and a temperature of zero.

[0045] FIG. 10H is a plot of the reduced voltage responsivity h _n of bilayer graphene versus frequency for carrier densities ranging from 0.5 X 10 ¹² cm ^-2 to 5 X 10 ¹² cm ^-2 and a temperature of zero.

[0046] FIG. 11A is a linear-scale plot of the linear conductivity s of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0047] FIG. 11B is a linear-scale plot of the second-order conductivity c of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0048] FIG. l lC is a linear-scale plot of the reduced current responsivity 77 _/ of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0049] FIG. 1 ID is a linear-scale plot of the reduced voltage responsivity h _n of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0050] FIG. 11E is a linear-scale plot of the linear conductivity s of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0051] FIG. 1 IF is a linear-scale plot of the second-order conductivity c of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0052] FIG. 11G is a linear-scale plot of the reduced current responsivity 77 _/ of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0053] FIG. 11H is a linear-scale plot of the reduced voltage responsivity h _n of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0054] FIG. 12A is a logarithmic-scale plot of the linear conductivity s of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0055] FIG. 12B is a logarithmic-scale plot of the second-order conductivity / of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 . [0056] FIG. 12C is a logarithmic-scale plot of the reduced current responsivity 77 _/ of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0057] FIG. 12D is a logarithmic- scale plot of the reduced voltage responsivity h _n of monolayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0058] FIG. 12E is a logarithmic-scale plot of the linear conductivity s of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0059] FIG. 12F is a logarithmic- scale plot of the second-order conductivity c of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0060] FIG. 12G is a logarithmic-scale plot of the reduced current responsivity 77 _/ of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

[0061] FIG. 12H is a logarithmic-scale plot of the reduced voltage responsivity h _n of bilayer graphene versus frequency and temperature for a carrier density n = 1 X 10 ¹² cm ^-2 .

DETAILED DESCRIPTION

[0062] FIGS. 1A and IB illustrate a current rectifier 100 that can rectify radio-frequency and terahertz electric fields (e.g., alternating current (AC) electric fields oscillating at 0.1 THz to 10 Hz). The current rectifier 100 includes an antenna 100 that is coupled to a thin patch of noncentrosymmetric material 120, which in turn is coupled to a pair of electrodes 130a and 130b (collectively, electrodes 130). The electrodes 130 connect the rectifier 100 to a load 140, which may have an impedance of ohms to kiloohms. In principle, there is no high-frequency coupling to these electrodes 130, but they may be DC-coupled (e.g., using a low-pass filter) to remove high-frequency electric field components and pass DC current.

[0063] In this example, the patch of noncentrosymmetric material 120 is rectangular, with two halves of the antenna 110 are coupled to one pair of opposing edges and the electrodes 130 coupled to the other pair of opposing edges. The antenna 110, noncentrosymmetric material 120, and electrodes 130 can be deposited on any suitable dielectric substrate (e.g., a silicon oxide substrate). Notably, this current rectifier 100 does not include or form a diode.

[0064] The antenna 110 is shown in FIGS. 1 A and IB as a bowtie antenna but could take any other suitable shape. It is made of metal or another conductive material and is configured to detect linearly polarized incident electric fields. For example, the bowtie antenna 110 has a long axis stretching horizontally in FIG. IB and is configured to detect a horizontally polarized electric field, as indicated by the double-headed arrow representing the electric field oscillation direction in FIG. IB. If the incident electric field is not linearly polarized along the bowtie antenna’s long axis— for example, if the incident electric field is unpolarized or elliptically polarized— then the bowtie antenna 110 couples the horizontal projection of the incident electric field into the noncentrosymmetric material 120.

[0065] The noncentrosymmetric material 120 can include a two-dimensional material, such as monolayer graphene, bilayer graphene, multi-layer graphene, or a transition metal dichalcogenide, or a three-dimensional material, such as a bulk transition metal dichalcogenide or a Weyl semimetal without inversion (e.g., tantalum arsenide (TaAs)). It is formed in a roughly planar rectangular patch whose dimensions can be about 1-10 pm wide to about 1-10 pm long to about 10-100 nm thick. The patch may be formed in other shapes (e.g., another polygon) so long as the noncentrosymmetric material 120 is coupled to the antenna 110 and electrodes 130. As explained in greater detail below, this patch may include multiple layers, including dielectric layers and/or additional electrodes that break the inversion symmetry of a two-dimensional material.

[0066] The noncentrosymmetric material 120 lacks inversion symmetry along an axis parallel to the DC current. This axis is generally, but not always, perpendicular to the noncentrosymmetric material’s largest surface (i.e., along the propagation direction of a normally incident electric field). For the transverse response illustrated in FIG. IB, the noncentrosymmetric material 120 lacks inversion symmetry along an axis perpendicular to the electric field oscillation direction.

[0067] In operation, the antenna 110 couples an incident electric field E _w into the noncentrosymmetric material 120, which exhibits a second-order nonlinear response to the incident electric field E _w as explained in greater detail below. More specifically, skew scattering by the lattice of the noncentrosymmetric material 120 produces a net DC current / _DC . In the rectifier 100 shown in FIGS. 1A and IB, this scattering is transverse— that is, the DC current / _DC propagates in a direction that is perpendicular to the oscillation direction (indicated by the double-headed arrow in FIG. IB) of the electric field E _{w .} The current direction is set by the polarization of the electric field E _w and the crystal asymmetry. The DC current / _DC flows within the layer of noncentrosymmetric material 120 to the load 140 via the electrodes 130. [0068] FIG. 2 shows a current rectifier 200 that generates a“longitudinally” propagating rectified current that / _DC from a high-frequency (e.g., 0.1-10 THz) electric field E _{w .} Like the current rectifier 100 in FIGS. 1A and IB, this current rectifier 200 includes an antenna 210; a thin, rectangular patch of noncentrosymmetric material 220; and electrodes 230a and 230b (collectively, electrodes 230) on a substrate 202. The electrodes 230 are electrically connected to a load 240, which may also be on the same substrate 202 as the current rectifier 200. In this case, however, the antenna elements and electrodes 230 are coupled to the same sides of the noncentrosymmetric material 220.

[0069] In operation, the antenna 210 couples the electric field E _w into the noncentrosymmetric material 220 such that the electric field E _w oscillates in the direction indicated by the double headed arrow in FIG. 2. Again, the noncentrosymmetric material 220 generates a net DC current / _DC in response to the electric field E _w via a second-order skew scattering response. In this case, however, the DC current / _DC propagates in a direction parallel to the electric field E _w oscillation direction. The DC current / _DC flows within the plane of noncentrosymmetric material 220 to the load 240 via the electrodes 230.

[0070] FIGS. 3 A, 3B, and 4 show different examples of noncentrosymmetric heterostructures suitable for use as the noncentrosymmetric materials in the current rectifiers shown in FIGS. 1 A, IB, and 2. More specifically, FIGS. 3A and 3B show side and top views, respectively, of a noncentrosymmetric heterostructure 320 formed of monolayer graphene 322 on a layer of dielectric boron nitride (BN) 324 on a substrate 302. This heterostructure 320 may be up to about 100 nm thick. By itself, the monolayer graphene 322 has inversion symmetry. The lattice mismatch between the monolayer graphene 322 and the BN layer 324 breaks the monolayer graphene’s inversion symmetry, causing the heterostructure 320 to be noncentrosymmetric. (A layer of silicon carbide (SiC) can be used instead of the BN layer 324 to break the monolayer graphene’s inversion symmetry.)

[0071] FIG. 4 is a side view of a noncentrosymmetric heterostructure 420 with bilayer graphene 422 sandwiched between a pair of BN layers 424a and 424b (collectively, BN layers 424), which are themselves sandwiched between a pair of electrodes 426a and 426b (collectively, electrodes 426). Each BN layer 424a, 424b may be about 10-30 nm thick (e.g., 20 nm) and can have different thicknesses (i.e., the upper BN layer 424a can be thicker or thinner than the lower BN layer 424b). Together, the bilayer graphene 422, BN layers 424, and electrodes 426 form a double capacitive structure on a substrate 402. [0072] Like monolayer graphene, the bilayer graphene 422 has inversion symmetry. This inversion symmetry is broken by applying bias fields V _BG and V _TG to the upper and lower layers of the bilayer graphene 422 with the electrodes 426. These bias fields are orthogonal to the surface of the substrate 402, to the polarization of the incident (AC) electric field, and to the rectified current. The bias fields generate a displacement field with the graphene bilayer 422 that are roughly proportional to (V _BG /d ₁ — V _TG /d ₂ ), where ci ₁ is the thickness of the upper BN layer 424a and d ₂ is the thickness of the lower BN layer 424b. The displacement field amplitude can range from about -1.5 V/nm to about +1.5 V/nm for applied voltages on the order of 2 V. The applied voltages and BN layer thicknesses are also selected based on the desired carrier density within the bilayer graphene 422. This carrier density is proportional to (V _BG /di + V _TG /d ₂ ) and may be about ±10 ¹³ cm ^-2 .

[0073] FIGS. 5A-5D and 6A-6D illustrate the second-order response by skew scattering on a honeycomb lattice for the longitudinal (jc _yxx ) and transverse (x _yyy ) cases, respectively. FIGS. 5A and 6A show the schematics of skew scattering. When a self-rotating wave packet (large circles) is scattered by an inversion-symmetric potential (small circle), its motion is deflected like the Magnus effect. Two wave packets moving in the opposite directions produces zero net current. However, if the two self-rotate in different directions, skew scattering produces net current in the perpendicular direction.

[0074] FIGS. 5B and 6B show the electric field and rectified current on a honeycomb lattice. There are two schemes for measuring a strong rectification effect: the AC oscillating field can be along either the zigzag direction (FIG. 5B) or armchair direction (FIG. 6B) and the DC rectification signal is collected along the armchair direction.

[0075] FIGS. 5C and 6C shows Fermi surface displacement at frequency w (solid outline). The oscillating electric field E(w) forces the Fermi surface to swing back and forth from its equilibrium position (dashed outline). Owing to the Fermi surface anisotropy, each valley yields finite velocity V along the k _x direction after time averaging. This velocity is canceled with the two valleys, and there is no DC current generated as a linear response.

[0076] FIGS. 5D and 6D show the stationary Fermi surface displacement. The electric field and skew scattering produce the stationary Fermi surface displacement (solid outline) from the equilibrium state (dashed outline) as a second-order response. Finite rectified current is observed when the contributions from the two valleys do not cancel.

Results Second-Order Response

[0077] We study the DC current in a homogeneous material generated at second order by an external electric field E of frequency w . The second-order response also involves a 2w component, corresponding to the second-harmonic generation, which we do not focus on in this work. We write down the second-order DC response as

Ja— abc^b^c 0)

where x _abc is the rank-three tensor for the second-order conductivity, E _a = E _a (w) is the (complex) amplitude of the external electric field of frequency w, and the indices a, b, c label the coordinates. x _abc satisfies x _abc = x _a ^* _{cb .} Since the current j _a is odd under inversion while the electric field E _a is even, finite x _abc is possible only in noncentrosymmetric media. The second-order response can be decomposed into two parts as j _a = X _abc 'E _b E _c + ix _ab c"E _b E _c , with x _abc ' = x _acb ' and X _abc " = ^~ Xacb " · The former describes the response to a linearly- polarized field and the latter to a circularly-polarized field. It is important to note that in an isotropic medium second-order response to a linearly-polarized field must vanish. In other words, a nonzero c' _abc requires the presence of crystal anisotropy.

[0078] Nonlinear responses are studied extensively in the optical frequency regime. The classical approach to nonlinear optics considers electrons bound to a nucleus by an anharmonic potential, so that the restoring force to the displacement of electrons becomes nonlinear. In quantum-mechanical theory, the nonlinear optical response at frequencies larger than the band gap is usually dominated by electric dipole transitions between different bands. On the other hand, with energy harvesting and infrared detection in mind, here we consider intra-band nonlinear response at frequencies below the inter-band transition threshold.

Semiclassical Boltzmann Analysis

[0079] We analyze the second-order electrical transport using the semiclassical Boltzmann theory. For a homogeneous system, the Boltzmann equation is given by (h = 1):

¾ + k -¾ = -e . (¾

where /( k, e) is the distribution function at energy e and C[f] denotes the collision integral. The time derivative of the wavevector is equal to the force felt by an electron k =— eE under the external electric field E. The distribution function / deviates from the equilibrium Fermi- Dirac distribution when electrons are accelerated by the external field. A nonequilibrium steady state is obtained by the balance of the acceleration and the relaxation due to scattering, described by the collision integral C[f\, whose explicit form is shown below in a general form. It includes the scattering rate w _k , _k , the probability of the transition per unit time from a Bloch state with a wavevector k to that with k'. We consider scattering by impurities. The impurity density should be sufficiently lower than the electron density for the semiclassical Boltzmann analysis to be valid; see below for further discussions.

[0080] To obtain a nonzero second-order conductivity j _a¾c , we should capture the effect of inversion breaking. For time-reversal-invariant systems, the energy dispersion e _k and the band velocity v ₀ (k) = V _k e _k satisfy the conditions e _k = e_ _k and v ₀ (k) =— v ₀ (— k) , which hold with or without inversion symmetry; that is, the absence of inversion is encoded in the wavefunction but not in the energy dispersion.

[0081] One consequence of inversion asymmetry is manifested in a scattering process, as the transition rate depends not only on the dispersion but also on the Bloch wavefunction ip _k . For and only for noncentrosymmetric crystals, the transition rate from k to k' and from— k to— k' can be different: w _k , _k ¹ w_ _k , _ _k . A scattering process with such asymmetry is referred to as skew scattering. Skew scattering is a source of the second-order DC response, i.e., rectification. The strength of skew scattering is characterized by = (w _k , _k — w_ _k , _ _k )/2, whereas the symmetric component w _{k k} = (w _k , _k + w_ _k , _ _k ) /2 does not rely on inversion breaking. The presence of time-reversal symmetry guarantees w _k , _k = w_ _k _ _k , , which allows us to write

[0082] The transition rate w _k , _k is usually calculated to the lowest order in scattering potential by Fermi’s golden rule. For elastic impurity scattering, it is given by = where V _k , _k is the matrix element of a single scatterer and ( ) denotes the impurity average. However, this formula is symmetric under the exchange k <-> k' and does not capture skew scattering. The latter arises at the next-leading order from the interference of a direct transition k ® k' and an indirect process k ® k' via an intermediate state q : which satisfies w[^ _k =— w _kk ^. Here we use the notation J _q = / d ^d q/(2n) ^d ( d is the spatial dimension).

[0083] Skew scattering provides an extrinsic contribution to the anomalous Hall effect and the spin Hall effect. These are linear response phenomena in systems with broken time-reversal and spin-rotational symmetries, respectively. Skew scattering exists and contributes to the second-order response even in materials with both time-reversal and spin-rotational symmetry. With regard to its microscopic origin, skew scattering can be caused by nonmagnetic impurities when the single-impurity potential V (r) is inversion asymmetric, such as a dipole or octupole potential. Even for a symmetric single-impurity potential V(r) = V(— r) , such as a delta function potential, skew scattering can also appear. In such cases, skew scattering is due to the chirality of the electron wavefunction in a noncentrosymmetric crystal. In the simplest case of a delta function impurity potential, we have F _k , _k oc (it _k , |ii _k ) where |u _k ) is the Bloch wavefunction within a unit cell. Then the expression given in Eq. (3) involves the wavefunction overlap at different k points on the Fermi surface, similar to the Berry phase formula, and also skew scattering is known to be related to Berry curvature in a certain limit.

[0084] While there is no net chirality in nonmagnetic systems, the Bloch wavefunction |u _k ) of an electron at finite momentum can be complex and carry k-dependent orbital magnetization, spin polarization, or Berry curvature. In other words, mobile electrons in inversion-breaking nonmagnetic crystals can exhibit momentum-dependent chirality that is opposite at k and— k. Then, skew scattering occurs due to electron chirality in a similar way as the classical Magnus effect, where a spinning object is deflected when moving in a viscous medium. The spinning motion is associated with the self-rotation of a wave packet and the viscosity results from scattering with impurities; see FIGS. 5 A and 6 A.

[0085] To obtain the second-order conductivity, the distribution function /(k) should be calculated to second order in the external electric field E(w). Moreover, since skew scattering is parametrically smaller than w^ ^s we expand the distribution function up to first order in w^ ^A Solving the Boltzmann equation self-consistently at each order yields the second-order DC current response j =— e J _k v ₀ (k)/(k) due to skew scattering in time- reversal-invariant systems. The second-order conductivity x _abc can be formally expressed in terms of various scattering times t _h , t' _{h .}

where t _h (w) is a shorthand for t _h (w) = t _h /(1— ίwt _h ) and likewise for t' _h (w) . The scattering times t _h , t' _h with n = 1,2 are associated with the dominant symmetric scattering see below for details. They determine the steady-state distribution function up to second order in the electric field, when skew scattering is neglected. A detailed discussion of a self- consistent solution of the Boltzmann equation and Joule heating appears below.

[0086] A rough order-of-magnitude estimate ofc can be obtained from Eq. (4) by using two scattering times t and f for the symmetric and skew scattering, respectively. Then, the second- order DC current of a metal or degenerate semiconductor under a linearly-polarized electric field is given by

with At = t for wt « 1 and At = l/w for wt » 1. Here v _F is the Fermi velocity, p _F = hk _F is the Fermi momentum, and n is the electron density.

[0087] This estimate provides a heuristic understanding of the second-order response. In the low-frequency limit w 0, the ratio of the second-order currentj ₂ and the linear response currentj = s _w E is a product of two dimensionless quantities eEx/p _F and t/t. The first term is the change of electron’s momentum under the external field during the mean free time, divided by the Fermi momentum. The second term t/t characterizes the strength of skew scattering which is responsible for second-order electrical response, relative to the symmetric scattering.

[0088] The short-circuit current responsivity ¾ _/ is a metric for a rectifier, defined as the ratio of the generated DC current to the power dissipation. With a sample dimension L ² , the current responsivity is given by ¾ _/ º j ₂ L/(j ₁ EL ² ). Similarly, the voltage responsivity 3i _v is defined for the generated DC voltage as 3i _v º (y ₂ f /s ₀ )/(ji^L ² ). Both ¾ _/ and 3i _v are independent of the magnitude of the external electric field. Using the linear response conductivity s _w ~ env _F (eX _FT /h)Re[(l— ίwt) ^-1 ] obtained from the same Boltzmann equation gives the short-circuit current responsivity and voltage responsivity:

[0089] Here we define the reduced current and voltage responsivities 77 _/ and h _n , respectively, which are independent of the sample size and the ratio t/t. The approximate relations for h _I and h _n hold in both low- and high-frequency limits. We also consider the ratio of the generated DC power and input power h _R º (j ₂ / ^s o) /( _& F ² ), which characterizes the power conversion efficiency (assuming the load resistance and internal resistance are comparable). We find

These figures of merit show that a material with low carrier density n or small Fermi momentum p _F is desirable for efficient rectification. Moreover, a longer mean free time leads to a higher current responsivity. In this respect, high-mobility semiconductors and semimetals are promising as rectifiers for infrared detection. In the following, we analyze graphene systems and estimate their efficiency of rectification.

Graphene Multilayers and Transition Metal Dichalcogenides

[0090] Pristine graphene has high mobility and low carrier density, but it is centrosymmetric. Nevertheless, as a van der Waals material, we can assemble multilayer stacks to break inversion. Realizations of noncentrosymmetric structures include monolayer graphene on a hexagonal boron nitride (hBN) substrate, electrically-biased bilayer graphene, and multilayer graphene, such as ABA-stacked trilayer graphene. Anisotropy of energy dispersions, also used for second-order response, naturally arises from trigonal crystal structures.

[0091] Skew scattering arises from the chirality of quantum wavefunctions and contributes to the second-order conductivity in graphene heterostructures with gaps induced by inversion symmetry breaking. The k p Hamiltonian considered here is

where s = ±1 denotes two valleys at K and K', respectively, and we define k+ = k _x ± ik _{y .} This Hamiltonian has a band gap of 2D in the energy spectrum and describes two graphene systems: (1) monolayer graphene on hBN, where the spatially varying atomic registries break the carbon sublattice symmetry and opens up gaps at Dirac points, and (2) bilayer graphene with an out-of-plane electric field that breaks the layer symmetry and opens up gaps. The difference between the monolayer and bilayer cases lies in the relative strength of v and A. For the monolayer (bilayer) case, Ak _F « vk _F (vk _F « Ak _F ) is responsible for the trigonal warping of the linear (quadratic) energy dispersion.

[0092] The Bloch wavefunction of H _s ( k) has two components associated with sublattice (layer) degrees of freedom in monolayer (bilayer) graphene. When the band gap is present ( D ¹ 0), the wavefunction in each valley is chiral and carries finite Berry curvature, leading to skew scattering from nonmagnetic impurities. The chirality and Berry curvature are opposite for valley K and K' due to time reversal symmetry. Because of the three-fold rotation symmetry of graphene, Berry curvature dipole vanishes, leaving skew scattering as the only mechanism for rectification.

[0093] We calculate the second-order conductivity from Eq. (4), and find nonvanishing elements c _ccg = c _cgc = c _gcc =— x _yyy º c, consistent with the point group symmetry. The relation among the electric field, induced current, and underlying crystalline lattice is depicted in FIGS. 5B and 6B, along with the schematic pictures of the Fermi surface displacement in FIGS. 5C, 5D, 6C, and 6D. The function c has the form / = e ³ v _F (x ³ /c)z(e _R , w), with the dimensionless function z given below. The ratio x/x is proportional to the gap D and the Fermi surface trigonal warping. The gap D describes the effect of inversion breaking on the electronic structure, while the Fermi surface trigonal warping is responsible for Fermi surface asymmetry within a valley E _s ( k) ¹ E _s (— k).

[0094] FIGS. 7A-7D show the second-order conductivity c and the reduced voltage responsivity h _n for both monolayer and bilayer gapped graphene as functions of frequency v = w/(2p). FIGS. 7A and 7C show the second-order conductivity c for monolayer graphene and bilayer graphene, respectively. FIGS. 7B and 7D show the reduced voltage responsivity h _n for monolayer graphene and bilayer graphene, respectively. The carrier density n is changed from 0.5 x 10 ¹² cm ^-2 to 5 x 10 ¹² cm ^-2 . We use the values v = 0.94 x 10 ⁶ m/s, D = 15 meV, and x = 1.13 ps for the monolayer case; D = 50 meV and x = 0.96 ps for the bilayer case, from transport, infrared spectroscopy, and scanning tunneling microscopy/spectroscopy measurements. For the bilayer case l = (2 m) ^-1 is determined by the effective mass m » 0.033m _e (m _e : electron mass) and v « 1 x 10 ⁵ m/s.

[0095] The frequency dependence of the second-order conductivity is qualitatively understood with Eq. (5). We can further simplify the relation in two dimensions for hl ² _r = constant. The second-order conductivity then behaves as c ~ e ³ v _F x ² /h ² · ( / ) for wc « 1 and c ~ e ³ v _F /(ho)) ² (t/t) for wc » 1.

[0096] The responsivity is affected by the linear conductivity s as well as the second-order conductivity c. At high frequencies with the carrier density fixed, we observe a decrease of the rectified current while the energy dissipation by the linear response also decreases. Since both s and c decrease as w ^-2 , the responsivity saturates in the high frequency limit. With a sample size of 10 mth , carrier density n = 0.5 x 10 ¹² cm ^-2 , and ratio x/x = 0.01 , reasonable estimates of mean free time and impurity concentration give a saturated voltage responsivity ¾iz ~ 30 V/W for monolayer graphene and 40 V/W for bilayer graphene (see below for the frequency dependence). As expected from the qualitative estimate Eq. (6), high responsivity is realized with low carrier density. The semiclassical approach is valid when the mean free path is much longer than the Fermi wavelength.

[0097] Both monolayer and bilayer graphene show broadband response. The saturation of responsivities lasts up to the onset frequency of an inter-band transition (not included in our study). For example, at the lowest density of n = 0.5 X 10 ¹² cm ^-2 shown in FIGS. 7A-7D, the Fermi energy for the monolayer case is e _F ~ 80 meV, so that due to the Pauli blocking the inter-band transition does not occur at frequencies below 2e _F ~ 40 THz . The inter-band threshold frequency also increases with the band gap, which is tunable in bilayer graphene by the displacement field.

[0098] Finite temperature affects the intra-band second-order response and responsivity through thermal smearing of the distribution function and the change of scattering times. Here, thermal smearing does not change the result much as the Fermi energy is much higher than room temperature. Since it utilizes material’s intrinsic nonlinearity, rectification from the intra band process considered here does not suffer from the noise associated with thermally excited electron-hole pairs in photodiodes.

[0099] Our analysis also applies to 2H-TMD monolayers, where transition metal and chalcogen atoms form trigonal crystal structures. Similar to gapped graphene, their band structures can also be described as massive Dirac fermions with trigonal warping, but with much larger band gaps. Hence, they should exhibit large intra-band second-order response as well. Topological insulator surface states also show this effect (see below). In addition to 2D systems, three-dimensional bulk materials without inversion are also worth considering. As discussed, they should have large responsivities in low carrier densities. From this respect, inversion-breaking Weyl semimetals could also be used in place of graphene for current rectification.

Discussion

[00100] So far, we have considered the skew scattering contribution to the second-order DC response by neglecting the effects from the Berry curvature W directly. It contributes to the distribution function through the collision integral and the electron’s velocity as anomalous and side-jump velocities. Hence the second-order conductivity is also modified. The Berry curvature is odd (even) under time reversal (inversion), and thus is allowed to exist in time- reversal-invariant noncentrosymmetric materials. The effect of Berry curvature to linear transport properties has been well studied in the context of the anomalous Hall effect. Finite Berry curvature enforces a coordinate shift <5r _kk , after scattering, giving a displacement of the center of an electron wave packet. It is accompanied by a potential energy shift <5e _k , _k in the presence of the external field, which results in the collision integral

e \f] = J _k , K, _k /(k, e _k ) - w _kk ,/(k', e _k + d e _k , _k )] . (9)

[00101] One may decompose the distribution function in nonequilibrium as / = / ₀ + / ^scatt + / ^adlst , where / ₀ is the equilibrium distribution function, / ^scatt describes the effect of scattering to the distribution (with _< 5e _k , _k = 0), and / ^adlst is the anomalous distribution due to finite <5e _k , _k. Because of the different origins, these distributions have distinct dependence on the scattering time t. For low frequencies with wt « 1, the t dependence of each term can be estimated by a simple power counting. Noting that the scattering rate w amounts to t ^_1 , we can expand the distribution function in powers of the electric field E , with coefficients depending on t:

in the weak impurity limit.

[00102] The electric current density is obtained by j =— e J _kv (k)/(k), where v(k) is the electron’s velocity,

[00103] Besides the velocity given by the band dispersion v ₀ , v contains two additional terms: The anomalous velocity v _av is associated with the Berry curvature and the last term v _sj describes the side-jump contribution, which arises from a combined effect of scattering and Berry curvature. Now, we can see several contributions to the second-order conductivity c from combinations in the product of the out-of-equilibrium distribution function (/ ^scatt + / ^adlst ) and the electron velocity (v ₀ + v _av + v _sj ).

[00104] The dominant contribution in the weak impurity limit (t ® oo) emerges from the skew scattering with v ₀ / ₂ ^scatt , resulting in c oc t ² . The Berry curvature dipole contribution corresponds to v _av fi ^catt , which amounts to c oc t. It shows the frequency dependence as c' oc e ³ r D , c" » 0 for low frequencies (wt « 1) and c' oc , c" oc b ³ ΰ/w for high frequencies (wt » 1), where D is the Berry curvature dipole. The Berry curvature dipole is a quantity with dimensions of length in a two-dimensional system; it is dimensionless in a three- dimensional system. The Berry curvature dipole and skew scattering have the same frequency dependence although they are distinct in the t dependence.

[00105] Inversion breaking results in finite skew scattering . Along with the accompanying anisotropic Fermi surface, the second-order electronic response induced by skew scattering should persist in any noncentrosymmetric material. This is in marked contrast to the Berry curvature dipole mechanism, which imposes more symmetry constraints in addition to inversion breaking. Some symmetry classes without inversion, e.g., those containing three-fold rotation axis, do not show nonlinear response from a Berry curvature dipole.

[00106] In summary, the analysis above is a systematic study of second-order electrical response due to an intra-band process. Skew scattering makes the dominant contribution in the weak impurity limit. Skew scattering is responsible for current rectification and opens a new way to low-power energy harvesters and terahertz detection.

Semiclassical Boltzmann Theory

[00107] We describe in detail the analysis of the semiclassical Boltzmann transport theory for the calculation of second-order response in noncentrosymmetric materials. The following analysis is not meant to be limiting and to some extent is parallel with that for the anomalous Hall effect. However, a fundamental difference is found in symmetry; for a second- order response, inversion is broken whereas time-reversal symmetry is broken for the anomalous Hall effect.

[00108] The total Hamiltonian is given by

H = H ₀ + V + U,

where H ₀ is the Hamiltonian for electrons in the clean limit with the eigenvalue e _k , V describes electron scattering, and U includes the external electric field, which drives the system into nonequilibrium. In a semiclassical description, U results in the force acting on an electron’s wave packet: V _r U =— F. In the following, we consider an external electric field of frequency w, namely, F =— eE(u>).

[00109] A formally rigorous analysis involves a quantum-mechanical calculation of density matrices, which usually requires considerable effort. The semiclassical analysis is simpler and offers intuitive understanding; the downside is that it may have limitations owing to the uncertainty principle. The semiclassical description relies on wave packets of Bloch states. A wave packet is localized both in the real space and the momentum space so that the mean position r and momentum k can be designated. However, the uncertainty principle imposes the limitation ArAk > 1. We desire position and momentum resolutions set by the mean free path - and the Fermi wavevector k _F ; i.e., k _F £ » 1 is assumed for the semiclassical analysis. In the following, we deal with electron scattering by impurities. Since f and k _F depend on the impurity and electron densities, respectively, the condition is satisfied when the impurity density is sufficiently lower than the electron density.

[00110] We deal with spatially homogeneous systems and hence the electron distribution function /(k, e ) does not depend on the spatial position, but on the wavevector k and the energy e. We assume that the frequency of the external field w is lower than the inter band spacing, so that inter-band transitions are suppressed and negligible. We calculate the distribution function in the presence of the external field by the semiclassical Boltzmann equation ( h = 1) i ^{+ k} £= -wl· < ¹² >

The time derivative of the momentum is equal to the force felt by an electron wave packet: k = F. The collision integral C[f] is given by

C[f] = f [ w _k , _k /( k, e) - w _kk ,f(k', e + 5e _k , _k )],

'k

where w _k , _k is the scattering rate from a state with momentum k to one with k' . In the following, we consider elastic scattering. Then, the scattering rate is defined by

w _{\ k} — 2TT| T _klk |<S(e _k ,— e _k ),

with the scattering T matrix T _k , _k. It is given by T _k , _k = (k'\V\ip _k ), where |fc) is the eigenstate of H _Q , the Hamiltonian in the clean limit, and | p _k ) is the eigenstate of H ₀ + V , including the scattering. | p _k ) can be obtained as the solution to the Lippman-Schwinger equation \ip _k ) =

[00111] The energy shift _< 5e _k , _k occurs after a scattering owing to the Berry curvature W. The Berry curvature is defined from the Berry connection A as

where | p _k ) is the lattice periodic part of the Bloch function at k and e _abc is the Levi-Civita symbol. The Berry curvature transforms under inversion T and time reversal T as : /2 _a (k) ® 12 _a (— k), T : W _a (k) ® -/2 _a (-k).

[00112] It is deduced that the Berry curvature vanishes everywhere in the Brillouin zone when PT symmetry exists. Finite Berry curvature causes the coordinate shift at a scattering process, denoted by <5r _k , _k , which describes the displacement of the center of a wave packet by a scattering k ® k'. The energy shift <5e _k , _k is given by using the coordinate shift <5r _k , _k as

<^ ^e _k , _k ^=— F <^r , .

For a weak scattering with a small momentum change, the coordinate shift can be approximated as <5r _k , _{k «} (k'— k) x W( k). The coordinate shift is independent of the scattering time t even though it is related to scattering; it describes the displacement after a single impurity scattering.

[00113] We decompose the distribution function / as

f = f ₀ + / ^scatt + / ^adist .

The first term / ₀ is the distribution function in equilibrium, i.e., the Fermi-Dirac distribution. The second term / ^scatt describes the scattering contribution without the Berry curvature and the last term / ^adlst is the anomalous distribution due to the Berry curvature and the energy shift

Scattering Rate in Noncentrosymmetric Media

[00114] In noncentrosymmetric media, the probabilities of a scattering process k ® k' and the inverted process— k ®— k' can be different. For time-reversal systems, this imbalance is captured by decomposing the scattering rate into two parts:

where the symmetric and antisymmetric parts and w^ ^A respectively, are defined by

When time reversal is preserved, the scattering rate satisfies the reversibility w _kk , = w_ _k , _ _k , which leads to the relations

Those equalities do not hold in general when time reversal is broken since states at k and— k usually have different energies. The optical theorem for elastic scattering guarantees the relation / _k , ^{) =} 0. (13)

When the elastic scattering is due to impurities, are obtained to the lowest order in the impurity scattering potential by

where ( ) denotes the average over the impurity distribution and F _k , _k is the matrix element for the single impurity scattering. The symmetric part of the scattering rate is obtained at the lowest order of the Bom approximation, whereas the antisymmetric part is found at the next leading order. The detailed balance is broken when the antisymmetric part is finite.

[00115] V. I. Belinicher et al.,“The photogalvanic effect in media lacking a center of symmetry,” Sov. Phys. Usp. 23, 199-223 (1980), reviews the photogalvanic effect in noncentrosymmetric materials, considering the asymmetry of scattering. The photovoltaic effect is studied with the Boltzmann equation and the antisymmetric component of the scattering rate w,( _;k , similarly to our analysis. A difference can be found in the origin of finite w _k , _k The asymmetry may be imposed on scattering potentials for impurity scattering, ionization, photoexcitation, and recombination. The asymmetry of the scattering rate can be finite due to an asymmetric scattering potential and also a wavefunction of a noncentrosymmetric medium. For the latter case, a scattering potential does not need to be inversion asymmetric to have finite w _k , _k , but even an isotropic scattering potential can generate w _{k k} through a wavefunction.

Formal Solutions

[00116] The semiclassical Boltzmann equation (12) should be solved self-consistently to obtain the distribution function /( k). Here, we expand the distribution function f = f ₀ + f ^scatt + / ^adlst with respect to the electric field E(w) and the asymmetric part of the scattering rate

The subscript n and the superscript m correspond to the orders of the electric field E and the antisymmetric scattering rate w^ ^A respectively. As we have seen in Eqs. is smaller than so that we treat the former as a perturbation. Then, the Boltzmann equation is decomposed for each to become

A prime symbol ' added to the distribution functions stands for the derivative with respect to the energy: /'(k, e ) = d _e /( k, e).

[00117] Each equation should be solved self-consistently; however, this is usually not easy. Here, we solve a part of the collision integral involving self-consistently, while the other part with is treated as a perturbation. The collision integrals with define the scattering times as follows:

do not satisfy the self-consistency, one may regard the equations above as the definitions for the relaxation time approximation. The notation of the scattering times t _h and t' _h used in the main text; they correspond to t _h º [00118] We assume t and t' are independent of momentum k in solving the equations. Then, we obtain the formal solutions for the distribution functions and g^. The second- order response contains different frequencies, 0 and ±2w. Since now we focus on rectification, i.e., zero-frequency response, we write down only the second-order solutions with w = 0. Now, we list below the formal solutions up to the second order in the electric field (n < 2) and to the first order i

(17)

(18)

are the shorthand notations for

Those formal solutions find the scattering time dependence of the distribution functions for low frequencies (wt « 1) :

Higher-order corrections for the anomalous distribution (n > 2) have a nested structure of f _n -i and f _n -f- In the weak disorder limit of r ® ¥ we obtain

adist _{K T} n-l _E n We can define another scattering time f from the antisymmetric part w^ ^A It is roughly speaking given by

[00119] The definition of f is accompanied with the distribution function because the optical theorem for elastic scattering concludes Eq. (13).

Current Response

[00120] From the distribution function, the electric current response is obtained by j = -e fv( k)f(k).

'k

v(k) is the electron’s group velocity, given by

[00121] The group velocity from the energy band dispersion v ₀ is independent of the Berry curvature, while the anomalous velocity v _av =— k X W = eE X W and the side-jump velocity v _sj = / w _k , _k <5r _k , _k are induced by the Berry curvature. In terms of the electric field E, v _av is linear in E whereas v ₀ and v _sj do not depend on E. Lastly, only the side-jump velocity depends on the scattering time: v _sj oc t ^_1 . Since there is no current in equilibrium, the electric current is described by

j = -e / _k ( ^v o + ^V av + v _sj ) (/ ^scatt + / ^adist ), (20)

which allows us to decompose the current into contributions of different origins. The contribution from the anomalous distribution / ^adlst could be recognized as a part of the side- jump effect since it is originated both from the Berry curvature and scattering.

[00122] Before calculating the second-order conductivity, consider the problem of energy dissipation or Joule heating because of j · E ¹ 0. Previous works indicates that a stationary state solution cannot be obtained for a closed system with energy conservation. To circumvent this issue, we neglect parts of the distribution function which are isotropically coupled to the electric field. To be more explicit, the distribution function at the second order in the electric field involves a term proportional to |E| ² /" _0. Technically, such an excess term arises when we substitute a distribution function f _n (n > 1) into the collision integral C[f], which involves only elastic scattering. Since this term is isotropic, it cannot be relaxed by elastic scattering, which conserves energy, so that a stationary state solution does not exist. Inelastic scattering resolves the problem as it does not conserve energy and excess energy is dissipated as heat. An inelastic scattering time is typically much longer than an elastic scattering time. In calculating the second-order conductivity, we simply subtract and neglect such isotropic terms although they potentially change the temperature for a closed system. The isotropic terms do not contribute to current from the band velocity since J _k v ₀ /" ₀ = J _k

F _k /'o = 0.

[00123] We can see the six possible combinations in Eq. (20). For the second-order current response, the scattering time dependence of each contribution for low frequencies wt « 1 is

In the weak disorder limit ofr ® ¥ the first term, the skew scattering contribution, v ₀ / ^scatt predominates, which we focus on here. This is equivalent to neglecting effects related to Berry curvature. The second-order current response arising from this term is

We define the second-order conductivity by the relation

72, a ^~ Xabc^b^c·

and we obtain from Eqs. (17) and (18)

X _abc vanishes identically in time-reversal-invariant systems because of e _k = e_ _k ; it follows from v ₀ (k) =— v ₀ (— k), and we see that the integrand is odd in k for Therefore, the second-order contribution from the skew-scattering contribution arises from X^ _c , i.e., x _abc = y ('·i)

abc

[00124] To calculate the second-order response, it is useful to express Eqs. (21) and (22) in the following forms:

[00125] We also show the expression for for reference although it vanishes. These expressions show that results at finite temperature can be obtained by considering the thermal broadening of the Fermi-Dirac distribution, where the Fermi energy e _F is replaced by the chemical potential /i _{c cm} with the carrier density fixed. More explicitly, from the following relation

-/' _o (k) = f de[-f _o (e^ _chem )]5(e _k - e _F ),

the second-order conductivity at finite temperature is obtained by

X _abc ^ _chem> F) = / de[— f o( ^e J / _em )] t _abc ( ^e _> F = 0).

Linear Response

[00126] We also calculate the linear response to obtain responsivities (defined below). From the Boltzmann equation, the linear current response is given by

leading to the linear conductivity

Graphene-Based Models with Trigonal Lattice Structures

[00127] FIGS. 8A-8C illustrate explicit models for evaluating the current response. We consider two graphene-based models, utilizing monolayer and bilayer graphene. FIG. 8A shows the lattice structure of the monolayer graphene along with the crystal axes. FIG. 8B shows the lattice structure of bilayer graphene. The shaded circles depict distinct sublattice sites with different potential energies, and the filled and empty circles correspond to atoms on different layers. Both lattices have reflection planes parallel to the y axis and its symmetry partners under the three-fold rotation. FIG. 8A shows Brillouin zone and energy contours for the monolayer model. With small doping, there are trigonally-warped Fermi surfaces around K and K' points, facing the opposite directions. The bilayer model exhibits the similar Fermi surfaces as the monolayer ones.

[00128] The crystalline lattices for the monolayer and bilayer models belong to a crystallographic point group D _3h and C _3v , respectively. Inversion is equivalent to the in-plane two-fold rotation in 2D models, which is absent in the present models. Around K and K' points,

^K = (i¾ , °) a„d K' the k p Hamiltonian to second order is

where k _x and k _y are measured from K or K' point, labeled by s = ±1, and k _s is defined by k _s = k _x + isk _{y .} The energy dispersion of the positive branch is given by e _k = jv ² k ² + A ² k ⁴ + D ² — 2svAk ³ cos39,

and the corresponding normalized wavefunction is

The Berry curvature is also calculated as

sA (v ² — 4 A ² k ⁴ )

^{Ώz =} 2

[00129] The Hamiltonian Eq. (24) is available to both monolayer and bilayer models, which are shown in FIGS. 8 A and 8B, respectively. The Pauli matrix s _a describes the two sublattices for monolayer graphene and two layers for AB (Bernal) stacked bilayer graphene. Another difference is found in which of v and l is dominant, which is explained below. A describes the potential energy difference at the two sublattice sites or layers, which opens a gap to the energy spectrum. We consider the two cases separately below.

Monolayer Graphene

[00130] For monolayer graphene, the terms with v dominate those with l ( vk _F » kp), where l characterizes the trigonal warping of the Fermi surface s = +1 corresponds to K valley and s =— 1 to K' valley. The Dirac mass can be induced for example by placing a graphene sheet on a hexagonal boron nitride (hBN) substrate. For a tight-binding model on a honeycomb lattice with the nearest-neighbor hopping amplitude t and the distance between sites a, the velocity v and the trigonal warping l are given by v = 3ta/2 and l = 3ta ² /8, respectively.

[00131] We assume that the trigonal warping is small and obtain the energy dispersion e _k to first order in l as

where the Fermi wavevector k _F is given by

The density of states (DOS) at the Fermi energy is

^€ F

Do 2pn ² '

and 6 _k denotes the polar angle of k. The Fermi surface can be parameterized by 0 _k as

k ²

k(0 _k ) = k _F + sX— cos30 _k .

v

[00132] We use the wavefunction on the Fermi surface to evaluate the scattering rate. The normalized wavefunction on the Fermi surface is given by

Also, it is useful to see the Berry curvature on the Fermi surface, which is given by

SV ² A

W + 0(X ² ).

2e _F

The second-order conductivity is considered in the presence of an elastic impurity scattering. Here, we assume short-ranged random scalar impurities

but restrict scattering within each valley. V _t is the impurity potential strength, 71; is the impurity density, and I}· denotes an impurity position, where the summation is taken over all impurity sites. Then, we obtain the scattering rates for the valley s to first order in l as

+1^ [sin(0 - 46f) + sin(40 - Q’) - sin(0 + 2Q’) + sin(20 + <9')]}5(e _k - _<¾ ,), with Q = 0 _k and 0' = 0 _k , parametrizing the Fermi surface.

[00133] In order to obtain the distribution function, we determine the scattering times that satisfy the relation Eq. (16). For the present Hamiltonian, the scattering times are obtained as

[00134] These solutions are obtained to 0(l°) , which are enough to evaluate conductivities to order 0(/ 1 ¹ ) (see the next paragraph). Also, there is no directional dependence in the scattering times at order A ⁰ . For clarity, we define

and the dimensionless quantity

which are defined to satisfy = gt.

[00135] We calculate the second-order contribution from each valley separately, namely

Xabc = Xabc + Xabc ^with bC

The contributions from each valley s are obtained from Eqs. (21) and (22) at zero temperature as

where the factor of two from spin is multiplied. D _3h symmetry means that (i) the imaginary part of the second-order DC conductivity vanishes and (ii) the system does not have second- order response to a circularly-polarized field. We define another scattering time f, related to

This quantity measures the amounts of inversion breaking (m) and the Fermi surface warping ( l ). The ratio t/t becomes

[00136] Roughly speaking, the ratio is determined by the product of the following dimensionless quantities: inversion breaking D /e _F , the Fermi surface warping Ak _F /( k _F ), and the impurity strength Vik _F /e _{F .} We can see from Eq. (25) that the antisymmetric scattering becomes finite when inversion is broken by D ¹ 0 but it does not require finite Fermi surface warping l. However, the second-order response vanishes when the Fermi surface warping is absent. The ratio t/t measures the amount of the contribution to the second-order response from skew scattering, but not the magnitude of skew scattering itself.

[00137] is originated from the absence of time-reversal or inversion symmetry within each valley, resulting from the Fermi surface anisotropy. However, time reversal is preserved with the two valleys, and thus contributions from different valleys cancel. This effect of c ⁽⁰⁾ can be seen in FIGS. 5C and 6C, where the two valleys generates the opposite static velocity V. In contrast, c ⁽¹⁾ adds up with the two valleys, since this contribution captures the opposite chirality of the wavefunctions around each valley, which is not visible from the Fermi surface shape. Finally, we obtain the nonvanishing elements of the second-order conductivity

with the average Fermi velocity on the Fermi surface

and the dimensionless function Here the factor of four in Eq. (26) corresponds to spin and valley degrees of freedom. The finite elements of x _abc are consistent with D _3h symmetry of the lattice.

Linear Conductivity and Responsivity of Monolayer Graphene

[00138] To estimate the magnitude of the second-order response and to calculate responsivities, we also calculate the linear conductivity, particularly _. From Eq. (23), we obtain

for each valley with the spin degrees of freedom included. The linear conductivity is diagonal, and hence we write

Since the linear conductivity can be written as s = ne/r with the electron density n and the mobility m , in low frequencies we find the relation

V

m— bt— .

e _F

The electron density is obtained by n = k /p.

[00139] The short-circuit current responsivity ¾ _/ is defined as the ratio of the generated DC current j ₂ to the power of the incident electric field absorbed by the sample. When a sample has a dimension of L ² and the incident field is uniform on the sample, the current responsivity is given by

[00140] The dimension L corresponds to either the width or length of a sample, depending on whether the second-order response of interest is longitudinal or transverse to the incident electric field. We also define the voltage responsivity 3i _v as the ratio of the voltage generated by the second-order response to the incident power:

[00141] Both current responsivity ¾ _/ and voltage responsivity 3i _v depend on the sample dimension L and the ratio t/t . To remove those factors, we define the reduced current responsivity 77 _/ and the reduced voltage responsivity h _n as follows:

For the present model, 77 _/ and h _n are obtained as

Bilayer Graphene

[00142] Bemal stacked bilayer graphene has C _3v symmetry when the two layers have different potential energies. Because of the symmetry, we can use the same Hamiltonian Eq. (24); however, v is treated as a perturbation instead of l for the bilayer graphene case (. k _p » vk _F ), as we can see from a tight-binging Hamiltonian. Treating v as a perturbation, we obtain the energy dispersion e _k to first order in v as

The DOS at the Fermi energy D _Q is given by

The normalized wavefunction on the Fermi surface to order v is

where k _F is defined by

The Berry curvature on the Fermi surface is

Using the wavefunction, we obtain the symmetric and antisymmetric scattering rates and _w O ⁾ , respectively:

+— [sin(0 + 2Q’) - sin(20 + 0') + sin(20 - 50') + sin(50 - 20')]}5(e _k - <¾,).

2 A-k-p

For the bilayer graphene case, the scattering times become

We define t as

resulting in t^ ⁰ ·* = t and = g ^_1 t , with the same g for the monolayer graphene case.

[00143] The contributions to the second-order conductivity from valley s are obtained from Eqs. (21) and (22) at zero temperature as

For the bilayer graphene model, the ratio t/t is given by

T _ k _F V _t

ΐ ^Dn 8l(€ ² + D ² ) ^'

As in the monolayer graphene case, x °^ _bc from each valley cancels but X^ _bc adds, leading to the second-order conductivity

where for the bilayer model the average velocity on the Fermi surface v _F is and the dimensionless function z (w) becomes

Linear Conductivity and Responsivity of Bilayer Graphene

[00144] We also evaluate the linear conductivity for each valley s with the spin degrees of freedom included [Eq. (23)]:

As in the monolayer graphene case, we define

For low frequencies, the mobility is calculated from the relation s = hbm by

where the electron density is n = k _p /p.

[00145] For the bilayer case, the reduced current responsivity becomes

and the reduced voltage responsivity is

Estimate of Parameters

[00146] Parameters for a realistic evaluation are taken and estimated from experiments on monolayer and bilayer graphene. For monolayer graphene on an hBN substrate, the imbalance of the potential energy within a sublattice is induced by the substrate and the lattice mismatch of graphene and hBN creates a superlattice with a periodicity of about 13 nm. In addition to the original (primary) Dirac points, the Moire superlattice produces new (secondary) superlattice Dirac points, located at a different energy. Those two have different velocities and energy gaps: at ambient pressure, v = 0.94 X 10 ⁶ m/s and 2/3 = 30 meV for the original Dirac point, and v = 0.5 X 10 ⁶ m/s and 2/3 = 20 meV for the superlattice Dirac points. The size of the gap depends on the superlattice periodicity and it could be enhanced by interaction. Even at the original Dirac point, the velocity is slightly modified. We consider the original Dirac points instead of the superlattice Dirac points. The latter have the valley degeneracy of six and anisotropic Dirac cones whereas the valley degeneracy of the former is two.

[00147] The scattering time is estimated from a transport measurement with a high- quality sample. At 300 K, the mobility m » 90,000 cm ² /V s with the carrier density n = 1 X 10 ¹² cm ^-2 , which amounts to a scattering time t ~ 1.13 ps. The corresponding frequency is (2pt) ^-1 » 141 GHz and the mean free path is £ ~ n _R t » 1.6 mih. The Fermi wavelength at n = 0.5 X 10 ¹² cm ^-2 is l _{R »} 0.05 mih, which is still shorter than the mean free path: l _R « £. To estimate the impurity strength, we need to know the impurity concentration. From transport and scanning tunneling microscopy/spectroscopy, we approximate the impurity density as n _L * 1 X 10 ⁹ cm ^-2 , which results in V _{t »} 2.77 X 10 ^-13 eV cm ² ( nV _{t »} 277 meV at n = 1 X 10 ¹² cm ^-2 ). Since the velocity is only slightly modified from its value without a superlattice structure, we employ the value l/n = a/4 obtained from a tight-binding model, where a is the carbon atom spacing, given by a = 1.42 A. Those values yield the ratio t/t » 0.003 at n = 1 x 10 ¹² cm ^-2 . In the present calculations, we do not consider the temperature dependence of the parameters or the origin of scattering and inhomogeneity of samples for simplicity, although the scattering time becomes longer in lower temperatures.

[00148] For bilayer graphene, high-quality bilayer graphene grown by chemical vapor deposition (CVD) and detached on hBN has a mobility m » 30,000 cm ² /V s at 300 K. An out- of-plane electric field (displacement field) D opens a band gap of about 2D = 100 meV with D » 1.1 V/nm . The coefficient l = (2m) ^-1 is determined by the effective mass m » 0.033m _e ( m _e : electron mass) and v » 1 x 10 ⁵ m/s. Evaluating at n = 1 x 10 ¹² cm ^-2 , we find the scattering time t ~ 0.96 ps, corresponding to a frequency (27tt) ^-1 ~ 166 GHz and a mean free path £ ~ n _R t » 0.35 mih . We assume an impurity concentration » 1 X 10 ⁹ cm ^-2 , considering its high mobility, to obtain an impurity strength V _L * 1.06 X 10 ^-13 eV cm ² ( nV _{t »} 106 meV at n = 1 x 10 ¹² cm ^-2 ). At a carrier density n = 1 x 10 ¹² cm ^-2 , we have the ratio t/t » 0.018.

[00149] FIGS. 9A, 9B, and 9C show the carrier density dependence of e _P , v _F , and s(0) at zero temperature, respectively, for the parameter values immediately above.

[00150] FIGS. 10A-10H show the frequency dependence of s , c , and m _/ for monolayer graphene (FIGS. 10A-10D) and bilayer graphene (FIGS. 10E-10H) models. Different traces represent different carrier densities, which range from 0.5 X 10 ¹² cm ^-2 to 5 x 10 ¹² cm ² , and the temperature is set to be zero. See above for the values of the parameters for each model.

[00151] FIGS. 11A-11H and 12A-12H show the frequency and temperature dependence of the response on linear and logarithmic scales, respectively. The scattering times used for the calculations are those at 300 K. The scattering time becomes longer in lower temperature, which makes the linear conductivity s, the second-order conductivity c, and the reduced current responsivity 77 _/ larger, whereas the low- and high-frequency values of the reduced voltage responsivity h _n barely change. The carrier density is fixed at n = 1 X 10 ¹² cm ^-2 . FIGS. 11A-11D and 12A-12D correspond to the monolayer graphene model and FIGS. 1 IE-11H and 12E-12H correspond to the bilayer graphene model. We assume that the scattering time t is a constant independent of temperature; see above for the values of parameters.

Surface State of a Topological Insulator

[00152] A surface state of a topological insulator also lacks inversion and thus hosts the rectification effect. Here, we take an example of the surface state of the topological insulator Bi ₂ Te ₃ , where the Hamiltonian is

leading to the energy spectrum

The origin k = 0 is located at G point and the k _x axis is parallel to GK line. Unlike the previous two cases based on graphene, a topological insulator has a single Dirac spectrum centered at G point. Furthermore, the topology of the wavefunction in the bulk prohibits a mass gap of the surface states. In this case, a surface breaks inversion by itself and the symmetry is C3V ·

[00153] For the present Hamiltonian, short-ranged impurities without momentum dependence do not produce finite w^ ^A nor does the second-order response. This issue can be circumvented instead by considering Coulomb impurities, where the impurity potential is inversely proportional to the distance from the impurity position: 1 r— ry | ^_1 . Assuming that the hexagonal warping of the Fermi surface is small, the matrix element of the Coulomb impurity potential on the Fermi surface can be expressed as

Then, we obtain the scattering rates

with the numerical factor c _» 0.93 obtained from the numerical integration over momenta.

[00154] Since there is a single Fermi surface and time-reversal symmetry relates the Fermi surface with itself, c ⁽⁰⁾ vanishes with the single Fermi surface. A finite second-order conductivity can be found in as follows:

Here, we assume for simplicity the constant scattering time t ~ (27rD ₀ n _i F _c ² ) ¹ obtained from w^ ^s and the scattering time originated from is defined by

2n _t V _c ³ e

f ^® = l- pn 7

The two numerical factors c _x ~ 0.87 and c ₂ ~ 0.14 are evaluated by numerical integrations.

Conclusion

[00155] While various inventive embodiments have been described and illustrated herein, those of ordinary skill in the art will readily envision a variety of other means and/or structures for performing the function and/or obtaining the results and/or one or more of the advantages described herein, and each of such variations and/or modifications is deemed to be within the scope of the inventive embodiments described herein. More generally, those skilled in the art will readily appreciate that all parameters, dimensions, materials, and configurations described herein are meant to be exemplary and that the actual parameters, dimensions, materials, and/or configurations will depend upon the specific application or applications for which the inventive teachings is/are used. Those skilled in the art will recognize or be able to ascertain, using no more than routine experimentation, many equivalents to the specific inventive embodiments described herein. It is, therefore, to be understood that the foregoing embodiments are presented by way of example only and that, within the scope of the appended claims and equivalents thereto, inventive embodiments may be practiced otherwise than as specifically described and claimed. Inventive embodiments of the present disclosure are directed to each individual feature, system, article, material, kit, and/or method described herein. In addition, any combination of two or more such features, systems, articles, materials, kits, and/or methods, if such features, systems, articles, materials, kits, and/or methods are not mutually inconsistent, is included within the inventive scope of the present disclosure.

[00156] Also, various inventive concepts may be embodied as one or more methods, of which an example has been provided. The acts performed as part of the method may be ordered in any suitable way. Accordingly, embodiments may be constructed in which acts are performed in an order different than illustrated, which may include performing some acts simultaneously, even though shown as sequential acts in illustrative embodiments.

[00157] All definitions, as defined and used herein, should be understood to control over dictionary definitions, definitions in documents incorporated by reference, and/or ordinary meanings of the defined terms.

[00158] The indefinite articles“a” and“an,” as used herein in the specification and in the claims, unless clearly indicated to the contrary, should be understood to mean“at least one.”

[00159] The phrase“and/or,” as used herein in the specification and in the claims, should be understood to mean“either or both” of the elements so conjoined, i.e., elements that are conjunctively present in some cases and disjunctively present in other cases. Multiple elements listed with“and/or” should be construed in the same fashion, i.e.,“one or more” of the elements so conjoined. Other elements may optionally be present other than the elements specifically identified by the“and/or” clause, whether related or unrelated to those elements specifically identified. Thus, as a non-limiting example, a reference to“A and/or B”, when used in conjunction with open-ended language such as“comprising” can refer, in one embodiment, to A only (optionally including elements other than B); in another embodiment, to B only (optionally including elements other than A); in yet another embodiment, to both A and B (optionally including other elements); etc.

[00160] As used herein in the specification and in the claims,“or” should be understood to have the same meaning as“and/or” as defined above. For example, when separating items in a list,“or” or“and/or” shall be interpreted as being inclusive, i.e., the inclusion of at least one, but also including more than one, of a number or list of elements, and, optionally, additional unlisted items. Only terms clearly indicated to the contrary, such as“only one of’ or“exactly one of,” or, when used in the claims,“consisting of,” will refer to the inclusion of exactly one element of a number or list of elements. In general, the term“or” as used herein shall only be interpreted as indicating exclusive alternatives (i.e.,“one or the other but not both”) when preceded by terms of exclusivity, such as“either,”“one of,”“only one of,” or “exactly one of.” “Consisting essentially of,” when used in the claims, shall have its ordinary meaning as used in the field of patent law.

[00161] As used herein in the specification and in the claims, the phrase“at least one,” in reference to a list of one or more elements, should be understood to mean at least one element selected from any one or more of the elements in the list of elements, but not necessarily including at least one of each and every element specifically listed within the list of elements and not excluding any combinations of elements in the list of elements. This definition also allows that elements may optionally be present other than the elements specifically identified within the list of elements to which the phrase“at least one” refers, whether related or unrelated to those elements specifically identified. Thus, as a non-limiting example,“at least one of A and B” (or, equivalently,“at least one of A or B,” or, equivalently“at least one of A and/or B”) can refer, in one embodiment, to at least one, optionally including more than one, A, with no B present (and optionally including elements other than B); in another embodiment, to at least one, optionally including more than one, B, with no A present (and optionally including elements other than A); in yet another embodiment, to at least one, optionally including more than one, A, and at least one, optionally including more than one, B (and optionally including other elements); etc.

[00162] In the claims, as well as in the specification above, all transitional phrases such as “comprising,” “including,” “carrying,” “having,” “containing,” “involving,” “holding,” “composed of,” and the like are to be understood to be open-ended, i.e., to mean including but not limited to. Only the transitional phrases“consisting of’ and“consisting essentially of’ shall be closed or semi-closed transitional phrases, respectively, as set forth in the United States Patent Office Manual of Patent Examining Procedures, Section 2111.03.