PTEROSOR/QuAcK
10
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-09-27 12:00:50 +02:00
Code Releases Activity

remove SRG-qsGW methods

Browse Source
This commit is contained in:
Pierre-Francois Loos 2024-09-11 17:58:36 +02:00
parent c2d92e2299
commit 9f98f7b856
22 changed files with 486 additions and 631 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
input/methods.default
View File

@ -12,8 +12,8 @@
F F F F F F F F
# G0F2 evGF2 qsGF2 ufGF2 G0F3 evGF3 # G0F2 evGF2 qsGF2 ufGF2 G0F3 evGF3
F F F F F F F F F F F F
# G0W0 evGW qsGW SRG-qsGW ufG0W0 ufGW # G0W0 evGW qsGW ufG0W0 ufGW
F F F F F F F F F F F
# G0T0pp evGTpp qsGTpp ufG0T0pp # G0T0pp evGTpp qsGTpp ufG0T0pp
F F F F F F F F
# G0T0eh evGTeh qsGTeh # G0T0eh evGTeh qsGTeh

2
input/options.default
View File

@ -4,7 +4,7 @@
F F
# CC: maxSCF thresh DIIS # CC: maxSCF thresh DIIS
64 0.00001 5 64 0.00001 5
# spin: TDA singlet triplet # LR: TDA singlet triplet
F T T F T T
# GF: maxSCF thresh DIIS lin eta renorm reg # GF: maxSCF thresh DIIS lin eta renorm reg
256 0.00001 5 F 0.0 0 F 256 0.00001 5 F 0.0 0 F

18
src/GW/GG0W0.f90
View File

@ -1,5 +1,5 @@
subroutine GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & subroutine GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF) linearize,eta,doSRG,nBas,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF)
! Perform G0W0 calculation ! Perform G0W0 calculation
implicit none implicit none
@ -22,7 +22,7 @@ subroutine GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: linearize logical,intent(in) :: linearize
double precision,intent(in) :: eta double precision,intent(in) :: eta
logical,intent(in) :: regularize logical,intent(in) :: doSRG
integer,intent(in) :: nBas integer,intent(in) :: nBas
integer,intent(in) :: nC integer,intent(in) :: nC
@ -77,11 +77,13 @@ subroutine GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
write(*,*) write(*,*)
end if end if
! TDA ! SRG regularization
if(TDA) then if(doSRG) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*) write(*,*)
end if end if
! Memory allocation ! Memory allocation
@ -110,9 +112,11 @@ subroutine GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
! Compute GW self-energy ! ! Compute GW self-energy !
!------------------------! !------------------------!
if(regularize) call GW_regularization(nBas,nC,nO,nV,nR,nS,eHF,Om,rho) if(doSRG) then
call GGW_SRG_self_energy_diag(nBas,nC,nO,nV,nR,nS,eHF,Om,rho,EcGM,SigC,Z)
else
call GGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,EcGM,SigC,Z) call GGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,EcGM,SigC,Z)
end if
!-----------------------------------! !-----------------------------------!
! Solve the quasi-particle equation ! ! Solve the quasi-particle equation !

10
src/GW/GGW.f90
View File

@ -1,5 +1,5 @@
subroutine GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE, & subroutine GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE, &
TDA_W,TDA,dBSE,dTDA,linearize,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nBas2,nC,nO,nV,nR,nS,EGHF,S,X,T,V,Hc, & TDA_W,TDA,dBSE,dTDA,linearize,eta,doSRG,nNuc,ZNuc,rNuc,ENuc,nBas,nBas2,nC,nO,nV,nR,nS,EGHF,S,X,T,V,Hc, &
ERI_AO,ERI,dipole_int_AO,dipole_int,PHF,cHF,eHF) ERI_AO,ERI,dipole_int_AO,dipole_int,PHF,cHF,eHF)
! GW module ! GW module
@ -30,7 +30,7 @@ subroutine GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF,thresh,max_diis,doACFDT,exchan
logical,intent(in) :: doppBSE logical,intent(in) :: doppBSE
logical,intent(in) :: linearize logical,intent(in) :: linearize
double precision,intent(in) :: eta double precision,intent(in) :: eta
logical,intent(in) :: regularize logical,intent(in) :: doSRG
integer,intent(in) :: nNuc integer,intent(in) :: nNuc
double precision,intent(in) :: ZNuc(nNuc) double precision,intent(in) :: ZNuc(nNuc)
@ -71,7 +71,7 @@ subroutine GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF,thresh,max_diis,doACFDT,exchan
call wall_time(start_GW) call wall_time(start_GW)
call GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & call GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
linearize,eta,regularize,nBas2,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF) linearize,eta,doSRG,nBas2,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF)
call wall_time(end_GW) call wall_time(end_GW)
t_GW = end_GW - start_GW t_GW = end_GW - start_GW
@ -88,7 +88,7 @@ subroutine GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF,thresh,max_diis,doACFDT,exchan
call wall_time(start_GW) call wall_time(start_GW)
call evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & call evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
linearize,eta,regularize,nBas2,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF) linearize,eta,doSRG,nBas2,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF)
call wall_time(end_GW) call wall_time(end_GW)
t_GW = end_GW - start_GW t_GW = end_GW - start_GW
@ -105,7 +105,7 @@ subroutine GGW(dotest,doG0W0,doevGW,doqsGW,maxSCF,thresh,max_diis,doACFDT,exchan
call wall_time(start_GW) call wall_time(start_GW)
call qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & call qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nBas2,nC,nO,nV,nR,nS,EGHF,S,X,T,V,Hc,ERI_AO,ERI, & eta,doSRG,nNuc,ZNuc,rNuc,ENuc,nBas,nBas2,nC,nO,nV,nR,nS,EGHF,S,X,T,V,Hc,ERI_AO,ERI, &
dipole_int_AO,dipole_int,PHF,cHF,eHF) dipole_int_AO,dipole_int,PHF,cHF,eHF)
call wall_time(end_GW) call wall_time(end_GW)

144
src/GW/GGW_SRG_self_energy.f90 Normal file
View File

