Numerical investigation of the fiber/matrix inter-phase damage of a PPS composite considering temperature and cooling liquid ageing

This paper presents a numerical study of fiber/matrix inter-phase damage of a PolyPhenylene Sulfide with 40%wt of short glass fibers (PPS GF40). Homogenization methods are used to estimate the macroscopic mechanical behavior of composites like glass fiber reinforced PPS. Starting from the study of the mechanical behavior of the PPS matrix, an analytical homogenization was performed. Comparison with composite’s experimental curves showed a good prediction of the elastic behavior. However, differences between experimental and predicted plastic behaviors were highlighted. Finite element analyses considering inter-phase damage were done at different temperatures and for several fiber orientations so to explain differences arising between analytical approach and experimental results. This work allowed a study of the evolution of the impact of this damage on mechanical properties as a function of temperature and fiber orientation. Finally, this work showed the development of a weakening of the fiber/matrix interface for a liquid aged composite and to quantify the decrease of the interface properties.


Introduction
Environmental issues are representing one of the main problematic to take into account for new industrial developments. For car manufacturers, one of the main thematic of CO 2 reduction strategy is to reduce weight of structures, especially by replacing metallic parts by plastic composites. The use of short fiber reinforced thermoplastic is thus widely increasing even in under-the-hood applications. Those applications represent technical challenges as the parts are subjected to temperatures up to 130 °C and are fully immersed in cooling liquid when activated. The understanding of the mechanical behavior of those materials is therefore essential in order to design mechanical parts. It is in this context that a first work was conducted in order to build a reliable analytical model, able to predict the evolution of the mechanical properties of a neat PPS as a function of temperature and cooling liquid ageing [1].
The purpose of this study is to extrapolate the model presented in this precedent paper to the composite level, by considering the adjunction of fibers. To do so, an analytical second-order Mori-Tanaka homogenization was first performed, in a revised form, in order to avoid mathematical fluctuations. The obtained behavior was confronted to experimental data, obtained with the same methodology that the one described in the first study [1]. Even if this homogenization scheme is known to be reliable [2], a noticeable difference was observed between plastic parts of the analytical and experimental curves. This difference was also highlighted in studies using other homogenization schemes [3,4]. In order to find the origin of this difference, recent studies tried to explain this phenomenon by the occurring of damage in the material.
The presented work tries to numerically validate this theory with an implementation of a fiber/matrix interphase damaging methodology. A finite elements equivalent homogenization model was thus built with a cohesive inter-phase between fibers and PPS matrix. Performing the analyses at different temperatures and for several preferential fiber orientations showed the evolution of the impact of inter-phase damage on mechanical properties, as a function of those different conditions. Finally, the developed inter-phase damaging methodology allowed a new way to understand the effects of liquid ageing on plastic composites, with an estimation of the induced decrease of fiber/matrix interface properties.
To do so, tensile tests were performed on an aged composite, for several preferential fiber orientations. The results of the experimental campaign are presented in Sect. 2. In order to validate the proposed multi-scale modelling by confrontation with experimental results, an estimation of the fiber orientation in the used samples was first performed and presented in Sect. 3. The development of the analytical homogenization and the confrontation to the experimental data are showed in Sect. 4. Finally, the numerical methods used in order to investigate the fiber/ matrix cohesion are described in Sect. 5.

