I have recently started doing DFT calculations as part of my PhD. The general approach to DFT in my department appears to be that it is regarded as a tool for doing calculations, which need not be understood in detail. While I understand this to some degree, I however still feel the need to at least understand the basic theoretical ideas behind a DFT calculation. I am therefore using this post to summarize a few of my questions about DFT starting from the basic theory.
As I understand it, the formalism of DFT rests upon the Hohenberg-Kohn theorems, stating that the ground state properties of a many-body system are unique functionals of the ground state density, which furthermore minimizes the total energy functional $E[n(\mathbf{r})]$ of the system:
$$ E[n(\mathbf{r})] = T[n(\mathbf{r})] + E_H[n(\mathbf{r})] + E_{ext}[n(\mathbf{r})] + E_{xc}[n(\mathbf{r})] $$
Here $T[n(\mathbf{r})]$ is the kinetic energy functional, $E_H[n]$ the Hartree functional and $E_{ext}[n(\mathbf{r})]$ the energy due to external potentials while $E_{xc}[n]$ is the unknown exchange-correlation energy functional, which encompasses the exchange and correlation effects of the many-body problem and must be suitably approximated by e.g. the LDA.
In practice DFT calculations are based on the Kohn-Sham approach, which makes the ansatz that there is a non-interacting system, which reproduces the density of the many-body system one wishes to describe. Consequently, instead of solving a many-body Hamiltonian, one solves a set of single-particle equations, the Kohn-Sham equations:
$$ (\nabla^2 + V^{KS})\phi_i(\mathbf{r}) = \epsilon_{i}\phi_{i}(\mathbf{r}) $$
with the Kohn-Sham potential given by:
$$ V^{KS}(\mathbf{r}) = V_{ext}(\mathbf{r}) + V_{H}(\mathbf{r}) + \frac{\delta E_{xc}}{\delta n(\mathbf{r})} $$
The density is then readily evaluated as a sum over the occupied orbitals:
$$ \rho(\mathbf{r}) = \sum_{occ}\lvert \phi_{i}(\mathbf{r}) \rvert^2 $$
My main question is essentially, if anyone can walk me through how the above theoretical ideas are used in a concrete DFT calculation. Say, e.g., I want to calculate the ground state energy of some molecule. How does a DFT calculation proceed using the equations outlined above?
My own (partial) answer: The first step is specifying the energy functional, which relies on choosing a suitable approximation for the exchange-correlation functional. My guess is then that the next step is to solve the Kohn-Sham single-particle equations. To do one needs the Kohn-Sham potential, which must be found from the functional derivative of the total energy functional. I am honestly not sure if this is how a DFT calculator proceeds, but it seems like the only way to get the Kohn-Sham potential. Having solved the KS system, one can then calculate the electron density, which should then equal the density of the interacting system. However, I speculate that the procedure may involve some sort of self-consistency such that one should in practice make an iterative procedure to find the correct interacting density.
Assuming everything is correct so far, how does one then get the ground state energy? So far we have only calculated the ground state density and the energies of the Kohn-Sham orbitals are not physical energies of the system.
