A novel methodology for predicting the high rate mechanical response of polymers from low rate data: application to (plasticised) poly(vinyl chloride)

Plasticised and unplasticised poly(vinyl chloride) (PVC) are used as engineering materials in many applications where they may be subjected to impact loading leading to high strain rate deformation at a variety of temperatures. It is therefore necessary to study the mechanical responses of these and similar materials over a range of loading conditions, especially as they exhibit strong rate and temperature dependence, and could include a low temperature brittle transition. In this paper, a model of the mechanical response of a PVC with 20 wt% plasticiser and one with no plasticiser is applied over a wide range of strain rates and strains and shown to have excellent agreement with experiments conducted in a previous study. As it is challenging to obtain high rate data on rubbery materials using conventional apparatus, such as the split-Hopkinson pressure bar (SHPB), an alternative approach is presented based on a novel modelling framework, which uses the time–temperature superposition principle and is fully calibrated using quasi-static experiments at different temperatures.

sitivity increases with increasing rate and decreasing temperature due to inhibitions in motions associated with the secondary (β) transitions (Bauwens-Crowet et al. 1969;Walley and Field 1994). However, at high strain rates (above c. 10 2 s −1 ), experimental characterisation of the behaviour becomes increasingly difficult. Conventional techniques such as the split-Hopkinson pressure bar (SHPB), which were developed for the characterisation of metals, do not give accurate measurements of the modulus due to experimental artefacts. Static equilibrium is required for a valid analysis; however, the low modulus and resulting wavespeed mean that the time taken to achieve this is similar to, or even greater than, the experimental duration. The low specimen strength can also lead to a poor signal-to-noise ratio (Siviour and Jordan 2016;Siviour 2017). Several modifications have been proposed to make the SHPB more appropriate for testing soft materials. These include lowering the impedence of the bars to match that of the specimen, reducing the bar cross-sectional area by hollowing it to increase the transmitted strain signal, and using piezoelectric lead zirconium titanate (PZT) gauges to increase sensitivity to the transmitted force signal (Chen et al. 1999;. Until a more robust experimental technique for characterising the high strain rate properties of soft polymers is developed, one must also explore other alternatives. Timetemperature superposition has been used in the past, as one such alternative, to validate high rate data (Siviour et al. 2005). In this work, it was shown that the inhibitions to micromolecular movement at high rates was analogous to the inhibitions due to low temperature. However, one of the key differences between high and low strain rate loading is that the former is adiabatic, and the specimen temperature increases because of the conversion of plastic work into heat. Previous work has also focussed on experimentally simulating a high rate test by either increasing the temperature during a low temperature isothermal test Siviour 2013, 2015) or through a series of strain rate jump tests (Furmanski et al. 2013). Much of the constitutive modelling has also tended to focus on predicting single values such as the yield stress (Bernard et al. 2018;Kendall and Siviour 2012) or made use of a large number of parameters characterised at high rates (Wang et al. 2017;Österlöf et al. 2014). For this reason, this paper looks to a novel modelling framework that has been developed by building on previous research Siviour 2017, 2018a,b), to predict the full high rate response of the PVC and PPVC via constitutive models calibrated entirely with data obtained accurately at low strain rates and a range of temperatures. These predictions are then compared to the high rate results presented within a previous study ).

Experimental characterisation 2.1 Material
In order to be able to validate the predictive capabilities of the modelling framework presented in this paper, the materials selected were from the same batch as that in the previous study . They are manufactured by Solvay under the trade name SolVin 264PC with and without 20 wt% of DINP plasticiser.
Not only does this material selection allow validation against an existing dataset, but relatively soft, plasticised PVC is one of the materials that would be challenging to test at high strain rates to obtain reliable and reproducible data. Another benefit of using this material is that the mechanical response is dominated by hyperelasticity at the lowest strain rates, or elevated temperatures, as well as at large strains.
Having shown this modelling framework to be promising with a more plasticised 60 wt% PPVC (Trivedi and Siviour 2018b), this less plasticised PVC provides an interesting trial of the modelling framework as the post-yield softening behaviour is different. The PVC has been selected due to the fact there appear to be fewer complex interactions between the primary glass transition α and the secondary β transitions at the temperatures and rates of interest. It will be shown that the total mechanical behaviour can be considered to be the additive sum of the α and β responses.

Compression experiments
New uniaxial compression experiments were not performed as a full regime of testing was conducted in a previous study . Right circular cylindrical specimens of diameter and length 5 mm were used. The low rate experiments were carried out on an Instron 5980 electromechanical static testing machine, intermediate rate experiments on a hydraulic press, and high rate experiments on an SHPB. These experiments were all performed at a temperature of 25 • C. Varying temperature experiments were also conducted at a lower rate of 10 −2 s −1 using the Instron coupled with a temperature controlled environmental chamber. The results of these uniaxial compression experiments for both PVC and PPVC can be seen in Fig. 1. In this paper, the modelling process is described for the PVC, but the same process is applied for the PPVC to obtain the results presented later.

