# Converting a periodic lattice into a trilayer geometry (ASA)

This tutorial shows how render a self-consistent calculation for crystal lattice into a trilayer geometry suitable for Landauer-Buttiker transport through an active region sanwiched between semi-infinite leads. This tutorial works with the cubic compound CsPbI_{3}, but lattice and basis vectors are deformed slightly to simulate the effect of finite temperature.

As a first step a self-consistent calculation is performed with the band code **lm**, and also with the crystal Green’s function code **lmgf**, in a periodic structure containing one formula unit. It is shown that the two codes generate the same densities.

Next, the single formula unit is doubled, using the superlattice editor. The system is still treated as periodic, but with a doubled unit cell. It is shown that the density made self-consistent for the single formula unit remains so in the doubled cell, for both **lm** and **lmgf**.

Finally two-unit superlattice is converted into a trilayer geometry suitable for the layer Green’s function code **lmpg**. The first unit becomes both the first layer of the active region and the the semi-infinite left lead, and the second becomes a second layer in the active region and the semi-infinite right lead. It is shown that the entire trilayer (consisting of the same formula unit in each layer) continues to be self-consistent.

### Table of Contents

- Preliminaries
- Single Formula Unit
- Doubled Formula Unit
- Conversion from Crystal to Trilayer Geometry
- Footnotes and references

**Single Formula Unit**

```
cp sitein.cspbi
blm --asa --gf --nk=8,8,8 --nl=3 --omax=0.16,0.21,0.25 --wsrmax=3.3 --rdsite --noshorten --nfile --wsite cspbi
lmscell ctrl.cspbi --stack~sort@h3~show@planes~wsite@short@fn=site1
lm -vfile=1 ctrl.cspbi -vnit=0
lmstr -vfile=1 ctrl.cspbi
mpix -n 16 lm -vnit=100 -vfile=1 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --rs=0,1 > out
cp rsta.cspbi rsta1.cspbi [for backup only]
lmgf -vgamma=t -vgf=1 -vnit=-1 -vfile=1 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
```

**Doubled Formula Unit**

```
lmscell -vfile=1 ctrl.cspbi cspbi --stack~rsasa@fn=rsta.cspbi@rdim~file@site1.cspbi@rsasa@rdim@fn=rsta.cspbi~wsite:short:fn=site2:rsasa=rsta2
cp rsta2 rsta.cspbi
lmstr -vfile=2 ctrl.cspbi
lm -vfile=2 ctrl.cspbi -vnit=0 --rs=1,0
mpix -n 16 lm -vnit=1 -vfile=2 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --quit=rho
mpix -n 16 lmgf -vgamma=t -vgf=1 -vnit=-1 -vfile=2 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
```

**Conversion from Crystal to Trilayer Geometry**

```
lmscell -vfile=2 -vpgf=2 cspbi --stack~rsasa@fn=rsta2~setpl:sort:range=9:shorps=0,0,1~wsite:fn=site3:rsasa=rsta3 >> out.lmscell.cspbi3
cp rsta.cspbi rsta2.cspbi [for backup only]
cp rsta3 rsta.cspbi
lmpg -vgamma=t -vpgf=1 -vnit=0 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
lmstr -vgamma=t -vpgf=2 -vnit=1 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
lmpg -vgamma=t -vpgf=2 -vnit=1 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
lmpg -vgamma=t -vpgf=1 -vnit=-1 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
```

### 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, or run serially by removing the prefix `mpix -n 16`

where it appears.

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

### Single Formula Unit

In this section the *ctrl* and *site* files are prepared.

Copy the contents of the box below to file *site.cspbi*. It is a site file for a pseudocubic structure of CsPbI_{3}. (The cubic structure has higher symmetry, but this structure provides a reasonable representation of CsPbI_{3} at room temperature.) Empty spheres have already been added to satisfy the space-filling requirement of the ASA.

