4
1
mirror of https://github.com/pfloos/quack synced 2024-10-20 06:48:21 +02:00

renaming in GW module

This commit is contained in:
Pierre-Francois Loos 2024-09-03 16:06:09 +02:00
parent c05432416e
commit ede79d9ff4
41 changed files with 187 additions and 145 deletions

View File

@ -109,7 +109,7 @@ subroutine GG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
! Compute spectral weights !
!--------------------------!
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
call GGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
!------------------------!
! Compute GW self-energy !

View File

@ -1,4 +1,4 @@
subroutine GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
subroutine GGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
! Compute excitation densities

View File

@ -72,7 +72,7 @@ subroutine GGW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,eta,nBas,nC,nO,nV,nR,nS,ERI,di
if(.not.TDA_W) call phLR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call GGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
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)

View File

@ -114,7 +114,7 @@ subroutine RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
! Compute spectral weights !
!--------------------------!
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
!------------------------!
! Compute GW self-energy !
@ -122,7 +122,7 @@ subroutine RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
if(regularize) call GW_regularization(nBas,nC,nO,nV,nR,nS,eHF,Om,rho)
call GW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,EcGM,SigC,Z)
call RGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,EcGM,SigC,Z)
!-----------------------------------!
! Solve the quasi-particle equation !
@ -144,13 +144,13 @@ subroutine RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
write(*,*) ' *** Quasiparticle energies obtained by root search *** '
write(*,*)
call GW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eGWlin,eHF,eGW,Z)
call RGW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eGWlin,eHF,eGW,Z)
end if
! Plot self-energy, renormalization factor, and spectral function
! call GW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eHF,Om,rho)
! call RGW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eHF,Om,rho)
!--------------------!
! Cumulant expansion !
@ -176,7 +176,7 @@ subroutine RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
if(dophBSE) then
call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE)
call RGW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE)
if(exchange_kernel) then
@ -210,7 +210,7 @@ subroutine RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
end if
call GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eHF,eGW,EcBSE)
call RGW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eHF,eGW,EcBSE)
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
@ -227,7 +227,7 @@ subroutine RG0W0(dotest,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA
if(doppBSE) then
call GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE)
call RGW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE)
EcBSE(2) = 3d0*EcBSE(2)

View File

