A shield tunneling pressure control method and device

Abstract

The present invention belongs to the technical field related to shield construction safety control, and discloses a shield propulsion pressure control method and equipment, comprising the following steps: (1) constructing a three-dimensional rotating body failure structure model of the composite stratum excavation surface based on the influence of the uncertainty of the soil parameters of the composite stratum on the stability of the shield excavation surface; (2) constructing a reliability analysis index system for the stability of the excavation surface, and calculating the value of the composite stratum ultimate support pressure under the soil parameters under the uncertainty condition under the reliability index; (3) constructing a load calculation module of the shield propulsion system and establishing a mathematical model of the shield hydraulic propulsion system, and then constructing a control model of the shield propulsion system, wherein the control model adopts a fuzzy PID control algorithm for propulsion force control; wherein the control model uses the input voltage of the proportional overflow valve of the shield hydraulic propulsion system as input and uses the propulsion pressure of the shield propulsion system as output. The present invention can accurately control the propulsion pressure.

Description

Technical Field

[0001] The present invention belongs to the technical field related to shield construction safety control, and more specifically, relates to a shield propulsion pressure control method and equipment. Background Art

[0002] Accidents such as surface collapse during shield tunnel construction are primarily caused by instability in the shield tunnel excavation face. In the case of ultra-large diameter composite strata, the risk of shield excavation face instability is even greater due to the diversity of strata and uncertainty in the soil. Whether using an earth pressure balance shield or a slurry balance shield, during shield tunnel construction, the excavation face is stabilized by applying a support stress equal to the original stratum stress to the soil layer ahead of the excavation face to prevent ground deformation caused by excavation. Too little support pressure on the excavation face can lead to soil collapse, while excessive support pressure can cause soil uplift. Therefore, determining and controlling the support pressure during shield construction to control excavation face stability is particularly important.

[0003] Therefore, the research on the stability and ultimate support pressure of shield excavation face in composite strata and the control and adjustment of the ultimate support pressure applied by the shield machine on the soil layer in front based on the soil pressure in front of the excavation face have important academic value and engineering application value. Summary of the Invention

[0004] In response to the above-mentioned defects or improvement needs of the prior art, the present invention provides a shield thrust pressure control method and equipment. The control method combines the calculation process of the ultimate support pressure solution of the spiral failure model of the composite stratum excavation face obtained under the condition of soil integrity as an external load module with the shield thrust control simulation model. By collecting, analyzing and calculating the soil parameters of the shield excavation face under graded reliability indicators, a suitable thrust force is obtained to balance the external soil pressure, ensuring that the excavation process is safe, efficient and effectively controlled, forming a complete set of shield excavation face stability control processes.

[0005] To achieve the above object, according to one aspect of the present invention, a shield propulsion pressure control method is provided, the method comprising the following steps:

[0006] (1) Based on the influence of uncertainty in soil parameters of composite strata on the stability of shield excavation face, a three-dimensional rotating body failure structure model of composite stratum excavation face is constructed. The ultimate support force calculation formula when the soil structure is destroyed under the critical failure state of the excavation face is derived by equating the work done by the external force of the system with the energy dissipated by friction in the soil.

[0007] (2) Construct a reliability analysis index system for excavation face stability and apply the Monte Carlo method to perform simulation calculations to obtain the values ​​of the ultimate support pressure of composite strata under uncertainty conditions of soil parameters under different reliability indicators;

[0008] (3) Construct a load calculation module for the shield propulsion system and establish a mathematical model of the shield hydraulic propulsion system, and then construct a control model for the shield propulsion system. The control model uses a fuzzy PID control algorithm to control the propulsion force; among them, the control model takes the input voltage of the proportional overflow valve of the shield hydraulic propulsion system as input and the propulsion pressure of the shield propulsion system as output.

