PTEROSOR/QuAcK
10
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-07-24 11:47:37 +02:00
Code Releases Activity

fix bug in GF2

Browse Source
This commit is contained in:
Pierre-Francois Loos 2019-04-29 09:52:21 +02:00
parent a963b22293
commit b93ee21573
2 changed files with 41 additions and 34 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

69
src/QuAcK/GF2_diag.f90
View File

@ -10,8 +10,13 @@ subroutine GF2_diag(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,V,e0)
integer,intent(in) :: maxSCF integer,intent(in) :: maxSCF
double precision,intent(in) :: thresh double precision,intent(in) :: thresh
integer,intent(in) :: max_diis integer,intent(in) :: max_diis
integer,intent(in) :: nBas,nC,nO,nV,nR integer,intent(in) :: nBas
double precision,intent(in) :: e0(nBas),V(nBas,nBas,nBas,nBas) integer,intent(in) :: nO
integer,intent(in) :: nC
integer,intent(in) :: nV
integer,intent(in) :: nR
double precision,intent(in) :: e0(nBas)
double precision,intent(in) :: V(nBas,nBas,nBas,nBas)
! Local variables ! Local variables
@ -20,7 +25,11 @@ subroutine GF2_diag(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,V,e0)
double precision :: eps double precision :: eps
double precision :: Conv double precision :: Conv
double precision :: rcond double precision :: rcond
double precision,allocatable :: eGF2(:),eOld(:),Bpp(:,:),error_diis(:,:),e_diis(:,:) double precision,allocatable :: eGF2(:)
double precision,allocatable :: eOld(:)
double precision,allocatable :: Bpp(:,:)
double precision,allocatable :: error_diis(:,:)
double precision,allocatable :: e_diis(:,:)
integer :: i,j,a,b,p integer :: i,j,a,b,p
@ -34,8 +43,7 @@ subroutine GF2_diag(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,V,e0)
! Memory allocation ! Memory allocation
allocate(Bpp(nBas,2),eGF2(nBas),eOld(nBas), & allocate(Bpp(nBas,2),eGF2(nBas),eOld(nBas),error_diis(nBas,max_diis),e_diis(nBas,max_diis))
error_diis(nBas,max_diis),e_diis(nBas,max_diis))
! Initialization ! Initialization
@ -59,39 +67,42 @@ subroutine GF2_diag(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,V,e0)
do p=nC+1,nBas-nR do p=nC+1,nBas-nR
do i=nC+1,nO do i=nC+1,nO
do j=nC+1,nO do j=nC+1,nO
do a=nO+1,nBas-nR do a=nO+1,nBas-nR
eps = eGF2(p) + e0(a) - e0(i) - e0(j) eps = eGF2(p) + e0(a) - e0(i) - e0(j)
Bpp(p,1) = Bpp(p,1) & Bpp(p,1) = Bpp(p,1) &
+ (2d0*V(p,a,i,j) - V(p,a,j,i))*V(p,a,i,j)/eps + (2d0*V(p,a,i,j) - V(p,a,j,i))*V(p,a,i,j)/eps
enddo end do
enddo end do
enddo end do
enddo end do
do p=nC+1,nBas-nR do p=nC+1,nBas-nR
do i=nC+1,nO do i=nC+1,nO
do a=nO+1,nBas-nR do a=nO+1,nBas-nR
do b=nO+1,nBas-nR do b=nO+1,nBas-nR
eps = eGF2(p) + e0(i) - e0(a) - e0(b) eps = eGF2(p) + e0(i) - e0(a) - e0(b)
Bpp(p,2) = Bpp(p,2) & Bpp(p,2) = Bpp(p,2) &
+ (2d0*V(p,i,a,b) - V(p,i,b,a))*V(p,i,a,b)/eps + (2d0*V(p,i,a,b) - V(p,i,b,a))*V(p,i,a,b)/eps
enddo end do
enddo end do
enddo end do
enddo end do
eGF2(:) = e0(:) & eGF2(:) = e0(:) + Bpp(:,1) + Bpp(:,2)
+ Bpp(:,1) + Bpp(:,2)
Conv = maxval(abs(eGF2 - eOld)) Conv = maxval(abs(eGF2 - eOld))
! Print results
call print_GF2(nBas,nO,nSCF,Conv,e0,eGF2)
! DIIS extrapolation ! DIIS extrapolation
n_diis = min(n_diis+1,max_diis) n_diis = min(n_diis+1,max_diis)
@ -99,17 +110,13 @@ subroutine GF2_diag(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,V,e0)
if(abs(rcond) < 1d-15) n_diis = 0 if(abs(rcond) < 1d-15) n_diis = 0
eOld = eGF2 eOld(:) = eGF2(:)
! Print results
call print_GF2(nBas,nO,nSCF,Conv,e0,eGF2)
! Increment ! Increment
nSCF = nSCF + 1 nSCF = nSCF + 1
enddo end do
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! End main SCF loop ! End main SCF loop
!------------------------------------------------------------------------ !------------------------------------------------------------------------
@ -124,6 +131,6 @@ subroutine GF2_diag(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,V,e0)
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!' write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*) write(*,*)
endif end if
end subroutine GF2_diag end subroutine GF2_diag

6
src/QuAcK/QuAcK.f90
View File

@ -359,7 +359,7 @@ program QuAcK
if(doADC) then if(doADC) then
call cpu_time(start_ADC) call cpu_time(start_ADC)
call ADC(singlet_manifold,triplet_manifold,maxSCF_GF,thresh_GF,n_diis_GF,nBas,nC,nO,nV,nR,eHF,ERI_MO_basis) call ADC(singlet_manifold,triplet_manifold,maxSCF_GF,thresh_GF,n_diis_GF,nBas,nC(1),nO(1),nV(1),nR(1),eHF,ERI_MO_basis)
call cpu_time(end_ADC) call cpu_time(end_ADC)
t_ADC = end_ADC - start_ADC t_ADC = end_ADC - start_ADC
@ -375,7 +375,7 @@ program QuAcK
if(doGF2) then if(doGF2) then
call cpu_time(start_GF2) call cpu_time(start_GF2)
call GF2_diag(maxSCF_GF,thresh_GF,n_diis_GF,nBas,nC,nO,nV,nR,ERI_MO_basis,eHF) call GF2_diag(maxSCF_GF,thresh_GF,n_diis_GF,nBas,nC(1),nO(1),nV(1),nR(1),ERI_MO_basis,eHF)
call cpu_time(end_GF2) call cpu_time(end_GF2)
t_GF2 = end_GF2 - start_GF2 t_GF2 = end_GF2 - start_GF2
@ -391,7 +391,7 @@ program QuAcK
if(doGF3) then if(doGF3) then
call cpu_time(start_GF3) call cpu_time(start_GF3)
call GF3_diag(maxSCF_GF,thresh_GF,n_diis_GF,renormalization,nBas,nC,nO,nV,nR,ERI_MO_basis,eHF) call GF3_diag(maxSCF_GF,thresh_GF,n_diis_GF,renormalization,nBas,nC(1),nO(1),nV(1),nR(1),ERI_MO_basis,eHF)
call cpu_time(end_GF3) call cpu_time(end_GF3)
t_GF3 = end_GF3 - start_GF3 t_GF3 = end_GF3 - start_GF3