Numerical simulation and shear strength research of reinforced concrete beam based on peridynamics

Abstract

To explore the applicability of three-dimensional (3D) peridynamics (PD) in complex stress, 36 reinforced concrete (RC) beams without web reinforcement were designed and tested, and investigated the effects of shear span ratio, longitudinal reinforcement, and flange width on the shear strength of beams. A 3D discretization model of all specimens in the test was established, and the specimens were simulated using the PD method. To consider the heterogeneity of concrete, a non-homogeneous PD model considering aggregate size was established, and the simulation results were compared with the original model. The results indicate that the shear span ratio, longitudinal reinforcement, and flange width have a significant impact on the shear strength of RC beams without web reinforcement. The 3D PD model has a good applicability for RC beams under complex stress. Without considering computational costs, heterogeneous models can obtain more accurate results than homogeneous models and better reflect the process of concrete beam failure.

Graphical abstract

Keywords

peridynamics / numerical simulation / heterogeneous model / shear strength / flange width

1 Introduction

As one of the most commonly used building materials, concrete is widely used in building structures due to its economy and good plasticity. However, due to the brittleness of concrete itself, most of concrete structural components are in a working state with cracks, resulting in a decrease in the durability of the components [1].

As an important component of reinforced concrete (RC) structures, the strength of concrete beam components has a significant impact on the safety and serviceability of RC structures. Meanwhile, as a flexural member, the bending moment resisted by the RC beams will lead to a more pronounced problem of cracking in the lower tensile zone. The failure modes of RC beams can be roughly divided into shear failure and bending failure. Although scholars already have theories that can perfectly explain bending failure, there is still no universally agreed theory of shear behavior due to its complexity [2]. The shear strength of RC beams is provided by the aggregate interlock, residual tensile stress in the fracture process zone, uncracked compression zone, and dowel action [3,4]. However, these mechanisms are not independent of each other, and their contributions to the shear strength are also influenced by such factors as beam height, reinforcement ratio, and load position. In addition, in previous studies, the influence of flanges on shear strength has often been ignored or oversimplified [5], resulting in a serious underestimation of the shear strength of T-beams. For example, in GB 50010-2010 [6], no enhancement has been made to the shear strength of T-section beams in calculation, which can lead to errors in the evaluation of structural bearing capacity and material waste.

Based on the traditional continuum mechanics, scholars have proposed various methods to simulate the continuous deformation of concrete, such as finite element method (FEM) [7,8], finite difference method (FDM) [9], etc., which have been successfully applied in the simulation of concrete structures. However, due to the limit of solving partial differential equations based on continuum mechanics, the above methods are constrained by the required spatial derivatives when dealing with discontinuous problems, resulting in difficulties in simulation. Therefore, scholars have improved the method by proposing extended finite element [10], cohesive zone model [11], smoothed particle hydrodynamics method [12], and so on. However, the above methods require the additional fracture criteria or limitations on crack paths, and some methods require smoothing the normal of the crack surface to avoid instability of the crack path [13], which imposes certain limitations on the simulation of crack initiation and development.

The peridynamic (PD) method proposed by Silling et al. [14–16] in 2000 as an integral form of equilibrium equation avoids the use of partial derivatives. Compared to the linear elastic fracture mechanics method, PD can directly model and analyze the problem of crack initiation, and use integral form to handle discontinuities directly. Meanwhile, PD can handle the bifurcation of cracks without special assumptions, and has a natural advantage in describing discontinuous problems such as cracking and failure [17–19]. For example, Ren et al. present a dual-horizon PD formulation which is less sensitive to the spatial than the original PD formulation [20,21], and generalize it to solve nonlocal problem [22,23]. In recent years, PD theory has been proven to be applicable to various concrete structures. For example, Friedrich et al. [24] proposed a PD model suitable for simulating the failure process of quasi brittle materials under mixed loads, followed by the verification in four point bending tests, L-shaped specimen tests, and fracture toughness tests. The simulation results show that the maximum error of the simulation results for the average peak force is always < 5.0%, and the simulated failure mode is consistent with the experimental results. Liu et al. [25] established a multi field PD model to simulated the expansion and cracking behavior of concrete under sulfate attack. The results show that the coupled PD model’s computational results are in good agreement with analytical solutions, FDM, and experimental results. Lu et al. [26] proposed a three-dimensional (3D) PD method incorporating non-conservative force for brittle material and validated the theory using dynamic impact tests. The results show that this method can reasonably reflect the influence of loading rate on the deformation, strength, and cracking of brittle materials.

However, the nonlocal nature of the PD method, in contrast to local theories, leads to increased computational demands, particularly evident in 3D models [27–29]. This computational cost can pose significant challenges when simulating large-scale or complex 3D structures, limiting the practical application of PD in such scenarios. Additionally, concrete, as a heterogeneous material, exhibits inherent randomness in its mechanical properties, particularly in the spatial distribution of strength. This variability can influence failure modes and lead to discrepancies between simulation results and experimental observations [30].

