TDFreezeOrbitals and Hartree energy
If we freeze orbitals, the Hartree energy is currently wrong. The code implements
E_H = \frac{1}{2}\int d^3 r n_{\rm tot} (r) v_H[n_{\rm tot}](r)
while the formula should be
E_H = \frac{1}{2}\int d^3 r n_{\rm KS} (r) v_H[n_{\rm KS}](r) + \int d^3 r n_{\rm KS} (r) v_H[n_{\rm frozen}](r)
where n_{\rm tot} = n_{\rm KS} + n_{\rm frozen}
, as the frozen Hartree term acts as an external potential.
Note the fact 1/2
difference for the last term. This comes from the fact that we use the sum of eigenvalues as a proxy to the kinetic energy.
However, as in the total energy we have -\int d^3 r n_{\rm KS} (r) v_{\rm ext}(r)
, I am suggesting to put in the code Hartree energy
E_H = \frac{1}{2}\int d^3 r n_{\rm KS} (r) v_H[n_{\rm tot}](r) - \frac{1}{2}\int d^3 r n_{\rm KS} (r) v_H[n_{\rm frozen}](r)
Alternatively, we could simply do
E_H = \frac{1}{2}\int d^3 r n_{\rm KS} (r) v_H[n_{\rm KS}](r)\,,
i.e. not add any frozen density in the Hartree energy.
After a kick (taken from the test real_time/13.xxx) the energy strongly oscillates with the current main
while it should remain constant.
This can be tested by setting XCFunctional=0
, which shows that the problem comes from the Hartree energy.
A quick fix implemented the equation above fixes the fluctuating energy.