← Back to branch summary

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

« back to all changes in this revision

Viewing changes to src/bin/cis/denom.c

  • Committer: Bazaar Package Importer
  • Author(s): Michael Banck
  • Date: 2006-09-10 14:01:33 UTC
  • Revision ID: james.westby@ubuntu.com-20060910140133-ib2j86trekykfsfv
Tags: upstream-3.2.3
Import upstream version 3.2.3

Show diffs side-by-side

added added

removed removed

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