Why is dipole moment equal to derivative of the energy? Maybe the following is obvious, but I haven't found an answer after thinking about it for some while, so I give it a try here.
Given some physical system in an external homogeneous electric field $\mathbf E$, we write the energy of the system as a function $E(\mathbf E)$ (I guess we consider the system to be in equilibrium such that no hysteresis occurs).
The electrostatic interaction energy between the charge distribution $\rho(\mathbf r)$ of the system and the external field is given by $E_\text{el}(\mathbf E) = - \mathbf \mu(\mathbf E) \cdot \mathbf E$, where $\mathbf \mu$ is the dipole moment. This energy is obviously not equal to $E(\mathbf E)$ if the field does work on the system, e.g. by polarizing it.
My question can (I think) alternatively be asked in one of the following ways:


*

*Why can the dipole moment be written as $\mu(\mathbf E) = - \frac{dE}{d\mathbf E}$ (note I'm writing $E$ and not $E_\text{el}$).

*Why is the work done on the system quadratic in the external field to lowest nonvanishing order?

*Why can I compute the total energy $E(\mathbf E)$ (including work done on the system) by integrating the dipole moment from zero to finite field? In other words: Why is the energy needed to polarize the system encoded in the dipole moment?
I know that an answer in the field of quantum mechanics is provided by the Hellmann-Feynman theorem if one substitutes dipole moment by expectation value of the dipole moment, but I'm looking for a more general argument.
 A: Let's think about the energy in terms of the potential $V$ instead of the field $\newcommand{\E}{\mathbf{E}}\E$. Then energy is given by $\int \rho V$, where $\rho$ is the charge density, which will adjust itself as the potential field $V$ is changed. The charge density $\rho$ adjusts itself to be in a minimum of energy. Thus is a change $\delta \rho$ is applied, the energy of the change will be zero: $\int V\delta \rho =0$. 
Now let's compute the change in energy that occurs when the potential field $V$ is changed. But first remember that when $V$ is changed, $\rho$ will react to adjust itself to a minimum of the new $V$. Thus we get $E+ \delta E=\int (\rho + \delta \rho)(V+\delta V) = \int \rho V + \int V\delta \rho  + \int \rho \delta V$. And now $\rho V$ is just the original energy $E$ and we have already argued that $\delta \rho V$ is zero. Thus we get that $\delta E = \int \rho \delta V$. Of course, this is exactly what $\delta E$ would have been if the charge density we fixed. But in the case of a fixed charge density, we know that $\dfrac{\partial{E}}{\partial{\E}}$ is $\boldsymbol{\mu}$. Therefore, we just also have that $\dfrac{\partial{E}}{\partial{\E}}=\boldsymbol{\mu}$ even when the charge density is allowed to relax in response to the field.
With this is place, it is easy to see that $E(\E)=\int_\mathbf{0}^\E \dfrac{\partial{E}}{\partial{\E}} d\E = \int_\mathbf{0}^\E \boldsymbol{\mu}(\E) d\E$.
This answers your first and third questions. The answer to the second question has two parts. First we assume the energy has a power series in $\E$ simply because the underlying laws of electrostatics are smooth, so we don't expect any non-analytic behavior. However, the energy may still be linear in $\E$ to lowest order. From what we have shown, this means the electric dipole would have to be zero with no applied field. This is almost always not the case, since a non-zero dipole produces a field which costs energy. Therefore, if there are no other energetics to consider, the ground state is the one with zero dipole moment, and when there is no applied field, $\dfrac{\partial{E}}{\partial{\E}}=\boldsymbol{\mu}=\mathbf{0}$
