4
1
mirror of https://github.com/pfloos/quack synced 2024-11-07 06:33:55 +01:00

commit b4 indiv energy

This commit is contained in:
Pierre-Francois Loos 2021-11-29 13:11:01 +01:00
parent 73f1ffaad3
commit ed05965399
6 changed files with 51 additions and 44 deletions

View File

@ -13,17 +13,17 @@
# GGA = 2: LYP,PBE # GGA = 2: LYP,PBE
# MGGA = 3: # MGGA = 3:
# Hybrid = 4: HF,B3LYP,PBE # Hybrid = 4: HF,B3LYP,PBE
0 H 1 VWN5
# quadrature grid SG-n # quadrature grid SG-n
1 0
# Number of states in ensemble (nEns) # Number of states in ensemble (nEns)
4 2
# occupation numbers # 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 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 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 1 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
@ -31,13 +31,13 @@
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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) # Ensemble weights: wEns(1),...,wEns(nEns-1)
1 0.0 0.0 0.75 0.0 0.0
# Ncentered ? # Ncentered ?
F F
# Parameters for CC weight-dependent exchange functional # Parameters for CC weight-dependent exchange functional
4 4
0.60601 -0.0631565 -0.0289751 0.00244785 0.0 0.0 0.0 0.0
-1.28839 -0.179261 0.105627 -0.0215269 0.0 0.0 0.0 0.0
0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
# choice of UCC exchange coefficient : 1 for Cx1, 2 for Cx2, 3 for Cx1*Cx2 # choice of UCC exchange coefficient : 1 for Cx1, 2 for Cx2, 3 for Cx1*Cx2
1 1

View File

@ -73,6 +73,8 @@ subroutine UCC_lda_exchange_potential(nEns,wEns,nCC,aCC,nGrid,weight,nBas,AO,rho
w1 = wEns(2) w1 = wEns(2)
Fx1 = 1d0 - w1*(1d0 - w1)*(a1 + b1*(w1 - 0.5d0) + c1*(w1 - 0.5d0)**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) Fx2 = 1d0 - w2*(1d0 - w2)*(a2 + b2*(w2 - 0.5d0) + c2*(w2 - 0.5d0)**2)
endif endif

View File

@ -177,7 +177,7 @@ subroutine UVWN5_lda_correlation_individual_energy(nGrid,weight,rhow,rho,doNcent
decdr = 0d0 decdr = 0d0
if(ra > threshold) decdr = decdr + decdra if(ra > threshold) decdr = decdr + decdra
if(rb > threshold) decdr = decdr + decdrb ! if(rb > threshold) decdr = decdr + decdrb
Ecrr(2) = Ecrr(2) - weight(iG)*decdr*r*r Ecrr(2) = Ecrr(2) - weight(iG)*decdr*r*r

View File

@ -149,39 +149,38 @@ subroutine print_unrestricted_individual_energy(nEns,ENuc,Ew,ET,EV,EJ,Ex,Ec,Exc,
! Total Energy and IP and EA ! Total Energy and IP and EA
!------------------------------------------------------------------------ !------------------------------------------------------------------------
! write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! write(*,'(A60)') ' IP AND EA FROM AUXILIARY ENERGIES ' write(*,'(A60)') ' ENERGY DIFFERENCES FROM AUXILIARY ENERGIES '
! write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! do iEns=2,nEns do iEns=2,nEns
! write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Omaux(iEns)+OmxcDD(iEns),' au' write(*,'(A40,I2,A1,F16.10,A3)') ' Energy difference 1 -> ',iEns,':',Omaux(iEns)+OmxcDD(iEns),' au'
! write(*,*) write(*,*)
! write(*,'(A44, F16.10,A3)') ' auxiliary energy contribution : ',Omaux(iEns), ' au' 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)') ' x ensemble derivative : ',OmxDD(iEns), ' au'
! write(*,'(A44, F16.10,A3)') ' c ensemble derivative : ',OmcDD(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(*,'(A44, F16.10,A3)') ' xc ensemble derivative : ',OmxcDD(iEns),' au'
! write(*,*) 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(*,'(A60)') '-------------------------------------------------'
! write(*,*)
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
write(*,'(A60)') ' IP and EA FROM INDIVIDUAL ENERGIES '
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(*,'(A60)') ' ENERGY DIFFERENCES FROM INDIVIDUAL ENERGIES '
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
do iEns=1,nEns 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 end do
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'

View File

@ -400,7 +400,7 @@ subroutine read_options_dft(nBas,method,x_rung,x_DFA,c_rung,c_DFA,SGn,nEns,wEns,
! Read parameters for weight-dependent functional ! Read parameters for weight-dependent functional
read(1,*) read(1,*)
read(1,*) nCC read(1,*) nCC
do iEns=2,nEns do iEns=2,maxEns
read(1,*) (aCC(iCC,iEns-1),iCC=1,nCC) read(1,*) (aCC(iCC,iEns-1),iCC=1,nCC)
end do end do
@ -410,12 +410,13 @@ subroutine read_options_dft(nBas,method,x_rung,x_DFA,c_rung,c_DFA,SGn,nEns,wEns,
write(*,*)'----------------------------------------------------------' write(*,*)'----------------------------------------------------------'
write(*,*)' Parameters for weight-dependent exchange functional ' write(*,*)' Parameters for weight-dependent exchange functional '
do iEns=2,maxEns
write(*,*)'----------------------------------------------------------' write(*,*)'----------------------------------------------------------'
do iEns=2,nEns write(*,'(A10,I2,A2)') ' State ',iEns,':'
write(*,'(A6,I2,A2)') 'State ',iEns,':' write(*,*)'----------------------------------------------------------'
do iCC=1,nCC write(*,*)
write(*,'(I2,F10.6)') iCC,aCC(iCC,iEns-1) call matout(nCC,1,acc(:,iEns-1))
end do write(*,*)
end do end do
write(*,*) write(*,*)

View File

@ -26,7 +26,12 @@ subroutine unrestricted_auxiliary_energy(nBas,nEns,eps,occnum,doNcentered,Eaux)
double precision,intent(out) :: Eaux(nspin,nEns) double precision,intent(out) :: Eaux(nspin,nEns)
! Memory allocation
allocate(nEl(nEns)) allocate(nEl(nEns))
! Compute the number of electrons
nEl(:) = 0d0 nEl(:) = 0d0
do iEns=1,nEns do iEns=1,nEns
do iBas=1,nBas do iBas=1,nBas