Universal Model of Optical-Field Electron Tunneling from Two-Dimensional Materials

Introduction.—Laser-matter interaction offers the capability for the manipulation of electron excitation and dynamics at ultrashort timescale, such as high-order harmonic generation [1, 2], carrier interband transition [3], spontaneous radiation [4], photoelectron emission [5, 6, 7, 8, 9, 10, 11, 12], quantum coherent control of excitation states [13, 14, 15, 16] and many others. Among them, laser-triggered photoemission from solids has gained considerable current attention, because of its crucial role in the development of high-resolution electron microscopy and diffraction [17, 18], free electron lasers [19], tabletop laser accelerators [20], coherent electron sources [21, 22] and quantum nano-vacuum electronics [23, 24, 25, 26].

At high laser intensity, the photoemission is due to the optical-field electron tunneling process (or optical-field emission) [27, 28], where the strong time-varying light field greatly suppresses the material interface potential barrier enabling electron tunneling in the vicinity of Fermi level (cf. Fig. 1). In this non-perturbative field-driven near-transient tunneling regime, incident optical electric field can achieve the steering of electron motion at a sub-cycle temporal scale, such as the subfemtosecond tunneling [29, 30] and attosecond electron pulse generation [31], which makes optical-field tunneling emission is important for ultrafast optical-field-driven electronics [32] and on-chip attosecond lightwave science [33, 34].

Other than metallic nanostrucutres [8, 10, 27, 35, 36, 37, 38], emerging two-dimensional (2D) materials [39, 40, 41, 42, 43] are employed to study their optical-field photoemission as well, due to the distinctive electronic band structures and nanoscale sharpness which provide the huge localized field enhancement factor allowing sufficient electron emission in using low-power lasers. For instance, the most common 2D material graphene has the atomic-scale thickness resulting in significantly high field enhancement near the edge [44], large electron mobility [45] and high thermal damage threshold [46], which explains graphene being used effectively to realize ultrafast optical-field electron tunneling emission [40]. Despite tremendous ongoing efforts to implement 2D materials as emitters for ultrafast applications, the theoretical description of 2D material electron emission relies heavily on the Fowler-Nordheim (FN) formalism (which was developed for the traditional bulk materials) that is fundamentally incompatible with the reduced dimensionality and electronic properties of 2D materials [47, 48, 49, 50, 51, 52, 53, 54]. The optical-field tunneling model of 2D materials thus remains an urgent and open question that has yet to be addressed.

In this work, we develop universal models for optical-field tunneling emission from the edge and surface of a wide class of 2D materials, which explicitly considers the reduction of dimensionality, non-parabolic energy dispersion of 2D material, quasi-time dependent tunneling and electron emission dynamics driven by ultrafast laser. We reveal a universal current-field scaling: In($J/|F|^{\beta})\propto 1/|F|$ with $\beta=3/2$ for the edge emission and $\beta$ = 1 for the vertical surface emission respectively, which are different from the conventional FN law of $\beta=2$ [55]. A peculiar saturation of photocurrent at ultrahigh laser intensity is identified due to the reduced dimensionality. Importantly, our model is in an excellent agreement to a recent experiment [40] of optical-field emission from monolayer graphene over different combinations of dc and laser fields. Our model provides a universal description of optical-field tunnelling phenomena for a wide variety of 2D materials and shall form an important foundation for the development of novel 2D-material based ultrafast optical-field-driven tunneling nano-vacuum optoelectronics and device engineering.

Theoretical formalism.—We consider a 2D material lying in the $x-y$ plane with a dc field $F_{0}$ under the illumination of a laser field $F_{1}(t)$ perpendicular to the edge. The optical-field-induced tunneling current density via the edge along $x$ direction is [50, 49]

$$J_{e}(t)=\frac{eg}{(2\pi)^{2}}\sum_{k_{\perp}^{(i)}}\int[v_{x}(\varepsilon_{\parallel})f_{FD}(\varepsilon_{\parallel})\Gamma(\varepsilon_{x},t)]d^{2}\emph{{k}}_{\parallel},$$

(1)

