I always get confused with this...
Let $S$ be the surface $$z=x^2+y^2, z\leq 1,$$ oriented so that the normal vector has positive $z$-coordinate, and let $F$ be the vector field $(yz,-xz+\sin(z), e^{x^2+y^2})$. Compute the surface integral $$\int_SF.n\, dS \,\,$$
Of course this would be much easier using Gauss divergence theorem, but I wanted to calculate it directly, because it is what the OP asked
Today I feel generous:
So let's go:
Let's use the following coordinates $\Phi:(\theta,z)\in\Omega:=[0,2\pi]\times[0,1]\mapsto\begin{pmatrix}z\cos(\theta) \\ z\sin(\theta) \\ z^2\end{pmatrix}$
We first compute one of the possible normal vectors $$\tilde{\vec{n}}=\frac{\Phi_\theta\times \Phi_z}{\|\Phi_\theta\times \Phi_z\|}=\frac{\begin{pmatrix}-z\sin(\theta) \\ z\cos(\theta) \\ 0\end{pmatrix}\times \begin{pmatrix}\cos(\theta) \\ \sin(\theta) \\ 2z\end{pmatrix}}{\|\Phi_\theta\times \Phi_z\|}= \frac{\begin{pmatrix}2z^2\cos(\theta) \\ 2z^2\sin(\theta) \\ -z\end{pmatrix}}{\|\Phi_\theta\times \Phi_z\| }.$$ As the normal vector has to be chosen with positive $z$ coordinate, we thus set $$\vec{n}=-\tilde{\vec{n}}$$
\begin{align*}\int_SF\cdot dA& =\int_\Omega F\circ\Phi(\theta,z)\cdot \vec{n}\|\Phi_\theta\times \Phi_z\|dzd\theta = \int_\Omega F\circ\Phi(\theta,z)\cdot-(\Phi_\theta\times \Phi_z) \ dzd\theta\\ & =\int_0^{2\pi}\int_0^1 \begin{pmatrix}\sin(\theta) z^3\\ -\cos(\theta)z^3+\sin(z^2) \\ e^{z^2} \end{pmatrix}\cdot \begin{pmatrix}-2z^2\cos(\theta) \\ -2z^2\sin(\theta) \\ z\end{pmatrix} dzd\theta \end{align*}
you can take it from here