An Analytical Solution for Stress and Displacement Field for Arbitrarily Positioned Twin-Tunnel Excavation at Great Depth

In this study, an analytical solution for the stress and displacement field for arbitrarily positioned twin-tunnel excavation at great depth is proposed. The twin-tunnel is subjected to equivalent loads on the tunnel periphery and uniform loads at infinity. The analytical continuation of the complex variable method is used for stress “stair”, owing to the equivalent loads of the tunnel support. Subsequently, the Schwarz alternatingmethod andCauchy integral formula are used.The redundant surface forces on the periphery are expanded into Laurent series in the sequence of iterations, and a reasonable series truncation is adopted. Comparing with several existing solutions covering special cases from different perspectives, it is found that the proposed solution meets accuracy requirements for tunnel center distance greater than 2.2 times tunnel radius after 1.5 iterations. The proposed solution is fast and yields reference value to tunnel engineering.


Introduction
Tunnel excavation at great depth is common for highway, railway, metro, and water passage.For instance, the recently constructed Xi'an-Chengdu high-speed railway in China passes through the Qinling Mountains at maximum tunnelling depth of more than 1000 m.The typical cross section is that of a twin-tunnel.The excavation of twin-tunnel at great depth is described as openings in an infinite plane.In the early 20th century, the bipolar-coordinate method (BCM) was popular and many models were proposed for it.Multiple analytical solutions for single, double, and multiple openings in an infinite plane using the BCM were conducted [1][2][3][4].The BCM is suitable for circular tunnels or openings, but can not account for tunnels with complicated shapes.Besides, the BCM does not yield solutions for displacement field, which is important in tunnel design and construction and has hindered the application of BCM.
The complex variable method (CVM) for tunnel excavations at great depth and openings in an infinite plane with complex boundary conditions were proposed [5].Using the CVM, solutions for twin-tunnel excavations at great depth are found based on single tunnel excavations at great depth and single openings in an infinite plane.Abundant researches of single tunnel excavation or single opening in an infinite plane are conducted.The stress and displacement in single tunnel excavation at great depth were studied [6,7].Analytical solutions for deep rectangular structures and a square tunnel were conducted analytically [8,9].A correction factor was adopted for improving the accuracy of the CVM [10].Single tunnels with liner by considering the stress continuity and deformation compatibility were investigated [11][12][13][14][15].The cover cracking on a bar at great depth due to corrosion and uniform stresses at infinity using the CVM was investigated [16].An infinite plane weakened by a hole with arbitrary shapes using Cauchy integral of CVM was conducted [17].
Combined with the Schwarz alternating method [18], the CVM analytically yields multitunnel excavation solutions as well as on multiopenings in an infinite plane.Tunnel excavation solutions using CVM were elaborated [19,20].Analytical solutions for the stress and displacement field in bicircular and multielliptical openings in an infinite plane were presented [21][22][23][24]."Pseudotraction" method was used to investigate multiopening problems [25].The stress concentration and distribution of bicircular openings in an infinite plane were studied [26][27][28].Moreover, particular conformal mapping functions exist to minimize iterations of the Schwarz alternating method and gain solution [29][30][31].The Cauchy integral of CVM is also used for an infinite plane weakened by finite holes [32,33].
The existing analytical solutions for twin-tunnel excavation and double openings in an infinite plane mentioned above can be divided into BCM and CVM.Because of its shortcomings, BCM can not deal with complicated-shape tunnel excavations and displacement field.CVM is more adaptable to complicated tunnel shapes and displacement field, but currently available solutions are complex.In general, the solutions of the CVM are derived by two methods.The first is by combining the Schwarz alternating method and Cauchy integral.Cauchy integral can not solve complex or mixed boundary problems [19]; therefore, its application is limited.Furthermore, the stress potential functions are long and computationally costly even for circular twintunnel [20].The second method is different in form, but based on the same theory.The twin-tunnel plane is mapped onto an annulus to minimize iterations, and a system of linear equations is solved; however, the derivation process is complex.Both methods are complex in practice; thus, simple analytical solutions are needed.

