Explicit Modelling of the Ignition Transient Structural Response of a Paraffin Wax Hybrid Rocket Motor Fuel Grain

Paraffin wax has been identified as a hybrid rocket motor fuel, which offers enhanced regression rates and improved combustion performance. While various investigations into the performance of this class of fuels are being conducted around the world, the consideration of its structural performance is often overlooked. The research presented here establishes a simplified, yet accurate method of defining the structural performance of a paraffin wax hybrid fuel grain to be introduced early in the design phase of a motor. The use of the Johnson–Cook (J–C) material model has been verified to work within the “low speed” ignition range experienced in paraffin wax/N2O hybrid motors, and therefore is used to predict failure in a variety of motors. The resultant stress profiles within the grains indicate that the grain outer to inner diameter (OD/ID) ratio, as well as the outer diameter (OD) itself, play an important role in the grain ability to withstand the loading conditions applied. Additionally, the grain structural properties, and the stiffness of the combustion chamber affect the severity of the internal stresses in the grain. The feasibility of large-scale pure paraffin wax grains without structural enhancement additives is thus found to be poor. Fuel additives should be considered for structural enhancement.


INTRODUCTION
The structural performance of a solid fuel grain is a necessary factor to consider for the propulsion of any solid rocket motor. In the past, solid rocket fuel structural performance has been investigated to reduce any risk of in-flight failure due to the development of cracks from flight loading and fuel storage (Fitzgerals and Hufferd 1971;Ho 2010;Kelley 1969;NASA 1973). A simple crack in the body of a solid rocket fuel is considered a structural failure due to the risk of exceeding the design grain surface exposure. The rapid exposure of combustible surface area, in this case, where the fuel and oxidiser are premixed, can result in a catastrophic failure of the vehicle due to over pressurisation (Gondouin 1993). The boundary layer combustion mechanism of hybrid rocket motors means that there is less risk associated with crack formation in terms of vehicle loss. However, there are still concerns attributed to fuel loss, which can result in an incomplete mission if fuel fragments were to dissociate from the grain body. Additionally, if the fuel fragments were to block the nozzle, over pressurisation of the combustion chamber may occur, and the vehicle may be lost. While a few small surface cracks may not cause failure, many small cracks or few large cracks can cause a potential failure situation in hybrid fuel grains (Sutton and Biblarz 2001).

Structural considerations
Throughout the life cycle of a motor, it undergoes numerous loads, which can be split into two categories. These are the specified loads and the induced loads (Fitzgerals and Hufferd 1971). These loads are fixed by the mission requirements and the propellant properties, respectively (Kumar et al. 2014;Simo 1987). Specified loads are caused by effects such as environmental temperature, gravity, vibration, shock, transportation and handling, as well as ageing conditions and humidity (Fitzgerals and Hufferd 1971). Induced loads are influenced by the propellant selection, manufacturing methods, cure shrinkage, chamber pressure, and the combined loads during flight (Fitzgerals and Hufferd 1971).
The loading conditions which are considered in the present research are the pressurisation loads. The chamber pressurisation is a result of the combustion process within the combustion chamber of the motor. This begins at the ignition of the motor and terminates at the end of the burn sequence, with compressive hydrostatic loading on the grain surface, as well as tensile hoop stresses (Fitzgerals and Hufferd 1971). Ignition pressurisation can be significantly higher than that of the mean operating pressure, and this can cause grain failure before launch; however, prolonged exposure to higher temperatures and high pressures can result in mid-flight grain failure too. Chamber pressure is generally considered a static load with a magnitude equal to that of the initial maximum pressure measured from a static firing, or a modelled theoretical prediction (Ho 2010). Depending on the propellant combination and ignition sequence, ignition pressurization times can vary from 0.005 to 0.5 s (NASA 1973). The speed of this transient pressure load can affect the grain structural response, and, thus, strain rate dependence must be investigated during the grain material characterisation.
In general, the fuel material is stronger but more brittle at lower temperatures, such as 20 °C, than it is at temperatures as low as 30 °C (Veale et al. 2018). Low-temperature firing could be a concern due to the rigidity of the propellant, residual cure stresses, and the high ignition pressures. If a fuel grain was to remain intact after the pressurisation load, any resulting flaws and cracks might fail after exposure to a temperature increase and the in-flight inertial loading. When analysing the effects of in-flight combined loads, the effect of thermal softening should be considered (Fitzgerals and Hufferd 1971). A number of paraffin structural enhancement techniques with the use of polymer additives have been researched by various institutions (Mengu and Kumar 2018) which generally show that these additives reduce the regression rate of the fuel, often by increasing the melt viscosity, resulting in reduced fuel entrainment. Other research has been conducted into the effect that structural discontinuities will have on the combustion performance of a motor (Andrianov et al. 2019). The preliminary analysis shows that induced cracks do not lead to immediate failure of the grain, but there is a clear regression inconsistency which can lead to chamber burn through if the burn time is long enough. Innovative structural matrices embedded within the grain are also being investigated in an attempt to minimise the inclusion of regression dampening materials. Initial tests indicate a reduced regression, yet the performance and structural properties are matrix material and geometry dependant (Hill et al. 2019).
Paraffin wax is considered a brittle material, which has a very low ultimate tensile strength (UTS). This will likely reduce the ability of the fuels to withstand the pressurisation and inertial loads, without cracking, during flight. The structural performance of a fuel grain is dependent on a number of factors, such as the size, length-to-diameter ratio and outer to inner diameter (OD/ID) ratio of the fuel grain. These design parameters of the fuel grain are often determined during the design phase to control the thrust, oxidiser to fuel ratio (O/F) ratio and mass flow rate of the motor (Karabeyoglu et al. 2001). Karabeyoglu (2011) has determined that the failure boundary is approached more rapidly with fuel grains with larger OD/ID ratios. In addition to the low strength of paraffin wax, the material has a low melting point and a large melting onset temperature range. This can be a problematic aspect to consider during analysis to ensure that the thermal penetration into the grain during combustion does not cause the fuel to prematurely melt or soften significantly.
Comprehensive material characterisations tests have been performed in an earlier publication (Veale et al. 2018). All test data used in this publication refer to this work. These tests were performed for pure paraffin wax (PW00) and paraffin wax with 40 Wt% aluminium powder (PW40AL). This publication details the tensile and compression testing conducted at various strain rates and temperatures. These results show a distinct strain rate dependency of paraffin wax. Additionally, small temperature changes in the material result in a significantly weaker, more ductile material. This is a result of the long transition temperature before the material begins to melt, which initiate at temperatures as low as 30 °C.

