LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-07-06 11:26:01 +02:00
Code Releases Activity

remove evDyn and cleanup dBSE@GW

Browse Source
This commit is contained in:
Pierre-Francois Loos 2023-07-21 10:21:54 +02:00
parent b53a1505be
commit 46ab457dc2
38 changed files with 195 additions and 661 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
input/options
View File

@ -13,6 +13,6 @@
# GT: maxSCF thresh DIIS n_diis lin eta TDA_T reg # 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 T 0.0 F F
# ACFDT: AC Kx XBS # ACFDT: AC Kx XBS
T T F F F T
# BSE: phBSE phBSE2 ppBSE dBSE dTDA evDyn # BSE: phBSE phBSE2 ppBSE dBSE dTDA
T T T T F F T T T T F

5
src/GF/G0F2.f90
View File

@ -1,4 +1,4 @@
subroutine G0F2(BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,linearize,eta,regularize, & subroutine G0F2(BSE,TDA,dBSE,dTDA,singlet,triplet,linearize,eta,regularize, &
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF) nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
! Perform a one-shot second-order Green function calculation ! Perform a one-shot second-order Green function calculation
@ -12,7 +12,6 @@ subroutine G0F2(BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,linearize,eta,regularize
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
logical,intent(in) :: linearize logical,intent(in) :: linearize
@ -95,7 +94,7 @@ subroutine G0F2(BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,linearize,eta,regularize
if(BSE) then if(BSE) then
call GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF2,EcBSE) call GF2_phBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF2,EcBSE)
end if end if

29
src/GF/GF2_phBSE2.f90
View File

@ -1,4 +1,4 @@
subroutine GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF,EcBSE) subroutine GF2_phBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF,EcBSE)
! Compute the second-order Bethe-Salpeter excitation energies ! Compute the second-order Bethe-Salpeter excitation energies
@ -10,7 +10,6 @@ subroutine GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
@ -75,22 +74,11 @@ subroutine GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
! Compute dynamic correction for BSE via perturbation theory ! Compute dynamic correction for BSE via perturbation theory
if(dBSE) then if(dBSE) &
if(evDyn) then
call GF2_phBSE2_dynamic_perturbation_iterative(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF, &
KA_sta,KB_sta,OmBSE,XpY,XmY)
else
call GF2_phBSE2_dynamic_perturbation(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF,KA_sta,KB_sta,OmBSE,XpY,XmY) call GF2_phBSE2_dynamic_perturbation(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF,KA_sta,KB_sta,OmBSE,XpY,XmY)
end if end if
end if
end if
!------------------- !-------------------
! Triplet manifold ! Triplet manifold
!------------------- !-------------------
@ -119,20 +107,9 @@ subroutine GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,n
! Compute dynamic correction for BSE via perturbation theory ! Compute dynamic correction for BSE via perturbation theory
if(dBSE) then if(dBSE) &
if(evDyn) then
call GF2_phBSE2_dynamic_perturbation_iterative(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF, &
KA_sta,KB_sta,OmBSE,XpY,XmY)
else
call GF2_phBSE2_dynamic_perturbation(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF,KA_sta,KB_sta,OmBSE,XpY,XmY) call GF2_phBSE2_dynamic_perturbation(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF,KA_sta,KB_sta,OmBSE,XpY,XmY)
end if end if
end if
end if
end subroutine end subroutine

4
src/GF/GF2_phBSE2_dynamic_perturbation.f90
View File

@ -47,8 +47,8 @@ subroutine GF2_phBSE2_dynamic_perturbation(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ER
! Memory allocation ! Memory allocation
allocate(OmDyn(nS),ZDyn(nS),X(nS),Y(nS),KAp_dyn(nS,nS),ZAp_dyn(nS,nS)) allocate(OmDyn(nS),ZDyn(nS),X(nS),Y(nS),KAp_dyn(nS,nS),ZAp_dyn(nS,nS), &
allocate(KAm_dyn(nS,nS),ZAm_dyn(nS,nS),KB_dyn(nS,nS)) KAm_dyn(nS,nS),ZAm_dyn(nS,nS),KB_dyn(nS,nS))
if(dTDA) then if(dTDA) then
write(*,*) write(*,*)

155
src/GF/GF2_phBSE2_dynamic_perturbation_iterative.f90
View File

@ -1,155 +0,0 @@
subroutine GF2_phBSE2_dynamic_perturbation_iterative(dTDA,ispin,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int, &
eHF,eGF,A_sta,B_sta,OmBSE,XpY,XmY)
! Compute self-consistently the dynamical effects via perturbation theory for BSE2
implicit none
include 'parameters.h'
! Input variables
logical,intent(in) :: dTDA
integer,intent(in) :: ispin
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) :: nS
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: dipole_int(nBas,nBas,ncart)
double precision,intent(in) :: eHF(nBas)
double precision,intent(in) :: eGF(nBas)
double precision,intent(in) :: A_sta(nS,nS)
double precision,intent(in) :: B_sta(nS,nS)
double precision,intent(in) :: OmBSE(nS)
double precision,intent(in) :: XpY(nS,nS)
double precision,intent(in) :: XmY(nS,nS)
! Local variables
integer :: ia
integer,parameter :: maxS = 10
double precision :: gapGF
integer :: nSCF
integer :: maxSCF = 10
double precision :: Conv
double precision :: thresh = 1d-3
double precision,allocatable :: OmDyn(:)
double precision,allocatable :: ZDyn(:)
double precision,allocatable :: OmOld(:)
double precision,allocatable :: X(:)
double precision,allocatable :: Y(:)
double precision,allocatable :: Ap_dyn(:,:)
double precision,allocatable :: Am_dyn(:,:)
double precision,allocatable :: ZAp_dyn(:,:)
double precision,allocatable :: ZAm_dyn(:,:)
double precision,allocatable :: B_dyn(:,:)
! Memory allocation
allocate(OmDyn(nS),OmOld(nS),ZDyn(nS),X(nS),Y(nS),Ap_dyn(nS,nS),ZAp_dyn(nS,nS))
if(.not.dTDA) allocate(Am_dyn(nS,nS),ZAm_dyn(nS,nS),B_dyn(nS,nS))
! Print main components of transition vectors
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,OmBSE,XpY,XmY)
if(dTDA) then
write(*,*)
write(*,*) '*** dynamical TDA activated ***'
write(*,*)
end if
gapGF = eGF(nO+1) - eGF(nO)
Conv = 1d0
nSCF = 0
OmOld(:) = OmBSE(:)
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*) ' First-order dynamical correction to static Bethe-Salpeter excitation energies '
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(A57,F10.6,A3)') ' BSE2 neutral excitation must be lower than the GF gap = ',gapGF*HaToeV,' eV'
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*)
do while(Conv > thresh .and. nSCF < maxSCF)
nSCF = nSCF + 1
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A15,I3)') 'Iteration n.',nSCF
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A5,1X,A20,1X,A20,1X,A20,A20)') '#','Static (eV)','Dynamic (eV)','Correction (eV)','Convergence (eV)'
write(*,*) '---------------------------------------------------------------------------------------------------'
do ia=1,min(nS,maxS)
X(:) = 0.5d0*(XpY(ia,:) + XmY(ia,:))
Y(:) = 0.5d0*(XpY(ia,:) - XmY(ia,:))
! Resonant part of the BSE correction
call GF2_phBSE2_dynamic_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,eGF,+OmOld(ia),Ap_dyn,ZAp_dyn)
if(dTDA) then
OmDyn(ia) = dot_product(X,matmul(Ap_dyn - A_sta,X))
ZDyn(ia) = dot_product(X,matmul(ZAp_dyn,X))
else
! Anti-resonant part of the BSE correction
call GF2_phBSE2_dynamic_kernel_A(ispin,eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,eGF,-OmOld(ia),Am_dyn,ZAm_dyn)
call GF2_phBSE2_dynamic_kernel_B(ispin,eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,eGF,B_dyn)
ZDyn(ia) = dot_product(X,matmul(ZAp_dyn,X)) &
+ dot_product(Y,matmul(ZAm_dyn,Y))
OmDyn(ia) = dot_product(X,matmul(Ap_dyn - A_sta,X)) &
- dot_product(Y,matmul(Am_dyn - A_sta,Y)) &
+ dot_product(X,matmul(B_dyn - B_sta,Y)) &
- dot_product(Y,matmul(B_dyn - B_sta,X))
end if
write(*,'(2X,I5,5X,F15.6,5X,F15.6,5X,F15.6,5X,F15.6)') &
ia,OmBSE(ia)*HaToeV,(OmBSE(ia)+OmDyn(ia))*HaToeV,OmDyn(ia)*HaToeV,(OmBSE(ia) + OmDyn(ia) - OmOld(ia))*HaToeV
end do
Conv = maxval(abs(OmBSE(:) + OmDyn(:) - OmOld(:)))*HaToeV
OmOld(:) = OmBSE(:) + OmDyn(:)
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A20,1X,F10.6)') ' Convergence = ',Conv
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*)
end do
! Did it actually converge?
if(nSCF == maxSCF) then
write(*,*)
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)' Convergence failed '
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)
endif
end subroutine