Assumptions and Problem Definition
This study focuses on deriving a simple analytical solution for the stress and displacement field for arbitrarily positioned circular twin-tunnel excavation at great depth.Rock mass and tunnel lining are simplified as follows: (1) The rock mass is homogeneous, isotropic with little deformation, and elastic.
(2) Both tunnels are at great depth, and gravity field is transformed into the initial stress field at infinity [5,19,20] (Figure 1).
(3) Both tunnels are circular and of the same size.
(4) Tunnel support is immediately applied after excavation and is simplified as equivalent loads uniformly and is radially distributed on the tunnel periphery, ignoring the loss of the initial stress field.
(5) Both tunnels are simultaneously excavated, but are arbitrarily positioned.
(6) Tunnelling face and three-dimensional effects are not considered in this study.
In Figure 1, tunnel peripheries are marked as "Tunnel 1" and "Tunnel 2" in blue and orange, respectively, and their centers are used to establish the complex plane  1 and  2 .In the complex plane  1 , the center distance of the tunnels is c, the modulus is ||, and argument is .The initial stress field at infinity consists of two uniformly distributed components: the horizontal   and the vertical   .The equivalent loads in red are uniform on the peripheries of both tunnels artificially, and normal load is .The complex plane  1 and  2 are conformally mapped onto planes  1 and  2 , respectively, and the tunnels are geometrically similar.
The proposed solution combines the CVM [5] and the Schwarz alternating method [18][19][20].At the initial phase of the solution, a Cauchy-type integral is used to solve the stress "stair" problem on the tunnel periphery, owing to the equivalent loads.The Schwarz alternating method is applied to the rest parts of the solution, and Laurent series expansion and truncation simplify the complex redundant surface traction to use Cauchy integral formula.Using Mathematica software, the iterations avoid complex formulae.The entire solution procedure is accurate and efficient.

Complex Variable Method
According to the CVM proposed by Muskhelishvili [5], the stresses vary with the stress potential functions on the abstract plane z, where   ,   , and  xy are the stress components in Cartesian coordinates and (z) and (z) are the stress potential functions.() is the conjugate of (z), and i= √ −1.
A similar equation exists for displacements, where  and V are the horizontal and vertical displacement, respectively, G stands for the shear modulus of the rock mass, and  is constant and =3-4, where  is the Poisson's ratio of the rock mass.
To use the properties of the unit circle, we transform the stress and displacement expression from Cartesian coordinates into polar coordinates.The stress components are where  is the direction angle of an arbitrary point of the rock mass in plane z and   ,   , and   are the stress components in polar coordinates.The rock mass in the exterior area of the tunnel on plane  is mapped onto the exterior area of the unit circle on the abstract plane , with the mapping function where  is the tunnel radius.Equation ( 4) shows that the tunnel on plane  and the mapping one on plane  share geometric similarity.So we have where  is the argument of the mapping point on plane  corresponding to the point on plane z, also known as the direction angle on plane .
Combining (1), (3), and (5), we obtain the stress components varying with the complex potential functions in polar coordinates where where  is the modulus.According to ( 2) and ( 4), the displacements is expressed in polar coordinates as where   and   are the radial and circumferential displacement components, respectively.The equivalent loads of tunnel support are uniformly and radially distributed on the circular periphery, with the resultant force being equal to zero [5,19,20].Thus, by eliminating the multivalued terms [5,19,20], the stress potential functions are where  and B' are real constants, Boundary condition of the surface force [5,19,20] is where  is the tunnel periphery on plane .

