CASE STUDY ON OPTIMIZING DISC CUTTER SPACING BASED ON VERTICAL FORCE THRESHOLD IN COLORADO RED GRANITE CUTTING USING LS-DYNA SOFTWARE

1. Introduction

Full-face tunnel boring machines (TBMs) have become a primary choice for modern tunnel construction across sectors like transportation and infrastructure. Among their components, the cutter head plays a pivotal role, where its design directly affects operational efficiency. Key design factors include the diameter of the cutters, the total number of disc cutters, applied vertical and torque forces, penetration depth, and the spacing arrangement between adjacent cutting discs.

Among these parameters, shear forces, spacing, and penetration depth are key factors affecting cutting efficiency (Cho et al., 2010). The cutting disc is the primary tool used for rock fragmentation in full-face tunnel boring machines. Therefore, analysing the forces acting on the cutting disc and examining the effects of penetration depth and spacing is essential for optimizing cutter head design and improving cutting efficiency (Ma et al., 2024; Wang et al., 2021).

To determine the aforementioned forces, full-scale linear cutting tests (LCM) can be employed, which is a reliable method that many researchers have utilized in their studies (Chang et al., 2006; Cho et al., 2013; Gertsch et al., 2007; Rostami and Ozdemir, 1996).

The most significant research conducted using this method pertains to Rostami in the years 1993 and 1997, who presented a set of relationships for estimating the forces applied on the shear disc of a constant cross-section. His primary relationship is based on the distribution of various pressures in the contact area between the shear disc and the rock for different cross-sectional profiles of the shear disc. This research was further extended through laboratory tests to develop a relationship with equivalent dimensions for estimating shear forces (Rostami, 2013; Rostami and Ozdemir, 1993).

Linear shear testing, despite its numerous advantages, requires substantial cost and time for each sample (Kazemi et al., 2024). In addition to theoretical modelling, recent experimental and numerical investigations have demonstrated the significant influence of operational and geomechanical parameters on excavation performance and tool wear. Using a custom-designed TBM laboratory simulator, it has been shown that factors such as cutter head rotation speed, penetration rate, and soil moisture content can drastically alter the wear rate of cutting tools and overall energy consumption during tunnelling (Mousapour et al., 2023). Integrated geomechanical analyses and slope stability modelling, in surface coal mining operations have further emphasized the role of accurate input parameters in ensuring design safety and compliance with geotechnical standards (Ahmeti and Maliqi, 2023). These studies highlight the growing importance of advanced simulation tools—such as LS-DYNA—in capturing complex interactions between excavation tools and geological media under varying site conditions.

To address this issue, various predictive models for the performance of full-face tunnel boring machines have been developed, including the Colorado School of Mines (CSM) model and the Norwegian University of Science and Technology (NTNU) model (Bruland, 1998; Rostami and Ozdemir, 1993). Previous research has produced numerous theoretical and semi-theoretical formulas; however, all these analytical models can only provide an approximate estimation of the forces acting on the disc. For instance, Foul, and Roxborough and Phillips examined shear force based on rock resistance and penetration, while Wang, Sanyo, and Farrokh focused solely on the rate of penetration (Farrokh et al., 2012; Fowell, 2016; Roxborough and Phillips, 1975; Sanio, 1985; Wang and Ozdemir, 1978). The shear process was optimized using the shear disc through Particle Swarm Optimization (PSO) methods, employing parameters such as RMR, RQD, and compressive and shear strength (Zhou et al., 2020).

Linear shear testing, in conjunction with numerical modelling, represents one of the best methods for predicting the performance of shear discs. In this approach, linear shear testing is utilized as a means of validating the results of modelling, subsequently leading to the development of models for other scenarios. Numerous researchers have employed numerical modelling to simulate rock shearing, with most opting for 2D modelling using the finite element method with Lagrangian algorithms to reduce computational costs and time. A challenge of this type of modelling, aside from its lack of three-dimensionality, is defining the erosion number for removing elements from the rock that have reached a certain strain level to prevent the divergence of the problem (Cardu et al., 2021; Fattahi et al., 2024).

