DOCI 3RDM implementation

[DanielCalero1]

In this PR I went address the issue #58 . I implemented the 3RDM computation for DOCI wave functions. Here I list the changes done:

1. I changed the source code in rdm.cpp to return two more 3d-matrices d3 and d4 in the case of DOCI wave functions. d3 stores the elements of the type p+ q+ r+ p q r. d4 stores the elements p+ q+ q+ p r r. Both matrices are pair preserving terms.
2. I changed the binding to python on rmd.cpp and pyci.h to recognize the d3, d4 matrices
3. I updated the function spinze_rdms in utility.py to fill the fill the 3rdm blocks in the DOCI case.
4. I added a new function called "spin_free_rdms". This is basically a wrapper of the function spinze_rdms, that returns the spin-traced rdms. It can be useful for people working with spin-free systems.

I tested all the functions with BH, H10 and H6, they all behaved correctly. The tested I did were the following:

• First, I tested the normalization of the 3rdm and its blocks (AAAAAA,AABAAB,BBABBA,BBBBBB) and it was correct.
• Then, I traced over the blocks of the 3rdm to obtain the 2rdm blocks. First tracing AAAAAA to obtain AAAA (the same for all beta), and then tracing over the mixed spin AABAAB to obtain both AAAA or ABAB depending on the indices that I was tracing over. At the end, the result was correct with the appropriate normalization factor (This is the same test they do in pyscf)
• Finally, I tested the the 3rdm blocks with the correspondent blocks of pyscf and the result was correct
• For the spin_free_rdms function, I tested that the result was the same to the 3rdm spin traced function in pyscf

I did not include any function in the module test because I am not sure where to put it, however if you give me some guidelines on how you want it implemented I can try to do it.
If not, I can pass the jupyter notebook where I was running the test such that anyone can test the implementation easily

[msricher]

This looks good.

Can you make your compute_rdms function separate from the old one? It's good to have one that only compute 1,2-RDMs, and then one that also computes the 3,4-RDMs.

You also need to add some tests. And make sure that you use very small test cases with your 3,4-RDM tests! Currently, the existing tests fail because it tries to allocate 400GB of memory to compute the 4-RDM of the test system. I would recommend adding a test FCIDUMP file for H2 or He in a minimal basis (sto-3g).

[DanielCalero1]

Sure, I can do the changes you mention. However let me ask you something about it

First, how should I call the new function compute_rdms, because since it will have the same inputs as the old one (the DOCI wave function and the coefficients) I think I can not just overload the function (am I wrong?), so I will have to create a new one.
Second, just to be sure, the function you want me to create separate from the old one is compute_rdms which is on the C++ source code right, not the spinize_rdms function in python.
Third, I can surely can add some tests, however I think H2 or He might not work, since those are two electron systems, the 3RDM will be full of zeros, I can use Be in the minimal basis sto-3g, or maybe H6 in sto-3g as well, is that ok?

[msricher]

Sure, I can do the changes you mention. However let me ask you something about it

First, how should I call the new function compute_rdms, because since it will have the same inputs as the old one (the DOCI wave function and the coefficients) I think I can not just overload the function (am I wrong?), so I will have to create a new one. Second, just to be sure, the function you want me to create separate from the old one is compute_rdms which is on the C++ source code right, not the spinize_rdms function in python. Third, I can surely can add some tests, however I think H2 or He might not work, since those are two electron systems, the 3RDM will be full of zeros, I can use Be in the minimal basis sto-3g, or maybe H6 in sto-3g as well, is that ok?

Yeah, make new C++ functions. Maybe compute_rdms_34?

Yep, you can maybe use Be and B in minimal basis set as test cases. My bad, I forgot it was 3/4 RDMS.

[DanielCalero1]

Hi, I already created the new function compute_rdms_34 and made all the bindings with python. However, the problems with the test persisted. Taking a look into the problems, its seems that the problems are coming from the function spinze_rdms in utility.py. In the line that I try to allocate memory for the 3rdm it crashes since is too much memory.
I tried to move that line inside the if statement such that only allocates that memory in the case of DOCI wfn, and now some of the test are working correctly, except the one that are for DOCI wfns.
I think that the simplest solution would be doing the same as in C++, just creating a new function spinze_rdms_34 in utility.py, and create new tests for that function with smaller systems. In that way, all the previous tests will still work with spinize_rdms.

[DanielCalero1]

I have implemented what I said in the last comment. I create a new function spinize_rdms_34 that only works for DOCI wave functions. I leave the function spinize_rdms as it was.

I also added a new test function. In this test I check the 3RDM trace, the traces of its blocks, and it also checks that tracing over two indices of the 3RDM gives the 2RDM (for all the blocks).

As far as I tested, all the tests were working perfectly.

SideNote: Just to be clear, this PR is only for the implementation of the 3RDM, I haven't implemented the 4RDM yet, although it would be pretty straightforward with this code, as the 4 RDM only has two new sets of elements, the rest would be obtained re using the 3rdm and 2rdm elements. The tricky part on the 4RDM (as it was for the 3RDM) is the symmetrization, the C++ code for the implementation I already have it (I believe). The 4RDM has a bunch of symmetries, and one has to be careful when adding all those symmetries.
I'll create another Issue describing the implementation. I hope I can work on this ASAP, however what I need right now for my project is the 3RDM, so, I'll leave it here for a time and then come back to implement the 4RDM

[msricher]

This looks OK to me!

I just have a few minor requests.

• Can you make the code style match a bit better? It's mostly missing spaces after commas and operators ('x, y' instead of 'x,y', 'a + b' instead of 'a+b', etc.). Just try to keep consistent with what already exists.

• I thought more about names, and can you rename the '_34' methods to '_1234'? Since it does compute all of them, it's more accurate.

Otherwise, the tests pass (under the Checks tab on this page), and the code is good. You can fix it up a bit and then squash merge with a descriptive commit msg it if the tests pass again.

[DanielCalero1]

Sure! I'll implement those changes right away