ondrejkrejci / PPSTM

OUTDATED: reposity moved to: github.com/Probe-Particle/PPSTM ; BEFORE: Code simulating variousSTM techniques, especially for tilting tip (dependent on ProkopHapala/ProbeParticleModel)
7 stars 2 forks source link

no sustained output #1

Closed Msunshinep closed 4 years ago

Msunshinep commented 5 years ago

Hello, I tried to simulate the STM image with the output of CP2K, I tested the example in "4N-coronene" which can perform and give the simulated image, but for my system, the code can't perform well , without new output after:

make: `IO' is up to date. reading FHI-AIMS LCAO coefficients for basis: sp

============ define atoms ============

forced skipped line :

Number of atoms: 216 atomic geometry read Reading CP2K MOs from:crazy_mol-cartesian-mos-10.MOLog Found 3288 MOs spanned by 3600 basis functions centered on 216 atoms. The Fermi Level: 0.8933 eV; in FHI-AIMS is the Fermi automatically 0. All coefficients read Grid has dimensios: (51, 151, 171, 3) Entering the dI/dV ( sp(d)-sp(d) ) procedure We're going to C++ calculating s orb. on a tip sp

also do not break after running 30 mininutes, could you tell me the reason for this problem? Is it common or just need more time?

ondrejkrejci commented 5 years ago

Hi,

I see that you are doing the things with your own system. I don't know what kind of system do you have, but it seems that it is rather kind big (3288 Molecular Orbital is a lot), and/or at the same time you are using a lot of states (MOs), maybe many of this is not useful. Therefore it is possible, that your system is simply too big. Generally speaking (all) STM simulations are much more complex and more expensive than PP-AFM simulations and they also require bigger inside, than in PP-AFM (it is not a black box, if you don't want to calculate really long time). For some big systems 1h in C++ is not a big deal. Did you tried to adjust the dIdV_test_4N-coronene_cp2k.py or did you used (adjusted) PPSTM_simple.py script in the folder? What I can see is that your system, has 216 atoms, but I don't think that all of them are in top-most layer (generally anything 3eV bellow the top-most layer is not important in most of cases). If you have your atoms in CP2K run adjusted upper layer first, then you can select only the top-most atoms via cut_at parameters. Then second thing you can try to cut out states, that are too far in energy and doesn't contribute to the tunnelling current. This can be done via adjusting cut_max & cut_min accordingly, to what you need. One of the possibilities is to use the PPSTM_simple.py, which cut out most of the states, but maybe if you have some metal states there, then it can be hard. Please have a look at this PPSTM_simple.py script, which is described here: https://github.com/ondrejkrejci/PPSTM/wiki#ppstm_simplepy Please, don't hesitate to ask me any questions about it. Other option is to use parallelization to speed up your calculation. The code has OMP parallelization and I can help you with that if you will want. I'm sorry, that my answer is so complicated and long, but I would say that this is a nature of STM. Cheers & good luck! Ondrej

Msunshinep commented 5 years ago

Thank you very much for the detailed reply! I tried something you pointed out above, it is pleased to tell you that it can give the STM simulation now. However, the simulated STM image looks very strange , I think it maybe arise from the wrong voltage that I used for the system. I just use the default dIdV_test_4N-coronene_cp2k.py with:

pbc=(0,0)
lvs = None
#lvs = np.array([[15., 0.,0.],[0.,15.,0],[0.,0.,15.]]
WorkFunction = 5.0 #more or less standart.
fermi=0.8933    # Normally is the Fermi Level read from the *.log file, but since CP2K gives the Fermi to the energy of HOMO and for this simulation we need it to be in the middle of the gap, we are changing it
orbs= 'sp'  # 'sp' works now, 'spd' works for fireball as well
cut_min=-0.894  # HOMO 0.893 below the Fermi Level, other orbitals cut
cut_max=+0.894  # LUMO 0.893 eV above the Fermi Level
cut_at=-1   # All atoms of the molecule
eta = 0.01  # very low, to pronounce the single orbitals only
WF_decay=1.0    # for STM only - how fast the exponential decay fall, with the applied bias ( if 1 - 1:1 correspondence with bias; if 0, it doesn't change)
nV = 9  # for STM only - number of STM integrational steps nV ~ V/eta
lower_atoms=[]  # No atoms has lowered hopping - be aware python numbering occurs here [0] - means lowering of the 1st atom
lower_coefs=[]  # Lowering of the hoppings

According to your reply,I know I should to change these parameter according to my system, but, I am sorry that it is my first time to use cp2k, could you tell me, how can I determine the suitable :”fermi“,“cut_min” and ”cut_max“ according to cp2k.

And another thing I want to ask that, can the simulation of STM image with a tip have mix orbital, like current0 = PS.dIdV( V, WorkFunction, eta, eigEn, tip_r2, Ratin, coefs, orbs=orbs, s=1.0, px=1.0, py=0.0, pz = 0.0) with the "s" and "pz"?

Looking forward for your reply!

Best wishes! Shine

ondrejkrejci commented 5 years ago

Hi! Sorry, for my late reply. I was busy in the beginning of the week and then I forgot about it, sorry for that. Yes, I guess that your system is on support, while in our (4N-coronene) system we tuned everything, so it would go with the calculated values of the system (freestanding molecule). Therefore the voltage and fermi values has more or less nothing connected with the real world in this system. As for the calculations and your parameters - I would need more informations about the system, to help you more concretely (you can email me - ondrej(dot)krejci(at)aalto(dot)fi -, if you don't want to share it here). Otherwise I would suggest using _PPSTMsimply.py script, which should have more or less or the constants somehow well defined. If you are using metallic/conducting support in you calculations, than I don't think, that there is any need to shift fermi, therefore fermi should be 0.0 (it is automatically taken from CP2K). On semiconductor, one could play with it, depending on experiment and doping, but there are more things, that can play role. As for cut_min & cut_max - they should be cut_min = minimal used voltage - 2eta and cut_max = maximal used voltage + 2eta. On a real system with substrate I would suggest eta = 0.1 eV. Since I don't know what you are looking for, but also considering into account that states can be shifted in DFT (kind a lot), I would suggest to do 'V-scan' from -2 to 2 V to see if you would see the pattern you want to look at. But of course, this means cut_min -2.2 and cut_max 2.2 which means huge amount of states in your example, therefore it would take really long time. Maybe run it overnight (and day). At least use the cut_at parameter (if possible) to speed things up. But with this, you will get all dI/dVs and STM's together within one (long) run.

And as for your last question - sure the tip can be mixed during the calculation. But I don't suggest it these days, since you don't know what kind of mixture one should use. I think it is better to calculate s, px & py, pz, .... contributions first store it as WSxM, XSF or NPY and later play with different contributions. Anyway, each contribution is calculated independently, therefore 2 contributions in one run take 2time, so it is better to separate it. maybe try to look for interesting voltage with "s"tip and then run the others with interesting voltages only. as for "s"and "pz" from my point of view they ends up with almost the same contrast. (maybe pz is a little bit "sharper", but the bright points stays bright and dark points dark). you can try to see difference in some low voltage, but otherwise, I would stick with "s"for most of cases - main difference between "s", "pz"or "dz2"would be on bare metals probably. The main difference (above molecules or molecules on substrate) is mainly with px & py tip, which basically do (d("s"-tip image)/dx) 2 + (d("s"-tip image)/dy) *2.

Hope this will help you! Regards, Ondrej

ondrejkrejci commented 4 years ago

Hello!

I added into the PPSTM_simple.py script a new scan feature - "states" -- it is especially good for dIdV above molecules, just at the energy of a states between <V,V_max>. I think, this is something it can help you. I tested it in different script, so please do not hesitate to write me, if it is not working properly. I'll try to write some deeper description on wiki, unfortunately, I'm really bad with time now. If you do not have any questions within a week, I will close this issue. Cheers, Ondrej