Shock and Vibration

Shock and Vibration / 2018 / Article

Research Article | Open Access

Volume 2018 |Article ID 3157457 |

Hrvoje Draganić, Damir Varevac, "Analysis of Blast Wave Parameters Depending on Air Mesh Size", Shock and Vibration, vol. 2018, Article ID 3157457, 18 pages, 2018.

Analysis of Blast Wave Parameters Depending on Air Mesh Size

Academic Editor: Isabelle Sochet
Received28 Feb 2018
Revised27 May 2018
Accepted11 Jun 2018
Published08 Jul 2018


Results of numerical simulations of explosion events greatly depend on the mesh size. Since these simulations demand large amounts of processing time, it is necessary to identify an optimal mesh size that will speed up the calculation and give adequate results. To obtain optimal mesh sizes for further large-scale numerical simulations of blast wave interactions with overpasses, mesh size convergence tests were conducted for incident and reflected blast waves for close range bursts (up to 5 m). Ansys Autodyn hydrocode software was used for blast modelling in axisymmetric environment for incident pressures and in a 3D environment for reflected pressures. In the axisymmetric environment only the blast wave propagation through the air was considered, and in 3D environment blast wave interaction and reflection of a rigid surface were considered. Analysis showed that numerical results greatly depend on the mesh size and Richardson extrapolation was used for extrapolating optimal mesh size for considered blast scenarios.

1. Introduction

Today's global society faces many challenges that have not been present before. Not only are military facilities exposed to terrorist threats, but also many civilian buildings and structures are in danger. Most of these potential targets cannot be protected since they are in public use and have to be accessible. This paper is part of a wider research of the impact of explosions on overpasses that cross over the public roads [1, 2]. The imagined scenario assumes that the vehicle, in which explosive is placed, is parked underneath the overpass and then detonated. Since it is impossible to carry out the tests in the 1:1 scale and handling explosive is very dangerous and requires special conditions, the most practical tool is computer simulation. In both cases, carrying out field tests or computer simulations, it is very demanding task. Due to the very nature of the explosion process and many parameters that affect measured effects, even the field tests cannot provide accurate and repeatable results. On the other hand, in the computer simulation some assumptions and idealizations must be adopted. Very important issues in these simulations are hourglass effect, material properties, fragmentation of material, and fireball effect [3] but the most important task is to determine adequate mesh size. Coarse mesh may result in useless results, but too fine mesh may occupy computer resources for hours and often days in case of complex problems. The purpose of this paper is to offer guidelines for selection of appropriate air mesh size. Two scenarios were considered for the purpose of further clarifying the influence of air mesh size on simulation results and enabling further development of blast loading analysis. The first scenario is ideal free-air explosion where incident pressures were calculated and the second is specific blast scenario where reflected pressures from the overpass superstructure were calculated. For this latter case, the convergence tests were carried out.

Several researchers have studied the mesh size effect, individually or as a part of broader research on using numerical tools, to study blast wave propagation and interaction with structures. Chapman et al. [4] conducted a comprehensive parametric study of model grid size to obtain optimal mesh size as a part of their research on the use of Autodyn 2D for simulation of blast wave interaction with structures. Their results were compared to those from existing simple experimental tests, and they concluded that Autodyn 2D is a suitable tool for blast wave investigation. Luccioni et al. [5] studied influence of mesh size on blast wave propagation as a part of their research on urban blast. They concluded that a finer air mesh (10 cm) gives better results, but this can be too expensive in the terms of computational time. Alternatively, a coarser air mesh (50 cm) may be used to obtain qualitative results for comparison. Shi et al. [6] studied the relationship between the scaled distance and blast load parameters in a series of numerical simulations with different mesh sizes and compared the results with those from TM5-1300 [7]. They found that the pressure and impulse are less sensitive to mesh size with an increase of the scaled distance; they also proposed a method for the adjustment of peak pressures for coarse meshes. Nam et al. [3] determined the maximum element size to ensure that analytical results are independent of mesh size; this eliminates load gradient differences and reduces the gap between explosive and internal energy. Deng and Jin [8] conducted a parametric study of free-air blast wave propagation as a part of their research on bridge damage induced by blast loading. Simulation results indicated that a 5 cm air mesh size guarantees effective blast loading and reliable computational results. Wilkinson et al. [9] in their study on modelling close-in detonations in the near field compared numerical results of overpressure obtained by Autodyn and Air3D. Mesh sensitivity analysis was conducted and results clearly indicated that mesh dependency of calculated results is existent. Shin et al. [1012] performed numerical 1D, 2D, and 3D simulation of close-in detonation of high explosives in order to validate CFD code for determining incident and reflected overpressure and impulses. It was concluded that mesh refinement studies are needed to establish cell sizes for air blast calculations and that a cell size of R/500 should be sufficient to predict overpressure histories in the near field, where R is the distance from the centre of the charge to the monitoring location. In addition, it was determined that the afterburning has a little effect on peak overpressure but can increase impulse if sufficient oxygen is present. Shin et al. [1214] as a part of their research on blast parameters of detonation of spherical high explosives in free-air performed mesh size analysis on a wide range of scaled distances. Analysis produced recommendation of new polynomials and design charts for detonation parameters of spherical TNT charges. It was concluded that cell size of R/500 is sufficient for calculating pressures in the scaled distance range of 8 m/kg1/3, similar to previous research [10].

Referred numerical simulations show that the accuracy of numerical results is strongly dependent on mesh size and blast scenario—a mesh size acceptable for one blast scenario might not be suitable for another. This is a clear indicator that a numerical mesh size convergence test needs to be carried out prior to large-scale numerical simulations.

With respect to these results and recommendations, for the purpose of this research, it was necessary to determine the amount of explosive which is likely to be placed under the overpass. Table 1 gives crude estimations of possible TNT quantities that can be placed in common vehicles and used for explosive deployment under a bridge [21, 22]. The first three quantities are assumed to be the most probable and may offer some hope of bridge survival; detonation of larger TNT quantities would for sure cause imminent bridge collapse. Physical characteristics of vehicles, which dictate the shape, and confinement of the charge, detonation point, and consequent flying debris (fragmentation) are not taken into consideration, and only transport capacity was used in order to estimate the amount of explosive.

Vehicle typeCharge mass [kg]

Compact car trunk115

Trunk of a large car230

Closed van680

