Viscoelastic Properties of Inert Solid Rocket Propellants Exposed to a Shock Wave

Inert solid rocket propellant samples were subjected to dynamic inflation experiments, to characterize the viscoelastic response at high strain rates. An oxyacetylene-driven shock tube created the shock wave, which was used to dynamically pressurize the surface of the samples during the inflation experiments. Two high-speed cameras captured the deforming samples, which were speckled to measure the full-field surface displacements using the digital image correlation (DIC) algorithm. An inverse finite element analysis (iFEA) was used to calibrate parameters of a generalized Maxwell model (i. e. Prony series), which was used to characterize the propellants’ viscoelastic response to shock wave exposure. The viscoelastic parameters calibrated using a Prony series with one Maxwell branch provided a better fit with the out-of-plane displacement data from DIC. At least 50 % of the energy dissipated, within the inert solid rocket propellant, occurred within 5 ms following shock wave exposure. The softening phenomenon, due to debonding of the particles embedded in the inert solid rocket propellant, occurred since there was a decrease in instantaneous elastic modulus with increased strain rate. The results of this study will add to the limited knowledge of the linear viscoelastic behavior of inert HTPB propellant at high strain rates and may improve the predictive capabilities of health-monitoring sensors that assess the solid rocket propellant's structural integrity.


