LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2024-12-22 12:23:50 +01:00
Code Releases Activity

Compare commits

base: LCPQ:1e428400c461c3eafa7d780bac6f5e66602852a5
Branches Tags
LCPQ:master
LCPQ:refact
LCPQ:trexio
LCPQ:old
LCPQ:add-license-1
LCPQ:0.9
...
compare: LCPQ:cd781a1e11332700fe790cad04ba56ffd264d861
Branches Tags
LCPQ:master
LCPQ:refact
LCPQ:trexio
LCPQ:old
LCPQ:add-license-1
LCPQ:0.9

2 Commits
1e428400c4 ... cd781a1e11

Author SHA1 Message Date
Pierre-Francois Loos
cd781a1e11 HF works 2020-07-06 21:27:28 +02:00
Pierre-Francois Loos
b70c4e1f7a more structure in the code 2020-07-06 20:57:27 +02:00
41 changed files with 939 additions and 162 deletions
Show Stats Expand all files Collapse all files

4
input/dft
View File

@ -6,7 +6,7 @@
# GGA = 2: RB88
# Hybrid = 4
# Hartree-Fock = 666
1 US51
666 HF
# correlation rung:
# Hartree = 0
# LDA = 1: RVWN5,RMFL20
@ -19,7 +19,7 @@
# Number of states in ensemble (nEns)
3
# Ensemble weights: wEns(1),...,wEns(nEns-1)
1 0.0
0.0 0.0
# Parameters for CC weight-dependent exchange functional
0.000000 0.0000000 0.000000
0.000000 0.0000000 0.0000000

6
src/eDFT/GOK_RKS.f90
View File

@ -236,8 +236,8 @@ subroutine GOK_RKS(restart,x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,aCC_
! Compute exchange potential
call exchange_potential(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas,Pw(:,:),ERI(:,:,:,:), &
AO(:,:),dAO(:,:,:),rhow(:),drhow(:,:),Fx(:,:),FxHF(:,:))
call restricted_exchange_potential(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas,Pw(:,:), &
ERI(:,:,:,:),AO(:,:),dAO(:,:,:),rhow(:),drhow(:,:),Fx(:,:),FxHF(:,:))
! Compute correlation potential
@ -294,7 +294,7 @@ subroutine GOK_RKS(restart,x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,aCC_
! Exchange energy
call exchange_energy(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
call restricted_exchange_energy(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
Pw(:,:),FxHF(:,:),rhow(:),drhow(:,:),Ex)
! Correlation energy

7
src/eDFT/GOK_UKS.f90
View File

@ -248,8 +248,9 @@ subroutine GOK_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGrid,weight,aCC_w1,aCC_w
! Compute exchange potential
do ispin=1,nspin
call exchange_potential(x_rung,x_DFA,nEns,wEns(:),nGrid,weight(:),aCC_w1,aCC_w2,nBas,Pw(:,:,ispin),ERI(:,:,:,:), &
AO(:,:),dAO(:,:,:),rhow(:,ispin),drhow(:,:,ispin),Fx(:,:,ispin),FxHF(:,:,ispin))
call unrestricted_exchange_potential(x_rung,x_DFA,nEns,wEns(:),nGrid,weight(:),aCC_w1,aCC_w2,nBas,Pw(:,:,ispin), &
ERI(:,:,:,:),AO(:,:),dAO(:,:,:),rhow(:,ispin),drhow(:,:,ispin),Fx(:,:,ispin), &
FxHF(:,:,ispin))
end do
! Compute correlation potential
@ -306,7 +307,7 @@ subroutine GOK_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGrid,weight,aCC_w1,aCC_w
! Exchange energy
do ispin=1,nspin
call exchange_energy(x_rung,x_DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
call unrestricted_exchange_energy(x_rung,x_DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
Pw(:,:,ispin),FxHF(:,:,ispin),rhow(:,ispin),drhow(:,:,ispin),Ex(ispin))
end do

4
src/eDFT/B88_gga_exchange_energy.f90 → src/eDFT/UB88_gga_exchange_energy.f90
View File

@ -1,4 +1,4 @@
subroutine B88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
subroutine UB88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
! Compute Becke's 88 GGA exchange energy
@ -47,4 +47,4 @@ subroutine B88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
end do
end subroutine B88_gga_exchange_energy
end subroutine UB88_gga_exchange_energy

4
src/eDFT/B88_gga_exchange_potential.f90 → src/eDFT/UB88_gga_exchange_potential.f90
View File

@ -1,4 +1,4 @@
subroutine B88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
subroutine UB88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
! Compute Becke's GGA exchange potential
@ -62,4 +62,4 @@ subroutine B88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
end do
end do
end subroutine B88_gga_exchange_potential
end subroutine UB88_gga_exchange_potential

4
src/eDFT/G96_gga_exchange_energy.f90 → src/eDFT/UG96_gga_exchange_energy.f90
View File

@ -1,4 +1,4 @@
subroutine G96_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
subroutine UG96_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
! Compute Gill's 96 GGA exchange energy
@ -46,4 +46,4 @@ subroutine G96_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
end do
end subroutine G96_gga_exchange_energy
end subroutine UG96_gga_exchange_energy

4
src/eDFT/G96_gga_exchange_potential.f90 → src/eDFT/UG96_gga_exchange_potential.f90
View File

@ -1,4 +1,4 @@
subroutine G96_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
subroutine UG96_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
! Compute Gill's GGA exchange poential
@ -63,4 +63,4 @@ subroutine G96_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
enddo
enddo
end subroutine G96_gga_exchange_potential
end subroutine UG96_gga_exchange_potential

4
src/eDFT/eDFT_UKS.f90
View File

@ -238,7 +238,7 @@ subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weig
! Compute exchange potential
do ispin=1,nspin
call exchange_potential(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
call unrestricted_exchange_potential(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
Pw(:,:,ispin),ERI(:,:,:,:),AO(:,:),dAO(:,:,:),rhow(:,ispin),drhow(:,:,ispin), &
Fx(:,:,ispin),FxHF(:,:,ispin))
end do
@ -317,7 +317,7 @@ subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weig
! Exchange energy
do ispin=1,nspin
call exchange_energy(x_rung,x_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, &
call unrestricted_exchange_energy(x_rung,x_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, &
Pw(:,:,ispin),FxHF(:,:,ispin),rhow(:,ispin),drhow(:,:,ispin),Ex(ispin))
end do

2
src/eDFT/huckel_guess.f90
View File

@ -36,7 +36,7 @@ subroutine huckel_guess(nBas,S,Hc,ERI,J,K,X,cp,F,Fp,e,c,P)
c(:,:) = matmul(X(:,:),cp(:,:))
call hartree_coulomb(nBas,P,ERI,J)
call fock_exchange_potential(nBas,P,ERI,K)
call restricted_fock_exchange_potential(nBas,P,ERI,K)
F(:,:) = Hc(:,:) + J(:,:) + 0.5d0*K(:,:)

8
src/eDFT/exchange_derivative_discontinuity.f90 → src/eDFT/restricted_exchange_derivative_discontinuity.f90
View File

@ -1,4 +1,4 @@
subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,drhow,ExDD)
subroutine restricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,drhow,ExDD)
! Compute the exchange part of the derivative discontinuity
@ -37,13 +37,13 @@ subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,aCC_w1,aCC_w2,nG
case(1)
call lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:))
call restricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:))
! GGA functionals
case(2)
call gga_exchange_derivative_discontinuity(DFA,nEns,wEns(:),nGrid,weight(:),rhow(:),drhow(:,:),ExDD(:))
call restricted_gga_exchange_derivative_discontinuity(DFA,nEns,wEns(:),nGrid,weight(:),rhow(:),drhow(:,:),ExDD(:))
! Hybrid functionals
@ -60,4 +60,4 @@ subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,aCC_w1,aCC_w2,nG
end select
end subroutine exchange_derivative_discontinuity
end subroutine restricted_exchange_derivative_discontinuity

