Shock-induced fracture of dolomite rock in small-scale blast tests

Pweł Brnowski,Michł Kucewicz,Mteusz Pytlik,Jerzy Młchowski

a Military University of Technology,Faculty of Mechanical Engineering,Institute of Mechanics & Computational Engineering,2 Gen.S.Kaliskiego Street,00-908,Warsaw,Poland

b Central Mining Institute,Conformity Assessment Body,40-166,Katowice,Poland

Keywords:Johnson-Holmquist II model (JH-2)Rock Blasting Fracture Small-scale test

ABSTRACT This paper attempts to study dolomite failure using small-scale blast tests.The experimental setup consisted of a cylindrical specimen with a central borehole fitted with a detonation cord inside a copper pipe.The specimen was confined using lead material.During the test,acceleration histories were recorded using sensors placed on the lead confinement.The results showed that heterogeneity and initial cracks significantly influenced the observed failure and cracking patterns.The tests were numerically represented using the previously validated Johnson-Holmquist II (JH-2) constitutive model.The properties of the detonation cord were first determined and verified in a special test with a lead specimen to compare the deformation in the test with that of numerical simulation.Then,the small-scale blast test was simulated,and the failure of the dolomite was compared with the test observations.Comparisons of acceleration histories,scabbing failure,and number of radial cracks and crack density confirmed the overall repeatability of the actual testing data.It is likely that the proposed model can be further used for numerical studies of blasting of dolomite rock.

1.Introduction

Rock fracturing via blasting has been widely studied by combining experimental tests with numerical simulations.Numerical simulations can efficiently model the expansion of high explosive (HE) materials and their interactions with rock or other brittle materials.However,the model scales of associated tests vary significantly in relevant studies (e.g.Yi,2013a,b; Hu et al.,2015,2017;Yi et al.,2016; Huo et al.,2020;Pytel et al.,2020,2021).The results of bench and cutting blasting are discussed (e.g.Cho and Kaneko,2005; Qu et al.,2008; Yi,2013b; Zheng et al.,2015).As the full-scale experiments and/or field tests are more costly and time-consuming than those performed at smaller scales,the latter is often adopted to investigate various materials subjected to blast loading (Katsabanis et al.,2006; Dehghan Banadaki and Mohanty,2012; Li and Shi,2015; Kukolj et al.,2018,2019; Chi et al.,2019;Gharehdash et al.,2020; Liu et al.,2020).

In early small-scale tests,Field and Ladegaard-Pedersen (1971)demonstrated the significance of controlling stress wave interactions through the position and time of HE detonation or the choice of free surface geometry because the reflection of the stress wave can induce spalling.Wilson and Holloway (1987) used concrete blocks with a size of up to 1 m3to study fracture and fragmentation in a laboratory setup with eight single-hole and two time-delayed two-hole tests.High-speed photographs were used to analyze the sequence and patterns of crack formation.Dehghan Banadaki and Mohanty (2012) presented accurate experimental tests and modeling based on single-hole blasting of granite rock specimens and analyzed the crack patterns resulting from the stress wave.The authors estimated parameters for the Johnson-Holmquist II (JH-2) model based on the measured mechanical data and blast-induced fracture and pressure values recorded from gauges inside the specimens.The JH-2 model was successfully calibrated and precisely reproduced the crack pattern.Nevertheless,two main drawbacks remained: experimental and numerical agreement was obtained only at the macro level,and plane-strain modeling was considered,thereby failing to analyze vertical fracture propagation in the specimen.

Several scientists have successfully used Dehghan Banadaki and Mohanty’s (2012) results to validate different models and simulation propositions.For instance,Gharehdash et al.(2020) utilized a fully meshless model with different particle resolutions and various patterns of particle arrangement and tensile instability.All parts of the model were represented using smoothed particle hydrodynamics (SPH),which was capable of reproducing blast-induced fractures.However,this approach requires the use of several million particles and adjustment of the numerical parameters of the SPH method.Then,an alternative,the three-dimensional (3D)multi-material arbitrary Lagrangian-Eulerian (MM-ALE) formulation,was effectively implemented for the numerical reproduction of Dehghan Banadaki and Mohanty’s(2012)small-scale laboratory test.Wang et al.(2018a) utilized the JH-2 model to simulate a granite specimen and obtained satisfactory agreement of the crack pattern and measured pressure values with experimental results.The MM-ALE method has been widely adopted in studies of rock blasting (e.g.Zhu et al.,2007; Xie et al.,2019; Liu et al.,2020; Tao et al.,2020; Wang et al.,2021).Zhu and Zhao (2021) used a numerical approach based on a peridynamic model to predict the small-scale blasting behavior of rock and compared the numerical outcomes with those obtained by Dehghan Banadaki and Mohanty(2012).