Given the objectives outlined, the aim of this project is to investigate the effects of penetration depth and spacing on the vertical and rotational forces applied to a shear disc using the Smoothed Particle Galerkin (SPG) method and LS-DYNA software, focusing on maximizing the shear disc’s resistance to the applied vertical force. The advantages of this approach include reduced costs and time, the elimination of the need for repeated experiments upon changes in any variable, and the potential for generalization to various disc and rock types. Furthermore, many researchers, such as Chen et al. and Gertsch et al., have utilized special energy to attain optimal spacing, while Yadong Geo has investigated the effects of slope and direction in addition to spacing (Gertsch et al., 2007; Xue et al., 2021; Xu et al., 2022). The influence of a shear disc traversing a new path or crossing a shear zone created by a previous disc, all employing special energy that is considerably more complex to calculate than vertical force, has been analysed (Li et al., 2022). It is important to note that special energy refers to the amount of force expended to shear a unit volume of rock.

To ensure the accuracy of the simulations, the rock sample used in this study is the red granite from Colorado, as utilized by Gertsch et al. (2007) The uniaxial compressive strength and shear resistance of this rock were determined through simulations of uniaxial compressive and Brazilian tests, subsequently compared with laboratory values. Following this, a simulation of the 17-inch A30581 disc was conducted in accordance with the tests performed by Gertsch et al. After achieving highly accurate results, further modelling of the same shear disc was carried out with varying penetration depths and spacings.

In this study, the Colorado Red Granite was selected as the reference rock type due to its well-documented mechanical behaviour and extensive prior use in experimental investigations. This granite originates from the Beulah region in the Wet Mountains of southern Colorado, United States (approximate coordinates: 38.08°N, 104.98°W), an area known for its homogeneous granitic formations of Precambrian age. The geological and mechanical consistency of this granite has made it standard material in linear cutting experiments, as previously employed in the study conducted by Gertsch et al. (2007), which forms the experimental basis for the validation of the present numerical simulations.

2. Structural Model of Rock

The Colorado Red Granite utilized in this study is a coarse-grained, porphyritic granite of Precambrian age, geologically located in the Wet Mountains region of southern Colorado, USA. The area is part of the southern Rocky Mountains and comprises crystalline basement rocks formed during the Proterozoic Eon, primarily as a result of magmatic intrusions associated with regional tectonic activity. These intrusions underwent slow cooling at depth, resulting in the formation of granitic bodies with well-developed crystalline structures. The granite exhibits a high degree of mineralogical homogeneity and structural integrity, making it suitable for controlled rock mechanics studies. Its evolutionary history is marked by multiple phases of tectonothermal activity, which contributed to the formation of its present-day mechanical and petrographic characteristics (Gertsch et al., 2007; Tweto, 1979; Wobus and Hedge, 1982).

The Johnson–Holmquist (JH) model was originally formulated to represent the dynamic response of brittle materials like concrete under conditions involving elevated pressures and strain rates (Holmquist and Johnson, 2011). A key feature of this model is its ability to differentiate between material behaviour in tension versus compression. By establishing a relationship between internal pressure and deviatoric stress in the damaged state, the JH model effectively estimates the material’s resistance characteristics under dynamic loading.

[TeX:] \frac{\tau}{f_{c}} = \left\lbrack A\left( 1 - D \right) + B\left( p^{*} \right)^{N} \right\rbrack\left\lbrack 1 + C\ln\left( {\dot{\varepsilon}}^{*} \right) \right\rbrack$\frac{\tau }{f_c}=[A(1-D)+B{\left(p^{*}\right)}^N][1+Cln\left({\dot{\epsilon }}^{*}\right)]$ (1)

Where:

Shear stress (τ) behaviour under tensile and compressive conditions can be analysed using the parameter Ks, which typically ranges from 0.778 to 1.0 to maintain the convex nature of the yield surface. Setting Ks equal to 1.0 removes the influence of the third invariant of shear stress, effectively simplifying the model back to the classical Drucker–Prager formulation. The corresponding yield surface for this condition is depicted by Equation 1 in Figure 1 (Hibbitt and Sorensen, 2005).

