Multi-scale modeling of blasting-induced fracture in polycrystalline granite with grain boundary effects

2.1. Microscopic simulation method of granite crystal minerals

The mineral composition, content, and grain morphology of granite were determined via petrographic analysis, X-ray diffraction (XRD), and scanning electron microscopy (SEM). In this study, Barre Granite was adopted as the representative material [10]. Its main mineral components are feldspar and quartz, accounting for about 65 % and 35 %, respectively. Barre Granite is characterized by a fine-grained structure, with average grain sizes of approximately 1.25 mm for quartz and 1.0 mm for feldspar.

Fine-grained and coarse-grained granites share highly similar mineral composition, structural characteristics, and geological formation history; their primary difference lies in grain size and crystallization conditions [11]. Therefore, the grain size of Barre Granite was scaled up by a factor of four to represent medium-grained granite, giving average grain sizes of 5.0 mm for quartz and 4.0 mm for feldspar. This scaling method is reasonable because the mineralogical and textural properties remain consistent among different grain-size granite types, and only the crystallization kinetics differs. The scaling factor is also supported by empirical grain-size distribution statistics in granitic rocks [12].

The microstructural anisotropy of granite can be effectively represented by the Voronoi tessellation. The Voronoi diagram discretizes the mineral grains of granite into discrete seed points, and the generated polygonal network reasonably reflects the spatial distribution and contact mode of grains. It can well characterize the crack propagation path, mechanical behavior, and multi-field coupling response, making it widely used in microstructural analysis, fracture simulation, and engineering applications [13-15].

In this work, Neper software [16] was used to generate the polycrystalline microstructure of granite. The grain growth morphology was assigned to the cells, which provides a wider grain size distribution and higher grain sphericity than the standard Poisson–Voronoi tessellation. A square model with a side length of 50 mm was generated, containing 168 grains in total, including 109 feldspar grains and 59 quartz grains, as shown in Fig. 1(a). The geometric shape and mosaic structure of the generated grains are consistent with the actual granite section [3], as shown in Fig. 1(b).

Fig. 1Voronoi model of granite grains and cross-section of granite sample: a) Voronoi model with dimensions of 50 mm × 50 mm; b) cross-section of granite sample

a)

b)

2.2. Finite-discrete simulation method with inserted cohesive elements

To accurately reproduce the dynamic fracture process of granite subjected to blasting loads, cohesive elements must be inserted both within mineral grains and along grain boundaries simultaneously. Traditional numerical methods struggle to realize synchronous insertion of intragranular and grain boundary cohesive elements and realize local mesh refinement near the borehole, which restricts the accurate characterization of dynamic crack propagation and blasting vibration responses.

To overcome this limitation, Neper2CAE software [17] was used to convert the two-dimensional granite microstructure into a three-dimensional geometric model in Abaqus, completing the reconstruction of the polycrystalline numerical model (Fig. 2). The established Abaqus model retains complete grain set information, which provides a prerequisite for the targeted insertion of intragranular and grain boundary cohesive elements.

On this basis, the polycrystalline model was further cut and trimmed to obtain a cylindrical granite blasting model, as shown in Fig. 3(a). The inner bore diameter is 3.225 mm, the outer diameter of the cylinder is 20 mm, and the radial thickness is 16.775 mm. The model was meshed with 2989 six-node linear triangular prism elements (C3D6), as shown in Fig. 3(b). The element size near the borehole wall is refined to 0.5 mm and gradually transitions to 1 mm toward the outer boundary, which ensures computational accuracy and efficiency for dynamic blasting simulation.

Fig. 2Abaqus polycrystalline model and grain assembly

Fig. 3Polycrystalline model and mesh division of cylindrical granite specimen: a) polycrystalline model of cylindrical granite specimen; b) mesh division of blasted granite specimen model

a)

b)

Fig. 4Cohesive sets of models: a) Sets of grain interface; b) whole sets of grains

a)

b)