Other noteworthy analyses were reported by Kukolj et al.(2018,2019),who investigated blast-induced fractures using experimental methods and numerical simulations.A cylindrical specimen was inserted in a special blast chamber that radially confined the cylinder while simultaneously capturing crack formation during blasting and reducing end-face spalling.Furthermore,internal crack patterns and fragmentation were analyzed by computed tomography(CT)and scanning electron microscopy(SEM).Numerical simulations realistically reproduced the crack patterns and fragment size distributions.Although the majority of studies have investigated rock,a very recent paper(Jeong et al.,2020)presented visual observations of crack propagation in polymethyl methacrylate (PMMA) cylindrical specimens.Experimental observations were conducted with the aid of a high-speed camera to capture fracturing processes and used to validate the developed numerical model.The authors highlighted the influence of stemming on shock-induced and gas-induced fractures.A detailed study of the impact of stemming on rock fracture behavior was presented by Zhang et al.(2020),who considered two constraint conditions and at least one specimen with and without stemming in each constraint type.It shows that stemmed blasts resulted in a better and finer fragmentation,whereas blasting without stemming wasted approximately 25% of the explosion energy.

Small-scale laboratory studies are typically used to validate numerical blasting models for further studies at different scales.However,such analyses require a reliable constitutive model.There are numerous constitutive models dedicated to simulating the behaviors of brittle rock materials.The models include the Johnson-Holmquist concrete model(JHC model)(Holmquist et al.,1993;Ren et al.,2017; Kucewicz et al.,2021),the Riedel-Hiermaier-Thoma(RHT) model (Mazurkiewicz et al.,2016; Yi et al.,2017),the brittle damage model (Wei et al.,2009; Morales et al.,2013),the Karagozian and case concrete (KCC) model (Malvar et al.,1997;Mardalizad et al.,2018; Kucewicz et al.,2020),the continuous surface cap model(Schwer and Murray,2002;Morales et al.,2013)and the JH-2 model (Holmquist et al.,1995,2001; Johnson and Holmquist,2008; Dehghan Banadaki and Mohanty,2012; Wang et al.,2018a,b; Baranowski et al.,2020a).

The present paper is part of a wider project aimed to simulate and optimize parallel cut-hole blasting.Several constitutive models have been tested in previous studies of the numerical reproduction of the fracture behavior of dolomite rock (Kucewicz et al.,2020,2021),including the KCC and JHC constitutive models.However,the JH-2 model was most effective in reproducing the cracking,fragmentation and large-scale failure of an entire rock face with several cut-holes (Baranowski et al.,2020a,b).Therefore,the JH-2 constitutive model was utilized in the present paper to further study the fracture of dolomite rock using small-scale blasting experiments and numerical simulations.This study provides a useful complement to the large number of rock blasting studies available for reproduction by other researchers and can serve as a reference for further studies of rock fracture and fragmentation.

The present paper makes several contributions to the literature on small-scale blast tests.First,a methodology for analyzing and determining cord properties based on analytical,experimental and numerical methods is proposed.Second,few of the numerous studies of rock blasting and fracture have generated experimental data.Therefore,the original laboratory setup presented here provides a valuable reference for future studies in the field.Finally,none of the numerous studies combining a constitutive modeling,simulations and experimental tests has examined dolomite rock,with the exception of the works published previously by the authors.

2.Small-scale blast tests

The dolomite specimens investigated in this paper have been previously described (Baranowski et al.,2020a; Kucewicz et al.,2020,2021).The specimens were obtained from Lower Silesia,Poland.Dolomite rock is generally strongly heterogeneous,with numerous imperfections,inclusions,pores and initial cracks.Consequently,specimens for experimental tests must be prepared carefully to maintain an intact structure.Typical dolomite specimens are shown in Fig.1.CT (Fig.1a) and microscopic analyses(Fig.1c),and macroscopic observations (Fig.1b) clearly revealed inclusions and initial cracks that will impact the fracture behavior of the material.The specimen presented in Fig.1a was investigated in the present paper,and the influence of the initial cracks was observed and considered.

