The Anisotropic Chemical Reaction Mechanism of 1,3,3-trinitroazetidine (TNAZ) under Different Shock Wave Directions by ReaxFF Reactive Molecular Dynamics Simulations

1,3,3-Trinitroazetidine (TNAZ) has good thermal stability and low shock sensitivity, among other properties, and it has broad prospects in insensitive ammunition applications. In this study, a molecular dynamics calculation based on the ReaxFF-lg force field and multiscale shock technique (MSST) was used to simulate the shock-induced chemical reaction of TNAZ with different shock wave directions. The results showed that the shock sensitivity of TNAZ was in the order of [100] > [010] > [001]. There were significant differences in molecular arrangements in different shock directions, which affected the reaction rate and reaction path in different directions. The molecular arrangement in the [010] and [001] directions formed a “buffer” effect. The formation and cleavage of bonds, formation of small molecules and growth of clusters were analyzed to show the effect of the “buffer”. The polymerization reactions in the [010] and [001] directions appeared later than that in the [100] direction, and the cluster growth in the [010] and [001] directions was slower than that in the [100] direction. In different shock loading directions, the formation and cleavage mechanisms of the N-O bonds of the TNAZ molecules were different, which resulted in differences in the initial reaction path and reaction rate in the three directions

Many theoretical and experimental studies have been conducted on the thermal decomposition of TNAZ. Archibald obtained the unit cell parameters of TNAZ crystals by X-ray diffraction [1]. Anderson used an STMBMS mass spectrometer to study the thermal decomposition of TNAZ. The results showed that the main initial decomposition products of TNAZ were NO 2 and NO, and the breakage of the N-NO 2 and C-NO 2 bonds was the main path of the initial thermal decomposition of TNAZ [7]. Olah also obtained similar conclusions by mass spectrometry [8]. Behrens studied the process of bond cleavage by simultaneous thermogravimetric modulated beam mass spectrometer. They found that the cleavage of the N-NO 2 group preceded that of the C-NO 2 cleavage in thermal decomposition of TNAZ [9]. Agrawal studied the melting point and isothermal compression properties of defect-free TNAZ crystals by molecular simulation. The unit cell parameters and thermodynamic melting point of TNAZ crystals obtained with the AMBERSTR force field were basically consistent with the experimental values [8]. Although TNAZ explosives have been the subject of some studies and have been used in some applications, few studies on shock reactions mechanism have been conducted to date, and no studies have been reported from the atomic and molecular microscopic levels. The mechanism and characteristics of the shock reactions of TNAZ explosives are not well understood. The reaction mechanism of TNAZ under shock conditions, which is studied using molecule dynamics, is of great significance to the shock response performance and safety evaluation of TNAZ explosives.
The multiscale shock technique (MSST) was proposed by Reed in 2003 [10]. Through the method, researchers can accurately obtain the evolution process of various thermodynamic parameters, which greatly reduces the computational cost. The MSST is an effective means to study the reaction of materials under shock [11][12][13][14][15][16][17][18][19][20][21]. In 2009, Manaa et al. used the MSST to study the shock response of 1,3,5-triamino-2,4,6-trinitrobenzene (TATB) crystals under different shock wave velocity loadings; they found that the formation of many clusters during the shock caused the low shock sensitivity of TATB [22]. In 2019, Liu et al. used the MSST with the ReaxFF molecular dynamics calculation method to study the initial reaction mechanism of ε-CL-20 crystals under shock. They found that the detonation reaction occurred when the cleavage rates of C-N and C-H were greater than 3.11 and 4.15%/ps, respectively [23]. In 2019, Liu et al. studied the initial decomposition mechanism of energetic cocrystal 2,4,6,8,10,12-hexantro-2,4,6,8,10,12-hexazaiso-wurzitane (CL-20)/1,3,5,7-tetranitro-1,3,5,7-tetrazacy-clooctane (HMX) under a steady shock wave. They found that after the shock wave application, the energetic cocrystals successively underwent an induction period, fast compression, slow compression and expansion processes. The fast and slow compression processes corresponded to the fast and slow decomposition of the reactants, respectively [24]. In 2019, Liu et al. studied ReaxFF molecular dynamics simulations of the shock-induced reaction initiation in TNT. They found that the formation of TNT dimers and decomposition to C7H5O5N3 were the dominant initial routes for shockinduced reaction initiation. In addition, they also analyzed the formation and evolution of carbon clusters [25]. In 2020, Huang et al. studied anisotropic hydrogen bond structures and the orientation dependence of shock sensitivity in 1,3,5-tri-amino-2,4,6-tri-nitrobenzene (TATB). They found that the temperature, stress, volume compressibility and decomposition rate of TATB were strongly dependent on the shockwave direction [26]. These studies show that by using the MSST method, researchers can explain the macroscopic phenomena at the atomic and molecular levels and obtain detailed physical properties (pressure and temperature) of the supercell under shock.
The results of previous studies did not discuss the reaction mechanism from the perspective of the system structure after shock compression. This paper discusses the similarities and differences between the reactions in various directions after shock from the perspective of the compressed crystal structure and intermolecular interaction. The shock in different directions leads to differences in the crystal structure and molecular structure. These differences lead to differences in the physical properties of supercells, which affect the subsequent chemical reaction pathways and products. In this paper, ReaxFF-lg [27,28] and MSST are used to study the molecular dynamics of the TNAZ shock reaction process. TNAZ in the [100], [010] and [001] directions at velocities of 8, 9, 10 and 11 km·s −1 is studied. TNAZ supercells at different shock loading rates are examined. The particle velocity, system pressure, shock Hugoniot relationship, the effect of the shock loading velocity on the reaction path of the explosive and the bond information are analyzed.

