3D Numerical Analysis Method for Simulating Collapse Behavior of RC Structures by Hybrid FEM/DEM

Abstract

Recent years have seen an increase in demand for the demolition of obsolete and potentially hazardous structures, including reinforced concrete (RC) structures, using blasting techniques. However, because the risk of failure is significantly higher when applying blasting to demolish RC structures than mechanical dismantling, it is critical to achieve the optimal demolition design and conditions using blasting by taking into account the major factors affecting a structure’s demolition. To this end, numerical analysis techniques have frequently been used to simulate the progressive failure resulting in the collapse of structures. In this study, the three-dimensional (3D) combined finite discrete element method (FDEM), which is accelerated by a parallel computation technique incorporating a general-purpose graphics processing unit (GPGPU), was coupled with the one-dimensional (1D) reinforcing bar (rebar) model as a numerical simulation tool for simulating the process of RC structure demolition by blasting. Three-point bending tests on the RC beams were simulated to validate the developed 3D FDEM code, including the calibration of 3D FDEM input parameters to simulate the concrete fracture in the RC beam accurately. The effect of the elements size for the concrete part on the RC beam’s fracture process was also discussed. Then, the developed 3D FDEM code was used to model the blasting demolition of a small-scale RC structure. The numerical simulation results for the progressive collapse of the RC structure were compared to the actual experimental results and found to be highly consistent.

1. Introduction

With the growth of civilian development, it is necessary to demolish either potentially hazardous or no longer used buildings. Most civil structures are reinforced concrete (RC) structures, and blasting-based demolition is preferred in their destruction applications due to the economical and time-efficient natures, compared with conventional mechanical deconstruction methods. However, the blasting-based demolition may involve numerous significant risks, including the potential of causing damage to surrounding buildings and the possibility of incomplete destruction, which may result in fatalities or significantly raise the cost. As a result, considerable effort is required to optimally design the demolition process.

Optimization of structural behavior when applying the blasting is essential for a successful demolition process, and a numerical method can be one of the most effective tools because it can provide a substantial amount of useful information (e.g., stress, strain, velocity, displacement, etc.) and enables the identifying the significant parameters affecting the progressive destruction of RC structures under various conditions. As a result, numerous studies have been conducted to simulate the destruction of the RC structure to determine the optimal demolition scenario in the event of an actual collapse.

A continuum-based approach has been widely adopted to reproduce RC structures’ destruction process using finite element method (FEM) software. Most RC structural models in FEM software use the constraint approach with a linear edge reinforcement model. The constraint method models the rebar and concrete separately, and the reinforcing effect of the RC structure is reproduced by constraining the movements of adjacent components when they overlap. In the FEM RC structural model, the failure process of the RC structure resulting in the collapse is expressed by deleting or fracturing the main load-bearing structural components. The frame structure-blasting demolition processes was simulated using FEM implemented in LS-DYNA to investigate the mechanism of collapse and its applicability to practical demolition modeling [1]. In [1], weakening of the structural part were compared for the cases with and without explosive loading (i.e., structural part deletion), and the results indicate that explosive blasting applied on the cut site is more comparable to the examination of the actual demolition. Meanwhile, FEM implemented in ANSYS software was also used to model the explosive demolition of old cylindrical silos [2]. This study provided guidance for pre-weakening cylindrical silos prior to blast destruction based on the findings obtained by structural analysis and numerical simulation. Due to the large number of reinforcement bars to represent in RC structures, efforts to simplify the reinforcement bars in the modeling of RC structure demolitions were also presented [3]. In [3], a simple structural analysis model was used to verify the proposed numerical approach in comparison to the experimental outcome. However, the simplified numerical modeling for demolition had limitations in simulating RC concrete fracture and fragmentation. The study for considering the effect of the size of FEM elements in structural blasting demolition simulation using FEM in LS-DYNA can be found [4]. To investigate the element size effect on the response of overall RC structure, the FEM models with element sizes of 20 cm and 40 cm were developed in this work. They demonstrated that the coarse element model reduced the structure’s collapse time due to an exaggeration of the structural impact effect, resulting in more severe structural damage than the fine element model. Most significantly, most continuum-based demolition simulations have used the element deletion approach to model RC concrete fracture processes. In [5], it was found that the simulation of material failures process via element erosion did result in the elimination of almost whole the structure, reducing the accuracy and reliability of the final result. They also demonstrated the disadvantages of the element erosion method and proposed a node split algorithm with a preprocessor for the FEM in LS-DYNA to solve the shortcoming. Using the node split technique, the destruction of the Sparkasse building in Hagen was reconstructed reasonably. Nonetheless, the findings obtained from the numerical simulations still had limits, since the abrupt separation of element connections caused by fracturing resulted in unrealistic and spurious behavior (e.g., shattering of structural parts).

On the other hand, the applied element method (AEM) [6], which is an effective numerical simulation method for the nonlinear analysis of framed and continuous structures, was recently developed as a feasible alternative to the conventional FEM for modeling the progressive failure process of the RC structure. In the AEM, the fracture process of a concrete structure is simulated using normal and shear springs between rigid elements, and the stress and strain in specific portion of the structure are also calculated using the two spring types. Numerous publications exist on modeling the demolition of RC buildings using the AEM. Especially, the AEM-based commercial software “Extreme Loading for Structures (ELS)” is commonly used to model the collapse of RC structures. Several researchers have simulated the demolition of small-scale RC structures to evaluate the capability of the AEM to analyze the progressive collapse of the RC structure. In [7], the AEM in ELS software was applied to model 1/5 scaled three and five stories of RC structures and compared the results to experimental results. The horizontal and vertical displacements obtained from these simulation results were similar to the experimental results. The AEM in the ELS software was also used to simulate the collapse of a prototype RC moment-resisting frame building induced by seismic loading [8]. The simulation results were compared to shake table test results and found to be in good agreement. The collapse mode was thoroughly investigated, and the findings were used to train a deep neural network (DNN) to predict the collapse mode of an RC structure. Several applications of the AEM for designing demolition sequences to actual large-scale RC structure collapse can also be found. Progressive collapse analysis of high-level and shear-wall apartment buildings subjected to blast loading was also simulated to demonstrate the capability of the AEM in ELS software [9]. The simulation results showed good agreement with observed demolition behavior. The mechanism of the progressive collapse of a structural component caused by local damage was addressed in [9] using several factors such as reinforcement, axial force and moment. The progressive structural collapse of the Perna Seca Hospital building in Brazil was also simulated to verify and optimize the demolition scenario using AEM in ELS software and a fully detailed model [10]. Ground vibration caused by debris impact during the demolition process was thoroughly investigated in order to minimize damage to adjacent buildings. Finally, they compared the AEM result to the actual observation, concluding that the AEM result closely matches the actual result in terms of the structure’s collapse pattern. However, in comparison to the FEM, the AEM is less accurate in the small displacement range and requires a large amount of CPU time [6,11]. Additionally, the AEM is limited in its ability to model Poisson’s effect because it is based on the assumption of rigid elements in which the deformation is not allowed.

