Decomposition Mechanisms and Kinetics of Novel Energetic Molecules BNFF-1 and ANFF-1: Quantum-Chemical Modeling

Decomposition mechanisms, activation barriers, Arrhenius parameters, and reaction kinetics of the novel explosive compounds, 3,4-bis(4-nitro-1,2,5-oxadiazol-3-yl)-1,2,5-oxadiazole (BNFF-1), and 3-(4-amino-1,2,5-oxadiazol-3-yl)-4-(4-nitro-1,2,5-oxadiazol-3-yl)-1,2,5-oxadiazole (ANFF-1) were explored by means of density functional theory with a range of functionals combined with variational transition state theory. BNFF-1 and ANFF-1 were recently suggested to be good candidates for insensitive high energy density materials. Our modeling reveals that the decomposition initiation in both BNFF-1 and ANFF-1 molecules is triggered by ring cleavage reactions while the further process is defined by a competition between two major pathways, the fast C-NO2 homolysis and slow nitro-nitrite isomerization releasing NO. We discuss insights on design of new energetic materials with targeted properties gained from our modeling.


Introduction
Explosive decomposition consists of three main stages: initiation, combustion, and detonation. The two latter aspects are fairly well understood especially on the macro-and meso-scales. By contrast, a detailed understanding of the initiation process has yet to be established. It is clear that the initiation of chemistry in energetic compounds is a complex process, and hence progress in searching for new and OPEN ACCESS improved energetic compounds has been rather slow. With recent advances in both computational capabilities [1,2] and ultra-fast time resolved spectroscopic (optical) diagnostics [3], attention has begun to focus on the characterization of the physicochemical responses of specific energetic ingredients such as localized mechanical stresses [4][5][6][7][8], crystal packing [9][10][11], morphology [12,13], and details of chemical decomposition reactions [14,15].
It is accepted that sensitivity to initiation of detonation (or a propensity of an energetic material to trigger a chemical reaction of explosive decomposition) depends upon many factors including chemical composition, bonding of functional groups, distribution of the charge density, stoichiometry, molecular ordering, grain size, presence of imperfections, nanoparticles, interfaces, etc. However, correlations between the structure, composition, and sensitivity have yet to be revealed. Fundamental insights into these relationships are essential so that energetic materials' properties can be tailored to influence macroscopic responses. Understanding factors that govern sensitivity and at the same time provide high performance represents an outstanding challenge. Making inroads into these questions is the goal of our study.

Computational Details
Calculations were performed using density functional theory (DFT) [28,29] with the PBE exchange correlation functional, [30] and hybrid B3LYP [31,32] and M06 [33] functionals, with double-zeta cc-pVDZ Dunning's correlation consistent basis set [34] as implemented in the GAUSSIAN 09 code package [35]. PBE functional is widely used for investigating solid state reaction mechanisms [2] (it is implemented in most periodic computer codes because of its efficiency) and hence gas-phase calculations are required before any studies of decomposition pathways in practical materials are attempted (see, for example, ref. [4][5][6]13]). Unlike PBE, both hybrid B3LYP and M06 functionals consume significant time and computer resources but are usually considered relatively reliable for modeling reactions in gaseous phase (see, for example, ref. [36]). Vibrational frequencies have been calculated for relevant atomistic configurations to distinguish energy minima and transition states and to determine corresponding zero-point energy (ZPE) corrections. The stationary points corresponding to the energy minimum have been positively identified by having no imaginary frequencies and the transition states had exactly one imaginary frequency. An intrinsic reaction coordinate (IRC) analysis was carried out by using the Hessian-based Predictor-Corrector integrator algorithm [37,38] for each transition state to ensure that it was the transition structure, connecting the desired reactants and products. Pre-exponential factors were calculated with conventional transition state theory (TST) [39] for the decomposition reactions that proceed through a formation of a transition state and with variational transition state theory (VTST) [40] for the homolytic cleavage pathways, as detailed elsewhere [13,[41][42][43].

