Dear all,
we are currently trying to simulate thermal oxidation of silicon nanowires using ReaxFF and LAMMPS 10 Aug 2015.
During our simulation we discovered that oxygen molecules far away from the nanowire become negatively charged.
This results in long-ranging repulsive behaviour with oxygen atoms attached to the wire, which leads to an unphysical slow-down of the oxidation process.
We were able to replicate this unusual charge redistribution using the following input file. Si.reaxFF is a file containing parameters by either Buehler et al. (PRL 96, 2006) or Kulkarni et al. (J. Phys. Chem. C 117, 2013). Both yield the abnormal charges even though the exact numbers differ from each other.
---------------------- Input File ----------------------
units real
atom_style full
boundary p p p
lattice diamond 5.431 origin 0.1 0.1 0.1 orient x 1 0 0 orient y 0 1 0 orient z 0 0 1
region box block 0 50 0 50 0 50 units lattice
create_box 2 box
region silicon block 0 0.1 0 0.1 0 0.1 units lattice
create_atoms 2 region silicon
region oxygen block 24 25 24 25 24 25 units lattice
create_atoms 1 region oxygen
mass 1 28.06
mass 2 15.999
pair_style reax/c NULL
pair_coeff * * Si.reaxFF Si O
neighbor 2 bin
neigh_modify every 10 delay 0 check no
compute type all property/atom type
timestep 0.05
fix qeq all qeq/reax 1 0.0 10.0 1.0e-6 reax/c
fix nve all nve
dump custom all custom 1 output.custom id c_type q x y z vx vy vz
run 0
---------------------- Output File ----------------------
ITEM: TIMESTEP
ITEM: NUMBER OF ATOMS
9
ITEM: BOX BOUNDS pp pp pp
0 271.55
ITEM: ATOMS id c_type q x y z vx vy vz
1 2 -0.328226 0.5431 0.5431 0.5431 0 0 0
2 1 0.0719583 130.887 130.887 130.887 0 0 0
6 1 -0.0070383 132.245 132.245 132.245 0 0 0
5 1 0.0370514 133.603 133.603 130.887 0 0 0
9 1 0.0507171 134.96 134.96 132.245 0 0 0
4 1 0.0370514 133.603 130.887 133.603 0 0 0
8 1 0.0507171 134.96 132.245 134.96 0 0 0
3 1 0.0370514 130.887 133.603 133.603 0 0 0
7 1 0.0507171 132.245 134.96 134.96 0 0 0
As you can see the charges are redistributed instantaneously. There seems to be no spatial cut-off during the charge equilibrium (QEq).
We already checked the available documentation and mailing list and found a similar problem described in http://lammps.sandia.gov/threads/msg29916.html. Unfortunately, no solution was available back in 2012. Is there a way to fix this problem in 2017?
Kind regards,
Georg Heinze & Florian Fuchs