Vibration mechanism study based on coupled dynamic model for TBM main drive system

Abstract

Tunnel boring machine (TBM) is seriously damaged during tunneling due to excessive vibration. Insufficient research on the TBM vibration mechanism has led to low accuracy in the calculation of the dynamic response of the TBM main drive system under uncertain load. In this paper, multi-degree-of-freedom coupled dynamics model of the main drive system containing the cutterhead, gear drive, shield, and propulsion system was established, taking bearing and drive gear nonlinear contact forces, uncertain load into account. A multi-position vibration test was carried out based on practical engineering. The vibration mechanism and vibration transmission law along the main structure are studied. The results show that the frequency response of the cutterhead presented a wide frequency domain bandwidth. The main frequency band of support shell is concentrated at 18–20 Hz, which is mainly caused by the excitation of the first three-order natural frequencies of the system under the action of the system’s broadband external excitation load. At the same time, it is found that the coupled vibration of radial overturning vibration and lateral translation vibration is easy to appear in the actual system vibration. The influence of the main mass and stiffness parameters of the system on the various modes of the system is analyzed. The overturning vibration of the shield body has a significant “amplification effect,” and its amplitude reaches 123.60% and 118.99% of the cutterhead, respectively. This study provides a basis for the system’s dynamic design and evaluation of the TBM performance in tunnel construction.

Introduction

In the process of tunneling, TBM bears impact loads. At the same time, the internal subsystems such as cutterhead, transmission system, propulsion system, and support system interact with each other, which makes the dynamic characteristic of the system extremely complex. It is difficult to predict the vibration law of each key part during the excavation process, which directly affects the overall stability and excavation performance of the main drive system. This easily leads to a series of failures of key components of the TBM main drive (such as large-scale cracking of the cutter head body, fracture of connecting bolts, failure of bearing seals, etc.), ¹ resulting in serious engineering problems. Therefore, studying the dynamic characteristic and mechanism of the TBM main drive system under the coupling action of multiple subsystems is the basis for effective dynamic performance design.

The rock breaking load of TBM cutterhead provides external excitation for studying the dynamic characteristics and transmission laws of the main drive system, which is the basis for studying the dynamic characteristics of TBM main drive system. It directly affects the accuracy of system dynamic characteristic calculation. ² The research on the rock breaking load of rolling cutters mainly focuses on three directions: obtaining the load from laboratory rolling cutter rock breaking experiments, ³ numerical model simulation, ⁴ and on-site measurement of rolling cutter rock breaking load. ⁵ In order to solve severe vibration problem of key components of the TBM main drive system and reveal the complex dynamic behavior of the system under impact load, scholars mainly study by establishing dynamic model of the main drive system and analyzing the relevant dynamic responses. Scholars^6,7 have established a dynamic model of the TBM main drive system based on gear dynamics and bearing dynamics in order to explore the dynamic characteristics of internal vibration loads and other dynamic features of the system. Gear dynamics is an important part of the dynamic modeling of TBM main drive systems to study their load sharing and meshing characteristics under strong sudden loads. Kahraman and Singh^8,9 established a 3-degree-of-freedom bending torsional coupling dynamic model for gear pairs, taking into account the time-varying mesh stiffness and tooth flank clearance of the gear pairs, as well as the coupling relationship between the bending and torsional vibration directions of the gears, and analyzed their nonlinear vibration response. Based on this gear model, Sun et al.^10,11 established a multi-degree-of-freedom coupled dynamic model of the TBM cutterhead system, considering the dynamic meshing stiffness of gear teeth and bearing stiffness. After that, Zhang et al.^12,13 based on the above dynamic model, used the finite element method to model the main drive system drive shaft subsystem, and finally used a hierarchical dynamic modeling method to model each subsystem of the TBM. The model is more accurate in analyzing the dynamic characteristics of the synchronous drive transmission system. Since the support propulsion system plays a major role in supporting the TBM, its dynamic characteristics significantly affect the vibration and attitude of the main drive system of the TBM during the excavation process. Therefore, based on the dynamic model of the TBM drive system, scholars further established the dynamic model of the TBM main drive support propulsion system. Liao et al. ¹⁴ considered the closed-chain and kinematic relationship of the TBM support system, and established a dynamic model of the support propulsion system. By using the CSM formula to calculate the equivalent external load, the model can study and analyze the attitude of the system under different surrounding rock conditions and the support state of each cylinder. Huo et al.^15,16 established a multi-degree-of-freedom coupling model of the TBM main drive system from the cutter head to the main beam, and measured the vibration of the cutter head of the main drive system. However, only the main bearing direction of the axial direction is considered. Zou et al. ¹⁷ established a dynamic model for the open TBM main drive system. Based on this, the vibration energy of the system under different penetration degrees and rock types is studied. However, the model mainly considers the three-dimensional translation degrees of freedom of the system, so it is difficult for the model to study and analyze the unbalanced overturning vibration that the system bears during the actual excavation process.

Based on the above-mentioned dynamic model of the TBM’s main drive system, scholars have conducted a lot of research on the inherent characteristics of the system and the dynamic characteristics under impact load. Ling et al. ¹⁸ obtained the modal energy distribution of each system on the basis of a multi-degree-of-freedom coupled dynamics model, further distinguished the modal vibrations of each order, identified the modal sensitive parameters. Wu et al. ¹⁹ analyzed the load characteristics of the main components of the system based on an equivalent dynamic model of the TBM main drive system. Also Huo et al. ²⁰ further used interval uncertainty theory to describe the uncertainty of strati-graphic parameters, combined with a dynamic model to analyze the vibration characteristics of TBMs in composite rock formation, and explored the influence of the ratio of soft and hard rock on its vibration characteristics. Xia et al. ²¹ studied the vibration frequency band of the support propulsion system based on a multi-body dynamic simulation platform. Then a field vibration of the TBM support propulsion system was measured. The results show that the vibration of its main beam part was dominated by axial vibration.

In summary, for the study of the dynamic behavior of the TBM main drive system, scholars have mainly studied and analyzed the inherent and dynamic characteristics of the system based on the drive system dynamic model. Existing dynamic models and vibration measurements are not systematic enough to determine the dynamics of the main drive system with reasonable accuracy. Existing studies can only roughly calculate the vibration level, and there is insufficient research on the vibration mechanism and variation law of the TBM system. This has led to difficulties in TBM systems’ dynamic design to effectively solve the severe system vibration.

The aim of this paper is to establish a multi-degree-of-freedom coupled dynamic model of the main drive system containing the cutterhead, gear drive, shield, and propulsion system based on the structural characteristics. Meanwhile, internal excitations such as bearing and drive gear nonlinear contact forces are also considered. The vibration characteristics of the main drive system were measured at multiple points of the whole machine to verify the accuracy of the dynamic model. The vibration transmission law and vibration mechanism of the dynamic system under typical external load are studied based on the dynamic model.

Methods and materials

Modeling of TBM main drive system dynamics

A typical TBM main drive system is shown in Figure 1. It consists of a cutter system, main drive system, shield system, and support propulsion system. The main drive system and the support propulsion system are cemented together by flanges at the front section of the main beam, and the support propulsion system provides axial support for the entire main drive system. Through the coordinated action of the main drive system and the support propulsion system, the TBM completes the rock-breaking action.^22,23

Figure 1.

The structure diagram of TBM main system.

Equivalent dynamic model of TBM main drive system

