††thanks: Corresponding author. Arun V. Kulkarni Tel.:0832 2580309; fax:
+91-8322557033 E-mail address<EMAIL_ADDRESS>
Also at ] Department of Physics, Parvatibai Chowgule College of Arts and
Science, Margao, Goa, 403602.

# Calculation of tunneling current across Trapezoidal potential barrier in a
Scanning Tunneling Microscope

Malati Dessai  [  Arun V Kulkarni<EMAIL_ADDRESS>Department of
Physics, Birla Institute of Technology and Science-Pilani,K K Birla Goa
campus, Goa, 403726.

###### Abstract

The Planar Model of the Electrode-Vacuum-Electrode configuration for STM in
which electrode surfaces are assumed to be infinite parallel planes, with
atomic size separation and vacuum between them, is used to calculate tunneling
current densities for both low and high bias voltages. Non WKB, Airy function
solutions for the Schrödinger Equation for the trapezoidal barrier in the
tunneling region are used to calculate the tunneling probability. Temperature
dependent Fermi Factors for each electrode are introduced and the calculation
involves integration over the electron energies. In order to convert the
current densities obtained in the planar model to tunneling currents the tip
and sample surfaces are modelled as confocal hyperboloids, and the tip sample
distance is replaced by the length of the line of force (field line). The
current is found by integrating the current density over a finite area of the
tip. The calculated tunnel currents for a few electrode pairs at room
temperature are plotted for several values of bias voltage, tip sample
distances, and tip radii of curvature. Pauli Effects are studied as a function
of bias voltage and tip-sample distance. Some estimate of lateral resolution
and its dependence on bias voltage and tip radius is also presented.

Keywords: Tunneling; Tunnel junctions; Trapezoidal barrier; confocal
hyperboloids; Electrode-vacuum-Electrode.

##  I Introduction

Scanning Tunneling Microscopes (STM’s) allow imaging of surface structures
with atomic scale resolution and have been the workhorses for most surface
science studies in recent decades. Scanning probe microscopes usually involve
a sharp tip brought close to a sample which is mostly flat except for atomic
scale surface features. The atoms of the tip interact with those of the sample
and in case of the STM this interaction causes tunelling of the electrons
between the tip and the sample. An experimentally measurable quantity such as
the tunneling current is measured as a function of the tip sample distance and
the bias potential between tip and sample. In typical STM applications,
extremely sharp tips are used and further a very close proximity between the
tip and sample is maintained. As a result, the tip-sample interaction is
highly localized and the tunneling current is sensitive only to the local
surface properties of the sample. This sensitivity is further enhanced due to
the STM current having an almost exponential dependence on the tip-sample
distance, thus achieving atomic scale resolutionbinnig5 . Since the tunneling
is a very sensitive function of the tip sample distance, the measured values
of the tunnel current almost resemble the sample surface topography when
raster scanning of the sample is carried out by the tip.

There are several other works in literature, that calculate tunneling current
density such as SimmonssimmonsI , Hartman and ChivianHartmanChi  in which the
WKB approximation is used for calculating the tunneling probability. While
tunneling current densities are not directly experimentally measurable
tunneling currents are directly measurable and these are of particular
interest in comparing theory with experiment.

Tersoff and HamannTersoffHamann  have developed a formalism starting from
Bardeen’s theoryBardeenI BardeenII  to calculate tunneling currents. However,
their calculations is likely to work best for very low tip-sample distances
and also for very low biases. In this paper the tip-sample distances exceed
this limit and range from a few $\AA$ to 20 $\AA$. The calculation of Chen
et.al.chen  uses Bessel functions of one-third order and therefore goes beyond
the WKB approximation. However only the current densities are reported and
effects of Pauli Blocking are not investigated. The self consistent
calculation of Banerjee et.al.SBPZ  uses the WKB approximation. They have not
studied Pauli Blocking effects and also report only current densities and not
currents.

The tip and sample here are initially treated as plane surfaces. In this so
called planar model, only tunneling current densities and not tunneling
currents can be calculated. The latter would be infinity in a planar model. In
the non planar model considered here, the tip and the sample are modelled as
confocal surfaces in a prolate spheroidal coordinate system. In this paper the
tunneling current densities are converted into currents through a method
described by Saenz and GarciaSG .

The energies of electrons inside the electrodes are distributed according to
the Fermi-Dirac distribution. With these assumptions the tunneling current
density and tunneling currents through the barrier are calculated for both
high and low bias voltages.

##  II THEORY

###  II.1 A. CURRENT DENSITY AND FERMI FACTOR

Consider two plane conducting electrodes with infinite plane surfaces,
separated by vacuum whose width ranges from 2 Å to 20 Å. The electrode
surfaces are initially assumed to be of infinite area. With this assumption,
both the Laplace equation for the electrostatic potential, and the Schrödinger
equation with this potential, reduce to one dimensional equations.

(a) zero bias

(b) Bias = 3.0 eV $(eV_{b}<\phi_{1})$

Figure 1: Dissimilar Electrodes with vacuum barrier at a) Zero Bias b) Bias =
3.0 eV

Let a bias potential $V_{b}$ be applied across the electrodes. Let $\eta_{1}$,
$\phi_{1}$, $\eta_{2}$, $\phi_{2}$, be the Fermi energies and work functions
of the left (first) electrode and the right (second) electrode respectively
(Fig.1(a)). The left electrode (electrode 1) is assumed to be at higher
potential than the right electrode (electrode 2) so that the electric field
points along the negative x axis, and if a free electron were introduced
between electrodes it would experience a force along the positive x axis as
shown in Fig.2. Tunneling of electrons between the two plane surfaces occurs
through the effective trapezoidal barrier in the vacuum region. The case for
zero bias between dissimilar electrodes is shown in Fig.1(a). The barrier
potential in this case is the difference $(\phi_{1}-\phi_{2})$ between the
work functions of the two conducting electrodes. When biased, the barrier
potential in the vacuum region, is the sum of the contact potential
$(\phi_{1}-\phi_{2})$ and the applied external bias $eV_{b}$ as shown in
Fig.1(b). Forward electron current density is along the positive x axis and
reverse electron current density is along the negative x axis.

