The instantaneous detonation model (IDM) is widely used in simulating underwater explosions due to its efficiency and ability to ignore the detonation reaction process. In this study, we propose a new IDM to predict the fluid structure in the detonation zone of an RS211 explosive charge. This model is based on high-order solutions provided by the detonation shock dynamics model, where the spatial term is discretized using fifth-order WENO reconstruction in characteristic space and Lax–Friedrichs’s splitting and the temporal terms are discretized using a third-order TVD Runge–Kutta scheme. The interface motion is captured using the level-set method combined with MGFM, and a programmed burn model is provided to describe the generation and propagation of the detonation wave. The self-similarity of detonation wave propagation is validated, and the quantitative calculation formula of the instantaneous detonation model is obtained by averaging or curve fitting the dimensionless results. Consequently, the IDM of the RS211 charge is established using high-order polynomial approximations of the Taylor rarefaction zone and a constant static zone for 1D planar, cylindrical, and spherical RS211 charges. The application of the IDM involves direct mapping from the radial direction to the spatial structured grid for 1D planar, 2D cylindrical, and 3D spherical charges. Numerical results demonstrate that the IDM proposed in this paper shows good accuracy and high computational efficiency.

## I. INTRODUCTION

Research on underwater explosions has spanned more than a century and has involved experimental, theoretical, and simulation investigations. Underwater explosions can be divided into three critical stages: the detonation wave, shock wave, and bubble pulsation stages.^{1–3} The generation and propagation of detonation waves in condensed-phase explosives are crucial to near-field underwater explosions. Discontinuities in shock waves and interfaces between multiphase materials with high density and pressure ratios can cause oscillations and computational collapse. When strong shock waves propagate into a solid energetic material, detonation occurs, and the surrounding media are severely compressed. The detonation wave is sustained and continues to propagate through the unreacted energetic material by the release of chemical energy.^{4–6} Subsequently, the detonation wave rapidly reaches the surface of the charge and produces a blast wave in the surrounding water or air. Baer and Nunziato^{7} proposed a two-phase model to study the deflagration-to-detonation (DTD) transition in solid heterogeneous materials, which has been used to model shock waves in compressible mixtures and detonation waves in granular materials. In this model, the solid and gas phases are not required to be in equilibrium, and the control equations and equations of state for each phase are independent. Saurel and Lemetayer^{8} introduced the concept of infinite relaxation parameters and instantaneous velocity and pressure equilibrium to study the DTD transition in solid energetic materials. Price *et al.*^{9} proposed a volume-fraction type approach to simulate explosive detonations and airblast applications in multi-material flows governed by the compressible Euler equations on three-dimensional unstructured grids in the five-equation system.

When the detonation wave reaches the surface of the explosive, the high-temperature and high-pressure detonation products not only produce a shock wave in the surrounding water medium but also induce topological motion at the explosion bubble interface. Experimental research on explosive detonations faces several challenges, such as high costs, difficulty in measurement, and limited data availability. With advancements in computational fluid dynamics (CFD) and computer science, numerical simulation methods have gained significant importance. The detonation phenomenon involves complex interactions between fluid mechanics and chemical reaction kinetics. Currently, two main methods are used to address the detonation process: the Detonation Shock Dynamics (DSD) model and the Instantaneous Detonation Model (IDM).^{4} The DSD method decouples the shock wave motion from the detonation product flow. This model employs detonation theory to calculate the motion of the detonation wave front, followed by numerical methods to determine the gaseous product flow behind the detonation wave front.^{10} The DSD model accurately describes the reaction growth from the shock wave and its transition to a detonation wave, providing high-precision solutions for the detonation wave structure. However, it is complex to simulate and requires significant computational resources.

