I am trying to implement a water model with drude oscillators that utilizes the Thole damping function, however, I am having trouble getting the system to stabilize. The errors I am getting are:
ERROR on proc 0: Bond atoms 971 975 missing on proc 0 at step 165
which happens when the bond energy becomes very high (due to the drude oscillator not converging?) and the temperature failing to thermostat to the correct temperature.
To try to understand what is happening, I am also outputing pair_write to plot the potential and force. I then also run the same input file without thole damping for comparison. I find that the potential barely changes at all, even at short distances, in the pair_write potential and force output for the input file with thole damping and without thole damping.
One possibility of the error is that I don't have the correct parameters (a damping parameter, alpha polarizability) for the thole damping pair coefficients. However, regardless of what I change these coefficients, the pair_write gives the same potential, and I think the force also does not change. I am hoping to be able to see the potential and force change so that I can tune the thole coefficients to proper values.
I have attached the input file and data file. Please let me know if there are any suggestions to make the model run on lammps.