[lammps-users] PPPM problem

Hi all,

I tried to run a MD simulation modeling a slab of liquid which contains a few charged atoms (the charges equal -1 or 1) on top of a solid charged surface (1% of the solid atoms are charged). I use the lj/cut/coul/long pair_style together with the PPPM option. I have periodic boundary condition in the x and y direction but not in the z direction (reflecting wall). I use use the kspace_modify slab 3.0 option
When the default value of the dielectric constant (1.0) is chosen, the simulation runs well with the PPPM initialisation giving:
G vector = 0.0965039
grid = 15 15 60
RMS precision = 5.61492*10-5
buck FFT buffer size/proc = 24800 13500 11160
So far, so good.

As I would like to model water at ambient temperature, I then fix the dielectric constant to 80.0 and launch the simulation.
I get the followig error message: PPPM stencil extends too far, reduce PPPM order.
Therefore, I add, in the input script, the line:
kspace_modify order 1
but it does not change the error message.
I then change the previous line by
kspace_modify mesh 15 15 60 order 1 gewald 0.0965 to mimic the values given by the PPPM initialisation when the dielectric constant equals 1.0
The simulation starts but crashed after somes steps. The error message is "Out of range atoms - cannot compute PPPM"

Do you have an idea of a way to make this simulation runs ?

Thank you ahead

PS: I join the input script.

DataSalt.inp (3.07 KB)


This problem should be easy to fix. First, don't change the PPPM
order, especially not to 1. This is a bad idea. Leave it at 5. I'll
admit that the first error message you saw is misleading. It would be
better if it encouraged you to use a denser mesh. What's happening is
that when you dramatically weaken the coulombic interactions by
setting dielectric to 80, the algorithm that chooses the mesh density
thinks that very few grid points are needed in order to meet the
accuracy that you've specified. If there are fewer grid points in any
dimension than the PPPM order, you'll get that error message because
LAMMPS simply can't do the PPPM in that case. So you'll need to either
force LAMMPS to use more grid points (good idea), reduce the PPPM
order (bad idea), decrease the cutoff (might work), or increase the
desired PPPM tolerance (probably a good idea in this case). You could
also simply switch to using regular Ewald (reasonably good idea,
especially if you have less than 1000 atoms in your box). I'd
recommend simply increasing the number of grid points. "kspace_modify
mesh 15 15 60" should be fine. But don't mess with the order or the
gewald parameter (another bad idea unless you really know what you're
doing). Also, it looks like your timestep is rather large. To get
better numerical stability, you might try reducing it. If all else
fails, try sending us your data file along with the input script, and
we'll see if we can get it running here.


Hi Paul,

thank for your help.
I tried the tips you gave me. Increasing the PPPM tolerance to 0.000001 and reducing the time step to 1 fs seems to work (I have just launched the calculation). I do not fix the order and gewald anymore. I am looking forward to see if it is still running tomorrow morning.


PS: I have been through mpich, mpiexec, torque, OpenPBS to be able to run the simulation on a cluster. Do you think it is worth saying some words about that for the community ?

Paul Crozier wrote:

1.0e-6 is a high tolerance for PPPM which means it will be
costly. Once you have your system setup correctly you should
be able to run at a lower tolerance and still get good answers.


Hi Steve,

I tried to relaunch the simulation from a restart file and to specified a value of 1.0e-4 for the PPPM tolerance but the PPPM solver was not able to initialize. I had to keep 1.0e-6. Hopefully, the simulation is running on 10 processors so that the time required to reach equilibrium is satisfactory.

All the best


Steve Plimpton wrote: