Suitability of Constitutive Models of the Structural Concrete Codes When Applied to Polyolefin Fibre Reinforced Concrete

The use of fibres as structural reinforcement in concrete is included in standards, providing guidelines to reproduce their behaviour, which have been proven adequate when steel fibres are used. Nevertheless, in recent years new materials, such as polyolefin fibres, have undergone significant development as concrete reinforcement. This work gives insight on how suitable the constitutive models proposed by the Model Code 2010 (MC2010) are in the case of such polymer fibres. A set of numerical models has been carried out to reproduce the material behaviour proposed by the MC2010 and the approach based on the softening function proposed by the authors. The results show remarkable differences between the experimental results and the numerical simulations when the constitutive models described in the MC2010 are employed for different polyolefin fibre reinforced concrete mixes, while the material behaviour can be reproduced with greater accuracy if the softening function proposed by the authors is employed when this type of macro-polymer fibres is used. Moreover, the relatively complex behaviour of polyolefin fibre reinforced concrete may be reproduced by using such constitutive model.


Introduction
Concrete has been one of the most employed construction material [1] in modern societies in spite of its impact on the environment [2]. Such use has been based on its resistance to natural weather conditions [3,4] and its remarkable compressive strength and modulus of elasticity [5,6]. However, plain concrete is also known as a quasi-brittle material which is not capable of sustaining high-tensile stresses. Adding reinforcing bars of concrete to form reinforced concrete was one of the possible solutions to this problem. Another possibility is adding fibres randomly distributed while mixing, forming fibre reinforced concrete (FRC). The effect of the addition of fibres in the mechanical properties of concrete depends on numerous parameters, such as the type of fibres used, dosage, shape and size, distribution, orientation and, among others, pull-out response [7,8]. Adding certain types of fibres in a determined volume fraction to concrete a material suitable for numerous applications can be obtained. For instance, steel fibre reinforced concrete (SFRC) has been successfully applied to pavements [9,10], tunnel linings [11,12], or even seismic structures [13,14]. However, the use of SFRC might have two main drawbacks, its durability under certain circumstances [15,16] and its environmental impact [17]. One possibility to overcome such durability issues is to use new polymer macro-structural fibres, which are chemically stable when applied in foundations, tunnel linings or marine environments [18][19][20][21]. Polyolefin fibre reinforced concrete (PFRC) is one of the materials available for these uses.
All the previous applications have been possible due to several national codes and recommendations [22][23][24][25][26][27]. These documents have set the conditions that a certain FRC should meet so that the improvements provided by the fibres added might enable a reduction in the amount of steel reinforcing bars used. Such requirements take into account the residual strength of the FRC studied at certain crack openings which are related to the service limit state (SLS) and the ultimate limit state (ULS). For checking if a FRC formulation boasts sufficient residual strength, bending tests are performed on FRC specimens in laboratory conditions. Once such tests have been performed, and if the FRC properties exceed the requirements set, the contribution of the fibres might be included in the structural design of the concrete element. Although it is true that both SFRC and PFRC have been shown as suitable for their application in structural concrete elements [7,28,29] it should not be overlooked that the extensive experience related with SFRC served as reference when the requirements of the material were set. Once such requirements have been checked, the aforementioned recommendations propose several constitutive models that should be applied to the structural design and that have been successfully applied when using SFRC. Similarly, the constitutive models proposed used several characteristics that were inspired by the mechanical behaviour of SFRC. Consequently, the accuracy of such approaches when dealing with PFRC is a matter that deserves being studied.
The behaviour of PFRC can meet the requirements in the standards [30] although the shapes of the fracture curves have different features compared with those of SFRC [29]. According to the Model Code 2010 (MC2010) [31], the first condition that FRC should meet is a minimum load value, proportional to the maximum strength in the limit of proportionality, at a crack mouth opening displacement of 0.5 mm; this limitation is intended to prevent brittle behaviour of the material. Although PFRC with an average fibre dosage might meet the cited requirement, polyolefin fibres need greater deformations to absorb the load that concrete is not capable of sustaining being the first unloading branch steeper than in SFRC. The second condition that MC2010 establishes is that at a crack mouth opening displacement of 2.5 mm the load must be over 20% of the strength at the limit of proportionality. This limitation for PFRC is met even for very low fibre dosages [30]. As there are remarkable differences in the post peak mechanical behaviours of SFRC and PFRC, the applicability of the constitutive models, which are built based on certain values of their residual strength, should be studied. Thus, this study seeks to address the main issues concerning the considerations in the standards regarding the mechanical behaviour of PFRC.
The main goals of this contribution are the following. First, an implementation of the linear and rigid-plastic models developed in MC2010 in a numerical code by means of a material subroutine is sought. Second, the differences between the material behaviour obtained with the constitutive models of MC2010 and those obtained using the trilinear approach proposed by the authors that delivers an accurate reproduction of the PFRC behaviour are assessed. This matter will be carried out in four mixes of PFRC with 3, 4.5, 6, and 10 kg/m 3 of polyolefin fibres. Lastly, an evaluation of the material behaviour at SLS and ULS and another regarding the fracture energy consumed in the tests are carried out for determining the reliability of using the studied approaches. Summarising, this manuscript supplies additional information for a possible future adaptation of MC2010, and other national concrete structural codes based on it when applied to PFRC.