@ -1,10 +1,7 @@
! ---
subroutine RGW(dotest, doG0W0, doevGW, doqsGW, doufG0W0, doufGW, doSRGqsGW, maxSCF, thresh, max_diis, doACFDT, &
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, &
S, X, T, V, Hc, ERI_AO, ERI_MO, dipole_int_AO, dipole_int_MO, PHF, cHF, eHF)
subroutine RGW(dotest,doG0W0,doevGW,doqsGW,doufG0W0,doufGW,doSRGqsGW,maxSCF,thresh,max_diis,doACFDT, &
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, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
! Restricted GW module
@ -46,7 +43,8 @@ subroutine RGW(dotest, doG0W0, doevGW, doqsGW, doufG0W0, doufGW, doSRGqsGW, maxS
double precision,intent(in) :: rNuc(nNuc,ncart)
double precision,intent(in) :: ENuc
integer,intent(in) :: nBas, nOrb
integer,intent(in) :: nBas
integer,intent(in) :: nOrb
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nV
@ -112,10 +110,10 @@ subroutine RGW(dotest, doG0W0, doevGW, doqsGW, doufG0W0, doufGW, doSRGqsGW, maxS
if(doqsGW) then
call wall_time(start_GW)
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, &
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 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, &
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
@ -131,11 +129,11 @@ subroutine RGW(dotest, doG0W0, doevGW, doqsGW, doufG0W0, doufGW, doSRGqsGW, maxS
if(doSRGqsGW) then
call wall_time(start_GW)
call SRG_qsGW(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 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
@ -169,7 +167,7 @@ subroutine RGW(dotest, doG0W0, doevGW, doqsGW, doufG0W0, doufGW, doSRGqsGW, maxS
call wall_time(start_GW)
! TODO
call ufGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,eHF)
call ufRGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,eHF)
call wall_time(end_GW)
t_GW = end_GW - start_GW

View File

@ -1,4 +1,4 @@
double precision function GW_ImSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
double precision function RGW_ImSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Compute diagonal of the correlation part of the self-energy
@ -27,7 +27,7 @@ double precision function GW_ImSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Initialize
GW_ImSigC = 0d0
RGW_ImSigC = 0d0
! Occupied part of the correlation self-energy
@ -35,7 +35,7 @@ double precision function GW_ImSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(i) + Om(m)
num = 2d0*rho(p,i,m)**2
GW_ImSigC = GW_ImSigC + num*eta/(eps**2 + eta**2)
RGW_ImSigC = RGW_ImSigC + num*eta/(eps**2 + eta**2)
end do
end do
@ -45,7 +45,7 @@ double precision function GW_ImSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(a) - Om(m)
num = 2d0*rho(p,a,m)**2
GW_ImSigC = GW_ImSigC - num*eta/(eps**2 + eta**2)
RGW_ImSigC = RGW_ImSigC - num*eta/(eps**2 + eta**2)
end do
end do

View File

@ -1,4 +1,4 @@
double precision function GW_ImdSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
double precision function RGW_ImdSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Compute the derivative of the correlation part of the self-energy
@ -27,7 +27,7 @@ double precision function GW_ImdSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Initialize
GW_ImdSigC = 0d0
RGW_ImdSigC = 0d0
! Occupied part of the correlation self-energy
@ -35,7 +35,7 @@ double precision function GW_ImdSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(i) + Om(m)
num = 2d0*rho(p,i,m)**2
GW_ImdSigC = GW_ImdSigC - 2d0*num*eps*eta/(eps**2 + eta**2)**2
RGW_ImdSigC = RGW_ImdSigC - 2d0*num*eps*eta/(eps**2 + eta**2)**2
end do
end do
@ -45,7 +45,7 @@ double precision function GW_ImdSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(a) - Om(m)
num = 2d0*rho(p,a,m)**2
GW_ImdSigC = GW_ImdSigC + 2d0*num*eps*eta/(eps**2 + eta**2)**2
RGW_ImdSigC = RGW_ImdSigC + 2d0*num*eps*eta/(eps**2 + eta**2)**2
end do
end do

View File

@ -1,4 +1,4 @@
subroutine GW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eGWlin,eOld,eGW,Z)
subroutine RGW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eGWlin,eOld,eGW,Z)
! Compute the graphical solution of the QP equation
@ -28,7 +28,7 @@ subroutine GW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eGWlin,eOld,eGW,Z)
integer :: nIt
integer,parameter :: maxIt = 64
double precision,parameter :: thresh = 1d-6
double precision,external :: GW_ReSigC,GW_RedSigC
double precision,external :: RGW_ReSigC,RGW_RedSigC
double precision :: SigC,dSigC
double precision :: f,df
double precision :: w
@ -54,8 +54,8 @@ subroutine GW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eGWlin,eOld,eGW,Z)
nIt = nIt + 1
SigC = GW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,eOld,Om,rho)
dSigC = GW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,eOld,Om,rho)
SigC = RGW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,eOld,Om,rho)
dSigC = RGW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,eOld,Om,rho)
f = w - eHF(p) - SigC
df = 1d0/(1d0 - dSigC)
w = w - df*f

View File

@ -1,4 +1,4 @@
double precision function GW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
double precision function RGW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Compute diagonal of the correlation part of the self-energy
@ -27,7 +27,7 @@ double precision function GW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Initialize
GW_ReSigC = 0d0
RGW_ReSigC = 0d0
! Occupied part of the correlation self-energy
@ -35,7 +35,7 @@ double precision function GW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(i) + Om(m)
num = 2d0*rho(p,i,m)**2
GW_ReSigC = GW_ReSigC + num*eps/(eps**2 + eta**2)
RGW_ReSigC = RGW_ReSigC + num*eps/(eps**2 + eta**2)
end do
end do
@ -45,7 +45,7 @@ double precision function GW_ReSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(a) - Om(m)
num = 2d0*rho(p,a,m)**2
GW_ReSigC = GW_ReSigC + num*eps/(eps**2 + eta**2)
RGW_ReSigC = RGW_ReSigC + num*eps/(eps**2 + eta**2)
end do
end do

View File

