Effect of Ewald/pppm splitting parameter on total energy of the system

Dear all,

I was trying to set up a MD run for a crystal with PBC in all three dimensions in NVE ensemble with Coulomb interaction. I’ve tried both Ewald sum and pppm method with different splitting parameters. For instance, for a system with 512 atoms in the periodic box with Ewald sum, the simulation results are as follows:

splitting parameter (A) total energy (eV)
2 -3279.79
4 -3279.96
6 -3279.7
9 -3278.55
10 -3278.9
15 -3279.75
18 -3279.75
23 -3277.28
The total energies of the system don’t seem well converged although in principle, the total energy should not vary with splitting parameters. Can anyone tell me what’s wrong with this? Or is it normal that the total energy fluctuates this way? Here is part of my input:

units metal
dimension 3
atom_style charge
boundary p p p

ATOM DEFINITION

read_data data.pos
neighbor 2 bin

---------- Define Interatomic Potential ---------------------

pair_style lj/cut/coul/long 10.0 18.0
kspace_style ewald 1.0e-4
pair_coeff 1 1 0.0026 3.29
pair_coeff 2 2 0.0039 4.36
pair_coeff 1 2 0.0032 3.825

timestep 0.001
thermo 100
velocity all create 300 1233445 dist gaussian mom yes rot yes

fix 3 all nve
run 100000
unfix 3

Thanks a lot!

Regards,
Jihong

how many standard deviations are the values apart? you should also check that, not only the mean!