When writing a new fix for LAMMPS in the source code, one sets a mask for the fix which defines at what point in the simulation it’s called, and then writes a method for the fix with the same name, e.g. setting the mask as
int mask = 0;
mask |= POST_FORCE;
then defining a method post_force which is called at the post_force stage.
I would like to know if there is any way to have a fix called after the forces from bonds, pairs, dihedrals and so on are set - but BEFORE contributions from any other fixes, such as fix NVT.
In other words, I would like to have LAMMPS call my method in my fix when the contents of the f[][] array has only the forces due to the gradient of the potential and nothing else such as thermostatting, langevin dynamics, and so on.
If not, is there any other way to get my hands on this information in each timestep?
When writing a new fix for LAMMPS in the source code, one sets a mask for
the fix which defines at what point in the simulation it's called, and then
writes a method for the fix with the same name, e.g. setting the mask as
int mask = 0;
mask |= POST_FORCE;
then defining a method post_force which is called at the post_force stage.
I would like to know if there is any way to have a fix called after the
forces from bonds, pairs, dihedrals and so on are set - but BEFORE
contributions from any other fixes, such as fix NVT.
In other words, I would like to have LAMMPS call my method in my fix when
the contents of the f[][] array has only the forces due to the gradient of
the potential and nothing else such as thermostatting, langevin dynamics,
and so on.
fixes are executed in the order they are defined. so all you'd have to
do is to make sure that your fix is defined first.
What Axel says is correct, within the categories of points at
which fixes can be called. The post-force catetgory is the point
in time at which only bond/pair/etc forces are present, so if your
fix is the first post-force fix defined, it will have the untouched forces.
You can also use fix store/force to store those untouched forces,
which your fix could then access, even if your fix operated later in
the timestep.
Thanks for your replies. I need to access and manipulate the forces that the velocity Verlet algorithm will use to advance, and to do this before the thermostat alters them. I have made my fix be set in the input script before the thermostat. Would post_force be the correct category for my fix to work in?
I have noticed that if I ask my fix to print out the forces from f[][] on the very first time step, they are different depending on if I set NVT or NVE. They are the same if I use initial_integrate instead, but this gives the wrong answer to a test case so it’s obviously not doing what I intend.
Thanks for your replies. I need to access and manipulate the forces that the
velocity Verlet algorithm will use to advance, and to do this before the
thermostat alters them. I have made my fix be set in the input script before
the thermostat. Would post_force be the correct category for my fix to work
in?
yes. have a look at Verlet::run()
the modify class calls all post_force fixes right after all forces are
computed and the reverse communication for newton's 3rd law across
domains is done.
I have noticed that if I ask my fix to print out the forces from f[][] on
the very first time step, they are different depending on if I set NVT or
NVE. They are the same if I use initial_integrate instead, but this gives
the wrong answer to a test case so it's obviously not doing what I intend.
please recall how the velocity verlet algorithm works and how it is
implemented in LAMMPS:
a) you start from initial positions and velocities and forces.
b) then update the velocities for a half step from the forces
c) then update the position for a full step from the velocities
d) then recompute the forces
e) then update the velocities for another half step based on the new forces
f) continue from step b) until done
step a) is done during the "setup" phase, where coordinates and
velocities are read in (or set) and then the initial forces computed.
step b) and c) is handled in the "initial integrate" phase of a fix
step d) after that (and after a potential neighborlist rebuild)
step e) at the "final integrate" phase
and "post force" is right between d) and e)
so the reason that you see a different between nve and nvt is not that
you didn't get to access the forces properly, but rather due to the
fact that fix nvt modified the forces of the *previous* step, i.e.
those computed during the setup phase. if you insert your fix before
the integrator in "initial integrate", then - of course - you will
still see the original initial forces.