Modelling crack propagation in RC beam-column joints

Accurate modelling is required to estimate crack propagation in a beam–column joint. In this study, a numerical method is developed to model crack propagation and failure loading in a beam–column joint under static load. To realize this objective, a four-node, thin-layer interface element is produced to model the fracture process zone and crack propagation. Moreover, the fracture criterion for determining the growth of a crack based on the release rate of strain energy is established. To validate the present model, ABAQUS software is used to simulate crack propagation by conventional cohesive elements. The numerical results obtained are extremely close to the experimental results within an accuracy level ranging from 4.3% to 6.7%. Meanwhile, the ABAQUS software data and the experimental data are predicted at a margin of error ranging from 12.4% to 16%.


INTRODUCTION
The beam-column joint is a structural element that is most vulnerable to sudden failure under static load, and thus, crack propagation in beam-column joints should be predicted [1].A crack in a joint starts at the point of intersection of the beam and the column on the tension face.Thus, fracture mechanics theory is the best method for predicting such cracks [2].A beam-column joint, which provides continuity to the structure, is a critical and significant part of concrete structures.Therefore, studying the behaviour, crack pattern, failure mode, and strengthening of beam-column joints is essential.
Presenting an accurate fracture mechanics model of a beam-column joint is significant because of the complex behaviour of joints.Taylor [3] was the first to study the behaviour of joints.The code developed by the American Concrete Institute committee [4] and the Uniform Building Code [5] were the first codes to recommend details for joint design.
Joint behaviour is complex because it involves different mechanisms such as gravity load, flexure, and shear.A crack that formed under static load may occur in the joint (that is, shear failure), beam, or column.A shear crack in the joint is generally more significant than other types of cracks because shear failure in joints can cause the entire structure to break down.A main crack begins at the conjunction point of the column and the beam because stress is concentrated at this point.When a shear crack occurs in the joint (that is, shear failure), the column and the beam rotate and lose their load capacity.
Many investigations [6][7][8][9][10][11][12][13][14][15] have been conducted to investigate the behaviour, performance, and design of several types of joints, such as exterior and interior joints, under static and dynamic loads [16].The use of fracture mechanics to investigate crack propagation in joints under static load has not yet been reported.Accordingly, the present study aims to use fracture mechanics to study cracks in joints.
To model a crack using fracture mechanics, the cohesive zone model (CZM) was first proposed by Dugdale [17] to analyse brittle fracture.Hillerborg et al. [18] were the first to use CZM, or the fictitious crack model, to calculate the softening fracture of concrete beams [19,20].
A fracture zone with a large and variable size was introduced by Hillerborg et al. [18] in front of the main crack.This fracture zone is referred to as the fracture process zone (FPZ).The FPZ can transfer normal and shear stresses close to the crack [21].In this zone, both stresses are reduced as crack opening increases.Several microcrack systems, such as microcracking in the matrix, fracture cement-matrix border, aggregate bridges, and crack branching, also exist in the FPZ.The length of the FPZ is 500 mm in normal concrete, 3000 mm in mass concrete (e.g., dam concrete), and 1000 mm in mass soil [22].The FPZ frequently occurs as a localized zone with a certain width.
Accurately predicting crack propagation in a concrete structure, such as a beam-column joint, is essential to improve reliability, durability, and serviceability of a structure [23].Therefore, accurate FPZ modelling is required.The width of the FPZ should be large and nonnegligible [24,25].Knowledge on the width of the FPZ in a cohesive model for the finite element method remains minimal [26].
Early researchers generally used conventional elements, such as a four-point isoparametric element for concrete, to study crack direction, as well as an interface element with two degrees of freedom per node to simulate crack propagation [27].A crack can deflect, and thus, the aforementioned elements do not consider the effect of drilling degrees of freedom on interface elements.In fact, in all previous numerical models, the effect of drilling degrees of freedom on interface elements was ignored.Inplane rotational degrees of freedom are also referred to as "drilling degrees of freedom".Two general types of problem are involved in plane analysis: plane stress and plane strain problems.These problems are defined by setting certain restrictions and assumptions on the strain in the stress and displacement fields.The restrictions can be overcome by including elements with drilling degrees of freedom, which can improve the limitations of the elements with two degrees of freedom.Using elements with drilling degrees of freedom to analyse cracks has not yet been reported in studies involving fracture mechanics.
The current study shows how an interface element can be developed to model Mode I fracture in a joint.A fournode, thin-layer interface element is developed to model the FPZ.Choosing such an element enables us to consider the width of the FPZ in modelling.To estimate crack propagation in a reinforced concrete (RC) joint, a crack is assumed to grow if the maximum tensile principal stress at the node in front of it reaches the maximum strength in the joint.

