Discrete Element Code Research Articles

Quantitative Assessment of Bed-Separation Dynamic Development Caused by Inclined Coal Seam Longwall Mining

Coal mining under the Quaternary thick loose layer affects key strata breakage, Bed-separations development, ground subsidence, and other studies. This paper presents a method for solving the deflection of a large-deflection inclined thin plate under a thick loose-layer cover with additional lateral loads and midplane forces. The methods presented are based on the principle of large-deflection of thin-plate, energy method, and fracture mechanics theory. The 7225 work face in Anhui Province, China, was studied. Combined with the large-deflection inclined thin plate model, the initial breakage distance within the main roof plate was calculated to be 33 m with the initial breakage angle of 61.2°, and the period breakage distance was calculated to be 21 m with the period breakage angle of 55.4°. The distribution range of “Vertical Three Zones” from 7225 working face to the ground, including the height of the caved zone is 38.07 m, the height of the fractured zone is 41.13 m, and the height of the curved zone with the thick loose layer removed is 187.56 m. During the dynamic development of the principal key strata (PKS), the deflection value develops from 0 mm to 2714 mm with 7225 working face mining, and the maximum value of the spatial volume is 56,485 m3, which is verified by Three-dimensional Discrete Element Code (3DEC) numerical simulation. The dynamic development of Bed-separation within the overlying strata, with a maximum development height of 545.2 mm and a maximum volume of 11,228.1 m3 of the Bed-separation cavity. The dynamic development of the Bed-separation height and the cavity under different mining length and width conditions of the working face are also discussed. The large-deflection inclined thin plate model proposed in this paper effectively explores the dynamic deflection and fragmentation law of the overlying strata induced by the inclined working face of Longwall mining and provides a theoretical basis and computational model for quantitatively evaluating the dynamic development of the Bed-separation cavity.

Particle fracture regimes from impact simulations.

We introduce an approach to particle breakage, wherein the particle is modeled as an aggregate of polyhedral cells with their common surfaces governed by the Griffith criterion of fracture. This model is implemented within a discrete element code to simulate and analyze the breakage behavior of a single particle impacting a rigid plane. We find that fracture dynamics involves three distinct regimes as a function of the normalized impact energy ω. At low values of ω, the particle undergoes elastic rebound and no cracks occur inside the particle. In the intermediate range, the particle is damaged by nucleation and propagation of cracks, and the effective restitution coefficient declines without breakup of the particle. Finally, for values of ω beyond a well-defined threshold, the particle breaks into fragments and the restitution coefficient increases with ω due to kinetic energy carried away by the fragments. We show that particle damage, restitution coefficient, and fracture efficiency (the amount of energy input consumed for particle fracture) collapse well as a function of dimensionless scaling parameters. Our data are also sufficiently accurate to scale fragment size and shape distributions. It is found that fragment masses (volumes) follow a power-law distribution with an exponent decreasing with fracture energy. Interestingly, the average elongation and flatness of fragments are very close to those observed in experiments and lunar samples at the optimal fracture efficiency.

Quantification of the crack propagation and global failure mechanism of single- and multi-ring masonry arch bridges

This paper presents the development of a novel crack-based damage quantification method to assess the damage accumulation in masonry arch bridges subjected to externally applied loading. Initially, a damage index was proposed, considering the initiation and extension of tensile and sliding cracks. Then, the proposed damage indexing was adopted into a two-dimensional discrete element code of analysis to evaluate single and multi-ring masonry arch bridges. The numerical model was validated by comparing the computational outputs against the experimental results obtained from the literature. Following the validation, masonry arch bridge models were developed, featuring the four most representative bond configurations in arch barrels. Also, different bond properties were assigned to investigate their impact on the global failure characteristics of single- and multi-ring masonry arch bridges. From the analyses of results of the multi-ring masonry arch bridges investigated in this study, it was shown that the rapid growth of tensile cracks in radial joints and the occurrence of shear slippage in circumferential joints between adjacent rings were observed to be early signs for multi-ring bridges before the ultimate load reached. By analysing the correlation between the normalised damage index with the normalised load, damage accumulation in different masonry arch bridges were identified. This study offers a novel perspective on quantitively evaluating the damage accumulation and failure mechanism of masonry arch bridges and the findings may provide valuable insights into the assessment and management of existing masonry arch bridges.

