LCPQ/quack
4
1
Fork 0
mirror of https://github.com/pfloos/quack synced 2025-01-03 01:56:09 +01:00
Code Releases Activity

Compare commits

base: LCPQ:372ccd2e4ea31f8cee4fa09e2c88ab58f98acbbd
Branches Tags
LCPQ:master
LCPQ:refact
LCPQ:trexio
LCPQ:old
LCPQ:add-license-1
LCPQ:0.9
...
compare: LCPQ:1e428400c461c3eafa7d780bac6f5e66602852a5
Branches Tags
LCPQ:master
LCPQ:refact
LCPQ:trexio
LCPQ:old
LCPQ:add-license-1
LCPQ:0.9

4 Commits
372ccd2e4e ... 1e428400c4

Author SHA1 Message Date
Pierre-Francois Loos
1e428400c4 scan w for cloclo 2020-07-06 16:07:22 +02:00
Pierre-Francois Loos
ec22768aa8 still cleaning up 2020-07-06 15:43:18 +02:00
Pierre-Francois Loos
2c2f02d701 clean up 2020-07-06 15:42:21 +02:00
Clotilde Marut
b5f7e5e6c9 testing some modifications 2020-07-06 12:19:49 +02:00
21 changed files with 119 additions and 99 deletions
Show Stats Expand all files Collapse all files

30
input/basis
View File

@ -1,25 +1,9 @@
1 10
S 4
1 528.5000000 0.0009400
2 79.3100000 0.0072140
3 18.0500000 0.0359750
4 5.0850000 0.1277820
1 3
S 3
1 38.3600000 0.0238090
2 5.7700000 0.1548910
3 1.2400000 0.4699870
S 1
1 1.6090000 1.0000000
S 1
1 0.5363000 1.0000000
S 1
1 0.1833000 1.0000000
1 0.2976000 1.0000000
P 1
1 5.9940000 1.0000000
P 1
1 1.7450000 1.0000000
P 1
1 0.5600000 1.0000000
D 1
1 4.2990000 1.0000000
D 1
1 1.2230000 1.0000000
F 1
1 2.6800000 1.0000000
1 1.2750000 1.0000000

12
input/dft
View File

@ -5,23 +5,23 @@
# LDA = 1: RS51,RMFL20
# GGA = 2: RB88
# Hybrid = 4
# Hartree-Fock = 666: RHF,UHF
1 US51
# Hartree-Fock = 666
1 US51
# correlation rung:
# Hartree = 0
# LDA = 1: RVWN5,RMFL20
# GGA = 2:
# Hybrid = 4:
# Hartree-Fock = 666
0 H
0 H
# quadrature grid SG-n
1
# Number of states in ensemble (nEns)
3
# Ensemble weights: wEns(1),...,wEns(nEns-1)
0.000 0.000
1 0.0
# Parameters for CC weight-dependent exchange functional
0.420243 0.0700561 -0.288301
0.135068 -0.00774769 -0.0278205
0.000000 0.0000000 0.000000
0.000000 0.0000000 0.0000000
# GOK-DFT: maxSCF thresh DIIS n_diis guess_type ortho_type
32 0.00001 T 5 1 1

30
input/weight
View File

@ -1,25 +1,9 @@
1 10
S 4
1 528.5000000 0.0009400
2 79.3100000 0.0072140
3 18.0500000 0.0359750
4 5.0850000 0.1277820
1 3
S 3
1 38.3600000 0.0238090
2 5.7700000 0.1548910
3 1.2400000 0.4699870
S 1
1 1.6090000 1.0000000
S 1
1 0.5363000 1.0000000
S 1
1 0.1833000 1.0000000
1 0.2976000 1.0000000
P 1
1 5.9940000 1.0000000
P 1
1 1.7450000 1.0000000
P 1
1 0.5600000 1.0000000
D 1
1 4.2990000 1.0000000
D 1
1 1.2230000 1.0000000
F 1
1 2.6800000 1.0000000
1 1.2750000 1.0000000

51
scripts/scan_w.sh Executable file
View File