In many cases, the explosive charge is considered to undergo ideal detonation, which can be simulated using a programmed burn model.^{9} The IDM simplifies the detonation process by ignoring details such as ignition initiation and detonation wave propagation within the explosive. It assumes that the detonation occurs instantaneously, forming a uniform gaseous product with high temperature and pressure, while the density of the detonation product remains the same as that of the unreacted explosive material. Baum *et al.*^{11} used the IDM to simulate the bomb explosion on the second floor of the World Trade Center in New York, solving the initial flow field that included hundreds of buildings and cars. Yu *et al.*^{12–14} investigated the dynamic behavior of cavitation induced by an underwater explosion near a free surface using an equivalent treatment of the detonation process. Daramizadeh and Ansari^{15} simulated the interaction between a shock wave and a free surface by simplifying the charge as an initial static high-pressure bubble. Jia *et al.*^{16} introduced self-similar variables into the symmetric Euler equations, obtaining state curves for the detonation product region under stable detonation of a one-dimensional planar TNT explosive through numerical solutions. These approaches based on the IDM can only qualitatively predict the state of the detonation zone due to their simplified assumptions and the use of an ideal gas equation of state for the detonation product. Generally, various methods define the initial state of detonation gaseous products in IDMs, but no universally accepted method has been proposed for predicting underwater explosions.

The DSD model demonstrates high precision because it accounts for the generation and propagation processes of the detonation wave within the explosive. However, it requires a very small grid size to accurately describe the flow in the reaction zone, leading to substantial computational demands and high costs, particularly for 2D and 3D problems. In contrast, the IDM simplifies the detonation process by ignoring specific details and treating the product region post-detonation as an equivalent initial physical domain. This approach offers high computational efficiency in engineering applications but at the cost of lower accuracy and limited generality, as it fails to accurately reflect the real flow structure in the detonation zone. To address these drawbacks, we have developed a high-precision IDM designed to provide foundational support for the numerical simulation of underwater explosions. This model operates using multicomponent fluid solvers without incorporating a detonation chemical reaction model, thereby achieving a balance between computational efficiency and accuracy. This new IDM aims to improve the representation of the detonation zone’s flow structure, enhancing the reliability of underwater explosion simulations in various engineering applications.

In this work, we first describe the numerical model for compressible multiphase flows using the DSD model. The complete system is solved using a simple fractional step approach based on finite difference schemes. The programmed burn model is integrated into the DSD model to capture the phase transition from unreacted energetic material to gaseous products and the release of detonation chemical energy. The Euler equation and level set equation are employed to control the flow of multicomponent fluids without detonation. Next, the DSD approach with the programmed burn model is used to reveal the detailed field structure in the detonation zone. From this, the high-precision IDM of the RS211 charge is derived through high-order approximation of the flow distribution characteristics. The IDM is suitable for various scenarios, including airblast and underwater explosions, as well as far-field and near-field explosions. This new model not only enhances the computational accuracy of simulations for underwater explosions and airblasts but also significantly reduces the required computational resources. Finally, we present a comparison between the IDM and the DSD model in underwater explosion simulations, demonstrating the effectiveness and efficiency of the IDM.

## II. NUMERICAL SIMULATION METHOD

^{17}

**U**denotes the conservative variables, F(

**U**) and G(

**U**) represent advective flux tensors, and S(

**U**) is the source term. The vector of conservative variables, advective flux tensor, and source term can be defined as

*ρ*,

*p*, and

*E*represent the density, pressure, and specific total energy of the fluid, respectively, and

*u*and

*v*are the velocity components in the coordinate directions

*x*and

*y*. The specific total energy is defined as

*E*=

*e*+ 0.5(

*u*

^{2}+

*v*

^{2}), where

*e*denotes the specific internal energy of the fluid.

*n*is the system parameter, which takes on a value of 0, 1, or 2. If

*n*= 0, it is for a 2D planar flow; if

*n*= 1, it is for a 2D axis-symmetric flow; if

*n*= 0 and

*v*= 0, it is for 1D plane flow; if

*n*= 1 and

*v*= 0, it is for 1D cylindrical flow; and if

*n*= 2 and

*v*= 0, it is for 1D spherical flow.

^{18}Therefore, the detonation wave is considered as a shock wave with constant velocity in the unreacted energetic material based on Chapman–Jouguet (CJ) theory. The programmed burn model requires that the burn factor

*λ*for the explosive cell should be found as

*t*

_{c}is the time when the detonation wave reaches the center of the element containing the explosive material,

*w*

_{z}is the equivalent size of the cell, and

*D*

_{CJ}is the detonation velocity. $V\u0304CJ$ is the relative density at the CJ state, $V\u0304CJ\u2009=\u2009\gamma /\gamma +1$, and $\gamma =1+DCJ\rho 0/2/e0$.

^{19}

