Dear LAMMPS developers and users!
I repeatedly encounter a problem using Berendsen press-conserving procedure with PETN described by reax/c forcefield. Here is an input file and log for a small test system, the one we are interested in (32x4x4 unit cells) shows the same behaviour. Basically, pressures and volume become NaN right after fix is applied.
Input file:
units real
atom_style charge
variable n_x equal 3
variable n_y equal 1
variable n_z equal 1
variable temperature equal 300
variable temperature_x2 equal ${temperature}*2
log log.petn.reaxc.T${temperature}K.{n_x}x{n_y}x${n_z}.txt
read_data data.reax
replicate 3 1 1
velocity all create ${temperature_x2} 9999
pair_style reax/c NULL
pair_coeff * * ffield.reax 1 2 3 4
timestep 0.1
fix 1 all nve
fix 2 all qeq/reax 1 0.0 10.0 1.0e-6 reax/c
fix 3 all reax/c/bonds 100 1 100 bonds.petn.reaxc.T${temperature}K.{n_x}x{n_y}x${n_z}.txt
dump d1 all cfg 1000 cfg.petn.T${temperature}K.{n_x}x{n_y}x${n_z}.*.cfg id type xs ys zs
dump_modify d1 element C H O N
dump d2 all custom 1000 dump.petn.T${temperature}K.{n_x}x{n_y}x${n_z}.*.dump id type x y z
dump_modify d2 element C H O N
thermo 100
thermo_style custom step temp ke pe pxx pyy pzz etotal vol
run 1000
fix 4 all press/berendsen aniso 0.0 0.0 1000.0
thermo 1
run 11000
unfix 4
thermo 100
compute my_ke all ke/atom
compute my_press all stress/atom
restart 10000 restart_petn_T${temperature}K.{n_x}x{n_y}x${n_z}.*
run 10000