diff --git a/flosic/frmorb_mpi.ftn b/flosic/frmorb_mpi.ftn index 7f820ad..42239e6 100644 --- a/flosic/frmorb_mpi.ftn +++ b/flosic/frmorb_mpi.ftn @@ -1,10 +1,11 @@ C UTEP Electronic Structure Lab (2020) SUBROUTINE SICLAGM(LSPN) use debug1 + use global_inputs,only : CALCTYPE1 use mpidat1,only : IRANK,NPROC,SHM_SIZE use mesh1,only : wmsh,rmsh,nmsh use common2,only : NSPN,IGGA !,RIDT,N_CON,LSYMMAX,N_POS,NFNCT - use common5,only : NWFS !,NWF,PSI + use common5,only : NWFS,CONVERGENCE!,NWF,PSI use mixpot1,only : POT=>POTOUT !POTIN ! use pot_dens,only : COULOMB,RHOG use for_diag1 @@ -262,15 +263,30 @@ C PRINT*,'SMCHG:',CHGUP,CHGDN,SMCHG,TMKIN,NWFS(1) 1000 FORMAT(' ',15F12.5) ! DEALLOCATE(SICP) ! DEALLOCATE(POTFILE) - DO LFM=1,NFRM(LSPN) - CLOSE(206+LFM+IRBOFS,STATUS='DELETE') - print *,206+LFM+IRBOFS," CLOSED AND DELETED" - END DO - IF(LIBXC1.OR.ISMGGA)THEN + IF(CONVERGENCE.AND. CALCTYPE1.NE.2) THEN + !Keep ZPOT FILES for post conv. proc. DO LFM=1,NFRM(LSPN) - CLOSE(106+LFM+IRBOFS,STATUS='DELETE') - print *,106+LFM+IRBOFS," CLOSED AND DELETED" + CLOSE(206+LFM+IRBOFS) + print *,206+LFM+IRBOFS," CLOSED" END DO + IF(LIBXC1.OR.ISMGGA)THEN + DO LFM=1,NFRM(LSPN) + CLOSE(106+LFM+IRBOFS) + print *,106+LFM+IRBOFS," CLOSED" + END DO + ENDIF + ELSE !.not. convergence + !Clean up ZPOT files for the next iter. step. + DO LFM=1,NFRM(LSPN) + CLOSE(206+LFM+IRBOFS,STATUS='DELETE') + print *,206+LFM+IRBOFS," CLOSED AND DELETED" + END DO + IF(LIBXC1.OR.ISMGGA)THEN + DO LFM=1,NFRM(LSPN) + CLOSE(106+LFM+IRBOFS,STATUS='DELETE') + print *,106+LFM+IRBOFS," CLOSED AND DELETED" + END DO + ENDIF ENDIF DEALLOCATE(ZMGGAS) DEALLOCATE(POTSMALL) diff --git a/flosic/main.ftn b/flosic/main.ftn index b46b278..49565d3 100644 --- a/flosic/main.ftn +++ b/flosic/main.ftn @@ -627,10 +627,8 @@ C IF(JSPNX.EQ.2)CALL FIXTMAT(1) CALL GTTIME(TIMET3) IF(CALCTYPE1.NE.2) THEN %ifdef ATOMFORCE - DO IORBX=1,NFRM(JSPNX) - WRITE(6,*) 'CALLING SICLAG_DER FOR FLO', IORBX - CALL SICLAG_DER(JSPNX,IORBX) - END DO + WRITE(6,*) 'CALLING SICLAG_DER FOR SPIN',JSPNX + CALL SICLAG_DER(JSPNX) %else WRITE(6,*) 'SKIPPING SICLAG_DER FOR FLO' WRITE(6,*) 'RECOMPILE WITH ATOMFORCE TO USE THE FEATURE' diff --git a/flosic/siclag_der.ftn b/flosic/siclag_der.ftn index 5fccd66..dd10247 100644 --- a/flosic/siclag_der.ftn +++ b/flosic/siclag_der.ftn @@ -1,5 +1,5 @@ C UTEP Electronic Structure Lab (2020) - SUBROUTINE SICLAG_DER(LSPN,LFM) + SUBROUTINE SICLAG_DER(LSPN) ! ! KAJ 1-8-2018 This computes ! for each atomic site, nu; compute for the orbital LFM @@ -8,7 +8,7 @@ C UTEP Electronic Structure Lab (2020) use common2,only : RIDT,N_CON,LSYMMAX,N_POS,NFNCT,NSPN use common5,only : NWF,NWFS,PSI use mixpot1,only : POTIN,POT=>POTOUT - use pot_dens,only : COULOMB,RHOG +! use pot_dens,only : COULOMB,RHOG !use V_SIC from ZPOT instead use for_diag1 !SIC modules use SICMAT,only : SIC,DERSIC @@ -23,14 +23,20 @@ C UTEP Electronic Structure Lab (2020) & IERR, IFNCT, IGR, ILOC, ISHDUM, ISHELLA, IWF, J_POS, JFM, JFN, & JPTS, KPTS, L_NUC, LFN, LI, LMAX1, LPV, M_NUC, MU, NDERV, NPV REAL*8 :: SYMBOL , APT1, APT2, CHGDN, CHGE, CHGUP, FMAT, - & RHI, RPTS, SICP, SMCHG, TMKIN, WMSA, ADD1, ADD2, - & ADD3,ADD4, AGRAD, FACTOR,over_test,testsic + & RHI, RPTS, SMCHG, TMKIN, WMSA, ADD1, ADD2, + & ADD3,ADD4, AGRAD, FACTOR,over_test,testsic, + & TIMET1,TIMET2 + !Additional vars + INTEGER :: IRBOFS + REAL(8),ALLOCATABLE :: POTSMALL(:,:) + CHARACTER*12 :: ZPOTSTR PARAMETER (NMAX=MPBLOCK) LOGICAL EXIST,FIRST,IUPDAT DIMENSION WMSA(NMAX),RPTS(3,NMAX) - DIMENSION SICP(NMAX) + !DIMENSION SICP(NMAX) !Changed the size + REAL(8) :: SICP(NMAX,MAX_OCC) !New size C DIMENSION ISIZE(3) DATA ISIZE/1,3,6/ @@ -38,21 +44,35 @@ C INTEGER MFRM SAVE + + CALL GTTIME(TIMET1) c c Print tests as needed... c c print *, 'in SICLAG_DER: LSPN LFM', LSPN, LFM c +! allocate array to store V_SIC temporarily + ALLOCATE(POTSMALL(NMAX,MAX_OCC)) + NGRAD=4 c c initialize test integrals c over_test = 0.0d0 testsic = 0.0d0 + IF(LSPN.EQ.1) dersic = 0.0d0 C C LOOP OVER ALL POINTS C CALL GTTIME(APT1) + IRBOFS=(LSPN-1)*NWFS(1) + + DO LFM=1,NFRM(LSPN) + WRITE(ZPOTSTR,'(A,I4.4)')'ZPOT',LFM+IRBOFS + print *,206+LFM+IRBOFS,"OPENING ",ZPOTSTR + OPEN(206+LFM+IRBOFS,FILE=ZPOTSTR,FORM='UNFORMATTED', + & STATUS='OLD') + END DO LPTS=0 10 CONTINUE IF(LPTS+NMAX.LT.NMSH)THEN @@ -60,9 +80,19 @@ C ELSE MPTS=NMSH-LPTS END IF +! Read potential from files + DO LFM=1,NFRM(LSPN) + DO IPTS=1,MPTS + READ(206+LFM+IRBOFS) POTSMALL(IPTS,LFM) + END DO + END DO + DO IPTS=1,MPTS WMSA(IPTS)=WMSH(LPTS+IPTS) - SICP(IPTS)=POT(LPTS+IPTS) +! SICP(IPTS)=POT(LPTS+IPTS) !commented + DO LFM=1,NFRM(LSPN) + SICP(IPTS,LFM)=POTSMALL(IPTS,LFM) + END DO DO J=1,3 RPTS(J,IPTS)=RMSH(J,LPTS+IPTS) END DO @@ -73,7 +103,9 @@ c written. It would follow the same format as in frcslv c c %ifndef MPI c - CALL SICLM_DERSLV(LSPN,LFM,MPTS,WMSA,SICP,RPTS, +! CALL SICLM_DERSLV(LSPN,LFM,MPTS,WMSA,SICP,RPTS, +! & over_test,testsic) + CALL SICLM_DERSLV(LSPN,MPTS,WMSA,SICP,RPTS, & over_test,testsic) C %else c CALL PAMLM_DERSIC(1,LSPN,LFM,MPTS,WMSA,SICP,RPTS) @@ -89,23 +121,37 @@ C PRINT*,'TIME FOR DERSIC:',LFM,APT2-APT1 C CALL PAMLMSIC(2,LSPN,LFM,MPTS,WMSA,SICP,RPTS) C %endif c - LFN = LFM + LFRM(1)*(LSPN-1) + !LFN = LFM + LFRM(1)*(LSPN-1) !commented c c print test integrals as needed... c c PRINT*, 'LFM OVER_TEST', LFN, over_test c PRINT*, 'LFM TESTSIC', LFN, testsic c PRINT*,'DERSIC MATRIX FOR ORBITAL LFM:',LFN -c do iid = 1,10 -c PRINT 1000,(DERSIC(IX,LFM,IID),IX=1,3) -c end do + DO LFM=1,LFRM(LSPN) + LFN=LFM+LFRM(1)*(LSPN-1) + PRINT*,'DERSIC MATRIX FOR ORBITAL LFM:',LFN + do iid = 1,10 + PRINT 1000,(DERSIC(IX,LFN,IID),IX=1,3) + end do + END DO 1000 FORMAT(' ',3F12.5) + + DO LFM=1,LFRM(LSPN) + CLOSE(206+LFM+IRBOFS,STATUS='DELETE') + print *,206+LFM+IRBOFS," CLOSED AND DELETED" + END DO + DEALLOCATE(POTSMALL) + + CALL GTTIME(TIMET2) + WRITE(6,*)'TIME FOR SICLAG_DER/SPIN:',TIMET2-TIMET1,LSPN + RETURN END !#################################################################### - SUBROUTINE SICLM_DERSLV(LSPN,LFM,MPTS,WMSA,SICP,RPTS, + SUBROUTINE SICLM_DERSLV(LSPN,MPTS,WMSA,SICP,RPTS, & over_test,testsic) c use mesh1,only : wmsh,rmsh,nmsh use common2,only : RIDT,N_CON,LSYMMAX,N_POS,NFNCT,NSPN @@ -117,7 +163,8 @@ c use mixpot1,only : POTIN,POT=>POTOUT use for_diag1 !SIC modules use SICMAT,only : SIC,DERSIC - use FRM,only : BFRM,RESULTS,LFRM,DEBDAX + use FRM,only : LFRM +! use FRM,only : BFRM,RESULTS,LFRM,DEBDAX use mpidat1,only : IRANK ! Converion to implicit none. Raja Zope Sun Aug 20 09:01:50 MDT 2017 @@ -135,10 +182,10 @@ c use mixpot1,only : POTIN,POT=>POTOUT PARAMETER (NMAX=MPBLOCK) LOGICAL EXIST,FIRST,IUPDAT DIMENSION WMSA(NMAX),RPTS(3,NMAX) - DIMENSION SICP(NMAX) + DIMENSION SICP(NMAX,MAX_OCC) DIMENSION FMAT(MAX_OCC,MAX_OCC,4,2),RHI(4) - REAL*8,allocatable :: PSIG(:,:), PTS(:,:), PSID(:,:,:) + REAL*8,allocatable :: PSIG(:,:,:), PTS(:,:), PSID(:,:,:) & ,GRAD(:,:,:,:,:),GRAD1(:,:,:,:),RVECA(:,:),GSUB(:,:) LOGICAL ,allocatable :: ICOUNT(:,:) C @@ -156,7 +203,7 @@ C c print *, 'in siclm_der_slv', LFM, LSPN c call flush(6) !YY allocate tmp2 arrays - ALLOCATE(PSIG(NMAX,10),STAT=IERR) + ALLOCATE(PSIG(NMAX,10,MAX_OCC),STAT=IERR) IF(IERR/=0)WRITE(6,*)'SICDERSLV:ERROR ALLOCATING PSIG' ALLOCATE(PSID(3,MAX_OCC,NMAX),STAT=IERR) IF(IERR/=0)WRITE(6,*)'SICDERSLV:ERROR ALLOCATING PSID' @@ -173,19 +220,19 @@ c call flush(6) ALLOCATE(ICOUNT(MAX_CON,3),STAT=IERR) IF(IERR/=0)WRITE(6,*)'SICDERSLV:ERROR ALLOCATING ICOUNT' - IF(LFM.EQ.0)THEN - PRINT*,IRANK, LFM,' FIRST LFM IS ZERO' - CALL STOPIT - END IF +! IF(LFM.EQ.0)THEN +! PRINT*,IRANK, LFM,' FIRST LFM IS ZERO' +! CALL STOPIT +! END IF ISPFAC = 2/NSPN NGRAD=4 - FORALL (IPT=1:MPTS,IGR=1:NGRAD) - PSIG(IPT,IGR) = 0.0D0 + FORALL (IWF=1:MAX_OCC,IPT=1:MPTS,IGR=1:NGRAD) + PSIG(IPT,IGR,IWF) = 0.0D0 END FORALL c c Get correct index for local orbital wave function c - LFN = LFM + (LSPN-1)*LFRM(1) + !LFN = LFM + (LSPN-1)*LFRM(1) IID=0 DO 86 IFNCT=1,NFNCT LMAX1=LSYMMAX(IFNCT)+1 @@ -243,16 +290,18 @@ C ILOC=ILOC+1 IF (ICOUNT(ICON,LI)) THEN c DO IWF=1,NWF !only need one orbital + DO LFM=1,LFRM(LSPN) !now get all orbital at once + LFN = LFM + (LSPN-1)*LFRM(1) FACTOR=PSI(ILOC,LFN,1) if(abs(FACTOR) .GT. 1.0d-10) then DO IGR=1,NGRAD DO LPV=1,NPV - PSIG(IPT+LPV,IGR)=PSIG(IPT+LPV,IGR) + PSIG(IPT+LPV,IGR,LFM)=PSIG(IPT+LPV,IGR,LFM) & +FACTOR*GRAD(LPV,IGR,MU,ICON,LI) END DO !LPV: fill in NSPEED points at a time END DO !IGR end if -c END DO !IWF + END DO !IWF END IF END DO !ICON END DO !MU @@ -265,11 +314,17 @@ c END DO !IWF C C CONSTRUCT TEST INTEGRALS: , and c + DO LFM=1,LFRM(LSPN) + !LFN = LFM + (LSPN-1)*LFRM(1) DO IPT=1,MPTS TESTSIC=TESTSIC + - & PSIG(IPT,1)*PSIG(IPT,1)*WMSA(IPT)*SICP(IPT) - OVER_test = over_test + psig(ipt,1)*psig(ipt,1)*wmsa(ipt) + & PSIG(IPT,1,LFM)*PSIG(IPT,1,LFM)*WMSA(IPT)*SICP(IPT,LFM) + OVER_test=over_test + psig(ipt,1,LFM) + & *psig(ipt,1,LFM)*wmsa(ipt) END DO + !PRINT *,"TEST INTEGRAL AND INDEX",over_test,LFM + !PRINT *,"TEST INTEGRAL ",TESTSIC,LFM + END DO 1000 FORMAT(' ',5F12.5) c C The following section computes dphi_i/dR = Sum_s {c^i_s * df_s/dR} @@ -376,16 +431,20 @@ C C C UPDATE MATRIX ELEMENTS C +!Orbital loop added + DO LFM=1,LFRM(LSPN) + LFN = LFM + (LSPN-1)*LFRM(1) DO 3100 IPT=1,MPTS C IWF=NWFS(1)*(ISPN-1) c DO 390 JWF=1,NWFS(ISPN) c IWF=IWF+1 DO 391 IX=1,3 DERSIC(IX,LFN,IID)=DERSIC(IX,LFN,IID) + ISPFAC* - & PSID(IX,LFN,IPT)*SICP(IPT)*PSIG(IPT,1)*WMSA(IPT) + & PSID(IX,LFN,IPT)*SICP(IPT,LFM)*PSIG(IPT,1,LFM)*WMSA(IPT) 391 CONTINUE 390 CONTINUE 3100 CONTINUE + END DO 30 CONTINUE 483 CONTINUE 485 CONTINUE @@ -403,11 +462,15 @@ c print 1000, (psid(1,LFN,ipt),ipt=1,5) C C CONSTRUCT TEST INTEGRALS: , and c +! Orbital loop added +! DO LFM=1,LFRM(LSPN) +! LFN = LFM + (LSPN-1)*LFRM(1) c DO IPT=1,MPTS c TESTSIC=TESTSIC + -c & PSIG(IPT,1)*PSIG(IPT,1)*WMSA(IPT)*SICP(IPT) +c & PSIG(IPT,1,LFM)*PSIG(IPT,1,LFM)*WMSA(IPT)*SICP(IPT,LFM) c OVER_test = over_test + psig(ipt,1)*psig(ipt,1)*wmsa(ipt) c END DO +! END DO c1000 FORMAT(' ',15F12.5) c call flush(6) c