Hint: Consider forming the constraints first. You have a sphere = 4/3*pi*r^3 volume cylinder volume is less than this.
An assumption that will make things easier is that the sphere is symmetrical everywhere so consider a cylinder with z axis corresponding to z axis of sphere and centre sphere at origin.
Cross section of cylinder will therefore be equal to cross section of sphere at that point. So if cross section at slice in xy plane for some z0 has radius r0 then the cylinder must also have at most a radius of r0.
Now if radius of cylinder is r_cyl then r_cyl < r (radius of sphere) and height can be calculated by considering cross section information of the sphere.
Volume of the cylinder is pi*r_cyl^2*h and h can be calculated using cross section information of the sphere.
So if we pick an r_cyl value we need to relate this to the maximum h.
If we have the equation of a circle as x^2 + y^2 = r_sphere^2 then we need to find y (which is r_cyl) in terms of x (which is the h* which is related to the height of cylinder).
So h* = SQRT(r_sphere^2 - r_cyl^2) and the heigh of the cylinder due to symmetry is 2*(r_sphere - h*).
So now you have constraints for r_cyl, and h* in terms of r_sphere and this is now just a normal optimization problem.