To address the limitations of homogeneous models, researchers have developed various heterogeneous models for concrete [31,32]. For instance, Wu et al. [33] established a stochastic homogeneous PD material model to simulate damage initiation and crack propagation during freeze–thaw cycles. Their results demonstrated that the simulated crack propagation patterns aligned well with experimental observations, indicating the model’s capability to accurately predict the macroscopic damage evolution in concrete. Similarly, Ren et al. [20] introduced a dual-horizon PD method, which extends the applicability of PD theory to multi-material problems with minimal modifications. Their findings suggest that this approach can effectively handle material heterogeneity and complex interactions.

Despite these advancements, the computational expense associated with heterogeneous models remains a significant barrier. As a result, most studies have focused on two-dimensional (2D) numerical simulations under quasi-static loading conditions, with relatively limited exploration of non-homogeneous 3D models. This gap highlights the need for further research into efficient computational strategies and high-performance computing (HPC) solutions to enable the practical application of 3D heterogeneous PD models in simulating real-world concrete structure.

To investigate the influence of cross-sectional shape, flange width, and other factors on the shear strength of RC beams without web reinforcement, 36 RC beams without web reinforcement were designed and test. Based on test results, the influence of shear span ratio, longitudinal reinforcement, and flange width on shear strength was explored. Using HPC platform, a 3D PD model was established and simulated for 36 beams in the test. The crack propagation paths of the 3D PD simulation and test were compared, and the reliability of the method in complex stress modes was verified. Additionally, to evaluate the impact of material heterogeneity on shear strength and crack path predictions, a heterogeneous PD model incorporating aggregate size was constructed for one selected beam. The results demonstrate that the 3D PD model effectively simulates RC beams under complex stress. Furthermore, the heterogeneous model, which accounts for the spatial variability of material properties, provides higher accuracy in predicting crack paths and failure modes compared with the homogeneous model in 3D simulations.

2 Introduction to peridynamics

PD is a non-local theory that replaces the differential equilibrium equations in continuum mechanics with integral forms, allowing for the generation of discontinuities. Due to the characteristics of this method, it is suitable for simulating crack initiation and development. Its advantages include: 1) PD theory can handle discontinuous problems because it does not rely on spatial derivatives of displacement; 2) PD theory supports meshless method; 3) numerical models based on PD theory are easy to simulate complex fracture processes, such as the branching and merging of cracks. In the past two decades, PD theory has been developed and applied in various engineering and scientific fields.

As shown in Fig.1, in PD, the region to be calculated is discretized into infinitely small material points, which interact with other material points within a field (which named horizon) through bonds

(1)\begin{array}{rl}\rho (\mathit{x})\ddot{\mathit{u}}(\mathit{x},t)=& \underset{H}{\int }[\mathit{t}\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x},t\right)\\ & -{\mathit{t}}^{\prime }\left(\mathit{u}-{\mathit{u}}^{\prime },\mathit{x}-{\mathit{x}}^{\prime },t\right)]\mathrm{d}H+\mathit{b}(\mathit{x},t),\end{array}

where \mathit{x}, {\mathit{x}}^{\prime } represents material points, t represents times, \rho ( \mathit{x}) represents mass density, \mathit{u} represents displacement, \mathit{b}(\mathit{x},t) represents body load vectors, t represents force densities, H represents a local region, it can be regarded as a circular area with \mathit{x} as the center and a radius of δ, called the horizon. Among them, δ is called the radius of horizon.

The properties of the horizon are as follows

(2)\mathrm{\forall }\left|\left({\mathit{x}}^{\prime }-\mathit{x}\right)\right|> \delta \Rightarrow t \left({\mathit{x}}^{\prime }-\mathit{x},{\mathit{u}}^{\prime }-\mathit{u}\right)=0.

That is, if the distance between any material point and other material points is greater than H, they will not impact each other [34].

As a special case, it is assumed that the magnitude of the force density vector is equal and parallel to the relative position vector in the deformed state, in order to satisfy the requirement of conservation of angular momentum. The calculation method for establishing force density vectors based on this assumption is called “bond based” peridynamics (BB-PD) [13]

(3)\rho (\mathit{x})\ddot{\mathit{u}}(\mathit{x},t)=\underset{H}{\int }\mathit{f}\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x},t\right)\mathrm{d}H+\mathit{b}(\mathit{x},t),

where \mathit{f} represents force density vector, referred to pairwise response function.

Ignore the influence of temperature, the force density vector can be assumed linearly dependent on the stretch between these material points in the form [14]

(4)\mathit{f}\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x}\right)=c\mathit{s}\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x}\right)\frac{{\mathit{y}}^{\prime }-\mathit{y}}{{\mathit{y}}^{\prime }-\mathit{y}},

