## CFD Modeling
### Flow Configuration  
The flow geometry practically reproduces a shock tube device with a square cross-section, where the spatial coordinates system \((x,y,z)\) is chosen with the first axis direction along the tube axis, corresponding to the traveling direction of the shock front. Based on previous numerical studies conducted for similar flow configurations, the extent of the computational domain is chosen such that \(-30<x/d_{0}<40\), \(|y|/d_{0}<30\). and \(|z|/d_{0}<30\), where \(d_{0}=4.8\) mm stands for the initial diameter of the water sphere. The current droplet size represents one of the most studied cases, corresponding to the classical experiments of [PERSON] and [PERSON] [7]. The above domain size serves as the necessary compromise between numerical accuracy and computational costs, as it could be determined by industrial researchers and developers. Following previous works, this size of the computational domain can be considered sufficiently large to avoid non-physical wave reflections on the boundaries [32]. In this study, differently from [19], the shock tube flow is not explicitly simulated, but the discontinuous air flow conditions are directly imposed. Practically, the experimental conditions are mimicked by initially setting two different thermodynamic states for the air occupying the two different subdomains that are separated by a virtual cross diaphragm perpendicular to the \(x\)-axis. The simulations are initiated with the shock front being positioned one diameter away from the droplet center, while the origin of the coordinate system corresponds to the initial position of the droplet leading edge. Practically, given the atmospheric pre-shock conditions in a section, the high-speed post-shock conditions in the other one uniquely correspond to the prescribed Mach number. The latter parameter was selected as \(\text{Ma}=1.47\), representing one of the cases mostly investigated in the engineering and scientific literature, both computationally [10; 15] and experimentally [7]. As it is empirically demonstrated, the proposed model can reproduce a plane shock wave traveling at the prescribed velocity that is \(V_{\text{S}}=505\) m/s. For brevity, the airflow parameters corresponding to the present simulations are summarized in Table 1. It is worth noting that owing to the actual post-shock thermodynamic conditions, the ambient air pressure results in being greater than the saturation pressure of water, and thus, the liquid body does not undergo a phase change. In practice, the two different fluids can be considered immiscible, which could be no more valid for increasing Mach number. The density and the dynamic viscosity of the liquid phase are initially set to \(\rho_{0}=998\,\text{kg/m}^{3}\) and \(\mu_{\text{0}}=1.00\times 10^{-3}\) Pa s, while the surface tension at the gas-liquid interface is assumed equal to \(\sigma=7.29\times 10^{-2}\) N/m. This way, the characteristic Ohnesorge and Weber numbers result in being \(\text{Oh}=\mu_{\text{0}}/\sqrt{\rho_{\text{0}}\sigma d_{0}}=1.64\times 10^{-4}\) and \(\text{We}=\rho_{\text{g}}\mu_{\text{g}}^{2}d_{0}/\sigma=7.3\times 10^{3}\), respectively. According to these values, the present configuration corresponds to the droplet aerobreakup under the SIE regime [8]. Moreover, based on the high Weber number, as well as the high Reynolds number, that is, \(\text{Re}=\rho_{\text{g}}u_{\text{g}}d_{0}/\mu_{\text{g}}=1.06\times 10^{5}\), both capillary and viscous effects might be neglected in practical calculations. Indeed, similar numerical studies have been mostly conducted by solving the compressible Euler equations without including the surface tension term [15]. For the sake of completeness, these effects are not ignored in the