where $e$ is the electron charge, $g$ is the spin-valley degeneracy factor, $\emph{{k}}_{\parallel}=(k_{x},k_{y})$ and $\varepsilon_{\parallel}$ are the electron wavevector and energy component in the $x-y$ plane, respectively. Here, $v_{x}(\varepsilon_{\parallel})=\hbar^{-1}\partial\varepsilon_{\parallel}/\partial k_{x}$ is the electron $x$-direction velocity ($\hbar$ is the reduced Planck constant), $f_{FD}(\varepsilon_{\parallel})$ is the zero-temperature Fermi-Dirac distribution function, $\Gamma(\varepsilon_{x},t)$ is the time-dependent electron transmission probability through the edge at energy component $\varepsilon_{x}$ [cf. Eq. (3)], $k_{\perp}^{(i)}\;(i=1,2,3,…)$ denotes the discrete bound states along the direction orthogonal to the $x-y$ plane, due to the confinement of electrons within the 2D plane, and $\sum_{k_{\perp}^{(i)}}$ denotes the summation of all these discrete bound states.

Consider a 2D material with a general form of anisotroppic energy dispersion,

$$\varepsilon_{\parallel}(k_{x},k_{y})=\left(\alpha_{l}k_{x}^{2}+\beta_{l}k_{y}^{2}\right)^{l/2},$$

(2)

where $\alpha_{l}$ and $\beta_{l}$ are the material-dependent parameters along $x$ and $y$ directions respectively, and $l\in\mathbb{Z}^{>0}$. Using Eq. (2), we derive $v_{x}(\varepsilon_{\parallel})=\hbar^{-1}l\alpha_{l}k_{x}(\alpha_{l}k_{x}^{2}+\beta_{l}k_{y}^{2})^{(l-2)/2}$, and $f_{FD}(\varepsilon_{\parallel})=H(\varepsilon_{F}-\varepsilon_{\parallel})=H(\varepsilon_{F}^{2/l}-\alpha_{l}k_{x}^{2}-\beta_{l}k_{y}^{2})$, where $H(x)$ denotes the Heaviside function.

In the optical high-field tunneling regime, incident optical field primarily behaves like an ac electric field that modulates the barrier, in which the effects of photoexcited nonthermal electrons can be ignored [40, 7]. Thus, the optical-field tunneling through the barrier resembles instantaneous field emission as a function of laser field, and the time-dependent transmission probability is approximated at quasi-static limit based on the Wentzel-Kramers-Brillouin (WKB) method, which gives

$$\Gamma(\varepsilon_{x},t)=H\left[-F(t)\right]\times D_{F}(t)\times\exp\left(\frac{\varepsilon_{x}-\varepsilon_{F}}{d_{F}(t)}\right).$$

(3)

Here, $F(t)$ is the summation of the dc field $F_{0}$ and the time-dependent optical field $F_{1}(t)$, and $H[-F(t)]$ denotes the Heaviside function to ensure electron emission occurs only during the negative field cycles. We omit $H[-F(t)]$ in the expressions below for simplicity. The tunneling factor is $D_{F}(t)=\exp\left[-2W/3d_{F}(t)\right]$, where $d_{F}(t)=e\hbar|F(t)|/2\sqrt{2mW}$, $W$ is the work function, $\varepsilon_{F}$ is the Fermi energy, and $m$ is the free-electron mass. Substituting Eq. (3) into Eq. (1) yields the optical-field edge tunneling current density, $J_{e}=A_{l}\times\alpha_{l}^{l/2}\times K(l)$, where $A_{l}=egl\varepsilon_{F}^{1/l}D_{F}/2\pi^{2}\hbar\sqrt{\beta_{l}}$, and $K(l)$ is a complex integral as shown in Eq. (S4) in the Supplemental Materials (SM). In optical-field tunneling, the emission electrons are mostly located near the Fermi level: $\alpha_{l}^{l/2}k_{x}^{l}\to\varepsilon_{F}$, thus $\sqrt{\varepsilon_{F}^{2/l}/\alpha_{l}k_{x}^{2}-1}\ll 1$. Under this approximation, the tunneling current $J_{e}$ can be analytically solved as (see Sec. I in the SM for details)

