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

« back to all changes in this revision

Viewing changes to src/bin/ccenergy/FT2.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 <libpsio/psio.h>
 
9
#include <libqt/qt.h>
 
10
#include "Params.h"
 
11
#include "MOInfo.h"
 
12
#define EXTERN
 
13
#include "globals.h"
 
14
 
 
15
namespace psi { namespace ccenergy {
 
16
 
 
17
void FT2(void)
 
18
{
 
19
  dpdfile2 tIA, tia, t1;
 
20
  dpdbuf4 newtIJAB, newtijab, newtIjAb, t2, t2a, t2b;
 
21
  dpdbuf4 F_anti, F;
 
22
  dpdbuf4 Z, X;
 
23
  int Gie, Gij, Gab, nrows, ncols, nlinks, Gi, Ge, Gj, i, I;
 
24
 
 
25
  if(params.ref == 0) { /** RHF **/
 
26
 
 
27
/*     dpd_buf4_init(&newtIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "New tIjAb"); */
 
28
 
 
29
/*     dpd_buf4_init(&Z, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z1(ij,ab)"); */
 
30
/*     dpd_buf4_init(&F, CC_FINTS, 0, 11, 5, 11, 5, 0, "F <ai|bc>"); */
 
31
/*     dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA"); */
 
32
/*     dpd_contract244(&tIA, &F, &Z, 1, 0, 0, 1, 0); */
 
33
/*     dpd_file2_close(&tIA);  */
 
34
/*     dpd_buf4_close(&F); */
 
35
 
 
36
/*     dpd_buf4_sort(&Z, CC_TMP0, qpsr, 0, 5, "Z2(ji,ba)"); */
 
37
/*     dpd_buf4_axpy(&Z, &newtIjAb, 1.0); */
 
38
/*     dpd_buf4_close(&Z); */
 
39
/*     dpd_buf4_init(&Z, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z2(ji,ba)"); */
 
40
/*     dpd_buf4_axpy(&Z, &newtIjAb, 1.0); */
 
41
/*     dpd_buf4_close(&Z); */
 
42
 
 
43
/*     dpd_buf4_close(&newtIjAb); */
 
44
 
 
45
    /* t(ij,ab) <-- t(j,e) * <ie|ab> + t(i,e) * <je|ba> */
 
46
    /* OOC code added 3/23/05, TDC */
 
47
    dpd_buf4_init(&X, CC_TMP0, 0, 0, 5, 0, 5, 0, "X(Ij,Ab)");
 
48
    dpd_buf4_init(&F, CC_FINTS, 0, 10, 5, 10, 5, 0, "F <ia|bc>");
 
49
    dpd_file2_init(&t1, CC_OEI, 0, 0, 1, "tIA");
 
50
    dpd_file2_mat_init(&t1);
 
51
    dpd_file2_mat_rd(&t1);
 
52
    for(Gie=0; Gie < moinfo.nirreps; Gie++) { 
 
53
      Gab = Gie; /* F is totally symmetric */
 
54
      Gij = Gab; /* T2 is totally symmetric */
 
55
      dpd_buf4_mat_irrep_init(&X, Gij);
 
56
      ncols = F.params->coltot[Gie];
 
57
 
 
58
      for(Gi=0; Gi < moinfo.nirreps; Gi++) {
 
59
        Gj = Ge = Gi^Gie; /* T1 is totally symmetric */
 
60
 
 
61
        nlinks = moinfo.virtpi[Ge];
 
62
        nrows = moinfo.occpi[Gj];
 
63
 
 
64
        dpd_buf4_mat_irrep_init_block(&F, Gie, nlinks);
 
65
 
 
66
        for(i=0; i < moinfo.occpi[Gi]; i++) {
 
67
          I = F.params->poff[Gi] + i;
 
68
          dpd_buf4_mat_irrep_rd_block(&F, Gie, F.row_offset[Gie][I], nlinks);
 
69
 
 
70
          if(nrows && ncols && nlinks)
 
71
            C_DGEMM('n','n',nrows,ncols,nlinks,1.0,t1.matrix[Gj][0],nlinks,F.matrix[Gie][0],ncols,
 
72
                    0.0,X.matrix[Gij][X.row_offset[Gij][I]],ncols);
 
73
        }
 
74
 
 
75
        dpd_buf4_mat_irrep_close_block(&F, Gie, nlinks);
 
76
      }
 
77
 
 
78
      dpd_buf4_mat_irrep_wrt(&X, Gij);
 
79
      dpd_buf4_mat_irrep_close(&X, Gij);
 
80
    }
 
81
    dpd_file2_mat_close(&t1);
 
82
    dpd_file2_close(&t1);
 
83
    dpd_buf4_close(&F);
 
84
    dpd_buf4_init(&t2, CC_TAMPS, 0, 0, 5, 0, 5, 0, "New tIjAb");
 
85
    dpd_buf4_axpy(&X, &t2, 1);
 
86
    dpd_buf4_close(&t2);
 
87
    dpd_buf4_sort_axpy(&X, CC_TAMPS, qpsr, 0, 5, "New tIjAb", 1);
 
88
    dpd_buf4_close(&X);
 
89
  }
 
90
  else if(params.ref == 1) { /** ROHF **/
 
91
 
 
92
    dpd_buf4_init(&newtIJAB, CC_TAMPS, 0, 0, 7, 2, 7, 0, "New tIJAB");
 
93
    dpd_buf4_init(&newtijab, CC_TAMPS, 0, 0, 7, 2, 7, 0, "New tijab");
 
94
    dpd_buf4_init(&newtIjAb, CC_TAMPS, 0, 0, 5, 0, 5, 0, "New tIjAb");
 
95
 
 
96
    dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA");
 
97
    dpd_file2_init(&tia, CC_OEI, 0, 0, 1, "tia");
 
98
 
 
99
    /*** AA ***/
 
100
 
 
101
    dpd_buf4_init(&F_anti, CC_FINTS, 0, 10, 7, 10, 5, 1, "F <ia|bc>");
 
102
    dpd_buf4_init(&t2, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (IJ,A>B)");
 
103
    dpd_contract424(&F_anti, &tIA, &t2, 1, 1, 1, 1, 0);
 
104
    dpd_buf4_sort(&t2, CC_TMP0, qprs, 0, 7, "T (JI,A>B)");
 
105
    dpd_buf4_close(&t2);
 
106
    dpd_buf4_init(&t2a, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (IJ,A>B)");
 
107
    dpd_buf4_init(&t2b, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (JI,A>B)");
 
108
    dpd_buf4_axpy(&t2b, &t2a, -1);
 
109
    dpd_buf4_axpy(&t2a, &newtIJAB, 1);
 
110
    dpd_buf4_close(&t2b);
 
111
    dpd_buf4_close(&t2a);
 
112
    dpd_buf4_close(&F_anti);
 
113
 
 
114
    /*** BB ***/
 
115
 
 
116
    dpd_buf4_init(&F_anti, CC_FINTS, 0, 10, 7, 10, 5, 1, "F <ia|bc>");
 
117
    dpd_buf4_init(&t2, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (IJ,A>B)");
 
118
    dpd_contract424(&F_anti, &tia, &t2, 1, 1, 1, 1, 0);
 
119
    dpd_buf4_sort(&t2, CC_TMP0, qprs, 0, 7, "T (JI,A>B)");
 
120
    dpd_buf4_close(&t2);
 
121
    dpd_buf4_init(&t2a, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (IJ,A>B)");
 
122
    dpd_buf4_init(&t2b, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (JI,A>B)");
 
123
    dpd_buf4_axpy(&t2b, &t2a, -1);
 
124
    dpd_buf4_axpy(&t2a, &newtijab, 1);
 
125
    dpd_buf4_close(&t2b);
 
126
    dpd_buf4_close(&t2a);
 
127
    dpd_buf4_close(&F_anti);
 
128
 
 
129
    /*** AB ***/
 
130
 
 
131
    dpd_buf4_init(&F, CC_FINTS, 0, 10, 5, 10, 5, 0, "F <ia|bc>");
 
132
    dpd_contract424(&F, &tia, &newtIjAb, 1, 1, 1, 1, 1);
 
133
    dpd_buf4_close(&F);
 
134
    dpd_buf4_init(&F, CC_FINTS, 0, 11, 5, 11, 5, 0, "F <ai|bc>");
 
135
    dpd_contract244(&tIA, &F, &newtIjAb, 1, 0, 0, 1, 1);
 
136
    dpd_buf4_close(&F);
 
137
 
 
138
    dpd_file2_close(&tIA); 
 
139
    dpd_file2_close(&tia);
 
140
 
 
141
    dpd_buf4_close(&newtIJAB);
 
142
    dpd_buf4_close(&newtijab);
 
143
    dpd_buf4_close(&newtIjAb);
 
144
  }
 
145
  else if(params.ref == 2) { /*** UHF ***/
 
146
 
 
147
    dpd_buf4_init(&newtIJAB, CC_TAMPS, 0, 0, 7, 2, 7, 0, "New tIJAB");
 
148
    dpd_buf4_init(&newtijab, CC_TAMPS, 0, 10, 17, 12, 17, 0, "New tijab");
 
149
    dpd_buf4_init(&newtIjAb, CC_TAMPS, 0, 22, 28, 22, 28, 0, "New tIjAb");
 
150
 
 
151
    dpd_file2_init(&tIA, CC_OEI, 0, 0, 1, "tIA");
 
152
    dpd_file2_init(&tia, CC_OEI, 0, 2, 3, "tia");
 
153
 
 
154
    /*** AA ***/
 
155
 
 
156
    dpd_buf4_init(&F, CC_FINTS, 0, 20, 7, 20, 5, 1, "F <IA|BC>");
 
157
    dpd_buf4_init(&t2, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (IJ,A>B)");
 
158
    dpd_contract424(&F, &tIA, &t2, 1, 1, 1, 1, 0);
 
159
    dpd_buf4_sort(&t2, CC_TMP0, qprs, 0, 7, "T (JI,A>B)");
 
160
    dpd_buf4_close(&t2);
 
161
    dpd_buf4_init(&t2a, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (IJ,A>B)");
 
162
    dpd_buf4_init(&t2b, CC_TMP0, 0, 0, 7, 0, 7, 0, "T (JI,A>B)");
 
163
    dpd_buf4_axpy(&t2b, &t2a, -1);
 
164
    dpd_buf4_axpy(&t2a, &newtIJAB, 1);
 
165
    dpd_buf4_close(&t2b);
 
166
    dpd_buf4_close(&t2a);
 
167
    dpd_buf4_close(&F);
 
168
 
 
169
    /*** BB ***/
 
170
 
 
171
    dpd_buf4_init(&F, CC_FINTS, 0, 30, 17, 30, 15, 1, "F <ia|bc>");
 
172
    dpd_buf4_init(&t2, CC_TMP0, 0, 10, 17, 10, 17, 0, "T (ij,a>b)");
 
173
    dpd_contract424(&F, &tia, &t2, 1, 1, 1, 1, 0);
 
174
    dpd_buf4_sort(&t2, CC_TMP0, qprs, 10, 17, "T (ji,a>b)");
 
175
    dpd_buf4_close(&t2);
 
176
    dpd_buf4_init(&t2a, CC_TMP0, 0, 10, 17, 10, 17, 0, "T (ij,a>b)");
 
177
    dpd_buf4_init(&t2b, CC_TMP0, 0, 10, 17, 10, 17, 0, "T (ji,a>b)");
 
178
    dpd_buf4_axpy(&t2b, &t2a, -1);
 
179
    dpd_buf4_axpy(&t2a, &newtijab, 1);
 
180
    dpd_buf4_close(&t2b);
 
181
    dpd_buf4_close(&t2a);
 
182
    dpd_buf4_close(&F);
 
183
 
 
184
    /*** AB ***/
 
185
 
 
186
    dpd_buf4_init(&F, CC_FINTS, 0, 24, 28, 24, 28, 0, "F <Ia|Bc>");
 
187
    dpd_contract424(&F, &tia, &newtIjAb, 1, 1, 1, 1, 1);
 
188
    dpd_buf4_close(&F);
 
189
    dpd_buf4_init(&F, CC_FINTS, 0, 28, 26, 28, 26, 0, "F <Ab|Ci>");
 
190
    dpd_contract244(&tIA, &F, &newtIjAb, 1, 2, 0, 1, 1);
 
191
    dpd_buf4_close(&F);
 
192
 
 
193
    dpd_file2_close(&tIA); 
 
194
    dpd_file2_close(&tia);
 
195
 
 
196
    dpd_buf4_close(&newtIJAB);
 
197
    dpd_buf4_close(&newtijab);
 
198
    dpd_buf4_close(&newtIjAb);
 
199
 
 
200
  }
 
201
 
 
202
}
 
203
}} // namespace psi::ccenergy