METHOD AND MATERIALS 2.1 Thin Layer Interface Element
A cohesive model assumes that the FPZ occurs along a line in 2D problems.Interface elements have been used in finite element methods to model the FPZ.The present study employs a four-node, thin-layer interface element to model the FPZ.This element with thickness t (Fig. 1) considers the width of the FPZ.By denoting the shape functions at the nodes as N i , the relative displacement vector u, v and θ can be expressed as follows: where q T is the displacement vector and N is shape function matrix.Shape functions based on ξ and η coordinate are as follows: ( )( ) where ξ i and η i are the coordinates of the i th node within the element in the x and y directions, respectively.If the relative displacements vector u, v and θ and the thickness of the interface element t are small, then the strain at any point can be defined as follows: where B is the strain-displacement transformation matrix.ε t and ε n are the shear and normal strains, respectively.The stiffness of the thin-layer interface element K in the softening zone is given by: where dA is the differential element of the crack surface area; T is the transpose; and D is a matrix that describes the relationship between stress and displacement, which is expressed as follows: Here, D n , D t and ν are the normal elastic stiffness, shear elastic stiffness, and Poison's ratio, respectively [28].The vector F of the cohesive forces in the nodes is given by: where σ ij = {τ σ σ} T is the stress vector caused by the constitutive model of the normal and shear stresses in the FPZ.By using Gaussian integration, the stiffness of the interface element and the nodal cohesive forces can be obtained.

Crack Propagation in the Joint
To implement CZM in its simplest form, two parameters are required: the energy release rate and cohesive strength.The former can be used to describe the crack propagation criterion in the fracture process in the crack tip.This approach states that the energy release rate, which is defined as the amount of energy stored in the FPZ that is required to form the crack must be sufficiently larger than the critical fracture energy.A suitable crack propagation criterion is required to evaluate crack propagation [29].The strain energy release rate, which is based on a small length variation, is as follows: Therefore, the Griffith criterion for crack propagation can be expressed as follows: where Δα and G F are the lengths of the interface element and the critical fracture energy of the concrete material, respectively.The stiffness matrix of the interface element in the crack tip can contribute to Eq. ( 8), which can be estimated using a finite difference procedure.In the present study, Eq. ( 8) is computed by assuming that the crack extension is equal to the length of the interface element.
To describe cohesive strength, assuming that the crack will grow if the maximum tensile principal stress at the node in front of it reaches the maximum strength in the joint is rational.
The relationship between the maximum strength in the joint and concrete compressive strength was proposed by Goto and Joh [30] based on experimental study.They found that maximum strength (that is, joint shear stress) is 5 c f ' from the scatter data, where f' c is the concrete compressive strength.Thus, the crack is assumed to propagate in the joint when the maximum tensile principal stress at the node in the crack tip reaches 5 c .f ' Moreover, they concluded that when the joint shear stress at the yielding beam bars reaches 4 c f ' , joint shear will not occur.In addition, when the maximum tensile principal stress at the node in the crack tip is between 5 c f ' and 5 c f ' , the crack will propagate in the beam.
To study other cracks in the beam or the column, the crack is assumed to grow if the maximum tensile principal stress attains the strength of the concrete as a result of the weak tension in the concrete materials.In the current investigation, the direction of a crack is modelled by a method where a crack follows existing inter-element borders [31].
This method has a simple algorithm and does not require remeshing.Crack propagation follows one of the inter-elements, that is, (AB) or (AC), in which the crack is assumed to stop and intersect with the main element (Fig. 2).Two cases are possible for the crack path.If the orientation angle θ is less than 45°, then the path of growth is AB; otherwise, it will be AC.Although crack paths are not smooth, those detected using this method exhibit good agreement with the correct crack path.The FEAPpv ® program code is developed to analyse 2D plane stress in concrete [32].A nonlinear element is adopted for the interface element in the User Subroutine FEAPpv® Fortran programming, whereas a nonlinear dynamic relaxation method is used for the interface element in the program [33].Four-node isoparametric elements are used for bulk concrete with linear elastic behaviour.
The major steps in the nonlinear crack propagation modelling process in a joint, which is based on the proposed method, are as follows.1) Conduct an analysis and identify the region to obtain the maximum tensile principal stress.
2) The nonlinear dynamic relaxation steps are repeated until the equilibrium equation is satisfied within an allowable tolerance.
3) When the maximum tensile principal stress at the node in the crack tip exceeds 5 c f ' , the crack will propagate in the joint in the identified direction.4) Once the FPZ propagates, the element stiffness is set to zero and the crack grows along the respective element by considering the direction in each step.5) Repeat Steps 2 to 4 until the structure fails.
A 2D plane stress finite element method is employed to predict crack propagation in a joint.An elastic-perfect plastic is considered to model the behavior of the steel bars.The steel bars are modelled by two-node beam-column elements with three degrees of freedom per node.Fournode elements with three degrees of freedom per node [34] are used for bulk concrete.The strategy is that cracks start from the surface and then are formed in the structure.