*ϕ*is the distance function between the current cell and the interface. This model can be extended to more fluids by including additional advection equations $\varphi ix,y,i=1,2,\u2026$ of other interfaces for more fluids. When the sharp interface is captured by the level-set model, a modified ghost fluid method (MGFM) with strong robustness and consistency is used to handle the strong shock wave impacting on the multiphase interface. This method solves the two-shock approximation to the Riemann problem at the interface, which has been provided in the Appendix.

^{20}

^{16,21}The JWL EOS can be described as

*p*,

*e*, and $V\u0304$ are the pressure, internal energy, and relative volume ($(V\u0304=\rho /\rho 0)$, where

*p*

_{0}is the initial density of the explosive), respectively, and

*A*,

*B*,

*R*

_{1},

*R*

_{2}, and ω are constants for a particular explosive. The JWL EOS properties for TNT and RS211 explosives are presented in Table I.

Types . | $AGPa$ . | $BGPa$ . | R_{1}
. | R_{2}
. | ω
. | $e0MJ/kg$ . | $\rho 0kg/m3$ . | D_{CJ}
. |
---|---|---|---|---|---|---|---|---|

TNT^{16} | 371.2 | 3.21 | 4.15 | 0.95 | 0.3 | 4.29 | 1630 | 6930 |

RS211^{21} | 758.0 | 8.51 | 4.9 | 1.1 | 0.2 | 6.20 | 1750 | 7521 |

^{22–24}

*γ*and

*p*

_{∞}are parameters characteristic of the material and

*p*

_{∞}is a pressure-like parameter determined from empirical fit to experimental results. For perfect gases,

*γ*is the ratio of specific heats and

*p*

_{∞}= 0; for water, the selection of these constants is related to the ambient pressure, which have not been completely determined. The selection of different combinations of parameters

*γ*and

*p*

_{∞}represents different dynamic characteristics of fluids. For water, the parameters

*γ*= 7.42, and

*p*

_{∞}= 296 MPa are recommended.

^{18}The sound speed for the SG EOS can be expressed as

*c*

^{2}=

*γ*(

*p*+

*p*

_{∞})/

*ρ*.

For the propagation of detonation waves and shock waves with a high density ratio and high-pressure ratio, it is necessary to introduce high-order schemes to capture the strong discontinuities. The WENO reconstruction scheme proposed by Shu and Osher^{25} has been widely adopted in computational fluid dynamics to resolve compressible fluids that involve both shock waves and complex smooth structures. The numerical fluxes are first locally decomposed onto the respective characteristic fields. Then, they are projected back into physical space after high-order reconstruction.^{26,27} The time-marching of the Euler equations is adopted by a third-order TVD Runge–Kutta scheme.

## III. RESULTS AND DISCUSSION

### A. One-dimensional detonation wave validation

This one-dimensional TNT slab detonation process is often considered an excellent benchmark for testing numerical methods in the detonation simulation field. In this section, we examined the detonation of a 0.1 m TNT slab explosive initiated at one end. Liu previously simulated this scenario using a SPH method to demonstrate the propagation of the detonation wave within the explosive.^{28} In Liu’s simulation, a symmetric condition is applied at the initiation end, which is equivalent to detonating a 0.2 m long slab ignited at the center point. The properties of the TNT charge are provided in Table I. The computational domain size is [0, 0.1] m, with a mesh size of *δ* = 0.025 mm. Figure 1 shows a comparison of the solutions for density, pressure, velocity, and internal energy profiles from 2 to 12 *µ*s, with intervals of 2 *µ*s. The solid black lines in the figure represent the solution of the current model, while the dotted color lines denote Liu’s solution. As shown in Fig. 1, there is good agreement between these two numerical methods, validating the accuracy and reliability of the current model.

### B. Detonation of a RS211 explosive and its self-similarity characteristics

In this section, we examine 1D, 2D, and 3D symmetrical RS211 charges with a central initiation pattern under ideal underwater detonation conditions using plane, cylindrical, and spherical models, respectively, as shown in Eq. (1). The computational domain size is [0, 0.1] m with a mesh size of *δ* = 0.01 mm. Numerical results for density, velocity, pressure, and internal energy profiles for the 1D plane, cylindrical, and spherical RS211 charges are presented from 2 to 10 *µ*s at intervals of 2 *µ*s in Figs. 2–4. From each case, it is observed that, except for a slightly lower peak value at *t* = 2 *µ*s, the peak values of density, velocity, pressure, and internal energy remain consistent at the other four moments (*t* = 4, 6, 8, and 10 *µ*s). These peak values are very close to the CJ state values of the RS211 charge. The region behind the detonation wave front is mainly composed of two parts: the static zone and the Taylor rarefaction zone. Although the range of the static zone increases with the propagation of the wave front, the state values remain unchanged.

