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

[stefano de gironcoli]

ciao Lorenzo,
   I think
      deband = delta_e() = - \int v_(r) rho_r(r) dr
   should be computed with the potential that generated the wfcs and not 
with the one obtained from the new density
   In this way i exactly cancel the corresponding term in eband

stefano

On 11/04/2013 12:13 PM, Lorenzo Paulatto wrote:
> 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.
>