In General Relativity Lecture Notes by Javier Rubio, Chapter 6-Einstein Equations we read
6.1 The energy-momentum tensor
$\cdots\cdots\cdots\cdots\cdots\cdots\cdots\cdots$
6.2 The microscopic description
The relation between $\;\rho\;$ and $\;p\;$ is usually characterized by an equation of state $\;p=p(\rho)\;$ which depends on the microscopic particles involved in the fluid. In order to get some insight about the possible equations of state, let me consider a macroscopic collection of N structureless point particles interacting through spatially localized collisions. The energy density associated to any of them is given by \begin{equation} T^{00}_{n}=E_{n}\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)=m_{n}\gamma_{n}\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)\,, \tag{6.19}\label{6.19} \end{equation} with $\;\gamma_n=1/\sqrt{1-\upsilon^2_n\,}\;$ and $\;n=1,...,N\;$ a label selecting the particular particle we are referring to. Taking into account the identity \begin{align} \int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau\,\delta^{(4)}\left(x-x\left(\tau\right)\vphantom{A^2}\right) & = \int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau\,\delta\left(t-t\left(\tau\right)\vphantom{A^2}\right)\delta^{(3)}\left(\mathbf x-\mathbf x\left(\tau\right)\vphantom{A^2}\right) \nonumber\\ & = \dfrac{\mathrm d\tau}{\mathrm dt}\delta^{(3)}\left(\mathbf x-\mathbf x\left(t\right)\vphantom{A^2}\right)=\dfrac{1}{\gamma}\delta^{(3)}\left(\mathbf x-\mathbf x\left(t\right)\vphantom{A^2}\right)\,, \tag{6.20}\label{6.20} \end{align} the non-Lorentz invariant 3-dimensional Dirac delta appearing in equation \eqref{6.21} can be transformed into a Lorentz invariant 4-dimensional Dirac delta5 \begin{equation} T^{00}_{n}=m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\,u^0_n u^0_n\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right) \,. \tag{6.21}\label{6.21} \end{equation} The same procedure can be applied to the spatial momentum density (or energy current) of the particle \begin{equation} \begin{split} T^{0i}_{n}&=p^i_n\,\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)\\ &=m_n\gamma_n\upsilon^i_n\,\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)=E_n\upsilon^i_n\,\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)\,,\\ \end{split} \tag{6.22}\label{6.22} \end{equation} and to the flux of the $\,i\,$ momentum component in the $\,j\,$ direction (or viceversa) \begin{equation} T^{i j}_{n}=p^i_n\upsilon^j_n\,\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)=p^j_n\upsilon^i_n\,\delta^{(3)}\left(\mathbf x-\mathbf x_n\left(t\right)\vphantom{A^2}\right)\,. \tag{6.23}\label{6.23} \end{equation} We obtain \begin{equation} \begin{split} T^{0i}_{n}&=m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\,u^0_n u^i_n\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right) \,,\\ T^{i j}_{n}&=m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\,u^i_n u^j_n\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right) \,.\\ \end{split} \tag{6.24}\label{6.24} \end{equation} Equations \eqref{6.21} and \eqref{6.24} can be rewritten in a very compact way in terms of the stress-energy-momentum tensor $\,T^{\mu\nu}_n\,$ \begin{equation} \begin{split} T^{\mu\nu}_n &=m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\,u^\mu_n u^\nu_n\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\\ &=\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\dfrac{p^\mu_n p^\nu_n}{m_{n}}\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right) \,,\\ \end{split} \tag{6.25}\label{6.25} \end{equation} which is manifestly symmetric and Lorentz invariant since $\,u^\mu_n u^\nu_n\,$ is a tensor under Lorentz transformations and both $\,m_{n}\,$ and $\,d\tau_n\,u^\mu_n u^\nu_n\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\,$ are Lorentz scalars. The total energy density of the whole system of particles can be written as the sum of the individual contributions, namely \begin{equation} T^{\mu\nu}=\sum\limits_{n=1}^{N}T^{\mu\nu}_n \,. \tag{6.26}\label{6.26} \end{equation}
6.2.1 Energy-momentum tensor conservation and geodesics
Let us see under which conditions the total energy momentum tensor (6.26) is conserved. Taking the derivative with respect to the coordinates we get \begin{equation} \partial_\mu T^{\mu\nu}=\sum\limits_{n=1}^{N}m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\,u^\mu_n\, u^\nu_n\,\partial_\mu\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\,, \tag{6.27}\label{6.27} \end{equation} which using \begin{equation} \begin{split} u^\mu_n\partial_\mu\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)&=\dfrac{\mathrm d x^\mu_n}{\mathrm d\tau_n}\dfrac{\partial}{\partial x^{\mu}}\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\\ &= -\dfrac{\mathrm d}{\mathrm d\tau_n}\left[\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\vphantom{\dfrac{a}{b}}\right] \,,\\ \end{split} \tag{6.28}\label{6.28} \end{equation} can be written as \begin{equation} \begin{split} \partial_\mu T^{\mu\nu}=&-\sum\limits_{n=1}^{N}m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\,\dfrac{\mathrm d}{\mathrm d\tau_n}\left[u^{\nu}_{n}\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\vphantom{\dfrac{a}{b}}\right]\\ &+\sum\limits_{n=1}^{N}m_{n}\!\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\, \dot u^\nu_n\,\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)\,.\\ \end{split} \tag{6.28a}\label{6.28a} \end{equation} The first term in the right hand side of the previous expression disappears in the particles are stable, i.e. if the orbits are closed or come from negative infinite time and disappear into positive infinite time. We are left then with the second term, which can be written as \begin{equation} \partial_\mu T^{\mu\nu}=\sum\limits_{n=1}^{N}\int\limits_{-\infty}^{+\infty}\!\!\mathrm d\tau_n\dfrac{\mathrm d p^\nu_n}{\mathrm d\tau_n}\delta^{(4)}\left(x-x_n\left(\tau_n\right)\vphantom{A^2}\right)=\sum\limits_{n=1}^{N}\dfrac{\mathrm d p^\nu_n}{\mathrm d t}\delta^{(3)}\left(\mathbf x-\mathbf x_n\vphantom{A^2}\right)\,, \tag{6.29}\label{6.29} \end{equation} with $\,p^\nu_n=m_n u^\nu_n\,$ the 4-momentum of the individual particles. The local energy momentum conservation $\,\partial_\mu T^{\mu\nu}=0\,$ requires the particles to be free. Or in other words, the condition $\,\partial_\mu T^{\mu\nu}=0\,$ is equivalent to the geodesic equation in Minkowski space-time, $\,\mathrm d p^\mu_n/\mathrm d\tau =0\,$. This will be also the case in curved spacetime.
$\text{------------------------------------------------------------------------------------}$
5 The fact that the 4-Dimensional Dirac delta $\;\delta^{(4)}\left(x\right)\;$ is Lorentz invariant follows directly from the definition $\;\int\mathrm d^4x\delta^{(4)}=1\;$ and the fact that the volume element $\;\mathrm d^4x\;$ is Lorentz invariant.
$\text{------------------------------------------------------------------------------------}$
I don't see how we got equation \eqref{6.27}. And how after using equation \eqref{6.28} we get equation \eqref{6.28a}?
Also for the \eqref{6.28} equation, can someone explain to me how this simplification can be done?