As the detonation waves propagate, the peak values at the wave front remain stable (close to the CJ state), and the ranges of both the Taylor rarefaction zone and the static zone increase correspondingly. By observing the shape of detonation waves at different times, a similarity is noted between them. Therefore, the dimensionless results of density profiles obtained by dividing the horizontal coordinates of the curves by the position of the wave front are presented in Fig. 5. It can be seen that the dimensionless curves at different times for the same case show good agreement with each other, indicating strong quantitative self-similarity. Similarly, the velocity, pressure, and internal energy profiles can be treated with the same approach, and self-similarity is also observed in these parameters. The quantitative calculation formula for the instantaneous detonation model can be derived by averaging or curve fitting the dimensionless results shown in Fig. 5. For underwater explosions and airblasts of spherical charges with arbitrary radii under central initiation conditions, the instantaneous detonation model can be directly mapped to the computation domain of the charge. This model does not require calculation from the initiation of the charge center, thereby saving significant computational resources.

### C. High precision equivalent IDM

In Sec. III B, it was observed that as the detonation waves propagate, the peak values at the wave front remain stable and very close to the CJ state. In addition, the range of the Taylor rarefaction zone and static zone increases correspondingly. To provide a clearer representation, the ordinate can be divided by the CJ state. The final dimensionless results of density, velocity, pressure, and internal energy profiles for the 1D plane, cylindrical, and spherical RS211 charges are shown in Fig. 6. A comprehensive formulation of the CJ state values for different charges can be obtained using the Rankine–Hugoniot jump relations. For the RS211 charge, the properties of the CJ state are as follows: *p* = 2277 kg/m^{3}, *u* = 1743 m/s, *p* = 22.9 GPa, and *e* = 7.7 MJ/kg. As depicted in Fig. 6, the demarcation point between the Taylor rarefaction zone and the static zone gradually shifts to the left for the three cases, and the value of the static zone decreases gradually. In addition, it is observed that the detonation zone of the spherical charge attenuates faster than that of the other two charges, which can be attributed to the size effect.

The Taylor rarefaction zone depicted in Fig. 6 can be accurately represented by high-order polynomials. The corresponding IDMs are presented in Tables II–IV for 1D plane, cylindrical, and spherical RS211 charges. Here, the static zone is described by a constant, while the Taylor rarefaction zone is fitted by fourth-order polynomials. To assess the accuracy of the fitting, a comparison between the accurate values and the fitting values of 1D plane, cylindrical, and spherical RS211 charges is illustrated in Fig. 7. Remarkably, the agreement between the two sets of results is excellent. It is worth noting that only two of the three variables for density, pressure, and internal energy in the detonation zone are independent, as the third can be determined by the JWL EOS of gaseous products. Although the peak values of density, velocity, pressure, and internal energy in the IDM are slightly smaller than those of the CJ state, this discrepancy has little effect on the overall simulation results in underwater explosion scenarios.

. | Taylor rarefaction zone $r\u0304\u22081,0.46$ . | |||||
---|---|---|---|---|---|---|

. | $fr\u0304=b0+b1r\u0304+b2r\u03042+b3r\u03043+b4r\u03044$ . | |||||

Detonation zone Physical quantity . | Static zone $r\u0304\u22080,0.46$ . | b_{0}
. | b_{1}
. | b_{2}
. | b_{3}
. | b_{4}
. |

Density | 0.690 | 1.433 | −4.778 | 10.28 | −8.660 | 2.737 |

Velocity | 0 | 2.291 | −14.64 | 31.20 | −25.88 | 8.067 |

Pressure | 0.298 | 1.053 | −4.543 | 8.680 | −6.009 | 1.852 |

Internal energy | 0.672 | 1.198 | −3.299 | 6.795 | 5.366 | 1.688 |

