From 2b58358681376054f5e7d89a176ca822a848bb2b Mon Sep 17 00:00:00 2001 From: Clotilde Marut Date: Thu, 16 Dec 2021 21:42:58 +0100 Subject: [PATCH 1/2] eDFT code OK --- input/dft | 18 ++++++------- ..._lda_exchange_derivative_discontinuity.f90 | 26 +++++++++++++++---- .../print_unrestricted_individual_energy.f90 | 4 +-- ...cted_exchange_derivative_discontinuity.f90 | 5 ++-- src/eDFT/unrestricted_individual_energy.f90 | 12 +-------- ..._lda_exchange_derivative_discontinuity.f90 | 5 ++-- 6 files changed, 39 insertions(+), 31 deletions(-) diff --git a/input/dft b/input/dft index 346489d..6065b05 100644 --- a/input/dft +++ b/input/dft @@ -13,31 +13,31 @@ # GGA = 2: LYP,PBE # MGGA = 3: # Hybrid = 4: HF,B3LYP,PBE -1 VWN5 +1 VWN5 # quadrature grid SG-n 1 # Number of states in ensemble (nEns) -2 +4 # 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 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 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 +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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 -0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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) - 0.55 0.0 0.0 + 0.3333 0.3333 0.0 # Ncentered ? T # Parameters for CC weight-dependent exchange functional 4 -0.642674 -0.07818 -0.0280307 0.00144198 -0.254939 -0.0893405 0.00765581 0. +0.502145 -0.0672146 0.00786214 -0.00202771 +-1.28842 -0.173117 0.0900511 -0.0118975 0.0 0.0 0.0 0.0 # choice of UCC exchange coefficient : 1 for Cx1, 2 for Cx2, 3 for Cx1*Cx2 -1 +2 diff --git a/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 b/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 index 8da35b4..78d6822 100644 --- a/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 @@ -1,4 +1,5 @@ -subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,weight,rhow,Cx_choice,doNcentered,ExDD) +subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,weight,rhow,Cx_choice,doNcentered,& + kappa,ExDD) ! Compute the unrestricted version of the curvature-corrected exchange ensemble derivative @@ -16,6 +17,7 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei double precision,intent(in) :: rhow(nGrid) integer,intent(in) :: Cx_choice logical,intent(in) :: doNcentered + double precision,intent(in) :: kappa(nEns) ! Local variables @@ -150,12 +152,26 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei ExDD(:) = 0d0 - do iEns=1,nEns - do jEns=2,nEns + if (doNcentered) then - ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*dExdw(jEns) + do iEns=1,nEns + do jEns=2,nEns + ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - kappa(iEns)* wEns(jEns))*dExdw(jEns) + + end do end do - end do + + else + + do iEns=1,nEns + do jEns=2,nEns + + ExDD(iEns) = ExDD(iEns) + (Kronecker_delta(iEns,jEns) - wEns(jEns))*dExdw(jEns) + + end do + end do + + endif end subroutine UCC_lda_exchange_derivative_discontinuity diff --git a/src/eDFT/print_unrestricted_individual_energy.f90 b/src/eDFT/print_unrestricted_individual_energy.f90 index f27e907..7274934 100644 --- a/src/eDFT/print_unrestricted_individual_energy.f90 +++ b/src/eDFT/print_unrestricted_individual_energy.f90 @@ -190,7 +190,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux write(*,'(A60)') '-------------------------------------------------' 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(*,'(A44, F16.10,A3)') ' H energy contribution : ',OmH(iEns), ' au' write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns), ' au' @@ -204,7 +204,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EH,Ex,Ec,Eaux write(*,'(A60)') '-------------------------------------------------' - 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(*,'(A44, F16.10,A3)') ' H energy contribution : ',OmH(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(iEns)*HaToeV, ' eV' diff --git a/src/eDFT/unrestricted_exchange_derivative_discontinuity.f90 b/src/eDFT/unrestricted_exchange_derivative_discontinuity.f90 index c465f8b..b5e3269 100644 --- a/src/eDFT/unrestricted_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/unrestricted_exchange_derivative_discontinuity.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nCC,aCC,nGrid,weight,rhow,drhow,& - Cx_choice,doNcentered,ExDD) + Cx_choice,doNcentered,kappa,ExDD) ! Compute the exchange part of the derivative discontinuity @@ -20,6 +20,7 @@ subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nCC double precision,intent(in) :: drhow(ncart,nGrid) integer,intent(in) :: Cx_choice logical,intent(in) :: doNcentered + double precision,intent(in) :: kappa(nEns) ! Local variables @@ -41,7 +42,7 @@ subroutine unrestricted_exchange_derivative_discontinuity(rung,DFA,nEns,wEns,nCC case(1) call unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns(:),nCC,aCC,nGrid,weight(:),& - rhow(:),Cx_choice,doNcentered,ExDD(:)) + rhow(:),Cx_choice,doNcentered,kappa,ExDD(:)) ! GGA functionals case(2) diff --git a/src/eDFT/unrestricted_individual_energy.f90 b/src/eDFT/unrestricted_individual_energy.f90 index 0227c75..3000e85 100644 --- a/src/eDFT/unrestricted_individual_energy.f90 +++ b/src/eDFT/unrestricted_individual_energy.f90 @@ -150,21 +150,11 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered do ispin=1,nspin call unrestricted_exchange_derivative_discontinuity(x_rung,x_DFA,nEns,wEns,nCC,aCC,nGrid,weight, & - rhow(:,ispin),drhow(:,:,ispin),Cx_choice,doNcentered,ExDD(ispin,:)) + rhow(:,ispin),drhow(:,:,ispin),Cx_choice,doNcentered,kappa,ExDD(ispin,:)) end do call unrestricted_correlation_derivative_discontinuity(c_rung,c_DFA,nEns,wEns,nGrid,weight,rhow,drhow,EcDD) -! Scaling derivative discontinuity for N-centered ensembles - - if(doNcentered) then - - do iEns=1,nEns - ExDD(:,iEns) = (1d0 - kappa(iEns))*ExDD(:,iEns) - EcDD(:,iEns) = (1d0 - kappa(iEns))*EcDD(:,iEns) - end do - - end if !------------------------------------------------------------------------ ! Total energy diff --git a/src/eDFT/unrestricted_lda_exchange_derivative_discontinuity.f90 b/src/eDFT/unrestricted_lda_exchange_derivative_discontinuity.f90 index 0617ab7..d33af22 100644 --- a/src/eDFT/unrestricted_lda_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/unrestricted_lda_exchange_derivative_discontinuity.f90 @@ -1,5 +1,5 @@ subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC,aCC,nGrid,weight,rhow,& - Cx_choice,doNcentered,ExDD) + Cx_choice,doNcentered,kappa,ExDD) ! Compute the exchange LDA part of the derivative discontinuity @@ -19,6 +19,7 @@ subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC, double precision,intent(in) :: rhow(nGrid) integer,intent(in) :: Cx_choice logical,intent(in) :: doNcentered + double precision,intent(in) :: kappa(nEns) ! Local variables @@ -38,7 +39,7 @@ subroutine unrestricted_lda_exchange_derivative_discontinuity(DFA,nEns,wEns,nCC, case (2) call UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,weight,rhow,& - Cx_choice,doNcentered,ExDD) + Cx_choice,doNcentered,kappa,ExDD) case default From 2c86a41d0dfc811da86b1a0d4279df11b6dd7b00 Mon Sep 17 00:00:00 2001 From: Clotilde Marut Date: Thu, 6 Jan 2022 11:18:07 +0100 Subject: [PATCH 2/2] safety commit --- input/dft | 14 +++++++------- input/methods | 4 ++-- src/eDFT/unrestricted_individual_energy.f90 | 4 +++- 3 files changed, 12 insertions(+), 10 deletions(-) diff --git a/input/dft b/input/dft index 6065b05..ed3a55d 100644 --- a/input/dft +++ b/input/dft @@ -19,25 +19,25 @@ # Number of states in ensemble (nEns) 4 # 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 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 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 -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 - 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 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 1 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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) - 0.3333 0.3333 0.0 + 1.0 0. 0.0 # Ncentered ? -T +F # Parameters for CC weight-dependent exchange functional 4 -0.502145 -0.0672146 0.00786214 -0.00202771 +0.60601 -0.0631565 -0.0289751 0.00244785 -1.28842 -0.173117 0.0900511 -0.0118975 0.0 0.0 0.0 0.0 # choice of UCC exchange coefficient : 1 for Cx1, 2 for Cx2, 3 for Cx1*Cx2 -2 +1 diff --git a/input/methods b/input/methods index 02ac24b..2e3d3d2 100644 --- a/input/methods +++ b/input/methods @@ -1,5 +1,5 @@ # RHF UHF KS MOM - T F F F + F F T F # MP2* MP3 MP2-F12 F F F # CCD pCCD DCD CCSD CCSD(T) @@ -13,7 +13,7 @@ # G0F2* evGF2* qsGF2* G0F3 evGF3 F F F F F # G0W0* evGW* qsGW* ufG0W0 ufGW - F F T F F + F F F F F # G0T0 evGT qsGT F F F # MCMP2 diff --git a/src/eDFT/unrestricted_individual_energy.f90 b/src/eDFT/unrestricted_individual_energy.f90 index 3000e85..d8cd907 100644 --- a/src/eDFT/unrestricted_individual_energy.f90 +++ b/src/eDFT/unrestricted_individual_energy.f90 @@ -177,7 +177,9 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered end do end if - + + print*,"LZ shift:",sum(LZH(:)) + sum(LZx(:)) + sum(LZc(:)),"au" + ! do iEns=1,nEns ! E(iEns) = sum(ET(:,iEns)) + sum(EV(:,iEns)) & ! + sum(EH(:,iEns)) + sum(Ex(:,iEns)) + sum(Ec(:,iEns)) &