Scaling for quantum tunneling current in nano- and subnano-scale plasmonic junctions

When two conductors are separated by a sufficiently thin insulator, electrical current can flow between them by quantum tunneling. This paper presents a self-consistent model of tunneling current in a nano- and subnano-meter metal-insulator-metal plasmonic junction, by including the effects of space charge and exchange correlation potential. It is found that the J-V curve of the junction may be divided into three regimes: direct tunneling, field emission, and space-charge-limited regime. In general, the space charge inside the insulator reduces current transfer across the junction, whereas the exchange-correlation potential promotes current transfer. It is shown that these effects may modify the current density by orders of magnitude from the widely used Simmons’ formula, which is only accurate for a limited parameter space (insulator thickness > 1 nm and barrier height > 3 eV) in the direct tunneling regime. The proposed self-consistent model may provide a more accurate evaluation of the tunneling current in the other regimes. The effects of anode emission and material properties (i.e. work function of the electrodes, electron affinity and permittivity of the insulator) are examined in detail in various regimes. Our simple model and the general scaling for tunneling current may provide insights to new regimes of quantum plasmonics.

It is found that Simmons' formula is only accurate in a limited parameter space in the direct tunneling regime (Figs. [2][3][4][5], when the insulating thin film is relatively thick (> 1 nm) and the barrier height is relatively large (> 3 eV). Its accuracy decreases when the effective barrier height decreases (Fig. 4), or when the permittivity of the insulator decreases (Fig. 5), where the self-consistent model would provide a more accurate prediction of the tunneling current. In the field emission regime and space-charge-limited regime, the self-consistent model may be used, as Simmons' formula becomes fairly unreliable. The proposed model reveals the general scaling for quantum tunneling current and its dependence on the bias voltage, the dimension and material properties of the tunneling junction. It can be applied to broad areas involving tunneling junctions. As an example, its application in quantum plasmonics will be briefly addressed in the Discussion Section.
Note that although the present model is developed for DC condition, it is applicable to plasmonics of up to Near Infrared frequency. The underlying reason is that the transit time for electron tunneling through a barrier of nm-scale thickness is typically less than 1 fs 4,29-31 , which is much shorter than the period of the driving fields (e.g. 10 fs for 0.4 eV optical energy). This transit time is even shorter for insulator of sub-nm thickness. Thus, the electron would see an almost constant barrier during its transit time, and the DC model applies. Such treatments have been extensively applied in quantum plasmonic modeling [2][3][4] . The DC calculation would not be valid if the driving field frequency is so high (e.g. Visible light frequency or higher) that its period is comparable or less than the electron transit time.

