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 cc3_gmatnew */ 20 SUBROUTINE CC3_GMATNEW(LISTL,IDLSTL,LISTB,IDLSTB, 21 * LISTC,IDLSTC, 22 * OMEGA1,OMEGA2, 23 * ISYRES, 24 * WORK,LWORK, 25 * LUTOC,FNTOC, 26 * LU3VI,FN3VI,LUDKBC3,FNDKBC3, 27 * LU3FOPX,FN3FOPX,LU3FOP2X,FN3FOP2X) 28* 29********************************************************************** 30* 31* Calculate following contributions to gmat: 32* 33* <L3|[H^B,T^{C}_{2}],tau_{1}]|HF> 34* <L3|[H^C,T^{B}_{2}],tau_{1}]|HF> 35* <L3|[H^BC,T^{0}_{2}],tau_{1}]|HF> 36* <L3|[H^BC,\tau_nu_2]|HF> 37* 38* L0 or L1 list may be used for the multipliers; thus gmat or cmat 39* contributions may be calculated, respectively. 40* 41* If we have L0 then WB3X = .false. and SKIPGEI = .true. 42* 43* If we have L1 (or other similar lists---like LE,M1,N2---that require 44* W intermediate) then WB3X = .true. and SKIPGEI = .false. 45* 46* 47* F. Pawlowski and P. Jorgensen, Spring 2003. 48* 49* (based on the cc3_bfmat routine) 50* 51* 52* April-2004, Aarhus, FP: VVVV integrals removed, flags LVVVV and SKIPGEI 53* introduced. 54* 55* 56********************************************************************** 57* 58 IMPLICIT NONE 59C 60#include "priunit.h" 61#include "ccorb.h" 62#include "ccl1rsp.h" 63#include "ccsdsym.h" 64#include "ccr1rsp.h" 65#include "dummy.h" 66#include "iratdef.h" 67#include "ccinftap.h" 68#include "ccexci.h" 69#include "ccsdinp.h" 70C 71 LOGICAL SKIPGEI 72C 73 INTEGER ISYM0 74 PARAMETER(ISYM0 = 1) 75C 76 CHARACTER LISTL0*3, LISTL*3,LISTB*3,LISTC*3,LISTL1R*3,LABELL1*8 77 CHARACTER*(*) FNTOC, FN3VI, FNDKBC3, FN3FOPX, FN3FOP2X 78C 79 CHARACTER*(8) FN3VI2 80 CHARACTER*(5) FN3FOP 81 CHARACTER*(6) FN3FOP2 82C 83 PARAMETER (FN3VI2 = 'CC3_VI12', FN3FOP='PTFOP', FN3FOP2='PTFOP2') 84C 85 CHARACTER*13 FNDELDR,FNDKBCR,FNDKBCR4,FNCKJDR,FN4VB,FN4VC,FN4VBC, 86 * FN3SRTR 87C 88 CHARACTER*6 FNGEI,FNFEI 89 CHARACTER*5 FNN1 90 PARAMETER( FNGEI = 'N1_GEI' , FNFEI = 'N1_FEI' , FNN1 = 'N1MAT' ) 91 INTEGER LUGEI,LUFEI,LUN1 92C 93 CHARACTER*7 FNGEIB,FNFEIB 94 CHARACTER*6 FNN1B 95 PARAMETER(FNGEIB='N1_GEIB' , FNFEIB='N1_FEIB' , FNN1B='N1MATB' ) 96 INTEGER LUGEIB,LUFEIB,LUN1B 97C 98 CHARACTER*7 FNGEIC,FNFEIC 99 CHARACTER*6 FNN1C 100 PARAMETER(FNGEIC='N1_GEIC' , FNFEIC='N1_FEIC' , FNN1C='N1MATC' ) 101 INTEGER LUGEIC,LUFEIC,LUN1C 102C 103 INTEGER IDLSTL0,IDLSTL,IDLSTB,IDLSTC,IDLSTL1R 104 INTEGER LUTOC, LU3VI, LUDKBC3, LU3FOPX, LU3FOP2X 105 INTEGER LU3VI2,LU3FOP,LU3FOP2 106 INTEGER LUDELDR,LUDKBCR,LUDKBCR4,LUCKJDR,LU4VB,LU4VC,LU4VBC 107 INTEGER LU3SRTR 108C 109 CHARACTER CDUMMY*1 110C 111 LOGICAL LOCDBG,LORXL1 112 PARAMETER (LOCDBG = .FALSE.) 113C 114 INTEGER ISYRES,LWORK 115 INTEGER ISYML1,ISYML1R,ISYMRB 116 INTEGER ISINT2,KRBJIA,KLAMP0,KLAMH0,KFOCKD,KFOCK0CK,KT2TP,KL1AM 117 INTEGER KL2TP,KEND1,LWRK1 118 INTEGER KL1L1,KL2L1,KFOCKL1 119 INTEGER KT2RB,KT2RC,KFOCK0 120 INTEGER IOPT 121 INTEGER ISINT1RB,ISINT2RB,ISINT1BC 122 INTEGER KBIOOOO,KBIOVVO,KBIOOVV,KBCIOOOO,KBCIOVVO,KBCIOOVV 123 INTEGER ISINT2L1R,ISYFCKL1R,KXIAJB,KT3BOG1,KT3BOL1,KT3BOG2 124 INTEGER KT3BOL2,KLAMPL1R,KLAMHL1R 125 INTEGER KW3XOGX1,KFOCKL1RCK,KW3BXOGX1,KW3BXOLX1 126 INTEGER KW3BXOG1,KW3BXOL1,KT1L1R 127 INTEGER LENGTH 128 INTEGER ISINT1,ISINT1L1R 129 INTEGER IOFF 130 INTEGER ISYMD,ISYCKBD0,ISYCKBDL1R,ISYCKBDRB,ISYCKBDRBC 131 INTEGER KVVVVB,KVVVVBC 132 INTEGER KT3BVDL1,KT3BVDL2,KT3BVDL3,KEND3,LWRK3 133 INTEGER KT3BVDG1,KT3BVDG2,KT3BVDG3 134 INTEGER KW3BXVDG1,KW3BXVDG2,KW3BXVDL1,KW3BXVDL2,KW3BXVDGX1 135 INTEGER KW3BXVDGX2,KW3BXVDLX1,KW3BXVDLX2 136 INTEGER KTRVIR,KTRVIR1,KEND4,LWRK4 137 INTEGER ISYMB,ISYALJB0,ISYALJD0,ISYALJBL1,ISYALJDL1,ISYMBD 138 INTEGER ISCKIJ,ISWBMAT,ISYCKD 139 INTEGER KSMAT2,KUMAT2,KDIAG,KINDSQ 140 INTEGER KDIAGWB,KINDSQWB 141 INTEGER KINDEX,KINDEX2 142 INTEGER KINDEXBL1,KINDEXDL1,KTMAT,KW3BMAT 143 INTEGER KT3BVBG1,KT3BVBG2,KT3BVBG3,KSMAT4,KUMAT4,KT3BVBL1 144 INTEGER KT3BVBL2,KT3BVBL3,KEND5,LWRK5 145 INTEGER KINTOC,KTROCR,KTROCR1 146 147 INTEGER LENSQ,LENSQWB 148 INTEGER ISTBD,KTB,LENSQTB,KINDSQTB,ISTB 149 150 INTEGER ISYMRC,ISYMBC,ISINT2RC,KCIOOOO,KCIOVVO,KCIOOVV 151 INTEGER ISYCKBDRC,KVVVVC,ISINT1RC 152C 153 INTEGER IR1TAMP 154 INTEGER ILSTSYM 155C 156 INTEGER KEND0,LWRK0 157C 158 INTEGER ISYMN1,ISYMN2,ISYMN1B,ISYMN2B,ISYMN1C,ISYMN2C 159 INTEGER KN2MAT,KN2MATB,KN2MATC,KINDSQN,KINDSQNB,KINDSQNC 160 INTEGER LENSQN,LENSQNB,LENSQNC,ISGEI,ISFEI,ISGEIB,ISFEIB 161 INTEGER ISGEIC,ISFEIC,LENGTHGEI,LENGTHGEIB,LENGTHGEIC 162 INTEGER KGEI,KFEI,KGEIB,KFEIB,KGEIC,KFEIC,IADR,KLAMPB 163 INTEGER KLAMHB,KLAMPC,KLAMHC 164c 165 integer kx3am 166c 167 LOGICAL WB3X 168C 169#if defined (SYS_CRAY) 170 REAL OMEGA1(*),OMEGA2(*) 171 REAL WORK(LWORK) 172 REAL FREQL1,FREQL1R 173 REAL DDOT,XNORMVAL 174#else 175 DOUBLE PRECISION OMEGA1(*),OMEGA2(*) 176 DOUBLE PRECISION WORK(LWORK) 177 DOUBLE PRECISION FREQL1,FREQL1R 178 DOUBLE PRECISION DDOT,XNORMVAL 179#endif 180C 181 CALL QENTER('GNEW') 182C 183c write(lupri,*)'BEFORE ' 184c write(lupri,*)'omega1 before CC3_GMATNEW, LISTL ', LISTL 185c call PRINT_MATAI(OMEGA1,ISYRES) 186c xnormval = ddot(NT2AM(ISYRES),OMEGA2,1,OMEGA2,1) 187c write(lupri,*)'norm omega2 before CC3_GMATNEW ', xnormval 188c CALL OUTPAK(OMEGA2,NT1AMX,1,LUPRI) 189C 190C---------------------------------------------------- 191C Initialise character strings and open files 192C---------------------------------------------------- 193C 194 CDUMMY = ' ' 195C 196 SKIPGEI = .FALSE. 197C 198 LU3VI2 = -1 199 LU3FOP = -1 200 LU3FOP2 = -1 201C 202 CALL WOPEN2(LU3VI2,FN3VI2,64,0) 203 CALL WOPEN2(LU3FOP,FN3FOP,64,0) 204 CALL WOPEN2(LU3FOP2,FN3FOP2,64,0) 205 206C 207C------------------------- 208C Open temporary files 209C------------------------- 210C 211 LU4VB = -1 212 LU4VC = -1 213 LU4VBC = -1 214 LU3SRTR = -1 215 LUCKJDR = -1 216 LUDELDR = -1 217 LUDKBCR = -1 218 LUDKBCR4 = -1 219C 220 FN4VB = 'CC3_FMAT_TMP1' 221 FN4VC = 'CC3_FMAT_TMP2' 222 FN4VBC = 'CC3_FMAT_TMP3' 223 FN3SRTR = 'CC3_FMAT_TMP4' 224 FNCKJDR = 'CC3_FMAT_TMP5' 225 FNDELDR = 'CC3_FMAT_TMP6' 226 FNDKBCR = 'CC3_FMAT_TMP7' 227 FNDKBCR4 = 'CC3_FMAT_TMP8' 228C 229 IF (LVVVV) THEN 230 CALL WOPEN2(LU4VB,FN4VB,64,0) 231 CALL WOPEN2(LU4VC,FN4VC,64,0) 232 CALL WOPEN2(LU4VBC,FN4VBC,64,0) 233 END IF 234 CALL WOPEN2(LU3SRTR,FN3SRTR,64,0) 235 CALL WOPEN2(LUCKJDR,FNCKJDR,64,0) 236 CALL WOPEN2(LUDELDR,FNDELDR,64,0) 237 CALL WOPEN2(LUDKBCR,FNDKBCR,64,0) 238 CALL WOPEN2(LUDKBCR4,FNDKBCR4,64,0) 239C 240C------------------------------------------------------------ 241C lists handling 242C------------------------------------------------------------ 243C 244 LISTL0 = 'L0 ' 245 IDLSTL0 = 0 246C 247 IF (LISTL(1:3).EQ.'L1 ') THEN 248 249 ! get symmetry, frequency and integral label from common blocks 250 ! defined in ccl1rsp.h 251 ISYML1 = ISYLRZ(IDLSTL) 252 FREQL1 = FRQLRZ(IDLSTL) 253 LABELL1 = LRZLBL(IDLSTL) 254 LORXL1 = LORXLRZ(IDLSTL) 255c 256c write(lupri,*)'LABELL1 in CC3_GMATNEW ', LABELL1 257c 258 259 IF (LORXL1) CALL QUIT('NO ORBITAL RELAX. IN CC3_GMATNEW') 260 261 LISTL1R = 'R1 ' 262 IDLSTL1R = IR1TAMP(LABELL1,LORXL1,FREQL1,ISYML1) 263 ! get symmetry and frequency from common blocks 264 ! defined in ccl1rsp.h 265 ISYML1R = ISYLRT(IDLSTL1R) 266 FREQL1R = FRQLRT(IDLSTL1R) 267C 268 IF (ISYML1 .NE. ISYML1R) THEN 269 WRITE(LUPRI,*)'ISYML1: ', ISYML1 270 WRITE(LUPRI,*)'ISYML1R: ', ISYML1R 271 CALL QUIT('Symmetry mismatch in CC3_GMATNEW') 272 END IF 273C 274 IF (FREQL1R .NE. FREQL1) THEN 275 WRITE(LUPRI,*)'FREQL1R: ', FREQL1R 276 WRITE(LUPRI,*)'FREQL1: ', FREQL1 277 CALL QUIT('Frequency mismatch in CC3_GMATNEW') 278 END IF 279 280C 281 ELSE IF (LISTL(1:3).EQ.'LE ') THEN 282 ! get symmetry, frequency and integral label from common blocks 283 ! defined in ccl1rsp.h 284 ISYML1 = ILSTSYM(LISTL,IDLSTL) 285 FREQL1 = -EIGVAL(IDLSTL) 286 LABELL1 = '- none -' 287C 288 !we don't have any "right" vector entering a right hand side 289 LISTL1R = '---' 290 IDLSTL1R = -99 291C 292 ELSE IF (LISTL(1:3).EQ.'L0 ') THEN 293* CONTINUE 294 SKIPGEI = .TRUE. 295 ELSE 296 CALL QUIT('Unknown left list in CC3_GMATNEW') 297 END IF 298C 299 IF (.NOT.LVVVV) THEN 300 !Open files for N1MAT intermediates 301 LUGEI = -1 302 LUFEI = -1 303 LUN1 = -1 304 CALL WOPEN2(LUFEI,FNFEI,64,0) 305 CALL WOPEN2(LUN1,FNN1,64,0) 306 IF (.NOT.SKIPGEI) THEN 307 CALL WOPEN2(LUGEI,FNGEI,64,0) 308 END IF 309 !Open files for N1MAT^B intermediates 310 LUGEIB = -1 311 LUFEIB = -1 312 LUN1B = -1 313 CALL WOPEN2(LUFEIB,FNFEIB,64,0) 314 CALL WOPEN2(LUN1B,FNN1B,64,0) 315 IF (.NOT.SKIPGEI) THEN 316 CALL WOPEN2(LUGEIB,FNGEIB,64,0) 317 END IF 318 !Open files for N1MAT^C intermediates 319 LUGEIC = -1 320 LUFEIC = -1 321 LUN1C = -1 322 CALL WOPEN2(LUFEIC,FNFEIC,64,0) 323 CALL WOPEN2(LUN1C,FNN1C,64,0) 324 IF (.NOT.SKIPGEI) THEN 325 CALL WOPEN2(LUGEIC,FNGEIC,64,0) 326 END IF 327 END IF 328 329 330c IF (LISTB(1:3).EQ.'R1 ') THEN 331c ! get symmetry, frequency and integral label for right list 332c ! from common blocks defined in ccr1rsp.h 333c ISYMRB = ISYLRT(IDLSTB) 334c ELSE 335c CALL QUIT('Unknown right list in CC3_GMATNEW') 336c END IF 337 ISYMRB = ILSTSYM(LISTB,IDLSTB) 338 ISYMRC = ILSTSYM(LISTC,IDLSTC) 339C 340 ISYMBC = MULD2H(ISYMRB,ISYMRC) 341C 342 ISINT1 = ISYM0 343 ISINT2 = ISYM0 344C 345 IF (.NOT.LVVVV) THEN 346 IF ( (LISTL(1:3).EQ.'L1 ') .OR. (LISTL(1:3).EQ.'LE ')) THEN 347 !Symmetries for N1 and N2 intermediates 348 ISYMN1 = MULD2H(ISYML1,ISYM0) 349 ISYMN2 = MULD2H(ISYML1,ISYM0) 350 !Symmetries for N1^B and N2^B intermediates 351 ISYMN1B = MULD2H(ISYML1,ISYMRB) 352 ISYMN2B = MULD2H(ISYML1,ISYMRB) 353 !Symmetries for N1^C and N2^C intermediates 354 ISYMN1C = MULD2H(ISYML1,ISYMRC) 355 ISYMN2C = MULD2H(ISYML1,ISYMRC) 356 ELSE IF (LISTL(1:3).EQ.'L0 ') THEN 357 !Symmetries for N1 and N2 intermediates 358 ISYMN1 = MULD2H(ISYM0,ISYM0) 359 ISYMN2 = MULD2H(ISYM0,ISYM0) 360 !Symmetries for N1^B and N2^B intermediates 361 ISYMN1B = MULD2H(ISYM0,ISYMRB) 362 ISYMN2B = MULD2H(ISYM0,ISYMRB) 363 !Symmetries for N1^C and N2^C intermediates 364 ISYMN1C = MULD2H(ISYM0,ISYMRC) 365 ISYMN2C = MULD2H(ISYM0,ISYMRC) 366 ELSE 367 CALL QUIT('Unknown left list in CC3_GMATNEW(x)') 368 END IF 369 END IF 370C 371 IF (LVVVV) THEN 372 KRBJIA = 1 373 ELSE 374 KN2MAT = 1 375 KN2MATB = KN2MAT + NCKIJ(ISYMN2) 376 KN2MATC = KN2MATB + NCKIJ(ISYMN2B) 377 KRBJIA = KN2MATC + NCKIJ(ISYMN2C) 378 END IF 379 KLAMP0 = KRBJIA + NT2SQ(ISYRES) 380 KLAMH0 = KLAMP0 + NLAMDT 381 KFOCKD = KLAMH0 + NLAMDT 382 KFOCK0CK = KFOCKD + NORBTS 383 KT2TP = KFOCK0CK + NT1AMX 384 KL1AM = KT2TP + NT2SQ(ISYM0) 385 KL2TP = KL1AM + NT1AM(ISYM0) 386 KEND0 = KL2TP + NT2SQ(ISYM0) 387 LWRK0 = LWORK - KEND0 388C 389 IF ( (LISTL(1:3).EQ.'L1 ') .OR. (LISTL(1:3).EQ.'LE ') ) THEN 390 KL1L1 = KEND0 391 KL2L1 = KL1L1 + NT1AM(ISYML1) 392 KFOCKL1 = KL2L1 + NT2SQ(ISYML1) 393 KEND0 = KFOCKL1 + N2BST(ISYML1) 394 LWRK0 = LWORK - KEND0 395 END IF 396C 397 KT2RB = KEND0 398 KFOCK0 = KT2RB + NT2SQ(ISYMRB) 399 KEND0 = KFOCK0 + N2BST(ISYM0) 400 LWRK0 = LWORK - KEND0 401C 402 KT2RC = KEND0 403 KEND0 = KT2RC + NT2SQ(ISYMRC) 404 LWRK0 = LWORK - KEND0 405C 406 IF (.NOT.LVVVV) THEN 407 KINDSQN = KEND0 408 KINDSQNB = KINDSQN + (6*NCKIJ(ISYMN2) - 1)/IRAT + 1 409 KINDSQNC = KINDSQNB + (6*NCKIJ(ISYMN2B) - 1)/IRAT + 1 410 KEND0 = KINDSQNC + (6*NCKIJ(ISYMN2C) - 1)/IRAT + 1 411 LWRK0 = LWORK - KEND0 412 END IF 413C 414 IF (LWRK0 .LT. 0) THEN 415 WRITE(LUPRI,*) 'Memory available : ',LWORK 416 WRITE(LUPRI,*) 'Memory needed : ',KEND0 417 CALL QUIT('Insufficient space in CC3_GMATNEW') 418 END IF 419C 420 CALL DZERO(WORK(KRBJIA),NT2SQ(ISYRES)) 421C 422 IF (.NOT.LVVVV) THEN 423 CALL DZERO(WORK(KN2MAT),NCKIJ(ISYMN2)) 424 CALL DZERO(WORK(KN2MATB),NCKIJ(ISYMN2B)) 425 CALL DZERO(WORK(KN2MATC),NCKIJ(ISYMN2C)) 426C 427 !index array for N2 428 LENSQN = NCKIJ(ISYMN2) 429 CALL CC3_INDSQ(WORK(KINDSQN),LENSQN,ISYMN2) 430 !index array for N2^B 431 LENSQNB = NCKIJ(ISYMN2B) 432 CALL CC3_INDSQ(WORK(KINDSQNB),LENSQNB,ISYMN2B) 433 !index array for N2^C 434 LENSQNC = NCKIJ(ISYMN2C) 435 CALL CC3_INDSQ(WORK(KINDSQNC),LENSQNC,ISYMN2C) 436 END IF 437C 438C------------------------------------- 439C Read in lamdap and lamdh 440C------------------------------------- 441C 442 CALL GET_LAMBDA0(WORK(KLAMP0),WORK(KLAMH0),WORK(KEND0),LWRK0) 443C 444C--------------------------------------------------------------------- 445C Read zeroth-order AO Fock matrix from file and trasform it to 446C lambda basis 447C--------------------------------------------------------------------- 448C 449 CALL GET_FOCK0(WORK(KFOCK0),WORK(KLAMP0),WORK(KLAMH0),WORK(KEND0), 450 * LWRK0) 451C 452C--------------------------------------------------------------------- 453C Read the matrix the property integrals and trasform it to lambda 454C basis for L1 list 455C--------------------------------------------------------------------- 456C 457 IF (LISTL(1:3).EQ.'L1 ') THEN 458 CALL GET_FOCKX(WORK(KFOCKL1),LABELL1,IDLSTL,ISYML1, 459 * WORK(KLAMP0),ISYM0, 460 * WORK(KLAMH0),ISYM0,WORK(KEND0),LWRK0) 461 END IF 462C 463C------------------------------------- 464C Read T2 zeroth-order amplitudes 465C------------------------------------- 466C 467 IOPT = 2 468 CALL GET_T1_T2(IOPT,.FALSE.,DUMMY,WORK(KT2TP),'R0',0,ISYM0, 469 * WORK(KEND0),LWRK0) 470C 471 IF (LOCDBG) WRITE(LUPRI,*) 'Norm of T2TP ', 472 * DDOT(NT2SQ(ISYM0),WORK(KT2TP),1,WORK(KT2TP),1) 473C 474C-------------------------------------------- 475C Read L1 and L2 zeroth-order multipliers 476C-------------------------------------------- 477C 478 IOPT = 3 479 CALL GET_T1_T2(IOPT,.FALSE.,WORK(KL1AM),WORK(KL2TP),LISTL0, 480 * IDLSTL0,ISYM0,WORK(KEND0),LWRK0) 481C 482 IF (LOCDBG) WRITE(LUPRI,*) 'Norm of L2TP ', 483 * DDOT(NT2SQ(ISYM0),WORK(KL2TP),1,WORK(KL2TP),1) 484C 485C--------------------------------------------- 486C Read L1L1 and L2L1 multipliers (L1 list) 487C--------------------------------------------- 488C 489 IF ( (LISTL(1:3).EQ.'L1 ') .OR. (LISTL(1:3).EQ.'LE ') ) THEN 490 IOPT = 3 491 CALL GET_T1_T2(IOPT,.FALSE.,WORK(KL1L1),WORK(KL2L1),LISTL, 492 * IDLSTL,ISYML1,WORK(KEND0),LWRK0) 493 END IF 494C 495C------------------------------------- 496C Read T1R1 and L2R1 amplitudes 497C------------------------------------- 498C 499 IOPT = 2 500 CALL GET_T1_T2(IOPT,.TRUE.,DUMMY,WORK(KT2RB),LISTB, 501 * IDLSTB,ISYMRB,WORK(KEND0),LWRK0) 502C 503C------------------------------------- 504C Read T1R1 and L2R1 amplitudes 505C------------------------------------- 506C 507 IOPT = 2 508 CALL GET_T1_T2(IOPT,.TRUE.,DUMMY,WORK(KT2RC),LISTC, 509 * IDLSTC,ISYMRC,WORK(KEND0),LWRK0) 510C 511C---------------------------------------- 512C Integrals [H,T1Y] where Y is LISTB 513C---------------------------------------- 514C 515 ISINT1RB = MULD2H(ISINT1,ISYMRB) 516 ISINT1RC = MULD2H(ISINT1,ISYMRC) 517 ISINT2RB = MULD2H(ISINT2,ISYMRB) 518 ISINT2RC = MULD2H(ISINT2,ISYMRC) 519C 520 ISINT1BC = MULD2H(ISINT1,ISYMBC) 521C 522 CALL CC3_3BARINT(ISYMRB,LISTB,IDLSTB,ISYMRC,LISTC,IDLSTC, 523 * IDUMMY,CDUMMY,IDUMMY,.FALSE., 524 * WORK(KLAMP0),WORK(KLAMH0),WORK(KEND0),LWRK0, 525 * LU3SRTR,FN3SRTR,LUCKJDR,FNCKJDR) 526C 527 CALL CC3_SORT1(WORK(KEND0),LWRK0,2,ISINT1BC,LU3SRTR,FN3SRTR, 528 * LUDELDR,FNDELDR,IDUMMY,CDUMMY,IDUMMY,CDUMMY, 529 * IDUMMY,CDUMMY) 530C 531 CALL CC3_SINT(WORK(KLAMH0),WORK(KEND0),LWRK0,ISINT1BC, 532 * LUDELDR,FNDELDR,LUDKBCR,FNDKBCR) 533C 534 !with option IOPT = 2 it just calculates LUDKBCR4 needed 535 !in contractions 536 IOPT = 2 537 CALL CC3_TCME(DUMMY,ISINT1BC,WORK(KEND0),LWRK0,IDUMMY,CDUMMY, 538 * LUDKBCR,FNDKBCR,IDUMMY,CDUMMY,IDUMMY,CDUMMY, 539 * IDUMMY,CDUMMY,LUDKBCR4,FNDKBCR4,IOPT) 540C 541C--------------------------------------------------------------- 542C Read canonical orbital energies and delete frozen orbitals 543C in Fock diagonal, if required 544C--------------------------------------------------------------- 545C 546 CALL GET_ORBEN(WORK(KFOCKD),WORK(KEND0),LWRK0) 547C 548C-------------------------------------------- 549C Sort the Fock matrix to get F(ck) block 550C-------------------------------------------- 551C 552 CALL SORT_FOCKCK(WORK(KFOCK0CK),WORK(KFOCK0),ISYM0) 553C 554C---------------------------------------- 555C If we want to sum the T3 amplitudes 556C---------------------------------------- 557C 558 if (.false.) then 559 kx3am = kend0 560 kend0 = kx3am + nrhft*nrhft*nrhft*nvirt*nvirt*nvirt 561 call dzero(work(kx3am),nrhft*nrhft*nrhft*nvirt*nvirt*nvirt) 562 lwrk0 = lwork - kend0 563 if (lwrk0 .lt. 0) then 564 write(lupri,*) 'Memory available : ',lwork 565 write(lupri,*) 'Memory needed : ',kend0 566 call quit('Insufficient space (T3) in CC3_GMATNEW') 567 END IF 568 endif 569C 570C write(lupri,*) 'WBMAT after dzero' 571C call print_pt3(work(kx3am),ISYML1,4) 572C 573C-------------------------------------------------------------- 574C Calculate the normal g^0 integrals for 575C OOOO, OVVO and OOVV integrals used in contraction 576C (VVVV is stored on file LU4VB and read in in D-loop) 577C-------------------------------------------------------------- 578C 579 KBIOOOO = KEND0 580 KBIOVVO = KBIOOOO + N3ORHF(ISINT1RB) 581 KBIOOVV = KBIOVVO + NT2SQ(ISINT1RB) 582 KEND0 = KBIOOVV + NT2SQ(ISINT1RB) 583 LWRK0 = LWORK - KEND0 584C 585 IF (LWRK0 .LT. 0) THEN 586 CALL QUIT('Out of memory in CC3_GMATNEW (g^0[2o2v] kind)') 587 ENDIF 588C 589 CALL DZERO(WORK(KBIOVVO),NT2SQ(ISINT1RB)) 590 CALL DZERO(WORK(KBIOOVV),NT2SQ(ISINT1RB)) 591C 592 CALL CC3_INT(WORK(KBIOOOO),WORK(KBIOOVV),WORK(KBIOVVO), 593 * ISINT1RB,LU4VB,FN4VB, 594 * .TRUE.,LISTB,IDLSTB,ISYMRB, 595 * .FALSE.,'DUMMY',IDUMMY,IDUMMY, 596 * WORK(KEND0),LWRK0) 597C 598C 599C-------------------------------------------------------------- 600C Calculate the normal g^0 integrals for 601C OOOO, OVVO and OOVV integrals used in contraction 602C (VVVV is stored on file LU4VC and read in in D-loop) 603C-------------------------------------------------------------- 604C 605 KCIOOOO = KEND0 606 KCIOVVO = KCIOOOO + N3ORHF(ISINT1RC) 607 KCIOOVV = KCIOVVO + NT2SQ(ISINT1RC) 608 KEND0 = KCIOOVV + NT2SQ(ISINT1RC) 609 LWRK0 = LWORK - KEND0 610C 611 IF (LWRK0 .LT. 0) THEN 612 CALL QUIT('Out of memory in CC3_GMATNEW (g^0[2o2v] kind)') 613 ENDIF 614C 615 CALL DZERO(WORK(KCIOVVO),NT2SQ(ISINT1RC)) 616 CALL DZERO(WORK(KCIOOVV),NT2SQ(ISINT1RC)) 617C 618 CALL CC3_INT(WORK(KCIOOOO),WORK(KCIOOVV),WORK(KCIOVVO), 619 * ISINT1RC,LU4VC,FN4VC, 620 * .TRUE.,LISTC,IDLSTC,ISYMRC, 621 * .FALSE.,'DUMMY',IDUMMY,IDUMMY, 622 * WORK(KEND0),LWRK0) 623C 624C--------------------------------------------------------------- 625C Calculate the g^B integrals for 626C OOOO^B, OVVO^B and OOVV^B integrals used in contraction 627C (VVVV^B is stored on file LU4VBC and read in in D-loop) 628C--------------------------------------------------------------- 629C 630 KBCIOOOO = KEND0 631 KBCIOVVO = KBCIOOOO + N3ORHF(ISINT1BC) 632 KBCIOOVV = KBCIOVVO + NT2SQ(ISINT1BC) 633 KEND0 = KBCIOOVV + NT2SQ(ISINT1BC) 634 LWRK0 = LWORK - KEND0 635C 636 IF (LWRK0 .LT. 0) THEN 637 CALL QUIT('Out of memory in CC3_GMATNEW (g^B[2o2v] kind)') 638 ENDIF 639C 640 CALL DZERO(WORK(KBCIOOOO),N3ORHF(ISINT1BC)) 641 CALL DZERO(WORK(KBCIOVVO),NT2SQ(ISINT1BC)) 642 CALL DZERO(WORK(KBCIOOVV),NT2SQ(ISINT1BC)) 643C 644 CALL CC3_INT(WORK(KBCIOOOO),WORK(KBCIOOVV),WORK(KBCIOVVO), 645 * ISINT1BC,LU4VBC,FN4VBC, 646 * .TRUE.,LISTB,IDLSTB,ISYMRB, 647 * .TRUE.,LISTC,IDLSTC,ISYMRC, 648 * WORK(KEND0),LWRK0) 649C 650C----------------------------- 651C Memory allocation. 652C----------------------------- 653C 654 IF (LISTL(1:3).EQ.'L1 ') THEN 655 ISINT2L1R = MULD2H(ISYML1R,ISINT2) 656 ISYFCKL1R = MULD2H(ISYMOP,ISYML1R) 657 END IF 658C 659 KXIAJB = KEND0 660 KEND0 = KXIAJB + NT2AM(ISYM0) 661 662 KT3BOG1 = KEND0 663 KT3BOL1 = KT3BOG1 + NTRAOC(ISYM0) 664 KT3BOG2 = KT3BOL1 + NTRAOC(ISYM0) 665 KT3BOL2 = KT3BOG2 + NTRAOC(ISYM0) 666 KEND0 = KT3BOL2 + NTRAOC(ISYM0) 667 LWRK0 = LWORK - KEND0 668C 669 IF (LISTL(1:3).EQ.'L1 ') THEN 670 KW3XOGX1 = KEND0 671 KFOCKL1RCK = KW3XOGX1 + NTRAOC(ISINT2L1R) 672 KW3BXOGX1 = KFOCKL1RCK + NT1AM(ISYFCKL1R) 673 KW3BXOLX1 = KW3BXOGX1 + NTRAOC(ISINT2L1R) 674 KEND0 = KW3BXOLX1 + NTRAOC(ISINT2L1R) 675 LWRK0 = LWORK - KEND0 676C 677 KLAMPL1R = KEND0 678 KLAMHL1R = KLAMPL1R + NLAMDT 679 KEND0 = KLAMHL1R + NLAMDT 680 LWRK0 = LWORK - KEND0 681C 682 END IF 683C 684 KW3BXOG1 = KEND0 685 KW3BXOL1 = KW3BXOG1 + NTRAOC(ISYM0) 686 KEND0 = KW3BXOL1 + NTRAOC(ISYM0) 687 LWRK0 = LWORK - KEND0 688C 689C 690 IF (LISTL(1:3).EQ.'L1 ') THEN 691 KT1L1R = KEND0 692 KEND0 = KT1L1R + NT1AM(ISYML1R) 693 LWRK0 = LWORK - KEND0 694 END IF 695C 696 KINTOC = KEND0 697 KTROCR = KINTOC + NTOTOC(ISINT1BC) 698 KTROCR1 = KTROCR + NTRAOC(ISINT1BC) 699 KEND0 = KTROCR1 + NTRAOC(ISINT1BC) 700 LWRK0 = LWORK - KEND0 701C 702 IF (LWRK0 .LT. 0) THEN 703 WRITE(LUPRI,*) 'Memory available : ',LWORK 704 WRITE(LUPRI,*) 'Memory needed : ',KEND0 705 CALL QUIT('Insufficient space in CC3_GMATNEW(x)') 706 END IF 707 708C 709C------------------------ 710C Construct L(ia,jb). 711C------------------------ 712C 713 LENGTH = IRAT*NT2AM(ISYM0) 714 715 REWIND(LUIAJB) 716 CALL READI(LUIAJB,LENGTH,WORK(KXIAJB)) 717 718 CALL CCSD_TCMEPK(WORK(KXIAJB),1.0D0,ISYM0,1) 719C 720C--------------------------------------------------- 721C Prepare to construct the occupied integrals... 722C--------------------------------------------------- 723C 724C isint1 - symmetry of integrals in standard H, transformed 725C with LambdaH_0 726C ISINT1L1R - symmetry of integrals in standard H, transformed 727C with LambdaH_L1R 728 729 ISINT1 = 1 730C 731C-------------------------- 732C Get Lambda for right list depended on left LISTL list 733C-------------------------- 734C 735 IF (LISTL(1:3).EQ.'L1 ') THEN 736C 737 ISINT1L1R = MULD2H(ISINT1,ISYML1R) 738C 739 CALL GET_LAMBDAX(WORK(KLAMPL1R),WORK(KLAMHL1R),LISTL1R, 740 * IDLSTL1R, 741 * ISYML1R, 742 * WORK(KLAMP0),WORK(KLAMH0),WORK(KEND0),LWRK0) 743 744C 745C------------------------------------------------------------------ 746C Calculate the F^L1R matrix (kc elements evaluated and stored 747C as ck) 748C------------------------------------------------------------------ 749C 750 IOPT = 1 751 CALL GET_T1_T2(IOPT,.FALSE.,WORK(KT1L1R),DUMMY,LISTL1R, 752 * IDLSTL1R,ISYML1R,WORK(KEND0),LWRK0) 753 CALL CC3LR_MFOCK(WORK(KFOCKL1RCK),WORK(KT1L1R),WORK(KXIAJB), 754 * ISYFCKL1R) 755C 756 END IF 757C 758C----------------------------------------------------------------- 759C Construct occupied integrals which are required to calculate 760C t3bar_0 multipliers 761C----------------------------------------------------------------- 762C 763 CALL INTOCC_T3BAR0(LUTOC,FNTOC,WORK(KLAMH0),ISYM0,WORK(KT3BOG1), 764 * WORK(KT3BOL1),WORK(KT3BOG2),WORK(KT3BOL2), 765 * WORK(KEND0),LWRK0) 766C 767C----------------------------------------------------------------- 768C Construct occupied integrals which are required to calculate 769C t3bar_Y multipliers 770C----------------------------------------------------------------- 771C 772 IF (LISTL(1:3).EQ.'L1 ') THEN 773 CALL INTOCC_T3BARX(.FALSE., 774 * LUTOC,FNTOC,ISYMOP,WORK(KLAMH0),ISYM0, 775 * ISINT1, 776 * WORK(KLAMHL1R),ISYML1R,ISINT1L1R, 777 * WORK(KW3BXOG1), 778 * WORK(KW3BXOL1),WORK(KW3BXOGX1), 779 * WORK(KW3BXOLX1), 780 * WORK(KEND0),LWRK0) 781 ELSE IF (LISTL(1:3).EQ.'LE ') THEN 782 CALL INTOCC_T3BARX(.TRUE., 783 * LUTOC,FNTOC,ISYMOP,WORK(KLAMH0),ISYM0, 784 * ISINT1, 785 * DUMMY,IDUMMY,IDUMMY, 786 * WORK(KW3BXOG1), 787 * WORK(KW3BXOL1),DUMMY, 788 * DUMMY, 789 * WORK(KEND0),LWRK0) 790 END IF 791C 792C-------------------------------------------------------------------- 793C Read in R1-transformed occupied integrals used in contractions, 794C transform with lambda_P^0 and sort 795C-------------------------------------------------------------------- 796C 797 IOFF = 1 798 IF (NTOTOC(ISINT1BC) .GT. 0) THEN 799 CALL GETWA2(LUCKJDR,FNCKJDR,WORK(KINTOC),IOFF,NTOTOC(ISINT1BC)) 800 ENDIF 801C 802 CALL CC3_TROCC(WORK(KINTOC),WORK(KTROCR),WORK(KLAMP0), 803 * WORK(KEND0),LWRK0,ISINT1BC) 804C 805 CALL CCFOP_SORT(WORK(KTROCR),WORK(KTROCR1),ISINT1BC,1) 806C 807 CALL CC3_LSORT1(WORK(KTROCR),ISINT1BC,WORK(KEND0),LWRK0,5) 808C 809 DO ISYMD = 1,NSYM 810C 811 ISYCKBD0 = MULD2H(ISYMD,ISYM0) 812C 813 ISYCKBDRB = MULD2H(ISYMD,ISINT2RB) 814 ISYCKBDRC = MULD2H(ISYMD,ISINT2RC) 815 ISYCKBDRBC = MULD2H(ISYMD,ISINT1BC) 816C 817 IF (.NOT.LVVVV) THEN 818 !Symmetry of arrays needed to construct N1MAT 819 ISGEI = MULD2H(ISYMN1,ISYMD) 820 ISFEI = MULD2H(ISYMN1,ISYMD) 821 !Symmetry of arrays needed to construct N1MAT^B 822 ISGEIB = MULD2H(ISYMN1B,ISYMD) 823 ISFEIB = MULD2H(ISYMN1B,ISYMD) 824 !Symmetry of arrays needed to construct N1MAT^C 825 ISGEIC = MULD2H(ISYMN1C,ISYMD) 826 ISFEIC = MULD2H(ISYMN1C,ISYMD) 827 END IF 828C 829 IF (LISTL(1:3).EQ.'L1 ') THEN 830 ISYCKBDL1R = MULD2H(ISYMD,ISYML1R) 831 END IF 832C 833 IF (LVVVV) THEN 834 KVVVVB = KEND0 835 KVVVVC = KVVVVB + NMAABC(ISYCKBDRB) 836 KVVVVBC = KVVVVC + NMAABC(ISYCKBDRC) 837 KEND1 = KVVVVBC + NMAABC(ISYCKBDRBC) 838 LWRK1 = LWORK - KEND1 839 ELSE 840 IF (SKIPGEI) THEN !we want to have the dummy allocations 841 !to skip GEI, but keep the code simple. 842 LENGTHGEI = 1 843 LENGTHGEIB = 1 844 LENGTHGEIC = 1 845 ELSE 846 LENGTHGEI = NCKATR(ISGEI) 847 LENGTHGEIB = NCKATR(ISGEIB) 848 LENGTHGEIC = NCKATR(ISGEIC) 849 END IF 850 !Arrays needed to construct N1MAT 851 KGEI = KEND0 852 KFEI = KGEI + LENGTHGEI 853 KEND1 = KFEI + NCKATR(ISFEI) 854 LWRK1 = LWORK - KEND1 855C 856 !Arrays needed to construct N1MAT^B 857 KGEIB = KEND1 858 KFEIB = KGEIB + LENGTHGEIB 859 KEND1 = KFEIB + NCKATR(ISFEIB) 860 LWRK1 = LWORK - KEND1 861C 862 !Arrays needed to construct N1MAT^C 863 KGEIC = KEND1 864 KFEIC = KGEIC + LENGTHGEIC 865 KEND1 = KFEIC + NCKATR(ISFEIC) 866 LWRK1 = LWORK - KEND1 867C 868 END IF 869C 870 IF (LWRK1 .LT. 0) THEN 871 WRITE(LUPRI,*) 'Memory available : ',LWORK 872 WRITE(LUPRI,*) 'Memory needed : ',KEND1 873 CALL QUIT('Insufficient space (isymd-loop) in CC3_GMATNEW') 874 END IF 875C 876 DO D = 1,NVIR(ISYMD) 877C 878 KT3BVDL1 = KEND1 879 KT3BVDL2 = KT3BVDL1 + NCKATR(ISYCKBD0) 880 KT3BVDL3 = KT3BVDL2 + NCKATR(ISYCKBD0) 881 KEND3 = KT3BVDL3 + NCKATR(ISYCKBD0) 882 LWRK3 = LWORK - KEND3 883 884 KT3BVDG1 = KEND3 885 KT3BVDG2 = KT3BVDG1 + NCKATR(ISYCKBD0) 886 KT3BVDG3 = KT3BVDG2 + NCKATR(ISYCKBD0) 887 KEND3 = KT3BVDG3 + NCKATR(ISYCKBD0) 888 LWRK3 = LWORK - KEND3 889 890 IF ( (LISTL(1:3).EQ.'L1 ') .OR. (LISTL(1:3).EQ.'LE ') ) THEN 891 KW3BXVDG1 = KEND3 892 KW3BXVDG2 = KW3BXVDG1 + NCKATR(ISYCKBD0) 893 KW3BXVDL1 = KW3BXVDG2 + NCKATR(ISYCKBD0) 894 KW3BXVDL2 = KW3BXVDL1 + NCKATR(ISYCKBD0) 895 KEND3 = KW3BXVDL2 + NCKATR(ISYCKBD0) 896 LWRK3 = LWORK - KEND3 897 END IF 898C 899 IF (LISTL(1:3).EQ.'L1 ') THEN 900 KW3BXVDGX1 = KEND3 901 KW3BXVDGX2 = KW3BXVDGX1 + NCKATR(ISYCKBDL1R) 902 KW3BXVDLX1 = KW3BXVDGX2 + NCKATR(ISYCKBDL1R) 903 KW3BXVDLX2 = KW3BXVDLX1 + NCKATR(ISYCKBDL1R) 904 KEND3 = KW3BXVDLX2 + NCKATR(ISYCKBDL1R) 905 LWRK3 = LWORK - KEND3 906 END IF 907C 908 KTRVIR = KEND3 909 KTRVIR1 = KTRVIR + NCKATR(ISYCKBDRBC) 910 KEND4 = KTRVIR1 + NCKATR(ISYCKBDRBC) 911 LWRK4 = LWORK - KEND4 912 913 IF (LWRK4 .LT. 0) THEN 914 WRITE(LUPRI,*) 'Memory available : ',LWORK 915 WRITE(LUPRI,*) 'Memory needed : ',KEND4 916 CALL QUIT('Insufficient space (d-loop) in CC3_GMATNEW') 917 END IF 918 919 IF (.NOT.LVVVV) THEN 920 CALL DZERO(WORK(KFEI),NCKATR(ISFEI)) 921 CALL DZERO(WORK(KFEIB),NCKATR(ISFEIB)) 922 CALL DZERO(WORK(KFEIC),NCKATR(ISFEIC)) 923 CALL DZERO(WORK(KGEI),LENGTHGEI) 924 CALL DZERO(WORK(KGEIB),LENGTHGEIB) 925 CALL DZERO(WORK(KGEIC),LENGTHGEIC) 926 END IF 927 928 929 IF (LVVVV) THEN 930C 931C-------------------------------------------------------------------- 932C Read in g^B_{vvvv} (used in contraction) for a given D 933C-------------------------------------------------------------------- 934C 935 IF (NMAABC(ISYCKBDRB) .GT. 0) THEN 936 IOFF = I3VVIR(ISYCKBDRB,ISYMD) 937 * + NMAABC(ISYCKBDRB)*(D-1) 938 * + 1 939 CALL GETWA2(LU4VB,FN4VB,WORK(KVVVVB),IOFF, 940 * NMAABC(ISYCKBDRB)) 941 ENDIF 942C 943C-------------------------------------------------------------------- 944C Read in g^C_{vvvv} (used in contraction) for a given D 945C-------------------------------------------------------------------- 946C 947 IF (NMAABC(ISYCKBDRC) .GT. 0) THEN 948 IOFF = I3VVIR(ISYCKBDRC,ISYMD) 949 * + NMAABC(ISYCKBDRC)*(D-1) 950 * + 1 951 CALL GETWA2(LU4VC,FN4VC,WORK(KVVVVC),IOFF, 952 * NMAABC(ISYCKBDRC)) 953 ENDIF 954C 955C--------------------------------------------------------------------- 956C Read in g^BC_{vvvv} (used in contraction) for a given D 957C--------------------------------------------------------------------- 958C 959 IF (NMAABC(ISYCKBDRBC) .GT. 0) THEN 960 IOFF = I3VVIR(ISYCKBDRBC,ISYMD) 961 * + NMAABC(ISYCKBDRBC)*(D-1) 962 * + 1 963 CALL GETWA2(LU4VBC,FN4VBC,WORK(KVVVVBC),IOFF, 964 * NMAABC(ISYCKBDRBC)) 965 ENDIF 966C 967 END IF !LVVVV 968C 969C----------------------------------------------------------------------- 970C Construct virtual integrals (for fixed D) which are required 971C to calculate t3bar_0 multipliers 972C----------------------------------------------------------------------- 973C 974 CALL INTVIR_T3BAR0_D(LU3FOPX,FN3FOPX,LU3FOP2X,FN3FOP2X, 975 * LUDKBC3,FNDKBC3,LU3VI,FN3VI,ISYM0, 976 * WORK(KT3BVDL1),WORK(KT3BVDG1), 977 * WORK(KT3BVDG2),WORK(KT3BVDL2), 978 * WORK(KT3BVDG3),WORK(KT3BVDL3), 979 * WORK(KLAMP0),ISYMD,D,WORK(KEND4),LWRK4) 980C 981C----------------------------------------------------------------------- 982C Construct virtual integrals (for fixed D) which are required 983C to calculate t3bar_X multipliers 984C----------------------------------------------------------------------- 985C 986 IF (LISTL(1:3).EQ.'L1 ') THEN 987 CALL INTVIR_T3BARX_D(.FALSE., 988 * ISYMOP,LU3VI,FN3VI,LU3VI2,FN3VI2, 989 * LU3FOP,FN3FOP,LU3FOP2,FN3FOP2, 990 * WORK(KW3BXVDGX1),WORK(KW3BXVDG1), 991 * WORK(KW3BXVDGX2),WORK(KW3BXVDG2), 992 * WORK(KW3BXVDLX1),WORK(KW3BXVDL1), 993 * WORK(KW3BXVDLX2),WORK(KW3BXVDL2), 994 * WORK(KLAMPL1R),ISYML1R,WORK(KLAMP0), 995 * ISYM0,ISYMD,D,WORK(KEND4),LWRK4) 996 ELSE IF (LISTL(1:3).EQ.'LE ') THEN 997 CALL INTVIR_T3BARX_D(.TRUE., 998 * ISYMOP,LU3VI,FN3VI,LU3VI2,FN3VI2, 999 * LU3FOP,FN3FOP,LU3FOP2,FN3FOP2, 1000 * DUMMY,WORK(KW3BXVDG1), 1001 * DUMMY,WORK(KW3BXVDG2), 1002 * DUMMY,WORK(KW3BXVDL1), 1003 * DUMMY,WORK(KW3BXVDL2), 1004 * DUMMY,IDUMMY,WORK(KLAMP0), 1005 * ISYM0,ISYMD,D,WORK(KEND4),LWRK4) 1006 END IF 1007C 1008C--------------------------------------------------------------------- 1009C Read and sort R1-transformed integrals used in contraction. 1010C--------------------------------------------------------------------- 1011C 1012 IF (NCKATR(ISYCKBDRBC) .GT. 0) THEN 1013 IOFF = ICKBD(ISYCKBDRBC,ISYMD)+NCKATR(ISYCKBDRBC)*(D - 1) 1014 * + 1 1015 CALL GETWA2(LUDKBCR,FNDKBCR,WORK(KTRVIR),IOFF, 1016 & NCKATR(ISYCKBDRBC)) 1017 ENDIF 1018C 1019 IF (LWRK4 .LT. NCKATR(ISYCKBDRBC)) THEN 1020 CALL QUIT('Insufficient space for allocation in '// 1021 & 'CC3_GMATNEW (TRVI)') 1022 END IF 1023C 1024 CALL CCSDT_SRVIR3(WORK(KTRVIR),WORK(KEND4),ISYMD,D,ISINT1BC) 1025C 1026C 1027 IF (NCKATR(ISYCKBDRBC) .GT. 0) THEN 1028 IOFF = ICKBD(ISYCKBDRBC,ISYMD)+NCKATR(ISYCKBDRBC)*(D - 1) 1029 * + 1 1030 CALL GETWA2(LUDKBCR4,FNDKBCR4,WORK(KTRVIR1),IOFF, 1031 & NCKATR(ISYCKBDRBC)) 1032 ENDIF 1033C 1034 IF (LWRK4 .LT. NCKATR(ISYCKBDRBC)) THEN 1035 CALL QUIT('Insufficient space for allocation in '// 1036 & 'CC3_GMATNEW (TRVI1)') 1037 END IF 1038C 1039 CALL CCSDT_SRVIR3(WORK(KTRVIR1),WORK(KEND4),ISYMD,D, 1040 * ISINT1BC) 1041C 1042 1043C 1044 DO ISYMB = 1,NSYM 1045C 1046 ISYALJB0 = MULD2H(ISYMB,ISYM0) 1047 ISYALJD0 = MULD2H(ISYMD,ISYM0) 1048 ISYMBD = MULD2H(ISYMD,ISYMB) 1049 ISCKIJ = MULD2H(ISYMBD,ISYM0) 1050 ISYCKD = MULD2H(ISYM0,ISYMB) 1051 IF ( (LISTL(1:3).EQ.'L1 ') 1052 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1053 ISYALJBL1 = MULD2H(ISYMB,ISYML1) 1054 ISYALJDL1 = MULD2H(ISYMD,ISYML1) 1055 ISWBMAT = MULD2H(ISCKIJ,ISYML1) 1056 END IF 1057C 1058 KSMAT2 = KEND4 1059 KUMAT2 = KSMAT2 + NCKIJ(ISCKIJ) 1060 KDIAG = KUMAT2 + NCKIJ(ISCKIJ) 1061 KINDSQ = KDIAG + NCKIJ(ISCKIJ) 1062 KEND5 = KINDSQ + (6*NCKIJ(ISCKIJ) - 1)/IRAT + 1 1063 1064 IF ( (LISTL(1:3).EQ.'L1 ') 1065 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1066 KDIAGWB = KEND5 1067 KINDSQWB = KDIAGWB + NCKIJ(ISWBMAT) 1068 KEND5 = KINDSQWB + (6*NCKIJ(ISWBMAT) - 1)/IRAT + 1 1069 END IF 1070 1071 KINDEX = KEND5 1072 KINDEX2 = KINDEX + (NCKI(ISYALJB0) - 1)/IRAT + 1 1073 KEND5 = KINDEX2 + (NCKI(ISYALJD0) - 1)/IRAT + 1 1074C 1075 IF ( (LISTL(1:3).EQ.'L1 ') 1076 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1077 KINDEXBL1 = KEND5 1078 KINDEXDL1 = KINDEXBL1 + (NCKI(ISYALJBL1)-1)/IRAT + 1 1079 KTMAT = KINDEXDL1 + (NCKI(ISYALJDL1)-1)/IRAT + 1 1080 KW3BMAT = KTMAT + MAX(NCKIJ(ISCKIJ),NCKIJ(ISWBMAT)) 1081 KEND5 = KW3BMAT + NCKIJ(ISWBMAT) 1082 LWRK5 = LWORK - KEND5 1083 ELSE 1084 KTMAT = KEND5 1085 KEND5 = KTMAT + NCKIJ(ISCKIJ) 1086 END IF 1087 1088 KT3BVBG1 = KEND5 1089 KT3BVBG2 = KT3BVBG1 + NCKATR(ISYCKD) 1090 KT3BVBG3 = KT3BVBG2 + NCKATR(ISYCKD) 1091 KEND5 = KT3BVBG3 + NCKATR(ISYCKD) 1092 LWRK5 = LWORK - KEND5 1093 1094 KSMAT4 = KEND5 1095 KUMAT4 = KSMAT4 + NCKIJ(ISCKIJ) 1096 KT3BVBL1 = KUMAT4 + NCKIJ(ISCKIJ) 1097 KT3BVBL2 = KT3BVBL1 + NCKATR(ISYCKD) 1098 KT3BVBL3 = KT3BVBL2 + NCKATR(ISYCKD) 1099 KEND5 = KT3BVBL3 + NCKATR(ISYCKD) 1100 LWRK5 = LWORK - KEND5 1101C 1102 IF (LWRK5 .LT. 0) THEN 1103 WRITE(LUPRI,*) 'Memory available : ',LWORK 1104 WRITE(LUPRI,*) 'Memory needed : ',KEND5 1105 CALL QUIT('Insufficient space(isymb-loop) 1106 * in CC3_GMATNEW') 1107 END IF 1108C 1109C-------------------------------------------- 1110C Construct part of the diagonal 1111C-------------------------------------------- 1112C 1113 CALL CC3_DIAG(WORK(KDIAG),WORK(KFOCKD),ISCKIJ) 1114 IF ( (LISTL(1:3).EQ.'L1 ') 1115 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1116 CALL CC3_DIAG(WORK(KDIAGWB),WORK(KFOCKD),ISWBMAT) 1117 END IF 1118C 1119C------------------------------------ 1120C Construct index arrays 1121C------------------------------------ 1122C 1123 LENSQ = NCKIJ(ISCKIJ) 1124 CALL CC3_INDSQ(WORK(KINDSQ),LENSQ,ISCKIJ) 1125 IF ( (LISTL(1:3).EQ.'L1 ') 1126 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1127 LENSQWB = NCKIJ(ISWBMAT) 1128 CALL CC3_INDSQ(WORK(KINDSQWB),LENSQWB,ISWBMAT) 1129 END IF 1130C 1131 CALL CC3_INDEX(WORK(KINDEX),ISYALJB0) 1132 CALL CC3_INDEX(WORK(KINDEX2),ISYALJD0) 1133 IF ( (LISTL(1:3).EQ.'L1 ') 1134 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1135 CALL CC3_INDEX(WORK(KINDEXBL1),ISYALJBL1) 1136 CALL CC3_INDEX(WORK(KINDEXDL1),ISYALJDL1) 1137 END IF 1138C 1139 DO B = 1,NVIR(ISYMB) 1140C 1141 IF ( (LISTL(1:3).EQ.'L1 ') 1142 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1143 CALL DZERO(WORK(KW3BMAT),NCKIJ(ISWBMAT)) 1144 END IF 1145C 1146C----------------------------------------------------------------------- 1147C Construct virtual integrals (for fixed B) which are 1148C required to calculate t3bar_0 multipliers 1149C (the same routine as in d-loop is used) 1150C----------------------------------------------------------------------- 1151C 1152 CALL INTVIR_T3BAR0_D(LU3FOPX,FN3FOPX,LU3FOP2X, 1153 * FN3FOP2X,LUDKBC3,FNDKBC3, 1154 * LU3VI,FN3VI,ISYM0,WORK(KT3BVBL1), 1155 * WORK(KT3BVBG1),WORK(KT3BVBG2), 1156 * WORK(KT3BVBL2),WORK(KT3BVBG3), 1157 * WORK(KT3BVBL3),WORK(KLAMP0), 1158 * ISYMB,B,WORK(KEND5),LWRK5) 1159C 1160C--------------------------------------------------------- 1161C Get T3bar_BD multipliers (using S and U) 1162C--------------------------------------------------------- 1163C 1164C 1165 IF ( (LISTL(1:3).EQ.'L0 ') 1166 * .OR. (LISTL(1:3).EQ.'L1 ') ) THEN 1167 CALL GET_T3BAR0_BD(ISYM0,WORK(KL1AM),ISYM0, 1168 * WORK(KL2TP),ISYM0,WORK(KTMAT), 1169 * WORK(KFOCK0CK),WORK(KFOCKD), 1170 * WORK(KDIAG),WORK(KXIAJB),ISYM0, 1171 * ISYM0,WORK(KINDSQ),LENSQ, 1172 * WORK(KSMAT2),WORK(KT3BVDG1), 1173 * WORK(KT3BVDG2),WORK(KT3BVDL1), 1174 * WORK(KT3BVDL2),WORK(KT3BOG1), 1175 * WORK(KT3BOL1),WORK(KINDEX), 1176 * WORK(KSMAT4),WORK(KT3BVBG1), 1177 * WORK(KT3BVBG2),WORK(KT3BVBL1), 1178 * WORK(KT3BVBL2),WORK(KINDEX2), 1179 * WORK(KUMAT2),WORK(KT3BVDG3), 1180 * WORK(KT3BVDL3),WORK(KT3BOG2), 1181 * WORK(KT3BOL2),WORK(KUMAT4), 1182 * WORK(KT3BVBG3),WORK(KT3BVBL3), 1183 * ISYMB,B,ISYMD,D,ISCKIJ, 1184 * WORK(KEND5),LWRK5) 1185c 1186c call sum_pt3(work(KTMAT),isymb,b,isymd,d, 1187c * ISYM0,work(kx3am),5) 1188 1189 END IF 1190C 1191C---------------------------------------------------- 1192C Get T3barX_BD multipliers (using W) 1193C---------------------------------------------------- 1194C 1195 IF (LISTL(1:3).EQ.'L1 ') THEN 1196 CALL GET_T3BARX_BD(.FALSE., 1197 * WORK(KTMAT),ISCKIJ,WORK(KFOCKL1), 1198 * ISYML1,WORK(KW3BMAT),ISWBMAT, 1199 * WORK(KL2TP),ISYM0,WORK(KFOCKL1RCK), 1200 * ISYFCKL1R,WORK(KW3BXVDLX2), 1201 * WORK(KW3BXVDLX1),WORK(KW3BXVDGX2), 1202 * WORK(KW3BXVDGX1),WORK(KW3BXOLX1), 1203 * WORK(KW3BXOGX1),ISINT2L1R, 1204 * WORK(KINDEX),WORK(KINDEX2), 1205 * WORK(KINDSQWB),LENSQWB,WORK(KL2L1), 1206 * ISYML1,WORK(KFOCK0CK),ISYM0, 1207 * WORK(KW3BXVDL2),WORK(KW3BXVDL1), 1208 * WORK(KW3BXVDG2),WORK(KW3BXVDG1), 1209 * WORK(KW3BXOL1),WORK(KW3BXOG1), 1210 * ISINT2,WORK(KINDEXBL1), 1211 * WORK(KINDEXDL1),WORK(KL1L1),ISYML1, 1212 * WORK(KXIAJB),ISINT1,-FREQL1, 1213 * WORK(KDIAGWB),WORK(KFOCKD),B,ISYMB, 1214 * D,ISYMD,ISYML1,WORK(KEND5),LWRK5) 1215c 1216c call sum_pt3(work(KW3BMAT),isymb,b,isymd,d, 1217c * ISWBMAT,work(kx3am),4) 1218 END IF 1219C 1220C 1221C--------------------------------------------------------------- 1222C Get E3bar_BD left eigenvectors (using S and U) 1223C--------------------------------------------------------------- 1224C 1225C 1226 IF (LISTL(1:3).EQ.'LE ') THEN 1227 CALL GET_L3BAR_BD(WORK(KTMAT), 1228 * WORK(KW3BMAT),ISWBMAT, 1229 * WORK(KINDSQWB),LENSQWB,WORK(KL2L1), 1230 * ISYML1,WORK(KFOCK0CK),ISYM0, 1231 * WORK(KW3BXVDL2),WORK(KW3BXVDL1), 1232 * WORK(KW3BXVDG2),WORK(KW3BXVDG1), 1233 * WORK(KW3BXOL1),WORK(KW3BXOG1), 1234 * ISINT2,WORK(KINDEXBL1), 1235 * WORK(KINDEXDL1),WORK(KL1L1),ISYML1, 1236 * WORK(KXIAJB),ISINT1,-FREQL1, 1237 * WORK(KDIAGWB),WORK(KFOCKD),B,ISYMB, 1238 * D,ISYMD,ISYML1,WORK(KEND5),LWRK5) 1239 1240c 1241c call sum_pt3(work(KTMAT),isymb,b,isymd,d, 1242c * ISYM0,work(kx3am),5) 1243 1244 END IF 1245 1246C 1247C------------------------------------------------------------------- 1248C Set up variables depending on the LISTL (L1 or L0) 1249C------------------------------------------------------------------- 1250C 1251C If WB3X = .TRUE. (i.e, LISTL = L1), then KTB contains wbar_X 1252C else 1253C we get tbar_0 in KTB 1254C 1255C Note Wbar_X and Tbar_0 is stored in the same way for fixed B and D 1256C 1257C Inside the routines CC3_W3_OMEGA1, CC3_W3_CY2V, CC3_W3_CY2O terms 1258C are selected depending on WB3X 1259 1260 IF ( (LISTL(1:3).EQ.'L1 ') 1261 * .OR. (LISTL(1:3).EQ.'LE ') ) THEN 1262 WB3X = .TRUE. 1263 ISTBD = ISWBMAT 1264 KTB = KW3BMAT 1265 LENSQTB = LENSQWB 1266 KINDSQTB = KINDSQWB 1267 ELSE 1268 WB3X = .FALSE. 1269 ISTBD = ISCKIJ 1270 !KTMAT will be destroyed so an extra array for 1271 !zeroth-order multipliers is needed 1272 KTB = KEND5 1273 KEND5 = KTB + NCKIJ(ISCKIJ) 1274 LWRK5 = LWORK - KEND5 1275 IF (LWRK5 .LT. 0) THEN 1276 WRITE(LUPRI,*) 'Memory available : ',LWORK 1277 WRITE(LUPRI,*) 'Memory needed : ',KEND5 1278 CALL QUIT('Too little space(setup) 1279 * in CC3_GMATNEW') 1280 END IF 1281 CALL DCOPY(NCKIJ(ISCKIJ),WORK(KTMAT),1,WORK(KTB),1) 1282 1283 LENSQTB = LENSQ 1284 KINDSQTB = KINDSQ 1285 END IF 1286 ISTB = MULD2H(ISTBD,ISYMBD) 1287C 1288 !Check the consistency of logical flags 1289 IF (WB3X.AND.(.NOT.SKIPGEI)) THEN 1290 CONTINUE 1291 ELSE IF ((.NOT.WB3X).AND.SKIPGEI) THEN 1292 CONTINUE 1293 ELSE 1294 WRITE(LUPRI,*)'WB3X = ',WB3X 1295 WRITE(LUPRI,*)'SKIPGEI = ',SKIPGEI 1296 WRITE(LUPRI,*)'WB3X and SKIPGEI should be opposite' 1297 CALL QUIT('Inconsistent flags in CC3_GMATNEW') 1298 END IF 1299C 1300C------------------------------------------------- 1301C Calculate <L3|[H^BC,\tau_nu_2]|HF> 1302C------------------------------------------------- 1303C 1304 CALL CC3_W3_CY2V(OMEGA2,ISYRES,WORK(KRBJIA), 1305 * WORK(KTB),ISTBD,WORK(KTMAT), 1306 * WORK(KTRVIR),WORK(KTRVIR1),ISINT1BC, 1307 * WORK(KEND5),LWRK5,WORK(KINDSQTB), 1308 * LENSQTB,ISYMB,B,ISYMD,D,WB3X) 1309C 1310 CALL CC3_W3_CY2O(OMEGA2,ISYRES,WORK(KTB), 1311 * ISTBD, 1312 * WORK(KTMAT),WORK(KTROCR), 1313 * WORK(KTROCR1),ISINT1BC, 1314 * WORK(KEND5),LWRK5,WORK(KINDSQTB), 1315 * LENSQTB,ISYMB,B,ISYMD,D,WB3X) 1316C 1317 IF (LVVVV) THEN 1318 !Calculate <L3|[H^B,T^{C}_{2}],tau_{1}]|HF> 1319 CALL CC3_W3_OMEGA1(OMEGA1,ISYRES,WORK(KTB), 1320 * WORK(KTMAT),ISTB, 1321 * WORK(KBIOOOO),WORK(KBIOVVO), 1322 * WORK(KBIOOVV),WORK(KVVVVB), 1323 * ISINT1RB,WORK(KT2RC),ISYMRC, 1324 * WORK(KEND5),LWRK5, 1325 * LENSQTB,WORK(KINDSQTB), 1326 * ISYMB,B,ISYMD,D,WB3X) 1327C 1328 !Calculate <L3|[H^C,T^{B}_{2}],tau_{1}]|HF> 1329 CALL CC3_W3_OMEGA1(OMEGA1,ISYRES,WORK(KTB), 1330 * WORK(KTMAT),ISTB, 1331 * WORK(KCIOOOO),WORK(KCIOVVO), 1332 * WORK(KCIOOVV),WORK(KVVVVC), 1333 * ISINT1RC,WORK(KT2RB),ISYMRB, 1334 * WORK(KEND5),LWRK5, 1335 * LENSQTB,WORK(KINDSQTB), 1336 * ISYMB,B,ISYMD,D,WB3X) 1337C 1338 !Calculate <L3|[H^BC,T^{0}_{2}],tau_{1}]|HF> 1339 CALL CC3_W3_OMEGA1(OMEGA1,ISYRES,WORK(KTB), 1340 * WORK(KTMAT),ISTB, 1341 * WORK(KBCIOOOO),WORK(KBCIOVVO), 1342 * WORK(KBCIOOVV),WORK(KVVVVBC), 1343 * ISINT1BC,WORK(KT2TP),ISYM0, 1344 * WORK(KEND5),LWRK5, 1345 * LENSQTB,WORK(KINDSQTB), 1346 * ISYMB,B,ISYMD,D,WB3X) 1347 ELSE 1348C 1349 CALL DSCAL(NCKIJ(ISTBD),-1.0D0,WORK(KTB),1) 1350 1351 !<L3|[H^B,T^{C}_{2}],tau_{1}]|HF> in terms of N1 and N2 1352 CALL WT2_N1N2(WORK(KTB),ISTB, 1353 * WORK(KT2RC),ISYMRC, 1354 * WORK(KGEIC),WORK(KFEIC), 1355 * ISYMN1C, 1356 * WORK(KN2MATC),ISYMN2C, 1357 * B,ISYMB,D,ISYMD, 1358 * WORK(KINDSQTB),LENSQTB, 1359 * WORK(KINDSQNC),LENSQNC, 1360 * WORK(KEND5),LWRK5, 1361 * WB3X) 1362C 1363 !<L3|[H^C,T^{B}_{2}],tau_{1}]|HF>in terms of N1 and N2 1364 CALL WT2_N1N2(WORK(KTB),ISTB, 1365 * WORK(KT2RB),ISYMRB, 1366 * WORK(KGEIB),WORK(KFEIB), 1367 * ISYMN1B, 1368 * WORK(KN2MATB),ISYMN2B, 1369 * B,ISYMB,D,ISYMD, 1370 * WORK(KINDSQTB),LENSQTB, 1371 * WORK(KINDSQNB),LENSQNB, 1372 * WORK(KEND5),LWRK5, 1373 * WB3X) 1374C 1375 !<L3|[H^BC,T^{0}_{2}],tau_{1}]|HF>in terms of N1 and N2 1376 CALL WT2_N1N2(WORK(KTB),ISTB, 1377 * WORK(KT2TP),ISYM0, 1378 * WORK(KGEI),WORK(KFEI), 1379 * ISYMN1, 1380 * WORK(KN2MAT),ISYMN2, 1381 * B,ISYMB,D,ISYMD, 1382 * WORK(KINDSQTB),LENSQTB, 1383 * WORK(KINDSQN),LENSQN, 1384 * WORK(KEND5),LWRK5, 1385 * WB3X) 1386 1387 END IF 1388 1389C 1390 ENDDO ! B 1391 ENDDO ! ISYMB 1392C 1393 IF (.NOT.LVVVV) THEN 1394C 1395C ---------------------------------------------------------- 1396C Put KGEI(ge,i)^F and KFEI(fe,i)^G (which are intermediates 1397C for N1MAT(fge,i) ) to files (for fixed F=D and G=D). 1398C ---------------------------------------------------------- 1399 1400 IF (.NOT.SKIPGEI) THEN 1401 !Put KGEI to file as (gei,F) (fixed F corresponds to D) 1402 IADR = ICKBD(ISGEI,ISYMD) + NCKATR(ISGEI)*(D-1) + 1 1403 CALL PUTWA2(LUGEI,FNGEI,WORK(KGEI),IADR,NCKATR(ISGEI)) 1404 END IF 1405C 1406 !Put KFEI to file as (fei,G) (fixed G corresponds to D) 1407 IADR = ICKBD(ISFEI,ISYMD) + NCKATR(ISFEI)*(D-1) + 1 1408 CALL PUTWA2(LUFEI,FNFEI,WORK(KFEI),IADR,NCKATR(ISFEI)) 1409C 1410 IF (.NOT.SKIPGEI) THEN 1411 !Put KGEIB to file as (gei,F) (fixed F corresponds to D) 1412 IADR = ICKBD(ISGEIB,ISYMD) + NCKATR(ISGEIB)*(D-1) + 1 1413 CALL PUTWA2(LUGEIB,FNGEIB,WORK(KGEIB),IADR,NCKATR(ISGEIB)) 1414 END IF 1415C 1416 !Put KFEIB to file as (fei,G) (fixed G corresponds to D) 1417 IADR = ICKBD(ISFEIB,ISYMD) + NCKATR(ISFEIB)*(D-1) + 1 1418 CALL PUTWA2(LUFEIB,FNFEIB,WORK(KFEIB),IADR,NCKATR(ISFEIB)) 1419C 1420 IF (.NOT.SKIPGEI) THEN 1421 !Put KGEIC to file as (gei,F) (fixed F corresponds to D) 1422 IADR = ICKBD(ISGEIC,ISYMD) + NCKATR(ISGEIC)*(D-1) + 1 1423 CALL PUTWA2(LUGEIC,FNGEIC,WORK(KGEIC),IADR,NCKATR(ISGEIC)) 1424 END IF 1425C 1426 !Put KFEIC to file as (fei,G) (fixed G corresponds to D) 1427 IADR = ICKBD(ISFEIC,ISYMD) + NCKATR(ISFEIC)*(D-1) + 1 1428 CALL PUTWA2(LUFEIC,FNFEIC,WORK(KFEIC),IADR,NCKATR(ISFEIC)) 1429C 1430 END IF 1431C 1432 ENDDO ! D 1433 ENDDO ! ISYMD 1434 1435C 1436C------------------------------------------------------ 1437C Accumulate RBJIA from the virtual contribution 1438C in OMEGA2 1439C------------------------------------------------------ 1440C 1441 CALL CC3_RBJIA(OMEGA2,ISYRES,WORK(KRBJIA)) 1442C 1443 IF (.NOT.LVVVV) THEN 1444C 1445 KLAMPB = KEND0 1446 KLAMHB = KLAMPB + NGLMDT(ISYMRB) 1447 KLAMPC = KLAMHB + NGLMDT(ISYMRB) 1448 KLAMHC = KLAMPC + NGLMDT(ISYMRC) 1449 KEND0 = KLAMHC + NGLMDT(ISYMRC) 1450 LWRK0 = LWORK - KEND0 1451 IF (LWRK0 .LT. 0) THEN 1452 WRITE(LUPRI,*) 'Memory available : ',LWORK 1453 WRITE(LUPRI,*) 'Memory needed : ',KEND0 1454 CALL QUIT('Insufficient space in CC3_GMATNEW(final)') 1455 END IF 1456C 1457 !Lambda^B (particle) needed in backtransformation of N1^C 1458 CALL GET_LAMBDAX(WORK(KLAMPB),WORK(KLAMHB),LISTB,IDLSTB,ISYMRB, 1459 * WORK(KLAMP0),WORK(KLAMH0), 1460 * WORK(KEND0),LWRK0) 1461C 1462 !Read (gei,F) and (fei,G) intermediates from files 1463 !add them and put the result to a file as (fge,I) 1464 CALL N1_RESORT(ISYMN1C,LUN1C,FNN1C,LUGEIC,FNGEIC,LUFEIC,FNFEIC, 1465 * WORK(KEND0),LWRK0,SKIPGEI) 1466C 1467 !Calculate <T3|[[H^B,T2^C],tau_ai]|HF> except VVVV contribution 1468 CALL N1N2_G(LUN1C,FNN1C, 1469 * ISYMN1C, 1470 * WORK(KN2MATC),ISYMN2C, 1471 * WORK(KBIOVVO),WORK(KBIOOVV), 1472 * WORK(KBIOOOO),ISINT1RB, 1473 * OMEGA1,ISYRES, 1474 * WORK(KINDSQNC),LENSQNC, 1475 * WORK(KEND0),LWRK0) 1476C 1477 !Calculate VVVV contribution to <T3|[[H^B,T2^C],tau_ai]|HF> 1478 IOPT = 1 !T1 one-index backtransformation 1479 CALL N1_GV4(IOPT, 1480 * LUN1C,FNN1C, 1481 * ISYMN1C, 1482 * WORK(KLAMPB),ISYMRB, 1483 * WORK(KLAMP0),1, 1484 * WORK(KLAMH0),1, 1485 * WORK(KLAMH0),1, 1486 * OMEGA1,ISYRES, 1487 * WORK(KEND0),LWRK0) 1488C 1489C REPEAT THINGS FOR N1MATB 1490C 1491 !Lambda^C (particle) needed in backtransformation of N1^B 1492 CALL GET_LAMBDAX(WORK(KLAMPC),WORK(KLAMHC),LISTC,IDLSTC,ISYMRC, 1493 * WORK(KLAMP0),WORK(KLAMH0), 1494 * WORK(KEND0),LWRK0) 1495C 1496 !Read (gei,F) and (fei,G) intermediates from files 1497 !add them and put the result to a file as (fge,I) 1498 CALL N1_RESORT(ISYMN1B,LUN1B,FNN1B,LUGEIB,FNGEIB,LUFEIB,FNFEIB, 1499 * WORK(KEND0),LWRK0,SKIPGEI) 1500C 1501 !Calculate <T3|[[H^C,T2^B],tau_ai]|HF> except VVVV contribution 1502 CALL N1N2_G(LUN1B,FNN1B, 1503 * ISYMN1B, 1504 * WORK(KN2MATB),ISYMN2B, 1505 * WORK(KCIOVVO),WORK(KCIOOVV), 1506 * WORK(KCIOOOO),ISINT1RC, 1507 * OMEGA1,ISYRES, 1508 * WORK(KINDSQNB),LENSQNB, 1509 * WORK(KEND0),LWRK0) 1510C 1511 !Calculate VVVV contribution to <T3|[[H^C,T2^B],tau_ai]|HF> 1512 IOPT = 1 !T1 one-index backtransformation 1513 CALL N1_GV4(IOPT, 1514 * LUN1B,FNN1B, 1515 * ISYMN1B, 1516 * WORK(KLAMPC),ISYMRC, 1517 * WORK(KLAMP0),1, 1518 * WORK(KLAMH0),1, 1519 * WORK(KLAMH0),1, 1520 * OMEGA1,ISYRES, 1521 * WORK(KEND0),LWRK0) 1522C 1523C REPEAT THINGS FOR N1MAT 1524C 1525C 1526 !Read (gei,F) and (fei,G) intermediates from files 1527 !add them and put the result to a file as (fge,I) 1528 CALL N1_RESORT(ISYMN1,LUN1,FNN1,LUGEI,FNGEI,LUFEI,FNFEI, 1529 * WORK(KEND0),LWRK0,SKIPGEI) 1530C 1531 !Calculate <T3|[[H^BC,T2^0],tau_ai]|HF> except VVVV contribution 1532 CALL N1N2_G(LUN1,FNN1, 1533 * ISYMN1, 1534 * WORK(KN2MAT),ISYMN2, 1535 * WORK(KBCIOVVO),WORK(KBCIOOVV), 1536 * WORK(KBCIOOOO),ISINT1BC, 1537 * OMEGA1,ISYRES, 1538 * WORK(KINDSQN),LENSQN, 1539 * WORK(KEND0),LWRK0) 1540C 1541 !Calculate VVVV contribution to <T3|[[H^BC,T2^0],tau_ai]|HF> 1542 IOPT = 1 !(double) T1 one-index backtransformation 1543 CALL N1_GV4(IOPT, 1544 * LUN1,FNN1, 1545 * ISYMN1, 1546 * WORK(KLAMPB),ISYMRB, 1547 * WORK(KLAMPC),ISYMRC, 1548 * WORK(KLAMH0),1, 1549 * WORK(KLAMH0),1, 1550 * OMEGA1,ISYRES, 1551 * WORK(KEND0),LWRK0) 1552C 1553 END IF 1554C 1555 1556C 1557C 1558c write(lupri,*)'AFTER ' 1559* write(lupri,*)'omega1 after CC3_GMATNEW,LISTL ', LISTL 1560* call PRINT_MATAI(OMEGA1,ISYRES) 1561* xnormval = ddot(NT1AM(ISYRES),OMEGA1,1,OMEGA1,1) 1562* write(lupri,*)'norm omega1 after CC3_GMATNEW ', xnormval 1563* xnormval = ddot(NT2AM(ISYRES),OMEGA2,1,OMEGA2,1) 1564* write(lupri,*)'omega2 after CC3_GMATNEW ' 1565* CALL OUTPAK(OMEGA2,NT1AMX,1,LUPRI) 1566COMMENT COMMENT 1567c write(lupri,*) 't3barx in CC3_GMATNEW' 1568c call print_pt3(work(kx3am),ISYM0,4) 1569COMMENT COMMENT 1570C 1571C 1572C----------------- 1573C Close files 1574C----------------- 1575C 1576 CALL WCLOSE2(LU3VI2,FN3VI2,'KEEP') 1577 CALL WCLOSE2(LU3FOP,FN3FOP,'KEEP') 1578 CALL WCLOSE2(LU3FOP2,FN3FOP2,'KEEP') 1579 1580C 1581C------------------------------------- 1582C Close and delete temporary files 1583C------------------------------------- 1584C 1585 IF (LVVVV) THEN 1586 CALL WCLOSE2(LU4VB,FN4VB,'DELETE') 1587 CALL WCLOSE2(LU4VC,FN4VC,'DELETE') 1588 CALL WCLOSE2(LU4VBC,FN4VBC,'DELETE') 1589 END IF 1590C 1591 IF (.NOT.LVVVV) THEN 1592 !Close files for N1MATB intermediates 1593 CALL WCLOSE2(LUFEIB,FNFEIB,'DELETE') 1594 CALL WCLOSE2(LUN1B,FNN1B,'DELETE') 1595 IF (.NOT.SKIPGEI) THEN 1596 CALL WCLOSE2(LUGEIB,FNGEIB,'DELETE') 1597 END IF 1598 !Close files for N1MATC intermediates 1599 CALL WCLOSE2(LUFEIC,FNFEIC,'DELETE') 1600 CALL WCLOSE2(LUN1C,FNN1C,'DELETE') 1601 IF (.NOT.SKIPGEI) THEN 1602 CALL WCLOSE2(LUGEIC,FNGEIC,'DELETE') 1603 END IF 1604 !Close files for N1MAT intermediates 1605 CALL WCLOSE2(LUFEI,FNFEI,'DELETE') 1606 CALL WCLOSE2(LUN1,FNN1,'DELETE') 1607 IF (.NOT.SKIPGEI) THEN 1608 CALL WCLOSE2(LUGEI,FNGEI,'DELETE') 1609 END IF 1610 END IF 1611C 1612 CALL WCLOSE2(LU3SRTR,FN3SRTR,'DELETE') 1613 CALL WCLOSE2(LUCKJDR,FNCKJDR,'DELETE') 1614 CALL WCLOSE2(LUDELDR,FNDELDR,'DELETE') 1615 CALL WCLOSE2(LUDKBCR,FNDKBCR,'DELETE') 1616 CALL WCLOSE2(LUDKBCR4,FNDKBCR4,'DELETE') 1617 1618C------------- 1619C End 1620C------------- 1621C 1622 CALL QEXIT('GNEW') 1623C 1624 RETURN 1625 END 1626 1627 1628