Experimental Study and Numerical Simulation Using Extended Finite Element Method (XFEM) Combined with Cohesive Zone Model (CZM), of Crack Growth Induced by Non-Explosive Expansive Material on Two Neighboring Circular Holes of A Gneiss Rock

Frank Ferry Kamgue Tiam^1,2*, Raïdandi Danwe^1,3, Noël Konaï¹, Lucien Meva’a¹

¹Laboratory of Civil and Mechanical Engineering, National Advanced School of Engineering, Yaoundé 1 University, Yaoundé, Cameroon.
²School of Geology and Mining Engineering, Ngaoundéré University, Meiganga, Cameroon.
³National Advanced School of Engineering of Maroua, Maroua, Cameroon.

DOI: 10.4236/ojapps.2020.1010042                

Abstract

The Non-explosive expansion material (NEEM) is a method more environmentally friendly than the harmful conventional rock fracturing techniques. However, it is slower and very costly. Thus, any means of economizing their use is very desirable. This paper investigates the crack growth between two neighboring holes of a gneiss rock internally pressurized by NEEM mixed with water with the aim to evaluate the influence of holes spacing (center-to-center distance), on the initiation and growth of cracks. Field experimental results reveal that crack starts earlier and grows faster with increasing ambient temperature. But when the ambient temperature is above 28°C, the NEEM is “blown out” of the holes. At these ambient temperatures, the surrounding rocks are hot and cannot dissipate efficiently the heat generated by the hydration reaction. The best filling time was found to be in the evening when the daily hot temperature has drooped. The time to first crack increases as hole diameter decreases. The 3D numerical modeling and simulation of crack growth between two neighboring holes internally pressurized by NEEM using ABAQUS (XFEM/CZM) software shows a good agreement with the theoretical and experimental results.

Keywords

Crack Growth, XFEM, Gneiss, Holes Spacing

1. Introduction

Rock fracturing process is one of the most important operations in quarry mining. Blasting is the most common method of rock splitting [1]. This method, that uses explosive energy to fragment rocks, generates shock waves and gas energy. These cause vibrations, rocks blocks projections, loud noise and huge dust. Although been considered to be the most economically viable rock fragmentation method [2], blasting tends thus to be quite harmful to the environment and the surrounding communities. Due to these environmental and social threats, new methods of rock splitting have been developed.

Nowadays, non-explosive rock fracturing methods are spreading in mining, especially in underground mining. Expansive cement, also known as Non-Explosive Expansive Material (NEEM), has been proved to be a safer, pollution free and silent rock splitting method [3]. In this method, rock fracturing is induced by drilling closely spaced circular holes, and filling them with a mixture of expansive cement and water at the adequate proportions (30% water by weight). The hydration reaction that follows creates a solid volume increase of the NEEM, thereby generating internal expansion pressure on the wall of the borehole, causing compressive stress in the radial direction and tensile stress in the tangential direction around each hole (Figure 1). Thus, fracture first occurs at the weakest (existing micro cracks or jointing) section along the inside surface of the hole [4], at a point where this surface intersects the hole.

When the stress intensity factor between two neighboring holes is equal to the rock fracture toughness, crack may be initiated [4]. Fracture toughness is one of the fundamental material properties that have been used in developing crack fragmentation models [3]. Rock fragmentation using NEEM predominantly causes tensile failure (mode I) due to the hoop stress developing around the walls of the filled holes [5].

As drawbacks, NEEM fracturing method is slower than blasting [6]. In fact, when the NEEM is mixed with water (at 28% - 33% water by weight) and poured

Figure 1. Rock fracturing mechanism between two neighboring holes filled with NEEM (Natanzi et al., 2016).

in the drilled holes, it takes a few hours (about 4 to 6 depending on the type of NEEM, the surrounding material, and ambient temperature) for the hydration peak to be reached [7]. Thus, it takes few hours for the internal pressure to be considerable and to start generating stress on the surrounding wall of the borehole. Pressure generated by the NEEM is the key factor of the rock fracturing process and needs to be well studied and characterized. Experimental investigation of the expansion pressure, done by using a steel thick walled cylinder subjected to inner pressure by NEEM, reveals that expansion pressure is nonlinearly related to three independent parameters: loading time, Young’s modulus and hole diameter [8] [9] [10]. But these researchers failed to consider the influence of ambient temperature on the pressure generation, and also to consider the influence of stress generated by the neighboring holes.

In 2016, Natanzi and Leafer studied the influence of ambient temperature on the expansion pressure developed in a filled bore hole. They experimentally demonstrated that, higher temperature leads to greater and earlier expansive pressure as well as solid volumetric expansion of the NEEM. They failed to develop a mathematical model to explain the evolution of expansion pressure with ambient temperature [7].