Introduction
Solid rocket propellants are energetic materials used as fuel in military applications to generate thrust for tactical or strategic rockets and missiles or generate gas in commercial applications (e. g. inflation of automobile airbags) [1,2,3,4]. The solid rocket propellants are comprised of a high-volume fraction (60 %-75 %) of microscopic solid oxidizer particles that are bonded to an elastomeric matrix such as hydroxyl-terminated polybutadiene (HTPB), but contain other chemical ingredients such as fuel, plasticizer, and stabilizer [2,5,6,7]. Since the matrix is elastomeric, HTPB is classified as a viscoelastic material; meaning that the stress response is strain-rate and temperature-dependent [3,8,9]. The ballistic performance of HTPB propellant is improved by increasing the volume fraction of the particulates, with a decrease in structural integrity as the trade-off [5,10]. The decreased structural integrity is attributed to the creation of an irregular network structure in the elastomeric matrix, due to increasing the volume fraction of the solid oxidizer particles [1]. Structural integrity is further compromised when the HTPB propellant is mechanically loaded multiaxially, during storage and/or the ignition and firing process of the solid rocket propellant [5,8,11]. Mechanical loading causes interfacial debonding (i. e. dewetting) of the particles that are bonded to the matrix, which leads to damage of the solid rocket propellant [12,13]. The damage manifests as cracks in the HTPB propellant, which affects the ballistic performance of the solid rocket propellant and the structural in-tegrity of the solid rocket motor [14]. In the case of solid rocket motors that experience bullet impact damage or fragment impact damage, this debonding takes place under high strain rate conditions, but conventional test methods do not allow for equivalent strain rates in the laboratory. As a result, designing health-monitoring sensors that can quantify the structural integrity of the solid rocket propellants would be invaluable [15]. However, the predictive capabilities of the health-monitoring sensors require additional research into the viscoelastic behavior of the solid rocket propellant at high strain rates [9,16] using non-conventional test methods such as dynamic inflation experiments using a shock tube as done by Bentil et al. [17].
Stress level, temperature, strain rate, and repeated loading are all factors that have been reported to influence dewetting of the viscoelastic HTPB propellant [6,9,16,18]. These factors were determined by first constructing the stress-strain relations for HTPB propellant, using commercial mechanical servohydraulic machines. HTPB propellant samples were clamped into the mechanical servohydraulic machines before uniaxial tensile experiments. The uniaxial tensile experiments, performed at strain rates less than 10 2 s À 1 , showed that both aging and increased cross-linking densities caused a decrease in the mechanical properties [19,20]. HTPB propellant samples were thermally aged in the study by Du et al. [19] and chemically aged in the study by Gligorijević et al. [20], to induce damage due to dewetting.
Stress-strain relations for HTPB propellant samples subjected to very high strain rates ( > 10 3 s À 1 ) have been obtained using the split-Hopkinson pressure bar (i. e. Kolsky bar) [21], even though these strain rates exceed the ignition rates (i. e. 1 m/s-6 m/s) of the solid rocket propellant [16]. In the Kolsky bar experiments, the stress-strain response of elastomers (e. g. HTPB propellant) under compressive loading can be measured as a function of strain rates that range from 10 2 s À 1 -10 4 s À 1 [22]. However, the literature is lacking high strain rate experiments of HTPB propellant under tensile loading, at rates that exceed the ignition rates. Thus, hindering the ability to completely characterize the mechanical behavior of HTPB propellant [16], which would be invaluable for increased confidence in the reliability of predicting the solid propellant's structural integrity [9]. Unfortunately, a commercially available mechanical testing machine that can construct the tensile stress-strain curves of soft materials, at high rates, does not exist [23,24]. Thus, making the complete characterization of HTPB propellant's dynamic mechanical behavior challenging. As a solution to the lack of a commercially available mechanical testing machine, this study combines an oxyacetylene-driven shock tube apparatus with the digital image correlation (DIC) method and inverse finite element analysis (iFEA) to characterize the viscoelastic behavior of inert HTPB propellant at high strain rates on the order of 10 4 s À 1 . The inert HTPB propellant samples were subjected to bending, using dynamic inflation, which produces both tensile and compressive loads for improved characterization of the high strain rate response.
Digital image correlation is a non-contact optical method that has been widely applied in many areas of science and engineering for over 30 years [25,26,27]. Contactless measurements of soft deformable materials are important since minor disturbances on the material surface, due to strain gauge attachment, can influence the mechanical behavior [17,28,29,30]. DIC is capable of quantifying both small and large deformation of soft materials at strain rates that range from quasistatic (10 À 6 s À 1 -10 0 s À 1 ) to shock wave (10 6 s À 1 -10 8 s À 1 ) [17,[31][32][33][34][35][36][37][38]. As a result, the DIC technique was applied in this study to measure the three-dimensional (3D) full-field surface displacement of inert solid rocket propellants at high strain rates. The results from DIC were used to validate the iFEA.
Inverse finite element analysis has been applied in the literature to calibrate constitutive model parameters describing the mechanical properties of viscoelastic materials [39,40]. During the iFEA, a comparison is made between the output from an experiment and simulated output following finite element simulations of the experiment. A match between the output implies that the mechanical properties of the material have been found. In this study, the iFEA was applied to determine the linear viscoelastic properties of inert solid rocket propellant samples exposed to a shock wave during a dynamic inflation experiment. Specifically, a comparison of calibrated constants using a Prony series, with either one, two, or three Maxwell branches, was considered for the linear viscoelastic model of the inert solid rocket propellant. A Prony series is generally used in the literature to model the linear viscoelastic behavior of HTPB propellants [9].

Sample Preparation
Inert solid rocket propellants, which were a simulant of the polymer matrix binder hydroxyl-terminated polybutadiene (HTPB) formulation, were fabricated for this study. Ammonium sulfate was used as the inert ingredient and served as the substitute for the oxidizer. The percent filler by volume for the inert propellant was 80 %, which is comparable with a real solid rocket propellant. A bimodal size distribution, similar to the energetic formulation, was used. The material was mixed in a one-gallon small-clearance batch mixer and was then vacuum cast into a large block. The block was then cured for approximately seven days at 60°C. Eight (8) solid rocket propellant samples were cut from the inert HTPB propellant block, using a horizontal milling machine and a band saw. Initially, 12.7-mm (0.5-in.) thick slabs were roughly cut from the block using the band saw. The length and width of the slab were 188-mm (7.4-in.) and 176-mm (7-in.), respectively. Opposite faces were then milled flat and parallel, using a fly cutter revolving at 1000 rpm. Approximately 1-mm (0.04-in.) of material was removed each time, until thin rectangles with a thickness of 5-mm (0.2-in.) was obtained. The rectangular HTPB propellant samples were shipped from Edwards Air Force Base to Iowa State University, for dynamic inflation experiments using the oxyacetylene-driven shock tube.
Prior to shock tube experiments, the rectangular HTPB propellant samples were cut into a circular plate (i. e. discshaped) geometry having a diameter of 152.4-mm (6-in.) with a thickness of 5-mm (0. 2-in.). The diameter of the HTPB propellant sample matched the outer diameter of the annuli used to clamp and attach the sample to the end of the oxyacetylene-driven shock tube. The inner diameter of the annuli plate was 7.44 cm (2.93-in.). The thickness of the HTPB propellant sample was chosen to facilitate characterization of the viscoelastic properties, due to shock wave exposure, using concepts from plate bending theory.

Shock Tube Experiments
An oxyacetylene-driven shock tube, with a 7.62 cm (3-in.) inner diameter, was used to conduct the dynamic inflation experiments (Figure 1). A 25.4 μm (1-mil) thick Mylar diaphragm was used to separate the shock tube's 30.48 cm (1ft) driver and 457.2 cm (15-ft) driven section. Rupturing the Mylar diaphragm, by igniting the oxyacetylene in the driver section, created a propagating shock wave through the driven section. The shock tube's driven section was instrumented with three (3) high frequency integrated circuit piezoelectric (ICP) pressure transducers (PCB Piezotronics, Model 102B15), to measure the speed of the shock wave and capture its pressure-time characteristics as a function of oxyacetylene volume in the driver section. A range of shock wave pressures was considered by varying the volume of oxyacetylene in the driver section, to study the influence on the material response.
Pressure transducers were placed 152.4 cm (5-ft) apart from each other. The input signals from the three pressure transducers were conditioned using a signal conditioner (PCB Piezotronics, Model 482C16), before being measured by an oscilloscope (Tektronix, MDO3024). Following diaphragm rupture, the oscilloscope was triggered to begin recording the pressure transducer signals when the shock wave reached pressure transducer 3. Initiation of pressure recording corresponded to when the time in this study was set to zero. Recorded pressure transducer data continued for either 70 ms, 140 ms, or 180 ms at a sampling rate of 100 kS/s, 50 kS/s, or 50 kS/s, respectively. Decreasing the sampling rate did not affect the resolution of the pressure waveform. The range of pressure-time signal duration was increased incrementally to facilitate investigation into the long-term deformation behavior of inert HTPB propellant.
Each HTPB propellant sample was speckled using white spray paint with a flat (i. e. matte) finish, to create the inherent speckle pattern needed for DIC analysis. After speckling, the HTPB propellant sample was clamped to the end of the shock tube ( Figure 2). Two synchronized high-speed digital cameras (Photron, FASTCAM SA-Z), with 105-mm macro lens (Nikon, AF-S VR Micro-NIKKOR 105-mm f/2.8G IF-ED) attached, were calibrated for stereo vision using the principle of stereo-triangulation [17]. Calibrating the cameras for stereo vision was needed to perform three-dimensional DIC analysis (3D-DIC) using the commercial software VIC-3D [41], to quantify the 3D full-field surface deformation (i. e. displacement and strain) of HTPB propellant. During the calibration, images were acquired using a standard rigid calibration target that had 14 dots × 10 dots, with a dot spacing of 14-mm. This calibration target was placed in front of the stereo camera pair and moved in three-dimensional space. Specifically, the motion of the calibration target included in-plane translation and/or rotation, along with out-of-plane translation and/or rotation. The size of the calibration target was selected such that the field of view was filled. The acquired images of the calibration target in motion were then analyzed, using the software VIC-3D [41], to complete the camera calibration. A calibration score of 0.036 � 0.002 was obtained each time, which was below the maximum value of 0.05 specified by the VIC-3D [41] software as a default to minimize measurement errors associated with 3D-DIC analysis. Following calibration, images of the deforming HTPB propellant sample due to shock wave exposure were acquired at a frame rate of 100,000 fps and a resolution of 640 × 280 for a period of Figure 1. (a) Computer-aided design (CAD) of the oxyacetylene shock tube, with a driver (red) and driven (blue) section that is separated by a Mylar diaphragm. (b) The propagating shock wave will dynamically inflate the clamped inert HTPB propellant sample and two high-speed cameras will record the deformation for 3D-DIC analysis. During the 3D-DIC analysis, VIC-3D's default subset size and step size of 29 pixels and 7 were applied, respectively. A preliminary study found that the lowest projection error of less than 10 % was generated if using VIC-3D's default settings. Subsets of speckles were analyzed in VIC-3D, to obtain measurements of the out-of-plane displacement W.
The measurements of W were used to validate the optimized viscoelastic properties for the inert HTPB propellant samples exposed to a shock wave. Since the inert HTPB propellant samples were dynamically inflated using a shock wave, the samples were therefore subjected to bending loading. Thus, during inflation, the speckled surface (Figure 2b) of the inert HTPB propellant sample was in tension while the unspeckled surface (i. e. opposite face of the discshaped sample) was under compression.

Inverse Finite Element Analysis (iFEA) to Extract Mechanical Properties
The iFEA was used to calibrate the linear viscoelastic properties of the disc-shaped HTPB propellant samples. Figure 3 shows the methodology for the iFEA, which required coupling the commercially available software MATLAB [42] and ABAQUS Explicit [43]. MATLAB's genetic algorithm function (i. e. ga) was used to perform the optimization of the constitutive model parameters. The ga function uses one of the following four criteria, to determine when to stop the solver: (i) the maximum number of generations is reached (default is 100 times the number of parameters/variables), (ii) if there is no change in the best fitness value for some time given in seconds (stall time limit), (iii) if there is no change in the best fitness value for some number of generations (maximum stall generations), or (iv) if the maximum time limit (in seconds) specified by the user is reached. The best fitness value was defined using the sum of the squared estimate of error (SSE), which is the sum of the squares of the residual. The residual is the difference between the out-ofplane displacement W at the apex obtained using iFEA and 3D-DIC (Equation 1).
where: n is a sampled point in time, which represents each frame from the acquired images of the deforming HTPB propellant sample. N p is the total number of sampled points in time. W iFEA t n ð Þ and W 3DÀ DIC t n ð Þ represent the out-of-plane displacement of the HTPB propellant sample's apex at the time t n , which is predicted by the iFEA and measured by 3D-DIC, respectively. W iFEA , at the apex of the sample in the iFEA, is dependent upon the estimated parameters of a Prony series describing the linear viscoelastic response of the HTPB propellant material. Figure 4 shows the boundary and loading conditions, along with the mesh applied to the disc-shaped HTPB propellant samples during the iFEA. The diameter of the samples (in the iFEA) was equivalent to the inner diameter of the annuli, which was used to clamp the HTPB propellant sample during the dynamic inflation experiments. The HTPB propellant sample was modeled as a deformable solid in 3D modeling space, using solid (or continuum) second-order tetrahedron elements for the mesh. Since both material nonlinearity (i. e. large deformations) and geometric nonlinearity (large displacements) are expected during the iFEA, the option NLGEOM was activated in ABAQUS to avoid convergence issues. A mesh convergence study was conducted to select the appropriate mesh size, which consisted of 3,724 elements. A generalized Maxwell model ( Figure 5), with either one (1), two (2), or three (3) Maxwell branches in parallel, was used to determine which linear viscoelastic constitutive model could best capture the behavior of the inert HTPB propellant exposed to a shock wave. The Maxwell branch consists of a spring and dashpot in series. A Prony series, in the time t domain, is used to mathematically represent the generalized Maxwell model using the relaxation modulus E(t).
Equation 2 defines the relaxation modulus E t ð Þ using the Prony series [44]. where: is the instantaneous elastic modulus and E 1 is the long-term (or equilibrium) elastic modulus. The relaxation time t i is obtained through the following relation: where h i is the viscosity of the dashpot in Figure 5. N represents the number of Maxwell branches in the Prony series, and t is time.   Pyrotech. 2021, 46, 1-16 www.pep.wiley-vch.de E t ð Þ contains parameters (i. e. E and t) that need to be optimized, to characterize the linear viscoelastic behavior of HTPB propellant using the generalized Maxwell models in Figure 5. However, implementation of E t ð Þ in the iFEA requires expression of the Prony series expansion (Equation 2) using both the shear relaxation modulus (Equation 3) and bulk relaxation modulus (Equation 4) [45]. As a result, the linear viscoelastic response of HTPB propellant was characterized using the following parameters: E 0 , n 0 , g, k, and t. where: Þ is the instantaneous shear modulus, g i is the shear relaxation modulus ratio, and t i is the relaxation time. n 0 is the instantaneous Poisson's ratio. N represents the number of Maxwell branches in the Prony series, and t is time. where: Þ is the instantaneous bulk modulus, k i is the bulk modulus ratio, and t i is the relaxation time. n 0 is the instantaneous Poisson's ratio. N represents the number of Maxwell branches in the Prony series, and t is time.
An isotropic response and the Boltzmann superposition principle were assumed for the linear viscoelastic HTPB propellant material, to obtain the constitutive relation describing the state of 3D stress and strain. This constitutive relation is defined using the decomposition of the 3D stress s ¼ and strain e ¼ tensors into their volumetric (i. e. hydrostatic) and deviatoric parts. Equation 5 is the constitutive relation for the volumetric state of stress, while Equation 6 describes the deviatoric state of stress for HTPB propellant. In applying the equations for stress, the assumption is that the stresses at one point in the HTPB propellant material depend only on the strains at that same point.
where: s vol and e vol are the volumetric stress and strain, respectively. t 0 is the time at which the deformation is imposed and K t ð Þ is the bulk relaxation modulus defined in Equation 4.
where: s ¼ dev and e ¼ dev are the deviatoric stress and strain tensors, respectively. G t ð Þ is the shear relaxation modulus defined in Equation 3.

Experiments: Dynamic Inflation
The oxyacetylene-driven shock tube was used to perform the dynamic inflation experiments by first sandwiching each inert HTPB propellant sample in between two annuli plates, before clamping to the end of the shock tube. Table 1 provides the density for each inert HTPB propellant sample, the maximum reflected shock wave amplitude that was applied on each sample's surface, along with the corresponding shock wave speed, strain rate, and projection error from 3D-DIC analysis. Each sample's density was calculated using the mass of the sample and the equation for the volume of a disc (i. e. cylinder). Reflected shock wave pressures, with a maximum amplitude below 313 kPa did not cause the HTPB propellant sample to pop out of the annuli plates. To study the influence of variable pressures on the material response, one of the samples (IP_3) was exposed to a lower reflected pressure. This lower pressure resulted in a wave speed that was subsonic. The supplementary material contains images and videos of HTPB propellant samples that popped out of the annulus, due to a reflected shock wave pressure above 313 kPa. The shock wave speed was calculated using the distance between pressure transducers 2 and 3 and the time interval for the shock wave to travel between these two transducers. The strain rates reported in Table 1 and elsewhere in this manuscript were calculated by dividing the shock wave speed by the thickness of the propellant. This strain rate can be considered as a normalized loading rate applied on the inert HTPB propellant material. From the 3D-DIC analysis, the inplane strain rates for the inert HTPB propellant material were calculated using the in-plane displacements and the applied frame rate. The in-plane strain rate ranges for the xdirection ( e xx dt ), y-direction ( e yy dt ), and xy-direction ( e xy dt ) were as follows: e xx dt ¼ 215=s À 1098=s, e yy dt ¼ 296=s À 1151=s, and e xy dt ¼ 141=s À 496=s. The xy-direction is the shear case. Figure 6 compares the time history for out-of-plane displacement, obtained from the eight HTPB propellant samples following 3D-DIC analysis. The maximum out-of-plane displacement for IP_2 was highest since the reflected shock wave pressure was the highest. IP_3 had the lowest maximum out-of-plane displacement since the reflected shock wave pressure was the lowest. Additionally, the density of IP_3 was lower than the other samples, which influenced the waveform for the time history of the out-of-plane displacement. Figure 7 compares the results of DIC and pressure transducer data, for the eight HTPB propellant samples. During the 3D-DIC experiments, the duration of the images acquired was longer than the duration of the pressure transducer data captured by the oscilloscope. As a result, the duration of W in Figure 7 was decreased to match the maximum duration of the pressure waveform, to facilitate comparison. Additionally, the pressure waveforms in Figure 7 were smoothed using a moving average to minimize noise.

iFEA: Dynamic Inflation
Three representative samples were considered for the iFEA to facilitate a comparison of the calibrated constants as a function of reflected shock wave pressures ( Table 1) that were classified in this study as low (IP_3), high (IP_2), and between low and high (IP_1). IP 4-IP_6 were all in the category between low and high reflected shock wave pressures and were not considered to avoid redundancy. The search space used during the iFEA for the instantaneous elastic modulus E 0 and relaxation time t i was 1 × 10 6 Pa-5 × 10 7 Pa and 0 s-0.01 s, respectively. A search space of 0-0.99 was used for the shear relaxation modulus ratio g i and bulk modulus ratio k i . The criteria used by the ga function, to stop the solver once the Prony series parameters were optimized, was when the average relative change in the penalty fitness value SSE was less than 1 × 10 À 6 over 50 generations (i. e. the default setting for the maximum stall generations). Table 1. Inert HTPB propellant material and shock wave properties, along with the projection error from 3D-DIC analysis. "IP" means "inert HTPB propellant". *Density and reflected pressure of IP_3 was lower than the other samples and the wave speed was subsonic. Figure 6. Out-of-plane displacement W at the apex from DIC, for the eight inert HTPB propellant samples. "IP" means "inert HTPB propellant".

Viscoelastic Properties of Inert Solid Rocket Propellants Exposed to a Shock Wave
These are not the final page numbers! �� Propellants Explos. Pyrotech. 2021, 46, 1-16 www.pep.wiley-vch.de Figure 8 compares the results of DIC and iFEA, for the IP_1 sample, using the optimized parameters from the linear viscoelastic constitutive models. The Prony series with a single Maxwell branch was better suited at capturing the linear viscoelastic behavior of HTPB propellant exposed to a shock wave. This is because the SSE for a Prony series with one Maxwell branch was the lowest when compared with the SSE for two and three Maxwell branches. The percent difference for SSE between one and two Maxwell branches, one and three Maxwell branches, and two and three Maxwell branches were 6.35 × 10 À 8 %, 4.04 × 10 À 9 %, and 5.94 × 10 À 8 %, respectively. Furthermore, the first three peaks of the out-of-plane displacement time history using one Maxwell branch (red curve) occurred around the same time as the initial three peaks from the DIC experiments (black curve) in Figure 8. Matching the initial three peaks implied that a Prony series with a single Maxwell branch was the optimal linear viscoelastic model. An additional advantage of using a Prony series with one Maxwell branch was that it reduced the computational time for the iFEA to optimize the parameters of the Prony series. Table 2 shows the parameters that were optimized using a Prony series having a single Maxwell branch in ABA-QUS Explicit, for the first three HTPB propellant samples in Table 1. The optimized parameters for a Prony series with two and three Maxwell branches, in ABAQUS Explicit, are shown in Table 3 for the sample labeled as IP_1 in Table 1. The instantaneous Poisson's ratio ν 0 was not optimized but was set to a value of 0.49 since the HTPB propellant sample was assumed to be incompressible.
The iFEA results showed that a Prony series with a single Maxwell branch yielded an instantaneous elastic modulus that was 77 % and 83 % different when comparing IP_3 with IP_1 and IP_3 and IP_2, respectively. One reason for the difference when comparing IP_3 with IP_1 and IP_2 was due to the lower reflected shock wave pressure ( Table 1) that was applied on the surface of the inert solid propellant sample. A lower applied reflected shock wave pressure would yield a lower out-of-plane displacement, which will manifest as a higher instantaneous elastic modulus (i. e. stiffer material). The percent difference for the maximum reflected shock wave pressure between IP_3 and IP_1 and IP_ 3 and IP_2 was 90.4 % and 99.3 %, respectively. The instantaneous elastic modulus for IP_1 and IP_2 was 7 % different from each other, with 11.5 % as the percent difference for the maximum reflected shock wave pressure. Figure 9 compares the results of DIC and iFEA, for the IP_2 and IP_3 samples, using each sample's pertinent opti- Figure 8. Out-of-plane displacement W at the apex from DIC (black) and iFEA, for IP_1. The iFEA is considered a Prony series with one Maxwell branch (red), two Maxwell branches (gold), and three Maxwell branches (blue). "IP" means "inert HTPB propellant". Table 2. Optimized parameters for a Prony series with one Maxwell branch, using the out-of-plane displacement and reflected pressure data for the samples labeled as IP_1, IP_2, and IP_3. Calculated parameters (i. e. G 0 and K 0 ), using the optimized E 0 and constant Poisson's ratio n 0 = 0.49, are also included. "IP" means "inert HTPB propellant". Table 3. Optimized parameters for a Prony series with two and three Maxwell branches, using the out-of-plane displacement and reflected pressure data for the HTPB Propellant 1 sample. Calculated parameters (i. e. G 0 and K 0 ), using the optimized E 0 and constant Poisson's ratio n 0 = 0.49, are also included. Additionally, the SSE is provided and rounded to two decimal places. "IP" means "inert HTPB propellant".

Viscoelastic Properties of Inert Solid Rocket Propellants Exposed to a Shock Wave
These are not the final page numbers! �� Propellants Explos. Pyrotech. 2021, 46, 1-16 www.pep.wiley-vch.de mized parameters for the Prony series with one Maxwell branch. Ideally, the optimized parameters in a constitutive model should be able to describe the mechanical behavior of all samples. Table 4 assesses the ability of the optimized parameters from three samples (i. e. IP_1, IP_2, and IP_3) at fitting the out-of-plane displacement data from 3D-DIC. SSE was used as a metric for goodness-of-fit. The expectation is that SSE will increase when using optimized parameters from one sample to predict the out-of-plane displacement time history of a different sample.
A small tolerance of 1 × 10 À 6 was used as the stopping criteria for the optimization, during the iFEA. Increasing this tolerance will provide a range of values for the Prony series parameters that may be used to yield a good fit when compared with the 3D-DIC experimental data. Table 5 highlights the range of parameters for IP_1, IP_2, and IP_3 that will yield the smallest SSE value in Table 2 when rounded down to the nearest integer. Specifically, the range of values (Table 5) for IP_1, IP_2, and IP_3 will have a calculated SSE of 1277 (mm 2 ), 1308 (mm 2 ), and 842 (mm 2 ) when using the Prony Series with one Maxwell branch, respectively. Figure 10 compares the contour plots from DIC and iFEA, for IP_1, using the optimized parameters for the Prony series with one Maxwell branch. Videos of the dynamic inflation for IP_1, IP_2, and IP_3 can be found in the supplementary material.
A logarithmic strain (i. e. true strain) definition was applied due to the large (i. e. finite) deformations of the samples, while the stress was defined using Cauchy stress (i. e. Figure 9. Out-of-plane displacement W at the apex from DIC (black) and iFEA (red), for (a) IP_2 and (b) IP_3. The iFEA is considered a Prony series with one Maxwell branch. "IP" means "inert HTPB propellant". Table 4. Sum of the squared estimate of error (SSE) for the optimized parameters in the Prony series with one Maxwell branch.* The SSE is rounded to two decimal places. "IP" means "inert HTPB propellant".
* SSE was obtained using the optimized parameters for sample IP_ X and the data (i. e. out-of-plane displacement from 3D-DIC and reflected pressure from the pressure transducer 3) for sample IP_Y. For example, when IP_X = IP_1 and IP_Y = IP_2 then SSE = 1510.68 mm 2 . This implied that the SSE was obtained by using the optimized parameters for the Prony Series with one Maxwell branch from IP_1 to fit the apex's out-of-plane displacement data from IP_2. † The SSE will be the lowest when IP_X = IP_Y, since the optimized parameters were obtained using the out-of-plane displacement and reflected pressure for the same sample. true stress). The logarithmic strain components were then used to calculate the maximum principal strain, which was used as a metric to compare the linear viscoelastic behavior of the HTPB propellant samples. Figure 11 shows the maximum principal true stress versus maximum principal true strain using the optimized parameters of the Prony series with one Maxwell branch, for IP_1, IP_2, and IP_3. Three hysteresis cycles were distinctly captured during the 25 ms time period, as illustrated in Figure 11. The energy dissipated during each cycle of inflation (loading) and de-flation (unloading) decreased, as evidenced by the decrease in the area within the hysteresis loop. The percent decrease in energy dissipated, from the first hysteresis cycle to the second hysteresis cycle, was 50 %, 60 %, and 63 % for IP_1, IP_2, and IP_3, respectively. This percent decrease, in energy dissipated from the first and second hysteresis cycles, occurred within the initial 5 ms following shock wave exposure. Additionally, the 5-mm thick HTPB propellant samples could sustain strains as high as 26 % before macro- Figure 10. Out-of-plane displacement W contour plots at the apex (denoted by the star) for IP_1, using Prony series with one Maxwell branch. 3D contour plot from (a) DIC and (b) iFEA. (c) Two-dimensional (2D) contour plot from DIC, with a color bar scaled using the W ranges obtained from the 3D-DIC analysis. The out-of-plane displacement color bar ranges in (a) and (b) have been made equivalent, to facilitate comparison of DIC and iFEA. Since the maximum W was higher in the iFEA, the red color in (b) covers a larger area than in (a). Figure 11. Maximum principal true stress versus maximum principal true strain at the apex for (a) IP_1, (b) IP_2, and (c) IP_3, using the optimized parameters for a Prony series with one Maxwell branch. Three hysteresis cycles were distinctly captured during the 25 ms duration, which is labeled in the legend. "IP" means "inert HTPB propellant". scopic failure due to the formation of cracks from the excessive deformation.

Discussion
HTPB propellant samples were exposed to shock waves during a dynamic inflation experiment, to characterize the mechanical properties at high strain rates. During the experiment, 3D-DIC was applied to capture the surface deformation of the HTPB propellant samples. An iFEA was applied to predict the material properties of HTPB propellant, by considering a linear viscoelastic constitutive model described using the generalized Maxwell model. Since the composition of the inert propellant was comparable with a real solid rocket propellant, the use of a linear viscoelastic constitutive model was appropriate. For instance, Shekhar et al. [46] described the linear viscoelastic behavior of solid rocket propellants using a Maxwell model comprised of a spring and dashpot in series. The solid rocket propellants were subjected to uniaxial tension at low strain rates (i. e. less than 10 À 1 s À 1 ), in the study by Shekhar et al. [46]. However, the Maxwell model could not capture the anelastic (i. e. time-dependent) recovery of the solid rocket propellant. To overcome this limitation of the Maxwell model, the generalized Maxwell model using the Prony series in the time domain was considered. A Prony series, with one Maxwell branch, provided a better fit with the 3D-DIC data when compared with a model containing two or three Maxwell branches since it yielded the smallest SSE. Additionally, a Prony series with one Maxwell branch required optimization of only four parameters (i. e. E 0 , g 1 , k 1 , and t 1 ), which reduced the computational time for the iFEA and made this model a good choice for HTPB propellant. Ideally, a good constitutive model is one that uses a relatively small number of fitted parameters to capture the mechanical response. Furthermore, the Prony series with the single Maxwell branch was also considered as good since this model was able to capture the time when the out-of-plane displacement peaks occurred when compared with the Prony series with two or three Maxwell branches.
In theory, introducing additional Maxwell branches into the generalized Maxwell model should improve the fit [47,48,49] of the iFEA results when compared with the 3D-DIC experimental data. As the fit improves, the calculated value for SSE should approach zero with a perfect match having an SSE ¼ 0. However, improved fit with increased Maxwell branches was not the case in this study. Rather, increasing the number of Maxwell branches yielded negligible changes in the calculated SSE. The main reason why the SSE value did not approach zero, with increased Maxwell branches, was attributed to the lack of a damage model in the iFEA and not the spatial resolution from DIC. The 3 pixels/mm spatial resolution from this study was good enough to have confidence in the DIC results since the 1megapixel image sensors for the high-speed digital cam-eras required larger features and spacing for the given fieldof-view and macro lens used [50]. Inclusion of a damage model (e. g. cohesive zone model [51,52,53]) will be considered in the future to facilitate the prediction of the debonding behavior of the particles embedded in the HTPB matrix caused by large deformations. The cohesive zone models will account for the damage in the bulk response that is caused by the local (microstructural) fracture due to each particle separating from the surrounding matrix [13], which would improve the fit of the iFEA and DIC and result in a lower SSE.
Yun et al. [54] examined the dewetting mechanism at the particle-binder interface attributed to debonding and found that the solid propellant will become less stiff as voids grow between the matrix and particles. Even though this study did not include a damage model, the results support Yun et al. [54]'s findings since a decrease in the stiffness were also obtained as evidenced by a decrease in the instantaneous elastic modulus E 0 with increased reflected shock wave pressure. Increasing the applied reflected shock wave pressure caused an increase in the strain and therefore an increase in debonding of the particles. Even though E 0 decreased with increased reflected shock wave pressure, the instantaneous elastic moduli values from this study were higher than the range reported in the literature for HTPB propellant, which is 1 MPa-10 MPa [6,20,21,46,54].
Higher moduli values are expected since the Prony series used during the iFEA did not include the softening phenomenon due to debonding of the particles [54]. The supplementary material contains a video that visually captures macroscopic failure due to large deformations, which manifests as cracks if the reflected shock wave pressure applied on the surface of the HTPB propellant sample is higher than 313 kPa. Future iFEA simulations will consider a Prony series with a single Maxwell branch that is coupled with a damage model [9], to better characterize the high strain rate viscoelastic response of HTPB propellant under bending loading. The ability of a nonlinear viscoelastic model that utilizes a nonlinear spring instead of a linear spring, as was done by Bentil et al. [49], will also be considered in the future as a means of characterizing damage due to particulate debonding.
Since the HTPB propellant samples were subjected to large displacements, a combination of plastic deformation and dewetting may have prevented the samples from returning to their original configuration. As a result, the time history plots of W from 3D-DIC did not have an apex that returned to the undeformed position of zero out-of-plane displacements. Additional evidence of plastic deformation and/or dewetting was found in the nonlinear waveform for the out-of-plane displacement over time, which became less defined if the initial reflected shock wave pressure magnitude was high. The contributing factor (i. e. plastic deformation and/or dewetting) causing the permanent change in the out-of-plane displacement time history will be examined in a future study. Energy dissipation within the first 5 ms following shock wave exposure contributed to the complex deformation time histories obtained from the dynamic inflation experiments.
Limitations of this study, which influenced the reported parameters for the Prony series, were attributed to the experimental setup and also the error from the 3D-DIC measurements. For instance, since inert HTPB propellant is a viscoelastic material, the mechanical behavior is strain rate dependent. Thus, the reported values for the optimized Prony series parameters pertain to 5-mm thick inert HTPB propellant samples subjected to high strain rates on the order of 10 4 s À 1 . Additionally, the experiments were conducted on an inert solid propellant using ammonium sulfate in place of the oxidizer ammonium perchlorate. Future dynamic inflation experiments using solid propellants with oxidizers are needed, for comparison of the linear viscoelastic behavior as a function of oxidizer.
The time history for the shock wave pressures that were measured by the pressure transducers was repeatable. However, the results presented are dependent upon the applied shock wave pressures on the surface of the inert HTPB propellant sample. For instance, the oxyacetylenedriven shock tube in this study produced shock wave pressures and also negative pressures, which both influenced the out-of-plane displacement behavior of the solid propellant. Negative pressure caused the HTPB propellant sample to inflate inside the shock tube, while the reflected shock wave pressure caused inflation outside of the shock tube. Therefore, dynamic inflation experiments using a shock tube that is not capable of generating negative pressures may influence the linear viscoelastic behavior by reducing debonding since there would be an absence of an additional tensile strain caused by the negative pressure phase.
During 3D-DIC analysis, the projection error was used as a metric to quantify the correlation accuracy (i. e. displacement measurement error). Factors that influence the projection error were optics (e. g. depth of field, lens focus, lighting, etc.) and the subset and step size. The default subset and step size values in VIC-3D generated an average projection error and standard deviation of 6.7 � 0.08 %, across all eight HTPB propellant samples. Decreasing the projection error would require increasing the depth of field to minimize blurred images caused by large deformations. Since optics influence, the calibration score, recalibrating the cameras using the rigid calibration target would also influence the projection error. However, the calibration score for this study was better than the suggested value of 0.05 by the VIC-3D [41] software. Thus, recalibrating is not expected to significantly decrease the projection error.

Conclusion
This study is one of the first to report the linear viscoelastic and tensile behavior of inert solid propellants with strain rates on the order of 10 4 s À 1 . Additionally, this study has provided data quantifying how well a generalized Maxwell model could describe the deformation of inert propellants exposed to shock waves. The experimental results presented will add to the limited knowledge of solid propellants at high strain rates, which may aid in improving the predictive capabilities of health-monitoring sensors and inform future work involving the development of a more elaborate constitutive model that can capture the complex micromechanics of the propellant through the inclusion of cohesive zone models. As such, the presented findings may be used to validate models that can assess the solid rocket propellant's structural integrity.

-16
Viscoelastic Properties of Inert Solid Rocket Propellants Exposed to a Shock Wave