# Advantages of a 3-parameter reduced constitutive model for the measurement of polymers elastic modulus using tensile tests

Arnaud Blaise, Stéphane André, Patrick Delobelle, Yves Meshaka, Christian Cunat

<figure><img src="image_1.png" /></figure>

## To cite this version:

Arnaud Blaise, Stéphane André, Patrick Delobelle, Yves Meshaka, Christian Cunat. Advantages of a 3-parameter reduced constitutive model for the measurement of polymers elastic modulus using tensile tests. Mechanics of Time-Dependent Materials, 2016, 20 (4), pp.553-577. 10.1007/s11043-016-9312-1. hal-04180369

HAL Id: hal-04180369

https://hal.science/hal-04180369

Submitted on 11 Aug 2023

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Advantages of a 3-parameter Reduced Constitutive Model for
the Measurement of Polymers Elastic Modulus using Tensile
Tests.

A. Blaise<sup>(1,2)</sup> . S. André<sup>(1)</sup> . P. Delobelle<sup>(3)</sup> . Y. Meshaka<sup>(4)</sup> , C. Cunat<sup>(1)</sup>

(1) LEMTA – CNRS UMR 7563, 2 avenue de la Forêt de Haye, 54504, Vandoeuvre-Lès-
Nancy, France

(2) Laboratoire 3SR – CNRS UMR 5521, 1270 rue de la piscine, 38400 Saint martin d’Hères, France

(3) Femto-ST, CNRS UMR 6174, 24 chemin de l'Epitaphe, 25000, Besançon, France

(4) Institut Jean Lamour – CNRS UMR 7198, Parc de Saurupt, 54042 Nancy Cedex, France

E-mail: stephane.andre@univ-lorraine.fr

**Abstract:** Exact measurements of the rheological parameters of time-dependent materials are crucial to improve our understanding of their intimate relation to the internal bulk microstructure. Concerning solid polymers and the apparently simple determination of Young's modulus in tensile tests, international standards rely on basic protocols that are known to lead to erroneous values. This paper describes an approach allowing a correct measurement of the instantaneous elastic modulus of polymers, by a tensile test. It is based on the use of an appropriate reduced model to describe the behavior of the material up to great strains, together with well-established principles of parameter estimation in engineering science. These principles are objective tools that are used to determine which parameters of a model can be correctly identified according to the informational content of a given data set. The assessment of the methodology and of the measurements is accomplished by comparing the results with those obtained from two other physical experiments, probing the material response at small temporal and length scales, namely: ultrasound measurements with excitation at 5MHz and modulated Nanoindentation tests over a few nanometers of amplitude.

**Keywords** inverse identification, parameter estimation, elastic modulus, tensile test, semi-crystalline polymer.
# 1 Introduction

This paper originates from the observation of a significant gap between recommended practices in the field of rheological parameter determination, which rely on international standards, and the efforts made by the scientific community to improve measurements. During recent decades, the scientific field of inverse methods, which deals precisely with the mathematical basis of parameter estimation through models (Model-Based Metrology, or MBM), emerged and grew intensively. MBM should nowadays be a common framework when considering the metrology of physical parameters through pertinent modeling of experiments. It is necessary for a real enhancement of the quality in measured properties of materials, which will in turn allow a deeper investigation of the link that can be made with the microstructural organization of materials and help to further develop the models. To illustrate this point of view, an example is given in this paper concerning the determination of the elastic modulus of polymers using the uniaxial tensile test. With respect to international standards (ISO 527-1 and ASTM D 638), it depends upon the characterization of a linear elastic regime, i.e. proportionality between stress and strain, at very small strains (short times). Although this method is efficient and physically well-founded for metallic materials, it is less satisfactory for polymers. Indeed, these materials quickly manifest a viscoelastic regime, which takes place in the beginning of the tensile test, and the response does not show any obvious linear part. A mathematical inverse approach of parameter estimation is required that may explain why this question, also basic and apparently simple, was not really considered earlier in detail and was consequently not well-documented. This latter point consists of the use of an optimization procedure based on uniaxial tensile test data and relying upon an adapted rheological constitutive model. The inverse problem is merely formulated as a least square optimization problem. It must be verified that the assessed parameters are not correlated and that their sensitivity is high enough to ensure a good quality of identification.

The aim of this paper is to present such an approach in view of the measurement of the elastic modulus of materials like HDPE (High Density PolyEthylene) one of the Semi-Crystalline Polymers (SCP). The paper starts with detailed information about the material we investigated and about the metrology we developed to produce curves of true strain-true stress relationships. In a second step, a rheological constitutive model is presented in its most reduced version. In view of a parameter estimation problem, a reduced model means a model that meets the identifiability criteria i.e. a good balance between the number of identifiable parameters and the relevance of the model to the experiments (parsimony principle). This model was shown to be able to describe different behavioral regimes (viscoelasticity, viscoplasticity, material hardening) with only 3 parameters. In a third part of the manuscript, the identification procedure and the sensitivity analysis of the estimated parameters are described. They rely on mathematical tools that have now spread throughout the engineering community (see the books of Aster et al., 2013; Beck, 1977; Walter and Pronzato, 1997 to cite a few references), but are not so frequently used in the field of material mechanics, where models are generally not questioned from the angle of the identifiability or measurability of their parameters. Essential works that address this question are generally based on Finite Element simulations at the structure scale (Ponthot, 2006; Cooreman, 2007) and make the inverse methods work on various observable features that can be local such as strain fields, or global like force, radius geometry, total elongation etc.. At the end, the results of the identification will be discussed and compared to identification as measured by two other scientific techniques. One of these is based on ultrasonic measurements (Piché, 1984; Legros, 1999). The other is based
on Nanoindentation tests (Oliver, 1992; Qasmi, 2006; Le Rouzic, 2009). Both prior approaches provide direct measurements of the Young modulus, meaning that no inverse approach is used.

# 2 Material and methods

## 2.1 Material

The material tested in this work is a “500 Natural” grade HDPE produced by Röchling Engineering Plastics KG. Two different products (A and B) were manufactured under the same reference at an interval of six years and were supplied as sheets (extrusion process). Information from the supplier indicates that the molecular weight and density are respectively 500,000 g/mol and 0.935 g/cm³. Differential scanning calorimetry gave a crystallinity index of 68 wt% for product A and 66 wt% for product B. The structural characteristics of this material in a non-deformed state and at room temperature were found to be a long period of the stack of 26.8 nm, a crystalline phase thickness of 18.5 nm and an amorphous layer thickness of 8.2 nm (from SAXS synchrotron measurements reported in Farge, 2013).