In 2018, Shang et al. studied the influence of stress created by neighboring holes loaded with NEEM on the fracturing process. They proposed a mathematical model and determined the optimum hole diameter and spacing, but no experimental study or numerical simulation were carried out [4].

Holes spacing is one of the important parameters related to rock fracturing with NEEM. This is because optimum spacing results in improving the fragmentation process and reducing the cost [3]. Among the investigations carried out during the last two decades to determine hole spacing, Arshadnejad et al. proposed a more accurate expression that considers the tensile stress (σ_t), the rock fracture toughness (K_1c) of mode 1, the hole diameter (d), and the expansion pressure (P) as shown in Equation (1) [3] [9].

$S=\left\{-0.0888{\left(\frac{P}{{\sigma }_t}\right)}^2+1.0824\left(\frac{P}{{\sigma }_t}\right)-2.1583\right\}\frac{P^2{d}^2}{K_{Ic}^2}$ (1)

Even though this model considers more parameters to calculate the holes spacing, it fails to consider the influence of the confining stress acting on the rock mass and also the evolution of expansion pressure (P) with loading time as predicted by many studies. According to Gholinejad et al. in 2012, the expansion pressure (P) is related to three parameters: time, holes diameter and the rock fracturing toughness of mode I as in Equation (2) [4],

$P=0.12r_0^{0.407}{t}^{0.933}{\left(37K_{Ic}-11.6\right)}^{0.493}$ (2)

In recent years, many numerical simulations have focused on the investigation of stress, crack initiation and propagation around NEEM loaded holes [7] [8] [9]. Numerical simulation tools not only allow for an accurate evaluation of local stress for a complex loading history but also act as a guide toward obtaining an optimal choice of parameters [7]. During the last decade, Finite Elements Methods (FEM) had become a very important numerical simulation method mostly used by engineers and researchers. However, despite its efficiency to model and simulate many engineering problems, FEM is unable to solve problems with discontinuities and cracks propagations [11] [12]. The main drawback is that the mesh must be updated at each propagation step. This task is computationally costly for very complex geometry. Moës et al. developed in 1999 a new approach where crack propagation is not meshed: The Extended Finite Element Method (XFEM) [11] [12].

Numerically, using cohesive zone model for brittle material with an assumption of some plasticity is found to be a good approach to predict the crack growth in rocks [13] [14]. Popularly used for fracture simulation in brittle materials, the cohesive zone model (CZM) uses traction-separation law. The traction separation law represents the material damage zone in front of the crack tip where the material elements are pulled apart.

In 2014, Zolenka et al. used combined pressure deformation cohesive zone Method (CZM) and Extended Finite Element Method (XFEM) to model and simulate cracks and fracture propagation during hydraulic fracturing of rocks [15]. The ABAQUS numerical solutions achieved from the couple modeling technique (CZM and XFEM) were compared with analytical solutions and they found that they accurately match up with the analytical solutions and converge as the mesh is refined. Two years later, G. Sivakumar et al. also combined XFEM approach with cohesive zone Method (CZM), to model crack growth in rocks during uniaxial compression test using finite element-based software ABAQUS [13]. Their results were compared to the experimental results achieved in laboratory test and those available in the literature on crack growth in rocks during uniaxial compression test. They concluded that numerical model predicting crack propagation using ABAQUS shows good agreements with theoretical and experimental results.

This paper combines the XFEM approach with cohesive zone model (CZM) to analyze the crack growth between two neighboring holes loaded with NEEM of a Gneiss hard brittle rock using ABAQUS software with the objective of determining the optimum spacing in order to minimize the rock fracturing project cost. The numerical results are compared with those achieved experimentally and with theoretical results available in the literature.

2. Experimental and Numerical Method

Experimental field works were carried out in a Gneiss quarry at “Nkolondom” in the Centre region of Cameroon during the colder season when the ambient temperature varies from 16˚C to 24˚C during the day and 14˚C to 19˚C during the night. During this period, the ambient temperature varies very little during the three sections of the day: the morning (from 5.30 am to 10.30 am); in the afternoon (from 12:30 p.m. to 3 p.m.); and in the evening (6.30 p.m. to 9 p.m.). Thus, we assumed a constant temperature during each of the slices of the day by considering the average temperature. Natanzi et al. (2016), concluded that hydration peak is reached after 4 to 6 hours after filling the predrilled holes with NEEM. We will therefore fill the holes at the launch of a day section as sliced previously, this will allow us to reach the hydration peak being in the same slice, before the ambient temperature varies considerably.

Gneiss is the predominant metamorphic rock of this area and is usually used as building stones and as ornamental stones for floor and wall tiling. Holes were drilled on an outcrop and rock blocks with a pneumatic rock drill with three drill bits: 30 mm, 40 mm and 50 mm. Measuring tools (electronic caliper, electronic meter, chronometer…) were used for the measurement of crack lengths and widths. Figure 2 presents the outcrop and two blocs of gneiss rocks of irregular form that were cracked; Block 1: length 1000 mm, width 675 mm and height 750 mm and Block 2: length 1200 mm, width 950 mm and height 775 mm. the drilling depth was 650 mm.