[0009] Furthermore, the spatial discretization technology is applied to generate the three-dimensional failure surface of the failure model point by point to obtain the coordinates of each point in the excavation face instability velocity field. After constructing the instability model of the excavation face in composite strata, the upper limit theorem of limit analysis is used to calculate the ultimate support pressure based on the power and energy consumption of the excavation working face.

[0010] Furthermore, the extension lines of OA and OE and the stratum boundary can divide the failure structure of the excavation surface of the composite stratum into four regions: I, II, III and IV; the instability velocity field of the excavation surface is completely determined by Z O and Y O Take the value and set Z O and Y O The two parameters are dimensionless through the polar coordinates of the tunnel center point G, with r G / D and θ G The two parameters represent the coordinates of point O in the rectangular coordinate system of the (Y, Z) plane as (Y O ,Z O ), the corresponding formula is:

[0011]

[0012]

[0013] The spatial discretization technology is applied to generate the three-dimensional failure surface of the failure model point by point to obtain the coordinates of each point in the excavation surface instability velocity field; among them, point A is the top of the tunnel, point B is the bottom of the tunnel, point A is the origin of the rectangular coordinate system of the (Y, Z) plane, and the excavation surface instability velocity field is composed of three logarithmic spiral curves AF, EF and BE, among which point E is the intersection of the logarithmic spiral curve and the stratum boundary line, and point F is the intersection of the two logarithmic spiral curves in the upper soil. The entire instability velocity field rotates with point O as the origin, point O is the origin of the polar coordinate system (r, θ), the polar axis is perpendicular to the surface plane and downward, and OA, OB, and OE are the lines connecting each point respectively.

[0014] Furthermore, the system power and energy consumption are related to the soil deadweight σ γ , support pressure σ T And the velocity discontinuity friction force f dissipation is calculated as:

[0015]

[0016]

[0017]

[0018] W f =W γ +W T

[0019] Among them, r A 、r B and r E are the lengths of OA, OB and OE respectively, θ A ,θ B and θ E are the angles between OA, OB and OE and the polar coordinate system (r,θ) relative to the polar axis direction, and are the cohesion, internal friction angle and soil weight of the upper and lower soil layers respectively, ω is the angular velocity, W f is the friction dissipation power, W γ is the soil gravity power, W T is the work power of support pressure.

[0020] Furthermore, the failure probability of the system is used as the basis for judging the risk state of the excavation face stability to confirm the risk state level of the excavation face stability structure. Each stratum change point is simulated using the Monte Carlo simulation method to simulate 1000 sets of corresponding soil parameter vectors based on the corresponding physical and mechanical properties of the soil layer and the characteristics of their respective uncertainties. Each set of soil parameter vectors can correspond to a composite stratum excavation face instability model, and the failure probability of the system is determined to obtain the failure curve of the model excavation face, and then the ultimate support pressure under soil uncertainty is calculated.

[0021] Furthermore, the formulas used in calculating the ultimate support pressure include:

[0022] p f =P(g(σ T )<0)

[0023]

[0024] β=Φ -1 (1-p f )

[0025] Where p f is the failure probability of the system, σ T is the simulated value of the ultimate support pressure of the shield face, represents the limit value of the ultimate support pressure of the shield surface, β is the reliability index, Φ -1 (·) is the inverse function of the standard normal distribution function.

[0026] Furthermore, the full dynamic range that needs to be controlled is obtained through the ultimate support pressure of the same level, that is, an algorithm is constructed to control the thrust value output by the shield thrust system within the full dynamic range.

[0027] Furthermore, based on the mechanical balance equation of the output force and load force of the propulsion hydraulic cylinder of the shield hydraulic propulsion system and the flow continuity equation of the propulsion hydraulic cylinder, the output hydraulic cylinder rodless chamber pressure P1 with respect to the input flow Q is obtained. L and external load force F L The transfer functions are:

[0028]

[0029]

[0030] Where ω1 and ζ1 are the comprehensive natural frequency and comprehensive damping ratio respectively, and their values ​​are:

