fix homo and lumo

input/dft

@@ -19,18 +19,18 @@
 # Number of states in ensemble (nEns)
   3
 # occupation numbers of orbitals
-1 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  
 1 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 
 
-1 1 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 
 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 
 # Ensemble weights: wEns(1),...,wEns(nEns-1)
   0.00 0.00 
 # Ncentered ? 0 for NO 
- 1 
+ 0 
 # Parameters for CC weight-dependent exchange functional
 0.420431 0.069097 -0.295049
 0.135075 -0.00770826 -0.028057

src/eDFT/eDFT.f90

@@ -239,7 +239,7 @@ program eDFT
 
     call cpu_time(start_KS)
     call eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight(:),maxSCF,thresh,max_diis,guess_type, & 
-                  nBas,AO,dAO,nO,nV,S,T,V,Hc,ERI,X,ENuc,Ew,occnum,Cx_choice)
+                  nBas,AO,dAO,S,T,V,Hc,ERI,X,ENuc,Ew,occnum,Cx_choice)
     call cpu_time(end_KS)
 
     t_KS = end_KS - start_KS

src/eDFT/eDFT_UKS.f90

@@ -1,5 +1,5 @@
 subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,maxSCF,thresh,max_diis,guess_type, &
-                   nBas,AO,dAO,nO,nV,S,T,V,Hc,ERI,X,ENuc,Ew,occnum,Cx_choice)
+                   nBas,AO,dAO,S,T,V,Hc,ERI,X,ENuc,Ew,occnum,Cx_choice)
 
 ! Perform unrestricted Kohn-Sham calculation for ensembles
 
@@ -22,7 +22,6 @@ subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weig
   double precision,intent(in)   :: AO(nBas,nGrid)
   double precision,intent(in)   :: dAO(ncart,nBas,nGrid)
 
-  integer,intent(in)            :: nO(nspin),nV(nspin)
   double precision,intent(in)   :: S(nBas,nBas)
   double precision,intent(in)   :: T(nBas,nBas)
   double precision,intent(in)   :: V(nBas,nBas)
@@ -372,17 +371,15 @@ subroutine eDFT_UKS(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,aCC_w1,aCC_w2,nGrid,weig
     end do
 !!!!!
 
-
-
 ! Compute final KS energy
 
-  call print_UKS(nBas,nO,eps,c,ENuc,ET,EV,EJ,Ex,Ec,Ew)
+  call print_UKS(nBas,nEns,occnum,wEns,eps,c,ENuc,ET,EV,EJ,Ex,Ec,Ew)
 
 !------------------------------------------------------------------------
 ! Compute individual energies from ensemble energy
 !------------------------------------------------------------------------
 
   call unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas, &
-                                      AO,dAO,nO,nV,T,V,ERI,ENuc,eps,Pw,rhow,drhow,J,Fx,FxHF,Fc,P,rho,drho,Ew,E,Om,occnum,Cx_choice)
+                                      AO,dAO,T,V,ERI,ENuc,eps,Pw,rhow,drhow,J,Fx,FxHF,Fc,P,rho,drho,Ew,E,Om,occnum,Cx_choice)
 
 end subroutine eDFT_UKS

src/eDFT/print_UKS.f90

@@ -1,12 +1,16 @@
-subroutine print_UKS(nBas,nO,eps,c,ENuc,ET,EV,EJ,Ex,Ec,Ew)
+subroutine print_UKS(nBas,nEns,occnum,wEns,eps,c,ENuc,ET,EV,EJ,Ex,Ec,Ew)
 
 ! Print one- and two-electron energies and other stuff for KS calculation
 
   implicit none
   include 'parameters.h'
 
+! Input variables
+
   integer,intent(in)                 :: nBas
-  integer,intent(in)                 :: nO(nspin)
+  integer,intent(in)                 :: nEns
+  double precision,intent(in)        :: occnum(nBas,nspin,nEns)
+  double precision,intent(in)        :: wEns(nEns)
   double precision,intent(in)        :: eps(nBas,nspin)
   double precision,intent(in)        :: c(nBas,nBas,nspin)
   double precision,intent(in)        :: ENuc