Simulating Chemical Decomposition Reactions
The geometric structures of the BNFF-1 and ANFF-1 molecules in equilibrium obtained with the M06 method are depicted in Figure 2. Both BNFF-1 and ANFF-1 molecules have fairly similar nonplanar structures. The calculated equilibrium bond distances of the BNFF-1 molecule are very close to those of ANFF-1. Figure 2 demonstrates that a substitution of the nitro group in BNFF-1 by the amino group in ANFF-1 affects only the fragment connected to the nitro-group (the C5-C6-N5-O4-N4 ring and C4-C5 bond distances) and introduces negligible changes in the rest of the molecule. The C5-C6, C6-N5, N5-O4 bond distances in ANFF-1 are slightly elongated and O4-N4 and C4-C5 bond lengths are slightly decreased as compared to BNFF-1. The torsion angles φ(C4C3C2C1) and φ(C3C4C5C6) corresponding to a rotation of the ring fragments about C2-C3 and C4-C5 bonds in BNFF-1 are 51.3° and −162.9°, respectively. A nearly planar fragment (φ(C3C4C5C6) = −177.8°) of the ANFF-1 molecule consists of the central and amino-furazan rings; the nitro-furazan fragment is rotated about C2-C3 bond by 53.9°. All geometry parameters obtained within different DFT functionals employed agree well with each other with the average discrepancy between bond distances in ANFF-1 calculated at M06 and B3LYP does not exceed 0.6%. Thermal stability of the BNFF-1 and ANFF-1 molecules was explored through modeling the five most plausible decomposition mechanisms, schematically illustrated in Figure 3. The homolytic cleavage of the C-N bond (reaction 1, Table 1), leading to NO 2 loss, is usually assumed to be the most favorable decomposition channel of nitro compounds. The two-step reaction of the nitro-nitrite rearrangement (CONO-isomer formation), proceeding via a pseudo rotation of the nitro group (reaction 2, Table 1), results in a release of NO and is often considered as an alternative primary dissociation reaction of nitro-arenes. The oxygen transfer involves the bond switching of the oxygen atom from the nitro group to the ring nitrogen (reaction 3, Table 1). We also probed the homolytic cleavage of the C-C bonds that connect two rings (reactions 4, Table 1) and a fragmentation of the molecules through the furazan ring cleavage (RC, reactions 5, Table 1).

The Homolytic C-NO 2 Break
The homolytic C-NO 2 break was simulated following the standard Equation (1): where X stands for a nitro-group in BNFF-1 and an amino-group in ANFF-1. According to our calculations, the reaction barrier in BNFF-1 requires 57.6-63.1 kcal/mol depending on the computational method used, indicating a slightly higher activation barrier for the nitro-group (II) than for the nitro-group (I), shown in Table 1. The difference in obtained energies is however very small and does not allow for a definite discrimination of one nitro-group vs another. All three DFT functionals, PBE, M06, and B3LYP give consistent energies, which are very close to each other, with the M06 value being the highest (63.1 kcal/mol for nitro-group (I) and 65.6 for the nitrogroup (II)), the B3LYP energy being the lowest (57.6 and 60.0 kcal/mol, respectively), and PBE energy falling right in the middle (61.1 and 62.2 kcal/mol, respectively). The similar energies obtained for the ANFF-1, 59.5-61.5 kcal/mol are even closer to each other, showing practically no dispersion. The estimated pre-exponential factors point towards the NO 2 loss as the fastest overall initiation reaction among all investigated mechanisms (Table 1).