Figure 1. the adopted failure criterion of the material model in two projections: (a) the meridional plane and (b) the deviatoric plane (Hibbitt and Sorensen, 2005).

The hydrostatic pressure within the model is determined through an equation of state, structurally similar to that introduced by Holmquist et al., which characterizes the relationship between pressure and volumetric changes (Hibbitt and Sorensen, 2005). Moreover, the volumetric strain (which is dimensionless) is defined based on the principal strain components and can be represented as follows:

[TeX:] \mu = \frac{V - V_{0}}{V_{0}} = \left( 1 + \varepsilon_{11} \right)\left( 1 + \varepsilon_{22} \right)\left( 1 + \varepsilon_{33} \right) - 1 \approx \varepsilon_{11} + \varepsilon_{22} + \varepsilon_{33}$\mu =\frac{V-V_0}{V_0}=(1+{\epsilon }_{11})(1+{\epsilon }_{22})(1+{\epsilon }_{33})-1\approx {\epsilon }_{11}+{\epsilon }_{22}+{\epsilon }_{33}$ (2)

Here, [TeX:] V_{0}$V_0$ and V denote the initial and deformed volumes of the material, respectively. When subjected to high hydrostatic pressures, the volumetric strain behaviour of rocks can be divided into three distinct zones, as illustrated in Figure 2. In this context, [TeX:] P_{\text{crus}h}(MPa)$P_{\mathrm{\text{crus}}h}(MPa)$ represents the pressure at which the pores are completely compacted, while [TeX:] \mu_{\text{crus}h}${\mu }_{\mathrm{\text{crus}}h}$ corresponds to the volumetric strain within the elastic threshold. During uniaxial compression experiments, [TeX:] P_{\text{lock}}(MPa)$P_{\mathrm{\text{lock}}}(MPa)$ and [TeX:] \mu_{\text{lock}}${\mu }_{\mathrm{\text{lock}}}$ mark the pressure and strain associated with the onset of material compaction. To quantify the evolution of plastic strain, a non-associated plastic flow rule is applied, with the plastic potential function g formulated as follows:

[TeX:] g = \sqrt{3J_{2}}$g=\sqrt{3J_2}$ (3)

Figure 2. The relationship between pressure and volumetric strain (Hibbitt and Sorensen, 2005).

The amount of plastic strain is described as follows:

[TeX:] d\varepsilon_{\text{ij}}^{p} = \dot{\lambda}\frac{\partial g}{\partial\sigma_{\text{ij}}}$d{\epsilon }_{\mathrm{\text{ij}}}^p=\dot{\lambda }\frac{\partial g}{\partial {\sigma }_{\mathrm{\text{ij}}}}$ (4)

Where [TeX:] \dot{\lambda}$\dot{\lambda }$ is referred to as plastic strain.

Fracture-induced damage is modelled through a cumulative approach analogous to that of the original formulation. This accumulated damage, derived from both equivalent plastic strain and volumetric plastic strain, is quantified using a compressive damage parameter as defined by Hibbitt and Sorensen (2005).

[TeX:] D = \sum_{}^{}\frac{\Delta\varepsilon_{\text{eff}}^{p} + \Delta\mu_{p}}{\varepsilon_{p}^{f} + \mu_{p}^{f}}$D={\sum }_{}^{}\frac{\Delta {\epsilon }_{\mathrm{\text{eff}}}^p+\Delta {\mu }_p}{{\epsilon }_p^f+{\mu }_p^f}$ (5)

Where [TeX:] \Delta\varepsilon_{p}$\Delta {\epsilon }_p$ and [TeX:] \Delta\mu_{p}$\Delta {\mu }_p$ represent the equivalent plastic strain and volumetric plastic strain, respectively, and denotes the plastic strain required for failure under constant pressure p, which is expressed as follows:

