In this project we performed some simple Coupled-Cluster Doubles (CCD) calculations for the electron gas. We used the basis and interactions described in Gustav Baardsens PhD-thesis, while the main part of the code was written by Ole Nordli as a part of his master thesis.
The first objective of the project was to reproduce some results from the literature, while the overall goal was to gradually extend the scope to the heterogenous neutron gas with interactions.
We wrote a class in C++ which sets up the basis and Fock-matrix for the electron gas. This class was successfully integrated into the prewritten code from Ole, and produced some nice results both for the Hartree-Fock energy and the CCD-calculations very quickly. (The results was in the same order of magnitude as those found in the literature) As we still lacked some of the insight into the theoretical foundation, we was not able to initialize the basis with the exact same parameters as found in the literature.
Pending work: * To alter the parameter $L$ ($L^3$ is the total volume of the box) so that it corresponds to the given $r_s$ from the literature.
To look into how to implement periodic boundary conditions.*
I have written a coupled cluster extension (class) for the restricted Hartree-Fock-Roothaan solver that I and Gøran created for the project in FYS4411. Despite extensive debugging, benchmarking and comparisons, I am unable to reproduce results from the literature.
I have beeen trying to reproduce the results from a RHF-CCSD calculation by NWChem, as described in the following document:
I am supposed to find
CCSD iterations
Iter | Residuum | Correlation | Cpu | Wall |
---|---|---|---|---|
1 | 0.0891232379551 | -0.0358672469179 | 0.4 | 0.4 |
2 | 0.0317596201320 | -0.0454068882657 | 0.4 | 0.4 |
3 | 0.0126828916023 | -0.0483870059027 | 0.4 | 0.4 |
4 | 0.0053832778844 | -0.0494370597647 | 0.4 | 0.4 |
5 | 0.0023954452285 | -0.0498391184890 | 0.4 | 0.4 |
6 | 0.0011108272683 | -0.0500021724029 | 0.4 | 0.4 |
7 | 0.0005330437725 | -0.0500711904756 | 0.4 | 0.4 |
8 | 0.0002625570400 | -0.0501014381364 | 0.4 | 0.4 |
9 | 0.0001317707641 | -0.0501150974135 | 0.4 | 0.4 |
10 | 0.0000669953658 | -0.0501214303300 | 0.4 | 0.4 |
11 | 0.0000343654836 | -0.0501244348663 | 0.4 | 0.4 |
12 | 0.0000177357555 | -0.0501258887096 | 0.4 | 0.4 |
13 | 0.0000091923943 | -0.0501266039080 | 0.4 | 0.4 |
14 | 0.0000047789221 | -0.0501269605251 | 0.4 | 0.4 |
15 | 0.0000024899940 | -0.0501271402835 | 0.4 | 0.4 |
16 | 0.0000012995442 | -0.0501272316751 | 0.4 | 0.4 |
17 | 0.0000006791083 | -0.0501272784536 | 0.4 | 0.4 |
18 | 0.0000003552400 | -0.0501273025229 | 0.4 | 0.4 |
19 | 0.0000001859739 | -0.0501273149581 | 0.4 | 0.4 |
20 | 0.0000000974237 | -0.0501273214031 | 0.4 | 0.4 |
My own RHS-CCD calculation (the result is comparable for CCSD) finds instead a value of -0.0875531 a.u., overestimating the energy by a factor of approximately 1.7. The table below shows the iterations and lists the separate contributions to the energy.
Iteration | CCDEnergy | 0th Cont. | Lin.Cont. | Quad.cont. |
---|---|---|---|---|
1 | -0.0358976 | 0 | 0 | 0 |
2 | -0.0562349 | -0.0358976 | -0.0203373 | 0 |
3 | -0.0682 | -0.0358976 | -0.0325562 | 0.000253831 |
4 | -0.0753594 | -0.0358976 | -0.040085 | 0.000623145 |
5 | -0.0797398 | -0.0358976 | -0.044754 | 0.000911875 |
6 | -0.082475 | -0.0358976 | -0.0476843 | 0.00110689 |
7 | -0.0842106 | -0.0358976 | -0.0495466 | 0.00123363 |
8 | -0.0853254 | -0.0358976 | -0.0507435 | 0.00131574 |
9 | -0.0860489 | -0.0358976 | -0.0515206 | 0.00136931 |
10 | -0.0865231 | -0.0358976 | -0.0520301 | 0.00140458 |
11 | -0.0868371 | -0.0358976 | -0.0523675 | 0.00142802 |
12 | -0.0870474 | -0.0358976 | -0.0525935 | 0.00144376 |
13 | -0.0871899 | -0.0358976 | -0.0527467 | 0.00145447 |
14 | -0.0872878 | -0.0358976 | -0.052852 | 0.00146186 |
15 | -0.087356 | -0.0358976 | -0.0529254 | 0.00146703 |
16 | -0.0874043 | -0.0358976 | -0.0529774 | 0.00147072 |
17 | -0.0874391 | -0.0358976 | -0.0530149 | 0.00147339 |
18 | -0.0874645 | -0.0358976 | -0.0530423 | 0.00147536 |
19 | -0.0874834 | -0.0358976 | -0.0530627 | 0.00147684 |
20 | -0.0874977 | -0.0358976 | -0.053078 | 0.00147797 |
21 | -0.0875086 | -0.0358976 | -0.0530898 | 0.00147885 |
22 | -0.0875171 | -0.0358976 | -0.053099 | 0.00147954 |
23 | -0.0875237 | -0.0358976 | -0.0531062 | 0.00148009 |
24 | -0.087529 | -0.0358976 | -0.0531119 | 0.00148053 |
25 | -0.0875333 | -0.0358976 | -0.0531165 | 0.00148089 |
26 | -0.0875367 | -0.0358976 | -0.0531202 | 0.00148118 |
27 | -0.0875395 | -0.0358976 | -0.0531233 | 0.00148142 |
28 | -0.0875418 | -0.0358976 | -0.0531258 | 0.00148162 |
29 | -0.0875436 | -0.0358976 | -0.0531278 | 0.00148178 |
30 | -0.0875452 | -0.0358976 | -0.0531295 | 0.00148192 |
31 | -0.0875465 | -0.0358976 | -0.0531309 | 0.00148204 |
32 | -0.0875475 | -0.0358976 | -0.0531321 | 0.00148213 |
33 | -0.0875484 | -0.0358976 | -0.053133 | 0.00148221 |
34 | -0.0875492 | -0.0358976 | -0.0531338 | 0.00148228 |
35 | -0.0875498 | -0.0358976 | -0.0531345 | 0.00148234 |
36 | -0.0875503 | -0.0358976 | -0.0531351 | 0.00148239 |
37 | -0.0875508 | -0.0358976 | -0.0531356 | 0.00148243 |
38 | -0.0875511 | -0.0358976 | -0.053136 | 0.00148246 |
39 | -0.0875515 | -0.0358976 | -0.0531363 | 0.00148249 |
40 | -0.0875517 | -0.0358976 | -0.0531366 | 0.00148251 |
41 | -0.0875519 | -0.0358976 | -0.0531368 | 0.00148253 |
42 | -0.0875521 | -0.0358976 | -0.0531371 | 0.00148255 |
43 | -0.0875523 | -0.0358976 | -0.0531372 | 0.00148256 |
44 | -0.0875524 | -0.0358976 | -0.0531374 | 0.00148258 |
45 | -0.0875525 | -0.0358976 | -0.0531375 | 0.00148259 |
46 | -0.0875526 | -0.0358976 | -0.0531376 | 0.0014826 |
47 | -0.0875527 | -0.0358976 | -0.0531377 | 0.0014826 |
48 | -0.0875528 | -0.0358976 | -0.0531378 | 0.00148261 |
49 | -0.0875528 | -0.0358976 | -0.0531378 | 0.00148262 |
50 | -0.0875529 | -0.0358976 | -0.0531379 | 0.00148262 |
51 | -0.0875529 | -0.0358976 | -0.0531379 | 0.00148262 |
52 | -0.087553 | -0.0358976 | -0.053138 | 0.00148263 |
53 | -0.087553 | -0.0358976 | -0.053138 | 0.00148263 |
54 | -0.087553 | -0.0358976 | -0.053138 | 0.00148263 |
55 | -0.087553 | -0.0358976 | -0.053138 | 0.00148263 |
56 | -0.0875531 | -0.0358976 | -0.0531381 | 0.00148264 |
57 | -0.0875531 | -0.0358976 | -0.0531381 | 0.00148264 |
58 | -0.0875531 | -0.0358976 | -0.0531381 | 0.00148264 |
59 | -0.0875531 | -0.0358976 | -0.0531381 | 0.00148264 |
As already noted, I am unable to find any errors in the implementation. There is, however, a number of aspects in the method that remains unclear to me. Some of these are:
(1) Shavitt-Bartlett derives the CCD equations, ending in the T2 amplitude equation (9.126) on page 288. In this equation they restrict the summations to $$ i>j, a>b $$ It is unclear to me if this is due to the symmetries in the matrix elementes $\langle ab || ij \rangle$ so that one may manually set the others for each iteration, or if it means that the T2 amplitude tensor will contain a lot of zeros. In my code I have interpreted this in the former manner.
(2) I have manually ensured that my coefficient matrix obtained from the RHF-procedure contains the same elements at the same locations as Ole Tobias Norlis elements. This means that the CC-iterations have a very similar initial condition as the ones described in the literature, as may be seen from the initial energy in the table above. However, upon initialization, Ole sets up his elements like
$\langle pr||qs\rangle = \langle pq|rs\rangle - \langle ps|rq\rangle$
(permuting the r and q index). I adopted this in my code to reproduce his basis, but I did not compensate for this in the CC iterations since I found no reason to do so.
The error in the code needs to be identified and fixed. As the initial energy is very comparable to the ones found in the literature, the error is most likely related to the iterative process itself.
I have made the output from the code more specific as to how much each term in the CCD T2 amplitudes contributes. Following Shavitt-Bartlett, this means breaking the full correlation energy down to the following parts:
My aim for today is to compare each constituent contribution to similar results from Ole Tobias Norlis. I have prepared calculations for both H2O and H2 using a STO-3g basis through a RHF procedure prior to the CCD calculation.
With help from Ole I finally identified the error in the code. It was actually a bug occuring in the for-loops for the linear L2b contribution.
I now have a working CCD and possibly a CCSD (needs more testing) implementation for my RHF solver.
With the new progress it is possible to perform benchmark calculations on certain molecules and atoms. I therefore plan to perform the following experiments today (if time permits it):
I observed some minor discrepancies in the H2O CCSD calculation performed yesterday, so I will probably spend some time reading through the singles contribution code today. I suspect there is a very minor indexing error in one contribution.
I performed some of the tests mentioned above. The results can be found in the notebook "Benchmarking.ipynb". The results was satisfactory compared to other similar software.
I also inspected all the contributions to the CCSD code and found no errors. I now feel a bit more confident that I have a working CCD and CCSD code. After a meeting with Morten, we agreed on some steps ahead:
In [ ]: