PTEROSOR/QuAcK
10
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2025-01-03 18:16:03 +01:00
Code Releases Activity

fix N-centered

Browse Source
This commit is contained in:
Pierre-Francois Loos 2021-12-06 15:27:39 +01:00
parent c14f6c9aa1
commit bbe1c4b2ff
18 changed files with 147 additions and 213 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

8
input/dft
View File

@ -6,18 +6,18 @@
# GGA = 2: B88,G96,PBE # GGA = 2: B88,G96,PBE
# MGGA = 3: # MGGA = 3:
# Hybrid = 4 HF,B3LYP,PBE # Hybrid = 4 HF,B3LYP,PBE
1 CC-S51 1 S51
# correlation rung: # correlation rung:
# Hartree = 0: H # Hartree = 0: H
# LDA = 1: PW92,VWN3,VWN5,eVWN5 # LDA = 1: PW92,VWN3,VWN5,eVWN5
# GGA = 2: LYP,PBE # GGA = 2: LYP,PBE
# MGGA = 3: # MGGA = 3:
# Hybrid = 4: HF,B3LYP,PBE # Hybrid = 4: HF,B3LYP,PBE
0 H 1 VWN5
# quadrature grid SG-n # quadrature grid SG-n
1 1
# Number of states in ensemble (nEns) # Number of states in ensemble (nEns)
4 2
# occupation numbers # 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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
@ -31,7 +31,7 @@
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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) # Ensemble weights: wEns(1),...,wEns(nEns-1)
1 0.0 0.0 0.55 0.0 0.0
# Ncentered ? # Ncentered ?
T T
# Parameters for CC weight-dependent exchange functional # Parameters for CC weight-dependent exchange functional

13
src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90
View File

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

48
src/eDFT/UCC_lda_exchange_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rhow,rho,Cx_choice,doNcentered,kappa,LZx,Ex) subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rhow,rho,Cx_choice,doNcentered,LZx,Ex)
! Compute the unrestricted version of the curvature-corrected exchange functional ! Compute the unrestricted version of the curvature-corrected exchange functional
@ -18,7 +18,6 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho
double precision,intent(in) :: rho(nGrid,nspin,nEns) double precision,intent(in) :: rho(nGrid,nspin,nEns)
integer,intent(in) :: Cx_choice integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Local variables ! Local variables
@ -32,7 +31,8 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho
! Output variables ! Output variables
double precision,intent(out) :: LZx(nspin,nEns), Ex(nspin,nEns) double precision,intent(out) :: LZx(nspin)
double precision,intent(out) :: Ex(nspin,nEns)
! Defining enhancements factor for weight-dependent functionals ! Defining enhancements factor for weight-dependent functionals
@ -52,7 +52,6 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho
c2 = aCC(3,2) c2 = aCC(3,2)
d2 = aCC(4,2) d2 = aCC(4,2)
w1 = wEns(2) w1 = wEns(2)
Fx1 = 1d0 + a1*w1 + b1*w1**2 + c1*w1**3 + d1*w1**4 Fx1 = 1d0 + a1*w1 + b1*w1**2 + c1*w1**3 + d1*w1**4
@ -100,7 +99,7 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho
! Compute LDA exchange matrix in the AO basis ! Compute LDA exchange matrix in the AO basis
Ex(:,:) = 0d0 Ex(:,:) = 0d0
LZx(:,:) = 0d0 LZx(:) = 0d0
do ispin=1,nspin do ispin=1,nspin
@ -108,43 +107,20 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho
r = max(0d0,rhow(iG,ispin)) r = max(0d0,rhow(iG,ispin))
if(doNcentered) then if(r > threshold) then
if(r > threshold) then e = Cx*r**(+1d0/3d0)
dedr = 1d0/3d0*Cx*r**(-2d0/3d0)
e = Cx*r**(+1d0/3d0) LZx(ispin) = LZx(ispin) - weight(iG)*dedr*r*r
dedr = 1d0/3d0*Cx*r**(-2d0/3d0)
do iEns=1,nEns do iEns=1,nEns
rI = max(0d0,rho(iG,ispin,iEns)) rI = max(0d0,rho(iG,ispin,iEns))
LZx(ispin,iEns) = LZx(ispin,iEns) - weight(iG)*kappa(iEns)*dedr*r*r if(rI > threshold) Ex(ispin,iEns) = Ex(ispin,iEns) + weight(iG)*(e+dedr*r)*rI
if(rI > threshold) Ex(ispin,iEns) = Ex(ispin,iEns) + weight(iG)*(kappa(iEns)*e+dedr*r)*rI end do
end do
endif
else
if(r > threshold) then
e = Cx*r**(+1d0/3d0)
dedr = 1d0/3d0*Cx*r**(-2d0/3d0)
LZx(ispin,:) = LZx(ispin,:) - weight(iG)*dedr*r*r
do iEns=1,nEns
rI = max(0d0,rho(iG,ispin,iEns))
if(rI > threshold) Ex(ispin,iEns) = Ex(ispin,iEns) + weight(iG)*(e+dedr*r)*rI
end do
endif
endif endif

