[lammps-users] Old vs New versions, Tersoff


I would have a question for the developers of the code about Tersoff potential. I am using an older version of LAMMPS, and running NVE run on 1 million Silicon atoms. I compared the parallel performances of the this version against the new on, on a quad-core machine. And it seems the older version runs much faster. The only difference I found are 2 coefficients for the 3 body problem, that are not present in my version. Would this only explain the difference in performances? And would the old version considered accurate?
Here are the tables (in seconds, for 100 steps):

LAMMPS, old version

N 1 cpu 2 cpu 3 cpu 4 cpu
tN 688 344 276 179
t1/tN 1 2 2.49 3.84

LAMMPS, latest version (26 Jul 08)

N 1 cpu 2 cpu 3 cpu 4 cpu
tN 981 628 368 260
t1/tN 1 1.56 2.67 3.77

Thanks in advance

I'll refer this Q to Aidan. Which older version are you using?
Did you get the same answers. Did things like the number
of neighbors or the number of reneighborings change between
the 2? Please post a simple-as-possible version of your
input script.



the version I am using is dated feb07, but patched up to nov07 I think. I don’t see no huge differences between the results. The pressure calculation at mid time is different, but seems to reach closer values at 100 steps. The number of neighbor is approximately the same (16402488 vs 16395990). Here is my input file (the “;” is for a newline):

variable x index 50 ; variable y index 50 ; variable z index 50
units metal ; atom_style atomic ; lattice diamond 5.43 ; region box block 0 $x 0 $y 0 $z
create_box 1 box ; create_atoms 1 box

pair_style tersoff
pair_coeff * * …(path)/potentials/Si.tersoff Si

mass * 28.0
velocity all create 3000.0 376847 loop geom
neighbor 1.0 bin ; neigh_modify every 1 delay 5 check yes
fix 1 all nve ; timestep 0.005
thermo 50

run 100

Please let me know if you want to see anything else.Thank you

2008/8/7 Steve Plimpton <[email protected]>


The addition of the two extra parameters added some extra floating-point
operations, and also required fetching some more values from memory, both of
which will cause some slow down to occur. However, there may be other
reasons for the differences in timing too. For example, did you use
identical compile and run platforms for both codes?

The agreement between the forces at step 0 should be very high. Can you
verify this? On later timesteps, the trajectories will diverge because of
finite machine precision.



thanks for your answer. I made sure that I used the same makefile and machine. But as you said, the agreement at step 0 is almost perfect. After 100 steps, the divergence in total energy is almost null, and about 1.5 % for the temperature, which is acceptable for a system of a million atoms I guess.

I would have one more question: I haven’t been through the Tersoff potential theory, but I am wondering about the origin of those parameters. How was the 3-body problem treated in LAMMPS before you add them? Does it follow some approximation already done in the theory, or was there another reason for them not to be there?

Thank you

2008/8/7 Thompson, Aidan <[email protected]>


    I have found using Berendsen thermostat will cause the system to have a
nonzero velocity. Is this due to the nature of this thermostat?
    I used nve together with Berendsen. Below is an example of my input file.
The nonzero velocity of center of mass appears in x and y direction, while in
the z direction I used a LJ93 potential.



units real
boundary s s s
atom_style molecular
read_data ./sim_10_473_.data
include mass.bs
include pair.bs
include ff.bs
special_bonds 0.0 0.0 1
neigh_modify delay 1
velocity all create 473 64289816 mom yes rot yes
fix 1 all temp/berendsen 473 473 200
fix 2 bcb rigid molecule
fix 3 ch2 nve
fix 4 all wall/lj93 zlo 0.0 0.44 3.7 7.65
timestep 4

thermo 10000
dump 1 all atom 10000 sim_10_473_.equili.dump
dump_modify 1 scale no
run 2500000

unfix 1
fix 5 all temp/berendsen 473 300 200
run 500000

unfix 5
fix 6 all temp/berendsen 300 300 400
undump 1
thermo 1000
restart 3200000 sim_10_473_.restart
dump 2 all atom 1000 sim_10_473_.dump
dump_modify 2 scale no
run 200000

Not sure what you are asking. Maybe it is that the system
has zero center-of-mass velocity in x,y initially but acquires it
later due to the Berendsen thermostat? I don't think this should
happen, although it might due to roundoff. It could be that
you are also using fix rigid. I would try thermostatting the
non-rigid atoms (e.g. the solvent) if possible. I've had more
success with using fix langevin as a thermostat on rigid bodies.


This turned out to be an interesting issue with the math libraries
different compilers use. I saw no difference in old vs current versions
of Tersoff with g++. But I saw something similar to your results with
the Intel compiler, which is also generates much faster Tersoff code
than g++.

The reason turned out to be that the pow() function was used as
pow(arg,m) in the
current version and as pow(arg,3.0) in the old version. M is now a variable
which can be 3 or 1. I think the Intel pow() must check for a constant
when compiling and avoid taking a log() or something like that when it
can shortcut
the answer.

Since all Tersoff variants we know of have m = 1 or 3, I just posted a patch
that uses m as a constant instead of a variable in the pow() and it
recovers nearly all of the old speed.


In the Berendsen thermstat, you scale all velocities by the same
factor at any given time step. Therefore if you start with zero
center of mass velocity, you should not be acquiring it later in
the simulation. I doubt if round-off errors can contribute to creating
any significant linear momentum in a randomized system - Steve,
has this been observed/reported to happen? However, I am not
sure what effect fix-rigid will have on the system.


Fix langevin can give random kicks to the overall
system causing the center-of-mass to undergo
a random walk. But its a small effect due to
being nearly zero in a big system, as you indicate.


Hi Steve,

   This is what I was asking. Sorry for the confusion.
   I 'm in fact running a simulation of polystyrene on substrate, where no
solvent is present. In this case, can I still use the fix langevin? From what
I read, it involves implicit solvent. So I'm concerned that it may not be
   Thank you!