Assessment of macro and micro mechanical properties of fresh and deteriorated ballast combining laboratory tests and 2D-discrete element methods

This research initially investigates the mechanical behavior of fresh and used ballast grain assemblies through a combination of laboratory experiments and discrete element modeling (DEM). The extended Los Angeles abrasion (LAA) test generated the used ballast grains defined at a fouling index (FI) of 40%, while fresh ballast grains serve as the reference material. Morphological analysis of angularity properties of fresh and used ballast grains was employed in producing 2D polygonal shapes implemented in a discrete element code. Numerical simulations using the contact dynamics method examined the bi-axial compression of the 2D-equivalent fresh and used ballast grain assemblies under varying confining pressures. Results show that both fresh and used ballast grain packings exhibit similar stress-strain curves, with fresh ballast packing displaying slightly higher shear strength (5%) at the residual state. Fresh ballast packing is denser than the used one at the peak shear state, while residual states show nearly identical volume fractions. Microscale analysis of contact and force networks reveals that the used ballast assembly is more anisotropic in contact at peak stress. At the same time, both packings exhibit similar contact anisotropy at the residual state. Fresh ballast packing exhibits higher force anisotropy, particularly in frictional forces, explaining the excess shear stress observed at the residual state. These findings suggest that used ballast grains mobilize local friction less efficiently than fresh ballast grains, potentially impacting the lateral strength properties of ballasted tracks. This study provides valuable insights into the behavior of used and fresh ballast grains, serving as a foundation for further investigations.

Investigation on the failure mechanism of a tunnel entrance high-steep slope with ultra-large wedge blocks based on discrete element method

An unstable wedge block with a volume of approximately 1,200,000 m3, intersected by an anti-dipped large fault and two outward-dipped faults, has been identified on the high-steep slope at the entrance of the tunnel spillway of the Shuangjiangkou hydropower station in China, a region characterized by high tectonic stresses. In this study, kinematic analysis was first conducted to preliminarily assess the possibility of block sliding, and then a three-dimensional discrete element code (3DEC) was utilized to analyze the mechanism of block instability during sequential excavation and reinforcement. The results showed that block sliding is caused by continuous excavation and deterioration of the mechanical properties of the bottom sliding surface. The failure behavior of the slope is primarily governed by the wedge block, which slides along the intersection line of two outward-dipped faults. After the excavation of Layers 1–4, the shear displacements of the faults exhibited accelerated and abrupt characteristics. Finally, the pretension cable proved effective in restraining block sliding, increasing the factor of safety for the slope from 1.16 to 1.28 after complete excavation and reinforcement. This study offers significant insights into the analysis of failure mechanisms in similar engineering challenges.

A macroscopic donor–acceptor-based discrete element model for contact electrification of insulating granular materials

Accurate modeling of contact electrification (CE) in insulating materials are essential for predicting and controlling their behavior. Here we propose a mesoscopic model based on the donor–acceptor mechanism, incorporating findings on surface patch features of identical dielectric materials from CE experiments. This model enables us to reproduce the charge transfer process and develop a macroscopic model that statistically captures this process. We implemented the macroscopic model in the discrete element code LIGGGHTS and validated its performance through comparisons with granular-vibration CE experiments. The results demonstrate the reliable performance of this model, as evidenced by the close agreement in overall particle charge distributions with experimental data. Furthermore, the model consistently reproduces the size-dependent bipolar charging feature. Our numerical tests also highlight the importance of stochastic discharging after breakdown and the influence of boundary conditions, which can reduce the mean charges and increase the variances, providing a better fit with experimental observations.

Transferring vision-based data to discontinuum analysis for the assessment of URM walls