Closed truck2270

Truck with a trailer13610

Truck with two trailers27220

Note: 1 pound = 0.45 kg.

2. Blast Wave Parameters

Explosion is an event of sudden and rapid energy release accompanied by a large pressure and temperature gradient on the wave front. It is very important to notice that blast parameters for free-air explosion are much more different from parameters which dictate the influence of explosion on engineering structure, in this case the overpass.

The blast wave is described by several parameters, but the main one is scaled distance, Z, which depends on charge weight, W [kg], and charge distance to target, R [m]. Scaled distance is calculated by the following expression:

Use of TNT as a referent explosive for determining Z is widely accepted. The first step in quantifying the blast waves from a source other than TNT is to transform the charge mass into the equivalent mass of TNT. In general, the TNT equivalent represents the mass of TNT that would result in an explosion with the same effect as the unit weight of explosive under consideration. It is defined as the ratio of the mass of TNT to the mass of the explosive resulting in the same magnitude of blast wave (or impulse, pressure, and energy) at the same radial distance for each charge, which assumes the scaling laws of Sachs and Hopkins. Common way to make conversion is to multiply the explosive mass by the conversion factor determined using one of several existing methods (the specific energy based concept, pressure based concept, impulse based concept, Chapman-Jouguet state based concept, etc.).

The primary reason for choosing TNT as the referent explosive is that most of the available experimental data regarding the characteristics of blast waves were collected using TNT. Several methods exist for determining the explosive characteristics of different explosives, but they do not yield the same values for the TNT equivalent. These values depend on the characteristic parameter of the blast wave, the geometry of the load, and the distance from the explosive charge [23]. The TNT equivalence of terrorist-manufactured explosive material is difficult to define precisely because of the variability of its formulation and the quality control of its manufacture. Table 2 gives the TNT equivalents for the most commonly used explosives [24].

ExplosiveSpecific energyTNT equivalent

RDX (Cyclonite)53601.185

Nitro-glycerine (liquid)67001.481




Note: 1 pound = 0.45 kg.

Pressure-time history for a blast wave is commonly described by the Friedlander equation [24]:where represents the blast wave overpressure, b is the waveform parameter, t0 is the positive phase duration, and t is considered time. Peak overpressure is the maximum blast pressure above normal atmospheric pressure (p0 = 101.325 kPa), which is dependent on the distance from the device. It can be calculated by numerous analytical formulations developed by researchers [2527] who were able to carry out experiments on blast wave propagation. Differences in results exist in peak overpressures for near- and far-field explosions, and they are incorporated in these formulations.

As the blast wave travels away from the detonation point, it undergoes reflection when the forward moving air molecules are brought to rest and further compressed upon meeting an obstacle. Combining peak overpressure and dynamic pressure, Rankine and Hugoniot derived the following equation for reflected overpressure:where represents the atmospheric pressure.

Further important parameters include specific impulse, time of arrival, and duration of positive and negative phase. Specific impulse represents the area beneath the pressure-time curve from arrival time to the end of the positive phase. Specific impulse is given by the following equation:where represents blast wave time of arrival, and is the duration of the positive phase. These specific blast times as well as other blast parameters can be obtained from US Army manual UFC 3-340-02 [7].

3. Numerical Simulations

Numerical models developed for this investigation were created in Ansys Autodyn hydrocode software [28], which is intended for fluid dynamic analysis. Autodyn is a fully three-dimensional explicit finite difference program which can utilize several different numerical techniques (Eulerian, Lagrangian, ALE, SPH, and Shell) to optimize the analysis of nonlinear dynamic problems. The structural analysis requires the ability to simulate both fluid/gas behaviour and structural response.

3.1. Numerical Techniques

The Lagrange processor is used for modelling solid continua and structures and operates on a structured (I-J-K) numerical mesh of quadrilateral (2D) or brick-type elements (3D) [29]. The numerical mesh moves and distorts with the material motion; no transport of material occurs from cell to cell. The advantage of such a scheme is that the motion of material is tracked very accurately, and the material interface and free surfaces are clearly defined. The primary disadvantage of the Lagrange formulation is that for severe material deformations or flow the numerical mesh will also become highly distorted, with attendant loss of accuracy and efficiency or outright failure of calculation [30].

The Euler processor is used for modelling fluids, gases, and large distortions. These processors include first-order and second-order accuracy schemes. A control volume method is used to solve the equations that govern conservation of mass, momentum, and energy. Numerical mesh is fixed in space, and material flows through it. The advantage of such a scheme is that large material flows and distortions can be easily treated. The disadvantage is that material interfaces and free surfaces are not naturally calculated, and sophisticated techniques must be utilized to track material interfaces. Additionally, history-dependent material behaviour is more difficult to track.

ALE (Arbitrary Lagrange Euler) is a hybrid processor wherein the numerical mesh moves and distorts according to user specification. An additional computational step is employed to move the grid and remap the solution onto a new grid. ALE is an extension of the Lagrange method and combines the best features of both Lagrange and Euler methods. Due to advantages of both Lagrange and Eulerian solver, ALE solver was used for model calculation.

Modelling with Ansys using explicit time integration is limited by the Courant-Friedrichs-Levy condition [31]. This condition implies that the time step is limited so that a disturbance (e.g., stress wave) cannot travel further than the smallest characteristic element dimension in the mesh in a single time step. Thus the time step criterion for solution stability iswhere dt is the time increment, dx and dy are characteristic dimensions of an element, and c is the local material sound speed in an element. Based on this condition one can conclude that smaller mesh size implies smaller time steps and consequently longer calculation time.

3.2. Material Models

The numerical models used in this study include three different materials: environment, explosive charge, and overpass superstructure. Environment was modelled as the air, which transmits pressure waves, presented as an ideal gas, and the pressure was related to energy using following expression:where ρ is air density, e is specific internal energy, and γ is a constant (i.e., ratio of specific heats). This expression is one of the simplest forms of equation of state for gases. Gaseous materials have no ability to support any kind of stresses, and consequently no strength or failure relations are associated with this material.

