10
1
mirror of https://github.com/pfloos/quack synced 2024-07-24 03:37:44 +02:00

individual energies for lda correlation

This commit is contained in:
Pierre-Francois Loos 2020-03-16 23:28:56 +01:00
parent b9d0e40407
commit 0667c8c5fd
10 changed files with 40 additions and 51 deletions

View File

@ -3,7 +3,7 @@
# exchange rung: Hartree = 0, LDA = 1 (RS51,S51), GGA = 2(G96,B88), hybrid = 4, Hartree-Fock = 666 # exchange rung: Hartree = 0, LDA = 1 (RS51,S51), GGA = 2(G96,B88), hybrid = 4, Hartree-Fock = 666
1 RS51 1 RS51
# correlation rung: Hartree = 0, LDA = 1 (W38,VWN5,C16,LF19), GGA = 2(LYP), hybrid = 4(B3LYP), Hartree-Fock = 666 # correlation rung: Hartree = 0, LDA = 1 (W38,VWN5,C16,LF19), GGA = 2(LYP), hybrid = 4(B3LYP), Hartree-Fock = 666
0 HF 1 RVWN5
# quadrature grid SG-n # quadrature grid SG-n
1 1
# Number of states in ensemble (nEns) # Number of states in ensemble (nEns)

View File

