Hybrid Finite-discrete Element Method Research Articles

A finite discrete element approach for modeling of desiccation fracturing around underground openings in Opalinus clay

Understanding and predicting potential failure mechanisms during the excavation and open drift stages of geological repository construction are among the crucial aspects of performance evaluation and safety assessment of nuclear waste storage facilities. The development of the Excavation Damage Zone (EDZ) and the generation of shrinkage-induced cracks during operational phases are prominent examples of failure mechanisms that can compromise the integrity of the repository systems. This study presents an integrated framework for investigating shrinkage-induced cracking of Opalinus Clay in niches and tunnels. To achieve this, the hybrid Finite Discrete Element Method (FDEM) is employed. The methodology incorporates a two-way staggered hydro-mechanical coupling scheme, where solid phase analysis relies on 2D FDEM and fluid flow is modeled using the nonlinear Richards’ equation and solved via Finite Volume discretization. To account for the effects of EDZ, characterized by a pronounced increase in hydraulic conductivity, a numerical simulation of tunnel excavation is first carried out. The resulting failure pattern around underground openings is then abstracted through the definition of an altered hydraulic conductivity field. Comparison of the numerical results with field observations demonstrates the framework’s ability to capture a wide range of failure mechanisms inherent in various stages of underground repository construction in Opalinus Clay.

A hybrid finite-discrete element method for modelling cracking processes in sandy mudstone containing a single edge-flaw under cyclic dynamic loading

Rock mass deformation and failure are macroscopic manifestations of crack initiation, propagation, and coalescence. However, simulating the transition of rocks from continuous to discontinuous media under cyclic dynamic loading remains challenging. This study proposes a hybrid finite-discrete element method (HFDEM) to model crack propagation, incorporating a frequency-dependent cohesive-zone model. The mechanical properties of standard sandy mudstone under quasi-static and cyclic dynamic loading were simulated using HFDEM, and the method's reliability was verified through experimental comparison. The comparative analysis demonstrates that HFDEM successfully captures crack interaction mechanisms and accurately simulates the overall failure behavior of specimens. Additionally, the effects of pre-existing flaw inclination angle and dynamic loading frequency on rock failure mechanisms were investigated. The numerical results reveal that rock samples exhibit significantly higher compressive strength under dynamic loading compared to quasi-static loading, with compressive strength increasing with higher cyclic dynamic load frequencies. Furthermore, by analyzing the strength characteristics, crack propagation, and failure modes of the samples, insights into the failure mechanisms of rocks under different frequency loads were obtained. This study provides valuable insights into crack development and failure of rocks under seismic loads, offering guidance for engineering practices.

Experimental and numerical analysis of NdFeB magnets under static and dynamic loading conditions

Neodymium iron boron (NdFeB) magnets, known for their superior magnetic properties, are pivotal in various high-performance applications, ranging from electric motors to wind turbines. The durability and reliability of these devices significantly depend on the mechanical integrity of NdFeB magnets under various stress conditions. In this work, the mechanical behavior and failure processes of NdFeB magnets are experimentally studied by in-door tests. We specifically examine the failure behavior of these magnets when subjected to dynamic loads, employing the Split Hopkinson Pressure Bar test both with and without initial loading perpendicular to the pressure bar direction. Furthermore, to complement the experimental findings, the experiments are numerically reproduced with a hybrid finite-discrete element method titled Continuous-Discontinuous Element Method (CDEM). Within the CDEM framework the Johnson-Holmquist 2 constitutive model is considered for the material under high strain rate loading conditions. Most material parameters are determined by experimental data and theoretical analysis. This research not only enhances our understanding of the properties of NdFeB magnets but also provides valuable insights for engineers and researchers aiming to improve the material’s resilience and durability.

Advanced tensile fracture analysis of alumina ceramics: Integrating hybrid finite-discrete element modeling with experimental insights

