Your equation is $Du = 0$ where
$$
D = \left(y^2\partial_{x}^2 - 2xy\partial_{x}\partial_{y} + x^2\partial_y^2 - \frac{y^2}x\partial_x - \frac{x^2}y\partial_y\right).
$$
With some work, we can factor $D$:
\begin{align*}
D & = (y\partial_x - x\partial_y)^2 + \left(\frac yx - \frac xy\right)x\partial_y + \left(\frac xy - \frac yx\right)y\partial_x \\
& = (y\partial_x - x\partial_y)^2 + \left(\frac xy - \frac yx\right)(y\partial_x - x\partial_y) \\
& = \left(y\partial_x - x\partial_y + \frac xy - \frac yx\right)(y\partial_x - x\partial_y).
\end{align*}
This suggests solving the PDE in two steps. The first one is to solve for $v$ from
$$
\left(y\partial_x - x\partial_y + \frac xy - \frac yx\right)v = 0.
$$
The second one is to solve for $u$ from
$$
\left(y\partial_x - x\partial_y\right)u = v.
$$
In both steps, the characteristics are defined by $\frac{\partial y}{\partial s} = -x$ and $\frac{\partial x}{\partial s} = y$ where $s$ is the curve parameter (or in your notation, $\frac{dy}{dx} = -\frac{x}y$, and the curve parameter is $x$). The equation for $v$ is $\frac{\partial v}{\partial s} = -\left(\frac xy - \frac yx\right)v$ while the equation for $u$ is $\frac{\partial u}{\partial s} = v$.
To solve for characteristics, differentiate equations for $y$ and $x$ once more to get $\frac{\partial^2 y}{\partial s^2} = -y$ and $\frac{\partial^2 x}{\partial s^2} = -x$. The solution is $x = c_1(r)\cos(s) + c_2(r)\sin(s)$, $y = -c_1(r)\sin(s) + c_2(r)\cos(s)$, where $r$ is the index of the characteristic curve. (Graphically, characteristics are circles, as you have seen from $x^2 + y^2 = \psi$. Your $\psi$ is the same as $c_1(r)^2 + c_2(r)^2$ in my notation.)
With two arbitrary functions of $r$, you do have some freedom here. You get to choose where each curve starts, which is the point on the circle that $s = 0$ corresponds to. Usually, you use the initial condition to decide, because making $s = 0$ correspond to the initial solution curve will simplify a lot of calculation. Anyway, since the initial condition is not given, I'll go ahead and do some computation.
Your lecturer's choice makes the positive $y$-axis correspond to $s = 0$ (so $c_1(r) = 0$), while your choice makes the positive $x$-axis correspond to $s = 0$ (so $c_2(r) = 0$). Your lecturer's choice leads to something like
\begin{align*}
x(r, s) & = r\sin s \\
y(r, s) & = r\cos s
\end{align*}
while your choice leads to something like
\begin{align*}
x(r, s) & = r\cos s \\
y(r, s) & = -r\sin s.
\end{align*}
It should not matter which one you choose. Note that with your lecturer's choice, you have
\begin{align*}
\frac yx - \frac xy & = \cot s - \tan s \\
& = 2\cot(2s)
\end{align*}
and with your choice, you have
\begin{align*}
\frac yx - \frac xy & = -\tan s + \cot s \\
& = 2\cot(2s).
\end{align*}
In both cases, $\frac yx - \frac xy$ becomes the same $2\cot(2s)$, and that is the expression you'll need to solve for $v$ and $u$:
\begin{align*}
\frac{\partial v}{\partial s} & = \left(\frac yx - \frac xy\right)v \\
& = 2\cot(2s)v\\
\therefore v & = C(r)\sin(2s) \\
\text{and from} \frac{\partial u}{\partial s} = v \quad \Longrightarrow \quad
u & = -\frac 12C(r)\cos(2s) + D(r) \\
\therefore u & = \tilde C(x^2 + y^2) \cdot (x^2 - y^2) + \tilde D(x^2 + y^2).
\end{align*}
Although $\cos(2s)$ becomes different expressions in different parametrizations, the final expression is the same because the sign difference is subsumed in $\tilde C$. I did get sloppy and ignored the case $\sin(2s) < 0$. However, the final solution works out if we keep in mind that $\tilde C$ and $\tilde D$, once solved from the initial condition, are valid only in the quadrant that contains the initial solution curve. It is natural to expect that the solution does not emanate past the $x$-axis or the $y$-axis because we have $x$ and $y$ in the denominators in the original PDE.