Computational Methods
The simulations were performed with LAMMPS using ReaxFF-lg and MSST. The single-crystal structure of TNAZ was obtained from the X-ray diffraction crystal data [1]. The space group of TNAZ is Pbca and a unit cell contains 8 molecules. Then, it was expanded along the a, b and c directions to construct a 6 × 3 × 1 TNAZ supercell. The 6 × 3 × 1 TNAZ supercell was constructed in Material Studio (MS) software and contained 144 molecules (2448 atoms). The parameters of the expanded supercell were a = 34.398 Å, b = 33.381 Å, c = 21.496 Å, α = 90.0 • , β = 90.0 • and γ = 90.0 • . A schematic diagram of the TNAZ molecule and supercell is shown in Figure 1. Before shock loading, the supercell was relaxed to obtain a stable structure under normal temperature and pressure. The structure was first relaxed by energy minimization, followed by thermalization at 298 K. A molecular dynamics (MD) simulation was first performed for 10 ps with the canonical ensemble (NVT, constant number of particles, volume and temperature) using the Berendsen thermostat. Then, a pressure relaxation at 0 GPa was performed for 5 ps using the isobaric-isothermal ensemble (NPT, constant number of particles, pressure and temperature) with the Nosé-Hoover thermostat and barostat. In the relaxation and formal calculations, the time step was 0.05 fs and the coupling coefficient was 100. We considered the use of different cut-off radii for different atom pairs proposed by reviewers. Here, different key types were set to different cut-off radii [14] (as shown in Table S1, and cut-off* is used in the following text). Finally, the crystal equilibrium structure at room temperature was obtained, and the density of the equilibrium structure was 1.84 g·cm −3 . Huang et al. verified the applicability of ReaxFF-lg to TNAZ, which is not discussed here [29]. After obtaining the crystal equilibrium structure, the TNAZ

