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

« back to all changes in this revision

Viewing changes to src/bin/cclambda/get_params.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 <stdlib.h>
 
3
#include <string.h>
3
4
#include <math.h>
4
5
#include <libipv1/ip_lib.h>
5
6
#include <libciomr/libciomr.h>
 
7
#include <psifiles.h>
 
8
#include <libqt/qt.h>
6
9
#define EXTERN
7
10
#include "globals.h"
8
11
 
9
12
void get_params(void)
10
13
{
11
 
  int errcod, iconv,i,j;
 
14
  int errcod, iconv,i,j,k,l,prop_sym,prop_root, excited_method=0;
 
15
        int *states_per_irrep, prop_all, lambda_and_Ls = 0;
12
16
  char lbl[32];
13
 
 
14
 
  params.maxiter = 50;
15
 
  errcod = ip_data("MAXITER","%d",&(params.maxiter),0);
 
17
  char *junk;
 
18
 
 
19
  /* check WFN keyword in input */
 
20
  errcod = ip_string("WFN", &(params.wfn), 0);
 
21
        excited_method = cc_excited(params.wfn);
 
22
 
 
23
  if(!strcmp(params.wfn,"CC2") || !strcmp(params.wfn,"EOM_CC2")) {
 
24
    psio_read_entry(CC_INFO, "CC2 Energy", (char *) &(moinfo.ecc),
 
25
                    sizeof(double));
 
26
    fprintf(outfile,  "\tCC2 energy          (CC_INFO) = %20.15f\n",moinfo.ecc);
 
27
    fprintf(outfile,  "\tTotal CC2 energy    (CC_INFO) = %20.15f\n",
 
28
            moinfo.eref+moinfo.ecc);
 
29
  }
 
30
  else if(!strcmp(params.wfn,"CCSD") || !strcmp(params.wfn,"EOM_CCSD")) {
 
31
    psio_read_entry(CC_INFO, "CCSD Energy", (char *) &(moinfo.ecc),
 
32
                    sizeof(double));
 
33
    fprintf(outfile,  "\tCCSD energy         (CC_INFO) = %20.15f\n",moinfo.ecc);
 
34
    fprintf(outfile,  "\tTotal CCSD energy   (CC_INFO) = %20.15f\n",
 
35
            moinfo.eref+moinfo.ecc);
 
36
  }
 
37
  else if(!strcmp(params.wfn,"CC3") || !strcmp(params.wfn,"EOM_CC3")) {
 
38
    psio_read_entry(CC_INFO, "CC3 Energy", (char *) &(moinfo.ecc),
 
39
                    sizeof(double));
 
40
    fprintf(outfile,  "\tCC3 energy          (CC_INFO) = %20.15f\n",moinfo.ecc);
 
41
    fprintf(outfile,  "\tTotal CC3 energy    (CC_INFO) = %20.15f\n",
 
42
            moinfo.eref+moinfo.ecc);
 
43
  }
 
44
 
 
45
  /* read in the easy-to-understand parameters */
 
46
 
16
47
  params.convergence = 1e-7;
17
48
  errcod = ip_data("CONVERGENCE","%d",&(iconv),0);
18
49
  if(errcod == IPE_OK) params.convergence = 1.0*pow(10.0,(double) -iconv);
 
50
 
19
51
  params.restart = 1;
20
52
  errcod = ip_boolean("RESTART", &(params.restart),0);
21
 
  /* If the MO orbital phases are screwed up, don't restart */
22
53
  if(!moinfo.phase) params.restart = 0;
23
54
 
24
55
  fndcor(&(params.memory),infile,outfile);
25
56
 
 
57
  params.print = 0;
 
58
  errcod = ip_data("PRINT", "%d", &(params.print),0);
 
59
 
26
60
  params.cachelev = 2;
27
61
  errcod = ip_data("CACHELEV", "%d", &(params.cachelev),0);
28
62
 
 
63
  params.sekino = 0;
 
64
  if(ip_exist("SEKINO",0)) {
 
65
    errcod = ip_boolean("SEKINO", &params.sekino, 0);
 
66
    if(errcod != IPE_OK) params.sekino = 0;
 
67
  }
 
68
 
 
69
  params.diis = 1;
 
70
  errcod = ip_boolean("DIIS", &params.diis, 0);
 
71
 
29
72
  params.aobasis = 0;
30
73
  errcod = ip_boolean("AO_BASIS", &(params.aobasis),0);
31
74
  params.aobasis = 0;  /* AO basis code not yet working for lambda */
32
75
 
33
 
  /* determine number of left-hand states per irrep */
34
 
  params.states_per_irrep = (int *) malloc(moinfo.nirreps * sizeof(int));
35
 
  for (i=0;i<moinfo.nirreps;++i) params.states_per_irrep[i] = 0;
36
 
 
37
 
  if (params.ground) { /* for ground state only find 1 A1 lambda */
38
 
    params.states_per_irrep[moinfo.sym] = 1;
39
 
  }
40
 
  else { /* excited state */
41
 
    if (ip_exist("STATES_PER_IRREP",0)) {
42
 
      ip_count("STATES_PER_IRREP", &i, 0);
43
 
      if (i != moinfo.nirreps) {
44
 
        fprintf(outfile,"Dim. of states_per_irrep vector must be %d\n", moinfo.nirreps) ;
45
 
        exit(0);
46
 
      }
47
 
      for (i=0;i<moinfo.nirreps;++i)
48
 
        errcod = ip_data("STATES_PER_IRREP","%d",&(params.states_per_irrep[i]),1,i);
49
 
    }
50
 
    else { fprintf(outfile,"Must have states_per_irrep vector in input.\n"); exit(0); }
51
 
  }
52
 
 
53
 
  /* determine Ls per irrep */
54
 
  params.Ls_per_irrep = (int *) malloc(moinfo.nirreps * sizeof(int));
55
 
  for (i=0;i<moinfo.nirreps;++i) {
56
 
    params.Ls_per_irrep[i^moinfo.sym] =  params.states_per_irrep[i];
57
 
  }
58
 
 
 
76
  if(ip_exist("ABCD",0)) {
 
77
    errcod = ip_string("ABCD", &(params.abcd), 0);
 
78
    if(strcmp(params.abcd,"NEW") && strcmp(params.abcd,"OLD")) {
 
79
      fprintf(outfile, "Invalid ABCD algorithm: %s\n", params.abcd);
 
80
      exit(PSI_RETURN_FAILURE);
 
81
    }
 
82
  }
 
83
  else params.abcd = strdup("NEW");
 
84
 
 
85
  params.num_amps = 10;
 
86
  if(ip_exist("NUM_AMPS",0)) {
 
87
    errcod = ip_data("NUM_AMPS", "%d", &(params.num_amps), 0);
 
88
  }
 
89
 
 
90
  /* Determine DERTYPE */
 
91
  params.dertype = 0;
 
92
  if(ip_exist("DERTYPE",0)) {
 
93
    errcod = ip_string("DERTYPE", &(junk),0);
 
94
    if(errcod != IPE_OK) {
 
95
      printf("Invalid value of input keyword DERTYPE: %s\n", junk);
 
96
      exit(PSI_RETURN_FAILURE); 
 
97
                }
 
98
    else if(!strcmp(junk,"NONE")) params.dertype = 0;
 
99
    else if(!strcmp(junk,"FIRST")) params.dertype = 1;
 
100
    else if(!strcmp(junk,"RESPONSE")) params.dertype = 3; /* linear response */
 
101
    else {
 
102
      printf("Invalid value of input keyword DERTYPE: %s\n", junk);
 
103
      exit(PSI_RETURN_FAILURE); 
 
104
    }
 
105
    free(junk);
 
106
  }
 
107
        else { /* DERTYPE is absent, assume 1 if jobtype=opt; 0 if jobtype=oeprop */
 
108
    ip_string("JOBTYPE", &(junk),0);
 
109
    if(!strcmp(junk,"OEPROP")) params.dertype = 0;
 
110
    else if(!strcmp(junk,"OPT")) params.dertype = 1;
 
111
    else {
 
112
      printf("Don't know what to do with DERTYPE missing and jobtype: %s\n", junk);
 
113
      exit(PSI_RETURN_FAILURE); 
 
114
    }
 
115
    free(junk);
 
116
        }
 
117
 
 
118
  /* begin local parameters */
59
119
  params.local = 0;
60
120
  errcod = ip_boolean("LOCAL", &(params.local),0);
61
121
  local.cutoff = 0.02;
62
122
  errcod = ip_data("LOCAL_CUTOFF", "%lf", &(local.cutoff), 0);
63
 
 
64
123
  if(ip_exist("LOCAL_METHOD",0)) {
65
124
    errcod = ip_string("LOCAL_METHOD", &(local.method), 0);
66
125
    if(strcmp(local.method,"AOBASIS") && strcmp(local.method,"WERNER")) {
82
141
  }
83
142
  else if(params.local) {
84
143
    local.weakp = (char *) malloc(4 * sizeof(char));
85
 
    sprintf(local.weakp, "%s", "MP2");
 
144
    sprintf(local.weakp, "%s", "NONE");
86
145
  }
87
146
  
88
 
  local.filter_singles = 1;
 
147
  if(params.dertype == 3)
 
148
    local.filter_singles = 0;
 
149
  else
 
150
    local.filter_singles = 1;
89
151
  ip_boolean("LOCAL_FILTER_SINGLES", &(local.filter_singles), 0);
90
152
 
 
153
  local.cphf_cutoff = 0.10;
 
154
  ip_data("LOCAL_CPHF_CUTOFF", "%lf", &(local.cphf_cutoff), 0);
 
155
 
 
156
  local.freeze_core = NULL;
 
157
  ip_string("FREEZE_CORE", &local.freeze_core, 0);
 
158
  if(local.freeze_core == NULL) local.freeze_core = strdup("FALSE");
 
159
 
 
160
  if(ip_exist("LOCAL_PAIRDEF",0)){
 
161
    errcod = ip_string("LOCAL_PAIRDEF", &(local.pairdef), 0);
 
162
    if(strcmp(local.pairdef,"BP") && strcmp(local.pairdef,"RESPONSE")) {
 
163
      fprintf(outfile, "Invalid keyword for strong/weak pair definition: %s\n", local.pairdef);
 
164
      exit(PSI_RETURN_FAILURE);
 
165
    }
 
166
  }
 
167
  else if(params.local && params.dertype == 3)
 
168
    local.pairdef = strdup("RESPONSE");
 
169
  else if(params.local)
 
170
    local.pairdef = strdup("BP");
 
171
 
 
172
        /* Now setup the structure which determines what will be solved */ 
 
173
        /* if --zeta, use Xi and solve for Zeta */
 
174
        /* if (DERTYPE == FIRST) determine ground vs. excited from wfn.
 
175
            if ground, do only lambda.
 
176
                        if excited, compute only one L chosen as described below.
 
177
                        */
 
178
        /* if (DERTYPE == RESPONSE), determine ground vs. excited from wfn.
 
179
            Compute lambda.
 
180
                        if excited, also do L(s) chosen as described below */
 
181
        /* if (DERTYPE == NONE) determine ground vs. excited from wfn.
 
182
            Compute lambda.
 
183
                        if excited, also do L(s) chosen as described below */
 
184
/* To determine which L(s) to compute for multiple L(s):
 
185
          Check PROP_ALL in input
 
186
                 - If (PROP_ALL == true), compute L for all excited states.
 
187
                 - If false, check PROP_SYM for irrep desired, and PROP_ROOT
 
188
                             for root desired, as in cceom. */
 
189
/* To determine which L(s) to compute for single L(s) 
 
190
                 - Check PROP_SYM for irrep desired, and PROP_ROOT
 
191
                             for root desired, as in cceom. */
 
192
 
 
193
  /* setup property variables for excited states */
 
194
  if (cc_excited(params.wfn)) {
 
195
    ip_count("STATES_PER_IRREP", &i, 0);
 
196
          states_per_irrep = (int *) malloc(moinfo.nirreps * sizeof(int));
 
197
    for (i=0;i<moinfo.nirreps;++i)
 
198
      errcod = ip_data("STATES_PER_IRREP","%d",&(states_per_irrep[i]),1,i);
 
199
 
 
200
          prop_all = 0;
 
201
          if (ip_exist("PROP_ALL",0)) ip_boolean("PROP_ALL",&prop_all,0);
 
202
                /* command-line overrides this keyword (at least for now) */
 
203
                if (params.all) prop_all = 1;
 
204
 
 
205
          if (ip_exist("PROP_SYM",0)) {  /* read symmetry of state for properties */
 
206
      ip_data("PROP_SYM","%d",&(prop_sym),0);
 
207
      prop_sym -= 1;
 
208
      prop_sym = moinfo.sym^prop_sym;
 
209
          }
 
210
    else { /* just use last irrep of states requested for symmetry of states */
 
211
      for (i=0;i<moinfo.nirreps;++i) {
 
212
                    if (states_per_irrep[i] > 0)
 
213
          prop_sym = i^moinfo.sym;
 
214
                  }
 
215
    }
 
216
 
 
217
    if (ip_exist("PROP_ROOT",0)) { /* read prop_root */
 
218
      ip_data("PROP_ROOT","%d",&(prop_root),0);
 
219
      prop_root -= 1;
 
220
          }
 
221
          else { /* just use highest root, if you need only one of them */
 
222
            prop_root = states_per_irrep[prop_sym^moinfo.sym];
 
223
                        prop_root -= 1;
 
224
          }
 
225
        }
 
226
        
 
227
 
 
228
  if (params.zeta) { /* only use Xi to solve for Zeta */
 
229
          params.nstates = 1;
 
230
    pL_params = (struct L_Params *) malloc(params.nstates * sizeof(struct L_Params));
 
231
    psio_read_entry(CC_INFO, "XI Irrep", (char *) &i,sizeof(int));
 
232
    fprintf(outfile,"\tIrrep of Zeta       (CC_INFO) = %d\n", i);
 
233
    pL_params[0].irrep = prop_sym = i; /* is this always A1? I forget */
 
234
    pL_params[0].root = prop_root = 0;
 
235
    pL_params[0].ground = 0;
 
236
    pL_params[0].cceom_energy = 0.0;
 
237
    pL_params[0].R0 = 0.0; /* <Zeta0|R0> = 0, since zeta_0 = 0 */
 
238
    sprintf(pL_params[0].L1A_lbl,"ZIA");
 
239
    sprintf(pL_params[0].L1B_lbl,"Zia");
 
240
    sprintf(pL_params[0].L2AA_lbl,"ZIJAB");
 
241
    sprintf(pL_params[0].L2BB_lbl,"Zijab");
 
242
    sprintf(pL_params[0].L2AB_lbl,"ZIjAb");
 
243
    sprintf(pL_params[0].L2RHF_lbl,"2ZIjAb - ZIjbA");
 
244
  }
 
245
        else if (params.dertype == 1) { /* analytic gradient, ignore prop_all */
 
246
          if (!cc_excited(params.wfn)) { /* do only lambda for ground state */
 
247
            params.nstates = 1;
 
248
      pL_params = (struct L_Params *) malloc(params.nstates * sizeof(struct L_Params));
 
249
      pL_params[0].irrep = 0;
 
250
      pL_params[0].root = -1;
 
251
      pL_params[0].ground = 1;
 
252
      pL_params[0].R0 = 1.0;
 
253
      pL_params[0].cceom_energy = 0.0;
 
254
      sprintf(pL_params[0].L1A_lbl,"LIA %d %d",0, -1);
 
255
      sprintf(pL_params[0].L1B_lbl,"Lia %d %d",0, -1);
 
256
      sprintf(pL_params[0].L2AA_lbl,"LIJAB %d %d",0, -1);
 
257
      sprintf(pL_params[0].L2BB_lbl,"Lijab %d %d",0, -1);
 
258
      sprintf(pL_params[0].L2AB_lbl,"LIjAb %d %d",0, -1);
 
259
      sprintf(pL_params[0].L2RHF_lbl,"2LIjAb - LIjbA %d %d",0, -1);
 
260
                }
 
261
                else { /* do only one L for excited state */
 
262
                  params.nstates = 1;
 
263
      pL_params = (struct L_Params *) malloc(params.nstates * sizeof(struct L_Params));
 
264
      pL_params[0].irrep = prop_sym;
 
265
      pL_params[0].root = prop_root;
 
266
      pL_params[0].ground = 0;
 
267
      if(!strcmp(params.wfn,"EOM_CC2")) {
 
268
        sprintf(lbl,"EOM CC2 Energy for root %d %d", prop_sym, prop_root);
 
269
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[0].cceom_energy),sizeof(double));
 
270
        sprintf(lbl,"EOM CC2 R0 for root %d %d", prop_sym, prop_root);
 
271
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[0].R0),sizeof(double));
 
272
      }
 
