3
\brief Enter brief description of file here
7
#include <libdpd/dpd.h>
13
namespace psi { namespace ccenergy {
17
int h, ij, ab, i, j, a, b, I, J, A, B;
18
int Isym, Jsym, Asym, Bsym;
20
dpdbuf4 tauIJAB, tauijab, tauIjAb, tauiJaB, tauIjbA;
21
dpdbuf4 tIJAB, tijab, tIjAb;
24
nirreps = moinfo.nirreps;
26
if(params.ref == 0) { /** RHF **/
28
dpd_buf4_init(&tIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tIjAb");
29
dpd_buf4_copy(&tIjAb, CC_TAMPS, "tauIjAb");
30
dpd_buf4_close(&tIjAb);
32
dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA");
33
dpd_file2_mat_init(&tIA);
34
dpd_file2_mat_rd(&tIA);
36
dpd_buf4_init(&tauIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tauIjAb");
38
for(h=0; h < nirreps; h++) {
40
dpd_buf4_mat_irrep_init(&tauIjAb, h);
41
dpd_buf4_mat_irrep_rd(&tauIjAb, h);
43
for(ij=0; ij < tauIjAb.params->rowtot[h]; ij++) {
44
i = tauIjAb.params->roworb[h][ij][0];
45
j = tauIjAb.params->roworb[h][ij][1];
46
I = tIA.params->rowidx[i];
47
J = tIA.params->rowidx[j];
48
Isym = tIA.params->psym[i];
49
Jsym = tIA.params->psym[j];
50
for(ab=0; ab < tauIjAb.params->coltot[h]; ab++) {
51
a = tauIjAb.params->colorb[h][ab][0];
52
b = tauIjAb.params->colorb[h][ab][1];
53
A = tIA.params->colidx[a];
54
B = tIA.params->colidx[b];
55
Asym = tIA.params->qsym[a];
56
Bsym = tIA.params->qsym[b];
58
if((Isym==Asym) && (Jsym==Bsym))
59
tauIjAb.matrix[h][ij][ab] +=
60
(tIA.matrix[Isym][I][A] * tIA.matrix[Jsym][J][B]);
65
dpd_buf4_mat_irrep_wrt(&tauIjAb, h);
66
dpd_buf4_mat_irrep_close(&tauIjAb, h);
69
/* This will generate the tauIjbA file from tauIjAb */
70
dpd_buf4_sort(&tauIjAb, CC_TAMPS, pqsr, 0, 5, "tauIjbA");
71
dpd_buf4_close(&tauIjAb);
73
dpd_file2_mat_close(&tIA);
74
dpd_file2_close(&tIA);
76
else if(params.ref == 1) { /** ROHF **/
78
dpd_buf4_init(&tIJAB, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tIJAB");
79
dpd_buf4_copy(&tIJAB, CC_TAMPS, "tauIJAB");
80
dpd_buf4_close(&tIJAB);
82
dpd_buf4_init(&tijab, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tijab");
83
dpd_buf4_copy(&tijab, CC_TAMPS, "tauijab");
84
dpd_buf4_close(&tijab);
86
dpd_buf4_init(&tIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tIjAb");
87
dpd_buf4_copy(&tIjAb, CC_TAMPS, "tauIjAb");
88
dpd_buf4_close(&tIjAb);
90
dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA");
91
dpd_file2_mat_init(&tIA);
92
dpd_file2_mat_rd(&tIA);
93
dpd_file2_init(&tia, CC_OEI, 0, 0, 1, "tia");
94
dpd_file2_mat_init(&tia);
95
dpd_file2_mat_rd(&tia);
97
dpd_buf4_init(&tauIJAB, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tauIJAB");
99
for(h=0; h < nirreps; h++) {
101
dpd_buf4_mat_irrep_init(&tauIJAB, h);
102
dpd_buf4_mat_irrep_rd(&tauIJAB, h);
104
for(ij=0; ij < tauIJAB.params->rowtot[h]; ij++) {
105
i = tauIJAB.params->roworb[h][ij][0];
106
j = tauIJAB.params->roworb[h][ij][1];
107
I = tIA.params->rowidx[i];
108
J = tIA.params->rowidx[j];
109
Isym = tIA.params->psym[i];
110
Jsym = tIA.params->psym[j];
111
for(ab=0; ab < tauIJAB.params->coltot[h]; ab++) {
112
a = tauIJAB.params->colorb[h][ab][0];
113
b = tauIJAB.params->colorb[h][ab][1];
114
A = tIA.params->colidx[a];
115
B = tIA.params->colidx[b];
116
Asym = tIA.params->qsym[a];
117
Bsym = tIA.params->qsym[b];
119
if((Isym==Asym) && (Jsym==Bsym))
120
tauIJAB.matrix[h][ij][ab] +=
121
(tIA.matrix[Isym][I][A] * tIA.matrix[Jsym][J][B]);
122
if((Isym==Bsym) && (Jsym==Asym))
123
tauIJAB.matrix[h][ij][ab] -=
124
(tIA.matrix[Isym][I][B] * tIA.matrix[Jsym][J][A]);
129
dpd_buf4_mat_irrep_wrt(&tauIJAB, h);
130
dpd_buf4_mat_irrep_close(&tauIJAB, h);
133
dpd_buf4_close(&tauIJAB);
135
dpd_buf4_init(&tauijab, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tauijab");
137
for(h=0; h < nirreps; h++) {
139
dpd_buf4_mat_irrep_init(&tauijab, h);
140
dpd_buf4_mat_irrep_rd(&tauijab, h);
142
for(ij=0; ij < tauijab.params->rowtot[h]; ij++) {
143
i = tauijab.params->roworb[h][ij][0];
144
j = tauijab.params->roworb[h][ij][1];
145
I = tia.params->rowidx[i];
146
J = tia.params->rowidx[j];
147
Isym = tia.params->psym[i];
148
Jsym = tia.params->psym[j];
149
for(ab=0; ab < tauijab.params->coltot[h]; ab++) {
150
a = tauijab.params->colorb[h][ab][0];
151
b = tauijab.params->colorb[h][ab][1];
152
A = tia.params->colidx[a];
153
B = tia.params->colidx[b];
154
Asym = tia.params->qsym[a];
155
Bsym = tia.params->qsym[b];
157
if((Isym==Asym) && (Jsym==Bsym))
158
tauijab.matrix[h][ij][ab] +=
159
(tia.matrix[Isym][I][A] * tia.matrix[Jsym][J][B]);
160
if((Isym==Bsym) && (Jsym==Asym))
161
tauijab.matrix[h][ij][ab] -=
162
(tia.matrix[Isym][I][B] * tia.matrix[Jsym][J][A]);
167
dpd_buf4_mat_irrep_wrt(&tauijab, h);
168
dpd_buf4_mat_irrep_close(&tauijab, h);
171
dpd_buf4_close(&tauijab);
173
dpd_buf4_init(&tauIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tauIjAb");
175
for(h=0; h < nirreps; h++) {
177
dpd_buf4_mat_irrep_init(&tauIjAb, h);
178
dpd_buf4_mat_irrep_rd(&tauIjAb, h);
180
for(ij=0; ij < tauIjAb.params->rowtot[h]; ij++) {
181
i = tauIjAb.params->roworb[h][ij][0];
182
j = tauIjAb.params->roworb[h][ij][1];
183
I = tIA.params->rowidx[i];
184
J = tia.params->rowidx[j];
185
Isym = tIA.params->psym[i];
186
Jsym = tia.params->psym[j];
187
for(ab=0; ab < tauIjAb.params->coltot[h]; ab++) {
188
a = tauIjAb.params->colorb[h][ab][0];
189
b = tauIjAb.params->colorb[h][ab][1];
190
A = tIA.params->colidx[a];
191
B = tia.params->colidx[b];
192
Asym = tIA.params->qsym[a];
193
Bsym = tia.params->qsym[b];
195
if((Isym==Asym) && (Jsym==Bsym))
196
tauIjAb.matrix[h][ij][ab] +=
197
(tIA.matrix[Isym][I][A] * tia.matrix[Jsym][J][B]);
202
dpd_buf4_mat_irrep_wrt(&tauIjAb, h);
203
dpd_buf4_mat_irrep_close(&tauIjAb, h);
206
/* This will generate the tauBA and tauIjbA files from tauIjAb */
207
dpd_buf4_sort(&tauIjAb, CC_TAMPS, pqsr, 0, 5, "tauIjbA");
208
dpd_buf4_close(&tauIjAb);
209
dpd_buf4_init(&tauIjbA, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tauIjbA");
210
dpd_buf4_sort(&tauIjbA, CC_TAMPS, qprs, 0, 5, "tauiJaB");
211
dpd_buf4_close(&tauIjbA);
213
dpd_file2_mat_close(&tIA);
214
dpd_file2_close(&tIA);
215
dpd_file2_mat_close(&tia);
216
dpd_file2_close(&tia);
219
else if(params.ref == 2) { /*** UHF ***/
221
dpd_buf4_init(&tIJAB, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tIJAB");
222
dpd_buf4_copy(&tIJAB, CC_TAMPS, "tauIJAB");
223
dpd_buf4_close(&tIJAB);
225
dpd_buf4_init(&tijab, CC_TAMPS, 0, 12, 17, 12, 17, 0, "tijab");
226
dpd_buf4_copy(&tijab, CC_TAMPS, "tauijab");
227
dpd_buf4_close(&tijab);
229
dpd_buf4_init(&tIjAb, CC_TAMPS, 0, 22, 28, 22, 28, 0, "tIjAb");
230
dpd_buf4_copy(&tIjAb, CC_TAMPS, "tauIjAb");
231
dpd_buf4_close(&tIjAb);
233
dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA");
234
dpd_file2_mat_init(&tIA);
235
dpd_file2_mat_rd(&tIA);
236
dpd_file2_init(&tia, CC_OEI, 0, 2, 3, "tia");
237
dpd_file2_mat_init(&tia);
238
dpd_file2_mat_rd(&tia);
240
dpd_buf4_init(&tauIJAB, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tauIJAB");
241
for(h=0; h < nirreps; h++) {
242
dpd_buf4_mat_irrep_init(&tauIJAB, h);
243
dpd_buf4_mat_irrep_rd(&tauIJAB, h);
244
for(ij=0; ij < tauIJAB.params->rowtot[h]; ij++) {
245
i = tauIJAB.params->roworb[h][ij][0];
246
j = tauIJAB.params->roworb[h][ij][1];
247
I = tIA.params->rowidx[i];
248
J = tIA.params->rowidx[j];
249
Isym = tIA.params->psym[i];
250
Jsym = tIA.params->psym[j];
251
for(ab=0; ab < tauIJAB.params->coltot[h]; ab++) {
252
a = tauIJAB.params->colorb[h][ab][0];
253
b = tauIJAB.params->colorb[h][ab][1];
254
A = tIA.params->colidx[a];
255
B = tIA.params->colidx[b];
256
Asym = tIA.params->qsym[a];
257
Bsym = tIA.params->qsym[b];
258
if((Isym==Asym) && (Jsym==Bsym))
259
tauIJAB.matrix[h][ij][ab] +=
260
(tIA.matrix[Isym][I][A] * tIA.matrix[Jsym][J][B]);
261
if((Isym==Bsym) && (Jsym==Asym))
262
tauIJAB.matrix[h][ij][ab] -=
263
(tIA.matrix[Isym][I][B] * tIA.matrix[Jsym][J][A]);
266
dpd_buf4_mat_irrep_wrt(&tauIJAB, h);
267
dpd_buf4_mat_irrep_close(&tauIJAB, h);
269
dpd_buf4_close(&tauIJAB);
271
dpd_buf4_init(&tauijab, CC_TAMPS, 0, 12, 17, 12, 17, 0, "tauijab");
272
for(h=0; h < nirreps; h++) {
273
dpd_buf4_mat_irrep_init(&tauijab, h);
274
dpd_buf4_mat_irrep_rd(&tauijab, h);
275
for(ij=0; ij < tauijab.params->rowtot[h]; ij++) {
276
i = tauijab.params->roworb[h][ij][0];
277
j = tauijab.params->roworb[h][ij][1];
278
I = tia.params->rowidx[i];
279
J = tia.params->rowidx[j];
280
Isym = tia.params->psym[i];
281
Jsym = tia.params->psym[j];
282
for(ab=0; ab < tauijab.params->coltot[h]; ab++) {
283
a = tauijab.params->colorb[h][ab][0];
284
b = tauijab.params->colorb[h][ab][1];
285
A = tia.params->colidx[a];
286
B = tia.params->colidx[b];
287
Asym = tia.params->qsym[a];
288
Bsym = tia.params->qsym[b];
289
if((Isym==Asym) && (Jsym==Bsym))
290
tauijab.matrix[h][ij][ab] +=
291
(tia.matrix[Isym][I][A] * tia.matrix[Jsym][J][B]);
292
if((Isym==Bsym) && (Jsym==Asym))
293
tauijab.matrix[h][ij][ab] -=
294
(tia.matrix[Isym][I][B] * tia.matrix[Jsym][J][A]);
297
dpd_buf4_mat_irrep_wrt(&tauijab, h);
298
dpd_buf4_mat_irrep_close(&tauijab, h);
300
dpd_buf4_close(&tauijab);
302
dpd_buf4_init(&tauIjAb, CC_TAMPS, 0, 22, 28, 22, 28, 0, "tauIjAb");
303
for(h=0; h < nirreps; h++) {
304
dpd_buf4_mat_irrep_init(&tauIjAb, h);
305
dpd_buf4_mat_irrep_rd(&tauIjAb, h);
306
for(ij=0; ij < tauIjAb.params->rowtot[h]; ij++) {
307
i = tauIjAb.params->roworb[h][ij][0];
308
j = tauIjAb.params->roworb[h][ij][1];
309
I = tIA.params->rowidx[i];
310
J = tia.params->rowidx[j];
311
Isym = tIA.params->psym[i];
312
Jsym = tia.params->psym[j];
313
for(ab=0; ab < tauIjAb.params->coltot[h]; ab++) {
314
a = tauIjAb.params->colorb[h][ab][0];
315
b = tauIjAb.params->colorb[h][ab][1];
316
A = tIA.params->colidx[a];
317
B = tia.params->colidx[b];
318
Asym = tIA.params->qsym[a];
319
Bsym = tia.params->qsym[b];
320
if((Isym==Asym) && (Jsym==Bsym))
321
tauIjAb.matrix[h][ij][ab] +=
322
(tIA.matrix[Isym][I][A] * tia.matrix[Jsym][J][B]);
325
dpd_buf4_mat_irrep_wrt(&tauIjAb, h);
326
dpd_buf4_mat_irrep_close(&tauIjAb, h);
328
dpd_buf4_close(&tauIjAb);
330
dpd_file2_mat_close(&tIA);
331
dpd_file2_close(&tIA);
332
dpd_file2_mat_close(&tia);
333
dpd_file2_close(&tia);
335
/* Resort IjAb to IjbA (used in Z.c) and iJaB */
336
dpd_buf4_init(&tauIjAb, CC_TAMPS, 0, 22, 28, 22, 28, 0, "tauIjAb");
337
dpd_buf4_sort(&tauIjAb, CC_TAMPS, pqsr, 22, 29, "tauIjbA");
338
dpd_buf4_close(&tauIjAb);
339
dpd_buf4_init(&tauIjbA, CC_TAMPS, 0, 22, 29, 22, 29, 0, "tauIjbA");
340
dpd_buf4_sort(&tauIjbA, CC_TAMPS, qprs, 23, 29, "tauiJaB");
341
dpd_buf4_close(&tauIjbA);
345
}} // namespace psi::ccenergy