Calculation Method of Loose Pressure in Surrounding Rock Mass

: With the implementation of the regional coordinated development strategy, trafﬁc ﬂow has grown explosively. The construction of a larger tunnel section becomes an effective way to solve the highway network’s insufﬁcient transport capability problem. Currently, there is little research on the factors inﬂuencing and methods of calculating the loose pressure in the surrounding rock mass for highway tunnels with super-large cross-sections. Based on the Bifurcation Tunnel, which is one of the sign projects in the past ﬁve years, this paper discusses the inﬂuencing factors for the range of loose zone in deeply buried tunnels using a combination of a numerical analysis and an orthogonal test. The weight of inﬂuencing factors is calculated via an efﬁciency evaluation method. This paper establishes a limit analysis model of the loose pressure in the surrounding rock mass under a non-linear failure criterion based on the ﬁtted boundary function and upper bound limit analysis method and deduces the correlations of the loose pressure. The distribution law of the loose pressure, obtained via the limit analysis method, is consistent with the pressure-monitoring results, verifying the correctness of the proposed calculation method. This study can provide a calculation basis for the design of a supporting structure and the selection of similar super-section tunnel projects.


Introduction
Owing to rapid economic development, logistics, and the vehicle industry, the construction of tunnels with large-span and extra-large-span sections has boomed in recent years and serves as the optimum scheme to improve road capacity and meet the existing travel demand.Compared to a single-hole project with a small span, an extra-large-spansection tunnel consists of a much larger span, a lower high-span ratio, and a greater construction risk.There are no corresponding design codes and construction standards for these projects.Although some construction experience has been obtained, theoretical research is relatively lacking.
To guarantee the rationality of the supporting structure design for super-span tunnels, the first problem to be solved is to clarify the method of predicting the loose pressure in the surrounding rock mass.At present, the empirical formulas commonly used for the loose pressure are divided into three categories: the equivalent thickness of the crown, statistics of collapse law, and surrounding rock classification, as shown in Figure 1.The equivalent thickness of the crown takes the gravity of the rock and soil mass above the vault as the primary influence for the loose pressure in the surrounding rock.Typical methods include Protodyakonov's theory [1], the Bierbaumer theory [2], the Terzaghi theory [3], and Xie Jiaxiao's theory [4].The statistics of collapse law method is based on practical engineering and is the leading engineering practice used in design codes in China [5][6][7][8].
The correlations are relatively simple and consider the influence of the surrounding rock classification.The surrounding rock classification method derives the surrounding rock loose pressure from the RQD, Q, and RMR systems [9][10][11][12].The correlations and applicable conditions of various calculation methods are presented in Table 1.
Appl.Sci.2023, 13, x FOR PEER REVIEW 2 of 22 theory [3], and Xie Jiaxiao's theory [4].The statistics of collapse law method is based on practical engineering and is the leading engineering practice used in design codes in China [5][6][7][8].The correlations are relatively simple and consider the influence of the surrounding rock classification.The surrounding rock classification method derives the surrounding rock loose pressure from the RQD, Q, and RMR systems [9][10][11][12].The correlations and applicable conditions of various calculation methods are presented in Table 1.
Small span Neglects the effect of geometrical dimension Bhasin and Grimstad correlation [10] Broken rock condition Neglects the effect of height
There are two deficiencies of the above-mentioned empirical formulas: (1) the calculations are derived or statistically obtained under the conditions of a small-span tunnel, and (2) the influencing factors are significantly different.The above two deficiencies restrict the applications of empirical formulas for loose pressure in the extra-large-span cross-section tunnels.
In addition to the empirical formulas, the upper bound theorem of limit analysis with a strictly theoretical basis has also been used recently to calculate the loose pressure in the surrounding rock mass.It focuses on the ultimate failure state of the rock-soil mass without considering the elastic-plastic deformation process.In addition, this method is not restricted by the span.Some scholars have published papers on the surrounding rock pressure of soil tunnels.Atkinson et al. [13] discussed the stability of a shallow circular tunnel in a cohesionless soil through the upper bound theorem of limit analysis together with results of tests performed in sand soil.Davis et al. [14] deduced the upper boundary stability solutions for a shallow tunnel excavation at three different positions under undrained conditions and proposed a new stability evaluation method.In order to analyze the stability of an excavation face, Mollon et al. [15,16] established the spatial discretization three-dimensional calculation method based on the limit analysis method, and the results of the calculation approached the real situation.Based on the Hoek-Brown failure criterion, Huang et al. [17] presented an upper bound solution for the collapse shape for a circular tunnel considered the influence of different geological conditions and pore water pressure.Lei et al. [18] established explicit expressions to calculate the loose pressure on a shallow, buried tunnel based on the upper bound theorem of limit analysis and the Terzaghi failure mode.Based on the spatial discretization model, Pam et al. [19][20][21] studied the impact mechanism of a heterogeneous stratum, pore water pressure, and a non-circular excavation face on the tunnel face stability.Zhang et al. [22] proposed an analytical method for calculating the ultimate support force of an excavation face in the clay layer and analyzed the local and whole stability of the excavation face for a largediameter shield tunnel in practical engineering.According to the tunnel face failure mode in purely cohesive soil proposed by Mollon [23], Huang et al. [24] established an analytical model considering the shear failure boundary which further improved the calculation accuracy of the upper bound solution for tunnel face stability in non-homogeneous and anisotropic undrained clay.To assess the stability of a tunnel face in inclined strata, Tu et al. [25] optimized the 3D rotating rigid body collapse mechanism by improving the spatial discretization technique.However, the failure mechanism constructed in the current research is highly subjective and lacks numerical simulation and field test verification.Such flaws can be solved by combining a numerical analysis and an upper bound limit analysis.With its lower cost, shorter period, and better repeatability, a numerical simulation is helpful for predicting the failure position [26,27].Wang et al. [28,29] verified the rationality of the analytical solution for the collapse curve via a vertical displacement cloud picture obtained from a simulation.
In view of the shortcomings of the existing research, this paper discussed the impact of rock conditions and cavern size on the range of the loose pressure zone by combining a numerical analysis and an orthogonal test.The investigation is based on the Bifurcation Tunnel.The limit analysis model of the loose pressure in the surrounding rock is established using the fitted boundary function of the loose zone and the upper bound limit analysis.Based on this model, this paper deduces the upper bound solution for a deep tunnel's loose pressure, which is verified by pressure-monitoring results obtained for the different sections.The results of this research can provide new ideas for the selection of a supporting structure for similar super-section tunnel projects.

