DFT Studies on Two Novel Explosives Based on the Guanidine-Fused Bicyclic Structure

Abstract

Density functional theory (DFT) calculations at the B3LYP/6-31G(d,p) theoretical level were performed for two novel explosives (compounds B and C) based on the guanidine-fused bicyclic skeleton $C_4N_6H_8$ (A). The heats of formation (HOFs) were calculated via isodesmic reaction. The detonation properties were evaluated by using the Kamlet-Jacobs equations. The bond dissociation energies (BDEs) for the thermolysis initiation bond were also analyzed to investigate the thermal stability. The results show that the compounds have high positive HOF values (B, 1064.68 $kJ{\cdot}mol^{-1}$; C, 724.02 $kJ{\cdot}mol^{-1}$), high detonation properties (${\rho}$, D and P values of 2.04 $g{\cdot}cm^{-3}$ and 2.21 $g{\cdot}cm^{-3}$, 9.98 $km{\cdot}s^{-1}$ and 10.99 $km{\cdot}s^{-1}$, 46.44 GPa and 59.91 Gpa, respectively) and meet the basic stability requirement. Additionally, feasible synthetic routes of the these high energy density compounds (HEDCs) were also proposed via retrosynthetic analysis.

Introduction

During the last few decades, high energy density compounds (HEDCs) have been receiving considerable attention because of their wide range of uses both in military and civilian applications.1-4 For example, from the original explosives 1,3,4,6-tetranitroglycouril (TNGU), hexahydro-1,3,5-trinitro- 1,3,5-triazine (RDX) and 1,3,5,7-tetranitro-1,3,5,7-tetraazacyclooctane (HMX) to the emerging explosives trans-1,4,5,8- tetranitro-1,4,5,8-tetraazadecalin (TNAD),5 2,4,6,8,10,12- hexanitro-2,4,6,8,10,12-hexaazaisowurtzitane (CL-20)6 and cis-2,4,6,8-tetranitro-1H,5H-2,4,6,8-tetraazabicyclo[3.3.0]octane (bicycle―HMX),7 they are all HEDCs with high positive heats of formation (HOFs) and excellent detonation properties. Although there are wide variety species of HEDCs, each has its own advantages and disadvantages. For example, TNGU, possesses superior detonation properties but poor stability in the moisture due to the carbonyl groups at the either end of the molecule. On the other hand, N-containing heterocyclic structures are traditional sources of energetic materials. Based on the above-described reasons, a new parent structure: guanidine-fused bicyclic skeleton C4N6H8 (Figure 1, A) was designed as an alternate of the glycoluril structure, in which the two ―C=O groups were instead of the ―C=N― groups. Coincidentally, C4N6H8 is the very structure of this type containing six nitrogen atoms in the ring Therefore, C4N6H8 (A) was selected as the parent structure for designing new HEDCs.

It is also fascinating that there are six N-hydrogen atoms located in the ring and can be substituted by additional functional groups. Among the functional groups, the ―NO2 and ―NF2 groups are two important groups in the synthesis of energetic materials because of their potential high density, energy and properties as solid propellant oxidizers.8 However, when attempting to synthetize the HEDCs with multiple ―NO2 and ―NF2 groups, it may be of great danger to humans and the environment both in synthesis and performance test. Fortunately, computer simulation, as an effective way in screening promising explosives without these shortcomings, provides us an efficient way to predict the detonation properties of them.

Figure 1.Structures of the title compounds and the guanidinefused bicyclic skeleton.

In the present study, two novel HEDCs (Figure 1, compounds B and C) were designed based on the structure of guanidinefused bicyclic skeleton C4N6H8. The molecular geometry, heat of formation, thermodynamic properties, thermal stability and detonation properties were investigated by using the density functional theory (DFT) method.9,10 The feasible synthetic routes of them were also proposed via retrosynthetic analysis. All these results may provide useful information for a better understanding of the two compounds.

 