Four cylindrical dolomite specimens with diameter and height of both 130 mm were used in the tests.A hole with a diameter of 8 mm was drilled at the central point of the specimen,and a copper pipe with a wall thickness of 1.5 mm was tightly installed in each specimen to prevent penetration of the explosion gas into the fractured dolomite (Dehghan Banadaki and Mohanty,2012).The copper pipe could easily expand inside the hole without failure,resulting in continuous interaction with the dolomite specimen.The detonation cord was inserted through the whole length of the copper pipe.The dolomite sample was placed in a steel pot with an inner diameter of 200 mm,which was then flooded with liquid lead to its height before placing the HE.The lead covering the specimen was used as a confinement and to prevent extensive failure of the specimen.Two piezoelectric sensors (PCB Piezotronics,Depew,NY),sensors Nos.1 and 2,were placed on the top of the lead confinement tube at a distance of approximately 84 mm from the center of the borehole.The sensors were connected to a data acquisition system(SIRIUS HS,Dewesoft,Slovenia),and the vertical acceleration history was recorded with a measurement range of up to 1×106m/s2and a sample frequency of 1 MHz.The experimental setup is illustrated in Fig.2.

3.Numerical model

The small-scale blast tests were reproduced using numerical simulations.To properly simulate the blast tests,additional tests were conducted to establish the parameters for the lead constitutive model and to determine the cord properties based on hole expansion of the lead material.

All simulations were performed using an explicit integration procedure with massively parallel processing(MPP)LS-Dyna code.The MM-ALE algorithm was used to simulate the detonation cord and its interaction with lead or dolomite,depending on the specific case.Previous studies used the finite element model (FEM)-SPH approach(Baranowski et al.,2020a,b),but the size of the HE model in the present paper is smaller than that considered in previous works,resulting in a very small time step in the calculations.Furthermore,it would be necessary to simulate the interactions of the SPH parts(HE,lead and copper pipe),which could influence the computational time.In addition,large deformation occurred in some parts,such as the copper pipe.Thus,in this paper,all parts of the detonation cord were treated as Eulerian parts,and the interactions between them were introduced by default.By contrast,brick finite elements with one integration point and hourglass control were adopted to model the other parts,depending on the analyzed case.The interactions between these parts were defined using a penalty-based contact procedure (Hallquist,2019).The same method was adopted for fluid-structure-interaction between the parts of the detonation cord and the Lagrangian parts.In the simulations,the second-order advection procedure was implemented(Hallquist,2019).

Fig.1.Tested dolomite specimen:(a)CT scan of the cylindrical specimen studied in the present paper;(b)Photograph of the specimen used in compression tests(Kucewicz et al.,2021); and (c) Microscopic image of the specimen used in the split Hopkinson pressure test (Kucewicz et al.,2021).

Fig.2.Experimental setup for the small-scale blast test with the dolomite specimen: (a) Close-up view of the detonation cord with detonator; (b) Steel pot with the specimen inside; (c) Close-up view of the detonation cord and its dimensions; (d) Top view of the experimental setup and (e) Photograph of the specimen and lead confinement with the sensors.

3.1.Model definition of the small-scale blast test

The conditions used in the numerical simulations matched with those in the experimental tests.A quarter of the model was considered,and the appropriate symmetric conditions were applied on the required faces of the model.The MM-ALE approach was adopted; all parts of the detonation cord were treated as Eulerian parts,whereas the dolomite specimen,lead confinement tube and steel pot were modeled using Lagrangian elements.The finite elements near the borehole were embedded in the Eulerian air domain,for which a non-reflecting boundary was introduced on its outer surfaces to consider the flow of the pressure outside the air domain.Air was considered as a simple ideal gas with vacuum properties.The JH-2 model was used with the parameters in Table 1 for the dolomite specimen.However,according to Baranowski et al.(2020a),the value of the maximum hydrostatic tensile pressure(T)must be adjusted depending on the finite element size.In the present paper,the dolomite specimen was discretized using a 0.8-mm mesh,which was also used by Baranowski et al.(2020a),and thus a value of T=43 MPa was used.During the tests,the steel pot remained intact,and consequently an elastic material model with the following properties was adopted for the steel pot:Esteel=210 GPa,ρsteel=7850 kg/m3,and νsteel=0.3.The lead sheath and lead confinement were modeled using the modified Johnson-Cook (MJC) constitutive model with the verified and correlated(see Section 4.1) parameters presented in Table 2.To simulate the blasting properties of the cord,the high explosive burn (HEB)constitutive model with the Jones-Wilkins-Lee (JWL) equation of state (EOS) was adopted (Table 3).The cord properties weredetermined from an additional test based on hole expansion of the lead material (see Section 4.2).Some of the JWL-EOS parameters were determined using an analytical approach.Additionally,the copper pipe inserted in the borehole was treated as a Eulerian part and modeled using the Johnson-Cook(JC)constitutive model with the parameters presented in Table 4.