3
src/GF/UG0F2.f90
View File

@ -1,4 +1,4 @@
subroutine UG0F2(BSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, & subroutine UG0F2(BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF, &
S,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,eHF) S,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,eHF)
! Perform unrestricted G0W0 calculation ! Perform unrestricted G0W0 calculation
@ -13,7 +13,6 @@ subroutine UG0F2(BSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,linearize,eta,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
logical,intent(in) :: linearize logical,intent(in) :: linearize

5
src/GF/evGF2.f90
View File

@ -1,4 +1,4 @@
subroutine evGF2(BSE,TDA,dBSE,dTDA,evDyn,maxSCF,thresh,max_diis,singlet,triplet, & subroutine evGF2(BSE,TDA,dBSE,dTDA,maxSCF,thresh,max_diis,singlet,triplet, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF) linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI,dipole_int,eHF)
! Perform eigenvalue self-consistent second-order Green function calculation ! Perform eigenvalue self-consistent second-order Green function calculation
@ -12,7 +12,6 @@ subroutine evGF2(BSE,TDA,dBSE,dTDA,evDyn,maxSCF,thresh,max_diis,singlet,triplet,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
integer,intent(in) :: maxSCF integer,intent(in) :: maxSCF
double precision,intent(in) :: thresh double precision,intent(in) :: thresh
integer,intent(in) :: max_diis integer,intent(in) :: max_diis
@ -143,7 +142,7 @@ subroutine evGF2(BSE,TDA,dBSE,dTDA,evDyn,maxSCF,thresh,max_diis,singlet,triplet,
if(BSE) then if(BSE) then
call GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF2,EcBSE) call GF2_phBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eHF,eGF2,EcBSE)
end if end if

3
src/GF/evUGF2.f90
View File

@ -1,4 +1,4 @@
subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip, & subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb, & eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,cHF,eHF) dipole_int_aa,dipole_int_bb,cHF,eHF)
@ -18,7 +18,6 @@ subroutine evUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,spin_conserved,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
double precision,intent(in) :: eta double precision,intent(in) :: eta

5
src/GF/qsGF2.f90
View File

@ -1,4 +1,4 @@
subroutine qsGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet, & subroutine qsGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,singlet,triplet, &
eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, & eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
@ -16,7 +16,6 @@ subroutine qsGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -262,7 +261,7 @@ subroutine qsGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,
if(BSE) then if(BSE) then
call GF2_phBSE2(TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO,eHF,eGF2,EcBSE) call GF2_phBSE2(TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO,eHF,eGF2,EcBSE)
write(*,*) write(*,*)
write(*,*)'-------------------------------------------------------------------------------' write(*,*)'-------------------------------------------------------------------------------'

3
src/GF/qsUGF2.f90
View File

@ -1,4 +1,4 @@
subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,regularize, & subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,regularize, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EUHF,S,X,T,V,Hc,ERI_AO, & nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EUHF,S,X,T,V,Hc,ERI_AO, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF)
@ -16,7 +16,6 @@ subroutine qsUGF2(maxSCF,thresh,max_diis,BSE,TDA,dBSE,dTDA,evDyn,spin_conserved,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
double precision,intent(in) :: eta double precision,intent(in) :: eta

5
src/GT/G0T0eh.f90
View File

@ -1,4 +1,4 @@
subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,evDyn,doppBSE, & subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE, &
singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF, & singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF, &
ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF,Vxc) ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF,Vxc)
@ -20,7 +20,6 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
logical,intent(in) :: linearize logical,intent(in) :: linearize
@ -186,7 +185,7 @@ subroutine G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,
! if(BSE) then ! if(BSE) then
! call Bethe_Salpeter(BSE2,TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE) ! call Bethe_Salpeter(BSE2,TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE)
! if(exchange_kernel) then ! if(exchange_kernel) then
! !

7
src/GT/G0T0pp.f90
View File

@ -1,4 +1,4 @@
subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,ppBSE,singlet,triplet, & subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,ppBSE,singlet,triplet, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF,Vxc) linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF,Vxc)
! Perform one-shot calculation with a T-matrix self-energy (G0T0) ! Perform one-shot calculation with a T-matrix self-energy (G0T0)
@ -17,7 +17,6 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
logical,intent(in) :: linearize logical,intent(in) :: linearize
@ -249,7 +248,7 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
if(BSE) then if(BSE) then
call GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, & 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, & Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,rho1aa,rho2aa, &
ERI_MO,dipole_int,eHF,eGT,EcBSE) ERI_MO,dipole_int,eHF,eGT,EcBSE)
@ -311,7 +310,7 @@ subroutine G0T0pp(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn,pp
if(ppBSE) then if(ppBSE) then
call GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nOOab,nVVab,nOOaa,nVVaa, & 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, & Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,rho1aa,rho2aa, &
ERI_MO,dipole_int,eHF,eGT,EcppBSE) ERI_MO,dipole_int,eHF,eGT,EcppBSE)

32
src/GT/GTpp_phBSE.f90
View File

@ -1,4 +1,4 @@
subroutine GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, & subroutine 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, & Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,rho1aa,rho2aa, &
ERI,dipole_int,eT,eGT,EcBSE) ERI,dipole_int,eT,eGT,EcBSE)
@ -13,7 +13,6 @@ subroutine GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,n
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
@ -143,22 +142,12 @@ subroutine GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,n
call print_excitation('phBSE@GTpp ',ispin,nS,OmBSE) call print_excitation('phBSE@GTpp ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE) call print_transition_vectors_ph(.true.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
if(dBSE) then ! Compute dynamic correction for BSE via renormalized perturbation theory
! Compute dynamic correction for BSE via perturbation theory (iterative or renormalized)
if(evDyn) then
print*, ' Iterative dynamical correction for BSE@GT NYI'
else
if(dBSE) &
call GTpp_phBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, & call GTpp_phBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, &
Om1ab,Om2ab,Om1aa,Om2aa,rho1ab,rho2ab,rho1aa,rho2aa,eT,eGT, & Om1ab,Om2ab,Om1aa,Om2aa,rho1ab,rho2ab,rho1aa,rho2aa,eT,eGT, &
dipole_int,OmBSE,XpY_BSE,XmY_BSE,TAab,TAaa) dipole_int,OmBSE,XpY_BSE,XmY_BSE,TAab,TAaa)
end if
end if
end if end if
@ -183,23 +172,12 @@ subroutine GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,n
call print_excitation('phBSE@GTpp ',ispin,nS,OmBSE) call print_excitation('phBSE@GTpp ',ispin,nS,OmBSE)
call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE) call print_transition_vectors_ph(.false.,nBas,nC,nO,nV,nR,nS,dipole_int,OmBSE,XpY_BSE,XmY_BSE)
if(dBSE) then ! Compute dynamic correction for BSE via renormalized perturbation theory
! Compute dynamic correction for BSE via perturbation theory (iterative or renormalized)
if(evDyn) then
print*, ' Iterative dynamical correction for BSE@GT NYI'
else
if(dBSE) &
call GTpp_phBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, & call GTpp_phBSE_dynamic_perturbation(ispin,dTDA,eta,nBas,nC,nO,nV,nR,nS,nOOab,nVVab,nOOaa,nVVaa, &
Om1ab,Om2ab,Om1aa,Om2aa,rho1ab,rho2ab,rho1aa,rho2aa,eT,eGT, & Om1ab,Om2ab,Om1aa,Om2aa,rho1ab,rho2ab,rho1aa,rho2aa,eT,eGT, &
dipole_int,OmBSE,XpY_BSE,XmY_BSE,TAab,TAaa) dipole_int,OmBSE,XpY_BSE,XmY_BSE,TAab,TAaa)
end if end if
end if
end if
end subroutine end subroutine

3
src/GT/GTpp_ppBSE.f90
View File

