What kind of free energy do we use for a superconductor in a magnetic field? My reasoning is as follows (using Gaussian units):
Start from the second law:
$$dU=TdS+dW,$$
where $dW$ is the work done by the magnetic field. To derive $dW$, we consider a solenoid with current $I$ and voltage $V$:
$$dW=I\cdot V\cdot dt,$$
with
$$\begin{align}
V&=\frac{N}{c}\cdot\frac{d(B\cdot A)}{dt} \qquad &&\text{(Faraday's law)} \\
I&=\frac{c}{4\pi}\cdot\frac{H}{n}. \qquad &&\text{(Ampere's law)}
\end{align}$$
Now plug in $I$ and $V$ in $dW$. We get
$$dW=\frac{1}{4\pi}\cdot HdB.$$
Since $B=H+4 \pi M$, if we exclude the energy of the magnetic field itself, we find
$$dW=HdM.$$
Thus the second law is
$$dU=TdS+HdM.$$
For the Helmholtz free energy we have
$$F=U-TS \Rightarrow dF=-SdT+HdM.$$
Similarly, for the Gibbs free energy:
$$G=F-HM \Rightarrow dG=-SdT-MdH.$$
In experiment, we control $T$ and $H$, therefore at any given $T$ and $H$, the system should minimize its Gibbs free energy $G$. Ergo the energy we talk about for superconductor in magnetic field should be Gibbs free energy. 
The questions is, why the famous superconducting textbooks (like Tinkham, Schmidt and de Gennes) use Helmholtz free energy $F$ instead of $G$? It doesn't make sense to me to minimize $F$ instead of $G$ (e.g. in deriving the GL equations using variation method).
 A: Ah, this is a classic and rather tragic issue in most statistical-mechanics books. I won't talk about superconductivity, but the general case of magnetic systems. The problem is thermodynamics as it is, was formulated for fluid systems for which pressure and volume were observables and could be manipulated. In that case, say I calculated the partition function of a given microscopic Hamiltonian of the system, then
$$
Z=\mathrm{Tr}\ e^{-\beta\mathcal{H}}
$$
and $F=-k_BT\ln Z$ is the Helmholtz free energy of the system as $(T,V)$ are kept constant (It is an (N,V,T) ensemble if we include the chemical potential $\mu$). Anyway, if we impose upon the system an external force, aka, pressure, then the ensemble changes to a $(p,T)$ one and we instead have the Gibbs potential to work with. So,
$$
Z[p]=\mathrm{Tr}\ e^{-\beta\left[\mathcal{H}+\int pV\right]}
$$
and $G(p,T)=-k_BT\ln Z[p]$. Now the analogue of pressure and volume in magnetic systems is total net magnetization ($M$) and the external field ($H$). They form a conjugate pair of generalized extensive variable and generalized intensive force. So, the correspondence is $M\leftrightarrow V$ and $H\leftrightarrow -p$.
$$
Z[H]=\mathrm{Tr}\ e^{-\beta\left[\mathcal{H}-\int MH\right]}
$$
Hence, in the presence of an external field, we actually have the Gibbs free energy $G(H,T)$ instead of the Helmholtz free energy $F(M,T)$, which will only be used in the null field case.
Most books gloss over this and continue to call the potential the Helmholtz free energy in all cases for magnetic systems and almost always get away with it. It doesn't change most of the calculations but it can be confusing conceptually, as  in your case.
A: You are actually correct. It should be Gibbs Free Energy. In some of the cases minimising Helmholtz free energy also works as it gives the same equations for order parameter as one minimizes Gibbs free energy. But in discussing the case, for example, of parallel critical fields of thin films, Gibbs free energy should be used.
In cases where minimization of Helmholtz free energy works, if you minimize the Gibbs free energy, you will find that when you minimize $G$ with respect to vector potential $A$, the difference between $G$ and $F$, $H_0.B = H_0. (\nabla \times A)$ comes out to be $0$. So it gives same equations as when $F$ is minimized.
Source: Superfluidity and Superconductivity by David R. Tilley and John Tilley, Section 6.3 and 6.6
