Solved it. The grand potential is $$ A =kT\int dE\,D(E)\log\left(1-\frac{1}{\exp\left(\frac{E-\mu}{kT}\right)+1}\right) =-kT\int dE\,D(E)\log\left(1+\exp\left(\frac{\mu-E}{kT}\right)\right)$$
from where we calculate $$ S=-\left.\frac{\partial A}{\partial T}\right)_{\mu,V}=-\frac{A}{T}-\frac{1}{T}\int dE\,D(E)\frac{1}{\exp\left(\frac{E-\mu}{kT}\right)+1}(\mu-E)$$$$ S=-\left.\frac{\partial A}{\partial T}\right|_{\mu,V}=-\frac{A}{T}-\frac{1}{T}\int dE\,D(E)\frac{(\mu-E)}{\exp\left(\frac{E-\mu}{kT}\right)+1}$$
and $$ p=-\left.\frac{\partial A}{\partial V}\right)_{T,V}=-\frac{A}{V}$$$$ p=-\left.\frac{\partial A}{\partial V}\right|_{T,V}=-\frac{A}{V}$$
and it is then straightforward to calculate $$G=U-TS+pV =\mu\int dE\,\frac{D(E)}{\exp\left(\frac{E-\mu}{kT}\right)+1} = \mu N$$
CQFD.