Experiments
Experimental investigations of the influence of temperature and fibers orientation on the mechanical behavior of the PPS GF40 and its matrix carried out in a precedent study [1] were extended to the ageing consequences on the behavior of the composite. Following ISO 527-2 standard, tensile tests were performed on aged samples with the same methodology and geometry of sample that it was first settled in a previous investigation [1]. The injection of short fiber reinforced composites induces a complex microstructure of the material. Then, one of the objective is to take into account the residual orientation in the samples. Indeed, the flow of the melted reinforced thermoplastic leads to a core-shell microstructure of the samples. Moreover, the geometry of the injected part can lead to side effects which can have consequences of the fiber distribution. The samples were thus extracted from injected plates in order to avoid any side effect and allowing the same good repeatability than the one observed in the reference study [1]. The plates, of size 100*100*3.2 mm, were molded by injection with a temperature of 136 °C and kindly supplied by Solvay. The locations of the specimens on the plates are shown in Fig. 1, in order to prevent end effects such as capillarity at the mold surface and obtain composite samples with different preferential fiber orientations. The matrix being more sensitive to end effects than the composite, matrix samples were cut following the transverse to mold flow direction location, allowing the extraction of three samples per plate. The highest stiffness of the composite is reached when the fibers are mostly aligned to the solicitation direction. The corresponding samples are then less sensitive to the end effects, even more as the eventual end effects are far from the working length L 0 . Consequently, those effects can be neglected. Thus, four samples aligned longitudinally with respect to the molflow direction per plate can be extracted, while only one sample with an orientation of 45° with respect to the moldflow can be extracted per plate. Shaded areas are unused areas.
The samples were then immersed in a cooling liquid composed of 50% of water and 50% of glycol in containers maintained at 90 °C for 400 h in order to reach saturation, as it was defined in the reference study [1]. Before the beginning of the tensile tests, the containers were cooled at ambient temperature and the samples were emerged from the cooling liquid right before the test in order to avoid liquid loss in the sample. The tensile tests were performed with an INSTRON™ 5585H, equipped with mechanical grips and controlled by a 10kN-capacity load cell. The local strain was controlled with a video extensometer. The monotonic tensile tests were strain controlled and for each test, the deformation rate was fixed at 1.10 −3 s −1 as specified in ISO-527. Each configuration was repeated three times. The presented results are the average of the three repetitions, as it was defined in the reference study [1].
The results of the tensile tests are presented in Fig. 2 and confronted to the results obtained for dry-as-molded samples, then defined in following figures as "DAM".
The consequences of cooling liquid immersion are evolving as a function of the fiber orientation. The ageing which the sample undergoes seems to decrease its mechanical properties when solicited in the longitudinal preferential orientation of fibers. However, the more the fibers are misaligned with the loading direction, the less the mechanical properties are affected by this process. Indeed, as shown in Fig. 2a, c slight increase of the Young modulus is even noticeable, traducing an opposite effect of ageing when the material is solicited in this direction than solicited in moldflow direction. This behavior, shown in Fig. 2c, can be explained by the fact that when the fibers are mostly transversally oriented with respect to the loading direction, the behavior of the composite is matrix dominant. The ageing consequence on the mechanical behavior of the composite is thus function of the fiber orientation and shows a different dynamic than the one observed for the neat PPS. Indeed, as it was shown in a precedent study [1], PPS matrix undergoes a physical ageing due to the imposed temperature during the ageing process inducing an increase of its properties, as shown in Fig. 3.
This dynamic is also observed for a transversally oriented composite with respect to the moldflow direction, as its behavior is matrix dominant for such configuration. The objective of the following sections is to take into account these experimental observations in the extrapolation of the model developed in a previous investigation of ageing effects on the mechanical behavior of a PPS [1] to the composite level, considering fiber fraction and fiber orientation.

Definition of the microstructure
The composite behavior is related to its microstructure. It is thus important to evaluate the repartition of fibers in the thickness of the sample, in order to estimate the macroscopic behavior of the composite. An injection simulation of the plates was performed in order to estimate the skin-core microstructure in the material, as presented in Fig. 4.
The simulation was performed with the software Digimat-RP® with the required common injection parameters for the injection of a PPS composite as the mold temperature was fixed at 136 °C. The geometry of the injected plates is the same as the geometry of the real plates supplied by Solvay, and the gate geometry is defined in order to be equivalent to the one used by the supplier. As the plate presents a thickness of 3.2 mm, a fine mesh was applied in order to be able to observe the core-shell microstructure. The injected plate was thus meshed with approximately 3 200 000 tetrahedral elements.
As presented in Fig. 5, it was then possible to estimate the fiber orientation in the samples, whatever the thickness and the position of the sample in the plate. The samples were numerically extracted from the simulated injected plate following the same extracting methodology that the one illustrated in Fig. 1. As it was shown in a previous study [1], the applied methodology of extraction of the samples allows a good repeatability of the tensile tests results for each studied orientation of the sample with respect to the moldflow direction. This can be explained by the fairly homogeneous microstructure in the characteristic working length L 0 in the sample, as shown in Fig. 5. Thus, a representative volume element (RVE) of the microstructure in this area can be extracted from those samples. This representative volume element is extracted at the middle of the working length, corresponding to the center of the plate, which is a location where samples were extracted for every studied orientation, as presented in Fig. 6.
This representative volume element was then used to perform an analytical homogenization of the model proposed in a previous study [1].