[TeX:] \varepsilon_{p}^{f} + \mu_{p}^{f} = D_{1}\left( p^{*} + T^{*} \right)^{D_{2}} \geq \varepsilon_{\text{fmin}}${\epsilon }_p^f+{\mu }_p^f=D_1{(p^{*}+T^{*})}^{D_2}\ge {\epsilon }_{\mathrm{\text{fmin}}}$ (6)

Where and are constants, and is the constant pressure as previously defined. This is designated to generate a limited amount of plastic strain for the failure of materials (Hibbitt and Sorensen, 2005). The parameters extracted by Hyun Lee et al., as shown in Table 1, have been utilized in the simulations for this structural model (Li and Du, 2016).

Table 1. Parameters of Colorado Red Granite Used in Numerical Analysis with the Johnson-Holmquist Material Model (Li and Du, 2016).

3. Validation of Structural Model Parameters

The calibration of the structural model for accurately simulating the behaviour of rocks requires the execution of a series of complex and costly experiments. These experiments include uniaxial and triaxial compression tests, Brazilian tests, point load tests, and even more specialized tests, such as the split-Hopkinson pressure bar test. To reduce costs and save time, all parameters utilized in this study have been sourced from the research conducted by Hyun Lee et al., which comprehensively addresses the calibration of numerical models and provides various parameters from compressive and tensile tests for different types of rocks under various conditions (Li and Du, 2016).

To ensure validity, two fundamental tests uniaxial and Brazilian were first simulated. These tests serve as key tools for evaluating and validating the numerical models used in this research. Initially, the uniaxial test was simulated to measure the compressive strength of Colorado red granite samples, and the results were then compared with experimental data. Subsequently, the Brazilian test was conducted to measure the indirect tensile strength of the samples.

3.1. Simulation of the Uniaxial Test

The uniaxial test is one of the primary and widely used experiments for assessing the compressive strength of materials, particularly in rock and concrete samples. This test is designed to measure a material's ability to withstand compressive loads without failure, and it is extensively employed in civil engineering and mining research. In the simulations conducted for this study, a Colorado red granite sample with a diameter of 54 mm and a height of 116 mm was selected. These dimensions conform to the general standards for uniaxial compression tests for construction materials (Bieniawski, 2020).

During the simulation process, the load is applied uniformly through the upper jaw of the testing apparatus to the rock sample, while the lower jaw remains fixed. This loading condition closely resembles laboratory settings and is essential for accurately assessing the compressive strength of the rock. Figure 3 illustrates the geometry designed in the LS-PREPOST software, in which the geometry of the rock sample and the loading conditions have been fully modelled and simulated.

Figure 3. View of the Geometry of Uniaxial Testing.

Based on numerous uniaxial tests conducted on Colorado red granite, the uniaxial compressive strength of this rock has been determined to be 158 megapascals. This value was calculated to be 162 megapascals in simulations carried out using calibrated coefficients referenced in the research conducted by Li and Du (2016) and Gertsch et al. (2007). This value corresponds to the peak point in the graph illustrated in Figure 4. The error obtained from comparing the simulation results with laboratory data is approximately 2.5%, which is quantitatively low and indicates a reliable accuracy of the model in the simulated model and its suitable alignment with experimental reality. Figure 5 presents a view of the simulated rock sample after the completion of the simulation process, effectively demonstrating the changes and failures that occurred within the sample. This simulated image clearly illustrates how the rock behaved under compressive loading in the numerical model and highlights the potential failures that occurred.

Figure 4. The stress (MPa)- displacement (mm) curve obtained from the simulated uniaxial test.

Figure 5. A view of the rock sample following the completion of the uniaxial test simulation.

3.2. Brazilian Test Simulation

The Brazilian test is a widely used method for indirectly measuring the tensile strength of rock and concrete specimens. This test is particularly relevant in rock engineering, as rock samples are prone to sudden failure when subjected to direct tensile forces. During the test, the specimen is subjected to compressive loading, with the standard requirement that the diameter of the sample should be twice its length.

In the simulation conducted for this research, the specimen's diameter was set at 54 mm, and its length was 27 mm. The geometry of the simulated Brazilian test is illustrated in Figure 6, which depicts the manner in which the load is applied and the stress distribution within the sample.

