Proble of reax/c with NVT and NVE

Dear All,
I am trying to simulate the graphene with reax/c implemented in Lammps. But I have encountered some difficulties. Two steps of simulations are wanted to be carried: NVT, then NVT. But I find that the temperature is rising unreasonably and the size of system is expanding unlimitedly. The input and output files are in the attachment. Anyone can give some suggestions?

Thanks in advance!
Wang Jian

The input script is as follows

units real
dimension 3
newton on
boundary s p s
atom_style charge
neighbor 0.3 nsq
neigh_modify check yes
read_data graph_reaxc.xyz
velocity all create 500 1212127641 mom yes rot yes dist gaussian units box

pair_style reax/c NULL checkqeq no
pair_coeff * * ffield.reax.cho 1

fix NVT all nvt temp 500 500 1
thermo_style custom step temp etotal vol
thermo_modify lost warn
thermo 10000

Run

timestep 0.5
dump data1 all custom 10000 outputnvt.xyz type x y z
run 100000

unfix NVT
fix nve all nve
dump data2 all custom 10000 outputnve.xyz type x y z
run 100000

The structure data is as follows.

2360 atoms
1 atom types

-11.541138 13.053984 xlo xhi
-125.433333 125.906667 ylo yhi
-1.675000 1.675000 zlo zhi

Masses

1 12

Atoms

1 1 0.000000 2.459512 -124.960000 0.000000
2 1 0.000000 4.919024 -124.960000 0.000000
3 1 0.000000 7.378536 -124.960000 0.000000
4 1 0.000000 9.838049 -124.960000 0.000000
5 1 0.000000 12.297561 -124.960000 0.000000
6 1 0.000000 -9.838049 -124.960000 0.000000
7 1 0.000000 -7.378536 -124.960000 0.000000
8 1 0.000000 -4.919024 -124.960000 0.000000

Part of log file is

LAMMPS (17 Oct 2012)
……
memory usage per processor = 1.90524 Mbytes
Step Temp TotEng Volume
0 500 -422493.8 3.9355046
10000 496.96668 -447001.31 34366.804
20000 511.68336 -446813.17 52925.553
30000 489.01261 -447035.91 52192.758
40000 494.57279 -447055.84 49979.958
50000 504.45062 -446935.76 52020.052
60000 516.99174 -446829.53 55296.869
70000 511.90547 -446981.81 60563.274
80000 493.29511 -447077.87 65691.772
90000 497.178 -446963.27 69837.58
100000 494.79137 -447110.48 67269.375
Loop time of 2077.01 on 4 procs for 100000 steps with 2360 atoms

Pair time () = 1987.14 (95.6728) Neigh time () = 12.6195 (0.607578)
Comm time () = 58.2898 (2.80642) Outpt time () = 0.036369 (0.00175102)
Other time (%) = 18.9314 (0.911474)

Nlocal: 590 ave 595 max 584 min
Histogram: 1 0 0 0 0 1 1 0 0 1
Nghost: 196.25 ave 200 max 190 min
Histogram: 1 0 0 0 0 0 1 0 0 2
Neighs: 32088.8 ave 32362 max 31777 min
Histogram: 1 0 0 0 0 2 0 0 0 1

Total # of neighbors = 128355
Ave neighs/atom = 54.3877
Neighbor list builds = 6838
Dangerous builds = 53

unfix NVT

fix nve all nve
dump data2 all custom 10000 outputnve.xyz type x y z

run 100000
Memory usage per processor = 1.90524 Mbytes
Step Temp TotEng Volume
100000 494.79137 -447110.48 67647.011
110000 1792.6493 -427406.65 63336.634
120000 19927.294 -13282.889 29380787
130000 20029.52 82340.806 7.0560577e+08
140000 20103.372 89890.917 2.5201139e+09
150000 20106.79 92315.111 5.4685969e+09
160000 20113.849 93225.42 9.5450081e+09
170000 20136.624 93850.547 1.4749347e+10
180000 20113.54 94312.333 2.1081615e+1

in_script.zip (360 KB)

Could be many things. Too big a timestep,
insufficient initial equilibration, a bad model or geometry,
etc.

Steve

Jian Wang,

0.5 is a pretty large time step for the cho parameter set. I’d try 0.1 fs.

I’ve recently published a paper detailing some of the properties of the cho parameter set with reax that might be helpful: http://dx.doi.org/10.1021/ct300491d

-Ben

Benjamin Jensen
PhD Candidate
Department of Mechanical Engineering - Engineering Mechanics
Michigan Technological University
Houghton, MI