The fnocc
and cc
problems are probably quite different, so let’s go one at a time.
My hunch for the problem in the cc
module is that your Fock matrix is incorrect. The cc
module obtains the Fock matrix when executing the relevant function here. Whether RHF or UHF, the Fock matrix (AO basis) is read from the reference wavefunction and transformed into the MO basis.
There’s no mechanism in your script to change the Fock matrix of the wavefunction you send in compared to the Fock matrix of your original Roothaan-Hartree-Fock solution, so the energy should only be correct if the AO basis Fock matrix is invariant to your change of orbitals. It isn’t.
The expression for the Fock matrix for not-necessarily-canonical orbitals is f^p_q = h^p_q + gbar^pi_qi
, where the double i
denotes summation over occupied orbitals. Different choices of occupied orbitals mean different two-electron integral terms, whatever your basis.
So, is your choice of occupied orbitals different from the SCF choice? Certainly. If your new orbitals preserved the occupied space, the Fock matrix should be block-diagonal in the occupied and virtual spaces, but I gather that it isn’t, if this is “pseudocanonicalized”. Therefore, your Fock matrix is incorrect in the AO basis.
As a demonstration, try the following
tricky_wfn = psi4.energy("hf/cc-pvdz", return_wfn=True)[1]
scf_wfn = psi4.energy('b3lyp/cc-pvdz', return_wfn=True)[1]
scf_wfn.Ca().copy(tricky_wfn.Ca())
scf_wfn.Cb().copy(tricky_wfn.Cb())
trick_ccsd_E = psi4.energy("ccsd/cc-pvdz", ref_wfn = scf_wfn)
Even though the orbitals you send in are perfect, the Fock matrix coming out of DFT isn’t the Fock matrix as defined by CC, so you don’t even get zero contribution to MP2 singles.
I leave it to @crawdad to suggest a fix. I know that the cc
mega-module has some orbital rotation machinery, so it may be possible to reconstruct the Fock matrix rather than trusting the wavefunction to be consistent.