where c represents PD material parameters, \mathit{y}= \mathit{y}(\mathit{x},t)=\mathit{x}+\mathit{u}, {\mathit{y}}^{\prime }=\mathit{y}\left({\mathit{x}}^{\prime },t\right)={\mathit{x}}^{\prime }+{\mathit{u}}^{\prime } represents the position vector of the two points after deformation, s represents stretch, s\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x}\right)={\displaystyle \frac{\left|{\mathit{y}}^{\prime }-\mathit{y}\right|-\left|{\mathit{x}}^{\prime }-\mathit{x}\right|}{\left|{\mathit{x}}^{\prime }-\mathit{x}\right|}}.

By making the energy density of PD equal to that of classical continuum theory, it can be determined that c in 3D problems is [15]

(5)\text{c}=\frac{12E}{\text{π}{\delta }^4}.

Due to the fundamental assumption of BB-PD, the Poisson’s ratio for isotropic, linear, and micro elastic materials is 1/4.

To consider material damage, a function related to the process \mu \left(\mathit{x}-{\mathit{x}}^{\prime },t\right) needs to be introduced here, and Eq. (5) needs to be rewritten as

(6)\mathit{f}\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x}\right)=\mu \left(\mathit{x}-{\mathit{x}}^{\prime },t\right)c\mathit{s}.

In the prototype microelastic brittle (PMB) model, only two states are included: bond unbreakable and bond broken, and bond breakage is irreversible. The process related function can be represented by the following equation

