4
1
mirror of https://github.com/pfloos/quack synced 2024-07-05 19:06:01 +02:00

modifs with antoine

This commit is contained in:
Pierre-Francois Loos 2024-07-02 14:41:12 +02:00
parent a2a63862c1
commit 5c5c1aa6d6
8 changed files with 44 additions and 25 deletions

View File

@ -25,6 +25,7 @@ subroutine GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambd
! Local variables
double precision,external :: Kronecker_delta
double precision :: dem,num
integer :: m
integer :: a,b,i,j
@ -62,7 +63,7 @@ subroutine GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambd
dem = eGW(j) - Om(m) - eGW(a)
num = rho(b,i,m)*rho(a,j,m)
KB_dyn(ab,ij) = KB_dyn(ab,ij) - num*dem/(dem**2 + eta**2)
KB_dyn(ab,ij) = KB_dyn(ab,ij) + num*dem/(dem**2 + eta**2)
dem = eGW(i) - Om(m) - eGW(a)
num = rho(a,i,m)*rho(b,j,m)
@ -72,7 +73,9 @@ subroutine GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambd
dem = eGW(i) - Om(m) - eGW(b)
num = rho(b,i,m)*rho(a,j,m)
KB_dyn(ab,ij) = KB_dyn(ab,ij) - num*dem/(dem**2 + eta**2)
KB_dyn(ab,ij) = KB_dyn(ab,ij) + num*dem/(dem**2 + eta**2)
KB_dyn(ab,ij) = 2d0*KB_dyn(ab,ij)/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(i,j)))
end do
@ -118,6 +121,8 @@ subroutine GW_ppBSE_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambd
KB_dyn(ab,ij) = KB_dyn(ab,ij) - num*dem/(dem**2 + eta**2)
KB_dyn(ab,ij) = 2d0*KB_dyn(ab,ij)
end do
end do

View File

@ -24,6 +24,7 @@ subroutine GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,eG
! Local variables
double precision,external :: Kronecker_delta
double precision :: dem,num
integer :: m
integer :: a,b,c,d
@ -64,8 +65,8 @@ subroutine GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,eG
dem = OmBSE - eGW(c) - Om(m) - eGW(a)
num = rho(b,c,m)*rho(a,d,m)
KC_dyn(ab,cd) = KC_dyn(ab,cd) - num*dem/(dem**2 + eta**2)
ZC_dyn(ab,cd) = ZC_dyn(ab,cd) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2
KC_dyn(ab,cd) = KC_dyn(ab,cd) + num*dem/(dem**2 + eta**2)
ZC_dyn(ab,cd) = ZC_dyn(ab,cd) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2
dem = OmBSE - eGW(d) - Om(m) - eGW(a)
num = rho(a,c,m)*rho(b,d,m)
@ -76,8 +77,11 @@ subroutine GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,eG
dem = OmBSE - eGW(d) - Om(m) - eGW(b)
num = rho(b,c,m)*rho(a,d,m)
KC_dyn(ab,cd) = KC_dyn(ab,cd) - num*dem/(dem**2 + eta**2)
ZC_dyn(ab,cd) = ZC_dyn(ab,cd) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2
KC_dyn(ab,cd) = KC_dyn(ab,cd) + num*dem/(dem**2 + eta**2)
ZC_dyn(ab,cd) = ZC_dyn(ab,cd) - num*(dem**2 - eta**2)/(dem**2 + eta**2)**2
KC_dyn(ab,cd) = 2d0*KC_dyn(ab,cd)/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))
ZC_dyn(ab,cd) = 2d0*ZC_dyn(ab,cd)/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))
end do
@ -127,6 +131,9 @@ subroutine GW_ppBSE_dynamic_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,eG
KC_dyn(ab,cd) = KC_dyn(ab,cd) - num*dem/(dem**2 + eta**2)
ZC_dyn(ab,cd) = ZC_dyn(ab,cd) + num*(dem**2 - eta**2)/(dem**2 + eta**2)**2
KC_dyn(ab,cd) = 2d0*KC_dyn(ab,cd)
ZC_dyn(ab,cd) = 2d0*ZC_dyn(ab,cd)
end do
end do

View File

@ -24,6 +24,7 @@ subroutine GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,eG
! Local variables
double precision,external :: Kronecker_delta
double precision :: dem,num
integer :: m
integer :: i,j,k,l
@ -64,8 +65,8 @@ subroutine GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,eG
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
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)
@ -76,8 +77,11 @@ subroutine GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,eG
dem = - OmBSE + eGW(l) - 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
KD_dyn(ij,kl) = 2d0*KD_dyn(ij,kl)/sqrt((1d0 + Kronecker_delta(i,j))*(1d0 + Kronecker_delta(k,l)))
ZD_dyn(ij,kl) = 2d0*ZD_dyn(ij,kl)/sqrt((1d0 + Kronecker_delta(i,j))*(1d0 + Kronecker_delta(k,l)))
end do
@ -127,6 +131,9 @@ subroutine GW_ppBSE_dynamic_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,eG
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) = 2d0*KD_dyn(ij,kl)
ZD_dyn(ij,kl) = 2d0*ZD_dyn(ij,kl)
end do
end do

View File

