Different state ordering with FERWE/FERDO y LDIAG = .FALSE between VASP 5.4.4 and 6.6.1

[wilver_muriel]

When constrained occupations are used, i.e., ISMEAR = -2 and LDIAG = .FALSE., together with FERWE and FERDO, the ordering of the states appears to be preserved only in the spin channel in which the excitation is performed.

To investigate this behavior, I performed the following test: I calculated the excitation of a carbon dimer in hBN using a (4\times4) supercell and 100 bands. The system has two defect-related levels corresponding to bands 64 and 65. Band 64 is predominantly composed of the \(p_z\) orbital of carbon atom 32, whereas band 65 is predominantly composed of the \(p_z\) orbital of carbon atom 31. Band 64 is initially occupied, while band 65 is unoccupied. The excitation was performed in the spin-up channel by occupying band 65 and leaving band 64 unoccupied.

I then performed a structural relaxation using VASP 5.4.4 and VASP 6.6.1. The results were as follows:

VASP 5.4.4: the ordering of the states is preserved in both spin channels. This behavior was verified by examining the orbital contributions to the bands in the PROCAR file (PROCAR_5.4.4):

spin up

band 64
...
29 0.000 0.000 0.009 0.000 0.000 0.000 0.000 0.000 0.000 0.009
30 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
31 0.000 0.000 0.037 0.000 0.000 0.000 0.000 0.000 0.000 0.037
32 0.000 0.000 0.192 0.000 0.000 0.000 0.000 0.000 0.000 0.192
tot 0.000 0.000 0.513 0.000 0.000 0.000 0.000 0.000 0.000 0.513

band 65
...
29 0.000 0.000 0.025 0.000 0.000 0.000 0.000 0.000 0.000 0.025
30 0.000 0.000 0.005 0.000 0.000 0.000 0.000 0.000 0.000 0.005
31 0.000 0.000 0.204 0.000 0.000 0.000 0.000 0.000 0.000 0.204
32 0.000 0.000 0.024 0.000 0.000 0.000 0.000 0.000 0.000 0.024
tot 0.000 0.000 0.512 0.000 0.000 0.000 0.000 0.000 0.000 0.512

spin down

band 64
...
30 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
31 0.000 0.000 0.029 0.000 0.000 0.000 0.000 0.000 0.000 0.029
32 0.000 0.000 0.226 0.000 0.000 0.000 0.000 0.000 0.000 0.226
tot 0.000 0.000 0.518 0.000 0.000 0.000 0.000 0.000 0.000 0.518

band 65
...
30 0.000 0.000 0.005 0.000 0.000 0.000 0.000 0.000 0.000 0.005
31 0.000 0.000 0.184 0.000 0.000 0.000 0.000 0.000 0.000 0.184
32 0.000 0.000 0.024 0.000 0.000 0.000 0.000 0.000 0.000 0.024
tot 0.000 0.000 0.499 0.000 0.000 0.000 0.000 0.000 0.000 0.499

VASP 6.6.1: the ordering of the states is preserved in the spin-up channel; however, in the spin-down channel, the ordering of the states changes. This behavior can be observed from the band contributions in the PROCAR file (PROCAR_6.6.1):

spin up

band 64
...
30 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
31 0.000 0.000 0.040 0.000 0.000 0.000 0.000 0.000 0.000 0.040
32 0.000 0.000 0.188 0.000 0.000 0.000 0.000 0.000 0.000 0.188
tot 0.000 0.000 0.505 0.000 0.000 0.000 0.000 0.000 0.000 0.505

band 65
...
30 0.000 0.000 0.005 0.000 0.000 0.000 0.000 0.000 0.000 0.005
31 0.000 0.000 0.195 0.000 0.000 0.000 0.000 0.000 0.000 0.195
32 0.000 0.000 0.028 0.000 0.000 0.000 0.000 0.000 0.000 0.028
tot 0.000 0.000 0.508 0.000 0.000 0.000 0.000 0.000 0.000 0.508

spin down

band 64
...
30 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
31 0.000 0.000 0.017 0.000 0.000 0.000 0.000 0.000 0.000 0.017
32 0.000 0.000 0.237 0.000 0.000 0.000 0.000 0.000 0.000 0.237
tot 0.000 0.000 0.521 0.000 0.000 0.000 0.000 0.000 0.000 0.521

band 65
...
30 0.000 0.000 0.002 0.000 0.000 0.000 0.000 0.000 0.000 0.002
31 0.000 0.000 0.002 0.000 0.000 0.000 0.000 0.000 0.000 0.002
32 0.000 0.000 0.006 0.000 0.000 0.000 0.000 0.000 0.000 0.006
tot 0.000 0.000 0.457 0.000 0.000 0.000 0.000 0.000 0.000 0.457

