Arab J Geosci (2016) 9: 145 DOI 10.1007/s12517-015-2196-6

ORIGINAL PAPER

Crack propagation in rock specimen under compressive loading using extended finite element method M. Eftekhari 1 & A. Baghbanan 1 & H. Hashemolhosseini 1

Received: 30 March 2015 / Accepted: 15 October 2015 / Published online: 24 February 2016 # Saudi Society for Geosciences 2016

Abstract Understanding of mechanisms of initiation and propagation of rock fracture is a key factor in many engineering problems in rock. Since most failures in rock occur under compression condition, in this study, the crack propagation in rock specimen with inclined crack under uniaxial compression has been investigated. In order to simulate the crack propagation under compressive loading, the maximum tangential stress (MTS) criterion was firstly improved and incorporated in developed extended finite element method (XFEM) code. The XFEM allows crack growing through the finite element mesh without re-meshing process. The developed code is validated in different bench mark tests under tension by MTS criterion which can precisely predict the direction of crack propagation under tension condition. Different specimens with variation of inclined cracks under uniaxial compression are examined using the developed code. The results of modeling show that the crack growth direction is deviated from the initial crack direction and propagated from the crack tip toward the loading plates. The results obtained by the XFEM are in good agreement with the experimental and numerical results in the literature. It shows that the developed improved stress-based criterion in this study could be effectively used in modeling crack initiation and propagation under both tension and compression conditions.

* M. Eftekhari [email protected]; [email protected] 1

Department of Mining Engineering, Isfahan University of Technology (IUT), Isfahan 84156-83111, Iran

Keywords Extended finite element method (XFEM) . Maximum tangential stress (MTS) criterion . Improved maximum tangential stress (IMTS) criterion . Crack propagation . Stress intensity factor

Introduction Generally, fracture mechanics is a science related to the occurrence of a crack and its propagation under applied external load. Since cracks may appear anywhere, the application domain of fracture mechanics is extensive. Discontinuity is a common formation phenomenon in rock mass. On the other hand, in brittle and semi-brittle material, the propagation of the crack is the main reason of failure in many cases. Therefore, the crack propagation in rock engineering is a great importance. The issue of crack propagation in rock is studied based on the theory of fracture mechanics and especially the theory of linear elastic fracture mechanics. Crack initiation and propagation mechanisms have been studied experimentally in laboratory (Bobet and Einstein 1998a; Brace and Bombolakis 1963; Jiefan et al. 1990; Li et al. 2005; Reyes 1991; Shen et al. 1995; Theocaris and Sakellariou 1991) or numerically (Bobet and Einstein 1998b; Chen et al. 1998; Ingraffea and Heuze 1980; Shen and Stephansson 1994; Wu and Wong 2012). The finite element method (FEM) is one of the most used methods, which has the ability of modeling any kind of boundary and geometry. The conventional finite element method can model the crack propagation by implementing singular elements to the crack tips. Nevertheless, it has limitations in modeling the crack propagation, such as the necessity of coinciding finite element mesh with the crack and also the need for re-meshing after each crack increment. Therefore, modeling of crack propagation for an arbitrary crack is almost impossible with

145 Page 2 of 10

FEM. This not only demands a large computation resource but also results in loss of accuracy while the data is mapped from old to the new meshes (Ahmed 2009). In order to overcome these limitations, novel methods based on classical ones have been developed and they are able to remove the limitations to an acceptable extent. One of the methods that benefits the advantages of the FEM and also has acceptably reduced the two recent problems is the extended finite element method (XFEM). This method is based on FEM and the partition of unity. The principle of this method is to solve the problems involving crack growth based on the enrichment of the polynomial approximation space of the classical finite element method. In this method, the crack is modeled independently from the finite element mesh due to the addition of enrichment functions to the approximation space and so increase of the nodal degrees of freedom. Enrichment in this method is only local and only subsets of nodes are enriched. In addition to having all the other important advantages of the FEM, the absence of dependence of the crack propagation process on the finite element mesh is the biggest and most important distinction of this method. Therefore, according to these advantages and capabilities of the XFEM, this numerical method is used in this study. In fracture mechanics, one of the main challenges is to specify the direction and path of the crack propagation. Although several different criteria have been proposed for determining direction of crack propagation under general mixed-mode loading conditions (Andrianopoulos and Theocaris 1985; Erdogan and Sih 1963; Hussain et al. 1974; Nuismer 1975; Palaniswamy 1972; Palaniswamy and Knauss 1972; Sih 1973; Sih 1974; Theocaris and Andrianopoulos 1982; Tirosh and Catz 1981), the most commonly used criteria are the maximum tangential stress (MTS) criterion (Erdogan and Sih 1963), maximum energy release rate criterion (Nuismer 1975), and also the minimum strain energy density (Sih 1974). However, the MTS criterion in which the crack is propagated from the crack tip along the direction of maximum tangential stress has been used more extensively in the XFEM for modeling crack growth. The direction of the crack in MTS criterion depends on the amount of the stress intensity factors. When a crack is located in a compressive stress field, the mode I stress intensity factor, KI, becomes negative. Since this criterion has been developed based on the tensile loading rather than the compressive, it holds good only in the positive region (crack opening) and cannot predict the direction of crack propagation in the negative zone (Al-Shayea 2005). Note that a negative value of mode I for closed crack does not have a correct physical sense in practice (Lin 2000). Since, unlike most engineering materials, in many

