Hi Lammps users,

I want to use reaxff potential to simulate 2d BCN. When I set dimension 2, bond start to break at the first step of minimization. When I set dimension 3, bond are not break, but after I give the initial velocity 10k, and run NVE, the temperature keep increasing over 1000k. Can anyone help me to check that where I got wrong?

Any help would be appreciated.

Sincerely,
Peng

Here is my dimension 3 input script

units real

processors 2 2 2

# variable

variable b equal 4.26

variable a equal \$b*sqrt(3.0)

variable 1_3 equal 1.0/3.0

variable 2_3 equal 2.0/3.0

variable 5_6 equal 5.0/6.0

variable 1_6 equal 1.0/6.0

# Setup Box

dimension 3

boundary p p p

atom_style charge

lattice custom 1.0 &

a1 \$a 0.0 0.0 &

a2 0.0 \$b 0.0 &

a3 0.0 0.0 10.0 &

basis 0.0 0.0 0.0 &

basis \${1_3} 0.0 0.0 &

basis \${2_3} 0.0 0.0 &

basis 0.0 \${2_3} 0.0 &

basis {1_6} {1_6} 0.0 &

basis 0.5 \${1_6} 0.0 &

basis {5_6} {1_6} 0.0 &

basis \${1_6} 0.5 0.0 &

basis 0.5 0.5 0.0 &

basis \${5_6} 0.5 0.0 &

basis {1_3} {2_3} 0.0 &

basis {2_3} {2_3} 0.0 &

region box block 0 10 0 10 -0.25 0.25

create_box 3 box

create_atoms 1 box &

basis 2 2 &

basis 3 3 &

basis 4 2 &

basis 5 3 &

basis 7 2 &

basis 8 3 &

basis 9 2 &

basis 12 3 &

group B type 1

group N type 2

group C type 3

mass 1 10.811

mass 2 14.0067

mass 3 12.0107

# Settings

pair_style reax/c NULL

pair_coeff * * ffield.reax.new B N C

fix 1 all qeq/reax 1 0.0 10.0 1.0e-6 reax/c

compute reax all pair reax/c

timestep 0.01

thermo 1

thermo_style custom step temp press pe ke etotal vol density epair lx ly lz

dump 99 all custom 1 dump.reax.atom id type mass xs ys zs

# Minimization

min_style cg

minimize 1.0e-10 1.0e-10 1000000 1000000

# Heating Stage

velocity all create 10.0 12345 rot yes dist gaussian

fix 2 all nve

run 50000

Here is my dimension 2 input script

units real

processors 2 2 2

# variable

variable b equal 4.26

variable a equal \$b*sqrt(3.0)

variable 1_3 equal 1.0/3.0

variable 2_3 equal 2.0/3.0

variable 5_6 equal 5.0/6.0

variable 1_6 equal 1.0/6.0

# Setup Box

dimension 2

boundary p p p

atom_style charge

lattice custom 1.0 &

a1 \$a 0.0 0.0 &

a2 0.0 \$b 0.0 &

basis 0.0 0.0 0.0 &

basis \${1_3} 0.0 0.0 &

basis \${2_3} 0.0 0.0 &

basis 0.0 \${2_3} 0.0 &

basis {1_6} {1_6} 0.0 &

basis 0.5 \${1_6} 0.0 &

basis {5_6} {1_6} 0.0 &

basis \${1_6} 0.5 0.0 &

basis 0.5 0.5 0.0 &

basis \${5_6} 0.5 0.0 &

basis {1_3} {2_3} 0.0 &

basis {2_3} {2_3} 0.0 &

region box block 0 10 0 10 -0.25 0.25

create_box 3 box

create_atoms 1 box &

basis 2 2 &

basis 3 3 &

basis 4 2 &

basis 5 3 &

basis 7 2 &

basis 8 3 &

basis 9 2 &

basis 12 3 &

group B type 1

group N type 2

group C type 3

mass 1 10.811

mass 2 14.0067

mass 3 12.0107

fix 88 all enforce2d

# Settings

pair_style reax/c NULL

pair_coeff * * ffield.reax.new B N C

fix 1 all qeq/reax 1 0.0 10.0 1.0e-6 reax/c

compute reax all pair reax/c

timestep 0.01

thermo 1

thermo_style custom step temp press pe ke etotal vol density epair lx ly lz

dump 99 all custom 1 dump.reax.atom id type mass xs ys zs

# Minimization

min_style cg

minimize 1.0e-10 1.0e-10 1000000 1000000

# Heating Stage

velocity all create 10.0 12345 rot yes dist gaussian

fix 2 all nve

run 50000

Please first visualize your structure. If your structure has bad geometry then there will be high energy, which will produce large forces, which result in significant accelerations, which impart high particle velocities, which in-turn give temperatures shooting up or boxes blowing apart. This type of behavior has been addressed many times and is one of the most common problems that hits the mailing list. Is there a particular reason you want to simulate in 2d? Just because a single BCN layer is considered 2d does not mean that it needs to be simulated in a 2d schemeā¦ What are you actually trying to study? What have other people done?

Best,
Dyaln