Figure 6. The geometry of the Brazilian test on Colorado Red Granite.

In the conducted experiments, the tensile strength of Colorado red granite was calculated to be 6.8 megapascals. This value was estimated to be 7.7 megapascals in the simulations performed. The tensile strength corresponds to the maximum value indicated in the graph shown in Figure 7. Additionally, a representation of the Colorado red granite sample is illustrated in Figure 8.

Figure 7. The stress (MPa)-displacement (mm) graph obtained from the Brazilian test simulation.

Figure 8. A view of the simulated rock specimen following the completion of the Brazilian test simulation.

Based on the results from the uniaxial and Brazilian test simulations, utilizing the calibrated values provided by Hyun Lee et al. the errors in the compressive and tensile strengths calculated through LS-DYNA simulations are less than 2.5% and 13%, respectively. This level of error indicates a high degree of accuracy in the calibrated structural model and demonstrates that the simulation results align satisfactorily with experimental values. This correlation enhances the credibility of numerical models and facilitates their application in more complex simulations, such as linear shear simulations of rocks using a cutting disk.

3.2. Petrophysical Properties of Colorado Red Granite

While Table 1 presents the numerical parameters used in the Johnson–Holmquist material model, the fundamental petrophysical properties of the Colorado Red Granite are summarized in Table 2. These values have been obtained from experimental investigations and literature sources and are used to validate the material behaviour simulated in this study. Presenting these properties independently offers a broader perspective on the physical characteristics of the rock, beyond model-specific parameters.

Table 2. Petrophysical Properties of Colorado Red Granite (Gertsch et al., 2007; Li and Du, 2016).

4. Characteristics of Linear Shear Modelling

LS-DYNA offers multiple finite element solvers tailored for dynamic simulations, including the Lagrangian method, Smoothed Particle Hydrodynamics (SPH), and the Smoothed Particle Galerkin (SPG) approach. The SPG method, an advanced evolution of SPH, is grounded in Galerkin-based formulations, and both techniques are categorized as mesh-free methods (Wu et al., 2015). According to Rokhy et al., in their investigation of rock shearing using LS-DYNA, mesh-free methods like SPH and SPG delivered more accurate outcomes compared to traditional Lagrangian approaches, primarily due to their independence from erosion-related limitations (Bhat and Maji, 2022; Rokhy et al., 2021).

In the present study, the Lagrangian solver was selected as the main computational framework for simulating uniaxial and Brazilian tests, with the intent of reducing overall computational demands. This method was also used for linear shear analyses involving the cutter disk and the surrounding rock segments not directly affected by cutting. Conversely, the SPG solver was selectively employed in the zone directly influenced by the cutter, defined by a depth equal to ten times the penetration depth and a width approximately three to four times the spacing of the cutter. The complete shear modelling arrangement is depicted in Figure 9.

Figure 9. Overview of linear shear modelling of rock in LS-DYNA.

In the SPG simulations, the particle spacing in all three directions (DX, DY, DZ) was set to 0.3 mm, which ensured sufficient resolution of crack initiation and propagation while maintaining computational efficiency. The time step was automatically controlled in LS-DYNA with a stability scale factor of 0.67, following the recommendations in the LS-DYNA manual. Boundary conditions were applied such that the movement of the rock specimen was restricted at the bottom and lateral faces, consistent with the configuration of the uniaxial compression test shown in Figure 10, while the upper surface was free to interact with the cutter. For the kernel function, the Semi-pseudo Lagrangian kernel was employed, as it provides stable particle interaction and accurate stress transfer within SPG formulations. These settings ensured stable and reproducible simulations across all penetration depths and cutter spacings investigated. All simulations were executed on a workstation with an Intel Core i9-12900K CPU and 64 GB RAM (no GPU). With 0.3 mm particle spacing and a stability scale factor of 0.67, each linear-cutting case completed in approximately 24–36 hours of wall-clock time, depending on penetration depth and cutter spacing, supporting the practical computational efficiency of the adopted setup.

