Subroutine DetCalcInit Use global_utilities Use IntegralsData, only: NFROZEN use SystemData, only: lms, lms2, nBasis, nBasisMax, nEl, SymRestrict, & Alat, arr, brr, boa, box, coa, ecore, g1, Beta, & tParity, tSpn, Symmetry, STot, NullBasisFn, tUHF, & tMolpro, tZeroDelimitedFCIDUMP use sym_mod use LoggingData, only: tLogDets use HElem use util_mod, only: get_free_unit, NECI_ICOPY Type(BasisFn) ISym integer i, j, ii, iunit integer ierr, norb integer nDetTot character(25), parameter :: this_routine = 'DetCalcInit' IF (.NOT. TCALCHMAT) THEN write(stdout, *) "Not storing the H matrix." IF (TENERGY .AND. .NOT. TBLOCK) THEN write(stdout, *) "Cannot calculate energies without blocking the Hamiltonian." TENERGY = .FALSE. end if IF (TENERGY .AND. NBLK /= 0) THEN !C.. We're doing a Lanczos without calculating the H mat write(stdout, *) "Cannot perform Lanczos without Hamiltonian" TENERGY = .FALSE. end if end if ! If we want to have UMat2D, and it isn't yet filled in, generate it ! here. All of the integrals setup/freezing etc is done... if (tDeferred_Umat2d .and. .not. tUMat2D) then ASSERT(.not. btest(nbasis, 0)) ! And fill in the array call SetupUMat2d_dense(nBasis) end if !Copied Specdet information from Calc.F, so if inspect is present, but no determinant/csf specified, it will still run. if (TSPECDET) then if (.not. associated(specdet)) then !specdet not allocated. Allocate it and copy fdet allocate(specdet(nel)) write(stdout, *) "TSPECDET set, but not allocated. using FDET" CALL NECI_ICOPY(NEL, FDET, 1, SPECDET, 1) else if (.not. ISVALIDDET(SPECDET, NEL)) then write(stdout, *) "TSPECDET set, but invalid. using FDET" CALL NECI_ICOPY(NEL, FDET, 1, SPECDET, 1) end if end if !C IF(TCALCHMAT.OR.NPATHS.NE.0.OR.DETINV.GT.0.OR.TBLOCK) THEN iExcitLevel = ICILEVEL IF (tFindDets) THEN !C..Need to determine the determinants IF (iExcitLevel /= 0) THEN write(stdout, *) "Performing truncated CI at level ", iExcitLevel IF (TSPECDET) THEN write(stdout, *) "Using SPECDET:" call write_det(stdout, SPECDET, .true.)! CALL NECI_ICOPY(NEL, SPECDET, 1, FDET, 1) ELSE write(stdout, *) "Using Fermi DET:" call write_det(stdout, FDET, .true.) end if !C.. if we're doing a truncated CI expansion CALL GENEXCIT(FDET, iExcitLevel, NBASIS, NEL, reshape([0], [1, 1]), & reshape([0], [1, 1]), NDET, 1, G1, .TRUE., NBASISMAX, .TRUE.) write(stdout, *) "NDET out of GENEXCIT ", NDET !C.. We need to add in the FDET NDET = NDET + 1 II = NDET NBLOCKS = 1 else if (TBLOCK) THEN write(stdout, *) "Determining determinants and blocks." IF (TPARITY) THEN write(stdout, *) "Using symmetry restriction:" CALL writeallsym(stdout, SymRestrict, .TRUE.) end if IF (TSPN) THEN write(stdout, *) "Using spin restriction:", LMS end if if (tZeroDelimitedFCIDUMP) then !When breaking spin symmetry in molpro, it is important to occupy alpha orbs preferentially CALL GNDTS_BLK(NEL, nBasis, BRR, NBASISMAX, NMRKS, .TRUE., & NDET, G1, II, NBLOCKSTARTS, NBLOCKS, TSPN, -LMS2, TPARITY, & SymRestrict, IFDET,.NOT. TREAD, NDETTOT, BLOCKSYM) else CALL GNDTS_BLK(NEL, nBasis, BRR, NBASISMAX, NMRKS, .TRUE., & NDET, G1, II, NBLOCKSTARTS, NBLOCKS, TSPN, LMS2, TPARITY, & SymRestrict, IFDET,.NOT. TREAD, NDETTOT, BLOCKSYM) end if write(stdout, *) "NBLOCKS:", NBLOCKS write(stdout, *) "Determining determinants." IF (TPARITY) THEN write(stdout, *) "Using symmetry restriction:" CALL writeallsym(stdout, SymRestrict, .TRUE.) end if IF (TSPN) THEN write(stdout, *) "Using spin restriction:", LMS end if if (tZeroDelimitedFCIDUMP) then !When breaking spin symmetry in molpro, it is important to occupy alpha orbs preferentially CALL GNDTS(NEL, nBasis, BRR, NBASISMAX, NMRKS, .TRUE., G1, TSPN, -LMS, TPARITY, SymRestrict, II, IFDET) else CALL GNDTS(NEL, nBasis, BRR, NBASISMAX, NMRKS, .TRUE., G1, TSPN, LMS, TPARITY, SymRestrict, II, IFDET) end if NBLOCKS = 1 NDET = II end if !C.. IF (II == 0) THEN write(stdout, *) "No determinants found. Cannot continue" call stop_all(this_routine, "No determinants found. Cannot continue") end if !C.. NEL now only includes active electrons write(stdout, *) "Number of determinants found to be: ", II write(stdout, *) "Allocating initial memory for calculation of energy..." CALL neci_flush(stdout) allocate(NMrks(nEl, II), stat=ierr) LogAlloc(ierr, 'NMRKS', NEL * II, 4, tagNMRKS) NMRKS(1:NEL, 1:II) = 0 allocate(NBLOCKSTARTS(NBLOCKS + 1), stat=ierr) call LogMemAlloc('NBLOCKSTARTS', NBLOCKS + 1, 4, this_routine, tagNBLOCKSTARTS, ierr) NBLOCKSTARTS(1:NBLOCKS + 1) = 0 allocate(BlockSym(NBLOCKS + 1), stat=ierr) LogAlloc(ierr, 'BLOCKSYM', NBLOCKS + 1, BasisFNSizeB, tagBlockSym) BLOCKSYM(1:NBLOCKS) = NullBasisFn !C.. NDET = II IF (iExcitLevel /= 0) THEN !C.. Use HAMIL to temporarily hold a list of excitation levels CALL NECI_ICOPY(NEL, FDET, 1, NMRKS, 1) allocate(Hamil(II), stat=ierr) LogAlloc(ierr, 'HAMIL', II, HElement_t_sizeB, tagHamil) NDET = 0 CALL GENEXCIT(& FDET, iExcitLevel, NBASIS, NEL, & NMRKS(:1, :2), int(HAMIL), NDET, 1, G1, .TRUE., NBASISMAX, .FALSE.) Deallocate(Hamil) LogDealloc(tagHamil) NDET = NDET + 1 NBLOCKSTARTS(1) = 1 NBLOCKSTARTS(2) = II + 1 IFDET = 1 else if (TBLOCK) THEN if (tZeroDelimitedFCIDUMP) then !When breaking spin symmetry in molpro, it is important to occupy alpha orbs preferentially CALL GNDTS_BLK(NEL, nBasis, BRR, NBASISMAX, NMRKS, .FALSE., NDET, G1, II, NBLOCKSTARTS, NBLOCKS, TSPN, -LMS2, TPARITY, & SymRestrict, IFDET,.NOT. TREAD, NDETTOT, BLOCKSYM) else CALL GNDTS_BLK(NEL, nBasis, BRR, NBASISMAX, NMRKS, .FALSE., NDET, G1, II, NBLOCKSTARTS, NBLOCKS, TSPN, LMS2, TPARITY, & SymRestrict, IFDET,.NOT. TREAD, NDETTOT, BLOCKSYM) end if ELSE if (tZeroDelimitedFCIDUMP) then !When breaking spin symmetry in molpro, it is important to occupy alpha orbs preferentially CALL GNDTS(NEL, nBasis, BRR, NBASISMAX, NMRKS, .FALSE., G1, TSPN, -LMS, TPARITY, SymRestrict, II, IFDET) else CALL GNDTS(NEL, nBasis, BRR, NBASISMAX, NMRKS, .FALSE., G1, TSPN, LMS, TPARITY, SymRestrict, II, IFDET) end if NBLOCKSTARTS(1) = 1 NBLOCKSTARTS(2) = II + 1 end if if (tLogDets) THEN iunit = get_free_unit() open(iunit, FILE='DETS', STATUS='UNKNOWN') DO I = 1, NDET call write_det(iunit, NMRKS(:, I), .false.) CALL GETSYM(NMRKS(:, I), NEL, G1, NBASISMAX, ISYM) CALL WRITESYM(iunit, ISym%Sym, .TRUE.) end do close(iunit) end if !Update: 14.03.2018, K.Ghanem !FDET has already been assigned in DetPreFreezeInit. !No idea why it is overwirtten here. !Therefore, I comment out the following code and hope for the best. !Instead, I look for FDET in the list of determinants NMRKS and assign the index to IFDET. !C.. Now generate the fermi determiant !C.. Work out the fermi det IFDET = 0 DO I = 1, NDET IF (ALL(NMRKS(:, I) == FDET)) THEN IFDET = I Exit END IF END DO IF (IFDET == 0) call stop_all("DetCalcInit", "Fermi determinant is not found in NMRKS!") write(stdout, *) ' NUMBER OF SYMMETRY UNIQUE DETS ', NDET !C write(stdout,*) ' TOTAL NUMBER OF DETS.' , NDETTOT IF (NEVAL == 0) THEN write(stdout, *) 'NEVAL=0. Setting NEVAL=NDET' NEVAL = NDET end if IF (NEVAL > NDET) THEN write(stdout, *) 'NEVAL>NDET. Setting NEVAL=NDET' NEVAL = NDET end if IF (ABS(DETINV) > NDET) THEN write(stdout, *) 'DETINV=', DETINV, '>NDET=', NDET write(stdout, *) 'Setting DETINV to 0' DETINV = 0 end if CALL neci_flush(stdout) !C ==----------------------------------------------------------------== !C..Set up memory for c's, nrow and the label IF (TCALCHMAT) THEN write(stdout, *) "CK Size", NDET * NEVAL * HElement_t_size allocate(CkN(nDet, nEval), stat=ierr) LogAlloc(ierr, 'CKN', nDet * nEval, HElement_t_sizeB, tagCKN) CKN = (0.0_dp) allocate(Ck(nDet, nEval), stat=ierr) LogAlloc(ierr, 'CK', nDet * nEval, HElement_t_sizeB, tagCK) CK = (0.0_dp) allocate(W(nEval), stat=ierr) LogAlloc(ierr, 'W', nEval, 8, tagW) W = 0.0_dp end if IF (TREAD) THEN #ifdef CMPLX_ call stop_all(this_routine, "not implemented for complex with tREAD") #else CALL READ_PSI(BOX, BOA, COA, NDET, NEVAL, NBASISMAX, NEL, CK, W) #endif end if end if End Subroutine DetCalcInit