Numerical Research on the Pullout Failure of GFRP

The fiber pullout is a main failure after the fiber is broken in the tension of glass fiber reinforced polymer (GFRP) bolt. In this paper, the numerical analysis is done on the distribution of both fiber normal force and interface shear stress. The results show that, on the ideal interface, the fiber pullout occurs from the lower end to the upper end of the matrix gradually, and both the normal stress of the fiber and the shear stress of the ideal interface gradually increase from the lower end to the upper end. With the increase of the interface layer thickness, the shear stress concentration area on the interface is enlarged while the stress applied is reduced, and the displacement of GFRP deformation is increasing sharply. This means that the capacity of GFRP deformation is enhanced. As a soft elastic body, the interface layer with a smaller elastic modulus can make the fiber stress and the interface shear stress sharply small and well dispersed. In addition, the load can be effectively transferred to the reinforced phase fibers in a bigger interfacial layer elastic modulus with a certain strength.


Introduction
The glass fiber reinforced polymer (GFRP) is widely used as structural composite material in the varied fields of engineering construction which is composed of the adhesive solution of resins (unsaturated polyester resin, epoxy resin, phenolic resin, etc.) as a matrix phase and the glass fiber as an enhancement phase.Compared with the traditional bolt, GFRP bolt has such characteristics as excellent corrosion resistance, high strength, low cost, and easy installation.
This paper focuses on the pullout failure at the fiber broken point which is one of main failure patterns of GFRP bolt [1].However, as the GFRP bolt is mainly used in underground tunnel and construction support, the mechanical research on GFRP bolt support is mainly concentrated in the field of macromechanics of bolt distribution and installation by the early literatures, and the researches on the microscopic damage of GFRP bolt are not enough.
Recently, there are two kinds of important theoretical research development in the microscopic mechanical field of GFRP composite.One is application and development of mesomechanics [2][3][4], and the other one is the proposed property study of interface phase structure [5][6][7].The fiber-matrix interface phase behavior plays an important role in determining the mechanical properties of GFRP composites such as elastic modulus, tensile strength, and fracture toughness, and it has received considerable attention and extensive investigations.Various interfacial interaction models based on classical shear-lag theory have been adopted in analysis of the stress transfer between enhancement phase and matrix phase, especially, it was a mature way to combine shear-lag model and Monte-Carlo simulation method [8].First of all, the van der Waals interfacial interaction between the carbon nanotubes (CNTs) and polymer matrix and the interfacial chemical bonding have been modelled with the Lennard-Jones potential [9], and the many-body bond-order potential [10,11], respectively, in analysis of the shear strength of CNTs-polymer interfaces.Modified classical shear-lag theory with introducing nonlinear interface laws for CNTspolymer interfaces has been established based on the van der Waals interfacial interaction [12,13] and chemical bonding [14], respectively.And cohesive zone finite element models with a nonlinear interface cohesive law have been established to model the pullout response of CNTs-coated fibers [15,16].However, linear interface law was also proposed for   CNTs pullout from a polymer matrix to investigate fracture toughness enhancement and bond breaking [17].Moreover, shear-lag behaviors widely exist in the following composite fields such as cementitious composites, ceramic matrix composites, metal matrix composites, and nanocomposites.A type of steel-concrete composite structure with double-box (the DBSCCS model) was proposed to illustrate and predict shear-lag characteristics by obtaining its strain functions and shear-lag coefficients.It provides more reference value for engineering design and structure optimization in some extent for the steel-concrete composite structures [18].Damage and failure process of the unidirectional fiber reinforced ceramic matrix composites (FRCMCs) under tensile stress was studied with consideration of evolving compliant interface (ECI) which can respond to the failure of debonding and sliding of interface efficiently [19].The strengthening behavior of SiC reinforced Al matrix composites was well explained by the modified shear-lag model [20].Short-fiber reinforced metal matrix composites were modelled by modified shear-lag theory based fatigue crack growth life prediction.It shows good agreement with the experimental data for the low cycle and high cycle fatigue applications [21].An analytical approach to predict initiation and growth of interfacial delamination in the through-silicon via structure was developed by combining a cohesive zone model with a shear-lag model [22].Interfacial shear stress transfer of a monolayer graphene on top of a polymer substrate subjected to uniaxial tension was investigated by a cohesive zone model integrated with a shearlag model.It provides valuable insight and design guidelines for enhancing interfacial shear stress transfer in nanocomposites, stretchable electronics, and other applications based on graphene and other 2D nanomaterials [23].
The abundant applications of modified classical shearlag theory with varying degrees of complexity in analyzing interfacial properties of composites above inspire us to establish an analytical model to discuss stress distribution and the pullout evolution for GFRP bolt.Of course, they can be determined accurately by the single fiber pullout experiment with less experimental material, and it is not affected by the strain performance of the matrix [24][25][26].However, there are a large number of experimental samples needed in the experiments which are complexly prepared in a relatively small size (microns).The critical length of the fiber embedded in the matrix can only be obtained by constant experiments (the length of the fiber embedded in each sample is changed based on the result of the previous experiment).So, these samples cannot be prepared at the same time in advance.In addition, no common interface model or standard specification is presently formulated for single fiber pullout experiment; there are too many complexes influencing factors in the experiment, and the experiment is not repeated and multiple [27,28].We emphatically derive the governing constitutive equations for a single fiber pullout model with a cohesive fiber-matrix interface by introducing an interfacial bilinear cohesive law into the classical shearlag model in Section 2.6.Then, we discuss shear stress in the interface, axial stress in the fiber, and the evolution of pullout with finite element calculations.

