I used the basic code below to generate a PE curve for a H2 molecule. The results are promising with a De of 0.174 Ha and a far asymptotic value of -1.0 Ha. However, when I compute the number of bound vibrational states, I find 21 rather than the 14 expected from experimental measurements.
To improve the accuracy of the PE curve, should I do a fragment calculation for large values of R? Would the total spin be zero? Do the spin states become uncorrelated at large distances? How would I blend the fragment calculation with a molecule calculation?
thank you very much for your advice,
import numpy as np
mol_tmpl = """H H 1 **R**""" # rvals = np.linspace(0.11,0.50,40) # rvals = np.append(rvals,0.55) rvals = np.linspace(0.6,5.0,45) fp = open('h2_pes.dat', 'w') for r in rvals: mol = psi4.geometry(mol_tmpl.replace("**R**", str(r))) energy = psi4.energy("ccsd/aug-cc-pVQZ", molecule=mol) print(r,energy, file = fp) fp.close()
The results are promising with a De = 0.174 Ha and a far asymptotic value of -1.0 Ha.