@@ -17,30 +21,35 @@ subroutine print_UKS(nBas,nO,eps,c,ENuc,ET,EV,EJ,Ex,Ec,Ew)
   double precision,intent(in)        :: Ec(nsp)
   double precision,intent(in)        :: Ew
 
+! Local variables
+
   integer                            :: ispin
+  integer                            :: iEns
+  integer                            :: iBas
   integer                            :: HOMO(nspin)
   integer                            :: LUMO(nspin)
   double precision                   :: Gap(nspin)
+  double precision                   :: nEl(nspin)
+
+! Number of electrons in the ensemble
+
+  nEl(:) = 0d0
+  do ispin=1,nspin
+    do iEns=1,nEns
+      do iBas=1,nBas
+        nEl(ispin) = nEl(ispin) + wEns(iEns)*occnum(iBas,ispin,iEns)
+      end do
+    end do
+  end do
 
 ! HOMO and LUMO
 
-
   do ispin=1,nspin
 
-    if(nO(ispin) > 0) then 
+      HOMO(ispin) = ceiling(nEl(ispin))
-
-      HOMO(ispin) = nO(ispin)
       LUMO(ispin) = HOMO(ispin) + 1
       Gap(ispin)  = eps(LUMO(ispin),ispin) - eps(HOMO(ispin),ispin)
 
-    else
-
-      HOMO(ispin) = 0
-      LUMO(ispin) = 0
-      Gap(ispin)  = 0d0
-
-    end if
-
   end do
 
 ! Dump results

src/eDFT/unrestricted_auxiliary_energy.f90

@@ -1,4 +1,4 @@
-subroutine unrestricted_auxiliary_energy(nBas,nEns,nO,eps,Eaux,occnum)
+subroutine unrestricted_auxiliary_energy(nBas,nEns,eps,Eaux,occnum)
 
 ! Compute the auxiliary KS energies 
 
@@ -10,7 +10,6 @@ subroutine unrestricted_auxiliary_energy(nBas,nEns,nO,eps,Eaux,occnum)
 
   integer,intent(in)            :: nBas
   integer,intent(in)            :: nEns
-  integer,intent(in)            :: nO(nspin)
   double precision,intent(in)   :: eps(nBas,nspin)
   double precision,intent(in)   :: occnum(nBas,nspin,nEns)
 

src/eDFT/unrestricted_individual_energy.f90

@@ -1,5 +1,5 @@
 subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered,nEns,wEns,aCC_w1,aCC_w2,nGrid,weight,nBas,AO,dAO, &
-                                          nO,nV,T,V,ERI,ENuc,eps,Pw,rhow,drhow,J,Fx,FxHF,Fc,P,rho,drho,Ew,E,Om,occnum,Cx_choice)
+                                          T,V,ERI,ENuc,eps,Pw,rhow,drhow,J,Fx,FxHF,Fc,P,rho,drho,Ew,E,Om,occnum,Cx_choice)
 
 ! Compute unrestricted individual energies as well as excitation energies
 
@@ -21,7 +21,6 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
   double precision,intent(in)   :: AO(nBas,nGrid)
   double precision,intent(in)   :: dAO(ncart,nBas,nGrid)
 
-  integer,intent(in)            :: nO(nspin),nV(nspin)
   double precision,intent(in)   :: T(nBas,nBas)
   double precision,intent(in)   :: V(nBas,nBas)
   double precision,intent(in)   :: ERI(nBas,nBas,nBas,nBas)
@@ -185,7 +184,7 @@ subroutine unrestricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,LDA_centered
 ! Compute auxiliary energies
 !------------------------------------------------------------------------
  
-  call unrestricted_auxiliary_energy(nBas,nEns,nO,eps,Eaux,occnum)
+  call unrestricted_auxiliary_energy(nBas,nEns,eps,Eaux,occnum)
 
 !------------------------------------------------------------------------
 ! Compute derivative discontinuities