@ -0,0 +1,144 @@
subroutine GGW_SRG_self_energy(nOrb,nC,nO,nV,nR,nS,e,Om,rho,EcGM,SigC,Z)
! Compute correlation part of the self-energy
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: nOrb
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nV
integer,intent(in) :: nR
integer,intent(in) :: nS
double precision,intent(in) :: e(nOrb)
double precision,intent(in) :: Om(nS)
double precision,intent(in) :: rho(nOrb,nOrb,nS)
! Local variables
integer :: i,j,a,b
integer :: p,q
integer :: m
double precision :: Dpim,Dqim,Dpam,Dqam,Diam
double precision :: s
! Output variables
double precision,intent(out) :: EcGM
double precision,intent(out) :: SigC(nOrb,nOrb)
double precision,intent(out) :: Z(nOrb)
!--------------------!
! SRG flow parameter !
!--------------------!
s = 500d0
!--------------------!
! SRG-GW self-energy !
!--------------------!
SigC(:,:) = 0d0
! Occupied part of the correlation self-energy
!$OMP PARALLEL &
!$OMP SHARED(SigC,rho,s,nS,nC,nO,nOrb,nR,e,Om) &
!$OMP PRIVATE(m,i,q,p,Dpim,Dqim) &
!$OMP DEFAULT(NONE)
!$OMP DO
do q=nC+1,nOrb-nR
do p=nC+1,nOrb-nR
do m=1,nS
do i=nC+1,nO
Dpim = e(p) - e(i) + Om(m)
Dqim = e(q) - e(i) + Om(m)
SigC(p,q) = SigC(p,q) &
+ rho(p,i,m)*rho(q,i,m)* &
(1d0-exp(-s*Dpim*Dpim)*exp(-s*Dqim*Dqim))*(Dpim + Dqim)/(Dpim*Dpim + Dqim*Dqim)
end do
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
! Virtual part of the correlation self-energy
!$OMP PARALLEL &
!$OMP SHARED(SigC,rho,s,nS,nC,nO,nR,nOrb,e,Om) &
!$OMP PRIVATE(m,a,q,p,Dpam,Dqam) &
!$OMP DEFAULT(NONE)
!$OMP DO
do q=nC+1,nOrb-nR
do p=nC+1,nOrb-nR
do m=1,nS
do a=nO+1,nOrb-nR
Dpam = e(p) - e(a) - Om(m)
Dqam = e(q) - e(a) - Om(m)
SigC(p,q) = SigC(p,q) &
+ rho(p,a,m)*rho(q,a,m)* &
(1d0-exp(-s*Dpam*Dpam)*exp(-s*Dqam*Dqam))*(Dpam + Dqam)/(Dpam*Dpam + Dqam*Dqam)
end do
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
!------------------------!
! Renormalization factor !
!------------------------!
Z(:) = 0d0
! Occupied part of the renormlization factor
do p=nC+1,nOrb-nR
do i=nC+1,nO
do m=1,nS
Dpim = e(p) - e(i) + Om(m)
Z(p) = Z(p) - rho(p,i,m)**2*(1d0-exp(-2d0*s*Dpim*Dpim))/Dpim**2
end do
end do
end do
! Virtual part of the renormlization factor
do p=nC+1,nOrb-nR
do a=nO+1,nOrb-nR
do m=1,nS
Dpam = e(p) - e(a) - Om(m)
Z(p) = Z(p) - rho(p,a,m)**2*(1d0-exp(-2d0*s*Dpam*Dpam))/Dpam**2
end do
end do
end do
Z(:) = 1d0/(1d0 - Z(:))
!-------------------------------------!
! Galitskii-Migdal correlation energy !
!-------------------------------------!
EcGM = 0d0
do i=nC+1,nO
do a=nO+1,nOrb-nR
do m=1,nS
Diam = e(a) - e(i) + Om(m)
EcGM = EcGM - rho(a,i,m)*rho(a,i,m)*(1d0-exp(-2d0*s*Diam*Diam))/Diam
end do
end do
end do
end subroutine

134
src/GW/GGW_SRG_self_energy_diag.f90 Normal file
View File

@ -0,0 +1,134 @@
subroutine GGW_SRG_self_energy_diag(nOrb,nC,nO,nV,nR,nS,e,Om,rho,EcGM,SigC,Z)
! Compute correlation part of the self-energy
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: nOrb
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nV
integer,intent(in) :: nR
integer,intent(in) :: nS
double precision,intent(in) :: e(nOrb)
double precision,intent(in) :: Om(nS)
double precision,intent(in) :: rho(nOrb,nOrb,nS)
! Local variables
integer :: i,j,a,b
integer :: p
integer :: m
double precision :: Dpim,Dpam,Diam
double precision :: s
! Output variables
double precision,intent(out) :: EcGM
double precision,intent(out) :: SigC(nOrb)
double precision,intent(out) :: Z(nOrb)
!--------------------!
! SRG flow parameter !
!--------------------!
s = 500d0
!--------------------!
! SRG-GW self-energy !
!--------------------!
SigC(:) = 0d0
! Occupied part of the correlation self-energy
!$OMP PARALLEL &
!$OMP SHARED(SigC,rho,s,nS,nC,nO,nOrb,nR,e,Om) &
!$OMP PRIVATE(m,i,p,Dpim) &
!$OMP DEFAULT(NONE)
!$OMP DO
do p=nC+1,nOrb-nR
do m=1,nS
do i=nC+1,nO
Dpim = e(p) - e(i) + Om(m)
SigC(p) = SigC(p) + rho(p,i,m)**2*(1d0-exp(-2d0*s*Dpim*Dpim))/Dpim
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
! Virtual part of the correlation self-energy
!$OMP PARALLEL &
!$OMP SHARED(SigC,rho,s,nS,nC,nO,nR,nOrb,e,Om) &
!$OMP PRIVATE(m,a,p,Dpam) &
!$OMP DEFAULT(NONE)
!$OMP DO
do p=nC+1,nOrb-nR
do m=1,nS
do a=nO+1,nOrb-nR
Dpam = e(p) - e(a) - Om(m)
SigC(p) = SigC(p) + rho(p,a,m)**2*(1d0-exp(-2d0*s*Dpam*Dpam))/Dpam
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
!------------------------!
! Renormalization factor !
!------------------------!
Z(:) = 0d0
! Occupied part of the renormlization factor
do p=nC+1,nOrb-nR
do i=nC+1,nO
do m=1,nS
Dpim = e(p) - e(i) + Om(m)
Z(p) = Z(p) - rho(p,i,m)**2*(1d0-exp(-2d0*s*Dpim*Dpim))/Dpim**2
end do
end do
end do
! Virtual part of the renormlization factor
do p=nC+1,nOrb-nR
do a=nO+1,nOrb-nR
do m=1,nS
Dpam = e(p) - e(a) - Om(m)
Z(p) = Z(p) - rho(p,a,m)**2*(1d0-exp(-2d0*s*Dpam*Dpam))/Dpam**2
end do
end do
end do
Z(:) = 1d0/(1d0 - Z(:))
!-------------------------------------!
! Galitskii-Migdal correlation energy !
!-------------------------------------!
EcGM = 0d0
do i=nC+1,nO
do a=nO+1,nOrb-nR
do m=1,nS
Diam = e(a) - e(i) + Om(m)
EcGM = EcGM - rho(a,i,m)*rho(a,i,m)*(1d0-exp(-2d0*s*Diam*Diam))/Diam
end do
end do
end do
end subroutine

16
src/GW/GGW_phBSE.f90
View File

@ -73,9 +73,21 @@ subroutine GGW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,eta,nBas,nC,nO,nV,nR,nS,ERI,di
call GGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta) call GGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta)
call GGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta) call GGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta)
EcBSE = 0d0
! Compute BSE excitation energies !-----!
! TDA !
!-----!
if(TDA) then
write(*,*) 'Tamm-Dancoff approximation activated in phBSE!'
write(*,*)
end if
!---------------------------------!
! Compute BSE excitation energies !
!---------------------------------!
EcBSE = 0d0
call phGLR_A(dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI,Aph) call phGLR_A(dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI,Aph)
if(.not.TDA) call phGLR_B(dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph) if(.not.TDA) call phGLR_B(dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)

115
src/GW/GGW_self_energy.f90
View File

