Cohesive Zone Model Based Numerical Analysis of Steel-Concrete Composite Structure Push-Out Tests

Push-out tests were widely used to determine the shear bearing capacity and shear stiffness of shear connectors in steel-concrete composite structures.Thefinite elementmethodwas one efficient alternative to push-out testing.This paper focused on a simulation analysis of the interface between concrete slabs and steel girder flanges as well as the interface of the shear connectors and the surrounding concrete. A cohesive zone model was used to simulate the tangential sliding and normal separation of the interfaces. Then, a zero-thickness cohesive element was implemented via the user-defined element subroutine UEL in the software ABAQUS, and a multiple broken line mode was used to define the constitutive relations of the cohesive zone. A three-dimensional numerical analysis model was established for push-out testing to analyze the load-displacement curves of the push-out test process, interface relative displacement, and interface stress distribution. This method was found to accurately calculate the shear capacity and shear stiffness of shear connectors.The numerical results showed that themultiple broken linesmode cohesive zonemodel could describe the nonlinear mechanical behavior of the interface between steel and concrete and that a discontinuous deformation numerical simulation could be implemented.


Introduction
The shear stiffness and shear bearing capacity of shear connectors in a steel-composite structure is usually assessed using push-out test [1][2][3][4][5][6][7].But push-out test is time consuming and expensive, and its results can be affected by interface bonding, boundary conditions, and other factors.Finite element methods can provide an efficient alternative to full-scale push-out tests.It can also be used to carry out parametrical analysis.During push-out tests, mechanical behavior of the interfaces between the concrete slab and the steel girder flange and between the shear connectors and the surrounding concrete is relatively complex.This complex interface mechanical behavior is one of the difficulties of nonlinear numerical analysis involving push-out tests.
One of the methods used in numerical analysis of pushout tests involves considering the elastic-plastic behavior of concrete and steel and neglecting the interface slip and separation.Oguejiofor and Hosain have developed a threedimensional numerical model using ANSYS software to analyze push-out specimens with perfobond rib connectors [8].The push-out test specimen was modeled therein using 3D reinforced concrete solid elements and shell elements for structural steel of beam flanges and perfobond rib connectors, by ignoring the interface mechanical behavior.Al-Darzi et al. also used this method to analyze similar perfobond rib connectors [9].Mirza and Uy developed a 3D nonlinear finite element model using ABAQUS to study the effects of combination of axial and shear loading on the behavior of headed stud steel connectors [10].Solid elements were used for the concrete slab, steel beam, and shear connectors.Concrete and steel element nodes at the interface were coupled and no interface elements were used.
Based on an experimental study, Johnson and Oehlers found that separation between the stud and the concrete on the surface of the stud shank opposite to the load can occur even at low load levels [4].To simulate this phenomenon, they assigned zero stiffness to the coincident concrete elements with stud shank surface where separation will occur.A similar method was used by Kalfas and Pavlidis [11] and by 2 Mathematical Problems in Engineering Kim et al. [12].Lam and El-Lobody developed a push-out test finite element model using ABAQUS [13].To simulate the separation between the stud root and its surrounding concrete, coincident stud nodes in the opposite direction of loading were detached from the surrounding concrete elements while nodes on the surface of the stud shank in the direction of loading were connected to the surrounding concrete nodes.Ellobody and Young have used a similar method to analyze push-out test of composite beams with profiled steel sheeting [14].
Guezouli and Lachal proposed a 2D nonlinear finite element model to study the influence of friction coefficients on push-out tests [15].This simplified 2D model showed strong convergence but did not consider local damage to the surrounding concrete of shear studs or spatial mechanical characteristics of the structure.Xu et al. developed a 3D finite element model of push-out testing with group studs [16].Contact interactions available in ABAQUS were used to simulate the interfaces between steel flanges and concrete slabs and between stud shafts and surrounding concrete [16].Okada et al. performed push-out tests on composite structures with grouped stud connectors and developed a 3D numerical analysis model which considered nonlinear properties of the material and interface bonding friction [17].The interface bonding model consisted of a linearly increasing curve rising curve and a peak platform line.A peak bonding stress of 0.9 MPa corresponding to a slip value of 0.06 mm was used based on experimental results.The interface bonding model did not consider the softening stage.
Nguyen and Kim have used the bilinear cohesive zone model available in ABAQUS to simulate the mechanical behavior of the interface between the concrete slab and the steel place of the push-out test specimen [18].Then, 0.1 cm and 0.1 cm were used for tangential stiffness and normal tensile stiffness, respectively, where  cm is the shear modulus and  cm is the elastic modulus of concrete.The critical relative displacement corresponding to peak cohesive stress and the maximum displacement at which the cohesive layer failed were determined by the authors to facilitate better agreement with experimental results.Then, tangential and normal critical relative displacements were assigned as 0.5 mm and 0.1 mm, respectively, and the displacement at failure was assigned as 0.8 mm.There were obvious differences between peak cohesive stresses adopted by the author and the actual cohesive stresses.For example, when C50 concrete was used, the values of its tangential and normal peak stress were 1.47 MPa and 3.68 MPa, respectively.In the finite element model, the shear stud nodes and its surrounding concrete nodes were tied together, and slip and separation of the interface between the shear stud and its surrounding concrete were not considered.
In summary, numerical simulation studies of pushout tests have been conducted by various researchers and documented in the literature as discussed above.However, detailed numerical analysis taking the complex mechanical behavior of the interfaces between the concrete slab and the steel girder flange and between the shear connectors and its surrounding concrete into account are not well documented.In this paper, a multiple broken lines mode cohesive zone model was used to describe the tangent slip and normal cracking at the interface of steel and concrete.Then a zerothickness cohesive element was implemented using a userdefined element subroutine UEL in ABAQUS [19].Finally, a three-dimensional numerical analysis model was presented simulating push-out testing.The load-displacement curve of the push-out test process, interface relative displacement, and interface stress distribution were analyzed.Numerical simulation of discontinuous deformation at the interface was achieved.

