# GAMESS-UK C2H2 tautomerisation DFT (BLYP) reaction path optimisation
# Step 1: minimise the endpoints

# Initialise ChemShell
from chemsh import *

# Read in reactant and product state geometries from disk
rs = Fragment(coords="reactant.xyz")
ps = Fragment(coords="product.xyz")

# Set up the QM calculation
rs_theory = GAMESS_UK(frag=rs, method="dft", functional="blyp", basis="3-21g")
ps_theory = GAMESS_UK(frag=ps, method="dft", functional="blyp", basis="3-21g")

# Optimise the endpoints using DL-FIND
optReactant = Opt(theory=rs_theory, tolerance=0.0003)
optReactant.run()

optProduct = Opt(theory=ps_theory, tolerance=0.0003)
optProduct.run() 

# Save the optimised structures for the NEB run
rs.save("reactant_opt.pun", "pun")
ps.save("product_opt.pun", "pun")