Explosive charge was modelled as a spherical TNT charge. Similar to the air model, no strength or failure relations are associated with explosive material. Detonation is a rapid process that converts explosive material into gaseous products, and it is usually completed at the start of a given simulation. TNT explosive was modelled with the Jones-Wilkins-Lee equation of state [32] using following expression:where p is hydrostatic pressure, χ is specific volume (1/ρ), and e is specific internal energy, and A, R1, B, R2, and ω are constants experimentally determined for several common explosives [33]. The shape of TNT charge in this study is spherical in order to ensure uniform propagation of the pressure waves and thus the clarity of results. Material properties of TNT and air are given in Table 3.



Equation of stateJWL

Reference density1.63[g/cm3]

Parameter A3.738E+08[kPa]

Parameter B3.747E+06[kPa]

Parameter R14.15

Parameter R20.90


C-J detonation velocity6.930E+03[m/s]

C-J energy/unit volume6.000E+06[kJ/m3]

C-J pressure2.100E+07[kPa]

Initial energy3.681E+06[mJ/mg]

Note: 1 inch = 2.54 cm; 1 pound = 0.45 kg; 1 psi = 6.89 kPa; 1 K = -272.15°C; †JWL = Jones-Wilkins-Lee.


Equation of stateIdeal gas

Reference density1.225E-03[g/cm3]


Reference temperature288.20[K]

Specific heat717.60[J/kgK]

Initial energy2.068E+05[mJ/mg]

Note: 1 inch = 2.54 cm; 1 pound = 453.59 g; 1 psi = 6.89 kPa; 1 K = -272.15°C.

In order to reduce the calculation time and get clear results, rigid material was chosen for the slab. This is a special feature of Autodyn software which may be used if the stresses and/or deformations of target object are not of interest. Insensitivity of the results on the rigid slab mesh size was confirmed by conducting several convergence tests and therefore is eliminated as an influential parameter in determining the dependency of blast wave parameters on air mesh size. Ground surface on which TNT charge was placed in 3D simulation before detonation was modelled as rigid boundary in order to provide full reflection of the blast wave and also eliminate ground mesh size as an influential parameter. Axisymmetric free-air simulations were performed using Eulerian solver, and 3D simulations were done using a coupled Lagrangian/Eulerian solver as it allows the combination of the best aspects of Lagrangian and Eulerian solvers.

3.3. Parametric Study

Modelling process of the impact of the explosion on overpass structure is divided into two parts. First, 2D FE model of free-air explosion is created. This model is composed of circular TNT charge, surrounding air, and outflow boundary, each with its properties (see Table 3) (free-air axisymmetric model). Radius of the TNT charge is determined by the specific reference density and chosen weight while radius of the air surrounding may be arbitrary defined (in this case 500 cm in order to shorten the calculation time). Using axial symmetry, only one circle sector may be modelled (Figure 2) and calculated which greatly reduces the calculation time. The results from this FE model are used to calculate free-air incident pressures. Second part is generation of the 3D model by revolving a 2D section 360° degrees and remapping the incident pressures to the overpass structure and rigid ground surface in order to compute reflected pressures (Figure 3).

Since the explosive charge is part of the model, it should be investigated whether its own mesh size (different from air mesh size) affects the results of the simulation. It was considered if the refinements of the charge mesh influence the level of energy release during detonation, time of blast wave arrival, or incident pressures, especially in 3D simulations, and should be taken into consideration. Due to this concern, additional simulations were conducted in order to study and determine the influence of charge mesh size. Simulations were conducted by varying air and charge mesh size for each of the three chosen quantities of explosive. Three charge weights, as mentioned earlier, 115 kg, 230 kg, and 680 kg of TNT (different line type), three air mesh sizes (different point type), and eight charge mesh sizes were considered. Numerical model is axisymmetric as shown in Figure 2, consisting of two separate parts, TNT charge and air surrounding. Simulation results are shown in Figure 1. Figure 1(a) shows normalized pressures while Figure 1(b) shows normalized total blast energy. Table 4 shows the values used for normalization, for each of the three air mesh sizes and three TNT quantities. These values were obtained in separate simulation.

Air mesh size  
Charge mesh size  
TNT (kg)
Total energy  
Total energy  
Total energy  




Note: 1 inch = 25.4 mm; 1 pound = 0.45 kg; 1 psi = 6.89 kPa.

Normalization is conducted in regard to the numerical model with the same charge and air mesh size; i.e., referent models were those with 5 mm, 10 mm, and 20 mm as both charge and air mesh sizes. Comparison of normalized total energies indicates that differences are minimal for all three observed parameters (under 5%), with somewhat larger but still acceptable differences in normalized pressures (under 10%). Due to minimal differences, charge meshes to coincide with air mesh in every simulation. This conclusion also has an impact on accelerating the modelling process because there is no need to mesh separately air and charge part of the model.

Free-air axisymmetric numerical simulation of explosion was performed using twelve different air mesh sizes (Table 5). Here, only the TNT in the air environment was modelled, and the wave propagation was observed. On the edge of the air environment, an outflow boundary condition was implemented to avoid the blast wave reflection and amplification; that is, to allow free gas outflow from designated air environment. Observed air environment is 5 m radius, with gauge points placed at equal distances (30.5 cm) in order to compare results for different scaled distances; air environment and gauge position are shown in Figure 2.

Mesh size [mm]0.390.781.563.136.2512.52550100200400800

Number of elements128216410320515988004002001005025136

Note: 1 inch = 25.4 mm.

The 3D model binds the air environment with a rigid slab, which simulates bridge superstructure and ground level for blast wave reflection and amplification. The Eulerian multimaterial solver was used for accurate simulation of hot gas expansion and interaction with structure primarily because there are no grid distortions during fluid flow through cells, and it allows large deformations, mixing of initially separate materials, and larger time steps. Its use requires more computations per cycle, finer zoning for similar accuracy, and extra cells for potential flow regions.

Parametric 3D numerical simulation of a ground explosion was performed using nine different air mesh sizes (Table 6). The mesh size of 800 mm was discarded because it was too coarse and gave unrealistically small pressures in comparison to AT-Blast [34]. Mesh sizes smaller than 3.13 mm were impractical because of high requirements for computer hardware components. The ground explosion model and gauge scheme are shown in Figure 3. The ideally rigid concrete slab was placed 5 m above the detonation point on the ground, which is typical clearing underneath the overpass. In addition, using the symmetry, only one-quarter of the numerical model presented in Figure 3 was modelled and simulated. Additional shortening of calculation time was achieved by remapping the final solution of axisymmetric analysis for selected charge quantities (115 kg, 230 kg, and 680 kg) as an input (load) to the 3D environment [35].

Mesh size [mm]3.136.2512.52550100200400800

Number of elements82944000010368000012960000162000020250025313316439649

Note: 1 inch = 25.4 mm.

4. Comparison of Results and Discussion

Numerical simulation results using coarse and fine air mesh sizes varied considerably. This indicates that mesh size convergence tests are needed to obtain the optimal mesh size for an accurate simulation of blast wave propagation. Detonation products influence blast wave pressures and consequently impulses; that is, they are within the blast fireball. The pressures depend on the type of explosive as each explosive has different detonation products. The exact influence of the fireball on the pressures could not be determined separately because additional analysis would be needed in which data are collected outside the fireball radius and compared to those obtained within the fireball.

4.1. Incident Pressures

Results obtained by free-air axisymmetric numerical simulation were compared with results acquired by AT-Blast, the software for analytical determination of the referent blast parameters. AT-Blast uses equations for blast estimation of overpressures at different ranges which were developed by Kingery and Bulmash [36]. These equations are widely accepted as engineering predictions for determining free-field pressures and loads on structures. Kingery and Bulmash have compiled results of explosive tests in which they detonated charges weighting from less than 1 kg to over 400,000 kg and used curve-fitting techniques to represent the data with high-order polynomial equations. Another reason why AT-Blast was chosen as the reference for the blast parameter comparison was that the software is publicly available, subjected to ITAR (International Traffic in Arms Regulations) controls. Unlike AT-Blast, ConWep [37], another widely used software for estimating blast properties, is not available to non-US citizens or organizations; therefore, it was not used as an additional control for the numerical results.

It can be seen that decreasing the mesh size allowed maximum incident pressure to converge to the referent value. With smaller mesh sizes (from 0.39 to 12.5 mm) maximum incident pressures exceeded values obtained by AT-Blast and can be accounted for by inaccuracies in the experimental pressure readings (the measuring instruments were influenced by very high pressures, and the blast wave duration and time of arrival were also very short, i.e., only a few milliseconds). The relative difference in maximum incident pressures measured at the last gauge point, 5 m from detonation point, between smaller mesh sizes and values generated by AT-Blast can clearly be seen from Table 7. The closest incident peak pressure was obtained for the 25 mm mesh size (marked by bold in Table 7), with about 2% maximum difference in relation to the AT-Blast value. Pressure histories for different mesh sizes relative to AT-Blast pressure values are plotted in Figure 4. It can be seen that the rate of pressure amplification from ambient to peak pressure becomes slower and the shape of the curve becomes flatter with the increase of the mesh size. The calculated incident peak pressure was also smaller with the larger mesh sizes.

Mesh sizeTNT mass
115 kg230 kg680 kg
Relative error %Pressure 
Relative error %Pressure 
error %


800 mm546.9857.20880.1257.121921.3751.41

400 mm777.4039.171039.2149.371738.8856.02

200 mm804.6737.041249.3939.132394.6839.44

100 mm1042.2918.441726.9115.863155.8320.19

50 mm1211.225.231888.427.993621.508.41

25 mm1288.310.812041.250.553877.961.93

12.5 mm1371.327.302128.353.704054.352.54

6.25 mm1429.7311.872185.076.464231.937.03

3.13 mm1454.8513.842224.998.414310.669.02

1.56 mm1476.0515.502258.5710.044386.3810.93

0.78 mm1495.1316.992290.7911.614537.9314.77

0.39 mm1494.0316.902299.1112.024565.7115.47

Note: 1 inch = 25.4 mm; 1 pound = 0.45 kg; 1 psi = 6.89 kPa.

If the positive peak incident pressure in relation to scaled distance is considered, its value decreases with the increase of scaled distance. This is in accordance with the initial assumption that with the increase of distance from detonation point to structure surface pressures exerted on structure surfaces are decreasing. Figure 5 shows that this is consistent for all air mesh sizes. With smaller air mesh size, pressures became closer to the referent values obtained by AT-Blast.

4.2. Reflected Pressures

Results obtained by 3D numerical simulation were also compared with the results acquired by AT-Blast. In terms of absolute values, the comparison gave the closest reflected pressures for 25 mm mesh size in an analysis for 115 kg and 230 kg of TNT; here, the minimum difference in relation to AT-Blast values was 8.57% for 115 kg and 0.08% for 230 kg. For 680 kg of TNT closest reflected pressure in comparison to AT-Blast was obtained for 12.5 mm mesh size with difference of 4.17% (Table 8). Reflected overpressure histories in relation to AT-Blast reflected overpressure values are plotted in Figure 6. Pressure amplification from ambient to peak reflected overpressure has a smaller rate of increase, and the shape of the curve becomes flatter with the increase of the mesh size as for the incident pressures. Pressure-time curve has additional pressure peak after the initial peak, which is caused by blast wave reflection and it can be clearly seen in Figure 6. Calculated reflected overpressure also becomes smaller with larger mesh sizes.

Mesh sizeTNT mass
115 kg230 kg680 kg
error %
Relative error %Pressure 
Relative error %


800 mm809.3689.291227.2991.071940.1493.68

400 mm1464.3780.612252.2283.614495.9785.35

200 mm2445.9067.624170.3369.658773.1971.42

100 mm5022.4233.518699.7636.6817584.1042.72

50 mm6890.838.7811566.1015.8124086.7021.53

25 mm8201.138.5713727.500.0829398.304.23

12.5 mm8875.4817.5015104.709.9431977.904.17

6.25 mm10161.6034.5215445.5012.4232894.807.16

3.125 mm10126.0034.0515568.2013.3232953,007,35

Note: 1 inch = 25.4 mm; 1 pound = 0.45 kg; 1 psi = 6.89 kPa.
4.3. Time of Arrival

Mesh size also influenced the blast wave time of arrival. Some authors [6] have concluded that the arrival time is independent of the mesh size, and if we consider that the time in blast problems is measured in milliseconds, differences arising from analyses with different mesh sizes exist. Tendency is that the smaller mesh sizes shorten the time required for the blast wave front to arrive to the designated measuring point. If the mesh size is sufficiently small, it no longer influences the arrival time, as can be seen in Figure 7. For larger mesh sizes for which pressure history is flatter, the time of maximum incident pressure is taken as the time of blast wave arrival.

