From 520203eb1c05793ce9bda85c56cdaa00e1b1f0e5 Mon Sep 17 00:00:00 2001 From: Akhil Shajan Date: Tue, 7 May 2024 12:21:53 -0400 Subject: [PATCH] Print QCint in binary files --- src/modules/quick_calculated_module.f90 | 52 +++++++++++++++++++++++++ src/modules/quick_method_module.f90 | 10 ++--- src/modules/quick_scf_module.f90 | 24 +++++++++++- src/subs/hrr.f90 | 11 +++++- 4 files changed, 90 insertions(+), 7 deletions(-) diff --git a/src/modules/quick_calculated_module.f90 b/src/modules/quick_calculated_module.f90 index 01b791203..496497d87 100644 --- a/src/modules/quick_calculated_module.f90 +++ b/src/modules/quick_calculated_module.f90 @@ -167,6 +167,11 @@ module quick_calculated_module ! Mulliken charge and Lowdin charge double precision,dimension(:),allocatable :: Mulliken,Lowdin + ! ERI index for 8 fold symmetry + integer, dimension(:,:,:,:), allocatable :: iarray2 + + ! ERI in 8 fold symmetry + double precision,dimension(:),allocatable :: aoint2e end type quick_qm_struct_type @@ -222,6 +227,10 @@ subroutine allocate_quick_qm_struct(self) integer nelec integer idimA integer nelecb + !For QC output + integer i,j,k,l + integer ij,kl,ijkl + integer npair,ns8 type (quick_qm_struct_type) self @@ -247,6 +256,39 @@ subroutine allocate_quick_qm_struct(self) if(.not. allocated(self%Mulliken)) allocate(self%Mulliken(natom)) if(.not. allocated(self%Lowdin)) allocate(self%Lowdin(natom)) + if (quick_method%QCint) then + if(.not. allocated(self%iarray2)) allocate(self%iarray2(nbasis,nbasis,nbasis,nbasis)) + self%iarray2=0.d0 + ijkl=0 + ij=0 + do i=1,nbasis + do j=1,i + ij=ij+1 + kl=0 + do k=1,i + do l=1,k + kl=kl+1 + if (ij .ge. kl) then + ijkl=ijkl+1 + self%iarray2(i,j,k,l)=ijkl + self%iarray2(i,j,l,k)=ijkl + self%iarray2(j,i,l,k)=ijkl + self%iarray2(j,i,k,l)=ijkl + self%iarray2(k,l,i,j)=ijkl + self%iarray2(k,l,j,i)=ijkl + self%iarray2(l,k,j,i)=ijkl + self%iarray2(l,k,i,j)=ijkl + endif + enddo + enddo + enddo + enddo + + npair=(nbasis*(nbasis+1))/2 + ns8=(ns8*(ns8+1))/2 + if(.not. allocated(self%aoint2e)) allocate(self%aoint2e(ns8)) + endif + ! if 1st order derivation, which is gradient calculation is requested if (quick_method%grad) then if(.not. allocated(self%gradient)) allocate(self%gradient(3*natom)) @@ -410,6 +452,11 @@ subroutine deallocate_quick_qm_struct(self) if (allocated(self%Mulliken)) deallocate(self%Mulliken) if (allocated(self%Lowdin)) deallocate(self%Lowdin) + if (quick_method%QCint) then + if (allocated(self%iarray2)) deallocate(self%iarray2) + if (allocated(self%aoint2e)) deallocate(self%aoint2e) + endif + ! if 1st order derivation, which is gradient calculation is requested if (quick_method%grad) then if (allocated(self%gradient)) deallocate(self%gradient) @@ -568,6 +615,11 @@ subroutine init_quick_qm_struct(self) call zeroMatrix(self%denseOld,nbasis) + ! For QC output + if (quick_method%QCint) then + self%aoint2e=0.0d0 + endif + ! if 1st order derivation, which is gradient calculation is requested if (quick_method%grad) then call zeroVec(self%gradient,3*natom) diff --git a/src/modules/quick_method_module.f90 b/src/modules/quick_method_module.f90 index 2499d8ff2..df25c5019 100644 --- a/src/modules/quick_method_module.f90 +++ b/src/modules/quick_method_module.f90 @@ -79,7 +79,7 @@ module quick_method_module ! this part is Div&Con options logical :: BEoff = .false. - logical :: Qint = .false. + logical :: QCint = .false. logical :: OWNfrag = .false. ! this is DFT grid @@ -253,7 +253,7 @@ subroutine broadcast_quick_method(self, ierr) call MPI_BCAST(self%DIVCON,1,mpi_logical,0,MPI_COMM_WORLD,mpierror) call MPI_BCAST(self%DCMP2only,1,mpi_logical,0,MPI_COMM_WORLD,mpierror) call MPI_BCAST(self%BEoff,1,mpi_logical,0,MPI_COMM_WORLD,mpierror) - call MPI_BCAST(self%Qint,1,mpi_logical,0,MPI_COMM_WORLD,mpierror) + call MPI_BCAST(self%QCint,1,mpi_logical,0,MPI_COMM_WORLD,mpierror) call MPI_BCAST(self%OWNfrag,1,mpi_logical,0,MPI_COMM_WORLD,mpierror) call MPI_BCAST(self%DNCRB,1,mpi_double_precision,0,MPI_COMM_WORLD,mpierror) call MPI_BCAST(self%DNCRB2,1,mpi_double_precision,0,MPI_COMM_WORLD,mpierror) @@ -503,7 +503,7 @@ subroutine print_quick_method(self,io,ierr) ! additional DIVCON options if (self%BEoff) write(io,'(" NO ELIMINATION OF EMBEDDED SUBSYSTEM IN DIVCON ")') - if (self%Qint) write(io,'(" PRINT DIVCON FRAGMENT INTEGRALS IN QISKIT FORMAT ")') + if (self%QCint) write(io,'(" PRINT DIVCON FRAGMENT INTEGRALS IN QISKIT FORMAT ")') ! computing cycles write(io,'(" MAX SCF CYCLES = ",i6)') self%iscf @@ -672,7 +672,7 @@ subroutine read_quick_method(self,keywd,ierr) if (index(keyWD,'EXTCHARGES').ne.0) self%EXTCHARGES=.true. if (index(keyWD,'FORCE').ne.0) self%grad=.true. if (index(keyWD,'BEOFF').ne.0) self%BEoff=.true. - if (index(keyWD,'QINT').ne.0) self%Qint=.true. + if (index(keyWD,'QCINT').ne.0) self%QCint=.true. ! The BEoff keyword is redundant when OWNFRAG is used. ! OWNFRAG sets radii of both buffer1 and buffer2 to zero, @@ -928,7 +928,7 @@ subroutine init_quick_method(self,ierr) self%DIVCON = .false. ! Div&Con self%DCMP2only = .false. ! Canonical HF and Div&Con for MP2 only self%BEoff = .false. ! Turns off bEliminate in Div&Con - self%Qint = .false. ! Prints 1e and 2e integrals in Qiskit format + self%QCint = .false. ! Prints 1e and 2e integrals in Qiskit format self%OWNfrag = .false. ! User defined fragmentation based on CORE.DNC and BUFFER.DNC files self%ifragbasis = 1 ! =2.residue basis,=1.atom basis(DEFUALT),=3 non-h atom basis diff --git a/src/modules/quick_scf_module.f90 b/src/modules/quick_scf_module.f90 index 15ca70e1d..60bfc2ed3 100644 --- a/src/modules/quick_scf_module.f90 +++ b/src/modules/quick_scf_module.f90 @@ -113,8 +113,12 @@ subroutine scf(ierr) logical :: done integer, intent(inout) :: ierr - integer :: jscf + integer :: jscf,ierr1 done=.false. + + if (quick_method%QCint) then + call quick_open(iIntFile,intFileName,'U','U','A',.false.,ierr1) + end if !----------------------------------------------------------------- ! The purpose of this subroutine is to perform scf cycles. At this @@ -206,6 +210,9 @@ subroutine electdiis(jscf,ierr) double precision :: c_coords(3),c_zeta,c_chg + ! For QC output + integer :: npair,ns8,ijkl, ijkl2 + !--------------------------------------------------------------------------- ! The purpose of this subroutine is to utilize Pulay's accelerated ! scf convergence as detailed in J. Comp. Chem, Vol 3, #4, pg 566-60, 1982. @@ -356,6 +363,21 @@ subroutine electdiis(jscf,ierr) call scf_operator(deltaO) + if (quick_method%QCint) then + if (jscf .eq. 2) then + write(iintfile)'1E INTEGRALS' + do ijkl=1, nbasis + write(iintfile)(quick_qm_struct%oneElecO(ijkl,ijkl2),ijkl2=1,nbasis) + enddo + write(iintfile)'2E INTEGRALS' + npair = (nbasis*(nbasis+1))/2 + ns8 = (npair*(npair+1))/2 + do ijkl=1,ns8 + write(iintfile)quick_qm_struct%aoint2e(ijkl) + enddo + endif + endif + quick_qm_struct%denseOld(:,:) = quick_qm_struct%dense(:,:) if (quick_method%debug) call debug_SCF(jscf) diff --git a/src/subs/hrr.f90 b/src/subs/hrr.f90 index 0554db469..4cf8e93f4 100644 --- a/src/subs/hrr.f90 +++ b/src/subs/hrr.f90 @@ -22,6 +22,9 @@ subroutine hrrwhole common /xiaostore/store common /hrrstore/II,JJ,KK,LL,NBI1,NBI2,NBJ1,NBJ2,NBK1,NBK2,NBL1,NBL2 + !For QC output + integer ijkl + select case (IJKLtype) case (0,10,1000,1010) @@ -282,7 +285,13 @@ subroutine hrrwhole end select 111 continue -! write(*,*) IJKLtype,mpirank, iii,jjj,kkk,lll,Y + + if (quick_method%QCint) then + ijkl=quick_qm_struct%iarray2(iii,jjj,kkk,lll) + quick_qm_struct%aoint2e(ijkl)=Y +print*,IJKL,quick_qm_struct%aoint2e(ijkl) + endif + end subroutine hrrwhole