In order to study the effects of anisotropy, drilled holes were spaced perpendicular and parallel to the foliation planes (weakness planes) of the gneiss rock.

The NEEM used during these experiments is CRACKAG that is manufactured in China. Expansive cement usually has a chemical composition as shown in Table 1.

NEEM is delivered from the manufacturer in a powder form and mixed with water to form slurry. The slurry is poured into drilled holes in rock. Water was poured in the powder expansive grout (in the adequate proportion) and mixed with a chemical electric liquid mixer. The mixing time was about five minutes to achieve a good homogeneity of the slurry. It was then poured in the drilledholes. Figure 3 present the different steps of the mixture process.

Figure 2. (a) Gneiss outcrop; (b) Gneiss blocks.

Table 1. Chemical composition of NEEM [2] [3].

3. Stress Distribution around Two Neighboring Circular Holes

Considering a volume element of the borehole located at a distance of r from the center of the first internally pressurized hole. The stress distribution (due to the internal pressure) is as shown in Figure 4, where σ_r is the radial stress, σ_θ the orthoradial stress, and τ_rθ the shear stress. The circular holes, internally pressurized by NEEM are of the same radius r₀ with the spacing between the two holes of S.

Based on the Theory of elasticity [4] [16], the equation of equilibrium of the volume element is given by

$\begin{array}{l}\left({\sigma }_r+\frac{\partial {\sigma }_r}{\partial r}\text{d}r\right)\left(r+\text{d}r\right)\text{d}\theta -{\sigma }_{r}r\text{d}\theta -\left({\sigma }_{\theta }+\frac{\partial {\sigma }_{\theta }}{\partial \theta }\text{d}\theta \right)\text{d}r\sin \left(\frac{\text{d}\theta }{2}\right)\\ -{\sigma }_{\theta }\text{d}r\sin \left(\frac{\text{d}\theta }{2}\right)+\left({\tau }_{r\theta }+\frac{\partial {\tau }_{r\theta }}{\partial \theta }\text{d}\theta \right)\text{d}r\cos \left(\frac{\text{d}\theta }{2}\right)-{\tau }_{r\theta }\text{d}r\cos \left(\frac{\text{d}\theta }{2}\right)=0\end{array}$ (3)

From Equation (3), we can obtain the following differential equations

