Dear LAMMPS users,
I try to model infinite film of ions interacting with 1 electron. I use kspace with slab option for this. I calculate force acting on electron. It does not go to zero for large distances, as it should be for 1D effective potential. It has value about 0.6. But my estimation gives for force=charge/square/2=(13ions)*1/10^2/2=0.065.
What is my misunderstanding? What is really calculated by LAMMPS when I use kspace with slab?
input script is below:
The slab option of the Kspace commands is for
2d slabs, not 1d, as the manual explains.
I don't really know what you're trying to do. Maybe
Paul can comment.
Steve
I try to calculate molecular dynamics for a system which is infinite in two directions (x,y), but thin in z. This is a charged film. I also have some separate charged particles. I take periodic boundaries for x and y, but fixed for z. Let separate charged particles also be periodically placed in x,y. Then I calculate force acting by film on particles, using “slab” option, because I have periodicity only in x & y.
The question is “What is that force?”. I found that it is not force acting on a separate charge by an infinite film, as described by a well-known simple formula. It can not also be the force which acts on all separate particles, because this should be infinite. Then, what is it?
yours
Phil
I try to calculate molecular dynamics for a system which is infinite in two
directions (x,y), but thin in z. This is a charged film. I also have some
separate charged particles. I take periodic boundaries for x and y, but
fixed for z. Let separate charged particles also be periodically placed in
x,y. Then I calculate force acting by film on particles, using "slab"
option, because I have periodicity only in x & y.
The question is "What is that force?". I found that it is not force acting
on a separate charge by an infinite film, as described by a well-known
simple formula. It can not also be the force which acts on all separate
particles, because this should be infinite. Then, what is it?
phil,
you are missing one important piece of information:
how large is the simulation cell in z-direction and
how far does the charge distribution extend?
this is crucial for a correct application of the slab option.
the size of the simulation cell in z direction has to be
_at least_ twice the extent of the charge distribution, assuming
that the particles are centered in the box.
so if you want to compute the force on a charge far away
from the surface, you have to make the slab (i.e. the z-direction)
large enough to have converged results. otherwise you don't
follow the boundary condition requirements of the poisson
equation solution that is applied and should get bogus forces.
also, what is the total charge of your system? is it overall neutral?
cheers,
axel.
Dear Axel, I had 500 for z boudaries of simulation box, and checked up to 20000. All charges are within 300 from center of the box. OK, 500 is small. But all the same questions exist for 20000. My system now is not neutral. I made it so to check how pppm with slab works. I did not find in manual that non-neutral system is bad. Is it so? Is it impossible to do what i want to do with LAMMPS?
phil
Dear Axel, I had 500 for z boudaries of simulation box, and checked up to
20000. All charges are within 300 from center of the box. OK, 500 is small.
But all the same questions exist for 20000. My system now is not neutral. I
made it so to check how pppm with slab works. I did not find in manual that
non-neutral system is bad. Is it so? Is it impossible to do what i want to
do with LAMMPS?
if your system is not neutral, then there may be a problem,
as the FFT will disregard the G=0 component in reciprocal
space which is equivalent to placing a neutralizing counter
charge in the environment around the simulation cell.
i've never tested this for the kind of system setup that you
describe, but there is a potential issue. my suggestion is
to build a minimal testbed environment with only a small
number of charges and test a neutral vs. a charged test
system and whether there is a visible difference.
it may be, that you require a real 2d ewald implementation,
but i am not 100% about that.
cheers,
axel.
Yes, i have just tested neutral vs charged system. The force per number of charged particles is the same in case of neutral system and charged one. So, it seems that the force which is calculated is stable and is not a mistake. But it is about 1 order larger than it should be for a charged film. I can’t understand this.
well, it has been too long since i worked on these issues to tell you off
the top of my head, where the problem is.
i have one more (simple?) recommendation. why don't you take your
system and compute the forces with run 0 but without periodicity and
then do a sequence of these calculations where you use the replicate
command to replicate the system in x and y and then see if you can
extrapolate to infinity before the calculations become too expensive?
unless you have a pathological case, that should get you to the right
ballpark and provide a hint where to look for the problem.
cheers,
axel.
Sounds to me like the calculation should be OK the way you’ve set it up. At this point, I’d recommend stepping back to a much more simple case that you can calculate analytically, like the force acting between a couple of isolated point charges at a fixed distance. Try it first with no periodicity and no kspace. Compare what you get analytically with what LAMMPS gives you. Then, you could try adding in a small number of charged particles and compare LAMMPS vs what you calculate manually. You could then try with 3D periodicity and kspace. Then the slab geometry option. I’d recommend maintaining net neutrality in all cases if possible.
Paul
yes, I see now that probably LAMMPS do it in correct way. I was wrong.