Scaling Turbulent Combustion Fields in Explosions

Featured Application

Scaling of Explosions.

Abstract

We considered the topic of explosions from spherical high-explosive (HE) charges. We studied how the turbulent combustion fields scale. On the basis of theories of dimensional analysis by Bridgman and similarity theories of Sedov and Barenblatt, we found that all fields scaled with the explosion length scale r₀. This included the blast wave, the mean and root mean squared (RMS) profiles of thermodynamic variables, combustion variables, velocities, vorticity, and turbulent Reynolds stresses. This was a consequence of the formulation of the problem and our numerical method, which both satisfied the similarity conditions of Sedov. We performed numerical simulations of 1 g charges and 1 kg charges; the solutions were identical (within roundoff error) when plotted in scaled variables. We also explored scaling laws related to three-phase pyrotechnic explosions. We show that although the scaling formally broke down, the fireball still essentially scaled with the explosion length scale r₀. However, the discrete Lagrange particles (DLP) (phase 2) and the heterogeneous continuum model (HCM) of the DLP wakes (phase 3) did not scale with r₀, and mean and RMS profiles could differ by a factor of 10 in some regions. This was because the DLP particles and wakes introduced an additional scale that broke the similarity conditions.

1. Introduction

Similarity theory dates back to Bridgman’s work [1] summarized in Dimensional Analysis, which he published in 1922. He points out that the differential equations describing physical phenomena can be rescaled into dimensionless units, resulting in equations that are invariant over dimensionless groups. Solutions are then invariant in settings where the dimensionless groups remain constant. This is codified in the “Buckingham $\mathsf{\Pi }$ Theorem” [2]. Sedov [3] applied such concepts to mechanics in Similarity and Dimensional Methods in Mechanics (1959). His work codified and summarized similarity methods in fluid mechanics. The similarity theory was used to derive the similarity solutions of the point explosion problem first published by G. I. Taylor [4] in 1941 and later by Sedov [5] in 1946. This same similarity solution method was used by Oppenheim and Kuhl to compute all possible similarity solutions bounded by a strong shock [6], a strong detonation wave [7], and self-similar explosion waves of variable energy at the front [8]. Scaling aspects of blast waves were summarized by G. I. Barenblatt in his book Scaling, Self-Similarity, and Intermediate Asymptotics [9] published in 1996.

Here, we consider the topic of explosions from spherical high-explosive (HE) charges. We are interested in the turbulent combustion fields in their fireballs, as illustrated in Figure 1. In this work, we derive scaling laws for all fields (including the thermodynamic variables, velocities, component mass-fractions, and turbulent Reynolds stresses) on the basis of the similarity theory for explosion fields, as described by Sedov.

Relevance of This Work

One of the reasons that such scaling laws are important is because we can study the characteristics of turbulent fireballs at the laboratory scale, such as 50 g charge explosions, as reported by Glumac and Kuhl [10], then scale them up to large scale explosions, such as the 15 kg charge explosion shown in Figure 1. Many more diagnostic techniques can be applied at the laboratory scale but are impractical at the field-test scale. Such laboratory measurements include schlieren/shadow photography; holographic interferometry of particles; spectral imaging in the IR, visible, and UV wavelengths; measurements of species composition via gas chronometry; collection of carbon particles that are the source of strong Planckian radiation [10]; and optical measurements inside the fireball, as demonstrated in [10], to name just a few. In other words, using the scaling laws demonstrated in this paper, we can measure intimate details of turbulent fireballs in laboratory experiments, then confidently apply them to large-scale explosion events like that shown in Figure 1.

In contrast, the particle and mist phases of pyrotechnic explosions do not obey Sedov’s scaling laws. Thus, calculations and experiments must be performed at each size of pyrotechnic charge under consideration. Even if the initial geometric charge size is scaled, the evolving flow fields of the particles and mist do not scale. Such facts have not been recognized in the past.

2. Gas Dynamics of Explosions

The turbulent combustion fields are governed by the inviscid conservation laws of gas dynamics. These may be expressed in strong conservation form as

Mass:

$${\partial }_t\rho +\nabla ·\rho \mathit{u}=0$$

(1)

Momentum:

$${\partial }_t\rho u+\nabla ·\left(\rho \mathit{u}\mathit{u}+p\right)=0$$

(2)

Energy:

$${\partial }_t\rho E+\nabla ·\left(\rho Eu+pu\right)=0$$

(3)

Component DP:

$${\partial }_t\rho Y_D+\nabla ·\left(\rho Y_{D}u\right)=0$$

(4)

In the equations above, $\rho$ is density, u represents velocity, E denotes total energy ($E=e+0.5 u·u$), and $Y_D$ is the mass fraction of detonation product species. Two components are recognized: detonation products (DP), which serve as the fuel, and air (A). Their mass fractions obey the conservation relation:

$$Y_D+Y_A=1$$

(5)

In pure computational cells, we have either $Y_D=1$ or $Y_A=1$; in mixed cells, $Y_D \mathrm{and} Y_A<1$. We assume that the components are in a state of thermodynamic equilibrium. This allows us to compute the pressure and temperature equation of states (EOS) according to

Pressure:

$$\mathrm{p}=f_1\left(\varrho ,e, Y_D\right)$$

(6)

Temperature:

$$\mathrm{T}=f_2\left(\varrho ,e, Y_D\right)$$

(7)

These functions are evaluated by using the thermodynamic equilibrium code Cheetah [11,12]. Results have been tabulated for computational efficiency.

We follow the evolution of component mass fractions of fuel (detonation products): $Y_F$, oxidizer; (air): $Y_O$; and products (combustion products): $Y_P$. They were computed according to the following conservation laws:

Fuel ($Y_F$):

$${\partial }_t\varrho Y_F+\nabla ·U\varrho Y_F=-\dot{\phi }$$

(8)

Oxidizer ($Y_O$):

$${\partial }_t\varrho Y_O+\nabla ·U\varrho Y_O=-\alpha \dot{\phi }$$

(9)

Products ($Y_P$):

$${\partial }_t\varrho Y_P+\nabla ·U\varrho Y_P=(1+\alpha \dot{)\phi }$$

(10)

Conservation:

$$Y_F+Y_O=Y_P$$

(11)

where $\alpha$ denotes the stoichiometric oxidizer/fuel ratio. We use the fast chemistry, infinite Damkohler number limit (i.e., gas dynamic limit), where turbulent mixing within the cell occurs according to the MILES concept of Boris [13]:

Combustion Rate:

