# Define looping variables
variable loopcount loop 42
variable latticeconst index 21.562500 20.117733 19.222324 18.571762 18.060498 17.639275 17.281077 16.969474 16.693729 16.446433 16.222262 16.017258 15.828400 15.653329 15.490169 15.337402 15.193782 15.058276 14.930016 14.808265 14.692396 14.581866 14.476204 14.375000 14.306103 14.233591 14.157065 14.076054 13.989999 13.898229 13.799933 13.694111 13.579511 13.454548 13.317157 13.164592 12.993079 12.797235 12.568999 12.295489 11.954184 11.500000 

# Define unit set and class of atomic model
units metal
atom_style atomic

# Periodic boundary conditions along all three dimensions
boundary p p p

# Create a fcc lattice using a single conventional (orthogonal)
# unit cell with a lattice constant from the 'latticeconst' variable defined on line 3 above
lattice fcc ${latticeconst}
region box block 0 1 0 1 0 1 units lattice
create_box 1 box
create_atoms 1 box
mass 1 15.9994

# Specify which KIM Model to use
pair_style kim Sim_LAMMPS_Vashishta_VashishtaKaliaRino_1990_SiO__SM_887826436433_000
pair_coeff * * O

# Set what thermodynamic information to print to log
thermo_style custom step atoms xlo xhi ylo yhi zlo zhi pe press pxx pyy pzz pxy pxz pyz
thermo 10 # Print every 10 steps

# Set what information to write to dump file
dump id all custom 10 output/lammps.dump id type x y z fx fy fz
dump_modify id format line "%d %d %16.7f %16.7f %16.7f %16.7f %16.7f %16.7f"

# Compute the energy and forces for this lattice spacing
run 0

# Define auxiliary variables to contain cohesive energy and equilibrium lattice constant
variable poteng    equal "c_thermo_pe"
variable natoms    equal "count(all)"
variable ecohesive equal "v_poteng/v_natoms"

# Output cohesive energy and equilibrium lattice constant
print "Cohesive energy = ${ecohesive} eV/atom"

# Queue next loop
clear # Clear existing atoms, variables, and allocated memory
next latticeconst # Increment latticeconst to next value
next loopcount # Increment loopcount to next value
jump SELF # Reload this input script with the new variable values