This research presents a new workflow to assess unreinforced masonry (URM) walls based on vision-based data. The morphological features (i.e. brick pattern) of URM walls are obtained through high-resolution orthogonal images and then imported into a three-dimensional discrete element code (3DEC) to perform structural analysis via a semi-automatic procedure. URM walls are represented using rigid blocks that can mechanically interact along their boundaries. The mechanical interaction is simulated based on the point-contact hypothesis considering spring frictional elements. Once the numerical approach is validated using the results available in the literature, the proposed data-driven modelling strategy is applied to different URM wall cross-sections, with or without openings. The results show that the adopted semi-automatic workflow can be a fast and robust solution to simulate the structural behaviour of existing URM buildings by considering their construction techniques efficiently.

Effect of chemical corrosion on propagation of complex fracture networks under different hydraulic pressures in enhanced geothermal systems

Chemical corrosion and hydraulic fracturing are important ways to improve the permeability of the enhanced geothermal systems. In this work, we developed a discrete element code for simulating the hydraulic fracturing of three different fracture network models after chemical corrosion. Mechanisms of the crack extension and coalescence in the field-scale models with fracture network were systematically revealed, in which the effects of corrosion time, inclination of the fracture (γ), pH value and hydraulic pressure were considered. The results show that the simulated cracking behaviors of granites agree well with those obtained through Brazilian tests and triaxial compression tests, verifying the validity of the developed discrete element code. During hydraulic fracturing, the impact of γ on the fracture network connectivity increases with the increment of hydraulic pressure (Pin). With the increase in hydraulic pressure, crack propagation accelerates whereas the connectivity of fracture network models decreases. The reasons are that the excessive hydraulic pressure causes the crack to spread in an unfavorable direction and increases the effect of γ on the connectivity of models of fracture networks. As the pH value increases from 2 to 12, the crack growth rate first decreases and subsequently raises. The numbers of cracks per time step of sparse mixed fracture networks with artificial fractures (SM), dense mixed fracture networks with artificial fractures (DM) and dense mixed fracture networks without artificial fractures (NDM) are in the ranges of 0.06601–0.1259, 0.05176–0.2072 and 0.1727–0.6846, respectively. The connectivity of SM is better than that of NDM under the same hydraulic fracturing conditions, indicating that artificial fractures can greatly improve the fracture network connectivity and the utilization efficiency of EGS.

Effect of the Number and Distribution of Circular Holes on the Damage and Fracture of Limestone Using PFC2D

The common occurrence of various number and distribution of hole flaws complicates the mechanical behavior and fracture mode of the rock masses that contain them. This study develops seven numerical models of limestone samples with different numbers and distributions of circular hole flaws using 2-D particle flow discrete element code (PFC2D) to investigate their impact on the mechanical properties of limestone while maintaining the same flaws area. In addition, it analyzes the effect of these factors on the mesomechanical features of each model, including the characteristic stress values (peak stress, crack initiation, and damage stress), crack evolution, stress, and displacement field. The results showed that the peak stress, crack initiation, and damage stress of the single-hole model are between those of multihole models. As the arrangement of dip angle increases, the peak stress, crack initiation, and damage stress of models with the same number of multihole flaws exhibit a V-shaped change. The characteristic stress values are the largest when the holes are vertically aligned. Model differences in crack development path, shape, and number, as well as stress concentration area and failure mode, are primarily due to the number and distribution of holes. The circular holes are arranged at approximately 45°, and the greater the number of defects, the more likely the model is to fail. The study’s findings can provide support and reference for the research system of deformation and failure of rock mass with hole flaws.

Numerical simulation of the effect of internal hole defect size on the mechanical properties of limestone.

To better understand the effect of the size of hole defects on the mechanical properties of a rock mass, the two-dimensional particle flow discrete element code (PFC2D) is applied to establish rock mass models with single circular hole defects of different diameters. Uniaxial compressive strength (UCS) tests are conducted on each model by only taking the defect size (area) as a variable. This study analyzes each model’s stress-strain, contact force chain, crack evolution, meso-damage and failure, and mechanical properties. The results showed that with the size enlargement of the circular hole defects, each model’s UCS and elastic modulus gradually decrease, and the defect size is negatively correlated with the mechanical strength of the rock samples. The size of the hole defects affects the entire process of contact force chain and crack evolution. The larger the aperture dimension of the circular hole defects in each model, the greater the concentration degree of the contact force chain, the earlier the crack initiation, and the higher the degree of crack coalescence in the post-peak stage. The number of cracks decreases as the hole size increases, and the model is more prone to failure. Rock models’ strength and failure characteristics with different numbers and arrangements of hole defects are discussed under the same defect area condition.