Table 1 Dolomite properties of the JH-2 model (Baranowski et al.,2020a,b).

Table 2 Lead material properties of the MJC model (Børvik et al.,2009; Kedzierski et al.,2019).

For the Eulerian parts,an element size of 0.15-0.3 mm was considered in the central area of the model,where major coupling between the Eulerian parts and Lagrangian dolomite occurred.By contrast,an average element size of 1.2 mm was used for the lead confinement tube and steel pot.A total of 997,445 elements were generated.To prevent non-physical leakage in the fluid-structureinteraction,several parameters were introduced,e.g.interface stiffness as the relation between pressure and penetration,a minimum volume fraction of detonation products in Eulerian elements to activate coupling (FRCMIN) of 0.01,and an increase in the number of coupling points distributed over each coupled Lagrangian surface segment (NQUAD) to 3 × 3.During the simulations,the vertical acceleration histories were measured in a node representing the gauge placement in the laboratory tests.The model with the appropriate initial boundary conditions is presented in Fig.3.

3.2.Constitutive model

3.2.1.Dolomite rock

The dolomite was described using the JH-2 constitutive model with the previously determined parameters presented in Table 1,and a detailed procedure can be found in Baranowski et al.(2020a).The model describes the relationships between normalized pressure and normalized equivalent stress for intact,damaged,and fractured surfaces (Fig.4).

Table 3 RDX material properties of the HEB model with the JWL-EOS (Cao et al.,2016).

Table 4 Copper material properties of the JC model(Johnson and Cook,1983;Kucewicz et al.,2019).

The normalized intact strength of the material presented in Fig.4 can be described using the following formula (Holmquist et al.,1995; Johnson and Holmquist,2008):

where σ*I=σI/σHELis the normalized intact equivalent stress,σIis the current equivalent stress,and σHELis the equivalent stress at the Hugoniot elastic limit (HEL); P*=P/PHELis the normalized hydrostatic pressure,P is the current hydrostatic pressure; T*=T/PHELis the normalized maximum tensile hydrostatic pressure; ˙ε*=is the dimensionless strain rate,is the current equivalent strain rate,and=1 s-1is the reference strain rate.The damaged state σ*Dis determined by(Holmquist et al.,1995;Johnson and Holmquist,2008):

Fig.3.(a) Close-up view of the detonation cord and borehole,and (b) Numerical model of the small-scale blast test of dolomite rock with initial boundary conditions.

Fig.4.Intact,damaged and fractured surfaces describing the JH-2 constitutive model(Baranowski et al.,2020a).

where σ*Fis the fractured surface,and D is a damage index and takes a value between 0 (undamaged) and 1 (fully damaged).

The fractured surface σ*Fof the material is described using the following formula(Holmquist et al.,1995;Johnson and Holmquist,2008):

where B and M are the fractured material constants.

3.2.2.Lead alloy used in the detonation cord sheath and confinement tube

The lead confinement tube and lead sheath covering the cord and hexogen with RDX were simulated using the MJC model with parameters adopted from previous studies (Børvik et al.,2009;Ke˛dzierski et al.,2019),see Table 2.

The MJC model includes the effects of damage,strain hardening,thermal softening and kinematic strengthening on yield stress,which can be described using the following equation(Borvik et al.,2001):

where εpis the equivalent plastic strain,and ˙εpis the equivalent plastic strain rate.

For the sake of credibility of the results,the MJC model with the parameters presented in Table 2 for the lead alloy used in the sheath and cylindrical specimen was verified and correlated in a quasi-static uniaxial compression test.Cylindrical samples with height and diameter of both 40 mm were used.In the numerical simulations,the specimen was discretized using eight-node solid elements with one integration point.The specimen was placed between two rigid walls:one was fixed,while the other was able to move in the axial direction.An explicit numerical solution was adopted,and the prescribed velocity function was applied for the moving surface to reduce inertial effects.The authors have effectively used this method in previous studies (Kucewicz et al.,2020,2021).To verify the robustness of the MJC model with the parameters for the alloy,the engineering stress and engineering plastic strain curves were compared with the experimental outcomes.Furthermore,to prevent non-physical kinematic hardening,the MJC JC strain rate sensitivity parameter (CMJC) was set to zero for the analyses.Fig.5 presents the quasi-static uniaxial compression model with corresponding initial boundary conditions.

