Hi,
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!
Camilo