BRMIX: very serious problems the old and the new charge density differ (for calculations with more than 272 atoms)

[reynaldo.putra]

Edit:
Hi, I tried uploading the input and output files, but I keep on getting the message "Sorry, the board attachment quota has been reached."...

Original post:
Hello,

I ran a calculation with the following parameters (the detailed input and output files are attached in the ZIP file in this post):
Type of calculation : ionic relaxation (PAW_PBE), surface calculation
No. of atoms : 352 atoms
VASP version : 6.3.1
K-Points : Monkhorst-Pack, 2x2x1
POTCAR files used : La, Ti_sv, O, N

I used an optimized bulk structure from the .cif file of a published study, and created a surface structure out of it.

Problem:
There appears the following warning starting from the 13th SCF iteration up until the 75th iteration in the slurm.out output file:

Code: Select all

DAV:  75    -0.108180352312E+07    0.16340E+07   -0.55867E+05 20016   0.411E+03BRMIX: very serious problems
 the old and the new charge density differ
 old charge density:  2066.91331 new 2050.34505
    0.240E+03

End result: job terminated with the following warning:

Code: Select all

|     Error EDDDAV: Call to ZHEGV failed. Returncode = 15 1 16                |
|                                                                             |
|       ---->  I REFUSE TO CONTINUE WITH THIS SICK JOB ... BYE!!! <---- 

Other notes:
Previous surface calculations containing 80-92 atoms were successful.
The goal of this test was to investigate the effect of increasing the surface area by 4x (axbxc to 2ax2bxc)

Possible questions:
1. Is there a limit to the number of atoms in a cell in VASP? Since this problem only exists when i increase the number of atoms to 272 and 352. I have also tried running a similar calculation with 192 atoms, and it has completed several ionic relaxation steps without any problems so far.
2. I briefly read that the Returncode with 15 and 1 in the first and second numbers indicates a problem in the matrix diagonalization routine and the initial phase of the diagonalization process respectively. I am currently trying ALGO = VeryFast, as opposed to ALGO = Normal which causes this error. Which one should I try, ALGO = VeryFast or ALGO = Fast? I have read that VeryFast applies the RMM-DIIS algorithm while Fast applies the RMM-DIIS combined with Block-Davidson. What does this practically mean?

Thank you for your attention.
Best regards,
Reynaldo

[manuel_engel1]

Hi Reynaldo,

Thanks for reaching out. We are aware of the upload problems and are trying to find a permanent solution. In the meantime, I have updated the upload quota. Could you please try again to attach the input and output files? Otherwise, it will be very difficult to understand what is happening.

Regarding your questions:

1. Is there a limit to the number of atoms in a cell in VASP? Since this problem only exists when i increase the number of atoms to 272 and 352. I have also tried running a similar calculation with 192 atoms, and it has completed several ionic relaxation steps without any problems so far.

In principle, there is no limit to the number of atoms in a VASP calculation. However in practice, there are limits due to numerical noise and computational complexity. A larger cell with more atoms is usually more difficult to converge to the same accuracy than a smaller cell with fewer atoms. It may be necessary to adjust certain INCAR tags in order to obtain reasonable results. But before I make any further assumptions, I would like to see your input/output files.

2. I briefly read that the Returncode with 15 and 1 in the first and second numbers indicates a problem in the matrix diagonalization routine and the initial phase of the diagonalization process respectively. I am currently trying ALGO = VeryFast, as opposed to ALGO = Normal which causes this error. Which one should I try, ALGO = VeryFast or ALGO = Fast? I have read that VeryFast applies the RMM-DIIS algorithm while Fast applies the RMM-DIIS combined with Block-Davidson. What does this practically mean?

This is indeed an error in the diagonalization routine for the Kohn-Sham Hamiltonian. You may try different values for ALGO, but I doubt that this will solve the underlying issue in this calculation. The corresponding VASP wiki page provides the list of different values and links that explain, for example, the advantages and disadvantages of RMM-DIIS.

[reynaldo.putra]

Hello Manuel,

Thank you for your assistance, I am indeed able to upload the input/output files now. Below are the attached files.

Thank you as well for addressing my questions. In relation to them, I have tried some calculations with ALGO=VeryFast. I have further questions regarding these methods.

