In short, I am trying to implement my own constraint to NPT dynamics in ASE (fixing center of mass, as well as no-overall-rotation. The latter is not implemented in the existing constrain classes, I think). My naive definition of my own constraint class didn’t work, so I’m trying to understand the how constraint is applied in NPT source code. However, I couldn’t find the places where the constraint is implement in the NPT source code. Please tell me which line the constrains are implemented.
In details, I can identify how other molecular dynamic handle constrains.
For instance, I found VelocityVerlet enforce constrain in the line “atoms.set_positions(r + self.dt * p / masses)”. Langevin has “atoms.set_positions(x + self.dt * self.v + self.rnd_pos)” where the set_positions from the atoms class will enforce the constrains.
The closest thing I found in NPT that apply constraint is
"
def _setbox_and_positions(self, h, q):
“”“Set the computational box and the positions.”“”
self.atoms.set_cell(h)
r = np.dot(q + 0.5, h)
self.atoms.set_positions(r)
"
but the r = np.dot(q + 0.5, h) doesn’t make sense to me. It doesn’t look like the expect q=q_old+ t * v.
I’m not very familiar with Nose-Hoover thermostats. Is the above code where the constrain is implemented?
Thanks,
Bohan