Superlattices in QSGW

This tutorial shows how to construct a superlattice using Questaal tools.

A short-period superlattice of InAs and GaSb is constructed on the [001] axis and a self-consistent QSGW potential is made for it. The superlattice is built with the superlattice editor.

The potential for longer period superlattices are assembled by using the sigma editor, embedding the QSGW self-energy into bulk regions away from the interface, calculated from QSGW calculations for bulk InAs and GaSb. This makes it possible to generate QSGW quality band structure for much large superlattices that would be possible by direct QSGW calculations.

A separate tutorial also shows how to build superlattices. Its purpose is to explain how to calculate transport in a metallic, magnetic trilayer using a Landauer-Buttiker framework.

Zincblende InAs

• Input file and site file
cp init.dat init.inas
blm -vinas=t --sort=no --shorten=no --wsitex,quad1 --findes,rmin=1.3,float,1spec --ehmx=-.3 --nit=20 --gw~gcutb=2.7~gcutx=2.3~pbtol=3e-4,3e-4,1e-3,1e-3,1e-2 --nk=6 --nkgw=4 --gmax=9.2 --eloc=-3.5 --ctrl=ctrl inas
mkdir inas-zb
cp ctrl.inas site.inas inas-zb

Zincblende GaSb

• Site file
cp init.dat init.gasb
blm -vinas=f --wsitex,quad1 --findes,rmin=1.3,float,1spec --ctrl=ctrl gasb
mkdir gasb-zb
cp ctrl.gasb site.gasb gasb-zb

Input files for superlattices

• Conventional cell(1) for InAs
echo .5 .5 0 .5 -.5 0 0 0 1 | lmscell inas --wsite~short~fn=site
lmscell inas '--stack~sort@h3-z/1000~wsite@short@fn=site'
mkdir inas2gasb0
cp site.inas inas2gasb0/site.inasgasb
• Conventional cell(1) for sheared GaSb
echo .5 .5 0 .5 -.5 0 0 0 1 | lmscell gasb --wsite~short~fn=site
lmscell ctrl.gasb -vz=11.5/11.43 '--stack~scale=1/z~stretch=z^3~sort@h3-z/1000~wsitex@short@fn=site'
mkdir inas0gasb2
cp site.gasb inas0gasb2/site.inasgasb

ctrl, basp, and site files for short-period [001] superlattices

