~ubuntu-branches/ubuntu/vivid/psicode/vivid

« back to all changes in this revision

Viewing changes to src/bin/ccresponse/LHX1Y2.c

  • Committer: Bazaar Package Importer
  • Author(s): Michael Banck
  • Date: 2008-06-07 16:49:57 UTC
  • mfrom: (2.1.2 hardy)
  • Revision ID: james.westby@ubuntu.com-20080607164957-8pifvb133yjlkagn
Tags: 3.3.0-3
* debian/rules (DEB_MAKE_CHECK_TARGET): Do not abort test suite on
  failures.
* debian/rules (DEB_CONFIGURE_EXTRA_FLAGS): Set ${bindir} to /usr/lib/psi.
* debian/rules (install/psi3): Move psi3 file to /usr/bin.
* debian/patches/07_464867_move_executables.dpatch: New patch, add
  /usr/lib/psi to the $PATH, so that the moved executables are found.
  (closes: #464867)
* debian/patches/00list: Adjusted.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
1
#include <stdio.h>
2
2
#include <string.h>
 
3
#include <libciomr/libciomr.h>
3
4
#include <libdpd/dpd.h>
 
5
#include <libqt/qt.h>
4
6
#define EXTERN
5
7
#include "globals.h"
6
8
 
11
13
  dpdbuf4 Z, Z1, Z2, I, Y2, L2, W;
12
14
  char lbl[32];
13
15
  double polar;
 
16
  int nirreps, Gbm, Gef, Gjf, Ge, Gf, Gj, bm, ef, jf;
 
17
  int Gfe, Gb, Gm, b, m, B, M, e, f, fe, nrows, ncols, nlinks;
 
18
  double *X;
14
19
 
15
20
  sprintf(lbl, "Z_%s_%1s_MI", pert_y, cart_y);
16
21
  dpd_file2_init(&z1, CC_TMP0, irrep_y, 0, 0, lbl);
65
70
  dpd_file2_close(&z);
66
71
 
67
72
  dpd_file2_init(&z, CC_TMP0, 0, 0, 1, "Z(I,A) Final");
68
 
  dpd_file2_init(&l1, CC_OEI, 0, 0, 1, "LIA");
 
73
  dpd_file2_init(&l1, CC_LAMPS, 0, 0, 1, "LIA 0 -1");
69
74
  polar = 2.0 * dpd_file2_dot(&z, &l1);
70
75
  dpd_file2_close(&l1);
71
76
  dpd_file2_close(&z);
130
135
  sprintf(lbl, "X_%s_%1s_IA (%5.3f)", pert_x, cart_x, omega_x);
131
136
  dpd_file2_init(&X1, CC_OEI, irrep_x, 0, 1, lbl);
132
137
 
133
 
  sprintf(lbl, "Z_%s_%1s_MbEj (Mb,jE)", pert_x, cart_x);
134
 
  dpd_buf4_init(&Z1, CC_TMP0, irrep_x, 10, 10, 10, 10, 0, lbl);
135
 
  dpd_buf4_init(&W, CC_HBAR, 0, 10, 5, 10, 5, 0, "WAmEf");
136
 
  dpd_contract244(&X1, &W, &Z1, 1, 2, 1, 1, 0);
 
138
  sprintf(lbl, "Z_%s_%1s_MbEj (bM,jE)", pert_x, cart_x);
 
139
  dpd_buf4_init(&Z1, CC_TMP0, irrep_x, 11, 10, 11, 10, 0, lbl);
 
140
  dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf");
 
141
  /*  dpd_contract244(&X1, &W, &Z1, 1, 2, 1, 1, 0); */
 
142
  /* add out-of-core algorithm here */
 
143
  /* Z(bm,jf) <-- X1(j,e) * W(bm,ef) */
 
144
  dpd_file2_mat_init(&X1);
 
145
  dpd_file2_mat_rd(&X1);
 
146
  for(Gbm=0; Gbm < moinfo.nirreps; Gbm++) {
 
147
    Gef = Gbm;  Gjf = Gbm^irrep_x;
 
148
 
 
149
    dpd_buf4_mat_irrep_init(&Z1, Gbm);
 
150
    dpd_buf4_mat_irrep_row_init(&W, Gbm);
 
151
    for(bm=0; bm < W.params->rowtot[Gbm]; bm++) {
 
152
      dpd_buf4_mat_irrep_row_rd(&W, Gbm, bm);
 
153
 
 
154
      for(Gj=0; Gj < moinfo.nirreps; Gj++) {
 
155
        Gf = Gjf^Gj;  Ge = Gef^Gf;
 
156
 
 
157
        ef = W.col_offset[Gbm][Ge];
 
158
        jf = Z1.col_offset[Gbm][Gj];
 
159
 
 
160
        nrows = moinfo.occpi[Gj];
 
161
        ncols = moinfo.virtpi[Gf];
 
162
        nlinks = moinfo.virtpi[Ge];
 
163
 
 
164
        if(nrows && ncols && nlinks)
 
165
          C_DGEMM('n','n', nrows, ncols, nlinks, 1.0, &(X1.matrix[Gj][0][0]), nlinks,
 
166
                  &(W.matrix[Gbm][0][ef]), ncols, 0.0, &(Z1.matrix[Gbm][bm][jf]), ncols);
 
167
      }
 
168
    }
 
169
    dpd_buf4_mat_irrep_row_close(&W, Gbm);
 
170
    dpd_buf4_mat_irrep_wrt(&Z1, Gbm);
 
171
    dpd_buf4_mat_irrep_close(&Z1, Gbm);
 
172
 
 
173
  }
 
174
  dpd_file2_mat_close(&X1);
 
175
  /* out-of-core algorithm done */
137
176
  dpd_buf4_close(&W);
 
177
 
138
178
  sprintf(lbl, "Z_%s_%1s_MbEj (ME,jb)", pert_x, cart_x);
139
 
  dpd_buf4_sort(&Z1, CC_TMP0, psrq, 10, 10, lbl);
 
179
  dpd_buf4_sort(&Z1, CC_TMP0, qsrp, 10, 10, lbl);
140
180
  dpd_buf4_close(&Z1);
141
181
 
142
182
  sprintf(lbl, "Z_%s_%1s_WMbEj (bM,Ej)", pert_x, cart_x);
143
183
  dpd_buf4_init(&Z1, CC_TMP0, irrep_x, 11, 11, 11, 11, 0, lbl);
144
 
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe");
 
184
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe (Mn,eI)");
145
185
  dpd_contract244(&X1, &W, &Z1, 0, 0, 0, -1, 0);
146
186
  dpd_buf4_close(&W);
147
187
  sprintf(lbl, "Z_%s_%1s_MbEj (ME,jb)", pert_x, cart_x);
150
190
 
151
191
  sprintf(lbl, "Z_%s_%1s_MbeJ (Mb,eJ)", pert_x, cart_x);
152
192
  dpd_buf4_init(&Z1, CC_TMP0, irrep_x, 10, 11, 10, 11, 0, lbl);
153
 
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe");
 
193
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe (Mn,eI)");
154
194
  dpd_contract424(&W, &X1, &Z1, 1, 0, 1, 1, 0);
155
195
  dpd_buf4_close(&W);
156
 
  dpd_buf4_init(&W, CC_HBAR, 0, 10, 5, 10, 5, 0, "WAmEf");
157
 
  dpd_contract424(&W, &X1, &Z1, 3, 1, 0, -1, 1);
158
 
  dpd_buf4_close(&W);
159
196
  sprintf(lbl, "Z_%s_%1s_MbeJ (Me,Jb)", pert_x, cart_x);
160
197
  dpd_buf4_sort(&Z1, CC_TMP0, prsq, 10, 10, lbl);
161
198
  dpd_buf4_close(&Z1);
162
199
 
 
200
  sprintf(lbl, "Z_%s_%1s_MbeJ (bM,eJ)", pert_x, cart_x);
 
201
  dpd_buf4_init(&Z1, CC_TMP0, irrep_x, 11, 11, 11, 11, 0, lbl);
 
202
  dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf");
 
203
  dpd_contract424(&W, &X1, &Z1, 3, 1, 0, -1, 0);
 
204
  dpd_buf4_close(&W);
 
205
  sprintf(lbl, "Z_%s_%1s_MbeJ (Me,Jb)", pert_x, cart_x);
 
206
  dpd_buf4_sort_axpy(&Z1, CC_TMP0, qrsp, 10, 10, lbl, 1);
 
207
  dpd_buf4_close(&Z1);
 
208
 
163
209
  dpd_file2_close(&X1);
164
210
 
165
211
  sprintf(lbl, "Z_%s_%1s_MbEj (ME,jb)", pert_x, cart_x);
230
276
 
231
277
  sprintf(lbl, "Z_%s_%1s_mj" , pert_x, cart_x);
232
278
  dpd_file2_init(&z, CC_TMP0, irrep_x, 0, 0, lbl);
233
 
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "2WMnIe - WnMIe");
 
