LCPQ
/
quack
mirror of https://github.com/pfloos/quack
4
1
Fork 0
Code Releases Activity

BSE@GT potentially correct

This commit is contained in:
Pierre-Francois Loos 2021-12-21 13:49:41 +01:00
parent de2c45ab70
commit 2a5d6892b8
10 changed files with 34 additions and 27 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

4
input/methods
View File

@ -11,11 +11,11 @@
# RPA* RPAx* crRPA ppRPA
F F F F
# G0F2* evGF2* qsGF2* G0F3 evGF3
F F F F F
T F F F F
# G0W0* evGW* qsGW* ufG0W0 ufGW
T F F F F
# G0T0 evGT qsGT
F F F
T F F
# MCMP2
F
# * unrestricted version available

4
input/options
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
@ -15,6 +15,6 @@
# ACFDT: AC Kx XBS
F F F
# BSE: BSE dBSE dTDA evDyn
F F F F
T T T F
# MCMP2: nMC nEq nWalk dt nPrint iSeed doDrift
1000000 100000 10 0.3 10000 1234 T

2
mol/h2.xyz
View File

@ -1,4 +1,4 @@
2
H 0. 0. 0.
H 0. 0. 0.7
H 0. 0. 1.5

1
src/GT/Bethe_Salpeter_Tmatrix.f90
View File