If a difference of ±15% was considered, all arrival times were within this range. Results presented in Figure 7 show that the numerical analysis provided earlier arrival times in comparison to AT-Blast. Similar to incident pressures, there was a tendency for smaller mesh sizes to shorten the arrival time; if the mesh size was sufficiently small, it no longer had an influence on the arrival time.

5. Result Convergence

Individual application of numerical methods differs in the level of accuracy with which a numerical model can reflect the physical phenomena. An analysis with the help of complex numerical models is particularly applicable if there is no closed analytical solution and if the experiments are impossible or insufficient. In all such cases, the question arises about the ability of such analysis to correctly predict the results of the physical process in question. The first recommended procedure within the verification process is the mesh refinement study. Its main objective is to evaluate the error of discretization and to check if the applied mesh is sufficiently refined.

5.1. Richardson Extrapolation

The mesh refinement study is conducted based on a comparison of the results for at least two but usually three meshes [29]. The Richardson extrapolation serves as a higher-order estimate of the evaluated quantity. A quantity f calculated for a mesh characterized by parameter (mesh size) h can be expressed using Taylor’s theorem. In practice, the Richardson extrapolation is generalized for any, also noninteger, p-th order approximations and the mesh refinement ratio r and is considered as p + 1 order approximation. The exact solution of extrapolation can be described as an asymptotic solution for a mesh with element dimension h approaching zero. Estimate of the asymptotic solution is given by where denotes asymptotic solution for h approaching 0 (extrapolated value); and second-order approximations of ; r is refinement ratio and p is order of convergence. Refinement ratio, r, can be set as arbitrary when using unstructured meshes or as a constant if regular meshes are considered. There is no recommendation about refinement ratio selection. In this study refinement ratio is set as a constant value of 2 which means that the finer meshes considered are generated by dividing node spacing in all directions into halves. Order of convergence, p, can be calculated no matter what refinement ratio is selected, constant or arbitrary. For structured meshes order of convergence can be calculated using (9) and for unstructured meshes using (10):where and . If constant refinement ratio is implemented in equation for calculating order of convergence for unstructured meshes (implemented in (10)), equal value as for equation for structured mesh should be obtained. All parameters used in extrapolation and obtained results are listed in Tables 9 and 10.

TNT mass [kPa]Min. [kPa]Max. [kPa]GCIrp

115 kg1493.961493.951494.

230 kg2302.022295.472302.750.162.001.95

680 kg4571.984557.874573.550.172.002.44

Note: 1 pound = 0.45 kg; 1 psi = 6.89 kPa.

TNT mass [kPa]Min. [kPa]Max. [kPa]GCIrp

115 kg10124.9910124.7310127.

230 kg15637.2315481.9115654.490.552.001.47

680 kg32956.9432948.0732957.930.012.003.97

Note: 1 pound = 0.45 kg; 1 psi = 6.89 kPa.

Given in a percentage manner, the grid convergence index (GCI) can be considered as a relative error bound showing how the solution calculated for the finest mesh is far from the asymptotic value. It predicts how much the solution would change with further refinement of the mesh. The smaller the value of the GCI, the better. This indicates that the computed solution is within the asymptotic range. The safety factors are arbitrarily set based on the accumulated experience on computation fluid dynamics (CFD) calculations [3840]. The safety factor represents 95% confidence for the estimated error bound. This assumption can be expressed as the following statement: There is 95% confidence that the converged solution is within the range [f1(1 − GCI12/100%), f1(1 + GCI12/100%)]. The GCI is defined aswhere is a safety factor. The recommended CFD values of the safety factor are

= 3.0 when two meshes are considered,

= 1.25 for three and more meshes.

Quantity ε defines relative difference between subsequent solutions:

Due to lack of experimental investigation and reliance on AT-Blast values, Richardson extrapolation was introduced to determine convergence of the simulations. Pressure values () were obtained through numerical simulations, using different mesh sizes with a constant refinement ratio (r) of 2. With known solutions, order of convergence (p) was determined using (9), based on three subsequent pressure values. Extrapolated pressures ( and ) were calculated using pressures from two subsequent meshes ( and ), order of convergence (p), and refinement ratio (r). Table 9 lists the pressure extrapolated values for axisymmetric simulations and Table 10 for 3D simulations together with estimated error bounds, GCI, refrainment ratio, and order of convergence. Extrapolated values were obtained by using three smallest air mesh sizes, 1.56 mm, 0.78 mm, and 0.39 mm for axisymmetric and 12.5 mm, 6.25 mm, and 3.13 mm for 3D simulations, respectively, following the previously stated calculation procedure given by [41].

With the reduction of the mesh size, pressures converge to a value dependent on the charge mass, as can be seen in Figure 8. This fact should be taken with caution because smaller mesh sizes substantially prolong calculation time and demand large quantities of computational memory. Mesh sizes which produce pressure difference in comparison to extrapolated pressure less than 10% are considered sufficiently accurate for engineering use. The analysis was based on close-in and intermediate scaled distances from 0 m/kg1/3 to 1 m/kg1/3 that coincide with the previously mentioned bridge superstructure distance from the road surface. If lower bound is considered, appropriate incident pressure values are obtained for 12.5 mm air mesh size for 115 kg, 230 kg, and 680 kg of TNT, respectively.

Gauge point of interest was located in the lower plane of bridge superstructure, at a distance of 5 m from ground plane, which is equivalent to scaled distances for each TNT quantity: 1.03 m/ for 115 kg, 0.82 m/ for 230 kg, and 0.57 m/ for 680 kg. Using Richardson extrapolation and accepting 10% upper and lower bound it was possible to obtain acceptable air mesh sizes for reflected pressures (Figure 9). Analysis resulted in air mesh sizes of 12.5 mm for 115 kg TNT and 25 mm for 230 kg and 680 kg TNT being adequate for pressure-structure interaction analysis. If it is considered that it is practically impossible to determine absolutely accurate blast pressures, using experimental measurements or numerical simulation, according to this analysis AT-Blast could provide an adequate approximation of blast wave pressures for both incident and reflected pressures. AT-Blast calculated pressures are very close or are within the value of lower bound (10%) of extrapolated pressures, clearly depicted in Figures 8 and 9, which is yet another argument for using AT-Blast for quick assessment of blast pressures.

5.2. Empirical and Analytical Expressions

