Volume charge density of H-atom

I have a problem where I am supposed to calculate the volume charge density of a neutral hydrogen atom. The potential is given to be $$\Phi = k \frac{e^{-ar}}{r} \left(1 + \frac{ar}{2}\right)$$ Now I tried to use the Poisson-Equation stating $$\Delta \Phi = \frac{\rho}{\varepsilon_0}$$ which leads me to $$\rho = \Delta \left( \underbrace{\frac{q}{4 \pi \epsilon_0}}_{= k} \frac{e^{-\alpha r}}{r} \left( 1 + \frac{\alpha r}{2} \right) \right) = k \Delta \Big( \frac{e^{-\alpha r}}{r} + \frac{\alpha e^{-\alpha r}}{2} \Big) = k \left( \Delta \left( \frac{e^{-\alpha r}}{r} \right) + \frac{\alpha}{2} \Delta \left( e^{-\alpha r} \right) \right)$$ Now I define $f = e^{-\alpha r}$ and $g = \frac{1}{r}$. The Laplacian of the product $fg$ is then $$\Delta(fg) = g \Delta(f) + f \Delta(g) + \nabla (f) \cdot \nabla (g)$$ and the derivatives are $$\nabla(f) = - \alpha e^{-\alpha r} \hat{ \mathbf r} \qquad \qquad \Delta(f) = \alpha^2 e^{-\alpha r}$$ $$\nabla(g) = - \frac{1}{r^2} \hat{ \mathbf r} \qquad \qquad \Delta(g) = - 4 \pi \delta(r)$$ $$\implies \nabla(f) \cdot \nabla(g) = \frac{\alpha e^{-\alpha r}}{r^2}$$ Inserting this back into the original equation yields $$\rho = k e^{-\alpha r}\Big( \frac{\alpha^2}{r} - 4 \pi \delta(r) + \frac{\alpha}{r^2} + \frac{\alpha^3}{2} \Big)$$ However, this seems somewhat wrong to me since I would have expected the expression to be increasing from the origin and then decreasing after some $r=R$ since the potential of the electron hull should take over.

Can anybody either confirm that this is correct or show me where I made the mistake?

Apart from taking the derivatives like in cartesian coordinates, I have tried calculating the Laplacian by calculation in spherical coordinates as well using the spherical Laplacian $$\Delta \Phi = \frac{1}{r^2} \frac{\partial}{\partial r} \left(r^2 \frac{\partial \Phi}{\partial r}\right)$$ but still got the same result.

• Always worth remembering that the $n = 0$ electron wave-function is at it's maximum at the center... – dmckee Sep 11 '13 at 16:51
• So you are saying that the calculation is sort of consistent since my electron charge has the highest density at $r=0$? @dmckee – user29498 Sep 11 '13 at 17:07

It's not clear to me exactly why you're unhappy with the answer you get. I would suggest phrasing your expectations in terms of the total charge contained in a sphere of radius $r$, $$4\pi\int_{0^-}^r \rho(r')r'^2\,\text dr'.$$ This should give the positive charge of the nucleus at $r\rightarrow 0^+$ (because the proton is point-sized in this model!) and then drop monotonically to zero as $r$ increases through $1/\alpha$ and beyond to $r\rightarrow\infty$, and your sphere includes more and more of the electron cloud that surrounds (and neutralizes) the nucleus. If this fails then you definitely need to check your calculations.
• Also, note that confusion between $\delta(r)$ and $\delta(\mathbf r)$ can be quite harmful here. – Emilio Pisanty Sep 11 '13 at 17:27
• Ok, thank you all, I have redone the calculations with first considering r>0 and afterwards letting r->0 and a Taylor expansion of the potential at $r=0$ Now I have $$\rho = - \varepsilon_0 \Phi = q(\delta(r) - \frac{\alpha^3}{8 \pi} e^{-\alpha r})$$ Which seems a bit better to me because it does no longer contain the $\frac{1}{r}$ and $\frac{1}{r^2}$ terms – user29498 Sep 11 '13 at 21:55
• No, the singular terms should be there (or at least they don't hurt). It doesn't matter that $\rho$ have a singularity as long as the singularity is integrable, which means that singularities of the type $r^{-s}$ all the way up to with $s<3$ are allowed. – Emilio Pisanty Sep 11 '13 at 22:23