← Back to branch summary

~ubuntu-branches/ubuntu/karmic/psicode/karmic

« back to all changes in this revision

Viewing changes to src/bin/ccdensity/Iij.cc

  • Committer: Bazaar Package Importer
  • Author(s): Michael Banck, Michael Banck, Daniel Leidert
  • Date: 2009-02-23 00:12:02 UTC
  • mfrom: (1.1.2 upstream)
  • Revision ID: james.westby@ubuntu.com-20090223001202-rutldoy3dimfpesc
Tags: 3.4.0-1
* New upstream release.

[ Michael Banck ]
* debian/patches/01_DESTDIR.dpatch: Refreshed.
* debian/patches/02_FHS.dpatch: Removed, applied upstream.
* debian/patches/03_debian_docdir: Likewise.
* debian/patches/04_man.dpatch: Likewise.
* debian/patches/06_466828_fix_gcc_43_ftbfs.dpatch: Likewise.
* debian/patches/07_464867_move_executables: Fixed and refreshed.
* debian/patches/00list: Adjusted.
* debian/control: Improved description.
* debian/patches-held: Removed.
* debian/rules (install/psi3): Do not ship the ruby bindings for now.

[ Daniel Leidert ]
* debian/rules: Fix txtdir via DEB_MAKE_INSTALL_TARGET.
* debian/patches/01_DESTDIR.dpatch: Refreshed.

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/bin/ccdensity/Iij.cc
 
1
/*! \file
 
2
    \ingroup CCDENSITY
 
3
    \brief Enter brief description of file here 
 
4
*/
 
5
#include <libdpd/dpd.h>
 
6
#include "MOInfo.h"
 
7
#include "Params.h"
 
8
#include "Frozen.h"
 
9
#define EXTERN
 
10
#include "globals.h"
 
11
 
 
12
namespace psi { namespace ccdensity {
 
13
 
 
14
    /* Iij(): Build the occupied-occupied block of the orbital Lagrangian
 
15
    ** using the expression given in lag.c.  Note that we include an
 
16
    ** addition term here referred to as the reference contribution,
 
17
    ** 2*fij.  This comes from the general spin-orbital SCF gradient
 
18
    ** expression and is present for all reference types (though for
 
19
    ** canonical unperturbed orbitals, only the diagonal elements
 
20
    ** contribute, of course).  
 
21
    **
 
22
    ** Note that the code currently produces only the unique I_IJ terms,
 
23
    ** but the actual contractions still need to be spin-adapted for
 
24
    ** greater efficiency.
 
25
    **
 
26
    ** TDC, 2/2008
 
27
    */
 
28
 
 
29
    void Iij(struct RHO_Params rho_params)
 
30
    {
 
31
      dpdfile2 I, F, D;
 
32
      dpdbuf4 G, Aints, Fints, Dints, Cints, Eints;
 
33
 
 
34
      if(params.ref == 0) { /** RHF **/
 
35
        /* I'IJ <-- sum_K fIK (DJK + DKJ) + sum_A fIA (DJA + DAJ) */
 
36
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
37
 
 
38
        dpd_file2_init(&F, CC_OEI, 0, 0, 0, "fIJ");
 
39
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
40
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 0.0);
 
41
        dpd_contract222(&F, &D, &I, 0, 1, 1.0, 1.0);
 
42
        dpd_file2_close(&D);
 
43
 
 
44
        /* Add reference contribution: I'IJ <-- 2 fIJ */
 
45
        dpd_file2_axpy(&F, &I, 2.0, 0);
 
46
        dpd_file2_close(&F);
 
47
 
 
48
        dpd_file2_init(&F, CC_OEI, 0, 0, 1, "fIA");
 
49
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
50
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
51
        dpd_file2_close(&D);
 
52
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
53
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
54
        dpd_file2_close(&D);
 
55
        dpd_file2_close(&F);
 
56
 
 
57
        dpd_file2_close(&I);
 
58
      }
 
59
      else if(params.ref == 1) { /** ROHF **/
 
60
 
 
61
        /* I'IJ <-- sum_K fIK (DJK + DKJ) + sum_A fIA (DJA + DAJ) */
 
62
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
63
 
 
64
        dpd_file2_init(&F, CC_OEI, 0, 0, 0, "fIJ");
 
65
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
66
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 0.0);
 
67
        dpd_contract222(&F, &D, &I, 0, 1, 1.0, 1.0);
 
68
        dpd_file2_close(&D);
 
69
 
 
70
        /* Add reference contribution: I'IJ <-- 2 fIJ */
 
71
        dpd_file2_axpy(&F, &I, 2.0, 0);
 
72
        dpd_file2_close(&F);
 
73
 
 
74
        dpd_file2_init(&F, CC_OEI, 0, 0, 1, "fIA");
 
75
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
76
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
77
        dpd_file2_close(&D);
 
78
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
79
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
80
        dpd_file2_close(&D);
 
81
        dpd_file2_close(&F);
 
82
 
 
83
        dpd_file2_close(&I);
 
84
 
 
85
        /* I'ij <-- sum_k fik (Djk + Dkj)  + sum_a fia (Dja + Daj) */
 
86
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
87
 
 
88
        dpd_file2_init(&F, CC_OEI, 0, 0, 0, "fij");
 
89
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.Dij_lbl);
 
90
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 0.0);
 
91
        dpd_contract222(&F, &D, &I, 0, 1, 1.0, 1.0);
 
92
        dpd_file2_close(&D);
 
93
 
 
94
        /* Add reference contribution: I'ij <-- 2 fij */
 
95
        dpd_file2_axpy(&F, &I, 2.0, 0);
 
96
        dpd_file2_close(&F);
 
97
 
 
98
        dpd_file2_init(&F, CC_OEI, 0, 0, 1, "fia");
 
99
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
100
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
101
        dpd_file2_close(&D);
 
102
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
103
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
104
        dpd_file2_close(&D);
 
105
        dpd_file2_close(&F);
 
106
 
 
107
        dpd_file2_close(&I);
 
108
      }
 
109
      else if(params.ref == 2) { /** UHF **/
 
110
 
 
111
        /* I'IJ <-- sum_K fIK (DJK + DKJ) + sum_A fIA (DJA + DAJ) */
 
112
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
113
 
 
114
        dpd_file2_init(&F, CC_OEI, 0, 0, 0, "fIJ");
 
115
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
116
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 0.0);
 
117
        dpd_contract222(&F, &D, &I, 0, 1, 1.0, 1.0);
 
118
        dpd_file2_close(&D);
 
119
 
 
120
        /* Add reference contribution: I'IJ <-- 2 fIJ */
 
121
        dpd_file2_axpy(&F, &I, 2.0, 0);
 
122
        dpd_file2_close(&F);
 
123
 
 
124
        dpd_file2_init(&F, CC_OEI, 0, 0, 1, "fIA");
 
125
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
126
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
127
        dpd_file2_close(&D);
 
128
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
129
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
130
        dpd_file2_close(&D);
 
131
        dpd_file2_close(&F);
 
132
 
 
133
        dpd_file2_close(&I);
 
134
 
 
135
        /* I'ij <-- sum_k fik (Djk + Dkj)  + sum_a fia (Dja + Daj) */
 
136
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
137
 
 
138
        dpd_file2_init(&F, CC_OEI, 0, 2, 2, "fij");
 
139
        dpd_file2_init(&D, CC_OEI, 0, 2, 2, rho_params.Dij_lbl);
 
140
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 0.0);
 
141
        dpd_contract222(&F, &D, &I, 0, 1, 1.0, 1.0);
 
142
        dpd_file2_close(&D);
 
143
 
 
144
        /* Add reference contribution: I'ij <-- 2 fij */
 
145
        dpd_file2_axpy(&F, &I, 2.0, 0);
 
146
        dpd_file2_close(&F);
 
147
 
 
148
        dpd_file2_init(&F, CC_OEI, 0, 2, 3, "fia");
 
149
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dia_lbl);
 
150
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
151
        dpd_file2_close(&D);
 
152
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dai_lbl);
 
153
        dpd_contract222(&F, &D, &I, 0, 0, 1.0, 1.0);
 
154
        dpd_file2_close(&D);
 
155
        dpd_file2_close(&F);
 
156
 
 
157
        dpd_file2_close(&I);
 
158
      }
 
159
 
 
160
      if(params.ref == 0) { /** RHF **/
 
161
        /* I'IJ <-- sum_KL <IK||JL> (D_KL + D_LK) + sum_kl <Ik|Jl> (D_kl + D_lk) */
 
162
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
163
 
 
164
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
165
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 1, "A <ij|kl>");
 
166
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
167
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
168
        dpd_buf4_close(&Aints);
 
169
        dpd_file2_close(&D);
 
170
 
 
171
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.Dij_lbl);
 
172
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 0, "A <ij|kl>");
 
173
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
174
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
175
        dpd_buf4_close(&Aints);
 
176
        dpd_file2_close(&D);
 
177
  
 
178
        dpd_file2_close(&I);
 
179
      }
 
180
      else if(params.ref == 1) { /** ROHF **/
 
181
 
 
182
        /* I'IJ <-- sum_KL <IK||JL> (D_KL + D_LK) + sum_kl <Ik|Jl> (D_kl + D_lk) */
 
183
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
184
 
 
185
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
186
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 1, "A <ij|kl>");
 
187
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
188
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
189
        dpd_buf4_close(&Aints);
 
190
        dpd_file2_close(&D);
 
191
 
 
192
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.Dij_lbl);
 
193
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 0, "A <ij|kl>");
 
194
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
195
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
196
        dpd_buf4_close(&Aints);
 
197
        dpd_file2_close(&D);
 
198
  
 
199
        dpd_file2_close(&I);
 
200
 
 
201
        /* I'ij <-- sum_kl <ik||jl> (D_kl + D_lk) + sum_KL <iK|jL> (D_KL + D_LK) */
 
202
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
203
 
 
204
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.Dij_lbl);
 
205
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 1, "A <ij|kl>");
 
206
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
207
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
208
        dpd_buf4_close(&Aints);
 
209
        dpd_file2_close(&D);
 
210
 
 
211
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
212
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 0, "A <ij|kl>");
 
213
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
214
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
215
        dpd_buf4_close(&Aints);
 
216
        dpd_file2_close(&D);
 
217
  
 
218
        dpd_file2_close(&I);
 
219
      }
 
220
      else if(params.ref == 2) { /** UHF **/
 
221
 
 
222
        /* I'IJ <-- sum_KL <IK||JL> (D_KL + D_LK) + sum_kl <Ik|Jl> (D_kl + D_lk) */
 
223
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
224
 
 
225
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
226
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 1, "A <IJ|KL>");
 
227
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
228
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
229
        dpd_buf4_close(&Aints);
 
230
        dpd_file2_close(&D);
 
231
 
 
232
        dpd_file2_init(&D, CC_OEI, 0, 2, 2, rho_params.Dij_lbl);
 
233
        dpd_buf4_init(&Aints, CC_AINTS, 0, 22, 22, 22, 22, 0, "A <Ij|Kl>");
 
234
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
235
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
236
        dpd_buf4_close(&Aints);
 
237
        dpd_file2_close(&D);
 
238
  
 
239
        dpd_file2_close(&I);
 
240
 
 
241
        /* I'ij <-- sum_kl <ik||jl> (D_kl + D_lk) + sum_KL <iK|jL> (D_KL + D_LK) */
 
242
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
243
 
 
244
        dpd_file2_init(&D, CC_OEI, 0, 2, 2, rho_params.Dij_lbl);
 
245
        dpd_buf4_init(&Aints, CC_AINTS, 0, 10, 10, 10, 10, 1, "A <ij|kl>");
 
246
        dpd_dot24(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
247
        dpd_dot24(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
248
        dpd_buf4_close(&Aints);
 
249
        dpd_file2_close(&D);
 
250
 
 
251
        dpd_file2_init(&D, CC_OEI, 0, 0, 0, rho_params.DIJ_lbl);
 
252
        dpd_buf4_init(&Aints, CC_AINTS, 0, 22, 22, 22, 22, 0, "A <Ij|Kl>");
 
253
        dpd_dot13(&D, &Aints, &I, 0, 0, 1.0, 1.0);
 
254
        dpd_dot13(&D, &Aints, &I, 1, 0, 1.0, 1.0);
 
255
        dpd_buf4_close(&Aints);
 
256
        dpd_file2_close(&D);
 
257
  
 
258
        dpd_file2_close(&I);
 
259
 
 
260
      }
 
261
 
 
262
      if(params.ref == 0) { /** RHF **/
 
263
        /* I'IJ <-- sum_KA <IK||JA> (D_KA + D_AK) + sum_ka <Ik|Ja> (D_ka + D_ak) */
 
264
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
265
 
 
266
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
267
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
268
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
269
        dpd_file2_close(&D);
 
270
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
271
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
272
        dpd_file2_close(&D);
 
273
        dpd_buf4_close(&Eints);
 
274
 
 
275
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
276
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
277
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
278
        dpd_file2_close(&D);
 
279
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
280
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
281
        dpd_file2_close(&D);
 
282
        dpd_buf4_close(&Eints);
 
283
  
 
284
        dpd_file2_close(&I);
 
285
      }
 
286
      else if(params.ref == 1) { /** ROHF **/
 
287
 
 
288
        /* I'IJ <-- sum_KA <IK||JA> (D_KA + D_AK) + sum_ka <Ik|Ja> (D_ka + D_ak) */
 
289
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
290
 
 
291
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
292
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
293
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
294
        dpd_file2_close(&D);
 
295
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
296
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
297
        dpd_file2_close(&D);
 
298
        dpd_buf4_close(&Eints);
 
299
 
 
300
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
301
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
302
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
303
        dpd_file2_close(&D);
 
304
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
305
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
306
        dpd_file2_close(&D);
 
307
        dpd_buf4_close(&Eints);
 
308
  
 
309
        dpd_file2_close(&I);
 
310
 
 
311
        /* I'ij <-- sum_ka <ik||ja> (D_ka + D_ak) + sum_KA <iK|jA> (D_KA + D_AK) */
 
312
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
313
 
 
314
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
315
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
316
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
317
        dpd_file2_close(&D);
 
318
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
319
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
320
        dpd_file2_close(&D);
 
321
        dpd_buf4_close(&Eints);
 
322
 
 
323
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
324
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
325
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
326
        dpd_file2_close(&D);
 
327
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
328
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
329
        dpd_file2_close(&D);
 
330
        dpd_buf4_close(&Eints);
 
331
  
 
332
        dpd_file2_close(&I);
 
333
      }
 
334
      else if(params.ref == 2) { /** UHF **/
 
335
 
 
336
        /* I'IJ <-- sum_KA <IK||JA> (D_KA + D_AK) + sum_ka <Ik|Ja> (D_ka + D_ak) */
 
337
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
338
 
 
339
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 20, 2, 20, 0, "E <IJ||KA> (I>J,KA)");
 
340
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
341
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
342
        dpd_file2_close(&D);
 
343
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
344
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
345
        dpd_file2_close(&D);
 
346
        dpd_buf4_close(&Eints);
 
347
 
 
348
        dpd_buf4_init(&Eints, CC_EINTS, 0, 22, 24, 22, 24, 0, "E <Ij|Ka>");
 
349
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dia_lbl);
 
350
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
351
        dpd_file2_close(&D);
 
352
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dai_lbl);
 
353
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
354
        dpd_file2_close(&D);
 
355
        dpd_buf4_close(&Eints);
 
356
  
 
357
        dpd_file2_close(&I);
 
358
 
 
359
        /* I'ij <-- sum_ka <ik||ja> (D_ka + D_ak) + sum_KA <iK|jA> (D_KA + D_AK) */
 
360
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
361
 
 
362
        dpd_buf4_init(&Eints, CC_EINTS, 0, 10, 30, 12, 30, 0, "E <ij||ka> (i>j,ka)");
 
363
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dia_lbl);
 
364
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
365
        dpd_file2_close(&D);
 
366
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dai_lbl);
 
367
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
368
        dpd_file2_close(&D);
 
369
        dpd_buf4_close(&Eints);
 
370
 
 
371
        dpd_buf4_init(&Eints, CC_EINTS, 0, 23, 27, 23, 27, 0, "E <iJ|kA>");
 
372
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
373
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
374
        dpd_file2_close(&D);
 
375
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
376
        dpd_dot24(&D, &Eints, &I, 0, 0, 1.0, 1.0);
 
377
        dpd_file2_close(&D);
 
378
        dpd_buf4_close(&Eints);
 
379
  
 
380
        dpd_file2_close(&I);
 
381
      }
 
382
 
 
383
      if(params.ref == 0) { /** RHF **/
 
384
        /* I'IJ <-- sum_AK <JK||IA> (D_AK + D_KA) + sum_ak <Jk|Ia> (D_ak + D_ka) */
 
385
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
386
 
 
387
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
388
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
389
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
390
        dpd_file2_close(&D);
 
391
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
392
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
393
        dpd_file2_close(&D);
 
394
        dpd_buf4_close(&Eints);
 
395
 
 
396
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
397
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
398
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
399
        dpd_file2_close(&D);
 
400
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
401
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
402
        dpd_file2_close(&D);
 
403
        dpd_buf4_close(&Eints);
 
404
  
 
405
        dpd_file2_close(&I);
 
406
      }
 
407
      else if(params.ref == 1) { /** ROHF **/
 
408
 
 
409
        /* I'IJ <-- sum_AK <JK||IA> (D_AK + D_KA) + sum_ak <Jk|Ia> (D_ak + D_ka) */
 
410
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
411
 
 
412
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
413
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
414
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
415
        dpd_file2_close(&D);
 
416
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
417
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
418
        dpd_file2_close(&D);
 
419
        dpd_buf4_close(&Eints);
 
420
 
 
421
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
422
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
423
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
424
        dpd_file2_close(&D);
 
425
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
426
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
427
        dpd_file2_close(&D);
 
428
        dpd_buf4_close(&Eints);
 
429
  
 
430
        dpd_file2_close(&I);
 
431
 
 
432
        /* I'ij <-- sum_ak <jk||ia> (D_ak + D_ka) + sum_AK <jK|iA> (D_AK + D_KA) */
 
433
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
434
 
 
435
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
436
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dia_lbl);
 
437
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
438
        dpd_file2_close(&D);
 
439
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.Dai_lbl);
 
440
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
441
        dpd_file2_close(&D);
 
442
        dpd_buf4_close(&Eints);
 
443
 
 
444
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
445
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
446
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
447
        dpd_file2_close(&D);
 
448
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
449
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
450
        dpd_file2_close(&D);
 
451
        dpd_buf4_close(&Eints);
 
452
  
 
453
        dpd_file2_close(&I);
 
454
      }
 
455
      else if(params.ref == 2) { /** UHF **/
 
456
 
 
457
        /* I'IJ <-- sum_AK <JK||IA> (D_AK + D_KA) + sum_ak <Jk|Ia> (D_ak + D_ka) */
 
458
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
459
 
 
460
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 20, 2, 20, 0, "E <IJ||KA> (I>J,KA)");
 
461
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
462
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
463
        dpd_file2_close(&D);
 
464
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
465
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
466
        dpd_file2_close(&D);
 
467
        dpd_buf4_close(&Eints);
 
468
 
 
469
        dpd_buf4_init(&Eints, CC_EINTS, 0, 22, 24, 22, 24, 0, "E <Ij|Ka>");
 
470
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dia_lbl);
 
471
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
472
        dpd_file2_close(&D);
 
473
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dai_lbl);
 
474
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
475
        dpd_file2_close(&D);
 
476
        dpd_buf4_close(&Eints);
 
477
  
 
478
        dpd_file2_close(&I);
 
479
 
 
480
        /* I'ij <-- sum_ak <jk||ia> (D_ak + D_ka) + sum_AK <jK|iA> (D_AK + D_KA) */
 
481
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
482
 
 
483
        dpd_buf4_init(&Eints, CC_EINTS, 0, 10, 30, 12, 30, 0, "E <ij||ka> (i>j,ka)");
 
484
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dia_lbl);
 
485
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
486
        dpd_file2_close(&D);
 
487
        dpd_file2_init(&D, CC_OEI, 0, 2, 3, rho_params.Dai_lbl);
 
488
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
489
        dpd_file2_close(&D);
 
490
        dpd_buf4_close(&Eints);
 
491
 
 
492
        dpd_buf4_init(&Eints, CC_EINTS, 0, 23, 27, 23, 27, 0, "E <iJ|kA>");
 
493
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DIA_lbl);
 
494
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
495
        dpd_file2_close(&D);
 
496
        dpd_file2_init(&D, CC_OEI, 0, 0, 1, rho_params.DAI_lbl);
 
497
        dpd_dot24(&D, &Eints, &I, 0, 1, 1.0, 1.0);
 
498
        dpd_file2_close(&D);
 
499
        dpd_buf4_close(&Eints);
 
500
  
 
501
        dpd_file2_close(&I);
 
502
 
 
503
      }
 
504
 
 
505
      if(params.ref == 0) { /** RHF **/
 
506
        /* I'IJ <-- sum_AB <IA||JB> (D_AB + D_BA) + sum_ab <Ia|Jb> (D_ab + D_ba) */
 
507
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
508
 
 
509
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.DAB_lbl);
 
510
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia||jb>");
 
511
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
512
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
513
        dpd_buf4_close(&Cints);
 
514
        dpd_file2_close(&D);
 
515
 
 
516
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.Dab_lbl);
 
517
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia|jb>");
 
518
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
519
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
520
        dpd_buf4_close(&Cints);
 
521
        dpd_file2_close(&D);
 
522
 
 
523
        dpd_file2_close(&I);
 
524
      }
 
525
      else if(params.ref == 1) { /** ROHF **/
 
526
 
 
527
        /* I'IJ <-- sum_AB <IA||JB> (D_AB + D_BA) + sum_ab <Ia|Jb> (D_ab + D_ba) */
 
528
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
529
 
 
530
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.DAB_lbl);
 
531
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia||jb>");
 
532
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
533
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
534
        dpd_buf4_close(&Cints);
 
535
        dpd_file2_close(&D);
 
536
 
 
537
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.Dab_lbl);
 
538
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia|jb>");
 
539
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
540
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
541
        dpd_buf4_close(&Cints);
 
542
        dpd_file2_close(&D);
 
543
 
 
544
        dpd_file2_close(&I);
 
545
 
 
546
        /* I'ij <-- sum_ab <ia||jb> (D_ab + D_ba) + sum_AB <iA|jB> (D_AB + D_BA) */
 
547
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
548
 
 
549
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.Dab_lbl);
 
550
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia||jb>");
 
551
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
552
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
553
        dpd_buf4_close(&Cints);
 
554
        dpd_file2_close(&D);
 
555
 
 
556
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.DAB_lbl);
 
557
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia|jb>");
 
558
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
559
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
560
        dpd_buf4_close(&Cints);
 
561
        dpd_file2_close(&D);
 
562
 
 
563
        dpd_file2_close(&I);
 
564
      }
 
565
      else if(params.ref == 2) { /** UHF **/
 
566
 
 
567
        /* I'IJ <-- sum_AB <IA||JB> (D_AB + D_BA) + sum_ab <Ia|Jb> (D_ab + D_ba) */
 
568
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
569
 
 
570
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.DAB_lbl);
 
571
        dpd_buf4_init(&Cints, CC_CINTS, 0, 20, 20, 20, 20, 0, "C <IA||JB>");
 
572
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
573
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
574
        dpd_buf4_close(&Cints);
 
575
        dpd_file2_close(&D);
 
576
 
 
577
        dpd_file2_init(&D, CC_OEI, 0, 3, 3, rho_params.Dab_lbl);
 
578
        dpd_buf4_init(&Cints, CC_CINTS, 0, 24, 24, 24, 24, 0, "C <Ia|Jb>");
 