$$J_{e}\approx\begin{cases}\vspace{5mm}B_{1}d_{F}[I_{1}(\frac{\varepsilon_{F}}{d_{F}})+L_{1}(\frac{\varepsilon_{F}}{d_{F}})]\text{exp}(-\frac{\varepsilon_{F}}{d_{F}}),&\text{$l=1$}\\ B_{2}[C_{1}\hskip 2.27626pt{}_{l-1}F_{l}(\zeta_{1},\xi_{1},\frac{\varepsilon_{F}^{2}}{4d_{F}^{2}})+C_{2}\frac{\varepsilon_{F}}{d_{F}}\hskip 1.42271pt_{l}F_{l+1}(\zeta_{2},\xi_{2},\frac{\varepsilon_{F}^{2}}{4d_{F}^{2}})]\text{exp}(-\frac{\varepsilon_{F}}{d_{F}})+\\ \vspace{5mm}B_{3}[C_{3}\frac{\varepsilon_{F}}{d_{F}}\hskip 1.42271pt_{l-1}F_{l}(\zeta_{3},\xi_{3},\frac{\varepsilon_{F}^{2}}{4d_{F}^{2}})+C_{4}\hskip 2.27626pt{}_{l}F_{l+1}(\zeta_{4},\xi_{4},\frac{\varepsilon_{F}^{2}}{4d_{F}^{2}})]\text{exp}(-\frac{\varepsilon_{F}}{d_{F}}),&\text{$l$ is odd (= 3, 5, 7, ...)}\\ B_{4}C_{5}\hskip 2.27626pt{}_{l/2}F_{l/2}(\zeta_{5},\xi_{5},\frac{\varepsilon_{F}}{d_{F}})\text{exp}(-\frac{\varepsilon_{F}}{d_{F}}),&\text{$l$ is even (= 2, 4, 6, ...)}\par\end{cases}$$

(4)

where $B_{1}=A_{1}\pi/2$, $A_{1}$ is $A_{l}$ with $l=1$, $I_{1}$ is the modified Bessel function of the first kind, $L_{1}$ is the modified Struve function, $B_{2}=A_{l}\varepsilon_{F}\frac{\Gamma(l/2)}{2\Gamma(3/2)\Gamma[(l+1)/2]}$, $B_{3}=A_{l}\varepsilon_{F}\{1-\frac{\Gamma(l/2)}{2\Gamma(3/2)\Gamma[(l+1)/2]}\}$, $\Gamma(x)$ is the complete gamma function, $B_{4}=A_{l}\varepsilon_{F}$, $C_{1},C_{2},C_{3},C_{4}$ and $C_{5}$ (including $C_{6}$, $C_{7}$ and $C_{8}$ below) all are $l$-dependent constants, and ${}_{p}F_{q}$ is the generalized hypergeometric function. The details of the constants and other parameters can be found in Sec. I in the SM.

The universal analytical $J_{e}$ [Eq. (4)] is characterized for three classes of the 2D materials [$l$ = 1, odd $l$ $(>1)$ and even $l$], which shows good agreements with the numerical results from Eq. (S4) (cf. Fig. S1). Equation (4) also indicates that $J_{e}$ is independent of the material parameter $\alpha_{l}$, and thus a possible universal edge emission scaling. At the practical field regime ($d_{F}\ll\varepsilon_{F}$) and the ultrahigh field regimes ($d_{F}\gg\varepsilon_{F}$), the analytical results for the universal model and three representative 2D systems are summarized in Table 1 (see Sec. II in the SM for detailed derivations). The validity of these analytical equations are verified in Fig. S2 in the SM.

Intriguingly, at practical intensity regime (column 2 of Table 1), the prefactor ($d_{F}\propto|F|$) of the field strength exhibits an unconventional universal scaling of $J_{e}\propto|F|^{\beta}$ with $\beta=3/2$ for $all$ $l$, which is in stark contrast to $\beta=2$ of the conventinoal FN law [55] derived for traditional 3D bulk materials. Interestingly, the same universal $\beta=3/2$ scaling has also been reported in the current-temperature scaling of the thermionic emission for lateral 2D-material-based Schottky contacts [49]. At ultrahigh field regime (column 3 of Table 1), $J_{e}$ saturates and becomes independent of field strength $|F|$ (or $d_{F}$). Such saturation reveals a dramatic consequence of the reduced dimensional of 2D material, in which the emission current saturates due to limited availability of electrons. Such a source-limited emission behavior [56] is completely different from the space-charge limitation of traditional bulk materials [24, 25], and shall offer a distinctive high-field transport signature of 2D-material-based electron optical-field emitters.

