The response of high explosives (HEs), due to mechanical and/or thermal insults, is of great importance for both safety and performance. A major component of how an HE responds to these stimuli stems from its reactant equation of state (EOS). Here, the tri-amino-tri-nitro-benzene based explosive PBX 9502 is investigated by examining recent experiments. Furthermore, a complete thermal EOS is calibrated based on the functional form devised by Wescott, Stewart, and Davis [J. Appl. Phys. 98, 053514 (2005)]. It is found, by comparing to earlier calibrations, that a variety of thermodynamic data are needed to sufficiently constrain the EOS response over a wide range of thermodynamic state space. Included in the calibration presented here is the specific heat as a function of temperature, isobaric thermal expansion, and shock Hugoniot response. As validation of the resulting model, isothermal compression and isentropic compression are compared with recent experiments.

The thermodynamic response of high explosives (HEs) is of paramount importance for determining the response to mechanical and thermal stimuli. During impact of projectiles/fragments, the pressures and temperatures achieved in the HE are determined solely by its reactant equation of state (EOS). Likewise, the EOS is responsible for determining the expansion/contraction during heating/cooling of the HE. This expansion/contraction leads to changes in the initial energy density of the HE, and thus the Chapman-Jouguet detonation velocity upon detonation. Other observable quantities such as the shock-polar matching between the HE and surrounding inert material rely on knowing the HE reactant EOS,2–5 as well as the process of shock transmission across the interface between the HE and a second material. The thermodynamic response is critical for future, higher fidelity, reactive flow modeling that includes effects due to changes of the initial HE temperature6,7 and in pseudo-reaction zone modeling of HE.8 

In a series of papers,1,9,10 a complete thermodynamic model of reactant EOS was developed. In those papers, a procedure for calibration and an initial set of the relevant thermodynamic parameters for PBX 9502 [95% tri-amino-tri-nitro-benzene (TATB), 5% polymeric binder Kel-F 800] were determined utilizing the data available at the time.1 The focus of the present study is to revisit the role of each of the model parameters to obtain a wider range of thermodynamic validity by incorporating experimental data not included in the above described calibration.

The outline of the remainder of the paper is as follows: Sec. II details the analytical model reactant equation of state along with a few newly derived relationships between model parameters, thermodynamic variables, and measurable quantities observable from experiments. This is then followed by Sec. III, which outlines a calibration procedure for systematically determining EOS parameters. Some derived properties of the model are presented in Sec. IV. Section V presents validation tests of the calibrated model via isentropic and isothermal compression experiments. Finally, Sec. VI gives conclusions and possible improvements for future work.

A thermally complete Mie-Grüneisen equation of state based on a reference isentrope has been devised, and has been utilized in reactive flow models of HE.1 Following Wescott et al., for a reactant EOS, the relations between specific internal energy, er, pressure, pr, and density, ρ, for the reactants (thus the subscript r) are given by the standard Mie-Grüneisen form

er(p,ρ)=ers(ρ)+pprs(ρ)ρΓr(ρ)
(1)

and

pr(e,ρ)=prs(ρ)+ρΓr(ρ)(eers(ρ)),
(2)

where Γr is the Grüneisen gamma, here taken to be only a function of ρ. The reference isentrope, denoted by a superscript s, of the reactants is given by the following piecewise functional form:

prs(ρ)={p̂[exp(4By)1],ρ<ρ0p̂[j=13(4By)jj!+C(4By)44!+y2(1y)4],ρρ0,
(3)

where y = 1 – ρ0/ρ and p̂=ρ0A2/4B. Here, A, B, and C are material dependent parameters and ρ0 is the reference ambient density. Note that the above two functional forms are continuous up to y3 at ρ = ρ0. Furthermore, the above prs(ρ) form is monotonically increasing with respect to ρ, assuming that the parameters A, B, and C are positive. Thus, along the principal isentrope dprsdρ>0, and a real sound speed will result, even when the material is distended. Note that even though the isentropic bulk modulus is always positive, the isothermal bulk modulus becomes negative at sufficiently low enough density, which is a common issue with solid EOS models.11 If no switch is employed in Eq. (3) (and C > 0), and the second term is used for all ρ values, then there would always come a point for low enough density where dprsdρ<0, implying an imaginary sound speed. Given this fact, the switch in the above equation is required to maintain reasonable mathematical behavior when ρ < ρ0. An example plot of both functions is shown in Fig. 1.

FIG. 1.

Isentropes from Eq. (3) with parameters from Table I. The red curve is for ρ < ρ0, blue curve for ρ ≥ ρ0. Note the non-monotonic behavior of the blue curve for ρ < ρ0.

FIG. 1.

Isentropes from Eq. (3) with parameters from Table I. The red curve is for ρ < ρ0, blue curve for ρ ≥ ρ0. Note the non-monotonic behavior of the blue curve for ρ < ρ0.

Close modal

Two other properties of the above fitting form are worth noting. First, by taking the derivative of prs(ρ) with respect to ρ, it is easy to verify that the parameter A is the ambient sound speed of the material at p = 0 and ρ = ρ0. The isentropic bulk modulus is defined as

K=ρdprsdρ,
(4)

along the principal isentrope. The bulk sound speed is related to the ambient isentropic bulk modulus, K0, via K0 = ρ0A2. Furthermore, by taking the derivative of the isentropic bulk modulus with respect to pressure along the principal isentrope, the following is obtained:

K=dKdprs=dKdρdρdprs=dKdρ/dprsdρ.
(5)

Evaluating this function at ρ = ρ0, it is found that

K0=dKdprs=4B+12B1.
(6)

The slope, s, of the Hugoniot in shock velocity, Us, and particle velocity, Up, space at Up = 0 is given by11 

s=14(K0+1)=B+18B.
(7)

Thus, when B ≫ 1, the initial slope of the Hugoniot in UsUp space is s  B. The C term, in Eq. (3), was added to modify the Hugoniot at high pressures, to ensure that the reactant Hugoniot does not cross the product Hugoniot.1 

From the thermodynamic identity de = Tdspdv, for an isentropic process, simply de = –pdv is attained. Upon substitution of specific volume, v, in favor of density, the following is obtained: de=pρ2dρ. The above reference pressure along the isentrope can then be integrated to obtain the reference energy along the isentrope

ers(ρ)=ρ0ρprsρ¯2dρ¯+E0,
(8)

where E0 is the integration constant. This integration constant is often taken to be the stored chemical potential energy of the explosive, but simply represents an offset in examining just the reactant EOS. For inert solid EOS, this integration constant is often taken to be zero, yielding zero energy at the reference state.

The Grüneisen parameter is taken to be

Γr(ρ)={Γr0,ρ<ρ0Γr0+Zy,ρρ0,
(9)

where Z is a constant used to describe the changes to Γr with respect to density. Note that the ρ < ρ0 branch is not explicitly mentioned in the original Ref. 1, but is often implemented in practice. Note also that y = 0 at ρ = ρ0, y → 1 as ρ and y → – as ρ → 0. Considering the fact that the sound speed11 contains derivatives of Γr(ρ), there will be a discontinuity in the sound speed at ρ = ρ0 for nonzero pressures and Z ≠ 0.

The reactant temperature, Tr, after a series of manipulations,1 can be obtained as a function of energy and density

Tr(e,ρ)=Trs(ρ)(1+αstCvr0Trs(ρ)[eers(ρ)]+1)11+αst,
(10)

where Trs(ρ) is the temperature along the reference isentrope

Trs(ρ)={T0(ρρ0)Γr0,ρ<ρ0T0exp(Zy)(ρρ0)(Γr0+Z),ρρ0.
(11)

Here, Cvr0 is the reactant specific heat at constant volume at the reference temperature, T0. The parameter αst determines how the specific heat changes with respect to temperature. In particular, Eq. (10) can be inverted to determine e(ρ, T), then holding ρ fixed and taking the derivative with respect to T yields the specific heat at constant volume

Cvr(ρ,T)=Cvr0(TTrs(ρ))αst.
(12)

Note that at the reference density, Cvr(ρ0,T)=Cvr0(TT0)αst. This functional form for the specific heat at constant volume is a power law, and is not fundamental, but rather a consequence of choices made in assuming the specific heat at constant volume was a linear function of entropy. Other choices, in general, can be made.12 

It is instructive to obtain pressure as an explicit function of density and temperature. Upon substituting energy in favor of pressure from Eq. (1) into Eq. (10), and solving for pressure, one obtains

pr(ρ,T)=prs(ρ)+ρΓr(ρ)Cvr0Trs(ρ)1+αst((TTrs(ρ))1+αst1).
(13)

From this formula, isochoric pressurization can be obtained while keeping the density fixed. This is the pressurization due to heating at a fixed volume.

Similarly, the pressure along a specified temperature may be obtained by varying the density. This is often accomplished in diamond anvil cells for HE constituents, and is discussed in Sec. V B.

For the reactants EOS, the following parameters are required to define the model: A, B, C, ρ0, Z, Γr0, E0, Cvr0, αst, and T0. Since E0 is associated with the stored chemical energy of the explosive, it is not investigated here. The specified initial pressing density range of PBX 9502 (Ref. 13) is from 1.885 g/cm3 to 1.895 g/cm3, so a nominal initial density will be taken as ρ0 = 1.890 g/cm3. Here, the reference temperature is taken to be T0 = 297 K.

It is difficult to directly measure the specific heat of a condensed phase material at constant volume over a wide temperature range. This is especially true for high explosives14 that decompose at temperatures from 465 to 600 K, yet there is interest in knowing thermal properties such as the shock temperatures under detonation conditions. These temperatures may reach 2000 K under shock loading, even prior to decomposition into detonation products. Considering that most high explosive molecules are rather large, and thus have many vibrational modes, the specific heat is expected to change significantly over these large temperature ranges. So, the best that can be done is to model15 the specific vibrational modes and frequencies16 to build a thermal model of Cvr. Here, the parameters Cvr0 and αst in Eq. (12) are fit to the model presented by Menikoff.15 The resulting parameters, determined by least squares fitting over 0 < T < 1800 K, are presented in Table I. A comparison between Menikoff's determined Cvr, the newly calibrated model, and the original model1 is shown in Fig. 2. Note that the original calibration matches Menikoff's work at relatively low temperatures, but the current calibration is reasonable over a wider range of temperatures; all models agree near ambient conditions. It should be pointed out that the functional form in Eq. (12) does not lend itself to simultaneously having the proper asymptotic behavior of dCvrdT=0 in the limits of T → 0 and T. The best fit of Eq. (12) typically will overshoot a Debye model at both high and low temperatures, while being too low at some intermediate temperature range. This is seen in Fig. 2.

TABLE I.

Davis reactant EOS parameters for PBX 9502.

ParameterValueUnit
A 1.80 mm/μ
B 4.6  
C 0.34  
ρ0 1.890 g/cm3 
Z 0.0  
Γr0 0.56  
Cvr0 0.001074 kJ/g K 
αST 0.4265  
T0 297 
ParameterValueUnit
A 1.80 mm/μ
B 4.6  
C 0.34  
ρ0 1.890 g/cm3 
Z 0.0  
Γr0 0.56  
Cvr0 0.001074 kJ/g K 
αST 0.4265  
T0 297 
FIG. 2.

Specific heat at constant volume, ρ = ρ0, versus temperature. Symbols are from Menikoff.15 The blue curve is the original calibration1 and the red curve is the current calibration.

FIG. 2.

Specific heat at constant volume, ρ = ρ0, versus temperature. Symbols are from Menikoff.15 The blue curve is the original calibration1 and the red curve is the current calibration.

Close modal

The Hugoniot of many HE, including PBX 9502,17,18 has nonlinear relations between the shock velocity, Us and particle velocity, Up. There is also a fair amount of uncertainty in the Hugoniot as shock velocities increase to the point where chemical reactions begin to take place. Furthermore, plastic bonded explosives are heterogenous and somewhat porous, which further complicates continuum modeling. As shown earlier, the parameters A and B can be related to the sound speed and the initial slope of the Hugoniot curve in UsUp space, respectively. From previous studies,1,17,19 there is a range of sound speeds from 1.75 mm/μs to above 2 mm/μs. Here, an intermediate value of A = 1.8 mm/μs was chosen. Likewise, the initial slope of the Hugoniot curve varies significantly from 3.15 to 5.2 in the references. Again, an intermediate value of B = 4.6 was chosen. Note that even in isothermal compression experiments, performed in diamond anvil cells, high uncertainty in the isothermal bulk modulus (30% variability) and its derivative with respect to pressure (100% variability)20 are observed.

The isobaric volumetric compression/expansion as a function of temperature is governed by a variety of thermodynamic parameters. Once the specific heat parameters, Cvr0 and αst, and the principal isentrope parameters, A and B, have been set, only the Grüneisen parameter is left to determine how the density varies with temperature at zero pressure (the C parameter plays a minor role, and only at the cold temperature extremes). Namely, the thermodynamic relation

Γr0=βA2Cp,
(14)

can be used to estimate thermal expansion and contraction. Here, β is the volumetric thermal expansion coefficient and Cp is the specific heat at constant pressure, which is typically only a few percent higher than Cvr for condensed phase materials. To go beyond estimation, one can compute the isobaric thermal expansion analytically by noticing that the term (eers(ρ)) in Eq. (10) becomes simply prs(ρ)/ρΓr(ρ) when p = 0, from examining Eq. (1). This leads to a rather simple piecewise expression (depending on whether in contraction or expansion) for the temperature as a function of density. This can be inverted to obtain density as a function of temperature, which can be compared to experimental observations. Equivalently, one can set the pressure on the left side of Eq. (2) to zero, and solve for T(ρ).

There have been several measurements of thermal expansion/contraction of PBX 9502 in the literature.21–26 TATB, the main component of PBX 9502, is quite anisotropic in its thermal expansion, which can lead to uncertainty in the bulk thermal expansion and contraction. And although the HE molding powder has random initial crystal orientation, once the material is pressed into PBX 9502, the crystals tend to reorient depending on the applied pressure direction during pressing. This leads to pressed parts which are then anisotropic. Another contributing factor to the uncertainty of the experimental measurements is the fact that PBX 9502 can experience “ratchet growth” upon cyclic heating and cooling,24 which leads to hysteresis. These complications are obviously beyond the scope of this fluid-like continuum EOS modeling, but should be kept in mind as possible sources of uncertainty and are worthy of future research.

To match the experimental thermal expansion, Γr0=0.56, is determined. The variation of the Grüneisen parameter with density is not used in this calibration, and so Z = 0 was taken. Taking Z = 0 also avoids the discontinuity of sound speed at the reference density. The isobaric thermal expansion for the original calibration, current calibration, and experimental data are shown in Fig. 3. For temperatures below 350 K, the current and original calibrations agree within about 1% in density. But, as the temperature increases beyond 350 K, significant differences between the original calibration and the experimental data are seen. For example, at 525 K, the measured data22 indicate a density of 1.7 g/cm3, whereas the original calibration1 is below 0.2 g/cm3 at this temperature. These differences will significantly affect the HE energy density, which will in turn affect the detonation velocity and pressure. Likewise, the Hugoniot response due to flyer impact will drastically be affected by such large density inconsistencies. Also, for safety modeling, the isochoric pressurization due to confined heating of HE charges will be rather different between these two models.

FIG. 3.

Isobaric thermal expansion. Experimental data are given by black symbols. Blue curve is the original calibration1 and the red fit is the current calibration.

FIG. 3.

Isobaric thermal expansion. Experimental data are given by black symbols. Blue curve is the original calibration1 and the red fit is the current calibration.

Close modal

The last parameter required for the EOS model is C. Once the other parameters are set, C determines the behavior at high pressures. In particular, it was originally introduced to keep the reactant Hugoniot from crossing the product Hugoniot in overdriven detonations; a crossing of Hugoniots implies that the products are more dense than the reactants, which is not expected in slightly overdriven detonations. From examination of the available Hugoniot data,17,18,27 a value of C = 0.34 is determined to adequately fit the observed Hugoniot. Figure 4 shows the original model, current model, and available experimentally measured Hugoniot points. Although the original and current models have quite different specific heats and thermal expansion properties, the Hugoniots are rather similar. This indicates that fitting to Hugoniot data, although very important, does not adequately constrain the complete thermodynamic EOS.

FIG. 4.

Ambient, T = T0, Hugoniot in the UsUp plane. The blue curve is the original calibration1 and the red curve is the current calibration.

FIG. 4.

Ambient, T = T0, Hugoniot in the UsUp plane. The blue curve is the original calibration1 and the red curve is the current calibration.

Close modal

The calibrated set of parameters for PBX 9502 are summarized in Table I.

It is worth exploring other thermodynamic properties of the current model, so as to compare with experimental data as it becomes available for validation.

From the isobaric thermal expansion, the densities at a variety of initial temperatures are obtained and listed in Table II. The reason for listing these specific temperatures is due to the variety of dynamic experiments performed at these temperatures.21–23,28–31 For three initial temperatures, the UsUp relations are shown in Fig. 5. As expected, the colder material has a stiffer response than the warmer material. The differences between the ambient and cold material in the UsUp space are hard to measure experimentally,28 and the results here are not in contradiction with the available data.

TABLE II.

Density versus temperature at p = 0.

Temperature (K)Temperature (°C)Density (g/cm3)
77 –196.15 1.942 
218.15 –55 1.914 
297 23.85 1.890 
348.15 75 1.869 
403.15 130 1.840 
525 251.85 1.697 
Temperature (K)Temperature (°C)Density (g/cm3)
77 –196.15 1.942 
218.15 –55 1.914 
297 23.85 1.890 
348.15 75 1.869 
403.15 130 1.840 
525 251.85 1.697 
FIG. 5.

Hugoniot in the UsUp plane. The blue curve is for material initially at T = 218.15 K, green curve is for T = 297 K and the red curve is for T = 348.15 K.

FIG. 5.

Hugoniot in the UsUp plane. The blue curve is for material initially at T = 218.15 K, green curve is for T = 297 K and the red curve is for T = 348.15 K.

Close modal

It is also worth examining the temperature along the Hugoniots for variations in the initial temperature. It is known that the HE sensitivity to shock depends quite strongly on the HE's initial temperature.22,28–31 Figure 6 shows the predicted shock temperature along the Hugoniot as a function of pressure for material at three different initial temperatures. It is observed that the temperature variation, between the initially cold material and hot material, increases with shock pressure. Future experimental efforts are being directed to measure the temperature of shocked PBX 9502.32 

FIG. 6.

Shock temperature as a function of shock pressure. The blue curve is for material initially at T = 218.15 K, green curve is for T = 297 K and the red curve is for T = 348.15 K.

FIG. 6.

Shock temperature as a function of shock pressure. The blue curve is for material initially at T = 218.15 K, green curve is for T = 297 K and the red curve is for T = 348.15 K.

Close modal

Another interesting property of the model is to compute the isochoric pressure increase due to thermal heating. HEs generally expand thermally much more than metals. As such, for the HE encased in metal without sufficient ullage, it is an important safety characteristic to know how the HE will pressurize internally. A rough approximation is to examine the isochoric pressurization of the HE from Eq. (2). This is shown in Fig. 7 for the current model. As seen, there can be significant pressurization of the HE, even prior to decomposition. The pressures reached after heating to 500 K are comparable to yield strengths in many metals.

FIG. 7.

Isochoric pressurization due to thermal heating starting with material initially at p = 0, T = T0.

FIG. 7.

Isochoric pressurization due to thermal heating starting with material initially at p = 0, T = T0.

Close modal

In this section, two validation tests are performed on the current model. Comparisons with isentropic and isothermal compression experiments are explored.

By loading an HE sample via a continuous pressure pulse, isentropic compression of the material can be achieved, at least for a limited thickness of the sample. The Sandia National Laboratory Z-machine33 is one such facility to isentropically compress solid materials to high pressures, and provides useful data for validating the current EOS model. Comparisons with reactive flow models of two such isentropic compression experiments have previously been made.34 The computational results presented here follow the Lagrangian methodology outlined previously.35,36

Shot Z1489 (Ref. 37) compressed PBX 9502 up to ∼5 GPa. A schematic of the experiment is shown in Fig. 8. The main diagnostic is the measured interfacial velocity between the sample HE and a LiF window. Figure 9 shows the results of the experiment and the present model. It is seen that the model “drive” matches the measured drive by construction.38 For the various PBX 9502 sample thicknesses, the model does a reasonable job predicting the experimentally observed traces, although the model is a bit too stiff at these low pressures, resulting in early interface motion. This stiffness at relatively low pressures is also observed in the Hugoniot in the range 0.2 mm/μs < Up < 0.6 mm/μs in Fig. 4. It is difficult for an analytical EOS to reproduce the “kink” observed in the Hugoniot from Fig. 4, without over- and under-shoots near the kink.

FIG. 8.

Schematic of the Z-machine compression shots. Only two of the five stations are shown in the schematic.

FIG. 8.

Schematic of the Z-machine compression shots. Only two of the five stations are shown in the schematic.

Close modal
FIG. 9.

Interface velocities for Z1489. Symbols are experimental results from various PBX 9502 sample thicknesses (for clarity only 5% of the experimental points are shown in the plot). Solid curves are computational results for the present model.

FIG. 9.

Interface velocities for Z1489. Symbols are experimental results from various PBX 9502 sample thicknesses (for clarity only 5% of the experimental points are shown in the plot). Solid curves are computational results for the present model.

Close modal

Shot Z1981 compressed PBX 9502 to ∼20 GPa. See Fig. 10 for comparison between the experimental and model interface velocities at various HE sample thicknesses. The current model matches the experiment quite well in both timing and amplitude of the waves. Again, it should be noted that the model results presented for these isentropic compression experiments are predictions and were not part of the calibration procedure.

FIG. 10.

Interface velocities for Z1981. Symbols are experimental results from various PBX 9502 sample thicknesses (for clarity only 5% of the experimental points are shown in the plot). Solid curves are computational results for the present model.

FIG. 10.

Interface velocities for Z1981. Symbols are experimental results from various PBX 9502 sample thicknesses (for clarity only 5% of the experimental points are shown in the plot). Solid curves are computational results for the present model.

Close modal

It should be noted in both Z1489 and Z1981 that there is early motion not well modeled by the simple fluid EOS utilized here. In particular, at very low pressures of less than 0.1 GPa, and particle velocities less than 0.025 mm/μs, the Hugoniot elastic limit17 has not been reached and strength effects of the HE are important. In reality, the composite and porous nature of the material at these low pressures are also important aspects and such modeling is outside the scope of the present study.

Isothermal compression of the PBX 9502 constituents has been investigated.20,39,40 Earlier works are cited in Refs. 20 and 40. Although discrete isothermal compression data points on TATB and Kel-F have been taken, there is no straightforward method to combine these discrete data to obtain the mixture's isothermal response. A simple methodology entails fitting the discrete data with an analytical fitting form for each of the constituents and determining the response to the overall mixture. One particularly convenient analytical form is from Murnaghan41 

p=KTKT((ρρ0)KT1),
(15)

where KT and KT are the isothermal bulk modulus at zero pressure and its derivative with respect to pressure, respectively. This particular fitting form allows analytical inversion to obtain ρ(p), which is convenient for obtaining the mixture response. From the data, for TATB, it is found that KT = 5.82 GPa and KT=1.85 with ρ0 = 1.937 g/cm3. For Kel-F 800, the parameters KT = 5.77 GPa and KT=3.26 with ρ0 = 1.998 g/cm3 are found. Note that the bulk modulus and its derivative are quite sensitive to the particular fitting form used,20 and these parameters represent a global fit to the available data. See the blue curve in Fig. 11 for the resulting mixture PBX 9502 isothermal compression based on the Murnaghan fit to the experimental constituent data. Note that at p = 0, the mixture is assumed to be at theoretical maximum density, which in this case is ∼1.94 g/cm3.

FIG. 11.

Density versus pressure for isothermal compression. The red curve is the multi-component isothermal compression curve for the mixture of 95% TATB and 5% Kel-F 800 assuming a Murnaghan fitting form for the constituents. The blue curve is the current model.

FIG. 11.

Density versus pressure for isothermal compression. The red curve is the multi-component isothermal compression curve for the mixture of 95% TATB and 5% Kel-F 800 assuming a Murnaghan fitting form for the constituents. The blue curve is the current model.

Close modal

For the current Davis EOS model calibrated in Sec. III, the isothermal compression curve, from Eq. (2), at T = T0 is shown as the red curve in Fig. 11. At p = 0, the density is 1.890 g/cm3. The discrepancy in initial density is due to the fact that PBX 9502 is porous, with ∼2.5% void space. Upon compression above the yield point, one expects the voids to be removed, and should be comparable to the initially void-free theoretical maximum density mixture. As can be seen in the figure, the isothermal compression of the current model quickly rises in density under applied compression, and then tapers off to a more gradual rise in density. The current model and Murnaghan mixture of constituents agree in density within ∼2% for the currently available data up to ∼45 GPa, which is not unreasonable given the uncertainties in the experimental data, Murnaghan fit to that data and the limitations of the Davis analytical EOS model.

Although the current model is quite adept at modeling a wide range of thermodynamics, there is room for further improvements. Specifically, this model (and many others) could benefit from a better underlying functional form for the specific heat at constant volume as a function of temperature, as well as the functional form of the Grüneisen parameter. Also, it would be beneficial to model porosity explicitly, as it is a key driver in shock initiation of HE. One may also wish to examine both from a modeling and experimental point of view any possible phase changes under applied pressure/temperature loading of TATB.

Experimentally, measuring shock temperatures, isentropic and isothermal compression at different initial temperatures, and thus effectively measuring isochoric pressurization are worthy pursuits for validation and testing of the model.

The author is grateful to a number of experimentalists who readily shared their results, specifically, Raja Chellappa, Dana Dattelbaum, Richard Gustavsen, Larry Hill, Shawn McGrane, and Darla Thompson. The author would also like to acknowledge John Bdzil, Carlos Chiquete, William Davis, Nicolas Desbiens, Ralph Menikoff, James Quirk, Mark Short, Remy Sorin, D. Scott Stewart, and Bradley Wescott for fruitful discussions on the modeling aspects of this work. John Bdzil, Ralph Menikoff, and Matthew Price are acknowledged for proofreading the manuscript and making many useful suggestions. This work was performed under the auspices of the U.S. Department of Energy National Nuclear Security Administration under Contract No. DE-AC52-06NA25396.

1.
B. L.
Wescott
,
D. S.
Stewart
, and
W. C.
Davis
, “
Equation of state and reaction rate for condensed-phase explosives
,”
J. Appl. Phys.
98
,
053514
(
2005
).
2.
T. D.
Aslam
and
J. B.
Bdzil
, “
Numerical and theoretical investigations on detonation-inert confinement interactions
,” in
Twelfth Symposium (International) on Detonation
(
2002
).
3.
G. J.
Sharpe
and
J. B.
Bdzil
, “
Interactions of inert confiners with explosives
,”
J. Eng. Math.
54
,
273
(
2006
).
4.
T. D.
Aslam
and
J. B.
Bdzil
, “
Numerical and theoretical investigations on detonation-inert confinement interactions
,” in
Thirteenth Symposium (International) on Detonation
(
2006
).
5.
L. G.
Hill
and
T. D.
Aslam
, “
The detonation confinement effect: Theory, observations, and experiments
,” in
Fifteenth Symposium (International) on Detonation
(
2014
).
6.
T. D.
Aslam
, “
Shock temperature dependent rate law for plastic bonded explosives
,”
J. Appl. Phys.
(unpublished).
7.
N.
Desbiens
,
C.
Matignon
, and
R.
Sorin
, “
Temperature-based model for condensed-phase explosive detonation
,” in
Journal of Physics: Conference Series
(
IOP Publishing
,
2014
), Vol. 500, No. 15, p.
152004
.
8.
C.
Chiquete
,
M.
Short
,
C. D.
Meyer
, and
J. J.
Quirk
, “
High explosive performance modeling using the Pseudo-Reaction Zone energy release model for Programmed Burn
,”
Combust. Flame
(unpublished).
9.
W. C.
Davis
, “
Complete equation of state for unreacted solid explosive
,”
Combust. Flame
120
,
399
(
2000
).
10.
D. S.
Stewart
,
W. C.
Davis
, and
S.
Yoo
, “
Equation of state for modeling the detonation reaction zone
,” in
Twelfth Symposium (International) on Detonation
(
2002
).
11.
R.
Menikoff
, “
Empirical equations of state for solids
,” in
Shock Wave Science and Technology Reference Library
(
Springer
,
Berlin, Heidelberg
,
2007
), pp.
143
188
.
12.
R.
Menikoff
, “
Complete Mie-Gruneisen equation of state
,”
Los Alamos National Laboratory Report LA-UR-16-21706
,
2016
.
13.
L. G.
Hill
and
T. D.
Aslam
, “
Detonation shock dynamics calibration for PBX 9502 with temperature, density, and material lot variations
,” in
Fourteenth International Detonation Symposium
(
2010
).
14.
T. R.
Gibbs
,
LASL Explosive Property Data
(
University of California Press
,
1980
), Vol.
4
.
15.
R.
Menikoff
, “
Complete EOS for PBX 9502
,”
Los Alamos National Laboratory Report LA-UR-09-06529
,
2009
.
16.
S. D.
McGrane
and
A. P.
Shreve
, “
Temperature-dependent Raman spectra of triaminotrinitrobenzene: Anharmonic mode couplings in an energetic material
,”
J. Chem. Phys.
119
,
5834
(
2003
).
17.
J. J.
Dick
,
C. A.
Forest
,
J. B.
Ramsay
, and
W. L.
Seitz
, “
The Hugoniot and shock sensitivity of a plastic-bonded TATB explosive PBX 9502
,”
J. Appl. Phys.
63
,
4881
(
1988
).
18.
R. L.
Gustavsen
,
S. A.
Sheffield
, and
R. R.
Alcon
, “
Measurements of shock initiation in the tri-amino-tri-nitro-benzene based explosive PBX 9502: Wave forms from embedded gauges and comparison of four different material lots
,”
J. Appl. Phys.
99
,
114907
(
2006
).
19.
B. D.
Lambourn
,
N. J.
Whitworth
,
C. A.
Handley
, and
H. R.
James
, “
A finite strain, non-reacted EOS FOR PBX9502
,”
AIP Conf. Proc.
955
,
137
(
2007
).
20.
L. L.
Stevens
,
N.
Velisavljevic
,
D. E.
Hooks
, and
D. M.
Dattelbaum
, “
Hydrostatic compression curve for triamino-trinitrobenzene determined to 13.0 GPa with powder x-ray diffraction
,”
Propellants, Explos., Pyrotech.
33
,
286
(
2008
).
21.
A. W.
Campbell
, “
Diameter effect and failure diameter of a TATB-based explosive
,”
Propellants, Explos., Pyrotech.
9
,
183
(
1984
).
22.
J. C.
Dallman
and
J.
Wackerle
, “
Temperature-dependent shock initiation of TATB-based high explosives
,”
Los Alamos National Laboratory Report LA-UR–93-2904, and Tenth Symposium
(International) on Detonation,
1993
.
23.
L. G.
Hill
,
J. B.
Bdzil
, and
T. D.
Aslam
, “
Front curvature rate stick measurements and detonation shock dynamics calibration for PBX 9502 over a wide temperature range
,” in
Eleventh International Detonation Symposium
(
1998
).
24.
C. B.
Skidmore
,
T. A.
Butler
, and
C. W.
Sandoval
, “
The elusive coefficients of thermal expansion in PBX 9502
,”
Los Alamos National Laboratory Report LA-14003
,
2003
.
25.
D. G.
Thompson
, private communication (
2015
).
26.
P. C.
Souers
,
P.
Lewis
,
M.
Hoffman
, and
B.
Cunningham
, “
Thermal expansion of LX-17, PBX 9502, and ultrafine TATB
,”
Propellants, Explos., Pyrotech.
36
,
335
(
2011
).
27.
S. P.
Marsh
,
LASL Shock Hugoniot Data
(
University of California Press
,
1980
).
28.
R. L.
Gustavsen
,
R. J.
Gehr
,
S. M.
Bucholtz
,
R. R.
Alcon
, and
B. D.
Bartram
, “
Shock initiation of the tri-amino-tri-nitro-benzene based explosive PBX 9502 cooled to –55°C
,”
J. Appl. Phys.
112
,
074909
(
2012
).
29.
B. C.
Hollowell
,
R. L.
Gustavsen
,
D. M.
Dattelbaum
, and
B. D.
Bartram
, “
Shock initiation of the TATB-based explosive PBX 9502 cooled to 77 Kelvin
,”
J. Phys.: Conf. Ser.
500
,
182014
(
2014
).
30.
R. L.
Gustavsen
,
R. J.
Gehr
,
S. M.
Bucholtz
,
A. H.
Pacheco
, and
B. D.
Bartram
, “
Shock initiation of the TATB-based explosive PBX-9502 heated to ∼ 76C
,”
AIP Conf. Proc.
1793
,
030017
(
2017
).
31.
R. L.
Gustavsen
,
B. D.
Bartram
,
L. L.
Gibson
,
A. H.
Pacheco
,
J. D.
Jones
, and
A. B.
Goodbody
, “
Shock initiation of the TATB-based explosive PBX 9502 heated to 130 degrees C
,” in
20th Biennial Conference of the APS Topical Group on Shock Compression of Condensed Matter
,
St. Louis, MO
(
2017
).
32.
S. D.
McGrane
,
T. D.
Aslam
,
T. H.
Pierce
, and
S. J.
Hare
, “
Temperature of shocked PBX 9502 measured with spontaneous Stokes/anti-Stokes Raman
,”
J. Appl. Phys.
(unpublished).
33.
J.-P.
Davis
,
C.
Deeney
,
M. D.
Knudson
,
R. W.
Lemke
,
T. D.
Pointon
, and
D. E.
Bliss
, “
Magnetically driven isentropic compression to multimegabar pressures using shaped current pulses on the Z accelerator
,”
Phys. Plasmas
12
,
056310
(
2005
).
34.
T. D.
Aslam
,
R. L.
Gustavsen
,
N. J.
Whitworth
,
H. J.
Lacy
,
B. D.
Lambourn
,
C. A.
Handley
,
H. R.
James
, and
S.
Root
, “
Multi-shock and isentropic compression of the tri-amino-tri-nitro-benzene based explosive PBX 9502: Evaluation of reactive flow models
,” in
Fifteenth Symposium (International) on Detonation
(
2014
).
35.
T. D.
Aslam
,
R.
Gustavsen
,
N.
Sanchez
, and
B. D.
Bartram
, “
An equation of state for polymethylpentene (TPX) including multi-shock response
,”
AIP Conf. Proc.
1426
,
767
(
2012
).
36.
T. D.
Aslam
,
R.
Gustavsen
, and
B. D.
Bartram
, “
An equation of state for polyurea aerogel based on multi-shock response
,”
J. Phys.: Conf. Ser.
500
,
032001
(
2014
).
37.
M. R.
Baer
,
M. L.
Hobbs
,
C. A.
Hall
,
D. E.
Hooks
,
R. L.
Gustavsen
,
D.
Dattelbaum
, and
S. A.
Sheffield
, “
Isentropic compression studies of energetic composite constituents
,”
AIP Conf. Proc.
955
,
1165
(
2007
).
38.
D. B.
Hayes
, “
Backward integration of the equations of motion to correct for free surface perturbations
,”
Sandia national Laboratory Report No. SAND2001-1440
,
2001
.
39.
R.
Chellappa
and
D. M.
Dattelbaum
, private communication (
2016
).
40.
L. L.
Stevens
,
D. M.
Dattelbaum
,
M.
Ahart
, and
R. J.
Hemley
, “
High-pressure elastic properties of a fluorinated copolymer: Poly (chlorotrifluoroethylene-co-vinylidene fluoride)(Kel-F 800)
,”
J. Appl. Phys.
112
,
023523
(
2012
).
41.
F. D.
Murnaghan
, “
The compressibility of media under extreme pressures
,”
Proc. Natl. Acad. Sci. U.S.A.
30
,
244
(
1944
).