Computational Method

All the computations were performed with the Gaussian 03 package11 at the (U) B3LYP12,13 method with 6-31G (d,p) basis set.14 The optimization was performed using the default convergence criteria in the programs. Vibrational frequencies for the optimized structure were calculated and the absence of any imaginary vibrational frequency confirms that the obtained structure corresponds to the true local energy minimum on the potential energy surface. Besides, based on the principle of statistical thermodynamics, standard molar heat capacity , standard molar entropy and standard molar ethalpy from 200 to 800 K were calculated to predict their thermodynamic properties.

The HOF for a HEDC is one of the most important parameter to calculate their detonation properties. In order to get the accurate HOFs, the isodesmic reactions, in which not only the bonds and electronic pairs of the reaction are kept equivalent both in products and reactants, but also could counterbalance the error of electronic correlation energies, were designed. Previous work also demonstrated its feasibility in estimating the accurate HOFs of HEDCs.15,16

The isodesmic reaction that was employed to calculate the accurate HOF value of the title compound at 298 K can be written as:

Now the most important task is to compute ΔH298 which can be calculated using the following expression:

where ΔHf,p and ΔHf,R are the HOFs of products and reactants at 298 K, respectively. Then the HOFs of the title compounds can be figured out when the heat of reaction ΔH298 is known. On the other hand, the HOFs at 298 K can be also written as the following equation:

where ΔE0 is the change in total energy between the products and the reactants at 0 K, ΔZPE is the difference between the zero-point energy (ZPE) of the products and the reactants, and ΔHT is the thermal correction from 0 to 298 K. Δ(PV) equals ΔnRT for the reaction in gas phase. For isodesmic reactions, Δn = 0.

In the isodesmic reaction (1), the experimental HOFs of the reference compound NH3 and NH2NO2 are available.17,18 As the experimental HOF of C4N6H8, NH2NF2 are unavailable, additional calculation were carried out at the G2 level19 from the atomization reaction: CaHbNcOd→aC(g) + bH(g) + cN(g) + dO(g) to predict its HOF accurately. Thus, the HOFs of the target molecules can be calculated out via Eqs. (1)-(3) in combination with the atomization reaction described.

However, for most energetic compounds, whose condensed phase are solid, the calculation of detonation properties require solid-phase HOFs (Δf,solid) instead of gas-phase HOFs (ΔHf,gas). According to Hess’s law of constant heat summation,20 the gas-phase heat of formation (ΔHf,gas) and heat of sublimation (ΔHsub) can be used to evaluate their solid-phase heats of formation (ΔHf,solid):

Additionally, Politzer et al.21 found that the heats of sublimation of energetic compounds can correlate well with the molecular surface area and electrostatic interaction index by the following expression:

where A is the surface area of the 0.001 e·bohr−3 isosurface of electronic density of the molecule, ν is the degree of balance between positive and negative potential on the iso-surface, and is a measure of variability of the electrostatic potential on the molecular surface. The coefficients a, b and c were determined by Rice et al.:a = 2.670 × 10−4 kcal·mol−1 A−4, b = 1.650 kcal·mol−1, and c = 2.966 kcal·mol−1.22 The descriptors A, ν and & were calculated using the computational procedures as described by Felipe et al.23 This approach has been demonstrated as a reliable way to predict the heats of sublimation of many energetic compounds. 24,25

As for the detonation velocity and detonation pressure, they were estimated by the empirical Kamlet-Jacobs equations: 26

where ρ, the loaded density of explosives; D, the detonation velocity; P, the detonation pressure; N, the moles of detonation gases per-gram explosive; The average molecular weight of these gases; and Q, the heat of detonation. For known explosives, their Q and ρ can be measured experimentally and then their D and P can be calculated according to Eqs. (6) and (7). However, for novel HEDCs whose Q and ρ cannot be obtained from experiments, these parameters should be calculated. ρ was obtained from the molecular weight divided by the average molecular volume. The crystal density can also be improved by the introduction of the interaction index