Figure 10. Boundary conditions governing the problem implemented in LS-DYNA software

The disk employed in all simulations is the same 17-inch A30581 disk used by the Colorado School of Mines for linear shear testing. A representation of this 17-inch shear disk is shown in Figure 11.

Figure 11. View of the 17-inch A30581 disk used at the Colorado School of Mines for linear shear tests on Colorado red granite (Gertsch et al., 2007).

5. Validation of Linear Cutting Simulations

The results obtained from validation simulations are of significant importance. The closer these simulations results are to the values derived from laboratory experiments and field results, the greater the accuracy and reliability of the simulations that lack corresponding experimental and field data. Consequently, every effort has been made to ensure that all relevant conditions governing the experimental setup, including boundary conditions, loading rates, disk movement velocities, and all other parameters, are meticulously adhered to during these simulations.

To validate the simulations and determine the error margin, twelve different simulations were conducted in accordance with the linear cutting tests performed by Gertsch et al. (2007), under the following conditions:

- Rock strength: 158 MPa,

- Disk rotational speed: 0.0015 radians per millisecond,

- Disk linear cutting speed: 0.33 millimetres per millisecond,

- Penetration depths: 2.5, 5.1, and 6.4 millimetres,

- Spacings: 25, 51, and 76 millimetres.

Given that three disks with specified spacings were utilized, the forces recorded in all sections pertain to the central disk. Due to the significant fluctuations in forces reported from the disks at both ends of the rock—attributable to the boundary conditions—it is preferable to consider the average forces in the middle of the rock, particularly during the time frame when the cutting disk is approximately positioned in the middle of the modelled rock (with a length of 304 millimetres), which occurs between 400 and 600 milliseconds. At this point, boundary conditions are minimal, and the forces are in equilibrium. This average is taken as the effective vertical and rotational forces exerted by the rock on the cutting disk. For instance, Figures 12 and 13 illustrate the vertical and rotational forces obtained using the Lagrangian-SPG coupled solution method and the Johnson-Holmquist concrete behaviour model, based on the simulations with a spacing of 76 millimetres and a penetration depth of 6.4 millimetres. Figure 14 presents a view of the simulated linear cutting of the rock.

Given that three disks with specified spacings were utilized, the forces recorded in all sections pertain to the central disk. Due to the significant fluctuations in forces reported from the disks at both ends of the rock—attributable to the boundary conditions—it is preferable to consider the average forces in the middle of the rock, particularly during the time frame when the cutting disk is approximately positioned in the middle of the modelled rock (with a length of 304 millimetres), which occurs between 400 and 600 milliseconds. At this point, boundary conditions are minimal, and the forces are in equilibrium. This average is taken as the effective vertical and rotational forces exerted by the rock on the cutting disk.

Figures 12 and 13 illustrate the vertical and rotational forces, based on the simulations with a spacing of 76 millimetres and a penetration depth of 6.4 millimetres. In both figures, the forces show a rapid increase during the initial engagement of the disks with the rock, followed by strong fluctuations before reaching a relatively stable trend between 400 and 600 milliseconds. These oscillations indicate progressive crack initiation and interaction within the crushed zone, which is consistent with laboratory observations. Figure 14 presents a view of the simulated linear cutting of the rock.

Figure 12. Simulated vertical force of the disk at a penetration depth of 6.4 millimetres and a spacing of 76 millimetres.

Figure 13. Simulated rotational force of the disk at a penetration depth of 6.4 millimetres and a spacing of 76 millimetres.

Figure 14. View of the simulated linear cutting of the rock.

A comparison of the results in Table 3 and Figures 15 and 16 shows that the simulated vertical forces deviate less than 10% from experiments, while rotational forces differ by about 13%. In both figures, vertical and rotational forces increase consistently with penetration depth and spacing. However, the error of the vertical force does not follow a fixed trend relative to the experimental results, indicating that local fracture mechanisms and disc–disc interactions influence the accuracy of the model. In contrast, the error of the simulated rotational forces generally decreases with penetration depth, except for the case of 51 mm spacing at 3.8 mm penetration depth, where the deviation is slightly higher than at 1.9 mm.

