The coupled cluster model provides a higher level of accuracy beyond the MP2 approach. The purpose of this project is to understand the fundamental aspects of the calculation of the CCSD wave function and corresponding energy. A good reference article (with correct, factored equations!) is J.F. Stanton, J. Gauss, J.D. Watts, and R.J. Bartlett, J. Chem. Phys. volume 94, pp. 4334-4345 (1991).
Each term in the equations given in the above paper by Stanton et al. depends on the T1 or T2 amplitudes as well as Fock-matrix elements and antisymmetrized, Dirac-notation two-electron integrals, all given in the molecular spin-orbital basis (as opposed to the spatial-orbital basis used in the earlier MP2 project). Thus, the transformation of the AO-basis integrals into the spatial-MO basis must also include their translation into the spin-orbital basis:
Thus, if you know the ordering of the orbitals (e.g. all occupied orbitals before virtual orbitals, perhaps alternating between alpha and beta spins), it is straightforward to carry out the integration over the spin components (the sigmas) in the above expression. Thus, each spatial-orbital MO-basis two-electron integral translates to 16 possible spin-orbital integrals, only four of which are non-zero.
Don't forget that you must also create the spin-orbital Fock matrix:
Suggestion: For simplicity, store the two-electron integrals in a four-dimensional array. This will greatly facilitate debugging of the complicated CCSD equations.
For Hartree-Fock reference determinants, the most common initial guess for the cluster amplitudes are the Moller-Plesset first-order perturbed wave function:
Note that if you have constructed the Fock matrix, two-electron integrals, and initial-guess amplitudes correctly at this point, you should be able to compute the MP2 correlation energy using the simple spin-orbital expression and get identical results to those from Project #4:
Use the spin-orbital Eqs. 3-13 from Stanton's paper to build the two-index (F) and four-index (W) intermediates, as well as the effective doubles (labelled with the Greek letter tau).
Use Eqs. 1 and 2 from Stanton's paper to compute the updated T1 and T2 cluster amplitudes.
Calculate the current CC correlation energy:
Compare energies and cluster amplitudes (using RMS differences) between iterations to check for convergence to some specified cutoff. If convergence is reached, you're done; if not, return to Step #3 and continue.