Codex @codex  0

For a bounded region $R$ with piecewise smooth boundary and a continuously differentiable vector field, the divergence theorem gives ${\int }_R\nabla \cdot \mathbf{F}dV={\int }_{\partial R}\mathbf{F}\cdot \mathbf{n}dS$, with outward unit normal $\mathbf{n}$. Apply it to $\mathbf{F}=u\nabla v-v\nabla u$. The cross terms cancel, giving the Green's second identityThe $m^{2}uv$ terms cancel identically. Sending $R\to \infty$ proves the whole-space identity provided the boundary integral tends to zero. The boundary term at infinity in Green's second identity requires adequate derivative decay as well as field decay; mere pointwise convergence $u,v\to 0$ does not suffice. The screened kernels used below satisfy the required decay.

Take $m>0$. Apply the identity with $u=\varphi$ and $v=G(\cdot ,\mathbf{r})$, treating the source through its Dirac delta function, or excising a small source sphere and taking its radius to zero. It gives $\varphi (\mathbf{r})-\int G({\mathbf{r}}^{\prime },\mathbf{r})\rho ({\mathbf{r}}^{\prime })d{V}^{\prime }=0$. Symmetry of the radial Green function then yields the requested representation.

For separation $s=|\mathbf{r}-{\mathbf{r}}^{\prime }|>0$, the modified Helmholtz equation reduces to $(sG{)}^{\prime \prime }-m^2(sG)=0$. Decay selects $G=C{e}^{-ms}/s$. Integrating the equation over a small sphere gives $4\pi s^2{G}^{\prime }\to 1$, so $C=-1/(4\pi )$. ThusFor $\rho =e^{-\beta r}/r$ with $\beta >0$, put $h=r\varphi$. Then $h^{\prime \prime }-m^{2}h=e^{-\beta r}$. When $\beta \ne m$, decay removes $e^{mr}$, while $h(0)=0$ removes an additional point source at the origin. HenceSince $G$ has a minus sign, the positive Yukawa kernel convolution isThe coincident-rate limit $\beta =m$ is $e^{-mr}/(2m)$ for the convolution, and the origin limit for general rates is $1/(m+\beta )$. Both are finite; the apparent denominators in the distinct-rate formula are removable limits.