[tex]\vec f(x,y,z)=x^4\,\vec\imath-x^3z^2\,\vec\jmath+4xy^2z\,\vec k\implies\nabla\cdot f=4x^3+4xy^2=4x(x^2+y^2)[/tex]
By the divergence theorem,
[tex]\displaystyle\iint_S\vec f\cdot\mathrm d\vec S=4\iiint_Rx(x^2+y^2)\,\mathrm dV[/tex]
where [tex]R[/tex] is the region with boundary [tex]S[/tex]. Convert to cylindrical coordinates, taking
[tex]x=u\cos v[/tex]
[tex]y=u\sin v[/tex]
[tex]z=z[/tex]
so that
[tex]\mathrm dV=u\,\mathrm du\,\mathrm dv\,\mathrm dz[/tex]
Then the integral is
[tex]\displaystyle4\iiint_Rx(x^2+y^2)\,\mathrm dV=4\int_0^{2\pi}\int_0^3\int_0^{u\cos v+3}u\cos v((u\cos v)^2+(u\sin v)^2)u\,\mathrm dz\,\mathrm du\,\mathrm dv[/tex]
[tex]\displaystyle=4\int_0^{2\pi}\int_0^3\int_0^{u\cos v+3}u^4\cos v\,\mathrm dz\,\mathrm du\,\mathrm dv[/tex]
[tex]=486\pi[/tex]
