# eq: init
units real
atom_style full

dimension 3
boundary p p p

read_data system.data &
  extra/bond/per/atom 30  &
  extra/angle/per/atom 30 &
  extra/dihedral/per/atom 30 &
  extra/improper/per/atom 30 &
  extra/special/per/atom 100
include system.ff
                
#inimize 1e-4 1e-4 1000 1000
# print "minimization done"

velocity all create 300 9953

timestep 1.0

fix sk all shake 0.001 20 0 m 1.008
# fix mom all momentum 1 linear 1 1 1

variable elapsed equal "elapsed"
variable tk equal "temp"
variable pr equal "press"
variable vl equal "vol"
variable dn equal "density"
variable et equal "etotal"
variable ep equal "pe"
variable ek equal "ke"
variable evdwl equal "evdwl"
variable ecoul equal "ecoul"
variable epair equal "epair"
variable eb equal "ebond"
variable eang equal "eangle"
variable edih equal "edihed"

molecule mol1 rxn1_stp1_unreacted.mol
molecule mol2 rxn1_stp1_reacted.mol
molecule mol3 rxn1_stp2_reacted.mol 

thermo 1

dump 1 all custom 1 test_vis.dump id type mol xu yu zu
# dump_modify 1 types labels

fix 1 all nve/limit .01
run 1000
unfix 1

fix myrxns all bond/react stabilization yes statted_grp .01 &
  react rxn1 all 1 0.0 4.0 mol1 mol2 rxn1_stp1.map 
  # react rxn2 all 1 0.0 5.0 mol2 mol3 rxn1_stp2.map

fix 1 statted_grp_REACT nvt temp 300 300 100

fix 4 bond_react_MASTER_group temp/rescale 1 300 300 10 1

thermo_style custom step temp press density f_myrxns[*]

run 10000