## 2 Smoothed Particle Hydrodynamics  
Smoothed particle hydrodynamics is a meshfree computational method initially developed by [PERSON] [11] and [PERSON] [12] for solving astrophysical problems. Later it was used to solve fluidynamics problems to overcome some of the limitations of the grid-based method in the case, for instance, of explosions and high velocity impact phenomena ([13; 14]). The method was also widely validated against shock waves in particular for the well known Riemann problem ([15; 16; 17]). The SPH approximation is based on the so-called integral representation of a function. Given the function \(f(\mathbf{r})\), defined in a volume \(V\), function of the three-dimensional position \(\mathbf{r}\), is defined by the identity  
\[f(\mathbf{r})=\iiint f(\mathbf{r}^{\prime})\delta(\mathbf{r}-\mathbf{r}^{ \prime})d\mathbf{r}^{\prime}, \tag{1}\]  
where \(\delta(\mathbf{r}-\mathbf{r}^{\prime})\) is the Dirac delta function defined as  
Figure 1: Geometry of the simulation box. Figure 2: Geometry of the shock tube with the cylindrical inhomogeneity. \[\delta(\mathbf{r}-\mathbf{r}^{\prime})= \begin{cases}+\infty&\mathbf{r}=\mathbf{r}^{ \prime}\\ 0&\mathbf{r}\
eq\mathbf{r}^{\prime}.\end{cases}  \tag{2}\]  
It is possible to approximate the integral representation by replacing the Dirac delta function with a bell-shaped function called smoothing function or kernel, \(W\), which depends on the position \(\mathbf{r}\) and on the so-called smoothing length, \(h\). When \(h\) approaches zero the kernel function \(W\) has the property  
\[\lim_{h\to 0}=W(\mathbf{r}-\mathbf{r}^{\prime},h)=\delta(\mathbf{r}-\mathbf{r}^{ \prime}), \tag{3}\]  
which approximates the integral representation, Equation (1), to the so-called kernel approximation:  
\[f(\mathbf{r})\approx\iiint f(\mathbf{r}^{\prime})W(\mathbf{r}-\mathbf{r}^{ \prime},h)d\mathbf{r}^{\prime}. \tag{4}\]  
In this work, the kernel used is Lucy kernel function  
\[W(R,h)= \begin{cases}\frac{1}{s}\left[1+3\frac{R}{h}\right]\left[1-\frac{R}{h} \right]^{3}&R\leq 1\\ 0&R>1,\end{cases}  \tag{5}\]  
where \(R=|\mathbf{r}-\mathbf{r}^{\prime}|\) and \(s\) is a parameter used to normalise the kernel function, which, for one, two and three dimensional space is, respectively, \(\frac{4h}{5}\), \(\frac{\pi h^{2}}{5}\) and \(\frac{16\pi h^{3}}{105}\). The last step is to approximate the infinitesimal volume \(d\mathbf{r}^{\prime}\) to a finite volume \(dr\) composed by computational particles with their own mass \(m=\rho dr\). With this approximation it is possible to discretise Equation (4)  
\[f(\mathbf{r})\approx\sum\frac{m_{i}}{\rho_{i}}f(\mathbf{r}_{i})W(\mathbf{r}- \mathbf{r}_{i},h), \tag{6}\]  
where \(m_{i}\), \(\rho_{i}\) and \(\mathbf{r}_{i}\) are mass, density and position of the \(i\)th particle. Only particles for which \(|\mathbf{r}-\mathbf{r}_{i}|<h\) are taken into account in the summation. With Equation (6) it is possible to discretise any set of equations such as the energy balance or the Navier-Stokes equation on an arbitrarily set of computational particles. Within the SPH framework, for instance, the momentum-conservation equation can be rewritten as  
\[m_{i}\frac{d\mathbf{v}_{i}}{dt}=\sum_{j}m_{i}m_{j}\left(\frac{P_{i}}{\rho_{i} }+\frac{P_{i}}{\rho_{i}}+\Pi_{ij}\right)\
abla_{j}W_{ij}, \tag{7}\]  
where \(W_{ij}=W(r_{j}-r_{i},h)\) and \(\
abla_{j}W_{ij}\) is the kernel gradient in the \(r_{j}\) direction. \(P\) is the pressure while \(\Pi_{ij}\) is the so-called artificial viscosity introduced by [PERSON] (1958) for simulating shock waves  
\[\Pi_{ij}=-ah\frac{c_{i}+c_{j}}{\rho_{i}+\rho_{j}}\frac{\mathbf{v}_{ij}\cdot \mathbf{r}_{ij}}{r_{ij}^{2}+\epsilon h^{2}}, \tag{8}\]where \(c_{i}\) and \(c_{j}\) are speed of sound of particles \(i\) and \(j\) and \(\alpha\) is dimensionless parameter that controls the strength of the viscous dissipation and \(\epsilon\approx 0.01\) is used to avoid singularities when particles are close to each other.