(7)\mu (\mathit{x}- {\mathit{x}}^{\prime },t)=\{\begin{array}{ll}1,& s⩽s_0,\\ 0,& s>s_0,\end{array}

where s_0 represents critical stretch, \mathit{\text{μ}}= 1 when the bond is not broken, and \mathit{\text{μ}}= 0 when the bond is broken.

s_0 can be determined from the critical energy release rate G_0 of the material. This method assumes the existence of a fracture surface in a homogeneous material and calculates the energy required to break the bond passing through the fracture surface. In the 3D situation, the S_0 determined based on the critical energy release rate G_0 is

(8)s_0=\sqrt{\frac{5G_0}{6E\delta }}.

The fracture model of the above bond only considers the brittle failure of the material and is not fully applicable to the simulation of concrete materials. To consider material damage [35] and balance computational efficiency, the method of Ref. [36] is adopted here to modify \mu. The modified damage model is shown in Fig.2

(9)

where s_{max} represents maximum critical stretch.

The s_0 in the softening damage model is not the same as s_0 in the PMB model [37]. In the 3D case, assuming that the stress component on the fracture surface before fracture is equal to the tensile strength f_{\mathrm{t}} of the material and that the bonds passing through the surface deform uniformly [36]

(10)f_{\text{t}}={\int }_0^{2\text{π}}{\int }_0^{\delta } {\int }_0^{\xi } {\int }_0^{{\mathrm{c}\mathrm{o}\mathrm{s}}^{-1}\frac{\mathit{\text{z}}}{\xi }}c{s}_0{\mathit{\xi }}^2\mathrm{c}\mathrm{o}\mathrm{s}\phi \mathrm{s}\mathrm{i}\mathrm{n}\phi \text{d}\phi \text{d}\mathit{\text{z}}\text{d}\xi \text{d}\theta ,

(11)s_0=\frac{f_{\mathrm{t}}}{2E} ,

where \mathit{\xi }= \mathit{x}- {\mathit{x}}^{\prime }.

The s_{\mathrm{c}} in the softening model can be calculated by equating G_0 with the energy dissipated during the degradation phase [14]

(12)w_0(\xi )=\underset{s_0}{\overset{s_{\mathrm{c}}}{\int }} \mathit{f}\left({\mathit{u}}^{\prime }-\mathit{u},{\mathit{x}}^{\prime }-\mathit{x}\right)\xi \text{d}s,

(13)s_{\mathrm{c}}=s_0+\frac{5G_{\mathrm{F}}}{6E\delta s_0},

where w_0 represents is the work required to break a stretched bond, and G_{\mathrm{F}} is the material fracture energy.

3 Shear test of reinforced concrete beams without web reinforcement

To investigate the applicability of PD in simulating large concrete specimens, this paper designed and fabricated 36 beams, of which 9 were rectangular section beams and the rest were T-shaped section beams. The design dimensions are 2100 (1850, 1500) mm × 125 mm × 250 mm (length × web width × height). In the design and production of the specimens, the effects of longitudinal reinforcement ratio, flange width, and shear span ratio were considered separately. The details of the test specimens are shown in Fig.3 and Fig.4.

To distinguish between the specimens and description test results conveniently, 36 beams were numbered in this experiment. Taking T-700-250-1 as an example, where T represents the cross-sectional shape of the specimen is T-shape, 700 represents the shear span, 250 represents the flange width (for rectangular section beams, this item is the cross-sectional width), 1 represents the reinforcement form, and the dimensions of the specimen are detailed in Appendix A in Electronic Supplementary Material. The yield strength of reinforcement in all specimens is 447 MPa, and the compressive strength of concrete is 25.52 MPa.

3.1 Process of specimen fabrication and testing

3.1.1 Specimen production process

The specimens used in this experiment were all processed in the same factory at the same time. The steel cages were tied according to the style specified in the test plan, and steel measuring points were arranged. Wait for all specimen templates to be supported and pour concrete simultaneously, and cure under the same conditions to ensure that the curing conditions for each concrete specimen are the same. The process of specimen production is shown in Fig.5.

3.1.2 Test process

In this test, the loading device is a 500 kN electro-hydraulic servo actuator, which distributes the load to the loading point through a distribution beam. After arranging the concrete strain gauges according to the test plan, connect the steel strain gauges and concrete strain gauges to the collection device separately.

To measure the load displacement curve of the beam and obtain the deflection curve under various levels of load, linear variable displacement transducers (LVDTs) were installed at the mid span and load point location of the beam in this test; Considering the importance of mid span deflection for subsequent analysis, an additional laser displacement gauge was installed at the mid span to more accurately measure the mid span deflection. At the same time, in order to detect the warping at both ends of the beam, LVDTs were arranged on the upper surface of the beam corresponding to the support position. The placement of LVDTs is shown in Fig.6.

This test was loaded with a first level load of 5 kN until the specimen failed. During the test, if no obvious shear cracks are observed after each level of loading, the crack trend under each level of load should be recorded on the back with a marker pen. After reaching 80% of the calculated load or observing obvious shear cracks, remove the displacement meter on the bottom surface of the beam to prevent the specimen from damaging the testing instrument after shear failure. The test loading process is shown in Fig.7.

3.2 Test results and analysis

3.2.1 Test results

36 beams, except for R700-125-1 occurred concrete crushing failure in the compression zone, while the rest of the beams failed with shear failure. The photos of the failure of each specimen are shown in Appendix B in Electronic Supplementary Material, and the test results are shown in Tab.1.

3.2.2 Test analysis

From the damage photos, it can be seen that the main cracks of most beams extend along the loading point to the supports, which conforms to the failure mode of shear failure. At the same time, it can be seen that except for the R-700-500-1 specimen, all other beams have very obvious shear cracks, and the majority of the cracks in the beams extend to the flange. However, some beams, especially those with wider flanges, reaches the ultimate strength, the crack width is significantly different from the crack width at the top of the web, and even occurred crack separating from the web plate along the bottom of the flange, showing strong discontinuity.

It can also be seen that the majority of the cracks on the left and right sides of the beam are not completely symmetric, especially the width of the shear cracks, and the vast majority of beams do not experience simultaneous damage on both sides. At the same time, it can also find that in cases where the design of the two sets of specimens is not significantly different, there are sometimes significant differences in their failure process and failure mode. Due to the presence of microcracks inside actual materials compared to ideal homogeneous materials. Microcracks gradually develop and penetrate with increasing load, leading to the appearance and development of macroscopic cracks. Among them, due to the heterogeneity and mesoscopic characteristics of concrete, the weak surfaces inside it exhibit stronger randomness, resulting in a certain degree of randomness in the distribution and path of cracks in concrete beams, which has a impact on the shear strength.

From the test results, it can be seen that with the increase of shear span ratio, reinforcement ratio, and flange width, the shear strength and cracking load of the beam have a significant upward trend, indicating that the selected parameters have a significant impact on the shear strength of the beam. Among them, the influence of shear span ratio on shear strength is particularly significant. This is because when the shear span ratio decreases, it will cause a change in the shear transmission mode of the component, from beam action as the dominant force to arch action as the dominant force, resulting in a significant difference in the strength at the time of failure.

Through the failure mode and loading process of the specimens, it can be seen that all beams have undergone brittle failure. However, although some beams with larger flange widths did not occurred diagonal compression failure, their deflection was still not significant. This phenomenon also indicates the restraining effect of beam flanges on deformation, and the contribution of flanges to the shear bearing capacity of concrete needs to be discussed in more detail. Meanwhile, due to the brittleness of shear failure and the heterogeneity of concrete materials, the deflection fluctuation range of test is relatively large.

The above test results indicate that for RC beams without web reinforcement, brittle shear failure occurs even when the shear span ratio is not large, and the selected parameters have a significant impact on their shear bearing capacity. This test can provide test basis for the applicability of PD simulation of RC beams.

4 Shear failure simulation of reinforced concrete beams based on peridynamics

4.1 Peridynamic modeling process

Most concrete materials exhibit quasi elastic brittle mechanical behavior under minimal influence of multiaxial stress [38,39]. The cracks formed by beam shear failure are mixed type I and II cracks, and the stress state is relatively complex. In PD theory, the macroscopic mechanical properties of a continuum are controlled by bond parameters. Considering the softening behavior of concrete during simulation can improve the accuracy of calculations and provide a basis for the process of shear failure and future simulations under multiaxial stress conditions. fracture model of bond in Eq. (9) would be applied.

For the simulation area, it is necessary to use discretization methods to discretize the calculation area. Unlike the traditional FEM of applying boundary conditions, in PD, boundary conditions usually need to be applied to a region, and the discretized result should be discrete points with position coordinates and volume. Therefore, the simplest approach during discretization is to convert the solid element into a point (position coordinates (x,y,z), volume, number), where the coordinates are the centroid of the solid element before transformation, the volume is the volume of the solid element before transformation, and the number is the number of the original element.

Because after discretization, it only represents the centroid coordinates of the volume and points, without distinguishing the shape of the points. To ensure the quality of the mesh, the reinforcement are divided into rectangular sections of equal area to divide the mesh. By using the above method to calculate the centroid positions and volumes of each element in the mesh file, the required input file can be obtained to avoid the problem of calculating the volume of boundary points. The schematic diagram of the discrete process is shown in Fig.8 and Fig.9.

4.2 Three-dimensional peridynamic simulation result

The simulation calculations were conducted using the HPC platform of Chang’an University. The specimen is discretized with \mathrm{\Delta } x = 20 mm and δ = 60 mm. The number of material points obtained ranges from 7590 to 18090. Among them, the test results in Ref. [40], the material fracture energy G_F of concrete is 107 J/m², the elastic modulus of concrete is 3 × 10⁴ N/mm², and the elastic modulus of steel bars is 2 × 10⁵ N/mm².

Discretize the 36 beams according to the method in Fig.10, and output the reaction force at the loading point. The simulated damage cloud map is shown in Fig.10, Comparison of test results and simulation results are shown in Appendix B in Electronic Supplementary Material.

From the simulation results, it can be seen that all beams in the simulation have undergone shear failure, and their main cracks extend along the loading end to the support. For some beams with flanges, there was a situation where the flanges separated from the web plate during failure. Comparing Fig.11 with the test result in Appendix B in Electronic Supplementary Material, it can be seen that the simulation results of PD are basically the same as those of the shear failure test of RC beams without web reinforcement. The above simulation results indicate that for 3D models, compared with traditional FEM, PD can spontaneously simulate the process of crack initiation and development in simulations; The characteristic of bond fracture in PD and the basic equation in integral form avoids the difficulties in solving large deformation and highly nonlinear problems caused by excessive twisting of elements in finite elements, and can accurately simulate the entire process of failure of concrete beam components.

To better explore the differences between simulation results and test, the simulated bearing capacity results were compared with the corresponding test results in Tab.1. The comparison results are shown in Tab.2.

From Tab.2, it can be seen that the simulation results of the shear strength of most beams are in good agreement with the test results. The simulation results of PD are consistent with the test results, indicating that as the flange width, reinforcement ratio, and shear span ratio increase, the shear strength of the beam shows an upward trend. The above simulation results indicate that the 3D PD method can accurately simulate RC components under complex stress states.

4.3 Three-dimensional peridynamic simulation result

In the simulation results, it can be observed that the first flexural crack occurred in the bending region for all beams, in which the flexural tensile stress was the highest, and the shear stress was zero. As the load increased, the observed flexural crack was extended to the neutral axis, which explained that there was no shear crack in the shear span due to the absence of stress in the shear region. With the increase of the load, more flexural cracks appeared in the shear span, which gradually became more inclined and propagated toward the point load due to the presence of the stress in the shear region. These flexural-shear cracks were extended progressively to form the failure crack called diagonal shear crack, which typically extended from the support to the loading point. The above simulation results are consistent with the test results.

1) a = 400 specimen group

