← Back to branch summary

~ubuntu-branches/ubuntu/precise/psicode/precise

« back to all changes in this revision

Viewing changes to src/lib/libdpd/dot13.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:
src/lib/libdpd/dot13.c
1
 
#include <stdio.h>
2
 
#include <libqt/qt.h>
3
 
#include "dpd.h"
4
 
 
5
 
/* the non-symmetric transpose cases have not been tested */
6
 
 
7
 
int dpd_dot13(dpdfile2 *T, dpdbuf4 *I, dpdfile2 *Z,
8
 
    int transt, int transz, double alpha, double beta)
9
 
{
10
 
  int h, Gp, Gq, Gr, Gs, GT, GI, GZ, Tblock, Zblock;
11
 
  int p, q, r, s;
12
 
  int P, Q, R, S;
13
 
  int row, col;
14
 
  int nirreps;
15
 
  double **X;
16
 
  double value;
17
 
 
18
 
  nirreps = T->params->nirreps;
19
 
  GT = T->my_irrep;
20
 
  GI = I->file.my_irrep;
21
 
  GZ = Z->my_irrep;
22
 
 
23
 
  /* Get the two-index quantities from disk */
24
 
  dpd_file2_mat_init(T);
25
 
  dpd_file2_mat_rd(T);
26
 
 
27
 
  dpd_file2_scm(Z, beta);
28
 
  dpd_file2_mat_init(Z);
29
 
  dpd_file2_mat_rd(Z);
30
 
 
31
 
#ifdef DPD_TIMER
32
 
  timer_on("dot13");
33
 
#endif
34
 
 
35
 
  /* loop over row irreps of bufI; h = Gpq = Grs^GI */
36
 
  for(h=0; h < nirreps; h++) {
37
 
 
38
 
    dpd_buf4_mat_irrep_init(I, h);
39
 
    dpd_buf4_mat_irrep_rd(I, h);
40
 
 
41
 
    /* Loop over row irreps of target Z, Gz = Gqs */
42
 
    for(Gq=0; Gq < nirreps; Gq++) {
43
 
      /* Gs = Gq;  Gp = Gr = h^Gq; */
44
 
      Gp = h^Gq; Gr = h^Gq^GT; Gs = Gq^GZ;
45
 
      if (!transt) Tblock = Gp; else Tblock = Gr;
46
 
      if (!transz) Zblock = Gq; else Zblock = Gs;
47
 
 
48
 
      /* Allocate space for the X buffer */
49
 
      if(T->params->ppi[Gp] && T->params->qpi[Gr])
50
 
        X = dpd_block_matrix(T->params->ppi[Gp],T->params->qpi[Gr]);
51
 
 
52
 
      /* Loop over orbitals of the target */
53
 
      for(q=0; q < Z->params->ppi[Gq]; q++) {
54
 
        Q = Z->params->poff[Gq] + q;
55
 
        for(s=0; s < Z->params->qpi[Gs]; s++) {
56
 
          S = Z->params->qoff[Gs] + s;
57
 
 
58
 
          /* Loop over orbitals of the two-index term */
59
 
          for(p=0; p < T->params->ppi[Gp]; p++) {
60
 
            P = T->params->poff[Gp] + p;
61
 
            for(r=0; r < T->params->qpi[Gr]; r++) {
62
 
              R = T->params->qoff[Gr] + r;
63
 
 
64
 
              /* Calculate row and column indices in I */
65
 
              if(!transt && !transz) {
66
 
                row = I->params->rowidx[P][Q];
67
 
                col = I->params->colidx[R][S];
68
 
              }
69
 
              else if(transt && !transz) {
70
 
                row = I->params->rowidx[R][Q];
71
 
                col = I->params->colidx[P][S];
72
 
              }
73
 
              else if(!transt && transz) {
74
 
                row = I->params->rowidx[P][S];
75
 
                col = I->params->colidx[R][Q];
76
 
              }
77
 
              else if(transt && transz) {
78
 
                row = I->params->rowidx[R][S];
79
 
                col = I->params->colidx[P][Q];
80
 
              }
81
 
 
82
 
              /* Build the X buffer */
83
 
              X[p][r] = I->matrix[h][row][col]; 
84
 
 
85
 
            }
86
 
          }
87
 
 
88
 
          value = dot_block(T->matrix[Tblock], X, T->params->ppi[Gp],
89
 
              T->params->qpi[Gr], alpha); 
90
 
 
91
 
          Z->matrix[Zblock][q][s] += value;
92
 
        }
93
 
      }
94
 
      if(T->params->ppi[Gp] && T->params->qpi[Gr])
95
 
        dpd_free_block(X, T->params->ppi[Gp],T->params->qpi[Gr]);
96
 
    }
97
 
    dpd_buf4_mat_irrep_close(I, h);
98
 
  }
99
 
 
100
 
#ifdef DPD_TIMER
101
 
  timer_off("dot13");
102
 
#endif
103
 
 
104
 
  /* Close the two-index quantities */
105
 
  dpd_file2_mat_close(T);
106
 
  dpd_file2_mat_wrt(Z);
107
 
  dpd_file2_mat_close(Z);
108
 
 
109
 
  return 0;
110
 
}