[0031]

[0032]

[0033] Where A1 is the area of ​​the rodless piston, p1 is the pressure of the rodless chamber, m1 is the total reduced mass of the system, and F L is the external load force acting on the piston, β1 is the effective bulk elastic modulus, λ1 is the leakage coefficient of the hydraulic cylinder, and V1 is the volume of the oil inlet chamber of the hydraulic cylinder.

[0034] Furthermore, the transfer functions of the current dynamic equation, flow equation, and valve core motion equation of the proportional relief valve are constructed as follows:

[0035]

[0036] U(s)=LsI(s)+(r c +r p )I(s)+K v sX1(s)

[0037]

[0038] Assume that the stable working point of the proportional relief valve is q0(x0,p L0 ), then linearization at q0 yields:

[0039] ΔQ=K q Δx1-K p Δp L

[0040] Q(s)=K q X1(s)

[0041]

[0042] A m P1(s)-F e (s)=(m2s 2 +Ds+K)X1(s)

[0043] F e (s) = K I I(s)

[0044] Where u is the valve solenoid coil input voltage, L is the coil inductance, r c and r p are the internal resistances of the coil and amplifier, respectively, K v is the back electromotive force coefficient of the coil, x1 is the valve core displacement, Q is the proportional relief valve flow, α d is the flow coefficient of the cone valve port, d1 is the aperture of the cone valve port, β is the semi-cone angle of the cone valve, ρ is the density of the hydraulic oil in the system, and p L is the equivalent load pressure, and Because K p is very small, so it can be ignored. m is the valve core area, F e is the electromagnetic output force, m2 is the mass of the armature assembly and the valve core, D is the damping coefficient, K is the spring stiffness coefficient, where F e =K I I(s)-K x x1(s), K I is the valve current gain, K x is the displacement force gain. Due to the characteristics of the proportional electromagnet, K x Very small and negligible.

[0045] Furthermore, based on the current dynamic equation, flow equation, and valve core motion equation of the proportional relief valve of the shield hydraulic propulsion system, the overall transfer function of the proportional relief valve is obtained as follows:

[0046] Where, Q is the proportional relief valve flow rate, x1 is the valve core displacement, p L is the equivalent load pressure, I is the current, ω2 and ζ2 are the natural frequency and damping ratio of the proportional relief valve respectively. ω e is the turning frequency of the proportional relief valve, ω e =(r c +r p )L.

[0047] According to another aspect of the present invention, a shield thrust pressure control device is provided, which includes a computer-readable storage medium and a processor. The computer-readable storage medium includes a stored computer program. When the computer program is executed by the processor, it controls the device to execute the shield thrust pressure control method as described above to control the shield thrust pressure.

[0048] In general, the above technical solutions conceived by the present invention, compared with the prior art, provide a shield tunneling pressure control method and device with the following advantages:

[0049] 1. The calculation process of the ultimate support pressure solution of the spiral failure model of the composite stratum excavation face under soil uncertainty conditions is combined with the shield thrust control simulation model as an external load module. By collecting, analyzing and calculating the soil parameters of the shield excavation face under graded reliability indicators, the appropriate thrust force is obtained to balance the external soil pressure, ensuring the safe, efficient and effective control of the excavation process.

[0050] 2. Based on the influence of uncertainty in soil parameters of composite strata on the stability of shield excavation face, a three-dimensional rotating body failure structure model of composite stratum excavation face is constructed. The model is accurate and advanced, and the failure form of the excavation face is analyzed.

[0051] 3. The control model adopts fuzzy PID control algorithm for propulsion control, which uses fuzzy PID online adjustment factor parameters to improve the control effect. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] Figure 1 This is a schematic flow chart of a shield tunneling propulsion pressure control method provided by the present invention;

[0053] Figure 2 It is the failure structure diagram of the excavation surface of the composite stratum;