@ -1,4 +1,4 @@
subroutine GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nOOab,nVVab,nOOaa,nVVaa, & subroutine 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, & Om1ab,X1ab,Y1ab,Om2ab,X2ab,Y2ab,rho1ab,rho2ab,Om1aa,X1aa,Y1aa,Om2aa,X2aa,Y2aa,rho1aa,rho2aa, &
ERI,dipole_int,eT,eGT,EcBSE) ERI,dipole_int,eT,eGT,EcBSE)
@ -13,7 +13,6 @@ subroutine GTpp_ppBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,n
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet

3
src/GT/UG0T0.f90
View File

@ -1,4 +1,4 @@
subroutine UG0T0(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn, & subroutine UG0T0(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA, &
spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV, & spin_conserved,spin_flip,linearize,eta,regularize,nBas,nC,nO,nV, &
nR,nS,ENuc,EUHF,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb, & nR,nS,ENuc,EUHF,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,Vxc,eG0T0) dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,Vxc,eG0T0)
@ -18,7 +18,6 @@ subroutine UG0T0(doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn, &
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
logical,intent(in) :: linearize logical,intent(in) :: linearize

5
src/GT/evGTeh.f90
View File

@ -1,4 +1,4 @@
subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,evDyn,doppBSE, & subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE, &
singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF, & singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF, &
cHF,eHF,Vxc) cHF,eHF,Vxc)
@ -23,7 +23,6 @@ subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,d
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: doppBSE logical,intent(in) :: doppBSE
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
@ -237,7 +236,7 @@ subroutine evGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,d
! if(BSE) then ! if(BSE) then
! call Bethe_Salpeter(BSE2,TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eGW,eGW,EcBSE) ! call Bethe_Salpeter(BSE2,TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eGW,eGW,EcBSE)
! if(exchange_kernel) then ! if(exchange_kernel) then

5
src/GT/evGTpp.f90
View File

@ -1,4 +1,4 @@
subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA,dBSE,dTDA,evDyn, & 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_AO,ERI_MO,dipole_int,PHF,cHF,eHF,Vxc) singlet,triplet,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF,Vxc)
! Perform eigenvalue self-consistent calculation with a T-matrix self-energy (evGT) ! Perform eigenvalue self-consistent calculation with a T-matrix self-energy (evGT)
@ -19,7 +19,6 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -233,7 +232,7 @@ subroutine evGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
if(BSE) then if(BSE) then
call GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,nOOt,nVVt, & 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, & Om1s,X1s,Y1s,Om2s,X2s,Y2s,rho1s,rho2s,Om1t,X1t,Y1t,Om2t,X2t,Y2t,rho1t,rho2t, &
ERI_MO,dipole_int,eGT,eGT,EcBSE) ERI_MO,dipole_int,eGT,eGT,EcBSE)

3
src/GT/evUGT.f90
View File

@ -1,5 +1,5 @@
subroutine evUGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, & subroutine evUGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, &
TDA_T,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,& TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip,&
eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,ERI_AO,ERI_aaaa, & eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,ERI_AO,ERI_aaaa, &
ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF, & ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF, &
Vxc,eG0T0) Vxc,eG0T0)
@ -21,7 +21,6 @@ subroutine evUGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, &
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
double precision,intent(in) :: eta double precision,intent(in) :: eta

5
src/GT/qsGTeh.f90
View File

@ -1,5 +1,5 @@
subroutine qsGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,TDA_T,TDA, & subroutine qsGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,TDA_T,TDA, &
dBSE,dTDA,evDyn,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, & dBSE,dTDA,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
! Perform a quasiparticle self-consistent GTeh calculation ! Perform a quasiparticle self-consistent GTeh calculation
@ -21,7 +21,6 @@ subroutine qsGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -292,7 +291,7 @@ subroutine qsGTeh(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,
! if(BSE) then ! if(BSE) then
! call Bethe_Salpeter(BSE2,TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO, & ! call Bethe_Salpeter(BSE2,TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO, &
! eGT,eGT,EcBSE) ! eGT,eGT,EcBSE)
! if(exchange_kernel) then ! if(exchange_kernel) then

5
src/GT/qsGTpp.f90
View File

@ -1,5 +1,5 @@
subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA, & subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T,TDA, &
dBSE,dTDA,evDyn,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, & dBSE,dTDA,singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
! Perform a quasiparticle self-consistent GT calculation ! Perform a quasiparticle self-consistent GT calculation
@ -20,7 +20,6 @@ subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -355,7 +354,7 @@ subroutine qsGTpp(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,TDA_T
if(BSE) then if(BSE) then
call GTpp_phBSE(TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,nOOs,nVVs,nOOt,nVVt, & 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, & Om1s,X1s,Y1s,Om2s,X2s,Y2s,rho1s,rho2s,Om1t,X1t,Y1t,Om2t,X2t,Y2t,rho1t,rho2t, &
ERI_MO,dipole_int_MO,eGT,eGT,EcBSE) ERI_MO,dipole_int_MO,eGT,eGT,EcBSE)

3
src/GT/qsUGT.f90
View File

@ -1,5 +1,5 @@
subroutine qsUGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, & subroutine qsUGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, &
TDA_T,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,& TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip,&
eta,regularize,nBas,nC,nO,nV,nR,nS,nNuc,ZNuc,rNuc,ENuc,EUHF,S,X,T,V,Hc,ERI_AO,ERI_aaaa,& eta,regularize,nBas,nC,nO,nV,nR,nS,nNuc,ZNuc,rNuc,ENuc,EUHF,S,X,T,V,Hc,ERI_AO,ERI_aaaa,&
ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF) ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF)
@ -20,7 +20,6 @@ subroutine qsUGT(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE, &
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
double precision,intent(in) :: eta double precision,intent(in) :: eta

5
src/GW/G0W0.f90
View File

@ -1,4 +1,4 @@
subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,doppBSE, & subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO, & singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO, &
dipole_int,PHF,cHF,eHF,Vxc) dipole_int,PHF,cHF,eHF,Vxc)
@ -20,7 +20,6 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
logical,intent(in) :: linearize logical,intent(in) :: linearize
@ -182,7 +181,7 @@ subroutine G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dT
if(dophBSE) then if(dophBSE) then
call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE) call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eHF,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

58
src/GW/GW_phBSE.f90
View File