The connection between the tunneling current density and the probability of
tuneling through the barrier, is well established and the earliest accounts of
it are by Fowler and NordheimFN  and SimmonssimmonsI simmonsII . In these
accounts, one of the principal assumptions is that the number of electrons
that can occupy a unit volume in the 6 dimensional phase space is
$\dfrac{2}{h^{3}}$, and therefore the number per unit volume (of coordinate
space) of electrons in the velocity space volume element
$d^{3}v=dv_{x}\,dv_{y}\,dv_{z}$ about the velocity vector ${\bf v}$ is
$\mathfrak{n}(v)d^{3}v=\frac{2m}{h^{3}}d^{3}v$ The energy of these electrons
is

$\displaystyle E=\frac{1}{2}m{\bf v}^{2}$ $\displaystyle=$
$\displaystyle\frac{1}{2}m(v_{x}^{2}+v_{y}^{2}+v_{z}^{2})$ (1)
$\displaystyle=$ $\displaystyle E_{x}+E_{y}+E_{z}$

where
$E_{x}=\frac{1}{2}mv_{x}^{2},E_{y}=\frac{1}{2}mv_{y}^{2},E_{z}=\frac{1}{2}mv_{z}^{2}$
are the kinetic energies of the electrons for translational motion along the
$x,\,y,\,\text{and}\,z$ directions respectively. The number density of
electrons in the electrode 1 that are available for tunneling would be
$n_{1}\,d^{3}v=\mathfrak{n}(v)f_{1}(E)\,d^{3}v$ where
$f_{1}(E)=[1+\exp({\beta(E-\eta_{1})})]^{-1}$ is the Fermi Dirac Distribution
function which can be regarded as the probability of occupation of the energy
level $E$ in the electrode 1. Note $\beta=\frac{1}{kT}$ where $k$ is the
Boltzmann Constant and $T$ is the temperature (in Kelvin) of both the
electrodes.

Figure 2: Electrode 1 is at higher potential ($U_{0}=\eta_{1}+\phi_{1}$) and
the elctric field points along the negative x axis

The number density of electrons that can tunnel through the barrier is

$n_{T}(v)d^{3}v=n_{1}(v)D(E_{x})\,d^{3}v$

where $D(E_{x})$ is the probability of tunneling through the barrier. In
theoretical calculation

$D(E_{x})=\rvert T(E_{x})\rvert^{2}\frac{j_{out}}{j_{inc}}$

where $T$ is the transmission amplitude for tunneling through the barrier and
if plane wave solutions are used within the electrodes then
${j_{inc}}=\frac{\hbar k_{1}}{m}$ and ${j_{out}}=\frac{\hbar k_{2}}{m}$ for
forward currents. where $\hbar k_{1}$ and $\hbar k_{2}$ are x-components of
the momenta of the electrons in the first and the second electrodes.

$\hbar k_{1}=\sqrt{2mE_{x}},\,\,\hbar k_{2}=\sqrt{2m(E_{x}-\Delta_{B})}$ (2)

where $E_{x}$ is the kinetic energy of motion of the electrons in the $x$\-
direction. For reverse currents ${j_{inc}}=\frac{\hbar k_{2}}{m}$ and
${j_{out}}=\frac{\hbar k_{1}}{m}$. Note that the transmission amplitude is
independent of the direction of tunneling. Hence it is the same for both
forward and reverse tunneling currents. Not all of the electrons emitted by
electrode 1 will necessarily be accepted by the electrode 2, since there will
be partial occupation of this level by electrons already present in the
electrode 2. This is due to the Pauli Exclusion Principle and may be called
the Pauli Effect. OtherssimmonsI Shu7 riza Ebeling  have referred to this
phenomena as Pauli Blocking.

To obtain the number density of electrons in the energy level $E$, that are
accepted by the electrode 2, one must multiply $n_{T}$ by the probability of
vacancy of this energy level, in the electrode 2. This is given by
$[1-f_{2}(E+\Delta_{b})]$ where

$f_{2}(E+\Delta_{b})=[1+\exp({\beta(E+\Delta_{b}-\eta_{2})})]^{-1}$

The energy argument in the Fermi Dirac function for electrode 2 will be
shifted by $+\Delta_{b}=\eta_{1}-\eta_{2}-eV_{b}$, because the zero level of
the energy in electrode 2 is shifted by $-\Delta_{b}$ compared to that of the
electrode 1. But as can be seen from Figures 2 and 6

$(E+\Delta_{b}-\eta_{2})=(E+eV_{b}-\eta_{1})$

Therefore

$f_{2}(E+\Delta_{b})=f_{1}(E+eV_{b})$

Thus the number density of electrons that will be found in the states of
electrode 2, and which contribute to the forward current are

$\displaystyle n_{\text{For}}d^{3}v$
$\displaystyle=\dfrac{2m^{3}}{h^{3}}\frac{k_{2}}{k_{1}}\lvert
T(E_{x})\rvert^{2}f_{1}(E)[1-f_{2}(E+\Delta_{b})]d^{3}v$ (3)
$\displaystyle=\dfrac{2m^{3}}{h^{3}}\frac{k_{2}}{k_{1}}\lvert
T(E_{x})\rvert^{2}f_{1}(E)[1-f_{1}(E^{\prime})]d^{3}v$ (4)

where $E^{\prime}=E+eV_{b}=E_{x}+E_{r}+eV_{b}$ The forward current density is
then

$J_{\text{For}}=e\int v_{x}\,n_{\text{For}}\,d^{3}v$

Now $v_{x}dv_{x}=(dE_{x}/m)$ and

$dv_{y}dv_{z}=v_{r}dv_{r}d\theta=(dE_{r}/m)d\theta$

where the $y,z$ velocity components are expressed in terms of their polar
components $(v_{r},\theta)$ in the $y-z$ plane. Also
$E_{r}=\frac{1}{2}mv_{r}^{2}$ and $E=E_{x}+E_{r}$. Putting it all together and
integrating over $\theta$ gives

