Micromechanical Modeling of the Biaxial Deformation-Induced Phase Transformation in Polyethylene Terephthalate

In this paper, a micromechanics-based constitutive representation of the deformation-induced phase transformation in polyethylene terephthalate is proposed and verified under biaxial loading paths. The model, formulated within the Eshelby inclusion theory and the micromechanics framework, considers the material system as a two-phase medium, in which the active interactions between the continuous amorphous phase and the discrete newly formed crystalline domains are explicitly considered. The Duvaut–Lions viscoplastic approach is employed in order to introduce the rate-dependency of the yielding behavior. The model parameters are identified from uniaxial data in terms of stress–strain curves and crystallization kinetics at two different strain rates and two different temperatures above glass transition temperature. Then, it is shown that the model predictions are in good agreement with available experimental results under equal biaxial and constant width conditions. The role of the crystallization on the intrinsic properties is emphasized thanks to the model considering the different loading parameters in terms of mechanical path, strain rate and temperature.


Introduction
Polyethylene terephthalate (PET) is one of the most used engineering polymers in applications where lightweight, chemical resistance, mechanical strength and thermal resistance are required. PET is a well-known example of a material system exhibiting a crystallization induced by mechanical loading. The phase transformation, highly dependent on the loading path, loading rate and temperature, implies a modification of the material response at the macroscale. Establishment of the structure-property relationship is a prerequisite for a reliable assessment of the material design whether under in-service or the manufacturing (e.g., hot-drawing near glass transition) process. Over the years, rheologically based models were formulated to reproduce the mechanical response of PET [1][2][3][4][5][6][7][8][9][10][11][12] and other crystallizable thermoplastics [13,14]. From the rheological viewpoint, these models combine different resistances to represent intermolecular and molecular network micromechanisms. The strain stiffening due to the appearance of the newly formed crystalline phase is introduced either implicitly via threshold conditions or explicitly by considering the crystallization kinetics via an Avrami expression. In the latter approach, the presence of crystallites is introduced thanks to the concept of volume fraction in which the polymer is seen as a two-phase composite. The concept was also extended to semicrystalline polymers [15][16][17][18][19]. Nonetheless, the active interactions between the continuous amorphous phase and the discrete crystalline domains are not considered, which constitutes a weakness from the physical viewpoint. Only micromechanical homogenization methods allow to incorporate local interactions and then lead to a physically convenient way to constitutively represent the heterogeneous problem. The latter is the basis of two main approaches to predict the effective mechanical properties of semicrystalline media, with different idealizations of the microstructure to derive the constitutive relations. In a first approach, the material can be defined as an aggregate of layered two-phase composite inclusions, each one represented by an amorphous layer and a crystalline lamella. The approach has been used in a viscoplastic framework [20,21], an elasto-viscoplastic framework [22][23][24][25] or a purely elastic framework [26,27]. In a second approach, the material system is treated as an Eshelby inclusion problem by seeing the crystals as reinforcing ellipsoidal inclusions embedded into a continuous amorphous matrix. The Eshelby-type inclusion approach was mainly used for the initial elastic behavior [26,[28][29][30][31][32][33], the linear viscoelastic behavior [34], the initial yield behavior [35] and the postyield behavior [36].
The main objective of this paper is to explore the relevance of a micromechanics-based elasto-viscoplastic modeling to represent the deformation-induced phase transformation in PET. The model arises from the Eshelby inclusion theory and the micromechanics framework. The intrinsic viscosity of the amorphous phase is considered using the Duvaut-Lions viscoplastic approach, and the crystallization kinetics is governed by an Avrami equation. The model parameters are calibrated using available uniaxial data of a PET stretched under two different strain rates at two different temperatures above glass transition temperature. The model capacities to predict the mechanical response along with the phase transformation are verified by comparison with the available experimental data of a PET stretched under equal biaxial and constant width conditions.
The outline of the present paper is as follows. The micromechanical model is presented in Section 2. Section 3 is devoted to the comparison of the model simulations with experimental observations. Concluding remarks are given in Section 4.

Model
The constitutive representation of the PET system is treated as an Eshelby-type inclusion problem in which the material volume element consists of an amorphous phase as the continuous phase and discrete newly born crystalline domains. The two constitutive phases are supposed to be isotropic and homogeneous media with elastic stiffness tensors C am and C cry .

Micromechanics-Based Theory for Deformation-Induced Phase Transformation
The constitutive relation between the macrostress tensor σ and the elastic part ε e of the macrostrain tensor ε is given by: in which C is the macroscopic elastic stiffness tensor of the semicrystalline material expressed as [37]: where I is the identity tensor, S is the Eshelby tensor and Y is a tensor expressed as: where φ cry is the volume fraction of the strain-induced crystalline phase for which the kinetics is governed by an Avrami-type formula that is specified below. The fourth-order isotropic elastic stiffness tensors, C am and C cry , are expressed, in Cartesian components, as follows: Polymers 2022, 14, 3028 3 of 13 C cry ijkl = E cry 2 1 + ν cry δ ik δ jl + δ il δ jk + 2ν cry 1 − 2ν cry δ ij δ kl (5) in which E am and E cry are the Young's moduli, and ν am and ν cry are the Poisson's ratios. The term δ ij denotes the Kronecker delta symbol. The plastic yielding is considered from the continuum plasticity theory. Regarding the associative plastic flow rule, the macroscopic plastic strain rate . ε p is expressed as: where . λ is the plastic multiplier and F is the macroscopic yield function. Considering the von Mises yield criterion with isotropic plastic hardening for the continuous amorphous phase, the yield function F can be expressed as [38]: in which e p is the macroscopic equivalent plastic strain, σ y is the initial yield strength of the amorphous phase and, h and q are the hardening parameters of the amorphous phase. The term T can be given by: where T 1 and T 2 are provided in Appendix A. The volume fraction of the strain-induced crystalline phase φ cry is given by: in which φ ∞_cry is the maximum crystal degree and κ is the total degree of transformation following the Avrami-type expression [4] modified by Ahzi et al. [5]: where α A is the Avrami exponent, K av is the phase transformation rate function, . ε is the applied strain rate and . ε re f is the reference strain rate. The transformation rate function K av takes the empirical form defined as follows: in which Nu is the number density of nuclei in the amorphous phase. This empirical formula first emerged for the study of spherulitic growth in thermallyinduced crystallization and may not be an optimized choice for all kinetics of newly formed crystals due to differences in morphology and in size. Nonetheless, it was also employed in previous studies in the context of strain-induced crystallization in PET [5,6] and in PLA [13,14] as a phenomenological description of the evolution of a newly formed phase.

Model Implementation
The model allows to relate elasto-viscoplastic macrobehavior to microstructure variations depending on the loading parameters in terms of mechanical path, strain rate and temperature. It is identified using uniaxial (UA) data and its predictability is verified using equal biaxial (EB) and constant width (CW) data. The components of the macrostress tensor σ under a general biaxial stretching are: σ 11 > 0, σ 22 = Rσ 11 and σ ij = 0 for all other components (12) where R = σ 22 /σ 11 is the stress biaxial ratio: UA : EB : CW : The CW condition is a particular biaxial loading in which a stretching is performed in one direction while the transversal one is kept constant. In this regard, the stress biaxial ratio is not constant but changes iteratively during the loading.
The terms T : σ and σ : T : σ in Equation (6) write: The macroscopic plastic strain rate . ε p becomes: The plastic multiplier . λ was computed from the plastic consistency condition: . λ . F = 0, the yield condition being formulated in a Kuhn-Tucker form by: . λ ≥ 0, λ F = 0. The model was coded in Matlab software using the flowchart provided in Figure 1.

Model Implementation
The model allows to relate elasto-viscoplastic macrobehavior to microstructure variations depending on the loading parameters in terms of mechanical path, strain rate and temperature. It is identified using uniaxial (UA) data and its predictability is verified using equal biaxial (EB) and constant width (CW) data.
The components of the macrostress tensor σ under a general biaxial stretching are: for all other components (12) where 22 11 R =  is the stress biaxial ratio: UA: EB: CW: 0 1 01 The CW condition is a particular biaxial loading in which a stretching is performed in one direction while the transversal one is kept constant. In this regard, the stress biaxial ratio is not constant but changes iteratively during the loading.
The terms : T σ and :: σ T σ in Equation (6) write: ( ) The macroscopic plastic strain rate p ε becomes: The plastic multiplier  was computed from the plastic consistency condition:   The model inputs are the material constants, in particular the components of the amorphous elastic stiffness tensor, the components of the crystalline elastic stiffness tensor and the maximum crystal degree. The loading parameters, in terms of strain rate and temperature, are also specified. Both stiffness and yield function are updated iteratively while, once the threshold is reached, the crystallization amount increases with the mechanical loading. The rate-dependency of the yielding behavior is also taken into account by using the Duvaut-Lions viscoplastic approach [39,40]: where η is a viscosity parameter, σ and by integrating the two previous equations using a backward Euler algorithm: where ∆t n+1 is the time step. When ∆t n+1 /η → ∞ , the inviscid solution is recovered, and when ∆t n+1 /η → 0 , the elastic solution is achieved.

Results and Discussion
In what follows, the model is quantitatively compared with experimental results of PET loaded at different loading conditions. The model parameters for the two constitutive phases are listed in Tables 1 and 2.   Figure 2a shows the variation with temperature θ of the amorphous stiffness E am around the glass transition temperature θ g . Like all amorphous polymers, the PET elastic modulus exhibits a plateau within the glassy state, a drastic drop within the glass transition region and a continued decrease within the rubbery region. This behavior is represented by the following function [8]: where E g is the amorphous modulus in the glassy region, E r is the amorphous modulus in the rubbery region, ∆θ is the interval of the temperature range across which the glass transition occurs and X g is the slope outside the glass transition region. The variation with temperature θ of the amorphous Poisson's ratio ν am around θ g is given by the following function: where ν g is the amorphous Poisson's ratio in the glassy region and ν r is the amorphous Poisson's ratio in the rubbery region.

Discrete Crystalline Phase Properties
In the identification exercise, the crystal features, in terms of shape factor and elastic properties, are determined using the overall elastic properties at ambient temperature taken from the work of Cosson et al. [11] on previously crystallized PET. The theoretical and experimental stiffening is shown in Figure 2b as a function of crystal amount. Note that the amorphous elastic properties are considered to be independent of crystal amount.

Overall Response
The other model parameters were identified using the UA data of Salem [41] in terms of stress-strain and crystallization curves. The identification exercise includes the rate effect using two available strain rates: 0.42/s and 2.1/s. Only the straining temperature of 90 °C is used for the identification. Figure 3 presents the model results in comparison with the experimental stress-strain and crystallization curves. The solid lines represent the model results while the symbols designate the experimental data.

Discrete Crystalline Phase Properties
In the identification exercise, the crystal features, in terms of shape factor and elastic properties, are determined using the overall elastic properties at ambient temperature taken from the work of Cosson et al. [11] on previously crystallized PET. The theoretical and experimental stiffening is shown in Figure 2b as a function of crystal amount. Note that the amorphous elastic properties are considered to be independent of crystal amount.

Overall Response
The other model parameters were identified using the UA data of Salem [41] in terms of stress-strain and crystallization curves. The identification exercise includes the rate effect using two available strain rates: 0.42/s and 2.1/s. Only the straining temperature of 90 • C is used for the identification. Figure 3 presents the model results in comparison with the experimental stress-strain and crystallization curves. The solid lines represent the model results while the symbols designate the experimental data. and experimental stiffening is shown in Figure 2b as a function of crystal amount. Note that the amorphous elastic properties are considered to be independent of crystal amount.

Overall Response
The other model parameters were identified using the UA data of Salem [41] in terms of stress-strain and crystallization curves. The identification exercise includes the rate effect using two available strain rates: 0.42/s and 2.1/s. Only the straining temperature of 90 °C is used for the identification. Figure 3 presents the model results in comparison with the experimental stress-strain and crystallization curves. The solid lines represent the model results while the symbols designate the experimental data.  It can be observed that the model is able to reproduce adequately the strain rate dependency on the overall UA response along with the crystallization. The crystallization kinetics and the intrinsic viscosity of the amorphous matrix are the two rate-dependent factors introduced in the model and affecting the overall response. The increase in strain rate slightly affects the initial yield region but significantly influences the progressive strain hardening and the dramatic strain hardening occurring at a strain of approximately 1. Moreover, the model is able to capture the acceleration of both the crystallization onset and the crystallization amount with the strain rate.  [2] are reported in the form of symbols, while the solid lines represent the model predictions. The crystallization prediction accompanies the stressstrain curves. Note that the papers of Buckley et al. [1] and Adams et al. [2] do not present crystallization data, except at a temperature of 86 • C under EB condition.

Comparison between Model and Experiments
A global view at these results shows that the model predictions are favorably compared with the experimental stress-strain curves for the two loading modes. Figures 4a and 5a show that the model adequately predicts the significant effect of the straining temperature on the stress-strain response including stiffness, yield strength and strain hardening region. The predicted crystallization is also presented in Figures 4b and 5b. Inversely to strain rate effects, the increase in loading temperature delays the onset of crystallization, which in turns affects the onset of strain hardening. The higher the loading temperature, the higher the strain level for which the dramatic strain hardening occurs. Figure 6 shows the theoretical and experimental CW stress-strain behavior at two strain rates: 1/s and 4/s. It can be seen that the model is able to capture the highly nonlinear mechanical response, including the dramatic strain hardening occurring at a strain of approximately 1. The higher the strain rate, the higher the onset of crystallization and the amount of crystallization. and the crystallization amount with the strain rate.

Comparison between Model and Experiments
Figures 4-6 present the model predictions under EB and CW loadings using the model parameters identified under UA loading. The data extracted from the works of Buckley et al. [1] and Adams et al. [2] are reported in the form of symbols, while the solid lines represent the model predictions. The crystallization prediction accompanies the stress-strain curves. Note that the papers of Buckley et al. [1] and Adams et al. [2] do not present crystallization data, except at a temperature of 86 °C under EB condition.     A global view at these results shows that the model predictions are favorably compared with the experimental stress-strain curves for the two loading modes. Figures 4a  and 5a show that the model adequately predicts the significant effect of the straining temperature on the stress-strain response including stiffness, yield strength and strain hardening region. The predicted crystallization is also presented in Figures 4b and 5b. Inversely to strain rate effects, the increase in loading temperature delays the onset of crystallization, which in turns affects the onset of strain hardening. The higher the loading

Implication of the Phase Changes
The phase transformation effect on the overall response can be analyzed thanks to the model considering the different loading parameters in terms of loading path, loading rate and loading temperature. When the presence of crystals is neglected in the model (dashed lines in Figures 3a, 4a, 5a and 6a), no dramatic strain hardening occurs and the rate-dependency becomes weaker. The amorphous phase stretching may also play a role in the dramatic strain-hardening region [13]. The actual model introduces the occurrence of the crystallization as the origin of the strain hardening but, acknowledging its effective contribution [42], the supplementary effect of the amorphous phase stretching may be considered by an adequate modification of the amorphous plastic hardening in the yield function [36].

Concluding Remarks
In this work, a micromechanics-based elasto-viscoplastic constitutive model is presented. Formulated within the Eshelby inclusion theory and the micromechanics framework, the constitutive representation considers the material system as a two-phase medium in which the crystal is the reinforcement element of a viscoplastic amorphous medium which increases during stretching. The model provides a quantitative relation between deformation-induced phase transformation and mechanical response under different loading conditions in terms of loading path, loading rate and loading temperature. The model capacities are verified using available experimental observations under uniaxial, equal biaxial and constant width conditions.
Although the model can quite well capture the elasto-viscoplastic response of PET along with the strain-induced crystallization, some improvements are still necessary. It is indeed necessary to extend the mathematical description of the polymer deformation behavior to finite deformation. Furthermore, both isotropic and kinematic plastic hardening deserve to be considered in the local yield function [43] in order to bring a better description of the different steps of the overall mechanical response.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.