Get two-particle density matrix from a CC calculation

[plafiosca]

Hello everybody,

I am new to Psi4, and I have just installed the 1.3.2 version. I was wondering if it is possible to retrieve the two-particle density matrix (in AOs) from a CC calculation. From the manual I have understood that there is a get_tpdm() method for the CIWavefunction, but I did not manage to find something similar for the CCWavefunction. In the CCDENSITY module seems to be available a kind of twopdm() function that calculate the two-particle density matrix contributions, but I don’t know how to use it.

Thank you very much, I appreciate every kind of help!

Piero

[jmisiewicz]

Can I ask what exactly your use case is? I can offer better advice based on that.

The Coupled Cluster code is old and not exactly well-maintained, but I can be convinced to go exploring it for the sake of reduced density matrices.

[plafiosca]

Hi,

thank you for you help! I need the TPDM because I want to compare the results with the one calculated by means of the HF method. I managed to recover the OPDM in different ways, such with the following input that I adapted from this forum:

---

import psi4
psi4.set_output_file(“2rdm.out”, False)
psi4.set_memory(‘20GB’)
psi4.set_num_threads(28)
mol_string = “”"
symmetry c1
0 1
O 0.009319 1.133156 0.000000
H 0.023452 0.185621 0.000000
H 0.906315 1.422088 0.000000
units angstrom
“”"
psimol = psi4.geometry(mol_string, name=“h2”)
psi4.set_options({
‘SCF_TYPE’: ‘PK’,
‘BASIS’: “3-21G”,
‘E_CONVERGENCE’: 1e-12,
‘D_CONVERGENCE’: 1e-12,
‘R_CONVERGENCE’: 1e-12,
‘REFERENCE’: ‘RHF’
})
E, wf = energy(‘scf’, return_wfn=True)
grad, ccwf = psi4.gradient(‘ccsd’, return_wfn=True, ref_wfn=wf)
print(ccwf.Da().to_array(copy=False, dense=True))
psi4.core.clean_options()
psi4.core.clean()

but if I don’t know how to extend it to the TPDM case.

Thank you very much.

Piero

[jmisiewicz]

After looking over the code, I cannot find an easy way to do this. This gets stored in some files, but I don’t know an easy way to convert those into the numpy arrays you’re looking for.

While I could figure it out with enough time, I don’t have the time right now. Some work in the code may be needed. If you need this badly enough that you’d be willing to work with the C++ code, I can see if I can get you in touch with somebody who’d be better able to help.

[crawdad]

I can give you a python code that will compute the CC twopdm. Would that help?

[plafiosca]

Hello everybody,

thank you very much for the kind answers!

@jmisiewicz I think I have enough time to work with the C++ code, even if I am not very used to this language. Surely some help will be useful, but still it would be really demanding.

However the python code of @crawdad seems to be the optimal choice, since I only need to read the twopdm without any further manipulation. If you could give me access to the script it would be amazing. Can I ask you with which level of CC the script is able to generate the twopdm? It would be very interesting if I manage to calculate the twopdm with all the level of CC implemented in Psi4 (such as CC2, CCSD, CC3 and so on), or if I can extend the code in order to do that.

Best,
Piero

[crawdad]

It can give you the twopdm for RHF-CCSD only, though other levels could be added easily. Psi4 doesn’t actually have CC3 twopdms.

Let me clean some things up, and then I’ll send you the link to the code.

[plafiosca]

Thank you very much, I will try to include other levels too.

[crawdad]

The code is available here: GitHub - lothian/pycc: PyCC is a simple, Python-based, reference implementation of the coupled cluster method of ab initio quantum chemistry.

I’m adding a lot to it right now, but the density code is in place with a test case and some documentation.

[plafiosca]

Thank you very much, I will go through the code as soon as possible.