LCPQ
/
quack
mirror of https://github.com/pfloos/quack
4
1
Fork 0
Code Releases Activity

include test in many different things... work on progress

This commit is contained in:
Pierre-Francois Loos 2023-11-13 16:15:15 +01:00
parent ae21a778c1
commit 8a5e094278
39 changed files with 477 additions and 213 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

14
input/methods
View File

@ -1,20 +1,20 @@
# RHF UHF GHF ROHF
T T T T
T F T F
# MP2 MP3
T T
# CCD pCCD DCD CCSD CCSD(T)
F F F F F
T T T T F
# drCCD rCCD crCCD lCCD
F F F F
T T T T
# CIS CIS(D) CID CISD FCI
F F F F F
T F F F F
# phRPA phRPAx crRPA ppRPA
T T F T
T T T T
# G0F2 evGF2 qsGF2 G0F3 evGF3
F F F F F
T T F F F
# G0W0 evGW qsGW SRG-qsGW ufG0W0 ufGW
F F F F F F
# G0T0pp evGTpp qsGTpp G0T0eh evGTeh qsGTeh
F F F F F F
# Rtest Utest Gtest
T T T
T F T

4
input/options
View File

@ -4,8 +4,8 @@
F
# CC: maxSCF thresh DIIS
64 0.0000001 5
# spin: TDA spin_conserved spin_flip
F T T
# spin: TDA singlet triplet
F T T
# GF: maxSCF thresh DIIS lin eta renorm reg
256 0.00001 5 F 0.0 0 F
# GW: maxSCF thresh DIIS lin eta TDA_W reg

24
src/CC/CCD.f90
View File

@ -1,4 +1,4 @@
subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
subroutine CCD(dotest,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! CCD module
@ -6,6 +6,8 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -29,7 +31,7 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e
integer :: nSCF
double precision :: Conv
double precision :: EcMP2,EcMP3,EcMP4
double precision :: ECCD,EcCCD
double precision :: ECC,EcCC
double precision,allocatable :: seHF(:)
double precision,allocatable :: sERI(:,:,:,:)
double precision,allocatable :: dbERI(:,:,:,:)
@ -191,13 +193,13 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t,EcCCD)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t,EcCC)
if(nSCF == 1) call MP3_correlation_energy(nC,nO,nV,nR,OOVV,t,v,delta_OOVV,EcMP3)
! Dump results
ECCD = ERHF + EcCCD
ECC = ERHF + EcCC
! DIIS extrapolation
@ -209,7 +211,7 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e
if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
@ -235,8 +237,8 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e
write(*,*)'----------------------------------------------------'
write(*,*)' CCD energy '
write(*,*)'----------------------------------------------------'
write(*,'(1X,A30,1X,F15.10)')' E(CCD) = ',ECCD
write(*,'(1X,A30,1X,F15.10)')' Ec(CCD) = ',EcCCD
write(*,'(1X,A30,1X,F15.10)')' E(CCD) = ',ECC
write(*,'(1X,A30,1X,F15.10)')' Ec(CCD) = ',EcCC
write(*,*)'----------------------------------------------------'
write(*,*)
@ -272,4 +274,12 @@ subroutine CCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,e
if(do_DIP_EOM_CC_2h) call DIP_EOM_CCD_2h(nC,nO,nV,nR,eO,OOVV,OOOO,t)
! Testing zone
if(dotest) then
call dump_test_value('R','CCD correlation energy',EcCC)
end if
end subroutine

12
src/CC/CCSD.f90
View File

@ -1,4 +1,4 @@
subroutine CCSD(maxSCF,thresh,max_diis,doCCSDT,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
subroutine CCSD(dotest,maxSCF,thresh,max_diis,doCCSDT,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! CCSD module
@ -6,6 +6,8 @@ subroutine CCSD(maxSCF,thresh,max_diis,doCCSDT,nBasin,nCin,nOin,nVin,nRin,ERI,EN
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -302,4 +304,12 @@ subroutine CCSD(maxSCF,thresh,max_diis,doCCSDT,nBasin,nCin,nOin,nVin,nRin,ERI,EN
end if
! Testing zone
if(dotest) then
call dump_test_value('R','CCSD correlation energy',EcCCSD)
end if
end subroutine

22
src/CC/DCD.f90
View File

@ -1,4 +1,4 @@
subroutine DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
subroutine DCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! DCD module
@ -6,6 +6,8 @@ subroutine DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -24,7 +26,7 @@ subroutine DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
integer :: nSCF
double precision :: Conv
double precision :: EcMP2,EcMP3,EcMP4
double precision :: EDCD,EcDCD
double precision :: ECC,EcCC
double precision,allocatable :: eO(:)
double precision,allocatable :: eV(:)
@ -123,13 +125,13 @@ subroutine DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! Compute correlation energy
EcDCD = 0d0
EcCC = 0d0
do i=1,nO-nC
do j=1,nO-nC
do a=1,nV-nR
do b=1,nV-nR
EcDCD = EcDCD + (2d0*OOVV(i,j,a,b) - OOVV(i,j,b,a))*t(i,j,a,b)
EcCC = EcCC + (2d0*OOVV(i,j,a,b) - OOVV(i,j,b,a))*t(i,j,a,b)
enddo
enddo
@ -138,10 +140,10 @@ subroutine DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! Dump results
EDCD = ERHF + EcDCD
ECC = ERHF + EcCC
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',EDCD+ENuc,'|',EcDCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
! Increment
@ -280,4 +282,12 @@ subroutine DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
endif
! Testing zone
if(dotest) then
call dump_test_value('R','DCD correlation energy',EcCC)
end if
end subroutine

20
src/CC/RCC.f90
View File

@ -1,4 +1,4 @@
subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD, &
subroutine RCC(dotest,doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD, &
maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
! Coupled-cluster module
@ -8,6 +8,8 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
! Input variables
logical :: dotest
logical :: doCCD
logical :: dopCCD
logical :: doDCD
@ -43,7 +45,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(doCCD) then
call wall_time(start_CC)
call CCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call CCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)
t_CC = end_CC - start_CC
@ -59,7 +61,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(doDCD) then
call wall_time(start_CC)
call DCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR, &
call DCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR, &
ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)
@ -78,7 +80,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(doCCSD) then
call wall_time(start_CC)
call CCSD(maxSCF,thresh,max_diis,doCCSDT,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call CCSD(dotest,maxSCF,thresh,max_diis,doCCSDT,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)
t_CC = end_CC - start_CC
@ -94,7 +96,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(dodrCCD) then
call wall_time(start_CC)
call drCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call drCCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)
t_CC = end_CC - start_CC
@ -110,7 +112,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(dorCCD) then
call wall_time(start_CC)
call rCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF,epsHF)
call rCCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF,epsHF)
call wall_time(end_CC)
t_CC = end_CC - start_CC
@ -126,7 +128,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(docrCCD) then
call wall_time(start_CC)
call crCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call crCCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)
t_CC = end_CC - start_CC
@ -142,7 +144,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(dolCCD) then
call wall_time(start_CC)
call lCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call lCCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)
t_CC = end_CC - start_CC
@ -158,7 +160,7 @@ subroutine RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD,
if(dopCCD) then
call wall_time(start_CC)
call pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call pCCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,epsHF)
call wall_time(end_CC)