@ -1,4 +1,4 @@
subroutine GGW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z) subroutine GGW_self_energy(eta,nOrb,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
! Compute correlation part of the self-energy and the renormalization factor ! Compute correlation part of the self-energy and the renormalization factor
@ -8,15 +8,15 @@ subroutine GGW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
! Input variables ! Input variables
double precision,intent(in) :: eta double precision,intent(in) :: eta
integer,intent(in) :: nBas integer,intent(in) :: nOrb
integer,intent(in) :: nC integer,intent(in) :: nC
integer,intent(in) :: nO integer,intent(in) :: nO
integer,intent(in) :: nV integer,intent(in) :: nV
integer,intent(in) :: nR integer,intent(in) :: nR
integer,intent(in) :: nS integer,intent(in) :: nS
double precision,intent(in) :: e(nBas) double precision,intent(in) :: e(nOrb)
double precision,intent(in) :: Om(nS) double precision,intent(in) :: Om(nS)
double precision,intent(in) :: rho(nBas,nBas,nS) double precision,intent(in) :: rho(nOrb,nOrb,nS)
! Local variables ! Local variables
@ -27,34 +27,30 @@ subroutine GGW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
! Output variables ! Output variables
double precision,intent(out) :: EcGM double precision,intent(out) :: EcGM
double precision,intent(out) :: Sig(nBas,nBas) double precision,intent(out) :: Sig(nOrb,nOrb)
double precision,intent(out) :: Z(nBas) double precision,intent(out) :: Z(nOrb)
! Initialize
Sig(:,:) = 0d0
Z(:) = 0d0
!----------------! !----------------!
! GW self-energy ! ! GW self-energy !
!----------------! !----------------!
Sig(:,:) = 0d0
! Occupied part of the correlation self-energy ! Occupied part of the correlation self-energy
!$OMP PARALLEL & !$OMP PARALLEL &
!$OMP SHARED(Sig,Z,rho,eta,nS,nC,nO,nBas,nR,e,Om) & !$OMP SHARED(Sig,Z,rho,eta,nS,nC,nO,nOrb,nR,e,Om) &
!$OMP PRIVATE(m,i,q,p,eps,num) & !$OMP PRIVATE(m,i,q,p,eps,num) &
!$OMP DEFAULT(NONE) !$OMP DEFAULT(NONE)
!$OMP DO !$OMP DO
do q=nC+1,nBas-nR do q=nC+1,nOrb-nR
do p=nC+1,nBas-nR do p=nC+1,nOrb-nR
do m=1,nS do m=1,nS
do i=nC+1,nO do i=nC+1,nO
eps = e(p) - e(i) + Om(m) eps = e(p) - e(i) + Om(m)
num = rho(p,i,m)*rho(q,i,m) num = rho(p,i,m)*rho(q,i,m)
Sig(p,q) = Sig(p,q) + num*eps/(eps**2 + eta**2) Sig(p,q) = Sig(p,q) + num*eps/(eps**2 + eta**2)
if(p == q) Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do end do
end do end do
@ -65,33 +61,84 @@ subroutine GGW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
! Virtual part of the correlation self-energy ! Virtual part of the correlation self-energy
!$OMP PARALLEL & !$OMP PARALLEL &
!$OMP SHARED(Sig,Z,rho,eta,nS,nC,nO,nBas,nR,e,Om) & !$OMP SHARED(Sig,Z,rho,eta,nS,nC,nO,nOrb,nR,e,Om) &
!$OMP PRIVATE(m,a,q,p,eps,num) & !$OMP PRIVATE(m,a,q,p,eps,num) &
!$OMP DEFAULT(NONE) !$OMP DEFAULT(NONE)
!$OMP DO !$OMP DO
do q=nC+1,nBas-nR do q=nC+1,nOrb-nR
do p=nC+1,nBas-nR do p=nC+1,nOrb-nR
do m=1,nS do m=1,nS
do a=nO+1,nBas-nR do a=nO+1,nOrb-nR
eps = e(p) - e(a) - Om(m) eps = e(p) - e(a) - Om(m)
num = rho(p,a,m)*rho(q,a,m) num = rho(p,a,m)*rho(q,a,m)
Sig(p,q) = Sig(p,q) + num*eps/(eps**2 + eta**2) Sig(p,q) = Sig(p,q) + num*eps/(eps**2 + eta**2)
if(p == q) Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do end do
end do end do
end do end do
end do end do
!$OMP END DO !$OMP END DO
!$OMP END PARALLEL !$OMP END PARALLEL
! Galitskii-Migdal correlation energy !------------------------!
! Renormalization factor !
!------------------------!
Z(:) = 0d0
! Occupied part of the renormalization factor
!$OMP PARALLEL &
!$OMP SHARED(Sig,Z,rho,eta,nS,nC,nO,nOrb,nR,e,Om) &
!$OMP PRIVATE(m,i,q,p,eps,num) &
!$OMP DEFAULT(NONE)
!$OMP DO
do p=nC+1,nOrb-nR
do m=1,nS
do i=nC+1,nO
eps = e(p) - e(i) + Om(m)
num = rho(p,i,m)*rho(q,i,m)
Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
! Virtual part of the renormalization factor
!$OMP PARALLEL &
!$OMP SHARED(Sig,Z,rho,eta,nS,nC,nO,nOrb,nR,e,Om) &
!$OMP PRIVATE(m,a,q,p,eps,num) &
!$OMP DEFAULT(NONE)
!$OMP DO
do p=nC+1,nOrb-nR
do m=1,nS
do a=nO+1,nOrb-nR
eps = e(p) - e(a) - Om(m)
num = rho(p,a,m)*rho(q,a,m)
Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
Z(:) = 1d0/(1d0 - Z(:))
!-------------------------------------!
! Galitskii-Migdal correlation energy !
!-------------------------------------!
EcGM = 0d0 EcGM = 0d0
do m=1,nS do m=1,nS
do a=nO+1,nBas-nR do a=nO+1,nOrb-nR
do i=nC+1,nO do i=nC+1,nO
eps = e(a) - e(i) + Om(m) eps = e(a) - e(i) + Om(m)
@ -102,8 +149,4 @@ subroutine GGW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
end do end do
end do end do
! Compute renormalization factor from derivative
Z(:) = 1d0/(1d0 - Z(:))
end subroutine end subroutine

53
src/GW/GW_regularization.f90
View File

@ -1,53 +0,0 @@
subroutine GW_regularization(nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Regularize GW excitation densities via SRG
implicit none
! Input variables
integer,intent(in) :: nBas
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nV
integer,intent(in) :: nR
integer,intent(in) :: nS
integer,intent(in) :: e(nBas)
integer,intent(in) :: Om(nS)
! Local variables
integer :: p,i,a,m
double precision :: s
double precision :: kappa
double precision :: Dpim,Dpam
! Output variables
double precision,intent(inout):: rho(nBas,nBas,nS)
! SRG flow parameter
s = 500d0
! Regularize excitation densities
do p=nC+1,nBas-nR
do m=1,nS
do i=nC+1,nO
Dpim = e(p) - e(i) + Om(m)
kappa = 1d0 - exp(-Dpim*Dpim*s)
rho(p,i,m) = kappa*rho(p,i,m)
end do
do a=nO+1,nBas-nR
Dpam = e(p) - e(a) - Om(m)
kappa = 1d0 - exp(-Dpam*Dpam*s)
rho(p,a,m) = kappa*rho(p,a,m)
end do
end do
end do
end subroutine

33
src/GW/RGW.f90
View File

@ -1,6 +1,6 @@
subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thresh,max_diis,doACFDT, & subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,maxSCF,thresh,max_diis,doACFDT, &
exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,singlet,triplet, & exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,singlet,triplet, &
linearize,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nOrb,nC,nO,nV,nR,nS,ERHF, & linearize,eta,doSRG,nNuc,ZNuc,rNuc,ENuc,nBas,nOrb,nC,nO,nV,nR,nS,ERHF, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
! Restricted GW module ! Restricted GW module
@ -17,7 +17,6 @@ subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
logical,intent(in) :: doqsGW logical,intent(in) :: doqsGW
logical,intent(in) :: doufG0W0 logical,intent(in) :: doufG0W0
logical,intent(in) :: doufGW logical,intent(in) :: doufGW
logical,intent(in) :: doSRGqsGW
integer,intent(in) :: maxSCF integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis integer,intent(in) :: max_diis
@ -36,7 +35,7 @@ subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
logical,intent(in) :: triplet logical,intent(in) :: triplet
logical,intent(in) :: linearize logical,intent(in) :: linearize
double precision,intent(in) :: eta double precision,intent(in) :: eta
logical,intent(in) :: regularize logical,intent(in) :: doSRG
integer,intent(in) :: nNuc integer,intent(in) :: nNuc
double precision,intent(in) :: ZNuc(nNuc) double precision,intent(in) :: ZNuc(nNuc)
@ -77,7 +76,7 @@ subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
call wall_time(start_GW) call wall_time(start_GW)
call RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet, & call RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet, &
linearize,eta,regularize,nBas,nOrb,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int_MO,eHF) linearize,eta,doSRG,nBas,nOrb,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int_MO,eHF)
call wall_time(end_GW) call wall_time(end_GW)
t_GW = end_GW - start_GW t_GW = end_GW - start_GW
@ -94,7 +93,7 @@ subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
call wall_time(start_GW) call wall_time(start_GW)
call evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & call evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
singlet,triplet,linearize,eta,regularize,nBas,nOrb,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int_MO,eHF) singlet,triplet,linearize,eta,doSRG,nBas,nOrb,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int_MO,eHF)
call wall_time(end_GW) call wall_time(end_GW)
t_GW = end_GW - start_GW t_GW = end_GW - start_GW
@ -111,7 +110,7 @@ subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
call wall_time(start_GW) call wall_time(start_GW)
call qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2, & call qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2, &
TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc, & TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet,eta,doSRG,nNuc,ZNuc,rNuc, &
ENuc,nBas,nOrb,nC,nO,nV,nR,nS,ERHF,S,X,T,V,Hc,ERI_AO,ERI_MO, & ENuc,nBas,nOrb,nC,nO,nV,nR,nS,ERHF,S,X,T,V,Hc,ERI_AO,ERI_MO, &
dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
call wall_time(end_GW) call wall_time(end_GW)
@ -122,26 +121,6 @@ subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
end if end if
!------------------------------------------------------------------------
! Perform SRG-qsGW calculation
!------------------------------------------------------------------------
if(doSRGqsGW) then
call wall_time(start_GW)
call SRG_qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta, &
nNuc,ZNuc,rNuc,ENuc,nBas,nOrb,nC,nO,nV,nR,nS, &
ERHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO, &
PHF,cHF,eHF)
call wall_time(end_GW)
t_GW = end_GW - start_GW
write(*,'(A65,1X,F9.3,A8)') 'Total wall time for qsGW = ',t_GW,' seconds'
write(*,*)
end if
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Perform ufG0W0 calculatiom ! Perform ufG0W0 calculatiom
!------------------------------------------------------------------------ !------------------------------------------------------------------------