[0054] Figure 3 (a), (b), (c), and (d) are discrete schematic diagrams of the composite stratum helical failure model;

[0055] Figure 4 It is a failure mode curve diagram based on the Monte Carlo method;

[0056] Figure 5 is the ultimate support pressure diagram under different reliability levels;

[0057] Figure 6 This is a simulation model diagram of the shield propulsion pressure control system based on fuzzy PID control;

[0058] Figure 7 It is a control simulation diagram of the propulsion working pressure. DETAILED DESCRIPTION

[0059] In order to make the objectives, technical solutions and advantages of the present invention more clearly understood, the present invention is further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely for the purpose of explaining the present invention and are not intended to limit the present invention. In addition, the technical features involved in the various embodiments of the present invention described below may be combined with each other as long as they do not conflict with each other.

[0060] See also Figure 1 and Figure 2 The present invention provides a shield propulsion pressure control method, which mainly includes the following steps:

[0061] (1) Based on the influence of the uncertainty of soil parameters of composite strata on the stability of shield excavation face, a three-dimensional rotating body failure structure model of composite stratum excavation face is constructed. The ultimate support force calculation formula when the soil structure is destroyed under the critical failure state of the excavation face is derived by equating the work done by the external force of the system with the energy dissipated by friction inside the soil.

[0062] Among them, based on the influence of the uncertainty of soil parameters of composite strata on the stability of shield excavation face, a three-dimensional rotating body destruction structure model of composite stratum excavation face is constructed. The model is accurate and advanced, and the destruction form of the excavation face is analyzed by this. At the same time, the calculation formula of the ultimate support force when the soil structure is destroyed under the critical destruction state of the excavation face is derived by equating the work done by the external force of the system with the energy dissipated by friction inside the soil. The spatial discretization technology is applied to generate the three-dimensional destruction surface of the failure model point by point to obtain the coordinates of each point in the instability velocity field of the excavation face; after constructing the instability model of the composite stratum excavation face, the ultimate support pressure can be calculated based on the power and energy consumption of the excavation working face using the upper limit theorem of limit analysis; the power and energy consumption of the system mainly include the weight of the soil, the support pressure and the friction dissipation of the velocity discontinuity surface. The extension lines of OA and OE and the stratum boundaries can divide the destruction structure of the composite stratum excavation face into four regions: I, II, III and IV, such as Figure 2 The calculation method is as follows:

[0063] The polar coordinate expressions of the three logarithmic spiral curves AF, EF and BE are:

[0064]

[0065]

[0066]

[0067] Among them, r A 、r B and r E are the lengths of OA, OB and OE respectively, θ A ,θB and θ E are the angles OA, OB, and OE included with the polar coordinate system (r, θ) relative to the polar axis direction. Their specific values ​​and calculation expressions are as follows:

[0068]

[0069]

[0070]

[0071]

[0072]

[0073] Since the elevation difference between point E and point O is Y O +H, so θ can be obtained from formula 9 E The value of .

[0074]

[0075] The θ obtained from Equation 9 is E Substituting into formula 6, we can get r E Since point F is the intersection of the two logarithmic spiral curves in the upper soil, we can obtain equation (10) by combining equations (1) and (2), θ F is obtained when θ takes a value such that equation (10) holds.

[0076]

[0077] In summary, the instability velocity field of the excavation face is completely determined by Z O and Y O To simplify the reasoning, Z O and Y O The two parameters are dimensionless through the polar coordinates of the tunnel center point G, with r G / D and θ G The two parameters represent the coordinates of point O in the rectangular coordinate system of the (Y, Z) plane as (Y O ,Z O )as follows:

[0078]

[0079]

[0080] Then, the spatial discretization technology is applied to generate the three-dimensional failure surface of the failure model point by point to obtain the coordinates of each point in the instability velocity field of the excavation surface, such as Figure 3 As shown, the calculation method is as follows:

[0081] The first radial plane contains O, A j and A′ j Three points, their angle θ with the Y axis j It can be expressed as:

[0082]

[0083] The second part of the radial plane can be discretized into m, and the final closed plane is Π jmax , two adjacent radial planes π j and Π j-1 The rotation angle δ around the OX axis θ for:

[0084]

[0085] Therefore, in the second partial plane Π j and Π j-1 The angle relationship between the polar coordinate system (r, θ) and the polar axis direction is:

[0086] θ j =θ j-1 +δ θ (16)

[0087] A local rectangular coordinate system (C j ,x j ,y j ), where C j is the corresponding plane Π j With point O as the center, r F The intersection point of the circles with the radius x j The direction of the axis is the same as the X-axis direction, j Axis direction and C j O direction is the same as C j The coordinates in the (X, Y, Z) coordinate system are:

[0088]

[0089] Each radial plane in the first part is as Figure 3 As shown, the local coordinate system inside the plane (C j ,x j ,y j ) j and A′ j The corresponding β′ j and β″ j Both vectors can be used and Obtain.

[0090]

[0091] in:

[0092]

[0093]

[0094] After constructing the instability model of the excavation face of the composite stratum, the upper limit theorem of limit analysis can be applied to calculate the ultimate support pressure based on the power and energy consumption on the excavation face. The power and energy consumption of the system are mainly composed of the soil deadweight σ γ , support pressure σ T And the velocity discontinuity friction force f dissipation, the calculation formula is as follows:

[0095]

[0096]

[0097]

[0098] W f =W γ +W T (twenty four)

[0099] Among them, r A 、r B and r E are the lengths of OA, OB and OE respectively, θ A ,θ B and θ E are the angles between OA, OB and OE and the polar coordinate system (r,θ) relative to the polar axis direction, and are the cohesion, internal friction angle and soil weight of the upper and lower soil layers respectively, ω is the angular velocity, W f is the friction dissipation power, W γ is the soil gravity power, W T is the work power of support pressure.

[0100] (2) A reliability analysis index system for excavation face stability is constructed, and the Monte Carlo method is used for simulation calculation to obtain the values ​​of the ultimate support pressure of composite strata under uncertainty conditions of soil parameters under different reliability indicators.

[0101] Among them, since the soil parameters are obtained through in-situ tests or indoor tests, and even the same soil type is by no means an "average" material, the uncertainty of the soil parameters will cause differences between the observed values ​​and the true values ​​of each link in the excavation construction process. Based on the actual shield excavation construction process, the failure probability of the system is used as the basis for judging the risk state of the excavation face stability to confirm the risk state level of the excavation face stability structure. Among them, each stratum change point can use the Monte Carlo simulation method to simulate 1000 groups of corresponding soil parameter vectors according to the corresponding physical and mechanical properties of the soil layer and the characteristics of their respective uncertainties. Each group of soil parameter vectors can correspond to a composite stratum excavation face instability model, determine the failure probability of the system, and obtain the failure curve of the model excavation face, such as Figure 4 Then the ultimate support pressure calculation under the relevant soil uncertainty is carried out. Figure 5 The calculation formula is as follows:

[0102] p f =P(g(σ T )<0) (25)

[0103]

[0104] β=Φ -1 (1-p f ) (27)

[0105] Among them, p f is the failure probability of the system, σ T is the simulated value of the ultimate support pressure of the shield face, represents the limit value of the ultimate support pressure of the shield surface, β is the reliability index, Φ -1 (·) is the inverse function of the standard normal distribution function.

[0106] The full dynamic range that needs to be controlled is obtained through the ultimate support pressure of the same level, that is, an algorithm is constructed to control the thrust value output by the shield thrust system within this range.