The shear span ratio of this group of RC beams without web reinforcement is about 2, and its failure process roughly goes through the initiation of bending cracks-crack propagation-shear crack initiation from the bottom of the beam or steel bars-shear crack propagation to the loading point-crack width propagation-specimen failure. The failure process is basically a typical shear compression failure, which is consistent with the simulation results.

In the experiment, we also observed that for T-section beams with flanges, their failure would result in multiple diagonal cracks in the web plate due to their flange constraint effect. This phenomenon was manifested in the simulation as the formation of a wide damage zone at the shear crack, which was successfully reflected in the simulation, indicating that the model can effectively consider the influence of the web plate on the shear failure process of this group of beams.

2) a = 575 specimen group

The shear span ratio of this group of beams is about 3, and its failure process is similar to that of the a = 400 specimen group. It also went through the process of bending crack initiation, shear crack penetration, and ultimately specimen failure. However, compared with the a = 400 specimen group, its bending cracks developed more completely and were closer to the form of shear compression failure. This is also reflected in the simulation.

In the test, it was able to observe that for T-section beams with flanges and beams with larger flange widths, cracks developed along the bottom of the flange, causing separation between the flange and the web plate during failure. This also occurred in the simulation. At the same time, we also noticed that compared with the a = 400 specimen group, the degree of damage to the wing edge during failure increased, which is consistent with the simulated failure process.

