Hello, all

Recently, I did a simple simulation, which was a tension
of a bar by fixing one end and exerting force on boundary
atoms of the other end (called hi_fixed hereafter).

There were three stages in my simulation:

first, relax the sytem to equilibrium,
then, exert tension gradually from 0 to the maximum,
finally, keep the maximum force for some time steps.

The corresponding statements I used in the last two stages are as follow:

dimension 3
boundary p p s
variable runstep_relax equal 40000
variable myscalor equal 1e0

# prepare boundary
# -----------------------
velocity lo-fixed set 0.0 0.0 0.0
fix 1 lo-fixed setforce 0.0 0.0 0.0
velocity hi-fixed set 0.0 0.0 0.0

# ------------
fix 5 all nve
unfix 4
unfix 5

# relax
# ------------
fix 4 hi-fixed aveforce 0.0 0.0 \{force\_every\_atom\} &nbsp;&nbsp;fix 5 all nve &nbsp;&nbsp;run {runstep_relax}
unfix 4
unfix 5

The center-of-mass of the boundary atoms group hi_fixed
is also computed :

compute hi_com hi-fixed com

However when I did post-processing, I found the center-of-mass
around the center-of-mass of hi_fixed when equilibrium is obtained
at the first relax stage. But when it entered the stage of keeping
the maximum load (the third stage), the oscillation center of the
center-of-mass of hi_fixed showed a jump and was discontinuous with
the second stage distinctly.

I also outputed the total force on the group of hi_fixed before
fix aveforce took action. It showed a mean value of 0 at the gradual

I think it is very strange. Could anyone point out the problem?

Wang Kai

I don't know. But you can monitor the COM of a group of
atoms during a simulation run (see the fix com command),
so you can diagnose when it does something unexpected.
And you can compare the COM values to your post-processing
to see if they agree.

Steve

Sorry, that should have been the compute com command.

Steve

Hi Kai,

In the loading stage the "fix 5 all nve" is applied to all so it
overwrites what you have done to hi-fixed, and this is why you see com
oscillate about its equilibrium position. The reason is the same for
the description of total forces, since group hi-fixed is just part of
NVE ensemble.

You should apply nve to group all-hi_fixed-lo_fixed, and probably
before "fix aveforce". Also make sure the force_every_atom value is
consistent with the system unit you are using. Otherwise the constant

Best,
Ray

Hi Plimpton,

I fix my problem accidentally.
I read the manual which is corresponding to
the lammps version (10Mar) I am using, and
found the equal-style variable hadn't mentioned
in the doc of fix aveforce. So I tried a
newer version 14Aug11, and it worked like
expect.

PS: is it possible to add a warning in the log file
for those command which do support equal-style variable
when users use equal-style variable in it?

Regards,
Wang Kai

Hi Shan,

I fix my problem accidentally.
I read the manual which is corresponding to
the lammps version (10Mar) I am using, and
found the equal-style variable hadn't mentioned
in the doc of fix aveforce. So I tried a
newer version 14Aug11, and it worked like
expect.

PS: I also tried using the group mobile
I found the center of mass of hi_fixed hadn't shown
variance from the equilibrium position, so I still
think using all in fix nve is right, for fix aveforce
don't perform time integration I guess.

Thank you all the same!

Regards,
Wang Kai

PS: is it possible to add a warning in the log file
for those command which do support equal-style variable
when users use equal-style variable in it?

I presume you mean, is it possible to generate an error
when the user specifies a variable in a command
like v_foo and LAMMPS is expecting a numeric value.

The answer is yes, it's possible, but it would require
adding lots of checks in lots of places in the code in
lots of commands. So it's not likely to happen soon.

A better thing to remember is that you should never
use the doc pages for the current version of LAMMPS
(e.g. by looking at them on the WWW site) when
using an old version of the code. Instead you
should always use the doc pages included in
that version, in its doc directory. Then the
commands will do exactly what the doc pages say.

Steve