diff --git a/src/eDFT/individual_energy.f90 b/src/eDFT/individual_energy.f90
index 025f21a..cefa131 100644
--- a/src/eDFT/individual_energy.f90
+++ b/src/eDFT/individual_energy.f90
@@ -44,7 +44,6 @@ subroutine individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGrid,weight,nB
   double precision              :: EJ(nsp,nEns)
   double precision              :: Ex(nspin,nEns)
   double precision              :: Ec(nsp,nEns)
-  double precision              :: EcLZ(nsp)
   double precision              :: EcDD(nsp,nEns) 
 
   double precision,external     :: trace_matrix
@@ -118,12 +117,6 @@ subroutine individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGrid,weight,nB
 
   end do
 
-!------------------------------------------------------------------------
-! Compute Levy-Zahariev shift
-!------------------------------------------------------------------------
-
-  call correlation_Levy_Zahariev_shift(c_rung,c_DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:,:),drho(:,:,:,:),EcLZ(:))
-
 !------------------------------------------------------------------------
 ! Compute derivative discontinuities
 !------------------------------------------------------------------------
@@ -136,7 +129,7 @@ subroutine individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGrid,weight,nB
 
   do iEns=1,nEns
     E(iEns) = ENuc + sum(ET(:,iEns)) + sum(EV(:,iEns)) + sum(EJ(:,iEns)) &
-                   + sum(Ex(:,iEns)) + sum(Ec(:,iEns)) + sum(EcLZ(:)) + sum(EcDD(:,iEns))
+                   + sum(Ex(:,iEns)) + sum(Ec(:,iEns)) + sum(EcDD(:,iEns))
   end do
 
 !------------------------------------------------------------------------
@@ -151,6 +144,6 @@ subroutine individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGrid,weight,nB
 ! Dump results
 !------------------------------------------------------------------------
 
-  call print_individual_energy(nEns,EJ,Ex,Ec,EcLZ,EcDD,E,Om)
+  call print_individual_energy(nEns,EJ,Ex,Ec,EcDD,E,Om)
 
 end subroutine individual_energy
diff --git a/src/eDFT/print_individual_energy.f90 b/src/eDFT/print_individual_energy.f90
index c82a44e..5dd4919 100644
--- a/src/eDFT/print_individual_energy.f90
+++ b/src/eDFT/print_individual_energy.f90
@@ -1,4 +1,4 @@
-subroutine print_individual_energy(nEns,EJ,Ex,Ec,EcLZ,EcDD,E,Om)
+subroutine print_individual_energy(nEns,EJ,Ex,Ec,EcDD,E,Om)
 
 ! Print individual energies for eDFT calculation
 
@@ -11,7 +11,6 @@ subroutine print_individual_energy(nEns,EJ,Ex,Ec,EcLZ,EcDD,E,Om)
   double precision,intent(in)        :: EJ(nsp,nEns)
   double precision,intent(in)        :: Ex(nspin,nEns)
   double precision,intent(in)        :: Ec(nsp,nEns)
-  double precision,intent(in)        :: EcLZ(nsp)
   double precision,intent(in)        :: EcDD(nsp,nEns)
   double precision,intent(in)        :: E(nEns)
   double precision,intent(in)        :: Om(nEns)
@@ -59,15 +58,6 @@ subroutine print_individual_energy(nEns,EJ,Ex,Ec,EcLZ,EcDD,E,Om)
   write(*,'(A60)')           '-------------------------------------------------'
   write(*,*)
 
-!------------------------------------------------------------------------
-! Compute Levy-Zahariev shift
-!------------------------------------------------------------------------
-
-  write(*,'(A60)')              '-------------------------------------------------'
-  write(*,'(A40,2X,2X,F16.10,A3)') ' Levy-Zahariev shifts: ',sum(EcLZ(:)),' au'
-  write(*,'(A60)')              '-------------------------------------------------'
-  write(*,*)
-
 !------------------------------------------------------------------------
 ! Compute derivative discontinuities
 !------------------------------------------------------------------------