```
% site-data vn=3.0 fast io=15 nbas=16 alat=11.82540057 plat= 1.00496473 0.0 0.0 0.00024653 0.9953719 0.0 0.0223297 -0.0000584 1.0183501
# pos
Cs 0.5238550 0.4975648 0.5289579
Pb -0.0336829 -0.0002003 0.9843686
I -0.1007007 0.4975195 0.0129042
I -0.0804758 -0.0001690 0.4726700
I 0.4549649 -0.0001821 0.9282127
E 0.4217061 -0.0000008 0.4664323
E -0.0160774 0.4975539 0.4836709
E 0.4451329 0.4976335 0.9997768
E1 -0.3508811 0.2180172 0.2044577
E1 -0.3509509 -0.2183507 0.2045600
E1 0.2004869 0.2553894 0.2279732
E1 0.2002804 -0.2556073 0.2282405
E1 0.1967869 -0.2731194 0.7477251
E1 0.1969437 0.2729484 0.7477949
E1 -0.2615046 -0.2467878 0.7555042
E1 -0.2614808 0.2464068 0.7555879
```

Copy the contents of the box below to file *ctrl.cspbi*. It is the input (*ctrl*) file for this system.

```
# Autogenerated from init.cspbi using:
# blm --asa --gf --nk=8,8,8 --nl=3 --omax=0.16,0.21,0.25 --wsrmax=3.3 --rdsite --noshorten --nfile --wsite cspbi
# Variables entering into expressions parsed by input
% const nit=10 asa=t
% const met=(asa?1:5)
% const so=0 nsp=so?2:1
% const lxcf=2 lxcf1=0 lxcf2=0 # for PBE use: lxcf=0 lxcf1=101 lxcf2=130
% const ccor=T twoc=F sx=0 gamma=sx scr=4 # ASA-specific
% const gf=1 nz=16 ef=-.238 c3=t pgf=0 sparse=0 # lmgf-specific variables
% const nk1=8 nk2=nk1 nk3=nk2 beta=.3
VERS LM:7 ASA:7 # FP:7
IO SHOW=f HELP=f IACTIV=f VERBOS=35,35 OUTPUT=*
EXPRESS
# Lattice vectors and site positions
%ifdef file>=1
file= site{file}
%endif
file= site
# Self-consistency
nit= {abs(nit)} # Maximum number of iterations
mix= B2,b={beta},k=7 # Charge density mixing parameters
conv= 1e-5 # Convergence tolerance (energy)
convc= 3e-5 # tolerance in RMS (output-input) density
# Brillouin zone
nkabc= {nk1},{nk2},{nk3} # 1 to 3 values
metal= {met} # Management of k-point integration weights in metals
# Potential
nspin= {nsp} # 2 for spin polarized calculations
so= {so} # 1 turns on spin-orbit coupling
xcfun= {lxcf},{lxcf1},{lxcf2} # set lxcf=0 for libxc functionals
SYMGRP find
OPTIONS
ASA[ CCOR={ccor} TWOC={twoc} ADNF=F GAMMA={gamma}] SCR={scr} SX={sx}
HAM PMIN=-1 # Constrain Pnu > P_free
# RDVEXT=2
GF MODE={gf} GFOPTS={?~c3~p3;~p2;}
BZ EMESH={nz},10,{ef}-1,{ef},.5,.3 # nz-pts;contour mode;emin;emax;ecc;bunching
PGF MODE={pgf} SPARSE={sparse} PLATL=0.02232969 -0.00005842 1.01835015 PLATR=0.02232969 -0.00005842 1.01835015
GFOPTS=p3;specfr;declead;padtol=1.001
SPEC
# SCLWSR=11 WSRMAX=3.3 OMAX1=0.16 0.21 0.25
ATOM=Cs Z= 55 R= 3.300000 LMX=2
ATOM=Pb Z= 82 R= 3.300000 LMX=2 IDMOD=0,0,1
ATOM=I Z= 53 R= 3.300000 LMX=2
ATOM=E Z= 0 R= 3.237117 LMX=2
ATOM=E1 Z= 0 R= 2.471028 LMX=1
START CNTROL={nit<=0} BEGMOM={nit<=0}
```