Mechanical Description of Discontinuous Deformation
Consider a discontinuous physical domain Ω as shown in Figure 1.The domain contains a cohesive crack, and the cohesive interfaces can be denoted by Γ +  and Γ −  .The prescribed tractions F are imposed on boundary Γ  and the prescribed displacement u on Γ  .The stress field inside the domain, , is related to the external loading F and the tractions t + and t − along the discontinuity through the equilibrium equations [20]: Here the traction t is a function of the relative displacement w between Γ +  and Γ −  , that is, t = t(w).The domain surrounding the discontinuity is assumed to be elastic.We further assume small strains and displacement condition.Thus, the constitutive law and geometric equation for the domain can be written as in which C denotes the material stiffness tensor.
The displacement u must be one of the set of kinematically admissible displacement, Using the principle of virtual work governing equations in integral form can be written as follows [21]:

Cohesive Zone Model
The interfaces between the concrete slabs and the steel girder flanges and between the shear connectors and the surrounding concrete of a typical push-out test are shown in Figure 2. Cohesive bonding stress exists at the interface of concrete and steel during the push-out test process.No slip is expected at the interface when longitudinal shear stresses are lower than the bonding resistance.As loads increase and longitudinal shear stresses exceed the bonding resistance, interface slippage occurs.If normal stress exceeds the tensile strength of the interface, crack initiation and propagation take place, which causes uplift forces on the shear connectors.To conduct continuous-discontinuous deformation analysis of a push-out test, a cohesive zone model was used here to describe the relationship between the interface shear stress and slip displacement and between the normal stress and the tensile displacement.
Dugdale proposed the cohesive zone model to describe the relationship between cohesive stress and cracking displacement during the material fracture process [22].Yang et al. developed a functional relationship between cohesive stress and the relative displacement to analyze mode-I and mode-II fracture using a criterion proposed by Wang and Suo [23][24][25][26].Cohesive zone models have been used to analyze the mechanical behavior of the bond interface between fiberreinforced polymers (FRP) and concrete [27,28].Ling et al. used a cohesive zone model based augmented finite element and analyzed progressive failure at the soil-structure interface [20,29].If parameters are properly selected, the cohesive zone model can indicate mechanical properties of the bond interface, such as modulus, strength, and toughness [30].
As shown in Figure 3, a multiple broken line mode cohesive zone model was used in this paper.In this model, Δ and Δ refer to normal displacement and slip displacement, respectively. and  represent the normal and shear stresses, respectively. 1 and  1 are peak stresses for mode-I and mode-II fractures, respectively.The multiple broken lines mode cohesive zone model can be written as follows: (5)

