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

« back to all changes in this revision

Viewing changes to src/bin/ccresponse/cc2_X1.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
#include <stdio.h>
 
2
#include <string.h>
 
3
#include <libdpd/dpd.h>
 
4
#include <libqt/qt.h>
 
5
#define EXTERN
 
6
#include "globals.h"
 
7
 
 
8
void denom1(dpdfile2 *X1, double omega);
 
9
void local_filter_T1(dpdfile2 *T1, double omega);
 
10
 
 
11
void cc2_X1_build(char *pert, char *cart, int irrep, double omega)
 
12
{
 
13
  int GX2, GX1, GW, Gam, Gim, Gef, Ga, Gi, Gm;
 
14
  int a, A, i, I, num_m, nlinks, length;
 
15
  dpdfile2 F, X1, X1new, Xme, Zia;
 
16
  dpdbuf4 W, X2, D, T2;
 
17
  char lbl[32];
 
18
 
 
19
  sprintf(lbl, "%sBAR_%1s_IA", pert, cart);
 
20
  dpd_file2_init(&X1new, CC_OEI, irrep, 0, 1, lbl);
 
21
  sprintf(lbl, "New X_%s_%1s_IA (%5.3f)", pert, cart, omega);
 
22
  dpd_file2_copy(&X1new, CC_OEI, lbl);
 
23
  dpd_file2_close(&X1new);
 
24
  dpd_file2_init(&X1new, CC_OEI, irrep, 0, 1, lbl);
 
25
 
 
26
  /*** S-S ***/
 
27
 
 
28
  sprintf(lbl, "X_%s_%1s_IA (%5.3f)", pert, cart, omega);
 
29
  dpd_file2_init(&X1, CC_OEI, irrep, 0, 1, lbl);
 
30
 
 
31
  dpd_file2_axpy(&X1, &X1new, -omega, 0);
 
32
 
 
33
  dpd_file2_init(&F, CC_OEI, 0, 1, 1, "FAE");
 
34
  dpd_contract222(&X1, &F, &X1new, 0, 0, 1, 1);
 
35
  dpd_file2_close(&F);
 
36
 
 
37
  dpd_file2_init(&F, CC_OEI, 0, 0, 0, "FMI");
 
38
  dpd_contract222(&F, &X1, &X1new, 1, 1, -1, 1);
 
39
  dpd_file2_close(&F);
 
40
 
 
41
  dpd_buf4_init(&W, CC2_HET1, 0, 10, 10, 10, 10, 0, "CC2 2 W(jb,ME) + W(Jb,Me)");
 
42
  dpd_contract422(&W, &X1, &X1new, 0, 0, 1, 1);
 
43
  dpd_buf4_close(&W);
 
44
 
 
45
  sprintf(lbl, "X_%s_%1s_ME", pert, cart);
 
46
  dpd_file2_init(&Xme, CC_OEI, irrep, 0, 1, lbl);
 
47
  dpd_buf4_init(&D, CC_DINTS, 0, 10, 10, 10, 10, 0, "D 2<ij|ab> - <ij|ba> (ia,jb)");
 
48
  dpd_contract422(&D, &X1, &Xme, 0, 0, 1, 0);
 
49
  dpd_buf4_close(&D);
 
50
 
 
51
  dpd_buf4_init(&T2, CC_TAMPS, 0, 10, 10, 10, 10, 0, "2 tIAjb - tIBja");
 
52
  dpd_contract422(&T2, &Xme, &X1new, 0, 0, 1, 1);
 
53
  dpd_buf4_close(&T2);
 
54
  dpd_file2_close(&Xme);
 
55
 
 
56
  dpd_file2_close(&X1);
 
57
 
 
58
 
 
59
  /*** S-D ***/
 
60
 
 
61
  dpd_file2_init(&F, CC_OEI, 0, 0, 1, "FME");
 
62
  sprintf(lbl, "X_%s_%1s_(2IjAb-IjbA) (%5.3f)", pert, cart, omega);
 
63
  dpd_buf4_init(&X2, CC_LR, irrep, 0, 5, 0, 5, 0, lbl);
 
64
  dpd_dot24(&F, &X2, &X1new, 0, 0, 1, 1);
 
65
  dpd_file2_close(&F);
 
66
 
 
67
  /** begin out of core contract442 **/
 
68
  dpd_buf4_init(&W, CC_HBAR, 0, 11, 5, 11, 5, 0, "WAmEf");
 
69
  GW = W.file.my_irrep;
 
70
  GX2 = X2.file.my_irrep;
 
71
  GX1 = X1new.my_irrep;
 
72
 
 
73
  dpd_file2_mat_init(&X1new);
 
74
  dpd_file2_mat_rd(&X1new);
 
75
 
 
76
  for(Gam=0; Gam < moinfo.nirreps; Gam++) {
 
77
    Gef = Gam^GW;
 
78
    Gim = Gef^GX2;
 
79
 
 
80
    for(Gi=0; Gi < moinfo.nirreps; Gi++) {
 
81
      Gm = Gef^Gi^GX2;
 
82
      Ga = Gi^GX1;
 
83
 
 
84
      num_m = W.params->qpi[Gm];
 
85
      dpd_buf4_mat_irrep_init_block(&W, Gam, num_m);
 
86
      dpd_buf4_mat_irrep_init_block(&X2, Gim, num_m);
 
87
 
 
88
      nlinks = W.params->coltot[Gam] * num_m;
 
89
      length = X1new.params->rowtot[Gi] * X1new.params->coltot[Ga];
 
90
 
 
91
      if(nlinks && length) {
 
92
        for(i=0; i < X1new.params->rowtot[Gi]; i++) {
 
93
 
 
94
          I = X2.params->poff[Gi] + i;
 
95
          dpd_buf4_mat_irrep_rd_block(&X2, Gim, X2.row_offset[Gim][I], num_m);
 
96
 
 
97
          for(a=0; a < X1new.params->coltot[Ga]; a++) {
 
98
 
 
99
            A = W.params->poff[Ga] + a;
 
100
            dpd_buf4_mat_irrep_rd_block(&W, Gam, W.row_offset[Gam][A], num_m);
 
101
 
 
102
            X1new.matrix[Gi][i][a] += C_DDOT(nlinks, X2.matrix[Gim][0], 1,
 
103
                                             W.matrix[Gam][0], 1);
 
104
          }
 
105
        }
 
106
      }
 
107
      dpd_buf4_mat_irrep_close_block(&X2, Gim, num_m);
 
108
      dpd_buf4_mat_irrep_close_block(&W, Gam, num_m);
 
109
    }
 
110
  }
 
111
 
 
112
  dpd_file2_mat_wrt(&X1new);
 
113
  dpd_file2_mat_close(&X1new);
 
114
 
 
115
  dpd_buf4_close(&W);
 
116
  dpd_buf4_close(&X2);
 
117
 
 
118
  sprintf(lbl, "X_%s_%1s_IjAb (%5.3f)", pert, cart, omega);
 
119
  dpd_buf4_init(&X2, CC_LR, irrep, 0, 5, 0, 5, 0, lbl);
 
120
  dpd_buf4_init(&W, CC_HBAR, 0, 0, 11, 0, 11, 0, "WMnIe - 2WnMIe (Mn,eI)");
 
121
  dpd_contract442(&W, &X2, &X1new, 3, 3, 1, 1);
 
122
  dpd_buf4_close(&W);
 
123
  dpd_buf4_close(&X2);
 
124
 
 
125
  if(params.local && local.filter_singles) local_filter_T1(&X1new, omega);
 
126
  else denom1(&X1new, omega);
 
127
  dpd_file2_close(&X1new);
 
128
}