# Computational Evaluation of Shock Wave Interaction with a Liquid Droplet
## 2 Computational Model
### CFD Analysis
#### 2.2.2 CFD Solution  
The mean turbulent flow in the virtual shock tube device was simulated by solving the unsteady compressible RANS equations. The shear-stress transport (SST) \(k\)-\(\omega\) two-equation eddy-viscosity model was used for turbulence closure, owing to its proven suitability for complex fluids engineering applications [27]. The flow governing equations, which are not reported here for brevity, can be found, for instance, in [28]. The numerical simulations were performed by using the industrial solver Ansys Fluent, which has been successfully employed in analogous works [13; 19], as well as in previous CFD studies using RANS models by the same research group [29; 30].  
The present pressure-based solver utilized the FV approach to discretize the governing equations and approximate the unsteady mean flow solution, where the conservation principles were applied over each computational cell [31; 32]. The FV grid was suitably refined in the flow region, just ahead of the diaphragm, that is interested by the SW formation and occupied by the deforming water body, including the near wake, as sketched in Figure 1. Moreover, in order to save computing power, the computational domain that was used corresponded to half the cylindrical shock tube domain, namely \(y\geq 0\), where a symmetry boundary condition was imposed at the centerplane. In fact, according to previous experimental visualizations, the flow can be assumed symmetric at the early stages of the interaction process [22]. Wall boundary conditions were employed for simulating the solid walls of the virtual shock tube device, as was done in similar studies, e.g., [12]. The second-order upwind scheme was used for spatial discretization, while the constant time-step of \(\Delta t=2\times 10^{-7}\) s was used for the explicit transient calculation.  
A grid sensitivity analysis was performed by examining the shock front and the postshock velocities predicted using three different FV grids with increasing resolution, as reported in Table 4. Therein, \(h_{min}\) stands for the minimum characteristic mesh size, which was evaluated as the cubic root of the cell volume. The total number of FV elements ranged between 6.5 and 15 million, where the finest local resolution resulted in 90 cells per original droplet diameter.  
The system of integral balance equations for mass, momentum, and energy was completed by employing two equations of state to model the density variations of the two different fluids. Specifically, along with the ideal gas law for the gaseous phase, the following equation of state for the liquid phase was considered [33],  
\[\frac{\rho_{l}}{\rho_{l0}}=\left(\frac{p+\mathcal{B}}{p_{0}+\mathcal{B}}\right) ^{\frac{1}{\kappa}}, \tag{7}\]  
where \(p_{0}\) and \(\rho_{l0}\) represent the reference pressure and density levels, while \(\mathcal{B}\) and \(\kappa\) are constant parameters. Following [12], the pressure-like parameter \(\mathcal{B}\) and the so-called adiabatic index \(\kappa\) were set to the values of 305 MPa and 6.68, respectively. However, given the maximum pressure level expected in the simulations, the liquid density can be considered practically constant for the present engineering analysis.  
\begin{table}
\begin{tabular}{c c c c c} \hline \hline
**Solution** & **\# of FV Cells** & \(h_{min}/d_{0}\) & \(V_{\mathbf{S}}\) (m/s)** & \(V_{\mathbf{2}}\) (m/s) \\ \hline CFD I & \(6.5\times 10^{6}\) & \(3.4\times 10^{-2}\) & \(504.1\) & \(225.6\) \\ CFD II & \(9\times 10^{6}\) & \(2.5\times 10^{-2}\) & \(504.3\) & \(225.7\) \\ CFD III & \(15\times 10^{6}\) & \(1.1\times 10^{-2}\) & \(504.4\) & \(225.8\) \\ Analytical & \(-\) & \(-\) & \(504.5\) & \(225.9\) \\ \hline \hline \end{tabular}
\end{table}
Table 4: Shock tube flow: traveling velocities for different mesh resolutions.