From a5163dbbcb29aecfe2a293fec0a30187ed3c7ffe Mon Sep 17 00:00:00 2001 From: Pierre-Francois Loos Date: Tue, 17 Mar 2020 22:39:04 +0100 Subject: [PATCH] Debug in progress --- input/basis | 46 +++++-------------- input/dft | 2 +- input/weight | 46 +++++-------------- ...a_correlation_derivative_discontinuity.f90 | 12 ++--- ..._lda_exchange_derivative_discontinuity.f90 | 36 +++++++-------- src/eDFT/restricted_elda_exchange_energy.f90 | 42 +++++++++++++++++ 6 files changed, 90 insertions(+), 94 deletions(-) create mode 100644 src/eDFT/restricted_elda_exchange_energy.f90 diff --git a/input/basis b/input/basis index 6d7be41..0e9b57e 100644 --- a/input/basis +++ b/input/basis @@ -1,35 +1,13 @@ -1 9 -S 8 - 1 2940.0000000 0.0006800 - 2 441.2000000 0.0052360 - 3 100.5000000 0.0266060 - 4 28.4300000 0.0999930 - 5 9.1690000 0.2697020 - 6 3.1960000 0.4514690 - 7 1.1590000 0.2950740 - 8 0.1811000 0.0125870 -S 8 - 1 2940.0000000 -0.0001230 - 2 441.2000000 -0.0009660 - 3 100.5000000 -0.0048310 - 4 28.4300000 -0.0193140 - 5 9.1690000 -0.0532800 - 6 3.1960000 -0.1207230 - 7 1.1590000 -0.1334350 - 8 0.1811000 0.5307670 -S 1 - 1 0.0589000 1.0000000 -S 1 - 1 0.0187700 1.0000000 +1 3 +S 3 + 1 30.1678710 0.15432897 + 2 5.4951153 0.53532814 + 3 1.4871927 0.44463454 +S 3 + 1 1.3148331 -0.0999672 + 2 0.3055389 0.3995128 + 3 0.0993707 0.7001155 P 3 - 1 3.6190000 0.0291110 - 2 0.7110000 0.1693650 - 3 0.1951000 0.5134580 -P 1 - 1 0.0601800 1.0000000 -P 1 - 1 0.0085000 1.0000000 -D 1 - 1 0.2380000 1.0000000 -D 1 - 1 0.0740000 1.0000000 + 1 1.3148331 0.1559163 + 2 0.3055389 0.6076837 + 3 0.0993707 0.3919574 diff --git a/input/dft b/input/dft index 53337d0..6a296de 100644 --- a/input/dft +++ b/input/dft @@ -19,6 +19,6 @@ # Number of states in ensemble (nEns) 2 # Ensemble weights: wEns(1),...,wEns(nEns-1) - 0.50000 0.00000 + 0.00000 0.00000 # GOK-DFT: maxSCF thresh DIIS n_diis guess_type ortho_type 64 0.0000001 T 5 1 1 diff --git a/input/weight b/input/weight index 6d7be41..0e9b57e 100644 --- a/input/weight +++ b/input/weight @@ -1,35 +1,13 @@ -1 9 -S 8 - 1 2940.0000000 0.0006800 - 2 441.2000000 0.0052360 - 3 100.5000000 0.0266060 - 4 28.4300000 0.0999930 - 5 9.1690000 0.2697020 - 6 3.1960000 0.4514690 - 7 1.1590000 0.2950740 - 8 0.1811000 0.0125870 -S 8 - 1 2940.0000000 -0.0001230 - 2 441.2000000 -0.0009660 - 3 100.5000000 -0.0048310 - 4 28.4300000 -0.0193140 - 5 9.1690000 -0.0532800 - 6 3.1960000 -0.1207230 - 7 1.1590000 -0.1334350 - 8 0.1811000 0.5307670 -S 1 - 1 0.0589000 1.0000000 -S 1 - 1 0.0187700 1.0000000 +1 3 +S 3 + 1 30.1678710 0.15432897 + 2 5.4951153 0.53532814 + 3 1.4871927 0.44463454 +S 3 + 1 1.3148331 -0.0999672 + 2 0.3055389 0.3995128 + 3 0.0993707 0.7001155 P 3 - 1 3.6190000 0.0291110 - 2 0.7110000 0.1693650 - 3 0.1951000 0.5134580 -P 1 - 1 0.0601800 1.0000000 -P 1 - 1 0.0085000 1.0000000 -D 1 - 1 0.2380000 1.0000000 -D 1 - 1 0.0740000 1.0000000 + 1 1.3148331 0.1559163 + 2 0.3055389 0.6076837 + 3 0.0993707 0.3919574 diff --git a/src/eDFT/RMFL20_lda_correlation_derivative_discontinuity.f90 b/src/eDFT/RMFL20_lda_correlation_derivative_discontinuity.f90 index 5b2b0e6..eb313aa 100644 --- a/src/eDFT/RMFL20_lda_correlation_derivative_discontinuity.f90 +++ b/src/eDFT/RMFL20_lda_correlation_derivative_discontinuity.f90 @@ -1,4 +1,4 @@ -subroutine RMFL20_lda_correlation_derivative_discontinuity(nEns,wEns,nGrid,weight,rhow,Ec) +subroutine RMFL20_lda_correlation_derivative_discontinuity(nEns,wEns,nGrid,weight,rhow,EcDD) ! Compute the restricted version of the eLDA correlation part of the derivative discontinuity @@ -17,12 +17,12 @@ subroutine RMFL20_lda_correlation_derivative_discontinuity(nEns,wEns,nGrid,weigh integer :: iEns,jEns double precision,allocatable :: aMFL(:,:) - double precision :: dEc(nEns) + double precision :: dEcdw(nEns) double precision,external :: Kronecker_delta ! Output variables - double precision,intent(out) :: Ec(nEns) + double precision,intent(out) :: EcDD(nEns) ! Allocation @@ -42,16 +42,16 @@ subroutine RMFL20_lda_correlation_derivative_discontinuity(nEns,wEns,nGrid,weigh do iEns=1,nEns - call restricted_elda_correlation_energy(nEns,aMFL(:,iEns),nGrid,weight(:),rhow(:),dEc(iEns)) + call restricted_elda_correlation_energy(nEns,aMFL(:,iEns),nGrid,weight(:),rhow(:),dEcdw(iEns)) end do - Ec(:) = 0d0 + EcDD(:) = 0d0 do iEns=1,nEns do jEns=1,nEns - Ec(iEns) = Ec(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*(dEc(jEns) - dEc(1)) + EcDD(iEns) = EcDD(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*(dEcdw(jEns) - dEcdw(1)) end do end do diff --git a/src/eDFT/RMFL20_lda_exchange_derivative_discontinuity.f90 b/src/eDFT/RMFL20_lda_exchange_derivative_discontinuity.f90 index ab6b337..d67a46e 100644 --- a/src/eDFT/RMFL20_lda_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/RMFL20_lda_exchange_derivative_discontinuity.f90 @@ -15,11 +15,10 @@ subroutine RMFL20_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight,r ! Local variables - integer :: iG - double precision :: Cx0 - double precision :: Cx1 - double precision :: rw - double precision :: dExdw + integer :: iEns,jEns + double precision :: Cx(nEns) + double precision :: dExdw(nEns) + double precision,external :: Kronecker_delta ! Output variables @@ -27,26 +26,25 @@ subroutine RMFL20_lda_exchange_derivative_discontinuity(nEns,wEns,nGrid,weight,r ! Weight-dependent Cx coefficient for RMFL20 exchange functional - Cx0 = -(4d0/3d0)*(1d0/pi)**(1d0/3d0) - Cx1 = -(176d0/105d0)*(1d0/pi)**(1d0/3d0) + Cx(1) = -(4d0/3d0)*(1d0/pi)**(1d0/3d0) + Cx(2) = -(176d0/105d0)*(1d0/pi)**(1d0/3d0) ! Compute correlation energy for ground, singly-excited and doubly-excited states - dExdw = 0d0 + do iEns=1,nEns - do iG=1,nGrid - - rw = max(0d0,rhow(iG)) - - if(rw > threshold) then - - dExdw = dExdw + weight(iG)*(Cx1 - Cx0)*rw**(4d0/3d0) - - end if + call restricted_elda_exchange_energy(nEns,Cx(iEns),nGrid,weight(:),rhow(:),dExdw(iEns)) end do - ExDD(1) = - wEns(2) *dExdw - ExDD(2) = (1d0 - wEns(2))*dExdw + ExDD(:) = 0d0 + + do iEns=1,nEns + do jEns=1,nEns + + ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*(dExdw(jEns) - dExdw(1)) + + end do + end do end subroutine RMFL20_lda_exchange_derivative_discontinuity diff --git a/src/eDFT/restricted_elda_exchange_energy.f90 b/src/eDFT/restricted_elda_exchange_energy.f90 new file mode 100644 index 0000000..f56ce0e --- /dev/null +++ b/src/eDFT/restricted_elda_exchange_energy.f90 @@ -0,0 +1,42 @@ +subroutine restricted_elda_exchange_energy(nEns,Cx,nGrid,weight,rho,Ex) + +! Compute the restricted LDA exchange energy of 2-glomium for various states + + implicit none + include 'parameters.h' + +! Input variables + + integer,intent(in) :: nEns + double precision,intent(in) :: Cx + integer,intent(in) :: nGrid + double precision,intent(in) :: weight(nGrid) + double precision,intent(in) :: rho(nGrid) + +! Local variables + + integer :: iG + double precision :: r + +! Output variables + + double precision,intent(out) :: Ex + + +! Compute eLDA exchange energy + + Ex = 0d0 + + do iG=1,nGrid + + r = max(0d0,rho(iG)) + + if(r > threshold) then + + Ex = Ex + weight(iG)*Cx*r**(4d0/3d0) + + end if + + end do + +end subroutine restricted_elda_exchange_energy