Cohesive Interface Element
A zero-thickness cohesive interface element was implemented using a user-defined subroutine UEL in ABAQUS [21,27,31,32].In the user-defined element, the element stiffness matrix (AMATRX), nodal residual force vector (RHS), and state variables (SVARS) must be defined.The eight-node cohesive interface element used in this paper is shown in Figure 4.
The nodal displacements of cohesive interface element in the global coordinate system are denoted by u; then, the relative displacement between the top and bottom nodes can be given as follows: Here N(, ) is the standard shape function.The matrix B(, ) is defined as follows: Then the relative displacement of the interface can be written as follows: Calculation of the transform matrix that describes the relationship between the local and global coordinates is shown below.During large deformations, initial configuration is given by x, and the reference surface state x  can be computed using a linear interpolation between the top and bottom nodes in their deformed state as follows: T 1 and T 2 indicate unit tangent vectors of the local coordinate element and T  is used to denote unit normal vector.The unit normal vector T  can be written as follows: Here ‖ ⋅ ‖ denotes the norm of a vector.Then the unit tangent vector can be given as follows: Then the transform matrix that describes the relationship between the local and global coordinates can be written as follows: Local displacements are then obtained as follows: Cohesive stresses can be calculated using the specified cohesive laws (Figure 3) and the relative displacement of the interface.Then node force vector can be obtained as follows: Here |J| is Jacobi matrix value of the transform matrix.
The tangent stiffness matrix of the cohesive interface element can be written as follows: Here D  = t loc / loc is tangent stiffness matrix of the cohesive zone model.A solution algorithm of cohesive interface element is shown in Figure 5.
In the finite element model, cohesive interface elements are utilized at the interface between concrete and steel to simulate initiation and propagation of cracks.Conventional solid elements can be used to model the concrete and the steel plate.

Geometry of Push-Out
Testing.The geometry of the pushout test specimen analyzed in this paper was the same as that used in an experimental study performed by Guezouli and Lachal [15].The geometry of the push-out test specimen is shown in Figure 6.The height and width of the steel beam were 260 mm, the thicknesses of the flange plate and web plate was 17.5 mm and 10 mm, respectively.The height, width, and thickness of the concrete slab were 620 mm, 600 mm, and 150 mm, respectively.The diameter of reinforcement in the concrete slab was 10 mm, the lengths of the transverse and longitudinal reinforcement were 520 mm and 550 mm, respectively.The height of the studs was 100 mm.The diameter of the stud shanks was 19 mm, and the diameter of the stud heads was 31.7 mm.

Material Parameters.
Constitutive relationship of the concrete used in this paper is shown in Figure 7(a).Young's modulus of the concrete slab   = 36,900 MPa and Poisson's ratio was equal to 0.2.The cylinder strength in compression   = 56 MPa and the one in tension   = 3.96 MPa were used in the model.Based on the information provided in literature, the proportional limit stress was set at 0.8   = 44.8MPa, and the corresponding strain was set at 0.0012 [14,17,18].The compressive strain associated with ultimate strength was equal to 0.0022.The ultimate strain of concrete at failure in compression and in tension was equal to 0.01 and 0.005, respectively.A damage plasticity model available in ABAQUS was utilized for the concrete element.
Young's moduli of the steel beam, shear stud, and reinforcement were all equal to 210,000 MPa.Poisson's ratio was equal to 0.3 for the steel.An ideal elastic-plastic model was used for the steel beam.The yield strength of the steel beam was equal to 355 MPa.The constitutive relationship of shear stud and reinforcement is shown in Figure 7(b).The yield stress and ultimate stress were 500 MPa and 550 MPa, respectively.Based on information available in the literature [10,17,33], strain before strain hardening and strain when ultimate stress is reached are set at 0.02 and 0.10, respectively.

Finite Element Model.
The whole geometric model of the push-out specimen is shown in Figure 8(a).Because of the symmetry, it was only necessary to model half of the actual structure using the ABAQUS program, as shown in Figure 8(b).Then 3D solid elements were used for concrete slabs, steel beams, and shear studs.Reinforcement was modeled using truss elements.The user-defined cohesive interface elements were implemented at the interfaces between the concrete slab and steel girder flange and between the shear connectors and the surrounding concrete.The finite element mesh is shown in Figure 8(c), in which the highlighted region is the position where the cohesive interface elements were implemented.
In this paper, two models with different boundary conditions were considered.In one of the simulation models, the concrete slab at the bottom was allowed to slide freely in lateral direction.Degree of freedom U2 was not constrained.This is hereafter referred to as the lateral free model.In the other simulation model, the concrete slab at the bottom was constrained in lateral direction.This is hereafter referred to as the lateral fixed model.For the actual push-out test experiment, the real boundary conditions of the concrete slab at the bottom involve contact with the base support.The load-bearing capacity of the shear stud in the experiment was found to be in between the values observed in the two simulation models.

