1. Introduction
Textile reinforced concrete (TRC) and fabric reinforced mortars (FRM) are cementitious composites reinforced with continuous, two-dimensional, or three-dimensional textiles/fabrics made of carbon, alkali-resistant glass, or polymer multifilament yarns [
1,
2,
3,
4]. These composites feature controlled multiple cracking under increasing deformation prior to reaching the tensile strength of the loaded yarns and are suitable as externally applied retrofitting and strengthening layers on concrete and masonry elements [
5]. The mechanical and physical properties of the reinforcing yarns allow the achievement of effective and durable strengthening properties in narrow thicknesses of up to 20 mm.
Although highly promising for applications involving ordinary service conditions, the energy dissipation capacity and damage tolerance of TRC/TRM require further enhancement for scenarios involving dynamic or repeated mechanical loading with extensive deformations and high energy input [
6]. This can be achieved through the targeted addition of short fibers, which can increase the cracking stress of the cementitious matrix, decrease the crack spacing and crack width, and mitigate matrix spalling [
7,
8].
Another promising material for such strengthening applications exists in the form of strain-hardening cement-based composites (SHCC), also known as engineered cementitious composites—ECC. SHCC are made of fine-grained cementitious matrices and short polymer micro-fibers, and they exhibit a strain-hardening tensile behavior with pronounced multiple cracking prior to failure localization [
9]. Despite their high strain capacity, high energy dissipation potential, and damage tolerance, SHCC demonstrate relatively low tensile strength, and their mechanical performance in real-scale applications can be negatively affected by production and application techniques [
10] and boundary conditions and structural size as well [
11,
12]. Furthermore, the multiple cracking in SHCC degrades stiffness, a development that can limit the degree of confinement offered by such strengthening layers without continuous reinforcement.
In this context, the combination of SHCC with textile reinforcement can effectively mitigate the disadvantages of each composite in particular and ensure superior mechanical features both under quasi-static [
13,
14] and impact loads [
15,
16]. However, aside from the mechanical features, the sustainability factors, especially the economical and the ecological ones, are essential with regard to the real-scale applicability of such materials. An economical material design requires the purposeful selection of the constituents and adjustment of their interaction with respect to specific performance requirements. Given that the experimental investigations of such composites imply high financial costs and require advanced technical equipment, high-fidelity numerical simulations are essential in enabling a wide and detailed study of various hypotheses and material parameters. Furthermore, these should serve as a reliable basis for assessing the strengthening performance of such composites at the structural level.
In order to analyze and predict accurately the influence of composition, application, and boundary conditions on the strength, deformation, and fracture properties of fiber-reinforced composites, the material models should imply reinforcement discretization both for textile [
17,
18] and short fibers [
19]. Finite element (FE) [
19] and discrete lattice models with discrete fiber representation [
20,
21,
22] imply the explicit consideration of the spatial inhomogeneity and anisotropy related to fiber distribution and orientation, but the substantial computational effort caused by the geometric properties of discrete fibers limits the applicability of this approach to small regions of interest [
23,
24].
For morphologically detailed but computationally efficient numerical simulations of hybrid fiber-reinforced composites, a continuum representation of the multiple cracking in SHCC can be adopted [
25,
26,
27,
28], while the inherent material variability can be accounted for within a probabilistic representation of the relevant material parameters, such as matrix strength and tensile strength [
29]. For applications requiring the accurate modeling of localized fracture, discrete crack models can be adopted or combined with the smeared crack model [
29,
30].
Commonly, the probabilistic approach is applied in phenomenological and numerical models of plain SHCC [
26,
27], while the continuous reinforcement is discretized in the context of steel-reinforced concrete elements [
25]. The paper at hand presents a FE model that combines these two approaches in the context of hybrid fiber-reinforced composites, i.e., the probabilistic approach for SHCC and the discrete representation of the continuous textile yarns. It enables a detailed analysis of the interaction mechanisms and of the effect of various material combinations on the macroscopic composite properties and facilitates an accurate prediction of the multiple cracking and fracture phenomena in the tensioned composites. The presented approach is applicable both to material design and optimization purposes as well as to assessing the structural behavior of elements strengthened with such composites.
The model and the numerical study presented in this paper complement a series of experimental investigations aimed at clarifying the influence of the constituents’ properties and their interaction on the mechanical behavior and strengthening performance of SHCC with textile reinforcement [
13,
14]. Given that the primary focus was laid on the pre-peak tensile properties of such composites, the numerical material model of SHCC implied a smeared crack formulation of both the strain-hardening and localization phases. The mesh-size sensitivity within the smeared crack model was counteracted by the probabilistic representation and spatial fluctuation of the material parameters equivalent to the matrix strength and collective crack-bridging strength in SHCC. The main goal of the probabilistic formulation was a realistic representation of the successive and spatially distributed multiple cracking in SHCC, depending on the inherent material variability in terms of flaw size distribution [
31,
32] and fiber orientation and distribution [
33,
34]. Moreover, the interaction with the continuous reinforcement was another important aspect analyzed with regard to the resulting extent of multiple cracking and crack width [
13,
14]. Following the multi-scale concept [
29,
35], the constitutive law of SHCC reflected the single-crack opening behavior. The adjustment of the probabilistic distribution of the matrix strength and tensile strength in SHCC allowed fitting the effective extent of multiple cracking of the modeled composite specimens.
The continuous reinforcement was discretized by truss elements embedded in the longitudinal (loading) direction of the composite specimens under tension. The influence of stiffness and elongation capacity of the continuous reinforcement, as well as of the bond strength between the continuous reinforcement and the SHCC, were analyzed in a parameter study. Besides the macroscopic tensile stress–strain relationships of the plain SHCC and of the hybrid fiber-reinforced composites, the multiple cracking pattern and the resulting crack width along the entire deformation range were analyzed and compared to the reference experimental results presented in [
13,
14].
2. Finite-Element Model
The finite-element simulations were performed using the software Atena Engineering 2D V5.6.1, by Cervenka Consulting s.r.o (Prague, Czech Republic). The constitutive nature and geometry of the modeled composite specimens allowed performing 2D simulations to limit the computational effort. The geometry and dimensions of the modeled composite specimens are presented in
Figure 1a according to [
13,
14]. The modeled region of interest corresponds to the gauge portion, having a width of 60 mm and length of 150 mm; see
Figure 1b. The 20 mm thickness of the experimental specimens defined the nominal thickness of the solid finite elements in the 2D simulations. The mesh size of 3 mm was defined according to the minimum average crack spacing observed in tensioned composite specimens with hybrid fiber reinforcement [
13,
14]. It was assumed that this crack spacing corresponds to the crack saturation achievable by the material compositions analyzed with the given specimen geometry and boundary conditions. The same applies to the defined stress parameters of the SHCC model.
The hybrid fiber-reinforced specimens accommodated five longitudinal yarns as presented in
Figure 1b. The truss elements were arranged over the nodes of the SHCC mesh, resulting in a spacing of 12 mm, slightly smaller than the spacing of the warp yarns in the actual textile of 12.7 mm [
13]. The degrees of freedom (DOFs) of the truss elements are independent of the DOFs of the underlying solid elements. Given the textile production technology of the reinforcing textiles by stitch bonding, the joints between the warp (longitudinal) and weft (transversal) yarns are neither rigid nor strong. In the context of uniaxial tensile loading, this allowed discretizing only the warp yarns. The weft yarns might contribute to the anchorage strength of the warp yarns, but this was implicitly considered by defining the bond properties based on warp–yarn pullout experiments with attached weft yarns [
14]. Furthermore, the weft yarns might weaken the corresponding cross-sections of the composite specimens, but this effect cannot be captured with transversal truss elements, and it was taken into account by the probabilistic representation of the matrix strength.
In the experimental specimens, the ends of the textile yarns embedded in the anchorage portions, aside from gauge length, were coated with epoxy-resin and sand in order to limit yarn slip. However, the additional coating did not ensure a perfect bond, and the deformation measurements on the gauge portion included the yarn delamination and partial pullout from the anchoring portions [
14]. Thus, disabling the relative slip of the truss elements at the ends of the gauge portion would not accurately reflect the real boundary conditions in the experiments. For this reason, additional anchorage portions with lengths of 15 mm, i.e., five elements, were defined at the ends of the gauge portion; see the darker solid elements in
Figure 1b. The additional solid elements had an identical Young’s modulus to the core SHCC elements but were designed to deform elastically throughout the entire load range. The slip of the truss elements relative to the solid elements was disabled at the boundary nodes.
The model specimen was loaded by a prescribed deformation of 0.1 mm with a step multiplier of 0.03 at the nodes of the top boundary elements, and the reaction forces were recorded and summed at the nodes of the bottom row, as presented in
Figure 1b. Only the deformation of the 150 mm long gauge portion was considered in deriving the macroscopic specimen deformation according to the relative displacement of the solid element nodes 1,2 and 3,4, as presented in
Figure 1b. The standard Newton Raphson method was used as the solution procedure with an iteration number limit of 40.
2.1. SHCC—Finite-Element Formulation
SHCC was discretized with plane quadrilateral elements of type CCIsoQuad, which are is oparametric elements featuring bilinear interpolation and four Gaussian integration points [
36]. The CC3DNonLinCementitious2User material model was applied to model SHCC due to the possibility of model strain-hardening. This material model is a fracture-plastic constitutive model, combining tensile (fracture) and compressive (plastic) behavior. The fracture model is based on the classical orthotropic, smeared crack formulation, and crack band model [
36,
37]. The hardening/softening plasticity model is based on the Menetrey–Willam failure surface. The tensile strain is decomposed in its elastic, plastic, and fracture components, as shown in Equation (1):
with the new stress state computed as in Equation (2):
is the increment of plastic strain, and
is the increment of fracturing strain. The Rankine failure criterion for concrete cracking is applied with stresses converted in the material directions:
where
is the trial stress, and
is the trial tensile strength in the material direction
i. The crack opening
is calculated from the total value of the calculated fracturing strain
in the direction
plus the current increment of the fracturing strain
, the sum being multiplied by the crack band size
as in Equation (4):
is calculated as the size of the element projected onto the crack direction. Given the user-defined, strain-hardening material law, the crack width is adjusted according to a predefined localization strain
and characteristic length, also described as band width limiter,
, as in Equation (5):
The constitutive law of SHCC represents the single-crack opening behavior as derived experimentally. The location of fracture localization in the model specimen resulted from the probabilistic and axial fluctuation of the tensile strength. Given that the probabilistic parameters only varied axially along the specimen, the fracture zone was limited to one transverse row of elements, thus allowing the characteristic length to be defined equal to the element size. In this way, the pre-peak and post-peak response of the SHCC elements strictly followed the corresponding user-defined constitutive laws as presented in
Section 2.2.
Given that the work at hand dealt exclusively with simulations of uniaxial tension experiments, the fixed crack model was applied with the default value of the shear retention factor [
36,
37,
38]. Similarly, the crack closure stiffness was set to its default value, implying unloading behavior parallel to the elastic stiffness. Refining this parameter would be necessary to model the partial crack closure accurately during the formation of new cracks, softening, or due to cyclic loading [
35].
2.2. SHCC—Constitutive Law
The reference SHCC in this study consisted of a high-strength matrix and 6 mm long ultra-high molecular weight polyethylene fibers (UHMWPE) in a volume fraction of 2% [
39,
40]. The composite yields an average compressive strength of 140 MPa and Young’s modulus of 29 GPa [
39,
41]. The average quasi-static tensile strength as derived on the specimen geometry modeled in the current work is approximately 6.7 MPa, and the average macroscopic strain at failure localization is approximately 1.3% [
13,
14]. Smaller SHCC specimens of identical compositions yield considerably lower crack spacing and higher strain capacity [
39,
40], this being traceable back to the strong influence of heterogeneity, size, and structural effects on the apparent material properties of SHCC.
The multilinear constitutive tensile law of SHCC in the post-elastic stage represents a simplified single-crack opening law and was derived by interpolating experimental data from different studies by the authors on notched and even specimens; see
Figure 2 [
13,
14,
42]. The matrix strength
and the tensile strength
are variable material parameters. The ranges indicated in
Figure 2 are partly empirical (lowest values) and partly hypothesized (highest values) due to the lack of experimental evidence regarding the strongest cross-sections in common SHCC specimens. The definition of the input laws of single SHCC elements imposed the condition of strain-hardening, i.e., tensile strength higher than the matrix strength. Given the partial overlap of the corresponding ranges of values, the violation of this condition was prevented by adjusting the axial fluctuation of the matrix strength and tensile strength along the model specimen, as explained in
Section 3. In terms of stress levels, the blue curve in
Figure 2 is representative of plain SHCC specimens according to the rule of the weakest link.
The constitutive strain relationships of the model are derived by dividing the presented deformation (crack opening) values to the element size of 3 mm. The crack opening at peak stress
has a constant value of 100 µm, and it corresponds to an equivalent strain at peak of 3.3%. The softening branch was defined based on
(stress at the kink point) and
(crack width at the kink point), which are ratios of
and
according to Equations (6) and (7). The values of
and
are given in
Table 1. The crack width at vanishing crack-bridging
was defined as half of the fiber length, i.e., 3 mm.
The strain values presented in
Table 1 are only valid for single elements and not for the composite specimens, in which the selective multiple cracking led to considerably lower macroscopic deformations; see
Section 3. Because the finite elements were designed to accommodate the elongation equivalent to one single crack, the term first-crack stress in this paper refers only to the model composite specimens, while the term matrix strength refers to the finite elements.
As facilitated by the axial fluctuation of the matrix strength and tensile strength, 50 monitoring points were set along the model specimen, i.e., one point for every row along the gauge portion, for deriving the crack widths under load. The corresponding crack widths were averaged over all the cracked elements and plotted next to the macroscopic stress–strain curves. In this way, the average crack width in the strain-hardening phase could be compared to that derived in the reference tension experiments.
2.3. Continuous Reinforcement
The reference continuous reinforcement in this study is a carbon textile produced by V. Fraas GmbH (Helmbrechts, Germany) under the brand name TUDALIT-BZT2 [
43] and carrying the index T1 in this work. Its mechanical and geometric properties in the warp, i.e., principal and direction, are summarized in
Table 2 based on the data by the producers and from experimental investigations by the authors [
13,
14,
15].
In the context of structural strengthening against impact loading, the strengthening layers should ensure a strong confinement of the strengthened element [
12] but at the same time exhibit a large capacity for energy dissipation in extreme loading cases. This implies that not only the SHCC but the continuous reinforcement as well should yield high inelastic deformability. As a model reinforcement of high elongation capacity, a UHMWPE textile was produced with identical geometric properties and impregnation material as the reference carbon textile [
14]. The properties of the UHMWPE filaments Dyneema SK75 according to the producer (DSM, Heerlen, The Netherlands) are given in [
44]. In both cases, the multi-filament yarns were impregnated with a watery styrene–butadiene based dispersion for ensuring adhesion and effective composite action among the single filaments within the yarns and with the surrounding cementitious material.
The continuous, unidirectional reinforcement was discretized by truss elements CCIsoTruss (type 2D) of the Reinforcement material model. These are isoparametric elements integrated by Gauss integration at one point with linear interpolation. The reinforcement truss elements have neither bending nor shear stiffness.
Figure 3 presents the tensile constitutive laws of the discrete reinforcement.
Whereas the numerical model T1 represents accurately the properties of the carbon warp yarns, the geometric features and the tensile strength of T2 differ from those presented in [
14]. For a systematic parameter study, in this work, they were set to be equal to the corresponding properties of T1. The interfacial equilibrium condition and stress-state in the embedded truss elements were derived based on the geometric cross-section parameters given in
Table 2.
2.4. Bond
The bond-slip was modeled with bond elements between the truss elements and the solid elements and defined by the Bond for Reinforcement material model [
45]. The defined bond-slip relationships were constant parameters along the entire length of the modeled element and were defined according to the pullout experiments presented in [
14,
15,
46]. These are presented in
Figure 4 as B1 for the carbon textile and B3 for the equivalent polymer textile.
The reference bond-slip relationships assume a simplified linear relationship to the slip strength followed by a stress drop after the critical slip and subsequently by a constant stress pullout. To avoid ill-conditioning, the input bond strength vs. slip law assumed an elastic range with no slip up to bond stress of 0.05 MPa; see
Figure 4.
In the previous study, the effect of bond strength was experimentally investigated by applying an additional coating of epoxy resin and sand to the yarns, which substantially increased the bond strength and altered both the tensile properties and fracture behavior of the fiber-reinforced composites [
14,
47]. In the case of the UHMWPE textile, the additional coating influenced the elongation behavior of the textile itself [
14]. However, the high bond strength caused premature yarn rupture in pullout experiments and did not allow the derivation of the bond-slip relationships. For this reason, the alternative bond properties analyzed in the frame of the numerical parameter study involved an adhesive bond strength twice higher than that of the reference bond B1 and B3, followed by a constant stress pullout with no stress drop after the critical slip displacement; see B2 and B4 in
Figure 4. Although not quantitatively comparable to the experimental data, the numerical results of B3 and B4 offered a qualitative evidence on the effect of bond strength on the tensile behavior of hybrid fiber-reinforced composites. The numerical parameters of the bond-slip relationships are summarized in
Table 3.
Besides the plain SHCC model, the parameter combinations, i.e., numerical material compositions, investigated in this work were: (1) SHCC-T1-B1, equivalent to M-PE-T in [
13]; (2) SHCC-T1-B2, equivalent to M1-PE-TE in [
47]; (3) SHCC-T2-B3; and (4) SHCC-T2-B4. The parameter combinations (3) and (4) are qualitatively equivalent to composites reinforced by a polymer textile, but they do not yield a quantitative comparison.
3. Probabilistic Material Parameters
The inherent material variability of SHCC in terms of flaw size and fiber distribution determines their effective tensile strength and ductility [
48]. The material inhomogeneity in real-scale applications is expected to be considerably more pronounced, and its negative effects can be partly mitigated in combination with continuous reinforcement. Proper material design in this respect should allow exploiting the potential ductility of SHCC and avoid premature crack localization [
12]. In numerical simulations, this effect can be highlighted by involving a probabilistic representation of the matrix strength and tensile strength of SHCC, as presented in
Figure 5.
The strength parameters of plain SHCC as shown by the blue curve in
Figure 5 are determined by the lowest values of the matrix strength and tensile strength and by the shape of the distribution functions of these parameters. Given predefined minimum values, distribution functions of the matrix strength with a lower median should result in a higher extent of multiple cracking than in the case of a larger median. This effect is presented in
Figure 5 by the green hatched area, which increases if the distribution function of the matrix strength yields a lower median.
It should be mentioned that the theoretical minimum crack spacing in plain SHCC is determined by the short fibers’ transmission length, which represents the tensile stress increment in the matrix with increasing distance from the crack flank [
49]. However, the analytically derived transmission length is an idealized parameter and hardly achievable in common SHCC samples.
The Weibull function [
50] was selected to define the probability distribution of the varying material parameters since its shape can be adjusted to match different distributions and demonstrate the above-formulated statements within a parameter study. Additionally, the Weibull distribution is commonly adopted to represent the mechanical strength of materials [
31,
51,
52]. The three-parameter Weibull distribution was analyzed in this work, with the corresponding probability density function
presented in Equation (8) and cumulative distribution function
in Equation (9). The following conditions apply:
The value of the shape parameter , also the slope parameter, controls the frequency distribution of the randomized parameter . The value of the scale parameter , characteristic life, is responsible for the width of the distribution range, given a constant . The location parameter locates the lower bound of the predefined range of .
The probabilistic material properties were generated in a preprocessing module. Firstly, the pseudo-random numbers
that follow a uniform distribution with values from 0 to 1 were generated. Subsequently, the cumulative probability function of the Weibull distribution
was equated to
according to [
53], resulting in the array of numbers
as shown in Equation (10):
The parameters to be defined besides the Weibull parameters are the size of the array with the minimum and the maximum values. In the presented study, the material inhomogeneity was only defined in the longitudinal direction of the model specimens, i.e., the properties did not vary transversally. This was intended to simulate the steady-state cracking in SHCC elements with small cross-sections, which implies instant crack propagation through the entire cross-section shortly after crack initiation. Given that the effective size of the random arrays generated in this work was equal to the number of elements along the gauge portion of the model specimen, i.e., 50, the probability of the smallest generated number to be equal to the predefined inferior boundary was extremely low, especially with larger β. For this reason, the location parameter was adjusted to shift the generated arrays to the predefined inferior boundaries.
The initially generated randomized values and their axially non-correlated fluctuation are schematically presented in
Figure 6a. These sets of numbers were transposed into predefined stress ranges depending on the material parameter to be represented, as shown in
Figure 6b. The next step consisted in the cross-correlation [
54] of the matrix-strength and tensile strength as well as in their axial correlation, as shown in
Figure 6b. The cross-correlation did not imply an explicit interdependence of the matrix strength and tensile strength but mainly targeted a congruent pattern of the green and purple fluctuations in
Figure 6b. In this way, strain-hardening was ensured in all finite elements along the model specimen.
The probabilistic parameters were axially correlated according to a predefined correlation length, which is a range that the algorithm uses to group the adjacent values spatially for steady variation. From the example in
Figure 6b, the correlation length is approximately equal to 45 mm or 15 elements. There is no experimental evidence on the spatial variability of the matrix strength and tensile strength in the reference SHCC samples. However, given that only axial variability was considered in this work, the representation in
Figure 6b is presumably more realistic than that in
Figure 6a, since the axial fluctuation must represent an average of the variability in the other directions, which levels out over the width and thickness of the specimens.
The probabilistic material parameters were assembled into readable input files for the FE software, which included the user-defined constitutive relationships attributed to specific elements according to their location in the model SHCC specimen.
Figure 6c shows the effective number of cracked elements at failure localization corresponding to exemplary probability distribution and axial fluctuation of matrix strength and tensile strength in
Figure 6b.
Note that the axial fluctuation of the probabilistic parameters is irrelevant to the effective extent of multiple cracking in plain SHCC. With no cross-correlation, the predefined regions of multiple cracking, the green hatched areas in
Figure 6b, would be distributed along the specimen as in
Figure 6a. However, their total area would not change, since that only depends on the number of elements having matrix strength lower than the tensile strength of the weakest element in the specimen.
Since the experimental data only provide evidence on the far-left tails of the distribution functions, the probability distributions presented in this section are not entirely empirical and are only intended to highlight the importance of probabilistic numerical modeling of such composites.
5. Conclusions and Outlook
The study presents a finite-element model of strain-hardening cement-based composites (SHCC) with discretized continuous reinforcement. Uniaxial tension experiments were simulated, with boundary conditions and material properties adopted from reference experimental studies on hybrid fiber-reinforced composites [
13,
14,
15,
47] in which the continuous reinforcement consisted of 2D textiles made of carbon and ultra-high molecular weight polyethylene (UHMWPE).
The SHCC material model was based on the smeared crack formulation and the constitutive law represented the single crack-opening behavior derived from experimental results. For an accurate assessment of the effect of continuous reinforcement and bond strength on the extent of multiple cracking and crack width in the hybrid fiber-reinforced composites, the matrix strength and tensile strength of SHCC were defined within a probabilistic framework following various Weibull probability distributions. This allowed for a selective multiple cracking in the model composites based on predefined distribution functions and axial fluctuation along the model specimens. The effect of the Weibull parameters on the composite response was investigated within a parameter, which demonstrated the strong influence of the probabilistic parameter definition and the importance of such an approach for achieving an accurate prediction of the composite response both with and without continuous reinforcement.
The simulations with continuous reinforcement could reliably reproduce the reference experimental results in terms of macroscopic stress–strain relationships, the extent of multiple cracking, and crack width, demonstrating that the model is suitable for extensive parameter studies on material composition and constituents’ properties.
Given the discrete representation of the continuous reinforcement, the influence of various mechanical (strain rate) and environmental (temperature) effects on the composite properties could be investigated by adjusting the properties of SHCC, continuous reinforcement, and bond based on existing experimental results. Furthermore, the model can be adapted and upscaled to simulate the strengthening performance of such composites on structural elements or to predict the mechanical performance of such composites with other types of continuous reinforcement, e.g., steel rebars.