#### Abstract

In order to obtain the shockwave load simplified algorithm model for the semiarmored projectile internal explosion in the cabin, this research made use of AUTODYN to provide a numerical modeling method for explosion in the cabin and verified the accuracy of the method via the experiment. Internal explosion simulation calculation was conducted on the operating condition numerical model with different cabin structural dimensions and different explosive loads. The cabin internal explosion space was divided into the noncorner central area, near-wall area, two-sided corner area, and three-sided corner area. Through regression of the abovementioned calculation results, an engineering model to calculate the shockwave load was obtained. It is hoped that the model can offer some references to the antiexplosion design for the ship cabin and for damage assessment of the internal explosion.

#### 1. Introduction

An explosion happening in an enclosed or semienclosed space under limited boundary conditions is defined as an internal explosion. Different from the explosion happening on the free field, the high-temperature, high-pressure products cannot be diffused externally and immediately [1–3]. The shock wave might go back and forth between the structural walls. Due to the superposition of the reflected waves, peak value of overpressure increases, coupling effect of the explosive shock wave in the corner area becomes increasingly obvious, and the time of reaction is lengthened, thus causing more damage to personnel, equipment, and so on, within the structure [4–6]. Recent years have witnessed vigorous research attempts made by experts and scholars from China and abroad to examine the damage caused by an explosive shockwave load in an enclosed space to the structure. Penneiter et al. [7] studied different reflection forms of shock waves caused by explosive explosion in an enclosed space and at typical position. The simulation calculation results showed a good agreement with the experiment results, and they also described distribution characteristics of the wall surface shock wave. Kurki et al. [8] carried out simulation calculation of the damage load caused by the semiarmored projectile internal explosion to obtain the damage mechanism of the damage load on the cabin structure. Yang et al. [9] obtained the overpressure time travel curve at typical measurement points inside an enclosed cabin, which generally coincided with experimental results. Their simulation results suggested that the overpressure at the cabin corner was several times higher than that at the bulkhead center. Geretto [10] studied explosion of certain warm-pressing explosive and calculated the shock wave of the warm-pressing explosive in AUTODYN, finding that the shock wave of the warm-pressing explosive was around twice as much as that of TNT explosive when the explosion happens on the free field and around 13 times as much as that of the latter when the explosion happens in the underground tunnel. This meant that the warm-pressing explosive was more suitable to create the internal explosion damage effect. Feldgun et al. [11] conducted a series of researches on internal explosive explosion in semienclosed spaces, such as structures and cabins. They studied the pressure curve of the structural wall surface through numerical simulation and analyzed its distribution rules and evolution process. Meanwhile, they conducted an explosive explosion experiment in a square enclosed space to study the evolution rules of the quasistatic pressure, compare them with the results given by the classical empirical formula, and put forward factors influencing the quasistatic pressure evolution, which laid a solid foundation for follow-up research. Peng [12] studied the dynamic response model of stiffened plates under impact load on the premise of considering the strain rate effect of the material itself. Jacob et al. [13] summarized the failure mode of the plate under explosive loading and pointed out that, in addition to the plastic deformation and tensile damage modes, there may also be local damage modes such as smashing discs and petal-shaped ruptures. Jiao et al. [14] deduced the dimensionless destruction number of multicompartment structure and gave a rapid damage prediction method of multicompartment structure. Hou et al. [15, 16] used MSC/DYTRAN simulation software to analyze the shockwave overpressure load characteristics when explosion happens in the cabin scale model. Results indicated that the superimposed effect of the corner could further increase the peak value of shockwave load overpressure, but experimental data could not provide solid evidences for this viewpoint. Shockwave data obtained by experiments showed that the largest and smallest shock waves were found at the measurement point at the two-sided corner and the three-sided corner, respectively. Li [17] analyzed convergence of the shock wave at the corner and preliminarily identified the high-pressure area of the corner through a series of simulation calculation. According to their finding, the cabin length-width ratio was a main factor that affected the high-pressure scope of the corner, and their research also casted light on the functional correlation between the high-pressure scope and the structural dimensions.

