Greens function of Laplace operator I am having trouble deriving the Greens function for the three dimensional Laplacian $\nabla^2$. I wish to solve the equation
$$ -\nabla_r^2 G(r,r') = \delta(r-r') \quad\quad\quad r,r'\in\mathbb{R}^3$$
I know that $\frac{1}{4\pi|r-r'|}$ is the solution, but I would like to derive this for myself. I believe that the correct approach is to apply a Fourier transform to go from $r$ space to $k$ space:
$$ -\int dr \; e^{-ikr}\nabla_r^2 G(r,r') = \int dr\; e^{-ikr}\delta(r-r') $$
$$ k^2 \tilde{G}(k,r') = e^{-ikr'} $$
Where $\tilde{G}(k,r') \equiv \int dr \; e^{-ikr} G(r,r')$ is the Fourier transform of $G$. We can solve for $\tilde{G}$ easily:
$$ \tilde{G}(k,r') = \frac{e^{-ikr'}}{k^2} $$ 
And now we transform back into $r$ space:
$$ G(r,r') = \int dk\; e^{+ikr}\tilde{G}(k,r') = \int dk\;\frac{e^{ik(r-r')}}{k^2}$$
But unfortunately this integral does not converge! We can see this by integrating by parts:
$$ G(r,r') = -\int dk\;e^{ik(r-r')}\frac{d}{dk}\frac{1}{k} $$
$$ G(r,r') = -\frac{e^{ik(r-r')}}{k} + i(r-r')\int dk\;\frac{e^{ik(r-r')}}{k}$$ 
$$ G(r,r') = -\frac{e^{ik(r-r')}}{k} - (r-r')^2 Ei\left(ik(r-r')\right)$$ 
Neither the boundary term nor the error function converge when we evaluate for $k\in\mathbb{R}^3$, and it is not at all obvious that they should cancel in any meaningful way. I don't see how to arrive at the correct solution $G(r,r')=\frac{1}{4\pi|r-r'|}$
 A: Let us define the Green's function by the equation,
$$\nabla^2 G(\mathbf{r},\mathbf{r}_0) = \delta(\mathbf{r}-\mathbf{r}_0).$$
Now let us define $S_\epsilon = \{\mathbf{r} : |\mathbf{r}-\mathbf{r}_0| \leq \epsilon \}$, from which we thus have that,
$$\int_{S_\epsilon} \nabla^2 G(\mathbf{r},\mathbf{r}_0) \, d\mathbf{r} = \int_{S_\epsilon} \delta(\mathbf{r}-\mathbf{r}_0) \, d\mathbf{r} = 1$$
by properties of the delta function. By Green's theorem (or analogously, Stokes') we have that,
$$\int_{\partial S_\epsilon} \nabla G(\mathbf{r},\mathbf{r}_0) \cdot \mathbf{n} ds = 1.$$ In two dimensions, $S_\epsilon$ is a circle, and we find,
$$\lim_{\epsilon \to 0} \int_0^{2\pi} \frac{\partial G}{\partial r} r \bigg\rvert_{r=\epsilon} d\theta = 1$$
from which one can deduce the leading order behaviour of the Green's function in two dimensions,
$$G \sim \frac{1}{2\pi} \ln r = \frac{1}{2\pi} \ln |\mathbf{r}-\mathbf{r}_0|.$$
Repeating this analysis for the case of a sphere instead, in three dimensions, you will find that,
$$G \sim -\frac{1}{4\pi r} = -\frac{1}{4\pi |\mathbf{r}-\mathbf{r}_0|}.$$
Note these are the behaviours of the Green's function near $\mathbf{r}_0$. However, if there are no boundary conditions and we consider all of $\mathbb R^3$, then this corresponds to the free space Green's function, as  we only require it to have the correct behaviour near the source.
If you choose the Laplace operator as $-\nabla^2$ so that it is elliptic, then this accounts for the difference in sign, compared to your cited result.
A: I prefer to start with the Green's function for a massive field (defined by: $\left[-\nabla^2 + m^2\right] G = \delta(\mathbf{x}-\mathbf{x}')$) and take the limit as $m\rightarrow 0$. That Green's function is:
$$G(r=|\mathbf{x}-\mathbf{x}'|) = \frac{1}{(2\pi)^{d/2}} \left(\frac{m}{r}\right)^{d/2-1} K_{d/2-1}(mr),$$ where $K_\nu(x)$ is the modified Bessel function of the second kind, $d$ is the number of dimensions of the space. Applying the limiting forms of the modified Bessel functions to take the limit gives:
$$\begin{align}
\lim_{m\rightarrow 0} G(r) &= \frac{1}{(2\pi)^{d/2}} \lim_{m\rightarrow 0}\left(\frac{m}{r}\right)^{d/2-1} \frac{\Gamma(d/2-1)}{2} \left(\frac{2}{mr}\right)^{d/2-1}\\
& = \frac{\Gamma\left(\frac{d}{2} - 1\right)}{4 \pi^{d/2} r^{d-2} }
\end{align}$$
when $d > 2$. For $d=2$ we get:
$$\lim_{m\rightarrow0} G(r) = \frac{1}{2\pi}\lim_{m\rightarrow 0} -\ln(mr).$$ The limit in this case diverges because the massive case requires that $\lim_{r\rightarrow \infty} G(r) = 0$, and that's not required for the massless case. In the massless case we can add a constant with respect to $r$, though, so the answer is:
$$G_2(r) = -\frac{\ln r}{2\pi}.$$
The 1-d case has a similar divergence as the 2-d case, but the pattern fits with:
$$G_1(x) = -\frac{|x|}{2}.$$
The pattern is:
$$-\frac{\partial G(r)}{\partial r} = \frac{1}{S_{d-1}(r)},$$ where $S_n(r)$ is the surface area of the $n$-sphere. The $1$-sphere is the circle, and the $0$ sphere is two isolated points (the measure of isolated points is a count of them).