279
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "2WMnIe - WnMIe (Mn,eI)");
234
280
  dpd_dot23(&X1, &W, &z, 0, 0, 1, 0);
235
281
  dpd_buf4_close(&W);
236
282
  dpd_buf4_init(&Z1, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z(Ij,Ab)");
243
289
  dpd_file2_close(&z);
244
290
 
245
291
 
246
 
 
247
 
  sprintf(lbl, "Z_%s_%1s_ae", pert_x, cart_x);
 
292
  sprintf(lbl, "Z_%s_%1s_AE", pert_x, cart_x);
248
293
  dpd_file2_init(&z, CC_TMP0, irrep_x, 1, 1, lbl);
249
294
 
250
 
  dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf 2(Am,Ef) - (Am,fE)");
251
 
  dpd_dot24(&X1, &W, &z, 0, 0, 1, 0);
 
295
/*   dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf 2(Am,Ef) - (Am,fE)"); */
 
296
/*   dpd_dot24(&X1, &W, &z, 0, 0, 1, 0); */
 
297
  dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf");
 
298
  dpd_file2_scm(&z, 0);
 
299
  dpd_file2_mat_init(&z);
 
300
  dpd_file2_mat_init(&X1);
 
301
  dpd_file2_mat_rd(&X1);
 
302
  for(Gbm=0; Gbm < moinfo.nirreps; Gbm++) {
 
303
    Gfe = Gbm; /* W is totally symmetric */
 
304
    dpd_buf4_mat_irrep_row_init(&W, Gbm);
 
305
    X = init_array(W.params->coltot[Gfe]);
 
306
    for(bm=0; bm < W.params->rowtot[Gbm]; bm++) {
 
307
      dpd_buf4_mat_irrep_row_rd(&W, Gbm, bm);
 
308
      b = W.params->roworb[Gbm][bm][0];
 
309
      m = W.params->roworb[Gbm][bm][1];
 
310
      Gb = W.params->psym[b];
 
311
      Gm = Gbm ^ Gb;
 
312
      Ge = Gm ^ irrep_x;
 
313
      Gf = Gfe ^ Ge;
 
314
      B = b - moinfo.vir_off[Gb];
 
315
      M = m - moinfo.occ_off[Gm];
 
316
      zero_arr(X, W.params->coltot[Gfe]);
 
317
      for(fe=0; fe < W.params->coltot[Gfe]; fe++) {
 
318
        f = W.params->colorb[Gfe][fe][0];
 
319
        e = W.params->colorb[Gfe][fe][1];
 
320
        ef = W.params->colidx[e][f];
 
321
        X[fe] = 2.0 * W.matrix[Gbm][0][fe] - W.matrix[Gbm][0][ef];
 
322
      }
 
323
      nrows = moinfo.virtpi[Gf];
 
324
      ncols = moinfo.virtpi[Ge];
 
325
      if(nrows & ncols) {
 
326
        C_DGEMV('n',nrows,ncols,1,&X[W.col_offset[Gfe][Gf]],ncols,
 
327
                X1.matrix[Gm][M],1,1,z.matrix[Gb][B],1);
 
328
      }
 
329
    }
 
330
    free(X);
 
331
    dpd_buf4_mat_irrep_row_close(&W, Gbm);
 
332
  }
 
333
  dpd_file2_mat_close(&X1);
 
334
  dpd_file2_mat_wrt(&z);
 
335
  dpd_file2_mat_close(&z);
252
336
  dpd_buf4_close(&W);
 
337
 
253
338
  dpd_buf4_init(&Z1, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z(Ij,Ab) temp");
254
339
  dpd_contract424(&Y2, &z, &Z1, 3, 1, 0, 1, 0);
255
340
  dpd_buf4_init(&Z, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z(Ij,Ab) Final");
267
352
  dpd_file2_init(&X1, CC_OEI, irrep_x, 0, 1, lbl);
268
353
  sprintf(lbl, "Z_%s_%1s_MnjI", pert_x, cart_x);
269
354
  dpd_buf4_init(&Z1, CC_TMP0, irrep_x, 0, 0, 0, 0, 0, lbl);
270
 
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe");
 
355
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe (Mn,eI)");
271
356
  dpd_contract244(&X1, &W, &Z1, 1, 2, 1, 1, 0);
272
357
  dpd_buf4_close(&W);
273
358
  sprintf(lbl, "Z_%s_%1s_MnIj", pert_x, cart_x);
286
371
  dpd_buf4_close(&Z1);
287
372
  dpd_buf4_close(&Z);
288
373
 
289
 
  sprintf(lbl, "Z_%s_%1s_AmIj (mA,Ij)", pert_y, cart_y);
290
 
  dpd_buf4_init(&Z1, CC_TMP0, irrep_y, 10, 0, 10, 0, 0, lbl);
 
374
  sprintf(lbl, "Z_%s_%1s_AmIj (Am,Ij)", pert_y, cart_y);
 
375
  dpd_buf4_init(&Z1, CC_TMP0, irrep_y, 11, 0, 11, 0, 0, lbl);
291
376
  sprintf(lbl, "X_%s_%1s_IjAb (%5.3f)", pert_y, cart_y, omega_y);
292
377
  dpd_buf4_init(&Y2, CC_LR, irrep_y, 0, 5, 0, 5, 0, lbl);
293
 
  dpd_buf4_init(&W, CC_HBAR, 0, 10, 5, 10, 5, 0, "WAmEf");
 
378
  dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf");
294
379
  dpd_contract444(&W, &Y2, &Z1, 0, 0, 1, 0);
295
380
  dpd_buf4_close(&W);
296
381
  dpd_buf4_close(&Y2);
297
382
 
298
 
  dpd_buf4_init(&Z, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z(Ij,bA)");
 
383
  dpd_buf4_init(&Z, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z(Ij,Ab) Final");
299
384
  sprintf(lbl, "X_%s_%1s_IA (%5.3f)", pert_x, cart_x, omega_x);
300
385
  dpd_file2_init(&X1, CC_OEI, irrep_x, 0, 1, lbl);
301
 
  dpd_contract244(&X1, &Z1, &Z, 0, 0, 1, -2, 0);
 
386
  dpd_contract424(&Z1, &X1, &Z, 1, 0, 0, -2, 1);
302
387
  dpd_file2_close(&X1);
303
388
  dpd_buf4_close(&Z1);
304
 
  dpd_buf4_sort_axpy(&Z, CC_TMP0, pqsr, 0, 5, "Z(Ij,Ab) Final", 1);
305
389
  dpd_buf4_close(&Z);
306
390
 
307
391
  dpd_buf4_init(&Z, CC_TMP0, 0, 0, 5, 0, 5, 0, "Z(Ij,Ab) Final");