Based on the multi-source drive mode and support structure of the TBM main drive system, the multi-degree-of-freedom coupling dynamic model of the TBM main drive system is established, as shown in Figure 2. The model mainly considers the translational and torsional degrees of freedom in the axial, horizontal, and vertical directions of the system. The model is established by the concentrated mass method, and the particle division mainly includes the cutter head-large ring gear, the support shield, the driving gears, the front section of the main beam, the middle section of the main beam, and the saddle. The mass points of the system are related to each other through nonlinear loads and spring damping elements. The large gear ring and the support shell are connected by bearing contact force. The driving gears and the large ring gear are connected by the meshing force of the gear teeth. The propulsion force between the main beam and the saddle is through the propulsion cylinder. The support shell and the surrounding rock are connected by contact force. The multi-source input load of the model is mainly reflected in two parts, one is the axial force, radial force, overturning moment and torque at the cutter head, and the other is the dynamic driving torque input by each drive motor. In the model, F_x, F_y, F_y are the horizontal, vertical, and axial equivalent loads on the cutter head. T_pi is the output torque of the ith motor. T_L is the equivalent total torque that the cutter head bears. M_x, M_y are the horizontal and vertical overturning moments of the cutter head. X_j, Y_j, Z_j, θ_j, θ_jx, θ_jy represent the horizontal displacement, vertical displacement, axial displacement, torsional displacement, horizontal overturning displacement, and vertical overturning displacement of the j component, respectively. m and I represent the equivalent mass and moment of inertia of each component. j = r, d, l1, l2, c, pi, mi represent the cutter head, the main drive support frame, the front section of the main beam, the middle section of the main beam, the ith drive gear of the saddle, and the ith drive motor, respectively.

Figure 2.

Dynamic model of TBM main drive system. (a) Main drive system dynamics model (horizontal direction). (b) Main drive system dynamics model (vertical direction). (c) Transmission system dynamics model (horizontal and axial directions). (d) Dynamic model of saddle and gripper shoe (axial direction).

The stiffness and damping of each component in the dynamic model are shown in Table 1. The dynamic parameters are calculated according to references. ¹¹

Table 1.

Main parameters of stiffness and damping in dynamic model.

Dynamic parameters	Symbols
Dynamic radial stiffness of main bearing	keqr
Radial stiffness of pinion	keqp
Radial stiffness of support shell	keqd
Axial connection stiffness between the support shell and the front section of the main beam	kadl 1i/cadl1i (i = 1–4)
Horizontal connection stiffness between the support shell and the front section of the main beam	khdl 1 /chdl 1
Vertical connection stiffness between the support shell and the front section of the main beam	kvdl 1 /cvdl 1
Left propulsion cylinder stiffness	kcy 1 /ccy 1
Stiffness of the left side gripper shoe cylinder	khcy 1 /chcy 1
Torsional stiffness of drive gear shaft	kmpQ
Dynamic meshing stiffness of gear pairs	kpri(t)
Axial connection stiffness between the cutterhead and the support shell	kardi/cardi (i = 1–4)
Axial connection stiffness of the front section of the main beam to the middle section of the main beam	kal 1l2i/cal1l2i (i = 1–4)
Horizontal connection stiffness of the front section of the main beam to the middle section of the main beam	khl 1l2/chl1l2
Vertical connection stiffness of the front section of the main beam to the middle section of the main beam	kvl 1l2 /cvl 1l2
Right propulsion cylinder stiffness	kcy 2 /ccy 2
Right propulsion cylinder stiffness	khcy 2 /chcy 2

The coupling relationship of each degree of freedom of the system

The relationship between the bending and torsional degrees of freedom of the large ring gear and the pinion

When the main drive system is working, each pinion synchronously drives the large ring gear to rotate. According to the structure of the main drive system, the drive mode can be divided into internal meshing form and external meshing form. When the gear drive is external meshing. The relative displacement relationship between the pinion gear and the large ring gear is shown in Figure 2(c).

The relative elastic deformation x_pri between the large ring gear and the pinion is:

$$\begin{array}{c}\multicolumn{1}{c}{x_{\mathrm{pr}j}=r_{b\mathrm{p}}{\theta }_{\mathrm{p}j}+r_{b\mathrm{r}}{\theta }_{\mathrm{r}}-V_{\mathrm{p}j}\sin \alpha +H_{\mathrm{p}j}\cos \alpha }\\ -X_{\mathrm{r}}\sin ({\phi }_j-\alpha )+Y_{\mathrm{r}}\cos ({\phi }_j-\alpha )-e_{\mathrm{pr}j}(t)\hfill \end{array}$$

(1)

Where u_r and u_pi are the displacements of the large ring gear and the pinion on the meshing line. X_r, Y_r, X_pi, and Y_pi are the horizontal and vertical displacements of the large ring gear and the pinion; α and φ_i are the meshing angle and phase angle of each small gear.

The equivalent meshing error of the large ring gear and the pinion on the meshing line is:

$$\begin{array}{c}\multicolumn{1}{c}{{\epsilon }_{\mathrm{pri}}(t)=E_r\sin ({\omega }_{r}t+\alpha +{\phi }_i-{\beta }_r)+A_r\sin (\alpha +{\phi }_i-{\gamma }_r)}\\ +E_{\mathrm{pi}}\sin ({\omega }_{p}t+\alpha -{\beta }_{\mathrm{pi}})+A_{\mathrm{pi}}\sin (\alpha -{\gamma }_{\mathrm{pi}})\hfill \\ +{\epsilon }_r+{\epsilon }_{\mathrm{pi}}+{\sigma }_r+{\sigma }_{\mathrm{pi}}(i=1-8)\hfill \end{array}$$

(2)

Where E_r, A_r, ε_r, σ_r are the manufacturing error, installation error, tooth thickness deviation, and tooth profile error of the large ring gear. E_pi, A_pi, ε_pi, σ_pi are manufacturing error, installation error, tooth thickness deviation, and tooth profile error of the pinion. β_r, γ_r, β_pi, γ_pi are the phase angles of each error. ω_r and ω_p are the angular velocity of the large ring gear and the pinion.

The time-varying meshing force F_pri and the time-varying meshing damping force D_pri of the gear are:

$$F_{\mathrm{pri}}=k_{\mathrm{pri}}(t)·f(x_{\mathrm{pri}},b)(i=1-8)$$

(3)

$$C_{\mathrm{pri}}=c_{\mathrm{pri}}(t)·x_{\mathrm{pri}}$$

(4)

where k_pri(t) and c_pri(t) are the time-varying meshing stiffness and damping of each pinion, f(x_pri, b) is the non-linear meshing clearance function.

