strange fluctuation in PBE0 convergence

[jasius]

One caveat. When I change from converged PBE0 ISMEAR=0 to ISMEAR=-5 (ALGO=NONE) I get that number of bands changed. Are the results still to be trusted? There is a warning on that:

found WAVECAR, reading the header
number of bands has changed, file: 944 present: 1168
trying to continue reading WAVECAR, but it might fail
WARNING: stress and forces are not correct
POSCAR, INCAR and KPOINTS ok, starting setup
FFT: planning ...
reading WAVECAR
random initialization beyond band 945
the WAVECAR file was read successfully
charge-density read from file: Chabazite_01Ns0
magnetization density read from file 1
entering main loop
N E dE d eps ncg rms rms(c)
NONE 1 -0.771861615614E+04 -0.77186E+04 0.00000E+00 0 0.000E+00
NONE 2 -0.771861615614E+04 0.00000E+00 0.00000E+00 0 0.000E+00
1 F= -.77186162E+04 E0= -.77186162E+04 d E =0.000000E+00 mag= 5.0000
writing wavefunctions

[alexey.tal]

Are you using the same number of MPI ranks in both calculations? If you change the number of rank the number of bands can also change (see NBANDS). Since the bands beyond 944 are random numbers, these should not be trusted and the DOS in this energy range is not going to be correct, but the bands below 944 should be fine.

[jasius]

Is it core number (quite possible) or changing from ISYM=1 for initial PBE to ISYM=3 for PBE0? I need to experiment a little but I do not believe I changed NCORE.

[alexey.tal]

Even if you don't change NCORE but change the number of MPI ranks, the number of bands may change so that NBANDS is divisible by the number of ranks. I would recommend using the same ISYM in both PBE and PBE0. Otherwise, the symmetry and the number of k-points may change and you won't be able to restart your calculation from WAVECAR.

[jasius]

Thank you! That led me to my last question - I have seen so many tutorials and suggestions online on calculating DOS using hybrid. Would band structure calculated using different ISYM, ALGO, directly calculated with ISMEAR=-5 and ALGO=Normal or first with ISMEAR=0, ALGO=D and then ISMEAR=-5 and ALGO=None would all be the same?

[jasius]

Sorry, still trying to understand. Ran simple PBE, switched to

LHFCALC = .TRUE.
ALGO = Damped
TIME = 0.8
ICHARG=11

still getting this nbands. All on the same HPC, no other changes in NCORE etc

found WAVECAR, reading the header
number of bands has changed, file: 944 present: 1184
trying to continue reading WAVECAR, but it might fail

[alexey.tal]

Thank you! That led me to my last question - I have seen so many tutorials and suggestions online on calculating DOS using hybrid. Would band structure calculated using different ISYM, ALGO, directly calculated with ISMEAR=-5 and ALGO=Normal or first with ISMEAR=0, ALGO=D and then ISMEAR=-5 and ALGO=None would all be the same?

One can use ALGO=None to recalculate quantities such as DOS, dielectric function etc. using the previously converged WAVECAR. We have a wiki page which explains in detail which smearing method to use. Since the xc in PBE0 depends on orbitals (not density) it is recommended to use ALGO=Damped, which doesn't use the density mixer and is usually more robust for the hybrid functional calculations. But in principle both algorithms Normal and Damped should converge to the same ground state.

Sorry, still trying to understand. Ran simple PBE, switched to

LHFCALC = .TRUE.
ALGO = Damped
TIME = 0.8
ICHARG=11

still getting this nbands. All on the same HPC, no other changes in NCORE etc

Could you please provide OUTCARs from these PBE and PBE0 calculations?

[jasius]