Computational Methods
The simulations were performed with LAMMPS using ReaxFF-lg and MSST. The single-crystal structure of TNAZ was obtained from the X-ray diffraction crystal data [1]. The space group of TNAZ is Pbca and a unit cell contains 8 molecules. Then, it was expanded along the a, b and c directions to construct a 6 × 3 × 1 TNAZ supercell. The 6 × 3 × 1 TNAZ supercell was constructed in Material Studio (MS) software and contained 144 molecules (2448 atoms). The parameters of the expanded supercell were a = 34.398 Å, b = 33.381 Å, c = 21.496 Å, α = 90.0°, β = 90.0° and γ = 90.0°. A schematic diagram of the TNAZ molecule and supercell is shown in Figure 1. Before shock loading, the supercell was relaxed to obtain a stable structure under normal temperature and pressure. The structure was first relaxed by energy minimization, followed by thermalization at 298 K. A molecular dynamics (MD) simulation was first performed for 10 ps with the canonical ensemble (NVT, constant number of particles, volume and temperature) using the Berendsen thermostat. Then, a pressure relaxation at 0 GPa was performed for 5 ps using the isobaric-isothermal ensemble (NPT, constant number of particles, pressure and temperature) with the Nosé-Hoover thermostat and barostat. In the relaxation and formal calculations, the time step was 0.05 fs and the coupling coefficient was 100. We considered the use of different cutoff radii for different atom pairs proposed by reviewers. Here, different key types were set to different cut-off radii [14] (as shown in Table S1, and cut-off* is used in the following text). Finally, the crystal equilibrium structure at room temperature was obtained, and the density of the equilibrium structure was 1.84 g•cm −3 . Huang et al. verified the applicability of ReaxFF-lg to TNAZ, which is not discussed here [29]. After obtaining the crystal equilibrium structure, the TNAZ supercell system was subjected to shock in the

Shock Wave Velocity and Particle Velocity
The selected particle velocity is the particle velocity at the moment when the system yields, which is used to describe the state of the system after shock compression. During the shock process, the shock velocity has a linear relationship with the postwave particle velocity, i.e., the Hugoniot relationship. The formula is where C is the sound velocity of the dielectric material, and S is the adiabatic bulk modulus.
The relationship between the TNAZ shock velocity and particle velocity and the relationship between the pressure and specific volume were obtained through calculations. As shown in Figure 2, the calculated data were consistent with the experimental results,

Shock Wave Velocity and Particle Velocity
The selected particle velocity is the particle velocity at the moment when the system yields, which is used to describe the state of the system after shock compression. During the shock process, the shock velocity has a linear relationship with the postwave particle velocity, i.e., the Hugoniot relationship. The formula is where C is the sound velocity of the dielectric material, and S is the adiabatic bulk modulus. The relationship between the TNAZ shock velocity and particle velocity and the relationship between the pressure and specific volume were obtained through calculations. As shown in Figure 2, the calculated data were consistent with the experimental results,   Figure 3 shows the hydrostatic pressure curve of the system with time i shock velocities. At the initial moment, when the system is subjected to shock sion, the pressure quickly jumps. With a greater shock wave velocity, a larger in wave pressure is generated, and the pressurization time is shorter. Then, the re gins to occur, and the pressure further rapidly increases. The rate of the pressu and its peak value are related to the shock wave velocity and direction. When wave velocities are identical, the pressures in different loading directions hav evolution trend. A greater shock wave velocity corresponds to a greater press system.  Table S2 shows the initial pressure value, stable pressure value and corr time in different directions. The initial pressure in this paper refers to the first sure when the shock wave was applied to the supercell. When the shock ve km•s −1 , the initial shock pressures are 33.55, 33.31 and 33.60 GPa in the [100], [001] directions, respectively. At a shock velocity of 8 km•s −1 , the reaction is slo pressure of the TNAZ system does not significantly increase within 150 ps. Ac the slope of the curve in Figure 4, the chemical reaction is divided into the fas reaction and the slow chemical reaction stage. Therefore, no statistical analysi formed at 8 km•s −1 . When the shock velocity is greater than or equal to 9 km•s − sure of the TNAZ system obviously increases and reaches a stable value. In Tabl