. | Taylor rarefaction zone $r\u0304\u22081,0.46$ . | |||||
---|---|---|---|---|---|---|

. | $fr\u0304=b0+b1r\u0304+b2r\u03042+b3r\u03043+b4r\u03044$ . | |||||

Detonation zone Physical quantity . | Static zone $r\u0304\u22080,0.46$ . | b_{0}
. | b_{1}
. | b_{2}
. | b_{3}
. | b_{4}
. |

Density | 0.690 | 1.433 | −4.778 | 10.28 | −8.660 | 2.737 |

Velocity | 0 | 2.291 | −14.64 | 31.20 | −25.88 | 8.067 |

Pressure | 0.298 | 1.053 | −4.543 | 8.680 | −6.009 | 1.852 |

Internal energy | 0.672 | 1.198 | −3.299 | 6.795 | 5.366 | 1.688 |

. | Taylor rarefaction zone $r\u0304\u22081,0.45$_{]}$fr\u0304=b0+b1r\u0304+b2r\u03042+b3r\u03043+b4r\u03044$
. | |||||
---|---|---|---|---|---|---|

Detonation zone Physical quantity . | Static zone $r\u0304\u22080,0.45$ . | b_{0}
. | b_{1}
. | b_{2}
. | b_{3}
. | b_{4}
. |

Density | 0.653 | 1.772 | −7.471 | 17.58 | −17.35 | 6.457 |

Velocity | 0 | 3.282 | −22.11 | 52.50 | −52.27 | 19.55 |

Pressure | 0.250 | 3.200 | −19.27 | 45.44 | −46.18 | 17.76 |

Internal energy | 0.642 | 1.907 | −8.378 | 19.83 | −20.01 | 7.630 |

. | Taylor rarefaction zone $r\u0304\u22081,0.45$_{]}$fr\u0304=b0+b1r\u0304+b2r\u03042+b3r\u03043+b4r\u03044$
. | |||||
---|---|---|---|---|---|---|

Detonation zone Physical quantity . | Static zone $r\u0304\u22080,0.45$ . | b_{0}
. | b_{1}
. | b_{2}
. | b_{3}
. | b_{4}
. |

Density | 0.653 | 1.772 | −7.471 | 17.58 | −17.35 | 6.457 |

Velocity | 0 | 3.282 | −22.11 | 52.50 | −52.27 | 19.55 |

Pressure | 0.250 | 3.200 | −19.27 | 45.44 | −46.18 | 17.76 |

Internal energy | 0.642 | 1.907 | −8.378 | 19.83 | −20.01 | 7.630 |

. | Taylor rarefaction zone $r\u0304\u22081,0.44fr\u0304=b0+b1r\u0304+b2r\u03042+b3r\u03043+b4r\u03044$ . | |||||
---|---|---|---|---|---|---|

Detonation zone Physical quantity . | Static zone $r\u0304\u22080,0.44$ . | b_{0}
. | b_{1}
. | b_{2}
. | b_{3}
. | b_{4}
. |

Density | 0.630 | 1.913 | −8.493 | 20.13 | −20.36 | 7.778 |

Velocity | 0 | 3.707 | −24.79 | 59.41 | −60.75 | 23.33 |

Pressure | 0.224 | 3.782 | −23.41 | 56.07 | −58.28 | 22.74 |

Internal energy | 0.626 | 2.218 | −10.47 | 24.92 | −25.60 | 9.896 |

. | Taylor rarefaction zone $r\u0304\u22081,0.44fr\u0304=b0+b1r\u0304+b2r\u03042+b3r\u03043+b4r\u03044$ . | |||||
---|---|---|---|---|---|---|

Detonation zone Physical quantity . | Static zone $r\u0304\u22080,0.44$ . | b_{0}
. | b_{1}
. | b_{2}
. | b_{3}
. | b_{4}
. |

Density | 0.630 | 1.913 | −8.493 | 20.13 | −20.36 | 7.778 |

Velocity | 0 | 3.707 | −24.79 | 59.41 | −60.75 | 23.33 |

Pressure | 0.224 | 3.782 | −23.41 | 56.07 | −58.28 | 22.74 |

Internal energy | 0.626 | 2.218 | −10.47 | 24.92 | −25.60 | 9.896 |

### D. Applications

