# Gas phase CO2 geometry optimisation

# Initialise ChemShell
from chemsh import *
from chemsh.utils import testutils

# Import initial CO2 geometry
co2 = Fragment(coords="co2.xyz")

# QM settings
qm_theory = NWChem(frag=co2,
                   method='dft',
                   functional='b3lyp',
                   basis='co2_nwchem.basis',
                   path="mpirun -np 4 nwchem")

# Run QM optimisation
opt = Opt(theory=qm_theory,
          algorithm='lbfgs', trust_radius='energy',
          maxcycle=200, maxene=300) 
opt.run()

ecalc = opt.result.energy

# Save optimised geometry
co2.save('co2_opt.xyz', 'xyz')
co2.save('co2_opt.pun', 'pun')

# Validate result against pre-calculated reference energy
eref = -188.446738574825
tolerance = 1.0e-7

testutils.validate(ecalc, eref, tol=tolerance)