Project Overview
The Shenzhen Transit Expressway is a major construction project in Shenzhen, Guangdong Province, China.It begins at Luosha Road, connects with the Hong Kong highway network through the port, and ends at the Huiyan expressway.The total length of the route is 32.5 km, while the total investment is about CNY 5.792 billion.It will become a critical traffic artery serving the development of the eastern coastal areas of South China after its completion.
The eastern transit passage of Shenzhen City starts from Luohu Port, connects with the Hong Kong highway network through the port, and ends at the Huiyan expressway.
The project is 32.5 km long in total, and it will become a critical traffic artery serving the development of the eastern coastal areas of South China after its completion.The passage planning is shown in Figure 2.
The Shenzhen Transit Expressway is a major construction project in Shenzhen, Guangdong Province, China.It begins at Luosha Road, connects with the Hong Kong highway network through the port, and ends at the Huiyan expressway.The total length of the route is 32.5 km, while the total investment is about CNY 5.792 billion.It will become a critical traffic artery serving the development of the eastern coastal areas of South China after its completion.
The eastern transit passage of Shenzhen City starts from Luohu Port, connects with the Hong Kong highway network through the port, and ends at the Huiyan expressway.The project is 32.5 km long in total, and it will become a critical traffic artery serving the development of the eastern coastal areas of South China after its completion.The passage planning is shown in Figure 2. The Bifurcation Tunnel is a key control project in the entire transit expressway line and has a total length of 1435 m [30].It is composed of three parts: the mainline section, the municipal section, and the port connecting section, of which the mainline section is divided into left and right lines.The Bifurcation Section is located at the entrance of the left line.The large-span section includes the maximum section, Gradient Section 1, Gradient Section 2, and Gradient Section 3. The Bifurcation Section is shown in Figure 3.The Bifurcation Tunnel is a key control project in the entire transit expressway line and has a total length of 1435 m [30].It is composed of three parts: the mainline section, the municipal section, and the port connecting section, of which the mainline section is divided into left and right lines.The Bifurcation Section is located at the entrance of the left line.The large-span section includes the maximum section, Gradient Section 1, Gradient Section 2, and Gradient Section 3. The Bifurcation Section is shown in Figure 3.
The Shenzhen Transit Expressway is a major construction project in Shenzhen, Guangdong Province, China.It begins at Luosha Road, connects with the Hong Kong highway network through the port, and ends at the Huiyan expressway.The total length of the route is 32.5 km, while the total investment is about CNY 5.792 billion.It will become a critical traffic artery serving the development of the eastern coastal areas of South China after its completion.
The eastern transit passage of Shenzhen City starts from Luohu Port, connects with the Hong Kong highway network through the port, and ends at the Huiyan expressway.The project is 32.5 km long in total, and it will become a critical traffic artery serving the development of the eastern coastal areas of South China after its completion.The passage planning is shown in Figure 2. The Bifurcation Tunnel is a key control project in the entire transit expressway line and has a total length of 1435 m [30].It is composed of three parts: the mainline section, the municipal section, and the port connecting section, of which the mainline section is divided into left and right lines.The Bifurcation Section is located at the entrance of the left line.The large-span section includes the maximum section, Gradient Section 1, Gradient Section 2, and Gradient Section 3. The Bifurcation Section is shown in Figure 3.The burial depth of the large-span bifurcation section is 78~106 m.The stratum at the tunnel site consists of artificial soil, sandstone, flood alluvial soil, and rhyolite.The rock mass is relatively complete, with only slightly developed fractures to which the degree of development increases with depth.The geological section is shown in Figure 4, and photos of the drill core are shown in Figure 5.The geological investigation, combined with the rock mass type, weathering degree, and joint development, provides the GSI (geological strength index) value of the Bifurcation Section, as shown in Table 2. Furthermore, the maximum excavation span is 30.01 m, and the maximum excavation area is 428.5 m 2 ; globally, this is a rare extra-large section highway tunnel.Therefore, the key to the construction project's structural design and safety is determining the influencing factors and distribution pattern of the loose pressure of the surrounding rock for the large-span section.
tunnel site consists of artificial soil, sandstone, flood alluvial soil, and rhyolite.The rock mass is relatively complete, with only slightly developed fractures to which the degree of development increases with depth.The geological section is shown in Figure 4, and photos of the drill core are shown in Figure 5.The geological investigation, combined with the rock mass type, weathering degree, and joint development, provides the GSI (geological strength index) value of the Bifurcation Section, as shown in Table 2. Furthermore, the maximum excavation span is 30.01 m, and the maximum excavation area is 428.5 m 2 ; globally, this is a rare extra-large section highway tunnel.Therefore, the key to the construction project's structural design and safety is determining the influencing factors and distribution pattern of the loose pressure of the surrounding rock for the large-span section.The surrounding rock is slightly weathered sandstone, and the rock mass is fractured with developed joint fissures.Generally, the groundwater is poor, and the permeability is low.

