Approximation method for yielding support analysis in high ground stress soft surrounding rock tunnels

In solving the whole process of interaction between soft rock and yielding support in high-stress environments in tunnels using mechanical analysis methods, it is challenging to simultaneously satisfy both displacement coordination and static equilibrium at the contact surface between the rock and the support structure. This paper, based on the mechanical analysis of rock and rigid support, considers the impact of the circumferential installation of yielding elements on radial displacement, and proposes displacement approximation and support force approximation methods using displacement coordination and static equilibrium as approximation conditions. The study fits curves of numerical simulation results and laboratory test results of yielding elements, and attempts to directly use the laboratory test data set of yielding elements as computational data. By calculating two circular tunnel examples and comparing the effects of the trisection method, bisection method, and substitution method on the convergence of the displacement approximation method, the effectiveness of the methods proposed in this paper is verified. The research results show that the two approximation algorithms proposed in this paper have good accuracy and reliability in calculating the relative displacement of rock and yielding support structure contact surfaces, and the support force of yielding support. The bisection method outperforms the trisection and substitution methods in terms of stability and convergence. However, there are certain limitations in this study, such as the effectiveness of the algorithm may be influenced by geological conditions; the complexity of actual geological conditions may exceed the assumptions of the current rock-support mechanical analysis model.


Introduction
Using yielding support instead of rigid support is an effective method to solve the problem of large deformation in soft rock tunnels under high ground stress conditions [1].Engineering practice has shown that increasing the stiffness of support to resist large deformations of soft rock by means of strong support and hard top often fails to meet the support needs and frequently requires multiple replacements or reinforcement of the support structure [2].This not only increases the project budget and extends the construction period but also puts the support structure under high stress, increasing the risk of failure [3].Yielding support is an improvement based on rigid support, such as installing yielding elements around the circumference of the initial rigid support [4].Yielding support allows for deformation of the support to release the energy of the surrounding rock and reduce the pressure of the surrounding rock, thereby improving the safety of the support [5], and is widely used in tunnels with high ground stress and weak surrounding rock.With the widespread application of yielding support in tunnels with high ground stress and weak surrounding rock, it is of great engineering significance to accurately consider the interaction mechanism between the surrounding rock and the yielding support system.
To reveal the interaction mechanism between the surrounding rock and the yielding support, scholars at home and abroad have conducted in-depth studies on this mechanism using methods such as numerical simulation and mechanical analysis.For example, numerical simulation methods have been used by researchers like Lian Chuanjie [6][7][8] to analyze the design parameters, simulation methods, and support effects of prestressed yielding anchor bolts; Wang Bo and others [9] constructed a yielding support system consisting of yielding anchor bolts, compressible steel arches, and deformation slots or yielding controllers to study the effects of the initial yielding point and yielding amount on the release of strain energy of the surrounding rock and support bearing; Tian and others [10] established a combined yielding support system of U-shaped steel frames and compressible foam concrete, analyzing the relationship between the thickness of the compressible layer and large deformations of soft rock; Qi Chun and others [11] used the Mohr-Coulomb criterion to establish a beam-spring model to simulate the contact effect between segmental lining and compressible layers, and a stiffness reduction model for compressible layers to analyze the effect of the presence or absence of compressible layers, their stiffness, and thickness on the internal force and deformation distribution of the support; Liu Fang and others [12] simulated the surrounding rock with the Cvisc creep model, using particle elements, solid elements, cable elements, and beam elements to simulate gravel filling layers, surrounding rock, grouting layers, anchor bolts, and linings, establishing a combined support model of segmental lining-compressible layer-anchor bolts.Hu and others [13] compared the yielding performance of two compressible materials (compressible concrete and expandable clay), illustrating the interaction mechanism between segmental lining and yield materials.However, the calculation accuracy of the numerical simulation method is highly related to the modeling quality of the yielding elements, susceptible to human factors, and the universality and predictive ability of numerical models need to be improved.
In the analysis of support displacement, yielding amount is added for calculation.Lei Sheng xiang and others [14] used the structural mechanics method to analyze the deformation pattern of the support structure and gave suggestions on the installation position of the circumferential yielding elements; Tian Yun and others [15] used the circular tunnel perimeter elastoplastic displacement formula, considering the aging weakening characteristics of the surrounding rock and the resistance of the buffer layer support, to propose a buffer layer parameter design method; Wu Kui and others [16][17][18][19] used an improved fractional order Burgers creep model to derive an analytical solution for tunnel displacement and support pressure considering the effects of the tunnel face and delayed support action; Dong Jianhua and others [20] derived the full process analytical solution for elastoplastic deformation of the surrounding rock under the action of yielding support.Hu Xiongyu and others [21] established a theoretical model of the interaction between the surrounding rock, the compressible layer of expanded clay, and the segmental lining, and determined the calculation method for the thickness design parameters of the expanded clay compressible layer.Yang Kai and others [22] calculated the radial shrinkage caused by the circumferential yielding amount in circular tunnels, adding this radial shrinkage to the radial support displacement and then calculating the new support force and surrounding rock displacement.However, the displacement of the surrounding rock around the tunnel calculated in a single calculation by superimposing the yielding amount is not equal to the displacement of the top of the yielding support, and the result only satisfies the static equilibrium condition, not the displacement coordination condition.
Therefore, to address the limitations of existing methods, this paper proposes a new approximation method for solving.This study is based on the interaction mechanism between rock and support, focusing on the tunnel steel arch yielding elements proposed by Dong Biao [23] and others as the research subject.In this study, we attempt to directly use the laboratory test data set of yielding elements as computational data and apply numerical approximation algorithms to calculate the relative displacement at the contact surface between the rock and yielding support, as well as the support force of the support structure.The core of the numerical approximation algorithm proposed in this paper lies in accurately calculating the relative displacement at the contact surface between the rock and yielding support and the support force of the support structure, analyzing the entire interaction process between rock and yielding support.
Additionally, through comparative analysis of numerical simulation results and laboratory test results of yielding elements, combined with numerical approximation algorithms, this study highlights the necessity of using the proposed approximation method in solving yielding support issues.This method not only satisfies the static equilibrium conditions but also meets the displacement coordination conditions, providing a solution closer to practical engineering needs.Finally, by analyzing two circular tunnel examples, the feasibility of the numerical approximation algorithm in solving yielding problems is further verified.