273
      else if(!strcmp(params.wfn,"EOM_CCSD")) {
 
274
        sprintf(lbl,"EOM CCSD Energy for root %d %d", prop_sym, prop_root);
 
275
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[0].cceom_energy),sizeof(double));
 
276
        sprintf(lbl,"EOM CCSD R0 for root %d %d", prop_sym, prop_root);
 
277
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[0].R0),sizeof(double));
 
278
      }
 
279
      else if(!strcmp(params.wfn,"EOM_CC3")) {
 
280
        sprintf(lbl,"EOM CC3 Energy for root %d %d", prop_sym, prop_root);
 
281
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[0].cceom_energy),sizeof(double));
 
282
        sprintf(lbl,"EOM CC3 R0 for root %d %d", prop_sym, prop_root);
 
283
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[0].R0),sizeof(double));
 
284
      }
 
285
      sprintf(pL_params[0].L1A_lbl,"LIA %d %d",prop_sym, prop_root);
 
286
      sprintf(pL_params[0].L1B_lbl,"Lia %d %d",prop_sym, prop_root);
 
287
      sprintf(pL_params[0].L2AA_lbl,"LIJAB %d %d",prop_sym, prop_root);
 
288
      sprintf(pL_params[0].L2BB_lbl,"Lijab %d %d",prop_sym, prop_root);
 