130~244 114 m
The surrounding rock is a greenish-gray, medium, and slightly weathered sandstone.The rock mass is relatively fractured with developed joint fissures, and the rock core is short, columnar, and massive.The  The burial depth of the large-span bifurcation section is 78~106 m.The stratum at the tunnel site consists of artificial soil, sandstone, flood alluvial soil, and rhyolite.The rock mass is relatively complete, with only slightly developed fractures to which the degree of development increases with depth.The geological section is shown in Figure 4, and photos of the drill core are shown in Figure 5.The geological investigation, combined with the rock mass type, weathering degree, and joint development, provides the GSI (geological strength index) value of the Bifurcation Section, as shown in Table 2. Furthermore, the maximum excavation span is 30.01 m, and the maximum excavation area is 428.5 m 2 ; globally, this is a rare extra-large section highway tunnel.Therefore, the key to the construction project's structural design and safety is determining the influencing factors and distribution pattern of the loose pressure of the surrounding rock for the large-span section.

Mileage Length Surrounding Rock Characteristics GSI 0~130 130 m
The surrounding rock is slightly weathered sandstone, and the rock mass is fractured with developed joint fissures.Generally, the groundwater is poor, and the permeability is low.

130~244 114 m
The surrounding rock is a greenish-gray, medium, and slightly weathered sandstone.The rock mass is relatively fractured with developed joint fissures, and the rock core is short, columnar, and massive.The   The surrounding rock is slightly weathered sandstone, and the rock mass is fractured with developed joint fissures.Generally, the groundwater is poor, and the permeability is low.

130~244 114 m
The surrounding rock is a greenish-gray, medium, and slightly weathered sandstone.The rock mass is relatively fractured with developed joint fissures, and the rock core is short, columnar, and massive.The groundwater is poor and does not have the conditions for water inrush.

Study of Factors Influencing the Loose Zone in the Surrounding Rock Mass
The formation of the loose zone due to excavation is an inherent characteristic of the surrounding rock mass.The scope of the loose zone directly determines the loose pressure in the surrounding rock mass.Based on the orthogonal test and the finite difference software FLAC3D, this paper designs multiple group model tests to determine the influence of various factors, such as the surrounding rock quality, cavern size, and tunnel burial depth, on the loose zone of an extra-large section tunnel.

Orthogonal Analysis
The orthogonal experimental design method, proposed by mathematical statistician R.A. Fisher in 1935, has been widely used in mathematics, economics, agriculture, aerospace, and other fields.Compared with the comprehensive test method, it has several advantages, such as fewer experiments, careful consideration of factors, etc.According to the orthogonality, representative test points are selected from all level combinations of all test factors.These test points have uniform dispersion and neat comparability.
For example, if an event has three factors (A, B, C), each factor has three levels (A i , B i , C i i = 1, 2, 3).A comprehensive test must be designed for 3 3 = 27 combination tests, while only nine tests are needed to analyze the degree of influence of each factor (A, B, C) via the orthogonal table L 9 (3 3 ).The distribution of the comprehensive and orthogonal test points is shown in Figure 6.
The formation of the loose zone due to excavation is an inherent characteristic of the surrounding rock mass.The scope of the loose zone directly determines the loose pressure in the surrounding rock mass.Based on the orthogonal test and the finite difference software FLAC3D, this paper designs multiple group model tests to determine the influence of various factors, such as the surrounding rock quality, cavern size, and tunnel burial depth, on the loose zone of an extra-large section tunnel.

Orthogonal Analysis
The orthogonal experimental design method, proposed by mathematical statistician R.A. Fisher in 1935, has been widely used in mathematics, economics, agriculture, aerospace, and other fields.Compared with the comprehensive test method, it has several advantages, such as fewer experiments, careful consideration of factors, etc.According to the orthogonality, representative test points are selected from all level combinations of all test factors.These test points have uniform dispersion and neat comparability.
For example, if an event has three factors ( , ,  ), each factor has three levels ( ,  ,   = 1, 2, 3).A comprehensive test must be designed for 3 3 = 27 combination tests, while only nine tests are needed to analyze the degree of influence of each factor (, , ) via the orthogonal table  (3 ).The distribution of the comprehensive and orthogonal test points is shown in Figure 6.

Hoek-Brown (H-B) Criterion
The H-B criterion, which perfectly reflects the non-linear failure characteristics of rocks and rock masses, was proposed by Hoek and Brown in 1980.After years of continuous development, it has become one of the most accepted criteria in the research on rock mass stability analysis and strength prediction [31].The expression of the criterion is as follows: where  and  are the maximum and minimum principal stress of the rock mass, respectively;  is the uniaxial compressive strength of the rock specimen;  , , and  are the material parameters of the rock mass related to lithology and the structural plane, which can be determined by the correlation:

Numerical Models and Simulation Scheme 3.2.1. Hoek-Brown (H-B) Criterion
The H-B criterion, which perfectly reflects the non-linear failure characteristics of rocks and rock masses, was proposed by Hoek and Brown in 1980.After years of continuous development, it has become one of the most accepted criteria in the research on rock mass stability analysis and strength prediction [31].The expression of the criterion is as follows: where σ 1 and σ 3 are the maximum and minimum principal stress of the rock mass, respectively; σ ci is the uniaxial compressive strength of the rock specimen; m b , s, and a are the material parameters of the rock mass related to lithology and the structural plane, which can be determined by the correlation: where m i is the material constant reflecting the softness and hardness of the rock mass.D is the parameter reflecting the disturbance of rock mass, and its value ranges from 0 to 1.

