Simulation of STM images¶

Author: Jorge Pilo (ICN2)

This tutorial provides some guided examples for the use of the tools available to simulate STM and STS images with Siesta. The important tools are two: the plstm executable and the stm executable (the latter one is also called ol-stm in this tutorial and it is also known with the more appropriate name wfs2ldos). Before starting it is important to read here to get an overview of the tools.

Note

For some background and a gallery of images, see the presentation by Jorge Pilo.

STM image of a benzene molecule¶

In this exercise we will obtain the STM image of a benzene molecule, We will compare the plstm and ol-stm ways and see how these tools can complement each other.

Options in the Siesta run¶

First of all, enter the folder STM_Benzene and explore the benzene.fdf file. The keywords we should focus on are:

The %block LocalDensityOfStates that asks Siesta to calculate the LDOS in an energy range.
write-denchar that asks Siesta to generate two files called SystemLabel.PLD and SystemLabel.DIM where information needed to plot the charge density and/or wavefunctions in real space is dumped.
COOP.Write the writes information on the wavefunctions at all k-points used to sample the BZ during self-consistency.
SaveElectrostaticPotential to save the Hartree potential
WFS.Energy.Min and WFS.Energy.Max to limit the number of wave functions plotted to the ones corresponding to an eigenvalue on this energy range.

Note

Some considerations about the structure: we can study only planes perpendicular to the z direction, therefore the structure must be prepared accordingly.

Run the calculation with:

siesta < benzene.fdf > benzene.out

STM image from direct LDOS calculation: plstm¶

We first use psltm to obtain a first STM image. For this tool only the %block LocalDensityOfStates was necessary in the siesta input and only the benzene.LDOS file is needed. This file contains LDOS(x,y,z) in the same grid used by Siesta in that calculation (i.e., it will depend on the MeshCutoff used). This LDOS function is periodic. Create a new folder and copy or link there the relevant file:

mkdir PLSTM
cd PLSTM
ln -sf ../benzene.LDOS ./benzene.LDOS

To get the STM image at a distance of 9.82 Bohr from z=0:

plstm -z 9.82 -o stm_9.82bohr benzene.LDOS

The “stm9.82bohr” is the output file that contains the STM data plottable with gnuplot. Open gnuplot (just type gnuplot) and enter the commands:

set pm3d
set pm3d map
set xrange [0:22.253958]
set yrange [0:20.845341]
set palette rgb 34,35,36
set terminal wxt size 440,440
splot 'stm_9.82bohr'

Congratulations! You got your first STM image! To exit gnuplot type quit.

More on the plstm will be cover in the next sections of this tutorial.

STM image from wavefunction information using ol-stm¶

The plstm program is very simple to use but it is limited in flexibility. For instance it does not allow to change the energy range or the grid for the calculation (unless you run again a new Siesta calculation). A more flaxible tool is the ol-stm program (called with the stm executable for historic reasons).

The ol-stm starts from the wavefunctions to first generate the LDOS, from where the STM/STS images can be obtained. Importantly, this code can optionally compute the values of the wavefunctions generated by Siesta on a reference plane and then extrapolate the value of these wavefunctions into the vacuum. In fact, for realistic experimental conditions the STM tip is relatively far from the surface. At this distance the accuracy of the Siesta calculation is not sufficiently good, because the basis set is short-ranged (and also because of other issues common to DFT approaches, including plane-wave simulations). To be able to extrapolate the value of the wavefunctions is a fairly important feature. More can be understood reading here.

Let’s go back to the folder of the siesta calculation cd ...

Located in the folder one can see the files generated by the Siesta run:

benzene.PLD
benzene.DIM
.ion
benzene.fullBZ.WFSX
benzene.VH

All these files are needed by the stm executable.

Warning

stm imposes to have the wavefunctions written in the SystemLabel.WFSX. Siesta returns this info in different files depending on the way used to obtain them. In our case we have the benzene.fullBZ.WFSX file name. This must be changed to benzene.WFSX.

ln -sf benzene.fullBZ.WFSX benzene.WFSX

The file benzene-stm.fdf was instead already present in the folder and it is the input for the stm executable. Explore this file. A part from the system information, some important keywords are:

STM.MinZ (real length): Defines the minimum value of the z-component of the 3D grid.
STM.MaxZ (real length): Defines the maximum value of the z-component of the 3D grid.
STM.RefZ (real length): Defines the z-component position of the reference plane from where the wavefunctions will be extrapolated into vacuum. The LDOS values for planes below the reference level are obtained from the original wavefunctions. Those above the reference level are obtained from projected wavefunction values.
STM.VacZ (real length): Defines the z-component position of the plane on which the zero of the Hartree potential in the vacuum region is taken. This should be a position at the center of the vacuum of the slab. You might want to check (better: plot) the V(z) values printed in the benzene.v_ave_z file. This file is produced by stm itself. Therefore you need a first guess and maybe rerun the calculation in case this value was wrong.
STM.Emin (real energy): Lower value of the energy window (eV).
STM.Emax (real energy): Upper value of the energy window (eV).
STM.NumberPointsX (integer): Number of subdivision of the grid in the direction of the first lattice vector.
STM.NumberPointsY (integer): Number of subdivision of the grid in the direction of the second lattice vector
STM.NumberPointsZ (integer): Number of subdivision of the grid in the z-direction (which must be the direction of the third lattice vector)

More keywords can be explored in the reference manual for the stm (or wfs2ldos program).

We can now run

stm < benzene-stm.fdf > benzene-stm.out

The file generated are:

benzene.CELL_STRUCT: Cell structure parameters
benzene.STM.LDOS: In STM mode, it contains the value of the LDOS (all relevant spin components) at the grid defined by the user, in Siesta grid format.
benzene.v_ave_z: An ASCII file containing the V(z) profile of the electrostatic potential. Plotting the information in this file is useful to set appropriately STM.Vacz and STM.RefZ.

To obtain the STM data we are going to use the rho2xsf util. We need only the benzene.XV file (produced by any Siesta run and containing the structure) and the file with the LDOS. Better create a separate folder, and take into account that the rho2xsf recognizes only files with single extension, therefore:

mkdir STM
cd STM
ln -sf ../benzene.XV ./benzene.XV
ln -sf ../benzene.STM.LDOS ./benzene.LDOS

Type rho2xsf. The program opens an interactive shell and you have to insert the following info: [Maybe better to prepare an input file to represent the standard input stream]

Specify SystemLabel (or ’siesta’ if none): benzene
Would you use Bohr (B) or Ang (A) ? B
Enter reference point in Bohr : 0.0 0.0 0.0
Enter 1st spanning vector in Bohr : 22.253958475 0.000000000 0.000000000
measured from reference point in Bohr : 0.000000000 22.253958475
Enter 2nd spanning vector in Bohr : 0.000000000 20.845341211 0.000000000
measured from reference point in Bohr : 0.000000000 20.845341211
Enter 3rd spanning vector in Bohr : 0.000000000 0.000000000 37.794537555
measured from reference point in Bohr : 0.000000000 37.794537555
Enter number of grid points along three vectors 40 40 40
Add grid property (LDOS, RHO, ...; or BYE if none): LDOS
Add grid property (LDOS, RHO, ...; or BYE if none): BYE

The benzene.XSF file has been created and can be opened with the xcrysden program. Type xcrysden. Go to File -> Open Structure -> Open XSF, select the XSF file just generated and “OK”. Then go to Tool -> Data Grid. The program should already select the LDOS file so just “OK”. A new window is opened to select the isosurface. Insert 0.01 as the value of the isosurface and untick the “Display isosurface” on the top left of the window. Go to the “Plain #1” tab. “Select color basis” to “BLACK-BROWN-WHITE”, tick “display color-plane” and “Submit”. Pictures of the process are in the presentation mentioned at the beginning of this tutorial. You have the STM of the benzene.

Note

An alternative would be to create a .cube file where the STM image is represented. This is done using g2c_ng.

g2c_ng -s benzene.CELL_STRUCT  -g benzene.STM.LDOS

And the .cube file can be visualized using VESTA (not by default in the Siesta Mobile)

A variety of images of the benzene STM can be found in the slide presentation mentioned at the top of this tutorial.