where M is the molecular mass (g·mol−1), and V is the volume defined as inside a contour of 0.001 e·bohr−3 density that was evaluated using a Monte Carlo integration. The coefficients β1, β2, and β3 are 0.9183, 0.0028, and 0.0443, respectively.27

Additionally, the strength of bonding, which could be evaluated by BDE, is a fundamental parameter for us to understand the thermal stability of a novel HEDC. In order to investigate the thermal stabilities of the titled compounds, the bond order and BDE of the weakest bond were calculated. The energy required for bond homolysis at 1.0 atm and 298.0 K, which defined as the bond dissociation enthalpy of A – B, corresponds to the enthalpy of reaction: A – B(g) → A·(g) + B·(g). Coincidentally, for most organic molecules, the terms “bond dissociation energy” and “bond dissociation enthalpy” often appear interchangeably.28 There-fore, the BDE of the homolytic bond can be calculated via the following equation:

The BDE with zero-point energy (ZPE) correction can be calculated by Eq. (10).

where ΔEZPEis the difference between the ZPEs of the products and the reactants.

 

Results and Discussion

Molecular Geometry and Electronic Structure. Because all the calculations were based on the optimized structure, it is necessary to examine the geometric structure of the title compounds before discussing the various properties. The results show that all of the optimized structures were characterized to be true local energy minima on the potential energy surfaces without imaginary frequencies. Besides, all the calculated bond lengths are between the normal C―N, N―N bond lengths (1.47 Å and 1.45 Å, respectively) and normal C=N, N=N bond lengths (1.28 Å and 1.25 Å, respectively). So, all the bonds tend to be average and form a big conjugative system which may facilitate the stability of the compounds. The optimized geometric structures of the title compounds can be expressed as shown in Figure 2.

Analysis of the molecular orbital can provide useful information on electronic structures.29 It is proposed that the larger energy gap of a compound, the lower reactivity in the chemical or photochemical processes with electron transfer. 30,31 The values were summarized as follows: HOMO (B, −8.59398 eV; C, −7.98907 eV), LUMO (B, −3.83820 eV; C, −2.20823 eV), and their gap ΔE (B, 4.75578 eV; C, 5.78084 eV). It is evident that the energies of HOMO, LUMO, and the gap ΔE of compound C are slightly higher than that of compound B, indicating that compound C may have a comparable chemical reactivity with compound B. Figure 3 illustrates HOMO, LUMO, and molecular electrostatic potentials (MEPs) of the two compounds. Compound B is from reference [2]. As is vividly seen, either C―N or N―N orbital participate in both the HOMO and the LUMO levels, indicating that the removal of an electron from the HOMO level or addition of an electron to the LUMO level could weaken the skeleton framework. Inspection of the MEPs for the two compounds, the positive potentials which increases some stabilization to a compound, ranges at the center of the bicyclic skeleton while the negative potentials appear to be distributed mostly on the oxygen or fluorine atoms.

Figure 2.The optimized geometric structures of the title compounds at the B3LYP/6-31G(d,p) level.

Thermodynamic Properties and Heat of Formation. Figure 4 presents the simulated IR spectrum of compounds B and C. There are several main absorption bands and the two compounds have the similar IR spectra. The modes in 3000-3300 cm−1 are associated with the C―H stretch vibration of the guanidine-fused bicyclic skeleton; the remarkable peak in 1600-1800 cm−1 are associated with the N=O or N―F asymmetric stretch of ―NO2 and –NF2 groups; band at 1000 cm−1 is composed of the N―N asymmetric stretch of heterocyclic skeleton together with C―H twisting out of plane; The bands less than 900 cm−1 which belong to the fingerprint spectrum, are mainly caused by the deformation of heterocyclic skeleton and the bending vibration of C―H and C―C bonds.