This study offers novel insights into the Mode I tensile response of an alumina ceramic through the use of computational modeling and the flattened Brazilian disk (FBD) experiments. A modified hybrid finite-discrete element method (HFDEM) is developed, integrating a coupled damage and friction cohesive model and a microscopic stochastic fracture model with a Weibull strength distribution by Monte Carlo simulation. The model is used to simulate direct tensile failure processes under quasi-static loading conditions, providing qualitative and quantitative predictions of direct tensile failure processes of an alumina ceramic. Concurrently, quasi-static flattened Brazilian disk tests (indirect tensile tests) are performed on a standard MTS machine coupled with a high-speed camera. The modified HFDEM model is also applied to reproduce the FBD experiments, and our simulated tensile strength is consistent with the experimental results. The results of the modified HFDEM model show three kinds of phenomena (i.e., “underestimation”, “reasonable estimation”, and “overestimation” of the indirect tensile strength) and four different types of associated fracture and fragment patterns of FBD testing. The integration of simulation and experimental results reveal relationships between fracture patterns, fragment geometry, tensile strength, and indirect tensile strength. The fracture and fragmentation patterns derived from our modified HFDEM model can be utilized to analyze the “tensile strength” measured in BD testing. Overall, this research offers important insights and direction for future Brazilian disk experiments and tensile strength assessments, enhancing our comprehension of ceramic failure mechanisms under tensile loadings.

Application of a Finite-Discrete Element Method Code for Modelling Rock Spalling in Tunnels: The Case of the Lyon-Turin Base Tunnel

Brittle failure, or spalling, occurs around openings excavated in hard rock masses with high in situ stresses. It takes place due to the nucleation and growth of cracks around the excavation boundary, induced by the redistribution of stresses following the excavation. Modelling this failure process is a tough challenge. The hybrid finite-discrete element method (FDEM) can overcome the boundary between continuum and discontinuum, capturing emergent discontinuities associated with brittle fracturing processes. In this study, FDEM is applied using a commercial code to show its applicability to model brittle behaviour around deep underground excavations in the case of the Torino-Lyon Base tunnel in three different stress conditions. Except for the hydrostatic condition, cracking is triggered immediately after the excavation. Spalling occurring around the tunnel is quite extended; therefore, an accurately designed support must be installed to prevent blocks from falling from the tunnel boundary. The obtained results are aligned with previous results existing in the literature. However, in this case, a deeper spalling is caused by the shape change due to the gradual stress redistribution. Such a phenomenon underlines the importance of using a code able to identify crack propagation, opening, and the formation of loose blocks that progressively modify the tunnel contour.

Empirical scaling of formation fracturing by high-energy impulsive mechanical loads

Recent technological advances to trigger high-energy seismic waves from within the wellbore have spurred interest in their potential application to induce fracturing within deep subsurface formations. Although a considerable body of experimental observations at the bench scale (on the order of 1 cubic foot) show promise, there remains considerable uncertainty in how the process may scale to the field. This work characterizes the relationship between the extent and intensity of fracturing and the duration and strength of the dynamic load. Our approach uses a direct numerical simulation of the elastodynamic equations while accounting for nonlinear fracture mechanics. A hybrid Finite-Discrete Element Method (FDEM) is applied whereby cohesive (elastoplastic) laws hold mesh elements together until complete failure. Beyond failure, elements act as deformable free bodies that can interact via contact constraints, including friction. The consistency of the model is first verified and it is validated with observations from previously reported physical experiments. Subsequently, we apply the model to conduct a numerical exploration of the dimensionless parameter space that defines plausible loading waveforms. We model an infinite domain that contains a spherical inclusion within which various impulsive loads are applied. The dynamic loading waveform is parameterized by its rise time to peak stress, followed by a decay period, and occurring within micro- to milliseconds. The dimensionless parameter space is sampled by varying loading characteristics (rise time, peak pressure, and impulse) to reveal the stimulated damaged bulk volume and crack intensity. Consistent with reported experimental observations, the model reproduces a nearly linear trend between the radius of damage and the maximum stress on the bench scale. Beyond that scale, however, the model predicts that this scaling slows considerably to a fractional power law between the damaged radius and the peak stress. This limitation coincides with a geometric increase in the intensity of damage within the stimulated volume. We present loci curves of the expected damage radii as functions of the loading waveform properties.

Impact of stress-driven crack growth on the emergence of anomalous transport in critically connected natural fracture networks