The surrounding rock-support interaction model
In analyzing the interaction mechanisms of elastic surrounding rock-rigid support and elastoplastic surrounding rock-rigid support, it is assumed that the initial ground stress of the rock mass is P 0 and the support resistance is P The tunnel excavation model is shown in

Support elastic solution
The displacement and stress of rigid support under uniform pressure P are calculated using the Lame ´solution, as follows: where σ ρ , σ φ , and u ρ represent radial stress, circumferential stress, and radial displacement, respectively.

Elastic solution for surrounding rock
The relative displacement of the elastic surrounding rock [24,25] can be expressed as When ρ = R, the relative displacement around the rock tunnel is The radial displacement of the support's top is When the displacement coordination condition is satisfied at the contact surface between the surrounding rock and the support, u R ¼ u 0 R the support force can be calculated using Eqs (5) and ( 6) If the displacement u 0 R around the rock tunnel is known, the support force can be calculated using Eq (5), which is

Elastic solution for surrounding rock
When not considering the volume expansion equation of elastoplastic bodies, the formula for calculating the elastoplastic displacement of general circular tunnels [23] is: where R p is the radius of the plastic zone, G is the shear modulus, C and θ are the cohesion and internal friction angle of the rock, and φ is the position angle.

Support yielding displacement
The displacement coordination relationship mentioned in this paper refers to the consistent relative radial displacement around the tunnel between the surrounding rock and the support.
After installing yielding elements circumferentially on the support, it is necessary to convert the circumferential yielding displacement of the support into radial yielding displacement based on the geometric relationship inside the tunnel.Before determining the circumferential yielding displacement, the concentrated force on the yielding elements should be calculated first, then search the circumferential yielding displacement u i corresponding to this concentrated force from the loading test data set.The formula for calculating the concentrated force is as follows: where A is the contact area between the yielding element and the rigid support.
Assuming that the circumferential yielding displacement is uniformly distributed along the midline of the circular ring, and the number of yielding elements is n the total radial yielding displacement u n is calculated based on the geometric relationship as By adding the total radial yielding displacement u n to the radial displacement u R of the rigid support, the total radial displacement of the yielding support can be obtained 3 The proximity algorithm Using Formula (10), the concentrated force on the pressure components is calculated from the rigid support ring's stress.This force is then used to search for experimental data on pressure components to obtain the displacement due to pressure.However, this method of obtaining displacement due to pressure does not take into account the impact of pressure action on the circumferential stress of the support ring.Compared to the actual displacement due to pressure on the support, its value is larger.The relative displacement of the pressure support calculated using Formula (12) at the interface between the support and the surrounding rock will struggle to meet the displacement coordination conditions.Furthermore, if only the displacement coordination condition between the surrounding rock and the pressure support is considered, without considering the static equilibrium condition, the calculated displacement value and support force will be inaccurate.To address the difficulty of simultaneously meeting the displacement coordination and static equilibrium conditions at the contact surface, this paper proposes the displacement proximity method and the support force proximity method.The numerical approximation algorithm technology roadmap is shown in Fig 3.

