Aloha :)
Als Parameterdarstellung für \(x\) und \(y\) würde ich die Standard-Form der Polarkoordinaten wählen, d.h. \(x\) mit dem Cosinus und \(y\) mit dem Sinus:$$\vec r=\begin{pmatrix}x\\y\\z\end{pmatrix}=\begin{pmatrix}x\\y\\x^2-y^2\end{pmatrix}=\begin{pmatrix}r\cos\varphi\\ r\sin\varphi\\r^2(\cos^2\varphi-\sin^2\varphi)\end{pmatrix}=\begin{pmatrix}r\cos\varphi\\ r\sin\varphi\\r^2\cos(2\varphi)\end{pmatrix}$$wobei \(r\in[0;1]\) und \(\varphi\in[0;2\pi]\) liegen.
Mit Hilfe des totalen Differentials wird klar, wie sich \(\vec r\) in Abhängigkeit von \(r\) und \(\varphi\) ändert:$$d\vec r=\frac{\partial \vec r}{\partial r}\,dr+\frac{\partial\vec r}{\partial\varphi}\,d\varphi$$Das von diesen beiden infinitesimalen Vektoren aufgespannte Flächenelement ist:
$$df=\left\|\left(\frac{\partial \vec r}{\partial r}\,dr\right)\times\left(\frac{\partial\vec r}{\partial\varphi}\,d\varphi\right)\right\|=\left\|\frac{\partial \vec r}{\partial r}\times\frac{\partial\vec r}{\partial\varphi}\right\|\,dr\,d\varphi$$$$\phantom{df}=\left\|\begin{pmatrix}\cos\varphi\\\sin\varphi\\2r\cos(2\varphi)\end{pmatrix}\times\begin{pmatrix}-r\sin\varphi\\r\cos\varphi\\-2r^2\sin(2\varphi)\end{pmatrix}\right\|\,dr\,d\varphi$$$$\phantom{df}=\left\|\begin{pmatrix}-2r^2\sin(2\varphi)\sin\varphi-2r^2\cos(2\varphi)\cos\varphi\\-2r^2\cos(2\varphi)\sin\varphi+2r^2\sin(2\varphi)\cos\varphi\\r\cos^2\varphi+r\sin^2\varphi\end{pmatrix}\right\|\,dr\,d\varphi$$
Das sieht schlimm aus, allerdings helfen uns die Additionstheoreme hier gut weiter:$$\sin(2\varphi)\sin\varphi+\cos(2\varphi)\cos\varphi=\cos(2\varphi-\varphi)=\cos\varphi$$$$\sin(2\varphi)\cos\varphi-\cos(2\varphi)\sin\varphi=\sin(2\varphi-\varphi)=\sin\varphi$$sodass sich der Betrag des Vektors in völliges Wohlgefallen auflöst:$$df=\left\|\begin{pmatrix}-2r^2\cos\varphi\\2r^2\sin\varphi\\r\end{pmatrix}\right\|\,dr\,d\varphi=\sqrt{4r^4+r^2}\,dr\,d\varphi$$
Damit formulieren wir das gesuchte Oberflächenintegral:$$I=\int\limits_{r=0}^1\;\int\limits_{\varphi=0}^{2\pi}\sqrt{4r^4+r^2}\,dr\,d\varphi=\int\limits_{\varphi=0}^{2\pi}d\varphi\cdot\int\limits_{r=0}^1\sqrt{4r^4+r^2}\,dr=2\pi\int\limits_0^1r\sqrt{4r^2+1}\,dr$$$$\phantom{I}=\pi\int\limits_0^1\sqrt{4r^2+1}\,2r\,dr=\pi\int\limits_0^1\sqrt{4r^2+1}\,d(r^2)=\pi\left[\frac{2}{3}\cdot\frac{1}{4}\left(4r^2+1\right)^{3/2}\right]_0^1$$$$\phantom{I}=\frac{\pi}{6}\left(5^{3/2}-1\right)=\frac{\pi}{6}\left(5\sqrt5-1\right)\approx5,3304$$