The CONO-Isomerization Precursor of the NO Loss
The nitro-nitrite isomerization proceeds via a pseudo rotation of a nitro group accompanied by breaking the C-NO 2 bond and creating the C-ONO bond [Equation (2)]: The formed nitrite isomers of BNFF-1 and ANFF-1 are shown in Figure 4a-c and the transitions states are displayed in Figure 5. The CONO-isomerization in BNFF-1 appears to require a somewhat lower activation barrier (by 4-15 kcal/mol depending on the method) than the NO 2 loss reaction pathway (Table 1), and the NO loss needs an even lower energy of 19.8-30.9 kcal/mol. Besides, the CONO-isomerization is an exothermic reaction, releasing heat [54] of 12.8-14.6 kcal/mol (Table 1).  Although the nitro-group (I) shows a slightly smaller activation barrier and a slightly higher reaction heat than the nitro-group (II), the observed difference is too small to claim a distinction. In the meantime, the dispersion of energies resulting from different functionals deserves some attention. PBE gives the lowest energy, 46.1 kcal/mol, for the CONO-isomerization versus 53.2 kcal/mol obtained with both M06 and B3LYP. PBE also yields the highest energy for the NO loss, 30.9 kcal/mol, versus 25.0 kcal/mol obtained with M06 and 19.8 kcal/mol obtained with B3LYP (see the reaction 2 in Table 1). Given that M06 and B3LYP produce the same value of 53.2 kcal/mol and taking into account that PBE tends to underestimate molecular rearrangements [36], the difference of 4 to 10 kcal/mol between the CONO-isomer formation and the NO 2 loss reactions appears more realistic than 4-15 kcal/mol. Similar results are obtained for the ANFF-1 molecule with close activation barriers of the CONO-isomerization stage ranging from 51.9 to 60.3 kcal/mol, an appreciable difference in energies of the NO-loss stage (18.6-29.8 kcal/mol), and the reaction heat of 8.2-10.7 kcal/mol.
An analysis of the kinetic parameters, generated by both PBE and M06, predicts a competition between the CONO formation and the C-NO 2 fission in BNFF-1 over a wide temperature range while kinetics obtained from B3LYP suggests that the C-NO 2 homolysis is favored over all other probed mechanisms (Table 1 and Figure 6). Unlike BNFF-1, only PBE supports the CONO-isomerization while both M06 and B3LYP prefer the C-NO 2 break in the initiation of decomposition of ANFF-1. Hence, the CONO preference achieved from the energetic considerations is practically negated once the kinetics is taken into account.
We emphasize here that energetic considerations should be coupled with kinetics to understand the thermal decomposition, especially, when several chemical pathways co-exist. Thus, the exothermic nitro-nitrite isomerization and the NO loss require a relatively low energy but the pre-exponential factor of this dissociation reaction is estimated to be much smaller than that of the NO 2 loss (Table 1). Hence, our analysis suggests that the C-NO 2 fission is the fastest dissociation process while the NO loss contributes to decomposition as a slow background reaction. Interestingly, neither of these two reactions is the dominant initiation mechanism in the decomposition of BNFF-1 and ANFF-1 (vide infra).

Oxygen Transfer
Oxygen transfer, another type of isomerization, results in a formation of nitroso-furoxan isomers of BNFF-1 and ANFF-1 displayed in Figures 4d-f [Equation (3)]. The O-transfer mechanisms in both molecules have significantly higher activation barriers (~70-80 kcal/mol in BNFF-1 and ~75-89 kcal/mol in ANFF-1) than the corresponding homolytic cleavage of the C-NO 2 bonds and the CONO-isomerization reactions and exhibit fairly low pre-exponential factors (reaction 3 in Table 1 and Figure 3). The difference obtained for the different nitro-groups is again small. Hence, this reaction should be ruled out as a possible candidate to initiate the decomposition of either BNFF-1 or ANFF-1. Figure 6. Reaction rates of the simulated decomposition channels of BNFF-1 (a-c) and ANFF-1 (d-f) molecules.
We do not consider a hydrogen transfer here for two main reasons: (1) the hydrogen migration is strongly dependent upon the molecular environment (and it is therefore unreliable to model it for molecules) and (2) the HONO reaction typically represents only a secondary isomerization pathway in nitro-arenes [52,62].