@ -1,4 +1,4 @@
subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eW,eGW,EcBSE) subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI,dipole_int,eW,eGW,EcBSE)
! Compute the Bethe-Salpeter excitation energies ! Compute the Bethe-Salpeter excitation energies
@ -12,7 +12,6 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
@ -53,8 +52,6 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
double precision,allocatable :: KB_sta(:,:) double precision,allocatable :: KB_sta(:,:)
double precision,allocatable :: W(:,:,:,:) double precision,allocatable :: W(:,:,:,:)
double precision,allocatable :: KA2_sta(:,:)
double precision,allocatable :: KB2_sta(:,:)
! Output variables ! Output variables
@ -63,8 +60,8 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
! Memory allocation ! Memory allocation
allocate(OmRPA(nS),XpY_RPA(nS,nS),XmY_RPA(nS,nS),rho_RPA(nBas,nBas,nS), & allocate(OmRPA(nS),XpY_RPA(nS,nS),XmY_RPA(nS,nS),rho_RPA(nBas,nBas,nS), &
A_sta(nS,nS),KA_sta(nS,nS),KA2_sta(nS,nS),OmBSE(nS),XpY_BSE(nS,nS),XmY_BSE(nS,nS)) A_sta(nS,nS),B_sta(nS,nS),KA_sta(nS,nS),KB_sta(nS,nS), &
allocate(B_sta(nS,nS),KB_sta(nS,nS),KB2_sta(nS,nS)) OmBSE(nS),XpY_BSE(nS,nS),XmY_BSE(nS,nS))
!--------------------------------- !---------------------------------
! Compute (singlet) RPA screening ! Compute (singlet) RPA screening
@ -96,9 +93,6 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI,A_sta) call phLR_A(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,eGW,ERI,A_sta)
if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,B_sta) if(.not.TDA) call phLR_B(ispin,dRPA,nBas,nC,nO,nV,nR,nS,1d0,ERI,B_sta)
A_sta(:,:) = A_sta(:,:) + KA_sta(:,:)
if(.not.TDA) B_sta(:,:) = B_sta(:,:) + KB_sta(:,:)
! Second-order BSE static kernel ! Second-order BSE static kernel
if(BSE2) then if(BSE2) then
@ -110,17 +104,17 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
write(*,*) write(*,*)
call static_kernel_W(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W) call static_kernel_W(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
call GW_phBSE2_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KA2_sta) call GW_phBSE2_static_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KA_sta)
if(.not.TDA) call GW_phBSE2_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KB2_sta) if(.not.TDA) call GW_phBSE2_static_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eW,W,KB_sta)
deallocate(W) deallocate(W)
A_sta(:,:) = A_sta(:,:) + KA2_sta(:,:)
if(.not.TDA) B_sta(:,:) = B_sta(:,:) + KB2_sta(:,:)
end if end if
A_sta(:,:) = A_sta(:,:) + KA_sta(:,:)
if(.not.TDA) B_sta(:,:) = B_sta(:,:) + KB_sta(:,:)
call phLR(TDA,nS,A_sta,B_sta,EcBSE(ispin),OmBSE,XpY_BSE,XmY_BSE) call phLR(TDA,nS,A_sta,B_sta,EcBSE(ispin),OmBSE,XpY_BSE,XmY_BSE)
call print_excitation('phBSE@GW ',ispin,nS,OmBSE) call print_excitation('phBSE@GW ',ispin,nS,OmBSE)
@ -130,21 +124,9 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
! Compute the dynamical screening at the BSE level ! Compute the dynamical screening at the BSE level
!------------------------------------------------- !-------------------------------------------------
if(dBSE) then if(dBSE) &
call GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
! Compute dynamic correction for BSE via perturbation theory (iterative or renormalized) OmBSE,XpY_BSE,XmY_BSE,KA_sta,KB_sta)
if(evDyn) then
call GW_phBSE_dynamic_perturbation_iterative(dTDA,eta,nBas,nC,nO,nV,nR,nS,eGW,dipole_int,OmRPA,rho_RPA, &
OmBSE,XpY_BSE,XmY_BSE)
else
call GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,dipole_int,OmRPA,rho_RPA, &
OmBSE,XpY_BSE,XmY_BSE,W,KA2_sta)
end if
end if
end if end if
@ -174,21 +156,9 @@ subroutine GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,n
! Compute the dynamical screening at the BSE level ! Compute the dynamical screening at the BSE level
!------------------------------------------------- !-------------------------------------------------
if(dBSE) then if(dBSE) &
call GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int,OmRPA,rho_RPA, &
! Compute dynamic correction for BSE via perturbation theory (iterative or renormalized) OmBSE,XpY_BSE,XmY_BSE,KA_sta,KB_sta)
if(evDyn) then
call GW_phBSE_dynamic_perturbation_iterative(dTDA,eta,nBas,nC,nO,nV,nR,nS,eGW,dipole_int,OmRPA,rho_RPA, &
OmBSE,XpY_BSE,XmY_BSE)
else
call GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,dipole_int,OmRPA,rho_RPA, &
OmBSE,XpY_BSE,XmY_BSE,W,KA2_sta)
end if
end if
end if end if

8
src/GW/GW_phBSE2_dynamic_kernel_A.f90
View File

@ -24,16 +24,10 @@ subroutine GW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE,A_dyn,
double precision,intent(out) :: A_dyn(nS,nS) double precision,intent(out) :: A_dyn(nS,nS)
double precision,intent(out) :: ZA_dyn(nS,nS) double precision,intent(out) :: ZA_dyn(nS,nS)
! Initialization
! A_dyn(:,:) = 0d0
! ZA_dyn(:,:) = 0d0
! Number of poles taken into account
! Build dynamic A matrix ! Build dynamic A matrix
jb = 0 jb = 0
!$omp parallel do default(private) shared(A_dyn,ZA_dyn,OmBSE,W,num,dem,eGW,nO,nBas,eta,nC,nR) !$omp parallel do default(private) shared(A_dyn,ZA_dyn,OmBSE,W,num,dem,eGW,nO,nBas,eta,nC,nR)
do j=nC+1,nO do j=nC+1,nO
do b=nO+1,nBas-nR do b=nO+1,nBas-nR

118
src/GW/GW_phBSE_dynamic_kernel.f90
View File

@ -1,118 +0,0 @@
subroutine GW_phBSE_dynamic_kernel(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rhO_RPA,OmBSE,Ap,Am,Bp,Bm)
! Compute the dynamic part of the Bethe-Salpeter equation matrices
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: nBas,nC,nO,nV,nR,nS
double precision,intent(in) :: eta
double precision,intent(in) :: lambda
double precision,intent(in) :: eGW(nBas)
double precision,intent(in) :: OmRPA(nS)
double precision,intent(in) :: rho_RPA(nBas,nBas,nS)
double precision,intent(in) :: OmBSE
! Local variables
integer :: maxS
double precision :: chi_A,chi_B
double precision :: chi_Ap,chi_Bp,eps_Ap,eps_Bp
double precision :: chi_Am,chi_Bm,eps_Am,eps_Bm
integer :: i,j,a,b,ia,jb,kc
! Output variables
double precision,intent(out) :: Ap(nS,nS)
double precision,intent(out) :: Am(nS,nS)
double precision,intent(out) :: Bp(nS,nS)
double precision,intent(out) :: Bm(nS,nS)
! Initialization
Ap(:,:) = 0d0
Am(:,:) = 0d0
Bp(:,:) = 0d0
Bm(:,:) = 0d0
! Number of poles taken into account
maxS = nS
! Build dynamic A matrix
ia = 0
do i=nC+1,nO
do a=nO+1,nBas-nR
ia = ia + 1
jb = 0
do j=nC+1,nO
do b=nO+1,nBas-nR
jb = jb + 1
chi_A = 0d0
chi_B = 0d0
do kc=1,maxS
chi_A = chi_A + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*OmRPA(kc)/(OmRPA(kc)**2 + eta**2)
chi_B = chi_B + rho_RPA(i,b,kc)*rho_RPA(a,j,kc)*OmRPA(kc)/(OmRPA(kc)**2 + eta**2)
enddo
Ap(ia,jb) = Ap(ia,jb) - 4d0*lambda*chi_A
Am(ia,jb) = Am(ia,jb) - 4d0*lambda*chi_A
Bp(ia,jb) = Bp(ia,jb) - 4d0*lambda*chi_B
Bm(ia,jb) = Bm(ia,jb) - 4d0*lambda*chi_B
chi_Ap = 0d0
chi_Am = 0d0
chi_Bp = 0d0
chi_Bm = 0d0
do kc=1,maxS
eps_Ap = + OmBSE - OmRPA(kc) - (eGW(a) - eGW(j))
chi_Ap = chi_Ap + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps_Ap/(eps_Ap**2 + eta**2)
eps_Ap = + OmBSE - OmRPA(kc) - (eGW(b) - eGW(i))
chi_Ap = chi_Ap + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps_Ap/(eps_Ap**2 + eta**2)
eps_Am = - OmBSE - OmRPA(kc) - (eGW(a) - eGW(j))
chi_Am = chi_Am + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps_Am/(eps_Am**2 + eta**2)
eps_Am = - OmBSE - OmRPA(kc) - (eGW(b) - eGW(i))
chi_Am = chi_Am + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps_Am/(eps_Am**2 + eta**2)
eps_Bp = - OmRPA(kc) - (eGW(a) - eGW(b))
chi_Bp = chi_Bp + rho_RPA(i,b,kc)*rho_RPA(a,j,kc)*eps_Bp/(eps_Bp**2 + eta**2)
eps_Bp = - OmRPA(kc) - (eGW(j) - eGW(i))
chi_Bp = chi_Bp + rho_RPA(i,b,kc)*rho_RPA(a,j,kc)*eps_Bp/(eps_Bp**2 + eta**2)
eps_Bm = - OmRPA(kc) - (eGW(a) - eGW(b))
chi_Bm = chi_Bm + rho_RPA(i,b,kc)*rho_RPA(a,j,kc)*eps_Bm/(eps_Bm**2 + eta**2)
eps_Bm = - OmRPA(kc) - (eGW(j) - eGW(i))
chi_Bm = chi_Bm + rho_RPA(i,b,kc)*rho_RPA(a,j,kc)*eps_Bm/(eps_Bm**2 + eta**2)
enddo
Ap(ia,jb) = Ap(ia,jb) - 2d0*lambda*chi_Ap
Am(ia,jb) = Am(ia,jb) - 2d0*lambda*chi_Am
Bp(ia,jb) = Bp(ia,jb) - 2d0*lambda*chi_Bp
Bm(ia,jb) = Bm(ia,jb) - 2d0*lambda*chi_Bm
enddo
enddo
enddo
enddo
end subroutine

13
src/GW/GW_phBSE_dynamic_kernel_A.f90
View File