Table 3. Results from the simulations and linear cutting experiments conducted by Gertsch et al. (2007)

Figure 15. The vertical force obtained from the simulation compared to the experimental results for three spacings of 25, 51, and 76 millimetres and three penetration depths of 1.9, 3.8, and 6.4 millimetres.

Figure 16. The rotational force obtained from the simulation alongside the experimental results for the same three spacings of 25, 51, and 76 millimetres and three penetration depths of 1.9, 3.8, and 6.4 millimetres.

6.1. Investigation of Optimal Spacing

There are various perspectives regarding the examination of optimal spacing; however, this research references the findings of Roby, which indicate that the maximum vertical force a 17-inch disk can withstand is 215 kN (Roby et al. 2008). According to Farokh, to enhance the service life of shear disks, the maximum vertical force applied to these disks should not exceed 80% of this value, equating to 180 kN (Farrokh et al., 2012). Consequently, based on the results obtained from preliminary modelling presented in Table 1, it was observed that the vertical force exerted on the shear disk did not reach this threshold at any of the penetration depths and spacings examined. Therefore, a penetration depth of 7.6 mm and a spacing of 100 mm were also investigated, with the results detailed in Table 4.

Table 4. Results from the Simulation of Shear Disk at Increased Penetration Depth and Spacing.

The results obtained from the simulations indicate that at a spacing of 100 mm, since the maximum vertical force exerted on the shear disk must equal 80% of the maximum vertical force the disk can withstand (180 kN), the penetration depth cannot exceed 5.1 mm. This limitation on increasing the penetration depth leads to a reduction in the penetration rate and efficiency of the tunnel boring machine. Additionally, as shown in Figures 12 and 13, at a spacing of 76 mm, the two side shear disks influence the central disk, resulting in lower vertical and rotational forces on this disk compared to the two side disks. However, at a spacing of 100 mm, for all penetration depths examined, the disks do not influence each other, resulting in a single-disk cutting scenario (where the cutting forces on the central disk equal those on the side disks), as depicted in Figures 17 and 18 for vertical and rotational forces, respectively. This condition leads to significantly higher forces exerted on the disks, which may prevent crack connectivity and the formation of a crushed zone in the rock. This behaviour highlights that the disks’ interaction is essential for crack coalescence and efficient rock fragmentation, which explains why 76 mm emerges as the optimal spacing.

Based on the results obtained and their comparison with the findings of Gertsch and colleagues, the optimal spacing of the disks is approximately 76 mm. This is because at spacings less than this, the force exerted on the rock, which is equal to the force exerted on the shear disk in the opposite direction, is significantly less than the shear disk's resistance. Conversely, at spacings greater than 76 mm, such as 100 mm with a penetration depth of 5.1 mm, the vertical force exerted on the shear disk reaches the critical region (80% of the maximum resistance of the disk against vertical force). Therefore, it is not feasible to increase the penetration depth of the disks to enhance the penetration rate.

Figure 17. Simulated Vertical Force on the Disk at a Penetration Depth of 6.4 mm and Spacing of 100 mm.

Figure 18. Simulated Rotational Force on the Disk at a Penetration Depth of 6.4 mm and Spacing of 100 mm.

7. Conclusions

In this study, a comprehensive three-dimensional numerical model was developed using the mesh-free Smoothed Particle Galerkin (SPG) method in LS-DYNA to simulate the linear cutting behaviour of hard rock by disc cutters. The Johnson-Holmquist (JH2) constitutive model was used to accurately capture the dynamic response of Colorado Red Granite under disc cutting conditions.

To validate the model, numerical simulations of uniaxial compressive strength (UCS) and Brazilian tensile strength (BTS) tests were performed. The differences between numerical and laboratory results were 2.5% and 13%, respectively, demonstrating acceptable accuracy. Furthermore, the model was calibrated and validated against experimental results from Gertsch et al. (2007), who conducted full-scale linear cutting tests on Colorado Red Granite at the Colorado School of Mines.