The displacement proximity method
The displacement approximation method involves using the bisection method to approximate the displacement around the surrounding rock tunnel to the displacement of the support top plate.The displacement of the surrounding rock and the support force are calculated and updated until the displacement coordination condition (convergence) is met.The specific steps are as follows: 1) Calculate the support force P using Formula (7); 2) Let ρ = R, calculate the circumferential stress σ φ,R and the radial displacement u R of the support using Formulas (2) and (3); 3) Substitute σ φ,R into Formula (10) to calculate the concentrated force F, and search for the corresponding displacement u i from the displacement-load data of the pressure components; 4) Substitute u i into Formulas (11) and (12) to calculate the total radial displacement of the support u * R ; 5) Calculate the radial displacement around the rock tunnel u 0 R using Formula (5); 6) Determine whether u 0 R equals u * R 6.1) If yes, the solution is complete; , calculate the support force P using Formula (8), and return to step 2) for further iteration.

Support force approximation method
For situations where it is inappropriate to pre-calculate the support force using the displacement coordination condition, this paper proposes a support force approximation algorithm.This algorithm requires an initial assumption of the support force range, suggested to be [0, 10P 0 ].Based on the positional relationship between the displacement around the rock tunnel and the radial displacement of the support, the support force range is gradually narrowed by comparing whether the displacement calculated with the preset support force meets the displacement coordination condition.The specific steps are as follows: 1) Assume P down = 0, P up = 10P 0 and preset the support force P = (P down +P up )/2; 2) Let ρ = R, calculate the circumferential stress σ φ,R and radial displacement u R of the support using Formulas (2) and (3); 3) Substitute σ φ,R into Formula (10) to calculate the concentrated force F, and search for the corresponding displacement u i from the displacement-load data of the pressure components; 4) Substitute u i into Formulas (11) and (12) to calculate the total radial displacement of the support u * R ; 5) Calculate the radial displacement around the rock tunnel u 0 R using Formula (5); 6) Determine whether u 0 R equals u * R : 6.1)If yes (u 0 R ¼ u * R ), the solution is complete; 6.2) If not, compare the sizes of u 0 R and u * R , update parameters, and return to step 1) for further iteration;

Construction of yielding elements
The yielding element is composed of an inner channel steel, an outer channel steel, and highstrength bolts, where the top of the inner channel steel and the bottom of the outer channel steel are connected to Q235 I-beam steel, respectively.The inner and outer channel steels are connected by GT/B 1228M16×60 high-strength bolts.The cross-sections of the inner channel steel, outer channel steel, and slide channel are shown in Figs 4-6, respectively, with detailed parameters of the channel steel and slide channel surface presented in Table 1 [23].

Radial force effect
A horizontal load of 10kN is applied to the top of the I-beam of the yielding component, and the bolt pre-tension is set to 60KN, resulting in a radial Mises stress contour map of the

Axial force effect
A vertical load of 200KN is applied to the top of the I-beam of the yielding component.The bolt pre-tension is the same as during the radial force effect, resulting in an axial Mises stress contour map of the yielding component (Fig 9).As can be seen from Fig 9, under this load, the slot displacement of the yielding component has been fully released.
The load-displacement curve of the yielding components is plotted ( Fig 10).As shown in Fig 10, the deformation process of the yielding components is divided into three stages: the elastic rise stage of bearing capacity (Ⅰ), during which the bearing capacity gradually increases from 0 to about 60KN with the increase of bolt friction; the constant bearing capacity stage (Ⅱ), where the bearing capacity is maintained within the range of 60kN to 62kN, and the bearing capacity remains stable as the slot smoothly slides; the compaction closing stage (Ⅲ), where the slot displacement of the yielding component has been fully released, the inner and outer slot steels are pressed together and closely combined, and the final bearing capacity stabilizes at 172KN.
The numerical simulation results show that once the yielding component enters the plastic state, the friction between the bolt and the slot steel can provide a constant bearing capacity.When the slot displacement of the yielding component is fully released, the compressed inner and outer slot steels provide a high bearing capacity, and the yielding component can still bear the load.
6 Laboratory test of yielding components