419
src/GW/SRG_qsUGW.f90
View File

@ -1,419 +0,0 @@
subroutine SRG_qsUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_W,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)
! Perform a quasiparticle self-consistent GW calculation
implicit none
include 'parameters.h'
! Input variables
logical,intent(in) :: dotest
integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis
double precision,intent(in) :: thresh
logical,intent(in) :: doACFDT
logical,intent(in) :: exchange_kernel
logical,intent(in) :: doXBS
logical,intent(in) :: BSE
logical,intent(in) :: TDA_W
logical,intent(in) :: TDA
logical,intent(in) :: dBSE
logical,intent(in) :: dTDA
logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip
double precision,intent(in) :: eta
logical,intent(in) :: regularize
integer,intent(in) :: nNuc
double precision,intent(in) :: ZNuc(nNuc)
double precision,intent(in) :: rNuc(nNuc,ncart)
double precision,intent(in) :: ENuc
integer,intent(in) :: nBas
integer,intent(in) :: nC(nspin)
integer,intent(in) :: nO(nspin)
integer,intent(in) :: nV(nspin)
integer,intent(in) :: nR(nspin)
integer,intent(in) :: nS(nspin)
double precision,intent(in) :: EUHF
double precision,intent(in) :: eHF(nBas,nspin)
double precision,intent(in) :: cHF(nBas,nBas,nspin)
double precision,intent(in) :: PHF(nBas,nBas,nspin)
double precision,intent(in) :: S(nBas,nBas)
double precision,intent(in) :: T(nBas,nBas)
double precision,intent(in) :: V(nBas,nBas)
double precision,intent(in) :: Hc(nBas,nBas)
double precision,intent(in) :: X(nBas,nBas)
double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas)
double precision,intent(inout):: ERI_aaaa(nBas,nBas,nBas,nBas)
double precision,intent(inout):: ERI_aabb(nBas,nBas,nBas,nBas)
double precision,intent(inout):: ERI_bbbb(nBas,nBas,nBas,nBas)
double precision,intent(in) :: dipole_int_AO(nBas,nBas,ncart)
double precision,intent(inout):: dipole_int_aa(nBas,nBas,ncart)
double precision,intent(inout):: dipole_int_bb(nBas,nBas,ncart)
! Local variables
logical :: dRPA
integer :: nSCF
integer :: nBasSq
integer :: ispin
integer :: ixyz
integer :: is
integer :: n_diis
integer :: nSa,nSb,nSt
double precision :: dipole(ncart)
double precision :: ET(nspin)
double precision :: EV(nspin)
double precision :: EJ(nsp)
double precision :: EK(nspin)
double precision :: EcRPA(nspin)
double precision :: EcGM(nspin)
double precision :: EqsGW
double precision :: EcBSE(nspin)
double precision :: Conv
double precision :: rcond(nspin)
double precision,external :: trace_matrix
double precision,allocatable :: err_diis(:,:,:)
double precision,allocatable :: F_diis(:,:,:)
double precision,allocatable :: Aph(:,:)
double precision,allocatable :: Bph(:,:)
double precision,allocatable :: Om(:)
double precision,allocatable :: XpY(:,:)
double precision,allocatable :: XmY(:,:)
double precision,allocatable :: rho(:,:,:,:)
double precision,allocatable :: c(:,:,:)
double precision,allocatable :: cp(:,:,:)
double precision,allocatable :: eGW(:,:)
double precision,allocatable :: P(:,:,:)
double precision,allocatable :: F(:,:,:)
double precision,allocatable :: Fp(:,:,:)
double precision,allocatable :: J(:,:,:)
double precision,allocatable :: K(:,:,:)
double precision,allocatable :: SigC(:,:,:)
double precision,allocatable :: SigCp(:,:,:)
double precision,allocatable :: Z(:,:)
double precision,allocatable :: err(:,:,:)
! Hello world
write(*,*)
write(*,*)'*************************************'
write(*,*)'* Unrestricted SRG-qsGW Calculation *'
write(*,*)'*************************************'
write(*,*)
! Warning
write(*,*) '!! ERIs in MO basis will be overwritten in qsUGW !!'
write(*,*)
! Stuff
nBasSq = nBas*nBas
dRPA = .true.
! TDA for W
if(TDA_W) then
write(*,*) 'Tamm-Dancoff approximation for dynamic screening!'
write(*,*)
end if
! TDA
if(TDA) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*)
end if
! Memory allocation
nSa = nS(1)
nSb = nS(2)
nSt = nSa + nSb
allocate(Aph(nSt,nSt),Bph(nSt,nSt),eGW(nBas,nspin),c(nBas,nBas,nspin),cp(nBas,nBas,nspin),P(nBas,nBas,nspin), &
F(nBas,nBas,nspin),Fp(nBas,nBas,nspin),J(nBas,nBas,nspin),K(nBas,nBas,nspin), &
SigC(nBas,nBas,nspin),SigCp(nBas,nBas,nspin),Z(nBas,nspin),Om(nSt),XpY(nSt,nSt),XmY(nSt,nSt), &
rho(nBas,nBas,nSt,nspin),err(nBas,nBas,nspin),err_diis(nBasSq,max_diis,nspin),F_diis(nBasSq,max_diis,nspin))
! Initialization
nSCF = -1
n_diis = 0
ispin = 1
Conv = 1d0
P(:,:,:) = PHF(:,:,:)
eGW(:,:) = eHF(:,:)
c(:,:,:) = cHF(:,:,:)
F_diis(:,:,:) = 0d0
err_diis(:,:,:) = 0d0
rcond(:) = 0d0
!------------------------------------------------------------------------
! Main loop
!------------------------------------------------------------------------
do while(Conv > thresh .and. nSCF < maxSCF)
! Increment
nSCF = nSCF + 1
! Buid Hartree matrix
do is=1,nspin
call Hartree_matrix_AO_basis(nBas,P(:,:,is),ERI_AO(:,:,:,:),J(:,:,is))
end do
! Compute exchange part of the self-energy
do is=1,nspin
call exchange_matrix_AO_basis(nBas,P(:,:,is),ERI_AO(:,:,:,:),K(:,:,is))
end do
!--------------------------------------------------
! AO to MO transformation of two-electron integrals
!--------------------------------------------------
do ixyz=1,ncart
call AOtoMO(nBas,nBas,c(:,:,1),dipole_int_AO(:,:,ixyz),dipole_int_aa(:,:,ixyz))
call AOtoMO(nBas,nBas,c(:,:,2),dipole_int_AO(:,:,ixyz),dipole_int_bb(:,:,ixyz))
end do
! 4-index transform for (aa|aa) block
call AOtoMO_ERI_UHF(1,1,nBas,c,ERI_AO,ERI_aaaa)
! 4-index transform for (aa|bb) block
call AOtoMO_ERI_UHF(1,2,nBas,c,ERI_AO,ERI_aabb)
! 4-index transform for (bb|bb) block
call AOtoMO_ERI_UHF(2,2,nBas,c,ERI_AO,ERI_bbbb)
! Compute linear response
call phULR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,1d0,eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,Aph)
if(.not.TDA) call phULR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nSa,nSb,nSt,1d0,ERI_aaaa,ERI_aabb,ERI_bbbb,Bph)
call phULR(TDA_W,nSa,nSb,nSt,Aph,Bph,EcRPA(ispin),Om,XpY,XmY)
!----------------------!
! Excitation densities !
!----------------------!
call UGW_excitation_density(nBas,nC,nO,nR,nSa,nSb,nSt,ERI_aaaa,ERI_aabb,ERI_bbbb,XpY,rho)
!------------------------------------------------!
! Compute self-energy and renormalization factor !
!------------------------------------------------!
call UGW_SRG_self_energy(nBas,nC,nO,nV,nR,nSt,eGW,Om,rho,EcGM,SigC,Z)
! Make correlation self-energy Hermitian and transform it back to AO basis
do is=1,nspin
call MOtoAO(nBas,nBas,S,c(:,:,is),SigC(:,:,is),SigCp(:,:,is))
end do
! Solve the quasi-particle equation
do is=1,nspin
F(:,:,is) = Hc(:,:) + J(:,:,is) + J(:,:,mod(is,2)+1) + K(:,:,is) + SigCp(:,:,is)
end do
! Check convergence
do is=1,nspin
err(:,:,is) = matmul(F(:,:,is),matmul(P(:,:,is),S(:,:))) - matmul(matmul(S(:,:),P(:,:,is)),F(:,:,is))
end do
if(nSCF > 1) Conv = maxval(abs(err))
! DIIS extrapolation
if(max_diis > 1) then
n_diis = min(n_diis+1,max_diis)
do is=1,nspin
if(nO(is) > 1) call DIIS_extrapolation(rcond(is),nBasSq,nBasSq,n_diis,err_diis(:,1:n_diis,is), &
F_diis(:,1:n_diis,is),err(:,:,is),F(:,:,is))
end do
end if
! Transform Fock matrix in orthogonal basis
do is=1,nspin
Fp(:,:,is) = matmul(transpose(X(:,:)),matmul(F(:,:,is),X(:,:)))
end do
! Diagonalize Fock matrix to get eigenvectors and eigenvalues
cp(:,:,:) = Fp(:,:,:)
do is=1,nspin
call diagonalize_matrix(nBas,cp(:,:,is),eGW(:,is))
end do
! Back-transform eigenvectors in non-orthogonal basis
do is=1,nspin
c(:,:,is) = matmul(X(:,:),cp(:,:,is))
end do
! Back-transform self-energy
do is=1,nspin
call AOtoMO(nBas,nBas,c(:,:,is),SigCp(:,:,is),SigC(:,:,is))
end do
! Compute density matrix
do is=1,nspin
P(:,:,is) = matmul(c(:,1:nO(is),is),transpose(c(:,1:nO(is),is)))
end do
!------------------------------------------------------------------------
! Compute total energy
!------------------------------------------------------------------------
! Kinetic energy
do is=1,nspin
ET(is) = trace_matrix(nBas,matmul(P(:,:,is),T(:,:)))
end do
! Potential energy
do is=1,nspin
EV(is) = trace_matrix(nBas,matmul(P(:,:,is),V(:,:)))
end do
! Hartree energy
EJ(1) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,1),J(:,:,1)))
EJ(2) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,1),J(:,:,2))) &
+ 0.5d0*trace_matrix(nBas,matmul(P(:,:,2),J(:,:,1)))
EJ(3) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,2),J(:,:,2)))
! Exchange energy
do is=1,nspin
EK(is) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,is),K(:,:,is)))
end do
! Total energy
EqsGW = sum(ET(:)) + sum(EV(:)) + sum(EJ(:)) + sum(EK(:))
!------------------------------------------------------------------------
! Print results
!------------------------------------------------------------------------
call dipole_moment(nBas,P(:,:,1)+P(:,:,2),nNuc,ZNuc,rNuc,dipole_int_AO,dipole)
call print_qsUGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,S,ENuc,ET,EV,EJ,EK,EcGM,EcRPA(ispin),EqsGW,SigCp,Z,dipole)
end do
!------------------------------------------------------------------------
! End main loop
!------------------------------------------------------------------------
! Did it actually converge?
if(nSCF == maxSCF) then
write(*,*)
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)' Convergence failed '
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)
stop
end if
! Deallocate memory
deallocate(cp,P,F,Fp,J,K,SigC,SigCp,Z,Om,XpY,XmY,rho,err,err_diis,F_diis)
! Perform BSE calculation
if(BSE) then
call UGW_phBSE(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, &
S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,c,eGW,eGW,EcBSE)
if(exchange_kernel) then
EcBSE(1) = 0.5d0*EcBSE(1)
EcBSE(2) = 0.5d0*EcBSE(2)
else
EcBSE(2) = 0.0d0
end if
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(2X,A50,F20.10,A3)') 'Tr@BSE@qsGW@UHF correlation energy (spin-conserved) = ',EcBSE(1),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@BSE@qsGW@UHF correlation energy (spin-flip) = ',EcBSE(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@BSE@qsGW@UHF correlation energy = ',sum(EcBSE),' au'
write(*,'(2X,A50,F20.10,A3)') 'Tr@BSE@qsGW@UHF total energy = ',ENuc + EqsGW + sum(EcBSE),' au'
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
! Compute the BSE correlation energy via the adiabatic connection
if(doACFDT) then
write(*,*) '--------------------------------------------------------------'
write(*,*) ' Adiabatic connection version of BSE@qsUGW correlation energy '
write(*,*) '--------------------------------------------------------------'
write(*,*)
if(doXBS) then
write(*,*) '*** scaled screening version (XBS) ***'
write(*,*)
end if
call UGW_phACFDT(exchange_kernel,doXBS,.true.,TDA_W,TDA,BSE,spin_conserved,spin_flip, &
eta,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb,eGW,eGW,EcRPA)
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
write(*,'(2X,A50,F20.10,A3)') 'AC@BSE@qsGW@UHF correlation energy (spin-conserved) = ',EcRPA(1),' au'
write(*,'(2X,A50,F20.10,A3)') 'AC@BSE@qsGW@UHF correlation energy (spin-flip) = ',EcRPA(2),' au'
write(*,'(2X,A50,F20.10,A3)') 'AC@BSE@qsGW@UHF correlation energy = ',sum(EcRPA),' au'
write(*,'(2X,A50,F20.10,A3)') 'AC@BSE@qsGW@UHF total energy = ',ENuc + EqsGW + sum(EcRPA),' au'
write(*,*)'-------------------------------------------------------------------------------'
write(*,*)
end if
end if
! Testing zone
if(dotest) then
call dump_test_value('U','qsGW correlation energy',EcRPA)
call dump_test_value('U','qsGW HOMOa energy',eGW(nO(1),1))
call dump_test_value('U','qsGW LUMOa energy',eGW(nO(1)+1,1))
call dump_test_value('U','qsGW HOMOa energy',eGW(nO(2),2))
call dump_test_value('U','qsGW LUMOa energy',eGW(nO(2)+1,2))
end if
end subroutine

8
src/GW/UG0W0.f90
View File

@ -86,11 +86,13 @@ subroutine UG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTD
write(*,*) write(*,*)
end if end if
! TDA ! SRG regularization
if(TDA) then if(doSRG) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*) write(*,*)
end if end if
! Memory allocation ! Memory allocation

