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