Figure 3.HOMO, LUMO, and electrostatic potentials mapped onto 0.001 electron·bohr−3 contour of the electronic density at the B3LYP/6-31G(d,p) level. Potential color range: from red (negative) to blue (positive).

Figure 4.Simulated IR spectra for the title compounds at the B3LYP/6-31G(d,p) level.

Based on the vibrational analysis results and statistical thermodynamic method, thermodynamic properties such as standard molar heat capacity , standard molar entropy and standard molar enthalpy from 200 to 800 K were calculated to predict their thermodynamic properties. Figure 5 presents the temperature-dependent relations for , and in the range of 200-800 K, and the correlation equations between the thermodynamic functions and different temperatures were also calculated as follows (where R2 is the correlation coefficients). It can be seen that the , and increase evidently with the temperature increasing. This is because the main contributions to the thermodynamic functions are from the translation and rotation of molecules when the temperature is low; however, at higher temperature, the vibrations are intensified and therefore make more contributions to the thermodynamic properties which lead to the increase in the thermodynamic functions. All of the data may provide useful information on the thermodynamic properties of them.

Figure 5.Relationships between the thermodynamic functions ( : J·mol.−1·K−1, : J·mol.−1·K.−1, : kJ·mol.−1) and temperature (T/K) for the title compounds.

Compound A:

= 64.1764 +1.1495T − 0.00061T2 R2 = 0.9997

= 317.92286 + 1.4596T − 0.00047T2 R2 = 0.9999

= −19.3536 + 0.2046T + 0.00027T2 R2 = 0.9998

Compound B:

= 63.9564 + 1.3072T − 0.00077T2 R2 = 0.9999

= 307.6264 + 1.6220T − 0.00055T2 R2 = 0.9999

= −25.7107 + 0.2435T + 0.00026T2 R2 = 0.9997

The heat of formation, which is a basic property of HEDCs, is usually taken as the indicator of the “energy content”. Thereby, heat of formation is frequently considered when the HEDCs are designed. Table 1shows the total energies (E0), thermal corrections (HT), zero-point energies (ZPE), and HOFs for the reference compounds being enlisted in the isodesmic reaction. Table 2 summarizes the total energies (E0), thermal corrections (HT), zero-point energies (ZPE), gas-phase HOF (ΔHf,gas) and solid-phase HOF (ΔHf,solid) of the two compounds at the B3LYP/6-31G(d,p) level. From Table 2, it is obvious that when the H atoms of the N−H in C4N6H8 (A) was replaced by −NO2 and −NF2 groups, the HOFs of the two compounds are up to 1064.68 kJ·mol−1 and 724.02 kJ·mol−1 respectively. On the other hand, the HOF reduced sharply from 1064.68 kJ·mol−1 to 724.02 kJ·mol−1 when the six −NO2 groups in compound B were replaced by the six −NF2 groups in compound C. This indicates that NO2 group has a superiority in increasing the HOF compared with NF2 group. However, the HOFs of compounds B and C are higher than that of un-substituted C4N6H8 (A) and thus, scientists devote their efforts to synthetize HEDCs with multi −NF2 and −NO2 groups nowadays.

Table 1.Calculated total energies (E0, au), thermal corrections (HT, kJ·mol−1), zero-point energies (ZPE, kJ·mol−1), and heats of formation (HOFs, kJ·mol−1) for the reference compounds

