BEGIN_PROVIDER [ double precision, kinetic_energy_alpha_p ]
 implicit none
 BEGIN_DOC  
! Alpha Kinetic Energy
 END_DOC

 kinetic_energy_alpha_p = 0.d0

 integer :: i, l
 do i=1,elec_alpha_num
  do l=1,3
   kinetic_energy_alpha_p = kinetic_energy_alpha_p+mo_grad_p(i,l)**2
  enddo
 enddo

END_PROVIDER

BEGIN_PROVIDER [ double precision, kinetic_energy_beta_p ]
 implicit none
 BEGIN_DOC  
! Beta Kinetic Energy
 END_DOC

 kinetic_energy_beta_p = 0.d0

 integer :: i,l
 do i=1,elec_beta_num
  do l=1,3
   kinetic_energy_beta_p = kinetic_energy_beta_p+mo_grad_p(i,l)**2
  enddo
 enddo

END_PROVIDER

BEGIN_PROVIDER [ double precision, kinetic_energy_p ]
 implicit none
 BEGIN_DOC  
! Kinetic energy
 END_DOC
 kinetic_energy_p = kinetic_energy_alpha_p + kinetic_energy_beta_p
END_PROVIDER


BEGIN_PROVIDER [ double precision, elf_value_p ]
 implicit none
 BEGIN_DOC  
! ELF value at point
 END_DOC
 double precision, parameter :: Cf = 2.871
 double precision :: D0, Ds
 double precision :: modulus_a2, modulus_b2

 integer :: l

 modulus_a2 = 0.d0
 modulus_b2 = 0.d0
 do l=1,3
   modulus_a2 = modulus_a2 + density_alpha_grad_p(l)**2
   modulus_b2 = modulus_b2 + density_beta_grad_p(l)**2
 enddo

 Ds = kinetic_energy_p - 0.125d0 *  &
  ( (modulus_a2/density_alpha_value_p) + &
    (modulus_b2/density_beta_value_p) )
 D0 = (2.d0*Cf*density_value_p**(5./3.))
 elf_value_p = 1.d0/ (1.d0 + (Ds/D0)**2)

END_PROVIDER