Shear Capacity and Shear Stiffness of the Shear Connector.
The load-slip curves of the push-out test process are shown in Figure 9(a).The ordinate value is the average force per stud defined as the total action load divided by the total number of studs.The abscissa is the average value of the slip at the top of the interface (point U in Figure 6) and the slip at the bottom of the interface (point D in Figure 6).The experimental results shown in Figure 9(a) were reported by Guezouli and Lachal [15].Results of the shear strength of shear connectors calculated by Eurocode-4 and AASHTO LRFD are also shown in Figure 9(a) [39,40].As shown, slip values calculated using the two different boundary models are similar when the applied load was relatively small (<60 kN).The shear capacity of lateral fixed model and lateral free model was 156 kN and 138 kN, respectively.The applied load when the slip value reached 5 mm was adopted as the shear capacity, because the applied load did not increase evidently when the slip values exceed 5 mm.The shear capacity of lateral fixed model was the same as the experimental result, while the shear capacity of lateral free model was 13.5% smaller than experimental value.Results indicated that the boundary conditions of the concrete slab at the bottom could influence the shear strength of the shear connectors, and friction at the concrete slab base increased the bearing capacity.This was consistent with test results reported by Johnson and Oehlers [4].Results of the secant shear stiffness are shown in Figure 9(b).These results indicate that the secant shear stiffness values calculated using the lateral fixed model and the lateral free model were both similar to the experimental results.At relatively small slip values, shear stiffness decreased rapidly as slip increased.

Separation between the Stud Root and Its Surrounding
Concrete.For the middle-row studs (as shown in Figure 6, stud A (hereafter referred to as the top stud), stud B (hereafter referred to as the middle stud), stud C (hereafter referred to as the bottom stud)), and the separation between the stud root and the concrete on the surface opposite to the load is shown in Figure 10.As shown, separation took place even at low load levels.This is consistent with the experimental results reported in the literature [3,4].Separation values calculated using the lateral fixed model and lateral free model were similar when the applied loads were relatively small.The differences among separation values calculated by the two simulation models became larger as the applied load increased.When the applied load reached 138 kN, separation values at the top stud, middle stud, and bottom stud calculated using the lateral free model were 2.1, 2.3, and 2.5 times the values calculated using the lateral fixed model, respectively.In the direction of loading, separation values at the three studs were not equal.Separation values were largest at the top stud and were lowest at the middle stud.The top stud was the closest to the applied load and carried the maximum shear force, thus separation values at the top stud were largest.Separation values at the middle stud were close to its values at the bottom stud.Although the bottom stud was far away from the applied load, it was closer than the middle stud to the bottom of concrete slab, and shear force of the bottom stud was slightly larger than the middle stud because of the influence of boundary.

Normal Separation between the Concrete Slab and Steel
Beam.Normal separation at the bottom centre of the interface between the concrete slab and the steel plate (point D in Figure 6) is shown in Figure 11.Results indicated that normal separation values calculated using the lateral free model were larger than those calculated using the lateral fixed model.When the applied load reached 138 kN, the normal separation value calculated using the lateral free model (2.34 mm) was about 35 times the value calculated using the lateral fixed model (0.066 mm).When the concrete slab at the bottom was able to slide freely in the lateral direction, the uplift phenomenon was more obvious and the shear studs were subjected to tensile forces.With this phenomenon, it reduced the shear bearing capacity of the shear connectors.The axial stress of the stud when the applied load reached 138 kN is shown in Figure 12.
Results of the normal separation distribution at the middle position (connection line between point U and point D in Figure 6) and at the side position of the interface in the push out direction are shown in Figure 13.It must be noted that the separation distribution is not uniform along the transverse direction, and the simplified 2D model is difficult to simulate this mechanical behavior accurately [15].It can be seen from Figure 13 that the separation value at the middle position is larger than at the side position around the top row shear studs.Around the bottom of the interface, when the concrete slab base can slide freely in the lateral direction, the normal separation value at the middle position is larger than at the side position.The normal separation values were relatively smaller when the concrete slab base was constrained in the lateral direction.Displacements in the lateral direction (U2) are shown in Figure 14 when the applied load reached 80 kN.

