This is an old revision of the document!
Table of Contents
Simulation of STM images for a graphene nanoribbon adsorbed on a metallic substrate
In this exercise we are going to see how we can generate Scanning Tunneling Microscope (STM) images from simulations. While this is a more elaborate example where we need additional scripts and techniques, the basic idea is simple: Imagine a STM cantilever scanning the surface in the constant current mode, the current is then proportional to the electron density, which is something we get from our simulation (remember the cube files?). Therefore, all we have to do is to find an isosurface of the electron density and plot this isosurfaces z-coordinate.
This exercise is based on ideas and material by Daniele Passerone.
Preparation
- On the server is a package for you to unpack (hohoho
), containing a number of input files. Run the following in a new and empty directory:
tar xf /users/tiziano/CHE437_ex7.tar.gz
- The scripts are contained in yet another python package:
pip install --user https://github.com/ltalirz/asetk/archive/master.zip
… and since you have setup the path variable in a previous exercise, you should now have the following new commands available:
stm.py
,cube-plot.py
,cp2k-sumbias.py
.
Geometry optimization
First take a look at the complete geometry in all.xyz
using vmd
. As you can see, this is a much larger geometry than what we had before and would therefore take a lot longer to optimize. This is why some tricks are used in geo_opt.inp
: For the nanoribbon (mol.xyz
), we use the Density Functional based Tight-Binding (DFTB) method instead of GPW and couple the systems together by Force Fields (FIST
). Look at the CP2K input file geo_opt.inp
to get an idea how that looks like and then run the simulation.
Like with the geometry optimization in the previous examples, you now should have a file nanoribbon-geo_opt-pos-1.xyz
.
From this we have to extract only the molecule into a new xyz
file since we are going to create a STM image only for the molecule, not the substrate:
echo "224" > nanoribbon.xyz # the molecule contains 224 atoms, first line in a XYZ file echo "optimized using DFTB/FF" >> nanoribbon.xyz # second line in a XYZ file, an optional comment tail -n 1064 nanoribbon-geo_opt-pos-1.xyz | head -n 224 >> nanoribbon.xyz # the complete system contains 1064, but we want only the first 224 (= the molecule)
Calculating the nanoribbon
To get the actual electron density, we are now going to run a full DFT calculation using a large basis set (TZV2P
) on the nanoribbon geometry. This time, the input nanoribbon.inp
is similar to what you are used to, except the fact that to to speed the calculation up, we use the wavefunctions from the previous calculation as a starting point (the RESTART_FILE_NAME
option).
This calculation will take a while to finish: run it in parallel using 4 processes (mpirun -np 4 …
) and in the background.
At the end of it, you can visualize the orbitals using either VMD as shown in previous exercises, or using one of the scripts provided in the newly installed Python package, which will generate a series of images (nanoribbon-WFN_*.png
), visualizing the different orbitals:
cube-plot.py --cubes *WFN*.cube --plot --position 22 # the position for the contour plot should be 2 angstrom above the surface
Generating the STM image
To get an actual STM image, we now have to combine the wavefunctions into a single one:
# use your output file of the full DFT calculation as your levelsfile! cp2k-sumbias.py --cubes *WFN*.cube --levelsfile nanoribbon.out --vmin -2.0 --vmax 2.0 --vstep 0.5 | tee sumbias.out # and pipe the output to the file sumbias.out and the screen simultaneously
The parameters –vmin
, –vmax
and –vstep
determine which bias voltages for the tip (the potential between the substrate/molecule and the tip) you want to simulate (in our case $-2.0$, $-1.5$, … $2.0$).
It is important to note that for a given bias voltage, for example $-2.0$ (current from the substrate/molecule towards the tip) all orbitals with an energy between $-2.0 eV$ and $0 eV$ have to be taken into account.
At this point you should have a new set of combined CUBE files: stm_-2.00V.cube
..stm_+0.00V.cube
..stm_+2.00V.cube
, one for each bias voltage, containing the respective electron density.
From these we can finally generate the actual STM images:
# zcut is the minim z-height stm.py --stmcubes stm_*.cube --isovalues 1.0e-7 --zcut 22 --plot
Which should give you a set of files stm_*V.cube.iso1e-07.png
.
Why are there no images for certain bias voltages? Would you expect the same for a metallic substrate?