← Back to branch summary

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

« back to all changes in this revision

Viewing changes to src/bin/psimrcc/index.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/psimrcc/index.cc
 
1
/***************************************************************************
 
2
 *  PSIMRCC : Copyright (C) 2007 by Francesco Evangelista and Andrew Simmonett
 
3
 *  frank@ccc.uga.edu   andysim@ccc.uga.edu
 
4
 *  A multireference coupled cluster code
 
5
 ***************************************************************************/
 
6
#include <iostream>
 
7
#include <algorithm>
 
8
 
 
9
#include <libmoinfo/libmoinfo.h>
 
10
 
 
11
#include "index.h"
 
12
#include "memory_manager.h"
 
13
 
 
14
extern FILE *infile, *outfile;
 
15
 
 
16
namespace psi{ namespace psimrcc{
 
17
 
 
18
using namespace std;
 
19
 
 
20
int                       CCIndex::nirreps=-1;
 
21
 
 
22
CCIndex::CCIndex(std::string str):
 
23
nelements(0),ntuples(0), greater_than_or_equal(false), greater_than(false), memory(0.0), label(str),tuples(0),
 
24
one_index_to_tuple(0),one_index_to_irrep(0),
 
25
two_index_to_tuple(0),two_index_to_irrep(0),
 
26
three_index_to_tuple(0),three_index_to_irrep(0)
 
27
{
 
28
  if(nirreps<0) nirreps = moinfo->get_nirreps();
 
29
  init();
 
30
}
 
31
 
 
32
CCIndex::~CCIndex()
 
33
{
 
34
  cleanup();
 
35
}
 
36
 
 
37
void CCIndex::init()
 
38
{
 
39
  // New orbital spaces must be added here
 
40
  for(int i =0;i<label.size();i++)
 
41
    if( label[i]=='o' || label[i]=='a' || label[i]=='v' || label[i]=='s' || label[i]=='n')
 
42
      nelements++;
 
43
 
 
44
  // Get the orbital spaces data pointers
 
45
  for(int i=0;i<label.size();i++)
 
46
  {
 
47
    if(label[i]=='o'){
 
48
      mospi.push_back(moinfo->get_occ());
 
49
      indices_to_pitzer.push_back(moinfo->get_occ_to_pitzer());
 
50
    }else if(label[i]=='v'){
 
51
      mospi.push_back(moinfo->get_vir());
 
52
      indices_to_pitzer.push_back(moinfo->get_vir_to_pitzer());
 
53
    }else if(label[i]=='a'){
 
54
      mospi.push_back(moinfo->get_actv());
 
55
      indices_to_pitzer.push_back(moinfo->get_act_to_pitzer());
 
56
    }else if(label[i]=='s'){
 
57
      mospi.push_back(moinfo->get_sopi());
 
58
      indices_to_pitzer.push_back(moinfo->get_so_to_pitzer());
 
59
    }else if(label[i]=='n'){
 
60
      mospi.push_back(moinfo->get_orbspi());
 
61
      indices_to_pitzer.push_back(moinfo->get_orbs_to_pitzer());
 
62
    }
 
63
  }
 
64
  for(int i=0;i<nelements;i++){
 
65
    first_mos.push_back(vector<int>(nirreps,0));
 
66
    dimension.push_back(0);
 
67
  }
 
68
  for(int i=0;i<nelements;i++){
 
69
    for(int h=0;h<nirreps;h++){
 
70
      first_mos[i][h]=dimension[i];
 
71
      dimension[i]+=mospi[i][h];
 
72
    }
 
73
  }
 
74
 
 
75
  switch(nelements){
 
76
    case 0:
 
77
      make_zero_index();
 
78
      break;
 
79
    case 1:
 
80
      make_one_index();
 
81
      break;
 
82
    case 2:
 
83
      make_two_index();
 
84
      break;
 
85
    case 3:
 
86
      make_three_index();
 
87
      break;
 
88
    default:{
 
89
      fprintf(outfile,"\n\n\tThe CCIndex class cannot handle %s because there are more than three indices!!!\n\n",label.c_str());
 
90
      fflush(outfile);
 
91
      exit(1);
 
92
    }
 
93
  }
 
94
  mem->add_allocated_memory(memory);
 
95
}
 
96
 
 
97
void CCIndex::cleanup()
 
98
{
 
99
  if(tuples!=0)
 
100
    free_smatrix(tuples,ntuples,dimension.size());
 
101
  if(one_index_to_tuple!=0)
 
102
    delete[] one_index_to_tuple;
 
103
  if(one_index_to_irrep!=0)
 
104
    delete[] one_index_to_irrep;
 
105
  if(two_index_to_tuple!=0)
 
106
    release2(two_index_to_tuple);
 
107
  if(two_index_to_irrep!=0)
 
108
    release2(two_index_to_irrep);
 
109
  if(three_index_to_tuple!=0)
 
110
    release3(three_index_to_tuple);
 
111
  if(three_index_to_irrep!=0)
 
112
    release3(three_index_to_irrep);
 
113
}
 
114
 
 
115
void CCIndex::make_zero_index()
 
116
{
 
117
  std::vector<std::vector<short> >  pairs;                  // The pairs ordered as a vector
 
118
  ntuples = 0;
 
119
  for(int h=0;h<nirreps;h++){
 
120
    first.push_back(ntuples);
 
121
    if(h==0){
 
122
      std::vector<short> pair;
 
123
      pairs.push_back(pair);
 
124
      ntuples++;
 
125
    }
 
126
    last.push_back(ntuples);
 
127
    pairpi.push_back(last[h]-first[h]);
 
128
  }
 
129
  // Allocate the memory for the tuples and store them
 
130
  memory+=(double)init_smatrix(tuples,1,1)/1048576.0;
 
131
  tuples[0][0] = 0;
 
132
}
 
133
 
 
134
void CCIndex::make_one_index()
 
135
{
 
136
  // The pairs ordered as a vector
 
137
  std::vector<std::vector<short> >  pairs;
 
138
 
 
139
  // Allocate the 1->tuple mapping array and set them to -1
 
140
  one_index_to_tuple = new size_t[dimension[0]];
 
141
  one_index_to_irrep = new int[dimension[0]];
 
142
  for(int i=0;i<dimension[0];i++){
 
143
    one_index_to_tuple[i] =  0;
 
144
    one_index_to_irrep[i] = -1;
 
145
  }
 
146
  memory+=(double)dimension[0]*2.0*sizeof(int)/1048576.0;
 
147
 
 
148
  ntuples = 0;
 
149
  for(int h=0;h<nirreps;h++){
 
150
    first.push_back(ntuples);
 
151
    for(int p=0;p<mospi[0][h];p++){
 
152
      one_index_to_tuple[ntuples]=p;
 
153
      one_index_to_irrep[ntuples]=h;
 
154
      std::vector<short> pair;
 
155
      pair.push_back(ntuples);
 
156
      pairs.push_back(pair);
 
157
      ntuples++;
 
158
    }
 
159
    last.push_back(ntuples);
 
160
    pairpi.push_back(last[h]-first[h]);
 
161
  }
 
162
 
 
163
  // Allocate the memory for the tuples and store them
 
164
  memory+=(double)init_smatrix(tuples,ntuples,1)/1048576.0;
 
165
  for(int n=0;n<pairs.size();n++)
 
166
    tuples[n][0] = pairs[n][0];
 
167
}
 
168
 
 
169
void CCIndex::make_two_index()
 
170
{
 
171
  std::vector<std::vector<short> >  pairs;                  // The pairs ordered as a vector
 
172
 
 
173
  // Allocate the 2->tuple mapping array and set them to -1
 
174
  allocate2(size_t,two_index_to_tuple,dimension[0],dimension[1]);
 
175
  allocate2(int,two_index_to_irrep,dimension[0],dimension[1]);
 
176
  
 
177
  for(int i=0;i<dimension[0];i++){
 
178
    for(int j=0;j<dimension[1];j++){
 
179
      two_index_to_tuple[i][j] = 0;
 
180
      two_index_to_irrep[i][j] = -1;
 
181
    }
 
182
  }
 
183
 
 
184
  // [X>=Y]
 
185
  if(label.find(">=")!=string::npos){
 
186
    greater_than_or_equal = true;
 
187
    ntuples               = 0;
 
188
    for(int h=0;h<nirreps;h++){
 
189
      first.push_back(ntuples);
 
190
      for(int p_sym=0;p_sym<nirreps;p_sym++){
 
191
        int q_sym = h^p_sym;
 
192
        int p = first_mos[0][p_sym];
 
193
        for(int p_rel=0;p_rel<mospi[0][p_sym];p_rel++){
 
194
          int q = first_mos[1][q_sym];
 
195
          for(int q_rel=0;q_rel<mospi[1][q_sym];q_rel++){
 
196
            if(p>=q){
 
197
              two_index_to_tuple[p][q]=ntuples-first[h];
 
198
              two_index_to_irrep[p][q]=h;
 
199
              std::vector<short> pair;
 
200
              pair.push_back(p);
 
201
              pair.push_back(q);
 
202
              pairs.push_back(pair);
 
203
              ntuples++;
 
204
            }
 
205
            q++;
 
206
          }
 
207
          p++;
 
208
        }
 
209
      }
 
210
      last.push_back(ntuples);
 
211
      pairpi.push_back(last[h]-first[h]);
 
212
    }
 
213
  }else if(label.find(">")!=string::npos){
 
214
    greater_than          = true;
 
215
    ntuples               = 0;
 
216
    for(int h=0;h<nirreps;h++){
 
217
      first.push_back(ntuples);
 
218
      for(int p_sym=0;p_sym<nirreps;p_sym++){
 
219
        int q_sym = h^p_sym;
 
220
        int p = first_mos[0][p_sym];
 
221
        for(int p_rel=0;p_rel<mospi[0][p_sym];p_rel++){
 
222
          int q = first_mos[1][q_sym];
 
223
          for(int q_rel=0;q_rel<mospi[1][q_sym];q_rel++){
 
224
            if(p>q){
 
225
              two_index_to_tuple[p][q]=ntuples-first[h];
 
226
              two_index_to_irrep[p][q]=h;
 
227
              std::vector<short> pair;
 
228
              pair.push_back(p);
 
229
              pair.push_back(q);
 
230
              pairs.push_back(pair);
 
231
              ntuples++;
 
232
            }
 
233
            q++;
 
234
          }
 
235
          p++;
 
236
        }
 
237
      }
 
238
      last.push_back(ntuples);
 
239
      pairpi.push_back(last[h]-first[h]);
 
240
    }
 
241
  }else{
 
242
    ntuples               = 0;
 
243
    for(int h=0;h<nirreps;h++){
 
244
      first.push_back(ntuples);
 
245
      for(int p_sym=0;p_sym<nirreps;p_sym++){
 
246
        int q_sym = h^p_sym;
 
247
        int p = first_mos[0][p_sym];
 
248
        for(int p_rel=0;p_rel<mospi[0][p_sym];p_rel++){
 
249
          int q = first_mos[1][q_sym];
 
250
          for(int q_rel=0;q_rel<mospi[1][q_sym];q_rel++){
 
251
            two_index_to_tuple[p][q]=ntuples-first[h];
 
252
            two_index_to_irrep[p][q]=h;
 
253
            std::vector<short> pair;
 
254
            pair.push_back(p);
 
255
            pair.push_back(q);
 
256
            pairs.push_back(pair);
 
257
            ntuples++;
 
258
            q++;
 
259
          }
 
260
          p++;
 
261
        }
 
262
      }
 
263
      last.push_back(ntuples);
 
264
      pairpi.push_back(last[h]-first[h]);
 
265
    }
 
266
  }
 
267
 
 
268
  // Allocate the memory for the tuples and store them
 
269
  memory+=(double)init_smatrix(tuples,ntuples,2)/1048576.0;
 
270
  for(int n=0;n<pairs.size();n++){
 
271
    tuples[n][0] = pairs[n][0];
 
272
    tuples[n][1] = pairs[n][1];
 
273
  }
 
274
}
 
275
 
 
276
 
 
277
void CCIndex::make_three_index()
 
278
{
 
279
  if(label.find(">")!=string::npos){
 
280
    fprintf(outfile,"\n\n\tThe CCIndex class cannot handle restricted loops for triplets!!!\n\n");
 
281
    fflush(outfile);
 
282
    exit(1);
 
283
  }
 
284
 
 
285
  std::vector<std::vector<short> >  pairs;                  // The pairs ordered as a vector
 
286
 
 
287
  // Allocate the 3->tuple mapping array and set them to -1
 
288
  allocate3(size_t,three_index_to_tuple,dimension[0],dimension[1],dimension[2]);
 
289
  allocate3(int,three_index_to_irrep,dimension[0],dimension[1],dimension[2]);
 
290
  for(int i=0;i<dimension[0];i++){
 
291
    for(int j=0;j<dimension[1];j++){
 
292
      for(int k=0;k<dimension[2];k++){
 
293
        three_index_to_tuple[i][j][k] =  0;
 
294
        three_index_to_irrep[i][j][k] = -1;
 
295
      }
 
296
    }
 
297
  }
 
298
 
 
299
 
 
300
  ntuples  = 0;
 
301
  for(int h=0;h<nirreps;h++){
 
302
    first.push_back(ntuples);
 
303
    for(int p_sym=0;p_sym<nirreps;p_sym++){
 
304
      for(int q_sym=0;q_sym<nirreps;q_sym++){
 
305
        int r_sym = h^p_sym^q_sym;
 
306
        int p = first_mos[0][p_sym];
 
307
        for(int p_rel=0;p_rel<mospi[0][p_sym];p_rel++){
 
308
          int q = first_mos[1][q_sym];
 
309
          for(int q_rel=0;q_rel<mospi[1][q_sym];q_rel++){
 
310
            int r = first_mos[2][r_sym];
 
311
            for(int r_rel=0;r_rel<mospi[2][r_sym];r_rel++){
 
312
              three_index_to_tuple[p][q][r]=ntuples-first[h];
 
313
              three_index_to_irrep[p][q][r]=h;
 
314
              std::vector<short> pair;
 
315
              pair.push_back(p);
 
316
              pair.push_back(q);
 
317
              pair.push_back(r);
 
318
              pairs.push_back(pair);
 
319
              ntuples++;
 
320
              r++;
 
321
            }
 
322
          q++;
 
323
          }
 
324
        p++;
 
325
        }
 
326
      }
 
327
    }
 
328
    last.push_back(ntuples);
 
329
    pairpi.push_back(last[h]-first[h]);
 
330
  }
 
331
 
 
332
  // Allocate the memory for the tuples and store them
 
333
  memory+=(double)init_smatrix(tuples,ntuples,3)/1048576.0;
 
334
  for(int n=0;n<pairs.size();n++){
 
335
    tuples[n][0] = pairs[n][0];
 
336
    tuples[n][1] = pairs[n][1];
 
337
    tuples[n][2] = pairs[n][2];
 
338
  }
 
339
}
 
340
 
 
341
void CCIndex::print()
 
342
{
 
343
  fprintf(outfile,"\n\n---------------------------------");
 
344
  fprintf(outfile,"\n\tPair Type %s has %d elements",label.c_str(),ntuples);
 
345
  fprintf(outfile,"\n---------------------------------");
 
346
  int index=0;
 
347
  for(int h=0;h<nirreps;h++){
 
348
    if(pairpi[h]>0)
 
349
      fprintf(outfile,"\n\t%s",moinfo->get_irr_labs(h));
 
350
    for(int tuple=0;tuple<pairpi[h];tuple++){
 
351
      fprintf(outfile,"\n\t\t( ");
 
352
      for(int k=0;k<nelements;k++)
 
353
        fprintf(outfile,"%d ",tuples[index][k]);
 
354
      fprintf(outfile,")");
 
355
      index++;
 
356
    }
 
357
  }
 
358
  fprintf(outfile,"\n---------------------------------");
 
359
}
 
360
 
 
361
}} /* End Namespaces */