A mesh convergence study was conducted with three density levels to verify numerical robustness. The baseline mesh (inner wall 0.5 mm, outer 1 mm) contains 2989 C3D6 elements. Compared with finer mesh (0.3/0.6 mm) and coarser mesh (0.7/1.4 mm), the baseline mesh yields stable results: the maximum principal stress difference is less than 1.8 % between medium and fine meshes, and crack initiation positions remain consistent. The coarse mesh causes an 8.3 % deviation in peak stress and distorts crack paths at grain boundaries. The minimum element size of 0.5 mm (about 15.5 % of the borehole diameter) effectively maintains the circular borehole shape and matches the grain size characteristics (quartz 5.0 mm, feldspar 4.0 mm). This mesh density (4-10 elements per grain) achieves the best balance between geometric fidelity, dynamic response resolution, and computational efficiency.

In the Abaqus Part module, two key sets were established: a grain interface set and a full grain set, as shown in Figs. 4(a-b). Using the Cohesive Insertion V4.1 plugin with the Special Set function, cohesive elements were generated separately for grain boundaries and intragranular interfaces, as shown in Figs. 5(a-b), respectively. This synchronous insertion strategy enables independent parameter assignment and accurate characterization of intragranular fracture and intergranular fracture under dynamic blasting loads.

Fig. 5Cohesive elements at grain boundaries and intragranular regions: a) grain boundary cohesive elements; b) intragranular cohesive elements

a)

b)

The traction-separation constitutive model was adopted to describe the mechanical behavior of cohesive elements, which is characterized by initial linear elasticity followed by damage initiation and evolution. for simulating blasting-induced progressive failure, the Maximum Nominal Stress Criterion (MAXS) [18] was used, as it can reliably capture tensile-shear mixed failure dominated by stress waves under blasting dynamics.

The MAXS criterion is defined as Eq. (1):

where ${\sigma }_n$, ${\sigma }_s$, and ${\sigma }_t$ represent normal and tangential stress components, and $t_n$, $t_s$, $t_t$ are corresponding critical strengths. Damage initiates when the ratio reaches 1. This criterion is suitable for dynamic rock fracture because it independently evaluates stress components, highlights tensile failure dominated by blast stress waves, and ensures high computational efficiency for dynamic simulation.

Cohesive element parameters and rock mechanical parameters were taken from Liu et al. (2025) [2], as listed in Table 1.

To apply dynamic blasting loads, PETN (pentaerythritol tetranitrate) explosive was adopted [19], with key Jones-Wilkins-Lee (JWL) equation-of-state parameters given in Table 2. The initial explosion pressure was set to 1.6 GPa. a pressure-time curve consistent with the JWL equation was applied as the borehole boundary condition to realistically reproduce the dynamic evolution of blasting pressure. The total simulation time was 0.032 ms, which sufficiently captures the dynamic response and early-stage crack propagation of granite under ultra-small diameter blasting.

3.1. Influence of crystal microstructure on blasting dynamic responses and fracture characteristics

The numerical results in this study show good qualitative agreement with laboratory blasting test observations. The simulated failure zones include a crushing zone near the explosive source, a radial fracture zone in the middle region, and an outer circumferential tensile fracture zone induced by reflected stress waves, as presented in Figs. 6(a-c). Although direct quantitative comparison is limited by experimental measurement accuracy and material variability, the overall failure morphology is consistent with physical test results [3].

Compared with the homogeneous continuum model, the polycrystalline granite model exhibits obvious non-uniformity and anisotropic characteristics in crack propagation, which can better reflect the real dynamic fracture behavior of granite under blasting. This qualitative verification method is consistent with previous studies [20, 21] that focus on fracture pattern consistency under dynamic loading.

Fig. 6Crack propagation and granite blasting test results: a) crack propagation diagram of the amorphous model; b) crack propagation diagram of the polycrystalline model; c) results of granite blasting test

a)

b)

c)

In the polycrystalline model, initial cracks can initiate both inside grains and at grain boundaries [22, 23]. Even if identical mechanical parameters are assigned to intragranular and grain-boundary cohesive elements, a large number of intergranular cracks still occur. This indicates that crack propagation is significantly controlled by grain morphology and shows strong dynamic anisotropy, as illustrated in Fig. 7(a).

In addition, blast-driven cracks pass through adjacent quartz and feldspar grains without obvious obstruction (Fig. 7(b)). This phenomenon reveals that the difference in mineral mechanical properties has limited influence on crack evolution under blasting dynamics. in other words, the mechanical parameters of cohesive elements play a more dominant role in controlling dynamic crack propagation paths than the bulk properties of mineral grains in the FDEM model.

