Dear all,
I am exploring the capabilities of LAMMPS for water simulations. For typical water models, the CPU requirements should scale linearly with the number of distances that need to be calculated for each pair of molecules, at least approximately. In SPC like models, there are 9 such distances (all pairs O-O, O-H and H-H), while for TIP4P models there are 10 distances (the previous 9 pairs + the distance between LJ interaction sites on each molecule). However, when I run short MD simulations using both models, TIP4P takes almost twice as long as SPC. Anybody knows why this is so? Apparently, Gromacs gets the right scaling.
The final lines of the log files are:
For SPC:
Dear all,
I am exploring the capabilities of LAMMPS for water simulations. For typical
water models, the CPU requirements should scale linearly with the number of
distances that need to be calculated for each pair of molecules, at least
approximately. In SPC like models, there are 9 such distances (all pairs
O-O, O-H and H-H), while for TIP4P models there are 10 distances (the
previous 9 pairs + the distance between LJ interaction sites on each
molecule). However, when I run short MD simulations using both models, TIP4P
takes almost twice as long as SPC. Anybody knows why this is so? Apparently,
Gromacs gets the right scaling.
the main difference is that gromacs supports dummy atoms on a general
level while LAMMPS doesn't.
there is not only the number of distances that matter but also the
effort to project the M-point forces back on O and H and do it
correctly so for flexible geometries and different M offsets. i don't
know the details of the gromacs implementation of tip4p, but i would
not be surprised if there are some corners that are being cut for
better performance at the expense of some generality and accuracy.
axel.
also, i cannot reproduce your timing differences when doing a real
apples to apples comparison. there is about 50% overhead in my tests.
please see the attached log files. this is the same density and the
same cutoff settings. TIP4P has to have a few more neighbors, because
it needs to have the coordinates of the hydrogen atoms within reach to
compute the position of M where with SPC/E they don't matter. this is
managed automatically by the tip4p pair style and you can see the
larger number of neighbors in the log.
log.tip4p-shake.gz (1.68 KB)
log.spce.gz (1.7 KB)
Dear Axel,
Thanks a lot for your detailed answer.
I was a little puzzled by the difference between your results and mine. It
took me some time to find out where the difference was. If I run in the NVT
ensemble I get:
SPC
Files.tgz (3.92 KB)
Dear Axel,
[...]
There is a huge difference (100% in pair time). I am attaching the full log
files.
Do you have any idea why this is happening? Are pair forces being calculated
yes, i have.
a second time for the TIP4P model in order to get the pressure?
no. what happens is that for SPC/E LAMMPS can use a much more
efficient way to compute the stress tensor since it has the
coordinates for all particles and then the virial stress can be
computed from f_i * r_i with a single loop over all atoms. however,
for TIP4P, this is not possible since the location of the M point is
computed on the fly and then the forces and virial stress
contributions have to be tallied for each pair of atoms.
i've implemented some optimizations to the process in the styles in
USER-OMP and it would be interesting to see how this works out. you
don't need to compile with OpenMP, just install the package and use
-suffix omp when you run. those optimizations are independent of using
OpenMP.
axel.