Two-Dimensional Numerical Model of the Fracture Process in Steel Fibre Reinforced Concrete with the Continuum Strong Discontinuity Approach and Functional Data Analysis

Abstract

This paper presents the formulation of a two-dimensional numeri-cal model able to describe the fracture process in structural mem-bers of steel fibre reinforced concrete (SFRC) from the volume ratio of the fibres and the mechanical properties of the compo-nents: a concrete matrix and a set of steel fibres with a random orientation. The relationship between the stress and the strain fields of the composite material is obtained using the mixture theory with a compatibility strain of its component materials. The concrete matrix is represented with a scalar damage constitutive model with a softening strain and a different strength in tension and compression. The mechanical strain of an insulated fibre and the slip between the fibre and the matrix are simultaneously de-scribed with a one-dimensional plasticity constitutive model. The cracking of the composite material indicates a jump in the dis-placement field and non-bounded values of the strain field, which are represented by the Continuum Strong Discontinuity Ap-proach. The model has been implemented in the framework of the nonlinear analysis with the Finite Element Method, using con-stant strain triangular elements. Moreover, the fibres distribution and orientation change randomly in each finite element and each simulation or observation. The structural responses of the simula-tions are treated as curves and analysed by tools from the Func-tional Data Analysis. Confidence intervals for the structural re-sponse are built using bootstrap methodology. Finally, experi-mental tests of SFRC members subjected to tension and bending are simulated. The structural response and the cracking patterns obtained from the numerical simulation are satisfactory.

Keywords:
Computational mechanics; fracture mechanics; strong discontinui-ty approach; mixture theory; steel fibres reinforced concrete (SFRC); functional data analysis (FDA)

1 INTRODUCTION

The mechanical behaviour of SFRC can be described by four types of approaches implemented with finite element methods, according to the analysis scale. The first represents SFRC by a unique constitutive model; the second establishes the behaviour of the composite material through the mixture theory and uses different constitutive models for the concrete matrix and the steel fibres. The third approach explicitly represents each fibre with lattice elements and the mesh of solid finite elements represents the plain concrete. Finally, in the fourth approach the mortar, the aggregates and the steel fibres are described as discrete entities with different constitutive models.

In the first approach type, Peng and Meyer (2000) characterise the composite material with a constitutive model of scalar continuum damage, where the damage parameter depends on the volumetric fraction of the steel fibres. Likewise, Oliveira de Barros (1995) implements a non-linear constitutive model which represents the fracture process by a rotating smeared crack model. This approach has a low computational cost but requires knowing the mechanical properties of the SFRC for a defined amount of steel as an input data.

In the second approach type, Nofal (1997) uses a damage-plasticity constitutive model which represents the degradation of the stiffness and the permanent strain after the unloading of the material. A set of randomly distributed elliptical micro-voids represents the cracking in the plane concrete matrix. The one-dimensional plasticity model describes the behaviour of the short fibres and their slipping. The stiffness matrix of the SFRC is obtained from the integral of a probability function of the random orientation of the fibres at the plane. A similar work, Luccioni et al. (2012) represents the plain concrete with a damage-plasticity constitutive model, the steel fibres and the slipping between them and the matrix through an orthotropic homogeneous plasticity model. Likewise, Li and Li (2001) apply the mixture theory for composite materials, representing the plain concrete, the steel fibres and their slipping with isotropic damage constitutive models. In another work Mora et al. (2013), the multi-field theory is applied to describe the axial behaviour of the fibres and the adherence between fibre and matrix. Each material point is formed by a concrete matrix, which is represented by a constitutive model of damage and by some packages of steel fibres with perpendicular orientations represented with one-dimensional plasticity constitutive models. This approach type needs the amount of fibre reinforcing and the mechanical properties of the concrete and the steel as input data.

In some works of the third approach type (Prós, 2011; Cunha et al. 2012; Kang et al. 2014; Yu et al. 2016), the ends of lattice element of the steel fibres do not match with the nodes of the finite element mesh of plain concrete. Contrarily, the displacement compatibility is defined by employing Immersed Boundary Methods. In another approach, in Radtke et al. (2010), the domain is represented by triangular elements of plain concrete, using an isotropic damage constitutive model, and the contribution of each fibre is defined using an axial force at its ends.

A sample of the fourth approach type for plane stress condition is developed by Bui (2007). He represents the cement matrix by six-node triangular elements, each aggregate grain by a set of six-node triangular elements and each steel fibre by a lattice element. These last elements are connected to the nodes of the triangular elements of the cement matrix. The Mohr-Coulomb criteria and a softening law establish the beginning and the evolution of the fracture process, directly applied at the nodes, preserving the elastic behaviour of triangular elements. Some authors present the extension of this model to the three-dimensional space (Leite et al. 2007; Gal and Kryvoruk 2011; Caggiano et al. 2012). This approach type allows describing the geometry of each material component, particularly, the shape of the aggregates and the length of fibres. However, the orientation, the distribution and the shape of the fibres and the size and the distribution of the aggregates are chosen randomly requiring a specific meshing technique and high computational cost.