22
src/CC/crCCD.f90
View File

@ -1,4 +1,4 @@
subroutine crCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
subroutine crCCD(dotest,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! Crossed-ring CCD module
@ -6,6 +6,8 @@ subroutine crCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -29,7 +31,7 @@ subroutine crCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
integer :: nSCF
double precision :: Conv
double precision :: EcMP2
double precision :: ECCD,EcCCD
double precision :: ECC,EcCC
double precision,allocatable :: seHF(:)
double precision,allocatable :: sERI(:,:,:,:)
double precision,allocatable :: dbERI(:,:,:,:)
@ -155,11 +157,11 @@ subroutine crCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCC)
! Dump results
ECCD = ERHF + EcCCD
ECC = ERHF + EcCC
! DIIS extrapolation
@ -171,7 +173,7 @@ subroutine crCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
@ -197,9 +199,15 @@ subroutine crCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
write(*,*)'----------------------------------------------------'
write(*,*)' crossed-ring CCD energy '
write(*,*)'----------------------------------------------------'
write(*,'(1X,A30,1X,F15.10)')' E(crCCD) = ',ECCD
write(*,'(1X,A30,1X,F15.10)')' Ec(crCCD) = ',EcCCD
write(*,'(1X,A30,1X,F15.10)')' E(crCCD) = ',ECC
write(*,'(1X,A30,1X,F15.10)')' Ec(crCCD) = ',EcCC
write(*,*)'----------------------------------------------------'
write(*,*)
if(dotest) then
call dump_test_value('R','crCCD correlation energy',EcCC)
end if
end subroutine

24
src/CC/drCCD.f90
View File

@ -1,4 +1,4 @@
subroutine drCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
subroutine drCCD(dotest,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! Direct ring CCD module
@ -6,6 +6,8 @@ subroutine drCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -29,7 +31,7 @@ subroutine drCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
integer :: nSCF
double precision :: Conv
double precision :: EcMP2
double precision :: ECCD,EcCCD
double precision :: ECC,EcCC
double precision,allocatable :: seHF(:)
double precision,allocatable :: sERI(:,:,:,:)
@ -147,12 +149,12 @@ subroutine drCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD)
EcCCD = 2d0*EcCCD
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCC)
EcCC = 2d0*EcCC
! Dump results
ECCD = ERHF + EcCCD
ECC = ERHF + EcCC
! DIIS extrapolation
@ -164,7 +166,7 @@ subroutine drCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
@ -190,9 +192,15 @@ subroutine drCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF
write(*,*)'----------------------------------------------------'
write(*,*)' direct ring CCD energy '
write(*,*)'----------------------------------------------------'
write(*,'(1X,A30,1X,F15.10)')' E(drCCD) = ',ECCD
write(*,'(1X,A30,1X,F15.10)')' Ec(drCCD) = ',EcCCD
write(*,'(1X,A30,1X,F15.10)')' E(drCCD) = ',ECC
write(*,'(1X,A30,1X,F15.10)')' Ec(drCCD) = ',EcCC
write(*,*)'----------------------------------------------------'
write(*,*)
if(dotest) then
call dump_test_value('R','drCCD correlation energy',EcCC)
end if
end subroutine

22
src/CC/lCCD.f90
View File

@ -1,4 +1,4 @@
subroutine lCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
subroutine lCCD(dotest,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF)
! Ladder CCD module
@ -6,6 +6,8 @@ subroutine lCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -29,7 +31,7 @@ subroutine lCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
integer :: nSCF
double precision :: Conv
double precision :: EcMP2
double precision :: ECCD,EcCCD
double precision :: ECC,EcCC
double precision,allocatable :: seHF(:)
double precision,allocatable :: sERI(:,:,:,:)
double precision,allocatable :: dbERI(:,:,:,:)
@ -168,11 +170,11 @@ subroutine lCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCCD)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t2,EcCC)
! Dump results
ECCD = ERHF + EcCCD
ECC = ERHF + EcCC
! DIIS extrapolation
@ -184,7 +186,7 @@ subroutine lCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
@ -210,9 +212,15 @@ subroutine lCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
write(*,*)'----------------------------------------------------'
write(*,*)' ladder CCD energy '
write(*,*)'----------------------------------------------------'
write(*,'(1X,A30,1X,F15.10)')' E(lCCD) = ',ECCD
write(*,'(1X,A30,1X,F15.10)')' Ec(lCCD) = ',EcCCD
write(*,'(1X,A30,1X,F15.10)')' E(lCCD) = ',ECC
write(*,'(1X,A30,1X,F15.10)')' Ec(lCCD) = ',EcCC
write(*,*)'----------------------------------------------------'
write(*,*)
if(dotest) then
call dump_test_value('R','lCCD correlation energy',EcCC)
end if
end subroutine

24
src/CC/pCCD.f90
View File