579
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
580
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
581
        dpd_buf4_close(&Cints);
 
582
        dpd_file2_close(&D);
 
583
 
 
584
        dpd_file2_close(&I);
 
585
 
 
586
        /* I'ij <-- sum_ab <ia||jb> (D_ab + D_ba) + sum_AB <iA|jB> (D_AB + D_BA) */
 
587
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
588
 
 
589
        dpd_file2_init(&D, CC_OEI, 0, 3, 3, rho_params.Dab_lbl);
 
590
        dpd_buf4_init(&Cints, CC_CINTS, 0, 30, 30, 30, 30, 0, "C <ia||jb>");
 
591
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
592
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
593
        dpd_buf4_close(&Cints);
 
594
        dpd_file2_close(&D);
 
595
 
 
596
        dpd_file2_init(&D, CC_OEI, 0, 1, 1, rho_params.DAB_lbl);
 
597
        dpd_buf4_init(&Cints, CC_CINTS, 0, 27, 27, 27, 27, 0, "C <iA|jB>");
 
598
        dpd_dot24(&D, &Cints, &I, 0, 0, 1.0, 1.0);
 
599
        dpd_dot24(&D, &Cints, &I, 1, 0, 1.0, 1.0);
 
600
        dpd_buf4_close(&Cints);
 
601
        dpd_file2_close(&D);
 
602
 
 
603
        dpd_file2_close(&I);
 
604
 
 
605
      }
 
606
 
 
607
      if(params.ref == 0) { /** RHF **/
 
608
        /* I'IJ <-- 2 [ 2 (Gjklm - Gjkml) <ik|lm>] */
 
609
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
610
 
 
611
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 0, 0, 0, 0, "GIjKl");
 
612
        dpd_buf4_scmcopy(&G, CC_GAMMA, "2 Gijkl - Gijlk", 2);
 
613
        dpd_buf4_sort_axpy(&G, CC_GAMMA, pqsr, 0, 0, "2 Gijkl - Gijlk", -1);
 
614
        dpd_buf4_close(&G);
 
615
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 0, 0, 0, 0, "2 Gijkl - Gijlk");
 
616
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 0, "A <ij|kl>");
 
617
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0); 
 
618
        dpd_buf4_close(&Aints);
 
619
        dpd_buf4_close(&G);
 
620
        dpd_file2_close(&I);
 
621
      }
 
622
      else if(params.ref == 1) { /** ROHF **/
 
623
 
 
624
        /* I'IJ <-- sum_KLM <IK||LM> G(JK,LM) + 2 sum_kLm <Ik|Lm> G(Jk,Lm) */
 
625
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
626
 
 
627
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 2, 0, 0, 1, "A <ij|kl>");
 
628
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 2, 2, 2, 0, "GIJKL");
 
629
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0);
 
630
        dpd_buf4_close(&G);
 
631
        dpd_buf4_close(&Aints);
 
632
 
 
633
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 0, "A <ij|kl>");
 
634
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 0, 0, 0, 0, "GIjKl");
 
635
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0);
 
636
        dpd_buf4_close(&G);
 
637
        dpd_buf4_close(&Aints);
 
638
 
 
639
        dpd_file2_close(&I);
 
640
 
 
641
        /* I'ij <-- sum_klm <ik||lm> G(jk,lm) + 2 sum_KlM <Ki|Ml> G(Kj,Ml) */
 
642
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
643
 
 
644
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 2, 0, 0, 1, "A <ij|kl>");
 
645
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 2, 2, 2, 0, "Gijkl");
 
646
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0);
 
647
        dpd_buf4_close(&G);
 
648
        dpd_buf4_close(&Aints);
 
649
 
 
650
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 0, 0, 0, 0, "A <ij|kl>");
 
651
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 0, 0, 0, 0, "GIjKl");
 
652
        dpd_contract442(&Aints, &G, &I, 1, 1, 2.0, 1.0);
 
653
        dpd_buf4_close(&G);
 
654
        dpd_buf4_close(&Aints);
 
655
 
 
656
        dpd_file2_close(&I);
 
657
 
 
658
      }
 
659
      else if(params.ref == 2) { /** UHF **/
 
660
 
 
661
        /* I'IJ <-- sum_KLM <IK||LM> G(JK,LM) + 2 sum_kLm <Ik|Lm> G(Jk,Lm) */
 
662
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
663
 
 
664
        dpd_buf4_init(&Aints, CC_AINTS, 0, 0, 2, 0, 0, 1, "A <IJ|KL>");
 
665
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 2, 2, 2, 0, "GIJKL");
 
666
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0);
 
667
        dpd_buf4_close(&G);
 
668
        dpd_buf4_close(&Aints);
 
669
 
 
670
        dpd_buf4_init(&Aints, CC_AINTS, 0, 22, 22, 22, 22, 0, "A <Ij|Kl>");
 
671
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 22, 22, 22, 0, "GIjKl");
 
672
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0);
 
673
        dpd_buf4_close(&G);
 
674
        dpd_buf4_close(&Aints);
 
675
 
 
676
        dpd_file2_close(&I);
 
677
 
 
678
        /* I'ij <-- sum_klm <ik||lm> G(jk,lm) + 2 sum_KlM <Ki|Ml> G(Kj,Ml) */
 
679
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
680
 
 
681
        dpd_buf4_init(&Aints, CC_AINTS, 0, 10, 12, 10, 10, 1, "A <ij|kl>");
 
682
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 12, 12, 12, 0, "Gijkl");
 
683
        dpd_contract442(&Aints, &G, &I, 0, 0, 2.0, 1.0);
 
684
        dpd_buf4_close(&G);
 
685
        dpd_buf4_close(&Aints);
 
686
 
 
687
        dpd_buf4_init(&Aints, CC_AINTS, 0, 22, 22, 22, 22, 0, "A <Ij|Kl>");
 
688
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 22, 22, 22, 0, "GIjKl");
 
689
        dpd_contract442(&Aints, &G, &I, 1, 1, 2.0, 1.0);
 
690
        dpd_buf4_close(&G);
 
691
        dpd_buf4_close(&Aints);
 
692
 
 
693
        dpd_file2_close(&I);
 
694
 
 
695
      }
 
696
 
 
697
      if(params.ref == 0) { /** RHF **/
 
698
        /* I'IJ <-- sum_ABC <IA||BC> G(JA,BC) + 2 sum_AbC <aI|bC> G(aJ,bC) */
 
699
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
700
 
 
701
        dpd_buf4_init(&G, CC_GAMMA, 0, 11, 5, 11, 5, 0, "GCiAb");
 
702
        dpd_buf4_scmcopy(&G, CC_GAMMA, "2 Gciab - Gciba", 2);
 
703
        dpd_buf4_sort_axpy(&G, CC_GAMMA, pqsr, 11, 5, "2 Gciab - Gciba", -1);
 
704
        dpd_buf4_close(&G);
 
705
 
 
706
        dpd_buf4_init(&G, CC_GAMMA, 0, 11, 5, 11, 5, 0, "2 Gciab - Gciba");
 
707
        dpd_buf4_sort(&G, CC_GAMMA, qpsr, 10, 5, "2 Giabc - Giacb");
 
708
        dpd_buf4_close(&G);
 
709
 
 
710
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 5, 10, 5, 0, "2 Giabc - Giacb");
 
711
        dpd_buf4_init(&Fints, CC_FINTS, 0, 10, 5, 10, 5, 0, "F <ia|bc>");
 
712
        dpd_contract442(&Fints, &G, &I, 0, 0, 2.0, 1.0);
 
713
        dpd_buf4_close(&Fints);
 
714
        dpd_buf4_close(&G);
 
715
 
 
716
        dpd_file2_close(&I);
 
717
      }
 