3.2.3.RDX used in the detonation cord

A detonation cord with RDX as the explosive material covered with a lead sheath was used.The outer diameter of the cord was 4.9 mm,and the diameter of the RDX core was 1.67 mm.The measured values of linear density and volume density were 2.806 g/m and 1281 kg/m3,respectively.Furthermore,the average detonation velocity measured by electronic probes was equal to DHE=4840 m/s.The RDX and its detonation products were described using the HEB model with the JWL-EOS.The HEB model includes a burn fraction (F) that controls the release of chemical energy and multiplies an EOS.The pressure in a single finite element of HE is calculated as follows (Hallquist,2019):

where pEOSis the pressure from the EOS,Vris the relative volume,and Edis the internal energy density per unit initial volume.

The JWL-EOS is described as follows (Hallquist,2019):

where VHE=ρ0HE/ρHE,ρ0HEis the initial density of the HE;ρHEis the actual density of the HE; E0is the detonation energy per unit volume and initial value of E of the HE;and AHE,BHE,R1,R2,and ω are empirical constants determined for a specific type of explosive material.

The basic parameters for the HE used in the present paper were taken from the literature(Cao et al.,2016).However,the analytical procedure proposed in Kamlet and Jacobs (1968) and effectively used in Cao et al.(2016) was adopted to determine some of the parameters of the JWL-EOS and HEB material model.Ultimately,the detonation heat Q,Chapman-Jouguet pressure pCJ,and initial energy E0were determined to be 3595 kJ/kg,8262 MPa and 4606 MPa,respectively.

Although the properties of the HE used in the present study are well studied,an additional test was conducted to correlate the cord properties,as these properties may differ within a certain range.Selected constants of the HEB model and JWL-EOS were iteratively changed until satisfactory agreement between lead specimen deformations in experimental tests and numerical simulations was achieved.In the experimental test,a cylindrical specimen made of lead with a diameter identical to that in the quasi-static uniaxial compression test was used.The specimen had a height of 37.5 mm,and a hole with a diameter of 5 mm was drilled in its center.A cord with a length of 137.5 mm was inserted in the hole of the cylindrical specimen such that the initiation point of detonation was 50 mm from the bottom of the specimen.The overall dimensions of the specimen were thoroughly measured before and after the test.To validate the numerical model of the detonation cord,the final deformation of the specimen and the diameter of the hole were compared.The simple experimental setup is illustrated in Fig.6a.

The abovementioned test was simulated to properly adjust the constitutive constants of the detonation cord.The modeling methodology used in this case was very similar to that mentioned earlier for simulating the small-scale blast test,i.e.quarter of the model,symmetry conditions,MM-ALE simulation,material models and their parameters.The lead specimen was modeled using Lagrangian elements,whereas all parts of the detonation cord were treated as Eulerian parts and were embedded in the air domain,which was considered as ideal gas with vacuum properties.On the outer surfaces of the Eulerian parts,a non-reflecting boundary was introduced to enable pressure flow outside the air domain.The lead specimen and lead sheath were modeled using the MJC constitutive model with the parameters presented in Table 2.To simulate the blasting properties of the cord,the HEB material model with the JWL-EOS was adopted with the final parameters presented in Table 3.For uniformity,the adopted element sizes were the same as those used in the numerical simulations of the small-scale blast tests.The model and initial boundary conditions are shown in Fig.6b.

Fig.5.Specimen of lead alloy for numerical simulations of the uniaxial quasi-static test with applied initial boundary conditions.

3.2.4.Copper used for pipe inserted in the borehole

The copper pipe was described using the JC model with parameters taken from Johnson and Cook (1983)and Kucewicz et al.(2019) (Table 4).The model provides the prediction of the flow stress σflowfor large strains and high strain rates and is given by the following relation (Hallquist,2019):

where AJC,BJC,CJC,NJC,and MJCare the material constants;TrJC,TmJC,and T are the reference,melting and absolute temperatures,respectively.The other constants have the same meaning as in the MJC model.

4.Results and discussion

4.1.Lead properties verification