33
src/eDFT/US51_lda_exchange_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,doNcentered,kappa,LZx,Ex) subroutine US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,LZx,Ex)
! Compute the restricted version of Slater's LDA exchange individual energy ! Compute the restricted version of Slater's LDA exchange individual energy
@ -12,8 +12,6 @@ subroutine US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,doNcen
double precision,intent(in) :: weight(nGrid) double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid,nspin) double precision,intent(in) :: rhow(nGrid,nspin)
double precision,intent(in) :: rho(nGrid,nspin,nEns) double precision,intent(in) :: rho(nGrid,nspin,nEns)
logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Local variables ! Local variables
@ -27,10 +25,10 @@ subroutine US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,doNcen
! Output variables ! Output variables
double precision,intent(out) :: LZx(nspin,nEns) double precision,intent(out) :: LZx(nspin)
double precision,intent(out) :: Ex(nspin,nEns) double precision,intent(out) :: Ex(nspin,nEns)
LZx(:,:) = 0d0 LZx(:) = 0d0
Ex(:,:) = 0d0 Ex(:,:) = 0d0
do ispin=1,nspin do ispin=1,nspin
@ -39,33 +37,12 @@ subroutine US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,doNcen
r = max(0d0,rhow(iG,ispin)) r = max(0d0,rhow(iG,ispin))
! if(doNcentered) then
! if(r > threshold) then
! e = CxLSDA*r**(+1d0/3d0)
! dedr = 1d0/3d0*CxLSDA*r**(-2d0/3d0)
! do iEns=1,nEns
! rI = max(0d0,rho(iG,ispin,iEns))
! LZx(ispin,iEns) = LZx(ispin,iEns) - weight(iG)*kappa(iEns)*dedr*r*r
! if(rI > threshold) Ex(ispin,iEns) = Ex(ispin,iEns) + weight(iG)*(kappa(iEns)*e+dedr*r)*rI
! end do
! endif
! else
if(r > threshold) then if(r > threshold) then
e = CxLSDA*r**(+1d0/3d0) e = CxLSDA*r**(+1d0/3d0)
dedr = 1d0/3d0*CxLSDA*r**(-2d0/3d0) dedr = 1d0/3d0*CxLSDA*r**(-2d0/3d0)
LZx(ispin,:) = LZx(ispin,:) - weight(iG)*dedr*r*r LZx(ispin) = LZx(ispin) - weight(iG)*dedr*r*r
do iEns=1,nEns do iEns=1,nEns
@ -77,8 +54,6 @@ subroutine US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,doNcen
endif endif
! endif
enddo enddo
enddo enddo

3
src/eDFT/UVWN5_lda_correlation_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine UVWN5_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,doNcentered,LZc,Ec) subroutine UVWN5_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,LZc,Ec)
! Compute VWN5 LDA correlation potential ! Compute VWN5 LDA correlation potential
@ -13,7 +13,6 @@ subroutine UVWN5_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,do
double precision,intent(in) :: weight(nGrid) double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid,nspin) double precision,intent(in) :: rhow(nGrid,nspin)
double precision,intent(in) :: rho(nGrid,nspin,nEns) double precision,intent(in) :: rho(nGrid,nspin,nEns)
logical,intent(in) :: doNcentered
! Local variables ! Local variables