$\{\begin{array}{l}\frac{\partial {\sigma }_r}{\partial r}+\frac{1}{r}\frac{\partial {\tau }_{r\theta }}{\partial \theta }+\frac{{\sigma }_r-{\sigma }_{\theta }}{r}=0\\ \frac{1}{r}\frac{\partial {\sigma }_{\theta }}{\partial \theta }+\frac{\partial {\tau }_{r\theta }}{\partial r}+\frac{2{\tau }_{r\theta }}{r}=0\end{array}$ (4)

Figure 3. Steps of the NEEM mixture process. (a) Adding water in NEEM powder; (b) Mixing; (c) Pouring into predrilled holes; (d) Holes filled with NEEM.

Figure 4. Stress distribution on an element of volume around the borehole filled with NEEM.

When the borehole is internally pressurized, the internal pressure is uniform on the hole walls, thus the radial deformation is uniform ( ${\tau }_{r\theta }=0$ ), and Equation (6) becomes

$\{\begin{array}{l}\frac{\partial {\sigma }_r}{\partial r}+\frac{{\sigma }_r-{\sigma }_{\theta }}{r}=0\\ \frac{\partial {\sigma }_{\theta }}{\partial \theta }=0\end{array}$ (5)

The solutions of Equation (5) is in the form of

${\sigma }_r=C{r}^n$, (6)

where C and n are constants.

The boundary conditions are as follows:

$\{\begin{array}{l}{\sigma }_r=P\text{when}r=r_0\\ {\sigma }_r=0\text{when}r=S\end{array}$ (7)

Equations (5) and (6) then give the following solutions:

${\sigma }_r=\frac{r_0^{2}P\left(1-\frac{S^2}{r}\right)}{S^2-r_0^2},{\sigma }_{\theta }=\frac{r_0^{2}P\left(1+\frac{S^2}{r}\right)}{S^2-r_0^2}$ (8)

With P the internal pressure generated by the NEEM, S the spacing between the two holes, r₀ the radius of the holes.

Substituting Equation (2) into Equation (8), we obtain:

$\begin{array}{l}{\sigma }_r=\frac{r_0^2\left(0.12r_0^{0.407}{t}^{0.933}{\left(37K_{Ic}-11.6\right)}^{0.493}\right)\left(1-\frac{S^2}{r}\right)}{S^2-r_0^2};\\ {\sigma }_{\theta }=\frac{r_0^2\left(0.12r_0^{0.407}{t}^{0.933}{\left(37K_{Ic}-11.6\right)}^{0.493}\right)\left(1+\frac{S^2}{r}\right)}{S^2-r_0^2}\end{array}$ (9)

According to Figure 3, the principal tensile stress is

${\sigma }_{yy}={\sigma }_r{\sin }^2\theta +{\sigma }_{\theta }{\cos }^2\theta$ (10)

4. Numerical Method

4.1. Cohesive Zone Models (CZM)

The viewpoint from which cohesive zone models originate regards fracture as a gradual phenomenon in which separation takes place across an extended crack tip, or cohesive zone, and is resisted by cohesive tractions [17]. Thus, cohesive zone elements do not represent any physical material, but describe the cohesive forces which occur when material elements (such as grains) are being pulled apart. Therefore, cohesive zone elements are placed between continuum (bulk) elements, as shown in Figure 5.

The idea of CZM is based on the assumption that the material’s failure process during fracture is limited to a narrow band in front of the main crack (Kuna et al.

Figure 5. Application of cohesive zone elements along bulk element boundaries.

2013). In CZM, the material follows the traction-separation law used for defining the shear traction and crack sliding displacement relations across the crack tip. Before the first principal stress reaches the tensile strength, the material behaves linearly elastic. As soon as the tensile strength is exceeded, the material begins to fail and the crack will get initiated. Crack initiation refers to the beginning of degradation of the cohesive response at an enriched element. The process of degradation begins when the stresses or the strains satisfy specified crack initiation criteria [18]. The degradation of the material is controlled by damage evolution law which describes the rate at which the cohesive stiffness will be degraded once the corresponding initiation criterion is reached.

When crack grows, the cohesive zone elements assigned in the mesh opens to simulate crack initiation. Since the crack path only follows the cohesive zone elements, crack propagation strongly depends on the mesh of the cohesive zone elements. This leads to the inclusion of Extended Finite Element Method (XFEM) approach, where the crack geometry is overlapped over the crack domain and their propagation happens without depending on the mesh.

4.2. Extended Finite Element Method (X-FEM)

The extended finite element method (XFEM) is a numerical technique which extends the classical finite element method approach focusing on crack propagation problems. The main idea behind this method is to deal with simple meshes and to take into account discontinuous displacements inside a finite element. Extended finite element methods (XFEM) allows simulation of crack growth without re meshing [13] [15]. The displacement approximation is enriched with discontinuous functions (Heaviside function).

Let’s consider 𝑥, a point in a finite element that is intersected by a crack. To calculate the displacement at point 𝑥 located within the domain, the approximation for a displacement vector function is,

$U\left(x\right)=\underset{\text{Classical FEM approximation}}{\underset{︸}{{\displaystyle \underset{i\in N}{\sum }{N}_i\left(x\right)u_i}}}+\underset{\text{Discontinuous enrichment}}{\underset{︸}{{\displaystyle \underset{i\in N_d}{\sum }{N}_i\left(x\right)H\left(x\right)a_i}}}+\underset{\text{Crack tip enrichment}}{\underset{︸}{{\displaystyle \underset{i\in N_p}{\sum }{N}_i\left(x\right)\left({\displaystyle \underset{j=1}{\overset{4}{\sum }}{F}_j\left(x\right)b_i^j}\right)}}}$ (11)

where N is the total nodes of the mesh,

$N_d\in N$ are all the enriched nodes by the discontinuity,

$N_p\in N$ are all the nodes near the crack tip,

N_i(x) are the usual nodal shape functions,

u_i are the displacement degrees of freedom at node i,

a_i are the nodal enriched degree of freedom vector,

H(x) is the Heaviside function associated to the discontinuous jump,

$b_i^j$ are the nodal enriched degree of freedom vector and

F_j(x) are the elastic asymptotic crack-tip functions and are given by the following expression:

$\left(F_j\left(x\right)\right)=\left(\sqrt{r}\sin \left(\frac{\theta }{2}\right);\sqrt{r}\cos \left(\frac{\theta }{2}\right);\sqrt{r}\sin \left(\frac{\theta }{2}\right)\sin \left(\theta \right);\sqrt{r}\cos \left(\frac{\theta }{2}\right)\sin \left(\theta \right)\right)$ (12)

(r, θ) are the polar co-ordinates related to the local axis of the crack tip and can be expressed in terms of the level sets as follows:

$r=\sqrt{ls{n}^2+ls{t}^2},\theta =\arctan \left(\frac{lst}{lsn}\right)$ (13)

These functions form the basis of the asymptotic field 1/r around the crack tip, and introduce additional degrees of freedom in each node, improving the solution accuracy near the crack tip. The first function $\sqrt{r}\sin \left(\theta /2\right)$ is discontinuous along the crack surfaces, giving the effect of required discontinuity in the approximation along the crack, as seen in Figure 6 and Figure 7.

With the use of the above mentioned near-tip enrichment functions an element partially cut by the crack can be modeled (Ahmed A., 2009).

Full XFEM enrichment is used only for the simulation of stationary cracks. The Near-tip asymptotic singularity is not considered for crack growth numerical analysis. Thus, only the displacement jump across a cracked element is considered in this study.

Figure 8 displays the meshing of the discontinuities and the tip crack by the XFEM method.

Figure 6. Near-tip enrichment functions (Ahmed A., 2009).

Figure 7. Enrichment function modeling the crack in a partially cut tip element (Ahmed A., 2009).

Figure 8. Crack propagation model with X-FEM.

4.3. Rock Numerical Modeling in ABAQUS Software

Equation (1) which is useful to determine the holes spacing is very complex and nonlinear. This is because the expansive pressure is one of the parameters of the equation, and it is well known that this pressure is dependent on parameters such as time, rock fracture toughness and holes radius as revealed by Equation (2), but also by ambient temperature [4]. The numerical model and simulation can be helpful to predict crack propagation on rocks and evaluate the optimum parameters for better field works. In this study, we used the extended finite element tool ABAQUS/Explicit to simulate the crack grow in a gneiss rock. ABAQUS add-ins, are written integrally with FORTRAN languages for calculation parts, in C++ languages for graphic display parts and in Python for scripts and parameterizations [14] [19]. Figure 9 displays the Solving Algorithm via ABAQUS 6.14 using XFEM.

The Gneiss rock parameters used to describe CZM in ABAQUS are: Young’s modulus E, the Poisson’s ratio ν, the shear modulus G, the rock fracture toughness of mode 1, the tensile strength σ_t and the rock density.

In Table 2 are the gneiss Rock parameters (properties) used in the numerical model [20] [21].

Damage initiation criterion used for this study is the Maximum principal stress criterion (MAXPS).Until the crack gets initiated, the material adopts the elastic properties. Once the material strength reaches its material limit, it will behave based on traction-separation law. Crack initiation occurs when the maximum

Figure 9. Solving Algorithm via ABAQUS 6.14 using XFEM.

Table 2. Rock parameters (properties) used in the numerical model.

principal stress reaches critical value. In regards to various calibrations that have been noted in different literatures (Chen, 2013) [22], calibration of the two main CZM parameters were made: The Cohesive Strength T₀ (yield stress (σ_𝑦)) and the Cohesive Energy Γ₀ (fracture energy (J_𝐼𝐶)).

The drilled hole was internally pressurized conferring to Equation (2). Crack growth simulation was carried out on the outcrop. Thus, the initial conditions are: three degree of freedom as shown in Figure 10(a) and in situ pressures are neglected. The meshing of the system is triangular and the mesh density is more concentrated around the holes as shown in Figure 10(b).

5. Results and Discussion

1) Influence of ambient temperature

Drilled Holes on the outcrop of diameter 50 mm were filled during three different hours of the day: In the morning with an ambient temperature of 20˚C, at midday with an ambient temperature of 28˚C, at the evening with an ambient temperature of 22˚C. The mixture temperature was 7˚C (NEEM was mixed with cold water). The experiments were carried out several times and the following results were achieved. Although the manufacturer recommended that the expansive cement could be used at temperatures of 25˚C to 40˚C, the NEEM was “blowout” of the holes filled at midday (ambient temperature of 28˚C). Figure 11 shows the crack evolution with time and ambient temperature. It shows that

Figure 10. (a) Numerical rock modeling in ABAQUS explicit with three degree of freedom; (b) meshing of the system in ABAQUS explicit.

Figure 11. Crack growth with time and ambient temperature.

the time for the first crack to appear was 8 hours for the morning filled holes and 10 hours 30 minutes for the evening filled holes.

From this figure, it appears that crack grows faster with increasing ambient temperature. Thus, the best filling time is in the evening when the daily hot temperature has dropped. Though the cracks started latter than of the holes filled in the morning, the propagations were faster, this may be because usually the ambient temperature drops in the evening, and night are colder.

The best time to first crack was achieved for holes filled at ambient temperature of 20˚C and at 22˚C, the “blown out» phenomenon did not happen as predicted by Natanzi et al. [4], who concluded that when the ambient temperature is above 21˚C, the NEEM will be blown out of holes several hours after filling (from 3 to 5 hours later). This contradiction may be caused by the fact that Natanzi’s experimental works were carried out in a thick-walled steel pipe, thus there were no surrounding rocks to dissipate the heat during the hydration reaction as in this paper. Hydration reaction is an exothermic reaction as shown in Equation (14).

$\text{CaO}+{\text{H}}_{\text{2}}\text{O}\to \text{Ca}{\left(\text{OH}\right)}_{\text{2}}+15.2\left(\text{kcal}/\text{mol}\right)\uparrow$ (14)

Nonetheless, the “blown out” phenomenon occurred after 4 hours for holes filled during midday when the ambient temperature was 28˚C, this may be because the surrounding rocks were hot (due to sun heating from morning to midday), and could not dissipate efficiently the heat during hydration reaction [3].

2) Influence of holes diameters

Holes of diameter 30 mm, 40 mm and 50 mm were drilled and filled on the outcrop. Figure 12 displays the crack growth with time per drilled holes diameter.

From Figure 12, it appears that the time to first crack are 8 hours 40 minutes for the 50 mm holes, 10 hours for the 40 mm holes and 10 hours 50 minutes for the 30 mm holes. Thus, the time to first crack increases as the holes diameter decreases. This result is the same as those obtained by previous experimental study on thick-walled steel pipe, aluminum pipes and concrete blocks [3] [7]. Figure 12 also reveals that the higher the holes diameter, the faster the crack grows. But higher diameter means more NEEM to be filled in the hole and this influences the cost of the rock splitting project.

3) Influence of rock Anisotropy.