The fracture process of the concrete in the previous approaches was represented using smeared crack models (Oliver et al. 1990; De Borst et al. 2004; Cunha et al. 2012), discrete models (De Borst et al. 2004; Caggiano et al. 2012; Manzoli et al. 2012), Continuum Strong Discontinuity Approach (CSDA) (Oliver 1996a), and strain injection approach ((Dias et al. 2011), among others.

On the other hand, a set of scalar data or curves can be treated employing the standard statistical or the functional data analysis (FDA) (Ramsay and Silverman 2005; Ferraty and Vieu 2006), respectively. Several authors have sustained efforts to adapt standard statistical methods for scalar values to this data type (Cardot et al. 1999, 2007; Locantore et al. 1999; Cuevas et al. 2002; Ferraty and Vieu, 2002; Dabo-Niang and Rhomari, 2009).

The structural response of an SFRC member is considered as a curve obtained from a continuous stochastic process, which represents the relationship between the displacement and the applied external loads. This structural response depends on random parameters of the composite material; therefore, a set of structural responses for different random parameters is treated employing tools belonging to functional data analysis. Some Authors has assumed that similar relationships in physical problems are functional data (Fletcher et al. 2004; Antoniadis and Sapatinas, 2007; Del Barrio et al. 2007; López- Pintado and Romo, 2007; Nerini and Ghattas, 2007; Tarrio, 2012).

In recent researches about the statistical control of processes, characterised by nonlinear relationships, Guevara and Vargas (2013) used the bootstrap method to provide a better estimation of the confidence intervals without requiring a specific distribution of probability.

Previous numerical approaches using statistical procedures to describe the functional relationship of the SFRC structural response were not found. Generally, the random parameters are represented by factors assigned to one simulation of the problem.

In essence, the objective of this work is to estimate the median of the mechanical response and its confidence interval, for structural members of SFRC by means of a stochastic procedure that is based on the functional data analysis, which has proven to be adequate for similar problems. The sample of continuous functions representing the mechanical response is obtained from a two-dimensional nonlinear finite element analysis using the formulation described in this paper. The relevance of this formulation is that the numerical model in the deterministic procedure has a low computational cost due to the composite material is homogenized by subdomains along the structural member. This allows obtaining a large number of simulations by applying a stochas-tic procedure in which the properties that depend on the distribution and direction of the fibres are randomly assigned to the subdomains. In the end, the result consists of a range of probable responses that have into account the variability that is proper of SFRC material.

The model is able to describe the fracture process in SFRC structural members, from the vol-ume ratio of the fibres and the mechanical properties of concrete and steel. A concrete matrix and a set of steel fibres with random orientation constitute this composite material. The relationship between the stress and strain fields of the composite material is obtained through the mixture theory with strain compatibility of its component materials (Truesdell and Toupin, 1960) as it is indicated in the second type of approaches. The concrete matrix is represented with a scalar damage constitutive model with a softening strain and a different strength in tension and compression (Oliver et al. 1990; Linero et al. 2007). The mechanical strain of discrete fibres and the slip between fibres and matrix are simultaneously described with a one-dimensional plasticity constitutive model (Simó and Hughes, 1998). A crack in the composite material indicates a jump displacement field and non-bounded values of the strain field, which are represented by the continuum stro,ng discontinuity approach (Oliver, 1996a). The model has been implemented into the finite element method using constant strain triangular elements. Moreover, the distribution and orientation of fibres change randomly in each finite element and each simulation or observation. he structural response of these observations is treated using tools from the Functional Data analysis (FDA) (Ramsay and Silverman 2005; Ferraty and Vieu, 2006). Mainly, the concept of the depth allows to sort the curves from the centre to the outer and to estimate the median curve. The confidence intervals for the structural response are found using bootstrap procedures for continuous functions defined in a compact interval $ℝ$, as is indicated in the reference Guevara and Vargas (2013).

2 MAIN PROCEDURE OF THE PROPOSED MODEL

SFRC is a composite material conformed by fibres of random direction and distribution, embedded in a plain concrete matrix. The proposed model represents the fracture process in structural members of SFRC, particularly considering plane stress condition, infinitesimal strain, and static loads.

The steel fibres can be oriented in any direction within the concrete, exhibiting a three-dimensional behaviour of SFRC. Nevertheless, when the stress components out-of-plane are null, the problem is simplified into a two-dimensional domain, through the statistical participation of steel fibres projected on the plane.

Each material point shows the macroscopic behaviour of the SFRC considering: the steel fibres uniformly distributed in the concrete matrix, the compatibility of the longitudinal strain between the materials, and the composite material stress equal to the weighted addition of the stress of each component material. Consequently, the numerical model does not depend on the relation between the fibres length and the finite element size.

This numerical model establishes an external procedure or stochastic procedure and an internal procedure or deterministic procedure as it is indicated in Table 1.

Table 1
Main procedure of proposed model.

input	Mechanical properties of the concrete and the steel, global volumetric ratio of the steel fibres, geometry, loads and boundary conditions of the solid.
1	Meshing of the geometry problem
2	Do		j-observation, from 1 to nj-observations
3			Random assignment of the orientation and the variability in the amount of steel fibres in each finite element (part of the stochastic procedure).
4			Deterministic Procedure: Non-linear finite element analysis with CSDA.
5			Partial output: Stress and strain state, structural response and cracking pattern.
7	End
8	Partial output: Set of nj-observations: structural responses or curves load vs displacement.
9	Stochastic procedure: Functional Data Analysis (FDA) of the set of nj-observations
Output	Median and confidence intervals of the set of structural responses.

In the deterministic procedure, the structural response and the cracking pattern of the structural member are obtained from a non-linear finite element analysis with CSDA and a set of parameters assigned to a simulation (Oliver and Huespe, 2004a; Oliver et al. 2008). An observation corresponds to a simulation where a random orientation and a distribution of fibres in each finite element are assigned. Likewise, in the stochastic procedure, the random parameters of each observation are computed, and an FDA is executed with the collection of deterministic results of the observations (Ramsay and Silverman, 2005; López-Pintado and Romo, 2009; Sun and Genton, 2011).

3 DETERMINISTIC PROCEDURE OF THE PROPOSED MODEL

The deterministic procedure of the proposed numerical model represents the fracture process of a SFRC structure, in which the orientation and distribution of the fibres are defined. The deterministic procedure is based on the Continuum Strong Discontinuity Approach (CSDA), the mixture theory for composite materials and the constitutive models of damage and of plasticity applied to each component material. Likewise, this procedure is implemented in a finite element method, considering infinitesimal strains, static loading, plane stress condition and the nonlinear behaviour of the material.

3.1. Constitutive model of the plain concrete

The two-dimensional behaviour of the concrete matrix is described with a constitutive model of isotropic scalar damage, in which the tensile strength ${\sigma }_t^c$ is different to the compressive strength ${\sigma }_c^c$ (Oliver et al. 1990) . The tangent constitutive equation of plain concrete ${\dot{\sigma }}^c=C_{tg}^c:{\dot{\epsilon }}^c$, relates the rate stress tensor ${\dot{\sigma }}^c$ and the rate strain tensor ${\dot{\epsilon }}^c$. The tangent constitutive tensor $C_{tg}^c$ and the factor $\varphi$ are:

$C_{tg}^c=\frac{q^c}{r^c}{C}_{}^c-\left(\frac{q^c-H^c{r}^c}{{\left(r^c\right)}^3}\right)\cdot \left[\frac{{\left(r^c\right)}^2}{\varphi }\left({\overline{\sigma }}^c\otimes C^c:{\partial }_{\overline{\sigma }}\varphi \right)+{\varphi }^2\left({\overline{\sigma }}^c\otimes {\overline{\sigma }}^c\right)\right]$ (1)

$\varphi =\frac{{\displaystyle \sum _{i=1}^3\langle {\overline{\sigma }}_i^c\rangle }}{{\displaystyle \sum _{i=1}^3\left|{\overline{\sigma }}_i^c\right|}}\left(1-\frac{{\sigma }_t^c}{{\sigma }_c^c}\right)+\frac{{\sigma }_t^c}{{\sigma }_c^c}$ (2)

In the previous equations, the effective stress tensor is equal to ${\overline{\sigma }}^c=C^c:{\epsilon }^c$, where ${\epsilon }^c$ is the strain tensor and $C^c$ is the elastic constitutive tensor, which is defined as $C^c=\left(E^c/1-{\left({\nu }^c\right)}_{}^2\right)\left[{\nu }^c(1\otimes 1)+(1-{\nu }^c)I\right]$ for plane stress conditions. The norm ${\tau }_{\epsilon }=\sqrt{{\overline{\sigma }}^c:{\left(C^c\right)}^{-1}:{\overline{\sigma }}^c}$ is a scalar function of the strain state and determines the elastic domain of the model. Mac-Auley parenthesis function $\langle {\overline{\sigma }}_i^c\rangle$ and the absolute value function $\left|{\overline{\sigma }}_i^c\right|$ are applied on the principal effective stresses ${\overline{\sigma }}_i^c$. The softening rule relates the stress-like internal variable $q^c$ and strain-like internal variable $r^c$, by the rate equation ${\dot{q}}^c=H^c(r^c){\dot{r}}^c$. The stress-like internal variable is obtained as $r^c\left(t\right)=\max \left(r^c\left(t^{\prime }\right),{\tau }_{\epsilon }\left(t^{\prime }\right)\right)$, where $t^{\prime }\in \left[\mathrm{0,}t\right]$ and this initial value is $r^c\left(t=0\right)={\sigma }_t^c/\sqrt{E^c}$. In this relation, the softening parameter is defined as $H^c(r^c)=d{q}^c/d{r}^c\le 0$ which is determined from the relation between the normal stress and the crack opening (Bazant and Planas 1998). The evolution law, the damage criteria, and the loading and unloading conditions of this constitutive model are explained by other authors (Oliver et al. 1990).

The constitutive model requires the Young’s modulus $E^c$, the Poisson’s ratio ${\nu }^c$, the fracture energy per surface unit $G_f^c$ and the tensile and compressive strengths. These mechanical properties can be obtained from experimental tests ASTM-C469, ASTM-C39, ASTM-C496 and RILEM TC-50 (RILEM 1985; ASTM 2004, 2014a, b) . Figure 1(a) and Figure 1(b) indicate the elastic domain and the uniaxial stress - strain curve in terms of principal stresses for a plane stress condition.

Figure 1
Constitutive models for the component materials: (a) elastic domain in function of the principal stresses for the plain concrete, (b) uniaxial relation between the normal stress and the longitudinal strain for the plain concrete, (c) relation between the normal stress and the longitudinal strain for the deformable-sliding fibre model, (d) relation between the shear stress and angular strain for dowel action model.

3.2. Constitutive model of the deformable-sliding fibre

The behaviour of the fibres is associated to three mechanisms: (i) the mechanical deformation due to axial stress, (ii) the slipping of the fibre relative to the concrete matrix and (iii) the shear transfer between the faces of a crack by the fibres. The pull-out tests of a steel fibre into concrete matrix describe the first and second mechanism simultaneously. A relation between the normal stress and the longitudinal strain is obtained from this test type (Naaman et al. 1991; Laranjeira et al. 2010; Zile and Zile, 2013).

This model considers two sets of parallel fibres with different orientation. The axial behaviour of each set of fibres is represented by one-dimensional plasticity isotropic constitutive model with a softening strain (Simó and Hughes, 1998; de Souza et al. 2008), named in this work as deformable-sliding fibre model. This model is defined by the elasticity modulus $E^f$, the maximum normal stress ${\sigma }_{\max }^f$ and the softening parameter $H^f$, as shown in Figure 1(c).

The tangent constitutive equation of this model ${\dot{\sigma }}^f=E_{tg}^f{\dot{\epsilon }}^f$, relates the normal stress rate ${\dot{\sigma }}^f$ and the longitudinal strain rate ${\dot{\epsilon }}^f$. The tangent constitutive modulus is $E_{tg}^f=E^f$ in elastic loading or unloading regimes and it is equal to $E_{tg}^f=E^f{H}^f/(E^f+H^f)$ in plasticity regime. The softening parameter $H^f$ is obtained from pull-out tests and its value is assumed as a constant.

The softening rule is defined as ${\dot{q}}^f=H^f{\dot{r}}^f$, where $q^f$ is the stress-like internal variable and $r^f$ is the strain-like internal variable, that is described by Simó and Hughes (Simó and Hughes, 1998) in the clasical one-dimensional plasticity constitutive model. Other elements of this constitutive model as the flow rule, the evolution law, the plastic yielding criterion, and the loading and unloading conditions are explained in several references (Simó and Hughes, 1998; de Souza et al. 2008).

3.3. Constitutive model of the dowel action for SFRC

When a crack is formed in the SFRC and its sliding is greater than its opening, the shear stress is mainly transferred to the steel fibres crossing the sides of the crack. This phenomenon is produced by two mechanisms: (i) bending of each fibre between the sides of the crack due to relative displacement perpendicular to its longitudinal axis, and (ii) direct shear of each fibre crossing the crack. Particularly in this work, the second mechanism is considered the most representative, because the crack opening is not substantially greater than the fibre diameter. A similar behaviour has been studied in concrete reinforced with steel bars under the name dowel action (Park and Paulay, 1975; Nawy, 2009).

The dowel action for SFRC can be represented with a one-dimensional plasticity isotropic model with a softening strain, which relates the shear stress rate ${\dot{\tau }}^d$ with the angular strain rate ${\dot{\gamma }}^d$ in the cross-section of the fibres, thus: ${\dot{\tau }}^d=G_{tg}^d{\dot{\gamma }}^d$. Figure 1(d) indicates the curve between the shear stress and the angular strain, where $G^d$ is the elasticity shear modulus, ${\tau }_{\max }^d$is the maximum shear stress and $H^d$ is the softening parameter.

The softening rule is ${\dot{q}}^d=H^d{\dot{r}}^d$, where $q^d$ is the stress-like internal variable and ${\dot{r}}^d$ is the strain-like internal variable. The flow rule, the evolution law, the plastic yielding criterion, and the loading and unloading conditions of this constitutive model are presented extensively in several references (Simó and Hughes, 1998; de Souza et al. 2008).

3.4. Mixture theory applied to SFRC

Three considerations of the deformation compatibility are assumed in the numerical model (Linero et al. 2007). First, the strain rate tensor of SFRC is equal to the strain rate tensor of the plain concrete. Second, the strain component rate of the plain concrete constitutive model in the fibre orientation is equal to the longitudinal strain rate of the deformable-sliding fibre constitutive model. Third, it is assumed that the angular strain rate of the dowel action model is equal to the angular strain rate in the plain concrete model after the crack formation.

A set of parallel fibres in the r direction and another set in the s direction are assumed at each material point, where r is perpendicular to s (Figure 2). Therefore, the compatibility conditions of the strain rate are:

${\dot{\epsilon }}^c=\dot{\epsilon }\left(plainconcrete\right)$ (3)

${\dot{\epsilon }}^{fr}=r\cdot \dot{\epsilon }\cdot r\left(deformable-slidingfibreinrdirection\right)$ (4)

${\dot{\epsilon }}^{fs}=s\cdot \dot{\epsilon }\cdot s\left(deformable-slidingfibreinsdirection\right)$ (5)

${\dot{\gamma }}^d=2r\cdot \dot{\epsilon }\cdot s\left(dowelactioninrsplane\right)$ (6)

where $\dot{\epsilon }$ is the tensor of strain rate of the SFRC, ${\dot{\epsilon }}^c$ is the tensor of strain rate of the plain concrete model, ${\dot{\gamma }}^d$is the angular strain rate of the dowel action model, ${\dot{\epsilon }}^{fr}$ and ${\dot{\epsilon }}^{fs}$ are the longitudinal strain rate of the deformable-sliding fibre constitutive model in r and s direction, respectively.

Figure 2
Representation of the cracking of SFRC: (a) structural member and (b) material point

The tensor of stress rate of the composite material is considered equal to the weighted addition of the stress rates of the constitutive models of the component materials (Truesdell and Toupin, 1960; Oller, 2003), thus:

$\dot{\sigma }=k^c{\dot{\sigma }}^c+k^{fr}(r\otimes r){\dot{\sigma }}_{}^{fr}+k^{fs}(s\otimes s){\dot{\sigma }}_{}^{fs}+2(k^{fr}+k^{fs}){(r\otimes s)}^{sym}{\dot{\tau }}^d$ (7)

The weight factors $k^c$,$k^{fr}$ and $k^{fs}$correspond to the volumetric participation of the plain concrete and the fibres oriented in r and s directions, respectively.

3.5. Description of fracture process by means of the CSDA

The Continuum Strong Discontinuity Approach (CSDA) can predict the formation and propagation of the cracks in a simple material (Oliver 1996a, b, Oliver and Huespe 2004a, b) or in a composite material (Linero et al. 2007; Manzoli et al. 2008; Oliver et al. 2008) . This approach establishes the existence of a jump of the displacement field $[\text{}[u]\text{}]$ on the discontinuity or fracture surface.

Let $\Omega$ be a solid exhibiting a strong discontinuity on the surface S of normal n, which splits the body into the domains ${\Omega }^{+}$ and ${\Omega }^{-}$ as is shown in Figure 2(b). The displacement rate in a material point x and at a time t is defined as (Oliver, 1996a):

$\dot{u}(x,t)=\dot{\overline{u}}(x,t)+M_S(x)\text{}[\text{}[\dot{u}]\text{}](x,t)$ (8)

Where $M_S(x)=H_S(x)-\phi (x)$, $H_S(x)=\{\begin{array}{c}0\forall x\in {\Omega }^{-}\\ 1\forall x\in {\Omega }^{+}\end{array}$ and $\phi (x)=\{\begin{array}{c}0\forall x\in {\Omega }^{-}/{\Omega }^h\\ 1\forall x\in {\Omega }^{+}/{\Omega }^h\end{array}$

The notation A/B indicates the set of all the material points of A that do not belong to B. Furthermore, $\dot{\overline{u}}$ and $\text{}[\text{}[\dot{u}]\text{}]$ correspond to the continuous part and to the jump of the displacement rate in the discontinuity, respectively.

In the continuum, the unit jump function $M_S(x)$ is computed as the subtraction between the Heaviside’s function $H_S(x)$ and an arbitrary continuous function $\varphi (x)$ limited by a small band ${\Omega }^h$ in $S$, as is indicated in Figure 3.

Figure 3
Special functions of displacement fields: (a) Heaviside’s function, (b) function $\phi$, and (c) unit jump function.

The strain rate field is evaluated with the symmetrical differential operator on the displacement rate field; so that, the strain is split in a compatible part in terms of continuous displacement and in an enhanced part in terms of displacement jump. Given that the gradient of Heaviside’s function produces an unbounded value equal to $\nabla H_S={\delta }_S^{}n$, the strain rate can be indicated as the sum of a regular or bounded part $\dot{\overline{\epsilon }}={(\nabla \otimes \dot{\overline{u}})}^{sym}-{(\nabla \phi \otimes [\text{}[\dot{u}]\text{}]\text{)}}^{sym}$ and a singular or unbounded part equal to ${({\delta }_S^{}n\otimes [\text{}[\dot{u}]\text{}]\text{)}}^{sym}$.

Dirac’s delta function ${\delta }_S^{}$ is approximately equal to a regularized function ${\delta }_S^k$ applied in a band $S^k$ inside the discontinuity $S$, whose regularization parameter $k$ is the minimum permitted by the computer precision. This regularized function is:

${\delta }_S^k(x)=\underset{k\to 0}{\lim }{\mu }_S(x)\frac{1}{k};{\mu }_S(x)=\{\begin{array}{c}1x\in S^k\\ 0x\notin S^k\end{array}$ (9)

where ${\mu }_S$ is a collocation function on $S^k$, therefore the strain rate $\dot{\epsilon }$ of the composite material is:

$\dot{\epsilon }={(\nabla \otimes \dot{u})}^{sym}=\dot{\overline{\epsilon }}+{({\delta }_S^{k}n\otimes [\text{}[\dot{u}]\text{}]\text{)}}^{sym}$ (10)

The previous equation shows an unbounded strain on the discontinuity surface. The constitutive models of the materials are modified to preserve bounded values for the stresses. The softening parameter is regularized in each material $•$, thus: $H_{•}=k{\overline{H}}_{•}$, where ${\overline{H}}_{•}$ is the intrinsic softening parameter.

Despite the discontinuity of the strain field, the discontinuity surface $S$ with normal vector $n$ exhibits the same traction vector between its sides, that is:

$\left[\left[\dot{\sigma }\right]\right]\cdot n={\dot{\sigma }}_{{\Omega }^{-}/S}\cdot n-{\dot{\sigma }}_{{\Omega }^{+}/S}\cdot n=0$ (11)

Likewise, the equilibrium equations in the volume and the boundary of the composite material are also fulfilled.

The activation time and the orientation of the discontinuity are determined by means of the material bifurcation analysis (Rice, 1976).

The deterministic procedure of the proposed model is based on a previous approach that was applied to structural members of concrete reinforced with steel continuum bars (Oliver et al. 2008). This approach has been modified to represent the short fibres with any orientation within SFRC.

4 STOCHASTIC PROCEDURE OF THE PROPOSED MODEL

The variability in the orientation and the amount of the steel fibres into the concrete matrix are random parameters related to the dispersion of the mechanical response of SFRC members. Consequently, the proposed model includes a stochastic procedure in which: the orientation and the variability of the amount of the fibres are randomly assigned to each finite element of the mesh for each observation, as is shown in Figure 4, moreover the results of the deterministic procedure of all observations are treated by means of Functional Data Analysis (FDA), as it is indicated in Line 3 and 9 of Table 1.

Figure 4
Orientation of the fibres into the finite element: (a) finite element mesh, (b) random direction r and s of three triangular elements, and (c) discrete values of the angle between x-axis and r-direction of a fibre.

4.1. Random assignment of parameters of the fibres

Previously, it was mentioned that the effect of the fibres on the SFRC mechanical behaviour depends on the amount and distribution of them into the concrete matrix. In practice, the average dosage of fibres is given by the weight of fibres per unit volume of SFRC $w^f$. Thus, the volume of fibres per unit volume of SFRC can be calculated as ${\rho }^f=w^f/{\gamma }_s$, where ${\gamma }_s$ is the steel specific weight of fibres.

For a two-dimensional problem, the effectiveness of fibres into the analysis plan can be accounted by the equivalent volumetric ratio $k_h^f$, that is the volume of fibres per unit volume of SFRC projected on the plane (Lamus, 2015) for an structural member with a defined thickness $t$and can be calculated through the form:

$k_h^f={\rho }^f{\vartheta }_t{c}_h^{dis}t$ (12)

Where ${\vartheta }_t$is the participation ratio on the plane and $c_h^{dis}$ is the dispersion coefficient.

The participation ratio on the plane ${\vartheta }_t$ is a function of the thickness of the structural member $t$ and the orientation range of the fibres outside of the plane $\alpha$. This ratio establishes the possible amount of fibres projected on the plane and is expressed as:

${\vartheta }_t={\displaystyle {\int }_0^t(2/z\alpha ){\displaystyle {\int }_0^{\alpha /2}\cos (\phi )d\phi }dz}$ (13)

where, $z$ is the distance between any fibre and the lateral boundary of the structural member and $\phi$ is the angle between any fibre and the plane under analysis.

The possible directions of any fibre whose centre distance to the structural member boundary is less than the half of its length are restrained as it is shown in the zone A of Figure 5(b). Because of that reason, ${\vartheta }_t$ is variable in this zone, taking a value of 1.0 when the fibre is near to the boundary and decreasing with the distance between the centre of the fibre and the boundary. Otherwise, a fibre at any other position can be oriented in any direction as it is shown in the zone B of Figure 5(a) and its participation ratio on the plane is equal to 0.636. In the examples presented in this work, the thickness of the structural member $t$ and the fibre length $l^f$ are constant and then the participation ratio on the plane is assumed to be constant.

Figure 5
Random assignment of parameters of the fibres: (a) and (b) participation ratio on the plane ${\vartheta }_t$ as a function of distance between the fibre center and the member boundary, and (c) normal probability distribution of the dispersion coefficient of the amount of fibres.

The dispersion coefficient $c_h^{dis}$ identifies the variability of the amount of fibres per unit volume along the structural member. This parameter exhibits a normal probability distribution which has a variation coefficient of 0.27 and a median associated to the global volumetric ratio of fibres (Soe, 2010). In this work, each finite element has a random dispersion coefficient choice of four values:$c_1^{dis}=\mathrm{0.688,}{c}_2^{dis}=\mathrm{0.912,}{c}_3^{dis}=\mathrm{1.084,}{c}_4^{dis}=1.308$, corresponding to the accumulated probabilities of 0.125, 0.375, 0.625 and 0.875, respectively, as shown in Figure 5(c). The four probable values of equivalent volumetric ratio of fibres $k_1^f,k_2^f,k_3^f$ and $k_4^f$ are obtained replacing the dispersion coefficients into Equation (12).

The equivalent volumetric ratio $k_{e,j}^f$ is the volume of fibres per volume unit of SFRC projected on the plane of a finite element $e$ for an observation $j$. The value of $k_{e,j}^f$ is randomly assigned among the values $k_1^f,k_2^f,k_3^f,k_4^f$, that is,

$k_{e,j}^f=\text{random}(k_1^f,k_2^f,k_3^f,k_4^f)$ (14)

The equivalent volumetric ratio of fibres $k_{e,j}^f$ includes the sum of the equivalent volumetric ratio of the two packages of fibres which are perpendicular to each other, that is, $k_{e,j}^f=k_{e,j}^{fr}+k_{e,j}^{fs}$, where $k_{e,j}^{fr}$ and $k_{e,j}^{fs}$ are the equivalent volumetric ratio in r direction and in s direction, respectively. The isotropy coefficient is $c_{\phi }=k_{e,j}^{fs}/k_{e,j}^{fr}$and indicates the relation of the amount of fibres between two perpendicular directions. Consequently, the equivalent volumetric ratio of the fibres in each direction can be expressed thus:

$k_{e,j}^{fr}=\frac{1}{c_{\phi }+1}{k}_{e,j}^f,k_{e,j}^{fs}=\frac{c_{\phi }}{c_{\phi }+1}{k}_{e,j}^f$ (15)

A methodology to obtain the isotropy coefficient from radiographies of SFRC specimens has been developed (Lamus, 2015) . This methodology depends on the size of the reference volume, the thickness of the structural member and the weight of the fibres per volume unit. In this paper, the isotropy coefficient is assumed as a constant equal to one, then the equivalent volumetric ratios $k_{e,j}^{fr}$ and $k_{e,j}^{fs}$ are equal to $0.5k_{e,j}^f$.

Particularly in this work, the angle between the x-axis and the direction r of a fibres package ${\theta }_l^{fr}$ can take one of ten values between ${\theta }_1^{fr}=90°$and ${\theta }_{10}^{fr}=-72°$, as shown in Figure 4(c). The angle between the x-axis and the r direction of a fibres package ${\theta }_{e,j}^{fr}$, at finite element e and in an observation j, is randomly chosen among the values of ${\theta }_l^{fr}$, that is,

${\theta }_{e,j}^{fr}=\text{random}({\theta }_1^{fr},\dots ,{\theta }_{10}^{fr})$ (16)

The s direction of the second fibres package is calculated as ${\theta }_{e,j}^{fs}={\theta }_{e,j}^{fr}+90°$.

According to the above, 40 material types are defined as the combination of 4 equivalent volume ratio and 10 r-direction of the fibres. A material type is randomly assigned to each finite element. This attempts to represent the real distribution of the fibres into the structural member and the discontinuity of the short fibres, keeping the global volumetric ratio of fibres constant.

In this work, observation is the deterministic procedure of the non-linear analysis of a finite element mesh with 40 types of materials. All the observations of a mechanical problem have the same mesh, but a different assignment of the material type to the finite elements, as shown in Figure 6. Here, the numbers from 1 to 40 identify the material type.

Figure 6
Randomly assignment of the materials type in the finite element mesh: (a) first observation and (b) second observation.

4.2. Overview of the functional data analysis of the set of the observations

The random distribution and orientation of the steel fibres into the SFRC are important factors in the dispersion of the mechanical response of the structural members. An experimental test is analysed by means of the proposed numerical model. Many simulations or observations are con-ducted to obtain representative results of the structural response. In the proposed numerical model, the structural response of each observation is different due to the change in the distribu-tion and orientation of the fibres in each finite element of the mesh. The estimation of the medi-an is made considering the nature of the structural response of the mechanical problem.

The structural response corresponds to a function between the applied load and a displacement component that represents the deformed shape. This relationship defines a continuum curve, which is studied by means of the FDA.

The functional data analysis done in this work has the following steps: (i) to generate a sample of continuous functions from the observations, (ii) to sort the sample and to measure the cen-trality of the curves, (iii) to detect and remove the outliers, and (iv) to build the confidence intervals.

4.3. Generation of continuous functions from the observations

In this step, the discretized data of each observation are used in order to fit a continuous function or curve using B-splines functions (Ramsay and Silverman, 2005). Figure 7 represents the relationship between the applied load and the displacement at the middle of the span of 1, 2 and $n_j$ observations of the bending test in a SFRC beam. Each relationship of an observation $j$ contains $n_i$ couples of the applied load $P_i^{(j)}$and of the displacement at the middle of the span ${\delta }_i^{(j)}$, where $i=\mathrm{1,2,}\dots ,n_i$. The curve $y^{(j)}$ is the one that minimize the least squares criterion with respect to the couples $\left\{({\delta }_1^{(j)},P_1^{(j)}),\cdots ,({\delta }_i^{(j)},P_i^{(j)}),\cdots ,({\delta }_{n_i}^{(j)},P_{n_i}^{(j)})\right\}$. The curves fitted of these $n_j$ observations can be grouped into the sample $X$, thus:

$X=\left\{y^{(1)},y^{(2)},\dots ,y^{(j)},\dots ,y^{(n_j)}\right\}$ (17)

Figure 7
Relationship between the applied load P and the displacement at the middle of the span $\delta$ as functional data.

4.4. Sorting of the sample and measurement of the centrality of the curves

In this step, the depth of each curve with respect to sample is computed using the graphical approach named Modified Band Depth (MBD) (López-Pintado and Romo, 2009) . Next, the curves are sorted decreasingly according to the depth value. The depth is a measure of the centrality of each curve with respect to the rest of the sample, where the curve with the greater depth is the more central curve. Following the MBD method is described.

Let $I$ be the compact interval where the functional relationships between load and displacement are defined. Likewise, let $R(I)$be the set of these continuous functions into the compact interval $I$. The plot of the function $y^{(j)}$is the subset into the space ${ℝ}^2$ defined by

$G(y^{(j)})=\left\{(t,y^{(j)}(t)):t\in I\right\}$ (18)

Particularly in this work, the band is the region delimited by two or three curves of the sample $X$. In a sample of $n_j$ curves, the number of possible bands is the combination of $n_j$ in 2 curves plus the combination of $n_j$ in 3 curves, that is $n_r=C_2^{n_j}+C_3^{n_j}$. A band delimited by two curves is defined as (López- Pintado and Romo, 2007):

$B_2^{(r)}(y(t))=\left\{{(t_i,y_i^{})}_{i=\mathrm{1,}{n}_i}:t_i\in I,\min \left[y_i^{(r_a)},y_i^{(r_b)}\right]\le y_i^{}\le \max \left[y_i^{(r_a)},y_i^{(r_b)}\right]\right\}$ (19)

and a band delimited by three curves is equal to:

$B_3^{(r)}(y(t))=\left\{{(t_i,y_i^{})}_{i=\mathrm{1,}{n}_i}:t_i\in I,\min \left[y_i^{(r_a)},y_i^{(r_b)},y_i^{(r_c)}\right]\le y_i^{}\le \max \left[y_i^{(r_a)},y_i^{(r_b)},y_i^{(r_c)}\right]\right\}$ (20)

where $r=\mathrm{1,}\dots ,n_r$ is the indicator of the combination of curves $y^{(r_a)},y^{(r_b)}$ for the band $B_2^{(r)}$ and curves $y^{(r_a)},y^{(r_b)},y^{(r_c)}$ for the band $B_3^{(r)}$. The super index $r_a,r_b,\dots$ indicates the number of the curve corresponding to the combination. Figure 8 shows the possible bands of a sample with $n_j=4$ curves. Particularly, the grey region in the subfigures (a) to (f) shows the band delimited by two curves, and in the subfigures (g) to (i) shows the band delimited by three curves. In Figure 8(a), the band is delimited by the curves $y^{(1)}$ and $y^{(2)}$.

For example, curve 4 is partially contained into the band delimited by curves 1 and 3 shown in Figure 8(b) and it is not contained into the band delimited by curves 1 and 2 shown in Figure 8(c). Curve 1 is fully contained into the band delimited by curves 2 and 3 shown in Figure 8(d) and into the band delimited by curves 2 and 4 shown in Figure 8(e).

The Modified Band Depth (MBD) of a curve in the sample is equal to the fraction of the curve contained into each possible band of the sample. This measure allows to rank each curve into the sample, depending on how central a curve is respect to the rest of the sample (Sun and Genton, 2011) .

Modified Band Depth (MBD) (López-Pintado and Romo, 2009) measures the depth of the curves partially contained into a band, by means of the ratio between Lebesgue’s measure of the curve segment contained into the band and Lebesgue’s measure of the whole curve. Considering that the bands are delimited by two or three curves, the modified band depth of a curve $j$ is equal to:

$MBD(y^{(j)})=\frac{1}{C_2^{n_j}}{\displaystyle \sum _{r=1}^{C_2^{n_j}}\frac{\lambda (A(y^{(j)},B_2^{(r)}))}{\lambda (G(y^{(j)}))}}+\frac{1}{C_3^{n_j}}{\displaystyle \sum _{r=1}^{C_3^{n_j}}\frac{\lambda (A(y^{(j)},B_3^{(r)}))}{\lambda (G(y^{(j)}))}}$ (21)

where $A(y^{(j)},B_2^{(r)})$is the subset of $y^{(j)}$ contained into the band $B_2^{(r)}$ and $A(y^{(j)},B_3^{(r)})$ is the subset of $y^{(j)}$ contained into the band $B_3^{(r)}$. The operator $\lambda (•)$indicates the Lebesgue’s measure.

Each curve of the sample has an associated value of modified band depth (MBD). The curve with maximum value of MBD is placed nearer to the center of the sample and the curve with minimum value of MBD is located further from this zone. In this procedure, the curves are sorted decreasingly into the sample according to the value of MBD. The curve with the maximum value of MBD corresponds to the sample median curve (López-Pintado and Romo, 2009). When two or more curves have the maximum MBD, the sample median curve is estimated as the average of the trajectories of these curves.

Figure 8
Possible bands of a sample with 4 curves: (a) to (f) bands delimited by 2 curves, and (g) to (i) bands delimited by 3 curves.

4.5. Detection and removal of outliers

In this step, the outliers are detected and removed by means of the functional boxplot (Sun and Genton, 2011). A functional outlier is a curve generated by a stochastic process with a different distribution than the one of the normal curves (Febrero et al. 2008). To identify and to remove outliers is important because any statistical technique may lead to misleading results and conclusions when the sample has statistical contaminated data.

The functional boxplot is a graphical tool that allows to identify and to remove outliers, when these are totally or partially on the outside of 1.5 times the limits of the 50% of the central region (Sun and Genton, 2011). The central region envelops the 50% of the deeper curves.

For a sample with 8 curves, Figure 9(a) shows the median curve and the central region of 50% of the curves. Figure 9(b) shows the median curve, the outliers and 1.5 times the central region of 50% of the curves.

Figure 9
Functional boxplot of a sample of 8 curves: (a) envelope of the 50% of the curves with greater MBD and median curve, and (b) 1.5 times the envelope of the 50% of the curves and the outliers.

4.6. Evaluation of the confidence intervals

In this step, a confidence interval of the median curve using bootstrap procedure is calculated (Guevara and Vargas, 2013) . First, the bootstrap method creates a random set of ghost samples from the original sample. Next, the curves into each ghost sample are sorted according to their MBD and a central region with 95% of the curves is defined (CR0.95). A set of upper limits that envelope the top of the (CR0.95) associated to the ghost samples is sorted using the MBD. A new central region with 95% of the upper limits is created and its bottom envelope is defined as the upper limit of the confidence interval. Likewise, a set of lower limits of all the ghost samples is sorted with respect to MBD. A new central region with 95% of the lower limits is created and its top envelope is defined as the lower limit of the confidence interval (Guevara and Vargas, 2013).

The bootstrap method consists on the construction of $n_k$ ghost samples ${\tilde{X}}^{(k)}$ from the random resampling with replacing of the original sample $X$of the form (Hesterberg, 2011):

$\begin{array}{l}{\tilde{X}}^{(1)}=\left\{{\tilde{y}}^{(1)(1)}(t),{\tilde{y}}^{(1)(2)}(t),\dots ,{\tilde{y}}^{(1)(j)}(t),\dots ,{\tilde{y}}^{(1)(n_j)}(t)\right\}\\ ⋮\\ {\tilde{X}}^{(k)}=\left\{{\tilde{y}}^{(k)(1)}(t),{\tilde{y}}^{(k)(2)}(t),\dots ,{\tilde{y}}^{(k)(j)}(t),\dots ,{\tilde{y}}^{(k)(n_j)}(t)\right\}\\ ⋮\\ {\tilde{X}}^{(n_k)}=\left\{{\tilde{y}}^{(n_k)(1)}(t),{\tilde{y}}^{(n_k)(2)}(t),\dots ,{\tilde{y}}^{(n_k)(j)}(t),\dots ,{\tilde{y}}^{(n_k)(n_j)}(t)\right\}\end{array}$ (22)

In each ghost sample, the MBD of the curves is calculated and these curves are sorted decreasingly with respect to MBD. Next, the envelope of the $100(1-\alpha )\%$ of the curves with greater MBD is made and the curves on the outside of this region are excluded. In this work $\alpha =0.05$and the confidence $100(1-\alpha )\%$ is equal to 95%. This envelope reduces the size of each ghost sample to $n_l$curves, which is equal to the integer part of $0.95n_j$, and renumbers the curves from greater to smaller value of MBD, that is $MBD({\tilde{y}}^{(1)})\ge MBD({\tilde{y}}^{(2)})\ge \dots \ge MBD({\tilde{y}}^{(l)})\ge \dots \ge MBD({\tilde{y}}^{(n_l)})$ with the index $l=\mathrm{1,}\dots ,n_l$.

The upper limit ${\tilde{g}}_u^{(k)}(t)$ and the lower limit ${\tilde{g}}_b^{(k)}(t)$ of the envelope are functions formed by the segments of the curves of the ghost sample ${\tilde{X}}^{(k)}$ with the maximum and minimum values, respectively. These functions can be expressed as:

${\tilde{g}}_u^{(k)}(t)=\left\{(t,\tilde{y}(t)):t\in I,{\max }_{l=\mathrm{1,}\dots ,n_l}({\tilde{y}}_{}^{(k)(l)}(t))\right\}$ (23)

${\tilde{g}}_b^{(k)}(t)=\left\{(t,\tilde{y}(t)):t\in I,{\min }_{l=\mathrm{1,}\dots ,n_l}({\tilde{y}}_{}^{(k)(l)}(t))\right\}$ (24)

The collection of the upper limits of the ghost samples $Z_u$and another collection of the lower limits of the ghost samples $Z_b$ are constructed, thus:

$Z_u=\left\{{\tilde{g}}_u^{(1)}(t),\dots ,{\tilde{g}}_u^{(k)}(t),\dots ,{\tilde{g}}_u^{(n_k)}(t)\right\}$ (25)

$Z_b=\left\{{\tilde{g}}_b^{(1)}(t),\dots ,{\tilde{g}}_b^{(k)}(t),\dots ,{\tilde{g}}_b^{(n_k)}(t)\right\}$ (26)

The curves of the collection $Z_u$ are sorted according to their MBD. Next, 95% of the curves with greater MBD of $Z_u$ is extracted and a new envelope is defined between the upper limit $h_{uu}^{}(t)$ and the lower limit $h_{ub}^{}(t)$, named upper envelope. The same procedure is applied to the collection $Z_b$, obtaining an envelope between the upper limit $h_{bu}^{}(t)$ and the lower limit $h_{bb}^{}(t)$, named lower envelope.

Finally, the confidence interval ${\tilde{C}}_{0.95}^{}$ is defined as the region delimited by the functions $h_{uu}^{}(t)$ as upper limit and $h_{bb}^{}(t)$ as the lower limit, that is:

${\tilde{C}}_{0.95}^{}=\left\{(t,\tilde{y}(t)):t\in I,h_{bb}^{}(t)\le \tilde{y}(t)\le h_{uu}^{}(t)\right\}$ (27)

5 EXAMPLES OF APPLICATION

SFRC specimens subjected to tension and bending are tested and simulated with the proposed numerical model. The estimated structural response and the representative crack pattern are obtained in the numerical model and compared to the experimental results. The numerical simulation considers plane stress condition, infinitesimal strains and applied static loads and represents the formation of a crack as the concentration of the field displacement isolines. In these examples of application, 200 observations with different orientations of fibres are simulated for each fibre volumetric ratio, to obtain the confidence interval for the median of the relationship between the applied load $P$ and the representative displacement $\delta$.

These sizes are the same that were used to estimate the dispersion coefficient $c_h^{dis}$ (Lamus, 2015). The element size has a low influence over the mean curve of the mechanical response Figure 10 (a).

To determine the influence of the finite element size in the results and to select the appropriate mesh size, a convergence study was performed over the meshes TM8, TM12, TM16 and TM20 with maximum element size of 8mm, 12mm, 16mm and 20mm, which correspond to 16%, 24%, 32% and 40% of the fibre length, respectively (Figure 10(b)). It is necessary to note that the geometry of the simulated specimens requires the use of small-sized elements mainly in the areas of the notch.

Figure 10
Influence of finite element size: (a) Curves load vs displacement for four finite element sizes (b) Four meshes used for the mesh analysis.

The percentages of outlier simulations presented in Table 2 allow identifying the efficiency of the model to obtain feasible results. For larger mesh sizes the number of outliers increased consider-ably. It was established that the ratio of admissible outliers with respect to the efficiency of the model is 5/100. The maximum element size used in this work was selected to be around 32% (TM16) of the length of the fibres used in the experimental tests. This limit corresponds to 16% of the structural member minor dimension in the notched region.

Table 2
Characteristics of the mesh used in estimating the influence of the finite element size.

Mesh	Maximum Element Size	Percentage of the fibre Length		Elements	Nodes	Outliers ratio
		Length=50mm	Length=60mm
TM8	8mm	16.00%	13.33%	2672	1426	2/200
TM12	12mm	24.00%	20.00%	1224	678	4/200
TM16	16mm	32.00%	26.67%	740	426	7/200
TM20	20mm	40.00%	33.33%	500	296	13/200

In each sample of 200 observations obtained for both, tension and bending specimens, less than four curves were considered as outliers using the functional boxplot and consequently, they were excluded. Due to the low probability of atypical events in the numerical simulation using the selected mesh size, the number of outliers does not exceed 2%. The behaviour of this curves can be associated to: (i) longitudinal bands of finite elements where the fibre orientation is parallel to the maximum principal stress direction, (ii) finite elements where the fibre orientation is perpendicular to the maximum principal stress direction at the cracking zone.

5.1. Direct tensile tests

Six cylindrical specimens of plain concrete tested according to NTC454 (ICONTEC, 2011), NTC673 (ICONTEC, 2010) and NTC722 (ICONTEC 2000), presented a compressive strength of 22MPa and an indirect tensile strength of 1.72MPa. The Young’s modulus estimated from the compressive strength is equal to 18GPa. The Poisson’s ratio and the fracture energy were assumed equal to 0.2 and 100N/m, respectively.

Steel fibres Dramix ® RL 45/50 BN were mixed with the plain concrete of the specimens. These fibres are 50mm long and have a shape factor of 45 and an ultimate tensile strength of 1115MPa, according to manufacturer. Five fibres are individually tested to tension, obtaining Young’s modulus of 200GPa and a yielding stress of 520MPa.

The slipping of the fibre embedded into the concrete matrix is defined by means of three repeats of the pull-out test Figure 11), which indicate the relationship between the applied load and the displacement. These curves allow to obtain the parameters of the deformable-sliding fibre model as pull-out strength of 494MPa, equivalent elasticity modulus of 24.56GPa and an equivalent softening parameter of 1564MPa.

Figure 11:
Three Load vs displacement curves for Pullout test using fibres 45/50 and concrete with a compressive strength of 22MPa.

Fifteen specimens of concrete reinforced with 0.5%, 1.0% and 1.5% of volumetric ratio of steel fibres were tested under direct tension as shown in Figure 12(a). This experiment is based on works of other authors (Lim et al. 1987), considering different dimensions of the specimen (Figure 12(b)).

This direct tensile test of the three volumetric ratios of fibres is simulated with the proposed numerical model (Figure 12(a) and Figure 12(b)). The finite element mesh has 568 nodes and 1016 linear triangular elements with a thickness of 150mm and a maximum size of 12mm, as shown in Figure 12(c). The applied loads and the boundary condition are shown in Figure 12(c).

Figure 12
Direct tensile test: (a) sketch of the test, (b) geometry of the specimen, (c) finite element mesh with boundary conditions.

The median and its confidence interval are calculated and compared to the experimental results. For a volumetric ratio of steel fibres equal to 0.5%, Figure 13(a) shows 200 estimated curves using B-splines and Figure 13(b) presents the five experimental curves, the median curve and its 95% confidence interval. The range of the maximum load, between 26kN and 33kN, in the numerical model and in the experimental test is the same. The confidence interval of the numerical result is like the envelope of the experimental curves. The estimation of the structural response of the SFRC with 1.0% and 1.5% is also satisfactory, as shown in Figure 13(c) and Figure 13(d). The maximum load increases slightly with the amount of the fibres, besides; the ductility of the specimens with 1.5% of fibres is significantly greater than that of the specimens with 0.5% of fibres.

The concentration of displacement isolines of a numerical simulation matches with the crack pattern observed in the experimental test as shown in Figure 14. The simulations of the other observations indicate the same result.

Figure 13
Direct tensile test: (a) load vs displacement curves for 0.5% volumetric ratio of fibres, (b) confidence interval, median curve and experimental curves for 0.5% volumetric ratio of fibres, (c) confidence interval, median curve and experimental curves for 1.0% volumetric ratio of fibres, and (d) confidence interval, median curve and experimental curves for 1.5% volumetric ratio of fibres.

Figure 14
Direct tensile test of SFRC with 0.5% of volumetric ratio of fibres in the last loading step: (a) cracking pattern in the test, (b) displacement isolines of the numerical simulation of an observation.

5.2. Three-point bending tests

Three points beams were tested to bending in previous works (Segura and Lamus, 2014), according with RILEM TC 162 standards (Vandewalle et al, 2002).. The concrete is reinforced with 1% and 2% of volumetric ratio of steel fibres. The dimensions, the supports and the applied load are shown in Figure 15(a). The beam has a notch of 25mm of depth at the mid-span, which induces the formation of the crack. This bending test is simulated with the proposed numerical model. The finite element mesh has 781 nodes and 1425 linear triangular elements with a thickness of 150mm and a maximum size of 12mm, as shown in Figure 15(b).

Figure 15
Three-point bending tests: (a) general sketch and (b) finite element mesh.

Six cylindrical specimens of plain concrete to NTC454 (ICONTEC, 2011), NTC673 (ICONTEC, 2010) and NTC722 (ICONTEC 2000), presented a compressive strength of 25.4MPa and an indirect tensile strength of 1.87MPa. The Young’s modulus estimated from the compressive strength was equal to 19GPa. The Poisson’s ratio and the fracture energy were assumed equal to 0.2 and 100N/m, respectively.

Steel fibres Dramix ® RC 80/60 BN were mixed with the plain concrete of the specimens. These fibres have a length of 60mm, a shape factor of 80 and an ultimate tensile strength of 1050MPa, according to manufacturer. Young’s modulus and the yielding stress are assumed equal to 200GPa and 590MPa, respectively.

The slipping of the fibre embedded into the concrete matrix is defined by means of three repeats of the pull-out test, which indicates the relationship between the applied load and the displacement. These curves allow to obtain the parameters of the deformable-sliding fibre model as pull-out strength of 825MPa, an equivalent elasticity modulus of 45GPa and an equivalent softening parameter of 1710MPa (Figure 16).

Figure 16
Three load vs displacement curves for Pullout test using fibres 80/60 and concrete with a compressive strength of 25.4MPa.

The median and confidence intervals of these curves are calculated and compared to the experimental results. Figure 17(a) shows the 200 obtained curves and Figure 17(b) presents the 95% confidence interval for the population median curve, the sample median curve and the three experimental curves for a volumetric ratio of steel fibres equal to 2.0%. The range of the maximum load is between 21kN and 38kN in the numerical result and between 29kN and 37kN in the experimental result. This indicates a suitable value of the upper limit, but a low approximation of the lower limit of the confidence interval. Likewise, the post-pick branch of the upper limit curve of the confidence interval is similar to the upper limit of the experimental envelope, but it is not alike at the lower limit. The estimation of the structural response of the SFRC with 1.0% is also acceptable, as shown in Figure 17[REMOVED REF FIELD](c).

Figure 17
Three-point bending tests: (a) load vs displacement curves for 2.0% volumetric ratio of fibres, (b) confidence interval, median curve and experimental curves for 2.0% volumetric ratio of fibres, (c) confidence interval, median curve and experimental curves for 1.0% volumetric ratio of fibres, (d) final cracking pattern in the test, (e) displacement isolines of the numerical simulation of an observation.

In the last loading step of the experimental test, a crack is formed at the notch and is vertically propagated, as shown in Figure 17(d). The concentration of displacement isolines of a numerical simulation matches with this cracking pattern observed in the experimental test as shown in Figure 17(e). The simulation of other observations indicates the same result.

6 CONCLUSIONS

The novelty of this work is the combination of a deterministic numerical procedure in order to obtain the structural response of SFRC members and a stochastic procedure in order to consider the random orientation and distribution of the steel fibres. This allows establishing the confidence interval of the statistical structural response from numerical simulations.

The orientation and distribution of the steel fibres in concrete are random. Consequently, the structural response of SFRC members is a statistical result that depends on the confidence interval. The proposed numerical model defines this statistical measure. Likewise, this idea could apply to mechanical problems with other random parameters.

The structural response of steel fibres reinforced concrete members is computed by means of a two-dimensional numerical model based on the continuum strong discontinuity approach and the mixture theory. A functional data analysis of the sample of structural response obtained for random orientations of the fibres is proposed. This model predicts the fracture process in members of concrete reinforced with short fibres of random orientation (SFRC), subjected to static load, plane stress condition and infinitesimal strain. The structural response and the cracking pattern of a specimen with defined orientation of fibres are obtained from a non-linear finite elements analysis. This procedure, named observation, is repeated with random orientations of fibres. The functional data analysis takes the structural response of the observations and estimates the median curve and its confidence interval.

The formulation of a constitutive model of SFRC at a material point integrates the behaviour of each component material and its interaction effects. This avoids the high computational cost associated to the discretization of the matrix and the fibres and conserves the random character of the fibres.

The structural response of each observation is a function between the applied load and the representative displacement. These functions are the statistical objects, which have been processed. The presented numerical model uses the multivariate statistics and particularly the functional data analysis, which allows estimating the median curve and its variability.

The numerical simulation of a tensile specimen and a bending beam, using the proposed model, shows satisfactory results with respect to the experimental tests. The gathering of displacement isolines of the numerical simulation is similar to the cracking pattern of the experimental test. Likewise, the estimation of the structural response of the numerical model is alike to the envelope of the experimental results.

The mixture theory allows to simulate mechanical problems where the SFRC has different global volumetric ratios of fibres. Consequently, this numerical model can show that the ductility of a structural member increases with the amount of fibres.

The Functional Data Analysis used in this work requires a high number of numerical analyses to obtain the median curve of the structural response and its confidence interval. It was found that the statistical results for samples with more than 100 observations (analysis) were convergent. In each test presented in this work, around 200 observations with different orientations of fibres were realized. In addition, each analysis demands a low computational cost: around 120 seconds using a 2.40 GHz Intel Core i7-4700M processor with 16GB of RAM for a model with 781 nodes and 1425 elements. However, the procedure has the advantage over less time-consuming methods that use deterministic approaches because its capability to estimate of the range of the experimental response.

REFERENCES

• Antoniadis, A. and Sapatinas, T. (2007). Estimation and inference in functional mixed-effects models. Computational Statistics and Data Analysis, 51: 4793-4813.
• ASTM (2014a). ASTM C39 / C39M-14a, Standard Test Method for Compressive Strength of Cylindrical Concrete Specimens. West Conshohocken, PA (USA): ASTM International.
• ASTM (2014b). ASTM C469 / C469M-14, Standard Test Method for Static Modulus of Elasticity and Poisson’s Ratio of Concrete in Compression. West Conshohocken, PA (USA): ASTM International.
• ASTM (2004). ASTM C496 / C496M-11, Standard Test Method for Splitting Tensile Strength of Cylindrical Concrete Specimens. West Conshohocken, PA (USA): ASTM International.
• Del Barrio, E., Cuesta - Albertos, J., Fraiman, R., Matran, C. (2007). The random projection method for goodness of fit for functional data. Computational Statistics & Data Analysis, this issue, doi: 10.1016/j.csda.2006.09.007.
  » https://doi.org/10.1016/j.csda.2006.09.007
• Bazant, Z. and Planas, J. (1998). Fracture and size effect in concrete and other quaisbrittle materials. (New York) CRC Press.
• De Borst, R., Remmers, J., Needleman, A., Abellan, M. A. (2004). Discrete vs smeared crack models for concrete fracture: bridging the gap. International Journal for Numerical and Analytical Methods in Geomechanics, 28: 583-607.
• Bui, T. T. (2007). Numerical simulation of fracture in plain and fibre-reinforced concrete. University of New South Wales. (Sydney)
• Caggiano, A., Etse, G., Martinelli, E. (2012). Zero-thickness interface model formulation for failure behavior of fiber-reinforced cementitious composites. Computers & Structures (on-line), 98-99: 23-32. http://www.sciencedirect.com/science/article/pii/S0045794912000284
  » http://www.sciencedirect.com/science/article/pii/S0045794912000284
• Cardot, H., Crambes, C., Kneip, A., Sarda, P. (2007). Smoothing spline estimators in functional linear regression with errors-in-variables. Computational Statisctics Data Analysis, 51: 4832-4848.
• Cardot, H., Ferraty, F., Sarda, P. (1999). Functional linear model. Stat. Probab. Lett., 45: 11-22.
• Cuevas, A., Febrero, M., Fraiman, R. (2002). Linear functional regression: the case of fixed design and functional response. Can. J. Stat., 30: 285-300.
• Cunha, V. M. C. F., Barros, J. A. O., Sena-Cruz, J. M. (2012). A finite element model with discrete embedded elements for fibre reinforced composites. Computers & Structures (on-line), 94-95: 22-33. http://www.sciencedirect.com/science/article/pii/S0045794911003038
  » http://www.sciencedirect.com/science/article/pii/S0045794911003038
• Dabo-Niang, S. and Rhomari, N. (2009). Kernel regression estimate in a banach space. J. Stat. Plan. Infer., 131.
• Dias, I. F., Oliver, J., Huespe, A. (2011). Strain injection, mixed formulation and strong discontinuities in fracture modeling of quasi-brittle materials. In: Congress in numerical methods in engineering. Coimbra.
• Febrero, M., Galeano, P., González-Manteiga, W. (2008). Outlier detection in functional data by depth measures, with application to identify abnormal NOx levels. Environmetrics, 19: 331-345.
• Ferraty, F. and Vieu, P. (2006). Nonparametric Functional Data Analysis. (Berlin): Springer-Verlag.
• Ferraty, F. and Vieu, P. (2002). The functional nonparametric model and application to spectrometric data. Computational Statistics, 17: 545-564.
• Fletcher, T., Lu, C., Pizer, S., Joshi, S. (2004). Principal Geodesic Analysis for the Study of Nonlinear Statistics of Shape. IEEE Transactions on medical imaging, 23: 995.
• Gal, E. and Kryvoruk, R. (2011). Meso-scale analysis of FRC using a two-step homogenization approach. Computers & Structures (on-line), 89: 921-929. http://www.sciencedirect.com/science/article/pii/S0045794911000393
  » http://www.sciencedirect.com/science/article/pii/S0045794911000393
• Guevara, R. D., Vargas, J. A. (2013). Process capability analysis for nonlinear profiles using depth functions. Quality and Rebiality Engineering International.
• Hesterberg, T. (2011). Bootstrap. WIREs Computational Statistics, 3: 497-526.
• ICONTEC (2011), NTC 454, Ingeniería y arquitectura. Concreto Fresco. Toma de muestras, (Bogotá).
• ICONTEC (2010), NTC 673, Ensayo de resistencia a la compresión de especímenes cilíndricos de concreto, (Bogotá).
• ICONTEC (2000), NTC 722, Concreto. Método de ensayo para determinar la resistencia a la tensión indirecta de especímenes cilíndricos de concreto, (Bogotá).
• Kang, J., Kim, K., Lim, Y. M., Bolander, J. E. (2014). Modeling of fiber-reinforced cement composites: Discrete representation of fiber pullout. International journal of solids and structures (on-line), 51: 1970-1979. http://www.sciencedirect.com/science/article/pii/S0020768314000547
  » http://www.sciencedirect.com/science/article/pii/S0020768314000547
• Lamus, F. A. (2015). Modelo numérico del comportamiento inelástico del concreto reforzado con fibras cortas de acero . Bogotá,: National University of Colombia.
• Laranjeira, F., Molins, C., Aguado, A. (2010). Predicting the pullout response of inclined hooked steel fibers. Cement and concrete research, 40: 1471-1487.
• Leite, J. P. B., Slowik, V., Apel, J. (2007). Computational model of mesoscopic structure of concrete for simulation of fracture processes. Computers & Structures (on-line), 85: 1293-1303. http://www.sciencedirect.com/science/article/pii/S0045794906003579
  » http://www.sciencedirect.com/science/article/pii/S0045794906003579
• Li, F. and Li, Z. (2001). Continuum damage mechanics based modeling of fiber reinforced concrete in tension. International journal of solids and structures, 38, 5: 777-793.
• Lim, T. Y., Paramasivam, P., Lee, S. L. (1987). Analytical model for tensile behaviour of steel - fibre concrete . ACI Materials Journal, Vol 84, No 4: 286-298.
• Linero, D. L., Oliver, X., Huespe, A. E. (2007). A model of material failure for reinforced concrete via continuum strong discontinuity approach and mixing theory. Barcelona: International Center for Numerical Methods in Engineering.
• Locantore, N., Marron, J. S., Simpson, D. G., Tripoli, N., Zhang, J. T., Cohen, K. L. (1999). Robust principal component analysis for functional. Test, 8: 1-74.
• López-Pintado, S. and Romo, J. (2007). Depth-based inference for functional data. Computational Statistics & Data Analysis 51: 4957-4968.
• López-Pintado, S. and Romo, J. (2009). On the Concept of Depth for Functional Data. Journal of the American Statistical Association, 104: 718-734.
• Luccioni, B., Ruano, G., Isla, F., Zerbino, R., Giaccio, G. (2012). A simple approach to model SFRC. Construction and building materials, 37: 111-124.
• Manzoli, O. L., Gamino, A. L., Rodrigues, E. A., Claro, G. K. S. (2012). Modeling of interfaces in two-dimensional problems using solid finite elements with high aspect ratio. Computers & Structures (on-line), 94-95: 70-82. http://www.sciencedirect.com/science/article/pii/S0045794911002938
  » http://www.sciencedirect.com/science/article/pii/S0045794911002938
• Manzoli, O. L., Oliver, X., Huespe, A., Diaz, G. (2008). A mixture theory based method for three dimensional modeling of reinforced concrete members with embedded crack finite elements. Computers and concrete, 5, 4: 401-416.
• Mora, D. F., Oliver, J., Huespe, A. (2013). Multifield-Based modeling of material failure in high performance reinforced cementitious composites. Monographs of the International Centre for Numerical Methods in Engineering (CIMNE).M139. (Barcelona). https://www.scipedia.com/public/Mora_et_al_2013a
  » https://www.scipedia.com/public/Mora_et_al_2013a
• Naaman, A., Namur, G., Alwan, J., Najm, H. (1991). Fiber pullout and bond slip II. Experimental validation. Journal of Structural Engineering ASCE, 117: 2791-2800.
• Nawy, E. (2009). Reinforced concrete: A fundamental approach (6th edition). London: Prentice Hall.
• Nerini, D. and Ghattas, B. (2007). Classifying densities using functional regression trees: Applications. Computational Statistics & Data Analysis 51: 4984-4993.
• Nofal, M. (1997). Continuum damage mechanics for plain, fibre-reinforced, and reinforced concrete materials and structures. Ottawa (Canada): Department of Civil and Environmental Engineering, Carleton University.
• Oliveira de Barros, J. A. (1995). Comportamento do betão reforçado com fibras, Análise experimental e simulação numérica. Porto: Faculdade de Engenharia, Universidade do Porto.
• Oliver, J. (1996a). Modelling strong discontinuities in solid mechanics via strain softening constitutive equations. Part I: Fundamentals. International Journal for Numerical Methods in Engineering, 39: 3575-3600.
• Oliver, J. (1996b). Modelling strong discontinuities in solid mechanics via strain softening constitutive equations. Part II: Numerical Simulation. International Journal for Numerical Methods in Engineering, 39: 3601-3623.
• Oliver, J., Cervera, M., Oller, S., Lubliner, J. (1990). Isotropic damage models and smeared crack analysis of concrete. In: N.B. et al. (ed.) SCI-C Computer Aided Analysis and Design of Concrete Structures. Swansea: N. Bicanic and H. Mang: 945-958.
• Oliver, J. and Huespe, A. (2004a). Continuum approach to material failure in strong discontinuity settings. Computer Methods in Applied Mechanics and Engineering, 193: 3195-3220.
• Oliver, J. and Huespe, A. (2004b). Theoretical and computational issues in modelling material failure in strong discontinuity scenarios. Computer Methods in Applied Mechanics and Engineering, 193: 2987-3014.
• Oliver, X., Linero, D. L., Huespe, A. E., Manzoli, O. L. (2008). Two-dimensional modeling of material failure in reinforced concrete by means of a continuum strong discontinuity approach. Computer Methods in Applied Mechanics and Engineering, 197: 332-348.
• Oller, S. (2003). Simulación numérica del comportamiento mecánico de materiales compuestos. (Barcelona): CIMNE.
• Park, R. and Paulay, T. (1975). Reinfored Concrete Structures. New York: John Wiley & Sons, Inc.
• Peng, X. and Meyer, C. (2000). A continuum damage mechanics model for concrete reinforced with randomly distributed short fibers. Computers and structures, 78, 4: 505-515.
• Prós, A. (2011). Numerical approach for modeling steel fiber reinforced concrete. (Barcelona): Universitat Politécnica de Catalunya.
• Radtke, F., Simone, A., Sluys, L. J. (2010). A computational model for failure analysis of fibre reinforced concrete with discrete treatment of fibres. Engineering fracture mechanics, 77, 4: 597-620.
• Ramsay, J. and Silverman, B. (2005). Functional data analysis. New York: Springer.
• Rice, J. R. (1976). The Localization of Plastic Deformation. In: ed W. T. Koiter (ed.) Proceedings of the 14th International Congress on Theoretical and Applied. Delft: NorthHolland Publishing Co: 207-220.
• RILEM (1985). Determination of fracture energy of mortar and concrete by means of three-point bend tests of notched beams. RILEM TC-50. Materials and Structures, 18, 4: 287-290.
• Segura, S. M. and Lamus, F. A. (2014). Energía de fractura en modo I del concreto de resistencia normal reforzado con fibras cortas de acero. Revista de la Escuela Colombiana de Ingeniería, 1: 47-56.
• Simó, J. and Hughes, T. H. R. (1998). Computational Inelasticity. New York: Springer-Verlag.
• Soe, T. N. (2010). Fracture processes in steel fibre reinforced concrete. Sydney: School of Civil and Environmental Engineering, The University of New South Wales.
• De Souza, E. A., Peric, D. and Owen, D. (2008). Computational Methods for Plasticity: Theory and Applications. Chichester: John Wiley & Sons.
• Sun, Y. and Genton, M. G. (2011). Functional Boxplots. Journal of Computational and Graphical Statistics, 20: 316-334.
• Tarrio, J. (2012). Evaluación y clasificación de materiales: un enfoque estadístico. Tesis Doctoral. La Coruña: Universidade da Coruña.
• Truesdell, C. and Toupin, R. (1960). The classical field theories. Springer. (Berlín).
• Vandewalle, L., Nemegeer, D., Balazs, L,. Barr, B., Barros, J., Bartos, P., Banthia, N., Criswell, M., Denarie, E., Di Prisco, M., Falkner, H., Gettu, R., Gopalaratnam, V., Groth, P., Hausler, V., Kooiman, A., Kovler, K., Massicotte, B., Mindess, S., Reinhardt, HW., Rossi, P., Schaerlaekens, S., Schumacher, P., Schnutgen, B., Shah, S., Skarendahl, A., Stang, Henrik., Stroeven, P., Swamy, R., Tatnall, P., Teutsch, M., Walraven, J. (2003), RILEM TC-162 Test and design methods for steel fibre reinforced concrete - sigma-epsilon-design method - Final recommendation. Materials and Structures, 36, 262: 560-567.
• Yu, R. C., Cifuentes, H., Rivero, I., Ruiz, G., Zhang, X. (2016). Dynamic fracture behaviour in fibre-reinforced cementitious composites. Journal of the Mechanics and Physics of Solids (on-line). http://www.sciencedirect.com/science/article/pii/S002250961530394X
  » http://www.sciencedirect.com/science/article/pii/S002250961530394X
• Zile, E. and Zile, O. (2013). Effect of the fiber geometry on the pullout response fo mechanically deformed steel fibers. Cement and concrete research, 44: 18-24.

• Available online: March 11, 2019.

Publication Dates

• Publication in this collection
  08 Apr 2019
• Date of issue
  2019

History

• Received
  24 Oct 2018
• Reviewed
  28 Feb 2019
• Accepted
  08 Mar 2019