Post-Peak Material Modelling
In this section, the numerical modelling of fracture used in the simulations is addressed. In the first part, the embedded cohesive crack model is briefly presented. Since the reader can find a detailed description in previous works published by the authors, this description has been kept to a minimum. In the second part, the model is adapted for PFRC following the guidelines of the MC2010.

Embedded Cohesive Crack Modeling
Fracture is modelled by means of a cohesive crack model, based on the work proposed by Hillerborg [32] and implemented as an embedded cohesive crack developed in a tri- angular finite element. This numerical modelling has been used for simulating fracture of concrete [33] and has been adapted for use with brick masonry [34] and, more recently, with PFRC [35].
In this model, fracture is assumed to develop mainly under mode I and with the cohesive stress vector t being parallel to the displacement vector w (thus, it is a central force model), as expressed by (1): with f w being the softening function, expressed in terms ofw, which corresponds to an equivalent crack opening that considers the maximum historical opening, necessary for correctly reproducing possible loading-unloading situations.
This formulation is implemented for triangular elements and, as cracks can only develop in lines parallel to the sides of the element, only three cracking orientations are allowed (see Figure 1). It follows a strong discontinuity approach, dividing the element into two parts, A + and A − , with t being parallel to the crack displacement and constant along the length of the crack L. It is obtained by (2): with A being the area of the element, h the height of the triangle over the side opposite to the solitary node and n the unit vector normal to the crack (see Figure 1).