289
      sprintf(pL_params[0].L2AB_lbl,"LIjAb %d %d",prop_sym, prop_root);
 
290
      sprintf(pL_params[0].L2RHF_lbl,"2LIjAb - LIjbA %d %d",prop_sym, prop_root);
 
291
                }
 
292
  }
 
293
        else if (params.dertype == 3) { /* response calculation */
 
294
          if (!cc_excited(params.wfn)) { /* ground state */
 
295
            params.nstates = 1;
 
296
      pL_params = (struct L_Params *) malloc(params.nstates * sizeof(struct L_Params));
 
297
      pL_params[0].irrep = 0;
 
298
      pL_params[0].root = -1;
 
299
      pL_params[0].ground = 1;
 
300
      pL_params[0].R0 = 1.0;
 
301
      pL_params[0].cceom_energy = 0.0;
 
302
      sprintf(pL_params[0].L1A_lbl,"LIA %d %d",0, -1);
 
303
      sprintf(pL_params[0].L1B_lbl,"Lia %d %d",0, -1);
 
304
      sprintf(pL_params[0].L2AA_lbl,"LIJAB %d %d",0, -1);
 
305
      sprintf(pL_params[0].L2BB_lbl,"Lijab %d %d",0, -1);
 
306
      sprintf(pL_params[0].L2AB_lbl,"LIjAb %d %d",0, -1);
 
307
      sprintf(pL_params[0].L2RHF_lbl,"2LIjAb - LIjbA %d %d",0, -1);
 
308
                }
 
309
                else { /* excited state */
 
310
                  lambda_and_Ls = 1; /* code is below */
 
311
                }
 
312
        }
 
