~ubuntu-branches/ubuntu/quantal/psicode/quantal

« back to all changes in this revision

Viewing changes to src/bin/cis/denom.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 CIS
 
3
    \brief Enter brief description of file here 
 
4
*/
 
5
#include <libdpd/dpd.h>
 
6
#include "MOInfo.h"
 
7
#include "Params.h"
 
8
#include "Local.h"
 
9
#define EXTERN
 
10
#include "globals.h"
 
11
 
 
12
namespace psi { namespace cis {
 
13
 
 
14
void denom(int irrep, double root)
 
15
{
 
16
  int Gij, Gab;
 
17
  int ij, ab, i, j, a, b, I, J, A, B, isym, jsym, asym, bsym;
 
18
  int nirreps;
 
19
  int *occpi, *virtpi, *occ_off, *vir_off;
 
20
  int *aoccpi, *avirtpi, *aocc_off, *avir_off;
 
21
  int *boccpi, *bvirtpi, *bocc_off, *bvir_off;
 
22
  double fii, fjj, faa, fbb;
 
23
  dpdfile2 fIJ, fij, fAB, fab;
 
24
  dpdbuf4 D;
 
25
  char lbl[32];
 
26
 
 
27
  nirreps = moinfo.nirreps;
 
28
 
 
29
  if(params.ref == 0) { /** RHF **/
 
30
 
 
31
    occpi = moinfo.occpi;
 
32
    virtpi = moinfo.virtpi;
 
33
    occ_off = moinfo.occ_off;
 
34
    vir_off = moinfo.vir_off;
 
35
 
 
36
    dpd_file2_init(&fIJ, CC_OEI, 0, 0, 0, "fIJ");
 
37
    dpd_file2_mat_init(&fIJ);
 
38
    dpd_file2_mat_rd(&fIJ);
 
39
  
 
40
    dpd_file2_init(&fij, CC_OEI, 0, 0, 0, "fij");
 
41
    dpd_file2_mat_init(&fij);
 
42
    dpd_file2_mat_rd(&fij);
 
43
  
 
44
    dpd_file2_init(&fAB, CC_OEI, 0, 1, 1, "fAB");
 
45
    dpd_file2_mat_init(&fAB);
 
46
    dpd_file2_mat_rd(&fAB);
 
47
  
 
48
    dpd_file2_init(&fab, CC_OEI, 0, 1, 1, "fab");
 
49
    dpd_file2_mat_init(&fab);
 
50
    dpd_file2_mat_rd(&fab);
 
51
 
 
52
    sprintf(lbl, "dIjAb[%d]", irrep);
 
53
    dpd_buf4_init(&D, CC_MISC, irrep, 0, 5, 0, 5, 0, lbl);
 
54
    for(Gij=0; Gij < nirreps; Gij++) {
 
55
 
 
56
      dpd_buf4_mat_irrep_init(&D, Gij);
 
57
 
 
58
      for(ij=0; ij < D.params->rowtot[Gij]; ij++) {
 
59
        i = D.params->roworb[Gij][ij][0];
 
60
        j = D.params->roworb[Gij][ij][1];
 
61
        isym = D.params->psym[i];
 
62
        jsym = D.params->qsym[j];
 
63
 
 
64
        I = i - occ_off[isym];
 
65
        J = j - occ_off[jsym];
 
66
 
 
67
        fii = fIJ.matrix[isym][I][I];
 
68
        fjj = fij.matrix[jsym][J][J];
 
69
          
 
70
        for(ab=0; ab < D.params->coltot[Gij^irrep]; ab++) {
 
71
          a = D.params->colorb[Gij^irrep][ab][0];
 
72
          b = D.params->colorb[Gij^irrep][ab][1];
 
73
          asym = D.params->rsym[a];
 
74
          bsym = D.params->ssym[b];
 
75
 
 
76
          A = a - vir_off[asym];
 
77
          B = b - vir_off[bsym];
 
78
 
 
79
          faa = fAB.matrix[asym][A][A];
 
80
          fbb = fab.matrix[bsym][B][B];
 
81
 
 
82
          D.matrix[Gij][ij][ab] = 1.0/(fii + fjj - faa - fbb + root);
 
83
      }
 
84
    }
 
85
 
 
86
    dpd_buf4_mat_irrep_wrt(&D, Gij);
 
87
    dpd_buf4_mat_irrep_close(&D, Gij);
 
88
 
 
89
  }
 
90
 
 
91
  dpd_buf4_close(&D);
 
92
 
 
93
}
 
94
  else if(params.ref == 2) { /** UHF **/
 
95
 
 
96
    aoccpi = moinfo.aoccpi;
 
97
    boccpi = moinfo.boccpi;
 
98
    avirtpi = moinfo.avirtpi;
 
99
    bvirtpi = moinfo.bvirtpi;
 
100
    aocc_off = moinfo.aocc_off;
 
101
    bocc_off = moinfo.bocc_off;
 
102
    avir_off = moinfo.avir_off;
 
103
    bvir_off = moinfo.bvir_off;
 
104
 
 
105
    dpd_file2_init(&fIJ, CC_OEI, 0, 0, 0, "fIJ");
 
106
    dpd_file2_mat_init(&fIJ);
 
107
    dpd_file2_mat_rd(&fIJ);
 
108
  
 
109
    dpd_file2_init(&fij, CC_OEI, 0, 2, 2, "fij");
 
110
    dpd_file2_mat_init(&fij);
 
111
    dpd_file2_mat_rd(&fij);
 
112
  
 
113
    dpd_file2_init(&fAB, CC_OEI, 0, 1, 1, "fAB");
 
114
    dpd_file2_mat_init(&fAB);
 
115
    dpd_file2_mat_rd(&fAB);
 
116
  
 
117
    dpd_file2_init(&fab, CC_OEI, 0, 3, 3, "fab");
 
118
    dpd_file2_mat_init(&fab);
 
119
    dpd_file2_mat_rd(&fab);
 
120
 
 
121
    sprintf(lbl, "dIJAB[%d]", irrep);
 
122
    dpd_buf4_init(&D, CC_MISC, irrep, 1, 6, 1, 6, 0, lbl);
 
123
    for(Gij=0; Gij < nirreps; Gij++) {
 
124
      dpd_buf4_mat_irrep_init(&D, Gij);
 
125
      for(ij=0; ij < D.params->rowtot[Gij]; ij++) {
 
126
        i = D.params->roworb[Gij][ij][0];
 
127
        j = D.params->roworb[Gij][ij][1];
 
128
        isym = D.params->psym[i];
 
129
        jsym = D.params->qsym[j];
 
130
        I = i - aocc_off[isym];
 
131
        J = j - aocc_off[jsym];
 
132
        fii = fIJ.matrix[isym][I][I];
 
133
        fjj = fIJ.matrix[jsym][J][J];
 
134
 
 
135
        for(ab=0; ab < D.params->coltot[Gij^irrep]; ab++) {
 
136
          a = D.params->colorb[Gij^irrep][ab][0];
 
137
          b = D.params->colorb[Gij^irrep][ab][1];
 
138
          asym = D.params->rsym[a];
 
139
          bsym = D.params->ssym[b];
 
140
          A = a - avir_off[asym];
 
141
          B = b - avir_off[bsym];
 
142
          faa = fAB.matrix[asym][A][A];
 
143
          fbb = fAB.matrix[bsym][B][B];
 
144
 
 
145
          D.matrix[Gij][ij][ab] = 1.0/(fii + fjj - faa - fbb + root);
 
146
        }
 
147
      }
 
148
      dpd_buf4_mat_irrep_wrt(&D, Gij);
 
149
      dpd_buf4_mat_irrep_close(&D, Gij);
 
150
    }
 
151
    dpd_buf4_close(&D);
 
152
 
 
153
    sprintf(lbl, "dijab[%d]", irrep);
 
154
    dpd_buf4_init(&D, CC_MISC, irrep, 11, 16, 11, 16, 0, lbl);
 
155
    for(Gij=0; Gij < nirreps; Gij++) {
 
156
      dpd_buf4_mat_irrep_init(&D, Gij);
 
157
      for(ij=0; ij < D.params->rowtot[Gij]; ij++) {
 
158
        i = D.params->roworb[Gij][ij][0];
 
159
        j = D.params->roworb[Gij][ij][1];
 
160
        isym = D.params->psym[i];
 
161
        jsym = D.params->qsym[j];
 
162
        I = i - bocc_off[isym];
 
163
        J = j - bocc_off[jsym];
 
164
        fii = fij.matrix[isym][I][I];
 
165
        fjj = fij.matrix[jsym][J][J];
 
166
 
 
167
        for(ab=0; ab < D.params->coltot[Gij^irrep]; ab++) {
 
168
          a = D.params->colorb[Gij^irrep][ab][0];
 
169
          b = D.params->colorb[Gij^irrep][ab][1];
 
170
          asym = D.params->rsym[a];
 
171
          bsym = D.params->ssym[b];
 
172
          A = a - bvir_off[asym];
 
173
          B = b - bvir_off[bsym];
 
174
          faa = fab.matrix[asym][A][A];
 
175
          fbb = fab.matrix[bsym][B][B];
 
176
 
 
177
          D.matrix[Gij][ij][ab] = 1.0/(fii + fjj - faa - fbb + root);
 
178
        }
 
179
      }
 
180
      dpd_buf4_mat_irrep_wrt(&D, Gij);
 
181
      dpd_buf4_mat_irrep_close(&D, Gij);
 
182
    }
 
183
    dpd_buf4_close(&D);
 
184
 
 
185
    sprintf(lbl, "dIjAb[%d]", irrep);
 
186
    dpd_buf4_init(&D, CC_MISC, irrep, 22, 28, 22, 28, 0, lbl);
 
187
    for(Gij=0; Gij < nirreps; Gij++) {
 
188
      dpd_buf4_mat_irrep_init(&D, Gij);
 
189
      for(ij=0; ij < D.params->rowtot[Gij]; ij++) {
 
190
        i = D.params->roworb[Gij][ij][0];
 
191
        j = D.params->roworb[Gij][ij][1];
 
192
        isym = D.params->psym[i];
 
193
        jsym = D.params->qsym[j];
 
194
        I = i - aocc_off[isym];
 
195
        J = j - bocc_off[jsym];
 
196
        fii = fIJ.matrix[isym][I][I];
 
197
        fjj = fij.matrix[jsym][J][J];
 
198
 
 
199
        for(ab=0; ab < D.params->coltot[Gij^irrep]; ab++) {
 
200
          a = D.params->colorb[Gij^irrep][ab][0];
 
201
          b = D.params->colorb[Gij^irrep][ab][1];
 
202
          asym = D.params->rsym[a];
 
203
          bsym = D.params->ssym[b];
 
204
          A = a - avir_off[asym];
 
205
          B = b - bvir_off[bsym];
 
206
          faa = fAB.matrix[asym][A][A];
 
207
          fbb = fab.matrix[bsym][B][B];
 
208
 
 
209
          D.matrix[Gij][ij][ab] = 1.0/(fii + fjj - faa - fbb + root);
 
210
        }
 
211
      }
 
212
      dpd_buf4_mat_irrep_wrt(&D, Gij);
 
213
      dpd_buf4_mat_irrep_close(&D, Gij);
 
214
    }
 
215
    dpd_buf4_close(&D);
 
216
 
 
217
    dpd_file2_mat_close(&fIJ);
 
218
    dpd_file2_mat_close(&fij);
 
219
    dpd_file2_mat_close(&fAB);
 
220
    dpd_file2_mat_close(&fab);
 
221
    dpd_file2_close(&fIJ);
 
222
    dpd_file2_close(&fij);
 
223
    dpd_file2_close(&fAB);
 
224
    dpd_file2_close(&fab);
 
225
  }
 
226
}
 
227
 
 
228
}} // namespace psi::cis