Please see attached. I am still trying to converge this case wit PBE0 using either ALGO=D (and then ISMEAR=0 followed by ISMEAR=-5 route or ALGO=NORMAL and ISMEAR=-5.
I am somewhat confused that after I run ISYM=2 PBE case and do ICGARGE=11 and ISYM=3 I get, in all cases, this below:

found WAVECAR, reading the header
number of bands has changed, file: 944 present: 1184
trying to continue reading WAVECAR, but it might fail
WARNING: stress and forces are not correct

[jasius]

When it rains it pours... What do I do if PBE0 ALGO=D does not converge? I attached files as well as PBE preconverged output

entering main loop
N E dE d eps ncg rms ort
gam= 0.000 g(H,U,f)= 0.693E+01 0.169E+01 0.340E-15 ort(H,U,f) = 0.000E+00 0.000E+00 0.000E+00
SDA: 1 -0.317521870433E+04 -0.31752E+04 -0.68954E+01 4736 0.862E+01 0.000E+00
gam= 0.225 g(H,U,f)= 0.166E+01 0.713E+00 0.708E-10 ort(H,U,f) = 0.529E+00 0.898E+00 0.102E-09
DMP: 2 -0.317917994753E+04 -0.39612E+01 -0.21522E+01 4736 0.237E+01 0.143E+01
gam= 0.225 g(H,U,f)= 0.545E+00 0.256E+00 0.203E-30 ort(H,U,f) =-0.264E+00 0.488E+00 0.460E-30
DMP: 3 -0.318032211585E+04 -0.11422E+01 -0.68177E+00 4736 0.802E+00 0.224E+00
gam= 0.225 g(H,U,f)= 0.395E+00 0.940E-01 0.724E-22 ort(H,U,f) =-0.261E+00 0.209E+00 0.209E-21
DMP: 4 -0.318063986564E+04 -0.31775E+00 -0.38195E+00 4736 0.489E+00-0.525E-01
gam= 0.225 g(H,U,f)= 0.289E+00 0.346E-01 0.467E-20 ort(H,U,f) =-0.227E+00 0.851E-01 0.161E-19
DMP: 5 -0.318077519985E+04 -0.13533E+00 -0.23357E+00 4736 0.324E+00-0.141E+00
gam= 0.225 g(H,U,f)= 0.169E+00 0.195E-01 0.224E-13 ort(H,U,f) =-0.164E+00 0.333E-01 0.907E-13
DMP: 6 -0.318083113023E+04 -0.55930E-01 -0.12743E+00 4736 0.189E+00-0.130E+00
gam= 0.225 g(H,U,f)= 0.832E-01 0.755E-02 0.519E-19 ort(H,U,f) =-0.909E-01 0.168E-01 0.137E-18
DMP: 7 -0.318086538918E+04 -0.34259E-01 -0.59270E-01 4736 0.908E-01-0.740E-01
gam= 0.225 g(H,U,f)= 0.301E-01 0.251E-02 0.548E-21 ort(H,U,f) =-0.360E-01 0.660E-02 0.181E-20
DMP: 8 -0.318088331547E+04 -0.17926E-01 -0.20799E-01 4736 0.326E-01-0.294E-01
gam= 0.225 g(H,U,f)= 0.100E-01 0.193E-01 0.110E-23 ort(H,U,f) =-0.123E-01 0.220E-02 0.768E-23
DMP: 9 -0.318086584314E+04 0.17472E-01 -0.21657E-01 4736 0.293E-01-0.101E-01
gam= 0.225 g(H,U,f)= 0.260E-02 0.902E-02 0.111E-17 ort(H,U,f) =-0.350E-02 0.135E-01 0.109E-16
DMP: 10 -0.318088066297E+04 -0.14820E-01 -0.11088E-01 4736 0.116E-01 0.997E-02
gam= 0.225 g(H,U,f)= 0.489E-03 0.308E-02 0.858E-26 ort(H,U,f) =-0.699E-03 0.703E-02 0.173E-24
DMP: 11 -0.318088874109E+04 -0.80781E-02 -0.39928E-02 4736 0.357E-02 0.633E-02
gam= 0.225 g(H,U,f)= 0.916E-04 0.136E-01 0.107E-28 ort(H,U,f) =-0.963E-04 0.238E-02-0.420E-27
DMP: 12 -0.318087538510E+04 0.13356E-01 -0.11335E-01 4736 0.137E-01 0.228E-02
gam= 0.225 g(H,U,f)= 0.341E-04 0.340E-01 0.106E-20 ort(H,U,f) =-0.244E-04 0.945E-02-0.867E-20
DMP: 13 -0.318084902388E+04 0.26361E-01 -0.28937E-01 4736 0.340E-01 0.943E-02
gam= 0.225 g(H,U,f)= 0.121E-04 0.153E-01 0.147E-28 ort(H,U,f) =-0.784E-05 0.241E-01 0.730E-28
DMP: 14 -0.318087314376E+04 -0.24120E-01 -0.16620E-01 4736 0.154E-01 0.241E-01
gam= 0.225 g(H,U,f)= 0.796E-05 0.255E-01 0.543E-22 ort(H,U,f) =-0.973E-05 0.120E-01 0.434E-22
DMP: 15 -0.318086004575E+04 0.13098E-01 -0.22541E-01 4736 0.255E-01 0.120E-01
gam= 0.225 g(H,U,f)= 0.130E-04 0.392E-01 0.957E-24 ort(H,U,f) =-0.543E-05 0.181E-01 0.527E-24
DMP: 16 -0.318084255022E+04 0.17496E-01 -0.34610E-01 4736 0.392E-01 0.181E-01
gam= 0.225 g(H,U,f)= 0.648E-05 0.666E-01 0.469E-29 ort(H,U,f) =-0.654E-05 0.280E-01 0.594E-29
DMP: 17 -0.318080730288E+04 0.35247E-01 -0.58322E-01 4736 0.666E-01 0.280E-01
gam= 0.225 g(H,U,f)= 0.111E-04 0.565E-01 0.112E-25 ort(H,U,f) =-0.101E-04 0.476E-01 0.185E-25
DMP: 18 -0.318082014161E+04 -0.12839E-01 -0.53792E-01 4736 0.565E-01 0.476E-01
gam= 0.225 g(H,U,f)= 0.643E-05 0.436E-01 0.285E-23 ort(H,U,f) =-0.114E-04 0.418E-01 0.651E-23
DMP: 19 -0.318083676660E+04 -0.16625E-01 -0.42420E-01 4736 0.436E-01 0.418E-01
gam= 0.225 g(H,U,f)= 0.755E-05 0.171E-01 0.202E-21 ort(H,U,f) =-0.429E-05 0.323E-01 0.515E-21
DMP: 20 -0.318087092746E+04 -0.34161E-01 -0.19489E-01 4736 0.171E-01 0.323E-01
gam= 0.225 g(H,U,f)= 0.649E-05 0.139E+00 0.180E-21 ort(H,U,f) =-0.978E-05 0.134E-01 0.421E-21
DMP: 21 -0.318071254439E+04 0.15838E+00 -0.11340E+00 4736 0.139E+00 0.134E-01
gam= 0.225 g(H,U,f)= 0.102E-04 0.658E-01 0.311E-23 ort(H,U,f) =-0.102E-04 0.979E-01 0.406E-23
DMP: 22 -0.318080867954E+04 -0.96135E-01 -0.70280E-01 4736 0.658E-01 0.979E-01

[alexey.tal]

I see now. So the reason NBANDS changes between PBE and PBE0 is that you specify MAGMOM in PBE. In the case of your PBE0 calculation, the MAGMOM values are not specified and the default value is 1.0. Hence, VASP needs to increase the number of bands to accommodate all these unpaired electrons.
Also, if you set ICHARG=11 in your PBE0 step, you select a non self-consistent calculations, i.e., the density will be kept constant. I am not sure that this is what you want.

[alexey.tal]

What do I do if PBE0 ALGO=D does not converge? I attached files as well as PBE preconverged output

The first thing to try when ALGO=D does not converge is to use a smaller TIME. I see that you have increased the time step by factor 2 wrt the default value of 0.4, which is likely too large.

[jasius]

Thank you for your help. I am following various pieces of information on VASP tutorials, advice provided here on the forum etc so perhaps it is hard for me to come up with something that works.
I am setting MAGMOM in PBE and not in PBE0 since I am hoping to read the guess, which is PBE. So I do not reset MAGMOM to some preset value but rather obtain it from the converged WAVECAR. Are you saying I should keep MAGMOM as is?

The second point is that ICHARG=11 is everywhere in tutorials and even in the beginning of this thread... Should I, from my now vast experience, reiterate the procedure, perhaps you can confirm?

First step: PBE, ALGO=NORMAL, ISYM=2, MAGMOM = 16*0 15*0 160*0 16*0 256*0 1*6, ISTART=1, ICHARG=1, LWAVE=TRUE, LCHARG=TRUE, ISMEAR=0, SIGMA = 0.05
Second step: LHFCALC = .TRUE., ALGO = Damped, TIME = 0.4, ISYM=3, MAGMOM = 16*0 15*0 160*0 16*0 256*0 1*6, ISTART=1, ICHARG=1(?), ISMEAR=0, SIGMA = 0.05
Third step: LHFCALC = .TRUE., ALGO = NONE, TIME = 0.4, ISYM=3, MAGMOM = 16*0 15*0 160*0 16*0 256*0 1*6, ISTART=1, ICHARG=1(?), ISMEAR=-5

[jasius]

What do I do if PBE0 ALGO=D does not converge? I attached files as well as PBE preconverged output

The first thing to try when ALGO=D does not converge is to use a smaller TIME. I see that you have increased the time step by factor 2 wrt the default value of 0.4, which is likely too large.

I am trying to follow what is available online, in this case Si tutorial online https://www.vasp.at/tutorials/latest/hybrids/part1/

[alexey.tal]

Thank you for your help. I am following various pieces of information on VASP tutorials, advice provided here on the forum etc so perhaps it is hard for me to come up with something that works.
I am setting MAGMOM in PBE and not in PBE0 since I am hoping to read the guess, which is PBE. So I do not reset MAGMOM to some preset value but rather obtain it from the converged WAVECAR. Are you saying I should keep MAGMOM as is?

I understand it can get quite confusing. I was just trying to explain why the number of bands changed when you restarted your calculation without setting MAGMOM. However, it shouldn't matter if extra conduction bands are added - it shouldn't affect the result, i.e., if you set MAGMOM in your PBE0 step it won't affect the calculation if the WAVECAR is found and read.

The second point is that ICHARG=11 is everywhere in tutorials and even in the beginning of this thread... Should I, from my now vast experience, reiterate the procedure, perhaps you can confirm?

ICHARG=11 is used in the tutorials for the step where the band structure is calculated. In this step, the calculation is restarted from WAVECAR and the eigenvalues are calculated at k-points along the trajectory, which does not require self-consistency because we use the converge orbitals, hence the ICHARG=11. However, it only works if you restart from a converged WAVECAR. In your case, you did a PBE0 calculation with ICHARG=11 starting from PBE WAVECAR.

I am trying to follow what is available online, in this case Si tutorial online https://www.vasp.at/tutorials/latest/hybrids/part1/

This tutorial is for non-magnetic silicon, where such a large time step works. In your case, the system is much more complicated and such a large step doesn't seem to work. Generally, I would recommend using the default values for tags if you are not sure what value to choose, as the VASP default values are quite a conservative choice and should work for most cases. Overall, the defaults in VASP are very helpful, so if you are not sure what a tag is doing and whether you need it - don't set it at all. For example, are you sure you need LPETIM? It doesn't exist in VASP. Or you set EDIFFG which controls the ionic steps convergence, but you don't do an ionic relaxation.

I don't know what exactly your ultimate goal is so it is difficult to advise on the workflow. Here is what I would suggest for a DOS calculation:
First step: PBE, ALGO=NORMAL, ISPIN=2, MAGMOM = [your values], ISMEAR=0, SIGMA = 0.05, LREAL=AUTO
Second step PBE0 restarting from PBE WAVECAR: ALGO=Damped, ISPIN=2, LHFCALC = .TRUE., TIME = 0.4 (or even 0.2), ISMEAR=0, SIGMA = 0.05, LREAL=Auto

If then you want to make a DOS plot with the tetrahedron method, i.e., no smearing:
Third step PBE0 restarting from PBE0 WAVECAR: ALGO=None, ISPIN=2, LHFCALC = .TRUE., ISMEAR=-5, LREAL=Auto

[jasius]

Thank you! Forgive me one more time, I do not understand your comment below. What you described is what I am doing, e.g. converge with PBE and then read WAVECAR for PBE0 but with ICHARG=11. So converged WAVECAR from PBE read for PBE0.

ICHARG=11 is used in the tutorials for the step where the band structure is calculated. In this step, the calculation is restarted from WAVECAR and the eigenvalues are calculated at k-points along the trajectory, which does not require self-consistency because we use the converge orbitals, hence the ICHARG=11. However, it only works if you restart from a converged WAVECAR. In your case, you did a PBE0 calculation with ICHARG=11 starting from PBE WAVECAR.