53
src/eDFT/print_unrestricted_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux,ExDD,EcDD,E,Om,Omx,Omc,Omaux,OmxDD,OmcDD) subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux,ExDD,EcDD,E,Om,OmH,Omx,Omc,Omaux,OmxDD,OmcDD)
! Print individual energies for eDFT calculation ! Print individual energies for eDFT calculation
@ -20,6 +20,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux
double precision,intent(in) :: EcDD(nsp,nEns) double precision,intent(in) :: EcDD(nsp,nEns)
double precision,intent(in) :: E(nEns) double precision,intent(in) :: E(nEns)
double precision,intent(in) :: OmH(nEns)
double precision,intent(in) :: Omx(nEns) double precision,intent(in) :: Omx(nEns)
double precision,intent(in) :: Omc(nEns) double precision,intent(in) :: Omc(nEns)
double precision,intent(in) :: Omaux(nEns) double precision,intent(in) :: Omaux(nEns)
@ -153,32 +154,32 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux
! Total Energy and IP and EA ! Total Energy and IP and EA
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! write(*,'(A60)') ' ENERGY DIFFERENCES FROM AUXILIARY ENERGIES ' write(*,'(A60)') ' ENERGY DIFFERENCES FROM AUXILIARY ENERGIES '
! write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! do iEns=2,nEns do iEns=2,nEns
! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Omaux(iEns)+OmxDD(iEns)+OmcDD(iEns),' au' write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Omaux(iEns)+OmxDD(iEns)+OmcDD(iEns),' au'
! write(*,*) write(*,*)
! write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns), ' au' write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns), ' au'
! write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns), ' au' write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns), ' au'
! write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns), ' au' write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns), ' au'
! write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxDD(iEns)+OmcDD(iEns), ' au' write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxDD(iEns)+OmcDD(iEns), ' au'
! write(*,*) write(*,*)
! write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',(Omaux(iEns)+OmxDD(iEns)+OmcDD(iEns))*HaToeV,' eV' write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',(Omaux(iEns)+OmxDD(iEns)+OmcDD(iEns))*HaToeV,' eV'
! write(*,*) write(*,*)
! write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns)*HaToeV, ' eV'
! write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV'
! write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV'
! write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',(OmxDD(iEns)+OmcDD(iEns))*HaToeV,' eV' write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',(OmxDD(iEns)+OmcDD(iEns))*HaToeV,' eV'
! write(*,*) write(*,*)
! end do end do
! write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! write(*,*) write(*,*)
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') ' ENERGY DIFFERENCES FROM INDIVIDUAL ENERGIES ' write(*,'(A60)') ' ENERGY DIFFERENCES FROM INDIVIDUAL ENERGIES '
@ -191,9 +192,10 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux
do iEns=2,nEns 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(*,*)
write(*,'(A44, F16.10,A3)') ' H energy contribution : ',OmH(iEns), ' au'
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns), ' au' write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns), ' au'
write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(iEns), ' au' write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(iEns), ' au'
write(*,'(A44, F16.10,A3)') ' xc energy contribution : ',Omx(iEns)+Omc(iEns), ' au' write(*,'(A44, F16.10,A3)') ' Hxc energy contribution : ',OmH(iEns)+Omx(iEns)+Omc(iEns), ' au'
write(*,*) write(*,*)
write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns), ' au' write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns), ' au'
write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns), ' au' write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns), ' au'
@ -204,9 +206,10 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux
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(*,*)
write(*,'(A44, F16.10,A3)') ' H energy contribution : ',OmH(iEns)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(iEns)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' xc energy contribution : ',(Omx(iEns)+Omc(iEns))*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' Hxc energy contribution : ',(OmH(iEns)+Omx(iEns)+Omc(iEns))*HaToeV, ' eV'
write(*,*) write(*,*)
write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV'

3
src/eDFT/unrestricted_auxiliary_energy.f90
View File

@ -1,4 +1,4 @@
subroutine unrestricted_auxiliary_energy(nBas,nEns,eps,occnum,doNcentered,Eaux) subroutine unrestricted_auxiliary_energy(nBas,nEns,eps,occnum,Eaux)
! Compute the auxiliary KS energies ! Compute the auxiliary KS energies
@ -12,7 +12,6 @@ subroutine unrestricted_auxiliary_energy(nBas,nEns,eps,occnum,doNcentered,Eaux)
integer,intent(in) :: nEns integer,intent(in) :: nEns
double precision,intent(in) :: eps(nBas,nspin) double precision,intent(in) :: eps(nBas,nspin)
double precision,intent(in) :: occnum(nBas,nspin,nEns) double precision,intent(in) :: occnum(nBas,nspin,nEns)
logical,intent(in) :: doNcentered
! Local variables ! Local variables

