← Back to branch summary

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

« back to all changes in this revision

Viewing changes to src/bin/ccdensity/dump_ROHF.cc

  • 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/bin/ccdensity/dump_ROHF.cc
 
1
/*! \file
 
2
    \ingroup CCDENSITY
 
3
    \brief Enter brief description of file here 
 
4
*/
 
5
#include <cstdio>
 
6
#include <libciomr/libciomr.h>
 
7
#include <libiwl/iwl.h>
 
8
#include <libdpd/dpd.h>
 
9
#include <psifiles.h>
 
10
#include "MOInfo.h"
 
11
#include "Params.h"
 
12
#include "Frozen.h"
 
13
#define EXTERN
 
14
#include "globals.h"
 
15
 
 
16
namespace psi { namespace ccdensity {
 
17
 
 
18
/* DUMP_ROHF(): Mulliken-order the ROHF-CCSD two-electron density and
 
19
** dump it to a file for subsequent backtransformation.  Basically all
 
20
** we have to do is swap indices two and three, e.g.
 
21
**
 
22
** G'(pr,qs) = G(pq,rs)
 
23
**
 
24
** In order for the Mulliken-ordered density to be valid for the
 
25
** backtransformation algorithm used in TRANSQT, the final density
 
26
** must have eight-fold permutational symmetry like the original
 
27
** integrals.  Unfortunately, there are a couple of complications
 
28
** introduced by the redundant storage I use for open-shell orbitals
 
29
** (useful for spin-restricted references --- see the CCSORT code). In
 
30
** particular, if the Mulliken-ordered density is not bra-ket
 
31
** symmetric, specific elements of the final density must be
 
32
** multiplied by two or they will not appear with the correct
 
33
** prefactor in the backtransformation.  This only affects the IJKA,
 
34
** IAJB, and ABCI Dirac-ordered densities, since the remaining three
 
35
** components are bra-ket symmetric in Mulliken order.
 
36
**
 
37
** I really need to give an example of this problem using specific
 
38
** elements of GIJKA so that the code below will be clearer.*/
 
39
 
 
40
void dump_ROHF(struct iwlbuf *OutBuf, struct RHO_Params rho_params)
 
41
{
 
42
  int nirreps, nmo, nfzv;
 
43
  int *qt_occ, *qt_vir;
 
44
  int h, row, col, p, q, r, s;
 
45
  dpdbuf4 G;
 
46
 
 
47
  qt_occ = moinfo.qt_occ;  qt_vir = moinfo.qt_vir;
 
48
  nirreps = moinfo.nirreps;
 
49
  nmo = moinfo.nmo;
 
50
  nfzv = moinfo.nfzv;
 
51
 
 
52
  psio_open(PSIF_MO_OPDM, PSIO_OPEN_OLD);
 
53
 /*  psio_write_entry(PSIF_MO_OPDM, "MO-basis OPDM", (char *) moinfo.opdm[0], */
 
54
  psio_write_entry(PSIF_MO_OPDM, rho_params.opdm_lbl, (char *) moinfo.opdm[0],
 
55
                   sizeof(double)*(nmo-nfzv)*(nmo-nfzv));
 
56
  psio_close(PSIF_MO_OPDM, 1);
 
57
 
 
58
if (!params.onepdm) {
 
59
  psio_open(PSIF_MO_LAG, PSIO_OPEN_OLD);
 
60
  psio_write_entry(PSIF_MO_LAG, "MO-basis Lagrangian", (char *) moinfo.I[0],
 
61
                   sizeof(double)*nmo*nmo);
 
62
  psio_close(PSIF_MO_LAG, 1);
 
63
 
 
64
  dpd_buf4_init(&G, CC_GAMMA, 0, 0, 0, 0, 0, 0, "GIjKl");
 
65
  dpd_buf4_sort(&G, CC_TMP0, prqs, 0, 0, "G(IK,JL)");
 
66
  dpd_buf4_close(&G);
 
67
  dpd_buf4_init(&G, CC_TMP0, 0, 0, 0, 0, 0, 0, "G(IK,JL)");
 
68
  dpd_buf4_dump(&G, OutBuf, qt_occ, qt_occ, qt_occ, qt_occ, 1, 0);
 
69
  dpd_buf4_close(&G);
 
70
 
 
71
  dpd_buf4_init(&G, CC_GAMMA, 0, 0, 10, 0, 10, 0, "GIjKa");
 
72
  dpd_buf4_sort(&G, CC_TMP0, prqs, 0, 10, "G(IK,JA)");
 
73
  dpd_buf4_close(&G);
 
74
  dpd_buf4_init(&G, CC_TMP0, 0, 0, 10, 0, 10, 0, "G(IK,JA)");
 
75
  
 
76
  for(h=0; h < nirreps; h++) {
 
77
    dpd_buf4_mat_irrep_init(&G, h);
 
78
    dpd_buf4_mat_irrep_rd(&G, h);
 
79
 
 
80
    for(row=0; row < G.params->rowtot[h]; row++) {
 
81
      p = G.params->roworb[h][row][0];
 
82
      q = G.params->roworb[h][row][1];
 
83
      for(col=0; col < G.params->coltot[h]; col++) {
 
84
        r = G.params->colorb[h][col][0];
 
85
        s = G.params->colorb[h][col][1];
 
86
 
 
87
        if((qt_occ[q] == qt_vir[s]) && (p == r))
 
88
          G.matrix[h][row][col] *= 2;
 
89
      }
 
90
    }
 
91
 
 
92
    dpd_buf4_mat_irrep_wrt(&G, h);
 
93
    dpd_buf4_mat_irrep_close(&G, h);
 
94
  }
 
95
 
 
96
  dpd_buf4_dump(&G, OutBuf, qt_occ, qt_occ, qt_occ, qt_vir, 0, 0);
 
97
  dpd_buf4_close(&G);
 
98
 
 
99
  dpd_buf4_init(&G, CC_GAMMA, 0, 0, 5, 0, 5, 0, "GIjAb");
 
100
  dpd_buf4_sort(&G, CC_TMP9, prqs, 10, 10, "G(IA,JB)");
 
101
  dpd_buf4_close(&G);
 
102
  dpd_buf4_init(&G, CC_TMP9, 0, 10, 10, 10, 10, 0, "G(IA,JB)");
 
103
  dpd_buf4_symm(&G);
 
104
  dpd_buf4_dump(&G, OutBuf, qt_occ, qt_vir, qt_occ, qt_vir, 1, 0);
 
105
  dpd_buf4_close(&G);
 
106
 
 
107
  dpd_buf4_init(&G, CC_GAMMA, 0, 10, 10, 10, 10, 0, "GIBJA");
 
108
  dpd_buf4_sort(&G, CC_TMP0, prqs, 0, 5, "G(IJ,AB)");
 
109
  dpd_buf4_close(&G);
 
110
  dpd_buf4_init(&G, CC_TMP0, 0, 0, 5, 0, 5, 0, "G(IJ,AB)");
 
111
  dpd_buf4_scm(&G, 0.5);
 
112
 
 
113
  for(h=0; h < nirreps; h++) {
 
114
    dpd_buf4_mat_irrep_init(&G, h);
 
115
    dpd_buf4_mat_irrep_rd(&G, h);
 
116
 
 
117
    for(row=0; row < G.params->rowtot[h]; row++) {
 
118
      p = G.params->roworb[h][row][0];
 
119
      q = G.params->roworb[h][row][1];
 
120
      for(col=0; col < G.params->coltot[h]; col++) {
 
121
        r = G.params->colorb[h][col][0];
 
122
        s = G.params->colorb[h][col][1];
 
123
 
 
124
        if((qt_occ[p] == qt_vir[r]) && (qt_occ[q] == qt_vir[s]))
 
125
          G.matrix[h][row][col] *= 2;
 
126
      }
 
127
    }
 
128
 
 
129
    dpd_buf4_mat_irrep_wrt(&G, h);
 
130
    dpd_buf4_mat_irrep_close(&G, h);
 
131
  }
 
132
  
 
133
  dpd_buf4_dump(&G, OutBuf, qt_occ, qt_occ, qt_vir, qt_vir, 0, 0);
 
134
  dpd_buf4_close(&G);
 
135
 
 
136
  dpd_buf4_init(&G, CC_GAMMA, 0, 11, 5, 11, 5, 0, "GCiAb");
 
137
  dpd_buf4_sort(&G, CC_TMP0, prqs, 5, 10, "G(ca,IB)");
 
138
  dpd_buf4_close(&G);
 
139
  dpd_buf4_init(&G, CC_TMP0, 0, 5, 10, 5, 10, 0, "G(ca,IB)");
 
140
 
 
141
  for(h=0; h < nirreps; h++) {
 
142
    dpd_buf4_mat_irrep_init(&G, h);
 
143
    dpd_buf4_mat_irrep_rd(&G, h);
 
144
 
 
145
    for(row=0; row < G.params->rowtot[h]; row++) {
 
146
      p = G.params->roworb[h][row][0];
 
147
      q = G.params->roworb[h][row][1];
 
148
      for(col=0; col < G.params->coltot[h]; col++) {
 
149
        r = G.params->colorb[h][col][0];
 
150
        s = G.params->colorb[h][col][1];
 
151
 
 
152
        if((qt_vir[p] == qt_occ[r]) && (q == s))
 
153
          G.matrix[h][row][col] *= 2;
 
154
      }
 
155
    }
 
156
 
 
157
    dpd_buf4_mat_irrep_wrt(&G, h);
 
158
    dpd_buf4_mat_irrep_close(&G, h);
 
159
  }
 
160
 
 
161
  dpd_buf4_dump(&G, OutBuf, qt_vir, qt_vir, qt_occ, qt_vir, 0, 0);
 
162
  dpd_buf4_close(&G);
 
163
  dpd_buf4_init(&G, CC_GAMMA, 0, 5, 5, 5, 5, 0, "GAbCd");
 
164
  dpd_buf4_sort(&G, CC_TMP0, prqs, 5, 5, "G(AC,BD)");
 
165
  dpd_buf4_close(&G);
 
166
  dpd_buf4_init(&G, CC_TMP0, 0, 5, 5, 5, 5, 0, "G(AC,BD)");
 
167
  dpd_buf4_dump(&G, OutBuf, qt_vir, qt_vir, qt_vir, qt_vir, 1, 0);
 
168
  dpd_buf4_close(&G);
 
169
 
 
170
  }
 
171
}
 
172
 
 
173
}} // namespace psi::ccdensity