[0107] (3) Construct a load calculation module for the shield propulsion system and establish a mathematical model of the shield hydraulic propulsion system, and then construct a control model for the shield propulsion system. The control model uses a fuzzy PID control algorithm to control the propulsion force; among them, the control model takes the input voltage of the proportional overflow valve of the shield hydraulic propulsion system as input and the propulsion pressure of the shield propulsion system as output.

[0108] In actual control, load changes or external interference will cause the parameters of the controlled object to change, and the nonlinear characteristics of the shield propulsion pressure make it impossible for conventional PID to achieve real-time adjustment under environmental changes, which has certain limitations. Therefore, fuzzy PID can be used to adjust the factor parameters online to improve the control effect. Construct a load calculation module for the shield propulsion system and establish a mathematical model of the shield hydraulic propulsion system. Establish a test model of the shield thrust control system in MATLAB / SIMULINK. The input is the input voltage of the proportional relief valve, and the output is the thrust pressure of the system. Figure 6 As shown, the input is the input voltage of the proportional relief valve, and the output is the propulsion pressure of the system. The relevant calculations are as follows:

[0109] 1. The mechanical balance equation of the output force and load force of the propulsion hydraulic cylinder is:

[0110]

[0111] Where A1 and A2 are the areas of the rodless and rod pistons, respectively; p1 and p2 are the pressures of the rodless and rod pistons, respectively; p2 is very small and can be ignored; m1 is the total mass of the system; x is the displacement of the piston rod; and F L is the external load force acting on the piston, B1 is the damping coefficient, and K is the stiffness coefficient.

[0112] 2. The flow continuity equation of the propulsion hydraulic cylinder is:

[0113]

[0114] Among them, Q L is the flow rate of the rodless cavity, λ1 is the leakage coefficient of the hydraulic cylinder, V1 is the volume of the hydraulic cylinder oil inlet cavity, and β1 is the effective volume elastic modulus.

[0115] After Laplace transformation of equations (28) and (29), we can obtain:

[0116] A1P1(s)=(m1s 2 +B1s+K)X(s)+F L (s) (30)

[0117]

[0118] By simplifying equations (30) and (31) and eliminating the intermediate variable X, we can obtain the output hydraulic cylinder rodless chamber pressure P1 with respect to the input flow Q L and external load force F L The transfer functions are:

[0119]

[0120]

[0121] Where ω1 and ζ1 are the comprehensive natural frequency and comprehensive damping ratio, respectively, and their values ​​are as follows:

[0122]

[0123]

[0124] The transfer function of the proportional relief valve is calculated as follows:

[0125] 1. The current dynamic equation of the proportional relief valve is:

[0126]

[0127] Where u is the valve solenoid coil input voltage, L is the coil inductance, r c and r p are the internal resistances of the coil and amplifier, respectively, K v is the back electromotive force coefficient of the coil, and x1 is the valve core displacement.

[0128] After Laplace transformation of Equation 36, we can get:

[0129] U(s)=LsI(s)+(r c +r p )I(s)+K v sX1(s) (37)

[0130] 2. The flow equation of the proportional relief valve is:

[0131]

[0132] Where Q is the proportional relief valve flow rate, α d is the flow coefficient of the cone valve port, d1 is the aperture of the cone valve port, β is the semi-cone angle of the cone valve, ρ is the density of the hydraulic oil in the system, and p L is the equivalent load pressure.

[0133] Assume that the stable operating point of the proportional relief valve is q0(x0,p L0 ), then Equation 38 can be linearized at q0 to obtain:

[0134] ΔQ=K q Δx1-K p Δp L (39)

[0135] Where, and Because K p It is very small and can be ignored.

[0136] After Laplace transformation of Equation 39, we can get:

[0137] Q(s)=K q X1(s) (40)

[0138] 3. The valve core motion equation of the proportional relief valve is:

[0139]

[0140] Where A m is the valve core area, F e is the electromagnetic output force, m2 is the mass of the armature assembly and the valve core, D is the damping coefficient, and K is the spring stiffness coefficient.