3
src/eDFT/unrestricted_correlation_derivative_discontinuity.f90
View File

@ -1,4 +1,4 @@
subroutine unrestricted_correlation_derivative_discontinuity(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow,kappa,Ec) subroutine unrestricted_correlation_derivative_discontinuity(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow,Ec)
! Compute the correlation part of the derivative discontinuity ! Compute the correlation part of the derivative discontinuity
@ -15,7 +15,6 @@ subroutine unrestricted_correlation_derivative_discontinuity(rung,DFA,nEns,wEns,
double precision,intent(in) :: weight(nGrid) double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid,nspin) double precision,intent(in) :: rhow(nGrid,nspin)
double precision,intent(in) :: drhow(ncart,nGrid,nspin) double precision,intent(in) :: drhow(ncart,nGrid,nspin)
double precision,intent(in) :: kappa(nEns)
! Local variables ! Local variables

7
src/eDFT/unrestricted_correlation_individual_energy.f90
View File

@ -1,5 +1,5 @@
subroutine unrestricted_correlation_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nGrid,weight, & subroutine unrestricted_correlation_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nGrid,weight, &
rhow,drhow,rho,drho,doNcentered,LZc,Ec) rhow,drhow,rho,drho,LZc,Ec)
! Compute the correlation energy of individual states ! Compute the correlation energy of individual states
@ -19,7 +19,6 @@ subroutine unrestricted_correlation_individual_energy(rung,DFA,LDA_centered,nEns
double precision,intent(in) :: drhow(ncart,nGrid,nspin) double precision,intent(in) :: drhow(ncart,nGrid,nspin)
double precision,intent(in) :: rho(nGrid,nspin,nEns) double precision,intent(in) :: rho(nGrid,nspin,nEns)
double precision,intent(in) :: drho(ncart,nGrid,nspin,nEns) double precision,intent(in) :: drho(ncart,nGrid,nspin,nEns)
logical,intent(in) :: doNcentered
! Output variables ! Output variables
@ -38,7 +37,7 @@ subroutine unrestricted_correlation_individual_energy(rung,DFA,LDA_centered,nEns
case(1) case(1)
call unrestricted_lda_correlation_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,doNcentered,LZc,Ec) call unrestricted_lda_correlation_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,LZc,Ec)
! GGA functionals ! GGA functionals
@ -56,7 +55,7 @@ subroutine unrestricted_correlation_individual_energy(rung,DFA,LDA_centered,nEns
case(4) case(4)
call unrestricted_hybrid_correlation_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,doNcentered,LZc,Ec) call unrestricted_hybrid_correlation_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,LZc,Ec)
end select end select

7
src/eDFT/unrestricted_exchange_derivative_discontinuity.f90
View File

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

7
src/eDFT/unrestricted_exchange_individual_energy.f90
View File

@ -1,5 +1,5 @@
subroutine unrestricted_exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,nBas, & subroutine unrestricted_exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,nBas, &
ERI,Pw,rhow,drhow,P,rho,drho,Cx_choice,doNcentered,kappa,LZx,Ex) ERI,Pw,rhow,drhow,P,rho,drho,Cx_choice,doNcentered,LZx,Ex)
! Compute the exchange individual energy ! Compute the exchange individual energy
@ -27,11 +27,10 @@ subroutine unrestricted_exchange_individual_energy(rung,DFA,LDA_centered,nEns,wE
double precision,intent(in) :: drho(ncart,nGrid,nspin,nEns) double precision,intent(in) :: drho(ncart,nGrid,nspin,nEns)
integer,intent(in) :: Cx_choice integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Output variables ! Output variables
double precision,intent(out) :: LZx(nspin,nEns) double precision,intent(out) :: LZx(nspin)
double precision,intent(out) :: Ex(nspin,nEns) double precision,intent(out) :: Ex(nspin,nEns)
select case (rung) select case (rung)
@ -47,7 +46,7 @@ subroutine unrestricted_exchange_individual_energy(rung,DFA,LDA_centered,nEns,wE
case(1) case(1)
call unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,& call unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,&
rhow,rho,Cx_choice,doNcentered,kappa,LZx,Ex) rhow,rho,Cx_choice,doNcentered,LZx,Ex)
! GGA functionals ! GGA functionals