718
      else if(params.ref == 1) { /** ROHF **/
 
719
 
 
720
        /* I'IJ <-- sum_ABC <IA||BC> G(JA,BC) + 2 sum_AbC <aI|bC> G(aJ,bC) */
 
721
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
722
 
 
723
        dpd_buf4_init(&Fints, CC_FINTS, 0, 10, 7, 10, 5, 1, "F <ia|bc>");
 
724
        dpd_buf4_init(&G, CC_GAMMA, 0, 11, 7, 11, 7, 0, "GCIAB");
 
725
        dpd_buf4_sort(&G, CC_TMP0, qprs, 10, 7, "GICAB");
 
726
        dpd_buf4_close(&G);
 
727
        dpd_buf4_init(&G, CC_TMP0, 0, 10, 7, 10, 7, 0, "GICAB");
 
728
        dpd_contract442(&Fints, &G, &I, 0, 0, -2.0, 1.0);
 
729
        dpd_buf4_close(&G);
 
730
        dpd_buf4_close(&Fints);
 
731
 
 
732
        dpd_buf4_init(&Fints, CC_FINTS, 0, 10, 5, 10, 5, 0, "F <ia|bc>");
 
733
        dpd_buf4_init(&G, CC_GAMMA, 0, 11, 5, 11, 5, 0, "GcIaB");
 
734
        dpd_buf4_sort(&G, CC_TMP0, qpsr, 10, 5, "GIcBa");
 
735
        dpd_buf4_close(&G);
 
736
        dpd_buf4_init(&G, CC_TMP0, 0, 10, 5, 10, 5, 0, "GIcBa");
 
737
        dpd_contract442(&Fints, &G, &I, 0, 0, 2.0, 1.0);
 
738
        dpd_buf4_close(&G);
 
739
        dpd_buf4_close(&Fints);
 
740
 
 
741
        dpd_file2_close(&I);
 
742
 
 
743
        /* I'ij <-- sum_abc <ia||bcC> G(ja,bc) + 2 sum_AbC <Ai|Bc> G(Aj,Bc) */
 
744
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
745
 
 
746
        dpd_buf4_init(&Fints, CC_FINTS, 0, 10, 7, 10, 5, 1, "F <ia|bc>");
 
747
        dpd_buf4_init(&G, CC_GAMMA, 0, 11, 7, 11, 7, 0, "Gciab");
 
748
        dpd_buf4_sort(&G, CC_TMP0, qprs, 10, 7, "Gicab");
 
749
        dpd_buf4_close(&G);
 
750
        dpd_buf4_init(&G, CC_TMP0, 0, 10, 7, 10, 7, 0, "Gicab");
 
751
        dpd_contract442(&Fints, &G, &I, 0, 0, -2.0, 1.0);
 
752
        dpd_buf4_close(&G);
 
753
        dpd_buf4_close(&Fints);
 
754
 
 
755
        dpd_buf4_init(&Fints, CC_FINTS, 0, 10, 5, 10, 5, 0, "F <ia|bc>");
 
756
        dpd_buf4_init(&G, CC_GAMMA, 0, 11, 5, 11, 5, 0, "GCiAb");
 
757
        dpd_buf4_sort(&G, CC_TMP0, qpsr, 10, 5, "GiCbA");
 
758
        dpd_buf4_close(&G);
 
759
        dpd_buf4_init(&G, CC_TMP0, 0, 10, 5, 10, 5, 0, "GiCbA");
 
760
        dpd_contract442(&Fints, &G, &I, 0, 0, 2.0, 1.0);
 
761
        dpd_buf4_close(&G);
 
762
        dpd_buf4_close(&Fints);
 
763
 
 
764
        dpd_file2_close(&I);
 
765
 
 
766
      }
 
767
      else if(params.ref == 2) { /** UHF **/
 
768
 
 
769
        /* I'IJ <-- sum_ABC <IA||BC> G(JA,BC) + 2 sum_AbC <aI|bC> G(aJ,bC) */
 
770
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
771
 
 
772
        dpd_buf4_init(&Fints, CC_FINTS, 0, 20, 7, 20, 5, 1, "F <IA|BC>");
 
773
        dpd_buf4_init(&G, CC_GAMMA, 0, 21, 7, 21, 7, 0, "GCIAB");
 
774
        dpd_buf4_sort(&G, CC_TMP0, qprs, 20, 7, "GICAB");
 
775
        dpd_buf4_close(&G);
 
776
        dpd_buf4_init(&G, CC_TMP0, 0, 20, 7, 20, 7, 0, "GICAB");
 
777
        dpd_contract442(&Fints, &G, &I, 0, 0, -2.0, 1.0);
 
778
        dpd_buf4_close(&G);
 
779
        dpd_buf4_close(&Fints);
 
780
 
 
781
        dpd_buf4_init(&Fints, CC_FINTS, 0, 24, 28, 24, 28, 0, "F <Ia|Bc>");
 
782
        dpd_buf4_init(&G, CC_GAMMA, 0, 25, 29, 25, 29, 0, "GcIaB");
 
783
        dpd_buf4_sort(&G, CC_TMP0, qpsr, 24, 28, "GIcBa");
 
784
        dpd_buf4_close(&G);
 
785
        dpd_buf4_init(&G, CC_TMP0, 0, 24, 28, 24, 28, 0, "GIcBa");
 
786
        dpd_contract442(&Fints, &G, &I, 0, 0, 2.0, 1.0);
 
787
        dpd_buf4_close(&G);
 
788
        dpd_buf4_close(&Fints);
 
789
 
 
790
        dpd_file2_close(&I);
 
791
 
 
792
        /* I'ij <-- sum_abc <ia||bc> G(ja,bc) + 2 sum_AbC <Ai|Bc> G(Aj,Bc) */
 
793
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
794
 
 
795
        dpd_buf4_init(&Fints, CC_FINTS, 0, 30, 17, 30, 15, 1, "F <ia|bc>");
 
796
        dpd_buf4_init(&G, CC_GAMMA, 0, 31, 17, 31, 17, 0, "Gciab");
 
797
        dpd_buf4_sort(&G, CC_TMP0, qprs, 30, 17, "Gicab");
 
798
        dpd_buf4_close(&G);
 
799
        dpd_buf4_init(&G, CC_TMP0, 0, 30, 17, 30, 17, 0, "Gicab");
 
800
        dpd_contract442(&Fints, &G, &I, 0, 0, -2.0, 1.0);
 
801
        dpd_buf4_close(&G);
 
802
        dpd_buf4_close(&Fints);
 
803
 
 
804
        dpd_buf4_init(&Fints, CC_FINTS, 0, 27, 29, 27, 29, 0, "F <iA|bC>");
 
805
        dpd_buf4_init(&G, CC_GAMMA, 0, 26, 28, 26, 28, 0, "GCiAb");
 
806
        dpd_buf4_sort(&G, CC_TMP0, qpsr, 27, 29, "GiCbA");
 
807
        dpd_buf4_close(&G);
 
808
        dpd_buf4_init(&G, CC_TMP0, 0, 27, 29, 27, 29, 0, "GiCbA");
 
809
        dpd_contract442(&Fints, &G, &I, 0, 0, 2.0, 1.0);
 
810
        dpd_buf4_close(&G);
 
811
        dpd_buf4_close(&Fints);
 
812
 
 
813
        dpd_file2_close(&I);
 
814
 
 
815
      }
 
816
 
 
817
      if(params.ref == 0) { /** RHF **/
 
818
        /* I'IJ <-- sum_KAB <IK||AB> G(JK,AB) + 2 sum_kAb <Ik|Ab> G(Jk,Ab) */
 
819
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
820
 
 
821
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 5, 0, 5, 0, "GIjAb");
 
822
        dpd_buf4_scmcopy(&G, CC_GAMMA, "2 Gijab - Gijba", 2);
 
823
        dpd_buf4_sort_axpy(&G, CC_GAMMA, pqsr, 0, 5, "2 Gijab - Gijba", -1);
 
824
        dpd_buf4_close(&G);
 
825
 
 
826
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 5, 0, 5, 0, "2 Gijab - Gijba");
 
827
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 5, 0, 5, 0, "D <ij|ab>");
 
828
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
829
        dpd_buf4_close(&Dints);
 
830
        dpd_buf4_close(&G);
 
831
 
 
832
        dpd_file2_close(&I);
 
833
      }
 