Specifically, the average results of the nine initial simulations reported in Table 3—each corresponding to cutting conditions that were experimentally tested by Gertsch et al.—showed strong agreement with laboratory data. The average errors in vertical and rotational forces obtained from the simulations were 4.67% and 7.99%, respectively, relative to experimental values. This close match confirms the high reliability and validity of the proposed model. Accordingly, simulations performed under different spacing and penetration conditions (as presented in Table 4), for which experimental data are not available, can also be considered valid and trustworthy.

The model also introduced a new approach for optimizing disc cutter spacing using 80% of the disc’s vertical force capacity instead of specific energy, which is commonly used in literature but requires more complex calculations. The optimal disc spacing determined through this method was 76 mm, which closely aligns with the 80 mm optimal value obtained from Gertsch’s experiments—showing less than a 5% deviation. This demonstrates that vertical force threshold is a viable and simpler alternative for practical design purposes.

It should be noted, however, that the 80% vertical force threshold is a deliberately simplified criterion aimed at providing a practical and field-deployable guideline. While it captures the fundamental mechanical limit of the cutter, it does not explicitly consider additional operational factors such as wear progression, rock anisotropy, cutter-head vibration, or crack propagation dynamics. These aspects are highly relevant to long-term performance and will be incorporated in future extensions of this framework through multi-objective optimisation approaches.

Parametric analysis confirmed that increases in both spacing and penetration depth lead to significant increases in vertical and rotational forces, which is consistent with the trends reported in prior experimental and analytical studies. The regression models developed for force prediction also achieved high accuracy, with R² values exceeding 0.9.

While classical models such as those of Rostami (1993, 2013) and Gertsch et al. (2007) have provided valuable frameworks for analysing disc cutter forces, they are limited by the reliance on complex specific energy calculations and simplified analytical assumptions. The present study aimed to overcome these limitations in several ways. First, the simulations were performed using the Smoothed Particle Galerkin (SPG) method in LS-DYNA. Unlike conventional Lagrangian approaches, which require an erosion number to delete highly deformed elements—often violating mass and energy conservation—SPG eliminates this need and fully preserves conservation laws, resulting in more accurate and stable rock-cutting simulations. Second, instead of adopting the conventional specific energy criterion, a simpler and more practical approach was introduced based on 80% of the maximum allowable vertical force on the cutter. This criterion not only shows strong agreement with experimental results but is also much easier to obtain in laboratory and even field conditions, making it directly applicable to TBM design and operation. Third, in addition to validation against UCS, Brazilian, and Gertsch’s linear cutting data, a statistical regression model was developed to predict cutting forces rapidly and at low computational cost using spacing and penetration depth as input parameters. Collectively, these contributions move the present study beyond a reiteration of previous models and provide a framework that combines numerical accuracy with practical applicability for TBM cutter design.

In conclusion, the combination of advanced numerical modelling, experimental validation, and a simplified optimization criterion offers a powerful and practical tool for cutter design and TBM performance prediction. The model's accuracy, validated against reliable laboratory data, supports its application in various tunnelling and rock excavation scenarios.

Future studies may explore the extension of this model to simulate rock cutting under more complex conditions, such as the presence of joints, anisotropy, and heterogeneity within the rock mass. Investigating the influence of wear and thermal effects on cutter-rock interaction can also enhance the realism of simulations. Additionally, coupling the SPG method with real-time data from TBM operation could pave the way for intelligent prediction systems. Finally, validating the model against additional rock types and full-scale TBM tunnelling data can broaden its applicability and reliability in field-scale engineering problems.

Funding

No funding was received for conducting this study.

Author’s contribution

Shahin Fattahi (MSc. Student): formal analysis, investigation, methodology, and writing – original draft. Ebrahim Farrokh (Associate Professor): writing – review & editing, and validation. Hesam Dehghani (Associate Professor): investigation, methodology, supervision, and validation. Aref Fayazi (PhD Student): formal analysis and software.

All authors have read and agreed to the published version of the manuscript.