Basic Assumptions
Numerical tests were conducted using the FLAC3D software.The following assumptions were made: (1) the models were established based on the plain strain hypothesis; (2) the surrounding rock mass was regarded as an ideal elastic-plastic body, following the H-B failure criterion; (3) the tectonic stresses were neglected in the numerical simulations; (4) the function of the supporting structure was not considered; (5) the model adopted the full-face method without considering the construction process and groundwater.

Simulation Scheme
Based on the traditional research method of loose zones in a surrounding rock mass, the following influencing factors were determined as the primary research objects: tunnel burial depth (H), tunnel span (B), tunnel height-span ratio (λ), lateral pressure coefficient (K), and GSI.There were five influencing factors in total, and each factor had four levels.This paper established 16 numerical models to carry out the orthogonal test in line with the orthogonal table L 16 4 5 .
According to the geological conditions and cross-sectional dimensions of the largespan section of the bifurcation, H was 80~110 m, B was 20~35 m, λ was 0.6~1.2,K is 0.5~2.0, and the GSI was 25~40.The orthogonal table and the rock mass mechanical parameters corresponding to GSI are shown in Tables 3 and 4, respectively.The bulk density and Poisson's ratio of the rock mass were identified by geological exploration reports.The σ ci value was mainly obtained from a laboratory test, and m b and s can be evaluated via Correlation 2. The GSI-based empirical equations [32,33] used to determine the elastic modulus E and m i are proposed as follows: The numerical calculation model of the unit thickness was set as: height × width = 180 m × 300 m (210 m × 300 m, 240 m × 300 m, 270 m × 300 m).The element mesh and model dimension are shown in Figure 7.  represents the vertical distance from the tunnel vault to the ground.The model's boundary conditions were set so that the upper boundary was a free surface boundary, and both sides and the bottom edges were constrained by the method of normal displacement.The initial in situ stress was the self-weight stress of the rock mass.

Method of Determining the Loose Zone in the Surrounding Rock Mass
The pressure arch in the surrounding rock mass was characterized by an arch structure.The pressure arch body bears its own load and the load of the rock mass above the arch so that the surrounding rock maintains a certain degree of self-stability after excavation.However, the rock mass below the pressure arch has poor stability because it is close to the tunnel excavation outline, which leads to the stress exceeding the rock mass strength and causing damage.Therefore, the loose zone in the surrounding rock mass is defined as the rock mass from the inner boundary of the pressure arch to the tunnel contour.The principal stress vector diagram after the excavation of the No. 16 model is in Figure 8, showing the pressure arch shape and the ring-shaped distribution of the maximum principal stress.H represents the vertical distance from the tunnel vault to the ground.The model's boundary conditions were set so that the upper boundary was a free surface boundary, and both sides and the bottom edges were constrained by the method of normal displacement.The initial in situ stress was the self-weight stress of the rock mass.

Method of Determining the Loose Zone in the Surrounding Rock Mass
The pressure arch in the surrounding rock mass was characterized by an arch structure.The pressure arch body bears its own load and the load of the rock mass above the arch so that the surrounding rock maintains a certain degree of self-stability after excavation.However, the rock mass below the pressure arch has poor stability because it is close to the tunnel excavation outline, which leads to the stress exceeding the rock mass strength and causing damage.Therefore, the loose zone in the surrounding rock mass is defined as the rock mass from the inner boundary of the pressure arch to the tunnel contour.The principal stress vector diagram after the excavation of the No. 16 model is in Figure 8, showing the pressure arch shape and the ring-shaped distribution of the maximum principal stress.Based on existing research results [34], this paper took the connection of the peak point of the surrounding rock's tangential stress after the excavation of the tunnel as the judgment basis for the pressure arch inner boundary.The inner boundary criterion was embedded into FLAC3D with the Fish programming language to conduct numerical simulations.The correlation of the tangential stress can be calculated as follows: The simulation was processed as follows: ① The model was established according to the orthogonal test scheme and the initial stress field was calculated.② The element stress component was extracted after the calculation was balanced and the tangential stress of each component was calculated using the correlation.③ Values of the tangential Based on existing research results [34], this paper took the connection of the peak point of the surrounding rock's tangential stress after the excavation of the tunnel as the judgment basis for the pressure arch inner boundary.The inner boundary criterion was embedded into FLAC3D with the Fish programming language to conduct numerical simulations.The correlation of the tangential stress can be calculated as follows: The simulation was processed as follows: 1 The model was established according to the orthogonal test scheme and the initial stress field was calculated. 2 The element stress component was extracted after the calculation was balanced and the tangential stress of each component was calculated using the correlation. 3Values of the tangential stress of each element on the search path and were compared, and the element with the peak tangential stress point was marked.Owing to the model's symmetry, half of the model was selected for the search.The specific shape of the loose zone was obtained by connecting the inner boundary position of the pressure arch on each search path. 4The loose zone was smoothed by the least square method to obtain the final scope of the loose zone in the surrounding rock mass, as shown in Figure 9.It should be noted that the influence of floor heave was not considered in this study.
Based on existing research results [34], this paper took the connection of the peak point of the surrounding rock's tangential stress after the excavation of the tunnel as the judgment basis for the pressure arch inner boundary.The inner boundary criterion was embedded into FLAC3D with the Fish programming language to conduct numerical simulations.The correlation of the tangential stress can be calculated as follows: The simulation was processed as follows: ① The model was established according to the orthogonal test scheme and the initial stress field was calculated.② The element stress component was extracted after the calculation was balanced and the tangential stress of each component was calculated using the correlation.③ Values of the tangential stress of each element on the search path and were compared, and the element with the peak tangential stress point was marked.Owing to the model's symmetry, half of the model was selected for the search.The specific shape of the loose zone was obtained by connecting the inner boundary position of the pressure arch on each search path.④ The loose zone was smoothed by the least square method to obtain the final scope of the loose zone in the surrounding rock mass, as shown in Figure 9.It should be noted that the influence of floor heave was not considered in this study.

