This is not an answer, as I don't know how to evaluate the integral in terms
of elementary functions. However,
just in case someone else does...
I am assuming that $r \le R$ and that $h$ is 'large' so that it completely
surrounds the bored hole. I am also assuming that the big cylinder has its centerline along the $y$ axis and that the drilled hole is along the $x$ axis.
If we can compute the volume of the material of the original cylinder that was removed by drilling then we can compute the volume of the cylinder with a hole.
The volume of removed material is given by the volume of
$D_{r,R} = \{ (x,y,z) | z^2+y^2 \le r^2, z^2+x^2 \le R^2 \}$. This can be written as
$D_{r,R} = \{ (x,y,z) | z \in [-r,r], |x| \le \sqrt{R^2-z^2}, |y| \le \sqrt{r^2-z^2} \}$
so we can write
$\operatorname{vol}D_{r,R} = \int_{z=-r}^r \int_{y=-{\sqrt{r^2-z^2}}}^{\sqrt{r^2-z^2}} \int_{x=-{\sqrt{R^2-z^2}}}^{\sqrt{R^2-z^2}} 1 \ dx dy dz = 8 \int_0^r \sqrt{r^2-z^2} \sqrt{R^2-z^2} dz$.
As a sanity check, if we fix $r$, then
$\lim_{R \to \infty} {D_{r,R} \over 2 R} = \pi r^2$.
However, I don't think that this has an elementary solution.
A numerical solution is straightforward, if that suffices.
Then the volume of the drilled cylinder is
$\pi R^2 h-D_{r,R}$.
Position the cube so that a vertex is at the origin and the cube lies in the first octant. The long diagonal has length $a\sqrt3$, so if the cylinder’s height is $h$, its near cap lies at a distance of $\frac12\left(a\sqrt3-h\right)$ from the origin. The normal to this cap’s plane is $(1,1,1)$, so an equation of this plane is $$x+y+z=\frac{\sqrt3}2\left(a\sqrt3-h\right).$$ The expression on the right-hand side is the axis-intercept of this plane with all three coordinate axes. The end cap touches the $x$-$y$ plane at the midpoint of the $x$- and $y$-intercepts, so using the Pythagorean theorem, we can obtain the square of its radius, namely $$\frac38\left(a\sqrt3-h\right)^2-\frac14(a\sqrt3-h)^2 = \frac18\left(a\sqrt3-h\right)^2.$$ The rest of the task, as you say, is a trivial calculus exercise.
Best Answer
Let $S$ be our shape. The picture below (made with Pov-Ray) illustrates $S$. It is the intersection of a torus (red) with a cylinder (blue, half-transparent).
Furthermore, let $h$ be its height and $r$ be the radius at the top (so the radius at the middle is $r+h$) (see image below, made with Paint :$).
Now we call the radius at a given height $x$ the function $R(x)$ (where $x$ goes from $x=-\frac12h$, bottom, to $x=\frac12h$, top. See the graph above. $x$ is on the x-axis). Luckily the shape of a circle is easy, so we can compute $R(x)$ by $$R(x)=r+\sqrt{\tfrac14h^2-x^2}$$ So now we can calculate the volume by rotating $R(x)$ around the x-axis (see this wikipedia article for more information about revolving an area between a curve and the x-axis):
\begin{align}\pi\int_{-h/2}^{h/2}R(x)^2 &=\pi\int_{-h/2}^{h/2}\left(r+\sqrt{\tfrac14h^2-x^2}\right)^2dx\\ &=\pi\int_{-h/2}^{h/2}\left(r^2+2r\sqrt{\tfrac14h^2-x^2}+\sqrt{\tfrac14h^2-x^2}^2\right)dx\\ &=\pi\left(\int_{-h/2}^{h/2}r^2dx+2r\int_{-h/2}^{h/2}\sqrt{\tfrac14h^2-x^2}dx+\int_{-h/2}^{h/2}|\tfrac14h^2-x^2|dx\right)\\ &=\pi\left(hr^2+2r(\tfrac12\pi\cdot \tfrac14h^2)+\int_{-h/2}^{h/2}\tfrac14h^2dx-\int_{-h/2}^{h/2}x^2dx\right)\\ &=\pi\left(hr^2+\tfrac14r\pi h^2+\tfrac14h^3-[\tfrac13x^3]_{-h/2}^{h/2}\right)\\ &=\pi\left(hr^2+\tfrac14r\pi h^2+\tfrac14h^3-\tfrac1{12}h^3\right)\\ &=h\pi\left(r^2+\tfrac14r\pi h+\tfrac16h^2\right)\\ \end{align}