The main part the *ctrl* was autogenerated using the **blm** utility. There are some modifications (added by hand) to the autogenerated template, as the box below explains.

**blm** was used to generate an initial template *ctrl* file from *site.cspbi*

```
$ blm --asa --gf --nk=8,8,8 --nl=3 --omax=0.16,0.21,0.25 --wsrmax=3.3 --rdsite --noshorten --nfile --wsite cspbi
```

**blm**’s command-line arguments are summarized here.

**–asa**and**–gf**target the*ctrl*for the ASA codes**lm**,**lmgf**,**lmpg****–nl=3**,**–omax=0.16,0.21,0.25**,**–wsrmax=3.3**affect parameters in augmentation spheres. These values were selected to allow spheres to fill space, while not letting any sphere get too large.**–rdsite**,**–noshorten**, and**–wsite**cause**blm**to read data from*sitein.cspbi*, and write the structural data without modifying it.

**blm** creates a template *actrl.cspbi*. Some modifications were necessary; the following **diff** command reveals the modifications of the template:

```
$ diff actrl.cspbi ctrl.cspbi
```

The output of **diff** is shown below, with annotations:

< % const ccor=T sx=0 gamma=sx scr=4 # ASA-specific < % const gfmode=1 nz=16 ef=0 c3=t # lmgf-specific variables --- > % const ccor=T twoc=F sx=0 gamma=sx scr=4 # ASA-specific > % const gf=1 nz=16 ef=-.238 c3=t pgf=0 sparse=0 # lmgf-specific variables < ASA[ CCOR={ccor} TWOC=F ADNF=F GAMMA={gamma}] SCR={scr} SX={sx} --- > ASA[ CCOR={ccor} TWOC={twoc} ADNF=F GAMMA={gamma}] SCR={scr} SX={sx}

Variable **twoc** is added to enable control of token **TWOC** Reasons for this are explained in the main text.

Variable **gfmod** was shortened to **gf**. It is use to control program flow in **lmgf** (**MODE={gf}**)

Variable **ef** is used by **lmgf** to set the energy integration window. An estimate of it was obtained by running **lm** to self-consistency (see below).

Variables **pgf** and **sparse** are used by **lmpg** to control program flow there;

< GF MODE={gfmode} GFOPTS={?~c3~p3;~p2;} < BZ EMESH={nz},10,-1,{ef},.5,.3 # nz-pts;contour mode;emin;emax;ecc;bunching --- > HAM PMIN=-1 # Constrain Pnu > P_free > # RDVEXT=2 > > GF MODE={gf} GFOPTS={?~c3~p3;~p2;} BZ> EMESH={nz},10,{ef}-1,{ef},.5,.3 # nz-pts;contour mode;emin;emax;ecc;bunching > PGF MODE={pgf} SPARSE={sparse} PLATL=0.02232969 -0.00005842 1.01835015 PLATR=0.02232969 -0.00005842 1.01835015 > GFOPTS=p3;specfr;declead;padtol=1.001

- A new
**HAM_PMIN=**was added for the same reason**IDMOD**was added (see below) - The energy window in
**BZ_EMESH**was modified **PLATL**and**PLATR**are tokens required by**lmpg**to define the left- and right- semi-infinite regions;**PGF_GFOPTS**carries switches specific to**lmpg**. Note that they are the same vectors as the third vector in**PLAT**(inspect the*site*file).

< ATOM=Pb Z= 82 R= 3.300000 LMX=2 --- > ATOM=Pb Z= 82 R= 3.300000 LMX=2 IDMOD=0,0,1