@ -1,6 +1,7 @@
subroutine exchange_individual_energy(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow,rho,drho,Ex) subroutine exchange_individual_energy(rung,DFA,nEns,wEns,nGrid,weight,nBas, &
ERI,P,FxHF,rhow,drhow,rho,drho,Ex)
! Compute the exchange energy of individual states ! Compute the exchange individual energy
implicit none implicit none
include 'parameters.h' include 'parameters.h'
@ -13,6 +14,10 @@ subroutine exchange_individual_energy(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow
double precision,intent(in) :: wEns(nEns) double precision,intent(in) :: wEns(nEns)
integer,intent(in) :: nGrid integer,intent(in) :: nGrid
double precision,intent(in) :: weight(nGrid) double precision,intent(in) :: weight(nGrid)
integer,intent(in) :: nBas
double precision,intent(in) :: ERI(nBas,nBas,nBas,nBas)
double precision,intent(in) :: P(nBas,nBas)
double precision,intent(in) :: FxHF(nBas,nBas)
double precision,intent(in) :: rhow(nGrid) double precision,intent(in) :: rhow(nGrid)
double precision,intent(in) :: drhow(ncart,nGrid) double precision,intent(in) :: drhow(ncart,nGrid)
double precision,intent(in) :: rho(nGrid) double precision,intent(in) :: rho(nGrid)
@ -21,7 +26,6 @@ subroutine exchange_individual_energy(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow
! Local variables ! Local variables
double precision :: ExLDA double precision :: ExLDA
double precision :: ExGGA
double precision :: ExHF double precision :: ExHF
! Output variables ! Output variables
@ -40,7 +44,9 @@ subroutine exchange_individual_energy(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow
case(1) case(1)
call lda_exchange_individual_energy(DFA,nEns,wEns(:),nGrid,weight(:),rhow(:),rho(:),Ex) call lda_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,rho,ExLDA)
Ex = ExLDA
! GGA functionals ! GGA functionals
@ -53,14 +59,16 @@ subroutine exchange_individual_energy(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow
case(4) case(4)
call print_warning('!!! Individual energies NYI for hybrids !!!') call print_warning('!!! Individual energies NYI for Hybrids !!!')
stop stop
! Hartree-Fock calculation ! Hartree-Fock calculation
case(666) case(666)
! call fock_exchange_individual_energy(nEns,wEns(:),nBas,P,FxHF,ExHF) call fock_exchange_potential(nBas,P,ERI,FxHF)
call fock_exchange_energy(nBas,P,FxHF,ExHF)
Ex = ExHF Ex = ExHF
end select end select

View File

@ -44,12 +44,6 @@ subroutine lda_correlation_energy(DFA,nEns,wEns,nGrid,weight,rho,Ec)
call C16_lda_correlation_energy(nGrid,weight(:),rho(:,:),Ec(:)) call C16_lda_correlation_energy(nGrid,weight(:),rho(:,:),Ec(:))
! Restricted version of the Marut-Fromager-Loos weight-dependent correlation functional
case ('RMFL20')
call RMFL20_lda_correlation_energy(nEns,wEns(:),nGrid,weight(:),rho(:,:),Ec(:))
case default case default
call print_warning('!!! LDA correlation functional not available !!!') call print_warning('!!! LDA correlation functional not available !!!')

View File

@ -48,12 +48,6 @@ include 'parameters.h'
call C16_lda_correlation_potential(nGrid,weight(:),nBas,AO(:,:),rho(:,:),Fc(:,:,:)) call C16_lda_correlation_potential(nGrid,weight(:),nBas,AO(:,:),rho(:,:),Fc(:,:,:))
! Restricted version of the Marut-Fromager-Loos weight-dependent correlation functional
case ('RMFL20')
call RMFL20_lda_correlation_potential(nEns,wEns(:),nGrid,weight(:),nBas,AO(:,:),rho(:,:),Fc(:,:,:))
case default case default
call print_warning('!!! LDA correlation functional not available !!!') call print_warning('!!! LDA correlation functional not available !!!')

View File

@ -29,6 +29,10 @@ subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow
ExDD(:) = 0d0 ExDD(:) = 0d0
case ('RS51')
ExDD(:) = 0d0
case ('RMFL20') case ('RMFL20')
! call MFL20_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight(:),rhow(:),ExDD(:)) ! call MFL20_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight(:),rhow(:),ExDD(:))

View File

@ -23,13 +23,13 @@ subroutine lda_exchange_individual_energy(DFA,nEns,wEns,nGrid,weight,rhow,rho,Ex
select case (DFA) select case (DFA)
case ('S51') case ('RS51')
! call S51_lda_exchange_individual_energy(nGrid,weight(:),rhow(:),rho(:),Ex) call RS51_lda_exchange_individual_energy(nGrid,weight(:),rhow(:),rho(:),Ex)
case ('MFL20') case ('RMFL20')
! call MFL20_lda_exchange_individual_energy(nEns,wEns,nGrid,weight(:),rhow(:),rho(:),Ex) ! call RMFL20_lda_exchange_individual_energy(nEns,wEns,nGrid,weight(:),rhow(:),rho(:),Ex)
case default case default

View File

@ -39,8 +39,8 @@ subroutine read_options(method,x_rung,x_DFA,c_rung,c_DFA,SGn,nEns,wEns,maxSCF,th
method = 'GOK-RKS' method = 'GOK-RKS'
x_rung = 1 x_rung = 1
c_rung = 1 c_rung = 1
x_DFA = 'S51' x_DFA = 'RS51'
c_DFA = 'W38' c_DFA = 'RVWN5'
SGn = 0 SGn = 0
wEns(:) = 0d0 wEns(:) = 0d0

View File

@ -38,24 +38,20 @@ subroutine restricted_correlation_energy(rung,DFA,nEns,wEns,nGrid,weight,rho,drh
case(1) case(1)
call lda_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:),Ec) call restricted_lda_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:),Ec)
! GGA functionals ! GGA functionals
case(2) case(2)
call gga_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:),drho(:,:),Ec) call print_warning('!!! restricted correlation energies NYI for GGAs !!!')
stop
! Hybrid functionals
case(4) case(4)
aC = 0.81d0 call print_warning('!!! restricted correlation energies NYI for Hybrids !!!')
stop
call lda_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:),EcLDA)
call gga_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:),drho(:,:),EcGGA)
Ec = EcLDA + aC*(EcGGA - EcLDA)
! Hartree-Fock calculation ! Hartree-Fock calculation