Time-temperature superposition (TTS)
Dynamic Mechanical Analysis (DMA) experiments were performed on a TA Instruments Q800 in order to quantify the frequency and temperature dependence. Isothermal frequency sweeps were performed over two decades of frequency and a temperature range spanning the two transitions. Rectangular specimens of dimension 5 mm × 10 mm × 60 mm were tested in dual cantilever configuration with a fixed strain amplitude of 0.1%, temperature range of −100 • C to 120 • C, and frequencies of 1, 10 and 100 Hz.
Results from the DMA experiments are shown in Fig. 2(a) for the PVC, where it is clear the α transition occurs around 80 • C with the transition temperature increasing with increasing frequency. This is also observed in the β transition, which can shift from −50 • C to ambient temperature at sufficiently high frequencies. In order to make a comparison between the frequency of the DMA experiment and the rate of an Instron experiment, a conversion ofε ≈ 4f ε 0 , where ε 0 is the strain amplitude of the DMA, was used.
Similar to plotting the temperature dependence of the storage modulus for constant frequencies, the variation of modulus for a constant temperature can be plotted ( Fig. 2(b)). The time-temperature superposition technique postulates that at every temperature, there exists a shift factor that can allow the isothermal frequency dependence to be translated laterally such that it overlaps the curve of a neighbouring temperature (Swallowe 1999;Ferry 1980). This procedure has been widely used to study the rate and temperature dependence of polymer properties (Bauwens-Crowet et al. 1969;Bauwens-Crowet 1973;Williamson et al. 2008). This allows a master curve of modulus as a function of frequency to be plotted at a given reference temperature. Figure 2(c) shows the isotherms shifted for a reference temperature of 25 • C and the resulting master curve. For more details on the master curve construction, refer to Appendix A. The inset shows the results in the frequencies of experimental relevance in the current study. One of the most common examples of an expression for the temperature dependence of shift factors is the Williams, Landel and Ferry (WLF) formula (Williams et al. 1955); however, a quadratic relationship is described showing the rate and temperature dependence of the modulus, (c) which are shifted laterally to create a master curve at a particular reference temperature, (d) using the shift factors as described (Colour figure online) here ( Fig. 2(d)) based on the empirically obtained shift factors used to construct the master curve. This relationship can be extended for different reference temperatures to inform the model of the temperature-rate dependent moduli, which is discussed later.

Characterisation of self-heating
At strain rates higher than 10 −1 s −1 , the compression of the specimen is no longer isothermal; the heat produced during plastic deformation does not conduct out of the specimen on the timescale of the experiment, and instead causes the temperature to rise, so-called adiabatic heating. One effect of this is observed in Fig. 1(a), the strain hardening is less at high strain rates due to thermal softening.
The reason for this phenomenon is that the timescale for heat diffusion out of the specimen is significant compared to the duration of the experiment itself. An analytical approximation to the rate above which this occurs may be found using Eqs. (1a)-(1c). Initially, the thermal diffusivity is obtained from the conductivity (k), density (ρ) and the heat capacity (C p ). Combining this with the characteristic length of the specimen (here 2.5 mm, half the specimen length), it is possible to find the characteristic time scale for diffusion out of the specimen. An approximate strain rate at which this timescale becomes significant is thus given by the reciprocal of the characteristic time scale.
Using the conductivity, density and heat capacity values for PVC from Matthews (1996), the lower and upper bounds of the diffusivity are 0.088 to 0.150 mm 2 s −1 . Using the characteristic length of the specimen, this gives the strain rate at which adiabatic effects start becoming significant to range approximately from 0.056 to 0.093 s −1 . Note that despite being calculated for PVC and not PPVC, the result is broadly the same as the 10 −1 s −1 observed experimentally. The exact value is not of importance as much as the order of magnitude, since the rate dependent material properties are often quoted in literature to depend on the logarithm of the strain rate.
Heat capacity is not constant with temperature, and as a result Differential Scanning Calorimetry (DSC) experiments were performed to obtain this temperature dependence. In order to separate complex transitions into the desired processes, modulated experiments were performed. Often, a glass transition in a conventional DSC may have similar characteristics to a melting transition due to an endothermic enthalpic relaxation. MDSC experiments allow the separation of the desired thermodynamic effects from the undesired kinetic effects. In such experiments, a sinusoidal temperature oscillation is superposed onto the conventional temperature ramp and the differences in heat flow required to maintain this temperature profile is monitored between the empty reference pan and the pan containing the specimen.
Parameters used for the MDSC experiments were an underlying ramp of 3 • C min −1 and a modulation of ±0.5 • C min −1 all over a temperature range from −90 • C to 200 • C. The results of these experiments are shown in Fig. 3(a). The results are split into reversing and non-reversing heat flow for the contribution from the instantaneous heating rate and the contribution that does not follow the modulation respectively. In finding the heat capacity, we are interested in the reversing heat flow since it relates to thermodynamically controlled changes. Two are provided, one on each side of the glass (α) transition (T g ), to describe the temperature dependence of the heat capacity. For the sub-T g fit, coefficient values are Once the temperature dependent heat capacity is known, the temperature rise due to this adiabatic self heating can be calculated using Eq. (2), where β is the Taylor-Quinney coefficient dictating the amount of mechanical work that is transferred to heating the specimen. Unless stated otherwise, it is assumed that β = 1 and all mechanical work is converted to heat. Figure 3(b) shows the theoretical temperature rises for the adiabatic compression experiments when the effect of the evolving heat capacity is included.