It is not difficult to see that the above research is only limited to the qualitative analysis of the internal explosion load propagation, and the understanding of the internal explosion mechanism is not perfect. Although the superimposed enhancement effect of the internal explosion load at the corners of the structure is found, the internal explosion shockwave propagation is lacking. The superimposed theoretical calculation model cannot quickly and accurately predict the structural failure mode caused by the internal explosion. It is difficult to meet the engineering needs for the damage assessment of the ship in the cabin explosion.

This research combined theoretical analysis, numerical simulation, and experimental verification to study a simplified algorithm model for the explosive shock wave when an explosion occurred inside the cabin. An explosion experiment in a simulation cabin was conducted, and the experimental data thus obtained were used to calibrate the numerical simulation model. Meanwhile, AUTODYN was used to obtain the data of the explosive explosion happening within the cabin model, based on which a numerical algorithm model was built. Next, the working conditions were simulated in accordance with cabin models of different sizes and explosives of different masses. The interior of the cabin was divided into zones in the light of the quantity of bulkheads whose shock wave was influenced by reflection. Through simulation calculation, shockwave load parameters and their distribution rules in different zones were obtained. On that basis, an engineering model which could be used to work out the shockwave load was established, which could provide data support and basis for the warhead design and the cabin target protective design.

#### 2. Explosion Experiment in the Cabin

##### 2.1. Experiment Goal and Plan

In order to study the distributional rules of explosion shockwave load in the cabin, a simulation cabin was built, with its internal dimensions being 900 mm × 600 mm × 600 mm. This model set up four measurement sites, including the bulkhead center, near-end two-sided corner, far-end two-sided corner, and three-sided corner, to install the high-frequency piezoelectric sensor. The pressure sensor frock geometric center at the corner was located 50 mm away from the bulkhead. Upon installation, the sensor manometry sensitive area and the cabin internal wall were ensured to be paralleled to each other. The spherical charge was lifted on the cabin top (a circular hole whose diameter is 70 mm) and located at the geometric center of the space. In order to measure the positive and negative reflection shockwave overpressure, the sensor P1 was installed at the planimetric center of the bulkhead that was close to the explosion center. The sensor P2 was installed at the center of the near-wall two-sided corner ridge. The sensor P3 was installed at the center of the far-end two-sided corner ridge. The sensor P4 was installed at the three-sided corner. The explosive for the experimental use was the TNT spherical bare charge. The explosive loads for the four measurement sites were 15 g, 20 g, 25 g, and 30 g, respectively. Each explosive load was experimented twice. The average of the valid data from the experiment under each explosive load was used to simulate the cabin sensor layout, whose schematic diagram is shown in Figure 1.

**(a)**

**(b)**

**(c)**

**(d)**

##### 2.2. Experimental Devices

Q235 steel whose wall thickness was 16 mm was adopted for the simulation cabin bulkhead. The explosive for experimental use was TNT bare charge. The data were collected using the nuxi-1008 high-speed data collector produced by Sichuan Tuopu Measurement and Control Technology. The high-frequency piezoelectric sensor was connected with the data collector. The sampling frequency of the shockwave load signal channel was 1 MHz, and its trigger mode featured external trigger. The CY-YD-205 sensor manufactured by Jiangsu Sinocera Piezotronics Inc., was adopted, with its range being 0∼10 MPa, output voltage being 0∼5 V, working temperature being −40∼150 °C, and response frequency being 100 kHz. Being highly elastic, insulative, and efficient in obstructing stress waves, nylon was used to design the nylon sleeve tooling for the sensor so as to reduce the experimental error caused by the interfering signals in the experiment process. The pressure-sensitive surface of the sensor and the inside of the bulkhead were maintained parallel to each other. The nylon tooling and the bulkhead interference were installed together compactly through interference fit [18]. The high-frequency piezoelectric sensor and nylon tooling for experimental use are shown in Figure 2. The layout of the simulation cabin and sensor is presented in Figure 3.