$$\dot{\mathsf{\phi }}=\{\begin{array}{c}\mathrm{burn} \mathrm{all} \mathrm{Air}: \mathsf{\Delta }{\varrho }_{\mathrm{A}}/\varrho \mathsf{\Delta }t \left(\mathrm{Fuel}-\mathrm{rich} \mathrm{cell}\right)\\ \mathrm{burn} \mathrm{all} \mathrm{Fuel}: \mathsf{\Delta }{\varrho }_{\mathrm{F}}/\varrho \mathsf{\Delta }t \left(\mathrm{Air}-\mathrm{rich} \mathrm{cell}\right)\end{array}$$

(12)

The problem is initialized by a similarity solution of Kuhl [14] for a constant velocity detonation wave propagating at the Chapman–Jouguet (CJ) speed when the detonation wave reaches the surface of the charge. See Figure 2 for the similarity profiles. (We used a perturbation on the density profile to break nonphysical numerical symmetries. We discuss HE density perturbation effects on the mean and root mean squared (RMS) profiles in the Section 6.)

The detonation products expand, doing $p·dV$ work on the surroundings, thereby creating a blast wave. The DP–air interface is unstable; perturbations grow due to vorticity creation by the inviscid baroclinic mechanism: $\dot{\omega }=\frac{1}{{\rho }^2}\nabla \rho \times \nabla p$, thereby creating a turbulent combustion layer near the edge of the fireball.

We solved the above conservation laws with an unsplit second-order Godunov method. We used the unsplit piecewise parabolic method (PPM) of Colella and Woodward [16,17] to advance the solution in time. PPM flattens slopes in cells near discontinuities to enforce monotonicity constraints (i.e., suppression of oscillations near discontinuities). This slope flattening induces a non-linear dissipation mechanism that acts on the cell level. It also reduces the scheme from second-order accurate in smooth regions of the flow, to first order near discontinuities. This is consistent with the MILES (monotone integrated large eddy simulations) concept of Jay Boris [13], whereby mixing/dissipation only occur on the smallest grid scale. We used adaptive mesh refinement (AMR) [18,19] to follow steep gradients and mixing structures.

3. Scaling Laws

In the inviscid gas dynamics formulation described above, variables such as density, velocity, pressure, and internal energy are related in ways controlled by their dimensional units. Sedov [3] and Barenblatt [9] show that chemical explosion fields are controlled by the non-dimensional independent variables of space, $x$, and time, $\tau$, according to

$$x=r/r_0$$

(13)

$$\tau =t/t_0$$

(14)

where the explosion length scale is

$${\mathrm{r}}_0={\left(m·\mathsf{\Delta }{H}_d/p_a\right)}^{1/\left(\mathrm{j}+1\right)}$$

(15)

and the explosion time scale is

$$t_0=r_0/a_a$$

(16)

In the above, m is the mass of the charge,$\mathsf{\Delta }{H}_d$ represents the heat of detonation of the explosive material, and $a_a$ denotes the sound speed of the ambient atmosphere. The variable j equals 0, 1, or 2 for planarly, cylindrically, or spherically symmetric flows. Similitude theory states [3] that the non-dimensional dependent variables of $U=\left\{\varrho , u,e, p, T\right\}$ of gas dynamics are functions of $x$ and $\tau$:

$$\mathrm{R}=\varrho /{\varrho }_a=f_1\left(x,\tau \right)$$

(17)

$$\mathrm{R}=\varrho /{\varrho }_a=f_1\left(x,\tau \right)$$

(18)

$$\epsilon =e/e_a=f_3\left(x,\tau \right)$$

(19)

$$P=p/p_a=f_4\left(x,\tau \right)$$

(20)

$$\Theta =T/T_a=f_5\left(x,\tau \right)$$

(21)

For spherical explosions in three-dimensional Cartesian co-ordinates, Equation (13) becomes

$$x=x/r_0 \mathrm{and} y=y/r_0 \mathrm{and} z=z/r_0$$

(22)

and the spherical (j = 2) explosion length scale becomes

$$r_0={\left(m·\mathsf{\Delta }{H}_d/p_a\right)}^{1/3}$$

(23)

The dependent variables then become

$$R=\varrho /{\varrho }_a=f_1\left(x,y, z, \tau \right)$$

(24)

$$U=u/a_a=f_2\left(x,y, z, \tau \right)$$

(25)

$$E=e/e_a=f_3\left(x,y, z, \tau \right)$$

(26)

$$P=p/p_a=f_4\left(x,y, z, \tau \right)$$

(27)

$$\Theta =T/T_a=f_6\left(x,y, z, \tau \right)$$

(28)

The above is the form of the solution of the inviscid conservation laws of gas dynamics for spherical explosions. Everything scales with the explosion length scale, $r_0$.

4. Results

4.1. Fireball Cross-Sections

Figure 3 depicts a cross-section of the temperature field of a fireball from a 1 g TNT (Trinitrotoluene explosive: $C_6{H}_2{\left(N{O}_2\right)}_{3}C{H}_3$) explosion in air at $100 \mathsf{\mu }\mathrm{s}/{\mathrm{g}}^{1/3}$. One can see a variety of scales in the mixing layer. Combustion in this inviscid formulation occurs in an exothermic sheet. Our numerical method satisfies same type of scaling laws as the fireball. Consequently, we expect the numerics to also be self-similar. We confirmed this observation by simulations at two scales.

Figure 4 and Figure 5 compare cross-sections of the flow fields from a 1 g charge at $t=0.5015 \mathrm{ms}/{\mathrm{g}}^{1/3}$ with flow fields from the 1 kg charge at the same scaled time, $t=5.015 \mathrm{ms}/{\mathrm{kg}}^{1/3}$. Figure 4 compares the density fields; it shows the outer shock (blast wave) and secondary expanding shock, as well as the turbulent mixing in the fireball. Figure 5a,b shows comparisons of the vorticity fields in the fireballs, showing mean fireball radii of $R_{FB}=30 \mathrm{cm}/{\mathrm{g}}^{1/3}$ and $R_{FB}=300\mathrm{cm}/{\mathrm{kg}}^{1/3}$, respectively. These are driven by the gas dynamic fields, which scale with the explosion length scale, $r_0$.

It is difficult to decern differences in the two flow fields from the above cross-section figures. Thus, we computed the maximum density difference between the two cases at each scaled time. Results are shown in Figure 6. Differences started off as roundoff errors $O~({10}^{-13})$. They grew exponentially at around $5 \mathsf{\mu }\mathrm{s}/{\mathrm{g}}^{1/3},$ and eventually saturated at values of tens of percent.

4.2. Mean and RMS Profiles

To quantitatively compare the flow field differences in Figure 4 and Figure 5, we azimuthally averaged the three dimensions fields in $\theta and \varphi$ to produce mean and RMS (root mean squared) radial profiles of the turbulent variables. The averaging method is presented in Appendix A.

Figure 7 and Figure 8 depict azimuthally averaged profiles of the thermodynamic profiles of temperature, pressure, and density from a 1 g charge at $t=0.5015 \mathrm{ms}$ (blue curves) compared with flow fields from the 1 kg charge at the same scaled time, $t=5.015 \mathrm{ms}$ (orange curves). The mean and RMS profiles from the two scales overlayed.

Figure 9 and Figure 10 present the azimuthally averaged combustion profiles of fuel (detonation products), oxidizer (air), and combustion products from the 1 g charge at $t=0.5015 ms$ (blue curves) compared with flow fields from the 1 kg charge at the same scaled time, $t=5.015 \mathrm{ms}$ (orange curves). The mean and RMS profiles from the two scales overlayed.

Figure 11 and Figure 12 depict azimuthally averaged velocity profiles ($u_r,u_{\theta }, \mathrm{and} u_{\varphi }$) from a 1 g charge at $t=0.5015 \mathrm{ms}$ (blue curves) compared with flow fields from the 1 kg charge at the same scaled time, $t=5.015 \mathrm{ms}$(orange curves). The mean and RMS profiles from the two scales overlayed. Profiles of mean values for $u_{\theta } \mathrm{and} u_{\varphi }$ were essentially zero.

Figure 13 shows azimuthally averaged off-diagonal turbulent Reynolds stresses (${\tau }_{r\theta },{\tau }_{r\varphi }, \mathrm{and} {\tau }_{\theta \varphi }$) from a 1 g charge at $t=0.5015 \mathrm{ms},$ compared with flow fields from the 1 kg charge at the same scaled time, $t=5.015 \mathrm{ms}$. The profiles were essentially zero, with random fluctuations due to the small ensemble size at small radii.

Figure 14 shows the maximum density histories for a 1 g charge (blue curve) and a 1 kg charge (orange curve); the curves overlayed. Maximum densities decayed because of the adiabatic expansion of the detonation products, with a spike at $t=0.1 \mathrm{ms}/{\mathrm{g}}^{1/3}$ due to the imploding shock [20].

Figure 15 presents the combustion histories for 1 g and (blue curve) and scaled 1 kg charge (orange curve); the curves overlayed. About 50% of the fuel mass was consumed by combustion at this time.

Let us azimuthally average the fireball temperature field in $\theta ,\varphi$ to find the evolution of the mean temperature field with time:

$${\overline{T}}_{\mathrm{FB}}\left(x,\tau \right)={{\displaystyle \iint }}^{ }T\left(x,\theta ,\varphi ,\tau \right)\mathrm{d}\mathsf{\theta }\mathrm{d}\varphi$$

(29)

Let us define the fireball radius as the radius where ${\overline{T}}_{FB}\left(x_{FB},\tau \right)=T_{max}/2$ (starting from outside the fireball and working inward). Then, the evolution of the fireball radius becomes

$$x_{FB}=\tilde{g}\left(\tau \right)$$

(30)

Dividing by the initial nondimensional charge radius $x_c=r_c/r_0$, one finds

$$R_{FB}/R_c=g\left(\tau \right)$$

(31)

From experiments and numerical simulations, when the spherical fireball has expanded to one atmosphere, one obtains

$$R_{FB}/R_c=g(\infty )\approx 30$$

(32)

The largest characteristic scale in the turbulent fireball is the fireball radius, and thus turbules in the fireball scale with $r_{FB}$.

4.3. Summary

We confirmed the similarity of turbulent combustion fields of fireballs from HE explosions. Figure 16 presents turbulent fireball structures from our gadynamic simulations (Figure 16a,c), which are compared with a photograph of a fireball (Figure 16b). They show similar turbulent structures that scaled with the explosion length scale, ${\mathrm{r}}_0$.

As expected, all gasdynamic profiles (e.g., thermodynamic profiles, combustion profiles, velocity profiles, and Reynolds stresses), both their mean and RMS values, scaled with the explosion scale $r_0={\left(m·\mathsf{\Delta }{H}_d/p_a\right)}^{1/3}$. The fireball radius scaled according to $R_{FB}/r_0=f\left(t/t_0\right)$. To sum up, all aspects of gaseous turbulent combustion in HE explosions scaled with the explosion radius $r_0$.

5. Three-Phase Pyrotechnic Explosions

Next, we considered pyrotechnic explosions. Figure 17 presents a cross-section of the three-phase temperature field of a pyrotechnic explosion from our paper in the 16th International Detonation Symposium [21]. It shows the turbulent mixing and combustion of the gas phase, which reached a temperature of 2500 (yellow regions) corresponding to the adiabatic flame temperature of the HE–air system. It also shows discrete Lagrange particles (black dots) that shed micro-mist wakes (red curves). This three-phase nature of such pyrotechnic explosions is what we explore in this section. Models are described next.

5.1. Gas Phase

The conservation laws of gas dynamics are used to model the expansion of the detonation products, mixing, and turbulent combustion of the detonation products with air. The gas phase is coupled to the particles and wakes phases with drag and heat transfer terms. The conservation laws of mass, momentum, total energy, and detonation products (DP) mass fraction are

$${\partial }_t\rho +\nabla ·\left(\rho u\right)=0$$

(33)

$${\partial }_t\rho u+\nabla ·\left(\rho uu+p\right)=+D_i+n_w{D}_w$$

(34)

$${\partial }_t\rho E+\nabla ·\left(\rho Eu+pu\right)=+\left(D_i+n_w{D}_W\right)·u+\dot{Q_i}+n_w{\dot{Q_w}}_{ }$$

(35)

$${\partial }_t\rho Y_D+\nabla ·\left(\rho Y_{D}u\right)=0$$

(36)

where $\rho , u, p, T, \mathrm{and} Y_D$ are the gas density, velocity, pressure, temperature, and mass-fraction of DP, respectively. Here, $E=e+u·u/2$ is the total energy and $\mathrm{e}$ is the internal energy. The equation of state is specified by tables $p, T, a=f_i\left(\rho ,e,Y_D\right)$ based upon the Cheetah code [11]. Drag and heat transfer effects with the particles and wakes are included in the source terms on the right-hand side of Equations (34) and (35).

5.2. Particle Phase

The dynamics of the particle phase is modeled by discrete Lagrange particles (DLPs) and their interactions with air, namely, mass, momentum, and energy [22]:

$$\frac{d{x}_i}{dt}=v_i\left(x_i\right)$$

(37)

$$\frac{d{m}_i}{dt}=-\dot{s_i}\left(x_i\right)$$

(38)

$$m_i\frac{d{v}_i}{dt}=-D_i\left(x_i\right)-\dot{s_i}{v}_i+m_i\mathrm{g}$$

(39)

$$m_i\frac{d{e}_i}{dt}=-\dot{Q_i}\left(x_i\right)-\dot{s_i}{e}_i$$

(40)

where $x_i,{ v}_i, m_i, \mathrm{and} e_i$ denote the position, velocity, mass, and internal energy of particle i, and $\mathrm{g}$ denotes gravity. Mass loss, drag, and heat transfer interactions with air are specified by the following relations:

$${\dot{s}}_i=s{d}_i^2$$

(41)

$$D_i=\left(1/8\right)\mathsf{\pi }{d}_i^2\rho \left(u-v_i\right)\left|u-v_i\right|C_D \mathrm{for} m_i$$

(42)

$$\dot{Q_i}=\pi d_{i}k\left(T_i-T\right)Nu \mathrm{for} m_i$$

(43)

where k represents thermal conductivity [k] = joules/(sec-m-K). The surface recession rate, $\dot{s_{\mathrm{i}}}$, is based on measurements [23] and empirical modeling [24,25] of Al–air combustion systems. The following drag and Nusselt correlations [26] (Appendix C) are assumed:

$$C_D=0.48+28R{e}^{-0.85}$$

(44)

$$Nu=2+1.4P{r}^{1/3}R{e}^{1/5}+0.13P{r}^{1/3}R{e}^{0.7}$$

(45)

5.3. Wake Phase

The wake phase is modeled by the heterogeneous continuum (HC) model of Nigmatulin [27] and Collins et al. [28]:

$${\partial }_t\sigma +\nabla ·\left(\sigma v\right)={\dot{s}}_i$$

(46)

$${\partial }_t\sigma v+\nabla ·\left(\sigma vv\right)={\dot{s}}_i{v}_i-n_w{D}_w$$

(47)

$${\partial }_t\sigma E_w+\nabla ·\left(\sigma E_{W}v\right)={\dot{s}}_i{e}_i-n_w\dot{Q_w}$$

(48)

Drag and heat transfer interactions between wakes and gas are modeled by

$$D_w=\left(1/8\right)\pi d_w^2\rho \left(u-v\right)\left|u-v\right|C_D$$

(49)

$$\dot{Q_w}=\pi d_{w}k\left(T_w-T\right)Nu$$

(50)

$$n_w=\sigma /m_w \mathrm{with} m_w=\left(1/6\right)\pi d_w^3{\rho }_w$$

(51)

The same drag and Nusselt correlations of Equations (44) and (45) are used.

5.4. Scaling of Pyrotechnic Explosions

The scaling laws for the gas phase were discussed in Section 2. For pyrotechnic explosions, we must deal with all three phases. One can say that the conservation laws for phases 2 and 3 describe inertial flow, since there are no pressure gradient forces. However, there are source/sink terms on the right-hand side of the equations, expressing mass transfer, drag, and heat transfer with the gas phase. They are functions of the Reynolds and Prandtl numbers. How these scale with change in problem scale is discussed below.

The initial conditions for a pyrotechnic charge are illustrated in Figure 2. It shows the detonation wave structure for the gas phase (as before) and the linear velocity profile for the droplet phase. This corresponds to similarity solution of inertial flow expansion into a vacuum, as published by Stanyukovich [15] in 1960. The ball radius of the droplet phase extends to 10% of the charge radius (i.e.,$R_B=0.1·R_c$). These initial conditions scale with the explosion scale, $r_0$, in other words $x_B=R_B/r_0=0.1·R_c/r_0$. Thus, we consider two pyrotechnic explosion scales:

•

A full-scale charge (the dynamics of this pyrotechnic charge system was published in the 16th Detonation Symposium [21]) where $R_c\equiv 100 \mathrm{cm}$ with a mass of 7175 kg;

•

1/3 scale charge where $R_c\equiv 100/3 \mathrm{cm}$ with a mass of 265.7 kg.

As mentioned above [21], the pyrotechnic charge contains a liquid ball of radius $R_B=0.1·R_c$, and thus a full-scale charge would have a ball radius of $R_B=10 \mathrm{cm}$, while a 1/3-scale charge would have a ball radius of $R_B=10/3 \mathrm{cm}$. The ball is made up of droplets; we estimated the droplet diameter as $d_i=3 \mathrm{mm}$, which will be the same for both scales. (We assumed the liquid sphere was put under strong radial strain, creating droplets. According to Tolman (1949) [29], the droplet size is controlled by surface tension, which is a material property. Thus, the droplet size will be the same for full scale and subscale.) The number of droplets scaled as $n_i={\left(2R_B/d_i\right)}^3$, which was equal to 296,296 droplets at full scale and 10,974 at 1/3 scale. The droplet wake diameter scaled as $d_w=d_i/100$, and thus was $d_w=3 \mathsf{\mu }\mathrm{m}$ in both cases. The wake equilibration length scaled as $l_w/d_w=10,000$, and thus the wake equilibration length was 30 cm in both cases. The droplet Reynolds number scaled as $R{e}_i=u·d_i/\upsilon$. Assuming a characteristic velocity of 1 km/s, one finds the droplet Reynolds number of $R{e}_i=2\times {10}^5$, and the wake Reynolds number of $R{e}_W=2\times {10}^3$. These have drag coefficients of $C_D=0.48$ and Nusselt numbers of $Nu~ 40$(see Appendix C for more Nusselt number correlations versus data). These results are summarized in

Table 1. Comparison of radii, drag, and heat transfer for pyrotechnic droplet explosions (d_i = 3 mm).

Variable	Scaling	Scale = 1	Scale = 1/3
———Initial Conditions———
HE charge radius	$R_c$	100 cm	100/3 cm
HE charge mass (Comp B)	$M_c={\rho }_{0}4\pi R_c^3/3$	$~ 7175$ kg	$~265.7$ kg
Liquid ball radius; diameter	$R_B=0.1·R_c; D_B=2R_B$	10 cm; 20 cm	10/3; 20/3 cm
Droplet diameter	$d_i=3 \mathrm{mm}$	3 mm	3 mm
Number of droplets	$n_i={\left(2R_B/d_i\right)}^3$	$296,296$	$10,974$
Droplets per unit volume $,{ \mathsf{\eta }}_{\mathrm{i}}$	${\eta }_i=6·{(2R_B/d_i)}^3/\pi {\left(2R_B\right)}^3$	37/cc	37/cc
Droplet equilibration length, $l_i$	$l_i/d_i=10,000$	30 m	$30$ m
Droplet wake diameter, $d_w$	$d_w=d_i/100$	$3 \mathsf{\mu }\mathrm{m}$	$3 \mathsf{\mu }\mathrm{m}$
Wake equilibration length, $l_w$	$l_w/d_w=10,000$	$30 \mathrm{cm}$	$30 \mathrm{cm}$
——— Drag & Heat Transfer Effects ———
Droplet Reynolds number *	$R{e}_i=u·d_i/\upsilon$	$2\times {10}^5$	$2\times {10}^5$
Droplet drag coefficient	$C_D=0.48+28/R{e}_i{}^{0.85}$	0.481	0.481
Droplet Nusselt number [23]	$Nu=f\left(Re,Pr\right)$	40	40
Wake Reynolds number	$R{e}_w=u·d_w/\upsilon$	$2\times {10}^3$	$2\times {10}^3$
Wake drag coefficient	$C_D=0.48+28/R{e}_w{}^{0.85}$	0.48	0.48
Wake Nusselt number [23]	$Nu=f\left(Re,Pr\right)$	40	40

* characteristic velocity u = 1 km/s is assumed.

. We note that here the particle sizes introduced an additional length scale, and thus we did not expect to maintain self-similarity.

5.5. Gas Phase Results

We investigated the scaling behavior by comparing simulations at two different scales. Pyrotechnic fireball temperature cross-sections at 1/3 scale and full scale are presented in Figure 18. The azimuthally averaged pyrotechnic fireball temperature profiles are presented in Figure 19. Curves of the mean and RMS profiles are shown; profiles at 1/3 scale (blue curves) and full scale (orange curves) overlayed, i.e., they scaled with ${\mathrm{r}}_0$. It is interesting to note that the maximum mean fireball temperature in the pyrotechnic explosion was $\overline{T}=1650 \mathrm{K}$, which was lower than the HE fireball case of $\overline{T}=2450 \mathrm{K}$ (see Figure 7a) because of heat losses to the DLPs.

Pyrotechnic fireball vorticity cross-sections are presented in Figure 20. Their azimuthally averaged mean and RMS profiles are depicted in Figure 21; results for 1/3 scale (blue curves) and full scale (orange curves) overlayed, i.e., the scale with $r_0$.

Pyrotechnic fireball combustion product cross-sections are presented in Figure 22 and Figure 23. Their azimuthally averaged mean and RMS combustion product profiles are shown. Results for 1/3 scale (blue curves) and full scale (orange curves) overlayed, i.e., the scale with $r_0$.

5.6. DLP Phase Velocity Results

Azimuthally averaged mean and RMS velocity profiles for the DLP phase of pyrotechnic explosions are presented in Figure 24 and Figure 25. The mean and RMS velocity profiles of the DLP phase did not overlay and thus the DLP phase results did not scale with ${\mathrm{r}}_0$.

5.7. Mist Phase Results

Pyrotechnic fireball mist temperature cross-sections are presented in Figure 26. Their azimuthally averaged mean and RMS mist temperature profiles are shown in Figure 27. Results for 1/3 scale (blue curves) and full scale (orange curves) did not overlay, i.e., their profiles did not scale with $r_0$.

Pyrotechnic fireball mist density cross-sections are presented in Figure 28. Their azimuthally averaged mean and RMS mist temperature profiles are shown in Figure 29. Results for 1/3 scale (blue curves) and full scale (orange curves) did not overlay, i.e., their profiles did not scale with $r_0$.

The azimuthally averaged mean and RMS mist velocity profiles are shown in Figure 30 and Figure 31. Results for 1/3 scale (blue curves) and full scale (orange curves) did not overlay, i.e., their profiles did not scale with $r_0$.

5.8. Summary

A summary of the scaling laws for pyrotechnic explosions is given in

Table 2. Scaling summary.

Variable	Scaled
HE charge radius, $R_C$	$x_C=R_C/r_0$
Ball radius, $R_B$	$x_B=R_B/r_0=R_C/3r_0=x_C/3$
Droplet diameter, $d_i$ *	$d_i=3 \mathrm{mm}=\mathrm{constant}$
Wake mist diameter, $d_w$	$d_w=d_i/100=30 \mathsf{\mu }\mathrm{m}=\mathrm{constant}$
Number droplets, $n_D$	$n_D={(2R_B/d_i)}^3$
Droplets per unit volume, ${\eta }_i$	${\eta }_i=6·{(2R_B/d_i)}^3/\mathsf{\pi }{\left(2R_B\right)}^3=\mathrm{constant}$
Droplet Reynolds number, $R{e}_i$	$R{e}_i=d_i·u/\upsilon =2\times {10}^5=\mathrm{constant}$
Wake Reynolds number, $R{e}_w$	$R{e}_w=d_w·u/\upsilon =2\times {10}^3=\mathrm{constant}$

* $d_i$ is controlled by surface tension; this value is estimated from experiments.

. When implemented in the three-phase model (Equations (33)–(51)), one finds that

•

Turbulent combustion in the gas phase fireball does scale with $r_0$;

•

The DLP flow field does not scale with $r_0$;

The mist flow field does not scale with $r_0$.

5.9. Conclusions

The gas phase of pyrotechnic explosions scaled with ${\mathrm{r}}_0$, indicating that the multi-phase effects did not disrupt the scaling behavior. The DLP particle phase behaved ballistically, being only subject to drag and heat transfer without pressure (blast wave) effects, and thus did not satisfy similarity conditions. The mist phase started with zero mass and was only created by a mass sink in the DLP equations. It also was inertial (i.e., devoid of pressure effects) and only subject to drag and heat transfer effects, and thus it also did not satisfy the similitude conditions.

6. Discussion

6.1. Initial HE Density Perturbations

When initializing the HE density profile shown in Figure 2, we introduced spatial sinusoidal variations with 1% maximum amplitude to break nonphysical numerical symmetries. This helped trigger the instabilities on the HE–air interface in order to make it less dependent on grid effects. To check whether these spatial perturbations might have any effect on scaling, we changed the sinusoidal base, and re-ran the HE fireball case shown in Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14 and Figure 15. The result was that the profiles of the two cases overlayed, and thus scaled with $r_0.$

6.2. Grid Refinement

We also explored the effect of grid refinement on scaling. We doubled the resolution and re-ran the HE fireball case shown in Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11, Figure 12, Figure 13, Figure 14 and Figure 15.

Table 3. Summary of grid resolution on mean and RMS profiles.

Variable	mean	RMS
$T$	OK *	OK
$\rho$	OK	OK
$m_F\left(t\right)$	OK	--
$\mathit{\omega }$	~ OK $(\mathsf{\delta }=9\%$)	~ OK ($\mathsf{\delta }=30\%)$
$\mathsf{\Delta }$	OK	OK
$v_r$	OK	OK
$v_{\theta }$	OK	~ OK ($\mathsf{\delta }=30\%)$
$v_{\varphi }$	OK	~ OK ($\mathsf{\delta }=30\%)$

* OK $\equiv$ scales with $r_0$, $\delta \equiv \mathrm{local} \mathrm{solution} \mathrm{variation}, \mathrm{resolution} 1 \mathrm{vs}. 2$.

summarizes the profile comparisons of 1× versus 2× zoning. We found that a factor of 2 increase in grid resolution had no effect on the mean and RMS profiles of thermodynamic variables, $T, \rho , m_F\left(t\right)$, and mean velocity-related profiles of $\omega ,v_r,v_{\theta },v_{\varphi }$. The RMS profiles also scaled with $r_0$, but had about 30% variation due to grid resolution.

7. Conclusions

We investigated scaling behavior for turbulent combustion fields in explosions. Sedov has derived a similarity theory for explosions, stating that the flow fields $U\left(x,y,z, \tau \right)$ scale as $x=r/r_0$ and $\tau =t/t_0$, where $r_0={\left(m·\mathsf{\Delta }{H}_d/p_a\right)}^{1/\left(j+1\right)}$ and $t_0=r_0/a_a$. This theory was derived for the blast wave (dilatational) velocity field. However, according to the Helmholtz decomposition theorem [30,31,32], the velocity vector field $\mathit{u}$ can be divided into its dilatational $u_{\mathsf{\Delta }}$ and rotational $u_{\omega }$ components (see Appendix B). The dilatational component corresponds to the blast wave solution of Sedov. The rotational component corresponds to the turbulent velocity field found in the fireball of explosions. Our algorithm numerically satisfies the same scaling behavior discretely, and thus we expect the numerical solution for inviscid gaseous explosions to show the same scaling relation. The two solutions were identical when displayed in scaled variables (i.e., the flow field profiles only differed by machine roundoff).

Numerical simulations of pyrotechnic fireballs performed at two different scales confirmed this observation; their azimuthally averaged scaled profiles overlayed, proving that they scaled with $r_0$. However, the DLP flow fields and the wake flow fields of pyrotechnic explosions did not scale with $r_0$ because particle sizes introduced another length scale that broke similarity. For example, scaled mist temperature, density, and velocity profiles can differ by a factor of 10 in certain regions. Thus, calculations and numerical simulations of pyrotechnic explosions must be performed at each scale of interest. This has not been recognized in the past.