@ -96,7 +96,7 @@ subroutine GW_ppBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOO,
write(*,'(2X,A5,1X,A20,1X,A20,1X,A20,1X,A20)') '#','Static (eV)','Dynamic (eV)','Correction (eV)','Renorm. (eV)'
write(*,*) '---------------------------------------------------------------------------------------------------'
do ij=1,min(nOO,maxOO)
do ij=max(1,nOO+1-maxOO),nOO
! if(.not.dTDA) call GW_ppBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,1d0,eGW,OmRPA,rho_RPA,OmBSE(ab),KB_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)

View File

@ -56,10 +56,10 @@ subroutine GW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda
do m=1,nS
eps = Om(m)**2 + eta**2
chi = chi - rho(i,a,m)*rho(j,b,m)*Om(m)/eps &
+ rho(i,b,m)*rho(a,j,m)*Om(m)/eps
- rho(i,b,m)*rho(a,j,m)*Om(m)/eps
end do
KB(ab,ij) = 2d0*lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(i,j)))
KB(ab,ij) = 4d0*lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(i,j)))
end do
end do
@ -90,7 +90,7 @@ subroutine GW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda
+ rho(i,b,m)*rho(a,j,m)*Om(m)/eps
end do
KB(ab,ij) = 2d0*lambda*chi
KB(ab,ij) = 4d0*lambda*chi
end do
end do
@ -121,7 +121,7 @@ subroutine GW_ppBSE_static_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda
+ rho(i,b,m)*rho(a,j,m)*Om(m)/eps
end do
KB(ab,ij) = lambda*chi
KB(ab,ij) = 2d0*lambda*chi
end do
end do

View File

@ -55,10 +55,10 @@ subroutine GW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,ERI
do m=1,nS
eps = Om(m)**2 + eta**2
chi = chi - rho(a,c,m)*rho(b,d,m)*Om(m)/eps &
+ rho(a,d,m)*rho(b,c,m)*Om(m)/eps
- rho(a,d,m)*rho(b,c,m)*Om(m)/eps
end do
KC(ab,cd) = 2d0*lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))
KC(ab,cd) = 4d0*lambda*chi/sqrt((1d0 + Kronecker_delta(a,b))*(1d0 + Kronecker_delta(c,d)))
end do
end do
@ -89,7 +89,7 @@ subroutine GW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,ERI
+ rho(a,d,m)*rho(b,c,m)*Om(m)/eps
end do
KC(ab,cd) = 2d0*lambda*chi
KC(ab,cd) = 4d0*lambda*chi
end do
end do
@ -120,7 +120,7 @@ subroutine GW_ppBSE_static_kernel_C(ispin,eta,nBas,nC,nO,nV,nR,nS,nVV,lambda,ERI
+ rho(a,d,m)*rho(b,c,m)*Om(m)/eps
end do
KC(ab,cd) = lambda*chi
KC(ab,cd) = 2d0*lambda*chi
end do
end do

View File

@ -55,10 +55,10 @@ subroutine GW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,ERI
do m=1,nS
eps = Om(m)**2 + eta**2
chi = chi - rho(i,k,m)*rho(j,l,m)*Om(m)/eps &
+ rho(i,l,m)*rho(j,k,m)*Om(m)/eps
- rho(i,l,m)*rho(j,k,m)*Om(m)/eps
end do
KD(ij,kl) = 2d0*lambda*chi/sqrt((1d0 + Kronecker_delta(i,j))*(1d0 + Kronecker_delta(k,l)))
KD(ij,kl) = 4d0*lambda*chi/sqrt((1d0 + Kronecker_delta(i,j))*(1d0 + Kronecker_delta(k,l)))
end do
end do
@ -89,7 +89,7 @@ subroutine GW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,ERI
+ rho(i,l,m)*rho(j,k,m)*Om(m)/eps
end do
KD(ij,kl) = 2d0*lambda*chi
KD(ij,kl) = 4d0*lambda*chi
end do
end do
@ -120,7 +120,7 @@ subroutine GW_ppBSE_static_kernel_D(ispin,eta,nBas,nC,nO,nV,nR,nS,nOO,lambda,ERI
+ rho(i,l,m)*rho(j,k,m)*Om(m)/eps
end do
KD(ij,kl) = lambda*chi
KD(ij,kl) = 2d0*lambda*chi
end do
end do

View File

@ -76,7 +76,7 @@ subroutine ppRRPA(dotest,TDA,doACFDT,singlet,triplet,nBas,nC,nO,nV,nR,ENuc,ERHF,
call ppLR(TDA,nOO,nVV,Bpp,Cpp,Dpp,Om1,X1,Y1,Om2,X2,Y2,EcRPA(ispin))
! call print_transition_vectors_pp(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call ppLR_transition_vectors(.true.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call print_excitation_energies('ppRPA@RHF','2p (singlet)',nVV,Om1)
call print_excitation_energies('ppRPA@RHF','2h (singlet)',nOO,Om2)
@ -108,7 +108,7 @@ subroutine ppRRPA(dotest,TDA,doACFDT,singlet,triplet,nBas,nC,nO,nV,nR,ENuc,ERHF,
call ppLR(TDA,nOO,nVV,Bpp,Cpp,Dpp,Om1,X1,Y1,Om2,X2,Y2,EcRPA(ispin))
! call print_transition_vectors_pp(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call ppLR_transition_vectors(.false.,nBas,nC,nO,nV,nR,nOO,nVV,dipole_int,Om1,X1,Y1,Om2,X2,Y2)
call print_excitation_energies('ppRPA@RHF','2p (triplet)',nVV,Om1)
call print_excitation_energies('ppRPA@RHF','2h (triplet)',nOO,Om2)