Gneiss is an anisotropic hard brittle rock. Holes of diameters 40 mm and 50 mm were drilled on the two blocks along the foliation lines. Figure 13 shows the

Figure 12. Crack evolution with holes diameter and time.

Figure 13. Cracks growth on foliation planes.

crack growth with time for each borehole diameter.

Figure 13 reveals that the time to first crack is 5 hours 30 minutes for the 50 mm holes and 6 hours 25 minutes for the 40 mm holes. Comparing these results with those obtained on the outcrop for the corresponding diameter (time to first crack were 8 hours 40 minutes for the 50 mm holes and 10 hours for the 40 mm holes on the outcrop), it appears that crack starts earlier and grows faster when the holes are drilled along the foliation line. In order to achieve fast rock splitting, it is then advisable to drill holes parallel to the weakness planes (foliation planes) of the rocks.

4) Numerical results

Aiming for a comparative analysis between the experimental results and the simulation, holes of diameters 30, 40 and 50 mm were drilled on the cohesive zone. Spacing S used during the study varies for each drilled hole from 120 to 250 mm. Figure 14(a) and Figure 14(b) represents a crack propagation between two neighboring holes after time t_a and t_b in hours respectively, with t_a < t_b. The red color in the figure means high stress concentration.

Figure 15 presents the numerical results of crack growth with time per diameter. It appears that the numerical crack propagation achieved by simulation with the coupled model of XFEM and CZM has the same evolution shape like the results obtained experimentally as shown in Figure 12. Figure 15 thus, reveals that cracks initiate earlier and grow faster for bigger holes diameters.