Temperature and Pressure in Different Shock Directions
The total energy curves of the TNAZ system under different shock wave loading conditions are shown in Figure S1. Figure 3 shows the hydrostatic pressure curve of the system with time in different shock velocities. At the initial moment, when the system is subjected to shock compression, the pressure quickly jumps. With a greater shock wave velocity, a larger initial shock wave pressure is generated, and the pressurization time is shorter. Then, the reaction begins to occur, and the pressure further rapidly increases. The rate of the pressure increase and its peak value are related to the shock wave velocity and direction. When the shock wave velocities are identical, the pressures in different loading directions have a similar evolution trend. A greater shock wave velocity corresponds to a greater pressure of the system.   Figure 3 shows the hydrostatic pressure curve of the system with time in di shock velocities. At the initial moment, when the system is subjected to shock com sion, the pressure quickly jumps. With a greater shock wave velocity, a larger initial wave pressure is generated, and the pressurization time is shorter. Then, the reacti gins to occur, and the pressure further rapidly increases. The rate of the pressure in and its peak value are related to the shock wave velocity and direction. When the wave velocities are identical, the pressures in different loading directions have a s evolution trend. A greater shock wave velocity corresponds to a greater pressure system.  Table S2 shows the initial pressure value, stable pressure value and correspo time in different directions. The initial pressure in this paper refers to the first peak sure when the shock wave was applied to the supercell. When the shock veloci km•s −1 , the initial shock pressures are 33.55, 33.31 and 33.60 GPa in the [100], [01 [001] directions, respectively. At a shock velocity of 8 km•s −1 , the reaction is slow, a pressure of the TNAZ system does not significantly increase within 150 ps. Accord the slope of the curve in Figure 4, the chemical reaction is divided into the fast ch reaction and the slow chemical reaction stage. Therefore, no statistical analysis wa formed at 8 km•s −1 . When the shock velocity is greater than or equal to 9 km•s −1 , the sure of the TNAZ system obviously increases and reaches a stable value. In Table S2 the shock velocity is identical, as the reaction progresses, the rate of pressure incre  Table S2 shows the initial pressure value, stable pressure value and corresponding time in different directions. The initial pressure in this paper refers to the first peak pressure when the shock wave was applied to the supercell. When the shock velocity is 8 km·s −1 , the initial shock pressures are 33.55, 33.31 and 33.60 GPa in the [100], [010] and [001] directions, respectively. At a shock velocity of 8 km·s −1 , the reaction is slow, and the pressure of the TNAZ system does not significantly increase within 150 ps. According to the slope of the curve in Figure 4, the chemical reaction is divided into the fast chemical reaction and the slow chemical reaction stage. Therefore, no statistical analysis was performed at 8 km·s −1 . When the shock velocity is greater than or equal to 9 km·s −1 , the pressure of the TNAZ system obviously increases and reaches a stable value. In Table S2, when the shock velocity is identical, as the reaction progresses, the rate of pressure increase in the [100] direction is in pregnancy [1,3,4]. Obesity has been identified as a significant risk factor for matern diabetes. A meta-analysis of 20 studies reported that women who were overweight (2 fold), obese (3.6-fold) or severely obese (8.6-fold) had a higher risk of developing diabet compared to normal-weight pregnant women [5].

Temperature and Pressure in Different Shock Directions
Maternal diabetes is associated with pregnancy complications and increased rates adverse maternal and neonatal outcomes [6,7]. Short-term complications include macr somia, large for gestational age (LGA), respiratory distress syndrome (RDS), neonatal h poglycaemia, neonatal intensive care unit (NICU) admission, intrauterine growth r striction, congenital anomalies, preterm birth, pre-eclampsia, caesarean section (C/S) a preterm birth while in the long-term both mothers and their babies have an increased ri of metabolic disease [8][9][10]. Women with GDM have a ~7-fold increased risk of developi T2DM [11] and a ~4-fold increased risk of developing cardiovascular and coronary arte disease after pregnancy [12], while pregestational diabetes predisposes women to dev oping diabetes-related complications such as retinopathy and nephropathy or may acc erate the course of these complications if they already exist [4,7,13].
It is widely reported that all types of maternal diabetes are associated with pregnan complications; although, adverse outcomes are more common in women with preges tional diabetes [14][15][16][17][18]. As adverse pregnancy outcomes are closely related to poor gl caemic control and the first trimester being a critical period for organogenesis, it is spec lated that preconception hyperglycaemia and the longer time of exposure to hyperglyca mia in utero may contribute to the complications associated with pregestational diabet [19].
Despite the large body of evidence that associates pregestational diabetes with mo frequent adverse pregnancy outcomes than GDM [20][21][22][23][24][25], conflicting results have be reported [17,22,[26][27][28]. This review aims to summarise and synthesise studies that ha compared adverse pregnancy outcomes in pregnancies complicated by pregestational d abetes and GDM. Three databases, Pubmed, Scopus and EBSCOhost were searched identify eligible studies, which were summarised and synthesised using systematic r view methods. Commonly reported adverse pregnancy outcomes in literature [29] we selected for inclusion in this review. These include congenital anomalies, pre-eclamps  Figure 4 shows the curve of the supercell temperature with time. First, the supercell is subjected to shock compression under shock, and kinetic energy is converted into internal energy. Behind the shock wave front, the temperature of the medium sharply increases; then, rapid chemical reactions occur, which release a large amount of heat. Here, the stage where the temperature rapidly increases is called the rapid chemical reaction stage, and the stage where temperature slowly increases is called the slow reaction stage. When the rapid chemical reactions are completed, there are many clusters in the products, and the clusters continue to decompose and produce the final small-molecule products. There are relatively slow chemical reactions in this process, and the released heat is much less than that released in the rapid chemical reaction stage. When the shock velocity is 8 km·s −1 , due to the low compression intensity of the shock, the rapid chemical reactions of different directions begin at 130 ps. When the shock velocity is 9 km·s −1 , the rapid chemical reactions of different directions end in 50-70 ps and enter the slow chemical reaction stage. When the shock velocity is 10 and 11 km·s −1 , the rapid chemical reactions of different directions are completed in 0-20 ps, and the slow chemical reaction stage begins.
At shock velocities of 8, 9, 10 and 11 km·s −1 , the temperatures in the [100] direction are 1243, 3987, 4903 and 6106 K at 50 ps, respectively. When the shock velocity is larger, the temperature of the system rises faster, the system starts to react quickly, the reaction rate is larger, and the system enters into the slow chemical reaction stage earlier. From the temperature curves of 9 and 10 km·s −1 , it can be seen that in the rapid chemical reaction stage, the temperature rise rate of the system in the [100] direction is significantly higher than that of the [010] and [001] directions. It shows that the rate of rapid chemical reaction of TNAZ under the shock loading in the [100] direction is higher than that of the other two directions. When the shock wave velocity is 8 km·s −1 , due to the low shock pressure, the chemical reaction is slow, and within the calculation time of 200 ps, only a small increase in the temperature of the system occurs in the [100] direction. For the shock wave velocity of 11 km·s −1 , due to the large shock loading intensity, the explosive completes a rapid chemical reaction in a very short time, and the temperature changes of the system are not significantly different under different shock loading directions. In conclusion, it can be seen from Figures 3 and 4 that the TNAZ exhibits obvious anisotropy characteristics for shock wave loading with different directions, and the difference is the most obvious for the velocity of 9 km·s −1 .

Study on the Causes of TNAZ Anisotropy
TNAZ has strong anisotropic pressure and temperature responses in different directions, which may be related to the molecular arrangements in different directions. Here, this paper takes the supercell with a shock velocity of 9 km·s −1 as an example for the analysis. Figure 5 shows the changes in internal structure of the TNAZ supercell shocked in the [100] direction. For the convenience of observation, this paper uses the "Y" figure to represent the approximate shape of the TNAZ single molecule in the [100] direction. From the perspective of molecular arrangement, TNAZ molecules are arranged in a straight line in the [100] direction, and the direction of the molecules on the same line is identical. When the supercell is shocked in the [100] direction, a vertex of the "Y" shape is squeezed with the adjacent "Y" shape. The squeezing easily breaks the molecular structure and makes the supercell react. As shown in Table 1

Study on the Causes of TNAZ Anisotropy
TNAZ has strong anisotropic pressure and temperature responses in different directions, which may be related to the molecular arrangements in different directions. Here, this paper takes the supercell with a shock velocity of 9 km•s −1 as an example for the analysis. Figure 5 shows the changes in internal structure of the TNAZ supercell shocked in the [100] direction. For the convenience of observation, this paper uses the "Y" figure to represent the approximate shape of the TNAZ single molecule in the [100] direction. From the perspective of molecular arrangement, TNAZ molecules are arranged in a straight line in the [100] direction, and the direction of the molecules on the same line is identical. When the supercell is shocked in the [100] direction, a vertex of the "Y" shape is squeezed with the adjacent "Y" shape. The squeezing easily breaks the molecular structure and makes the supercell react. As shown in Table 1, when t = 0 ps, the initial number of TNAZ molecules is 144. When t = 10 ps, the number of TNAZ molecules is reduced to 94 in the [100] direction. Compared to the [010] and [001] directions, the [100] direction has more small molecules. There are five NO2 molecules, two HNO2 molecules, one O2 molecule and two NO3 molecules, respectively. Along the [100] direction, the arrangement of TNAZ molecules cannot be maintained, and the molecular structure has been significantly deformed. Thus, in this direction, the reaction rate of the supercell is the fastest.    Figure 6 shows changes in the internal structure of the TNAZ supercell shocked in the [010] direction. From the perspective of the molecular arrangement, the TNAZ molecules are also arranged in a straight line in the [010] direction. However, the directions of adjacent TNAZ molecules on the same line are opposite. When the supercell is shocked in the [010] direction, a vertex of the "Y" shape (Y1) and two vertices of the adjacent "Y" shape (Y2 and Y3) are near each other. In the [010] direction, there is a large space between adjacent TNAZ molecules to form a "buffer" for the shock in this direction. The supercell remains relative stable in this direction. As shown in Table 2, when t = 0 ps, the initial number of TNAZ molecules is 144. When t = 10 ps, the TNAZ molecules are shocked in the [010] direction, and the number of TNAZ molecules is reduced to 121. The TNAZ molecules maintain a linear arrangement in the [010] direction, and the molecular structure remains intact. At 10 ps, the number of unreacted TNAZ molecules in the system after being shocked in the [100] direction is less than that in the [100] direction. There are not many small molecules in the [010] direction. The number of NO 2 , HNO 2 , O 2 and NO 3 is two, zero, one and two, respectively. Thus, in this direction, the rate of reactions is lower than in the [100] direction.   Figure 6 shows changes in the internal structure of the TNAZ supercell shocked in the [010] direction. From the perspective of the molecular arrangement, the TNAZ molecules are also arranged in a straight line in the [010] direction. However, the directions of adjacent TNAZ molecules on the same line are opposite. When the supercell is shocked in the [010] direction, a vertex of the "Y" shape (Y1) and two vertices of the adjacent "Y" shape (Y2 and Y3) are near each other. In the [010] direction, there is a large space between adjacent TNAZ molecules to form a "buffer" for the shock in this direction. The supercell remains relative stable in this direction. As shown in Table 2, when t = 0 ps, the initial number of TNAZ molecules is 144. When t = 10 ps, the TNAZ molecules are shocked in the [010] direction, and the number of TNAZ molecules is reduced to 121. The TNAZ molecules maintain a linear arrangement in the [010] direction, and the molecular structure remains intact. At 10 ps, the number of unreacted TNAZ molecules in the system after being shocked in the [100] direction is less than that in the [100] direction. There are not many small molecules in the [010] direction. The number of NO2, HNO2, O2 and NO3 is two, zero, one and two, respectively. Thus, in this direction, the rate of reactions is lower than in the [100] direction.   Figure 7 shows the changes in internal structure of the TNAZ supercell shocked in the [001] direction. From the perspective of the molecular arrangement, the TNAZ molecules are arranged as a wavelike line in the [001] direction. When the supercell is shocked in the [001] direction, a vertex of the "Y" shape (Y1) is also near two vertices of the adjacent "Y" shape (Y2 and Y3). When adjacent TNAZ molecules are close, there is a large space between the molecules, which form a type of "buffer" for the shock in this direction. Thus, the molecular structure is not easily destroyed in the [001] direction, and the supercell is stable in this direction. As shown in Table 3, when t = 0 ps, the initial number of TNAZ molecules is also 144. When t = 10 ps, the number of TNAZ molecules is reduced to 132 in the [001] direction. The TNAZ molecules maintain the arrangement of a wavelike line in the [001] direction, and many molecules are intact. There are not many small molecules in the [010] direction. The number of NO 2 , HNO 2 , O 2 and NO 3 is four, zero, one and two, respectively. Hence, in this direction, the rate of reactions is also lower than that in the [100] direction.    Based on the steric hindrance model proposed by Dick et al. [31], during the impact sliding process, molecules avoid close contact between each other through rotational deformation and other methods to cushion the impact compression. In the [100] direction, the molecules are arranged in a good order, the front-row molecules directly collide with the back-row molecules, and the impact compression effect cannot be buffered. In the [010] and [001] directions, due to the staggered arrangement of molecules, along the compression direction, there are some slip spaces between the molecules, so the molecules move into the intermolecular space, rotate and deform to avoid molecular collisions. Therefore, the impact sensitivity of TNAZ in the [100] direction is higher than that in the other two directions.

Initial Reactions and the Chemical Bonds
The high-frequency reactions of the first 10 ps at 9 km·s −1 is shown in Table 4. Bimolecular polymerization first occurs in the shock. By the time of bimolecular polymerization in different directions, the "buffering" effect in different directions can be revealed. In the  To study the formation of clusters in the early stage of the reactions, the bonds (O-O, N-O and H-N) with obvious changes in the early stage of the reactions were analyzed, as shown in Figure 8. In Figure 8, the positive value plotted on the Y-axis represents the number of positive net bonds, and the negative value represents the number of negative net bonds, that is, the number of net broken bonds. It can be seen from Table 4 that the polymerization of the TNAZ system frequ occurs before 10 ps. Therefore, in the first 10 ps reactions, the TNAZ molecules in the shock direction are mainly combined into clusters by The occurrence of the chemical reaction must be accompanied by the breaki bonds. The shock directions affect the formation and cleavage of bonds at the initial In other words, the shock directions affect the initial reaction path of the TNAZ syst

Main Intermediate Products and Cluster Analysis
In different shock directions, the small-molecule products of TNAZ are H2O, N2 NO2, O2, NO3, HNO2, etc. Here, typical small-molecule products such as O2, NO HNO2 were selected for the analysis. When the shock velocities are identical, the m It can be seen from Table 4 that the polymerization of the TNAZ system frequently occurs before 10 ps. Therefore, in the first 10 ps reactions, the TNAZ molecules in the In the [010] and [001] shock directions, the trend of net bond formation is not obvious, and almost all of them show net bond breaking. This phenomenon may be related to the shock direction. In the [100] direction, the number of N-O bonds has an upward peak during the first 5-10 ps. This is due to the close distance between the N and O atoms in adjacent TNAZ molecules in the [100] direction. Shock wave compression makes the distance between the N and O atoms shorter than the cut-off radius, resulting in an increase in the number of N-O bonds.
The occurrence of the chemical reaction must be accompanied by the breaking of bonds. The shock directions affect the formation and cleavage of bonds at the initial stage. In other words, the shock directions affect the initial reaction path of the TNAZ system.

Main Intermediate Products and Cluster Analysis
In different shock directions, the small-molecule products of TNAZ are H 2 O, N 2 , CO 2 , NO 2 , O 2 , NO 3 , HNO 2 , etc. Here, typical small-molecule products such as O 2 , NO 3 and HNO 2 were selected for the analysis. When the shock velocities are identical, the mechanisms of reactions in different directions to form small molecules are not exactly identical. When the supercell is shocked in a certain direction, the distance between the atoms in the direction of the shock decreases. When the atoms or groups on some molecules are detached, the decrease in distance makes it easy for these atoms or groups to combine, which causes the differences between the reaction paths of the small molecules in subsequent reactions. Table 5      The evolutions of the number of clusters and molecular mass of maximum clusters in the TNAZ system with different shocks are shown in Figure 9a,b. For the convenience of statistics and analysis, "clusters" refers to molecules whose relative molecular weight is more than that of a single TNAZ molecule. First, the most rapid increase in number of clusters occurs in the [100] direction, but the most rapid decrease in number of clusters also occurs in this direction. The numbers of clusters in the [010] and [001] directions reach the maximum value at very similar times, but they are later than that in the [100] direction. When the number of clusters decreases, the molecular weight of the largest cluster rapidly increases, which indicates that the reduced clusters are aggregated to the largest cluster.

Conclusion
In this paper, the ReaxFF-lg force field and MSST method were used to study anisotropic reaction mechanism of the TNAZ shock compression process. The chem reactions in different shock directions were calculated. When the shock velocities identical, the increasing rates of pressure and temperature in the [100] direction greater than those in the [010] and [001] directions. There were significant differenc molecular arrangements between different directions, which affected the chemical tion rate. Due to the differences in arrangements, the shock sensitivity of the TNAZ su cells in different directions was ranked as [001] directions can form a buffer for the shock wave. This paper provides not only m information about the properties of explosives for other scholars studying TNAZ ex sives but also new ideas and methods for studying the anisotropic reaction princip energetic materials.
Supplementary Materials: The following supporting information can be downloade www.mdpi.com/xxx/s1. Figure S1: Total energy curve of the TNAZ system under different s wave loading conditions. It can be seen from the figure that the energy curve is continuously r indicating that the total energy of the computing system is not conserved. In the slow che reaction stage, the nonconservation of energy may be caused by the explicit integration algo of the LAMMPS calculation program, Table S1: List of bond order minimum values used to d mine molecules, Table S2: Initial shock pressure and stable pressure.

Conclusions
In this paper, the ReaxFF-lg force field and MSST method were used to study the anisotropic reaction mechanism of the TNAZ shock compression process. The chemical reactions in different shock directions were calculated. When the shock velocities were identical, the increasing rates of pressure and temperature in the [100] direction were greater than those in the [010] and [001] directions. There were significant differences in molecular arrangements between different directions, which affected the chemical reaction rate. Due to the differences in arrangements, the shock sensitivity of the TNAZ supercells in different directions was ranked as [100] > [010] > [001]. The shock directions affected the formation and cleavage of the N-O bonds in the initial stage. In other words, the shock directions affected the initial reaction path of the TNAZ system. The shock directions also affected the reaction paths and appearance time of small molecules such as O 2 , NO 3 and HNO 2 . According to the reactions of the first 10 ps and the evolution of the clusters, the polymerization reactions in the [010] and [001] directions appeared later than that in the [100] direction, and the cluster growth in the [010] and [001] directions was slower than that in the [100] direction. Therefore, the TNAZ molecules in the [010] and [001] directions can form a buffer for the shock wave. This paper provides not only more information about the properties of explosives for other scholars studying TNAZ explosives but also new ideas and methods for studying the anisotropic reaction principle of energetic materials.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/molecules27185773/s1. Figure S1: Total energy curve of the TNAZ system under different shock wave loading conditions. It can be seen from the figure that the energy curve is continuously rising, indicating that the total energy of the computing system is not conserved. In the slow chemical reaction stage, the nonconservation of energy may be caused by the explicit integration algorithm of the LAMMPS calculation program, Table S1: List of bond order minimum values used to determine molecules, Table S2: Initial shock pressure and stable pressure.