Reactive molecular dynamics simulations were conducted to study the flame speed enhancement phenomenon of a solid mono-propellant, Pentaerythritol Tetranitrate (PETN), when coupled to highly conductive multi-walled carbon nanotubes (MWCNTs). The simulations were based on the first-principles derived reactive force field, ReaxFF, which includes both the physical changes such as thermal transport and the chemical changes such as bond breaking and forming. An annular deposition of a PETN layer around the MWCNTs was considered. The thickness of the PETN layer and the diameter of the MWCNT were varied to understand the effect of the MWCNT loading ratio on the flame propagation. Flame speed enhancements up to 3 times the bulk value were observed. An optimal MWCNT loading ratio was determined. The enhancement was attributed to the layering of the PETN molecules around the MWCNT, which increased the heat transport among the PETN molecules near the MWCNT surface, thus causing the flame to travel faster. Furthermore, a stronger ignition source was required for the MWCNT-PETN complex because of the higher thermal transport among the PETN molecules along the MWCNT, which makes the ignition energy dissipate more quickly. Lastly, the MWCNT remained unburned during the PETN combustion process.

Carbon-based structures, such as carbon nanotubes,1 2-D and 3-D grapheme,2,3 and graphite sheets,4 because of their high thermal conductivity (1600–4600 W/m K)5–10 and large surface-to-volume ratio, have been used as thermal interface materials,11–14 as fillers to enhance the thermal conductivity of various composites,3,12,15–31 and as heat exchangers in nano-electronic devices.32,33 For example, the thermal conductivity of epoxy, as shown by Biercuk et al.,15 was increased by 125% by using Single-walled carbon nanotubes (SWCNTs) at 1% volumetric loading. However, Huang et al.17 conducted experiments with vertically aligned multi-walled carbon nanotubes (MWCNTs) and found the thermal conductivity of the polymer/MWCNT complex to be increased by 280% at 0.3% weight fraction. Thus, the amount of enhancement that is obtained depends on whether the carbon nanotubes (CNTs) are vertically aligned (280%) or randomly orientated (125%). Kwon et al.34 studied the thermal conductivity of MWCNT/polymer composites as a function of the MWCNT volume fraction and found enhancements up to 390% for a MWCNT volume fraction of 1.4%. MWCNTs with very high aspect ratios were used (>2500). Bonnet et al.35 also studied the thermal conductivity enhancement of polymethylmethacrylate using SWCNTs but found enhancements up to only 55% for 7% SWCNT volume fraction. This could be due to the randomly orientated CNTs, which decreases the net thermal conductivity of the sample, and also because of the low aspect ratio of the CNTs used in their study. Liao et al.36 and Lee et al.37 performed MD simulations to predict the thermal conductivity enhancements using CNTs. Liao et al.36 investigated the thermal conductivity of aligned carbon nanotube-polyethylene composites (ACPCs) and found enhancements up to 2 times. They attributed such a considerable enhancement to the alignment of the CNTs and the polymer and to the high thermal conductivity of the CNTs. Lee et al.37 investigated the effect of the SWCNTs on the thermal conductivity of water and found enhancements up to 370% for a 7% volume fraction.

The performance of solid-propellant micro-thrusters and other micro-power devices largely depends on the control and improvement of the combustion wave propagation velocities of the solid monopropellants. Carbon-based materials have been shown to enhance the burning rate of various solid monopropellants.4,40–42 In most of the solid propellant micro-thruster systems, confinement is used for the anisotropic combustion release.38,39 Choi et al.40 proposed anisotropic flame speed enhancements in a TNA (trinitramine) solid-monopropellant by coupling the exothermic reactions to a thermally conductive base such as MWCNTs. Flame speed enhancements up to 104 times were achieved. Jain et al.4,41 used graphite sheets, graphene foam, and graphene nano-pellets as thermally conductive substrates for a nitrocellulose monopropellant and found flame speed enhancements up to 8 times. Jain et al.4 also performed simple 1-D modelling and revealed the mechanisms responsible for the flame speed enhancements. Similar to the observations made by Choi et al.,40 the thermal conductivity of the graphite sheet was regarded as the main parameter governing the flame speed enhancement. Zhang et al.42 also studied flame speed enhancement of nitrocellulose with the addition of graphene oxide pellets. Flame speed enhancements up to 7 times the bulk value were reported.

Motivated by the above discoveries, this work investigated the flame propagation process of Pentaerythritol Tetranitrate (PETN) coupled with MWCNTs using reactive molecular dynamics simulations to examine the physical and the chemical changes occurring at the atomic scale. This is in contrast to the previous works,4,40 where large-scale 1-D modeling was used. Although a number of experiments have been performed to study the flame propagation process, a complete atomic-level understanding of the thermal transport at nanoscale especially at the interface of the two materials is still missing. Thus, in this work, the first-principles based reactive force field ReaxFF was used for simulating the chemical reactions and to understand the coupling between the chemical and the physical processes. The average flame speeds were determined as a function of the MWCNT loading ratio (%), and an optimum loading ratio (%) was found corresponding to the maximum enhancement. The flame propagation process of PETN on CNTs was compared to that of pure PETN. The fundamental mechanism for the flame speed enhancement was identified on a molecular scale. Lastly, the effect of using highly conductive MWCNTs on the minimum ignition energy of PETN was also studied.