Figures 16-18 display the comparative study between experimental and numerical simulation results. They reveal that the numerical solution converges with experimental field results. These figures also illustrate that cracks initiate very much earlier numerically than experimentally (time to first crack for a hole of diameter 40 mm is 10 hours experimentally and 1 hour 30 min for numerical simulation). This is because ABAQUS software does not consider the hydration time of the expansive cement. In fact, when the NEEM is mixed with water and poured in the drilled holes, it takes few hours (about 4 to 6 depending on the type of NEEM powder, surrounding material and the ambient temperature) for

Figure 14. (a) Crack propagation between two neighboring holes after t_a (in hours); (b) Crack propagation between two neighboring holes after t_b (in hours).

Figure 15. Numerical results of crack growth with time and hole diameter.

Figure 16. Experimental and numerical crack growth with time for holes of diameter 30.

Figure 17. Experimental and numerical crack growth with time for holes of diameter 40.

Figure 18. Experimental and numerical crack growth with time for holes of diameter 50.

the hydration peak to be reached [4]. Thus, it takes few hours for the internal pressure to be considerable and to start generating stress on the surrounding wall of the borehole.

The hydration process depends on the type of NEEM powder (percentage of CaO) and ambient temperature (Natanzi et al. 2016). However, a quantitative relation/equation between ambient temperature, NEEM type and the performance of NEEM in terms of expansive pressure, is still not available to our knowledge. Such further work needs to be performed in this regard to further improve the prediction performance of Equation (2) (Shang et al. in 2018), and reproduce the nature of the expansive grout, that takes several hours for the internal pressure to become considerable.

Though the in-situ pressures were neglected, rock modeling and crack numerical simulation with XFEM method exhibit some similarities with the experimental field results. Figure 19 presents the crack propagation between two neighboring holes achieved by experimental field works (Figure 19(a)), and by numerical simulation with CZM approach combined with XFEM method (Figure 19(b)).

5) Optimum hole spacing determination