Note

One can now use plstm to postprocess the benzene.STM.LDOS file and get the appropriate slices (either constant-current or constant-height). The results can be plotted with gnuplot.

STM images for a monolayer of Cr with Spin-Orbit coupling¶

We will now explore the possibility to plot STM images for magnetic systems.

Enter the folder STM_Cr_3atms

You can see that the Cr.fdf file contains the keywords for the initialization of the spin and to activate spin orbit coupling (if you are not familiar with these concepts look at the Spin-orbit coupling tutorial).

Run siesta, it will take about 5 minutes:

siesta < Cr.fdf > Cr.out

Is important to check that you are in good agree with the final magnetic moments direction, these are shown because we specified the keyword WriteMullikenPop 1 in the input file. We show you the final directions for this example (can be found in the “mulliken” section of the Cr.out):

ATOM       Charge    Mag.Mom.     Sx      Sy      Sz
   1       5.99983   4.93634      2.466   4.276   0.000
   2       6.00020   4.93714      2.470  -4.275  -0.000
   3       5.99997   4.93661     -4.937   0.003  -0.000
Total      18.00000   0.00406     -0.001   0.004  -0.000

It is important to say that these values are in cartesian coordinates.

Get the stm image with plstm:

mkdir PLSTM
cd PLSTM
ln -sf ../Cr.RHO ./Cr.RHO
plstm -z 4.46 -s s -X 2 -Y 2 -o stm_s_4.46bohr Cr.RHO

To notice:

We now used the .RHO file, in fact the plstm recognizes both LDOS and RHO, so they can be used interchangeably.
We used the spin option -s set to s to plot the total spin magnitude [Spin(up)-Spin(down)] instead of the total charge (q that is the default).
We used -X and -Y to ask to plot more than just one unit cell.

You can use the School_STM.gplot script to plot the result with gnuplot:

ln -sf ../School_STM.gplot ./School_STM.gplot
gnuplot School_STM.gplot

Run with no spin:

plstm -z 4.46 -s q -X 2 -Y 2 -o stm_q_4.46bohr Cr.RHO

and selecting just one spin component:

plstm -z 4.46 -s x -X 2 -Y 2 -o stm_x_4.46bohr Cr.RHO

Remember to change the file name inside the gnoplot script in order to plot the new files.

Images with magnetic tips¶

The use of the spin option of the plstm (-s) is useful to explore the spin dependence of the LDOS, however a scenario that is more close to an experimental set-up is the situation when a magnetic tip is used.

We propose two different vectors for the magnetic moment of the tip: V1=(0.5, 0.866, 0) and V2=(0.866, -0.5, 0). Important to say that these vectors have to be normalized. Type:

plstm -z 4.46 -v '0.5 0.866 0' -X 3 -Y 3 -o stmV14.46bohr Cr.RHO

And then:

plstm -z 4.46 -v '0.866 -0.5 0' -X 3 -Y 3 -o stmV24.46bohr Cr.RHO

Use the script School_STM.gplot to post processing both files, after changing the file name inside the script.

Look at the presentation at the beginning of this page for more insights on the results.

If we are not happy with the grid or energy we can use the procedure involving the ol-stm program explained in the previous section. More on that on the slides at the beginning of this page.

Conclusions¶

We presented some simple examples to show how to use the STM utils from Siesta. In a real life calculation it is important to be sure about the convergence of all the parameters and inputs that the job needs.

As we could see, there are different ways to porecess a STM image with SIESTA, and this overview is not sufficient to cover all the aspects of the topic. For more info you can check the user’s manual for the codes (links from here).

There are some articles in which the authors use the STM Tool from SIESTA. We invite you to read:

E. Cisternas and J.D. Correa: Theoretical reproduction of superstructures revealed by STM on bilayer graphene, Chemical Physics 409, 74-78 (2012). DOI: 10.1016/j.chemphys.2012.09.021.
Ana Sanz-Matías, et al.: Computational insight into the origin of unexpected contrast in chiral markers as revealed by STM, Nanoscale 10, 1680-1694 (2018). DOI: 10.1039/C7NR07395J.