@ -119,6 +119,7 @@ subroutine Bethe_Salpeter_Tmatrix(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,
call linear_response_Tmatrix(ispin,.false.,TDA,eta,nBas,nC,nO,nV,nR,nS,1d0,eGT,ERI,TA,TB, &
EcBSE(ispin),OmBSE(:,ispin),XpY_BSE(:,:,ispin),XmY_BSE(:,:,ispin))
call print_excitation('BSE@GT ',ispin,nS,OmBSE(:,ispin))
call print_transition_vectors(.true.,nBas,nC,nO,nV,nR,nS,dipole_int, &
OmBSE(:,ispin),XpY_BSE(:,:,ispin),XmY_BSE(:,:,ispin))

3
src/GT/G0T0.f90
View File

@ -199,9 +199,6 @@ subroutine G0T0(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,sing
if(BSE) then
! eG0T0(1) = -0.5507952119d0
! eG0T0(2) = +1.540259769d0
call Bethe_Salpeter_Tmatrix(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,nOOt,nVVt, &
Omega1s,X1s,Y1s,Omega2s,X2s,Y2s,rho1s,rho2s,Omega1t,X1t,Y1t,Omega2t,X2t,Y2t,rho1t,rho2t, &
ERI_MO,dipole_int,eHF,eG0T0,EcBSE)

12
src/GT/dynamic_Tmatrix_A.f90
View File

@ -58,14 +58,16 @@ subroutine dynamic_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGT,Omega1,O
chi = 0d0
do cd=1,nVV
chi = chi + rho1(i,j,cd)*rho1(a,b,cd)*Omega1(cd)/(Omega1(cd)**2 + eta**2)
eps = + Omega1(cd)
chi = chi + rho1(i,j,cd)*rho1(a,b,cd)*eps/(eps**2 + eta**2)
end do
do kl=1,nOO
chi = chi + rho2(i,j,kl)*rho2(a,b,kl)*Omega2(kl)/(Omega2(kl)**2 + eta**2)
do kl=1,nOO
eps = - Omega2(kl)
chi = chi + rho2(i,j,kl)*rho2(a,b,kl)*eps/(eps**2 + eta**2)
end do
A_dyn(ia,jb) = A_dyn(ia,jb) + 1d0*lambda*chi
A_dyn(ia,jb) = A_dyn(ia,jb) - 1d0*lambda*chi
chi = 0d0
@ -75,7 +77,7 @@ subroutine dynamic_Tmatrix_A(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,eGT,Omega1,O
end do
do kl=1,nOO
eps = + OmBSE - Omega2(kl) - (eGT(a) + eGT(b))
eps = + OmBSE + Omega2(kl) - (eGT(a) + eGT(b))
chi = chi + rho2(i,j,kl)*rho2(a,b,kl)*eps/(eps**2 + eta**2)
end do

13
src/GT/static_Tmatrix_TA.f90
View File

@ -26,6 +26,7 @@ subroutine static_Tmatrix_TA(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Omega1,r
! Local variables
double precision :: chi
double precision :: eps
integer :: i,j,a,b,ia,jb,kl,cd
! Output variables
@ -44,16 +45,18 @@ subroutine static_Tmatrix_TA(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Omega1,r
chi = 0d0
do cd=1,nVV
! chi = chi + lambda*rho1(i,j,cd)*rho1(a,b,cd)*Omega1(cd)/(Omega1(cd)**2 + eta**2)
chi = chi + rho1(i,j,cd)*rho1(a,b,cd)*Omega1(cd)/(Omega1(cd)**2 + eta**2)
eps = + Omega1(cd)
! chi = chi + lambda*rho1(i,j,cd)*rho1(a,b,cd)*eps/(eps**2 + eta**2)
chi = chi + rho1(i,j,cd)*rho1(a,b,cd)*eps/(eps**2 + eta**2)
enddo
do kl=1,nOO
! chi = chi + lambda*rho2(i,j,kl)*rho2(a,b,kl)*Omega2(kl)/(Omega2(kl)**2 + eta**2)
chi = chi - rho2(i,j,kl)*rho2(a,b,kl)*Omega2(kl)/(Omega2(kl)**2 + eta**2)
eps = - Omega2(kl)
! chi = chi - lambda*rho2(i,j,kl)*rho2(a,b,kl)*eps/(eps**2 + eta**2)
chi = chi + rho2(i,j,kl)*rho2(a,b,kl)*eps/(eps**2 + eta**2)
enddo
TA(ia,jb) = TA(ia,jb) - 1d0*lambda*chi
TA(ia,jb) = TA(ia,jb) + 1d0*lambda*chi
enddo
enddo

9
src/GT/static_Tmatrix_TB.f90
View File

@ -26,6 +26,7 @@ subroutine static_Tmatrix_TB(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Omega1,r
! Local variables
double precision :: chi
double precision :: eps
integer :: i,j,a,b,ia,jb,kl,cd
! Output variables
@ -44,16 +45,18 @@ subroutine static_Tmatrix_TB(eta,nBas,nC,nO,nV,nR,nS,nOO,nVV,lambda,ERI,Omega1,r
chi = 0d0
do cd=1,nVV
eps = Omega1(cd)
! chi = chi + lambda*rho1(i,b,cd)*rho1(a,j,cd)*Omega1(cd)/Omega1(cd)**2 + eta**2
chi = chi + rho1(i,b,cd)*rho1(a,j,cd)*Omega1(cd)/Omega1(cd)**2 + eta**2
chi = chi + rho1(i,b,cd)*rho1(a,j,cd)*eps/(eps**2 + eta**2)
enddo
do kl=1,nOO
eps = - Omega2(kl)
! chi = chi + lambda*rho2(i,b,kl)*rho2(a,j,kl)*Omega2(kl)/Omega2(kl)**2 + eta**2
chi = chi - rho2(i,b,kl)*rho2(a,j,kl)*Omega2(kl)/Omega2(kl)**2 + eta**2
chi = chi + rho2(i,b,kl)*rho2(a,j,kl)*eps/(eps**2 + eta**2)
enddo
TB(ia,jb) = TB(ia,jb) - 1d0*lambda*chi
TB(ia,jb) = TB(ia,jb) + 1d0*lambda*chi
enddo
enddo

7
src/GW/Bethe_Salpeter_A_matrix_dynamic.f90
View File

@ -48,15 +48,16 @@ subroutine Bethe_Salpeter_A_matrix_dynamic(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,Om
jb = jb + 1
chi = 0d0
do kc=1,maxS
chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*OmRPA(kc)/(OmRPA(kc)**2 + eta**2)
eps = OmRPA(kc)
chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps/(eps**2 + eta**2)
enddo
A_dyn(ia,jb) = A_dyn(ia,jb) - 4d0*lambda*chi
chi = 0d0
do kc=1,maxS
eps = + OmBSE - OmRPA(kc) - (eGW(a) - eGW(j))

6
src/LR/linear_response_pp.f90
View File

@ -47,8 +47,8 @@ subroutine linear_response_pp(ispin,TDA,nBas,nC,nO,nV,nR,nOO,nVV,lambda,e,ERI,Om
! Memory allocation
allocate(B(nVV,nOO),C(nVV,nVV),D(nOO,nOO),M(nOO+nVV,nOO+nVV),Z(nOO+nVV,nOO+nVV),Omega(nOO+nVV))
write(*,*) 'nOO', nOO
write(*,*) 'nVV', nVV
!write(*,*) 'nOO', nOO
!write(*,*) 'nVV', nVV
!-------------------------------------------------!
! Solve the p-p eigenproblem !
!-------------------------------------------------!
@ -88,7 +88,7 @@ write(*,*) 'nVV', nVV
M( 1:nVV ,nVV+1:nOO+nVV) = - B(1:nVV,1:nOO)
M(nVV+1:nOO+nVV, 1:nVV) = + transpose(B(1:nVV,1:nOO))
call matout(nOO+nVV,nOO+nVV,M)
!call matout(nOO+nVV,nOO+nVV,M)
! Diagonalize the p-h matrix
if(nOO+nVV > 0) call diagonalize_general_matrix(nOO+nVV,M,Omega,Z)