From 73f1ffaad33c0dc96226f6ad93086a56734666fe Mon Sep 17 00:00:00 2001 From: Clotilde Marut Date: Thu, 25 Nov 2021 17:45:48 +0100 Subject: [PATCH] CC functional --- input/dft | 48 +++++++------- ..._lda_exchange_derivative_discontinuity.f90 | 65 ++++++++++++------- src/eDFT/UCC_lda_exchange_energy.f90 | 55 ++++++++++------ .../UCC_lda_exchange_individual_energy.f90 | 51 ++++++++++----- src/eDFT/UCC_lda_exchange_potential.f90 | 53 +++++++++------ .../print_unrestricted_individual_energy.f90 | 52 +++++++-------- 6 files changed, 194 insertions(+), 130 deletions(-) diff --git a/input/dft b/input/dft index 85cd470..7a4ae5a 100644 --- a/input/dft +++ b/input/dft @@ -4,38 +4,40 @@ # Hartree = 0: H # LDA = 1: S51,CC-S51 # GGA = 2: B88,G96,PBE -# MGGA = 3: -# Hybrid = 4: HF,B3LYP,PBE - 1 S51 +# MGGA = 3: +# Hybrid = 4 HF,B3LYP,PBE +1 S51 # correlation rung: # Hartree = 0: H # LDA = 1: PW92,VWN3,VWN5,eVWN5 # GGA = 2: LYP,PBE # MGGA = 3: -# Hybrid = 4: HF,B3LYP,PBE - 1 VWN5 +# Hybrid = 4: HF,B3LYP,PBE +0 H # quadrature grid SG-n - 0 + 1 # Number of states in ensemble (nEns) - 2 -# 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 - 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 +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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 - 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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 - 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 +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 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 - 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 +0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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.95 0.00 0.00 -# N-centered? - F + 1 0.0 0.0 +# Ncentered ? +F # Parameters for CC weight-dependent exchange functional -3 -0.0 0.0 0.0 +4 +0.60601 -0.0631565 -0.0289751 0.00244785 +-1.28839 -0.179261 0.105627 -0.0215269 +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/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 b/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 index a8b0d1a..27fe75f 100644 --- a/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 +++ b/src/eDFT/UCC_lda_exchange_derivative_discontinuity.f90 @@ -26,8 +26,8 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei double precision,allocatable :: dExdw(:) double precision,external :: Kronecker_delta - double precision :: a1,b1,c1,w1 - double precision :: a2,b2,c2,w2 + double precision :: a1,b1,c1,d1,w1 + double precision :: a2,b2,c2,d2,w2 double precision :: dCxdw1,dCxdw2 ! Output variables @@ -44,41 +44,43 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei allocate(dExdw(nEns)) -! Parameters for first state - - a1 = aCC(1,1) - b1 = aCC(2,1) - c1 = aCC(3,1) - -! Parameters for second state - - a2 = aCC(1,2) - b2 = aCC(2,2) - c2 = aCC(3,2) - - w1 = wEns(2) - w2 = wEns(3) - ! Defining enhancements factor for weight-dependent functionals if (doNcentered) then - + +! Parameters for first state + + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + d1 = aCC(4,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + d2 = aCC(4,2) + + w1 = wEns(2) + w2 = wEns(3) + select case (Cx_choice) case(1) - dCxdw1 = 2.d0*a1*(w1-1.d0)+(2.d0+3.d0*(w1-2.d0)*w1)*b1+2.d0*(w1-1.d0)*(1.d0+2.d0*(w1-2.d0)*w1)*c1 + dCxdw1 = a1 + 2.d0*b1*w1 + 3.d0*c1*w1**2 + 4.d0*d1*w1**3 dCxdw2 = 0.d0 case(2) dCxdw1 = 0.d0 - dCxdw2 = 2.d0*a2*(w2-1.d0)+(2.d0+3.d0*(w2-2.d0)*w2)*b2+2.d0*(w2-1.d0)*(1.d0+2.d0*(w2-2.d0)*w2)*c2 + dCxdw2 = a2 + 2.d0*b2*w2 + 3.d0*c2*w2**2 + 4.d0*d2*w2**3 case(3) - dCxdw1 = (2.d0*a1*(w1-1.d0)+(2.d0+3.d0*(w1-2.d0)*w1)*b1+2.d0*(w1-1.d0)*(1.d0+2.d0*(w1-2.d0)*w1)*c1) & - * (1d0 - w2*(2d0 - w2)*(a2 + b2*(w2 - 1d0) + c2*(w2 - 1d0)**2)) + dCxdw1 = (a1 + 2.d0*b1*w1 + 3.d0*c1*w1**2 + 4.d0*d1*w1**3) & + * (1d0 + a2*w2 + b2*w2**2 + c2*w2**3 + d2*w2**4) - dCxdw2 = (1d0 - w1*(2d0 - w1)*(a1 + b1*(w1 - 1.d0) + c1*(w1 - 1.d0)**2)) & - * (2.d0*a2*(w2-1.d0)+(2.d0+3.d0*(w2-2.d0)*w2)*b2+2.d0*(w2-1.d0)*(1.d0+2.d0*(w2-2.d0)*w2)*c2) + dCxdw2 = (1d0 + a1*w1 + b1*w1**2 + c1*w1**3 + d1*w1**4) & + * (a2 + 2.d0*b2*w2 + 3.d0*c2*w2**2 + 4.d0*d2*w2**3) case default dCxdw1 = 0d0 @@ -88,6 +90,21 @@ subroutine UCC_lda_exchange_derivative_discontinuity(nEns,wEns,nCC,aCC,nGrid,wei else +! Parameters for first state + + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + + w1 = wEns(2) + w2 = wEns(3) + select case (Cx_choice) case(1) diff --git a/src/eDFT/UCC_lda_exchange_energy.f90 b/src/eDFT/UCC_lda_exchange_energy.f90 index 5bc917d..0d8bc00 100644 --- a/src/eDFT/UCC_lda_exchange_energy.f90 +++ b/src/eDFT/UCC_lda_exchange_energy.f90 @@ -22,43 +22,60 @@ subroutine UCC_lda_exchange_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho,Cx_choice, integer :: iG double precision :: r - double precision :: a1,b1,c1,w1 - double precision :: a2,b2,c2,w2 + double precision :: a1,b1,c1,d1,w1 + double precision :: a2,b2,c2,d2,w2 double precision :: Fx1,Fx2,Cx ! Output variables double precision :: Ex -! Parameters for first state - - a1 = aCC(1,1) - b1 = aCC(2,1) - c1 = aCC(3,1) - -! Parameters for second state - - a2 = aCC(1,2) - b2 = aCC(2,2) - c2 = aCC(3,2) ! Defining enhancements factor for weight-dependent functionals if(doNcentered) then + +! Parameters for first state + + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + d1 = aCC(4,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + d2 = aCC(4,2) + w1 = wEns(2) - Fx1 = 1d0 - w1*(2d0 - w1)*(a1 + b1*(w1 - 1d0) + c1*(w1 - 1d0)**2) + Fx1 = 1d0 + a1*w1 + b1*w1**2 + c1*w1**3 + d1*w1**4 w2 = wEns(3) - Fx2 = 1d0 - w2*(2d0 - w2)*(a2 + b2*(w2 - 1d0) + c2*(w2 - 1d0)**2) + Fx2 = 1d0 + a2*w2 + b2*w2**2 + c2*w2**3 + d2*w2**4 else - w1 = wEns(2) - Fx1 = 1d0 - w1*(1d0 - w1)*(a1 + b1*(w1 - 0.5d0) + c1*(w1 - 0.5d0)**2) +! Parameters for first state - w2 = wEns(3) - Fx2 = 1d0 - w2*(1d0 - w2)*(a2 + b2*(w2 - 0.5d0) + c2*(w2 - 0.5d0)**2) + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + + + w1 = wEns(2) + Fx1 = 1d0 - w1*(1d0 - w1)*(a1 + b1*(w1 - 0.5d0) + c1*(w1 - 0.5d0)**2) + + w2 = wEns(3) + Fx2 = 1d0 - w2*(1d0 - w2)*(a2 + b2*(w2 - 0.5d0) + c2*(w2 - 0.5d0)**2) endif diff --git a/src/eDFT/UCC_lda_exchange_individual_energy.f90 b/src/eDFT/UCC_lda_exchange_individual_energy.f90 index cd0666e..9977288 100644 --- a/src/eDFT/UCC_lda_exchange_individual_energy.f90 +++ b/src/eDFT/UCC_lda_exchange_individual_energy.f90 @@ -26,8 +26,8 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho double precision :: e_p,dedr double precision :: Exrr,ExrI,ExrrI - double precision :: a1,b1,c1,w1 - double precision :: a2,b2,c2,w2 + double precision :: a1,b1,c1,d1,w1 + double precision :: a2,b2,c2,d2,w2 double precision :: Fx1,Fx2,Cx ! Output variables @@ -39,30 +39,45 @@ subroutine UCC_lda_exchange_individual_energy(nEns,wEns,nCC,aCC,nGrid,weight,rho double precision,external :: electron_number -! Parameters for first state - - a1 = aCC(1,1) - b1 = aCC(2,1) - c1 = aCC(3,1) - -! Parameters for second state - - a2 = aCC(1,2) - b2 = aCC(2,2) - c2 = aCC(3,2) - ! Defining enhancements factor for weight-dependent functionals if(doNcentered) then - w1 = wEns(2) - Fx1 = 1d0 - w1*(2d0 - w1)*(a1 + b1*(w1 - 1d0) + c1*(w1 - 1d0)**2) +! Parameters for first state - w2 = wEns(3) - Fx2 = 1d0 - w2*(2d0 - w2)*(a2 + b2*(w2 - 1d0) + c2*(w2 - 1d0)**2) + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + d1 = aCC(4,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + d2 = aCC(4,2) + + + w1 = wEns(2) + Fx1 = 1d0 + a1*w1 + b1*w1**2 + c1*w1**3 + d1*w1**4 + + w2 = wEns(3) + Fx2 = 1d0 + a2*w2 + b2*w2**2 + c2*w2**3 + d2*w2**4 else +! Parameters for first state + + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + w1 = wEns(2) Fx1 = 1d0 - w1*(1d0 - w1)*(a1 + b1*(w1 - 0.5d0) + c1*(w1 - 0.5d0)**2) diff --git a/src/eDFT/UCC_lda_exchange_potential.f90 b/src/eDFT/UCC_lda_exchange_potential.f90 index 114a204..a7ca60e 100644 --- a/src/eDFT/UCC_lda_exchange_potential.f90 +++ b/src/eDFT/UCC_lda_exchange_potential.f90 @@ -24,44 +24,57 @@ subroutine UCC_lda_exchange_potential(nEns,wEns,nCC,aCC,nGrid,weight,nBas,AO,rho integer :: mu,nu,iG double precision :: r,vAO - double precision :: a1,b1,c1,w1 - double precision :: a2,b2,c2,w2 + double precision :: a1,b1,c1,d1,w1 + double precision :: a2,b2,c2,d2,w2 double precision :: Fx1,Fx2,Cx ! Output variables double precision,intent(out) :: Fx(nBas,nBas) -! Parameters for first state - - a1 = aCC(1,1) - b1 = aCC(2,1) - c1 = aCC(3,1) - -! Parameters for second state - - a2 = aCC(1,2) - b2 = aCC(2,2) - c2 = aCC(3,2) ! Defining enhancements factor for weight-dependent functionals - if(doNcentered) then + if(doNcentered) then + +! Parameters for first state + + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) + d1 = aCC(4,1) + +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + d2 = aCC(4,2) w1 = wEns(2) - Fx1 = 1d0 - w1*(2d0 - w1)*(a1 + b1*(w1 - 1d0) + c1*(w1 - 1d0)**2) + Fx1 = 1d0 + a1*w1 + b1*w1**2 + c1*w1**3 + d1*w1**4 w2 = wEns(3) - Fx2 = 1d0 - w2*(2d0 - w2)*(a2 + b2*(w2 - 1d0) + c2*(w2 - 1d0)**2) + Fx2 = 1d0 + a2*w2 + b2*w2**2 + c2*w2**3 + d2*w2**4 else - w1 = wEns(2) - Fx1 = 1d0 - w1*(1d0 - w1)*(a1 + b1*(w1 - 0.5d0) + c1*(w1 - 0.5d0)**2) +! Parameters for first state - w2 = wEns(3) - Fx2 = 1d0 - w2*(1d0 - w2)*(a2 + b2*(w2 - 0.5d0) + c2*(w2 - 0.5d0)**2) + a1 = aCC(1,1) + b1 = aCC(2,1) + c1 = aCC(3,1) +! Parameters for second state + + a2 = aCC(1,2) + b2 = aCC(2,2) + c2 = aCC(3,2) + + w1 = wEns(2) + Fx1 = 1d0 - w1*(1d0 - w1)*(a1 + b1*(w1 - 0.5d0) + c1*(w1 - 0.5d0)**2) + Fx2 = 1d0 - w2*(1d0 - w2)*(a2 + b2*(w2 - 0.5d0) + c2*(w2 - 0.5d0)**2) + endif select case (Cx_choice) diff --git a/src/eDFT/print_unrestricted_individual_energy.f90 b/src/eDFT/print_unrestricted_individual_energy.f90 index 70b076a..f75eaa0 100644 --- a/src/eDFT/print_unrestricted_individual_energy.f90 +++ b/src/eDFT/print_unrestricted_individual_energy.f90 @@ -149,39 +149,39 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc, ! Total Energy and IP and EA !------------------------------------------------------------------------ - write(*,'(A60)') '-------------------------------------------------' - write(*,'(A60)') ' IP AND EA FROM AUXILIARY ENERGIES ' - write(*,'(A60)') '-------------------------------------------------' +! write(*,'(A60)') '-------------------------------------------------' +! write(*,'(A60)') ' IP AND EA FROM AUXILIARY ENERGIES ' +! write(*,'(A60)') '-------------------------------------------------' - do iEns=2,nEns - write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Omaux(iEns)+OmxcDD(iEns),' au' - write(*,*) - write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns), ' au' - write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns), ' au' - write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns), ' au' - write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(iEns),' au' - write(*,*) +! do iEns=2,nEns +! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Omaux(iEns)+OmxcDD(iEns),' au' +! write(*,*) +! write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns), ' au' +! write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns), ' au' +! write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns), ' au' +! write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(iEns),' au' +! write(*,*) - write(*,'(A60)') '-------------------------------------------------' - write(*,*) +! write(*,'(A60)') '-------------------------------------------------' +! write(*,*) - write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',(Omaux(iEns)+OmxcDD(iEns))*HaToeV,' eV' - write(*,*) - write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns)*HaToeV, ' eV' - write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV' - write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV' - write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(iEns)*HaToeV,' eV' - write(*,*) - end do +! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',(Omaux(iEns)+OmxcDD(iEns))*HaToeV,' eV' +! write(*,*) +! write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns)*HaToeV, ' eV' +! write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV' +! write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV' +! write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(iEns)*HaToeV,' eV' +! write(*,*) +! end do - write(*,'(A60)') '-------------------------------------------------' - write(*,*) +! write(*,'(A60)') '-------------------------------------------------' +! write(*,*) - write(*,'(A60)') '-------------------------------------------------' + write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') ' IP and EA FROM INDIVIDUAL ENERGIES ' write(*,'(A60)') '-------------------------------------------------' do iEns=1,nEns - write(*,'(A40,I2,A2,F16.10,A3)') ' Individual energy state ',iEns,': ',E(iEns) + ENuc,' au' +! write(*,'(A40,I2,A2,F16.10,A3)') ' Individual energy state ',iEns,': ',E(iEns) + ENuc,' au' end do write(*,'(A60)') '-------------------------------------------------' @@ -208,7 +208,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc, write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(iEns)*HaToeV, ' eV' write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(iEns)*HaToeV,' eV' - write(*,*) + write(*,*) end do write(*,'(A60)') '-------------------------------------------------' write(*,*)