memory 4 gb molecule monomer { 0 1 # charge and multiplicity. units angstrom Ar 0.0 0.0 0.0 # Unfortunately, Psi4 1.x does not write out the orbital energies # correctly if symmetry is used. So we disable it here by forcing C1. # This will have the consequence that the calculation will proceed # more slowly, and, for symmetric systems, the calculation may even # fail to converge. If symmetry is enabled then the e_HOMO and AC-shift # values may be incorrect and these should be checked and corrected manually. Symmetry C1 } # Define the XC functional here dft_func = 'pbe0' set basis aug-cc-pVTZ set freeze_core TRUE # get the charge and multiplicity of the monomer charge = monomer.molecular_charge() mult = monomer.multiplicity() set reference rks EN, wfn = energy(dft_func, return_wfn=True) # We need to find the HOMO eigenvalue of this system nelec_a = wfn.nalpha() # The epsilon_a list is not accessible like this in Psi4 1.3.2 and later versions: # epsilon_a = wfn.epsilon_a() # e_homo = epsilon_a[nelec_a-1] # So we set e_homo = 0. This means that the AC-shift must be computed manually e_homo = 0.0 # now remove and electron and increase the multiplicity. monomer.set_molecular_charge(charge+1) monomer.set_multiplicity(mult+1) set reference uks EN1 = energy(dft_func) IP = EN1 - EN AC_shift = IP + e_homo psi4.print_out("\n") psi4.print_out("%s energies and IP \n\n" % (dft_func)) psi4.print_out(" Num Electrons Charge E_int [a.u.] \n") psi4.print_out("-----------------------------------------------------\n") psi4.print_out(" N %3d %10.6f \n" % (charge,EN)) psi4.print_out(" N-1 %3d %10.6f \n" % (charge+1,EN1)) psi4.print_out("-----------------------------------------------------\n") psi4.print_out(" IP %10.6f a.u. \n" % (IP)) psi4.print_out(" e_HOMO %10.6f a.u. \n" % (e_homo)) psi4.print_out(" AC-shift %10.6f a.u. \n" % (AC_shift)) psi4.print_out("-----------------------------------------------------\n")