• ctrl file for (InAs)2(GaSb)2, which can be applied to any [001] superlattice. Also make the (InAs)2(GaSb)2 site file.
cp ./inas-zb/* .
cp inas2gasb0/site.inasgasb site.inas
cp inas0gasb2/site.inasgasb site.gasb
lmscell inas '--stack~show~file@site.gasb@dpos=0,0,-.032~show~wsite@short'
cp sites.inas sitein.inasgasb
blm --molstat --ewald=1d-10 --rdsite --omax=-.04 --addes,end --nosort --shorten=no --findes,rmin=1.3,float,1spec --ehmx=-.3 --nit=20 --gw~gcutb=2.7~gcutx=2.3~pbtol=3e-4,3e-4,1e-3,1e-3,1e-2 --nk=6,6,6 --nkgw=4,4,4 --gmax=9.2 --eloc=-3.5 --ctrl=ctrl inasgasb
• basp file for (InAs)2(GaSb)2, which can be applied to any [001] superlattice.
lmfa -vles=1 --basfile=basp ctrl.inasgasb
lmfa ctrl.inasgasb
echo -vnk2=4 -vnkgw2=2 > switches-for-lm
mkdir -p inas2gasb2/relax0
cp ctrl.inasgasb basp.inasgasb site.inasgasb switches-for-lm inas2gasb2/relax0
• Complete the input file setup for the two conventional cells
echo -vnk2=4 -vnkgw2=4 > switches-for-lm
cp ctrl.inasgasb basp.inasgasb switches-for-lm inas2gasb0
cp ctrl.inasgasb basp.inasgasb switches-for-lm inas0gasb2
• Site file and other essential inputs for the (4,4) and (6,6) superlattices
lmscell inas '--stack~show~file site.inas dup=1~file@site.gasb@dpos=0,0,-.032@dup=2~show~wsite@short'p
echo -vnk2=2 -vnkgw2=1 > switches-for-lm
cp sites.inas site.inasgasb
mkdir -p inas4gasb4/relax0
cp ctrl.inasgasb basp.inasgasb site.inasgasb switches-for-lm inas4gasb4/relax0
lmscell inas '--stack~show~file site.inas dup=2~file@site.gasb@dpos=0,0,-.032@dup=3~show~wsite@short'
echo -vnk2=2 -vnkgw2=1 > switches-for-lm
cp sites.inas site.inasgasb
mkdir -p inas6gasb6/relax0
cp ctrl.inasgasb basp.inasgasb site.inasgasb switches-for-lm inas6gasb6/relax0

LDA, QSGW self-consistency, energy bands in conventional cells

• InAs and GaSb conventional cell, LDA self-consistency + QSGW self-consistency
cp inas2gasb0/* .
lmfa ctrl.inasgasb
mpix -n 16 lmf ctrl.inasgasb

The preceding applies to InAs. Use the same instructions for GaSb, substituting inas2gasb0inas0gasb2.

lmfgwd --job=-1 ctrl.inasgasb cat switches-for-lm
lmgwsc --wt --insul=28 --mpi=16,16 --tol=2e-5 --sym --getsigp --save= inas2gasb0 inasgasb > out.job
mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm; mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm --rs=1,0 -vnit=1 >  out.lmf
cp rst.inasgasb  out.lmf inas2gasb0
• Energy bands in conventional cells, with SO coupling, with and without scaling Σ0
lmf ctrl.inasgasb --wsig:fbz; cp sigm2.inasgasb sigm.inasgasb
mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1
mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1 --band~fn=syml
echo -7,6 / | plbnds -nocol -fplot~sh~scl=.7 -ef=0 -scl=13.6 -lbl=L,G,X bnds.inasgasb
open fplot.ps
mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1 --mixsig=0.8
mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1 --band~fn=syml --mixsig=0.8
echo -7,6 / | plbnds -nocol -fplot~sh~scl=.7 -ef=0 -scl=13.6 -lbl=L,G,X bnds.inasgasb
open fplot.ps

LDA + QSGW self-consistency in superlattices

Set shell variable dir to be one of the following inas2gasb2/relax0, inas4gasb4/relax0, inas6gasb6/relax0

Initially make a self-consistent density within the LDA

cp $dir/* . rm -f rst.inasgasb mixm.inasgasb lmfa ctrl.inasgasb mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm QSGW self-consistency lmfgwd --job=-1 ctrl.inasgasb cat switches-for-lm lmgwsc --wt --insul=$nel --mpi=16,16 --tol=2e-5 --sym --getsigp --save= $dir inasgasb > out.job mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm; mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm --rs=1,0 -vnit=1 > out.lmf cp rst.inasgasb out.lmf$dir

Partial DOS and reference energies

The steps below apply to any of the superlattices. Set shell variables according to the table.

• Generate moms file for partial DOS
cp $dir/* . ; rm -f mixm.inasgasb; ln -s sigm sigm.inasgasb mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm --rs=1,0 -vnit=1 lmf ctrl.inasgasb --wsig:fbz; cp sigm2.inasgasb sigm mpix -n 16 lmf -vso=1 ctrl.inasgasb -vnk1=$nk1 -vnk2=$nk2 --rs=1,0 -vnit=1 --rsig:fbz --pdos~mode=0~sites=$sitelist --mixsig=.80 --pr35,40 > out.pdos
cp out.pdos $dir Partial DOS and integrated DOS lmdos -vso=1 ctrl.inasgasb -vnk1=$nk1 -vnk2=$nk2 --rs=1,0 -vnit=1 --rsig:fbz --pdos~mode=0~sites=$sitelist --mixsig=.80 --dos:npts=2001:window=-1,1
echo 0.999 4 -3 2 | pldos -fplot~sh~long~open~tmy=.25~dmin=0.40~xl=E -esclxy=13.6 -ef=0 -lst="$doschan" dos.inasgasb cp dosp.dat$dir/pdos.80
open fplot.ps

Integrated DOS

lmdos -vso=1 ctrl.inasgasb -vnk1=$nk1 -vnk2=$nk2 --rs=1,0 -vnit=1 --rsig:fbz --pdos~mode=0~sites=$sitelist --mixsig=.80 --dos:npts=40001:window=-3,1:rdm:idos cp dos.inasgasb$dir/idos.80
fplot -frme 0,.7,0,.4 -x -3,0 -y 0,7 -tmx .1:5 -colsy 1+1:1+4 -qr dos.inasgasb
open fplot.ps

In d DOS in the (6,6) superlattice

fplot -frme 0,.7,0,.4 -x -1.2,-1 -y 0,5 -tmx .05 -1p -map -i 'x1>-1.2&x1<-1' 1:nr  -colsy 4:24:4 -lt 1,col=0,0,0 -qr dos.inasgasb
open fplot.ps

Valence Band Offset

To be completed; for now follow the Tutorial

Relaxation of the interface

To be completed; for now follow the Tutorial

Longer Period superlattices

To be completed

Preliminaries

This tutorial assumes ~/lm is the top-level directory for the Questaal repository, and that Questaal executables are in your path.

This tutorial uses mpix to execute MPI jobs. Typically mpirun is used for MPI jobs. We assume your system has 16 processors on a core, but you can use any number.

For drawing postscript files, this tutorial assumes you are use the Apple open. Substitute your postscript viewer of choice for open.

Basic strategy

QSGW can make very high quality band structures. In principle it can be done for many atoms, but in practice QSGW is limited to relative small sizes (up to about 50 atoms, depending on the context).

By using the “near sightedness” principle and cutting and pasting self-energies in various regions, we can using small calculations to construct a high-quality potential for much larger systems.

Here we use QSGW to perform calculations on :

1. Bulk InAs in the zincblende structure
2. Bulk GaSb in the a sheared zincblende structure, tailored to match InAs
3. Short-period superlattice of (InAs)n(GaSb)n along [001]

Finally the three calculations are merged to construct an approximate self-energy for a longer period superlattice.

There are many steps required to accomplish this, but the result is a very high quality band structure for the superlattice.

This tutorial is somewhat involved, and you may wish to go through some basic tutorial such as the PbTe tutorial, the Bi2Te3 tutorial which describes the basis set, and the introductory QSGW tutorial before embarking on this one.

Getting Started

Three files are required to make any calculation. All are autogenerated once supplied structural information.

• A ctrl file, which is the main input file. We will autogenerate this file using the blm utility.
• A site file with structural data. (The ctrl file can also hold structural information, but we keep it separate so only the site changes as the structure changes.) blm also makes this file.
• A basp file, which defines the basis set. The free-atom code lmfa makes this file.
• Along [001], the zincblende lattice has a repeat length of two unit cells. It is simpler (but not essential) to use a doubled unit cell as the building block for the superlattice. We will call this the “conventional cell.” (1)

One more user-supplied file is needed, switches-for-lm. It contains algebraic variables, enabling control of parameters such as the number of k points without altering the ctrl file.

For lmfa to make parameters for the basp file, it needs information about each kind of atom and the augmentation sphere radii. We want to use a fixed basis set for all calculations, which means we need to set up a minimal superlattice which contains all the atoms before doing any electronic structure calculations. This will be (InAs)2(GaSb)2 [001] superlattice.

The prescription is then:

• Use blm to build site files for bulk InAs and GaSb in the zincblende structure.
• Use these files and the lmscell utility to assemble a site file for the (InAs)2(GaSb)2 [001] superlattice.
• Construct the ctrl and basp files from this structure. They can remain (more or less) fixed for different structures; Only files site and switches-for-lm change between structures.

Site files for zincblende InAs and GaSb

This init file contains the requisite structural information for bulk InAs and GaSb. Copy the box below into file init.dat.

LATTICE
% ifdef inas
% const a=11.43   # For InAs
% endif
% const a=11.50  # For GaSb
ALAT={a} PLAT=    0.0 0.5 0.5 0.5 0.0 0.5 0.5 0.5 0.0

% ifdef inas
SITE
ATOM=In    X= 0   0   0
ATOM=As    X= 1/4 1/4 1/4
% else
SITE
ATOM=Ga    X= 0   0   0
ATOM=Sb    X= 1/4 1/4 1/4
% endif

Notes:

• init.dat is designed to work with blm to make site files for both InAs and GaSb.
• Both form in the zincblende structure, but InAs and GaSb are slightly lattice-mismatched. We assume here that InAs is thick and that the global lattice constant is that of InAs. To make GaSb conform to that lattice, it must be sheared by elongating it on the [001] direction. We make the simple (and reasonably accurate) assumption that the amount of shear is fixed by keeping the volume constant.
• Reliable treatment of the In and As d states will require local orbitals. blm will automatically include them; however it will be necessary to set the tolerances for identifying them a bit tighter than the default; see --eloc below.

Make the site file for zincblende GaSb:

cp init.dat init.gasb
blm -vinas=f --wsitex,quad1 --findes,rmin=1.3,float,1spec --ctrl=ctrl gasb
mkdir gasb-zb
cp ctrl.gasb site.gasb gasb-zb

blm creates ctrl.inas and site.inas. Its switches a do the following:

 -vinas=f sets preprocessor variable inas (allowing init.inas to be used for both InAs and GaSb) –wsite,quad1 Translates all sites to be in the first quadrant –findes,rmin=1.3,float,1spec Identifies open interstitial regions in which to place floating orbitalsb –ctrl=ctrl Write output ctrl file directly to ctrl.gasb

a Many of these switches have no effect on the site files (which is all we are after in this context) but they are included here anyway since they will be used later.
b The lmf basis is slightly too small for converged GW calculations. This will be remedied when Jigsaw Puzzle orbitals become available.

Site files for conventional cells of InAs and GaSb

We define the conventional cell(1) to be a doubled unit cell along [001]. This will be the basic building block for (InAs)n(GaSb)m superlattices. The InAs and GaSb conventional cells are then the (InAs)2(GaSb)0 and The (InAs)0(GaSb)2 superlattices.

echo .5 .5 0 .5 -.5 0 0 0 1 | lmscell gasb --wsite~short~fn=site
lmscell ctrl.gasb -vz=11.5/11.43 '--stack~scale=1/z~stretch=z^3~sort@h3-z/1000~wsitex@short@fn=site'
mkdir inas0gasb2
cp site.gasb inas0gasb2/site.inasgasb

lmscell is first invoked in the supercell mode. This takes as input any desired primitive lattice vectors. (For the supercell to be well defined, these vectors should be integer multiples of the starting cell). In this case it generates a site file for the conventional cell.

lmscell is invoked a second time, invoking the superlattice editor. This editor is tailored to build supercells. In the present case, it does not stack any cells together but merely shears the lattice and orders the sites:

 -vz=11.5/11.43 declares variable z as the ratio of GaSb/InAs lattice constants, needed to determine the shear –stack~… invokes the superlattice mode –stack~scale=1/z scales the lattice constant by 1/z –stack~stretch=z^3 shears (stretches) the lattice along the third axisa by z3. –stack~sort@h3-z/1000 sorts the atoms in the supercellb. Note the @ is used as a delimiter separating subarguments belonging to ~sort. –stack~wsitex@short@fn=site Writes a site file to file site.gasb, overwriting the existing file.

a Taken together, the volume compression and the shear keep the cell volume constant.
b Sites are sorted by increasing height h3, which is the axis along the third lattice vector (the z axis in this case). The extra small term of -z/1000 is use to disambiguate sites of equal height: the physical site with z>0, takes precedence over the floating orbital site with z=0.

Make the corresponding site files for InAs. No shear is required in this case.

cp init.dat init.inas
blm -vinas=t --wsitex,quad1 --findes,rmin=1.3,float,1spec --ctrl=ctrl inas
mkdir inas-zb
cp ctrl.inas site.inas inas-zb

(The ctrl.ext file is retained only because it is needed to construct superlattices).

The conventional cell:

echo .5 .5 0 .5 -.5 0 0 0 1 | lmscell inas --wsite~short~fn=site
lmscell inas '--stack~sort@h3-z/1000~wsite@short@fn=site'
mkdir inas2gasb0
cp site.inas inas2gasb0/site.inasgasb

At the close of lmscell, the following table is printed:

Found 2 q points in original cell that map to q=0:
q (2*pi/a)                    q (multiples of qlat)
0.000000   0.000000   1.000000      0.500000   0.500000   0.000000
0.000000   0.000000   0.000000      0.000000   0.000000   0.000000

It tells you the other k point (an X point) that gets folded into the Γ point because of the artificial construction of the superlattice.

Site, ctrl, and basp files using the (2,2) superlattice as template

The following instructions build the ctrl file we will use for all calculations, and a site for the (InAs)2(GaSb)2 superlattice.

cp ./inas-zb/* .
cp inas2gasb0/site.inasgasb site.inas
cp inas0gasb2/site.inasgasb site.gasb
lmscell inas '--stack~show~file@site.gasb@dpos=0,0,-.032~show~wsite@short'
cp sites.inas sitein.inasgasb
blm --molstat --ewald=1d-10 --rdsite --omax=-.04 --addes,end --nosort --shorten=no --findes,rmin=1.3,float,1spec --ehmx=-.3 --nit=20 --gw~gcutb=2.7~gcutx=2.3~pbtol=3e-4,3e-4,1e-3,1e-3,1e-2 --nk=6,6,6 --nkgw=4,4,4 --gmax=9.2 --eloc=-3.5 --ctrl=ctrl inasgasb

This time lmscell stacks a superlattice of the InAs and GaSb conventional cells together.

 ~show Not required, but shows the structure at the current stage of construction ~file@site.gasb@dpos=0,0,-.032 stacks the GaSb conventional cell on top, rigidly shifting the cell by (0,0,−0.032) a ~wsite@short Writes a site file. The default name is sites

a The (InAs)m(GaSb)n superlattice has two interfaces. At one interface frontier atoms (In,Sb) are nearest neighbors, while (Ga,As) are neighbors at the other. Away from the interface the (InAs) and (GaSb) bulk part look like the bulk; however the entire (GaSb) block can shift relative to the (InAs), increasing one frontier bond and decreasing the other. To a first approximation to a full lattice relaxation, the entire (GaSb)n part of the lattice can be shifted, increasing one frontier bond length and decreasing the other. It turns out that residual forces are small if this constrained relaxation is allowed, and (0,0,−0.032) is the rigid shift that minimizes the energy. Relaxations are explained in this section.

blm is invoked with a somewhat complex combination of switches

 –molstat Adds a category for molecular statics (for doing lattice relaxations) –ewald=1d-10 Sets a stricter tolerance for Ewald sums. This is needed because long cells require it –rdsite Read structure information for a site file. By default it reads from file sitein.ext) –omax=-.04 Restrict the augmentation sphere overlap to −4 % of touching.a –addes,end Tells blm to put the E species last in the SPEC category (for aesthetics) –nosort Tells blm not to sort the sites (for aesthetics) –shorten=no Do not shorten the basis vectors –findes,rmin=1.3,float,1spec Identifies open interstitial regions in which to place floating orbitalsb –ehmx=-.3 Specifies the upper bound to the smooth Hankel energy –nit=20 maximum number of iterations (not required) –gw Tailor input file for a GW calculation –gw~gcutb=2.7~gcutx=2.3 Set G cutoffs for 1-particle and 2-particle quantities in the GW calculationc –gw~pbtol=3e-4,3e-4,1e-3,1e-3,1e-2 product basis tolerances in the GW calculationc –nk=6,6,6 Number of k divisions for lmf –nkgw=4,4,4 Number of k divisions for GWd –gmax=9.2 G cutoff for lmfe –eloc=-3.5 Set default energy cutoff when searching for local orbitalsf –ctrl=ctrl Write ctrl file to ctrl.inasgasb

a This is not necessary for the bulk, but gives some “wiggle room” to avoid spheres overlapping at an interface. It turns out it is not necessary for this interface either, but you can verify that the band structures depend little on the choice of overlap.
b The lmf basis is slightly too small for converged GW calculations. This will be remedied when Jigsaw Puzzle orbitals become available.
Meanwhile we augment the basis by adding smooth Hankels without augmentation spheres (“floating orbitals”).
All floating orbitals are constrained to be the same species (the modifiers to --findes are described in the command-line switches page).
c You do not need to supply these cutoffs but the default values are a bit more conservative. These slightly coarser values speed up the calculation without any significant change the result.
d The GW mesh of 4×4×4 divisions slightly overestimates the gap (by 0.05 eV). We keep this mesh for economy.
In any case the QSGW gap is generally too large, and we will scale the QSGW $\Sigma^0$ to reduce the gap.
e gmax is not known in advance; however lmfa will print out a recommended value.
f With the default cutoff, lmfa does not pick up the deep As 3d state, which reduces the gap by 0.05 eV.

Construct the basis set file, basp.inasgasb

lmfa -vles=1 --basfile=basp ctrl.inasgasb
lmfa ctrl.inasgasb
echo -vnk2=3 -vnkgw2=2 > switches-for-lm
mkdir -p inas2gasb2/relax0
cp ctrl.inasgasb basp.inasgasb site.inasgasb switches-for-lm inas2gasb2/relax0

lmfa creates the basis file basp.inas and the atomic density file atm.ext. It must be run twice, initially to determine what local orbitals are to be included. A second pass is necessary to properly partition the valence-core density. Arguments do the following:

 -vles=1 variable causes E “sites” to have an augmentation radius, like at “atom”: with atomic number 0a –basfile=basp write the auto-generated basis to basp.inas (default = basp0.inas)

a This temporarily turns the floating orbitals into “empty spheres” with augmentation radii. Empty spheres are not used in the electronic structure calculations, but the radius is used to make the basis set.

The essential files needed for the superlattice are saved in inas2gasb2/relax0.

Inspect basp.inas. Note that In 4d is the normally augmented valence state; the 5d is treated as a local orbital. The corresponding As 3d state is deeper than the In 4d, so for As, the 3d is treated as a local orbital, and the 4d is the normally augmented valence state.

Copy necessary files to complete the set of required input files for basic building blocks

echo -vnk2=4 -vnkgw2=4 > switches-for-lm
cp ctrl.inasgasb basp.inasgasb switches-for-lm inas2gasb0
cp ctrl.inasgasb basp.inasgasb switches-for-lm inas0gasb2

The rigid shift (0,0,-.032). optimized for the (2,2) superlattice is nearly optimal for any (InAs)n(GaSb)m superlattice.

Make various (InAs)n(GaSb)n superlattices (with and without shifts)

lmscell inas '--stack~show~file site.inas dup=1~file@site.gasb@dpos=0,0,-.032@dup=2~show~wsite@short'
echo -vnk2=2 -vnkgw2=1 > switches-for-lm
cp sites.inas site.inasgasb
mkdir -p inas4gasb4/relax0
cp ctrl.inasgasb basp.inasgasb site.inasgasb switches-for-lm inas4gasb4/relax0
lmscell inas '--stack~show~file site.inas dup=2~file@site.gasb@dpos=0,0,-.032@dup=3~show~wsite@short'
echo -vnk2=2 -vnkgw2=1 > switches-for-lm
cp sites.inas site.inasgasb
mkdir -p inas6gasb6/relax0
cp ctrl.inasgasb basp.inasgasb site.inasgasb switches-for-lm inas6gasb6/relax0

LDA + QSGW self-consistency, bulk systems (conventional cells)

Make a self-consistent density and potential in the LDA:

cp inas2gasb0/* .
lmfa ctrl.inasgasb
mpix -n 16 lmf ctrl.inasgasb

Note that at the LDA level, there is no bandgap. You can see this in these lines from the standard output:

BZWTS : --- Tetrahedron Integration ---
... only filled or empty bands encountered:  ev=0.100216  ec=0.100216
VBmax = 0.100216  CBmin = 0.100216  gap = 0.000000 Ry = 0.00000 eV
BZINTS: Fermi energy:      0.100216;  56.000000 electrons;  D(Ef):    0.000
Sum occ. bands:  -72.6386525  incl. Bloechl correction:    0.000000

The following makes InAs self-consistent at the QSGW level and saves essential files in inas2gasb0

lmfgwd --job=-1 ctrl.inasgasb cat switches-for-lm
lmgwsc --wt --insul=28 --mpi=16,16 --tol=2e-5 --sym --getsigp --save= inas2gasb0 inasgasb > out.job
mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm
mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm --rs=1,0 -vnit=1 >  out.lmf
cp rst.inasgasb  out.lmf inas2gasb0

You should find lmgwsc converges to a high degree of self-consistency in 4 iterations. The first invocation of lmf isn’t necessary. It touches up the restart file making a self-consistent density for the given QSGW Σ0.

Before drawing energy bands, check the bandgap:

grep gap out.lmf

You should see that the bandgap is 0.76 eV, significantly larger than the experimental value of 0.42 eV. The error can be traced to several sources:

1. incomplete k convergence
2. neglect of spin-orbit coupling
3. neglect of ladder diagrams in the RPA dielectric function
4. neglect of phonon contributions to the self-energy

All of them are relatively small; but each reduces the gap by 0.05 eV to 0.1 eV.

Energy bands in conventional cells with spin-orbit coupling

Spin-orbit coupling plays a non-negligible role. It is readily taken into account because  ctrl.inas is already set up for a calculation with SO coupling.

However if you run lmf straight off it will fail with the message

Exit -1 rdsigm unexpected value 21 for file sigm nqp ... expected 64

This is because the QSGW self-energy was constructed on the irreducible wedge of the zincblende lattice; however spin-orbit coupling reduces the symmetry. The simplest solution (though not most efficient) workaround is to remake the QSGW self energy file Σ0 at all k points

lmf ctrl.inasgasb --wsig:fbz
cp sigm2.inasgasb sigm.inasgasb

The following steps will:

• find the Fermi level with SO coupling included (only needed to tell the band plotter knows where to set the energy zero)
• generate the energy bands
• make a postscript file

We use a symmetry lines file for the fcc lattice from the Questaal distribution. since we have an artificial supercell, use only the first two symmetry lines :

mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1
mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1 --band~fn=syml
echo -7,6 / | plbnds -nocol -fplot~sh~scl=.7 -ef=0 -scl=13.6 -lbl=L,G,X bnds.inasgasb
open fplot.ps

lmf is invoked initially to establish the Fermi level. This is not essential but it is convenient when drawing the energy bands.

Inspect the first row of bnd2.dat (the Γ point is the first point in the second panel of the energy bands). You should find the following:

 k-point split off band light hole band heavy hole band conduction band 0.00000 -0.3618   -0.3618 -0.0299   -0.0299 0.0000   0.0001 0.6623   0.6623 default cutoffs

The split-off band, at −0.36 eV, is split slightly less than the experimental splitting (−0.43 eV). The light-hole and heavy-hold bands should be degenerate; there is a slight splitting because the cubic symmetry was artificially broken(1) (the c axis is not equivalent, e.g. the k mesh is different). The bandgap, 0.66 eV, still significantly larger than the experimental value of 0.42 eV. Band gaps are overestimated, and indeed they should be because the RPA screened coulomb interaction W is not sufficiently screened. And also zero-point motion has been left out. Improvement of W by the addition of ladder diagrams indeed does improve the gaps, as Kutepov has shown. However, these higher order diagrams are expensive, and we adopt a simpler approach here. The dielectric constant is about 80% of experiment (it seems always to be the case for semiconductors and insulators), because ladder diagrams are missing in the RPA. We have found that scaling ε by 0.8 (assuming the full $\epsilon(\mathbf{r},\mathbf{r'},\omega)$ gets scaled by a fixed number), or alternatively using a hybrid of 80% QSGW and 20% LDA we can mimic the effect of the ladders. This eliminates most of the error.

mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1 --mixsig=0.8
mpix -n 16 lmf ctrl.inasgasb -vso=1 cat switches-for-lm --rs=1,0 -vnit=1 --band~fn=syml --mixsig=0.8
echo -7,6 / | plbnds -nocol -fplot~sh~scl=.7 -ef=0 -scl=13.6 -lbl=L,G,X bnds.inasgasb

Click here to view the band structure of InAs and GaSb comparing QSGW with hybrid QSGW, $(\Sigma^0{\times}0.8{+}V_\mathrm{xc}^\mathrm{LDA}{\times}0.2)$. With the scaling, agreement with ellipsometry measurements is excellent (shown as circles). The fundamental gap drops to 0.42 eV in InAs.

GaSb proceeds in much the same way as InAs. Repeat the preceding LDA+QSGW instructions, substituting inas2gasb0inas0gasb2. The direct gap is 1.21 eV, larger than the experimental value 0.82 eV. With SO coupling and scaling it drops to 0.78 eV. Note that to get more perfect agreement with experiment, the scaling factor 0.8 can be modified slightly.

The InAs/GaSb [001] interface

The procedure for a superlattice is nearly identical to the bulk. We have already generated the necessary input files.

Set shell variable dir to inas2gasb2/relax0 or inas4gasb4/relax0, or another directory where a site file resides. The inas2gasb2 case is useful to construct the ctrl file and to determine the forces at the interface but the system is to small to be of relevance for electronic structure.

cp $dir/* . lmfa ctrl.inasgasb mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm You should find that the forces are all small. Forces on the frontier atoms will be modest but still less that 10 mRy/au. Lattice relaxations are discussed later. For the QSGW calculation, do the following. The number of electrons (set here through shell variable nel) isn’t essential. lmfgwd --job=-1 ctrl.inasgasb cat switches-for-lm lmgwsc --wt --insul=$nel --mpi=16,16 --tol=2e-5 --sym --getsigp --save= $dir inasgasb > out.job mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm; mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm --rs=1,0 -vnit=1 > out.lmf cp rst.inasgasb out.lmf$dir

Self-consistency should be reached in 4 iterations. (InAs)2(GaSb)2 takes about 30 minutes/iteration; (InAs)4(GaSb)4 about 12 hours/iteration; (InAs)6(GaSb)6 about two days.

Site-resolved Partial and Integrated Densities of States

If you are not familiar with Questaal’s implementation of Densities-of-States (DOS) resolved by site, it is recommended that you review the tutorial on partial DOS in Cr3Si6 before going through this section.

By projecting DOS onto sites, we can resolve the superlattice DOS layer by layer and compare to bulk InAs and GaSb. This is useful in determining the valence band offset. Here we set up the DOS for the reference bulk systems.

The instructions below refer to bulk InAs, bulk GaSb, and the (InAs)4(GaSb)4 superlattice. The quantities in the table below distinguish which system the instructions apply to

 Structure/variable directory k-divisions in-plane k-divisions normal site list dos channels dir nel nk1 nk2 sitelist doschan InAs inas2gasb0 28 12 9 1:8:2 1:4 GaSb inas0gasb2 28 12 9 1:8:2 1:4 (2,2) SL inas2gasb2/relax0 56 12 4 1:16:2 1:4;seq=5:16:4 (4,4) SL inas4gasb4/relax0 112 12 2 1:32:2 1:4;seq=5:32:4 (6,6) SL inas6gasb6/relax0 168 12 3 1:48:2 1:4;seq=5:48:4

Assign shell variables to each of these quantities, e.g. in bash, dir=inas2gasb0 nk1=12 nk2=9 sitelist=1:8:2 doschan="1:4" for the InAs case.

Variable sitelist is used in the argument --pdos~sites in lmf. It is a an integer list, e.g. 1:31:2 → 1,3,5,…31.

Variable doschan is used in the argument -lst= in pldos. It is a generalized integer list; the argument to seq= is another list, which creates a new panel for each element in the list. In particular, seq=1:4;seq=5:32:4 makes eight panels consisting of the following channels
1:4   5:8   9:12   13:16   17:20   21:24   25:28   29:32
while seq=1:4;seq=5:48:4 makes twelve panels consisting of the following channels
1:4   5:8   9:12   13:16   17:20   21:24   25:28   29:32   33:36   37:40   41:44   45:48

These correspond to the individual (InAs) or (GaSb) layers when making partial DOS in the superlattice

To make the DOS and iDOS, a moments file must be created.

cp $dir/* . ; rm -f mixm.inasgasb ln -s sigm sigm.inasgasb mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm --rs=1,0 -vnit=1 The last step is not needed, but run it anyway to confirm that the RMS DQ is small and that the bandgap is as it should be (0.76 eV in InAs). It is important to include spin-orbit coupling, so we need to render the QSGW Σ0 suitable for reduced symmetry operations. Here we use the 80% hybrid to well reproduce the experimental band structure. Do: lmf ctrl.inasgasb --wsig:fbz; cp sigm2.inasgasb sigm mpix -n 16 lmf -vso=1 ctrl.inasgasb -vnk1=$nk1 -vnk2=$nk2 --rs=1,0 -vnit=1 --rsig:fbz --pdos~mode=0~sites=$sitelist --mixsig=.80 --pr35,40 > out.pdos
cp out.pdos dir Note the large number of k-divisions used. This is not necessary but a fine k mesh generates a smooth DOS. The --pdos switch is explained in some detail in the partial densities of states tutorial. See also documentation for the command-line switch for –pdos. lmf generates file moms.inasgasb, which contains the weights necessary to generate site-resolved DOS. In what follows, DOS are generated in two forms: first to extract the partial DOS near the band edges, for drawing local DOS, and also the integrated DOS (iDOS) to extract the position of the deep d states, which will be used align reference levels in the bulk with the superlattice, to extract the valence band offset. Partial DOS resolved by site Partial DOS are made using the lmdos utility, once lmf has generated the appropriate file. This utility generates the projection of DOS into augmentation spheres2, given a moms file generated by a band program, as documented in the partial densities of states tutorial, and also documentation for the --dos command-line argument. You must invoke lmdos with the same switches that were supplied to lmf, and also an additional argument, --dos. lmdos -vso=1 ctrl.inasgasb -vnk1=nk1 -vnk2=$nk2 --rs=1,0 -vnit=1 --rsig:fbz --pdos~mode=0~sites=$sitelist --mixsig=.80 --dos:npts=2001:window=-1,1

On exit, lmdos should print the following table (inas2gasb0 case)

Channels in dos file generated by LMDOS:
site class label   spin-1                       spin-2
1    1   In     1                            2
3    2   As     3                            4
5    4   In2    5                            6
7    5   As2    7                            8

As the table shows, DOS are resolved for eight channels, interleaving (spin-up,spin-down). This system is nonmagnetic so the two are equal and there is no point in distinguishing them. Also the two In and two As should be equivalent. lmdos writes the DOS to file dos.inasgasb in a (not-particularly friendly) format. pldos reads files in this format to make figures, as documented in the pldos manual.

There are thus four channels per layer. Since we don’t need to resolve DOS by site yet, make a figure of the In+As DOS for a single layer by combining DOS for a single layer, i.e. channels 1:4. For the bulk InAs and GaSb cases, shell variable doschan should be set to “1:4”. Make a postscript file for the DOS with (instructions also apply to the superlattice, as described below)

echo 0.999 4 -3 2 | pldos -fplot~sh~long~open~tmy=.25~dmin=0.40~xl=E -esclxy=13.6 -ef=0 -lst="$doschan" dos.inasgasb cp dosp.dat$dir/pdos.80
open fplot.ps

fplot~... tells pldos to create a DOS data file dosp.dat, formatted in Standard Questaal format for xy data. The second command saves this file for future use, to draw reference and superlattice DOS on the same figure. pldos also makes a script plot.dos for the fplot utility, which will in turn generate a postscript file, fplot.ps. (You can edit plot.dos to customize the figure as we will do when aligning superlattice and reference DOS). Modifiers ~long~open~tmy=.125~dmin=0.40~xl=E make the figure more aesthetic, while ~sh tells pldos to run fplot in a subshell, making fplot.ps directly. -esclxy=13.6 converts ordinate and abscissa from atomic units to eV and eV−1, respectively, and the abscissa is shifted to place the Fermi level at zero (ef=0).

The arguments from standard input (0.999 5 -3 2) cap the DOS to just under 1, with a frame of about 5 cm, in an energy window (−3,2)  eV.

Superlattice DOS

For the (InAs)4(GaSb)4 or (InAs)6(GaSb)6 superlattice, follow essentially the same steps from the beginning of the pdos section, and it is instructive to resolve the DOS by layer. Now there are a multiplicity of layers (8 or 12) and 32 or 48 channels, four channels per layer; pldos will draw 8 or 12 panels of DOS, one panel for each layer. doschan is constructed to group the DOS for each layer into one panel (see sitelist and doschan shell variables explained).

In+As DOS in InAs. Valence band maximum shown by the blue dashed line.

Ga+Sb DOS in strained GaSb

Layer resolved DOS in (InAs)6(GaSb)6 superlattice, with bulk DOS shifted by local electrostatic shifts shown as dashed lines.

Integrated DOS to extract position of core levels

Invoke lmdos as before, but change the --dos switch to --dos:npts=40001:window=-3,1:rdm:idos. This generates the integrated DOS on a very fine energy mesh over a wide energy window, writing data in Standard Questaal format for xy data.

lmdos -vso=1 ctrl.inasgasb -vnk1=$nk1 -vnk2=$nk2 --rs=1,0 -vnit=1 --rsig:fbz --pdos~mode=0~sites=$sitelist --mixsig=.80 --dos:npts=40001:window=-3,1:rdm:idos cp dos.inasgasb$dir/idos.80

Inspect dos.inasgasb. Energy is stored in the first column and the  iDOS  starting in column 2. The second instruction saves the  iDOS  for later, when we align bands to extract the valence band maximum.

For the bulk cases, e.g. InAs, draw a picture of the integrated DOS over a broad energy window (−3,1) Ry:

fplot -frme 0,.7,0,.4 -x -3,0 -y 0,7 -tmx .1:5 -colsy 1+1:1+4 -qr dos.inasgasb
open fplot.ps

Four  iDOS  are drawn but spin-1 and spin-2 overlay each other. In one curve, at −2.8 Ry there are two closely spaced step-function changes in the  iDOS  to 2, and then 5 electrons. These are the five 5 As d states split into 2 + 3 levels by spin-orbit coupling. The other curve shows a similar shape near −1.1 . Above −0.8 Ry where the valence bands begin, the  iDOS  begin to rise again.

For the superlattices there are many levels. To see the evolution of In d levels in a (−1.2,−1) Ry energy window, try

fplot -frme 0,.7,0,.4 -x -1.2,-1 -y 0,5 -tmx .05 -1p -map -i 'x1>-1.2&x1<-1' 1:nr  -colsy 4:24:4 -lt 1,col=0,0,0 -qr dos.inasgasb
open fplot.ps

For the strained bulk GaSb or the superlattice, repeat the steps starting at Site-resolved Partial Densities of States, substituting the appropriate set of shell variables.

In+As iDOS in InAs

Ga+Sb iDOS in strained GaSb

In 4d iDOS in the (InAs)6(GaSb)6 superlattice on a finer energy scale

Valence band offset

Reference points needed for superlattice band offsets

Interfacial effects — apart from possible long range electrostatic shifts — die out quickly. Then the local valence band maximum (l-VBM), sufficiently far from an interface, will not depend on the interface. A good way to determine the alignment of the l-VBM on two sides of an interface, is to choose a reference energy or potential whose position relative to the l-VBM is constant far enough away from the interface; call it Eref. Two good choices for Eref are:

• the position of semicore states such as 3d or Sb 4d
• the potential on the augmentation sphere boundary

This construction enables us to determine a local valence band maximum (l-VBM) for any layer, but because of nonlocal quantum effects, its meaning is not precise close to the interface.

It is common practice to align bands from separate bulk calculations. This is how “natural band offsets” are determined, for example. Aligning bands for two independent materials relies on an assumption, namely that if an interface were formed between two materials, no change in the dipole contribution to the electrostatic potential would arise with respect superposition of the potentials from the independent materials. The “natural band offsets” scheme aligns bands from independent calculations of two surface structures. Other assumptions can be made, for example from the following. Imagine a lattice of spheres that fill space. If each sphere is neutral, a sensible definition of the electrostatic potential is to choose 0 at the sphere boundary. As a next step, transfer charge between spheres. The electrostatic potential (the Madelung potential) is now completely defined, since the reference was defined. In any case the assumption can be tested by comparing how the potential differs from the potential of an actual interface. Since in fact the offset depends on details of the relaxation of the interface, any assumption is likely to be suspect.

Absent additional dipoles forming at a junction relative to the bulk potentials as defined by the implementation, valence band maxima are merely those of the bulk systems.

This information is contained in the file out.pdos, of the bulk directory, e.g. for InAs (directory inas2gasb0) there should be this table:

BZWTS : --- Tetrahedron Integration ---
... only filled or empty bands encountered:  ev=0.120480  ec=0.151794
VBmax = 0.120480  CBmin = 0.151794  gap = 0.031314 Ry = 0.42587 eV

Extract the VBM for the two bulk systems:

grep -A 2 'VBmax' inas2gasb0/out.pdos | head -1 | awk '{print $3}' > vb1 grep -A 2 'VBmax' inas0gasb2/out.pdos | head -1 | awk '{print$3}' > vb2

If no additional dipole forms at the interface the band offset would be (0.120480-0.164545)×13.6 = −0.599 eV. That is, the VBM of InAs falls below that of GaSb by about 0.6 eV.

We need to add to this the dipole shift, which can be obtained by taking the difference in reference energies.

Dipole shifts by aligning d states in the superlattice to the bulk

We align the d states (16 of them in the (4,4) superlattice), against the corresponding reference states. Thus we need

1. Energies of the bulk InAs + bulk GaSb d states
2. Energies of the d states in the superlattice

For the former, figures were already drawn identifying the levels pictorially. There is some arbitrariness in how a numerical value is extracted, but the arbitrariness is not large except for a possible uniform constant which is irrelevant. For convenience let’s pick the energy where the electron count is exactly 4 for each lattice. This readily accomplished with the mcx calculator. For the bulk systems

mcx -qr inas2gasb0/idos.80 -at 4 x2  -w:nohead . > d1
mcx -qr inas2gasb0/idos.80 -at 4 x4  -w:nohead . >> d1
mcx -qr inas0gasb2/idos.80 -at 4 x2  -w:nohead . > d2
mcx -qr inas0gasb2/idos.80 -at 4 x4  -w:nohead . >> d2

-at val expr  finds the abscissa where the ordinate matches  expr. Note that the first and third channels (In spin-1 and As spin-1) lie in columns 2 and 4. Inspect files d1 and d2. You should get about −1.038 Ry for the In level and −2.673 Ry for the As level.

The d levels of a superlattice are a bit more complicated as there are many of them.
They are readily extracted in mcx using the looping construct [ ... ]. Set shell variable dir according to the table and invoke mcx for the (4,4) superlattice as follows:

mcx -qr dir/idos.80 -a idos [ ix=1:31:2 idos -at 4 x\{ix+1} -tog -pop '?ix>1' -rcat ] > d3 For the (6,6) superlattice substitute 1:31:2 → 1:47:2. mcx extracts the crossing point of 4 electrons for columns 2,4,6,…32 and concatenates them into a 16×2 array. 1:31:2 is the integer list 1,3,5,…31. Note however that these channels correspond to columns 1+1,1+3,1+5,…1+31 in file idos.80. Before entering the loop […] mcx loads idos.80 onto the stack and names it array idos. Within the loop arguments are expanded into idos -at 4 x2 -tog -pop idos -at 4 x4 -tog -pop -rcat .... Note that -rcat is executed only after the iteration because of the conditional ?ix>1. For each column xn in the list, mcx finds the abscissa where column xn crosses 4 and pushes on the stack a 1×2 array with the first column being the abscissa. idos, initially loaded onto the stack, is still present; thus -tog -pop toggles and pops the stack leaving only the desired result in as 1×2 array. In the second and subsequent iterations -rcat concatenates the top two arrays, thus eventually accumulating a 16×2 array. Inspect file d3. You should see the 16 d levels of the superlattice. The difference in the bulk and superlattice data is the local shift in VBM. Compute this difference and convert to eV. For the (4,4) superlattice: mcx [ ix=1:4 d1 '?ix>1' -rcat ] [ ix=1:4 d2 '?ix>1' -rcat ] -rcat d3 -ccat -e2 x3 x1 -e1 x1-x2 -s13.6 For the (6,6) superlattice, substitute 1:4 → 1:6. Inspect the output of this command. The (4,4) superlattice still has strong interface effects, so we discuss the (6,6) case instead. Around row 6 the shift is about 0.05 eV; this corresponds to the middle (InAs) layer. Around row 18 it is 0.01 eV. Note that all the elements can be shifted by a uniform constant, which arises from arbitrariness in the choice of the G=0 part of the potential. Thus net dipole shift is about 0.04 eV. Extracting electrostatic potentials from standard output An alternate way to obtain the band offset is to align bulk and superlattice electrostatic potentials on the augmentation sphere boundary. Apart from interfacial effects, it should yield the same result as the d level shifts. lmf calculates the potential at each augmentation sphere boundary and sets the average to zero. This information is also contained in out.pdos; for example for InAs there should be a table Average es pot at MT boundaries after vconst shift Site ves Site ves Site ves Site ves 1 -0.033546 | 3 0.036403 | 5 -0.033548 | 7 0.036403 2 0.714034 | 4 0.662523 | 6 0.714049 | 8 0.662523 Sites 2,4,6,8 are floating orbital sites which we will ignore for the purpose of aligning bulk and superlattice states. For the (InAs)4(GaSb)4 and (InAs)6(GaSb)6 superlattices there are 32 or 48 potentials, of which the odd numbered sites correspond to the real atoms. The following instructions extract these potentials from out.pdos for the two bulk systems and the (InAs)4(GaSb)4 superlattice, saving results in files v1, v2, v3. In each of the instructions grep, awk and sed are used to extract the table and render it into a form of xy data that mcx can read. Then mcx extracts the relevant potentials. grep -A 2 'Average es pot at MT boundaries' inas2gasb0/out.pdos | awk '1 ~ /^[0-9]+$/{print }' | sed 's/|//g' | mcx -nc=2 . -sort x1 -inc 'i%4==1|i%4==2' > v1 grep -A 2 'Average es pot at MT boundaries' inas0gasb2/out.pdos | awk '$1 ~ /^[0-9]+$/{print }' | sed 's/|//g' | mcx -nc=2 . -sort x1 -inc 'i%4==1|i%4==2' > v2 ... for the (4,4) superlattice: grep -A 9 'Average es pot at MT boundaries' inas4gasb4/relax0/out.pdos | awk '$1 ~ /^[0-9]+$/{print }' | sed 's/|//g' | mcx -nc=2 . -sort x1 -inc 'i%2==1' > v3 ... for the (6,6) superlattice: grep -A 13 'Average es pot at MT boundaries' inas6gasb6/relax0/out.pdos | awk '$1 ~ /^[0-9]+$/{print }' | sed 's/|//g' | mcx -nc=2 . -sort x1 -inc 'i%2==1' > v3 Aligning superlattice levels to reference We now have files d1d3… for the d levels and files v1v3… for the electrostatic potentials. It is convenient to tabulate the levels as a function of height from the basal plane. This is easily accomplished with the lmplan utility. For one of the superlattices, e.g. the (6,6) superlattice, set shell variables according to the table and do cp$dir/* .
echo q | lmplan ctrl.inasgasb

lmplan prints out a table of positions of each site its projection h on the axis perpendicular to the basal plane:

Gtplan: 24 different planes normal to (    0.000000    0.000000    1.000000)
Projection of normal onto P1, P2, P3:      0.000000    0.000000    6.055456
ib      Class   Plane  PL       x-x.h     y-y.h     z-z.h      h       del h

1       1:In      1    0      0.00000   0.00000   0.00000   0.00000   0.28662
2       5:E       1    0      0.50000   0.00000   0.00000   0.00000   0.00000
3       2:As      2    0      0.25000   0.25000   0.00000   0.25000   0.25000
4      10:E2      2    0      0.25000  -0.25000   0.00000   0.25000   0.00000
5       6:In2     3    0      0.50000   0.00000   0.00000   0.50000   0.25000
6      15:E3      3    0      0.00000   0.00000   0.00000   0.50000   0.00000
7       7:As2     4    0      0.25000  -0.25000   0.00000   0.75000   0.25000
8      20:E4      4    0      0.25000   0.25000   0.00000   0.75000   0.00000

Extract h from this table and store in file h, excluding floating sites:

echo q | lmplan ctrl.inasgasb | grep -A 100 Plane | grep -v :E | awk '$1 ~ /^[0-9]+$/{print $8,"#",$2}' > h

Make a single file with the superlattice and bulk electrostatic potentials, and another with the corresponding d levels. Note the conversion to electron volts. Also plot these changes as a function of h, solid lines for the electrostatic potential and dashed for the d levels. The following applies for the (4,4) case; for the (6,6) superlattice, substitute 1:4 → 1:6.

mcx v3 [ ix=1:4 v1 '?ix>1' -rcat ] [ ix=1:4 v2 '?ix>1' -rcat ] -rcat -ccat -e2 x2 x4 -s13.6 -w:nohead . | pr -m -t - h > estat
mcx d3 [ ix=1:4 d1 '?ix>1' -rcat ] [ ix=1:4 d2 '?ix>1' -rcat ] -rcat -ccat -e2 x1 x3 -s13.6 -w:nohead . | pr -m -t - h > dlevel
fplot -frme 0,.7,0,.7 -ab x3 -ord x1-x2 estat -ab x3 -ord x1-x2 -lt 2 dlevel

Shifts in d levels (dashed line) and potentials (solid line) at the augmentation boundary in the (InAs)6(GaSb)6 superlattice, relative to the bulk, as a function of height h from the basal plane.

In the figure, the vertical dot-dashed lines show the center point of the (InAs) and (GaSb) layers, respectively. Note the following points:

1. The two references are slightly inconsistent, especially at the boundary layers. These are consequences of finite-size effects.
Even where the potential is most bulk-like there is still a nearly constant 0.02 eV discrepancy.
2. A global constant can be added to either of these shifts, since they are determined by the arbitrary definition of the G=0 potential.
Taking this into account, the relevant dipole shift is the change at the midpoints (compare blue and red lines). The shift is approximately 0.05 eV.
3. The potential shift is not quite flat even in the middle regions. This is not an artifact of the small cell: if a still longer cell were calculated it would become clear that there is a linear term, which if many layers are present, would cause a dielectric breakdown. This is to be expected because the [001] interface is polar. In the real world the crystal will reconstruct somehow, possibly by introducing charged defects at the interface.
4. The dipole adds to the “natural” (or reference) valence band offset which is −0.599 eV. Thus the band offset is 0.55 eV.
Interpretation of (6,6) superlattice

With this information in hand let’s revisit the partial DOS of the (6,6) superlattice and observe the following

1. DOS in the frontier layers (InAs)1, (InAs)6, (GaSb)1 (GaSb)6 differs significantly from the bulk counterpart. GaSb valence bands just near the valence band maximum tunnel into the InAs (see (InAs)1 near the VBM); similarly InAs conduction band states near the conduction band minimum tunnel into the GaSb (see (GaSb)1 near the CBM). These states decay exponentially away from the interface.
2. Even while the band offset exceeds the bulk InAs bandgap (0.43 eV), the gap does not close, even in the (6,6) superlattice. Look for the gap in inas6gasb6/relax0/out.pdos: it is 0.28 eV. This appears to originate from confinement effects: the GaSb acts like an effective quantum well confining the (highly extended) states in the InAs conduction band.

To more clearly see this, make and inspect the energy bands along the c axis around a (−1,1) eV window around the Fermi energy

echo   11   0  0  0     0 0 .5  > symlz.inasgasb
mpix -n 11 lmf -vso=1 ctrl.inasgasb -vnk1=12 -vnk2=2 --rs=1,0 -vnit=1 --rsig:fbz --mixsig=.80 --band~mq~col=1:30,dup=936,seq=79:391:78~col2=40:69,dup=936,seq=118:430:78~fn=symlz
echo -1,1,5,10 | plbnds -fplot~scl=.7~sh -ef=0 -scl=13.6 -dat=green -lt=1,bold=3,col=0,0,1,colw=1,0,0,colw2=0,1,0 -lbl=G,X  bnds.inasgasb

Bands colored blue are of GaSb character. There is a nearly dispersionless state at the VBM (E=0). To the extent a band is red, it is of InAs character. The lowest conduction band has significant, but by no means all InAs character: InAs has a very small effective mass, and this state it is higher than it would be in a longer period superlattice because of confinement effects. The extended range of this state is also evident by its significant dispersion. States colored green are of As character. The highest valence-band state of nearly pure As character lies at around −0.8 eV. Those between −0.8 eV and −0.1 eV are of mixed GaSb/As character, with increasing GaSb approaching the E=0.

Atomic Relaxation of Superlattices

This section explains how a limited relaxation of the (2,2) superlattice is accomplished. As noted when it was constructed (see lmscell and blm arguments explained The (2,2) superlattice was used as a template to construct ctrl.inasgasb, and a the (GaSb) layers were stacked on with a particular offset (0,0,−0.032) that minimized the total energy.

Here we explain how that shift was determined by parameterizing the shift as (0,0,delta) and examining the total energy and forces as a function of delta

Set shell variable delta to (0,0,−0.032), e.g. in the bash shell, delta=-.032, and revisit the generation of the (InAs)2(GaSb)2 site file

cp ./inas-zb/* .
cp inas2gasb0/site.inasgasb site.inas
cp inas0gasb2/site.inasgasb site.gasb
lmscell inas "--stack~show~file@site.gasb@dpos=0,0,$delta~show~wsite@short" cp inas2gasb2/relax0/* . cp sites.inas site.inasgasb You should find, for delta=−0.032, that site.inasgasb and inas2gasb2/relax0/site.inasgasb are identical, but change with delta. Repeat the instructions above for delta=0 and make the calculation self-consistent rm -f mixm.inasgasb rst.inasgasb lmfa ctrl.inasgasb mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm -vdelta=$delta

Forces are printed in a table from stdout:

Forces, with eigenvalue correction
ib           estatic                  eigval                    total
1    0.00    0.00  -13.34     0.00    0.00   60.85     0.00    0.00   47.51
...
7    0.00    0.00   18.44     0.00    0.00   17.52     0.00    0.00   35.96
9    0.00    0.00   26.55     0.00    0.00  -56.90     0.00    0.00  -30.34
...
15    0.00    0.00   -7.94     0.00    0.00  -45.71     0.00    0.00  -53.65

Forces are large and opposite at the (1,15) interface, and similarly so at the (7,9) interface By rigidly translating the GaSb block all of the forces should drop.

Repeat the setup and relaxation for delta=−0.01 and the maximum force should drop from 53.7 to 38.9. Also you should find that the total energy drops by about 9 mRy. Repeat the steps for the sequence delta=0, −0.01, −0.02, −0.03, −0.04 and you should see the following emerge

 delta Mean abs force Maximum force Total energy 0 22.15 53.65 -66205.7300418 -0.01 15.56 38.89 -66205.7386342 -0.02 9.20 25.00 -66205.7442142 -0.03 3.51 12.08 -66205.7468116 -0.032 2.98 9.62 -66205.7469775 -0.034 2.48 7.20 -66205.7470269 -0.04 3.84 11.52 -66205.7464823

The mean absolute value of the force is minimum around delta=−0.034, as is the minimum in total energy. For the SL calculations we selected delta=−0.032 — not quite optimum because a slightly cruder minimization process was used. At the relaxed point the maximum force is 9.6 (delta=−0.032) and 7.2 (delta=−0.034).

The lattice should properly be fully relaxed. This is readily accomplished through the DYN category. Adding -vminx=t includes this category in the input:

mpix -n 16 lmf ctrl.inasgasb cat switches-for-lm -vminx=t > out &

but the forces are already small at the first stage and the additional relaxation hardly affects the band structure.

… To be finished

Footnotes and references

1 The zincblende lattice doubled along the z is a conventional cell in the sense that all primitive lattice vectors are orthogonal. It differs from the proper conventional for the zincblende lattice (a cube, with 8 atoms), and is not a true conventional cell because it reduces the symmetry from cubic to tetragonal symmetry.

2 The In+As DOS do not sum to the total DOS because there is a missing interstitial contribution. This is unimportant for the present purposes.