View File

@ -22,8 +22,6 @@ subroutine restricted_correlation_potential(rung,DFA,nEns,wEns,nGrid,weight,nBas
! Local variables ! Local variables
double precision,allocatable :: FcLDA(:,:) double precision,allocatable :: FcLDA(:,:)
double precision,allocatable :: FcGGA(:,:)
double precision :: aC
! Output variables ! Output variables
@ -43,26 +41,21 @@ subroutine restricted_correlation_potential(rung,DFA,nEns,wEns,nGrid,weight,nBas
case(1) case(1)
call lda_correlation_potential(DFA,nEns,wEns(:),nGrid,weight(:),nBas,AO(:,:),rho(:),Fc(:,:)) call restricted_lda_correlation_potential(DFA,nEns,wEns(:),nGrid,weight(:),nBas,AO(:,:),rho(:),Fc(:,:))
! GGA functionals ! GGA functionals
case(2) case(2)
call gga_correlation_potential(DFA,nEns,wEns(:),nGrid,weight(:),nBas,AO(:,:),dAO(:,:,:),rho(:),drho(:,:),Fc(:,:)) call print_warning('!!! restricted correlation potentials NYI for GGAs !!!')
stop
! Hybrid functionals ! Hybrid functionals
case(4) case(4)
allocate(FcLDA(nBas,nBas),FcGGA(nBas,nBas)) call print_warning('!!! restricted correlation potentials NYI for Hybrids !!!')
stop
aC = 0.81d0
call lda_correlation_potential(DFA,nEns,wEns(:),nGrid,weight(:),nBas,AO(:,:),rho(:),FcLDA(:,:))
call gga_correlation_potential(DFA,nEns,wEns(:),nGrid,weight(:),nBas,AO(:,:),dAO(:,:,:),rho(:),drho(:,:),FcGGA(:,:))
Fc(:,:) = FcLDA(:,:) + aC*(FcGGA(:,:) - FcLDA(:,:))
! Hartree-Fock calculation ! Hartree-Fock calculation

View File

@ -75,7 +75,7 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGri
!------------------------------------------------------------------------ !------------------------------------------------------------------------
do iEns=1,nEns do iEns=1,nEns
call hartree_coulomb(nBas,P(:,:,iEns),ERI,J(:,:)) call hartree_coulomb(nBas,P(:,:,iEns),ERI(:,:,:,:),J(:,:))
EJ(iEns) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,iEns),J(:,:))) EJ(iEns) = 0.5d0*trace_matrix(nBas,matmul(P(:,:,iEns),J(:,:)))
end do end do
@ -85,8 +85,8 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGri
do iEns=1,nEns do iEns=1,nEns
call exchange_energy_individual_energy(x_rung,x_DFA,nEns,wEns(:),nGrid,weight(:),nBas,P(:,:,iEns),FxHF(:,:), & call exchange_individual_energy(x_rung,x_DFA,nEns,wEns(:),nGrid,weight(:),nBas,ERI(:,:,:,:), &
rho(:,iEns),drho(:,:,iEns),Ex(iEns)) P(:,:,iEns),FxHF(:,:),rhow(:),drhow(:,:),rho(:,iEns),drho(:,:,iEns),Ex(iEns))
end do end do
@ -133,7 +133,7 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGri
! Dump results ! Dump results
!------------------------------------------------------------------------ !------------------------------------------------------------------------
call print_restricted_individual_energy(nEns,ET(:),EV(:),EJ(:),Ex(:),Ec(:),Exc(:),ExLZ,EcLZ,ExcLZ, & call print_restricted_individual_energy(nEns,ET(:),EV(:),EJ(:),Ex(:),Ec(:),Exc(:), &
ExDD(:),EcDD(:),ExcDD(:),E(:),Om(:)) ExDD(:),EcDD(:),ExcDD(:),E(:),Om(:))
end subroutine restricted_individual_energy end subroutine restricted_individual_energy