Analytical homogenization
The selected homogenization scheme is based on the well-known Mori-Tanaka method [2,5], with a two-steps pseudo-grain decomposition in order to take into account the complexity of the microstructure. As presented in the study of the matrix behavior [1], the PPS presents a hygrothermo-visco-elasto-pseudo-plastic behavior which can be simplified to (Eq. 1): This model is composed of several variables and parameters: • The Young's modulus C * 0 • The hardening exponent m • The hardening modulus R ∞ in MPa • The overall strain ε • The plastic strain p • An Arrhenius-type variable e , allowing the prediction of elastic behaviour through temperature  • An Arrhenius-type variable p , allowing the prediction of pseudo-viscoplastic behaviour through temperature • The reference temperature T 0 • The actual temperature T eq • The stress relaxation time τ In order to consider the non-linearity induced by viscosity and the pseudo-plasticity of the PPS, an incremental linearization method [3] was used on this model. In order to avoid fluctuations, the used version of Mori-Tanaka is the second-order one, in a revised form. This homogenization was carried out with Digimat-MF® software [5].
The material parameters of the matrix and of the reinforcement are listed in Table 1, for a temperature of solicitation fixed at 70 °C. Where The matrix model allowing the prediction of the evolution of its mechanical properties as a function of temperature and ageing effects, correlations were performed at different temperatures and on an aged material, for several preferential fiber orientations, as described in Fig. 7. The predictions, presented as "M-T" in Fig. 7, are confronted to experimental data obtained in a precedent study using the same experimental protocol [1], and describing the influence of temperature on mechanical properties on the PPS GF40 and also confronted to experimental data describing the influence of ageing process on this behavior presented in Sect. 2.
As shown in Fig. 7a, b, the homogenization prediction shows a very good correlation on the elastic part, whatever the temperature and the preferential fiber orientation on a DAM material. However, a slight difference between the analytical prediction and the experimental data is observed for the plastic part. This difference is even more pronounced when the behavior shows a major pseudoplastic behavior, i.e. when the test is performed at high temperature. This difference was already observed in several studies, even with FEA based homogenization [4,6].
Recently, some studies tried to explain this observed difference by dealing with the hypothesis of the existence of a fiber/matrix inter-phase damage occurring for monotonic loadings thanks to different methods. Some methods tried to take this effect into account with analytical methods [7][8][9][10][11][12], with the study of the fiber/matrix interface [13][14][15][16][17][18][19][20][21][22] or with cohesive inter-phases surrounding the fibers [22][23][24][25][26][27]. Indeed, SEM observations performed in precedent studies showed that micro-crack at fiber end, plasticization of the matrix and/or fiber/matrix debonding can occur with a quasi-static monotonic solicitation [28][29][30][31][32][33]. This degradation of the fiber/matrix cohesion can lead to a loss of reinforcement effect and so, to a decrease of the effective mechanical properties of the composite. The difference between the experimental curve and the analytical prediction observed at 70 °C can also be due   to the modelling process. A modelling cannot be perfect and a loss of information is possible. Thus, there is a loss of information linked to the injection simulation, another one from the linearization method and another from the homogenization scheme. However, a change of dynamic is clearly noticeable on an aged composite as a function of the preferential fiber orientation. Indeed, even if the behavior of the composite with mostly transversally oriented fibers is well predicted, the analytical homogenization does not allow the prediction of the mechanical behavior of the aged composite with a longitudinal preferential orientation of fibers, as shown in Fig. 7c. As the matrix is stiffer once aged than when dry-as-molded, the homogenization scheme predicts a macroscopic behavior of the composite stiffer than the behavior obtained for a dry-as-molded material. The consideration of a loss of information linked to the simulation workflow only cannot explain this change of dynamic.