@ -1,4 +1,4 @@
double precision function GW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
double precision function RGW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Compute the derivative of the correlation part of the self-energy
@ -27,7 +27,7 @@ double precision function GW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
! Initialize
GW_RedSigC = 0d0
RGW_RedSigC = 0d0
! Occupied part of the correlation self-energy
@ -35,7 +35,7 @@ double precision function GW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(i) + Om(m)
num = 2d0*rho(p,i,m)**2
GW_RedSigC = GW_RedSigC - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
RGW_RedSigC = RGW_RedSigC - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
end do
@ -45,7 +45,7 @@ double precision function GW_RedSigC(p,w,eta,nBas,nC,nO,nV,nR,nS,e,Om,rho)
do m=1,nS
eps = w - e(a) - Om(m)
num = 2d0*rho(p,a,m)**2
GW_RedSigC = GW_RedSigC - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
RGW_RedSigC = RGW_RedSigC - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
end do

View File

@ -0,0 +1,47 @@
subroutine RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
! Compute excitation densities
implicit none
! Input variables
integer,intent(in) :: nBas
integer,intent(in) :: nC
integer,intent(in) :: nO
integer,intent(in) :: nR
integer,intent(in) :: nS
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: XpY(nS,nS)
! Local variables
integer :: ia,jb,p,q,j,b
! Output variables
double precision,intent(out) :: rho(nBas,nBas,nS)
rho(:,:,:) = 0d0
!$OMP PARALLEL &
!$OMP SHARED(nC,nBas,nR,nO,nS,rho,ERI,XpY) &
!$OMP PRIVATE(q,p,jb,ia) &
!$OMP DEFAULT(NONE)
!$OMP DO
do q=nC+1,nBas-nR
do p=nC+1,nBas-nR
jb = 0
do j=nC+1,nO
do b=nO+1,nBas-nR
jb = jb + 1
do ia=1,nS
rho(p,q,ia) = rho(p,q,ia) + ERI(p,j,q,b)*XpY(ia,jb)
end do
end do
end do
end do
end do
!$OMP END DO
!$OMP END PARALLEL
end subroutine

View File

@ -1,4 +1,4 @@
subroutine GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,BSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eW,e,EcAC)
subroutine RGW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,BSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eW,e,EcAC)
! Compute the correlation energy via the adiabatic connection fluctuation dissipation theorem
@ -83,10 +83,10 @@ subroutine GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,BSE,singlet,triplet,e
if(.not.TDA_W) call phLR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call GW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA)
call GW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB)
call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA)
call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB)
! Singlet manifold
@ -113,10 +113,10 @@ subroutine GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,BSE,singlet,triplet,e
if(.not.TDA_W) call phLR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,lambda,ERI,Bph)
call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call GW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
call GW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
end if
@ -167,10 +167,10 @@ subroutine GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,BSE,singlet,triplet,e
if(.not.TDA_W) call phLR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,lambda,ERI,Bph)
call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call GW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
call GW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KA)
call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,OmRPA,rho_RPA,KB)
end if

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eW,eGW,EcBSE)
subroutine RGW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eW,eGW,EcBSE)
! Compute the Bethe-Salpeter excitation energies
@ -74,10 +74,10 @@ subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,
if(.not.TDA_W) call phLR_B(isp_W,dRPA_W,nBas,nC,nO,nV,nR,nS,1d0,ERI,Bph)
call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call GW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta)
call GW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta)
call RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KA_sta)
call RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,KB_sta)
!-------------------
! Singlet manifold
@ -103,10 +103,10 @@ subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,
write(*,*) '*** Second-order BSE static kernel activated! ***'
write(*,*)
call GW_phBSE_static_kernel(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
call GW_phBSE2_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KA_sta)
call RGW_phBSE_static_kernel(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
call RGW_phBSE2_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KA_sta)
if(.not.TDA) call GW_phBSE2_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KB_sta)
if(.not.TDA) call RGW_phBSE2_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KB_sta)
deallocate(W)
@ -131,7 +131,7 @@ subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,
!----------------------------------------------------!
if(dBSE) &
call GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
call RGW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
OmBSE,XpY_BSE,XmY_BSE,KA_sta,KB_sta)
end if
@ -163,7 +163,7 @@ subroutine GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,
!-------------------------------------------------
if(dBSE) &
call GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
call RGW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
OmBSE,XpY_BSE,XmY_BSE,KA_sta,KB_sta)
end if

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE,KA_dyn,ZA_dyn)
subroutine RGW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE,KA_dyn,ZA_dyn)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE2_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,eGW,W,KB_dyn)
subroutine RGW_phBSE2_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,eGW,W,KB_dyn)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE2_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eW,W,KA2_sta)
subroutine RGW_phBSE2_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eW,W,KA2_sta)
! Compute the second-order static BSE kernel for the resonant block (only for singlets!)

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE2_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,eW,W,KB2_sta)
subroutine RGW_phBSE2_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,eW,W,KB2_sta)
! Compute the second-order static BSE kernel for the coupling block (only for singlets!)

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rho_RPA,OmBSE,KA_dyn,ZA_dyn)
subroutine RGW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rho_RPA,OmBSE,KA_dyn,ZA_dyn)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rho,KB)
subroutine RGW_phBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rho,KB)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int, &
subroutine RGW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int, &
OmRPA,rho_RPA,OmBSE,XpY,XmY,KA_sta,KB_sta)
! Compute dynamical effects via perturbation theory for BSE
@ -69,7 +69,7 @@ subroutine GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,e
write(*,*)
allocate(W(nBas,nBas,nBas,nBas))
call GW_phBSE_static_kernel(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
call RGW_phBSE_static_kernel(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
end if
@ -90,9 +90,9 @@ subroutine GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,e
! Resonant part of the BSE correction for dynamical TDA
call GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,+OmBSE(ia),KAp_dyn,ZAp_dyn)
call RGW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,+OmBSE(ia),KAp_dyn,ZAp_dyn)
if(dophBSE2) call GW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE(ia),KAp_dyn,ZAp_dyn)
if(dophBSE2) call RGW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE(ia),KAp_dyn,ZAp_dyn)
if(dTDA) then
@ -103,10 +103,10 @@ subroutine GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,e
! Resonant and anti-resonant part of the BSE correction
call GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,-OmBSE(ia),KAm_dyn,ZAm_dyn)
call GW_phBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
call RGW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,-OmBSE(ia),KAm_dyn,ZAm_dyn)
call RGW_phBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
if(dophBSE2) call GW_phBSE2_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,eGW,W,KB_dyn)
if(dophBSE2) call RGW_phBSE2_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,eGW,W,KB_dyn)
! Renormalization factor of the resonant and anti-resonant parts

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE_static_kernel(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,W)
subroutine RGW_phBSE_static_kernel(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,W)
! Compute the second-order static BSE kernel for the resonant block (only for singlets!)

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,KA)
subroutine RGW_phBSE_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,KA)
! Compute the BSE static kernel for the resonant block