[0141] After Laplace transformation of Equation 41, we can get:

[0142] A m P1(s)-F e (s)=(m2s 2 +Ds+K)X1(s) (42)

[0143] Among them F e =K I I(s)-K x x1(s), K I is the valve current gain, K x is the displacement force gain. Due to the characteristics of the proportional electromagnet, K x is very small and can be ignored, so Equation 42 can be transformed to:

[0144] F e (s) = K I I(s) (43)

[0145] Since the influence of the coil induced back electromotive force can be ignored and the turning frequency of the proportional relief valve is large, the second-order oscillation link will play a dominant role. Therefore, the overall transfer function of the proportional relief valve can be simplified as follows:

[0146]

[0147] Substituting this calculation formula into the propulsion pressure control system, the required propulsion pressure can be calculated.

[0148] Among them, regarding the shield face thrust as the external load value range for thrust pressure control modeling, a shield thrust pressure control algorithm was designed. By analyzing the main components of the shield thrust system, the respective transfer functions were calculated according to their respective mechanical equilibrium equations, current dynamic equations and flow equations and built in the SIMULINK module of MATLAB software; then actual experimental verification was carried out, and finally the fuzzy PID control algorithm was applied to the shield thrust control to obtain its actual effect on shield thrust control, realizing the control of shield thrust pressure during shield tunneling, such as Figure 7 As shown in the results, it can be seen that fuzzy PID control can ensure that the shield propulsion working pressure is within the set working pressure range, so this shield propulsion pressure control method is feasible and effective.

[0149] It will be easily understood by those skilled in the art that the above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.

Claims

1. A shield tunneling pressure control method, characterized in that: The control method comprises the following steps: (1) Based on the influence of uncertainty in soil parameters of composite strata on the stability of shield excavation face, a three-dimensional rotating body failure structure model of composite stratum excavation face is constructed. The ultimate support force calculation formula when the soil structure is destroyed under the critical failure state of the excavation face is derived by equating the work done by the external force of the system with the energy dissipated by friction in the soil. (2) Construct a reliability analysis index system for excavation face stability and apply the Monte Carlo method to perform simulation calculations to obtain the values ​​of the ultimate support pressure of composite strata under uncertainty conditions of soil parameters under different reliability indicators; (3) Construct a load calculation module for the shield propulsion system and establish a mathematical model of the shield hydraulic propulsion system. Then, construct a control model for the shield propulsion system. The control model uses a fuzzy PID control algorithm to control the propulsion force. The control model uses the input voltage of the proportional relief valve of the shield hydraulic propulsion system as input and the propulsion pressure of the shield propulsion system as output. Using spatial discretization technology, a three-dimensional failure surface of the failure model is generated point by point to obtain the coordinates of each point in the excavation face instability velocity field. After constructing the composite stratum excavation face instability model, the ultimate support pressure is calculated based on the power and energy consumption of the excavation working face using the upper bound theorem of limit analysis. System power and energy consumption are determined by soil deadweight σ γ , support pressure σ T And the velocity discontinuity friction force f dissipation is calculated as: IN f =In γ +W T Among them, r A 、r B and r E are the lengths of OA, OB and OE respectively, θ A ,θ B and θ E are the angles between OA, OB and OE and the polar coordinate system (r,θ) relative to the polar axis direction, and are the cohesion, internal friction angle and soil weight of the upper and lower soil layers respectively, ω is the angular velocity, W f is the friction dissipation power, W γ is the soil gravity power, W T is the work power of support pressure.

