[Q-e-developers] recompute the total energy after scf is finished

[Lorenzo Paulatto]

Dear all,
I'm trying to recompute the total energy after the scf cycle is finished 
recycling the wavefunctions that are stored to file at the end.

The worflow is like this:
1. call weights to recompute bands occupations, and eventually demet 
(smearing contribution)
2. call sum_bands to recompute rho from the wfcs+occupations, this also 
computes eband
3. call v_of_rho to recompute the potential, this also computes ehart, 
etxc, etxcc, etc.
4. call deband = delta_e() = \int v_(r) rho_r(r) dr (function from 
electrons)
5. recompute ewald term
6. etot = eband + deband + (etxc - etxcc) + ewld + ehart + demet

All the terms are identical to the SCF calculation, except for deband 
which is slightly different.

I've noticed that the value deband computed during the scf cycle uses 
the value of v(r) from the second-last iteration. I.e. if the 
calculations converges in N iterations:
  deband_N = \int v_{N-1}(r) \rho_N(r) dr

I would expect the two values to differ by something of the order
  descf = \int (\rho_{N-1)-rho_N(r)) v_{N-1}(r) dr
but it turns out that the difference is more of order sqrt(descf).

 From this I deduce that the form I use for the total energy may not be 
the correct variational one. Finally, I'm confused. Anybody has the time 
to help me find out what I may be doing wrong?

kind regards

P.S. I can provide the modified files (electrons.f90 + one new file) to 
anybody on request, I'm not posting them here for size concerns.

-- 
Dr. Lorenzo Paulatto
IdR @ IMPMC -- CNRS & Université Paris 6
+33 (0)1 44 275 084 / skype: paulatz
http://www-int.impmc.upmc.fr/~paulatto/
23-24/4é16 Boîte courrier 115, 4 place Jussieu 75252 Paris Cédex 05