$$$J_{\text{For}}=\dfrac{4\pi me}{h^{3}}\int dE_{x}\dfrac{k_{2}}{k_{1}}\lvert
T(E_{x})\rvert^{2}\,\times\\\
\int\limits_{0}^{\infty}dE_{r}f_{1}(E)[1-f_{1}(E^{\prime})]$$$ (5)

In a like manner one can also calculate the reverse current density to be

$$$J_{\text{Rev}}=\dfrac{4\pi me}{h^{3}}\int dE_{x}\dfrac{k_{1}}{k_{2}}\lvert
T(E_{x})\rvert^{2}\,\times\\\
\int\limits_{0}^{\infty}dE_{r}[1-f_{1}(E)]f_{1}(E^{\prime})$$$ (6)

The net current density is given by
$J_{\text{Net}}=J_{\text{For}}-J_{\text{Rev}}$, which becomes

$$$J_{\text{Net}}=\dfrac{4\pi me}{h^{3}}\int dE_{x}\lvert
T(E_{x})\rvert^{2}\,\times\\\
\int\limits_{0}^{\infty}dE_{r}\Big{[}f_{1}(E)[1-f_{1}(E^{\prime})]\,\dfrac{k_{2}}{k_{1}}-\\\
[1-f_{1}(E)]f_{1}(E^{\prime})]\,\dfrac{k_{1}}{k_{2}}\Big{]}$$$ (7)

Define the Fermi Factor $\mathcal{F}(E_{x})$ as

$$$\mathcal{F}(E_{x})=\int\limits_{0}^{\infty}dE_{r}\Big{[}f_{1}(E)[1-f_{1}(E^{\prime})]\,\dfrac{k_{2}}{k_{1}}-\\\
[1-f_{1}(E)]f_{1}(E^{\prime})\Big{]}\,\dfrac{k_{1}}{k_{2}}$$$ (8)

Carrying out the integral over $E_{r}$

$$$\mathcal{F}(E_{x})=\dfrac{k_{2}}{k_{1}}F_{1}(E_{x})-\dfrac{k_{1}}{k_{2}}F_{1}(E_{x}+eV_{b})\,+\\\
\frac{\Big{[}\dfrac{k_{1}}{k_{2}}-\dfrac{k_{2}}{k_{1}}\Big{]}}{(1-e^{-\beta
eV_{b}})}\Big{[}F_{1}(E_{x}+eV_{b})-e^{-\beta eV_{b}}F_{1}(E_{x})\Big{]}$$$
(9)

where

$F_{1}(E_{x})=\dfrac{1}{\beta}\text{Log}[1+e^{-\beta(E_{x}-\eta_{1})}]$

Figure 3: Fermi Factor for Bias voltage of $3$ and $5$ V for Al-vac-Cu. The
solid red line refers to $\mathcal{F}(E_{x})$, the dotted line refers to Non-
Pauli term, and dot-dashed line refers to Pauli term.

The first two terms are the Non-Pauli contribution; these are the terms that
one would get if Pauli Blocking were completely ignored. The third term in
equation (9) explicitely introduces Pauli Effects. The behaviour of these
terms $viz$ Non-Pauli, $\mathcal{F}(E_{x})$, and the Pauli term (third term
only), as a function of $E_{x}$ is displayed in Fig.3. The Fermi Function
$\mathcal{F}(E_{x})$ (solid red line) is seen to peak at $8.7\,\text{eV}$
$\approx\eta_{1}-eV_{b}$ for $eV_{b}=3\,\text{eV}$; and it peaks at
$6.7\,\text{eV}$ $\approx\eta_{1}-eV_{b}$ for $eV_{b}=5\,\text{eV}$.
Subsequently $\mathcal{F}(E_{x})$ is found to decrease rapidly to zero as
$E_{x}$ exceeds $\eta_{1}$. This is found true for all bias voltages. The
dotted line refers to the Non-Pauli(NP) term, which may have negative values
for low bias voltages. For these voltages the reverse current exceeds the
forward current. This may be due to the different flux ratios for the forward
and the reverse currents. The dot-dashed line refers to the third term of
equation (9). The third term is seen to go to zero at about
$E_{x}=\eta_{1}-eV_{b}$. This term contains Pauli Effects and it is studied in
more detail in Figures 4 and 5.

Figure 4: Third Term as a function of $E_{x}$ for bias voltages $2$V,$5$V,and
$10$V for Al-vac-Cu. Figure 5: Third Term as a function of $V_{b}$ for
electron energies of $2$eV,$5$eV,and $10$eV for Al-vac-Cu.

Figure 4 plots the third term for fixed bias voltages ($V_{b}=2$V,
$V_{b}=5$V,$V_{b}=10$V) as a function of electron energies in the range ($1$
eV to $12$ eV). The graph shows that Pauli Effects are significant for low
energy electrons and especially for low bias voltages. For
$E_{x}>\eta_{1}-eV_{b}$ the third term is essentially zero implying that there
is negligible Pauli Effect. The inset in this figure shows the third term for
$2$V rising rapidly by at least $4$ orders of magnitude for very low values of
$E_{x}$. Figure 5 plots the third term for fixed electron energies
($E_{x}=2$eV, $E_{x}=5$eV,$E_{x}=10$eV) as a function of the bias voltage in
the range ( $1$V to $10$V). The inset in this figure shows the third term for
$E_{x}=2$eV rising rapidly by at least $3$ orders of magnitude for very low
values of bias voltage. All curves in Fig.5 pass through zero at
$eV_{b}=\eta_{1}-\eta_{2}$ which in case of electrodes Al and Cu is given by
$V_{b}=4.7$ V. This is because $k_{1}=k_{2}$ for this value of bias, and the
zero level of the energy in the two electrodes become identical. For bias
voltages lesser than $eV_{b}=\eta_{1}-\eta_{2}$ the third term is positive and
its effect is to increase the net tunneling current beyond the non-Pauli
contribution. For bias voltages larger than this value the third term becomes
negative and the Pauli Effect will be to reduce the net tunneling current from
its Non-Pauli value.

The net current density becomes