The results of the quasi-static uniaxial compression tests of the lead alloy are shown in Fig.7.Excellent correlation of the numerical results with the actual tests can be observed,and the post-yield curves are in good agreement.The curve obtained from numerical simulations starts at the defined yield point equal Re=21.2 MPa,which corresponds to the experimentally determined value under parameter AMJCin the MJC model(see Table 2).Thus,the literature data of the MJC model for lead alloy can be used to effectively simulate the lead material used in the described test.

4.2.Cord properties verification

Fig.6.Testing of the cord properties: (a) Experimental scheme and (b) Numerical model with initial boundary conditions.

Fig.8 presents the comparison of lead specimen deformation between the numerical simulation and experimental results.To better visualize the results,the specimen and its numerical representation are cut in half.The expansion of the detonation cord resulted in plastic deformation of the lead specimen.The diameter was measured in several places across the length of the specimen,and the average value was obtained.The diameter of the specimen changed from Dinitial=40 mm to Davg=41.42 mm in the actual test.In the numerical model of the lead specimen,Davg=40.93 mm was measured.The hole expanded significantly:the diameter along its length was not uniform,and larger deformation occurred at the end of the hole.The diameters were similar,with values of approximately Dend=12.82 mm and Dend=12.41 mm in the experiment and numerical simulation,respectively.In the central part of the specimen,diameters of Dcenter=10.95 mm in the experimental specimen and Dcenter=10.89 mm in the numerical simulation were observed.During the numerical simulations,the values of BHE,R1HEand R2HEin the JWL-EOS (Eq.(6)) were iteratively changed until satisfactory agreement was obtained.Ultimately,values of BHE=4690 MPa,R1=4 and R2=0.4 were selected.The actual tests were well reproduced in the numerical simulations,and the maximum deviation between results was approximately 3% and occurred in the diameters measured at the end of the specimen hole.Notably,at the middle height of the specimen,excellent agreement between the numerical and actual outcomes was achieved,with an error less than 1% (Table 5).Therefore,the HEB material model with the JWL-EOS correctly describes the actual detonation cord properties and can be used for further analyses.

Table 5 Comparison of lead specimen deformation obtained from simulation and experiment.

4.3.Small-scale blast test

4.3.1.Experimental results

Fig.7.Comparison of engineering stress with engineering plastic strain curves obtained from the experimental test and numerical simulation.

After verifying the properties of the detonation cord,small-scale blasting studies were performed.The vertical acceleration histories measured with sensors Nos.1 and 2 are presented in Fig.9 for the four dolomite specimens (No.1 to No.4).The acceleration history from sensor No.1 for specimen No.4 has an unphysical characteristic,and the results from this measurement were not further analyzed.All signals have a damped vibration characteristic,and some correspondence with a typical HE pressure vs.time can be observed.The sudden pressure growth responsible for the initial positive peaks produces a negative slope of the curve.The results are generally repeatable,although there are some discrepancies between specimens.The first maximum acceleration peaks have values between 4,509 m/s2and 12,804 m/s2.Furthermore,there are slight phase differences between the measurements.The observed differences can be attributed to the heterogeneous structure of the dolomite and pre-crack occurrence inside the specimens,which leads to reflection and weakening of stress waves in the sample.The CTof the selected specimen revealed the existence of pre-crack surfaces,and a previous study (Kucewicz et al.,2021) reported that the tested dolomite included pores and inclusions varying in size from 1 mm to 10 mm.Moreover,some initial cracks can be observed on the specimen face,and these cracks extend to the interior of the specimen(see Fig.1).Such imperfections may impact rock brittle behavior,crack propagation and final failure of the specimen,as evidenced by the failure of samples with similar inclusions under uniaxial and triaxial compression in previous laboratory tests by the authors(Kucewicz et al.,2020,2021).The induced stress state is complex,and it may be assumed that such imperfections can be crack initiators.Finally,all acceleration histories from the experimental tests were compared with the results of the numerical simulations.

An exemplary specimen failure is presented in Fig.10.Two characteristic zones in the vicinity of the blast hole can be distinguished: a crushing zone around the hole and a fracture zone,which consists of mainly radial cracks.Compared with the bottom surface,the intensity of fracture is clearly more visible on the top surface,where a larger number of cracks occurred.This discrepancy is due to the location of the detonation point at the bottom side of the specimen.Consequently,the effect of stress wave superposition plays a key role in the extent of fracture of the top surface.Similar phenomena have been reported and discussed in other studies in terms of small-scale rock blasting (Dehghan Banadaki,2010;Jayasinghe et al.,2019; Gharehdash et al.,2020).At a distance of approximately 30 mm from the top surface,a fracture surface can be seen.As the stress wave propagated parallel to the axis of the blast hole,it struck the top free specimen surface.Part of the stress wave was reflected as a tensile stress wave,resulting in spalling failure.Additionally,an extensive crushing zone due to scabbing and compressive damage can be observed in the top surface.

