Hello!
Interesting question, I actually also discussed this with my colleague Jonathan because we were both not sure how this is handled internally in VASP. But before diving into details let me mention the INCAR tag LATTICE_CONSTRAINTS which may help you to perform the special form of structure relaxation you mentioned without having to modify the VASP source code.
If I understood your argument correctly it boils down to this: there may be situations where the net forces acting on ions are arbitrarily close to zero while the lattice is not yet relaxed (stress is non-zero). A simple example would be the following: imagine a perfectly relaxed simple cubic lattice, i.e., a cubic box with a single atom inside (or even a supercell thereof). Here, the energy is at a minimum, all forces are zero and there is no stress, everything is relaxed. If I now scale one of the lattice vectors, then the situation changes in the way we seek to understand: the box is squeezed (or stretched) in one direction and hence the energy is out of the minimum and we expect a non-zero stress tensor. However, due to symmetry the total forces acting on the ions are still all zero!
Now, in the VASP Wiki we write that setting EDIFFG < 0 will check only for the norms of all force vectors. Hence, as you correctly pointed out, this criterion alone cannot suffice to properly relax the situation described above. We did a little research in the code and found that the Wiki is actually incomplete! VASP does handle the above situation by introducing another criterion that is in line with your expectation that there should be a "stress stopping criterion". Depending on your VASP version somewhere around line 4400 in main.F one can find two lines which deal with the stopping criterion:
Code: Select all
! This one is the force criterion:
IF (SQRT(TIFOR(1,NI)**2+TIFOR(2,NI)**2+TIFOR(3,NI)**2) &
& >ABS(DYN%EDIFFG)) LSTOP2=.FALSE.
! ... and a few lines below there is one for the stress components:
IF (ABS(D2SIF(J,I))/FACTSI/T_INFO%NIONS>ABS(DYN%EDIFFG)) LSTOP2=.FALSE.
While the first one captures the force magnitudes as expected, the second criterion is actually undocumented but active whenever ISIF > 2. It reuses the value of EDIFFG multiplied by a conversion factor (from force to stress units, incorporating the step size POTIM) to define the desired stress criterion. Hence, in the "EDIFFG < 0 mode" VASP will check the forces and the stress to determine whether the ionic loop should be stopped (if ISIF > 2).
Therefore, it is safe to assume that VASP will not incorrectly stop at an unrelaxed situation similar to the one described above. However, it is true that there is no separate INCAR tag which allows one to control the stress threshold. I can only speculate that the conversion factor has been set up in a robust way since it has not been changed in roughly 20 years now! Anyway, I will ask my colleagues if there is any further insight on this topic.
All the best,
Andreas Singraber