Hyperelasticity
In order to describe the large elasticity that many polymers exhibit, a linear elastic model like Hooke's law is often inadequate. To capture this type of non-linear elasticity, a hyperelastic model is needed. Many models have been proposed for this behaviour: the early invariant based continuum Mooney- Rivlin model (Mooney 1940;Rivlin 1948); the continuum stretch based Ogden model (Ogden 1972); and the macromolecular Langevin chain statistic based Arruda-Boyce model (Arruda and Boyce 1993;Arruda et al. 1995). Although these models are popular for capturing large stretch dominated elasticity, by themselves, their use is limited in describing plasticity or rate dependence. Suitable adaptations, such as the inclusion of rate dependence and thermal effects must be made in order to incorporate these features into the overall model.

Viscoelasticty
Despite the fact that the rate dependence of many materials has been known for some time (Tobolsky 1952;Ferry 1980), the modelling of this behaviour can still be quite challenging.
A common method to describe viscoelasticity has been to make use of linear combinations of springs and dashpots. The most basic of these is the Standard Linear Solid or Zener model (Zener and Siegel 1949), which uses either Maxwell or Kelvin-Voigt elements in parallel or series with a spring. This can be generalised with numerous such elements forming a Prony series (Fung 1965;Tschoegl 1989). However, this approach has many disadvantages, one of which is the fact that each decade of frequency requires an element in the Prony series to describe the viscoelastic behaviour appropriately, thereby adding numerous parameters requiring characterisation. A more satisfactory approach is shown below in which a more generalised fractional derivative model (FDM) (Koeller 1986;Torvik 1983, 1986) will be fitted to the results of the DMA experiment in order to describe the non-linear viscoelastic behaviour.

Plasticity
For large strain behaviour, plasticity must be included in the model. Rate dependent plasticity has been implemented previously with the Mulliken-Boyce (MB) model Jordan et al. 2008), and this is the model that forms the basis of this study. The model includes the effects of thermal softening through adiabatic heating by explicitly increasing the temperature based on the theoretical increases as calculated above ). However, the seamless integration of the effects of adiabatic heating into a model with predictive capabilities at high strain rates, entirely calibrated at low strain rates, has not been fully explored and is therefore the focus of this study.

Basics of the model
The Mulliken-Boyce model is based on the premise that since the molecular motions of the α and β processes are rate and temperature controlled, by using the results of a DMA experiment and its information on the rate dependence of α and β transitions, the modulus can be shifted accordingly. In the approach initially proposed by Mulliken, the frequency dependence of each transition is used to shift modulus-temperature data obtained from DMA experiments at low frequencies to reconstruct the expected modulus-temperature curves at high strain rates. Each transition is shifted independently to room temperature and then reconstructed to obtain the high rate behaviour at the temperature. The model then considers the contributions to the stress to be an additive sum of the α and β components and the hyperelastic limiting behaviour: σ = σ α + σ β + σ .
The α and β components are both represented by a spring and viscoplastic dashpot, which correspond to the intermolecular resistance for each of these processes. The hyperelasticity is represented by a non-linear, non-Gaussian Langevin spring; one that requires extra stress to stretch out an incrementally larger amount, and corresponds to the entropic resistance to chain alignment ( Fig. 4(a)).
The hyperelastic Langevin stress is calculated based on the network back-stress approach used in the Arruda-Boyce model (Arruda and Boyce 1993): where E 0 is the long-term rubbery modulus; √ N is the limiting chain extensibility based on random walk statistics; λ chain is the stretch on a chain given by √ I 1 /3 (I 1 is the first invariant where for an incompressible 1D case it reduces to U 2 + 2U −1 , where U = exp ε, the stretch based on strain); L −1 (x) is the inverse Langevin function: L −1 (x) = coth x − x −1 ; and the final term is based on the stress expression that would be obtained with the chain stretch described above and a 1D neo-Hookean free-energy expression: ψ = 3 2 (λ 2 − 1). To work out the stress in either α or β component, one can consider the stress in the linear spring portion. In this manner, the α and β stresses would be given by (i = α, β) where E i is the modulus and ε e is the elastic strain. Since the strain across each of the elements would be identical for compatibility reasons, assuming a constant strain rate, ε =ε e i +ε p i . Taking time derivatives of Eq. (4) and combining with the previous equation, it is possible to obtain the total stress rate asσ i = E i (ε e −ε p i ). Mulliken and Boyce used the shifting method described above to obtain E i ; in this paper the FDM will instead be used.
The shear strain rate in the viscoplastic dashpots can be given by the following equation: whereγ p 0 is proportional to the attempt frequency; Q is the activation energy; k is Boltzmann's constant; θ is the absolute temperature; τ is the shear stress in the dashpot; s is the athermal shear strength, an internal variable that models the post-yield strain softening behaviour; α p is the pressure coefficient dictating the pressure (p) dependence. In order to calculate the temperature, the following equation is used, where only plastic strain gives a temperature rise since stored elastic strain energy is recovered: The rate at which the athermal shear strength evolves is given bẏ where h is the softening slope and s ss is the preferred steady state value. The initial value is defined as where μ is the shear modulus; ν is the Poisson's ratio, assumed to be 0.5 for both elements in the model due to the incompressibility assumption.