Embedded Cohesive Crack Modeling
Fracture is modelled by means of a cohesive crack model, based on posed by Hillerborg [32] and implemented as an embedded cohesive crack a triangular finite element. This numerical modelling has been used for si ture of concrete [33] and has been adapted for use with brick masonry [ recently, with PFRC [35].
In this model, fracture is assumed to develop mainly under mode I an hesive stress vector being parallel to the displacement vector (thus, force model), as expressed by (1): with ( ) being the softening function, expressed in terms of , which c an equivalent crack opening that considers the maximum historical open for correctly reproducing possible loading-unloading situations. This formulation is implemented for triangular elements and, as crack velop in lines parallel to the sides of the element, only three cracking orien lowed (see Figure 1). It follows a strong discontinuity approach, dividing th two parts, A + and A − , with being parallel to the crack displacement and the length of the crack . It is obtained by (2)  The material around the crack is elastic and the stress tensor is obtain strain vector that results from subtracting the strain due to the crack disp the apparent strain , obtained under elastic considerations and by using placements (Equation (3)): where stands for the elastic tangent tensor, for the apparent strain v with the nodal displacements, for the gradient vector of the shape func itary node, which can be obtained as: Superscript S makes reference to the symmetric part of the tensor, the dot product (( : ) = ) and ⨂ the usual direct product (( ⨂ ) The material around the crack is elastic and the stress tensor is obtained by using a strain vector that results from subtracting the strain due to the crack displacement from the apparent strain ε a , obtained under elastic considerations and by using the nodal displacements (Equation (3)): where E stands for the elastic tangent tensor, ε a for the apparent strain vector obtained with the nodal displacements, b + for the gradient vector of the shape function of the solitary node, which can be obtained as: Superscript S makes reference to the symmetric part of the tensor, the usual double-dot product ((A : B) ij = A ijkl b kl ) and the usual direct product ((a b) ij = a i b j ). Lastly, by using the expression of the stress tensor (3), the stress vector can be obtained as: This can be rewritten as: where 1 identifies the second-order identity tensor. By using this expression, the crack displacement can be solved for a certain set of nodal displacements by means of the Newton-Raphson method. This model has been implemented by using a user subroutine UMAT of ABAQUS ® and, since it needs information about the geometry of the model (the nodes coordinates of each element), it reads a previously prepared external file.

Model Adopted by the MC2010
The MC2010 provides two possible constitutive laws to reproduce the behaviour of FRC, with one describing a rigid-plastic behaviour and another one a linear behaviour. These models are defined by parameters that make reference to the values of the residual flexural tensile strength and the crack mouth opening displacement (CMOD) at specific instants during a three-point bending test. In particular, the expression of the residual flexural tensile strength for a specific value of the CMOD may be obtained by expression (7): where f R,j is the residual flexural tensile strength in MPa corresponding to a value of the CMOD = CMOD j in mm, F j is the load in N applied when CMOD = CMOD j (in mm), l is the span length, b the specimen width and h sp the height of the notched section (ligament) where the crack propagates, with all these dimensions being in mm (see Figure 2b). Consequently, f R1 stands for the residual strength at a CMOD of 0.5 mm and f R3 stands for the residual strength at a CMOD of 2.5 mm. Figure 2a shows a typical load-CMOD curve for FRC where some values of F and CMOD are identified according to the previous description, as described in the MC2010. Expression (7) can be easily obtained by applying the basic principles of the strength of materials at the notched section of the specimen. As mentioned before, MC2010 proposes two possible constitutive laws for FRC. The rigid-plastic model simply describes a constant stress behaviour from the point at which cracking begins (w = 0) until it reaches a limit value (w u ). In this model, the constant stress is obtained by considering the ultimate state when CMOD = CMOD 3 and the ultimate value of w is considered as w u = CMOD 3 .
The linear model is defined by two reference values, f Fts and f Ftu , which can be obtained by (9) and (10): These values define a line which, as can be observed in these equations, is based on the experimental values for and . Neverheless, the ultimate crack opening wu is not specifically defined and is kept as the maximum crack opening accepted in structural design, which depends on the ductility required. Figure 3 shows the schemes for both models: the rigid-plastic and the linear model. It should be noted that while in this figure the linear model shows a load-descending behaviour, depending on the performance of the PFRC, the linear model may present a positive slope which in the MC2010 is referred to as post-crack hardening behaviour.

Trilinear Model Proposed for PFRC
In past papers published by the authors, the behaviour of hardened specimens of self-compacting concrete (SCC) and vibrated conventional concrete (VCC) was studied, with particular examination of how the increase of fibre dosages affects the flexural behaviour of specimens under three-point bending tests [35]. The authors proposed the use of a trilinear softening curve able to reproduce the classical three-point bending tests that develops under pure mode I conditions, as well as a modified version of the test proposed in [36], where cracking is initiated under a combination of modes (I) and (II) [37].
This trilinear model can be described by three linear parts, as can be seen in Figure 4. The first linear part reproduces the load decrease due to the failure of the concrete matrix up to the instant at which fibres make a relevant contribution to the specimen flexural strength. The second linear part reproduces the contribution of the fibres up to their These values define a line which, as can be observed in these equations, is based on the experimental values for CMOD 1 and CMOD 3 . Neverheless, the ultimate crack opening w u is not specifically defined and is kept as the maximum crack opening accepted in structural design, which depends on the ductility required. Figure 3 shows the schemes for both models: the rigid-plastic and the linear model. It should be noted that while in this figure the linear model shows a load-descending behaviour, depending on the performance of the PFRC, the linear model may present a positive slope which in the MC2010 is referred to as post-crack hardening behaviour. These values define a line which, as can be observed in these equations, is b the experimental values for and . Neverheless, the ultimate crack wu is not specifically defined and is kept as the maximum crack opening accepted tural design, which depends on the ductility required. Figure 3 shows the schemes for both models: the rigid-plastic and the linea It should be noted that while in this figure the linear model shows a load-descen haviour, depending on the performance of the PFRC, the linear model may prese itive slope which in the MC2010 is referred to as post-crack hardening behaviour

Trilinear Model Proposed for PFRC
In past papers published by the authors, the behaviour of hardened speci self-compacting concrete (SCC) and vibrated conventional concrete (VCC) was with particular examination of how the increase of fibre dosages affects the flex haviour of specimens under three-point bending tests [35]. The authors proposed of a trilinear softening curve able to reproduce the classical three-point bending t develops under pure mode I conditions, as well as a modified version of the test p in [36], where cracking is initiated under a combination of modes (I) and (II) [37] This trilinear model can be described by three linear parts, as can be seen in F The first linear part reproduces the load decrease due to the failure of the concret up to the instant at which fibres make a relevant contribution to the specimen strength. The second linear part reproduces the contribution of the fibres up

Trilinear Model Proposed for PFRC
In past papers published by the authors, the behaviour of hardened specimens of selfcompacting concrete (SCC) and vibrated conventional concrete (VCC) was studied, with particular examination of how the increase of fibre dosages affects the flexural behaviour of specimens under three-point bending tests [35]. The authors proposed the use of a trilinear softening curve able to reproduce the classical three-point bending tests that develops under pure mode I conditions, as well as a modified version of the test proposed in [36], where cracking is initiated under a combination of modes (I) and (II) [37].
This trilinear model can be described by three linear parts, as can be seen in Figure 4. The first linear part reproduces the load decrease due to the failure of the concrete matrix up to the instant at which fibres make a relevant contribution to the specimen flexural strength. The second linear part reproduces the contribution of the fibres up to their maximum tensile capacity and in the third part the load decreases due to the final failure of fibres, which can be due to fibre tensile failure or because of the slip of fibres inside the concrete. rials 2022, 15, x FOR PEER REVIEW maximum tensile capacity and in the third part the load decrease of fibres, which can be due to fibre tensile failure or because of th concrete.

Stress-Crack Opening Diagrams Used in This Study
In order to compare the suitability of each of the three appr this study examines the experimental behaviour of four PFRC m predictions of each of them.
The reference experimental results correspond to four PFRC fibre dosages: 3, 4.5, 6 and 10 kg/m 3 . These are selected results mental campaign that studied several PFRC mixes and analysed ties through a series of two prismatic and nine cylindrical specim tails regarding these experimental works may be consulted in [2 Regarding the constitutive laws defined by the MC2010, wh before, Table 1 shows the parameters for each fibre dosage and f bilities described in the standard: rigid-plastic model and linear ear model is defined by two points: one for w = 0 (defined by) by the value of w = CMOD . These diagrams are shown in Figure  following the MC2010, thus a limit value of the crack opening ployed and Equations (8)-(10) used with the values of fR1 and fR3 e which can be consulted in [28,29].

Stress-Crack Opening Diagrams Used in This Study
In order to compare the suitability of each of the three approaches described above, this study examines the experimental behaviour of four PFRC mixes with the numerical predictions of each of them.
The reference experimental results correspond to four PFRC mixes of VCC with four fibre dosages: 3, 4.5, 6 and 10 kg/m 3 . These are selected results from a broader experimental campaign that studied several PFRC mixes and analysed their mechanical properties through a series of two prismatic and nine cylindrical specimens of each mix; the details regarding these experimental works may be consulted in [28,29].
Regarding the constitutive laws defined by the MC2010, which have been described before, Table 1 shows the parameters for each fibre dosage and for each of the two possibilities described in the standard: rigid-plastic model and linear model. Note that the linear model is defined by two points: one for w = 0 (defined by) and another one defined by the value of w = CMOD 3 . These diagrams are shown in Figure 5 and have been defined following the MC2010, thus a limit value of the crack opening of 2.5 mm has been employed and Equations (8)-(10) used with the values of f R1 and f R3 experimentally obtained, which can be consulted in [28,29].

Numerical Simulations
Numerical simulations were carried out by means of the finite element method with the commercial software ABAQUS ® (version 6.18). The softening behaviour was reproduced by the embedded cohesive crack formulation described above, implemented through a user subroutine UMAT in ABAQUS ® . Since, as required by the mathematical formulation described before, only linear triangular elements have been used, and the mesh is refined along the line where cracking develops (see Figure 6a); this refinement was designed based on previous works (see [37] and [38]) where the mesh dependence was already analysed.
In the case of models where fracture is reproduced by a trilinear model and by the linear model of the MC2010 with softening (negative slope of the diagram), all the elements have been modelled by using the same cohesive formulation as that presented before. In the case of those constitutive models with abrupt decrease of strength, like the rigid-plastic model and that termed "linear model with limited development", convergence issues are usual when employing an implicit integration scheme. The convergence problem in the models that uses "linear softening with unlimited development" with hardening (positive slope of the diagram) provides a different explanation. In these cases, As for the trilinear model, Table 2 shows the coordinates of the points that define them for each fibre dosage. The inverse analysis performed to obtain these values may be consulted in [35]. Table 2. Coordinates of the points that define the trilinear model for each fibre dosage according to the point nomenclature provided in Figure 3.

Numerical Simulations
Numerical simulations were carried out by means of the finite element method with the commercial software ABAQUS ® (version 6.18). The softening behaviour was reproduced by the embedded cohesive crack formulation described above, implemented through a user subroutine UMAT in ABAQUS ® . Since, as required by the mathematical formulation described before, only linear triangular elements have been used, and the mesh is refined along the line where cracking develops (see Figure 6a); this refinement was designed based on previous works (see [37,38]) where the mesh dependence was already analysed.
gressing in the element that is already damaged). This leads to an increasing wider region of the specimen that develops, with it being more difficult for the model to find the crack path since, as load displacement increases, there are more possible elements to develop damage which means that convergence is not reached. Figure 6b shows one of these models where damage started developing along the expected crack path but, at some point started to develop around it and finally convergence was not reached when there were many possible elements to develop damage. In such cases with convergence problems, the path of fracture was tracked by defining elements that could be damaged only in the crack path, with the rest of elements behaving elastically. Figure 7 shows the load-deflection curves for each of the approaches described for various fibre dosages: 3, 4.5, 6 and 10 kg/m³, therefore reproducing the experimental results of [28,29].

Results
It is interesting to observe the differences among these approaches. In the case of the rigid-plastic model, diagrams always show a plateau and a subsequent load decrease that describes an exponential-like curve with small load jumps. As these drops are due to the abrupt load decrease described by the rectangular diagram, when an element reaches , there is a sudden load decrease. In the case of the "linear model with unlimited development", the load-deflection diagrams exhibit a linear behaviour after the concrete matrix fails and, depending on the fibre dosage the specimen, may show a softening behaviour (for 3 and 4.5 kg/m 3 ) or a hardening behaviour (for 6 and 10 kg/m 3 ). Obviously, in the case of diagrams with hardening, this behaviour cannot be realistic, since the load may increase without a limit. Nevertheless, it is interesting to compare these results with those of the linear model with limited development that display an analogous behaviour until the most damaged element reaches , when the load-deflection diagrams show a highly similar behaviour to that observed with the rigid-plastic model. In all cases, the load drop follows an exponential-like shape with small jumps that can be explained as for the rigid-plastic model.
Lastly, if the results obtained with a trilinear model are observed, remarkable differences can be highlighted compared with any of the results obtained with the other approaches. First, the maximum load reached that identifies the beginning of the fracture process is higher than in any of the other cases and the shape shows an initial load decay In the case of models where fracture is reproduced by a trilinear model and by the linear model of the MC2010 with softening (negative slope of the diagram), all the elements have been modelled by using the same cohesive formulation as that presented before. In the case of those constitutive models with abrupt decrease of strength, like the rigid-plastic model and that termed "linear model with limited development", convergence issues are usual when employing an implicit integration scheme. The convergence problem in the models that uses "linear softening with unlimited development" with hardening (positive slope of the diagram) provides a different explanation. In these cases, the material becomes harder with increasing values of w, (if an element is damaged it becomes harder, meaning that damage is more likely to develop around it instead of progressing in the element that is already damaged). This leads to an increasing wider region of the specimen that develops, with it being more difficult for the model to find the crack path since, as load displacement increases, there are more possible elements to develop damage which means that convergence is not reached. Figure 6b shows one of these models where damage started developing along the expected crack path but, at some point, started to develop around it and finally convergence was not reached when there were many possible elements to develop damage.
In such cases with convergence problems, the path of fracture was tracked by defining elements that could be damaged only in the crack path, with the rest of elements behaving elastically. Figure 7 shows the load-deflection curves for each of the approaches described for various fibre dosages: 3, 4.5, 6 and 10 kg/m 3 , therefore reproducing the experimental results of [28,29].

Results
It is interesting to observe the differences among these approaches. In the case of the rigid-plastic model, diagrams always show a plateau and a subsequent load decrease that describes an exponential-like curve with small load jumps. As these drops are due to the abrupt load decrease described by the rectangular diagram, when an element reaches w = w u , there is a sudden load decrease.
In the case of the "linear model with unlimited development", the load-deflection diagrams exhibit a linear behaviour after the concrete matrix fails and, depending on the fibre dosage the specimen, may show a softening behaviour (for 3 and 4.5 kg/m 3 ) or a hardening behaviour (for 6 and 10 kg/m 3 ). Obviously, in the case of diagrams with hardening, this behaviour cannot be realistic, since the load may increase without a limit. Nevertheless, it is interesting to compare these results with those of the linear model with limited development that display an analogous behaviour until the most damaged element reaches w = w u , when the load-deflection diagrams show a highly similar behaviour to that observed with the rigid-plastic model. In all cases, the load drop follows an exponential-like shape with small jumps that can be explained as for the rigid-plastic model. followed by a load recovery and another load decay later. Second, the last load decay, which identifies the last part of the test until failure, shows a curved shape that has the opposite curvature if compared with results obtained with the rigid-plastic model and the linear model with limited development.

Discussion
In Figure 8 it can be seen that in the case of PFRC3 there are remarkable changes in the shape of the curves obtained through using the post-peak constitutive functions analysed. Regarding the rigid-plastic model and the linear model proposed by the MC2010, none is capable of reproducing the peak load registered in the experimental tests (it should be taken into account that usually 3 kg/m 3 of polyolefin fibres do not allow PFRC to be classified as structural concrete). However, the trilinear approach accurately captures such a peak load. This phenomenon might be explained because the trilinear approach is the only one that establishes a progressive unloading after reaching peak load. Similarly, the trilinear approach is able to simulate the unloading process of the PFRC3 formulation. In the case of the rigid-plastic approach, when the deflection increases the load-deflection curve features a plateau. This stretch ends when a deflection of 3 mm is reached. From that deflection onwards, a progressive unloading takes place which does not fit the experimental behaviour. In the case of the "linear approach (unlimited)", as the slope is negative Lastly, if the results obtained with a trilinear model are observed, remarkable differences can be highlighted compared with any of the results obtained with the other approaches. First, the maximum load reached that identifies the beginning of the fracture process is higher than in any of the other cases and the shape shows an initial load decay, followed by a load recovery and another load decay later. Second, the last load decay, which identifies the last part of the test until failure, shows a curved shape that has the opposite curvature if compared with results obtained with the rigid-plastic model and the linear model with limited development.

Discussion
In Figure 8 it can be seen that in the case of PFRC3 there are remarkable changes in the shape of the curves obtained through using the post-peak constitutive functions analysed. Regarding the rigid-plastic model and the linear model proposed by the MC2010, none is capable of reproducing the peak load registered in the experimental tests (it should be taken into account that usually 3 kg/m 3 of polyolefin fibres do not allow PFRC to be classified as structural concrete). However, the trilinear approach accurately captures such a peak load. This phenomenon might be explained because the trilinear approach is the only one that establishes a progressive unloading after reaching peak load. Similarly, the trilinear approach is able to simulate the unloading process of the PFRC3 formulation. In the case of the rigid-plastic approach, when the deflection increases the load-deflection curve features a plateau. This stretch ends when a deflection of 3 mm is reached. From that deflection onwards, a progressive unloading takes place which does not fit the experimental behaviour. In the case of the "linear approach (unlimited)", as the slope is negative the load decreases as the deflection increases. Consequently, such simulation is able to maintain its loadbearing capacity beyond 6 mm of deflection. The results obtained with the "linear model (limited)" softening curve also boast a descending slope until a progressive unloading starts at 3 mm of deflection. From that point onwards, the curve is not capable of reproducing the experimental behaviour. It should be underlined that the rigid-plastic approach and both linear approaches obtain load-deflection curves that are above the experimental one until a deflection of 2.5 mm is reached. For greater deflections, the experimental curve is above that predicted by both linear approaches. In the case of the rigid-plastic model, a similar comment could be made. However, in the case of the trilinear approach a close reproduction of the experimental behaviour PFRC3 is achieved throughout the test. the load decreases as the deflection increases. Consequently, such simulation is able to maintain its load-bearing capacity beyond 6 mm of deflection. The results obtained with the "linear model (limited)" softening curve also boast a descending slope until a progressive unloading starts at 3 mm of deflection. From that point onwards, the curve is not capable of reproducing the experimental behaviour. It should be underlined that the rigidplastic approach and both linear approaches obtain load-deflection curves that are above the experimental one until a deflection of 2.5 mm is reached. For greater deflections, the experimental curve is above that predicted by both linear approaches. In the case of the rigid-plastic model, a similar comment could be made. However, in the case of the trilinear approach a close reproduction of the experimental behaviour PFRC3 is achieved throughout the test. In the part of Figure 8 that corresponds to PFRC4.5, it may be observed that the differences perceived in the curves related with PFRC3 still appear. The peak load, as in the previous case, is not captured by any approach proposed by the MC2010 or even the modification proposed for the linear model. However, the trilinear constitutive model is able to simulate the flexural tensile behaviour of the material prior to the development of the In the part of Figure 8 that corresponds to PFRC4.5, it may be observed that the differences perceived in the curves related with PFRC3 still appear. The peak load, as in the previous case, is not captured by any approach proposed by the MC2010 or even the modification proposed for the linear model. However, the trilinear constitutive model is able to simulate the flexural tensile behaviour of the material prior to the development of the failure surface. As in the case of PFRC3, the trilinear softening function was able to simulate the experimental behaviour of the PFRC4.5 in all its stages. Regarding the behaviour of the rigid-plastic model, it may be seen that its main feature is a constant load-bearing capacity until a 3 mm deflection is reached. From that point onwards the rigid-plastic model is below the experimental curve. As regards to both linear approaches are above the experimental curve between 0.5 and 2.5 mm of deflection. There are few differences between the rigid-plastic approach and the linear approach because the slope of the linear model is close to zero. One feature that could be highlighted is that the unlimited linear model behaviour is close to the experimental one for deflections between 2.5 and 6 mm.
Regarding the curves shown for PFRC6, a similar trend may be perceived if the experimental curve is compared with the one obtained with the rigid-plastic approach. For deflections ranging from 0.3 mm to 2.5 mm, the simulated curve is above the experimental one as the plateau is close to the experimental maximum post-peak load. From 2.5 mm of deflection, the unloading branch of the rigid-plastic model features a steep degree of unloading that is below the experimental behaviour. Regarding the linear model with the unlimited option, it should be pointed out that as the slope of the load-deflection curve is positive, the model is not able to capture the material failure as its load bearing capacity increases with deflection. Thus, beyond 3.5 mm of deflection the linear model is above the experimental curve. This also occurs between 0.3 mm and 2.5 mm of deflection. If the limited linear model is studied, similar comments to those previously mentioned regarding its unlimited version could be made until 2.5 mm of deflection is reached. From greater deflection values there are remarkable differences between the behaviour of the model and the real one because steep unloading appears in the simulated curve. As regards the simulation performed with the trilinear approach, it could be said that it is able to reproduce the shape and values of the experimental curves although there are certain differences from 3.5 mm onwards. In this case, the experimental behaviour is slightly above the simulated one.
Regarding PFRC10, similar features to those observed in PFRC6 can be seen. The rigid-plastic approach has similar characteristics than those previously mentioned in the case of PFRC6. A similar trend may be mentioned when the differences between the experimental curve and both linear approaches are evaluated. On this occasion, the slope of the linear model is positive, and no failure of the material occurs. If the model is limited a under evaluation of the behaviour of PFRC10 is observed. Regarding the trilinear approach, for deflection values greater than the maximum post-peak load such a curve marginally underestimates the load-bearing capacity of the material. Figure 9 provides a comparison among the loads registered at different deflection values in the simulations performed and takes into account the three constitutive models studied. The values shown are relative to the value obtained in the experimental test which is considered as unity. Consequently, if a load value obtained in the simulation is greater than the one obtained in the test, the column that appears in Figure 9 is greater than one. On the contrary, values lower than the unity reflect that the simulation underestimates the real load-bearing capacity of the material. Such load analysis has been performed at three deflections: 0.5 mm, 2.5 mm and 5 mm. The values of deflection studied are based on the almost linear relation that exists between the crack width and the deflection. Consequently, it may be considered that a deflection value of 0.5 mm is close to the crack width used for checking the SLS of the structural element. Correspondingly, the deflection value of 2.5 mm has been chosen because such a value is close to the crack width chosen when analysing the ultimate limit state of the structural element. Lastly, regarding the 5 mm of deflection, such a value has been chosen as it may reflect the suitability of the constitutive models to detect the ultimate bearing capacity of the material at great deflections. When the analysis of the SLS is performed for the PFRC formulations, several relevant features are evident. In the case of PFRC3, it may be seen that while the trilinear approach underestimates the capacity of the material, the two proposals of the MC2010 boast loads greater than the experimental one. In the case of the linear model, both versions obtain a load approximately 1.6 times the experimental one. In the case of the rigid-plastic model, although it also overestimates the experimental value registered, it is much closer to the real one. In the case of PFRC4.5, all the values obtained in the simulations are greater than the experimental ones (except the trilinear one which is much closer to the real value). On this occasion, the load values of the rigid-plastic model are close to those obtained with the linear models. Similar results were obtained for the PFRC6 and PFRC10 formulations.
These features seem to be a consequence of the load-crack mouth opening displacement curve taken as reference to set the requirements needed for considering the contribution of fibres in FRC (see Figure 2a). As such a reference curve boasted an 3 value lower than fR1 if the design values are obtained by reducing the load obtained in 3 , such a constitutive model might have an accurate correlation with material behaviour. Nevertheless, in the case of PFRC, given that the value of 3 is greater than 1 , the design value obtained for reduced strains might not be an accurate reproduction of the experimental behaviour.
When analysing the ULS, similar to the behaviour that appears at 2.5 mm of deflection, it may be seen that the trilinear approach and those proposed by the MC2010 reflect a notable degree of accuracy with the experimental behaviour of the material.
Regarding the loads obtained with the trilinear approach, all load values are within a 20% margin from the experimental at 5 mm of deflection. In the case of the rigid-plastic model and the limited linear model, both are far from the experimental results for any fibre dosage. In the case of the unlimited linear model, the load values are below the experimental ones when the slope is negative. On the contrary, when the slope is positive such values are slightly above those registered in the tests.
Another material property that might be of interest is the fracture energy required to generate certain deflections of the concrete element. Consequently, the fracture energy at the aforementioned deflections has been obtained for all the concrete formulations and constitutive models analysed. Those results can be seen in Table 3, where the variation of the fracture energy of the curves obtained in the simulations performed by using the constitutive models implemented has been related to that obtained in the experimental tests. When the analysis of the SLS is performed for the PFRC formulations, several relevant features are evident. In the case of PFRC3, it may be seen that while the trilinear approach underestimates the capacity of the material, the two proposals of the MC2010 boast loads greater than the experimental one. In the case of the linear model, both versions obtain a load approximately 1.6 times the experimental one. In the case of the rigid-plastic model, although it also overestimates the experimental value registered, it is much closer to the real one. In the case of PFRC4.5, all the values obtained in the simulations are greater than the experimental ones (except the trilinear one which is much closer to the real value). On this occasion, the load values of the rigid-plastic model are close to those obtained with the linear models. Similar results were obtained for the PFRC6 and PFRC10 formulations.
These features seem to be a consequence of the load-crack mouth opening displacement curve taken as reference to set the requirements needed for considering the contribution of fibres in FRC (see Figure 2a). As such a reference curve boasted an f R3 value lower than fR1 if the design values are obtained by reducing the load obtained in f R3 , such a constitutive model might have an accurate correlation with material behaviour. Nevertheless, in the case of PFRC, given that the value of f R3 is greater than f R1 , the design value obtained for reduced strains might not be an accurate reproduction of the experimental behaviour.
When analysing the ULS, similar to the behaviour that appears at 2.5 mm of deflection, it may be seen that the trilinear approach and those proposed by the MC2010 reflect a notable degree of accuracy with the experimental behaviour of the material.
Regarding the loads obtained with the trilinear approach, all load values are within a 20% margin from the experimental at 5 mm of deflection. In the case of the rigid-plastic model and the limited linear model, both are far from the experimental results for any fibre dosage. In the case of the unlimited linear model, the load values are below the experimental ones when the slope is negative. On the contrary, when the slope is positive such values are slightly above those registered in the tests.
Another material property that might be of interest is the fracture energy required to generate certain deflections of the concrete element. Consequently, the fracture energy at the aforementioned deflections has been obtained for all the concrete formulations and constitutive models analysed. Those results can be seen in Table 3, where the variation of the fracture energy of the curves obtained in the simulations performed by using the constitutive models implemented has been related to that obtained in the experimental tests.    In the case of the trilinear model, it can be seen that a notable precision was obtained at all the deflections analysed. Regarding the rigid-plastic model, there was a certain degree of scattering between the values obtained at 0.5 mm between the experimental values and the simulated ones. Such scattering was reduced at 2.5 mm of deflection and increased at 5 mm of deflection. In the case of the linear model, the accuracy of the prediction of the fracture energy absorbed by the material increased as the deflection grew. In Table 3 it may be seen that the proposed limitation of the linear approach did not enhance the accuracy of the prediction of the fracture energy absorbed.

Conclusions
The present contribution has implemented the constitutive relations proposed by the MC2010 for FRC in a user material subroutine seeking to assess the adaptability of these constitutive relations to PFRC. In addition, a trilinear model and an unlimited version of the linear model proposed by MC2010 have been employed as comparisons.
The simulations carried out with the implementation of the trilinear softening function proposed by the authors showed a remarkable resemblance with the experimental behaviour of PFRC. Regarding the simulations performed through using the constitutive models proposed by MC2010, it was observed that they might not replicate the flexural behaviour of PFRC precisely. Neither the peak load nor the unloading process appeared in the simulations performed with the rigid-plastic and linear constitutive relations.
When analysing the load values obtained at 0.5 mm of deflection, it is seen that the constitutive relations proposed by the MC2010 are above the experimental results. At ULS, the relations proposed by the MC2010 and that proposed by the authors are able to reproduce with notable exactness the experimental load values. Regarding the ultimate load-bearing capacity of the material, which was assumed at a deflection of 5 mm, both the trilinear and unlimited linear model were close to the load-bearing capacity of the material. On the contrary, the rigid-plastic model and the limited version of the linear model clearly underestimated the mechanical properties of PFRC.
Based on the aforementioned arguments, it seems that the rigid-plastic and the linear models proposed by the MC2010 might be employed with caution at SLS when applied to PFRC. This might be a result of the limited degree of accuracy that the idealisation considered in MC2010 may offer when applied to FRC formulations where f R3 is notably greater than f R1 . The authors consider that if the results of the EN-14651 fracture tests boast a clear three-stretch behaviour after reaching the peak load (unloading-reloadingfinal unloading), a trilinear constitutive relation might be adopted as an alternative to the models proposed by MC2010. This is relevant given that the latter models were developed considering that steel fibres were the main option for reinforcing concrete with fibres.
The limitation of this study might be that the differences detected were determined for certain constitutive models stated in MC2010 and consequently if such models are modified the present conclusions should be re-evaluated.

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