What is a mean field? Consider an interacting electron gas in a box. The Hamiltonian will have an interaction term
$$H = \sum_{i,j}u\:c_j^{\dagger}c_jc_i^{\dagger}c_i$$
$u$ is somehow dependent on length such that only particles that are close interact.
Doing a mean field approximation we can get a quadratic Hamiltonian instead.
$$H = \sum_{i,j}u\left< c_j^{\dagger}c_j\right>c_i^{\dagger}c_i + u\left< c_i^{\dagger}c_i\right>c_j^{\dagger}c_j - u\left< c_i^{\dagger}c_i\right>\left< c_j^{\dagger}c_j\right>$$
In order to solve the problem I will minimize the free energy of the problem with respect to the mean field parameters and use those values in my Hamiltonian.
However, I am now confused because to me it seems like there are two ways to think of these mean fields, and I don't know which one are correct.


*

*I can think of them as the mean density of the particles. Meaning wherever there are lots of particles these mean fields will also be large.

*I can put choose their value however I want in order to get the lowest energy state. I would then put all my particles at the boundary of the box, and then all the mean field in the middle, which to me seems like the lowest energy state.
These two interpretations give different results, and I don't know how to think about it. So which one is correct, or are there another way to think of it that I haven't considered?
Edit: I am still not sure I understand whether minimization of free energy is equivalent of solving the self consistent equation in my situation because I have a constraint coming from the wall of the box.
I mean a minimal free energy does not imply $\frac{\partial F}{\partial \left<{a^\dagger}a\right>} = 0$ anymore, and I don't understand how we then get the self consistent equation anymore.
 A: The mean field approximation doesn't refer to some external field, but to the assumption that we can replace interactions between the particles with some effective average behavior. Under this assumption, when we have an interaction of the sort you gave, we treat the average densities as numbers, and then we indeed have
$$H_I^{MF} = u\sum_{i,j}(n_j c^{\dagger}_i c_i + n_i c^{\dagger}_j c_j - n_i n_j )$$
with $n_j$ the average density at point $j$ (I assume that $i$ and $j$ are position dependent). (By the way - this is only the interaction, I assume you also have a free part to the Hamiltonian?)
Now what you should do is solve the set of equations self-consistently. These $n_j$ are not parameters that you can choose freely! They are the true densities of the particles. So basically you have a free effective Hamiltonian, with some unknown parameters, but this parameters can be deduced from the Hamiltonian itself, by requiring a minimization of the free energy.
The mean field approximation basically assumes that fluctuations and correlations are weak, and we can use the averages to get a good approximation to the behavior on the microscopic level.