Analysis of the Factors Influencing the Loose Zone in the Surrounding Rock Mass
Three geometric quantities were selected as evaluation indexes reflecting the scope of the closed annular loose zone, namely  (the span of the loose zone at the tunnel vault),   (the height of the loose zone at the tunnel vault), and  (the area of the loose zone).The geometric quantity diagram of the loose zone is shown in Figure 10.The data were processed for the 16 groups of working conditions, and the test number of the closed annular loose zone and the evaluation indexes of the loose zone scope were then summarized.The evaluation indexes of the loose area are shown in Table 5.

Analysis of the Factors Influencing the Loose Zone in the Surrounding Rock Mass
Three geometric quantities were selected as evaluation indexes reflecting the scope of the closed annular loose zone, namely L (the span of the loose zone at the tunnel vault), H 1 (the height of the loose zone at the tunnel vault), and S (the area of the loose zone).The geometric quantity diagram of the loose zone is shown in Figure 10.The data were processed for the 16 groups of working conditions, and the test number of the closed annular loose zone and the evaluation indexes of the loose zone scope were then summarized.The evaluation indexes of the loose area are shown in Table 5.The CCR model quantitatively evaluates the factors influencing the loose zone in the surrounding rock mass using a fuzzy restriction on the weight [35].This model can realize the evaluation and prediction of multi-index input and multi-index output problems with the help of the linear programming method.The model overcomes the shortcomings of the traditional CCR model, which discards influencing factors (the weight value is 0) when solving for optimal efficiency.The maximum fuzzy membership function constrains the upper and lower limits of the weight values to obtain the model's optimal weight solution.The CCR model with a fuzzy restriction on the weight can be expressed as: where µ is the maximum fuzzy membership function; y rj is the output value of the r-th term; x ij is the input value of the i-th term; µ r is the weight value of the r-th output; v i is the weight value of the i-th input; u U r and v U i are the upper limit values of the output weight and input weight, respectively.
When analyzing the factors influencing the loose pressure in the surrounding rock mass of a super-span tunnel, the input factors are H, B, λ, K, and GSI.The values of each impact factor are shown in Table 2.The output factors are L, H 1 , and S, and their geometric statistics are shown in Table 6.The weight values of the influencing factors are summarized in Table 1.It can be concluded that the order of the weight values of the influencing factors in the loose zone of a super-large section tunnel is B > GSI > K > λ > H, which indicates that two factors (B for the tunnel span and GSI for surrounding rock quality) have a high degree of influence on the surrounding rock loose zone, while the other factors have less of an impact.

Boundary Function of the Surrounding Rock's Loose Zone
The boundary curve of the loose zone was fitted by a non-linear function based on the scope of the loose zone in the surrounding rock mass.Owing to the symmetry of the loose zone in the surrounding rock, half of the boundary curve was taken as the research object.The Cartesian rectangular coordinate system was adopted, taking the arch vertex of the loose zone as the coordinate origin.The positive direction on the X-axis was vertically downward, and the positive direction on the Y-axis was horizontal and to the left.Ten characteristic points were selected from the coordinate information on the boundary curve and were fitted by the built-in function in MATLAB according to the changing trend of data points to determine the boundary function of the loose zone.After data sorting and analysis, the quadratic function had the highest degree of fit for the dataset.Therefore, this paper took it as the boundary function with a, b, and c being undetermined coefficients.Taking Model 7 as an example, the fitting effect and relevant information are represented in Figure 11 and Table 7.
The boundary curve of the loose zone was fitted by a non-linear function based on the scope of the loose zone in the surrounding rock mass.Owing to the symmetry of the loose zone in the surrounding rock, half of the boundary curve was taken as the research object.The Cartesian rectangular coordinate system was adopted, taking the arch vertex of the loose zone as the coordinate origin.The positive direction on the X-axis was vertically downward, and the positive direction on the Y-axis was horizontal and to the left.Ten characteristic points were selected from the coordinate information on the boundary curve and were fitted by the built-in function in MATLAB according to the changing trend of data points to determine the boundary function of the loose zone.After data sorting and analysis, the quadratic function had the highest degree of fit for the dataset.Therefore, this paper took it as the boundary function with , , and  being undetermined coefficients.Taking Model 7 as an example, the fitting effect and relevant information are represented in Figure 11 and Table 7.The following can be concluded from the correlation between the evaluation indexes of the loose zone and the boundary function: Taking Correlation 5 into the boundary function  =  +  + , the undetermined coefficients , , and  can be approximated:  The following can be concluded from the correlation between the evaluation indexes of the loose zone and the boundary function: Taking Correlation 5 into the boundary function y = ax 2 + bx + c, the undetermined coefficients a, b, and c can be approximated: Using the data in Tables 3-5, the multiple regression equations of L and H 1 can be concluded via the process described below: (1) Based on the principle of the least squares method, the univariate regression equation between the influencing factors of the loose zone with L and H 1 are established, as shown in Table 8.According to the manual, the correlation coefficient threshold value is 0.83.As illustrated in Table 8, the correlation coefficients of the five influencing factors all exceed 0.83, indicating that the univariate regression equations of L and H 1 are correct.
The correctness of Correlation 7 and Correlation 8 are verified by the variance analysis, as shown in Table 9.The F distribution value is F(6, 4) = 6.59 in Correlation 7 and F(6, 4) = 7.02 in Correlation 8.Both are greater than F 0.05 (6, 4) = 6.16 at the 0.05 significance level, indicating that the accuracy of the regression equation meets the requirements.
The parameter values must be defined within the domain.