Probabilistic approach to assess URM walls with openings using discrete rigid block analysis (D-RBA)

This study aims to improve our current understanding of the seismic assessment of load-bearing unreinforced masonry (URM) systems by proposing a probabilistic computational modeling framework using the discrete element method (DEM). The main objective is to predict the structural behavior and capacity of URM walls with openings subjected to lateral loading, considering uncertainties in material properties. The proposed modeling strategy represents masonry as an assembly of rigid blocks interacting along their boundaries by adopting the point-contact hypothesis. Fracture energy-based softening contact models are implemented into a commercial discrete element code (3DEC) to better simulate both the pre- and post-peak behavior of masonry. The results highlight the influence of material properties on the force capacity, displacement capacity (drift limits), and collapse mechanisms of walls with openings. Based on the applied non-spatial probabilistic analyses, the most commonly observed failure mechanisms are further assessed using a simplified macro-block formulation. As a result, practical, yet necessary, inferences are made, providing valuable contributions. Furthermore, the validated discontinuum analysis framework is demonstrated as an accurate structural analysis strategy and a useful approach to simulating the potential collapse mechanism of load-bearing URM structures.

Phase transition in monodisperse granular materials: How to model it by using a strain hardening visco‐elastic‐plastic constitutive relationship

AbstractDuring heating (loading characterised by a progressive increase in strain rate) and cooling (loading characterised by a progressive decrease in strain rate) numerical tests, performed by using Discrete Element codes, granular materials experience phase transition phenomena, named in this paper ‘dry liquefaction’ (when from solid the material starts behaving like a fluid) and ‘dry resolidification’ (freezing, when the material from fluid starts behaving like a solid). The aim of this paper consists in reproducing phase transition phenomena by using a strain hardening visco‐elastic‐plastic model based on the critical state concept and kinetic theories of granular gases. The authors demonstrate that crucial is the role of isotropic softening/hardening, which describes the size of the elastic domain and the capability of the solid skeleton of storing elastic energy according to permanent force chains. The main ingredients of the model are: (i) the additivity of quasi‐static and collisional stresses, (ii) the energy balance equation governing the evolution of the granular temperature, interpreted this latter as an additional internal variable for the system for the collisional contribution, (iii) the mixed isotropic and kinematic hardening characterising the quasi‐static incremental constitutive relationship. The model has been both calibrated and validated on Discrete Element Method (DEM) numerical results, obtained by testing dry assemblies of monodisperse spheres under true triaxial loading conditions.

An experiment-based cohesive-frictional constitutive model for cemented materials

Experiment-based damage laws at the contact bond level are often absent in the computational modelling of cemented materials using cohesive contact-bond models. In this paper, in-situ X-Ray Computed Tomography on element-scale cement bonds subjected to pure compression, shear and tensile loading was conducted to investigate their failure mechanisms and associated damage laws. A new cohesive-frictional constitutive model was subsequently developed using results obtained from X-Ray CT tests. The new constitutive model was capable of describing the behaviour of contact bonds subjected to complex loading conditions, including mixed-mode and compression damage. The new model was then implemented in an open-access Discrete Element Code (YADE) to simulate several challenging laboratory-scale tests (e.g., triaxial, Brazilian and Splitting tests) and very good agreements with the experimental data were achieved.

A Vertical Joint Spacing Calculation Method for UDEC Modeling of Large-Scale Strata and Its Influence on Mining-Induced Surface Subsidence

While universal discrete element code (UDEC) is widely used for understanding the mechanism of large-scale strata movement and the effects of mining subsidence on the environment, the fundamental knowledge of how to set vertical joint spacing (VJS) in UDEC is still not fully understood. To address the knowledge gap, we first present a novel VJS calculation method, then conduct UDEC experiments, and finally compare the predictions of UDEC models with field subsidence observation. The results suggest the following: (1) when compared to the conventional VJS setting (1× to 3× bed thickness), the maximum surface subsidence (MSS) prediction via UDEC models based on our proposed VJS setting method is closer to field observation; (2) a smaller but varying VJS setting can also have the effect of a larger VJS setting; and (3) with the increase in VJS, MSS first drops, then rises, and reaches the minimum when VJS is set at approximately 7× bed thickness. This paper provides an explanation of the VJS setting in UDEC and establishes a bridge between the KS theory and VJS, which is helpful for the sustainable development of such an UDEC modeling strategy and for a better understanding of the influences of mining subsidence on the environment in mining-affected areas.