MATERIAL MODELS
The strain-rate dependency of a material generally results in an increased yield stress at a lower strain at high strain rate.
Generally, standard finite element analyses (FEA) packages make use of quasi-static material properties when considering standard use failures, meaning that the internal forces of the material can be neglected. In the case of impact analysis, where there would be a rapid deformation of the material, such as an explosion, these internal forces significantly change the response of the material. Simulations on ignition pressurisation would be more appropriate if the strain rate dependency of the material can be accurately modelled.
Due to the low melting point, and low melting onset temperature, the thermal effects on the structural performance of paraffin wax are evident at temperatures as low as 30 °C. Paraffin wax hybrid motors burn within a boundary layer formed between a melted fuel layer and the oxidiser stream. The liquid layer is constantly being stripped away throughout the burn and replenished with newly melted fuel. The thermal penetration into the fuel grain has been shown to be quite slow, considering the temperatures involved in combustion (Santi et al. 2017). This is primarily due to the low thermal conductivity of the paraffin wax and the repeated formation of the insulating liquid layer. The tests conducted by Santi et al. (2017) were only considered for pure paraffin wax, while aluminised paraffin wax has a thermal conductivity coefficient nearing a factor of three higher (Veale et al. 2018). While some early thermal simulations indicate that the expected thermal wave penetration into an aluminised grain is still not faster than the fuel regression, thus not causing slumping, this property should still be considered in the case of high ambient temperature launches.
The Johnson-Cook (J-C) constitutive model (Johnson and Cook 1983) was developed to determine an equivalent stress (σ ep ) as a function of plastic strain, strain rate and temperature in metals. The use of this model is designed specifically for a range of metals, but is known to have been used in some plastic applications within a verified tested range. The model is represented as follows (Eq. 1): (1) where A, B, n, c, and m are the empirically determined material constants, ε p is the accumulated plastic strain, ε • * is the dimensionless strain rate and T* is the homologous temperature. The material constants are defined as: A: quasi-static yield stress, B: hardening constant, c: strain rate dependency constant, n: hardening exponent and m: thermal softening constant. Also (Eqs. 2 and 3), where ε • p is the plastic strain rate, ε • 0 is the user-defined reference strain rate, T is the working temperature, T 0 is the room temperature, and m is the melting temperature.
The five constants presented in the J-C equation (Eq. 1) are to be determined empirically. The constant A is the yield stress of the material determined from a quasi-static tensile test. B and n are the hardening constants and are defined as the power fit curve of the true stress plastic strain curve after yield has occurred. The constant c represents the strain rate dependence of the material and m is the thermal softening constant.
There are several published methods to determine the constants necessary for the J-C constitutive model, such as optical strain rate measurements (Stopel and Skibicki 2016), ballistic impact (Burley et al. 2018), Taylor anvil testing (Šlais et al. 2012) and tensile testing (Banerjee et al. 2015;Dehgolan et al. 2016). In general, the agreement between the methods requires material testing to be conducted at both quasi-static strain-rates, high strain-rates, and elevated temperatures. The J-C model was originally developed for use on a specific range of materials, such as copper, brass, iron, aluminium and steel. Often, FEA modelling for plastics makes use of a bilinear plastic material model which does not take into consideration the strain-rate and temperature dependency of the material. In certain cases, the J-C material model can be used on other materials, such as polymers, with verification for each case (Louche et al. 2009). To determine the J-C model constants, the true stress and strain curves from torsion and tensile tests at various strain rates and temperatures were modified with the Bridgeman correction and presented in equivalent tensile flow form. The Bridgman correction (Bridgman 1952) is used to determine the equivalent stress of a round tensile specimen after necking occurs due to the tri-axial nature of the loading condition at that time. In the case of brittle materials, which demonstrate negligible necking, this correction is not necessary as it is a function of the necking diameter. Other researchers, such as Banerjee et al. (2014Banerjee et al. ( , 2015, make use of only tensile testing to determine the constants required. This simplified approach for determining the constants is utilised in this work. The validity of this assumption is verified through FEA modelling, as discussed later in this work. In addition to the J-C material model, there is also the J-C failure model (Johnson and Cook 1985). This is a modified version of the criterion proposed by Hancock and Mackenzie (1976) to include factors such as stress triaxiality, strain rate and temperature effects. This model defines material failure at the point where the accumulated damage (D) during plastic straining reaches a critical value. This is defined as follows (Eq. 4): (4) where the fracture strain is ε f , the damage threshold is ε p,d and the critical damage is D c .
The J-C failure model define the fracture strain as Eq. 5: where the D 1 is the strain at the initiation of necking for the quasi-static test, D 2 and D 3 define the variation in the failure strain in relation to the triaxiality (σ*) of the material, D 4 is the strain rate dependent failure parameter, and D 5 defines the effect of temperature on the failure strain. For the purposes of this analysis, the triaxiality constants are considered negligible (Hancock and Mackenzie 1976) as the tensile testing conducted showed effectively no necking.

Determining the J-C constants
The determination of the J-C material model constants required a range of tensile tests to be conducted at a range of strain rates and temperatures. In the case of paraffin wax, which is a brittle microcrystalline wax, the temperature range necessary for testing was quite small when compared to other materials, which are typically described by the J-C material model. In the case of the quasi-static testing, the strain rate was selected where the acceleration effects on the load measuring devices were negligible, and where it was not necessary to consider the propagation of stress waves through the specimen.
The quasi-static testing was conducted at a strain rate of 2 × 10 -2 s -1 , which was taken as the reference strain rate (ε • 0 ) for determining the other constants. The constant A was defined as the yield strength of the material at this reference strain-rate. In the case of paraffin wax, this is not a clearly defined point and was thus determined using a 0.1% offset yield due to the low elastic modulus and plastic strain. This low offset yield was applied due to the shape of the stress-strain curve, and smooth transition into plastic deformation. The constants B and n were defined by the plastic region of the stress-strain curve. The stress-strain curve was trimmed to only reveal the plastic region, after the yield point. A power law trend line was fitted to these data on the basis of the least fit square approximation and the constants of the resulting equation in the form of Bε n p were used to define the hardening constants. The strain rate dependence constant (c) required testing at a variety of strain rates. In the case of this testing regime, the tests were conducted at the reference strain-rate of 2 × 10 -2 s -1 , as well as at 2 × 10 -1 and 2 s -1 . The constant c was established by plotting the true stress at the yield strain for each strain rate, against ln(ε • *). The gradient of the linear fit that can be derived from these data was used to estimate the constant c. This plot can be seen in Fig. 1. Finally, the thermal softening constant, m, was determined by conducting tensile tests at a variety of elevated temperatures.
Due to the low melting point of paraffin wax, these tests were only conducted at 30 and 40 °C. Equation 1 can be rearranged to give Eq. 6, which shows a linear relationship between ln(T*) and ln(k -σ eq ), where K is the remaining constant from Eq. 1.
The constant m is the gradient of the curve obtained when plotting ln(k -σ eq ) against ln(T*). The stress values at a strain of 0.5% and at a strain rate of 2 × 10 -1 s -1 for each temperature are substituted into this equation, along with the other constants which have already been defined to produce a best fit linear curve. This can be seen in Fig. 2.
All experimental tests used to find these material constants were conducted a minimum of five times with the average points utilised in this work. The experimental aspects can be found in a previous publication (Veale et al. 2018). The J-C constants for both the pure and aluminised version of the paraffin wax and were determined using the above method. These constants, defined in Table 1, were used as the material input criteria for the ANSYS (2018a) explicit finite element model detailed later. The constants utilized in this model required verification prior to the application into the full grain. This was achieved by replicating the tensile tests in ANSYS at various temperatures and strain rates and comparing these results with the experimental results.
The constants required for the failure criteria are determined similarly to that of the material model. The constant D 1 is defined by the point in the quasi-static test where necking is initiated. In the case of paraffin wax, no necking is evident, and therefore, the plastic strain at failure is considered for this variable. For a similar reason, constants D 2 and D 3 are negligible. Constant D 4 , which defines the strain-rate dependence failure, is defined by the gradient of a linear curve fitted to data for plastic strain at failure versus ln(ε • *) at different strain rates. Finally, the temperature-dependent failure is defined by the gradient of a linear curve fitted to data for plastic strain at failure versus T*. The simulated stress-strain data were validated against the experimental data, as described in the next section.

Material model validation
The material constants, which were determined from the experimental data, were then validated through FEA modelling. The tensile test scenario was represented in a 1/8 th symmetry model of a tensile test specimen using the ANSYS Mechanical Explicit Solver (AUTODYN).
The size of the model could be reduced as a result of geometrical and boundary condition symmetry. The tensile specimen was discretised with a finer mesh in the active region (region of the extensometer), and a coarse mesh in the grip region. Due to the simplicity of the geometry, and loading, a full hex mesh was applied to the model in order to reduce computational expense. First order hexahedral elements are susceptible to hour-glassing affects when experiencing large deflections, and thus hourglass dampening controls were required in this validation model to prevent unrealistic deformation of the elements. The hourglass dampening control utilized in this model was the Flanagan Belytschko method (Flanagan and Belytschko 1981), with a stiffness coefficient of 0.03 and a viscous coefficient of 0.1.
Finite element models employing explicit solution regimes are computationally expensive, thus, the model duration is normally minimised and utilised for impact or shock loading situations. In the case of modelling experimental tensile tests, the modelled duration, especially in the case of quasi-static testing, is significantly longer than what would normally be represented. The time taken to run a simulation is determined by the number, size and type of elements employed, as well as the time step applied. In this case, the shockwave propagation within the material needs to be considered so that the resultant shockwave front is not able to travel further than the length of an element for each time step (ANSYS 2018b). This ensures that the reaction of the material to loading and deformation is considered. The time step (Δt) size is therefore calculated according to the ratio of the material speed of sound (C) and the length of the element (l), as described in Eq. 7 (ANSYS 2018b). (7) where C is defined using (Eq. 8) and the Modulus of Elasticity (E) and the density (ρ)   (8) The meshed tensile test specimen model can be seen in Fig. 3. The mesh is split into two parts, which are joined by connection nodes. These parts define the grip area and the test area between the grips. This allows the velocity load to be applied to the whole grip area. A 50 mm extensometer was used to measure the strain of the region of interest. A comparison of simulated stressstrain data and experimental data is presented in Fig. 4, and provides an indication to the level of agreement between data sets. The plots in Figs. 4 and 5 show the difference between the empirical tensile data and the simulated data. The shapes of the stressstrain curved between the experimental and simulated results differ, with a discrepancy being particularly evident near the yield point, while the elastic modulus, plastic region, as well as the failure point, correlate closely. The reason behind this discrepancy is the development of the J-C model based on materials with more clearly defined transitions between the elastic and plastic regions, while the transition in paraffin wax is comparatively extended over a strain range. As described, the definition of the first three constants in the J-C constitutive model (Eq. 1) represent the yield point (A), and the power law that defines the plastic region (Bε n n ). The results presented indicate the plausibility of using the J-C model within the regime presented in the tensile tests, and may therefore be applied for the purpose of estimating the structural integrity of a paraffin wax fuel grain as the point of interest is the failure.  Table 2. The various geometrical and performance properties of these vehicles were determined using the internally developed Hybrid Rocket Performance Simulator (HYROPS). More information on the design process and the use of HYROPS can be found in the references Balmogim et al. (2015), Genevieve et al. (2015) and Leverone et al. (2013). The ratio of the outer diameter (OD) (b) to the inner diameter (ID) (a) (OD/ID ratio) has been used as a general design rule to reduce internal stresses of the grain as per Fig. 6, with an understanding that the material approaches the failure boundary as the 9 OD/ID ratio (or b/a as depicted in the figure) increases (Majdalani 2011). The radial distance, r/a, is the ratio of any point (r) within the grain thickness to the port diameter (a). In Fig. 6, the simulated stress profile within an unspecified grain geometry is similar to that of a thick-walled pressure vessel with an internal pressure, P i , with some reaction given by the motor casing. The stress distribution over the thickness of the grain, without a reaction from the motor casing, is graphically presented in Fig. 7a. It can be seen that the tangential (or hoop) stresses in this schematic are positive, while the radial stresses are negative, with the maximum in both cases occurring closest to the port. When comparing this result to that presented in Fig. 6, it can be seen that, in some cases, the hoop stresses become compressive as a reaction to the motor casing. The extent of this reaction force (P o ) depends greatly on the material of the fuel and combustion chamber, as well as the internal pressure. When there is a large material stiffness difference between the grain and combustion chamber, it is estimated, based on multi-layered thick-walled pressure vessels design calculations, that the stress distribution will look more similar to that presented in Fig. 7b, where both the tangential and radial stresses are compressive, while still exhibiting a similar trend. This means that a smaller OD/ID ratio may result in primarily compressive loads within the motor, reducing the likelihood of failure. It is also known that an increased fuel grain OD further increases the maximum hoop stresses, and therefore, the OD/ID ratio should be considered in conjunction with the maximum OD of the grain.

Model preparation
The simulations considered apply the J-C material and J-C failure models onto the paraffin wax fuel grain encased in an aluminium combustion chamber. The use of the J-C model was verified to account for the nonlinear properties and strain-rate dependence of paraffin wax within the tested range. The modelling regimes, at this stage, only consider the ignition transient pressure loading, which, in effect, allows the thermal effects, regression of the grain and the inertial loads to be neglected.
The motors listed above were modelled as full bodies, with only the grain, post-and pre-combustion chambers, and the motor casing present. An explicit model implemented in ANSYS Explicit (AUTODYN) was utilised due to the ability of the software to model material responses to time dependant, high impact loads over very short time frames. In all the cases above, the combustion chamber was longer than the grain, as it would be in reality. The pre-and post-combustion chambers were included as their influence on the grain would need to be considered. The thermal liner was removed from the model at this early stage of the modelling regime, as the complex nature of composite materials would affect the interpretation and understanding of the grain response. It should, however, be included for any quantitative structural modelling after adequate material properties are obtained.
The fore-end of the motor housing was fully constrained in all degrees of freedom, while the aft-end was allowed freedom in the axial direction. This allows the chamber length to expand and contract axially, as it would be able to in reality, while restraining the ends in the radial direction due to the stiff bulkheads. The fuel grain and pre-and post-combustion chambers were set to move freely within the motor housing, however, not able to pass the axial limit of the combustion chamber. This is a representative limit enforced by the bulkheads on either side of the motor. This constraint was applied as the parts are typically loosely fitted within the combustion chamber during manufacture and assembly. Of course, this constraint is based on the manufacturing method applied in these specific motors and may vary depending on the manufacturing and assembly methods used in other designs. An effective pressure load was applied as a stepwise, time dependent load to all internal surfaces of the combustion chamber. The pressure data applied to the model was determined from pressure transducers placed at the fore-end of the motor during experimental hot fire testing as well as the average chamber pressure modelled from HYROPS. The application of a constant pressure profile is an over simplification of the complex pressure distribution within an igniting hybrid rocket motor. The simplification allows the analyses of the resultant pressure profile and application of the J-C models before applying complexity. Body interactions between the grain, pre-and post-combustion chambers and the motor housing were activated to ensure that the varying material deformation and reaction forces were adequately accounted for. The load application was similar to that presented by Hsiao (2013) in his work concerning solid propellant grains under ignition pressurisation.
A hex mesh, with mapped surfaces, was applied to all parts for computational simplicity, which means that hourglass dampening, similar to that described for the tensile test models, needed to be applied. Mesh size independency and quality were verified for each model before final simulation results were documented. A portion of the mesh used for the P-1A simulation can be seen in Fig. 8.
A variety of scenarios were simulated for each motor depending on what pressure data was available for each motor. These scenarios included actual pressure data from available hot-fire tests, the modelled theoretical pressure data based on the HYROPS software estimations (including an estimated pressure spike), and a reduced pressure profile where the chamber pressure does not exceed the design pressure (excluding the pressure spike). All pressure profiles are available in references Balmogim et al. (2015), Genevieve et al. (2012Genevieve et al. ( , 2015 and Leverone et al. (2013Leverone et al. ( , 2019. It is important to note that the pressure profiles and start-up transients greatly affect the ability of paraffin wax fuel grains to withstand the loads applied, as the material is noticeably strain-rate dependent.

Case study
Two sets of analyses were considered in this study. The primary analyses considered the existing Phoenix sounding rocket motors and their susceptibility to failure given their modelled, and, if available, tested ignition pressure profiles. The vast difference in the modelled and experimentally measured pressure profiles, as well as the difference in material response between the pure and aluminised waxes, affected the structural response of each grain. Table 3 shows a summary of grain failure predictions in the respective simulations based on the J-C failure model. Due to the method of load application, if failure of the surface was to occur, and the subsequent elements began to erode, the loading profile would no longer be applied as intended, and thus the result would be at the point where failure occurs. This means that the simulation method utilised here can be used to predict if failure was to occur, but not its extent. Having said this, the prediction of failure initiation is an important exercise, which can lay the foundation for the development of a computationally intensive explosive model to more completely model the expected pressure application. In addition, it should be noted that based on the results discussed below, the maximum strain-rates for all the scenarios simulated fall within the range of the test data, meaning that the "slow" rate of applied load in hybrid ignition makes the use of the J-C model suitable. The aim of performing the secondary set of analyses was to investigate the effect that fuel grain geometry has on the stress profile within each of the grains of the motor when subjected to the same pressure loading profile. This was achieved by applying the same pressure profile on each of the described grains, and by applying the same pressure profile to a grain of constant OD, but with a modified OD/ID ratio. In this case, secondary influencing bodies, such as the pre-and post-combustion chambers, were removed so that even pressure could be applied to all exposed grain surfaces, and to ensure there were no external influences when interpreting the stress profiles. The results obtained show the influences grain geometries have on the stress profile of the fuel grain, which would influence the design considerations at an early stage.

P-1A simulation
A series of images are presented in Fig. 9 detailing the simulated von Mises equivalent stress profiles within the pure wax grain of the P-1A motor, at the end of the ignition time during the hot-fire motor test. Due to the lower maximum pressure in this test, the resultant stresses are relatively low. A peak equivalent stress of 2.4 MPa is predicted near the aft-end of the grain where the post-combustion chamber ring interacts with the fuel grain and casing. There was no evidence of any form of failure in this motor after this hot-fire test nor in this ignition simulation, lending validation to the proposed model.   In the modelled loading scenario of the P-1A, the artificial pressure spike of 15% (47 bar) (Genevieve et al. 2012) over the start-up transient of 0.027 s resulted in grain failure at 0.0243 s. At this time, the material exceeded its maximum allowable stress at the effective strain rate. The maximum equivalent stress at the moment before failure is 3.3 MPa. The simulation of the modelled pressure transient without the artificial pressure spike did not result in grain failure, but resulted in a maximum equivalent stress of 3.2 MPa. When comparing these two equivalent stresses, it can be seen that the simulation that did not exhibit failure reached an equivalent stress very close to that of the simulation that did fail. This reinforces concerns regarding the design pressure and the attributed risk associated with ignition spikes.
As discussed, structural response predictions, after the initiation of failure occurs, are not likely to be accurate, as some nodes through which loads are transferred fail, and are removed from the active model, affecting the intended load application. While there are various ways to account for this, the purpose of this work is to determine the effectiveness of a simplified analysis for the purposes of verifying a fuel grain design. It is, however, interesting to visualise where the grain failure is initiated. Figure 10 shows the initial predicted failure, and the propagation of the failure of the grain as determined by the J-C failure model, which eliminates elements where the damage threshold is exceeded. The failure initiates at a few points at mid-length on the OD of the grain and then propagates diagonally around the OD. At some time after this initial failure, as the pressure continues to increase, the grain fails on the surfaces in direct contact with the applied ignition pressure. At this time, it is apparent that the grain is failing as a direct response to the surface pressure on the port surface and the fore and aft-ends.

P-1B simulation
The P-1B motor grain was predicted to react similarly to that of the P-1A motor. This is attributed to the similar grain outer dimensions. The OD/ID ratio of the P-1B grain is only slightly more than that of the P-1A grain, which implies that this may result in a higher stress profile; however, the reduced OD also acts to reduce this stress profile slightly. Additionally, the actual hot-fire start-up transient time of the P-1B motor is more than double that of the P-1A motor, which affects the material response due to its strain rate dependency. Figure 10 shows a cross-section of the equivalent stress profile for the P-1B grain under the actual loading conditions measured during its hot-fire test. A maximum equivalent stress of 2.5 MPa is predicted. The simulation with the pressure spike (43 bar)  resulted in grain failure with a maximum equivalent stress of 3.05 MPa measured just prior to failure at 0.0165 s. When the modelled pressure scenario was simulated without the pressure spike, the grain did not fail and resulted in a maximum equivalent stress of 2.95 MPa at 0.02 s. These values are not directly comparable to the P-1A results due to the different pressure profiles.  Figure 11. P-1B actual pressure profile, equivalent stress, side cross-section.

P-2A simulation
The P-2A motor has a larger OD when compared to the previous two motors. While the fuel thicknesses are similar, the tangential stresses within this motor are expected to be larger due to the increased diameter. On the other hand, however, the lower OD/ID ratio should counteract this effect. With the modelled loading condition determined during the design of the P-2A motor (Leverone et al. 2013), with no ignition pressure spike, failure was predicted to occur by the associated simulation at 0.0086 s into the start-up transient.
The result presented in Fig. 11 shows that the maximum equivalent stress prior to failure was predicted to be 3.38 MPa, which is significantly lower than the failure thresholds of the smaller motors. While the rate at which the load is applied does affect the vehicle ability to withstand the stresses, the low toughness of this material indicates that the low OD/ID ratio of 1.44 is not enough to reduce the increased tangential stresses which are present due to the increase in diameter. This suggests that larger diameter grains should likely be strengthened with additives in order to resist premature failure.

P-2A II simulation
The P-2A II (aluminised) motor has similar dimensions to that of the P-2A motor, except that it is shorter as a result of the increase in the density specific impulse, and lower optimal O/F ratio of the fuel material. The material properties of the aluminised wax were considered in this simulation. Based on the material testing described above, it is clear that the aluminised wax considered in this study is notably more brittle than the pure wax. A simulation of the P-2A II fuel grain response to the modelled pressure profile, with no pressure spike, predicted the initiation of failure at 0.00682 s at a peak equivalent stress of 2.99 MPa. The cross section of the equivalent stress profile can be seen in Fig. 12. At higher strain-rates, the aluminised wax is noticeably more brittle than the pure paraffin wax, having a higher strength with less allowable deformation. The rapid pressure load therefore results in a rapidly increasing stress throughout. This is evident in the image showing the point of grain failure, Fig. 14, where the wax grain is predicted to experience many internal failure initiations throughout its volume, as opposed to the localised failure, crack propagation predicted in the pure wax grains. While this image does not demonstrate the actual failure expected, it is interesting to note the predicted differences between the failure profiles of pure wax and aluminised wax indicating an increased brittleness of the aluminised grain.

Summary of results for existing motor designs
From the simulations presented above, it is difficult to define a direct correlation between loading and geometry due to the number of variables associated with failure. Table 4 lists the various peak stresses and strains for each motor at either the end of the simulation or just prior to failure. It is noted that all principal stresses are in compression. When considering the peak equivalent strains associated with each of the grains considered, it is clear that the higher the OD/ID ratio, the higher the strain experienced at the end of the simulation. This correlation is not evident in respect of the tangential stresses, which decrease as the OD/ID ratio is lowered, and increase as the OD is raised. It appears that the increasing OD has more of an effect on the stress profile changes in the OD/ID ratio.

Effect of geometry on stress profile
It is clear from the results above that the radial and tangential stresses are, in fact, fully negative. However, to fully understand the estimated stress profile presented in Fig. 7 (b), these solutions cannot be directly utilised as numerous variables affect the result, including, primarily, the maximum pressure and the rate of loading. To have comparable results based solely on the geometry of the motor, a simplified loading scenario was simulated for each motor of different geometry by applying the same fictional loading condition to each motor. The maximum pressure simulated was 40 bar, reached at a start-up time of 0.02 s. The preand post-combustion chambers were removed so that the entire system was simplified and only predicted the direct response between the applied pressure and grain dimensions. However, the fore and aft-end pressures were still applied. The radial and tangential stresses predicted for each scenario are plotted in Fig. 14, showing the stress profiles through the thickness of each grain. A similar trend to that presented by Karabeyoglu (2011) was observed, except that all the predicted resultant stresses in this case are in compression. The maximum equivalent stress is predicted to be 2.72 MPa for P-1A grain, 3 MPa for the P1-B grain and 2.81MPa for the P-2A grain, while the associated maximum equivalent strain is predicted to be 0.004, 0.0052 and 0.0043 mm·mm -1 , respectively.
Radial Distance, r/a P-1B Tangential P-1B Radial P-1A Tangential P-1A Radial P-2A Tangential  The work presented by Karabeyoglu (2011) in Fig. 6 shows that for grains with the same OD, the grains with higher OD/ID ratio have more positive tangential and radial stresses, while the stresses predicted in Fig. 14 are not ordered in the same fashion.
This is attributed to the effect of the different OD of each grain in the presented work, indicating that both the OD/ID ratio and the absolute OD have an effect on the stress profile within the fuel grain. The predicted stresses in this graph indicate that increasing OD negatively affects the structural integrity of the grain design, and reducing the OD/ID ratio may not be sufficient to prevent failure.
To verify this consideration, an additional simulation was conducted on fictional versions of the P-1A fuel grain where the OD remained the same, and the port diameter was modified to either reduce or increase the OD/ID ratio. The predicted tangential and radial stresses of these simulations can be seen in Fig. 15. There is a distinct similarity to the result presented in Fig. 6, also showing that the higher the OD/ID ratio, the closer the grain approaches failure. The change in the tangential stress profiles for each grain with varying ratios is quite large, while the effect of an increasing OD has an opposite effect.  Figure 16. Tangential and radial stress vs. radial distance for P-1A grain with various OD/ID ratios.

CONCLUSION
While it is difficult to ascertain a direct correlation between the results obtained and the geometry of a motor grain given the current data set, these results do lead to the conclusion that that the J-C material and failure model and modelling method are suitable in the grain structural integrity validation during the initial design phases. Fuel grain dimensions, port size and the chamber pressure are determined through the mission objectives and design constraints. The usefulness of applying the J-C material and failure models to capture the strain rate dependency of a material that is not metallic, such as paraffin wax, is limited. However, the data and testing presented in this work show that, within the range specified, these models can be utilised for determining the grain material response and in particular the failure point.
Previous fuel grain stress estimations have considered the fuel grain to be a constrained, thick-walled pressure vessel, with a different material encapsulating it and restricting its deformation, resulting in compressive stress profiles. In the models presented in this paper, the addition of a fore-and aft-end loading condition is necessary for an accurate determination of the stress profile. While it could be said that the pressure loading imposed on the fuel grain should not be applied evenly to all wetted surfaces, as the chamber pressure tends to decrease from the fore-to aft-end, this effect was neglected at this early stage in the analysis. Unsymmetrical loading would add significant complexity for a very small increase in accuracy with the interpretation of the resultant stresses and applicability of the J-C material models at this early stage of development. This being said, the next phase of analysis and simulation should allow for this attribute to be included. When analysing the principal radial, tangential and longitudinal stresses in each design case there is no distinct response pattern, as the same loading condition is not applied to each associated geometry. These results do, however, highlight how these simulation methods offer more insight into the grain response to unintentional pressure spikes, or geometric attributes.
A hypothetical set of analyses was conducted with the same pressure profile applied to each geometry to visualise the direct effect that the geometry has on the stress profile. It is clear from these results that, while the OD/ID ratio does have an effect on the point of failure, the maximum OD of the motor is likely to have a larger effect.
The estimated tangential and radial stress profiles presented in this work show a similar trend to what would be expected in a thick-walled pressure vessel. The tangential stress profile is, however, in compression due to the constraining force imposed by the combustion chamber, which has vastly different material properties. This constraining force will depend on the material and geometric characteristics of the combustion chamber and any thermal liner, if present.
The results presented here aim to validate a procedure for determining the structural response of paraffin wax hybrid rocket motor fuel grains by means of simulations utilising the J-C material and failure models within a verified range of testing. While it is impossible to determine the validity of each grain simulation without extensive testing, the implication from the stress profile shows the plausibility of utilising the J-C models for more grain analyses. Full validation of these results can only be achieved with extensive motor testing, which is beyond the scope of this paper; however, it should be noted that the simulations which utilised data from actual hot-fire tests did not predict failure, and no failure was discovered after the tests. Additionally, should the presented modelling methods be utilized, not only for the initial design stages, but additionally for the final design verification, the ignition pressure profile should be applied with an expanding gas model or explosive model to more realistically apply the loads inside the combustion chamber after failure occurs. This would allow the pressure profile to be applied to the exposed grain surface independent of surface failure, and where complex pressure gradients throughout the grain length can be applied more realistically. Further studies into this topic should tackle the material characterisation of a wider range of mechanically strengthened paraffin wax grains, and the improvements in geometrical properties that can be obtained with these structural enhancements with their associated effects on rocket motor performance.

DATA AVAILABILITY STATEMENT
All data sets were generated or analysed in the current study.