Fig. 7Crack distribution and mineral distribution in polycrystalline model: a) crack distribution in the polycrystalline model during blast wave propagation; b) mineral distribution in the polycrystalline fracture model

a)

b)

3.2. Influence of grain boundary strength on dynamic fracture type under blasting

A parametric study was carried out to investigate the dynamic fracture mode transition regulated by grain boundary strength in granite under blasting load. Three groups of cohesive element parameters were designed to represent different grain boundary strengths, while the intragranular parameters remained unchanged.

In the baseline case, identical strength parameters were applied to both intragranular and grain boundary cohesive elements ($t_n=$30 MPa, $t_s=$60 MPa, $t_t=$60 MPa). the results show that transgranular fracture dominates the failure pattern, with cracks propagating directly through mineral grains and only local intergranular cracking occurring at individual boundaries (Fig. 8(a)).

When grain boundary strength is moderately reduced ($t_n=$25 MPa, $t_s=$ 55 MPa, $t_t=$55 MPa), intergranular fracture increases noticeably. Cracks tend to propagate along relatively weak grain boundaries, forming a mixed transgranular-intergranular failure mode (Fig. 8(b)).

When grain boundary strength is further decreased to approximately two-thirds of the original level ($t_n=$20 MPa, $t_s=$40 MPa, $t_t=$40 MPa), the failure mode changes significantly. The rock mass is dominated by intergranular fracture, with cracks extensively developing and extending along grain boundaries (Fig. 8(c)).

These results demonstrate that grain boundary strength serves as a key controlling parameter for determining the dominant fracture type under dynamic blasting load. With decreasing grain boundary strength, the failure mode gradually transforms from transgranular fracture to mixed fracture and finally to intergranular fracture, which directly affects the dynamic energy dissipation, crack propagation path and blasting-induced fracture characteristics of granite.

Fig. 8Crack distribution under different grain boundary strengths: a) crack distribution with grain boundary strength tn= 30 MPa, ts= 60 MPa, tt= 60 MPa; b) crack distribution with grain boundary strength tn= 25 MPa, ts= 55, tt= 55; c) crack distribution with grain boundary strength tn= 20 MPa, ts= 40 MPa, tt= 40 MPa

a)

b)

c)

3.3. Dynamic stress transfer and crack initiation in the near-borehole region

Dynamic stress distribution and energy transfer near the ultra-small borehole are critical to understanding the blasting-induced fracture mechanism. at the early stage of blasting ($t=$ 3.4 µs), the ultra-small borehole causes severe stress concentration within the mineral grains immediately surrounding the explosive. To accurately characterize the complex dynamic stress state under blast loading, the Von Mises equivalent stress (${\sigma }_{e}q$) is adopted as the evaluation indicator, since it comprehensively reflects the combined effect of normal and shear stresses under dynamic conditions.

As shown in the Von Mises contour (Fig. 9(a)), widespread intragranular cohesive cracks initiate inside the highly stressed grains near the borehole wall. at the later moment $t=$ 7.3 µs, significant stress concentration zones appear at partial grain boundaries, where the Von Mises stress can reach up to 200 MPa (Fig. 9(b)). These high-stress regions coincide exactly with the initiation sites of new grain-boundary cracks and the tips of propagating cracks, revealing the key role of dynamic stress redistribution in controlling crack initiation and propagation.

This phenomenon indicates that under moderate grain boundary strength, the blast energy is first transmitted intragranularly in the near-field, and then released and dissipated at the weaker grain boundaries. Such a mechanism forms a typical chain-type dynamic failure path: intragranular energy transfer to grain boundary fracture, which dominates the anisotropic fracture evolution and energy dissipation characteristics of granite under ultra-small diameter blasting.

Fig. 9Von Mises stress distribution in near-explosion zone: a): Von Mises stress distribution contour plot of the near-explosion zone at t= 3.4 µs, illustrating element-averaged stress magnitudes within the range of 0-200 MPa. the color gradient from blue (low stress) to red (high stress) indicates the intensity of the equivalent stress state; b): Von Mises stress distribution contour plot of the near-explosion zone at t= 7.3 µs, also illustrating element-averaged stress magnitudes within the same range. High-stress zones (red regions) align with crack initiation sites at grain boundaries

a)

b)

4.1. Mechanisms of anisotropic dynamic crack propagation

