# Applying the constant number of electrons in DFT density optimisation

(Previous related question: Finding mathematically the ground state density in DFT)

I am studying the density optimisation procedure (in particular for Orbital-Free DFT) this thesis.

The derivative of the energy is given as:

$$\frac{\delta E[n(r)]}{\delta n(r)} = v(r) + v_i(r) + \frac{1}{2}(3\pi^2)^{2/3}\left(n(r)\right)^{2/3},$$

where $$v(r)$$ is the external potential and $$v_i(r)$$ is the Hartree potential, or electron-electron repulsion.

The chemical potential can thus be obtained after integration: $$\mu = \frac{1}{V} \int_{\Omega} \frac{\delta E[n(r)]}{\delta n(r)}dr.$$

What confuses me is the next equation, which allows to find the density of the next step $$k+1$$: $$n(r)^{(k+1)} = n(r)^{(k)} - t\left(\frac{\delta E[n(r)]}{\delta n(r)} - \mu\right),$$ where $$t$$ is a step size to help convergence.

If I understand this correctly, $$\mu$$ is an average variation of energy. The density at a given $$r$$ is thus varied by the difference between the derivative at that position and the average derivative.

$$v(r)$$ is a constant for a given position, and $$v_i(r)$$ depends only on the density, just like the kinetic term above (Thomas-Fermi in this case). It thus seems like the variation of density will never quite get to zero except if the density is exactly the right value to balance the external potential. This seems simplistic and wouldn't explain orbitals beyond the $$s$$ orbital.

In parallel, I am implementing a simple code for the hydrogen atom, but I can't get anything to converge nor give reasonable values; I must be doing something wrong somewhere.