The Homolytic C-C Break
The homolytic fission of the C-C bonds in the molecules [Equation 4(a,b)] are also energetically unfavorable reaction pathways with the calculated dissociation energies of ~121-125 kcal/mol in BNFF-1 and ~123-128 kcal/mol in ANFF-1, which are much higher than the corresponding energies required for breaking the C-N bonds (reactions 4, Table 1 and Figure 3). These estimates agree with the C-C dissociation energy of a DADNE molecule (120.2 kcal/mol) [52].
The difference in energy between splitting off of the peripheral nitro-furzan rings in BNFF-1 does not exceed 2 kcal/mol. The distinction between the nitro-furazan ring's loss and amino-furazan ring loss is comparable, ~4 kcal/mol, signifying a measure of the substitution effect (reaction 4, Table 1). Hence, the splitting off the rings would not contribute to the initial stages of the overall decomposition of either molecule.

A Heterocyclic Ring Cleavage (RC)
Opening of the heterocyclic ring may be in principle achieved through breaking one of the ring's bonds (see, for example, ref. [47]). We explored concerted fragmentation mechanisms of BNFF-1 and ANFF-1 molecules (path 5 in Figure 3 and reaction 5 in Table 1 The transition state structures, shown in Figure 7, indicate an elongation of N-O (~1 Å) and C-C (~0.6 Å) bonds in the furazan rings. The calculations reveal that the ring cleavage reactions in BNFF-1 require 42.9-51.1 kcal/mol of energy. The endothermic cleavage of the central ring ( Figure 3) in BNFF-1 requires 42.9-47.6 kcal/mol, with the lowest barrier (42.9 kcal/mol) obtained with PBE and very close barriers of 47.6 and 47.5 kcal/mol obtained with M06 and B3LYP, respectively. This energy is slightly (~3 kcal/mol) lower than the barriers for the opening of peripheral rings (reaction 5 in Table 1). An energy dispersion obtained with different functionals is ~5 kcal/mol, with a systematic difference observed between PBE and hybrid methods and almost identical energies received from M06 and B3LYP. The energy difference between outer nitro-furazan rings does not exceed 1 kcal/mol. Consistent activation barriers (42.7-55.6 kcal/mol) and reaction energies are obtained for ANFF-1. The same trends are seen: the central ring opening needs a lower energy than the peripheral rings, and all three methods produce corresponding energies with the lowest barrier (42.7 kcal/mol) delivered by PBE and nearly identical barriers (47.6 and 47.1 kcal/mol obtained from M06 and B3LYP). The difference in energy between cleavage of the nitro-furazan ring and the amino-furazan ring is negligible.
Most importantly, the concerted ring opening reactions exhibit the lowest activation barriers among all mechanisms explored in both molecules, BNFF-1 and in ANFF-1. These values agree with the activation energy of 42.3 kcal/mol obtained from differential scanning calorimetry measurements of BNFF [21]. Pre-exponential factors (log(A/ s−1 ) ~15, Table 1) of the ring cleavage mechanisms are still lower than the NO 2 loss but visibly higher than the CONO isomerization. The estimated reaction constant is consistent with the log(A/ s−1 ) = 13.68, reported for BNFF [21]. Kinetic behavior is illustrated in Figure 6, which shows that the fragmentation through the central ring opening dominates the decomposition initiation in both molecules at low and moderate temperatures (regardless of the computational method used). The NO 2 loss, while being the fastest reaction, prevails only at high temperatures. The exothermic CONO-isomerization pathway makes a modest contribution to decomposition at the earliest stages.