16
src/eDFT/exchange_energy.f90 → src/eDFT/restricted_exchange_energy.f90
View File

@ -1,4 +1,4 @@
subroutine exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,P,FxHF,rho,drho,Ex)
subroutine restricted_exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,P,FxHF,rho,drho,Ex)
! Compute the exchange energy
@ -43,7 +43,7 @@ subroutine exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,w
case(1)
call lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,ExLDA)
call restricted_lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,ExLDA)
Ex = ExLDA
@ -51,7 +51,7 @@ subroutine exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,w
case(2)
call gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,ExGGA)
call restricted_gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,ExGGA)
Ex = ExGGA
@ -63,9 +63,9 @@ subroutine exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,w
aX = 0.72d0
aC = 0.81d0
call lda_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,ExLDA)
call gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,ExGGA)
call fock_exchange_energy(nBas,P,FxHF,ExHF)
call restricted_lda_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,ExLDA)
call restricted_gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,ExGGA)
call restricted_fock_exchange_energy(nBas,P,FxHF,ExHF)
Ex = ExLDA &
+ cX*(ExHF - ExLDA) &
@ -75,10 +75,10 @@ subroutine exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,w
case(666)
call fock_exchange_energy(nBas,P,FxHF,ExHF)
call restricted_fock_exchange_energy(nBas,P,FxHF,ExHF)
Ex = ExHF
end select
end subroutine exchange_energy
end subroutine restricted_exchange_energy

10
src/eDFT/exchange_individual_energy.f90 → src/eDFT/restricted_exchange_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, &
subroutine restricted_exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, &
ERI,Pw,P,rhow,drhow,rho,drho,Ex)
! Compute the exchange individual energy
@ -48,7 +48,7 @@ subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC
case(1)
call lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,ExLDA)
call restricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,ExLDA)
Ex = ExLDA
@ -56,7 +56,7 @@ subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC
case(2)
call gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,ExGGA)
call restricted_gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,ExGGA)
Ex = ExGGA
@ -71,10 +71,10 @@ subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC
case(666)
call fock_exchange_individual_energy(nBas,Pw,P,ERI,ExHF)
call restricted_fock_exchange_individual_energy(nBas,Pw,P,ERI,ExHF)
Ex = ExHF
end select
end subroutine exchange_individual_energy
end subroutine restricted_exchange_individual_energy

17
src/eDFT/exchange_potential.f90 → src/eDFT/restricted_exchange_potential.f90
View File

