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======================================================================== 20C Revision 1.2 2000/05/05 11:19:33 hjj 21C bug fix: now always calculate ci diag. in CIDAG, previously it was only 22C done if LUIT2 .ge. 0 (must be a leftover, the code LUIT2 .lt. 0 was not 23C used anywhere, rather UPDGRD expected diagonal to be calculated with 24C LUIT2 .lt. 0 !!!). 25C======================================================================== 26C /* Deck tpblm2 */ 27 SUBROUTINE TPBLM2(A,AT,LBLR,NBLR,LBLC,NBLC,IFLAG,IORD) 28C 29C A BLOCKED MATRIX A IS GIVEN . 30C 31C THE BLOCK STRUCTURE CAN BE OF THE FOLLOWING TYPES 32C IORD = 1 : 33C LOOP OVER BLOCK OF ROWS 34C LOOP OVER BLOCK OF COLUMNS ALLOWED FOR GIVEN ROW BLOCK 35C LOOP OVER COLUMNS OF THIS BLOCK 36C LOOP OVER ROWS OF THIS BLOCK 37C 38C IORD = 2 : 39C LOOP OVER BLOCK OF ROWS 40C LOOP OVER BLOCK OF COLUMNS ALLOWED FOR GIVEN ROW BLOCK 41C LOOP OVER ROWS OF THIS BLOCK 42C LOOP OVER COLUMNS OF THIS BLOCK 43C 44C FOR IORD = 2 ARE THE INDIVIDUAL BLOCKS THUS TRANSPOSED 45C 46C THE COMBINATION OF TWO BLOCKS IABL AND IBBL ARE ALLOWED 47C IF IFLAG(IABL,IBBL) = 1 48C 49C TRANSPOSE THE INDIVIDUAL BLOCKS OF THIS MATRIX TO GIVE AT 50C THE ORDER OF THE BLOCKS ARE NOT CHANGED 51C 52C JEPPE OLSEN , NOVEMBER 1988 53C 54#include "implicit.h" 55 DIMENSION A(*),AT(*) 56 DIMENSION LBLR(NBLR),LBLC(NBLC) 57 DIMENSION IFLAG(NBLR,NBLC) 58C 59 IOFF = 1 60 DO 200 IBLR = 1, NBLR 61 DO 100 IBLC = 1, NBLC 62 IF(IFLAG(IBLR,IBLC) .EQ. 1 ) THEN 63 LR = LBLR(IBLR) 64 LC = LBLC(IBLC) 65 IF( IORD .EQ. 1 ) THEN 66 CALL TRPMAT(A(IOFF),LR,LC,AT(IOFF)) 67 ELSE IF( IORD .EQ. 2 ) THEN 68 CALL TRPMAT(A(IOFF),LC,LR,AT(IOFF)) 69 END IF 70 IOFF = IOFF + LR * LC 71 END IF 72 100 CONTINUE 73 200 CONTINUE 74C 75 RETURN 76 END 77C /* Deck prblm2 */ 78 SUBROUTINE PRBLM2(A,LBLR,NBLR,LBLC,NBLC,IFLAG,IORD) 79C 80C A BLOCKED MATRIX A IS GIVEN . 81C 82C THE BLOCK STRUCTURE CAN BE OF THE FOLLOWING TYPES 83C IORD = 1 : 84C LOOP OVER BLOCK OF ROWS 85C LOOP OVER BLOCK OF COLUMNS ALLOWED FOR GIVEN ROW BLOCK 86C LOOP OVER COLUMNS OF THIS BLOCK 87C LOOP OVER ROWS OF THIS BLOCK 88C 89C IORD = 2 : 90C LOOP OVER BLOCK OF ROWS 91C LOOP OVER BLOCK OF COLUMNS ALLOWED FOR GIVEN ROW BLOCK 92C LOOP OVER ROWS OF THIS BLOCK 93C LOOP OVER COLUMNS OF THIS BLOCK 94C 95C FOR IORD = 2 ARE THE INDIVIDUAL BLOCKS THUS TRANSPOSED 96C 97C THE COMBINATION OF TWO BLOCKS IABL AND IBBL ARE ALLOWED 98C IF IFLAG(IABL,IBBL) = 1 99C 100C PRINT THIS MATRIX ! 101C 102C JEPPE OLSEN , NOVEMBER 1988 103C 104#include "implicit.h" 105#include "priunit.h" 106 DIMENSION A(*) 107 DIMENSION LBLR(NBLR),LBLC(NBLC) 108 DIMENSION IFLAG(NBLR,NBLC) 109C 110 IOFF = 1 111 DO 200 IBLR = 1, NBLR 112 DO 100 IBLC = 1, NBLC 113 IF(IFLAG(IBLR,IBLC) .EQ. 1 ) THEN 114 WRITE(LUPRI,*) ' BLOCK INDICES ',IBLR,IBLC 115 LR = LBLR(IBLR) 116 LC = LBLC(IBLC) 117 IF(IORD.EQ.1) THEN 118 CALL WRTMAT(A(IOFF),LR,LC,LR,LC,1) 119 ELSE IF( IORD .EQ. 2 ) THEN 120 CALL WRTMAT(A(IOFF),LC,LR,LC,LR,1) 121 END IF 122 IOFF = IOFF + LR * LC 123 END IF 124 100 CONTINUE 125 200 CONTINUE 126C 127 RETURN 128 END 129C /* Deck traci */ 130 SUBROUTINE TRACI(NCVEC,CVEC,LCVEC,UMO,XNDXCI,WRK,LFREE,IPRTCI) 131C 132C Driver for TRACI using determinant routines 133C Version with CSF business 134C 135C JULY 15 1988 Jeppe Olsen 136C 137C 138C MOTECC-90: The algorithms used in this module, TRACI, are 139C described in Chapter 8 Section D.5 of MOTECC-90 140C "Counter Rotations of CI coefficients" 141C 142! module dependencies 143 use lucita_mcscf_ci_cfg 144 use lucita_mcscf_ci_interface_procedures 145 use mcscf_or_gasci_2_define_cfg, only : 146 & define_lucita_cfg_dynamic 147 148#include "implicit.h" 149 DIMENSION CVEC(LCVEC,NCVEC), UMO(*), XNDXCI(*), WRK(LFREE) 150C 151C Used from common blocks: 152C PRIUNIT : IPRSTAT 153C INFINP : LSYM 154C INFORB : MULD2H(8,8), NASHT, N2ASHX 155C INFVAR : NCONF 156C DETBAS : KIASTR,KIBSTR, ... 157C STRNUM : NAEL,NBEL, ... 158C CIINFO : NDTASM 159C CSFBAS : Pointers for CSF information 160C CBREOR : SLREOR 161C 162#include "maxorb.h" 163#include "mxpdim.h" 164#include "priunit.h" 165#include "infinp.h" 166#include "inforb.h" 167#include "infvar.h" 168#include "detbas.h" 169#include "strnum.h" 170#include "csfbas.h" 171#include "ciinfo.h" 172#include "cbreor.h" 173#include "dummy.h" 174 LOGICAL CSFEXP 175 integer :: xctype1, xctype2 176C 177C 178 CSFEXP = .NOT.FLAG(27) 179 IF (CSFEXP) THEN 180 NCDET = NDTASM(LSYM) 181 IF (NCONF .NE. NCSASM(LSYM)) THEN 182 WRITE (LUPRI,*) 'TRACI ERROR, NCONF .ne. NCSASM(LSYM):', 183 * NCONF, NCSASM(LSYM) 184 CALL QUIT('TRACI ERROR, NCONF .ne. NCSASM(LSYM)') 185 END IF 186 ELSE 187 NCDET = NCONF 188 IF (NCONF .NE. NDTASM(LSYM)) THEN 189 WRITE (LUPRI,*) 'TRACI ERROR, NCONF .ne. NDTASM(LSYM):', 190 * NCONF, NDTASM(LSYM) 191 CALL QUIT('TRACI ERROR, NCONF .ne. NDTASM(LSYM)') 192 END IF 193 END IF 194c Local memory 195 KLSOT = 1 196 KLC2 = KLSOT + N2ASHX 197 KLFREE = KLC2 + MAX(N2ASHX, NCDET ) 198 IF(CSFEXP) THEN 199 KLDET1 = KLFREE 200 KLDET2 = KLDET1 + NCDET 201 KLFREE = KLDET2 + NCDET 202 END IF 203 lfree_local = lfree - KLFREE 204 205! note: C2 only needs max of nia*nib (see TRACI1) 206 IF ( KLFREE-1 .GT. LFREE ) CALL ERRWRK('TRACI',KLFREE-1,LFREE) 207 208! 1p-matrix 209! write(lupri,*) ' 1p-mat ' 210! call wrtmatmn(UMO,1,nasht**2,1,nasht**2,lupri) 211 212 if(ci_program .eq. 'SIRIUS-CI')then 213 IF (.NOT. SLREOR)THEN 214 CALL DCOPY(N2ASHX,UMO,1,WRK(KLSOT),1) 215 ELSE 216! reorder from Sirius order to Lunar order 217 CALL REORMT(UMO,WRK(KLSOT),NASHT,NASHT,XNDXCI(KLTSOB), 218 & XNDXCI(KLTSOB) ) 219 END IF 220! get single orbital transformation matrix 221 CALL SOTMAT(NASHT,WRK(KLSOT),IFAIL) 222 223C.. TRANSPOSE SOT MATRIX IN ACCORDANCE WITH LUNAR DESCRIPTION 224 CALL DGETRN(WRK(KLSOT),NASHT,NASHT) 225c Store CI offsets in C arrays 226 IPRXXX = MAX(0,IPRTCI - 10) 227 CALL CIOFF(LSYM,1,XNDXCI,IPRXXX) 228c CALL CIOFF(IREFSM,ICORHC,XNDXCI,NTEST) 229 230 else if(ci_program .eq. 'LUCITA ')then 231 call dcopy(n2ashx,umo,1,wrk(klsot),1) 232 call sotmat(nasht,wrk(klsot),ifail) 233 end if 234 235 IF ( IFAIL .NE. 0 ) THEN 236 WRITE(LUPRI,'(///,A,I5,/)') 237 & ' TRACI : ERROR IN SOTMAT, IFAIL = ',IFAIL 238 CALL QUIT(' TRACI : ERROR IN SOTMAT ') 239 END IF 240 241 IF(IPRTCI .GT. 10 ) THEN 242 WRITE(LUPRI,'(/A)') ' SINGLE ORBITAL TRANSFORMATION MATRIX ' 243 CALL WRTMAT(WRK(KLSOT),NASHT,NASHT,NASHT,NASHT,0) 244 END IF 245c 246 DO 450 IVEC = 1, NCVEC 247 IF(CSFEXP) THEN 248c 249c CSF to DET transformation, rotate in DET basis, DET to CSF transformation 250c 251 CALL COPVEC(CVEC(1,IVEC),WRK(KLDET1),NCONF) 252 CALL CSDTVC(WRK(KLDET1),WRK(KLDET2),1,XNDXCI(KDTOC), 253 & XNDXCI(KICTS(1)),LSYM,0,IPRXXX) 254C CSDTVC(CSFVEC,DETVEC,IWAY,DTOCMT,ICTSDT,IREFSM,ICOPY,NTEST) 255 CALL TRACI2(WRK(KLSOT),NASHT,LSYM,NAEL,XNDXCI(KIASTR), 256 & XNDXCI(KTAIJ),XNDXCI(KTATO),XNDXCI(KTASYM), 257 & XNDXCI(KSTASA),XNDXCI(KSTBAA),NBEL, 258 & XNDXCI(KIBSTR),XNDXCI(KTBIJ), 259 & XNDXCI(KTBTO),XNDXCI(KTBSYM),XNDXCI(KSTASB), 260 & XNDXCI(KSTBAB),NASTR,NBSTR,XNDXCI(KIANNI), 261 & WRK(KLDET2),XNDXCI(KCOFF),WRK(KLC2), 262 & NOCTPA,XNDXCI(KNSSOA),XNDXCI(KISSOA), 263 & XNDXCI(KTPFSA),NOCTPB,XNDXCI(KNSSOB), 264 & XNDXCI(KISSOB), 265 & XNDXCI(KTPFSB),XNDXCI(KIOCOC),XNDXCI(KICOOS), 266 & XNDXCI(KCDTAS),MULD2H,MAXSYM,IPRTCI) 267 CALL CSDTVC(WRK(KLDET1),WRK(KLDET2),2,XNDXCI(KDTOC), 268 & XNDXCI(KICTS(1)),LSYM,0,IPRXXX) 269C CSDTVC(CSFVEC,DETVEC,IWAY,DTOCMT,ICTSDT,IREFSM,ICOPY,NTEST) 270 CALL COPVEC(WRK(KLDET1),CVEC(1,IVEC),NCONF) 271 ELSE 272c 273c Rotate in DET expansion 274c 275! default: transform cref but in case we have a multi-root MCSCF switch to xtype1 == 2 276 xctype2 = -1 277 xctype1 = 1 278 if(ncvec > 1) xctype1 = 2 279 280!#define LUCI_DEBUG 281#ifdef LUCI_DEBUG 282 write(lupri,*) ' before tra-ci run: ci-vec # =',ivec 283 call wrtmatmn(CVEC(1,IVEC),1,NCDET,1,NCDET,lupri) 284 write(lupri,*) ' before tra-ci run: sot-mat' 285 call wrtmatmn(WRK(KLSOT),1,n2ashx,1,n2ashx,lupri) 286#endif 287 288 if(ci_program .eq. 'LUCITA ')then 289 call define_lucita_cfg_dynamic(lsym, ! rhs symmetry ! swap see subroutine below 290 & lsym, ! lhs symmetry ! swap see subroutine below 291 & 0, ! no istaci 292 & 0, ! spin for operator 293 & 0, ! spin for operator 294 & 1, ! one root 295 & ispin, ! singlet, doublet, triplet, ... 296 & -1, ! # of davidson iterations 297 & -1, ! 1/2-particle density calcs 298 & iprtci, ! print level 299 & 1.0d-10, ! screening factor 300 & 0.0d0, ! davidson ci convergence 301 & 0.0d0, ! davidson ci convergence (auxiliary) 302 & .false., ! do not calculate nat. orb. occ. num. (internally inside LUCITA) 303 & .false., ! wave function analysis 304 & .true., ! no 2-el operators active 305 & .true., ! integrals provided by / return density matrices to the MCSCF environment 306 & .false., ! calculate 2-particle density matrix 307 & .false., ! calculate transition density matrix 308 & xctype1, ! vector_exchange_type1 309 & xctype2, ! vector_exchange_type2 310 & .false., ! vector exchange active in parallel runs in mc2lu interface (both mc-->lu and lu-->mc) 311 & .false.) ! io2io vector exchange between mcscf and lucita: default mem2io/io2mem 312 313 314 call mcscf_lucita_interface(cvec(1,ivec),wrk(klc2), 315 & wrk(klsot),vdummy, 316 & 'rotate Cvec', 317 & wrk(klfree),lfree_local,iprtci) 318 319 else if(ci_program .eq. 'SIRIUS-CI')then 320 321 CALL TRACI2(WRK(KLSOT),NASHT,LSYM,NAEL,XNDXCI(KIASTR), 322 & XNDXCI(KTAIJ),XNDXCI(KTATO),XNDXCI(KTASYM), 323 & XNDXCI(KSTASA),XNDXCI(KSTBAA),NBEL, 324 & XNDXCI(KIBSTR),XNDXCI(KTBIJ), 325 & XNDXCI(KTBTO),XNDXCI(KTBSYM),XNDXCI(KSTASB), 326 & XNDXCI(KSTBAB),NASTR,NBSTR,XNDXCI(KIANNI), 327 & CVEC(1,IVEC),XNDXCI(KCOFF),WRK(KLC2), 328 & NOCTPA,XNDXCI(KNSSOA),XNDXCI(KISSOA), 329 & XNDXCI(KTPFSA),NOCTPB,XNDXCI(KNSSOB), 330 & XNDXCI(KISSOB), 331 & XNDXCI(KTPFSB),XNDXCI(KIOCOC),XNDXCI(KICOOS), 332 & XNDXCI(KCDTAS),MULD2H,MAXSYM,IPRTCI) 333! & XNDXCI(KCDTAS),MULD2H,MAXSYM,100) 334 335 end if 336#ifdef LUCI_DEBUG 337 write(lupri,*) ' after tra-ci run: ci-vec # =',ivec 338 call wrtmatmn(CVEC(1,IVEC),1,NCDET,1,NCDET,lupri) 339 write(lupri,*) ' after tra-ci run: sot-mat' 340 call wrtmatmn(WRK(KLSOT),1,n2ashx,1,n2ashx,lupri) 341#undef LUCI_DEBUG 342#endif 343 END IF 344 450 CONTINUE 345C 346 RETURN 347 END 348C /* Deck traci2 */ 349 SUBROUTINE TRACI2(T,NORB,ICSM, 350 & NAEL,IASTR,TAIJ,TATO,TASYM,NSTASA,ISTBAA, 351 & NBEL,IBSTR,TBIJ,TBTO,TBSYM,NSTASB,ISTBAB, 352 & NASTR,NBSTR, IANNI,C,ICOFF,C2, 353 & NOCTPA,NSTAOS,ISTAOS,ITPAST, 354 & NOCTPB,NSTBOS,ISTBOS,ITPBST, 355 & IOCOC,ICOOS,NCDTAS,SYMPRO,MAXSYM,NTEST ) 356C 357C COUNTER ROTATE CI - COEFFICIENTS CORRESPONDING 358C TO ORBITAL ROTATIONS DEFINED IN T 359C INPUT CI VECTOR IS C AND OUTPUT CI VECTOR OVERWRITERS C 360C C2 IS SCRATCH ABLE TO HOLD LARGEST SYMMETRY BLOCK 361C SYMMETRY OF CIN AND COUT ARE ASSUMED TO BE IDENTICAL AND EQUAL TO 362C ICSM 363C 364C d2h RAS Version , Jeppe Olsen November 1988 365C ( Debugged ! Error corrected 890809-hjaaj) 366c 367c New format of single excitations, October 1990 368c 369#include "implicit.h" 370#include "priunit.h" 371 INTEGER SYMPRO(8,8) 372 INTEGER TAIJ(*),TATO(*) 373 INTEGER TBIJ(*),TBTO(*) 374 INTEGER TASYM(MAXSYM+1,*),TBSYM(MAXSYM+1,*) 375 INTEGER NSTASA(MAXSYM),NSTASB(MAXSYM) 376 INTEGER ISTBAA(MAXSYM),ISTBAB(MAXSYM) 377 INTEGER ICOFF(MAXSYM),IANNI(NORB ** 2 ) 378 INTEGER IASTR(NAEL,NASTR),IBSTR(NBEL,NBSTR) 379 INTEGER NSTAOS(NOCTPA,*),ISTAOS(NOCTPA,*) 380 INTEGER NSTBOS(NOCTPB,*),ISTBOS(NOCTPB,*) 381 INTEGER IOCOC(NOCTPA,NOCTPB),ICOOS(NOCTPB,NOCTPA,MAXSYM) 382 INTEGER NCDTAS(*), ITPAST(*),ITPBST(*) 383C 384 DIMENSION T(*), C(*),C2(*) 385C LENGTH OF C2 MUST BE AT LEAST 386C 387 IF( NTEST .GE. 20 ) WRITE(LUPRI,*) ' --- INFORMATION FROM TRACI2' 388C 389C 390 DO 2000 IASM = 1, MAXSYM 391 IF( NTEST .GE. 25 ) WRITE(LUPRI,*) ' INFO FOR SYMMETRY... ',IASM 392 KDET = NCDTAS(IASM) 393 IBSM = SYMPRO(IASM,ICSM) 394 IOFF11 = ICOOS(1,1,IASM) 395C.. TRANSFORM THE BETA STRINGS 396 IF(NTEST.GE.25) WRITE(LUPRI,*) ' TRANSFORMATION OF BETA STRINGS' 397C.. TRANSFORM ORBITAL IORB 398 DO 900 IORB = 1, NORB 399 CALL SETVEC(C2,0.0D0,KDET) 400 IF(NTEST .GE. 25)WRITE(LUPRI,*)' *** ORBITAL TO BE ROTATED ', 401 & IORB 402 DO 890 IATP = 1, NOCTPA 403 DO 880 IBTP = 1, NOCTPB 404 IF( IOCOC(IATP,IBTP) .NE. 1 ) GOTO 880 405 IOFF = ICOOS (IBTP,IATP,IASM) 406C 407 NIA = NSTAOS(IATP,IASM) 408 NIB = NSTBOS(IBTP,IBSM) 409 IF ( NTEST .GE. 25 ) THEN 410 WRITE(LUPRI,*) ' INITIAL CI BLOCK ' 411 CALL WRTMAT(C(IOFF),NIA,NIB,NIA,NIB,0) 412 END IF 413C 414 IBSTRT = ISTBOS(IBTP,IBSM) 415 IASTRT = ISTAOS(IATP,IASM) 416 DO 800 IB = IBSTRT, IBSTRT+NIB-1 417 IOFFI = IOFF + (IB-IBSTRT)*NIA 418C IS ORBITAL IORB OCCUPIED IN IB ? 419 IOCC = 0 420 DO 750 IEL = 1, NBEL 421 IF(IBSTR(IEL,IB) .EQ. IORB ) IOCC = 1 422 750 CONTINUE 423C 424 IF ( IOCC .EQ. 0 ) THEN 425 CALL VECSUM(C2(IOFFI-IOFF11+1),C2(IOFFI-IOFF11+1), 426 & C(IOFFI),1.0D0,1.0D0,NIA) 427 ELSE 428 DO 700 IEX = TBSYM(1,IB),TBSYM(2,IB)-1 429 IJEX = TBIJ(IEX) 430 IF(IANNI(IJEX) .EQ. IORB ) THEN 431 JB = TBTO(IEX) 432 IF( JB .GT. NBSTR ) THEN 433 JB = JB - NBSTR 434 FACTOR = -T(IJEX) 435 ELSE 436 FACTOR = T(IJEX) 437 END IF 438C 439 JBTP = ITPBST(JB) 440 IF( IOCOC(IATP,JBTP) .NE. 1 ) GOTO 700 441 JBEFF = JB - ISTBOS(JBTP,IBSM)+1 442 IOFFJ = ICOOS(JBTP,IATP,IASM) - IOFF11 443 & + (JBEFF-1)*NIA + 1 444 CALL VECSUM(C2(IOFFJ),C2(IOFFJ), 445 & C(IOFFI),1.0D0,FACTOR,NIA) 446 END IF 447 700 CONTINUE 448 END IF 449 800 CONTINUE 450C 451 880 CONTINUE 452 890 CONTINUE 453 CALL COPVEC(C2,C(IOFF11),KDET) 454 IF( NTEST .GE. 50 ) THEN 455 WRITE(LUPRI,*) ' ROTATED CI BLOCK ' 456 CALL PRBLM2(C(IOFF11),NSTAOS(1,IASM),NOCTPA, 457 & NSTBOS(1,IBSM),NOCTPB,IOCOC,1) 458 END IF 459 900 CONTINUE 460 IF( NTEST .GE. 25 .AND. NTEST .LT. 50) THEN 461 WRITE(LUPRI,*) ' ROTATED CI BLOCK ' 462 CALL PRBLM2(C(IOFF11),NSTAOS(1,IASM),NOCTPA, 463 & NSTBOS(1,IBSM),NOCTPB,IOCOC,1) 464 END IF 465C 466C TPBLM2(A,AT,LBLR,NBLR,LBLC,NBLC,IFLAG,IORD) 467C PRBLM2(A,LBLR,NBLR,LBLC,NBLC,IFLAG,IORD) 468C.. BETA STRINGS HAVE NOW BEEN TRANSFORMED , TRANSFORM ALPHA QQ 469C STRINGS 470C.. TRANSPOSE TO EASE INDEXING 471 IF(NTEST .GE. 25)WRITE(LUPRI,*)' TRANSFORMATION OF ALPHA '// 472 & 'STRINGS' 473C CALL TRPMAT(C2,NIA,NIB,C(IOFF)) 474 CALL TPBLM2(C2,C(IOFF11),NSTAOS(1,IASM),NOCTPA, 475 & NSTBOS(1,IBSM),NOCTPB,IOCOC,1) 476 DO 1900 IORB = 1, NORB 477 IF( NTEST .GE. 25 )WRITE(LUPRI,*)' ORBITAL TO BE ROTATED ', 478 & IORB 479 CALL SETVEC(C2,0.0D0,KDET) 480 DO 1890 IATP = 1, NOCTPA 481 DO 1880 IBTP = 1, NOCTPB 482 IF( IOCOC(IATP,IBTP) .NE. 1 ) GOTO 1880 483 NIA = NSTAOS(IATP,IASM) 484 NIB = NSTBOS(IBTP,IBSM) 485 IF ( NTEST .GE. 25 ) THEN 486 WRITE(LUPRI,*) ' INITIAL CI BLOCK - alpha strings' 487 CALL WRTMAT(C(IOFF11),NIA,NIB,NIA,NIB,0) 488 END IF 489 IASTRT = ISTAOS(IATP,IASM) 490 IBSTRT = ISTBOS(IBTP,IBSM) 491 IOFF = ICOOS(IBTP,IATP,IASM) 492 IOFFR = IOFF - IOFF11 + 1 493 DO 1800 IA = IASTRT, IASTRT+NIA-1 494 IOFFI = IOFFR + (IA-IASTRT)*NIB 495C IS ORBITAL IORB OCCUPIED IN IA ? 496 IOCC = 0 497 DO 1750 IEL = 1, NAEL 498 IF(IASTR(IEL,IA) .EQ. IORB ) IOCC = 1 499 1750 CONTINUE 500C 501 IF ( IOCC .EQ. 0 ) THEN 502 CALL VECSUM(C2(IOFFI),C2(IOFFI), 503 & C(IOFFI+IOFF11-1),1.0D0,1.0D0,NIB) 504 ELSE 505 DO 1700 IEX = TASYM(1,IA),TASYM(2,IA)-1 506 IJEX = TAIJ(IEX) 507 IF(IANNI(IJEX) .EQ. IORB ) THEN 508 JA = TATO(IEX) 509 IF( JA .GT. NASTR ) THEN 510 JA = JA - NASTR 511 FACTOR = -T(IJEX) 512 ELSE 513 FACTOR = T(IJEX) 514 END IF 515 JATP = ITPAST(JA) 516C 517 IF(IOCOC(JATP,IBTP) .NE. 1 ) GOTO 1700 518C 519 JASTRT = ISTAOS(JATP,IASM) 520 IOFFJ = ICOOS(IBTP,JATP,IASM)-IOFF11 521 & + (JA-JASTRT)*NIB + 1 522 CALL VECSUM(C2(IOFFJ),C2(IOFFJ),C(IOFFI+IOFF11-1), 523 & 1.0D0,FACTOR,NIB) 524C 525 END IF 526 1700 CONTINUE 527 END IF 528 1800 CONTINUE 529 1880 CONTINUE 530 1890 CONTINUE 531 CALL COPVEC(C2,C(IOFF11),KDET) 532 1900 CONTINUE 533C CALL TRPMAT(C2,NIB,NIA,C(IOFF) ) 534 CALL TPBLM2(C2,C(IOFF11),NSTAOS(1,IASM),NOCTPA, 535 & NSTBOS(1,IBSM),NOCTPB,IOCOC,2) 536 IF( NTEST .GE. 20 ) THEN 537 WRITE(LUPRI,*) ' ROTATED CI BLOCK FOR IASM ...', IASM 538 CALL PRBLM2(C(IOFF11),NSTBOS(1,IBSM),NOCTPB, 539 & NSTAOS(1,IBSM),NOCTPB,IOCOC,1) 540 END IF 541 2000 CONTINUE 542C 543 RETURN 544 END 545C /* Deck typstr */ 546 SUBROUTINE TYPSTR(STRING,NEL,NSTRIN,ITYPE,IAB,NTEST) 547C 548C OCCUPATION TYPE OF STRINGS 549C 550#include "implicit.h" 551#include "priunit.h" 552 INTEGER STRING(NEL,NSTRIN),ITYPE(NSTRIN) 553C 554 DO 100 ISTRIN = 1, NSTRIN 555 ITYPE(ISTRIN) = IOCTYP(STRING(1,ISTRIN),IAB,1) 556 100 CONTINUE 557C 558 IF ( NTEST .GT. 8 ) THEN 559 WRITE(LUPRI,*) ' TYPSTR: OCCUPATION TYPES OF STRINGS ' 560 CALL IWRTMA(ITYPE,1,NSTRIN,1,NSTRIN) 561 END IF 562C 563 RETURN 564 END 565C /* Deck cidiag */ 566 SUBROUTINE CIDIAG(ICSYM,NOH2,FCAC,H2AC,XNDXCI,DIAGC,WRK,LFREE) 567C 568C Written 18-jan-1988 J.O. ( after cidiag of hjaaj ) 569C NOH2 parameter 900717-hjaaj 570C 571C Part of routines for determinant based ci 572C 573C Purpose: 574C 575C Calculate the diagonal of the CI matrix and, 576C if LUIT2 .gt. 0, save it on LUIT2. 577C 578C 579! module dependencies 580 use lucita_mcscf_ci_interface_procedures 581 use mcscf_or_gasci_2_define_cfg, only : 582 & define_lucita_cfg_dynamic 583 584#include "implicit.h" 585#include "priunit.h" 586 DIMENSION FCAC(*),H2AC(*),DIAGC(*) 587 DIMENSION XNDXCI(*), WRK(LFREE) 588 LOGICAL NOH2, FNDLAB 589C 590C Used from common blocks: 591C INFINP : FLAG(*) 592C INFORB : N2ASHX 593C INFVAR : NCONF 594C INFTAP : LUIT2 595C INFPRI : IPRDIA 596C CIINFO : NDTASM(*) 597C SPINFO : ? 598C CSFBAS : CSF core allocation for XNDXCI 599C 600#include "maxorb.h" 601#include "dummy.h" 602#include "infinp.h" 603#include "inforb.h" 604#include "infvar.h" 605#include "inftap.h" 606#include "infpri.h" 607#include "ciinfo.h" 608#include "spinfo.h" 609#include "csfbas.h" 610C 611C *** local variables: 612 LOGICAL CSFEXP 613 CHARACTER*8 TABLE(4) 614 DATA TABLE/'********','CIDIAG1 ','CIDIAG2 ', 615 & 'EODATA '/ 616C 617 CSFEXP = (.NOT.FLAG(27) .AND. ICSYM .EQ. LSYM ) 618C ... 980820-hjaaj: only CSF for ref.sym. LSYM, 619C determinants for other symmetries 620C 621C *** Construct diagonal of CI matrix 622 if(ci_program .eq. 'LUCITA ')then 623 KW1 = 1 624 LW1 = LFREE - KW1 625 call define_lucita_cfg_dynamic(icsym, ! rhs symmetry 626 & icsym, ! lhs symmetry 627 & 0, ! no istaci 628 & 0, ! spin for operator 629 & 0, ! spin for operator 630 & 1, ! one root 631 & ispin, ! singlet, doublet, triplet, ... 632 & -1, ! # of davidson iterations 633 & -1, ! 1/2-particle density calcs 634 & iprdia, ! print level 635 & 1.0d-10, ! screening factor 636 & 0.0d0, ! davidson ci convergence 637 & 0.0d0, ! davidson ci convergence (auxiliary) 638 & .false., ! do not calculate nat. orb. occ. num. (internally inside LUCITA) 639 & .false., ! wave function analysis 640 & NOH2, ! no 2-el operators active 641 & .true., ! integrals provided by / return density matrices to the MCSCF environment 642 & .false., ! calculate 2-particle density matrix 643 & .false., ! calculate transition density matrix 644 & -1, ! vector exchange type 1 645 & -1, ! vector exchange type 2 646 & .false., ! vector exchange active in parallel runs in mc2lu interface (both mc-->lu and lu-->mc) 647 & .false.) ! io2io vector exchange between mcscf and lucita: default mem2io/io2mem 648 649 call dzero(diagc,nconf) 650 call mcscf_lucita_interface(diagc,vdummy,fcac,h2ac, 651 & 'return CIdia',wrk(kw1),lw1, 652 & iprdia) 653!#define LUCI_DEBUG 654#ifdef LUCI_DEBUG 655 write(lupri,*) ' after diag run: hdiag' 656 call wrtmatmn(diagc,1,nconf,1,nconf,lupri) 657#undef LUCI_DEBUG 658#endif 659 else if(ci_program .eq. 'SIRIUS-CI')then 660 KFIJ = 1 661 KFJI = KFIJ + N2ASHX 662 KW1 = KFJI + N2ASHX 663 LW1 = LFREE - KW1 664 IF (NOH2) THEN 665 CALL DZERO(WRK(KFIJ),2*N2ASHX) 666 ELSE 667 CALL GETFIJ(WRK(KFIJ),WRK(KFJI),H2AC) 668 END IF 669 IF (.NOT.CSFEXP) THEN 670 CALL DZERO (DIAGC,NCONF) 671 CALL CIDIAD(DIAGC,FCAC,WRK(KFIJ),WRK(KFJI),XNDXCI,ICSYM, 672 * WRK(KW1),LW1) 673#ifdef LUCI_DEBUG 674 write(lupri,*) ' after diag run: hdiag' 675 call wrtmatmn(diagc,1,nconf,1,nconf,lupri) 676#undef LUCI_DEBUG 677#endif 678 ELSE 679C 680C.. Change to CSF format 681C 682 KDETDG = KW1 683 KW1 = KDETDG + NDTASM(ICSYM) 684 LW1 = LFREE - KW1 685 CALL CIDIAD(WRK(KDETDG),FCAC,WRK(KFIJ),WRK(KFJI),XNDXCI, 686 & ICSYM,WRK(KW1),LW1) 687 CALL CSDIAG(DIAGC,WRK(KDETDG),NCNFTP(1,ICSYM),NTYP, 688 & XNDXCI(KICTS(1)),NDTFTP,NCSFTP,0,IPRDIA) 689C CSDIAG(CSFDIA,DETDIA,NCNFTP,NTYP, 690C & ICTSDT,NDTFTP,NCSFTP,IFLAG,NTEST) 691 END IF 692 end if ! ci_program switch 693C 694 IF (LUIT2 .GT. 0) THEN 695C 696C *** Write diagonal of CI matrix on LUIT2 (label: 'CIDIAG2 ' ) 697C 698 REWIND LUIT2 699 IF( FNDLAB(TABLE(4),LUIT2) ) THEN 700 BACKSPACE LUIT2 701 ELSE 702 REWIND LUIT2 703 END IF 704 WRITE (LUIT2) TABLE(1),TABLE(1),TABLE(1),TABLE(3) 705 NC4 = MAX(4,NCONF) 706 CALL WRITT(LUIT2,NC4,DIAGC) 707 WRITE (LUIT2) TABLE(1),TABLE(1),TABLE(1),TABLE(4) 708 END IF 709 IF (IPRDIA.GT.0 .OR. IPRSIR .GT. 20) THEN 710 write(lupri,*) 'CIDIAG LUJIT2',LUIT2 711 !call qdump(lupri) 712 write(lupri,'(/A)') 'CIDIAG: first 20 elements of CI diagonal' 713 NPRINT = MIN(20,NCONF) 714 write(lupri,'(I5,F20.10)') (I, DIAGC(I),I=1,NPRINT) 715 END IF 716C 717C *** End of subroutine CIDIAG 718C 719 RETURN 720 END 721C /* Deck cidiad */ 722 SUBROUTINE CIDIAD(DIAG,FCAC,FIJ,FJI,XNDXCI,ICSYM,WRK,LFREE) 723C 724C *** CONSTRUCT DIAGONAL PART OF CI MATRIX 725C 726C SOME INPUT 727C FCAC : inactive Fock matrix 728C (modified one body integrals) 729C FIJ : (II/JJ) integrals 730C FJI : (IJ/IJ) integrals 731C 732#include "implicit.h" 733 DIMENSION DIAG(*), FCAC(*), FIJ(*), FJI(*), WRK(LFREE) 734 DIMENSION XNDXCI(*) 735#include "priunit.h" 736C 737C 738C Used from common blocks: 739C INFINP : 740C INFORB : MULD2H(8,8),NASHT,N2ASHX 741C INFIND : NSM(NASHT) 742C INFPRI : IPRDIA 743C 744C DETBAS : core allocation 745C STRNUM : NAEL,NBEL,NASTR,NBSTR,? 746C CIINFO : NDTASM(*),ICOMBI 747C 748#include "maxash.h" 749#include "maxorb.h" 750#include "infinp.h" 751#include "inforb.h" 752#include "infind.h" 753#include "infpri.h" 754C 755#include "mxpdim.h" 756#include "detbas.h" 757#include "strnum.h" 758#include "ciinfo.h" 759C 760 CALL GETTIM(T0,W0) 761 ECORE = 0.0D0 762c Store CI offsets in C arrays 763 CALL CIOFF(ICSYM,1,XNDXCI,IPRDIA) 764C CALL CIOFF(IREFSM,ICORHC,XNDXCI,NTEST) 765C 766C 767 KFREE = 1 768 CALL MEMADD(KXA, NASHT, KFREE,2) 769 CALL MEMADD(KXB, NASHT, KFREE,2) 770 CALL MEMADD(KSCR, 2*NASHT,KFREE,2) 771 CALL MEMADD(KH1DIA,NASHT, KFREE,2) 772 CALL MEMADD(KRJ, N2ASHX, KFREE,2) 773 CALL MEMADD(KRK, N2ASHX, KFREE,2) 774 CALL MEMADD(KLSCR, N2ASHX, KFREE,2) 775 IF (IPRDIA .GT. 15) THEN 776 WRITE (LUPRI,'(/2A/8I8)') ' CIDIAD : local pointers', 777 * ' KXA KXB KSCR KH1DIA KRJ KRK KFREE LFREE', 778 * KXA,KXB,KSCR,KH1DIA,KRJ,KRK,KFREE,LFREE 779 WRITE (LUPRI,*) ' ICSYM IN CIDIAD ', ICSYM 780 END IF 781 IF (KFREE .GT. LFREE) CALL ERRWRK('CIDIAD',KFREE,LFREE) 782C 783C Pack information as wanted by CIDIA4 784C 785 DO 100 I = 1,NASHT 786 WRK(KXA-1+I) = FCAC( I*(I+1)/2 ) 787 100 CONTINUE 788 CALL REORMT(WRK(KXA),WRK(KH1DIA),NASHT,1,XNDXCI(KLTSOB),1) 789C ... diagonal of modified one-body integrals 790C 791 CALL REORMT(FIJ,WRK(KRJ),NASHT,NASHT,XNDXCI(KLTSOB), 792 & XNDXCI(KLTSOB) ) 793 CALL REORMT(FJI,WRK(KRK),NASHT,NASHT,XNDXCI(KLTSOB), 794 & XNDXCI(KLTSOB) ) 795C ... Coulomb and exchange integrals, FIJ and FJI 796C 797C 798C 799 IF ( IPRDIA .GE. 15 ) THEN 800 WRITE(LUPRI,'(/A)') 801 * ' CIDIAD : Diagonal of FCAC, the modified one body matrix' 802 WRITE(LUPRI,'(5F15.8)') (WRK(KH1DIA-1+I),I=1,NASHT) 803 WRITE(LUPRI,'(/A)') 804 * ' CIDIAD : Coulomb integrals RJ(u,v) = (uu/vv)' 805 CALL OUTPUT(WRK(KRJ),1,NASHT,1,NASHT,NASHT,NASHT,1,LUPRI) 806 WRITE(LUPRI,'(/A)') 807 * ' CIDIAD : Exchange integrals RK(u,v) = (uv/uv)' 808 CALL OUTPUT(WRK(KRK),1,NASHT,1,NASHT,NASHT,NASHT,1,LUPRI) 809 END IF 810C 811C 812 NCDET = NDTASM(ICSYM) 813 CALL GETTIM(T1,W1) 814 CALL CIDIA4(NAEL,NASTR,XNDXCI(KIASTR), 815 & NBEL,NBSTR,XNDXCI(KIBSTR), 816 & NASHT,DIAG,NCDET,XNDXCI(KSTASA),XNDXCI(KSTASB), 817 & MAXSYM,WRK(KH1DIA), 818 & XNDXCI(KSTBAA),XNDXCI(KSTBAB), 819 & ICSYM,WRK(KXA),WRK(KXB),WRK(KSCR),WRK(KRJ), 820 & WRK(KRK),MULD2H,XNDXCI(KNSSOA),XNDXCI(KNSSOB), 821 & XNDXCI(KIOCOC),NOCTPA,NOCTPB,XNDXCI(KISSOA), 822 & XNDXCI(KISSOB),ECORE,XNDXCI(KICOOS),IPRDIA ) 823C CALL CIDIA4(NAEL,NASTR,IASTR,NBEL,NBSTR,IBSTR, 824C & NORB,DIAG,NDET,NSTASA,NSTASB, 825C & MAXSYM,H,ISTBAA,ISTBAB, 826C & ICSYM,XA,XB,SCR,RJ,RK, 827C & SYMPRO,NSSOA,NSSOB,IOCOC,NOCTPA,NOCTPB, 828C & ISSOA,ISSOB,ECORE,ICOOS,NTEST) 829 CALL GETTIM(T2,W2) 830 IF (IPRDIA .GE. 1) WRITE (LUPRI,'(/A,2I10/A,2F10.3)') 831 * ' CIDIAD.cidia4 : KFREE, LFREE = ',KFREE,LFREE, 832 * ' CPU and wall time (sec) :',T2-T1,W2-W1 833C 834C... end of cidiad so : 835 RETURN 836 END 837C /* Deck cisigc */ 838 SUBROUTINE CISIGC(NSIM,BCVECS,SCVECS,LSCVEC, 839 * FCAC,H2AC,XNDXCI,WRK,LFREE) 840C 841C Written 18-JAN-1988 by J.O. 842C ( from hjaajs cisigc ) 843C 844C Purpose: 845C Compute CI sigma vector(s) for direct CI, 846C the inactive energy EMY is not included. 847C calls lunar determinant routines 848C 849C Output: 850C SCVECS; CI sigma vector(s) 851C 852C Scratch: 853C WRK(LFREE) 854C 855C 856! module dependencies 857 use lucita_mcscf_ci_cfg 858#include "implicit.h" 859 DIMENSION BCVECS(NCONST,*),FCAC(*),H2AC(*),SCVECS(LSCVEC,*) 860 DIMENSION XNDXCI(*), WRK(LFREE) 861 LOGICAL NOH2, IH8SM 862 PARAMETER ( NOH2 = .FALSE. , IH8SM = .TRUE. , D0 = 0.0D0 ) 863C NOH2 ... 2-electron part to be included 864C IH8SM ... integrals have 8-fold symmetry 865C 866C Used from common blocks: 867C INFORB : NASHT, N2ASHX 868C INFLIN : LSYMST,NCONST (set in sirset.F) 869C INFTIM : NCALLS,TIMCPU,TIMWAL ! IDTIM is index for these 870C 871#include "maxorb.h" 872#include "infinp.h" 873#include "inforb.h" 874#include "inflin.h" 875#include "infpri.h" 876C --- 877 PARAMETER (IDTIM = 2) 878#include "inftim.h" 879#include "priunit.h" 880C 881!#define LUCI_DEBUG 882#ifdef LUCI_DEBUG 883 real(8), allocatable :: btmp(:) 884#endif 885C 886C *** SCVECS(I,K) = SUM(J) OF L_CI(I,J)*BCVECS(J,K) 887C 888C 889 CALL GETTIM(T0,W0) 890C 891C ** One- and two-electron contributions 892C 893 CALL GETTIM(T1,W1) 894 CALL GETTIM(T2,W2) 895 896 ISPIN1 = 0 897 ISPIN2 = 0 898C ... singlet-singlet coupling of 2-electron integrals 899 900 KUFCAC = 1 901 KW1 = KUFCAC + N2ASHX 902 LW1 = LFREE - KW1 903! Unpack FCAC for SIRIUS-CI 904 if(ci_program .eq. 'SIRIUS-CI')then 905 CALL DSPTSI(NASHT,FCAC,WRK(KUFCAC)) 906! ... unpack FCAC using CALL DSPTSI(N,ASP,ASI) 907 else if(ci_program .eq. 'LUCITA ')then 908 call dcopy(NNASHX,FCAC,1,WRK(KUFCAC),1) 909! 910! set flag for CI trial vector in MCSCF (needed for LUCITA interface) 911 mcscf_ci_trial_vector = .true. 912 913! set vector exchange type: cref 914 if(cref_is_active_bvec_for_sigma)then 915 vector_exchange_type1 = 2 916 cref_is_active_bvec_for_sigma = .false. 917 else 918! set vector exchange type: bvec trial vector 919 vector_exchange_type1 = 2 920 end if 921 end if 922 923#ifdef LUCI_DEBUG 924 write(lupri,*) ' mcscf_ci_trial_vector, vector_exchange_type1 '// 925 & ' cref_is_active_bvec_for_sigma , nsim , noh2', 926 & mcscf_ci_trial_vector, vector_exchange_type1, 927 & cref_is_active_bvec_for_sigma , nsim, noh2 928 allocate (btmp(nconst)) 929#endif 930 931 DO isim = 1,nsim ! loop over sigma vectors to calculate 932 933 call dzero(scvecs(1,isim),nconst) 934 935#ifdef LUCI_DEBUG 936 btmp(1:nconst) = bcvecs(1:nconst,isim) 937 do i = 1, nconst 938 if (abs(bcvecs(i,isim)) .gt. 0.2d0 ) then 939 write(lupri,*) 'big B in',i, bcvecs(i,isim) 940 end if 941 end do 942#endif 943 944 CALL CISIGD(LSYMST,LSYMST,NCONST,NCONST, 945 & BCVECS(1,ISIM),SCVECS(1,ISIM), WRK(KUFCAC),H2AC, 946 & NOH2,IH8SM,XNDXCI,ISPIN1,ISPIN2,WRK(KW1),LW1) 947C CALL CISIGD(ICSYM,IHCSYM,NCDET,NHCDET, C,HC, FCAC,H2AC, 948C & NOH2,IH8SM,XNDXCI,ISPIN1,ISPIN2,WRK,LFREE) 949 950#ifdef LUCI_DEBUG 951 btmp(1:nconst) = btmp(1:nconst) - bcvecs(1:nconst,isim) 952 do i = 1, nconst 953 if (abs(btmp(i)) .gt. 1.0D-14) then 954 write(lupri,*) 'B changed ',i, bcvecs(i,isim), btmp(i) 955 end if 956 if (abs(bcvecs(i,isim)) .gt. 0.2d0 ) then 957 write(lupri,*) 'big B ',i, bcvecs(i,isim), scvecs(i,isim) 958 end if 959 end do 960#endif 961 end do 962 963#ifdef LUCI_DEBUG 964 deallocate(btmp) 965#undef LUCI_DEBUG 966#endif 967C 968! reset flag for CI trial vector in MCSCF (needed for LUCITA interface) 969 mcscf_ci_trial_vector = .false. 970 971 CALL GETTIM(T3,W3) 972C 973C 974C Acuumulate timing for TIMOPT 975C 976 NCALLS( IDTIM) = NCALLS( IDTIM) + 1 977 NVECS ( IDTIM) = NVECS ( IDTIM) + NSIM 978 TIMCPU(1,IDTIM) = TIMCPU(1,IDTIM) + T3 - T0 979 TIMCPU(2,IDTIM) = TIMCPU(2,IDTIM) + T1 - T0 980 TIMCPU(3,IDTIM) = TIMCPU(3,IDTIM) + T2 - T1 981 TIMCPU(4,IDTIM) = TIMCPU(4,IDTIM) + T3 - T2 982 TIMWAL(1,IDTIM) = TIMWAL(1,IDTIM) + W3 - W0 983 TIMWAL(2,IDTIM) = TIMWAL(2,IDTIM) + W1 - W0 984 TIMWAL(3,IDTIM) = TIMWAL(3,IDTIM) + W2 - W1 985 TIMWAL(4,IDTIM) = TIMWAL(4,IDTIM) + W3 - W2 986C 987C 988C 989C *** End of subroutine CISIGC 990C 991 RETURN 992 END 993C /* Deck ciprp */ 994 SUBROUTINE CIPRP(NSIM,CREF,SCVECS,LSCVEC,PRPAC,XNDXCI,WRK,LFREE) 995C 996C Written Dec 1989 hjaaj 997C 998C Purpose: 999C Compute CI sigma vector(s) for NSIM one-electron operators 1000C SCVECS(I) = SUM(J) <I|PRP|J> * CREF(J) 1001C 1002C calls lunar determinant routines 1003C 1004C Output: 1005C SCVECS; CI sigma vector(s) 1006C 1007C Scratch: 1008C WRK(LFREE) 1009C 1010C 1011 use lucita_mcscf_ci_cfg 1012#include "implicit.h" 1013 DIMENSION CREF(NCONRF),PRPAC(N2ASHX,NSIM),SCVECS(LSCVEC,NSIM) 1014 DIMENSION XNDXCI(*), WRK(LFREE) 1015 LOGICAL NOH2, IH8SM 1016 PARAMETER ( NOH2 = .TRUE. , IH8SM = .FALSE. , D2 = 2.0D0 ) 1017#include "thrzer.h" 1018C NOH2 ... No 2-electron part 1019C IH8SM ... integrals do not have 8-fold symmetry 1020C 1021C Used from common blocks: 1022C INFORB : N2ASHX 1023C INFLIN : LSYMRF,NCONRF,LSYMST,NCONST 1024C 1025#include "maxorb.h" 1026#include "inforb.h" 1027#include "inflin.h" 1028#include "priunit.h" 1029C 1030C 1031C *** SCVECS(I) = SUM(J) <I|PRP|J> * CREF(J) 1032C 1033 KW1 = 1 1034 LW1 = LFREE - KW1 1035C 1036C Set spin couplings 1037C 1038 ISPIN1 = 0 1039 ISPIN2 = 0 1040! 1041! set flag for CI trial vector in MCSCF (needed for LUCITA interface) 1042 mcscf_ci_trial_vector = .true. 1043 1044! set vector exchange type: cref 1045 if(cref_is_active_bvec_for_sigma)then 1046 vector_exchange_type1 = 2 1047 cref_is_active_bvec_for_sigma = .false. 1048 else 1049! set vector exchange type: bvec trial vector 1050 vector_exchange_type1 = 2 1051 end if 1052 1053C ... singlet-singlet coupling of 2-electron integrals 1054 DO 100 ISIM = 1, NSIM 1055 CALL DZERO(SCVECS(1,ISIM),NCONST) 1056 CALL CISIGD(LSYMRF,LSYMST,NCONRF,NCONST, CREF,SCVECS(1,ISIM), 1057 & PRPAC(1,ISIM),DUMMY,NOH2,IH8SM, 1058 & XNDXCI,ISPIN1,ISPIN2,WRK(KW1),LW1) 1059C CALL CISIGD(ICSYM,IHCSYM,NCDET,NHCDET, C,HC, FCAC,H2AC, 1060C & NOH2,IH8SM,XNDXCI,ISPIN1,ISPIN2,WRK,LFREE) 1061 100 CONTINUE 1062C 1063! set flag for CI trial vector in MCSCF (needed for LUCITA interface) 1064 mcscf_ci_trial_vector = .false. 1065C 1066C *** End of subroutine CIPRP 1067C 1068 RETURN 1069 END 1070C /* Deck cirpr_triplet */ 1071 SUBROUTINE CIPRP_TRIPLET(NSIM,CREF,SCVECS,LSCVEC, 1072 & PRPAC,XNDXCI,WRK,LFREE) 1073C 1074C hjaaj Nov 2015: triplet version of subroutine CIPRP 1075C 1076C Purpose: 1077C Compute CI sigma vector(s) for NSIM one-electron operators 1078C SCVECS(I) = SUM(J) <I|PRP(triplet)|J> * CREF(J) 1079C 1080C calls lunar determinant routines 1081C 1082C Output: 1083C SCVECS; CI sigma vector(s) 1084C 1085C Scratch: 1086C WRK(LFREE) 1087C 1088C 1089 use lucita_mcscf_ci_cfg 1090#include "implicit.h" 1091 DIMENSION CREF(NCONRF),PRPAC(N2ASHX,NSIM),SCVECS(LSCVEC,NSIM) 1092 DIMENSION XNDXCI(*), WRK(LFREE) 1093 LOGICAL NOH2, IH8SM 1094 PARAMETER ( NOH2 = .TRUE. , IH8SM = .FALSE. , D2 = 2.0D0 ) 1095#include "thrzer.h" 1096C NOH2 ... No 2-electron part 1097C IH8SM ... integrals do not have 8-fold symmetry 1098C 1099C Used from common blocks: 1100C INFORB : N2ASHX 1101C INFLIN : LSYMRF,NCONRF,LSYMST,NCONST 1102C 1103#include "maxorb.h" 1104#include "inforb.h" 1105#include "inflin.h" 1106#include "priunit.h" 1107C 1108C 1109C *** SCVECS(I) = SUM(J) <I|PRP|J> * CREF(J) 1110C 1111 KW1 = 1 1112 LW1 = LFREE - KW1 1113C 1114C Set spin couplings 1115C 1116 ISPIN1 = 1 1117 ISPIN2 = 0 1118! 1119! set flag for CI trial vector in MCSCF (needed for LUCITA interface) 1120 mcscf_ci_trial_vector = .true. 1121 1122! set vector exchange type: cref 1123 if(cref_is_active_bvec_for_sigma)then 1124 vector_exchange_type1 = 2 1125 cref_is_active_bvec_for_sigma = .false. 1126 else 1127! set vector exchange type: bvec trial vector 1128 vector_exchange_type1 = 2 1129 end if 1130 1131C ... singlet-singlet coupling of 2-electron integrals 1132 DO 100 ISIM = 1, NSIM 1133 CALL DZERO(SCVECS(1,ISIM),NCONST) 1134 CALL CISIGD(LSYMRF,LSYMST,NCONRF,NCONST, CREF,SCVECS(1,ISIM), 1135 & PRPAC(1,ISIM),DUMMY,NOH2,IH8SM, 1136 & XNDXCI,ISPIN1,ISPIN2,WRK(KW1),LW1) 1137C CALL CISIGD(ICSYM,IHCSYM,NCDET,NHCDET, C,HC, FCAC,H2AC, 1138C & NOH2,IH8SM,XNDXCI,ISPIN1,ISPIN2,WRK,LFREE) 1139 100 CONTINUE 1140C 1141! set flag for CI trial vector in MCSCF (needed for LUCITA interface) 1142 mcscf_ci_trial_vector = .false. 1143C 1144C *** End of subroutine CIPRP_TRIPLET 1145C 1146 RETURN 1147 END 1148C /* Deck cisigd */ 1149 SUBROUTINE CISIGD(ICSYM,IHCSYM,NCVEC,NHCVEC,C,HC, 1150 * FCAC,H2AC,ONLYH1,IH8SM, 1151 * XNDXCI,JSPIN1,JSPIN2,WRK,LFREE) 1152C 1153C MASTER ROUTINE FOR DIRECT CI CALCULATION 1154C 1155C MOTECC-90: The algorithms used in this module, CISIGD, 1156C are described in Chapter 8 Section D.3 of 1157C MOTECC-90 "Direct CI for RAS Expansions" 1158C 1159C 1160C SOME INPUT : 1161C C : CI vector, of length NCVEC 1162C FCAC : ONE-ELECTRON HAMILTONIAN MODIFIED FOR CORE ELECTRONS 1163C *** FULL MATRIX *** 1164C 1165C XNDXCI: ARRAY CONTAINING STRING INFORMATION 1166C AS OBTAINED IN DETINF 1167C OUTPUT : 1168C HC : H TIMES C, OF SYMMETRY IHCSYM AND LENGTH NHCVEC 1169C 1170C 1171! module dependencies 1172 use lucita_mcscf_ci_cfg 1173 use lucita_mcscf_ci_interface_procedures 1174 use mcscf_or_gasci_2_define_cfg, only : 1175 & define_lucita_cfg_dynamic 1176#include "implicit.h" 1177 DIMENSION C(*), HC(*) 1178 DIMENSION FCAC(*), H2AC(*), WRK(LFREE), XNDXCI(*) 1179 LOGICAL ONLYH1,IH8SM 1180C 1181C THRSML : THRESHOLD FOR INTEGRALS AND CI COEFFICIENTS TO BE ZAPPED 1182C 1183 PARAMETER ( D1 = 1.0D0, THRSML = 1.0D-14 ) 1184C 1185#include "priunit.h" 1186C 1187C Used from common blocks: 1188C INFINP : FLAG(27), LSYM,ISPIN 1189C INFPRI : IPRSIG 1190C 1191C CIINFO : NDTASM(*) 1192C CSFBAS : CSF core allocation for XNDXCI 1193C CBESPN : ISPIN1, ISPIN2 1194C 1195#include "maxorb.h" 1196#include "infinp.h" 1197#include "infpri.h" 1198C 1199#include "ciinfo.h" 1200#include "csfbas.h" 1201#include "cbespn.h" 1202C 1203C 1204 LOGICAL CSFEXP 1205 integer :: xctype1, xctype2 1206C 1207C Transfer spin-couplings to CBESPN 1208C 1209 ISPIN1 = JSPIN1 1210 ISPIN2 = JSPIN2 1211C 1212C 1213C 980820-hjaaj: now use CSF for LSYM and det.s for other sym.s 1214C 1215 CSFEXP = .NOT.FLAG(27) .AND. 1216 & (ICSYM .EQ. LSYM .OR. IHCSYM .EQ. LSYM) 1217C 1218 IF (CSFEXP) THEN 1219 IERR = 0 1220 ICCSF = 0 1221 IF ( ICSYM.EQ.LSYM ) THEN 1222C ABACUS and RESPONSE use CSF for singlet and det for triplet, 1223C we find what is the case for this call by comparing 1224C to NCSASM and NDTASM: 1225 IF ( NCVEC .EQ. NCSASM( ICSYM) ) THEN 1226 ICCSF = 1 1227 ELSE IF ( NCVEC .NE. NDTASM( ICSYM) ) THEN 1228 IERR=IERR+1 1229 END IF 1230 ELSE 1231 IF ( NCVEC .NE. NDTASM( ICSYM) ) IERR=IERR+1 1232 END IF 1233 IHCCSF = 0 1234 IF (IHCSYM.EQ.LSYM ) THEN 1235C ABACUS and RESPONSE use CSF for singlet and det for triplet, 1236C we find what is the case for this call by comparing 1237C to NCSASM and NDTASM: 1238 IF (NHCVEC .EQ. NCSASM(IHCSYM) ) THEN 1239 IHCCSF = 1 1240 ELSE IF (NHCVEC .NE. NDTASM(IHCSYM) ) THEN 1241 IERR=IERR+1 1242 END IF 1243 ELSE 1244 IF (NHCVEC .NE. NDTASM(IHCSYM) ) IERR=IERR+1 1245 END IF 1246 IF (IERR .GT. 0) THEN 1247 WRITE (LUPRI,*) 'CSF ERROR in CISIGD, LSYM =',LSYM 1248 WRITE (LUPRI,*)' ISPIN, ISPIN1, ISPIN2:', 1249 & ISPIN,ISPIN1,ISPIN2 1250 WRITE (LUPRI,*) 1251 * ' ICSYM, NCVEC, NCSASM( ICSYM), NDTASM( ICSYM):', 1252 * ICSYM, NCVEC, NCSASM( ICSYM), NDTASM( ICSYM) 1253 WRITE (LUPRI,*) 1254 * 'IHCSYM, NHCVEC, NCSASM(IHCSYM), NDTASM(IHCSYM):', 1255 * IHCSYM, NHCVEC, NCSASM(IHCSYM), NDTASM(IHCSYM) 1256 CALL QUIT('NCVEC/NHCVEC disagree with NCSASM(:) in cisigd') 1257 END IF 1258 IF (ICCSF .EQ. 0 .AND. IHCCSF .EQ. 0) CSFEXP = .FALSE. 1259 END IF 1260C 1261 IF (CSFEXP) THEN 1262C .. change from CSf basis to determinant basis 1263 NCDET = NDTASM(ICSYM) 1264 NHCDET = NDTASM(IHCSYM) 1265 NDET = MAX(NCDET,NHCDET) 1266 KFREE = 1 1267 CALL MEMADD(KCDET,NDET,KFREE,2) 1268 CALL MEMADD(KHCDET,NDET,KFREE,2) 1269 IF (ICCSF .EQ. 1) THEN 1270 CALL COPVEC(C,WRK(KHCDET),NCVEC) 1271 CALL CSDTVC(WRK(KHCDET),WRK(KCDET),1,XNDXCI(KDTOC), 1272 & XNDXCI(KICTS(1)),ICSYM,0,IPRSIG) 1273 ELSE 1274 CALL COPVEC(C,WRK(KCDET),NCVEC) 1275 END IF 1276 CALL SETVEC(WRK(KHCDET),0.0D0,NHCDET) 1277 LW1 = LFREE - KFREE 1278 CALL CISGD2(ICSYM,IHCSYM,NCDET,NHCDET,WRK(KCDET),WRK(KHCDET), 1279 * FCAC,H2AC,ONLYH1,IH8SM,XNDXCI,WRK(KFREE),LW1) 1280 IF (IHCCSF .EQ. 1) THEN 1281 CALL CSDTVC(WRK(KCDET),WRK(KHCDET),2,XNDXCI(KDTOC), 1282 & XNDXCI(KICTS(1)),IHCSYM,0,IPRSIG) 1283 CALL DAXPY(NHCVEC,D1,WRK(KCDET),1,HC,1) 1284 ELSE 1285 CALL DAXPY(NHCVEC,D1,WRK(KHCDET),1,HC,1) 1286 END IF 1287 ELSE 1288 IF (NCVEC .NE. NDTASM( ICSYM) .OR. 1289 * NHCVEC .NE. NDTASM(IHCSYM)) THEN 1290 WRITE (LUPRI,*) 'ERROR in CISIGD' 1291 WRITE (LUPRI,*) ' ICSYM, NCVEC, NDTASM( ICSYM):', 1292 * ICSYM, NCVEC, NDTASM( ICSYM) 1293 WRITE (LUPRI,*) 'IHCSYM, NHCVEC, NDTASM(IHCSYM):', 1294 * IHCSYM, NHCVEC, NDTASM(IHCSYM) 1295 CALL QUIT('NCVEC/NHCVEC disagree with NDTASM(:) in cisigd') 1296 END IF 1297 1298 if(ci_program .eq. 'LUCITA ')then 1299 1300! vector_exchange_type1 is set in calling upper level routine 1301 xctype1 = vector_exchange_type1 1302 xctype2 = -1 1303 1304 call define_lucita_cfg_dynamic(icsym, ! rhs symmetry 1305 & ihcsym, ! lhs symmetry 1306 & 0, ! no istaci 1307 & ispin1, ! spin for operator 1308 & ispin2, ! spin for operator 1309 & 1, ! one root 1310 & ispin, ! singlet, doublet, triplet, ... 1311 & -1, ! # of davidson iterations 1312 & -1, ! 1/2-particle density calcs 1313 & iprsig, ! print level 1314 & 1.0d-10, ! screening factor 1315 & 0.0d0, ! davidson ci convergence 1316 & 0.0d0, ! davidson ci convergence (auxiliary) 1317 & .false., ! do not calculate nat. orb. occ. num. (internally inside LUCITA) 1318 & .false., ! wave function analysis 1319 & onlyh1, ! no 2-el operators active 1320 & .true., ! integrals provided by / return density matrices to the MCSCF environment 1321 & .false., ! calculate 2-particle density matrix 1322 & .false., ! calculate transition density matrix 1323 & xctype1, ! vector exchange type 1 1324 & xctype2, ! vector exchange type 2 1325 & .false., ! vector exchange active in parallel runs in mc2lu interface (both mc-->lu and lu-->mc) 1326 & .false.) ! io2io vector exchange between mcscf and lucita: default mem2io/io2mem 1327 1328!#define LUCI_DEBUG 1329#ifdef LUCI_DEBUG 1330 write(lupri,*) ' before sigma vec run: c' 1331 call wrtmatmn(c,1,NCVEC,1,NCVEC,lupri) 1332 write(lupri,*) ' before sigma vec run: hc' 1333 call wrtmatmn(hc,1,NHCVEC,1,NHCVEC,lupri) 1334#endif 1335 1336 call mcscf_lucita_interface(c,hc,fcac,h2ac,'sigma vec ', 1337 & wrk,lfree,iprsig) 1338 1339#ifdef LUCI_DEBUG 1340 write(lupri,*) ' after sigma vec run: c' 1341 call wrtmatmn(c,1,NCVEC,1,NCVEC,lupri) 1342 write(lupri,*) ' after sigma vec run: hc' 1343 call wrtmatmn(hc,1,NHCVEC,1,NHCVEC,lupri) 1344#endif 1345 1346 else if(ci_program .eq. 'SIRIUS-CI')then 1347 1348#ifdef LUCI_DEBUG 1349 write(lupri,*) ' before sigma vec run: c' 1350 call wrtmatmn(c,1,NCVEC,1,NCVEC,lupri) 1351 write(lupri,*) ' before sigma vec run: hc' 1352 call wrtmatmn(hc,1,NHCVEC,1,NHCVEC,lupri) 1353#endif 1354 1355 CALL CISGD2(ICSYM,IHCSYM,NCVEC,NHCVEC,C,HC, 1356 * FCAC,H2AC,ONLYH1,IH8SM,XNDXCI,WRK,LFREE) 1357 1358#ifdef LUCI_DEBUG 1359 write(lupri,*) ' after sigma vec run: c' 1360 call wrtmatmn(c,1,NCVEC,1,NCVEC,lupri) 1361 write(lupri,*) ' after sigma vec run: hc' 1362 call wrtmatmn(hc,1,NHCVEC,1,NHCVEC,lupri) 1363#undef LUCI_DEBUG 1364#endif 1365 end if 1366 1367 END IF 1368 1369 RETURN 1370 END 1371C /* Deck cisgd2 */ 1372 SUBROUTINE CISGD2(ICSYM,IHCSYM,NCDET,NHCDET,C,HC, 1373 * FCAC,H2AC,ONLYH1,IH8SM,XNDXCI,WRK,LFREE) 1374C 1375C MASTER ROUTINE FOR DIRECT CI CALCULATION 1376C 1377C SOME INPUT : 1378C C : VECTOR TO BE MULTIPLIED WITH H, 1379C SYMMETRY OF C IS ICSYM, LENGTH IS NCDET 1380C FCAC : ONE-ELECTRON HAMILTONIAN MODIFIED FOR CORE ELECTRONS 1381C *** FULL MATRIX *** 1382C 1383C XNDXCI: ARRAY CONTAINING STRING INFORMATION 1384C AS OBTAINED IN DETINF 1385C OUTPUT : 1386C HC : H TIMES C, OF SYMMETRY IHCSYM AND LENGTH NHCDET 1387C 1388C 1389#include "implicit.h" 1390 DIMENSION C(*), HC(*), FCAC(*) 1391 DIMENSION H2AC(*), WRK(LFREE), XNDXCI(*) 1392 LOGICAL ONLYH1,IH8SM,PERMSM 1393C 1394C THRSML : THRESHOLD FOR INTEGRALS AND CI COEFFICIENTS TO BE ZAPPED 1395C 1396 PARAMETER ( THRSML = 1.0D-14 ) 1397C 1398#include "priunit.h" 1399#include "mxsmob.h" 1400C 1401C Used from common blocks: 1402C INFINP : FLAG(66) 1403C INFORB : MULD2H(8,8), NASHT, N2ASHX,NNASHX 1404C INFPRI : IPRSIG 1405C 1406C STRNUM : EQUAL,NAEL,NASTR,NAEXCI, NBEL,... 1407C DETBAS : core allocation for XNDXCI 1408C MXBLK : MXSASM,MXVBLK 1409C CBESPN : ISPIN1,ISPIN2 1410C SPINFO : MS2,? 1411C CBGETDIS : DISTYP, IADINT 1412C INFSPI : ISGN1,ISGN2 1413C 1414#include "maxorb.h" 1415#include "infinp.h" 1416#include "inforb.h" 1417#include "infpri.h" 1418C 1419#include "mxpdim.h" 1420#include "strnum.h" 1421#include "detbas.h" 1422#include "mxblk.h" 1423#include "cbespn.h" 1424#include "spinfo.h" 1425#include "cbgetdis.h" 1426#include "infspi.h" 1427C 1428C 1429 CALL GETTIM(T0,W0) 1430C 1431C *** Set some control variables *** 1432C 1433C IDOH2 : <> 0 ' NORMAL ' DIRECT CI WITH ONE- AND TWO- BODY TERMS 1434C : = 0 DIRECT CI WITH ONE- BODY TERMS ONLY 1435C 1436C.Spin permutation simplifications 1437C 980826-hjaaj: added ICSYM.eq.IHCSYM test, PERMSM only works 1438C for totally symmetric operators. 1439C 990427-hjaaj: tests show PERMSM does not always work in this version, 1440C we therefore never use PERMSM. 1441Chj IF(MS2.EQ.0.AND..NOT.FLAG(66) .AND. ICSYM.EQ.IHCSYM 1442Chj & .AND.ISPIN1.EQ.0.AND.ISPIN2.EQ.0) THEN 1443Chj PERMSM = .TRUE. 1444Chj IF(MOD((MULTS+1)/2,2).EQ.1) THEN 1445Chj PSIGN = 1.0D0 1446Chj ELSE 1447Chj PSIGN = -1.0D0 1448Chj END IF 1449Chj ELSE 1450 PERMSM = .FALSE. 1451 PSIGN = 999999999.D0 1452Chj END IF 1453C 1454 IF (ONLYH1) THEN 1455 IDOH2 = 0 1456 ELSE 1457 IDOH2 = 1 1458 END IF 1459C 1460C INTFRM = 1 : use GETIN2 1461C = 3 : get integrals directly from H2AC 1462C 1463 IF (DISTYP .EQ. 1 .AND. IADINT .LT. 0) THEN 1464 INTFRM = 3 1465 ELSE 1466 INTFRM = 1 1467 END IF 1468C 1469C Make sure ISGN1/ISGN2 is defined for the DISTYPs where 1470C they are needed in GETIN2. They are now only defined in 1471C RSP2GR /950524-hjaaj 1472C 1473 ISGN1 = (-1)**ISPIN1 1474 ISGN2 = (-1)**ISPIN2 1475C 1476C core energy neglected 1477 ECORE = 0.0D0 1478c. Store determinant CI off-sets in C arrays and HC arrays 1479 CALL CIOFF(ICSYM ,1,XNDXCI,IPRSIG) 1480 CALL CIOFF(IHCSYM,2,XNDXCI,IPRSIG) 1481C CALL CIOFF(IREFSM,ICORHC,XNDXCI,NTEST) 1482C 1483C** 2 : DO H TIMES C 1484c 1485 KFREE = 1 1486 CALL MEMADD(KUFCAC,N2ASHX, KFREE,2) 1487 CALL MEMADD(KRIJKL,N2ASHX*MXSMOB, KFREE,2) 1488 CALL MEMADD(KVEC1,MXSASM,KFREE,2) 1489 CALL MEMADD(KVEC2,MXSASM*MXSMOB,KFREE,2) 1490 CALL MEMADD(KVEC3,MXSASM*MXSMOB,KFREE,2) 1491 CALL MEMADD(KINDX1,MXSASM,KFREE,1) 1492 CALL MEMADD(KINDX2,MXSASM,KFREE,1) 1493 CALL MEMADD(KL,MXSASM*MXSMOB,KFREE,1) 1494 CALL MEMADD(KR,MXSASM*MXSMOB,KFREE,1) 1495 CALL MEMADD(KC2,MXVBLK,KFREE,2) 1496 CALL MEMADD(KWIJKL,N2ASHX,KFREE,2) 1497 CALL MEMADD(KINDE2,MXPST*MXPTP,KFREE,1) 1498 CALL MEMADD(KF3,MXPST*MXPTP,KFREE,2) 1499 IF (KFREE .GT. LFREE) CALL ERRWRK('CISIGD',KFREE,LFREE) 1500c 1501c Reorder FCAC from Sirius order to Lunar order. 1502c 1503 CALL REORMT(FCAC,WRK(KUFCAC),NASHT,NASHT,XNDXCI(KLTSOB), 1504 & XNDXCI(KLTSOB)) 1505C 1506 CALL GETTIM(T1,W1) 1507C 1508 CALL CISIG9(XNDXCI(KTAIJ),XNDXCI(KTATO),NAEL,NASTR,NAEXCI, 1509 & XNDXCI(KTBIJ),XNDXCI(KTBTO),NBEL,NBSTR,NBEXCI, 1510 & NASHT,C,HC,NCDET,NHCDET,XNDXCI(KSTASA),XNDXCI(KSTASB), 1511 & XNDXCI(KCOFF),XNDXCI(KHCOFF),MAXSYM,XNDXCI(KISSYM), 1512 & XNDXCI(KSTBAA),XNDXCI(KSTBAB),XNDXCI(KTASYM),XNDXCI(KTBSYM), 1513 & WRK(KUFCAC),WRK(KRIJKL),ICSYM,IHCSYM,PERMSM, 1514 & XNDXCI(KIOCOC),NOCTPA,NOCTPB, 1515 & XNDXCI(KICSO), XNDXCI(KIHCSO),XNDXCI(KICOOS),XNDXCI(KIHOOS), 1516 & XNDXCI(KNSSOA),XNDXCI(KISSOA),XNDXCI(KNSSOB),XNDXCI(KISSOB), 1517 & XNDXCI(KKLTP), 1518 & XNDXCI(KICREA),XNDXCI(KIANNI),WRK(KVEC1),WRK(KVEC2), 1519 & WRK(KC2),WRK(KINDX1),WRK(KINDX2),WRK(KL), 1520 & WRK(KR),PSIGN,WRK(KVEC3),XNDXCI(KTPFSA),XNDXCI(KTPFSB), 1521 & XNDXCI(KCDTAS),XNDXCI(KHDTAS),XNDXCI(KKLCAN),XNDXCI(KTPFOB), 1522 & HC,IDOH2,ECORE,H2AC,WRK(KWIJKL),IPRSIG,XNDXCI(KIASTR), 1523 & XNDXCI(KIBSTR),MXSASM,XNDXCI(KLTSOB),XNDXCI(KSTLOB), 1524 & ISPIN1,ISPIN2,IH8SM,INTFRM,WRK(KINDE2),WRK(KF3)) 1525 CALL GETTIM(T2,W2) 1526 IF (IPRSIG .GE. 3) WRITE (LUPRI,'(/A,2I10)') 1527 & ' CISIGD : KFREE, LFREE = ',KFREE,LFREE 1528 IF (IPRSIG .GE. 1) WRITE (LUPRI,'(A,2F10.3)') 1529 & ' CISIGD : CPU and wall time (sec) :',T2-T1,W2-W1 1530C 1531 RETURN 1532 END 1533C /* Deck cisigo */ 1534 SUBROUTINE CISIGO(NOSIM,SOVECS,CREF,EMYX,FXCAC,H2XAC, 1535 * XNDXCI,WRK,LFREE) 1536C 1537C Written 18-Jan-1988 J.O. 1538C (After cisigo for casguga of hjaaj ) 1539C 1540C Parameter list: 1541C NOSIM number of simultaneous orbital expansion vectors 1542C CREF reference CI-vector 1543C SOVECS NOSIM sigma vectors of orbital trial vectors 1544C EMYX the inactive "energy" from the one-index transformed 1545C "Hamiltonian" 1546C FXCAC 1-ind. transf. inactive Fock matrix 1547C H2XAC 1-ind. transf. active 2-el. integrals 1548C XNDXCI CI information 1549C 1550C Output: 1551C SOVECS(I) = SUM(sr) of K(I,sr)*BOVECS(sr) 1552C 1553C Scratch: 1554C WRK work area containing : 1555C 1556C 1557! module dependencies 1558 use lucita_mcscf_ci_cfg 1559#include "implicit.h" 1560C 1561 DIMENSION CREF(NCONRF), SOVECS(NVARPT,NOSIM),EMYX(NOSIM),XNDXCI(*) 1562 DIMENSION FXCAC(NNASHX,*),H2XAC(NNASHX,NNASHX,*) 1563 DIMENSION WRK(LFREE) 1564 PARAMETER (IDTIM = 10) 1565 PARAMETER ( D2 = 2.0D0 ) 1566 LOGICAL NOH2, IH8SM 1567 PARAMETER ( NOH2 = .FALSE. , IH8SM = .TRUE. ) 1568C NOH2 ... 2-electron part to be included 1569C IH8SM ... integrals have 8-fold symmetry 1570C 1571C 1572C Used from common blocks: 1573C 1574C INFORB : NASHT,NNASHX,N2ASHX 1575C INFLIN : LSYMRF,LSYMPT,LSYMST,NCONRF,NCONST 1576C INFTIM : NCALLS,TIMCPU,TIMWAL ! IDTIM is index for these 1577C 1578#include "priunit.h" 1579#include "maxorb.h" 1580#include "inforb.h" 1581#include "inflin.h" 1582#include "inftim.h" 1583#include "infinp.h" 1584C 1585C 1586 KUFCAC = 1 1587 KW1 = KUFCAC + N2ASHX 1588 LW1 = LFREE - KW1 1589C 1590 ISPIN1 = 0 1591 ISPIN2 = 0 1592C ... singlet-singlet coupling of 2-electron integrals 1593C 1594C 1595C *** Calculate sigma vectors with modified integrals 1596C 1597 CALL GETTIM(T0,W0) 1598 CALL GETTIM(T1,W1) 1599 CALL GETTIM(T2,W2) 1600 CALL GETTIM(T3,W3) 1601C 1602! 1603! set flag for orbital trial vector in MCSCF (needed for LUCITA interface) 1604 mcscf_orbital_trial_vector = .true. 1605C 1606 DO 300 IOSIM = 1, NOSIM 1607C 1608 if(ci_program .eq. 'SIRIUS-CI')then 1609 CALL DSPTSI(NASHT,FXCAC(1,IOSIM),WRK(KUFCAC)) 1610C ... unpack FXCAC using CALL DSPTSI(N,ASP,ASI) 1611 else if(ci_program .eq. 'LUCITA ')then 1612 CALL DCOPY(NNASHX,FXCAC(1,IOSIM),1,WRK(KUFCAC),1) 1613! set vector exchange type: cref 1614 vector_exchange_type1 = 1 1615 end if 1616 CALL DZERO(SOVECS(1,IOSIM),NCONST) 1617 CALL CISIGD(LSYMRF,LSYMST,NCONRF,NCONST, 1618 & CREF,SOVECS(1,IOSIM),WRK(KUFCAC),H2XAC(1,1,IOSIM), 1619 & NOH2,IH8SM,XNDXCI,ISPIN1,ISPIN2,WRK(KW1),LW1) 1620C CALL CISIGD(ICSYM,IHCSYM,NCDET,NHCDET, C,HC, FCAC,H2AC, 1621C & NOH2,IH8SM,XNDXCI,ISPIN1,ISPIN2,WRK,LFREE) 1622C 1623C 1624C --- Add inactive energy contributions 1625C --- multiply orbital sigma vectors by 2 to get final SOVECS 1626C 1627 IF (LSYMPT .EQ. 1) THEN 1628 DO 200 I = 1,NCONST 1629 SOVECS(I,IOSIM) = D2*(SOVECS(I,IOSIM) + CREF(I)*EMYX(IOSIM)) 1630 200 CONTINUE 1631 ELSE 1632 CALL DSCAL(NCONST,D2,SOVECS(1,IOSIM),1) 1633 END IF 1634C 1635 300 CONTINUE 1636C 1637! 1638! re-set flag for orbital trial vector in MCSCF (needed for LUCITA interface) 1639 mcscf_orbital_trial_vector = .false. 1640C 1641 CALL GETTIM(T4,W4) 1642C 1643 NCALLS( IDTIM) = NCALLS( IDTIM) + 1 1644C!!! NVECS ( IDTIM) = NVECS ( IDTIM) + NOSIM 1645 TIMCPU(1,IDTIM) = TIMCPU(1,IDTIM) + T4 - T0 1646 TIMCPU(2,IDTIM) = TIMCPU(2,IDTIM) + T1 - T0 1647 TIMCPU(3,IDTIM) = TIMCPU(3,IDTIM) + T2 - T1 1648 TIMCPU(4,IDTIM) = TIMCPU(4,IDTIM) + T3 - T2 1649 TIMCPU(5,IDTIM) = TIMCPU(5,IDTIM) + T4 - T3 1650 TIMWAL(1,IDTIM) = TIMWAL(1,IDTIM) + W4 - W0 1651 TIMWAL(2,IDTIM) = TIMWAL(2,IDTIM) + W1 - W0 1652 TIMWAL(3,IDTIM) = TIMWAL(3,IDTIM) + W2 - W1 1653 TIMWAL(4,IDTIM) = TIMWAL(4,IDTIM) + W3 - W2 1654 TIMWAL(5,IDTIM) = TIMWAL(5,IDTIM) + W4 - W3 1655C 1656C 1657C *** End of subroutine CISIGO 1658C 1659 RETURN 1660 END 1661