Finite element simulation of a motorway bridge collapse using the concrete damage plasticity model

. The aim of this work is to show how the concrete damage plasticity model developed by Lubliner et al. can be applied for calculation of a motorway bridge collapse occurred in the Amur region, Russia. The concrete structural behaviour is highly complex. Being a quasi-brittle material, concrete demonstrates softening behaviour that is numerically complex due to the loss of positive definiteness of the tangent rigidity matrix of the material, and hence the loss of the ellipticity of the equilibrium rate equation. This eventually leads to the loss of well-posedness of the rate boundary value problem. Besides that, concrete behaviour in compression differs from that in tension. There are a few different failure modes of concrete material: tension cracking, compression crushing, spalling of concrete, etc.


Introduction
The problem of finite element simulation of concrete structures is related to complexities that one faces when describing non-linear material behavior that leads to the development of cracks. In order to include those non-linear effects one needs to specify non-linear parameters of material, which is not a straightforward process since data dispersion is present when it comes to experimental testing of concrete. In this paper, finite element simulation is performed regarding one of bridge supports after a long period of operation. The goal of this work is modeling a part of the bridge support using concrete material model that includes non-linear behavior of the material.
As a sample for calculations we took the motorway bridge crossing railway tracks in the town of Svobodnyi, Amur region, Russia. Two spans of this bridge collapsed on October 9, 2018. At the moment of collapse a lorry truck was on the bridge and a freight train was passing under. Figure 1 contains photo showing the collapsed bridge.
The bridge superstructure consists of several spans with the following lengths: 1×18.0 m + 3×33.0 m + 1×18.0 m. In the present paper, one bridge support is considered which consists of an L-shaped beam lying on six columns (Fig. 2). The support is made of reinforced concrete. Rebars are introduced only in the L-shaped beam since the crack appeared in the beam and was followed by bridge collapse.

Concrete Damage Plasticity (CDP) model
The Concrete Damage Plasticity (CDP) model is selected for this study as this model has the potential to represent complete inelastic behavior of concrete both in compression and in tension, including damage evolution in material. This model assumes two main failure mechanisms in concrete -compressive crushing and tensile cracking. The complete behavior of material is captured by uniaxial compression (Fig. 3) and tension curves (Fig. 4).  Since concrete is exposed to softening under compression and tension, the user should define modulus of elasticity as well as plasticity parameters and properties of compression/tension behavior. Compression and tension behavior is determined by inelastic strains and by corresponding damage variables: , , , . Inelastic strains are calculated by formulae: where / is total strain, , 0 are elastic strains in compression and tension.
Knowing the values of inelastic strains, one can calculate plastic strains using damage variables: where 0 is non-damaged modulus of elasticity.
Damage variables are assumed to take the values of 0 ≤ , ≤ 1, where 0 corresponds to non-damaged material and 1 corresponds to complete failure of material at a point of calculation (integration point). There are a lot of different damage evolution equations [1,2,3]. In this work exponential evolution is assumed: where b and g are material constants: b defines how fast the damage (d) approaches 1 and g governs the shape of the curve; 0 represents threshold value.
Simulation is based on real material data tested in laboratory. Since only compression results reaching the peak value of compression strength were tested, the rest of the material data are taken from reference papers [4,5].

Experimental researches
In order to obtain the strain and stress properties of concrete some experimental research was undertaken in the Mechanical Laboratory named after Prof. N.A. Belelubsky, Emperor Alexander I St. Petersburg State Transport University. The research was performed on the core samples with 53 mm diameter that were extracted from a piece of the bridge concrete structure brought directly from the site of the bridge collapse. Tests were made in December 2018 in accordance with the acting standard GOST 28570-90: Concrete. Methods of determining strength using structural samples. To prepare the stress-strain diagram and to determine the Poisson ratio 4 strain-gauge sensors (2 longitudinal and 2 transverse ones) were attached to each sample; and the MIC-036 measuring and computing complex was used to record signals (Fig. 5). As a result, the concrete compression and damage diagrams were obtained (Fig. 6, 7), which were then used for the purpose of finite element calculation.

