The teacher solved an exercise in class which required you to prove that, if $\Omega$ is a bounded domain in $\mathbb R^n$ and $G$ its Green function, then $G$ is symmetric, i.e. $G(x,y)=G(y,x)$ for all $x,y\in\Omega$.
He starts by arguing that, by Green's representation formula, we can obtain:
$$G(x,y)=\int_{\partial(\Omega\smallsetminus B_\epsilon(x))} G(x,Y)\partial_\nu G(y,Y)d\sigma(Y),$$
or, if we set $G_x(y):=G(x,y)$, that:
$$G_x(y)=\int_{\partial(\Omega\smallsetminus B_\epsilon(x))} G_x(Y)\partial_\nu G_y(Y)d\sigma(Y).$$
Trouble is, this formula doesn't convince me. The formula supposedly used here is that:
$$u(y)=\int_{\partial\Omega\smallsetminus B_\epsilon(x))} u(Y)\partial_\nu G_y(Y)d\sigma(Y),$$
for any $u$ harmonic in $\Omega\smallsetminus B_\epsilon(x)$. But in my opinion, this requires to use the Green function not of $\Omega$, which is $G$, but of $\Omega\smallsetminus B_\epsilon(x)$, which is a different domain and thus should have, if it has one, a different Green function, say $G^\epsilon$. So I would see that formula as:
$$G_x(y)=\int_{\partial(\Omega\smallsetminus B_\epsilon(x))}G_x(Y)\partial_\nu G_y(Y)d\sigma(Y).$$
Now I know that $G_x(y)=N_x(y)+h_x(y)$, where $N_x(y)=N(x-y)$ is the Newton potential and $h_x(y)$ is a function harmonic in $\Omega$, continuous on $\overline\Omega$ and having value $-N_x$ at $\partial\Omega$. For $G_\epsilon$, we will then have $h_x^\epsilon$ analogously. If I add and subtract $h_y(Y)$ from $G^\epsilon$ in the above integral, I will get:
$$G_x(y)=\int_{\partial(\Omega\smallsetminus B_\epsilon(x))}G_x(Y)\partial_\nu G_y(Y)d\sigma(Y)+\int_{\partial(\Omega\smallsetminus B_\epsilon(x))}G_x(Y)\partial_\nu(h_y^\epsilon(Y)-h_y(Y))d\sigma(Y).$$
An application of Green's second identity gives me:
\begin{multline*} \int_{\partial(\Omega\smallsetminus B_\epsilon(x))}G_x(Y)\partial_\nu(h_y^\epsilon(Y)-h_y(Y))d\sigma(Y)=\int_{\partial(\Omega\smallsetminus B_\epsilon(x))}\partial_\nu G_x(Y)\cdot(h_y^\epsilon(Y)-h_y(Y))d\sigma(Y)+{} \\ {}+\int_{\Omega\smallsetminus B_\epsilon(x)}[G_x(Y)\Delta(h_y^\epsilon-h_y)(Y)-\Delta G_x(Y)\cdot(h_y^\epsilon(Y)-h_y(Y))]dY. \end{multline*}
Both functions are harmonic on the integration domain, since $x$ is excluded, so the volume integral is zero. The surface integral will split into integral on $\partial\Omega$ minus integral on ball boundary. The integral on $\partial\Omega$ sees $h_y=h_y^\epsilon=-N_y$, so it is zero. So I am left with:
$$\int_{\partial(\Omega\smallsetminus B_\epsilon(x))}G_x(Y)\partial_\nu(h_y^\epsilon(Y)-h_y(Y))d\sigma(Y)=-\int_{\partial B_\epsilon(x)}(h_y^\epsilon(Y)-h_y(Y))\partial_\nu G_x(Y)d\sigma(Y).$$
An analogous manipulation can be done for the integral in the teacher's formula, the $\partial\Omega$ integral vanishing since $G_x(Y)=0$ there, so I'd be left with:
$$G_x(y)=-\int_{\partial B_\epsilon(x)}G_x(Y)\partial_\nu G_y(Y)d\sigma(Y)-\int_{\partial B_\epsilon(x)}(h_y^\epsilon(Y)-h_y(Y))\partial_\nu G_x(Y)d\sigma(Y).$$
I see no reason whatsoever for the second integral to be zero, which would make me agree with what the teacher said. This probably makes proving the symmetry way harder. Am I missing something?
More calcs
Starting from the above formula, I got:
\begin{align*} G_x(y)-G_y(x)={}&-\int_{\partial B_\epsilon(x)}G_x(Y)\partial_\nu G_y(Y)d\sigma(Y)+{} \\ &{}+\int_{\partial B_\epsilon(y)}G_y(Y)\partial_\nu G_x(Y)d\sigma(Y)+{} \\ &{}-\int_{\partial B_\epsilon(x)}(h_y^\epsilon(Y)-h_y(Y))\partial_\nu G_x(Y)d\sigma(Y)+{} \\ &{}+\int_{\partial B_\epsilon(y)}(h_x^\epsilon(Y)-h_x(Y))\partial_\nu G_y(Y)d\sigma(Y). \end{align*}
Substituting $G_x=N_x+h_x$ and $G_y=N_y+h_y$ in the first two terms gives four types of terms:
- The ones with $h_x\partial_\nu h_y$ or $x,y$ swapped, which are continuous (harmonic, in fact) functions integrated over sets of infinitesimal measure (as $\epsilon\to0^+$ the measure of $\partial B_\epsilon(x),\partial B_\epsilon(y)$ is infinitesimal), so they can be estimated with $C\epsilon^k$ where $C$ is the supremum of the integrand in a big solid ball (to let the estimate hold for all the boundaries) and $k$ is an exponent (a positive one) depending on the dimension; so these go to zero as $\epsilon\to0^+$;
- The ones with $N_x\partial\nu N_y$ and $x,y$ swapped, which can be shown to be equal and thus cancel out;
- The ones with $N_x\partial_\nu h_y$ or $x,y$ swapped, where the Newton potential is constant $C\epsilon^{2-n}$ ($n$ being the dimension of the space) on the ball boundary, then we estimate those normal derivatives by the sup norm of $\nabla h_x$ or $\nabla h_y$ depending on the term, and the measure of the ball boundary is $\omega_n\epsilon^{n-1}$ so in the end these are estimated by $K\epsilon$, and thus are $o(\epsilon)$;
- Those with $h_x\partial N_y$ or $x,y$ swapped, which are infinitesimal like the first type, since $h_x\partial N_y$ is integrated on $\partial B_\epsilon(x)$ where it is nonsingular, and analogously for the other term.
So since the equality holds for all $\epsilon$, we can send $\epsilon\to0^+$ and neglect the contribution of those terms.
Let us see what happens to the other two terms as $\epsilon\to0^+$. Those can be split by separating $h_x$ and $h_x^\epsilon$, $h_y$ and $h_y^\epsilon$. $h_x^\epsilon$ is harmonic in $\Omega\smallsetminus B_\epsilon(y)$, continuous on the closure of that, and coincides with $-N_x$ on the boundary, in particular on $\partial B_\epsilon(y)$. So the first of our two terms contains the following integral:
\begin{align*} \int_{\partial B_\epsilon(y)}h_x^\epsilon(Y)\partial_\nu G_y(Y)d\sigma(Y)={}&-\int_{\partial B_\epsilon(y)}N_x^\epsilon(Y)\partial_\nu G_y(Y)d\sigma(Y)={} \\ {}={}&-\int_{\partial B_\epsilon(y)}N_x(Y)\partial_\nu N_y(Y)d\sigma(Y)-\int_{\partial B_\epsilon(y)}N_x(Y)\partial_\nu h_y(Y)d\sigma(Y). \end{align*}
The second integral here is infinitesimal just like we have said twice in the above analysis of the other terms. So the first one is left. But then, when we analyse the term with $h_y^\epsilon$, we will end up with this same term with $x,y$ swapped, which, as said above, cancels this one, and another term like the one we have said to be infinitesimal. So these contributions are infinitesimal. Hence, we are left with:
\begin{align*} G(x,y)-G(y,x)={}&\int_{\partial B_\epsilon(x)}h_y(Y)\partial_\nu G_x(Y)d\sigma(Y)+{} \\ &{}-\int_{\partial B_\epsilon(y)}h_x(Y)\partial_\nu G_y(Y)d\sigma(Y)+o(\epsilon), \end{align*}
$o(\epsilon)$ summing up all the above considerations. Again, we expand the $G$'s and get terms with $N$'s and terms with $h$'s in the normal derivative. The latter terms are estimated as we have said a thousand times now, and hence are chucked into the $o(\epsilon)$. So:
\begin{align*} G(x,y)-G(y,x)={}&\int_{\partial B_\epsilon(x)}h_y(Y)\partial_\nu N_x(Y)d\sigma(Y)+{} \\ &{}-\int_{\partial B_\epsilon(y)}h_x(Y)\partial_\nu N_y(Y)d\sigma(Y)+o(\epsilon). \end{align*}
I mistakenly thought the ones with $G_x$ and $G_y$ which we just split would be $h_y(x)-h_x(y)$ because of the integrand, but the integration surface is not the right one. If we add in the appropriate $h_x^{B_\epsilon}$ and $h_y^{B_\epsilon}$, these will indeed become $h_y(x)-h_x(y)$, and the added terms will then have to be subtracted, but I'm sure you can guess where they will end up: in the $o(\epsilon)$, since those functions are harmonic in the balls and thus can be estimated etc etc. Actually: no. Dammit. Those functions depend on $\epsilon$, and on the boundaries they are the Newton potential, so no way I cannot estimate them. Then again, IIRC, in the proof of Green's representation formula one step is precisely proving these integrals tend to $h_y(x)$ and $h_x(y)$ respectively as $\epsilon\to0^+$. So finally:
$$h_y(x)-h_x(y)=\lim_{\epsilon\to0^+}\dots=G(x,y)-G(y,x)=N_x(y)+h_x(y)-N_y(x)-h_y(x)=h_x(y)-h_y(x),$$
hence:
$$h_x(y)=h_y(x),$$
and we are done. Is there anything wrong in the above? And again, are those extra terms that weren't in the teacher's argument just zero?
Update
I wrote to the professor, and he suggested I use:
$$u(y)=\int_{\partial\Omega}(u(Y)\partial_\nu N_y(Y)-N_y(Y)\partial_\nu u(Y))d\sigma(Y),$$
for $u$ harmonic in $\Omega$. Use that with $G_x$ in $\Omega\smallsetminus B_\epsilon(x)$. The integral is then split into integral over $\partial\Omega$ minus integral over $\partial B_\epsilon(x)$. The $u(Y)$-$\partial\Omega$ part vanishes since $u=G_x$ is zero on that border. So we get:
$$G_x(y)=-\int_{\partial\Omega}N_y(Y)\partial_\nu G_x(Y)d\sigma(Y)-\int_{\partial B_\epsilon(x)}G_x(Y)\partial_\nu N_y(Y)d\sigma(Y)+\int_{\partial B_\epsilon(x)}N_y(Y)\partial_\nu G_x(Y)d\sigma(Y).$$
Let us now recall what we are aiming at:
$$G_x(y)=\int_{\partial(\Omega\smallsetminus B_\epsilon(x))}G_x(Y)\partial_\nu G_y(Y)d\sigma(Y). \tag{TO BE PROVEN}$$
The integral here is also split into two pieces, of which only the $\partial B_\epsilon$ one is left:
$$G_x(y)=-\int_{\partial B_\epsilon(x)}G_x(Y)\partial_\nu G_y(Y)d\sigma(Y). \tag{TO BE PROVEN B}$$
We see we almost have that above: we are missing an $h_y$ term for the middle term above to be our term. So our TO BE PROVEN reduces to:
$$-\int_{\partial\Omega}N_y(Y)\partial_\nu G_x(Y)d\sigma(Y)+\int_{\partial B_\epsilon(x)}G_x(Y)\partial_\nu h_y(Y)d\sigma(Y)+\int_{\partial B_\epsilon(x)}N_y(Y)\partial_\nu G_x(Y)d\sigma(Y)=0.$$
Using $G_x=N_x+h_x$ on the last integral, the $N$-only part reduces to the spherical average of $N_y$ over $\partial B_\epsilon$, and $N_y$ is harmonic on $B_\epsilon$, hence that gives me $N_y(x)$, and I am left with proving:
$$N_y(x)=\int_{\partial\Omega}N_y(Y)\partial_\nu G_x(Y)d\sigma(Y)-\int_{\partial B_\epsilon(x)}G_x(Y)\partial_\nu h_y(Y)d\sigma(Y)-\int_{\partial B_\epsilon(x)}N_y(Y)\partial_\nu h_x(Y)d\sigma(Y).$$
And I'm stuck again. How can I proceed?
Update 2
If I use Green's formula above for $\Omega\smallsetminus B_\epsilon(y)$ and $N_y$, I should get:
$$N_y(x)=\int_{\partial\Omega}(N_y\partial_\nu N_x-N_x\partial_\nu N_y)d\sigma+\int_{\Omega\smallsetminus B_\epsilon(y)}N_x\Delta N_ydy-\int_{\partial B_\epsilon(y)}(N_y\partial_\nu N_x-N_x\partial_\nu N_y)d\sigma.$$
$N_y$ is harmonic where the middle term integrates, hence that term is zero. We have shown in general that the last integral tends to $N_x(y)$, so since that is $N_y(x)$ we cancel it with the LHS and get:
$$\int_{\partial\Omega}(N_y\partial_\nu N_x-N_x\partial_\nu N_y)d\sigma=0,$$
surprising though it may seem. But does that help? Seems not to me…