mirror of
https://github.com/pfloos/quack
synced 2024-12-22 20:34:46 +01:00
clean up GTpp
This commit is contained in:
parent
ee0cb0fd6a
commit
4fec6e88f5
@ -15,5 +15,5 @@
|
||||
# G0W0* evGW* qsGW* SRG-qsGW ufG0W0 ufGW
|
||||
F F F F F F
|
||||
# G0T0pp evGTpp qsGTpp G0T0eh evGTeh qsGTeh
|
||||
F F F T F F
|
||||
T T F F F F
|
||||
# * unrestricted version available
|
||||
|
@ -11,7 +11,7 @@
|
||||
# GW: maxSCF thresh DIIS n_diis lin eta TDA_W reg
|
||||
256 0.00001 T 5 T 0.0 F F
|
||||
# GT: maxSCF thresh DIIS n_diis lin eta TDA_T reg
|
||||
256 0.00001 T 5 T 0.0 F F
|
||||
256 0.00001 T 5 F 0.0 F F
|
||||
# ACFDT: AC Kx XBS
|
||||
F F T
|
||||
# BSE: phBSE phBSE2 ppBSE dBSE dTDA
|
||||
|
@ -41,22 +41,22 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
|
||||
integer :: ispin
|
||||
integer :: iblock
|
||||
integer :: nOOab,nOOaa
|
||||
integer :: nVVab,nVVaa
|
||||
integer :: nOOs,nOOt
|
||||
integer :: nVVs,nVVt
|
||||
double precision :: EcRPA(nspin)
|
||||
double precision :: EcBSE(nspin)
|
||||
double precision :: EcGM
|
||||
double precision,allocatable :: Bpp(:,:)
|
||||
double precision,allocatable :: Cpp(:,:)
|
||||
double precision,allocatable :: Dpp(:,:)
|
||||
double precision,allocatable :: Om1ab(:),Om1aa(:)
|
||||
double precision,allocatable :: X1ab(:,:),X1aa(:,:)
|
||||
double precision,allocatable :: Y1ab(:,:),Y1aa(:,:)
|
||||
double precision,allocatable :: rho1ab(:,:,:),rho1aa(:,:,:)
|
||||
double precision,allocatable :: Om2ab(:),Om2aa(:)
|
||||
double precision,allocatable :: X2ab(:,:),X2aa(:,:)
|
||||
double precision,allocatable :: Y2ab(:,:),Y2aa(:,:)
|
||||
double precision,allocatable :: rho2ab(:,:,:),rho2aa(:,:,:)
|
||||
double precision,allocatable :: Om1s(:),Om1t(:)
|
||||
double precision,allocatable :: X1s(:,:),X1t(:,:)
|
||||
double precision,allocatable :: Y1s(:,:),Y1t(:,:)
|
||||
double precision,allocatable :: rho1s(:,:,:),rho1t(:,:,:)
|
||||
double precision,allocatable :: Om2s(:),Om2t(:)
|
||||
double precision,allocatable :: X2s(:,:),X2t(:,:)
|
||||
double precision,allocatable :: Y2s(:,:),Y2t(:,:)
|
||||
double precision,allocatable :: rho2s(:,:,:),rho2t(:,:,:)
|
||||
double precision,allocatable :: Sig(:)
|
||||
double precision,allocatable :: Z(:)
|
||||
double precision,allocatable :: eGT(:)
|
||||
@ -74,20 +74,20 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
|
||||
! Dimensions of the pp-RPA linear reponse matrices
|
||||
|
||||
nOOab = nO*nO
|
||||
nVVab = nV*nV
|
||||
nOOs = nO*nO
|
||||
nVVs = nV*nV
|
||||
|
||||
nOOaa = nO*(nO - 1)/2
|
||||
nVVaa = nV*(nV - 1)/2
|
||||
nOOt = nO*(nO - 1)/2
|
||||
nVVt = nV*(nV - 1)/2
|
||||
|
||||
! Memory allocation
|
||||
|
||||
allocate(Om1ab(nVVab),X1ab(nVVab,nVVab),Y1ab(nOOab,nVVab), &
|
||||
Om2ab(nOOab),X2ab(nVVab,nOOab),Y2ab(nOOab,nOOab), &
|
||||
rho1ab(nBas,nBas,nVVab),rho2ab(nBas,nBas,nOOab), &
|
||||
Om1aa(nVVaa),X1aa(nVVaa,nVVaa),Y1aa(nOOaa,nVVaa), &
|
||||
Om2aa(nOOaa),X2aa(nVVaa,nOOaa),Y2aa(nOOaa,nOOaa), &
|
||||
rho1aa(nBas,nBas,nVVaa),rho2aa(nBas,nBas,nOOaa), &
|
||||
allocate(Om1s(nVVs),X1s(nVVs,nVVs),Y1s(nOOs,nVVs), &
|
||||
Om2s(nOOs),X2s(nVVs,nOOs),Y2s(nOOs,nOOs), &
|
||||
rho1s(nBas,nBas,nVVs),rho2s(nBas,nBas,nOOs), &
|
||||
Om1t(nVVt),X1t(nVVt,nVVt),Y1t(nOOt,nVVt), &
|
||||
Om2t(nOOt),X2t(nVVt,nOOt),Y2t(nOOt,nOOt), &
|
||||
rho1t(nBas,nBas,nVVt),rho2t(nBas,nBas,nOOt), &
|
||||
Sig(nBas),Z(nBas),eGT(nBas),eGTlin(nBas))
|
||||
|
||||
!----------------------------------------------
|
||||
@ -99,18 +99,18 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
|
||||
! Compute linear response
|
||||
|
||||
allocate(Bpp(nVVab,nOOab),Cpp(nVVab,nVVab),Dpp(nOOab,nOOab))
|
||||
allocate(Bpp(nVVs,nOOs),Cpp(nVVs,nVVs),Dpp(nOOs,nOOs))
|
||||
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVab,1d0,eHF,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOab,1d0,eHF,ERI_MO,Dpp)
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVs,1d0,eHF,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOs,1d0,eHF,ERI_MO,Dpp)
|
||||
|
||||
call ppLR(TDA_T,nOOab,nVVab,Bpp,Cpp,Dpp,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,EcRPA(ispin))
|
||||
call ppLR(TDA_T,nOOs,nVVs,Bpp,Cpp,Dpp,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcRPA(ispin))
|
||||
|
||||
deallocate(Bpp,Cpp,Dpp)
|
||||
|
||||
call print_excitation('pp-RPA (N+2)',iblock,nVVab,Om1ab(:))
|
||||
call print_excitation('pp-RPA (N-2)',iblock,nOOab,Om2ab(:))
|
||||
call print_excitation('pp-RPA (N+2)',iblock,nVVs,Om1s(:))
|
||||
call print_excitation('pp-RPA (N-2)',iblock,nOOs,Om2s(:))
|
||||
|
||||
!----------------------------------------------
|
||||
! alpha-alpha block
|
||||
@ -121,78 +121,50 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
|
||||
! Compute linear response
|
||||
|
||||
allocate(Bpp(nVVaa,nOOaa),Cpp(nVVaa,nVVaa),Dpp(nOOaa,nOOaa))
|
||||
allocate(Bpp(nVVt,nOOt),Cpp(nVVt,nVVt),Dpp(nOOt,nOOt))
|
||||
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVaa,1d0,eHF,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOaa,1d0,eHF,ERI_MO,Dpp)
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVt,1d0,eHF,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOt,1d0,eHF,ERI_MO,Dpp)
|
||||
|
||||
call ppLR(TDA_T,nOOaa,nVVaa,Bpp,Cpp,Dpp,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,EcRPA(ispin))
|
||||
call ppLR(TDA_T,nOOt,nVVt,Bpp,Cpp,Dpp,Om1t,X1t,Y1t,Om2t,X2t,Y2t,EcRPA(ispin))
|
||||
|
||||
deallocate(Bpp,Cpp,Dpp)
|
||||
|
||||
call print_excitation('pp-RPA (N+2)',iblock,nVVaa,Om1aa(:))
|
||||
call print_excitation('pp-RPA (N-2)',iblock,nOOaa,Om2aa(:))
|
||||
call print_excitation('ppRPA (N+2) ',iblock,nVVt,Om1t)
|
||||
call print_excitation('ppRPA (N-2) ',iblock,nOOt,Om2t)
|
||||
|
||||
!----------------------------------------------
|
||||
! Compute T-matrix version of the self-energy
|
||||
!----------------------------------------------
|
||||
|
||||
EcGM = 0d0
|
||||
Sig(:) = 0d0
|
||||
Z(:) = 0d0
|
||||
|
||||
iblock = 3
|
||||
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,ERI_MO,X1ab,Y1ab,rho1ab,X2ab,Y2ab,rho2ab)
|
||||
|
||||
if(regularize) then
|
||||
|
||||
call regularized_self_energy_Tmatrix_diag(eta,nBas,nC,nO,nV,nR,nOOab,nVVab,eHF,Om1ab,rho1ab,Om2ab,rho2ab,EcGM,Sig)
|
||||
call regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOab,nVVab,eHF,Om1ab,rho1ab,Om2ab,rho2ab,Z)
|
||||
|
||||
else
|
||||
|
||||
call GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOOab,nVVab,eHF,Om1ab,rho1ab,Om2ab,rho2ab,EcGM,Sig,Z)
|
||||
|
||||
end if
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,ERI_MO,X1s,Y1s,rho1s,X2s,Y2s,rho2s)
|
||||
|
||||
iblock = 4
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI_MO,X1t,Y1t,rho1t,X2t,Y2t,rho2t)
|
||||
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,ERI_MO,X1aa,Y1aa,rho1aa,X2aa,Y2aa,rho2aa)
|
||||
|
||||
if(regularize) then
|
||||
|
||||
call regularized_self_energy_Tmatrix_diag(eta,nBas,nC,nO,nV,nR,nOOaa,nVVaa,eHF,Om1aa,rho1aa,Om2aa,rho2aa,EcGM,Sig)
|
||||
call regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOaa,nVVaa,eHF,Om1aa,rho1aa,Om2aa,rho2aa,Z)
|
||||
|
||||
else
|
||||
|
||||
call GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOOaa,nVVaa,eHF,Om1aa,rho1aa,Om2aa,rho2aa,EcGM,Sig,Z)
|
||||
|
||||
end if
|
||||
|
||||
Z(:) = 1d0/(1d0 - Z(:))
|
||||
call GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eHF,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t,EcGM,Sig,Z)
|
||||
|
||||
!----------------------------------------------
|
||||
! Solve the quasi-particle equation
|
||||
!----------------------------------------------
|
||||
|
||||
eGTlin(:) = eHF(:) + Z(:)*Sig(:)
|
||||
|
||||
if(linearize) then
|
||||
|
||||
write(*,*) ' *** Quasiparticle energies obtained by linearization *** '
|
||||
write(*,*)
|
||||
|
||||
eGT(:) = eGTlin(:)
|
||||
eGT(:) = eHF(:) + Z(:)*Sig(:)
|
||||
|
||||
else
|
||||
|
||||
write(*,*) ' *** Quasiparticle energies obtained by root search (experimental) *** '
|
||||
write(*,*)
|
||||
|
||||
call GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOOaa,nVVaa,eHF,Om1aa,rho1aa,Om2aa,rho2aa,eGTlin,eGT)
|
||||
call GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eHF,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t,eHF,eGT)
|
||||
|
||||
end if
|
||||
|
||||
@ -205,26 +177,26 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
ispin = 1
|
||||
iblock = 3
|
||||
|
||||
allocate(Bpp(nVVab,nOOab),Cpp(nVVab,nVVab),Dpp(nOOab,nOOab))
|
||||
allocate(Bpp(nVVs,nOOs),Cpp(nVVs,nVVs),Dpp(nOOs,nOOs))
|
||||
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOab,nVVab,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVab,1d0,eGT,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOab,1d0,eGT,ERI_MO,Dpp)
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVs,1d0,eGT,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOs,1d0,eGT,ERI_MO,Dpp)
|
||||
|
||||
call ppLR(TDA_T,nOOab,nVVab,Bpp,Cpp,Dpp,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,EcRPA(ispin))
|
||||
call ppLR(TDA_T,nOOs,nVVs,Bpp,Cpp,Dpp,Om1s,X1s,Y1s,Om2s,X2s,Y2s,EcRPA(ispin))
|
||||
|
||||
deallocate(Bpp,Cpp,Dpp)
|
||||
|
||||
ispin = 2
|
||||
iblock = 4
|
||||
|
||||
allocate(Bpp(nVVaa,nOOaa),Cpp(nVVaa,nVVaa),Dpp(nOOaa,nOOaa))
|
||||
allocate(Bpp(nVVt,nOOt),Cpp(nVVt,nVVt),Dpp(nOOt,nOOt))
|
||||
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOaa,nVVaa,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVaa,1d0,eGT,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOaa,1d0,eGT,ERI_MO,Dpp)
|
||||
if(.not.TDA_T) call ppLR_B(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,1d0,ERI_MO,Bpp)
|
||||
call ppLR_C(iblock,nBas,nC,nO,nV,nR,nVVt,1d0,eGT,ERI_MO,Cpp)
|
||||
call ppLR_D(iblock,nBas,nC,nO,nV,nR,nOOt,1d0,eGT,ERI_MO,Dpp)
|
||||
|
||||
call ppLR(TDA_T,nOOaa,nVVaa,Bpp,Cpp,Dpp,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,EcRPA(ispin))
|
||||
call ppLR(TDA_T,nOOt,nVVt,Bpp,Cpp,Dpp,Om1t,X1t,Y1t,Om2t,X2t,Y2t,EcRPA(ispin))
|
||||
|
||||
deallocate(Bpp,Cpp,Dpp)
|
||||
|
||||
@ -237,8 +209,8 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
|
||||
if(dophBSE) then
|
||||
|
||||
call GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, &
|
||||
Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,rho1aa,rho2aa, &
|
||||
call GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,nOOt,nVVt, &
|
||||
Om1s,X1s,Y1s,Om2s,X2s,Y2s,rho1s,rho2s,Om1t,X1t,Y1t,Om2t,X2t,Y2t,rho1t,rho2t, &
|
||||
ERI_MO,dipole_int,eHF,eGT,EcBSE)
|
||||
|
||||
if(exchange_kernel) then
|
||||
@ -274,8 +246,8 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
end if
|
||||
|
||||
call GTpp_phACFDT(exchange_kernel,doXBS,.false.,TDA_T,TDA,dophBSE,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS, &
|
||||
nOOab,nVVab,nOOaa,nVVaa,Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa, &
|
||||
Om2aa,X2aa,Y2aa,rho1aa,rho2aa,ERI_MO,eHF,eGT,EcBSE)
|
||||
nOOs,nVVs,nOOt,nVVt,Om1s,X1s,Y1s,Om2s,X2s,Y2s,rho1s,rho2s,Om1t,X1t,Y1t, &
|
||||
Om2t,X2t,Y2t,rho1t,rho2t,ERI_MO,eHF,eGT,EcBSE)
|
||||
|
||||
if(exchange_kernel) then
|
||||
|
||||
@ -299,8 +271,8 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp
|
||||
|
||||
if(doppBSE) then
|
||||
|
||||
call GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nOOab,nVVab,nOOaa,nVVaa, &
|
||||
Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,rho1aa,rho2aa, &
|
||||
call GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt, &
|
||||
Om1s,X1s,Y1s,Om2s,X2s,Y2s,rho1s,rho2s,Om1t,X1t,Y1t,Om2t,X2t,Y2t,rho1t,rho2t, &
|
||||
ERI_MO,dipole_int,eHF,eGT,EcBSE)
|
||||
|
||||
write(*,*)
|
||||
|
@ -108,11 +108,11 @@ subroutine GT(doG0T0pp,doevGTpp,doqsGTpp,doG0T0eh,doevGTeh,doqsGTeh,unrestricted
|
||||
call wall_time(start_GT)
|
||||
if(unrestricted) then
|
||||
call evUGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
|
||||
eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, &
|
||||
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb, &
|
||||
PHF,cHF,epsHF)
|
||||
else
|
||||
call evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,singlet,triplet, &
|
||||
eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF)
|
||||
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI,dipole_int,PHF,cHF,epsHF)
|
||||
end if
|
||||
call wall_time(end_GT)
|
||||
|
||||
|
@ -1,4 +1,5 @@
|
||||
subroutine GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOO,nVV,eHF,Om1,rho1,Om2,rho2,eGTlin,eGT)
|
||||
subroutine GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eHF,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t,eGTlin,eGT)
|
||||
|
||||
implicit none
|
||||
include 'parameters.h'
|
||||
@ -9,15 +10,15 @@ subroutine GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOO,nVV,eHF,Om1,rho1,Om2,rho2,eGTl
|
||||
integer,intent(in) :: nO
|
||||
integer,intent(in) :: nV
|
||||
integer,intent(in) :: nR
|
||||
integer,intent(in) :: nOO
|
||||
integer,intent(in) :: nVV
|
||||
integer,intent(in) :: nOOs,nOOt
|
||||
integer,intent(in) :: nVVs,nVVt
|
||||
|
||||
double precision,intent(in) :: eta
|
||||
double precision,intent(in) :: eHF(nBas)
|
||||
double precision,intent(in) :: Om1(nVV)
|
||||
double precision,intent(in) :: rho1(nBas,nBas,nVV)
|
||||
double precision,intent(in) :: Om2(nOO)
|
||||
double precision,intent(in) :: rho2(nBas,nBas,nOO)
|
||||
double precision,intent(in) :: Om1s(nVVs),Om1t(nVVt)
|
||||
double precision,intent(in) :: rho1s(nBas,nBas,nVVs),rho1t(nBas,nBas,nVVt)
|
||||
double precision,intent(in) :: Om2s(nOOs),Om2t(nOOt)
|
||||
double precision,intent(in) :: rho2s(nBas,nBas,nOOs),rho2t(nBas,nBas,nOOt)
|
||||
|
||||
double precision,intent(in) :: eGTlin(nBas)
|
||||
|
||||
@ -54,8 +55,8 @@ subroutine GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOO,nVV,eHF,Om1,rho1,Om2,rho2,eGTl
|
||||
|
||||
nIt = nIt + 1
|
||||
|
||||
sigC = GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,eHF,Om1,rho1,Om2,rho2)
|
||||
dsigC = GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,eHF,Om1,rho1,Om2,rho2)
|
||||
sigC = GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eHF,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t)
|
||||
dsigC = GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eHF,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t)
|
||||
write (*,*) sigC
|
||||
f = w - eHF(p) - sigC
|
||||
df = 1d0 - dsigC
|
||||
|
@ -1,4 +1,5 @@
|
||||
double precision function GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2)
|
||||
double precision function GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,e,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t)
|
||||
|
||||
! Compute diagonal of the correlation part of the self-energy
|
||||
|
||||
@ -15,12 +16,14 @@ double precision function GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,
|
||||
integer,intent(in) :: nO
|
||||
integer,intent(in) :: nV
|
||||
integer,intent(in) :: nR
|
||||
integer,intent(in) :: nOO,nVV
|
||||
integer,intent(in) :: nOOs,nOOt
|
||||
integer,intent(in) :: nVVs,nVVt
|
||||
|
||||
double precision,intent(in) :: e(nBas)
|
||||
double precision,intent(in) :: Om1(nVV)
|
||||
double precision,intent(in) :: rho1(nBas,nBas,nVV)
|
||||
double precision,intent(in) :: Om2(nOO)
|
||||
double precision,intent(in) :: rho2(nBas,nBas,nOO)
|
||||
double precision,intent(in) :: Om1s(nVVs),Om1t(nVVt)
|
||||
double precision,intent(in) :: rho1s(nBas,nBas,nVVs),rho1t(nBas,nBas,nVVt)
|
||||
double precision,intent(in) :: Om2s(nOOs),Om2t(nOOt)
|
||||
double precision,intent(in) :: rho2s(nBas,nBas,nOOs),rho2t(nBas,nBas,nOOt)
|
||||
|
||||
! Local variables
|
||||
|
||||
@ -36,10 +39,17 @@ double precision function GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,
|
||||
!----------------------------------------------
|
||||
|
||||
do i=nC+1,nO
|
||||
do cd=1,nVV
|
||||
eps = w + e(i) - Om1(cd)
|
||||
GTpp_SigC = GTpp_SigC + rho1(p,i,cd)**2*eps/(eps**2 + eta**2)
|
||||
|
||||
do cd=1,nVVs
|
||||
eps = w + e(i) - Om1s(cd)
|
||||
GTpp_SigC = GTpp_SigC + rho1s(p,i,cd)**2*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
do cd=1,nVVt
|
||||
eps = w + e(i) - Om1t(cd)
|
||||
GTpp_SigC = GTpp_SigC + rho1t(p,i,cd)**2*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
enddo
|
||||
|
||||
!----------------------------------------------
|
||||
@ -47,10 +57,17 @@ double precision function GTpp_SigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,
|
||||
!----------------------------------------------
|
||||
|
||||
do a=nO+1,nBas-nR
|
||||
do kl=1,nOO
|
||||
eps = w + e(a) - Om2(kl)
|
||||
GTpp_SigC = GTpp_SigC + rho2(p,a,kl)**2*eps/(eps**2 + eta**2)
|
||||
|
||||
do kl=1,nOOs
|
||||
eps = w + e(a) - Om2s(kl)
|
||||
GTpp_SigC = GTpp_SigC + rho2s(p,a,kl)**2*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
do kl=1,nOOt
|
||||
eps = w + e(a) - Om2t(kl)
|
||||
GTpp_SigC = GTpp_SigC + rho2t(p,a,kl)**2*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
enddo
|
||||
|
||||
end function
|
||||
|
@ -1,4 +1,5 @@
|
||||
double precision function GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2)
|
||||
double precision function GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,e,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t)
|
||||
|
||||
! Compute diagonal of the correlation part of the self-energy
|
||||
|
||||
@ -15,13 +16,14 @@ double precision function GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1
|
||||
integer,intent(in) :: nO
|
||||
integer,intent(in) :: nV
|
||||
integer,intent(in) :: nR
|
||||
integer,intent(in) :: nOO
|
||||
integer,intent(in) :: nVV
|
||||
integer,intent(in) :: nOOs,nOOt
|
||||
integer,intent(in) :: nVVs,nVVt
|
||||
|
||||
double precision,intent(in) :: e(nBas)
|
||||
double precision,intent(in) :: Om1(nVV)
|
||||
double precision,intent(in) :: rho1(nBas,nBas,nVV)
|
||||
double precision,intent(in) :: Om2(nOO)
|
||||
double precision,intent(in) :: rho2(nBas,nBas,nOO)
|
||||
double precision,intent(in) :: Om1s(nVVs),Om1t(nVVt)
|
||||
double precision,intent(in) :: rho1s(nBas,nBas,nVVs),rho1t(nBas,nBas,nVVt)
|
||||
double precision,intent(in) :: Om2s(nOOs),Om2t(nOOt)
|
||||
double precision,intent(in) :: rho2s(nBas,nBas,nOOs),rho2t(nBas,nBas,nOOt)
|
||||
|
||||
! Local variables
|
||||
|
||||
@ -37,12 +39,17 @@ double precision function GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1
|
||||
!----------------------------------------------
|
||||
|
||||
do i=nC+1,nO
|
||||
do cd=1,nVV
|
||||
|
||||
eps = w + e(i) - Om1(cd)
|
||||
GTpp_dSigC= GTpp_dSigC - rho1(p,i,cd)**2*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
|
||||
do cd=1,nVVs
|
||||
eps = w + e(i) - Om1s(cd)
|
||||
GTpp_dSigC = GTpp_dSigC - rho1s(p,i,cd)**2*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
|
||||
do cd=1,nVVt
|
||||
eps = w + e(i) - Om1t(cd)
|
||||
GTpp_dSigC = GTpp_dSigC - rho1t(p,i,cd)**2*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
|
||||
enddo
|
||||
|
||||
!----------------------------------------------
|
||||
@ -51,10 +58,17 @@ double precision function GTpp_dSigC(p,w,eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1
|
||||
!----------------------------------------------
|
||||
|
||||
do a=nO+1,nBas-nR
|
||||
do kl=1,nOO
|
||||
eps = w + e(a) - Om2(kl)
|
||||
GTpp_dSigC = GTpp_dSigC - rho2(p,a,kl)**2*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
|
||||
do kl=1,nOOs
|
||||
eps = w + e(a) - Om2s(kl)
|
||||
GTpp_dSigC = GTpp_dSigC - rho2s(p,a,kl)**2*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
|
||||
do kl=1,nOOt
|
||||
eps = w + e(a) - Om2t(kl)
|
||||
GTpp_dSigC = GTpp_dSigC - rho2t(p,a,kl)**2*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
|
||||
enddo
|
||||
|
||||
end function
|
||||
|
@ -1,4 +1,4 @@
|
||||
subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcGM,Sig,Z)
|
||||
subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,e,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t,EcGM,Sig,Z)
|
||||
|
||||
! Compute the correlation part of the T-matrix self-energy and the renormalization factor
|
||||
|
||||
@ -13,13 +13,13 @@ subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcG
|
||||
integer,intent(in) :: nO
|
||||
integer,intent(in) :: nV
|
||||
integer,intent(in) :: nR
|
||||
integer,intent(in) :: nOO
|
||||
integer,intent(in) :: nVV
|
||||
integer,intent(in) :: nOOs,nOOt
|
||||
integer,intent(in) :: nVVs,nVVt
|
||||
double precision,intent(in) :: e(nBas)
|
||||
double precision,intent(in) :: Om1(nVV)
|
||||
double precision,intent(in) :: rho1(nBas,nBas,nVV)
|
||||
double precision,intent(in) :: Om2(nOO)
|
||||
double precision,intent(in) :: rho2(nBas,nBas,nOO)
|
||||
double precision,intent(in) :: Om1s(nVVs),Om1t(nVVt)
|
||||
double precision,intent(in) :: rho1s(nBas,nBas,nVVs),rho1t(nBas,nBas,nVVt)
|
||||
double precision,intent(in) :: Om2s(nOOs),Om2t(nOOt)
|
||||
double precision,intent(in) :: rho2s(nBas,nBas,nOOs),rho2t(nBas,nBas,nOOt)
|
||||
|
||||
! Local variables
|
||||
|
||||
@ -32,6 +32,12 @@ subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcG
|
||||
double precision,intent(inout):: Sig(nBas,nBas)
|
||||
double precision,intent(inout):: Z(nBas)
|
||||
|
||||
! Initialization
|
||||
|
||||
Sig(:,:) = 0d0
|
||||
Z(:) = 0d0
|
||||
EcGM = 0d0
|
||||
|
||||
!----------------------------------------------
|
||||
! Occupied part of the T-matrix self-energy
|
||||
!----------------------------------------------
|
||||
@ -39,14 +45,21 @@ subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcG
|
||||
do p=nC+1,nBas-nR
|
||||
do q=nC+1,nBas-nR
|
||||
do i=nC+1,nO
|
||||
do cd=1,nVV
|
||||
|
||||
eps = e(p) + e(i) - Om1(cd)
|
||||
num = rho1(p,i,cd)*rho1(q,i,cd)
|
||||
do cd=1,nVVs
|
||||
eps = e(p) + e(i) - Om1s(cd)
|
||||
num = rho1s(p,i,cd)*rho1s(q,i,cd)
|
||||
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
|
||||
|
||||
enddo
|
||||
|
||||
do cd=1,nVVt
|
||||
eps = e(p) + e(i) - Om1t(cd)
|
||||
num = rho1t(p,i,cd)*rho1t(q,i,cd)
|
||||
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
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
@ -58,14 +71,21 @@ subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcG
|
||||
do p=nC+1,nBas-nR
|
||||
do q=nC+1,nBas-nR
|
||||
do a=nO+1,nBas-nR
|
||||
do kl=1,nOO
|
||||
|
||||
eps = e(p) + e(a) - Om2(kl)
|
||||
num = rho2(p,a,kl)*rho2(q,a,kl)
|
||||
do kl=1,nOOs
|
||||
eps = e(p) + e(a) - Om2s(kl)
|
||||
num = rho2s(p,a,kl)*rho2s(q,a,kl)
|
||||
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
|
||||
|
||||
enddo
|
||||
|
||||
do kl=1,nOOt
|
||||
eps = e(p) + e(a) - Om2t(kl)
|
||||
num = rho2t(p,a,kl)*rho2t(q,a,kl)
|
||||
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
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
@ -76,26 +96,40 @@ subroutine GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcG
|
||||
|
||||
do i=nC+1,nO
|
||||
do j=nC+1,nO
|
||||
do cd=1,nVV
|
||||
|
||||
eps = e(i) + e(j) - Om1(cd)
|
||||
num = rho1(i,j,cd)*rho1(i,j,cd)
|
||||
do cd=1,nVVs
|
||||
eps = e(i) + e(j) - Om1s(cd)
|
||||
num = rho1s(i,j,cd)*rho1s(i,j,cd)
|
||||
EcGM = EcGM + num*eps/(eps**2 + eta**2)
|
||||
|
||||
enddo
|
||||
|
||||
do cd=1,nVVt
|
||||
eps = e(i) + e(j) - Om1t(cd)
|
||||
num = rho1t(i,j,cd)*rho1t(i,j,cd)
|
||||
EcGM = EcGM + num*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
|
||||
do a=nO+1,nBas-nR
|
||||
do b=nO+1,nBas-nR
|
||||
do kl=1,nOO
|
||||
|
||||
eps = e(a) + e(b) - Om2(kl)
|
||||
num = rho2(a,b,kl)*rho2(a,b,kl)
|
||||
do kl=1,nOOs
|
||||
eps = e(a) + e(b) - Om2s(kl)
|
||||
num = rho2s(a,b,kl)*rho2s(a,b,kl)
|
||||
EcGM = EcGM - num*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
do kl=1,nOOt
|
||||
eps = e(a) + e(b) - Om2t(kl)
|
||||
num = rho2t(a,b,kl)*rho2t(a,b,kl)
|
||||
EcGM = EcGM - num*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
|
||||
Z(:) = 1d0/(1d0 - Z(:))
|
||||
|
||||
end subroutine
|
||||
|
@ -1,4 +1,5 @@
|
||||
subroutine GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho2,EcGM,Sig,Z)
|
||||
subroutine GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,e,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t, &
|
||||
EcGM,Sig,Z)
|
||||
|
||||
! Compute diagonal of the correlation part of the T-matrix self-energy
|
||||
|
||||
@ -13,13 +14,13 @@ subroutine GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho
|
||||
integer,intent(in) :: nO
|
||||
integer,intent(in) :: nV
|
||||
integer,intent(in) :: nR
|
||||
integer,intent(in) :: nOO
|
||||
integer,intent(in) :: nVV
|
||||
integer,intent(in) :: nOOs,nOOt
|
||||
integer,intent(in) :: nVVs,nVVt
|
||||
double precision,intent(in) :: e(nBas)
|
||||
double precision,intent(in) :: Om1(nVV)
|
||||
double precision,intent(in) :: rho1(nBas,nBas,nVV)
|
||||
double precision,intent(in) :: Om2(nOO)
|
||||
double precision,intent(in) :: rho2(nBas,nBas,nOO)
|
||||
double precision,intent(in) :: Om1s(nVVs),Om1t(nVVt)
|
||||
double precision,intent(in) :: rho1s(nBas,nBas,nVVs),rho1t(nBas,nBas,nVVt)
|
||||
double precision,intent(in) :: Om2s(nOOs),Om2t(nOOt)
|
||||
double precision,intent(in) :: rho2s(nBas,nBas,nOOs),rho2t(nBas,nBas,nOOt)
|
||||
|
||||
! Local variables
|
||||
|
||||
@ -32,19 +33,31 @@ subroutine GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho
|
||||
double precision,intent(inout) :: Sig(nBas)
|
||||
double precision,intent(inout) :: Z(nBas)
|
||||
|
||||
!----------------------------------------------
|
||||
! Occupied part of the T-matrix self-energy
|
||||
!----------------------------------------------
|
||||
! Initialization
|
||||
|
||||
Sig(:) = 0d0
|
||||
Z(:) = 0d0
|
||||
EcGM = 0d0
|
||||
|
||||
!--------------------------------------!
|
||||
! Occupied part of the Tpp self-energy !
|
||||
!--------------------------------------!
|
||||
|
||||
do p=nC+1,nBas-nR
|
||||
do i=nC+1,nO
|
||||
do cd=1,nVV
|
||||
|
||||
eps = e(p) + e(i) - Om1(cd)
|
||||
num = rho1(p,i,cd)**2
|
||||
do cd=1,nVVs
|
||||
eps = e(p) + e(i) - Om1s(cd)
|
||||
num = rho1s(p,i,cd)**2
|
||||
Sig(p) = Sig(p) + num*eps/(eps**2 + eta**2)
|
||||
Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
|
||||
do cd=1,nVVt
|
||||
eps = e(p) + e(i) - Om1t(cd)
|
||||
num = rho1t(p,i,cd)**2
|
||||
Sig(p) = Sig(p) + num*eps/(eps**2 + eta**2)
|
||||
Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
@ -55,14 +68,21 @@ subroutine GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho
|
||||
|
||||
do p=nC+1,nBas-nR
|
||||
do a=nO+1,nBas-nR
|
||||
do kl=1,nOO
|
||||
|
||||
eps = e(p) + e(a) - Om2(kl)
|
||||
num = rho2(p,a,kl)**2
|
||||
do kl=1,nOOs
|
||||
eps = e(p) + e(a) - Om2s(kl)
|
||||
num = rho2s(p,a,kl)**2
|
||||
Sig(p) = Sig(p) + num*eps/(eps**2 + eta**2)
|
||||
Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
|
||||
enddo
|
||||
|
||||
do kl=1,nOOt
|
||||
eps = e(p) + e(a) - Om2t(kl)
|
||||
num = rho2t(p,a,kl)**2
|
||||
Sig(p) = Sig(p) + num*eps/(eps**2 + eta**2)
|
||||
Z(p) = Z(p) - num*(eps**2 - eta**2)/(eps**2 + eta**2)**2
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
|
||||
@ -72,26 +92,40 @@ subroutine GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Om1,rho1,Om2,rho
|
||||
|
||||
do i=nC+1,nO
|
||||
do j=nC+1,nO
|
||||
do cd=1,nVV
|
||||
|
||||
eps = e(i) + e(j) - Om1(cd)
|
||||
num = rho1(i,j,cd)**2
|
||||
do cd=1,nVVs
|
||||
eps = e(i) + e(j) - Om1s(cd)
|
||||
num = rho1s(i,j,cd)**2
|
||||
EcGM = EcGM + num*eps/(eps**2 + eta**2)
|
||||
|
||||
enddo
|
||||
|
||||
do cd=1,nVVt
|
||||
eps = e(i) + e(j) - Om1t(cd)
|
||||
num = rho1t(i,j,cd)**2
|
||||
EcGM = EcGM + num*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
|
||||
do a=nO+1,nBas-nR
|
||||
do b=nO+1,nBas-nR
|
||||
do kl=1,nOO
|
||||
|
||||
eps = e(a) + e(b) - Om2(kl)
|
||||
num = rho2(a,b,kl)**2
|
||||
do kl=1,nOOs
|
||||
eps = e(a) + e(b) - Om2s(kl)
|
||||
num = rho2s(a,b,kl)**2
|
||||
EcGM = EcGM - num*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
do kl=1,nOOt
|
||||
eps = e(a) + e(b) - Om2t(kl)
|
||||
num = rho2t(a,b,kl)**2
|
||||
EcGM = EcGM - num*eps/(eps**2 + eta**2)
|
||||
enddo
|
||||
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
|
||||
Z(:) = 1d0/(1d0 - Z(:))
|
||||
|
||||
end subroutine
|
||||
|
@ -1,5 +1,5 @@
|
||||
subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, &
|
||||
singlet,triplet,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int,PHF,cHF,eHF)
|
||||
subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,singlet,triplet, &
|
||||
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_MO,dipole_int,PHF,cHF,eHF)
|
||||
|
||||
! Perform eigenvalue self-consistent calculation with a T-matrix self-energy (evGT)
|
||||
|
||||
@ -21,6 +21,7 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
|
||||
logical,intent(in) :: dTDA
|
||||
logical,intent(in) :: singlet
|
||||
logical,intent(in) :: triplet
|
||||
logical,intent(in) :: linearize
|
||||
double precision,intent(in) :: eta
|
||||
logical,intent(in) :: regularize
|
||||
|
||||
@ -162,23 +163,14 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
|
||||
! Compute T-matrix version of the self-energy
|
||||
!----------------------------------------------
|
||||
|
||||
EcGM = 0d0
|
||||
Sig(:) = 0d0
|
||||
Z(:) = 0d0
|
||||
|
||||
iblock = 3
|
||||
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,ERI_MO,X1s,Y1s,rho1s,X2s,Y2s,rho2s)
|
||||
|
||||
call GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT,Om1s,rho1s,Om2s,rho2s,EcGM,Sig,Z)
|
||||
|
||||
iblock = 4
|
||||
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI_MO,X1t,Y1t,rho1t,X2t,Y2t,rho2t)
|
||||
|
||||
call GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT,Om1t,rho1t,Om2t,rho2t,EcGM,Sig,Z)
|
||||
|
||||
Z(:) = 1d0/(1d0 - Z(:))
|
||||
call GTpp_self_energy_diag(eta,nBas,nC,nO,nV,nR,nooS,nVVt,nOOt,nVVt,eGT,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t,EcGM,Sig,Z)
|
||||
|
||||
! Solve the quasi-particle equation
|
||||
|
||||
@ -186,8 +178,23 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
|
||||
! Solve the quasi-particle equation
|
||||
!----------------------------------------------
|
||||
|
||||
if(linearize) then
|
||||
|
||||
write(*,*) ' *** Quasiparticle energies obtained by linearization *** '
|
||||
write(*,*)
|
||||
|
||||
eGT(:) = eHF(:) + Sig(:)
|
||||
|
||||
else
|
||||
|
||||
write(*,*) ' *** Quasiparticle energies obtained by root search (experimental) *** '
|
||||
write(*,*)
|
||||
|
||||
call GTpp_QP_graph(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eHF,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t,eOld,eGT)
|
||||
|
||||
end if
|
||||
|
||||
! Convergence criteria
|
||||
|
||||
Conv = maxval(abs(eGT(:) - eOld(:)))
|
||||
|
@ -217,44 +217,14 @@ subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,T
|
||||
|
||||
! Compute correlation part of the self-energy
|
||||
|
||||
EcGM = 0d0
|
||||
Sig(:,:) = 0d0
|
||||
Z(:) = 0d0
|
||||
|
||||
iblock = 3
|
||||
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOs,nVVs,ERI_MO,X1s,Y1s,rho1s,X2s,Y2s,rho2s)
|
||||
|
||||
if(regularize) then
|
||||
|
||||
call regularized_self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT,Om1s,rho1s,Om2s,rho2s,EcGM,Sig)
|
||||
|
||||
call regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT,Om1s,rho1s,Om2s,rho2s,Z)
|
||||
|
||||
else
|
||||
|
||||
call GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,eGT,Om1s,rho1s,Om2s,rho2s,EcGM,Sig,Z)
|
||||
|
||||
end if
|
||||
|
||||
iblock = 4
|
||||
|
||||
call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI_MO,X1t,Y1t,rho1t,X2t,Y2t,rho2t)
|
||||
|
||||
if(regularize) then
|
||||
|
||||
call regularized_self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT,Om1t,rho1t,Om2t,rho2t,EcGM,Sig)
|
||||
|
||||
call regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT, &
|
||||
Om1t,rho1t,Om2t,rho2t,Z)
|
||||
|
||||
else
|
||||
|
||||
call GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOOt,nVVt,eGT,Om1t,rho1t,Om2t,rho2t,EcGM,Sig,Z)
|
||||
|
||||
end if
|
||||
|
||||
Z(:) = 1d0/(1d0 - Z(:))
|
||||
call GTpp_self_energy(eta,nBas,nC,nO,nV,nR,nOOs,nVVs,nOOt,nVVt,eGT,Om1s,rho1s,Om2s,rho2s, &
|
||||
Om1t,rho1t,Om2t,rho2t,EcGM,Sig,Z)
|
||||
|
||||
! Make correlation self-energy Hermitian and transform it back to AO basis
|
||||
|
||||
|
Loading…
Reference in New Issue
Block a user