(a)
\left\langle\frac{1}{\left|r_{1}-r_{2}\right|}\right\rangle=\left(\frac{8}{\pi a^{3}}\right)^{2} \int \underbrace{\left[\int \frac{e^{-4\left(r_{1}+r_{2}\right) / a}}{\sqrt{r_{1}^{2}+r_{2}^{2}-2 r_{1} r_{2} \cos \theta_{2}}} d^{3} r _{2}\right]}_{\blacklozenge} d^{3} r _{1} .
\blacklozenge=2 \pi \int_{0}^{\infty} e^{-4\left(r_{1}+r_{2}\right) / a} \underbrace{\left[\int_{0}^{\pi} \frac{\sin \theta_{2}}{\sqrt{r_{1}^{2}+r_{2}^{2}-2 r_{1} r_{2} \cos \theta_{2}}} d \theta_{2}\right]}_{\star} r_{2}^{2} d r_{2} .
\star=\left.\frac{1}{r_{1} r_{2}} \sqrt{r_{1}^{2}+r_{2}^{2}-2 r_{1} r_{2} \cos \theta_{2}}\right|_{0} ^{\pi}=\frac{1}{r_{1} r_{2}}\left[\sqrt{r_{1}^{2}+r_{2}^{2}+2 r_{1} r_{2}}-\sqrt{r_{1}^{2}+r_{2}^{2}-2 r_{1} r_{2}}\right] .
=\frac{1}{r_{1} r_{2}}\left[\left(r_{1}+r_{2}\right)-\left|r_{1}-r_{2}\right|\right]=\left\{\begin{array}{l} 2 / r_{1}\left(r_{2}<r_{1}\right) \\ 2 / r_{2}\left(r_{2}>r_{1}\right) \end{array}\right. .
\blacklozenge=4 \pi e^{-4 r_{1} / a}\left[\frac{1}{r_{1}} \int_{0}^{r_{1}} r_{2}^{2} e^{-4 r_{2} / a} d r_{2}+\int_{r_{1}}^{\infty} r_{2} e^{-4 r_{2} / a} d r_{2}\right] .
\frac{1}{r_{1}} \int_{0}^{r_{1}} r_{2}^{2} e^{-4 r_{2} / a} d r_{2}=\left.\frac{1}{r_{1}}\left[-\frac{a}{4} r_{2}^{2} e^{-4 r_{2} / a}+\frac{a}{2}\left(\frac{a}{4}\right)^{2} e^{-4 r_{2} / a}\left(-\frac{4 r_{2}}{a}-1\right)\right]\right|_{0} ^{r_{1}} .
=-\frac{a}{4 r_{1}}\left[r_{1}^{2} e^{-4 r_{1} / a}+\frac{a r_{1}}{2} e^{-4 r_{1} / a}+\frac{a^{2}}{8} e^{-4 r_{1} / a}-\frac{a^{2}}{8}\right] .
\int_{r_{1}}^{\infty} r_{2} e^{-4 r_{2} / a} d r_{2}=\left.\left(\frac{a}{4}\right)^{2} e^{-4 r_{2} / a}\left(-\frac{4 r_{2}}{a}-1\right)\right|_{r_{1}} ^{\infty}=\frac{a r_{1}}{4} e^{-4 r_{1} / a}+\frac{a^{2}}{16} e^{-4 r_{1} / a} .
\blacklozenge=4 \pi\left\{\frac{a^{3}}{32 r_{1}} e^{-4 r_{1} / a}+\left[-\frac{a r_{1}}{4}-\frac{a^{2}}{8}-\frac{a^{3}}{32 r_{1}}+\frac{a r_{1}}{4}+\frac{a^{2}}{16}\right] e^{-8 r_{1} / a}\right\} .
=\frac{\pi a^{2}}{8}\left\{\frac{a}{r_{1}} e^{-4 r_{1} / a}-\left(2+\frac{a}{r_{1}}\right) e^{-8 r_{1} / a}\right\} .
\left\langle\frac{1}{\left|r_{1}-r_{2}\right|}\right\rangle=\frac{8}{\pi a^{4}} \cdot 4 \pi \int_{0}^{\infty}\left[\frac{a}{r_{1}} e^{-4 r_{1} / a}-\left(2+\frac{a}{r_{1}}\right) e^{-8 r_{1} / a}\right] r_{1}^{2} d r_{1} .
=\frac{32}{a^{4}}\left\{a \int_{0}^{\infty} r_{1} e^{-4 r_{1} / a} d r_{1}-2 \int_{0}^{\infty} r_{1}^{2} e^{-8 r_{1} / a} d r_{1}-a \int_{0}^{\infty} r_{1} e^{-8 r_{1} / a} d r_{1}\right\}.
=\frac{32}{a^{4}}\left\{a \cdot\left(\frac{a}{4}\right)^{2}-2 \cdot 2\left(\frac{a}{8}\right)^{3}-a \cdot\left(\frac{a}{8}\right)^{2}\right\}=\frac{32}{a}\left(\frac{1}{16}-\frac{1}{128}-\frac{1}{64}\right)=\frac{5}{4 a} .
(b)
V_{e e} \approx \frac{e^{2}}{4 \pi \epsilon_{0}}\left\langle\frac{1}{\left|r_{1}-r_{2}\right|}\right\rangle=\frac{5}{4} \frac{e^{2}}{4 \pi \epsilon_{0}} \frac{1}{a} =\frac{5}{4} \frac{m}{\hbar^{2}}\left(\frac{e^{2}}{4 \pi \epsilon_{0}}\right)^{2}=\frac{5}{2}\left(-E_{1}\right)=\frac{5}{2}(13.6 eV )=34 eV .
E_{0}+V_{e e}=(-109+34) eV =-75 eV , which is pretty close to the experimental value (-79 eV).