Arab J Geosci (2016) 9: 145

important engineering applications in rocks such as cutting tools, mining, and fragmentation, failure occurs through cracking under compression, the study about initiation and propagation of the crack under compressive condition becomes necessary. The main aim of this research work is to implement a criterion based on the condition of the stress near the crack tip in the XFEM to determine the crack propagation path which could be used in both tension and compression loading conditions. This is in line with the research work where the author started the development of an objectoriented code called MEX-FEM, based on XFEM, and employed in simulating fracture propagation in some specimens of rock such as cracked Brazilian disc (Eftekhari et al. 2014; Eftekhari et al. 2015a) and semi-circular bend specimen (Eftekhari et al. 2015b). This paper consists of three parts. Firstly, the basic logic of XFEM, developed numerical code, MEX-FEM, and also the criterion of crack propagation are described. Then in the next part, the suggested criterion has been validated in tension mode and mixed modes I and II. Finally, crack propagation path in specimen with crack under uniaxial compressive loading is determined.

Extended finite element method The XFEM was firstly developed by Belytschko and Black (1999) and further extended by Moës et al. (1999) and Dolbow (1999). It is based on FEM for modeling the crack propagation without any remeshing stages. In this method, the finite element approximation is enriched through the enrichment functions which are the modified Heaviside function, H, (Eq. 1) and the near-tip asymptotic functions (Eq. 2). All nodes of the elements that are completely cut by the crack and shown by a circle symbol in Fig. 1 are enriched by Heaviside function while the nodes of element containing the crack tip are enriched by near-tip asymptotic functions which are shown by a square symbol in Fig. 1.
H ðxÞ ¼



F j ðr; θÞ

þ1 −1

f or x > 0 f or x < 0


 j¼1;2;3;4

¼

ð1Þ

pffiffi pffiffi θ pffiffi θ pffiffi θ θ rsin ; rcos ; rsin sinθ; rcos sinθ 2 2 2 2



ð2Þ where r and θ are polar coordinates in the local crack-tip coordinate system. Among the near-tip asymptotic functions,

Arab J Geosci (2016) 9: 145

Page 3 of 10 145 a)

a ag

Crack tip ck

crack

cra

