Hello LAMMPS users,
I am running into issues when using the distance-dependent dielectric (coul/diel) correction in LAMMPS. Any guidance you might provide would be much appreciated.
Specifically, I am seeking simulate explicit ions models in LAMMPS. I’ve been attempting to reproduce the results of Lenart et al (dx.doi.org/10.1063/1.2431169 ), using the parameters supplied in their paper and pair_style hybrid/overlay with the coul/long, coul/diel, lj/cut, and gauss/cut. When I use the distance-dependent dielectric, the system of 1 M NaCl explodes. Without the distance-dependent dielectric, the salt behaves as a liquid.
I’ve also seen odd behavior also when using the parameters of Freeman et al.(dx.doi.org/10.1063/1.3652956). With their parameters, use of the distance-dependent dielectric causes the ions to begin crystallize. Without the distance-dependent dielectric, the system behaves as a liquid.
I’ve compared the calculated forces from LAMMPS with the equation for the dielectric correction and found them to be in agreement. To me this suggests that I am making a mistake in my input file. I’m reasonably sure that the problem is not in the pair_coeff’s. I’ve reproduced the relevant parts of the input file below for the Lenart forcefield (Na=1, Cl=2).
Does anyone out there have experience with the pair_style coul/diel? If so, did you see similar behavior?
Am I using the pair_style hybrid/overlay correctly?
Thank you for your time and help,
Dan
dielectric 78
pair_style hybrid/overlay coul/long 50.0 gauss/cut 18.0 lj/cut 18.0 coul/diel 12
atom_modify map array sort 0 0.0
#Coulombic interactions
pair_coeff * * coul/long
kspace_style pppm 1.0e-4
Lennard-Jones attraction
pair_coeff 1 1 lj/cut 0.02847 2.443 # Na±Na+
pair_coeff 1 2 lj/cut 0.08427 2.796 # Na±Cl-
pair_coeff 2 2 lj/cut 0.23401 3.4872 # Cl–Cl-
Hydration barrier
pair_coeff 1 1 gauss/cut 8.14 3.1 0.57
pair_coeff 1 2 gauss/cut 8.14 3.1 0.57
pair_coeff 2 2 gauss/cut 8.14 3.1 0.57
Distance-dependent dielectric
pair_coeff 1 1 coul/diel 78 2.7 0.57
pair_coeff 1 2 coul/diel 78 2.7 0.57
pair_coeff 2 2 coul/diel 78 2.7 0.57
neighbor 2.0 multi
neigh_modify check yes
timestep 20.0
velocity all create {T} {random} rot yes mom yes dist gaussian
fix 1 all langevin {T} {T} 1000 ${random}
fix 2 all nve
</input file snippet>