Fractional derivative model
In order to define the modulus in Eq. (4) over a wide range of rates, conventionally a Prony series as described in Eq. (8) would be used, where E ∞ is the long term modulus, E i is the modulus for an individual element, and t i is the relaxation time for that element:  Fractional model The disadvantage with using this definition is that in order to get an accurate depiction of the modulus, a very large number of elements is needed. In this application, over the frequency ranges described in the master curve in Fig. 2(c), at least 25 elements would be required each with an unknown requiring parameterisation. Meanwhile, with an FDM as described in Eq. (9) (and Fig. 4(b)), the number of parameters can be reduced to a minimum of six for the same frequency range: Hence, the FDM is a preferable choice for representing the modulus over a wide frequency range. The new parameter β lies between zero and one, and dictates the slope of the transition in modulus, which is not a rapid change, as is the case with a single Maxwell element. Note the similar forms of Eqs. (8) and (9). The physical representation is that when all of the β terms are one, the model represents a Newtonian fluid and the FDM reduces to the Prony series, and when they are zero, the model represents a linear elastic body. Most materials exhibit a mixture of these behaviours, hence the FDM offers a more meaningful model to fit to this type of non-linear viscoelastic behaviour. In this application, both Prony series and FDM representations are applied to the DMA results for comparison, the parameters are shown in Table 1 and the fits in Fig. 5. Using the three element version of the FDM as in Fig. 4(b), the additive sum of three elements along with the long-term modulus gives the full response. The three elements used represent the alpha and beta transitions, and also a third transition representing the change in modulus between the first two.
The biggest advantage of the FDM representation is that it complements the Mulliken-Boyce modelling framework. It is for this reason that it is preferred to the conventional Prony series. A two component fractional derivative model would map directly to the two terms (α and β) used for viscoplasticity, making each FDM component physically meaningful. In order to maintain this advantage when the model was implemented in this study, the modulus associated with the third transition was split into two and allocated equally to the alpha and beta components (Fig. 4(c)).

Thermal relaxation of the modulus
The modulus is depicted to change temporally with the approximation described as follows. The storage modulus from the DMA experiments describes the variation with frequency. The frequency is approximately equal to the reciprocal of the time. For a given constant strain rateε and a maximum strain ε max , the test is conducted over a time t max = ε max /ε. This can then be related to the change in storage modulus over that time to give an approximation as to how the modulus changes with time.
In characterising the self-heating phenomenon, the temperature rise for each simulation was obtained using Eq. (6). Although the MB model includes the change in θ in Eq. (5) due to adiabatic effects, the changing modulus is not captured. Not only does the modulus relax with time, but there is a thermal relaxation, too. In this modelling framework, we propose that the thermal effects to the modulus can be captured by considering the shifting of the master curve to the temperature resulting from the adiabatic heating.
In order to obtain the modulus a certain temperature above the reference temperature, E (log f, T 0 + T ), there are two possible options as seen in Fig. 6(a). One is to fully reconstruct the master curve at the new temperature and read off the value of the modulus at the frequency of interest. This, however, is cumbersome when we require the modulus at each incremental temperature increase due to the adiabatic self-heating effect. The alternative method is to use the same master curve for the reference temperature, but read off the value at a frequency lower than that of interest, E (log f − log a T , T 0 ), where log a T is also a function of the temperature.
Since the temperature dependence of the shift factors is known (Fig. 2(d)), the additional shift required to compensate for the higher temperature and take the master curve from the reference temperature (which is the same as that of the uniaxial experiments for the room temperature experiments, T 0 = 25 • C) to the increased temperature is simply obtained by substituting T − T for T . At each timestep, this calculation is performed to obtain the temporally and thermally relaxed modulus to be used in the modelling framework as shown in Fig. 6(b).