31
src/GW/UGW.f90
View File

@ -1,6 +1,6 @@
subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thresh,max_diis,doACFDT, & subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,maxSCF,thresh,max_diis,doACFDT, &
exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
linearize,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EUHF,S,X,T,V,Hc, & linearize,eta,doSRG,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) ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF)
! GW module ! GW module
@ -17,7 +17,6 @@ subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
logical,intent(in) :: doqsGW logical,intent(in) :: doqsGW
logical,intent(in) :: doufG0W0 logical,intent(in) :: doufG0W0
logical,intent(in) :: doufGW logical,intent(in) :: doufGW
logical,intent(in) :: doSRGqsGW
integer,intent(in) :: maxSCF integer,intent(in) :: maxSCF
integer,intent(in) :: max_diis integer,intent(in) :: max_diis
@ -36,7 +35,7 @@ subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
logical,intent(in) :: linearize logical,intent(in) :: linearize
double precision,intent(in) :: eta double precision,intent(in) :: eta
logical,intent(in) :: regularize logical,intent(in) :: doSRG
integer,intent(in) :: nNuc integer,intent(in) :: nNuc
double precision,intent(in) :: ZNuc(nNuc) double precision,intent(in) :: ZNuc(nNuc)
@ -79,7 +78,7 @@ subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
call wall_time(start_GW) call wall_time(start_GW)
call UG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, & call UG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI_aaaa,ERI_aabb,ERI_bbbb, & linearize,eta,doSRG,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,cHF,eHF) dipole_int_aa,dipole_int_bb,cHF,eHF)
call wall_time(end_GW) call wall_time(end_GW)
@ -97,7 +96,7 @@ subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
call wall_time(start_GW) call wall_time(start_GW)
call evUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA, & call evUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA, &
spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S, & spin_conserved,spin_flip,linearize,eta,doSRG,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eHF) ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eHF)
call wall_time(end_GW) call wall_time(end_GW)
@ -115,25 +114,7 @@ subroutine UGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thre
call wall_time(start_GW) call wall_time(start_GW)
call qsUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA, & call qsUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,TDA,dBSE,dTDA, &
spin_conserved,spin_flip,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EUHF, & spin_conserved,spin_flip,eta,doSRG,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)
call wall_time(end_GW)
t_GW = end_GW - start_GW
write(*,'(A65,1X,F9.3,A8)') 'Total CPU time for qsGW = ',t_GW,' seconds'
write(*,*)
end if
!------------------------------------------------------------------------
! Perform SRG-qsGW calculation
!------------------------------------------------------------------------
if(doSRGqsGW) then
call wall_time(start_GW)
call SRG_qsUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_W,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) 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)
call wall_time(end_GW) call wall_time(end_GW)