On the role of fibre bonds on the elasticity of low-density papers: a micro-mechanical approach

Fine prediction of the elastic properties of paper materials can now be obtained using sophisticated fibre scale numerical approaches. However, there is still a need, in particular for low-density papers, for simple and compact analytical models that enable the elastic properties of these papers to be estimated from the knowledge of various structural information about their fibres and their fibrous networks. For that purpose, we pursued the analysis carried out in Marulier et al. (Cellulose 22:1517–1539, 2015. https://doi.org/10.1007/s10570-015-0610-6) with low-density papers that were fabricated with planar random and orientated fibrous microstructures and different fibre contents. The fibrous microstructures of these papers were imaged using X-ray synchrotron microtomography. The corresponding 3D images revealed highly connected fibrous networks with small fibre bond areas. Furthermore, the evolutions of their Young’s moduli were non-linear and evolved as power-laws with the fibre content. Current analytical models of the literature do not capture these trends. In light of these experimental data, we developed a fibre network model for the in-plane elasticity of papers in which the main deformation mechanisms of the micromechanical model is the shear of the numerous fibre bonds and their vicinity, whereas the fibre parts far from these zones were considered as rigid bodies. The stiffness tensor of papers was then estimated both numerically using a discrete element code and analytically using additional assumptions. Both approaches nicely fit the experimental trends by adjusting a unique unknown micromechanical parameter, which is the shear stiffness of bonding zones. The estimate of this parameter is relevant in light of several recently reported experimental results.

A coupled experimental and numerical simulation of concretejoints' behaviors in tunnel support using concrete specimens

The experimentally tested modelled specimens were simulated by a two-dimensional particle flow code to study the behavior of rock mass surrounding a tunnel interacted with a nearby rock joint or discontinuity. The specially prepared specimens are tested in the laboratory and the measured results are provided. Then a numerical modelling of these tests is accomplished by a calibrated two-dimensional particle flow code to study the rock tunnel behavior while interacting with a neighboring joint. The two-dimensional discrete element code was calibrated using Brazilian tensile test. Then the modelled specimens are provided so that various configurations of tunnel cross sections and the neighboring joints were tested under uniaxial compression condition. This study showed that the tensile cracks are the most dominant mode occurred in these modelled samples. The wing cracks initiated from the joint tips when the joints interacted in a position less than that of the tunnel height. These cracks are then propagated and interacted with the tunnel ceiling. When the joint interacting the tunnel head and the interaction angle is negative the tunnel can be in a stable position considering the joints effect on its instability situation. But for positive interaction angles and the case of joint existing near the tunnel head, the wing cracks may initiate and propagate towards the tunnel ceiling. As the distance of joint from the tunnel ceiling is increased its effect on the tunnel instability is decreased because the failure stress is increased. The number of joints and their distance with the tunnel boundary (ceiling) have also a profound effect on the stability condition of the tunnel. The failure stress reached its maximum value for the increase from -30 to -60 degrees or increase from 30 to 60 degrees. The failure stress also decreased as the number of notches and their lengths increased. In all these interaction scenarios, the corresponding numerical and experimental values compared and it is concluded that the failure stresses are very close to each other which verified the accuracy and applicability of the proposed modeling technic.

Behaviors and Overlying Strata Failure Law for Underground Filling of a Gently Inclined Medium-Thick Phosphate Deposit

In this study, the No. 6 pit in the eastern mining area of the Jinning phosphate mine in China was taken as the research background. In order to reduce the cost of filling, an improved pillarless sublevel caving method is proposed. This method greatly improves the ore recovery rate by adding a recovery route. In addition, the combination of similar material simulation experiments and numerical simulations (discrete element and universal distinct element code) revealed the deformation and failure laws of the surrounding roof rock and the characteristics of the surface subsidence. The results indicate the following. (1) The similar simulation experimental results indicate that the deformation of the overlying rock layer originated from the direct roof of the goaf and gradually developed into the deep part of the rock layer. An irregular stepped caving zone formed in the goaf. The maximum surface subsidence was located above the phosphorus orebody, and it gradually decreased toward both sides. As the stope approached propulsion, the location of the maximum subsidence gradually moved toward the propulsion direction. (2) The numerical results revealed that the displacement of the overlying strata was nonlinear, and it decreased with increasing roof height. A support pressure concentration area was formed within a certain range of the stope roof. The numerical simulation results are basically consistent with the similar simulation experimental results.

Seismic vulnerability assessment of masonry arch bridges

A seismic assessment study of two historical masonry arch bridges, with different geometries and spans, is presented. The seismicity of the two locations was analyzed in order to define the seismic action to be applied. Artificial ground motion records were generated accordingly, and real ground motions were also used. The two numerical models were created with a 3D discrete element code. The accuracy of the deformability of the block models was verified by comparison with a finite element model. In situ measurements of the natural frequencies were available for one of the bridges, providing further means for model calibration. Time domain dynamic analyses were performed with each set of ground motions, allowing an evaluation of damage and failure modes under seismic action. The shorter span bridge was found to withstand the regulatory records with minimum damage. In contrast, the longer span bridge displayed a significant vulnerability to earthquake action. The study shows that rigid block models have the ability to represent effectively the dynamic elastic response of the bridges and can be calibrated by ambient vibration measurements. Furthermore, the discrete element models provide a powerful tool to investigate the nonlinear dynamic response and failure of these structures, which has to be based on a comprehensive set of dynamic records.

Numerical study and field performance of rockbolt support schemes in TBM-excavated coal mine roadways: A case study

There is an increasing number of TBMs applied in coal mine roadway excavation projects due to the higher safety, faster tunnelling speed and fewer labor works. Nevertheless, the experiences gained from TBM tunnelling projects of transportation or water conveyance tunnels are not well applicable in underground coal mines because of their special engineering and geological conditions. This paper proposed a case study of the excavation disturbance behaviors and support design of TBM-excavated roadway in Zhangji coal mine in Anhui Province, China. Mechanical properties of roadway surrounding rocks were obtained based on parameters of intact rocks and Rock Quality Designation (RQD) values of rock masses. The Universal Discrete Element Code (UDEC) software was used to establish numerical model of the roadway and micro-parameters are also calibrated. The excavation disturbance behaviors of unsupported roadway were firstly simulated. Moreover, the excavation disturbance behaviors and rockbolts loading characters under two different support schemes were studied and compared based on simulation and monitoring results. In this case, simulation results are in good agreement with monitoring results. Research results and experiences obtained from the working site suggest that installing rockbolts perpendicularly to rock surface can better satisfy the requirements of TBM tunnelling in weak grounds of coal mine roadways.

In-plane structural performance of dry-joint stone masonry Walls: A spatial and non-spatial stochastic discontinuum analysis

In this study, the in-plane structural behavior, capacity, and performance of dry-joint stone masonry walls (DJ-SMWs) and the effects of the vertical stress level on these factors are investigated via a stochastic discontinuum analysis that considers the material uncertainty. A discontinuum type of analysis is performed based on the discrete element method (DEM), where each stone masonry unit is explicitly represented in the computational model. To better simulate the cracking and shear failure modes within the stone units, a coupled fracture energy-based contact constitutive model is implemented into a commercial discrete element code, 3DEC. First, the proposed modeling approach is validated by comparing to experimental findings in literature. Then, the approach is used to explore the failure mechanism and the force–displacement behavior of DJ-SMWs, considering different vertical stress levels and material properties. The results of the novel modeling strategy provide a better understanding of the progressive collapse mechanism of DJ-SMWs and the influence of the vertical stress level. Furthermore, the outcomes of this research indicate the major role of the frictional resistance at the joints in the safety and performance assessment of the dry-joint load-bearing masonry walls. Finally, important inferences are made regarding the non-spatial and spatial stochastic discontinuum analysis.