evangelistalab / forte

http://www.evangelistalab.org
GNU Lesser General Public License v3.0
54 stars 30 forks source link

Remix of PR #420 #423

Closed fevangelista closed 1 month ago

fevangelista commented 1 month ago

This is a reworked version of PR #420. Here is the original description:

Description

This PR implements new orbital orthonormalization procedure for initial orbitals in Forte. The new procedure keeps the relative ordering of the guess orbitals, that is, no canonicalization procedure is performed. Now initial orbitals can come from a different basis set (see added test case df-casscf-4-basis).

To actually check these initial guess orbitals, some Python codes can be added in the input:

molecule Cr2{
Cr
Cr 1 3.5
symmetry c1
}

set globals{
scf_type             df
reference            rhf
guess                sadno
basis                cc-pwcvtz
df_basis_scf         def2-universal-jkfit
df_basis_mp2         def2-universal-jkfit
restricted_docc       [18]
active                [22]
}
escf, wfn = energy('scf', return_wfn=True)

from forte._forte import make_mo_space_info_from_map
mo_space_info = make_mo_space_info_from_map(wfn.nmopi(),
                                            Cr2.point_group().symbol(),
                                            {"ACTIVE": [22], "RESTRICTED_DOCC": [18]}, [])

# same basis, different geometry
wfn_old = wfn.from_file('wfn1.npy')  # cc-pwCVTZ @ 3.25 Angstrom
from forte.proc.orbital_helpers import ortho_orbs_new
wfn.Ca().copy(ortho_orbs_new(wfn.Ca(), wfn.S(), wfn_old.Ca(), mo_space_info))

# different basis, different geometry
wfn_old = wfn.from_file('wfn2.npy')  # cc-pVDZ @ 3.25 Angstrom
from forte.proc.orbital_helpers import basis_projection
Ca = basis_projection(wfn_old, wfn, mo_space_info)
wfn.Ca().copy(Ca)

# plot orbitals
orbs = [19 + i for i in range(22)]
set cubeprop_orbitals $orbs
cubeprop(wfn)

User Notes

Checklist

lcyyork commented 1 month ago

@fevangelista Thank you for clearing up the logic. I only have a few minor comments.

fevangelista commented 1 month ago

Sorry this took a while. The newer code I just added can handle symmetry and works fine. Note that I changed the test since projection and frozen MCSCF orbitals do not make sense together.

lcyyork commented 2 weeks ago

@fevangelista This new procedure fails if the target basis set is smaller than the one in the reference. For example, energy('forte', ref_wfn=wfn, basis='cc-pvdz') fails if wfn is obtained from cc-pVTZ. The dimension of Ca is incorrect.

I also have a concern of doing general basis set projection without separating occupied and unoccupied. Is the ordering of occupied orbitals preserved between the two sets of orbitals? This seems important for maintaining a good guess of active space for CASSCF.