Enhanced Fractional Thermo-Visco-Hyperelasto-Plastic (FTVHEP) model
In this section, the procedure will be described for the PVC in order to obtain the relevant parameters for the model.

Langevin component
For the hyperelastic Langevin stress component, there is only one parameter that needs identifying -the limiting chain extensibility √ N . This can be obtained by applying Eq.
(3) to fit to the lowest rate quasi-static compression experiment. This is as the large strain response is dominated by the hyperelastic component.
As seen in Fig. 7, the experimental results for the uniaxial compression conducted at 10 −3 s −1 is shifted downwards to isolate the large strain hyperelastic component. In this way, the fitting can made appropriately to obtain the value of √ N = 1.3.

α component
The α component relates to the glass transition and this element dominates at the low strain rates. For this reason, the parameters associated with the α component were obtained by fitting the σ α expression to the low rate experiments conducted from 10 −3 s −1 to 10 −1 s −1 .
The results of the fitting can be seen in Fig. 8. Both the yield stress and the large strain response are captured extremely well. There is, however, some discrepancy between the initial modulus as obtained from the experiments and that obtained from the FDM. This is likely to be a result of the difficulty in calculating an accurate modulus in compression experiments on small specimens. The smooth pre-yield behaviour is hard to capture with the MB model since for the α element, there is only one activation energy. In polymers, however, there is a spectrum of energies associated with the varying length chains in the microstructure. This feature of the experimental response can be captured in a future model iteration by considering a distribution of activation energies as is done in the Hasan-Boyce model (Hasan and Boyce 1995).
The parameters that are fitted in this step are the pre-exponential coefficient related to the attempt frequency (γ p 0,α ), the activation energy (Q α ), the steady state athermal shear strength (s ss,α ), the softening slope (h α ) and the pressure coefficient (α p,α ).

α and β components
Once the α components are set and fixed, they can be used whilst allowing the β components to fit to the low rate compression experiments at varying temperatures. Since the β compo-   (2014) Pre-exponential factor,γ 0 6.33 × 10 19 s −1 9.64 × 10 16 s −1 Curve fits, low rate data in  Activation energy, Q 3.28 × 10 −19 J 1 .53 × 10 −19 J Curve fits, low rate data in  Softening slope, h 75.3 MPa N/A Curve fits, low rate data in  Steady state athermal shear strength, S ss 41.0 MPa N/A Curve fits, low rate data in  nent to the mechanical response is only active at either low temperatures or high rates, the rate used was 10 −2 s −1 and the temperatures were varied from −60 • C to 0 • C. The results of the fitting are seen in Fig. 9. Once again, the yield stress is captured very well. The initial modulus varies for the different temperatures. Broadly speaking, the large strain response is well captured for all of the experiments, however could be better for the −60 • C and the −20 • C experimental fits. Further experiments at these temperatures could help to understand how much of this discrepancy is due to experimental variation.
The parameters that are fitted in this step are the pre-exponential coefficient related to the attempt frequency (γ p 0,β ), the activation energy (Q β ), and the pressure coefficient (α p,β ). Since the β component is not expected to soften, the steady state athermal shear strengh (s ss,β ) and the softening slope (h β ) do not exist in this model. A summary of the parameters required for the model, their values, and how they are obtained is provided in Table 2.

High rate predictions
The purpose of this modelling framework is to be able to predict the mechanical response of the PVC at high strain rates. None of the parameters in the modelling framework have been obtained from high rate experiments.
Presented in Fig. 10(a) are the predictions for the mechanical response of PVC at various strain rates up to the order of 10 3 s −1 and their comparisons to the previous study. Figure 10(b) shows the predictions and experimental comparisons for the PPVC obtained by following the same modelling procedure as outlined above.
For the PVC, the predicted high rate mechanical response matches the experimental results exceptionally well. This fact gives more confidence in the predictions for the intermediate strain rates, since those experiments are challenging to perform with a constant true strain rate. The incorporation of the adiabatic effects can also be seen where the large strain response is not dominated by hyperelasticity as it is for the lower strain rates. On the whole, the predictions for the PPVC are also very good and closely match experimental results from the previous study. There is some discrepancy at the intermediate strain rates; however, this could also be due to the challenge with performing experiments at these rates. Careful experiments can be done at these rates to interrogate this further.
To the best of our knowledge, the ability to predict the full mechanical response of PVC and PPVC at these strain rates is something that has not previously been achieved.

Peak stresses
In the literature, it is common to plot the yield or peak stress variation with the logarithm of the strain rate. For this reason, Fig. 11(a) shows the predicted peak stresses of PVC and the experimental comparisons, and Fig. 11(b) shows the same for the PPVC.
For the PVC, the peak stress predictions match very closely to the experimental results, as is expected since the whole stress-strain profile matched well for the PVC. For the PPVC, there is a consistent over-prediction for the peak stresses, but the magnitude of this is not that large and overall predicts the peak stresses well for this material too.

