1!
2!  Dalton, a molecular electronic structure program
3!  Copyright (C) by the authors of Dalton.
4!
5!  This program is free software; you can redistribute it and/or
6!  modify it under the terms of the GNU Lesser General Public
7!  License version 2.1 as published by the Free Software Foundation.
8!
9!  This program is distributed in the hope that it will be useful,
10!  but WITHOUT ANY WARRANTY; without even the implied warranty of
11!  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
12!  Lesser General Public License for more details.
13!
14!  If a copy of the GNU LGPL v2.1 was not distributed with this
15!  code, you can obtain one at https://www.gnu.org/licenses/old-licenses/lgpl-2.1.en.html.
16!
17!
18C
19c /* deck CC_ETADRV1 */
20*=====================================================================*
21      SUBROUTINE CC_ETADRV1(TYPE,LABEL,ISYMS,ISTAT,EIGV,
22     &                      ISYMO,FREQS,LORX,ICAU,NVEC,MAXVEC,
23     &                      WORK,LWORK)
24*---------------------------------------------------------------------*
25*
26*    Purpose: calculate response eta vectors, used to build the
27*             right-hand-side vectors for the lagrangian multipliers
28*             and as intermediates in the hyperpolarizability
29*             and n-photon-transition matrix calculations
30*
31*             for excited states the X vectors are identical to the
32*             rhs vectors for the left eigenvector response equations
33*
34*     implemented:  L:   ORDER = 2, 3
35*                   LE:  ORDER = 1, 2
36*                   CL:  ORDER = 2
37*
38*     Written by Christof Haettig april/june/july 1997.
39*     extensions for Cauchy eta vectors in March 1998.
40*
41*     TENTATIVE MODIFIED VERSION FOR RELAXED RHS for EL1   Sonia
42*     one perturbation with orb-rel
43*     Orbital relaxation for EX1 introduce. Sonia Coriani, April 2000
44*=====================================================================*
45#if defined (IMPLICIT_NONE)
46      IMPLICIT NONE
47#else
48#include "implicit.h"
49#endif
50#include "priunit.h"
51#include "ccsdinp.h"
52#include "ccsdsym.h"
53#include "ccorb.h"
54#include "dummy.h"
55Csonia
56#include "cclists.h"
57
58* local parameters:
59      CHARACTER*(20) MSGDBG
60      PARAMETER (MSGDBG = '[debug] CC_ETADRV1> ')
61      LOGICAL LOCDBG
62      PARAMETER (LOCDBG = .false.)
63
64      CHARACTER TYPE*(*)
65
66      INTEGER NVEC, MAXVEC, LWORK
67      INTEGER ISYMO(MAXVEC,*), ICAU(MAXVEC,*)
68      INTEGER ISYMS(MAXVEC,*), ISTAT(MAXVEC,*)
69      LOGICAL LORX(MAXVEC,*)
70
71      CHARACTER*8 LABEL(MAXVEC,*)
72
73#if defined (SYS_CRAY)
74      REAL FREQS(MAXVEC,*), EIGV(MAXVEC,*)
75      REAL WORK(LWORK)
76      REAL ZERO
77      REAL DDOT, XNORM, RDUM
78#else
79      DOUBLE PRECISION FREQS(MAXVEC,*), EIGV(MAXVEC,*)
80      DOUBLE PRECISION WORK(LWORK)
81      DOUBLE PRECISION ZERO
82      DOUBLE PRECISION DDOT, XNORM, RDUM
83#endif
84      PARAMETER (ZERO = 0.0d0)
85
86      CHARACTER MODEL*(10), CDUM*3
87      INTEGER MX0KTRAN, MX1GTRAN, MX2FTRAN, MX1FATRAN
88      INTEGER MXTRAN, MX0GTRAN, MX1FTRAN, MX0FATRAN, MXEATRAN
89      INTEGER K0KTRAN, K1GTRAN, K2FTRAN, K1FATRAN
90      INTEGER K0GTRAN, K1FTRAN, K0FATRAN, KEATRAN
91      INTEGER N0KTRAN, N1GTRAN, N2FTRAN, N1FATRAN
92      INTEGER N0GTRAN, N1FTRAN, N0FATRAN, NEATRAN
93      INTEGER IOPT, ISYM, IVEC, ORDER, MPERM, NSTAT, IOPT1
94      INTEGER KEND0, LEND0, KEND1, LEND1, LMAX1, LMAX2, KCHI1, KCHI2
95      INTEGER KEND2, LEND2, KRHS1, KRHS2, IDUM
96
97* external functions
98      INTEGER ILSTSYM
99
100*---------------------------------------------------------------------*
101* check number of required eta/rhs vectors, if zero return immediatly:
102*---------------------------------------------------------------------*
103      IF (NVEC.EQ.0) RETURN
104
105*---------------------------------------------------------------------*
106* print header for eta/rhs vector section
107*---------------------------------------------------------------------*
108      WRITE (LUPRI,'(7(/1X,2A),/)')
109     & '------------------------------------',
110     &                               '-------------------------------',
111     & '|                 OUTPUT FROM ETA/RH',
112     &                               'S VECTOR SECTION  (ETADRV1)   |',
113     & '------------------------------------',
114     &                               '-------------------------------'
115      CALL FLSHFO(LUPRI)
116
117*---------------------------------------------------------------------*
118      IF (.NOT. (CCS .OR. CC2 .OR. CCSD) ) THEN
119         CALL QUIT(
120     *      'CC_ETADRV called for unknown Coupled Cluster.')
121      END IF
122
123      IF (TYPE(1:2).EQ.'X1') THEN
124        WRITE(LUPRI,*) 'X1 vectors not implemented in CC_ETADRV,'
125        WRITE(LUPRI,*) 'routine CCRHSVEC should be used instead.'
126        CALL QUIT('X1 vectors not implemented in CC_ETADRV.')
127      ELSE IF (TYPE(1:2).EQ.'X2') THEN
128        ORDER = 2
129        NSTAT = 0
130        MPERM = 2
131      ELSE IF (TYPE(1:2).EQ.'X3') THEN
132        ORDER = 3
133        NSTAT = 0
134        MPERM = 6
135C     ELSE IF (TYPE(1:2).EQ.'X4') THEN
136C       ORDER = 4
137C       NSTAT = 0
138C       MPERM = ??
139      ELSE IF (TYPE(1:3).EQ.'EX1') THEN
140        ORDER = 1
141        NSTAT = 1
142        MPERM = 1
143      ELSE IF (TYPE(1:3).EQ.'EX2') THEN
144        ORDER = 2
145        NSTAT = 1
146        MPERM = 2
147        WRITE(LUPRI,*) 'warning: X vectors ',TYPE(1:3),' not tested!!!.'
148      ELSE IF (TYPE(1:3).EQ.'CX2') THEN
149        ORDER = 2
150        NSTAT = 0
151        MPERM = 2
152      ELSE
153        WRITE(LUPRI,*) 'rhs vectors ',TYPE(1:2),' not implemented.'
154        CALL QUIT('required rhs vectors not implemented.')
155      END IF
156
157
158* print some debug/info output
159      IF (IPRINT .GT. 10)
160     *       WRITE(LUPRI,*) 'CC_ETADRV Workspace:',LWORK
161
162*---------------------------------------------------------------------*
163* allocate & initialize work space for lists
164*---------------------------------------------------------------------*
165
166      MXTRAN  = MPERM * NVEC
167
168      MX0KTRAN  = 5 * MXTRAN
169      MX0GTRAN  = MXDIM_GTRAN  * MXTRAN
170      MX1GTRAN  = MXDIM_GTRAN  * MXTRAN
171      MX1FTRAN  = MXDIM_FTRAN  * MXTRAN
172      MX2FTRAN  = MXDIM_FTRAN  * MXTRAN
173      MX0FATRAN = MXDIM_FATRAN * MXTRAN
174      MX1FATRAN = MXDIM_FATRAN * MXTRAN
175      MXEATRAN  = MXDIM_XEVEC  * MXTRAN
176
177      K0KTRAN  = 1
178      K0GTRAN  = K0KTRAN  + MX0KTRAN
179      K1GTRAN  = K0GTRAN  + MX0GTRAN
180      K1FTRAN  = K1GTRAN  + MX1GTRAN
181      K2FTRAN  = K1FTRAN  + MX1FTRAN
182      K0FATRAN = K2FTRAN  + MX2FTRAN
183      K1FATRAN = K0FATRAN + MX0FATRAN
184      KEATRAN  = K1FATRAN + MX1FATRAN
185      KEND0    = KEATRAN  + MXEATRAN
186      LEND0    = LWORK - KEND0
187
188      IF (LEND0 .LT. 0 ) THEN
189        WRITE(LUPRI,*) 'Insufficient work space in CC_ETADRV.'
190      END IF
191
192*---------------------------------------------------------------------*
193* set up lists for G, F and F{A} transformations and ETA{O} vectors:
194*---------------------------------------------------------------------*
195      CALL CC_ETA_SETUP1(TYPE,NSTAT,ORDER,LABEL,ISTAT,EIGV,
196     &                   ISYMO,FREQS,LORX,ICAU,
197     &                   NVEC, MAXVEC,  MXTRAN,
198     &                   WORK(K0KTRAN), N0KTRAN,
199     &                   WORK(K0GTRAN), N0GTRAN,
200     &                   WORK(K1GTRAN), N1GTRAN,
201     &                   WORK(K1FTRAN), N1FTRAN,
202     &                   WORK(K2FTRAN), N2FTRAN,
203     &                   WORK(K0FATRAN),N0FATRAN,
204     &                   WORK(K1FATRAN),N1FATRAN,
205     &                   WORK(KEATRAN), NEATRAN  )
206
207*---------------------------------------------------------------------*
208* initialize ETA vector files:
209* open files and fill with zeros. Later we always add to vecs on file.
210* sonia
211*---------------------------------------------------------------------*
212      LMAX1 = 0
213      LMAX2 = 0
214      DO ISYM = 1, NSYM
215        LMAX1 = MAX(LMAX1,NT1AM(ISYM))
216        LMAX2 = MAX(LMAX2,NT2AM(ISYM))
217      END DO
218
219      KCHI1 = KEND0
220      KCHI2 = KCHI1 + LMAX1
221      KEND1 = KCHI2 + LMAX2
222      LEND1 = LWORK - KEND1
223
224      IF (LEND1 .LT. 0 ) THEN
225        CALL QUIT('Insufficient work space in CC_ETADRV.')
226      END IF
227
228      CALL DZERO(WORK(KCHI1),LMAX1)
229      IF (.NOT.CCS) CALL DZERO(WORK(KCHI2),LMAX2)
230
231      IF (CCS) THEN
232         MODEL = 'CCS       '
233         IOPT  = 1
234      ELSE IF (CC2) THEN
235         MODEL = 'CC2       '
236         IOPT  = 3
237      ELSE IF (CCSD) THEN
238         MODEL = 'CCSD      '
239         IOPT  = 3
240      ELSE
241         CALL QUIT('Unknown coupled cluster model in CC_ETADRV.')
242      END IF
243
244      DO IVEC = 1, NVEC
245        ISYM = ILSTSYM(TYPE,IVEC)
246        CALL CC_WRRSP(TYPE,IVEC,ISYM,IOPT,MODEL,IDUMMY,
247     &                WORK(KCHI1),WORK(KCHI2),WORK(KEND1),LEND1)
248      END DO
249
250*---------------------------------------------------------------------*
251* calculate H matrix contributions:
252*---------------------------------------------------------------------*
253      IF (TYPE(1:2).EQ.'X3') THEN
254        IOPT1 = 4
255        CALL CC_HMAT('L0','R1','R1','R1',TYPE,N0KTRAN, 0,
256     &               WORK(K0KTRAN),IDUMMY,IDUMMY,
257     &               WORK(KEND0), LEND0, IOPT1 )
258      END IF
259
260      IF (LOCDBG) THEN
261        WRITE(LUPRI,*)
262     &      MSGDBG, 'NORM^2 of ETA vectors after H matrix terms:'
263        DO IVEC = 1, NVEC
264          ISYM = ILSTSYM(TYPE,IVEC)
265          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
266     &                  WORK(KCHI1),WORK(KCHI2))
267          XNORM = DDOT(NT1AM(ISYM),WORK(KCHI1),1,WORK(KCHI1),1)
268          IF (.NOT.CCS)
269     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KCHI2),1,WORK(KCHI2),1)
270          WRITE(LUPRI,*) MSGDBG, IVEC,XNORM
271        END DO
272      END IF
273
274*---------------------------------------------------------------------*
275* calculate G matrix contributions:
276*---------------------------------------------------------------------*
277      IF (TYPE(1:2).EQ.'X2') THEN
278        IOPT1 = 4
279        CALL CC_GMATRIX('L0 ','R1 ','R1 ',TYPE,N0GTRAN, 0,
280     &                 WORK(K0GTRAN),IDUM,RDUM,WORK(KEND0),LEND0,IOPT1)
281      ELSE IF (TYPE(1:2).EQ.'X3') THEN
282        IOPT1 = 4
283        CALL CC_GMATRIX('L0 ','R2 ','R1 ',TYPE,N0GTRAN, 0,
284     &                 WORK(K0GTRAN),IDUM,RDUM,WORK(KEND0),LEND0,IOPT1)
285        IOPT1 = 4
286        CALL CC_GMATRIX('L1 ','R1 ','R1 ',TYPE,N1GTRAN, 0,
287     &                 WORK(K1GTRAN),IDUM,RDUM,WORK(KEND0),LEND0,IOPT1)
288      ELSE IF (TYPE(1:3).EQ.'EX2') THEN
289        IOPT1 = 4
290        CALL CC_GMATRIX('LE ','R1 ','R1 ',TYPE,N0GTRAN, 0,
291     &                 WORK(K0GTRAN),IDUM,RDUM,WORK(KEND0),LEND0,IOPT1)
292      ELSE IF (TYPE(1:3).EQ.'CX2') THEN
293        IOPT1 = 4
294        CALL CC_GMATRIX('L0 ','RC ','RC ',TYPE,N0GTRAN, 0,
295     &                 WORK(K0GTRAN),IDUM,RDUM,WORK(KEND0),LEND0,IOPT1)
296      END IF
297
298      IF (LOCDBG) THEN
299        WRITE(LUPRI,*)
300     &       MSGDBG, 'NORM^2 of ETA vectors after G matrix terms:'
301        DO IVEC = 1, NVEC
302          ISYM = ILSTSYM(TYPE,IVEC)
303          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
304     &                  WORK(KCHI1),WORK(KCHI2))
305          XNORM = DDOT(NT1AM(ISYM),WORK(KCHI1),1,WORK(KCHI1),1)
306          IF (.NOT.CCS)
307     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KCHI2),1,WORK(KCHI2),1)
308          WRITE(LUPRI,*)
309     &       MSGDBG, IVEC,XNORM
310        END DO
311      END IF
312
313*---------------------------------------------------------------------*
314* calculate F matrix contributions:
315*---------------------------------------------------------------------*
316      IF (TYPE(1:2).EQ.'X2') THEN
317        IOPT1 = 4
318        CALL CC_FMATRIX(WORK(K1FTRAN),N1FTRAN,'L1 ','R1 ',IOPT1,TYPE,
319     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
320      ELSE IF (TYPE(1:3).EQ.'CX2') THEN
321        IOPT1 = 4
322        CALL CC_FMATRIX(WORK(K1FTRAN),N1FTRAN,'LC ','RC ',IOPT1,TYPE,
323     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
324      ELSE IF (TYPE(1:2).EQ.'X3') THEN
325        IOPT1 = 4
326        CALL CC_FMATRIX(WORK(K1FTRAN),N1FTRAN,'L1 ','R2 ',IOPT1,TYPE,
327     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
328        IOPT1 = 4
329        CALL CC_FMATRIX(WORK(K2FTRAN),N2FTRAN,'L2 ','R1 ',IOPT1,TYPE,
330     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
331      ELSE IF (TYPE(1:3).EQ.'EX2') THEN
332        IOPT1 = 4
333        CALL CC_FMATRIX(WORK(K1FTRAN),N1FTRAN,'LE ','R2 ',IOPT1,TYPE,
334     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
335        IOPT1 = 4
336        CALL CC_FMATRIX(WORK(K2FTRAN),N2FTRAN,'EL1','R1 ',IOPT1,TYPE,
337     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
338      ELSE IF (TYPE(1:3).EQ.'EX1') THEN
339        IOPT1 = 4    !note IOPT = 4  added to vector on file!!!!!!!!
340        CALL CC_FMATRIX(WORK(K1FTRAN),N1FTRAN,'LE ','R1 ',IOPT1,TYPE,
341     &                  IDUM, RDUM, 0, WORK(KEND0), LEND0)
342      END IF
343
344      IF (LOCDBG) THEN
345        WRITE(LUPRI,*)
346     &       MSGDBG, 'NORM^2 of ETA vectors after F matrix terms:'
347        DO IVEC = 1, NVEC
348          ISYM = ILSTSYM(TYPE,IVEC)
349          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
350     &                  WORK(KCHI1),WORK(KCHI2))
351          XNORM = DDOT(NT1AM(ISYM),WORK(KCHI1),1,WORK(KCHI1),1)
352          IF (.NOT.CCS)
353     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KCHI2),1,WORK(KCHI2),1)
354          WRITE(LUPRI,*)
355     &       MSGDBG, IVEC,XNORM
356        END DO
357      END IF
358
359*---------------------------------------------------------------------*
360* calculate F{O} matrix contributions:
361*---------------------------------------------------------------------*
362      IF (TYPE(1:2).EQ.'X2') THEN
363        CALL CCQR_FADRV('L0 ','o1 ','R1 ',TYPE,N0FATRAN, 0,
364     &                   WORK(K0FATRAN),IDUMMY,IDUMMY,
365     &                   WORK(KEND0), LEND0, 'FILE' )
366      ELSE IF (TYPE(1:3).EQ.'CX2') THEN
367        CALL CCQR_FADRV('L0 ','o1 ','RC ',TYPE,N0FATRAN, 0,
368     &                   WORK(K0FATRAN),IDUMMY,IDUMMY,
369     &                   WORK(KEND0), LEND0, 'FILE' )
370      ELSE IF (TYPE(1:2).EQ.'X3') THEN
371        CALL CCQR_FADRV('L0 ','o1 ','R2 ',TYPE,N0FATRAN, 0,
372     &                   WORK(K0FATRAN),IDUMMY,IDUMMY,
373     &                   WORK(KEND0), LEND0, 'FILE' )
374        CALL CCQR_FADRV('L1 ','o1 ','R1 ',TYPE,N1FATRAN, 0,
375     &                   WORK(K1FATRAN),IDUMMY,IDUMMY,
376     &                   WORK(KEND0), LEND0, 'FILE' )
377      ELSE IF (TYPE(1:3).EQ.'EX2') THEN
378        CALL CCQR_FADRV('LE ','o1 ','R1 ',TYPE,N0FATRAN, 0,
379     &                   WORK(K0FATRAN),IDUMMY,IDUMMY,
380     &                   WORK(KEND0), LEND0, 'FILE' )
381      END IF
382
383      IF (LOCDBG) THEN
384        WRITE(LUPRI,*)
385     &       MSGDBG,'NORM^2 of ETA vectors after F{O} matrix terms:'
386        DO IVEC = 1, NVEC
387          ISYM = ILSTSYM(TYPE,IVEC)
388          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
389     &                  WORK(KCHI1),WORK(KCHI2))
390          XNORM = DDOT(NT1AM(ISYM),WORK(KCHI1),1,WORK(KCHI1),1)
391          IF (.NOT.CCS)
392     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KCHI2),1,WORK(KCHI2),1)
393          WRITE(LUPRI,*)
394     &       MSGDBG, IVEC,XNORM
395        END DO
396      END IF
397
398*---------------------------------------------------------------------*
399* calculate ETA{O} vector contributions:
400*---------------------------------------------------------------------*
401      IF (TYPE(1:2).EQ.'X2') THEN
402        IOPT1 = 4
403        CALL CC_XIETA(WORK(KEATRAN),NEATRAN,IOPT1, ORDER, 'L1 ',
404     &                CDUM, IDUM, RDUM, TYPE, IDUM, RDUM,
405     &                .FALSE.,0, WORK(KEND0),LEND0)
406      ELSE IF (TYPE(1:3).EQ.'CX2') THEN
407        IOPT1 = 4
408        CALL CC_XIETA(WORK(KEATRAN),NEATRAN,IOPT1, ORDER, 'LC ',
409     &                CDUM, IDUM, RDUM, TYPE, IDUM, RDUM,
410     &                .FALSE.,0, WORK(KEND0),LEND0)
411      ELSE IF (TYPE(1:2).EQ.'X3') THEN
412        IOPT1 = 4
413        CALL CC_XIETA(WORK(KEATRAN),NEATRAN,IOPT1, ORDER, 'L2 ',
414     &                CDUM, IDUM, RDUM, TYPE, IDUM, RDUM,
415     &                .FALSE.,0, WORK(KEND0),LEND0)
416      ELSE IF (TYPE(1:3).EQ.'EX2') THEN
417        IOPT1 = 4
418        CALL CC_XIETA(WORK(KEATRAN),NEATRAN,IOPT1, ORDER, 'EL1',
419     &                CDUM, IDUM, RDUM, TYPE, IDUM, RDUM,
420     &                .FALSE.,0, WORK(KEND0),LEND0)
421      ELSE IF (TYPE(1:3).EQ.'EX1') THEN
422        !sonia
423        WRITE(LUPRI,*) 'CC_ETADRV1: I am entering CC_XIETA'
424        IOPT1 = 4
425        CALL CC_XIETA(WORK(KEATRAN),NEATRAN,IOPT1, ORDER, 'LE ',
426     &                CDUM, IDUM, RDUM, TYPE, IDUM, RDUM,
427     &                .FALSE.,0, WORK(KEND0),LEND0)
428
429      END IF
430
431      IF (LOCDBG) THEN
432        WRITE(LUPRI,*)
433     &     MSGDBG,'NORM^2 of ETA vectors after ETA{O} vec. terms:'
434        DO IVEC = 1, NVEC
435          ISYM = ILSTSYM(TYPE,IVEC)
436          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
437     &                  WORK(KCHI1),WORK(KCHI2))
438          XNORM = DDOT(NT1AM(ISYM),WORK(KCHI1),1,WORK(KCHI1),1)
439          IF (.NOT.CCS)
440     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KCHI2),1,WORK(KCHI2),1)
441          WRITE(LUPRI,*) MSGDBG, IVEC,XNORM
442        END DO
443      END IF
444*---------------------------------------------------------------------*
445* test (static) EX1 vectors by calculating the excited state FOP's
446*---------------------------------------------------------------------*
447      IF (LOCDBG .AND. TYPE(1:3).EQ.'EX1') THEN
448        KRHS1 = KEND1
449        KRHS2 = KRHS1 + LMAX1
450        KEND2 = KRHS2 + LMAX2
451        LEND2 = LWORK - KEND2
452
453        IF (LEND2 .LT. 0 ) THEN
454          CALL QUIT('Insufficient work space in CC_ETADRV.')
455        END IF
456
457        WRITE(LUPRI,*) MSGDBG,'excited state first-order properties:'
458        DO IVEC = 1, NVEC
459        IF (ISYMO(IVEC,1).EQ.1 .AND. FREQS(IVEC,1).EQ.ZERO) THEN
460          ISYM = ILSTSYM(TYPE,IVEC)
461          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
462     &                  WORK(KCHI1),WORK(KCHI2))
463          CALL CC_RDRSP('RE',ISTAT(IVEC,1),ISYMS(IVEC,1),IOPT,MODEL,
464     &                  WORK(KRHS1),WORK(KRHS2))
465          XNORM = DDOT(NT1AM(ISYM),WORK(KCHI1),1,WORK(KRHS1),1)
466          IF (.NOT. CCS)
467     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KCHI2),1,WORK(KRHS2),1)
468          WRITE(LUPRI,'(A,I3,2X,F12.8,2X,A,2X,L2,2X,F12.8)') MSGDBG,
469     &              ISTAT(IVEC,1),EIGV(IVEC,1),
470     &              LABEL(IVEC,1),LORX(IVEC,1),XNORM
471        ELSE
472          WRITE(LUPRI,'(A,I3,2X,F12.8,2X,A,2X,L2,2X,F12.8)') MSGDBG,
473     &              ISTAT(IVEC,1),EIGV(IVEC,1),
474     &              LABEL(IVEC,1),LORX(IVEC,1),ZERO
475        END IF
476        END DO
477
478      END IF
479*---------------------------------------------------------------------*
480* that's it:
481*---------------------------------------------------------------------*
482
483      RETURN
484      END
485
486*=====================================================================*
487*              END OF SUBROUTINE CC_ETADRV                            *
488*=====================================================================*
489
490c /* deck CC_ETA_SETUP1 */
491*=====================================================================*
492      SUBROUTINE CC_ETA_SETUP1(TYPE,NSTAT,ORDER,LAB,ISTAT,
493     &                        EIGV,ISYMO,FREQ,LORX,ICAU,
494     &                        NVEC, MAXVEC, MXTRAN,
495     &                        I0KTRAN, N0KTRAN,
496     &                        I0GTRAN, N0GTRAN,
497     &                        I1GTRAN, N1GTRAN,
498     &                        I1FTRAN, N1FTRAN,
499     &                        I2FTRAN, N2FTRAN,
500     &                        I0FATRAN,N0FATRAN,
501     &                        I1FATRAN,N1FATRAN,
502     &                        IEATRAN, NEATRAN  )
503
504*---------------------------------------------------------------------*
505*
506*    Purpose: set up for CC_ETA section
507*                - list of G matrix transformations
508*                - list of F matrix transformations
509*                - list of F{O} matrix transformations
510*                - list of ETA{O} vector calculations
511*
512*     Written by Christof Haettig, april/june/july 1997.
513*     extensions for Cauchy eta vectors in march 1998.
514*     orb.-relax. or derivatives by Sonia Coriani, Apr. 2000
515*     TENTATIVE NEW VERSION TO USE XIETA...
516*     DO we need Second order operators??? SOP
517*=====================================================================*
518#if defined (IMPLICIT_NONE)
519      IMPLICIT NONE
520#else
521#  include "implicit.h"
522#endif
523#include "priunit.h"
524#include "cclists.h"
525
526* local parameters:
527      CHARACTER*(23) MSGDBG
528      PARAMETER (MSGDBG = '[debug] CC_ETA_SETUP1> ')
529      LOGICAL LOCDBG
530      PARAMETER (LOCDBG = .false.)
531
532
533      INTEGER MXORD, MXORD2, MXORD3, MXSTAT
534      PARAMETER (MXORD  = 4, MXSTAT = 2)
535      PARAMETER (MXORD2 = MXORD *(MXORD-1)/2 )
536      PARAMETER (MXORD3 = MXORD2*(MXORD-2)/3 )
537
538
539      INTEGER MXTRAN, NSTAT, ORDER, MAXVEC, NVEC
540
541      CHARACTER*(*) TYPE
542
543      CHARACTER*(8) LAB(MAXVEC,*)
544      INTEGER ISTAT(MAXVEC,*), ICAU(MAXVEC,*), ISYMO(MAXVEC,*)
545      LOGICAL LORX(MAXVEC,*)
546
547#if defined (SYS_CRAY)
548      REAL FREQ(MAXVEC,*), EIGV(MAXVEC,*)
549#else
550      DOUBLE PRECISION FREQ(MAXVEC,*), EIGV(MAXVEC,*)
551#endif
552
553      INTEGER I0KTRAN(5,MXTRAN)
554      INTEGER I0GTRAN(MXDIM_GTRAN,MXTRAN)
555      INTEGER I1GTRAN(MXDIM_GTRAN,MXTRAN)
556      INTEGER I1FTRAN(MXDIM_FTRAN,MXTRAN)
557      INTEGER I2FTRAN(MXDIM_FTRAN,MXTRAN)
558      INTEGER I0FATRAN(MXDIM_FATRAN,MXTRAN)
559      INTEGER I1FATRAN(MXDIM_FATRAN,MXTRAN)
560      INTEGER IEATRAN(MXDIM_XEVEC,MXTRAN)
561
562      INTEGER N0KTRAN, N1GTRAN, N2FTRAN, N1FATRAN
563      INTEGER          N0GTRAN, N1FTRAN, N0FATRAN, NEATRAN
564
565      INTEGER IVEC, ISYML, ITRAN, I, IDX, IDXA, IDXB, IDXAB, IDXS
566
567      INTEGER A, B, C, D
568      PARAMETER (A = 1, B = 2, C = 3, D = 4)
569      INTEGER AB, AC, AD, BC, BD, CD
570      PARAMETER (AB = 1, AC = 2, BC = 3, AD = 4, BD = 5, CD = 6)
571      INTEGER ABC, ABD, ACD, BCD
572      PARAMETER (ABC = 1, ABD = 2, ACD = 3, BCD = 4)
573
574      INTEGER NS2A, NS3A, NP3AB, NP4AB, NT4ABC
575      PARAMETER (NS2A = 2, NS3A = 3, NP3AB = 3, NP4AB = 6, NT4ABC = 4)
576
577      INTEGER ISA(NS3A), ISB(NS3A), ISC(NS3A)
578      INTEGER IPAB(NP4AB), IPC(NP4AB), IPD(NP4AB), IPCD(NP4AB)
579      INTEGER ITABC(NT4ABC), ITD(NT4ABC)
580
581      DATA ISA  / A, B, C/
582      DATA ISB  / B, A, A/
583      DATA ISC  / C, C, B/
584
585      DATA IPAB / AB, AC, BC, AD, BD, CD /
586      DATA IPC  / C,  B,  A,  B,  A,  A  /
587      DATA IPD  / D,  D,  D,  C,  C,  B  /
588      DATA IPCD / CD, BD, AD, BC, AC, AB /
589
590      DATA ITABC / ABC, ABD, ACD, BCD /
591      DATA ITD   / D,   C,   B,   A   /
592
593
594      INTEGER IL0
595      PARAMETER (IL0 = 0)  ! index for zeroth-order zeta vector
596      INTEGER IL1(MXORD),  IR1(MXORD), IOP(MXORD), ISYM(MXORD)
597      INTEGER IL2(MXORD2), IR2(MXORD2)
598      INTEGER IE0(MXSTAT), IE1(MXORD,MXSTAT), ISYMS(MXSTAT)
599      INTEGER IRELAX(MXORD)
600      INTEGER LEN
601      INTEGER NRELAX
602
603      CHARACTER CLASS*(5)
604
605
606
607* external functions:
608      INTEGER IROPER
609      INTEGER IR1TAMP
610      INTEGER IR1KAPPA
611      INTEGER IR2TAMP
612      INTEGER IL1ZETA
613      INTEGER IL2ZETA
614      INTEGER IEL1AMP
615      INTEGER IEL2AMP
616      INTEGER ILRCAMP
617      INTEGER ILC1AMP
618
619      CALL QENTER('CC_ETA_SETUP1')
620
621*---------------------------------------------------------------------*
622* initializations:
623*---------------------------------------------------------------------*
624      N0KTRAN  = 0
625      N0GTRAN  = 0
626      N1GTRAN  = 0
627      N1FTRAN  = 0
628      N2FTRAN  = 0
629      N0FATRAN = 0
630      N1FATRAN = 0
631      NEATRAN  = 0
632
633*---------------------------------------------------------------------*
634* start loop over all requested ETA-vectors:
635*---------------------------------------------------------------------*
636
637      DO IVEC = 1, NVEC
638
639* eigenvectors that contribute:
640        IF (NSTAT.EQ.1) THEN
641          DO IDXS = 1, NSTAT
642            IE0(IDXS) = ISTAT(IVEC,IDXS)
643          END DO
644        END IF
645
646* operators:
647        IF (ORDER.GE.1) THEN
648          DO IDXA = 1, ORDER
649            IOP(IDXA) = IROPER(LAB(IVEC,IDXA),ISYM(IDXA))
650          END DO
651        END IF
652
653* relaxation flags:
654      IF (TYPE(1:1).EQ.'X' .OR. TYPE(1:2).EQ.'EX') THEN
655        NRELAX = 0
656        DO IDXA = 1, ORDER
657         IF ( LORX(IVEC,IDXA) ) THEN
658           IRELAX(IDXA) = IR1KAPPA(LAB(IVEC,IDXA),
659     &                            FREQ(IVEC,IDXA),ISYM(IDXA))
660           NRELAX = NRELAX + 1
661         ELSE
662           IRELAX(IDXA) = 0
663         END IF
664        END DO
665      ELSE
666        NRELAX = 0
667        DO IDXA = 1, ORDER
668         IRELAX(IDXA) = 0
669        END DO
670      END IF
671
672      IF (NRELAX.GT.1) THEN
673         CALL QUIT('NRELAX TOO LARGE IN CC_ETA_SETUP1.')
674      END IF
675
676
677* operators and first-order vectors that contribute:
678
679        IF (TYPE(1:1).EQ.'X' .AND. ORDER.GT.1) THEN
680          DO IDXA = 1, ORDER
681            IL1(IDXA) = IL1ZETA(LAB(IVEC,IDXA),LORX(IVEC,IDXA),
682     &                          FREQ(IVEC,IDXA),ISYM(IDXA))
683            IR1(IDXA) = IR1TAMP(LAB(IVEC,IDXA),LORX(IVEC,IDXA),
684     &                          FREQ(IVEC,IDXA),ISYM(IDXA))
685          END DO
686        END IF
687        IF (TYPE(1:2).EQ.'CX' .AND. ORDER.GT.1) THEN
688          DO IDXA = 1, ORDER
689            IR1(IDXA) = ILRCAMP(LAB(IVEC,IDXA),ICAU(IVEC,IDXA),ISYML)
690            IL1(IDXA) = ILC1AMP(LAB(IVEC,IDXA),ICAU(IVEC,IDXA),ISYML)
691          END DO
692        END IF
693Csonia
694        IF (TYPE(1:2).EQ.'EX' .AND. ORDER.GE.1) THEN
695          !I would need relaxation here for case EX1
696          DO IDXA = 1, ORDER
697            IR1(IDXA) = IR1TAMP(LAB(IVEC,IDXA),LORX(IVEC,IDXA),
698     &                          FREQ(IVEC,IDXA),ISYM(IDXA))
699          END DO
700          IF (ORDER.GT.1) THEN
701            call quit(' Sonia: please insert LORXA and LPROJ in call')
702            IE1(IDXA,1) =
703     &           IEL1AMP(ISTAT(IVEC,1), EIGV(IVEC,1),ISYMS(1),
704     &                   LAB(IVEC,IDXA),FREQ(IVEC,IDXA),ISYM(IDXA))
705                         !what about lprj and lorxa here?
706Cdef  INTEGER FUNCTION IEL1AMP(IEXCI,  EIGVNEW,ISYMS,
707Cdef *                         NEWLBLA,FRQANEW,ISYMA,LORXA,LPROJ )
708          END IF
709        END IF
710
711* second-order vectors that contribute:
712      IF (ORDER.GT.2 .OR. (ORDER.GE.2 .AND. NSTAT.GE.1)) THEN
713
714        IDXAB  = 0
715        DO IDXB = 2, ORDER
716        DO IDXA = 1, IDXB-1
717         IDXAB = IDXAB + 1
718         IR2(IDXAB) =
719     &       IR2TAMP(LAB(IVEC,IDXA),.FALSE.,FREQ(IVEC,IDXA),ISYM(IDXA),
720     &               LAB(IVEC,IDXB),.FALSE.,FREQ(IVEC,IDXB),ISYM(IDXB))
721        END DO
722        END DO
723
724       IF (TYPE(1:2).NE.'EX') THEN
725        IDXAB  = 0
726        DO IDXB = 2, ORDER
727        DO IDXA = 1, IDXB-1
728         IDXAB = IDXAB + 1
729         IL2(IDXAB) =
730     &           IL2ZETA(LAB(IVEC,IDXA),FREQ(IVEC,IDXA),ISYM(IDXA),
731     &                   LAB(IVEC,IDXB),FREQ(IVEC,IDXB),ISYM(IDXB))
732        END DO
733        END DO
734       END IF
735
736      END IF
737
738
739*---------------------------------------------------------------------*
740* set up list of H matrix transformations
741*---------------------------------------------------------------------*
742        IF (TYPE(1:2).EQ.'X3') THEN
743          N0KTRAN = N0KTRAN + 1
744          I0KTRAN(1,N0KTRAN) = IL0
745          I0KTRAN(2,N0KTRAN) = IR1(A)
746          I0KTRAN(3,N0KTRAN) = IR1(B)
747          I0KTRAN(4,N0KTRAN) = IR1(C)
748          I0KTRAN(5,N0KTRAN) = IVEC
749        END IF
750*---------------------------------------------------------------------*
751* set up list of G matrix transformations
752*---------------------------------------------------------------------*
753        IF (TYPE(1:2).EQ.'X2') THEN
754          N0GTRAN = N0GTRAN + 1
755          I0GTRAN(1,N0GTRAN) = IL0
756          I0GTRAN(2,N0GTRAN) = IR1(A)
757          I0GTRAN(3,N0GTRAN) = IR1(B)
758          I0GTRAN(4,N0GTRAN) = IVEC
759        ELSE IF (TYPE(1:3).EQ.'CX2') THEN
760          N0GTRAN = N0GTRAN + 1
761          I0GTRAN(1,N0GTRAN) = IL0
762          I0GTRAN(2,N0GTRAN) = IR1(A)
763          I0GTRAN(3,N0GTRAN) = IR1(B)
764          I0GTRAN(4,N0GTRAN) = IVEC
765        ELSE IF (TYPE(1:2).EQ.'X3') THEN
766          DO IDX = 1, NP3AB
767            N0GTRAN = N0GTRAN + 1
768            I0GTRAN(1,N0GTRAN) = IL0
769            I0GTRAN(2,N0GTRAN) = IR2(IPAB(IDX))
770            I0GTRAN(3,N0GTRAN) = IR1(IPC(IDX))
771            I0GTRAN(4,N0GTRAN) = IVEC
772          END DO
773
774          DO IDX = 1, NS3A
775            N1GTRAN = N1GTRAN + 1
776            I1GTRAN(1,N1GTRAN) = IL1(ISA(IDX))
777            I1GTRAN(2,N1GTRAN) = IR1(ISB(IDX))
778            I1GTRAN(3,N1GTRAN) = IR1(ISC(IDX))
779            I1GTRAN(4,N1GTRAN) = IVEC
780          END DO
781        ELSE IF (TYPE(1:3).EQ.'EX2') THEN
782          N0GTRAN = N0GTRAN + 1
783          I0GTRAN(1,N0GTRAN) = IE0(1)
784          I0GTRAN(2,N0GTRAN) = IR1(A)
785          I0GTRAN(3,N0GTRAN) = IR1(B)
786          I0GTRAN(4,N0GTRAN) = IVEC
787        END IF
788
789*---------------------------------------------------------------------*
790* set up list of F matrix transformations
791*---------------------------------------------------------------------*
792        IF (TYPE(1:2).EQ.'X2') THEN
793          N1FTRAN = N1FTRAN + 1
794          I1FTRAN(1,N1FTRAN) = IL1(A)
795          I1FTRAN(2,N1FTRAN) = IR1(B)
796          I1FTRAN(3,N1FTRAN) = IVEC
797
798          N1FTRAN = N1FTRAN + 1
799          I1FTRAN(1,N1FTRAN) = IL1(B)
800          I1FTRAN(2,N1FTRAN) = IR1(A)
801          I1FTRAN(3,N1FTRAN) = IVEC
802        ELSE IF (TYPE(1:3).EQ.'CX2') THEN
803          N1FTRAN = N1FTRAN + 1
804          I1FTRAN(1,N1FTRAN) = IL1(A)
805          I1FTRAN(2,N1FTRAN) = IR1(B)
806          I1FTRAN(3,N1FTRAN) = IVEC
807
808          N1FTRAN = N1FTRAN + 1
809          I1FTRAN(1,N1FTRAN) = IL1(B)
810          I1FTRAN(2,N1FTRAN) = IR1(A)
811          I1FTRAN(3,N1FTRAN) = IVEC
812        ELSE IF (TYPE(1:2).EQ.'X3') THEN
813          DO IDX = 1, NP3AB
814            N1FTRAN = N1FTRAN + 1
815            I1FTRAN(1,N1FTRAN) = IL1(IPC(IDX))
816            I1FTRAN(2,N1FTRAN) = IR2(IPAB(IDX))
817            I1FTRAN(3,N1FTRAN) = IVEC
818          END DO
819
820          DO IDX = 1, NP3AB
821            N2FTRAN = N2FTRAN + 1
822            I2FTRAN(1,N2FTRAN) = IL2(IPAB(IDX))
823            I2FTRAN(2,N2FTRAN) = IR1(IPC(IDX))
824            I2FTRAN(3,N2FTRAN) = IVEC
825          END DO
826        ELSE IF (TYPE(1:3).EQ.'EX2') THEN
827          N1FTRAN = N1FTRAN + 1
828          I1FTRAN(1,N1FTRAN) = IE0(1)
829          I1FTRAN(2,N1FTRAN) = IR2(AB)
830          I1FTRAN(3,N1FTRAN) = IVEC
831
832          N2FTRAN = N2FTRAN + 1
833          I2FTRAN(1,N2FTRAN) = IE1(A,1)
834          I2FTRAN(2,N2FTRAN) = IR1(B)
835          I2FTRAN(3,N2FTRAN) = IVEC
836
837          N2FTRAN = N2FTRAN + 1
838          I2FTRAN(1,N2FTRAN) = IE1(B,1)
839          I2FTRAN(2,N2FTRAN) = IR1(A)
840          I2FTRAN(3,N2FTRAN) = IVEC
841        ELSE IF (TYPE(1:3).EQ.'EX1') THEN
842          !sonia
843          N1FTRAN = N1FTRAN + 1
844          I1FTRAN(1,N1FTRAN) = IE0(1)
845          I1FTRAN(2,N1FTRAN) = IR1(A)
846          I1FTRAN(3,N1FTRAN) = IVEC
847        END IF
848
849*---------------------------------------------------------------------*
850* set up list of F{O} matrix transformations
851*---------------------------------------------------------------------*
852        IF (TYPE(1:2).EQ.'X2') THEN
853          N0FATRAN = N0FATRAN + 1
854          I0FATRAN(1,N0FATRAN) = IL0
855          I0FATRAN(2,N0FATRAN) = IOP(A)
856          I0FATRAN(3,N0FATRAN) = IR1(B)
857          I0FATRAN(4,N0FATRAN) = IVEC
858          I0FATRAN(5,N0FATRAN) = 0
859
860          N0FATRAN = N0FATRAN + 1
861          I0FATRAN(1,N0FATRAN) = IL0
862          I0FATRAN(2,N0FATRAN) = IOP(B)
863          I0FATRAN(3,N0FATRAN) = IR1(A)
864          I0FATRAN(4,N0FATRAN) = IVEC
865          I0FATRAN(5,N0FATRAN) = 0
866        ELSE IF (TYPE(1:3).EQ.'CX2') THEN
867          IF (ICAU(IVEC,A).EQ.0) THEN
868            N0FATRAN = N0FATRAN + 1
869            I0FATRAN(1,N0FATRAN) = IL0
870            I0FATRAN(2,N0FATRAN) = IOP(A)
871            I0FATRAN(3,N0FATRAN) = IR1(B)
872            I0FATRAN(4,N0FATRAN) = IVEC
873            I0FATRAN(5,N0FATRAN) = 0
874          END IF
875
876          IF (ICAU(IVEC,B).EQ.0) THEN
877            N0FATRAN = N0FATRAN + 1
878            I0FATRAN(1,N0FATRAN) = IL0
879            I0FATRAN(2,N0FATRAN) = IOP(B)
880            I0FATRAN(3,N0FATRAN) = IR1(A)
881            I0FATRAN(4,N0FATRAN) = IVEC
882            I0FATRAN(5,N0FATRAN) = 0
883          END IF
884        ELSE IF (TYPE(1:2).EQ.'X3') THEN
885          DO IDX = 1, NP3AB
886            N0FATRAN = N0FATRAN + 1
887            I0FATRAN(1,N0FATRAN) = IL0
888            I0FATRAN(2,N0FATRAN) = IOP(IPC(IDX))
889            I0FATRAN(3,N0FATRAN) = IR2(IPAB(IDX))
890            I0FATRAN(4,N0FATRAN) = IVEC
891            I0FATRAN(5,N0FATRAN) = 0
892          END DO
893
894          DO IDX = 1, NP3AB
895            N1FATRAN = N1FATRAN + 1
896            I1FATRAN(1,N1FATRAN) = IL1(ISA(IDX))
897            I1FATRAN(2,N1FATRAN) = IOP(ISB(IDX))
898            I1FATRAN(3,N1FATRAN) = IR1(ISC(IDX))
899            I1FATRAN(4,N1FATRAN) = IVEC
900            I1FATRAN(5,N1FATRAN) = 0
901            N1FATRAN = N1FATRAN + 1
902            I1FATRAN(1,N1FATRAN) = IL1(ISA(IDX))
903            I1FATRAN(2,N1FATRAN) = IOP(ISC(IDX))
904            I1FATRAN(3,N1FATRAN) = IR1(ISB(IDX))
905            I1FATRAN(4,N1FATRAN) = IVEC
906            I1FATRAN(5,N1FATRAN) = 0
907          END DO
908        ELSE IF (TYPE(1:3).EQ.'EX2') THEN
909          N0FATRAN = N0FATRAN + 1
910          I0FATRAN(1,N0FATRAN) = IE0(1)
911          I0FATRAN(2,N0FATRAN) = IOP(A)
912          I0FATRAN(3,N0FATRAN) = IR1(B)
913          I0FATRAN(4,N0FATRAN) = IVEC
914          I0FATRAN(5,N0FATRAN) = 0
915
916          N0FATRAN = N0FATRAN + 1
917          I0FATRAN(1,N0FATRAN) = IE0(1)
918          I0FATRAN(2,N0FATRAN) = IOP(B)
919          I0FATRAN(3,N0FATRAN) = IR1(A)
920          I0FATRAN(4,N0FATRAN) = IVEC
921          I0FATRAN(5,N0FATRAN) = 0
922        END IF
923
924*---------------------------------------------------------------------*
925* set up list of ETA{O} vector calculations:
926* Use the IXETRAN list
927* For the time being only relax flag in the EX1 case...
928*---------------------------------------------------------------------*
929        IF (TYPE(1:2).EQ.'X2') THEN
930          NEATRAN = NEATRAN + 1
931          IEATRAN(1,NEATRAN) = IOP(B)
932          IEATRAN(2,NEATRAN) = IL1(A)
933          IEATRAN(3,NEATRAN) = -1
934          IEATRAN(4,NEATRAN) = IVEC
935          IEATRAN(5,NEATRAN) = 0
936          IEATRAN(6,NEATRAN) = 0
937          IEATRAN(7,NEATRAN) = 0
938          IEATRAN(8,NEATRAN) = 0
939C
940          NEATRAN = NEATRAN + 1
941          IEATRAN(1,NEATRAN) = IOP(A)
942          IEATRAN(2,NEATRAN) = IL1(B)
943          IEATRAN(3,NEATRAN) = -1
944          IEATRAN(4,NEATRAN) = IVEC
945          IEATRAN(5,NEATRAN) = 0
946          IEATRAN(6,NEATRAN) = 0
947          IEATRAN(7,NEATRAN) = 0
948          IEATRAN(8,NEATRAN) = 0
949
950        ELSE IF (TYPE(1:3).EQ.'CX2') THEN
951          IF (ICAU(IVEC,B).EQ.0) THEN
952            NEATRAN = NEATRAN + 1
953            IEATRAN(1,NEATRAN) = IOP(B)
954            IEATRAN(2,NEATRAN) = IL1(A)
955            IEATRAN(3,NEATRAN) = -1
956            IEATRAN(4,NEATRAN) = IVEC
957            IEATRAN(5,NEATRAN) = 0
958            IEATRAN(6,NEATRAN) = 0
959            IEATRAN(7,NEATRAN) = 0
960            IEATRAN(8,NEATRAN) = 0
961          END IF
962
963          IF (ICAU(IVEC,A).EQ.0) THEN
964            NEATRAN = NEATRAN + 1
965            IEATRAN(1,NEATRAN) = IOP(A)
966            IEATRAN(2,NEATRAN) = IL1(B)
967            IEATRAN(3,NEATRAN) = -1
968            IEATRAN(4,NEATRAN) = IVEC
969            IEATRAN(5,NEATRAN) = 0
970            IEATRAN(6,NEATRAN) = 0
971            IEATRAN(7,NEATRAN) = 0
972            IEATRAN(8,NEATRAN) = 0
973          END IF
974        ELSE IF (TYPE(1:2).EQ.'X3') THEN
975          DO IDX = 1, NP3AB
976            NEATRAN = NEATRAN + 1
977            IEATRAN(1,NEATRAN) = IOP(IPC(IDX))
978            IEATRAN(2,NEATRAN) = IL2(IPAB(IDX))
979            IEATRAN(3,NEATRAN) = -1
980            IEATRAN(4,NEATRAN) = IVEC
981            IEATRAN(5,NEATRAN) = 0
982            IEATRAN(6,NEATRAN) = 0
983            IEATRAN(7,NEATRAN) = 0
984            IEATRAN(8,NEATRAN) = 0
985          END DO
986        ELSE IF (TYPE(1:3).EQ.'EX2') THEN
987          NEATRAN = NEATRAN + 1
988          IEATRAN(1,NEATRAN) = IOP(B)
989          IEATRAN(2,NEATRAN) = IE1(A,1)
990          IEATRAN(3,NEATRAN) = -1
991          IEATRAN(4,NEATRAN) = IVEC
992          IEATRAN(5,NEATRAN) = 0
993          IEATRAN(6,NEATRAN) = 0
994          IEATRAN(7,NEATRAN) = 0
995          IEATRAN(8,NEATRAN) = 0
996
997          NEATRAN = NEATRAN + 1
998          IEATRAN(1,NEATRAN) = IOP(A)
999          IEATRAN(2,NEATRAN) = IE1(B,1)
1000          IEATRAN(3,NEATRAN) = -1
1001          IEATRAN(4,NEATRAN) = IVEC
1002          IEATRAN(5,NEATRAN) = 0
1003          IEATRAN(6,NEATRAN) = 0
1004          IEATRAN(7,NEATRAN) = 0
1005          IEATRAN(8,NEATRAN) = 0
1006        ELSE IF (TYPE(1:3).EQ.'EX1') THEN
1007          !sonia
1008          NEATRAN = NEATRAN + 1
1009          IEATRAN(1,NEATRAN) = IOP(A)
1010          IEATRAN(2,NEATRAN) = IE0(1)
1011          IEATRAN(3,NEATRAN) = -1
1012          IEATRAN(4,NEATRAN) = IVEC
1013          IEATRAN(5,NEATRAN) = IRELAX(A)
1014          IEATRAN(6,NEATRAN) = 0
1015          IEATRAN(7,NEATRAN) = 0
1016          IEATRAN(8,NEATRAN) = 0
1017        END IF
1018
1019*---------------------------------------------------------------------*
1020* end loop over all requested ETA vectors
1021*---------------------------------------------------------------------*
1022      END DO
1023
1024*---------------------------------------------------------------------*
1025* print the lists:
1026*---------------------------------------------------------------------*
1027* general statistics:
1028      IF (TYPE(1:1).EQ.'X') THEN
1029          LEN = 2
1030          CLASS = ' eta '
1031      ELSE IF (TYPE(1:2).EQ.'CX') THEN
1032          LEN = 3
1033          CLASS = ' eta '
1034      ELSE IF (TYPE(1:2).EQ.'EX') THEN
1035          LEN = 3
1036          CLASS = ' rhs '
1037      ELSE
1038          LEN = 2
1039          CLASS = '     '
1040      END IF
1041      WRITE(LUPRI,'(/,/3X,A,I3,1X,3A)') 'For the requested',NVEC,
1042     &      TYPE(1:LEN),CLASS,' vectors'
1043      WRITE(LUPRI,'((8X,A,I3,A))')
1044     &   ' - ',N0KTRAN,  ' H matrix transformations ',
1045     &   ' - ',N0GTRAN,  ' G matrix transformations ',
1046     &   ' - ',N1GTRAN,  ' generalized G matrix transformations ',
1047     &   ' - ',(N1FTRAN+N2FTRAN),
1048     &                   ' generalized F matrix transformations ',
1049     &   ' - ',N0FATRAN, ' F{O} matrix transformations ',
1050     &   ' - ',N1FATRAN, ' generalized F{O} matrix transformations ',
1051     &   ' - ',NEATRAN,  ' generalized ETA{O} vector calculations '
1052      WRITE(LUPRI,'(3X,A,/,/)') 'will be performed.'
1053
1054
1055      IF (LOCDBG) THEN
1056
1057* H matrix transformations:
1058        IF (N0KTRAN.GT.0)
1059     &     WRITE(LUPRI,*)'List of H matrix transformations:'
1060        DO ITRAN = 1, N0KTRAN
1061          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1062     &     (I0KTRAN(I,ITRAN),I=1,4)
1063        END DO
1064        WRITE(LUPRI,*)
1065
1066* G matrix transformations:
1067        IF (N0GTRAN.GT.0)
1068     &       WRITE(LUPRI,*)'List of G matrix transformations:'
1069
1070        DO ITRAN = 1, N0GTRAN
1071          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1072     &     (I0GTRAN(I,ITRAN),I=1,4)
1073        END DO
1074        WRITE(LUPRI,*)
1075
1076        IF (N1GTRAN.GT.0)
1077     &      WRITE(LUPRI,*) 'List of (T^1 C) matrix transformations:'
1078        DO ITRAN = 1, N1GTRAN
1079          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1080     &     (I1GTRAN(I,ITRAN),I=1,4)
1081        END DO
1082        WRITE(LUPRI,*)
1083
1084* F matrix transformations:
1085        IF (N1FTRAN.GT.0)
1086     &      WRITE(LUPRI,*) 'List of (T^1 B) matrix transformations:'
1087        DO ITRAN = 1, N1FTRAN
1088          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1089     &     (I1FTRAN(I,ITRAN),I=1,3)
1090        END DO
1091        WRITE(LUPRI,*)
1092
1093        IF (N2FTRAN.GT.0)
1094     &    WRITE(LUPRI,*) 'List of (T^2 B) matrix transformations:'
1095        DO ITRAN = 1, N2FTRAN
1096          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1097     &     (I2FTRAN(I,ITRAN),I=1,3)
1098        END DO
1099        WRITE(LUPRI,*)
1100
1101* F{O} matrix transformations:
1102        IF (N0FATRAN.GT.0)
1103     &     WRITE(LUPRI,*) 'List of F{O} matrix transformations:'
1104        DO ITRAN = 1, N0FATRAN
1105          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1106     &     (I0FATRAN(I,ITRAN),I=1,4)
1107        END DO
1108        WRITE(*,*)
1109
1110        IF (N1FATRAN.GT.0)
1111     &      WRITE(LUPRI,*) 'List of (T^1 B{O}) matrix transformations:'
1112        DO ITRAN = 1, N1FATRAN
1113          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1114     &     (I1FATRAN(I,ITRAN),I=1,4)
1115        END DO
1116        WRITE(LUPRI,*)
1117
1118* ETA{O} vector calculations:
1119        IF (NEATRAN.GT.0)
1120     &    WRITE(LUPRI,*) 'List of (T^n A{O}) matrix transformations:'
1121        DO ITRAN = 1, NEATRAN
1122          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1123     &     (IEATRAN(I,ITRAN),I=1,5)
1124        END DO
1125        WRITE(LUPRI,*)
1126
1127      END IF
1128
1129      CALL QEXIT('CC_ETA_SETUP1')
1130      RETURN
1131      END
1132
1133*---------------------------------------------------------------------*
1134*              END OF SUBROUTINE CC_ETA_SETUP                         *
1135*---------------------------------------------------------------------*
1136