$J_{\text{Net}}=\dfrac{4\pi
me}{h^{3}}\int\limits_{E^{\text{min}}_{x}}^{E^{\text{max}}_{x}}dE_{x}\lvert
T(E_{x})\rvert^{2}\mathcal{F}(E_{x})$ (10)

where the limits of the integral in equation (10) are given by

$E_{x}^{\text{min}}=\text{Max}[0,\Delta_{b}]\,\,\text{and}\,\,E_{x}^{\text{max}}=\eta_{1}+\phi_{1}$

$E_{x}^{\text{min}}$ is the minimum of energy in the energy Stage I and
$E_{x}^{\text{max}}$ is the maximum of energy in the energy Stage II. The
energy stages I, II, and III are described in the next section.

###  II.2 B. TUNNEL AMPLITUDE ${\bm{T(E_{x})}}$

Let the two planar electrodes described in the previous section be parallel to
each other and parallel to the y-z plane. Let their surfaces be at $x=0$ and
$x=d$ respectively, with the intervening space being occupied by vacuum. The
entire extent of the x-axis can be divided up into three spatial regions and
the potential energy of an electron in these three regions is given by

$U(x)=\left\\{\begin{array}[]{@{}ll@{}}0&x<0\\\ U_{II}(x)&0\leqslant
x\leqslant d\\\ \Delta_{b}&x>d\\\ \end{array}\right.$ (11)

where

$U_{II}(x)=(\eta_{1}+\phi_{1})-(\phi_{1}-\phi_{2}+eV_{b})\dfrac{x}{d}$ (12)

The wavefunctions in two regions $(x<0)$ and $(x>d)$ are given by

$\psi_{1}(x)=e^{ik_{1}x}+Re^{-ik_{1}x},\hskip 14.45377ptx<0$ (13)
$\psi_{3}(x)=Te^{ik_{2}x},\hskip 14.45377ptx>d$ (14)

where, $R$, $T$ are the reflection and transmission amplitudes. The
wavefunction $\psi_{b}$ in the barrier region $0\leqslant x\leqslant d$ is
given by the Schrödinger equation for the linear (trapezoidal) potential

$\frac{d^{2}\Psi_{b}}{dx^{2}}-(A-Bx)\Psi_{b}=0$ (15)

where

$A=\dfrac{2m}{\hbar^{2}}(\eta_{1}+\phi_{1}-E_{x})$ (16)
$B=\dfrac{2m}{\hbar^{2}d}(\phi_{1}-\phi_{2}+eV_{b})$ (17)

If the electrodes were at absolute zero temperature, only electrons at the
Fermi levels of the electrodes would undergo tunneling. At higher
temperatures, electrons with a wider range of energies will be involved in
tunneling. The energy range of such electrons can be divided into three Energy
Stages. Fig.6 shows the three Energy Stages for Al-vac-Cu system for which
$\eta_{1}=11.7$ eV, $\phi_{1}=4.16$ eV,$\eta_{2}=7.0$ eV, $\phi_{2}=4.81$ eV,
$V_{b}=7$ volts.

Figure 6: Energy Stages I ($0$ eV to $9.51$ eV), II ($9.51$ eV to $15.86$ eV),
and III ($>15.86$ eV) for Al-vac-Cu, and $eV_{b}$ = 7.0 eV
$(eV_{b}>\phi_{1})$. The upper limits are indicated by blue dotted line (Stage
I), pink dotted line (Stage II).

Energy Stage I: The barrier region in this stage does not contain a turning
point. The energy in this region ranges from
Max$[0,\Delta_{b}]<E_{x}<(\eta_{1}+\phi_{1}-eV_{b})$. Therefore, the term
$(A-Bx)$ will never be zero $\forall\,x\,\in\,[0,d]$.
Energy Stage II: In this stage, there is a turning point inside the barrier
region at $x=x_{T}=\frac{A}{B}$. The energy range in this stage is
$(\eta_{1}+\phi_{2}-eV_{b})<E_{x}<(\eta_{1}+\phi_{1})$ as shown in Fig.6. The
barrier is triangular and electrons tunneling out of this triangular barrier
are real (not virtual) electrons which have not entered the the second
electrode and those that do so enter the second electrode, contribute to the
tunneling current.
Energy Stage III: The energy of the electrons in this stage exceeds
$(\eta_{1}+\phi_{1})$ and these electrons should be emitted without having to
tunnel through any barrier suggesting that the transmission probability would
be very high. However for most values of $\phi_{1}$ the electron energy is so
far above the Fermi energy $\eta_{1}$ that the Fermi factor (for room
temperatures) is essentially zero and the contribution of electrons in this
energy stage at these temperatures to the current density is vanishingly
small. Therefore unless the electrode temperatures are very high (i.e.
comparable to the Fermi temperature of the electrodes), the contribution of
electrons in this Energy Stage will be neglected alltogether.

In Energy Stages I and II the Schrödinger equation (4) reduces to the Airy
equation.

$\frac{d^{2}\Psi_{b}}{dh^{2}}-h\Psi_{b}=0$ (18)

where $h=h(x)=\dfrac{A}{B^{2/3}}-B^{1/3}x$. The general solution of the Airy
Differential equation can be written as AS

$\psi_{b}(x)=C\phi^{(1)}(x)+D\phi^{(2)}(x)$ (19)

where

$\phi^{(1)}(x)=\dfrac{Ai[h(x)]}{Ai[h(\bar{x})]}\,\,\phi^{(2)}(x)=\dfrac{Bi[h(x)]}{Bi[h(\bar{x})]}$
(20)

The denominator of the functions $\phi^{(1)}$ and $\phi^{(2)}$ contain Airy
functions evaluated at $\bar{x}$ (a fixed value of x near the centre of the
barrier region). These denominators are sometimes required to avoid overflow
and underflow errors in the calculation of the wavefunctions for different
values of $x$ in the barrier. The constants $C$ and $D$ are determined from
matching of the wavefunction and its derivative at the relevant spatial
boundaries. In Energy Stage I these boundaries are at $x=0$ and $x=d$. In
Energy Stage II, two separate wavefunctions $\psi_{b1}$ and $\psi_{b2}$ have
to be determined for $x<x_{T}$, and $x>x_{T}$ respectively. The spatial
boundaries for the former are at $x=0$ and $x=x_{T}$ and that for the latter,
are at $x=x_{T}$ and $x=d$. The Tunneling amplitude is given by

$T(E_{x})=e^{-ik_{2}d}\Big{\\{}\begin{array}[]{@{}ll@{}}\psi_{b}(d)&\text{in
Stage I}\\\ \psi_{b2}(d)&\text{in Stage II}\end{array}$ (21)

Fig.7 shows that the $\rvert T(E_{x})\rvert^{2}$ increases with $E_{x}$ which
is to be expected.

Figure 7: Plot of Transmission Probability versus $E_{x}$ for Bias Voltage = 3
V $\&$ $d=5,10,15,20$ $\AA$

Fig.8, shows a plot of the product of calculated (Airy) $\rvert
T(E_{x})\rvert^{2}$ and the Fermi Factor $\mathcal{F}$ for metal electrodes
Al(tip) and Cu(sample). The product $\rvert
T(E_{x})\rvert^{2}\times\mathcal{F}$ increases exponentially upto Fermi level
$(\eta_{1}=11.7\,eV)$ of the left electrode and drops sharply thereafter, well
before reaching the upper limit of Energy Stage II
$(\eta_{1}+\phi_{1}=15.86\,eV)$. Thus significant contribution to the
tunneling current density arises from a narrow energy band (for $E_{x}$)
around $\eta_{1}$. The Fermi Factor is thus seen to damp out the contribution
to the tunneling current density for larger values of the kinetic energy
$E_{x}$ for the translational motion of the electrons in the x-direction.

Figure 8:  Plot of product of Transmission Probability and Fermi Factor versus
$E_{x}$ for Bias Voltage = 3 V and $d=5,\,10,\,15,\,20$ $\AA$ Figure 9: Plot
of Transmission Probability ( Airy compared with WKB ) versus $E_{x}$ for Bias
Voltage = 3 V and $d=5,10$ $\AA$ (Al-vac-Cu).

Fig.9 compares $\rvert T(E_{x})\rvert^{2}\times\mathcal{F}$ calculated using
Airy functions with that using WKB approximation for an intermediate bias of
$3$ V. For this bias voltage, this product $\rvert
T(E_{x})\rvert^{2}\times\mathcal{F}$ using the WKB approximation is found to
exceed the Airy function determined $\rvert
T(E_{x})\rvert^{2}\times\mathcal{F}$ for low tip-sample distances such as
$d=5\AA$. The two are nearly equal for $d=10\AA$, however its Airy value
exceeds its WKB value for $d=15\AA$. Thus the WKB approximation does not
provide an accurate dependence of the tunnel current densities on the tip-
sample distances. This feature suggests that determination of surface
topogrophies of samples would be inaccurate, if the image processing software
uses WKB approximation.

###  II.3 C. ${\bm{T(E_{x})}}$ FOR VERY LOW BIAS

In the low bias regime (0.1V to 1V) the Airy Function $Ai$ becomes very small
causing underflow problems and the Airy function $Bi$ becomes very large
causing overflow problems as in the calculation. Even introduction of
denominators as in equation (20), doesn’t seem to help. Instead it is more
useful to convert the Schrödinger equation to the following inhomogeneous
equation

$\mathscr{L}\psi(x)=f(x)\psi(x)$ (22)

where

$\mathscr{L}=\dfrac{d^{2}}{dx^{2}}-\kappa^{2}$ (23)
$\dfrac{\hbar^{2}\,\kappa^{2}}{2m}=\eta_{1}+\phi_{1}-E_{x}$ (24)
$\dfrac{\hbar^{2}}{2m}f(x)=-(\phi_{1}-\phi_{2}+eV_{b})\dfrac{x}{d}$ (25)

The term $f(x)$ in the equation (22), is usually very small for low bias
voltages and low contact potentials, permitting a perturbative solution of
equation (22).

Fig.8, shows that the maximum contribution to the tunneling occurs for $E_{x}$
lying in a very narrow energy band whose upper limit is $\eta_{1}$. For
$E_{x}\leqslant\eta_{1}$ the minimum value of
$\frac{\hbar^{2}\,\kappa^{2}}{2m}$ will be $\phi_{1}$. The condition that
$\frac{\hbar^{2}}{2m}|f(x)|$ be smaller than this minimum value of
$\frac{\hbar^{2}\,\kappa^{2}}{2m}$ is $eV_{b}<\phi_{2}$. Hence for low bias
voltages, $|f(x)|<\kappa$ is true and perturbative solutions are valid. For
$eV_{b}>\phi_{2}$ the Airy function methods will do fine.

Define the Green function for the operator $\mathscr{L}$ of equation (23) as

$\mathscr{L}G(x,x^{\prime})=\delta(x-x^{\prime})$ (26)

Let the wavefunction in the barrier region be

$\psi(x)=\phi_{0}(x)+\phi(x)$ (27)

where $\phi_{0}(x)$ is the solution to the homogeneous equation

$\mathscr{L}\phi_{0}(x)=0$ (28)

From equations (27) and (28), $\phi(x)$ will satisfy

$\mathscr{L}\phi(x)=f(x)\,[\phi_{0}(x)+\phi(x)]$ (29)
$\displaystyle\phi(x^{\prime})$
$\displaystyle=\int_{0}^{d}G(x,x^{\prime})f(x)[\phi_{0}(x)+\phi(x)]\,dx$
$\displaystyle+\,\text{end point terms}$ (30)

The end point terms ($viz$ at $x=0$ and $x=d$) can be made to vanish if the
Green function and the function $\phi$ are both made to vanish at $x=0$ and
$x=d$. Equation (II.3) can be iterated to obtain terms containing higher and
higher powers of $f(x)$. Here $f(x)$ will be regarded as small enough to
neglect higher than first order terms. Keeping only the first order term,
gives

$\phi(x^{\prime})=\int_{0}^{d}G(x,x^{\prime})f(x)\,\phi_{0}(x)dx$ (31)

and therefore

$\psi(x)=\phi_{0}(x)+\int_{0}^{d}G(x,x^{\prime})f(x^{\prime})\,\phi_{0}(x^{\prime})dx^{\prime}$
(32)

The tunneling amplitude $T^{Gr}$ is calculated using the wavefunction in
equation (32) and joining it smoothly with the wavefunctions given in
equations (13) and (14).

For very low bias, the energy of the electrons in the Energy Stage II is well
above the Fermi level of the first electrode and the Fermi Factor is extremely
small for these erergies. Therefore these electrons will make a negligible
contribution to tunneling and hence the Green function method outlined above
is not extended to the second stage. Instead, the second stage contribution
for very low bias voltages is zeroed out.

###  II.4 D. CURRENT FROM CURVED TIPS USING PLANAR MODEL CURRENT DENSITIES

As mentioned in the introduction it is necessary to go beyond the planar model
and treat the tip as a sharp pointed surface and let the sample be flat. For a
sharp tip, its radius of curvature at the point closest to the sample must be
small, so that the length of the electric field lines increases rapidly as one
goes away from the centre. Since the tunneling current densities fall off
exponentially with increasing electric field line length therefore these
current densities when integrated over the tip area will lead to a finite
value for the current.

A good model in which the tip is sharp and the sample is flat is provided by
using confocal hyperboloids for the tip and sample shapes (Fig.10). These
surfaces are coordinate surfaces in the prolate spheroidal coordinate
systemcutler Patil8 date russel SG jun moon morse Aizpurua . This is described
by

$x=\rho\cos\phi,\quad y=\rho\sin\phi,\quad z=a\xi\eta$ (33)

where $\rho=a\sqrt{\xi^{2}-1}\sqrt{1-\eta^{2}}$, and $a$ is a constant whose
value is determined by the tip radius of curvature $R$ and the tip sample
distance $d$ as $a=\sqrt{d(d+R)}$. The range of coordinate values are given by

$1\leqslant\xi\leqslant\infty,\,\,-1\leqslant\eta\leqslant
1,\,\,0\leqslant\phi\leqslant 2\pi$ (34)

The Laplace equation for the electrostatic potential in this coordinate system
becomes one dimensional and equipotential surfaces are confocal hyperboloids
which are a one parameter family of surfaces, whose equation in cartesian
coordinates is

$\frac{z^{2}}{a^{2}\eta^{2}}-\frac{(x^{2}+y^{2})}{a^{2}(1-\eta^{2})}=1$ (35)

The tip surface (treated here as an emitter) is represented by a hyperboloid
described by equation (35) in which $\eta=\eta_{tip}=d/a$ The sample surface
is characterised by $\eta=\eta_{sample}=0$ which is the $x-y$ plane.

Figure 10: Hyperboloid Tip and Spheroidal shaped field line.

The equation of the field line that passes through the point $(\xi,\phi)$ is
independent of $\phi$ and is given by

$\frac{\rho^{2}}{a^{2}(\xi^{2}-1)}+\frac{z^{2}}{a^{2}\xi^{2}}=1$ (36)

and the length of this field line stretching from $\eta=\eta_{sample}$ to
$\eta=\eta_{tip}$ is

$d_{fl}(\xi)=a\int_{\eta_{sample}}^{\eta_{tip}}\frac{\sqrt{\xi^{2}-\eta^{2}}}{\sqrt{1-\eta^{2}}}d\eta$
(37)

In Fig.10, a small shaded area on the tip is shown along with, a couple of
electric field lines, both lines, parametrized by $\xi$. This small shaded
area on the tip, and the corresponding area on the sample define a planar
tunnel junction of infinitesimal area, in which the length of the field line
$d_{fl}(\xi)$ plays the role of the tip sample distance. Thus the net current
density $J_{Net}(\xi)$ for this infinitesimal area planar tunnel junction is
given by equation (10). The product of the area of the shaded region $dS(\xi)$
with the current density of this infinitesimal tunnel junction gives the
current carried by it. The sum of these currents for all such junctions
covering a reasonably wide area on the tip will give the total net current for
the entire non planar tip. Thus planar model current densities can be used to
calculate the total net current for a given non planar shape of the tip. This
procedure has been well described by Saenz and Garcia SG .

###  II.5 THE CALCULATION

The total current density $J_{Net}(\xi)$ is obtained by integrating the
product of the tunneling probability and the Fermi Factor over the energy
$E_{x}$ of the tunneling electrons for the translational motion along the
x-direction. The tunneling amplitude $T=T[E_{x},d_{fl}(\xi)]$ depends not only
upon the energy $E_{x}$ but also on the location of the point on the tip
(dictated by $\xi$) from where tunneling occurs. The current density
$J_{Net}(\xi)$, given by an equation similar to equation (10) is a function of
$\xi$ through its dependence on the length of the field line $d_{fl}(\xi)$,
and since $J(\xi)$ decreases rapidly with $d_{fl}(\xi)$ it will also decrease
equally rapidly with $\xi$. The range of energies $E_{x}$ in the integration
spans those in Energy Stages I and II described earlier in section B.

Figure 11: Plot of Current Density versus Length of field line with
$R=20A^{0}$, $d=10A^{0}$ and Bias $=3$ V for Pt-vac-Ag

The net current is given by

$I=\int_{1}^{\xi_{max}}J(\xi)\,dS(\xi)$ (38)

where the upper limit $\xi_{max}$ is chosen to have a value of $5$, because
the current density $J_{Net}(\xi)$ is essentially zero beyond this value as
shown in Fig.11. In this figure it can be seen that the current density drops
by roughly $20$ orders of magnitude as $\xi$ increases from $\xi=1$ to
$\xi=5$. The surface element $dS(\xi)$ in equation (38) is given by

$dS(\xi)=2\pi
a^{2}\sqrt{(\xi^{2}-\eta_{tip}^{2})\,\,(1-\eta_{tip}^{2})}\,d\xi$ (39)

is the area of the shaded region on the tip as shown in Fig.10. The curved
line in this Fig.10 is the field line correponding to the parameter $\xi$
$(1\leqslant\xi\leqslant\xi_{max})$.

###  II.6 RESULTS AND DISCUSSION

Figure 12 shows calculated I-V characteristics for identical electrodes
(W-vac-W) in which the tip is negatively biased, for several tip-sample
distances. The current increases almost exponentially with bias voltage at all
distances.

Figure 12: Plot of Current versus Bias Potential for $R=20\,\,A^{0}$ and
$d=5,\,6,\,8,\,10\,\,A^{0}$ for similar electrodes Figure 13: Plot of Current
versus Bias Potential for $R=20\,\,A^{0}$ and $d=5,\,6,\,8,\,10\,\,A^{0}$ for
dissimilar electrodes

Figure 13 shows calculated I-V characteristics for dissimilar electrodes (Pt-
vac-Ag) and the I-V curves behave similar to the case of identical electrodes,
$viz$ the current increases almost exponentially with bias voltage for all
distances.

Figure 14: Plot of Current versus Bias Potential $(<0.1\,V)$ for
$R=20,30,50A^{0}$ and $d=5A^{0}$ for W-vac-W

Ting et. al.ting  have reported values of tunneling currents in the W-vac-W
electrode system for very low bias voltages and for several cross currents
$I_{\text{cross}}$, including for $I_{\text{cross}}=0$. Here
$I_{\text{cross}}$ is a current made to flow in the surface of the sample. For
$I_{\text{cross}}=0$, their experimental setup corresponds to the tip-vacuum-
sample system for which tunneling currents are calculated in this paper.
Tunneling currents of $0.1\text{nA}$ are reported by Ting et. al.ting  for low
bias voltages in the range $[0\,\text{to}\,70\,\text{mV}]$. Fig.14 plots the
calculated tunneling current for (W-vac-W) for very low bias for a fixed tip
sample distance $d=5\,\AA$ and several tip radii of curvature
$(R=20,30,50\,\AA)$. These plots show the tunneling current increasing with
the bias voltage as expected. Also the tunneling currents are seen to increase
with increasing tip radii. Further for $R=30$ $\AA$, a good numerical
agreement with results for $I_{\text{cross}}=0$ of Ting et. al.ting  is
obtained.

Figure 15: Plot of Current versus Tip - Sample distance for $R=20A^{0}$ and
Bias Voltage = 3, 4 and 5 V for similar electrodes Figure 16: Plot of Current
versus Tip - Sample distance for $R=20A^{0}$ and Bias Voltage = 3, 4 and 5 V
for Dissimilar electrodes

Figures 15 and 16 show the plots of calculated tunneling currents $I$ versus
tip sample distance $d$ in similar (Al-vac-Al) and dissimilar (Pt-Vac-Ag)
pairs of electrodes. These plots are for bias voltages of $3$ V, $4$ V and $5$
V, and tip radius of curvature of $R=20$ $\AA$. Therefore $I$ is found to
decrease exponentially with increasing $d$. This behaviour is also
qualitatively reproduced by almost all calculations of tunneling current
density including those that use the WKB approximation simmonsI simmonsII SBPZ
Hofer .

For fixed bias voltage $V$ and tip sample distance $d$, increasing the tip
radius $R$, leads to a flatter tip, which means that more area in the vicinity
of the center lies closer to the sample and the current density for these
areas is greatly increased. Since the current $I$ is an integral of the
current density over tip area, $I$ is expected to increase with increasing $R$
for all $V$ and $d$. Fig.17 indeed confirms this behaviour by showing how the
$I-V$ curves correspond to higher currents for larger $R$. The general
behaviour of $I$ verses $V$ is not much influenced by $R$.

Figure 17: Plot of Current versus Bias Voltage for $d=10$ $A^{0}$ and
$R=5,20,40,100\,A^{0}$ Figure 18: Log Plot of Current Ratio (CR) as a function
of Bias Voltage ($2$V to $5$V) for $R=50\AA$ and $d=5,7,10\AA$. Figure 19:
Plot of Current Ratio (CR) as a function of Bias Voltage ($4$V to $8$V) for
$R=50\AA$ and $d=5,7,10\AA$

The Pauli Effects have been discussed in some detail following Figures 3,4,and
5. The effect of the third term on the tunneling currents is studied by
considering the Current Ratio $\text{CR}=1-(I_{\text{Non-Pauli}}/I)$. Figures
18,19 show that the ratio of the third term contribution to the net tunneling
current which is labelled as the CR is quite small for all bias voltages. The
value of the CR and hence the strength of the Pauli Effect decreases very fast
with increasing tip-sample distance. The CR is positive for bias voltages less
than $eV_{b}=\eta_{1}-\eta_{2}=4.7$eV, for which voltage range
$I_{\text{Net}}>I_{\text{Non-Pauli}}$ and for bias voltages greater than this
value the CR is negative, implying that $I_{\text{Net}}<I_{\text{Non-Pauli}}$.
The Pauli Effect has a larger magnitude at low bias voltages as compared to
higher bias voltages.

Figure 20:  Plot of r ( Measure of lateral Resolution ) versus R ( Radius of
curvature of the tip ) for Bias Voltage = 5 V and Tip sample distances = 5, 7,
10, 15 $\AA$

The ratio $\dfrac{I}{J_{Net}(\xi=1)}$ gives an area $\sigma$, which can be
regarded as the area of a circle of lateral resolution; i.e., sample profile
features that lie within this circle would be poorly resolved. The radius
$r=\sqrt{\sigma/\pi}$ of this circle would therefore be a good measure of the
lateral resolution limit, and can be named as the Lateral Resolution Parameter
(LRP).

Figure 21:  Plot of r ( Measure of lateral Resolution ) versus R ( Radius of
curvature of the tip ) for Tip sample distances = 10 $\AA$ and Bias Voltage=1,
3, 5 V

Fig.20 and Fig.21 show the dependence of the $r$ (LRP) with the tip radius $R$
which is varied from $(5\,\AA\,\text{to}\,100\,\AA)$. Both figures show that
sharper the tip ($viz$ lower the value of $R$), the greater is the lateral
resolution ($viz$ smaller the value of $r$). Fig.20 shows that for fixed bias
the lateral resolution is degraded with increasing tip sample distance. In
Fig.21 the tip sample distance is kept constant, while the bias voltage is
changed. This figure shows that the lateral resolution is degraded as the bias
voltage increases; which feature is found to be more pronounced at higher,
rather than at lower bias voltages.

###  II.7 SUMMARY

The model used for calculating tunnelling current for the biased M-V-M tunnel
junctions treats the effective barrier potential as a trapezoidal potential.
Exact analytic expressions for the tunneling probability and current are
obtained by solving the Schrödinger’s equation for the trapezoidal potential,
for which the wavefunctions can be expressed in terms of Airy functions.
Electrons with x-directed translational energies ranging from $0$ eV to
$\eta_{1}+\phi_{1}$ contribute to tunneling.

The current is calculated for a wide range of input parameters such as tip-
sample distance, bias voltage, and radius of curvature of the tip, both
electrodes being at room temperature. Pauli Effects are explicitly introduced
and their behaviour as a function of bias voltage and tip-sample distance is
studied. A measure of lateral resolution of the tunnel junction is introduced
and its behaviour as a function of tip radius of curvature for several values
of bias voltage, tip sample distance is also presented.

## References

  * (1) G. Binnig, H. Rohrer, Ch. Gerber and E. Weibel ”Surface Studies by Scanning Tunneling Microscopy” Phys. Rev. Lett. 1982, 49:57
  * (2) J. G. Simmons ”Generalized Formula for the Electric Tunnel Effect between Similar Electrodes Separated by a Thin Insulating Film” Journal of applied Physics 1963, 34:1793-1803.
  * (3) Thomas E. Hartman and Jay S. Chivian ”Electron Tunneling Through Thin Aluminum Oxide Films” Physical Review 1964, 134:A1094-A1101.
  * (4) J. Tersoff, D. Hamann ”Theory of the scanning tunneling microscope” Phys. Rev. B 1985, 31:805-813.
  * (5) J. Bardeen ”Tunneling From a Many-Particle Point Of View” Phys. Rev. Lett. 1961,6:57-59.
  * (6) J. Bardeen ”Surface States and Rectification at a Metal SemiConductor Contact” Physical Review 1947, 71:717-727.
  * (7) C. Julian Chen ”Role of Tip Material in Scanning Tunneling Microscopy” Mat. Res. Soc. Symp. Proc. Vol. 159. 1990.
  * (8) S. Banerjee and P. Zhang ”A generalized self-consistent model for quantum tunneling current in dissimilar metal-insulator-metal junction” AIP Advances 9, 085302 (2019)
  * (9) J. J. Saenz, R. Garcia ”Near field emission Scanning tunelling Microscopy” Appl. Phys. Lett. 1994, 65:3022-3024.
  * (10) R. H. Fowler and L. Nordheim ”Proc. Roy. Sot. (London) A119 (1928) 173.
  * (11) J. G. Simmons ”Electric Tunnel Effect between Dissimilar Electrodes Separated by a Thin Insulating Film” Journal of applied Physics 1963, 34:2581-2590.
  * (12) Q.Q. Shu a, Y. Jiang b, S. Meng b, G. Lin b, W.G. Ma ”Exactly solvable model for metal-insulator-metal stepped boundary tunnel junctions” Thin Solid Films 2002, 414:136–142.
  * (13) C. Bracher, M. Riza, and M. Kleber ”Propagator theory of scanning tunneling microscopy” Physical Review B Volume 56, Number 12
  * (14) W. Ebeling, D Blaschke, R Redmer, H Reinholz, and G Ropke ”The influence of Pauli blocking effects on the properties of dense hydrogen” J. Phys. A Math. Theor. 2009, 42:214033(1-8pp).
  * (15) M. Abramowitz and I. A. Stegun ”Handbook of Mathematical Functions”: National Bureau of Standards Applied Mathematics Series-55, June 1964, Tenth Printing, December 1972, p. 446.
  * (16) P.H. Cutler, Jun He, J. Miller, N. M. Miskovsky, B. Weiss and T. E. Sullivan ”Theory of electron emission in high fields from atomically sharp emitters: Validity of the Fowler-Nordheim equation” Progress in Surface Science, Vol. 42, pp. 169-185 1993
  * (17) S. Patil, A. Kulkarni, C.V Dharmadhikari ”Study of the electrostatic force between a conducting tip in proximity with a metallic surface: Theory and experiment” Journal of Applied Phys. 2000, 88:6940-6942.
  * (18) K. Date, A. Kulkarni, C.V Dharmadhikari ”Force on a Conducting Tip near a Metallic Surface Coated with a Polarizable Dielectric Layer: Theory and Experiment” e-J. Surf. Sci. Nanotech. Vol. 9 (2011) 206-209
  * (19) A. M. Russell ”Electron Trajectories in a Field Emission Microscope” Journal of Applied Science 1962, 33:970
  * (20) J. He, P.H. Cutler, N.M. Miskovsky,T.E. Feuchtwang ”Derivation of the image interaction for non-planar pointed emitter geometries: application to field emission I-V characteristics” Surface Science 1991, 246: 348-364.
  * (21) P. Moon and D. E. Spencer ”Field Theory Handbook” 2nd Edition, Springer-Verlag Berlin Heildelberg New York (1971) Table 1.06, pp-28
  * (22) P. Morse and H. Feshbach ”Methods of Theoretical Physics” McGraw Hill Book Company, Inc (1953) pp-661, Volume-I
  * (23) J. Aizpurua, S. P. Apell, R. Berndt ”Role of tip shape in light emission from the scanning tunneling microscope” Physical Review B Volume 62, Number 3, (2000-I)
  * (24) T. Xie,M. Dreyer,D. Bowen, C. Krafft ”On local sensing of spin Hall effect in tungsten films by using STM-based measurements” IEEE Transactions on Nanotechnology PP(99) July 2017
  * (25) W. A. Hofer,A. S. Foster,A. L. Shluger ”Theories of scanning probe microscopes at the atomic scale” Reviews of Modern Physics, Volume 75,1287, October 2003.