4
1
mirror of https://github.com/pfloos/quack synced 2024-12-22 20:35:36 +01:00

fix probleme HOMO/LUMO

This commit is contained in:
Pierre-Francois Loos 2022-01-06 11:18:34 +01:00
parent 310f4ab132
commit 94a4a2937c
2 changed files with 47 additions and 28 deletions

View File

@ -13,16 +13,16 @@
# GGA = 2: LYP,PBE # GGA = 2: LYP,PBE
# MGGA = 3: # MGGA = 3:
# Hybrid = 4: HF,B3LYP,PBE # Hybrid = 4: HF,B3LYP,PBE
1 VWN5 0 H
# quadrature grid SG-n # quadrature grid SG-n
1 1
# Number of states in ensemble (nEns) # Number of states in ensemble (nEns)
2 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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 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
@ -31,9 +31,9 @@
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 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)
0.55 0.0 0.0 1.0 0.0 0.0
# Ncentered ? # Ncentered ?
T F
# Parameters for CC weight-dependent exchange functional # Parameters for CC weight-dependent exchange functional
4 4
0.642674 -0.07818 -0.0280307 0.00144198 0.642674 -0.07818 -0.0280307 0.00144198

View File

@ -29,11 +29,13 @@ subroutine print_UKS(nBas,nEns,occnum,Ov,wEns,eps,c,ENuc,ET,EV,EH,Ex,Ec,Ew,dipol
integer :: ispin integer :: ispin
integer :: iEns integer :: iEns
integer :: iBas integer :: iBas
integer :: HOMO(nspin) integer :: iHOMOa,iHOMOb
integer :: LUMO(nspin) integer :: iLUMOa,iLUMOb
double precision :: Gap(nspin) double precision :: HOMOa,HOMOb,HOMO
double precision :: S_exact,S2_exact double precision :: LUMOa,LUMOb,LUMO
double precision :: S,S2 double precision :: Gapa,Gapb,Gap
! double precision :: S_exact,S2_exact
! double precision :: S,S2
double precision :: nO(nspin) double precision :: nO(nspin)
@ -48,21 +50,34 @@ subroutine print_UKS(nBas,nEns,occnum,Ov,wEns,eps,c,ENuc,ET,EV,EH,Ex,Ec,Ew,dipol
end do end do
end do end do
print*,'nO = ',nO
! HOMO and LUMO ! HOMO and LUMO
do ispin=1,nspin iHOMOa = ceiling(nO(1))
iLUMOa = iHOMOa + 1
HOMO(ispin) = ceiling(nO(ispin)) iHOMOb = ceiling(nO(2))
LUMO(ispin) = HOMO(ispin) + 1 iLUMOb = iHOMOb + 1
Gap(ispin) = eps(LUMO(ispin),ispin) - eps(HOMO(ispin),ispin)
end do HOMOa = -huge(0d0)
if(iHOMOa > 0) HOMOa = eps(iHOMOa,1)
LUMOa = +huge(0d0)
if(iLUMOa <= nBas) LUMOa = eps(iLUMOa,1)
HOMOb = -huge(0d0)
if(iHOMOb > 0) HOMOb = eps(iHOMOb,2)
LUMOb = +huge(0d0)
if(iLUMOb <= nBas) LUMOb = eps(iLUMOb,1)
HOMO = max(HOMOa,HOMOb)
LUMO = min(LUMOa,LUMOb)
Gapa = LUMOa - HOMOa
Gapb = LUMOb - HOMOb
Gap = LUMO - HOMO
! Spin comtamination ! Spin comtamination
! S2_exact = dble(nO(1) - nO(2))/2d0*(dble(nO(1) - nO(2))/2d0 + 1d0) ! S2_exact = (nO(1) - nO(2))/2d0*(nO(1) - nO(2))/2d0 + 1d0
! S2 = S2_exact + nO(2) - sum(matmul(transpose(c(:,1:nO(1),1)),matmul(Ov,c(:,1:nO(2),2)))**2) ! S2 = S2_exact + nO(2) - sum(matmul(transpose(c(:,1:nO(1),1)),matmul(Ov,c(:,1:nO(2),2)))**2)
! S_exact = 0.5d0*dble(nO(1) - nO(2)) ! S_exact = 0.5d0*dble(nO(1) - nO(2))
@ -104,19 +119,23 @@ subroutine print_UKS(nBas,nEns,occnum,Ov,wEns,eps,c,ENuc,ET,EV,EH,Ex,Ec,Ew,dipol
write(*,'(A40,1X,F16.10,A3)') ' Nuclear repulsion: ',ENuc,' au' write(*,'(A40,1X,F16.10,A3)') ' Nuclear repulsion: ',ENuc,' au'
write(*,'(A40,1X,F16.10,A3)') ' Kohn-Sham energy: ',Ew + ENuc,' au' write(*,'(A40,1X,F16.10,A3)') ' Kohn-Sham energy: ',Ew + ENuc,' au'
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
write(*,'(A40,F13.6,A3)') ' KS HOMO a energy:',eps(HOMO(1),1)*HatoeV,' eV' write(*,'(A40,F13.6,A3)') ' KS HOMO a energy:',HOMOa*HatoeV,' eV'
write(*,'(A40,F13.6,A3)') ' KS LUMO a energy:',eps(LUMO(1),1)*HatoeV,' eV' write(*,'(A40,F13.6,A3)') ' KS LUMO a energy:',LUMOa*HatoeV,' eV'
write(*,'(A40,F13.6,A3)') ' KS HOMOa-LUMOa gap:',Gap(1)*HatoeV,' eV' write(*,'(A40,F13.6,A3)') ' KS HOMOa-LUMOa gap:',Gapa*HatoeV,' eV'
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
write(*,'(A40,F13.6,A3)') ' KS HOMO b energy:',eps(HOMO(2),2)*HatoeV,' eV' write(*,'(A40,F13.6,A3)') ' KS HOMO b energy:',HOMOb*HatoeV,' eV'
write(*,'(A40,F13.6,A3)') ' KS LUMO b energy:',eps(LUMO(2),2)*HatoeV,' eV' write(*,'(A40,F13.6,A3)') ' KS LUMO b energy:',LUMOb*HatoeV,' eV'
write(*,'(A40,F13.6,A3)') ' KS HOMOb-LUMOb gap :',Gap(2)*HatoeV,' eV' write(*,'(A40,F13.6,A3)') ' KS HOMOb-LUMOb gap :',Gapb*HatoeV,' eV'
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
write(*,'(A40,1X,F16.6)') ' S (exact) :',2d0*S_exact + 1d0 write(*,'(A40,F13.6,A3)') ' KS HOMO energy:',HOMO*HatoeV,' eV'
write(*,'(A40,1X,F16.6)') ' S :',2d0*S + 1d0 write(*,'(A40,F13.6,A3)') ' KS LUMO energy:',LUMO*HatoeV,' eV'
write(*,'(A40,1X,F16.6)') ' <S**2> (exact) :',S2_exact write(*,'(A40,F13.6,A3)') ' KS HOMO -LUMO gap :',Gap*HatoeV,' eV'
write(*,'(A40,1X,F16.6)') ' <S**2> :',S2
write(*,'(A60)') '-------------------------------------------------' write(*,'(A60)') '-------------------------------------------------'
! write(*,'(A40,1X,F16.6)') ' S (exact) :',2d0*S_exact + 1d0
! write(*,'(A40,1X,F16.6)') ' S :',2d0*S + 1d0
! write(*,'(A40,1X,F16.6)') ' <S**2> (exact) :',S2_exact
! write(*,'(A40,1X,F16.6)') ' <S**2> :',S2
! write(*,'(A60)') '-------------------------------------------------'
write(*,'(A45)') ' Dipole moment (Debye) ' write(*,'(A45)') ' Dipole moment (Debye) '
write(*,'(19X,4A10)') 'X','Y','Z','Tot.' write(*,'(19X,4A10)') 'X','Y','Z','Tot.'
write(*,'(19X,4F10.6)') (dipole(ixyz)*auToD,ixyz=1,ncart),norm2(dipole)*auToD write(*,'(19X,4F10.6)') (dipole(ixyz)*auToD,ixyz=1,ncart),norm2(dipole)*auToD