Upper Bound Solution for
Loose Pressure in the Surrounding Rock Mass 5.1.Upper Boundary Solution for Loose Pressure for Extra-Large Section Tunnels 5.1.1.Basic Assumptions Four basic assumptions must be made when the upper bound theory is applied to analyze and calculate the loose pressure of surrounding rock: 1 The rock-soil mass is regarded as a perfectly elastic-plastic material that obeys the H-B yield criteria. 2 The tunnel arch is under a vertical even load q, and the side walls are under a horizontal even load e.Assuming that q is proportional to e, the proportional coefficient is the side pressure coefficient K, that is, e = Kq 3 The influences of the construction process and floor heave are not considered.

Description of Failure Mechanism
The basic shape of the loose zone obtained from the orthogonal test shows that the failure forms of the deeply buried tunnel include the tunnel vault subsidence and side wall shrinkage.Based on the determination of the failure zone, the mass of rock and soil in the failure mechanism is divided into rigid blocks that meet the conditions.The velocity and direction of each rigid block are obtained according to the geometric relationship and the sine theorem corresponding to the failure mechanism.
In summary, the failure mechanism of a deeply buried tunnel is represented in Figure 12.The tunnel section is approximately rectangular, ABCD, with height h and span B. M and N are located at the midpoint of the tunnel vault and the midpoint of the loose zone boundary curve, respectively.O is located on the straight line MN and is the auxiliary point.H represents the burial depth of the cavity H 1 represents the distance from O to the ground plane, and the support reaction forces of the tunnel arch and the side walls are the uniform loads q and e, respectively.
Owing to the symmetry of the failure mechanism, this study takes the left half of the failure mechanism for the corresponding description.AB is extended to intersect the boundary curve at point A m .∠A m ON is divided into m + 1 parts with angles for β 0 , β 1 , . . ., β m .∠AA m C is divided into n parts with angles β m+1 , β m+2 , . . ., β m+n .The intersections of each angle and the boundary of the loose zone are recorded as A 0 , A 1 , . . ., A m , . . ., A m+n+1 .Connecting N A 0 , A 0 A 1 , A 1 A 2 , . . ., A m A m+1 , . . ., A m+n−1 C, when m and n are large enough, the error between each connecting line and the boundary curve is therefore small enough.
GSI ≤ 40.The parameter values must be defined within the domain.

Basic Assumptions
Four basic assumptions must be made when the upper bound theory is applied to analyze and calculate the loose pressure of surrounding rock: ① The rock-soil mass is regarded as a perfectly elastic-plastic material that obeys the H-B yield criteria.
② The tunnel arch is under a vertical even load , and the side walls are under a horizontal even load .Assuming that  is proportional to , the proportional coefficient is the side pressure coefficient , that is,  =  ③ The influences of the construction process and floor heave are not considered.

Description of Failure Mechanism
The basic shape of the loose zone obtained from the orthogonal test shows that the failure forms of the deeply buried tunnel include the tunnel vault subsidence and side wall shrinkage.Based on the determination of the failure zone, the mass of rock and soil in the failure mechanism is divided into rigid blocks that meet the conditions.The velocity and direction of each rigid block are obtained according to the geometric relationship and the sine theorem corresponding to the failure mechanism.
In summary, the failure mechanism of a deeply buried tunnel is represented in Figure 12.The tunnel section is approximately rectangular, , with height ℎ and span . and  are located at the midpoint of the tunnel vault and the midpoint of the loose zone boundary curve, respectively. is located on the straight line  and is the auxiliary point. represents the burial depth of the cavity  represents the distance from  to the ground plane, and the support reaction forces of the tunnel arch and the side walls are the uniform loads  and , respectively.According to the associated flow rule, the velocity vector sum of the rigid slider is zero, and the angle between the relative velocity vector and the discontinuity line should be ϕ .Figure 13 shows the kinematically admissible velocity field corresponding to the failure mechanism.The velocity of block ∆A 0 ON is V 0 in the vertically downward direction.The velocity of block ∆AOM is V m+n+1 in the vertically downward direction.The velocities of the other blocks, from top to bottom, are V 1 , . . .V m+n , and the relative velocities of each rigid slider on the discontinuity line are V 0,1 , . . ., V m+n−1,m+n , V m+n,m+n+1 .The shape of the blocks is determined by α 1 , . . ., α m and β 1 , . . ., β m .With reference to Protodyakonov's theory and the Terzaghi theory, the shapes can be assumed to be β 0 = π/4 − ϕ'/4.Based on the assumptions mentioned above, the loose zone is regarded as an ideal elastic-plastic body which obeys the H-B criterion.The relation between  and  can be obtained as below using the tangent method (Figure 14) [36]:  Based on the assumptions mentioned above, the loose zone is regarded as an ideal elastic-plastic body which obeys the H-B criterion.The relation between c t and ϕ t can be obtained as below using the tangent method (Figure 14) [36]:  Based on the assumptions mentioned above, the loose zone is regarded as an ideal elastic-plastic body which obeys the H-B criterion.The relation between  and  can be obtained as below using the tangent method (Figure 14) [36]:  The non-linear shear strength index ϕ t in the above correlations serves as an unknown parameter in the upper-bound limit analysis, reflecting the instantaneous variability of the tangent, which is obtained by the optimization algorithm.c t is obtained by Correlation 10 after ϕ t is determined.
(2) Calculation of geometric quantities Owing to the symmetry of the failure mechanism, this paper takes the left part of MN as an example for geometric quantity calculation.The point N is taken as the coordinate origin and a rectangular coordinate system is established.The vertical is the X-axis, where the down is positive, and the horizontal direction is the Y-axis, where the left is positive.According to the loose zone boundary equation, the expression of the curve NC is y = ax 2 + bx + c.
According to the relevant description of failure mechanism and the research results of the loose zone in the surrounding rock mass, it can be known that |AM| obtained in which the maximum value is the upper bound solution for the loose pressure in the surrounding rock mass.
The horizontal support reaction can be calculated as: 5.2.Upper Boundary Solution for Loose Pressure for Extra-Large Section Tunnels Based on the above research results, the process of calculating the loose pressure for a super-span tunnel is represented in Figure 15.The loose pressure of the large-span section (the maximum Section, Gradient Section 1, and Gradient Section 2) was calculated according to the calculation process outlined above, with the geological conditions of the tunnel taken into consideration.The calculation parameters and the calculation results are represented in Tables 10 and 11, respectively.
With the commands of the MATLAB optimization toolbox utilized and the constraints met, a certain number of the loose pressure values in the surrounding rock mass are obtained in which the maximum value is the upper bound solution for the loose pressure in the surrounding rock mass.
The horizontal support reaction can be calculated as:

