Charge differences restarting qeq run

Hello LAMMPS users,

I am running simulations using qeq/point which requires stopping and starting simulations. At first, when running these simulations, I was noticing extended periods of time in the restarted simulations where the density was decreasing ~20% before recovering to the expected value. This was confirmed to be caused by the restart by running two sequential short simulations, and one longer simulation with the same initial configuration. When I looked at the charges for these simulations, I noticed that the charges for the larger simulations matched those of the first simulation perfectly, but once the second short simulation started, the charges were drastically different, at times over one e different. Relevant files can be found at the link below Am I doing something I shouldn’t be doing or is this being caused by a bug in LAMMPS?

Thanks for your help,
Michael Humbert


thanks for providing a complete and illustrative input de​ck for testing.
there are multiple things at play here:

- fix qeq uses initial guesses extrapolated from the past 5 QEQ
calculations. these are not restarted, but reinitialized at the beginning
of a run. for systems with a delicate balance between charges, this can
lead to different minima. this effect applies to restarting with a restart
file and using 2 run statements with half the steps in the same input.

- you are using fix npt, which will result in volume changes and thus can
trigger changes in the ewald sum parameter estimator. this can enhance
differences. when using a restart file, fix npt will restore its internal
state, when restarting from a data file, it cannot. if at all possible,
consider switching to fix nvt. the fact, that you use a rather coarse
kspace convergence, enhances the impact of the volume change. for cells
with changing volume, the kspace convergence should be selected tighter, so
that there is only a limited impact on the accuracy from volume changes.

Hi Axel,

Thank you for your suggestions. I have tried running simulations with the change in the code, as well as a kspace convergence of 1e-8, but I am still running into the issue with the density (see attached image). I strongly believe that this large density deviation is a result of only performing qeq on a subset of the simulation. I have tried running a simulation with all atoms in the qeq group (without too much parameterization on the other atoms so take the results with a grain of salt), and the fluctuations in the density for the second short run are similar to those of the longer run.

Michael Humbert



there isn’t currently much that i can do about this situation, that to acknowledge, that this doesn’t look good and needs some careful testing and review.

there are two more suggestions:

  1. can you check whether the same or a similar issue appears when using the other qeq variants?

  2. can you please report the situation with the multiple input decks included at so that this doesn’t get lost and forgotten. there is no time frame for resolving issues posted there, but at least it will remain visible and somebody running into the same issue can see it as well as perhaps somebody may feel compelled to correct the problem and, for example, make fix qeq properly restartable.



Maybe Ray has comments on this.


In response to checking if other methods work, it seems like point will converge to an incorrect value, fire does not converge at all no matter the tolerance and number of cycles but it appears that the dynamic method will converge to the correct value

Michael Humbert