Unphysical energy fluctuations vs. forces during BFGS relaxation (depends on pseudopotentials)

[Nejc Flajšar]

Dear all,

I am studying copper surfaces with various adsorbed species, primarily chlorine atoms. During my calculations, the BFGS relaxation occasionally gets stuck. Deleting the BFGS history in Quantum ESPRESSO does not resolve the issue. From the output files, it appears that the optimizer attempts to move atoms along the force direction, but any displacement results in an increase in total energy.

To investigate this, I performed additional calculations on a structure where the relaxation was stuck. It seems that atoms can become trapped in a local energy minimum that is not reflected in the forces. I have included example QE and environ input/output files for one such case (stuck_relax*). For one such stuck structure, I varied only the x-coordinate of the top Cl atom and computed the energy and forces at each point. The total energy and the force in the x-direction on the Cl atom are plotted below.

Plots are calculated with norm-conserving pseudopotentials, i have included the input/ouput files for calculation at x=0.0 (environ.in, scf.in, scf.out)

Clear fluctuations in the energy can be observed. They are very small, but enough to confuse BFGS. Interestingly, the severity of this issue depends on the choice of pseudopotentials (PPs). I tested several combinations, including ultrasoft (US) and PAW, but norm-conserving pseudopotentials appear to give the most stable results. The deriv_low_pass filter can also help, although I have not found a consistent set of parameters (p1 and p2) that works in all cases.

I further examined the debug files for the local potentials and found that the fluctuations are seen in vsoftcavity.cube. When projecting these potentials onto a 2D plane, both the fluctuations and the differences between pseudopotentials become clearly visible (see attached pictures). Setting environ_type = 'input' and enabling only the dielectric contribution of water—while disabling the surface tension term—seems to reduce the fluctuations. The pictures show vsoftcavity.cube, first is calculated for ecutrho = 1000, second for 300. Used pseudopotentials are PAW type.

Since the differences in energy and the positions of the minima are relatively small, I am currently using relaxation algorithms that rely only on forces (e.g., damped dynamics), which appear to work reasonably well. However, I would like to understand whether I might be missing something, whether this issue is known, or if there are better ways to avoid it.

If I have not provided sufficient context or explanation, please let me know—I would be happy to share additional input/output examples.

Thank you for your time and any suggestions.

Best regards,
Nejc