Solution Procedure
We iteratively combine the CVM and the Schwarz alternating method to solve the stress potential functions.Thus, the procedure can be stated as follows: (1) Subjected to the initial stress field in Figure 1, Tunnel 1 is excavated and immediately supported.Equivalent loads are applied on the periphery of Tunnel 1, and corresponding stress potential functions are obtained using analytical continuation.This is the initial stage of the solution process and is named as "the 0.5th iteration".
(2) The surface force on the periphery of Tunnel 2 is named as "redundant surface traction", and is acquired from the stress field of the Tunnel 1 excavation.Then, opposite loads of the same magnitude that are named as "antisurface traction" are artificially imposed on the periphery with the equivalent loads of the Tunnel 2 support.Both forces generate the stress potential functions of the Tunnel 2 excavation, and this is "the first iteration".
(3) As above, a round of iterations is complete.Henceforth, we eliminate the redundant surface traction on the peripheries of Tunnels 1 and 2 and obtain the corresponding stress potential functions repeatedly, until the accuracy requirements are met.It is noteworthy that the equivalent loads and initial stress field should not be recounted in.Based on the superposition principle, stress potential functions accumulate after the coordinates transformation [20] and the total stress potential functions are obtained.
The details of the solution process are given below, where subscripts 1 and 2 are used to denote Tunnels 1 and 2, respectively, and the superscripts denote twice of iterations.
4.1.The 0.5th Iteration.The application of equivalent loads on the periphery of Tunnel 1 causes "stair" phenomenon; thus, analytical continuation [19] is appropriate.Equivalent loads, with the same direction of the corresponding polar axes, have positive sign, and the analytical partition function is where  1 ( 1 ) is the analytical partition function of the 0.5th iteration, p is the value of equivalent loads, and  1 stands for the periphery of Tunnel 1 in coordinates  1 .

Mathematical Problems in Engineering
Solving (12) yields where  1+ ( 1 ) and  1− ( 1 ) stand for the analytical partition functions in the interior and exterior areas of the unit circle on plane  1 respectively, and  1 1 ( 1 ) serves as the intermediate stress potential function of the 0.5th iteration.The three stress potential functions  1  1 ( 1 ),  1 1 ( 1 ), and The combination of ( 13) and ( 14) yields Integrating (15) and ignoring the constant terms, we obtain the following: The 0.5th iteration is complete.

The First Iteration.
In the first iteration, we should not recount the initial stress field but note that the surface forces on the periphery of Tunnel 2 have two components: the equivalent loads of Tunnel 2 and the antisurface traction.For convenience, we can separately calculate the stress potential functions of each part and add them later.For equivalent loads, the analytical continuation mentioned above suits; meanwhile, Cauchy integral formula [5,19,20] serves well for the antisurface traction.Similar to the 0.5th iteration, analytical continuation yields where  2 2(1) ( 2 ) and  2 2(1) ( 2 ) are the stress potential functions of the equivalent loads on Tunnel 2 periphery.By subtracting (17) from ( 16), the solutions obtained are the same as those proposed by Chen [19] and Lu [20].
The stress potential functions of the equivalent loads are obtained by the process mentioned above.Next, we use the Cauchy integral formula to obtain the stress potential functions for  2 ( 2 ), also known as the antisurface traction on the periphery of Tunnel 2 owing to the excavation of Tunnel 1.
Figure 1 shows that the center distance between Tunnels 1 and 2 on plane  1 is c, which transforms to c/R by the conformal mapping of (4).Therefore, based on the coordinates transformation [20],  2 ( 2 ) could be obtained by combining (11) and ( 16) Clearly, for each term of ( 18), || > 2 and | 2 |=1.Based on the incomplete derivation and power series expansion, we obtain the following: Substituting ( 19) into (18), we can get the Laurent series where   is the coefficient of the Laurent series.Equation ( 20) is similar to Lu's solution [20].Because || > 2, (/||) 10 < 0.001, we can truncate (19) and reserve terms with absolute value of power less than 10; thus, ( 20) is written as Each coefficient of ( 21) is implicitly obtained using Mathematica, and   (2 ≤  ≤ 9) equals to zero, and we get the following: To neutralize the effect of redundant surface traction owing to the excavation of Tunnel 1, we artificially impose the antisurface traction  2 2 ( 2 )=− 1 2 ( 2 ) on the periphery of Tunnel 2, and we obtain the following: where  2 2(2) ( 2 ) and  2 2(2) ( 2 ) are the stress potential functions on the periphery of Tunnel 2 owing to  2 2 ( 2 ).By combining ( 21) and ( 22), we obtain By taking conjugate of ( 24), we obtain By applying the Harnack theorem [5,19,20] and Cauchy integral formula to ( 24) and ( 25), we obtain We add ( 17) and ( 26) to obtain Applying the coordinates transformation formula [20] to (27) yields The first iteration is complete.Thus, the 0.5th and the first iterations compose the first round of iterations.