Moldflow direcƟon
The only difference between the PPS GF40 and its matrix is the presence of glass fibers in the composite. As glass is supposed to be insensitive to ageing effects, a way to explain this duality of behaviors is to focus on the fiber/ matrix interface degradation due to liquid ageing, as it was already experimentally observed in other studies [28][29][30][31][32][33].
The aim of the next section is thus to build a methodology of inter-phase damage simulation with finite elements analysis in order to quantify this damage under monotonic tests and study its evolution for different configurations of temperature and fiber orientation.

Geometry
The first step was to build a finite element model, equivalent to the analytical prediction. The representative volume element (RVE) was created with the help of the injection simulation described in Sect. 3. This RVE is thus the same that the one used for the analytical homogenization. The implementation of inter-phase damage requires the introduction of another material associated to a damage law. Thus, a thin inter-phase surrounding all the fibers was modelled between the fibers and the matrix with a 5% relative thickness. The geometry was generated with the help of the software Digimat-FE®, as well as the meshing process and the simulation, and is presented in Fig. 8.
The aspect ratio of fibers and the material parameters are the same that the ones used to perform the analytical homogenization.
A fine-level meshing was chosen for the inter-phase, as this area is subjected to a damage law, and a mediumlevel meshing was chosen for the fibers and the matrix, as presented in Fig. 9.

Boundary conditions
In order to simulate a tensile test, the simulation is strain controlled at a quasi-static speed, fixed at 1.10 −3 s −1 . Then, a periodic boundary condition was applied to unsure that the flux of the field variable is periodic with respect to the faces of the volume element [4], as presented in Fig. 10.
This type of boundary condition is also known to lead to best predictions when compared to the Dirichlet and mixed boundary condition type [4].
The tetrahedral mesh was defined in order to find a compromise between calculation time and accuracy of the results. As the RVE used for the numerical simulation is the same that the one used for analytical homogenization, the tetrahedral mesh was validated by a confrontation to the analytical result. As shown in Fig. 14, for a RVE without any damage law, the "FEA" result is equivalent to the "M-T" prediction, validating the mesh applied on the matrix and on the fibers and also validating the applied boundary conditions.

Damage implementation and model validation
The modeled inter-phase area is defined by a cohesive behavior, with a damaging evolution following a traction-separation law, described in Fig. 11. The damage initiation criterion used is defined in (Eq. 3).
where n , s and t are respectively the normal, first shear and second shear directions.
In order to be as representative as possible of the real microstructure of the composite, the cohesive zone parameters were set as far as possible on the mechanical parameters of the PPS matrix. As the matrix behavior assumes a plastic part in its behavior, the initial damage criterion was fixed at the maximum mechanical stress of the matrix, which can be predicted at any temperature by the model developed in a previous investigation [1], in order to consider an inter-phase which starts to damage when the matrix fails in this area. This is a hypothesis    made by the authors, considering the additives present in the PPS Ryton R-4-200 in order to strengthen the fiber/matric cohesion and described by the supplier. The traction-separation law used has a linear shape. Its evolution is based on a normalized displacement and its parameters are listed in Table 2, for a temperature fixed at 70 °C.
The normalized displacement at failure was estimated with respect to the inter-phase thickness. The influence of the normalized displacement at failure parameter δ on the mechanical properties of the composite is presented in Fig. 12. Figure 12 clearly shows that, for a normalized displacement lower than 0.1, the maximum stress is the most affected. However, as shown in Fig. 7b, the ultimate tensile stress of the composite solicited in the moldflow direction shows an ultimate tensile stress corresponding to a normalized displacement at failure superior to 0,1. Therefore, in order to take into account the experimental results, the normalized displacement at failure was fixed at 0.1.
As it was also shown in this previous study [1], the higher the temperature is the more viscosity affects the mechanical behavior. Thus, the more the temperature rises, the more the ultimate strain of the material rises. This phenomenon also has an influence on the normalized displacement at failure δ. The evolution of this parameter with respect to temperature is presented in Fig. 13.
A simulation was then performed at 70 °C, for a longitudinal preferential orientation of fibers, in order to validate the implemented method. Figure 14 shows that, without the definition of the traction-separation law, the finite element analysis results are equivalent to the ones given by analytical homogenization, validating the FE equivalent model. As also shown in Fig. 14, a good correlation between the finite elements analysis with inter-phase damage and the experimental curve is now observed. This result is in line with the hypothesis that inter-phase damage occurs even in monotonic tensile tests [28][29][30][31][32][33].
The influence of the interface thickness on the mechanical properties of the composite solicited in the moldflow direction was finally studied. The results, presented in Fig. 15, showed that this relative thickness induce a relatively small dispersion of the results. The relative thickness of the cohesive zone was thus arbitrary fixed to 5% of the thickness of the fiber, in order to reduce the number of elements and to avoid convergence issues. As the method was validated, a study of the evolution of inter-phase damage at different temperatures and for different fiber orientations was then performed with the defined mesh, which was generated with approximately 1 300 000 nodes in order to provide accurate results.