Optical-field tunenling from 2D material surface.—We now consider the case of optical-field tunneling of electrons that occurs vertically from the plane of the 2D materials. In this case, the dc field $F_{0}$ and high-intensity optical field $F_{1}(t)$ are considered to be perpendicular to the 2D surface. During this out-of-plane tunneling process, carrier scattering effects, such as electron-electron [57] and defects scatterings [58], violates the conservation of in-plane momentum $\emph{{k}}_{\parallel}$, which leads to the coupling between $\emph{{k}}_{\parallel}$ and $k_{\perp}^{(i)}$. Thus, the vertical optical-field tunneling current density from 2D material surface is expressed as [49],

$$J_{s}(t)=\frac{eg}{(2\pi)^{2}}\sum_{k_{\perp}^{(i)}}\frac{v_{\perp}^{(i)}(k_{\perp}^{(i)})}{L_{\perp}}\int[f_{FD}(\varepsilon_{\parallel})\Gamma(\varepsilon_{\parallel},\varepsilon_{\perp}^{(i)},t)]d^{2}\emph{{k}}_{\parallel},$$

(5)

where $L_{\perp}$ is the thickness of 2D material and $v_{\perp}^{(i)}(k_{\perp}^{(i)})$ is the cross-plane group velocity of the electron with discrete energy state $\varepsilon_{\perp}^{(i)}(k_{\perp}^{(i)})$. Consider the surface tunneling from one subband near the Fermi level, Eq. (5) becomes

$$J_{s}(t)=\frac{ev_{\perp}}{L_{\perp}}D_{F}(t)\int_{0}^{\varepsilon_{F}}\{D(\varepsilon_{\parallel})\text{exp}[(\varepsilon_{\parallel}-\varepsilon_{F})/d_{F}(t)]\}d{\varepsilon_{\parallel}},$$

(6)

where $D(\varepsilon_{\parallel})$ is the density of state (DOS), and $D_{F}(t)$ and $d_{F}(t)$ are defined in Eq. (3).

Employing a general 2D anisotropic DOS $D(\varepsilon_{\parallel})=g\times C_{k,q}\times\varepsilon_{\parallel}^{k}$ with $C_{k,q}$ being material-dependent coefficients and $k,q\in\mathbb{Z}^{\geq 0}$ [54], Eq. (6) yields a generalized analytical form

$$J_{s}(t)=\frac{egv_{\perp}}{L_{\perp}}D_{F}C_{k,q}(-d_{F})^{k+1}U\text{exp}(-\frac{\varepsilon_{F}}{d_{F}}),$$

(7)

where $U=\Gamma(1+k)-\Gamma(1+k,-\varepsilon_{F}/{d_{F}})$, and $\Gamma(s)$ and $\Gamma(s,x)$ are the complete gamma function and the upper incomplete gamma function respectively. By taking the Taylor series expansion of $\Gamma(1+k,-\varepsilon_{F}/{d_{F}})$ for $d_{F}\ll\varepsilon_{F}$ and $d_{F}\gg\varepsilon_{F}$, $J_{s}$ is solved under these two limits for the universal formula and three different materials in Table 2 (see Sec. III of the SM for more details).

The vertical surface optical-field emission current all converge to an universal scaling of In($J_{s}/|F|^{\beta})\propto 1/|F|$ with $\beta=1$ at the practical field regime, which is different from $\beta=3/2$ of edge emission presented in Table 1 and $\beta=2$ of the classical FN law. The current also saturates at the ultrahigh intensity due to the reduced dimensionality. The validity of these analytical solutions is shown in Fig. S3 in the SM. By comparing Figs. S2 and S3, it is found under the same illumination condition, vertical surface optical-field tunneling generates significantly higher emission current than the edge tunneling. Notably the analytical formulas in this work can be directly employed to study dc field emission by replacing $F(t)$ with dc field $F_{0}$.

Current-field characteristics.—We verify our model in comparison with a recently measured optical-field tunneling emission current from the edge of monolayer graphene [40] as shown in Fig. 2(a). Due to the significantly low backscattering efficiency on the carbon surface [59], the contributions from rescattered electrons are excluded in the calculations of Fig. 2 (see Sec. IV of the SM for the method). Our model produces a good agreement with the experiment by using the dc and laser field enhancement factors of about 146 and 22 respectively, which are smaller than the reported values [40] based on the traditional bulk material FN law.