11
src/GW/UGW_phBSE.f90
View File

@ -57,7 +57,7 @@ subroutine UGW_phBSE(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,nBas,nC,nO
double precision,intent(out) :: EcBSE(nspin) double precision,intent(out) :: EcBSE(nspin)
! Memory allocation ! Memory allocation
nS_aa = nS(1) nS_aa = nS(1)
nS_bb = nS(2) nS_bb = nS(2)
@ -69,6 +69,15 @@ subroutine UGW_phBSE(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,nBas,nC,nO
allocate(OmRPA(nS_sc),XpY_RPA(nS_sc,nS_sc),XmY_RPA(nS_sc,nS_sc),rho_RPA(nBas,nBas,nS_sc,nspin)) allocate(OmRPA(nS_sc),XpY_RPA(nS_sc,nS_sc),XmY_RPA(nS_sc,nS_sc),rho_RPA(nBas,nBas,nS_sc,nspin))
!-----!
! TDA !
!-----!
if(TDA) then
write(*,*) 'Tamm-Dancoff approximation activated in phBSE!'
write(*,*)
end if
!--------------------------! !--------------------------!
! Spin-conserved screening ! ! Spin-conserved screening !
!--------------------------! !--------------------------!

18
src/GW/evGGW.f90
View File

@ -1,5 +1,5 @@
subroutine evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & subroutine evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF) linearize,eta,doSRG,nBas,nC,nO,nV,nR,nS,ENuc,EGHF,ERI,dipole_int,eHF)
! Perform self-consistent eigenvalue-only GW calculation ! Perform self-consistent eigenvalue-only GW calculation
@ -27,7 +27,7 @@ subroutine evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
logical,intent(in) :: doppBSE logical,intent(in) :: doppBSE
logical,intent(in) :: linearize logical,intent(in) :: linearize
double precision,intent(in) :: eta double precision,intent(in) :: eta
logical,intent(in) :: regularize logical,intent(in) :: doSRG
integer,intent(in) :: nBas integer,intent(in) :: nBas
integer,intent(in) :: nC integer,intent(in) :: nC
@ -79,11 +79,13 @@ subroutine evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
write(*,*) write(*,*)
end if end if
! TDA ! SRG regularization
if(TDA) then if(doSRG) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*) write(*,*)
end if end if
! Linear mixing ! Linear mixing
@ -127,9 +129,11 @@ subroutine evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
! Compute correlation part of the self-energy ! Compute correlation part of the self-energy
if(regularize) call GW_regularization(nBas,nC,nO,nV,nR,nS,eGW,Om,rho) if(doSRG) then
call GGW_SRG_self_energy_diag(nBas,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z)
else
call GGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z) call GGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z)
end if
! Solve the quasi-particle equation ! Solve the quasi-particle equation