3) a = 700 specimen group

The shear span ratio of this group of unreinforced beams is about 4, and its failure process is the initiation of bending cracks-the propagation of bending cracks-the initiation of shear cracks-the penetration of shear cracks-specimen failure. Compared with the above two types of failure processes, its brittle nature is more significant, and it is a typical diagonal tensile failure. In the simulation, it can also be seen that the degree of damage in the pure bending section is significantly improved compared with the first two groups, which is consistent with the failure characteristics of the experiment.

Through the simulated failure process, we can see that the degree of damage to both the belly plate and the flange is significantly higher compared with the first two groups of specimens. The bending shear section of the belly plate is basically completely destroyed, which is consistent with the failure mode in the experiment. However, for this group of specimens, when the web and flange were separated, there was no significant damage area on the flange as in the test, which may also be the reason for the higher simulation results of the strength in this group.

4.4 Numerical simulation and comparative analysis of three-dimensional heterogeneous model reinforced concrete beam based on peridynamics

4.4.1 Heterogeneous modeling process

From the simulation results in the previous section, it can be seen that the homogenization model has shortcomings in simulating RC beams. Therefore, it is necessary to establish a heterogeneous model to conduct a more detailed analysis of the shear failure of RC beams. Taking the R-400-125-2 specimen as an example, explores the influence of heterogeneous models on the simulation results of shear bearing capacity of RC beams without web reinforcement.

This section will establish a heterogeneous geometric model based on the file divided into mesh in Subsection 4.1. First, a list of aggregate diameters is generated based on particle size distribution. Then, aggregates are generated based on their diameter information and placed in the geometric center of the mesh file in descending order. After checking the intersection, the aggregate placement area is marked and the mesh file is generated. Among them, completely intersecting and aggregate center element are marked as aggregate areas, partially intersecting areas are marked as interface areas, and completely non intersecting areas are marked as cement areas. The process of generating heterogeneous models using this method is shown in Fig.11. Finally, based on the generated heterogeneous model, the process method in Fig.10 is used to discretize each region.

In the process of establishing a heterogeneous model, considering the computational cost and in order to simplify the model, the influence of porosity was not taken into account in the heterogeneous model. At the same time, in order to consider the interface transition zone and eliminate the influence of excessive grid size on aggregate size changes, the strength of the bonds connected to the interface material points will be reduced.

Based on the above method, a heterogeneous model is established consisting of six types of bonds, namely aggregate and aggregate A-A, mortar and mortar C-C, transition zone and transition zone I-I, aggregate and mortar A-C, aggregate and transition zone A-I, and transition zone and mortar C-I. The schematic diagram of the discrete process is shown in Fig.12.

4.4.2 Heterogeneous model simulation results

In the process of establishing this model, the parameters of the material are referenced to the heterogeneous material parameters of materials similar to those in this article [41,42], and the volume parameter of concrete aggregates is taken as 40%. The parameters of the material are shown in Tab.3.

To include the geometric shape of the aggregate in the heterogeneous model while considering computational costs, only aggregates with a size greater than or equal to 10 mm are considered here. When discretized, \mathrm{\Delta }x = 5 mm and δ = 20 mm. The discretization results of the model are shown in Fig.13(a), and the simulation results are shown in Fig.13(b)–Fig.13(e).

It can be seen that the simulated failure mode is basically the same as the test failure mode. The simulated failure load is 145.866 kN, and the relative error of the failure load is 10%, indicating that this simulation method is suitable for simulating the shear failure of concrete beams. Meanwhile, compared to the simulation results of the homogeneous model, this simulation result can more clearly see the formation process of the main crack and more realistically reflect the experimental results.

The comparison between the load displacement curves of heterogeneous and homogeneous models and the experimental curves is shown in Fig.14. It can be seen from Fig.14 that the load displacement curves of homogeneous and heterogeneous models are in good agreement with the experimental curves. Meanwhile, the heterogeneous model can better reflect the impact of aggregate interlocking in experiments, and its simulation results are more accurate.

5 Conclusions

36 concrete beams without web reinforcement were designed and manufactured in this paper, and explores the effects of shear span ratio, longitudinal reinforcement ratio, and flange width on shear strength through four point loading tests. Meanwhile, based on this test, a 3D PD model was established and numerical simulations were conducted to verify the applicability of PD in simulating concrete beams under complex stress states. Finally, a heterogeneous model of RC beams without web reinforcement was established, and the differences between the simulation results of heterogeneous and homogeneous models, as well as the application of heterogeneous models in the simulation of RC structures, were explored. The conclusion are as follows.

