diff --git a/src/GT/G0T0eh.f90 b/src/GT/G0T0eh.f90 index b2d4900..25b0168 100644 --- a/src/GT/G0T0eh.f90 +++ b/src/GT/G0T0eh.f90 @@ -124,16 +124,9 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE, ! Compute GW self-energy ! !------------------------! - if(regularize) then - - write(*,*) 'Regularization not yet implemented at the G0T0eh level!' - stop + if(regularize) call GTeh_regularization(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rhoL,rhoR) - else - - call GTeh_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rhoL,rhoR,EcGM,Sig,Z) - - end if + call GTeh_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eHF,Om,rhoL,rhoR,EcGM,Sig,Z) !-----------------------------------! ! Solve the quasi-particle equation ! diff --git a/src/GT/G0T0pp.f90 b/src/GT/G0T0pp.f90 index c33dcbb..40684b8 100644 --- a/src/GT/G0T0pp.f90 +++ b/src/GT/G0T0pp.f90 @@ -153,6 +153,8 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,dopp ! Compute T-matrix version of the self-energy !---------------------------------------------- +! if(regularize) call GTpp_regularization(eta,nBas,nC,nO,nV,nR,nOO,nVV,eHF,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t) + 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) diff --git a/src/GT/evGTeh.f90 b/src/GT/evGTeh.f90 index c7f37bf..4e5e582 100644 --- a/src/GT/evGTeh.f90 +++ b/src/GT/evGTeh.f90 @@ -130,16 +130,9 @@ subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,d ! Compute correlation part of the self-energy - if(regularize) then + if(regularize) call GTeh_regularization(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR) -! call regularized_self_energy_correlation_diag(COHSEX,eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rho,EcGM,Sig) -! call renormalization_factor_SRG(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rho,Z) - - else - - call GTeh_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR,EcGM,Sig,Z) - - end if + call GTeh_self_energy_diag(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR,EcGM,Sig,Z) ! Solve the quasi-particle equation diff --git a/src/GT/evGTpp.f90 b/src/GT/evGTpp.f90 index d167f14..c1e4292 100644 --- a/src/GT/evGTpp.f90 +++ b/src/GT/evGTpp.f90 @@ -167,11 +167,15 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T iblock = 4 call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI,X1t,Y1t,rho1t,X2t,Y2t,rho2t) + !---------------------------------------------- + ! Compute T-matrix version of the self-energy + !---------------------------------------------- + + ! if(regularize) call GTpp_regularization(eta,nBas,nC,nO,nV,nR,nOO,nVV,eGT,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t) + 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 - !---------------------------------------------- ! Solve the quasi-particle equation !---------------------------------------------- diff --git a/src/GT/qsGTeh.f90 b/src/GT/qsGTeh.f90 index 3c6ae93..d299884 100644 --- a/src/GT/qsGTeh.f90 +++ b/src/GT/qsGTeh.f90 @@ -183,16 +183,9 @@ subroutine qsGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2, call GTeh_excitation_density(nBas,nC,nO,nR,nS,ERI_MO,XpY,XmY,rhoL,rhoR) - if(regularize) then + if(regularize) call GTeh_regularization(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR) -! call regularized_self_energy_correlation(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR,EcGM,Sig) -! call regularized_renormalization_factor(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR,Z) - - else - - call GTeh_self_energy(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR,EcGM,Sig,Z) - - endif + call GTeh_self_energy(eta,nBas,nC,nO,nV,nR,nS,eGT,Om,rhoL,rhoR,EcGM,Sig,Z) ! Make correlation self-energy Hermitian and transform it back to AO basis diff --git a/src/GT/qsGTpp.f90 b/src/GT/qsGTpp.f90 index 0974eb0..8bf51d8 100644 --- a/src/GT/qsGTpp.f90 +++ b/src/GT/qsGTpp.f90 @@ -223,6 +223,8 @@ subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,T iblock = 4 call GTpp_excitation_density(iblock,nBas,nC,nO,nV,nR,nOOt,nVVt,ERI_MO,X1t,Y1t,rho1t,X2t,Y2t,rho2t) + ! if(regularize) call GTpp_regularization(eta,nBas,nC,nO,nV,nR,nOO,nVV,eGT,Om1s,rho1s,Om2s,rho2s,Om1t,rho1t,Om2t,rho2t) + 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) diff --git a/src/GT/regularized_renormalization_factor_Tmatrix.f90 b/src/GT/regularized_renormalization_factor_Tmatrix.f90 deleted file mode 100644 index 68e311f..0000000 --- a/src/GT/regularized_renormalization_factor_Tmatrix.f90 +++ /dev/null @@ -1,71 +0,0 @@ -subroutine regularized_renormalization_factor_Tmatrix(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Omega1,rho1,Omega2,rho2,Z) - -! Compute renormalization factor of the T-matrix self-energy - - implicit none - include 'parameters.h' - -! Input variables - - double precision,intent(in) :: eta - integer,intent(in) :: nBas,nC,nO,nV,nR - integer,intent(in) :: nOO - integer,intent(in) :: nVV - double precision,intent(in) :: e(nBas) - double precision,intent(in) :: Omega1(nVV) - double precision,intent(in) :: rho1(nBas,nBas,nVV) - double precision,intent(in) :: Omega2(nOO) - double precision,intent(in) :: rho2(nBas,nBas,nOO) - -! Local variables - - integer :: i,a,p,cd,kl - double precision :: eps - - double precision :: kappa - double precision :: fk,dfk - -! Output variables - - double precision,intent(out) :: Z(nBas) - -!-----------------------------------------! -! Parameters for regularized calculations ! -!-----------------------------------------! - - kappa = 1d0 - -! Occupied part of the T-matrix self-energy - - do p=nC+1,nBas-nR - do i=nC+1,nO - do cd=1,nVV - - eps = e(p) + e(i) - Omega1(cd) - - fk = (1d0 - exp(-2d0*eps**2/kappa**2))/eps - dfk = - fk/eps + 4d0*kappa**2*exp(-2d0*eps**2/kappa**2) - Z(p) = Z(p) - rho1(p,i,cd)**2*dfk - - enddo - enddo - enddo - -! Virtual part of the T-matrix self-energy - - do p=nC+1,nBas-nR - do a=1,nV-nR - do kl=1,nOO - - eps = e(p) + e(nO+a) - Omega2(kl) - - fk = (1d0 - exp(-2d0*eps**2/kappa**2))/eps - dfk = - fk/eps + 4d0*kappa**2*exp(-2d0*eps**2/kappa**2) - - Z(p) = Z(p) - rho2(p,nO+a,kl)**2*dfk - - enddo - enddo - enddo - -end subroutine regularized_renormalization_factor_Tmatrix diff --git a/src/GT/regularized_self_energy_Tmatrix.f90 b/src/GT/regularized_self_energy_Tmatrix.f90 deleted file mode 100644 index ee807d9..0000000 --- a/src/GT/regularized_self_energy_Tmatrix.f90 +++ /dev/null @@ -1,80 +0,0 @@ -subroutine regularized_self_energy_Tmatrix(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Omega1,rho1,Omega2,rho2,SigT) - -! Compute the correlation part of the T-matrix self-energy - - 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) :: nOO - integer,intent(in) :: nVV - double precision,intent(in) :: e(nBas) - double precision,intent(in) :: Omega1(nVV) - double precision,intent(in) :: rho1(nBas,nBas,nVV) - double precision,intent(in) :: Omega2(nOO) - double precision,intent(in) :: rho2(nBas,nBas,nOO) - -! Local variables - - integer :: i,a,p,q,cd,kl - double precision :: eps - - double precision :: kappa - double precision :: fk - -! Output variables - - double precision,intent(inout) :: SigT(nBas,nBas) - -!-----------------------------------------! -! Parameters for regularized calculations ! -!-----------------------------------------! - - kappa = 1d0 - -!---------------------------------------------- -! Occupied part of the T-matrix self-energy -!---------------------------------------------- - - 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) - Omega1(cd) - fk = (1d0 - exp(-2d0*eps**2/kappa**2))/eps - - SigT(p,q) = SigT(p,q) + rho1(p,i,cd)*rho1(q,i,cd)*fk - - enddo - enddo - enddo - enddo - -!---------------------------------------------- - ! Virtual part of the T-matrix self-energy -!---------------------------------------------- - - 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) - Omega2(kl) - fk = (1d0 - exp(-2d0*eps**2/kappa**2))/eps - - SigT(p,q) = SigT(p,q) + rho2(p,a,kl)*rho2(q,a,kl)*fk - - enddo - enddo - enddo - enddo - -end subroutine regularized_self_energy_Tmatrix diff --git a/src/GT/regularized_self_energy_Tmatrix_diag.f90 b/src/GT/regularized_self_energy_Tmatrix_diag.f90 deleted file mode 100644 index 978bd2a..0000000 --- a/src/GT/regularized_self_energy_Tmatrix_diag.f90 +++ /dev/null @@ -1,76 +0,0 @@ -subroutine regularized_self_energy_Tmatrix_diag(eta,nBas,nC,nO,nV,nR,nOO,nVV,e,Omega1,rho1,Omega2,rho2,SigT) - -! Compute diagonal of the correlation part of the T-matrix self-energy - - 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) :: nOO - integer,intent(in) :: nVV - double precision,intent(in) :: e(nBas) - double precision,intent(in) :: Omega1(nVV) - double precision,intent(in) :: rho1(nBas,nBas,nVV) - double precision,intent(in) :: Omega2(nOO) - double precision,intent(in) :: rho2(nBas,nBas,nOO) - -! Local variables - - integer :: i,a,p,cd,kl - double precision :: eps - - double precision :: kappa - double precision :: fk - -! Output variables - - double precision,intent(inout) :: SigT(nBas) - -!-----------------------------------------! -! Parameters for regularized calculations ! -!-----------------------------------------! - - kappa = 1d0 - -!---------------------------------------------- -! Occupied part of the T-matrix self-energy -!---------------------------------------------- - - do p=nC+1,nBas-nR - do i=nC+1,nO - do cd=1,nVV - - eps = e(p) + e(i) - Omega1(cd) - fk = (1d0 - exp(-2d0*eps**2/kappa**2))/eps - - SigT(p) = SigT(p) + rho1(p,i,cd)**2*fk - - enddo - enddo - enddo - -!---------------------------------------------- -! Virtual part of the T-matrix self-energy -!---------------------------------------------- - - do p=nC+1,nBas-nR - do a=nO+1,nBas-nR - do kl=1,nOO - - eps = e(p) + e(a) - Omega2(kl) - fk = (1d0 - exp(-2d0*eps**2/kappa**2))/eps - - SigT(p) = SigT(p) + rho2(p,a,kl)**2*fk - - enddo - enddo - enddo - -end subroutine regularized_self_energy_Tmatrix_diag diff --git a/src/GW/GW_regularization.f90 b/src/GW/GW_regularization.f90 index b3bda9a..3919db5 100644 --- a/src/GW/GW_regularization.f90 +++ b/src/GW/GW_regularization.f90 @@ -35,7 +35,7 @@ subroutine GW_regularization(nBas,nC,nO,nR,nS,e,Om,rho) do m=1,nS do i=nC+1,nO - Dpim = e(p) - e(i) - Om(m) + Dpim = e(p) - e(i) + Om(m) kappa = 1d0 - exp(-Dpim*Dpim*s) rho(p,i,m) = kappa*rho(p,i,m) enddo