9
src/GW/evRGW.f90
View File

@ -83,6 +83,15 @@ subroutine evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
write(*,*) write(*,*)
end if end if
! SRG regularization
if(doSRG) then
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*)
end if
! Linear mixing ! Linear mixing
linear_mixing = .false. linear_mixing = .false.

8
src/GW/evUGW.f90
View File

@ -90,11 +90,13 @@ subroutine evUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE
write(*,*) write(*,*)
end if end if
! TDA ! SRG regularization
if(TDA) then if(doSRG) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*) write(*,*)
end if end if
! Initialization ! Initialization

18
src/GW/qsGGW.f90
View File

@ -1,5 +1,5 @@
subroutine qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, & subroutine qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nBas2,nC,nO,nV,nR,nS,EGHF,Ov,Or,T,V,Hc,ERI_AO, & eta,doSRG,nNuc,ZNuc,rNuc,ENuc,nBas,nBas2,nC,nO,nV,nR,nS,EGHF,Ov,Or,T,V,Hc,ERI_AO, &
ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
! Generalized version of quasiparticle self-consistent GW ! Generalized version of quasiparticle self-consistent GW
@ -25,7 +25,7 @@ subroutine qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: doppBSE logical,intent(in) :: doppBSE
double precision,intent(in) :: eta double precision,intent(in) :: eta
logical,intent(in) :: regularize logical,intent(in) :: doSRG
integer,intent(in) :: nNuc integer,intent(in) :: nNuc
double precision,intent(in) :: ZNuc(nNuc) double precision,intent(in) :: ZNuc(nNuc)
@ -130,11 +130,13 @@ subroutine qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
write(*,*) write(*,*)
end if end if
! TDA ! SRG regularization
if(TDA) then if(doSRG) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*) write(*,*)
end if end if
! Memory allocation ! Memory allocation
@ -261,9 +263,11 @@ subroutine qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
call GGW_excitation_density(nBas2,nC,nO,nR,nS,ERI_MO,XpY,rho) call GGW_excitation_density(nBas2,nC,nO,nR,nS,ERI_MO,XpY,rho)
if(regularize) call GW_regularization(nBas2,nC,nO,nV,nR,nS,eGW,Om,rho) if(doSRG) then
call GGW_SRG_self_energy(nBas2,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z)
else
call GGW_self_energy(eta,nBas2,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z) call GGW_self_energy(eta,nBas2,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z)
end if
! Make correlation self-energy Hermitian and transform it back to AO basis ! Make correlation self-energy Hermitian and transform it back to AO basis

8
src/GW/qsUGW.f90
View File

@ -126,11 +126,13 @@ subroutine qsUGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE
write(*,*) write(*,*)
end if end if
! TDA ! SRG regularization
if(TDA) then if(doSRG) then
write(*,*) 'Tamm-Dancoff approximation activated!'
write(*,*) '*** SRG regularized qsGW scheme ***'
write(*,*) write(*,*)
end if end if
! Memory allocation ! Memory allocation

8
src/QuAcK/QuAcK.f90
View File