For the solid fuel, PETN (C5H8N4O12) was selected because of its high enthalpy of combustion and wide use as a powerful secondary explosive material. Moreover, a number of experiments and simulations of pure PETN as a function of the initial hydrostatic pressure have been performed,43–45 thus providing a base case check for the MD simulations performed in this study. A high pressure of 3 GPa was chosen as the initial hydrostatic pressure in order to simulate the combustion process under shock compression conditions.43 

The simulations were performed using Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS),46 an open source MD simulation code developed by Sandia National Lab. The interactions between the atoms were calculated using ReaxFF interaction potential, which was initially developed by Duin et al.47 and was implemented within LAMMPS using the USER-REAXC module.48 The bond order was determined, as a first approximation, between a pair of atoms from their interatomic distance, using which the potential energy of the whole system was calculated. The particular ReaxFF force-field used in this study was developed by Budzien et al.49 for PETN using DFT (density function theory) calculations for specific decomposition reactions in PETN and cold compression curve.

The unit cell of a PETN crystal50 (space group symmetry P421C) contains two molecules and has a size of 9.2759 Å × 9.2759 Å × 6.6127 Å. The center of one of the molecules is at the center of the unit cell, whereas the center of the other molecule coincides with front lower left vertex. The resulting unit cell, as shown in Fig. 1(a), was then repeated in x, y, and z directions, giving the final simulation domain. The pure PETN case consist of 1075 unit cells (43 × 5 × 5) with the long side of the simulation cell being 40 nm and the total number of atoms in the computational domain being 62 350. Figs. 1(b) and 1(c) show the isotropic and front views of the final simulation domain, respectively. The PETN crystal orientation has a major effect on its detonation51 and combustion properties.52 In the present work, flame propagation was carried out along the [1,0,0] direction.

FIG. 1.

(a) PETN unit cell. (b) Isotropic view and (c) front view of the pure PETN simulation domain.

FIG. 1.

(a) PETN unit cell. (b) Isotropic view and (c) front view of the pure PETN simulation domain.

Close modal

For the PETN-MWCNT complex, a cylinder was cut in the center of the PETN computational domain along the [1,0,0] direction and a MWCNT was then inserted into the hollow space. A Matlab code was written to generate the simulation data file, which made sure that only the full PETN molecules were deleted and not the individual atoms. The initial thickness of the annular PETN layer was varied from 6.6 Å to 13.5 Å, whereas the initial MWCNT's outer diameter was varied from 16 Å (20,0) to 27 Å (35,0) with its initial inner diameter being 3.914 Å (5,0). This was done to vary the MWCNT loading ratio (%), which ranged from 24% to 64%. The initial wall-to-wall distance between two tubes of a MWCNT was 1.957 Å with the individual tubes being achiral zig-zag type. However, as explained in Section II D, the MWCNT's diameter, the wall-to-wall thickness, and the PETN layer thickness will change after performing the pressure and temperature equilibration. The length of the PETN-MWCNT complex was kept fixed at around 40 nm to be consistent with the pure PETN case. Figs. 2(a) and 2(b) show the front and the isotropic view of the simulation domain of a typical PETN-MWCNT complex, respectively.

FIG. 2.

(a) Front and (b) isotropic views of the PETN-MWCNT complex. The particular case shown corresponds to the MWCNT with 15.6 Å (initial) outer diameter and 3.914 Å (initial) inner diameter. The PETN layer thickness was 9.4 Å (initial). The PETN-MWCNT complex length was 40 nm.

FIG. 2.

(a) Front and (b) isotropic views of the PETN-MWCNT complex. The particular case shown corresponds to the MWCNT with 15.6 Å (initial) outer diameter and 3.914 Å (initial) inner diameter. The PETN layer thickness was 9.4 Å (initial). The PETN-MWCNT complex length was 40 nm.

Close modal

After the initial set-up of the computational domain, as shown in Figs. 1 and 2, the temperature equilibration under the NVT (constant number of particles (N), temperature (T), and volume (V)) conditions to 300 K using the Nosé-Hoover temperature thermostat53,54 was performed. The NVT equilibration was run for 70 ps with the time step and the relaxation time being 0.1 fs and 10 fs, respectively. After the NVT equilibration, the system was equilibrated to the desired pressure of 3 GPa under the NPT (constant number of particles (N), temperature (T), and pressure (P)) conditions. The relaxation time for the Nosé-Hoover thermostat and the barostat55 was set to 10 fs with the time step again being 0.1 fs. Only the y and the z dimensions of the simulation box were altered under the NPT equilibration, while the x dimension was held constant. The NPT equilibration was run for 50 ps. The final MWCNT diameter and the PETN layer thickness for the 6 different cases that were run, corresponding to different MWCNT loading ratios (%), are listed in Table I.

TABLE I.

MD simulation matrix.

Case No.SystemMWCNT loading ratio (%)No. of tubes in the MWCNTFinal PETN thicknessFinal MWCNT diameter
Pure PETN … … … 
PETN-MWCNT 24 11.35 Å 18 Å 
PETN-MWCNT 35 8.8 Å 18 Å 
PETN-MWCNT 50 5 Å 18 Å 
PETN-MWCNT 57 5.55 Å 26 Å 
PETN-MWCNT 64 4.8 Å 31 Å 
Case No.SystemMWCNT loading ratio (%)No. of tubes in the MWCNTFinal PETN thicknessFinal MWCNT diameter
Pure PETN … … … 
PETN-MWCNT 24 11.35 Å 18 Å 
PETN-MWCNT 35 8.8 Å 18 Å 
PETN-MWCNT 50 5 Å 18 Å 
PETN-MWCNT 57 5.55 Å 26 Å 
PETN-MWCNT 64 4.8 Å 31 Å 

After the NPT equilibration, the NVE ensemble was used to simulate the flame propagation process. Periodic boundary conditions were employed in all the three directions, and a time step of 0.1 fs was used as suggested by Sergeev.43 To simulate an ignition source, only the temperature of the PETN molecules in the middle of the computational domain, being 5-unit cell thick in the X direction, was set to 4000 K for both the PETN-MWCNT case and the pure PETN, while the temperature of the rest of the molecules remained at 300 K. The ignition temperature was based on the minimum value needed for the flame propagation in the PETN-MWCNT complexes. Since the minimum ignition temperature was found to increase for the PETN molecules when coupled to MWCNTs as compared to the pure PETN, the minimum ignition temperature corresponding to that of the PETN-MWCNT complexes was used for all the cases simulated, in order to have the same ignition energy input. To identify the species produced in the combustion process, the bond order minimum values as listed by Budzien et al.49 were used. A pair of atoms that had the bond order value greater than the threshold bond order value (Table II) was considered to be bonded.

TABLE II.

Minimum bond order values.49 

Atom typeAtom typeBond order
0.55 
0.40 
0.30 
0.65 
0.55 
0.55 
0.40 
0.55 
0.40 
0.65 
Atom typeAtom typeBond order
0.55 
0.40 
0.30 
0.65 
0.55 
0.55 
0.40 
0.55 
0.40 
0.65 

Figure 3 shows the spatial temperature profile along the x-direction at different times during a typical pure PETN combustion. The flame propagates from the ignition zone in the middle towards the two ends. The initial ignition zone, as stated above, was set to 5-unit cell thick in the x-direction, which corresponds to approximately 5 nm. The reaction front, which is required to calculate the flame speed, was identified using 2 different criteria. The first criterion is based on the location of the peak NO2 concentration. NO2 is one of the key intermediate species formed during the PETN decomposition; thus, the position of the reaction zone could be identified from the NO2 formation. For this criterion, the simulation domain was divided into 15 slabs in the X direction with each slab being 2.7 nm or 3-unit cell thick. For the second criterion, the reaction zone was identified by tracking the temperature changes along the PETN sample. The location of three different temperatures (1000, 1500, and 2000 K, respectively) along the x-direction as a function of time was used as an indicator of the reaction front. The local temperature at a specific location was computed by considering all the atoms in the local layer that had a thickness in the X direction corresponding to the 1-unit cell. The temperature was determined by the average kinetic energy of the particles in the local layer.

FIG. 3.

Spatial temperature profile of pure PETN during the flame propagation process at different times.

FIG. 3.

Spatial temperature profile of pure PETN during the flame propagation process at different times.

Close modal

Figure 4(a) plots the location of the flame front as a function of time using three different temperatures. The three curves are nearly parallel to each other, giving an average flame propagation speed of 106 ± 10 m/s. This indicates that the selection of the temperature as an indicator of the flame front does not influence the calculated flame speed. Because of the transient nature of the ignition process, the slabs placed closed to the ignition zone were ignored and the temperature data were only gathered from the slabs with at least 2 slab thickness away from the ignition zone, after which a uniform flame propagation was obtained with every point on the temperature profile traveling at the same speed. Figure 4(b) shows the location of the peak NO2 concentration as a function of time. Linear fitting to data gave an average flame propagation speed of 110 ± 5 m/s, similar to the speeds obtained using the temperature criterion. In other words, both criteria produce nearly identical flame speeds. The computed flame speed for pure PETN was compared to the experimental data obtained by Andreev44 and Foltz45 along with the simulated values obtained by Sergeev.43 As shown in Figure 5, good agreement is obtained.

FIG. 4.

Flame speed determined from (a) three different temperature profiles and (b) peak NO2 concentration.

FIG. 4.

Flame speed determined from (a) three different temperature profiles and (b) peak NO2 concentration.

Close modal
FIG. 5.

Comparison of the flame speed of pure PETN with other computational and experimental data.

FIG. 5.

Comparison of the flame speed of pure PETN with other computational and experimental data.

Close modal

Figure 6 shows the species profiles of reactants, intermediate species, and final products as a function of time for the pure PETN. In the figure, XM represents the number of molecules of a species per initial PETN molecules in a given slab (3-unit cell thick) located at a distance of 4 nm. The results show that the NO2 is the dominant initial product formed during the PETN decomposition. This indicates that the NO2-O bond in the PETN molecule, which has the weakest bond energy, breaks first during initiation. HNO2 and NO are the other major intermediate products formed after NO2. The water formation begins after the appearance of OH and NO, which is followed by an increase in the concentrations of CO2 and CO. Thus, the temporal species evolution from the present MD simulation is consistent with that shown by Sergeev,43 except that the peak NO concentration obtained in the present study was greater than the NO2 peak concentration. This could be attributed to the different minimum bond order values used, which could redistribute the identification of molecules between NO2 and NO. For example, Sergeev43 used a minimum bond order value of 0.5 between O and O atoms, comparing this to the present study where a minimum bond order value of 0.65 was used. Thus, a molecule identified as NO or CO in the present study could have been identified as NO2 or CO2 in the work performed by Sergeev.43 

FIG. 6.

Species distribution as a function of time, at a location = 4 nm for PETN combustion.

FIG. 6.

Species distribution as a function of time, at a location = 4 nm for PETN combustion.

Close modal

Figure 7 shows the propagation of the combustion wave along a PETN-MWCNT complex (case 2). As can be seen, the MWCNTs remain intact during the combustion process although their structure in the burned zone is modified from twisting and bending.

FIG. 7.

Combustion wave propagation along the PETN-MWCNT. Case 2 is shown. The ignition zone is 5 nm.

FIG. 7.

Combustion wave propagation along the PETN-MWCNT. Case 2 is shown. The ignition zone is 5 nm.

Close modal

Figure 8 compares the temperature profile of the PETN-MWCNT complex with that of the pure PETN at different times after ignition. The equilibrium temperature of the PETN-MWCNT complex was lower than that of the pure PETN because some of the energy released during the exothermic reaction was used in heating the CNTs, which acted as heat sinks and thus reduced the equilibrium temperature of the PETN molecules. At t = 5 ps, the temperature of the unburned PETN close to the reaction zone remained at 300 K for pure PETN, whereas it increased to about 700 K for the PETN-MWCNT. At t = 10 ps, the unburned PETN temperature (near the reaction zone) was raised to around 900 K for the PETN-MWCNT, whereas it still remained at 300 K for the pure PETN. Thus, comparing the temperature profiles of the two cases, the PETN-MWCNT has a much wider reaction and a pre-heat zone. This can be attributed to the high thermal conductivity of the PETN-MWCNT composite,5–10 which leads to faster heat propagation and thus more unburned portions of the fuel are heated ahead of the reaction front.

FIG. 8.

Comparing the spatial temperature profiles of pure PETN and PETN-MWCNT at different times.

FIG. 8.

Comparing the spatial temperature profiles of pure PETN and PETN-MWCNT at different times.

Close modal

Similar to the pure PETN case, the reaction front was identified using 2 different criteria, i.e., the peak NO2 concentration and fixed temperatures. The thickness of the local layer used to calculate the local temperature and the local peak NO2 concentration was also kept the same as that of the pure PETN case, i.e., 1-unit cell and 3-unit cell thick, respectively. Figure 9(a) plots the location of the flame front as a function of time using two different temperatures. The two curves are nearly parallel to each other, giving an average flame propagation speed of 320 ± 10 m/s. Figure 9(b) shows the NO2 peak location as a function of time giving an average flame speed of 330 ± 10 m/s. Thus, nearly identical flame speeds are obtained from either criterion. The particular simulation shown in Figure 9 corresponds to the case 4 (Table I), which was the optimum case as will be shown later.

FIG. 9.

Flame speed determined from (a) 2 different temperature profiles and (b) peak NO2 concentration (case 4).

FIG. 9.

Flame speed determined from (a) 2 different temperature profiles and (b) peak NO2 concentration (case 4).

Close modal

Figure 10 shows the species profiles for the PETN-MWCNT (case 4) as a function of time at a particular location of 4 nm. The location was kept the same as that of the pure PETN in order to facilitate one-to-one comparison. Again, in Fig. 10, XM represents the number of molecules of a species per initial PETN molecules, in a given slab. Comparing Figs. 6 and 10, the species distribution curve is shifted to the left for the PETN-MWCNT, but a similar reaction path was obtained. Again, NO2 was the dominant initial product formed during the PETN decomposition. HNO2 and NO were other major intermediate products formed after NO2 with the water molecules being formed after OH and NO appearance. Thus, the species distribution of the PETN decomposition remains unchanged after the addition of MWCNTs. However, the rate of production of CO2 and CO was slightly less in PETN-MWCNTs as compared to that in pure PETN, which could be attributed to the lower flame temperature. Some of the heat released during the exothermic reactions of the PETN molecules was used in heating the MWCNTs' carbon atoms, thus decreasing the flame temperature. The fact that the MWCNT remains unburned during the PETN combustion could be confirmed by looking at the mole fractions of CO2 and CO (Fig. 10). Since CO2 and CO do not appear until after H2O formation, the temporal evolution for these oxides is consistent with that of the pure PETN. Moreover, their peak mole fractions are on the same order of magnitude as that observed in pure PETN. Thus, the MWCNT is not consumed during the PETN combustion and only acts to increase the layering of the PETN molecules along its surface that facilitates the transfer of heat from the burned to the unburned portions of the fuel much faster, which in turn causes the species distribution curves to shift to the left. The deformation of the MWCNTs, as shown in Fig. 7, was from their increased temperature resulting from the PETN combustion.

FIG. 10.

Species distribution as a function of time, location = 4 nm for PETN-MWCNT (case 4).

FIG. 10.

Species distribution as a function of time, location = 4 nm for PETN-MWCNT (case 4).

Close modal

Five different PETN-MWCNT combinations were simulated (cases 2 to 6), as shown in Table I, in order to determine the effect of the MWCNT loading ratio (%) on the flame speed enhancement. The results are shown in Fig. 11. The MWCNT loading ratio (%) is defined as the mass of the MWCNT per total mass of the system (PETN + MWCNT). For very thick PETN layers, a large amount of energy was released, but heat transfer among the PETN molecules near the MWCNT surface was not high enough to conduct the heat efficiently from the exothermic reactions to aid in reaction propagation. On the contrary for very thin PETN layers, although the thermal transport among the PETN molecules was enhanced, the amount of heat reaching the unburned portions of the fuel was substantially reduced because some of the energy released during the exothermic reaction was used in heating the MWCNTs, which acted as heat sinks and thus lowered the reaction propagation speeds. Consequently, an optimum loading ratio (%) exists. An optimal loading ratio of around 55% was found for which the flame speed enhancement was around 3 times the bulk speed of 110 cm/s.