@ -1,4 +1,5 @@
subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,P,ERI,AO,dAO,rho,drho,Fx,FxHF)
subroutine restricted_exchange_potential(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,P, &
ERI,AO,dAO,rho,drho,Fx,FxHF)
! Compute the exchange potential
@ -47,13 +48,13 @@ subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGri
case(1)
call lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
call restricted_lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
! GGA functionals
case(2)
call gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
call restricted_gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
! Hybrid functionals
@ -64,9 +65,9 @@ subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGri
cX = 0.20d0
aX = 0.72d0
call lda_exchange_potential(DFA,nGrid,weight,nBas,AO,rho,FxLDA)
call gga_exchange_potential(DFA,nGrid,weight,nBas,AO,dAO,rho,drho,FxGGA)
call fock_exchange_potential(nBas,P,ERI,FxHF)
call restricted_lda_exchange_potential(DFA,nGrid,weight,nBas,AO,rho,FxLDA)
call restricted_gga_exchange_potential(DFA,nGrid,weight,nBas,AO,dAO,rho,drho,FxGGA)
call restricted_fock_exchange_potential(nBas,P,ERI,FxHF)
Fx(:,:) = FxLDA(:,:) &
+ cX*(FxHF(:,:) - FxLDA(:,:)) &
@ -76,10 +77,10 @@ subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGri
case(666)
call fock_exchange_potential(nBas,P,ERI,FxHF)
call restricted_fock_exchange_potential(nBas,P,ERI,FxHF)
Fx(:,:) = FxHF(:,:)
end select
end subroutine exchange_potential
end subroutine restricted_exchange_potential

25
src/eDFT/restricted_fock_exchange_energy.f90 Normal file
View File

@ -0,0 +1,25 @@
subroutine restricted_fock_exchange_energy(nBas,P,Fx,Ex)
! Compute the (exact) Fock exchange energy
implicit none
! Input variables
integer,intent(in) :: nBas
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: Fx(nBas,nBas)
! Local variables
double precision,external :: trace_matrix
! Output variables
double precision,intent(out) :: Ex
! Compute HF exchange energy
Ex = 0.5d0*trace_matrix(nBas,matmul(P,Fx))
end subroutine restricted_fock_exchange_energy

6
src/eDFT/fock_exchange_individual_energy.f90 → src/eDFT/restricted_fock_exchange_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine fock_exchange_individual_energy(nBas,Pw,P,ERI,Ex)
subroutine restricted_fock_exchange_individual_energy(nBas,Pw,P,ERI,Ex)
! Compute the Fock exchange potential
@ -24,8 +24,8 @@ subroutine fock_exchange_individual_energy(nBas,Pw,P,ERI,Ex)
allocate(Fx(nBas,nBas))
call fock_exchange_potential(nBas,Pw(:,:),ERI(:,:,:,:),Fx(:,:))
call restricted_fock_exchange_potential(nBas,Pw(:,:),ERI(:,:,:,:),Fx(:,:))
Ex = trace_matrix(nBas,matmul(P(:,:),Fx(:,:))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:),Fx(:,:)))
end subroutine fock_exchange_individual_energy
end subroutine restricted_fock_exchange_individual_energy

4
src/eDFT/fock_exchange_potential.f90 → src/eDFT/restricted_fock_exchange_potential.f90
View File

@ -1,4 +1,4 @@
subroutine fock_exchange_potential(nBas,P,ERI,Fx)
subroutine restricted_fock_exchange_potential(nBas,P,ERI,Fx)
! Compute the Fock exchange potential
@ -33,4 +33,4 @@ subroutine fock_exchange_potential(nBas,P,ERI,Fx)
Fx(:,:) = 0.5d0*Fx(:,:)
end subroutine fock_exchange_potential
end subroutine restricted_fock_exchange_potential

6
src/eDFT/gga_exchange_derivative_discontinuity.f90 → src/eDFT/restricted_gga_exchange_derivative_discontinuity.f90
View File

@ -1,4 +1,4 @@
subroutine gga_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow,drhow,ExDD)
subroutine restricted_gga_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow,drhow,ExDD)
! Compute the exchange GGA part of the derivative discontinuity
@ -26,7 +26,7 @@ subroutine gga_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow
select case (DFA)
case ('RB88')
case ('B88')
ExDD(:) = 0d0
@ -37,4 +37,4 @@ subroutine gga_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow
end select
end subroutine gga_exchange_derivative_discontinuity
end subroutine restricted_gga_exchange_derivative_discontinuity

14
src/eDFT/gga_exchange_energy.f90 → src/eDFT/restricted_gga_exchange_energy.f90
View File

@ -1,4 +1,4 @@
subroutine gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,Ex)
subroutine restricted_gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,Ex)
! Select GGA exchange functional for energy calculation
@ -22,17 +22,9 @@ subroutine gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,Ex)
select case (DFA)
case ('G96')
call G96_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
case ('RB88')
call RB88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
case ('B88')
call B88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
call RB88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
case default
@ -41,4 +33,4 @@ subroutine gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,Ex)
end select
end subroutine gga_exchange_energy
end subroutine restricted_gga_exchange_energy