The ASA can be a bit tricky when floating the band center-of-gravity for localized orbitals. The Pb *d* orbital is such a case, and **IDMOD=0,0,1** freezes its logarithmic derivative at its starting point (which is a safe value). **HAM_PMIN=-1** (above) also restricts logarithmic derivatives in all channels to remain above the free-electron value

51c56 < ATOM=E1 Z= 0 R= 2.471028 LMX=2 --- > ATOM=E1 Z= 0 R= 2.471028 LMX=1

The radius of the second empty sphere (**E1**) is small enough to reduce the basis to *sp* orbitals on that site. This significantly reduces the rank of the hamiltonaian at minimal loss of accuracy.

When later we convert to a trilayer geometry, the trilayer axis will be parallel to **PLATL** and **PLATR** (approximately along *z*). Definition of the trilayer requires partitioning of the sites into principal layers.

**lmscell** can make this construction automatically using its superlattice editor, but it requires that the sites be ordered by increasing distance along the trilayer axis. It is advisable at this stage to rearrange the sites to simplify synchronization of the crystal and trilayer geometries.

The following instruction

```
$ lmscell ctrl.cspbi --stack~sort@h3~show@planes~wsite@short@fn=site1
```

- sorts the basis by
**h3**(**h3**is the trilayer axis), - displays the projection of (now reordered) site positions onto
**h3**(it is called**h**in the printout) - writes the reordered structural information to file
*site1.cspbi*.

*site.cspbi* and *site1.cspbi* contain the same structural information but order the basis differently. The rest of the tutorial will use *site1.cspbi* (or later *site2.cspbi*); which file is used will be controlled from the command line.

Editor instructions are documented here.

#### Self-consistent density through the band code

In the following we use the band code **lm** to a self-consistent ASA density

```
$ lmstr -vfile=1 ctrl.cspbi
$ lm -vfile=1 ctrl.cspbi -vnit=0
$ mpix -n 16 lm -vnit=100 -vfile=1 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --rs=0,1 > out
$ cp rsta.cspbi rsta1.cspbi
```

The first two steps create the structure matrix and a (crude) trial density, and the next step converges the density to self-consistency.

The (last) copy command is not necessary but it is useful if you want to redo the tutorial starting from this point. (Just copy *rsta1.cspbi* to *rsta.cspbi*).

In the iterations to reach self-consistency, there are several points to note:

At the end of the file you should see

`it 17 of100 ehf= -99965.003750 ehk= -99965.003749 From last iter ehf= -99965.003748 ehk= -99965.003748 diffe(q)= -0.000001 (0.000027) tol= 0.000010 (0.000030) more=F c nit=100 file=1 met=2 ccor=0 twoc=1 ehf=-99965.0037495 ehk=-99965.0037494 cpudel ... Time this iter: time(s): 1.58 total: 28.1s Jolly good show ! You converged to rms DQ= 0.000061`

Self-consistency was reached in 17 iterations. Note that the HK and HF functionals are nearly identical.

- See what Fermi level
**lm**finds:`$ grep Fermi out`

The Fermi energy should be about −0.238 Ry. This is basis for selecting variable

**ef=−.238**in the*ctrl*file.**lm**doesn’t need this tag, but an energy window (part of**BZ_EMESH**) are parameters**lmgf**and**lmpg**require to perform a contour integration. - Note the Hamitonian rank
`$ grep Makidx out`

The Hamiltonian consists of a basis of 104 orbitals. There are 40 orbitals omitted because of the basis for empty sphere

**E1**was reduced (**LMX=1**). - Note the command line argument
**-vfile=1**. When the command-line switch`blm ... --nfile ...`

is used, the following lines are inserted into the*ctrl*file:`%ifdef file>=1 file= site{file} %endif`

If variable

**file**is defined and its value exceeds unity, the preprocessed*ctrl*file picks up an extra line, which if**file=1**, reads`file= site1`

, thus causing structural data to read from this file. Variables

`... -vccor=f -vtwoc=t`