RESULTS AND DISCUSSION
A beam-column joint under static load is considered to study the behaviour of a joint based on the proposed model.The same scenario was analysed by Parvin and Granata [35].The geometry and boundary conditions of the joint are shown in Fig. 3.The compressive strength, tensile strength, and elastic modulus of the concrete are set as 27.6 MPa, 2.8 MPa, and 25.8 GPa, respectively.The elastic modulus and yield stress of the steel are 200 GPa and 413.7 MPa, respectively.The axial load on the column (P 1 ) is 26.7 kN.The cross sections of the column and the beam are illustrated in Figs.4a and 4b, respectively.The joint is reinforced as a strong column/weak beam, as designed by Parvin and Granata [35].To validate the present model, ABAQUS software is used in this study to simulate crack propagation through conventional cohesive elements (COH2D4P).The initial mesh is illustrated in Fig. 5.The mesh of the joint in the present model has 1036 elements.The average element size is 16 × 16 mm, with finer mesh (14 × 12 mm) and 305 interface elements.The joint is modelled using ABAQUS software with 9874 C3DBR S4R (average size: 15 × 15 × 12 mm) elements.Fig. 7 shows the plot of the load versus the relative rotation on the column face for different concrete compressive strengths.As shown in Fig. 7, the concrete compressive strength of the beam and the column considerably affects the load-rotation curve of the joint.
The load increases with linear behaviour during the initial phase as the crack grows in the joint, which is understandable because no microcrack occurs in the crack tip.The load is nearly constant and exhibits plastic behaviour because of the significant effects of steel under a certain load.That is, only the FPZ fully increases at the beginning, and then, the effective crack grows in the concrete without steel yielding.When the steel yields, rotation increases as the crack in concrete grows.Moreover, Fig. 7 compares the load-displacement in the beam in the present model with that in the model obtained using ABAQUS software (concrete compressive strength = 35 MPa).As shown in Fig. 7, the curve of the present model is slightly below that of the model obtained using ABAQUS software.Fig. 8 illustrates the crack paths in the beam-column joint using the present model with 10 scales at a failure load (P 2 = 9 kN).Fig. 9 shows the crack paths for the same joint using ABAQUS software.In both cases, the predicted crack paths occur in the beam, as expected, because the joint has been designed to have a strong column/weak beam.The cracks initially grow in a straight line, and then gradually follow an inclined slope, in the center of the joint.The initiation and location of some cracks may change slightly because of the size of the meshes.The number of cracks is ten, with two flexural-shear cracks predicated in the model.The main crack occurs in the beam.The two flexural-shear cracks, which appear near the column face, grow in the higher half of the beam depth.
Figs. 10 and 11 show the crack pattern in the joint obtained using the present model and ABAQUS software, respectively, with a concrete compressive strength of 27.6 MPa.Both models show that the shear crack occurs in the joint.The length and the number of cracks increase with the decrease in concrete compressive strength in both models.Thus, the decrease in concrete compressive strength moves the crack to the joint.Another example is the beam-column joint tested by Goto and Joh [30].The geometry and boundary conditions of this joint are shown in Fig. 12.The compressive strength, tensile strength, and elastic modulus of the concrete are set to 32.3 MPa, 3.2 MPa, and 30.0 GPa, respectively.The elastic modulus of the steel is 210 GPa.The joint, which is not reinforced, is modelled using ABAQUS software with 1248 C3DBR S4R (average size: 14 × 14 × 12 mm) elements.
Fig. 13 compares the results of the present numerical model, the ABAQUS software data, and the experimental data [30] of the principal compressive stress versus the diagonal tensile strain curve on the joint face.The present numerical results are extremely close to the experimental results within the accuracy level of approximately 4.3% to 6.7%, whereas the ABAQUS software data and the experimental data are predicted within a 12.4% to 16% margin of error.As shown in Fig. 13, the principal compressive stress initially increases with the increase in diagonal tensile strain, and then decreases.The joint suddenly breaks at the shear crack in the concrete because of the absence of reinforcement.