Fig.8.Comparison of lead specimen deformation after detonation cord expansion: (a) Numerical simulation and (b) Experimental results.

4.3.2.Numerical simulation results

Fig.11 presents the vertical acceleration obtained from the numerical simulations and the experimental tests,and the latter is shown as a shaded zone representing the area between the limiting curves.The value of the first acceleration peak obtained from the simulations is in good correlation with the actual measurement of the maximum peak,with a value of 12,100 m/s2.The numerical curve is approximately in the middle of the boundaries of the shaded area.However,its length is smaller than the measured length.Unfortunately,the later portions of the curve do not match the experimentally determined trend.This divergence may be due to the factors mentioned earlier,including the heterogeneity of the dolomite and the presence of initial cracks and pores.Moreover,in the numerical simulations,the physical loss of rock continuity was not reproduced since erosion of elements was not introduced.The conditions at the boundaries of the copper pipe,dolomite and lead confinement were not fully known,which could further affect the observed discrepancies.The differences in the results can also be attributed to lead confinement shrinkage,which introduces slight triaxiality and is not included in the numerical simulations.These issues need to be further studied.

Fig.12 presents the blasting process and fractures of the dolomite.The fringe represents the damage index (D).Once it reachesthe critical value of 1,a fully damaged material is observed (the fracture surface is reached (see Fig.4).In the first stage,the detonation of the RDX causes the shock wave to travel into the specimen,resulting in compressive damage in the close vicinity of the blast hole and generating the crushing zone.Radial cracks form because the tensile stress exceeds tensile strength of the material.These cracking characteristics are generated within the perpendicular cross-section of the specimen along the axis of the detonation cord and blast hole.Spalling fracture at some distance from the top surface is also observed in the numerical model,which is the result of the stress wave reflection from the top free surface.Similar to the actual laboratory tests,higher cracking intensity can be observed on the top surface compared with that on the bottom surface of the specimen.

The blast-induced crack patterns observed at the top and bottom surfaces are presented in Fig.13.The photographs of the actual specimen were carefully analyzed,and the observed cracks were redrawn to visualize them more clearly.For the bottom surface,the numerical results match quite well with the actual outcomes.The compressive zones near the borehole are similar,with diameters of Dexp=23.3 mm and Dnum=22.8 mm from the experimental and numerical studies,respectively.However,in the numerical simulations,no transverse cracks were obtained due to the real-world properties of the specimen mentioned earlier,particularly the presence of initial cracks,heterogeneity,pores,and irregularities(see Fig.1).This also shows that dolomite study is difficult,as its fragility makes obtaining fully intact specimens for testing extremely challenging.

To validate the proposed model,the experimentally and numerically obtained radial cracks were counted and compared.Furthermore,the number of cracks was predicted based on fracture mechanics theory (Hajibagherpour et al.,2020) according to the following formula:

Fig.9.Vertical acceleration obtained from the experimental tests.

Fig.10.Exemplary dolomite specimen failure after the experimental small-scale blast test:(a) Radial cracks and the crushing zone on the top surface; (b) Scabbing failure on the bottom surface and (c) Spalling failure at approximately 30 mm from the top surface.

Fig.11.Vertical acceleration obtained from the experimental tests and numerical simulations.

where n is the number of radial cracks,R is the radius of the blast hole,c is the sound speed in the rock,KICis the rock fracture toughness and ˙ε is obtained from ˙ε=Pmax(1+ν)/(Etr)(Pmaxis the peak blast hole pressure,tris the time corresponding to Pmax).

The number of cracks calculated from Eq.(8) was compared with that in the middle cross-section of the numerical model and with that in the bottom surface in the experimental and numerical tests.Due to the intensive scabbing fracture of the rock in the top surface of the real specimen,it was extremely difficult to determine the number of cracks.Thus,the scabbing area was measured and compared with the numerical outcomes for this case.For this purpose,the area of fully damaged elements in the scabbing area of the top surface was measured(Fig.14) after blanking the elements in the outer layers of the top surface.The computationally determined scabbing failure was very close to the experimental measurements,with a discrepancy of less than 2%.

Fig.12.Failure process of the dolomite specimen obtained from numerical simulations.

For the bottom surface,the numerical simulations satisfactorily reproduced the number of radial cracks.Nevertheless,the simulations slightly underestimated the number of cracks,most likely due to the presence of initial cracks and inclusions in the dolomite specimen.The model also showed good correlation with fracture mechanics theory,and a similar number of cracks were obtained for the surface in the middle of the specimen.However,the numerical model slightly overestimated the number of cracks,20 compared to 19 calculated using Eq.(8).

Fig.13.Comparison of blast-induced crack patterns obtained from (a) Experimental tests and (b) Numerical simulations.

The crack density,defined as the crack length to area ratio,was also calculated for the bottom surface in the laboratory tests and numerical simulations.The transversal cracks and pre-cracking of the specimen were not considered during the calculations.Consequently,three zones with identical widths were distinguished:zone 1 close to the borehole,zone 2 in the middle,and zone 3 close to the specimen edge.Similar analyses were conducted previously by others (Dehghan Banadaki,2010; Gharehdash et al.,2020;Hajibagherpour et al.,2020).The crack density on the top surface was only calculated for zone 3 due to the scabbing as mentioned previously.Table 6 summarizes the measured and calculated crack densities for both surfaces (numerical and experimental),which are presented along with the three zones in Fig.15.On the bottom surface,the crack densities estimated based on the numerical outcomes are slightly smaller than the experimentally determined values.When moving from the center to the surface edge,the difference between the actual and numerical outcomes increases.The results are closest for zone 1 in the vicinity of the borehole.The crack density calculated for zone 3 is satisfactory compared with the experimental data,with an underestimation in the numerical simulation.

Fig.14.Surface used for estimating the scabbing failure area,number of radial cracks and crack density: (a) Top surface with blanked face elements and (b) Cross-section in the middle of the specimen length.

Table 6 Comparison of the number of crack,crack density and scabbing failure obtained from the simulations and experiments.

Fig.15.Top and bottom surfaces for which the crack density was estimated with marked zones: (a) Experiment and (b) Numerical simulation.

5.Conclusions

Experimental and numerical studies of small-scale blast tests were presented in this paper.Based on the results,the following conclusions can be drawn:

(1) The JH-2 material model adopted to simulate the behavior of dolomite was further validated in terms of its ability to reproduce the cracking and fragmentation of rock,and satisfactory results were obtained.The investigations will be continued and expanded to study cracking and fragmentation in blast tests of different types of specimens.

(2) Although literature data are available for the RDX material,several constants of the JWL-EOS model were determined based on the analytical procedure and hole expansion of the lead material.

(3) The experimental small-scale blast tests revealed a strong influence of the dolomite structure on the observed material response and failure.However,it was possible to distinguish typical radial cracks in the bottom and top surfaces of the specimen and scabbing failure of the top part of the specimen.Spalling failure was also observed during tests.

(4) The overall numerical reproduction of the results was sufficiently effective according to quantitative and qualitative comparisons of the results.Nevertheless,the acceleration history recorded in the numerical simulations did not match the experimental characteristics in the latter parts of the curve.There may be several reasons for these discrepancies,including the heterogeneous structure of the dolomite or the lack of consideration of physical loss of continuity in the simulations.Satisfactory similarity in terms of the peak value and the beginning of the peak was observed in the first peak.

(5) The proposed model satisfactorily reproduced the failure of the dolomite specimen.The results were compared in terms of crack pattern,crack number,scabbing failure and crack density.The analytically calculated number of cracks was also in good agreement with the number of cracks in the middle surface of the specimen model.

In future work,we will examine the possibility of simulating the loss of continuity of dolomite due to failure using erosion criteria coupled with meshless methods and transformation of finite elements into SPH particles.In addition,we will implement the proposed method of modeling and the JH-2 constitutive model for blast tests of dolomite,including parallel cut-hole blasting.The present study provides a foundation for the design of simple laboratory setups for conducting controlled blast tests with rock specimens and a methodology for analyzing fracture and cracking by coupling experiments with simulations.Further optimization of the cut-hole pattern could increase the effectiveness of the blasting method and rock fragmentation and minimize the costs of the procedure in the mining industry.

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgments

This research was supported by the Interdisciplinary Center for Mathematical and Computational Modeling (ICM),University of Warsaw (Grant No.GA73-19).The article was written as part of the implementation of the Military University of Technology(Grant No.22-876).This support is gratefully acknowledged.The numerical models were prepared using Altair® HyperMesh software.