Dear LAMMPS Developers and Users,

I have a 3d system of 500 charged particles interacting via Yukawa Pair potential. I am using the NVT thermostat (using ‘fix nvt’ command) to equilibrate the system at desire temperature. I am calculating the total potential energy and pair energy (both are equal for this case, as expected) using the following commands:

pair_style yukawa {kappa} {cutf}

pair_coeff * * 3.2891e-20

thermo_style custom step temp pe epair

thermo 500

Beside that, I have also dumped all the positions and velocities by following command:

dump dp1 all custom 500 dump_nvt.lammpstrj id type x y z vx vy vz

By knowing the positions of all the particle at a given time step, I can calculate total pair energy of the system at that time step by calculating

sum_(i,j) [(A./r_{ij}).*exp(-{kappa}.*r_{ij})] ,

With the condition " r_c(cutf) >= r_{ij} > 0 "

where r_i = sqrt((x_i).^2 + (y_i).^2 + (z_i).^2) and r_{ij} = abs(r_j - r_i)

I expected potential energies from these two different methods should be same or might be different by a constant normalization factor, which will be indipendent of ‘kappa’. But I am not getting this. Am I missing any conceptual issue by doing the calculation of pair potential energy?

Any help or suggestion will be appreciable.

I am using ‘lammps-17Nov16’ version of LAMMPS.

Thank you in advance.

With regards

Srimanta Maity

Institute for Plasma Research, Gandhinagar, India

Dear LAMMPS Developers and Users,

I have a 3d system of 500 charged particles interacting via Yukawa Pair

potential. I am using the NVT thermostat (using 'fix nvt' command) to

equilibrate the system at desire temperature. I am calculating the total

potential energy and pair energy (both are equal for this case, as expected)

using the following commands:

pair_style yukawa \{kappa\} {cutf}

pair_coeff * * 3.2891e-20

thermo_style custom step temp pe epair

thermo 500

Beside that, I have also dumped all the positions and velocities by

following command:

dump dp1 all custom 500 dump_nvt.lammpstrj id type x y z vx

vy vz

By knowing the positions of all the particle at a given time step, I can

calculate total pair energy of the system at that time step by calculating

sum_(i,j) [(A./r_{ij}).*exp(-{kappa}.*r_{ij})] ,

With the condition " r_c(cutf) >= r_{ij} > 0 "

where r_i = sqrt((x_i).^2 + (y_i).^2 + (z_i).^2) and r_{ij} = abs(r_j -

r_i)

1) Is this a periodic system? Did you remembering to consider

periodic images in the summation above?

I expected potential energies from these two different methods should be

same or might be different by a constant normalization factor, which will be

indipendent of 'kappa'. But I am not getting this. Am I missing any

conceptual issue by doing the calculation of pair potential energy?

2) One suggestion would be to compare differences between the energies

of two different conformations, rather than the absolute energy. You

can add a constant to the energy without effecting the physics of the

simulation.

3) Are you using pair_modify shift? (this will add a constant to the

energy of the system)

4) You can ask LAMMPS to generate a table of energy (and force) as a

function of inter-particle spacing (r) using the "pair_write" command.

5) This does not effect you, but (generally speaking) if you are using

a pair_style with a long range electrostatic contribution, then you

have to consider the contributions from periodic images carefully.

That's what comes to my mind.

I don't know if this helped at all.

Cheers

Andrew

Thank you for your suggestions. Yes, I have missed to point out that I am using periodic boundary condition in my simulation. Here I am using Yukawa pair potential, which is a short-range interaction limited up to the distance r = r_c (force cutf). So in this particular case, I wonder is it necessary to consider periodic images in the summation?

No, I am not using “pair_modify shift” command. But I have used the “pair_write” command to get the profile ® of pair energy which is exactly matching with the analytical expression [(A./r).*exp(-{kappa}.*r)]. But it does not give me the total pair potential energy of the entire system.

Thank you

With regards

Srimanta

Thank you for your suggestions. Yes, I have missed to point out that I am using periodic boundary condition in my simulation. Here I am using Yukawa pair potential, which is a short-range interaction limited up to the distance r = r_c (force cutf). So in this particular case, I wonder is it necessary to consider periodic images in the summation?

yes, absolutely. please grab your MD text book and look it up.

axel.