molecule h2 {
  symmetry c1
  H
  H 1 R
  }

range=[0.3+0.1*x for x in range(25)]

set {
  basis 6-31G
  num_roots 3
}

for R in range:
  h2.R=R
  en=energy('fci')
  en1=psi4.get_variable("CI ROOT 1 TOTAL ENERGY")
  en2=psi4.get_variable("CI ROOT 2 TOTAL ENERGY")
  print(R, en, en1, en2)