FIG. 11.

The effect of the MWCNT loading ratio (%) on the average flame speeds.

FIG. 11.

The effect of the MWCNT loading ratio (%) on the average flame speeds.

Close modal

In addition to the reactive MD simulations, two additional non-reactive molecular dynamics (MD) simulations were conducted to better understand the mechanisms contributing to the thermal conductivity enhancement of the composite and in turn the flame speed enhancement.

First, a non-reactive reverse non-equilibrium MD simulation (RNEMD) was conducted using LAMMPS to investigate the interfacial heat transfer in the PETN-MWCNT composite. The MD study conducted was based on the procedure outlined in the studies performed previously by Zahedi et al.56 and Alaghemandi et al.57 In the RNEMD approach, a constant heat flux is imposed on the simulation box by performing velocity exchanges between the coldest particle (from the hot layer) and the hottest particle (from the cold layer) in a given direction. If the masses of the particles being exchanged are different, then an exchange of velocities relative to the center of mass motion of the two atoms is performed, to conserve the total kinetic energy of the system. The RNEMD simulation was performed under the NVE (constant volume and energy) conditions at a chosen temperature of 330 K. The system was first equilibrated to a temperature and a pressure of 330 K and 3 GPa, respectively. The relaxation time for the Nosé-Hoover thermostat and the barostat was set to 10 fs with the time step being 0.2 fs. After the equilibration, a constant heat flux was applied using the Muller-Plathe algorithm58 under the NVE conditions. Sufficient energy and temperature conservation were obtained using the timestep of 0.2 fs. At higher timesteps, deviations in the total energy were observed. In this study, the velocity exchanges were performed between the CNT atoms (located in the slab in the middle of the simulation box) and the PETN atoms (located in the slab at the maximum separation in the y direction from the CNT atoms). The simulation box was divided into 16 slabs, 0.271 nm thick, in the direction of the heat flux (y-axis) and 13 slabs, 0.33 nm thick, in the direction perpendicular to the heat flux (z-axis). Moreover, the velocity exchanges were performed between the two atoms every 20 fs. The heat flux as computed by the LAMMPS fix thermal conductivity command is given by57 

jy=12tA(mhotVhot2mcoldVcold2)2.
(1)

In the above equation, mhot and mcold are the masses of the hot and the cold particle, respectively, whose velocity is being exchanged, A is the cross-sectional area perpendicular to the heat flux direction (z-x), vhot and vcold are the velocities of the hot and the cold particle, respectively, and t is the total simulation time. A factor of 2 is needed in Eq. (1) because of the periodic boundary conditions used in the direction of the heat flux.58 From the imposed heat flux, the thermal conductivity value can be obtained as follows:57 

ky=jydTdy.
(2)

In the above equation, ky is the thermal conductivity value in the y-direction and dT/dy is the temperature gradient due to the imposed heat flux. The z-direction was divided into 13 slabs of 0.33 nm thickness, and the ky value was obtained by looking at the temperature gradient (dT/dy) for a z-slab located at (y,0).

Figure 12 shows a typical temperature profile in the y-direction for the z-slab located at (y,0). As can be seen, the temperature profile is linear in the individual regions belonging to PETN, interface, and MWCNTs. From the linear temperature profiles, the thermal conductivity values for each region can be calculated using Eq. (2). The heat flux value is the same for all the regions. An effective thermal conductivity of 0.172 ± 7% (W/m K), 0.045 ± 5% (W/m K), and 0.7 ± 10% (W/m K) was obtained for PETN, interface, and MWCNTs, respectively. The interface thermal conductivity value obtained was 4 times lower than that of the PETN, which could be attributed to the mismatch of the thermal transport regimes in PETN and MWCNTs.57,59 In MWCNTs, the heat is transferred through the ballistic regime, whereas, in PETN, the thermal transport occurs in the diffusive regime.60 This sudden transition from the ballistic to the diffusive regime limits the net thermal conductivity enhancement of the composite.57,59

FIG. 12.

(a) The slabs in the y and z directions. (b) Temperature profile in the Y-direction (case 2).

FIG. 12.

(a) The slabs in the y and z directions. (b) Temperature profile in the Y-direction (case 2).

Close modal

Alaghemandi et al.57 investigated the thermal conductivity of composites of single-walled carbon nanotubes and polymamide-6,6 (PA) using reverse non-equilibrium MD simulations and found the interface thermal conductivity value to be only 0.003 W/m K, which was 1–2 orders of magnitude lower than the thermal conductivity of pure PA (0.24 W/m K). The interface thermal conductivity value of 0.045 W/m K obtained in this work is an order of magnitude higher than the interface thermal conductivity value of 0.003 W/m K obtained by Alaghemandi et al.57 The difference could be attributed to different materials and simulation conditions used. The present simulations were conducted at an extremely high pressure of 3 GPa, as opposed to the atmospheric pressures in the simulations performed by Alaghemandi et al.57 Because of such a high thermal interface resistance, there must be a different mechanism responsible for the increased global thermal conductivity of the composite.