affect tags**OPTIONS_ASA_CCOR**and**OPTIONS_ASA_TWOC**. They turn off**lm**’s “combined correction term”, and also reduce**lm**’s standardd 3-center hamiltonian with a simpler two-center one. Neither of these approximations is necessary (both reduce the accuracy of the bands) but in making them the hamiltonian generated by**lm**is formally identical to the Green’s function**lmgf**makes.**lmgf**has no analogous combined correction term, and while it does have a third order hamiltonian, the third order addition is slightly different from that of**lm**. Thus we employ these approximations to demonstrate that**lm**,**lmgf**, and**lmpg**are essentially equivalent.- Switch
`--rs=0,1`

causes**lm**to write a restart file*rsta.cspbi*(ASA format).

#### Self-consistent density through the Green’s function code

Starting from a self-consistent density generated by **lm**, we should find that the density generated by **lmgf** is the same. This is not trivial, since **lmgf** requires a contour integration that must encompass all the poles of the occupied valence bands. The integration must be performed on the real axis below the deepest pole to the Fermi level.

To confirm that **lmgf** indeed generates a self-consistent density, try

```
$ rm -f mixm.cspbi
$ lmgf -vgamma=t -vgf=1 -vfile=1 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi -vnit=-1 --rs=1,0 --quit=rho > out
```

`-vgamma=t`

tells**lmgf**to work in the “gamma representation.” On the real axis, the poles are equivalent, but there can be false poles as well. They are pushed far from the real axis in this representation. Using this switch can be especially important when third-order (i.e. three-center) potential functions are used.`-vgf=1`

selects the mode**lmgf**operates in (self-consistency cycle)`vnit=-1`

makes the**START**category appear as`START CNTROL=1 BEGMOM=1`

which causes

**lmgf**to begin with energy moments , , and of the density, which are sufficient to generate a density and corresponding potential. The potential is used in turn to construct the hamiltonian for the Green’s function.`--rs=1,0 --quit=rho`

causes**lmgf**to read , , and from*rsta.cspbi*, and also to stop after the charge density is generated.*rsta.cspbi*stores this information by site.

You should find something the following at the end of the output:

```
PQMIX: read 0 iter from file mixm. RMS DQ=1.45e-5
AMIX: nmix=0 mmix=8 nelts=160 beta=0.3 tm=10 rmsdel=1.45e-5
h gamma=1 gf=1 file=1 met=2 ccor=0 twoc=1 nit=-1 ehf=-99965.0036905 ehk=-99965.0036905
```

Deviation from exact self-consistency (**1.45e-5**) is small, establishing that the same density is self-consistent in **lm** and **lmgf**.

Midway through the output you should see

```
--- GFASA (gamma) : make qnu, idos
GFASA: integrated quantities to efermi = -0.238
PL D(Ef) N(Ef) E_band 2nd mom Q-Z
1.923116 31.999960 -20.311026 15.739134 -0.000040
deviation from charge neutrality: -0.00004
```

It shows that the Fermi level (−0.238 Ry) found by **lm** also corresponds to the charge-neutrality point in **lmgf**.

Also inspect the potential shift file *vshft.cspbi*. The first line should read something similar to

```
ef=-0.2380000 vconst=-0.0033597
```

**lmgf** shifted the crystal potential to keep the charge neutrality point fixed at −0.238 Ry. The shift is small (about 3 mRy) and causes the Fermi level to fall just above the valence band maximum. The system is an insulator without a density-of-states at the Fermi level. Thus, the Fermi level (or the shift) is very “soft” and is sensitive to details, e.g. the number of points on the contour energy integration.

### Doubled Formula Unit

The trilayer geometry requires both left semi-infinite and right semi-infinite regions. The original unit cell can serve as both, if we duplicate the unit and double the cell. Then the first unit will correspond to the left region and the second to the right.

To prepare for the trilayer geometry, construct a doubled cell with periodic boundary conditions using the superlattice editor

