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

« back to all changes in this revision

Viewing changes to src/bin/detci/s1.c

  • 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
 
/*
2
 
** S1.C
3
 
** 
4
 
** File contains code to calculate sigma1 in various ways, all
5
 
** block-at-a-time now.
6
 
**
7
 
** C. David Sherrill
8
 
** Center for Computational Quantum Chemistry
9
 
** University of Georgia
10
 
** 
11
 
*/
12
 
 
13
 
#include <stdio.h>
14
 
#include <libciomr/libciomr.h>
15
 
#include "structs.h"
16
 
#define EXTERN
17
 
#include "globals.h"
18
 
 
19
 
extern unsigned char ***Occs;
20
 
 
21
 
#define MIN0(a,b) (((a)<(b)) ? (a) : (b))
22
 
#define MAX0(a,b) (((a)>(b)) ? (a) : (b))
23
 
 
24
 
extern void b2brepl(unsigned char **occs, int *Jcnt, int **Jij, int **Joij, 
25
 
      int **Jridx, signed char **Jsgn, struct olsen_graph *Graph,
26
 
      int Ilist, int Jlist, int len);
27
 
 
28
 
 
29
 
/*
30
 
** S1_BLOCK_FCI(): 
31
 
**
32
 
** Calculate the sigma_1 vector as described by
33
 
** equation (20) of RAS Paper (Olsen, Roos, Jorgensen, Aa. Jensen JCP 1988)
34
 
**
35
 
** This sigma1 routine is for Full CI's only.
36
 
** currently assumes that (ij|ij)'s have not been halved!! 
37
 
** 
38
 
** David Sherrill, 21 June 1995
39
 
** Based on previous code by David Sherrill, 1994
40
 
**
41
 
** Updated 3/27/94 to include g matrix for RAS
42
 
** Modified 4/8/94 to make C and s one-dimensional
43
 
** Modified 4/10/94 to make FCI-only (for now) and use new string structs
44
 
** Modified 6/21/95 for use in new RAS program
45
 
**
46
 
*/
47
 
void s1_block_fci(struct stringwr **alplist, struct stringwr **betlist, 
48
 
      double **C, double **S, double *oei, double *tei, double *F,
49
 
      int nlists, int nas, int nbs, int Ib_list, int Jb_list, 
50
 
      int Jb_list_nbs)
51
 
{
52
 
   struct stringwr *Ib, *Kb;
53
 
   unsigned int Ia_idx, Ib_idx, Kb_idx, Jb_idx;
54
 
   unsigned int Ibcnt, Kbcnt, Kb_list, Ib_ex, Kb_ex;
55
 
   unsigned int *Ibridx, *Kbridx;
56
 
   int nirreps, *Ibij, *Kbij;
57
 
   signed char *Ibsgn, *Kbsgn;
58
 
   int ij,kl,ijkl;
59
 
   double Kb_sgn, Jb_sgn;
60
 
   double tval;
61
 
 
62
 
   nirreps = CalcInfo.nirreps;
63
 
 
64
 
   /* loop over I_b */
65
 
   for (Ib=betlist[Ib_list], Ib_idx=0; Ib_idx < nbs; Ib_idx++, Ib++) {
66
 
      zero_arr(F, Jb_list_nbs);
67
 
 
68
 
      /* loop over excitations E^b_{kl} from |B(I_b)> */
69
 
      for (Kb_list=0; Kb_list < nlists; Kb_list++) {
70
 
         Ibcnt = Ib->cnt[Kb_list];
71
 
         Ibridx = Ib->ridx[Kb_list];
72
 
         Ibsgn = Ib->sgn[Kb_list];
73
 
         Ibij = Ib->ij[Kb_list];
74
 
         for (Ib_ex=0; Ib_ex < Ibcnt; Ib_ex++) {
75
 
            kl = *Ibij++;
76
 
            Kb_idx = *Ibridx++;
77
 
            Kb_sgn = (double) *Ibsgn++;
78
 
 
79
 
            /* B(K_b) = sgn(kl) * E^b_{kl} |B(I_b)> */
80
 
            Kb = betlist[Kb_list] + Kb_idx;
81
 
            if (Kb_list == Jb_list) F[Kb_idx] += Kb_sgn * oei[kl];
82
 
 
83
 
            /* loop over excitations E^b_{ij} from |B(K_b)> */
84
 
            /* Jb_list pre-determined because of C blocking */
85
 
            Kbcnt = Kb->cnt[Jb_list];
86
 
            Kbridx = Kb->ridx[Jb_list];
87
 
            Kbsgn = Kb->sgn[Jb_list];
88
 
            Kbij = Kb->ij[Jb_list];
89
 
            for (Kb_ex=0; Kb_ex < Kbcnt; Kb_ex++) {
90
 
               Jb_idx = *Kbridx++;
91
 
               Jb_sgn = (double) *Kbsgn++;
92
 
               ij = *Kbij++;
93
 
               ijkl = ioff[MAX0(ij,kl)] + MIN0(ij,kl) ;
94
 
               F[Jb_idx] += 0.5 * Kb_sgn * Jb_sgn * tei[ijkl] ;
95
 
               }
96
 
            } /* end loop over Ib excitations */
97
 
         } /* end loop over Kb_list */
98
 
 
99
 
      
100
 
      for (Ia_idx=0; Ia_idx < nas; Ia_idx++) {
101
 
         tval = 0.0;
102
 
         for (Jb_idx=0; Jb_idx < Jb_list_nbs; Jb_idx++) {
103
 
            tval += C[Ia_idx][Jb_idx] * F[Jb_idx];
104
 
            }
105
 
         S[Ia_idx][Ib_idx] += tval;
106
 
         }
107
 
      } /* end loop over Ib */
108
 
 
109
 
}
110
 
 
111
 
 
112
 