44
src/eDFT/unrestricted_hartree_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine unrestricted_hartree_individual_energy(nBas,nEns,Pw,P,ERI,doNcentered,kappa,LZH,EH) subroutine unrestricted_hartree_individual_energy(nBas,nEns,Pw,P,ERI,LZH,EH)
! Compute the hartree contribution to the individual energies ! Compute the hartree contribution to the individual energies
@ -12,8 +12,6 @@ subroutine unrestricted_hartree_individual_energy(nBas,nEns,Pw,P,ERI,doNcentered
double precision,intent(in) :: Pw(nBas,nBas,nspin) double precision,intent(in) :: Pw(nBas,nBas,nspin)
double precision,intent(in) :: P(nBas,nBas,nspin,nEns) double precision,intent(in) :: P(nBas,nBas,nspin,nEns)
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas) double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: kappa(nEns)
logical,intent(in) :: doNcentered
! Local variables ! Local variables
@ -26,47 +24,31 @@ subroutine unrestricted_hartree_individual_energy(nBas,nEns,Pw,P,ERI,doNcentered
! Output variables ! Output variables
double precision,intent(out) :: LZH(nsp,nEns) double precision,intent(out) :: LZH(nsp)
double precision,intent(out) :: EH(nsp,nEns) double precision,intent(out) :: EH(nsp,nEns)
! Compute HF exchange matrix ! Compute HF exchange matrix
allocate(J(nBas,nBas,nspin)) allocate(J(nBas,nBas,nspin))
LZH(:,:) = 0.d0 LZH(:) = 0.d0
EH(:,:) = 0.d0 EH(:,:) = 0.d0
do ispin=1,nspin do ispin=1,nspin
call unrestricted_hartree_potential(nBas,Pw(:,:,ispin),ERI,J(:,:,ispin)) call unrestricted_hartree_potential(nBas,Pw(:,:,ispin),ERI,J(:,:,ispin))
end do end do
LZH(1) = - 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,1)))
LZH(2) = - 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,2))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,1)))
LZH(3) = - 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,2)))
do iEns=1,nEns do iEns=1,nEns
! if(doNcentered) then EH(1,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,1)))
EH(2,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,2))) &
! LZH(1,iEns) = - 0.5d0*kappa(iEns)*kappa(iEns)*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,1))) + trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,1)))
! LZH(2,iEns) = - 0.5d0*kappa(iEns)*kappa(iEns)*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,2))) & EH(3,iEns) = trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,2)))
! - 0.5d0*kappa(iEns)*kappa(iEns)*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,1)))
! LZH(3,iEns) = - 0.5d0*kappa(iEns)*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,2)))
! EH(1,iEns) = kappa(iEns)*trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,1)))
! EH(2,iEns) = kappa(iEns)*trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,2))) &
! + kappa(iEns)*trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,1)))
! EH(3,iEns) = kappa(iEns)*trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,2)))
! else
LZH(1,iEns) = - 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,1)))
LZH(2,iEns) = - 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,2))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,1)))
LZH(3,iEns) = - 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,2)))
EH(1,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,1)))
EH(2,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,2))) &
+ trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,1)))
EH(3,iEns) = trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,2)))
! endif
end do end do

3
src/eDFT/unrestricted_hybrid_correlation_individual_energy.f90
View File