6
src/eDFT/gga_exchange_individual_energy.f90 → src/eDFT/restricted_gga_exchange_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,Ex)
subroutine restricted_gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,Ex)
! Compute GGA exchange energy for individual states
@ -25,7 +25,7 @@ subroutine gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,
select case (DFA)
case ('RB88')
case ('B88')
call RB88_gga_exchange_individual_energy(nGrid,weight(:),rhow(:),drhow(:,:),rho(:),drho(:,:),Ex)
@ -36,4 +36,4 @@ subroutine gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,
end select
end subroutine gga_exchange_individual_energy
end subroutine restricted_gga_exchange_individual_energy

14
src/eDFT/gga_exchange_potential.f90 → src/eDFT/restricted_gga_exchange_potential.f90
View File

@ -1,4 +1,4 @@
subroutine gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
subroutine restricted_gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
! Select GGA exchange functional for potential calculation
@ -26,17 +26,9 @@ subroutine gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drh
select case (DFA)
case ('G96')
call G96_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
case ('RB88')
call RB88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
case ('B88')
call B88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
call RB88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
case default
@ -45,4 +37,4 @@ subroutine gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drh
end select
end subroutine gga_exchange_potential
end subroutine restricted_gga_exchange_potential

7
src/eDFT/restricted_individual_energy.f90
View File

@ -91,8 +91,8 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,n
!------------------------------------------------------------------------
do iEns=1,nEns
call exchange_individual_energy(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas,ERI(:,:,:,:), &
Pw(:,:),P(:,:,iEns),rhow(:),drhow(:,:),rho(:,iEns),drho(:,:,iEns),Ex(iEns))
call restricted_exchange_individual_energy(x_rung,x_DFA,LDA_centered,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),nBas, &
ERI(:,:,:,:),Pw(:,:),P(:,:,iEns),rhow(:),drhow(:,:),rho(:,iEns),drho(:,:,iEns),Ex(iEns))
end do
!------------------------------------------------------------------------
@ -114,7 +114,8 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,n
! Compute derivative discontinuities
!------------------------------------------------------------------------
call exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),drhow(:,:),ExDD(:))
call restricted_exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:), &
rhow(:),drhow(:,:),ExDD(:))
call restricted_correlation_derivative_discontinuity(c_rung,c_DFA,nEns,wEns(:),nGrid,weight(:),rhow(:),drhow(:,:),EcDD(:))

18
src/eDFT/lda_exchange_derivative_discontinuity.f90 → src/eDFT/restricted_lda_exchange_derivative_discontinuity.f90
View File

@ -1,4 +1,4 @@
subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,ExDD)
subroutine restricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,ExDD)
! Compute the exchange LDA part of the derivative discontinuity
@ -28,26 +28,18 @@ subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,aCC_w1,aCC_w2,nGr
select case (DFA)
case ('US51')
case ('S51')
ExDD(:) = 0d0
case ('RS51')
ExDD(:) = 0d0
case ('RMFL20')
case ('MFL20')
call RMFL20_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight(:),rhow(:),ExDD(:))
case ('RCC')
case ('CC')
call RCC_lda_exchange_derivative_discontinuity(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:))
case ('UCC')
call UCC_lda_exchange_derivative_discontinuity(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:))
case default
call print_warning('!!! LDA exchange derivative discontinuity not available !!!')
@ -55,4 +47,4 @@ subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,aCC_w1,aCC_w2,nGr
end select
end subroutine lda_exchange_derivative_discontinuity
end subroutine restricted_lda_exchange_derivative_discontinuity

20
src/eDFT/lda_exchange_energy.f90 → src/eDFT/restricted_lda_exchange_energy.f90
View File

@ -1,4 +1,4 @@
subroutine lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,Ex)
subroutine restricted_lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,Ex)
! Select LDA exchange functional
@ -25,31 +25,23 @@ subroutine lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,we
select case (DFA)
case ('US51')
call US51_lda_exchange_energy(nGrid,weight,rho,Ex)
case ('RS51')
case ('S51')
call RS51_lda_exchange_energy(nGrid,weight,rho,Ex)
case ('RMFL20')
case ('MFL20')
call RMFL20_lda_exchange_energy(LDA_centered,nEns,wEns,nGrid,weight,rho,Ex)
case ('RCC')
case ('CC')
call RCC_lda_exchange_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,Ex)
case ('UCC')
call UCC_lda_exchange_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,Ex)
case default
call print_warning('!!! LDA exchange functional not available !!!')
call print_warning('!!! LDA restricted exchange functional not available !!!')
stop
end select
end subroutine lda_exchange_energy
end subroutine restricted_lda_exchange_energy

18
src/eDFT/lda_exchange_individual_energy.f90 → src/eDFT/restricted_lda_exchange_individual_energy.f90
View File

