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 SUBROUTINE C3FOCK(IBEQC,FDTI,VECB,VECC, 20 * ZYMB,ZYMC,KZYVB,KZYVC, 21 * ISYMA,ISYMB,ISYMC, 22 * CMO,FC,MJWOP,WRK,LWRK) 23C 24C PURPOSE: 25C To calculate Fock matrix (FDTI) with doubly one-index transformed 26C integrals to be used in direct cubic RPA. 27C March 1997: tsaue Just one call to SIRFCK, big speedup ! 28C 29#include "implicit.h" 30C 31 DIMENSION FDTI(N2ORBX) 32 DIMENSION VECB(*),VECC(*) 33 DIMENSION ZYMB(NORBT,NORBT),ZYMC(NORBT,NORBT) 34 DIMENSION CMO(*),FC(*),WRK(*),ISYMDM(3) 35C 36C INFDIM : NASHDI 37C INFINP : DIRFCK 38C WRKRSP : KSYMOP 39C 40#include "maxorb.h" 41#include "maxash.h" 42#include "priunit.h" 43#include "inforb.h" 44#include "infinp.h" 45#include "infvar.h" 46 DIMENSION MJWOP(2,MAXWOP,8) 47#include "wrkrsp.h" 48#include "infrsp.h" 49C 50 CALL QENTER('C3FOCK') 51 IF (8*N2BASX .GT. LWRK) THEN 52 WRITE (LUPRI,'(A)') ' Due to lack of available memory, I '// 53 & 'use a slower C3FOCK routine' 54 CALL C3FCKM(IBEQC,FDTI,VECB,VECC, 55 * ZYMB,ZYMC,KZYVB,KZYVC, 56 * ISYMA,ISYMB,ISYMC, 57 * CMO,FC,MJWOP,WRK,LWRK) 58 ELSE 59 CALL C3FCKO(IBEQC,FDTI,VECB,VECC, 60 * ZYMB,ZYMC,KZYVB,KZYVC, 61 * ISYMA,ISYMB,ISYMC, 62 * CMO,FC,MJWOP,WRK(1),LWRK) 63 END IF 64 CALL QEXIT('C3FOCK') 65 RETURN 66 END 67C 68 SUBROUTINE C3FCKO(IBEQC,FDTI,VECB,VECC, 69 * ZYMB,ZYMC,KZYVB,KZYVC, 70 * ISYMA,ISYMB,ISYMC, 71 * CMO,FC,MJWOP,WRK,LWRK) 72C 73C PURPOSE: 74C To calculate Fock matrix (FDTI) with doubly one-index transformed 75C integrals to be used in direct cubic RPA. 76C March 1997: tsaue Just one call to SIRFCK, big speedup ! 77C 78#include "implicit.h" 79#include "dummy.h" 80#include "thrzer.h" 81C 82 DIMENSION FDTI(N2ORBX) 83 DIMENSION VECB(*),VECC(*) 84 DIMENSION ZYMB(NORBT,NORBT),ZYMC(NORBT,NORBT) 85 DIMENSION CMO(*),FC(*),WRK(*),ISYMDM(3),IFCTYP(3) 86C 87C INFDIM : NASHDI 88C INFINP : DIRFCK 89C WRKRSP : KSYMOP 90C 91#include "maxorb.h" 92#include "maxash.h" 93#include "priunit.h" 94#include "inforb.h" 95#include "infinp.h" 96#include "infvar.h" 97 DIMENSION MJWOP(2,MAXWOP,8) 98#include "wrkrsp.h" 99#include "infrsp.h" 100#include "infpri.h" 101C 102 PARAMETER ( D1 = 1.0D0, D2 = 2.0D0, D4 = 4.0D0 ) 103C 104C Allocate work space for Fock matrices 105C Corresponding Fock matrices in equations 106C 107C F-1 <--> F12 + F21 108C F-2 <--> F1 109C F-3 <--> F2 110C 111 CALL QENTER('C3FCKO') 112 KF1AO = 1 113 KF2AO = KF1AO + N2BASX 114 KF3AO = KF2AO + N2BASX 115 KD1AO = KF3AO + N2BASX 116 KD2AO = KD1AO + N2BASX 117 KD3AO = KD2AO + N2BASX 118 KFREE = KD3AO + N2BASX 119 LFREE = LWRK - KFREE 120 IF (LFREE.LT.0) CALL ERRWRK('C3FCKO',LFREE,LWRK) 121C 122C Unpack the response vectors 123C 124 CALL GTZYMT(1,VECB,KZYVB,ISYMB,ZYMB,MJWOP) 125 CALL GTZYMT(1,VECC,KZYVC,ISYMC,ZYMC,MJWOP) 126C 127C 128C Construct density matrices in AO basis 129C 130C 131 CALL DZERO(WRK(KF1AO),6*N2BASX) 132C 133C Alternative construction need to assign KTMP, KBCZYM, ISYMBC 134C Further scale D1AO with factor 2 135C 136 ISYMDM(1) = MULD2H(ISYMB,ISYMC) 137 ISYMDM(2) = ISYMB 138 ISYMDM(3) = ISYMC 139C 140 KBCZYM = KF1AO 141 KTMP = KF2AO 142C 143 CALL ZYTRA1(ISYMB,ZYMB,CMO,WRK(KTMP),WRK(KD2AO)) 144 CALL ZYTRA1(ISYMC,ZYMC,CMO,WRK(KTMP),WRK(KD3AO)) 145C 146 CALL DZERO(WRK(KTMP),N2BASX) 147 CALL ZYMUL2(ISYMB,ISYMC,ZYMB,ZYMC,WRK(KBCZYM)) 148 CALL ZYMUL2(ISYMC,ISYMB,ZYMC,ZYMB,WRK(KBCZYM)) 149 CALL ZYTRA2(ISYMDM(1),WRK(KBCZYM),CMO,WRK(KTMP),WRK(KD1AO)) 150C 151C 152C Transpose and scale by 2 to agree with convention used by SIRFCK 153C Further scale by 2 due to permutations. D2AO and D3AO scaled in one call. 154C 155 CALL DSCAL(3*N2BASX,D4,WRK(KD1AO),1) 156C 157 IF (IPRRSP.GT.200) THEN 158 WRITE(LUPRI,'(//A)')'D12 in C3FOCK' 159 CALL OUTPUT(WRK(KD1AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 160 WRITE(LUPRI,'(//A)')'D1 in C3FOCK' 161 CALL OUTPUT(WRK(KD2AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 162 WRITE(LUPRI,'(//A)')'D2 in C3FOCK' 163 CALL OUTPUT(WRK(KD3AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 164 END IF 165C 166C Construct Fock matrices in AO basis 167C 168 NFMATX = 1 169C 170C IFCTYP = XY 171C X indicates symmetry about diagonal 172C X = 0 No symmetry 173C X = 1 Symmetric 174C X = 2 Anti-symmetric 175C Y indicates contributions 176C Y = 0 no contribution ! 177C Y = 1 Coulomb 178C Y = 2 Exchange 179C Y = 3 Coulomb + Exchange 180C 181C Check if density matrix is unsymmetric (IX=0), 182C symmetric (IX=10), antisymmetric (IX=20), or zero matrix (IX=30) 183C to threshold THRZER 184C 185 JD1AO = KD1AO 186 DO I = 1,3 187 IX = 10 * MATSYM(NBAST,NBAST,WRK(JD1AO),THRZER) 188C INTEGER FUNCTION MATSYM(N,NDIM,AMAT,THRZER) 189C 190 IF (IX .EQ. 30) THEN 191C zero density matrix, do nothing ! 192 IFCTYP(I) = 0 193CCCC ELSE IF ("triplet Fock matrix") THEN 194C only exchange ! hjaaj 990505: triplet not implemented! 195CCCC IFCTYP(I) = IX + 2 196 ELSE IF (IX .EQ. 20) THEN 197C Only exchange if antisymmetric density matrix 198 IFCTYP(I) = IX + 2 199 ELSE 200C Coulomb+exchange 201 IFCTYP(I) = IX + 3 202 END IF 203 JD1AO = JD1AO + N2BASX 204 END DO 205C 206 CALL DZERO(WRK(KF1AO),3*N2BASX) 207 IF (IBEQC.EQ.1) THEN 208 NFMATX = 2 209 CALL SIRFCK(WRK(KF1AO),WRK(KD1AO),NFMATX,ISYMDM,IFCTYP, 210 * DIRFCK,WRK(KFREE),LFREE) 211 CALL DCOPY(N2BASX,WRK(KF2AO),1,WRK(KF3AO),1) 212 ELSE 213 NFMATX = 3 214 CALL SIRFCK(WRK(KF1AO),WRK(KD1AO),NFMATX,ISYMDM,IFCTYP, 215 * DIRFCK,WRK(KFREE),LFREE) 216 ENDIF 217C 218 IF (IPRRSP.GT.200) THEN 219 WRITE(LUPRI,'(//A)')'F12 in AO basis in C3FOCK' 220 CALL OUTPUT(WRK(KF1AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 221 WRITE(LUPRI,'(//A)')'F1 in AO basis in C3FOCK' 222 CALL OUTPUT(WRK(KF2AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 223 WRITE(LUPRI,'(//A)')'F2 in AO basis in C3FOCK' 224 CALL OUTPUT(WRK(KF3AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 225 END IF 226C 227C Transform Fock matrices to MO basis 228C Reuse memory used for density matrices 229C 230 KF2MO = KD1AO 231 KF3MO = KF2MO + N2ORBX 232C 233 CALL FAOMO(ISYMDM(1),WRK(KF1AO),FDTI,CMO,WRK(KFREE),LFREE) 234 CALL FAOMO(ISYMDM(2),WRK(KF2AO),WRK(KF2MO),CMO,WRK(KFREE),LFREE) 235 CALL FAOMO(ISYMDM(3),WRK(KF3AO),WRK(KF3MO),CMO,WRK(KFREE),LFREE) 236C 237 IF (IPRRSP.GT.200) THEN 238 WRITE(LUPRI,'(//A)')'F12 in MO basis in C3FOCK' 239 CALL OUTPUT(FDTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 240 WRITE(LUPRI,'(//A)')'F1 in MO basis in C3FOCK' 241 CALL OUTPUT(WRK(KF2MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 242 WRITE(LUPRI,'(//A)')'F2 in MO basis in C3FOCK' 243 CALL OUTPUT(WRK(KF3MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 244 END IF 245C 246C Evaluate FDTI see formula 247C 248C 2*F2 + [k2,Fc] --> F-3 249C 2*F1 + [k1,Fc] --> F-2 250C 251 KSYMOP = ISYMC 252 CALL FCKOIN(1,FC,DUMMY,ZYMC,WRK(KF3MO),DUMMY) 253 KSYMOP = ISYMB 254 CALL FCKOIN(1,FC,DUMMY,ZYMB,WRK(KF2MO),DUMMY) 255C 256C F12 + [k2,F-3] --> FDTI 257C FDTI + [k1,F-2] --> FDTI 258C 259 CALL OITH1(ISYMC,ZYMC,WRK(KF2MO),FDTI,ISYMDM(2)) 260 CALL OITH1(ISYMB,ZYMB,WRK(KF3MO),FDTI,ISYMDM(3)) 261C 262C Scale with 1/2 ( from definition of E3 ) 263C 264 CALL DSCAL(NORBT*NORBT,1/D2,FDTI,1) 265C 266 IF (IPRRSP.GT.100) THEN 267 WRITE(LUPRI,'(//A)')'Final result in C3FCKO' 268 CALL OUTPUT(FDTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 269 END IF 270C 271 CALL QEXIT('C3FCKO') 272 RETURN 273 END 274C 275 SUBROUTINE C3FCKM(IBEQC,FDTI,VECB,VECC, 276 * ZYMB,ZYMC,KZYVB,KZYVC, 277 * ISYMA,ISYMB,ISYMC, 278 * CMO,FC,MJWOP,WRK,LWRK) 279C 280C PURPOSE: 281C To calculate Fock matrix (FDTI) with doubly one-index transformed 282C integrals to be used in direct cubic RPA. 283C 284#include "implicit.h" 285#include "dummy.h" 286#include "thrzer.h" 287C 288 DIMENSION FDTI(N2ORBX) 289 DIMENSION VECB(*),VECC(*) 290 DIMENSION ZYMB(NORBT,NORBT),ZYMC(NORBT,NORBT) 291 DIMENSION CMO(*),FC(*),WRK(*) 292C 293C INFDIM : NASHDI 294C INFINP : DIRFCK 295C WRKRSP : KSYMOP 296C 297#include "maxorb.h" 298#include "maxash.h" 299#include "priunit.h" 300#include "inforb.h" 301#include "infinp.h" 302#include "infvar.h" 303 DIMENSION MJWOP(2,MAXWOP,8) 304#include "wrkrsp.h" 305#include "infrsp.h" 306#include "infpri.h" 307C 308 PARAMETER ( D1 = 1.0D0, D2 = 2.0D0, D4 = 4.0D0 ) 309C 310 NFMATX = 1 311C 312C Allocate work space for Fock matrices 313C Corresponding Fock matrices in equations 314C 315C F-1 <--> F12 + F21 316C F-2 <--> F1 317C F-3 <--> F2 318C 319 CALL QENTER('C3FCKM') 320 KFAO = 1 321 KDAO = KFAO + N2BASX 322 KFRE1 = KDAO + N2BASX 323 KFREE = KFRE1 + N2BASX 324 LFREE = LWRK - KFREE 325 IF (LFREE.LT.0) CALL ERRWRK('C3FCKM',LFREE,LWRK) 326C 327 KFREE = KFRE1 328 LFREE = LWRK - KFREE 329C 330C Unpack the response vectors 331C 332 CALL GTZYMT(1,VECB,KZYVB,ISYMB,ZYMB,MJWOP) 333 CALL GTZYMT(1,VECC,KZYVC,ISYMC,ZYMC,MJWOP) 334C 335C 336C Construct density matrices in AO basis 337C 338C 339 CALL DZERO(WRK(KFAO),3*N2BASX) 340 CALL DZERO(FDTI,N2ORBX) 341 CALL CDENS2(ISYMB,ISYMC,CMO,ZYMB,ZYMC, 342 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),FDTI) 343 CALL CDENS2(ISYMC,ISYMB,CMO,ZYMC,ZYMB, 344 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),FDTI) 345 CALL DSCAL(N2BASX,D2,WRK(KDAO),1) 346 IF (IPRRSP.GT.200) THEN 347 WRITE(LUPRI,'(//A)')'D12 in C3FOCK' 348 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 349 END IF 350C 351 ISYMDM = MULD2H(ISYMB,ISYMC) 352 CALL DZERO(WRK(KFAO),N2BASX) 353Chj 354C Check if density matrix is unsymmetric (IX=0), 355C symmetric (IX=10), antisymmetric (IX=20), or zero matrix (IX=30) 356C to threshold THRZER 357C 358 IX = 10 * MATSYM(NBAST,NBAST,WRK(KDAO),THRZER) 359C 360 IF (IX .EQ. 30) THEN 361C zero density matrix, do nothing ! 362 IFCTYP = 0 363 ELSE IF (IX .EQ. 20) THEN 364C Only exchange if antisymmetric density matrix 365 IFCTYP = IX + 2 366 ELSE 367C Coulomb+exchange 368 IFCTYP = IX + 3 369 END IF 370 IF (IFCTYP .NE. 0) THEN 371 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 372 * DIRFCK,WRK(KFREE),LFREE) 373 END IF 374 IF (IPRRSP.GT.200) THEN 375 WRITE(LUPRI,'(//A)')'F12 in AO basis in C3FOCK' 376 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 377 END IF 378C 379 CALL FAOMO(ISYMDM,WRK(KFAO),FDTI,CMO,WRK(KFREE),LFREE) 380 IF (IPRRSP.GT.200) THEN 381 WRITE(LUPRI,'(//A)')'F12 in MO basis in C3FOCK' 382 CALL OUTPUT(FDTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 383 END IF 384C 385C Second contribution 386C 387 CALL DZERO(WRK(KFAO),2*N2BASX) 388 CALL CDENS1(ISYMB,CMO,ZYMB,WRK(KDAO),WRK(KFREE),LFREE) 389 CALL DSCAL(N2BASX,D4,WRK(KDAO),1) 390 IF (IPRRSP .GT. 200) THEN 391 WRITE(LUPRI,'(//A)')'D1 in C3FOCK' 392 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 393 END IF 394C 395 ISYMDM = ISYMB 396Chj 397C Check if density matrix is unsymmetric (IX=0), 398C symmetric (IX=10), antisymmetric (IX=20), or zero matrix (IX=30) 399C to threshold THRZER 400C 401 IX = 10 * MATSYM(NBAST,NBAST,WRK(KDAO),THRZER) 402C 403 IF (IX .EQ. 30) THEN 404C zero density matrix, do nothing ! 405 IFCTYP = 0 406 ELSE IF (IX .EQ. 20) THEN 407C Only exchange if antisymmetric density matrix 408 IFCTYP = IX + 2 409 ELSE 410C Coulomb+exchange 411 IFCTYP = IX + 3 412 END IF 413 IF (IFCTYP .NE. 0) THEN 414 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 415 * DIRFCK,WRK(KFREE),LFREE) 416 END IF 417 IF (IPRRSP .GT. 200) THEN 418 WRITE(LUPRI,'(//A)')'F1 in AO basis in C3FOCK' 419 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 420 END IF 421C 422 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 423 IF (IPRRSP .GT. 200) THEN 424 WRITE(LUPRI,'(//A)')'F1 in MO basis in C3FOCK' 425 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 426 END IF 427C 428 KSYMOP = ISYMB 429 CALL FCKOIN(1,FC,DUMMY,ZYMB,WRK(KDAO),DUMMY) 430 IF (IBEQC .EQ. 1) THEN 431C 432C If B and C vector identical, all that remains to do is to scale 433C this term by a factor of two. 434C 435 CALL DSCAL(N2ORBX,D2,WRK(KDAO),1) 436 ELSE 437 CALL OITH1(ISYMC,ZYMC,WRK(KDAO),FDTI,ISYMDM) 438C 439C Third contribution 440C 441 CALL DZERO(WRK(KFAO),2*N2BASX) 442 CALL CDENS1(ISYMC,CMO,ZYMC,WRK(KDAO),WRK(KFREE),LFREE) 443 CALL DSCAL(N2BASX,D4,WRK(KDAO),1) 444 IF (IPRRSP .GT. 200) THEN 445 WRITE(LUPRI,'(//A)')'D2 in C3FOCK' 446 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 447 END IF 448C 449 ISYMDM = ISYMC 450Chj 451C Check if density matrix is unsymmetric (IX=0), 452C symmetric (IX=10), antisymmetric (IX=20), or zero matrix (IX=30) 453C to threshold THRZER 454C 455 IX = 10 * MATSYM(NBAST,NBAST,WRK(KDAO),THRZER) 456C 457 IF (IX .EQ. 30) THEN 458C zero density matrix, do nothing ! 459 IFCTYP = 0 460 ELSE IF (IX .EQ. 20) THEN 461C Only exchange if antisymmetric density matrix 462 IFCTYP = IX + 2 463 ELSE 464C Coulomb+exchange 465 IFCTYP = IX + 3 466 END IF 467 IF (IFCTYP .NE. 0) THEN 468 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 469 * DIRFCK,WRK(KFREE),LFREE) 470 END IF 471 IF (IPRRSP .GT. 200) THEN 472 WRITE(LUPRI,'(//A)')'F2 in AO basis in C3FOCK' 473 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 474 END IF 475C 476 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 477 IF (IPRRSP .GT. 200) THEN 478 WRITE(LUPRI,'(//A)')'F2 in MO basis in C3FOCK' 479 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 480 END IF 481C 482 KSYMOP = ISYMC 483 CALL FCKOIN(1,FC,DUMMY,ZYMC,WRK(KDAO),DUMMY) 484 CALL OITH1(ISYMB,ZYMB,WRK(KDAO),FDTI,ISYMDM) 485 END IF 486C 487C All contributions to FDTI now finished, we 488C scale with 1/2 ( from definition of E3 ) 489C 490 CALL DSCAL(NORBT*NORBT,1/D2,FDTI,1) 491C 492 IF (IPRRSP.GT.100) THEN 493 WRITE(LUPRI,'(//A)')'Final result in C3FCKM' 494 CALL OUTPUT(FDTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 495 END IF 496C 497 CALL QEXIT('C3FCKM') 498 RETURN 499 END 500C 501 SUBROUTINE C4FOCK(IBCDEQ,FTTI,VECB,VECC,VECD, 502 * ZYMB,ZYMC,ZYMD,KZYVB,KZYVC,KZYVD, 503 * ISYMA,ISYMB,ISYMC,ISYMD, 504 * CMO,FC,MJWOP,WRK,LWRK) 505C 506C PURPOSE: 507C To calculate Fock matrix (FTTI) with triply one-index transformed 508C integrals to be used in direct cubic RPA. 509C March 1997: tsaue Just one call to SIRFCK, big speedup ! 510C July 1997: kruud. Interface routine; further path depends on available 511C memory 512C 513#include "implicit.h" 514C 515 DIMENSION FTTI(N2ORBX) 516 DIMENSION VECB(*),VECC(*),VECD(*) 517 DIMENSION ZYMB(NORBT,NORBT),ZYMC(NORBT,NORBT),ZYMD(NORBT,NORBT) 518 DIMENSION CMO(*),FC(*),WRK(*) 519C 520C INFDIM : NASHDI 521C INFINP : DIRFCK 522C WRKRSP : KSYMOP 523C 524#include "maxorb.h" 525#include "maxash.h" 526#include "priunit.h" 527#include "inforb.h" 528#include "infinp.h" 529#include "infvar.h" 530 DIMENSION MJWOP(2,MAXWOP,8) 531#include "wrkrsp.h" 532#include "infrsp.h" 533C 534 IF (14*N2BASX .GT. LWRK) THEN 535 WRITE (LUPRI,'(A)') ' Due to potential lack of available '// 536 & 'memory, I use a slower C4FOCK routine' 537 CALL C4FCKM(IBCDEQ,FTTI,VECB,VECC,VECD, 538 * ZYMB,ZYMC,ZYMD,KZYVB,KZYVC,KZYVD, 539 * ISYMA,ISYMB,ISYMC,ISYMD, 540 * CMO,FC,MJWOP,WRK,LWRK) 541 ELSE 542 CALL C4FCKO(IBCDEQ,FTTI,VECB,VECC,VECD, 543 * ZYMB,ZYMC,ZYMD,KZYVB,KZYVC,KZYVD, 544 * ISYMA,ISYMB,ISYMC,ISYMD, 545 * CMO,FC,MJWOP,WRK,LWRK) 546 END IF 547C 548 RETURN 549 END 550C 551 SUBROUTINE C4FCKO(IBCDEQ,FTTI,VECB,VECC,VECD, 552 * ZYMB,ZYMC,ZYMD,KZYVB,KZYVC,KZYVD, 553 * ISYMA,ISYMB,ISYMC,ISYMD, 554 * CMO,FC,MJWOP,WRK,LWRK) 555#include "implicit.h" 556#include "dummy.h" 557C 558 DIMENSION FTTI(N2ORBX) 559 DIMENSION VECB(*),VECC(*),VECD(*) 560 DIMENSION ZYMB(NORBT,NORBT),ZYMC(NORBT,NORBT),ZYMD(NORBT,NORBT) 561 DIMENSION CMO(*),FC(*),WRK(*),ISYMDM(7),IFCTYP(7) 562C 563C INFDIM : NASHDI 564C INFINP : DIRFCK 565C WRKRSP : KSYMOP 566C 567#include "maxorb.h" 568#include "maxash.h" 569#include "priunit.h" 570#include "inforb.h" 571#include "infinp.h" 572#include "infvar.h" 573 DIMENSION MJWOP(2,MAXWOP,8) 574#include "wrkrsp.h" 575#include "infrsp.h" 576#include "infpri.h" 577C 578 PARAMETER ( D1 = 1.0D0, D2 = 2.0D0, D6 = 6.0D0 , D8 = 8.0D0 , 579 & D12 = 12.0D0) 580C Allocate work space for Fock matrices 581C Corresponding Fock matrices in equations 582C 583C F-1 <--> F123 + permutations 584C F-2 <--> F12 + F21 585C F-3 <--> F13 + F31 586C F-4 <--> F23 + F32 587C F-5 <--> F1 588C F-6 <--> F2 589C F-7 <--> F3 590C 591 KF1AO = 1 592 KF2AO = KF1AO + N2BASX 593 KF3AO = KF2AO + N2BASX 594 KF4AO = KF3AO + N2BASX 595 KF5AO = KF4AO + N2BASX 596 KF6AO = KF5AO + N2BASX 597 KF7AO = KF6AO + N2BASX 598 KD1AO = KF7AO + N2BASX 599 KD2AO = KD1AO + N2BASX 600 KD3AO = KD2AO + N2BASX 601 KD4AO = KD3AO + N2BASX 602 KD5AO = KD4AO + N2BASX 603 KD6AO = KD5AO + N2BASX 604 KD7AO = KD6AO + N2BASX 605 KFREE = KD7AO + N2BASX 606 LFREE = LWRK - KFREE 607C 608 IF (IPRRSP.GT.7) WRITE(LUPRI,'(//A,2(/A,I10))') 609 *' Memory allocation in C4FCKO', 610 *' Available workspace: ', LWRK, 611 *' Allocated : ', KFREE 612C 613C Unpack the response vectors 614C 615 CALL GTZYMT(1,VECB,KZYVB,ISYMB,ZYMB,MJWOP) 616 CALL GTZYMT(1,VECC,KZYVC,ISYMC,ZYMC,MJWOP) 617 CALL GTZYMT(1,VECD,KZYVD,ISYMD,ZYMD,MJWOP) 618C 619 CALL DZERO(WRK(KF1AO),14*N2BASX) 620C 621C Alternativ construction need to assign KTMP, K..ZYM, K...ZYM 622C 623 KBCZYM = KF1AO 624 KBCDZYM = KF2AO 625 KTMP = KF3AO 626 KBDZYM = KBCZYM 627 KCDZYM = KBCZYM 628C 629 ISYMDM(1) = ISYMA 630 ISYMDM(2) = MULD2H(ISYMB,ISYMC) 631 ISYMDM(3) = MULD2H(ISYMB,ISYMD) 632 ISYMDM(4) = MULD2H(ISYMC,ISYMD) 633 ISYMDM(5) = ISYMB 634 ISYMDM(6) = ISYMC 635 ISYMDM(7) = ISYMD 636 ISYMBCD = MULD2H(ISYMDM(2),ISYMD) 637C 638 CALL ZYTRA1(ISYMB,ZYMB,CMO,WRK(KTMP),WRK(KD5AO)) 639 CALL ZYTRA1(ISYMC,ZYMC,CMO,WRK(KTMP),WRK(KD6AO)) 640 CALL ZYTRA1(ISYMD,ZYMD,CMO,WRK(KTMP),WRK(KD7AO)) 641C 642 CALL DZERO(WRK(KBCZYM),3*N2ORBX) 643C 644 CALL ZYMUL2(ISYMB,ISYMD,ZYMB,ZYMD,WRK(KBDZYM)) 645 CALL ZYMUL2(ISYMD,ISYMB,ZYMD,ZYMB,WRK(KBDZYM)) 646 CALL ZYMUL3(ISYMDM(3),ISYMC,WRK(KBDZYM),ZYMC,WRK(KBCDZYM)) 647 CALL ZYTRA2(ISYMDM(3),WRK(KBDZYM),CMO,WRK(KTMP),WRK(KD3AO)) 648C 649 CALL DZERO(WRK(KBDZYM),N2ORBX) 650 CALL DZERO(WRK(KTMP),N2BASX) 651 CALL ZYMUL2(ISYMB,ISYMC,ZYMB,ZYMC,WRK(KBCZYM)) 652 CALL ZYMUL2(ISYMC,ISYMB,ZYMC,ZYMB,WRK(KBCZYM)) 653 CALL ZYMUL3(ISYMDM(2),ISYMD,WRK(KBCZYM),ZYMD,WRK(KBCDZYM)) 654 CALL ZYTRA2(ISYMDM(2),WRK(KBCZYM),CMO,WRK(KTMP),WRK(KD2AO)) 655C 656 CALL DZERO(WRK(KCDZYM),N2ORBX) 657 CALL DZERO(WRK(KCDZYM),N2ORBX) 658 CALL ZYMUL2(ISYMC,ISYMD,ZYMC,ZYMD,WRK(KCDZYM)) 659 CALL ZYMUL2(ISYMD,ISYMC,ZYMD,ZYMC,WRK(KCDZYM)) 660 CALL ZYMUL3(ISYMDM(4),ISYMB,WRK(KCDZYM),ZYMB,WRK(KBCDZYM)) 661 CALL ZYTRA2(ISYMDM(4),WRK(KCDZYM),CMO,WRK(KTMP),WRK(KD4AO)) 662C 663 CALL DZERO(WRK(KTMP),N2BASX) 664 CALL ZYTRA1(ISYMBCD,WRK(KBCDZYM),CMO,WRK(KTMP),WRK(KD1AO)) 665C 666C Transpose and scale by 2 to agree with convention used by SIRFCK 667C Further scale by 3 due to permutations. D2AO -> D7AO scaled in one call 668C 669 CALL DSCAL(N2BASX,D8,WRK(KD1AO),1) 670 CALL DSCAL(3*N2BASX,D12,WRK(KD2AO),1) 671 CALL DSCAL(3*N2BASX,D6,WRK(KD5AO),1) 672C 673 IF (IPRRSP.GT.200) THEN 674 WRITE(LUPRI,'(//A)')'D123 in C4FOCK' 675 CALL OUTPUT(WRK(KD1AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 676 WRITE(LUPRI,'(//A)')'D12 in C4FOCK' 677 CALL OUTPUT(WRK(KD2AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 678 WRITE(LUPRI,'(//A)')'D13 in C4FOCK' 679 CALL OUTPUT(WRK(KD3AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 680 WRITE(LUPRI,'(//A)')'D23 in C4FOCK' 681 CALL OUTPUT(WRK(KD4AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 682 WRITE(LUPRI,'(//A)')'D1 in C4FOCK' 683 CALL OUTPUT(WRK(KD5AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 684 WRITE(LUPRI,'(//A)')'D2 in C4FOCK' 685 CALL OUTPUT(WRK(KD6AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 686 WRITE(LUPRI,'(//A)')'D3 in C4FOCK' 687 CALL OUTPUT(WRK(KD7AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 688 END IF 689C 690C Construct Fock matrices in AO basis 691C 692 DO I = 1,7 693 IFCTYP(I) = 3 694 ENDDO 695 NFMATX = 7 696 IF (IBCDEQ.EQ.24) THEN 697 NFMATX = 5 698 IFCTYP(3) = 0 699 IFCTYP(4) = 0 700 IFCTYP(6) = 0 701 IFCTYP(7) = 0 702 ELSE IF (IBCDEQ.EQ.2) THEN 703 IFCTYP(4) = 0 704 IFCTYP(6) = 0 705 ELSE IF (IBCDEQ.EQ.3) THEN 706 NFMATX = 6 707 IFCTYP(4) = 0 708 IFCTYP(7) = 0 709 ELSE IF (IBCDEQ.EQ.4) THEN 710 NFMATX = 6 711 IFCTYP(3) = 0 712 IFCTYP(7) = 0 713 ENDIF 714C 715 CALL DZERO(WRK(KF1AO),7*N2BASX) 716 CALL SIRFCK(WRK(KF1AO),WRK(KD1AO),NFMATX,ISYMDM,IFCTYP, 717 * DIRFCK,WRK(KFREE),LFREE) 718C 719 IF (IBCDEQ.EQ.24) THEN 720 CALL DCOPY(N2BASX,WRK(KF2AO),1,WRK(KF3AO),1) 721 CALL DCOPY(N2BASX,WRK(KF2AO),1,WRK(KF4AO),1) 722 CALL DCOPY(N2BASX,WRK(KF5AO),1,WRK(KF6AO),1) 723 CALL DCOPY(N2BASX,WRK(KF5AO),1,WRK(KF7AO),1) 724 ELSE IF (IBCDEQ.EQ.2) THEN 725 CALL DCOPY(N2BASX,WRK(KF3AO),1,WRK(KF4AO),1) 726 CALL DCOPY(N2BASX,WRK(KF5AO),1,WRK(KF6AO),1) 727 ELSE IF (IBCDEQ.EQ.3) THEN 728 CALL DCOPY(N2BASX,WRK(KF2AO),1,WRK(KF4AO),1) 729 CALL DCOPY(N2BASX,WRK(KF5AO),1,WRK(KF7AO),1) 730 ELSE IF (IBCDEQ.EQ.4) THEN 731 CALL DCOPY(N2BASX,WRK(KF2AO),1,WRK(KF3AO),1) 732 CALL DCOPY(N2BASX,WRK(KF6AO),1,WRK(KF7AO),1) 733 END IF 734C 735 IF (IPRRSP.GT.200) THEN 736 WRITE(LUPRI,'(//A)')'F123 in AO basis in C4FOCK' 737 CALL OUTPUT(WRK(KF1AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 738 WRITE(LUPRI,'(//A)')'F12 in AO basis in C4FOCK' 739 CALL OUTPUT(WRK(KF2AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 740 WRITE(LUPRI,'(//A)')'F13 in AO basis in C4FOCK' 741 CALL OUTPUT(WRK(KF3AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 742 WRITE(LUPRI,'(//A)')'F23 in AO basis in C4FOCK' 743 CALL OUTPUT(WRK(KF4AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 744 WRITE(LUPRI,'(//A)')'F1 in AO basis in C4FOCK' 745 CALL OUTPUT(WRK(KF5AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 746 WRITE(LUPRI,'(//A)')'F2 in AO basis in C4FOCK' 747 CALL OUTPUT(WRK(KF6AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 748 WRITE(LUPRI,'(//A)')'F3 in AO basis in C4FOCK' 749 CALL OUTPUT(WRK(KF7AO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 750 END IF 751C 752C Transform Fock matrices to MO basis 753C Reuse memory used for density matrices 754C 755 KF2MO = KD1AO 756 KF3MO = KF2MO + N2ORBX 757 KF4MO = KF3MO + N2ORBX 758 KF5MO = KF4MO + N2ORBX 759 KF6MO = KF5MO + N2ORBX 760 KF7MO = KF6MO + N2ORBX 761C 762 CALL FAOMO(ISYMDM(1),WRK(KF1AO),FTTI,CMO,WRK(KFREE),LFREE) 763 CALL FAOMO(ISYMDM(2),WRK(KF2AO),WRK(KF2MO),CMO,WRK(KFREE),LFREE) 764 CALL FAOMO(ISYMDM(3),WRK(KF3AO),WRK(KF3MO),CMO,WRK(KFREE),LFREE) 765 CALL FAOMO(ISYMDM(4),WRK(KF4AO),WRK(KF4MO),CMO,WRK(KFREE),LFREE) 766 CALL FAOMO(ISYMDM(5),WRK(KF5AO),WRK(KF5MO),CMO,WRK(KFREE),LFREE) 767 CALL FAOMO(ISYMDM(6),WRK(KF6AO),WRK(KF6MO),CMO,WRK(KFREE),LFREE) 768 CALL FAOMO(ISYMDM(7),WRK(KF7AO),WRK(KF7MO),CMO,WRK(KFREE),LFREE) 769C 770 IF (IPRRSP.GT.200) THEN 771 WRITE(LUPRI,'(//A)')'F123 in MO basis in C4FOCK' 772 CALL OUTPUT(FTTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 773 WRITE(LUPRI,'(//A)')'F12 in MO basis in C4FOCK' 774 CALL OUTPUT(WRK(KF2MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 775 WRITE(LUPRI,'(//A)')'F13 in MO basis in C4FOCK' 776 CALL OUTPUT(WRK(KF3MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 777 WRITE(LUPRI,'(//A)')'F23 in MO basis in C4FOCK' 778 CALL OUTPUT(WRK(KF4MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 779 WRITE(LUPRI,'(//A)')'F1 in MO basis in C4FOCK' 780 CALL OUTPUT(WRK(KF5MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 781 WRITE(LUPRI,'(//A)')'F2 in MO basis in C4FOCK' 782 CALL OUTPUT(WRK(KF6MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 783 WRITE(LUPRI,'(//A)')'F3 in MO basis in C4FOCK' 784 CALL OUTPUT(WRK(KF7MO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 785 END IF 786C 787C Evaluate FTTI see formula 788C 789C 3*F3 + [k3,Fc] --> F-7 790C 3*F2 + [k2,Fc] --> F-6 791C 3*F1 + [k1,Fc] --> F-5 792C 793 KSYMOP = ISYMD 794 CALL FCKOIN(1,FC,DUMMY,ZYMD,WRK(KF7MO),DUMMY) 795 KSYMOP = ISYMC 796 CALL FCKOIN(1,FC,DUMMY,ZYMC,WRK(KF6MO),DUMMY) 797 KSYMOP = ISYMB 798 CALL FCKOIN(1,FC,DUMMY,ZYMB,WRK(KF5MO),DUMMY) 799C 800C F23 + [k2,F-7] --> F-4 801C F-4 + [k3,F-6] --> F-4 802C 803 CALL OITH1(ISYMC,ZYMC,WRK(KF7MO),WRK(KF4MO),ISYMDM(7)) 804 CALL OITH1(ISYMD,ZYMD,WRK(KF6MO),WRK(KF4MO),ISYMDM(6)) 805C 806C F13 + [k1,F-7] --> F-3 807C F-3 + [k3,F-5] --> F-3 808C 809 CALL OITH1(ISYMB,ZYMB,WRK(KF7MO),WRK(KF3MO),ISYMDM(7)) 810 CALL OITH1(ISYMD,ZYMD,WRK(KF5MO),WRK(KF3MO),ISYMDM(5)) 811C 812C F12 + [k2,F-5] --> F-2 813C F-2 + [k1,F-6] --> F-2 814C 815 CALL OITH1(ISYMC,ZYMC,WRK(KF5MO),WRK(KF2MO),ISYMDM(5)) 816 CALL OITH1(ISYMB,ZYMB,WRK(KF6MO),WRK(KF2MO),ISYMDM(6)) 817C 818C F123 + [k1,F-4] --> FTTI 819C FTTI + [k2,F-3] --> FTTI 820C FTTI + [k3,F-2] --> FTTI 821C 822 CALL OITH1(ISYMB,ZYMB,WRK(KF4MO),FTTI,ISYMDM(4)) 823 CALL OITH1(ISYMC,ZYMC,WRK(KF3MO),FTTI,ISYMDM(3)) 824 CALL OITH1(ISYMD,ZYMD,WRK(KF2MO),FTTI,ISYMDM(2)) 825C 826C Scale with 1/6 ( from definition of E4, without minus sign ) 827C 828 CALL DSCAL(NORBT*NORBT,1/D6,FTTI,1) 829C 830 IF (IPRRSP.GT.100) THEN 831 WRITE(LUPRI,'(//A)')'Final result in C4FOCK' 832 CALL OUTPUT(FTTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 833 END IF 834C 835 RETURN 836 END 837C 838 SUBROUTINE C4FCKM(IBCDEQ,FTTI,VECB,VECC,VECD, 839 * ZYMB,ZYMC,ZYMD,KZYVB,KZYVC,KZYVD, 840 * ISYMA,ISYMB,ISYMC,ISYMD, 841 * CMO,FC,MJWOP,WRK,LWRK) 842C 843C PURPOSE: 844C To calculate Fock matrix (FTTI) with triply one-index transformed 845C integrals to be used in direct cubic RPA. 846C July 1997: kruud: Memory efficient version of C4FOCK 847C 848#include "implicit.h" 849#include "dummy.h" 850C 851 DIMENSION FTTI(N2ORBX) 852 DIMENSION VECB(*),VECC(*),VECD(*) 853 DIMENSION ZYMB(NORBT,NORBT),ZYMC(NORBT,NORBT),ZYMD(NORBT,NORBT) 854 DIMENSION CMO(*),FC(*),WRK(*) 855C 856C INFDIM : NASHDI 857C INFINP : DIRFCK 858C WRKRSP : KSYMOP 859C 860#include "maxorb.h" 861#include "maxash.h" 862#include "priunit.h" 863#include "inforb.h" 864#include "infinp.h" 865#include "infvar.h" 866 DIMENSION MJWOP(2,MAXWOP,8) 867#include "wrkrsp.h" 868#include "infrsp.h" 869#include "infpri.h" 870C 871 PARAMETER ( D1 = 1.0D0, D2 = 2.0D0, D6 = 6.0D0, D3 = 3.0D0) 872C 873C Allocate work space for Fock matrices 874C Corresponding Fock matrices in equations 875C 876C F-1 <--> F123 + permutations 877C F-2 <--> F12 + F21 878C F-3 <--> F13 + F31 879C F-4 <--> F23 + F32 880C F-5 <--> F1 881C F-6 <--> F2 882C F-7 <--> F3 883C 884 IFCTYP = 3 885 NFMATX = 1 886C 887 KFAO = 1 888 KDAO = KFAO + N2BASX 889 KFRE1 = KDAO + N2BASX 890 KFRE2 = KFRE1 + N2BASX 891 KFREE = KFRE2 + N2BASX 892 LFREE = LWRK - KFREE 893 IF (LFREE.LT.0) CALL ERRWRK('C4FCKM',LFREE,LWRK) 894C 895 IF (IPRRSP.GT.7) WRITE(LUPRI,'(//A,2(/A,I10))') 896 *' Memory allocation in C4FCKM', 897 *' Available workspace: ', LWRK, 898 *' Allocated : ', KFREE 899C 900 KFREE = KFRE1 901 LFREE = LWRK - KFREE 902C 903C Unpack the response vectors 904C 905 CALL GTZYMT(1,VECB,KZYVB,ISYMB,ZYMB,MJWOP) 906 CALL GTZYMT(1,VECC,KZYVC,ISYMC,ZYMC,MJWOP) 907 CALL GTZYMT(1,VECD,KZYVD,ISYMD,ZYMD,MJWOP) 908C 909C F123 contribution 910C 911 CALL DZERO(WRK(KFAO),4*N2BASX) 912 CALL CDENS3(ISYMB,ISYMC,ISYMD,CMO,ZYMB,ZYMC,ZYMD, 913 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 914 CALL CDENS3(ISYMB,ISYMD,ISYMC,CMO,ZYMB,ZYMD,ZYMC, 915 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 916 CALL CDENS3(ISYMC,ISYMB,ISYMD,CMO,ZYMC,ZYMB,ZYMD, 917 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 918 CALL CDENS3(ISYMC,ISYMD,ISYMB,CMO,ZYMC,ZYMD,ZYMB, 919 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 920 CALL CDENS3(ISYMD,ISYMB,ISYMC,CMO,ZYMD,ZYMB,ZYMC, 921 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 922 CALL CDENS3(ISYMD,ISYMC,ISYMB,CMO,ZYMD,ZYMC,ZYMB, 923 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 924C 925 CALL DSCAL(N2BASX,D2,WRK(KDAO),1) 926 IF (IPRRSP.GT.200) THEN 927 WRITE(LUPRI,'(//A)')'D123 in C4FOCK' 928 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 929 END IF 930C 931 ISYMDM = ISYMA 932 CALL DZERO(WRK(KFAO),N2BASX) 933 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 934 * DIRFCK,WRK(KFREE),LFREE) 935 IF (IPRRSP.GT.200) THEN 936 WRITE(LUPRI,'(//A)')'F123 in AO basis in C4FOCK' 937 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 938 END IF 939C 940 CALL FAOMO(ISYMDM,WRK(KFAO),FTTI,CMO,WRK(KFREE),LFREE) 941 IF (IPRRSP.GT.200) THEN 942 WRITE(LUPRI,'(//A)')'F123 in MO basis in C4FOCK' 943 CALL OUTPUT(FTTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 944 END IF 945C 946C F12 matrix 947C 948 CALL DZERO(WRK(KFAO),4*N2BASX) 949 CALL CDENS2(ISYMB,ISYMC,CMO,ZYMB,ZYMC, 950 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 951 CALL CDENS2(ISYMC,ISYMB,CMO,ZYMC,ZYMB, 952 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 953 CALL DSCAL(N2BASX,D6,WRK(KDAO),1) 954 IF (IPRRSP .GT. 200) THEN 955 WRITE(LUPRI,'(//A)')'D12 in C4FOCK' 956 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 957 END IF 958C 959 ISYMDM = MULD2H(ISYMB,ISYMC) 960 CALL DZERO(WRK(KFAO),N2BASX) 961 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 962 * DIRFCK,WRK(KFREE),LFREE) 963 IF (IPRRSP .GT. 200) THEN 964 WRITE(LUPRI,'(//A)')'F12 in AO basis in C4FOCK' 965 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 966 END IF 967C 968 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 969 IF (IPRRSP .GT. 200) THEN 970 WRITE(LUPRI,'(//A)')'F12 in MO basis in C4FOCK' 971 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 972 END IF 973 IF (IBCDEQ .EQ. 24) THEN 974 CALL DSCAL(N2BASX,D3,WRK(KDAO),1) 975 END IF 976 CALL OITH1(ISYMD,ZYMD,WRK(KDAO),FTTI,ISYMDM) 977C 978C F13 matrix 979C 980 IF (IBCDEQ .NE. 24) THEN 981 CALL DZERO(WRK(KFAO),4*N2BASX) 982 CALL CDENS2(ISYMB,ISYMD,CMO,ZYMB,ZYMD, 983 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 984 CALL CDENS2(ISYMD,ISYMB,CMO,ZYMD,ZYMB, 985 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 986 CALL DSCAL(N2BASX,D6,WRK(KDAO),1) 987 IF (IPRRSP .GT. 200) THEN 988 WRITE(LUPRI,'(//A)')'D13 in C4FOCK' 989 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 990 END IF 991C 992 ISYMDM = MULD2H(ISYMB,ISYMD) 993 CALL DZERO(WRK(KFAO),N2BASX) 994 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 995 * DIRFCK,WRK(KFREE),LFREE) 996 IF (IPRRSP .GT. 200) THEN 997 WRITE(LUPRI,'(//A)')'F13 in AO basis in C4FOCK' 998 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 999 END IF 1000C 1001 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 1002 IF (IPRRSP .GT. 200) THEN 1003 WRITE(LUPRI,'(//A)')'F13 in MO basis in C4FOCK' 1004 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1005 END IF 1006 CALL OITH1(ISYMC,ZYMC,WRK(KDAO),FTTI,ISYMDM) 1007C 1008C F23 contribution 1009C 1010 CALL DZERO(WRK(KFAO),4*N2BASX) 1011 CALL CDENS2(ISYMC,ISYMD,CMO,ZYMC,ZYMD, 1012 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 1013 CALL CDENS2(ISYMD,ISYMC,CMO,ZYMD,ZYMC, 1014 * WRK(KDAO),WRK(KFAO),WRK(KFRE1),WRK(KFRE2)) 1015 CALL DSCAL(N2BASX,D6,WRK(KDAO),1) 1016 IF (IPRRSP .GT. 20) THEN 1017 WRITE(LUPRI,'(//A)')'D23 in C4FOCK' 1018 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1019 END IF 1020C 1021 ISYMDM = MULD2H(ISYMC,ISYMD) 1022 CALL DZERO(WRK(KFAO),N2BASX) 1023 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 1024 * DIRFCK,WRK(KFREE),LFREE) 1025 IF (IPRRSP .GT. 200) THEN 1026 WRITE(LUPRI,'(//A)')'F23 in AO basis in C4FOCK' 1027 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1028 END IF 1029C 1030 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 1031 IF (IPRRSP .GT. 200) THEN 1032 WRITE(LUPRI,'(//A)')'F23 in MO basis in C4FOCK' 1033 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1034 END IF 1035 CALL OITH1(ISYMB,ZYMB,WRK(KDAO),FTTI,ISYMDM) 1036 END IF 1037C 1038C F1 contribution 1039C 1040 CALL DZERO(WRK(KFAO),2*N2BASX) 1041 CALL CDENS1(ISYMB,CMO,ZYMB,WRK(KDAO), 1042 * WRK(KFREE),LFREE) 1043 CALL DSCAL(N2BASX,D6,WRK(KDAO),1) 1044 IF (IPRRSP .GT. 200) THEN 1045 WRITE(LUPRI,'(//A)')'D1 in C4FOCK' 1046 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1047 END IF 1048C 1049 ISYMDM = ISYMB 1050 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 1051 * DIRFCK,WRK(KFREE),LFREE) 1052 IF (IPRRSP .GT. 200) THEN 1053 WRITE(LUPRI,'(//A)')'F1 in AO basis in C4FOCK' 1054 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1055 END IF 1056C 1057 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 1058 IF (IPRRSP .GT. 200) THEN 1059 WRITE(LUPRI,'(//A)')'F1 in MO basis in C4FOCK' 1060 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1061 END IF 1062C 1063 KSYMOP = ISYMB 1064 CALL FCKOIN(1,FC,DUMMY,ZYMB,WRK(KDAO),DUMMY) 1065 CALL DZERO(WRK(KFAO),N2BASX) 1066 CALL OITH1(ISYMD,ZYMD,WRK(KDAO),WRK(KFAO),ISYMDM) 1067 IF (IBCDEQ .EQ. 24) THEN 1068 CALL DSCAL(N2ORBX,D3,WRK(KFAO),1) 1069 END IF 1070 CALL OITH1(ISYMC,ZYMC,WRK(KFAO),FTTI,MULD2H(ISYMB,ISYMD)) 1071C 1072 CALL DZERO(WRK(KFAO),N2BASX) 1073 CALL OITH1(ISYMC,ZYMC,WRK(KDAO),WRK(KFAO),ISYMDM) 1074 IF (IBCDEQ .EQ. 24) THEN 1075 CALL DSCAL(N2ORBX,D3,WRK(KFAO),1) 1076 END IF 1077 CALL OITH1(ISYMD,ZYMD,WRK(KFAO),FTTI,MULD2H(ISYMB,ISYMC)) 1078C 1079C F2 contribution 1080C 1081 IF (IBCDEQ .NE. 24) THEN 1082 CALL DZERO(WRK(KFAO),2*N2BASX) 1083 CALL CDENS1(ISYMC,CMO,ZYMC,WRK(KDAO), 1084 * WRK(KFREE),LFREE) 1085 CALL DSCAL(N2BASX,D6,WRK(KDAO),1) 1086 IF (IPRRSP .GT. 200) THEN 1087 WRITE(LUPRI,'(//A)')'D2 in C4FOCK' 1088 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1089 END IF 1090C 1091 ISYMDM = ISYMC 1092 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 1093 * DIRFCK,WRK(KFREE),LFREE) 1094 IF (IPRRSP .GT. 200) THEN 1095 WRITE(LUPRI,'(//A)')'F2 in AO basis in C4FOCK' 1096 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1097 END IF 1098C 1099 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 1100 IF (IPRRSP .GT. 200) THEN 1101 WRITE(LUPRI,'(//A)')'F2 in MO basis in C4FOCK' 1102 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1103 END IF 1104C 1105 KSYMOP = ISYMC 1106 CALL FCKOIN(1,FC,DUMMY,ZYMC,WRK(KDAO),DUMMY) 1107 CALL DZERO(WRK(KFAO),N2BASX) 1108 CALL OITH1(ISYMD,ZYMD,WRK(KDAO),WRK(KFAO),ISYMDM) 1109 CALL OITH1(ISYMB,ZYMB,WRK(KFAO),FTTI,MULD2H(ISYMC,ISYMD)) 1110C 1111 CALL DZERO(WRK(KFAO),N2BASX) 1112 CALL OITH1(ISYMB,ZYMB,WRK(KDAO),WRK(KFAO),ISYMDM) 1113 CALL OITH1(ISYMD,ZYMD,WRK(KFAO),FTTI,MULD2H(ISYMB,ISYMC)) 1114C 1115C F3 contribution 1116C 1117 CALL CDENS1(ISYMD,CMO,ZYMD,WRK(KDAO),WRK(KFREE),LFREE) 1118 CALL DSCAL(N2BASX,D6,WRK(KDAO),1) 1119 IF (IPRRSP .GT. 200) THEN 1120 WRITE(LUPRI,'(//A)')'D3 in C4FOCK' 1121 CALL OUTPUT(WRK(KDAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1122 END IF 1123C 1124 ISYMDM = ISYMD 1125 CALL SIRFCK(WRK(KFAO),WRK(KDAO),NFMATX,ISYMDM,IFCTYP, 1126 * DIRFCK,WRK(KFREE),LFREE) 1127 IF (IPRRSP .GT. 200) THEN 1128 WRITE(LUPRI,'(//A)')'F3 in AO basis in C4FOCK' 1129 CALL OUTPUT(WRK(KFAO),1,NBAST,1,NBAST,NBAST,NBAST,1,LUPRI) 1130 END IF 1131C 1132 CALL FAOMO(ISYMDM,WRK(KFAO),WRK(KDAO),CMO,WRK(KFREE),LFREE) 1133 IF (IPRRSP .GT. 200) THEN 1134 WRITE(LUPRI,'(//A)')'F3 in MO basis in C4FOCK' 1135 CALL OUTPUT(WRK(KDAO),1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1136 END IF 1137C 1138 KSYMOP = ISYMD 1139 CALL FCKOIN(1,FC,DUMMY,ZYMD,WRK(KDAO),DUMMY) 1140 CALL DZERO(WRK(KFAO),N2BASX) 1141 CALL OITH1(ISYMC,ZYMC,WRK(KDAO),WRK(KFAO),ISYMDM) 1142 CALL OITH1(ISYMB,ZYMB,WRK(KFAO),FTTI,MULD2H(ISYMC,ISYMD)) 1143C 1144 CALL DZERO(WRK(KFAO),N2BASX) 1145 CALL OITH1(ISYMB,ZYMB,WRK(KDAO),WRK(KFAO),ISYMDM) 1146 CALL OITH1(ISYMC,ZYMC,WRK(KFAO),FTTI,MULD2H(ISYMB,ISYMD)) 1147 END IF 1148C 1149C All contributions now gathered in FTTI, and we finish by 1150C scaling with 1/6 ( from definition of E4, without minus sign ) 1151C 1152 CALL DSCAL(NORBT*NORBT,1/D6,FTTI,1) 1153C 1154 IF (IPRRSP.GT.100) THEN 1155 WRITE(LUPRI,'(//A)')'Final result in C4FOCK' 1156 CALL OUTPUT(FTTI,1,NORBT,1,NORBT,NORBT,NORBT,1,LUPRI) 1157 END IF 1158C 1159 RETURN 1160 END 1161 SUBROUTINE FAOMO(ISYMF,FAO,FMO,CMO,WRK,LWRK) 1162C 1163#include "implicit.h" 1164C 1165 DIMENSION FAO(N2BASX),FMO(N2ORBX),CMO(*),WRK(*) 1166C 1167#include "maxorb.h" 1168#include "maxash.h" 1169#include "inforb.h" 1170#include "infinp.h" 1171#include "wrkrsp.h" 1172C 1173 CALL DZERO(FMO,N2ORBX) 1174 DO 100 ISYM=1,NSYM 1175 JSYM = MULD2H(ISYM,ISYMF) 1176 NORBI = NORB(ISYM) 1177 NORBJ = NORB(JSYM) 1178 IF (NORBI.EQ.0 .OR. NORBJ.EQ.0) GO TO 100 1179 CALL AUTPV(ISYM,JSYM,CMO(ICMO(ISYM)+1),CMO(ICMO(JSYM)+1), 1180 * FAO,NBAS,NBAST,FMO,NORB,NORBT,WRK,LWRK) 1181 100 CONTINUE 1182 RETURN 1183 END 1184