[lammps-users] RDF calculation for a ReaxFF simulation


I am currently testing the O-H parameterizations, and I would like to
produce a G(r) for my simulation.

My input file is the following:

# Water

units real

atom_style charge
read_restart restart.C

pair_style reax 2.12 1.0e-6
pair_coeff * * ffield.reax 3 2

neighbor 2.0 bin
neigh_modify every 1 delay 0 check no

compute OOrdf all rdf 250 1 1
compute OHrdf all rdf 250 1 2
compute HHrdf all rdf 250 2 2

fix 1 all temp/berendsen 300.0 300.0 100.0
fix 2 all nve
fix 3 all reax/bonds 400 h2o.bonds
fix 4 all ave/time 100 1 100 c_OOrdf file OOrdf mode vector
fix 5 all ave/time 100 1 100 c_OHrdf file OHrdf mode vector
fix 6 all ave/time 100 1 100 c_HHrdf file HHrdf mode vector

thermo 40
thermo_style custom step cpu temp press etotal
thermo_modify format float %19.8f

dump 1 all custom 40 h2o.dump id type q x y z
dump_modify 1 format "%d %d %10.5f %10.5f %10.5f %10.5f"

restart 400 restart.A restart.B

timestep 0.25
run 40000

Where entry labeled 3 next to 'ffield.reax' is Oxygen, and 2 is Hydrogen.

The problem that I have is that the function is not converging to 1,
but to 1.5. The LJ-ish parameter in ReaxFF that determines long
distance cutoffs, is the bond order cutoff term, which is set to
0.100. THis should correspond to 10 Angstroms, from how this
parameter is treated in the code, as it gets multiplied by 100. The
simulation I am running is in a cubic box with 1029 atoms (343 waters,
7^3, i.e.). A side is 21.729 Ang.

Thanks for any help / suggestions!


I'm not sure what you're asking. The compute rdf command
knows nothing about ReaxFF. It can just compute an rdf out
to the cutoff distance. If you're saying that isn't at the value
you expect, maybe the distance isn't long enough. But you'll
have to post-process an RDF if you want it out to very long