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 CC_GATHEROO */
21*=====================================================================*
22      SUBROUTINE CC_GATHEROO(XMO,XOO,ISYINT)
23*---------------------------------------------------------------------*
24*     Purpose: gather occupied/occupied integrals from one-electron
25*              MO integral matrix and store them according to IMATIJ
26*
27*     Christof Haettig, January 1997
28*=====================================================================*
29#if defined (IMPLICIT_NONE)
30      IMPLICIT NONE
31#else
32# include "implicit.h"
33#endif
34#include "priunit.h"
35#include "ccsdsym.h"
36#include "ccorb.h"
37
38      LOGICAL LOCDBG
39      PARAMETER (LOCDBG = .FALSE.)
40
41      INTEGER ISYINT, KOFF1, KOFF2, ISYMJ, ISYMI
42
43#if defined (SYS_CRAY)
44      REAL XMO(*)
45      REAL XOO(*)      ! dimension (NMATIJ(ISYINT))
46#else
47      DOUBLE PRECISION XMO(*)
48      DOUBLE PRECISION XOO(*)      ! dimension (NMATIJ(ISYINT))
49#endif
50
51      DO ISYMJ = 1, NSYM
52        ISYMI = MULD2H(ISYMJ,ISYINT)
53        DO J = 1, NRHF(ISYMJ)
54          KOFF1 = IFCRHF(ISYMI,ISYMJ) + NORB(ISYMI)*(J-1) + 1
55          KOFF2 = IMATIJ(ISYMI,ISYMJ) + NRHF(ISYMI)*(J-1) + 1
56          CALL DCOPY(NRHF(ISYMI),XMO(KOFF1),1,XOO(KOFF2),1)
57        END DO
58      END DO
59
60      RETURN
61      END
62*---------------------------------------------------------------------*
63c/* Deck CC_GATHEROV */
64*=====================================================================*
65      SUBROUTINE CC_GATHEROV(XMO,XOV,ISYINT)
66*---------------------------------------------------------------------*
67*     Purpose: gather occupied/virtual integrals from one-electron
68*              MO integral matrix and store them according to IT1AM
69*
70*     Christof Haettig, January 1997
71*=====================================================================*
72#if defined (IMPLICIT_NONE)
73      IMPLICIT NONE
74#else
75# include "implicit.h"
76#endif
77#include "priunit.h"
78#include "ccsdsym.h"
79#include "ccorb.h"
80
81      LOGICAL LOCDBG
82      PARAMETER (LOCDBG = .FALSE.)
83
84      INTEGER ISYINT, KOFF1, KOFF2, ISYMJ, ISYMB
85
86#if defined (SYS_CRAY)
87      REAL XMO(*)
88      REAL XOO(*)     ! dimension (NT1AM(ISYINT))
89#else
90      DOUBLE PRECISION XMO(*)
91      DOUBLE PRECISION XOV(*)     ! dimension (NT1AM(ISYINT))
92#endif
93
94      DO ISYMJ = 1, NSYM
95        ISYMB = MULD2H(ISYMJ,ISYINT)
96        DO J = 1, NRHF(ISYMJ)
97          KOFF1 = IFCVIR(ISYMJ,ISYMB) + J
98          KOFF2 = IT1AM(ISYMB,ISYMJ)  + NVIR(ISYMB)*(J-1) + 1
99          CALL DCOPY(NVIR(ISYMB),XMO(KOFF1),NORB(ISYMJ),XOV(KOFF2),1)
100        END DO
101      END DO
102
103      RETURN
104      END
105*---------------------------------------------------------------------*
106c/* Deck CC_GATHERVV */
107*=====================================================================*
108      SUBROUTINE CC_GATHERVV(XMO,XVV,ISYINT)
109*---------------------------------------------------------------------*
110*     Purpose: gather occupied/occupied integrals from one-electron
111*              MO integral matrix and store them according to IMATAB
112*
113*     Christof Haettig, January 1997
114*=====================================================================*
115#if defined (IMPLICIT_NONE)
116      IMPLICIT NONE
117#else
118# include "implicit.h"
119#endif
120#include "priunit.h"
121#include "ccsdsym.h"
122#include "ccorb.h"
123
124      LOGICAL LOCDBG
125      PARAMETER (LOCDBG = .FALSE.)
126
127      INTEGER ISYINT, KOFF1, KOFF2, ISYMB, ISYMC
128
129#if defined (SYS_CRAY)
130      REAL XMO(*)
131      REAL XOO(*)  ! dimension (NMATAB(ISYINT))
132#else
133      DOUBLE PRECISION XMO(*)
134      DOUBLE PRECISION XVV(*)  ! dimension (NMATAB(ISYINT))
135#endif
136
137      DO ISYMC = 1, NSYM
138        ISYMB = MULD2H(ISYMC,ISYINT)
139        DO C = 1, NVIR(ISYMC)
140          KOFF1 = IFCVIR(ISYMB,ISYMC) + NORB(ISYMB)*(C-1)+NRHF(ISYMB)+1
141          KOFF2 = IMATAB(ISYMB,ISYMC) + NVIR(ISYMB)*(C-1)+1
142          CALL DCOPY(NVIR(ISYMB),XMO(KOFF1),1,XVV(KOFF2),1)
143        END DO
144      END DO
145
146
147      RETURN
148      END
149*---------------------------------------------------------------------*
150c/* Deck CC_XBAR */
151*=====================================================================*
152      SUBROUTINE CC_XBAR(XBAR,XLIAJB,ISYOVOV,T2AMP,ISYAMP,WORK,LWORK)
153*---------------------------------------------------------------------*
154*     Purpose: calculate XBAR intermediate needed for EMAT2
155*              the XBAR intermediate in similar to the X intermediate
156*              in the left transformation, but the two indeces
157*              transposed and the Zeta vector exchanged by L(ia|jb).
158*
159*     symmetries:   ISYOVOV -- XLIAJB
160*                   ISYAMP  -- T2AMP
161*
162*     Christof Haettig, January 1997
163*=====================================================================*
164#if defined (IMPLICIT_NONE)
165      IMPLICIT NONE
166#else
167# include "implicit.h"
168#endif
169#include "priunit.h"
170#include "ccsdsym.h"
171#include "ccorb.h"
172
173      LOGICAL LOCDBG
174      PARAMETER (LOCDBG = .FALSE.)
175
176      INTEGER ISYMX, ISYMI, ISYMJ, MAXJ, NIJ, NJI
177      INTEGER ISYOVOV, ISYAMP, LWORK
178
179#if defined (SYS_CRAY)
180      REAL XBAR(*)      ! dimension (NMATIJ(MULD2H(ISYOVOV,ISYAMP)))
181      REAL XLIAJB(*)    ! dimension (NT2SQ(ISYOVOV))
182      REAL T2AMP(*)     ! dimension (NT2AM(ISYAMP))
183      REAL WORK(LWORK)
184      REAL SWAP
185#else
186      DOUBLE PRECISION XBAR(*)   ! dim. (NMATIJ(MULD2H(ISYOVOV,ISYAMP)))
187      DOUBLE PRECISION XLIAJB(*) ! dim. (NT2SQ(ISYOVOV))
188      DOUBLE PRECISION T2AMP(*)  ! dim. (NT2AM(ISYAMP))
189      DOUBLE PRECISION WORK(LWORK)
190      DOUBLE PRECISION SWAP
191#endif
192
193
194* call CC_XI for the expansive part:
195      CALL CC_XI(XBAR,XLIAJB,ISYOVOV,T2AMP,ISYAMP,WORK,LWORK)
196
197* transpose the result:
198      ISYMX = MULD2H(ISYOVOV,ISYAMP)
199      DO ISYMI = 1, NSYM
200        ISYMJ = MULD2H(ISYMI,ISYMX)
201        IF (ISYMJ .LE. ISYMI) THEN
202          DO I = 1, NRHF(ISYMI)
203            MAXJ =  NRHF(ISYMJ)
204            IF (ISYMJ .EQ. ISYMI) MAXJ = I-1
205          DO J = 1, MAXJ
206            NIJ = IMATIJ(ISYMI,ISYMJ) + NRHF(ISYMI)*(J-1) + I
207            NJI = IMATIJ(ISYMJ,ISYMI) + NRHF(ISYMJ)*(I-1) + J
208            SWAP      = XBAR(NIJ)
209            XBAR(NIJ) = XBAR(NJI)
210            XBAR(NJI) = SWAP
211          END DO
212          END DO
213        END IF
214      END DO
215
216
217      RETURN
218      END
219*---------------------------------------------------------------------*
220c/* Deck CC_YBAR */
221*=====================================================================*
222      SUBROUTINE CC_YBAR(YBAR,XLIAJB,ISYOVOV,T2AMP,ISYAMP,WORK,LWORK)
223*---------------------------------------------------------------------*
224*     Purpose: calculate YBAR intermediate needed for EMAT1
225*              the YBAR intermediate in similar to the X intermediate
226*              in the left transformation, but has the opposite sign
227*              and the Zeta is vector exchanged by L(ia|jb).
228*
229*     symmetries:   ISYOVOV -- XLIAJB
230*                   ISYAMP  -- T2AMP
231*
232*     Christof Haettig, January 1997
233*=====================================================================*
234#if defined (IMPLICIT_NONE)
235      IMPLICIT NONE
236#else
237# include "implicit.h"
238#endif
239#include "priunit.h"
240#include "ccsdsym.h"
241#include "ccorb.h"
242
243      LOGICAL LOCDBG
244      PARAMETER (LOCDBG = .FALSE.)
245
246      INTEGER ISYMY, LEN
247      INTEGER ISYOVOV, ISYAMP, LWORK
248
249#if defined (SYS_CRAY)
250      REAL YBAR(*)       ! dimension (NMATAB(MULD2H(ISYOVOV,ISYAMP)))
251      REAL XLIAJB(*)     ! dimension (NT2SQ(ISYOVOV))
252      REAL T2AMP(*)      ! dimension (NT2AM(ISYAMP))
253      REAL WORK(LWORK)
254#else
255      DOUBLE PRECISION YBAR(*)   ! dim. (NMATAB(MULD2H(ISYOVOV,ISYAMP)))
256      DOUBLE PRECISION XLIAJB(*) ! dim. (NT2SQ(ISYOVOV))
257      DOUBLE PRECISION T2AMP(*)  ! dim. (NT2AM(ISYAMP))
258      DOUBLE PRECISION WORK(LWORK)
259#endif
260
261
262* call CC_YI for the expansive part:
263      CALL CC_YI(YBAR,XLIAJB,ISYOVOV,T2AMP,ISYAMP,WORK,LWORK)
264
265* invert the sign of the result:
266      ISYMY = MULD2H(ISYOVOV,ISYAMP)
267      LEN   = NMATAB(ISYMY)
268      CALL DSCAL(LEN,-1.0d0,YBAR,1)
269
270
271      RETURN
272      END
273*---------------------------------------------------------------------*
274c /* deck iccset1 */
275*=====================================================================*
276       INTEGER FUNCTION ICCSET1(IARRAY,LIST,IDLST,NARR,MAXARR,LAPPEND)
277*---------------------------------------------------------------------*
278*
279*    Purpose: set up and maintain a nonredundant array of response
280*             vectors characterized by LIST and IDLST
281*
282*             IARRAY - array of vectors
283*             LIST   - vector type, given as a list name 'R0','R1',...
284*             IDLST  - index of the vector in LIST
285*
286*             NARR   - used length of IARRAY
287*             MAXARR - maximum dimension of IARRAY
288*
289*    Written by Christof Haettig, january 1997.
290*    PL1 added march 2000, Sonia
291*    QL (Lanczos) added 2010, Sonia
292*=====================================================================*
293      IMPLICIT NONE
294
295#include "priunit.h"
296#include "cciccset.h"
297
298      LOGICAL LOCDBG
299      PARAMETER (LOCDBG = .FALSE.)
300
301      LOGICAL LAPPEND
302      CHARACTER*(*) LIST
303      INTEGER NARR,MAXARR, IDLST, ILIST, I, INDEX, LEN
304      INTEGER IARRAY(2,MAXARR)
305
306      LOGICAL LTST, LEND
307      INTEGER IT,IV,IDX
308
309* statement  functions:
310      LTST(IDX,IT,IV) = IARRAY(1,IDX).EQ.IV .AND. IARRAY(2,IDX).EQ.IT
311      LEND(IDX) = IDX.GT.NARR
312
313*---------------------------------------------------------------------*
314* if LOCDBG echo input:
315*---------------------------------------------------------------------*
316      IF (LOCDBG) THEN
317        WRITE (LUPRI,*) 'LIST, IDLST:',LIST, IDLST
318        WRITE (LUPRI,*) 'NARR,MAXARR:',NARR,MAXARR
319        WRITE (LUPRI,*) 'LAPPEND:',LAPPEND
320        WRITE (LUPRI,*) 'IARRAY:'
321        WRITE (LUPRI,'((/5X,2I5))')
322     &        ((IARRAY(I,INDEX),I=1,2),INDEX=1,NARR)
323      END IF
324
325*---------------------------------------------------------------------*
326* maintain list of nonredundant vectors:
327*---------------------------------------------------------------------*
328
329* convert LIST to an integer:
330      IF (      LIST(1:1).EQ.'R' .OR. LIST(1:1).EQ.'O'
331     &     .OR. LIST(1:1).EQ.'L' .OR. LIST(1:1).EQ.'X'
332     &     .OR. LIST(1:1).EQ.'N' .OR. LIST(1:1).EQ.'M'
333     &     .OR. LIST(1:1).EQ.'Q'                       !Lanczos
334     &     .OR. LIST(1:2).EQ.'D0'                      ) THEN
335         LEN = 2
336      ELSE IF  (LIST(1:1).EQ.'E' .OR. LIST(1:1).EQ.'C'
337     &                           .OR. LIST(1:1).EQ.'P') THEN
338         LEN = 3
339      END IF
340      ILIST = 0
341      DO I = 1, MAXTAB
342       IF (VTABLE(I)(1:LEN).EQ.LIST(1:LEN)) ILIST = I
343      END DO
344      IF (ILIST .EQ.0) CALL QUIT('Unknown list in ICCSET1.')
345
346      INDEX = 1
347      DO WHILE ( .NOT. (LTST(INDEX,ILIST,IDLST) .OR. LEND(INDEX)) )
348        INDEX = INDEX + 1
349      END DO
350
351      IF (INDEX.GT.MAXARR) THEN
352        CALL QUIT('ERROR> list overflow in CCCR_SET1.')
353      END IF
354
355      IF (.NOT.LTST(INDEX,ILIST,IDLST) .OR. LEND(INDEX)) THEN ! append:
356        IF (LAPPEND) THEN
357          IARRAY(1,INDEX) = IDLST
358          IARRAY(2,INDEX) = ILIST
359          NARR = NARR + 1
360        ELSE
361          CALL QUIT('ICCSET1 failed: requested element was '//
362     &          'not registered.')
363        END IF
364      END IF
365
366      ICCSET1 = INDEX
367
368      RETURN
369      END
370
371*---------------------------------------------------------------------*
372*              END OF SUBROUTINE ICCSET1                              *
373*---------------------------------------------------------------------*
374c /* deck iccset2 */
375*=====================================================================*
376      INTEGER FUNCTION ICCSET2(IARRAY,LISTA,IDLSTA,LISTB,IDLSTB,
377     &                         NARR,MAXARR,LAPPEND)
378*---------------------------------------------------------------------*
379*
380*    Purpose: set up and maintain a nonredundant array of pairs of
381*             response vectors characterized by LIST and IDLST
382*
383*            IARRAY   - array of vectors
384*            LISTA/B  - vector type, given as a list name 'R0','R1',...
385*            IDLSTA/B - index of the vector in LIST
386*
387*            MAXARR - maximum dimension of IARRAY
388*
389*    Written by Christof Haettig, january 1997.
390*    PL1 vectors, Sonia 2000
391*    Lanczos QL added in august 2010, Sonia
392*=====================================================================*
393      IMPLICIT NONE
394#include "priunit.h"
395
396#include "cciccset.h"
397
398      LOGICAL LAPPEND
399      CHARACTER*(*) LISTA, LISTB
400      INTEGER NARR,MAXARR, IDLSTA, IDLSTB, I
401      INTEGER IARRAY(4,MAXARR)
402
403      LOGICAL LTST, LEND
404      INTEGER ITA,ITB,IVA,IVB,IDX,ILISTA,ILISTB, INDEX, LENA, LENB
405
406* statement  functions:
407      LTST(IDX,ITA,IVA,ITB,IVB) =
408     &  ( IARRAY(1,IDX).EQ.IVA .AND. IARRAY(2,IDX).EQ.ITA .AND.
409     &    IARRAY(3,IDX).EQ.IVB .AND. IARRAY(4,IDX).EQ.ITB      ) .OR.
410     &  ( IARRAY(3,IDX).EQ.IVA .AND. IARRAY(4,IDX).EQ.ITA .AND.
411     &    IARRAY(1,IDX).EQ.IVB .AND. IARRAY(2,IDX).EQ.ITB      )
412
413      LEND(IDX) = IDX.GT.NARR
414
415*---------------------------------------------------------------------*
416* maintain list of nonredundant vectors:
417*---------------------------------------------------------------------*
418
419* convert LISTA & LISTB to an integer:
420      IF (      LISTA(1:1).EQ.'R' .OR. LISTA(1:1).EQ.'O'
421     &     .OR. LISTA(1:1).EQ.'L' .OR. LISTA(1:1).EQ.'X'
422     &     .OR. LISTA(1:1).EQ.'N' .OR. LISTA(1:1).EQ.'M'
423     &     .OR. LISTA(1:1).EQ.'Q'
424     &     .OR. LISTA(1:2).EQ.'D0'                       ) THEN
425         LENA = 2
426      ELSE IF  (LISTA(1:1).EQ.'E' .OR. LISTA(1:1).EQ.'C'
427     &     .OR. LISTA(1:1).EQ.'P'                        ) THEN
428         LENA = 3
429      END IF
430      IF (      LISTB(1:1).EQ.'R' .OR. LISTB(1:1).EQ.'O'
431     &     .OR. LISTB(1:1).EQ.'L' .OR. LISTB(1:1).EQ.'X'
432     &     .OR. LISTB(1:1).EQ.'N' .OR. LISTB(1:1).EQ.'M'
433     &     .OR. LISTB(1:1).EQ.'Q'
434     &     .OR. LISTB(1:2).EQ.'D0'                       ) THEN
435         LENB = 2
436      ELSE IF  (LISTB(1:1).EQ.'E' .OR. LISTB(1:1).EQ.'C'
437     &     .OR. LISTB(1:1).EQ.'P'                        ) THEN
438         LENB = 3
439      END IF
440      ILISTA = 0
441      ILISTB = 0
442      DO I = 1, MAXTAB
443       IF (VTABLE(I)(1:LENA).EQ.LISTA(1:LENA)) ILISTA = I
444       IF (VTABLE(I)(1:LENB).EQ.LISTB(1:LENB)) ILISTB = I
445      END DO
446      IF (ILISTA .EQ.0) CALL QUIT('Unknown list in ICCSET2.')
447      IF (ILISTB .EQ.0) CALL QUIT('Unknown list in ICCSET2.')
448
449
450      INDEX = 1
451      DO WHILE ( .NOT. (LTST(INDEX,ILISTA,IDLSTA,ILISTB,IDLSTB)
452     &                  .OR. LEND(INDEX)) )
453        INDEX = INDEX + 1
454      END DO
455
456      IF (INDEX.GT.MAXARR) THEN
457        CALL QUIT('ERROR> list overflow in CCCR_SET2.')
458      END IF
459
460      IF (.NOT.LTST(INDEX,ILISTA,IDLSTA,ILISTB,IDLSTB)
461     &    .OR.LEND(INDEX)) THEN
462        IF (LAPPEND) THEN
463          IARRAY(1,INDEX) = IDLSTA
464          IARRAY(2,INDEX) = ILISTA
465          IARRAY(3,INDEX) = IDLSTB
466          IARRAY(4,INDEX) = ILISTB
467          NARR = NARR + 1
468        ELSE
469          WRITE (LUPRI,*) 'error in ICCSET2:'
470          WRITE (LUPRI,*) 'requested pair was not registered.'
471          WRITE (LUPRI,*) IDLSTA, ILISTA, IDLSTB, ILISTB
472          WRITE (LUPRI,*) 'LIST:'
473          WRITE(*,'((/5X,4I5))') ((IARRAY(I,IDX),I=1,4),IDX=1,NARR)
474          CALL QUIT(
475     *         'ICCSET2 failed: requested pair was not registered.')
476        END IF
477      END IF
478
479      ICCSET2 = INDEX
480
481      RETURN
482      END
483
484*---------------------------------------------------------------------*
485*              END OF SUBROUTINE ICCSET2                              *
486*---------------------------------------------------------------------*
487