Zahedi et al.56 investigated the structural properties of the polymer matrix around the CNTs. A highly ordered polymer matrix structure was observed in the interphase region as quantified using the normalized density profiles around the CNTs. They concluded that because of this wrapping of the PA molecules around the CNTs, the PA molecules are predominantly tangential to the CNT surface, which increases the heat transport along the CNTs but decreases the heat transport in the perpendicular direction. Motivated by this, the layering of the PETN molecules as a result of their interactions with MWCNTs was also examined. An equilibrium non-reactive MD simulation was conducted under the NVE conditions. Again, the system was first equilibrated to a temperature and a pressure of 300 K and 3 GPa, respectively. The relaxation time for the Nosé-Hoover thermostat and the barostat was set to 10 fs along with a timestep of 0.2 fs. After the equilibration, the density profile calculations were performed under the NVE conditions. The simulation box was divided into cylindrical bins having a length of 4 nm and a radial thickness of 0.07 nm.

Figure 13 shows a typical normalized density profile of the PETN molecules around the MWCNT (case 2). As can be seen, the PETN molecules are indeed ordered around the MWCNT. This organized interface structure increases the thermal transport in the direction parallel to the CNT surface but decreases the thermal transport in the direction perpendicular to it57 and thus contributes to the net thermal conductivity enhancement of the PETN-MWCNT composite.

FIG. 13.

Normalized density of the PETN molecules around the MWCNT (case 2).

FIG. 13.

Normalized density of the PETN molecules around the MWCNT (case 2).

Close modal

In this section, the effect of adding MWCNTs to PETN on the minimum ignition energy required to initiate successful flame propagation along the PETN sample was examined. To achieve this goal, the temperature of the ignition zone was varied with its length unchanged for both pure PETN and PETN-MWCNT cases. Figure 14 plots the average flame speeds as a function of various ignition temperatures in the range of 3000–5000 K. The minimum ignition temperature is defined as the temperature below which the flame propagation could not be sustained and the system eventually cools down. The minimum ignition temperature for the PETN molecules was found to increase from 3000 K to 4000 K when coupled to MWCNTs. This was again attributed to the high thermal transport among the PETN molecules near the MWCNT surface, which resulted in a faster heat dissipation (or heat loss) and thus a higher minimum ignition temperature was required. Nevertheless, above the minimum ignition temperature, the flame speed values remain unchanged and no over-driven ignition characteristic was observed. Atwood et al.61 suggested that the overdriven condition occurs most often at lower pressures. Since the present simulations were conducted at extremely high pressures (3 GPa), the overdriven phenomenon may not have occurred. Another reason for not observing the over-driven ignition could be that the applied ignition energy was simply not high enough. Atwood et al.61 observed over-driven ignition in gun propellants (at 1.72 MPa) only when the heat flux was increased 3 times.

FIG. 14.

The effect of the ignition temperature on the flame speeds for (a) PETN-MWCNT and (b) pure PETN.

FIG. 14.

The effect of the ignition temperature on the flame speeds for (a) PETN-MWCNT and (b) pure PETN.

Close modal

Reactive MD simulations of flame propagation of a monopropellant (PETN) coupled with a MWCNT were conducted. The thickness of the PETN layer and the MWCNT's diameter were varied to study the effect of the MWCNT loading ratio (%) on the amount of the flame speed enhancement. Flame speed enhancements up to 3 times the bulk value were observed, and an optimal MWCNT loading ratio (%) of around 55% was found. In addition to the reactive MD simulations, two additional non-reactive molecular dynamics (MD) simulations were conducted to better understand the mechanism contributing to the thermal conductivity enhancement of the composite and in turn the flame speed enhancement. The enhancement was attributed to the layering of the PETN molecules along the MWCNT surface, which resulted in the faster heat conduction in the PETN molecules, thus causing the flame to travel faster. Moreover, the PETN-MWCNT complex requires higher minimum ignition energy than pure PETN to initiate successful flame propagation, where the minimum ignition temperature for the PETN molecules was found to increase from 3000 K to 4000 K when coupled to MWCNTs. Lastly, the temporal distribution of the species was also studied, which confirmed that the MWCNT remained unburned during the PETN combustion.

This research was supported by the Air Force Office of Scientific Research (AFOSR) with Dr. Chiping Li as the technical monitor.