Results
Self-consistent model for tunneling current. Consider two metallic electrodes separated by a thin insulating film, as shown in Fig. 1. Since the insulating film is assumed to be sufficiently thin (in the subnano-and nano-meter scale), charge trapping may be ignored 17,32 . The electrons in the electrodes would see a potential barrier formed between the two electrodes, where E F and W are the Fermi level and the work function of the metal electrodes respectively; X is electron affinity of the insulator; is the image charge potential energy including the effect of anode screening 4,14,33 , where e is the electron charge, є 0 is the permittivity of free space, є r is the relative permittivity of the insulator, and D is the gap distance; ( ) eV x is the electric potential, which is the sum of the potential due to the external applied voltage V g and the potential due to the electron space charge; and Φ ( ) x xc is the electron exchange-correlation potential, where the exchange potential is related to the Pauli Exclusion Principle, and the correlation potential denotes the quantum-mechanical part of the Coulomb interaction between electrons. The term  35 and the correlation energy 34 respectively, Figure 1. Metal-insulator-metal tunneling junction. The metal electrodes have Fermi level E F and work function W. The insulator thin film has electron affinity X, relative permittivity є r , and thickness D. The applied voltage bias is V g , the effective potential between the electrodes is Φ( ) x . The current densities emitted from the cathode and the anode into the gap are J 1 and J 2 , respectively. , and c, A, a 1 , b 1 , b 2 , b 3 , and b 4 are parametrized constants obtained using the random phase approximation 34 .
Following Simmons 14 , the probability ( ) D E x that an electron with longitudinal energy E x (normal to the surface) can penetrate the potential barrier of height Φ (x) is given by the WKBJ approximation 37 , where x 1 and x 2 are the two roots of Φ − ( ) = E x 0 x . The current density tunneling through the barrier from electrode 1 to the right is calculated by 4,14,28,38,39 , Figure 2. Current density as a function of applied gap voltage V g , for two gold (Au) electrodes ( = . W 5 1 eV) 39 separated by a vacuum gap ( = , = є X 1 0 r eV) of width = D 1 nm, at = T 300 K, (a) in normalized form in terms of CL law, J CL , (b) in unit of A/cm 2 . The calculations in (a) are from three methods: direction integration of eqs 3 and 5 (or 12 and 13), SCM with no V xc included, and full SCM with both space charge and V xc included. Simmons is for eq. 14, Simmons (V ~ 0) is for eq. 15.
is the total number of electrons inside electrode 2 with longitudinal energy between E x and + E dE x x impinging on the surface of electrode 1 across a unit area per unit time, calculated by free-electron theory of metal 40 .
Inside the insulator between the two electrodes, < < x D 0 , we use the mean-field theory [23][24][25] to solve the electric potential ( ) eV x and the exchange-correlation potential Φ ( ) x xc , as appeared in eq. (1). Thus, we solve the coupled Shrodinger equation and the Poisson equation [23][24][25] , where ψ is the complex electron wave function, ψψ = ⁎ n is the electron density, and E 0 is the electron emission energy (with respect to the Fermi energy E F ). Note that ψ is the superposition of two streams of electrons, one travelling from electrode 1 to electrode 2, and the other from electrode 2 to electrode 1 ( Fig. 1), both with emission energy of E 0 . We assume E 0 = 0 in the calculation 4,23,24,28 . For a given gap bias voltage V g , we have ( ) = V 0 0, and ( ) = V D V g . The boundary conditions on the wave function ψ are derived from the conditions that ψ and ψ/ d dx are continuous at = x 0, and = x D. Charge conservation requires that the current density J net = e(iħ/2m) (ψψ * ′ − ψ * ψ ′) = J 1 − J 2 be constant for all x, where a prime denotes a derivative with respect to x, and = − i 1 . It is convenient to introduce nondimensional quantities, , and E H is the Hartree energy. The wave function may be expressed in the normalized form to read 23 x are respectively the nondimensional amplitude and phase, both assumed real. Thus, the coupled Shrodinger equation and the Poisson equation, eqs 7 and 8, are expressed in their normalized form as, where eqs. 11c and 11d are derived by matching the wave function and its derivative at = x 1. The normalized emission current density J 1 and J 2 in eqs. (3) and (5) Note that the integrations in eqs. (12) and (13) are independent of the Fermi level E F . By solving eqs. (9)-(13) iteratively, we are able to self-consistently obtain the numerically converged results of the complete potential barrier profile of Φ (x) [eq (1)], the current density emitted from both electrodes J 1 , J 2 , and therefore the net current density J net , for any given materials of the electrodes (E F , W), thin film insulator (є r , X), film thickness (D), and external applied bias voltage (V g ). This is referred as the self-consistent method (SCM) thereafter.
Main results. Figure 2a shows the normalized current density γ (in terms of CL law) as a function of applied gap voltage V g , for two gold (Au) electrodes ( = .
W 5 1 eV) 39 separated by a = D 1 nm vacuum gap ( = , = є X 1 0 r eV). The current density in A/cm 2 is shown in Fig. 2b. The current densities are calculated from three methods: (1) direct integration using eqs 3 and 5, where space charge potential and exchange correlation potential V xc are not included in eq 1, (2) SCM without V xc , i.e. only space charge potential is included, and (3) complete SCM with both space charge potential and exchange correlation potential V xc included. As shown in Fig. 2, the − J V curves may be roughly divided into three regimes: direct tunneling regime ( In the direct tunneling regime, the tunneling current density from cathode γ 1 and that from anode γ 2 are comparable, where the latter was ignored in Refs. 4,25,28 . The net current density, which is the difference between γ 1 and γ 2 , γ γ γ = − net 1 2 , may therefore be orders of magnitude lower than both γ 1 and γ 2 . Thus, in the direct tunneling regime, both anode emission and cathode emission have to be considered to give an accurate evaluation in the tunneling current of the junction. The difference between γ 1 and γ 2 increases as V g increases. The three methods mentioned above give almost identical results for the current densities γ 1 , γ 2 , and γ net when < V 1 g V, which implies that the space charge potential and exchange-correlation potential are not important in the direct tunneling regime, for the given Au-vacuum-Au junction with 1 nm gap spacing. The − J V characteristic in the direct tunneling regime is linear, which indicates that the tunneling junction acts like an ohmic resister. The results are compared with the Simmons formula 14,19,21 for general V g , x T  (14) and (15) give a fairly good estimation in the − J V behavior of the given Au-Vacuum-Au structure when the applied bias < V 1 g V. It has been checked that the − J V curves for the Au-Vacuum-Au structure in Fig. 2 is very insensitive to temperature: only with an increase of < 2% from = T 0 K to 600 K, which is consistent with the relative small T dependence in eqs. 14 and 15. Physically, this is because the apparent barrier height of the Au-Vacuum-Au structure ( − W X = 5.1 eV) is much higher than the width change of the Fermi function (~0.5 eV), so that the majority electrons would still see an almost unchanged tunneling barrier for = T 0 K to 600 K. The temperature dependence would become important for junctions with small barrier heights (e.g. − W X ~1 eV). In the regime of 1 V < < V 10 g V, the tunneling current from anode γ 2 is much smaller compared to the cathode current γ 1 . By applying an appreciable bias voltage V g , the effective barrier height for the cathode is reduced, indicating an increase of current J 1 or γ 1 with V g . However, due to the down shift of the "effective" Fermi level (Fig. 1), the effective barrier height seen by electrons in the anode is increased, leading to a dramatic drop of current J 2 or γ 2 with V g . The tunneling behavior of the junction resembles field emission, thus we denote this regime the field emission regime. Field emission is most widely modeled by Fowler-Nordheim (FN) law 38,44,45  AeVV −2 and = . × B 6 8308 10 is the applied electric field. FN law is derived by assuming no anode screening. As shown in Figure 2, although the net current density γ net is approaching the FN law as V g increases, in general FN law is not sufficiently accurate to model the tunneling current in such a nano-scale junction 28 . In this regime, Simmons formula, eq. (14), gives a more accurate fit to the self-consistent SCM result. Note the breakdown of eq. (14) around Ṽ 10 g V, where the effective barrier height is depressed by V g below the Fermi level of the cathode. When V g is approaching 10 V, the current from direct integration (eqs. (3) and (5)) is closely fitted by Simmons formula, eq. (14). The current calculated from SCM by including only the space charge effect is slightly reduced. However, when both exchange-correlation V xc and space charge effects are included in the SCM, the resulting current is enhanced by one order of magnitude, indicating the profound effect of exchange-correlation energy in the field emission regime.
In the space-charge-limited regime of > V 10 g V, the direct integration method, which ignores both the space charge effect and the exchange-correlation effect, cannot provide a reliable estimate of the current. When only the space charge potential is included in the SCM calculation, the resulting current is reduced and is approaching classical Child-Langmuir (CL) law, J CL . However, when exchange-correlation potential is also included in the SCM calculation, the emitted current is enhanced in general. When V g reaches 100 V, the cathode current γ 1 (and therefore the net current γ net ) approaches the quantum CL law (QCL) 23,24 , which gives the maximum current density that can be transported across a vacuum nano-gap for a given V g and D, with quantum corrections. Figure 3a shows the net current density J net as a function of V g , for various gap width D for the Au-Vacuum-Au tunneling junction. Similar to Figure 2, the − J V curve may be roughly divided into three regimes: direct tunneling regime, field emission regime, and space-charge-limited regime. As gap width D decreases, the voltage range for both the direct tunneling regime and the space-charge-limited regime expands towards the field emission regime, whose voltage range decreases with D. In the direct tunneling regime, when > D 1 nm, the direct integration method and the SCM give almost identical results, where the Simmons formula (eq. (14)), which fits the direct integration well, is a very good approximation. However, when the gap width is in the sub nanometer range, < D 1 nm, the direct integration method (and therefore Simmons formula) underestimates the net current, thus the SCM including the effects of both space charge and exchange-correlation needs to be used to give more accurate calculation. In general, direct integration method would not be accurate in the field emission regime and space-charge-limited regime, where the SCM has to be applied. In the space-charge-limited regime, J net approaches QCL limit as V g increases. Figure 3b shows the net current density J net as a function of insulator thin film thickness D, for various V g for the Au-Vacuum-Au tunneling junction. It is important to see that the tunneling current, therefore the tunneling conductivity, is extremely sensitive to the thickness of the insulating thin film in MIM tunnel junctions. It is clear that for the limited parameter space, e.g. > D 1 nm and < V 10 g V, the direct integration calculation is accurate. Note that the values of gap voltage V g and gap spacing D in Fig. 3 are within the typical range of quantum plasmonic applications 3,4,39 .
The − J V characteristics of a MIM junction (Figure 1) is very sensitive to its apparent barrier height, − W X. Figure 4 shows J net as a function of V g for MIM junctions formed by various metal electrodes separated by a 1 nm wide vacuum gap. When the work function of the electrodes increases from = W 2 eV (Cs) to 5.1 eV (Au), J net in the direct tunneling regime ( < V 1 V g ) decreases by 6 orders of magnitude for a given bias. Simmons formula (eq. (14)) and the direct integration method are only accurate when − > W X 4 eV for a junction with vacuum gap = D 1 nm. When V g approaches 100 V, the current density J net converges to the same asymptotic value of QCL, since the space-charge-limited current density depends only on V g and D, but not on W. The effect of the electron affinity X of the insulating thin film (Fig. 1) on − J V would be similar, that is, increasing X would be equivalent to decreasing W, provided the relative permittivity є r of the insulator is unchanged.
It is interesting to note the nonmonotonic behavior of some curves in Fig. 3a (D = 0.5 nm when V g < 1 V) and Fig. 4 ). This is due to the profound effects of the nonlinear exchange-correlation potential, where the normalized insulating gap space , and the normalized gap voltage φ = / eV E g g H << 1 so that the space charge potential is not important compared to the exchange-correlation potential 24,25,28,36,46 .
The effect of relative permittivity є r of the insulating thin film is shown in Fig. 5. In the direct tunneling regime ( < V 1 g V), J net decreases as є r increases for a given V g . This is due to the fact that the image charge potential Φ image decreases as є r increases, as seen from the second line after eq. (1). Thus, the overall potential barrier will increase, leading to smaller tunneling current. In contrast, in the space-charge-limited regime ( > V 10 g V), J net calculated by SCM (solid lines) increases with є r , as clearly Scientific RepoRts | 5:09826 | DOi: 10.1038/srep09826 seen from Figure 5. This is because a larger є r reduces the effect of space charge, as seen from eq. (8) or (10), thus resulting in a larger SCL current. Note that J net calculated by direct integration (dashed lines) shows very different trends from that of SCM, indicating the dominant effects of space charge, which have to be included to give reliable predictions in the space-charge-limited regime. Thus, Simmons formula and the direct integration method are only accurate in the direct tunneling regime, when > є 2 r for junctions with 1 nm thickness and − = W X 3 eV. It is important to note that if є r is temperature dependent, the − J V behavior would also be temperature dependent 14 , even for tunneling junctions with relative big barrier height.