Effect of heating
In Sect. 2.4, the self heating at high strain rates was characterised and it was found that at rates above ≈ 10 −1 s −1 , self heating must be considered. To highlight the effect of self heating, Fig. 12 shows the difference in predictions for PVC with a model incorporating self heating through adiabatic effects, and one that is purely isothermal.

Simulating adiabatic experiments
In a previous study involving the PPVC (Kendall and Siviour 2013), an experimental method was proposed that would effectively simulate high rate experiments with ones performed at lower rates and lower temperatures, but with a gradually increasing temperature during the test. The selection of the low temperature was based on the temperature at which the initial modulus and yield stress behaved similarly to that of the high rate experiment. This technique is possible within the modelling framework presented in this paper too. Figure 13 shows the model high rate predictions (red dashed) and the experimental validation for a compression experiment conducted at 3300 s −1 (red solid). A quasi-static experiment conducted at 10 −2 s −1 and −40 • C (blue solid) shows the most similar modulus and yield response to the high rate one. The differences are that the yield is slightly lower and that the large strain response is obviously that which is expected for an isothermal experiment -one largely dominated by hyperelasticity and showing a gradual upturn in stress. In addressing these issues, it is simple to use the model here with the same strain rate, a decreased temperature of −55 • C, and the self-heating 'switched on'.
The decreased temperature increases the yield strength, however in order to match the large strain experimental results, a lower temperature than ideal was used. Note that the isothermal model prediction for −55 • C (black solid) matches quite closely to the isothermally run model simulation for 3300 s −1 (red solid). By using this 10 −2 s −1 simulation at −55 • C and turning the self-heating on, the response to the high rate experiment can be simulated (black dashed) as they are in Kendall and Siviour (2013). It may seem as if the simulated adiabatic response is better at capturing the high rate response of the PPVC, but it should be noted that the adiabatic simulation run at low temperature does not take into account the fact that the modulus for a model run at 3300 s −1 would be larger at all strain values considered than that for the 10 −2 s −1 simulation at −55 • C leading to higher calculated stresses. Furthermore, the viscoplastic effect would also influence the predicted stress for the high rate model prediction, leading to a higher response than the adiabatic simulation of the response despite the isothermal curves for both looking similar. Further work is needed to link the temperature-frequency dependence of the modulus and the yield stress.
The key benefit in using this modelling framework is that no further real experiments need to be conducted other than those required to parameterise the model in the first instance, and those needed for comparison at high rate.

Conclusions
In this paper, we have shown that it is possible to use a combination of various modelling techniques to predict the full stress-strain response of high rate compression experiments for both PVC and PPVC.
Starting with the MB model as a base, enhancements were made to allow this predictive capability. A relaxation modulus that not only shows its temporal decay via an approximation, but the drop associated with self heating of the specimen is included by considering the appropriate shifting of the master curve based on the TTS technique. By using a FDM representation for the modulus, the number of parameters required was greatly reduced when compared to the conventional Prony series representation.
After incorporating these changes, the enhanced model was parameterised using solely quasi-static experiments at varying low rates and varying temperatures. Predictions for the high rate mechanical response were then made using the model and validated against carefully conducted experiments from a previous study ) and shown to match very closely. In this paper, the pressure dependence is obtained through compression experiments solely, which is not ideal since the yield stress in compression is more than that in tension and so at the very least both sets of experiments should be used.
The effect of self heating in the model was also investigated to show that the temperature increase during the higher rate experiments cannot be neglected in the modelling process and must be included. Furthermore, it was shown how this modelling framework can assist in the simulation of adiabatic experiments. Upon calibrating the model at varying low strain rates at room temperature, and a low rate at varying low temperatures, this model is able to provide insight into various isothermal and adiabatic experiments.
In future, this novel modelling framework will be able to form the basis of predicting the high strain rate compression response for numerous materials where the behaviour is strongly rate and temperature dependent. The predictive capabilities of this modelling framework will therefore reduce the number of high rate experiments needing to be conducted as they would only be required for validation purposes and to interrogate particular rates of interest.