Empirical and analytical expressions are given in literature for estimating blast overpressure. Some of the expressions are based on the analysis of experimentally acquired pressure data while other are product of purely mathematical algorithms. Expressions are given for calculating overpressures generated from spherical or hemispherical airbursts. Based on the parameters used for numerical simulations (stand of distance of 5 m and charge weight of 115 kg, 230 kg, and 640 kg) calculation of analytical and empirical expressions for incident overpressure estimation by different authors was conducted. Calculated results were used for numerical and AT-Blast overpressure evaluation. The list of analytical expressions was given by Ullah et al. (2017) [42]. Table 11 gives some perspective on calculated values and differences of expressions with AT-Blast and FEM simulations. Based on the results given in Table 11 calculations provide a wide variety of overpressure values with differences ranging from 0.1% to 166%.

TNT mass, W115 kg230 kg640 kg

Range, R5 m

Scaled distance, Z1.0280.8160.580


Richardson extrapolation1493.962302.024571.98


Authors115 kg230 kg640 kgAT-BlastFEM

Sadovskyi (1952)1120.852063.435236.09-12.30.532.4-

Brode (1955)716.401332.803530.40-43.9-35.1-10.7-44.4-34.7-9.0

Naumyenko and Petrovskyi (1956)894.981833.905277.90-30.0-10.633.5-30.5-10.236.1

Adushkin and Korotkov (1961)46.43-74.01-678.98-96.4-103.6-117.2-96.4-103.6-117.5

Henrych and Major (1979)734.891081.542022.94-42.5-47.3-48.8-43.0-47.0-47.8

Held (1983)1891.863003.145941.2348.046.350.346.847.153.2

Kinney and Graham (1985)949.711551.233044.87-25.7-24.4-23.0-26.3-24.0-21.5

Mills (1987)1627.443221.648920.1327.357.0125.626.357.8130.0

Hopkins-Brown and Bailey (1998)849.741459.842734.52-33.5-28.9-30.8-34.0-28.5-29.5

Gelfand and Silnikov (2004)886.731423.762714.82-30.6-30.6-31.3-31.2-30.3-30.0

Bajić (2007)1799.633355.678638.9940.863.5118.539.764.4122.8



Newmark and Hansen (1961)906.121647.064138.65-29.1-19.84.7-29.7-19.36.7

Wu and Hao (2005)953.231730.864216.17-25.4-15.76.6-26.0-15.28.7

Ahmad et al. (2012)2284.054232.7510524.1178.7106.2166.277.3107.4171.4

Siddiqui and Ahmad (2007)964.421499.412876.67-24.5-26.9-27.2-25.1-26.5-25.8

Ahmad et al. (2012a)1473.352858.187603.0215.339.392.314.440.096.1

Swisdak (1994)1276.392048.123827.12-0.1-0.2-3.2-0.90.3-1.3


Note: 1 ft = 0.3048 m; 1 pound = 0.45 kg; 1 psi = 6.89 kPa.

If results from Table 11 are analysed it can be seen that for some equations calculated overpressure is higher than from AT-Blast, and in other cases overpressure is lower than from AT-Blast. The same can be said for FEM results. Exact value of overpressures (ps) for analysed charge quantities is obtained by Swisdak’s equation (see (13)) but that is because it is based on the Kingery-Bulmash equation and parameters, and AT-Blast internal equations are exactly those defined by Kingery-Bulmash:where Z is the scaled distance and A1, B1, C1, D1, E1, F1, and G1 are simplified Kingery air blast coefficients given in Table 12.

Z [m/kg1/3]A1B1C1D1E1F1G1

0.2 - 2.97.211-2.107-0.3230.1120.0690.0000.000

2.9 - 23.87.594-3.0520.4100.026-0.0130.0000.000

23.8 - 198.56.054-1.4070.0000.0000.0000.0000.000

Due to large discrepancy between overpressures calculated using different expressions and AT-Blast (and numerical) results it is hard to evaluate overall performance of AT-Blast. Expressions are mathematical equations obtained through correlation with experimental data and most of these approaches are limited by extent of overall experimental data set. The tendency of empirical or analytical expressions is that the range of the charge (scaled distance) is decreasing so is the accuracy reduced. This is due to uncertainties in measured experimental data in near field and contact explosions.

6. Validation of Numerical Models

Validation of numerical models was conducted based on comparison of numerical simulation results and experimental measurements. Purpose of this comparison is to establish the degree of error of parametrically determined acceptable air mesh size in obtaining field measured blast pressures. Results are compared, not only with measured but also with estimated pressures given in research listed in Tables 13 and 14. Researchers estimated pressures using CONWEP and/or analytical calculation. Incident pressures were measured using free-field pressure gauges positioned on some distance from explosive charge while reflected pressures were measured by pressure gauges positioned in such a way to reduce influence of material deformation and deterioration to pressure values. Sensors were placed on a rigid containment structure that provided planed support and boundary conditions, for test specimens during experiment. In experimental testing performed by [15, 20] sensors were positioned on the surface of tested specimen and they were not able to measure pressures in all tests. From eight sensors, positioned by [20] just three to four were able to successfully measure blast pressures; the rest of the data was corrupted due to sensor malfunctioning. Table 13 provides comparison of measured and simulated incident pressures, while Table 14 gives comparison of reflected pressures. Pressures were additionally estimated by AT-Blast in order to validate the software as a tool for quick pressure and impulse assessment.

Free-field Pressure [MPa]Relative error [%]
12,5 mmM/EM/FEMM/ATB





Note: 1 kg = 0.453592 kg; 1 ft = 0.3048 m; 1 psi = 6.89 kPa; M: measured; E: estimated; FEM: finite element model; ATB: AT-Blast; N/A: not available.

Reflected Pressure [MPa]Relative error [%]






Note: 1 kg = 0.453592 kg; 1 ft = 0.3048 m; 1 psi = 6.89 kPa; M: measured; E: estimated; FEM: finite element model; ATB: AT-Blast.

In Table 13 it can be seen that the numerical results overestimate measured incident pressures from 1% to 25%. From Table 14 it is impossible to draw general conclusion about reflected pressure variations; nevertheless, pressures are within acceptable error margin from which it can be concluded that numerical models provide sufficiently accurate results. If compared to estimated values by authors and AT-Blast, numerical results provide smaller errors, in previously stated range. Estimated values are generally larger than those measured which provides pressure overestimation and accounts for all possible measurement errors. This pressure overestimation can lead to conservative structural design but from the aspect of human life protection in structures of higher importance this is acceptable.