In recent years, the combined finite-discrete element method (FDEM) [12] has gained popularity as a useful numerical method for evaluating civil and rock engineering issues involving fracturing and contact processes between fragments, i.e., discrete bodies. The FDEM is able to model continuum behavior, fracturing (i.e., crack initiation and propagation) and contact between material surfaces including newly generated fracture surface. These features indicate that the FDEM may be used to simulate the demolition of reinforced concrete buildings. There are numerous research publications available to simulate the fracturing of solid structure with reinforcement based on the FDEM. Most publications on reinforcement in RC structures using the FDEM have incorporated the linear edge reinforcement model into the finite concrete model via the embedded method [13,14,15,16,17,18,19,20,21]. In other words, the reinforcement model is conceptually separated by interfaces of the solid finite elements corresponding to the concrete (hereafter, concrete elements) part and remains fixed within the concrete elements. The reinforcement model can be deformed only through the deformation of the concrete element. Previous research [13,14,15,16,17,18,19] have actively developed the reinforcement model by integrating FDEM and solid elements containing rebar, i.e., embedded rebar elements, and interface elements connecting the rebars, i.e., rebar joint elements. In the embedded rebar finite elements, a linear-elastic model is used, whereas the constitutive model in the rebar joint elements is applied by classifying the before and after the occurrence of a crack in the concrete. When no concrete separation occurs, the penalty function assures the continuity of the rebar finite components in this model [13]. After concrete cracking, on the other hand, the rebar model in the joint element adopts a path-dependent mechanical model [22]. The developed rebar model was validated through simulations of several simple experimental tests, including three- or four-point flexural loading and tension loading under monotonous or cyclic static and dynamic loading conditions. Good agreement between the FDEM and experimental results was obtained for the simulation of the fracture processes in the RC structures such as RC walls, Prothyron structures and exterior RC beam–column joints. Meanwhile, in [20,21], a reinforcement model for rock structure reinforcement via linear rock support (i.e., rockbolts, dowels, cables, etc.) was developed in the 1D reinforcement model. The reinforcement model includes a reinforcement finite element and reinforcement joint. In these works, the axial stress of the reinforcement finite element is calculated under the assumption of linear elasticity. The axial stress of the reinforcement joint is computed by selecting between the linear elasticity and elasto-plastic model, and tangential stress of the reinforcement joint is calculated by the linear elasticity model. Performance of the 1D support model was compared to the explicit reinforcement bars model (i.e., FDEM rebar model) by simulating the pull-out load and double shear test demonstrating excellent agreement. The results of laboratory double shear tests [23] were also compared to numerical simulations and were found to be in good agreement in terms of the damage of the concrete and the rebar deformation. Recently, the so-called two-node bolt element model with a spring-slider shear connector has been further developed into the GPGPU-based FDEM code in order to rapidly simulate the reinforcement of rock mass around underground excavations by fully-grouted rockbolts [24]. In the work, numerical instability caused by specific intersections between triangular mesh and rockbolts is resolved by adopting only the equivalent length of the bolt finite element without separation by the finite rock element. The spring-slider shear connection is also used to simulate the existence of a grout annulus between the rockbolt and the surrounding rock. To validate rockbolt logic, the pull-out tests were simulated, and rock bolting in underground excavations was successfully reproduced. However, existing research are mostly limited to two-dimensional (2D) study. Although the previous work [25] described the fracturing of the three-dimensional RC structure model, a detailed description was not provided in the paper. Therefore, this work attempts to expand to the analysis of 3D reinforced structures. Additionally, because most demolition with blasting for RC structures occurs in 3D manners rather than 2D, the development of a numerical simulation tool for the demolition of 3D RC structures is significantly important.

Based on the above research background, this study develops and applies the GPGPU-accelerated 3D FDEM to simulate the collapse of the RC structure. Here, as a fundamental study, we developed and implemented the GPGPU-accelerated 3D FDEM code with a 1D reinforcement bar (rebar) model. The three-point bending test for RC beam was simulated to validate the developed 3D FDEM code, in which the input parameters for the concrete material model were first calibrated appropriately by simulating the uniaxial compressive strength (UCS) and Brazilian tensile strength (BTS) tests to ensure an accurate simulation of fracturing process in RC. In addition, to examine the effect of element size on the fracture process of the RC structure, the three-point bending test was simulated for the RC beam with different finite concrete elements. The fracture process of the RC beam was discussed by comprehensively investigating the fracture process of the reinforced and plain concrete beams in the three-point bending test. Finally, the developed code was also applied to simulate the progressive collapse of a small-scale RC structure to show its applicability for blast-induced demolition RC structure.

2. Development of the GPGPU-Accelerated 3D FDEM Code for the RC Structure

2.1. GPGPU-Accelerated 3D FDEM for Modeling Concrete Fracturing

A 3D FDEM code has been developed to model the fracture of rock-like materials such as rocks and concrete. This code was originally developed in a C++ programming language and is based on the open-source Y2D/Y3D libraries. This code has been successfully used to deal with a variety of rock engineering issues by modeling the fracture of rock and rock-like material [26,27,28]. However, since the applied algorithms in the original Y-code were sequential, its prior applications were restricted to small-scale 2D issues with relatively coarse meshes. To further improve computation processing such as speed and capabilities, the FDEM code has incorporated a parallel computing using a GPGPU device controlled by computing unified device architecture (CUDA) C++ [29,30].

The FDEM is primarily based on continuum mechanics, the cohesive zone model (CZM), and contact mechanics within an explicit FEM framework. In this section, the only necessary component applied in this paper is addressed, while other comprehensive concepts of the FDEM code can be found in literature [29,30].

