I am simulating a liquid-vapor interface and I want to use fix
recenter so that the middle of the liquid phase is always centered in
the simulation cell. However, when I used fix recenter (fix 1 all
recenter 0.5 NULL NULL units fraction) to constrain the center of mass
to the middle of the simulation cell I got unexpected results. The
middle of the liquid phase wasn't even close to the center of the
simulation cell. I could only get this to work if I modified
Group::xcm() to not unwrap atom coordinates. I don't care about the
"correct" center of mass, I just want the liquid phase to be in the
middle. Am I missing something and is there a way to do this without
modifying the code? If not, a "no unwrapping" option in fix recenter
shown below might be nice to include in LAMMPS.
Thanks,
Stan
Here is the original code:
/* ----------------------------------------------------------------------
compute the center-of-mass coords of group of atoms
masstotal = total mass
return center-of-mass coords in cm[]
must unwrap atoms to compute center-of-mass correctly
------------------------------------------------------------------------- */
void Group::xcm(int igroup, double masstotal, double *cm)
if (rmass) {
for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) {
xbox = (image[i] & 1023) - 512;
ybox = (image[i] >> 10 & 1023) - 512;
zbox = (image[i] >> 20) - 512;
massone = rmass[i];
cmone[0] += (x[i][0] + xbox*xprd) * massone;
cmone[1] += (x[i][1] + ybox*yprd) * massone;
cmone[2] += (x[i][2] + zbox*zprd) * massone;
}
}
I changed the code to this to get what I wanted:
if (rmass) {
for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) {
massone = rmass[i];
cmone[0] += x[i][0] * massone;
cmone[1] += x[i][1] * massone;
cmone[2] += x[i][2] * massone;
}
The reason fix recenter uses unwrapped coords, is that
if it didn't then when an atom crossed a periodic boundary,
you would get a delta change in the COM, and thus
in the position of all atoms which is unphysical. With unwrapped coords,
you get only an epsilon change in the group COM, even
when an atom crosses a PBC.
I am guessing the problem you are having initially,
is that you have image flags with big values (e.g. you've
run for a while already), so that the COM is not where you
think it is. If you use the INIT option,
it should choose the COM to be whatever it currently is,
and running from there with unwrapped coords should be fine.
Steve
The reason fix recenter uses unwrapped coords, is that
if it didn't then when an atom crossed a periodic boundary,
you would get a delta change in the COM, and thus
in the position of all atoms which is unphysical. With unwrapped coords,
you get only an epsilon change in the group COM, even
when an atom crosses a PBC.
I am guessing the problem you are having initially,
is that you have image flags with big values (e.g. you've
run for a while already), so that the COM is not where you
think it is. If you use the INIT option,
it should choose the COM to be whatever it currently is,
and running from there with unwrapped coords should be fine.
wouldn't it be better to use fix spring or fix momentum
on this kind of system? if i am not mistaken, fix recenter
can easily hide the "flying icecube syndrom" and thus
it should be important to monitor the total center of mass
momentum.
axel.
Steve,
I had been running the simulation for a while and when I used fix
recenter to constrain the COM to the middle of the cell I got the
attached results, which don't make much sense to me. Using the INIT
option could work, but it doesn't allow me to specify the location of
the COM, only to use the current value, right? So is the only way to
get the center of mass right in the middle of the simulation cell to
use fix recenter from the very beginning of the simulation?
The reason fix recenter uses unwrapped coords, is that
if it didn't then when an atom crossed a periodic boundary,
you would get a delta change in the COM, and thus
in the position of all atoms which is unphysical. With unwrapped coords,
you get only an epsilon change in the group COM, even
when an atom crosses a PBC.
I don't understand the problem with changing the position of all of
the atoms together--basically now you are just using a relative
coordinate system instead of an absolute one. You could think this as
taking a bunch of snapshots of a wandering interface and then shifting
the coordinates of the snapshots (using PBCs) so that the interface is
always in the same spot. I don't see how this would affect the
dynamics of the simulation--can you explain more?
Thanks,
Stan
Dear Stan,
I was wondering if the issue could be the following: the liquid/vapor
interface results from a dynamical equilibrium, with atoms constantly
evaporating and condensating. As a result, each atom in the system
could have crossed a different number of PBC and have very different
unwrapped coordinates, resulting in this surprising center of mass.
Anyway Axel's suggestion to use fix momentum should work better in this case.
Best,
Laurent
Comments below.
Steve
Steve,
I had been running the simulation for a while and when I used fix
recenter to constrain the COM to the middle of the cell I got the
attached results, which don't make much sense to me. Using the INIT
option could work, but it doesn't allow me to specify the location of
the COM, only to use the current value, right? So is the only way to
get the center of mass right in the middle of the simulation cell to
use fix recenter from the very beginning of the simulation?
All LAMMPS knows how to do is displace the atoms so
their COM is at some location. If you run for a while and
the COM drifts off to some arbitrary location, then you ask
it to move it (at the start of the next run) to some radically different
location, then you will likely break your model. If you
want to move atoms between simulations you should look
at the displace_atoms command.
I would only use fix recenter to keep the COM where it already is.
In your case, maybe that means you need to start that at
the very beginning.
The reason fix recenter uses unwrapped coords, is that
if it didn't then when an atom crossed a periodic boundary,
you would get a delta change in the COM, and thus
in the position of all atoms which is unphysical. With unwrapped coords,
you get only an epsilon change in the group COM, even
when an atom crosses a PBC.
I don't understand the problem with changing the position of all of
the atoms together--basically now you are just using a relative
coordinate system instead of an absolute one. You could think this as
taking a bunch of snapshots of a wandering interface and then shifting
the coordinates of the snapshots (using PBCs) so that the interface is
always in the same spot. I don't see how this would affect the
dynamics of the simulation--can you explain more?
If you are using fix recenter on all the atoms, this should be fine.
It's when you use it on a subset of the atoms and they displace
a large distance (in one step) relative to other atoms, that you
have problems. Some people do this, b/c they want their
substrate or solute molecule to stay in the center of the box,
and let water molecules move freely around it. Big displacements
then cause overlaps.