View File

@ -1,4 +1,4 @@
subroutine GW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,KB)
subroutine RGW_phBSE_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,ERI,Om,rho,KB)
! Compute the BSE static kernel for the coupling block

View File

@ -1,4 +1,4 @@
subroutine GW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eGW,Om,rho)
subroutine RGW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eGW,Om,rho)
! Dump several GW quantities for external plotting
@ -24,7 +24,7 @@ subroutine GW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eGW,Om,rho)
integer :: p,g
integer :: nGrid
double precision :: wmin,wmax,dw
double precision,external :: GW_ReSigC,GW_ImSigC,GW_RedSigC
double precision,external :: RGW_ReSigC,RGW_ImSigC,RGW_RedSigC
double precision,allocatable :: w(:)
double precision,allocatable :: ReSigC(:,:),ImSigC(:,:)
double precision,allocatable :: Z(:,:)
@ -56,9 +56,9 @@ subroutine GW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eGW,Om,rho)
do g=1,nGrid
do p=nC+1,nBas-nR
ReSigC(p,g) = GW_ReSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
ImSigC(p,g) = GW_ImSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
Z(p,g) = GW_RedSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
ReSigC(p,g) = RGW_ReSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
ImSigC(p,g) = RGW_ImSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
Z(p,g) = RGW_RedSigC(p,w(g),eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
end do
end do
@ -77,10 +77,10 @@ subroutine GW_plot_self_energy(nBas,eta,nC,nO,nV,nR,nS,eHF,eGW,Om,rho)
! Dump quantities in files as a function of w
open(unit=8 ,file='GW_SigC.dat')
open(unit=9 ,file='GW_freq.dat')
open(unit=10 ,file='GW_Z.dat')
open(unit=11 ,file='GW_A.dat')
open(unit=8 ,file='RGW_SigC.dat')
open(unit=9 ,file='RGW_freq.dat')
open(unit=10 ,file='RGW_Z.dat')
open(unit=11 ,file='RGW_A.dat')
do g=1,nGrid
write(8 ,*) w(g)*HaToeV,(ReSigC(p,g)*HaToeV,p=nC+1,nBas-nR)

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eW,eGW,EcBSE)
subroutine RGW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eW,eGW,EcBSE)
! Compute the Bethe-Salpeter excitation energies at the pp level
@ -83,7 +83,7 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
call phLR(TDA_W,nS,Aph,Bph,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY_RPA,rho_RPA)
deallocate(XpY_RPA,XmY_RPA,Aph,Bph)
@ -111,9 +111,9 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
! Compute BSE excitation energies
call GW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,ERI,OmRPA,rho_RPA,KC_sta)
call GW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,ERI,OmRPA,rho_RPA,KD_sta)
if(.not.TDA) call GW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,ERI,OmRPA,rho_RPA,KB_sta)
call RGW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,ERI,OmRPA,rho_RPA,KC_sta)
call RGW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,ERI,OmRPA,rho_RPA,KD_sta)
if(.not.TDA) call RGW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,ERI,OmRPA,rho_RPA,KB_sta)
call ppLR_C(ispin,nBas,nC,nO,nV,nR,nVV,1d0,eGW,ERI,Cpp)
call ppLR_D(ispin,nBas,nC,nO,nV,nR,nOO,1d0,eGW,ERI,Dpp)
@ -132,7 +132,7 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
!----------------------------------------------------!
if(dBSE) &
call GW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
call RGW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
Om1,X1,Y1,Om2,X2,Y2,KB_sta,KC_sta,KD_sta)
@ -163,9 +163,9 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
! Compute BSE excitation energies
call GW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,ERI,OmRPA,rho_RPA,KC_sta)
call GW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,ERI,OmRPA,rho_RPA,KD_sta)
if(.not.TDA) call GW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,ERI,OmRPA,rho_RPA,KB_sta)
call RGW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,ERI,OmRPA,rho_RPA,KC_sta)
call RGW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,ERI,OmRPA,rho_RPA,KD_sta)
if(.not.TDA) call RGW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,ERI,OmRPA,rho_RPA,KB_sta)
call ppLR_C(ispin,nBas,nC,nO,nV,nR,nVV,1d0,eGW,ERI,Cpp)
call ppLR_D(ispin,nBas,nC,nO,nV,nR,nOO,1d0,eGW,ERI,Dpp)
@ -184,7 +184,7 @@ subroutine GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,
!----------------------------------------------------!
if(dBSE) &
call GW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
call RGW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
Om1,X1,Y1,Om2,X2,Y2,KB_sta,KC_sta,KD_sta)
deallocate(Om1,X1,Y1,Om2,X2,Y2,Bpp,Cpp,Dpp,KB_sta,KC_sta,KD_sta)

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGW,Om,rho,KB_dyn)
subroutine RGW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGW,Om,rho,KB_dyn)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,eGW,Om,rho,OmBSE,KC_dyn,ZC_dyn)
subroutine RGW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,eGW,Om,rho,OmBSE,KC_dyn,ZC_dyn)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,eGW,Om,rho,OmBSE,KD_dyn,ZD_dyn)
subroutine RGW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,eGW,Om,rho,OmBSE,KD_dyn,ZD_dyn)
! Compute the dynamic part of the Bethe-Salpeter equation matrices

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,eW,eGW,ERI,dipole_int, &
subroutine RGW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,eW,eGW,ERI,dipole_int, &
OmRPA,rho_RPA,Om1,X1,Y1,Om2,X2,Y2,KB_sta,KC_sta,KD_sta)
! Compute dynamical effects via perturbation theory for BSE
@ -76,16 +76,16 @@ subroutine GW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,
if(dTDA) then
call GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,eGW,OmRPA,rho_RPA,Om1(ab),KC_dyn,ZC_dyn)
call RGW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,eGW,OmRPA,rho_RPA,Om1(ab),KC_dyn,ZC_dyn)
Z1_dyn(ab) = + dot_product(X1(:,ab),matmul(ZC_dyn,X1(:,ab)))
Om1_dyn(ab) = + dot_product(X1(:,ab),matmul(KC_dyn - KC_sta,X1(:,ab)))
else
call GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
call GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,eGW,OmRPA,rho_RPA,Om1(ab),KC_dyn,ZC_dyn)
call GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,eGW,OmRPA,rho_RPA,Om1(ab),KD_dyn,ZD_dyn)
call RGW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
call RGW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,eGW,OmRPA,rho_RPA,Om1(ab),KC_dyn,ZC_dyn)
call RGW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,eGW,OmRPA,rho_RPA,Om1(ab),KD_dyn,ZD_dyn)
Z1_dyn(ab) = dot_product(X1(:,ab),matmul(ZC_dyn,X1(:,ab))) &
+ dot_product(Y1(:,ab),matmul(ZD_dyn,Y1(:,ab)))
@ -119,16 +119,16 @@ subroutine GW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,
if(dTDA) then
call GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,eGW,OmRPA,rho_RPA,Om2(ij),KD_dyn,ZD_dyn)
call RGW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,eGW,OmRPA,rho_RPA,Om2(ij),KD_dyn,ZD_dyn)
Z2_dyn(kl) = dot_product(Y2(:,ij),matmul(ZD_dyn,Y2(:,ij)))
Om2_dyn(kl) = - dot_product(Y2(:,ij),matmul(KD_dyn - KD_sta,Y2(:,ij)))
else
call GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
call GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,eGW,OmRPA,rho_RPA,Om2(ij),KC_dyn,ZC_dyn)
call GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,eGW,OmRPA,rho_RPA,Om2(ij),KD_dyn,ZD_dyn)
call RGW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
call RGW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,1d0,eGW,OmRPA,rho_RPA,Om2(ij),KC_dyn,ZC_dyn)
call RGW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,1d0,eGW,OmRPA,rho_RPA,Om2(ij),KD_dyn,ZD_dyn)
Z2_dyn(kl) = dot_product(X2(:,ij),matmul(ZC_dyn,X2(:,ij))) &
+ dot_product(Y2(:,ij),matmul(ZD_dyn,Y2(:,ij)))

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Om,rho,KB)
subroutine RGW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Om,rho,KB)
! Compute the VVOO block of the static screening W for the pp-BSE

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,ERI,Om,rho,KC)
subroutine RGW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,ERI,Om,rho,KC)
! Compute the VVVV block of the static screening W for the pp-BSE

