4
1
mirror of https://github.com/pfloos/quack synced 2024-11-07 06:33:55 +01:00

B part of dBSE@GT

This commit is contained in:
Pierre-Francois Loos 2022-01-12 09:17:39 +01:00
parent fc0f9bb5f0
commit 96e4f15aeb
9 changed files with 143 additions and 33 deletions

View File

@ -13,9 +13,9 @@
# G0F2* evGF2* qsGF2* G0F3 evGF3
F F F F F
# G0W0* evGW* qsGW* ufG0W0 ufGW
T F F F F
F F F F F
# G0T0 evGT qsGT
F F F
T F F
# MCMP2
F
# * unrestricted version available

View File

@ -5,7 +5,7 @@
# CC: maxSCF thresh DIIS n_diis
64 0.00001 T 5
# spin: TDA singlet triplet spin_conserved spin_flip
F T T T T
F T F T T
# GF: maxSCF thresh DIIS n_diis lin eta renorm reg
256 0.00001 T 5 T 0.0 3 F
# GW: maxSCF thresh DIIS n_diis lin eta COHSEX SOSEX TDA_W G0W GW0 reg

View File

@ -85,8 +85,8 @@ subroutine Bethe_Salpeter_Tmatrix(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,
! call excitation_density_Tmatrix(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,ERI,X1s,Y1s,rho1s,X2s,Y2s,rho2s)
call static_Tmatrix_A(ispin,eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,1d0,ERI,Omega1s,rho1s,Omega2s,rho2s,TAs)
if(.not.TDA) call static_Tmatrix_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,1d0,ERI,Omega1s,rho1s,Omega2s,rho2s,TBs)
call static_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,1d0,Omega1s,rho1s,Omega2s,rho2s,TAs)
if(.not.TDA) call static_Tmatrix_B(eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,1d0,Omega1s,rho1s,Omega2s,rho2s,TBs)
! print*,'ab block of TA'
! call matout(nS,nS,TAs)
@ -105,8 +105,8 @@ subroutine Bethe_Salpeter_Tmatrix(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,
! call excitation_density_Tmatrix(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI,X1t,Y1t,rho1t,X2t,Y2t,rho2t)
call static_Tmatrix_A(ispin,eta,nBas,nC,nO,nV,nR,nS,nOOt,nVVt,1d0,ERI,Omega1t,rho1t,Omega2t,rho2t,TAt)
if(.not.TDA) call static_Tmatrix_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOOt,nVVt,1d0,ERI,Omega1t,rho1t,Omega2t,rho2t,TBt)
call static_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOOt,nVVt,1d0,Omega1t,rho1t,Omega2t,rho2t,TAt)
if(.not.TDA) call static_Tmatrix_B(eta,nBas,nC,nO,nV,nR,nS,nOOt,nVVt,1d0,Omega1t,rho1t,Omega2t,rho2t,TBt)
! print*,'aa block of TA'
! call matout(nS,nS,TAt)

View File

@ -62,12 +62,12 @@ subroutine dynamic_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGT,Omega1,O
chi = chi + rho1(i,b,cd)*rho1(j,a,cd)*eps/(eps**2 + eta**2)
end do
! do kl=1,nOO
! eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
! chi = chi + rho2(i,b,kl)*rho2(j,a,kl)*eps/(eps**2 + eta**2)
! end do
do kl=1,nOO
eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
chi = chi + rho2(i,b,kl)*rho2(j,a,kl)*eps/(eps**2 + eta**2)
end do
TA(ia,jb) = TA(ia,jb) + 1d0*lambda*chi
TA(ia,jb) = TA(ia,jb) + lambda*chi
chi = 0d0
@ -76,12 +76,12 @@ subroutine dynamic_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGT,Omega1,O
chi = chi + rho1(i,b,cd)*rho1(j,a,cd)*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
! do kl=1,nOO
! eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
! chi = chi + rho2(i,b,kl)*rho2(j,a,kl)*(eps**2 - eta**2)/(eps**2 + eta**2)**2
! end do
do kl=1,nOO
eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
chi = chi + rho2(i,b,kl)*rho2(j,a,kl)*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
ZA(ia,jb) = ZA(ia,jb) - 1d0*lambda*chi
ZA(ia,jb) = ZA(ia,jb) - lambda*chi
end do
end do

View File