Overview of the test
A microcomputer-controlled electro-hydraulic servo press (YAW-5000F) is used to load the yielding components.The loading method is program-controlled, with a load size of 200 kN and a speed of 1 kN/s.The sample diagram of the yielding component is shown in Fig 11 [23].

Test results
As can be seen from Fig 12, local failure occurred at the slot node, damage occurred on the contact surface between the bolts and the inner and outer slot steels, and the stiffening plate of the inner slot steel bent.After loading, the failure form of the yielding components is shown in Fig 12 [23].https://doi.org/10.1371/journal.pone.0299426.g007

Comparison of results
The load-displacement curve of the yielding components is plotted and compared with the curve from the numerical simulation results (Fig 13).The comparison shows that the numerical simulation results of the yielding components are highly consistent with the experimental results, but there is a certain numerical difference between them.The main reason for this difference is the manufacturing errors in the slots and bolt holes of the yielding components, such as different roughness of different slot sections.
In the numerical simulation, the deformation of the yielding component when entering the constant bearing capacity stage is 5mm, while in the laboratory test, the deformation is 12mm; however, both stabilize at a load of 60KN when entering the constant bearing capacity stage.After yielding, the numerical simulation and experimental results show maximum load values of 172kN and 169kN, and maximum displacement values of 46mm and 56mm, respectively.

Case analysis
This paper compares the numerical simulation results of yielding elements with the results of laboratory experiments and uses the interaction models between elastic yielding support and elastic surrounding rock, as well as between elastic yielding support and elasto-plastic surrounding rock, as examples to analyze the entire interaction process between the surrounding rock and yielding support.The model parameters are shown in Table 2.

