My doubt is regarding the last part. How does one find the differential equation and then prove that they lie on the given surface?
I am aware of the method, by writing dx/u = dy/v (or its equivalent in cylindrical). It's just that I am unable to get any answer.
Reference: Part of this question has been asked here - To find the velocity potential
The velocity field is given by $\mathbf{v} = - \nabla \phi$. Using spherical coordinates $(r,\theta,\varphi)$ where the polar angle $\theta$ is measured with respect to the $x$-axis we have
$$\mathbf{v} = -\left(\frac{\partial \phi}{\partial r}\mathbf{e}_r + \frac{1}{r} \frac{\partial \phi}{\partial \theta}\mathbf{e}_\theta + \frac{1}{r \sin \theta}\frac{\partial \phi}{\partial \varphi}\right)$$
The potential, $\phi = mr^{-1} - Ur \cos \theta$, does not depend on the azimuthal angle $\varphi$, and the flow is axisymmetric with velocity components
$$v_r = - \frac{\partial }{\partial r}\left(\frac{m}{r}- U r \cos \theta \right) = \frac{m}{r^2}+ U \cos \theta, \\v_\theta = -\frac{1}{r} \frac{\partial }{\partial \theta}\left(\frac{m}{r}- U r \cos \theta \right) = -U \sin \theta$$
For three-dimensional axisymmetric flow, the velocity components are related to the streamfunction, $\psi$, by
$$v_r = \frac{1}{r^2\sin \theta} \frac{\partial \psi}{\partial \theta}, \quad v_\theta = -\frac{1}{r\sin \theta} \frac{\partial \psi}{\partial r}$$
Hence,
$$\tag{1}\frac{\partial \psi}{\partial \theta} = r^2 \sin \theta \,v_r = m \sin \theta + U r^2 \cos \theta \sin \theta,$$ $$\tag{2}\frac{\partial \psi}{\partial r} = -r \sin \theta \,v_\theta = Ur \sin^2 \theta$$
Integrating (1) (with respect to $\theta$) and (2) (with respect to $r$) we get we get $$\psi = -m \cos \theta + \frac{1}{2} U r^2 \sin^2 \theta + F(r), \\ \psi = \frac{1}{2} U r^2 \sin^2 \theta + G(\theta)$$
Comparing we see that $F(r) = 0$ and $G(\theta ) = - m \cos \theta$ and the stream function (up to an arbitrary constant) is
$$\psi = -m \cos \theta + \frac{1}{2} U r^2 \sin^2 \theta,$$
and streamlines lie on surfaces where $\psi = \text{constant}$, or equivalently, $ U r^2 \sin^2 \theta -2m \cos \theta= \text{constant}$.