834
      else if(params.ref == 1) { /** ROHF **/
 
835
 
 
836
        /* I'IJ <-- sum_KAB <IK||AB> G(JK,AB) + 2 sum_kAb <Ik|Ab> G(Jk,Ab) */
 
837
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
838
 
 
839
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 7, 0, 7, 0, "D <ij||ab> (ij,a>b)");
 
840
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 7, 2, 7, 0, "GIJAB");
 
841
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
842
        dpd_buf4_close(&G);
 
843
        dpd_buf4_close(&Dints);
 
844
 
 
845
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 5, 0, 5, 0, "D <ij|ab>");
 
846
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 5, 0, 5, 0, "GIjAb");
 
847
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
848
        dpd_buf4_close(&G);
 
849
        dpd_buf4_close(&Dints);
 
850
 
 
851
        dpd_file2_close(&I);
 
852
 
 
853
        /* I'ij <-- sum_kab <ik||ab> G(jk,ab) + 2 sum_KaB <iK|aB> G(jK,aB) */
 
854
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
855
 
 
856
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 7, 0, 7, 0, "D <ij||ab> (ij,a>b)");
 
857
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 7, 2, 7, 0, "Gijab");
 
858
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
859
        dpd_buf4_close(&G);
 
860
        dpd_buf4_close(&Dints);
 
861
 
 
862
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 5, 0, 5, 0, "D <ij|ab>");
 
863
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 5, 0, 5, 0, "GIjAb");
 
864
        dpd_contract442(&Dints, &G, &I, 1, 1, 2.0, 1.0);
 
865
        dpd_buf4_close(&G);
 
866
        dpd_buf4_close(&Dints);
 
867
 
 
868
        dpd_file2_close(&I);
 
869
 
 
870
      }
 
871
      else if(params.ref == 2) { /** UHF **/
 
872
 
 
873
        /* I'IJ <-- sum_KAB <IK||AB> G(JK,AB) + 2 sum_kAb <Ik|Ab> G(Jk,Ab) */
 
874
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
875
 
 
876
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 7, 0, 7, 0, "D <IJ||AB> (IJ,A>B)");
 
877
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 7, 2, 7, 0, "GIJAB");
 
878
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
879
        dpd_buf4_close(&G);
 
880
        dpd_buf4_close(&Dints);
 
881
 
 
882
        dpd_buf4_init(&Dints, CC_DINTS, 0, 22, 28, 22, 28, 0, "D <Ij|Ab>");
 
883
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 28, 22, 28, 0, "GIjAb");
 
884
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
885
        dpd_buf4_close(&G);
 
886
        dpd_buf4_close(&Dints);
 
887
 
 
888
        dpd_file2_close(&I);
 
889
 
 
890
        /* I'ij <-- sum_kab <ik||ab> G(jk,ab) + 2 sum_KaB <iK|aB> G(jK,aB) */
 
891
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
892
 
 
893
        dpd_buf4_init(&Dints, CC_DINTS, 0, 10, 17, 10, 17, 0, "D <ij||ab> (ij,a>b)");
 
894
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 17, 12, 17, 0, "Gijab");
 
895
        dpd_contract442(&Dints, &G, &I, 0, 0, 2.0, 1.0);
 
896
        dpd_buf4_close(&G);
 
897
        dpd_buf4_close(&Dints);
 
898
 
 
899
        dpd_buf4_init(&Dints, CC_DINTS, 0, 22, 28, 22, 28, 0, "D <Ij|Ab>");
 
900
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 28, 22, 28, 0, "GIjAb");
 
901
        dpd_contract442(&Dints, &G, &I, 1, 1, 2.0, 1.0);
 
902
        dpd_buf4_close(&G);
 
903
        dpd_buf4_close(&Dints);
 
904
 
 
905
        dpd_file2_close(&I);
 
906
 
 
907
      }
 
908
 
 
909
      if(params.ref == 0) { /** RHF **/
 
910
        /* I'IJ <-- 2 sum_AKB <IA||KB> G(JA,KB) + 2 sum_aKb <Ia|Kb> G(Ja,Kb) -
 
911
           2 sum_akB <Ik|Ba> GJakB(Jk,Ba) */
 
912
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
913
 
 
914
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia||jb>");
 
915
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIBJA");
 
916
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
917
        dpd_buf4_close(&G);
 
918
        dpd_buf4_close(&Cints);
 
919
 
 
920
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia|jb>");
 
921
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIbJa");
 
922
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
923
        dpd_buf4_close(&G);
 
924
        dpd_buf4_close(&Cints);
 
925
 
 
926
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 5, 0, 5, 0, "D <ij|ab>");
 
927
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIbjA");
 
928
        dpd_buf4_sort(&G, CC_TMP0, prsq, 0, 5, "GIbjA (Ij,Ab)");
 
929
        dpd_buf4_close(&G);
 
930
        dpd_buf4_init(&G, CC_TMP0, 0, 0, 5, 0, 5, 0, "GIbjA (Ij,Ab)");
 
931
        dpd_contract442(&Dints, &G, &I, 0, 0, -2.0, 1.0);
 
932
        dpd_buf4_close(&G);
 
933
        dpd_buf4_close(&Dints);
 
934
 
 
935
        dpd_file2_close(&I);
 
936
      }
 
937
      else if(params.ref == 1) { /** ROHF **/
 
938
 
 
939
        /* I'IJ <-- 2 sum_AKB <IA||KB> G(JA,KB) + 2 sum_aKb <Ia|Kb> G(Ja,Kb) -
 
940
           2 sum_akB <Ik|Ba> GJakB(Jk,Ba) */
 
941
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
942
 
 
943
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia||jb>");
 
944
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIBJA");
 
945
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
946
        dpd_buf4_close(&G);
 
947
        dpd_buf4_close(&Cints);
 
948
 
 
949
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia|jb>");
 
950
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIbJa");
 
951
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
952
        dpd_buf4_close(&G);
 
953
        dpd_buf4_close(&Cints);
 
954
 
 
955
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 5, 0, 5, 0, "D <ij|ab>");
 
956
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIbjA");
 
957
        dpd_buf4_sort(&G, CC_TMP0, prsq, 0, 5, "GIbjA (Ij,Ab)");
 
958
        dpd_buf4_close(&G);
 
959
        dpd_buf4_init(&G, CC_TMP0, 0, 0, 5, 0, 5, 0, "GIbjA (Ij,Ab)");
 
960
        dpd_contract442(&Dints, &G, &I, 0, 0, -2.0, 1.0);
 
961
        dpd_buf4_close(&G);
 
962
        dpd_buf4_close(&Dints);
 
963
 
 
964
        dpd_file2_close(&I);
 
965
 
 
966
        /* I'ij <-- 2 sum_akb <ia||kb> G(ja,kb) + 2 sum_AkB <iA|kB> G(jA,kB) +
 
967
           2 sum_AKb <iK|bA> GjAKb(jK,bA) */
 
968
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
969
 
 
970
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia||jb>");
 
971
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "Gibja");
 
972
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
973
        dpd_buf4_close(&G);
 
974
        dpd_buf4_close(&Cints);
 
975
 
 
976
        dpd_buf4_init(&Cints, CC_CINTS, 0, 10, 10, 10, 10, 0, "C <ia|jb>");
 
977
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GiBjA");
 
978
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
979
        dpd_buf4_close(&G);
 
980
        dpd_buf4_close(&Cints);
 
981
 
 
982
        dpd_buf4_init(&Dints, CC_DINTS, 0, 0, 5, 0, 5, 0, "D <ij|ab>");
 
983
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GiBJa");
 
984
        dpd_buf4_sort(&G, CC_TMP0, prsq, 0, 5, "GiBJa (iJ,aB)");
 
985
        dpd_buf4_close(&G);
 
986
        dpd_buf4_init(&G, CC_TMP0, 0, 0, 5, 0, 5, 0, "GiBJa (iJ,aB)");
 
987
        dpd_contract442(&Dints, &G, &I, 0, 0, -2.0, 1.0);
 
988
        dpd_buf4_close(&G);
 
989
        dpd_buf4_close(&Dints);
 
990
 
 
991
        dpd_file2_close(&I);
 
992
 
 
993
      }
 