We study the role of in-situ stresses in fluid flow and mass transport through a two-dimensional fractured rock. The fracture network is based on a real outcrop with a connectivity state around the percolation threshold. We simulate stress-dependent fracture propagation and deformation using a geomechanics model based on the hybrid finite-discrete element method. The hydrodynamic transport through the deformed rock mass is then modeled based on particle tracking simulations. By exploring various scenarios of different stress ratios and orientations, we systematically analyzed how the superimposed multiple geomechanical processes (i.e., compression-induced closure, shear-induced dilation, and stress-induced propagation) control the aperture patterns, flow fields, and transport behaviors of the fracture network. Our results show that an increased stress ratio tends to promote fracture shear dilation and thus enhance flow channelization, resulting in an anomalously early arrival of solute transport. Furthermore, the increased differential stress attempts to drive the growth of wing cracks that bridge pre-existing fractures and form new flow pathways. We observe a qualitative transformation of the flow structure from a disconnected to a connected regime accompanied by a significant alteration of its overall transport properties. This phenomenon is attributed to the superimposed effect of stress-driven crack propagation and shear-induced fracture dilation. By quantitative parameter analyses of global flow and transport properties, we further found that the breakthrough time median of particles is well correlated with the equivalent permeability of the fractured rock when the effect of new cracks is removed. However, this correlation is greatly reduced when fracture propagation is considered; instead, the breakthrough time shows a good correlation with mean tortuosity. Our findings have important implications for the strategy design of geothermal production and waste disposal in critically connected fractured geological media under various stress conditions.

Vulnerability of Pointed Masonry Barrel Vaults Subjected to Differential Settlement Simulated with a GPGPU-Parallelized FDEM

