Convergence issue with NEB

[anirudhnk]

Dear Prof. Henkelman,

I have been trying to run a NEB calculation for a structure but it has never converged despite running around 200 ionic steps. Recently, it ran for 12 days (6 days first and then re-ran it) and still reported an error of "JOB ON cg11-1 CANCELLED AT 2022-09-02T18:44:58 DUE TO TIME LIMIT." The structure im trying to study is Na23Ge136. The symmetric structure is Na24Ge136. I have removed one sodium atom to see the migration of a sodium vacancy from one site to another. Please find attached the files.

Thank you.

[graeme]

Check if you have attached the files - they would help me see what's going on.

[anirudhnk]

Sorry, i forgot to attach them. Have attached them now.

[graeme]

So I will almost always recommend using cheaper setting to find minima and saddles. You can always easily reconverge them with more precise settings. If you look at other questions on this forum, you can see many examples where people struggle getting convergence just because they focus all their effort on precision of the electronic structure and not on geometric convergence.

In what I'm posting, I got rid of all your extra electrons and just used the standard Na and Ge potentials. Then I could use a cutoff of 180 eV. Starting with the gamma point, I was able to converge the endpoints and the band in an hour or so (min1, min2, and neb). Note that a major problem was the poor geometric convergence of your endpoints. There was a significant geometric relaxation of the Na atoms that you missed. Again, I blame it on putting too much weight on ENCUT, core electrons and KPOINTS, so that you just get tired of waiting for a relaxation. My reduction in simulation parameters probably saves two orders of magnitude in computational time and I was able to see the important relaxation in a hour on my old cluster with 24 core nodes.

Then, again, it is very easy to turn up the settings and check for convergence. Going from gamma point to a 2x2x2 kpoint mesh did change the barrier from 0.7 to 1.0 (see min1_2k, min2_2k, and neb_2k), but the geometries are very similar. You can follow along these lines and check for convergence in other degrees of freedom.

On a positive note, the minimum energy paths and barriers look smooth.

[anirudhnk]

Thank you very much for your response. So do you recommend always starting with cheaper settings to find the minima?
Because i've always been using my current settings and its worked. Only this time i faced this error..
Also, i see that you used NSW = 500 and had to run it thrice for it to converge, can i instead set it to 1500 or 2000 so that it finishes in one run itself?

[graeme]

Yes, you will save computer time starting with cheaper settings. You will also get a sense of the uncertainty and importance of parameters when you make them more accurate.

The IBRION=3 POTIM=0.1 setting is nice and safe but also pretty inefficient. You can switch to IBRION=1 to save time, especially if you are away from high-force regions.

You can use whatever NSW makes sense. I generally prefer to have something reasonable to the jobs finish and I can assess how they are doing.

[anirudhnk]

Dear Prof. Henkelman,

I did try running the neb with the above settings and was surprised to see that it converged within an hour itself! However i get an energy barrier like this - please see thr attached figure. Is it ever possible that the initial and final states have higher energy than the transition state itself?

[graeme]

Check that the settings on the endpoints are the same as the band. You will (almost certainly) find that they are different. Note: I put OUTCAR files in 00 and 08, which need to be consistent with the settings of your NEB calculation.

[anirudhnk]

I see, but i have another question. Wont the cell shape and volume to change when doing geometry relaxation? By using ISIF = 2 aren't you constraining the cell shape? I generally allow the volume to change by using ISIF=3 and then in the 2nd step, i use ISIF =2...

[graeme]

That is more of a scientific question rather than a technical one. The nice thing about a constant volume is that you benefit from cancelation of errors in terms of the plane wave basis and kpoint sampling. However, when changes in the cell volume are important, then they should be considered.

[anirudhnk]

Got it, thank you very much

[anirudhnk]

Dear Professor,

I have a general question to ask. Why do some of the atoms, particularly on the edges, disappear after a relaxation or an neb calculation? Like in the CONTCAR file, the no. of atoms still remains the same, but when i visualize it in VESTA, some of the atoms have disappeared. This was after using the looser settings for relaxation (not with my earlier stringent conditions). Please see the attached images for comparison.

Thank you.

[graeme]

I expect that atoms that are precisely on the boundary (such as the yellow atom on the corners) are shown on all of the corners. If that atom relaxes, it will only be shown in one location (top-right) instead of on all 8 corners.

[anirudhnk]

Is it the case with atoms within the cell as well? As two yellow atoms that are inside have also disappeared after relaxation..

[graeme]

Yes, those are on the bottom face and being duplicated on the top face before relaxation.