##### 2.3. Experimental Results and Analysis

The typical shockwave load experimental data of the 25 g TNT explosive load are displayed in Figure 4. The measurement point P1 was located at the center of the cabin wall. The shockwave generated by the explosion of the explosive in the cabin center first arrived at the geometric center, the bulkhead, and caused positive reflection. The measuring point P2 was located at the center of the near-end two-sided corner ridge. The measuring point P3 was located at the center of the far-end cabin two-sided corner ridge. According to the data curve obtained by the experiment, obvious macroscopic pulsation phenomenon of shockwaves appeared at the cabin corner. The shockwave overpressure peak value of the measuring points P3 and P2 obtained by the sensor was widely apart by around twice, which suggested a significant convergence effect of the cabin side ratio on the shockwave at the two-sided corner. The measuring point P4 was located at the three-sided corner of the cabin. The overpressure time travel line demonstrated obvious multipeak characteristics. Therefore, when studying the damage effect of the shockwave load of explosion inside the cabin on the cabin, one should not only consider the overpressure peak value of the pulse at the cabin corner but also pay attention to the follow-up damage effect of the pulse. The superimposed effect at the corner could further strengthen the capability of the shockwave to damage. Therefore, the corner structure should be a focus of the antiexplosion design for the cabin.

**(a)**

**(b)**

**(c)**

**(d)**

#### 3. Establishment and Verification of Numerical Simulation Model

##### 3.1. Settings of Model Parameters

The simulation model was established in this paper by the following numerical simulation approach. Multimaterial Euler algorithm was adopted to calculate the air and explosive parameters. The model’s axis of symmetry was defined as a symmetric boundary. The flow-out boundary conditions were added to the air boundary. The simulation model resorted to explosion at the center. The three-dimensional symmetric one-eighth model was adopted. The explosive was detonated at the geometric center of the cabin. The explosive load was set to be 25 g for numerical simulation calculation. The grid size was 2 mm. The multimaterial Euler algorithm was adopted. The materials were directly selected from the AUTODYN material bank. The air parameters were worked out by the ideal-gas state equation. Table 1 demonstrates the material parameters:

The experimental explosive was the spherical TNT bare charge, whose parameters are presented in Table 2. TNT is filled in the air domain. Euler/Lagrange flow-solid interaction was adopted for the coupling between the air and the cabin wall surface.where *E* is the internal energy per unit mass, *V* is the relative specific volume of the explosive red product, and *A*, *B*, *R*_{1}, *R*_{2}, and are the constants and the equation of state and parameter values and TNT explosive-related parameters: density *ρ*_{E}, detonation pressure *P*_{CJ}, and detonation speed *D*.

Q235 steel was adopted as the bulkhead material, whose wall thickness was 16 mm; the cabin numerical model was maintained the same to that of the experiment, with its internal space being 0.9 m × 0.6 m × 0.6 m; the chosen strength model was Johnson-Cook. The material parameters are shown in Table 3.

##### 3.2. Verification of Numerical Simulation Model

In order to verify the validity of the numerical simulation algorithm model, this paper conducted simulation calculation of the experimental working conditions. The relative error between the data acquired by the experiment and by the simulation calculation was used to verify the accuracy of the numerical model established by this paper. Figure 5 demonstrates the propagation process of shockwave formed by the explosive explosion in the cabin. The red area was the outmost wave front of the shockwave. Through simulation, the evolution process of the explosion shockwave inside the cabin is shown as follows.

**(a)**

**(b)**

**(c)**

**(d)**

