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*=====================================================================* 20 SUBROUTINE CCSDT_O32_NODDY(LISTX,IDLSTX,LISTY,IDLSTY, 21 * LISTXY,IDLSTXY, 22 * XLAMDP0,XLAMDH0,FOCK0, 23 * WORK,LWORK) 24*---------------------------------------------------------------------* 25* Purpose: compute triples contributions to the rhs vectors for * 26* the second-order amplitude response equations (O2) * 27* * 28* B t^x t^y + A{x} t^y + A{y} t^x * 29* * 30* Note: preliminary version, only the contributions of * 31* A{x} t^y + A{y} t^x to the triples rhs vector * 32* are implemented, contributions from B t^x t^y terms, * 33* correction to singles/doubles terms are not implemented * 34* partitioning step is not done... * 35* * 36* Written by Christof Haettig, Mai 2003, Friedrichstal * 37*---------------------------------------------------------------------* 38 IMPLICIT NONE 39#include "priunit.h" 40#include "ccorb.h" 41#include "ccsdsym.h" 42#include "dummy.h" 43#include "ccr1rsp.h" 44 45 INTEGER ISYM0 46 PARAMETER (ISYM0 = 1) 47 48 LOGICAL LOCDBG, HAVET31 49 PARAMETER (LOCDBG=.FALSE., HAVET31=.FALSE.) 50 51 52 CHARACTER*3 LISTXY, LISTX, LISTY 53 INTEGER LUDELD,LUCKJD,LUDKBC 54 INTEGER LWORK, IDLSTX, IDLSTY, IDLSTXY 55 56 CHARACTER FNDELD*6, FNCKJD*6, FNDKBC*4 57 PARAMETER (FNDELD = 'CKDELD' , FNCKJD = 'CKJDEL') 58 PARAMETER (FNDKBC = 'DKBC') 59 60#if defined (SYS_CRAY) 61 REAL XLAMDP0(*), XLAMDH0(*), FOCK0(*) 62 REAL WORK(LWORK) 63 REAL HALF, ONE, ZERO, TWO, DDOT, XNORM 64#else 65 DOUBLE PRECISION XLAMDP0(*), XLAMDH0(*), FOCK0(*) 66 DOUBLE PRECISION WORK(LWORK) 67 DOUBLE PRECISION HALF, ONE, TWO, ZERO, DDOT, XNORM 68#endif 69 PARAMETER( ZERO=0.0D0, HALF=0.5D0, ONE=1.0D0, TWO = 2.0D0 ) 70 71 CHARACTER MODEL*10, LABELX*8, LABELY*8 72 LOGICAL HAVEX, HAVEY, HAVEXY 73 INTEGER KEND1, LWRK1, ISYMX, ISYMY, KXPROP, KYPROP, KXYMAT, 74 & KTAMPX3, KTAMPY3, KWYINT, KTHETAY, KWXINT, KTHETAX, 75 & KTAMP03, NCK, NBJ, NAI, KAIBJCK, KCKAIBJ, KBJCKAI, 76 & KXPROP_AO, KYPROP_AO, KFCKBUF, KEND2, KLAMPX, KLAMPY, 77 & KLAMHX, KLAMHY, KTAMX1, KTAMY1, IOPT, LWRK2, IRREP, IERR, 78 & KT3AMPXY, KTAMP02, KTAMPX2, KTAMPY2, KTHEOCC, KWXYV, 79 & KWXYO, KTHVIRT 80 81 INTEGER ILSTSYM 82 83*---------------------------------------------------------------------* 84* Begin: 85*---------------------------------------------------------------------* 86 IF (NSYM.GT.1) CALL QUIT('No symmetry in CCSDT_O32_NODDY!') 87 88 ISYMX = ILSTSYM(LISTX,IDLSTX) 89 ISYMY = ILSTSYM(LISTY,IDLSTY) 90 91 KEND1 = 1 92 93 KXPROP = KEND1 94 KYPROP = KXPROP + NBAST*NBAST 95 KXYMAT = KYPROP + NBAST*NBAST 96 KEND1 = KXYMAT + NBAST*NBAST 97 98 IF (HAVET31) THEN 99 KTAMPX3 = KEND1 100 KTAMPY3 = KTAMPX3 + NT1AMX*NT1AMX*NT1AMX 101 KEND1 = KTAMPY3 + NT1AMX*NT1AMX*NT1AMX 102 END IF 103 104 KWYINT = KEND1 105 KTHETAY = KWYINT + NT1AMX*NT1AMX*NT1AMX 106 KWXINT = KTHETAY + NT1AMX*NT1AMX*NT1AMX 107 KTHETAX = KWXINT + NT1AMX*NT1AMX*NT1AMX 108 KTAMP03 = KTHETAX + NT1AMX*NT1AMX*NT1AMX 109 KEND1 = KTAMP03 + NT1AMX*NT1AMX*NT1AMX 110 LWRK1 = LWORK - KEND1 111 112 IF (LWRK1.LT.0)CALL QUIT('Out of memory in CCSDT_O32_NODDY (2)') 113 114*---------------------------------------------------------------------* 115* Open some files needed: 116*---------------------------------------------------------------------* 117 LUDELD = -1 118 CALL WOPEN2(LUDELD,FNDELD,64,0) 119 120 LUCKJD = -1 121 CALL WOPEN2(LUCKJD,FNCKJD,64,0) 122 123 LUDKBC = -1 124 CALL WOPEN2(LUDKBC,FNDKBC,64,0) 125 126*---------------------------------------------------------------------* 127* Precalculate a number of intermediates needed for the 128* contributions from A{x} t^y + A{y} t^x, i.e. the w and theta 129* intemediates and some one-electron matrices: 130*---------------------------------------------------------------------* 131 KLAMPX = KEND1 132 KLAMHX = KLAMPX + NLAMDT 133 KLAMPY = KLAMHX + NLAMDT 134 KLAMHY = KLAMPY + NLAMDT 135 KEND2 = KLAMHY + NLAMDT 136 LWRK2 = LWORK - KEND2 137C arnfinn: remove compiler warning 138 IF (.NOT.HAVET31)THEN 139 KTAMPX3 = KEND2 140 KTAMPY3 = KEND2 141 ENDIF 142 IF (LWRK2.LT.0)CALL QUIT('Out of memory in CCSDT_O32_NODDY (2a)') 143 144 CALL CCSDT_T32_AAPREP_NODDY( 145 * LISTX,IDLSTX,ISYMX,WORK(KWXINT),WORK(KTHETAX), 146 * WORK(KXPROP),WORK(KLAMPX),WORK(KLAMHX),WORK(KTAMPX3), 147 * LISTY,IDLSTY,ISYMY,WORK(KWYINT),WORK(KTHETAY), 148 * WORK(KYPROP),WORK(KLAMPY),WORK(KLAMHY),WORK(KTAMPY3), 149 * WORK(KTAMP03),WORK(KXYMAT),XLAMDP0,XLAMDH0,FOCK0, 150 * LUDELD,FNDELD,LUDKBC,FNDKBC,LUCKJD,FNCKJD, 151 * LOCDBG,HAVET31,WORK(KEND2),LWRK2) 152 153*---------------------------------------------------------------------* 154* Compute contributions: 155*---------------------------------------------------------------------* 156 HAVEX = (LISTX(1:3).EQ.'R1 ') 157 HAVEY = (LISTY(1:3).EQ.'R1 ') 158 HAVEXY = HAVEX .OR. HAVEY 159 160 KT3AMPXY = KEND1 161 KEND1 = KT3AMPXY + NT1AMX*NT1AMX*NT1AMX 162 163 KTAMP02 = KEND1 164 KTAMPX2 = KTAMP02 + NT1AMX*NT1AMX 165 KTAMPY2 = KTAMPX2 + NT1AMX*NT1AMX 166 KTHEOCC = KTAMPY2 + NT1AMX*NT1AMX 167 KWXYO = KTHEOCC + NT1AMX*NRHFT*NRHFT 168 KWXYV = KWXYO + NT1AMX*NRHFT*NRHFT 169 KEND2 = KWXYV + NT1AMX*NRHFT*NRHFT 170 171 LWRK2 = LWORK - KEND2 172 IF (LWRK2.LT.NT2AMX) 173 & CALL QUIT('Out of memory in CCSDT_O32_NODDY (3)') 174 175 IOPT = 2 176 CALL CC_RDRSP('R0',0,ISYMOP,IOPT,MODEL,DUMMY,WORK(KEND2)) 177 CALL CC_T2SQ(WORK(KEND2),WORK(KTAMP02),ISYM0) 178 179 IOPT = 2 180 CALL CC_RDRSP(LISTX,IDLSTX,ISYMX,IOPT,MODEL,DUMMY,WORK(KEND2)) 181 CALL CCLR_DIASCL(WORK(KEND2),TWO,ISYMX) 182 CALL CC_T2SQ(WORK(KEND2),WORK(KTAMPX2),ISYMX) 183 184 IOPT = 2 185 CALL CC_RDRSP(LISTY,IDLSTY,ISYMY,IOPT,MODEL,DUMMY,WORK(KEND2)) 186 CALL CCLR_DIASCL(WORK(KEND2),TWO,ISYMY) 187 CALL CC_T2SQ(WORK(KEND2),WORK(KTAMPY2),ISYMY) 188 189 190 CALL DZERO(WORK(KT3AMPXY),NT1AMX*NT1AMX*NT1AMX) 191 192 CALL CCSDT_O32VIR_NODDY(WORK(KWXINT),WORK(KXPROP), 193 & WORK(KTAMPX2),HAVEX, 194 & WORK(KWYINT),WORK(KYPROP), 195 & WORK(KTAMPY2),HAVEY, 196 & WORK(KTAMP03),WORK(KTAMP02), 197 & WORK(KXYMAT),HAVEXY, 198 & WORK(KTHEOCC),WORK(KWXYV),WORK(KWXYO), 199 & WORK(KT3AMPXY),.TRUE., 200 & DUMMY,.FALSE.,DUMMY,.FALSE. ) 201 202 KTHVIRT = KTAMPY2 + NT1AMX*NT1AMX 203 KWXYV = KTHVIRT + NT1AMX*NVIRT*NVIRT 204 KWXYO = KWXYV + NT1AMX*NVIRT*NVIRT 205 KEND2 = KWXYO + NT1AMX*NVIRT*NVIRT 206 207 LWRK2 = LWORK - KEND2 208 IF (LWRK2.LT.0) 209 & CALL QUIT('Out of memory in CCSDT_O32_NODDY (4)') 210 211 CALL CCSDT_O32OCC_NODDY(WORK(KTHETAX),WORK(KXPROP),HAVEX, 212 & WORK(KTHETAY),WORK(KYPROP),HAVEY, 213 & WORK(KT3AMPXY),WORK(KTHVIRT), 214 & WORK(KWXYV),WORK(KWXYO),.TRUE.,.TRUE., 215 & DUMMY,.FALSE.,DUMMY,.FALSE.) 216 217 IF (LOCDBG) THEN 218 WRITE(LUPRI,*) 'CCSDT_O32_NODDY> (incomplete) 2.ord. rhs O^XY:' 219C CALL DSCAL(NT1AMX*NT1AMX*NT1AMX,HALF,WORK(KT3AMPXY),1) 220 CALL CCSDT_CLEAN_T3(WORK(KT3AMPXY),NT1AMX,NVIRT,NRHFT) 221 CALL PRINT_PT3_NODDY(WORK(KT3AMPXY)) 222 END IF 223 224*---------------------------------------------------------------------* 225* print debug output, close files and return: 226*---------------------------------------------------------------------* 227 CALL WCLOSE2(LUDELD,FNDELD,'KEEP') 228 CALL WCLOSE2(LUCKJD,FNCKJD,'KEEP') 229 CALL WCLOSE2(LUDKBC,FNDKBC,'KEEP') 230 231 RETURN 232 END 233*---------------------------------------------------------------------* 234* END OF SUBROUTINE CCSDT_O32_NODDY * 235*---------------------------------------------------------------------* 236*=====================================================================* 237 SUBROUTINE CCSDT_O32VIR_NODDY(WXINT,XPROP,TAMPX2,HAVEX, 238 & WYINT,YPROP,TAMPY2,HAVEY, 239 & TAMP03,TAMP02,XYMAT,HAVEXY, 240 & THEOCC,WXYV,WXYO, 241 & T3AMPXY,LT3AMPXY, 242 & WXYOINT,LWXYOINT, 243 & WXYVINT,LWXYVINT) 244*---------------------------------------------------------------------* 245* Purpose: set up loop over virtual, sort Wbar, W, theta as they * 246* would be sorted in the real code and compute the * 247* contributions to second-order rhs vectors * 248* * 249* Written by Christof Haettig, May 2003, Friedrichstal * 250*---------------------------------------------------------------------* 251 IMPLICIT NONE 252#include "priunit.h" 253#include "ccorb.h" 254#include "ccsdsym.h" 255 256 LOGICAL HAVEX, HAVEY, HAVEXY, LT3AMPXY, LWXYOINT, LWXYVINT 257 258#if defined (SYS_CRAY) 259 REAL XPROP(NORBT,NORBT), YPROP(NORBT,NORBT) 260 REAL XYMAT(NORBT,NORBT) 261 REAL WYINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 262 REAL WXINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 263 REAL TAMP03(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 264 REAL WXYOINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 265 REAL WXYVINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 266 REAL T3AMPXY(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 267 REAL THEOCC(NVIRT,NRHFT,NRHFT,NRHFT) 268 REAL WXYV(NVIRT,NRHFT,NRHFT,NRHFT) 269 REAL WXYO(NVIRT,NRHFT,NRHFT,NRHFT) 270 REAL TAMP02(NVIRT,NRHFT,NVIRT,NRHFT) 271 REAL TAMPX2(NVIRT,NRHFT,NVIRT,NRHFT) 272 REAL TAMPY2(NVIRT,NRHFT,NVIRT,NRHFT) 273 REAL CONTRIB 274#else 275 DOUBLE PRECISION XPROP(NORBT,NORBT), YPROP(NORBT,NORBT) 276 DOUBLE PRECISION XYMAT(NORBT,NORBT) 277 DOUBLE PRECISION WYINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 278 DOUBLE PRECISION WXINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 279 DOUBLE PRECISION TAMP03(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 280 DOUBLE PRECISION WXYOINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 281 DOUBLE PRECISION WXYVINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 282 DOUBLE PRECISION T3AMPXY(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 283 DOUBLE PRECISION THEOCC(NVIRT,NRHFT,NRHFT,NRHFT) 284 DOUBLE PRECISION WXYV(NVIRT,NRHFT,NRHFT,NRHFT) 285 DOUBLE PRECISION WXYO(NVIRT,NRHFT,NRHFT,NRHFT) 286 DOUBLE PRECISION TAMP02(NVIRT,NRHFT,NVIRT,NRHFT) 287 DOUBLE PRECISION TAMPX2(NVIRT,NRHFT,NVIRT,NRHFT) 288 DOUBLE PRECISION TAMPY2(NVIRT,NRHFT,NVIRT,NRHFT) 289 DOUBLE PRECISION CONTRIB 290#endif 291 292 293 DO D = 1, NVIRT 294 DO B = 1, NVIRT 295 296 CALL DZERO(WXYV,NVIRT*NRHFT*NRHFT*NRHFT) 297 CALL DZERO(WXYO,NVIRT*NRHFT*NRHFT*NRHFT) 298 299 IF (HAVEX) THEN 300c ----------------------------------------------------------- 301c build THEOCC^DB(ai,lm) = w(bm,ai,dl)+w(dl,bm,ai)+w(ai,dl,bm) 302c with the w intermediate for Y 303c (should be available from A density routines !) 304c ----------------------------------------------------------- 305 DO A = 1, NVIRT 306 DO I = 1, NRHFT 307 DO L = 1, NRHFT 308 DO M = 1, NRHFT 309 THEOCC(A,I,L,M) = WYINT(B,M,A,I,D,L)+ 310 & WYINT(D,L,B,M,A,I)+WYINT(A,I,D,L,B,M) 311 END DO 312 END DO 313 END DO 314 END DO 315 END IF 316 317 IF (HAVEX) THEN 318c ----------------------------------------------------------- 319c Do one-index transformation: 320c WXYV(ai,lm) -= THEOCC^DB(ci,lm) X_ac 321c WXYO(ai,lm) += THEOCC^DB(ak,lm) X_ki 322c (should be available by modification of routines for W^X !) 323c ----------------------------------------------------------- 324 DO A = 1, NVIRT 325 DO I = 1, NRHFT 326 DO L = 1, NRHFT 327 DO M = 1, NRHFT 328 DO C = 1, NVIRT 329 WXYV(A,I,L,M) = WXYV(A,I,L,M) - 330 & THEOCC(C,I,L,M) * XPROP(NRHFT+A,NRHFT+C) 331 END DO 332 DO K = 1, NRHFT 333 WXYO(A,I,L,M) = WXYO(A,I,L,M) + 334 & THEOCC(A,K,L,M) * XPROP(K,I) 335 END DO 336 END DO 337 END DO 338 END DO 339 END DO 340 END IF 341 342 IF (HAVEY) THEN 343c ----------------------------------------------------------- 344c build THEOCC^DB(ai,lm) = w(bm,ai,dl)+w(dl,bm,ai)+w(ai,dl,bm) 345c with the w intermediate for X 346c (should be available from A density routines !) 347c ----------------------------------------------------------- 348 DO A = 1, NVIRT 349 DO I = 1, NRHFT 350 DO L = 1, NRHFT 351 DO M = 1, NRHFT 352 THEOCC(A,I,L,M) = WXINT(B,M,A,I,D,L)+ 353 & WXINT(D,L,B,M,A,I)+WXINT(A,I,D,L,B,M) 354 END DO 355 END DO 356 END DO 357 END DO 358 END IF 359 360 IF (HAVEY) THEN 361c ----------------------------------------------------------- 362c Do one-index transformation: 363c WXYV(ai,lm) -= THEOCC^DB(ci,lm) Y_ac 364c WXYO(ai,lm) += THEOCC^DB(ak,lm) Y_ki 365c (should be available by modification of routines for W^X !) 366c ----------------------------------------------------------- 367 DO A = 1, NVIRT 368 DO I = 1, NRHFT 369 DO L = 1, NRHFT 370 DO M = 1, NRHFT 371 DO C = 1, NVIRT 372 WXYV(A,I,L,M) = WXYV(A,I,L,M) - 373 & THEOCC(C,I,L,M) * YPROP(NRHFT+A,NRHFT+C) 374 END DO 375 DO K = 1, NRHFT 376 WXYO(A,I,L,M) = WXYO(A,I,L,M) + 377 & THEOCC(A,K,L,M) * YPROP(K,I) 378 END DO 379 END DO 380 END DO 381 END DO 382 END DO 383 END IF 384 385 IF (HAVEX) THEN 386c ----------------------------------------------------------- 387c Add the comutator [[X,T^y_2],T^0_2] 388c (should be (almost) be available from routines for W^X !) 389c ----------------------------------------------------------- 390 DO A = 1, NVIRT 391 DO I = 1, NRHFT 392 DO L = 1, NRHFT 393 DO M = 1, NRHFT 394 DO C = 1, NVIRT 395 DO K = 1, NRHFT 396 WXYO(A,I,L,M) = WXYO(A,I,L,M) - ( 397 & TAMP02(D,L,A,K) * TAMPY2(B,M,C,I) * XPROP(K,NRHFT+C)+ 398 & TAMPY2(D,L,A,K) * TAMP02(B,M,C,I) * XPROP(K,NRHFT+C)+ 399 & TAMP02(B,M,A,K) * TAMPY2(D,L,C,I) * XPROP(K,NRHFT+C)+ 400 & TAMPY2(B,M,A,K) * TAMP02(D,L,C,I) * XPROP(K,NRHFT+C)) 401 END DO 402 END DO 403 END DO 404 END DO 405 END DO 406 END DO 407 END IF 408 409 IF (HAVEY) THEN 410c ----------------------------------------------------------- 411c Add the comutator [[Y,T^x],T^0_2] 412c (should be (almost) available from routines for W^X !) 413c ----------------------------------------------------------- 414 DO A = 1, NVIRT 415 DO I = 1, NRHFT 416 DO L = 1, NRHFT 417 DO M = 1, NRHFT 418 DO C = 1, NVIRT 419 DO K = 1, NRHFT 420 WXYO(A,I,L,M) = WXYO(A,I,L,M) - ( 421 & TAMP02(D,L,A,K) * TAMPX2(B,M,C,I) * YPROP(K,NRHFT+C)+ 422 & TAMPX2(D,L,A,K) * TAMP02(B,M,C,I) * YPROP(K,NRHFT+C)+ 423 & TAMP02(B,M,A,K) * TAMPX2(D,L,C,I) * YPROP(K,NRHFT+C)+ 424 & TAMPX2(B,M,A,K) * TAMP02(D,L,C,I) * YPROP(K,NRHFT+C)) 425 END DO 426 END DO 427 END DO 428 END DO 429 END DO 430 END DO 431 END IF 432 433 IF (HAVEXY) THEN 434c ----------------------------------------------------------- 435c Add comutator [(X^Y+Y^X),T^0]: 436c (should be available from routines for W^X !) 437c ----------------------------------------------------------- 438 DO A = 1, NVIRT 439 DO I = 1, NRHFT 440 DO L = 1, NRHFT 441 DO M = 1, NRHFT 442 DO C = 1, NVIRT 443 WXYV(A,I,L,M) = WXYV(A,I,L,M) - 444 & TAMP03(C,I,D,L,B,M) * XYMAT(NRHFT+A,NRHFT+C) 445 END DO 446 DO K = 1, NRHFT 447 WXYO(A,I,L,M) = WXYO(A,I,L,M) + 448 & TAMP03(A,K,D,L,B,M) * XYMAT(K,I) 449 END DO 450 END DO 451 END DO 452 END DO 453 END DO 454 END IF 455 456c ----------------------------------------------------------- 457c Assuming that the partioning problem with W intermediates 458c has been solved in the real code, the noddy code uses now 459c a short cut and and applies C^abd_iml while storing in an 460c N^6 array 461c ----------------------------------------------------------- 462 IF (LT3AMPXY) THEN 463 DO A = 1, NVIRT 464 DO I = 1, NRHFT 465 DO L = 1, NRHFT 466 DO M = 1, NRHFT 467 CONTRIB = WXYO(A,I,L,M) + WXYV(A,I,L,M) 468 T3AMPXY(A,I,B,M,D,L) = T3AMPXY(A,I,B,M,D,L) + CONTRIB 469 T3AMPXY(B,M,D,L,A,I) = T3AMPXY(B,M,D,L,A,I) + CONTRIB 470 T3AMPXY(D,L,A,I,B,M) = T3AMPXY(D,L,A,I,B,M) + CONTRIB 471 END DO 472 END DO 473 END DO 474 END DO 475 END IF 476 477 IF (LWXYOINT) THEN 478 DO A = 1, NVIRT 479 DO I = 1, NRHFT 480 DO L = 1, NRHFT 481 DO M = 1, NRHFT 482 WXYOINT(A,I,B,M,D,L) = WXYO(A,I,L,M) 483 END DO 484 END DO 485 END DO 486 END DO 487 END IF 488 489 IF (LWXYVINT) THEN 490 DO A = 1, NVIRT 491 DO I = 1, NRHFT 492 DO L = 1, NRHFT 493 DO M = 1, NRHFT 494 WXYVINT(A,I,B,M,D,L) = WXYV(A,I,L,M) 495 END DO 496 END DO 497 END DO 498 END DO 499 END IF 500 501 END DO 502 END DO 503 504C IF (LWXYOINT) THEN 505C WRITE(LUPRI,*) 'WXYOINT in CCSDT_O32VIR_NODDY:' 506C CALL PRINT_PT3_NODDY(WXYOINT) 507C END IF 508 509C IF (LWXYVINT) THEN 510C WRITE(LUPRI,*) 'WXYVINT in CCSDT_O32VIR_NODDY:' 511C CALL PRINT_PT3_NODDY(WXYVINT) 512C END IF 513 514 RETURN 515 END 516*---------------------------------------------------------------------* 517* END OF SUBROUTINE CCSDT_O32VIR_NODDY * 518*---------------------------------------------------------------------* 519*=====================================================================* 520 SUBROUTINE CCSDT_O32OCC_NODDY(THETAX,XPROP,HAVEX, 521 & THETAY,YPROP,HAVEY, 522 & T3AMPXY,THVIRT,WXYV,WXYO, 523 & LOTRN,LVTRN, 524 & TXYOINT,LTXYOINT, 525 & TXYVINT,LTXYVINT) 526*---------------------------------------------------------------------* 527* Purpose: set up loop over occupied, sort Wbar, W, theta as they * 528* would be sorted in the real code and compute the * 529* contributions to second-order rhs vectors * 530* * 531* Written by Christof Haettig, May 2003, Friedrichstal * 532*---------------------------------------------------------------------* 533 IMPLICIT NONE 534#include "priunit.h" 535#include "ccorb.h" 536#include "ccsdsym.h" 537 538 LOGICAL HAVEX, HAVEY, HAVEXY, LOTRN, LVTRN, LTXYOINT, LTXYVINT 539 540#if defined (SYS_CRAY) 541 REAL XPROP(NORBT,NORBT), YPROP(NORBT,NORBT) 542 REAL THETAY(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 543 REAL THETAX(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 544 REAL TXYOINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 545 REAL TXYVINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 546 REAL T3AMPXY(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 547 REAL THVIRT(NVIRT,NVIRT,NVIRT,NRHFT) 548 REAL WXYO(NVIRT,NVIRT,NVIRT,NRHFT) 549 REAL WXYV(NVIRT,NVIRT,NVIRT,NRHFT) 550 REAL CONTRIB 551#else 552 DOUBLE PRECISION XPROP(NORBT,NORBT), YPROP(NORBT,NORBT) 553 DOUBLE PRECISION THETAY(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 554 DOUBLE PRECISION THETAX(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 555 DOUBLE PRECISION TXYOINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 556 DOUBLE PRECISION TXYVINT(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 557 DOUBLE PRECISION T3AMPXY(NVIRT,NRHFT,NVIRT,NRHFT,NVIRT,NRHFT) 558 DOUBLE PRECISION THVIRT(NVIRT,NVIRT,NVIRT,NRHFT) 559 DOUBLE PRECISION WXYO(NVIRT,NVIRT,NVIRT,NRHFT) 560 DOUBLE PRECISION WXYV(NVIRT,NVIRT,NVIRT,NRHFT) 561 DOUBLE PRECISION CONTRIB 562#endif 563 564 565 DO L = 1, NRHFT 566 DO M = 1, NRHFT 567 568 CALL DZERO(WXYO,NT1AMX*NVIRT*NVIRT) 569 CALL DZERO(WXYV,NT1AMX*NVIRT*NVIRT) 570 571 IF (HAVEX) THEN 572c ----------------------------------------------------------- 573c build THVIRT^DB(ai,lm) = theta(b-bar m,ai,dl)+ 574c theta(d-bar l,bm,ai) + theta(a-bar i,dl,bm) 575c with the theta intermediate for Y 576c (should be available from A density routines !) 577c ----------------------------------------------------------- 578 DO I = 1, NRHFT 579 DO A = 1, NVIRT 580 DO B = 1, NVIRT 581 DO D = 1, NVIRT 582 THVIRT(B,D,A,I) = THETAY(B,M,A,I,D,L)+ 583 & THETAY(D,L,B,M,A,I)+THETAY(A,I,D,L,B,M) 584 END DO 585 END DO 586 END DO 587 END DO 588 END IF 589 590 IF (HAVEX) THEN 591c ----------------------------------------------------------- 592c Do one-index transformation: WXY(bdai) = 593c + THVIRT^LM(bdci) X_ac - THVIRT^LM(bdak) X_ki 594c ----------------------------------------------------------- 595 DO I = 1, NRHFT 596 DO A = 1, NVIRT 597 DO B = 1, NVIRT 598 DO D = 1, NVIRT 599 IF (LVTRN) THEN 600 DO C = 1, NVIRT 601 WXYV(B,D,A,I) = WXYV(B,D,A,I) - 602 & THVIRT(B,D,C,I) * XPROP(NRHFT+A,NRHFT+C) 603 END DO 604 END IF 605 IF (LOTRN) THEN 606 DO K = 1, NRHFT 607 WXYO(B,D,A,I) = WXYO(B,D,A,I) + 608 & THVIRT(B,D,A,K) * XPROP(K,I) 609 END DO 610 END IF 611 END DO 612 END DO 613 END DO 614 END DO 615 END IF 616 617 IF (HAVEY) THEN 618c ----------------------------------------------------------- 619c build THVIRT^DB(ai,lm) = theta(b-bar m,ai,dl)+ 620c theta(d-bar l,bm,ai) + theta(a-bar i,dl,bm) 621c with the theta intermediate for X 622c (should be available from A density routines !) 623c ----------------------------------------------------------- 624 DO I = 1, NRHFT 625 DO A = 1, NVIRT 626 DO B = 1, NVIRT 627 DO D = 1, NVIRT 628 THVIRT(B,D,A,I) = THETAX(B,M,A,I,D,L)+ 629 & THETAX(D,L,B,M,A,I)+THETAX(A,I,D,L,B,M) 630 END DO 631 END DO 632 END DO 633 END DO 634 END IF 635 636 IF (HAVEX) THEN 637c ----------------------------------------------------------- 638c Do one-index transformation: WXY(bdai) = 639c + THVIRT^LM(bdci) X_ac - THVIRT^LM(bdak) X_ki 640c ----------------------------------------------------------- 641 DO I = 1, NRHFT 642 DO A = 1, NVIRT 643 DO B = 1, NVIRT 644 DO D = 1, NVIRT 645 IF (LVTRN) THEN 646 DO C = 1, NVIRT 647 WXYV(B,D,A,I) = WXYV(B,D,A,I) - 648 & THVIRT(B,D,C,I) * YPROP(NRHFT+A,NRHFT+C) 649 END DO 650 END IF 651 IF (LOTRN) THEN 652 DO K = 1, NRHFT 653 WXYO(B,D,A,I) = WXYO(B,D,A,I) + 654 & THVIRT(B,D,A,K) * YPROP(K,I) 655 END DO 656 END IF 657 END DO 658 END DO 659 END DO 660 END DO 661 END IF 662 663c ----------------------------------------------------------- 664c Assuming that the partioning problem with W intermediates 665c has been solved in the real code, the noddy code uses now 666c a short cut and and applies C^abd_iml while storing in an 667c N^6 array 668c ----------------------------------------------------------- 669 DO I = 1, NRHFT 670 DO A = 1, NVIRT 671 DO B = 1, NVIRT 672 DO D = 1, NVIRT 673 CONTRIB = WXYO(B,D,A,I) + WXYV(B,D,A,I) 674 T3AMPXY(A,I,B,M,D,L) = T3AMPXY(A,I,B,M,D,L) + CONTRIB 675 T3AMPXY(B,M,D,L,A,I) = T3AMPXY(B,M,D,L,A,I) + CONTRIB 676 T3AMPXY(D,L,A,I,B,M) = T3AMPXY(D,L,A,I,B,M) + CONTRIB 677 END DO 678 END DO 679 END DO 680 END DO 681 682 IF (LTXYVINT) THEN 683 DO I = 1, NRHFT 684 DO A = 1, NVIRT 685 DO B = 1, NVIRT 686 DO D = 1, NVIRT 687 TXYVINT(A,I,B,M,D,L)=TXYVINT(A,I,B,M,D,L) + WXYV(B,D,A,I) 688 TXYVINT(B,M,D,L,A,I)=TXYVINT(B,M,D,L,A,I) + WXYV(B,D,A,I) 689 TXYVINT(D,L,A,I,B,M)=TXYVINT(D,L,A,I,B,M) + WXYV(B,D,A,I) 690 END DO 691 END DO 692 END DO 693 END DO 694 ENDIF 695 696 IF (LTXYOINT) THEN 697 DO I = 1, NRHFT 698 DO A = 1, NVIRT 699 DO B = 1, NVIRT 700 DO D = 1, NVIRT 701 TXYOINT(A,I,B,M,D,L)=TXYOINT(A,I,B,M,D,L) + WXYO(B,D,A,I) 702 TXYOINT(B,M,D,L,A,I)=TXYOINT(B,M,D,L,A,I) + WXYO(B,D,A,I) 703 TXYOINT(D,L,A,I,B,M)=TXYOINT(D,L,A,I,B,M) + WXYO(B,D,A,I) 704 END DO 705 END DO 706 END DO 707 END DO 708 ENDIF 709 710 END DO 711 END DO 712 713C IF (LTXYVINT) THEN 714C WRITE(LUPRI,*) 'TXYINT in CCSDT_O32OCC_NODDY:' 715C CALL PRINT_PT3_NODDY(TXYVINT) 716C END IF 717 718 RETURN 719 END 720*---------------------------------------------------------------------* 721* END OF SUBROUTINE CCSDT_O32OCC_NODDY * 722*---------------------------------------------------------------------* 723*=====================================================================* 724 SUBROUTINE CCSDT_T32_AAPREP_NODDY( 725 * LISTX,IDLSTX,ISYMX,WXINT,THETAX,XPROP, 726 * XLAMDPX,XLAMDHX,TAMPX3, 727 * LISTY,IDLSTY,ISYMY,WYINT,THETAY,YPROP, 728 * XLAMDPY,XLAMDHY,TAMPY3, 729 * T3AMP0,XYMAT,XLAMDP0,XLAMDH0,FOCK0, 730 * LUDELD,FNDELD,LUDKBC,FNDKBC,LUCKJD,FNCKJD, 731 * LOCDBG,HAVET31,WORK,LWORK) 732*---------------------------------------------------------------------* 733* Purpose: prepare intermediates needed for the A{x} t^y + A{y} t^x * 734* contribution to the triples part of the second-order * 735* amplitude response * 736* * 737* Written by Christof Haettig, Jun 2003, Friedrichstal * 738*---------------------------------------------------------------------* 739 IMPLICIT NONE 740#include "priunit.h" 741#include "ccorb.h" 742#include "ccsdsym.h" 743#include "dummy.h" 744#include "ccr1rsp.h" 745#include "ccnoddy.h" 746 747 INTEGER ISYM0 748 PARAMETER (ISYM0 = 1) 749 750 CHARACTER*(3) LISTX, LISTY 751 CHARACTER*(*) FNDELD, FNDKBC, FNCKJD 752 INTEGER IDLSTX, IDLSTY, LWORK, LUDELD, LUDKBC, LUCKJD, 753 & ISYMX, ISYMY 754 LOGICAL LOCDBG, HAVET31 755 756#if defined (SYS_CRAY) 757 REAL WXINT(*), THETAX(*), XPROP(*), TAMPX3(*) 758 REAL WYINT(*), THETAY(*), YPROP(*), TAMPY3(*) 759 REAL T3AMP0(*), XYMAT(*) 760 REAL XLAMDP0(*), XLAMDH0(*), FOCK0(*), WORK(*) 761 REAL XLAMDPX(*), XLAMDHX(*), XLAMDPY(*), XLAMDHY(*) 762 REAL ONE, THREE, FREQX, FREQY 763#else 764 DOUBLE PRECISION WXINT(*), THETAX(*), XPROP(*), TAMPX3(*) 765 DOUBLE PRECISION WYINT(*), THETAY(*), YPROP(*), TAMPY3(*) 766 DOUBLE PRECISION T3AMP0(*), XYMAT(*) 767 DOUBLE PRECISION XLAMDP0(*), XLAMDH0(*), FOCK0(*), WORK(*) 768 DOUBLE PRECISION XLAMDPX(*), XLAMDHX(*), XLAMDPY(*), XLAMDHY(*) 769 DOUBLE PRECISION ONE, THREE, FREQX, FREQY 770#endif 771 PARAMETER(ONE = 1.0D0, THREE = 3.0D0) 772 773 CHARACTER*10 MODEL 774 CHARACTER*8 LABELX, LABELY 775 INTEGER KXPROP_AO, KTAMX1, KYPROP_AO, KTAMY1, 776 & KFCKBUF, KEND2, LWRK2, IRREP, IERR, IOPT, 777 & NCK, NBJ, NAI, KAIBJCK, KCKAIBJ, KBJCKAI, 778 & KFOCKD, KINT1S0, KINT2S0, KLAMPB, KLAMHB, KFOCKB, 779 & KINT1SB, KINT2SB, LUTEMP, LUSIFC 780 781*---------------------------------------------------------------------* 782* Get w^y and theta^y intermediate into memory: 783*---------------------------------------------------------------------* 784 IF (LOCDBG) THEN 785 write(lupri,*) 'entered CCSDT_T32_AAPREP_NODDY...' 786 write(lupri,*) 'listx,idlstx:',listx,idlstx 787 write(lupri,*) 'listy,idlsty:',listy,idlsty 788 END IF 789 790*---------------------------------------------------------------------* 791* Get w^y and theta^y intermediate into memory: 792*---------------------------------------------------------------------* 793 IF (LISTY(1:3).EQ.'R1 ') THEN 794 795 CALL CCSDT_GWTIC(LISTY,IDLSTY,WYINT,THETAY, 796 & T3AMP0,HAVET31,TAMPY3,LOCDBG, 797 & XLAMDP0,XLAMDH0,FOCK0, 798 & LUDELD,FNDELD,LUDKBC,FNDKBC,LUCKJD,FNCKJD, 799 & WORK,LWORK) 800 801 ELSE IF (LISTY(1:3).EQ.'RE ' .OR. LISTY(1:3).EQ.'RC ') THEN 802 803 KEND2 = 1 804 805 KFOCKD = KEND2 806 KEND2 = KFOCKD + NORBT 807 808 KINT1S0 = KEND2 809 KINT2S0 = KINT1S0 + NT1AMX*NVIRT*NVIRT 810 KEND2 = KINT2S0 + NRHFT*NRHFT*NT1AMX 811 812 KLAMPB = KEND2 813 KLAMHB = KLAMPB + NLAMDT 814 KFOCKB = KLAMHB + NLAMDT 815 KEND2 = KFOCKB + N2BST(1) 816 817 KINT1SB = KEND2 818 KINT2SB = KINT1SB + NT1AMX*NVIRT*NVIRT 819 KEND2 = KINT2SB + NRHFT*NRHFT*NT1AMX 820 821 LWRK2 = LWORK - KEND2 822 IF (LWRK2 .LT. 0) 823 & CALL QUIT('Out of memory in CCSDT_T32_AAPREP_NODDY...') 824 825C --------------------------------------- 826C Read precalculated integrals from file: 827C XINT1S0 = (CK|BD) 828C XINT2S0 = (CK|LJ) 829C --------------------------------------- 830 CALL CCSDT_READ_NODDY(.FALSE.,DUMMY,DUMMY,DUMMY,DUMMY, 831 & .FALSE.,DUMMY,DUMMY, 832 & .TRUE.,WORK(KINT1S0),WORK(KINT2S0), 833 & .FALSE.,DUMMY,DUMMY, 834 & .FALSE.,DUMMY,DUMMY,DUMMY,DUMMY, 835 & NORBT,NLAMDT,NRHFT,NVIRT,NT1AMX) 836 837 LUSIFC = -1 838 CALL GPOPEN(LUSIFC,'SIRIFC','OLD',' ','UNFORMATTED',IDUMMY, 839 & .FALSE.) 840 REWIND LUSIFC 841 CALL MOLLAB('TRCCINT ',LUSIFC,LUPRI) 842 READ (LUSIFC) 843 READ (LUSIFC) (WORK(KFOCKD-1+I), I=1,NORBT) 844 CALL GPCLOSE(LUSIFC,'KEEP') 845 846 847 CALL DZERO(WYINT,NT1AMX*NT1AMX*NT1AMX) 848 849 CALL CCSDT_T31_NODDY(WYINT,LISTY,IDLSTY,FREQY,.FALSE., 850 & .FALSE.,WORK(KINT1S0),WORK(KINT2S0), 851 & .FALSE.,DUMMY,DUMMY, 852 & .FALSE.,DUMMY,DUMMY, 853 & WORK(KINT1SB),WORK(KINT2SB), 854 & WORK(KLAMPB),WORK(KLAMHB),WORK(KFOCKB), 855 & XLAMDP0,XLAMDH0,FOCK0,DUMMY,WORK(KFOCKD), 856 & WORK(KEND2),LWRK2) 857 858 CALL DSCAL(NT1AMX*NT1AMX*NT1AMX,-ONE,WYINT,1) 859 860 IF (HAVET31) THEN 861 CALL DCOPY(NT1AMX*NT1AMX*NT1AMX,WYINT,1,TAMPY3,1) 862 END IF 863 864 ! Theta-Y is set to zero for this case 865 CALL DZERO(THETAY,NT1AMX*NT1AMX*NT1AMX) 866 867 ! account for the fact that in the real code a combination 868 ! of three W intermediates give the complete Tbar3 869 CALL DSCAL(NT1AMX*NT1AMX*NT1AMX,ONE/THREE,WYINT,1) 870 871 LUTEMP = -1 872 CALL GPOPEN(LUTEMP,FILNODT30,'UNKNOWN',' ','UNFORMATTED', 873 & IDUMMY,.FALSE.) 874 READ(LUTEMP) (T3AMP0(I),I=1,NT1AMX*NT1AMX*NT1AMX) 875 CALL GPCLOSE(LUTEMP,'KEEP') 876 877 ELSE 878 CALL QUIT('Unknown or illegal LISTY="'//LISTY(1:3)// 879 & '" in CCSDT_T32_AAPREP_NODDY') 880 END IF 881 882C IF (LOCDBG) THEN 883C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> W^Y:' 884C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> LISTY:',LISTY,IDLSTY 885C CALL PRINT_PT3_NODDY(WYINT) 886C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> Theta^Y:' 887C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> LISTY:',LISTY,IDLSTY 888C CALL PRINT_PT3_NODDY(THETAY) 889C END IF 890 891*---------------------------------------------------------------------* 892* Get w^x and theta^x intermediate into memory: 893*---------------------------------------------------------------------* 894 IF (LISTX(1:3).EQ.'R1 ') THEN 895 896 CALL CCSDT_GWTIC(LISTX,IDLSTX,WXINT,THETAX, 897 & T3AMP0,HAVET31,TAMPX3,LOCDBG, 898 & XLAMDP0,XLAMDH0,FOCK0, 899 & LUDELD,FNDELD,LUDKBC,FNDKBC,LUCKJD,FNCKJD, 900 & WORK,LWORK) 901 902 ELSE IF (LISTX(1:3).EQ.'RE ' .OR. LISTX(1:3).EQ.'RC ') THEN 903 904 KEND2 = 1 905 906 KFOCKD = KEND2 907 KEND2 = KFOCKD + NORBT 908 909 KINT1S0 = KEND2 910 KINT2S0 = KINT1S0 + NT1AMX*NVIRT*NVIRT 911 KEND2 = KINT2S0 + NRHFT*NRHFT*NT1AMX 912 913 KLAMPB = KEND2 914 KLAMHB = KLAMPB + NLAMDT 915 KFOCKB = KLAMHB + NLAMDT 916 KEND2 = KFOCKB + N2BST(1) 917 918 KINT1SB = KEND2 919 KINT2SB = KINT1SB + NT1AMX*NVIRT*NVIRT 920 KEND2 = KINT2SB + NRHFT*NRHFT*NT1AMX 921 922 LWRK2 = LWORK - KEND2 923 IF (LWRK2 .LT. 0) 924 & CALL QUIT('Out of memory in CCSDT_T32_AAPREP_NODDY...') 925 926C --------------------------------------- 927C Read precalculated integrals from file: 928C XINT1S0 = (CK|BD) 929C XINT2S0 = (CK|LJ) 930C --------------------------------------- 931 CALL CCSDT_READ_NODDY(.FALSE.,DUMMY,DUMMY,DUMMY,DUMMY, 932 & .FALSE.,DUMMY,DUMMY, 933 & .TRUE.,WORK(KINT1S0),WORK(KINT2S0), 934 & .FALSE.,DUMMY,DUMMY, 935 & .FALSE.,DUMMY,DUMMY,DUMMY,DUMMY, 936 & NORBT,NLAMDT,NRHFT,NVIRT,NT1AMX) 937 938 LUSIFC = -1 939 CALL GPOPEN(LUSIFC,'SIRIFC','OLD',' ','UNFORMATTED',IDUMMY, 940 & .FALSE.) 941 REWIND LUSIFC 942 CALL MOLLAB('TRCCINT ',LUSIFC,LUPRI) 943 READ (LUSIFC) 944 READ (LUSIFC) (WORK(KFOCKD-1+I), I=1,NORBT) 945 CALL GPCLOSE(LUSIFC,'KEEP') 946 947 948 CALL DZERO(WXINT,NT1AMX*NT1AMX*NT1AMX) 949 950 CALL CCSDT_T31_NODDY(WXINT,LISTX,IDLSTX,FREQX,.FALSE., 951 & .FALSE.,WORK(KINT1S0),WORK(KINT2S0), 952 & .FALSE.,DUMMY,DUMMY, 953 & .FALSE.,DUMMY,DUMMY, 954 & WORK(KINT1SB),WORK(KINT2SB), 955 & WORK(KLAMPB),WORK(KLAMHB),WORK(KFOCKB), 956 & XLAMDP0,XLAMDH0,FOCK0,DUMMY,WORK(KFOCKD), 957 & WORK(KEND2),LWRK2) 958 959 CALL DSCAL(NT1AMX*NT1AMX*NT1AMX,-ONE,WXINT,1) 960 961 IF (HAVET31) THEN 962 CALL DCOPY(NT1AMX*NT1AMX*NT1AMX,WXINT,1,TAMPX3,1) 963 END IF 964 965 ! Theta-Y is set to zero for this case 966 CALL DZERO(THETAX,NT1AMX*NT1AMX*NT1AMX) 967 968 ! account for the fact that in the real code a combination 969 ! of three W intermediates give the complete Tbar3 970 CALL DSCAL(NT1AMX*NT1AMX*NT1AMX,ONE/THREE,WXINT,1) 971 972 LUTEMP = -1 973 CALL GPOPEN(LUTEMP,FILNODT30,'UNKNOWN',' ','UNFORMATTED', 974 & IDUMMY,.FALSE.) 975 READ(LUTEMP) (T3AMP0(I),I=1,NT1AMX*NT1AMX*NT1AMX) 976 CALL GPCLOSE(LUTEMP,'KEEP') 977 978 ELSE 979 CALL QUIT('Unknown or illegal LISTX="'//LISTY(1:3)// 980 & '" in CCSDT_T32_AAPREP_NODDY') 981 END IF 982 983C IF (LOCDBG) THEN 984C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> W^X:' 985C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> LISTX:',LISTX,IDLSTX 986C CALL PRINT_PT3_NODDY(WXINT) 987C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> Theta^X:' 988C WRITE(LUPRI,*)'CCSDT_T32_AAPREP_NODDY> LISTX:',LISTX,IDLSTX 989C CALL PRINT_PT3_NODDY(THETAX) 990C END IF 991 992*---------------------------------------------------------------------* 993* Allocate some intermediates needed for the next sections: 994*---------------------------------------------------------------------* 995 KXPROP_AO = 1 996 KYPROP_AO = KXPROP_AO + NBAST*NBAST 997 KFCKBUF = KYPROP_AO + NBAST*NBAST 998 KEND2 = KFCKBUF + NBAST*NBAST 999 1000 KTAMX1 = KEND2 1001 KTAMY1 = KTAMX1 + NT1AMX 1002 KEND2 = KTAMY1 + NT1AMX 1003 1004 LWRK2 = LWORK - KEND2 1005 IF (LWRK2.LT.0) CALL 1006 & QUIT('Out of memory in CCSDT_T32_AAPREP_NODDY') 1007 1008*---------------------------------------------------------------------* 1009* Get property matrices in Lambda-MO basis: 1010*---------------------------------------------------------------------* 1011 CALL DZERO(WORK(KXPROP_AO),NORBT*NORBT) 1012 CALL DZERO(WORK(KYPROP_AO),NORBT*NORBT) 1013 1014 CALL DZERO(XPROP,NORBT*NORBT) 1015 CALL DZERO(YPROP,NORBT*NORBT) 1016 1017 IF (LISTX(1:3).EQ.'R1 ') THEN 1018 LABELX = LRTLBL(IDLSTX) 1019 1020 ! read property integrals from file: 1021 CALL CCPRPAO(LABELX,.TRUE.,WORK(KXPROP_AO),IRREP,ISYMX,IERR, 1022 & WORK(KEND2),LWRK2) 1023 IF ((IERR.GT.0) .OR. (IERR.EQ.0 .AND. IRREP.NE.ISYMX)) THEN 1024 CALL QUIT('CCSDT_T32_NODDY: error reading operator '//LABELX) 1025 ELSE IF (IERR.LT.0) THEN 1026 CALL DZERO(WORK(KXPROP_AO),N2BST(ISYMX)) 1027 END IF 1028 CALL DCOPY(NORBT*NORBT,WORK(KXPROP_AO),1,XPROP,1) 1029 1030 ! transform property integrals to Lambda-MO basis 1031 CALL CC_FCKMO(XPROP,XLAMDP0,XLAMDH0, 1032 & WORK(KEND2),LWRK2,ISYMX,1,1) 1033 END IF 1034 1035 IF (LISTY(1:3).EQ.'R1 ') THEN 1036 LABELY = LRTLBL(IDLSTY) 1037 ! read property integrals from file: 1038 CALL CCPRPAO(LABELY,.TRUE.,WORK(KYPROP_AO),IRREP,ISYMY,IERR, 1039 & WORK(KEND2),LWRK2) 1040 IF ((IERR.GT.0) .OR. (IERR.EQ.0 .AND. IRREP.NE.ISYMY)) THEN 1041 CALL QUIT('CCSDT_T32_NODDY: error reading operator '//LABELY) 1042 ELSE IF (IERR.LT.0) THEN 1043 CALL DZERO(WORK(KYPROP_AO),N2BST(ISYMY)) 1044 END IF 1045 CALL DCOPY(NORBT*NORBT,WORK(KYPROP_AO),1,YPROP,1) 1046 1047 ! transform property integrals to Lambda-MO basis 1048 CALL CC_FCKMO(YPROP,XLAMDP0,XLAMDH0, 1049 & WORK(KEND2),LWRK2,ISYMY,1,1) 1050 END IF 1051 1052*---------------------------------------------------------------------* 1053* Construct sum of one-index transformed property matrices: 1054* XYMAT = [X,T^Y_1] + [Y,T^X_1] 1055*---------------------------------------------------------------------* 1056 IOPT = 1 1057 CALL CC_RDRSP(LISTX,IDLSTX,ISYMX,IOPT,MODEL,WORK(KTAMX1),DUMMY) 1058 1059 CALL CCLR_LAMTRA(XLAMDP0,XLAMDPX,XLAMDH0,XLAMDHX, 1060 & WORK(KTAMX1),ISYMX) 1061 1062 IOPT = 1 1063 CALL CC_RDRSP(LISTY,IDLSTY,ISYMY,IOPT,MODEL,WORK(KTAMY1),DUMMY) 1064 1065 CALL CCLR_LAMTRA(XLAMDP0,XLAMDPY,XLAMDH0,XLAMDHY, 1066 & WORK(KTAMY1),ISYMY) 1067 1068 1069 CALL DZERO(XYMAT,NORBT*NORBT) 1070 1071 IF (LISTX(1:3).EQ.'R1 ') THEN 1072 ! add [X,T^Y_1] 1073 CALL DCOPY(NORBT*NORBT,WORK(KXPROP_AO),1,WORK(KFCKBUF),1) 1074 CALL CC_FCKMO(WORK(KFCKBUF),XLAMDPY,XLAMDH0, 1075 & WORK(KEND2),LWRK2,ISYMX,ISYMY,ISYM0) 1076 CALL DAXPY(NORBT*NORBT,ONE,WORK(KFCKBUF),1,XYMAT,1) 1077 1078 CALL DCOPY(NORBT*NORBT,WORK(KXPROP_AO),1,WORK(KFCKBUF),1) 1079 CALL CC_FCKMO(WORK(KFCKBUF),XLAMDP0,XLAMDHY, 1080 & WORK(KEND2),LWRK2,ISYMX,ISYM0,ISYMY) 1081 CALL DAXPY(NORBT*NORBT,ONE,WORK(KFCKBUF),1,XYMAT,1) 1082 END IF 1083 1084 IF (LISTY(1:3).EQ.'R1 ') THEN 1085 ! add [Y,T^X_1] 1086 CALL DCOPY(NORBT*NORBT,WORK(KYPROP_AO),1,WORK(KFCKBUF),1) 1087 CALL CC_FCKMO(WORK(KFCKBUF),XLAMDPX,XLAMDH0, 1088 & WORK(KEND2),LWRK2,ISYMY,ISYMX,ISYM0) 1089 CALL DAXPY(NORBT*NORBT,ONE,WORK(KFCKBUF),1,XYMAT,1) 1090 1091 CALL DCOPY(NORBT*NORBT,WORK(KYPROP_AO),1,WORK(KFCKBUF),1) 1092 CALL CC_FCKMO(WORK(KFCKBUF),XLAMDP0,XLAMDHX, 1093 & WORK(KEND2),LWRK2,ISYMY,ISYM0,ISYMX) 1094 CALL DAXPY(NORBT*NORBT,ONE,WORK(KFCKBUF),1,XYMAT,1) 1095 END IF 1096 1097 RETURN 1098 END 1099*=====================================================================* 1100