@ -17,7 +17,6 @@ subroutine GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rh
! Local variables ! Local variables
integer :: maxS
double precision :: chi double precision :: chi
double precision :: eps double precision :: eps
integer :: i,j,a,b,ia,jb,kc integer :: i,j,a,b,ia,jb,kc
@ -32,14 +31,10 @@ subroutine GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rh
A_dyn(:,:) = 0d0 A_dyn(:,:) = 0d0
ZA_dyn(:,:) = 0d0 ZA_dyn(:,:) = 0d0
! Number of poles taken into account
maxS = nS
! Build dynamic A matrix ! Build dynamic A matrix
jb = 0 jb = 0
!$omp parallel do default(private) shared(A_dyn,ZA_dyn,OmRPA,OmBSE,eGW,rho_RPA,nO,nBas,maxS,chi,eps,eta,nC,nR,lambda) !$omp parallel do default(private) shared(A_dyn,ZA_dyn,OmRPA,OmBSE,eGW,rho_RPA,nO,nBas,nS,chi,eps,eta,nC,nR,lambda)
do j=nC+1,nO do j=nC+1,nO
do b=nO+1,nBas-nR do b=nO+1,nBas-nR
jb = (b-nO) + (j-1)*(nBas-nO) jb = (b-nO) + (j-1)*(nBas-nO)
@ -51,7 +46,7 @@ subroutine GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rh
chi = 0d0 chi = 0d0
do kc=1,maxS do kc=1,nS
eps = OmRPA(kc) eps = OmRPA(kc)
chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps/(eps**2 + eta**2) chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps/(eps**2 + eta**2)
enddo enddo
@ -60,7 +55,7 @@ subroutine GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rh
chi = 0d0 chi = 0d0
do kc=1,maxS do kc=1,nS
eps = + OmBSE - OmRPA(kc) - (eGW(a) - eGW(j)) eps = + OmBSE - OmRPA(kc) - (eGW(a) - eGW(j))
chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps/(eps**2 + eta**2) chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*eps/(eps**2 + eta**2)
@ -73,7 +68,7 @@ subroutine GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rh
A_dyn(ia,jb) = A_dyn(ia,jb) - 2d0*lambda*chi A_dyn(ia,jb) = A_dyn(ia,jb) - 2d0*lambda*chi
chi = 0d0 chi = 0d0
do kc=1,maxS do kc=1,nS
eps = + OmBSE - OmRPA(kc) - (eGW(a) - eGW(j)) eps = + OmBSE - OmRPA(kc) - (eGW(a) - eGW(j))
chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*(eps**2 - eta**2)/(eps**2 + eta**2)**2 chi = chi + rho_RPA(i,j,kc)*rho_RPA(a,b,kc)*(eps**2 - eta**2)/(eps**2 + eta**2)**2

60
src/GW/GW_phBSE_dynamic_kernel_B.f90 Normal file
View File

@ -0,0 +1,60 @@
subroutine GW_phBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,lambda,eGW,OmRPA,rhO,OmBSE,KB)
! Compute the dynamic part of the Bethe-Salpeter equation matrices
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: nBas,nC,nO,nV,nR,nS
double precision,intent(in) :: eta
double precision,intent(in) :: lambda
double precision,intent(in) :: eGW(nBas)
double precision,intent(in) :: OmRPA(nS)
double precision,intent(in) :: rho(nBas,nBas,nS)
double precision,intent(in) :: OmBSE
! Local variables
double precision :: chi,eps
integer :: i,j,a,b,ia,jb,kc
! Output variables
double precision,intent(out) :: KB(nS,nS)
! Initialization
KB(:,:) = 0d0
! Build dynamic A matrix
ia = 0
do i=nC+1,nO
do a=nO+1,nBas-nR
ia = ia + 1
jb = 0
do j=nC+1,nO
do b=nO+1,nBas-nR
jb = jb + 1
chi = 0d0
do kc=1,nS
eps = - OmRPA(kc) - (eGW(a) - eGW(b))
chi = chi + rho(i,b,kc)*rho(a,j,kc)*eps/(eps**2 + eta**2)
eps = - OmRPA(kc) - (eGW(j) - eGW(i))
chi = chi + rho(i,b,kc)*rho(a,j,kc)*eps/(eps**2 + eta**2)
enddo
KB(ia,jb) = KB(ia,jb) - 2d0*lambda*chi
enddo
enddo
enddo
enddo
end subroutine

65
src/GW/GW_phBSE_dynamic_perturbation.f90
View File

@ -1,5 +1,5 @@
subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,dipole_int, & subroutine GW_phBSE_dynamic_perturbation(dophBSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eGW,ERI,dipole_int, &
OmRPA,rho_RPA,OmBSE,XpY,XmY,W,A_stat) OmRPA,rho_RPA,OmBSE,XpY,XmY,KA_sta,KB_sta)
! Compute dynamical effects via perturbation theory for BSE ! Compute dynamical effects via perturbation theory for BSE
@ -8,7 +8,7 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
! Input variables ! Input variables
logical,intent(in) :: BSE2 logical,intent(in) :: dophBSE2
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
double precision,intent(in) :: eta double precision,intent(in) :: eta
integer,intent(in) :: nBas integer,intent(in) :: nBas
@ -18,6 +18,7 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
integer,intent(in) :: nR integer,intent(in) :: nR
integer,intent(in) :: nS integer,intent(in) :: nS
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: eW(nBas) double precision,intent(in) :: eW(nBas)
double precision,intent(in) :: eGW(nBas) double precision,intent(in) :: eGW(nBas)
double precision,intent(in) :: dipole_int(nBas,nBas,ncart) double precision,intent(in) :: dipole_int(nBas,nBas,ncart)
@ -26,8 +27,8 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
double precision,intent(in) :: OmBSE(nS) double precision,intent(in) :: OmBSE(nS)
double precision,intent(in) :: XpY(nS,nS) double precision,intent(in) :: XpY(nS,nS)
double precision,intent(in) :: XmY(nS,nS) double precision,intent(in) :: XmY(nS,nS)
double precision,intent(in) :: W(nBas,nBas,nBas,nBas) double precision,intent(in) :: KA_sta(nS,nS)
double precision,intent(in) :: A_stat(nS,nS) double precision,intent(in) :: KB_sta(nS,nS)
! Local variables ! Local variables
@ -41,23 +42,19 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
double precision,allocatable :: X(:) double precision,allocatable :: X(:)
double precision,allocatable :: Y(:) double precision,allocatable :: Y(:)
double precision,allocatable :: Ap_dyn(:,:) double precision,allocatable :: KAp_dyn(:,:)
double precision,allocatable :: KAm_dyn(:,:)
double precision,allocatable :: ZAp_dyn(:,:) double precision,allocatable :: ZAp_dyn(:,:)
double precision,allocatable :: Bp_dyn(:,:)
double precision,allocatable :: ZBp_dyn(:,:)
double precision,allocatable :: Am_dyn(:,:)
double precision,allocatable :: ZAm_dyn(:,:) double precision,allocatable :: ZAm_dyn(:,:)
double precision,allocatable :: Bm_dyn(:,:) double precision,allocatable :: KB_dyn(:,:)
double precision,allocatable :: ZBm_dyn(:,:)
double precision,allocatable :: W(:,:,:,:)
! Memory allocation ! Memory allocation
maxS = min(nS,maxS) allocate(OmDyn(maxS),ZDyn(maxS),X(nS),Y(nS),KAp_dyn(nS,nS),ZAp_dyn(nS,nS), &
allocate(OmDyn(maxS),ZDyn(maxS),X(nS),Y(nS),Ap_dyn(nS,nS),ZAp_dyn(nS,nS)) KAm_dyn(nS,nS),ZAm_dyn(nS,nS),KB_dyn(nS,nS))
allocate(Am_dyn(nS,nS),ZAm_dyn(nS,nS),Bp_dyn(nS,nS),ZBp_dyn(nS,nS),Bm_dyn(nS,nS),ZBm_dyn(nS,nS))
if(dTDA) then if(dTDA) then
write(*,*) write(*,*)
@ -65,6 +62,17 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
write(*,*) write(*,*)
end if end if
if(dophBSE2) then
write(*,*)
write(*,*) '*** Second-order BSE static kernel activated! ***'
write(*,*)
allocate(W(nBas,nBas,nBas,nBas))
call static_kernel_W(eta,nBas,nC,nO,nV,nR,nS,1d0,ERI,OmRPA,rho_RPA,W)
end if
gapGW = eGW(nO+1) - eGW(nO) gapGW = eGW(nO+1) - eGW(nO)
write(*,*) '---------------------------------------------------------------------------------------------------' write(*,*) '---------------------------------------------------------------------------------------------------'
@ -75,7 +83,7 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
write(*,'(2X,A5,1X,A20,1X,A20,1X,A20,1X,A20)') '#','Static (eV)','Dynamic (eV)','Correction (eV)','Renorm. (eV)' write(*,'(2X,A5,1X,A20,1X,A20,1X,A20,1X,A20)') '#','Static (eV)','Dynamic (eV)','Correction (eV)','Renorm. (eV)'
write(*,*) '---------------------------------------------------------------------------------------------------' write(*,*) '---------------------------------------------------------------------------------------------------'
do ia=1,maxS do ia=1,min(nS,maxS)
X(:) = 0.5d0*(XpY(ia,:) + XmY(ia,:)) X(:) = 0.5d0*(XpY(ia,:) + XmY(ia,:))
Y(:) = 0.5d0*(XpY(ia,:) - XmY(ia,:)) Y(:) = 0.5d0*(XpY(ia,:) - XmY(ia,:))
@ -86,32 +94,29 @@ subroutine GW_phBSE_dynamic_perturbation(BSE2,dTDA,eta,nBas,nC,nO,nV,nR,nS,eW,eG
! Resonant part of the BSE correction for dynamical TDA ! Resonant part of the BSE correction for dynamical TDA
call GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,OmBSE(ia),Ap_dyn,ZAp_dyn) call GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,+OmBSE(ia),KAp_dyn,ZAp_dyn)
if(BSE2) call GW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE(ia),Ap_dyn,ZAp_dyn,W) if(dophBSE2) call GW_phBSE2_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,eGW,W,OmBSE(ia),KAp_dyn,ZAp_dyn)
ZDyn(ia) = dot_product(X,matmul(ZAp_dyn,X)) ZDyn(ia) = dot_product(X,matmul(ZAp_dyn,X))
OmDyn(ia) = dot_product(X,matmul( Ap_dyn - A_stat,X)) OmDyn(ia) = dot_product(X,matmul(KAp_dyn - KA_sta,X))
else else
! Resonant and anti-resonant part of the BSE correction ! Resonant and anti-resonant part of the BSE correction
call GW_phBSE_dynamic_kernel(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,OmBSE(ia),Ap_dyn,Am_dyn,Bp_dyn,Bm_dyn) call GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,-OmBSE(ia),KAm_dyn,ZAm_dyn)
call GW_phBSE_dynamic_kernel_B(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,KB_dyn)
! Renormalization factor of the resonant and anti-resonant parts ! Renormalization factor of the resonant and anti-resonant parts
call GW_phBSE_dynamic_kernel_Z(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,OmBSE(ia),ZAp_dyn,ZAm_dyn,ZBp_dyn,ZBm_dyn)
ZDyn(ia) = dot_product(X,matmul(ZAp_dyn,X)) & ZDyn(ia) = dot_product(X,matmul(ZAp_dyn,X)) &
- dot_product(Y,matmul(ZAm_dyn,Y)) & + dot_product(Y,matmul(ZAm_dyn,Y))
+ dot_product(X,matmul(ZBp_dyn,Y)) &
- dot_product(Y,matmul(ZBm_dyn,X))
OmDyn(ia) = dot_product(X,matmul(Ap_dyn,X)) & OmDyn(ia) = dot_product(X,matmul(KAp_dyn - KA_sta,X)) &
- dot_product(Y,matmul(Am_dyn,Y)) & - dot_product(Y,matmul(KAm_dyn - KA_sta,Y)) &
+ dot_product(X,matmul(Bp_dyn,Y)) & + dot_product(X,matmul(KB_dyn - KB_sta,Y)) &
- dot_product(Y,matmul(Bm_dyn,X)) - dot_product(Y,matmul(KB_dyn - KB_sta,X))
end if end if