Material parameters
In order to describe compression softening, approximation of compression curve was obtained using the formula Hsu L.S. and Hsu, C.-T.T. (1994) (Fig. 5, 6): where Е 0 represents the initial modulus of elasticity, 0 represents the peak strain ( ). The maximum value of strain is taken to be at the softening part of the curve where stress equals to 0.3 of the maximum stress.
When it comes to tension, different forms of tension stiffening models are presented in the literature. In this work a modified version of the model developed by Nayal and Rasheed (2006) [6] was used. This model was selected since it indicates similarity to the tension stiffening model needed for Abaqus CDP model. This model was originally based on homogenized stress-strain relationship developed by Gilbert and Warner (1978), which accounts for tension stiffening, local bond slip effects and tension softening. This model describes two stages of crack formation, primary and secondary. Nayal and Rasheed model and its modification [6] are presented below (Fig. 8). The modification is introduced to avoid run time errors in Abaqus. The sudden drop from maximum tensile stress 0 to 0.8 0 at critical tensile strain , as used by Nayal and Rasheed (2006) [6] and Gilbert and Warner (1978), is slanted from ( , 0 ) to (1.25 , 0.77 0 ) to avoid run time errors in Abaqus material model. After this point, the stress-strain curve follows exactly the Nayal and Rasheed (2006) tension stiffening model in both primary and secondary cracking areas but stops at (8.7 , 0.10 0 ) to avoid Abaqus run time errors. This procedure matches the procedure used in [7]. Fig. 9 and 10 shows concrete strain diagrams for 5% and 10% ratio between concrete tensile strength and concrete compression strength. Table 1 contains parameters needed for the concrete model description: β is the angle in the plane of the first and the second stress tensor invariants, f b0 /f c0 is the ratio between concrete uniaxial and biaxial stress strength, e is the eccentricity accounting for the approach speed of plasticity potential function asymptote [8]   Parameters describing the concrete damage curve are both obtained in the experiment and borrowed from [1]. They are shown in Table 2.

Finite element model
In this work, the concrete damage plasticity model was used to investigate the damage development in the bridge support that had led to collapse of the superstructure. Finite element calculations were performed in the Abaqus finite element software. In order to reduce time of calculation a submodeling technique was used. To perform this kind of analysis a global model with coarse mesh was created, whose results were used for setting boundary conditions for the submodel. The mesh of the submodel was then made finer to achieve better results. Globally, a linear material model was used with non-linear displacements and deformations. Further, a smaller part was cut from the general model, and boundary conditions were set for it. The node-based approach was applied to the submodel boundary conditions that were presented by displacements calculated within the global model. The material model was then changed to the damage plasticity model for concrete. A cubic element C3D8R (linear cubic element with 8 nodes and reduced integration) was used for concrete and a truss element T3D2 (truss element with linear displacements) was used for rebars. Fig. 11 shows the finite element model and its components. Table 3 describes load applied to the bridge support, which equals to the total of bridge structures dead load and the live load of a dump truck near the support.

Conclusion
As far as the results are concerned, it can be clearly seen that damage areas differ for materials with two varying tension properties. In the material with the higher tensile strength the damage tends to spread more compared to the material with the lower yield value where the damage is concentrated in the corner of the L-shaped beam. This damage concentration can result from numerical inaccuracy that is typical of all models with softening behavior, where results are highly dependent on the mesh. Besides that, it can be seen that stress in rebars differs for two models, especially around the area with bigger damage evolution. In the 3.3 MPa tensile strength model the stress in rebars is around 2.5 times bigger than in the 6.6 MPa tensile strength model. This can be explained by damage localization and loss of concrete strength at those locations. Further investigations are needed to check the mesh sensitivity in the model with lower tensile strength.