7. Conclusion

No recommendations exist for mesh size selection in blast load numerical simulations since the size of the finite element strongly depends on the blast scenario. For this reason the mesh size convergence tests were carried out for ground explosion scenarios. Tests comprised axisymmetric and 3D blast numerical simulations in which different air mesh sizes were used to compare the results with AT-Blast and extrapolated values. A 10% difference of blast parameters in relation to extrapolated values was considered sufficiently accurate for engineering use, and based on this assumption air mesh size recommendations were made for two different blast scenarios.

In both cases (free-air and reflected scenario), coarse mesh led to stretching of the pressure-time profile without a distinct arrival time of the blast front. If the rate of initial pressure rise is not very steep (close to 90°) but rather gradual, it is clear sign that the mesh size is too coarse. Another indication of inadequacy of mesh size is lack of clear secondary pressure peak (see Figure 6). Finer mesh shortened the time required for blast wave front to arrive to the structure, and if the mesh size was sufficiently small, it no longer influenced the arrival time. In the general case of free-air explosion, the largest mesh size, which is adequate to capture all the blast parameters, is 12.5 mm. For the specific blast scenario (detonation of the charge beneath the bridge structure), the largest adequate mesh sizes are 12.5 mm for 115 kg of TNT and 25 mm for 230 and 680 kg of TNT.

Results supported the initial presumption that, with the increase of distance between the detonation point and structure, pressures are decreasing but the sensitivity on the mesh size and scattering of the results are greater for smaller scaled distances. The pressure values for larger scaled distances were less sensitive to air mesh size.

It should be emphasized that this research was carried out for ideal conditions (spherical shape of charge, rigid bridge deck, and rigid ground). Further research should aim to investigate pressure variations in less ideal cases: boxlike charges, nonrigid concrete slab, compressible ground, and explosive confined in car trunk (with offset from ground level).

Data Availability

Authors state that the majority of data acquired from research and input parameters needed for numerical simulations (repeatability of simulations) is listed in tables implemented in the submitted article. If additional data is requested, authors will be more than glad to share it with interested researchers. Authors think that publicly available data contributes to research better quality; it is an opportunity for discussion and exchange of new ideas among researchers.

Conflicts of Interest

The authors declare that they have no conflicts of interest.


