#this is the file to calculate vacancy formation energy of carbon fiber
#-----------------INITIALIZATION-----------------
clear
units metal
dimension 3
boundary p p p
atom_style charge
#-----------------ATOM DEFINITION----------------
read_data data
#-----------------FORCE FIELD--------------------
pair_style comb3 polar_off
pair_coeff * * ffield.comb3 H
neighbor 2.0 bin 
neigh_modify delay 10 check yes
mass 1 1
#---------------------------Settings----------------------------
compute eng all pe/atom 
compute eatoms all reduce sum c_eng
#----------------------Run Minimization-------------------------
reset_timestep	0

thermo 10
thermo_style custom step pe lx ly lz press pxx pyy pzz c_eatoms 

dump 1 all custom 200 dump.relax.1.* id type xs ys zs c_eng 
fix 1 all qeq/comb 10 1e-4
min_style cg
minimize 1e-15 1e-15 5000 5000

run 0
undump 1

variable N equal count(all)
variable No equal $N
variable E equal "c_eatoms"
variable Ei equal $E

#---------------------------------------------------------------

######################################
# SIMULATION DONE
print "All done"
print "Total number of atoms = ${No}"
print "Initial energy of atoms = ${Ei}"