Skip to content

CABS, rework build_cabs_space #2982

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 20 commits into from
Sep 11, 2023
Merged

Conversation

EricaCMitchell
Copy link
Contributor

Description

Provides users the ability to form a basisset composed of two combined basissets, the two can be combined simply or through use of the complementary auxiliary basis set (CABS) approach of Valeev et al.

Dev notes & details

  • [X ] Uncomments pyconstruct_combined in qcdb with minor changes to fit updated python
  • [X ] Changes input to build_ri_space in OrbitalSpace to require a prebuilt combined BasisSet object
  • [X ] Changes build_cabs_space in OrbitalSpace to use a full SVD
  • [X ] Adds pytest test_orbitalspace.py to check orthonormality between orbital basis set and CABS
  • [X ] Adds updated F12 basissets from the Basis Set Exchange
  • [X ] Use of the CABS created works with my F12 plugin giving correct energies

Questions

  • Shoud I generalize _pybuild_basis aka BasisSet.build to work with a list of keys, targets, roles, and other? Or create a new function to more easily access the building of a combined basisset?

Checklist

Status

  • [X ] Ready for review
  • Ready for merge

@susilehtola
Copy link
Member

The fix to one-electron integrals should be its own pull request. It has nothing to do with the other changes, which are pretty extensive by themselves to review.

@EricaCMitchell EricaCMitchell force-pushed the CABS branch 3 times, most recently from 484314c to be7ab1b Compare June 14, 2023 13:48
Copy link
Member

@loriab loriab left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for all the basis set work!

key=lambda shell: shell.l)

# Molecule and parser prepped, call the constructor
mol.set_basis_all_atoms(name, "CABS")
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This isn't allowing basis { ... } style input, but that's fine for now.

@loriab loriab added this to the Psi4 1.9 milestone Jul 1, 2023
@loriab loriab added the feature Extends an existing Psi feature or develops a new one. label Jul 1, 2023
@JonathonMisiewicz
Copy link
Contributor

@konpat, do you think this PR would benefit from you looking it over? I know you have some F12 work in mind.

@konpat
Copy link
Contributor

konpat commented Jul 11, 2023

Yup, @JonathonMisiewicz, I do have this PR on my radar and it's overall a super useful contribution! Unfortunately, I'm very short on time for the next couple weeks and it might still be a while until I properly test it out.

Copy link
Contributor

@konpat konpat left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great contribution! It took some work, but I made my development F12 code to work with it, and the results agree with my old CABS code (which I'm happy to ditch).

Sorry it took so long.

S_mo = np.linalg.multi_dot([C1.T, S_ao, C2])

if o1 != o2:
np.testing.assert_allclose(np.dot(S_mo.T, S_mo), np.zeros((C2.shape[1], C2.shape[1])), rtol=1e-05, atol=1e-07)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this test really need to verify that STS is identity or zero? I think we just need to check the matrix S itself for orthonormality (or orthogonality between different spaces):

    if o1 != o2:
        np.testing.assert_allclose(S_mo, np.zeros((C1.shape[1], C2.shape[1])), rtol=1e-05, atol=1e-07)
    else:
        np.testing.assert_allclose(S_mo, np.eye(C1.shape[1]), rtol=1e-05, atol=1e-07)

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Another (very picky) issue: shouldn't we prefer to check if o1 and o2 are approximately equal, to within some numerical threshold?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@konpat You're right! I've fixed the test to just look at the S_mo matrix.

@JonathonMisiewicz The two spaces are going to be unequal just by construction with the CABS space consisting of the basis functions of the OBS. What I'm checking for here is that we are fulfilling the condition that S12 C = 0, as described in the paper by Valeev. I guess I do check more cases than is necessary in the test though.

Copy link
Member

@jturney jturney left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I've been meaning to approve this for a while. Looks good!

Copy link
Contributor

@JonathonMisiewicz JonathonMisiewicz left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

A few nitpicky additional comments.

S_ao = np.array(mints.ao_overlap(bs1, bs2))
S_mo = np.linalg.multi_dot([C1.T, S_ao, C2])

if o1 != o2:
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What if o1 and o2 are ever so slightly off? I understand that in practice, we won't encounter that, but please either don't check this case at all, or have a more robuster check to see if those two matrices are similar.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The function spaces returns the orbital basis set (OBS) and the complementary auxiliary basis set (CABS). By definition, the CABS is the orthogonal complement to the OBS so the overlap of the MOs between the two should be rigorously zero.

Copy link
Contributor

@JonathonMisiewicz JonathonMisiewicz Sep 11, 2023

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We're talking past each other, so let's back up. Is the condition o1 != o2 ever going to be false? The code tells me both true and false are possible. What you are telling me tells me that o1 != o2 should always be true.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@pytest.mark.parametrize("o1,o2", [(0, 0), (0, 1), (1, 0), (1, 1)])

The test looks at the orthonormality of the OBS with itself (o1 = 0 and o2 =0), the OBS with the CABS (o1 = 0 and o2 = 1, or vice versa), and the CABS with itself (o1 = 1 and o2 =1). So o1 != o2 is sometimes false.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Okay, that makes much more sense. I'm used to different parameter values being entirely different computations, not different checks on a single computation. Can you please add a docstring explaining what this test does, so that doesn't bite anybody else? Example: "Test that the CABS and OBS spaces are orthonormal." Renaming "o1" and "o2" to "idx1" and "idx2" would also make a lot of difference.

@JonathonMisiewicz
Copy link
Contributor

I'll happily approve after these final comments. Unless anybody else objects, I can merge this in last. Thanks again for your patience!

@JonathonMisiewicz
Copy link
Contributor

Okay, that makes much more sense. I'm used to different parameter values being entirely different computations, not different checks on a single computation. Can you please add a docstring explaining what this test does, so that doesn't bite anybody else? Example: "Test that the CABS and OBS spaces are orthonormal." Renaming "o1" and "o2" to "idx1" and "idx2" would also make a lot of difference.

Just want to make sure this one didn't get lost. 😅

Copy link
Contributor

@JonathonMisiewicz JonathonMisiewicz left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Devs have had ample time to review this one if interested, so I'll add this to merge queue.

Thanks for the code!

@JonathonMisiewicz
Copy link
Contributor

The reason for docs build failure is unrelated to this PR. I've asked Lori to take a look.

@loriab
Copy link
Member

loriab commented Sep 11, 2023

Almost certainly because Ben took down some old qcel/qcng docs in preparation for the qcfractal release. I'll hunt down the replacement. Probably this can still join the merge queue b/c only Azure is required to pass.

@JonathonMisiewicz JonathonMisiewicz added this pull request to the merge queue Sep 11, 2023
Merged via the queue into psi4:master with commit 61e1550 Sep 11, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
feature Extends an existing Psi feature or develops a new one.
Projects
None yet
Development

Successfully merging this pull request may close these issues.

6 participants