View File

@ -1,4 +1,4 @@
subroutine GW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,ERI,Om,rho,KD)
subroutine RGW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,ERI,Om,rho,KD)
! Compute the OOOO block of the static screening W for the pp-BSE

View File

@ -1,4 +1,4 @@
subroutine GW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
subroutine RGW_self_energy(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
! Compute correlation part of the self-energy and the renormalization factor

View File

@ -1,4 +1,4 @@
subroutine GW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
subroutine RGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,e,Om,rho,EcGM,Sig,Z)
! Compute diagonal of the correlation part of the self-energy and the renormalization factor

View File

@ -1,10 +1,7 @@
! ---
subroutine SRG_qsGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel, doXBS, &
BSE, BSE2, 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)
subroutine SRG_qsRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS, &
BSE,BSE2,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)
! Perform a quasiparticle self-consistent GW calculation
@ -230,7 +227,7 @@ subroutine SRG_qsGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel,
call wall_time(tex1)
call GW_excitation_density(nOrb,nC,nO,nR,nS,ERI_MO,XpY,rho)
call RGW_excitation_density(nOrb,nC,nO,nR,nS,ERI_MO,XpY,rho)
call wall_time(tex2)
tex=tex+tex2-tex1
@ -353,7 +350,7 @@ subroutine SRG_qsGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel,
if(BSE) then
call GW_phBSE(BSE2, TDA_W, TDA, dBSE, dTDA, singlet, triplet, eta, nOrb, &
call RGW_phBSE(BSE2, TDA_W, TDA, dBSE, dTDA, singlet, triplet, eta, nOrb, &
nC, nO, nV, nR, nS, ERI_MO, dipole_int_MO, eGW, eGW, EcBSE)
if(exchange_kernel) then
@ -388,7 +385,7 @@ subroutine SRG_qsGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel,
end if
call GW_phACFDT(exchange_kernel, doXBS, .true., TDA_W, TDA, BSE, singlet, triplet, &
call RGW_phACFDT(exchange_kernel, doXBS, .true., TDA_W, TDA, BSE, singlet, triplet, &
eta, nOrb, nC, nO, nV, nR, nS, ERI_MO, eGW, eGW, EcBSE)
write(*,*)

View File

@ -125,7 +125,7 @@ subroutine evGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
! Compute spectral weights
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
call GGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
! Compute correlation part of the self-energy

View File

@ -127,13 +127,13 @@ subroutine evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
! Compute spectral weights
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
! Compute correlation part of the self-energy
if(regularize) call GW_regularization(nBas,nC,nO,nV,nR,nS,eGW,Om,rho)
call GW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z)
call RGW_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eGW,Om,rho,EcGM,SigC,Z)
! Solve the quasi-particle equation
@ -149,7 +149,7 @@ subroutine evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
write(*,*) ' *** Quasiparticle energies obtained by root search *** '
write(*,*)
call GW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eOld,eOld,eGW,Z)
call RGW_QP_graph(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rho,eOld,eOld,eGW,Z)
end if
@ -221,7 +221,7 @@ subroutine evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
if(dophBSE) then
call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eGW,eGW,EcBSE)
call RGW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eGW,eGW,EcBSE)
if(exchange_kernel) then
@ -255,7 +255,7 @@ subroutine evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
end if
call GW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eGW,eGW,EcBSE)
call RGW_phACFDT(exchange_kernel,doXBS,dRPA,TDA_W,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,eGW,eGW,EcBSE)
write(*,*)
write(*,*)'-------------------------------------------------------------------------------'
@ -272,7 +272,7 @@ subroutine evRGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
if(doppBSE) then
call GW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE)
call RGW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGW,EcBSE)
EcBSE(2) = 3d0*EcBSE(2)

