The important point to be taken here is that we do not need to work in the coordinates we have been given from the exercise.
Observe that if $\log p(x,y)$ is a solution then $\log p(x-X, y-Y)$ is also a solution for arbitrary constants $X, Y$. This corresponds to moving the origin. But note also that transformation that rotates the $xy$ plane around the origin also creates a new solution. In other words, we can assume that the polynomial is in the form $p(x,y) = Ax^2 + By^2 + C,$ linear terms being destroyed by moving the origin appropriately and the crossterm by a suitable rotation. Finding the solutions is now simple since $ 0 = u_{xx} + u_{yy} = {p(p_{xx} + p_{yy}) - p_x^2 - p_y^2 \over p^2}$ and so $ (Ax^2 + By^2 + C)(2A + 2B) = 4A^2x^2 + 4B^2y^2.$ If $A = -B$ then the LHS vanishes and $A = B = 0$ giving us a constant solution with $C$ arbitrary. Otherwise we must have $A = B$. Therefore the only nonconstant solutions in these coordinates are $u = \log(A(x^2 + y^2)) = \log(x^2 + y^2) + \log A$. Moreover, since this solution is invariant to rotations (being only a function of the radial distance $r^2 = x^2 + y^2$) the full set of solutions is recovered from this by translations. Therefore, the only nonconstant solutions must be of the form $u = \log ((x-X)^2 + (y-Y)^2) + Z$ with $X, Y, Z$ arbitrary. You can check that all of these are indeed solutions.