Abstract
A numerical simulation has been carried out to examine the response of steel plates with different arrangement of stiffeners and subjected to noncontact underwater explosion (UNDEX) with different shock loads. Numerical analysis of the underwater explosion phenomena is implemented in the nonlinear finite element code ABAQUS/Explicit. The aim of this work is to enhance the dynamic response to resist UNDEX. Special emphasis is focused on the evolution of mid-point displacements. Further investigations have been performed to study the effects of including material damping and the rate-dependant material properties at different shock loads. The results indicate that stiffeners configurations and shock loads affect greatly the overall performance of steel plates and sensitive to the materials data. Also, the numerical results can be used to obtain design guidelines of floating structures to enhance resistance of underwater shock damage, since explosive tests are costly and dangerous.
1. Introduction
The dynamic responses of floating and submerged structures subjected to noncontact underwater explosion have received attention as a concerned topic for a national defense since the 1950s. When a naval ship is attacked by an underwater explosion (UNDEX), the ship can be severely damaged by shock waves and gas bubble pulse. Predicting the shock response of ships to noncontact UNDEX from underwater weapon is of great importance for the warhead design of underwater weapon.
The design and analysis of structure subjected to UNDEX require a detailed understanding of explosion phenomena and the dynamic response of various structural elements [1]. UNDEX structural damage depends on standoff distance and explosion weight. The most important way to reduce the damages due to UNDEX is to provide sufficient standoff distance between the structure and explosion source in order to decrease the effect of the UNDEX wave so that the structure is not highly damaged. To accomplish these objectives, it is necessary to do various scenarios to evaluate the behavior of the ship structure to blast loading. These scenarios should include studying such aspects (explosive magnitude, distance from source of explosion, structure scantling, complex fluid-structure interaction phenomena, structure geometry, etc.) [2].
Rajendran [3] carried a numerical simulation of underwater explosion experiments on plane plates using ANSYS/LS-DYNA for elastic and inelastic response of circular and rectangular plates. Gupta et al. [4] predicted the failure modes of stiffened and unstiffened square plates subjected to underwater explosion using an elastoplastic model with isotropic hardening, strain rate effects, and fracture criterion. Qiankun and Gangyi [5] predict the shock response of a ship section to noncontact UNDEX using ABAQUS. The results showed that the size of fluid mesh and the fluid thickness from the wetted surface to the outer boundary of fluid is of great importance for improving numerical accuracy. Rajendran and Narasimhan [6] studied the damage prediction of clamped circular plates subjected to contact underwater explosion. Jacinto et al. [7] studied a linear dynamic analysis of the plate models under explosive loading using ABAQUS.
Zhang et al. [8] developed a procedure which employs the finite element method coupled with the DAA method to study the transient dynamic response of a surface ship model subjected to an underwater explosion bubble. And they found that the ship model’s vertical response subjected to an underwater explosion bubble is mainly the global response, which is composed of two parts: rigid-body motion and elastic deformation.
When a charge of high explosive is detonated below the surface of water, the sudden release of energy from underwater explosions generates a shockwave propagate in a spherical pattern at very high speed. The explosive also forms highly compressed gas bubble in the surrounding water, which are very long compared with the times for shock passage. Of the total energy released from a 1500-lb TNT underwater explosion approximately 53% goes into the shockwave and 47% goes into bubble pulsation. Most cases demonstrate that the damage done to marine structures, such as the surface of ships and submarines, occurs early and is due to the primary shockwaves. This investigation only considers the effects of the shockwaves as demonstrated in [9].
Zhang et al. [10] calculated the dynamic bending moment of bubble acting on hulls based on the underwater explosion bubble theory. And he found that underwater explosion bubble will cause general longitudinal bending of ship and lead to hogging and sagging damage of ship.
Forghani et al. [11] studied the modeling of damage development in blast loaded composite laminates to blast loads. It was shown that the tie-break interface option in LS-DYNA can be used successfully in simulating cohesive cracks. Chirica and Boazu [2] studied the response of ship hull laminated plates to blast loads using finite-element computer code COSMOS/M. Various scenarios to evaluate the behavior of the ship structure laminated plate to blast loading are studied like explosive magnitude, distance from source of explosion, and plate thickness.
Jen [12] used finite element software ABAQUS to simulate and analyze the transient dynamic response of a midget submersible vehicle (MSV) pressure hull that experiences loading by an acoustic pressure shock wave resulting from an underwater explosion (UNDEX). The results showed that, the local cavitations formed on the fluid–structure interface, after the shockwave hit the midget submersible vehicle. The aim of this work is to enhance the dynamic response to resist UNDEX and obtain the optimal configuration for stiffened plates to resist underwater shock loading
2. The UNDEX Loading
Noncontact underwater explosion is the major source of threat to ships and submarines. The responses and damages of submerged structures are divided into two categories: near-field explosion and far-field explosion depending on the distance between the explosive charge and the target (standoff distance) [12].
Figure 1 shows the different events occurring during the UNDEX event in a pressure against time history plot [13]. The loading mechanisms due to underwater explosion include incident shock wave, free surface reflection wave, bottom reflection wave, gas bubble oscillation, bubble-pulse loading, and bulk and hull cavitations [14]. The reflection of the shock wave from the bottom of the ocean is a compression wave that adds additional load to the structure. The reflection of the shock wave from free ocean surface causes a reduction in the pressure produced by the shock wave [11, 15]. The incident wave is the shock wave produced by the UNDEX charge. The scattered wave is the acoustic field generated by the interaction of the incident wave and the submerged structure. The initial shock wave is modeled as a spherical incident shock wave applied as a transient load active on both the acoustic and structural meshes at their common surfaces (the wetted interface). The distribution of this shock wave onto the plate is obtained using the incident pressure wave equations as demonstrated in [16]. The pressure history of the shockwave at a fixed location starts with an instantaneous pressure peak, , followed by a decline which is usually approximated by an exponential function. The empirically determined equation of the pressure profile [9, 12, 17] has the following form: where is the time elapsed after the detonation of charge in (msec); is the arrival time of shock wave to the target after the detonation of charge in (msec); is the peak pressure magnitude (MPa) at the shock wave front, and is the time decay constant which describes the exponential decay in (msec) and is equal to the time taken by the peak pressure to fall to 1/ times its initial value. The peak pressure and decay constant depend upon the size of the explosion and the standoff distance from this charge at which pressure is measured. The peak pressure, , and decay constant, , can be expressed as follows [18]: In case of TNT charge, the maximum bubble radii of explosive gas,, and first pulsation periods, , are expressed as [19]: where , and , are constants that depend on explosive charge type when different explosives are used. These input constants are as stated in Table 1 [20]. is the weight of the explosive charge in (Kg); is a charge depth in (m), and is the distance between the explosive charge and target in (m). Also Cole [9] provided further information regarding the systematic presentation of physical effects associated with underwater explosions.
(a)
(b)
When pressure from an underwater explosion impinges upon a flexible surface, such as the hull of surface ship, the reflected pressure on the fluid-structure interaction surface can be reasonably and accurately predicted using Taylor’s plate theory [12, 15]. For an air backed plate of mass per unit area (m) subjected to an incident plane shockwave, , a reflection wave of pressure, , leaves the plate which is moving at velocity. Using Newton’s second law of motion,
The fluid particle velocities behind the incident and reflected shockwave are and , respectively; thus, the velocity of the plate becomes
The incidence and reflective shockwave pressures are and , respectively, where is fluid density and is the acoustic velocity. By substituting the pressure into (4) and solving with (1), becomes Equation (4) can then be rewritten as
Differentiating (7) yields the following expression for plate velocity, where and . The total pressure on the plate is
A shock factor (SF), which is proportional to the energy density of the shockwave arriving at a structure, due to various combinations of charge weight and standoff distance is derived by
3. Models Description
Six different steel plates are considered using three-dimensional parts with an extruded shell base feature available in ABAQUS/Explicit version 6.10. All the plates are 20 mm in thickness and 3000 × 3000 mm2 with rectangular stiffeners 30 mm in thickness and 100 mm in height. Stiffeners are created by adding an extruded shell feature which implies perfect connection without any additional constraint. Care is taken not to overlap the material of the stiffener by offsetting its reference surface from mid-surface, thus avoiding the possibility of additional stiffness at the junction in real fabrication. Figure 2 shows the different stiffener configurations used in the numerical study.
4. Numerical Model of the Plates
The nonlinear finite element program ABAQUS/Explicit is used to undertake a three-dimensional (3D) analysis of the problem.
4.1. Models Geometry
ABAQUS/Explicit offers an element library for a wide range of geometric models. In the present study, the fourth node shell element (S4R-6 DOF at each node) with reduced integration and hourglass control was used to model the geometry of the plates and stiffeners. Six different models consisting of grids of shell elements of size 0.075 m were used as shown in the Figure 2. The fluid region of the model is represented by an assemblage of 4-node acoustic tetrahedral elements (AC3D4). In explicit dynamic analysis, the shortest wavelength can be estimated using the highest frequency in the UNDEX load or prescribed boundary conditions [21]. Since the acoustic wave length decreases with increasing frequency, there exists a highest frequency for a given acoustic domain mesh. Based on inequality: where is the upper frequency range of the shockwave; represent the maximum internodal interval of an element in a mesh; is the number of internodal intervals per acoustic wavelength ( is recommended); is the acoustic speed of the fluid; and are the bulk modulus of the acoustic medium and its density.
Artificial bulk viscosity is activated to properly represent propagation of the induced compressive stress wave by employing quadratic and linear functions of volumetric strain rates with default values of 1.2 and 0.06, respectively [21]. The outer boundary of the external fluid is represented by half cylindrical surface as shown in Figure 3.
4.2. Boundary Conditions and Fluid-Structure Coupling
The panel on the ship’s frame is typically stiffened by beams or stringers; thus, the panel can be divided into many small panels. The restraining moment of the borders of these panels is the torsion rigidity of a girder of stringer. During analysis, fully clamped boundary conditions are imposed on the four sides of the panels.
The boundaries of the fluid may cause shockwave refraction or reflection, resulting in its superposition or cancellation by the incident wave [21, 22]. To prevent this phenomenon, the boundary condition of the fluid element is set as a nonreflective boundary during analysis except for the free surface where zero pressure boundary condition was applied to it as shown in Figure 3. The acoustic-structural interaction between the wet surfaces of the plat and the acoustic interaction surfaces (the wetted interface) was implemented by use of a surface-based “tie” constraint. The location of the charge and the stand-off point defined as reference points, prior to the interaction, and the effect of the UNDEX event are transferred to the structure by means of the incident wave loading feature [16].
4.3. Material Properties
The stiffened panel is made of mild steel. The numerical model uses the constitutive law for elastic/plastic materials to model the stiffened panel. Isotropic hardening rules are adopted in the hardening model. The parameters of steel used in the numerical model are as follows: Poisson ratio is equal to 0.3 and mass density is equal to 7800.0 kg/m3. The initial yield stress is 300 MPa and the yield stress increases to 400 MPa at a plastic strain of 35%. Table 2 shows the plastic material properties for steel used in this study.
When the material sustains momentary dynamic loading the effect of the strain rate causes the material’s dynamic strength to exceed the strength during a static experiment. Thus, the effect of strain rate must be considered during the analysis to match actual situations. As recommended by Jones [23], this study adopts the Cowper-Symonds strain rate mode as follows: where is the material’s dynamic yielding stress; is the material’s static yielding stress; is strain rate, and and are material parameters whose values are normally s−1 and accepted values for steel, respectively.
The fluid region of the model is represented by the acoustic fluid domain. Its properties are the bulk modulus and density. In this numerical investigation, commonly accepted values for the sea water were stated in [16]. The bulk modulus is 2140.4 MPa and the density of the sea water is 1000 kg/m3.
4.4. Convergence Study
To verify the accuracy of analytical findings a convergence study has been carried out to choose the optimum mesh for the model. The finite element mesh contains one variable () (number of divisions along the plate side) that affects the number of elements in the model. Varying the number () affects the accuracy of the results. The elements throughout the convergence study have aspect ratio 1 : 1.
The steel plate of model (1) has been used to perform the convergence study on it. This plate is a square steel plate with fixed sides of side length 3000 mm and thickness 20 mm. Three different models consisting of grids of shell elements of size 0.0375 m, 0.075 m, and 0.15 m representing fine, medium, and coarse meshes, respectively, were used to verify the accuracy of the finite element models of the plates. The plate is loaded due to the UNDEX with charge mass 7 kg TNT at standoff distance of 5 m from the central point of the plate. For each mesh considered, the maximum displacement in the middle is compared.
The results of the convergence study show that the maximum displacement for the model consisting of grids of shell elements size of 0.15 m is 11.42 cm. For the model consisting of grids of shell elements size of 0.075 m, the maximum displacement is 11.72 cm. It can be noted that the difference is 2.56%. Also, the maximum displacement for model consisting of grids of shell elements size of 0.0375 m is 11.8 cm and the difference is only 0.67%. Consequently, the mesh with elements size of 0.075 m is chosen to perform the whole analysis since all displacement are less than 1% different from the results for the case of elements size 0.0375 m and give minimum possible element size and analysis time. Figure 4 shows the displacement distribution at the central point of the plate of the three different models with time.
5. Results and Discussions
This study is based on the plate’s performance improvement ratio , which can be expressed as the ratio between the performances of unstiffened to stiffened plate for a specified case. All models are subjected to UNDEX with different amounts of explosives (5, 7, and 9 kg) at standoff distance that is equal to ( m) representing SF (0.2, 0.238, and 0.27). The parameters that play an effective influence on improving the plate’s performance include the stiffeners configurations, damping effect, strain rate sensitivity and the structural integrity.
5.1. Effect of Stiffeners Configurations
The submerged or floating structure within the vicinity of an UNDEX will be subjected to loading from both the shock wave and the bubble pulse pressure waves and its performance according to the strength of these waves and the resilience of the structure.
Figure 5(a) shows the maximum central-point displacement under SF of (0.2, 0.238, and 0.27). Generally, it can be noticed that the maximum central-point displacement for the models increase with the increase of SF for each model. This is axiomatic because of increasing the explosive charge weight with fixed standoff distance.
(a) Max displacement at the center of the plates
(b) Improvement ratio
As shown in Figure 5(a), at shock factor equal to 0.2, the maximum displacement for model 1 is 99.69 mm. The introduction of one simple stiffener results in decreasing the peak central point displacement to 69.37 as in model 2, while for models 3, 4, 5, and 6 the maximum central-point displacement is 70 mm, 57.82 mm, 55.8 mm, and 62.93 mm, respectively. Similar trend of results is observed for all other models at SF equal to 0.238. The maximum displacement for all models is greater than the maximum displacement at shock factor equal to 0.2. The maximum displacement for model 1 is 11.72 mm while for models 2, 3, 4, 5, and 6 the maximum central-point displacement is 8.65 mm, 8.6 mm, 7.4 mm, 7.02 mm, and 7.62 mm, respectively.
At SF equal to 0.27, the maximum displacement for model 1 is 13.27 mm while for models 2, 3, 4, 5, and 6 the maximum central-point displacement is 10.6 mm, 10.05 mm, 8.77 mm, 8.24 mm, and 8.73 mm, respectively. From the above result, it can be found that at low SF the maximum central-point displacement of model 2 is better than model 3 and, by increasing the SF until 0.27, the maximum central-point displacement of model 3 becomes better than model 2. Also, the minimum displacement occurs at model 5. Therefore, the configurations of stiffeners have an important influence on the response of the stiffened plates and leads to better response of the plates to different shock waves; this result is matching to the results reported in [24–26].
Figure 5(b) shows the improvement ratio () for all models. As shock factor (SF) increases the decreases, at SF equal to 0.2, for model 2 is 30.41%, while for models 3, 4, 5, and 6 the is 29.78%, 42%, 44.03%, and 36.87%, respectively.
The best is for model 5 and the lowest is for model 3. By increasing the SF to 0.238, it can be found that the best is 40.07% for model 5 and the lowest is 26.19% for model 2. When the SF increases to 0.27, is 20.12%, 24.27%, 33.93%, 37.94%, and 34.21% for models 2, 3, 4, 5, and 6, respectively. This indicates that the stiffened steel plate with one cross stiffener configuration (model 5) exhibits the highest response reduction. Also, the lowest is 20.12% for model 2.
For better understanding and evaluating the performance of the plates to the shock loading, the displacement of the center point of all plate’s models is monitored with time under the effect of different shock loads as shown in Figure 6. It can be further noted that the time it takes for each model to reach maximum displacement at the same position is roughly the same.
Also, the transverse center-line deformation of all models is recorded as shown in Figure 7. It can be noted that the maximum deformation occurs at the center of the plates and reduced gradually towards the fixed boundaries. From Figure 7(a), model 2 has the maximum deformation until the distances are equal to 1.1 m, after that it decreases towards the center of the plates. Model 6 has the minimum displacement from distance equaling 0.9 m to 1.1 m after that the displacement increases rapidly towards the center of the plate. Also, model 5 has the minimum displacement.
(a) All models at SF = 0.2 (12.5 msec)
(b) All models at SF = 0.237 (12.5 msec)
(c) All models at SF = 0.27 (12.5 msec)
(d) All models at SF = 0.2 (15 msec)
(e) All models at SF = 0.238 (15 msec)
(f) All models at SF = 0.27 (15 msec)
Figures 7(b), 7(c), 7(d), and 7(f) show the transverse center-line deformation. Model 1 has the maximum deformation and this deformation increases linearly as the distance from the edge of the plate increases and reaches the maximum value at the center of the plate. Also model 5 has the best performance at all SFs. Figure 8 shows the transverse center-line deformation for all models at different SFs. The displacement distribution contours with time are presented in Figure 9.
5.2. Strain Rate Sensitivity
In order to investigate the strain rate sensitivity on the plate’s performance under the effect of different shock loading, the analysis has been first carried out ignoring this effect, after that the strain-rate effect is included by adjusting the material dynamic yield stress at each Gauss point according to the Cowper-Symonds strain rate mode as mentioned in (11). In this investigation the material parameters defined as s−1 and being the commonly accepted values for steel have been used as recommended by Jones [23].
Figure 10 shows the displacement history of the central node for both with and without rate dependant for all models at different SFs. It shows the effect of inclusion of the strain rate on the performance of the different plates. The material with rate dependant analysis leads to decreasing the maximum displacement, which assures that the results are very sensitive to the material data. Figure 11 shows the maximum displacement values of the central node for both with/without rate dependant analysis and the improving ratio at SF equal to 0.2, 0.238, and 0.27, respectively.
(a) Max displacement at the center of the plates
(b) Improvement ratio
It can be concluded that, for model 1, the maximum displacement decrease from 99.69 mm, 117.2 mm, and 132.7 mm to 84.88 mm, 101 mm, and 114.7 mm at SF equal to 0.2, 0.238, and 0.27, respectively. For model 2, the maximum mid-point displacement without strain rate is 69.37 mm, 86.5 mm, and 106 mm and with strain rate are 60.8 mm, 74.67 mm, and 87.06 mm with equal to 39%, 36.28%, and 34.39% at SF equal to 0.2, 0.238, and 0.27, respectively.
For model 3, the maximum mid-point displacement with strain rate is 60.88 mm, 73.16 mm, and 84.73 mm with equal to 38.9%, 37.57%, and 36.15% at SF equal to 0.2, 0.238, and 0.27, respectively.
For model 4, the maximum mid-point displacement with strain rate is 49.89 mm, 62.71 mm, and 73.93 mm with equal to 49.95%, 46.49%, and 44.28% at SF equal to 0.2, 0.238, and 0.27, respectively.
Model 5 has the least deformation and the material with rate dependant analysis leads to decreasing the maximum displacement to 47.4 mm, 59.67 mm, and 69.61 mm with equal to 52.45%, 49.08%, and 47.54% at SF equal to 0.2, 0.238, and 0.27, respectively.
For model 6, the material with rate dependant analysis leads to decreasing the maximum displacement to 52.7 mm, 64.4 mm, and 74.1 mm with equal to 47.13%, 45.05%, and 44.16% at SF equal to 0.2, 0.238, and 0.27, respectively.
Finally it can be concluded that the best is 52.45% for model 5 and the lowest is 38.93% for model 3 at SF equal to 0.2. By increasing the SF the decreases, the best is 49.08% also for model 5 and the lowest is 36.28% for model 2 at SF equal to 0.238. At SF equal to 0.27 the best is 47.54% for model 5 and the lowest is 34.39% for model 2. The results show that at low SF the performance of model 2 is better than model 3 and by increasing the SF the performance of model 3 is better than model 2 which refers to the importance of different SFs in this analysis.
The results show that when the strain-rate effect is taken into account, the yield stress increases because the elastic modulus is higher than the plastic modulus. It is noted that the analysis with strain rate will be much stiffer, resulting in a decrease in the mid-point displacement as in [26–28].
5.3. Damping Effect
Under the effect of load applying, the undamped structures continue to vibrate with constant amplitude. A constant amplitude vibration is not the response that would be expected in practice since the vibrations in this type of structure would tend to die out over time and effectively disappear. The energy loss typically occurs by a variety of mechanisms, so it is needed to consider the presence of damping in the analysis to model this energy loss. The results of the damped analysis clearly show the effect of mass proportional damping.
Figure 12 shows the displacement history of the central node for both the damped and undamped material for all models at different SFs. It show that the material with damping analysis leads to decreasing the maximum displacement. For model 1, the maximum displacement decrease from 99.69 mm, 117.2 mm, and 132.7 mm to 93.7 mm, 110.1 mm, and 125.4 mm at SF equal to 0.2, 0.238, and 0.27, respectively. Also Figure 12 shows that the time it takes for all the models to reach maximum displacement at the same position is roughly the same.
Figure 13(a) shows the values of maximum central displacement in case of damped and undamped material and, at SF equal to 0.2, the best performance is model 5 which has the least deformation. The material with damping analysis leads to decreasing the maximum displacement up to 47.97% (model 5), 46.35% (model 4), 40.31% (model 6), 35.45% (model 2), and 34% in model 3. These values are based on the comparison to the value obtained from model 1 without damping. Figure 13(b) shows the improving ratio in case of damped material analysis based on the comparison to the value obtained from model 1 with damping. It can be noted that as the SF increases the decreases and the best performance is model 5 with of 44.6%, 40.53%, and 38.7% at SF equal to 0.2, 0.238, and 0.27, respectively. This confirms that the proposal with stiffeners can help the structure to sustain shock loads resulting from an underwater explosion and this result is matching the one in [27]. So, it can be concluded that the damping properties play a very important role in the response prediction due to underwater shock load.
(a) Max displacement at the center of the plates
(b) Improvement ratio
5.4. Structural Integrity
Structural integrity is one of the major concerns of such analysis. The equivalent plastic strain in a material (PEEQ) is a scalar variable that is used to represent the material’s inelastic deformation. If this variable is greater than zero, the material has yielded [17]. Figure 14 shows a contour plot of accumulated PEEQ for all models and permanent deformations are monitored which emphasize that the largest permanent deformation is obtained on model 1 and model 6. The FE analysis predicts the central peak observed and permanent deformations which first occurred at the middle of the sides of the plates and then progressed towards the corners according to experimental results reported in [24–26]. The lowest permanent deformation is obtained on model 2, model 4, and model 5. Figure 15 shows a contour plot of von Mises stress at the end of the analysis. The regions of high stress concentration are analyzed which confirm that model 5 offered greater resistance and its is 16.13%. It is found that for models 4 and 6 is 13.05% and 4.86%, respectively. The lowest is for model 3, which matches the experimental results reported in [24, 25, 29].
6. Conclusions
From the numerical simulation of stiffened plates subjected to noncontact underwater explosion the following can be concluded.
The displacement-time histories under different shock loadings are presented which will be used in designing stiffened panels to enhance resistance to underwater shock damage. The effect of stiffener configurations is very important, since it can drastically affect the overall behavior of the plates as indicated. The proposal with special damping system can help the structure to sustain shock loads resulting from an underwater explosion. Model 5 offered greater resistance to UNDEX load than all models which assure that stiffener location influenced the plate response. The inclusion of strain-rate effect and the proposal with special damping system results in lower mid-point displacement. Thus, the damping and the strain-rate effect should be taken into account when analyzing structures subjected to underwater shock loading. The inclusion of strain-rate effect is more effective than the proposal with special damping system to reduce the displacement of the plate.
Conflict of Interests
The authors declare that there is no conflict of interest is regarding the publication of this paper.
Authors’ Contribution
Professor Tong Lili and Dr. Mahmoud Helal have contributed equally to this work.