The anisotropic and non-uniform dynamic crack propagation observed in the polycrystalline model (Section 3.1) originates from the inherent meso-heterogeneity of mineral grains and grain boundaries. Grain boundaries naturally act as weak planes with distinct mechanical properties (strength, fracture energy) compared with the intragranular crystal lattice [22, 23]. When a blast-induced stress wave or crack front impinges on a grain boundary under dynamic loading, several typical responses may occur: crack arrest, deflection, transmission across the boundary, or new crack initiation along the boundary itself.

The final propagation behavior is jointly governed by the boundary orientation relative to the principal stress, the strength contrast between grain interior and boundary, and the local dynamic stress state. The present numerical results confirm that the properties of cohesive elements assigned at grain interfaces play a more decisive role in dynamic crack propagation paths than the bulk mechanical properties of mineral grains, which is consistent with the fundamental understandings in computational dynamic fracture mechanics.

This mechanism explains the significant anisotropic fracture characteristics in polycrystalline granite under ultra-small borehole blasting, and provides a micro-scale theoretical basis for analyzing the dynamic stress wave propagation, crack evolution, and blast-induced vibration responses from a mesoscopic viewpoint.

4.2. The role of grain boundary strength in dynamic energy dissipation and fracture path

The parametric analysis on grain boundary strength (Section 3.2) clearly demonstrates that this parameter dominates the dynamic fracture mode and energy dissipation behavior of granite under blasting. The transition from transgranular to intergranular failure with decreasing grain boundary strength is physically consistent and mechanically reasonable.

Stronger grain boundaries force dynamic cracks to propagate through mineral grains, while weaker boundaries provide preferential paths for crack growth and induce interfacial debonding along grain boundaries.

This difference significantly affects dynamic energy dissipation during blasting. Transgranular fracture typically dissipates more energy due to the breaking of strong intragranular crystal bonds. in contrast, intergranular fracture along weak boundaries represents a more energy-efficient failure mode under dynamic loading.

The chain-type dynamic failure path revealed in Section 3.3 further highlights this energy transfer mechanism: blast energy is first absorbed and transferred inside grains near the explosion source, then released efficiently by fracturing weaker grain boundaries.

This energy transfer and dissipation routine directly determines the dynamic fragmentation efficiency, crack distribution, and blast-induced vibration characteristics of granite under ultra-small diameter blasting. It also provides a mesoscopic mechanism for understanding and regulating the dynamic fracture behavior and vibration response in polycrystalline rock blasting.

4.3. Implications for dynamic blasting mechanics and engineering practice of ultra-small diameter boreholes

The multi-scale FDEM model established in this study, which considers explicit grain microstructure and distinguished grain boundary properties, provides an effective tool for revealing the dynamic micromechanics of ultra-small diameter blasting in crystalline granite. The sensitive response of fracture mode to grain boundary strength indicates that the mineral composition, microstructure, and especially grain boundary characteristics of actual granite will significantly affect its dynamic blasting response and vibration generation mechanism. This understanding breaks through the traditional assumption that rock is a homogeneous continuum.

For engineering applications including rock fragmentation, blasting vibration control, and dynamic disaster prevention, the main implications are as follows:

1) Blasting parameter optimization. Clarifying the dominant fracture mode controlled by grain boundary strength can help optimize explosive type, charge concentration, and initiation timing. by promoting a favorable energy dissipation path, it is possible to achieve expected fragmentation quality and simultaneously reduce unfavorable blast-induced vibration.

2) Blast-induced vibration control. The anisotropic crack propagation and the interaction between stress waves and microstructure directly affect the amplitude, frequency distribution and propagation law of blasting vibration. the multi-scale dynamic model proposed in this study can improve the prediction accuracy of blasting vibration and provide a theoretical basis for vibration reduction design.

3) Dynamic rock characterization. The key role of grain boundary properties highlights the necessity of using advanced microscopic testing methods to obtain key meso-parameters for specific rock masses, so as to realize more accurate dynamic fracture and vibration response prediction.

This study shows that incorporating grain boundary effects into multi-scale dynamic modeling is crucial for revealing the micro-mechanism of granite blasting, and provides a reliable method for dynamic response analysis, vibration prediction and parameter optimization of ultra-small diameter blasting in engineering practice.