For inviscid, rotational flow a weaker form of the Bernoulli equation $p + \frac{1}{2} \rho | \mathbf{u}|^2 + \rho g z = C(\psi)$ is valid, where $C(\psi)$ is constant along a streamline. You can work with this or solve for the pressure field directly from the governing equations for inviscid flow, which reveals the Bernoulli relationship as well.
Using cylindrical coordinates we have $x = r \cos \theta, \, y = r \sin \theta,$ and the basis vectors
$$\mathbf{e}_r = \cos \theta \, \mathbf{e}_x + \sin \theta \, \mathbf{e}_y, \\\mathbf{e}_\theta = -\sin \theta \, \mathbf{e}_x + \cos \theta \, \mathbf{e}_y $$
The velocity field is, for $r^2 = x^2 + y^2 \leqslant a^2 $,
$$u_r = 0, \,\, u_\theta = - \frac{\Gamma r}{a^2}, \,\, u_z = 0,$$
and for $x^2 + y^2 > a^2$,
$$u_r = 0, \,\, u_\theta = - \frac{\Gamma }{r}, \,\, u_z = 0$$
The components of the Euler equation for steady, inviscid flow, $\rho\mathbf{u} \cdot \nabla \mathbf{u} \,\, = - \nabla p + \rho \mathbf{g}$, reduce to
$$\frac{\partial p}{\partial r} = \frac{\rho u_\theta^2}{r} = \begin{cases}\frac{\rho \Gamma^2r}{a^4}, \,\, x^2+y^2 \leqslant a^2 \\ \frac{\rho \Gamma^2}{r^3}, \,\,\,\, x^2+y^2 > a^2\end{cases}$$
and
$$\frac{\partial p}{\partial z}= -\rho g$$
Integrating, we get for $x^2 + y^2 \leqslant a^2$,
$$p = \alpha - \rho g z + \frac{\rho \Gamma^2 r^2}{2a^4} = \alpha - \rho g z + \frac{\rho \Gamma^2 (x^2 + y^2)}{2a^4},$$
and, for $x^2 + y^2 > a^2$,
$$p = \beta - \rho g z - \frac{\rho \Gamma^2}{2r^2} = \beta - \rho g z - \frac{\rho \Gamma^2}{2(x^2 + y^2)} $$
The integration constant $\alpha$ can be evaluated in terms of $\beta \,$ by matching the inner pressure field with the outer pressure field at $x^2 + y^2 = a$, where we have
$$\alpha - \rho g z + \frac{\rho \Gamma^2 (a^2)}{2a^4} = \beta - \rho g z - \frac{\rho \Gamma^2}{2a^2} $$
Solving for $\alpha$ we get
$$\alpha = \beta - \frac{\rho \Gamma^2}{a^2}, $$
and the inner pressure field is
$$p = \beta - \frac{\rho \Gamma^2}{a^2} - \rho g z + \frac{\rho \Gamma^2 (x^2 + y^2)}{2a^4}$$
The free-surface where $p = C$ is given by
$$z = \frac{\beta-C}{\rho g} - \frac{\Gamma^2}{g a^2} + \frac{ \Gamma^2 (x^2 + y^2)}{2ga^4}$$
Best Answer
The field is always directed along the first coordinate axis. Therefore, its streamlines are lines parallel to that axis. Very easy to sketch. If you wish to connect this with the slope-from-the-ratio method, consider that the ratio $0/(x_1x_2)$ is zero. So, the slope of streamlines is everywhere zero.
Sure, the points with $x_1x_2=0$ are exceptional: they are stationary points of the field. You may want to mark them somehow to indicate this fact.