The quantitative calculation formula of the IDM illustrated in Tables II–IV is utilized to predict the distribution of the detonation wave for 1D, 2D, and 3D RS211 charges, as depicted in Figs. 8–10. In Fig. 8, the prediction results of the 1D plane RS211 charge are compared to Liu’s results, as shown in Fig. 1. The position of the detonation wave front can be determined by the detonation velocity D_{CJ} under ideal detonation theory. Furthermore, the results of the 2D cylindrical RS211 charge are obtained using the IDM illustrated in Table III on the radial direction. Similarly, in Fig. 10, the distribution of density, velocity, pressure, and internal energy is conducted by the IDM illustrated in Table IV on the radial direction. Thus, the simulation of the detonation process for 1D plane, 2D cylindrical, and 3D spherical RS211 charges with center initiation under the ideal detonation hypothesis can be simplified. The final state of the detonation zone can be directly obtained from the IDM, providing a high-precision flow structure for the detonation zone.

### E. Comparison between the IDM and DSD models in underwater explosions

In the early stages of underwater explosion, two distinct phases are observed: the detonation wave stage and the subsequent shock wave propagation stage. To compare the performance of the IDM and DSD models, we focus on the detonation wave stage. For this comparative analysis, we examine the underwater explosion of a RS211 charge with a radius of 20 mm using both the DSD models and IDMs. The simulations are conducted on a uniform grid with a grid size of δ = 1 mm. To conserve computational resources, a 2D axisymmetric system is employed instead of a full 3D system. The computational domain size is set to [0, 0.6] × [−0.06, 0.06] m^{2}, with the y-axis serving as the symmetry axis. The explosive center is positioned at the origin of the coordinate system. The surrounding water is simulated at atmospheric pressure to replicate a shallow water explosion environment. The interface between the charge and water is accurately captured using a combination of the level-set method and MGFM treatment.

According to the CJ theory, the detonation wave is expected to reach the surface of the spherical charge at ∼2.66 *µ*s after central initiation. Consequently, the distribution results obtained from the IDM are derived through direct mapping. This setup allows for a direct comparison between the IDM and DSD models, providing insights into their respective capabilities and performance in accurately simulating the detonation wave stage of underwater explosions.

#### 1. Detonation wave stage

The density and pressure profiles obtained from the IDM at the detonation stage are compared with those from the DSD model solution from 0.6 to 2.4 *µ*s, with a time interval of 0.6 *µ*s. This comparison is illustrated in Figs. 11 and 12. It is observed that the IDM exhibits higher accuracy in capturing the peak density and pressure of the detonation wave front compared to the DSD model. This enhanced accuracy can be attributed to the fact that the IDM utilizes solutions provided by the DSD model, albeit on a finer grid with a grid size of *δ* = 0.01 mm.

During the mapping process from the radial grid to the spatial structured grid, a zigzag distribution pattern may emerge in the detonation wave front due to size mismatches. However, this zigzag shape tends to disappear shortly after the detonation wave reaches the surface of the charge, and its influence on the overall calculation accuracy can be disregarded. Moreover, the pressure histories at two points located at different distances from the center of the charge are compared between the IDM and DSD models, as shown in Fig. 13. It is evident that the IDM offers higher accuracy in predicting the pressure history compared to the DSD model solutions.

#### 2. Shock wave propagation stage

Numerical results of density and pressure profiles from 60 to 240 *µ*s under the two models are depicted in Figs. 14 and 15. Once the detonation wave reaches the charge surface at 2.6 *µ*s, the detonation process concludes, and a high-pressure, high-velocity gaseous product bubble forms in the water, causing a strong shock wave propagation. This interaction triggers a rarefaction wave that propagates into the bubble, resulting in a noticeable low-pressure zone inside the explosion bubble. Comparatively, there is little difference in the calculation accuracy between the solutions of the two models.

The pressure history at four points with different distance-to-radius (r/R_{0}) ratios is illustrated in Fig. 16. It is observed that the peak pressure of the IDM is slightly lower than that of the DSD model, with the maximum relative error between the two models being less than 5%. However, there is no noticeable difference in the pressure histories between the two models. Generally, the accuracy of the IDM is equivalent to that of the DSD model under the condition of the same grid size. Numerical results demonstrate that the IDM proposed in this paper exhibits good accuracy and high computational efficiency.

## IV. CONCLUSIONS

