Solution (source code)

= Solution

For $k=1$, the <vector field> is $(x^2y,xy^2,z)^T$. On the cone,
$$
\mathbf F\cdot d\mathbf S_+
=\left(\frac{x^3y+xy^3}{\rho}+1-\rho\right)dx\,dy
=(xy\rho+1-\rho)dx\,dy.
$$
The first term has zero <integral> by symmetry. Using <polar coordinates> for the remaining <double integral> gives
$$
\Phi_+=2\pi\int_0^1(1-\rho)\rho\,d\rho=\frac\pi3.
$$
On the paraboloid $z=\rho^2-1$,
$$
\mathbf F\cdot d\mathbf S_-
=\bigl(2x^3y+2xy^3-(\rho^2-1)\bigr)dx\,dy
=(2xy\rho^2+1-\rho^2)dx\,dy.
$$
Again the $xy$ term cancels by symmetry, giving
$$
\Phi_-=2\pi\int_0^1(1-\rho^2)\rho\,d\rho=\frac\pi2.
$$
Therefore direct <surface integral> yields
$$
\boxed{\Phi=\frac\pi3+\frac\pi2=\frac{5\pi}{6}}.
$$
No separate disk at $z=0$ is included: that rim joins the two surfaces and is not an additional face of the closed boundary.