313
        else if (params.dertype == 0) {
 
314
          if (!cc_excited(params.wfn)) { /* ground state */
 
315
            params.nstates = 1;
 
316
      pL_params = (struct L_Params *) malloc(params.nstates * sizeof(struct L_Params));
 
317
      pL_params[0].irrep = 0;
 
318
      pL_params[0].root = -1;
 
319
      pL_params[0].ground = 1;
 
320
      pL_params[0].R0 = 1.0;
 
321
      pL_params[0].cceom_energy = 0.0;
 
322
      sprintf(pL_params[0].L1A_lbl,"LIA %d %d",0, -1);
 
323
      sprintf(pL_params[0].L1B_lbl,"Lia %d %d",0, -1);
 
324
      sprintf(pL_params[0].L2AA_lbl,"LIJAB %d %d",0, -1);
 
325
      sprintf(pL_params[0].L2BB_lbl,"Lijab %d %d",0, -1);
 
326
      sprintf(pL_params[0].L2AB_lbl,"LIjAb %d %d",0, -1);
 
327
      sprintf(pL_params[0].L2RHF_lbl,"2LIjAb - LIjbA %d %d",0, -1);
 
328
                }
 
329
                else { /* excited state */
 
330
                  lambda_and_Ls = 1; /* code is below */
 
331
                }
 
332
        }
 
