Why is the Schwarzschild black hole a vacuum solution? It is always said that Schwarzschild black holes are solutions of the vacuum Einstein equation, but I don't quite understand why. I understand that $R_{{\mu}{\nu}}=0$ at each and every spacetime point except for $r=0$ and thus, clearly, $T_{{\mu}{\nu}}=0$ at each and every point except for $r=0$. But since the $r=0$ set of points can not be said to be excluded from the considered spacetime (because a particle at a finite $r$ would reach there in finite proper-time under a free-fall, we must care about $r=0$ as well. And, certainly, $R_{{\mu}{\nu}}$ doesn't vanish there but rather it blows up. So, shouldn't Schwarzschild solution be regarded as a solution of Einstein equations in a spacetime with some peculiar infinite density at a point and vacuum everywhere else? From some of the references I have browsed, I have a hint that I might be confusing the astrophysical blackholes with the vacuum solution blackholes, but I don't exactly understand where and how. 
 A: General relativity is a theory grounded in a fairly rigorous mathematical model. In this model, the Ricci tensor is a function of the form
$$\text{Ric} : T\mathcal{M} \times T\mathcal{M} \to \Bbb R$$
mapping two elements of the tangent space of the manifold to the real numbers. The Ricci tensor is defined in term of the Riemann tensor. For the Schwarschild metric, this will involve components of the form $\approx r^{-n}$, which is only defined on $\Bbb R \setminus \{0\}$. Hence within the framework of general relativity, it is impossible to add this point to the manifold. It is what is referred to as a singular boundary point, and more specifically a scalar singularity. For those, there exists no spacetimes $(\mathcal M', g')$ that can extend the existing spacetime $(\mathcal M, g)$ and remove those singularities. 
It is possible to find extensions of general relativity where the manifold is defined at all points but in which case you have to forgo tensors being defined at points, using a more distributional approach, in which case the Ricci tensor for instance would be of the form
$$\text{Ric} : \mathcal G(\mathcal M, T\mathcal{M} \times T\mathcal{M}) \to \Bbb R$$
where $\mathcal G$ is a generalized section of a vector bundle by the way of some Colombeau generalized functions on the manifold. In this case, the Ricci tensor can be calculated to yield something of the form (for the $tt$ component)
$$R^t_t = -4\pi \delta_\varepsilon(r)$$
where $\delta_\varepsilon$ is a mollified delta function in the sense of Colombeau algebras. While it does have a value at $r = 0$ in a general sense, that value isn't actually in $\Bbb R$, and it should be more understood as a distribution. 
While you can do all that, there isn't much point in doing it as it makes things mostly more complicated without actually bringing much illumination to the whole affair, outside of the fact that you can describe the Schwarzschild solution as the distribution of a point particle. 
