10
1
mirror of https://github.com/pfloos/quack synced 2024-11-03 20:53:53 +01:00

eDFT code OK

This commit is contained in:
Clotilde Marut 2021-12-16 21:42:58 +01:00
parent 1de6de072e
commit 2b58358681
6 changed files with 39 additions and 31 deletions

View File

@ -13,31 +13,31 @@
# GGA = 2: LYP,PBE
# MGGA = 3:
# Hybrid = 4: HF,B3LYP,PBE
1 VWN5
1 VWN5
# quadrature grid SG-n
1
# Number of states in ensemble (nEns)
2
4
# occupation numbers
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
# Ensemble weights: wEns(1),...,wEns(nEns-1)
0.55 0.0 0.0
0.3333 0.3333 0.0
# Ncentered ?
T
# Parameters for CC weight-dependent exchange functional
4
0.642674 -0.07818 -0.0280307 0.00144198
0.254939 -0.0893405 0.00765581 0.
0.502145 -0.0672146 0.00786214 -0.00202771
-1.28842 -0.173117 0.0900511 -0.0118975
0.0 0.0 0.0 0.0
# choice of UCC exchange coefficient : 1 for Cx1, 2 for Cx2, 3 for Cx1*Cx2
1
2

View File

@ -1,4 +1,5 @@
subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,weight,rhow,Cx_choice,doNcentered,ExDD)
subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,weight,rhow,Cx_choice,doNcentered,&
kappa,ExDD)
! Compute the unrestricted version of the curvature-corrected exchange ensemble derivative
@ -16,6 +17,7 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei
double precision,intent(in) :: rhow(nGrid)
integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Local variables
@ -150,12 +152,26 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei
ExDD(:) = 0d0
do iEns=1,nEns
do jEns=2,nEns
if (doNcentered) then
ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*dExdw(jEns)
do iEns=1,nEns
do jEns=2,nEns
ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - kappa(iEns)* wEns(jEns))*dExdw(jEns)
end do
end do
end do
else
do iEns=1,nEns
do jEns=2,nEns
ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*dExdw(jEns)
end do
end do
endif
end subroutine UCC_lda_exchange_derivative_discontinuity

View File

@ -190,7 +190,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux
write(*,'(A60)') '-------------------------------------------------'
do iEns=2,nEns
write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Om(iEns), ' au'
! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Om(iEns), ' au'
write(*,*)
write(*,'(A44, F16.10,A3)') ' H energy contribution : ',OmH(iEns), ' au'
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns), ' au'
@ -204,7 +204,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Om(iEns)*HaToeV, ' eV'
! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Om(iEns)*HaToeV, ' eV'
write(*,*)
write(*,'(A44, F16.10,A3)') ' H energy contribution : ',OmH(iEns)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns)*HaToeV, ' eV'

View File

@ -1,5 +1,5 @@
subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nCC,aCC,nGrid,weight,rhow,drhow,&
Cx_choice,doNcentered,ExDD)
Cx_choice,doNcentered,kappa,ExDD)
! Compute the exchange part of the derivative discontinuity
@ -20,6 +20,7 @@ subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nCC
double precision,intent(in) :: drhow(ncart,nGrid)
integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Local variables
@ -41,7 +42,7 @@ subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nCC
case(1)
call unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),nCC,aCC,nGrid,weight(:),&
rhow(:),Cx_choice,doNcentered,ExDD(:))
rhow(:),Cx_choice,doNcentered,kappa,ExDD(:))
! GGA functionals
case(2)

View File

@ -150,21 +150,11 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
do ispin=1,nspin
call unrestricted_exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,nCC,aCC,nGrid,weight, &
rhow(:,ispin),drhow(:,:,ispin),Cx_choice,doNcentered,ExDD(ispin,:))
rhow(:,ispin),drhow(:,:,ispin),Cx_choice,doNcentered,kappa,ExDD(ispin,:))
end do
call unrestricted_correlation_derivative_discontinuity(c_rung,c_DFA,nEns,wEns,nGrid,weight,rhow,drhow,EcDD)
! Scaling derivative discontinuity for N-centered ensembles
if(doNcentered) then
do iEns=1,nEns
ExDD(:,iEns) = (1d0 - kappa(iEns))*ExDD(:,iEns)
EcDD(:,iEns) = (1d0 - kappa(iEns))*EcDD(:,iEns)
end do
end if
!------------------------------------------------------------------------
! Total energy

View File

@ -1,5 +1,5 @@
subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC,aCC,nGrid,weight,rhow,&
Cx_choice,doNcentered,ExDD)
Cx_choice,doNcentered,kappa,ExDD)
! Compute the exchange LDA part of the derivative discontinuity
@ -19,6 +19,7 @@ subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC,
double precision,intent(in) :: rhow(nGrid)
integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Local variables
@ -38,7 +39,7 @@ subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC,
case (2)
call UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,weight,rhow,&
Cx_choice,doNcentered,ExDD)
Cx_choice,doNcentered,kappa,ExDD)
case default