From the simulation of crack growth as shown in Figure 14(a) and Figure 14(b), it is possible to determine the stress distribution around the filled borehole and thus study its evolution with time. In this paper, tensile stress evolution was evaluated at the midpoint between the center points of two neighboring holes: ( $r=S/2;\theta =0$ ).

Figures 20-22 illustrate the evolution of the tensile stress at the midpoint for holes of diameter 30 mm, 40 mm and 50 mm respectively. The time after NEEM loading considered is: 4, 8, 12, 16, 20, and 25 h.

When the maximum principal stress reaches critical value (greater than the tensile strength of the rock), crack initiation occurs and propagate. Figures 20-22,

Figure 19. Experimental and numerical crack propagation between two neighboring hoes.

Figure 20. Tensile stress evolution at midpoint for hole of diameter 30 mm.

Figure 21. Tensile stress evolution at midpoint for hole of diameter 40 mm.

Figure 22. Tensile stress evolution at midpoint for hole of diameter 50 mm.

thus give us the optimum spacing for a particular diameter, and also reveal the corresponding fragmentation time. For example, Figure 22 reveal that for a diameter of 50 mm, when the spacing is equal to 146.5 mm, fragmentation time is equal to 25 hours.

From Figures 20-22, optimum holes spacing evolution with fragmentation time and radius can be deduced as displayed in Figure 23.

From Figure 23, it is then possible, using the polynomial regression, to determine the optimal hole spacing S_r that corresponds to a particular radius and the desired fragmentation time as given by Equation (15).

Figure 23. Evolution of Spacing with fragmentation time.

