Calculating Simulated STM

Message

sjensen313 · #1 Post by **sjensen313** » Mon Dec 19, 2011 7:42 pm

Hi Everyone,
I am trying to create simulated STM images for TiO2. How I understand it is that you do a normal run first. You then run a second run where you only electronically relax, not geometrically and create a WAVECAR, CHG, and CHGCAR. This is done with the following INCAR:

SYSTEM = slab - testPAWs: 1st run
Starting parameters for this run:
NWRITE = 2 write-flag
ISTART = 1 job : 0-new, 1-cont, 2-samecut
ICHARG = 0 charge: 0-wave, 1-file, 2-atom, >10-const
INIWAV = 1 electr: 0-lowe 1-rand 2-diag
LWAVE = .TRUE. wavefuc. is written out
LCHARG = .TRUE. charge density is written out
GGA = 91
ENCUT = 400.000
LREAL = Auto

Electronic Relaxation:
PREC = Normal low | medium | high
NELM = 500 number of ELM steps
EDIFF = 1E-06 energy stopping-criterion for ELM
LREAL = Auto real-space projection (.FALSE., .TRUE., On, Auto)
ALGO = Normal algorithm (8=CG for small, 48=RMM for big systems)
NSIM = 4 NSIM bands are optimized simultaneously (IALGO=48 only)

Ionic Relaxation:
NSW = 1 max number of geometry steps
IBRION = 2 ionic relax: 0-MD, 1-quasi-Newton, 2-CG, 3-Damped MD
EDIFFG = -0.020 force (eV/A) stopping-criterion for geometry steps
ISIF = 2 (1:force=y stress=y ions=n shape=n volume=y)
POTIM = 0.10 initial time step for geo-opt (increase for soft sys)

PAW tags:
LMAXMIX = 2
ADDGRID = .FALSE.

Dipole correction:
LDIPOL=.TRUE.

DOS related values:
ISMEAR = 0 (-5:tetrahedron/Blochl, -1:Fermi, 0:gaussian, 1:Methfessel/Paxton,)
SIGMA = 0.2 broadening in eV

"Parallelization flags:
NPAR = 1
LPLANE = .TRUE."

ISYM = 0

The parts in blue are what I changed from a normal calculation. Then, to create the PARCHG file, I create a new folder, copy the POSCAR, POTCAR, INCAR, KPOINTS, CHG, CHGCAR, and WAVECAR into it and then change the following in the INCAR:

ICHARG = 11 charge: 0-wave, 1-file, 2-atom, >10-const
NSW = 0 max number of geometry steps

STM:
LPARD = .TRUE. !Partial charge densities
NBMOD = -3 !Type of energy interval, here with respect to Fermi if I recall -3 correctly.
EINT = 0.0 1.0 !Energy interval from Fermi to +1 eV
LSEPB = .FALSE. !Whether you want band-decomposed partial CDs

The weird thing is when I run this, the outcar claims to finish calculating the partial charge density in a few minutes but it continues to run in the queue until it times out and when I load the PARCHG it has all 0s. Am I doing something wrong here?

Thanks!
<span class='smallblacktext'>[ Edited ]</span>

#2 Post by **admin** » Tue Dec 20, 2011 1:44 pm

Everything seems OK. The calculation of PARCHG takes only several minutes.

sjensen313 · #3 Post by **sjensen313** » Wed Dec 21, 2011 12:23 am

Thanks for the response! If my approach for computing the STM seems okay, what do you think is wrong and is there a good approach to use to try and track down the issue? The last few lines of my OUTCAR are:

k-point no. 1; bands inside the range 0

k-point no. 2; bands inside the range 0

k-point no. 3; bands inside the range 0

k-point no. 4; bands inside the range 0

Finished calculating partial charge density.
VASP will stop now.

I think that there should be numbers besides 0 here. It's so weird because a colleague got it to work using the same parameters so I am not sure what can be going on.

gopalkr · #4 Post by **gopalkr** » Wed Dec 21, 2011 1:13 am

I hate to hijack your thread but I'm trying to simulate STM images as well. But the problem I'm facing is that I get the error:
" WARNING: chargedensity file is incomplete
ERROR: charge density could not be read from file CHGCAR for ICHARG>10
FORTRAN STOP "

I've seen 3-4 threads here in the forum where people had posted similar issues but I didn't find a clear-cut solution. I've tried everything that has been suggested in previous threads but I've not been successful! Then, I went in to main.f90 and added print statements to see what is happening; because according to the code, the above error should occur if (ICHARG - 10) equals 0 which is not the case as the ICHARG used is 11. A 'READCH' routine is called in line 1926 of the code after which (ICHARG - 10) becomes zero (I don't know how). It is 1 before READCH is called and 0 in the next line. I tried to find where readch has been defined but I couldn't! And yes, the CHGCAR is in the same folder where I'm running the code..

It would be great if someone could explain what I'm doing wrong. My INCAR was the same as above (sjensen313); as a start, I'm just trying with diamond-cubic Ge system.

Thanks in advance.
(Admin: My question is relevant to the title of this thread; so I didn't create a new thread!)

Danny · #5 Post by **Danny** » Wed Dec 21, 2011 3:27 pm

There is a very small error in your setup:

EINT = 0.0 1.0 !Energy interval from Fermi to +1 eV
If you want to use NBMOD=-3 you only give one value at EINT--> EINT = 1.0
this will take all charge going from the fermi level up to the fermilevel +1eV

Also make sure that you get out of the bandgap

, VASP actually uses the WAVECAR you give it, and doesn't start a random one (different number of cpu's that the initial calculation can coase this,...)

A code to visualise/simulate STM based on this VASP data you can find in my sig

@gopalkr: are you sure that the CHGCAR you have is complete! accidents at copying, jobkill during write of CHGCAR,...
--> so check if the tail of your CHGCAR contains the Augmentation part, else rerun the calculation to generate a new CHGCAR.

Cheers
Danny

gopalkr · #6 Post by **gopalkr** » Wed Dec 21, 2011 5:53 pm

Thanks for the heads-up Danny. My CHGCAR contains the Augmentation part. So there are two ions in my system and there are 138 occupancies for each ion. I'm not sure if the file is supposed to contain more details.. The excerpts of the file are:
"Ge Bulk
5.857500000000000
0.484473 0.484473 0.000000
0.000000 0.484473 0.484473
0.484473 0.000000 0.484473
2
Direct
0.000000 0.000000 0.000000
0.250000 0.250000 0.250000

augmentation occupancies 1 138
<<< 138 numbers>>>

augmentation occupancies 2 138
<<< 138 numbers>>>

End of File

Thanks,
Gopal

Danny · #7 Post by **Danny** » Fri Dec 23, 2011 3:02 pm

Have you tried running both calculation in a single submit/job? (to make sure you use exactly the same setup(same machine/environment/version of VASP)--> ie run a single job script which first generates the SCF calculation, and then immediatly the partial charge calculation, using the same number of CPU's, same number of NBANDS, NSIM, NPAR=1 (IIRC--> NPAR needs to be a specific value for ICHARG>=10 calculations)

Cheers
Danny

#8 Post by **admin** » Wed Jan 04, 2012 5:16 pm

dear gopalkr
have you checked that the number of grid points of CHGCAR (NGX,...)
are the same in the run that has written CHGCAR and the run that reads CHGCAR?