- Email: [email protected]
The modelling of rock breakage process by TBM rolling cutters using 3D FEM-SPH coupled method
Tunnelling and Underground Space Technology 61 (2017) 90–103
Contents lists available at ScienceDirect
Tunnelling and Underground Space Technology journal homepage: www.elsevier.com/locate/tust
The modelling of rock breakage process by TBM rolling cutters using 3D FEM-SPH coupled method Nan Xiao a,b, Xiao-Ping Zhou a,b,⇑, Qiu-Ming Gong c a
State Key Laboratory of Coal Mine Disaster Dynamics and Control, Chongqing University, Chongqing 400045, China School of Civil Engineering, Chongqing University, Chongqing 400045, China c Key Laboratory of Urban Security and Disaster Engineering of Ministry of Education, Beijing University of Technology, Beijing 100124, China b
a r t i c l e
i n f o
Article history: Received 9 May 2016 Received in revised form 19 September 2016 Accepted 12 October 2016
Keywords: FEM-SPH Rock breakage process TBM Rock fragment The ejection speed
a b s t r a c t The rock breakage process by a single disc cutter and the ejection speed of rock fragments are simulated using three-dimensional Finite Element Method (FEM) coupled with Smoothed Particle Hydrodynamics (SPH) method in this paper. Firstly, the rock penetration process can be divided into three stages including the formation stage of crushed zone, the formation stage of micro-cracks and the propagation stage of major cracks. Then, during the rock penetration process, three different zones beneath the cutter, namely crushed zone, cracked zone and intact rock zone, can be distinguished based on the degree of cracking. Finally, the ejection speed of fragments during rock penetration process is obtained using FEM-SPH method, which is a more convenient way to study the effect of the harmful ejection of rock fragments on TBM machine. Moreover, the three-dimensional numerical results are in good agreement with experimental observations. Ó 2016 Elsevier Ltd. All rights reserved.
1. Introduction The tunnel boring machine (TBM) has been extensively adopted in tunneling constructions, due to its rapid advance rate, high efficiency, nice tunnel formation and little impact on the surrounding environment and security (Zhao, 2007). The modern era of tunnel boring machines was born in the early 1950s. Since Robbins (1987) summarized the development and application of tunnel boring machines from the 1950s to the 1980s, various applications and models developed for tunnel boring machines have been received extensive attention in the last 30 years, such as application of shield method to urban tunneling (Grandori, 1989); application of multi-micro shield tunneling method to large rectangular cross-section tunnels (Mori and Abe, 2005), application of compact shield tunneling method to urban underground construction (Maeda and Kushiyama, 2005), application of a new hard rock TBM performance prediction model to project planning (Hassanpour et al., 2011), and to blocky rock conditions (Delisio and Zhao, 2014), application of a new method to predict the TBM performance in mixed-face ground for project planning and optimization (Tóth et al., 2013), analysis of TBM performance in
⇑ Corresponding author at: State Key Laboratory of Coal Mine Disaster Dynamics and Control, Chongqing University, Chongqing 400045, China. E-mail address: [email protected] (X.-P. Zhou). http://dx.doi.org/10.1016/j.tust.2016.10.004 0886-7798/Ó 2016 Elsevier Ltd. All rights reserved.
highly jointed rock masses and fault zones (Paltrinieri et al., 2016), and so on. The performance prediction and rock breakage mechanism by TBM cutters are becoming considerably important issues. Many experimental model were designed to predict the performance and to study the rock breakage mechanism of TBM cutters. Chang et al. (2006) used the linear cutting experiment system to derive the optimum cutting condition for Hwangdeung Granite, and suggested that a designer should pay more attention to thrust capacity of a TBM in which the normal force is used as an input data. Zhao et al. (2015) conducted a series of linear cutting experiments to study the crack pattern of Beishan granite under one single normal and inclined disc cutter, and found that the crack patterns for normal and inclined cuttings at the same penetration depth are greatly different due to different inclination angles and cutting forces in the two scenarios. Ma et al. (2016) used fullscale cutting tests with large intact granite specimens to investigate the effect of different confining stress conditions on TBM performance, and found that the confining stress has significant impact on the normal force and rolling force. Gong et al. (2016) developed a new rock breakage experimental platform, and the experimental results showed that this platform can satisfy its design requirements to guide the design and operation of excavation machine. Ejection of rock fragments during rock breakage process by cutters may be a important factor causing wear or damage of TBM and TBM shield jamming disaster. For the ejection of rock
91
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
fragments, Zhang et al. (2004) provided a direct, explicit and quantitative method to determine the average fragment size generated by crack coalescence and the fragment ejection velocity during the dynamic fragmentation process, and concluded that the value of the ejection velocity increases with an increase in the strain rate. Jiang et al. (2015) developed a fully designed observation scheme and an algorithm for velocity estimation, and found that the relationship between the ejection speeds and the elastic modulus of rocks is approximately linear. These above experiments have great
(a) The initial establishment of the finite element model of rock mass
significance for understanding rock fragmentation process and rock fragmentation mechanisms by an indenter. However, there are several shortcomings for experiments: (1) It is very difficult to trace fracture development in experiments in detail. (2) The structure of rock materials used in tests is usually unknown and it is difficult to investigate the effects of structural parameters, such as cracks, joints and faults, on rock fragmentation. (3) The experiments are generally expensive, time consuming and often need a long design phase to be arranged. During the past few decades, many numerical models have been developed based on site survey and data from tunnel project using TBM, which played an increasingly important role in tunnel construction. Moreover, with the rapid development of computer technology, there are many kinds of numerical methods which are able to study rock breakage mechanism, such as Finite Element Method (FEM), and the discrete element method (DEM). Franzius and Potts (2005) investigated the effects of the geometry and the dimension of a 3D FEM model on tunnel-induced surface settlement predictions. Menezes et al. (2014) simulated a rock fragmentation process during the mechanical cutting of rocks using an explicit finite element code, LS-DYNA, and actually predicted the seperation of rock fragments from the base rock slab. Cho et al. (2010) used a finite element code, AUTODYN-3D, to simulate
(b) Generation of SPH particles in the element based on the value of stress or deformation.
(a) The cracking condition for Mode I fracture
(c) The generated SPH particle inherits the distance, quality, density, velocity, poisson’s ratio and constitutive model from the original finite element.
(d) The original finite elementis deleted and no longer involved in the calculation after the physical and mechanical properties of the materials are converted into SPH particles.
(b) The cracking condition for Mode II fracture Fig. 2. The cracking condition.
Table 1 The parameters of Beshan granite (Zhao et al., 2015).
(e) The generatedSPH particles repalce the original finite element to compute the stress and strain. Fig. 1. The coupling and transformation in the coupling FEM-SPH model.
Mass density Young’s modulus Poisson’s ratio Tensile failure stress Fracture toughness Fracture energy
2.6 g/cm3 23.01 GPa 0.188 6.4 MPa 25.7 MPa mm1/2 28.71 N/m
92
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
Fig. 3. The dimension of the numerical cutter (mm).
(a) The numerical model and the constraint condition
(b) The detailed constraint condition of the cutters in the numerical model Fig. 4. The numerical model and the constraint condition of the cutters.
(a) The linear cutting machine
(b) The detailed constraint condition of the cutter in the machine Fig. 5. Linear cutting machine and the constraint condition of the physical cutters.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
three-dimensional dynamic failure observed in linear cutting experiments, and found that the numerical simulations by AUTODYN-3D could be possibly substituted for linear cutting experiments used for performance assessment of TBM. However, FEM remains limiting because it requires that the finite element meshes coincide with cracks. This limitation makes the task of generating the mesh very difficult (Bouchard et al., 2000; Liang et al., 2012). The discrete element method (DEM) was introduced to simulate the rock fragmentation by TBM cutters. Gong et al. (2005, 2006) used Universal Distinct Element Code (UDEC) to study the correlation between joint orientation and rock fragmentation by one single cutter, the numerical results show that the joint orientation and joint spacing can obviously affect the initiation and propagation of crack as well as the fragmentation pattern. Li et al. (2016) used Distinct Element Method (DEM) to simulate the process of TBMs indentation, and concluded that tensile cracks are resulted from the chipping process and shear cracks are induced by the crushing force. However, DEM still has the following limitations (Donze et al., 2009): (1) Fracture is closely related to the size of elements, and that is so called size effect. (2) Cross effect exists because of the difference between the numerical elements and real grains. (3) In order to establish the relationship between the local and macroscopic constitutive laws, data obtained from classical geomechanical tests which may be impractical are used. Smoothed particle hydrodynamics (SPH) is a mesh-free method (Lucy, 1977), which is applied to solve discontinuous problems through continuum mechanics. Libersky and Petschek (1991) and Benz and Asphaug (1995) successfully demonstrated that SPH was applicable for simulation of the fracture process in brittle materials. Ma et al. (2011) applied SPH to simulate the failure of
93
brittle heterogeneous materials by tracing the propagation of the microscopic cracks as well as the macro-mechanical behaviors. Pramanik and Deb (2015) studied dynamic failure mechanism in surrounding brittle rock media under blast-induced stress waves using SPH. Each of these methods has its own merits and shortcomings, and combining them would allow the best performance to be obtained. The Finite Element Method (FEM) based on grid has been widely used to analyze the stability of underground engineering, and it is very efficient and accurate for small deformation, but it has some disadvantages in the large deformation region. Smoothed Particle Hydrodynamics (SPH) without grid is more accurate for large deformation region, but it takes a longer time compared with FEM. FEM-SPH coupled method inherits the advantages of both, and the model has been defined with Lagrangian finite elements at the beginning of the analysis, then converted to SPH particles after the deformation becomes significant, so that FEM-SPH coupled method can be quick and accurate in solving dynamic problem of rock fragmentation by TBM cutters. The coupling FEMSPH model was first proposed by Attaway and Heinstein (1994) and incorporated into PRONTO. Sauer (2000) used the gap function to achieve the coupling FEM-SPH model. De Vuyst et al. (2005) applied the SPH–FEM coupled model to study fluid–structure impact problems. Fourey et al. (2010) simulated the water entry of elastic wedge with high speed and the impact of the water column on elastic plate using the coupling FEM-SPH method. However, to the best of my knowledge, the coupling FEM-SPH model was not applied to investigate the rock breakage process and the ejection of rock fragments by disc cutters. Most studies on rock breakage by TBM cutters mentioned above mainly concentrated on two-dimensional rock penetration
(a) The global meshing strategy of the numerical model
(b) The partial enlarged detail of the numerical model Fig. 6. The global and local meshing strategy of the numerical model.
94
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
process, while three-dimensional rock breakage process still need to be studied. Moreover, the ejection speed of fragments during rock penetration process is not investigated. Therefore, in this paper, a series of numerical simulations are performed using a 3D FEM-SPH coupled method to study the mechanism of the rock breakage and the ejection speed of fragments. Meanwhile, the numerical results are discussed and compared with the physical experiment results. 2. The coupling FEM-SPH model 2.1. The theory of coupling FEM-SPH Lagrangian approach is the basic method to illustrate the physical governing equation. The typical one is FEM which is based on Lagrangian grid (Zienkiewicz and Taylor, 2000; Liu and Quek, 2003). FEM is a very efficient and mature numerical tool to analyze the stability of underground engineering. The accuracy and efficiency of FEM has been proved in small deformation problems. However, due to its formulas derived from grid, FEM encounters a problem when the relative movement becomes too large. This problem is that time history maybe lost, even the calculation will fail if the grid is too small. The main reason is that the time step is controlled by the minimum cell size. One way to solve this problem is to repartition the grid or re-mesh in the problem domain so as to calculate on the new grid without distortion. However, process of re-meshing is very tedious and time-consuming in Lagrangian calculation. Moreover, the diffusion of mass or loss of history will happen during every re-meshing process (Benson, 1992; Charles and Anderson, 1987; Hans, 1999). These defects cause
difficulty to calculate in simulating large deformation although FEM has many advantages. Smoothed particle hydrodynamics method (SPH) is a mesh-free method first proposed in astrophysics, and also based on Lagrangian method (Gingold and Monaghan, 1977; Lucy, 1977), but there is no need to define nodes or elements in this method, instead, we just need to present a given body by a series of points called as particles or pseudo-particles. The basic concept of SPH is to reproduce the physical field variables in an integral representation. Most importantly, the adaptability of SPH has been obtained in the early stage of approximation for field variables, and the arbitrary particle distribution with the adaptability of SPH has no effect on the construction of SPH formulas, so SPH has the ability to handle large deformation problems easily. The other attractive property is the perfect combination of Lagrangian formulas and particle approximation. SPH particles can move under the interaction of external and internal force carrying the properties of material, so it is quite convenient for SPH particles to handle large deformation problems including cracking of brittle rock, rockburst, etc. (Qing et al., 2012). Tracking fragments from primary impacts through a large volume until secondary impact occurs can be very expensive in a FEM analysis but comes at no additional cost in a smoothed particle hydrodynamic analysis. However, the disadvantage of SPH is obvious because the computational efficiency of SPH is lower than that of FEM due to each particle carrying its own mass, velocity and position. In addition, the computational accuracy of SPH has little superiority over that of FEM in small deformation zone. Therefore, the coupling FEM-SPH scheme is selected to reveal the mechanism of rock fragmentation by disc cutters. FEM is used
(a) The normal section of the formation stage of crushed zone
(b) The profile of the formation stage of crushed zone Fig. 7. The formation stage of crushed zone.
95
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
to calculate the stress and strain in the zone undergoing small deformation, while the SPH is employed to simulate large deformation zone of the rock mass. The finite elements within the large deformation zone are removed based on the threshold of conversion to particles, and then be replaced by SPH particles located in the center of the deleted element. The generated particles inherit all the physical and mechanical properties from the deleted element including distance, quality, density, velocity, Poisson’s ratio and constitutive model. The calculation and conversion mode are described as shown in Fig. 1. By adopting the coupling FEM-SPH scheme, these finite elements with excessive deformation are converted into SPH particles, and the excessive deformed finite elements continue to be involved in the calculation after they are converted to SPH particles. Meanwhile, the crack propagation can be traced based on the stress state of SPH particles. The coupling calculation also works for the neighboring finite elements. The coupling FEM-SPH model overcomes the excessive element deformation, which is very common in the conventional FEM analysis. The coupling FEM-SPH method in this paper is numerically stable, so the rounding error caused by numerical instabilities does not need to be considered.
among cracks is ignored. By introducing the consistency condition, Rankine criterion can be expressed in following tensor form.
C ¼ C t;
rI;II ¼ 0
ð1Þ
where rI;II represents a tension softening model (Mode I fracture) in the case of the direct tensile components of stress, and a shear softening/retention model (Mode II fracture) in the case of the shear components of stress, and t represents r1 in the local coordinate system. The cracking condition of brittle materials is more complex than the classical plastic yield condition of the elastic-plastic material. The brittle materials may be in an active opening crack state or a closing/reopening crack state, while the elastic-plastic material will only be in a single plastic yield state. For Mode I fracture, Eq. (1) can be further expressed as
C nn ¼ C nn tnn ; rIt ¼ t nn rIt ðeck nn Þ ¼ 0
where rIt ðeck nn Þ is the tension softening evolution for an active opening crack. From Eq. (1), Rankine criterion for a closing/reopening crack can be obtained as
C nn ¼ C nn tnn ; rIc ¼ t nn rIc eck nn eopen ¼ 0 nn
2.2. Consitituve model of rock breakage by cutters
where
The modelling of crack initiation and propagation is based on the Maximum Stress Criteria (Rankine criterion), which has been commonly used to predict the failure of brittle material. It states that failure occurs when the maximum tensile stress reaches the tensile strength of rocks. In this study, cracking is set to be irrecoverable throughout the whole simulation, and the interaction
ð2Þ
ð3Þ
rIc eck is the crack closing/reopening evolution, which nn eopen nn
depends on the maximum crack opening strain, eopen ¼ maxov erhistory ðeck nn nn Þ. The cracking condition for Mode I fracture is plotted in Fig. 2(a), which represents the tension softening model adopted for the cracking behavior normal to crack surfaces. It should be noted that, although the cracking condition has been written for all
(a) The normal section of the formation stage of microcracks
(b) The profile of the formation stage of microcracks Fig. 8. The formation stage of microcracks.
96
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
(a) The normal section of the propagation stage of major cracks
(b) The profile of the propagation stage of major cracks Fig. 9. The propagation stage of major cracks.
Fig. 10. The stress vector of the failure element.
Fig. 11. Numerical cracking modes of rock under penetration of a single cutter.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
97
Table 2 Sizes and areas of different zones obtained from the proposed method. Partition zone
Area of a single failure Number of the element/mm2 failure element
Failure area/mm2
Crushed zone Intensive cracked zone Sparse cracked zone Total cracked zone
49 49 49
98 875 931 1904
2 17 19
pre-existing cracks in general, only the components of C that refer to the pre-existing cracks are considered in the computations. For Mode II fracture, Eq. (1) can be rewritten as
ck ck C nt ¼ C nt t nt ; rIIs ¼ t nt rIIs g ck nt ; enn ; ett ¼ 0
ð4Þ
ck ck where rIIs g ck is the shear evolution that linearly depends nt ; enn ; ett on the shear strain and the crack opening strain. The cracking condition for Mode II fracture is plotted in Fig. 2(b).
Fig. 13. The displacement of SPH particles 1 along the direction of the X, Y, and Z axis.
(a) The state of the marked particle 1 at the penetration depth of 0.5mm at 0.01s
(b) The state of the marked particle 1 at the penetration depth of 0.5mm at 0.015s
(c) The state of the marked particle 1 at the penetration depth of 0.5mm at 0.02s Fig. 12. The state of the marked SPH particle 1.
98
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
2.3. The numerical model of rock fragmentation The dimension of the numerical sample is 1000 mm 1000 mm 600 mm. The mechanical properties of the numerical model are taken from Beshan Jijicao granite in Yumen, Gansu, China. A series of tests were performed to obtain
the parameters of Beshan granite (Zhao et al., 2015). The testing results are listed in Table 1. In the numerical model, the cutter is set to be the same size as the TBM cutter with 17 in. (432 mm) in diameter and 21 mm in width, as shown in Fig. 3. The material properties of the cutters are obtained from the Robins cutters, which are made of AISI4340
(a) The state of the marked particle 2 at the penetration depth of 1mm at 0.01s
(b) The state of the marked particle 2 at the penetration depth of 1mm at 0.015s
(c) The state of the marked particle 2 at the penetration depth of 0.5mm at 0.02s Fig. 14. The state of the marked SPH particle 2.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
steel. In order to model the cutting process, the constraint conditions of the cutter and rock mass are set to be the same as that of the linear cutting tests, as shown in Figs. 4 and 5. The numerical cutter is fixed with the left and right direction or the direction of the Z-axis because the cutter rolls along the horizontal direction and then intrudes downwards into the rock, as shown in Fig. 4 (a). There are three different penetration depths of one single cutter corresponding to three different stages in this paper (i.e. 0.5 mm, 1.5 mm, 6 mm). The same material properties are employed in all these simulations. The numerical simulation is controlled by displacement. The total cutting time of the numerical experiment is 50 s. The linear cutting machine or rock cutting machine is plotted in Fig. 5(a). A part of the numerical model with the dimension of 250 mm 400 mm 1000 mm below the cutter is meshed with the element size of 7 mm, and the other part of the numerical model is meshed with the approximate element size of 50 mm for the limitation of the computer calculation capacity, as shown in Fig. 6(a). Three FEM elements are just accurately covered by thin sharp side of the cutter, as shown in Fig. 6(b).
3. Numerical simulation of rock penetration by one single disc cutter using the coupling FEM-SPH method Numerical results of the rock penetration by one single disc cutter obtained from FEM-SPH coupled method are plotted in Figs. 7–9. Note that the penetration of the cutter maintains constant during the invasion process, and ranges from 0.5 mm to 1.5 mm during the whole simulation. Cracks propagate as the penetration of cutter increases. The process of rock penetration can be divided into three stages based on the different penentration ranging from 0.5 mm to 1.5 mm: (1) the formation stage of crushed zone, as shown in Fig. 7(a); (2) the formation stage of microcracks as shown in Fig. 8(a); and (3) the propagation stage of major cracks, as shown in Fig. 9(a). The stress state of the rock sample is depicted in Figs. 7(b), 8(b) and 9(b) at the beginning of three stages during rock penetration process. Firstly, as mentioned before, the FEM elements are converted into SPH particles with a critical degree of deformation when the rock is fractured to form cracks, so the area of crack propagation can be obtained from the FEM elements converted into the SPH particles. According to works by Gong et al. (2006), the cracks including the lateral cracks and the median cracks are mainly initiated from the tensile failure element zone, and propagate along the tensile failure element. The crack propagation is caused by the tensile damage at the tip of the cracks. Moreover, it is found from Figs. 9(a) and 10 that the cracks propagate along the orthogonal direction of the tensile stress vector of SPH particles. Therefore, it can be more likely to understand the crack propagation in the rock fragmentation by one single disc cutter by using FEM-SPH coupled model. As the cutter intrudes into the rock, a number of tensile and compressive failure elements are observed on the tip edge of the cracks during the formation stage of the crushed zone. Hertzian crack is initiated at the edge of the cutter, and a symmetric crushed zone forms underneath the cutter, as shown in Fig. 7(a) and (b). The rest of the rock relatively remains intact due to the high hydrostatic pressure applied (Cook et al., 1984; Chiaia, 2001; Liu et al., 2002). In the crushed zone, there are a number of tensile and compressive failure elements. The part of the rock is crushed into powder or extremely small particles (Pang and Goldsmith, 1990), which burst outwards from the bottom of the cutter, as shown in Fig. 7(a).
99
With the increase of penentration, the stress field undergoes a dynamic expansion process, as shown in Figs. 7(b), 8(b) and 9(b), and the microcrack is formed under the crushed zone, as shown in Fig. 8(a) and (b). The rock penetration is a three-dimensional process. The cutter intrudes downward into the rock to form a vertical fan-shaped failure zone under the normal force, and rolls forward to form a horizontal fan-shaped failure zone under the rolling force. As the cutter continue to intrude into the rock, both the crushed zone and lateral cracks remain unchanged, while the median cracks continue to propagate along the loading direction. As a result, the growth length of the median cracks is larger than that of the lateral cracks, as shown in Fig. 9(a) and (b). It is found from Fig. 11 that there are three zones beneath the indenter in the rock indentation test. They are crushed zone, cracked zone and intact rock zone. After the rock is partially cracked into rock powder or extreme small particles, some of them burst outward carrying the kinetic energy, and the others form the crushed zone. The fractured zone which locates between the crushed zone and the intact rock zone can be divided into the fragmentation zone and the crack zone. The crack becomes more and more intensive from outside to inside in the crack zone so that the crack zone can be divided into the intensive crack zone and the sparse crack zone based on the intensive degree of the crack. The intact rock zone remains unchanged after rock has been intruded by the cutter. The sizes and areas of crushed and cracked zone are listed in Table 2. The failure elements within the crushed zone, the microcrack zone and major crack zone are converted into SPH particles based on the strain criteria. Depending on the failure states (compression and tension) of the finite elements, the generated SPH particles are moved outwards and inwards, respectively. During this process, all the physical and mechanical properties of the FEM element are transferred to SPH particles. The displacement and ejection speeds of the element and rock fragment can be easily tracked by this proposed method. For example, during the rock penetration process at the penetration depth of 0.5 mm, 1 mm and 1.5 mm, some SPH particles were marked to monitor its displacement and ejection speed, as shown in Figs. 12, 14 and 16. The value of its displacement and ejection speeds can be acquired from the numerical results, as shown in Figs. 13, 15 and 17. During the rock penetration process at the penetration depth of 0.5 mm, the parabolic property of the monitored SPH particle is obvious over a long distance under the consideration of gravity, as shown in Figs. 12 and 13. Moreover, the ejection speeds of
Fig. 15. The displacement of SPH particles 2 along the direction of the X, Y, and Z axis.
100
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
SPH particles are not identical to each other. A portion of SPH particles demonstrate high speeds and the other portion of SPH particles exhibit low speeds. The initial speed of the ejected SPH particles varies from 24.7 m/s to 51.8 m/s approximately, the ejection speed of the marked SPH particle, namely, the separated rock fragment can reach about 45.2 m/s. The initial throwing angle of the ejected SPH particles is random within a certain degree. The data collected from the numerical simulation shows that the initial throwing angle of the ejected SPH particles varies from 30° to 30°. During the rock penetration process at the penetration depth of 1 mm, a interesting SPH particle was marked, as shown in Figs. 14 and 15. The reason for this interesting SPH particle is that this SPH particle ejects along the approximate horizontal trajectory, then strikes the surface of the rock and bounces out of it. Meanwhile, this SPH particle loses a part of kinetic energy, which leads to local imbalance of stress field. As the penetration increases, the range of the ejection speed and throwing angle of the SPH particles became large. The initial speed of the ejected SPH particles vary from 17.4 m/s to 55.2 m/s approximately. The ejection speed of the
marked SPH particle reaches about 18.3 m/s. The initial throwing angle of the ejected SPH particles is random within a certain degree. The data obtained from the proposed method show that the initial throwing angle of the ejected SPH particles varies from 80° to 80°. During the rock penetration process at the penetration depth of 1.5 mm, two SPH particles are marked to investigate the seperation behavior of the fragment, as shown in Figs. 16 and 17. The data obtained from the proposed method show that although the two particles have the same ejection speed along the direction of the Z axis, the two SPH particles has a different ejection speed along the direction of the X axis and Y axis with a initial integrated state. Moreover, their throwing angles are negative, which is different from the above two cases. The initial speed of the ejected SPH particles varies from 14.2 m/s to 64.7 m/s approximately. The ejection speed of the marked SPH particle reaches about 16.1 m/s. The data obtained from the proposed method show that the initial throwing angle of the ejected SPH particles varies from 140° to 140°. In addition, the present analysis indicates the relationship between the penetration and ejection speed, and the relationship
(a) The state of the marked particle 3, 4 at the penetration depth of 0.5mm at 0.01s
(b) The state of the marked particle 3, 4 at the penetration depth of 0.5mm at 0.015s
(c) The state of the marked particle 3, 4 at the penetration depth of 0.5mm at 0.02s Fig. 16. The state of the marked SPH particle.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
101
of the rolling and normal force is different at the different penetration depth. For example, the peak values of the rolling and normal force at the penetration depth of 0.5 mm are about 310 kN and 24 kN, respectively. While the peak values of the rolling and normal force at the penetration depth of 1 mm are about 370 kN and 27 kN, respectively.
4. The comparison of the numerical results and experimental data
Fig. 17. The displacement of SPH particles 3 and 4 along the direction of the X, Y, and Z axis.
between the penetration and throwing angle. The deeper the penetration, the higher the ejection speed and the larger the throwing angle. This result well agrees with the performance of TBM on rock underground engineering, where rocks that have a higher penetration are determined to have a higher ejection speed and a larger throwing angle. The rolling and normal force at the penetration depth of 1.5 mm are respectively plotted in Fig. 18(a) and (b). Fig. 18 shows that both the rolling and normal force have a elastic-brittle characteristics during the rock penetration process, and the tendency of the rolling and normal force are almost the same, but the peak value
In this paper, the parameters of the rock are mainly derived from the experiments performed by Zhao et al. (2015). The numerical model is the same as experimental one conducted by Zhao et al. (2015), so that it can be naturally and meaningfully compared with the experiment performed by Zhao et al. (2015). The numerical result of the normal force during rock breakage process is in good agreement with the experimental data, as shown in Fig. 19. The numerical result of cracking mode under one single cutter is also in good agreement with the experimental observation, as shown in Fig. 20. The comparison of the crushed and cracked zone are listed in Table 3 (Zhao et al., 2015). Besides, the crack propagation pattern obtained from the proposed model is in good agreement with that obtained from Colorado School of Mines (CSM, 2007). The present numerical result is also similar to the previous numerical result obtained from DEM (Gong et al., 2006). All of the above numerical results show that the coupling FEM-SPH method has a good ability to deal with the rock penetration problem.
(a) Normal force obtained from the proposed method (a) Normal force vs time during rock fragmentation process
(b) Normal force obtained from experiment (Zhao et., 2015) (b) Rolling force vs time during rock fragmentation process Fig. 18. Force vs time in different direction during rock fragmentation process.
Fig. 19. Comparison of the numerical normal force and the experimental normal force.
102
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
(a) Cracking mode under one single cutter obtained from the proposed method
(b) Cracking mode under one single cutter obtained from experiment(Zhao, et al. 2015) Fig. 20. Comparison of cracking mode obtained from the proposed method and experiment.
Table 3 The comparison of the numerical results and experimental data (Zhao et al., 2015). Partition zone
Failure area obtained Failure area obtained from from experiments/mm2 the numerical results/mm2
Crushed zone Intensive cracked zone Sparse cracked zone Total cracked zone
80 891 1131 2102
98 875 931 1904
5. Discussion Although many researchers investigated the rock breakage mechanism by TBM cutters, the ejection of rock fragments during the rock breakage process is not studied. During the rock breakage process, rock fragments may strike the TBM cutters or eject into the machine by the release of internal energy. As a result, wear or damage of TBM and TBM shield jamming disaster occur. In this paper, the 3D coupling FEM-SPH method is first proposed to investigate the mechanism of rock breakage by TBM cutters and ejection phenomenon during rock breakage process. The ejection response of SPH particles or fragments mainly depends on Young’s modulus, Poisson’s ratio, fracture energy and strength of rocks. Besides, the released kinetic energy during rock penetration process is embodied by the harmful ejection of rock fragments at high speed. Therefore, the analysis of the ejection speed of fragments during rock
penetration process obtained by SPH method is a more convenient way to study the effect of the harmful ejection of rock fragments on TBM machine. Then, we carefully compared the numerical results with the previous experimental data so as to verify the numerical rock breakage process. Although the numerical results are in good agreement with the experimental data in many aspects, there still exist some differences. For example, although the tendency of cutter forces obtained from the proposed method is the same as that obtained from experiments, the numerical results seem to be higher than the experimental data. The main reasons are that (1) the cutter force gathered from the groove closest to the edge is eliminated in experiments, while it is taken into account in numerical model; (2) the internal structure of the actual rock is different from that of the numerical specimen, there exist many micro-cracks in the actual rock, while there do not exist micro-cracks in the numerical specimen in this paper. 6. Conclusions The mechanism of rock fragmentation by one single disc cutters is revealed using three-dimensional coupling FEM-SPH method. The main conclusions can be summarized as follows: (1) In numerical simulation, the rock penetration process by one single disc cutter can be divided into three stages, namely, the formation of the crushed zone, the formation of the
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
micro cracks, the propagation of major cracks. During the rock penetration process, three different zones beneath the cutter can be distinguished, namely, crushed zone, cracked zone and intact rock zone based on the degree of cracking. Moreover, the cracked zone can be further subdivided into intensive cracked zone and sparse cracked zone. The numerical results show that the coupling FEM-SPH method has a good ability to deal with the rock penetration problem. (2) The penetration depth affects the ejection speed and throwing angle of the rock fragments. With the increase of the penetration depth, the range of the initial ejected speed and initial throwing angle become wide. The study of rock breakage mechanism and ejection speed of rock fragments during rock breakage process are helpful to understand TBM performance in rock engineering.
Acknowledgements The work is supported by project 973 (No. 2014CB046903), the National Natural Science Foundation of China (No. 51325903, 51679017) and National Natural Science Foundation Project of CQ CSTC (No. CSTC, cstc2013kjrc-ljrccj0001, cstc2015jcyjys30001, cstc2015jcyjys30006 and cstc2016jcyjys0005), and Research fund by the Doctoral Program of Higher Education of China (No. 20130191110037). References Attaway, S., Heinstein, J., 1994. Coupling of smooth particle hydrodynamics with the finite element method. Nucl. Eng. Des. 150, 199–205. Benz, W., Asphaug, E., 1995. Simulations of brittle solids using smooth particle hydrodynamics. Comput. Phys. Commun. 87, 253–265. Benson, D.J., 1992. Computational methods in Lagrangian and Eulerian hydrocodes. Comput. Methods Appl. Mech. Engrg. 99, 235–394. Bouchard, P.O., Baya, F., Chastela, Y., Tovena, I., 2000. Crack propagation modelling using an advanced remeshing technique. Comput. Methods Appl. Mech. Engrg. 189 (3), 723–742. Chang, S.H., Choi, S.W., Bae, G.J., Jeon, S., 2006. Performance prediction of tbm disc cutting on granitic rock by the linear cutting test. Tunn. Undergr. Space Technol. 21 (3), 271. Charles, E., Anderson Jr., 1987. An overview of the theory of hydrocodes. Int. J. Impact Eng 5, 33–59. Chiaia, B., 2001. Fracture mechanisms induced in a brittle material by a hard cutting indenter. Int. J. Solids Struct. 38, 7747–7768. Cho, J.W., Jeon, S., Yu, S.H., Chang, S.H., 2010. Optimum spacing of tbm disc cutters: a numerical simulation using the three-dimensional dynamic fracturing method. Tunn. Undergr. Space Technol. 25 (3), 230–244. Cook, N.G.W., Hood, M., Tsai, F., 1984. Observations of crack growth in hard rock loaded by an indenter. Int. J. Rock Mech. Min. Sci. Geomech. Abstr. 21 (2), 97– 107. CSM, 2007. CSM Computer Model for TBM Performance Prediction
103
Grandori, C., 1989. Development and current experiences with double-shield TBMs in Italy. Tunn. Undergr. Space Technol. 4 (3), 333–336. Gong, Q.M., Zhao, J., Jiao, Y.Y., 2005. Numerical modeling of the effects of joint orientation on rock fragmentation by TBM cutters. Tunn. Undergr. Space Technol. 20 (2), 183–191. Gong, Q.M., Zhao, J., Jiao, Y.Y., 2006. Numerical modeling of the effects of joint spacing on rock fragmentation by TBM cutters. Tunn. Undergr. Space Technol. 21 (1), 46–55. Gong, Q.M., Du, X., Li, Z., Wang, Q.X., 2016. Development of a mechanical rock breakage experimental platform. Tunn. Undergr. Space Technol. 57, 129–136. Hans, U.M., 1999. Review: hydrocodes for structure response to underwater explosions. Shock Vib. 6 (2), 81–96. Hassanpour, J., Rostami, J., Zhao, J., 2011. A new hard rock tbm performance prediction model for project planning. Tunn. Undergr. Space Technol. 26 (5), 595–603. Jiang, Q., Su, G.S., Feng, X.T., Cui, J., Pan, P.Z., Jiang, J.Q., 2015. Observation of rock fragment ejection in post-failure response. Int. J. Rock Mech. Min. Sci. 74, 30–37. Li, X.F., Li, H.B., Liu, Y.Q., Zhou, Q.C., Xia, X., 2016. Numerical simulation of rock fragmentation mechanisms subject to wedge penetration for TBMs. Tunn. Undergr. Space Technol. 53 (5), 96–108. Liang, Z.Z., Xing, H., Wang, S.Y., Williams, D.J., Tang, C.A., 2012. A three-dimensional numerical investigation of the fracture of rock specimens containing a preexisting surface flaw. Comput. Geotech. 45, 19–33. Libersky, L.D., Petschek, A.G., 1991. Smoothed particle hydrodynamics with strength of materials. In: Trease, H., Friits, J., Crowley, W. (Eds.), Proceedings of the Next Free Lagrange Conference, vol. 395. Springer, New York, pp. 248–257. Liu, G.R., Quek, 2003. The Finite Element Method: A Practical Course. Butterworth Heinemann. Liu, H.Y., Kou, S.Q., Lindqvist, P.A., Tang, C.A., 2002. Numerical simulation of the rock fragmentation process induced by indenters. Int. J. Rock Mech. Min. Sci. 39, 491–505. Lucy, L.B., 1977. Numerical approach to testing the fission hypothesis. Astronom. J. 82, 1013–1024. Ma, G.W., Wang, X.J., Ren, F., 2011. Numerical simulation of compressive failure of heterogeneous rock-like materials using SPH method. Int. J. Rock Mech. Min. Sci. 48 (3), 353–363. Maeda, M., Kushiyama, K., 2005. Use of compact shield tunneling method in urban underground construction. Tunn. Undergr. Space Technol. 20 (2), 159–166. Ma, H.S., Gong, Q.M., Wang, J., Yin, L.J., Zhao, X.B., 2016. Study on the influence of confining stress on tbm performance in granite rock by linear cutting test. Tunn. Undergr. Space Technol. 57, 145–150. Menezes, P.L., Lovell, M.R., Avdeev, I.V., Fred Higgs III, C., 2014. Studies on the formation of discontinuous rock fragments during cutting operation. Int. J. Rock Mech. Min. Sci. 71, 131–142. Mori, K., Abe, Y., 2005. Large rectangular cross-section tunneling by the multi-micro shield tunneling (MMST) method. Tunn. Undergr. Space Technol. 20 (2), 129– 141. Paltrinieri, E., Sandrone, F., Zhao, J., 2016. Analysis and estimation of gripper tbm performances in highly fractured and faulted rocks. Tunn. Undergr. Space Technol. 52, 44–61. Pang, S.S., Goldsmith, W., 1990. Investigation of crack formation during loading of brittle rock. Rock Mech. Rock Eng. 23, 53–63. Pramanik, R., Deb, D., 2015. Implementation of smoothed particle hydrodynamics for detonation of explosive with application to rock fragmentation. Rock Mech. Rock Eng. 48 (4), 1683–1698. Qing, Y., Van, J., Leigh, M., 2012. Free-surface flow interactions with deformable structures using an SPH-FEM model. Ocean Eng. 55, 136–147. Robbins, R.J., 1987. Machine tunnelling in the twenty-first century. Tunn. Undergr. Space Technol. 2 (2), 147–154. Sauer, M., 2000. Adaptive Koppling Des Netzfreien SPH-Verfahrens Mit Fniten Elementen Zur Berechnung Von Impaktvorgaengen. University of the Federal Armed Forces, Munich. Tóth, A., Gong, Q.M., Zhao, J., 2013. Case studies of tbm tunneling performance in rock–soil interface mixed ground. Tunn. Undergr. Space Technol. 38 (9), 140– 150. Zhang, Y.Q., Lu, Y., Hao, H., 2004. Analysis of fragment size and ejection velocity at high strain rate. Int. J. Mech. Sci. 46 (1), 27–34. Zhao, J., 2007. Tunnelling in rock-present technology and future challenges. In: Proceedings of World Tunnel Congress, Prague, pp. 22–32. Zhao, X.B., Yao, X.H., Gong, Q.M., Ma, H.S., Li, X.Z., 2015. Comparison study on rock crack pattern under a single normal and inclined disc cutter by linear cutting experiments. Tunn. Undergr. Space Technol. 50, 479–489. Zienkiewicz, O.C., Taylor, R.L., 2000. The Finite Element Method. Butterworth Heinemann.
Contents lists available at ScienceDirect
Tunnelling and Underground Space Technology journal homepage: www.elsevier.com/locate/tust
The modelling of rock breakage process by TBM rolling cutters using 3D FEM-SPH coupled method Nan Xiao a,b, Xiao-Ping Zhou a,b,⇑, Qiu-Ming Gong c a
State Key Laboratory of Coal Mine Disaster Dynamics and Control, Chongqing University, Chongqing 400045, China School of Civil Engineering, Chongqing University, Chongqing 400045, China c Key Laboratory of Urban Security and Disaster Engineering of Ministry of Education, Beijing University of Technology, Beijing 100124, China b
a r t i c l e
i n f o
Article history: Received 9 May 2016 Received in revised form 19 September 2016 Accepted 12 October 2016
Keywords: FEM-SPH Rock breakage process TBM Rock fragment The ejection speed
a b s t r a c t The rock breakage process by a single disc cutter and the ejection speed of rock fragments are simulated using three-dimensional Finite Element Method (FEM) coupled with Smoothed Particle Hydrodynamics (SPH) method in this paper. Firstly, the rock penetration process can be divided into three stages including the formation stage of crushed zone, the formation stage of micro-cracks and the propagation stage of major cracks. Then, during the rock penetration process, three different zones beneath the cutter, namely crushed zone, cracked zone and intact rock zone, can be distinguished based on the degree of cracking. Finally, the ejection speed of fragments during rock penetration process is obtained using FEM-SPH method, which is a more convenient way to study the effect of the harmful ejection of rock fragments on TBM machine. Moreover, the three-dimensional numerical results are in good agreement with experimental observations. Ó 2016 Elsevier Ltd. All rights reserved.
1. Introduction The tunnel boring machine (TBM) has been extensively adopted in tunneling constructions, due to its rapid advance rate, high efficiency, nice tunnel formation and little impact on the surrounding environment and security (Zhao, 2007). The modern era of tunnel boring machines was born in the early 1950s. Since Robbins (1987) summarized the development and application of tunnel boring machines from the 1950s to the 1980s, various applications and models developed for tunnel boring machines have been received extensive attention in the last 30 years, such as application of shield method to urban tunneling (Grandori, 1989); application of multi-micro shield tunneling method to large rectangular cross-section tunnels (Mori and Abe, 2005), application of compact shield tunneling method to urban underground construction (Maeda and Kushiyama, 2005), application of a new hard rock TBM performance prediction model to project planning (Hassanpour et al., 2011), and to blocky rock conditions (Delisio and Zhao, 2014), application of a new method to predict the TBM performance in mixed-face ground for project planning and optimization (Tóth et al., 2013), analysis of TBM performance in
⇑ Corresponding author at: State Key Laboratory of Coal Mine Disaster Dynamics and Control, Chongqing University, Chongqing 400045, China. E-mail address: [email protected] (X.-P. Zhou). http://dx.doi.org/10.1016/j.tust.2016.10.004 0886-7798/Ó 2016 Elsevier Ltd. All rights reserved.
highly jointed rock masses and fault zones (Paltrinieri et al., 2016), and so on. The performance prediction and rock breakage mechanism by TBM cutters are becoming considerably important issues. Many experimental model were designed to predict the performance and to study the rock breakage mechanism of TBM cutters. Chang et al. (2006) used the linear cutting experiment system to derive the optimum cutting condition for Hwangdeung Granite, and suggested that a designer should pay more attention to thrust capacity of a TBM in which the normal force is used as an input data. Zhao et al. (2015) conducted a series of linear cutting experiments to study the crack pattern of Beishan granite under one single normal and inclined disc cutter, and found that the crack patterns for normal and inclined cuttings at the same penetration depth are greatly different due to different inclination angles and cutting forces in the two scenarios. Ma et al. (2016) used fullscale cutting tests with large intact granite specimens to investigate the effect of different confining stress conditions on TBM performance, and found that the confining stress has significant impact on the normal force and rolling force. Gong et al. (2016) developed a new rock breakage experimental platform, and the experimental results showed that this platform can satisfy its design requirements to guide the design and operation of excavation machine. Ejection of rock fragments during rock breakage process by cutters may be a important factor causing wear or damage of TBM and TBM shield jamming disaster. For the ejection of rock
91
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
fragments, Zhang et al. (2004) provided a direct, explicit and quantitative method to determine the average fragment size generated by crack coalescence and the fragment ejection velocity during the dynamic fragmentation process, and concluded that the value of the ejection velocity increases with an increase in the strain rate. Jiang et al. (2015) developed a fully designed observation scheme and an algorithm for velocity estimation, and found that the relationship between the ejection speeds and the elastic modulus of rocks is approximately linear. These above experiments have great
(a) The initial establishment of the finite element model of rock mass
significance for understanding rock fragmentation process and rock fragmentation mechanisms by an indenter. However, there are several shortcomings for experiments: (1) It is very difficult to trace fracture development in experiments in detail. (2) The structure of rock materials used in tests is usually unknown and it is difficult to investigate the effects of structural parameters, such as cracks, joints and faults, on rock fragmentation. (3) The experiments are generally expensive, time consuming and often need a long design phase to be arranged. During the past few decades, many numerical models have been developed based on site survey and data from tunnel project using TBM, which played an increasingly important role in tunnel construction. Moreover, with the rapid development of computer technology, there are many kinds of numerical methods which are able to study rock breakage mechanism, such as Finite Element Method (FEM), and the discrete element method (DEM). Franzius and Potts (2005) investigated the effects of the geometry and the dimension of a 3D FEM model on tunnel-induced surface settlement predictions. Menezes et al. (2014) simulated a rock fragmentation process during the mechanical cutting of rocks using an explicit finite element code, LS-DYNA, and actually predicted the seperation of rock fragments from the base rock slab. Cho et al. (2010) used a finite element code, AUTODYN-3D, to simulate
(b) Generation of SPH particles in the element based on the value of stress or deformation.
(a) The cracking condition for Mode I fracture
(c) The generated SPH particle inherits the distance, quality, density, velocity, poisson’s ratio and constitutive model from the original finite element.
(d) The original finite elementis deleted and no longer involved in the calculation after the physical and mechanical properties of the materials are converted into SPH particles.
(b) The cracking condition for Mode II fracture Fig. 2. The cracking condition.
Table 1 The parameters of Beshan granite (Zhao et al., 2015).
(e) The generatedSPH particles repalce the original finite element to compute the stress and strain. Fig. 1. The coupling and transformation in the coupling FEM-SPH model.
Mass density Young’s modulus Poisson’s ratio Tensile failure stress Fracture toughness Fracture energy
2.6 g/cm3 23.01 GPa 0.188 6.4 MPa 25.7 MPa mm1/2 28.71 N/m
92
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
Fig. 3. The dimension of the numerical cutter (mm).
(a) The numerical model and the constraint condition
(b) The detailed constraint condition of the cutters in the numerical model Fig. 4. The numerical model and the constraint condition of the cutters.
(a) The linear cutting machine
(b) The detailed constraint condition of the cutter in the machine Fig. 5. Linear cutting machine and the constraint condition of the physical cutters.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
three-dimensional dynamic failure observed in linear cutting experiments, and found that the numerical simulations by AUTODYN-3D could be possibly substituted for linear cutting experiments used for performance assessment of TBM. However, FEM remains limiting because it requires that the finite element meshes coincide with cracks. This limitation makes the task of generating the mesh very difficult (Bouchard et al., 2000; Liang et al., 2012). The discrete element method (DEM) was introduced to simulate the rock fragmentation by TBM cutters. Gong et al. (2005, 2006) used Universal Distinct Element Code (UDEC) to study the correlation between joint orientation and rock fragmentation by one single cutter, the numerical results show that the joint orientation and joint spacing can obviously affect the initiation and propagation of crack as well as the fragmentation pattern. Li et al. (2016) used Distinct Element Method (DEM) to simulate the process of TBMs indentation, and concluded that tensile cracks are resulted from the chipping process and shear cracks are induced by the crushing force. However, DEM still has the following limitations (Donze et al., 2009): (1) Fracture is closely related to the size of elements, and that is so called size effect. (2) Cross effect exists because of the difference between the numerical elements and real grains. (3) In order to establish the relationship between the local and macroscopic constitutive laws, data obtained from classical geomechanical tests which may be impractical are used. Smoothed particle hydrodynamics (SPH) is a mesh-free method (Lucy, 1977), which is applied to solve discontinuous problems through continuum mechanics. Libersky and Petschek (1991) and Benz and Asphaug (1995) successfully demonstrated that SPH was applicable for simulation of the fracture process in brittle materials. Ma et al. (2011) applied SPH to simulate the failure of
93
brittle heterogeneous materials by tracing the propagation of the microscopic cracks as well as the macro-mechanical behaviors. Pramanik and Deb (2015) studied dynamic failure mechanism in surrounding brittle rock media under blast-induced stress waves using SPH. Each of these methods has its own merits and shortcomings, and combining them would allow the best performance to be obtained. The Finite Element Method (FEM) based on grid has been widely used to analyze the stability of underground engineering, and it is very efficient and accurate for small deformation, but it has some disadvantages in the large deformation region. Smoothed Particle Hydrodynamics (SPH) without grid is more accurate for large deformation region, but it takes a longer time compared with FEM. FEM-SPH coupled method inherits the advantages of both, and the model has been defined with Lagrangian finite elements at the beginning of the analysis, then converted to SPH particles after the deformation becomes significant, so that FEM-SPH coupled method can be quick and accurate in solving dynamic problem of rock fragmentation by TBM cutters. The coupling FEMSPH model was first proposed by Attaway and Heinstein (1994) and incorporated into PRONTO. Sauer (2000) used the gap function to achieve the coupling FEM-SPH model. De Vuyst et al. (2005) applied the SPH–FEM coupled model to study fluid–structure impact problems. Fourey et al. (2010) simulated the water entry of elastic wedge with high speed and the impact of the water column on elastic plate using the coupling FEM-SPH method. However, to the best of my knowledge, the coupling FEM-SPH model was not applied to investigate the rock breakage process and the ejection of rock fragments by disc cutters. Most studies on rock breakage by TBM cutters mentioned above mainly concentrated on two-dimensional rock penetration
(a) The global meshing strategy of the numerical model
(b) The partial enlarged detail of the numerical model Fig. 6. The global and local meshing strategy of the numerical model.
94
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
process, while three-dimensional rock breakage process still need to be studied. Moreover, the ejection speed of fragments during rock penetration process is not investigated. Therefore, in this paper, a series of numerical simulations are performed using a 3D FEM-SPH coupled method to study the mechanism of the rock breakage and the ejection speed of fragments. Meanwhile, the numerical results are discussed and compared with the physical experiment results. 2. The coupling FEM-SPH model 2.1. The theory of coupling FEM-SPH Lagrangian approach is the basic method to illustrate the physical governing equation. The typical one is FEM which is based on Lagrangian grid (Zienkiewicz and Taylor, 2000; Liu and Quek, 2003). FEM is a very efficient and mature numerical tool to analyze the stability of underground engineering. The accuracy and efficiency of FEM has been proved in small deformation problems. However, due to its formulas derived from grid, FEM encounters a problem when the relative movement becomes too large. This problem is that time history maybe lost, even the calculation will fail if the grid is too small. The main reason is that the time step is controlled by the minimum cell size. One way to solve this problem is to repartition the grid or re-mesh in the problem domain so as to calculate on the new grid without distortion. However, process of re-meshing is very tedious and time-consuming in Lagrangian calculation. Moreover, the diffusion of mass or loss of history will happen during every re-meshing process (Benson, 1992; Charles and Anderson, 1987; Hans, 1999). These defects cause
difficulty to calculate in simulating large deformation although FEM has many advantages. Smoothed particle hydrodynamics method (SPH) is a mesh-free method first proposed in astrophysics, and also based on Lagrangian method (Gingold and Monaghan, 1977; Lucy, 1977), but there is no need to define nodes or elements in this method, instead, we just need to present a given body by a series of points called as particles or pseudo-particles. The basic concept of SPH is to reproduce the physical field variables in an integral representation. Most importantly, the adaptability of SPH has been obtained in the early stage of approximation for field variables, and the arbitrary particle distribution with the adaptability of SPH has no effect on the construction of SPH formulas, so SPH has the ability to handle large deformation problems easily. The other attractive property is the perfect combination of Lagrangian formulas and particle approximation. SPH particles can move under the interaction of external and internal force carrying the properties of material, so it is quite convenient for SPH particles to handle large deformation problems including cracking of brittle rock, rockburst, etc. (Qing et al., 2012). Tracking fragments from primary impacts through a large volume until secondary impact occurs can be very expensive in a FEM analysis but comes at no additional cost in a smoothed particle hydrodynamic analysis. However, the disadvantage of SPH is obvious because the computational efficiency of SPH is lower than that of FEM due to each particle carrying its own mass, velocity and position. In addition, the computational accuracy of SPH has little superiority over that of FEM in small deformation zone. Therefore, the coupling FEM-SPH scheme is selected to reveal the mechanism of rock fragmentation by disc cutters. FEM is used
(a) The normal section of the formation stage of crushed zone
(b) The profile of the formation stage of crushed zone Fig. 7. The formation stage of crushed zone.
95
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
to calculate the stress and strain in the zone undergoing small deformation, while the SPH is employed to simulate large deformation zone of the rock mass. The finite elements within the large deformation zone are removed based on the threshold of conversion to particles, and then be replaced by SPH particles located in the center of the deleted element. The generated particles inherit all the physical and mechanical properties from the deleted element including distance, quality, density, velocity, Poisson’s ratio and constitutive model. The calculation and conversion mode are described as shown in Fig. 1. By adopting the coupling FEM-SPH scheme, these finite elements with excessive deformation are converted into SPH particles, and the excessive deformed finite elements continue to be involved in the calculation after they are converted to SPH particles. Meanwhile, the crack propagation can be traced based on the stress state of SPH particles. The coupling calculation also works for the neighboring finite elements. The coupling FEM-SPH model overcomes the excessive element deformation, which is very common in the conventional FEM analysis. The coupling FEM-SPH method in this paper is numerically stable, so the rounding error caused by numerical instabilities does not need to be considered.
among cracks is ignored. By introducing the consistency condition, Rankine criterion can be expressed in following tensor form.
C ¼ C t;
rI;II ¼ 0
ð1Þ
where rI;II represents a tension softening model (Mode I fracture) in the case of the direct tensile components of stress, and a shear softening/retention model (Mode II fracture) in the case of the shear components of stress, and t represents r1 in the local coordinate system. The cracking condition of brittle materials is more complex than the classical plastic yield condition of the elastic-plastic material. The brittle materials may be in an active opening crack state or a closing/reopening crack state, while the elastic-plastic material will only be in a single plastic yield state. For Mode I fracture, Eq. (1) can be further expressed as
C nn ¼ C nn tnn ; rIt ¼ t nn rIt ðeck nn Þ ¼ 0
where rIt ðeck nn Þ is the tension softening evolution for an active opening crack. From Eq. (1), Rankine criterion for a closing/reopening crack can be obtained as
C nn ¼ C nn tnn ; rIc ¼ t nn rIc eck nn eopen ¼ 0 nn
2.2. Consitituve model of rock breakage by cutters
where
The modelling of crack initiation and propagation is based on the Maximum Stress Criteria (Rankine criterion), which has been commonly used to predict the failure of brittle material. It states that failure occurs when the maximum tensile stress reaches the tensile strength of rocks. In this study, cracking is set to be irrecoverable throughout the whole simulation, and the interaction
ð2Þ
ð3Þ
rIc eck is the crack closing/reopening evolution, which nn eopen nn
depends on the maximum crack opening strain, eopen ¼ maxov erhistory ðeck nn nn Þ. The cracking condition for Mode I fracture is plotted in Fig. 2(a), which represents the tension softening model adopted for the cracking behavior normal to crack surfaces. It should be noted that, although the cracking condition has been written for all
(a) The normal section of the formation stage of microcracks
(b) The profile of the formation stage of microcracks Fig. 8. The formation stage of microcracks.
96
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
(a) The normal section of the propagation stage of major cracks
(b) The profile of the propagation stage of major cracks Fig. 9. The propagation stage of major cracks.
Fig. 10. The stress vector of the failure element.
Fig. 11. Numerical cracking modes of rock under penetration of a single cutter.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
97
Table 2 Sizes and areas of different zones obtained from the proposed method. Partition zone
Area of a single failure Number of the element/mm2 failure element
Failure area/mm2
Crushed zone Intensive cracked zone Sparse cracked zone Total cracked zone
49 49 49
98 875 931 1904
2 17 19
pre-existing cracks in general, only the components of C that refer to the pre-existing cracks are considered in the computations. For Mode II fracture, Eq. (1) can be rewritten as
ck ck C nt ¼ C nt t nt ; rIIs ¼ t nt rIIs g ck nt ; enn ; ett ¼ 0
ð4Þ
ck ck where rIIs g ck is the shear evolution that linearly depends nt ; enn ; ett on the shear strain and the crack opening strain. The cracking condition for Mode II fracture is plotted in Fig. 2(b).
Fig. 13. The displacement of SPH particles 1 along the direction of the X, Y, and Z axis.
(a) The state of the marked particle 1 at the penetration depth of 0.5mm at 0.01s
(b) The state of the marked particle 1 at the penetration depth of 0.5mm at 0.015s
(c) The state of the marked particle 1 at the penetration depth of 0.5mm at 0.02s Fig. 12. The state of the marked SPH particle 1.
98
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
2.3. The numerical model of rock fragmentation The dimension of the numerical sample is 1000 mm 1000 mm 600 mm. The mechanical properties of the numerical model are taken from Beshan Jijicao granite in Yumen, Gansu, China. A series of tests were performed to obtain
the parameters of Beshan granite (Zhao et al., 2015). The testing results are listed in Table 1. In the numerical model, the cutter is set to be the same size as the TBM cutter with 17 in. (432 mm) in diameter and 21 mm in width, as shown in Fig. 3. The material properties of the cutters are obtained from the Robins cutters, which are made of AISI4340
(a) The state of the marked particle 2 at the penetration depth of 1mm at 0.01s
(b) The state of the marked particle 2 at the penetration depth of 1mm at 0.015s
(c) The state of the marked particle 2 at the penetration depth of 0.5mm at 0.02s Fig. 14. The state of the marked SPH particle 2.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
steel. In order to model the cutting process, the constraint conditions of the cutter and rock mass are set to be the same as that of the linear cutting tests, as shown in Figs. 4 and 5. The numerical cutter is fixed with the left and right direction or the direction of the Z-axis because the cutter rolls along the horizontal direction and then intrudes downwards into the rock, as shown in Fig. 4 (a). There are three different penetration depths of one single cutter corresponding to three different stages in this paper (i.e. 0.5 mm, 1.5 mm, 6 mm). The same material properties are employed in all these simulations. The numerical simulation is controlled by displacement. The total cutting time of the numerical experiment is 50 s. The linear cutting machine or rock cutting machine is plotted in Fig. 5(a). A part of the numerical model with the dimension of 250 mm 400 mm 1000 mm below the cutter is meshed with the element size of 7 mm, and the other part of the numerical model is meshed with the approximate element size of 50 mm for the limitation of the computer calculation capacity, as shown in Fig. 6(a). Three FEM elements are just accurately covered by thin sharp side of the cutter, as shown in Fig. 6(b).
3. Numerical simulation of rock penetration by one single disc cutter using the coupling FEM-SPH method Numerical results of the rock penetration by one single disc cutter obtained from FEM-SPH coupled method are plotted in Figs. 7–9. Note that the penetration of the cutter maintains constant during the invasion process, and ranges from 0.5 mm to 1.5 mm during the whole simulation. Cracks propagate as the penetration of cutter increases. The process of rock penetration can be divided into three stages based on the different penentration ranging from 0.5 mm to 1.5 mm: (1) the formation stage of crushed zone, as shown in Fig. 7(a); (2) the formation stage of microcracks as shown in Fig. 8(a); and (3) the propagation stage of major cracks, as shown in Fig. 9(a). The stress state of the rock sample is depicted in Figs. 7(b), 8(b) and 9(b) at the beginning of three stages during rock penetration process. Firstly, as mentioned before, the FEM elements are converted into SPH particles with a critical degree of deformation when the rock is fractured to form cracks, so the area of crack propagation can be obtained from the FEM elements converted into the SPH particles. According to works by Gong et al. (2006), the cracks including the lateral cracks and the median cracks are mainly initiated from the tensile failure element zone, and propagate along the tensile failure element. The crack propagation is caused by the tensile damage at the tip of the cracks. Moreover, it is found from Figs. 9(a) and 10 that the cracks propagate along the orthogonal direction of the tensile stress vector of SPH particles. Therefore, it can be more likely to understand the crack propagation in the rock fragmentation by one single disc cutter by using FEM-SPH coupled model. As the cutter intrudes into the rock, a number of tensile and compressive failure elements are observed on the tip edge of the cracks during the formation stage of the crushed zone. Hertzian crack is initiated at the edge of the cutter, and a symmetric crushed zone forms underneath the cutter, as shown in Fig. 7(a) and (b). The rest of the rock relatively remains intact due to the high hydrostatic pressure applied (Cook et al., 1984; Chiaia, 2001; Liu et al., 2002). In the crushed zone, there are a number of tensile and compressive failure elements. The part of the rock is crushed into powder or extremely small particles (Pang and Goldsmith, 1990), which burst outwards from the bottom of the cutter, as shown in Fig. 7(a).
99
With the increase of penentration, the stress field undergoes a dynamic expansion process, as shown in Figs. 7(b), 8(b) and 9(b), and the microcrack is formed under the crushed zone, as shown in Fig. 8(a) and (b). The rock penetration is a three-dimensional process. The cutter intrudes downward into the rock to form a vertical fan-shaped failure zone under the normal force, and rolls forward to form a horizontal fan-shaped failure zone under the rolling force. As the cutter continue to intrude into the rock, both the crushed zone and lateral cracks remain unchanged, while the median cracks continue to propagate along the loading direction. As a result, the growth length of the median cracks is larger than that of the lateral cracks, as shown in Fig. 9(a) and (b). It is found from Fig. 11 that there are three zones beneath the indenter in the rock indentation test. They are crushed zone, cracked zone and intact rock zone. After the rock is partially cracked into rock powder or extreme small particles, some of them burst outward carrying the kinetic energy, and the others form the crushed zone. The fractured zone which locates between the crushed zone and the intact rock zone can be divided into the fragmentation zone and the crack zone. The crack becomes more and more intensive from outside to inside in the crack zone so that the crack zone can be divided into the intensive crack zone and the sparse crack zone based on the intensive degree of the crack. The intact rock zone remains unchanged after rock has been intruded by the cutter. The sizes and areas of crushed and cracked zone are listed in Table 2. The failure elements within the crushed zone, the microcrack zone and major crack zone are converted into SPH particles based on the strain criteria. Depending on the failure states (compression and tension) of the finite elements, the generated SPH particles are moved outwards and inwards, respectively. During this process, all the physical and mechanical properties of the FEM element are transferred to SPH particles. The displacement and ejection speeds of the element and rock fragment can be easily tracked by this proposed method. For example, during the rock penetration process at the penetration depth of 0.5 mm, 1 mm and 1.5 mm, some SPH particles were marked to monitor its displacement and ejection speed, as shown in Figs. 12, 14 and 16. The value of its displacement and ejection speeds can be acquired from the numerical results, as shown in Figs. 13, 15 and 17. During the rock penetration process at the penetration depth of 0.5 mm, the parabolic property of the monitored SPH particle is obvious over a long distance under the consideration of gravity, as shown in Figs. 12 and 13. Moreover, the ejection speeds of
Fig. 15. The displacement of SPH particles 2 along the direction of the X, Y, and Z axis.
100
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
SPH particles are not identical to each other. A portion of SPH particles demonstrate high speeds and the other portion of SPH particles exhibit low speeds. The initial speed of the ejected SPH particles varies from 24.7 m/s to 51.8 m/s approximately, the ejection speed of the marked SPH particle, namely, the separated rock fragment can reach about 45.2 m/s. The initial throwing angle of the ejected SPH particles is random within a certain degree. The data collected from the numerical simulation shows that the initial throwing angle of the ejected SPH particles varies from 30° to 30°. During the rock penetration process at the penetration depth of 1 mm, a interesting SPH particle was marked, as shown in Figs. 14 and 15. The reason for this interesting SPH particle is that this SPH particle ejects along the approximate horizontal trajectory, then strikes the surface of the rock and bounces out of it. Meanwhile, this SPH particle loses a part of kinetic energy, which leads to local imbalance of stress field. As the penetration increases, the range of the ejection speed and throwing angle of the SPH particles became large. The initial speed of the ejected SPH particles vary from 17.4 m/s to 55.2 m/s approximately. The ejection speed of the
marked SPH particle reaches about 18.3 m/s. The initial throwing angle of the ejected SPH particles is random within a certain degree. The data obtained from the proposed method show that the initial throwing angle of the ejected SPH particles varies from 80° to 80°. During the rock penetration process at the penetration depth of 1.5 mm, two SPH particles are marked to investigate the seperation behavior of the fragment, as shown in Figs. 16 and 17. The data obtained from the proposed method show that although the two particles have the same ejection speed along the direction of the Z axis, the two SPH particles has a different ejection speed along the direction of the X axis and Y axis with a initial integrated state. Moreover, their throwing angles are negative, which is different from the above two cases. The initial speed of the ejected SPH particles varies from 14.2 m/s to 64.7 m/s approximately. The ejection speed of the marked SPH particle reaches about 16.1 m/s. The data obtained from the proposed method show that the initial throwing angle of the ejected SPH particles varies from 140° to 140°. In addition, the present analysis indicates the relationship between the penetration and ejection speed, and the relationship
(a) The state of the marked particle 3, 4 at the penetration depth of 0.5mm at 0.01s
(b) The state of the marked particle 3, 4 at the penetration depth of 0.5mm at 0.015s
(c) The state of the marked particle 3, 4 at the penetration depth of 0.5mm at 0.02s Fig. 16. The state of the marked SPH particle.
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
101
of the rolling and normal force is different at the different penetration depth. For example, the peak values of the rolling and normal force at the penetration depth of 0.5 mm are about 310 kN and 24 kN, respectively. While the peak values of the rolling and normal force at the penetration depth of 1 mm are about 370 kN and 27 kN, respectively.
4. The comparison of the numerical results and experimental data
Fig. 17. The displacement of SPH particles 3 and 4 along the direction of the X, Y, and Z axis.
between the penetration and throwing angle. The deeper the penetration, the higher the ejection speed and the larger the throwing angle. This result well agrees with the performance of TBM on rock underground engineering, where rocks that have a higher penetration are determined to have a higher ejection speed and a larger throwing angle. The rolling and normal force at the penetration depth of 1.5 mm are respectively plotted in Fig. 18(a) and (b). Fig. 18 shows that both the rolling and normal force have a elastic-brittle characteristics during the rock penetration process, and the tendency of the rolling and normal force are almost the same, but the peak value
In this paper, the parameters of the rock are mainly derived from the experiments performed by Zhao et al. (2015). The numerical model is the same as experimental one conducted by Zhao et al. (2015), so that it can be naturally and meaningfully compared with the experiment performed by Zhao et al. (2015). The numerical result of the normal force during rock breakage process is in good agreement with the experimental data, as shown in Fig. 19. The numerical result of cracking mode under one single cutter is also in good agreement with the experimental observation, as shown in Fig. 20. The comparison of the crushed and cracked zone are listed in Table 3 (Zhao et al., 2015). Besides, the crack propagation pattern obtained from the proposed model is in good agreement with that obtained from Colorado School of Mines (CSM, 2007). The present numerical result is also similar to the previous numerical result obtained from DEM (Gong et al., 2006). All of the above numerical results show that the coupling FEM-SPH method has a good ability to deal with the rock penetration problem.
(a) Normal force obtained from the proposed method (a) Normal force vs time during rock fragmentation process
(b) Normal force obtained from experiment (Zhao et., 2015) (b) Rolling force vs time during rock fragmentation process Fig. 18. Force vs time in different direction during rock fragmentation process.
Fig. 19. Comparison of the numerical normal force and the experimental normal force.
102
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
(a) Cracking mode under one single cutter obtained from the proposed method
(b) Cracking mode under one single cutter obtained from experiment(Zhao, et al. 2015) Fig. 20. Comparison of cracking mode obtained from the proposed method and experiment.
Table 3 The comparison of the numerical results and experimental data (Zhao et al., 2015). Partition zone
Failure area obtained Failure area obtained from from experiments/mm2 the numerical results/mm2
Crushed zone Intensive cracked zone Sparse cracked zone Total cracked zone
80 891 1131 2102
98 875 931 1904
5. Discussion Although many researchers investigated the rock breakage mechanism by TBM cutters, the ejection of rock fragments during the rock breakage process is not studied. During the rock breakage process, rock fragments may strike the TBM cutters or eject into the machine by the release of internal energy. As a result, wear or damage of TBM and TBM shield jamming disaster occur. In this paper, the 3D coupling FEM-SPH method is first proposed to investigate the mechanism of rock breakage by TBM cutters and ejection phenomenon during rock breakage process. The ejection response of SPH particles or fragments mainly depends on Young’s modulus, Poisson’s ratio, fracture energy and strength of rocks. Besides, the released kinetic energy during rock penetration process is embodied by the harmful ejection of rock fragments at high speed. Therefore, the analysis of the ejection speed of fragments during rock
penetration process obtained by SPH method is a more convenient way to study the effect of the harmful ejection of rock fragments on TBM machine. Then, we carefully compared the numerical results with the previous experimental data so as to verify the numerical rock breakage process. Although the numerical results are in good agreement with the experimental data in many aspects, there still exist some differences. For example, although the tendency of cutter forces obtained from the proposed method is the same as that obtained from experiments, the numerical results seem to be higher than the experimental data. The main reasons are that (1) the cutter force gathered from the groove closest to the edge is eliminated in experiments, while it is taken into account in numerical model; (2) the internal structure of the actual rock is different from that of the numerical specimen, there exist many micro-cracks in the actual rock, while there do not exist micro-cracks in the numerical specimen in this paper. 6. Conclusions The mechanism of rock fragmentation by one single disc cutters is revealed using three-dimensional coupling FEM-SPH method. The main conclusions can be summarized as follows: (1) In numerical simulation, the rock penetration process by one single disc cutter can be divided into three stages, namely, the formation of the crushed zone, the formation of the
N. Xiao et al. / Tunnelling and Underground Space Technology 61 (2017) 90–103
micro cracks, the propagation of major cracks. During the rock penetration process, three different zones beneath the cutter can be distinguished, namely, crushed zone, cracked zone and intact rock zone based on the degree of cracking. Moreover, the cracked zone can be further subdivided into intensive cracked zone and sparse cracked zone. The numerical results show that the coupling FEM-SPH method has a good ability to deal with the rock penetration problem. (2) The penetration depth affects the ejection speed and throwing angle of the rock fragments. With the increase of the penetration depth, the range of the initial ejected speed and initial throwing angle become wide. The study of rock breakage mechanism and ejection speed of rock fragments during rock breakage process are helpful to understand TBM performance in rock engineering.
Acknowledgements The work is supported by project 973 (No. 2014CB046903), the National Natural Science Foundation of China (No. 51325903, 51679017) and National Natural Science Foundation Project of CQ CSTC (No. CSTC, cstc2013kjrc-ljrccj0001, cstc2015jcyjys30001, cstc2015jcyjys30006 and cstc2016jcyjys0005), and Research fund by the Doctoral Program of Higher Education of China (No. 20130191110037). References Attaway, S., Heinstein, J., 1994. Coupling of smooth particle hydrodynamics with the finite element method. Nucl. Eng. Des. 150, 199–205. Benz, W., Asphaug, E., 1995. Simulations of brittle solids using smooth particle hydrodynamics. Comput. Phys. Commun. 87, 253–265. Benson, D.J., 1992. Computational methods in Lagrangian and Eulerian hydrocodes. Comput. Methods Appl. Mech. Engrg. 99, 235–394. Bouchard, P.O., Baya, F., Chastela, Y., Tovena, I., 2000. Crack propagation modelling using an advanced remeshing technique. Comput. Methods Appl. Mech. Engrg. 189 (3), 723–742. Chang, S.H., Choi, S.W., Bae, G.J., Jeon, S., 2006. Performance prediction of tbm disc cutting on granitic rock by the linear cutting test. Tunn. Undergr. Space Technol. 21 (3), 271. Charles, E., Anderson Jr., 1987. An overview of the theory of hydrocodes. Int. J. Impact Eng 5, 33–59. Chiaia, B., 2001. Fracture mechanisms induced in a brittle material by a hard cutting indenter. Int. J. Solids Struct. 38, 7747–7768. Cho, J.W., Jeon, S., Yu, S.H., Chang, S.H., 2010. Optimum spacing of tbm disc cutters: a numerical simulation using the three-dimensional dynamic fracturing method. Tunn. Undergr. Space Technol. 25 (3), 230–244. Cook, N.G.W., Hood, M., Tsai, F., 1984. Observations of crack growth in hard rock loaded by an indenter. Int. J. Rock Mech. Min. Sci. Geomech. Abstr. 21 (2), 97– 107. CSM, 2007. CSM Computer Model for TBM Performance Prediction
103
Grandori, C., 1989. Development and current experiences with double-shield TBMs in Italy. Tunn. Undergr. Space Technol. 4 (3), 333–336. Gong, Q.M., Zhao, J., Jiao, Y.Y., 2005. Numerical modeling of the effects of joint orientation on rock fragmentation by TBM cutters. Tunn. Undergr. Space Technol. 20 (2), 183–191. Gong, Q.M., Zhao, J., Jiao, Y.Y., 2006. Numerical modeling of the effects of joint spacing on rock fragmentation by TBM cutters. Tunn. Undergr. Space Technol. 21 (1), 46–55. Gong, Q.M., Du, X., Li, Z., Wang, Q.X., 2016. Development of a mechanical rock breakage experimental platform. Tunn. Undergr. Space Technol. 57, 129–136. Hans, U.M., 1999. Review: hydrocodes for structure response to underwater explosions. Shock Vib. 6 (2), 81–96. Hassanpour, J., Rostami, J., Zhao, J., 2011. A new hard rock tbm performance prediction model for project planning. Tunn. Undergr. Space Technol. 26 (5), 595–603. Jiang, Q., Su, G.S., Feng, X.T., Cui, J., Pan, P.Z., Jiang, J.Q., 2015. Observation of rock fragment ejection in post-failure response. Int. J. Rock Mech. Min. Sci. 74, 30–37. Li, X.F., Li, H.B., Liu, Y.Q., Zhou, Q.C., Xia, X., 2016. Numerical simulation of rock fragmentation mechanisms subject to wedge penetration for TBMs. Tunn. Undergr. Space Technol. 53 (5), 96–108. Liang, Z.Z., Xing, H., Wang, S.Y., Williams, D.J., Tang, C.A., 2012. A three-dimensional numerical investigation of the fracture of rock specimens containing a preexisting surface flaw. Comput. Geotech. 45, 19–33. Libersky, L.D., Petschek, A.G., 1991. Smoothed particle hydrodynamics with strength of materials. In: Trease, H., Friits, J., Crowley, W. (Eds.), Proceedings of the Next Free Lagrange Conference, vol. 395. Springer, New York, pp. 248–257. Liu, G.R., Quek, 2003. The Finite Element Method: A Practical Course. Butterworth Heinemann. Liu, H.Y., Kou, S.Q., Lindqvist, P.A., Tang, C.A., 2002. Numerical simulation of the rock fragmentation process induced by indenters. Int. J. Rock Mech. Min. Sci. 39, 491–505. Lucy, L.B., 1977. Numerical approach to testing the fission hypothesis. Astronom. J. 82, 1013–1024. Ma, G.W., Wang, X.J., Ren, F., 2011. Numerical simulation of compressive failure of heterogeneous rock-like materials using SPH method. Int. J. Rock Mech. Min. Sci. 48 (3), 353–363. Maeda, M., Kushiyama, K., 2005. Use of compact shield tunneling method in urban underground construction. Tunn. Undergr. Space Technol. 20 (2), 159–166. Ma, H.S., Gong, Q.M., Wang, J., Yin, L.J., Zhao, X.B., 2016. Study on the influence of confining stress on tbm performance in granite rock by linear cutting test. Tunn. Undergr. Space Technol. 57, 145–150. Menezes, P.L., Lovell, M.R., Avdeev, I.V., Fred Higgs III, C., 2014. Studies on the formation of discontinuous rock fragments during cutting operation. Int. J. Rock Mech. Min. Sci. 71, 131–142. Mori, K., Abe, Y., 2005. Large rectangular cross-section tunneling by the multi-micro shield tunneling (MMST) method. Tunn. Undergr. Space Technol. 20 (2), 129– 141. Paltrinieri, E., Sandrone, F., Zhao, J., 2016. Analysis and estimation of gripper tbm performances in highly fractured and faulted rocks. Tunn. Undergr. Space Technol. 52, 44–61. Pang, S.S., Goldsmith, W., 1990. Investigation of crack formation during loading of brittle rock. Rock Mech. Rock Eng. 23, 53–63. Pramanik, R., Deb, D., 2015. Implementation of smoothed particle hydrodynamics for detonation of explosive with application to rock fragmentation. Rock Mech. Rock Eng. 48 (4), 1683–1698. Qing, Y., Van, J., Leigh, M., 2012. Free-surface flow interactions with deformable structures using an SPH-FEM model. Ocean Eng. 55, 136–147. Robbins, R.J., 1987. Machine tunnelling in the twenty-first century. Tunn. Undergr. Space Technol. 2 (2), 147–154. Sauer, M., 2000. Adaptive Koppling Des Netzfreien SPH-Verfahrens Mit Fniten Elementen Zur Berechnung Von Impaktvorgaengen. University of the Federal Armed Forces, Munich. Tóth, A., Gong, Q.M., Zhao, J., 2013. Case studies of tbm tunneling performance in rock–soil interface mixed ground. Tunn. Undergr. Space Technol. 38 (9), 140– 150. Zhang, Y.Q., Lu, Y., Hao, H., 2004. Analysis of fragment size and ejection velocity at high strain rate. Int. J. Mech. Sci. 46 (1), 27–34. Zhao, J., 2007. Tunnelling in rock-present technology and future challenges. In: Proceedings of World Tunnel Congress, Prague, pp. 22–32. Zhao, X.B., Yao, X.H., Gong, Q.M., Ma, H.S., Li, X.Z., 2015. Comparison study on rock crack pattern under a single normal and inclined disc cutter by linear cutting experiments. Tunn. Undergr. Space Technol. 50, 479–489. Zienkiewicz, O.C., Taylor, R.L., 2000. The Finite Element Method. Butterworth Heinemann.