Discussion
Recently, the quantum-corrected model (QCM) 2,4 has been introduced to study charge transfer plasmon (CTP) 47,48 due to quantum tunneling, by accounting for the tunneling current across the gap via the insertion of an effective conductive medium in the gap. With the classical description, the permittivity ω ( ) є g of the effective medium is related to its DC conductivity σ g as ω σ ω where є 0 is the free space permittivity, and ω is the oscillating frequency. In the Drude model, the dielectric response of the effective conducting medium in the gap is characterized by where ω g is the plasmon frequency (typically set to the bulk plasma frequency of the surrounding resonators), and γ g is the tunneling damping parameter, which can thus be calculated as γ ω σ , under the assumption that γ g >>ω. The optical responses and the induced local fields of the quantum plasmon system are then obtained by standard classical approaches solving Maxwell's equations 2,47 . The validity of the calculation is crucially dependent on the two key parameters σ g and γ g , which describe the quantum tunneling resistance introduced by the presence of the gap.
As an example, in Fig. 6, σ g and γ g obtained from the proposed self-consistent model (SCM) are compared to those by direction integration (eqs. 3 and 5), and by the SCM but switching off the emission from anode (similar to Refs. 4,28 ), for a tunneling junction with = D 1 nm vacuum gap, and electrode work function = . W 2 9 eV. For simplicity, we estimate the DC quantum gap conductivity as σ = / J E g n et is the applied electric field across the tunneling gap. In direct tunneling regime ( < E 10 9 V/m) and the field emission regime ( < < × E 10 3 10 9 9 V/m), direct integration method (or Simmons formula) underestimates the gap conductivity and overestimates the tunneling damping. In the space-charge-limited regime ( > × E 3 10 9 V/m), direct integration method is generally not reliable. Ignoring the current emission form anode (i.e. set = J 0 2 in eq 5) would result in a much higher σ g and much lower γ g in the direct tunneling regime. The relative large damping γ g calculated from SCM in the direct tunneling regime suggests that CTP via quantum tunneling in this regime would be very difficult to observe experimentally. Instead, for a given junction, by simply increasing the driving field E to reach the field emission or space-charge-limited regime, the damping γ g can be significantly reduced so that the experimental realization of CTP via tunneling could be relatively easier. Figure 6. Quantum tunneling gap DC conductivity σ g and the tunneling damping parameter γ g for a MIM plasmonic tunneling junction, as a function of applied electric field E, with a vacuum gap of D = 1 nm, X = 0, and work function of the electrodes W = 2.9 eV at 300 K. The plasmon frequency is assumed to be ω π = × ( × ) 2 2 10 g 15 rad/s for the calculation of γ g . The solid lines are for SCM, dashed lines for direct integration of eqs 3 and 5, and dotted lines for SCM with anode emission being switched off (i.e. set J 2 = 0 in eq 5 and γ 2 = 0 in eq 13).
Scientific RepoRts | 5:09826 | DOi: 10.1038/srep09826 In summary, we have developed a self-consistent model to characterize the tunneling current of nanoand subnano-scale plasmonic junctions, by taking into account of the effects of both space charge and exchange-correlation potential. The effects of material properties, including the work function of the electrodes W, the permittivity є r and the electron affinity X of the insulator, are examined in detail. In general, the − J V curves may be divided into three regimes: direct tunneling regime, field emission regime, and space-charge-limited regime. It is found that Simmons formula (eqs. (14) and (15)) are good approximations of the tunneling current for a limited parameter space in the direct tunneling regime only. Their accuracy decreases when the effective barrier height decreases, i.e. W decreases or X increases, or when the permittivity of the insulator є r decreases. They become unreliable when the insulator thickness is in the sub-nanometer scale, < D 1 nm, where the self-consistent model would give a more accurate evaluation.
In this formulation, we have made the following widely used assumptions: 1) the electron transmission probability during the emission process is approximated by the WKBJ solution, where the metal electrodes are based on the free electron gas model; 2) the surfaces of the electrodes are flat and the problem is assumed one-dimensional; 3) the image potential is approximated by the classical image charge methods. The effects of electrodes geometry, nature of the ion lattice of the electrodes, possible charge trapping inside the insulator film, frequency dependence, and dissimilar electrodes will be subjects of future studies.