@ -1,5 +1,5 @@
subroutine unrestricted_hybrid_correlation_individual_energy(DFA,nEns,wEns,nGrid,weight, & subroutine unrestricted_hybrid_correlation_individual_energy(DFA,nEns,wEns,nGrid,weight, &
rhow,drhow,rho,drho,doNcentered,LZc,Ec) rhow,drhow,rho,drho,LZc,Ec)
! Compute the hybrid correlation energy for individual states ! Compute the hybrid correlation energy for individual states
@ -17,7 +17,6 @@ subroutine unrestricted_hybrid_correlation_individual_energy(DFA,nEns,wEns,nGrid
double precision,intent(in) :: drhow(ncart,nGrid) double precision,intent(in) :: drhow(ncart,nGrid)
double precision,intent(in) :: rho(nGrid,nEns) double precision,intent(in) :: rho(nGrid,nEns)
double precision,intent(in) :: drho(ncart,nGrid,nEns) double precision,intent(in) :: drho(ncart,nGrid,nEns)
logical,intent(in) :: doNcentered
! Output variables ! Output variables

3
src/eDFT/unrestricted_hybrid_exchange_derivative_discontinuity.f90
View File

@ -1,5 +1,5 @@
subroutine unrestricted_hybrid_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC,aCC,nGrid,weight,rhow,& subroutine unrestricted_hybrid_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC,aCC,nGrid,weight,rhow,&
Cx_choice,doNcentered,kappa,ExDD) Cx_choice,doNcentered,ExDD)
! Compute the exchange part of the derivative discontinuity for hybrid functionals ! Compute the exchange part of the derivative discontinuity for hybrid functionals
@ -19,7 +19,6 @@ subroutine unrestricted_hybrid_exchange_derivative_discontinuity(DFA,nEns,wEns,n
double precision,intent(in) :: rhow(nGrid) double precision,intent(in) :: rhow(nGrid)
integer,intent(in) :: Cx_choice integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Local variables ! Local variables

94
src/eDFT/unrestricted_individual_energy.f90
View File

@ -52,8 +52,8 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
double precision :: EH(nsp,nEns) double precision :: EH(nsp,nEns)
double precision :: Ex(nspin,nEns) double precision :: Ex(nspin,nEns)
double precision :: Ec(nsp,nEns) double precision :: Ec(nsp,nEns)
double precision :: LZH(nsp,nEns) double precision :: LZH(nsp)
double precision :: LZx(nspin,nEns) double precision :: LZx(nspin)
double precision :: LZc(nsp) double precision :: LZc(nsp)
double precision :: Eaux(nspin,nEns) double precision :: Eaux(nspin,nEns)
@ -116,18 +116,18 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
! Individual Hartree energy ! Individual Hartree energy
!------------------------------------------------------------------------ !------------------------------------------------------------------------
LZH(:,:) = 0d0 LZH(:) = 0d0
EH(:,:) = 0d0 EH(:,:) = 0d0
call unrestricted_hartree_individual_energy(nBas,nEns,Pw,P,ERI,doNcentered,kappa,LZH,EH) call unrestricted_hartree_individual_energy(nBas,nEns,Pw,P,ERI,LZH,EH)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Individual exchange energy ! Individual exchange energy
!------------------------------------------------------------------------ !------------------------------------------------------------------------
LZx(:,:) = 0d0 LZx(:) = 0d0
Ex(:,:) = 0d0 Ex(:,:) = 0d0
call unrestricted_exchange_individual_energy(x_rung,x_DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,nBas,ERI, & call unrestricted_exchange_individual_energy(x_rung,x_DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,nBas,ERI, &
Pw,rhow,drhow,P,rho,drho,Cx_choice,doNcentered,kappa,LZx,Ex) Pw,rhow,drhow,P,rho,drho,Cx_choice,doNcentered,LZx,Ex)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Individual correlation energy ! Individual correlation energy
@ -136,13 +136,13 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
LZc(:) = 0d0 LZc(:) = 0d0
Ec(:,:) = 0d0 Ec(:,:) = 0d0
call unrestricted_correlation_individual_energy(c_rung,c_DFA,LDA_centered,nEns,wEns,nGrid,weight, & call unrestricted_correlation_individual_energy(c_rung,c_DFA,LDA_centered,nEns,wEns,nGrid,weight, &
rhow,drhow,rho,drho,doNcentered,LZc,Ec) rhow,drhow,rho,drho,LZc,Ec)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Compute auxiliary energies ! Compute auxiliary energies
!------------------------------------------------------------------------ !------------------------------------------------------------------------
call unrestricted_auxiliary_energy(nBas,nEns,eKS,occnum,doNcentered,Eaux) call unrestricted_auxiliary_energy(nBas,nEns,eKS,occnum,Eaux)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Compute derivative discontinuities ! Compute derivative discontinuities
@ -150,64 +150,84 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
do ispin=1,nspin do ispin=1,nspin
call unrestricted_exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,nCC,aCC,nGrid,weight, & call unrestricted_exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,nCC,aCC,nGrid,weight, &
rhow(:,ispin),drhow(:,:,ispin),Cx_choice,doNcentered,kappa,ExDD(ispin,:)) rhow(:,ispin),drhow(:,:,ispin),Cx_choice,doNcentered,ExDD(ispin,:))
end do end do
call unrestricted_correlation_derivative_discontinuity(c_rung,c_DFA,nEns,wEns,nGrid,weight,rhow,drhow,kappa,EcDD) call unrestricted_correlation_derivative_discontinuity(c_rung,c_DFA,nEns,wEns,nGrid,weight,rhow,drhow,EcDD)
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Total energy ! Total energy
!------------------------------------------------------------------------ !------------------------------------------------------------------------
if(doNcentered) then
do iEns=1,nEns
E(iEns) = sum(Eaux(:,iEns)) &
+ kappa(iEns)*(sum(LZH(:)) + sum(LZx(:)) + sum(LZc(:))) &
+ sum(ExDD(:,iEns)) + sum(EcDD(:,iEns))
end do
print*,E
else
do iEns=1,nEns
E(iEns) = sum(Eaux(:,iEns)) &
+ sum(LZH(:)) + sum(LZx(:)) + sum(LZc(:)) &
+ sum(ExDD(:,iEns)) + sum(EcDD(:,iEns))
end do
end if
do iEns=1,nEns do iEns=1,nEns
E(iEns) = sum(Eaux(:,iEns)) & E(iEns) = sum(ET(:,iEns)) + sum(EV(:,iEns)) &
+ sum(LZH(:,iEns)) + sum(LZx(:,iEns)) + sum(LZc(:)) & + sum(EH(:,iEns)) + sum(Ex(:,iEns)) + sum(Ec(:,iEns)) &
+ sum(ExDD(:,iEns)) + sum(EcDD(:,iEns)) + kappa(iEns)*(sum(LZH(:)) + sum(LZx(:)) + sum(LZc(:))) &
+ sum(ExDD(:,iEns)) + sum(EcDD(:,iEns))
end do end do
!E(2) = (1.d0/2.d0)*(E(2)) print*,E
! print*,'test shift =',(1.d0+wEns(2)/2.d0)*(sum(Eaux(:,2)) + sum(LZH(:,2)) &
! + sum(LZx(:,2)) + sum(LZc(:)))+(1.d0-kappa(2)*wEns(2))*(sum(Eaux(:,1))&
! +sum(LZH(:,1)) + sum(LZx(:,1)) + sum(LZc(:)))
! print*,'test=',(1.d0+wEns(2)/2.d0)*(sum(Eaux(:,2)))+(1.d0-kappa(2)*wEns(2))*(sum(Eaux(:,1))&
! )
! print*, 'ensemble energy=',Ew
! Alternative way of calculating individual energies
! do iEns=1,nEns
! E(iEns) = sum(ET(:,iEns)) + sum(EV(:,iEns)) &
! + sum(EH(:,iEns)) + sum(Ex(:,iEns)) + sum(Ec(:,iEns)) &
! + sum(LZH(:)) + sum(LZx(:)) + sum(LZc(:)) &
! + sum(ExDD(:,iEns)) + sum(EcDD(:,iEns))
! end do
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Excitation energies ! Excitation energies
!------------------------------------------------------------------------ !------------------------------------------------------------------------
do iEns=1,nEns do iEns=1,nEns
Om(iEns) = E(iEns) - E(1) Om(iEns) = E(iEns) - E(1)
OmH(iEns) = sum(EH(:,iEns)) - sum(EH(:,1)) OmH(iEns) = sum(EH(:,iEns)) - sum(EH(:,1))
Omx(iEns) = sum(Ex(:,iEns)) - sum(Ex(:,1)) Omx(iEns) = sum(Ex(:,iEns)) - sum(Ex(:,1))
Omc(iEns) = sum(Ec(:,iEns)) - sum(Ec(:,1)) Omc(iEns) = sum(Ec(:,iEns)) - sum(Ec(:,1))
Omaux(iEns) = sum(Eaux(:,iEns)) - sum(Eaux(:,1)) Omaux(iEns) = sum(Eaux(:,iEns)) - sum(Eaux(:,1))
OmxDD(iEns) = sum(ExDD(:,iEns)) - sum(ExDD(:,1)) OmxDD(iEns) = sum(ExDD(:,iEns)) - sum(ExDD(:,1))
OmcDD(iEns) = sum(EcDD(:,iEns)) - sum(EcDD(:,1)) OmcDD(iEns) = sum(EcDD(:,iEns)) - sum(EcDD(:,1))
end do end do
if(doNcentered) then
do iEns=1,nEns
OmH(iEns) = OmH(iEns) + (kappa(iEns) - kappa(1))*sum(LZH(:))
Omx(iEns) = Omx(iEns) + (kappa(iEns) - kappa(1))*sum(LZx(:))
Omc(iEns) = Omc(iEns) + (kappa(iEns) - kappa(1))*sum(LZc(:))
Omaux(iEns) = Omaux(iEns) &
+ (kappa(iEns) - kappa(1))*(sum(LZH(:)) + sum(LZx(:)) + sum(LZc(:)))
OmxDD(iEns) = kappa(iEns)*sum(ExDD(:,iEns)) - kappa(1)*sum(ExDD(:,1))
OmcDD(iEns) = kappa(iEns)*sum(EcDD(:,iEns)) - kappa(1)*sum(EcDD(:,1))
end do
end if
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! Dump results ! Dump results
!------------------------------------------------------------------------ !------------------------------------------------------------------------
call print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux,ExDD,EcDD,E,Om,Omx,Omc,Omaux,OmxDD,OmcDD) call print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux,ExDD,EcDD,E,Om,OmH,Omx,Omc,Omaux,OmxDD,OmcDD)
end subroutine unrestricted_individual_energy end subroutine unrestricted_individual_energy

