Solution (source code)

= Solution

Use the diagram's labels: $a,b$ are horizontal side lengths and $c$ is the height of the <rectangular parallelepiped>. Its four top vertices are at horizontal distance $\sqrt{a^2+b^2}/2$ from the axis. Thus the contact <constraint> and <volume> are
$$
c+\frac{H(a^2+b^2)}{4R^2}=H,\qquad V_p=abc.
$$
Apply a <Lagrange multiplier>:
$$
\mathcal L=abc-\lambda\left(c+\frac{H(a^2+b^2)}{4R^2}-H\right).
$$
For positive dimensions, the <derivative> equations are
$$
bc=\lambda\frac{Ha}{2R^2},\qquad
ac=\lambda\frac{Hb}{2R^2},\qquad ab=\lambda.
$$
Substituting the last equation into the first two gives
$$
c=\frac{Ha^2}{2R^2}=\frac{Hb^2}{2R^2}.
$$
Hence $a=b$, and the <constraint> becomes $2c=H$. Therefore
$$
\boxed{a=b=R,\qquad c=\frac H2,\qquad V_{p,\max}=\frac{R^2H}{2}}.
$$
These are the <maximal cylinder and box in an elliptic paraboloid> dimensions for the rectangular case. A horizontal rotation of the square footprint gives the same <volume> by rotational symmetry.