Printing density matrix, error in tce dipole

[Stephen H. Yuwono]

Dear NWChem devs,

I have a couple of questions related to the TCE module. My first question is: is there a way to print the 1-body reduced density matrix (or transition density matrix) from a coupled-cluster calculation?

My second question is somewhat related. I was trying to run a property calculation at the CCSD and CCSDT level (using the "dipole" keyword inside tce), but I was faced with an error that says (4 copies of this error message, since I was using 4 processors):

3:3:ngai_put_common:cannot locate region: ./h2_ccsd.d1x [1:0 ,1:1 ]:: -999

(rank:3 hostname:dev-intel18 pid:10079):ARMCI DASSERT fail. ../../ga-5.7.1/armci/src/common/armci.c:ARMCI_Error():208 cond:0

In contrast, running an energy calculation (for both ground and excited states) works just fine.

I would appreciate your help with this matter.

Thank you,

Stephen

[nicholas.bauman]

Dear Stephen,

Would you be so kind as to attach your input for comparison. We need to figure out what is causing that message.

As for printing out the density, When you run a calculation with 'DIPOLE' keyword, it will print out the largest elements of the P-P, P-H, H-P, and H-H blocks. To print the entire density matrix simply add the line 

set tce:printtol 0.0

after the TCE block in your input.

Best,

Nick

[Stephen H. Yuwono]

Dear Nick,

I believe I have replied to this message last week, but for some reason it seems that it did not go through (probably some error with posting). Thus, please allow me to re-reply your message (and apologies if it is a repetition!).

First, I am sorry for not including the input file from the get go - it is a habit from using a different ticket submission that does not allow file attachment.

In any case, please find attached four sets of input and output files run at EOMCCSDT/6-31G level of theory:

• he2_ccsdt_6-31g.* > energy only, successful run
• he2_ccsdt_6-31g_dipole.* > energy and dipole, failed
• hf_ccsdt_6-31g.* > energy only, successful run
• hf_ccsdt_6-31g+dipole.* > energy and dipole, successful run

Thank you very much for the pointer for printing 1-RDM through TCE. Unfortunately, I was not able to find the parts of the output related to the 1-RDM. I do see that now it prints out all T, R, and L amplitudes (up to double excitation) with the threshold set to 0.0. Could you please help point me to the right direction? I am also attaching the input and output files for this run (hf_ccsdt_6-31g_rdm.*).

All of the jobs were run using NWChem-7.0.0, compiled with Intel Parallel Studio XE 2019.

Thank you very much for your help!

Best,

Stephen

[nicholas.bauman]

Dear Stephen,

The 1-PDM printing was only implemented for CCSD. So if you change the method to CCSD you will see this printing. It has separate P-H and H-P sections. Keep in mind that some other codes may use symmetrized matrix elements in which these blocks are the same. 

As for CCSDT, there is no printing and these routines should be taken with some careful consideration. The routines for the CCSDT density matrix were implemented a while back and tested on only a few molecules. A different approach compared to CCSD was implemented that symmetrizes the density matrices. While the NWChem package does contain routines for separate P-H and H-P calculations, they are not connected to the DIPOLE calculation. I can describe how to implement printing for this case, but without an extensive investigation into all of the equations and testing, I cannot give my guarantee that everything in the pipeline is correct to that point. 

Let me know how you want to proceed.

The failed calculation is a separate problem. I think it has to do with such a small basis set. If you use a larger one, it completes fine. I get an error in communication related to the dipole file. 

Best,

Nick

[Stephen H. Yuwono]

Dear Nick,

Thank you very much for the detailed response. I have tried printing the 1-PDM using CCSD and it works as expected. I do like the fact that it is not symmetrized, since that is exactly what I need to compare with my own implementation of 1-PDM.

I would really appreciate it if you could point me to where the CCSDT density matrix is coded and how I could print the CCSDT 1-PDM. As I said above, the routines for separate P-H and H-P calculations will be ideal, but the symmetrization process would be helpful as well. I think I have the working equations for the CCSDT 1-PDM in spin integrated forms, so I can double check what NWChem produces.

I have also checked the failed calculations and I agree that it might have to do with such a small basis set (for example, in the case of HF molecule the job is running just fine). I used it to make sure that my own implementation of 1-PDM is correct, but I guess I will have to use something else for it.

Thank you again!

Best,

Stephen

[nicholas.bauman]

Dear Stephen,

The separate PP, PH, HP, and HH routines for the ccsdt density matrix are in the src/tce/gradients directory.

There is a file called tce_energy.F in src/tce directory. You can search for a line with "! CCSDT   left" . This is where the lambda is computed. Either after this loop or the part that that is evaluated after convergence "if (residual .lt. thresh)" you can call the PP, PH, HP, and HH routines for CCSDT. Just follow the prescription in "ccsd_lambda.F". The only thing you will be changing is the call to the separate P/H block routines. Just make sure the variables match up with the call and routines. That prescription also contains the call to print out the blocks. After that, recompile and run.

Some suggestions in checking your density matrices:

1. While the reduced density matrix is not generally Hermitian due to the biorthogonality nature of the wave function, it is symmetric for cases in which the theory is exact, hence 3 or less electron systems or systems with only 3 or less valence orbitals for CCSDT. A symmetric density matrix in this case does not imply that the equations are implemented correctly, but can help identify problems in your equations.

2. Since you are dealing with spin-integrated equations, a check that is always good to perform is to switch the number of alpha and beta electrons for an open-shell system. You should end up with the same results. The only thing to be careful with is to make sure that if routines for different spin combinations are called  and updated separately, that their order is reversed to reflect to switch in the number of alpha and beta electrons.

All the best,

Nick