11
src/eDFT/unrestricted_lda_correlation_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine unrestricted_lda_correlation_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,doNcentered,LZc,Ec) subroutine unrestricted_lda_correlation_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,LZc,Ec)
! Compute LDA correlation energy for individual states ! Compute LDA correlation energy for individual states
@ -15,7 +15,6 @@ subroutine unrestricted_lda_correlation_individual_energy(DFA,LDA_centered,nEns,
double precision,intent(in) :: weight(nGrid) double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid,nspin) double precision,intent(in) :: rhow(nGrid,nspin)
double precision,intent(in) :: rho(nGrid,nspin,nEns) double precision,intent(in) :: rho(nGrid,nspin,nEns)
logical,intent(in) :: doNcentered
! Output variables ! Output variables
@ -28,19 +27,19 @@ subroutine unrestricted_lda_correlation_individual_energy(DFA,LDA_centered,nEns,
case (1) case (1)
! call UW38_lda_correlation_individual_energy(nGrid,weight,rhow,rho,doNcentered,kappa,Ec) ! call UW38_lda_correlation_individual_energy(nGrid,weight,rhow,rho,LZc,Ec)
case (2) case (2)
! call UPW92_lda_correlation_individual_energy(nGrid,weight,rhow,rho,doNcentered,kappa,Ec) ! call UPW92_lda_correlation_individual_energy(nGrid,weight,rhow,rho,LZc,Ec)
case (3) case (3)
! call UVWN3_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,doNcentered,LZc,Ec) ! call UVWN3_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,LZc,Ec)
case (4) case (4)
call UVWN5_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,doNcentered,LZc,Ec) call UVWN5_lda_correlation_individual_energy(nEns,nGrid,weight,rhow,rho,LZc,Ec)
case default case default