Establishment of Finite Element Model
2.1.Geometric Model.In the body of the GFRP bolt for the numerical research in this paper, the long and evenly distributed glass fiber and the bolt were isotropic and perpendicular to the cross-section, and the matrix was made of unsaturated polyester resin.The abstract model built for the calculations could be described as one fiber at a diameter of 8 m that was embedded in a cylindrical matrix at a diameter of 18 m, and each fiber had the same length of 40 m.The lower end face of each fiber was fixed, and its upper one was pulled with a gradual increase of tensile.As the model was axisymmetric, only a quarter of the whole geometric model was used in order to reduce amount of calculation, as shown in Figure 1.

Material Properties.
In the numerical analysis on GFRP bolt in this paper, the fiber and the matrix were homogeneous, that is, the impact of inhomogeneity on the fiber pullout was ignored.The elastic modulus and Poisson's ratio of the fiber and the matrix especially presented a negative change with temperature.Table 1 shows the physical and mechanical property parameters of the fiber and the matrix.
In the initial state, it was assumed that the interface between the fiber and the matrix was fine without debonding or sliding, and there were no such defects as bubble, impurity, and residual stress in the fiber and the matrix.

Definition of Boundary Conditions and
Load.In the GFRP bolt numerical analysis in this paper, the upper direction parallel to the fiber axis was the positive direction of -axis, and the right direction perpendicular to the fiber axis was the positive direction of -axis; then the direction of -axis was determined by Descartes coordinate system.The lower end face of the XOY plane was fixed.Only unidirectional stretching was loaded in the positive direction of -axis at a loading rate of 0.04 mm/s, and no compression, shear failure, or rotation occurred.

Unit Selection and Meshing
. 163 explicit structural shell units were selected as the calculating nits, and they supported all nonlinear characteristics in the explicit dynamics numerical analysis [29,30].In addition, the fastest Belytschko-Tsay algorithm was applied [31].The regular hexahedral element was easily used to divide the regular geometry model [31].Therefore, the mapped mesh method was used which could reduce the total number of elements to improve the computational efficiency.The model after meshing is shown in Figure 2.

Material Model.
In the constitutive relation of elastic plastic material, the bilinear isotropic (BISO) was represented to build a stress-strain curve [32].In the BISO model, the stress-strain relationship curve was linear according to the bigger specific elastic modulus before the material yielding.However, it was linear according to the smaller specific elastic modulus after the material yielding.Moreover, the yield strength of the material increased in one direction, and so did the yield strength in the other direction accordingly.

Mechanical Model and Constitutive Equation.
The mechanical model numerical built in this paper is shown in Figure 3.
The average shear strength of the interface   is where  denotes load at a complete pullout of the fiber,  stands for the length of the fiber embedded in the matrix, and  is the diameter of the fiber.Assuming the length of the fiber embedded in the matrix is , and the stress of the fiber before pullout is , then, the following equation is satisfied [33]: where  is the radius of the matrix;  is the radius of the fiber;   is elastic modulus of matrix; and   is elastic modulus of the fiber.Formula ( 2) can be integrated as follows: When  = 0, When  = , The equilibrium equation of shear-lag model is Formula ( 3) is substituted into Formula (6);   is got as follows: By Formula (7), the shear stresses of the upper and lower ends of the glass fiber embedded in the matrix are  1 and  2 , respectively: According to Formulas ( 8) and ( 9),  1 is smaller than  2 in the ideal interface.2.7.Failure Criterion.The numerical simulation for the pullout of single fiber was carried out in two steps.The first step was to judge if the fibers were broken or not.If not, the pullout of single fiber likely occurred; if yes, the length of fiber embedded in the matrix was reduced and the process above was repeated until the fiber was not broken.The second step was to repeat the numerical simulation according the new adjusted length of the fiber embedded in the matrix.
For the first step the max tensile stress criterion was applied; that is, the axial tensile stress at any point on the fiber was greater than or equal to the yield strength of the fiber, as expressed by the following Formula (10); the fiber was determined to be broken.
where   is the axial tensile stress at any point on the fiber and [   ] is the yield strength of the fiber.For the second step the max shear stress criterion was applied; that is, the shear stress at any point on the interface was greater than or equal to the yield shear strength of the interface, as expressed by the following Formula (11); the fiber pullout occurred.
where   is the shear stress at any point on the interface and [   ] is the yield shear strength of the interface.Based on the above said, the numerical research on the pullout failure of GFRP bolt was done finally in the ANSYS14.0LS-DYNA secondary development platform, using the numerical model built in the APDL-An ANSYS parametric design language, the extension program written by Visual C ++ language, and the calculating units chosen as 163 explicit structural shell units with the fastest Belytschko-Tsay algorithm and element life & death technology.