Results of the displacement approximation method
The relative displacement values of elastic surrounding rock-rigid support interaction in a stable state are taken as the initial values for iterative calculation (0th iteration).Subsequently, based on the displacement coordination condition and static equilibrium condition, the entire process of the interaction between the surrounding rock and pressure-relief support is iteratively solved.The displacement and support force of the elastic surrounding rock interaction model are calculated using the displacement approximation method, and the relative displacement graph of elastic surrounding rock and support (Fig 14  relative displacement between the rock and the support tends to converge.As shown in Fig 15, the initial support force value is 100KN, which decreases with the number of iterations; after 5 iterations, the support force value tends to converge.The calculation results show that after using pressure components in rigid support, the surrounding rock and support transition from a rigid support stable state to a pressure-relief support stable state.

Results of the support force approximation method
The lower and upper limits of the support force value are denoted as P down and P up , respectively, with the lower limit being 0 and the upper limit being 10 times the initial geostress.The relative displacement and support force of the elastic-plastic surrounding rock and pressurerelief support interaction model are calculated using the support force approximation method, and the displacement graph of elastic-plastic surrounding rock and support (

Convergence of different numerical approximation methods
The relative displacement of elastic surrounding rock and pressure-relief support interaction is calculated using the trisection method (u 0 ), and substitution method (u 0 R ¼ u * R ), and the impact of these three it can be seen that the trisection method has the slowest convergence speed, the displacement curve of the substitution method has sudden changes, and the bisection method has faster convergence speed and higher stability.

Stability analysis
The stability of two tunnel cases was analyzed by solving with a numerical approximation algorithm, verifying the accuracy and reliability of the numerical approximation algorithm in calculating the relative displacement and support force at the contact surface between the surrounding rock and yielding support.The results of the numerical approximation algorithm for the two tunnel cases are shown in Table 3.
Table 1 shows that the initial displacement changes of the elastic surrounding rock with support are minor, then gradually increase but the growth rate slows down.Compared to the elasto-plastic type, the stability of elastic surrounding rock and support may be achieved more quickly.The initial displacement changes of elasto-plastic surrounding rock with support are larger and then gradually stabilize.As the number of iterations increases, the displacement changes of the surrounding rock and support gradually decrease, tending towards a fixed value.

Comparison of methods
The method in this paper shows significant advantages over traditional methods such as linear superposition in dealing with the interaction between surrounding rock and support structures in tunnels under high ground stress conditions with soft surrounding rock.The displacement approximation method and support force approximation method proposed in this paper better handle the complex issues of displacement coordination and static equilibrium.Although these methods are relatively complex in the calculation process, they provide higher accuracy, which means that in practical applications, it can reduce the repeated calculations and design modifications caused by inaccurate models, thereby saving time and resources.

Existing limitations and challenges
This research also reveals several limitations when proposing new calculation methods: First, the complexity and uncertainty of actual geological conditions may become the main challenges in implementing these methods.Due to the variability of geological conditions, the required models must be flexible enough to adapt to different engineering environments.Second, although these methods perform well in research, their universality in different geological conditions and engineering environments needs further verification, especially in some extreme or uncommon geological conditions, where these methods may need to be adjusted to better fit the actual situation, beyond the basic assumptions of the current model.In addition, the effectiveness of these methods may vary in different practical application scenarios.
Despite these challenges, the methods in this paper provide a computational approach for calculating the interaction between surrounding rock and support structures.Although

Conclusion
To calculate the relative displacement at the contact surface between surrounding rock and pressure-relief support and the support force of pressure-relief support, this paper starts from the interaction mechanism between surrounding rock and rigid support.Based on the displacement coordination condition and static equilibrium condition, the displacement approximation method and support force approximation method are presented.Through the analysis of case studies, the effectiveness of the methods proposed in this paper is verified, and the main conclusions are as follows: (1) The displacement approximation method and support force approximation method can efficiently and accurately calculate the relative displacement at the contact surface between the surrounding rock and yielding support and the support force of the yielding support.The calculation results not only meet the static equilibrium condition but also the displacement coordination condition, proving the accuracy and reliability of these two approximation methods.In addition, these two methods have fewer iterations and faster convergence.
(2) The convergence of the bisection method, the trisection method, and the substitution method, three numerical approximation algorithms, is comparatively analyzed.The bisection method performs better in terms of stability and convergence speed when solving the interaction problem between surrounding rock and yielding support, showing superiority over the trisection method and substitution method in terms of convergence speed and stability.
(3) This paper adopts a method combining theoretical analysis, numerical approximation algorithms, and laboratory experimental data sets, providing an important reference for the design of yielding support structures in tunnels with soft surrounding rock under high ground stress conditions.
Fig 1, and the rigid support mechanical model is depicted in Fig 2.

Fig 8 .
Fig 8. Radial mises stress contour map of yielding components.https://doi.org/10.1371/journal.pone.0299426.g008 ), and the support force curve of elastic surrounding rock (Fig 15) are drawn.As shown in Fig 14, the relative displacement value at the beginning of the iteration is small, and the relative displacement between the surrounding rock and the support increases with the number of iterations; after 5 iterations, the

Fig 9 .
Fig 9. Axial mises stress contour map of yielding components.https://doi.org/10.1371/journal.pone.0299426.g009 Fig 16) and the support force and its range change graph (Fig 17) are drawn.As shown in Fig 16, at the 10th iteration, the relative displacement of the surrounding rock and support is equal (meeting the displacement coordination condition).As shown in Fig 17, the range of the support force value is halved with each iteration, and from the 5th iteration, the change in the support force value decreases, converging within no more than 10 iterations.

Fig 12 .
Fig 12. Diagram of the failure form of yielding components.(a) Failure form of slot and bolt hole.(b) Failure form of inner and outer slot steel.https://doi.org/10.1371/journal.pone.0299426.g012

Table 3 . Displacement of surrounding rock and support in tunnel cases-numerical approximation algorithm solution results. Elastic surrounding rock and support Unit: 10 -3 m Elasto-plastic surrounding rock and support Unit: 10 -3 m Iteration number Surrounding rock displacement Support displacement Surrounding rock displacement Support displacement
://doi.org/10.1371/journal.pone.0299426.t003practical applications may encounter challenges of geological complexity and technical requirements, in the long run, these methods still have room for development in increasing the accuracy of engineering design, ensuring structural safety, and improving economic efficiency. https