1) The parameters selected in this test, especially the flange width, have a significant impact on the shear strength of RC beams without web reinforcement, which has often been overlooked in previous research.

2) By discretizing the mesh file and considering the softening behavior of quasi brittle materials, a PD model of RC beams without web reinforcement considering softening was established; a PD model considering softening behavior was established for 36 specimens, and simulations were conducted based on the HPC platform of Chang’an University. The simulation results indicate that the model has high accuracy for 3D simulation of concrete structures under complex stress conditions.

3) A non-homogeneous model of RC was established using the pick-place algorithm and discretized material point properties. Based on the HPC platform of Chang’an University, the non-homogeneous model was simulated using the R-400-125-2 specimen as an example. The simulation results show that, without considering the computational cost, the heterogeneous model outperforms the homogeneous model in simulating the failure mode and shear strength, and can effectively consider the heterogeneity of concrete materials.

References

[1] Fayyad T M, Lees J M. Experimental investigation of crack propagation and crack branching in lightly reinforced concrete beams using digital image correlation. Engineering Fracture Mechanics, 2017, 182: 487–505

[2] Filiagi PastoreM. Influence of loading arrangement on shear enhancement in reinforced concrete beams. Dissertation for the Doctoral Degree. London: Imperial College London, 2021

[3] Joint ASCE-ACI Task Committee 426 on Shear, Diagonal Tension of the Committee on Masonry, Reinforced Concrete of the Structural Division. The shear strength of reinforced concrete members. Journal of the Structural Division, 1973, 99(6): 1091–1187

[4] ASCE-ACI Committee 445 on Shear. Recent approaches to shear design of structural concrete. Journal of Structural Engineering, 1998, 124(12): 1375–1417

[5] Hawileh R A, Mhanna H H, Abdalla J A. Effect of flange geometry on the shear capacity of RC T-beams. Procedia Structural Integrity, 2022, 42: 1198–1205

[6] GB50010-2010. Code for Design of Concrete Structures. Beijing: Ministry of Housing and Urban-Rural Development of the People’s Republic of China, 2010

[7] Firouzi N, Alzaidi A S. Non-linear elastic beam deformations with four-parameter Timoshenko beam element considering through-the-thickness stretch parameter and reduced integration. Symmetry, 2024, 16(8): 984

[8] Živaljić N, NikolićŽ, Smoljanović H, Munjiza A. Numerical simulation of reinforced concrete structures under impact loading. Materials Science and Engineering Technology, 2019, 50(5): 599–610

[9] Jiang P, Zheng H, Xiong J, Wen P. Nonlinear elastic–plastic analysis of reinforced concrete column-steel beam connection by RBF-FD method. Engineering Analysis with Boundary Elements, 2021, 128: 188–194

[10] Kozlowski M, Kadela M, Gwozdz-Lason M. Numerical fracture analysis of foamed concrete beam using XFEM method. Applied Mechanics and Materials, 2016, 837: 183–186

[11] Bouziadi F, Haddi A, Tahenni T, Boulekbache B, Hamrat M, Naser M Z, Amziane S. Development of a local bond shear stress-slip model of RC beams externally strengthened with FRP materials. Journal of Composite Materials, 2023, 57(14): 2261–2285

[12] Mokhatar S N, Sonoda Y, Kueh A B H, Jaini Z M. Quantitative impact response analysis of reinforced concrete beam using the Smoothed Particle Hydrodynamics (SPH) method. Structural Engineering and Mechanics, 2015, 56(6): 917–938

[13] Wang Y T, Zhou X P, Kou M M. Three-dimensional numerical study on the failure characteristics of intermittent fissures under compressive-shear loads. Acta Geotechnica, 2019, 14(4): 1161–1193

[14] Silling S A. Reformulation of elasticity theory for discontinuities and long-range forces. Journal of the Mechanics and Physics of Solids, 2000, 48(1): 175–209

[15] Silling S A, Askari E. A meshfree method based on the peridynamic model of solid mechanics. Computers & Structures, 2005, 83(17–18): 1526–1535

[16] Silling S A, Epton M, Weckner O, Xu J, Askari E. Peridynamic states and constitutive modeling. Journal of Elasticity, 2007, 88(2): 151–184

[17] Wu L, Huang D. Peridynamic modeling and simulations on concrete dynamic failure and penetration subjected to impact loadings. Engineering Fracture Mechanics, 2022, 259: 108135

[18] Shi C, Shi Q, Tong Q, Li S. Peridynamics modeling and simulation of meso-scale fracture in recycled coarse aggregate (RCA) concretes. Theoretical and Applied Fracture Mechanics, 2021, 114: 102949

