=====Simulation of STM and AFM images for two short graphene nanoribbons with different chemical termination=====
<note warning>
In case you do not want to use the quantum-mobile VM, you will need to install the asetk and ProbeParticle packages:
<code>
git clone https://github.com/ltalirz/asetk
pip install -e asetk
</code>
and
<code>
git clone https://github.com/ProkopHapala/ProbeParticleModel.git
cd ProbeParticleModel/
git checkout dev
</code>

</note>
download from [[https://polybox.ethz.ch/index.php/s/CH5VdcI40YdELez|here]] the tar file exercise_10.tar, move the file to your exercise directory, and extract the content


<code>
tar -xvf exercise_10.tar
cd exercise_10 
</code>




We consider two possible chemical terminations for a finite size 7-AGNR.
In TASK_1 the ribbon is terminated with a C-H2 bonding while in TASK_2 the termination is C-H
The additional H atom present at the termini of the ribbon of TASK_1 will suppress the spin polarized
edge states that are evident in the ribbon of TASK_2


===TASK_1===
<code>
cd TASK_1
</code>

Have a look to the cp2k input file cp2k.inp
used to obtain quickly the optimized geometry of a ribbon adsorbed on a Au substrate.
The ribbon is modelled within DFTB (similar to tight binding) while the substrate is modelled
via Embedded Atom Model.
An empirical potential in the form of C6/R^6 plus a Pauli repulsion term
are added to couple the adsorbate/substrate systems.


Two geometry fiels are present: mol.xyz and all.xyz
The cp2k program will need both of them.

Have a look at the geometry of the system using ASE or vmd
for both all.xyz and mol.xyz:

<code>
ipython
In [1]: from ase.io import read

In [2]: from ase.visualize import view

In [3]: s=read("all.xyz")

In [4]: view(s)

In [5]: exit()
</code>

<note important>
submit the geometry optimisation run
<code>
./run
</code>


Have a look at the final geometry optained (you can extract the last frame from
the file PROJ-pos-1.xyz)
After completion of the optimization you should extract the final coordinates of the molecule
(first 80 atoms) and copy them in the STM directory (call them p.xyz) to compute the KS orbitals 
and to compute the STM images

<code>
tail -442 PROJ-pos-1.xyz | head -82 > p.xyz
mv p.xyz STM
</code>


Now go to the directorySTM
<code>
cd STM
</code>
and have a look to the input file cp2k.inp used to converge the
wavefunction of the system and to plot the cube files for the KS orbitals.
Execute the program
<code>
cd STM
./run
</code>
The program will compute the 4 highest occupied and 4 lowest unoccupied KS orbitals.
Visualize the orbitals with VMD (remember **+** and **-**)



To obtain the STM images you have to combine different KS orbitals (depending on the bias voltage applied)
into a single cube file:

<code>
./run_sumbias
</code>
you will then obtain a cube file for each desired bias voltage (see the script run_sumbias)

Now you can compute a constant current STM image running the script

<code>
./run_stm
</code>

Please note that we are simulating a molecule, we do not include the electrons of the substrate
thus we have a discrete spectrum of energies 

<note warning>
why some of the STM images look empty?
</note>
</note>

Now we can simulate for the same ribbon a nc-AFM image:
<note important>
Go the the AFM directory of TASK_1 (and have a look to the parameter file params.ini)
copy there the p.xyz file that you havein the STM directory
and execute:

<code>
./run_PP
</code>
It will take ~ 5 minutes, then you will find a dir containing the AFM simulated image.
</note>
===TASK_2===
Modify the geometry of TASK_1 removing one H atom from each C-H2 at the termini of the ribbon (remove two H atoms in total).
Create the corresponding mol.xyz and all.xyz files, optimize the geometry, compute STM and nc-AFM images
repeating all the instructions of TASK_1 for the scripts present in the dir TASK_2
<note warning>
Be careful: here we do a spin polarised simulation. **When doing the STM simulation (ONLY for the STM)**
we have to distinguish the three C atoms of one terminus of the ribbon from the 
three of the opposite terminus calling them C1 and C2. For these atoms
we will define a guess electronic configuration with spin up on one side and spin down on the opposite side.
This is achieved defining a occupation unbalance in the alpha and beta orbitals (try to identify this section of the input
and note that the calculation is performed for a spin multiplicity of 1)

The file p.xyz in teh STM directory should look similar to:
<code>
     78 
 i =       49, E =      -140.2738100175
  H         4.2914607718       10.2130614763       21.2815435017
  H         4.2778729017        7.7509218987       21.2954986738
  H         4.2782723704       12.6751364096       21.2955091278
  C         7.4704758534        6.5236639413       21.2352922076
  .
  .
  .
  C1        5.3788157746        7.7465647443       21.2687198580
  .
  .
  C1        5.3936844253       10.2129317839       21.2797918647
  .
  .
  C1        5.3792136407       12.6792819903       21.2687263656
  .
  .
  .
  C2       21.1530397078        7.7456205579       21.2687376504
  .
  C2       21.1385072480       10.2118877383       21.2797955201
  .
  C2       21.1533012965       12.6781430186       21.2687326678
  .
  .
  
</code>

</note>

<note important>
Look at the KS orbitals (especially HOMO and LUMO) for both spin UP and DOWN
</note>
<note important>
Notice the difference between the images in TASK_2 and the images in TASK_1
In TASK_2 we have KS states localised at the termini of the ribbon.
These states are suppressed by the addiitonal H atoms in TASK_1
</note>
<note important>
why some STM images are remarkably asymmetric? Is this correct?
</note>