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

« back to all changes in this revision

Viewing changes to src/bin/ccenergy/tau.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:
 
1
/*! \file
 
2
    \ingroup CCENERGY
 
3
    \brief Enter brief description of file here 
 
4
*/
 
5
#include <cstdio>
 
6
#include <cstdlib>
 
7
#include <libdpd/dpd.h>
 
8
#include "MOInfo.h"
 
9
#include "Params.h"
 
10
#define EXTERN
 
11
#include "globals.h"
 
12
 
 
13
namespace psi { namespace ccenergy {
 
14
 
 
15
void tau_build(void)
 
16
{
 
17
  int h, ij, ab, i, j, a, b, I, J, A, B;
 
18
  int Isym, Jsym, Asym, Bsym;
 
19
  int nirreps;
 
20
  dpdbuf4 tauIJAB, tauijab, tauIjAb, tauiJaB, tauIjbA;
 
21
  dpdbuf4 tIJAB, tijab, tIjAb;
 
22
  dpdfile2 tIA, tia;
 
23
 
 
24
  nirreps = moinfo.nirreps;
 
25
 
 
26
  if(params.ref == 0) { /** RHF **/
 
27
 
 
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);
 
31
 
 
32
    dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA");
 
33
    dpd_file2_mat_init(&tIA);
 
34
    dpd_file2_mat_rd(&tIA);
 
35
 
 
36
    dpd_buf4_init(&tauIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tauIjAb");
 
37
 
 
38
    for(h=0; h < nirreps; h++) {
 
39
 
 
40
      dpd_buf4_mat_irrep_init(&tauIjAb, h);
 
41
      dpd_buf4_mat_irrep_rd(&tauIjAb, h);
 
42
 
 
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];
 
57
 
 
58
          if((Isym==Asym) && (Jsym==Bsym))
 
59
            tauIjAb.matrix[h][ij][ab] +=
 
60
              (tIA.matrix[Isym][I][A] * tIA.matrix[Jsym][J][B]);
 
61
 
 
62
        }
 
63
      }
 
64
 
 
65
      dpd_buf4_mat_irrep_wrt(&tauIjAb, h);
 
66
      dpd_buf4_mat_irrep_close(&tauIjAb, h);
 
67
    }
 
68
 
 
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);
 
72
 
 
73
    dpd_file2_mat_close(&tIA);
 
74
    dpd_file2_close(&tIA);
 
75
  }
 
76
  else if(params.ref == 1) { /** ROHF **/
 
77
 
 
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);
 
81
 
 
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);
 
85
 
 
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);
 
89
 
 
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);
 
96
 
 
97
    dpd_buf4_init(&tauIJAB, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tauIJAB");
 
98
 
 
99
    for(h=0; h < nirreps; h++) {
 
100
 
 
101
      dpd_buf4_mat_irrep_init(&tauIJAB, h);
 
102
      dpd_buf4_mat_irrep_rd(&tauIJAB, h);
 
103
 
 
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];
 
118
 
 
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]);
 
125
 
 
126
        }
 
127
      }
 
128
 
 
129
      dpd_buf4_mat_irrep_wrt(&tauIJAB, h);
 
130
      dpd_buf4_mat_irrep_close(&tauIJAB, h);
 
131
    }
 
132
 
 
133
    dpd_buf4_close(&tauIJAB);
 
134
 
 
135
    dpd_buf4_init(&tauijab, CC_TAMPS, 0, 2, 7, 2, 7, 0, "tauijab");
 
136
 
 
137
    for(h=0; h < nirreps; h++) {
 
138
 
 
139
      dpd_buf4_mat_irrep_init(&tauijab, h);
 
140
      dpd_buf4_mat_irrep_rd(&tauijab, h);
 
141
 
 
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];
 
156
 
 
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]);
 
163
 
 
164
        }
 
165
      }
 
166
 
 
167
      dpd_buf4_mat_irrep_wrt(&tauijab, h);
 
168
      dpd_buf4_mat_irrep_close(&tauijab, h);
 
169
    }
 
170
 
 
171
    dpd_buf4_close(&tauijab);
 
172
 
 
173
    dpd_buf4_init(&tauIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "tauIjAb");
 
174
 
 
175
    for(h=0; h < nirreps; h++) {
 
176
 
 
177
      dpd_buf4_mat_irrep_init(&tauIjAb, h);
 
178
      dpd_buf4_mat_irrep_rd(&tauIjAb, h);
 
179
 
 
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];
 
194
 
 
195
          if((Isym==Asym) && (Jsym==Bsym))
 
196
            tauIjAb.matrix[h][ij][ab] +=
 
197
              (tIA.matrix[Isym][I][A] * tia.matrix[Jsym][J][B]);
 
198
 
 
199
        }
 
200
      }
 
201
 
 
202
      dpd_buf4_mat_irrep_wrt(&tauIjAb, h);
 
203
      dpd_buf4_mat_irrep_close(&tauIjAb, h);
 
204
    }
 
205
 
 
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);
 
212
 
 
213
    dpd_file2_mat_close(&tIA);
 
214
    dpd_file2_close(&tIA);
 
215
    dpd_file2_mat_close(&tia);
 
216
    dpd_file2_close(&tia);
 
217
 
 
218
  }
 
219
  else if(params.ref == 2) { /*** UHF ***/
 
220
 
 
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);
 
224
 
 
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);
 
228
 
 
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);
 
232
 
 
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);
 
239
 
 
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]);
 
264
        }
 
265
      }
 
266
      dpd_buf4_mat_irrep_wrt(&tauIJAB, h);
 
267
      dpd_buf4_mat_irrep_close(&tauIJAB, h);
 
268
    }
 
269
    dpd_buf4_close(&tauIJAB);
 
270
 
 
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]);
 
295
        }
 
296
      }
 
297
      dpd_buf4_mat_irrep_wrt(&tauijab, h);
 
298
      dpd_buf4_mat_irrep_close(&tauijab, h);
 
299
    }
 
300
    dpd_buf4_close(&tauijab);
 
301
 
 
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]);
 
323
        }
 
324
      }
 
325
      dpd_buf4_mat_irrep_wrt(&tauIjAb, h);
 
326
      dpd_buf4_mat_irrep_close(&tauIjAb, h);
 
327
    }
 
328
    dpd_buf4_close(&tauIjAb);
 
329
 
 
330
    dpd_file2_mat_close(&tIA);
 
331
    dpd_file2_close(&tIA);
 
332
    dpd_file2_mat_close(&tia);
 
333
    dpd_file2_close(&tia);
 
334
 
 
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);
 
342
 
 
343
  } /*** UHF ***/
 
344
}
 
345
}} // namespace psi::ccenergy