$$\begin{array}{c}f(x_{\mathrm{pri}},b)=\{\begin{array}{c}{x}_{\mathrm{pri}}-b\begin{array}{c}{x}_{\mathrm{pri}}>b\hfill \end{array}\hfill \\ \multicolumn{1}{c}{0\begin{array}{c}-b<x_{\mathrm{pri}}<b\hfill \end{array}}\\ \multicolumn{1}{c}{x_{\mathrm{pri}}+b\begin{array}{c}{x}_{\mathrm{pri}}<-b\hfill \end{array}}\end{array}\hfill \end{array}$$

(5)

Where the backlash of the gear pair is 2b, and x_pri is the dynamic relative displacement of each gear pair along the meshing line.

In the dynamic model, δ_pdxi, δ_pdyi, δ_dpxi, δ_dpyi are the horizontal and vertical displacements of the i-th drive gear and the support frame, respectively. The calculation method is shown in the equation (6).

$$\begin{array}{c}\{\begin{array}{c}{\delta }_{\mathrm{pdxi}}=X_d+Y_{\mathrm{pi}}\cos {\phi }_i+X_{\mathrm{pi}}\sin {\phi }_i\hfill \\ \multicolumn{1}{c}{{\delta }_{\mathrm{dpxi}}=X_{\mathrm{pi}}-Y_d\cos {\phi }_i+X_d\sin {\phi }_i}\end{array}\hfill \\ \{\begin{array}{c}{\delta }_{\mathrm{pdyi}}=Y_d+Y_{\mathrm{pi}}\sin {\phi }_i-X_{\mathrm{pi}}\cos {\phi }_i\hfill \\ \multicolumn{1}{c}{{\delta }_{\mathrm{dpyi}}=Y_{\mathrm{pi}}+Y_d\sin {\phi }_i+X_d\cos {\phi }_i}\end{array}\hfill \end{array}$$

(6)

Relative displacement relationship of propulsion system

The bending, axial, and overturning coupling relationship of the support system is mainly used to describe the relative positional relationship between the front section of the main beam and the gripper shoe, and the relative displacement relationship is shown in Figure 3. Since the TBM gripper shoe can be equivalent to a stationary state along the driving direction, only the horizontal direction of the gripper shoe cylinder is considered in the dynamic model. Its relative relationship with the middle section of the main beam is shown in Figure 3.Where the intermediate variables x_l1cy1 and x_l1cy2 are the relative elastic deformation of the main beam’s front section and the gripper shoe, and the calculation method is shown in equation (7).

$$\{\begin{array}{c}{x}_{l1\mathrm{cy}1}=X_{l1}\cos {\beta }_1-z_{l1}\sin {\beta }_1+a_{l1m}{\theta }_{l1x}\hfill \\ x_{l1\mathrm{cy}2}=X_{l1}\cos {\beta }_1+z_{l1}\sin {\beta }_1+a_{l1m}{\theta }_{l1x}\hfill \end{array}$$

(7)

Where β₁ is the support angle of the main propulsion cylinder, and a_l1m is the projected distance from the center of mass of the front section of the main beam to the support direction of the main propulsion cylinder.

Figure 3.

The relative displacement relationship.

The intermediate variables δ_l2c1 and δ_l2c2 are the relative displacements of the main beam’s middle section and the saddle, respectively, as shown in formula (8).

$$\{\begin{array}{c}{\delta }_{l2c1}=X_{l2}\cos {\beta }_2+Y_{l2}\sin {\beta }_2-X_c\cos {\beta }_2\hfill \\ {\delta }_{l2c2}=-X_{l2}\cos {\beta }_2+Y_{l2}\sin {\beta }_2+X_c\cos {\beta }_2\hfill \end{array}$$

(8)

Where X_l2, Y_l2, X_c, Y_c are the horizontal and vertical displacements of the middle section of the main beam and the saddle, respectively, and β₂ is the support angle of the anti-twist cylinder.

Overturning-axial displacement relationship of main drive system

The relative relationship between the overturning displacement and the axial displacement of the cutter head, support shell, and main beam is shown in Figure 4.

Figure 4.

Axial-overturning relationship of the main system.

The intermediate variable Z_ri (i = 1–4) is the axial displacement of the lower, left, upper, and right positions of the cutter head. Z_di (i = 1–4) is the axial displacement of the lower, left, upper, and right positions of the shield. Z_l1i (i = 1–4) is the axial displacement of the lower, left, upper, and right positions of the front section of the main beam.

$$\{\begin{array}{c}{z}_{r1}=z_r-a_r{\theta }_{\mathrm{rx}}\hfill \\ z_{r2}=z_r+a_r{\theta }_{\mathrm{ry}}\hfill \\ z_{r3}=z_r+a_r{\theta }_{\mathrm{rx}}\hfill \\ z_{r4}=z_r-a_r{\theta }_{\mathrm{ry}}\hfill \end{array}\{\begin{array}{c}{z}_{d1}=z_d-a_{l1}{\theta }_{\mathrm{dx}}\hfill \\ z_{d2}=z_d+a_{\mathrm{lq}}{\theta }_{\mathrm{dy}}\hfill \\ z_{d3}=z_d+a_{l1}{\theta }_{\mathrm{dx}}\hfill \\ z_{d4}=z_d-a_{\mathrm{lq}}{\theta }_{\mathrm{dy}}\hfill \end{array}\{\begin{array}{c}{z}_{l11}=z_{l1}-a_{\mathrm{lq}}{\theta }_{l1x}\hfill \\ z_{l12}=z_{l1}+a_{\mathrm{lq}}{\theta }_{l1y}\hfill \\ z_{l13}=z_{l1}+a_{\mathrm{lq}}{\theta }_{l1x}\hfill \\ z_{l14}=z_{l1}-a_{\mathrm{lq}}{\theta }_{l1y}\hfill \end{array}$$

(9)

Where a_r is the radius where the axial thrust roller of the bearing is located, a_l1 is half the length of the shield body’s flange and the front main beam, and a_lq is half the width of shield body’s flange and the front main beam.

Establishment of system dynamics differential equations

According to the dynamic equivalent model and the relative displacement relationship of each part, the dynamic equations corresponding to the dynamic model of the TBM main drive system are shown in equations (10) to (13).

For cutterhead:

$$\{\begin{array}{c}{m}_r{\stackrel{··}{X}}_r+\sum _{i=1}^N(F_{pri}+D_{pri})\sin ({\phi }_i+\alpha )+C_{\mathrm{eqr}}({\stackrel{·}{X}}_r-{\stackrel{·}{X}}_d)+k_{\mathrm{eqr}}(X_r-X_d)=F_X\hfill \\ m_r{\stackrel{··}{Y}}_r-\sum _{i=1}^N(F_{pri}+D_{pri})\cos ({\phi }_i+\alpha )+C_{\mathrm{eqr}}({\stackrel{·}{Y}}_r-{\stackrel{·}{Y}}_d)+k_{\mathrm{eqr}}(Y_r-Y_d)=F_Y\hfill \\ m_r{\stackrel{··}{z}}_r=F_z-\sum _{i=1}^4[k_{\mathrm{ardi}}(z_{\mathrm{ri}}-z_{\mathrm{di}})+c_{\mathrm{adl}1i}({\stackrel{·}{z}}_{\mathrm{ri}}-{\stackrel{·}{z}}_{\mathrm{di}})]\hfill \\ I_{\mathrm{rx}}{\stackrel{··}{\theta }}_{\mathrm{rx}}=M_x+a_r[k_{\mathrm{ard}1}(z_{r1}-z_{d1})+c_{\mathrm{ard}1}({\stackrel{·}{z}}_{r1}-{\stackrel{·}{z}}_{d1})-k_{\mathrm{ard}3}(z_{r3}-z_{d3})-c_{\mathrm{ard}3}({\stackrel{·}{z}}_{r3}-{\stackrel{·}{z}}_{d3})]\hfill \\ I_{\mathrm{ry}}{\stackrel{··}{\theta }}_{\mathrm{ry}}=M_y+a_r[k_{\mathrm{ard}4}(z_{r4}-z_{d4})+c_{\mathrm{ard}4}({\stackrel{·}{z}}_{r4}-{\stackrel{·}{z}}_{d4})-k_{\mathrm{ard}2}(z_{r2}-z_{d2})-c_{\mathrm{ard}2}({\stackrel{·}{z}}_{r2}-{\stackrel{·}{z}}_{d2})]\hfill \\ I_r{\stackrel{··}{\theta }}_r-\sum _{i=1}^N(F_{\mathrm{pri}}+D_{\mathrm{pri}})r_{\mathrm{br}}=-T_L\hfill \end{array}$$

(10)

For drive gears and motors:

$$\{\begin{array}{c}{m}_{pi}{\stackrel{··}{X}}_{\mathrm{pi}}+(F_{pri}+D_{pri})\cos \alpha +C_{\mathrm{eqpi}}{\stackrel{·}{\delta }}_{\mathrm{dpxi}}+k_{\mathrm{eqp}i}{\delta }_{\mathrm{dpxi}}=0\hfill \\ m_{pi}{\stackrel{··}{Y}}_{\mathrm{pi}}+(F_{pri}+D_{pri})\sin \alpha +C_{\mathrm{eqpi}}{\stackrel{·}{\delta }}_{\mathrm{dpyi}}+k_{\mathrm{eqpi}}{\delta }_{\mathrm{dpyi}}=0\hfill \\ I_{\mathrm{mi}}{\stackrel{··}{\theta }}_{\mathrm{mi}}+C_{\mathrm{mpQ}}({\stackrel{·}{\theta }}_{\mathrm{mi}}-{\stackrel{·}{\theta }}_{\mathrm{pi}})+k_{\mathrm{mpQ}}({\theta }_{\mathrm{mi}}-{\theta }_{\mathrm{pi}})=T_{\mathrm{pi}}\hfill \\ I_{\mathrm{pi}}{\stackrel{··}{\theta }}_{\mathrm{pi}}+(F_{\mathrm{pri}}+D_{\mathrm{pri}})r_{\mathrm{bp}}+C_{\mathrm{mpQ}}({\stackrel{·}{\theta }}_{\mathrm{pi}}-{\stackrel{·}{\theta }}_{\mathrm{mi}})+k_{\mathrm{mpQ}}({\theta }_{\mathrm{pi}}-{\theta }_{\mathrm{mi}})=0\hfill \end{array}$$

(11)

For shield body:

$$\{\begin{array}{c}{m}_d{\stackrel{··}{X}}_d+C_{\mathrm{eqr}}({\stackrel{·}{X}}_d-{\stackrel{·}{X}}_r)+C_{\mathrm{eqd}}{\stackrel{·}{X}}_d+k_{\mathrm{eqr}}(X_d-X_r)+k_{\mathrm{eqd}}{X}_d+\hfill \\ \sum _{i=1}^8(C_{\mathrm{eqp}}{\stackrel{·}{\delta }}_{\mathrm{pdxi}}+k_{\mathrm{eqp}}{\delta }_{\mathrm{pdxi}})+[k_{\mathrm{hdl}1}(X_d-X_{l1})+c_{\mathrm{hdl}1}({\stackrel{·}{X}}_d-{\stackrel{·}{X}}_{l1})]=0\hfill \\ m_d{\stackrel{··}{Y}}_d+C_{\mathrm{eqr}}({\stackrel{·}{Y}}_d-{\stackrel{·}{Y}}_r)+C_{\mathrm{eqd}}{\stackrel{·}{Y}}_d+k_{\mathrm{eqr}}(Y_d-Y_r)+k_{\mathrm{eqd}}{Y}_d+\hfill \\ \sum _{i=1}^8(C_{\mathrm{eqp}}{\stackrel{·}{\delta }}_{\mathrm{pdyi}}+k_{\mathrm{eqp}}{\delta }_{\mathrm{pdyi}})+k_{\mathrm{vdl}1}(Y_d-Y_{l1})+c_{\mathrm{vdl}1}({\stackrel{·}{Y}}_d-{\stackrel{·}{Y}}_{l1})]=0\hfill \\ m_d{\stackrel{··}{z}}_d=\sum _{i=1}^4[k_{\mathrm{ardi}}(z_{\mathrm{ri}}-z_{\mathrm{di}})+c_{\mathrm{adl}1i}({\stackrel{·}{z}}_{\mathrm{ri}}-{\stackrel{·}{z}}_{\mathrm{di}})]-\hfill \\ \sum _{i=1}^4[k_{\mathrm{adl}1i}(z_{\mathrm{di}}-z_{l1i})+c_{\mathrm{adl}1i}({\stackrel{·}{z}}_{\mathrm{di}}-{\stackrel{·}{z}}_{l1i})]\hfill \\ I_{\mathrm{dx}}{\stackrel{··}{\theta }}_{\mathrm{dx}}=-a_r[k_{\mathrm{ard}1}(z_{r1}-z_{d1})+c_{\mathrm{ard}1}({\stackrel{·}{z}}_{r1}-{\stackrel{·}{z}}_{d1})-k_{\mathrm{ard}3}(z_{r3}-z_{d3})\hfill \\ -c_{\mathrm{ard}3}({\stackrel{·}{z}}_{r3}-{\stackrel{·}{z}}_{d3})]+a_{l1}[k_{\mathrm{adl}11}(z_{d1}-z_{l11})+c_{\mathrm{adl}11}({\stackrel{·}{z}}_{d1}-{\stackrel{·}{z}}_{l11})\hfill \\ -k_{\mathrm{adl}13}(z_{d3}-z_{l13})-c_{\mathrm{adl}13}({\stackrel{·}{z}}_{d3}-{\stackrel{·}{z}}_{l13})]-f_l{r}_d-f_r{r}_d\hfill \\ I_{\mathrm{dy}}{\stackrel{··}{\theta }}_{\mathrm{dy}}=-a_r[k_{\mathrm{ard}4}(z_{r4}-z_{d4})+c_{\mathrm{ard}4}({\stackrel{·}{z}}_{r4}-{\stackrel{·}{z}}_{d4})-k_{\mathrm{ard}2}(z_{r2}-z_{d2})\hfill \\ -c_{\mathrm{ard}2}({\stackrel{·}{z}}_{r2}-{\stackrel{·}{z}}_{d2})]+a_{\mathrm{lq}}[k_{\mathrm{adl}14}(z_{d4}-z_{l14})+c_{\mathrm{adl}14}({\stackrel{·}{z}}_{d4}-{\stackrel{·}{z}}_{l14})\hfill \\ -k_{\mathrm{adl}12}(z_{d2}-z_{l12})-c_{\mathrm{adl}12}({\stackrel{·}{z}}_{d2}-{\stackrel{·}{z}}_{l12})]-f_b{r}_d-f_a{r}_d\hfill \end{array}$$

(12)

For main beam:

$$\{\begin{array}{c}{m}_{l1}{\stackrel{··}{X}}_{l1}=[k_{\mathrm{hdl}1}(X_d-X_{l1})+c_{\mathrm{hdl}1}({\stackrel{·}{X}}_d-{\stackrel{·}{X}}_{l1})]-[k_{\mathrm{hl}1l2}(X_{l1}-X_{l2})+c_{\mathrm{hl}1l2}({\stackrel{·}{X}}_{l1}-{\stackrel{·}{X}}_{l2})]\hfill \\ -[X_{l1\mathrm{cy}1}{k}_{\mathrm{cy}1}+{\stackrel{·}{X}}_{l1\mathrm{cy}1}{c}_{\mathrm{cy}1}]\cos {\beta }_1-[X_{l1\mathrm{cy}2}{k}_{\mathrm{cy}2}+{\stackrel{·}{X}}_{l1\mathrm{cy}2}{c}_{\mathrm{cy}2}]\cos {\beta }_1\hfill \\ m_{l1}{\stackrel{··}{Y}}_{l1}=k_{\mathrm{vdl}1}(Y_d-Y_{l1})+c_{\mathrm{vdl}1}({\stackrel{·}{Y}}_d-{\stackrel{·}{Y}}_{l1})]-[k_{\mathrm{vl}1l2}(Y_{l1}-Y_{l2})+c_{\mathrm{vl}1l2}({\stackrel{·}{Y}}_{l1}-{\stackrel{·}{Y}}_{l2})]\hfill \\ m_{l2}{\stackrel{··}{X}}_{l2}=[k_{\mathrm{hl}1l2}(X_{l1}-X_{l2})+c_{\mathrm{hl}1l2}({\stackrel{·}{X}}_{l1}-{\stackrel{·}{X}}_{l2})]-(k_{\mathrm{vcy}1}{\delta }_{l2c1}+c_{\mathrm{vcy}1}{\stackrel{·}{\delta }}_{l2c1})\cos {\beta }_2\hfill \\ +(k_{\mathrm{vcy}2}{\delta }_{l2c2}+c_{\mathrm{vcy}2}{\stackrel{·}{\delta }}_{l2c2})\cos {\beta }_2-[k_{h12c}(X_{l2}-X_c)+c_{h12c}({\stackrel{·}{X}}_{l2}-{\stackrel{·}{X}}_c)]\hfill \\ m_{l2}{\stackrel{··}{Y}}_{l2}=[k_{\mathrm{vl}1l2}(Y_{l1}-Y_{l2})+c_{\mathrm{vl}1l2}({\stackrel{·}{Y}}_{l1}-{\stackrel{·}{Y}}_{l2})]-(k_{\mathrm{vcy}1}{\delta }_{l2c1}+c_{\mathrm{vcy}1}{\stackrel{·}{\delta }}_{l2c1})\sin {\beta }_2\hfill \\ -(k_{\mathrm{vcy}2}{\delta }_{l2c2}+c_{\mathrm{vcy}2}{\stackrel{·}{\delta }}_{l2c2})\sin {\beta }_2\hfill \\ m_c{\stackrel{··}{X}}_c=[k_{h12c}(X_{l2}-X_c)+c_{h12c}({\stackrel{·}{X}}_{l2}-{\stackrel{·}{X}}_c)]-(k_{\mathrm{hcy}1}{X}_c+c_{\mathrm{hcy}1}{\stackrel{·}{X}}_c)-(k_{\mathrm{hcy}2}{X}_c+c_{\mathrm{hcy}2}{\stackrel{·}{X}}_c)\hfill \\ +(k_{\mathrm{vcy}1}{\delta }_{l2c1}+c_{\mathrm{vcy}1}{\stackrel{·}{\delta }}_{l2c1})\cos {\beta }_2-(k_{\mathrm{vcy}2}{\delta }_{l2c2}+c_{\mathrm{vcy}2}{\stackrel{·}{\delta }}_{l2c2})\cos {\beta }_2\hfill \\ m_{l1}{\stackrel{··}{z}}_{l1}=\sum _{i=1}^4[k_{\mathrm{adl}1i}(z_{\mathrm{di}}-z_{l1i})+c_{\mathrm{adl}1i}({\stackrel{·}{z}}_{\mathrm{di}}-{\stackrel{·}{z}}_{l1i})]-\sum _{i=1}^4[k_{\mathrm{al}1l2i}{z}_{l1i}+c_{\mathrm{al}1l2i}{\stackrel{·}{z}}_{l1i}]\hfill \\ +[x_{l1\mathrm{cy}1}{k}_{\mathrm{cy}1}+{\stackrel{·}{x}}_{l1\mathrm{cy}1}{c}_{\mathrm{cy}1}]\sin \alpha -[x_{l1\mathrm{cy}2}{k}_{\mathrm{cy}2}+{\stackrel{·}{x}}_{l1\mathrm{cy}2}{c}_{\mathrm{cy}2}]\sin \alpha \hfill \\ I_{l1x}{\stackrel{··}{\theta }}_{l1x}=-a_{l1}[k_{\mathrm{adl}11}(z_{d1}-z_{l11})+c_{\mathrm{adl}11}({\stackrel{·}{z}}_{d1}-{\stackrel{·}{z}}_{l11})-k_{\mathrm{adl}13}(z_{d3}-z_{l13})-c_{\mathrm{adl}13}({\stackrel{·}{z}}_{d3}-{\stackrel{·}{z}}_{l13})]\hfill \\ +a_{l1}[k_{\mathrm{al}1l21}{z}_{l11}+c_{\mathrm{al}1l21}{\stackrel{·}{z}}_{l11}-k_{\mathrm{al}1l23}{z}_{l13}-c_{\mathrm{al}1l23}{\stackrel{·}{z}}_{l13}]\hfill \\ +a_{l1m}[-x_{l1\mathrm{cy}1}{k}_{\mathrm{cy}1}-{\stackrel{·}{x}}_{l1\mathrm{cy}1}{c}_{\mathrm{cy}1}-x_{l1\mathrm{cy}2}{k}_{\mathrm{cy}2}-{\stackrel{·}{x}}_{l1\mathrm{cy}2}{c}_{\mathrm{cy}2}]\hfill \\ I_{l1y}{\stackrel{··}{\theta }}_{l1y}=-a_{\mathrm{lq}}[k_{\mathrm{adl}14}(z_{d4}-z_{l14})+c_{\mathrm{adl}14}({\stackrel{·}{z}}_{d4}-{\stackrel{·}{z}}_{l14})-k_{\mathrm{adl}12}(z_{d2}-z_{l12})-c_{\mathrm{adl}12}({\stackrel{·}{z}}_{d2}-{\stackrel{·}{z}}_{l12})]\hfill \\ +a_{l\mathrm{q}}[k_{\mathrm{al}1l24}{z}_{l14}+c_{\mathrm{al}1l24}{\stackrel{·}{z}}_{l14}-k_{\mathrm{al}1l22}{z}_{l12}-c_{\mathrm{al}1l22}{\stackrel{·}{z}}_{l12}]\hfill \end{array}$$

(13)

Multi-point vibration measurement of TBM main drive system

As the vibration of the cutter head part is the strongest. The lower part of the cutter head often has the presence of water. The installation position of the sensor will obvious influence its normal operation. The location of the sensor on the TBM’s cutterhead is shown in Figure 5(a) and (b). The acceleration sensor in the shield part is placed behind the bearing ring as shown in Figure 5(c), and the sensor in the front part of the main beam is placed in the place where the main beam and the support shield are connected. The environmental rock parameters, excavation parameters, and main support state parameters of the system during the vibration measurement process as shown in Table 2.

Figure 5.

Layout of cutterhead acceleration sensor and test diagram: (a) cutterhead T1, (b) support shell T2, (c) main beam T3, and (d) tunneling parameters.

Table 2.

Main driving parameters in vibration measurement.

Rock type	Class II	Rock compressive strength	90–100 MPa
Cutterhead rotation speed	6 rpm	Penetration rate	5.8 mm/r
Total motor torque	2925 kNm	Thrust cylinder pressure Rodless chamber/rodded Chamber	182/2 bar
Excavation inclination	0.6%	Roof shield cylinder pressure	82 bar
Right wedge block cylinder Pressure	82 bar	Left wedge block cylinder Pressure	92 bar
Right shield support cylinder pressure Rodless Chamber/Rodded Chamber	54/104 bar	Left shield support cylinder pressure Rodless Chamber/Rodded Chamber	64/88 bar
Right support shoe cylinder Pressure	4 bar	Left support shoe cylinder Pressure	6 bar
Right torque cylinder pressure Rodless chamber/rodded chamber	78/92 bar	Left torque cylinder pressure Rodless chamber/rodded Chamber	102/68 bar
Right rear support cylinder Pressure	56 bar	Left rear support cylinder Pressure	56 bar

Results

Equivalent load of cutterhead

The load of the cutter at each position of the cutterhead is synthesized to obtain the equivalent load of the cutterhead as shown in Figure 6.

Figure 6.

Equivalent load of cutterhead.

From Figure 5, it can be seen that the dynamic loads of the cutterhead exhibit significant abrupt changes in all directions, among which the transverse and longitudinal overturning moments and the transverse and longitudinal unbalanced loads also show obvious periodic variation characteristics, with the variation period consistent with the rotation period of the cutterhead. After statistical analysis, the RMS values of axial load and torque were 5610.9 kN and 2872.2 kNm, respectively, with maximum amplitudes of 1448.1 kN and 190.2 kNm, respectively; The RMS values of its horizontal and vertical unbalanced loads are 192.1 and 164.2 kN, respectively, with maximum amplitudes of 637.7 and 480.1 kN, respectively; The RMS values of its transverse and longitudinal overturning moments are 1723.0 and 1763.0 kN, respectively, with maximum amplitudes of 3250.9 and 3388.3 kNm, respectively.

Vibration characteristics of main drive system

In the dynamic model of the main drive system, the structural stiffness of the front and middle sections of the main beam is calculated by the finite element method. The dynamic parameters to be modified mainly include the structural stiffness of the front and middle sections of the main beam. The estimated range of each stiffness is shown in equation:

$$\begin{array}{c}5\times {10}^8\le k_{\mathrm{adl}1i}\le 1\times {10}^9\hfill \\ 4\times {10}^8\le k_{\mathrm{al}1l2i}\le 1\times {10}^9\hfill \\ 1\times {10}^9\le k_{\mathrm{hdl}1}=k_{\mathrm{vdl}1}\le 2\times {10}^9\hfill \\ 1\times {10}^9\le k_{\mathrm{hl}1l2}=k_{\mathrm{vl}1l2}\le 2\times {10}^9\hfill \end{array}$$

(14)

In order to improve the calculation accuracy of the coupled model for the system vibration response, the main structural stiffness of the TBM main drive was modified based on the excavation data of the field test. The optimization process takes the minimum error between the vibration response RMS value of the dynamic model and the test data as the optimization goal. The optimization method adopts the Hooke-Jeeves method (K and M, 1997; Table 3).

Table 3.

Modified TBM structure stiffness parameters.

Axial stiffness of front section of main beam N/m	Axial stiffness of middle section of main beam N/m	Radial stiffness of front section of main beam N/m	Radial stiffness of the middle section of the main beam N/m
2.34e8	2.08e8	1.25e9	1.05e9

After the model is modified, the three-directional vibration of the cutter head, support shield, and main beam of the system is calculated based on the coupled dynamics model, and the comparison with the measured vibration data of each component is shown in Figure 7 and Table 4.

Figure 7.

Comparison between vibration test data and calculation value. (a) Axial acceleration of the cutterhead. (b) Horizontal acceleration of the cutterhead. (c) Vertical acceleration of the cutterhead. (d) Axial acceleration of support shell. (e) Horizontal acceleration of support shell. (f) Vertical acceleration of support shell. (g) Axial acceleration of main beam. (h) Horizontal acceleration of main beam. (i) Vertical acceleration of main beam.

Table 4.

Test data and calculated value statistics.

Measured position	Axial direction			Vertical direction			Horizontal direction
	Calculated	Measured	Error%	Calculated	Measured	Error%	Calculated	Measured	Error%
Cutterhead	0.733	0.827	11.30	0.510	0.491	3.85	0.527	0.642	17.81
Support shell	0.401	0.331	20.99	0.188	0.221	14.64	0.276	0.272	1.47
Main beam	0.241	0.246	2.11	0.181	0.148	22.00	0.154	0.129	19.60

As shown in Figure 7, the maximum errors of the RMS values of the three-directional vibration at the specific position of the TBM cutter head, the support shield, and the front section of the main beam are within 18%, 21% and 22%, respectively.

The frequency response data corresponding to the main vibration directions of the cutter head, the interface between the support shell and the main beam are extracted. It was compared with the frequency domain response calculated by the dynamic model, as shown in Figure 8.

Figure 8.

Frequency domain of vibration test data and calculation value: (a) axial acceleration of the cutterhead and (b) axial acceleration of support shell.

It can be seen that the axial acceleration frequency domain response of the cutter head is a typical broadband frequency domain vibration, which is caused by the uncertain load of the surrounding rock in the front section during the field excavation process. The vibration distribution of the measured frequency domain response is relatively uniform in the frequency domain of 0–80 Hz, and shows a slightly larger proportion in the frequency domain of 20–30 Hz and 50–60 Hz. But the calculated frequency domain response of the dynamic model shows a slightly larger proportion of vibration at frequencies of 10–15 Hz, 25, 42, 50, and 62 Hz. This frequency domain response is very close to the frequency domain response of the axial equivalent load of the cutterhead, ¹⁹ so the vibration of the cutter head is mainly reflected as the forced vibration of the system under the external excitation load.

The vibration test part of the support shield is located at the connecting flange between the support shield and the front section of the main beam. The frequency-domain response of the axial vibration acceleration is obviously concentrated around 18–20 Hz, and it also shows larger vibration at 25–30 Hz. The main frequency band of the frequency domain response of vibration acceleration at this position calculated by the dynamic model is concentrated around 17.4 Hz. From the analysis of the natural frequency of the system in section “Inherent characteristics of main drive system and multi-parameter influence analysis” of this paper, it can be seen that it is the first order natural frequency of the system, and its corresponding mode shape is horizontal overturning vibration. The main frequency band of the system is relatively consistent with the measured frequency domain response.

Inherent characteristics of main drive system and multi-parameter influence analysis

Analysis of inherent characteristics of main drive system

The first 10 natural frequencies of the TBM main drive system and the corresponding modes are shown in Figure 9. It can be seen that the first-order mode of the system is 17.44 Hz, which is dominated by the horizontal overturning vibration of the system, accompanied by the horizontal translation vibration of the system. The second-order and third-order modes are 19.23 and 21.39 Hz, respectively, and are dominated by overturning vibration. The first three modes of the system are similar, and the coupled vibration of the first to third modes is easy to occur in the actual vibration. The mode is the coupled vibration of the overturning vibration in each direction and the horizontal translation vibration. The natural frequency corresponding to the axial mode is the highest, followed by the vertical translation mode, and the natural frequency corresponding to the horizontal mode is the smallest. Therefore, the rigidity of the system is the highest in the axial direction, followed by the vertical direction, and the smallest in the horizontal direction.

Figure 9

Mode corresponding to natural frequency of the system.

Multi-parameter analysis of natural frequency

Equivalent support stiffness of main propulsion cylinder

The main propulsion cylinder plays a major role in supporting the axial direction of the system and provides the power required for the axial movement of the system. The relationship between the equivalent support stiffness of the cylinder and the natural frequency of the system is shown in Figure 10(a).

Figure 10.

Relationship between propulsion cylinder stiffness and natural frequency.

In the figure, X_i, Y_i, and Z_i represent the natural frequencies corresponding to the i-th horizontal, vertical, and axial translational vibration modes of the system, respectively. θ_xi and θ_yi represent the natural frequencies corresponding to the i-th horizontal and vertical overturning vibration modes of the system, respectively.

It can be seen that with the increase of the equivalent support stiffness of the main propulsion cylinder, the corresponding natural frequencies of the first-order and second-order horizontal overturning modes of the system show a significant increasing trend. When the equivalent support stiffness is about 1e9 N/m, the phenomenon of modal superposition occurs, and the coupled vibration of horizontal and vertical overturning vibration occurs. Within the overall range of equivalent stiffness changes, the natural frequency of the first-order axial vibration of the system shows a relatively obvious increasing trend. At the same time, the natural frequency corresponding to the first-order vertical translation vibration of the system does not change significantly.

Equivalent support stiffness of shield body

It can be seen from Figure11(b) that the equivalent support stiffness of the shield body has a significant influence on the second-order horizontal overturning mode and the first-order horizontal and vertical translation mode of the system. When the equivalent support stiffness of the shield is greater than 8e8 N/m, the second-order horizontal overturning mode increases significantly.

Figure 11.

Relationship between main bearing support stiffness and natural frequency: (a) main bearing axial stiffness and (b) main bearing radial stiffness.

Equivalent support stiffness of main bearing

The main bearing of TBM is the main bearing structure of the system. The relationship between the equivalent axial and radial equivalent stiffness of the main bearing and the system mode is shown in the Figure 11.

It can be seen that the natural frequencies of the first-order overturning vibration modes show a significant increasing trend. In the range of 1e8–3e8 N/m, the two modes show the phenomenon of mode separation, and in the range of 3e8–1e9 N/m, the two modes show the trend of mode overlap. The equivalent axial stiffness of the main bearing has a great influence on the various order modes of the system, especially on the axial translation vibration and overturning vibration modes of the system.

The modal energy analysis of the system

On the basis of system modal analysis, we have analyzed the modal energy of the system. The specific analysis is as follows. The kth mode kinetic energy W_i and mode deformation energy E_i of the main drive system are expressed as follows:

$$E_i=\frac{1}{2}{{\left\{\mathit{\phi }\right\}}_k}^T\overline{\mathit{K}}{\left\{\mathit{\phi }\right\}}_k$$

(15)

$$W_i=\frac{1}{2}{\omega }_k^2{{\left\{\mathit{\phi }\right\}}_k}^T\mathit{M}{\left\{\mathit{\phi }\right\}}_k$$

(16)

Where ω_k is the kth natural frequency, K and M are the stiffness matrix and mass matrix of the system, respectively, and {φ} _k is the principal mode vector corresponding to the kth natural frequency.

The modal deformation energy of each component is shown in Figure 12. From Figure 12, it can be seen that the total deformation energy of the TBM main drive system is mainly concentrated in the first 10 modes of the system. When the modes are greater than 10, the deformation energy of the system is very low, indicating that the deformation degree of the system is relatively large during low-frequency vibration. In the first 10 natural frequencies, the natural frequencies corresponding to the transverse and longitudinal overturning vibration of the first and second orders also exhibit higher total system deformation energy. The modal deformation energy corresponding to the lateral overturning vibration mode is the maximum supported by the shell. At the same time, the natural frequencies corresponding to axial vibration and longitudinal overturning vibration in the 7th to 10th order also exhibit high modal deformation energy. The deformation energy of the cutterhead is the highest, followed by the supporting shell, and the main beam is the smallest.

Figure 12.

Modal deformation energy distribution. (a) The modal deformation energy of cutterhead. (b)The modal deformation energy of shield. (c) The modal deformation energy of main beam. (d) Proportion of modal deformation energy.

The modal kinetic energy of each component is shown in Figure 13. From Figure 13, it can be concluded that the total modal kinetic energy of the system is mainly concentrated in the first 10 modes of vibration of the TBM main drive system. The total kinetic energy is mainly concentrated in the transverse and longitudinal overturning modes of order 1–2 and the axial translation and longitudinal overturning modes of order 8–10. The total kinetic energy of the torsional vibration modes of the third and fifth order systems along the excavation direction is not significant. In its first-order lateral overturning mode, the kinetic energy of the supporting shell accounts for the largest proportion. In the first order axial translational vibration mode, the kinetic energy of the cutterhead accounts for the largest proportion.

Figure 13.

System modal kinetic energy distribution. (a) The kinetic energy energy of cutterhead. (b) The kinetic energy energy of shield. (c) The kinetic energy energy of main beam. (d) Proportion of modal kinetic energy.

Sensitivity analysis of system dynamic parameters

Sensitivity analysis of cutterhead mass

This paper analyzes the sensitivity of the mass of the cutterhead to system vibration. The mass range varies from 80% of the original mass to 120%, and the analysis results of the three-dimensional vibration amplitude and RMS value of each component are shown in the following Figure 14.

Figure 14.

Relationship between cutterhead quality and cutterhead vibration. (a) The relationship between cutterhead mass and vibration RMS value. (b) The relationship between cutterhead mass and vibration amplitude.

It can be seen that during the process of the mass of the cutterhead changing from 80% of the original mass to 120%, the RMS value of the cutterhead’s translational vibration has undergone significant changes. As the mass coefficient of the cutterhead increases, the RMS value of its translational vibration shows a significant downward trend in the mass coefficient, and the rate of decrease is similar. During the process of increasing the mass coefficient of the cutterhead from 0.8 to 1.2, the RMS value of its axial vibration decreased by 29.5%, and the RMS values of its transverse and longitudinal vibrations decreased by 36.9% and 36.2%, respectively. The amplitude of its isotropic translational vibration also shows a significant decreasing trend, with the axial vibration amplitude of the cutterhead decreasing at a similar rate within the overall range of mass coefficient variation. When the mass coefficient reaches 1.2, the vibration amplitude decreases by 33.0%. The lateral vibration amplitude shows a continuous decreasing trend during the process of changing the mass coefficient, and the decreasing speed gradually decreases, with a decrease amplitude of 31.7%. The longitudinal vibration amplitude decreases rapidly within the range of 0.85–0.95 and 1.15–1.2 of the mass coefficient variation, and the longitudinal vibration amplitude decreases by 45.5% within the overall range of mass coefficient variation.

Sensitivity analysis of hydraulic cylinder stiffness

The main propulsion cylinder plays a major axial support role in the main drive system. In order to study the influence of the axial stiffness of the main propulsion cylinder on the dynamic characteristics of the system, the relationship between its stiffness and the overturning and translational vibration of the main beam was extracted as shown in Figure 15.

Figure 15.

Relationship between the stiffness of the propulsion cylinder and the RMS value of the main beam vibration. (a) The relationship between propulsion cylinder stiffness and translational vibration RMS value. (b) The relationship between propulsion cylinder stiffness and overturning vibration RMS value. (c)The relationship between propulsion cylinder stiffness and translational vibration amplitude. (d)The relationship between propulsion cylinder stiffness and overturning vibration amplitude.

It can be seen that the angle of the overturning vibration of the main beam and the stiffness of the propulsion cylinder have a significant impact on the lateral overturning vibration of the main beam, while they have almost no effect on the longitudinal overturning vibration of the main beam. When the stiffness of the propulsion cylinder is less than 4e8 N/m, the RMS value and amplitude of the transverse overturning vibration in the front section of the main beam do not change significantly, at the levels of 0.6–0.7 rad/s ² and 1.8–1.9 rad/s ² , respectively. When the stiffness of the propulsion cylinder is greater than 4e8 N/m, the RMS value and amplitude of the transverse overturning vibration in the front section of the main beam show a significant downward trend. When the stiffness of the propulsion cylinder increases to 1e9 N/m, the RMS value of its lateral overturning vibration decreases by 59.93% and the amplitude decreases by 57.66%. From the perspective of translational vibration of the main beam, the stiffness of the main propulsion cylinder has the most significant impact on the axial vibration of the main beam, followed by its lateral vibration, and has little effect on its longitudinal vibration. The axial vibration amplitude and RMS value significantly decrease with the increase of its equivalent support stiffness. During the process of increasing its support stiffness from 1e8 to 10E8, the axial vibration amplitude of the main beam decreased by 33.44%, and its axial vibration RMS decreased by 43.60%. The RMS value of lateral vibration shows a slight decrease with the increase of the stiffness of the supporting oil cylinder, which is 2.62% lower than that when the stiffness is 1e8 N/m. The amplitude and RMS value of its longitudinal vibration are not significantly related to the equivalent support stiffness of the main propulsion cylinder. The lateral vibration amplitude shows a significant change with the increase of equivalent support stiffness. When the lateral vibration amplitude of the main beam is below 5e8 N/m, it shows a negative correlation with the equivalent support stiffness of the main propulsion cylinder. When the lateral vibration amplitude is above 5e8 N/m, it shows a positive correlation with the equivalent support stiffness of the main propulsion cylinder.

Sensitivity analysis of radial stiffness of main bearing

In order to study the influence of the radial stiffness of the main bearing on its dynamic characteristics, the relationship between the radial stiffness of the main bearing and the radial translational vibration of the cutterhead and support housing was extracted as shown in Figure 16.

Figure 16.

Relationship between main bearing radial stiffness and system acceleration. (a) The relationship between main bearing stiffness and vibration RMS value of cutterhead. (b) The relationship between main bearing stiffness and vibration RMS value of shield. (c)The relationship between main bearing stiffness and vibration amplitude of cutterhead. (d)The relationship between main bearing stiffness and vibration amplitude of shield.

It can be seen that the RMS value of the transverse and longitudinal translational vibration of the cutterhead fluctuates slightly with the increase of the radial stiffness of the main bearing. The RMS value of its lateral translational vibration is in the range of 5.03–5.45 m/s ² , and the RMS value of its lateral translational vibration is in the range of 4.89–5.30 m/s ² . The lateral translational vibration amplitude of the cutterhead significantly decreases within the range of 5–7e8 N/m of the main bearing radial stiffness, with a reduction of 25.00%. The amplitude of its longitudinal translational vibration shows a slight overall increase trend during the increase of the radial stiffness of the main bearing, with an increase of 13.66%. The RMS value and amplitude of the lateral and longitudinal translational vibration of the supporting shell significantly increase with the increase of the radial stiffness of the main bearing. The RMS values of its lateral and longitudinal translational vibrations increased by 1 and 0.66 times, respectively, while the amplitude of its lateral and longitudinal translational vibrations increased by 1.72 and 0.56 times, respectively.

Vibration transmission law of TBM main drive system

Figure 17 shows the vibration data of acceleration at the position of the cutter head, shield body, main beam and center mass of the gripper shoe. The variation law of the RMS value and amplitude of the vibration acceleration of the system along the main drive driving direction is shown in the Figure 14 and Table 5. It can be seen that the axial translation vibration of the system presents a continuous attenuation trend along the driving direction. The RMS value of the axial translation vibration at the mass center of the shield body is reduced to 59.01% of the cutterhead, and the maximum amplitude is reduced to 55.97% of the cutterhead. The RMS value is reduced to 0.2411 g at the front section of the main beam, which is 35.51% of the cutter head, and the maximum amplitude is reduced to 0.8217 g, which is 44.79% of the cutter head. The horizontal and vertical translation vibration of the system is significantly reduced when the vibration is transmitted to the shield body, the RMS value is reduced to 25.94% and 27.04% of the cutterhead, and the amplitude is reduced to 21.23% and 18.50% of the cutter head. In general, the horizontal and vertical translation vibration first showed a significant decrease to a lower level at the shield body, then increased slightly at the front of the main beam, and finally gradually decreased along the back of the main beam. This shows that the shield body plays a major radial support role for the TBM main drive system.

Figure 17.

Vibration transmission law of TBM main drive system. (a) Translation vibration RMS value transmission law of TBM main drive system. (b) Overturning vibration RMS value transmission law of TBM main drive system. (c) Translation vibration amplitude transmission law of TBM main drive system. (d) Overturning vibration amplitude transmission law of TBM main drive system.

Table 5.

Acceleration response statistics.

Position	Vibration direction	Cutterhead	Support shell	Main beam	Gripper shoe cylinder
RMS (g)	Axial translation	0.6788	0.4006	0.2411	–
	Horizontal translation	0.5274	0.1368	0.1544	0.028
	Vertical translation	0.5103	0.1380	0.1808	–
	Horizontal overturning	0.6703	1.0404	0.6703	–
	Vertical overturning	0.4806	0.6253	0.4680	–
Max amplitude (g)	Axial translation	1.8347	1.0269	0.8217	–
	Horizontal translation	1.9672	0.4176	0.4208	0.070
	Vertical translation	2.2719	0.4204	0.6100	–
	Horizontal overturning	2.1792	2.6936	1.9116	–
	Vertical overturning	1.6295	1.9390	1.4610	–

The horizontal and vertical overturning vibration of the system shows an obvious increasing trend at the shield body. Its RMS value reaches 155.21% and 130.11% of the cutter head, and its amplitude reaches 123.60% and 118.99% of the cutter head. This trend may be due to the small size of the connecting flange between the main beam and the shield body, which is prone to overturning vibration. Figure 18(a)–(c) show the frequency domain diagrams of the displacement vibration of the shield in different directions.

Figure 18.

Frequency response of support shell displacement vibration: (a) axial vibration, (b) radial vibration, and (c) overturning vibration.

It can be seen from Figure 18 that the main vibration form of the shield body is the coupled vibration around 15–20 Hz. Radial translation vibration and overturning vibration also show significant peaks at low frequencies (below 1 Hz), which are mainly due to the periodic change with a frequency of 0.093 Hz under the rotation period of the cutter head.

In order to better compare the vibration energy ratio of the system, after normalizing the vibration energy, the statistical value of the vibration energy ratio of each part is shown in Table 6.

Table 6.

Statistical value of vibration energy proportion.

Components	Horizontal translation (%)	Vertical translation (%)	Axial translation (%)	Horizontal overturning (%)	Vertical overturning (%)
Cutterhead	16.92	17.53	20.47	23.38	21.71
Support shell	7.83	9.44	34.17	30.78	17.78
Main beam	23.80	19.08	21.36	17.18	18.60

It can be seen from Figure 19 that at the position of the cutter head, the radial vibration, axial vibration, and overturning vibration energy of the system are close. At the shield position, the proportion of radial vibration is significantly reduced, and the proportion of vibration energy is reduced to 7.83% and 9.44%, respectively, and the radial vibration energy only accounts for less than 20% of the total system energy. The proportion of overturning vibration energy increased significantly, and the proportion of its overturning vibration energy reached 30.78% and 17.78%, respectively. When in the position of the main beam, the vibration energy in all directions is in an average condition again. Its horizontal translation vibration energy is the largest, reaching 23.8%. According to the statistical value of the vibration energy ratio in each direction in Table 6, the variation law of the vibration energy ratio is shown in Figure 16.

Figure 19.

Law of vibration energy proportion change.

Conclusion

Based on the distributed support of the main drive system and the structure of multi-source input, this paper analyzed the coupling relationship between the degrees of freedom of the system. Considering the nonlinear internal excitation of bearing dynamic stiffness, gear meshing error, tooth side clearance, this paper established a multi DOF coupling dynamic model of the main drive system, including the cutter head, transmission system, support, and propulsion system. The vibration of the main part of the TBM is measured. Based on the measured data, the dynamic model of the main drive system is corrected, and the maximum errors in the calculation of the three-dimensional vibration of the TBM cutterhead, the shield body, and the front section of the main beam are within 18%, 21%, and 22%, respectively.

The analysis shows that the axial acceleration frequency domain response of the cutter head is the vibration of a wide frequency domain bandwidth. At the same time, it shows a larger proportion in the frequency domain of 20–30 Hz and 50–60 Hz. The vibration of the cutter head is mainly reflected in the forced vibration of the system under the external excitation load, which is consistent with the characteristic that the cutter head is in direct contact with the surrounding rock at the front end of the system. For the connection flange between the shield body and the front section of the main beam, the frequency domain response of the axial vibration acceleration is obviously concentrated at 18–20 Hz. The main frequency band is relatively consistent with the measured frequency domain response, which is mainly caused by the excitation of the first three-order natural frequencies of the system under the action of the system’s broadband external excitation load.

The coupled vibration of radial overturning vibration and lateral translation vibration is easy to appear in the actual vibration of the system. The influence of the main mass parameters and stiffness parameters of the system on the various modes of the system is analyzed, and the modal transition points under each dynamic parameter are obtained. By studying the nonlinear variation law of the main components of the system along the main drive tunneling direction, it is found that the overturning vibration of the shield body has a significant “amplification effect,” and its amplitude reaches 123.60% and 118.99% of the cutter head, respectively.

In terms of dynamic modeling research, the coupling relationship between hydraulic and mechanical systems should be considered in the future, which may have an influence on the vibration and load characteristics of the main load-bearing area of the system. At the same time, the clearance and assembly parameters between key components of the TBM main drive system also have an influence on its dynamic behavior. In subsequent research, the dynamic model should be further improved to address this issue. In terms of validating the dynamic model, future research can combine vibration big data to calibrate the dynamic response of the dynamic model under different complex rock conditions, thereby improving the applicability and universality of the dynamic model.