@ -1,4 +1,4 @@
subroutine lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex)
subroutine restricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex)
! Compute LDA exchange energy for individual states
@ -26,26 +26,18 @@ subroutine lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_
select case (DFA)
case ('US51')
call US51_lda_exchange_individual_energy(nGrid,weight,rhow,rho,Ex)
case ('RS51')
case ('S51')
call RS51_lda_exchange_individual_energy(nGrid,weight,rhow,rho,Ex)
case ('RMFL20')
case ('MFL20')
call RMFL20_lda_exchange_individual_energy(LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,Ex)
case ('RCC')
case ('CC')
call RCC_lda_exchange_individual_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex)
case ('UCC')
call UCC_lda_exchange_individual_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex)
case default
call print_warning('!!! LDA exchange individual energy not available !!!')
@ -53,4 +45,4 @@ subroutine lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_
end select
end subroutine lda_exchange_individual_energy
end subroutine restricted_lda_exchange_individual_energy

50
src/eDFT/restricted_lda_exchange_potential.f90 Normal file
View File

@ -0,0 +1,50 @@
subroutine restricted_lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
! Select LDA correlation potential
implicit none
include 'parameters.h'
! Input variables
logical,intent(in) :: LDA_centered
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: AO(nBas,nGrid)
double precision,intent(in) :: rho(nGrid)
! Output variables
double precision,intent(out) :: Fx(nBas,nBas)
! Select exchange functional
select case (DFA)
case ('S51')
call RS51_lda_exchange_potential(nGrid,weight,nBas,AO,rho,Fx)
case ('MFL20')
call RMFL20_lda_exchange_potential(LDA_centered,nEns,wEns,nGrid,weight,nBas,AO,rho,Fx)
case ('CC')
call RCC_lda_exchange_potential(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
case default
call print_warning('!!! LDA exchange functional not available !!!')
stop
end select
end subroutine restricted_lda_exchange_potential

63
src/eDFT/unrestricted_exchange_derivative_discontinuity.f90 Normal file
View File

@ -0,0 +1,63 @@
subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,drhow,ExDD)
! Compute the exchange part of the derivative discontinuity
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: rung
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid)
double precision,intent(in) :: drhow(ncart,nGrid)
! Local variables
! Output variables
double precision,intent(out) :: ExDD(nEns)
select case (rung)
! Hartree calculation
case(0)
ExDD(:) = 0d0
! LDA functionals
case(1)
call unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:))
! GGA functionals
case(2)
call unrestricted_gga_exchange_derivative_discontinuity(DFA,nEns,wEns(:),nGrid,weight(:),rhow(:),drhow(:,:),ExDD(:))
! Hybrid functionals
case(4)
call print_warning('!!! exchange part of derivative discontinuity NYI for hybrids !!!')
stop
! Hartree-Fock calculation
case(666)
ExDD(:) = 0d0
end select
end subroutine unrestricted_exchange_derivative_discontinuity

84
src/eDFT/unrestricted_exchange_energy.f90 Normal file
View File

@ -0,0 +1,84 @@
subroutine unrestricted_exchange_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,P,FxHF,rho,drho,Ex)
! Compute the exchange energy
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: rung
character(len=12),intent(in) :: DFA
logical,intent(in) :: LDA_centered
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: FxHF(nBas,nBas)
double precision,intent(in) :: rho(nGrid)
double precision,intent(in) :: drho(ncart,nGrid)
! Local variables
double precision :: ExLDA,ExGGA,ExHF
double precision :: cX,aX,aC
! Output variables
double precision,intent(out) :: Ex
select case (rung)
! Hartree calculation
case(0)
Ex = 0d0
! LDA functionals
case(1)
call unrestricted_lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,ExLDA)
Ex = ExLDA
! GGA functionals
case(2)
call unrestricted_gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,ExGGA)
Ex = ExGGA
! Hybrid functionals
case(4)
cX = 0.20d0
aX = 0.72d0
aC = 0.81d0
call unrestricted_lda_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,ExLDA)
call unrestricted_gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,ExGGA)
call unrestricted_fock_exchange_energy(nBas,P,FxHF,ExHF)
Ex = ExLDA &
+ cX*(ExHF - ExLDA) &
+ aX*(ExGGA - ExLDA)
! Hartree-Fock calculation
case(666)
call unrestricted_fock_exchange_energy(nBas,P,FxHF,ExHF)
Ex = ExHF
end select
end subroutine unrestricted_exchange_energy

80
src/eDFT/unrestricted_exchange_individual_energy.f90 Normal file
View File

@ -0,0 +1,80 @@
subroutine unrestricted_exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, &
ERI,Pw,P,rhow,drhow,rho,drho,Ex)
! Compute the exchange individual energy
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: rung
character(len=12),intent(in) :: DFA
logical,intent(in) :: LDA_centered
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: Pw(nBas,nBas)
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: rhow(nGrid)
double precision,intent(in) :: drhow(ncart,nGrid)
double precision,intent(in) :: rho(nGrid)
double precision,intent(in) :: drho(ncart,nGrid)
! Local variables
double precision :: ExLDA
double precision :: ExGGA
double precision :: ExHF
! Output variables
double precision,intent(out) :: Ex
select case (rung)
! Hartree calculation
case(0)
Ex = 0d0
! LDA functionals
case(1)
call unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,ExLDA)
Ex = ExLDA
! GGA functionals
case(2)
call unrestricted_gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,ExGGA)
Ex = ExGGA
! Hybrid functionals
case(4)
call print_warning('!!! Individual energies NYI for Hybrids !!!')
stop
! Hartree-Fock calculation
case(666)
call unrestricted_fock_exchange_individual_energy(nBas,Pw,P,ERI,ExHF)
Ex = ExHF
end select
end subroutine unrestricted_exchange_individual_energy

