Calculating the point charge potential via Gauss's law I want to calculate the potential of a point charge via Gauss's law
$$\vec{\nabla}\cdot\vec{E}=4\pi\rho(\vec{r})$$
My Idea:
we can find an expression for $4\pi$ that has the form of an integral,
by using spherical coordinates we can derive:
$$\int_{V}\vec{\nabla}^2 \left( \frac{1}{r} \right)dV =\int_V \vec{\nabla}\cdot\left( \vec{\nabla}  \frac{1}{r} \right)dV$$
by using the Divergence theorem and applying spherical coordinates we obtain:
$$\oint_S \vec{\nabla} \left(\frac{1}{r}\right) \cdot \vec{n}dS = \oint_S \partial_r \left( \frac{1}{r} \right)r^2\sin(\theta) d\theta dr = -4\pi$$
$$\Rightarrow \vec{\nabla}\cdot \vec{E}=-\int_V \vec{\nabla}^2 \frac{\rho(\vec{r})}{|\vec{r}-\vec{r}'|}dV = -\vec{\nabla}^2\int_V  \frac{\rho(\vec{r})}{|\vec{r}-\vec{r}'|}dV$$
since following expression is known:
$$\vec{\nabla}\cdot \vec{E}=-\vec{\nabla}^2\phi$$
we can deduce that:
$$\Rightarrow \phi(\vec{r})=\int_V \frac{\rho(\vec{r}')}{|\vec{r}-\vec{r}'|}dV$$
My problem is that I'm not sure how to solve this integral. I'm aware that the charge density for discrete charges is:
$$\rho(\vec{r})=\sum_i^nq_i \delta(\vec{r}-\vec{r}')$$
and since we are talking about a single charge at the origin ($\vec{r}'=0$) this sum will equal $q$.
$$\Rightarrow \phi(\vec{r})=\int_V \frac{q}{r}dV=q\int_0^{2\pi}\int_0^\pi \int_0^R \frac{1}{r}r^2\sin(\theta) \space d\theta d\varphi dr$$
However, this delivers:
$$\phi(\vec{r})=2q\pi r^2$$
Which is of course wrong.
Any help is appreciated
 A: Note that $\mathbf{\nabla}\cdot\mathbf{E}=4\pi\rho(\mathbf{r})$ is equivalent to $$\begin{align}
\int_V \mathbf{\nabla}\cdot\mathbf{E}\, d\tau & = 4\pi\int_V \rho\,d\tau \\
\oint_S \mathbf{E}\cdot d\mathbf{a}&=4\pi \,Q_{enc}
\end{align}$$ 
by the divergence theorem. Define a sphere around the point. Then, $\mathbf{E}$ is spherically symmetric, constant, and comes out of the integral. If we say the point has charge $q$, then we have 
$$
\begin{align}
\oint_S |\mathbf{E}|\,da &= 4\pi q \\
|\mathbf{E}|\oint da &= 4\pi q \\
|\mathbf{E}|\cdot 4\pi r^2 &= 4\pi q \\
|\mathbf{E}| & = \frac{q}{r^2} \\
\mathbf{E} & = \frac{q}{r^2}\hat{\mathbf{r}}
\end{align}
$$
Now, we have $$V(\mathbf{r})=-\int^\mathbf{r}_\infty \mathbf{E}\cdot d\mathbf{l}$$
Therefore, $$\begin{align}
V(\mathbf{r})&=-\int^\mathbf{r}_\infty \frac{q}{r'^2}\, dr' \\
&= \frac{q}{r'}\bigg|_\infty^r \\
&= \frac{q}{r}
\end{align}$$
A: Your reasoning is correct until your equation
$$\rho(\vec{r})=\sum_i^nq_i \delta(\vec{r}-\vec{r}'_i)$$
(I've added an index $i$ here for clarity,
because $\vec{r}'_i$ denotes the place of charge $q_i$)
But then, for the special case of a single charge $q$
at the origin ($\vec{r}' = \vec{0}$) this sum will be equal
$$\rho(\vec{r})=q \delta(\vec{r})$$
From this it follows by integration
$$\begin{align}
\phi(\vec{r})&=\int_V \frac{\rho(\vec{r}')}{|\vec{r}-\vec{r}'|}dV \\
&=\int_V \frac{q\delta(\vec{r}')}{|\vec{r}-\vec{r}'|}dV \\
&=\frac{q}{|\vec{r}|}
\end{align}$$
(The last step here used the definition of the $\delta$ function)