Comparing BNFF-1 and ANFF-1
Next, let us turn to the difference between the BNFF-1 and ANFF-1 decomposition processes, or discuss the substitution effect here. In both molecules, the cleavage of the central ring seems to be the most favorable dissociation process from an energetic point of view (42.9-47.6 kcal/mol in BNFF-1 vs. 42.7-47.6 kcal/mol in ANFF-1), with a negligible difference between ANFF-1 and BNFF-1. The energy of the C-NO 2 bond dissociation in ANFF-1 (61.1 delivered by PBE, 63.1 by M06, and 57.6 kcal/mol by B3LYP) also shows only a small change as compared to BNFF-1 (61.3, 61.5, and 59.5 kcal/mol, respectively). The C-NO 2 bond pre-exponential factors evidence that the detachment of NO 2 from BNFF-1 and ANFF-1 will proceed at comparable rates and faster than other reactions. At the same time, dissociation energies coupled with kinetic parameters point towards the ring opening as a dominating dissociation process at low and moderate temperatures while the C-NO 2 homolysis becomes important at high temperatures.
The most notable effect of the BNFF-1's nitro group replacement with the ANFF-1's amino group appears in the barrier heights and reaction energies of the CONO decomposition channel. The CONOisomerization of ANFF-1 (51.9 kcal/mol found by PBE, 59.1 kcal/mol by B3LYP, and 60.3 kcal/mol by M06) requires ~6 kcal/mol higher energy in comparison to BNFF-1 (46.1 kcal/mol obtained by PBE and 53.2 kcal/mol by M06 and B3LYP) while the exothermicity of this pathway is reduced by ~4 kcal/mol (Table 1). Due to a somewhat higher activation barrier and somewhat lower reaction heat of ANFF-1 than those of BNFF-1, the contribution of the slow exothermic CONO-isomerization reaction in BNFF-1 is expected to be more pronounced as compared to ANFF-1. However, we suggest that this process has only a minor effect on the overall thermal dissociation and hence the presence of nitro-or amino-groups in these heterocyclic molecules does not make a significant difference.

Summary and Conclusions
Modeling of candidate decomposition channels of new energetic compounds BNFF-1 and ANFF-1 was performed by means of density functional theory combined with transition state theory. The obtained activation barriers, reaction energies, nature of the transition states, pre-exponential factors, and reaction rates were analyzed for five plausible decomposition channels, including the homolytic NO 2 loss, step-wise NO loss, switching of N-O bonds through the oxygen transfer, splitting off of the heterocyclic rings by cleavage of the C-C bond, and concerted opening of the heterocyclic rings. A careful evaluation revealed that the cleavage of the central furazan ring dominates the overall decomposition initiation of both BNFF-1 and ANFF-1 molecules, and breaking the peripheral rings requires only slightly higher energy. The obtained decomposition activation barriers of both BNFF-1 and ANFF-1 agree with DSC measurements of BNFF. The next important decomposition mechanism is the homolytic C-NO 2 bond cleavage. Despite being the fastest reaction, the C-NO 2 homolysis prevails only at high temperatures due to a notably higher activation barrier than that of the ring opening process. The two-step NO loss in both molecules, while it exhibits a somewhat higher activation barrier than the ring opening reactions and may co-exists with the primary dissociation channel, proceeds at a much lower rate.
Thus, our calculations established that the decomposition scenario in heterocyclic nitro-furazanbased molecules differs from known nitro-arenes. Unlike nitro-arenes, the thermal decomposition in BNFF analogs is strongly affected by the presence of heterocyclic rings and the ring's opening reactions trigger the overall bond dissociation. A further development of the dissociation processes is defined by the interplay of the slow exothermic CONO-isomerization and the fast C-NO 2 homolysis reactions. While the activation barrier of the NO loss requires a lower energy than the NO 2 detachment, the reaction kinetics happened to be orders of magnitude slower for the CONO-isomerization than for C-NO 2 homolysis.
The substitution of one nitro group in BNFF-1 with the amino group in ANFF-1 does not affect the ring opening reactions, negligibly decreases the C-NO 2 bond strength, tends to increase the CONO-isomerization activation barrier and to reduce the reaction heat.
The analyzed decomposition pathways provide valuable insight on the dissociation of bonds in these complex molecules. However the derived conclusions cannot be simply extrapolated to predict materials' behavior. The next step in this research is sophisticated solid-state periodic calculations involving modeling of similar chemical reactions in the crystalline environment [4,5,13,14] and on surfaces [13,42,43,63]. In addition, the effect of crystallographic and electronic defects [15,52,62] has to be taken into account as, at least for DADNE, it was discovered that the shear-strain [4,5,14] and other imperfections [12] severely affect the decomposition process by reducing the decomposition barriers and/or enhancing exothermicity [52] and hence decrease sensitivity of the material to initiation of detonation.