A.2 Master curve construction process
Instead of plotting the temperature dependence of the moduli for each of the three frequencies as in Fig. 14, it is possible to plot instead the frequency dependence at a given temperature. When plotted against the frequency on a log scale, each of these isotherms is linear as shown in Fig. 15(a). A horizontal shifting algorithm is then employed to construct the master curve as shown diagrammatically in Fig. 15(b). Other isotherms can be shifted to a reference temperature, T ref (here, 25 • C to correspond with the temperature of the compression experiments).
If T 1 < T ref , the isotherm is shifted in the positive (to the right) direction to overlap with the reference isotherm. Similarly, if T 2 > T ref , the isotherm is shifted in the negative direction. The amount of shift required is called the shift factor and is described as log [a T ].
Although there exists an ISO standard for time-temperature superposition (ISO 18437 2012), this was not used as the existing DMA data used only had three data points per temperature and lacked the detail prescribed . For this reason, an simplified alternative method for TTS was used, which has been previously published (Bauwens-Crowet 1973;Drodge et al. 2010;Yoon et al. 2014).
Each isotherm is represented by a first order polynomial fitted through the three data points: E = m log [f (Hz)] + c. In order to shift an isotherm to lay on top of the reference one, we consider the fact that the y-value (here the storage modulus) must be the same to achieve this. The shift factor is hence the difference in x-values (here the logarithm of the frequency) as described in Eq. (10): Once all of the individual isotherms (532 in this study) are shifted using this shifting algorithm, a master curve will be produced at the desired reference temperature as seen in Fig. 15(c). Although relatively short isotherms consisting of three data points were used to construct the master curve, the fact that a smooth master curve has been produced as a result of the shifting process provides confidence in the method employed.

A.3 Error mitigation
In attempting to understand the extent to which isotherms characterised by only three data points per temperature (PPT) can be used to accurately construct the master curve, the following modelling based investigation is undertaken. Start with the DMA based experimental master curve modelled with the Prony series (Eq. (8), Table 1). By making use of the temperature dependent shift factors (Fig. 2(d)), a series of smooth isotherms can be defined. Each isotherm shows the dependence of the storage modulus on the frequency. By multiplying the frequencies of interest with the shift factor, a T at a particular temperature from Fig. 2(d), reduced frequencies are obtained. Substituting these reduced frequencies into Eq. (8) gives the values of the storage moduli as required for the reconstructed isotherm at that temperature. The temperature range was kept Fig. 16 Comparison of master curves constructed with either three PPT or ten PPT. For completeness, the original master curve is also provided (Colour figure online) the same as the original DMA experiment; from −100 • C to 120 • C in 2.5 • C steps. For the isotherms characterised by three PPT, the same frequencies as the DMA experiments were chosen; 1, 10 and 100 Hz. For the series of smoother isotherms, the frequencies described in the ISO standard (ISO 18437 2012) were used; 0.1, 0.2, 0.3, 0.5, 1, 2, 3, 5, 10 and 20 Hz. These smooth isotherms are each characterised by ten PPT as opposed to the experimental isotherms, which were only three PPT.
The shifting algorithm as defined above is then employed to construct the master curve. For the sake of completeness, reconstructed isotherms characterised by only three PPT are also used to produce a master curve, in order to compare directly with the original master curve presented in this paper.
By comparing the master curves produced with isotherms characterised by either three or ten PPT, it is possible to determine whether a master curve can be produced with short segments. In Fig. 16, these master curves are shown for comparison. As is evident, there is no significant inaccuracy resulting from using only three points to characterise each isotherm as the master curves effectively overlap each other perfectly. There is some discrepancy between the master curves produced using the Prony series and temperature dependent shift factors versus the original master curve using empirical shift factors. This is due to slight inaccuracies in using the quadratic shift factor definition to model the empirical shift factors measured, especially at 25 • C, as can be seen in Fig. 2(d). This results in a 0.8 decade mismatch between the empirical master curve and the two Prony derived reconstructions.
In order to mitigate the additive nature of the error in constructing the master curve, an additional shift factor has been introduced. The purpose of this 'zero-shift' factor is to take the whole master curve once it has been constructed and translate it such that the value of the shift factor at the reference temperature is zero. Further confidence in the accuracy of the constructed master curve is provided by the fact that they can be used in the modelling strategy presented in this paper to accurately predict the mechanical response of (P)PVC. If there were significant errors, then the viscoelastic definition at a given strain rate would not be accurate and the model would not have been able to predict the behaviour with such fidelity. Of course, the overall modelling approach proposed in this paper can be used with any suitable TTS methodology.

Appendix B: Modelling framework overview
In this appendix, the modelling framework used in this paper is presented as a flowchart (Fig. 17) and detailed with step-by-step information to guide the user.

B.1 Quasi-static and thermomechanical experiments
The first step is to conduct the various quasi-static and thermomechanical experiments required in the modelling framework. DMA experiments are conducted over a range of temperatures covering the transitions of interest, and at discrete frequencies. This combination allows the temperature and frequency (or rate) dependent properties of the material to be obtained. Quasi-static compression experiments are also performed both at an ambient temperature and a variety of low strain rates, and at sub-ambient temperatures at a fixed low strain rate. DSC experiments allow for the temperature dependent heat capacity to be deduced, alongside a corroboration of the transition temperatures for the material.