In both calculations, I used only the (\Gamma) point and the same general calculation settings. I have attached the corresponding PROCAR, POSCAR, and INCAR files for both calculations.

I would like to know whether this behavior could be caused by an issue with the input files or by a change in the implementation of these options between the different VASP versions. In particular, I would like to determine whether this is expected behavior or whether it could indicate a possible bug in VASP 6.6.1.

[christopher_sheldon1]

Dear Wilver Muriel,

Thank you for your question. Could you upload the OUTCAR files and I will look into this? Could you also tell me which POTCAR files you used?

In principle, it shouldn't matter, as both states that you are looking at are unoccupied, so they have no impact on the energy, density, or forces. Setting LDIAG = .FALSE. is likely the cause of this. In the LDIAG wiki page, it says:

Generally using LDIAG = .FALSE. is only advised, if one wants to maintain a certain orbital order, for instance, when reading the orbitals from an existing WAVECAR file.

Without LDIAG = .TRUE., the order will not be consistent between runs.

I'd recommend trying the calculation with LDIAG = .TRUE. and see if this fixes the issue.

Best wishes,

Chris

[wilver_muriel]

Dear Chris,

Thank you for your response, and I apologize for the delayed reply.

I have attached the OUTCAR files from calculations performed with VASP 6.6.1 and VASP 5.4.4 (OUTCARs.zip).

I used the PAW pseudopotentials from potpaw_PBE.64.

Regarding LDIAG = .TRUE., unfortunately, this is not a suitable option for my calculation because I am specifically interested in calculating excited states of defects using constrained occupations. Therefore, I need to use LDIAG = .FALSE. together with FERWE and FERDO.

Thank you again for your help and for looking into this issue.

Best regards,

Wilver Muriel

[christopher_sheldon1]

Dear Wilver,

Thank you for the files. I'm afraid only the 6.6.1 OUTCAR was uploaded; the other is the PROCAR from 5.4.4.

I understand your problem a little better. In the 6.6.1 version, the PROCAR shows that spin-down band 65 is not localised on C31. Instead, it is delocalised across B5, B7, B12, and B14. It has switched with band 68, which has the expected C31 localisation. Since they are both unoccupied, it will not affect the excited state.

The two bands crossed after the first ionic step. There is a big jump in the energy there:

free energy TOTEN = -276.22762964 eV
free energy TOTEN = -275.43861716 eV
free energy TOTEN = -276.49043521 eV

Once they have crossed, they stay that way for the rest of the structure relaxation. The ionic step may be too large, so reducing POTIM = 0.2/0.3, or using IBRION=1 might help. I suspect that you have two different final structures, which might explain the differences in PROCAR. OUTCAR_5.4.4 should clarify this.

Best,

Chris

[wilver_muriel]

Dear Chris,

I apologize for the mistake when uploading the files.

I ran a calculation using VASP 6.6.1 with POTIM = 0.2, and obtained the expected behavior. In any case, I have attached the OUTCAR_5.4.4, PROCAR_6.6.1_POTIM_0.2, OUTCAR_6.6.1_POTIM_0.2 files for your reference.

Thank you for your suggestion.

Best regards

Wilver Muriel

[christopher_sheldon1]

Dear Wilver,

The issue was due to a single line in main.F in VASP 5.4.4, which gave inconsistent energies for fixed electron occupancies. This is related to a known issue (https://vasp.at/wiki/index.php/Known_is ... ownIssue76), where the total minimiser stops at too high an energy for ALGO = Damped and ALGO = All. This was fixed in version 6.3.0, so I would avoid using VASP 5.4.4 for these calculations.

Which allowed 5.4.4 to correctly match the energies:

Code: Select all

                    Unpatched 5.4.4          5.4.4 + fix              6.6.1
Single point (eV)   -276.0997                -276.22762984            -276.22762951
Relaxation (eV)     -276.394 (23 steps)      -276.53848734 (15 steps) -276.53848704 (15 steps)

I'd recommend using 6.6.1, as this has the bug fix, as well as a few others (https://vasp.at/wiki/index.php/Known_is ... ownIssue58 - fixed in 6.4.0, and https://vasp.at/wiki/index.php/Known_is ... ownIssue60 - fixed in 6.6.0), and so will give the correct energies and orderings. In VASP 6.6.0, the default for LDIAG = .TRUE. was updated so that Loewdin diagonalisation (IFLAG = 4) is used by default - which preserves order - rather than full diagonalisation (IFLAG = 3) - which can reorder; this was not an issue in VASP 5.4.4 as IFLAG was not used in the default routine. Let me know if there is anything else that needs clarifying; otherwise I'll close this thread next week.

Best wishes,

Chris

PS: there is a section on a wiki page about this: https://vasp.at/wiki/Delta_self-consist ... occupation