Result and Analysis
3.1.The Fiber Pullout on the Ideal Interface.Figure 4 shows the development of the single fiber pullout on the ideal interface.As can be seen from the figure, the fiber pullout occurs from the lower end of the matrix, gradually extending to the upper end, and the complete pullout of a single fiber can be realized finally.This is because there is a maximum shear stress at the lower end of the fiber embedded in the matrix.This agrees with the previous theoretical analyses of stress concentration at the end of the fiber and also presents the characteristics of unidirectional fiber reinforced composite [34][35][36].
The axial stress distribution of the fiber and the shear stress distribution of the ideal interface are shown in Figures 5 and 6, respectively.The axial stress of the fiber gradually increases from the lower end to the upper end, and the axial stress at the upper end is an order of magnitude compared with that of the lower end.It is consistent with the experimental result that the fiber is always broken at the upper end in the single fiber pullout experiments.However, the maximum shear stress of the ideal interface appears at the lower end of the matrix.Once the fiber pullout comes into being, the shear stress of the ideal interface gradually becomes bigger and bigger from the broken point to the upper end of the matrix because of the energy release and shear stress reduction at every broken point.
Figure 7 shows the stress-strain curves of ideal interface during the pullout of single fiber.It can be divided into the four sections of linear-elastic deformation, partial pullout, complete pullout, and interface friction.This agrees with the conclusions drawn in literature [37,38].However, there is an obvious stress oscillation in the interface friction section in our study in a full consideration of the friction between the fiber and the matrix and the use of material model BISO.
The tension force-time curve of the fiber and the shear stress-time curve of the ideal interface are, respectively, shown in Figures 8 and 9.While  = 0.6 s, the tension force of the fiber and the shear stress of the ideal interface are 0.97 mN and 3.8 MPa, respectively.According to Formula (1) and considering calculation of a quarter of the whole model, the tension of the fiber force () is available as a theoretical value, that is,  = 0.954 mN, which is close to the numerical simulation value of 0.97 mN.Therefore, the accuracy of the numerical simulation is proved.

Effect of the Interface Layer Thickness on Single Fiber
Pullout.The elastic modulus of the interface is set as 50 GPa, while the interface layer thickness values are set as 1/8, 2/8, and 3/8 of the fiber radius, and then the corresponding values are, respectively, 0.5 m, 1 m, and 1.5 m.It should be noted that the interface is a separate layer and should be divided into fine grids to ensure the accuracy of the numerical simulation.
Figure 10 shows the development of single fiber pullout with the interface thickness of 3/8 of fiber radius.It is apparent that the fiber pullout starts from the lower end and then gradually extends to the upper end, and the complete pullout of a single fiber is available finally.The single fiber pullout with the other two interface thickness values develops in the same process.However, a bigger gap of pullout can be obviously observed with thickening of the bigger interface layer.
Figures 11 and 12 show, respectively, both the axial stress distribution of the fiber and the shear stress distribution of the interface at a moment of  = 0.38 s at different interface layer thicknesses of 1/8, 2/8, and 3/8 of the fiber radius.Without obvious effect from the interface layer thickness, no matter what the interface layer thickness value is, the axial normal stress concentration of the fiber has the nearly same distribution with the bigger concentration at the upper end, as shown in Figure 11.However, with thickening of the interface layer, the shear stress concentration area of the interface is enlarged but the degree of shear stress concentration is reduced.This means that the interface layer thickness can effectively reduce the shear stress concentration, as shown in Figure 12.
Figures 13 and 14 show, respectively, the axial stressstrain curve of the fiber and the shear stress-strain curve of the interface with different interface layer thicknesses of 1/8, 2/8, and 3/8 of the fiber radius.For each stressstrain curve and shear stress-strain curve based on thickness values of interface layers, there is nearly the same varying curve which includes three stages: linear-elastic deformation, plastic deformation, and pullout.Because the elastic modulus (50 GPa) of interface is unchanged, the GFRP bolt has nearly the same linear-elastic deformation and is not affected by the interface layer thickness.However, both the axial stress of the fiber and the shear stress of the interface increase a little bit and the displacement in the plastic deformation increases sharply, and the plastic deformation is enhanced as the interface layer thickness increases.There is obvious concussion observed in the pullout stage because of the friction of interface.