1. Is it true that in general ALGO=VeryFast requires more NELM than ALGO=Normal? I tested with a run which I have successfully converged using ALGO=Normal, and at almost every ionic step the electronic self-consistency was not achieved.

2. Using ALGO=VeryFast, some warnings as shown below would appear, but the calculations are still running as if nothing is wrong. What do they signify?

Code: Select all

WARNING in EDDRMM: call to ZHEGV failed, returncode =   8  4     4
WARNING in EDDRMM: call to ZHEGV failed, returncode =   8  4    10
WARNING in EDDRMM: call to ZHEGV failed, returncode =   8  4    11
WARNING in EDDRMM: call to ZHEGV failed, returncode =   8  4    17

[manuel_engel1]

1. Is it true that in general ALGO=VeryFast requires more NELM than ALGO=Normal? I tested with a run which I have successfully converged using ALGO=Normal, and at almost every ionic step the electronic self-consistency was not achieved.

This is hard to say in general, since RMM-DIIS and Davidson are too entirely different algorithms. However, what I can say is that RMM-DIIS is usually less robust than (blocked) Davidson, albeit much faster. This means that RMM-DIIS can more easily run into issues like the one you are experiencing. Using pure RMM-DIIS is a bit risky as the final result depends on the quality of the original orbitals. Therefore, we recommend to combine ALGO=VeryFast with a large NELMDL (the wiki, for example, recommends NELMDL=12).

2. Using ALGO=VeryFast, some warnings as shown below would appear, but the calculations are still running as if nothing is wrong. What do they signify?

I'm not 100% sure. It is likely related to what I mentioned in my response to your first question.

I checked your INCAR file and there a few things I can comment on

• You run a spin-polarized calculation (ISPIN=2) but do not specify the initial magnetic moments (MAGMOM). This could explain some of the convergence problems.
• The INCAR tag LPETIM is unknown to me. Is this maybe a typo?
• For ionic relaxations with volume changes (ISIF=3), we recommend an even larger energy cutoff (in this case, at least ENCUT=520) and PREC=Accurate.
• Sometimes, real-space projections (LREAL=Auto) can be inaccurate and it pays off to switch back to LREAL=False.

I hope that some of these comments help resolve your problem.

[reynaldo.putra]

Dear Manuel,

Thank you very much for your inputs! I will give an update once I have newer results.

Best,
Reynaldo

[reynaldo.putra]

Hello Manuel,

I have some follow-up questions in relation to the ISPIN and MAGMOM tags. The system I am studying seems to not be reported to have any significant magnetic configurations. How can I know how to specify the MAGMOM tag?
For a large/many-atom system should I initially run ISPIN=1, and then run ISPIN=2 from the optimized structure?

As for the ISIF tag, I keep ISIF=3 even for an 88-atom slab. I initally used ISIF=3 for a bulk structure, and has been keeping it the same ever since. Since I am going to a many-atom system, should I keep ISIF=2 or 4?
What things should I consder when deciding whether to allow the cell to change shapes or not?

[manuel_engel1]

The system I am studying seems to not be reported to have any significant magnetic configurations. How can I know how to specify the MAGMOM tag?

If the system is non-magnetic, then there is no need to run with ISPIN=2. If the magnetic configuration is unknown, I would try to switch off symmetry (ISYM=-1) and let the code run with the default MAGMOM. You might be able to learn about the magnetic configuration from the final results, but it might also be stuck in a local minimum.

For a large/many-atom system should I initially run ISPIN=1, and then run ISPIN=2 from the optimized structure?

Yes, this can be a good strategy to speed up convergence.

As for the ISIF tag, I keep ISIF=3 even for an 88-atom slab. I initally used ISIF=3 for a bulk structure, and has been keeping it the same ever since. Since I am going to a many-atom system, should I keep ISIF=2 or 4?
What things should I consder when deciding whether to allow the cell to change shapes or not?

The safest strategy is always to do things separately. You can focus first on volume/cell-shape changes and afterwards on the change of atomic positions. So after setting ISIF=3, it's always a good idea to rerun the calculation with ISIF=1 or ISIF=2. It might even be good to restart the calculation after the initial volume relaxation without any changes to account for Pulay Stress.