333
 
 
334
 
 
335
  /* do lambda for ground state AND do L(s) for excited states */
 
336
  if (lambda_and_Ls) {
 
337
    /* determine number of states to converge */
 
338
          params.nstates = 1; /* for ground state */
 
339
                if (prop_all) {
 
340
                  for (i=0; i<moinfo.nirreps; ++i) 
 
341
                          params.nstates += states_per_irrep[i]; /* do all L(s) */
 
342
                }
 
343
                else {
 
344
                  params.nstates += 1; /* do only one L */
 
345
                }
 
346
 
 
347
    pL_params = (struct L_Params *) malloc(params.nstates * sizeof(struct L_Params));
 
348
 
 
349
                /* ground state */
 
350
    pL_params[0].irrep = 0;
 
351
    pL_params[0].root = -1;
 
352
    pL_params[0].ground = 1;
 
353
    pL_params[0].R0 = 1.0;
 
354
    pL_params[0].cceom_energy = 0.0;
 
355
    sprintf(pL_params[0].L1A_lbl,"LIA %d %d",0, -1);
 
356
    sprintf(pL_params[0].L1B_lbl,"Lia %d %d",0, -1);
 
357
    sprintf(pL_params[0].L2AA_lbl,"LIJAB %d %d",0, -1);
 
358
    sprintf(pL_params[0].L2BB_lbl,"Lijab %d %d",0, -1);
 
359
    sprintf(pL_params[0].L2AB_lbl,"LIjAb %d %d",0, -1);
 
360
    sprintf(pL_params[0].L2RHF_lbl,"2LIjAb - LIjbA %d %d",0, -1);
 
361
 
 
362
                if (prop_all) { /* do all L(s) */
 
363
                  k=1;
 
364
                    for (i=0; i<moinfo.nirreps; ++i)  { /* look over irrep of L(s) */
 
365
                                  for (j=0; j < states_per_irrep[i^moinfo.sym]; ++j) {
 
366
            pL_params[k].irrep = i;
 
367
            pL_params[k].root = j;
 
368
            pL_params[k].ground = 0;
 
369
 
 
370
            if(!strcmp(params.wfn,"EOM_CC2")) {
 
371
              sprintf(lbl,"EOM CC2 Energy for root %d %d", i, j);
 
372
              psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[k].cceom_energy),sizeof(double));
 
373
              sprintf(lbl,"EOM CC2 R0 for root %d %d", i, j);
 
374
              psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[k].R0),sizeof(double));
 
375
            }
 
