One of the implicit assumptions that the question is using is that the mass of the solid is constant, and that the density of the solid is constant. If the mass isn't constant, you could set it to zero. If the density wasn't constant, you could indeed crush the mass into a thin axis of infinite density, but that wouldn't be characteristic of many solids!
Constant mass means
$$\int_0^h \rho \pi r^2 \, \mathrm{d}z = M = \text{const.}$$
or, equivalently, as a constraint equation
$$G(r) = \int_0^h \rho \pi r^2 - \frac{M}{h} \, \mathrm{d}z = 0$$
(Constant density is easier to enforce: simply let the variable $\rho$ be constant)
Therefore, the problem is now:
Minimise $I(r) = \int_0^h \frac{1}{2}\rho \pi r^4 \, \mathrm{d}z$ such that $G(r) = \int_0^h \rho \pi r^2 - \frac{M}{h} \, \mathrm{d}z = 0$
That is, you now have a constrained variational problem, which requires the use of a Lagrangian multiplier $\lambda$. So the cost function is instead
$$J(r,\lambda) = I(r) + \lambda G(r) = \int_0^h \frac{1}{2} \rho \pi r^4 + \lambda \left(\rho \pi r^2 - \frac{M}{h}\right) \, \mathrm{d}z = \int_0^h F(r,\lambda) \, \mathrm{d}z$$
The Lagrangian multiplier $\lambda$ now behaves as an addition variable to minimise over, so the solution is obtained by evaluating
$$\frac{\mathrm{d}}{\mathrm{d}z}\left(\frac{\partial F}{\partial r'}\right) - \frac{\partial F}{\partial r} = 0$$
$$\frac{\mathrm{d}}{\mathrm{d}z}\left(\frac{\partial F}{\partial \lambda'}\right) - \frac{\partial F}{\partial \lambda} = 0$$
where the prime notation denotes differentiation with respect to $z$. (Since no $r'$ or $\lambda'$ terms are present, the first terms of the above equations are equal to zero).
Solving this will then yield the cylinder as the optimal solid.