Identification of damaging system
Analysis of the FE results allows the identification of the leading damaging mechanism. As shown in Fig. 16, the maximum damage occurs at the tip of the fiber, in the inter-phase.
This clearly shows that the major damaging mechanism involved in the decrease of the mechanical properties is the micro-cracks at fiber ends, which was already experimentally observed in precedent works [28][29][30][31][32][33].
The geometry of the RVE and the geometry of the coating were thus validated, as well as the implemented material properties for the matrix, fibers and cohesive material. This section also highlighted and confirmed the nature of the principal damaging mechanism inducing a decrease of the mechanical properties, in accordance with SEM observations, allowing other simulations for other configurations.

Damage evolution
Other simulations were then performed at 50 °C and 70 °C, for a 0°, a 45° as well as a 90° oriented sample with respect to the moldflow direction, in order to study the evolution of inter-phase damaging impact on the mechanical properties of the PPS GF40 as a function of temperature and fiber orientation. The simulations for different preferential fiber orientations are performed with a rotation of the RVE with respect to the loading direction. The simulations showed that the inter-phase damage occurs when the macroscopic plasticity of composite starts, as presented in Fig. 17. The exponential evolution of the damage indicates that the more important the plastic strain is, the more weakened the fiber/matrix cohesion is.
As shown in Fig. 18, the developed methodology provides good predictions for other conditions. The comparison between Figs. 7b, 18 clearly shows that the implemented damaging process corrects the plastic part of the behavior for different configurations by following the exponential law showed in Fig. 17, validating the methodology.
The normalized impact of the inter-phase damage on the mechanical properties is resumed in Fig. 19.
The evolution of the inter-phase damage seems to increase following an exponential law described in (Eq. 4) as the temperature increases.
where T is the temperature in (°C) and a and b are coefficients.
= ae bT Figure 19a shows that the damaging impact increases exponentially as a function of temperature, following (Eq. 4). Figure 19b shows that the damaging impact on mechanical properties seems to decrease linearly the more the fibers are misaligned with respect to the loading direction.
The fact that the damaging sensitivity is more pronounced at higher temperatures can be explained by the decrease of the matrix properties as a function of temperature. This leads to a damaging process which starts sooner, as presented in Fig. 20a. Moreover, as previously shown in Fig. 17, damage process starts when the deformation become plastic. As it was shown in a previous work [1], the more the temperature increase, the more the inelastic behavior of the PPS becomes important, explaining the interphase damage to be more pronounced and the dynamic of evolution shown in Fig. 19a.
It is also interesting to note that, as the principal damaging system is the micro-crack at fiber-end, the damage as an influence on the mechanical behavior of the composite which can be detected sooner for a material solicited in the moldflow direction as shown in Fig. 20b. The appearing time is also a normalized value. This result must be relativized as the maximum strain reached in the tensile tests evolves with respect to the testing conditions. Indeed, as the maximum strain is the smallest in the longitudinal direction, the damage artificially appears sooner.