After detonation of the explosive, the explosion shockwave, before reaching the internal wall of the device, first underwent three-dimensional axis-symmetric inflation. The outmost part of the red area is the wave front of the shockwave. Because of the same distance between the explosive load and the lower and upper wall, the shockwave inflation was symmetric. The pressure cloud chart of the wall on two sides is completely the same. As shown in Figure 5(a), the explosion shockwave reached the wall the first time at 0.15 ms, and the positive reflection occurred there. After that, the remaining points of the bulkhead started receiving the shockwave load, which underwent the regular oblique reflection and Mach reflection in turn. As shown in Figure 5, with the passage of time, shockwave reflected multiple times in an enclosed space, because of the restriction of the wall in four directions during the propagation process. In Figure 5(b), obvious Mach wave was formed because of the irregular reflection of the wall at 0.25 ms. The initial shockwave and the Mach wave from the wall converged, and the action time between the two lasted long. One can observe from Figure 5(c) that the explosion shockwave generated the local pressure convergence at the interaction between two walls at 0.45 ms and that the shockwave kept on propagating along the ridge. This can also explain why overpressure soared up. In Figure 5(d), when the shockwave arrived at the ridge destination at 0.65 ms, the convergence pressure could be formed again at the three-wall corner area. The shockwave overpressure curve is characterized by multiple peaks and macroscopic pulsation. Numerical simulation results could well reflect the practical situations of this experiment. The experimental data of the bulkhead positive reflection were compared with simulation results for calibration and verification of the validity of the numerical simulation model. The pressure-time curve of the bulkhead positive reflection measuring point at 7 ms and with 25 g explosion load was compared with the experimental data, as shown in Figure 6.

Numerical simulation calculation results can directly show the phenomenon of multiple pulsation of pressure load. The shockwave load peak overpressure of the operating condition curve was 2.62 MPa. The shockwave load peak overpressure was read and taken down by this approach and comparatively analyzed with experimental data. The analysis results are presented in Table 4.

Table 4 shows that the deviation between the calculated value and the experimental value ranges from −7.487%∼8.177%, which suggests favorable applicability of the numerical simulation approach to calculate the explosion shockwave load in the cabin.

#### 4. Simplified Algorithm Model for Shockwave Load

##### 4.1. Theoretical Analysis of Shockwave Convergence at the Corner

The overpressure propagation rules of the shockwave generated by the explosive after internal explosion in an enclosed space are extremely complex. Restricted by the cabin wall, the shockwave will undergo reflection, superposition, convergence, and long-lasting quasistatic pressure many times, which can cause serious damage to personnel and devices within the cabin. Propagation of the shockwave in the cabin adheres to the wall surface reflection principle, which includes the positive reflection and oblique reflection. The oblique reflection can be further divided into the regular reflection and Mach reflection. When the angle of incidence is smaller than the critical angle, *α*_{e} of the Mach reflection, the shockwave reflection will be the regular oblique reflection. On the contrary, the shockwave reflection will be the Mach reflection. The critical angle, *α*_{e}, of the Mach reflection depends on the intensity of the incident wave. The theoretical calculation formulas of various reflection behaviors are shown as follows:where denotes the incident overpressure; denotes the reflection overpressure; denotes the atmospheric pressure; denotes the incident angle; denotes the math reflection peak overpressure; denotes the peak overpressure formed by TNT of the equivalent weight under rigid wall explosion; and denotes the corresponding incident angle.