86
src/eDFT/unrestricted_exchange_potential.f90 Normal file
View File

@ -0,0 +1,86 @@
subroutine unrestricted_exchange_potential(rung,DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,P, &
ERI,AO,dAO,rho,drho,Fx,FxHF)
! Compute the exchange potential
implicit none
include 'parameters.h'
! Input variables
integer,intent(in) :: rung
character(len=12),intent(in) :: DFA
logical,intent(in) :: LDA_centered
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: AO(nBas,nGrid)
double precision,intent(in) :: dAO(ncart,nBas,nGrid)
double precision,intent(in) :: rho(nGrid)
double precision,intent(in) :: drho(ncart,nGrid)
! Local variables
double precision,allocatable :: FxLDA(:,:),FxGGA(:,:)
double precision :: cX,aX
! Output variables
double precision,intent(out) :: Fx(nBas,nBas),FxHF(nBas,nBas)
! Memory allocation
select case (rung)
! Hartree calculation
case(0)
Fx(:,:) = 0d0
! LDA functionals
case(1)
call unrestricted_lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
! GGA functionals
case(2)
call unrestricted_gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
! Hybrid functionals
case(4)
allocate(FxLDA(nBas,nBas),FxGGA(nBas,nBas))
cX = 0.20d0
aX = 0.72d0
call unrestricted_lda_exchange_potential(DFA,nGrid,weight,nBas,AO,rho,FxLDA)
call unrestricted_gga_exchange_potential(DFA,nGrid,weight,nBas,AO,dAO,rho,drho,FxGGA)
call unrestricted_fock_exchange_potential(nBas,P,ERI,FxHF)
Fx(:,:) = FxLDA(:,:) &
+ cX*(FxHF(:,:) - FxLDA(:,:)) &
+ aX*(FxGGA(:,:) - FxLDA(:,:))
! Hartree-Fock calculation
case(666)
call unrestricted_fock_exchange_potential(nBas,P,ERI,FxHF)
Fx(:,:) = FxHF(:,:)
end select
end subroutine unrestricted_exchange_potential

25
src/eDFT/unrestricted_fock_exchange_energy.f90 Normal file
View File

@ -0,0 +1,25 @@
subroutine unrestricted_fock_exchange_energy(nBas,P,Fx,Ex)
! Compute the (exact) Fock exchange energy
implicit none
! Input variables
integer,intent(in) :: nBas
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: Fx(nBas,nBas)
! Local variables
double precision,external :: trace_matrix
! Output variables
double precision,intent(out) :: Ex
! Compute HF exchange energy
Ex = trace_matrix(nBas,matmul(P,Fx))
end subroutine unrestricted_fock_exchange_energy

31
src/eDFT/unrestricted_fock_exchange_individual_energy.f90 Normal file
View File

@ -0,0 +1,31 @@
subroutine unrestricted_fock_exchange_individual_energy(nBas,Pw,P,ERI,Ex)
! Compute the Fock exchange potential
implicit none
! Input variables
integer,intent(in) :: nBas
double precision,intent(in) :: Pw(nBas,nBas)
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
! Local variables
double precision,allocatable :: Fx(:,:)
double precision,external :: trace_matrix
! Output variables
double precision,intent(out) :: Ex
! Compute HF exchange matrix
allocate(Fx(nBas,nBas))
call unrestricted_fock_exchange_potential(nBas,Pw(:,:),ERI(:,:,:,:),Fx(:,:))
Ex = trace_matrix(nBas,matmul(P(:,:),Fx(:,:))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:),Fx(:,:)))
end subroutine unrestricted_fock_exchange_individual_energy

34
src/eDFT/unrestricted_fock_exchange_potential.f90 Normal file
View File

@ -0,0 +1,34 @@
subroutine unrestricted_fock_exchange_potential(nBas,P,ERI,Fx)
! Compute the Fock exchange potential
implicit none
! Input variables
integer,intent(in) :: nBas
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
! Local variables
integer :: mu,nu,la,si
! Output variables
double precision,intent(out) :: Fx(nBas,nBas)
! Compute HF exchange matrix
Fx(:,:) = 0d0
do si=1,nBas
do la=1,nBas
do nu=1,nBas
do mu=1,nBas
Fx(mu,nu) = Fx(mu,nu) - P(la,si)*ERI(mu,la,si,nu)
enddo
enddo
enddo
enddo
end subroutine unrestricted_fock_exchange_potential

40
src/eDFT/unrestricted_gga_exchange_derivative_discontinuity.f90 Normal file
View File

@ -0,0 +1,40 @@
subroutine unrestricted_gga_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow,drhow,ExDD)
! Compute the exchange GGA part of the derivative discontinuity
implicit none
include 'parameters.h'
! Input variables
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid)
double precision,intent(in) :: drhow(ncart,nGrid)
! Local variables
! Output variables
double precision,intent(out) :: ExDD(nEns)
! Select correlation functional
select case (DFA)
case ('B88')
ExDD(:) = 0d0
case default
call print_warning('!!! GGA exchange derivative discontinuity not available !!!')
stop
end select
end subroutine unrestricted_gga_exchange_derivative_discontinuity

