ReaxFF simulation of water

LAMMPS LAMMPS Mailing List Mirror
1

Hello,
I am using reax force field for water molecules. I have taken parameters for water from the reference , 2017, 121 (24), pp 6021–6032.

I am facing the following problems.
once the simulation enters 1fs it’s showing the following warning and error message
WARNING: Fix qeq/reax CG convergence failed after 200 iterations at 151119 step (…/fix_qeq_reax.cpp:741)

ERROR: Lost atoms: original 12600 current 3862 (…/thermo.cpp:441)
Last command: run 50000

My system consists of 4120 water molecules which I calculated using the number density and box dimension is 505050 Angstrom cube
The following are my input script,qeq parameters are given below

Input script

units real

newton on
boundary p p p

atom_style charge
read_data data.waterchrg

pair_style reax/c lmp_control2
pair_coeff * * ffield_reaxh2o H O

compute reax all pair reax/c

variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]

neighbor 2 bin
neigh_modify every 10 delay 0 check no

fix 2 all qeq/reax 1 0.0 10.0 1e-6 param.qeq

thermo 1000
thermo_style custom step temp epair etotal press v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq

velocity all create 300.0 4712398

dump 1 all cfg 1000 *.cfg mass type xs ys zs
dump_modify 1 element H O

fix 3 all reax/c/species 1 10 1000 speciesRT.out element H O

dump 5 all xyz 1000 dump_reaxwatRT

timestep 0.1
fix 11 all nve
run 50000
unfix 11

timestep 0.25
fix 11 all nve
run 50000
unfix 11

timestep 0.50
fix 11 all nve
run 50000
unfix 11

timestep 1.00
fix 11 all nve
run 50000
unfix 11

qeq parameters are given below

1 5.8678 14.0000 0.9000 for Hydrogen
2 8.5000 17.9978 1.0500 for Oxygen

Water molecule was packed using packmol and data file was generated using vmd.
I was facing the same issue even after using xyz coordinates taken by equilibrating water molecules.

Expecting a solution to this issue.

Thank you

2

I am using lammps-22august-2018 version

3

Why are you increasing time step size? Use no larger than 0.2 fs for ReaxFF.

Ray

4

Thank You for your reply Ray Shan,

Is it because ReaxFF have charge equilibration, we have to use smaller time steps? Is there any option to find the optimum time step for the system other than trial and error?

Regarding qeq parameters, how will I choose those parameters as in most of the journals supporting information, only force field parameters are given.

Thank you

5

All non-bonded potentials with light atoms such as H should always use a smaller time step. For ReaxFF the largest timestep size is recommended to be 0.2 fs. Consult a text book (e.g. Allen & Tildesley) on how to determine the optimum time step size.

In the fix qeq/reax command just use “reax/c” instead of assigning a param.qeq file – this is explained on the LAMMPS ReaxFF doc page.

Ray

6

Thank you for the clarification.