As shown in Figure 7, the explosive was located at the center of the cabin. Half of the length of the longer and shorter side on the horizontal plane of the cabin was set to be *R* and *H*, respectively. Then, the side ratio of the cabin can be written as *λ* = *R*/*H*. When the incident angle, *α*, is larger than the Mach reflection critical angle, *α*_{e} (usually being *α*_{e} = 40°), the shockwave reflection will be the Mach reflection. In other words, when *λ* = *R*/*H* > tan*α*_{e}, the initial shockwave formed by the explosive explosion will undergo the Mach reflection on the wall surface A. Similarly, when *H*/*R* > tan*α*_{e}, the initial shockwave formed by explosive explosion on the wall surface B will undergo the math reflection too. So, when 1/tan*α*_{e} > *λ* > tan*α*, the initial shockwave will undergo the Mach reflection on the wall surfaces A and B. The Mach wave is pushed to the corner *M* along the wall surface horizontally to have superposition and convergence with the initial incident wave, thus forming the local high pressure. When *λ* > 1/tan*α*_{e,} the initial shockwave formed by the explosive explosion will experience the Mach experience on the wall surface A, but the regular oblique reflection on the wall surface B. Since the reflected wave propagates in the air with a rising temperature, after being compressed by the incident wave, the wave front speed of the reflected wave is faster than the propagation speed of the incident wave front. The Mach wave front becomes wider and wider in the propagation process to form an obvious Mach stem. The Mach stem progresses along the bulkhead with the longer side. As the value of *λ* keeps on increasing, the Mach wave front arrives at the corner before the incident wave. The incident shockwave directly propagates to the wall surface B of the cabin along the triple-point wave track without converging at the corner.

##### 4.2. Division of the Cabin Space

Under the effect of internal explosion, the shockwave experiences reflection and supervision constantly in the cabin, thus resulting in an uneven distribution of the shockwave load distribution in the cabin and convergence of the shockwave at the corner of the structure.

Therefore, it is necessary to divide the cabin space before the calculation of the shockwave overpressure peak in the cabin. Different methods were adopted to work out the explosion shockwave load in different areas. By comparing data obtained by the experiment and simulation, this paper verified the accuracy of the numerical simulation approach. Under the prerequisite of not changing various parameters of numerical simulation, the numerical simulation approach was used to examine the internal explosion shockwave propagation process and study rules governing the impact of different cabin structural dimensions on the distributional characteristics of the shockwave load at the corner. The statistical analysis results of the numerical simulation results suggest that the area division under different explosive loads was generally the same [19]. The simulation model maintained the square with the width and height both being 0.6 m and designed by the experiment to simulate the cabin unchanged. The TNT spherical bare charge was 100 g. The dimensions of the longer side were changed to study the distribution rules and convergence effect of the shockwave load distribution. Investigation of the side ratio of the ship’s typical cabin suggested that the range of the side ratio, *λ*, 1 < *λ* < 3, was universally applicable. Therefore, 12 cabin simulation operating conditions with different side ratios were designed within the above range. One-eighth model of the cabin was designed, and measuring points were evenly distributed on the central plane of the cabin crossing the explosion point. Since the scaled distance of the measuring point which was within 0.1 m away from the explosion center was defined as the initial position of the observation point, the pressure peak value of various measuring points was extracted to draw the pressure contour line. The pressure distribution cloud chart of the cabin internal explosion central plane with different side ratios under typical operating conditions is shown in Figure 8.

The short side length of the cabin two-dimensional plane was *H*. According to the cabin central plane overpressure distribution chart, the position where the wall surface reflection pressure of the incident wave is equal to the shockwave overpressure peak of the measuring point at the corner is defined as the demarcation boundary. The boundary 1 was set to be 1 and away from the bulkhead with the longer and shorter side, respectively. Through translation of the boundary between the cabin horizontal plane and vertical plane, the internal space of the cabin could be divided into four areas. First is the central area. Similar to the free field propagation, the peak value of the shockwave overpressure time travel curve appears at the pulse arriving first. Second is the near-wall area. The shockwave experiences reflection at the cabin wall. The peak of the overpressure time travel curve appears at the follow-up pulse, and the peak keeps increasing. Third is the two-sided corner area. This area is subject to the influence of the two-sided bulkhead reflection, which can reflect in the superposition effect of the shockwave. Fourth is the three-sided corner area. As to the area subject to the influence of the three-sided bulkhead reflection, the overpressure peak of the corner position will be further increased because of the convergence effect of the shockwave. Under 12 simulation operating conditions, the shockwave distribution rules at the cabin plane corner were analyzed. When the cabin side ratio was 1, the high-pressure area scope of the corner was the same to the bulkhead length on two sides of the cabin. Take the cabin side ratio, *λ*, as the *x*-coordinate, and the ratio of l/*H* to /*H* as the *y*-coordinate. The correlation between the cabin corner high-pressure area scope and the side ratio can be expressed as that in Figure 9. Polynomial fitting was conducted of data points to obtain the functional relation between the high-pressure area scope at the cabin two-sided corner and the side ratio of the structure.

After polynomial fitting of the function correlation between l/*H* and side ratio, *λ* (1 <*λ* < 3), the following equation can be obtained; after polynomial fitting of the function correlation between /*H* and side ratio, *λ*, the following equation can be obtained:

##### 4.3. Establishment of Simplified Algorithm Model for Shockwave Load

The experimental data and the simulation data were compared to verify the accuracy of the numerical simulation approach. Without changing the parameters of numerical simulation, the shockwave load parameter algorithm model for different areas in the cabin space was built. Five typical ship cabin models of different sizes but of the same shape were designed. The reference points were evenly laid out on the vertical and horizontal planes at the center of the cabin. The bulkhead thickness was both 16 mm. In each cabin, explosives of five different masses were denoted, respectively, and simulation calculation was conducted on 25 operating conditions. The numerical simulation operating conditions are presented in Table 5, and the numerical calculation models in Figure 10.

According to the corner boundary identification approach, the simplified algorithm model of shockwave overpressure peak value in different areas of the cabin shockwave space was studied. Before the establishment of the shockwave peak overpressure, *p*_{r}, algorithm model, it was necessary to identify the model form. Based on the similarity theory, the dimensional analysis was conducted, suggesting that the model form featured a function of the scaled distance, Z. It can then be expanded into a polynomial form, where the calculation model of shockwave load parameters is

After the algorithm model form was pinpointed, the numerical simulation results were combined for a fitting of shockwave parameters at different points. The fitting results of the central area data are shown in Figure 11.

Bring the proposed determined coefficient into equation (6), and then derive the shockwave load calculation model in the central area as follows:

According to the above fitting method, the shock wave load parameters of each working condition observation point are calculated, and fitting of the load peak overpressure in different reflection areas is conducted, thus obtaining the shockwave load model of each area:where ; TNT denotes the explosive; *Q* denotes the quantity, kg; *R* denotes the explosion distance, *m*; and *p*_{r} denotes the peak overpressure, MPa.

#### 5. Conclusions

This research discussed the convergence effect of the shockwave load at the corner under the cabin internal explosion, revealed the formation mechanism of the shockwave high-pressure area at the corner of the ship cabin, and analyzed rules influencing the convergence effect of the cabin structural dimensions on the shockwave load at the corner. Based on the above discussions and analysis, this research came to the following conclusions:(1)If the explosive exploded in the cabin, the shockwave load at the cabin corner would generate an obvious convergence effect. The overpressure time travel curve featured multiple peaks. The overpressure peak of the pulse arriving first was not the maximum. To study the damage effect of the internal explosion shockwave on the corner requires taking the effect of the follow-up pulse into consideration.(2)As the cabin side ratio kept on increasing, the high-pressure area at the shorter side of the corner gradually expanded as well. The high-pressure area at the longer side of the corner first narrowed and then expanded. When the side ratio was 2, the convergence effect of the two-sided corner on the shockwave reached minimum. This finding can provide some references for antiexplosion design for the cabin.(3)This research provided the way to identify the boundary of various areas in the cabin. Based on numerical calculation of cabin models of different dimensions and under different operating conditions with different explosive loads, the simplified algorithm model for the peak overpressure load of different areas in the cabin was built.

#### Data Availability

The datasets supporting the conclusions of this article are available from the corresponding authors upon request.

#### Conflicts of Interest

The authors declare that they have no conflicts of interest to report regarding the present study.