In 3D FDEM, the continuous deformation including the stress wave propagation in the intact regime of materials such as concrete is modeled using a continuum 4-node tetrahedral finite element model (TET4s) (Figure 1) [29]. The CZM is applied with the concept of smeared crack to describe the transition from continuum to discontinuum [31]. The tensile and shear softening curves are applied using the assembly of 6-node cohesive elements (CE6) (Figure 1), which are initially zero-thickness, and inserted into all boundaries of the TET4s at the start of the simulation to model the behavior of the fracture process zone (FPZ) in front of the crack tip. Tensile and shear softening curves as a function of crack opening (o) and sliding displacements (s) control the opening and sliding of the faces in each CE6 (Figure 2). The normal and shear cohesive tractions (${\sigma }_{coh}$ and ${\tau }_{coh}$) acting on each side of CE6 are computed using Equations (1) and (2), respectively.

(1)$\sigma =\left\{\begin{array}{cc}\frac{2o}{o_{overlap}}{T}_s,& \mathrm{if} o<0\hfill \\ \left[\frac{2o}{o_p}-{\left(\frac{o}{o_p}\right)}^2\right]f(D)T_s,& \mathrm{if} 0\le o\le o_p\hfill \\ f(D)T_s,& \mathrm{if} o_p<0\hfill \end{array}\right.$

(2)${\tau }_{coh}\left\{\begin{array}{cc}\left[\frac{2\left|s\right|}{s_p}-{\left(\frac{\left|s\right|}{s_p}\right)}^2\right](-{\sigma }_{coh}\tan \varphi +f(D)c),& \mathrm{if} o\le \left|s\right|\le s_p\hfill \\ -{\sigma }_{coh}\tan \varphi +f(D)c,& \mathrm{if} s_p<\left|s\right|\hfill \end{array}\right.$

where $o_p$ and $s_p$ denote the corresponding artificial elastic limits of $o$ and $s$. $o_{overlap}$ is the representative overlap, which controls the artificial elastic behavior of CE6s when $o$ is negative, i.e., overlapping. $T_s$, $c$ and $\varphi$ are the microscopic tensile strength, the microscopic cohesion and the microscopic internal friction angle of the CE6, respectively. The $o_p$, $s_p$ and $o_{overlap}$ are specified as $o_p=2h{T}_s/P_{overlap}$, $s_p=2hc/P_{tan}$ and $o_{overlap}=2h{T}_s/P_{overlap}$, respectively. $h$ denotes the size of the CE6, and $P_{open}$, $P_{tan}$ and $P_{overlap}$ indicate the CE6’s artificial penalty factors for opening, sliding and overlapping. Equation (2) is equivalent to the Mohr–Coulomb (MC) shear strength model with a tension cut-off. The function $f\left(D\right)$ is the characteristic function for the tensile and shear softening curves (Figure 2), and it varies with respect to the damage value D of CE6. D not only reflects the opening and sliding of the faces of CE6, but also mixed-mode crack behavior according to the following equation [31]:

(3)$f(D)=\begin{array}{c}\left[1-\frac{A+B-1}{A+B}\exp \left(D\frac{A+CB}{(A+B)(1-A-B)}\right)\right]\hfill \\ \left[A(1-D)+B{(1-D)}^C\right] (0\le \mathrm{D}\le 1)\hfill \end{array}$

(4)$$\begin{gathered} D=\min \left[1,\sqrt{{\left(\frac{o-o_p}{o_t}\right)}^2+{\left(\frac{\left|s\right|-s_p}{s_t}\right)}^2}\right] \\ \mathrm{if} o>o_p \mathrm{or} \left|s\right|>s_p, \mathrm{otherwise} 0 \end{gathered}$$

where $A$, $B$ and $C$ in Equation (3) are the intrinsic rock properties that define the shape of the softening curves, and $o_t$ and $s_t$ in Equation (4) denote the critical $o$ and $s$ values where CE6 fractures and forms macroscopic fracture surfaces. Meanwhile, the $o_t$ and $s_t$ in Equation (3) are computed based on the mode I and mode II fracture energies $G_{fI}$ and $G_{fII}$, respectively, based on the following equations [31,32]:

(5)$G_{fI}={\displaystyle \underset{o_p}{\overset{o_t}{\int }}{\sigma }_{coh}(o)do}$

(6)$G_{fII}+W_{res}={\displaystyle \underset{s_p}{\overset{s_t}{\int }}{\tau }_{coh}(\left|s\right|)d\left|s\right|}$

where $W_{res}$ is the work performed by the residual stress component in the MC shear strength model per the area of CE6. Note that the current formulation does not distinguish between mode II and III fracture modes. Mode II and mode III reactions to microcracks (i.e., CE6) are considered to be simply described by the parameters $G_{fII}$, since defining the crack tips and planes in the three-dimensionally complex fracture networks generated by FDEM is challenging.

2.2. Embedded Reinforcement Bar Model in 3D FDEM Concrete Model

A reinforcement bar (rebar) model is integrated into TET4s and CE6s and implemented in 3D FDEM code to simulate rebar reinforcement in the RC structure. A one-dimensional rebar model is adopted, and it is placed within the TET4s. Figure 3 depicts the RC structural model overlapped with the rebar model. Each rebar object is defined by its start and end points, and the entire rebar is further discretized into two-node rebar elements (REs) and two-node joint rebar elements (RJEs) based on the intersection points between the TET4s and rebar model interfaces.

The RE is embedded within the TET4 because the intersected RE’s two points with the TET4 are fixed. Thus, the deformation of TET4 influences the deformation of RE, while the axial strain of RE, ${\epsilon }_{axial}^{RE}$, is calculated using Equation (7) while the RE’s axial force, $F_{axial}^{RE}$, is derived by applying the linear-elastic constitutive model.

(7)${\epsilon }_{axial}^{RE}=\frac{l_c^{RE}-l_i^{RE}}{l_i^{RE}}$

(8)$F_{axial}^{RE}=A^{\mathrm{Re}bar}{\sigma }_{axial}^{RE}=A^{\mathrm{Re}bar}{E}^{\mathrm{Re}bar}{\epsilon }_{axial}^{RE}$

where $A^{Rebar}$ and $E^{Rebar}$ denote the cross-sectional area and Young’s modulus of the each rebar. $F_{axial}^{RE}$ is assembled as the internal force of the corresponding TET4. It is worth noting that each RE is embedded in its host TET4 during the FDEM simulation.

When cracks initiate between REs, the initially zero-length two-node RJE begins to influence the mechanical behavior of concrete. The induced forces by RJEs consist of two components, i.e., axial force ($F_n^{RJE})$ and shear force ($F_s^{RJE}$). The $F_n^{RJE}$ and $F_s^{RJE}$ are derived by the axial stress (${\sigma }_n^{RJE})$and shear stress (${\tau }_s^{RJE})$, which are determined as a function of the relative opening ($o$) and slip($s$) between two nodes of RJE (Figure 4), respectively, as presented in Equations (8) and (9).

(9)$F_n^{RJE}=A^{\mathrm{Re}bar}{\sigma }_n^{RJE}=A^{\mathrm{Re}bar}{K}_n^{\mathrm{Re}bar}o$

(10)$F_s^{RJE}=A^{\mathrm{Re}bar}{\tau }_s^{RJE}=A^{\mathrm{Re}bar}{K}_s^{\mathrm{Re}bar}s$

where $K_n^{Rebar}$ and $K_s^{Rebar}$ denote normal and tangential stiffness of each rebar, respectively. Equivalent constitutive models in [20,21] were used for computing the ${\sigma }_n^{RJE}$ and ${\tau }_s^{RJE}$. ${\sigma }_n^{RJE}$follows the elastic-perfect plastic model, while ${\tau }_s^{RJE}$ follows the elastic model. As can be seen in Figure 5, the ${\sigma }_n^{RJE}$and ${\tau }_s^{RJE}$ in the RJE initially behave like elastic material, while the ${\sigma }_n^{RJE}$ remains constant when it reaches the yield stress of the rebar. The yield stress is determined by a user as the mechanical properties of the rebar. Finally, the $F_n^{RJE}$ and $F_s^{RJE}$ are assembled into the cohesive force in the CE6.

3. Verification of Developed 3D FDEM for Modeling of RC Structure Collapse

3.1. Determination of the 3D FDEM Input Parameters for Concrete Model

Three-dimensional FDEM input parameters of the concrete model should be established to reproduce its fracture behavior reasonably. In our previous research [33], the UCS and BTS tests, which are often used to evaluate the mechanical characteristics of concrete, were simulated for determining numerical parameters using the 3D FDEM. The FDEM simulation results were compared to the experimental results to calibrate the FDEM input parameters. Concrete specimens with a typical UCS of 35 MPa were used in the calibrating process, and the physical and mechanical properties of the concrete sample evaluated through experiments are summarized in Table 1. The elastic parameters for the 3D FDEM models were determined using the values in Table 1. The calibrated parameters for the FDEM simulation of the concrete structure in [33] are listed in Table 2 and Table 3, and these are used in the next section.

3.2. Modeling Procedure of Three-Point Bending Test for RC Beam Fracturing

The three-point bending test of a supported RC beam subjected to monotonically increasing three-point loads [34] was selected (Figure 6) to validate the rebar model. The following FDEM simulations use the calibrated values listed in Table 2 and Table 3, because the mechanical properties of 35 MPa-Concrete used in [34] are similar to those of 35 MPa-concrete used in our previous research [33]. Physical and mechanical properties of the rebar such as diameter (d = 8 mm), elastic modulus ($E_{steel}=200GPa)$ and yield stress (${\sigma }_y=300MPa$) are used as the input parameters in the rebar model. Figure 7 depicts the 3D FDEM model, including concrete, and this model is used to simulate the three-point bending testing for the RC beam. Since FDEM is a mesh-dependent method in which crack paths are affected by element size, several average mesh sizes, $h_{ave}$ (=100 mm, 50 mm, 25 mm and 12.5 mm), are chosen to investigate the effect of the element size on the fracture process of the RC structure. The load was applied with a critical damping scheme to achieve a quasi-static loading state [29] using a loading rod whose velocity is kept to 1 mm/s. The simulation is continued until the RC beam is collapsed. The FDEM simulation results analyses were finally compared to the experimental results [34].

3.3. Simulation Results of the Three-Point Bending Test

Figure 8 shows the numerical simulation results of the three-point bending test for the RC beam in which the RC model’s stress distribution and crack propagation are presented. At first, the point load is applied to the RC beam through contact with the top-loading rod with the constant velocity, while the bottom-loading rod remains fixed. After further displacement of the loading rod, the vertical loading to the RC specimen linearly increases due to the specimen’s elastic deformation, as shown in Figure 8a. The tensile stress is concentrated at the bottom part of the beam, while the microscopic crack is not presented at this stage. With the additional displacement of the loading rod, the axial loading to the reinforced concrete specimen continues to increase linearly, and tensile fracturing finally occurs at the bottom of the beam, as presented in Figure 8b. With the further loading, the tensile crack propagates upward (Figure 8c). After that, additional fracture initiation and propagation occurs along with the rebar due to the further deflection of the rebar around the center (Figure 8d). Then, the reinforced concrete beam collapses completely, resulting in significant concrete fracturing and separation. Here, the crack is propagated along with the rebar model due to severe plastic deformation of the rebar (Figure 8e). The resultant crack patterns for four different element sizes ($h_{ave}$ = 100 mm, 50 mm, 25 mm and 12.5 mm) are shown in Figure 9. The case of $h_{ave}$ = 100 mm shows that the RC beam is split into two halves (Figure 9a), and fragments move almost parallel to the applied loading, which is not reasonable compared with the experimental observation [34]. On the other hand, the remaining cases with $h_{ave}\le$ 50 mm exhibit a similar crack pattern, with numerous cracks distributed throughout the lower part in the model (Figure 9b,c), which shows better agreement with the experiment [34]. Hence, it is indicated that the value of $h_{ave}\le$ 50 mm is at least required to obtain an acceptable fracture pattern in the case of this three-point bending simulation.

Figure 10 illustrates the applied load (P)—displacement (d) curves for the RC beam model with various $h_{ave}s$ during three-point bending testing simulation. As shown in Figure 10, when $h_{ave}$ = 100 mm, the RC beam loses its bearing capacity at d = 1 mm, while the remaining cases exhibit nearly identical behavior. The sudden drop in the P in the case of $h_{ave}$ = 100 mm is owing to the crack propagation mode identified as splitting mode (Figure 9a). On the other hand, for the cases with $h_{ave}$ = 50 mm, 25 mm and 12.5 mm, the shapes of the P-d curves are almost similar both for the elastic regime and plastic regime for d > 1 mm. Thus, at least from this problem, it may be concluded that the FDEM simulations incorporating the rebar element are less sensitive to $h_{ave}$ as long as $h_{ave}$ is taken small enough in the way that the fracture process and pattern are accurate.

The fracture behavior of the plain concrete without rebar is also analyzed to investigate the reinforcement effect of the applied rebar model (i.e., RE and RJE). Figure 11 presents the mid-span loading–displacement curves for plain and reinforced concrete structures obtained from FDEM simulations. As can be seen the specific point “A” to “B” in Figure 11, the mid-span loading is linearly raised up to 15.12 kN in the elastic deformation regime for both cases. During the elastic deformation, the tensile stress is concentrated at the bottom part of the beam. With the further displacement, the mid-span load of the plain concrete beam is dropped significantly after the peak stress, indicating typical brittle concrete fracturing (“C” in Figure 11). On the other hand, the stiffness of the RC beam is slightly decreased due to the crack initiation and propagation of concrete. When the displacement is increased further, the loading shows a continuous increment due to the reinforcement effect of the rebar and a small variation in loading due to the concrete crack initiation and propagation (“D” in Figure 11). As the concrete fragments, it no longer supports axial loading, and only the yield stress of the rebar resists axial loading (“E” in Figure 11)

3.4. Comparison of Fracture Behavior in FDEM with Experimental Results

Validation of the developed 3D FDEM for simulating the collapse of RC structures requires comparison to actual fracture behavior. The mid-span loading–displacement curves and resultant fracture pattern for the RC beam between FDEM and the experiment are compared in Figure 12 and Figure 13, respectively. From these comparisons, it can be found that the FDEM can capture both the mid-span loading-displacement curve and fracture pattern in experiments. It is worth noting that the crack pattern presents higher consistency with the experiment result if more refined mesh with $h_{ave}$ $\le$ 25 mm is used (see Figure 9b,c and Figure 13a). Since concrete includes numerous aggregates and cement in general, fracturing easily occurs near the interfaces between cement and aggregate. Besides, because the concrete commonly uses a maximum coarse aggregate size as 25 mm [35], the 35 MPa-concrete material modeled in this study may be independent to the mesh dependency when simulating its fracture process with a finer mesh of $h_{ave}\le$ 25 mm. Therefore, the careful selection of the element size for the concrete considering the internal structure of target concrete is also significant to properly simulate the fracture process of the RC structure. Based on the extensive analysis and comparison with the experimental results of the fracturing process of the RC beams, the developed 3D FDEM is successfully validated, and further applications for simulating the collapse of the RC structure using the proposed 3D FDEM may be possible.

4. 3D FDEM Modeling of Demolition of a Small-scale RC Structure

4.1. Modeling Procedure of Small-Scale RC Structure Collapse

In this section, we show the blasting demolition of a small-scale RC structure model to investigate the applicability of the developed 3D FDEM code with a rebar model to the RC structure demolition. As shown in Figure 14a, the RC structure is modeled by piling up to three stories as a 1 × 1 bay structure, including one slab, four beams, two sub-columns and two main columns. The ground was also modeled as a fixed thin rigid plate consisting of TET4s to sustain RC structure. The rebar is embedded in each part of the RC structure, as shown in Figure 14b. The RC structure is discretized as 22,624 of TET4s for concrete, and the RE and RJE rebar elements are discretized as 841 and 1680, respectively. Although the developed code can deal with a massive number of elements using GPGPU, here, the main focus is not on the computational performance but rather on verification and validation using small-scale models, and the FDEM simulation of demolition of massive RC structure is considered as our future work. The input parameters listed in Table 2 and Table 3 are again used to simulate the demolition of RC structure. The numerical simulation is performed with two stages (i.e., initial static analysis and dynamic analysis). At the first stage, to realize the initial static stress state in the RC structure under gravity, the gravity acceleration (=9.81 $\mathrm{m}/{\mathrm{s}}^2$) is applied to the RC structure along the negative y-direction in Figure 14 until the stress equilibrium statement is achieved. Here, to accelerate the convergence to the static stress state, the local damping scheme coupled with an optimized mass scaling technique is applied [30]. Then, at the second stage for dynamic analysis, to mimic the real blasting demolition of the RC structure, its two-main columns in the first floor are instantaneously removed at the beginning of the simulation, which triggers the collapse of the RC structure due to the loss of stress equilibrium. During the dynamic analysis, gravity is remained to be applied to the RC structure.

4.2. Simulation Results of the Demolition of a Small-Scale RC Structure

Figure 15 compares the demolition process of an RC structure between 3D FDEM simulation and experiment observation [7] at equivalent elapsed times after removing the two-main columns in 1st floor due to blast load. Note that, although the two columns were removed by a detonator in the experiment, 3D FDEM simulation just removed the corresponding parts at once, which should be acceptable since the time-scale of interest is several hundred to thousand milliseconds (ms), while the detonation-induced removal completes roughly within several hundred microseconds. At 400 ms, the RC structures both in FDEM and the experiment have clearly started to tilt, and it continuously inclines until 600 ms due to gravity. With the further inclination, the newly created surface due to the removal of the columns on the 1st floor dynamically contacts with the ground. At 1000 ms, the RC structure overturns using the initial contact part as a rotation axis. After that, one of the original RC structure surfaces finally makes contact with the ground at 1200 ms. Then occurs the sliding of the contact face against the ground at 1400 ms in both the FDEM and experiment. Up to this stage, the FDEM and experiment show very similar behavior. However, after this stage, the RC structure is broken into multiple segments in the experiment but is nearly preserved in the numerical simulation at 1600 ms. A possible source of difference between experiment and simulation at this stage is the reinforcement bar’s formulation. The beam and column were connected in the experiment using epoxy resin. However, in FDEM, the beam and column are reinforced, which makes the nominal strength of the RC structure in the FDEM higher than that in the experiment.

The time profile of x- and y-displacement at the specific point “A” (see Figure 14a) on the RC structure model is also compared as shown in Figure 16 and Figure 17. It can be found that the FDEM and experiment are observed after 800 ms, which could be explained by the aforementioned reason (i.e., experiment used the epoxy resin to connect the beams and columns). As illustrated in Figure 16, the differences in the direction of the y-axis between numerical simulation and experiment begin at 800 ms. From 800 ms on, the y-axis displacement in numerical simulations moves more easily toward the ground as the equal time elapses. At 1000 ms, Figure 15 illustrates this behavior. Both the experimental and numerical models are collapsing after contact with the ground at 1000 ms. However, the experimental model exhibits more brittle deformation than the numerical model. The significant deformation in the bottom part of the numerical model, in particular, results in differences in the y-axis displacement between them. On the other hand, the differences in the direction of the x-axis in Figure 17 demonstrate a similar behavior between experiment and numerical simulation. However, after 1600 ms, the experimental model continues to move due to structure collapse, as illustrated in Figure 15.

Our future study includes conducting the blasting demolition experiment that realizes the almost same condition as in the FDEM simulation in Figure 14. Thus, this paper does not discuss this demolition any further in detail.

5. Conclusions

In this study, the GPGPU-accelerated 3D FDEM has been developed with the introduction of the 1D rebar model to analyze the demolition process of the RC structure in detail. The fracturing of the concrete is modeled on the framework of the 3D FDEM and can reproduce the fracture process of the concrete by modeling continuum behavior, crack initiation and propagation and contact between discrete elements. The reinforcement of the RC structure is also modeled using the 1D rebar model, which includes two-node REs and two-node RJEs. As the constitutive model, linear elastic and elastic-perfect plastic models were applied to describe stress statements in the RE and RJE induced by its deformation. The 3D FDEM parameter was firstly determined based on our previous study [33] to properly simulate the fracture process of the plain concrete by simulating fundamental mechanical tests, i.e., UCS and BTS test. On the other hand, RE was modelled as 1D linear elastic material while RJE was modeled as 1D elastic-plastic material. The three-point bending test for the RC beam was simulated to verify the proposed 3D FDEM. Different element sizes in the FDEM mesh for the concrete model were constructed and investigated to consider the size effect of element size on the fracture behavior. The mid-span loading displacement curves and resultant fracture pattern for the RC beam were also discussed based on the results of the three-point load testing for the reinforced and plain concrete beams. In addition, the blasting demolition of the small-scale RC structure, composed of the three stories as a 1 × 1 bay structure, was also simulated to verify the applicability of the proposed FDEM for simulating the progressive collapse of the RC structure. The following conclusions can be taken from this study:

• (1). The fracture process of the concrete part in RC structure is highly dependent on the FDEM element size and, depending on the selected element size, exhibits different mechanical behavior when subjected to three-point loading. The RC beam model with an average element size of 100 mm showed a splitting mode whose direction of main fracture was parallel to the three-point loading direction, which is not observed in reality. On the other hand, with the decrease of the element size (50 mm, 25 mm and 12.5 mm in this study), the obtained crack propagation pattern became more and more realistic compared with the experiment. Therefore, the optimal size of the element size for concrete must be determined by carefully checking the effect of the element size on the fracture process of concrete in the RC structure. Once the proper element size is determined, the rebar model incorporated in this paper was found to be less sensitive with respect to the element size. In addition, by comparing the cases with and without the reinforcement, brittle fracturing occurred in the former case, while the latter case showed the clear effect of reinforcement characterized by its yielding. The loading-displacement curves and resultant fracture pattern for the RC beam obtained by FDEM were compared with the experimental results, and these exhibit a high degree of agreement. It may be concluded that the developed 3D FDEM can be used to simulate the collapse process of the RC structure effectively.

• (2). The FDEM simulation of collapse process in the small-scale RC structure including the beams, slabs, and columns for blasting demolition was successful in terms of the reinforcing effect of each component of the RC structure. When compared to the history of displacement at a specific point in the RC structure, the numerical simulation results agreed well with experimental observation to some extent. However, the differences between the FDEM simulation and the experiment were also noticeable after the falling RC structure dynamically contacts with the ground. This may be owing to the fact that the beam and column were connected in the experiment using epoxy resin while the FDEM simulation did not consider this effect. As a future work, additional blasting demolition tests for the small-scale RC structure intended for the comparison of 3D FDEM simulation should be conducted.

• (3). Further verification of the 3D FDEM with the 1D rebar model should also be performed to ensure accurate analysis of the collapse RC structure by simulating the mechanical load test for RC structure components such as the beam–column joint and slabs. Because only the three-point bending test for RC beams was simulated in this study for validation, the result may not be sufficient to ensure that our methodology can fully simulate the collapse of complex RC structures.

Author Contributions

Conceptualization, G.M., D.F. and S.C.; methodology, G.M., D.F. and S.C.; investigation, G.M., D.F. and S.C.; writing—original draft preparation, G.M., D.F. and S.C.; writing—review and editing, G.M., D.F. and S.C.; visualization, G.M. and S.C.; supervision, S.C.; project administration, S.C. All authors have read and agreed to the published version of the manuscript.

Funding

This Work was supported by National University Development Project in 2019.

Acknowledgments

The authors like to express their gratitude to Hongyuan Liu at the University of Tasmania in Australia for establishing a framework for the development of 3D FDEM and offering invaluable advice.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures and Tables

Figure 1. Two types of finite elements (TET4 and CE6) in 3D FDEM.

Figure 2. Tensile (left) and shear (right) softening curves in the ICZM.

Figure 3. The RC structure as a numerical representation (separation between TET4s is exaggerated for visualization).

Figure 4. The opening and sliding definitions for the RJE (point B-C in the figure).

Figure 5. RJE constitutive model in (a) axial and (b) tangential directions (arrows indicate loading, unloading and reloading paths).

Figure 6. Schematic diagrams showing (a) the RC beam and (b) the center cross-sectional view of the beam for the three-point bending test.

Figure 7. Numerical simulation model for the three-point bending test of the RC beam: (a) the specimen geometry and (b) the rebar inserted in the specimen.

Figure 8. Spatial distribution of maximum principal stress (left column) and macroscopic crack (right column) in the RC beam subjected to the three-point loading.

Figure 9. The resultant crack pattern of RC beam as a result of numerical simulations with various finite concrete element sizes: (a) [Forumla omitted. See PDF.] = 100 mm, (b) [Forumla omitted. See PDF.] = 50 mm, (c) [Forumla omitted. See PDF.] = 25 mm and (d) [Forumla omitted. See PDF.] = 12.5 mm.

Figure 10. Applied loading (P)–displacement (d) curves in the three-point bending test for the RC beam with different element size [Forumla omitted. See PDF.].

Figure 11. Loading–displacement curves of the simulation in three-point bending testing for the RC and plain concrete beam.

Figure 12. Loading—displacement curves of simulation and experiment in three-point bending testing (This figure has been reproduced from [34] with permission from American Concrete Institute).

Figure 13. Crack patterns as a result of (a) the numerical simulation and (b) experiment (This figure has been reproduced from [34] with permission from American Concrete Institute).

Figure 14. A numerical model for simulating the demolition of a small-scale RC structure: (a) numerical model of the RC structure and (b) embedded rebars in the RC structure.

Figure 15. Comparison of 3D-FDEM (left-side) and experiment (right-side) results from the demolition of a small-scale RC structure (This figure has been reproduced from [7] under the terms and conditions of the CC BY-NC license).

Figure 16. Comparison of the y-axis displacement in experiment and numerical simulation (This figure has been reproduced from [7] under the terms and conditions of the CC BY-NC license).

Figure 17. Comparison of the x-axis displacement in experiment and numerical simulation (This figure has been reproduced from [7] under the terms and conditions of the CC BY-NC license).

Table 1

Experimentally investigated physical and mechanical parameters of 35 MPa concrete in previous research [33].

Property	Unit	Value
Density, $\rho$	$\mathrm{kg}/{\mathrm{m}}^3$	2200
Young’s modulus, E	$\mathrm{GPa}$	30.0
Poisson’s ratio, v	-	0.20
Indirect tensile strength, $S_T$	$\mathrm{MPa}$	2.62
Cohesion, c	$\mathrm{MPa}$	7.87
Internal friction angle, $\varphi$	Degrees	35.0
Uniaxial compressive strength, $S_c$	$\mathrm{MPa}$	35.0

Table 2

Calibrated input parameters for the FDEM simulation of fracture behaviors in the concrete structures [33].

FDEM Parameter	Unit	Unit
		Concrete
Density ($\rho$)	$\mathrm{kg}/{\mathrm{m}}^3$	2200
Young’s modulus ($E$)	$\mathrm{GPa}$	30.0
Poisson’s ratio ($v$)	-	0.20
Microscopic tensile strength ($S_T$)	$\mathrm{MPa}$	2.10
Microscopic cohesion ($c$)	$\mathrm{MPa}$	6.70
Microscopic internal friction angle ($\varphi$)	Degrees	35.0
Microscopic mode I fracture energy ($G_{If}$)	$\mathrm{N}/\mathrm{m}$	30.0
Microscopic mode II fracture energy ($G_{IIf}$)	$\mathrm{N}/\mathrm{m}$	90.0
Artificial cohesive penalty $(P_{open}$, $P_{tan}$)	$\mathrm{GPa}$	3000
Artificial cohesive penalty ($P_{overlap}$)	$\mathrm{GPa}$	30,000

Table 3

Calibrating FDEM input parameters for modeling contact behavior in concrete structure [33].

FDEM Parameter	Unit	Value
Normal contact penalty between concrete—concrete interface ($P_{n\_con}^{CC}$)	$\mathrm{GPa}$	300.0
Normal contact penalty between concrete—platen interface ($P_{n\_con}^{CP}$)	$\mathrm{GPa}$	2000.0
Friction coefficient between concrete—concrete interface (${\mu }_{\mathrm{fric}}^{\mathrm{CC}}$)	-	0.6
Friction coefficient between concrete—platen interface (${\mu }_{\mathrm{fric}}^{\mathrm{CP}}$)	-	0.1

References

1. Xie, C.M.; Tao, T.J.; Chi, E.A.; Yang, J.; Zhao, M.S. Numerical Simulation of Frame Structure Blasting Demolition under Explosive Loading. Appl. Mech. Mater.; 2014; 553, pp. 164-169. [DOI: https://dx.doi.org/10.4028/www.scientific.net/AMM.533.164]

2. Choi, H.; Kim, H.-J.; Park, H.; Yoon, S.-J. Structural analysis of the pre-weakening of a cylindrical concrete silo for the application of overturning explosive demolition method. Explos. Blasting; 2009; 27, pp. 12-18.

3. Choi, J.-H.; Jung, J.-W.; Jeon, S.-W. Simulation of Blasting Demolition of Reinforced Concrete Structures and Ground Vibration using Finite Element Method. Tunn. Undergr. Sp.; 2009; 19, pp. 190-202.

4. Yang, G.-l.; Jiang, L.-l.; Che, Y.-l. Size Effect of Element in Structural Blasting Demolition. Proceedings of the 2011 2nd International Conference on Digital Manufacturing & Automation; Zhangjiajie, China, 5–7 August 2011; pp. 1094-1097.

5. Michaloudis, G.; Mattern, S.; Schweizerhof, K. Computer Simulation for Building Implosion Using LS-DYNA. High Performance Computing in Science and Engineering’10; Springer: Berlin/Heidelberg, Germany, 2011; pp. 519-528.

6. Meguro, K.; Tagel-Din, H. Applied element method for structural analysis theory and application for linear materials. Doboku Gakkai Ronbunshu; 2000; 2000, pp. 31-45. [DOI: https://dx.doi.org/10.2208/jscej.2000.647_31]

7. Park, H.; Suk, C.-G.; Kim, S.-K. Collapse modeling of model RC structure using applied element method. Tunn. Undergr. Space; 2009; 19, pp. 43-51.

8. Sediek, O.A.; El-Tawil, S.; McCormick, J. Seismic Debris Field for Collapsed RC Moment Resisting Frame Buildings. J. Struct. Eng.; 2021; 147, 04021045. [DOI: https://dx.doi.org/10.1061/(ASCE)ST.1943-541X.0002985]

9. Kim, H. Numerical Simulation of Demolition of the Shear Wall Apartment and Progressive Collapse-Local Damage Transition of RC Building Based on AEM; Jeonbuk National University: Jeonju, Korea, 2012.

10. Helmy, H.M.; Elfouly, A.K.; Salem, H.M. Numerical simulation of demolition of Perna Seca hospital using the applied element method. Proceedings of the Structures Congress 2012; Chicago, IL, USA, 29–31 March 2012; pp. 279-294.

11. Özmen, H.; Soyluk, K.; Anil, Ö. Analysis of RC structures with different design mistakes under explosive based demolition. Struct. Concr.; 2021; 22, pp. 1462-1486. [DOI: https://dx.doi.org/10.1002/suco.201900367]

12. Munjiza, A.A. The Combined Finite-Discrete Element Method; John Wiley & Sons: Hoboken, NJ, USA, 2004.

13. Zivaljic, N.; Smoljanovic, H.; Nikolic, Z. A combined finite-discrete element model for RC structures under dynamic loading. Eng. Comput.; 2013; 30, pp. 982-1010. [DOI: https://dx.doi.org/10.1108/EC-03-2012-0066]

14. Zivaljic, N.; Nikolic, Z.; Munjiza, A. A combined finite-discrete element model for reinforced concrete under seismic load. Proceedings of the XI International Conference on Computational Plasticity: Fundamentals and Applications (COMPLAS 2011); Barcelona, Spain, 7–9 September 2011; pp. 788-795.

15. Živaljić, N.; Nikolić, Ž.; Smoljanović, H. Computational aspects of the combined finite-discrete element method in modelling of plane reinforced concrete structures. Eng. Fract. Mech.; 2014; 131, pp. 669-686. [DOI: https://dx.doi.org/10.1016/j.engfracmech.2014.10.017]

16. Smoljanović, H.; Živaljić, N.; Nikolić, Ž.; Munjiza, A. Numerical model for confined masonry structures based on finite discrete element method. Int. J. Eng. Model.; 2017; 30, pp. 19-35.

17. Smoljanović, H.; Zivaljic, N.; Nikolić, Ž.; Munjiza, A. Numerical simulations of impacts using combined finite-discrete element method. Proceedings of the European Congress on Computational Methods in Applied Sciences and Engineering (ECCOMAS 2012); Viena, Austria, 10–14 September 2012; pp. 5019-5028.

18. Živaljić, N.; Smoljanović, H.; Nikolić, Ž. Sensitivity analysis of numerical parameters in FEM/DEM model for RC structures. Int. J. Eng. Model.; 2012; 25, pp. 7-17.

19. Nikolić, Ž.; Živaljić, N.; Smoljanović, H. The application of fem/dem method in analysis of rc beam-column joints. Proceedings of the 4th International Conference on Computational Methods in Structural Dynamics and Earthquake Engineering; Kos Island, Greece, 12–14 June 2013.

20. Tatone, B.S.A.; Lisjak, A.; Mahabadi, O.K.; Vlachopoulos, N. Verification of the Implementation of Rock-Reinforcement Elements in Numerical Analyses Based on the Hybrid Combined Finitediscrete Element Method (FDEM). Proceedings of the 49th US Rock Mechanics/Geomechanics Symposium; San Francisco, CA, USA, 28 June 2015.

21. Tatone, B.S.A.; Lisjak, A.; Mahabadi, O.K.; Vlachopoulos, N. Incorporating rock reinforcement elements into numerical analyses based on the hybrid finite-discrete element method (FDEM). Proceedings of the 13th ISRM International Congress of Rock Mechanics; Montreal, QC, Canada, 10–13 May 2015.

22. Soltani, M.; Maekawa, K. Path-dependent mechanical model for deformed reinforcing bars at RC interface under coupled cyclic shear and pullout tension. Eng. Struct.; 2008; 30, pp. 1079-1091. [DOI: https://dx.doi.org/10.1016/j.engstruct.2007.06.013]

23. Forbes, B.; Vlachopoulos, N.; Oke, J.; Hyett, A.J. The application of distributed optical sensing for monitoring temporary support schemes. Proceedings of the 67th Canadian Geotechnical Conference (GeoRegina 2014); Regina, SK, Canada, 28 September–1 October 2014.

24. Lisjak, A.; Young-Schultz, T.; Li, B.; He, L.; Tatone, B.S.A.; Mahabadi, O.K. A novel rockbolt formulation for a GPU-accelerated, finite-discrete element method code and its application to underground excavations. Int. J. Rock Mech. Min. Sci.; 2020; 134, 104410. [DOI: https://dx.doi.org/10.1016/j.ijrmms.2020.104410]

25. Munjiza, A.; Smoljanović, H.; Živaljić, N.; Mihanovic, A.; Divić, V.; Uzelac, I.; Nikolić, Ž.; Balić, I.; Trogrlić, B. Structural applications of the combined finite-discrete element method. Comput. Part. Mech.; 2020; 7, pp. 1029-1046. [DOI: https://dx.doi.org/10.1007/s40571-019-00286-5]

26. An, H.M.; Liu, H.Y.; Han, H.; Zheng, X.; Wang, X.G. Hybrid finite-discrete element modelling of dynamic fracture and resultant fragment casting and muck-piling by rock blast. Comput. Geotech.; 2017; 81, pp. 322-345. [DOI: https://dx.doi.org/10.1016/j.compgeo.2016.09.007]

27. Liu, H.Y.; Kang, Y.M.; Lin, P. Hybrid finite-discrete element modeling of geomaterials fracture and fragment muck-piling. Int. J. Geotech. Eng.; 2015; 9, pp. 115-131. [DOI: https://dx.doi.org/10.1179/1939787913Y.0000000035]

28. Liu, H.Y.; Han, H.; An, H.M.; Shi, J.J. Hybrid finite-discrete element modelling of asperity degradation and gouge grinding during direct shearing of rough rock joints. Int. J. Coal Sci. Technol.; 2016; 3, pp. 295-310. [DOI: https://dx.doi.org/10.1007/s40789-016-0142-1]

29. Fukuda, D.; Mohammadnejad, M.; Liu, H.; Zhang, Q.; Zhao, J.; Dehkhoda, S.; Chan, A.; Kodama, J.-i.; Fujii, Y. Development of a 3D hybrid finite-discrete element simulator based on GPGPU-parallelized computation for modelling rock fracturing under quasi-static and dynamic loading conditions. Rock Mech. Rock Eng.; 2020; 53, pp. 1079-1112. [DOI: https://dx.doi.org/10.1007/s00603-019-01960-z]

30. Fukuda, D.; Mohammadnejad, M.; Liu, H.; Dehkhoda, S.; Chan, A.; Cho, S.H.; Min, G.J.; Han, H.; Kodama, J.-i.; Fujii, Y. Development of a GPGPU-parallelized hybrid finite-discrete element method for modeling rock fracture. Int J. Numer. Anal. Met.; 2019; 43, pp. 1797-1824. [DOI: https://dx.doi.org/10.1002/nag.2934]

31. Munjiza, A.; Andrews, K.R.F.; White, J.K. Combined single and smeared crack model in combined finite-discrete element analysis. Int. J. Numer. Meth. Eng; 1999; 44, pp. 41-57. [DOI: https://dx.doi.org/10.1002/(SICI)1097-0207(19990110)44:1<41::AID-NME487>3.0.CO;2-A]

32. Mahabadi, O.K.; Lisjak, A.; Munjiza, A.; Grasselli, G. Y-Geo: New Combined Finite-Discrete Element Numerical Code for Geomechanical Applications. Int. J. Geomech.; 2012; 12, pp. 676-688. [DOI: https://dx.doi.org/10.1061/(ASCE)GM.1943-5622.0000216]

33. Min, G.; Fukuda, D.; Oh, S.; Kim, G.; Ko, Y.; Liu, H.; Chung, M.; Cho, S. Three-Dimensional Combined Finite-Discrete Element Modeling of Shear Fracture Process in Direct Shearing of Rough Concrete–Rock Joints. Appl. Sci.; 2020; 10, 8033. [DOI: https://dx.doi.org/10.3390/app10228033]

34. Carpinteri, A.; Carmona, J.R.; Ventura, G. Failure Mode Transitions in Reinforced Concrete Beams—Part 2: Experimental Tests. ACI Struct. J.; 2011; 108, pp. 286-293.

35. Mamlouk, M.S.; Zaniewski, J.P. Materials for Civil and Construction Engineers; Pearson: London, UK, 2014.