coul/cut and coul/long gives the same force/energy for like charge and opposite charge interactions

Hi everyone!

I’m a new user of Lammps, but I’ve searched the mail list and the documentation extensively for my problem so hopefully my problem hasn’t been asked before (although I suspect it might still be a RTFM-issue).

Trying simulate liquid NaCl I noticed that even below the boiling temperature (p = 1 bar) the radial distribution function was that of a gas. The same thing happened for both

pair_style born/coul/long 7.6 15.0 (I used a previous answer from the mailing list that treated NaCl and modified it)

or using a hybrid/overlay pair_style with the repulsive part being a table taken from Kittels Introduction to Solid State Physics (page 62 of 8th edition paperback) and a coul/long term.

I used pair_write from the below test code to see what kind of potential was being calculated (version: LAMMPS 64-bit 15Apr2020​).

Setup

dimension 3
units metal
boundary p p p
atom_style charge

Na atoms

lattice fcc 5.64 origin 0 0 0

region box block 0 10 0 10 0 10
create_box 2 box
create_atoms 1 box

Cl atoms

lattice fcc 5.64 origin 0.5 0 0

create_atoms 2 box

Other settings

mass 1 28.990
mass 2 35.453
set type 1 charge +1.0
set type 2 charge -1.0

Potentials

pair_style coul/cut 10.0

pair_coeff 1 1
pair_coeff 1 2
pair_coeff 2 2

pair_write 1 1 1000 r 0.2 10.2 table1.txt coul11

pair_write 1 2 1000 r 0.2 10.2 table2.txt coul12

Plotting the energy and force values for the tables gave the same exact curves and both were always positive. I also tried the same with born/coul/long, and coul/long but the plotted values were still positive.

In the documentation for pair_write it states “If the pair potential is computed between charged atoms, the charges of the pair of interacting atoms can optionally be specified. If not specified, values of Qi = Qj = 1.0 are used.” I’m not sure if the output from pair_write only takes into account Qi=Qj = 1.0 interactions even though I set the charges in the code or if in fact the potential is always positive, no matter the charges in the pair, which might explain why the system does not form a liquid even below the boiling point. Am I missing something?

Best regards,
William Green

William,

the pair write command does not look at the assigned charges but only at atom types. while in your case those coincide with atom type assignments, that need not be, so the pair_write command cannot know and thus uses the charges of +1 for both atom types. the outcome is different, if you do a real computation via “run 0”

to get the expected outcome you would have to use (BTW no need to use separate files, new tables get appended):
pair_write 1 1 1000 r 0.2 10.2 table.txt coul11 1.0 1.0
pair_write 1 2 1000 r 0.2 10.2 table.txt coul12 1.0 -1.0

pair_write 1 2 1000 r 0.2 10.2 table.txt coul22 -1.0 -1.0

axel.

Hi Axel,

Thanks a lot for the quick reply! Using the arguments you gave in the pair_write command gave correct behavior for the coulomb plots and the rapid expansion of the system upon changing

set type 1 charge +1.0
set type 2 charge -1.0

to

set type 1 charge +1.0
set type 2 charge +1.0

​showed that my previous code most certainly respects the charge of the atoms.

Thanks again.

Best regards,

William Green