(Δ

e

nl

tio

Crack interior

th ng

rop

p

new crack tip

θc

crack tip

Fig. 3 Illustration of crack propagation direction

Fig. 1 Crack tip and crack interior enrichment nodes

pffiffi only the function rsin 2θ is discontinuous that takes into account the discontinuity across the crack face. The XFEM approximation of displacement field is explained as follows (Moës et al. 1999): uh ðxÞ ¼

X i∈I

N i ðxÞui þ

X

N j ðxÞH j ðxÞa j þ

j∈ J

X k∈K

N k ðx Þ

4 X

F l ðxÞblk

l¼1

ð3Þ where N(x)’s are the shape functions, ui are the nodal displacements (standard degrees of freedom), aj are vectors of additional degrees of nodal freedom associated with the Heaviside function, and blk are vectors of additional degrees of nodal freedom associated with the elastic asymptotic crack-tip functions. In Eq. 3, I is the set of all nodes in the mesh, J is the set of enriched nodes with discontinuous enrichment, and K is the set of nodes enriched with asymptotic enrichment. The equations and associated matrix in the XFEM are formed and solved as similar as FEM procedure. From the proposed approximation in Eq. 3, the standard discrete equation Kd = f is obtained, where f is the vector of external nodal

forces and K is the stiffness matrix. In order to simulate the crack propagation, a numerical code called MEX-FEM has been developed based on XFEM formulation with C++ programming language.

Crack propagation criterion The direction of crack propagation in MTS criterion is determined according to the amount of the SIF as follows: 0 1 sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  2 1 K K I I ð4Þ θc ¼ 2tan−1 @  þ 8A ; K II 4 K II where θc is the crack growth angle and KI and KII are stress intensity factors of modes I and II, respectively. If KII = 0, then

σ

Initial Crack

a

rd

W Crack tip

σ (a) Fig. 2 Effective elements in the interaction integral (the hachured elements)

(b)

Fig. 4 Plate with edge crack under tension. a Sschematic view. b Discretized finite element meshes

145 Page 4 of 10

Arab J Geosci (2016) 9: 145 Pa

Initial Crack Tip

where W(1,2) is the interaction strain energy and is formulated as follows: ð1Þ

ð2Þ

ð2Þ

ð1Þ

W ð1;2Þ ¼ σi j ϵi j ¼ σi j ϵ i j

Crack Propagation

ð7Þ

(1) (1) (2) (2) (2) (σ(1) ij , ϵij , ui ) and (σij , ϵij , ui ) are the actual and auxiliary state, respectively, of stress, strain, and displacement matrixes. The relationship between the mode I and II stress intensity factors and interaction integral is as follows (Moës et al. 1999):

2 ð1Þ ð2Þ ð1Þ ð2Þ I ð1;2Þ ¼ 0 K I K I þ K I I K II ; ð8Þ E

Fig. 5 Tangential stress concentration near the crack tip in plate with edge crack under tension

θc = 0 and pure extension, mode I occurs while when KII > 0 and KII < 0, the crack growth angles are determined as θc < 0 and θc > 0, respectively. A more efficient expression for θc is developed (Sukumar and Prévost 2003) as follows: 3 2 6 θc ¼ 2tan−1 4

−2K II =K I 7 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 5 2 1 þ 1 þ 8ðK II =K I Þ

ð5Þ

The value of SIF is calculated by the interaction integral, I (Eq. 6), on the elements around the crack tip. These elements are hachured in the Fig. 2. # Z " ð2Þ ð1Þ ∂q ð1Þ ∂ui ð2Þ ∂ui ð1;2Þ ð1;2Þ I ¼ σi j þ σi j −W δ1 j dA ð6Þ ∂x1 ∂x1 ∂x j A

Fig. 6 Comparison of crack propagation path obtained from MTS (dashed line) and IMTS criterion (bold line) in plate with edge crack under tension

where E′ is defined as follows: ( E 0 f or plane strain condition E ¼ 1−υ2 E f or plane stress condition

ð9Þ

Making suitable choice of state 2 as the pure mode I as(2) ymptotic fields with K(2) I = 1 and KII = 0 gives the SIF of mode I in terms of the interaction integral while if K(2) I = 0 and (2) KII = 1, mode II is obtained (Moës et al. 1999). When a specimen with a crack experiences uniaxial compression loading, the SIF mode I becomes negative in which the crack faces are penetrated into each other. Since there is not any physical sense for this phenomenon in a closed crack, the MTS criterion cannot determine the direction of the crack propagation correctly. This issue was addressed in the literature when Bobet (1997) compared the results of the initiation angle using MTS criterion with the results of the uniaxial compressive test in the laboratory and stated that MTS criterion does not provide a good prediction of closed crack. Also due to the same reason, Shen and Stephansson (1994) developed the so-called F criterion or the modified BG^

Arab J Geosci (2016) 9: 145

Page 5 of 10 145

σ

criterion to investigate the propagation of closed crack in compressive loading. It is notable that crack propagation under compressive conditions is caused by tensile stresses that act near the tips of pre-existing cracks (Sharafisafa and Nazem 2014). It should be noted that in all abovementioned methods, the deciding factor in crack propagation process is stress condition only near the crack tip not in crack surfaces. Therefore, in this paper, a criterion based on the condition of the stress near the crack tip is employed to determine the crack propagation direction. The logic of this criterion is similar to Bobet’s criterion for crack propagation in Boundary Element Method (Bobet and Einstein 1998b), and a tension crack is propagated from the tip of pre-existing crack in the direction of the maximum tensile tangential stress (as shown in Fig. 3) and follows: ∂σθ ∂2 σ θ ¼0; <0 ∂θ ∂θ2

a

α=60

L

ð10Þ

This criterion may be considered as improved MTS criterion, called IMTS, which the direction of crack propagation in this criterion is similar to the MTS criterion, while the initiation angle is determined regarding the values of stress near the crack tip and it is completely independent of the value of the SIF. In this case, SIF value does not compute by interaction integral procedure. In order to find the direction of maximum tangential stress, the stress state in near the crack tip is processed in two steps. Firstly, the values of the stresses in the numerical integration points near the crack tip are transformed from Cartesian system to polar coordinate system in order to calculate the value of tangential stress. Since the numerical integration points are irregularly located in near the crack tip, in the second step, the values of the tangential stress in the proposed points with regular intervals in different angles are calculated by interpolating process. By this way, the direction of the maximum tangential stress is identified.

Crack propagation under different loading paths using XFEM Since the MTS criterion precisely predicts the direction of the crack propagation under tension conditions, i.e., mode I and mixed mode I/II, therefore, to validate the proposed criterion in the XFEM, some bench mark tests are conducted and the obtained results are compared to the MTS criterion results.

W

σ Fig. 7 Schematic view of plate with angled center crack under tension

edges, as shown in Fig. 4a. The plate has a length L = 11 cm, width W = 5.4 cm, and crack length a = 1.5 cm. Plane strain conditions are assumed and the domain is partitioned with a uniform mesh of 18 × 36 4-node quadrilateral elements, as shown in Fig. 4b. In this case, the exact mode I stress intensity factor is given by the following: pffiffiffiffiffiffi K I ¼ Cσ aπ ð11Þ where C is the finite geometry correction factor (Richardson et al. 2011) and is calculated as follows: a

a 2 a 3 a 4 þ 10:55 C ¼ 1:12−0:231 −21:72 þ 30:39 W W W W

ð12Þ

Table 1

Test 1, plate with edge crack under tension The first example involves a straight edge crack in a rectangular plate which is loaded by tension over the top and bottom

Analytical XFEM

Analytical and computed values of the stress intensity factor KI (KPa.m0.5)

KII (KPa.m0.5)

0.52 0.55

0.90 0.92

145 Page 6 of 10

Arab J Geosci (2016) 9: 145 Pa

Initial Crack Tip

Crack Propagation

are shown in Fig. 6. It can be seen that the crack propagation path predicted by the IMTS criterion matches well to the MTS criterion. Test 2, plate with angled center crack under tension The second example has the same geometry of plate as the first example, but the crack is angled center crack (Fig. 7). The length and angle of crack (i.e., a and α) are 1.5 cm and 60°, respectively. In this geometry of crack, the values of the stress intensity factor are determined as follows (Richardson et al. 2011): pffiffiffiffiffiffi K I ¼ σ aπcos2 ðαÞ

ð13Þ

pffiffiffiffiffiffi K II ¼ σ aπsinðαÞcosðαÞ

ð14Þ

Fig. 8 Tangential stress concentration near the crack tip in plate with angled center crack under tension

The analytical and computed values of mode I stress intensity factor for this example are 4.6 and 4.5 KPa.m0.5 , respectively. The condition of the tangential stress near the crack tip and also the predicted direction of the crack propagation in this specimen are shown in Fig. 5. As can be seen in Fig. 5, the crack is propagated under pure mode I. The calculated initiation angle in the proposed and the MTS criteria, respectively, are 0.9° and 0.02° which are very small and close to zero. The crack propagation path using the IMTS criterion and also the MTS criterion

Fig. 9 Comparison of crack propagation path obtained from MTS (dashed line) and IMTS criterion (bold line) in plate with angled center crack under tension

The analytical and computed values of the stress intensity factor for this example are summarized in Table 1. Due to the crack geometry, the combination of modes I and II occurred. The condition of the tangential stress near the crack tip in polar coordinates in this specimen is shown in Fig. 8. The initiation angle in the suggested criterion and the MTS criterion are −60° and −59.7°, respectively. The crack propagation path predicted by the IMTS criterion and also the MTS criterion are shown in Fig. 9. There is a good agreement between the IMTS criterion and the MTS criterion in this example. Therefore, this criterion could be employed in the rock

Arab J Geosci (2016) 9: 145

Page 7 of 10 145

P

Ri

β

a

Ro

a

P Fig. 10 Geometry of a HCCD specimen

specimens containing cracks which are subjected to compressive loading.

Crack propagation in rock specimen under compressive loading In order to validate the IMTS criterion in the rock specimen under mixed mode I/II, Hollow Centre Cracked Disc (HCCD) specimen is employed. This specimen is one of the suggested specimens to determine the toughness of the rock (Shiryaev and Kotkis 1983). Schematic view and geometrical dimensions of HCCD specimen is presented in Fig. 10. As can be seen in Fig. 10, HCCD

Fig. 11 Comparison of crack propagation path obtained from MTS (dashed line) and IMTS criterion (bold line) in HCCD specimen with crack angle of 20°

is a disc with a radius of Ro in which a central hole with radius of Ri is drilled. Two straight central cracks with length of a are created from the surface of the hole. Simple geometry, easy test setup procedure, and convenient specimen preparation are the advantages of this specimen for determining fracture toughness (Amrollahi et al. 2011). In addition, it is possible to measure pure mode I, pure mode II, and mixed mode I–II fracture toughness by changing the pre-exiting crack angle with respect to the diametrical loading direction. The crack propagation path in the HCCD specimen with crack angle of 20° using the IMTS criterion and the MTS criterion are shown in Fig. 11. It can be concluded that this criterion has the ability to predict the crack propagation. The IMTS criterion is employed in the cracked specimen under uniaxial compression. The size of the specimen is similar to the tensile examples (i.e., L = 0.11 m, W = 0.054 m, and a = 0.015 m). The geometry and boundary condition of the model are shown in Fig. 12. The different crack angle directions (α) with 30, 45, 60, and 80° are considered as sensitivity analysis of crack propagation under compressive loading conditions using IMTS criterion. The numerical modeling in our research work is based on an assumption that the value of friction coefficient during the crack propagation could be set as zero value as also considered in the reported numerical works in Haeri et al. (2014) and Sharafisafa and Nazem (2014). The condition of stress in polar coordinates near the crack tip and direction of the crack propagation in the specimen with an initial crack angle of 60 are shown in Fig. 13. Analyses of the stress concentration around the tips of the crack have shown that, in polar coordinates, the tangential

145 Page 8 of 10

Arab J Geosci (2016) 9: 145

σ

a

In order to investigate the effect of the crack propagation increment on the crack propagation path, the different crack length increments (i.e., 2, 3, and 4 mm) have been examined in the model with different initial crack angles. The result of the crack propagation path for different crack length increments in a model with initial crack angle of 60° is shown in Fig. 14. It can be seen that the crack initiation angle and general direction of crack propagation in all models with different crack length increments are identical. It means that the initial crack propagates toward the direction of the applied loading. Therefore, the crack length increment of Δa = 4 mm has been selected to simulate the crack propagation in all models with different initial crack angles. The path of the crack propagation in different crack conditions is shown in Fig. 15. As can be seen, the crack is propagated from its tip toward the direction of the applied loading (wing cracks). In Table 2, the obtained results from this study, numerical modeling results from Wu and Wong (2012), and also experimental test (Bobet 2000) on the same specimen and loading condition are summarized. As can be seen that with increasing crack inclination angle, the initiation of crack propagation angle is decreased and there are no ssignificant differences between XFEM results and other laboratory and numerical modeling results (Bobet 2000; Wu and Wong 2012).

α

W

σ Fig. 12 Schematic view of cracked specimen under uniaxial compression

a=0.004

a=0.003

a=0.002

0.09

stresses have two maxima, one with tensile stresses and the other one with compressive stresses (Bobet 1997).

0.085

Pa

0.08

Crack Propagation

Y (m)

0.075

0.07

Initial Crack Tip 0.065

0.06

Initial Crack Tip 0.055 0.02

0.025

0.03

0.035

0.04

X (m)

Fig. 13 Tangential stress concentration near the crack tip in cracked specimen under uniaxial compression with crack angle of 60°

Fig. 14 The effect of crack increment length on crack propagation path in the specimen with crack angle of 60°

Arab J Geosci (2016) 9: 145

Page 9 of 10 145

(a)

(b)

Therefore, the IMTS criterion in XFEM can precisely predict the crack propagation direction in cracked specimen under tension and also uniaxial compression.

Conclusion and discussion

(c)

(d)

Fig. 15 Crack propagation path predicted by proposed criterion in specimen under uniaxial compression with crack angle of a 30, b 45, c 60, and d 80

Table 2 Comparison between results of XFEM, laboratory test, and numerical modeling under uniaxial compression Crack inclination Initiation angle, θ (degree) angle, α (degree) XFEM result laboratory test Numerical modeling (Bobet 2000) (Wu and Wong 2012) 30 45 60 80

81 66 50 20

85 67 54 –

79 65 52 –

Although crack propagation in rock specimen under compressive loading is very important, study about this phenomenon in the literature is very limited and is the main objective of this research work. A numerical code based on XFEM has been developed in which a stress-based criterion (IMTS criterion) has been implemented to determine the direction of crack propagation under different loading condition and different specimen shape and sizes. In this criterion, the tension crack is propagated from the tip of pre-existing crack in direction of the maximum tensile tangential stress, and unlike the MTS criterion, there is no need to calculate the values of stress intensity factors. Since the MTS criterion correctly predicts the direction of crack propagation in tension condition, this approach is validated by MTS criterion in the specimens under tension. The results of crack propagation under compressive loading show that the cracks are propagated from its tip toward direction of compressive loading. The crack initiation angle decreases by increasing the crack inclination angle. The values of crack initiation angle for various angles of pre-existing crack have been compared with the experimental and numerical results reported in the literature. The comparisons indicate that the IMTS criterion could correctly predict the direction of crack propagation under compressive loading. This criterion could be a generalized form of the MTS criterion when it can predict the direction of crack propagation under both compression and tension loading conditions. Since most failure modes in rocks occur under compressive loading in the nature, the IMTS criterion and developed code in this research work may be adopted in different rock engineering problems from laboratory to field scale which are in plan for future research works. The numerical modeling in our research work is subjected to an assumption that the value of friction coefficient during the crack propagation could be set as zero value. Although the variation of friction angle of crack influences the stress condition near the crack surface and finally the magnitude of stresses at which cracking phenomena occur (Park and Bobet 2009), however, according to the literature, the frictional coefficient may (Wong et al. 2001; Wong and Chau 1998) or may not (Park and Bobet 2009) violate the crack propagation angle. The sensitivity analysis about the effect of mechanical properties of crack and material is planned for further research in our ongoing study.

145 Page 10 of 10

References Ahmed A (2009) Extended finite element method (XFEM)—modeling arbitrary discontinuities and failure analysis. Dissertation, Istituto Universitario di Studi Superiori di Pavia Al-Shayea NA (2005) Crack propagation trajectories for rocks under mixed mode I–II fracture. Eng Geol 81:84–97 Amrollahi H, Baghbanan A, Hashemolhosseini H (2011) Measuring fracture toughness of crystalline marbles under modes I and II and mixed mode I–II loading conditions using CCNBD and HCCD specimens. Int J Rock Mech Min Sci 48:1123–1134 Andrianopoulos N, Theocaris P (1985) The Griffith-Orowan fracture theory revisited: the T-criterion. Int J Mech Sci 27:793–801 Belytschko T, Black T (1999) Elastic crack growth in finite elements with minimal remeshing. Int J Numer Methods Eng 45:601–620 Bobet A (1997) Fracture coalescence in rock materials: experimental observations and numerical predictions. Dissertation, Massachusetts Institute of Technology Bobet A (2000) The initiation of secondary cracks in compression. Eng Fract Mech 66:187–219 Bobet A, Einstein H (1998a) Fracture coalescence in rock-type materials under uniaxial and biaxial compression. Int J Rock Mech Min Sci 35:863–888 Bobet A, Einstein HH (1998b) Numerical modeling of fracture coalescence in a model rock material. Int J Fract 92:221–252 Brace WF, Bombolakis E (1963) A note on brittle crack growth in compression. J Geophys Res 68:3709–3713 Chen C-S, Pan E, Amadei B (1998) Fracture mechanics analysis of cracked discs of anisotropic rock using the boundary element method. Int J Rock Mech Min Sci 35:195–218 Dolbow JE (1999) An extended finite element method with discontinuous enrichment for applied mechanics. Dissertation, Northwestern university Eftekhari M, Baghbanan A, Hashemolhosseini H (2014) Determining stress intensity factor for cracked Brazilian disc using extended finite element method. Int J Sci Eng Technol 3:890–893 Eftekhari M, Baghbanan A, Hashemolhosseini H (2015a) Extended finite element simulation of crack propagation in Cracked Brazilian Disc. J Min Environ 6:95–102 Eftekhari M, Baghbanan A, Hashemolhosseini H (2015b) Fracture propagation in a cracked semicircular bend specimen under mixed mode loading using extended finite element method. Arab J Geosci:1–12, DOI 10.1007/s12517-12015-11906-12514 Erdogan F, Sih G (1963) On the crack extension in plates under plane loading and transverse shear. J Basic Eng 85:519 Haeri H, Shahriar K, Marji MF, Moarefvand P (2014) Cracks coalescence mechanism and cracks propagation paths in rock-like specimens containing pre-existing random cracks under compression. J Cent South Univ 21:2404–2414 Hussain M, Pu S, Underwood J (1974) Strain energy release rate for a crack under combined mode I and mode II. Fracture analysis, ASTM STP 560: American Society for Testing and Materials, pp. 2–28 Ingraffea AR, Heuze FE (1980) Finite element models for rock fracture mechanics. Int J Numer Anal Methods Geomech 4:25–43 Jiefan H, Ganglin C, Yonghong Z, Ren W (1990) An experimental study of the strain field development prior to failure of a marble plate under compression. Tectonophysics 175:269–284

Arab J Geosci (2016) 9: 145 Li Y-P, Chen L-Z, Wang Y-H (2005) Experimental research on precracked marble under compression. Int J Solids Struct 42:2505– 2516 Lin X (2000) Estimation of stress intensity factors with considering crack surface contact. In: Thirteen European Conference on Fracture (ECF 13), San Sebastian, Spain Moës N, Dolbow J, Belytschko T (1999) A finite element method for crack growth without remeshing. Int J Numer Methods Eng 46:131– 150 Nuismer R (1975) An energy release rate criterion for mixed mode fracture. Int J Fract 11:245–250 Palaniswamy K (1972) Crack propagation under general in-plane loading. Dissertation, California Institute of Technology Palaniswamy K, Knauss W (1972) Propagation of a crack under general, in-plane tension. Int J Fract 8:114–117 Park C, Bobet A (2009) Crack coalescence in specimens with open and closed flaws: a comparison. Int J Rock Mech Min Sci 46:819–829 Reyes OML (1991) Experimental study and analytical modelling of compressive fracture in brittle materials. Dissertation, Massachusetts Institute of Technology Richardson CL, Hegemann J, Sifakis E, Hellrung J, Teran JM (2011) An XFEM method for modeling geometrically elaborate crack propagation in brittle materials. Int J Numer Methods Eng 88:1042–1065 Sharafisafa M, Nazem M (2014) Application of the distinct element method and the extended finite element method in modelling cracks and coalescence in brittle materials. Comput Mater Sci 91:102–121 Shen B, Stephansson O (1994) Modification of the G-criterion for crack propagation subjected to compression. Eng Fract Mech 47:177–189 Shen B, Stephansson O, Einstein HH, Ghahreman B (1995) Coalescence of fractures under shear stresses in experiments. J Geophys Res Solid Earth 100:5975–5990 Shiryaev A, Kotkis A (1983) Methods for determining fracture toughness of brittle porous materials. Indust Lab 48:917–918 Sih GC (1973) Methods of analysis and solutions of crack problems. Noordhoff International Publishing, The Netherlands Sih GC (1974) Strain-energy-density factor applied to mixed mode crack problems. Int J Fract 10:305–321 Sukumar N, Prévost J-H (2003) Modeling quasi-static crack growth with the extended finite element method part I: computer implementation. Int J Solids Struct 40:7513–7537 Theocaris P, Andrianopoulos N (1982) A modified strain-energy density criterion applied to crack propagation. J Appl Mech 49:81–86 Theocaris P, Sakellariou M (1991) Crack propagation in brittle materials under compressive stresses studied by caustics. J Mater Sci 26: 1640–1646 Tirosh J, Catz E (1981) Mixed-mode fracture angle and fracture locus of materials subjected to compressive loading. Eng Fract Mech 14:27– 38 Wong R, Chau K, Tang C, Lin P (2001) Analysis of crack coalescence in rock-like materials containing three flaws—part I: experimental approach. Int J Rock Mech Min Sci 38:909–924 Wong RH, Chau K (1998) Crack coalescence in a rock-like material containing two cracks. Int J Rock Mech Min Sci 35:147–164 Wu Z, Wong LNY (2012) Frictional crack initiation and propagation analysis using the numerical manifold method. Comput Geotech 39:38–53