I'm trying to calculate the Laplacian $\triangle\Phi_W^{IJ}$ of the following "generalized Whitehead potential":
$$ \Phi_W^{IJ}(\vec{x}) = \int\rho^I(\vec{x}')\rho^J(\vec{x}'')F(\vec{x},\vec{x}',\vec{x}'')\,d^3x'\,d^3x''\,, $$
where
$$ F(\vec{x},\vec{x}',\vec{x}'') = \frac{\vec{x} - \vec{x}'}{|\vec{x} - \vec{x}'|^3} \cdot \left(\frac{\vec{x}' - \vec{x}''}{|\vec{x} - \vec{x}''|} - \frac{\vec{x} - \vec{x}''}{|\vec{x}' - \vec{x}''|}\right) $$
and $\rho^I$, $\rho^J$ are smooth functions with compact support. The difficulty that arises here is that I cannot simply take the pointwise Laplacian of $F$, since it does not agree with the distributional Laplacian. Apparently I am missing terms of the form
$$ c_1\int\rho^I(\vec{x}')\rho^J(\vec{x}'')\frac{\delta(\vec{x} - \vec{x}')}{|\vec{x} - \vec{x}''|}\,d^3x'\,d^3x''\,, \quad c_2\int\rho^I(\vec{x}')\rho^J(\vec{x}'')\frac{\delta(\vec{x} - \vec{x}'')}{|\vec{x} - \vec{x}'|}\,d^3x'\,d^3x''\,, $$
and I need to determine the constants $c_1$ and $c_2$. I know how to do this in the "simple" case of a 1/r potential (i.e., solving the Poisson equation), but here it seems to be more complicated.
I'd be thankful for any suggestions.
EDIT: I made some progress and got a result, but it seems to be wrong. So here's what I did:
First of all, since F depends only on coordinate differences such as $\vec{x} - \vec{x}'$, we have
$$ (\nabla + \nabla' + \nabla'')F = 0 $$
and thus
$$ \triangle F = (\nabla' + \nabla'')^2F\,. $$
We then introduce new coordinates via
$$ \vec{x}^{\pm} = \frac{1}{2}(\vec{x}' \pm \vec{x}'')\,, \quad \nabla^{\pm} = \nabla' \pm \nabla'' $$
and express $F$ in the new coordinates. The integral then reads
$$ \Phi_W^{IJ}(\vec{x}) = 8\int_{\mathbb{R}^3}\int_{\mathbb{R}^3}\rho^I(\vec{x}^+ + \vec{x}^-)\rho^J(\vec{x}^+ - \vec{x}^-)\triangle^+F\,d^3x^+\,d^3x^-\,, $$
where I explicitly wrote down the domains of integration, since these will become important in the next step. Now the inner integral over $\vec{x}^+$ contains a Laplacian $\triangle^+$ with respect to $\vec{x}^+$, which is much nicer. We can thus determine at which values of $\vec{x}^+$ the function $F$ has singularities, cut out a ball around them and split the integration domain into these balls plus the rest. For the balls we use the surface theorem, the rest we can simply integrate as usual. Let $B_{\epsilon}'$ and $B_{\epsilon}''$ be balls of radius $\epsilon$ around the singularities $\vec{x}^+ = \vec{x} + \vec{x}^-$ and $\vec{x}^+ = \vec{x} - \vec{x}^-$. Partial integration (and omitting terms that vanish in the limit $\epsilon \to 0$) then yields
$$ 8\lim_{\epsilon \to 0}\int_{\mathbb{R}^3}\int_{\mathbb{R}^3 \setminus (B_{\epsilon}' \cup B_{\epsilon}'')}\rho^I(\vec{x}^+ + \vec{x}^-)\rho^J(\vec{x}^+ - \vec{x}^-)\triangle^+F\,d^3x^+\,d^3x^- = 2\partial_{\alpha}\partial_{\beta}U^I\partial_{\alpha}\partial_{\beta}\chi^J - 2\triangle(U^IU^J) - 16\pi\rho^IU^J - 8\pi\rho^JU^I\,, $$
$$ 8\lim_{\epsilon \to 0}\int_{\mathbb{R}^3}\int_{B_{\epsilon}'}\rho^I(\vec{x}^+ + \vec{x}^-)\rho^J(\vec{x}^+ - \vec{x}^-)\triangle^+F\,d^3x^+\,d^3x^- = 8\lim_{\epsilon \to 0}\int_{\mathbb{R}^3}\int_{\partial B_{\epsilon}'}\rho^I(\vec{x}^+ + \vec{x}^-)\rho^J(\vec{x}^+ - \vec{x}^-)\nabla^+F \cdot \vec{n}\,d^2x^+\,d^3x^- = \frac{16\pi}{3}\rho^IU^J\,, $$
$$ 8\lim_{\epsilon \to 0}\int_{\mathbb{R}^3}\int_{B_{\epsilon}''}\rho^I(\vec{x}^+ + \vec{x}^-)\rho^J(\vec{x}^+ - \vec{x}^-)\triangle^+F\,d^3x^+\,d^3x^- = 8\lim_{\epsilon \to 0}\int_{\mathbb{R}^3}\int_{\partial B_{\epsilon}''}\rho^I(\vec{x}^+ + \vec{x}^-)\rho^J(\vec{x}^+ - \vec{x}^-)\nabla^+F \cdot \vec{n}\,d^2x^+\,d^3x^- = 4\pi\rho^JU^I\,, $$
where $\vec{n}$ is the unit normal to the surface of the ball and the potentials $U^I$ and $\chi^I$ are
$$ U^I(\vec{x}) = \int\frac{\rho^I(\vec{x}')}{|\vec{x} - \vec{x}'|}\,d^3x'\,, $$
$$ \chi^I(\vec{x}) = -\int\rho^I(\vec{x}')|\vec{x} - \vec{x}'|\,d^3x'\,, $$
so that in total I get
$$ \triangle\Phi_W^{IJ} = 2\partial_{\alpha}\partial_{\beta}U^I\partial_{\alpha}\partial_{\beta}\chi^J - 2\triangle(U^IU^J) - \frac{32\pi}{3}\rho^IU^J - 4\pi\rho^JU^I\,. $$
This has more or less the expected form. Unfortunately it does not reproduce the literature result
$$ \triangle\Phi_W = 2\partial_{\alpha}\partial_{\beta}U\partial_{\alpha}\partial_{\beta}\chi - 2\triangle(U^2) - 12\pi\rho U $$
for the simpler Whitehead potential $\rho^I = \rho^J = \rho$ (Will, Theory and Experiment in Gravitational Physics, Cambridge, 1993, eq. (4.36)), so I guess my result must be wrong...
Any help would be highly appreciated.
EDIT²: Following an "experimental" approach and simply inserting some test functions $\rho^I$, $\rho^J$ showed that the correct result should be:
$$ \triangle\Phi_W^{IJ} = 2\partial_{\alpha}\partial_{\beta}U^I\partial_{\alpha}\partial_{\beta}\chi^J - 2\triangle(U^IU^J) - 8\pi\rho^IU^J - 4\pi\rho^JU^I\,, $$
which is consistent with the literature result. Of course this is not a rigorous proof of this formula, so I am not really satisfied, but at least I have something I can use in my calculations... I will leave this question open until I found a rigorous proof.
