# Nb/Ni Superlattice

The electronic structure of a Nb(110)/Ni(111) superlattice is created using Questaal’s superlattice editor *lmscell*, and the interface relaxed with *lmf*, using a constrained form of molecular statics.

JMRAM devices typically consist of films of Nb, oriented on the (110) direction, followed by a Cu spacer layer, followed by Ni, with another Cu spacer layer and finally capped with Nb. Cu and Ni are oriented along (111). Ni and Cu are both fcc with room temperature lattice constants 3.524 Å and 3.615 Å, respectively. Nb on the other hand has a lattice constant of 3.300 Å. To make an A/B superlattice unit cells of the A and B regions in the plane normal to the interface must be coincident. This is a challenging requirement, in part because the lattices are not matched, but also the (110) and (111) are not compatible.

The tutorial begins by rotating Nb [110] axis and the Ni [111] axis respectively to [001]. Then superlattice vectors are identified that are nearly coincident. By applying a small shear to the Ni, the lattice vectors are made exactly coincident.

At this point superlattices can be constructed. Some are made and the site positions at the frontier layers relaxed.

A separate tutorial takes the structural information made in this tutorial to investigate Landauer-Buttiker transport through a Nb/Ni/Nb trilayer.

*Note: at present that tutorial is written for Nb[110]/Ni(110)/Nb([110] trilayer*

### Table of Contents

- Nb/Ni(3
*n*)/Nb superlattices - Preliminaries
- Introduction
- 1. Reconstructed unit cells for elemental Nb and Ni
- 2. Lattice structure of the Nb(
*m*)/Ni(*n*)/Nb(*m*) superlattice - 3. Relaxing the Nb(2)/Ni(3)/Nb(2) superlattice
- Footnotes and references

Sections build on prior sections, though any section can be treated independently by copying all the files indicated in the Setup instructions.

** Section 1. Coincident site Nb and Ni site lattices**

Creates files *site20.nb* and *site42.ni*, and *spec.nb* and *spec.ni*, and also *site10a.nb*, *site10b.nb* and *site14a.ni*, *site14b.ni* and *site14c.ni*.

Requires files *init.nb* and *init.nb* Copy and paste *init.nb* and *init.ni* file to your working directory, or do:

```
touch init.nb; rm *.nb *.ni;
cp ~/lmf/fp/test/{init.nb,init.ni} .
```

Make *site20.nb*,*spec.nb*, *site10a.nb* and *site10b.nb*

```
blm --nfile --ctrl=ctrl --omax=-.01 nb
lmchk ctrl.nb --wspec
cp spec.nb fpspec.nb
blm --nfile --ctrl=ctrl --omax=.14 nb
lmscell --plx~m~-1,-4,0,2,-2,0,1,1,2 '--rot=(1,-1,0)-pi/2,z:pi/4,z:acos(5/sqrt(43))' nb --wsitex~fn=site2 --sort:"x3 x2"
cp site2.nb site20.nb
lmchk ctrl.nb --wspec > out.lmchk.nb
lmscell -vfile=20 ctrl.nb '--stack~rmsite@targ=11:20~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,1/sqrt(2)~wsite@fn=site10a'
lmscell -vfile=20 ctrl.nb '--stack~rmsite@targ=1:10~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,1/sqrt(2)~wsite@fn=site10b'
```

Make *site42.ni*,*spec.ni*, *site[abc].ni*

```
blm --nfile --ctrl=ctrl --omax=-.01 ni
lmchk ctrl.ni --wspec > out.lmchk.ni
cp spec.ni fpspec.ni
blm --nfile --ctrl=ctrl --omax=.14 ni
lmchk ctrl.ni --wspec > out.lmchk.ni
lmscell '--rot=z:pi/4,x:acos(sqrt(1/3)),z:acos(4/sqrt(19))' --plx~m~3,0,2,-1,0,4,-1,3,-1 ni --wsitex~fn=site2 --sort:"x3 x2"
echo 3.27871926 0.0 0.0 1.82998284 2.15665546 0.0 0.0 0.0 1.41421356 | mcx -f9f12.6 . -s6.226 | sed 's/rows 1 cols 9 real/rows 3 cols 3 real/' > pnb
echo 3.082207 0.0 0.0 1.62221421 1.9668302 0.0 0.0 0.0 1.73205081 | mcx -f9f12.6 . -s6.643 | sed 's/rows 1 cols 9 real/rows 3 cols 3 real/' > pni
cp ctrl.ni ctrl.ni.bak
echo STRUC >> ctrl.ni
echo ' DEFGRD= ' `mcx -f9f12.8 pni -i pnb -x -t -a:nr=1 s s -w:nohead . | sed 's:[0-9.]*$:1:'` >> ctrl.ni
lmchk -vfile=2 ni --wsitex~fn=site3
cp ctrl.ni.bak ctrl.ni
lmchk -vfile=2 ni | grep vol; lmchk -vfile=3 ni | grep vol
sed -i.bak 's:\(plat.*[1-9.]*\)$:\1/1.00983:' site3.ni
lmchk -vfile=2 ni | grep vol; lmchk -vfile=3 ni | grep vol
lmchk -vfile=3 ni --wsitex~scala=6.226/6.643~fn=site42
sed -i.bak 's:\(plat.*1.83\).*$:\1:' site42.ni
lmchk ctrl.ni --wspec > out.lmchk.ni
lmscell -vfile=42 ctrl.ni '--stack~rmsite@targ=15:42~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,.61~wsite@fn=site14a'
lmscell -vfile=42 ctrl.ni '--stack~rmsite@targ=1:14,29:42~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,.61~wsite@fn=site14b'
lmscell -vfile=42 ctrl.ni '--stack~rmsite@targ=1:28~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,.61~wsite@fn=site14c'
```

** Section 2. Initial estimate for Nb(110)/Ni(111)**

Creates files *site20.nb* and *site42.ni*.

#### Nb/Ni(3*n*)/Nb superlattices

This section requires the files *ctrl.nb*, *site20.nb*, *spec.nb*, *site42.ni*, *spec.ni*.

Set initial conditions with

```
touch spec.nb; rm *.nb *.ni
cp ../input.fp.Nb110.Ni111/{init.nb,site20.nb,spec.nb,site42.ni,spec.ni} .
cp ~/lmf/fp/test/nbni/{init.nb,site20.nb,spec.nb,site42.ni,spec.ni} .
```

Build an ideal superlattice and confirm Nb/Ni and Ni/Nb have equal spacing

```
blm --nfile --ctrl=ctrl --omax=-.01 nb
lmscell -vfile=20 ctrl.nb -vdz='1/2/sqrt(2)-.61/2' --stack~file@site42.ni@dpos=0,0,-dz~file@site20.nb~wsite@fn=sitex
cp sitex.nb sitein.dat
blm --nfile --rdsite=sitein --rdspec~fn=spec.nb~fn=spec.ni --usefiler --nosort --shorten=no --ctrl=ctrl dat
lmplan dat --pledit:q
```

Make *site82.dat* and get shift in positions from `--mino`

relaxation.

```
lmchk dat --wsite~short~fn=site82 --wpos=pos0
lmchk dat --mino~dxmx=.001~xtol=.0001~nkill=10~maxit=200~11:34,49:72 --wsite~short~fn=site82 --wpos=pos
mcx pos.dat pos0.dat -- > dpos.dat
```

** Section 3. Relaxing the Nb(2)/Ni(3)/Nb(2) superlattice**

OLD:

```
cp site82.dat sitein.nbni
blm --autobas~loc=5 --rsham~el=-.3,-1.2~rscut=1.6 --dhftol=45 --molstat~xtol=0~gtol=0.01 --nfile --rdsite --omax=-.01 --nosort --shorten=no --nit=20 --nk=4,5,2 --gmax=10 --mixb=.1 --frzwf --wsite~ib=1:10,35:48,73:82~rlx=000~long~fn=site82 nbni
cp actrl.nbni ctrl.nbni
```

NEW … allow Nb/Ni spacing to vary (variable dp3)

```
cp ../input.fp.Nb110.Ni111/ctrl.nb ../input.fp.Nb110.Ni111/site20.nb ../input.fp.Nb110.Ni111/site82.dat ../input.fp.Nb110.Ni111/site42.ni ../input.fp.Nb110.Ni111/dpos.dat ../input.fp.Nb110.Ni111/fpspec.nb ../input.fp.Nb110.Ni111/fpspec.ni .
cp site82.dat sitein.nbni
lmscell -vfile=20 ctrl.nb -vdz0='1/2/sqrt(2)' -vdz=dz0-.61/2 -vdp3=-.00 --stack~file@site42.ni@dpos=0,0,-dz+dp3/2~file@site20.nb@dpos=0,0,dp3~setlbl@spid=Nbx@targ=11:20,nbas-19:nbas-10~setlbl@spid=Nix@targ=21:34,nbas-33:nbas-20~pad=0,0,dp3~wsite@fn=sitex
cp sitex.nb sitein.nbni
blm --rdspec~fn=fpspec.nb~fn=fpspec.ni --usefiler --rdpos=dpos.dat --autobas~loc=5 --rsham~el=-.3,-1.2~rscut=1.6 --dhftol=45 --molstat~xtol=0~gtol=0.01 --nfile --rdsite --omax=-.01 --nosort --shorten=no --nit=20 --nk=4,5,2 --gmax=10 --mixb=.2 --frzwf --wsite~ib=1:10,35:48,73:82~rlx=000~long~fn=site82 nbni
cat actrl.nbni | sed s/NFORCE=0/NFORCE=20/ > ctrl.nbni
```

Repeat for dp3=-.04 and -.08; use same ctrl file

Extract pos0 = ideal positions, pos5 = relaxed positions and dpos5 the difference

```
lmf ctrl.nbni -vfile=3 -vnk1=4 -vnk2=4 -vnk3=1 -vfrzwf=t -vconv=1e-4 -vconvc=1e-3 --rpos=pos36 --rs=11,0,0 -vnit=0 --wpos=pos5
lmchk ctrl.nbni --rdpos~scale=-1~fn=dpos.dat --wpos=pos0
mcx -f3f13.7 pos5.nbni pos0.nbni -- > dpos.2+3+2
```

### Preliminaries

This tutorial makes extensive use of many of the Questaal packages. including *lmscell*, *lmplan* *lmf*, and miscellaneous utilities such as **vextract**. They are assumed to be in your path, or in the directory where you cloned the Questaal package (assumed here to be *~/lm*).

It invokes *mpix* to execute MPI jobs. *mpix* is a proxy for *mpirun*, the default command that launches MPI jobs. We assume your system has 16 processors on a core, but you can use any number, or run in serial.

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

### Introduction

This tutorial builds a Nb(110)/Ni(111) superlattice, using Questaal tools *blm*, *lmscell* and *lmchk*. Since Nb[110] and Ni[111] are poorly lattice matched, unit cells in the plane normal to [110] and [111] respectively must be found that are nearly coincident. Within the coincident lattice, bonds across the interface are very irregular, and there is a large attendant lattice relaxation which must be dealt with. This is done in two steps: a trial relaxation that minimizes some measure of overlaps between atomic density, and subsequently using the LDA code *lmf* to relax the structure by minimizing the forces.

#### Basic strategy

To stack unit cells of Nb and Ni on top of one another, the two lattice vectors normal the stacking plane must match (they must be coincident). This entails:

- a reconstruction of bcc Nb with 10 atoms in a plane normal to [110]. The Nb lattice is also rotated orienting the [110] axis along
*z*. - a reconstruction of fcc Ni with 14 atoms in a plane normal to [111]. The Ni lattice is also rotated orienting [111] axis along
*z*. - The Ni lattice is sheared slightly to make Nb and Ni exactly coincident in the plane normal to
*z*. It is also strained along*z*to partially accommodate the volume change from the in-plane shear, according to Poisson’s ratio.

The prescription is in outline:

- Use
*blm*to build input files for elemental Nb and Ni in reconstructed lattices. - Use the supercell maker
*lmscell*in combination with*blm*to generate different reconstructions of the elemental lattices that can be joined together. - Use the superlattice editor
`lmscell --stack`

to stack reconstructed Nb and Ni lattices into superlattices. - Use
*lmchk*to make a rough estimate of the atomic relaxation of frontier atoms at the Nb/Ni interfaces, followed by the LDA code*lmf*to perform a constrained relaxation the interface, allowing atoms in the frontier planes to relax until forces on them are small.

### 1. Reconstructed unit cells for elemental Nb and Ni

This tuturial requires only two *init* files corresponding to element Nb and Ni to complete.

This tutorial concludes with a self-consistent *finish*

Three files are required to use *lmf* to make the self-consistent electronic structure for a specified lattice. 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.

We need to assemble files in several contexts:

- For the end leads, Nb
- For the central region, Ni
- Ni/Nb superlattices, built so that the interface can be relaxed in a constrained manner. The constraint is necessary for later use, when superlattices are converged into trilayer (… Nb Nb / Ni / Nb Nb …) suitable for study of transport of an active region. The semi-infinite left and right leads are required to have a periodic lattice.

#### 1a. Structural information for elemental Nb and Ni

Copy the boxes below into files *init.nb* and *init.ni*. They contain the requisite structural information for elemental bcc Nb and fcc Ni. Here we use the 0K lattice constants.

```
# Init file for Nb. Lattice constant at 0K
LATTICE
ALAT=6.226 PLAT= -.5 .5 .5 .5 -.5 .5 .5 .5 -.5
SITE
ATOM=Nb X=0 0 0
```

and

```
# Init file for Ni. Lattice constant at 0K
LATTICE
ALAT=6.643 PLAT= .5 .5 .0 .5 .0 .5 .0 .5 .5
SPEC
ATOM=Ni MMOM=0,0,0.6 PZ=0,0,0
SITE
ATOM=Ni X=0 0 0
```

#### 1b. Reconstructed Nb, and Ni, coincident site lattices

##### Initial reconstruction of Nb

Nb is bcc, and it will be rotated so the [110] direction points along *z* which will be axis of the superlattice.

Build a template *ctrl* and *site* file for bcc Nb. This file is incomplete but it sufficient for steps needed for the reconstruction.

```
blm --nfile --ctrl=ctrl nb
```

`--nfile`

tells *blm* to insert the following lines into *ctrl.nb* :

%ifdef file | file==0 file= site{file} %endif

The other switches supply information *blm* cannot determine on its own: the *k* mesh and the plane wave cutoff for the interstitial density.

This allows Questaal codes to read structural data from *site*n*.ext*, when preprocessor variable **file** is assigned to *n*. These lines have no effect if **file** is not assigned, in which case Questaal codes will read from the usual file *site.ext*.

Next we rotate the crystal with the following (admittedly odd) construction. The total rotation consists of a sequence of three separate ones: first, a rotation of the [110] axis to the [001], then a 45 degree rotation the new *z* axis, finally another rotation around the same axis. `--rot=(1,-1,0)-pi/2`

rotates [110] to *z* which will be the superlattice axis. The latter rotations will be made clear when we align Ni and Nb.

```
lmscell --noplx '--rot=(1,-1,0)-pi/2,z:pi/4,z:acos(5/sqrt(43))' nb --wsitex~fn=site1
```

Next, make a superlattice of this rotated cell. It reads site information from *site1.nb* which *lmscell* generated in the previous instruction.

```
lmscell -vfile=1 --plx~m~-1,-4,0,2,-2,0,1,1,2 nb --wsitex~fn=site2 --sort:"x3 x2"
```

This superlattice is artfully constructed from a linear combination of integer multiples bcc lattice vectors. These particular multiples yield the following lattice vectors:

```
plx plx/plat qlx
3.27872 0.00000 0.00000 -1.00000 -4.00000 0.00000 0.30500 -0.25880 0.00000 ← PLAT(1)
1.82998 2.15666 0.00000 2.00000 -2.00000 0.00000 0.00000 0.46368 0.00000 ← PLAT(2)
0.00000 0.00000 1.41421 1.00000 1.00000 2.00000 0.00000 0.00000 0.70711 ← PLAT(3)
```

This superlattice consists of two planes 10 atoms each. The first plane passes through the origin, and the second through , which is the spacing between (110) planes in the bcc lattice. To confirm this, do

```
lmplan -vfile=2 nb --pledit~q
```

bcc planes are stacked in an ABAB… pattern. A and B planes are equivalent but B is translated in the lateral direction. The third plane at a height , coincides with the first without a lateral translation. Thus this rotation/superlattice combination puts **PLAT(3)** along [001] with a length .

##### Initial reconstruction of Ni

Ni is fcc, and it will be rotated so the [110] direction points along *z*.

Make a template *ctrl* and *site* file for fcc Ni:

```
blm --nfile --ctrl=ctrl ni
```

In the Ni case we will rotate the [111] to *z* using `--rot=z:pi/4,x:acos(sqrt(1/3))`

. An additional rotation `z:acos(4/sqrt(19))`

is made for reasons that will soon become clear. We also need a superlattice of this rotated cell, with different superlattice vectors than the Nb case. We can accomplish both the rotation and the superlattice in one step as follows:

```
lmscell '--rot=z:pi/4,x:acos(sqrt(1/3)),z:acos(4/sqrt(19))' --plx~m~3,0,2,-1,0,4,1,1,1 ni --wsitex~fn=site1 --sort:"x3 x2"
```

Compare the first lines of *site1.ni* and *site2.nb*. The lattice vectors **PLAT(1)** and **PLAT(2)** (when scaled by their respective lattice constants) are very similar.

This particular superlattice consists of 14 atoms all in a single plane normal to the new *z*, which is [111] in the unrotated frame. To confirm this, do

```
lmplan -vfile=1 ni --pledit~q
```

Ni is fcc whose planes are stacked ABC. We will use as a basic unit a superlattice with three planes (42 atoms) that include all the ABC planes. As needed we will split the ABC layers into single planes.

```
lmscell '--rot=z:pi/4,x:acos(sqrt(1/3)),z:acos(4/sqrt(19))' --plx~m~3,0,2,-1,0,4,-1,3,-1 ni --wsitex~fn=site2 --sort:"x3 x2"
```

You can confirm that **PLAT(1)** and **PLAT(2)** in files *site1.ni* and *site2.ni* are the same; only the third axis now points along *z*.

*lmplan* confirms that the new cell of 42 atoms consists of 3 planes, normal to [111] in the unrotated frame.

```
lmplan -vfile=2 ni --pledit~q
```

You can confirm that sites 15:28 are uniformly translated relative to sites 1:14, and sites 29:42 are translated by twice that amount. Each plane is spaced by in units of the lattice constant.

##### Nb and Ni, coincident site lattices

Note that *lmscell* created the superlattices using `--wsitex`

instead of `--wsite`

. This saves the basis vectors in the form of fractional multiples of the lattice vectors rather than Cartesian coordinates. This will be important for the shear which we need to do next, because to stack Nb and Ni into a superlattice **PLAT(1)** and **PLAT(2)** have to be exactly coincident.

From each of the files *site2.nb* and *site2.ni* , extract the lattice vectors and scale by their respective lattice constants.

```
echo 3.27871926 0.0 0.0 1.82998284 2.15665546 0.0 0.0 0.0 1.41421356 | mcx -f9f12.6 . -s6.226 | sed 's/rows 1 cols 9 real/rows 3 cols 3 real/' > pnb
echo 3.082207 0.0 0.0 1.62221421 1.9668302 0.0 0.0 0.0 1.73205081 | mcx -f9f12.6 . -s6.643 | sed 's/rows 1 cols 9 real/rows 3 cols 3 real/' > pni
```

The `sed`

instruction is appended to render *pnb* and *pni* into the 3×3 matrices they actually are. Inspect them and you should see:

```
% rows 3 cols 3 real
20.413306 0.000000 0.000000 11.393473 13.427337 0.000000 0.000000 0.000000 8.804894
% rows 3 cols 3 real
20.475101 0.000000 0.000000 10.776369 13.065653 0.000000 0.000000 0.000000 11.506014
```

This shows that the superlattice vectors are nearly coincident, with the first vector in each case pointing along *x*. This explains why these particular superlattice vectors and rotations were chosen.

We must shear one or the other lattice to render the first two vectors **PLAT(1)** and **PLAT(2)** exactly coincident. (**PLAT(3)** has no such constraint and we will deal with it later.) We elect to shear the Ni lattice, since real Ni layers are disordered anyway.

Do the following:

```
mcx -f3f12.8 pni -i pnb -x
```

You should get

```
0.99698194 0.00000000 0.00000000
0.04972026 1.02768205 0.00000000
0.00000000 0.00000000 0.76524277
```

This matrix tells us how much we have to shear Ni to match Nb. For now we are interested only in the 2×2 subblock since the third axis is not dictated by matching planes of the two lattices. For now let’s restrict the deformation to the (**PLAT(1)**, **PLAT(2)**) in the plane normal to *z*.

Temporarily add these lines to *ctrl.ni*:

```
STRUC
DEFGRD= 0.99698194 0.00000000 0.00000000 0.04972026 1.02768205 0.00000000 0.00000000 0.00000000 1
```

and create a new file *site3.ni* with the transformed structural data.

```
lmchk -vfile=2 ni --wsitex~fn=site3
```

Remove or comment out **DEFGRD** before continuing.

Token **DEFGRD** supplies a general linear transformation of the lattice: it scales both lattice and basis vectors by a 3×3 matrix. Compare lattice vectors of *site2.nb* and *site3.ni*: you should see that (**PLAT(1)**, **PLAT(2)**) are coincident (remember that **alat** is different in the two cases)

```
head -1 site2.nb; head -1 site3.ni
```

Finally we must decide how much to scale **PLAT(3)**. Were Ni incompressible we would know how much we should strain Ni along *z*: any volume change generated by **DEFGRD** would be compensated by a tensile strain along *z*. To determine the volume change induced by **DEFGRD** do the following (be sure you have removed or commented out the **DEFGRD** line first!):

```
lmchk -vfile=2 ni | grep vol; lmchk -vfile=3 ni | grep vol
```

You should find that the volume changed by this ratio:

```
dval 3153.75620/3078.09526 = 1.02458
```

As a sanity check, copy *site3.ni* to *site3.ni.bak* and scale the 9^{th} element of **plat** by 1/1.02458, or simply use *sed*:

```
sed -i.bak 's:\(plat.*[1-9.]*\)$:\1/1.02458:' site3.ni
```

Run *lmchk* again and confirm the unit cell volumes of *site2.ni* and *site3.ni* are the same. This is the required strain for the incompressible case, which corresponds to a Poisson ratio of 0.5. The opposite limit, Poisson ratio of 0, means that change in lateral area has no effect on tensile strain along *z* (cork is approximately like that). You can look up Poisson ratios for elements here : https://periodictable.com/Properties/A/PoissonRatio.v.html; it is 0.31 for Ni, close to 3/5.

Since Ni isn’t incompressible, this strain is tool large. This volume increase should be reduced by 2/5, which corrresponds to a strain `dval 1+.02458*2/5`

= 1.00983 . Replace the 9^{th} element of **plat** with `dval 1.73205081/1.00983`

= 1.71519. To do it on the command-line, invoke

```
cp site3.ni.bak site3.ni
sed -i.bak 's:\(plat.*[1-9.]*\)$:\1/1.00983:' site3.ni
```

Run *lmchk* as before and confirm the ratios of cell volumes are `dval 3123.05566/3078.09526`

= 1.01461 . Thus the ~2.5% strain is reduced to a ~1.5% strain.

To make the final reconstructed Ni *site* file, do the following

```
lmchk -vfile=3 ni --wsitex~scala=6.226/6.643~fn=site42
```

Command line argument `scala`

causes *lmchk* to scale the lattice constant **ALAT** by a factor and the lattice and basis vectors by its reciprocal. The structure is unchanged, only the distribution between **ALAT** and **PLAT** is changed. In this context **ALAT** becomes the same for Ni and Nb. The name *site42.ni* was used simply to mark how many atoms are in the unit cell.

Note **plat(3,3)** is 1.83006861, which is very close to 1.83, which is a multiple of 0.61. For convenience replace **plat(3,3)** with 1.83, since we will need to divide this cell into three subcells corresponding to each of the three planes. Use your text editor or use *sed*:

```
sed -i.bak 's:\(plat.*1.83\).*$:\1:' site42.ni
```

Before closing this section, it is convenient to construct partial site files with individual Nb or Ni planes

```
lmscell -vfile=20 ctrl.nb '--stack~rmsite@targ=11:20~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,1/sqrt(2)~wsite@fn=site10a'
lmscell -vfile=20 ctrl.nb '--stack~rmsite@targ=1:10~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,1/sqrt(2)~wsite@fn=site10b'
lmscell -vfile=42 ctrl.ni '--stack~rmsite@targ=15:42~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,.61~wsite@fn=site14a'
lmscell -vfile=42 ctrl.ni '--stack~rmsite@targ=1:14,29:42~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,.61~wsite@fn=site14b'
lmscell -vfile=42 ctrl.ni '--stack~rmsite@targ=1:28~assign@z=pos(3,1)~addpos@dx=0,0,-z~p3@0,0,.61~wsite@fn=site14c'
```

### 2. Lattice structure of the Nb(*m*)/Ni(*n*)/Nb(*m*) superlattice

The Nb(*m*)/Ni(*n*)/Nb(*m*) superlattice (*m* and *n* refer to numbers of planes in each region) has an Nb/Ni and an Ni/Nb interface; also periodic boundary conditions apply and there is a Nb/Nb interface as well. Provided the right layer ends in a B plane and the left begins with an A plane, the Nb/Nb interface will just be a continuum of 4 bcc Nb planes along *z* ([110] in the unrotated lattice). Periodic boundary conditions requires that cannot be translated with respect to each other, so we impose the constraint that the left Nb lead begins with an untranslated A plane and the right B lead ends with an untranslated B plane.

The Nb/Ni and Ni/Nb interfaces are more complicated. The main difficulty is that stacking ideal Nb and Ni planes will bring some Ni-Nb bonds much too close together while leaving void regions. Atoms can be relaxed, but there are some constraints we want to follow. We would like to structure the stacks so the superlattice is built in a lego-like fashion, layer by layer. That way construction can be automatic and easily generated for any number of layers. This is greatly facilitated by restricting atomic relaxation to atoms on frontier planes. This is an approximation, but it is a reasonable one. In any case the real material will be highly non-ideal, which is very difficult to model, so there is little point in trying to allow more general relaxations. The main thing is that the frontier atoms are not pressed too close together or kept too far apart. The fact that the displacements are limited only to the interface greatly simplifies the adoption the lego-like, layer-by-layer construction, since each layer is essentially independent of the others.

To stack Ni on Nb, we must address these issues:

- What is the spacing between Nb and Ni planes?

As a first guess we will take the average of the Nb-Nb and Ni-Ni spacing. - Should the Ni lattice be translated relative to the Nb?

In principle this question would be resolved by lattice relaxations. However, the interfaces are complicated and there may be many metastable configurations where forces vanish.

To handle any number of Ni layers we need to consider three distinct configurations:

- the
**Nb/**(**A’B’C’**)_{n}**/Nb**family with 3*n*Ni planes - the
**Nb/**(**A’B’C’**)_{n}**A’/Nb**family with 3*n*+1 Ni planes - the
**Nb/**(**A’B’C’**)_{n}**A’B’/Nb**family with 3*n*+2 Ni planes

Within any family the interfaces will be the same, but between families they will not, and each family must be dealt with independently of the others.

#### 2a. Nb/Ni(3*n*)/Nb superlattices

We will distinguish the three kinds of lattices by denoting them as 3*n*, 3*n*+1 and 3*n*+2 superlattices, respectively. This section builds the simplest case, the Nb(2)/Ni(3)/Nb(2) superlattice. If the files below are not already present in your working directory set up the directory with

```
cp ~/lmf/fp/test/nbni/{ctrl.nb,site20.nb,spec.nb,site42.ni,spec.ni,site14[abc].ni,site10[ab].nb} .
```

Even without any translations or reconstructions, there are several possible stacking arrangements. Consider the following stacking possibilities: AB123AB, AB231AB, AB312AB. Here A and B refer to the two Nb planes, and 1,2,3 refer to the three Ni planes. The simplest choice: **AB|123|AB** turns out to be the best among the three possibilities. Build the stack with the superlattice editor *lmscell* starting with the reconstructed Nb lattice as a base.

```
lmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.00 --stack~file@site42.ni@dpos=0,0,-dz~file@site20.nb~wsite@fn=sitex
```

This makes the **AB|123|AB** superlattice, using the average of Nb-Nb and Ni-Ni spacings as the spacing between Nb and Ni layers. The superlattice structure is written to file *sitex.nb*.

`-vfile=20 ctrl.nb`

tells*lmscell*to build the stack starting from the recontructed 20-atom Nb cell as a base.The Nb-Nb spacing is and the Ni-Ni spacing is . is half this difference.

**dp3**isn’t used now, but later it will become an ajustment to the spacings between frontier layers.`--stack`

invokes the superlattice editor. It is invoked in batch mode with**~**separating commands and does the following:`~file@site42.ni@dpos=0,0,-dz`

stacks the 42-atom reconstructed Ni layer on top of Nb. Absent the modifier**dpos**, the first Ni plane would be placed at the boundary of the 20-atom reconstructed Nb layer, which lies at a height above the topmost Nb layer. Rigidly shifting the 42 Ni atoms by causes the Nb-Ni spacing to take the average of the Nb-Nb and Ni-Ni spacings.`~file@site20.nb`

stacks the 20atom recontructed Nb layer on top of the stack. Since the Ni was shifted by , the second Ni-Nb is automatically at the average of Ni-Ni and Nb-Nb spacings.`~wsite@fn=sitex`

writes the file*sitex.nb*, which is a*site*file of the superlattice.

Use *blm* to construct a temporary *ctrl* and *site* file with extension *.dat*. At this stage we only want a *ctrl* file with enough information to examine the lattice structure. *blm* can read structural information from a *site* file, a feature which we use here:

```
cp sitex.nb sitein.dat
blm --nfile --rdsite=sitein --rdspec~fn=spec.nb~fn=spec.ni --usefiler --nosort --shorten=no --ctrl=ctrl dat
```

Command-line arguments to *blm* are documented on this page.

`--nfile`

was explained above`--rdsite=sitein`

tells*blm*to read structural data from*sitein.dat*`--rdspec~fn=spec.nb~fn=spec.ni`

supplies*blm*with species-specific information gleaned from*spec.nb*and*spec.ni*, notably in this case sphere radii.`--usefiler`

tells*blm*to accept without change the muffin-tin radius (**R**) it read from*spec.nb*and*spec.ni*. It is useful in this context because will obtain a crude estimate for reconsructiong by minimizing overlaps between spheres on neighboring sites.`--nosort`

suppresses re-ordering of site positions`--shorten=no`

suppresses automatic shortening of site positions by adding multiples of lattice translation vectors.`--ctrl=ctrl`

writes the template input file directly to*ctrl.dat*(instead of the default*actrl.dat*)

This lattice has 20+42+20 = 82 sites. Use *lmscell* to see the structure of stack

```
lmplan dat --pledit:q
```

Inspect the table below this line

```
ib Class Plane PL x-x.h y-y.h z-z.h h del h
```

**delh** indicates the change in projection along z compared to the preceding site. **delh** is always zero except for the first atom in a new plane:

```
1 1:Nb 1 0 0.00000 0.00000 0.00000 0.00000 0.70711
...
11 21:Nb11 2 0 2.05873 0.10783 0.00000 0.70711 0.70711
...
21 2:Ni 3 0 0.00000 0.00000 0.00000 1.36566 0.65855
...
35 30:Ni15 4 0 0.43389 0.05135 0.00000 1.97566 0.61000
...
49 58:Ni29 5 0 0.86779 0.10270 0.00000 2.58566 0.61000
...
63 41:Nb21 6 0 0.00000 0.00000 0.00000 3.24421 0.65855
...
73 61:Nb31 7 0 2.05873 0.10783 0.00000 3.95132 0.70711
```

This confirms that planes are stacked (Nb, Nb, Ni, Ni, Ni, Nb, Nb) with Nb-Nb spacings **0.70711**, Ni-Ni spacings **0.61** and Nb-Ni spacings **0.65855** (the average of **0.70711** and **0.61**). Note also the very first site has a spacing **0.70711**. This is the spacing between the last plane, sites 73:82, and the first plane, sites 1:10. It is , as it should be because of periodic boundary conditions.

Next we need to examine Ni-Nb bonds, and relax frontier sites to space bonds more evenly. The simplest way to do this is to minimize sphere overlaps. *blm* automatically determined sphere radii in the elemental compounds by enlarging the sphere radius to match the Wigner Seitz cell volume, and thus making them overlap: 14% in the bcc case, 11% in the close-packed fcc case. (*blm* imported these radii into *ctrl.dat* via `--rdspec`

and `--usefiler`

.) Such geometry violation is always present in ASA calculations, but in this context we want to use the sphere overlaps as a measure of how bonds are spaced.

Do the following

```
lmchk dat --wpos=pos0
```

Near the end you should see this line

```
OVMIN, 1426 pairs: fovl = 3.78727e-5 <ovlp> = 18.3% max ovlp = 36%
```

The maximum overlap is large, but slightly less than what it would be if we had chosen interfaces other than **AB|A’B’C’|AB**, e.g. **ABA|A’B’C’|AB** or **AB|A’B’C’|BAB**.

*lmchk* has a facility to move atoms that minimize sphere overlaps. It does this by minimizing an objective function (e.g. the sixth power of the sum of all positive overlaps). We will use this facility to get an initial estimate of the reconstruction. As noted above, we do not want to move all sites, only the frontier ones so we can build any 3*n* lattice lego-style. Do the following

```
lmchk dat --mino~dxmx=.001~xtol=.0001~nkill=10~maxit=200~11:34,49:72 --wsite~short~fn=site82 --wpos=pos
```

Note that only sites 11:34,49:72 move (compare to table above).

Command-line arguments to *lmchk* are documented on this page. Note that **~** separates arguments.

`--mino~dxmx=.001~xtol=.0001~nkill=10~maxit=200~11:34,49:72`

tells*lmchk*to move a selected list of sites to minimize the objective function.

The last argument selects the list (see here for syntax of integer lists).

**11:34,49:72**encompasses all the sites in the frontier layers (see table above).

**dxmx**,**xtol**,**nkill**and**maxit**are explained on the command line arguments web page.`--wsite~short~fn=site82`

tells*lmchk*to write*site*file*site82.dat*.`--wpos=pos`

tells*lmchk*to write site positions to file*pos.dat*.

Near the end you should find the following:

```
minimized: fovl = 6.96742e-6 <ovlp> = 13.8% max ovlp = 19.6%
```

The maximum overlap is reduced to 19.6% from an initial 36%. Since there is a heavy energy penalty for atoms being too close together, this is a much better starting point than the unrelaxed lattice.

*lmchk* wrote the unrelaxed and relaxed positions to files *pos0.dat* and *pos.dat* respectively. To see how much atoms moved, do

```
mcx pos.dat pos0.dat -- > dpos3.dat
```

Displacements are significant, on the order of 10% of the bond length. Note that if *all* the sites were included in the list, the maximum overlap would fall to 17.3% — not much lower than the constrained case. This provides a rough justification for constraining the relaxation to frontier atoms.

#### 2b. Nb/Ni(3*n*)+1/Nb superlattices

Consider the following possibibilities stacking the Nb/Ni4/Nb superlattice. (Here A and B refer to the two Nb planes, and 1,2,3 refer to the three Ni planes). AB1231AB, AB2312AB, AB3123AB. It turns out the AB1231AB has the largest overlap (38.1%), while AB2312AB and AB3123AB both have 36% overlaps.

To build the AB3123AB stack, copy the following box into file *batch3123ab*

```
file@site14c.ni@dpos=0,0,-dz+dp3/2
file@site14a.ni@dpos=0,0,-dz+dp3/2
file@site14b.ni@dpos=0,0,-dz+dp3/2
file@site14c.ni@dpos=0,0,-dz+dp3/2
file@site20.nb@dpos=0,0,dp3
```

Then do

```
lmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.00 --stack~batch=batch3123ab~wsite@fn=sitex
cp sitex.nb sitein.dat
blm --nfile --rdsite=sitein --rdspec~fn=spec.nb~fn=spec.ni --usefiler --nosort --shorten=no --ctrl=ctrl dat
lmchk dat --pr20 --wpos=pos0
```

You should find the maximum sphere overlap is 36%. We will use the ideal AB3123AB as a starting point.

Do the crude relaxation by minimizing sphere overlaps

```
lmchk dat --mino~dxmx=.001~xtol=.0001~nkill=10~maxit=200~11:34,63:86 --wsite~short~fn=site82 --wpos=pos
mcx pos.dat pos0.dat -- > dpos4.dat
```

Near the end you should find the following:

```
minimized: fovl = 6.05759e-6 <ovlp> = 13.5% max ovlp = 19.8%
```

The initial maximum overlap is a bit larger than the AB123AB of Section 2a, but the minimizer nearly equalized the two.

#### 2c. Nb/Ni(3*n*)+2/Nb superlattices

Use the same patterning for the Nb/Ni5/Nb superlattice as for the Nb/Ni4/Nb superlattice. To build the AB12312AB stack, copy the following box into file *batch12312ab*

```
file@site14a.ni@dpos=0,0,-dz+dp3/2
file@site14b.ni@dpos=0,0,-dz+dp3/2
file@site14c.ni@dpos=0,0,-dz+dp3/2
file@site14a.ni@dpos=0,0,-dz+dp3/2
file@site14b.ni@dpos=0,0,-dz+dp3/2
file@site20.nb@dpos=0,0,dp3
```

Then do

```
lmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.00 --stack~batch=batch12312ab~wsite@fn=sitex
cp sitex.nb sitein.dat
blm --nfile --rdsite=sitein --rdspec~fn=spec.nb~fn=spec.ni --usefiler --nosort --shorten=no --ctrl=ctrl dat
lmchk dat --pr20 --wpos=pos0
```

You should find the maximum sphere overlap is 36%. We will use the ideal AB3123AB as a starting point, as the overlaps are better than AB23123AB (37.5%) or AB31231AB (38.1%).

Minimize sphere overlaps

```
lmchk dat --mino~dxmx=.001~xtol=.0001~nkill=10~maxit=200~11:34,77:100 --wsite~short~fn=site82 --wpos=pos
mcx -f3f14.8 pos.dat pos0.dat -- > dpos5.dat
```

### 3. Relaxing the Nb(2)/Ni(3)/Nb(2) superlattice

If the files below are not already present in your working directory, set up this section with

```
cp ~/lmf/fp/test/nbni/{ctrl.nb,site20.nb,spec.nb,site42.ni,spec.ni,fpspec.nb,fpspec.ni,dpos3.dat,site14[abc].ni,site10[ab].nb,nb2ni4nb2/dpos4.dat,nb2ni4nb2/batch3123ab} .
```

#### 3a. Relaxing the Nb(2)/Ni(3)/Nb(2) superlattice

The previous section supplied a suitable structure to perform a density-functional calculation of the Nb(2)/Ni(3)/Nb(2) superlattice, assuming for the Nb/Ni spacing an average of the Nb-Nb spacing and Ni-Ni spacing. It also provided a rough estimate for how sites relax around their ideal geometries (*dpos3.dat*).

Here we will use total energy minimization to properly relax the structure by minimizing forces. At the same time the Nb/Ni interplanar spacing should not be assumed, but calculated. We also will differentiate the frontier Nb and Ni species from their bulk counterparts. It is not crucial here but it will be when transport is calculated in the ASA.

```
lmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.00 --stack~file@site42.ni@dpos=0,0,-dz+dp3/2~file@site20.nb@dpos=0,0,dp3~setlbl@spid=Nbx@targ=11:20,nbas-19:nbas-10~setlbl@spid=Nix@targ=21:34,nbas-33:nbas-20~pad=0,0,dp3~wpos@fn=pos0~wsite@fn=sitex
cp sitex.nb sitein.nbni
blm --rdspec~fn=fpspec.nb~fn=fpspec.ni --usefiler --rdpos=dpos3.dat --autobas~loc=5 --rsham~el=-.1,-0.9~rscut=1.6 --dhftol=45 --molstat~xtol=0~gtol=0.01 --nfile --rdsite --omax=-.01 --nosort --shorten=no --nit=20 --nk=4,5,2 --gmax=10 --mixb=.2 --frzwf --wsite~ib=1:10,35:48,73:82~rlx=000~long~fn=site82 nbni
cat actrl.nbni | sed s/NFORCE=0/NFORCE=20/ > ctrl.nbni
```

The initial flags are similar to Section 2, with some addition

`-vfile=20 ctrl.nb`

tells*lmscell*to build the stack starting from the reconstructed 20-atom Nb cell as a base.The Nb-Nb spacing is and the Ni-Ni spacing is . is half this difference.

**dp3**isn’t used now, but later it will become an ajustment to the spacings between frontier layers.`--stack`

invokes the superlattice editor. It is invoked in batch mode with**~**separating commands and does the following:`~file@site42.ni@dpos=0,0,-dz+dp3/2`

stacks the 42-atom reconstructed Ni layer on top of Nb. Absent the modifier**dpos**, the first Ni plane would be placed at the boundary of the 20-atom reconstructed Nb layer, which lies at a height above the topmost Nb layer. Rigidly shifting the 42 Ni atoms by causes the Nb-Ni spacing to take the average of the Nb-Nb and Ni-Ni spacings, plus a shift**p3/2**.`file@site20.nb@dpos=0,0,dp3`

stacks the 20 atom reconstructed Nb layer on top of the stack. Since the Ni was shifted by , the second Ni-Nb is automatically at the average of Ni-Ni and Nb-Nb spacings.**p3**takes into account the shift from changing the entire stack size.`~wsite@fn=sitex`

writes the file*sitex.nb*, which is a*site*file of the superlattice.`setlbl@spid=Nbx@targ=11:20,nbas-19:nbas-10`

reassigns labels for frontier Nb sites to Nbx;`setlbl@spid=Nix@targ=21:34,nbas-33:nbas-20`

does the same for Ni sites.`pad=0,0,dp3`

adds**dp3**to the stack size.`wpos@fn=pos0`

writes the site positions to file*pos0.nb*. These are site positions for the ideal lattice.`wsite@fn=sitex`

writes file*sitex.nb*

The entire stack height can be modified by **p3**. This is distributed in the cell by shifting sites after the first interface by **p3/2** and sites after the second by **p3**.

Command-line arguments to *blm* are documented on this page and more details are provided there.

`--rdspec~fn=fpspec.nb~fn=fpspec.ni`

read species data from files*fpspec.nb*and*fpspec.ni*)`--usefiler`

tells*blm*to retain any sphere radii it reads from species files (in this case*fpspec.nb*and*fpspec.ni*)`--rdpos=dpos3.dat`

shifts positions of each site by the contents of*dpos3.dat*(Section 2’s simple estimate for relaxation by minimizing sphere overlaps)`--autobas~loc=5`

Modifies AUTOBAS_LOC to designate low-lying local orbitals as of the conventional type.`--rsham~el=-.3,-1.2`

Adds category STR for real-space strux. This will play a role when the JPO basis is available.`--dhftol=45`

Adds convergence criteria (ITER_DHFTOL) for self-consistency. This allows charge convergence tolerance to be set loosely (see below)`--molstat~xtol=0~gtol=0.01`

Adds DYN category for molecular statics, using gradients only as a check for convergence`--nfile`

was explained above`--rdsite`

tells*blm*to read structural data from*sitein.nbni*`--omax=-.01`

select sphere radii so that the maximum sphere overlap is −1%`--nosort`

suppress re-ordering of site positions`--shorten=no`

suppress automatic shortening of site positions by adding multiples of lattice translation vectors.`--nit=20`

variable controlling the number of iterations stop at if convergence is not reached. If modifies the input file via tag**EXPRESS_nit**, an alias for ITER_NIT`--nk=4,5,2`

*k*mesh for BZ_NKABC Brillouin zone integration. This choice spaces*k*along the three axes approximately equally.`--gmax=10`

specify mesh density cutoff HAM_GMAX.*lmfa*can find it for you (which is how the value was arrived at).`--mixb=.1`

Specify the charge mixing parameter**beta**.**beta=0.2**was found empirically to be a more efficient choice than the default**beta=0.3**.`--frzwf`

Add a variable**frzwf**which if set, suppresses updating partial waves through the tag HAM_FRZWF`--wsite~ib=1:10,35:48,73:82~rlx=000~long~fn=site82`

writes a file*site82.nbni*sets the rlx flag to 000 for all sites that are not frontier layers. This restricts relaxation to frontier atoms.

*blm* makes *actrl.nbni* and *site82.nbni* in a form suitable for *lmf*. The last command `cat actrl.nbni | sed s/NFORCE=0/NFORCE=20/ > ctrl.nbni`

creates the *ctrl* we will use.

We will use this to relax the lattice via self-consistent LDA calculations. Self-consistency is not the purpose, but rather to obtain lattice displacments needed by the ASA. It will have bo be repeated for transport in the ASA, but the ASA cannot be used to relax the lattice. The bands generated here can also serve as a benchmark to validate the quality of the ASA.

It is always a good idea to confirm that sphere radii and overlaps are sensible. The input file has one radius Nb and one for Ni and two other radii for the frontier species, all of which are reasonably similar. Click on the line below for further details.

- This invocation of
*blm*makes a*ctrl.nbni*to be used in self-consistent LDA calculations. Some care now has to be taken concerning chemical attributes in addition to structural ones, e.g. choice of augmentation sphere radii.*blm*found good choices for the sphere radii automatically: they were made as large as possible subject to the constraint that overlap between any pair not exceed −1%.

To confirm this, try`lmchk ctrl.nbni -vfile=82`

These lines should appear in the output:`Cell volume= 7949.71676 Sum of sphere volumes= 4542.64918 (0.57142) ... OVMIN, 1404 pairs: fovl = 5.77296e-14 <ovlp> = 0.6% max ovlp = 0.6%`

See this page for general considerations in how to select sphere radii. Here we would like to make them as large as we can while avoiding inaccuracies from geometry violation when spheres overlap. As that page notes

*lmf*can tolerate a limited amount of overlap. By starting with and initial −1%, we may expect some positive overlaps to appear in the relaxation process. (The maximum overlap, 0.6%, comes from Ni-Ni pairs whose radii were both constrained. Overlaps between pairs where at least one radius can vary do not exceed −1%.) As we will see the maximum overlap will ultimately be about 2% when the system is finally relaxed, which is acceptable. Also the average sphere packing is 57%, which is quite healthy, so the basis should be accurate. **--molstat**adds provision in the input file for molecular statics, to relax the lattice to equilibrium positions. There is a tutorial for a simple case, which you might want to consult before continuing with this one. This is a restricted relaxation for reasons explained above Only the interfacial atoms will be free to move (information is fed to*lmf*via*site82.ni*, the**rlx**column).*blm*sets**rlx**to zero for all sites*not*on a frontier layer.- The
**k**mesh,**4,5,2**specifies the division of the reciprocal lattice into microcells. The requisite mesh density cannot be reliably predicted in advance: some trial and error is needed. This mesh was found to accurate be fairly accurate, certainly enough that the forces will be reliable. We are using tetrahedron integration and as the calculation prints out the “Bloechl correction,” which gives a measure of the mesh quality. (It is about 2 mRy/au in this case, which for 82 atoms is quite good.) The cell has 82 atoms: it is a bit expensive so we want keep the mesh relatively light.**4,5,2**was selected because it was found empirically to be dense enough to be accurate, and also space lengths of the edges of parallelpipeds as evenly as possible. You can do determine this by hand, but*lmf*can do it for you. Temporarily insert a line`nkabc= -50`

near the beginning of the EXPRESS category, at any rate before the existing tag

**nkabc**. (Always the first tag is read; later ones are ignored.) Then do`lmf ctrl.nbni -vfile=82 --quit=ham`

and look for this line:`BZMESH: 22 irreducible QP from 40 ( 4 5 2 ) shift= F F F`

If you did this check, be sure to remove the temporary line

`nkabc= -50`

. `--rsham`

was included anticipating that we will use Pashov’s redesign of*lmf*.*lmf*will use alternative set of critical blocks that are significantly faster than the ones they replace, and moreover*lmf*can be run on multicore machines with both OMP and BLAS threading as well as MPI parallelization over**k**. This instruction is put there in anticipation of the Jigsaw Puzzle Orbital basis. Even though it does not work as of this writing, when it does work no change to the input file will be needed. In particular these lines supply information for the JPO basis. At present they supply Hankel energies EL. (The JPO basis does not permit an*l*dependent EH.)`% ifdef rsham>1 STR RMAXA={rscut} ENV[MODE=1 NEL=2 EL={el1} {el2}] % endif`

We use

**el1=−0.3**and**el2=−-1.2**.

We need a starting density and basis set both supplied by *lmfa*:

```
cp site82.nbni site1.nbni
lmfa -vfile=1 ctrl.nbni
cp basp0.nbni basp.nbni
lmfa -vfile=1 ctrl.nbni
```

`-vfile=1`

tells *lmfa* to read site file *site1.nbni*, as can seen by inspecting *ctrl.nbni*. (We could have used `-vfile=82`

and read from *site82.nbni*, but the calculations will need to be repeated for different **p3**, *site1.nbni* will be designated as the file for **p3=0**.) *lmfa* writes the basis file to *basp0.nbni* and information about the atomic density to *atm.nbni*. It has to be run a second time once information about the local orbitals is found from the first pass (this can redistribute the weight between core and valence density).

We are ready to run *lmf*. Initially we want to make the system self-consistent without relaxing any atoms. Then we can see how large the forces are with our simple but rude estimate for lattice displacements, and to monitor the precision of the **k** point integration.

With 82 atoms, running in serial mode is not a viable option. Exactly how you run *lmf* will depend on the hardware available to you. *lmf* will parallelize over **k** points, but we need to know how many. Extract that information this way:

```
lmf ctrl.nbni -vfile=1 --quit=ham | grep BZMESH
```

There are 22 irreducible QP. As the calculation is a bit expensive, we will make the initial runs with a smaller *k*-mesh using 4×4×1 divisions in place of 4×5×2. The difference is not all that large, and we can tune the *k* mesh more finely near convergence, if needed.

```
lmf ctrl.nbni -vnk1=4 -vnk2=4 -vnk3=1 -vfile=1 --quit=ham | grep BZMESH
```

About parallelization: suppose you have a single node machine with 16 cores, and *lmf* is compiled with Intel’s MKL library (this tutorial was written using such a machine). You could do either one of these:

```
mpix -n 10 lmf -vbeta=.1 -vnit=30 -vnk1=4 -vnk2=4 -vnk3=1 -vfile=1 --v8 --tbeh -vrsham=2 ctrl.nbni > out
env OMP_NUM_THREADS=2 MPI_NUM_THREADS=2 | mpix -n 10 lmf -vbeta=.1 -vnit=30 -vnk1=4 -vnk2=4 -vnk3=1 -vfile=1 --v8 --tbeh -vrsham=2 ctrl.nbni > out
```

The latter uses 20 threads which is more than 16, but the extra cost of swapping threads is made up for making better use of the 16 available cores. If your machine has 32 cores you can do `env OMP_NUM_THREADS=3 MPI_NUM_THREADS=3`

, and so on. You can also invoke MPI across nodes.

*Caution:* When this tutorial was written the `--v8`

option did not compute the forces correctly, so traditional algorithms are used:

```
mpix -n 10 lmf -vbeta=.1 -vnit=30 -vnk1=4 -vnk2=4 -vnk3=1 -vfile=1 ctrl.nbni > out
```

Before relaxing the structure, inspect the output file. Note in particular:

Do

`grep more= out`

to see the change**RMS DQ**in output − input charge. As the iterations proceeded it decreased to a relatively small value, less than 10^{−4}.Also inspect

*save.nbni*and observe that the Harris-Foulkes and Kohn-Sham total energies are near each other in the final iteration.These results are both good indications that the density is close to self-consistent.

- The forces are not small, but they are not massive either:
`grep 'Maximum Harris force' out`

The maximum force should be about 67 mRy/a.u. This is a measure of the reliability of the crude relaxation by minimizing sphere overlaps.

- The relatively crude 4×4×1
*k*mesh isn’t bad either. Check out the change to the band sum evaluated with Bloechl weights instead of standard tetrahedron weights:`grep 'Bloechl correction' out`

1 mRy for 82 atoms is quite modest; it means that forces will not be strongly affected by a finer

*k*mesh.

Run *lmf* again using the (nearly) self-consistent *rst.nbni* file, this time allowing frontier atoms to relax.

```
rm -f mixm.nbni
mpix -n 10 lmf -vfile=1 -vnk1=4 -vnk2=4 -vnk3=1 -vfrzwf=t -vconv=1e-4 -vconvc=1e-3 -vminx=5 --wpos=pos36 ctrl.nbni > out1
```

This uses a Broyden algorithm to relax the structure for five iterations. It writes the positions file for what would be iteration six into *pos36.nbni*. ‘6’ designates iteration 6, and ‘3’ designates three ML of Ni. Later relaxations of of the

`-vconv=1e-4`

: Set a tolerance rougher than the default for convergence in the energy before shifting atoms`-vconvc=1e-3`

: Set a tolerance rougher than the default for convergence in the charge density before shifting atoms. This is further explained below.`-vfrzwf=t`

: suppresses updating of the shape of the partial waves. The constraint is not essential, but it stabilizes convergence and most of the change is already realized in the initial self-consistent cycle which did allow the potential to float. (Keeping it frozen to the shpae of the free atom is essentially what a pseudopotential does.) Since the atomic relaxation will only slightly modify the this potential, it is safer and more stable with little loss of accuracy. Also without this switch the forces are not quite exact. Finally the error it can always be checked after the relaxation by letting the potential float once more.`-vminx=5`

: Perform at most 5 relaxation steps. It gets pretty close to the minimum energy point; a second pass will be taken to refine the positions.`-wpos=pos36`

: Writes the position file after every move, including the last one.

Try `grep RELAX out1`

to see how the gradient evolved in successive relaxation steps:

```
RELAX: (warning) failed to read Hessian; set to unity
RELAX: completed Broyden iter 1; max shift=-0.0100 |g|=0.253
RELAX: completed Broyden iter 2; max shift=0.01000 |g|=0.169
RELAX: completed Broyden iter 3; max shift=-0.0100 |g|=0.129
RELAX: completed Broyden iter 4; max shift=0.01000 |g|=0.0589
RELAX: completed Broyden iter 5; max shift=0.00493 |g|=0.0235
```

**|g|** is the length of the 144 component force “vector”; it should drop from an initial value **0.237** to **0.0275** in the 5^{th} iteration.

The total energy also drops significantly as `egrep -E '^[cx]' out1`

shows:

```
c file=1 nk1=4 nk2=4 nk3=1 frzwf=1 conv=1e-4 convc=.001 minx=5 ehf=-432804.2830013 ehk=-432804.2832668
x file=1 nk1=4 nk2=4 nk3=1 frzwf=1 conv=1e-4 convc=.001 minx=5 ehf=-432804.3377062 ehk=-432804.3382388
x file=1 nk1=4 nk2=4 nk3=1 frzwf=1 conv=1e-4 convc=.001 minx=5 ehf=-432804.3724307 ehk=-432804.3727466
x file=1 nk1=4 nk2=4 nk3=1 frzwf=1 conv=1e-4 convc=.001 minx=5 ehf=-432804.4007881 ehk=-432804.4011704
x file=1 nk1=4 nk2=4 nk3=1 frzwf=1 conv=1e-4 convc=.001 minx=5 ehf=-432804.4154372 ehk=-432804.4152356
```

- Concerning the relaxation algorithm
- Parameters in
**DYN_MSTAT**govern*lmf*’s relaxation algorithm. Tokens are explained in the input file; you can also get a quick summary with

`lmf --input | grep -A 3 DYN_MSTAT`

.**DYN_MSTAT_MODE=6**invokes a Broyden method.- Convergence is checked only on the forces (gradients) since
**XTOL=0**. **STEP=.010**limits the steps size to a relatively conservative displacement. It could be made larger, but the disruption to the potential becomes greater and it is safer not keep displacements fairly small.- Because of the constraints (see discussion of
**rlx**above) only frontier atoms are allowed to move. There are 10 and 14 Nb and Ni atoms per layer; thus there are 48×3 = 144 degrees of freedom in the relaxation.

- Iterations to convergence
*lmf*has two nested loops: the inner loop converged the charge at the given atomic positions, the outer loop relaxes the nuclei. This disrupts the potential, so that it has to be converged again at the updated positions. You can see this by doing`egrep -E 'DQ|Maximum Harris|more=F' out5`

.

This extracts lines with RMS change in charge**DQ**, the maximum value of the Harris force, and summary lines terminating The RMS change in the charge**DQ**drops steadily until the last iteration at given positions (**more=F**). the inner cycle. Atoms are then moved, the potential is disrupted and it the density has to be converged again.`Maximum Harris force = 63.2 mRy/au (site 64) Max eval correction = 0.9 ← iter 1 rlx 0 Maximum Harris force = 63.2 mRy/au (site 64) Max eval correction = 0.8 ← iter 2 rlx 0 Maximum Harris force = 233 mRy/au (site 64) Max eval correction = 355.5 ← iter 1 rlx 1 ... Maximum Harris force = 36.8 mRy/au (site 64) Max eval correction = 111.8 ← iter 20 rlx 1 ...`

Lines show maximum force on any site and an estimate for the correction to this force.

- When that correction is on the order of the force itself, usually the density is converged enough that there is little advantage to driving self-consistency much farther. Correction is small in first two lines above because the density was already nearly self-consistent from the initial run of
*lmf*. - In the third line the atoms had been moved; it is the first iteration of the inner loop. The correction is huge and the forces not very meaningful.
- The next line [iter 20 rlx 1] is the last step in the inner loop (iterations are limited to 20, even if self-consistency has not been reached). The density is nevertheless reasonably converged because the Harris-Foulkes and Kohn-Sham forces are nearly equal. But when the convergence criteria were not met, the summary line e.g.
`x file=1 frzwf=1 conv=1e-4 convc=.001 minx=5 beta=.2 ehf=-432804.3830901 ehk=-432804.3833868`

has an

**x**in the first column, rather than**c**. - The energy and charge convergence tolerances (
**conv=1e-4 convc=.001**) are set fairly loosely. This is safe because an extra condition is imposed (note**dhftol**in the*ctrl*file). The correction term to the Harris-Foulkes forces must fall below this value. We don’t care so much whether the energy is well converged, but that the forces are well described. The forces are good enough for a rough relaxation, which since it will be refined is all we need at this stage. As it happens in this case the relaxation proceeds smoothly anyway.

- When that correction is on the order of the force itself, usually the density is converged enough that there is little advantage to driving self-consistency much farther. Correction is small in first two lines above because the density was already nearly self-consistent from the initial run of

The maximum force on an individual atom can be seen from `grep 'Maximum Harris force' out1`

```
Maximum Harris force = 67.1 mRy/au (site 64) Max eval correction = 0.3 ← before relaxation
...
Maximum Harris force = 26.7 mRy/au (site 70) Max eval correction = 144.1 ← after relaxation, iteration 5
```

The largest forces by far are on atoms constrained not to relax (largest force is on site 47, middle Ni layer), which will not be improved by further iteration. You can confirm this by running the *lmf* once again at the positions for what would have been iteration 6:

```
rm -f mixm.nbni hssn.nbni
mpix -n 10 lmf ctrl.nbni -vfile=$file -vnk1=4 -vnk2=4 -vnk3=1 -vfrzwf=t -vconv=1e-4 -vconvc=1e-3 --rpos=pos36 --rs=11,2,1 -vnit=20
```

The energy drops another 2 mRy, small compared to prior changes in energy. Thus the positions in file _pos36.nbni{:.path} are nearly optimum within the constraints. You can also check the **k** convergence (one iteration is enough)

```
mpix -n 10 lmf ctrl.nbni -vfile=$file -vfrzwf=t -vconv=1e-4 -vconvc=1e-3 --rpos=pos36 --rs=12,0 -vnit=1
```

and you should see that the charge density and total energy are minimally affected.

As noted, frontier atoms continue to have significant forces. If need be we could reduce this force by including more layers in the relaxation. But, this complicates our ability to assemble stacks “lego” style, and since the real interface is probably rather non-ideal and almost certainly different from this one, it isn’t worth the trouble to reduce the residual forces.

##### Optimizing the Nb/Ni layer spacing

There is no good justification to assume that the Nb-Ni spacing is . The instructions above can vary spacing via **p3**: spacing at each interface shrinks by **p3/2**.

Repeat the `lmscell ... / cp ... / blm ...`

block of Section 3a but assign **dp3** as **-vdp3=-.04** and designate *site2.nbni* as the site file for **p3=−0.04** : `cp site82.nbni site2.nbni`

. Do not repeat the last step (`cat actrl.nbni | sed s/NFORCE=0/NFORCE=20/ > ctrl.nbni`

), but preserve *ctrl.nbni* from Section 3a to minimize cancellation of errors. Sphere radii could have changed because of the different relaxation process, though they don’t in practice in this case. (If you no longer have this file, copy it from *~/lmf/fp/test/nbni/ctrl.nbni* .)

Perform the same procedure as before (substituting **-vfile=2** for **-vfile=1**) to relax the *site2.nbni* setup.

Make a new *site* file for **p3=−0.08**, now designating *site3.nbni* as the site file, and relax the *site3.nbni* setup; to the same for **p3=−0.12**, designating *site4.nbni* as the site file

For the four relaxations you should find the following total energies

```
0.00 -432804.4154372 ← site1
-.04 -432804.4578181 ← site2
-.08 -432804.4852606 ← site3
-.12 -432804.4949157 ← site4
```

The energy gained by allowing **p3** to vary is comparable to the internal site relaxation. To estimate the minimum, copy the first two columns into a file *dat* and fit a second order polynomial to it:

```
echo q | pfit 2 dat
```

The fit puts the minimum at &minus.125, which is close enough to the *site4* structure.

*Sanity checks:* the following files should be nearly identical to those in the distribution, *~/lmf/pgf/test/nbni-111* :

```
diff pos36.nbni ~/lmf/pgf/test/nbni-111/pos36.nbni
```

##### Preparing a site file for the ASA

```
cp ~/lmf/fp/test/nbni/{ctrl.nb,site20.nb,fpspec.nb,site42.ni,fpspec.ni} .
```

The ASA calculation of transmission will need site information to proceed. It requires a site file without any relaxations (created below as *site30.nbni*), and information about the relaxed positions. This is file *pos36.nbni* made in the previous section.

The steps below follow the setup of Section 3a but produce *site30.nbni* (the unrelaxed site file corresponding to **p3=−.12**). ‘3’ marks three ML of Ni, and ‘0’ designates no relaxation.

```
lmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.12 --stack~file@site42.ni@dpos=0,0,-dz+dp3/2~file@site20.nb@dpos=0,0,dp3~setlbl@spid=Nbx@targ=11:20,nbas-19:nbas-10~setlbl@spid=Nix@targ=21:34,nbas-33:nbas-20~pad=0,0,dp3~wpos@fn=pos0~wsite@fn=sitex
cp sitex.nb sitein.nbni
blm --rdspec~fn=fpspec.nb~fn=fpspec.ni --usefiler --autobas~loc=5 --rsham~el=-.1,-0.9~rscut=1.6 --dhftol=45 --molstat~xtol=0~gtol=0.01 --nfile --rdsite --omax=-.01 --nosort --shorten=no --nit=20 --nk=4,5,2 --gmax=10 --mixb=.2 --frzwf --wsite~ib=1:10,35:48,73:82~rlx=000~long~fn=site30 nbni
```

*Sanity checks:* the following files should be nearly identical to those in the distribution, *~/lmf/pgf/test/nbni-111* :

```
```

#### 3b. Relaxing the Nb(2)/Ni(4)/Nb(2) superlattice

The files below are needed to set up this section

```
cp ~/lmf/fp/test/nbni/{ctrl.nb,ctrl.nbni,site20.nb,spec.nb,site42.ni,spec.ni,fpspec.nb,fpspec.ni,site14[abc].ni,site10[ab].nb,nb2ni4nb2/dpos4.dat,nb2ni4nb2/batch3123ab} .
```

We can reasonably assume **p3** will be the same for the Nb(2)/Ni(4)/Nb(2) superlattice as for the Nb(2)/Ni(3)/Nb(2) superlattice. Make a new *site* file for **p3=−0.08**, now designating *site43.nbni* as the site file

```
clmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.12 --stack~batch=batch3123ab~setlbl@spid=Nbx@targ=11:20,nbas-19:nbas-10~setlbl@spid=Nix@targ=21:34,nbas-33:nbas-20~pad=0,0,dp3~wpos@fn=pos0~wsite@fn=sitex
cp sitex.nb sitein.nbni
blm --rdspec~fn=fpspec.nb~fn=fpspec.ni --usefiler --rdpos=dpos4.dat --autobas~loc=5 --rsham~el=-.1,-0.9~rscut=1.6 --dhftol=45 --molstat~xtol=0~gtol=0.01 --nfile --rdsite --omax=-.01 --nosort --shorten=no --nit=20 --nk=4,5,2 --gmax=10 --mixb=.2 --frzwf --wsite~ib=1:10,35,96-9:96,35:62~rlx=000~long~fn=site82 nbni
cp site82.nbni site43.nbni
```

We can use the same *ctrl.nbni* as for the Nb(2)/Ni(3)/Nb(2) superlattice. But because the cell is becoming long, it is safer to tighten tolerances on the Ewald sums. Starting from the existing *ctrl* file, do

```
cat ctrl.nbni | awk '{ if ($3 == "rsham>1") {print " TOL=1d-12"; print "EWALD TOL=1d-12"}; print}' > ctrl.nbni2
cp ctrl.nbni2 ctrl.nbni
```

Following the relaxation procedure once again, relax the *site43.nbni* setup. This time write the positions file into *pos46.nbni*. ‘6’ designates iteration 6, and ‘4’ designates four ML of Ni.

#### 3c. Relaxing the Nb(2)/Ni(5)/Nb(2) superlattice

The files below are needed to set up this section

```
cp ~/lmf/fp/test/nbni/{ctrl.nb,ctrl.nbni,site20.nb,spec.nb,site42.ni,spec.ni,fpspec.nb,fpspec.ni,site14[abc].ni,site10[ab].nb,nb2ni5nb2/dpos5.dat,nb2ni5nb2/batch12312ab}
```

Use the same **p3** for the Nb(2)/Ni(5)/Nb(2) superlattice. Make a new *site* file for **p3=−0.08**, now designating *site53.nbni* as the site file

```
lmscell -vfile=20 ctrl.nb -vdz0='1/sqrt(2)' -vdz=dz0/2-.61/2 -vdp3=-.12 --stack~batch=batch12312ab~setlbl@spid=Nbx@targ=11:20,nbas-19:nbas-10~setlbl@spid=Nix@targ=21:34,nbas-33:nbas-20~pad=0,0,dp3~wpos@fn=pos0~wsite@fn=sitex
cp sitex.nb sitein.nbni
blm --rdspec~fn=fpspec.nb~fn=fpspec.ni --usefiler --rdpos=dpos5.dat --autobas~loc=5 --rsham~el=-.1,-0.9~rscut=1.6 --dhftol=45 --molstat~xtol=0~gtol=0.01 --nfile --rdsite --omax=-.01 --nosort --shorten=no --nit=20 --nk=4,5,2 --gmax=10 --mixb=.2 --frzwf --wsite~ib=1:10,35:76,110-9:110~rlx=000~long~fn=site82 nbni
cp site82.nbni site54.nbni
```

Use the same *ctrl.nbni* as for the Nb(2)/Ni(4)/Nb(2) superlattice.

```
cat ctrl.nbni | awk '{ if ($3 == "rsham>1") {print " TOL=1d-12"; print "EWALD TOL=1d-12"}; print}' > ctrl.nbni2
cp ctrl.nbni2 ctrl.nbni
```

Following the relaxation procedure once again, relax the *site53.nbni* setup.

### Footnotes and references

^{1} Dimitar Pashov, Swagata Acharya, Walter R. L. Lambrecht, Jerome Jackson, Kirill D. Belashchenko, Athanasios Chantis, Francois Jamet, Mark van Schilfgaarde, *Questaal: a package of electronic structure methods based on the linear muffin-tin orbital technique*, Comp. Phys. Comm. **249**, 107065 (2020).