# Gaussian HF single point energy

from chemsh       import *
from chemsh.utils import testutils

water = Fragment(coords="water.cjson")

# default root state: 1 (the first excited state)
gaussian = Gaussian(frag=water, method="tdhf", basis="3-21g", nstates=8)

# TODO: sp = SP(theory=gaussian, frag=water)
sp = SP(theory=gaussian, gradients=True)

sp.run()


ecalc_abs = sp.result.energy
# NB: Python starts with 0!
ecalc_exc = sp.result.excitations['energies'][0]
fcalc_exc = sp.result.excitations['f'][0]

# Test result
# Py-ChemShell/Gaussian16 as of 11/01/2021

# absolute energy (in a.u.) of the S1 (B1) state
eref_abs = -75.232969774400

# excitation energy (in eV) of the S1 (B1) state
eref_exc = 9.5871

# oscillator strength of the S1 (B1) state
fref_exc = 0.0067


testutils.validate(ecalc_abs, eref_abs)
testutils.validate(ecalc_exc, eref_exc)
testutils.validate(fcalc_exc, fref_exc)