Upper Boundary Solution for Loose Pressure for Extra-Large Section Tunnels
Based on the above research results, the process of calculating the loose pressure for a super-span tunnel is represented in Figure 15.The loose pressure of the large-span section (the maximum Section, Gradient Section 1, and Gradient Section 2) was calculated according to the calculation process outlined above, with the geological conditions of the tunnel taken into consideration.The calculation parameters and the calculation results are represented in Table 10 and Table 11, respectively.

Comparison between the Correlation Results and In Situ Measured Values
Several difficulties and risks existed in the construction of the super-large section.In order to ensure the construction quality and grasp the distribution pattern of the surrounding rock pressure, our research group conducted a half-year on-site monitoring of the large-span section from April 2018 to October 2018.The monitoring items included the displacement of the deep surrounding rock, the contact pressure between the surrounding rock and the structure's steel, and the stress of steel arch.The Maximum Section, the Gradient Section 1 and the Gradient Section 2 were selected as surrounding rock pressure monitoring sections.The site photos of the earth pressure cell are represented in Figure 16.
rounding rock pressure, our research group conducted a half-year on-site monitoring of the large-span section from April 2018 to October 2018.The monitoring items included the displacement of the deep surrounding rock, the contact pressure between the surrounding rock and the structure's steel, and the stress of steel arch.The Maximum Section, the Gradient Section 1 and the Gradient Section 2 were selected as surrounding rock pressure monitoring sections.The site photos of the earth pressure cell are represented in The stable values of the loose pressure of each monitoring section's surrounding rock were collected to determine the circumferential distribution.They were then compared with the value calculated for the loose pressure in the surrounding rock mass based on the upper bound limit analysis method, as shown in Figure 17.The figure shows that the distribution law of the surrounding rock's calculated loose pressure is consistent with the measured results except for the right side of the maximum section.The main reason for monitoring data distortion at the right side of the maximum section was the failure to install pressure sensors in a timely manner after excavation.As the surrounding rock pressure is monitored after the completion of the supporting structure, the monitoring results are the variation relative to the initial stress, and the monitoring values are usually less than the theoretical calculation values.As illustrated in Figure 17, the calculated results have a good envelope effect on the measured data, which proves the correctness and practicability of the calculation method.The stable values of the loose pressure of each monitoring section's surrounding rock were collected to determine the circumferential distribution.They were then compared with the value calculated for the loose pressure in the surrounding rock mass based on the upper bound limit analysis method, as shown in Figure 17.The figure shows that the distribution law of the surrounding rock's calculated loose pressure is consistent with the measured results except for the right side of the maximum section.The main reason for monitoring data distortion at the right side of the maximum section was the failure to install pressure sensors in a timely manner after excavation.As the surrounding rock pressure is monitored after the completion of the supporting structure, the monitoring results are the variation relative to the initial stress, and the monitoring values are usually less than the theoretical calculation values.As illustrated in Figure 17, the calculated results have a good envelope effect on the measured data, which proves the correctness and practicability of the calculation method.

Conclusions
(1) With the combination of a numerical analysis and an orthogonal test, the factors affecting the loose zone range of the surrounding rock were studied.According to CCR model with fuzzy restrictions on the weight in statistics, the order of factors influencing the weight value from large to small is  > GSI >  >  > .
(2) The loose boundary zone expression of a deeply buried tunnel was determined to be  =  +  + , and a calculation method for the undetermined coefficients of the boundary expression was proposed.The correlations of typical geometric quantities  and   , representing the loose zone range, were established based on the unitary regression equations between the factors affecting the loose zone and typical geometric quantities.
(3) Based on the boundary function of the loose zone and the upper bound limit analysis method, a limit analysis model of the loose pressure in the surrounding rock mass under a non-linear failure criterion was established via the tangent method and virtual power principle, and the correlations of the loose pressure were deduced.The calculated values of the loose pressure were consistent with the pressure-monitoring results in the different sections of the large-span section, verifying the correctness of the proposed calculation method.The research results can provide new ideas for the selection of a supporting structure for similar super-section tunnel projects.