40
src/eDFT/unrestricted_gga_exchange_energy.f90 Normal file
View File

@ -0,0 +1,40 @@
subroutine unrestricted_gga_exchange_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,Ex)
! Select GGA exchange functional for energy calculation
implicit none
include 'parameters.h'
! Input variables
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rho(nGrid)
double precision,intent(in) :: drho(3,nGrid)
! Output variables
double precision :: Ex
select case (DFA)
case ('G96')
call UG96_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
case ('B88')
call UB88_gga_exchange_energy(nGrid,weight,rho,drho,Ex)
case default
call print_warning('!!! GGA exchange energy not available !!!')
stop
end select
end subroutine unrestricted_gga_exchange_energy

35
src/eDFT/unrestricted_gga_exchange_individual_energy.f90 Normal file
View File

@ -0,0 +1,35 @@
subroutine unrestricted_gga_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,Ex)
! Compute GGA exchange energy for individual states
implicit none
include 'parameters.h'
! Input variables
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid)
double precision,intent(in) :: drhow(ncart,nGrid)
double precision,intent(in) :: rho(nGrid)
double precision,intent(in) :: drho(ncart,nGrid)
! Output variables
double precision :: Ex
! Select correlation functional
select case (DFA)
case default
call print_warning('!!! GGA exchange individual energy not available !!!')
stop
end select
end subroutine unrestricted_gga_exchange_individual_energy

44
src/eDFT/unrestricted_gga_exchange_potential.f90 Normal file
View File

@ -0,0 +1,44 @@
subroutine unrestricted_gga_exchange_potential(DFA,nEns,wEns,nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
! Select GGA exchange functional for potential calculation
implicit none
include 'parameters.h'
! Input variables
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: AO(nBas,nGrid)
double precision,intent(in) :: dAO(3,nBas,nGrid)
double precision,intent(in) :: rho(nGrid)
double precision,intent(in) :: drho(3,nGrid)
! Output variables
double precision,intent(out) :: Fx(nBas,nBas)
! Select GGA exchange functional
select case (DFA)
case ('G96')
call UG96_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
case ('B88')
call UB88_gga_exchange_potential(nGrid,weight,nBas,AO,dAO,rho,drho,Fx)
case default
call print_warning('!!! GGA exchange potential not available !!!')
stop
end select
end subroutine unrestricted_gga_exchange_potential

64
src/eDFT/unrestricted_individual_energy.f90
View File

@ -118,15 +118,15 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
! Checking Hartree contributions for each individual states
!------------------------------------------------------------------------
print*,'Hartree contributions for each individual states'
print*,''
print*,''
print*,'EJ(aa,1)=',EJ(1,1),'EJ(ab,1)=',EJ(2,1),'EJ(bb,1)=',EJ(3,1)
print*,''
print*,'EJ(aa,2)=',EJ(1,2),'EJ(ab,2)=',EJ(2,2),'EJ(bb,2)=',EJ(3,2)
print*,''
print*,'EJ(aa,3)=',EJ(1,3),'EJ(ab,3)=',EJ(2,3),'EJ(bb,3)=',EJ(3,3)
print*,''
! print*,'Hartree contributions for each individual states'
! print*,''
! print*,''
! print*,'EJ(aa,1)=',EJ(1,1),'EJ(ab,1)=',EJ(2,1),'EJ(bb,1)=',EJ(3,1)
! print*,''
! print*,'EJ(aa,2)=',EJ(1,2),'EJ(ab,2)=',EJ(2,2),'EJ(bb,2)=',EJ(3,2)
! print*,''
! print*,'EJ(aa,3)=',EJ(1,3),'EJ(ab,3)=',EJ(2,3),'EJ(bb,3)=',EJ(3,3)
! print*,''
!------------------------------------------------------------------------
@ -135,7 +135,7 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
do iEns=1,nEns
do ispin=1,nspin
call exchange_individual_energy(x_rung,x_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,ERI, &
call unrestricted_exchange_individual_energy(x_rung,x_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,ERI, &
Pw(:,:,ispin),P(:,:,ispin,iEns),rhow(:,ispin),drhow(:,:,ispin), &
rho(:,ispin,iEns),drho(:,:,ispin,iEns),Ex(ispin,iEns))
end do
@ -144,32 +144,30 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
!------------------------------------------------------------------------
! Checking exchange contributions for each individual states
!------------------------------------------------------------------------
print*,''
print*,''
print*,'Exchange contributions for each individual states'
print*,''
print*,''
print*,'Ex(aa,1) =' ,Ex(1,1),'Ex(bb,1) =' ,Ex(2,1)
print*,''
print*,'Ex(aa,2) =' ,Ex(1,2),'Ex(bb,2) =' ,Ex(2,2)
print*,''
print*,'Ex(aa,3) =' ,Ex(1,3),'Ex(bb,3) =' ,Ex(2,3)
! print*,''
! print*,''
! print*,'Exchange contributions for each individual states'
! print*,''
! print*,''
! print*,'Ex(aa,1) =' ,Ex(1,1),'Ex(bb,1) =' ,Ex(2,1)
! print*,''
! print*,'Ex(aa,2) =' ,Ex(1,2),'Ex(bb,2) =' ,Ex(2,2)
! print*,''
! print*,'Ex(aa,3) =' ,Ex(1,3),'Ex(bb,3) =' ,Ex(2,3)
!------------------------------------------------------------------------
! Checking number of alpha and beta electrons for each individual states
!------------------------------------------------------------------------
print*,''
print*,''
print*,'Checking number of alpha and beta electrons for each individual states'
print*,''
print*,''
print*,'nEl(a,1) = ',electron_number(nGrid,weight,rho(:,1,1)),'nEl(b,1) = ',electron_number(nGrid,weight,rho(:,2,1))
print*,''
print*,'nEl(a,2) = ',electron_number(nGrid,weight,rho(:,1,2)),'nEl(b,2) = ',electron_number(nGrid,weight,rho(:,2,2))
print*,''
print*,'nEl(a,3) = ',electron_number(nGrid,weight,rho(:,1,3)),'nEl(b,3) = ',electron_number(nGrid,weight,rho(:,2,3))
! print*,''
! print*,''
! print*,'Checking number of alpha and beta electrons for each individual states'
! print*,''
! print*,''
! print*,'nEl(a,1) = ',electron_number(nGrid,weight,rho(:,1,1)),'nEl(b,1) = ',electron_number(nGrid,weight,rho(:,2,1))
! print*,''
! print*,'nEl(a,2) = ',electron_number(nGrid,weight,rho(:,1,2)),'nEl(b,2) = ',electron_number(nGrid,weight,rho(:,2,2))
! print*,''
! print*,'nEl(a,3) = ',electron_number(nGrid,weight,rho(:,1,3)),'nEl(b,3) = ',electron_number(nGrid,weight,rho(:,2,3))
!------------------------------------------------------------------------
! Individual correlation energy
@ -192,7 +190,7 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
do ispin=1,nspin
call exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight, &
call unrestricted_exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight, &
rhow(:,ispin),drhow(:,:,ispin),ExDD(ispin,:))
end do

46
src/eDFT/unrestricted_lda_exchange_derivative_discontinuity.f90 Normal file
View File

@ -0,0 +1,46 @@
subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,ExDD)
! Compute the exchange LDA part of the derivative discontinuity
implicit none
include 'parameters.h'
! Input variables
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid)
! Local variables
! Output variables
double precision,intent(out) :: ExDD(nEns)
! Select correlation functional
select case (DFA)
case ('S51')
ExDD(:) = 0d0
case ('CC')
call UCC_lda_exchange_derivative_discontinuity(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:))
case default
call print_warning('!!! LDA exchange derivative discontinuity not available !!!')
stop
end select
end subroutine unrestricted_lda_exchange_derivative_discontinuity