Effect of the Interface Layer Elastic Modulus on Single Fiber
Pullout.The thickness of the interface layer is set as 1/8 of the fiber radius, that is, 0.5 m, while the interface layer elastic moduli are set as 70 GPa, 50 GPa, and 2.5 GPa, respectively.interface at the moment () of 0.38 s with different interface layer elastic moduli of 70 GPa, 50 GPa, and 2.5 GPa.No matter what the interface layer elastic modulus is, the axial normal stress of the fiber becoming bigger from the lower end to the upper end has the nearly same distribution with a concentration at the upper end, as shown in Figure 15.For the shear stress of the interface at the broken point of the lower end, there is a concentration on the interface layer with the three different elastic moduli of 70 GPa, 50 GPa, and 2.5 GPa.However, at the upper end, the shear stress concentration cannot be observed on the interface layer with an elastic modulus of 2.5 GPa but on the interface layer with bigger elastic moduli of 70 GPa and 50 GPa, as shown in Figure 16.This is because the shear stress of the interface can be well dispersed by the bigger plastic deformation on the interface layer with a smaller elastic modulus of 2.5 GPa.It corresponds with the viewpoints from literatures [39].Figures 17 and 18 show, respectively, both the axial stressstrain curve of the fiber and the shear stress-strain curve of the interface with different interface layer elastic moduli of 70 GPa, 50 GPa, and 2.5 GPa.They have three stages of linear-elastic deformation, plastic deformation, and pullout and present nearly the same distribution with relative bigger interface layer elastic moduli of 70 GPa and 50 GPa.In contrast, the stress of the fiber and shear stress of the interface become less sharp and are obtained, respectively, with the smaller interface layer elastic modulus of 2.5 GPa.To carry the stress, the interface with the smaller elastic modulus, as a soft elastic body, can transfer the stress evenly on the entire interface to avoid the stress concentration.This is conducive to improving the macroductility of GFRP.Of course, the interface should have a certain elastic modulus to ensure that the load can be effectively transferred to the reinforce phase fibers at a certain strength.This is conducive to improving the macrostrength of GFRP [40].Therefore, there is an optimal interfacial layer elastic modulus worthy of exploring.

Conclusions
The numerical mechanical research has been carried out for the pullout failure of the GFRP bolt in different interface conditions.There are three deformation stages included in the pullout process: linear-elastic deformation, plastic deformation, and pullout.The interface structure plays a key role for the mechanical behavior in the pullout of GFRP bolt.Firstly, there is a stress concentration at the fiber broken point at the lower and the upper ends of the matrix and the concentration at the lower end is distributed more widely and rigidly.Secondly, the interface with bigger layer thickness but less elastic modulus can make the interface shear stress distributed more evenly with less concentration.Moreover, the plastic deformation capacity of the interface increases, which means that GFRP bolt should present properties of high strength and ductility with the nonoccurrence of brittleness failure.It should be noted that it is not a fully predictive model because compressive stresses on the fiber due to matrix shrinkage, stress transfer in the interface, matrix cracking, and fiber breaking are not considered.Consequently, the limitations of the simplified models are unavoidable, and it can only essentially add to the knowledge base of numerical investigation for single fiber pullout.

Figure 1 :
Figure 1: A quarter of the geometric model.

Figure 2 :
Figure 2: Mesh generation of the model.

Figure 3 :
Figure 3: Analysis of stress concentration in the model of single fiber pullout.

Figure 5 :
Figure 5: The axial stress distribution of the fiber with the ideal interface.

Figure 6 :Figure 7 :
Figure 6: The shear stress distribution of the ideal interface.

Figure 8 :Figure 9 :
Figure 8: The tension force-time curve of the fiber under the condition of ideal interface.

Figure 10 :Figure 11 :
Figure 10: Development of the single fiber pullout with the interface layer thickness of 3/8 of the fiber radius at (a) 10%, (b) 30%, (c) 50%, and (d) 100% of the whole interface layer.

8 Figure 14 :
Figure 14: The shear stress-strain curve of the interface in the single fiber pullout with the different thickness of interface layer.

Figure 15 :Figure 16 :
Figure 15: The axial stress distribution of the fiber at  = 0.38 s with different interface layer elastic modulus of (a) 70 GPa, (b) 50 GPa, and (c) 2.5 GPa.

Figure 17 :Figure 18 :
Figure 17: The axial stress-strain curve of the fiber with different elastic moduli of interface layer.

Table
: Physical and mechanical property parameters of the fiber and the matrix.