The optical-field edge emission current density $J$ from the graphene and its power dependence on the laser field are calculated as a function of laser field $F_{1}$ and dc field $F_{0}$ in Fig. 2(b) and 2(c), respectively. A Gaussian laser pulse is chosen, $F_{1}(t)=F_{1}e^{-4\textnormal{In}2\,t^{2}/\tau^{2}}\cos(\omega t)$, where $F_{1}$ is the peak of field strength, $\tau$ is the pulse width, and $\omega$ is the angular frequency. The laser wavelength is 800 nm with a pulse width of $\tau$ = 8 fs (3 cycles) and the integration time [cf. Eq. (S40)] is from -30 to 30 fs. Increasing either laser intensity $F_{1}$ or dc bias $F_{0}$ will increase the optical-field emission [cf. Fig. 2(b)]. Interestingly, increasing $F_{0}$ changes the order of nonlinear power dependence of $J$ on $F_{1}$ from the decreasing to increasing trends [cf. Fig. 2(c) and dotted lines in Fig.2(d)]. This transition occurs at around $F_{0}=10$ V/nm and exists in the case without the consideration of emitted electron dynamics in the vacuum. Such behavior arises from the having higher electron tunneling probability through the barrier with stronger dc bias, thus allowing relatively more electrons to contribute in the high-intensity regime. When the relative number of emitted electrons is sufficiently large, the unconventional increasing exponential growth rate of current with laser intensity takes place [cf. the case with $F_{0}>10$ V/nm as shown in Fig. 2(c)]. This nonlinear order increment is also in agreement with experimental observations [40].

Sub-cycle electron emission characteristics.—The time-dependent optical-field emission current $J_{e,\text{Gr}}(t)\delta(t)$ [cf. Eq. (S40)] in Fig. 2(e) shows that most electrons are emitted during the negative first half cycles due to the attainment of sufficient kinetic energy to escape the confinement of positive optical cycle. Beacsue of the reduced kinetic energy, parts of the electrons emitted in the negative second half cycles can be driven back to the edge by the deceleration of following positive cycle. This confirms the sub-optical-cycle emission dynamics in the optical-field emission regime. Increasing the dc field enlarges the rescattering time (cf. the sharp edge of emission yield), which implies that under the stronger dc acceleration, more electrons could escape the laser field and contribute to the net photocurrent. Finally, we calculate the kinetic energy spectra at the position of around 25 nm away from the graphene edge for different laser field $F_{1}$ at zero dc field ($F_{0}$ = 0). The extracted main spectrum peaks are plotted in Fig. 2(f). The cutoff kinetic energy representing the maximum kinetic energy of collected photoelectron as a function of laser intensity is depicted in the inset of Fig. 2(f). The sublinear dependence of the cutoff energy on the laser intensity indicates a sub-cycle-timescale electron emission characteristics [60] and the electron dynamic motion can be controlled by the optical field strength.

Conclusion.—In summary, we have presented analytical and universal models for ultrafast optical-field tunneling emission of electrons from the edge and surface of general anisotropic 2D materials. Our model reports a universal scaling between the optical tunneling current density $J$ and laser electric field strength $F$ given by $\textnormal{In}(J/|F|^{\beta})\propto 1/|F|$ with $\beta$ = 3/2 in the edge emission and $\beta$ = 1 in the vertical surface emission, which both are different from the traditional scaling of $\beta=2$ based on the FN law. At ultrahigh laser field, we observe a saturation effect that the tunneling current no longer depends on the field, due to the reduction of dimensionality. These universal scaling and saturation are explicitly exhibited in three selected materials: monolayer grapehene, highly doped multilayer black phosphorus and $ABC$-stacked trilayer grapehene. Our calculation also well reproduces prior experimental measurement, and uncovers the dc induced modulation to the optical-field emission current and the power-dependence of photocurrent on laser field. The developed analytical formulation provides a theoretical foundation for understanding a broad range of 2D materials in applications related to ultrafast optical-field-driven optoelectronics and field emission in vaccum nanoelectronics. The model shall pave a way towards a generalized model combining thermal, field and photoemission [61] for 2D materials which is highly desirable for the design of 2D-material-based electron emitter.

This work is supported by the Singapore A*STAR IRG grant (A2083c0057). Y.S.A. is supported by the Singapore Ministry of Education (MOE) Academic Research Fund (AcRF) Tier 2 under Grant No. MOE-T2EP50221-0019.