From 09d3b39e6378b6710d893b524e5a0630553c6204 Mon Sep 17 00:00:00 2001 From: Clotilde Marut Date: Thu, 2 Jul 2020 18:49:42 +0200 Subject: [PATCH] modifications for the parameters of the weight-dependent functionnals --- ...C_lda_exchange_derivative_discontinuity.f90 | 18 +++++++++++------- .../UCC_lda_exchange_individual_energy.f90 | 15 ++++++++------- src/eDFT/UCC_lda_exchange_potential.f90 | 18 ++++++++++-------- src/eDFT/eDFT_UKS.f90 | 4 ++-- src/eDFT/exchange_derivative_discontinuity.f90 | 6 ++++-- src/eDFT/exchange_individual_energy.f90 | 6 ++++-- src/eDFT/exchange_potential.f90 | 6 ++++-- .../lda_exchange_derivative_discontinuity.f90 | 9 ++++++--- src/eDFT/lda_exchange_individual_energy.f90 | 6 ++++-- src/eDFT/lda_exchange_potential.f90 | 8 +++++--- src/eDFT/unrestricted_individual_energy.f90 | 8 +++++--- 11 files changed, 63 insertions(+), 41 deletions(-) diff --git a/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 b/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 index 9941626..17a92ae 100644 --- a/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 @@ -1,4 +1,4 @@ -subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight,rhow,ExDD) +subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,ExDD) ! Compute the unrestricted version of the curvature-corrected exchange ensemble derivative @@ -9,6 +9,8 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight,rhow 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) @@ -59,15 +61,17 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight,rhow ! Parameters for He N -> N-1 - a1 = 0.420243d0 - b1 = 0.0700561d0 - c1 = -0.288301d0 + a1 = aCC_w1(1) + b1 = aCC_w1(2) + c1 = aCC_w1(3) + ! Parameters for He N -> N+1 - a2 = 0.135068d0 - b2 = -0.00774769d0 - c2 = -0.0278205d0 + a2 = aCC_w2(1) + b2 = aCC_w2(2) + c2 = aCC_w2(3) + ! Cx coefficient for unrestricted Slater LDA exchange diff --git a/src/eDFT/UCC_lda_exchange_individual_energy.f90 b/src/eDFT/UCC_lda_exchange_individual_energy.f90 index 6c52e6b..65a4ff2 100644 --- a/src/eDFT/UCC_lda_exchange_individual_energy.f90 +++ b/src/eDFT/UCC_lda_exchange_individual_energy.f90 @@ -1,4 +1,4 @@ -subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nGrid,weight,rhow,rho,Ex) +subroutine UCC_lda_exchange_individual_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,Ex) ! Compute the unrestricted version of the curvature-corrected exchange functional @@ -9,6 +9,7 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nGrid,weight,rhow,rho,Ex 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) @@ -54,15 +55,15 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nGrid,weight,rhow,rho,Ex ! Parameters for He N -> N-1 - a1 = 0.420243d0 - b1 = 0.0700561d0 - c1 = -0.288301d0 + a1 = aCC_w1(1) + b1 = aCC_w1(2) + c1 = aCC_w1(3) ! Parameters for He N -> N+1 - a2 = 0.135068d0 - b2 = -0.00774769d0 - c2 = -0.0278205d0 + a2 = aCC_w2(1) + b2 = aCC_w2(2) + c2 = aCC_w2(3) ! Cx coefficient for unrestricted Slater LDA exchange diff --git a/src/eDFT/UCC_lda_exchange_potential.f90 b/src/eDFT/UCC_lda_exchange_potential.f90 index 4e73d74..7c853e8 100644 --- a/src/eDFT/UCC_lda_exchange_potential.f90 +++ b/src/eDFT/UCC_lda_exchange_potential.f90 @@ -1,4 +1,4 @@ -subroutine UCC_lda_exchange_potential(nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) +subroutine UCC_lda_exchange_potential(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx) ! Compute the unrestricted version of the curvature-corrected exchange potential @@ -9,6 +9,8 @@ subroutine UCC_lda_exchange_potential(nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) 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 @@ -54,15 +56,15 @@ subroutine UCC_lda_exchange_potential(nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) ! Parameters for He N -> N-1 - a1 = 0.420243d0 - b1 = 0.0700561d0 - c1 = -0.288301d0 - + a1 = aCC_w1(1) + b1 = aCC_w1(2) + c1 = aCC_w1(3) + ! Parameters for He N -> N+1 - a2 = 0.135068d0 - b2 = -0.00774769d0 - c2 = -0.0278205d0 + a2 = aCC_w2(1) + b2 = aCC_w2(2) + c2 = aCC_w2(3) ! Cx coefficient for unrestricted Slater LDA exchange diff --git a/src/eDFT/eDFT_UKS.f90 b/src/eDFT/eDFT_UKS.f90 index b4f0604..e025da1 100644 --- a/src/eDFT/eDFT_UKS.f90 +++ b/src/eDFT/eDFT_UKS.f90 @@ -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(:),nGrid,weight(:),nBas,Pw(:,:,ispin),ERI(:,:,:,:), & + call 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 @@ -366,7 +366,7 @@ subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weig ! Compute individual energies from ensemble energy !------------------------------------------------------------------------ - call unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas, & + call unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, & AO,dAO,nO,nV,T,V,ERI,ENuc,eps,Pw,rhow,drhow,J,Fx,FxHF,Fc,P,rho,drho,Ew,E,Om) end subroutine eDFT_UKS diff --git a/src/eDFT/exchange_derivative_discontinuity.f90 b/src/eDFT/exchange_derivative_discontinuity.f90 index 088e616..9f822de 100644 --- a/src/eDFT/exchange_derivative_discontinuity.f90 +++ b/src/eDFT/exchange_derivative_discontinuity.f90 @@ -1,4 +1,4 @@ -subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nGrid,weight,rhow,drhow,ExDD) +subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,drhow,ExDD) ! Compute the exchange part of the derivative discontinuity @@ -11,6 +11,8 @@ subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nGrid,weight,rho 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) @@ -35,7 +37,7 @@ subroutine exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nGrid,weight,rho case(1) - call lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),nGrid,weight(:),rhow(:),ExDD(:)) + call lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:)) ! GGA functionals diff --git a/src/eDFT/exchange_individual_energy.f90 b/src/eDFT/exchange_individual_energy.f90 index ee81b97..9ee0bbc 100644 --- a/src/eDFT/exchange_individual_energy.f90 +++ b/src/eDFT/exchange_individual_energy.f90 @@ -1,4 +1,4 @@ -subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas, & +subroutine 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 @@ -13,6 +13,8 @@ subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nGrid,weig 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 @@ -46,7 +48,7 @@ subroutine exchange_individual_energy(rung,DFA,LDA_centered,nEns,wEns,nGrid,weig case(1) - call lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,ExLDA) + call lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,rho,ExLDA) Ex = ExLDA diff --git a/src/eDFT/exchange_potential.f90 b/src/eDFT/exchange_potential.f90 index 758be66..b7fd6cd 100644 --- a/src/eDFT/exchange_potential.f90 +++ b/src/eDFT/exchange_potential.f90 @@ -1,4 +1,4 @@ -subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas,P,ERI,AO,dAO,rho,drho,Fx,FxHF) +subroutine 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 @@ -12,6 +12,8 @@ subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas, 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 @@ -45,7 +47,7 @@ subroutine exchange_potential(rung,DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas, case(1) - call lda_exchange_potential(DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) + call lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx) ! GGA functionals diff --git a/src/eDFT/lda_exchange_derivative_discontinuity.f90 b/src/eDFT/lda_exchange_derivative_discontinuity.f90 index ed0ce0c..45fd060 100644 --- a/src/eDFT/lda_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/lda_exchange_derivative_discontinuity.f90 @@ -1,4 +1,4 @@ -subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow,ExDD) +subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,rhow,ExDD) ! Compute the exchange LDA part of the derivative discontinuity @@ -10,6 +10,9 @@ subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow 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) @@ -39,11 +42,11 @@ subroutine lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nGrid,weight,rhow case ('RCC') - call RCC_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight(:),rhow(:),ExDD(:)) + 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,nGrid,weight(:),rhow(:),ExDD(:)) + call UCC_lda_exchange_derivative_discontinuity(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),ExDD(:)) case default diff --git a/src/eDFT/lda_exchange_individual_energy.f90 b/src/eDFT/lda_exchange_individual_energy.f90 index e9ea3cc..b8e90fd 100644 --- a/src/eDFT/lda_exchange_individual_energy.f90 +++ b/src/eDFT/lda_exchange_individual_energy.f90 @@ -1,4 +1,4 @@ -subroutine lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weight,rhow,rho,Ex) +subroutine 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 @@ -11,6 +11,8 @@ subroutine lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weigh 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) @@ -42,7 +44,7 @@ subroutine lda_exchange_individual_energy(DFA,LDA_centered,nEns,wEns,nGrid,weigh case ('UCC') - call UCC_lda_exchange_individual_energy(nEns,wEns,nGrid,weight(:),rhow(:),rho(:),Ex) + call UCC_lda_exchange_individual_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight(:),rhow(:),rho(:),Ex) case default diff --git a/src/eDFT/lda_exchange_potential.f90 b/src/eDFT/lda_exchange_potential.f90 index 14d46b8..c486267 100644 --- a/src/eDFT/lda_exchange_potential.f90 +++ b/src/eDFT/lda_exchange_potential.f90 @@ -1,4 +1,4 @@ -subroutine lda_exchange_potential(DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) +subroutine lda_exchange_potential(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx) ! Select LDA correlation potential @@ -12,6 +12,8 @@ subroutine lda_exchange_potential(DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas,A 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 @@ -40,11 +42,11 @@ subroutine lda_exchange_potential(DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas,A case ('RCC') - call RCC_lda_exchange_potential(nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) + call RCC_lda_exchange_potential(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx) case ('UCC') - call UCC_lda_exchange_potential(nEns,wEns,nGrid,weight,nBas,AO,rho,Fx) + call UCC_lda_exchange_potential(nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,rho,Fx) case default diff --git a/src/eDFT/unrestricted_individual_energy.f90 b/src/eDFT/unrestricted_individual_energy.f90 index 61b33f9..f6f0b8e 100644 --- a/src/eDFT/unrestricted_individual_energy.f90 +++ b/src/eDFT/unrestricted_individual_energy.f90 @@ -1,4 +1,4 @@ -subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,nGrid,weight,nBas,AO,dAO, & +subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,dAO, & nO,nV,T,V,ERI,ENuc,eps,Pw,rhow,drhow,J,Fx,FxHF,Fc,P,rho,drho,Ew,E,Om) ! Compute unrestricted individual energies as well as excitation energies @@ -13,6 +13,8 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered 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 @@ -116,7 +118,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,nGrid,weight,nBas,ERI, & + call 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 @@ -143,7 +145,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,nGrid,weight, & + call exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight, & rhow(:,ispin),drhow(:,:,ispin),ExDD(ispin,:)) end do