SEMITIP V6, program MultPlane3
Introduction
This program computes the electrostatic potential and the resulting tunnel current between a metallic tip and a nonuniform (inhomogeneous) semiconducting sample, for a fully 3-dimensional geometry. The nonuniformity in the semiconductor is of the type that affects the electrostatic potential, e.g. layers of differing doping in a heterostructure. Additionally, for the tunnel current, a region with differing band edge positions can be defined (e.g. as would be encountered for a quantum dot at the surface). The tunnel current is computed using a plane wave expansion to solve Schrödinger equation, for a restricted spatial region centered around the point on the surface opposite the tip apex.
The plane wave basis is formed by matching plane waves in the semiconductor with decaying exponentials in the vacuum, as described in Ref. [1]. In the z-direction a periodic system is formed consisting of the semiconductor slab and a vacuum region, with the basis states having either odd or even parity relative to the center of the slab.
In the x and y directions a periodic system as also formed, with sine and cosine functions used for the basis in the semiconductor. The solution to Schrödinger's equation within planecurr3.f assumes a mirror plane lying along the x-axis (i.e. only cosine basis states in the y direction).
The most critical parameters in the plane wave solution to Schrödinger's equation are the number of plane waves (k-points) used. This is set in the PARAMETER statement within planecurr3.f by the nkx, nky, and nkz parameters. A minimum value of 4 for these parameters will allow the program to run quickly, but with quite unreliable results. Values of 8 or 10 are much more reliable, but time consuming. To provide some reduction in computation time, the parameter on line 53 of the FORT.9 input file can be used. A setting of 1 for that parameter means that all the valence bands and the conduction band are evaluated for each voltage, whereas a setting of 0 means that only the valence bands are used for V<0 and the conduction band for V>=0. The latter restriction on the bands used is sufficient for many applications (unless accumulation or inversion occurs, in which case the plane wave method is not applicable anyway since it is not self-consistent, as further discussed in the
SEMITIP V6 Technical Manual).
Larger values for the number of plane waves are very time consuming and require parallel processing for routine evaluation of the tunnel current.
The spatial size of the periodic supercell used in the solution to Schrödinger's equation is determined by the number of plane waves and by the energy limits for the plane waves as specified in the FORT.9 input file (with different limits used for each band due to their different effective masses). For number of k-points NK, effective EFFM, and energy limit EMAX, the extent of the supercell is NK*2*PI/SQRT(C*EFFM*EMAX) where C equals 2m/hbar^2 where m is the free electron mass. Here, SQRT(C*EFFM*EMAX)/NK is the minimum nonzero wavevector, so the extent of the supercell (in each direction) is just the size corresponding to maximum wavelength of the problem. EMAX is generally kept fixed, corresponding to a scale on the order of 1 eV such that states of that energy (as needed for the voltage range of a spectrum) are included. Hence, increasing NK produces an increase in the size of the supercell. The sampling of spatial points within the supercell is set in the program to occur on a grid with 4*NK points in each dimension. For the maximum wavevector, which has NK wavelengths in the supercell, it will thus be sampled at 4 points along its wavelength.
Parameters internal to the planecurr3.f routine are:
-
IPOT - has value of 1 to include the computed potential from semitip3.f in the eigenvalue solution to Schrödinger's equation, or a value of 0 to not include it (this parameter is used mainly for debugging purpose). Default value is 1.
-
IVAC - has value of 1 to include the vacuum region in the potential used for the eigenvalue solution to Schrödinger's equation, or a value of 0 to not include it (without it, the problem solved consists only of a periodic semiconductor region, with sine and cosine functions used as basis states). Default value is 1.
-
Several parameters that affect the output of intermediate values such as wavefunction values or potential maps, as described in the Output section below.
Version information
Version 6.3; see top of MultPlane-6.3.f source code in example(s) below for prior version information.
Usage
Source code files for the program are available in the example(s) provided below. In order to specify the inhomogeneities in the semiconductor, it is necessary to modify the source code as follows:
-
The routine IGETREG (at the bottom of the source code) defines the various regions in the semiconductor. The input parameters x, y, and S are the X, Y, and z positions in the semiconductor, respectively. The routine should return an integer value giving the number of the region for those specified positions (i.e. return a value of 1 to use the first set of parameters in the FORT.9 input file defining the semiconductor region, a value of 2 for the second set of parameters, etc.).
-
The maximum number of region types in the semiconductor is specified by the parameter NREGDIM and the maximum number of area types on the surface by the parameter NARDIM. The values of these parameters are set in the main programs and its associated routines by PARAMETER statements. By default these values are set to 2, but they can be increased by changing all of the PARAMETER statements.
Additional details of these source code modifications can be found in the
SEMITIP V6 Technical Manual. With the appropriate modifications, the source code of the main program, together with the subroutines listed below, must by re-compiled and linked in order to produce the executable code.
All input for the program (other than that modifications to the source code detailed above) comes from the file FORT.9, which is also included in the example(s) provided below.
Output
Output from the program is contained in the following files
(output depends on the value of the output parameter IWRIT as specified
in the input file FORT.9):
- FORT.10 - gives the numerical results for the following quantities:
- tip radius of curvature (nm)
- tip-sample separation (nm)
- sample-tip bias voltage (V)
- contact potential (eV)
- Pot0 - the surface potential at the semiconductor surface (eV).
- FORT.11 - provides the electrostatic potential energy (eV) along the z axis, as a
function of distance (output for IWRIT>=1). Also, the electrostatic potential plus the
vacuum barrier energy is output to FORT.95, and the energy of the valence and conduction
band edges (as used in computing the tunnel current) are output to FORT.96 and FORT.97, respectively.
- FORT.12 - provides the potential (eV) along the surface, for specified angle from x-axis, as a function
of the radial distance from the central axis (output for IWRIT>=1)
- FORT.13 - gives the entire array of potential values (eV) (output for IWRIT>=3); see
VERSION 6 Technical Manual
for more details.
- FORT.14 - provides the current (A/nm^2) (column 2) as a function of sample voltage (V) (column 1). Also, columns 3 and 4 provide the contributions to the current of extended states and
localized states, respectively. Separate contributions from the valence band and conduction band
go in to FORT.91 and FORT.92, respectively.
- FORT.15 - provides the conductance dI/dV (A/(V nm^2)) (column 2) as a function of sample voltage (V) (column 1). Also, columns 3 and 4 provide the contributions to the conductance of extended states and
localized states, respectively. Separate contributions from the valence band and conduction band
go in to FORT.93 and FORT.94, respectively.
- FORT.16 - gives an exact copy of the output to the console
- FORT.17 - provides the charge densities on the central axis (column 2) as a function of z-distance along the central axis (column 1). Also, columns 3 and 4 provide the contributions to the charge densities of extended states and localized states, respectively.
- FORT.18 - provides the surface charge densities (column 2) as a function of radial distance away from the central axis (column 1). Also, columns 3 and 4 provide the contributions to the charge densities of extended states and localized states, respectively.
- FORT.19 - provides the surface charge density (column 2) at the point opposite the tip apex (i.e. on the central axis), as a function of bias voltage. Also, columns 3 and 4 provide the contributions to the charge densities of extended states and localized states, respectively.
- FORT.31, FORT.41 (odd and even parity, respectively) - z values (column 1) and band edge potential (column 2) (output for IWRIT>=5).
- FORT.32, FORT.42 (odd and even parity, respectively) - z values (column 1) and wavefunction values for the first five states (columns 2-7). The states that are output here can be changed by changing the values of IWF1, ... IWF5 in the program (output for IWRIT>=5).
- FORT.33, FORT.43 (odd and even parity, respectively) - x values (column 1) and band edge potential (column 2) (output for IWRIT>=5).
- FORT.34, FORT.44 (odd and even parity, respectively) - x values (column 1) and wavefunction values for the first five states (columns 2-7). The states that are output here can be changed by changing the values of IWF1, ... IWF5 in the program (output for IWRIT>=5).
- FORT.35, FORT.45 (odd and even parity, respectively) - eigenvalues of Hamiltonian matrix (output for IWRIT>=4).
- FORT.40 - ordered list of eigenvalues (column 1), magnitude squared of wavefunction at surface (column 2), z-derivative at surface of log of wavefunction (column 3), expectation value of parallel wavevector (column 4). The values in column 4 are used for computing the decay constant for states in the vacuum as detailed in Ref. [1] (output for IWRIT>=4).
- FORT.80, FORT.81 (odd and even parity, respectively) - energy eigenvalues (column 1) and local density of states at the point on the surface opposite the tip apex (column 2) (output for IWRIT>=6).
- FORT.91, FORT.92,... - listing of surface charge density (column 2) vs. Fermi energy (column 1), for various areas (output for IWRITE>=3)
- img1.PGM, img2.PGM, img3.PGM, ..., img9.PGM provide images (PGM format) of the potential along planes located at z-values of 0 (i.e. the surface), first z-value, second z-value, etc. For a setting of IPOT=0 in the program, these images provide a convenient means of examining a localized potential such as that produced by a quantum dot (i.e. to ensure that the code defining the quantum dot is correct). Parameters for the gray-scale plotting of the images can be modified in the program code (output for IWRIT>=7).
All of the parameters in the program can be varied using the input file FORT.9, with the exception of the array sizes, the specification of a surface state density other than a uniform or Gaussian shaped one, the specification of spatial arrangement of bulk or surface charge density, and the specification of local variations in the band edge potential. The latter is accomplished with the routine DELVBEDGE and DELCBEDGE, which provide the change in band edge locations relative to the overall band gap defined in the FORT.9 input file. See
SEMITIP V6 Technical Manual
for additional information on these user-defined functions. Modification of those functions
can be accomplished by changing the source code of the program. The source code is available, in the following files (version numbers follow the dash in the names):
-
MultPlane3-6.3.f -
main program, available in example(s) below.
-
BoundWell3E-6.0.f -
routine for forming plane wave basis states for conduction band (including extended states).
-
BoundWell3V-6.0.f -
routine for forming plane wave basis states for valence band.
-
contr3-6.0.f -
routine for outputting contour plot.
-
dgsect-6.0.f -
general purpose double precision Golden Section search routine, used in forming plane wave basis.
-
GetCurr-6.0.f -
sums up states, to form the tunnel current.
-
gsect-6.0.f -
general purpose Golden Section search routine, for handling nonlinear aspects of the electrostatic solution.
-
planecurr3-6.1.f -
performs the eigenvalue solution of Schrödinger's equation with a plane wave basis.
-
PlotGray-6.0.f -
make gray scale plots of 2D functions (cuts of the potentials).
-
potcut3-6.0.f -
takes a cut along the the central axis (r=0) of the potential from semitip3, used to form the potential in the vacuum.
-
potperiod3-6.0.f -
expands the potential evaluated on the grid from semitip3.f into a finer grid with uniform spacing, as needed by planecurr3.f.
-
rseispackALL.f -
EISPACK routines for solving the eigenvalue problem of a real, symmetric, double precision matrix.
-
semirhomult-6.0.f -
routines for computing semiconductor charge densities.
-
semitip3-6.1.f -
performs the detailed finite element solution of Poisson's equation.
-
surfrhomult-6.2.f -
routines for handling surface charge densities.
All routines are written in Fortran. The source code can be downloaded
directly from the above locations, and it can be compiled and linked
on any platform. Sample input and output from the program is shown in the examples below.
Illustrative Examples of Running the Code
-
example 1
References
1. S. Gaan, G. He, R. M. Feenstra, J. Walker, and E. Towe, Size, shape, composition, and electronic properties of InAs/GaAs quantum dots by scanning tunneling microscopy and spectroscopy, J. Appl. Phys. 108, 114315 (2010). For preprint, see
http://www.cmu.edu/physics/stm/publ/93/.