145
src/GW/GW_phBSE_dynamic_perturbation_iterative.f90
View File

@ -1,145 +0,0 @@
subroutine GW_phBSE_dynamic_perturbation_iterative(dTDA,eta,nBas,nC,nO,nV,nR,nS,eGW,dipole_int,OmRPA,rho_RPA,OmBSE,XpY,XmY)
! Compute self-consistently the dynamical effects via perturbation theory for BSE
implicit none
include 'parameters.h'
! Input variables
logical,intent(in) :: dTDA
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) :: nS
double precision,intent(in) :: eGW(nBas)
double precision,intent(in) :: dipole_int(nBas,nBas,ncart)
double precision,intent(in) :: OmRPA(nS)
double precision,intent(in) :: rho_RPA(nBas,nBas,nS)
double precision,intent(in) :: OmBSE(nS)
double precision,intent(in) :: XpY(nS,nS)
double precision,intent(in) :: XmY(nS,nS)
! Local variables
integer :: ia
integer :: maxS = 10
double precision :: gapGW
integer :: nSCF
integer :: maxSCF = 10
double precision :: Conv
double precision :: thresh = 1d-3
double precision,allocatable :: OmDyn(:)
double precision,allocatable :: OmOld(:)
double precision,allocatable :: X(:)
double precision,allocatable :: Y(:)
double precision,allocatable :: Ap_dyn(:,:)
double precision,allocatable :: ZAp_dyn(:,:)
double precision,allocatable :: Bp_dyn(:,:)
double precision,allocatable :: Am_dyn(:,:)
double precision,allocatable :: ZAm_dyn(:,:)
double precision,allocatable :: Bm_dyn(:,:)
! Memory allocation
maxS = min(nS,maxS)
allocate(OmDyn(maxS),OmOld(maxS),X(nS),Y(nS),Ap_dyn(nS,nS),ZAp_dyn(nS,nS))
allocate(Am_dyn(nS,nS),ZAm_dyn(nS,nS),Bp_dyn(nS,nS),Bm_dyn(nS,nS))
if(dTDA) then
write(*,*)
write(*,*) '*** dynamical TDA activated ***'
write(*,*)
end if
gapGW = eGW(nO+1) - eGW(nO)
Conv = 1d0
nSCF = 0
OmOld(1:maxS) = OmBSE(1:maxS)
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*) ' First-order dynamical correction to static Bethe-Salpeter excitation energies '
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(A57,F10.6,A3)') ' BSE neutral excitation must be lower than the GW gap = ',gapGW*HaToeV,' eV'
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*)
do while(Conv > thresh .and. nSCF < maxSCF)
nSCF = nSCF + 1
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A15,I3)') 'Iteration n.',nSCF
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A5,1X,A20,1X,A20,1X,A20,A20)') '#','Static (eV)','Dynamic (eV)','Correction (eV)','Convergence (eV)'
write(*,*) '---------------------------------------------------------------------------------------------------'
do ia=1,maxS
X(:) = 0.5d0*(XpY(ia,:) + XmY(ia,:))
Y(:) = 0.5d0*(XpY(ia,:) - XmY(ia,:))
! First-order correction
if(dTDA) then
! Resonant part of the BSE correction
call GW_phBSE_dynamic_kernel_A(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,OmOld(ia),Ap_dyn,ZAp_dyn)
OmDyn(ia) = dot_product(X(:),matmul(Ap_dyn(:,:),X(:)))
else
! Anti-resonant part of the BSE correction
call GW_phBSE_dynamic_kernel(eta,nBas,nC,nO,nV,nR,nS,1d0,eGW,OmRPA,rho_RPA,OmOld(ia),Ap_dyn,Am_dyn,Bp_dyn,Bm_dyn)
OmDyn(ia) = dot_product(X(:),matmul(Ap_dyn(:,:),X(:))) &
- dot_product(Y(:),matmul(Am_dyn(:,:),Y(:))) &
+ dot_product(X(:),matmul(Bp_dyn(:,:),Y(:))) &
- dot_product(Y(:),matmul(Bm_dyn(:,:),X(:)))
end if
write(*,'(2X,I5,5X,F15.6,5X,F15.6,5X,F15.6,5X,F15.6)') &
ia,OmBSE(ia)*HaToeV,(OmBSE(ia)+OmDyn(ia))*HaToeV,OmDyn(ia)*HaToeV,(OmBSE(ia) + OmDyn(ia) - OmOld(ia))*HaToeV
end do
Conv = maxval(abs(OmBSE(1:maxS) + OmDyn(:) - OmOld(:)))*HaToeV
OmOld(:) = OmBSE(1:maxS) + OmDyn(:)
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,'(2X,A20,1X,F10.6)') ' Convergence = ',Conv
write(*,*) '---------------------------------------------------------------------------------------------------'
write(*,*)
end do
! Did it actually converge?
if(nSCF == maxSCF) then
write(*,*)
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)' Convergence failed '
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
write(*,*)
endif
end subroutine

5
src/GW/SRG_qsGW.f90
View File