2. The shield propulsion pressure control method according to claim 1, characterized in that: The extension lines of OA and OE and the stratum boundary can divide the failure structure of the excavation surface of the composite stratum into four regions: I, II, III and IV; the instability velocity field of the excavation surface is completely determined by Z O and Y O Take the value and set Z O and Y O The two parameters are dimensionless through the polar coordinates of the tunnel center point G, with r G / D and θ G The two parameters represent the coordinates of point O in the rectangular coordinate system of the (Y, Z) plane as (Y O ,Z O ), the corresponding formula is: The spatial discretization technology is applied to generate the three-dimensional failure surface of the failure model point by point to obtain the coordinates of each point in the excavation surface instability velocity field; among them, point A is the top of the tunnel, point B is the bottom of the tunnel, point A is the origin of the rectangular coordinate system of the (Y, Z) plane, and the excavation surface instability velocity field is composed of three logarithmic spiral curves AF, EF and BE, among which point E is the intersection of the logarithmic spiral curve and the stratum boundary line, and point F is the intersection of the two logarithmic spiral curves in the upper soil. The entire instability velocity field rotates with point O as the origin, point O is the origin of the polar coordinate system (r, θ), the polar axis is perpendicular to the surface plane and downward, and OA, OB, and OE are the lines connecting each point respectively.

3. The shield propulsion pressure control method according to claim 1, characterized in that: The failure probability of the system is used as the basis for judging the risk state of the excavation face stability to confirm the risk state level of the excavation face stability structure. Each stratum change point is simulated using the Monte Carlo simulation method to simulate 1000 sets of corresponding soil parameter vectors based on the physical and mechanical properties of the corresponding soil layer and the characteristics of their respective uncertainties. Each set of soil parameter vectors can correspond to a composite stratum excavation face instability model, and the failure probability of the system is determined to obtain the failure curve of the model excavation face. Then, the ultimate support pressure under soil uncertainty is calculated.

4. The shield propulsion pressure control method according to claim 3, characterized in that: Formula used to calculate the ultimate support pressure include: p f JP(g(σ). T )<0) β=Φ -1 (1-p f ) Where p f is the failure probability of the system, σ T is the simulated value of the ultimate support pressure of the shield face, represents the limit value of the ultimate support pressure of the shield surface, β is the reliability index, Φ -1 (·) is the inverse function of the standard normal distribution function.

5. The shield propulsion pressure control method according to claim 4, characterized in that: The full dynamic range that needs to be controlled is obtained through the ultimate support pressure of the same level, that is, an algorithm is constructed to control the thrust value output by the shield thrust system within the full dynamic range.

6. The shield propulsion pressure control method according to any one of claims 1 to 5, characterized in that: Based on the mechanical balance equation of the output force and load force of the propulsion hydraulic cylinder of the shield hydraulic propulsion system and the flow continuity equation of the propulsion hydraulic cylinder, the output hydraulic cylinder rodless chamber pressure P1 with respect to the input flow Q is obtained. L and external load force F L The transfer functions are: Where ω1 and ζ1 are the comprehensive natural frequency and comprehensive damping ratio respectively, and their values ​​are: Where A1 is the area of ​​the rodless piston, p1 is the pressure of the rodless chamber, m1 is the total reduced mass of the system, and F L is the external load force acting on the piston, β1 is the effective bulk elastic modulus, λ1 is the leakage coefficient of the hydraulic cylinder, and V1 is the volume of the oil inlet chamber of the hydraulic cylinder.

7. The shield propulsion pressure control method according to any one of claims 1 to 5, characterized in that: Based on the current dynamic equation, flow equation, and valve core motion equation of the proportional relief valve of the shield hydraulic propulsion system, the overall transfer function of the proportional relief valve is obtained as follows: Where, Q is the proportional relief valve flow rate, x1 is the valve core displacement, p L is the equivalent load pressure, I is the current, ω2 and ζ2 are the natural frequency and damping ratio of the proportional relief valve respectively. .

8. A shield thrust pressure control device, characterized by: The device includes a computer-readable storage medium and a processor, wherein the computer-readable storage medium includes a stored computer program. When the computer program is run by the processor, the device is controlled to execute the shield propulsion pressure control method according to any one of claims 1 to 7 to control the shield propulsion pressure.

Images