7
src/eDFT/unrestricted_lda_exchange_derivative_discontinuity.f90
View File

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

9
src/eDFT/unrestricted_lda_exchange_individual_energy.f90
View File

@ -1,5 +1,5 @@
subroutine unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,rhow,& subroutine unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nCC,aCC,nGrid,weight,rhow,&
rho,Cx_choice,doNcentered,kappa,LZx,Ex) rho,Cx_choice,doNcentered,LZx,Ex)
! Compute LDA exchange energy for individual states ! Compute LDA exchange energy for individual states
@ -20,11 +20,10 @@ subroutine unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEn
double precision,intent(in) :: rho(nGrid,nspin,nEns) double precision,intent(in) :: rho(nGrid,nspin,nEns)
integer,intent(in) :: Cx_choice integer,intent(in) :: Cx_choice
logical,intent(in) :: doNcentered logical,intent(in) :: doNcentered
double precision,intent(in) :: kappa(nEns)
! Output variables ! Output variables
double precision :: LZx(nspin,nEns) double precision :: LZx(nspin)
double precision :: Ex(nspin,nEns) double precision :: Ex(nspin,nEns)
! Select correlation functional ! Select correlation functional
@ -33,12 +32,12 @@ subroutine unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEn
case (1) case (1)
call US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,doNcentered,kappa,LZx,Ex) call US51_lda_exchange_individual_energy(nEns,nGrid,weight,rhow,rho,LZx,Ex)
case (2) case (2)
call UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rhow,rho, & call UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rhow,rho, &
Cx_choice,doNcentered,kappa,LZx,Ex) Cx_choice,doNcentered,LZx,Ex)
case default case default