diff --git a/src/eDFT/print_restricted_individual_energy.f90 b/src/eDFT/print_restricted_individual_energy.f90
index eb1657f..45b8a3e 100644
--- a/src/eDFT/print_restricted_individual_energy.f90
+++ b/src/eDFT/print_restricted_individual_energy.f90
@@ -1,4 +1,4 @@
-subroutine print_restricted_individual_energy(nEns,ET,EV,EJ,Ex,Ec,Exc,ExLZ,EcLZ,ExcLZ,ExDD,EcDD,ExcDD,E,Om)
+subroutine print_restricted_individual_energy(nEns,ET,EV,EJ,Ex,Ec,Exc,ExDD,EcDD,ExcDD,E,Om)
 
 ! Print individual energies for eDFT calculation
 
@@ -12,7 +12,6 @@ subroutine print_restricted_individual_energy(nEns,ET,EV,EJ,Ex,Ec,Exc,ExLZ,EcLZ,
   double precision,intent(in)        :: EV(nEns)
   double precision,intent(in)        :: EJ(nEns)
   double precision,intent(in)        :: Ex(nEns),Ec(nEns),Exc(nEns)
-  double precision,intent(in)        :: ExLZ,EcLZ,ExcLZ
   double precision,intent(in)        :: ExDD(nEns),EcDD(nEns),ExcDD(nEns)
   double precision,intent(in)        :: E(nEns)
   double precision,intent(in)        :: Om(nEns)
@@ -86,17 +85,6 @@ subroutine print_restricted_individual_energy(nEns,ET,EV,EJ,Ex,Ec,Exc,ExLZ,EcLZ,
   write(*,'(A60)')           '-------------------------------------------------'
   write(*,*)
 
-!------------------------------------------------------------------------
-! Compute Levy-Zahariev shift
-!------------------------------------------------------------------------
-
-  write(*,'(A60)')              '-------------------------------------------------'
-  write(*,'(A40,2X,2X,F16.10,A3)') '  x Levy-Zahariev shifts: ',ExLZ, ' au'
-  write(*,'(A40,2X,2X,F16.10,A3)') '  c Levy-Zahariev shifts: ',EcLZ, ' au'
-  write(*,'(A40,2X,2X,F16.10,A3)') ' xc Levy-Zahariev shifts: ',ExcLZ,' au'
-  write(*,'(A60)')              '-------------------------------------------------'
-  write(*,*)
-
 !------------------------------------------------------------------------
 ! Compute derivative discontinuities
 !------------------------------------------------------------------------
diff --git a/src/eDFT/restricted_individual_energy.f90 b/src/eDFT/restricted_individual_energy.f90
index 7145834..7e58b82 100644
--- a/src/eDFT/restricted_individual_energy.f90
+++ b/src/eDFT/restricted_individual_energy.f90
@@ -43,7 +43,6 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGri
   double precision              :: EV(nEns)
   double precision              :: EJ(nEns)
   double precision              :: Ex(nEns),Ec(nEns),Exc(nEns)
-  double precision              :: ExLZ,EcLZ,ExcLZ
   double precision              :: ExDD(nEns),EcDD(nEns),ExcDD(nEns)
 
   double precision,external     :: trace_matrix
@@ -102,16 +101,6 @@ subroutine restricted_individual_energy(x_rung,x_DFA,c_rung,c_DFA,nEns,wEns,nGri
 
   end do
 
-!------------------------------------------------------------------------
-! Compute Levy-Zahariev shift
-!------------------------------------------------------------------------
-
-! call exchange_Levy_Zahariev_shift(x_rung,x_DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:),drho(:,:,:),ExLZ)
-
-! call restricted_correlation_Levy_Zahariev_shift(c_rung,c_DFA,nEns,wEns(:),nGrid,weight(:),rho(:,:),drho(:,:,:),EcLZ)
-
-! ExcLZ = ExLZ + ExLZ
-
 !------------------------------------------------------------------------
 ! Compute derivative discontinuities
 !------------------------------------------------------------------------