994
      else if(params.ref == 2) { /** UHF **/
 
995
 
 
996
        /* I'IJ <-- 2 sum_AKB <IA||KB> G(JA,KB) + 2 sum_aKb <Ia|Kb> G(Ja,Kb) -
 
997
           2 sum_akB <Ik|Ba> GJakB(Jk,Ba) */
 
998
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
999
 
 
1000
        dpd_buf4_init(&Cints, CC_CINTS, 0, 20, 20, 20, 20, 0, "C <IA||JB>");
 
1001
        dpd_buf4_init(&G, CC_GAMMA, 0, 20, 20, 20, 20, 0, "GIBJA");
 
1002
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
1003
        dpd_buf4_close(&G);
 
1004
        dpd_buf4_close(&Cints);
 
1005
 
 
1006
        dpd_buf4_init(&Cints, CC_CINTS, 0, 24, 24, 24, 24, 0, "C <Ia|Jb>");
 
1007
        dpd_buf4_init(&G, CC_GAMMA, 0, 24, 24, 24, 24, 0, "GIbJa");
 
1008
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
1009
        dpd_buf4_close(&G);
 
1010
        dpd_buf4_close(&Cints);
 
1011
 
 
1012
        dpd_buf4_init(&Dints, CC_DINTS, 0, 22, 28, 22, 28, 0, "D <Ij|Ab>");
 
1013
        dpd_buf4_init(&G, CC_GAMMA, 0, 24, 27, 24, 27, 0, "GIbjA");
 
1014
        dpd_buf4_sort(&G, CC_TMP0, prsq, 22, 28, "GIbjA (Ij,Ab)");
 
1015
        dpd_buf4_close(&G);
 
1016
        dpd_buf4_init(&G, CC_TMP0, 0, 22, 28, 22, 28, 0, "GIbjA (Ij,Ab)");
 
1017
        dpd_contract442(&Dints, &G, &I, 0, 0, -2.0, 1.0);
 
1018
        dpd_buf4_close(&G);
 
1019
        dpd_buf4_close(&Dints);
 
1020
 
 
1021
        dpd_file2_close(&I);
 
1022
 
 
1023
        /* I'ij <-- 2 sum_akb <ia||kb> G(ja,kb) + 2 sum_AkB <iA|kB> G(jA,kB) +
 
1024
           2 sum_AKb <iK|bA> GjAKb(jK,bA) */
 
1025
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
1026
 
 
1027
        dpd_buf4_init(&Cints, CC_CINTS, 0, 30, 30, 30, 30, 0, "C <ia||jb>");
 
1028
        dpd_buf4_init(&G, CC_GAMMA, 0, 30, 30, 30, 30, 0, "Gibja");
 
1029
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
1030
        dpd_buf4_close(&G);
 
1031
        dpd_buf4_close(&Cints);
 
1032
 
 
1033
        dpd_buf4_init(&Cints, CC_CINTS, 0, 27, 27, 27, 27, 0, "C <iA|jB>");
 
1034
        dpd_buf4_init(&G, CC_GAMMA, 0, 27, 27, 27, 27, 0, "GiBjA");
 
1035
        dpd_contract442(&Cints, &G, &I, 0, 0, 2.0, 1.0);
 
1036
        dpd_buf4_close(&G);
 
1037
        dpd_buf4_close(&Cints);
 
1038
 
 
1039
        dpd_buf4_init(&Dints, CC_DINTS, 0, 23, 29, 23, 29, 0, "D <iJ|aB>");
 
1040
        dpd_buf4_init(&G, CC_GAMMA, 0, 27, 24, 27, 24, 0, "GiBJa");
 
1041
        dpd_buf4_sort(&G, CC_TMP0, prsq, 23, 29, "GiBJa (iJ,aB)");
 
1042
        dpd_buf4_close(&G);
 
1043
        dpd_buf4_init(&G, CC_TMP0, 0, 23, 29, 23, 29, 0, "GiBJa (iJ,aB)");
 
1044
        dpd_contract442(&Dints, &G, &I, 0, 0, -2.0, 1.0);
 
1045
        dpd_buf4_close(&G);
 
1046
        dpd_buf4_close(&Dints);
 
1047
 
 
1048
        dpd_file2_close(&I);
 
1049
 
 
1050
      }
 
1051
 
 
1052
      if(params.ref == 0) { /** RHF **/
 
1053
        /* I'IJ <-- 2 sum_KLA <IK||LA> G(JK,LA) + 2 sum_kLa <Ik|La> G(Jk,La)
 
1054
           + 2 sum_kAl <kI|lA> G(kJ,lA) */
 
1055
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
1056
 
 
1057
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GIjKa");
 
1058
        dpd_buf4_scmcopy(&G, CC_GAMMA, "2 Gijka - Gjika", 2);
 
1059
        dpd_buf4_sort_axpy(&G, CC_GAMMA, qprs, 0, 10, "2 Gijka - Gjika", -1);
 
1060
        dpd_buf4_close(&G);
 
1061
 
 
1062
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "2 Gijka - Gjika");
 
1063
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1064
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1065
        dpd_contract442(&Eints, &G, &I, 1, 1, 2.0, 1.0);
 
1066
        dpd_buf4_close(&Eints);
 
1067
        dpd_buf4_close(&G);
 
1068
 
 
1069
        dpd_file2_close(&I);
 
1070
      }
 
1071
      else if(params.ref == 1) { /** ROHF **/
 
1072
 
 
1073
        /* I'IJ <-- 2 sum_KLA <IK||LA> G(JK,LA) + 2 sum_kLa <Ik|La> G(Jk,La)
 
1074
           + 2 sum_kAl <kI|lA> G(kJ,lA) */
 
1075
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
1076
 
 
1077
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
1078
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 2, 10, 0, "GIJKA");
 
1079
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1080
        dpd_buf4_close(&G);
 
1081
        dpd_buf4_close(&Eints);
 
1082
 
 
1083
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1084
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GIjKa");
 
1085
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1086
        dpd_buf4_close(&G);
 
1087
        dpd_buf4_close(&Eints);
 
1088
 
 
1089
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1090
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GiJkA");
 
1091
        dpd_contract442(&Eints, &G, &I, 1, 1, 2.0, 1.0);
 
1092
        dpd_buf4_close(&G);
 
1093
        dpd_buf4_close(&Eints);
 
1094
 
 
1095
        dpd_file2_close(&I);
 
1096
 
 
1097
        /* I'ij <-- 2 sum_kla <ik||la> G(jk,la) + 2 sum_KlA <iK|lA> G(jK,lA)
 
1098
           + 2 sum_KaL <Ki|La> G(Kj,La) */
 
1099
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
1100
 
 
1101
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
1102
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 2, 10, 0, "Gijka");
 
1103
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1104
        dpd_buf4_close(&G);
 
1105
        dpd_buf4_close(&Eints);
 
1106
 
 
1107
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1108
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GiJkA");
 
1109
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1110
        dpd_buf4_close(&G);
 
1111
        dpd_buf4_close(&Eints);
 
1112
 
 
1113
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1114
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GIjKa");
 
1115
        dpd_contract442(&Eints, &G, &I, 1, 1, 2.0, 1.0);
 
1116
        dpd_buf4_close(&G);
 
1117
        dpd_buf4_close(&Eints);
 
1118
 
 
1119
        dpd_file2_close(&I);
 
1120
      }
 
1121
      else if(params.ref == 2) { /** UHF **/
 
1122
 
 
1123
        /* I'IJ <-- 2 sum_KLA <IK||LA> G(JK,LA) + 2 sum_kLa <Ik|La> G(Jk,La)
 
1124
           + 2 sum_kAl <kI|lA> G(kJ,lA) */
 
1125
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
1126
 
 
1127
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 20, 2, 20, 0, "E <IJ||KA> (I>J,KA)");
 
1128
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 20, 2, 20, 0, "GIJKA");
 
1129
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1130
        dpd_buf4_close(&G);
 
1131
        dpd_buf4_close(&Eints);
 
1132
 
 
1133
        dpd_buf4_init(&Eints, CC_EINTS, 0, 22, 24, 22, 24, 0, "E <Ij|Ka>");
 
