4
$\begingroup$

I am trying to use Mathematica to check an integral from a book, and cannot get it done. The book talks about a spherically symmetric mass distribution in $R^3$ of the form

$\rho(r)=r^{-\gamma}, \gamma<3$

and it goes on by stating that the Newtonian potential to that distribution is given by

$\phi(r) = 4π\ln r + \mathrm{const.}, \gamma=2$

and

$\phi(r) = \frac{4π}{(2-\gamma)(3-\gamma)}r^{2-\gamma}+\mathrm{const.}, \gamma\neq2$ .

(There are a couple more constants in the book which I set to 1 here for simplicity.)

I wanted to verify this by directly evaluating the integral for the potential

$\phi(\mathbf x) = -\int_{R^3}\frac{\rho(\mathbf x')}{|\mathbf x-\mathbf x'|}\,d^3x' =-4π\int_0^\infty\frac{\rho(r')}{|r-r'|}r'^2\,dr'$ ,

and I am failing.

Here is my code to check the $\gamma=2$ case:

Assuming[{r > 0, g == 2}, -4 π Integrate[ s^(-g) / Abs[r - s] s^2, {s, 0, ∞}]]

Mathematica gives me the message: "Integrate: Integral of 1/Abs[r-s] does not converge on {0,[Infinity]}"

My code for the other case is

Assuming[{r > 0, g != 2, g < 3}, -4 π Integrate[ s^(-g) / Abs[r - s] s^2, {s, 0, ∞}]]

which also yields "Integrate: Integral of s^(2-g)/Abs[r-s] does not converge on {0, ∞}."

What am I missing?

Improve this question
$\endgroup$
7
  • 1
    $\begingroup$ One does problems of this sort with the help of the Gauss theorem. The way you try it, is much too difficult, and directly it will not work. $\endgroup$ Commented Mar 3, 2020 at 16:15
  • $\begingroup$ Hi @AlexeiBoulbitch , thanks for the answer! Sure, there might be simpler ways than directly evaluating the integral. I wanted to try this, since later I want to generalise this to more complex distributions and integrate numerically. But I found a mistake now: I cannot simply replace $|\mathbf x-\mathbf x'|$ in the 3D integrand by $|r-r'|$ in the 1D integrand! I will correct this and try again tomorrow and then report back! $\endgroup$ Commented Mar 3, 2020 at 22:42
  • $\begingroup$ FullSimplify[Integrate[s^(2-g)/Sqrt[(r-s)^2],s]] does this help? You can now try Limit and Series with Assuming. $\endgroup$ Commented Mar 4, 2020 at 5:44
  • $\begingroup$ @Britzel If your (more complex) distributions are spherically- or cylindrically symmetric, the Gauss theorem approach is the way to implement. Do you know it? If you have an asymmetric problem you will have a bad time. I am just struggling against a comparable problem, where I need to integrate a Green function over a distribution. [Here][1] you can see the level of mathematical difficulties arising there. Have fun! [1]: mathematica.stackexchange.com/questions/215614/… $\endgroup$ Commented Mar 4, 2020 at 8:58
  • $\begingroup$ @AlexeiBoulbitch For the start they will be spherically symmetric, yes. Perhaps it is a good idea to go with a Gauß theorem approach for now, as you suggest. Eventually though they will become axisymmetric. But not for now. $\endgroup$ Commented Mar 4, 2020 at 13:16

1 Answer 1

Reset to default
6
$\begingroup$

The problem is with the expansion: to integrate $\frac{1}{\lvert \mathbf{x'} - \mathbf{x}\rvert}$ over a volume with a spherically invariant factor (the charge density), the standard (and the most sensible trick) is to apply spherical harmonic decomposition:

$$\frac{1}{\lvert \mathbf{x'} - \mathbf{x}\rvert}=\sum\limits_{l=0}^\infty\frac{r_<^l}{r_>^{l+1}}\left(\frac{4\pi}{2l+1}\right)\sum\limits_{m=-l}^{l}Y_{lm}(\theta,\phi)Y_{lm}^*(\theta',\phi')$$

With this expansion, the integral will be much simpler; in particular, the angular integral will immediately fix radial dependence. In fact, for spherically symmetric overall factor, the angular integration simply reads

$$\begin{align} \int d\Omega'\rho(r')\frac{1}{\lvert \mathbf{x'} - \mathbf{x}\rvert}=&\sum\limits_{l=0}^\infty\frac{r_<^l}{r_>^{l+1}}\left(\frac{4\pi}{2l+1}\right)\sum\limits_{m=-l}^{l}Y_{lm}(\theta,\phi)\rho(r')\int d\Omega'Y_{lm}^*(\theta',\phi')\\ =&\sum\limits_{l=0}^\infty\frac{r_<^l}{r_>^{l+1}}\left(\frac{4\pi}{2l+1}\right)\sum\limits_{m=-l}^{l}Y_{lm}(\theta,\phi)\rho(r')\sqrt{4\pi}\delta_l^0\delta_m^0\\ =&\frac{4\pi}{r_>}\rho(r') \end{align}$$

The radial integral is now straightforward:

$$\begin{align} \int d^3x'\frac{\rho(r')}{\lvert \mathbf{x'} - \mathbf{x}\rvert}=&\int r'^2 dr'\int d\Omega'\rho(r')\frac{1}{\lvert \mathbf{x'} - \mathbf{x}\rvert}\\ =&4\pi\int \rho(r')r' dr'\\ =&4\pi\int r'^{1-\gamma} dr' \end{align}$$

which yields your expected result.

Improve this answer
$\endgroup$
3
  • 1
    $\begingroup$ Fantastic, thanks! I will have to look back into spherical harmonics in order to understand the angular integration you performed. Will do so now! $\endgroup$ Commented Mar 4, 2020 at 13:18
  • 1
    $\begingroup$ Soner, can you show how you did this with the Wolfram Language? $\endgroup$ Commented Mar 4, 2020 at 18:47
  • 1
    $\begingroup$ @CATrevillian I am unaware of how to do intermediate steps I wrote above with the Wolfram Language. $\endgroup$ Commented Mar 5, 2020 at 4:47

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.