376
            else if(!strcmp(params.wfn,"EOM_CCSD")) {
 
377
              sprintf(lbl,"EOM CCSD Energy for root %d %d", i, j);
 
378
              psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[k].cceom_energy),sizeof(double));
 
379
              sprintf(lbl,"EOM CCSD R0 for root %d %d", i, j);
 
380
              psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[k].R0),sizeof(double));
 
381
            }
 
382
            else if(!strcmp(params.wfn,"EOM_CC3")) {
 
383
              sprintf(lbl,"EOM CC3 Energy for root %d %d", i, j);
 
384
              psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[k].cceom_energy),sizeof(double));
 
385
              sprintf(lbl,"EOM CC3 R0 for root %d %d", i, j);
 
386
              psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[k].R0),sizeof(double));
 
387
            }
 
388
 
 
389
            sprintf(pL_params[k].L1A_lbl,"LIA %d %d",i, j);
 
390
            sprintf(pL_params[k].L1B_lbl,"Lia %d %d",i, j);
 
391
            sprintf(pL_params[k].L2AA_lbl,"LIJAB %d %d",i, j);
 
392
            sprintf(pL_params[k].L2BB_lbl,"Lijab %d %d",i, j);
 
393
            sprintf(pL_params[k].L2AB_lbl,"LIjAb %d %d",i, j);
 
394
            sprintf(pL_params[k].L2RHF_lbl,"2LIjAb - LIjbA %d %d",i, j);
 
395
                                                k++;
 
396
                                        }
 
397
                          }
 
398
                }
 
399
                else { /* use prop_sym and prop_root determined above from input or inferrence */
 
400
      pL_params[1].irrep = prop_sym;
 
401
      pL_params[1].root = prop_root;
 
402
      pL_params[1].ground = 0;
 
403
 
 
404
      if(!strcmp(params.wfn,"EOM_CC2")) {
 
405
        sprintf(lbl,"EOM CC2 Energy for root %d %d", prop_sym, prop_root);
 
406
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[1].cceom_energy),sizeof(double));
 