[19] Torrence C E, Trageser J E, Jones R E, Rimsza J M. Sensitivity of the strength and toughness of concrete to the properties of the interfacial transition zone. Construction & Building Materials, 2022, 336: 126875

[20] Ren H, Zhuang X, Rabczuk T. Dual-horizon peridynamics: A stable solution to varying horizons. Computational Methods in Applied Mathematics, 2017, 318: 762–782

[21] Ren H, Zhuang X, Cai Y, Rabczuk T. Dual-horizon peridynamics. International Journal for Numerical Methods in Engineering, 2016, 108(12): 1451–1476

[22] Rabczuk T, Ren H, Zhuang X. A nonlocal operator method for partial differential equations with application to electromagnetic waveguide problem. Computers, Materials & Continua, 2019, 59(1): 31–55

[23] Ren H, Zhuang X, Rabczuk T. A nonlocal operator method for solving partial differential equations. Computational Methods in Applied Mathematics, 2020, 358: 112621

[24] Friedrich L F, Iturrioz I, Colpo A B, Vantadori S. Fracture failure of quasi-brittle materials by a novel peridynamic model. Composite Structures, 2023, 323: 117402

[25] Liu Y, Li W, Guan J, Zhou X, Guo L. Fully coupled peridynamic model for analyzing the chemo-diffusion-mechanical behavior of sulfate attack in concrete. Construction & Building Materials, 2023, 409: 133874

[26] Lu D, Song Z, Wang G, Zhou X, Du X. A developed 3D peridynamic method incorporating non-conservative force for brittle materials. Engineering Fracture Mechanics, 2022, 274: 108772

[27] Liu Z, Cheng A, Wang H. An hp-Galerkin method with fast solution for linear peridynamic models in one dimension. Computers & Mathematics with Applications, 2017, 73(7): 1546–1565

[28] Hermann A, Shojaei A, Steglich D, Höche D, Zeller-Plumhoff B, Cyron C J. Combining peridynamic and finite element simulations to capture the corrosion of degradable bone implants and to predict their residual strength. International Journal of Mechanical Sciences, 2022, 220: 107143

[29] Mossaiby F, Sheikhbahaei P, Shojaei A. Multi-adaptive coupling of finite element meshes with peridynamic grids: robust implementation and potential applications. Engineering with Computers, 2023, 39(4): 2807–2828

[30] Wu P, Zhao J, Chen Z, Bobaru F. Validation of a stochastically homogenized peridynamic model for quasi-static fracture in concrete. Engineering Fracture Mechanics, 2020, 237: 107293

[31] Yaghoobi A, Chorzepa M G, Kim S S, A S. Mesoscale fracture analysis of multiphase cementitious composites using peridynamics. Materials, 2017, 10(2): 162

[32] Li W, Guo L. Meso-fracture simulation of cracking process in concrete incorporating three-phase characteristics by peridynamic method. Construction & Building Materials, 2018, 161: 665–675

[33] Wu P, Liu Y, Peng X, Chen Z. Peridynamic modeling of freeze-thaw damage in concrete structures. Mechanics of Advanced Materials and Structures, 2023, 30(14): 2826–2837

[34] MadenciEOterkus E. Peridynamic Theory and Its Applications. New York: Springer New York, 2013, 19–43

[35] Gerstle W, Sau N, Silling S. Peridynamic modeling of concrete structures. Nuclear Engineering and Design, 2007, 237(12–13): 1250–1258

[36] SillingS A. Introduction to Peridynamics. Handbook of Peridynamic Modeling. Florida: Chapman and Hall/CRC, 2016, 63–98

[37] Sheikhbahaei P, Mossaiby F, Shojaei A. An efficient peridynamic framework based on the arc-length method for fracture modeling of brittle and quasi-brittle problems with snapping instabilities. Computers & Mathematics with Applications, 2023, 136: 165–190

[38] Sfer D, Carol I, Gettu R, Etse G. Study of the behavior of concrete under triaxial compression. Journal of Engineering Mechanics, 2002, 128(2): 156–163

[39] Jin Y, Li L, Jia Y, Shao J, Rougelot T, Burlion N. Numerical study of shrinkage and heating induced cracking in concrete materials and influence of inclusion stiffness with peridynamics method. Computers and Geotechnics, 2021, 133: 103998

[40] Hobbs M, Hattori G, Orr J. Predicting shear failure in reinforced concrete members using a three-dimensional peridynamic framework. Computers & Structures, 2022, 258: 106682

[41] Li G, Zhao Y, Pang S S, Li Y. Effective Young’s modulus estimation of concrete. Cement and Concrete Research, 1999, 29(9): 1455–1462

[42] Zhang N L, Guo X M, Zhu B B, Guo L. A mesoscale model based on Monte-Carlo method for concrete fracture behavior study. Science China Technological Sciences, 2012, 55(12): 3278–3284

RIGHTS & PERMISSIONS

Higher Education Press