@ -0,0 +1,51 @@
#! /bin/bash
MOL=$1
BASIS=$2
w_start=0.0
w_end=1.05
dw=0.05
w2=0.0
XF=$3
CF=$4
aw1="0.000000 0.0000000 0.000000"
aw2="0.000000 0.0000000 0.0000000"
for w1 in $(seq $w_start $dw $w_end)
do
### w2=${w1}
echo "Weights = " $w1 $w2
echo "# Restricted or unrestricted KS calculation" > input/dft
echo " eDFT-UKS" >> input/dft
echo "# exchange rung:" >> input/dft
echo "# Hartree = 0" >> input/dft
echo "# LDA = 1: RS51,RMFL20" >> input/dft
echo "# GGA = 2: RB88" >> input/dft
echo "# Hybrid = 4" >> input/dft
echo "# Hartree-Fock = 666" >> input/dft
echo " 1 $XF " >> input/dft
echo "# correlation rung: " >> input/dft
echo "# Hartree = 0" >> input/dft
echo "# LDA = 1: RVWN5,RMFL20" >> input/dft
echo "# GGA = 2: " >> input/dft
echo "# Hybrid = 4: " >> input/dft
echo "# Hartree-Fock = 666" >> input/dft
echo " 0 $CF " >> input/dft
echo "# quadrature grid SG-n" >> input/dft
echo " 1" >> input/dft
echo "# Number of states in ensemble (nEns)" >> input/dft
echo " 3" >> input/dft
echo "# Ensemble weights: wEns(1),...,wEns(nEns-1)" >> input/dft
echo " ${w1} ${w2} " >> input/dft
echo "# Parameters for CC weight-dependent exchange functional" >> input/dft
echo ${aw1} >> input/dft
echo ${aw2} >> input/dft
echo "# GOK-DFT: maxSCF thresh DIIS n_diis guess_type ortho_type" >> input/dft
echo " 32 0.00001 T 5 1 1" >> input/dft
./GoXC $MOL $BASIS > ${MOL}_${BASIS}_${XF}_${CF}_${w1}.out
done

2
src/eDFT/RB88_gga_exchange_individual_energy.f90
View File

@ -40,7 +40,7 @@ subroutine RB88_gga_exchange_individual_energy(nGrid,weight,rhow,drhow,rho,drho,
r = max(0d0,0.5d0*rhow(iG))
rI = max(0d0,0.5d0*rho(iG))
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
g = 0.25d0*(drho(1,iG)**2 + drho(2,iG)**2 + drho(3,iG)**2)
x = sqrt(g)/r**(4d0/3d0)

2
src/eDFT/RCC_lda_exchange_individual_energy.f90
View File

@ -83,7 +83,7 @@ subroutine RCC_lda_exchange_individual_energy(nEns,wEns,aCC_w1,aCC_w2,nGrid,weig
r = max(0d0,rhow(iG))
rI = max(0d0,rho(iG))
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
e_p = Cx*r**(1d0/3d0)
dedr = 1d0/3d0*Cx*r**(-2d0/3d0)

2
src/eDFT/RMFL20_lda_exchange_individual_energy.f90
View File

@ -42,7 +42,7 @@ subroutine RMFL20_lda_exchange_individual_energy(LDA_centered,nEns,wEns,nGrid,we
r = max(0d0,rhow(iG))
rI = max(0d0,rho(iG))
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
e_p = Cxw*r**(1d0/3d0)
dedr = 1d0/3d0*Cxw*r**(-2d0/3d0)

2
src/eDFT/RVWN5_lda_correlation_individual_energy.f90
View File

@ -42,7 +42,7 @@ subroutine RVWN5_lda_correlation_individual_energy(nGrid,weight,rhow,rho,Ec)
r = max(0d0,rhow(iG))
rI = max(0d0,rho(iG))
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
rs = (4d0*pi*r/3d0)**(-1d0/3d0)
x = sqrt(rs)

5
src/eDFT/US51_lda_exchange_energy.f90
View File

@ -1,4 +1,4 @@
subroutine US51_lda_exchange_energy(nGrid,weight,rho,Ex,wEns,nEns)
subroutine US51_lda_exchange_energy(nGrid,weight,rho,Ex)
! Compute Slater's LDA exchange energy
@ -11,9 +11,6 @@ subroutine US51_lda_exchange_energy(nGrid,weight,rho,Ex,wEns,nEns)
double precision,intent(in) :: weight(nGrid)
double precision,intent(in) :: rho(nGrid)
integer,intent(in) :: nEns
double precision,intent(in) :: wEns(nEns)
! Local variables
integer :: iG

10
src/eDFT/US51_lda_exchange_individual_energy.f90
View File

@ -32,12 +32,18 @@ subroutine US51_lda_exchange_individual_energy(nGrid,weight,rhow,rho,Ex)
r = max(0d0,rhow(iG))
rI = max(0d0,rho(iG))
if(r > threshold .or. rI > threshold) then
if(r > threshold) then
e = alpha*r**(1d0/3d0)
dedr = 1d0/3d0*alpha*r**(-2d0/3d0)
Ex = Ex + weight(iG)*(e*rI + dedr*r*rI - dedr*r*r)
Ex = Ex - weight(iG)*dedr*r*r
if(rI > threshold) then
Ex = Ex + weight(iG)*(e*rI + dedr*r*rI)
endif
endif

2
src/eDFT/UVWN5_lda_correlation_energy.f90
View File

@ -73,7 +73,7 @@ subroutine UVWN5_lda_correlation_energy(nGrid,weight,rho,Ec)
! alpha-beta contribution
if(ra > threshold .and. rb > threshold) then
if(ra > threshold .or. rb > threshold) then
r = ra + rb
rs = (4d0*pi*r/3d0)**(-1d0/3d0)

2
src/eDFT/UWN5_lda_correlation_individual_energy.f90 → src/eDFT/UVWN5_lda_correlation_individual_energy.f90
View File

@ -28,7 +28,7 @@ subroutine UVWN5_lda_correlation_individual_energy(nGrid,weight,rhow,rho,Ec)
! Output variables
double precision :: Ec(nspin)
double precision :: Ec(nsp)
! Parameters of the functional

2
src/eDFT/UW38_lda_correlation_energy.f90
View File

@ -35,7 +35,7 @@ subroutine UW38_lda_correlation_energy(nGrid,weight,rho,Ec)
ra = max(0d0,rho(iG,1))
rb = max(0d0,rho(iG,2))
if(ra > threshold .and. rb > threshold) then
if(ra > threshold .or. rb > threshold) then
r = ra + rb

2
src/eDFT/UW38_lda_correlation_individual_energy.f90
View File

@ -45,7 +45,7 @@ subroutine UW38_lda_correlation_individual_energy(nGrid,weight,rhow,rho,Ec)
r = ra + rb
rI = raI + rbI
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
epsc = ra*rb/(r + d*r**(2d0/3d0))
dFcdra = epsc*(d/(3d0*r**(4d0/3d0)*(1d0 + d*r**(-1d0/3d0))) - 1d0/r + 1d0/ra)

2
src/eDFT/eDFT_UKS.f90
View File

@ -311,7 +311,7 @@ subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weig
EJ(1) = 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,1)))
EJ(2) = 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,2))) &
+ 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,1))) !!!!!!
+ 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,1)))
EJ(3) = 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,2)))
! Exchange energy