$\{\begin{array}{l}{S}_{15}=-0.0833t^2+4.25t+75.833\\ S_{20}=-0.001t^2+1.0068t+152.22\\ S_{25}=-0.0017t^2+1.2173t+195.61\end{array}$ (15)

Hole spacing is an important parameter which influences the rock fragmentation process, in fact when the holes are too close (small hole spacing), cracks will occur and grow rapidly but the project cost will be very high because more holes will need to be drilled. In the other hand, large spacing will result in few holes to drill, but crack growth will be delayed or may not even happen. Optimal spacing as given by Equating (15) therefore helps to predict the fragmentation time for a radius at the adequate spacing. This study focused only on three diameters, but more simulation with XFEM can be done to obtain S_r relation with fragmentation time for other desired radii.

Figures 20-22 are similar to those obtained at the same position ( $r=S/2;\theta =0$ ) with analytical method by Shang et al. in 2018. The analytical equation used by Shang et al. to determine the maximum principal stress is very complex and is nonlinear when the azimuth angle, θ is not equal to zero. With this coupled simulation method (CZM and XFEM) using ABAQUS software, it is therefore possible to determine the principal stress at any point around the internally pressurized holes with the desired azimuth.

6. Conclusions

In this paper, experimental field works and numerical simulation were carried out to investigate the crack growth between two neighboring holes of a gneiss rock internally pressurized by NEEM. Experimental results reveal that ambient temperature, hole diameter, rock anisotropy and holes spacing influence the crack growth. It appears that crack grows faster with increasing temperature, but when the ambient temperature is above 22˚C, NEEM will be blown out of the holes. The numerical models simulated with ABAQUS software (XFEM coupled with CZM) of crack growth between two neighboring holes internally pressurized by NEEM show a good agreement with the theoretical/analytical and experimental results.

Further study will be carried out to investigate on the heat dissipation on surrounding rocks during the hydration reaction.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Author Contributions

The initial draft of the article was done by Frank Ferry KAMGUE TIAM. The drafted manuscript was reviewed and corrected by Noël KONAÏ and Lucien MEVA’A. The concept was conceived by Raïdandi DANWE and he also carried out the final check of the paper.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Saharan, M.R. and Mitri, H.S. (2008) Numerical Procedure for Dynamic Simulation of Discrete Fractures Due to Blasting. Rock Mechanics and Rock Engineering, 41, 641-670. https://doi.org/10.1007/s00603-007-0136-9
[2]	Huynh, M.P. and Laefer, D.F. (2009) Expansive Cements and Soundless Chemical Demolition Agents: State of Technology Review. Proceedings of 11th Conference on Science and Technology, Ho Chi Minh City, 21-23 October 2009, 209-214.
[3]	De Silva, R., Ranjith, P.G. and Mandadige, A.P. (2016) An Alternative to Conventional Rock Fragmentation Methods Using SCDA: A Review. Energies, 9, 958. https://doi.org/10.3390/en9110958
[4]	Shang, J., Zhao, Z. and Aliyu, M.M. (2018) Stresses Induced by a Demolition Agent in Non-Explosive Rock Fracturing. International Journal of Rock Mechanics and Mining Sciences, 107, 172-180. https://doi.org/10.1016/j.ijrmms.2018.04.049
[5]	Tang, S.B., et al. (2017) Study of the Fracture Process in Heterogeneous Materials around Boreholes Filled with Expansion Cement. International Journal of Solids and Structures, 112, 1-15. https://doi.org/10.1016/j.ijsolstr.2017.03.002
[6]	Arshadnejad, S., Goshtasbi, K. and Aghazadeh, J. (2009) Stress Concentration Analysis between Two Neighboring Circular Holes under Internal Pressure of a Non-Explosive Expansion Material. Yerbilimleri, 30, 259-270.
[7]	Natanzi, A.S., Laefer, D.F., et al. (2016) Cold and Moderate Ambient Temperatures Effects on Expansive Pressure Development in Soundless Chemical Demolition Agent. Construction and Building Materials, 110, 117-127. https://doi.org/10.1016/j.conbuildmat.2016.02.016
[8]	Laefer, D.F., Ceribasi, S., Wortman, J., Abrozevitch-Cooper, N., Huynh, M.P. and Midgette, J. (2010) Expansive Fracture Agent Behaviour for Concrete Cracking. Magazine of Concrete Research, 62, 443-452. https://doi.org/10.1680/macr.2010.62.6.443
[9]	Arshadnejad, S., Goshtasbi, K. and Aghazadeh, J. (2011) A Model to Determine Hole Spacing in the Rock Fracture Process by Non-Explosive Expansion Material. International Journal of Minerals Metallurgy and Materials, 18, 509-514. https://doi.org/10.1007/s12613-011-0470-5
[10]	El Dessouki, S. and Mitri, A.H. (2011) Rock Breakage Using Expansive Cement. Engineering, 3, 168. https://doi.org/10.4236/eng.2011.32020
[11]	Moes, N. and Belytschko, T. (2002) Extended Finite Element Method for Cohesive Crack Growth. Engineering Fracture Mechanics, 69, 831-833. https://doi.org/10.1016/S0013-7944(01)00128-X
[12]	Dolbow, J., Moës, N. and Belytschko, T. (2001) An Extended Finite Element Method for Modeling Crack Growth with Frictional Contact. Computer Methods in Applied Mechanics and Engineering, 190, 6825-6846. https://doi.org/10.1016/S0045-7825(01)00260-2
[13]	Sivakumar, G. and Maji, V.B. (2016) Simulation of Crack Propagation in Rocks by XFEM. Recent Advances in Rock Engineering (RARE 2016), Bengaluru, 16-18 November 2016. https://doi.org/10.2991/rare-16.2016.46
[14]	Haddad, M. and Sepehrnoori, K. (2016) XFEM-Based CZM for the Simulation of 3D Multiple-Cluster Hydraulic Fracturing in Quasi-Brittle Shale Formations. Rock Mechanics and Rock Engineering, 49, 4731-4748. https://doi.org/10.1007/s00603-016-1057-2
[15]	Zielonka, M., et al. (2014) Development and Validation of Fully-Coupled Hydraulic Fracturing Simulation Capabilities. 2014 SIMULIA Community Conference, SCC2014, Providence, 19-21 May 2014.
[16]	Murakami, Y. (2016) Theory of Elasticity and Stress Concentration. Wiley, West Sussex, 445. https://doi.org/10.1002/9781119274063
[17]	Blal, N., Daridon, L., Monerie, Y. and Pagano, S. (2011) On Mesh Size to Cohesive Zone Parameters Relationships. ACOMEN 2011, Liège, November 2011.
[18]	Daridon, L., Monerie, Y., Pagano, S. and Blal, N. (2012) Artificial Compliance Inherent to the Intrinsic Cohesive Zone Models: Criteria and Application to Planar meshes. International Journal of Fracture, 178, 71-83.
[19]	Giner, E., Sukumar, N., Tarancón, J. and Fuenmayor, F. (2008). An Abaqus Implementation of the Extended Finite Element Method. Engineering Fracture Mechanics, 76, 347-368. https://doi.org/10.1016/j.engfracmech.2008.10.015
[20]	Gasc-Barbier, M. and Wassermann, J. (2013) Etude de l’anisotropie des roches par la méthode ultrasonique—Application au gneiss de Valabres. Bulletin des Laboratories des Ponts et Chaussees, No. 280-281, 139-153.
[21]	Kawa, Z. and Alshkane, Y. (2018) Physical and Mechanical Properties of Metamorphic Rocks. Journal of Garmian University, 5, 194-209. https://doi.org/10.24271/garmian.334
[22]	Chen, X., Deng, X.M. and Sutton, M.A. (2013) Simulation of Stable Tearing Crack Growth Events Using the Cohesive Zone Model Approach. Engineering Fracture Mechanics, 99, 223-238. https://doi.org/10.1016/j.engfracmech.2012.12.017