Pointed masonry barrel vaults are widely used in classical historic structures, such as cathedrals and aqueducts, and they are very sensitive to differential settlement. These vaults are assemblages of masonry units and mortar. Since the bonding strength of mortar degrades over ages, dry-joint assumption is widely accepted. Failure behavior of dry-joint pointed masonry barrel vaults subjected to differential settlement is highly complex, discontinuous, and nonlinear. In this study, a 3D GPGPU-parallelized hybrid finite-discrete element method (FDEM), which is an advanced extension of finite element method (FEM) and discrete element method (DEM), is employed to investigate the capacity of pointed masonry barrel vaults subjected to differential settlement. When modeling barrel vaults with 3D FDEM, each masonry unit is discretized into a couple of four-node tetrahedral elements whose deformability is characterized by standard finite element formulation. Thus, structural deformation and interaction forces can be obtained in an accurate manner. Numerical examples are presented and validated with results from literatures. A base case is selected, and the influence of embrace angle ([Formula: see text], sharpness (Sh), stockiness (St), and out-of-plane length ([Formula: see text] on the failure behavior is parametrically investigated. The larger the [Formula: see text] or Sh, the smaller the ultimate settlement. The same applies to St in general, while an excessively large St results in small ultimate settlement due to sliding. The influence of [Formula: see text] can be mitigated should it is large enough compared with the span. It is demonstrated that the 3D GPGPU-parallelized FDEM is a robust tool for analyzing the vulnerability of pointed masonry barrel vaults subjected to differential settlement.

Numerical and Experimental Simulation of Hydraulic Fracture Propagation Mechanism in Conglomerate Formation Based on Hybrid Finite-Discrete Element Method

Hydraulic fracturing was the main technology to achieve the economic development of conglomerate reservoirs, knowing that the hydraulic fracture propagation mode was of great significance for improving the development of conglomerate reservoirs. This paper proposed a new method to understand the hydraulic fracture behavior based on a hybrid finite-discrete element method. The simulation indicated that a complex fracture network was created near the wellbore in the studied conglomerate reservoir, and hydraulic fracture propagation around the gravel layer was the main failure mode when the hydraulic fracture reached the gravel layer. From the simulations, it was shown that under small differences in horizontal stress and tensile strength, the hydraulic fracture propagated more easily around the gravel layer, while it could cross the gravel under large differences in horizontal stress and tensile strength. Greater tensile strength differences can reduce the complexity of the fracture network. In addition, higher pumping rates and viscosities of fracturing fluid contribute to the complex fracture network and also can produce more gravel crosses when the hydraulic fracture is met. The main reason was that a higher pumping rate and higher viscosity of fracturing fluid can obtain a higher net pressure, which can ensure the hydraulic fracture crosses the gravel layer.

Numerical investigation of progressive damage and associated seismicity on a laboratory fault

Understanding rock shear failure behavior is crucial to gain insights into slip-related geohazards such as rock avalanches, landslides, and earthquakes. However, descriptions of the progressive damage on the shear surface are still incomplete or ambiguous. In this study, we use the hybrid finite-discrete element method (FDEM) to simulate a shear experiment and obtain a detailed comprehension of shear induced progressive damage and the associated seismic activity. We built a laboratory fault model from high resolution surface scans and micro-CT imaging. Our results show that under quasi-static shear loading, the fault surface experiences local dynamic seismic activities. We found that the seismic activity is related to the stress concentration on interlocking asperities. This interlocking behavior (i) causes stress concentration at the region of contact that could reach the compressive strength, and (ii) produces tensile stress up to the tensile strength in the region adjacent to the contact area. Thus, different failure mechanisms and damage patterns including crushing and sub-vertical fracturing are observed on the rough surface. Asperity failure creates rapid local slips resulting in significant stress perturbations that alter the overall stress condition and may trigger the slip of adjacent critically stressed asperities. We found that the spatial distribution of the damaged asperities and the seismic activity is highly heterogeneous; regions with intense asperity interactions formed gouge material, while others exhibit minimal to no damage. These results emphasize the important role of surface roughness in controlling the overall shear behavior and the local dynamic seismic activities on faults.

Fracture Process and Failure Mode of Brazilian Discs with Cracks of Different Angles: A Numerical Study

In order to determine the effect of internal cracks on the tensile failure of materials, a hybrid finite–discrete element method was used to analyze the Brazilian disc test with cracks of different angles. When the pre-crack angle is between 0° and 60°, the wing crack is initiated from the pre-crack end. When the pre-crack is 90°, the crack initiated from the pre-crack center. When the pre-crack angle is between 0° and 60°, the maximum principal stress and plastic strain are concentrated at the pre-crack end. When the pre-crack angle is 90°, the maximum principal stress and plastic strain are concentrated in the pre-crack center. As the crack angle increased from 0° to 90°, the failure mode of Brazilian discs with cracks transits from splitting into two parts to splitting into four parts. The influence of crack length is further studied. When the crack length is less than 5 mm, the crack angle has little influence on the disc failure mode; Brazilian discs with cracks of different angles undergoes splitting failure along the loading axis. When the crack length is larger than 5 mm, the crack angle has a great effect on the disc failure mode.

Three-dimensional modelling on the impact fracture of glass using a GPGPU-parallelised FDEM

Due to the brittleness and the wide use of glass in modern engineering applications, its vulnerability to impact actions and the corresponding fracture behaviour attracted growing attentions from academics and engineers. In this study, impact fracture responses of glass have been modelled and simulated using a 3D GPGPU-parallelised hybrid finite-discrete element method, i.e., the FDEM. Glass is discretised into discrete elements where finite element formulation is incorporated, enabling accurate predictions on contact forces and structural deformation. A cohesive fracture model accounting for the rupture of glass is implemented, and numerical examples are presented and validated with results from literatures. The influence of impact velocity, boundary condition and projectile nose shape on the fracture of glass has been investigated. It is found that: (i) fracture pattern changes with the change of velocity; (ii) a rigid boundary support can be used should no damage occur in the edge of glass; (iii) under the same circumstance, a larger contact surface results in more severe damage. The GPGPU-parallelised FDEM provides a practical, efficient and robust computational approach in analysing the impact transient dynamic behaviour of glass in 3D.

Cracking behavior of concrete/rock bi-material specimens containing a parallel flaw pair under compression

Crack propagation and coalescence of pre-existing flaws in concrete/rock bi-material are critical for stability assessments of rock/concrete structures. A hybrid finite-discrete element method (FDEM) was used to simulate the cracking behavior of the concrete/rock bi-material containing a parallel flaw pair under compression. The model was validated by physical tests. The results indicate that the crack evolution does not occur synchronously in the concrete and rock layers. The cracks in the different material layers penetrate the interface and coalesce at lower interface dip angles, whereas an interface crack always forms at higher interface dip angles, preventing crack coalescence in different material layers. A threshold of the strength ratio exists and determines whether the cracks in different material layers coalesce. If the strength ratio is lower than the threshold, the crack propagations are limited to one material due to the interface crack. If the strength ratio exceeds the threshold, the cracks in the different materials penetrate the interface and coalesce. Two critical values of the interface reduction coefficient exist for a given interface dip angle. The cracking pattern of the bi-material can be classified into cracking along the interface, local cracking in the different materials, and crack coalescence in different materials.

Development of a three-dimensional grain-based combined finite-discrete element method to model the failure process of fine-grained sandstones

This paper introduces a novel three-dimensional grain-based hybrid finite-discrete element method parallelized based on general-purpose graphic processing units (GPGPUs) and applies it to investigate the failure process of sedimentary fine-grained sandstones. The grain-based method considers the actual microstructures of rocks with Voronoi or grain-growth tessellations to model their failures, including transgranular, intergranular and intragranular crack propagations. The novel semi-adaptive contact activation approach (semi-ACAA) proposed by the authors and the efficient tetrahedron-to-point (TtoP) contact interaction algorithm developed by HOSS are implemented to speed up grain-based modelling in addition to GPGPU parallelization. Semi-ACAA and TtoP are approximately 2–20 times and 1.5 times faster than the brute-force contact activation approach and tetrahedron-to-triangle contact interaction algorithm, respectively, which are prevalent in the FDEM community. The grain-based modelling elucidates that most intragranular cracks are observed in the medium- and high-strength grains of the sandstones, while most intergranular cracks occur at the grain boundaries between low- and high-strength grains. The transgranular cracks do not discriminate any grains on their propagating paths but propagate and coalesce with intragranular and intergranular cracks, which results in the final failure of the sandstones.

Investigation of the Dynamic Pure-Mode-II Fracture Initiation and Propagation of Rock during Four-Point Bending Test Using Hybrid Finite–Discrete Element Method

A hybrid finite–discrete element method (FDEM) is proposed to investigate dynamic pure-mode-II fracture behaviors. The transition of continuum to discontinuum was applied to the FDEM through the use of three fracture modes, so that the whole fracture process could be modeled naturally. The FDEM was then employed to model the dynamic pure-mode-II fracture behavior of rock during a four-point bending test with a prefabricated notch. The results showed that the fracture initiated from the tip of the prefabricated notch under a relatively lower loading rate, i.e., 1 m/s and 5 m/s. However, when the loading rate reached higher levels, i.e., 10 m/s and 50 m/s, the prefabricated notch played a small role in the fracture patterns. Under these conditions, the fracture initiated from the center of the beam bottom or the stress concentration vicinity, instead of the tip of the prefabricated notch. Regardless of the loading rate, the obtained force-loading displacement curves showed a typical brittle material failure process. Additionally, by incorporating the empirical correlation between the static and dynamic strengths obtained from the dynamic rock fracture tests, the hybrid finite–discrete element method could effectively reflect the impact of the loading rate on the strength of the rock. To conclude, the hybrid finite–discrete element method is an effective instrument to investigate the fracture initiation and propagation of rock, since it can both naturally simulate the process of rock fracture and capture the effect of the loading rate on the rock behaviors.

Study on the failure mechanism in shale-sand formation based on hybrid finite-discrete element method

During the hydraulic fracturing in shale formation, the mechanical properties and stress field caused by sand and shale interbedding were the major factors affecting the hydraulic fracture propagation. In this paper, a hydraulic fracturing simulation model was proposed based on the hybrid finite-discrete element method, and the factors affecting the fracture propagation were analyzed. The simulation results indicated that the fracture width experienced a sudden decrease when the hydraulic fracture interacted with a bedding plane. The reservoir with high stress difference and Young’s modulus was favorable for the opening of bedding plane, and high injection rate was favorable for the longitudinal expansion of hydraulic fracture. This was mainly because the lower stress difference and Young’s modulus was, the opening degree of the bedding plane was much smaller, and the less the energy consumption of opening the bedding surface was.

Hybrid Finite-Discrete Element Modelling of Various Rock Fracture Modes during Three Conventional Bending Tests

The numerical techniques for modelling the rock fracture have been briefly reviewed. A hybrid finite-discrete element method (HFDEM) is proposed to simulate various fracture types of rock. A fracture model is implemented in the HFDEM for simulation of the three main fracture types. In addition, the influence of the strain rate is considered during the HFDEM modelling rock behavior. Then, two typical rock mechanism tests are employed to calibrate the HFDEM. The proposed method has well modelled the rock fracture processes and can obtain reasonable stress distribution and force–displacement curves. After that, the HFDEM is used to model three convention bending tests. The obtained rock fracture processes indicates that the HFDEM can simulate various fracture types. The obtained rock strengths and fracture toughness indicate that the HFDEM can reflect the influence of the strain rate. It is concluded that the HFDEM can model the entire and complete rock fracture process during the three convention bending tests, and it also can capture the rock’s behavior on the strain rate.

Hydro-mechanical modelling of gas transport in clayey host rocks for nuclear waste repositories

Host rocks, as the final impediment to waste migration, play a significant role in the nuclear waste repositories. Modelling of gas migration in saturated host rocks as well as its coupled hydro-mechanical (HM) behaviors is of importance to the assessment of geological disposal facilities. A comprehensive literature review is carried out on the models for simulating gas migration in saturated rock materials. Several aspects have been provided and discussed, including the material properties and experimental interpretations, governing equations, constitutive models for hydraulic and mechanical processes, fracture propagation models. Specifically, these models are discussed in detail with respect to their performance in describing the recorded experimental behaviors. It is found that the visco-capillary two phase flow model with enriched intrinsic permeability is commonly used to describe the advective gas transport in saturated host rocks. The embedded fracture model (EFM) or enriched EFM seems to be the most favored model as it accounts for the fracturing mechanism which is more representative to implicitly simulate the preferential gas pathways. To describe the mechanical behavior of rocks during the gas migration process, linear/nonlinear elastic, elastoplastic and damage models have been included in the mechanical processes. The HM models incorporating plasticity or damage may not be applicable in most experimental studies, which are not competent to simulate the key experimental behavior associated with the development of gas preferential pathways. Current models that can explicitly describe the gas induced micro-fracturing are seldomly reported, the existed ones are not able to represent all the key experimental behaviors related to preferential gas flow. Advanced approaches, i.e., phase field (PF) method, extended finite element method (XFEM), discrete element method (DEM), hybrid finite discrete element method (FDEM) can be integrated with cohesive zone model (CZM) to provide some promising aspects in representing the development of preferential gas pathways. Lastly, the conclusions and recommendations for future modelling are given.

Development and application of a three-dimensional GPGPU-parallelized FDEM for modelling rock fragmentation by blast

Numerical modelling of rock fragmentation by blast is an important step to the optimum fragmentation. Combined finite-discrete element method (FDEM) has recently been proven to be one of the most promising methods for modelling the fracture and fragmentation of rocks under various loading conditions including rock blasting. However, almost all FDEM modellings of rock blasts are conducted in two-dimension. Correspondingly, this study implements a three-dimensional (3D) hybrid finite-discrete element method (HFDEM) parallelized on the basis of the general-purpose graphic-processing-unit (GPGPU) to model the fracture and fragmentation process of rock by blast. The main components of HFDEM3D include robust fracturing algorithm, efficient contact activation approach and various implementations of gas and rock interactions, which are firstly briefly reviewed here and introduced in detail in authors’ various former publications. The GPGPU-parallelized HFDEM3D is then applied to model the rock fracture and fragmentation process by the detonation of a single borehole and simultaneous and consecutive detonations of multiple boreholes.

Experimental and numerical studies on failure behaviours of sandstones subject to freeze-thaw cycles

The freeze-thaw induced damage of rock affects the durability and serviceability of geo-structures, especially those constructed in the regions frequently impacted by climatic changes. A series of laboratory tests including P-wave velocity tests, freeze-thaw tests, uniaxial compression strength (UCS) tests and Brazilian tensile strength (BTS) tests are conducted to investigate the physical-mechanical properties and failure behaviours of Tasmanian sandstones subjected to various freeze-thaw cycles. It is observed that the P-wave velocity, BTS and UCS of the sandstone decrease as the number of freeze-thaw cycles increases, in which the decreasing rate from 0 to 20 freeze-thaw cycles was more pronounced than that from 20 to 40 and 40 to 60 freeze-thaw cycles. Moreover, it is found that the main failure mode of the sandstone changes from axial splitting to shearing along a single plane in the UCS tests and from central smooth fractures to a central zigzag fracture in the BTS tests with the number of freeze-thaw cycles increasing. Three-dimensional (3D) numerical modellings are then conducted using a self-developed 3D hybrid finite-discrete element method (HFDEM) parallelized on the basis of general-purpose graphic processing units (GPGPU) to further investigate the failure mechanisms of Tasmanian sandstones subjected to various freeze-thaw cycles in the UCS and BTS tests. The 3D numerical modellings agree very well with the experimental observations that the physical-mechanical parameters of the sandstone degrade with the increasing number of the freeze-thaw cycles. Moreover, the 3D numerical modellings reveal the deterioration and failure mechanisms of sandstones subjected to various freeze-thaw cycles. For the sandstone specimens without subjecting to freeze-thaw cycles, axial splitting is the main failure pattern while tensile and mixed-mode damages are the dominant failure mechanism in the UCS tests. For the sandstone subjecting to various freeze-thaw cycles, the increasing number of freeze-thaw cycles causes the macroscopic cracks to propagate, interact and coalesce in the shear behaviour resulting in the final shear fracture pattern in the UCS test. The 3D numerical modellings of the BTS test show that, although, for both the models with and without subjecting to freezing and thawing cycles, a central fracture is the eventual failure pattern, the failure surface becomes more zigzag as the number of freeze-thaw cycles increases.