mirror of https://github.com/pfloos/quack
working on qsUGW
This commit is contained in:
parent
efaf67b649
commit
2e6662056c
|
|
@ -1,9 +1,9 @@
|
|||
# RHF UHF KS MOM
|
||||
T F F F
|
||||
F T F F
|
||||
# MP2* MP3 MP2-F12
|
||||
F F F
|
||||
# CCD CCSD CCSD(T)
|
||||
T F F
|
||||
F F F
|
||||
# drCCD rCCD lCCD pCCD
|
||||
F F F F
|
||||
# CIS* CIS(D) CID CISD
|
||||
|
|
@ -12,8 +12,8 @@
|
|||
F F F
|
||||
# G0F2 evGF2 G0F3 evGF3
|
||||
F F F F
|
||||
# G0W0* evGW* qsGW
|
||||
F F F
|
||||
# G0W0* evGW* qsGW*
|
||||
F F T
|
||||
# G0T0 evGT qsGT
|
||||
F F F
|
||||
# MCMP2
|
||||
|
|
|
|||
|
|
@ -5,11 +5,11 @@
|
|||
# CC: maxSCF thresh DIIS n_diis
|
||||
64 0.0000001 T 5
|
||||
# spin: TDA singlet triplet spin_conserved spin_flip
|
||||
F T F T T
|
||||
F T T T T
|
||||
# GF: maxSCF thresh DIIS n_diis lin eta renorm
|
||||
256 0.00001 T 5 T 0.0 3
|
||||
# GW/GT: maxSCF thresh DIIS n_diis lin eta COHSEX SOSEX TDA_W G0W GW0
|
||||
256 0.00001 T 5 T 0.0 F F F F F
|
||||
128 0.00001 T 5 T 0.00367493 F F T F F
|
||||
# ACFDT: AC Kx XBS
|
||||
F F T
|
||||
# BSE: BSE dBSE dTDA evDyn
|
||||
|
|
|
|||
|
|
@ -1,4 +1,4 @@
|
|||
2
|
||||
\ce{CH},^2\Pi,CC3,aug-cc-pVTZ
|
||||
|
||||
C 0.00000000 0.00000000 -0.08596945
|
||||
H 0.00000000 0.00000000 1.02362355
|
||||
H 0.00000000 0.00000000 1.02362355
|
||||
|
|
|
|||
|
|
@ -1,4 +1,4 @@
|
|||
2
|
||||
|
||||
H 0.0 0.0 0.0
|
||||
H 0.0 0.0 5.0
|
||||
H 0.0 0.0 1.4
|
||||
|
|
|
|||
|
|
@ -174,8 +174,8 @@ subroutine UHF(maxSCF,thresh,max_diis,guess_type,nNuc,ZNuc,rNuc,ENuc,nBas,nO,S,T
|
|||
|
||||
n_diis = min(n_diis+1,max_diis)
|
||||
do ispin=1,nspin
|
||||
if(nO(ispin) > 1) call DIIS_extrapolation(rcond(ispin),nBasSq,nBasSq,n_diis,err_diis(:,:,ispin),F_diis(:,:,ispin), &
|
||||
err(:,:,ispin),F(:,:,ispin))
|
||||
if(nO(ispin) > 1) call DIIS_extrapolation(rcond(ispin),nBasSq,nBasSq,n_diis,err_diis(:,1:n_diis,ispin), &
|
||||
F_diis(:,1:n_diis,ispin),err(:,:,ispin),F(:,:,ispin))
|
||||
end do
|
||||
|
||||
! Reset DIIS if required
|
||||
|
|
|
|||
|
|
@ -153,8 +153,6 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE
|
|||
|
||||
call unrestricted_linear_response(ispin,.true.,TDA_W,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, &
|
||||
eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
|
||||
|
||||
|
||||
endif
|
||||
|
||||
!----------------------!
|
||||
|
|
@ -203,8 +201,8 @@ subroutine evUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,BSE
|
|||
|
||||
n_diis = min(n_diis+1,max_diis)
|
||||
do is=1,nspin
|
||||
call DIIS_extrapolation(rcond(ispin),nBas,nBas,n_diis,error_diis(:,:,is), &
|
||||
e_diis(:,:,is),eGW(:,is)-eOld(:,is),eGW(:,is))
|
||||
call DIIS_extrapolation(rcond(ispin),nBas,nBas,n_diis,error_diis(:,1:n_diis,is), &
|
||||
e_diis(:,1:n_diis,is),eGW(:,is)-eOld(:,is),eGW(:,is))
|
||||
end do
|
||||
|
||||
! Reset DIIS if required
|
||||
|
|
|
|||
|
|
@ -1,4 +1,4 @@
|
|||
subroutine print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigC,Z,EcRPA,EcGM,EqsGW)
|
||||
subroutine print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigC,Z,EcRPA,EqsGW)
|
||||
|
||||
|
||||
! Print one-electron energies and other stuff for qsGW
|
||||
|
|
@ -9,7 +9,7 @@ subroutine print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigC,Z
|
|||
! Input variables
|
||||
|
||||
integer,intent(in) :: nBas,nO,nSCF
|
||||
double precision,intent(in) :: ENuc,EcRPA,EcGM,Conv,thresh
|
||||
double precision,intent(in) :: ENuc,EcRPA,Conv,thresh
|
||||
double precision,intent(in) :: eHF(nBas),eGW(nBas),c(nBas),P(nBas,nBas)
|
||||
double precision,intent(in) :: T(nBas,nBas),V(nBas,nBas)
|
||||
double precision,intent(in) :: J(nBas,nBas),K(nBas,nBas),F(nBas,nBas)
|
||||
|
|
@ -35,8 +35,7 @@ subroutine print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigC,Z
|
|||
EV = trace_matrix(nBas,matmul(P,V))
|
||||
EJ = 0.5d0*trace_matrix(nBas,matmul(P,J))
|
||||
Ex = 0.25d0*trace_matrix(nBas,matmul(P,K))
|
||||
! Ec = 0.5d0*trace_matrix(nBas,matmul(P,SigC))
|
||||
Ec = 0d0
|
||||
Ec = 0.5d0*trace_matrix(nBas,matmul(P,SigC))
|
||||
EqsGW = ET + EV + EJ + Ex + Ec
|
||||
|
||||
! Dump results
|
||||
|
|
@ -57,7 +56,6 @@ subroutine print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigC,Z
|
|||
'|',x,'|',eHF(x)*HaToeV,'|',(eGW(x)-eHF(x))*HaToeV,'|',Z(x),'|',eGW(x)*HaToeV,'|'
|
||||
enddo
|
||||
|
||||
|
||||
write(*,*)'-------------------------------------------------------------------------------'
|
||||
write(*,'(2X,A10,I3)') 'Iteration ',nSCF
|
||||
write(*,'(2X,A19,F15.5)')'max(|FPS - SPF|) = ',Conv
|
||||
|
|
@ -67,11 +65,9 @@ subroutine print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigC,Z
|
|||
write(*,'(2X,A30,F15.6)') 'qsGW HOMO-LUMO gap (eV):',Gap*HaToeV
|
||||
write(*,*)'-------------------------------------------'
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW total energy =',EqsGW + ENuc
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW GM total energy =',EqsGW + ENuc + EcGM
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW exchange energy =',Ex
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW correlation energy =',Ec
|
||||
write(*,'(2X,A30,F15.6)') 'RPA@qsGW correlation energy =',EcRPA
|
||||
write(*,'(2X,A30,F15.6)') 'GM@qsGW correlation energy =',EcGM
|
||||
write(*,*)'-------------------------------------------'
|
||||
write(*,*)
|
||||
|
||||
|
|
|
|||
|
|
@ -0,0 +1,204 @@
|
|||
subroutine print_qsUGW(nBas,nO,Ov,nSCF,Conv,thresh,eGW,cGW,PGW,T,V,J,K,ENuc,EHF,SigC,Z,EcRPA,dipole)
|
||||
|
||||
! Print one-electron energies and other stuff for qsUGW
|
||||
|
||||
implicit none
|
||||
include 'parameters.h'
|
||||
|
||||
! Input variables
|
||||
|
||||
integer,intent(in) :: nBas
|
||||
integer,intent(in) :: nO(nspin)
|
||||
double precision,intent(in) :: Ov(nBas,nBas)
|
||||
integer,intent(in) :: nSCF
|
||||
double precision,intent(in) :: ENuc
|
||||
double precision,intent(in) :: EHF
|
||||
double precision,intent(in) :: EcRPA
|
||||
double precision,intent(in) :: Conv
|
||||
double precision,intent(in) :: thresh
|
||||
double precision,intent(in) :: eGW(nBas,nspin)
|
||||
double precision,intent(in) :: cGW(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: PGW(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: T(nBas,nBas)
|
||||
double precision,intent(in) :: V(nBas,nBas)
|
||||
double precision,intent(in) :: J(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: K(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: SigC(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: Z(nBas,nspin)
|
||||
double precision,intent(in) :: dipole(ncart)
|
||||
|
||||
! Local variables
|
||||
|
||||
integer :: p
|
||||
integer :: ispin,ixyz
|
||||
double precision :: HOMO(nspin)
|
||||
double precision :: LUMO(nspin)
|
||||
double precision :: Gap(nspin)
|
||||
double precision :: ET(nspin)
|
||||
double precision :: EV(nspin)
|
||||
double precision :: EJ(nsp)
|
||||
double precision :: Ex(nspin)
|
||||
double precision :: Ec(nsp)
|
||||
double precision :: EqsGW
|
||||
double precision :: S_exact,S2_exact
|
||||
double precision :: S,S2
|
||||
double precision,external :: trace_matrix
|
||||
|
||||
! HOMO and LUMO
|
||||
|
||||
do ispin=1,nspin
|
||||
if(nO(ispin) > 0) then
|
||||
HOMO(ispin) = eGW(nO(ispin),ispin)
|
||||
LUMO(ispin) = eGW(nO(ispin)+1,ispin)
|
||||
Gap(ispin) = LUMO(ispin) - HOMO(ispin)
|
||||
else
|
||||
HOMO(ispin) = 0d0
|
||||
LUMO(ispin) = eGW(1,ispin)
|
||||
Gap(ispin) = 0d0
|
||||
end if
|
||||
end do
|
||||
|
||||
S2_exact = dble(nO(1) - nO(2))/2d0*(dble(nO(1) - nO(2))/2d0 + 1d0)
|
||||
S2 = S2_exact + nO(2) - sum(matmul(transpose(cGW(:,1:nO(1),1)),matmul(Ov,cGW(:,1:nO(2),2)))**2)
|
||||
|
||||
S_exact = 0.5d0*dble(nO(1) - nO(2))
|
||||
S = -0.5d0 + 0.5d0*sqrt(1d0 + 4d0*S2)
|
||||
|
||||
!------------------------------------------------------------------------
|
||||
! Compute total energy
|
||||
!------------------------------------------------------------------------
|
||||
|
||||
! Kinetic energy
|
||||
|
||||
do ispin=1,nspin
|
||||
ET(ispin) = trace_matrix(nBas,matmul(PGW(:,:,ispin),T(:,:)))
|
||||
end do
|
||||
|
||||
! Potential energy
|
||||
|
||||
do ispin=1,nspin
|
||||
EV(ispin) = trace_matrix(nBas,matmul(PGW(:,:,ispin),V(:,:)))
|
||||
end do
|
||||
|
||||
! Coulomb energy
|
||||
|
||||
EJ(1) = 0.5d0*trace_matrix(nBas,matmul(PGW(:,:,1),J(:,:,1)))
|
||||
EJ(2) = trace_matrix(nBas,matmul(PGW(:,:,1),J(:,:,2)))
|
||||
EJ(3) = 0.5d0*trace_matrix(nBas,matmul(PGW(:,:,2),J(:,:,2)))
|
||||
|
||||
! Exchange energy
|
||||
|
||||
do ispin=1,nspin
|
||||
Ex(ispin) = 0.5d0*trace_matrix(nBas,matmul(PGW(:,:,ispin),K(:,:,ispin)))
|
||||
end do
|
||||
|
||||
! Correlation energy
|
||||
|
||||
Ec(1) = 0.5d0*trace_matrix(nBas,matmul(PGW(:,:,1),SigC(:,:,1)))
|
||||
Ec(2) = trace_matrix(nBas,matmul(PGW(:,:,1),SigC(:,:,2)))
|
||||
Ec(3) = 0.5d0*trace_matrix(nBas,matmul(PGW(:,:,2),SigC(:,:,2)))
|
||||
|
||||
! Total energy
|
||||
|
||||
EqsGW = sum(ET(:)) + sum(EV(:)) + sum(EJ(:)) + sum(Ex(:)) + sum(Ec(:))
|
||||
|
||||
! Dump results
|
||||
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
if(nSCF < 10) then
|
||||
write(*,'(1X,A21,I1,A1,I1,A12)')' Self-consistent qsG',nSCF,'W',nSCF,' calculation'
|
||||
else
|
||||
write(*,'(1X,A21,I2,A1,I2,A12)')' Self-consistent qsG',nSCF,'W',nSCF,' calculation'
|
||||
endif
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
write(*,'(A1,A3,A1,A30,A1,A30,A1,A30,A1,A30,A1)') &
|
||||
'|',' ','|','e_HF ','|','Sig_c ','|','Z ','|','e_QP ','|'
|
||||
write(*,'(A1,A3,A1,2A15,A1,2A15,A1,2A15,A1,2A15,A1)') &
|
||||
'|','#','|','up ','dw ','|','up ','dw ','|','up ','dw ','|','up ','dw ','|'
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
|
||||
do p=1,nBas
|
||||
write(*,'(A1,I3,A1,2F15.6,A1,2F15.6,A1,2F15.6,A1,2F15.6,A1)') &
|
||||
'|',p,'|',eGW(p,1)*HaToeV,eGW(p,2)*HaToeV,'|',SigC(p,p,1)*HaToeV,SigC(p,p,2)*HaToeV,'|', &
|
||||
Z(p,1),Z(p,2),'|',eGW(p,1)*HaToeV,eGW(p,2)*HaToeV,'|'
|
||||
enddo
|
||||
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
write(*,'(2X,A10,I3)') 'Iteration ',nSCF
|
||||
write(*,'(2X,A19,F15.5)')'max(|FPS - SPF|) = ',Conv
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW HOMO energy (eV):',maxval(HOMO(:))*HaToeV
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW LUMO energy (eV):',minval(LUMO(:))*HaToeV
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW HOMO-LUMO gap (eV):',(minval(LUMO(:))-maxval(HOMO(:)))*HaToeV
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW total energy =',EqsGW + ENuc
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW exchange energy =',sum(Ex(:))
|
||||
write(*,'(2X,A30,F15.6)') 'qsGW correlation energy =',sum(Ec(:))
|
||||
write(*,'(2X,A30,F15.6)') 'RPA@qsGW correlation energy =',EcRPA
|
||||
write(*,*)'-------------------------------------------------------------------------------&
|
||||
-------------------------------------------------'
|
||||
write(*,*)
|
||||
|
||||
! Dump results for final iteration
|
||||
|
||||
if(Conv < thresh) then
|
||||
|
||||
write(*,*)
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,'(A40)') ' Summary '
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' One-electron energy: ',sum(ET(:)) + sum(EV(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' One-electron a energy: ',ET(1) + EV(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' One-electron b energy: ',ET(2) + EV(2),' au'
|
||||
write(*,*)
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Kinetic energy: ',sum(ET(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Kinetic a energy: ',ET(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Kinetic b energy: ',ET(2),' au'
|
||||
write(*,*)
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Potential energy: ',sum(EV(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Potential a energy: ',EV(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Potential b energy: ',EV(2),' au'
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Two-electron energy: ',sum(EJ(:)) + sum(Ex(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Two-electron aa energy: ',EJ(1) + Ex(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Two-electron ab energy: ',EJ(2),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Two-electron bb energy: ',EJ(3) + Ex(2),' au'
|
||||
write(*,*)
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Coulomb energy: ',sum(EJ(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Coulomb aa energy: ',EJ(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Coulomb ab energy: ',EJ(2),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Coulomb bb energy: ',EJ(3),' au'
|
||||
write(*,*)
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Exchange energy: ',sum(Ex(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Exchange a energy: ',Ex(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Exchange b energy: ',Ex(2),' au'
|
||||
write(*,*)
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Correlation energy: ',sum(Ec(:)),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Correlation aa energy: ',Ec(1),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Correlation ab energy: ',Ec(2),' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Correlation bb energy: ',Ec(3),' au'
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Electronic energy: ',EqsGW,' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' Nuclear repulsion: ',ENuc,' au'
|
||||
write(*,'(A40,1X,F16.10,A3)') ' qsUGW energy: ',EqsGW + ENuc,' au'
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,'(A40,F13.6)') ' S (exact) :',2d0*S_exact + 1d0
|
||||
write(*,'(A40,F13.6)') ' S :',2d0*S + 1d0
|
||||
write(*,'(A40,F13.6)') ' <S**2> (exact) :',S2_exact
|
||||
write(*,'(A40,F13.6)') ' <S**2> :',S2
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,'(A45)') ' Dipole moment (Debye) '
|
||||
write(*,'(19X,4A10)') 'X','Y','Z','Tot.'
|
||||
write(*,'(19X,4F10.6)') (dipole(ixyz)*auToD,ixyz=1,ncart),norm2(dipole)*auToD
|
||||
write(*,'(A60)') '-------------------------------------------------'
|
||||
write(*,*)
|
||||
|
||||
endif
|
||||
|
||||
end subroutine print_qsUGW
|
||||
|
|
@ -35,7 +35,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,SOSE
|
|||
double precision,intent(in) :: cHF(nBas,nBas)
|
||||
double precision,intent(in) :: PHF(nBas,nBas)
|
||||
double precision,intent(in) :: S(nBas,nBas)
|
||||
double precision,intent(in) :: T(nBas,nBAs)
|
||||
double precision,intent(in) :: T(nBas,nBas)
|
||||
double precision,intent(in) :: V(nBas,nBas)
|
||||
double precision,intent(in) :: Hc(nBas,nBas)
|
||||
double precision,intent(in) :: X(nBas,nBas)
|
||||
|
|
@ -131,7 +131,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,SOSE
|
|||
|
||||
! Initialization
|
||||
|
||||
nSCF = 0
|
||||
nSCF = -1
|
||||
n_diis = 0
|
||||
ispin = 1
|
||||
Conv = 1d0
|
||||
|
|
@ -147,6 +147,10 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,SOSE
|
|||
|
||||
do while(Conv > thresh .and. nSCF <= maxSCF)
|
||||
|
||||
! Increment
|
||||
|
||||
nSCF = nSCF + 1
|
||||
|
||||
! Buid Coulomb matrix
|
||||
|
||||
call Coulomb_matrix_AO_basis(nBas,P,ERI_AO_basis,J)
|
||||
|
|
@ -223,11 +227,7 @@ subroutine qsGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,SOSE
|
|||
! Print results
|
||||
|
||||
! call print_excitation('RPA ',ispin,nS,Omega(:,ispin))
|
||||
call print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigCp,Z,EcRPA,EcGM,EqsGW)
|
||||
|
||||
! Increment
|
||||
|
||||
nSCF = nSCF + 1
|
||||
call print_qsGW(nBas,nO,nSCF,Conv,thresh,eHF,eGW,c,ENuc,P,T,V,J,K,F,SigCp,Z,EcRPA,EqsGW)
|
||||
|
||||
enddo
|
||||
!------------------------------------------------------------------------
|
||||
|
|
|
|||
|
|
@ -0,0 +1,386 @@
|
|||
subroutine qsUGW(maxSCF,thresh,max_diis,doACFDT,exchange_kernel,doXBS,COHSEX,SOSEX,BSE,TDA_W,TDA, &
|
||||
G0W,GW0,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta,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,dipole_int_aa, &
|
||||
dipole_int_bb,PHF,cHF,eHF)
|
||||
|
||||
! Perform a quasiparticle self-consistent GW calculation
|
||||
|
||||
implicit none
|
||||
include 'parameters.h'
|
||||
|
||||
! Input variables
|
||||
|
||||
integer,intent(in) :: maxSCF
|
||||
integer,intent(in) :: max_diis
|
||||
double precision,intent(in) :: thresh
|
||||
logical,intent(in) :: doACFDT
|
||||
logical,intent(in) :: exchange_kernel
|
||||
logical,intent(in) :: doXBS
|
||||
logical,intent(in) :: COHSEX
|
||||
logical,intent(in) :: SOSEX
|
||||
logical,intent(in) :: BSE
|
||||
logical,intent(in) :: TDA_W
|
||||
logical,intent(in) :: TDA
|
||||
logical,intent(in) :: dBSE
|
||||
logical,intent(in) :: dTDA
|
||||
logical,intent(in) :: evDyn
|
||||
logical,intent(in) :: G0W
|
||||
logical,intent(in) :: GW0
|
||||
logical,intent(in) :: spin_conserved
|
||||
logical,intent(in) :: spin_flip
|
||||
double precision,intent(in) :: eta
|
||||
|
||||
integer,intent(in) :: nNuc
|
||||
double precision,intent(in) :: ZNuc(nNuc)
|
||||
double precision,intent(in) :: rNuc(nNuc,ncart)
|
||||
double precision,intent(in) :: ENuc
|
||||
|
||||
integer,intent(in) :: nBas
|
||||
integer,intent(in) :: nC(nspin)
|
||||
integer,intent(in) :: nO(nspin)
|
||||
integer,intent(in) :: nV(nspin)
|
||||
integer,intent(in) :: nR(nspin)
|
||||
integer,intent(in) :: nS(nspin)
|
||||
|
||||
double precision,intent(in) :: EUHF
|
||||
double precision,intent(in) :: eHF(nBas,nspin)
|
||||
double precision,intent(in) :: cHF(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: PHF(nBas,nBas,nspin)
|
||||
double precision,intent(in) :: S(nBas,nBas)
|
||||
double precision,intent(in) :: T(nBas,nBas)
|
||||
double precision,intent(in) :: V(nBas,nBas)
|
||||
double precision,intent(in) :: Hc(nBas,nBas)
|
||||
double precision,intent(in) :: X(nBas,nBas)
|
||||
double precision,intent(in) :: ERI_AO(nBas,nBas,nBas,nBas)
|
||||
double precision,intent(inout):: ERI_aaaa(nBas,nBas,nBas,nBas)
|
||||
double precision,intent(inout):: ERI_aabb(nBas,nBas,nBas,nBas)
|
||||
double precision,intent(inout):: ERI_bbbb(nBas,nBas,nBas,nBas)
|
||||
double precision,intent(in) :: dipole_int(nBas,nBas,ncart)
|
||||
double precision,intent(in) :: dipole_int_aa(nBas,nBas,ncart)
|
||||
double precision,intent(in) :: dipole_int_bb(nBas,nBas,ncart)
|
||||
|
||||
! Local variables
|
||||
|
||||
logical :: doGWPT = .false.
|
||||
|
||||
integer :: nSCF
|
||||
integer :: nBasSq
|
||||
integer :: ispin
|
||||
integer :: is
|
||||
integer :: n_diis
|
||||
integer :: nS_aa,nS_bb,nS_sc
|
||||
double precision :: dipole(ncart)
|
||||
|
||||
double precision :: EcRPA
|
||||
double precision :: EcBSE(nspin)
|
||||
double precision :: EcAC(nspin)
|
||||
double precision :: Conv
|
||||
double precision :: rcond(nspin)
|
||||
double precision,external :: trace_matrix
|
||||
double precision,allocatable :: error_diis(:,:,:)
|
||||
double precision,allocatable :: F_diis(:,:,:)
|
||||
double precision,allocatable :: OmRPA(:)
|
||||
double precision,allocatable :: XpY_RPA(:,:)
|
||||
double precision,allocatable :: XmY_RPA(:,:)
|
||||
double precision,allocatable :: rho_RPA(:,:,:,:)
|
||||
double precision,allocatable :: c(:,:,:)
|
||||
double precision,allocatable :: cp(:,:,:)
|
||||
double precision,allocatable :: eGW(:,:)
|
||||
double precision,allocatable :: P(:,:,:)
|
||||
double precision,allocatable :: F(:,:,:)
|
||||
double precision,allocatable :: Fp(:,:,:)
|
||||
double precision,allocatable :: J(:,:,:)
|
||||
double precision,allocatable :: K(:,:,:)
|
||||
double precision,allocatable :: SigC(:,:,:)
|
||||
double precision,allocatable :: SigCp(:,:,:)
|
||||
double precision,allocatable :: SigCm(:,:,:)
|
||||
double precision,allocatable :: Z(:,:)
|
||||
double precision,allocatable :: error(:,:,:)
|
||||
|
||||
! Hello world
|
||||
|
||||
write(*,*)
|
||||
write(*,*)'*************************************************'
|
||||
write(*,*)'| Self-consistent unrestricted qsGW calculation |'
|
||||
write(*,*)'*************************************************'
|
||||
write(*,*)
|
||||
|
||||
! Warning
|
||||
|
||||
write(*,*) '!! ERIs in MO basis will be overwritten in qsUGW !!'
|
||||
write(*,*)
|
||||
|
||||
! Stuff
|
||||
|
||||
nBasSq = nBas*nBas
|
||||
|
||||
! SOSEX correction
|
||||
|
||||
if(SOSEX) then
|
||||
write(*,*) 'SOSEX correction activated but BUG!'
|
||||
stop
|
||||
end if
|
||||
|
||||
! COHSEX approximation
|
||||
|
||||
if(COHSEX) then
|
||||
write(*,*) 'COHSEX approximation activated!'
|
||||
write(*,*)
|
||||
end if
|
||||
|
||||
! TDA for W
|
||||
|
||||
if(TDA_W) then
|
||||
write(*,*) 'Tamm-Dancoff approximation for dynamic screening!'
|
||||
write(*,*)
|
||||
end if
|
||||
|
||||
! TDA
|
||||
|
||||
if(TDA) then
|
||||
write(*,*) 'Tamm-Dancoff approximation activated!'
|
||||
write(*,*)
|
||||
end if
|
||||
|
||||
! Memory allocation
|
||||
|
||||
nS_aa = nS(1)
|
||||
nS_bb = nS(2)
|
||||
nS_sc = nS_aa + nS_bb
|
||||
|
||||
allocate(eGW(nBas,nspin),c(nBas,nBas,nspin),cp(nBas,nBas,nspin),P(nBas,nBas,nspin),F(nBas,nBas,nspin),Fp(nBas,nBas,nspin), &
|
||||
J(nBas,nBas,nspin),K(nBas,nBas,nspin),SigC(nBas,nBas,nspin),SigCp(nBas,nBas,nspin),SigCm(nBas,nBas,nspin), &
|
||||
Z(nBas,nspin),OmRPA(nS_sc),XpY_RPA(nS_sc,nS_sc),XmY_RPA(nS_sc,nS_sc),rho_RPA(nBas,nBas,nS_sc,nspin), &
|
||||
error(nBas,nBas,nspin),error_diis(nBasSq,max_diis,nspin),F_diis(nBasSq,max_diis,nspin))
|
||||
|
||||
! Initialization
|
||||
|
||||
nSCF = -1
|
||||
n_diis = 0
|
||||
ispin = 1
|
||||
Conv = 1d0
|
||||
P(:,:,:) = PHF(:,:,:)
|
||||
eGW(:,:) = eHF(:,:)
|
||||
c(:,:,:) = cHF(:,:,:)
|
||||
F_diis(:,:,:) = 0d0
|
||||
error_diis(:,:,:) = 0d0
|
||||
|
||||
!------------------------------------------------------------------------
|
||||
! Main loop
|
||||
!------------------------------------------------------------------------
|
||||
|
||||
do while(Conv > thresh .and. nSCF < maxSCF)
|
||||
|
||||
! Increment
|
||||
|
||||
nSCF = nSCF + 1
|
||||
|
||||
! Buid Coulomb matrix
|
||||
|
||||
do is=1,nspin
|
||||
call Coulomb_matrix_AO_basis(nBas,P(:,:,is),ERI_AO(:,:,:,:),J(:,:,is))
|
||||
end do
|
||||
|
||||
! Compute exchange part of the self-energy
|
||||
|
||||
do is=1,nspin
|
||||
call exchange_matrix_AO_basis(nBas,P(:,:,is),ERI_AO(:,:,:,:),K(:,:,is))
|
||||
end do
|
||||
|
||||
!--------------------------------------------------
|
||||
! AO to MO transformation of two-electron integrals
|
||||
!--------------------------------------------------
|
||||
|
||||
! 4-index transform for (aa|aa) block
|
||||
|
||||
call AOtoMO_integral_transform(1,1,1,1,nBas,c,ERI_AO,ERI_aaaa)
|
||||
|
||||
! 4-index transform for (aa|bb) block
|
||||
|
||||
call AOtoMO_integral_transform(1,1,2,2,nBas,c,ERI_AO,ERI_aabb)
|
||||
|
||||
! 4-index transform for (bb|bb) block
|
||||
|
||||
call AOtoMO_integral_transform(2,2,2,2,nBas,c,ERI_AO,ERI_bbbb)
|
||||
|
||||
! Compute linear response
|
||||
|
||||
if(.not. GW0 .or. nSCF == 0) then
|
||||
|
||||
call unrestricted_linear_response(ispin,.true.,TDA_W,.false.,eta,nBas,nC,nO,nV,nR,nS_aa,nS_bb,nS_sc,nS_sc,1d0, &
|
||||
eGW,ERI_aaaa,ERI_aabb,ERI_bbbb,OmRPA,rho_RPA,EcRPA,OmRPA,XpY_RPA,XmY_RPA)
|
||||
|
||||
endif
|
||||
|
||||
!----------------------!
|
||||
! Excitation densities !
|
||||
!----------------------!
|
||||
|
||||
call unrestricted_excitation_density(nBas,nC,nO,nR,nS_aa,nS_bb,nS_sc,ERI_aaaa,ERI_aabb,ERI_bbbb,XpY_RPA,rho_RPA)
|
||||
|
||||
!------------------------------------------------!
|
||||
! Compute self-energy and renormalization factor !
|
||||
!------------------------------------------------!
|
||||
|
||||
if(G0W) then
|
||||
|
||||
call unrestricted_self_energy_correlation_diag(eta,nBas,nC,nO,nV,nR,nS_sc,eHF,OmRPA,rho_RPA,SigC)
|
||||
call unrestricted_renormalization_factor(eta,nBas,nC,nO,nV,nR,nS_sc,eHF,OmRPA,rho_RPA,Z)
|
||||
|
||||
else
|
||||
|
||||
call unrestricted_self_energy_correlation_diag(eta,nBas,nC,nO,nV,nR,nS_sc,eGW,OmRPA,rho_RPA,SigC)
|
||||
call unrestricted_renormalization_factor(eta,nBas,nC,nO,nV,nR,nS_sc,eGW,OmRPA,rho_RPA,Z)
|
||||
|
||||
endif
|
||||
|
||||
! Make correlation self-energy Hermitian and transform it back to AO basis
|
||||
|
||||
do is=1,nspin
|
||||
SigCp(:,:,is) = 0.5d0*(SigC(:,:,is) + transpose(SigC(:,:,is)))
|
||||
SigCm(:,:,is) = 0.5d0*(SigC(:,:,is) - transpose(SigC(:,:,is)))
|
||||
end do
|
||||
|
||||
do is=1,nspin
|
||||
call MOtoAO_transform(nBas,S,c(:,:,is),SigCp(:,:,is))
|
||||
end do
|
||||
|
||||
! Solve the quasi-particle equation
|
||||
|
||||
do is=1,nspin
|
||||
F(:,:,is) = Hc(:,:) + J(:,:,is) + J(:,:,mod(is,2)+1) + K(:,:,is) + SigCp(:,:,is)
|
||||
end do
|
||||
|
||||
! Check convergence
|
||||
|
||||
do is=1,nspin
|
||||
error(:,:,is) = matmul(F(:,:,is),matmul(P(:,:,is),S(:,:))) - matmul(matmul(S(:,:),P(:,:,is)),F(:,:,is))
|
||||
end do
|
||||
|
||||
if(nSCF > 1) conv = maxval(abs(error(:,:,:)))
|
||||
|
||||
! DIIS extrapolation
|
||||
|
||||
n_diis = min(n_diis+1,max_diis)
|
||||
do is=1,nspin
|
||||
if(nO(is) > 1) call DIIS_extrapolation(rcond(is),nBasSq,nBasSq,n_diis,error_diis(:,1:n_diis,is), &
|
||||
F_diis(:,1:n_diis,is),error(:,:,is),F(:,:,is))
|
||||
end do
|
||||
|
||||
! Reset DIIS if required
|
||||
|
||||
if(minval(rcond(:)) < 1d-15) n_diis = 0
|
||||
|
||||
! Transform Fock matrix in orthogonal basis
|
||||
|
||||
do is=1,nspin
|
||||
Fp(:,:,is) = matmul(transpose(X(:,:)),matmul(F(:,:,is),X(:,:)))
|
||||
end do
|
||||
|
||||
! Diagonalize Fock matrix to get eigenvectors and eigenvalues
|
||||
|
||||
cp(:,:,:) = Fp(:,:,:)
|
||||
do is=1,nspin
|
||||
call diagonalize_matrix(nBas,cp(:,:,is),eGW(:,is))
|
||||
end do
|
||||
|
||||
! Back-transform eigenvectors in non-orthogonal basis
|
||||
|
||||
do is=1,nspin
|
||||
c(:,:,is) = matmul(X(:,:),cp(:,:,is))
|
||||
end do
|
||||
|
||||
! Compute density matrix
|
||||
|
||||
do is=1,nspin
|
||||
P(:,:,is) = matmul(c(:,1:nO(is),is),transpose(c(:,1:nO(is),is)))
|
||||
end do
|
||||
|
||||
! Print results
|
||||
|
||||
call dipole_moment(nBas,P(:,:,1)+P(:,:,2),nNuc,ZNuc,rNuc,dipole_int,dipole)
|
||||
call print_qsUGW(nBas,nO,S,nSCF,Conv,thresh,eGW,c,P,T,V,J,K,ENuc,EHF,SigC,Z,EcRPA,dipole)
|
||||
|
||||
enddo
|
||||
!------------------------------------------------------------------------
|
||||
! End main loop
|
||||
!------------------------------------------------------------------------
|
||||
|
||||
! Compute second-order correction of the Hermitization error
|
||||
|
||||
!if(doGWPT) call qsGW_PT(nBas,nC,nO,nV,nR,nS,eGW,SigCm)
|
||||
|
||||
! Did it actually converge?
|
||||
|
||||
if(nSCF == maxSCF) then
|
||||
|
||||
write(*,*)
|
||||
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
|
||||
write(*,*)' Convergence failed '
|
||||
write(*,*)'!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!'
|
||||
write(*,*)
|
||||
|
||||
stop
|
||||
|
||||
endif
|
||||
|
||||
! Deallocate memory
|
||||
|
||||
deallocate(c,cp,P,F,Fp,J,K,SigC,SigCp,SigCm,Z,OmRPA,XpY_RPA,XmY_RPA,rho_RPA,error,error_diis,F_diis)
|
||||
|
||||
! Perform BSE calculation
|
||||
|
||||
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,ERI_aaaa,ERI_aabb,ERI_bbbb,dipole_int_aa,dipole_int_bb,cHF,eGW,eGW,EcBSE)
|
||||
|
||||
if(exchange_kernel) then
|
||||
|
||||
EcBSE(1) = 0.5d0*EcBSE(1)
|
||||
EcBSE(2) = 1.5d0*EcBSE(2)
|
||||
|
||||
end if
|
||||
|
||||
write(*,*)
|
||||
write(*,*)'-------------------------------------------------------------------------------'
|
||||
write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsUGW correlation energy (spin-conserved) =',EcBSE(1)
|
||||
write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsUGW correlation energy (spin-flip) =',EcBSE(2)
|
||||
write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsUGW correlation energy =',EcBSE(1) + EcBSE(2)
|
||||
write(*,'(2X,A50,F20.10)') 'Tr@BSE@qsUGW total energy =',ENuc + EUHF + EcBSE(1) + EcBSE(2)
|
||||
write(*,*)'-------------------------------------------------------------------------------'
|
||||
write(*,*)
|
||||
|
||||
! Compute the BSE correlation energy via the adiabatic connection
|
||||
|
||||
if(doACFDT) then
|
||||
|
||||
write(*,*) '--------------------------------------------------------------'
|
||||
write(*,*) ' Adiabatic connection version of BSE@qsUGW correlation energy '
|
||||
write(*,*) '--------------------------------------------------------------'
|
||||
write(*,*)
|
||||
|
||||
if(doXBS) then
|
||||
|
||||
write(*,*) '*** scaled screening version (XBS) ***'
|
||||
write(*,*)
|
||||
|
||||
end if
|
||||
|
||||
call unrestricted_ACFDT(exchange_kernel,doXBS,.true.,TDA_W,TDA,BSE,spin_conserved,spin_flip, &
|
||||
eta,nBas,nC,nO,nV,nR,nS,ERI_aaaa,ERI_aabb,ERI_bbbb,eGW,eGW,EcAC)
|
||||
|
||||
write(*,*)
|
||||
write(*,*)'-------------------------------------------------------------------------------'
|
||||
write(*,'(2X,A50,F20.10)') 'AC@BSE@qsUGW correlation energy (spin-conserved) =',EcAC(1)
|
||||
write(*,'(2X,A50,F20.10)') 'AC@BSE@qsUGW correlation energy (spin-flip) =',EcAC(2)
|
||||
write(*,'(2X,A50,F20.10)') 'AC@BSE@qsUGW correlation energy =',EcAC(1) + EcAC(2)
|
||||
write(*,'(2X,A50,F20.10)') 'AC@BSE@qsUGW total energy =',ENuc + EUHF + EcAC(1) + EcAC(2)
|
||||
write(*,*)'-------------------------------------------------------------------------------'
|
||||
write(*,*)
|
||||
|
||||
end if
|
||||
|
||||
end if
|
||||
|
||||
end subroutine qsUGW
|
||||
|
|
@ -0,0 +1,94 @@
|
|||
subroutine unrestricted_self_energy_correlation(eta,nBas,nC,nO,nV,nR,nSt,e,Omega,rho,SigC)
|
||||
|
||||
! Compute diagonal of the correlation part of the self-energy
|
||||
|
||||
implicit none
|
||||
include 'parameters.h'
|
||||
|
||||
! Input variables
|
||||
|
||||
double precision,intent(in) :: eta
|
||||
integer,intent(in) :: nBas
|
||||
integer,intent(in) :: nC(nspin)
|
||||
integer,intent(in) :: nO(nspin)
|
||||
integer,intent(in) :: nV(nspin)
|
||||
integer,intent(in) :: nR(nspin)
|
||||
integer,intent(in) :: nSt
|
||||
double precision,intent(in) :: e(nBas,nspin)
|
||||
double precision,intent(in) :: Omega(nSt)
|
||||
double precision,intent(in) :: rho(nBas,nBas,nSt,nspin)
|
||||
|
||||
! Local variables
|
||||
|
||||
integer :: i,j,a,b,p,q,jb
|
||||
double precision :: eps
|
||||
|
||||
! Output variables
|
||||
|
||||
double precision,intent(out) :: SigC(nBas,nBas,nspin)
|
||||
|
||||
! Initialize
|
||||
|
||||
SigC(:,:,:) = 0d0
|
||||
|
||||
!--------------!
|
||||
! Spin-up part !
|
||||
!--------------!
|
||||
|
||||
! Occupied part of the correlation self-energy
|
||||
|
||||
do p=nC(1)+1,nBas-nR(1)
|
||||
do q=nC(1)+1,nBas-nR(1)
|
||||
do i=nC(1)+1,nO(1)
|
||||
do jb=1,nSt
|
||||
eps = e(p,1) - e(i,1) + Omega(jb)
|
||||
SigC(p,q,1) = SigC(p,q,1) + rho(p,i,jb,1)*rho(q,i,jb,1)*eps/(eps**2 + eta**2)
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
|
||||
! Virtual part of the correlation self-energy
|
||||
|
||||
do p=nC(1)+1,nBas-nR(1)
|
||||
do q=nC(1)+1,nBas-nR(1)
|
||||
do a=nO(1)+1,nBas-nR(1)
|
||||
do jb=1,nSt
|
||||
eps = e(p,1) - e(a,1) - Omega(jb)
|
||||
SigC(p,q,1) = SigC(p,q,1) + rho(p,a,jb,1)*rho(q,a,jb,1)*eps/(eps**2 + eta**2)
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
|
||||
!----------------!
|
||||
! Spin-down part !
|
||||
!----------------!
|
||||
|
||||
! Occupied part of the correlation self-energy
|
||||
|
||||
do p=nC(2)+1,nBas-nR(2)
|
||||
do q=nC(2)+1,nBas-nR(2)
|
||||
do i=nC(2)+1,nO(2)
|
||||
do jb=1,nSt
|
||||
eps = e(p,2) - e(i,2) + Omega(jb)
|
||||
SigC(p,q,2) = SigC(p,q,2) + rho(p,i,jb,2)*rho(q,i,jb,2)*eps/(eps**2 + eta**2)
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
|
||||
! Virtual part of the correlation self-energy
|
||||
|
||||
do p=nC(2)+1,nBas-nR(2)
|
||||
do q=nC(2)+1,nBas-nR(2)
|
||||
do a=nO(2)+1,nBas-nR(2)
|
||||
do jb=1,nSt
|
||||
eps = e(p,2) - e(a,2) - Omega(jb)
|
||||
SigC(p,q,2) = SigC(p,q,2) + rho(p,a,jb,2)*rho(q,a,jb,2)*eps/(eps**2 + eta**2)
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
end do
|
||||
|
||||
end subroutine unrestricted_self_energy_correlation
|
||||
|
|
@ -884,9 +884,22 @@ program QuAcK
|
|||
if(doqsGW) then
|
||||
|
||||
call cpu_time(start_qsGW)
|
||||
call qsGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,SOSEX, &
|
||||
BSE,TDA_W,TDA,G0W,GW0,dBSE,dTDA,evDyn,singlet,triplet,eta_GW, &
|
||||
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF)
|
||||
|
||||
if(unrestricted) then
|
||||
|
||||
call qsUGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,SOSEX,BSE,TDA_W,TDA, &
|
||||
G0W,GW0,dBSE,dTDA,evDyn,spin_conserved,spin_flip,eta_GW,nNuc,ZNuc,rNuc,ENuc,nBas,nC,nO, &
|
||||
nV,nR,nS,EUHF,S,X,T,V,Hc,ERI_AO,ERI_MO_aaaa,ERI_MO_aabb,ERI_MO_bbbb,dipole_int,dipole_int_aa, &
|
||||
dipole_int_bb,PHF,cHF,eHF)
|
||||
|
||||
else
|
||||
|
||||
call qsGW(maxSCF_GW,thresh_GW,n_diis_GW,doACFDT,exchange_kernel,doXBS,COHSEX,SOSEX, &
|
||||
BSE,TDA_W,TDA,G0W,GW0,dBSE,dTDA,evDyn,singlet,triplet,eta_GW, &
|
||||
nBas,nC,nO,nV,nR,nS,ENuc,ERHF,S,X,T,V,Hc,ERI_AO,ERI_MO,dipole_int,PHF,cHF,eHF)
|
||||
|
||||
end if
|
||||
|
||||
call cpu_time(end_qsGW)
|
||||
|
||||
t_qsGW = end_qsGW - start_qsGW
|
||||
|
|
|
|||
Loading…
Reference in New Issue