407
        sprintf(lbl,"EOM CC2 R0 for root %d %d", prop_sym, prop_root);
 
408
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[1].R0),sizeof(double));
 
409
      }
 
410
      else if(!strcmp(params.wfn,"EOM_CCSD")) {
 
411
        sprintf(lbl,"EOM CCSD Energy for root %d %d", prop_sym, prop_root);
 
412
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[1].cceom_energy),sizeof(double));
 
413
        sprintf(lbl,"EOM CCSD R0 for root %d %d", prop_sym, prop_root);
 
414
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[1].R0),sizeof(double));
 
415
      }
 
416
      else if(!strcmp(params.wfn,"EOM_CC3")) {
 
417
        sprintf(lbl,"EOM CC3 Energy for root %d %d", prop_sym, prop_root);
 
418
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[1].cceom_energy),sizeof(double));
 
419
        sprintf(lbl,"EOM CC3 R0 for root %d %d", prop_sym, prop_root);
 
420
        psio_read_entry(CC_INFO, lbl, (char *) &(pL_params[1].R0),sizeof(double));
 
421
      }
 
422
 
 
423
      sprintf(pL_params[1].L1A_lbl,"LIA %d %d", prop_sym, prop_root);
 
424
      sprintf(pL_params[1].L1B_lbl,"Lia %d %d", prop_sym, prop_root);
 
425
      sprintf(pL_params[1].L2AA_lbl,"LIJAB %d %d", prop_sym, prop_root);
 
426
      sprintf(pL_params[1].L2BB_lbl,"Lijab %d %d", prop_sym, prop_root);
 
427
      sprintf(pL_params[1].L2AB_lbl,"LIjAb %d %d", prop_sym, prop_root);
 
428
      sprintf(pL_params[1].L2RHF_lbl,"2LIjAb - LIjbA %d %d", prop_sym, prop_root);
 
429
                }
 
430
  }
 
431
 
 
432
 
 
433
  params.maxiter = 50 * params.nstates;
 
434
  errcod = ip_data("MAXITER","%d",&(params.maxiter),0);
 
435
 
91
436
  fprintf(outfile, "\n\tInput parameters:\n");
92
437
  fprintf(outfile, "\t-----------------\n");
93
 
  fprintf(outfile, "\tMaxiter     =    %4d\n", params.maxiter);
94
 
  fprintf(outfile, "\tConvergence = %3.1e\n", params.convergence);
95
 
  fprintf(outfile, "\tRestart     =     %s\n", params.restart ? "Yes" : "No");
96
 
  fprintf(outfile, "\tCache Level =     %1d\n", params.cachelev);