The research presented in this paper is a part of the research project “Behaviour of Medium and Short Length Bridges Subjected to Blast Load” supported by the Faculty of Civil Engineering Osijek, and support for this research is gratefully acknowledged.


  1. H. Draganić, Overpasses Subjected to Blast Load, Faculty of Civil Engineering Osijek, Josip Juraj Strossmayer University of Osijek, Osijek, 2014. View at: Publisher Site
  2. H. Draganić and D. Varevac, “Numerical simulation of effect of explosive action on overpasses,” Gradjevinar, vol. 69, no. 6, pp. 437–451, 2017. View at: Google Scholar
  3. J. W. Nam, J. J. Kim, S. B. Kim, N. H. Yi, and K. J. Byun, “A study on mesh size dependency of finite element blast structural analysis induced by non-uniform pressure distribution from high explosive blast wave,” KSCE Journal of Civil Engineering, vol. 12, no. 4, pp. 259–265, 2008. View at: Publisher Site | Google Scholar
  4. T. C. Chapman, T. A. Rose, and P. D. Smith, “Blast wave simulation using AUTODYN2D: a parametric study,” International Journal of Impact Engineering, vol. 16, no. 5-6, pp. 777–787, 1995. View at: Publisher Site | Google Scholar
  5. B. Luccioni, D. Ambrosini, and R. Danesi, “Blast load assessment using hydrocodes,” Engineering Structures, vol. 28, no. 12, pp. 1736–1744, 2006. View at: Publisher Site | Google Scholar
  6. Y. Shi, Z. Li, and H. Hao, “Mesh size effect in numerical simulation of blast wave propagation and interaction with structures,” Transactions of Tianjin University, vol. 14, no. 6, pp. 396–402, 2008. View at: Publisher Site | Google Scholar
  7. UFC, Structures to Resist The Effects of Accidental Explosions, 2008.
  8. R.-B. Deng and X.-L. Jin, “Numerical Simulation of Bridge Damage under Blast Loads,” WSES Transactions on Computers, vol. 8, p. 10, 2009. View at: Google Scholar
  9. W. Wilkinson, D. Cormie, J. Shin, and A. Whittaker, “Modelling close-in detonations in the near-field,” in Proceedings of the 15th International Symposium on Interaction of the Effects of Munitions with Structures, Potsdam, Germany, 2013. View at: Google Scholar
  10. J. Shin, A. S. Whittaker, D. Cormie, and W. Wilkinson, “Numerical modeling of close-in detonations of high explosives,” Engineering Structures, vol. 81, pp. 88–97, 2014. View at: Publisher Site | Google Scholar
  11. J. Shin, A. Whittaker, D. Cormie, and A. Aref, Verification and validation of a CFD code for modeling detonations of high explosives, 2015.
  12. J. Shin, Air-blast effects on civil structures, State University of New York, Buffalo, 2014.
  13. J. Shin, A. S. Whittaker, and D. Cormie, “Incident and Normally Reflected Overpressure and Impulse for Detonations of Spherical High Explosives in Free Air,” Journal of Structural Engineering (United States), vol. 141, no. 12, 2015. View at: Google Scholar
  14. J. Shin, A. Whittaker, and D. Cormie, “Estimating incident and reflected airblast parameters: updated design charts,” in Proceedings of the 16th international symposium for the interaction of the effects of munitions with structures (ISIEMS 16), 2015. View at: Google Scholar
  15. N. Yi, J. J. Kim, T. Han, Y. Cho, and J. H. Lee, “Blast-resistant characteristics of ultra-high strength concrete and reactive powder concrete,” Construction and Building Materials, vol. 28, no. 1, pp. 694–707, 2012. View at: Publisher Site | Google Scholar
  16. C. Q. Wu, L. Huang, and D. J. Oehlers, “Blast testing of aluminum foam-protected reinforced concrete slabs,” Journal of Performance of Constructed Facilities, vol. 25, no. 5, pp. 464–474, 2011. View at: Publisher Site | Google Scholar
  17. A. Ghani Razaqpur, A. Tolba, and E. Contestabile, “Blast loading response of reinforced concrete panels reinforced with externally bonded GFRP laminates,” Composites Part B: Engineering, vol. 38, no. 5-6, pp. 535–546, 2007. View at: Publisher Site | Google Scholar
  18. T. S. Lok and J. R. Xiao, “Steel-fibre-reinforced concrete panels exposed to air blast loading,” Proceedings of the Institution of Civil Engineers: Structures and Buildings, vol. 134, no. 4, pp. 319–331, 1999. View at: Publisher Site | Google Scholar
  19. J. Li, C. Wu, and H. Hao, “Residual loading capacity of ultra-high performance concrete columns after blast loads,” International Journal of Protective Structures, vol. 6, no. 4, pp. 649–669, 2015. View at: Publisher Site | Google Scholar
  20. L. Chen, Q. Fang, J. Fan, Y. Zhang, H. Hao, and J. Liu, “Responses of masonry infill walls retrofitted with CFRP, steel wire mesh and laminated bars to blast loadings,” Advances in Structural Engineering, vol. 17, no. 6, pp. 817–836, 2014. View at: Publisher Site | Google Scholar
  21. A. H. Hameed, “Dynamic Behaviour Of Reinforced Concrete Structures Subjected To External Explosion, A Thesis,” Jumada Elawal, vol. 1428, 2007. View at: Google Scholar
  22. J. S. Humphreys, “Plastic deformation of impulsively loaded straight clamped beams,” Journal of Applied Mechanics, vol. 32, no. 1, pp. 7–10, 1964. View at: Google Scholar
  23. A. Freidenberg, A. Aviram, L. K. Stewart, D. Whisler, H. Kim, and G. A. Hegemier, “Demonstration of tailored impact to achieve blast-like loading,” International Journal of Impact Engineering, vol. 71, pp. 97–105, 2014. View at: Publisher Site | Google Scholar
  24. G. C. Mays and P. D. Smith, “Blast effects on buildings. Design of buildings to optimize resistance to blast loading,” T. Telford, 1995. View at: Google Scholar
  25. N. Newmark and R. Hansen, Design of blast resistant structures, Shock and vibration Handbook, vol. 3, 1961, Design of blast resistant structures, Shock and vibration Handbook.
  26. J. Drake, L. Twisdale, R. Frank et al., “Protective Construction Design Manual: Resistance of Structural Elements (Section IX),” Report ESL-TR-87-57, Air Force Engineering and Services Center, Tyndall Air Force, Base, 1989. View at: Google Scholar
  27. P. D. Smith, J. G. Hetherington, and P. Smith, Blast and ballistic loading of structures, Butterworth-Heinemann Oxford, 1994.
  28. ANSYS, “ANSYS AUTODYN User manual,” in ANSYS Release 14.0, pp. 1–464, Canonsburg, PA, USA, 2010. View at: Google Scholar
  29. L. Kwasniewski, “Application of grid convergence index in FE computation,” Bulletin of the Polish Academy of Sciences—Technical Sciences, vol. 61, no. 1, pp. 123–128, 2013. View at: Google Scholar
  30. N. K. Bimbaum, J. Tancreto, and K. Hager, “Calculation of blast loading in the high performance magazine with AUTODYN-3D,” DTIC Document, 1994. View at: Google Scholar
  31. R. Courant, K. Friedrichs, and H. Lewy, “Über die partiellen differenzengleichungen der mathematischen physik,” Mathematische Annalen, vol. 100, no. 1, pp. 32–74, 1928. View at: Publisher Site | Google Scholar | MathSciNet
  32. G. Baudin and R. Serradeill, “Review of Jones-Wilkins-Lee equation of state,” in Proceedings of the EPJ Web of Conferences, vol. 10, 00021, 2010. View at: Google Scholar
  33. E. L. Lee, H. C. Hornig, and J. W. Kury, “Adiabatic Expansion Of High Explosive Detonation Products,” Tech. Rep. UCRL--50422, California Univ., Livermore. Lawrence Radiation Lab., 1968. View at: Publisher Site | Google Scholar
  34. ARA, A.T.-Blast, in, Applied Research Associates, Inc., Vicksburg, Mississippi, USA, 2013.
  35. ANSYS, “Remapping Tutorial Revision 4.3,” in Autodyn; Explicite Software for Nonlinear Dynamics, pp. 1–32, Canonsburg, PA, USA, 2005. View at: Google Scholar
  36. G. K. Schleyer, M. J. Lowak, M. A. Polcyn, and G. S. Langdon, “Experimental investigation of blast wall panels under shock pressure loading,” International Journal of Impact Engineering, vol. 34, no. 6, pp. 1095–1118, 2007. View at: Publisher Site | Google Scholar
  37. J.-A. Calgaro, M. Tschumi, and H. Gulvanessian, Designers' Guide to Eurocode 1: Actions on Bridges, 2010.
  38. I. Celik and O. Karatekin, “Numerical experiments on application of richardson extrapolation with nonuniform grids,” Journal of Fluids Engineering, vol. 119, no. 3, pp. 584–590, 1997. View at: Publisher Site | Google Scholar
  39. J. Slater, Examining spatial (grid) convergence, Public tutorial on CFD verification and validation, vol. MS, 86, NASA Glenn Research Centre, 2006.
  40. L. E. Schwer, “Is your mesh refined enough? Estimating discretization error using GCI,” in LS-DYNA Forum, Germany, Bamberg, 2008. View at: Google Scholar
  41. I. B. Celik, “Procedure for estimation and reporting of uncertainty due to discretization in CFD applications,” Journal of Fluids Engineering, vol. 130, no. 7, pp. 0780011–0780014, 2008. View at: Publisher Site | Google Scholar
  42. A. Ullah, F. Ahmad, H.-W. Jang, S.-W. Kim, and J.-W. Hong, “Review of analytical and empirical estimations for incident blast pressure,” KSCE Journal of Civil Engineering, pp. 2211–2225, 2017. View at: Publisher Site | Google Scholar

Copyright © 2018 Hrvoje Draganić and Damir Varevac. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

More related articles

 PDF Download Citation Citation
 Download other formatsMore
 Order printed copiesOrder

Related articles