@ -1,4 +1,4 @@
subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
subroutine pCCD(dotest,maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! pair CCD module
@ -6,6 +6,8 @@ subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -21,7 +23,7 @@ subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
integer :: nSCF
double precision :: Conv
double precision :: ECCD,EcCCD
double precision :: ECC,EcCC
double precision,allocatable :: eO(:)
double precision,allocatable :: eV(:)
@ -91,10 +93,10 @@ subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
t(:,:) = -0.5d0*OOVV(:,:)/delta_OV(:,:)
EcCCD = 0d0
EcCC = 0d0
do i=1,nO-nC
do a=1,nV-nR
EcCCD = EcCCD + OOVV(i,a)*t(i,a)
EcCC = EcCC + OOVV(i,a)*t(i,a)
end do
end do
@ -170,16 +172,16 @@ subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! Compute correlation energy
EcCCD = 0d0
EcCC = 0d0
do i=1,nO-nC
do a=1,nV-nR
EcCCD = EcCCD + OOVV(i,a)*t(i,a)
EcCC = EcCC + OOVV(i,a)*t(i,a)
end do
end do
! Dump results
ECCD = ERHF + EcCCD
ECC = ERHF + EcCC
! DIIS extrapolation
@ -191,7 +193,7 @@ subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
! if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
@ -212,4 +214,10 @@ subroutine pCCD(maxSCF,thresh,max_diis,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eHF)
endif
if(dotest) then
call dump_test_value('R','pCCD correlation energy',EcCC)
end if
end subroutine

23
src/CC/rCCD.f90
View File

@ -1,4 +1,4 @@
subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF,eGW)
subroutine rCCD(dotest,maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,eHF,eGW)
! Ring CCD module
@ -6,6 +6,8 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -31,7 +33,7 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
integer :: nSCF
double precision :: Conv
double precision :: EcMP2
double precision :: ECCD,EcCCD
double precision :: ECC,EcCC
double precision,allocatable :: seHF(:)
double precision,allocatable :: seGW(:)
double precision,allocatable :: sERI(:,:,:,:)
@ -161,11 +163,11 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
! Compute correlation energy
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t,EcCCD)
call CCD_correlation_energy(nC,nO,nV,nR,OOVV,t,EcCC)
! Dump results
ECCD = ERHF + EcCCD
ECC = ERHF + EcCC
! DIIS extrapolation
@ -177,7 +179,7 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
if(abs(rcond) < 1d-15) n_diis = 0
write(*,'(1X,A1,1X,I3,1X,A1,1X,F16.10,1X,A1,1X,F10.6,1X,A1,1X,F10.6,1X,A1,1X)') &
'|',nSCF,'|',ECCD+ENuc,'|',EcCCD,'|',Conv,'|'
'|',nSCF,'|',ECC+ENuc,'|',EcCC,'|',Conv,'|'
enddo
write(*,*)'----------------------------------------------------'
@ -203,8 +205,8 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
write(*,*)'----------------------------------------------------'
write(*,*)' ring CCD energy '
write(*,*)'----------------------------------------------------'
write(*,'(1X,A30,1X,F15.10)')' E(rCCD) = ',ECCD
write(*,'(1X,A30,1X,F15.10)')' Ec(rCCD) = ',EcCCD
write(*,'(1X,A30,1X,F15.10)')' E(rCCD) = ',ECC
write(*,'(1X,A30,1X,F15.10)')' Ec(rCCD) = ',EcCC
write(*,*)'----------------------------------------------------'
write(*,*)
@ -223,4 +225,11 @@ subroutine rCCD(maxSCF,thresh,max_diis,nBasin,nCin,nOin,nVin,nRin,ERI,ENuc,ERHF,
if(do_EE_EOM_CC_1h1p) call EE_EOM_CCD_1h1p(nC,nO,nV,nR,eO,eV,OOVV,OVVO,t)
if(dotest) then
call dump_test_value('R','rCCD correlation energy',EcCC)
end if
end subroutine

8
src/CI/CID.f90
View File

@ -1,4 +1,4 @@
subroutine CID(singlet_manifold,triplet_manifold,nBasin,nCin,nOin,nVin,nRin,ERIin,eHFin,E0)
subroutine CID(dotest,singlet,triplet,nBasin,nCin,nOin,nVin,nRin,ERIin,eHFin,E0)
! Perform configuration interaction with doubles
@ -7,8 +7,10 @@ subroutine CID(singlet_manifold,triplet_manifold,nBasin,nCin,nOin,nVin,nRin,ERIi
! Input variables
logical,intent(in) :: singlet_manifold
logical,intent(in) :: triplet_manifold
logical,intent(in) :: dotest
logical,intent(in) :: singlet
logical,intent(in) :: triplet
integer,intent(in) :: nBasin
integer,intent(in) :: nCin
integer,intent(in) :: nOin

8
src/CI/CISD.f90
View File

@ -1,4 +1,4 @@
subroutine CISD(singlet_manifold,triplet_manifold,nBasin,nCin,nOin,nVin,nRin,ERIin,eHFin,E0)
subroutine CISD(dotest,singlet,triplet,nBasin,nCin,nOin,nVin,nRin,ERIin,eHFin,E0)
! Perform configuration interaction with singles and doubles
@ -7,8 +7,10 @@ subroutine CISD(singlet_manifold,triplet_manifold,nBasin,nCin,nOin,nVin,nRin,ERI
! Input variables
logical,intent(in) :: singlet_manifold
logical,intent(in) :: triplet_manifold
logical,intent(in) :: dotest
logical,intent(in) :: singlet
logical,intent(in) :: triplet
integer,intent(in) :: nBasin
integer,intent(in) :: nCin
integer,intent(in) :: nOin

20
src/CI/RCI.f90
View File

@ -1,4 +1,4 @@
subroutine RCI(doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR,nS,ERI,dipole_int, &
subroutine RCI(dotest,doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR,nS,ERI,dipole_int, &
epsHF,EHF,cHF,S)
! Configuration interaction module
@ -8,11 +8,13 @@ subroutine RCI(doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR
! Input variables
logical :: doCIS
logical :: doCIS_D
logical :: doCID
logical :: doCISD
logical :: doFCI
logical,intent(in) :: dotest
logical,intent(in) :: doCIS
logical,intent(in) :: doCIS_D
logical,intent(in) :: doCID
logical,intent(in) :: doCISD
logical,intent(in) :: doFCI
logical,intent(in) :: singlet
logical,intent(in) :: triplet
@ -40,7 +42,7 @@ subroutine RCI(doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR
if(doCIS) then
call wall_time(start_CI)
call CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,epsHF)
call RCIS(dotest,singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,epsHF)
call wall_time(end_CI)
t_CI = end_CI - start_CI
@ -56,7 +58,7 @@ subroutine RCI(doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR
if(doCID) then
call wall_time(start_CI)
call CID(singlet,triplet,nBas,nC,nO,nV,nR,ERI,epsHF,EHF)
call CID(dotest,singlet,triplet,nBas,nC,nO,nV,nR,ERI,epsHF,EHF)
call wall_time(end_CI)
t_CI = end_CI - start_CI
@ -72,7 +74,7 @@ subroutine RCI(doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR
if(doCISD) then
call wall_time(start_CI)
call CISD(singlet,triplet,nBas,nC,nO,nV,nR,ERI,epsHF,EHF)
call CISD(dotest,singlet,triplet,nBas,nC,nO,nV,nR,ERI,epsHF,EHF)
call wall_time(end_CI)
t_CI = end_CI - start_CI

20
src/CI/CIS.f90 → src/CI/RCIS.f90
View File

@ -1,4 +1,4 @@
subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
subroutine RCIS(dotest,singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
! Perform configuration interaction single calculation`
@ -7,6 +7,8 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: singlet
logical,intent(in) :: triplet
logical,intent(in) :: doCIS_D
@ -69,6 +71,14 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
maxS = min(maxS,nS)
if(doCIS_D) call CIS_D(ispin,nBas,nC,nO,nV,nR,nS,maxS,eHF,ERI,Om(1:maxS),A(:,1:maxS))
! Testing zone
if(dotest) then
call dump_test_value('R','RCIS singlet excitation energy',Om(1))
end if
endif
if(triplet) then
@ -97,6 +107,14 @@ subroutine CIS(singlet,triplet,doCIS_D,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF)
maxS = min(maxS,nS)
if(doCIS_D) call CIS_D(ispin,nBas,nC,nO,nV,nR,nS,maxS,eHF,ERI,Om(1:maxS),A(:,1:maxS))
! Testing zone
if(dotest) then
call dump_test_value('R','RCIS triplet excitation energy',Om(1))
end if
endif
end subroutine

16
src/CI/UCI.f90
View File

@ -1,4 +1,4 @@
subroutine UCI(doCIS,doCIS_D,doCID,doCISD,doFCI,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS, &
subroutine UCI(dotest,doCIS,doCIS_D,doCID,doCISD,doFCI,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,epsHF,EHF,cHF,S,F)
! Configuration interaction module
@ -8,11 +8,13 @@ subroutine UCI(doCIS,doCIS_D,doCID,doCISD,doFCI,spin_conserved,spin_flip,nBas,nC
! Input variables
logical :: doCIS
logical :: doCIS_D
logical :: doCID
logical :: doCISD
logical :: doFCI
logical,intent(in) :: dotest
logical,intent(in) :: doCIS
logical,intent(in) :: doCIS_D
logical,intent(in) :: doCID
logical,intent(in) :: doCISD
logical,intent(in) :: doFCI
logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip
@ -44,7 +46,7 @@ subroutine UCI(doCIS,doCIS_D,doCID,doCISD,doFCI,spin_conserved,spin_flip,nBas,nC
if(doCIS) then
call wall_time(start_CI)
call UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb, &
call UCIS(dotest,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb, &
ERI_bbbb,dipole_int_aa,dipole_int_bb,epsHF,cHF,S)
call wall_time(end_CI)

20
src/CI/UCIS.f90
View File

@ -1,4 +1,4 @@
subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb, &
subroutine UCIS(dotest,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,eHF,cHF,S)
! Perform configuration interaction single calculation`
@ -8,6 +8,8 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip
integer,intent(in) :: nBas
@ -88,6 +90,14 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
write(*,*)
endif
! Testing zone
if(dotest) then
call dump_test_value('U','UCIS singlet excitation energy',Om_sc(1))
end if
deallocate(A_sc,Om_sc)
endif
@ -128,6 +138,14 @@ subroutine UCIS(spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,E
write(*,*)
endif
! Testing zone
if(dotest) then
call dump_test_value('U','UCIS triplet excitation energy',Om_sf(1))
end if
deallocate(A_sf,Om_sf)
endif

23
src/GF/GG0F2.f90
View File

@ -1,4 +1,4 @@
subroutine GG0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
subroutine GG0F2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
! Perform a one-shot second-order Green function calculation
@ -8,6 +8,8 @@ subroutine GG0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: dophBSE
logical,intent(in) :: doppBSE
logical,intent(in) :: TDA
@ -39,10 +41,11 @@ subroutine GG0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
! Hello world
write(*,*)
write(*,*)'************************************************'
write(*,*)'| One-shot second-order Green function |'
write(*,*)'************************************************'
write(*,*)'********************************'
write(*,*)'* Generalized G0F2 Calculation *'
write(*,*)'********************************'
write(*,*)
! Memory allocation
@ -80,7 +83,7 @@ subroutine GG0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
! Print results
call GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,eGF,Ec)
call GMP2(.false.,regularize,nBas,nC,nO,nV,nR,ERI,ENuc,EHF,eGF,Ec)
call print_G0F2(nBas,nO,eHF,SigC,eGF,Z,ENuc,ERHF,Ec)
! Perform BSE2 calculation
@ -113,4 +116,14 @@ subroutine GG0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
! end if
! Testing zone
if(dotest) then
call dump_test_value('G','GG0F2 correlation energy',Ec)
call dump_test_value('G','GG0F2 HOMO energy',eGF(nO))
call dump_test_value('G','GG0F2 LUMO energy',eGF(nO+1))
end if
end subroutine

16
src/GF/GGF.f90
View File

@ -1,4 +1,4 @@
subroutine GGF(doG0F2,doevGF2,doqsGF2,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
subroutine GGF(dotest,doG0F2,doevGF2,doqsGF2,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI,dipole_int_AO,dipole_int,PHF,cHF,epsHF)
! Green's function module
@ -8,9 +8,11 @@ subroutine GGF(doG0F2,doevGF2,doqsGF2,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA
! Input variables
logical :: doG0F2
logical :: doevGF2
logical :: doqsGF2
logical,intent(in) :: dotest
logical,intent(in) :: doG0F2
logical,intent(in) :: doevGF2
logical,intent(in) :: doqsGF2
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
@ -61,7 +63,7 @@ subroutine GGF(doG0F2,doevGF2,doqsGF2,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA
if(doG0F2) then
call wall_time(start_GF)
call GG0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
call GG0F2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,epsHF)
call wall_time(end_GF)
@ -78,7 +80,7 @@ subroutine GGF(doG0F2,doevGF2,doqsGF2,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA
if(doevGF2) then
call wall_time(start_GF)
call evGGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,linearize,eta,regularize, &
call evGGF2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,epsHF)
call wall_time(end_GF)
@ -95,7 +97,7 @@ subroutine GGF(doG0F2,doevGF2,doqsGF2,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA
if(doqsGF2) then
call wall_time(start_GF)
! call qsGGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,eta,regularize,nNuc,ZNuc,rNuc,ENuc, &
! call qsGGF2(dotest,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,eta,regularize,nNuc,ZNuc,rNuc,ENuc, &
! nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI,dipole_int_AO,dipole_int,PHF,cHF,epsHF)
call wall_time(end_GF)

36
src/GF/G0F2.f90 → src/GF/RG0F2.f90
View File

@ -1,5 +1,5 @@
subroutine G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
subroutine RG0F2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
! Perform a one-shot second-order Green function calculation
@ -8,6 +8,8 @@ subroutine G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regu
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: dophBSE
logical,intent(in) :: doppBSE
logical,intent(in) :: TDA
@ -42,9 +44,9 @@ subroutine G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regu
! Hello world
write(*,*)
write(*,*)'************************************************'
write(*,*)'| One-shot second-order Green function |'
write(*,*)'************************************************'
write(*,*)'*******************************'
write(*,*)'* Restricted G0F2 Calculation *'
write(*,*)'*******************************'
write(*,*)
! Memory allocation
@ -82,7 +84,7 @@ subroutine G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regu
! Print results
call RMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eGF,Ec)
call RMP2(.false.,regularize,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eGF,Ec)
call print_G0F2(nBas,nO,eHF,SigC,eGF,Z,ENuc,ERHF,Ec)
! Perform BSE2 calculation
@ -95,8 +97,8 @@ subroutine G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regu
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(2X,A50,F20.10)') 'Tr@phBSE@G0F2 correlation energy (singlet) =',EcBSE(1)
write(*,'(2X,A50,F20.10)') 'Tr@phBSE@G0F2 correlation energy (triplet) =',EcBSE(2)
write(*,'(2X,A50,F20.10)') 'Tr@phBSE@G0F2 correlation energy =',sum(EcBSE(:))
write(*,'(2X,A50,F20.10)') 'Tr@phBSE@G0F2 total energy =',ENuc + ERHF + sum(EcBSE(:))
write(*,'(2X,A50,F20.10)') 'Tr@phBSE@G0F2 correlation energy =',sum(EcBSE)
write(*,'(2X,A50,F20.10)') 'Tr@phBSE@G0F2 total energy =',ENuc + ERHF + sum(EcBSE)
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
@ -108,15 +110,27 @@ subroutine G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regu
call GF2_ppBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,ERI,dipole_int,eGF,EcBSE)
EcBSE(2) = 3d0*EcBSE(2)
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 correlation energy (singlet) =',EcBSE(1),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 correlation energy (triplet) =',3d0*EcBSE(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 correlation energy =',EcBSE(1) + 3d0*EcBSE(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 total energy =',ENuc + ERHF + EcBSE(1) + 3d0*EcBSE(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 correlation energy (triplet) =',EcBSE(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 correlation energy =',sum(EcBSE),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@ppBSE@G0F2 total energy =',ENuc + ERHF + sum(EcBSE),' au'
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
end if
! Testing zone
if(dotest) then
call dump_test_value('R','RG0F2 correlation energy',Ec)
call dump_test_value('R','RG0F2 HOMO energy',eGF(nO))
call dump_test_value('R','RG0F2 LUMO energy',eGF(nO+1))
end if
end subroutine

3
src/GF/G0F3.f90 → src/GF/RG0F3.f90
View File

@ -1,4 +1,4 @@
subroutine G0F3(renormalization,nBas,nC,nO,nV,nR,V,e0)
subroutine RG0F3(dotest,renormalization,nBas,nC,nO,nV,nR,V,e0)
! Perform third-order Green function calculation in diagonal approximation
@ -7,6 +7,7 @@
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: renormalization
integer,intent(in) :: nBas,nC,nO,nV,nR
double precision,intent(in) :: e0(nBas),V(nBas,nBas,nBas,nBas)

32
src/GF/RGF.f90
View File

@ -1,4 +1,4 @@
subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_diis, &
subroutine RGF(dotest,doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_diis, &
dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regularize, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI, &
dipole_int_AO,dipole_int,PHF,cHF,epsHF)
@ -10,11 +10,13 @@ subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
! Input variables
logical :: doG0F2
logical :: doevGF2
logical :: doqsGF2
logical :: doG0F3
logical :: doevGF3
logical,intent(in) :: dotest
logical,intent(in) :: doG0F2
logical,intent(in) :: doevGF2
logical,intent(in) :: doqsGF2
logical,intent(in) :: doG0F3
logical,intent(in) :: doevGF3
integer :: renorm
integer,intent(in) :: maxSCF
@ -68,8 +70,8 @@ subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doG0F2) then
call wall_time(start_GF)
call G0F2(dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,epsHF)
call RG0F2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,epsHF)
call wall_time(end_GF)
t_GF = end_GF - start_GF
@ -85,9 +87,9 @@ subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doevGF2) then
call wall_time(start_GF)
call evGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF, &
ERI,dipole_int,epsHF)
call evRGF2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF, &
ERI,dipole_int,epsHF)
call wall_time(end_GF)
t_GF = end_GF - start_GF
@ -103,8 +105,8 @@ subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doqsGF2) then
call wall_time(start_GF)
call qsGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc, &
nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI,dipole_int_AO,dipole_int,PHF,cHF,epsHF)
call qsRGF2(dotest,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc, &
nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI,dipole_int_AO,dipole_int,PHF,cHF,epsHF)
call wall_time(end_GF)
t_GF = end_GF - start_GF
@ -120,7 +122,7 @@ subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doG0F3) then
call wall_time(start_GF)
call G0F3(renorm,nBas,nC,nO,nV,nR,ERI,epsHF)
call RG0F3(dotest,renorm,nBas,nC,nO,nV,nR,ERI,epsHF)
call wall_time(end_GF)
t_GF = end_GF - start_GF
@ -136,7 +138,7 @@ subroutine RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doevGF3) then
call wall_time(start_GF)
call evGF3(maxSCF,thresh,max_diis,renorm,nBas,nC,nO,nV,nR,ERI,epsHF)
call evRGF3(dotest,maxSCF,thresh,max_diis,renorm,nBas,nC,nO,nV,nR,ERI,epsHF)
call wall_time(end_GF)
t_GF = end_GF - start_GF

37
src/GF/UG0F2.f90
View File

@ -1,4 +1,4 @@
subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, &
subroutine UG0F2(dotest,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,eHF)
! Perform unrestricted G0W0 calculation
@ -9,6 +9,8 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regula
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: BSE
logical,intent(in) :: TDA
logical,intent(in) :: dBSE
@ -44,18 +46,17 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regula
double precision,allocatable :: Z(:,:)
integer :: nS_aa,nS_bb,nS_sc
double precision,allocatable :: eGF2lin(:,:)
double precision,allocatable :: eGF2(:,:)
double precision,allocatable :: eGFlin(:,:)
double precision,allocatable :: eGF(:,:)
! Output variables
! Hello world
write(*,*)
write(*,*)'************************************************'
write(*,*)'| One-shot G0F2 calculation |'
write(*,*)'| *** Unrestricted version *** |'
write(*,*)'************************************************'
write(*,*)'*********************************'
write(*,*)'* Unrestricted G0F2 Calculation *'
write(*,*)'*********************************'
write(*,*)
! TDA
@ -71,7 +72,7 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regula
nS_bb = nS(2)
nS_sc = nS_aa + nS_bb
allocate(SigC(nBas,nspin),Z(nBas,nspin),eGF2(nBas,nspin),eGF2lin(nBas,nspin))
allocate(SigC(nBas,nspin),Z(nBas,nspin),eGF(nBas,nspin),eGFlin(nBas,nspin))
!---------------------!
! Compute self-energy !
@ -91,14 +92,14 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regula
! Solve the quasi-particle equation !
!-----------------------------------!
eGF2lin(:,:) = eHF(:,:) + Z(:,:)*SigC(:,:)
eGFlin(:,:) = eHF(:,:) + Z(:,:)*SigC(:,:)
if(linearize) then
write(*,*) ' *** Quasiparticle energies obtained by linearization *** '
write(*,*)
eGF2(:,:) = eGF2lin(:,:)
eGF(:,:) = eGFlin(:,:)
else
@ -111,11 +112,11 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regula
! Compute MP2 correlation energy
call UMP2(nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EUHF,eGF2,Ec)
call UMP2(.false.,nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EUHF,eGF,Ec)
! Dump results
call print_UG0F2(nBas,nO,eHF,ENuc,EUHF,SigC,Z,eGF2,Ec)
call print_UG0F2(nBas,nO,eHF,ENuc,EUHF,SigC,Z,eGF,Ec)
! Perform BSE calculation
@ -125,4 +126,16 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regula
end if
! Testing zone
if(dotest) then
call dump_test_value('U','UG0F2 correlation energy',Ec)
call dump_test_value('U','UG0F2 HOMOa energy',eGF(nO(1),1))
call dump_test_value('U','UG0F2 LUMOa energy',eGF(nO(1)+1,1))
call dump_test_value('U','UG0F2 HOMOa energy',eGF(nO(2),2))
call dump_test_value('U','UG0F2 LUMOa energy',eGF(nO(2)+1,2))
end if
end subroutine

20
src/GF/UGF.f90
View File

@ -1,4 +1,4 @@
subroutine UGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_diis, &
subroutine UGF(dotest,doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_diis, &
dophBSE,doppBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regularize, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
@ -10,11 +10,13 @@ subroutine UGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
! Input variables
logical :: doG0F2
logical :: doevGF2
logical :: doqsGF2
logical :: doG0F3
logical :: doevGF3
logical,intent(in) :: dotest
logical,intent(in) :: doG0F2
logical,intent(in) :: doevGF2
logical,intent(in) :: doqsGF2
logical,intent(in) :: doG0F3
logical,intent(in) :: doevGF3
integer :: renorm
integer,intent(in) :: maxSCF
@ -71,7 +73,7 @@ subroutine UGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doG0F2) then
call wall_time(start_GF)
call UG0F2(dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regularize, &
call UG0F2(dotest,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,epsHF)
call wall_time(end_GF)
@ -89,7 +91,7 @@ subroutine UGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doevGF2) then
call wall_time(start_GF)
call evUGF2(maxSCF,thresh,max_diis,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
call evUGF2(dotest,maxSCF,thresh,max_diis,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,cHF,epsHF)
call wall_time(end_GF)
@ -107,7 +109,7 @@ subroutine UGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm,maxSCF,thresh,max_di
if(doqsGF2) then
call wall_time(start_GF)
call qsUGF2(maxSCF,thresh,max_diis,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,regularize, &
call qsUGF2(dotest,maxSCF,thresh,max_diis,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,regularize, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
call wall_time(end_GF)

22
src/GF/evGGF2.f90
View File

@ -1,4 +1,4 @@
subroutine evGGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
subroutine evGGF2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
! Perform eigenvalue self-consistent second-order Green function calculation
@ -8,6 +8,8 @@ subroutine evGGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: dophBSE
logical,intent(in) :: doppBSE
logical,intent(in) :: TDA
@ -50,9 +52,9 @@ subroutine evGGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
! Hello world
write(*,*)
write(*,*)'************************************************'
write(*,*)'| Second-order Green function calculation |'
write(*,*)'************************************************'
write(*,*)'*********************************'
write(*,*)'* Generalized evGF2 Calculation *'
write(*,*)'*********************************'
write(*,*)
! Memory allocation
@ -105,7 +107,7 @@ subroutine evGGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
! Print results
call GMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eGF,Ec)
call GMP2(.false.,regularize,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eGF,Ec)
call print_evGF2(nBas,nO,nSCF,Conv,eHF,SigC,Z,eGF,ENuc,ERHF,Ec)
! DIIS extrapolation
@ -170,4 +172,14 @@ subroutine evGGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis, &
! end if
! Testing zone
if(dotest) then
call dump_test_value('G','evGGF2 correlation energy',Ec)
call dump_test_value('G','evGGF2 HOMO energy',eGF(nO))
call dump_test_value('G','evGGF2 LUMO energy',eGF(nO+1))
end if
end subroutine

23
src/GF/evGF2.f90 → src/GF/evRGF2.f90
View File

@ -1,4 +1,4 @@
subroutine evGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,triplet, &
subroutine evRGF2(dotest,dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,triplet, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
! Perform eigenvalue self-consistent second-order Green function calculation
@ -8,6 +8,8 @@ subroutine evGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,tr
! Input variables
logical,intent(in) :: dotest
logical,intent(in) :: dophBSE
logical,intent(in) :: doppBSE
logical,intent(in) :: TDA
@ -51,10 +53,11 @@ subroutine evGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,tr
! Hello world
write(*,*)
write(*,*)'************************************************'
write(*,*)'| Second-order Green function calculation |'
write(*,*)'************************************************'
write(*,*)'********************************'
write(*,*)'* Restricted evGF2 Calculation *'
write(*,*)'********************************'
write(*,*)
! Memory allocation
@ -107,7 +110,7 @@ subroutine evGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,tr
! Print results
call RMP2(regularize,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eGF,Ec)
call RMP2(.false.,regularize,nBas,nC,nO,nV,nR,ERI,ENuc,ERHF,eGF,Ec)
call print_evGF2(nBas,nO,nSCF,Conv,eHF,SigC,Z,eGF,ENuc,ERHF,Ec)
! DIIS extrapolation
@ -176,4 +179,14 @@ subroutine evGF2(dophBSE,doppBSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,tr
end if
! Testing zone
if(dotest) then
call dump_test_value('R','evRGF2 correlation energy',Ec)
call dump_test_value('R','evRGF2 HOMO energy',eGF(nO))
call dump_test_value('R','evRGF2 LUMO energy',eGF(nO+1))
end if
end subroutine

4
src/GF/evGF3.f90 → src/GF/evRGF3.f90
View File

@ -1,4 +1,4 @@
subroutine evGF3(maxSCF,thresh,max_diis,renormalization,nBas,nC,nO,nV,nR,V,e0)
subroutine evRGF3(dotest,maxSCF,thresh,max_diis,renormalization,nBas,nC,nO,nV,nR,V,e0)
! Perform third-order Green function calculation in diagonal approximation
@ -7,6 +7,8 @@
! Input variables
logical,intent(in) :: dotest
double precision,intent(in) :: thresh
integer,intent(in) :: maxSCF,max_diis,renormalization
integer,intent(in) :: nBas,nC,nO,nV,nR

13
src/GF/evUGF2.f90
View File

@ -1,4 +1,4 @@
subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
subroutine evUGF2(dotest,maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,cHF,eHF)
@ -9,6 +9,8 @@ subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_f
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -61,10 +63,11 @@ subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_f
! Hello world
write(*,*)
write(*,*)'**************************************************'
write(*,*)'| Self-consistent unrestricted evGF2 calculation |'
write(*,*)'**************************************************'
write(*,*)'*********************************'
write(*,*)'* Unrestricted G0F2 Calculation *'
write(*,*)'*********************************'
write(*,*)
! TDA
@ -133,7 +136,7 @@ subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_f
! Compute MP2 correlation energy
call UMP2(nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EUHF,eGF2,Ec)
call UMP2(.false.,nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EUHF,eGF2,Ec)
! Print results

23
src/GF/qsGF2.f90 → src/GF/qsRGF2.f90
View File

@ -1,4 +1,4 @@
subroutine qsGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet, &
subroutine qsRGF2(dotest,maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet, &
eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
@ -9,6 +9,8 @@ subroutine qsGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,tr
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -78,10 +80,11 @@ subroutine qsGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,tr
! Hello world
write(*,*)
write(*,*)'************************************************'
write(*,*)'| Self-consistent qsGF2 calculation |'
write(*,*)'************************************************'
write(*,*)'********************************'
write(*,*)'* Restricted qsGF2 Calculation *'
write(*,*)'********************************'
write(*,*)
! Warning
@ -216,7 +219,7 @@ subroutine qsGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,tr
! Correlation energy
call RMP2(regularize,nBas,nC,nO,nV,nR,ERI_MO,ENuc,EqsGF2,eGF,Ec)
call RMP2(.false.,regularize,nBas,nC,nO,nV,nR,ERI_MO,ENuc,EqsGF2,eGF,Ec)
! Total energy
@ -288,4 +291,14 @@ subroutine qsGF2(maxSCF,thresh,max_diis,dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,tr
end if
! Testing zone
if(dotest) then
call dump_test_value('R','qsRGF2 correlation energy',Ec)
call dump_test_value('R','qsRGF2 HOMO energy',eGF(nO))
call dump_test_value('R','qsRGF2 LUMO energy',eGF(nO+1))
end if
end subroutine

13
src/GF/qsUGF2.f90
View File

@ -1,4 +1,4 @@
subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,regularize, &
subroutine qsUGF2(dotest,maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,regularize, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EUHF,S,X,T,V,Hc,ERI_AO, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF)
@ -9,6 +9,8 @@ subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_f
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
@ -89,10 +91,11 @@ subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_f
! Hello world
write(*,*)
write(*,*)'**************************************************'
write(*,*)'| Self-consistent unrestricted qsGF2 calculation |'
write(*,*)'**************************************************'
write(*,*)'**********************************'
write(*,*)'* Unrestricted qsGF2 Calculation *'
write(*,*)'**********************************'
write(*,*)
! Warning
@ -290,7 +293,7 @@ subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_f
! Correlation energy
call UMP2(nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EqsGF2,eGF2,Ec)
call UMP2(.false.,nBas,nC,nO,nV,nR,ERI_aaaa,ERI_aabb,ERI_bbbb,ENuc,EqsGF2,eGF2,Ec)
! Total energy

4
src/QuAcK/GQuAcK.f90
View File

@ -233,7 +233,7 @@ subroutine GQuAcK(dotest,doGHF,dostab,dosearch,doMP2,doMP3,dophRPA,dophRPAx,dopp
if(doGF) then
call wall_time(start_GF)
call GGF(doG0F2,doevGF2,doqsGF2,maxSCF_GF,thresh_GF,max_diis_GF,dophBSE,doppBSE,TDA,dBSE,dTDA,lin_GF,eta_GF,reg_GF, &
call GGF(dotest,doG0F2,doevGF2,doqsGF2,maxSCF_GF,thresh_GF,max_diis_GF,dophBSE,doppBSE,TDA,dBSE,dTDA,lin_GF,eta_GF,reg_GF, &
nNuc,ZNuc,rNuc,ENuc,nBas2,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
call wall_time(end_GF)
@ -252,7 +252,7 @@ subroutine GQuAcK(dotest,doGHF,dostab,dosearch,doMP2,doMP3,dophRPA,dophRPAx,dopp
if(doGW) then
call wall_time(start_GW)
call GGW(doG0W0,doevGW,doqsGW,maxSCF_GW,thresh_GW,max_diis_GW,doACFDT,exchange_kernel,doXBS, &
call GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF_GW,thresh_GW,max_diis_GW,doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,lin_GW,eta_GW,reg_GW,nNuc,ZNuc,rNuc,ENuc, &
nBas,nBas2,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
call wall_time(end_GW)

10
src/QuAcK/RQuAcK.f90
View File

@ -225,7 +225,7 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
if(doCC) then
call wall_time(start_CC)
call RCC(doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD, &
call RCC(dotest,doCCD,dopCCD,doDCD,doCCSD,doCCSDT,dodrCCD,dorCCD,docrCCD,dolCCD, &
maxSCF_CC,thresh_CC,max_diis_CC,nBas,nC,nO,nV,nR,ERI_MO,ENuc,EHF,epsHF)
call wall_time(end_CC)
@ -244,7 +244,7 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
if(doCI) then
call wall_time(start_CI)
call RCI(doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO, &
call RCI(dotest,doCIS,doCIS_D,doCID,doCISD,doFCI,singlet,triplet,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO, &
epsHF,EHF,cHF,S)
call wall_time(end_CI)
@ -282,7 +282,7 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
if(doGF) then
call wall_time(start_GF)
call RGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm_GF,maxSCF_GF,thresh_GF,max_diis_GF, &
call RGF(dotest,doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm_GF,maxSCF_GF,thresh_GF,max_diis_GF, &
dophBSE,doppBSE,TDA,dBSE,dTDA,singlet,triplet,lin_GF,eta_GF,reg_GF, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO, &
dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
@ -303,7 +303,7 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
if(doGW) then
call wall_time(start_GW)
call RGW(doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF_GW,thresh_GW,max_diis_GW,doACFDT, &
call RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF_GW,thresh_GW,max_diis_GW,doACFDT, &
exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,singlet,triplet, &
lin_GW,eta_GW,reg_GW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc, &
ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
@ -324,7 +324,7 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
if(doGT) then
call wall_time(start_GT)
call RGT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,maxSCF_GT,thresh_GT,max_diis_GT,doACFDT, &
call RGT(dotest,doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,maxSCF_GT,thresh_GT,max_diis_GT,doACFDT, &
exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_T,TDA,dBSE,dTDA,singlet,triplet, &
lin_GT,eta_GT,reg_GT,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc, &
ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)

8
src/QuAcK/UQuAcK.f90
View File

@ -259,7 +259,7 @@ subroutine UQuAcK(dotest,doUHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,do
if(doCI) then
call wall_time(start_CI)
call UCI(doCIS,doCIS_D,doCID,doCISD,doFCI,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS, &
call UCI(dotest,doCIS,doCIS_D,doCID,doCISD,doFCI,spin_conserved,spin_flip,nBas,nC,nO,nV,nR,nS, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,epsHF,EHF,cHF,S)
call wall_time(end_CI)
@ -297,7 +297,7 @@ subroutine UQuAcK(dotest,doUHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,do
if(doGF) then
call wall_time(start_GF)
call UGF(doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm_GF,maxSCF_GF,thresh_GF,max_diis_GF, &
call UGF(dotest,doG0F2,doevGF2,doqsGF2,doG0F3,doevGF3,renorm_GF,maxSCF_GF,thresh_GF,max_diis_GF, &
dophBSE,doppBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,lin_GF,eta_GF,reg_GF, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
@ -318,7 +318,7 @@ subroutine UQuAcK(dotest,doUHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,do
if(doGW) then
call wall_time(start_GW)
call UGW(doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF_GW,thresh_GW,max_diis_GW,doACFDT, &
call UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF_GW,thresh_GW,max_diis_GW,doACFDT, &
exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
lin_GW,eta_GW,reg_GW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc, &
ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
@ -339,7 +339,7 @@ subroutine UQuAcK(dotest,doUHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,do
if(doGT) then
call wall_time(start_GT)
call UGT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,maxSCF_GT,thresh_GT,max_diis_GT,doACFDT, &
call UGT(dotest,doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,maxSCF_GT,thresh_GT,max_diis_GT,doACFDT, &
exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
lin_GT,eta_GT,reg_GT,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc, &
ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)

22
src/RPA/crRRPA.f90
View File

@ -105,10 +105,10 @@ subroutine crRRPA(dotest,TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(2X,A50,F20.10)') 'Tr@crRPA correlation energy (singlet) =',EcRPA(1)
write(*,'(2X,A50,F20.10)') 'Tr@crRPA correlation energy (triplet) =',EcRPA(2)
write(*,'(2X,A50,F20.10)') 'Tr@crRPA correlation energy =',EcRPA(1) + EcRPA(2)
write(*,'(2X,A50,F20.10)') 'Tr@crRPA total energy =',ENuc + EHF + EcRPA(1) + EcRPA(2)
write(*,'(2X,A50,F20.10,A3)') 'Tr@crRRPA correlation energy (singlet) =',EcRPA(1),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@crRRPA correlation energy (triplet) =',EcRPA(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@crRRPA correlation energy =',sum(EcRPA),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@crRRPA total energy =',ENuc + EHF + sum(EcRPA),' au'
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
@ -125,13 +125,19 @@ subroutine crRRPA(dotest,TDA,doACFDT,exchange_kernel,singlet,triplet,nBas,nC,nO,
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(2X,A50,F20.10)') 'AC@crRPA correlation energy (singlet) =',EcRPA(1)
write(*,'(2X,A50,F20.10)') 'AC@crRPA correlation energy (triplet) =',EcRPA(2)
write(*,'(2X,A50,F20.10)') 'AC@crRPA correlation energy =',EcRPA(1) + EcRPA(2)
write(*,'(2X,A50,F20.10)') 'AC@crRPA total energy =',ENuc + EHF + EcRPA(1) + EcRPA(2)
write(*,'(2X,A50,F20.10,A3)') 'AC@crRRPA correlation energy (singlet) =',EcRPA(1),' au'
write(*,'(2X,A50,F20.10,A3)') 'AC@crRRPA correlation energy (triplet) =',EcRPA(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'AC@crRRPA correlation energy =',sum(EcRPA),' au'
write(*,'(2X,A50,F20.10,A3)') 'AC@crRRPA total energy =',ENuc + EHF + sum(EcRPA),' au'
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
end if
if(dotest) then
call dump_test_value('R','crRRPA correlation energy',sum(EcRPA))
end if
end subroutine

12
test/Gtest_ref.dat
View File

@ -14,3 +14,15 @@
-0.368057033788489
ppGRPA correlation energy
-0.092561239023951
GG0F2 correlation energy
-0.139599095383592
GG0F2 HOMO energy
-0.398113441401304
GG0F2 LUMO energy
0.189499108558060
evGGF2 correlation energy
-0.143423802624327
evGGF2 HOMO energy
-0.374883361361466
evGGF2 LUMO energy
0.188029526841004

44
test/Rtest_ref.dat
View File

@ -6,13 +6,45 @@
0.203278954950938
RHF dipole moment
0.611349538338893
ROHF energy
-85.160473714509976
RMP2 correlation energy
-0.128988144386404
-0.128988144318865
CCD correlation energy
-0.134831644514626
DCD correlation energy
-0.137864343293459
CCSD correlation energy
-0.135521483429938
drCCD correlation energy
-0.138552794834488
rCCD correlation energy
-0.184028508433731
crCCD correlation energy
0.116371767570026
lCCD correlation energy
-0.092561235738883
pCCD correlation energy
-0.033110922532151
RCIS singlet excitation energy
0.345631869593838
RCIS triplet excitation energy
0.310356023571962
phRRPA correlation energy
-0.138552809856833
-0.138552809810619
phRRPAx correlation energy
-0.197284981952336
-0.197284981858331
crRRPA correlation energy
-0.161690198064065
ppRRPA correlation energy
-0.092561239071529
-0.092561239024299
RG0F2 correlation energy
-0.139599095383591
RG0F2 HOMO energy
-0.398113441401301
RG0F2 LUMO energy
0.189499108558074
evRGF2 correlation energy
-0.143423802624325
evRGF2 HOMO energy
-0.374883361361465
evRGF2 LUMO energy
0.188029526841018

4
test/Utest_ref.dat
View File

@ -12,6 +12,10 @@
0.611349538338891
UMP2 correlation energy
-0.128988144318865
UCIS singlet excitation energy
0.310356023571988
UCIS triplet excitation energy
0.310356023571958
phURPA correlation energy
-0.138552809810790
phURPAx correlation energy

14
test/methods.test
View File

@ -1,20 +1,20 @@
# RHF UHF GHF ROHF
T T T T
T F T F
# MP2 MP3
T T
# CCD pCCD DCD CCSD CCSD(T)
F F F F F
T T T T F
# drCCD rCCD crCCD lCCD
F F F F
T T T T
# CIS CIS(D) CID CISD FCI
F F F F F
T F F F F
# phRPA phRPAx crRPA ppRPA
T T F T
T T T T
# G0F2 evGF2 qsGF2 G0F3 evGF3
F F F F F
T T F F F
# G0W0 evGW qsGW SRG-qsGW ufG0W0 ufGW
F F F F F F
# G0T0pp evGTpp qsGTpp G0T0eh evGTeh qsGTeh
F F F F F F
# Rtest Utest Gtest
T T T
T F T

4
test/options.test
View File

@ -4,8 +4,8 @@
F
# CC: maxSCF thresh DIIS
64 0.0000001 5
# spin: TDA spin_conserved spin_flip
F T T
# spin: TDA singlet triplet
F T T
# GF: maxSCF thresh DIIS lin eta renorm reg
256 0.00001 5 F 0.0 0 F
# GW: maxSCF thresh DIIS lin eta TDA_W reg