```
$ lmscell -vfile=1 ctrl.cspbi cspbi --stack~rsasa@fn=rsta.cspbi@rdim~file@site1.cspbi@rsasa@rdim@fn=rsta.cspbi~wsite:short:fn=site2:rsasa=rsta2
```

See the user’s manual for the editor instructions.

This invocation of **lmscell** doubles the unit cell along the third axis, creating both a *site* file for structural data (file *site2.cspbi*) and an *rsta* file for sphere moments *rsta2*. Thus it automatically sets up files for a self-consistent potential in the doubled cell.

To confirm that this is the case, do the following.

- copy
*rsta2*to the appropriately named file - make the structure constants for the doubled cell
- make the potential from the moments in
*rsta.cspbi* - run
**lm**one pass to confirm it is self-consistent in the doubled cell - run
**lmgf**one pass to confirm it is self-consistent in the doubled cell

```
$ cp rsta2 rsta.cspbi
$ lmstr -vfile=2 ctrl.cspbi
$ lm -vfile=2 ctrl.cspbi -vnit=0 --rs=1,0
$ mpix -n 16 lm -vnit=1 -vfile=2 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --quit=rho -vnk3=4
$ mpix -n 16 lmgf -vgamma=t -vgf=1 -vnit=-1 -vfile=2 -vmet=2 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho -vnk3=4
```

A new argument `-vnk3=4`

halves the number of *k* divisions on the third axis, while keeping fixed the first two axes (see the *ctrl* file). Since the unit cell was doubled along the third axis, the number of *k* divisions can be halved on that axis to recover the same mesh.

Search for **RMS DQ** in the output and note that it remains less than 10^{−5} for both **lm** and **lmgf**. Note also that *vshft.cspbi* is essentially the same as before.

### Conversion from Crystal to Trilayer Geometry

Having established a self-consistent density with two formula units for the periodic case, we can proceed to converting the geometry to a trilayer form. In this tutorial we make the active region two layers: layer 0 will be taken from the first unit of the preceding section (as will the left semi-infinite layers); layer 1 will be taken from the second unit of the preceding section (as will the right semi-infinite layers).

**lmpg** requires that the active region be divided into principal layers. **lmscell** has the ability to make the conversion to a trilayer geometry, including partitioning the active region into principal layers. It requires as input **PLATL** and **PLATR**. They are already included in the *ctrl* file; note that are identical to the third axis in the original *site* file (*site1.cspbi*).

The active region will consist of the same lattice vectors as in the periodic system, only now the third lattice vector does not repeat periodically. Instead, attached on the left are an infinite stack of repeating layers each with the third vector **PLATL**; similarly stacked are repeating layers on the right with the third vector **PLATR**.

The assignment of principal layers proceeds as follows.

Recall that atoms in

*site1.cspbi*are ordered by projection along the third axis. All sites with projection between 0 and the height corresponding to**PLATL**will form the one layer in the left semi-infinite region; similarly all sites between the last site and sites below it with projection**PLATR**form a layer in the right semi-infinite region.Once the end layers are established (the are taken from the active region and shifted left and right by

**PLATL**and**PLATR**respectively),**lmscell**will start from the left layer, form the first layer (called layer 0) by finding all sites with a connecting vector that touch any site in the left semi-infinite layer, the layer 1 by finding all sites that touch any site in layer 0, and so on.

This is accomplished with the following instruction:

```
$ lmscell -vfile=2 -vpgf=2 cspbi --stack~rsasa@fn=rsta2~setpl:sort:range=9:shorps=0,0,1~wsite:fn=site3:rsasa=rsta3 >> out.lmscell.cspbi3
```

Switch `-vfile=2`

tells **lmscell** (and later **lmpg**) to read structural information from *site2.cspbi*. In this tutorial, the active region has only these two layers, so the *site2.cspbi* can be used both for the periodic case and the trilayer geometry; only the third lattice vector has different meanings in two two contexts, as noted above.

