Background
In PySCF, one can use from_chkfile() to Read SCF results from PySCF chkfile and transform 1e, 2e integrals using the SCF orbitals. The transformed integrals are then written to FCIDUMP file. It seems that the code supports writing the FCIDUMP file for the case where the calculation is RHF, but not for UHF. For instance, in the following input of a UHF calculation on O++:
mol = pyscf.M(
atom = '''
O
''',
unit = 'angstrom',
basis = {
'O' : parse_gaussian.load('STO-3G.gbs', 'O')
},
charge = 2,
spin = 2,
verbose = 9,
symmetry = True,
output = name +'.txt',
symmetry_subgroup = 'D2h',
max_memory = 4000,
)
mf = mol.UHF().set(conv_tol=1e-10,max_cycle=999,direct_scf_tol=1e-14,chkfile=name+'.chk',init_guess='atom',irrep_nelec={'Ag': 4, 'B3u':1 , 'B2u':1 ,'B1u':0 })
mf.kernel()
pyscf.tools.fcidump.from_chkfile('fcidump', name+'.chk', tol=1e-18, float_format=' %.16g')
calling from_chkfile() with UHF calculation would give the following error:
s = reduce(numpy.dot, (mo_coeff.conj().T, mol.intor_symmetric('int1e_ovlp'), mo_coeff))
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
ValueError: shapes (5,5,2) and (5,5) not aligned: 2 (dim 2) != 5 (dim 0)
Obviously that is because in this case, mo_coeff is a tuple of two ndarrays, which I think corresponds to Alpha and Beta spins.
My Approach
Now I am trying to directly use from_integrals() to convert the given 1e, 2e integrals to FCIDUMP. I think that there would be two sets of integrals each coming from a spin channel. I try to access that with the following:
mol, scf_rec = scf.chkfile.load_scf(name+'.chk')
mo_coeff = numpy.array(scf_rec['mo_coeff'])
orbsym_alpha = symm.label_orb_symm(mol, mol.irrep_id,mol.symm_orb, mo_coeff[0], check=False)
orbsym_beta = symm.label_orb_symm(mol, mol.irrep_id,mol.symm_orb, mo_coeff[1], check=False)
h1ao = scf.hf.get_hcore(mol)
h1e_alpha = reduce(numpy.dot, (mo_coeff[0].T, h1ao, mo_coeff[0]))
h1e_beta = reduce(numpy.dot, (mo_coeff[1].T, h1ao, mo_coeff[1]))
eri_alpha = ao2mo.full(mol, mo_coeff[0], verbose=0)
eri_beta = ao2mo.full(mol, mo_coeff[1], verbose=0)
nuc = mol.energy_nuc()
ms=mol.spin
from_integrals('fcidump_alpha', h1e_alpha, eri_alpha, h1e_alpha.shape[0], mol.nelec, nuc, ms, orbsym_alpha, tol=1e-18, float_format=' %.16g')
from_integrals('fcidump_beta', h1e_beta, eri_beta, h1e_beta.shape[0], mol.nelec, nuc, ms, orbsym_beta, tol=1e-18, float_format=' %.16g')
The full input and output files can be found in Input/Output
Question
I should mention that I am not sure about the approach I am following, I suspect that orbsym should be coming from only one of the spin channels as I have seen in some of the example files found in FCIdump. Thus, my question: How should I handle the orbital symmetry and the converted integrals in the case of UHF calculation so that the final integrals are written into a single FCIDUMP file?
Any suggestions or corrections are appreciated. Thanks!