The 1.5th
Iteration.The 1.5th iteration means that when using the Schwarz alternating method in the second round of iterations, the redundant surface traction on Tunnel 1 is eliminated, whereas the one on Tunnel 2 is not.Lu [20] proposed that, after 1.5 iterations, the stress potential functions yield sufficiently accurate results.
Substituting ( 28) into (11) yields Similar to (20), we expand (29) using the incomplete derivation and the power series expansions, and obtain the following: We substitute (30) and ( 31) into ( 29) and obtain the redundant surface traction of the 1.5th iteration where   is the coefficient of the Laurent series.We truncate (30) and ( 31) and reserve terms with absolute value of power less than 10, as in the first iteration.Equation (32) then is written as Similar to the first iteration, the coefficients of ( 33) are also obtained implicitly using Mathematica.After applying the antisurface traction  3 1 ( 1 )=− 2 1 ( 1 ) on the periphery of Tunnel 1, we obtain As in the first iteration, we apply Harnarck theorem and Cauchy integral formula to (34) and its conjugate and obtain the following: The 1.5th iteration is complete.Thus, we acquire the stress potential functions after 1.5 iterations, and ( 16), (28), and (35) are the results.

Stress and Displacement Field.
Based on the superposition principle, we add ( 16), (28), and (35) to obtain the total stress potential functions on plane  1 after 1.5 iterations Using the coordinates transformation [20], we obtain the total stress potential functions in plane  2 We substitute (36) and (37) into ( 6) and ( 8), respectively, and obtain the stress and displacement field of the rock mass surrounding Tunnels 1 and 2 where Re and Im are the real and imaginary parts of the expressions in the brackets.Similarly, we solve (39) to obtain the radial and circumferential displacements on the tunnel periphery as

Results and Discussion
Ling [4] proposed analytical solutions for an infinite plane containing two circular openings of the same size and at the same horizontal level subjected to all-around, longitudinal, and transverse tension of unit force at infinity (Figure 2) and analyzed the circumferential stresses on the periphery of the right opening.Tables 1, 2, and 3 list the circumferential stresses calculated by the proposed solution after 1.5 iterations and the ones provided by Ling.Tables 1, 2, and 3 present that the errors greatly decrease with the increasing center distance of the two openings.As the center distance increases, R/|c| decreases, and the accuracy improves for the same truncation terms in (19), (30), and (31).Tables 1, 2, and 3 show that the errors reach the 5% tolerance interval for the center distance thrice larger than the opening radius.We believe that the proposed solution would meet accuracy requirements for Whereas the accuracy of the proposed solution involving the center distance between twice and thrice of tunnel radius remains uncertain, further verification is needed.Circumferential stresses at particular positions for circular openings of the same size and horizontal level within an infinite plane subjected to longitudinal and transverse tension of the unit force were studied [25,29] (Figures 2(b) and 2(c)).The center distance covers twice and thrice the opening radius.They reported consistent results, which suggest agreement with Ling [4].Tables 4 and 5 give Haddon's circumferential stress results on the opening periphery and the angles to obtain them [29].Results by the proposed solution are given for comparison.
Tables 4 and 5 show that most errors are within 5% tolerance interval, when center distance is 2.2 times larger than the opening radius.Similar results were reported [28].Several errors at 180 ∘ exceed 5% tolerance interval, and we have to notice that the proposed solution has defects on depicting interaction between two openings for truncation is employed for the proposed solution.
For further accuracy verification of present solution, we should dense calculating results of circumferential stress on the periphery of the openings.Analytical solutions were given [3] for an infinite plane containing two circular openings subjected to longitudinal tension (Figure 2(b)).The circumferential stress on the periphery was with divisions of 10 ∘ for the polar angle, while the ratio of the center distance to the opening radius is 6.67.Comparisons between the proposed solution and Green's solution are given in Table 6.
Table 6 shows that the results of the proposed solution are approximately similar to Green's, which suggests that the

Periphery of Tunnel 1 PeripheryFigure 2 :
Figure 2: Schematic diagram of the unit force tension for all-around tension (a), longitudinal tension (b), and transverse tension (c).

Table 1 :
Circumferential stress owing to the all-around tension of unit force.

Table 2 :
Circumferential stress owing to the longitudinal tension of unit force.

Table 3 :
Circumferential stress owing to the transverse tension of unit force

Table 4 :
Circumferential stress owing to the transverse tension of unit force.