1#define FCKTRA_DEBUG -1 2! 3! Dalton, a molecular electronic structure program 4! Copyright (C) by the authors of Dalton. 5! 6! This program is free software; you can redistribute it and/or 7! modify it under the terms of the GNU Lesser General Public 8! License version 2.1 as published by the Free Software Foundation. 9! 10! This program is distributed in the hope that it will be useful, 11! but WITHOUT ANY WARRANTY; without even the implied warranty of 12! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 13! Lesser General Public License for more details. 14! 15! If a copy of the GNU LGPL v2.1 was not distributed with this 16! code, you can obtain one at https://www.gnu.org/licenses/old-licenses/lgpl-2.1.en.html. 17! 18! 19 20! FILE: DALTON/sirius/sirfcktra.F 21! (c) Copyright Hans Joergen Aa. Jensen, hjj@sdu.dk (2016) 22 23! =================================================================== 24! sirfcktra.F section 0: Sirius interface 25! =================================================================== 26 subroutine SIR_INTOPEN_FCKTRA() 27! 28! Written by Hans Joergen Aa. Jensen October 2014 29! 30 implicit none 31#include "priunit.h" 32 33#if FCKTRA_DEBUG > 0 34 write(LUPRI,*) 'SIR_INTOPEN_FCKTRA called' 35#endif 36 37 end 38 39 subroutine SIR_FCKTRA_CTL(TYPE_TEXT, ITRSIR, THR_MOTWOINT, CMO, 40 & WORK, LWORK) 41! 42! Written by Hans Joergen Aa. Jensen Apr 2017 43! 44 implicit none 45 integer ITRSIR, LWORK 46 real*8 THR_MOTWOINT, CMO(*), WORK(LWORK) 47#include "priunit.h" 48 49! Used from include files: 50! gnrinf.h : PARCAL 51! inforb.h : NNORBX, NORBT, ??? 52! inftra.h : IPRTRA, LBUF 53! inftap.h : LUINTM, LBINTM 54! cbihr2.h : IFTHRS 55! ccom.h : THRS 56#include "iratdef.h" 57#include "maxaqn.h" 58#include "thrzer.h" 59#include "gnrinf.h" 60#include "inforb.h" 61#include "inftra.h" 62#include "inftap.h" 63#include "cbihr2.h" 64#include "ccom.h" 65 66 character*4 TYPE_TEXT 67 integer JTRSIR, NTRLVL, IOSVAL, MX_NDMAT 68 integer KFREE, LFREE, KFRSAV, KUCMO, KDMAT, KFMAT, I, J, K 69 integer MISH_TEST(8), MASH_TEST(8), IDAMAX, L_H2CD 70 integer, allocatable :: ICDTRA(:), ITRTYP(:) 71 character*8 TABLE(3), LAB123(3), INT_LABEL 72 data TABLE /'MUL_2INT','DRC_2INT','END_2INT'/ 73 data LAB123/'********','********','********'/ 74 logical FEXIST, OLDDX, FNDLAB 75 76C --- start of executable code 77 78 call QENTER('SIR_FCKTRA_CTL') 79 80 JTRSIR = abs(ITRSIR) 81 82#if FCKTRA_DEBUG > 0 83 IPRTRA = MAX(IPRTRA,FCKTRA_DEBUG) 84 write(lupri,*) 'DEBUG: IPRTRA set to',IPRTRA 85#endif 86 87 IF (IPRTRA .GT. 0) THEN 88 CALL TITLER('.FCKTRA AO->MO integral transformation','*',200) 89 WRITE(LUPRI,'(A,I4,1P,D10.2)') 90 & ' Sirius integral transformation level & screening threshold:', 91 & JTRSIR,THR_MOTWOINT 92 FLUSH(LUPRI) 93 END IF 94 95 ! call SIR_INTOPEN_FCKTRA() 96 97C Test if MO integrals already available 98 99 if (TYPE_TEXT .EQ. 'COUL') then 100 INT_LABEL = TABLE(1) 101 L_H2CD = NNORBT 102 else if (TYPE_TEXT .EQ. 'EXCH') then 103 INT_LABEL = TABLE(2) 104 L_H2CD = N2ORBT 105 else 106 call QUIT('Unknown TYPE_TEXT : '//TYPE_TEXT) 107 end if 108 LBUF = IRAT*L_H2CD 109 110 REWIND(LUINTM) 111 IF ( FNDLAB(INT_LABEL,LUINTM) ) THEN ! true iff existing MOTWOINT is created by .FCKTRA 112 REWIND(LUINTM) 113 READ (LUINTM) 114 READ (LUINTM) J,NTRLVL ! originally LBINTM, JTRLVL 115 IF (NTRLVL .EQ. 3 .AND. JTRSIR .EQ. 2) NTRLVL = -1 ! (aa|ii) and (ai|ai) missing 116 IF (NTRLVL .EQ. 2 .AND. JTRSIR .EQ. 3) NTRLVL = 3 117 IF (NTRLVL .EQ. 5 .AND. JTRSIR .EQ. 6) NTRLVL = 6 118 IF (NTRLVL .EQ. 6 .AND. JTRSIR .EQ. 4) NTRLVL = -1 ! (gg|ii) and (gi|gi) missing 119 IF (NTRLVL .EQ. 6 .AND. JTRSIR .EQ. 5) NTRLVL = -1 ! (gg|gi) missing 120 IF (NTRLVL .GE. JTRSIR) THEN 121C level OK, now read CMO matrix from LUINTM 122C and subtract from input CMO matrix 123 READ (LUINTM, IOSTAT=IOSVAL) WORK(1:NCMOT) 124 IF (IOSVAL .EQ. 0) THEN 125 WORK(1:NCMOT) = WORK(1:NCMOT) - CMO(1:NCMOT) 126 I = IDAMAX(NCMOT,WORK,1) 127 READ (LUINTM) MISH_TEST(1:8), MASH_TEST(1:8) 128 MISH_TEST(1:8) = MISH_TEST(1:8) - NISH(1:8) 129 MASH_TEST(1:8) = MASH_TEST(1:8) - NASH(1:8) 130 J = 0 131 DO K = 1,8 132 J = J + ABS(MISH_TEST(K)) + ABS(MASH_TEST(K)) 133 END DO 134 IF (ABS(WORK(I)) .LE. THRZER .AND. J .EQ. 0) THEN 135 IF (IPRTRA.GE.0) WRITE (LUPRI,'(/A)') 136 & ' SIR_FCKTRA_CTL: Integral transformation abandoned,' 137 & //' the required MO integrals are already available.' 138 GO TO 9999 139 END IF 140 END IF ! IOSVAL test 141 END IF ! NTRLVL test 142 END IF 143 144 IF (THR_MOTWOINT .LT. THRS) THEN 145 WRITE(LUPRI,'(/A,1P,D10.2,A,D10.2)') 146 & '.FCKTRA INFO: AO screening threshold raised from', 147 & THR_MOTWOINT,' to the general integral threshold',THRS 148 THR_MOTWOINT = THRS 149 FLUSH(LUPRI) 150 END IF 151 152 allocate(ICDTRA(NNORBX)) 153 154 KFREE = 1 155 LFREE = LWORK 156 KFRSAV = KFREE 157 158! translate sirius integral level code to FCKTRA and NEWTRA integral level code 159 call N_TRALVL(JTRSIR, NTRLVL) 160 IF (JTRSIR .GE. 1 .AND. JTRSIR .LE. 4) NTRLVL = NTRLVL + 1 161 ! we need third order, because USEDRC not defined for FCKTRA 162#if FCKTRA_DEBUG > 0 163 write(LUPRI,*) 'DEBUG: SIR_FCKTRA_CTL called, debug level', 164 & FCKTRA_DEBUG 165 write(LUPRI,*) 'DEBUG: ITRSIR, NTRLVL, LWORK:',ITRSIR,NTRLVL,LWORK 166 flush(LUPRI) 167#endif 168 169! call N_TRASET(sirntra.F) in order to be able to call N_TRALIM(sirntra.F) 170! which sets : 171! ICDTRA: index array for C,D distributions (**|CD) to calculate 172! ITRTYP(1:NORBT) = number of integral indices in which this orbital 173! enters (i.e. 0,1,2,3, or 4) 174 call N_TRASET(-1,LWORK) 175#if FCKTRA_DEBUG > 0 176 IPRTRA = MAX(IPRTRA,FCKTRA_DEBUG) 177 write(lupri,*) 'DEBUG: IPRTRA set to',IPRTRA 178#endif 179 allocate(ITRTYP(NORBT)) 180 call N_TRALIM(NTRLVL,ICDTRA,ITRTYP) 181 deallocate ( ITRTYP ) 182 183 REWIND(LUINTM) 184 WRITE (LUINTM) (0, I=1,8) ! info not used for FCKTRA 185 WRITE (LUINTM) LBINTM,JTRSIR,NSYM,NORB,NBAS 186 WRITE (LUINTM) CMO(1:NCMOT) 187 WRITE (LUINTM) NISH(1:8), NASH(1:8) 188 189 call GETDAT(LAB123(2),LAB123(3)) ! place date in LAB123(2) and time in LAB123(3) 190 WRITE (LUINTM) LAB123,INT_LABEL 191 WRITE (LUINTM) NNORBX, ICDTRA 192! 193!... Open direct access file for final MO integrals 194! (the ICDTRA array which links a C,D to a record has been saved on LUINTM) 195! 196 if (TYPE_TEXT .eq. 'COUL' ) then 197 198 CALL GPINQ('MO2INT_FCKTRA','EXIST',FEXIST) 199 IF (FEXIST) THEN 200 IF (LUMINT .LE. 0) CALL GPOPEN(LUMINT, 201 & 'MO2INT_FCKTRA','OLD','DIRECT',' ',LBUF,OLDDX) 202 CALL GPCLOSE(LUMINT,'DELETE') 203 END IF 204 LUMINT = -1 205 CALL GPOPEN(LUMINT,'MO2INT_FCKTRA','NEW', 206 & 'DIRECT','UNFORMATTED',LBUF,OLDDX) 207 else if (TYPE_TEXT .eq. 'EXCH') then 208 209 CALL GPINQ('DRC2INT_FCKTRA','EXIST',FEXIST) 210 IF (FEXIST) THEN 211 IF (LUMINT .LE. 0) CALL GPOPEN(LUMINT, 212 & 'DRC2INT_FCKTRA','OLD','DIRECT',' ',LBUF,OLDDX) 213 CALL GPCLOSE(LUMINT,'DELETE') 214 END IF 215 LUMINT = -1 216 CALL GPOPEN(LUMINT,'DRC2INT_FCKTRA','NEW', 217 & 'DIRECT','UNFORMATTED',LBUF,OLDDX) 218 else 219 call QUIT('Unknown TYPE_TEXT : '//TYPE_TEXT) 220 end if 221 222 223 CALL FCKTRA_CTL(TYPE_TEXT, ICDTRA, THR_MOTWOINT, CMO, WORK, LWORK) 224 225 226 WRITE (LUINTM) LAB123,TABLE(3) 227 REWIND(LUINTM) 228 CALL GPCLOSE(LUINTM, 'KEEP') 229 230 deallocate ( ICDTRA ) 231 9999 call QEXIT('SIR_FCKTRA_CTL') 232 return 233 end 234 235! =================================================================== 236! sirfcktra.F section 1: calculate MO integrals 237! =================================================================== 238 subroutine FCKTRA_CTL( TYPE_TEXT, ICDTRA, THR_MOTWOINT, CMO, 239 & WORK, LWORK) 240! 241! Written by Hans Joergen Aa. Jensen 2016 242! 243 implicit none 244 character*4 TYPE_TEXT 245 integer ICDTRA(1:NNORBX), LWORK 246 real*8 THR_MOTWOINT, CMO(*), WORK(LWORK) 247#include "priunit.h" 248 249! Used from include files: 250! gnrinf.h : PARCAL 251! inforb.h : NNORBX, NORBT, ??? 252! inftra.h : IPRTRA, LBUF 253! inftap.h : LUSUPM 254! cbihr2.h : IFTHRS 255! ccom.h : THRS 256#include "iratdef.h" 257#include "maxaqn.h" 258#include "gnrinf.h" 259#include "inforb.h" 260#include "inftra.h" 261#include "inftap.h" 262#include "cbihr2.h" 263#include "ccom.h" 264 265 integer MX_NDMAT, LUSUPM_save, IFTHRS_save 266 integer KFREE, LFREE, KFRSAV, KUCMO, KDMAT, KFMAT, I, J 267 integer MISH_TEST(8), MASH_TEST(8), IDAMAX 268 integer, allocatable :: ICD_NODE(:) 269 character*8 TABLE(2),LAB123(3) 270 data TABLE /'MUL_2INT','END INTM'/ 271 data LAB123/'********','********','********'/ 272 logical FEXIST, OLDDX, FNDLAB 273 274C --- start of executable code 275 276 call QENTER('FCKTRA_CTL') 277 278#if FCKTRA_DEBUG > 0 279 IPRTRA = MAX(IPRTRA,FCKTRA_DEBUG) 280 write(lupri,*) 'FCKTRA_DEBUG: IPRTRA set to',IPRTRA 281 write(lupri,*) 'FCKTRA_DEBUG: THR_MOTWOINT ',THR_MOTWOINT 282 flush(lupri) 283#endif 284 285 IF (IPRTRA .GT. 0) THEN 286! WRITE(LUPRI,'(A,1P,D10.2)') 287! & ' Transformation screening threshold:', THR_MOTWOINT 288 FLUSH(LUPRI) 289 END IF 290 291 IFTHRS_save = IFTHRS 292 IFTHRS = NINT(-LOG10(THR_MOTWOINT)) 293 294 295 KFREE = 1 296 LFREE = LWORK 297 KFRSAV = KFREE 298 299 300 if (FCKTRA_TYPE .eq. 1 .AND. PARCAL) then 301 302#ifdef VAR_MPI 303 CALL FCKTRA_DISTRIBUTED_MASTER( TYPE_TEXT, ICDTRA, IFTHRS, 304 & CMO, WORK, LWORK) 305 306#else 307 CALL QUIT('Parallel FCKTRA only possible if compiled with MPI') 308#endif 309 else ! not PARCAL 310 311 call MEMGET2('REAL','UCMO', KUCMO,NBAST*NORBT, 312 & WORK,KFREE,LFREE) 313 314 MX_NDMAT = LFREE/(4*N2BASX) ! use half of free memory for DMAT and FMAT 315 MX_NDMAT = MIN(NNORBX,MX_NDMAT) 316 call MEMGET2('REAL','DMAT', KDMAT,MX_NDMAT*N2BASX, 317 & WORK,KFREE,LFREE) 318 call MEMGET2('REAL','FMAT', KFMAT,MX_NDMAT*N2BASX, 319 & WORK,KFREE,LFREE) 320 321 call UPKCMO(CMO,WORK(KUCMO)) 322 323 allocate(ICD_NODE(NNORBX)) 324 ICD_NODE(:) = 0 325 326 LUSUPM_save = LUSUPM 327 if (FCKTRA_TYPE .eq. 3) LUSUPM = -9999 ! use RDSUPM 328 if (FCKTRA_TYPE .eq. 4) LUSUPM = -1 ! do not use RDSUPM 329 330 call FCKTRA_FCK_DISTRIBUTED(TYPE_TEXT, ICDTRA, ICD_NODE, 331 & WORK(KUCMO), MX_NDMAT, 332 & WORK(KDMAT), WORK(KFMAT), WORK(KFREE), LFREE) 333 334 LUSUPM = LUSUPM_save 335 deallocate ( ICD_NODE ) 336 337 call MEMREL('after FCKTRA_FCK_DISTRIBUTED',WORK, 338 & 1,KFRSAV,KFREE,LFREE) 339 340 end if ! if (FCKTRA_TYPE .eq. 1 .AND. PARCAL) then ... else ... 341 342 IFTHRS = IFTHRS_save 343 call QEXIT('FCKTRA_CTL') 344 return 345 end 346 347 subroutine FCKTRA_FCK_DISTRIBUTED( TYPE_TEXT, ICDTRA, 348 & ICD_NODE, UCMO, MX_NDMAT, DMAT, FMAT, WORK, LWORK) 349! 350! Written by Hans Joergen Aa. Jensen 2016 351! 352 implicit none 353 character*4 TYPE_TEXT 354 integer ICDTRA(NNORBX), ICD_NODE(NNORBX) 355 integer MX_NDMAT, LWORK 356 real*8 DMAT(N2BASX,MX_NDMAT), FMAT(N2BASX,MX_NDMAT) 357 real*8 UCMO(NBAST,NORBT), WORK(LWORK) 358#include "priunit.h" 359 360! Used from include files: 361! infpri.h : IPRFCK 362! inftra.h : IPRTRA 363! infinp.h : DIRFCK, ADDSRI 364! inforb.h : N2BASX, NORBT,NNORBX, ... 365! infind.h : ISMO(:) 366! inftap.h : LUMINT 367! dftcom.h : HFXFAC 368! cbihrs.h : NOSSUP, ONLY_J 369! infpar.h : MYNUM, MASTER 370! mtags.h : MTAG7 371#include "maxorb.h" 372#include "maxash.h" 373#include "infpri.h" 374#include "inftra.h" 375#include "infinp.h" 376#include "inforb.h" 377#include "infind.h" 378#include "inftap.h" 379#include "dftcom.h" 380#include "cbihrs.h" 381#include "infpar.h" 382#include "mtags.h" 383 384 integer ICD, IC, ID, NDMAT, L_H2CD 385 integer, allocatable :: ICD_REC(:), IFCTYP(:), ISYMDM(:) 386 logical FEXIST, OLDDX, FCKTRA_TIMING 387 real*8 TCPU0,TCPU1,TCPU2,TWAL0,TWAL1,TWAL2 388 real*8 TCPU_FCK,TCPU_WR,TWAL_FCK,TWAL_WR 389 390 logical NOSSUP_SAVE, ONLY_J_SAVE, ADDSRI_save 391 logical ABCD_EQ_BACD 392 real*8 HFXFAC_SAVE 393 394 allocate ( ICD_REC(MX_NDMAT) ) 395 allocate ( IFCTYP(MX_NDMAT) ) 396 allocate ( ISYMDM(MX_NDMAT) ) 397 398 FCKTRA_TIMING = IPRTRA .gt. 0 .or. FCKTRA_DEBUG .ge. 0 399 IF (FCKTRA_TIMING) THEN 400 call GETTIM(TCPU0,TWAL0) 401 TCPU_FCK = 0.0D0 402 TCPU_WR = 0.0D0 403 TWAL_FCK = 0.0D0 404 TWAL_WR = 0.0D0 405 END IF 406#if FCKTRA_DEBUG > 0 407 write(LUPRI,*) 'FCKTRA_FCK_DISTRIBUTED called' 408 write(LUPRI,*) 'MYNUM, MX_NDMAT, LWORK:',MYNUM,MX_NDMAT,LWORK 409#if FCKTRA_DEBUG > 5 410 write(LUPRI,*) 'Test output of UCMO:' 411 call OUTPUT(UCMO,1,NBAST,1,NORBT,NBAST,NORBT,-1,LUPRI) 412#endif 413 flush(lupri) 414#endif 415 416! if srdft, only long range integrals 417 ADDSRI_save = ADDSRI 418 ADDSRI = .FALSE. 419 420! for RDSUPM calls 421 NOSSUP_SAVE = NOSSUP 422 NOSSUP = .true. 423 ONLY_J_SAVE = ONLY_J 424 425 HFXFAC_SAVE = HFXFAC 426 427 if (TYPE_TEXT(1:4) .eq. 'COUL') then ! Coulomb/Mulliken integrals 428 HFXFAC = 0.0D0 429 ONLY_J = .true. 430 ABCD_EQ_BACD= .true. 431 IFCTYP(1:MX_NDMAT) = 11 ! symmetric density matrix, only Coulomb 432 else if (TYPE_TEXT(1:4) .eq. 'EXCH') then ! Exchange/Dirac integrals 433 HFXFAC = 2.0D0 434 ABCD_EQ_BACD= .false. 435 IFCTYP(1:MX_NDMAT) = 2 ! non-symmetric density matrix, only exchange 436 if (LUSUPM .ne. -1) then 437 call quit( 438 & 'ERROR: FCKTRA "EXCH" integrals not yet for "SUPMAT"') 439 end if 440 else if (TYPE_TEXT(1:4) .eq. 'ANTI') then ! antisymmetrized exchange integrals 441 HFXFAC = 2.0D0 442 ABCD_EQ_BACD= .false. 443 IFCTYP(1:MX_NDMAT) = 3 ! non-symmetric density matrix, Coulomb and exchange 444 CALL QUIT('FCKTRA ERROR: '// 445 & 'Antisymmetrized Dirac integrals not implemented yet.') 446 if (LUSUPM .ne. -1) then 447 call quit( 448 & 'ERROR: FCKTRA "EXCH" integrals not yet for "SUPMAT"') 449 end if 450 else 451 write(LUPRI,'(///2A)') 'FCKTRA ERROR: '// 452 & 'Unknown MO integral type: ',TYPE_TEXT 453 CALL QUIT('FCKTRA ERROR: '// 454 & 'Requested MO integral type not implemented.') 455 end if 456 457 458 NDMAT = 0 459 ICD = 0 460 do IC = 1, NORBT 461 do ID = 1, IC 462 ICD = ICD + 1 463#if FCKTRA_DEBUG > 4 464 write (lupri,*) 'mynum, icd, icdtra, icd_node', 465 & mynum, icd, icdtra(icd), icd_node(icd) ; flush(lupri) 466#endif 467 468 if (ICDTRA(ICD) .gt. 0 .and. 469 & ICD_NODE(ICD) .eq. MYNUM) then 470 471 NDMAT = NDMAT + 1 472 ISYMDM(NDMAT) = MULD2H(ISMO(IC), ISMO(ID)) 473 ICD_REC(NDMAT) = ICDTRA(ICD) 474 call FCKTRA_CD_DMAT(UCMO(1,IC),UCMO(1,ID), 475 & DMAT(1,NDMAT),ABCD_EQ_BACD) 476#if FCKTRA_DEBUG > 2 477 write(lupri,*) 478 & 'mynum, NDMAT, IC, ID, ICD, ICDTRA, ISYMDM:', 479 & mynum,NDMAT,IC,ID,ICD,ICDTRA(ICD), ISYMDM(NDMAT) 480#if FCKTRA_DEBUG > 5 481 write(lupri,*) 'DMAT no. NDMAT' 482 call output(DMAT(1,NDMAT),1,NBAST,1,NBAST, 483 & NBAST,NBAST,-1,LUPRI) 484#endif 485 flush(lupri) 486#endif 487 end if 488 489 if (NDMAT .eq. MX_NDMAT .or. 490 & (NDMAT .gt. 0 .and. ICD .eq. NNORBX) ) then 491 call DZERO(FMAT,NDMAT*N2BASX) 492 if (FCKTRA_TIMING) call GETTIM(TCPU1,TWAL1) 493 494 if (LUSUPM .ne. -1) then 495#if FCKTRA_DEBUG > 1 496 write (lupri,*) 'calling RDSUPM, NDMAT',NDMAT 497 flush(lupri) 498#endif 499 call RDSUPM(NDMAT,FMAT,DMAT,ISYMDM,WORK,LWORK) 500 else 501#if FCKTRA_DEBUG > 2 502 IPRFCK = MAX(IPRFCK,FCKTRA_DEBUG) 503 write (lupri,*)'calling SIRFCK; NDMAT, IPRFCK, DIRFCK' 504 & ,NDMAT, IPRFCK, DIRFCK 505 flush(lupri) 506#endif 507 call SIRFCK(FMAT,DMAT,NDMAT,ISYMDM,IFCTYP,DIRFCK, 508 & WORK,LWORK) 509 end if 510 if (FCKTRA_TIMING) then 511 call GETTIM(TCPU2,TWAL2) 512#if FCKTRA_DEBUG > 2 513 write (lupri,*) 'CPU and wall time', 514 & TCPU2-TCPU1,TWAL2-TWAL1 515 flush(lupri) 516#endif 517 TCPU_FCK = TCPU_FCK + TCPU2 - TCPU1 518 TWAL_FCK = TWAL_FCK + TWAL2 - TWAL1 519 end if 520 521 ! (**|CD) in AO basis now in FMAT 522 if (ABCD_EQ_BACD) then 523 call FCKTRA_AB_TO_MO(NDMAT,UCMO,FMAT,DMAT) 524 ! (**|CD) in MO basis now in DMAT 525 CALL FCKTRA_AB_PACK(NDMAT,DMAT,FMAT,ISYMDM) 526 ! packed (**|CD) in MO basis now in FMAT 527 L_H2CD = NNORBT 528 else 529 call FCKTRA_AB_TO_MO(NDMAT,UCMO,FMAT,FMAT) 530 ! (**|CD) in MO basis now in FMAT 531 call quit('packing abcd .ne. bacd not implemented') 532 L_H2CD = N2ORBT 533 end if 534#if FCKTRA_DEBUG > 2 535 write (lupri,*) 'Back from FCKTRA_AB_TO_MO'; flush(lupri) 536#endif 537 if (MYNUM .ne. 0) then 538 call MPIXSEND(MYNUM,1,'INTEGER',MASTER,MTAG7) 539 call MPIXSEND(NDMAT,1,'INTEGER',MASTER,MTAG7) 540 call MPIXSEND(ICD_REC,NDMAT,'INTEGER',MASTER,MTAG7) 541 call MPIXSEND(FMAT,NDMAT*L_H2CD,'DOUBLE',MASTER,MTAG7) 542#if FCKTRA_DEBUG > 1 543 write(lupri,*)'done to master, mtag7',master,mtag7 544 flush(lupri) 545#endif 546 else 547#if FCKTRA_DEBUG > 1 548 write(lupri,*) 'writing to MOTWOINT, ndmat',ndmat 549 flush(lupri) 550#endif 551 call FCKTRA_WR_LUMINT(LUMINT,NDMAT,FMAT, 552 & L_H2CD,ICD_REC) 553 end if 554 if (FCKTRA_TIMING) then 555 call GETTIM(TCPU1,TWAL1) 556 TCPU_WR = TCPU_WR + TCPU1 - TCPU2 557 TWAL_WR = TWAL_WR + TWAL1 - TWAL2 558 end if 559 NDMAT = 0 560 end if 561 end do ! ID = 1, IC 562 end do ! IC = 1, NORBT 563 564 if (FCKTRA_TIMING) then 565 call GETTIM(TCPU1,TWAL1) 566 write(lupri,'(A,2F20.3)') 'CPU and WALL time (s) AO Fock', 567 & TCPU_FCK,TWAL_FCK 568 write(lupri,'(A,2F20.3)') 'CPU and WALL time (s) AO->MO ', 569 & TCPU_WR,TWAL_WR 570 TCPU2 = TCPU1 - TCPU0 - TCPU_FCK - TCPU_WR 571 TWAL2 = TWAL1 - TWAL0 - TWAL_FCK - TWAL_WR 572 write(lupri,'(A,2F20.3)') 'CPU and WALL time (s) rest ', 573 & TCPU2,TWAL2 574 flush(lupri) 575 end if 576 577 HFXFAC = HFXFAC_SAVE 578 NOSSUP = NOSSUP_SAVE 579 ONLY_J = ONLY_J_SAVE 580 ADDSRI = ADDSRI_save 581 582 deallocate ( ICD_REC ) 583 deallocate ( IFCTYP ) 584 deallocate ( ISYMDM ) 585 586 return 587 end 588! =================================================================== 589 subroutine FCKTRA_CD_DMAT(CMO_C,CMO_D,DMAT_CD,SYM_CD) 590! 591! Written by Hans Joergen Aa. Jensen 2016 592! 593 implicit none 594 real*8 CMO_C(NBAST), CMO_D(NBAST), DMAT_CD(NBAST,NBAST) 595 logical SYM_CD 596 597! Used from include files: 598! inforb.h : NBAST 599#include "inforb.h" 600 601 integer JA, JB 602 603 do JB = 1, NBAST 604 do JA = 1,NBAST 605 DMAT_CD(JA,JB) = CMO_C(JA)*CMO_D(JB) 606 end do 607 end do 608 609 if (SYM_CD) then ! symmetrize 610 call DGETSI(NBAST,DMAT_CD,DMAT_CD) 611 end if 612 613 end 614! =================================================================== 615 subroutine FCKTRA_AB_TO_MO(NDMAT,UCMO,H2CD_AO,H2CD_MO) 616! 617! Written by Hans Joergen Aa. Jensen 2016-2017 618! 619! Transform (ab|CD) to (AB|CD); ABCD MO indices, ab AO indices 620! 621! H2CD_AO contains NDMAT "CD" records of (ab|CD) on input, 622! H2CD_MO contains NDMAT "CD" records of (AB|CD) on output. 623! 624! Note: H2CD_AO and H2CD_MO may be equivalenced 625! (i.e. share memory), as N2ORBX .le. N2BASX 626! 627! 628 implicit none 629 integer NDMAT 630 real*8 UCMO(NBAST,NORBT) 631 real*8 H2CD_AO(N2BASX,NDMAT), H2CD_MO(N2ORBX,NDMAT) 632 real*8, allocatable :: H2CD_TMP(:) 633#include "priunit.h" 634 635! Used from include files: 636! inforb.h : NBAST,NORBT,N2BASX,N2ORBX, ... 637#include "inforb.h" 638 639 integer IDMAT 640 641 allocate (H2CD_TMP(N2BASX)) 642 643! TO DO: this should be done on GPU if available /hjaaj Aug. 2016 644 do IDMAT = 1, NDMAT 645 call DGEMM('T','N',NORBT,NBAST,NBAST, 1.0D0, 646 & UCMO,NBAST, 647 & H2CD_AO(1,IDMAT),NBAST, 0.0D0, 648 & H2CD_TMP,NORBT) 649 call DGEMM('N','N',NORBT,NORBT,NBAST, 1.0D0, 650 & H2CD_TMP,NORBT, 651 & UCMO,NBAST, 0.0D0, 652 & H2CD_MO(1,IDMAT),NORBT) 653#if FCKTRA_DEBUG > 5 654 write (LUPRI,*) 'H2CD_MO, record',IDMAT 655 call OUTPUT(H2CD_MO(1,IDMAT),1,NORBT,1,NORBT, 656 & NORBT,NORBT,-1,LUPRI) 657#endif 658 end do 659 660 deallocate(H2CD_TMP) 661 662 end 663! =================================================================== 664 subroutine FCKTRA_AB_PACK(NDMAT,H2CD_MO,H2CD_MO_PK,ICDSYM) 665! 666! Written by Hans Joergen Aa. Jensen 2016-2017 667! 668! Purpose: Pack Mulliken MO integrals before writing them to disk, 669! using (AB|CD) = (BA|CD) and (AB|CD)=0 if sym(AB) .ne. sym(CD). 670! 671! H2CD_MO contains NDMAT "CD" records of (AB|CD) on input, 672! packed into H2CD_MO_PK on output. 673! 674! 675 implicit none 676 integer NDMAT, ICDSYM(NDMAT) 677 real*8 H2CD_MO(N2ORBX,NDMAT), H2CD_MO_PK(NNORBT,NDMAT) 678 679#include "priunit.h" 680 681! Used from include files: 682! inforb.h : N2ORBX,NNORBT, NSYM, ... 683#include "inforb.h" 684 685 integer IDMAT 686 real*8, allocatable :: H2CD_TMP(:) 687 688 if (NSYM > 1) allocate (H2CD_TMP(NNORBX)) 689 690 do IDMAT = 1, NDMAT 691 if (NSYM > 1) then 692 call DGETSP(NORBT, H2CD_MO(1,IDMAT), H2CD_TMP) 693 call TRDPAK( H2CD_TMP, H2CD_MO_PK(1,IDMAT), 694 & NORB,IORB,NORBT,ICDSYM(IDMAT),1) 695 else 696 call DGETSP(NORBT, H2CD_MO(1,IDMAT), H2CD_MO_PK(1,IDMAT)) 697 end if 698#if FCKTRA_DEBUG > 5 699 write (LUPRI,*) 'IDMAT, sym', 700 & IDMAT,ICDSYM(IDMAT) 701 if (ICDSYM(idmat) == 1) then 702 write (LUPRI,*) 'Packed H2CD_MO' 703 call OUTPKB(H2CD_MO_PK(1,IDMAT),NORB,NSYM,-1,LUPRI) 704 end if 705#endif 706 end do 707 708 if (NSYM > 1) deallocate(H2CD_TMP) 709 710 end 711 712! =================================================================== 713 subroutine FCKTRA_WR_LUMINT(LUMINT,NDMAT,H2CD,L_H2CD,ICD_REC) 714! 715! Written by Hans Joergen Aa. Jensen 2016 716! 717 implicit none 718 integer LUMINT, NDMAT, L_H2CD, ICD_REC(NDMAT), IOS 719 real*8 H2CD(L_H2CD,NDMAT) 720 721#include "priunit.h" 722 723 integer IDMAT 724 725 do IDMAT = 1, NDMAT 726 ! write NDMAT (**|CD) MO integral distributions to file 727 write(LUMINT, rec=ICD_REC(IDMAT), iostat=IOS) 728 & H2CD(1:L_H2CD,IDMAT) 729 if (IOS .ne. 0) then 730 call QENTER('FCKTRA_WR_LUMINT') 731 write(lupri,*) 'IDMAT, rec, IOS',IDMAT,ICD_REC(IDMAT),IOS 732 write(lupri,*) 'LUMINT, L_H2CD',LUMINT,L_H2CD 733 call QUIT('IOS .ne. 0') 734 end if 735 end do 736 737 end 738 739! =================================================================== 740! sirfcktra.F section 2: FCKTRA type 1 routines (if VAR_MPI) 741! =================================================================== 742#ifdef VAR_MPI 743 subroutine FCKTRA_DISTRIBUTED_MASTER( TYPE_TEXT, ICDTRA, IFTHRS, 744 & CMO, WORK, LWORK) 745! 746! Written by Hans Joergen Aa. Jensen 2016 747! 748 implicit none 749 character*4 TYPE_TEXT 750 integer ICDTRA(NNORBX), IFTHRS, LWORK 751 real*8 CMO(NCMOT), WORK(LWORK) 752 753#include "priunit.h" 754 755#include "mpif.h" 756 757! Used from include files: 758! gnrinf.h : CHIVAL 759! inforb.h : N2BASX, NORBT,NNORBX, ... 760! infind.h : ISMO(:) 761! inftra.h : IPRTRA 762! inftap.h : LUMINT 763! infpar.h : MYNUM 764! iprtyp.h : defined parallel calculation types 765! mtags.h : MTAG7 766 767#include "maxorb.h" 768#include "maxash.h" 769#include "gnrinf.h" 770#include "inforb.h" 771#include "infind.h" 772#include "inftra.h" 773#include "inftap.h" 774#include "infpar.h" 775#include "iprtyp.h" 776#include "mtags.h" 777 778 integer, allocatable :: ICD_NODE(:) 779 integer IPRTYP, IPRINT_slaves, ICD, NCDTRA, IWHO, NWHO 780 integer KFREE, LFREE 781 integer NDMAT, IDIST, N2GAB, L_H2CD 782 783 call qenter('FCKTRA_DISTRIBUTED_MASTER') 784 785 if (TYPE_TEXT .eq. 'COUL') then 786 L_H2CD = NNORBT 787 else if (TYPE_TEXT .eq. 'EXCH') then 788 L_H2CD = N2ORBT 789 else 790 call QUIT('Unknown TYPE_TEXT '//TYPE_TEXT) 791 end if 792 793! start nodes on the distributed 2-electron integral transformation task 794 795 IPRTYP = CALL_FCKTRA_DISTRIBUTED ! get code number from iprtyp.h 796 call MPIXBCAST(IPRTYP,1,'INTEGER',MASTER) 797 IPRINT_slaves = IPRTRA ! TODO set print level in input 798 call MPIXBCAST(IPRINT_slaves,1,'INTEGER',MASTER) 799 800! Now nodes are in subroutine FCKTRA_DISTRIBUTED_NODE. 801! First make sure basis set information is known to nodes 802 803 KFREE = 1 804 LFREE = LWORK 805 CALL HER_sendreceive_molinfo(IPRINT_slaves,WORK,KFREE,LFREE) 806 807! Transfer integral transformation information to nodes 808 809 allocate( ICD_NODE(NNORBX) ) 810 ICD_NODE(:) = 0 811 NCDTRA = 0 812 DO ICD = 1,NNORBX 813 IF (ICDTRA(ICD) .EQ. 0) CYCLE 814 ICD_NODE(ICD) = MOD(NCDTRA,NODTOT) + 1 815 NCDTRA = NCDTRA + 1 816 END DO 817 818#if FCKTRA_DEBUG > 2 819 write (*,*) 'FCKTRA_DISTRIBUTED_MASTER has started slaves' 820 write (lupri,*) 'FCKTRA_DISTRIBUTED_MASTER has started slaves' 821 write (lupri,*) 'IPRTYP,IPRINT_slaves=',IPRTYP,IPRINT_slaves 822 write (lupri,*) 'ICD_NODE is' 823 write (lupri,'(10I4,5X,10I4)') ICD_NODE(1:NNORBX) 824 flush(lupri) 825#endif 826 call MPIXBCAST(TYPE_TEXT,4,'CHARACTER',MASTER) 827 call MPIXBCAST(IFTHRS,1,'INTEGER',MASTER) 828 call MPIXBCAST(CHIVAL,1,'DOUBLE',MASTER) 829 call MPIXBCAST(ICDTRA,NNORBX,'INTEGER',MASTER) 830 call MPIXBCAST(ISMO,NORBT,'INTEGER',MASTER) 831 call MPIXBCAST(ICD_NODE,NNORBX,'INTEGER',MASTER) 832 call MPIXBCAST(CMO,NCMOT,'DOUBLE',MASTER) 833 834 IDIST = 0 835 100 CONTINUE 836 837 IWHO = -1 838 call MPIXRECV(NWHO,1,'INTEGER',IWHO,MTAG7) 839 call MPIXRECV(NDMAT,1,'INTEGER',NWHO,MTAG7) 840 IF (NDMAT .LE. 0) THEN 841 write(LUPRI,*) 842 & 'ERROR: master received ndmat.le.0, ndmat, nwho',ndmat,nwho 843 call quit('NDMAT .le. 0') 844 END IF 845 ! we reuse ICD_NODE to store ICD_REC 846 ! CALL MPIXRECV(ICD_REC,NDMAT,'INTEGER',NWHO,MTAG7) 847 CALL MPIXRECV(ICD_NODE,NDMAT,'INTEGER',NWHO,MTAG7) 848 IF (NDMAT*L_H2CD .GT. LWORK) CALL QUIT('NDMAT*L_H2CD .gt. lwork') 849 CALL MPIXRECV(WORK,NDMAT*L_H2CD,'DOUBLE',NWHO,MTAG7) 850 851 call FCKTRA_WR_LUMINT(LUMINT,NDMAT,WORK,L_H2CD,ICD_NODE) 852 853#if FCKTRA_DEBUG > 1 854 write (lupri,*) 'Received from nwho =',nwho 855 write (lupri,*) 'Received and saved ndmat =',ndmat 856 write(lupri,'(A,(20I5))') 'The CD records:',ICD_NODE(1:NDMAT) 857#endif 858 859 IDIST = IDIST + NDMAT 860 IF (IDIST .GT. NCDTRA) THEN 861 write(lupri,*) 'ERROR: IDIST .gt. NCDTRA',IDIST,NCDTRA 862 call quit('IDIST .gt. NCDTRA, uha') 863 else if (IDIST .LT. NCDTRA) then 864 go to 100 865 end if 866 867 deallocate( ICD_NODE ) 868 call qexit('FCKTRA_DISTRIBUTED_MASTER') 869 return 870 end subroutine FCKTRA_DISTRIBUTED_MASTER 871! =================================================================== 872 subroutine FCKTRA_DISTRIBUTED_NODE(WORK,LWORK,IPRINT) 873! 874! Written by Hans Joergen Aa. Jensen 2016 875! 876 implicit none 877 integer LWORK, IPRINT 878 real*8 WORK(LWORK) 879 880 integer, allocatable :: ICDTRA(:), ICD_NODE(:) 881 integer ICD, NCDTRA, MX_NDMAT 882 integer KFREE, LFREE, KDMAT, KFMAT, KUCMO 883 logical PARHER_save, SLAVE_save, DIRFCK_save 884 character*4 TYPE_TEXT 885#include "priunit.h" 886 887! Used from include files: 888! gnrinf.h : PARCAL, CHIVAL 889! infinp.h : DIRFCK 890! infind.h : ISMO(:) 891! inforb.h : NORBT, NNORBX, NCMOT 892! infpar.h : MYNUM,SLAVE,PARHER,? 893! inftap.h : LUSUPM 894! cbihr2.h : IFTHRS 895#include "maxorb.h" 896#include "maxash.h" 897#include "gnrinf.h" 898#include "infinp.h" 899#include "infind.h" 900#include "inforb.h" 901#include "infpar.h" 902#include "inftap.h" 903#include "cbihr2.h" 904 905! mpif.h : for MPI 906#include "mpif.h" 907 908! First make sure basis set information is known to this node 909 910 KFREE = 1 911 LFREE = LWORK 912 CALL HER_sendreceive_molinfo(IPRINT,WORK,KFREE,LFREE) 913 914#if FCKTRA_DEBUG > 1 915 write (lupri,*) 916 &'FCKTRA_DISTRIBUTED_NODE has started, MYNYM =',MYNUM 917 write (lupri,*) 'NORBT, NNORBX, NCMOT', 918 & NORBT, NNORBX, NCMOT 919 flush(lupri) 920#endif 921 922! We check that orbital information has been transferred to slaves, 923! and then we call SETORB to define inforb.h and infind.h (we need ISMO(:)) 924 925 IF (NORBT .LE. 0) THEN 926 CALL QUIT( 927 & 'ERROR: Orbital information not transferred to slaves.') 928 END IF 929 CALL SETORB 930 931! Some set up for this task (AO to MO transformation) 932 933 KFREE = 1 934 LFREE = LWORK 935 936 SLAVE_save = SLAVE 937 SLAVE = .FALSE. ! Be sure the Fock matrices are NOT constructed in parallel, 938 939 PARHER_save = PARHER 940 PARHER = .FALSE. ! be sure the Fock matrices are NOT constructed in parallel 941 ! we are parallelizing on the distributions 942 943 DIRFCK_save = DIRFCK 944 DIRFCK = .TRUE. ! NEVER read from AOTWOINT when parallel 945 946 LUSUPM = -1 ! Do not use AO integrals from file AO2_JINT 947 948 949 allocate( ICDTRA(NNORBX) ) 950 allocate( ICD_NODE(NNORBX) ) 951 952 call MEMGET2('REAL','UCMO', KUCMO,NBAST*NORBT, 953 & WORK,KFREE,LFREE) 954 955 call MPIXBCAST(TYPE_TEXT,4,'CHARACTER',MASTER) 956 call MPIXBCAST(IFTHRS,1,'INTEGER',MASTER) ! accuracy of MO integrals (screening) 957 call MPIXBCAST(CHIVAL,1,'DOUBLE',MASTER) 958 call MPIXBCAST(ICDTRA,NNORBX,'INTEGER',MASTER) 959 call MPIXBCAST(ISMO,NORBT,'INTEGER',MASTER) 960 call MPIXBCAST(ICD_NODE,NNORBX,'INTEGER',MASTER) 961 962 call MPIXBCAST(WORK(KFREE),NCMOT,'DOUBLE',MASTER) 963 call UPKCMO(WORK(KFREE),WORK(KUCMO)) 964 965 NCDTRA = 0 966 DO ICD = 1,NNORBX 967 IF (ICD_NODE(ICD) .EQ. MYNUM) NCDTRA = NCDTRA + 1 968 END DO 969 970 IF (NCDTRA .EQ. 0) THEN 971 ! nothing to do for me this time ... 972 GO TO 1000 973 END IF 974 975 MX_NDMAT = LFREE/(4*N2BASX) ! use up to half of free memory for DMAT and FMAT 976 MX_NDMAT = MIN(NCDTRA,MX_NDMAT) 977 978 IF (IPRINT .GT. 0) THEN 979 write(lupri,'(A,2I12)') 'FCKTRA # of distributions and'// 980 & ' # of distributions per load',NCDTRA, MX_NDMAT 981 flush(lupri) 982 END IF 983 984 IF (MX_NDMAT .EQ. 0) THEN 985 CALL QUIT('MX_NDMAT = 0, too little memory') 986 END IF 987 call MEMGET2('REAL','DMAT', KDMAT,MX_NDMAT*N2BASX, 988 & WORK,KFREE,LFREE) 989 call MEMGET2('REAL','FMAT', KFMAT,MX_NDMAT*N2BASX, 990 & WORK,KFREE,LFREE) 991 992#if FCKTRA_DEBUG > 2 993 write(lupri,*) 'NCDTRA,MX_NDMAT,NNORBX =',NCDTRA,MX_NDMAT,NNORBX 994 write (lupri,'(10I4,5X,10I4)') ICD_NODE(1:NNORBX) 995 write(lupri,*) 'chival',chival 996 write (lupri,*) mynum,'node; ICD_NODE is' 997 write (lupri,*) MYNUM,' is calling fcktra_fck_distributed' 998 flush(lupri) 999#endif 1000 1001 CALL FCKTRA_FCK_DISTRIBUTED( TYPE_TEXT, 1002 & ICDTRA, ICD_NODE, WORK(KUCMO), MX_NDMAT, 1003 & WORK(KDMAT), WORK(KFMAT), WORK(KFREE), LFREE) 1004 1005 1006 1000 CONTINUE 1007 1008 deallocate( ICDTRA ) 1009 deallocate( ICD_NODE ) 1010 1011 PARHER = PARHER_save 1012 SLAVE = SLAVE_save 1013 DIRFCK = DIRFCK_save 1014 1015 return 1016 end subroutine FCKTRA_DISTRIBUTED_NODE 1017#endif /* VAR_MPI */ 1018 1019! =================================================================== 1020! sirfcktra.F section 3: read MO integrals 1021! =================================================================== 1022 subroutine FCKTRA_NXTH2M(IC,ID,H2CD,NEEDTP,WORK,IDIST) 1023! 1024! Written by Hans Joergen Aa. Jensen 2016 1025! This version is interface routine for .FCKTRA integral transformation. 1026! 1027! Purpose: 1028! Read next Mulliken two-electron integral distribution (**|cd) 1029! where |cd) distribution is needed according to NEEDTP(ITYPCD) 1030! 1031! Usage: 1032! 1) Set IDIST = 0 before first call. 1033! Do **NOT** change IDIST in calling routine 1034! until last distribution has been read (signalled by IDIST .eq. -1) 1035! 2) if IDIST<0 on input, then IC, ID must be set to correspond to 1036! the distribution wanted in H2CD 1037! 1038 1039 implicit none 1040#include "priunit.h" 1041 1042 integer :: IC, ID, NEEDTP(-4:6), IDIST 1043 real*8 :: H2CD(NORBT,NORBT), WORK(*) 1044 1045! Used from include files: 1046! inforb.h : NORBT, NNORBT, NNORBX, ??? 1047! infind.h : ISMO(:) 1048! inftra.h : LBUF 1049! inftap.h : LUINTM, LUMINT 1050#include "maxorb.h" 1051#include "maxash.h" 1052#include "dummy.h" 1053#include "iratdef.h" 1054#include "inforb.h" 1055#include "infind.h" 1056#include "inftra.h" 1057#include "inftap.h" 1058#include "orbtypdef.h" 1059 1060 integer :: MMORBX, ICD, ICD_1, ICDSYM 1061 integer :: ITYPC, ITYPD, ITYPCD 1062 integer, save :: ICD_SAVE 1063 integer, save, allocatable :: ICDTRA(:) 1064 logical :: OLDDX 1065 1066#if FCKTRA_DEBUG > 3 1067 write(LUPRI,*) 'FCKTRA_NXTH2M called, IC, ID, IDIST =', 1068 & IC, ID, IDIST 1069 flush(LUPRI) 1070#endif 1071 1072 if (IDIST .le. 0) then 1073 if (allocated(ICDTRA)) deallocate(ICDTRA) 1074 allocate(ICDTRA(1:NNORBX)) 1075 if (LUINTM .le. 0) THEN 1076 CALL GPOPEN(LUINTM,FNINTM,'OLD',' ', 1077 & 'UNFORMATTED',IDUMMY,.FALSE.) 1078 END IF 1079 LBUF = IRAT*NNORBT 1080 IF (LUMINT .LE. 0) THEN 1081 CALL GPOPEN(LUMINT,'MO2INT_FCKTRA','OLD','DIRECT',' ', 1082 & LBUF,OLDDX) 1083 END IF 1084 rewind (LUINTM) 1085 call MOLLAB('MUL_2INT',LUINTM,LUPRI) 1086 read (LUINTM) MMORBX, ICDTRA(1:MMORBX) 1087 if (NNORBX .ne. MMORBX) then 1088 call QENTER('FCKTRA_NXTH2M') 1089 write (LUPRI,'(//A,2I10)') 'FATAL ERROR '// 1090 & 'NNORBX on MO2INT_FCKTRA .ne. NNORBX',MMORBX,NNORBX 1091 call QUIT( 1092 & 'NNORBX on MO2INT_FCKTRA .ne. NNORBX in FCKTRA_NXTH2M') 1093 end if 1094 ICD_SAVE = 0 1095 end if 1096 1097 if (IDIST .lt. 0) then 1098 1099 if (IC .le. 0 .or. IC .gt. NORBT .and. 1100 & ID .le. 0 .or. ID .gt. NORBT) then 1101 write (LUPRI,*) 'FCKTRA_NXTH2M: invalid input IC, ID',IC,ID 1102 call quit('FCKTRA_NXTH2M: invalid input IC, ID') 1103 end if 1104 if (IC .ge. ID) then 1105 ICD = IC*(IC-1)/2 + ID 1106 else 1107 ICD = ID*(ID-1)/2 + IC 1108 end if 1109 IDIST = ICDTRA(ICD) 1110 if (IDIST .le. 0) then 1111 call QUIT('IC,ID record no. .le. 0 in FCKTRA_NXTH2M') 1112 end if 1113 deallocate(ICDTRA) 1114 1115 else ! IDIST .ge. 0 1116 1117 100 continue 1118 IDIST = -1 1119 1120 ICD_1 = ICD_SAVE + 1 1121 do ICD = ICD_1,NNORBX 1122 if (ICDTRA(ICD) .gt. 0) then 1123 IDIST = ICDTRA(ICD) 1124 ICD_SAVE = ICD 1125 exit 1126 end if 1127 end do 1128 if (IDIST .eq. -1) then 1129 ! finished 1130 deallocate(ICDTRA) 1131 return 1132 end if 1133 1134 do IC = 1, NORBT 1135 ICD = IC*(IC+1)/2 1136 if (ICD .lt. ICD_SAVE) cycle 1137 ID = ICD_SAVE - ICD + IC 1138 exit 1139 end do 1140 ITYPC = IOBTYP(IC) 1141 ITYPD = IOBTYP(ID) 1142 ITYPCD = IDBTYP(ITYPC,ITYPD) 1143 if ( NEEDTP(ITYPCD) .eq. 0 ) go to 100 1144 1145 end if 1146 1147 ICDSYM = MULD2H( ISMO(IC), ISMO(ID) ) 1148 read ( LUMINT, rec=IDIST ) WORK(1:NNORBT) 1149 call TRDPAK(WORK(1+NNORBT),WORK,NORB,IORB,NORBT, 1150 & ICDSYM,-1) 1151 call DSPTGE(NORBT,WORK(1+NNORBT),H2CD) 1152 1153#if FCKTRA_DEBUG > 9 1154 write (LUPRI,*) 'FCKTRA_NXTH2M, rec,IC,ID',IDIST,IC,ID 1155 write (LUPRI,*) 'Unpacked H2CD, ICDSYM',ICDSYM 1156 call output(H2CD,1,NORBT,1,NORBT,NORBT,NORBT,-1,LUPRI) 1157 if (ICDSYM.eq.1) then ! only print some of them, of ICDSYM.eq.1 1158 write (LUPRI,*) 'Packed H2CDPK' 1159 call outpkb(WORK,NORB,NSYM,-1,LUPRI) 1160 end if 1161#endif 1162 1163 return 1164 end 1165 1166 subroutine FCKTRA_NXTH2D(IC,ID,H2CD,NEEDTP,IDIST) 1167! 1168! Written by Hans Joergen Aa. Jensen 2016 1169! This version is interface routine for .FCKTRA integral transformation. 1170! 1171! Purpose: 1172! Read next Dirac two-electron integral distribution <**|cd> = (*c|*d) 1173! where (cd) distribution is needed according to NEEDTP(ITYPCD) 1174! 1175! 1176 1177 implicit none 1178#include "priunit.h" 1179 1180 integer IC, ID, NEEDTP(*), IDIST 1181 real*8 H2CD(NORBT,NORBT) 1182 1183! Used from include files: 1184! inforb.h : NORBT, ??? 1185#include "inforb.h" 1186 1187 call quit('ERROR: FCKTRA_NXTH2D is not implemented') 1188 1189 1190 end 1191! -- end of DALTON/sirius/sirfcktra.F -- 1192