A particular geometry of dog-bone shaped specimens was used for the tensile tests. A narrow material “gauge” cube of 6×6×6 mm³ was deliberately arranged between the specimen shoulders so that necking was always initiated at the center of the specimen. A-type specimens were cut from a 6 mm-thick sheet of polymer along the extrusion direction “$A_{//}$” (reference samples throughout this paper) and perpendicularly to it “$A_{⊥}$”. Some specimens of the $A_{//}$ type were machine-finished symmetrically to produce 3 mm thick core samples referred to as “$A_c$”.B samples were cut along the extrusion direction only, from a 4 mm sheet produced six years earlier. Even if specimens A and B are supposedly of the same chemical material they exhibit a slightly different mechanical behavior, leading to different values of the parameters identified by the model. Experiments performed on B samples at delivery and 6 years later gave exactly the same results. Thus ageing or the time interval between the dates of their production cannot be considered to be responsible for the difference between samples A and B. Instead it is more likely to be due to a different manufacturing process. The microstructure and mechanical behavior of a polymer are clearly dependent on the cooling process which rules the development of a crystalline phase and its distribution and organization in the bulk among the amorphous phase. This material has been studied for more than 15 years and the interested reader can refer to (Farge, 2013a; 2013b; Blaise, 2011; 2010) to gain insight into all its microstructural aspects that are not covered in the present study.

The various types of specimen ($A_{//}$, $A_{⊥}$, $A_c$, B) for the same material will be very useful in view of the discussion about the results obtained through tensile tests and the validation of the methodology we used.
## 2.2 Video-controlled tensile tests and experimental tensile curves

All mechanical tests were performed on a servo-hydraulic MTS 810 load frame with a Flextest SE electronic controller. A video-extensometer (VidéoTraction®, G'sell et al., 2002) gives access to the real-time measurement of local true strain in the central part of the specimen (Fig. 1). Seven dot markers were printed on the front face of the sample prior to deformation. The markers are black, nearly round and have a diameter of about 0.4 mm. A CCD camera mounted on a telescopic drive records images during the test and follows the elementary frame during deformation. Local transversal and longitudinal deformations are measured at the barycenter of two neighboring dots, according to Hencky's definition of the true strain:

$$ \varepsilon(t) = \ln\left(\frac{l(t)}{l_0}\right) \qquad (1) $$

Longitudinal strains are interpolated through Lagrangian polynomials to obtain exactly the longitudinal strain value $\varepsilon_{11}$ in the section where the transversal strains are measured (averaged FCG line in Fig. 1a). An important feature of the system is that this latter measurement is used in real-time to control the machine's servo-valve through a feedback loop. Any desired input path for the true longitudinal strain $\varepsilon_{11}$ can be imposed on the gauge volume. A standard deviation of nearly $10^{-4}$ is commonly achieved for the noise that corrupts this strain signal (it increases logically as the marker deformation proceeds). VidéoTraction measurements were assessed with a 3D-DIC system (ARAMIS 3D Digital Image Stereo-Correlation from GOM Instruments), relying on the correlation of successive images because of a black and white pattern sprayed on the reverse-side of the specimen (Fig. 1b). For example, the longitudinal strain as measured by VidéoTraction has been proved to remain within an error of less than 4%, which accounts for the 2D plane measurement bias, since 3D-DIC takes into account the out-of plane displacement when necking occurs.

In addition to the longitudinal strain, the applied force $F$ is measured directly by a 5 kN load cell. These two measurements enable the axial true stress – true strain behavior of the material to be obtained. True stress takes into account the reduction of the cross-sectional area $S(t) < S_0$ undergone by the sample while it is stretched according to:

$$ \sigma_{11}(\varepsilon_{22}, \varepsilon_{33}) = \frac{F}{S} = \frac{F}{S_0}\exp(-\varepsilon_{22} - \varepsilon_{33}) \qquad (2) $$
<figure><img src="image_3.png" /><figcaption>Figure 1: a) Principle of the true strain measurements as put into action in the VidéoTraction system. b) Pictures of the specimen in non-deformed (left) and highly-deformed (right) states.</figcaption></figure>

Because the measurements of the transversal strains $\varepsilon_{22}$ and/or $\varepsilon_{33}$ are not very reliable with this system, or in the case where measurement of the transversal strain is simply not available, a constant volume strain hypothesis can be used. This hypothesis is nearly met for our HDPE samples, as proved with 3D-DIC measurements of the strain field on both lateral sides of a specimen. Fig. 2 shows that even if the transverse strain measurements exhibit some lack of reproducibility and precision, especially at high strains, the curves are nearly linear with approximately a -0.5 slope. Hence, we will consider the case where the Poisson coefficient $\nu$ is set to 0.5.

This criterion is summarized as follows:

$$ \varepsilon_{22}^* = \varepsilon_{33}^* = -\nu \varepsilon_{11} = -0.5 \varepsilon_{11} \quad (3) $$

A star superscript denotes the calculated transversal strain obtained under a constant volume strain assumption. This leads to a null volume strain:

$$ \varepsilon_v = \varepsilon_{11} + \varepsilon_{22}^* + \varepsilon_{33}^* = 0 \quad (4) $$

and to a true stress calculated with:

$$ \sigma_{11}^* = \frac{F}{S_0} \exp(-\varepsilon_{22}^* - \varepsilon_{33}^*) = \frac{F}{S_0} \exp(\varepsilon_{11}) \quad (5) $$

where the observed $\sigma_{11}^*$ is obtained directly from the measurements of $F$ and $\varepsilon_{11}$. In Fig. 2 the curves obtained for $\sigma_{11}^*$ are plotted for three repeated experiments and the curves obtained for $\sigma_{11}$ (equation 2), where the stress is calculated using the measurement of $\varepsilon_{22}$ and the transversal isotropy assumption $\varepsilon_{22} = \varepsilon_{33}$, are shown. The limits of this assumption were clearly established in a recent paper (Farge et al.,
2013). For the three distinct experiments shown here, it is clear that the $\sigma_{11}(\varepsilon_{11})$ curves present some discrepancies due to uncertainties in the transversal strain measurements, whereas the $\sigma^*(\varepsilon_{11})$ curves are superimposed. As both options introduce some bias, we chose to work with $\sigma^*(\varepsilon_{11})$ which ensures a high reproducibility of experiments. In the approach developed in section 4 of this paper, the analysis of post-identification residuals allows the extent of the overall metrological bias to be quantified.

Fig. 3 represents true stress $\sigma^*$ - true strain $\varepsilon_{11}$ curves for five different strain rates $\dot{\varepsilon}_{11}$ (covering a 3 decades range). Because only the longitudinal behavior in direction 1 is under investigation, subscript (11) will be omitted from now on in the notations of $\varepsilon$, $\sigma$ or other mechanical variables.

<figure><img src="image_4.png" /><figcaption>Figure 2: True stress $\sigma$ or $\sigma^*$ (left axis - the latter in solid black) and transversal true strain $\varepsilon_{22}$ (right axis), plotted versus longitudinal true strain $\varepsilon_{11}$ for 3 repeated experiments –</figcaption></figure>

Specimen $A_{//}$ at $\dot{\varepsilon}_{11} = 0.005 \text{ s}^{-1}$.
<figure><img src="image_5.png" /><figcaption>Figure 3: True stress σ* versus true strain ε₁₁ curves for different strain rates (Specimen A//).</figcaption></figure>

According to information that is well-known for SCP's, a plot is shown in Fig. 4 of the true stress as a function of the variable strain εHT defined as:

$$ \epsilon_{HT} = \lambda^2 - \lambda^{-1} = \exp(2\epsilon) - \exp(-\epsilon) \quad (6) $$

with $\lambda = \exp(\epsilon) = l/l_0$, the extension ratio.

This “Gaussian” strain variable arises naturally from microscopic modeling of the entropic elasticity of a molecular network (Treloar, 1975; Haward, 1993). It favors the observation of the mechanical regime at high strains which in SCP's generally appear to follow an almost linear relation

$$ \sigma^{hard} = G \epsilon_{HT} \quad (7) $$

where G stands for the rubbery (or hardening or hyper-elastic) modulus (in MPa) and the superscript *hard* is for the *hardening* stress, or stress at high strains. As can be seen in Fig. 4, this is so for our HDPE material.
<figure><img src="image_6.png" /><figcaption>Figure 4: True stress σ* versus strain variable ε<sub>HT</sub> (specimen A// - Different strain rates).</figcaption></figure>

Various domains can be clearly identified on the true stress σ* - true strain ε curves of HDPE (Fig. 3).
At very low strains, the mechanical behavior is viscoelastic, i.e. there is no existence of a pure elastic regime.
The yield point identifies the initiation of a necking phenomenon (Fig. 1b), which takes place for a strain
range of 0.05 < ε<sub>yield</sub> < 0.1, depending on the applied strain rate. At this point, the force reaches a
maximum in agreement with the first Considere condition (Hiss et al., 1999). After the yield point, the
mechanical response reveals the presence of a pseudo-plateau corresponding to polymer softening followed
by a plastic or flow regime. The point where the strain hardening coefficient is greater than 1 (Gaucher-Miri,
1996) occurs at a value of about 1.2. This corresponds more or less to the second Considere condition
(Haward, 2007) and marks the onset of hardening. This strain level corresponds to the propagation of necking
over the whole specimen. The hardening stage initiates in the viscoplastic flow regime and is followed in this
study up to strains of about ε ≈ 1.9. More details on the kinematic behavior of HDPE necking and its
underlying microstructural aspects are discussed in Ye (Ye et al., 2015) and will not be considered in this
study which is methodologically oriented.

3 Behavioral model for HDPE (an SCP)

3.1 General Formulation
To characterize accurately and correctly the HDPE mechanical behavior, a necessary condition is to develop a general constitutive model that can capture the physics involved. The model used in this study relies on Thermodynamics of Irreversible Processes (TIP.) (Callen, 1985; Kuiken, 1994). In irreversible thermodynamics, State laws are classically decomposed into an instantaneous, or unrelaxed, component and a delayed one. In the following considerations, superscripts *u* and *d* refer respectively to the *unrelaxed* and *delayed* components. This decomposition is classical in other thermodynamic approaches that make use of internal variables (e.g. Maugin and Muschik, 1994). In a tensile test where the strain is imposed, its dual thermodynamic variable is the stress, and this corresponds to the system response. The delayed response $\sigma^d$ is described using a modal approach (j subscripts in the equations will denote the jth mode). Meixner (1949) proved the existence of a modal basis to describe dissipative mechanisms in thermodynamics. In short, this means that the matrix that would represent the effect of a relaxation kernel in a viscoelastic integral formulation in algebraic form can always be diagonalized in an appropriate basis. As a result, any vector representing variables that govern internal processes can be transformed into vector modes. This classical way of thinking in mechanics leads to a simple description of relaxation through a first order kinetic model. In this model, the stress ‘fluctuation’ or difference between the modal true stress $\sigma_j$ and the modal relaxed stress $\sigma_j^r$ regresses according to the modal time $\tau_j$, which represents the relaxation time governing the kinetics of each dissipation mode j. The overall true stress $\sigma$ is the sum of all the modal components, $\sigma_j$. All this information is summarized in the following rate-type constitutive equation, which is given for a single component of the stress-strain tensor in its principal framework (Cunat, 1991; 2001a; André et al., 2012):

$$ \dot{\sigma} = \dot{\sigma}^u + \dot{\sigma}^d = \sum_{j=1}^{N} \dot{\sigma}_j = \sum_{j=1}^{N} \left( E_j^u \dot{\epsilon} - \frac{\sigma_j - \sigma_j^r}{\tau_j} \right) \quad (8) $$

$E_j^u$ corresponds to the modal *unrelaxed* modulus. $E_j^u$ and $\sigma_j^r$ can be defined by: $E_j^u = p_j^0 E^u$ and $\sigma_j^r = p_j^0 \sigma^r$ where $p_j^0$ is a coefficient that weights the modal component j with respect to the overall quantity. $\sigma^r$ refers to the stress in the *relaxed* state. This latter corresponds to the steady-state regime of internal mechanisms, where the non-equilibrium forces do not evolve with time. This thermodynamic state is strictly defined by $\dot{A} = 0$ and is named the "iso-affin" state according to the recommendation of Prigogine (1946). A is the affinity state variable introduced by De Donder (1936) to take into account chemical reactions. It is used by mechanical or material scientists in solid rheology to describe internal reorganizations that take place in matter under deformation processes (Nowick and Berry, 1972; Kuiken, 1994; Gutman, 1998). $E^u$ stands for the common elastic (Young’s) modulus and hence the modal weights must fulfill the normalized condition:
$$ \sum_{j=1}^{N} p_{j}^{0} = 1 \qquad (9) $$

From equation 8, we can infer that for a tensile test at a very low strain rate or when all relaxation times are depleted i.e. viscoelastic effects are finished and the material enters the hardening regime, the true stress is expected to behave as does the relaxed one (the non-equilibrium forces are constant). According to this, the relaxed state can be seen as a pseudo-equilibrium that is still dynamic. It must not be mistaken for the equilibrium state, which is rigorously defined in TIP by a state of zero affinity $A=0$ obtained when all external forces are stopped. As an example, this distinction is clearly made in a same type of modeling using the viscoplasticity theory based on overstress (VBO) approach of Krempl, 2001, although the same stress notation ($g$) and name ("equilibrium") are used for both. $g$ can be considered as the equilibrium stress "when all time rates are zero" and must be equal to 0 in the absence of applied stress. At the same time the equilibrium stress has to evolve under the application of an external force. It is then said to correspond to the stress "that must be overcome to generate inelastic deformation" and corresponds to the variable $\sigma^r$ in the present model.

The constitutive equation (8) is generally referred to as the DNLR approach (Distribution of Non-Linear Relaxations) and has been used in the past to reproduce different experiments (various polymers and loading conditions, Rahouadj et al., 2003; Mrabet et al., 2005). It was recently applied to analyze thermomechanical experiments with strain and infrared imaging (André et al., 2012). It offers two entry points for modeling: the spectrum of relaxation times and the relaxed state.

## 3.2 Spectrum of relaxation times

The full linear spectrum of relaxation times is defined by the following equation

$$ \tau_j = \tau_{\max} \cdot 10^{-\left(\frac{N-j}{N-1}\right)d} \qquad (10) $$

which requires three values to be assigned to parameters: $\tau_{\max}, d, N$. In practice and inherent to the modal thermodynamic approach one has to consider vanishingly small relaxation times. As equation 10 distributes relaxation times along a logarithmic scale ($\tau_{\max} = \tau_{\min} \cdot 10^d$) in the interval $[\tau_{\min}, \tau_{\max}]$, the number of decades $d$ has to be as high as possible. In practice, it is set to the minimum value that preserves the solution from being dependent on, or sensitive to, larger values of this parameter. In other words, the solutions of the Estimation Problem will not depend on its value once it is set to a maximum value. A number $d=6$ decades is generally required. This is so for the second (hyper)parameter of the model $N$, the number of dissipative modes or the number of Ordinary Differential Equations (ODE's) considered in the system of equation (8). This behaves for example like the number of nodes of a finite element calculation, the solution should not depend on its value. $N=50$ dissipative modes for example, build a quasi-continuous spectrum of relaxation times.
As a result, only one parameter is needed to describe this spectrum in linear viscoelasticity, the longest relaxation time $\tau_{max}$. This has the advantage of being closer to the physical perception of an experimenter. It is named $\tau_{max}^T$ where the “T” superscript denotes the Tensile phase of the test. In the case where a non-linear spectrum is needed, a multiplying (shift-)factor ($a(T,t,...)$) is classically used to drive the $\tau_{max}$ value (Facio-Toussaint et al., 2001; Mrabet et al., 2005).

Following the postulate of an equipartition of the entropy created that accompanies regression of the fluctuation around an equilibrium state, a 'universal' spectrum can be defined that links each dissipative modal weight $p_j^0$ to its corresponding relaxation time $\tau_j$ (Cunat, 1991; 2001a; 2001b):

$$P_j^0 = \frac{\sqrt{\tau_j}}{\sum_{j=1}^{N} \sqrt{\tau_j}} \qquad (11)$$

### 3.3 Relaxed state model adopted for HDPE under tensile tests

The second entry point for modeling lies in the model used for the relaxed stress $\sigma^r$, supposed to describe the plastic-hardening phase. The most basic form assumes the linear form $\sigma^r = E^r \varepsilon$, but this can only produce linear hardening in a plastic regime. As elaborated in section 2, a constitutive law for HDPE must comply with a more complex continuous modeling, from the viscoplastic to hardening behaviors. Hence, different modeling must be considered for the relaxed state. A SCP in a highly deformed state is generally considered as a rubber-type tridimensional network with entanglements. A model for elastomer-hardening first considered by Wang and Guth (1952) and modified by Arruda and Boyce (1993) put forward cells with eight sub-chains per node which is close to a random distribution of the chains in a real material and does not favor any spatial direction. Starting from the works of Haward & Thackray (Haward, 1993) and Treolar (1975), Arruda and Boyce (1993) showed that during hardening, the “back stress”, which is understood as the stress necessary to overcome the anisotropic resistance to chain alignment, can be described by the following formula:

$$\sigma^{\text{backstress}} = \frac{N_0 k_B T}{3 \lambda_c} \sqrt{n} L^{-1} \left( \frac{\lambda_c}{\sqrt{n}} \right) \varepsilon_{\text{HT}} \qquad (12)$$

with $\lambda_c = \sqrt{(\exp(2\varepsilon) + 2\exp(\varepsilon))/3}$ and where $N_0$ stands for the density of chains per unit volume, $n$ is the number of segments per chain, which controls the behavior at large strains until the maximum extensibility of the network is reached, $k_B$ is the Boltzmann constant and $T$ the temperature. $L^{-1}$ represents the inverse Langevin function.

According to Krempl (2001), the term “back stress” is often used in the literature but does not express clearly the subtle differences that exist between an equilibrium stress and overstress, called a relaxed state in our approach, already discussed in the "General Formulation" section (chapter 3.1). Furthermore, Negahban
(2006) considers that the “back stress” concept is responsible for the delayed, or inelastic, response. In the framework of Irreversible Thermodynamics based on the concept of Affinity, the distinction between true, relaxed, and equilibrium states is self-contained.

Returning to the experimental results presented above (Fig. 4), it is clear that the true stress plotted as function of the Haward-Thackray strain variable $\varepsilon_{HT}$ exhibits an almost linear relation, irrespective of the strain rate. This fact has already been highlighted in other works (Van Melick et al., 2003). In our paper and in view of proposing the most reduced model possible in terms of the number of parameters, the relaxed state of our model will be considered to follow the simplified version of equation 12 given below:

$$ \sigma^r = G \varepsilon_{HT} \qquad (13) $$

Only the single parameter G (called the hardening modulus) is expected to be correctly measured from this relaxed state description.

Finally, the reduced model that will be considered in this study contains only three parameters: $E^u$, G and $\tau_{max}^T$ and is now fully described by equation 14:

$$ \dot{\sigma} = \sum_{j=1}^{N} \dot{\sigma}_j = \sum_{j=1}^{N} \left( p_j^0 E^u \dot{\varepsilon} - \frac{\sigma_j - p_j^0 G \varepsilon_{HT}}{\tau_j} \right) \qquad (14) $$

This constitutive relationship corresponds to a set of very simple ODE’s that can be solved analytically, either directly or using the Laplace transformation. The analytical solution for a tensile test at a constant strain rate ($\varepsilon(t) = \dot{\varepsilon}_0 t$) is:

$$ \sigma_j(t) = \dot{\varepsilon}_0 \left[ p_j^0 \tau_j E^u (1 - e^{-t/\tau_j}) + G \left( \frac{3p_j^0 \dot{\varepsilon}_0 e^{-t/\tau_j} - p_j^0 (2\dot{\varepsilon}_0 + 1/\tau_j) e^{-\dot{\varepsilon}_0 t} + p_j^0 (1/\tau_j - \dot{\varepsilon}_0) e^{2\dot{\varepsilon}_0 t}}{\tau_j \dot{\varepsilon}_0 (1/\tau_j - \dot{\varepsilon}_0) (2\dot{\varepsilon}_0 + 1/\tau_j)} \right) \right] \qquad (15) $$
for each modal "branch"

$$ \sigma(t) = \sum_{j=1}^{N} \sigma_j(t) \quad \text{for the overall response} $$

The mathematical structure of this model (equation 14) rigorously corresponds to the Biot model for viscoelasticity (Biot, 1958) or to a generalized Zener model (Tschoegl, 1989), except that a recursive-type relation links both the relaxation times and weights (equations 10, 11; André et al., 2003). In the case where the relaxed stress is modeled through a linear relationship $\sigma^r = E^r \varepsilon$, this model has an exact equivalent analogic representation in terms of springs and dashpots, which is given in Fig. 5. The following correspondence between the springs and dashpot constants and the parameters of the model is required:
$$E_j^u = p_j^0 E^u = E_1$$

$$E_j^r = p_j^0 E^r = \frac{E_1 E_2}{E_1 + E_2}$$

$$\tau_j = \eta_j / (E_1 + E_2)$$

<figure><img src="image_7.png" /><figcaption>Figure 5: Approximate analogical representation of the model described by equations 14-15.</figcaption></figure>

In the case where the relaxed stress is modeled through the more complicated relation of equation 13, this schema no longer holds (there is no analogical schema available for the Kelvin unit) but the mathematical structure of equation 14 remains identical.

The relevance of this 3-parameter model will now be proven through the sensitivity analysis and parameter identification procedures.
# 4. Identification Procedure

All fundamental mathematical concepts of the Parameter Estimation Theory and major definitions have been summarized briefly in the Appendix to recreate a sequential presentation of the results. Nevertheless, they are also at the heart of the rational study presented below.

## 4.1 Adjustment Results

The model described above was applied to experimental data $\sigma^*(\varepsilon)$ to identify the following parameter vector: $\beta = [E^u, G, \tau_{max}^T]$. We made computations with Matlab Software and used the Levenberg-Marquardt algorithm to minimize the least square criterion (equation A-2). Once the convergence was reached and the parameters determined, we could calculate the sensitivity coefficients at this point of the parameter space (equations A-4a, A-4b) along with the resulting optimal confidence bounds (equation A-15). We roughly determined the initial $\beta_{init}$ values used as starting approximations to calculate the parameter vector, according to various strategies. We derived an initial instantaneous modulus from the approximate slope at the origin of the stress-strain curve. An initial value for the maximum relaxation time can be derived by assuming a pure exponential-type behavior of the rising part of the curve as a result of a step input. For the experiment performed at a strain rate of $5 \times 10^{-3} s^{-1}$ for example (Fig. 6), a characteristic strain of the order of 0.04 can be deduced from the rule of tangents, or logarithmic plots, and this leads to a typical time of about 8 seconds. Finally, a "visual" linear regression applied to the experimental data of Fig. 4 gives a first rough estimate of the order of magnitude of $G = 110/50 = 2.2 \text{ MPa}$. We intentionally performed the identification of the 3 parameters for each experiment carried out at each strain rate. That is to say, there is no parameter in the model that accounts for a strain rate dependency. This latter may be analyzed from the results brought about by the identification for the two parameters that can physically carry this dependency, namely $G$ and $\tau_{max}^T$.

Two different identification intervals can be considered. In cases where small strain levels are examined ($\varepsilon \le \varepsilon_{yield} \approx 0.13$ – Interval I), the way the relaxed state is modeled has no influence on the solution. Only $E^u$ and $\tau_{max}^T$ parameters can be identified. The modeling cannot be sensitive to the hardening phase since the polymer is still in the viscoelastic regime. Identifications performed on either $\sigma$ or $\sigma^*$ lead to the same results. In Fig. 6, both the experimental data and the model for the optimal identified values (reported in the insert in the figure) are plotted. Setting aside the initial points of the curve, which are highly dependent on the feedback-loop control efficiency, the identification residuals lie in the [-1, 1] MPa range, which represents a maximal discrepancy of 3% with respect to the yield stress value (of the order of 33 MPa). If the identification interval covers a high range of strains ($0 < \varepsilon < 2$ – Interval II), then the complete model can be used to account for the hardening.
<figure><img src="image_8.png" /><figcaption>Figure 6: Experimental data, fitted curve (equation 15) and residuals obtained when identifying parameters $E^u$ and $\tau_{\max}^T$ on interval I (Specimen A$_{II}$ - Strain rate of $5 \times 10^{-3} s^{-1}$).</figcaption></figure>

Figs. 7 and 8 show the experimental and fitted curves obtained for specimen A$_{II}$ for two different strain rates. Good agreement between the model and the experimental tensile curves was obtained. The residuals, that is the difference between the modeled data for the optimally identified parameters and the experimental data, are also plotted in the same figures (right axis MPa scale). In the ideal case, the model would perfectly represent the true behavior of the material and the experimental set-up would be in total agreement with the assumptions of the model (e.g., perfect ramp excitation, exact sensors etc.). Then these residuals would be unsigned and distributed around the zero-value. They should in fact represent the only measurement noise carried by signals $\sigma^*$ and $\epsilon_{11}$. If and only if the identifiability of the parameters of this model is proved, null residuals mean that the parameters are perfectly identified. The confidence interval is then mainly due to the contribution of the noise (equations A-15 or A-17). In the two experiments reported in Figs. 7 and 8, the residuals show a signed ("non-ideal") character, which means that some bias exists (incompatibility between model and data). This bias means an additional error in the estimated parameters (biased estimations). If the bias remains limited (of the order of the noise level as is the case in Figs. 7 and 8) then the parameters can still be considered to be well-identified. In the present case and according to the sensitivity analysis tools presented in the Appendix, this means that the estimation process produces the unique optimal set of parameters that minimizes the residuals in the least-square sense. Assuming no bias in either the data or the model but the same amount of noise, the variances on the unbiased parameter estimates can be shown to be lower-bounded by those calculated using equation A-15. The existence of a bias generally helps the experimentalist to improve the model-experiment couplet. For example in this metrological application, a possible origin of the bias may be the non-ideal input command. In our case, the input ramp $\epsilon(t) = \dot{\epsilon}t$ is
very good due to the special care taken in the selection of the PID settings used by the control feed-back
loop. The bias here could be due to the choice of working with the isochoric assumption to produce our
observable stresses and strains (equations 3 to 5), which is the price to pay to avoid using (bad) transversal
strain measurements to calculate the cross-section variations of the specimen with respect to time. This bias
can also be due to the approximation used in the relaxed state modeling (equation 14) as evidenced by the
curves shown in Fig. 4 since the linear behavior is not perfect. If this bias is too great, trying to reduce it by
tracking defects with respect to the idealized experimental conditions is the first thing to do. A refinement of
the model can be a legitimate second step. Otherwise, a refinement of the model stimulated by biased data
can be dramatic in view of the Parameter Estimation Problem (PEP). For example, trying to use here the
Arruda & Boyce model of equation 13 will introduce additional unnecessary parameters (N₀ and n in
isothermal conditions) which will distort the PEP and introduce poor identifiability conditions. In so far as
the residuals remain small enough with respect to the signal to noise ratio, which is the case here, the model
is judged efficient in terms of adjustment ability and one now has to investigate it in terms of identifiability
of its parameters, especially the instantaneous Young modulus.

<figure><img src="image_9.png" /><figcaption>Figure 7: Experimental data, fitted curve (Eq. 15) and residuals obtained when identifying parameters $E^u$, $G$ and $\tau_{max}^T$ on interval II (Specimen $A_{//}$ - Strain rate of $2.5 \times 10^{-3} s^{-1}$).</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>Figure 8: Experimental data, fitted curve (Eq. 15) and residuals obtained when identifying parameters $E^u$, G and $\tau_{max}^T$ on interval II (Specimen A$_{II}$ - Strain rate of $5 \times 10^{-3} s^{-1}$).</figcaption></figure>

It is also interesting to show the results of the identification procedure when the tensile test is followed by a relaxation (initiated around a strain of 1.9). The reduced model can still be used but with two variations. When relaxation starts, the material is highly deformed. This new excitation applies to a material which is totally different and implies a change in the modeling associated with the relaxed state and the time spectrum. For the relaxed state a simple constant stress value noted $\sigma^r = \sigma_\infty$, "equilibrium" state at rest, is assumed to be reached at long times. The spectrum conserves its properties of linearity, number of decades $d$ and number of modes $N$, but is simply shifted according to a new maximum relaxation time noted $\tau_{max}^R$ that is a new property of the material in a complete fibrillary state for such a high deformation. Note that in view of our constitutive model, the instantaneous, or unrelaxed, modulus has not to be considered different. For a whole tensile-relaxation test, the parameter vector is now $\beta = [E^u, G, \tau_{max}^T, \sigma_\infty, \tau_{max}^R]$. The results of the identification are shown in Fig. 9 for specimen B. The agreement is very good, both in the tensile and relaxation stages with the same $E^u$ value. Physically, this means that the SCP elastic, or instantaneous, properties are the same in the initial bi-phasic microstructure and in the fibrillary state, when unloading is triggered.

For both tensile and tensile-relaxation tests the fit is good, although a simple linear viscoelastic approach is considered for the model; the spectrum of relaxation times is fixed once with a single parameter and it does not change during the test.
<figure><img src="image_11.png" /><figcaption>Figure 9: Tensile test followed by a relaxation for specimen B. Experimental data, fitted curve (Eq. 15) and residuals plotted as a function of time. Identified parameters are $E^u$, $G$, $\tau_{max}^T$, $\sigma_\infty$ and $\tau_{max}^R$ - Strain rate of $5 \times 10^{-3} s^{-1}$ (400 s correspond to a true strain of 2).</figcaption></figure>

## 4.2 Sensitivity Analysis

Fig. 10 plots the normalized sensitivity coefficients (equation A-4b) of the model parameters and Table 1 presents the results of the stochastic analysis (Matrix $\tilde{\Delta}$ of equation A-17).

<table><thead><tr><th>Interval I</th><th>E<sup>u</sup></th><th>τ<sub>max</sub><sup>T</sup></th></tr></thead><tbody><tr><td>E<sup>u</sup></td><td><strong>2.15 %</strong></td><td>-0.9667</td></tr><tr><td>τ<sub>max</sub><sup>T</sup></td><td>-0.9667</td><td><strong>3.61 %</strong></td></tr></tbody></table>

a)

<table><thead><tr><th>Interval II</th><th>E<sup>u</sup></th><th>G</th><th>τ<sub>max</sub><sup>T</sup></th></tr></thead><tbody><tr><td>E<sup>u</sup></td><td><strong>1.37 %</strong></td><td>-0.18</td><td>-0.9988</td></tr><tr><td>G</td><td>-0.18</td><td><strong>0.06 %</strong></td><td>0.147</td></tr><tr><td>τ<sub>max</sub><sup>T</sup></td><td>-0.9988</td><td>0.147</td><td><strong>1.38 %</strong></td></tr></tbody></table>

b)

<caption>Table 1. Variance-correlation matrix $\tilde{\Delta}$ of parameters $E^u$ and $\tau_{max}^T$ on interval I (a) and $E^u$, $G$, $\tau_{max}^T$ on interval II (b).</caption>

Parameters $E^u$ and $\tau_{max}^T$ apparently seem to be strongly correlated over a large part of the curve (between $\epsilon \approx 0.1$ and $\epsilon \approx 0.8$). If the identification was only carried out on the 'plateau', it would be
impossible to gain simultaneous information about $E^u$ and $\tau_{\max}^T$. Indeed, the correlation coefficient between $E^u$ and $\tau_{\max}^T$ has a value very close to 1 (Tables 1a and 1b). Therefore further investigation is required. Figure 11 plots the sensitivity of the elastic modulus $E^u$ versus that of the relaxation time $\tau_{\max}^T$. One can see that at small strain values ($\varepsilon \le 0.13$), the sensitivities of $E^u$ and $\tau_{\max}^T$ show a nearly linear behavior with a zero intercept, which explains in part the high degree of correlation found. Nevertheless, there is sufficient divergence in their behavior to have good identification conditions. Indeed, the singular point around $\varepsilon_{11} = 0.13$ that marks the transition between the different regimes can be seen in the insert of figure 11. Furthermore, the relative errors calculated for both $E^u$ and $\tau_{\max}^T$ parameters are of the order of a few percent. The parameter $G$ can be estimated with a very low confidence bound (error of 0.06 %) as a consequence of its very weak or inexistent correlation with the other parameters ($\rho_{12} = -0.18$ and $\rho_{23} = 0.14$). The robustness of the identification has been checked by changing the initial values of the optimization algorithm. If two parameters were correlated, then the estimated values would have changed for each different run. Converging towards a single identified parameter vector, whatever its initialization, is a good verification that the PEP is well-conditioned.

Concerning the use of interval I to perform the identification of parameters $E^u$ and $\tau_{\max}^T$ only, the sensitivity analysis helps in their definition. For strain values below 0.13, the sensitivity to $G$ is only 5% of the maximum sensitivities of $E^u$ and $\tau_{\max}^T$ and therefore this parameter can be omitted. When Tables 1a and 1b are compared, it can be seen that the variances of parameters $E^u$ and $\tau_{\max}^T$ are lowered by a factor of 2 when considering interval II instead of interval I. Therefore using the whole curve is recommended. Additionally, the sensitivities of both $E^u$ and $\tau_{\max}^T$ remain at their maximum value after the yield. This is physically sound as they actually determine the value of the yield stress plateau.
<figure><img src="image_12.png" /><figcaption>Figure 10: Normalized sensitivities to $E^u$, G and $\tau_{max}^T$ parameters on interval II.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Figure 11: Normalized Sensitivities $\tilde{X}(E^u)$ as a function of $\tilde{X}(\tau_{max}^T)$.</figcaption></figure>
## 5. Identification results and validation of the approach

The MBM approach can be validated by analyzing the results obtained for tensile tests performed at different strain rates. Table 2 gives the values estimated for specimen A on both interval I ($E^u, \tau_{max}^T$) and interval II ($E^u, \tau_{max}^T, G$).

<table><caption>Table 2. Identified values of parameters $E^u, \tau_{max}^T$ (Identification interval I) and $E^u, \tau_{max}^T, G$ (Identification interval II) for Specimen A<sub>II</sub> at different strain rates</caption><thead><tr><th>Strain rate</th><th>5.10<sup>-5</sup> s<sup>-1</sup></th><th>4.10<sup>-4</sup> s<sup>-1</sup></th><th>2.5.10<sup>-3</sup> s<sup>-1</sup></th><th>5.10<sup>-3</sup> s<sup>-1</sup></th><th>10<sup>-2</sup> s<sup>-1</sup></th></tr></thead><tbody><tr><td><strong>Eu (MPa) Interval I</strong></td><td>2777</td><td>2783</td><td>2787</td><td>2828</td><td>2846</td></tr><tr><td>Eu (MPa) Interval II</td><td>2811</td><td>2809</td><td>2940</td><td>2812</td><td>2772</td></tr><tr><td><strong>τ<sub>max</sub><sup>T</sup> (s) Interval I</strong></td><td>433.64</td><td>64.17</td><td>11.40</td><td>6.01</td><td>3.17</td></tr><tr><td>τ<sub>max</sub><sup>T</sup> (s) Interval II</td><td>419.72</td><td>63.23</td><td>11.02</td><td>6.01</td><td>3.21</td></tr><tr><td>G (MPa)</td><td>1.53</td><td>1.85</td><td>2.31</td><td>2.31</td><td>2.46</td></tr></tbody></table>

Firstly, we compared the measurements obtained for the elastic modulus to those given by two other physical techniques and for the two specimens, A and B. Secondly, we analyzed the consistency of the values obtained for the 3 parameters of this rheological model.

### Elastic moduli measurements - Comparison with other techniques.

As can be seen from Table 2, the elastic moduli identified for sample A, specimen $A_{II}$, lie in the range from 2720 to 2950 MPa for all experiments including the results obtained from both identification intervals and more than 3 experiments were reproduced at each different strain rate. This represents a variation of ±4% around the central value, which corresponds to the order of magnitude of the estimated variances yielded by the stochastic analysis (Table 1b). The uncertainty of these measurements has been calculated by changing the value of the identified $E^u$ within a ±Δ$E^u$ interval so that the direct model produces theoretical curves that flank the experimental ones. This allows the bias effects evidenced earlier by the residual plotting to be taken into account. For specimen $A_{II}$, at a strain rate of $5 \times 10^{-3}$ s<sup>-1</sup>, this uncertainty
is ±180MPa that corresponds to a variation of ±6,5% around the nominal value. The results for all specimens tested are: 2830±180 MPa, 2940±200 MPa, 2400±260 MPa, and 2270±150 MPa for specimens A$_{‖}$, A$_{c}$, A$_{⊥}$ and B, respectively.

To check the validity of the results obtained by this MBM approach, two independent direct measurements of the elastic modulus of the material were carried out. One relies upon the Pulse-Echo technique, which is based on ultrasound measurements (Optel OPBOX system). Direct measurements of the velocity of a sonic wave at 5MHz were made in the longitudinal and transversal directions of the sample which allowed an elastic modulus to be calculated. The other measurement depends upon the Nanoindentation Continuous Stiffness Method. Indentation tests were performed with Nanoindenter II$^{S}$ at the L.M.A Femto-ST laboratory (France). Harmonic excitations of the pyramidal Berkovich probe tip at 45 Hz allowed the measurement of a complex modulus (and hardness) of the material with good precision (Olivier et al., 1992; Le Rouzic et al., 2009).

Table 3 gathers together the values of the elastic modulus obtained for all 500 Natural HDPE A (A$_{‖}$, A$_{⊥}$ and A$_{c}$) and B specimens by the three metrological techniques (MBM, PE and Nano subscripts, are for respectively the Model Based Metrology applied to tensile curves, Pulse-Echo and Nanoindentation techniques).

<table><caption>Table 3. Comparison between elastic modulus values given by standards applied to tensile curves, the present model-based parameter estimation, the Pulse-Echo technique and the Nanoindentation technique, for all A (A$_{‖}$, A$_{⊥}$, A$_{c}$) and B samples (in MPa). Discrepancy percentages in brackets were calculated with reference to E<sup>u</sup><sub>MBM</sub>.</caption><thead><tr><th>Specimen</th><th>Manufacturer</th><th>ISO 527-1<br>ε̇<sub>N</sub> = 5<br>mm/min<br>(values for<br>ASTMD638)</th><th>ISO 527-1<br>ε̇<sub>N</sub> = 6<br>mm/min<br>(values for<br>ASTMD638)</th><th>Present Model-<br>Based parameter<br>estimation<br>E<sup>u</sup><sub>MBM</sub></th><th>Ultrasonic<br>technique<br>E<sup>u</sup><sub>PE</sub></th><th>Nanoindentation<br>E<sup>u</sup><sub>Nano</sub></th></tr></thead><tbody><tr><td rowspan="3">A</td><td>A<sub>‖</sub></td><td>1138<br>(1141)</td><td>1475<br>(1479)</td><td>2830</td><td>2780<sup>±40</sup><br>(-2%)</td><td>2380<sup>±40</sup><br>(-16%)</td></tr><tr><td>A<sub>c</sub></td><td>-</td><td>-</td><td>2940</td><td>2800<sup>±230</sup><br>(-5%)</td><td>2730<sup>±60</sup><br>(-8%)</td></tr><tr><td>A<sub>⊥</sub></td><td>-</td><td>-</td><td>2400</td><td>2780<sup>±40</sup><br>(= E<sup>u</sup><sub>A<sub>‖</sub></sub>)</td><td>2320<sup>±40</sup><br>(-3%)</td></tr><tr><td>B</td><td></td><td>772<br>(763)</td><td>-</td><td>2270</td><td>2220<sup>±230</sup><br>(-2%)</td><td>2180<sup>±110</sup><br>(-5%)</td></tr></tbody></table>
The data obtained by applying standards on tensile curves are also given. These last values show that standards are far from producing (i) reliable data and (ii) physically founded estimations. The discrepancy between standards and experimental techniques can be as high as 300% and the values obtained depend on the applied displacement rate. Nevertheless, using tensile curves to precisely identify the elastic modulus of polymers is possible if a correct MBM is used. This can be asserted from the comparison made in Table 3 where the values yielded by the two other physical techniques have all the corresponding uncertainty intervals.

For all specimens, the three techniques gave very similar results. The agreement was very good between MBM and PE values, but local differences from the Nanoindentation test occurred. These differences must be analyzed carefully. For specimen A//, the value given by the “Nano” technique is 16% lower than by the MBM and PE techniques. This can be explained since this former technique probes a small volume just behind the surface. The value identified is then sensitive to a skin effect which is generally observed on thick extruded samples (such as for the 6mm A// specimen). In contrast, both the MBM and PE methods probe the entire volume. As specimen B was manufactured directly with a 4 mm thickness, this effect is not visible within the confidence intervals of the three measurements. Looking at specimen A_c, this effect should disappear. Higher values of $E_{Nano}^u$ are indeed obtained and the three techniques give more convergent results. This proves that (i) from a physical point of view core samples are slightly more rigid than full samples, with less differences in microstructure and (ii) from a metrological point of view, the MBM technique is sensitive to such small variations.

Examining the incidence of these results on the material property, it is clear that even though this commercial HDPE is produced under the same reference for specimens A and B, it can exhibit major changes in its elastic behaviour (around 600 MPa). The increase in elasticity of A//, A_c samples compared to sample B indicates an oriented texture due to the extrusion process. The presence of such anisotropy is confirmed by a positive anisotropy index measured for the nondeformed HDPE A// specimens, thanks to X-ray microtomography (Blaise et al., 2010). This result conforms to those reported by Legros et al. (1999), where un-oriented and oriented HDPE are investigated by the Pulse-Echo technique. Oriented HDPE produces a higher longitudinal and shear wave signal and hence a higher elastic modulus than does un-oriented HDPE.

Lastly, as expected from an instantaneous modulus, the identified values can be considered as independent from each other.

## Analysis of the results for the other parameters

The results were next analyzed with respect to viscosity, hyper-elasticity (or hardening) and time-dependency. Table 4 brings together the estimated maximal relaxation times $\tau_{max}^T$, for specimens A// and B
for different strain rates. The corresponding Weissenberg number (similar to the Deborah number) is also shown, as is the identified hardening modulus G.

<table><caption>Table 4. Average estimated maximal relaxation time $\tau_{max}^T$, corresponding Weissenberg number $We_{max}$ and hardening modulus G for specimens A// and B at different strain rates.</caption><thead><tr><th>Specimen</th><th>Strain rate (s<sup>-1</sup>)</th><th>τ<sub>max</sub><sup>T</sup> (s)</th><th>We<sub>max</sub></th><th>G (MPa)</th></tr></thead><tbody><tr><td rowspan="5">A<sub>//</sub></td><td>5.10<sup>-5</sup></td><td>419.72</td><td>0.0209</td><td>1.53</td></tr><tr><td>4.10<sup>-4</sup></td><td>63.23</td><td>0.0252</td><td>1.85</td></tr><tr><td>2,5.10<sup>-3</sup></td><td>11.02</td><td>0.0276</td><td>2.31</td></tr><tr><td>5.10<sup>-3</sup></td><td>6.01</td><td>0.0301</td><td>2.31</td></tr><tr><td>10<sup>-2</sup></td><td>3.21</td><td>0.0321</td><td>2.46</td></tr><tr><td rowspan="2">B</td><td>5.10<sup>-3</sup></td><td>7.16</td><td>0.0358</td><td>2.30</td></tr><tr><td>10<sup>-2</sup></td><td>3.79</td><td>0.0379</td><td>2.39</td></tr></tbody></table>

The Weissenberg number corresponds to the ratio between the time constant $t_{mat}$ characterizing the intrinsic ‘fluidity’ of the material and the time scale of the experiment or of the observer $t_{exp}$. The fluidity of the material is inversely proportional to the Weissenberg number. In the present case, its maximal value can be calculated by the following formula:

$$We_{max} = \frac{t_{mat}}{t_{exp}} = \tau_{max} \dot{\epsilon} \qquad (16)$$

where $\tau_{max}$, the maximum relaxation time of the spectrum, is used for the material characteristic time and $1/\dot{\epsilon}$ for the experimental characteristic excitation time.

The Weissenberg number varies logarithmically from 0.02 to 0.032 (+50%) over nearly 3 decades of strain rates. A regression gives $We_{max} \propto 0.005 \log(\dot{\epsilon})$. The increase in the Weissenberg number agrees
well with the observation that plasticity of the material is enhanced when strain rates are diminished. The low values of $We_{max}$ also indicate that the time response of the material always remains two logs below the time scale imposed by the experiment. This argues in favor of a similar scenario for the succession of internal mechanisms taking place at the microstructure levels. Finally, it has been shown by André et al., (2012) that a characteristic time of the order of 1s can be derived from energetic measurement for a strain rate of $\dot{\epsilon}=10^{-2} s^{-1}$, which validates the order of magnitude found here for $\tau_{max}^T$. Additionally, a good agreement of the values found for $E^u$ and $\tau_{max}^T$ with respect to the strain rate dependency of the yield stress and the identification intervals was found. Re-examining Table 2 shows two things:

- When identified values of $E^u$ increase slightly between the two identification intervals I and II, then the corresponding $\tau_{max}^T$ value decreases, and the contrary is true.

- The $E^u \tau_{max}^T$ product increases with strain rate.

Both effects can be explained as follows. At short times, the model (Eq. 15) acts like
$$ \sigma(t) = \sum_j \dot{\epsilon}_0 p_j^0 \tau_j E^u (1 - e^{-t/\tau_j}) $$
The behavior of the material is roughly considered to be described by a simple exponential relation like
$$ \sigma(t) = \sigma_{yield} (1 - e^{-t/\tau}) $$
As a result the pre-factors can be compared through $k \dot{\epsilon}_0 \tau_{max} E^u \square \sigma_{yield}$ (k is a constant which weights $\tau_{max}$ to have a barycentric relaxation time corresponding to $\langle \tau \rangle$). The yield stress increases with strain rate in the experimental curves and for a given strain rate a constant $\sigma_{yield}$ imposes that $E^u$ varies inversely proportionally to $\tau_{max}$.

Concerning the value of the hardening modulus $G$ that describes the behavior of HDPE at a high strain, the identified values for all strain rates are also in a logarithmic dependency with respect to the strain rates. The values identified for $G$ are close to those found in the experimental results of uniaxial tensile tests reported by Haward for HDPE (Haward, 1993; Haward, 2007) and the value found by Bartczak and Kozanecki (2005) for linear polyethylene in plane-strain compression. The estimated hardening modulus appears to be the same for both specimens A and B, which is in line with the idea that the microstructure transformation in uniaxial drawing is highly reproducible for a given polymer even if the initial state may be different.

Considering the consistency of all estimated parameters, one can remember that the sensitivity analysis showed that the identification of parameters $E^u$, $\tau_{max}^T$ and $G$ of the reduced model (Eq. 15) was possible with good variances (confidence bounds) and that the order of magnitude was confirmed by the uncertainties obtained from the residuals analysis. For each strain rate tested, the same value was identified for the elastic modulus, as expected from the instantaneous character of this physical parameter. Young's modulus results obtained with the MBM approach have been corroborated by two other independent and direct physical measurements that probe the material at high frequencies (MHz) or at small spatial scales (nm-µm). This is clear evidence of the consistency of the approach with the underlying physics, because this agreement is the
corollary of an identified maximum relaxation time of the order of a few seconds. As a spectrum of 6 decades is considered in the model to arrive at convergence of estimates, this means that relaxation times as small as 1 µs are necessary to allow for a good description of the data and a determination of the 'truly' instantaneous modulus. The pulse-echo technique relies on excitation at 5 MHz, which means that the time scale required by the model is in perfect agreement with this technique, and supports the identification of an "instantaneous" modulus. For the Nanoindentation technique, the spatial size of the representative elementary volume is less than a micrometer, which means that small time scale events are also probed.

## 6. Conclusion

In the light of the characterization of SCP trough tensile tests, both ideas of using an appropriate reduced metrological model in the relationship of sensitivity analysis principles have been shown to be key concepts for a correct estimation of the parameters of the model through optimization procedures. The 3 parameter rheological model we have developed is sufficient to apprehend correctly the physics of the tensile behavior and to give reliable estimates of the parameters. Based on such measurements, much precise information can be gained concerning the behavior of the polymer. An illustration has been given of the apparently simple problem of measuring an instantaneous Young modulus, where contemporary recommended standards lead to measurements twice as low as the experimental values. If good confidence bounds are obtained for the estimated parameters, associated phenomena like strain rate effects in our case, can be more securely quantified. Additionally, any deviation in the residuals can be considered to result from a bias (in the model or in the experimental conditions) which in all cases, give the modeler-experimentalist new trails to improve the model and experimental conditions.

Like various other works in different fields of physics, we hope that the research described here demonstrates that major improvements can be expected if simple identifiability principles are used more frequently in material science studies.
Appendix: Sensitivity Analysis Fundamentals

Let $y^{\text{meas}}(t)$ be the measured stress output of a system (our material specimen), and $y(\beta, t)$ the theoretical stress output of the behavioral model with a parameter vector $\beta$ of dimension $n$, representing in our case the $n$ constitutive material parameters. The output error $e(t)$ or the residuals can be defined as:

$$e(t) = y^{\text{meas}}(t) - y(\beta, t) \qquad (A-1)$$

The role of an estimator is to minimize this output error. A least square criterion can be used classically and is written as follows:

$$E_{LS} = \sum_{i=1}^{m} (y^{\text{meas}}(t_i) - y(\beta, t_i))^2 \qquad (A-2)$$

where the sum is made over the successive acquisition times $t_i$, up to a total number of measurement data points $m$.

The minimization of the criterion is done when its derivatives with respect to parameters $\beta_j$ ($j=1...n$) are null:

$$\begin{align}
\forall j \in [1, n], \quad \frac{\partial E_{LS}}{\partial \beta_j} = 0 \quad &\rightarrow \sum_{i=1}^{m} \left( \frac{\partial y(\beta, t_i)}{\partial \beta_j} [y^{\text{meas}}(t_i) - y(\beta, t_i)] \right) = 0 \tag{A-3} \\
&\rightarrow \sum_{i=1}^{m} \left( X_j(\beta, t_i) [y^{\text{meas}}(t_i) - y(\beta, t_i)] \right) = 0 \nonumber
\end{align}$$

From this equation, the sensitivity coefficient vector component $X_j$ associated to the parameter $\beta_j$ can be clearly recognized as:

$$X_j(\beta, t_i) = \frac{\partial y(\beta, t_i)}{\partial \beta_j} \qquad (A-4a)$$

Normalized sensitivity coefficients:

$$\tilde{X}_j = \beta_j \frac{\partial y(\beta, t_j)}{\partial \beta_j} \qquad (A-4b)$$
are introduced to check graphically the level of sensitivity and possible correlations between parameters: different intervals of the time-independent variable *t* may be advantageously considered. These coefficients are dimensionally homogeneous to the signal itself, the stress in MPa.

Sensitivity coefficients express how much a model reacts to some small variation of the parameters. Sensitivity coefficients play a fundamental role in the conditioning of the inverse parameter estimation process and, as a consequence, in the errors made in the estimates of confidence bounds. It is obvious that major sensitivities are sought for, when designing an experiment for metrological purposes. In almost all cases, sensitivity coefficients are non-linear functions of the parameters themselves when the model *y* is a non-linear function of *β*<sub>j</sub>.

Switching to a matrix formulation and defining vectors **Y**<sup>meas</sup> and **Y** as:

$$ 
\mathbf{Y}^{\text{meas}} = \begin{bmatrix} y^{\text{meas}}(t_1) \\ y^{\text{meas}}(t_2) \\ \vdots \\ y^{\text{meas}}(t_m) \end{bmatrix} \quad \text{and} \quad \mathbf{Y} = \begin{bmatrix} y(\mathbf{\beta}, t_1) \\ y(\mathbf{\beta}, t_2) \\ \vdots \\ y(\mathbf{\beta}, t_m) \end{bmatrix} \qquad (\text{A-5}) $$

makes the minimization process expressed as

$$ \mathbf{X}^T (\mathbf{Y}^{\text{meas}} - \mathbf{Y}) = 0 \qquad (\text{A-6}) $$

with **X**, the *m*×*n* sensitivity rectangular matrix, where *m* is the number of rows (dimension of the experimental observation time vector **t** ) and *n* is the number of columns (dimension of the parameter vector **β** ):

In the case of a linear model with respect to the parameters (linear estimation problem) for which the matrix of the sensitivity coefficients does not depend on the parameters, we have:

$$ \mathbf{Y} = \mathbf{X}\hat{\mathbf{\beta}} \qquad (\text{A-7}) $$

The estimated parameter vector, denoted by $\hat{\mathbf{\beta}}$, corresponds to the value reached by **β** when the criterion is minimized. Therefore using (A-7) we can rewrite equation (A-6) as:

$$ \mathbf{X}^T (\mathbf{Y}^{\text{meas}} - \mathbf{X}\hat{\mathbf{\beta}}) = \mathbf{0} \qquad (\text{A-8}) $$

The (A-8) relation can be inverted to obtain the expression of $\hat{\mathbf{\beta}}$ in the case of a linear estimation problem:

$$ \hat{\mathbf{\beta}} = (\mathbf{X}^T \mathbf{X})^{-1} \mathbf{X}^T \mathbf{Y}^{\text{meas}} \qquad (\text{A-9}) $$
However, most parameter estimation problems are not linear and require an iterative linearizing procedure. This is obtained by developing the solution at rank $n+1$ in the neighborhood of the solution obtained for the prior iteration at rank $n$:

$$ \mathbf{Y}^{(n+1)} = \mathbf{Y}^{(n)} + \mathbf{X}^{(n)} (\hat{\mathbf{\beta}}^{(n+1)} - \hat{\mathbf{\beta}}^{(n)}) \quad (A-10) $$

If the relation (A-6) is written at rank $n+1$ for parameters estimated at rank $n$, we obtain: $\mathbf{X}^{T(n)}(\mathbf{Y}^{\text{meas}} - \mathbf{Y}^{(n+1)}) = \mathbf{0}$ and when this is combined with the relation (A-10) we obtain the following relation of recurrence between estimated parameters at rank $n+1$ and rank $n$:

$$ \hat{\mathbf{\beta}}^{(n+1)} = \hat{\mathbf{\beta}}^{(n)} + (\mathbf{X}^{T(n)}\mathbf{X}^{(n)})^{-1} \mathbf{X}^{T(n)} (\mathbf{Y}^{\text{meas}} - \mathbf{Y}^{(n)}) \quad (A-11) $$

which defines the iterative procedure that can be used to estimate the parameters (Gauss-Newton algorithm).

The statistical estimator's properties depend on the noise $\varepsilon(t)$ of the signal. If the theoretical model is assumed to be unbiased (perfect) then we have:

$$ Y^{\text{meas}}(t) = Y(t, \mathbf{\beta}) + \varepsilon(t) \quad (A-12) $$

If classical statistical assumptions are made concerning the experimental noise $\varepsilon(t)$ in the measured stress signal (Beck et al., 1977), it is possible to obtain an approximation of the errors that can be made in the estimation process for the different parameters. These assumptions are:

a. zero mean value of the signal in the absence of excitation, which corresponds to a zero expectancy for noise (expected value $E(\varepsilon) = 0$);

b. constant variance or standard deviation of the noise: $V(\varepsilon) = \sigma_0^2$.

In the case of a 1st order linearized estimation, its expected value can be proved to be:

$$ E(\hat{\mathbf{\beta}}) = \mathbf{\beta} \quad (A-13) $$

This means that there is no error or bias made on the identified parameters.

The variance-covariance matrix $\Delta$ on the estimated parameters (generalization of the scalar-valued variance to a higher dimension) involves the sensitivity coefficients. It is calculated as

$$ \Delta = E \left[ (\hat{\mathbf{\beta}} - E(\hat{\mathbf{\beta}})) (\hat{\mathbf{\beta}} - E(\hat{\mathbf{\beta}}))^{T} \right] = \sigma_{0}^{2} (\mathbf{X}^{T} \mathbf{X})^{-1} \quad \text{which in the expanded form gives:} $$
$$ \Delta = \begin{bmatrix} V(\beta_1) & \operatorname{cov}(\beta_1, \beta_2) & \dots & \operatorname{cov}(\beta_1, \beta_p) \\ \operatorname{cov}(\beta_1, \beta_2) & V(\beta_2) & \dots & \operatorname{cov}(\beta_2, \beta_p) \\ \vdots & \vdots & \ddots & \vdots \\ \operatorname{cov}(\beta_1, \beta_p) & \operatorname{cov}(\beta_2, \beta_p) & \dots & V(\beta_p) \end{bmatrix} \quad (A-14) $$

A stochastic analysis can be made that consists of calculating the variance-covariance matrix through equation (A-14) theoretically according to a given noise and set of parameters. This symmetric square matrix has a dimension equal to the number of parameters. The diagonal terms correspond directly to the variance of each parameter $V(\beta_j)$. They can be used to determine the error made for each parameter. This error can be expressed as a %:

$$ Err(\beta_j) = \frac{\sqrt{V(\beta_j)}}{\beta_j} \quad (A-15) $$

The off-diagonal terms can be used to calculate the correlation coefficients $\rho_{mn}$ which express the degree of correlation of the parameters:

$$ \rho_{rs} = \frac{\operatorname{cov}(\beta_r, \beta_s)}{\sqrt{V(\beta_r)} \sqrt{V(\beta_s)}} \quad (A-16) $$

The $|\rho_{rs}|$ values lie between 0 and 1. In the case of a model with strongly correlated parameters, the correlation coefficients are close to 1, which means that two columns of the sensitivity matrix $X$ are nearly proportional to each other. The resulting confidence bounds interval for two correlated parameters are therefore generally very high. This means that a large number of solutions exist for these two parameters to allow for a good fit to the experimental curve. A deterministic algorithm like the steepest gradient technique, used for the minimization process, is consequently very sensitive to the initial estimate made for the parameters. A strategy to produce first approximations using the physical background is highly recommended. But even so, the estimation problem is ill-conditioned and indicates to the experimentalist that the physical description involved is probably not appropriate and must be changed.

In the following considerations, the identifiability of the model parameters will be analyzed through the $\tilde{\Delta}$ array, which combines the variance-covariance matrix and the correlation matrix. $\tilde{\Delta}$ puts the main diagonal of $\Delta$ in place and the correlation coefficients on the off-diagonal terms.

$$ \tilde{\Delta} = \begin{bmatrix} Err(\beta_1) & \rho_{12} & \dots & \rho_{1p} \\ \rho_{12} & Err(\beta_2) & \dots & \rho_{2p} \\ \vdots & \vdots & \ddots & \vdots \\ \rho_{1p} & \rho_{2p} & \dots & Err(\beta_p) \end{bmatrix} \quad (A-17) $$
When there is no model bias, or perfect agreement between the conditions of the experiment and the model based on it, the relation between the estimated parameter vector for the non-linear estimation problem and its ‘exact’ value can be given at convergence by the following formula:

$$ \hat{\beta} = \beta + (X^T X)^{-1} X^T \varepsilon \qquad (A-18) $$

In this case, the residuals $e(t)$ that are the difference between the model and the data correspond exactly to the noise. Their standard deviation corresponds to the experimental standard deviation of the noise and the residuals remain unsigned with no large fluctuation around the zero level.

The matrix $\tilde{\Delta}$ is a good tool to investigate the identification conditions of a parameter estimation problem.
References

André, S., Meshaka, Y., Cunat, C., Rheological constitutive equation of solids: a link between models based on irreversible thermodynamics and on fractional order derivative equations. Rheologica Acta. 42, 500-515 (2003)

André, S., Renault, N., Meshaka, Y., Cunat, C., From the thermodynamics of constitutive laws to thermomechanical experimental characterization of materials: An assessment based on inversion of thermal images, Continuum Mechanics and Thermodynamics, 24, 1-20 (2012) doi: 10.1007/s00161-011-0205-x

Arruda, E.M., Boyce, M.C., A three-dimensional constitutive model for the large stretch behavior of rubber elastic materials, Journal of the Mechanics and Physics of Solids, 41-2, 389. (1993)

Aster, R.C., Borchers, B., Thurber, C.H, Parameter Estimation and Inverse Problems, S²d ed., Academic Press. (2013)

Baravian, C., André, S., Renault, N., Moumini, N., Cunat C., Optical techniques for in situ dynamical investigation of plastic damage, Rheologica Acta, 47, 555-564 (2008)

Bartczak, Z., Kozanecki, M., Influence of molecular parameters on high-strain deformation of polyethylene in the plane-strain compression, Part I. Stress-strain behaviour. Polymer. 46, 8210-8221 (2005)

Beck, J.V., Arnold, K.J., Parameter Estimation in Engineering and Science, John Wiley & Sons, New York (1977)

Biot, M.A., Linear thermodynamics and the mechanics of solids. Proceeding of the third US National Congress of Applied Mechanics, ASME. 1, 1-18 (1958)

Blaise, A., Baravian, C., André, S., Dillet, J., Michot, L.J., Mokso, R., Investigation of the mesostructure of a mechanically deformed HDPE by synchrotron microtomography, Macromolecules, 43, 8143-8152 (2010). (DOI: 10.1021/ma101033b)

Blaise, A., Baravian, C., Dillet, J., Michot, L.J., André, S., Characterization of the mesostructure of HDPE under "in-situ" uniaxial tensile test by incoherent polarized steady-light transport, Journal of Polymer ScienceB: Polymer Physics, 50(5), 328-337, (2011). DOI: 10.1002/polb.23020

Callen, H.B., Thermodynamics and an Introduction to Thermostatics, Wiley, New York (1985)

Cooreman, S.; Lecompte, D.; Sol, H., Vantomme, J., Debruyne, D., Elasto-plastic material parameter identification by inverse methods: Calculation of the sensitivity matrix, Int. J. of solids and Structures, 44, 13, 4329-4341 (2007)

Cunat, C., A thermodynamic theory of relaxation based on a distribution of non-linear processes. J. Non-Crystalline Solids. 131/133, 196-199 (1991)

Cunat, C., The DNLR approach and relaxation phenomena : part I : Historical account and DNLR formalism. Mech. Of Time-Depend. Mater. 5, 39-65 (2001)

Faccio-Toussaint, E., Ayadi, Z., Pilvin, P., Cunat, C., Modeling of the Mechanical Behavior of a Nickel alloy by Using a Time-Dependent Thermodynamic Approach to Relaxations of Continuous Media. Mech. Of Time-Depend. Mater. 5, 1-25 (2001)

De Donder, T., Thermodynamic theory of affinity: A book of principle, Oxford University Press. (1936)
Farge, L., André, S., Meneau, F., Dillet, J., Cunat, C., A common multiscale feature of the deformation mechanisms of a semi-crystalline polymer, Macromolecules, 46(24), 9659-9668 (2013). dx.doi.org/10.1021/ma4019747/

Farge, L., André, S., Pawlak, A., Baravian, C., Irvine, S.C., Philippe, A.M., A study of the deformation-induced whitening phenomenon for cavitating and non cavitating semi-crystalline polymers, Journal of Polymer ScienceB: Polymer Physics, 51(10), 826-841, (2013). DOI: 10.1002/polb.23267

Gaucher-Miri, V., François, P., Séguéla, R., 1996. On the mechanisms of initiation and propagation of plastic instability in polyethylene under tensile drawing, J. Polymer Science, PartB: Polymer Physics, 34, 1113-1125.

Gutman, E., Mechanochemistry of Materials. Cambridge International Science Publishing, Cambridge (1988)

G'sell, C., Hiver, J.M., Dahoun, A., Experimental characterization of deformation damage in solid polymers under tension, and its interrelation with necking. Int. J. Solids and Structures. 39, 3857-3872 (2002)

Haward, R.N., Strain Hardening of Thermoplastics. Macromolecules. 26, 5860-5869 (1993)

Haward, R.N., Strain hardening of High Density Polyethylene. J. Polymer Sciences, Part B: Polymer Physics. 45, 1090-1099 (2007)

Hiss, R., Hobeika, S., Lynn, C. and Strobl, G., Network Stretching, Slip Processes, and Fragmentation of Crystallites during Uniaxial Drawing of Polyethylene and Related Copolymers. A Comparative Study. Macromolecules. 32, 4390-4403 (1999)

Krempl, E., Relaxation behaviour and modeling. International Journal of Plasticity, 17, 1419-1436 (2001)

Kuiken, G.D.C., Thermodynamics of Irreversible Processes: Applications to Diffusion and Rheology, Wiley (1994)

Legros, N., Jen, C.-K., Ihara, I., Ultrasonic evaluation and application of oriented polymer rods. Ultrasonics. 37, 291-297 (1999)

Le Rouzic, J., Delobelle, P., Vairac, P., Cretin, B., Comparaison of three different scales techniques for the dynamic mechanical characterization of two polymers (PDMS and SU8). The European Physical Journal Applied Physics. 48, 11201 (2009)

Maugin, G., Muschink, W., Thermodynamics with Internal Variables. Part I General Concepts. Journal of Non-Equilibrium Thermodynamics. 19, 217 (1994)

Meixner, J.Z., Thermodynamik und Relaxationserscheinungen. Naturforsch. 4a, 504-600 (1949)

Mrabet, K., Rahouadj, R., Cunat, C., An irreversible model for semicristalline polymers submitted to multisequence loading at large strain. Polymer Engineering and Science. 45 (1), 42 (2005)

Negahban, M., Experimentally Evaluating the Equilibrium Stress in Shear of Glassy Polycarbonate. Journal of Engineering Materials and Technology, ASME. 128, 537-542 (2006)

Nowick, A. S. and Berry, B. S., Anelastic relaxation in crystalline solids, Academic Press, New York and London, (1972).

Oliver, W.C., Pharr, G.M., An improved technique for determining hardness and elastic modulus using load and displacement sensing indentation experiments. J. Mat. Res.. 7, 6, 1563 (1992)

Piché, S., Ultrasonic velocity measurement for the determination of density in polyethylene. Polymer Engineering & Science. 24, 1354-1358 (1984)
Ponthot, J-P.; Kleinermann, J-P., A cascade optimization methodology for automatic parameter identification and shape/process optimization in metal forming simulation, Computer Methods in applied Mechanics and Engineering, 195, 41-43, 5472-5508 (2006)

Prigogine, I., Defay, R.: Thermodynamique chimique conformément aux méthodes de Gibbs et De Donder, Tomes I-II. Gauthier-Villars (1944-1946).

Qasmi, M., Delobelle, P., Influence of the average roughness Rms on the precision of the Young's modulus and hardness determination using nanoindentation technique with a Berkovich indenter' Surf. Coat. Techn.. 201, 1191 (2006)

Rahouadj, R., Ganghoffer, J.F., Cunat, C., A thermodynamic approach with internal variables using Lagrange formalism. Part I: General framework. Mechanics Research Communications. 30 (2), 109 (2003)

Treloar, L.R.G., The Physics of Rubber Elasticity. Clarendon Press, Oxford, UK. (1975)

Tschoegl, N.W., The Phenomenological Theory of Linear Viscoelastic Behaviour, An Introduction, Springer-Verlag. (1989)

Van Melick, H.G.H., Govaert, L.E., Meijer, H.E.H., On the origin of strain hardening in glassy polymers. Polymer. 44, 2493 (2003).

Walter, E, Pronzato, L, Identification of Parametric Models from Experimental Data Springer-Verlag Heidelberg (1997).

Wang, M.C., Güth, E., Statistical theory of networks of non-gaussian flexible chains, Journal of chemical physics, 20(7), 1144 (1952).

Ye, J, Andre, S., Farge, L., Kinematic study of necking in a semi-crystalline polymer through 3D digital image correlation, Int. J. Solids & Structures, 59, 58-72 (2015).