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 CCRHSVEC */
20*=====================================================================*
21      SUBROUTINE CCRHSVEC(TYPE,LABEL,ISYMS,ISTAT,EIGV,ISYMO,
22     &                    FREQS,LORX,ICAU,NVEC,MAXVEC,IOFFV,
23     &                    WORK,LWORK)
24*---------------------------------------------------------------------*
25*
26*    Purpose: calculate right hand side vectors for higher-order
27*             coupled cluster amplitude response equations,
28*             left and right excited state response equations
29*
30*             if called for ORDER=n, the solutions for ORDER=n-1
31*             must be available on file
32*
33*     implemented:  T:   ORDER = 1, 2, 3, 4
34*                   RE:  ORDER = 1, 2
35*                   CR:  ORDER = 1, 2
36*
37*    Written by Christof Haettig maj 1997, extension to RE july '97
38*                                          extension to CR march '98
39*                                          extension to O1 jan '99
40*    orb.-relax. or derivatives by Christof Haettig, Aug '99.
41*
42*=====================================================================*
43#if defined (IMPLICIT_NONE)
44      IMPLICIT NONE
45#else
46#include "implicit.h"
47#endif
48#include "priunit.h"
49#include "ccsdinp.h"
50#include "ccsdsym.h"
51#include "ccorb.h"
52#include "cclists.h"
53#include "dummy.h"
54Cholesky
55#include "maxorb.h"
56#include "ccdeco.h"
57Cholesky
58
59* local parameters:
60      CHARACTER*(18) MSGDBG
61      PARAMETER (MSGDBG = '[debug] CCRHSVEC> ')
62      LOGICAL LOCDBG
63      PARAMETER (LOCDBG = .FALSE. )
64
65
66
67      CHARACTER TYPE*(*), LISTR*3
68
69      INTEGER NVEC, MAXVEC, IOFFV, LWORK
70      INTEGER ISYMS(MAXVEC,*), ISYMO(MAXVEC,*)
71      INTEGER ISTAT(MAXVEC,*), ICAU(MAXVEC,*)
72      LOGICAL LORX(MAXVEC,*)
73
74      CHARACTER*8 LABEL(MAXVEC,*)
75
76#if defined (SYS_CRAY)
77      REAL FREQS(MAXVEC,*), EIGV(MAXVEC,*)
78      REAL WORK(LWORK)
79      REAL ZERO, RDUM
80      REAL XNORM, DDOT
81#else
82      DOUBLE PRECISION FREQS(MAXVEC,*), EIGV(MAXVEC,*)
83      DOUBLE PRECISION WORK(LWORK)
84      DOUBLE PRECISION ZERO, RDUM
85      DOUBLE PRECISION XNORM, DDOT
86#endif
87      PARAMETER (ZERO = 0.0d0)
88
89      CHARACTER MODEL*(10), MODELW*(10)
90      CHARACTER APROXR12*(3)
91      LOGICAL NEW_RHS
92      INTEGER IOPT, ISYM, IVEC, MPERM, NSTAT, ORDER, IDUM, IOPTE
93      INTEGER MXTRAN,MXDTRAN,MXCTRAN,MXBTRAN,MXBATRAN,MXAATRAN,MXXETRAN
94      INTEGER MXCATRAN
95      INTEGER KDTRAN,KCTRAN,KB1TRAN,KB2TRAN,KBA1TRAN,KAA1TRAN,KXETRAN
96      INTEGER NDTRAN,NCTRAN,NB1TRAN,NB2TRAN,NBA1TRAN,NAA1TRAN,NXETRAN
97      INTEGER KCATRAN, KBA2TRAN, KAA2TRAN
98      INTEGER NCATRAN, NBA2TRAN, NAA2TRAN
99      INTEGER KEND0, LEND0, LMAX1, LMAX2, KRHS1, KRHS2, KEND1, LEND1
100      INTEGER KLHS1, KLHS2, KEND2, LEND2, IDXVEC
101      INTEGER KRHSR12, LMAXR12, IOPTR12, MODLEN
102
103* external functions:
104      INTEGER ILSTSYM, ILRCAMP
105
106
107*---------------------------------------------------------------------*
108* check number of required rhs vectors, if zero return immediatly:
109*---------------------------------------------------------------------*
110      IF (NVEC.EQ.0) RETURN
111
112*---------------------------------------------------------------------*
113* print header for rhs vector section
114*---------------------------------------------------------------------*
115      WRITE (LUPRI,'(7(/1X,2A),/)')
116     & '------------------------------------',
117     &                               '-------------------------------',
118     & '|          OUTPUT FROM AMPLITUDE RHS',
119     &                               ' VECTOR SECTION               |',
120     & '------------------------------------',
121     &                               '-------------------------------'
122      CALL FLSHFO(LUPRI)
123
124*---------------------------------------------------------------------*
125      IF (.NOT. (CCS .OR. CC2 .OR. CCSD .OR. CC3) ) THEN
126         CALL QUIT('CCRHSVEC called for unknown Coupled Cluster model.')
127      END IF
128
129      NEW_RHS = .FALSE.
130
131      IF (TYPE(1:3).EQ.'O1 ') THEN
132        ORDER = 1
133        NSTAT = 0
134        MPERM = 1
135      ELSE IF (TYPE(1:2).EQ.'O2') THEN
136        ORDER = 2
137        NSTAT = 0
138        MPERM = 2
139        ! compute complete O2 vector in B matrix module
140        NEW_RHS = .TRUE.
141      ELSE IF (TYPE(1:2).EQ.'O3') THEN
142        ORDER = 3
143        NSTAT = 0
144        MPERM = 3
145      ELSE IF (TYPE(1:2).EQ.'O4') THEN
146        ORDER = 4
147        NSTAT = 0
148        MPERM = 12
149        WRITE (LUPRI,*) 'warning: rhs vectors ',TYPE(1:2),
150     &       ' not tested!!!.'
151      ELSE IF (TYPE(1:3).EQ.'EO1') THEN
152        ORDER = 1
153        NSTAT = 1
154        MPERM = 1
155        ! compute complete EO1 vector in B matrix module
156        NEW_RHS = .TRUE.
157      ELSE IF (TYPE(1:3).EQ.'EO2') THEN
158        ORDER = 2
159        NSTAT = 1
160        MPERM = 2
161        WRITE (LUPRI,*) 'warning: rhs vectors ',TYPE(1:3),
162     &       ' not tested!!!.'
163      ELSE IF (TYPE(1:3).EQ.'CO1') THEN
164        ORDER = 1
165        NSTAT = 0
166        MPERM = 1
167      ELSE IF (TYPE(1:3).EQ.'CO2') THEN
168        ORDER = 2
169        NSTAT = 0
170        MPERM = 2
171      ELSE
172        WRITE (LUPRI,*) 'rhs vectors ',TYPE(1:2),' not implemented.'
173        CALL QUIT('required rhs vectors not implemented.')
174      END IF
175
176* Cholesky check: only CC2 O1 has been implemented
177      IF (CHOINT .AND. CC2) THEN
178         IF (TYPE(1:2).NE.'O1') THEN
179            WRITE (LUPRI,*)
180     &      'rhs vectors ',TYPE(1:2),' not implemented for Cholesky.'
181            CALL QUIT('required rhs vectors not implemented.')
182         ENDIF
183      ENDIF
184
185* print some debug/info output
186      IF (IPRINT .GT. 10 .OR. LOCDBG) THEN
187        WRITE(LUPRI,*) 'CCRHSVEC Workspace:',LWORK
188      END IF
189
190*---------------------------------------------------------------------*
191* allocate & initialize work space for lists
192*---------------------------------------------------------------------*
193
194      MXTRAN   = MPERM * NVEC
195
196      MXDTRAN  = MXDIM_DTRAN  * MXTRAN
197      MXCTRAN  = MXDIM_CTRAN  * MXTRAN
198      MXBTRAN  = MXDIM_BTRAN  * MXTRAN
199      MXBATRAN = MXDIM_BATRAN * MXTRAN
200      MXCATRAN = MXDIM_CATRAN * MXTRAN
201      MXAATRAN = MXDIM_AATRAN * MXTRAN
202      MXXETRAN = MXDIM_XEVEC  * MXTRAN
203
204      KDTRAN   = 1
205      KCTRAN   = KDTRAN   + MXDTRAN
206      KB1TRAN  = KCTRAN   + MXCTRAN
207      KB2TRAN  = KB1TRAN  + MXBTRAN
208      KCATRAN  = KB2TRAN  + MXBTRAN
209      KBA1TRAN = KCATRAN  + MXCATRAN
210      KBA2TRAN = KBA1TRAN + MXBATRAN
211      KAA1TRAN = KBA2TRAN + MXBATRAN
212      KAA2TRAN = KAA1TRAN + MXAATRAN
213      KXETRAN  = KAA2TRAN + MXAATRAN
214      KEND0    = KXETRAN  + MXXETRAN
215      LEND0    = LWORK    - KEND0
216
217
218      IF (LEND0 .LT. 0 ) THEN
219        WRITE (LUPRI,*) 'Insufficient work space in CCRHSVEC.'
220        WRITE (LUPRI,*) 'KEND0, LEND0, LWORK:',KEND0,LEND0,LWORK
221        WRITE (LUPRI,*) 'MXTRAN:',MXTRAN
222        CALL QUIT('Insufficient work space in CCRHSVEC.')
223      END IF
224
225*---------------------------------------------------------------------*
226* set up lists for D, C, B, B{O} and A{O} transformations:
227*---------------------------------------------------------------------*
228      CALL CC_RHS_SETUP(TYPE,NSTAT,ORDER,LABEL,ISTAT,EIGV,ISYMO,FREQS,
229     &                  LORX, ICAU, NVEC, MAXVEC, IOFFV, MXTRAN,
230     &                  NEW_RHS,
231     &                  WORK(KDTRAN),  NDTRAN,
232     &                  WORK(KCTRAN),  NCTRAN,
233     &                  WORK(KB1TRAN), NB1TRAN,
234     &                  WORK(KB2TRAN), NB2TRAN,
235     &                  WORK(KCATRAN), NCATRAN,
236     &                  WORK(KBA1TRAN),NBA1TRAN,
237     &                  WORK(KBA2TRAN),NBA2TRAN,
238     &                  WORK(KAA1TRAN),NAA1TRAN,
239     &                  WORK(KAA2TRAN),NAA2TRAN,
240     &                  WORK(KXETRAN), NXETRAN )
241
242*---------------------------------------------------------------------*
243* initialize rhs vector files:
244*---------------------------------------------------------------------*
245      IF( TYPE(1:3).NE.'O1 ' .AND. TYPE(1:3).NE.'X1 ' .and.
246     &    TYPE(1:3).NE.'CO1'                                ) THEN
247
248        LMAX1 = 0
249        LMAX2 = 0
250        LMAXR12 = 0
251        DO ISYM = 1, NSYM
252          LMAX1 = MAX(LMAX1,NT1AM(ISYM))
253          LMAX2 = MAX(LMAX2,NT2AM(ISYM))
254          IF (CCR12) LMAXR12 = MAX(LMAXR12,NTR12AM(ISYM))
255        END DO
256
257        KRHS1 = KEND0
258        KRHS2 = KRHS1 + LMAX1
259        KRHSR12 = KRHS2 + LMAX2
260        KEND1 = KRHSR12 + LMAXR12
261        LEND1 = LWORK - KEND1
262
263        IF (LEND1 .LT. 0 ) THEN
264          WRITE (LUPRI,*) 'Insufficient work space in CCRHSVEC. (2)'
265          WRITE (LUPRI,*) 'KEND1, LEND1, LWORK:',KEND1,LEND1,LWORK
266          CALL QUIT('Insufficient work space in CCRHSVEC. (2)')
267        END IF
268
269        CALL DZERO(WORK(KRHS1),LMAX1)
270        IF (.NOT.CCS) CALL DZERO(WORK(KRHS2),LMAX2)
271        IF (CCR12) CALL DZERO(WORK(KRHSR12),LMAXR12)
272
273        IF (CCS) THEN
274           MODEL = 'CCS       '
275           IOPT  = 1
276        ELSE IF (CC2) THEN
277           MODEL = 'CC2       '
278           IOPT  = 3
279        ELSE IF (CCSD) THEN
280           MODEL = 'CCSD      '
281           IOPT  = 3
282        ELSE IF (CC3) THEN
283           MODEL = 'CC3       '
284           ! intialize usual and effective rhs vector
285           IOPT  = 3
286           IOPTE = 24
287        ELSE
288           CALL QUIT('Unknown coupled cluster model in CCRHSVEC.')
289        END IF
290        IF (CCR12) THEN
291          APROXR12 = '   '
292          IOPTR12 = 32
293        END IF
294        CALL CCSD_MODEL(MODELW,MODLEN,10,MODEL,10,APROXR12)
295
296        DO IVEC = IOFFV+1, IOFFV+NVEC
297         ISYM = ILSTSYM(TYPE,IVEC)
298         CALL CC_WRRSP(TYPE,IVEC,ISYM,IOPT,MODELW,IDUMMY,
299     &                 WORK(KRHS1),WORK(KRHS2),WORK(KEND1),LEND1)
300         IF (CCR12) THEN
301           CALL CC_WRRSP(TYPE,IVEC,ISYM,IOPTR12,MODELW,IDUMMY,
302     &                   IDUMMY,WORK(KRHSR12),WORK(KEND1),LEND1)
303         END IF
304         IF (CCSDT) THEN
305           CALL CC_WRRSP(TYPE,IVEC,ISYM,IOPTE,MODELW,IDUMMY,
306     &                   WORK(KRHS1),WORK(KRHS2),WORK(KEND1),LEND1)
307           CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPTE,MODEL,
308     &                   WORK(KRHS1),WORK(KRHS2))
309         END IF
310        END DO
311
312      END IF
313
314*---------------------------------------------------------------------*
315* calculate D matrix contributions:
316*---------------------------------------------------------------------*
317      IF (TYPE(1:2).EQ.'O4') THEN
318        IOPT = 4
319        CALL CC_DMAT(WORK(KDTRAN),NDTRAN,
320     &               'R1 ','R1 ','R1 ','R1 ',IOPT,TYPE,
321     &               IDUM, RDUM, 0, WORK(KEND0), LEND0 )
322      END IF
323
324      IF (LOCDBG .AND. TYPE(1:3).NE.'O1 '.AND. TYPE(1:3).NE.'X1 '
325     &           .AND. TYPE(1:3).NE.'CO1'                        ) THEN
326        WRITE (LUPRI,*) MSGDBG,
327     &        'NORM^2 of RHS vectors after D matrix terms:'
328        DO IVEC = IOFFV+1, IOFFV+NVEC
329          IOPT = 3
330          ISYM = ILSTSYM(TYPE,IVEC)
331          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
332     &                  WORK(KRHS1),WORK(KRHS2))
333          XNORM = DDOT(NT1AM(ISYM),WORK(KRHS1),1,WORK(KRHS1),1)
334          IF (.NOT. CCS)
335     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KRHS2),1,WORK(KRHS2),1)
336          WRITE (LUPRI,*) MSGDBG, IVEC,XNORM
337        END DO
338      END IF
339
340
341*---------------------------------------------------------------------*
342* calculate C matrix contributions:
343*---------------------------------------------------------------------*
344      IF      (TYPE(1:2).EQ.'O4') THEN
345        IOPT = 4
346        CALL CC_CMAT(WORK(KCTRAN),NCTRAN,'R2 ','R1 ','R1 ',IOPT,TYPE,
347     &               IDUM, RDUM, 0, WORK(KEND0), LEND0 )
348      ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
349        IOPT = 4
350        CALL CC_CMAT(WORK(KCTRAN),NCTRAN,'R1 ','R1 ','R1 ',IOPT,TYPE,
351     &               IDUM, RDUM, 0, WORK(KEND0), LEND0 )
352      ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
353        IOPT = 4
354        CALL CC_CMAT(WORK(KCTRAN),NCTRAN,'R1 ','R1 ','RE ',IOPT,TYPE,
355     &               IDUM, RDUM, 0, WORK(KEND0), LEND0 )
356      END IF
357
358      IF (LOCDBG .AND. TYPE(1:3).NE.'O1 '.AND. TYPE(1:3).NE.'X1 '
359     &           .AND. TYPE(1:3).NE.'CO1'                        ) THEN
360        WRITE (LUPRI,*) MSGDBG,
361     &        'NORM^2 of RHS vectors after C matrix terms:'
362        DO IVEC = IOFFV+1, IOFFV+NVEC
363          IOPT = 3
364          IF (CC3) IOPT = 24
365          ISYM = ILSTSYM(TYPE,IVEC)
366          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
367     &                  WORK(KRHS1),WORK(KRHS2))
368          XNORM = DDOT(NT1AM(ISYM),WORK(KRHS1),1,WORK(KRHS1),1)
369          IF (.NOT. CCS)
370     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KRHS2),1,WORK(KRHS2),1)
371          WRITE (LUPRI,*) MSGDBG, IVEC,XNORM
372        END DO
373      END IF
374
375*---------------------------------------------------------------------*
376* calculate B matrix contributions:
377*---------------------------------------------------------------------*
378      IF      ( TYPE(1:2).EQ.'O4' ) THEN
379        IOPT = 4
380        CALL CC_BMAT(WORK(KB1TRAN), NB1TRAN,'R3 ','R1 ',IOPT,TYPE,
381     &               IDUM, RDUM, 0, .FALSE.,WORK(KEND0), LEND0  )
382        IOPT = 4
383        CALL CC_BMAT(WORK(KB2TRAN), NB2TRAN,'R2 ','R2 ',IOPT,TYPE,
384     &               IDUM, RDUM, 0, .FALSE.,WORK(KEND0), LEND0 )
385      ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
386        IOPT = 4
387        CALL CC_BMAT(WORK(KB1TRAN), NB1TRAN,'R2 ','R1 ',IOPT,TYPE,
388     &               IDUM, RDUM, 0, .FALSE.,WORK(KEND0), LEND0 )
389      ELSE IF ( TYPE(1:2).EQ.'O2' ) THEN
390        IOPT = 4
391        CALL CC_BMAT(WORK(KB1TRAN), NB1TRAN,'R1 ','R1 ',IOPT,TYPE,
392     &               IDUM, RDUM, 0, NEW_RHS,WORK(KEND0), LEND0 )
393      ELSE IF ( TYPE(1:3).EQ.'EO2') THEN
394        IOPT = 4
395        CALL CC_BMAT(WORK(KB1TRAN), NB1TRAN,'R2 ','RE ',IOPT,TYPE,
396     &               IDUM, RDUM, 0, .FALSE.,WORK(KEND0), LEND0 )
397        IOPT = 4
398        CALL CC_BMAT(WORK(KB2TRAN), NB2TRAN,'R1 ','ER1',IOPT,TYPE,
399     &               IDUM, RDUM, 0, .FALSE.,WORK(KEND0), LEND0 )
400      ELSE IF ( TYPE(1:3).EQ.'EO1') THEN
401        IOPT = 4
402        CALL CC_BMAT(WORK(KB1TRAN), NB1TRAN,'R1 ','RE ',IOPT,TYPE,
403     &               IDUM, RDUM, 0, NEW_RHS,WORK(KEND0), LEND0 )
404      ELSE IF ( TYPE(1:3).EQ.'CO2' ) THEN
405        IOPT = 4
406        CALL CC_BMATRIX(WORK(KB1TRAN), NB1TRAN,'RC ','RC ',IOPT,TYPE,
407     &                  IDUM, RDUM, 0, .FALSE.,WORK(KEND0), LEND0 )
408      END IF
409
410      IF (LOCDBG .AND. TYPE(1:3).NE.'O1 '.AND. TYPE(1:3).NE.'X1 '
411     &           .AND. TYPE(1:3).NE.'CO1'                        ) THEN
412        WRITE (LUPRI,*) MSGDBG,
413     &        'NORM^2 of RHS vectors after B matrix terms:'
414        DO IVEC = IOFFV+1, IOFFV+NVEC
415          IOPT = 3
416          IF (CC3) IOPT = 24
417          ISYM = ILSTSYM(TYPE,IVEC)
418          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
419     &                  WORK(KRHS1),WORK(KRHS2))
420          XNORM = DDOT(NT1AM(ISYM),WORK(KRHS1),1,WORK(KRHS1),1)
421          IF (.NOT. CCS)
422     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KRHS2),1,WORK(KRHS2),1)
423          WRITE (LUPRI,*) MSGDBG, IVEC,XNORM
424        END DO
425      END IF
426
427*---------------------------------------------------------------------*
428* calculate C{O} matrix contributions:
429*---------------------------------------------------------------------*
430      IF      ( TYPE(1:2).EQ.'O4' .AND. NCATRAN.NE.0) THEN
431        IOPT = 4
432c       CALL CC_CAMAT(WORK(KCATRAN),NCATRAN,'o1 ','R1 ','R1 ','R1 ',
433c    &                IOPT, TYPE, IDUM, RDUM, 0, WORK(KEND0), LEND0 )
434        CALL QUIT('cc_camat routine not yet implememted.')
435      END IF
436
437*---------------------------------------------------------------------*
438* calculate B{O} matrix contributions:
439*---------------------------------------------------------------------*
440      IF      ( TYPE(1:2).EQ.'O4' ) THEN
441        IOPT = 4
442        CALL CC_BAMAT(WORK(KBA1TRAN),NBA1TRAN,'o1 ','R2 ','R1 ',IOPT,
443     &                TYPE, IDUM, RDUM, 0,WORK(KEND0), LEND0 )
444        IOPT = 4
445        CALL CC_BAMAT(WORK(KBA2TRAN),NBA2TRAN,'o2 ','R1 ','R1 ',IOPT,
446     &                TYPE, IDUM, RDUM, 0,WORK(KEND0), LEND0 )
447      ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
448        IOPT = 4
449        CALL CC_BAMAT(WORK(KBA1TRAN),NBA1TRAN,'o1 ','R1 ','R1 ',IOPT,
450     &                TYPE, IDUM, RDUM, 0,WORK(KEND0), LEND0 )
451      ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
452        IOPT = 4
453        CALL CC_BAMAT(WORK(KBA1TRAN),NBA1TRAN,'o1 ','R1 ','RE ',IOPT,
454     &                TYPE, IDUM, RDUM, 0,WORK(KEND0), LEND0 )
455      END IF
456
457      IF (LOCDBG .AND. TYPE(1:3).NE.'O1 '.AND. TYPE(1:3).NE.'X1 '
458     &           .AND. TYPE(1:3).NE.'CO1'                        ) THEN
459        WRITE (LUPRI,*) MSGDBG,
460     &        'NORM^2 of RHS vectors after B{O} matrix terms:'
461        DO IVEC = IOFFV+1, IOFFV+NVEC
462          IOPT = 3
463          IF (CC3) IOPT = 24
464          ISYM = ILSTSYM(TYPE,IVEC)
465          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
466     &                  WORK(KRHS1),WORK(KRHS2))
467          XNORM = DDOT(NT1AM(ISYM),WORK(KRHS1),1,WORK(KRHS1),1)
468          IF (.NOT. CCS)
469     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KRHS2),1,WORK(KRHS2),1)
470          WRITE (LUPRI,*) MSGDBG, IVEC,XNORM
471        END DO
472      END IF
473
474*---------------------------------------------------------------------*
475* calculate A{O} matrix contributions:
476*---------------------------------------------------------------------*
477      IF      ( TYPE(1:2).EQ.'O4' ) THEN
478        IOPT = 4
479        CALL CC_AAMAT(WORK(KAA1TRAN),NAA1TRAN,'o1 ','R3 ',IOPT,TYPE,
480     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
481        IOPT = 4
482        CALL CC_AAMAT(WORK(KAA2TRAN),NAA2TRAN,'o2 ','R2 ',IOPT,TYPE,
483     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
484      ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
485        IOPT = 4
486        CALL CC_AAMAT(WORK(KAA1TRAN),NAA1TRAN,'o1 ','R2 ',IOPT,TYPE,
487     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
488        IOPT = 4
489        CALL CC_AAMAT(WORK(KAA2TRAN),NAA2TRAN,'o2 ','R1 ',IOPT,TYPE,
490     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
491      ELSE IF ( TYPE(1:2).EQ.'O2' .AND. (.NOT.NEW_RHS)) THEN
492        IOPT = 4
493        CALL CC_AAMAT(WORK(KAA1TRAN),NAA1TRAN,'o1 ','R1 ',IOPT,TYPE,
494     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
495      ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
496        IOPT = 4
497        CALL CC_AAMAT(WORK(KAA1TRAN),NAA1TRAN,'o1 ','ER1',IOPT,TYPE,
498     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
499        IOPT = 4
500        CALL CC_AAMAT(WORK(KAA2TRAN),NAA2TRAN,'o2 ','RE ',IOPT,TYPE,
501     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
502      ELSE IF ( TYPE(1:3).EQ.'EO1' .AND. (.NOT.NEW_RHS)) THEN
503        IOPT = 4
504        CALL CC_AAMAT(WORK(KAA1TRAN),NAA1TRAN,'o1 ','RE ',IOPT,TYPE,
505     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
506      ELSE IF ( TYPE(1:3).EQ.'CO2' ) THEN
507        IOPT = 4
508        CALL CC_AAMAT(WORK(KAA1TRAN),NAA1TRAN,'o1 ','RC ',IOPT,TYPE,
509     &                IDUMMY,DUMMY,1,WORK(KEND0), LEND0 )
510      END IF
511
512      IF (LOCDBG .AND. TYPE(1:3).NE.'O1 '.AND. TYPE(1:3).NE.'X1 '
513     &           .AND. TYPE(1:3).NE.'CO1'                        ) THEN
514        WRITE (LUPRI,*) MSGDBG,
515     &        'NORM^2 of RHS vectors after A{O} matrix terms:'
516        DO IVEC = IOFFV+1, IOFFV+NVEC
517          IOPT = 3
518          IF (CC3) IOPT = 24
519          ISYM = ILSTSYM(TYPE,IVEC)
520          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
521     &                  WORK(KRHS1),WORK(KRHS2))
522          XNORM = DDOT(NT1AM(ISYM),WORK(KRHS1),1,WORK(KRHS1),1)
523          IF (.NOT. CCS)
524     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KRHS2),1,WORK(KRHS2),1)
525          WRITE (LUPRI,*) MSGDBG, IVEC,XNORM
526        END DO
527      END IF
528
529*---------------------------------------------------------------------*
530* calculate Xi{O} vector contributions:
531*---------------------------------------------------------------------*
532
533Cholesky
534*
535*   the Cholesky routine is *far* less general, hence most of the
536*   input is implicit: all that's needed are the operator indices
537*   in KXETRAN.
538
539      IF (CHOINT .AND. TYPE(1:2).EQ.'O1') THEN
540         CALL CC_CHOXI0(WORK(KXETRAN),NXETRAN,WORK(KEND0),LEND0)
541         GOTO 1234
542      END IF
543Cholesky
544
545      IF ( TYPE(1:3).EQ.'O1 '.OR. TYPE(1:3).EQ.'X1 ') THEN
546        IOPT = 3
547        CALL CC_XIETA(WORK(KXETRAN),NXETRAN,IOPT, ORDER, 'L0 ',
548     &                'O1 ', IDUM, RDUM, 'X1 ', IDUM, RDUM,
549     &                .FALSE.,0, WORK(KEND0),LEND0)
550      ELSE IF ( TYPE(1:2).EQ.'O2' .OR. TYPE(1:2).EQ.'X2' ) THEN
551        IOPT = 3
552        CALL CC_XIETA(WORK(KXETRAN),NXETRAN,IOPT, ORDER, 'L0 ',
553     &                'O2 ', IDUM, RDUM, 'X2 ', IDUM, RDUM,
554     &                .FALSE.,0, WORK(KEND0),LEND0)
555      ELSE IF ( TYPE(1:3).EQ.'CO1' ) THEN
556        IOPT = 3
557        CALL CC_XIETA(WORK(KXETRAN),NXETRAN,IOPT, ORDER, 'L0 ',
558     &                'RC ', IDUM, RDUM, '---', IDUM, RDUM,
559     &                .TRUE.,0, WORK(KEND0),LEND0)
560      END IF
561
562 1234 CONTINUE     ! From Cholesky
563
564      IF (LOCDBG) THEN
565        LMAX1 = 0
566        LMAX2 = 0
567        LMAXR12 = 0
568        DO ISYM = 1, NSYM
569          LMAX1 = MAX(LMAX1,NT1AM(ISYM))
570          LMAX2 = MAX(LMAX2,NT2AM(ISYM))
571          IF (CCR12) THEN
572            LMAXR12 = MAX(LMAXR12,NTR12AM(ISYM))
573          END IF
574        END DO
575
576        KRHS1 = KEND0
577        KRHS2 = KRHS1 + LMAX1
578        KRHSR12 = KRHS2 + LMAX2
579        KEND1 = KRHSR12 + LMAXR12
580        LEND1 = LWORK - KEND1
581
582        IF (LEND1 .LT. 0 ) THEN
583          WRITE (LUPRI,*) 'Insufficient work space in CCRHSVEC. (3)'
584          WRITE (LUPRI,*) 'KEND1, LEND1, LWORK:',KEND1,LEND1,LWORK
585          CALL QUIT('Insufficient work space in CCRHSVEC. (3)')
586        END IF
587
588        WRITE (LUPRI,*) MSGDBG,
589     &       'NORM^2 of RHS vectors after Xi{O} matrix terms:'
590        LISTR = TYPE(1:3)
591        IF (TYPE.EQ.'CO1') LISTR ='RC '
592        DO IVEC = IOFFV+1, IOFFV+NVEC
593          IOPT = 3
594          IF (CC3) IOPT = 24
595          ISYM = ILSTSYM(LISTR,IVEC)
596          IF (TYPE.EQ.'CO1') THEN
597           WRITE(LUPRI,*) 'Cauchy order:',ICAU(IVEC,1)
598           IDXVEC=ILRCAMP(LABEL(IVEC,1),ICAU(IVEC,1)-1,ISYM)
599          ELSE
600           IDXVEC= IVEC
601          END IF
602          CALL CC_RDRSP(LISTR,IDXVEC,ISYM,IOPT,MODEL,
603     &                  WORK(KRHS1),WORK(KRHS2))
604          XNORM = DDOT(NT1AM(ISYM),WORK(KRHS1),1,WORK(KRHS1),1)
605Chol      IF (.NOT. CCS)
606          IF ((.NOT. CCS) .AND. (.NOT. (CHOINT.AND.CC2)))
607     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KRHS2),1,WORK(KRHS2),1)
608          IF (CCR12) THEN
609            CALL CC_RDRSP(LISTR,IDXVEC,ISYM,IOPTR12,MODEL,
610     &                  DUMMY,WORK(KRHSR12))
611            XNORM = XNORM + DDOT(NTR12AM(ISYM),WORK(KRHSR12),1,
612     &                           WORK(KRHSR12),1)
613          END IF
614          WRITE (LUPRI,*) MSGDBG, IVEC,XNORM
615        END DO
616      END IF
617
618*---------------------------------------------------------------------*
619* test (static) EO1 vectors by calculating the excited state FOP's
620*---------------------------------------------------------------------*
621      IF (LOCDBG .AND. TYPE(1:3).EQ.'EO1') THEN
622        LMAX1 = 0
623        LMAX2 = 0
624        DO ISYM = 1, NSYM
625          LMAX1 = MAX(LMAX1,NT1AM(ISYM))
626          LMAX2 = MAX(LMAX2,NT2AM(ISYM))
627        END DO
628
629        KLHS1 = KEND1
630        KLHS2 = KLHS1 + LMAX1
631        KEND2 = KLHS2 + LMAX2
632        LEND2 = LWORK - KEND2
633
634        IF (LEND2 .LT. 0 ) THEN
635          WRITE (LUPRI,*) 'Insufficient work space in CCRHSVEC. (4)'
636          WRITE (LUPRI,*) 'KEND2, LEND2, LWORK:',KEND2,LEND2,LWORK
637          CALL QUIT('Insufficient work space in CCRHSVEC. (4)')
638        END IF
639
640        WRITE (LUPRI,*) MSGDBG, 'excited state first order properties:'
641        DO IVEC = IOFFV+1, IOFFV+NVEC
642        IF (ISYMO(IVEC,1).EQ.1 .AND. FREQS(IVEC,1).EQ.ZERO) THEN
643          IOPT = 3
644          ISYM = ILSTSYM(TYPE,IVEC)
645          CALL CC_RDRSP(TYPE,IVEC,ISYM,IOPT,MODEL,
646     &                  WORK(KRHS1),WORK(KRHS2))
647          CALL CCLR_DIASCL(WORK(KRHS2),0.5d0,ISYM)
648          CALL CC_RDRSP('LE',ISTAT(IVEC,1),ISYMS(IVEC,1),IOPT,MODEL,
649     &                  WORK(KLHS1),WORK(KLHS2))
650          XNORM = DDOT(NT1AM(ISYM),WORK(KLHS1),1,WORK(KRHS1),1)
651          IF (.NOT. CCS)
652     &     XNORM = XNORM+DDOT(NT2AM(ISYM),WORK(KLHS2),1,WORK(KRHS2),1)
653          WRITE (LUPRI,'(A,I3,2X,F12.8,2X,A,2X,F12.8)') MSGDBG,
654     &              ISTAT(IVEC,1),EIGV(IVEC,1),LABEL(IVEC,1),XNORM
655        ELSE
656          WRITE (LUPRI,'(A,I3,2X,F12.8,2X,A,2X,F12.8)') MSGDBG,
657     &              ISTAT(IVEC,1),EIGV(IVEC,1),LABEL(IVEC,1),ZERO
658        END IF
659        END DO
660      END IF
661
662*---------------------------------------------------------------------*
663* that's it:
664*---------------------------------------------------------------------*
665
666      RETURN
667      END
668
669*=====================================================================*
670*              END OF SUBROUTINE CCRHSVEC                             *
671*=====================================================================*
672c /* deck CC_RHS_SETUP */
673*=====================================================================*
674      SUBROUTINE CC_RHS_SETUP(TYPE,NSTAT,ORDER,LAB,ISTAT,
675     &                        EIGV,ISYMO,FREQ,LORX,ICAU,
676     &                        NVEC,MAXVEC,IOFFV,MXTRAN,
677     &                        NEW_RHS,
678     &                        IDTRAN,  NDTRAN,
679     &                        ICTRAN,  NCTRAN,
680     &                        IB1TRAN, NB1TRAN,
681     &                        IB2TRAN, NB2TRAN,
682     &                        ICATRAN, NCATRAN,
683     &                        IBA1TRAN,NBA1TRAN,
684     &                        IBA2TRAN,NBA2TRAN,
685     &                        IAA1TRAN,NAA1TRAN,
686     &                        IAA2TRAN,NAA2TRAN,
687     &                        IXETRAN, NXETRAN )
688*---------------------------------------------------------------------*
689*
690*    Purpose: set up for CCRHSVEC section
691*                - list of D matrix transformations
692*                - list of C matrix transformations
693*                - list of B matrix transformations
694*                - list of B{O} matrix transformations
695*                - list of A{O} matrix transformations
696*                - list of Xi{O} vector calculations
697*
698*     Written by Christof Haettig, maj 1997.
699*     O1, O2, O3, O4, EO1, EO2 with one perturbation including
700*     orb.-relax. or derivatives by Christof Haettig, Aug '99.
701*
702*=====================================================================*
703#if defined (IMPLICIT_NONE)
704      IMPLICIT NONE
705#else
706#  include "implicit.h"
707#endif
708#include "priunit.h"
709#include "ccroper.h"
710#include "cclists.h"
711
712* local parameters:
713      CHARACTER*(22) MSGDBG
714      PARAMETER (MSGDBG = '[debug] CC_RHS_SETUP> ')
715      LOGICAL LOCDBG
716      PARAMETER (LOCDBG = .FALSE.)
717      INTEGER MXORD, MXORD2, MXORD3
718      PARAMETER (MXORD  = 4)
719      PARAMETER (MXORD2 = MXORD *(MXORD-1)/2 )
720      PARAMETER (MXORD3 = MXORD2*(MXORD-2)/3 )
721
722      INTEGER MXTRAN, NSTAT, ORDER, MAXVEC, NVEC, IOFFV
723
724      CHARACTER*(*) TYPE
725
726      CHARACTER*(8) LAB(MAXVEC,*)
727      INTEGER ISTAT(MAXVEC,*), ICAU(MAXVEC,*), ISYMO(MAXVEC,*)
728      LOGICAL LORX(MAXVEC,*), NEW_RHS
729
730#if defined (SYS_CRAY)
731      REAL FREQ(MAXVEC,*), EIGV(MAXVEC)
732#else
733      DOUBLE PRECISION FREQ(MAXVEC,*), EIGV(MAXVEC)
734#endif
735
736      INTEGER IDTRAN( MXDIM_DTRAN, MXTRAN)
737      INTEGER ICTRAN( MXDIM_CTRAN, MXTRAN)
738      INTEGER IB1TRAN(MXDIM_BTRAN, MXTRAN)
739      INTEGER IB2TRAN(MXDIM_BTRAN, MXTRAN)
740      INTEGER ICATRAN(MXDIM_CATRAN,MXTRAN)
741      INTEGER IBA1TRAN(MXDIM_BATRAN,MXTRAN)
742      INTEGER IBA2TRAN(MXDIM_BATRAN,MXTRAN)
743      INTEGER IAA1TRAN(MXDIM_AATRAN,MXTRAN)
744      INTEGER IAA2TRAN(MXDIM_AATRAN,MXTRAN)
745      INTEGER IXETRAN(MXDIM_XEVEC,MXTRAN)
746
747      INTEGER NDTRAN,NCTRAN,NB1TRAN,NB2TRAN,NBA1TRAN,NAA1TRAN,NXETRAN
748      INTEGER NBA2TRAN,NAA2TRAN,NCATRAN
749
750      INTEGER IOP(MXORD), IOP2(MXORD2)
751      INTEGER IR1(MXORD), IR2(MXORD2), IR3(MXORD3)
752      INTEGER IEX, IE1(MXORD), IET1(MXORD), IET2(MXORD2)
753      INTEGER ISYMS, ISYM(MXORD), IRELAX(MXORD)
754
755      INTEGER A, B, C, D
756      PARAMETER (A = 1, B = 2, C = 3, D = 4)
757      INTEGER AB, AC, AD, BC, BD, CD
758      PARAMETER (AB = 1, AC = 2, BC = 3, AD = 4, BD = 5, CD = 6)
759      INTEGER ABC, ABD, ACD, BCD
760      PARAMETER (ABC = 1, ABD = 2, ACD = 3, BCD = 4)
761
762
763      INTEGER NS2A, NS3A, NS4A, NP3AB, NP4AB, NT4ABC
764      PARAMETER (NS2A = 2, NS3A = 3, NS4A = 4)
765      PARAMETER (NP3AB = 3, NP4AB = 6, NT4ABC = 4)
766
767      INTEGER ISA(NS4A), ISB(NS4A), ISC(NS4A), ISD(NS4A)
768      INTEGER IPAB(NP4AB), IPCD(NP4AB)
769      INTEGER IPA(NP4AB), IPB(NP4AB), IPC(NP4AB), IPD(NP4AB)
770      INTEGER ITABC(NT4ABC), ITD(NT4ABC)
771
772      DATA ISA  / A, B, C, D/
773      DATA ISB  / B, A, A, A/
774      DATA ISC  / C, C, B, B/
775      DATA ISD  / D, D, D, C/
776
777      DATA IPAB / AB, AC, BC, AD, BD, CD /
778      DATA IPA  / A,  A,  B,  A,  B,  C  /
779      DATA IPB  / B,  C,  C,  D,  D,  D  /
780      DATA IPCD / CD, BD, AD, BC, AC, AB /
781      DATA IPC  / C,  B,  A,  B,  A,  A  /
782      DATA IPD  / D,  D,  D,  C,  C,  B  /
783
784      DATA ITABC / ABC, ABD, ACD, BCD /
785      DATA ITD   / D,   C,   B,   A   /
786
787      CHARACTER*8 LABSOP
788      INTEGER IDXA, IDXB, IDXC, ITRAN, IDX, IDXAB, IDXABC, IVEC
789      INTEGER NRELAX, ISGNSOP, ISYSOP, INUM
790
791* external functions:
792      INTEGER IROPER
793      INTEGER IROPER2
794      INTEGER IETA1
795      INTEGER ICHI2
796      INTEGER IR1KAPPA
797      INTEGER IR1TAMP
798      INTEGER IR2TAMP
799      INTEGER IR3TAMP
800      INTEGER IER1AMP
801      INTEGER IER2AMP
802      INTEGER ILRCAMP
803      INTEGER ICR2AMP
804
805
806*---------------------------------------------------------------------*
807* initializations:
808*---------------------------------------------------------------------*
809      NDTRAN   = 0
810      NCTRAN   = 0
811      NB1TRAN  = 0
812      NB2TRAN  = 0
813      NCATRAN  = 0
814      NBA1TRAN = 0
815      NBA2TRAN = 0
816      NAA1TRAN = 0
817      NAA2TRAN = 0
818      NXETRAN  = 0
819
820*---------------------------------------------------------------------*
821* start loop over all requested rhs-vectors:
822*---------------------------------------------------------------------*
823
824      DO IVEC = IOFFV+1, IOFFV+NVEC
825
826* eigenvectors that contribute:
827      IF (NSTAT.EQ.1) THEN
828        IEX = ISTAT(IVEC,1)
829      END IF
830
831* first-order operators:
832      DO IDXA = 1, ORDER
833        IOP(IDXA)=IROPER(LAB(IVEC,IDXA),ISYM(IDXA))
834      END DO
835
836* relaxation flags:
837      IF (TYPE(1:1).EQ.'O'  .OR. TYPE(1:1).EQ.'X' .OR.
838     &    TYPE(1:2).EQ.'EO' .OR. TYPE(1:2).EQ.'EX'     ) THEN
839        NRELAX = 0
840        DO IDXA = 1, ORDER
841         IF ( LORX(IVEC,IDXA) ) THEN
842           IRELAX(IDXA) = IR1KAPPA(LAB(IVEC,IDXA),
843     &                            FREQ(IVEC,IDXA),ISYM(IDXA))
844           NRELAX = NRELAX + 1
845         ELSE
846           IRELAX(IDXA) = 0
847         END IF
848        END DO
849      ELSE
850        NRELAX = 0
851        DO IDXA = 1, ORDER
852         IRELAX(IDXA) = 0
853        END DO
854      END IF
855
856      IF (NRELAX.GT.1) THEN
857         CALL QUIT('NRELAX TOO LARGE IN CC_RHS_SETUP.')
858      END IF
859
860* second-order operators that contribute:
861      IF (     (TYPE(1:1).EQ.'O'  .AND. ORDER.GE.2)
862     &    .OR. (TYPE(1:2).EQ.'EO' .AND. ORDER.GE.2) ) THEN
863        IDXAB  = 0
864        DO IDXB = 2, ORDER
865        DO IDXA = 1, IDXB-1
866         IDXAB = IDXAB + 1
867         IF (IRELAX(IDXA).GT.1 .OR. LPDBSOP(IOP(IDXA)) .OR.
868     &       IRELAX(IDXB).GT.1 .OR. LPDBSOP(IOP(IDXB))      ) THEN
869           INUM        = IROPER2(LAB(IVEC,IDXA),LAB(IVEC,IDXB),
870     &                           LABSOP,ISGNSOP,ISYSOP)
871           IOP2(IDXAB) = IROPER(LABSOP,ISYSOP)
872         ELSE
873           IOP2(IDXAB) = -1
874         END IF
875        END DO
876        END DO
877      END IF
878
879* first-order vectors that contribute:
880      IF (     (TYPE(1:1).EQ.'O'  .AND. ORDER.GT.1)
881     &    .OR. (TYPE(1:2).EQ.'EO' .AND. ORDER.GE.1) ) THEN
882        DO IDXA = 1, ORDER
883         IR1(IDXA)=IR1TAMP(LAB(IVEC,IDXA),LORX(IVEC,IDXA),
884     &                     FREQ(IVEC,IDXA),ISYM(IDXA))
885        END DO
886      END IF
887      IF (TYPE(1:3).EQ.'O1 ') THEN
888        DO IDXA = 1, ORDER
889         IET1(IDXA) = IETA1(LAB(IVEC,IDXA),LORX(IVEC,IDXA),
890     &                      FREQ(IVEC,IDXA),ISYM(IDXA))
891        END DO
892      END IF
893      IF (TYPE(1:2).EQ.'EO' .AND. ORDER.GT.1) THEN
894        call quit('Sonia: please define LPROJ in IER1AMP call')
895        DO IDXA = 1, ORDER
896         IE1(IDXA)=IER1AMP(ISTAT(IVEC,1),EIGV(IVEC),ISYMS,
897     &                     LAB(IVEC,IDXA),FREQ(IVEC,IDXA),ISYM(IDXA))
898        END DO
899      END IF
900      IF (TYPE(1:2).EQ.'CO' .AND. ORDER.GT.1) THEN
901        DO IDXA = 1, ORDER
902         IR1(IDXA)=ILRCAMP(LAB(IVEC,IDXA),ICAU(IVEC,IDXA),ISYM(IDXA))
903        END DO
904      END IF
905
906* second-order vectors that contribute:
907      IF (     (TYPE(1:1).EQ.'O'  .AND. ORDER.GT.2)
908     &    .OR. (TYPE(1:2).EQ.'EO' .AND. ORDER.GE.2) ) THEN
909        IDXAB  = 0
910        DO IDXB = 2, ORDER
911        DO IDXA = 1, IDXB-1
912         IDXAB = IDXAB + 1
913         IR2(IDXAB) =
914     &        IR2TAMP(LAB(IVEC,IDXA),LORX(IVEC,IDXA),
915     &                   FREQ(IVEC,IDXA),ISYM(IDXA),
916     &                LAB(IVEC,IDXB),LORX(IVEC,IDXB),
917     &                   FREQ(IVEC,IDXB),ISYM(IDXB))
918        END DO
919        END DO
920      END IF
921      IF (TYPE(1:2).EQ.'O2' .AND. IOP2(AB).GT.0) THEN
922         IET2(IDXAB) = ICHI2(LAB(IVEC,1),LORX(IVEC,1),
923     &                              FREQ(IVEC,1),ISYM(1),
924     &                       LAB(IVEC,2),LORX(IVEC,2),
925     &                              FREQ(IVEC,2),ISYM(2))
926      END IF
927
928* third-order vectors that contribute:
929      IF (ORDER .GT. 3) THEN
930       IDXABC = 0
931       DO IDXC = 3, ORDER
932       DO IDXB = 2, IDXC-1
933       DO IDXA = 1, IDXB-1
934        IDXABC = IDXABC + 1
935        IR3(IDXABC) =
936     &            IR3TAMP(LAB(IVEC,IDXA),FREQ(IVEC,IDXA),ISYM(IDXA),
937     &                    LAB(IVEC,IDXB),FREQ(IVEC,IDXB),ISYM(IDXB),
938     &                    LAB(IVEC,IDXC),FREQ(IVEC,IDXC),ISYM(IDXC))
939       END DO
940       END DO
941       END DO
942      END IF
943
944*---------------------------------------------------------------------*
945* set up list of D matrix transformations:
946*---------------------------------------------------------------------*
947        IF ( TYPE(1:2).EQ.'O4' ) THEN
948          NDTRAN = NDTRAN + 1
949          IDTRAN(1,NDTRAN) = IR1(A)
950          IDTRAN(2,NDTRAN) = IR1(B)
951          IDTRAN(3,NDTRAN) = IR1(C)
952          IDTRAN(4,NDTRAN) = IR1(D)
953          IDTRAN(5,NDTRAN) = IVEC
954        END IF
955
956*---------------------------------------------------------------------*
957* set up list of C matrix transformations:
958*---------------------------------------------------------------------*
959        IF      ( TYPE(1:2).EQ.'O4' ) THEN
960          DO IDX = 1, NP4AB
961            NCTRAN = NCTRAN + 1
962            ICTRAN(1,NCTRAN) = IR2(IPAB(IDX))
963            ICTRAN(2,NCTRAN) = IR1(IPC(IDX))
964            ICTRAN(3,NCTRAN) = IR1(IPD(IDX))
965            ICTRAN(4,NCTRAN) = IVEC
966          END DO
967        ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
968          NCTRAN = NCTRAN + 1
969          ICTRAN(1,NCTRAN) = IR1(1)
970          ICTRAN(2,NCTRAN) = IR1(2)
971          ICTRAN(3,NCTRAN) = IR1(3)
972          ICTRAN(4,NCTRAN) = IVEC
973        ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
974          NCTRAN = NCTRAN + 1
975          ICTRAN(1,NCTRAN) = IR1(1)
976          ICTRAN(2,NCTRAN) = IR1(2)
977          ICTRAN(3,NCTRAN) = IEX
978          ICTRAN(4,NCTRAN) = IVEC
979        END IF
980
981*---------------------------------------------------------------------*
982* set up list of B matrix transformations
983*---------------------------------------------------------------------*
984        IF      ( TYPE(1:2).EQ.'O4' ) THEN
985          DO IDX = 1, NT4ABC
986            NB1TRAN = NB1TRAN + 1
987            IB1TRAN(1,NB1TRAN) = IR3(ITABC(IDX))
988            IB1TRAN(2,NB1TRAN) = IR1(ITD(IDX))
989            IB1TRAN(3,NB1TRAN) = IVEC
990          END DO
991
992          DO IDX = 1, NP4AB
993            NB2TRAN = NB2TRAN + 1
994            IB2TRAN(1,NB2TRAN) = IR2(IPAB(IDX))
995            IB2TRAN(2,NB2TRAN) = IR2(IPCD(IDX))
996            IB2TRAN(3,NB2TRAN) = IVEC
997          END DO
998        ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
999          DO IDX = 1, NP3AB
1000            NB1TRAN = NB1TRAN + 1
1001            IB1TRAN(1,NB1TRAN) = IR2(IPAB(IDX))
1002            IB1TRAN(2,NB1TRAN) = IR1(IPC(IDX))
1003            IB1TRAN(3,NB1TRAN) = IVEC
1004          END DO
1005        ELSE IF ( TYPE(1:2).EQ.'O2' ) THEN
1006          NB1TRAN = NB1TRAN + 1
1007          IB1TRAN(1,NB1TRAN) = IR1(1)
1008          IB1TRAN(2,NB1TRAN) = IR1(2)
1009          IB1TRAN(3,NB1TRAN) = IVEC
1010        ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
1011          NB1TRAN = NB1TRAN + 1
1012          IB1TRAN(1,NB1TRAN) = IR2(1)
1013          IB1TRAN(2,NB1TRAN) = IEX
1014          IB1TRAN(3,NB1TRAN) = IVEC
1015
1016          DO IDX = 1, NS2A
1017            NB2TRAN = NB2TRAN + 1
1018            IB2TRAN(1,NB2TRAN) = IR1(ISA(IDX))
1019            IB2TRAN(2,NB2TRAN) = IE1(ISB(IDX))
1020            IB2TRAN(3,NB2TRAN) = IVEC
1021          END DO
1022        ELSE IF ( TYPE(1:3).EQ.'EO1' ) THEN
1023          NB1TRAN = NB1TRAN + 1
1024          IB1TRAN(1,NB1TRAN) = IR1(1)
1025          IB1TRAN(2,NB1TRAN) = IEX
1026          IB1TRAN(3,NB1TRAN) = IVEC
1027        ELSE IF ( TYPE(1:3).EQ.'CO2' ) THEN
1028          NB1TRAN = NB1TRAN + 1
1029          IB1TRAN(1,NB1TRAN) = IR1(1)
1030          IB1TRAN(2,NB1TRAN) = IR1(2)
1031          IB1TRAN(3,NB1TRAN) = IVEC
1032        END IF
1033
1034*---------------------------------------------------------------------*
1035* set up list of C{O} matrix transformations:
1036*---------------------------------------------------------------------*
1037        IF      ( TYPE(1:2).EQ.'O4' ) THEN
1038          DO IDX = 1, NS4A
1039            IF (IRELAX(ISA(IDX)).GT.0) THEN
1040              NCATRAN = NCATRAN + 1
1041              ICATRAN(1,NCATRAN) = IOP(ISA(IDX))
1042              ICATRAN(2,NCATRAN) = IR1(ISB(IDX))
1043              ICATRAN(3,NCATRAN) = IR1(ISC(IDX))
1044              ICATRAN(4,NCATRAN) = IR1(ISD(IDX))
1045              ICATRAN(5,NCATRAN) = IVEC
1046              ICATRAN(6,NCATRAN) = IRELAX(ISA(IDX))
1047              ICATRAN(7,NCATRAN) = 0
1048              ICATRAN(8,NCATRAN) = 0
1049              ICATRAN(9,NCATRAN) = 0
1050            END IF
1051          END DO
1052        END IF
1053
1054*---------------------------------------------------------------------*
1055* set up list of B{O} matrix transformations:
1056*---------------------------------------------------------------------*
1057        IF      ( TYPE(1:2).EQ.'O4' ) THEN
1058          DO IDX = 1, NP4AB
1059            NBA1TRAN = NBA1TRAN + 1
1060            IBA1TRAN(1,NBA1TRAN) = IOP(IPC(IDX))
1061            IBA1TRAN(2,NBA1TRAN) = IR2(IPAB(IDX))
1062            IBA1TRAN(3,NBA1TRAN) = IR1(IPD(IDX))
1063            IBA1TRAN(4,NBA1TRAN) = IVEC
1064            IBA1TRAN(5,NBA1TRAN) = IRELAX(IPC(IDX))
1065            IBA1TRAN(6,NBA1TRAN) = 0
1066            IBA1TRAN(7,NBA1TRAN) = 0
1067            IBA1TRAN(8,NBA1TRAN) = 0
1068
1069            NBA1TRAN = NBA1TRAN + 1
1070            IBA1TRAN(1,NBA1TRAN) = IOP(IPD(IDX))
1071            IBA1TRAN(2,NBA1TRAN) = IR2(IPAB(IDX))
1072            IBA1TRAN(3,NBA1TRAN) = IR1(IPC(IDX))
1073            IBA1TRAN(4,NBA1TRAN) = IVEC
1074            IBA1TRAN(5,NBA1TRAN) = IRELAX(IPD(IDX))
1075            IBA1TRAN(6,NBA1TRAN) = 0
1076            IBA1TRAN(7,NBA1TRAN) = 0
1077            IBA1TRAN(8,NBA1TRAN) = 0
1078
1079            IF (IOP2(IPAB(IDX)).GT.0) THEN
1080              NBA2TRAN = NBA2TRAN + 1
1081              IBA2TRAN(1,NBA2TRAN) = IOP2(IPAB(IDX))
1082              IBA2TRAN(2,NBA2TRAN) = IR1(IPC(IDX))
1083              IBA2TRAN(3,NBA2TRAN) = IR1(IPD(IDX))
1084              IBA2TRAN(4,NBA2TRAN) = IVEC
1085              IBA2TRAN(5,NBA2TRAN) = IRELAX(IPA(IDX))
1086              IBA2TRAN(6,NBA2TRAN) = IRELAX(IPB(IDX))
1087              IBA2TRAN(7,NBA2TRAN) = 0
1088              IBA2TRAN(8,NBA2TRAN) = 0
1089            END IF
1090          END DO
1091        ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
1092          DO IDX = 1, NS3A
1093            NBA1TRAN = NBA1TRAN + 1
1094            IBA1TRAN(1,NBA1TRAN) = IOP(ISA(IDX))
1095            IBA1TRAN(2,NBA1TRAN) = IR1(ISB(IDX))
1096            IBA1TRAN(3,NBA1TRAN) = IR1(ISC(IDX))
1097            IBA1TRAN(4,NBA1TRAN) = IVEC
1098            IBA1TRAN(5,NBA1TRAN) = IRELAX(ISA(IDX))
1099            IBA1TRAN(6,NBA1TRAN) = 0
1100            IBA1TRAN(7,NBA1TRAN) = 0
1101            IBA1TRAN(8,NBA1TRAN) = 0
1102          END DO
1103        ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
1104          DO IDX = 1, NS2A
1105            NBA1TRAN = NBA1TRAN + 1
1106            IBA1TRAN(1,NBA1TRAN) = IOP(ISA(IDX))
1107            IBA1TRAN(2,NBA1TRAN) = IR1(ISB(IDX))
1108            IBA1TRAN(3,NBA1TRAN) = IEX
1109            IBA1TRAN(4,NBA1TRAN) = IVEC
1110            IBA1TRAN(5,NBA1TRAN) = IRELAX(ISA(IDX))
1111            IBA1TRAN(6,NBA1TRAN) = 0
1112            IBA1TRAN(7,NBA1TRAN) = 0
1113            IBA1TRAN(8,NBA1TRAN) = 0
1114          END DO
1115        END IF
1116
1117*---------------------------------------------------------------------*
1118* set up list of A{O} vector calculations:
1119*---------------------------------------------------------------------*
1120        IF      ( TYPE(1:2).EQ.'O4' ) THEN
1121          DO IDX = 1, NT4ABC
1122            NAA1TRAN = NAA1TRAN + 1
1123            IAA1TRAN(1,NAA1TRAN) = IOP(ITD(IDX))
1124            IAA1TRAN(2,NAA1TRAN) = IR3(ITABC(IDX))
1125            IAA1TRAN(3,NAA1TRAN) = IVEC
1126            IAA1TRAN(4,NAA1TRAN) = IRELAX(ITD(IDX))
1127            IAA1TRAN(5,NAA1TRAN) = 0
1128            IAA1TRAN(6,NAA1TRAN) = 0
1129            IAA1TRAN(7,NAA1TRAN) = 0
1130          END DO
1131          DO IDX = 1, NP4AB
1132            IF (IOP2(IPAB(IDX)).GT.0) THEN
1133              NAA2TRAN = NAA2TRAN + 1
1134              IAA2TRAN(1,NAA2TRAN) = IOP2(IPAB(IDX))
1135              IAA2TRAN(2,NAA2TRAN) = IR2(IPCD(IDX))
1136              IAA2TRAN(3,NAA2TRAN) = IVEC
1137              IAA2TRAN(4,NAA2TRAN) = IRELAX(IPA(IDX))
1138              IAA2TRAN(5,NAA2TRAN) = IRELAX(IPB(IDX))
1139              IAA2TRAN(6,NAA2TRAN) = 0
1140              IAA2TRAN(7,NAA2TRAN) = 0
1141            END IF
1142          END DO
1143        ELSE IF ( TYPE(1:2).EQ.'O3' ) THEN
1144          DO IDX = 1, NP3AB
1145            NAA1TRAN = NAA1TRAN + 1
1146            IAA1TRAN(1,NAA1TRAN) = IOP(IPC(IDX))
1147            IAA1TRAN(2,NAA1TRAN) = IR2(IPAB(IDX))
1148            IAA1TRAN(3,NAA1TRAN) = IVEC
1149            IAA1TRAN(4,NAA1TRAN) = IRELAX(IPC(IDX))
1150            IAA1TRAN(5,NAA1TRAN) = 0
1151            IAA1TRAN(6,NAA1TRAN) = 0
1152            IAA1TRAN(7,NAA1TRAN) = 0
1153
1154            IF (IOP2(IPAB(IDX)).GT.0) THEN
1155              NAA2TRAN = NAA2TRAN + 1
1156              IAA2TRAN(1,NAA2TRAN) = IOP2(IPAB(IDX))
1157              IAA2TRAN(2,NAA2TRAN) = IR1(IPC(IDX))
1158              IAA2TRAN(3,NAA2TRAN) = IVEC
1159              IAA2TRAN(4,NAA2TRAN) = IRELAX(IPA(IDX))
1160              IAA2TRAN(5,NAA2TRAN) = IRELAX(IPB(IDX))
1161              IAA2TRAN(6,NAA2TRAN) = 0
1162              IAA2TRAN(7,NAA2TRAN) = 0
1163            END IF
1164          END DO
1165        ELSE IF ( TYPE(1:2).EQ.'O2' ) THEN
1166          DO IDX = 1, NS2A
1167            NAA1TRAN = NAA1TRAN + 1
1168            IAA1TRAN(1,NAA1TRAN) = IOP(ISB(IDX))
1169            IAA1TRAN(2,NAA1TRAN) = IR1(ISA(IDX))
1170            IAA1TRAN(3,NAA1TRAN) = IVEC
1171            IAA1TRAN(4,NAA1TRAN) = IRELAX(ISB(IDX))
1172            IAA1TRAN(5,NAA1TRAN) = 0
1173            IAA1TRAN(6,NAA1TRAN) = 0
1174            IAA1TRAN(7,NAA1TRAN) = 0
1175          END DO
1176        ELSE IF ( TYPE(1:3).EQ.'EO2' ) THEN
1177          DO IDX = 1, NS2A
1178            NAA1TRAN = NAA1TRAN + 1
1179            IAA1TRAN(1,NAA1TRAN) = IOP(ISB(IDX))
1180            IAA1TRAN(2,NAA1TRAN) = IE1(ISA(IDX))
1181            IAA1TRAN(3,NAA1TRAN) = IVEC
1182            IAA1TRAN(4,NAA1TRAN) = IRELAX(ISB(IDX))
1183            IAA1TRAN(5,NAA1TRAN) = 0
1184            IAA1TRAN(6,NAA1TRAN) = 0
1185            IAA1TRAN(7,NAA1TRAN) = 0
1186          END DO
1187          IF (IOP2(AB).GT.0) THEN
1188            NAA2TRAN = NAA2TRAN + 1
1189            IAA2TRAN(1,NAA2TRAN) = IOP2(AB)
1190            IAA2TRAN(2,NAA2TRAN) = IEX
1191            IAA2TRAN(3,NAA2TRAN) = IVEC
1192            IAA2TRAN(4,NAA2TRAN) = IRELAX(A)
1193            IAA2TRAN(5,NAA2TRAN) = IRELAX(B)
1194            IAA2TRAN(6,NAA2TRAN) = 0
1195            IAA2TRAN(7,NAA2TRAN) = 0
1196          END IF
1197        ELSE IF ( TYPE(1:3).EQ.'EO1' ) THEN
1198          NAA1TRAN = NAA1TRAN + 1
1199          IAA1TRAN(1,NAA1TRAN) = IOP(1)
1200          IAA1TRAN(2,NAA1TRAN) = IEX
1201          IAA1TRAN(3,NAA1TRAN) = IVEC
1202          IAA1TRAN(4,NAA1TRAN) = IRELAX(1)
1203          IAA1TRAN(5,NAA1TRAN) = 0
1204          IAA1TRAN(6,NAA1TRAN) = 0
1205          IAA1TRAN(7,NAA1TRAN) = 0
1206        ELSE IF ( TYPE(1:3).EQ.'CO2' ) THEN
1207          DO IDX = 1, NS2A
1208            IF (ICAU(IVEC,ISB(IDX)).EQ.0) THEN
1209              NAA1TRAN = NAA1TRAN + 1
1210              IAA1TRAN(1,NAA1TRAN) = IOP(ISB(IDX))
1211              IAA1TRAN(2,NAA1TRAN) = IR1(ISA(IDX))
1212              IAA1TRAN(3,NAA1TRAN) = IVEC
1213              IAA1TRAN(4,NAA1TRAN) = IRELAX(ISB(IDX))
1214              IAA1TRAN(5,NAA1TRAN) = 0
1215              IAA1TRAN(6,NAA1TRAN) = 0
1216              IAA1TRAN(7,NAA1TRAN) = 0
1217            END IF
1218          END DO
1219        END IF
1220
1221*---------------------------------------------------------------------*
1222* set up list of Xi{O} vector calculations:
1223* Note, that we set up here a list for the simultaneous calculation
1224* of the first-order xi "O1" and the first-order eta "X1" vectors.
1225* Xi and eta vectors are only precalculated for orbital relaxed
1226* "operators" or for field-dependent basis sets. For simple unrelaxed
1227* one-electron perturbations they are calculated on the fly when needed
1228*---------------------------------------------------------------------*
1229        IF ( TYPE(1:3).EQ.'O1 ') THEN
1230C         IF ( IRELAX(A).EQ.1 .OR. LPDBSOP(IOP(A)) ) THEN
1231            NXETRAN = NXETRAN + 1
1232            IXETRAN(1,NXETRAN) = IOP(A)
1233            IXETRAN(2,NXETRAN) = 0      ! L0 for first-order ETA vec.
1234            IXETRAN(3,NXETRAN) = IVEC
1235            IXETRAN(4,NXETRAN) = IET1(A)
1236            IXETRAN(5,NXETRAN) = IRELAX(A)
1237            IXETRAN(6,NXETRAN) = 0
1238            IXETRAN(7,NXETRAN) = 0
1239            IXETRAN(8,NXETRAN) = 0
1240C         END IF
1241        ELSE IF ( TYPE(1:2).EQ.'O2' ) THEN
1242          IF ( IOP2(AB).GT.0 ) THEN
1243            NXETRAN = NXETRAN + 1
1244            IXETRAN(1,NXETRAN) = IOP2(AB)
1245            IXETRAN(2,NXETRAN) = 0      ! L0 for second-order ETA vec.
1246            IXETRAN(3,NXETRAN) = IVEC
1247            IXETRAN(4,NXETRAN) = IET2(AB)
1248            IXETRAN(5,NXETRAN) = IRELAX(A)
1249            IXETRAN(6,NXETRAN) = IRELAX(B)
1250            IXETRAN(7,NXETRAN) = 0
1251            IXETRAN(8,NXETRAN) = 0
1252          END IF
1253        ELSE IF ( TYPE(1:3).EQ.'CO1') THEN
1254            NXETRAN = NXETRAN + 1
1255            IXETRAN(1,NXETRAN) = IOP(A)
1256            IXETRAN(2,NXETRAN) = 0      ! L0 for first-order ETA vec.
1257            IXETRAN(3,NXETRAN) = IVEC
1258            IXETRAN(4,NXETRAN) = -1
1259            IXETRAN(5,NXETRAN) = IRELAX(A)
1260            IXETRAN(6,NXETRAN) = 0
1261            IXETRAN(7,NXETRAN) = 0
1262            IXETRAN(8,NXETRAN) = 0
1263        END IF
1264
1265*---------------------------------------------------------------------*
1266* end loop over all requested rhs vectors
1267*---------------------------------------------------------------------*
1268      END DO
1269
1270*---------------------------------------------------------------------*
1271* print the lists:
1272*---------------------------------------------------------------------*
1273* general statistics:
1274      WRITE(LUPRI,'(/,/3X,A,I3,I2,3A)') 'For the requested',NVEC,ORDER,
1275     &      'th.-order amplitude rhs vectors "',TYPE,'".'
1276      WRITE(LUPRI,'((8X,A,I3,A))')
1277     &   ' - ',NDTRAN,            ' D matrix transformations ',
1278     &   ' - ',NCTRAN,            ' C matrix transformations ',
1279     &   ' - ',NB1TRAN+NB2TRAN,   ' B matrix transformations ',
1280     &   ' - ',NCATRAN,           ' C{O} matrix transformations ',
1281     &   ' - ',NBA1TRAN+NBA2TRAN, ' B{O} matrix transformations ',
1282     &   ' - ',NAA1TRAN+NAA2TRAN, ' A{O} matrix transformations ',
1283     &   ' - ',NXETRAN,           'Xi{O} vector calculations    '
1284      IF (NEW_RHS) WRITE(LUPRI,'(14X,A)')
1285     &   '(A{O} matrix included in B matrix)'
1286      WRITE(LUPRI,'(3X,A/,/)') 'will be performed.'
1287
1288
1289* D matrix transformations:
1290      IF (LOCDBG) THEN
1291        WRITE (LUPRI,*) 'List of D matrix transformations:'
1292        DO ITRAN = 1, NDTRAN
1293          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1294     &     (IDTRAN(IDX,ITRAN),IDX=1,5)
1295        END DO
1296        WRITE (LUPRI,*)
1297      END IF
1298
1299* C matrix transformations:
1300      IF (LOCDBG) THEN
1301        WRITE (LUPRI,*) 'List of C matrix transformations:'
1302        DO ITRAN = 1, NCTRAN
1303          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1304     &     (ICTRAN(IDX,ITRAN),IDX=1,4)
1305        END DO
1306        WRITE (LUPRI,*)
1307      END IF
1308
1309* B matrix transformations:
1310      IF (LOCDBG) THEN
1311        WRITE (LUPRI,*) 'List of B matrix transformations (type1):'
1312        DO ITRAN = 1, NB1TRAN
1313          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1314     &     (IB1TRAN(IDX,ITRAN),IDX=1,3)
1315        END DO
1316        WRITE (LUPRI,*) 'List of B matrix transformations (type2):'
1317        DO ITRAN = 1, NB2TRAN
1318          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1319     &     (IB2TRAN(IDX,ITRAN),IDX=1,3)
1320        END DO
1321        WRITE (LUPRI,*)
1322      END IF
1323
1324* C{O} matrix transformations:
1325      IF (LOCDBG) THEN
1326        WRITE (LUPRI,*) 'List of C{O} matrix transformations:'
1327        DO ITRAN = 1, NCATRAN
1328          WRITE(LUPRI,'(A,5I5,5X,(12I5,20X))') MSGDBG,
1329     &     (ICATRAN(IDX,ITRAN),IDX=1,5)
1330        END DO
1331        WRITE (LUPRI,*)
1332      END IF
1333
1334* B{O} matrix transformations:
1335      IF (LOCDBG) THEN
1336        WRITE (LUPRI,*) 'List of B{O} matrix transformations (type1):'
1337        DO ITRAN = 1, NBA1TRAN
1338          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1339     &     (IBA1TRAN(IDX,ITRAN),IDX=1,4)
1340        END DO
1341        WRITE (LUPRI,*) 'List of B{O} matrix transformations (type 2):'
1342        DO ITRAN = 1, NBA2TRAN
1343          WRITE(LUPRI,'(A,4I5,5X,(12I5,20X))') MSGDBG,
1344     &     (IBA2TRAN(IDX,ITRAN),IDX=1,4)
1345        END DO
1346        WRITE (LUPRI,*)
1347      END IF
1348
1349* A{O} matrix calculations:
1350      IF (LOCDBG) THEN
1351        WRITE (LUPRI,*) 'List of A{O} matrix transformations (type1):'
1352        DO ITRAN = 1, NAA1TRAN
1353          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1354     &     (IAA1TRAN(IDX,ITRAN),IDX=1,3)
1355        END DO
1356        WRITE (LUPRI,*) 'List of A{O} matrix transformations (type 2):'
1357        DO ITRAN = 1, NAA2TRAN
1358          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1359     &     (IAA2TRAN(IDX,ITRAN),IDX=1,3)
1360        END DO
1361        WRITE (LUPRI,*)
1362        CALL FLSHFO(LUPRI)
1363      END IF
1364
1365* Xi{O} vector calculations:
1366      IF (LOCDBG) THEN
1367        WRITE (LUPRI,*) 'List of Xi{O} vector calculations:'
1368        DO ITRAN = 1, NXETRAN
1369          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
1370     &     (IXETRAN(IDX,ITRAN),IDX=1,2)
1371        END DO
1372        WRITE (LUPRI,*)
1373        CALL FLSHFO(LUPRI)
1374      END IF
1375
1376
1377      RETURN
1378      END
1379
1380*---------------------------------------------------------------------*
1381*              END OF SUBROUTINE CC_RHS_SETUP                         *
1382*---------------------------------------------------------------------*
1383