/*
113
 
** S1_BLOCK_RAS.C: 
114
 
** 
115
 
** Calculate the sigma_1 vector as described by
116
 
** equation (20) of RAS Paper (Olsen, Roos, Jorgensen, Aa. Jensen JCP 1988)
117
 
**
118
 
** This sigma1 routine is for RAS CI's.
119
 
** currently assumes that (ij|ij)'s have not been halved!!
120
 
** 
121
 
** David Sherrill, 2 August 1995
122
 
** Based on previous code by David Sherrill, 1994
123
 
**
124
 
** Updated 3/27/94 to include g matrix for RAS
125
 
** Modified 4/8/94 to make C and s one-dimensional
126
 
** Modified 4/10/94 to make FCI-only (for now) and use new string structs
127
 
** Modified 6/21/95 for use in new RAS program (C, s now 2D again!)
128
 
** Modified 8/2/95 to make RAS again
129
 
**
130
 
*/
131
 
void s1_block_ras(struct stringwr **alplist, struct stringwr **betlist, 
132
 
      double **C, double **S, double *oei, double *tei, double *F,
133
 
      int nlists, int nas, int nbs, int Ib_list, int Jb_list, 
134
 
      int Jb_list_nbs)
135
 
{
136
 
   struct stringwr *Ib, *Kb;
137
 
   unsigned int Ia_idx, Ib_idx, Kb_idx, Jb_idx;
138
 
   unsigned int Ibcnt, Kbcnt, Kb_list, Ib_ex, Kb_ex;
139
 
   unsigned int *Ibridx, *Kbridx;
140
 
   int nirreps,  *Ibij, *Kbij, *Iboij, *Kboij;
141
 
   signed char *Ibsgn, *Kbsgn;
142
 
   int ij,kl,ijkl,oij,okl;
143
 
   double Kb_sgn, Jb_sgn;
144
 
   double tval;
145
 
 
146
 
   nirreps = CalcInfo.nirreps;
147
 
 
148
 
   /* loop over I_b */
149
 
   for (Ib=betlist[Ib_list], Ib_idx=0; Ib_idx < nbs; Ib_idx++, Ib++) {
150
 
      zero_arr(F, Jb_list_nbs);
151
 
 
152
 
      /* loop over excitations E^b_{kl} from |B(I_b)> */
153
 
      for (Kb_list=0; Kb_list < nlists; Kb_list++) {
154
 
         Ibcnt = Ib->cnt[Kb_list];
155
 
         Ibridx = Ib->ridx[Kb_list];
156
 
         Ibsgn = Ib->sgn[Kb_list];
157
 
         Ibij = Ib->ij[Kb_list];
158
 
         Iboij = Ib->oij[Kb_list];
159
 
         for (Ib_ex=0; Ib_ex < Ibcnt; Ib_ex++) {
160
 
            kl = *Ibij++;
161
 
            okl = *Iboij++;
162
 
            Kb_idx = *Ibridx++;
163
 
            Kb_sgn = (double) *Ibsgn++;
164
 
 
165
 
            /* B(K_b) = sgn(kl) * E^b_{kl} |B(I_b)> */
166
 
            Kb = betlist[Kb_list] + Kb_idx;
167
 
            /* note okl on next line, not kl */
168
 
            if (Kb_list == Jb_list) F[Kb_idx] += Kb_sgn * oei[okl];
169
 
 
170
 
            /* loop over excitations E^b_{ij} from |B(K_b)> */
171
 
            /* Jb_list pre-determined because of C blocking */
172
 
            Kbcnt = Kb->cnt[Jb_list];
173
 
            Kbridx = Kb->ridx[Jb_list];
174
 
            Kbsgn = Kb->sgn[Jb_list];
175
 
            Kbij = Kb->ij[Jb_list];
176
 
            Kboij = Kb->oij[Jb_list];
177
 
            for (Kb_ex=0; Kb_ex < Kbcnt; Kb_ex++) {
178
 
               Jb_idx = *Kbridx++;
179
 
               Jb_sgn = (double) *Kbsgn++;
180
 
               ij = *Kbij++;
181
 
               oij = *Kboij++;
182
 
               ijkl = ioff[MAX0(ij,kl)] + MIN0(ij,kl) ;
183
 
               if (oij > okl) 
184
 
                  F[Jb_idx] += Kb_sgn * Jb_sgn * tei[ijkl] ;
185
 
               else if (oij == okl) 
186
 
                  F[Jb_idx] += 0.5 * Kb_sgn * Jb_sgn * tei[ijkl] ;
187
 
               }
188
 
            } /* end loop over Ib excitations */
189
 
      } /* end loop over Kb_list */
190
 
 
191
 
      
192
 
   for (Ia_idx=0; Ia_idx < nas; Ia_idx++) {
193
 
      tval = 0.0;
194
 
      for (Jb_idx=0; Jb_idx < Jb_list_nbs; Jb_idx++) {
195
 
         tval += C[Ia_idx][Jb_idx] * F[Jb_idx];
196
 
         }
197
 
      S[Ia_idx][Ib_idx] += tval;
198
 
      }
199
 
   } /* end loop over Ib */
200
 
 
201
 
}
202
 
 
203
 
 
204
 
 
205
 
