From 37d312c789de0da8666ec6d1d836c35580dc01aa Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Mon, 16 Sep 2024 11:00:19 +0200 Subject: [PATCH] ppBSE upfolded --- src/GW/RGW_phBSE.f90 | 6 +-- src/GW/RGW_ppBSE.f90 | 2 +- src/GW/RGW_ppBSE_dynamic_kernel_D.f90 | 47 ++++++++++++----------- src/GW/RGW_ppBSE_dynamic_perturbation.f90 | 6 +-- src/GW/RGW_ppBSE_upfolded.f90 | 26 ++++++------- 5 files changed, 44 insertions(+), 43 deletions(-) diff --git a/src/GW/RGW_phBSE.f90 b/src/GW/RGW_phBSE.f90 index 24b7eed..29bed89 100644 --- a/src/GW/RGW_phBSE.f90 +++ b/src/GW/RGW_phBSE.f90 @@ -146,9 +146,9 @@ subroutine RGW_phBSE(dophBSE2,exchange_kernel,TDA_W,TDA,dBSE,dTDA,singlet,triple ! Upfolded phBSE ! !----------------! - call RGW_phBSE_upfolded_sym(ispin,nOrb,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eW) +! call RGW_phBSE_upfolded_sym(ispin,nOrb,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eW) - call RGW_phBSE_upfolded(ispin,nOrb,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eGW) +! call RGW_phBSE_upfolded(ispin,nOrb,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eGW) end if @@ -185,7 +185,7 @@ subroutine RGW_phBSE(dophBSE2,exchange_kernel,TDA_W,TDA,dBSE,dTDA,singlet,triple ! Upfolded phBSE ! !----------------! - call RGW_phBSE_upfolded(ispin,nOrb,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eGW) +! call RGW_phBSE_upfolded(ispin,nOrb,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eGW) end if diff --git a/src/GW/RGW_ppBSE.f90 b/src/GW/RGW_ppBSE.f90 index b64fd70..8dd35a9 100644 --- a/src/GW/RGW_ppBSE.f90 +++ b/src/GW/RGW_ppBSE.f90 @@ -369,7 +369,7 @@ subroutine RGW_ppBSE(TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nOrb,nC,nO,nV,nR,nS ! Upfolded ppBSE ! !----------------! - call RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eGW) +! call RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho_RPA,OmRPA,eGW) deallocate(KB_sta,KC_sta,KD_sta) deallocate(Om1,X1,Y1,Om2,X2,Y2) diff --git a/src/GW/RGW_ppBSE_dynamic_kernel_D.f90 b/src/GW/RGW_ppBSE_dynamic_kernel_D.f90 index e56fa80..75024bb 100644 --- a/src/GW/RGW_ppBSE_dynamic_kernel_D.f90 +++ b/src/GW/RGW_ppBSE_dynamic_kernel_D.f90 @@ -56,29 +56,29 @@ subroutine RGW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,e do m=1,nS - dem = - OmBSE + eGW(k) - Om(m) + eGW(j) + dem = OmBSE - eGW(k) + Om(m) - eGW(i) num = rho(i,k,m)*rho(j,l,m) - KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) - ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 + KD_dyn(ij,kl) = KD_dyn(ij,kl) - num*dem/(dem**2 + eta**2) + ZD_dyn(ij,kl) = ZD_dyn(ij,kl) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 - dem = - OmBSE + eGW(k) - Om(m) + eGW(i) + dem = OmBSE - eGW(k) + Om(m) - eGW(j) num = rho(j,k,m)*rho(i,l,m) - KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) - ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 + KD_dyn(ij,kl) = KD_dyn(ij,kl) - num*dem/(dem**2 + eta**2) + ZD_dyn(ij,kl) = ZD_dyn(ij,kl) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 - dem = - OmBSE + eGW(l) - Om(m) + eGW(i) + dem = OmBSE - eGW(l) + Om(m) - eGW(j) num = rho(i,k,m)*rho(j,l,m) - KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) - ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 + KD_dyn(ij,kl) = KD_dyn(ij,kl) - num*dem/(dem**2 + eta**2) + ZD_dyn(ij,kl) = ZD_dyn(ij,kl) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 - dem = - OmBSE + eGW(l) - Om(m) + eGW(j) + dem = OmBSE - eGW(l) + Om(m) - eGW(i) num = rho(j,k,m)*rho(i,l,m) - KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) - ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 + KD_dyn(ij,kl) = KD_dyn(ij,kl) - num*dem/(dem**2 + eta**2) + ZD_dyn(ij,kl) = ZD_dyn(ij,kl) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 end do @@ -106,30 +106,31 @@ subroutine RGW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,e kl = kl + 1 do m=1,nS - dem = - OmBSE + eGW(k) - Om(m) + eGW(j) + + dem = OmBSE - eGW(k) + Om(m) - eGW(i) num = rho(i,k,m)*rho(j,l,m) - KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) - ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 - - dem = - OmBSE + eGW(k) - Om(m) + eGW(i) - num = rho(j,k,m)*rho(i,l,m) - KD_dyn(ij,kl) = KD_dyn(ij,kl) - num*dem/(dem**2 + eta**2) ZD_dyn(ij,kl) = ZD_dyn(ij,kl) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 - dem = - OmBSE + eGW(l) - Om(m) + eGW(i) - num = rho(i,k,m)*rho(j,l,m) + dem = OmBSE - eGW(k) + Om(m) - eGW(j) + num = rho(j,k,m)*rho(i,l,m) KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 - dem = - OmBSE + eGW(l) - Om(m) + eGW(j) - num = rho(j,k,m)*rho(i,l,m) + dem = OmBSE - eGW(l) + Om(m) - eGW(j) + num = rho(i,k,m)*rho(j,l,m) KD_dyn(ij,kl) = KD_dyn(ij,kl) - num*dem/(dem**2 + eta**2) ZD_dyn(ij,kl) = ZD_dyn(ij,kl) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 + dem = OmBSE - eGW(l) + Om(m) - eGW(i) + num = rho(j,k,m)*rho(i,l,m) + + KD_dyn(ij,kl) = KD_dyn(ij,kl) + num*dem/(dem**2 + eta**2) + ZD_dyn(ij,kl) = ZD_dyn(ij,kl) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2 + end do KD_dyn(ij,kl) = 2d0*KD_dyn(ij,kl) diff --git a/src/GW/RGW_ppBSE_dynamic_perturbation.f90 b/src/GW/RGW_ppBSE_dynamic_perturbation.f90 index 255cd85..e6c9393 100644 --- a/src/GW/RGW_ppBSE_dynamic_perturbation.f90 +++ b/src/GW/RGW_ppBSE_dynamic_perturbation.f90 @@ -114,14 +114,14 @@ subroutine RGW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nOrb,nC,nO,nV,nR,nS,nOO write(*,*) '---------------------------------------------------------------------------------------------------' kl = 0 - do ij=max(1,nOO+1-maxOO),nOO + do ij=nOO,max(1,nOO+1-maxOO),-1 kl = kl + 1 if(dTDA) then call RGW_ppBSE_dynamic_kernel_D(ispin,eta,nOrb,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))) + 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 @@ -144,7 +144,7 @@ subroutine RGW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nOrb,nC,nO,nV,nR,nS,nOO Om2_dyn(kl) = Z2_dyn(kl)*Om2_dyn(kl) write(*,'(2X,I5,5X,F15.6,5X,F15.6,5X,F15.6,5X,F15.6)') & - ij,Om2(ij)*HaToeV,(Om2(ij)+Om2_dyn(kl))*HaToeV,Om2_dyn(kl)*HaToeV,Z2_dyn(kl) + kl,Om2(ij)*HaToeV,(Om2(ij)+Om2_dyn(kl))*HaToeV,Om2_dyn(kl)*HaToeV,Z2_dyn(kl) end do write(*,*) '---------------------------------------------------------------------------------------------------' diff --git a/src/GW/RGW_ppBSE_upfolded.f90 b/src/GW/RGW_ppBSE_upfolded.f90 index 5e96d9d..445e509 100644 --- a/src/GW/RGW_ppBSE_upfolded.f90 +++ b/src/GW/RGW_ppBSE_upfolded.f90 @@ -25,7 +25,7 @@ subroutine RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho,Om,eGW) integer :: i,j,k,l integer :: a,b,c,d integer :: m,ij,kl,ijm - integer,parameter :: maxH = 20 + integer,parameter :: maxH = 100 double precision :: tmp,tmp1,tmp2,tmp3,tmp4 integer :: n1h,n1p,n2h,n2p,n1h1p,n3h1p,n3p1h,n2h2p,nH @@ -42,9 +42,9 @@ subroutine RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho,Om,eGW) ! Hello world write(*,*) - write(*,*)'**********************************************' - write(*,*)'| Unfolded ppBSE@GW calculation |' - write(*,*)'**********************************************' + write(*,*)'*********************************' + write(*,*)'* Upfolded ppBSE@GW Calculation *' + write(*,*)'*********************************' write(*,*) ! TDA for W @@ -118,8 +118,8 @@ subroutine RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho,Om,eGW) do l=k,nO kl = kl + 1 - H(ij,kl) = (eGW(i) + eGW(j))*Kronecker_delta(i,k)*Kronecker_delta(j,l) & - - (ERI(i,j,k,l) + ERI(i,j,l,k))/sqrt((1d0 + Kronecker_delta(i,j))*(1d0 + Kronecker_delta(k,l))) + H(ij,kl) = - (eGW(i) + eGW(j))*Kronecker_delta(i,k)*Kronecker_delta(j,l) & + + (ERI(i,j,k,l) + ERI(i,j,l,k))/sqrt((1d0 + Kronecker_delta(i,j))*(1d0 + Kronecker_delta(k,l))) end do end do @@ -145,8 +145,8 @@ subroutine RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho,Om,eGW) do l=k+1,nO kl = kl + 1 - H(ij,kl) = (eGW(i) + eGW(j))*Kronecker_delta(i,k)*Kronecker_delta(j,l) & - - (ERI(i,j,k,l) - ERI(i,j,l,k)) + H(ij,kl) = - (eGW(i) + eGW(j))*Kronecker_delta(i,k)*Kronecker_delta(j,l) & + + (ERI(i,j,k,l) - ERI(i,j,l,k)) end do end do @@ -172,14 +172,14 @@ subroutine RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho,Om,eGW) kl = kl + 1 tmp1 = sqrt(2d0)*Kronecker_delta(j,l)*rho(i,k,m) - tmp2 = sqrt(2d0)*Kronecker_delta(j,l)*rho(i,l,m) + tmp2 = sqrt(2d0)*Kronecker_delta(j,k)*rho(i,l,m) tmp3 = sqrt(2d0)*Kronecker_delta(i,l)*rho(j,k,m) tmp4 = sqrt(2d0)*Kronecker_delta(i,k)*rho(j,l,m) - H(n2h+0*n3h1p+ijm,kl ) = -tmp1 - H(kl ,n2h+0*n3h1p+ijm) = +tmp3 + H(n2h+0*n3h1p+ijm,kl ) = +tmp1 + H(kl ,n2h+0*n3h1p+ijm) = +tmp4 - H(n2h+1*n3h1p+ijm,kl ) = -tmp1 + H(n2h+1*n3h1p+ijm,kl ) = +tmp1 H(kl ,n2h+1*n3h1p+ijm) = +tmp4 H(n2h+2*n3h1p+ijm,kl ) = -tmp2 @@ -254,7 +254,7 @@ subroutine RGW_ppBSE_upfolded(ispin,nOrb,nC,nO,nV,nR,nS,ERI,rho,Om,eGW) do s=1,min(nH,maxH) if(Z(s) > 1d-7) & write(*,'(1X,A1,1X,I3,1X,A1,1X,F15.6,1X,A1,1X,F15.6,1X,A1,1X)') & - '|',s,'|',OmBSE(s)*HaToeV,'|',Z(s),'|' + '|',s,'|',-OmBSE(s)*HaToeV,'|',Z(s),'|' end do write(*,*)'-------------------------------------------'