# Import the KohnSham module from ATK
from ATK.KohnSham import *
import datetime
# Load molecule from VNLFile
h2o = VNLFile("h2o.vnl").readAtomicConfigurations()["h2o"]

print '    MeshCutoff (Rydberg)     TotalEnergy (eV)'
print '------------------------------------------------'
use_old=False
for mcut in range(30,200,10)*Rydberg:
    el = electronDensityParameters(mesh_cutoff=mcut)
    method = KohnShamMethod(electron_density_parameters=el)
    if use_old:
        density = executeSelfConsistentCalculation(
            method=method,
            atomic_configuration=h2o,
            initial_calculation=old_scf
            )
    else:
        density = executeSelfConsistentCalculation(
            method=method,
            atomic_configuration=h2o,
            )
        use_old=True
    old_scf=density
    print '\t', mcut.inUnitsOf(Rydberg), '\t\t\t', \
          calculateTotalEnergy(density).inUnitsOf(eV)