/*
206
 
** S1_BLOCK_RAS_ROTF
207
 
** 
208
 
** String replacements on-the-fly version
209
 
**
210
 
** This sigma1 routine is for RAS CI's.
211
 
** currently assumes that (ij|ij)'s have not been halved!! 
212
 
** 
213
 
** David Sherrill, 13 August 1995
214
 
** Based on previous code by David Sherrill, 1994
215
 
**
216
 
** Updated 3/27/94 to include g matrix for RAS
217
 
** Modified 4/8/94 to make C and s one-dimensional
218
 
** Modified 4/10/94 to make FCI-only (for now) and use new string structs
219
 
** Modified 6/21/95 for use in new RAS program (C, s now 2D again!)
220
 
** Modified 8/2/95 to make RAS again
221
 
**
222
 
*/
223
 
void s1_block_ras_rotf(int *Cnt[2], int **Ij[2], int **Oij[2],
224
 
      int **Ridx[2], signed char **Sgn[2], unsigned char **Toccs,
225
 
      double **C, double **S,
226
 
      double *oei, double *tei, double *F, int nlists, int nas, int nbs,
227
 
      int Ib_list, int Jb_list, int Jb_list_nbs)
228
 
{
229
 
   int Ia_idx, Ib_idx, Kb_idx, Jb_idx;
230
 
   int Ibcnt, Kbcnt, Kb_list, Ib_ex, Kb_ex;
231
 
   int *Ibridx, *Kbridx;
232
 
   int nirreps,  *Ibij, *Kbij, *Iboij, *Kboij;
233
 
   signed char *Ibsgn, *Kbsgn;
234
 
   int i,ij,kl,ijkl,oij,okl;
235
 
   double Kb_sgn, Jb_sgn;
236
 
   double tval;
237
 
 
238
 
   nirreps = CalcInfo.nirreps;
239
 
 
240
 
   for (Kb_list=0; Kb_list < nlists; Kb_list++) {
241
 
      b2brepl(Occs[Ib_list], Cnt[0], Ij[0], Oij[0], Ridx[0],
242
 
         Sgn[0], BetaG, Ib_list, Kb_list, nbs);
243
 
 
244
 
      /* loop over I_b */
245
 
      for (Ib_idx=0; Ib_idx < nbs; Ib_idx++) {
246
 
 
247
 
         if ((Ibcnt = Cnt[0][Ib_idx]) < 0) continue;
248
 
         zero_arr(F, Jb_list_nbs);
249
 
 
250
 
         /* loop over excitations E^b_{kl} from |B(I_b)> */
251
 
         Ibridx = Ridx[0][Ib_idx];
252
 
         Ibsgn = Sgn[0][Ib_idx];
253
 
         Ibij = Ij[0][Ib_idx];
254
 
         Iboij = Oij[0][Ib_idx];
255
 
 
256
 
         for (i=0; i<Ibcnt; i++) 
257
 
            Toccs[i] = Occs[Kb_list][Ibridx[i]];
258
 
 
259
 
         b2brepl(Toccs, Cnt[1], Ij[1], Oij[1], Ridx[1], Sgn[1],
260
 
            BetaG, Kb_list, Jb_list, Ibcnt);
261
 
 
262
 
         for (Ib_ex=0; Ib_ex < Ibcnt; Ib_ex++) {
263
 
            kl = *Ibij++;
264
 
            okl = *Iboij++;
265
 
            Kb_idx = *Ibridx++;
266
 
            Kb_sgn = (double) *Ibsgn++;
267
 
 
268
 
            /* B(K_b) = sgn(kl) * E^b_{kl} |B(I_b)> */
269
 
            /* note okl on next line, not kl */
270
 
            if (Kb_list == Jb_list) F[Kb_idx] += Kb_sgn * oei[okl];
271
 
 
272
 
            /* loop over excitations E^b_{ij} from |B(K_b)> */
273
 
            /* Jb_list pre-determined because of C blocking */
274
 
            Kbcnt = Cnt[1][Ib_ex];
275
 
            Kbridx = Ridx[1][Ib_ex];
276
 
            Kbsgn = Sgn[1][Ib_ex];
277
 
            Kbij = Ij[1][Ib_ex];
278
 
            Kboij = Oij[1][Ib_ex];
279
 
            for (Kb_ex=0; Kb_ex < Kbcnt; Kb_ex++) {
280
 
               Jb_idx = *Kbridx++;
281
 
               Jb_sgn = (double) *Kbsgn++;
282
 
               ij = *Kbij++;
283
 
               oij = *Kboij++;
284
 
               ijkl = ioff[MAX0(ij,kl)] + MIN0(ij,kl) ;
285
 
               if (oij > okl) 
286
 
                  F[Jb_idx] += Kb_sgn * Jb_sgn * tei[ijkl] ;
287
 
               else if (oij == okl) 
288
 
                  F[Jb_idx] += 0.5 * Kb_sgn * Jb_sgn * tei[ijkl] ;
289
 
               }
290
 
            } /* end loop over Ib excitations */
291
 
 
292
 
      for (Ia_idx=0; Ia_idx < nas; Ia_idx++) {
293
 
         tval = 0.0;
294
 
         for (Jb_idx=0; Jb_idx < Jb_list_nbs; Jb_idx++) {
295
 
            tval += C[Ia_idx][Jb_idx] * F[Jb_idx];
296
 
            }
297
 
         S[Ia_idx][Ib_idx] += tval;
298
 
         }
299
 
      } /* end loop over Ib */
300
 
   } /* end loop over Kb_list */
301
 
 
302
 
}
303