B.2 Obtaining E(t)
Starting with the results of the DMA experiments, one constructs a master curve to obtain the frequency dependence of the storage modulus, E (f ), for a fixed reference temperature (usually set to room temperature). Once the master curve has been constructed, a function can be fit to it in order to describe its form. Traditionally, a Prony series of Maxwell elements has been used as a fitting function. This requires an additional element per decade of frequency covered by the master curve. By using an FDM representation, the number of fitting parameters is drastically reduced from two per element per decade of frequency for the Prony series (plus one to describe the long-term modulus) to three per element (plus one). The key difference is that with the FDM, one element can accurately describe the transition behaviour across many decades of frequency. In PVC and PPVC, there are two transitions and so a minimum of two elements of the FDM are required. Since there is a complex interaction between the two transitions, an extra element is introduced to describe the behaviour in between the two distinct transitions. A big advantage of using a two component FDM over the Prony series is that the components map directly to the two terms used for viscoplasticity in this modelling framework. In order to maintain this advantage when the model was implemented in this study, this third element's contribution to the modulus is split evenly between the two elements describing the α and β transitions.
The task now is to convert from E (f ) to the relaxation modulus, E(t). When using the FDM, analytical conversions only exist for special cases of the fractional power, and so for arbitrary fractional derivatives numerical interconversion methods may be used (Schapery and Park 1999;Park 2001). A further approximation to this numerical conversion can simply state E(t) ≈ E (f ), where f = 1/t, this was deemed sufficient for obtaining the temporally relaxing modulus for the purposes of this modelling framework.

B.3 Obtaining E(t, T )
From the DSC experiments, functions can be fit to the results of the heat capacity versus temperature graphs to describe this dependence and enable its implementation into the modelling framework. Here two fits were required, one to describe the dependence before the transition temperature, and one to describe the dependence post-transition.
For high strain rate simulations, since the temperature rise due to adiabatic effects and its effect on the modulus must be considered, information from the DMA and DSC experiments have to be combined to obtain E(t, T ). Since the FDM expression can lead to E(t), it forms the foundation. Initially, the modulus for the simulation is taken as the room temperature value, which is also what the FDM uses as the reference temperature. Until the plastic strain is non-zero, this remains true. Once the plastic strain increases, Eq. (6) can be used calculate the incremental temperature rise at each integration step. By using the same master curve for the reference temperature, one can read off the value at a frequency lower than that of interest, E (log f − log a T , T 0 ), where log a T is also a function of the temperature. Since the temperature dependence of the shift factors is known, the additional shift required to compensate for the higher temperature and take the master curve from the reference temperature to the increased temperature is simply obtained by substituting T − T for T . At each timestep, this calculation is performed to obtain the temporally and thermally relaxed modulus to be used in the modelling framework.

B.4 Obtaining Langevin parameters
In order to obtain the Langevin (L) parameters, the results to the lowest strain rate uniaxial compression experiments are required. Since the mechanical response of the lowest strain rate ambient compression is dominated by the large-strain hyperelastic response, by fitting to it with the Langevin stress expression, it is possible to obtain the two parameters: E 0 , the long-term rubbery modulus and √ N , the limiting chain extensibility based on random walk statistics.

B.5 Obtaining α parameters
The α component dominates the response at low strain rates, and for this reason its parameters are obtained by fitting to ambient uniaxial compression experiments at a variety of low strain rates. The modulus used in these fittings are obtained directly from the FDM α element, which is supplemented by half of the transitory third element. Since these experiments are at low strain rates, there is no need to worry about the thermal softening of the modulus. Fits are made simultaneously to all relevant curves using a least squares minimisation scheme. Once obtained from these fits, the α parameters are fixed.

B.6 Obtaining β parameters
By noting the time-temperature equivalence in that the observed mechanical response at low temperatures is akin to the response at high strain rates and an ambient temperature, the remaining model parameters of the β viscoplastic element are fitted to low rate uniaxial compression experiments conducted at low temperatures. This is possible since the molecular motions associated with the β transition are inhibited by a low temperature, and also by the lack of time duration for molecular motions during a high strain rate experiment.
Here, the modulus includes the full FDM description including contributions from all of the elements. For implementation into this modelling framework, the modulus is split into the α component comprised of the FDM α element and half of the transitory third element, and the β component comprising the FDM β element and the other half of the transitory element. As with the α parameters, fits are made simultaneously to all relevant curves and once obtained from these fits, the β parameters are also fixed.

B.7 Predicting and validating high strain rate response
The combination of the previously obtained α, β, L parameters and temporally and thermally relaxed modulus, E(t, T ) can be used in combination in this modelling framework to predict the high strain rate response of the PVC and PPVC. Until now, the thermal relaxation of the modulus has not been required, but in predicting the high strain rate response, this cannot be neglected. In the simulation process, at each integration step, the temperature increase is calculated and used to obtain the new modulus value to be used in the next integration step.
The model is run with the strain rate now set to the value for which the response is desired to be predicted, whilst retaining the fixed model parameters as obtained above. These simulations can then be validated against carefully conducted high strain rate experiments at the same strain rates. In this paper, these experiments were conducted by Kendall for a previous study .