# [Pw_forum] Electrons kinetic energy

Pietro Bonfa' pietro.bonfa at fis.unipr.it
Thu Oct 17 19:03:00 CEST 2013

Dear Prof. Giannozzi,

I believe that this is exactly what you suggested

DO ik = 1, nk
DO ibnd = 1, nbnd
DO j = 1, npw
IF(gamma_only)THEN !.and.j>1)then
ek = ek +  2*conjg(evc(j,ibnd)) * evc(j,ibnd) * &
g2kin(j) * wg(ibnd,ikk)
ELSE
ek = ek +  conjg(evc(j,ibnd)) * evc(j,ibnd) * &
g2kin(j) * wg(ibnd,ikk)
ENDIF
ENDDO
ENDDO
ENDDO

I also started from the local kinetic energy density just to
double-check the code since I'm not familiar with the QE structure and
variables.

I encountered problems obtaining the total kinetic energy by summing
over the real grid the kinetic energy density.
I found that the problem is solved if I change the from 'Smooth' to
'Dense' the fft grid in the elf.f90 code (line 121) and I was wondering
if I did something wrong or not.

Anyway, it was almost only a curiosity. I'm happy with the code reported
above.

Thanks again,
Pietro Bonfa'

On 10/17/13 18:24, Paolo Giannozzi wrote:
> What do you want exactly to calculate? The "average electron
> kinetic energy" for norm-conserving pseudopotentials is a
> rather simple object. In atomic (Ry) units:
>   \sum_k \sum_i f_{i,k} <\psi_{i,k}|\nabla^2|\psi_{i,k}> =
>   \sum_k \sum_i f_{i,k} \sum_G |\psi_i(k+G)|^2 (k+G)^2
> where f_{i,k} is the occupancy of Kohn-Sham state \psi_{i,k}
> whose plane-wave components are \psi_i(k+G). YOu do not need
> the "local kinetic energy density":
>   \sum_k \sum_i f_{i,k} |\nabla psi_{k,i}(r)|**2
> P.
>
> On Wed, 2013-10-16 at 17:04 +0200, Pietro Bonfa' wrote:
>> Dear All,
>>
>> following the recipes given in this thread:
>>
>> http://qe-forge.org/pipermail/pw_forum/2009-May/087347.html
>>
>> I added two routines to the pp.x code to calculate the average
>> electronic kinetic energy. I started from both the elf.f90 subroutine in
>> PP code and from the pw2casino implementation.
>>
>> It's almost only a matter of copy/paste but there is a thing I cannot
>> understand.
>>
>> Let "kkin" be the variable with the kinetic energy density, following
>> the elf.f90 code I have:
>>
>>      0    ! symmetrize the local kinetic energy if needed
>>      1	  IF ( .not. gamma_only) THEN
>>      2	     !
>>      3	     CALL sym_rho_init ( gamma_only )
>>      4	     !
>>      5	     aux(:) =  cmplx ( kkin (:), 0.0_dp, kind=dp)
>>      6	     CALL fwfft ('Dense', aux, dfftp)
>>      7	     ALLOCATE (aux2(ngm))
>>      8	     aux2(:) = aux(nl(:))
>>      9	     !
>>     10	     ! aux2 contains the local kinetic energy in G-space to be
>> symmetrized
>>     11	     !
>>     12	     CALL sym_rho ( 1, aux2 )
>>     13	     !
>>     14	     aux(:) = (0.0_dp, 0.0_dp)
>>     15	     aux(nl(:)) = aux2(:)
>>     16	     DEALLOCATE (aux2)
>>     17
>>     18	     CALL invfft ('Dense', aux, dfftp)
>>     19	     kkin (:) = dble(aux(:))
>>     20	     !
>>     21	  ENDIF
>>     22	  T = sum( kkin(:) )*omega / ( dfftp%nr1*dfftp%nr2*dfftp%nr3 )
>>
>> To get the correct kinetic energy (T), with respect to elf.f90, I had to
>> change
>> CALL fwfft ('Smooth', aux, dfftp)
>> to
>> CALL fwfft ('Dense', aux, dfftp)
>> on line 6 because otherwise I couldn't perform the sum correctly (even
>> if using dffts%nr1 nr2 nr3).
>>
>> With this trick the elf-derived and the pw2casino-derived routines give
>> similar results in all the cases I used for testing.
>> Nonetheless, I'm completely unfamiliar with the QE implementation, and
>> probably there is no need to change Smooth to Dense.
>> Moreover, I don't really know if performing sym_rho on the Dense grid is
>> fine or not.
>>
>>
>> Thanks,
>> Pietro
>>
>> ps: the subroutines are attached.
>>
>>
>> _______________________________________________
>> Pw_forum mailing list
>> Pw_forum at pwscf.org
>> http://pwscf.org/mailman/listinfo/pw_forum
>

--
Pietro Bonfa' - PhD student
Dipartimento di Fisica e Scienze della Terra "Macedonio Melloni"
Viale delle Scienze 7A
43124 Parma - Italy