In this article, we propose an efficient and highly precise IDM for an RS211 charge to simulate the detonation process using a direct mapping approach. This model yields accurate solutions for the fluid structure within the detonation zone of the charge. The IDM is derived from the high-order solutions provided by the detonation shock dynamics model. In the DSD model, the spatial term is discretized using fifth-order WENO reconstruction in characteristic space and Lax–Friedrichs’s splitting, while the temporal terms are discretized using a third-order TVD Runge–Kutta scheme. In addition, interface motion is captured using a level-set method combined with the multiphase Godunov front-tracking method (MGFM), and a new programmed burn model is introduced to describe the generation and propagation of the detonation wave.

The self-similarity of detonation wave propagation is validated using the DSD model, and the quantitative calculation formula of the instantaneous detonation model is derived by averaging or curve fitting the dimensionless results. The IDM is established by approximating the Taylor rarefaction zone and constant static zone using high-order polynomials for 1D plane, cylindrical, and spherical spaces.

The application of the IDM involves direct mapping from the radial direction to the spatial structured grid for 1D plane, 2D cylindrical, and 3D spherical charges. A comparison between the IDM and the DSD model in underwater explosion simulation reveals that the accuracy of the IDM is equivalent to that of the DSD model under the condition of the same grid size. Numerical results demonstrate that the IDM proposed in this paper exhibits good accuracy and high computational efficiency. Furthermore, the approach adopted in this study can be easily extended to other types of charges.

## AUTHOR DECLARATIONS

### Conflict of Interest

The authors have no conflicts to disclose.

### Author Contributions

**Yi Hao (郝轶)**: Conceptualization (equal); Formal analysis (equal); Investigation (equal); Resources (equal); Validation (equal); Writing – original draft (equal). **Jun Wang (汪俊)**: Data curation (equal); Investigation (equal); Methodology (equal). **Jiping Chen (陈继平)**: Data curation (equal); Formal analysis (equal); Resources (equal); Software (equal). **Zhenxin Sheng (盛振新)**: Data curation (equal); Investigation (equal); Methodology (equal). **Guozhen Liu (刘国振)**: Formal analysis (equal); Software (equal); Validation (equal); Visualization (equal). **Jun Yu (余俊)**: Conceptualization (equal); Formal analysis (equal); Methodology (equal); Resources (equal); Writing – review & editing (equal).

## DATA AVAILABILITY

The data that support the finding of this study are available from the corresponding author upon reasonable request.

### APPENDIX: MODIFIED GHOST FLUID METHOD (MGFM)

*et al.*

^{20}to enhance the treatment of two-phase fluids involving strong shocks at multi-material interfaces. This method offers a consistent solution at the interface by applying a two-shock approximation to the underlying Riemann problem. The solution of the two-shock Riemann problem entails formulating two nonlinear scalar equations based on pressure and velocity balance across the contact discontinuity,

*L*and

*R*represent the fluid at the left and right of discontinuity and the variables $\rho *L,\rho *R,u\u2217,p\u2217$ denote the unknown interfacial density, velocity, and pressure, respectively. The functions

*V*

_{L},

*V*

_{R},

*P*

_{L}, and

*P*

_{R}are obtained from the Rankine–Hugoniot relations when the two states are separated by two-shock waves.

*h*denotes the specific enthalpy,

*h*=

*e*+

*p*/

*ρ*, while the subscript

*k*represents the left or right state, and the subscript * indicates the state behind the shock front.

*p*

_{∗}and the density

*ρ*

_{*k}can be obtained. For SG EOS, these relations have algebraic expressions and can be defined as

*f*and

*g*denote the Mie–Grüneisen function. When the relationship between pressure and density has been confirmed by solving Eq. (A5), Eqs. (A3) and (A4) can be simplified to a nonlinear equation containing only interfacial pressure

*p*

_{∗}, which can be solved by an iterative algorithm (such as the Newton iterative method). Therefore, the Riemann solutions on both sides of the discontinuity are completely obtained, and they can be extended to the WENO reconstruction template, thereby realizing two-phase flow interface treatment based on MGFM.

## REFERENCES

*Underwater Explosion*

*Numerical Simulation Method of Explosion Flow Field Including Moving Interface and its Application*

*Smoothed Particle Hydrodynamics: A Meshfree Particle Method*