2
src/eDFT/elda_correlation_individual_energy.f90
View File

@ -39,7 +39,7 @@ subroutine elda_correlation_individual_energy(aLF,nGrid,weight,rhow,rho,Ec)
r = ra + rb
rI = raI + rbI
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
ec_p = aLF(1)/(1d0 + aLF(2)*r**(-1d0/6d0) + aLF(3)*r**(-1d0/3d0))

2
src/eDFT/lda_exchange_energy.f90
View File

@ -27,7 +27,7 @@ subroutine lda_exchange_energy(DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,we
case ('US51')
call US51_lda_exchange_energy(nGrid,weight,rho,Ex,wEns,nEns)
call US51_lda_exchange_energy(nGrid,weight,rho,Ex)
case ('RS51')

38
src/eDFT/print_unrestricted_individual_energy.f90
View File

@ -136,18 +136,18 @@ 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)') '-------------------------------------------------'
write(*,'(A40,F16.10,A3)') ' Ionization Potential 1 -> 2 :',Omaux(2)+OmxcDD(2),' au'
write(*,'(A43,F16.10,A4)') ' Ionization Potential 1 -> 2:',Omaux(2)+OmxcDD(2),' au'
write(*,*)
write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(2), ' au'
write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(2), ' au'
write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(2), ' au'
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(2),' au'
write(*,*)
write(*,'(A40,F16.10,A3)') ' Electronic Affinity 1 -> 3 :',Omaux(3)+OmxcDD(3),' au'
write(*,'(A43,F16.10,A4)') ' Electronic Affinity 1 -> 3:',Omaux(3)+OmxcDD(3),' au'
write(*,*)
write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(3), ' au'
write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(3), ' au'
@ -155,17 +155,17 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(3),' au'
write(*,*)
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') '-------------------------------------------------'
write(*,*)
write(*,'(A40,F16.10,A3)') ' Ionization Potential 1 -> 2 :',(Omaux(2)+OmxcDD(2))*HaToeV,' eV'
write(*,'(A40,F16.10,A3)') ' Ionization Potential 1 -> 2:',(Omaux(2)+OmxcDD(2))*HaToeV,' eV'
write(*,*)
write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(2)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(2)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(2)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(2)*HaToeV,' eV'
write(*,*)
write(*,'(A40,F16.10,A3)') ' Electronic Affinity 1 -> 3 :',(Omaux(3)+OmxcDD(3))*HaToeV,' eV'
write(*,'(A40,F16.10,A3)') ' Electronic Affinity 1 -> 3:',(Omaux(3)+OmxcDD(3))*HaToeV,' eV'
write(*,*)
write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(3)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' x ensemble derivative : ',OmxDD(3)*HaToeV, ' eV'
@ -173,18 +173,18 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(3)*HaToeV,' eV'
write(*,*)
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') '-------------------------------------------------'
write(*,*)
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') ' IP and EA FROM INDIVIDUAL ENERGIES '
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'
end do
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A40,F16.10,A3)') ' Ionization Potential 1 -> 2 :',Om(2), ' au'
write(*,'(A43,F16.10,A4)') ' Ionization Potential 1 -> 2:',Om(2), ' au'
write(*,*)
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(2), ' au'
write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(2), ' au'
@ -194,7 +194,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(2), ' au'
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(2),' au'
write(*,*)
write(*,'(A40,F16.10,A3)') ' Electronic Affinity 1 -> 3 :',Om(3), ' au'
write(*,'(A43,F16.10,A4)') ' Electronic Affinity 1 -> 3:',Om(3), ' au'
write(*,*)
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(3), ' au'
write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(3), ' au'
@ -205,9 +205,9 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(3),' au'
write(*,*)
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A40,F16.10,A3)') ' Ionization Potential 1 -> 2 :',Om(2)*HaToeV, ' eV'
write(*,'(A43,F16.10,A4)') ' Ionization Potential 1 -> 2:',Om(2)*HaToeV, ' eV'
write(*,*)
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(2)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(2)*HaToeV, ' eV'
@ -217,7 +217,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(2)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(2)*HaToeV,' eV'
write(*,*)
write(*,'(A40,F16.10,A3)') ' Electronic Affinity 1 -> 3 :',Om(3)*HaToeV, ' eV'
write(*,'(A43,F16.10,A4)') ' Electronic Affinity 1 -> 3:',Om(3)*HaToeV, ' eV'
write(*,*)
write(*,'(A44, F16.10,A3)') ' x energy contribution : ',Omx(3)*HaToeV, ' eV'
write(*,'(A44, F16.10,A3)') ' c energy contribution : ',Omc(3)*HaToeV, ' eV'
@ -228,7 +228,7 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
write(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(3)*HaToeV,' eV'
write(*,*)
write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') '-------------------------------------------------'
write(*,*)

2
src/eDFT/restricted_elda_correlation_individual_energy.f90
View File

@ -33,7 +33,7 @@ subroutine restricted_elda_correlation_individual_energy(aMFL,nGrid,weight,rhow,
r = max(0d0,rhow(iG))
rI = max(0d0,rho(iG))
if(r > threshold .and. rI > threshold) then
if(r > threshold .or. rI > threshold) then
ec_p = aMFL(1)/(1d0 + aMFL(2)*r**(-1d0/6d0) + aMFL(3)*r**(-1d0/3d0))

8
src/eDFT/unrestricted_correlation_energy.f90
View File

@ -38,13 +38,13 @@ subroutine unrestricted_correlation_energy(rung,DFA,nEns,wEns,nGrid,weight,rho,d
case(1)
call unrestricted_lda_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:),Ec(:))
call unrestricted_lda_correlation_energy(DFA,nEns,wEns,nGrid,weight,rho,Ec)
! GGA functionals
case(2)
call unrestricted_gga_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:),drho(:,:,:),Ec(:))
call unrestricted_gga_correlation_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,Ec)
! Hybrid functionals
@ -52,8 +52,8 @@ subroutine unrestricted_correlation_energy(rung,DFA,nEns,wEns,nGrid,weight,rho,d
aC = 0.81d0
call unrestricted_lda_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:),EcLDA(:))
call unrestricted_gga_correlation_energy(DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:),drho(:,:,:),EcGGA(:))
call unrestricted_lda_correlation_energy(DFA,nEns,wEns,nGrid,weight,rho,EcLDA)
call unrestricted_gga_correlation_energy(DFA,nEns,wEns,nGrid,weight,rho,drho,EcGGA)
Ec(:) = EcLDA(:) + aC*(EcGGA(:) - EcLDA(:))

10
src/eDFT/unrestricted_individual_energy.f90
View File

@ -42,7 +42,6 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
double precision,intent(in) :: Fc(nBas,nBas,nspin)
double precision :: Ew
! Local variables
double precision :: ET(nspin,nEns)
@ -105,9 +104,8 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
EJ(1,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,1))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,1)))
EJ(2,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,2))) &
+ trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,1))) &
! 2.0d0*trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,2))) &
EJ(2,iEns) = trace_matrix(nBas,matmul(P(:,:,1,iEns),J(:,:,2))) &
+ trace_matrix(nBas,matmul(P(:,:,2,iEns),J(:,:,1))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,1),J(:,:,2))) &
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,1)))
@ -115,7 +113,7 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
- 0.5d0*trace_matrix(nBas,matmul(Pw(:,:,2),J(:,:,2)))
end do
!------------------------------------------------------------------------
! Checking Hartree contributions for each individual states
!------------------------------------------------------------------------
@ -142,7 +140,7 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
rho(:,ispin,iEns),drho(:,:,ispin,iEns),Ex(ispin,iEns))
end do
end do
!------------------------------------------------------------------------
! Checking exchange contributions for each individual states
!------------------------------------------------------------------------