View File

@ -261,7 +261,7 @@ subroutine qsGGW(dotest,maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dop
! Compute correlation part of the self-energy
call GW_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)

View File

@ -206,11 +206,11 @@ subroutine qsRGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel, doX
! Compute correlation part of the self-energy
call GW_excitation_density(nOrb, nC, nO, nR, nS, ERI_MO, XpY, rho)
call RGW_excitation_density(nOrb, nC, nO, nR, nS, ERI_MO, XpY, rho)
if(regularize) call GW_regularization(nOrb, nC, nO, nV, nR, nS, eGW, Om, rho)
call GW_self_energy(eta, nOrb, nC, nO, nV, nR, nS, eGW, Om, rho, EcGM, SigC, Z)
call RGW_self_energy(eta, nOrb, nC, nO, nV, nR, nS, eGW, Om, rho, EcGM, SigC, Z)
! Make correlation self-energy Hermitian and transform it back to AO basis
@ -316,7 +316,7 @@ subroutine qsRGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel, doX
if(dophBSE) then
call GW_phBSE(dophBSE2, TDA_W, TDA, dBSE, dTDA, singlet, triplet, eta, &
call RGW_phBSE(dophBSE2, TDA_W, TDA, dBSE, dTDA, singlet, triplet, eta, &
nOrb, nC, nO, nV, nR, nS, ERI_MO, dipole_int_MO, eGW, eGW, EcBSE)
if(exchange_kernel) then
@ -351,7 +351,7 @@ subroutine qsRGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel, doX
end if
call GW_phACFDT(exchange_kernel, doXBS, .true., TDA_W, TDA, dophBSE, singlet, triplet, &
call RGW_phACFDT(exchange_kernel, doXBS, .true., TDA_W, TDA, dophBSE, singlet, triplet, &
eta, nOrb, nC, nO, nV, nR, nS, ERI_MO, eGW, eGW, EcBSE)
write(*,*)
@ -369,7 +369,7 @@ subroutine qsRGW(dotest, maxSCF, thresh, max_diis, doACFDT, exchange_kernel, doX
if(doppBSE) then
call GW_ppBSE(TDA_W, TDA, dBSE, dTDA, singlet, triplet, eta, nOrb, &
call RGW_ppBSE(TDA_W, TDA, dBSE, dTDA, singlet, triplet, eta, nOrb, &
nC, nO, nV, nR, nS, ERI_MO, dipole_int_MO, eHF, eGW, EcBSE)
EcBSE(2) = 3d0*EcBSE(2)

View File

@ -106,7 +106,7 @@ subroutine ufG0W0(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
! Compute spectral weights !
!--------------------------!
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
deallocate(Aph,Bph,XpY,XmY)

View File

@ -1,4 +1,4 @@
subroutine ufGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
subroutine ufRGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
! Unfold GW equations
@ -106,7 +106,7 @@ subroutine ufGW(dotest,TDA_W,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,eHF)
! Compute spectral weights !
!--------------------------!
call GW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
call RGW_excitation_density(nBas,nC,nO,nR,nS,ERI,XpY,rho)
deallocate(Aph,Bph,XpY,XmY)