@ -0,0 +1,91 @@
subroutine dynamic_Tmatrix_B(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGT,Omega1,Omega2,rho1,rho2,OmBSE,TB,ZB)
! Compute the off-diagonal dynamic part of the Bethe-Salpeter equation matrices for GT
implicit none
include 'parameters.h'
! Input variables
double precision,intent(in) :: eta
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) :: nOO
integer,intent(in) :: nVV
double precision,intent(in) :: lambda
double precision,intent(in) :: eGT(nBas)
double precision,intent(in) :: OmBSE
double precision,intent(in) :: Omega1(nVV)
double precision,intent(in) :: Omega2(nOO)
double precision,intent(in) :: rho1(nBas,nBas,nVV)
double precision,intent(in) :: rho2(nBas,nBas,nOO)
! Local variables
double precision :: chi
double precision :: eps
integer :: i,j,a,b,ia,jb,cd,kl
! Output variables
double precision,intent(out) :: TB(nS,nS)
double precision,intent(out) :: ZB(nS,nS)
! Initialization
TB(:,:) = 0d0
ZB(:,:) = 0d0
! Build dynamic A matrix
ia = 0
do i=nC+1,nO
do a=nO+1,nBas-nR
ia = ia + 1
jb = 0
do j=nC+1,nO
do b=nO+1,nBas-nR
jb = jb + 1
chi = 0d0
do cd=1,nVV
eps = + OmBSE - Omega1(cd) + (eGT(i) + eGT(j))
chi = chi + rho1(i,j,cd)*rho1(b,a,cd)*eps/(eps**2 + eta**2)
end do
do kl=1,nOO
eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
chi = chi + rho2(i,j,kl)*rho2(b,a,kl)*eps/(eps**2 + eta**2)
end do
TB(ia,jb) = TB(ia,jb) + lambda*chi
chi = 0d0
do cd=1,nVV
eps = + OmBSE - Omega1(cd) + (eGT(i) + eGT(j))
chi = chi + rho1(i,j,cd)*rho1(b,a,cd)*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
do kl=1,nOO
eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
chi = chi + rho2(i,a,kl)*rho2(b,a,kl)*(eps**2 - eta**2)/(eps**2 + eta**2)**2
end do
ZB(ia,jb) = ZB(ia,jb) - lambda*chi
end do
end do
end do
end do
end subroutine dynamic_Tmatrix_B

View File

@ -46,7 +46,6 @@ subroutine evGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS, &
! Local variables
logical :: linear_mixing
integer :: nSCF
integer :: n_diis
double precision :: rcond

View File

@ -70,7 +70,6 @@ subroutine qsGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,T
integer :: nOOs,nOOt
integer :: nVVs,nVVt
logical :: print_W = .false.
double precision,allocatable :: error_diis(:,:)
double precision,allocatable :: F_diis(:,:)
double precision,allocatable :: c(:,:)
@ -212,22 +211,47 @@ subroutine qsGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,T
call excitation_density_Tmatrix(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,ERI_MO, &
X1s,Y1s,rho1s,X2s,Y2s,rho2s)
call self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT, &
Omega1s,rho1s,Omega2s,rho2s,EcGM,SigT)
if(regularize) then
call renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT, &
Omega1s,rho1s,Omega2s,rho2s,Z)
call regularized_self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT, &
Omega1s,rho1s,Omega2s,rho2s,EcGM,SigT)
call regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT, &
Omega1s,rho1s,Omega2s,rho2s,Z)
else
call self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT, &
Omega1s,rho1s,Omega2s,rho2s,EcGM,SigT)
call renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT, &
Omega1s,rho1s,Omega2s,rho2s,Z)
end if
iblock = 4
call excitation_density_Tmatrix(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI_MO, &
X1t,Y1t,rho1t,X2t,Y2t,rho2t)
call self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
Omega1t,rho1t,Omega2t,rho2t,EcGM,SigT)
if(regularize) then
call renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
Omega1t,rho1t,Omega2t,rho2t,Z)
call self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
Omega1t,rho1t,Omega2t,rho2t,EcGM,SigT)
call renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
Omega1t,rho1t,Omega2t,rho2t,Z)
else
call regularized_self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
Omega1t,rho1t,Omega2t,rho2t,EcGM,SigT)
call regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
Omega1t,rho1t,Omega2t,rho2t,Z)
end if
Z(:) = 1d0/(1d0 - Z(:))

View File

@ -1,4 +1,4 @@
subroutine static_Tmatrix_A(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Omega1,rho1,Omega2,rho2,TA)
subroutine static_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,Omega1,rho1,Omega2,rho2,TA)
! Compute the OOVV block of the static T-matrix for the resonant block
@ -7,7 +7,6 @@ subroutine static_Tmatrix_A(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Ome
! Input variables
integer,intent(in) :: ispin
double precision,intent(in) :: eta
integer,intent(in) :: nBas
integer,intent(in) :: nC
@ -18,7 +17,6 @@ subroutine static_Tmatrix_A(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Ome
integer,intent(in) :: nOO
integer,intent(in) :: nVV
double precision,intent(in) :: lambda
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: Omega1(nVV)
double precision,intent(in) :: rho1(nBas,nBas,nVV)
double precision,intent(in) :: Omega2(nOO)

View File

@ -1,4 +1,4 @@
subroutine static_Tmatrix_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Omega1,rho1,Omega2,rho2,TB)
subroutine static_Tmatrix_B(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,Omega1,rho1,Omega2,rho2,TB)
! Compute the OVVO block of the static T-matrix for the coupling block
@ -7,7 +7,6 @@ subroutine static_Tmatrix_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Ome
! Input variables
integer,intent(in) :: ispin
double precision,intent(in) :: eta
integer,intent(in) :: nBas
integer,intent(in) :: nC
@ -18,7 +17,6 @@ subroutine static_Tmatrix_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Ome
integer,intent(in) :: nOO
integer,intent(in) :: nVV
double precision,intent(in) :: lambda
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: Omega1(nVV)
double precision,intent(in) :: rho1(nBas,nBas,nVV)
double precision,intent(in) :: Omega2(nOO)