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
19*=====================================================================*
20c /* deck ccmcd_setup */
21*=====================================================================*
22      SUBROUTINE CCMCD_SETUP(MXTRAN, MXVEC,
23     &                       IGTRAN, IGDOTS, NGTRAN,
24     &                       IFATRAN,IFADOTS,NFATRAN,
25     &                       IEATRAN,IEADOTS,NEATRAN,
26     &                       IEATRA1,IEADOT1,NEATRA1,
27     &                       IEATRA2,IEADOT2,NEATRA2,
28     &                       IFTRAN, IFDOTS, NFTRAN,
29     &                       IFTRA1, IFDOT1, NFTRA1,
30     &                       IXE2TRA,IX2DOTS,IE2DOTS,NXE2TRA,
31     &                       IX2TRAN,IX2DOT1,NX2TRAN,
32     &                       WORK, LWORK  )
33*---------------------------------------------------------------------*
34*
35*    Purpose: set up for MCD section
36*                - list of G matrix transformations
37*                - list of F{Op} matrix transformations
38*                - list of ETA{Op} vector calculations
39*    if (LUSEPL1) also
40*                - list of B transformations
41*
42*    Written by:       Sonia Coriani 1997-98
43*    Restructured by:  Sonia Coriani 26/01-2000
44*    Orbital relaxation for B operator introduced March 2000. Sonia
45*=====================================================================*
46#if defined (IMPLICIT_NONE)
47      IMPLICIT NONE
48#else
49#  include "implicit.h"
50#endif
51#include "priunit.h"
52#include "ccorb.h"
53#include "ccmcdinf.h"
54#include "ccroper.h"
55#include "ccexci.h"
56#include "cclists.h"
57
58* local parameters:
59      CHARACTER*(23) MSGDBG
60      PARAMETER (MSGDBG = '[debug] CCMCD_SETUP> ')
61      LOGICAL LOCDBG
62      PARAMETER (LOCDBG = .FALSE.)
63
64      INTEGER MXVEC, MXTRAN, LWORK
65
66#if defined (SYS_CRAY)
67      REAL ZERO, D05, TWO
68      REAL EIGVAF, WORK(LWORK)
69#else
70      DOUBLE PRECISION ZERO, D05, TWO
71      DOUBLE PRECISION EIGVAF, WORK(LWORK)
72#endif
73      PARAMETER (ZERO = 0.0D00, D05  = 0.5D00, TWO = 2.0D00)
74
75      INTEGER IGTRAN(MXDIM_GTRAN,MXTRAN)
76      INTEGER IGDOTS(MXVEC,MXTRAN)
77
78      INTEGER IFATRAN(MXDIM_FATRAN,MXTRAN)
79      INTEGER IFADOTS(MXVEC,MXTRAN)
80
81      INTEGER IEATRAN(MXDIM_XEVEC,MXTRAN)
82      INTEGER IEADOTS(MXVEC,MXTRAN)
83
84      INTEGER IEATRA1(MXDIM_XEVEC,MXTRAN)
85      INTEGER IEADOT1(MXVEC,MXTRAN)
86
87      INTEGER IEATRA2(MXDIM_XEVEC,MXTRAN)
88      INTEGER IEADOT2(MXVEC,MXTRAN)
89
90      INTEGER IFTRAN(MXDIM_FTRAN,MXTRAN)
91      INTEGER IFDOTS(MXVEC,MXTRAN)
92
93      INTEGER IFTRA1(MXDIM_FTRAN,MXTRAN)
94      INTEGER IFDOT1(MXVEC,MXTRAN)
95
96      INTEGER IXE2TRA(MXDIM_XEVEC,MXTRAN)
97      INTEGER IX2DOTS(MXVEC,MXTRAN), IE2DOTS(MXVEC,MXTRAN)
98      INTEGER IX2TRAN(MXDIM_XEVEC,MXTRAN)
99      INTEGER IX2DOT1(MXVEC,MXTRAN)
100
101      INTEGER NGTRAN, NFATRAN, NEATRAN
102      INTEGER NEATRA1,NEATRA2,NFTRAN,NFTRA1,NXE2TRA,NX2TRAN
103      INTEGER NBTERMS
104
105      INTEGER ISYMA, ISYMB, ISYMC, ISYMAB, ISYMSF, ISYSOP, ISGNSOP
106      INTEGER IFREQ, IOPER
107
108      CHARACTER*8 LABELA,LABELB,LABELC,LABSOP
109      INTEGER IOPERA, IOPERB, IOPERC, IOPSOP
110      INTEGER ITAMPA, ITAMPB, ITAMPC
111      INTEGER IKAPA,IKAPB
112      INTEGER IVEC, ITRAN, I, ISTATE, INUM, K
113      INTEGER ISTATF, IEXCIF
114      INTEGER IM1F,IZETAA,IPL1A,IZETA0
115
116      LOGICAL LORXA,LORXB,LORXC,LPDBSA,LPDBSB,LPDBSC,LRELAX,LPROJ
117
118* external functions:
119      INTEGER IR1TAMP
120      INTEGER IR1KAPPA
121      INTEGER IPL1ZETA
122      INTEGER ILRMAMP
123      INTEGER IROPER
124      INTEGER IRHSR2
125
126*---------------------------------------------------------------------*
127* initializations:
128*---------------------------------------------------------------------*
129      DO ITRAN = 1, MXTRAN
130        DO I = 1, MXDIM_GTRAN
131          IGTRAN(I,ITRAN) = 0
132        END DO
133        DO I = 1, MXDIM_FATRAN
134          IFATRAN(I,ITRAN) = 0
135        END DO
136        DO I = 1, MXDIM_XEVEC
137          IEATRAN(I,ITRAN) = 0
138          IEATRA1(I,ITRAN) = 0
139          IEATRA2(I,ITRAN) = 0
140          IXE2TRA(I,ITRAN) = 0
141          IX2TRAN(I,ITRAN) = 0
142        END DO
143        IEATRAN(3,ITRAN) = -1
144        IEATRA1(3,ITRAN) = -1
145        IEATRA2(3,ITRAN) = -1
146        IXE2TRA(3,ITRAN) = -1
147        IX2TRAN(3,ITRAN) = -1
148        IX2TRAN(4,ITRAN) = -1
149
150        DO I = 1, MXDIM_FTRAN
151          IFTRAN(I,ITRAN) = 0
152          IFTRA1(I,ITRAN) = 0
153        END DO
154
155        DO IVEC  = 1, MXVEC
156          IGDOTS(IVEC, ITRAN) = 0
157          IFADOTS(IVEC,ITRAN) = 0
158          IEADOTS(IVEC,ITRAN) = 0
159          IEADOT1(IVEC,ITRAN) = 0
160          IEADOT2(IVEC,ITRAN) = 0
161          IFDOTS(IVEC,ITRAN)  = 0
162          IFDOT1(IVEC,ITRAN)  = 0
163          IX2DOTS(IVEC,ITRAN)  = 0
164          IE2DOTS(IVEC,ITRAN)  = 0
165          IX2DOT1(IVEC,ITRAN)  = 0
166        END DO
167      END DO
168
169      NGTRAN  = 0
170      NFATRAN = 0
171      NEATRAN = 0
172      NEATRA1 = 0
173      NEATRA2 = 0
174      NFTRAN  = 0
175      NFTRA1  = 0
176      NXE2TRA = 0
177      NX2TRAN = 0
178
179      NBTERMS = 0   !nr operator triplets matching selected exc. state.
180
181      LPROJ = .FALSE.
182
183*---------------------------------------------------------------------*
184* loop # triples A,B,C
185*     loop # states
186*---------------------------------------------------------------------*
187
188      DO IOPER = 1, NMCDOPER
189
190         IOPERA = IAMCDOP(IOPER)
191         IOPERB = IBMCDOP(IOPER)
192         IOPERC = ICMCDOP(IOPER)
193
194         ISYMA  = ISYOPR(IOPERA)
195         ISYMB  = ISYOPR(IOPERB)
196         ISYMC  = ISYOPR(IOPERC)
197
198         ISYMAB = MULD2H(ISYMA,ISYMB)
199
200         IF (ISYMAB.EQ.ISYMC) THEN
201
202           LABELA = LBLOPR(IOPERA)
203           LABELB = LBLOPR(IOPERB)
204           LABELC = LBLOPR(IOPERC)
205
206           LORXA  = LAMCDRX(IOPER)
207           LORXB  = LBMCDRX(IOPER)
208           LORXC  = LCMCDRX(IOPER)
209
210           LPDBSA = LPDBSOP(IOPERA)
211           LPDBSB = LPDBSOP(IOPERB)
212           LPDBSC = LPDBSOP(IOPERC)
213
214           LRELAX = (LORXA.OR.LPDBSA.OR.
215     &               LORXB.OR.LPDBSB.OR.LORXC.OR.LPDBSC)
216
217           DO ISTATE = 1, NMCDST
218              ISYMSF = IMCDSTSY(ISTATE)
219              ISTATF = IMCDSTNR(ISTATE)
220              IEXCIF = ISYOFE(ISYMSF) + ISTATF
221              EIGVAF = EIGVAL(IEXCIF)
222
223              IF (ISYMSF.EQ.ISYMC) THEN
224
225                 NBTERMS = NBTERMS + 1
226c
227                 ITAMPA = IR1TAMP(LABELA,LORXA,-EIGVAF,ISYMA)
228                 ITAMPB = IR1TAMP(LABELB,LORXB,ZERO,ISYMB)
229                 IKAPA  = 0
230                 IKAPB  = 0
231                 IF (LORXA) IKAPA = IR1KAPPA(LABELA,-EIGVAF,ISYMA)
232                 IF (LORXB) IKAPB = IR1KAPPA(LABELB,ZERO,ISYMB)
233
234                 IF (LOCDBG) THEN
235                    WRITE (LUPRI,*) 'CCMCD_SETUP> LABELA: ', LABELA,
236     &                      '  ISYMA: ',ISYMA,
237     &                      '  LORXA: ', LORXA, ' LPDBSA: ', LPDBSA
238                    WRITE (LUPRI,*) 'CCMCD_SETUP> LABELB: ', LABELB,
239     &                      '  ISYMB: ',ISYMB,
240     &                      '  LORXB: ', LORXB, ' LPDBSB: ', LPDBSB
241                    WRITE (LUPRI,*) 'LUSEPL1: ', LUSEPL1
242                    CALL FLSHFO(LUPRI)
243                 ENDIF
244
245                 IZETA0 = 0
246                 IM1F   = ILRMAMP(IEXCIF,EIGVAF,ISYMC)
247
248                 IF (LUSEPL1) THEN
249                    IF (ISYMB.EQ.1) LPROJ = .TRUE.
250                    IPL1A  = IPL1ZETA(LABELA,LORXA,-EIGVAF,ISYMA,
251     &                                LPROJ,IEXCIF,EIGVAF,ISYMC)
252                 END IF
253*---------------------------------------------------------------------*
254* set up (asymmetric) list of G matrix transformations:(G*Ta*Tb)*Ef
255*---------------------------------------------------------------------*
256                  CALL CC_SETG112(IGTRAN,IGDOTS,MXTRAN,MXVEC,IZETA0,
257     &                            ITAMPA,ITAMPB,IEXCIF,ITRAN,IVEC)
258                  NGTRAN = MAX(NGTRAN,ITRAN)
259*---------------------------------------------------------------------*
260* set up (asymmetric) list of F{O} matrix transformations
261*---------------------------------------------------------------------*
262* (F{A}*TB)*Ef
263                  CALL CC_SETFB12(IFATRAN,IFADOTS,MXTRAN,MXVEC,IZETA0,
264     &                           IOPERA,IKAPA,ITAMPB,IEXCIF,ITRAN,IVEC)
265                  NFATRAN = MAX(NFATRAN,ITRAN)
266* (F{B}*TA)*Ef
267                  CALL CC_SETFB12(IFATRAN,IFADOTS,MXTRAN,MXVEC,IZETA0,
268     &                           IOPERB,IKAPB,ITAMPA,IEXCIF,ITRAN,IVEC)
269                  NFATRAN = MAX(NFATRAN,ITRAN)
270
271*---------------------------------------------------------------------*
272* set up list of ETA{O} vector calculations
273*---------------------------------------------------------------------*
274                  CALL CC_SETXE('Eta',IEATRAN,IEADOTS,
275     &                           MXTRAN,MXVEC,
276     &                           IEXCIF,IOPERA,IKAPA,0,0,0,ITAMPB,
277     &                           ITRAN,IVEC)
278                  NEATRAN = MAX(NEATRAN,ITRAN)
279
280                  IF (LUSEPL1) THEN
281                     ! two permutations
282                     CALL CC_SETXE('Eta',IEATRA1,IEADOT1,
283     &                              MXTRAN,MXVEC,
284     &                              IM1F,IOPERA,IKAPA,0,0,0,ITAMPB,
285     &                              ITRAN,IVEC)
286                     NEATRA1 = MAX(NEATRA1,ITRAN)
287                     CALL CC_SETXE('Eta',IEATRA1,IEADOT1,
288     &                              MXTRAN,MXVEC,
289     &                              IM1F,IOPERB,IKAPB,0,0,0,ITAMPA,
290     &                              ITRAN,IVEC)
291                     NEATRA1 = MAX(NEATRA1,ITRAN)
292                     !
293                     CALL CC_SETXE('Eta',IEATRA2,IEADOT2,
294     &                             MXTRAN,MXVEC,
295     &                             IPL1A,IOPERB,IKAPB,0,0,0,IEXCIF,
296     &                             ITRAN,IVEC)
297                     NEATRA2 = MAX(NEATRA2,ITRAN)
298                  END IF
299*---------------------------------------------------------------------*
300* set up (asymmetric) list of generalized F matrix transformations
301*---------------------------------------------------------------------*
302                  IF (LUSEPL1) THEN
303
304* M^f*(B*T^A*T^B) = (M^f*B*T^A)*T^B = (F[M^f]*T^A)*T^B
305
306                     CALL CCQR_SETF(IFTRAN,IFDOTS,MXTRAN,MXVEC,
307     &                              IM1F,ITAMPA,ITAMPB,ITRAN,IVEC)
308                     NFTRAN = MAX(NFTRAN,ITRAN)
309
310* Z^A*(B*T^B*E^f) = (Z^A*B*T^B)*E^f = (F[Z^A]*T^B)*E^f
311
312                     CALL CC_SETF12(IFTRA1,IFDOT1,MXTRAN,MXVEC,
313     &                              IPL1A,ITAMPB,IEXCIF,ITRAN,IVEC)
314                     NFTRA1 = MAX(NFTRA1,ITRAN)
315
316                  END IF
317*---------------------------------------------------------------------*
318* set up list of Xksi{O2} and Eta{O2} vector calculations
319* Please note: we only use RELAXED B.
320*---------------------------------------------------------------------*
321                  IF (LUSEPL1.AND.(LPDBSB.OR.LORXB)) THEN
322                     CALL CC_FIND_SO_OP(LBLOPR(IOPERA),LBLOPR(IOPERB),
323     &                                  LABSOP,ISYSOP,ISGNSOP,INUM,
324     &                                  WORK,LWORK)
325                     IOPSOP = IROPER(LABSOP,ISYSOP)
326                     !Xksi{O2} dot-multiplied by M^f
327                     CALL CC_SETXE('Xi ',IXE2TRA,IX2DOTS,MXTRAN,MXVEC,
328     &                            IZETA0,IOPSOP,IKAPA,IKAPB,0,0,IM1F,
329     &                             ITRAN,IVEC)
330                     NXE2TRA = MAX(NXE2TRA,ITRAN)
331                     !Eta{O2} dot-multiplied by E^f
332                     CALL CC_SETXE('Eta',IXE2TRA,IE2DOTS,MXTRAN,MXVEC,
333     &                            IZETA0,IOPSOP,IKAPA,IKAPB,0,0,IEXCIF,
334     &                             ITRAN,IVEC)
335                     NXE2TRA = MAX(NXE2TRA,ITRAN)
336                     !Xksi{O2} dot-multiplied by Ebar^f
337                     CALL CC_SETXE('Xi ',IX2TRAN,IX2DOT1,MXTRAN,MXVEC,
338     &                            IZETA0,IOPSOP,IKAPA,IKAPB,0,0,IEXCIF,
339     &                             ITRAN,IVEC)
340                     NX2TRAN = MAX(NX2TRAN,ITRAN)
341                  ELSE IF (LPDBSB.OR.LORXB) THEN
342                    CALL QUIT(
343     &                 'PDBS only programmed for B and LUSEPL1 = true')
344                  END IF
345*---------------------------------------------------------------------*
346* end loops
347*---------------------------------------------------------------------*
348               END IF                ! <F| = C
349            END DO                   ! do # selected MCD states
350         END IF                      ! endif matching sym AB = C
351      END DO                         !do loop on # triples ABC
352
353*---------------------------------------------------------------------*
354* print the lists:
355*---------------------------------------------------------------------*
356* general statistics:
357      WRITE(LUPRI,'(/,/3X,A,I3,A)') 'For the requested ',NBTERMS,
358     &      ' contributions to B terms of magn.circ.dichr. '
359      WRITE(LUPRI,'((8X,A,I3,A))')
360     &   ' - ',NGTRAN,  ' G matrix transformations ',
361     &   ' - ',NFATRAN, ' F{O} matrix transformations ',
362     &   ' - ',NEATRAN, ' ETA{O} vector calculations ',
363     &   ' - ',NEATRA1, ' extra ETA{O} vector calculations (1)',
364     &   ' - ',NEATRA2, ' extra ETA{O} vector calculations (2)',
365     &   ' - ',NFTRAN,  ' extra B matrix transformation (1)',
366     &   ' - ',NFTRA1,  ' extra B matrix transformation (2)(1)'
367      WRITE(LUPRI,'(3X,A,/,/)') 'will be performed.'
368
369
370* G matrix transformations:
371      IF (LOCDBG) THEN
372        WRITE (LUPRI,*) 'List of G matrix transformations:'
373        DO ITRAN = 1, NGTRAN
374          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
375     &     (IGTRAN(I,ITRAN),I=1,3),(IGDOTS(I,ITRAN),I=1,MXVEC)
376        END DO
377        WRITE (LUPRI,*)
378      END IF
379
380* F{O} matrix transformations:
381      IF (LOCDBG) THEN
382        WRITE (LUPRI,*) 'List of F{O} matrix transformations:'
383        DO ITRAN = 1, NFATRAN
384          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
385     &     (IFATRAN(I,ITRAN),I=1,3),(IFADOTS(I,ITRAN),I=1,MXVEC)
386        END DO
387        WRITE (LUPRI,*)
388      END IF
389
390* ETA{O} vector calculations:
391      IF (LOCDBG) THEN
392        WRITE (LUPRI,*) 'List of ETA{O} vector calculations:'
393        DO ITRAN = 1, NEATRAN
394          WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG,
395     &     (IEATRAN(I,ITRAN),I=1,5),(IEADOTS(I,ITRAN),I=1,MXVEC)
396        END DO
397        WRITE (LUPRI,*)
398        IF (LUSEPL1) THEN
399* extra ETA{O} vector calculations:
400          WRITE (LUPRI,*)
401     &         'List of additional ETA{O} vector calculations:'
402          DO ITRAN = 1, NEATRA1
403            WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG,
404     &     (IEATRA1(I,ITRAN),I=1,5),(IEADOT1(I,ITRAN),I=1,MXVEC)
405          END DO
406          WRITE (LUPRI,*)
407          DO ITRAN = 1, NEATRA2
408            WRITE(LUPRI,'(A,2I5,5X,(12I5,20X))') MSGDBG,
409     &     (IEATRA2(I,ITRAN),I=1,5),(IEADOT2(I,ITRAN),I=1,MXVEC)
410          END DO
411          WRITE (LUPRI,*)
412        END IF
413      END IF
414
415* Additional B (generalized F) matrix transformations:
416      IF (LOCDBG.AND.LUSEPL1) THEN
417        WRITE (LUPRI,*) 'List of B/F matrix transformations:'
418        DO ITRAN = 1, NFTRAN
419          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
420     &     (IFTRAN(I,ITRAN),I=1,3),(IFDOTS(I,ITRAN),I=1,MXVEC)
421        END DO
422        WRITE (LUPRI,*)
423        DO ITRAN = 1, NFTRA1
424          WRITE(LUPRI,'(A,3I5,5X,(12I5,20X))') MSGDBG,
425     &     (IFTRA1(I,ITRAN),I=1,3),(IFDOT1(I,ITRAN),I=1,MXVEC)
426        END DO
427        WRITE (LUPRI,*)
428      END IF
429
430      RETURN
431      END
432
433*---------------------------------------------------------------------*
434*              END OF SUBROUTINE CCMCD_SETUP                          *
435*---------------------------------------------------------------------*
436