43
src/eDFT/unrestricted_lda_exchange_energy.f90 Normal file
View File

@ -0,0 +1,43 @@
subroutine unrestricted_lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,Ex)
! Select LDA exchange functional
implicit none
include 'parameters.h'
! Input variables
character(len=12),intent(in) :: DFA
logical,intent(in) :: LDA_centered
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rho(nGrid)
! Output variables
double precision,intent(out) :: Ex
! Select correlation functional
select case (DFA)
case ('S51')
call US51_lda_exchange_energy(nGrid,weight,rho,Ex)
case ('CC')
call UCC_lda_exchange_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rho,Ex)
case default
call print_warning('!!! LDA exchange functional not available !!!')
stop
end select
end subroutine unrestricted_lda_exchange_energy

44
src/eDFT/unrestricted_lda_exchange_individual_energy.f90 Normal file
View File

@ -0,0 +1,44 @@
subroutine unrestricted_lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex)
! Compute LDA exchange energy for individual states
implicit none
include 'parameters.h'
! Input variables
logical,intent(in) :: LDA_centered
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rhow(nGrid)
double precision,intent(in) :: rho(nGrid)
! Output variables
double precision :: Ex
! Select correlation functional
select case (DFA)
case ('S51')
call US51_lda_exchange_individual_energy(nGrid,weight,rhow,rho,Ex)
case ('CC')
call UCC_lda_exchange_individual_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex)
case default
call print_warning('!!! LDA exchange individual energy not available !!!')
stop
end select
end subroutine unrestricted_lda_exchange_individual_energy

46
src/eDFT/unrestricted_lda_exchange_potential.f90 Normal file
View File

@ -0,0 +1,46 @@
subroutine unrestricted_lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
! Select LDA correlation potential
implicit none
include 'parameters.h'
! Input variables
logical,intent(in) :: LDA_centered
character(len=12),intent(in) :: DFA
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
double precision,intent(in) :: aCC_w1(3)
double precision,intent(in) :: aCC_w2(3)
integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: AO(nBas,nGrid)
double precision,intent(in) :: rho(nGrid)
! Output variables
double precision,intent(out) :: Fx(nBas,nBas)
! Select exchange functional
select case (DFA)
case ('S51')
call US51_lda_exchange_potential(nGrid,weight,nBas,AO,rho,Fx)
case ('CC')
call UCC_lda_exchange_potential(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx)
case default
call print_warning('!!! LDA exchange functional not available !!!')
stop
end select
end subroutine unrestricted_lda_exchange_potential