97
 
  fprintf(outfile, "\tAO Basis        =     %s\n", 
 
438
  fprintf(outfile, "\tMaxiter       =    %4d\n", params.maxiter);
 
439
  fprintf(outfile, "\tConvergence   = %3.1e\n", params.convergence);
 
440
  fprintf(outfile, "\tRestart       =     %s\n", params.restart ? "Yes" : "No");
 
441
  fprintf(outfile, "\tCache Level   =     %1d\n", params.cachelev);
 
442
  fprintf(outfile, "\tModel III     =     %s\n", params.sekino ? "Yes" : "No");
 
443
  fprintf(outfile, "\tDIIS          =     %s\n", params.diis ? "Yes" : "No");
 
444
  fprintf(outfile, "\tAO Basis      =     %s\n", 
98
445
          params.aobasis ? "Yes" : "No");
99
 
  fprintf(outfile, "\tExcited State Computation =     %s\n", 
100
 
          params.ground ? "No" : "Yes");
 
446
  fprintf(outfile, "\tABCD            =     %s\n", params.abcd);
 
447
  fprintf(outfile, "\tLocal CC        =     %s\n", params.local ? "Yes" : "No");
101
448
  if(params.local) {
102
449
    fprintf(outfile, "\tLocal Cutoff    = %3.1e\n", local.cutoff);
103
450
    fprintf(outfile, "\tLocal Method    =    %s\n", local.method);
104
451
    fprintf(outfile, "\tWeak pairs      =    %s\n", local.weakp);
105
452
    fprintf(outfile, "\tFilter singles  =    %s\n", local.filter_singles ? "Yes" : "No");
106
 
  }
107
 
  fprintf(outfile, "\tStates sought per irrep     =");
108
 
  for (i=0;i<moinfo.nirreps;++i)
109
 
    fprintf(outfile, " %s %d,", moinfo.labels[i],
110
 
        params.states_per_irrep[i]);
111
 
  fprintf(outfile,"\n");
112
 
  fprintf(outfile, "\tLs sought per irrep         =");
113
 
  for (i=0;i<moinfo.nirreps;++i)
114
 
    fprintf(outfile, " %s %d,", moinfo.labels[i],
115
 
        params.Ls_per_irrep[i]);
116
 
  fprintf(outfile,"\n");
117
 
  fprintf(outfile, "\tLocal CC        =     %s\n", params.local ? "Yes" : "No");
118
 
 
119
 
  /* get EOM energies and R0s from CC_INFO */
120
 
  /* compute total number of states */
121
 
  params.cceom_energy = (double **) malloc(moinfo.nirreps*sizeof(double *));
122
 
  params.R0 = (double **) malloc(moinfo.nirreps*sizeof(double *));
123
 
  for (i=0;i<moinfo.nirreps;++i) {
124
 
    params.cceom_energy[i] = malloc(params.Ls_per_irrep[i]*sizeof(double)); 
125
 
    params.R0[i] = malloc(params.Ls_per_irrep[i]*sizeof(double)); 
126
 
  }
127
 
  if (params.ground) { /* just a ground state calculation */
128
 
    for (i=0;i<moinfo.nirreps;++i) {
129
 
      if (params.Ls_per_irrep[i]) {
130
 
        params.cceom_energy[i][0] = 0.0;
131
 
        params.R0[i][0] = 1.0;
132
 
      }
133
 
    }
134
 
  }
135
 
  else { /* excited state */
136
 
    for (i=0;i<moinfo.nirreps;++i) {
137
 
      for (j=0; j<params.Ls_per_irrep[i]; ++j) {
138
 
        sprintf(lbl,"EOM CCSD Energy for root %d %d", i, j);
139
 
        psio_read_entry(CC_INFO, lbl, (char *) &(params.cceom_energy[i][j]),sizeof(double));
140
 
        sprintf(lbl,"EOM CCSD R0 for root %d %d", i, j);
141
 
        psio_read_entry(CC_INFO, lbl, (char *) &(params.R0[i][j]),sizeof(double));
142
 
      }
143
 
    }
144
 
  }
145
 
 
146
 
  for (i=0;i<moinfo.nirreps;++i)
147
 
    for (j=0; j<params.Ls_per_irrep[i]; ++j) {
148
 
      fprintf(outfile,"\tparams.cceom_energy[%d][%d] = %15.10lf\n",i,j,params.cceom_energy[i][j]);
149
 
      fprintf(outfile,"\tparams.R0[%d][%d] = %15.10lf\n",i,j,params.R0[i][j]);
150
 
    }
151
 
 
152
 
  /* determine L0 value */
153
 
  params.L0 = (params.ground ? 1.0 : 0.0);
154
 
 
 
453
    fprintf(outfile, "\tLocal pairs       =    %s\n", local.pairdef);
 
454
    fprintf(outfile, "\tLocal CPHF cutoff =  %3.1e\n", local.cphf_cutoff);
 
455
  }
 
456
 
 
457
  fprintf(outfile,"\tParamaters for left-handed eigenvectors:\n");
 
458
  fprintf(outfile,"\t    Irr   Root  Ground-State?    EOM energy        R0\n");
 
459
  for (i=0; i<params.nstates; ++i) {
 
460
    fprintf(outfile,"\t%3d %3d %5d %10s %18.10lf %14.10lf\n", i+1, pL_params[i].irrep, pL_params[i].root+1,
 
461
            (pL_params[i].ground ? "Yes":"No"), pL_params[i].cceom_energy, pL_params[i].R0);
 
462
  }
 
463
 
 
464
  for (i=0; i<params.nstates; ++i) {
 
465
    fprintf(outfile,"\tLabels for eigenvector %d:\n\t%s, %s, %s, %s, %s, %s\n",
 
466
            i+1,pL_params[i].L1A_lbl,pL_params[i].L1B_lbl,pL_params[i].L2AA_lbl,pL_params[i].L2BB_lbl,
 
467
            pL_params[i].L2AB_lbl, pL_params[i].L2RHF_lbl);
 
468
  }
 
469
 
 
470
  fflush(outfile);
155
471
  return;
156
472
}
157