@ -1,5 +1,5 @@
subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,TDA_W,TDA, & subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE2,TDA_W,TDA, &
dBSE,dTDA,evDyn,singlet,triplet,eta,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, & dBSE,dTDA,singlet,triplet,eta,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF, &
S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
! Perform a quasiparticle self-consistent GW calculation ! Perform a quasiparticle self-consistent GW calculation
@ -21,7 +21,6 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -308,7 +307,7 @@ subroutine SRG_qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,BSE,BSE
if(BSE) then if(BSE) then
call GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO,eGW,eGW,EcBSE) call GW_phBSE(BSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO,eGW,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

5
src/GW/UG0W0.f90
View File

@ -1,4 +1,4 @@
subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip, & subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb, & linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,EUHF,S,ERI,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,Vxc,eGW) dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,Vxc,eGW)
@ -19,7 +19,6 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
logical,intent(in) :: linearize logical,intent(in) :: linearize
@ -193,7 +192,7 @@ subroutine UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA,dBSE,dTDA,ev
if(BSE) then if(BSE) then
call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS,S, & call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS,S, &
ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eHF,eGW,EcBSE) ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eHF,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

5
src/GW/evGW.f90
View File

@ -1,4 +1,4 @@
subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,doppBSE, & subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF, & singlet,triplet,linearize,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc,ERHF,ERI_AO,ERI_MO,dipole_int,PHF, &
cHF,eHF,Vxc) cHF,eHF,Vxc)
@ -23,7 +23,6 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: doppBSE logical,intent(in) :: doppBSE
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
@ -234,7 +233,7 @@ subroutine evGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop
if(dophBSE) then if(dophBSE) then
call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eGW,eGW,EcBSE) call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int,eGW,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

5
src/GW/evUGW.f90
View File

@ -1,5 +1,5 @@
subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA, & subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA, &
dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc, & dBSE,dTDA,spin_conserved,spin_flip,eta,regularize,nBas,nC,nO,nV,nR,nS,ENuc, &
EUHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,Vxc,eG0W0) EUHF,S,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,PHF,cHF,eHF,Vxc,eG0W0)
! Perform self-consistent eigenvalue-only GW calculation ! Perform self-consistent eigenvalue-only GW calculation
@ -23,7 +23,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -241,7 +240,7 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE
if(BSE) then if(BSE) then
call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, & call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, &
S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eGW,eGW,EcBSE) S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eGW,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

5
src/GW/qsGW.f90
View File

@ -1,4 +1,4 @@
subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,doppBSE, & subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE, &
singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF,S,X,T,V,Hc,ERI_AO, & singlet,triplet,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,ERHF,S,X,T,V,Hc,ERI_AO, &
ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF) ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,eHF)
@ -21,7 +21,6 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: doppBSE logical,intent(in) :: doppBSE
logical,intent(in) :: singlet logical,intent(in) :: singlet
logical,intent(in) :: triplet logical,intent(in) :: triplet
@ -292,7 +291,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,dophBSE,dop
if(dophBSE) then if(dophBSE) then
call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO,eGW,eGW,EcBSE) call GW_phBSE(dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta,nBas,nC,nO,nV,nR,nS,ERI_MO,dipole_int_MO,eGW,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

5
src/GW/qsUGW.f90
View File

@ -1,5 +1,5 @@
subroutine qsUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA, & subroutine qsUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE,TDA_W,TDA, &
dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO, & dBSE,dTDA,spin_conserved,spin_flip,eta,regularize,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO, &
nV,nR,nS,EUHF,S,X,T,V,Hc,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa, & nV,nR,nS,EUHF,S,X,T,V,Hc,ERI_AO,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_AO,dipole_int_aa, &
dipole_int_bb,PHF,cHF,eHF) dipole_int_bb,PHF,cHF,eHF)
@ -22,7 +22,6 @@ subroutine qsUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip
double precision,intent(in) :: eta double precision,intent(in) :: eta
@ -371,7 +370,7 @@ subroutine qsUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE
if(BSE) then if(BSE) then
call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, & call unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta,nBas,nC,nO,nV,nR,nS, &
S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,c,eGW,eGW,EcBSE) S,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,c,eGW,eGW,EcBSE)
if(exchange_kernel) then if(exchange_kernel) then

3
src/GW/unrestricted_Bethe_Salpeter.f90
View File

