See my answer at this question. A somehow more abstract explanation can also be found here.
I took a look at K. Huang's book Statistical Mechanics (chap. 17) to be sure about this: yes, we put it by hand, and we assume that the external field $h$ is weak so that only $h$ enters and no higher powers.
- This is the most cumbersome point. I have always been quite surprised by that term, to be honest, and I suspect that Landau originally put it by hand using is incredible physical intuition.
In Huang's book, I found no justification for it, but one
can be found here:
We expect that the Gizburg Landau functional should be of
$$f = f_0 + \alpha_2 M^2 + \alpha_4 M^4 + ... + \kappa | \nabla M|^2 + ...$$
The first line is familiar from (spatially independent) Landau theory in zero external field that we have previously discussed. The second term gives the functional some spatial dependence.
Let us try to justify this form. As above, we want to expand the free energy in powers of $M$. In addition here we want to expand the free energy in increasing powers of derivatives. The justification for this is that we expect spatial changes to be rather long length scales, so increasing number of derivatives will be increasingly small. In more detail our functional needs to be a scalar, so no terms linear in $\nabla M$ can enter
. We also expect the ground state should be the case where
$M$ is spatially uniform — and this is true so long as $\kappa$ is positive.
One way to understand this $|\nabla M|^2$ term is to see that it is
disfavoring any domain walls in the system. If $M$ is pointing in some
particular direction, then it should want to point the same direction
nearby. (Perhaps this is even more obviously correct in the case of
Heisenberg magnets where we can imagine $M$ rotating just a little bit
from place to place).
Another justification can be found here, where the author basically says that we include the gradient term to include the effect of fluctuations.
In conclusion, the gradient term takes into account the local spatial variations of the order parameter, and since we assume that fluctuations have a long length scale (we are close to the critical point) we can neglect higher order derivatives. Anyway, if we were to consider other gradient terms, they would all be even powers, because $M$ must be a scalar (see for example here).
As for the coefficient $1/2$, Huang says that it is chosen this way "to fix the scale of $m(x)$" (the order parameter), and I think he means that we arbitrarily fix it and then redefine the order parameter accordingly.