1.
P.
Kim
,
L.
Shi
,
A.
Majumdar
, and
P. L.
McEuen
,
Phys. Rev. Lett.
87
,
215502
(
2001
).
2.
L. M.
Viculis
,
J. J.
Mack
,
O. M.
Mayer
,
H. T.
Hahn
, and
R. B.
Kaner
,
J. Mater. Chem.
15
,
974
(
2005
).
3.
Z.
Chen
,
W.
Ren
,
L.
Gao
,
B.
Liu
,
S.
Pei
, and
H.-M.
Cheng
,
Nat. Mater.
10
,
424
(
2011
).
4.
S.
Jain
,
O.
Yehia
, and
L.
Qiao
,
J. Appl. Phys.
119
,
094904
(
2016
).
5.
W.
Cai
,
A. L.
Moore
,
Y.
Zhu
,
X.
Li
,
S.
Chen
,
L.
Shi
, and
R. S.
Ruoff
,
Nano Lett.
10
,
1645
(
2010
).
6.
S.
Chen
,
A. L.
Moore
,
W.
Cai
,
J. W.
Suk
,
J.
An
,
C.
Mishra
,
C.
Amos
,
C. W.
Magnuson
,
J.
Kang
,
L.
Shi
, and
R. S.
Ruoff
,
ACS Nano
5
,
321
(
2011
).
7.
S.
Ghosh
,
W.
Bao
,
D. L.
Nika
,
S.
Subrina
,
E. P.
Pokatilov
,
C. N.
Lau
, and
A. A.
Balandin
,
Nat. Mater.
9
,
555
(
2010
).
8.
K.
Sun
,
M. A.
Stroscio
, and
M.
Dutta
,
J. Appl. Phys.
105
,
074316
(
2009
).
9.
M. K.
Samani
,
N.
Khosravian
,
G. C. K.
Chen
,
M.
Shakerzadeh
,
D.
Baillargeat
, and
B. K.
Tay
,
Int. J. Therm. Sci.
62
,
40
(
2012
).
10.
S.
Berber
,
Y.-K.
Kwon
, and
D.
Tománek
,
Phys. Rev. Lett.
84
,
4613
(
2000
).
11.
T.
Tong
,
Y.
Zhao
,
L.
Delzeit
,
A.
Kashani
,
M.
Meyyappan
, and
A.
Majumdar
,
IEEE Trans. Compon. Packag. Technol.
30
,
92
(
2007
).
12.
J.
Xu
and
T. S.
Fisher
,
Int. J. Heat Mass Transfer
49
,
1658
(
2006
).
13.
M. A.
Panzer
,
G.
Zhang
,
D.
Mann
,
X.
Hu
,
E.
Pop
,
H.
Dai
, and
K. E.
Goodson
,
J. Heat Transfer
130
,
052401
(
2008
).
14.
W.
Park
,
J.
Hu
,
L. A.
Jauregui
,
X.
Ruan
, and
Y. P.
Chen
,
Appl. Phys. Lett.
104
,
113101
(
2014
).
15.
M. J.
Biercuk
,
M. C.
Llaguno
,
M.
Radosavljevic
,
J. K.
Hyun
,
A. T.
Johnson
, and
J. E.
Fischer
,
Appl. Phys. Lett.
80
,
2767
(
2002
).
16.
V.
Goyal
and
A. A.
Balandin
,
Appl. Phys. Lett.
100
,
073113
(
2012
).
17.
H.
Huang
,
C. H.
Liu
,
Y.
Wu
, and
S.
Fan
,
Adv. Mater.
17
,
1652
(
2005
).
18.
A.
Yu
,
P.
Ramesh
,
M. E.
Itkis
,
E.
Bekyarova
, and
R. C.
Haddon
,
J. Phys. Chem. C
111
,
7565
(
2007
).
19.
F.
Yavari
,
H. R.
Fard
,
K.
Pashayi
,
M. A.
Rafiee
,
A.
Zamiri
,
Z.
Yu
,
R.
Ozisik
,
T.
Borca-Tasciuc
, and
N.
Koratkar
,
J. Phys. Chem. C
115
,
8753
(
2011
).
20.
K. M. F.
Shahil
and
A. A.
Balandin
,
Nano Lett.
12
,
861
(
2012
).
21.
A.
Yu
,
P.
Ramesh
,
X.
Sun
,
E.
Bekyarova
,
M. E.
Itkis
, and
R. C.
Haddon
,
Adv. Mater.
20
,
4740
(
2008
).
22.
H.
Ji
,
D. P.
Sellan
,
M. T.
Pettes
,
X.
Kong
,
J.
Ji
,
L.
Shi
, and
R. S.
Ruoff
,
Energy Environ. Sci.
7
,
1185
(
2014
).
23.
X.
Zhang
,
K. K.
Yeung
,
Z.
Gao
,
J.
Li
,
H.
Sun
,
H.
Xu
,
K.
Zhang
,
M.
Zhang
,
Z.
Chen
,
M. M. F.
Yuen
, and
S.
Yang
,
Carbon
66
,
201
(
2014
).
24.
Z.
Liu
,
D.
Shen
,
J.
Yu
,
W.
Dai
,
C.
Li
,
S.
Du
,
N.
Jiang
,
H.
Li
, and
C.-T.
Lin
,
RSC Adv.
6
,
22364
(
2016
).
25.
M.-T.
Hung
,
O.
Choi
,
Y. S.
Ju
, and
H. T.
Hahn
,
Appl. Phys. Lett.
89
,
023117
(
2006
).
26.
Y.-H.
Zhao
,
Z.-K.
Wu
, and
S.-L.
Bai
,
Composites, Part A
72
,
200
(
2015
).
27.
M. B.
Bryning
,
D. E.
Milkie
,
M. F.
Islam
,
J. M.
Kikkawa
, and
A. G.
Yodh
,
Appl. Phys. Lett.
87
,
161909
(
2005
).
28.
S.
Harish
,
D.
Orejon
,
Y.
Takata
, and
M.
Kohno
,
Appl. Therm. Eng.
80
,
205
(
2015
).
29.
M. T.
Pettes
,
H.
Ji
,
R. S.
Ruoff
, and
L.
Shi
,
Nano Lett.
12
,
2959
(
2012
).
30.
L.
Chen
,
R.
Zou
,
W.
Xia
,
Z.
Liu
,
Y.
Shang
,
J.
Zhu
,
Y.
Wang
,
J.
Lin
,
D.
Xia
, and
A.
Cao
,
ACS Nano
6
,
10884
(
2012
).
31.
J.-N.
Shi
,
M.-D.
Ger
,
Y.-M.
Liu
,
Y.-C.
Fan
,
N.-T.
Wen
,
C.-K.
Lin
, and
N.-W.
Pu
,
Carbon
51
,
365
(
2013
).
32.
K. H.
Baloch
,
N.
Voskanian
,
M.
Bronsgeest
, and
J.
Cumings
,
Nat. Nanotechnol.
7
,
316
(
2012
).
33.
Z.
Yan
,
G.
Liu
,
J. M.
Khan
, and
A. A.
Balandin
,
Nat. Commun.
3
,
827
(
2012
).
34.
S. Y.
Kwon
,
I. M.
Kwon
,
Y.-G.
Kim
,
S.
Lee
, and
Y.-S.
Seo
,
Carbon
55
,
285
(
2013
).
35.
P.
Bonnet
,
D.
Sireude
,
B.
Garnier
, and
O.
Chauvet
,
Appl. Phys. Lett.
91
,
201910
(
2007
).
36.
Q.
Liao
,
Z.
Liu
,
W.
Liu
,
C.
Deng
, and
N.
Yang
,
Sci. Rep.
5
,
16543
(
2015
).
37.
J. W.
Lee
,
J. A. J.
Meade
,
E. V.
Barrera
, and
J. A.
Templeton
,
J. Heat Transfer
137
,
072401
(
2015
).
38.
E.
Rudnyi
,
T.
Bechtold
,
J.
Korvink
, and
C.
Rossi
, in
NanoTech 2002 - “At the Edge of Revolution”
(
American Institute of Aeronautics and Astronautics
,
Houston, Texas
,
2002
).
39.
Z.
Kaili
,
S. K.
Chou
, and
S. S.
Ang
,
J. Microelectromech. Syst.
13
,
165
(
2004
).
40.
W.
Choi
,
S.
Hong
,
J. T.
Abrahamson
,
J.-H.
Han
,
C.
Song
,
N.
Nair
,
S.
Baik
, and
M. S.
Strano
,
Nat. Mater.
9
,
423
(
2010
).
41.
S.
Jain
,
W.
Park
,
Y. P.
Chen
, and
L.
Qiao
,
J. Appl. Phys.
120
,
174902
(
2016
).
42.
X.
Zhang
,
W. M.
Hikal
,
Y.
Zhang
,
S. K.
Bhattacharia
,
L.
Li
,
S.
Panditrao
,
S.
Wang
, and
B. L.
Weeks
,
Appl. Phys. Lett.
102
,
141905
(
2013
).
43.
O. V.
Sergeev
and
A. V.
Yanilkin
,
Combust., Explos. Shock Waves
50
,
323
(
2014
).
44.
K. K.
Andreev
,
Thermal Decompistion and Combustion of Explosives
, revised ed. (
Foreign Technology Div, Wright-Patterson AFB
,
OH/Nauka/Moscow
,
1966
), Vol.
2
.
45.
M. F.
Foltz
,
Pressure dependence on the reaction propagation rate of PETN at high pressure
, Boston, MA, UCRL-JC-111316, Lawrence Livermore National Laboratory,
1993
.
46.
S.
Plimpton
,
J. Comput. Phys.
117
,
1
(
1995
).
47.
A. C. T.
van Duin
,
S.
Dasgupta
,
F.
Lorant
, and
W. A.
Goddard
,
J. Phys. Chem., A
105
,
9396
(
2001
).
48.
H. M.
Aktulga
,
J. C.
Fogarty
,
S. A.
Pandit
, and
A. Y.
Grama
,
Parallel Comput.
38
,
245
(
2012
).
49.
J.
Budzien
,
A. P.
Thompson
, and
S. V.
Zybin
,
J. Phys. Chem. B
113
,
13142
(
2009
).
50.
G. R.
Miller
and
A. N.
Garroway
,
A Review of the Crystal Structures of Common Explosives Part I: RDX, HMX, TNT, PETN, and Tetryl
(
Naval Research Laboratory
,
Washington, DC
,
2001
).
51.
J. J.
Dick
,
Appl. Phys. Lett.
44
,
859
(
1984
).
52.
S.
Tzu-Ray
and
P. T.
Aidan
,
J. Phys.: Conf. Ser.
500
,
172009
(
2014
).
53.
W. G.
Hoover
,
Phys. Rev. A
31
,
1695
(
1985
).
54.
S.
Nosé
,
J. Chem. Phys.
81
,
511
(
1984
).
55.
S.
Melchionna
,
G.
Ciccotti
, and
B.
Lee Holian
,
Mol. Phys.
78
,
533
(
1993
).
56.
M. R.
Gharib-Zahedi
,
M.
Tafazzoli
,
M. C.
Bohm
, and
M.
Alaghemandi
,
Phys. Chem. Chem. Phys.
17
,
14502
(
2015
).
57.
M.
Alaghemandi
,
F.
Müller-Plathe
, and
M. C.
Böhm
,
J. Chem. Phys.
135
,
184905
(
2011
).
58.
F.
Müller-Plathe
,
J. Chem. Phys.
106
,
6082
(
1997
).
59.
S. U. S.
Choi
,
Z. G.
Zhang
,
W.
Yu
,
F. E.
Lockwood
, and
E. A.
Grulke
,
Appl. Phys. Lett.
79
,
2252
(
2001
).
60.
J.
Wang
and
J.-S.
Wang
,
Appl. Phys. Lett.
88
,
111909
(
2006
).
61.
A. I.
Atwood
,
K. P.
Ford
,
A. L.
Daniels
,
C. J.
Wheeler
,
P. O.
Curran
,
T. L.
Boggs
, and
J.
Covino
,
Ignition and Combustion Studies of Hazard Division 1.1 and 1.3 Substances
(
Naval Air Warfare Center
,
China Lake, California
,
2010
).