1134
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 24, 22, 24, 0, "GIjKa");
 
1135
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1136
        dpd_buf4_close(&G);
 
1137
        dpd_buf4_close(&Eints);
 
1138
 
 
1139
        dpd_buf4_init(&Eints, CC_EINTS, 0, 23, 27, 23, 27, 0, "E <iJ|kA>");
 
1140
        dpd_buf4_init(&G, CC_GAMMA, 0, 23, 27, 23, 27, 0, "GiJkA");
 
1141
        dpd_contract442(&Eints, &G, &I, 1, 1, 2.0, 1.0);
 
1142
        dpd_buf4_close(&G);
 
1143
        dpd_buf4_close(&Eints);
 
1144
 
 
1145
        dpd_file2_close(&I);
 
1146
 
 
1147
        /* I'ij <-- 2 sum_kla <ik||la> G(jk,la) + 2 sum_KlA <iK|lA> G(jK,lA)
 
1148
           + 2 sum_KaL <Ki|La> G(Kj,La) */
 
1149
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
1150
 
 
1151
        dpd_buf4_init(&Eints, CC_EINTS, 0, 10, 30, 12, 30, 0, "E <ij||ka> (i>j,ka)");
 
1152
        dpd_buf4_init(&G, CC_GAMMA, 0, 10, 30, 12, 30, 0, "Gijka");
 
1153
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1154
        dpd_buf4_close(&G);
 
1155
        dpd_buf4_close(&Eints);
 
1156
 
 
1157
        dpd_buf4_init(&Eints, CC_EINTS, 0, 23, 27, 23, 27, 0, "E <iJ|kA>");
 
1158
        dpd_buf4_init(&G, CC_GAMMA, 0, 23, 27, 23, 27, 0, "GiJkA");
 
1159
        dpd_contract442(&Eints, &G, &I, 0, 0, 2.0, 1.0);
 
1160
        dpd_buf4_close(&G);
 
1161
        dpd_buf4_close(&Eints);
 
1162
 
 
1163
        dpd_buf4_init(&Eints, CC_EINTS, 0, 22, 24, 22, 24, 0, "E <Ij|Ka>");
 
1164
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 24, 22, 24, 0, "GIjKa");
 
1165
        dpd_contract442(&Eints, &G, &I, 1, 1, 2.0, 1.0);
 
1166
        dpd_buf4_close(&G);
 
1167
        dpd_buf4_close(&Eints);
 
1168
 
 
1169
        dpd_file2_close(&I);
 
1170
 
 
1171
      }
 
1172
 
 
1173
      if(params.ref == 0) { /** RHF **/
 
1174
        /* I'IJ <-- 2 sum_AKL <K>L||IA> G(K>L,JA) + 2 sum_aKl <Kl|Ia> G(Kl,Ja) */
 
1175
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
1176
 
 
1177
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GIjKa");
 
1178
        dpd_buf4_scmcopy(&G, CC_GAMMA, "2 Gijka - Gjika", 2);
 
1179
        dpd_buf4_sort_axpy(&G, CC_GAMMA, qprs, 0, 10, "2 Gijka - Gjika", -1);
 
1180
        dpd_buf4_close(&G);
 
1181
 
 
1182
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "2 Gijka - Gjika");
 
1183
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1184
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1185
        dpd_buf4_close(&Eints);
 
1186
        dpd_buf4_close(&G);
 
1187
 
 
1188
        dpd_file2_close(&I);
 
1189
      }
 
1190
      else if(params.ref == 1) { /** ROHF **/
 
1191
 
 
1192
        /* I'IJ <-- 2 sum_AKL <K>L||IA> G(K>L,JA) + 2 sum_aKl <Kl|Ia> G(Kl,Ja) */
 
1193
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
1194
 
 
1195
        dpd_buf4_init(&Eints, CC_EINTS, 0, 2, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
1196
        dpd_buf4_init(&G, CC_GAMMA, 0, 2, 10, 2, 10, 0, "GIJKA");
 
1197
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1198
        dpd_buf4_close(&G);
 
1199
        dpd_buf4_close(&Eints);
 
1200
 
 
1201
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1202
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GIjKa");
 
1203
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1204
        dpd_buf4_close(&G);
 
1205
        dpd_buf4_close(&Eints);
 
1206
 
 
1207
        dpd_file2_close(&I);
 
1208
 
 
1209
        /* I'ij <-- 2 sum_akl <k>l||ia> G(k>l,ja) + 2 sum_AkL <kL|iA> G(kL,jA) */
 
1210
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'ij");
 
1211
 
 
1212
        dpd_buf4_init(&Eints, CC_EINTS, 0, 2, 10, 2, 10, 0, "E <ij||ka> (i>j,ka)");
 
1213
        dpd_buf4_init(&G, CC_GAMMA, 0, 2, 10, 2, 10, 0, "Gijka");
 
1214
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1215
        dpd_buf4_close(&G);
 
1216
        dpd_buf4_close(&Eints);
 
1217
 
 
1218
        dpd_buf4_init(&Eints, CC_EINTS, 0, 0, 10, 0, 10, 0, "E <ij|ka>");
 
1219
        dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GiJkA");
 
1220
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1221
        dpd_buf4_close(&G);
 
1222
        dpd_buf4_close(&Eints);
 
1223
 
 
1224
        dpd_file2_close(&I);
 
1225
 
 
1226
      }
 
1227
      else if(params.ref == 2) { /** UHF **/
 
1228
 
 
1229
        /* I'IJ <-- 2 sum_AKL <K>L||IA> G(K>L,JA) + 2 sum_aKl <Kl|Ia> G(Kl,Ja) */
 
1230
        dpd_file2_init(&I, CC_OEI, 0, 0, 0, "I'IJ");
 
1231
 
 
1232
        dpd_buf4_init(&Eints, CC_EINTS, 0, 2, 20, 2, 20, 0, "E <IJ||KA> (I>J,KA)");
 
1233
        dpd_buf4_init(&G, CC_GAMMA, 0, 2, 20, 2, 20, 0, "GIJKA");
 
1234
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1235
        dpd_buf4_close(&G);
 
1236
        dpd_buf4_close(&Eints);
 
1237
 
 
1238
        dpd_buf4_init(&Eints, CC_EINTS, 0, 22, 24, 22, 24, 0, "E <Ij|Ka>");
 
1239
        dpd_buf4_init(&G, CC_GAMMA, 0, 22, 24, 22, 24, 0, "GIjKa");
 
1240
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1241
        dpd_buf4_close(&G);
 
1242
        dpd_buf4_close(&Eints);
 
1243
 
 
1244
        dpd_file2_close(&I);
 
1245
 
 
1246
        /* I'ij <-- 2 sum_akl <k>l||ia> G(k>l,ja) + 2 sum_AkL <kL|iA> G(kL,jA) */
 
1247
        dpd_file2_init(&I, CC_OEI, 0, 2, 2, "I'ij");
 
1248
 
 
1249
        dpd_buf4_init(&Eints, CC_EINTS, 0, 12, 30, 12, 30, 0, "E <ij||ka> (i>j,ka)");
 
1250
        dpd_buf4_init(&G, CC_GAMMA, 0, 12, 30, 12, 30, 0, "Gijka");
 
1251
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1252
        dpd_buf4_close(&G);
 
1253
        dpd_buf4_close(&Eints);
 
1254
 
 
1255
        dpd_buf4_init(&Eints, CC_EINTS, 0, 23, 27, 23, 27, 0, "E <iJ|kA>");
 
1256
        dpd_buf4_init(&G, CC_GAMMA, 0, 23, 27, 23, 27, 0, "GiJkA");
 
1257
        dpd_contract442(&Eints, &G, &I, 2, 2, 2.0, 1.0);
 
1258
        dpd_buf4_close(&G);
 
1259
        dpd_buf4_close(&Eints);
 
1260
 
 
1261
        dpd_file2_close(&I);
 
1262
 
 
1263
      }
 
1264
    }
 
1265
 
 
1266
  }} // namespace psi::ccdensity