Detonation Properties. Detonation velocity and detonation pressure are two important parameters for energetic materials. The Kamlet–Jacobs equations show that ρ is a key factor to influence D and P. Thus, density is one of the most important physical properties for all energetic materials. The values of ρ, D and P of compounds B and C are listed in Table 3. It is seen that all the compounds have excellent detonation properties and compound C is calculated to have higher ρ, D and P values than that of compound B. This indicates that the −NO2, and −NF2 groups are effective units for improving the detonation properties of a compound, especially the −NF2 group. This supports the conclusion that although −NF2 decreases the HOF a little, it increases ρ most greatly and thus makes compound C possesses the highest D and P; −NO2 increases HOF and ρ moderately and thus makes D and P of compound B stand in the second. For a comparison, the detonation properties of the well-known explosives RDX, HMX and CL-20 are also listed in Table 3. Clearly, the designed compounds possess higher ρ, D and P than that of RDX, HMX and CL-20 and meets the quantitative criteria of a HEDC (that is, ρ ≈ 1.9 g·cm−3, D ≈ 9.0 km·s−1, and P ≈ 40.0 GPa). If these cyclic derivatives can be synthesized, they will have more exploitable values in the future.

Thermal Stability. Bond dissociation energy (BDE) of the trigger bond is another key parameter which can provide useful information for understanding the stability and sensitivity of HEDCs. Generally, the smaller the energy is needed for breaking a bond, the weaker the bond is, and becomes a trigger bond when heated or assaulted.

Table 2.Calculated total energies (E0, au), thermal corrections (HT, au), zero-point energies (ZPE, au), ΔHf,gas, ΔHf,sub, ΔHf,solid (kJ·mol−1) for the title compounds

Table 3.aExperimental value from Ref..33 bCalculated value from Ref.16

On the other hand, people nowadays have reached a consensus that N−NO2 or N−NF2 bond often represents the primary cause of initiation reactivity of organic polynitro or difluoroamino compounds.34,35 Therefore, the weakest N−N bonds, which were screened according to the “principle of smallest bond order (PSBO)” (the smallest bond order is, the least stable is),36 were selected as the breaking bond (B, N(7)−N(22); C, N(1)−N(17)) to calculate bond dissociation energy (BDE) at UB3LYP/6-31G(d,p) level. The values of BDEs of the relatively weaker bonds for compounds B and C were 91.47 kJ·mol−1 and 91.39 kJ·mol−1, respectively. Though the BDE of compound B is slightly higher than that of compound C, they are very close to each other. In other words, it is to say compound C has a lower stability (higher sensitivity) than compound B and may be explode firstly when heated or assaulted. This is in harmonic agreement with the fact that −NF2 group is more sensitive than −NO2 group for an explosive because of the strong electron attracting. And thus lower the stability of the molecule when −NO2 is substituted by −NF2. Furthermore, take the practical requirements into consideration, a quantitative criteria associated stability (BDE of the trigger bond) requirement, 𝑖.𝑒., BDE80∼120 kJ·mol−1, is proposed and employed by Chung et al. to filtrate a potential HEDC.37 Based on the law, it is found that both of compounds B and C satisfy the basic requirements and seem to be the potential HEDCs.

Feasible Synthetic Routes. Based on the above-calculated data, it is clear to see that both of the compounds are potential HEDCs with high detonation properties and exploitable values. Thus, the feasible synthetic rotes for the two compounds were proposed in Figure 6.

From the scheme, it is found that all of the crude materials are easily obtained with low price and it provides a necessary condition for the mass industrial production. Besides, the synthetic rotes are simple. All the above-described reasons prove the feasibility of the synthetic rotes.

Figure 6.Synthetic rotes of the title compounds B and C.

 

Conclusions

In the present work, the electronic structure, thermodynamic properties, heats of formation, detonation properties, and thermal stability were studied for two novel explosives based on guanidine-fused bicyclic skeleton C4N6H8 by using the B3LYP/6-31G(d,p) method of DFT theory. Calculation results show that the two compounds possess high positive heats of formation; the density and detonation properties of them are larger than that of RDX, HMX and CL-20; the bond dissociation energies of the pyrolysis initiation reaction calculated at the B3LYP/6-31G(d,p) level indicate that the title compounds meet basic requirements and seem to be the potential HEDCs. Besides, the synthetic rotes proposed are simple and feasible. All the data may provide useful information for a better understanding of physical and chemical properties of compounds B and C.