@ -1,4 +1,4 @@
subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta, & subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta, &
nBas,nC,nO,nV,nR,nS,S,ERI_aaaa,ERI_aabb,ERI_bbbb, & nBas,nC,nO,nV,nR,nS,S,ERI_aaaa,ERI_aabb,ERI_bbbb, &
dipole_int_aa,dipole_int_bb,cW,eW,eGW,EcBSE) dipole_int_aa,dipole_int_bb,cW,eW,eGW,EcBSE)
@ -13,7 +13,6 @@ subroutine unrestricted_Bethe_Salpeter(TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,
logical,intent(in) :: TDA logical,intent(in) :: TDA
logical,intent(in) :: dBSE logical,intent(in) :: dBSE
logical,intent(in) :: dTDA logical,intent(in) :: dTDA
logical,intent(in) :: evDyn
logical,intent(in) :: spin_conserved logical,intent(in) :: spin_conserved
logical,intent(in) :: spin_flip logical,intent(in) :: spin_flip

48
src/QuAcK/QuAcK.f90
View File

@ -101,7 +101,7 @@ program QuAcK
logical :: DIIS_GT,TDA_T,linGT,regGT logical :: DIIS_GT,TDA_T,linGT,regGT
double precision :: eta_GT double precision :: eta_GT
logical :: dophBSE,dophBSE2,doppBSE,dBSE,dTDA,evDyn logical :: dophBSE,dophBSE2,doppBSE,dBSE,dTDA
! Hello World ! Hello World
@ -145,7 +145,7 @@ program QuAcK
maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW,regGW,COHSEX,TDA_W, & maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW,regGW,COHSEX,TDA_W, &
maxSCF_GT,thresh_GT,DIIS_GT,n_diis_GT,linGT,eta_GT,regGT,TDA_T, & maxSCF_GT,thresh_GT,DIIS_GT,n_diis_GT,linGT,eta_GT,regGT,TDA_T, &
doACFDT,exchange_kernel,doXBS, & doACFDT,exchange_kernel,doXBS, &
dophBSE,dBSE,dTDA,evDyn,doppBSE,dophBSE2) dophBSE,dBSE,dTDA,doppBSE,dophBSE2)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Read input information ! Read input information
@ -667,13 +667,13 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call UG0F2(dophBSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,linGF,eta_GF,regGF, & call UG0F2(dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,linGF,eta_GF,regGF, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, & nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, &
dipole_int_aa,dipole_int_bb,epsHF) dipole_int_aa,dipole_int_bb,epsHF)
else else
call G0F2(dophBSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,linGF,eta_GF,regGF, & call G0F2(dophBSE,TDA,dBSE,dTDA,singlet,triplet,linGF,eta_GF,regGF, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_MO,dipole_int_MO,epsHF) nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_MO,dipole_int_MO,epsHF)
end if end if
@ -696,13 +696,13 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call evUGF2(maxSCF_GF,thresh_GF,n_diis_GF,dophBSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip, & call evUGF2(maxSCF_GF,thresh_GF,n_diis_GF,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
eta_GF,regGF,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, & eta_GF,regGF,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, &
dipole_int_aa,dipole_int_bb,cHF,epsHF) dipole_int_aa,dipole_int_bb,cHF,epsHF)
else else
call evGF2(dophBSE,TDA,dBSE,dTDA,evDyn,maxSCF_GF,thresh_GF,n_diis_GF, & call evGF2(dophBSE,TDA,dBSE,dTDA,maxSCF_GF,thresh_GF,n_diis_GF, &
singlet,triplet,linGF,eta_GF,regGF,nBas,nC,nO,nV,nR,nS,ENuc,EHF, & singlet,triplet,linGF,eta_GF,regGF,nBas,nC,nO,nV,nR,nS,ENuc,EHF, &
ERI_MO,dipole_int_MO,epsHF) ERI_MO,dipole_int_MO,epsHF)
@ -726,13 +726,13 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call qsUGF2(maxSCF_GF,thresh_GF,n_diis_GF,dophBSE,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta_GF,regGF, & call qsUGF2(maxSCF_GF,thresh_GF,n_diis_GF,dophBSE,TDA,dBSE,dTDA,spin_conserved,spin_flip,eta_GF,regGF, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO, & nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO, &
ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF) ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
else else
call qsGF2(maxSCF_GF,thresh_GF,n_diis_GF,dophBSE,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta_GF,regGF,nNuc,ZNuc,rNuc,ENuc, & call qsGF2(maxSCF_GF,thresh_GF,n_diis_GF,dophBSE,TDA,dBSE,dTDA,singlet,triplet,eta_GF,regGF,nNuc,ZNuc,rNuc,ENuc, &
nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF) nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
end if end if
@ -806,12 +806,12 @@ program QuAcK
call cpu_time(start_GW) call cpu_time(start_GW)
if(unrestricted) then if(unrestricted) then
call UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,dophBSE,TDA_W,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip, & call UG0W0(doACFDT,exchange_kernel,doXBS,COHSEX,dophBSE,TDA_W,TDA,dBSE,dTDA,spin_conserved,spin_flip, &
linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, & linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, &
dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF,Vxc) dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF,Vxc)
else else
call G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,doppBSE,singlet,triplet, & call G0W0(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet, &
linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc) linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc)
end if end if
@ -833,14 +833,14 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call evUGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,dophBSE,TDA_W,TDA, & call evUGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,dophBSE,TDA_W,TDA, &
dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc, & dBSE,dTDA,spin_conserved,spin_flip,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc, &
EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa,dipole_int_bb, & EHF,S,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa,dipole_int_bb, &
PHF,cHF,epsHF,Vxc) PHF,cHF,epsHF,Vxc)
else else
call evGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS, & call evGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,doppBSE,singlet,triplet,linGW,eta_GW,regGW, & dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet,linGW,eta_GW,regGW, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc) nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc)
end if end if
call cpu_time(end_GW) call cpu_time(end_GW)
@ -862,14 +862,14 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call qsUGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,dophBSE,TDA_W,TDA, & call qsUGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,dophBSE,TDA_W,TDA, &
dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta_GW,regGW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO, & dBSE,dTDA,spin_conserved,spin_flip,eta_GW,regGW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO, &
nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_AO, & nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_AO, &
dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF) dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
else else
call qsGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS, & call qsGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta_GW,regGW,nNuc,ZNuc,rNuc,ENuc, & dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,singlet,triplet,eta_GW,regGW,nNuc,ZNuc,rNuc,ENuc, &
nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF) nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
end if end if
@ -896,7 +896,7 @@ program QuAcK
else else
call SRG_qsGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn, & call SRG_qsGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA, &
singlet,triplet,eta_GW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO, & singlet,triplet,eta_GW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO, &
dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF) dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
@ -952,14 +952,14 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call UG0T0(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,evDyn, & call UG0T0(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA, &
spin_conserved,spin_flip,linGT,eta_GT,regGT,nBas,nC,nO,nV, & spin_conserved,spin_flip,linGT,eta_GT,regGT,nBas,nC,nO,nV, &
nR,nS,ENuc,EHF,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, & nR,nS,ENuc,EHF,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb, &
dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF,Vxc) dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF,Vxc)
else else
call G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,evDyn,doppBSE,singlet,triplet, & call G0T0pp(doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T,TDA,dBSE,dTDA,doppBSE,singlet,triplet, &
linGT,eta_GT,regGT,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc) linGT,eta_GT,regGT,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc)
end if end if
@ -983,7 +983,7 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call evUGT(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, & call evUGT(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, &
dophBSE,TDA_T,TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,& dophBSE,TDA_T,TDA,dBSE,dTDA,spin_conserved,spin_flip,&
eta_GT,regGT,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO, & eta_GT,regGT,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO, &
ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa, & ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int_aa, &
dipole_int_bb,PHF,cHF,epsHF,Vxc) dipole_int_bb,PHF,cHF,epsHF,Vxc)
@ -991,7 +991,7 @@ program QuAcK
else else
call evGTpp(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, & call evGTpp(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, &
dophBSE,TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta_GT,regGT, & dophBSE,TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta_GT,regGT, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO, & nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO, &
PHF,cHF,epsHF,Vxc) PHF,cHF,epsHF,Vxc)
@ -1016,13 +1016,13 @@ program QuAcK
if(unrestricted) then if(unrestricted) then
call qsUGT(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T, & call qsUGT(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS,dophBSE,TDA_T, &
TDA,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta_GT,regGT,nBas,nC,nO,nV,& TDA,dBSE,dTDA,spin_conserved,spin_flip,eta_GT,regGT,nBas,nC,nO,nV,&
nR,nS,nNuc,ZNuc,rNuc,ENuc,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,& nR,nS,nNuc,ZNuc,rNuc,ENuc,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,&
ERI_MO_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF) ERI_MO_bbbb,dipole_int_AO,dipole_int_aa,dipole_int_bb,PHF,cHF,epsHF)
else else
call qsGTpp(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, & call qsGTpp(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, &
dophBSE,TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta_GT,regGT, & dophBSE,TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta_GT,regGT, &
nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc, & nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc, &
ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF) ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
@ -1050,7 +1050,7 @@ program QuAcK
else else
call G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,evDyn,doppBSE,singlet,triplet, & call G0T0eh(doACFDT,exchange_kernel,doXBS,dophBSE,dophBSE2,TDA_W,TDA,dBSE,dTDA,doppBSE,singlet,triplet, &
linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc) linGW,eta_GW,regGW,nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc)
end if end if
@ -1075,7 +1075,7 @@ program QuAcK
else else
call evGTeh(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, & call evGTeh(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,evDyn,doppBSE,singlet,triplet,linGT,eta_GT,regGT, & dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,doppBSE,singlet,triplet,linGT,eta_GT,regGT, &
nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc) nBas,nC,nO,nV,nR,nS,ENuc,EHF,ERI_AO,ERI_MO,dipole_int_MO,PHF,cHF,epsHF,Vxc)
end if end if
call cpu_time(end_GT) call cpu_time(end_GT)
@ -1099,7 +1099,7 @@ program QuAcK
else else
call qsGTeh(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, & call qsGTeh(maxSCF_GT,thresh_GT,n_diis_GT,doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,evDyn,singlet,triplet,eta_GT,regGT,nNuc,ZNuc,rNuc,ENuc, & dophBSE,dophBSE2,TDA_T,TDA,dBSE,dTDA,singlet,triplet,eta_GT,regGT,nNuc,ZNuc,rNuc,ENuc, &
nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF) nBas,nC,nO,nV,nR,nS,EHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int_AO,dipole_int_MO,PHF,cHF,epsHF)
end if end if

7
src/QuAcK/read_options.f90
View File

@ -6,7 +6,7 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t
maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW,regGW,COHSEX,TDA_W, & maxSCF_GW,thresh_GW,DIIS_GW,n_diis_GW,linGW,eta_GW,regGW,COHSEX,TDA_W, &
maxSCF_GT,thresh_GT,DIIS_GT,n_diis_GT,linGT,eta_GT,regGT,TDA_T, & maxSCF_GT,thresh_GT,DIIS_GT,n_diis_GT,linGT,eta_GT,regGT,TDA_T, &
doACFDT,exchange_kernel,doXBS, & doACFDT,exchange_kernel,doXBS, &
dophBSE,dophBSE2,doppBSE,dBSE,dTDA,evDyn) dophBSE,dophBSE2,doppBSE,dBSE,dTDA)
! Read desired methods ! Read desired methods
@ -74,7 +74,6 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t
logical,intent(out) :: doppBSE logical,intent(out) :: doppBSE
logical,intent(out) :: dBSE logical,intent(out) :: dBSE
logical,intent(out) :: dTDA logical,intent(out) :: dTDA
logical,intent(out) :: evDyn
! Local variables ! Local variables
@ -227,17 +226,15 @@ subroutine read_options(maxSCF_HF,thresh_HF,DIIS_HF,n_diis_HF,guess_type,ortho_t
doppBSE = .false. doppBSE = .false.
dBSE = .false. dBSE = .false.
dTDA = .true. dTDA = .true.
evDyn = .false.
read(1,*) read(1,*)
read(1,*) answer1,answer2,answer3,answer4,answer5,answer6 read(1,*) answer1,answer2,answer3,answer4,answer5
if(answer1 == 'T') dophBSE = .true. if(answer1 == 'T') dophBSE = .true.
if(answer2 == 'T') dophBSE2 = .true. if(answer2 == 'T') dophBSE2 = .true.
if(answer3 == 'T') doppBSE = .true. if(answer3 == 'T') doppBSE = .true.
if(answer4 == 'T') dBSE = .true. if(answer4 == 'T') dBSE = .true.
if(answer5 == 'F') dTDA = .false. if(answer5 == 'F') dTDA = .false.
if(answer6 == 'T') evDyn = .true.
! Close file with options ! Close file with options