Slip Distribution of the Interface between the Concrete
Slab and Steel Beam.Results of slip distribution of the interface between the concrete slab and steel beam in push out direction are shown in Figure 15.The slip value of the interface calculated using the lateral free model was larger than that of the lateral fixed model.When the applied load reached 138 kN, the slip value at the top of the interface (point U in Figure 6) calculated using the lateral free model was 1.92 times the value calculated using the lateral fixed model.
Results indicated that the interface slip was not evenly distributed.For the region between the top row studs and the bottom row studs, the slip value at the side position was larger than at the middle position.Slip values were higher in the regions above the top row studs and below the bottom row studs because of the compression deformation of the steel flange and concrete slab.For regions around the three row studs in the push-out direction, slip values of the interface around the top row studs were larger than in the other two rows.This is because the top row studs carry a higher share of the total shear force.When the applied load reached 80 kN, the slip at the top of the interface (point U in Figure 6) and the slip at the bottom of the interface (point D in Figure 6) calculated using the lateral free model were 0.58 mm and 0.38 mm, respectively.The former value is about 53% larger than the latter.The slip values calculated using the lateral fixed model were 0.41 mm (point U in Figure 6) and 0.31 mm (point U in Figure 6), respectively.The former value is about 52% larger than the latter.When the applied load reached  138 kN, the slip values at the top of the interface calculated using the lateral free model and lateral fixed model were about 9% and 8% larger than the slip values at the bottom of the interface, respectively.Results of displacement in the pushout direction (U3) are shown in Figure 16, when the applied load reached 80 kN.due to compressive stress.As calculated using the lateral fixed model, the load values corresponding to the concrete deformation values of 2.0/1000 and 3.5/1000 were 23 kN and 32 kN, respectively.While its values calculated using the lateral free model were 24 kN and 34 kN, respectively.Concrete around the stud head can undergo plastic strain in tension because of excessive principal tensile stresses.As calculated using the lateral free model, the most serious crack state appears on the concrete around the top stud head.As calculated using the lateral fixed model, the most serious crack state appears on the concrete around the bottom stud head.Failure modes calculated using the proposed models  can agree with the experimental results reported by Ollgaard et al. [5].

Conclusions
A multiple broken lines mode cohesive zone model was used in this study to describe the tangent slip and normal crack at the interfaces between concrete slab and steel girder flange and between shear connectors and surrounding concrete.A zero-thickness cohesive element was incorporated into a finite element model using the user-defined element subroutine UEL in ABAQUS.A three-dimensional numerical analysis model was established for push-out testing, and a load-displacement curve of the push-out test process, interface relative displacement, and interface stress distribution were analyzed.The following conclusions are drawn from this study.
(1) The method proposed in this paper can accurately calculate (a) the shear strength and shear stiffness of the shear connectors, (b) the normal separation and tangential slip of the interfaces between concrete slab and steel girder flange, and (c) the normal separation and tangential slop of the interfaces between shear connectors and its surrounding concrete.Separation between the stud root and the concrete on the surface opposite to the load took place even at low load levels.
For the push-out test model analyzed in this paper, separation values were largest at the top stud, lowest at the middle stud.Results indicated that normal separation and tangential slip of the interface between concrete slab and steel girder flange are not evenly distributed.
(2) The boundary conditions of the concrete slab at the bottom can influence the shear strength of the shear connectors.For the push-out test model analyzed in this paper, shear capacity when the concrete slab base was constrained (156 kN) was about 14% larger than and when the concrete slab base was free in the lateral direction (138 kN).The separation values between the stud root and the concrete on the surface opposite to the load calculated using the lateral fixed model and lateral free model were similar when the applied loads were relatively small, and the difference in the separation values calculated using the two simulation models became lagrer as the applied load increased.Constraints of the concrete slab base were found to markedly influence the normal separation of the interface between concrete slab and steel girder     flange.When the concrete slab at the bottom was able to slide freely in the lateral direction, normal separation values were much larger and the uplift phenomenon was more obvious than when it was not.Because of this phenomenon, a free concrete slab reduces the shear bearing capacity of shear connectors.The tangential slip of the interface between the concrete slab and the steel girder flange calculated using the lateral free model was larger than that of the lateral fixed model.
(3) The multiple broken lines mode cohesive zone model used in this paper was found to effectively describe the nonlinear mechanical properties of the interface between the concrete and the steel.Thus, the discontinuous deformation numerical simulation of the interface was achieved.Shear strength and shear stiffness of the shear connectors were calculated accurately.

Figure 2 :
Figure 2: Interface between steel and concrete of typical push-out test.

Figure 5 :
Figure 5: Solution algorithm of cohesive interface element.

Figure 10 :
Figure 10: Separation between the stud root and the surrounding concrete.

6. 5 .Figure 11 :
Figure 11: Normal separation at the bottom central portion of the interface.

Figure 12 :Figure 13 :
Figure 12: Axial stress of the stud when the applied load reached 138 kN.

Figure 14 :Figure 15 :
Figure 14: Displacement in lateral direction when the applied load reached 80 kN.

Figure 16 :
Figure 16: Displacement in push out direction when the applied load reached 80 kN.

Figure 17 :
Figure 17: Plastic strain magnitude when the applied load reached 138 kN.

Figure 18 :
Figure 18: Equivalent plastic strain in tension when the applied load reached 138 kN.