Ageing consequence on inter-phase strength
The numerical methodology developed in the precedent sections was then used in order to explain the observed decrease of the mechanical properties of the composite once aged. For a longitudinal solicitation, the experimental aged tensile curve is reached by the FE model with an inter-phase behavior considering 50% of the stiffness and the maximum stress of the matrix (and so 50% of the initial damage initiation criterion for a DAM material), as presented in Fig. 21.
In order to validate the decrease of 50% of the strength of the inter-phase, simulations were performed for other fiber orientations with exactly the same material parameters.
The predictions, presented in Fig. 22, show very good correlations. The accuracy obtained for the extrapolations for other fiber orientations tends to confirm the implemented rate of degradation.
The prediction presented in Fig. 22a for a preferential orientation of 45° of fibers shows a behavior which is quite less sensitive to fiber/matrix partial debonding than the one obtained for a longitudinal preferential orientation of fibers, which agrees with the dynamic observed in Fig. 19. Figure 22b also shows good prediction capabilities of the numerical model, even for an orientation less impacted by damaging effects than others. Indeed, as shown in Fig. 19, an ambient temperature and a transversal orientation of the sample with respect to the molflow direction represent the configuration which the partial fiber/matrix debonding has the least influence on the macroscopic behavior of the composite. This result perfectly agrees with the experimental observations and analytical predictions presented in the precedent sections. Indeed, as shown in Fig. 2, when oriented transversally with the loading direction, the composite seems to be affected by ageing in the same way than the matrix, as its Young modulus shows a slight increase. Thus, the Mori-Tanaka prediction showed in Fig. 7 is reliable for this preferential orientation of fibers.
The good predictions of the numerical model and the fact that it reports correctly the influence of the fibers orientation on the influence of the damage of the interphase on the mechanical behavior of the composite tends to confirm the hypothesis that the immersion in cooling liquid maintained at 90 °C for 400 h halves the strength of the fiber/matrix inter-phase, leading to the observed decrease of the macroscopic behavior of the PPS GF40. This result is perfectly in line with the precedent experimental SEM observations showing that ageing effects decrease the fiber/matrix interface cohesion and enhance the damage mechanisms and matrix plasticization [28][29][30][31][32][33]. This study completes those observations by allowing a quantification of the decrease of the interface strength due to such ageing process.
A synthesis of damaging impact on the mechanical behavior is presented in Fig. 23 and clearly shows that the inter-phase damage of an aged material can be more than twice more important than the maximum damaging coefficient obtained for a Dry-As-Molded material. This influence is less sensitive for a transverse orientation of the sample, as the reinforcement is less pronounced, and shows the same dynamic that the influence of temperature on the composites mechanical properties studied in a precedent work on the neat PPS [1].

Conclusion
The aim of the study was to develop a numerical methodology of fiber/matrix interface damage in order to understand the difference of behaviors observed for the plastic parts of the tensile curves obtained experimentally and analytically. The analytical homogenization was performed with a second-order Mori-Tanaka method, in a revised form, in order to avoid mathematical fluctuations. The numerical process was validated thanks to correlations on tensile curves and simulations were performed at other temperatures and for several preferential fiber orientations. The principal damaging system was clearly identified as micro-crack at fiber-end and a comparative study of the damaging impact on the mechanical properties was then performed with respect to the different solicitation conditions. This comparative study showed that this impact increases exponentially with respect to temperature and linearly with respect to preferential fiber orientation.
This work also allows a quantification of the consequences of liquid ageing on composites behavior. Indeed, this study shows that inter-phase damage drastically increases for a cooling liquid aged sample, as cooling liquid immersion decreases the fiber/matrix interphase properties by approximately 50%.
The estimation of the decrease of fiber/matrix interphase is fully based on successive modelling methods. The interest of this workflow is thus also its weakness. Indeed, as this estimation is based on modelling considering injection simulation, analytical homogenization and finite elements methods, mathematical fluctuations and loss of information can happen at each step of the modelling process. The error observed between the analytical prediction and experimental results at 70 °C can thus be explained partially or fully by mathematical fluctuation of the homogenization method or by inherent loss of information due to the precision of the injection simulation. Microtomography could be performed on the used samples in order to get the real core-shell microstructure, allowing an improvement of the homogenization results. An improved homogenization considering real RVE would lead to a refinement of the traction-separation law parameters and thus a refinement

Compliance with ethical standards
Conflict of interest On behalf of all the authors, the corresponding author states that there is no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.