Conclusions
(1) With the combination of a numerical analysis and an orthogonal test, the factors affecting the loose zone range of the surrounding rock were studied.According to CCR model with fuzzy restrictions on the weight in statistics, the order of factors influencing the weight value from large to small is B > GSI > K > λ > H.
(2) The loose boundary zone expression of a deeply buried tunnel was determined to be y = ax 2 + bx + c, and a calculation method for the undetermined coefficients of the boundary expression was proposed.The correlations of typical geometric quantities L and H 1 , representing the loose zone range, were established based on the unitary regression equations between the factors affecting the loose zone and typical geometric quantities.
(3) Based on the boundary function of the loose zone and the upper bound limit analysis method, a limit analysis model of the loose pressure in the surrounding rock mass under a non-linear failure criterion was established via the tangent method and virtual power principle, and the correlations of the loose pressure were deduced.The calculated values of the loose pressure were consistent with the pressure-monitoring results in the different sections of the large-span section, verifying the correctness of the proposed

Figure 2 .
Figure 2. Route planning of the Greater Bay Area.

Figure 3 .
Figure 3. Plan view of bifurcation section.

Figure 2 .
Figure 2. Route planning of the Greater Bay Area.

Figure 2 .
Figure 2. Route planning of the Greater Bay Area.

Figure 3 .
Figure 3. Plan view of bifurcation section.Figure 3. Plan view of bifurcation section.

Figure 3 .
Figure 3. Plan view of bifurcation section.Figure 3. Plan view of bifurcation section.

Figure 4 .Figure 5 .
Figure 4. Engineering geologic profile of the project.

Figure 4 .
Figure 4. Engineering geologic profile of the project.

Figure 4 .Figure 5 .
Figure 4. Engineering geologic profile of the project.

Figure 5 .
Figure 5. Site photos of the drill core.(a) Small clear-distance section; (b) large-span section.

Figure 9 .
Figure 9. Search path and loose zone boundary.

Figure 9 .
Figure 9. Search path and loose zone boundary.

Figure 10 .
Figure 10.Loose zone of surrounding rock mass.

Figure 10 .
Figure 10.Loose zone of surrounding rock mass.

Figure 11 .
Figure 11.The fitting curve of the surrounding rock's loose zone.

Figure 11 .
Figure 11.The fitting curve of the surrounding rock's loose zone.

Figure 12 .Figure 12 .
Figure 12.Failure mechanisms.Owing to the symmetry of the failure mechanism, this study takes the left half of the failure mechanism for the corresponding description. is extended to intersect the boundary curve at point  .∠  is divided into  + 1 parts with angles for  ,  , … ,  .∠  is divided into n parts with angles  ,  , … ,  .The Figure 12.Failure mechanisms.

Figure 13 .
Figure 13 shows the kinematically admissible velocity field corresponding to the failure mechanism.The velocity of block ∆  is  in the vertically downward direction.The velocity of block ∆ is  in the vertically downward direction.The velocities of the other blocks, from top to bottom, are  , …  , and the relative velocities of each rigid slider on the discontinuity line are  , , … ,

Figure 13 .
Figure 13.The velocity field schematic diagram.5.1.3.Derivation of the Upper Bound Solution for the Loose Pressure (1) Limit Analysis Based on H-B Strength CriterionBased on the assumptions mentioned above, the loose zone is regarded as an ideal elastic-plastic body which obeys the H-B criterion.The relation between  and  can be obtained as below using the tangent method (Figure14)[36]:

Figure 13 .
Figure 13.The velocity field schematic diagram.5.1.3.Derivation of the Upper Bound Solution for the Loose Pressure (1) Limit Analysis Based on H-B Strength CriterionBased on the assumptions mentioned above, the loose zone is regarded as an ideal elastic-plastic body which obeys the H-B criterion.The relation between c t and ϕ t can be obtained as below using the tangent method (Figure14)[36]:

Figure 13 .
Figure 13 shows the kinematically admissible velocity field corresponding to the failure mechanism.The velocity of block ∆  is  in the vertically downward direction.The velocity of block ∆ is  in the vertically downward direction.The velocities of the other blocks, from top to bottom, are  , …  , and the relative velocities of each rigid slider on the discontinuity line are  , , … ,

Figure 13 .
Figure 13.The velocity field schematic diagram.5.1.3.Derivation of the Upper Bound Solution for the Loose Pressure (1) Limit Analysis Based on H-B Strength CriterionBased on the assumptions mentioned above, the loose zone is regarded as an ideal elastic-plastic body which obeys the H-B criterion.The relation between  and  can be obtained as below using the tangent method (Figure14)[36]:

Figure 15 .
Figure 15.flow of calculation for the loose pressure in the surrounding rock mass.

Figure 15 .
Figure 15.Flow of calculation for the loose pressure in the surrounding rock mass.

22 Figure 17 .
Figure 17.Comparison between the calculated and measured values.

Figure 17 .
Figure 17.Comparison between the calculated and measured values.

Table 1 .
Different correlations for calculation methods.

Table 2 .
Characteristics of the surrounding rock of the project.

Table 2 .
Characteristics of the surrounding rock of the project.

Table 4 .
Mechanical parameters of the rock mass corresponding to different GSI values.× 300 m (210 m × 300 m, 240 m × 300 m, 270 m × 300 m).The element mesh and model dimension are shown in Figure 7.
The numerical calculation model of the unit thickness was set as: height × width = 180 m

Table 5 .
Summary of evaluation indexes of the loose area.

Table 5 .
Summary of evaluation indexes of the loose area.

Table 6 .
Weight value of influencing factors.

Table 7 .
Results for the fitting curve.

Table 7 .
Results for the fitting curve.

Table 9 .
Variance analysis of multiple regression equation.

Table 11 .
Loose pressure in the surrounding rock mass based on the upper bound analysis.