Starting with Green's formula $\iint_{D} \frac{\partial Q}{\partial x}-\frac{\partial P}{\partial y} dxdy=\oint_{bD}Pdx+Qdy$, take $Q=u_x$, $P=-u_y$, and $D$ a small circle of radious r centered at $(x_0,y_0)$. We get:
$$ \iint_{D} \Delta u \space dA = \oint_{\partial D} \nabla u \cdot (dy,-dx) .$$
If we parametrized the boundary of $D$ as:
$$ \begin{equation}
x(\theta) = x_0 + r\cos(\theta) \\
y(\theta) = y_0 + r\sin(\theta) \\
\end{equation} $$
then
$$ (dy,-dx) = r (\cos(\theta), \sin(\theta)) d\theta = r \nu d\theta $$
where $\nu$ is the exterior normal to $\partial D$, so the Green's formula you mentioned implies the formula in the problem suggestion.
Since $u$ is harmonic, the left side is $0$ and we get:
$$ 0 = \oint_{\partial D} \nabla u \cdot (dy,-dx)=\\
\int_{0}^{2\pi} \{u_x(x_0 +r\cos(\theta), y_0 + r\sin(\theta)) r\cos(\theta) + \\
\space u_y(x_0 +r\cos(\theta), y_0 + r\sin(\theta)) r\sin(\theta)\} \space d\theta $$
Now we compare what we have and what we want to get. Let me point some obvious facts (hind side is 20/20...) :
* the last formula involves derivatives of $u$, but the problem has only $u$, without derivatives, which suggest we will need to differentiate something
* the right side of the formula involves $r$, but the left side has no $r$, so its derivative with respect to $r$ is $0$
You can give it the final kick.
$$\dots\dots\dots$$
Now we observed that
$$ 0 = r \space \int_{0}^{2\pi} \{u_x(x_0 +r\cos(\theta), y_0 + r\sin(\theta)) \cos(\theta) + \\
\space u_y(x_0 +r\cos(\theta), y_0 + r\sin(\theta)) \sin(\theta)\} \space d\theta =\\
= r \space \frac{\partial}{\partial r} \int_{0}^{2\pi} u(x_0 +r\cos(\theta), y_0 + r\sin(\theta)) \space d\theta .$$
Therefor the integral over the boundary of the little circle is constant, and the value can be found by taking limit as $r \rightarrow 0^+$