Dear LAMMPS Developers,
I am interested in carrying out a NVE simulation of SPC/E water with a custom end_of_step()-Fix. That fix makes a tiny modification to velocities and coordinates and I ran into a problem when trying to extend FixShake in such a way that it can take the changes into account. Let’s assume that the fix has the following form:
void FixSomefix::end_of_step() {
// modify x and v a tiny bit
domain->pbc();
}
In FixShake::post_force() the forces are modified such that the coordinates will satisfy the constraints after the next integration step exactly (w.r.t. the specified tolerance). However, if FixSomefix::end_of_step() modifies velocities and/or coordinates after that, the forces need to be adjusted again to take the modifications into account. If user-defined fixes modify the forces before FixShake::post_force() is invoked, everything should be fine. However, in my case I want to do that at the end of the timestep to avoid other complications.
I was thinking of the following procedure:
Steve, thank you for your comments and your help! Actually, I got it to work right now after I removed domain->pbc() and another bug (not in the methods of the previous posts, they were OK).
With regard to the motivation of doing this, I completely agree with you and I realized that even without calling shake twice, the geometry of the molecules was still preserved up to a tolerance of 1.e-5 (max. error). I am more than happy with that, but I needed to know that I can control the error, if it was necessary. It was more like a proof of concept. Now everything is fine!
Please allow me to ask you two more questions:
-
When you said “better do that check in fix_shake”, were you referring to testing that the constraints are satisfied or to checking whether there are any other end_of_step()-fixes after fix_shake?
-
How can I make sure that fix_shake::end_of_step() is the last end_of_step()-fix without messing around in modify, i.e. fix[list_end_of_step[i]]->end_of_step()? Should I just change the order in the input file or is there another way?
Best wishes,
Peter