@ -12,7 +12,7 @@ program QuAcK
logical :: doCIS,doCIS_D,doCID,doCISD,doFCI logical :: doCIS,doCIS_D,doCID,doCISD,doFCI
logical :: dophRPA,dophRPAx,docrRPA,doppRPA logical :: dophRPA,dophRPAx,docrRPA,doppRPA
logical :: doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3 logical :: doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3
logical :: doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW logical :: doG0W0,doevGW,doqsGW,doufG0W0,doufGW
logical :: doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh logical :: doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh
integer :: nNuc integer :: nNuc
@ -103,7 +103,7 @@ program QuAcK
dophRPA,dophRPAx,docrRPA,doppRPA, & dophRPA,dophRPAx,docrRPA,doppRPA, &
doG0F2,doevGF2,doqsGF2,doufG0F02, & doG0F2,doevGF2,doqsGF2,doufG0F02, &
doG0F3,doevGF3, & doG0F3,doevGF3, &
doG0W0,doevGW,doqsGW,doSRGqsGW, & doG0W0,doevGW,doqsGW, &
doufG0W0,doufGW, & doufG0W0,doufGW, &
doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp, & doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp, &
doG0T0eh,doevGTeh,doqsGTeh, & doG0T0eh,doevGTeh,doqsGTeh, &
@ -210,7 +210,7 @@ program QuAcK
if(doRQuAcK) & if(doRQuAcK) &
call RQuAcK(doRtest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,doCCSD,doCCSDT, & call RQuAcK(doRtest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,doCCSD,doCCSDT, &
dodrCCD,dorCCD,docrCCD,dolCCD,doCIS,doCIS_D,doCID,doCISD,doFCI,dophRPA,dophRPAx,docrRPA,doppRPA, & dodrCCD,dorCCD,docrCCD,dolCCD,doCIS,doCIS_D,doCID,doCISD,doFCI,dophRPA,dophRPAx,docrRPA,doppRPA, &
doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW, & doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3,doG0W0,doevGW,doqsGW,doufG0W0,doufGW, &
doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh, & doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh, &
nNuc,nBas,nOrb,nC,nO,nV,nR,ENuc,ZNuc,rNuc, & nNuc,nBas,nOrb,nC,nO,nV,nR,ENuc,ZNuc,rNuc, &
S,T,V,Hc,X,dipole_int_AO,ERI_AO,maxSCF_HF,max_diis_HF,thresh_HF,level_shift, & S,T,V,Hc,X,dipole_int_AO,ERI_AO,maxSCF_HF,max_diis_HF,thresh_HF,level_shift, &
@ -226,7 +226,7 @@ program QuAcK
if(doUQuAcK) & if(doUQuAcK) &
call UQuAcK(doUtest,doUHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,doCCSD,doCCSDT, & call UQuAcK(doUtest,doUHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,doCCSD,doCCSDT, &
dodrCCD,dorCCD,docrCCD,dolCCD,doCIS,doCIS_D,doCID,doCISD,doFCI,dophRPA,dophRPAx,docrRPA,doppRPA, & dodrCCD,dorCCD,docrCCD,dolCCD,doCIS,doCIS_D,doCID,doCISD,doFCI,dophRPA,dophRPAx,docrRPA,doppRPA, &
doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW, & doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3,doG0W0,doevGW,doqsGW,doufG0W0,doufGW, &
doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh, & doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh, &
nNuc,nBas,nC,nO,nV,nR,ENuc,ZNuc,rNuc, & nNuc,nBas,nC,nO,nV,nR,ENuc,ZNuc,rNuc, &
S,T,V,Hc,X,dipole_int_AO,ERI_AO,maxSCF_HF,max_diis_HF,thresh_HF,level_shift, & S,T,V,Hc,X,dipole_int_AO,ERI_AO,maxSCF_HF,max_diis_HF,thresh_HF,level_shift, &

8
src/QuAcK/RQuAcK.f90
View File

@ -1,6 +1,6 @@
subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,doCCSD,doCCSDT, & subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,doDCD,doCCSD,doCCSDT, &
dodrCCD,dorCCD,docrCCD,dolCCD,doCIS,doCIS_D,doCID,doCISD,doFCI,dophRPA,dophRPAx,docrRPA,doppRPA, & dodrCCD,dorCCD,docrCCD,dolCCD,doCIS,doCIS_D,doCID,doCISD,doFCI,dophRPA,dophRPAx,docrRPA,doppRPA, &
doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW, & doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3,doG0W0,doevGW,doqsGW,doufG0W0,doufGW, &
doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh, & doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp,doG0T0eh,doevGTeh,doqsGTeh, &
nNuc,nBas,nOrb,nC,nO,nV,nR,ENuc,ZNuc,rNuc, & nNuc,nBas,nOrb,nC,nO,nV,nR,ENuc,ZNuc,rNuc, &
S,T,V,Hc,X,dipole_int_AO,ERI_AO,maxSCF_HF,max_diis_HF,thresh_HF,level_shift, & S,T,V,Hc,X,dipole_int_AO,ERI_AO,maxSCF_HF,max_diis_HF,thresh_HF,level_shift, &
@ -25,7 +25,7 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
logical,intent(in) :: doCIS,doCIS_D,doCID,doCISD,doFCI logical,intent(in) :: doCIS,doCIS_D,doCID,doCISD,doFCI
logical,intent(in) :: dophRPA,dophRPAx,docrRPA,doppRPA logical,intent(in) :: dophRPA,dophRPAx,docrRPA,doppRPA
logical,intent(in) :: doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3 logical,intent(in) :: doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3
logical,intent(in) :: doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW logical,intent(in) :: doG0W0,doevGW,doqsGW,doufG0W0,doufGW
logical,intent(in) :: doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp logical,intent(in) :: doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp
logical,intent(in) :: doG0T0eh,doevGTeh,doqsGTeh logical,intent(in) :: doG0T0eh,doevGTeh,doqsGTeh
@ -309,12 +309,12 @@ subroutine RQuAcK(dotest,doRHF,doROHF,dostab,dosearch,doMP2,doMP3,doCCD,dopCCD,d
! GW module ! ! GW module !
!-----------! !-----------!
doGW = doG0W0 .or. doevGW .or. doqsGW .or. doufG0W0 .or. doufGW .or. doSRGqsGW doGW = doG0W0 .or. doevGW .or. doqsGW .or. doufG0W0 .or. doufGW
if(doGW) then if(doGW) then
call wall_time(start_GW) call wall_time(start_GW)
call RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF_GW,thresh_GW,max_diis_GW, & call RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,maxSCF_GW,thresh_GW,max_diis_GW, &
doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,doppBSE,TDA_W,TDA,dBSE,dTDA,singlet,triplet, & 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,nOrb,nC,nO,nV,nR,nS,ERHF,S,X,T, & lin_GW,eta_GW,reg_GW,nNuc,ZNuc,rNuc,ENuc,nBas,nOrb,nC,nO,nV,nR,nS,ERHF,S,X,T, &
V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)

12
src/QuAcK/read_methods.f90
View File

@ -6,7 +6,7 @@ subroutine read_methods(doRHF,doUHF,doGHF,doROHF, &
dophRPA,dophRPAx,docrRPA,doppRPA, & dophRPA,dophRPAx,docrRPA,doppRPA, &
doG0F2,doevGF2,doqsGF2,doufG0F02, & doG0F2,doevGF2,doqsGF2,doufG0F02, &
doG0F3,doevGF3, & doG0F3,doevGF3, &
doG0W0,doevGW,doqsGW,doSRGqsGW, & doG0W0,doevGW,doqsGW, &
doufG0W0,doufGW, & doufG0W0,doufGW, &
doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp, & doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp, &
doG0T0eh,doevGTeh,doqsGTeh, & doG0T0eh,doevGTeh,doqsGTeh, &
@ -25,7 +25,7 @@ subroutine read_methods(doRHF,doUHF,doGHF,doROHF, &
logical,intent(out) :: doCIS,doCIS_D,doCID,doCISD,doFCI logical,intent(out) :: doCIS,doCIS_D,doCID,doCISD,doFCI
logical,intent(out) :: dophRPA,dophRPAx,docrRPA,doppRPA logical,intent(out) :: dophRPA,dophRPAx,docrRPA,doppRPA
logical,intent(out) :: doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3 logical,intent(out) :: doG0F2,doevGF2,doqsGF2,doufG0F02,doG0F3,doevGF3
logical,intent(out) :: doG0W0,doevGW,doqsGW,doSRGqsGW,doufG0W0,doufGW logical,intent(out) :: doG0W0,doevGW,doqsGW,doufG0W0,doufGW
logical,intent(out) :: doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp logical,intent(out) :: doG0T0pp,doevGTpp,doqsGTpp,doufG0T0pp
logical,intent(out) :: doG0T0eh,doevGTeh,doqsGTeh logical,intent(out) :: doG0T0eh,doevGTeh,doqsGTeh
@ -147,18 +147,16 @@ subroutine read_methods(doRHF,doUHF,doGHF,doROHF, &
doG0W0 = .false. doG0W0 = .false.
doevGW = .false. doevGW = .false.
doqsGW = .false. doqsGW = .false.
doSRGqsGW = .false.
doufG0W0 = .false. doufG0W0 = .false.
doufGW = .false. doufGW = .false.
read(1,*) read(1,*)
read(1,*) ans1,ans2,ans3,ans4,ans5,ans6 read(1,*) ans1,ans2,ans3,ans4,ans5
if(ans1 == 'T') doG0W0 = .true. if(ans1 == 'T') doG0W0 = .true.
if(ans2 == 'T') doevGW = .true. if(ans2 == 'T') doevGW = .true.
if(ans3 == 'T') doqsGW = .true. if(ans3 == 'T') doqsGW = .true.
if(ans4 == 'T') doSRGqsGW = .true. if(ans4 == 'T') doufG0W0 = .true.
if(ans5 == 'T') doufG0W0 = .true. if(ans5 == 'T') doufGW = .true.
if(ans6 == 'T') doufGW = .true.
! Read GTpp methods ! Read GTpp methods