Switch `-pgf=2`

gives token **PGF_MODE** the value 2. When it is nonzero, **lmscell** will try to read **PLATL** and **PLATR**.

Breaking down the syntax of the `--stack`

command (see also the instruction manual):

**~rsasa@fn=rsta2~**tells**lmscell**to read energy moments of the current (2 formula-unit) structure into memory.**lmscell**has not yet been told to switch to a trilayer geometry.**~setpl**causes**lmscell**to convert the periodic geometry to a trilayer form, and assign a principal layer index to each site following the prescription noted above.

The modifiers serve the following functions:**:sort**re-orders sites by increasing height along the PL axis**:range=9**is required. This supplies to**lmscell**the length of the longest connecting vector It is defined by the range of the structure constants, and you should pick a number consistent with this range.**:shorps=0,0,1**places site in 1st quadrant on the third axis

**~wsite:fn=site3:rsasa=rsta3**causes**lmscell**to write a new*site*file, and a moments file*rsta3*which contains moments for all sites in the active region, and also sites in the left- and right- principal layers.

Compare *site3.cspbi* to *site2.cspbi*. They share in common the same site and lattice vectors, but *site3.cspbi* also includes information about principal layers.

Compare *rsta3* to *rsta2*. *rsta3* has moments for 64 sites, while *rsta3* has moments for 32 sites. The left and right layers contain 16 sites each (this information was copied from the first and last layer in the active region)

#### Self-consistent density in the trilayer geometry

This section demonstrates that the density in the trilayer geometry is self-consistent. This should be the case since all the layers are identical, but it is a nontrivial check because **lmpg** does not make use of periodic boundary conditions on the third axis.

There are in fact three regions that can be checked: the left- and right- semi-infinite regions, and the active region.

The following step is not needed, but you might want to save the existing *rsta* file so you can repeat earlier steps without starting from the beginning:

```
$ cp rsta.cspbi rsta2.cspbi
```

- copy
*rsta3*to the appropriately named file - make the potential from the moments in
*rsta.cspbi*for all atoms in the trilayer geometry - make the structure constants for the trilayer geometry
- run
**lmpg**in mode 2, to converge the left- and right- regions to self-consistency. By construction, the potential should already be self-consistent: running**lmpg**confirms that it is. - run
**lmpg**in mode 1, to converge the active region to self-consistency. Similarly in this case the potential should already be self-consistent.

```
$ cp rsta3 rsta.cspbi
$ lmpg -vgamma=t -vpgf=1 -vnit=0 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
$ lmstr -vgamma=t -vpgf=2 -vnit=1 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
$ lmpg -vgamma=t -vpgf=2 -vnit=1 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
$ lmpg -vgamma=t -vpgf=1 -vnit=-1 -vfile=3 -vccor=f -vtwoc=t ctrl.cspbi --rs=1,0 --quit=rho
```

You should find that the end layers are converged to a precision of order 10^{−5}, and the trilayer converged to a precision of order 10^{−4}. Both of these are acceptable.

Inspect *vshft.cspbi* now. You should see something similar to the following:

```
ef=-0.2380000 vconst=0.0152048 vconst(L)=-0.0017305 vconst(R)=-0.0017305
```

Potential shifts for the end layers are generated in the self-consistency cycle for the end layers (**PGF_MODE=2**), and shift for the active layer (**vconst=**) in the self-consistency cycle for the trilayer (**PGF_MODE=1**).

### Footnotes and references

Implementation of Green’s functions in the ASA context, and its connection to the LMTO-ASA hamiltonian is described in detail in O. K. Andersen, Z. Pawlowska, and O. Jepsen, Phys. Rev. B 34, 5253 (1986).

Implemention of the layer Green’s function technique, including the non-equilibrium capability can be found in S. V. Faleev, et al, Phys. Rev. B71, 195422 (2005).