CONCLUSION
In the present study, RC joints were analysed by applying the finite element method in fracture mechanics, with particular attention on modelling the FPZ and crack propagation.A thin-layer interface element was used to consider the width of the FPZ.Moreover, a criterion for crack propagation in joints and a new approach to identify the crack path were presented.The present numerical model and the ABAQUS software data were then compared.The present numerical results were close to the previous experimental test results (error: 4.3% to 6.7%).Meanwhile, the ABAQUS software results were predicted to have an accuracy ranging from 12.4% to 16% compared with the experimental data.Therefore, a conclusion can be drawn that the present model is acceptable by considering the width of the FPZ.The model has been verified computationally and has been shown to be capable of predicting the crack pattern in a joint.The compressive strength of the concrete used in the beam and in the column considerably affects the load-rotation curve and the shear crack in the joint.In addition, the shear stress of the unreinforced joint is inversely and linearly proportional to the beam/column depth ratio.

Figure 1
Figure 1 Four-node thin layer: a -global coordinate, b -local coordinate

Figure 2
Figure 2 Two possible cases for direction of propagation

Figure 3
Figure 3 The geometry and the boundary conditions of joint analysed by Parvin and Granata [35] (unit=mm)

Figure 4
Figure 4 The cross section (unit=mm): a -the column, b -the beam

Figure 5
Figure 5The initial mesh of the joint

Fig. 6
Fig. 6 shows the load-relative rotation on the column face in the present model with different meshes.Meshes (1), (2), and (3) have 305, 75, and 55 interface elements, respectively.The approximate matching of the three curves demonstrates model independence from mesh size; the model also exhibits rapid convergence.

Figure 6
Figure 6 Load-relative rotation in face of column with different meshes Tab. 1 shows the results for the beam-column joint in the current research and the previous model by Parvin and Granata [35].As shown in the table, the present results are close to the previous results.

Figure 7 Figure 8 Figure 9
Figure 7 Load-relative rotation in face of column with different concrete compressive strength

Figure 10 Figure 11
Figure 10 Shear crack and the crack pattern with 25 MPa for the concrete compressive strength by the present model

Figure 12 Figure 13 Figure 15
Figure 12The geometry and the boundary conditions of the joint (unit=mm)

Table 1
Results in the beam-column joint