Improved Mobilized Strength Design Method for Multi-Support Excavation Deformation Analysis

Abstract

The safe and reliable design of underground spaces ensures the safety of a structure itself and its surroundings. The traditional Mobilized Strength Design (MSD) method for a multi-support excavation deformation analysis ignores the effects of soil parameters and excavation boundary conditions. Therefore, to compensate for the shortcomings of the existing MSD method, this paper proposes an improved mobilized strength design (IMSD) method for a multi-support excavation deformation analysis. The improved incremental deformation mechanism further considers the effect of the soil friction angle, and the effect of excavation depth and the first support on deformation energy are also considered. Further, the excavation calculation process based on the IMSD method is given, and the effects of different calculation parameters on the IMSD solution of excavation deformation are discussed. The results show that the IMSD method can effectively consider the effect of boundary conditions and the excavated process on the excavation deformation. The traditional MSD method underestimates the excavation deformation and surface settlement by an average of 15–23%, while the IMSD solution is more consistent with the measured values. The study results can provide a theoretical reference for the design of multi-support excavation.

1. Introduction

With the implementation of the strategy of a country with a strong transportation network, urban infrastructure construction and rail transportation have been rapidly developed [1,2,3,4]. The construction of underground projects has been required due to the massive construction of high-rise buildings, commercial centers, and subway projects [5,6,7]. Underground projects have been mostly located in urban centers, where the environment was very complex, especially around subway excavation with dense buildings, roads, and lifelines [8,9,10]. Once an excavation failed, it would be a fatal disaster for the surrounding buildings and lifelines [11,12,13]. Therefore, to ensure the safety of an excavation, surrounding buildings, and lifelines, it is very significant to study a reliable theoretical system to accurately predict excavation deformation and force.

With the continuous development of urban construction, the building environment around urban excavation has been becoming increasingly complex, making the deformation control design of excavation engineering increasingly significant [14,15]. Therefore, the prediction method of excavation deformation has received extensive attention. The current methods for analyzing excavation deformation and forces have been divided into three main categories: (1) the classical method; (2) the elastic foundation beam method; and (3) the finite element method [16]. The classical method simplifies the excavation boundary conditions and differs from the actual force pattern [17]. The elastic foundation beam method treats the soil as a homogeneous elastic semi-infinite body, ignoring the plastic deformation of the soil [18]. The finite element method has a complex selection of constitutive models with uncertainties in the input parameters [7,9,19]. Therefore, it is of great significance to investigate new reliable methods for excavation design.

With the development of the excavation theory, a new method from the energy point of view has emerged: the Mobilized Strength Design (MSD) method [20,21,22]. The MSD method adopted the plastic deformation mechanism to predict the excavation deformation and surface settlement through the energy conservation relationship between the soil and the supporting system [23,24,25]. Overall, compared with traditional design methods, the MSD method has three major advantages: (1) tightly formulated derivations using energy conservation and clear physical concepts [26,27]; (2) the nonlinear relationship of soil stress–strain and the dynamic development of excavation deformation can be considered; and (3) the non-homogeneous soil condition, burial depth of soil, and the soil principal stress rotation caused by the excavation can be considered [28].

Currently, the MSD method has been successfully applied to cantilever retaining walls [29,30,31], ground surface settlement caused by tunnel excavation [32,33,34], excavation-induced deformation and ground surface settlement [35,36,37], and nonlinear pile group response [38,39]. Bolton et al. [20,23] first proposed the mobilizable strength design (MSD) method and successfully applied it to wall displacement and stability studies in over-consolidated clays. Osman and Bolton [40] proposed a design method for cantilevered excavation considering the nonlinear soil stress–strain relationship based on the MSD method and elastoplasticity theory. Then, the effect of anisotropic soil was considered, and the application of the MSD method in excavation was further extended and improved [37]. Lam and Bolton [21] proposed different deformation mechanisms for wide and narrow excavation, considering the effect of excavation width on the surrounding soil deformation pattern. Moreover, the accuracy of the excavation MSD method, considering the mobilized shear strength of soils with different burial depths, has been verified and improved by model experiments and engineering cases [28,35]. Further, Liu et al. [41] and Wang et al. [22,42] considered the effect of soil non-homogeneity and anisotropy and introduced the compression deformation energy and bending deformation energy of supporting structures in the energy conservation equation.

However, the traditional incremental deformation mechanism adopted for excavation is inconsistent with the actual situation. The current MSD method to calculate the excavation deformation does not consider the effect of the soil friction angle. The range of plastic deformation caused by excavation only considers the length effect of the supporting piles, ignoring the effect of the excavation depth on the deformation energy. In addition, the traditional energy conservation formula does not give enough consideration to the deformation energy of the excavation supporting structure, especially the deformation energy generated by the first support.

Therefore, to address the above research deficiencies, this paper proposes an improved mobilized strength design (IMSD) method for multi-support excavation deformation analysis. The improved incremental deformation mechanism further considers the effect of soil friction angle, and the effect of excavation depth and first support on deformation energy is also considered. The reasonableness and accuracy of the newly proposed IMSD method were verified using engineering examples. Finally, the effects of different parameters on the IMSD solution of excavation deformation were discussed. These study results improve the prediction accuracy of excavation deformation and provide a theoretical basis for the excavation design.

2. Application of the IMSD Method in Excavation

2.1. Traditional MSD Method

The selection of the incremental deformation mechanism is the key to the MSD method. The MSD method proposed by Lam and Bolton [21] assumed that the incremental deformation function should be a cosine function, as in Equation (1). However, this method only considered the soil deformation between the lowest support and the hard soil layer, which was not consistent with the actual excavation deformation. Liu et al. [41] further improved the incremental deformation function by considering the soil deformation between the pile top and the hard soil layer. The incremental deformation curve was assumed to be a segmented cosine function based on the excavation deformation in the layered strata, as in Equation (2).

$${\delta }_{\mathrm{x}}={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi y}{\lambda }}\right)\right]$$

(1)

$$\left\{\begin{array}{l}{\delta }_{\mathrm{x}}={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{\pi y}{k\lambda }}\right)\right],-k\lambda \le y\le 0\\ {\delta }_{\mathrm{x}}={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{\pi (y-(1-2k)\lambda }{(1-k)\lambda }}\right)\right],-\lambda \le y\le -k\lambda \\ \lambda =\alpha s\end{array}\right.$$

(2)

where $y$ is the distance from the lowest support of any point below the lowest support, ${\delta }_x$ is the incremental displacement of the wall at $y$, ${\delta }_{\max }$ is the maximum incremental displacement of the soil, $s$ is the length of the supporting structure below the lowest support, $\lambda$ is the deformation wavelength, k is the ratio of the distance from the peak to the origin to the wavelength, defined as the crest coefficient of variation, and $\alpha$ is the deformation influence coefficient, which takes the value of 0.5–2.

2.2. Improved MSD Method

2.2.1. Plastic Deformation Mechanism of Excavation

To compensate for the shortcomings of the traditional MSD method, the incremental deformation mechanism of the excavation was improved (see Figure 1). In this paper, the range of incremental deformation outside the excavation lay between O’rourke [43] and Liu et al. [41], while the influence range inside the excavation was reduced. From Figure 1, the plastic deformation flow is shown as the dashed line with arrows, including active zone ABCF, sectors CEF and EHF, and passive zone FHI. The ABCF zone considers the horizontal displacement of the area above the excavation surface and forms a quadrilateral zone where the angle between the flow line and the horizontal is the angle of internal friction. The CFE zone is a sector centered at point F. The angle between the sideline CF and the vertical is the angle of internal friction. The FH line is perpendicular to the CF line, dividing the sector FEH zone and the triangle FHI zone. The influence range of the improved incremental deformation mechanism varies with the excavation depth. In addition, a plastic deformation field, considering the displacement coordination and soil friction angle, is proposed.

It is assumed that the excavation width is greater than $2(1+{\alpha }_2)h\sec \theta$, and the deformation is not affected by the interaction of the two sides of the supporting structures. According to the excavation deformation in the actual layered stratigraphy, the virtual top deformation zero point of the deformation curve is assumed to be located above the excavation top. The improved plastic deformation field further considers the effect of the soil friction angle. Meanwhile, two parameters, the virtual fixed point location coefficient, ${\alpha }_1$, and the influence range coefficient of excavation depth, ${\alpha }_2$, are introduced into the incremental deformation function so as to reflect the influence of excavation depth and the first support on the deformation energy. If the incremental deformation function ignores the compression deformation energy of the first support, the calculated excavation top displacement is always zero. In addition, if the deformation increment function assumes that the deformation wavelength is a fixed value and fails to correlate with the excavation depth, the influence range of soil deformation is quite different from the actual situation when the excavation depth is small.

The improved incremental deformation cosine function is expressed as follows:

$${\delta }_{\mathrm{x}}={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(y+{\alpha }_{1}h\right)}{\lambda }}\right)\right]$$

(3)

$$\lambda =\left(1+{\alpha }_1+{\alpha }_2\right)h$$

(4)

$$k={\displaystyle \frac{h}{\lambda }}$$

(5)

where ${\alpha }_1$ is the coefficient of the virtual fixed point position, ${\alpha }_2$ is the influence range coefficient of excavation depth, $h$ is the excavation depth, $\lambda$ is the deformation wavelength, and $k$ is the crest coefficient of variation.

2.2.2. Total Loss of Potential Energy of the Soil

Excavation unloading disrupts the initial stress balance. Horizontal additional stresses caused by excavation act to deform the supporting structure, and the soil behind the supporting structure is also displaced horizontally and vertically. The tendency of soil bodies to move relative to each other results in shear stresses and, hence, shear displacements. This result is accompanied by the bending deformation of the supporting structure and the compression deformation of the inner support, which generates the corresponding deformation energy, respectively. Based on Lam and Bolton [21], the effect of the incremental compression deformation energy of internal supports was further considered in the energy conservation equation:

$$\Delta P=\Delta W+\Delta U+\Delta V$$

(6)

where $\Delta P$ is the total loss of potential energy of the soil, $\Delta W$ is the total work carried out in shearing the soil, $\Delta U$ is elastic strain energy, and $\Delta V$ is the incremental compression deformation energy of the internal support.

The increment of the total potential energy of the excavation system in the mth stage of construction can be expressed as follows:

$$\Delta P={\displaystyle \sum _{i=1}^n\left[{\displaystyle \int {\gamma }_{sat}\left(m,i\right){\delta }_y\left(m,i\right)dA}\right]}$$

(7)

where ${\gamma }_{sat}\left(m,i\right)$ is the saturated unit weight of soil in the ith layer for the mth construction, ${\delta }_y\left(m,i\right)$ is the vertical component of displacement of soil in the ith layer for the mth construction, and $A$ is the site soil area within the influence range.

According to the deformation mechanism shown in Figure 1, the excavation soil body can be divided into four deformation zones: ABCF, CEF, EHF, and FHI. The incremental displacement at the mth excavation condition can be expressed as follows:

$${\delta }_y=y\left({\delta }_{\max }\right)$$

(8)

$${\delta }_x=x\left({\delta }_{\max }\right)$$

(9)

The details are shown in the following Equations (10)–(17), where $r^2=x^2+y^2$, $r$ and $x$ with the same positive and negative sign.

(1) In the zone ABCF, with point A as the origin of the coordinates, then:

$${\delta }_{\mathrm{x}}=-{\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\sqrt{x^2+y^2}+\left(1+{\alpha }_1\right)h\right)\tan \theta }{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]\sin \theta$$

(10)

$${\delta }_y={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\sqrt{x^2+y^2}+\left(1+{\alpha }_1\right)h\right)\tan \theta }{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]\cos \theta$$

(11)

(2) In the zone CEF, with point F as the origin of the coordinates, then:

$${\delta }_{\mathrm{x}}=-{\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\sqrt{x^2+y^2}+\left(1+{\alpha }_1\right)h\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]{\displaystyle \frac{y}{\sqrt{x^2+y^2}}}$$

(12)

$${\delta }_y={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\sqrt{x^2+y^2}+\left(1+{\alpha }_1\right)h\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]{\displaystyle \frac{x}{\sqrt{x^2+y^2}}}$$

(13)

(3) In the zone FEH, with point F as the origin of the coordinates, then:

$${\delta }_{\mathrm{x}}={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\sqrt{x^2+y^2}+\left(1+{\alpha }_1\right)h\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]{\displaystyle \frac{y}{\sqrt{x^2+y^2}}}$$

(14)

$${\delta }_y=-{\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\sqrt{x^2+y^2}+\left(1+{\alpha }_1\right)h\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]{\displaystyle \frac{x}{\sqrt{x^2+y^2}}}$$

(15)

(4) In the zone FHI, with point F as the origin of the coordinates, then:

$${\delta }_{\mathrm{x}}={\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\frac{x}{\tan \theta }+y+\left(1+{\alpha }_1\right)h\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]\sin \theta$$

(16)

$${\delta }_y=-{\displaystyle \frac{{\delta }_{\max }}{2}}\left[1-\cos \left({\displaystyle \frac{2\pi \left(\frac{x}{\tan \theta }+y+\left(1+{\alpha }_1\right)h\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]\cos \theta$$

(17)

2.2.3. The Total Work Performed in Shearing the Soil

The calculation of the internal energy increment of the soil can be referred to Liu et al. [41]:

$$\Delta W={\displaystyle \sum _{i=1}^n\left[{\displaystyle \int C_{\mathrm{u},\mathrm{mob}}\left(m,i\right)\left|{\delta }_{\gamma }\left(m,i\right)\right|}dA\right]}$$

(18)

where $C_{\mathrm{u},\mathrm{mob}}\left(m,i\right)$ is the undrained shear strength of soil in the ith layer for the mth stage of construction, and ${\delta }_{\gamma }\left(m,i\right)$ is the shear strain increment of soil in the ith layer for the mth stage of construction.

${\delta }_{\gamma }\left(m,i\right)$ can be expressed as a partial differential function of the displacement:

$${\delta }_{\gamma }\left(m,i\right)=\sqrt{{\left({\displaystyle \frac{\partial {\delta }_x}{\partial y}}+{\displaystyle \frac{\partial {\delta }_y}{\partial x}}\right)}^2-4{\displaystyle \frac{\partial {\delta }_x}{\partial x}}\cdot {\displaystyle \frac{\partial {\delta }_y}{\partial y}}}$$

(19)

In the zone ABCF:

$${\delta }_{\gamma }={\displaystyle \frac{2\pi {\delta }_{\max }}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\left|\sin \theta \times \sin \left({\displaystyle \frac{2\pi \left(x+y\cot \theta \right)\tan \theta }{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right|$$

(20)

In the zones CEF and EHF:

$$\begin{array}{l}{\delta }_{\gamma }={\displaystyle \frac{{\delta }_{\max }}{\sqrt{2}\left(1+{\alpha }_1+{\alpha }_2\right)hr}}\left|\sin \left({\displaystyle \frac{\pi (h+r)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\times \left[4{\pi }^2{r}^2\left(1+\cos \left({\displaystyle \frac{2\pi (h+r)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right)+\right.\right.\\ {\left(1+{\alpha }_1+{\alpha }_2\right)}^2{h}^2\left(1-\cos \left({\displaystyle \frac{2\pi (h+r)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right)-4\pi r\left(1+{\alpha }_1+{\alpha }_2\right)h\times \sin \left.\left.\left({\displaystyle \frac{2\pi (h+r)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right]\right|\end{array}$$

(21)

In the zone EFI:

$${\delta }_{\gamma }={\displaystyle \frac{\pi {\delta }_{\max }}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\left|\csc \theta \times \sin \left({\displaystyle \frac{2\pi \left(h+y+x\cot \theta \right)}{\left(1+{\alpha }_1+{\alpha }_2\right)h}}\right)\right|$$

(22)

2.2.4. Elastic Strain Energy

Elastic strain energy can be obtained by integrating the bending stiffness with the incremental horizontal displacement within the burial wall depth:

$$\Delta U={\displaystyle \frac{EI}{2}}{{\displaystyle {\int }_0^h\left[\frac{d^2{\delta }_x}{d{y}^2}\right]}}^{2}dy$$

(23)

The excavation deformation curve conforms to the cosine function:

$$\Delta U={\displaystyle \frac{\pi EI{\delta }_{\max }}{2{\left(1+{\alpha }_1+{\alpha }_2\right)}^2{h}^2}}\left[\sin {\displaystyle \frac{2\pi \left(1+{\alpha }_1\right)}{\left(1+{\alpha }_1+{\alpha }_2\right)}}-\sin {\displaystyle \frac{2\pi {\alpha }_1}{\left(1+{\alpha }_1+{\alpha }_2\right)}}\right]$$

(24)

where $EI$ is the bending stiffness of the supporting structure.

2.2.5. Compression Deformation Energy of the Internal Support

The calculation of compression deformation energy of the internal support can be referred to Liu et al. [41]:

$$\Delta V={\displaystyle \sum \left[\frac{E_{\mathrm{P}}{A}_{\mathrm{P}}}{2l_{\mathrm{P}}}{\left({\delta }_{\mathrm{x},\mathrm{p}}\right)}^2\sin \theta \right]}$$

(25)

where $E_{\mathrm{P}}{A}_{\mathrm{P}}$ is the bending stiffness of internal support.

2.3. Determination of Calculation Parameters

Three significant calculation parameters affect the IMSD solution, including the following: (1) the deformation-influenced wavelength, (2) the soil shear strength, and (3) the soil stress–strain curve. To match the actual engineering, these significant calculation parameters in the IMSD method were improved.

2.3.1. Improved Deformation-Influenced Wavelength

The current incremental deformation function assumes the deformation-influenced wavelength as a fixed value and fails to account for the effect of excavation depth. When the excavation depth is small, the influence scope of soil deformation differs greatly from the actual deformation. Compared with the fixed deformation-influenced wavelength in the traditional MSD method, the deformation-influenced wavelength determined by the support load-sharing method introduced in this paper can more realistically reflect the influence range of the soil body [27,41]. Therefore, the support load-sharing method was used to determine the deformation-influenced wavelength in the MSD method.

It is assumed that the maximum soil deformation influence range of the excavation is two times the actual insertion depth of the soil. Therefore, the deformation-influenced wavelength can be expressed as a relation between the insertion depth below and above the reverse bend points of the supporting structure and the excavation depth (see Figure 2). The improved deformation-influenced wavelength is calculated as follows:

$$\lambda =2(x+y)+h$$

(26)

where $\lambda$ is the deformation wavelength, x and y are the insertion depth below and above the reverse bend points of the supporting structure, respectively, and $h$ is the excavation depth.

2.3.2. Soil Shear Strength

Stress changes in the soil caused by excavation in soft clay can be regarded as undrained behavior, and the soil strength eigenvalues for this process can be adopted as undrained shear strength. The undrained shear strength of soft clays is characterized by anisotropy due to the depositional history and initial consolidation state and is related to the overlying effective stress. Due to the soil stress spindle rotation caused by the excavation unloading, the soil stress state is different at different excavation locations, so the soil shear strength at different locations is different. Based on the model proposed by Casagrande and Carillo [44], formulas for the shear strength of soils considering different stress states and depth conditions are obtained:

$$s_{u\alpha }=\eta \xi +\xi (1-\eta ){\cos }^2(\theta -{\displaystyle \frac{\pi }{4}}+{\displaystyle \frac{{\phi }^{\prime }}{2}}){{\sigma }^{\prime }}_v$$

(27)

$${{\sigma }^{\prime }}_v={\displaystyle \sum _i{\gamma }_{sat,i}{z}_i}$$

(28)

where $s_{u\alpha }$ is the soil shear strength, $\eta$ is the anisotropy ratio, $\xi$ is the ratio of the vertical undrained shear strength to the overlying effective stress, $\theta$ is the deflection angle of the main axis, ${{\sigma }^{\prime }}_v$ is the overlying effective stress, ${\phi }^{\prime }$ is the effective angle of internal friction, ${\gamma }_{sat}$ is the saturated unit weight of the soil, and z is the depth.

2.3.3. Stress–Strain Curve of Soil

The theoretical core of the MSD method is the shear strength mobilization coefficient of the soil and the function of soil shear strain. Vardanega and Bolton [33] obtained undrained shear stress–strain curves for cohesive soils conforming to an exponential function considering the effects of soil nonlinearity:

$$\beta ={\displaystyle \frac{S_{\mathrm{u},\mathrm{mob}}}{S_{\mathrm{u}}}}=a{\left({\displaystyle \frac{\gamma }{{\gamma }_{0.5}}}\right)}^b$$

(29)

where $S_{\mathrm{u}}\left(i\right)$ is the peak strength of the stress–strain curve in the ith layer, ${\gamma }_{0.5}\left(i\right)$ is the shear strain at corresponding to $0.5S_{\mathrm{u}}\left(i\right)$ on the stress–strain curve, a and b are the fitting parameters, and $\beta$ is the undrained shear strength exertion coefficient.

3. Suggested Calculation Process

This paper proposes a multi-support excavation deformation analysis method based on the improved MSD method. The main calculation process is divided into eight steps (see Figure 3):

(1) Input the basic parameters: Input soil parameters, such as saturated unit weight, ${\gamma }_{\mathrm{sat}}$, shear stress, $S_{\mathrm{u}}$, and shear strain, ${\gamma }_{0.5}$. Input the supporting structure parameters, such as bending stiffness (EI), length ($l$), and compressive stiffness ($EA$) of the inner support, the angle between the support and the horizontal direction of the supporting structure ($\theta$), excavation depth ($h$), and the layered excavation depth.

(2) Calculate wavelength and crest variation coefficient: Calculate the wavelength, $\lambda$, and crest variation coefficient, $k$, using Equation (3).

(3) Calculate the increment of different systems: Equations (18), (23), and (25) are used to calculate the total work performed in shearing the soil, $\Delta W$, the elastic deformation energy, $\Delta U$, and the incremental compression deformation energy of the internal support, $\Delta V$, respectively.

(4) Calculate the total loss of potential energy of the soil: According to the energy conservation in the excavation system, the total loss of potential energy of the soil, $\Delta P$, is calculated by Equation (6).

(5) Determine the cosine function of the incremental deformation curve: The maximum horizontal displacement increment of the supporting structure is determined, and the cosine function of the incremental deformation curve is then determined.

(6) Determine whether the excavation step m is equal to the final excavation step M: If $m=M$, output the full deformation curve. If $m\ne M$, continue to input parameters for the superposition of incremental deformation.

(7) Judge whether to meet the control standard: When the obtained deformation does not meet the requirements of the control standards, the excavation construction program can be adjusted appropriately and repeatedly accounted for until the requirements are met.

(8) Output the results of excavation calculation: Output the results of the excavation displacement and the ground surface settlement.

4. Engineering Verification

4.1. Case 1

Liu et al. [41] analyzed the excavation deformation characteristics using the MSD method. The Jinzhongjie Station excavation of Tianjin Metro Line 6 was 499.9 m long × 16.6 m deep × 20.7 m wide. The support system adopted a combination of a diaphragm wall (0.8 m wide × 34.1 m deep) and mixed internal support with the 1st concrete strut (0.8 × 0.8 m in section size, $EI=1.28\times {10}^6 \mathrm{kN}·{\mathrm{m}}^2$) and 2nd, 3rd, and 4th steel pile struts (0.609 m or 0.8 m in diameter and 0.016 m in thickness, $E=210 \mathrm{GPa}$). The soil around the excavation was mainly silty clay and fine sand. Excavation resulted in a weakened resistance of the supporting structure to the soil. The horizontal resistance coefficients for the soils below and above the excavation surface were $2.4\times {10}^6 \mathrm{N}·{\mathrm{m}}^{-3}$ and $1.2\times {10}^6 \mathrm{N}·{\mathrm{m}}^{-3}$, respectively. The corresponding deformation coefficients for different soil resistance coefficients were 0.1195 and 0.1005, respectively. Refer to Liu et al. [41] for specific project overview details.

Figure 4 compares the pile deformation and surface settlement calculated by different methods. By comparing the measured values and the MSD solution proposed by Liu et al. [41], the IMSD solution proposed in this paper was more consistent with the measured pile deformation and ground surface settlement. This finding was because the proposed IMSD solution adopted a plastic incremental deformation mechanism considering the internal friction angle, which could more accurately reflect the soil deformation and significantly improve the accuracy of the maximum excavation displacement. In addition, the MSD solution proposed by Liu et al. [41] underestimated the measured pile deformation, which differs from the measured value. This result was because the same incremental deformation function was used for the surface settlement and pile deformation in the plastic incremental deformation mechanism. The incremental deformation mechanism adopted by Liu et al. [41] ignored the effect of excavation depth and first support on deformation energy. In addition, the proposed model solution had an error of 5% with the measured data, but the accuracy was improved by 20–30% compared with the traditional model solution, so the error was within the controllable range. The reason for the error was that there was some error in the selection of the soil pressure model in the deformation influence wavelength of the method in this paper. The determination of earth pressure parameters requires the improvement of existing earth pressure models. This is an important research component that will be further investigated in the future.

Figure 5 shows the relationship between the calculated parameters and excavation depth under different zones. From Figure 5, the potential energy increment and internal energy increment of the soil were positively correlated with the excavation depth. This positive correlation was due to the increase in the deformation influence range caused by the increase in excavation depth. The most sensitive of each energy increment to excavation depth was the potential energy increment in the ABCF zone, with the least significant change in the FEH zone. This finding was because excavation depth affected the area of each plastic deformation zone differently. The incremental change in internal energy for each zone was more similarly affected by the increase in excavation depth, and the shear strain in zone ABCF was the largest in all zones. However, considering that the soil shear strength was positively correlated with the burial depth, the internal energy increment in the ABCF zone was slightly smaller than that in the CEF zone.

4.2. Case 2

Tan et al. [45] studied the characterization of semi-top-down excavation for a subway station in Shanghai soft ground based on measured data. The excavation with a dimension of 139.6 × 21.8 × 24.8 m (length × width × depth) and the supporting structures consisting of the diaphragm wall (1 m wide) and mixed strut systems were adopted. The six internal support systems with the 1st concrete strut (0.8 × 1 m in section size and 6 m in horizontal spacing) and 2nd–6th steel pile struts (0.609 m in diameter, 2–3 m in horizontal spacing) were adopted. The soils around the excavation were mostly silty clay, muddy clay, clay and silty clay, and silty sand. For details of the excavation, refer to Tan et al. [45].

Figure 6 shows the comparison of measured and calculated wall displacement under different working conditions. From Figure 6, the wall displacement calculated by the proposed method was consistent with the measured values in the overall trend. The depth where the maximum wall deformation was basically the same as the measured values. When excavating to 12.3 m, 19.3 m, and 24.8 m, respectively, the average error between the measured and calculated wall displacement was 6.07%. This value verified the accuracy of the proposed calculation model in this paper. The main reasons for this error were (1) the influence of cyclic traffic loads on the excavation and (2) the temperature effect on excavation deformation. It is worth noting that there is some error between the measured and calculated wall bottom displacement. This outcome is due to the fact that the proposed model in this paper is mainly applied to special cases in which the basal strata is too soft, or the length of embedded retaining wall is insufficient.

4.3. Case 3

This paper relied on the subway excavation project in soft soil at Daliang Station (266 m long × 25.3 m deep × 19.9 m wide) reported by Han et al. [15]. The excavation support structure consisted of a diaphragm wall (0.8 m wide × 30 m deep) and four layers of hybrid internal support. The first layer was concrete support (0.7 × 0.9 m in cross-section size) and the second to fourth layers were steel supports (0.609 m in diameter, 0.016 m in thickness). The concrete and steel supports were spaced horizontally at 9 m and 4.5 m, respectively. The soils around the excavation were mainly silty clay, muddy soil, fully weathered muddy siltstone, and strongly weathered muddy siltstone. See Han et al. [15] for an overview of specific projects. The length of the diaphragm wall for subway excavation at Daliang Station was 30 m, and the embedment depth was 5 m. According to the geological condition and design data, the embedded depth of the excavation wall of Daliang station was insufficient, and the soil condition was soft soil, which belongs to the special excavation case to which the model of this paper applies.

Figure 7 shows a comparison of wall deformation when excavating to the third and fourth support positions. From Figure 7, the measured wall displacements were consistent with the overall trend of this paper’s solution, which demonstrates the reasonableness and accuracy of the proposed calculation model. As the excavation top support in the Foshan area was affected by the temperature, it led to the wall top deformation towards the outside of the excavation. And the proposed method in this paper does not consider the influence of temperature effect, so there is some error in the displacement value at the top of the wall. Therefore, the model proposed in this paper applies to the special case where the embedment depth of the wall of Daliang Station is insufficient and the soil condition situation is soft soil. In addition, the method proposed by Liu et al. [41] underestimated the wall displacement, while the method proposed by Lam and Bolton [21] overestimated the wall displacement, especially the large difference in displacement at the top and bottom of the wall. This finding was because these two methods proposed by Liu et al. [41] and Lam and Bolton [21] did not consider the actual boundary conditions of the excavation.

Figure 8 and

Table 1. Comparison of wall displacement.

Depth (m)	Wall Displacement When Excavating to 14.55 m (mm)				Depth (m)	Wall Displacement When Excavating to 19.7 m (mm)
	Measured Values	This Paper	Liu et al. [41]	Lam and Bolton [21]		Measured Values	This Paper	Liu et al. [41]	Lam and Bolton [21]
0	−3.21	0.15	0	0	0	−4.31	0	0.12	0
2	0.25	1.32	0.86	1.83	2	0	2	3.26	2
2.5	1.2	2.05	1.54	3	2.5	1.5	2.5	4.32	2.5
3.5	2.6	3.25	2	4	3.5	3.5	3.5	6.57	3.5
4.5	4.36	4.52	3.4	5.4	4.5	4.9	4.5	8.38	4.5
5.5	6.68	6.82	5.72	7.7	5.5	6	5.5	10.51	5.5
6.5	8.34	8.51	7.1	9.42	6.5	7.72	6.5	12.58	6.5
7.5	9.44	9.38	8.1	10.2	7.5	11.12	7.5	14.88	7.5
8.5	9.88	10.15	8.6	11.05	8.5	14.94	8.5	18.2	8.5
9.5	9.87	10.33	8.9	11.34	9.5	18.89	9.5	20.08	9.5
10.5	9.34	9.53	8	10.48	10.5	22.69	10.5	23.6	10.5
11.5	7.85	8.33	6.9	9.41	11.5	25.85	11.5	26.8	11.5
12	7.32	7.56	6	8.42	12	27.69	12	28.3	12
13.5	5.34	6.27	4.8	7.18	13.5	31.13	13.5	31.8	13.5
14	5.01	5.7	3.8	6.8	14.5	32.61	14.5	33.07	14.5
15.5	4.2	5.2	3.4	5.7	15.5	33.27	15.5	33.8	15.5
16	4	4.9	3	5.2	16	33.31	16	33.5	16
18.5	3.3	4.1	2.4	4.6	18.5	29.64	18.5	29.87	18.5
21	2.8	3.8	2	3.7	21	21.8	21	23.8	21
23	2.4	3.4	1.6	3.2	23	18.5	23	19.5	23
25	2	3.1	1.2	2.2	25	16.2	25	14.3	25
26.5	2	2.9	0.8	1.5	26.5	14.8	26.5	12.5	26.5
27.5	2	2.8	0.5	1	27.5	14	27.5	11.5	27.5
28	1.9	2.7	0.1	0.4	28	13.5	28	10.8	28
28.5	1.9	2.7	0	0	28.5	13	28.5	9.5	28.5

show the comparison of measured and predicted wall displacements. From Figure 8 and Table 1, the IMSD solution proposed in this paper can predict the wall displacements better and can improve the prediction accuracy of the MSD method, with the error basically in the range of 10–15%. Compared with the measured wall displacements, the method proposed by Liu et al. [41] underestimated the wall displacements, while the method proposed by Lam and Bolton [21] overestimated the wall displacements, with the error basically within 25%. When the wall displacement was small, the error between the predicted and measured values was relatively large. This result was because the methods proposed by Liu et al. [41] and Lam and Bolton [21] both ignored the effect of the actual boundary conditions of the excavation.

5. Parametric Analysis

The accuracy of the IMSD solution relies on the selection of parameters. Parameters affecting the IMSD solution can be divided into two categories: internal and external parameters. Various parameter combinations affect the deformation energy to different degrees. Therefore, the selection of parameter combinations, such as the size of ${\alpha }_1$ and ${\alpha }_2$, and the determination of the deformation wavelength may present some challenges and uncertainties. The internal parameters include two parameters used in the deformation mechanism: the coefficient of the virtual fixed point location, ${\alpha }_1$, and the influence range coefficient of excavation depth, ${\alpha }_2$. The external parameters include parameters such as the soil internal friction angle, $\phi$, the soil shear strength factor, $\xi$, the flexural stiffness of pile or wall, EI, and the compressive stiffness of support, EA. The internal parameters mainly change the deformation function form, and the external factors affect the maximum displacement increment of the MSD solution.

5.1. Effect of Virtual Fixed-Point Position Coefficients

Figure 9 shows the excavation displacements under different ${\alpha }_1$. When ${\alpha }_1=0$, the excavation top displacement was zero, and the maximum displacement was near the excavation surface. As ${\alpha }_1$ increased, the maximum excavation displacement and pile top displacement grew gradually, and the depth of the maximum excavation displacement shifted upward. When ${\alpha }_1$ was small, the depth of the maximum pile displacement was positively correlated with the excavation depth. When ${\alpha }_1$ increased to a certain extent, the depth of the maximum pile displacement gradually shifted upward with the excavation depth. When the maximum pile displacement occurred at the excavation top, it corresponded to the deformation pattern of the cantilevered excavation. Therefore, different values of ${\alpha }_1$ can reflect the excavation deformation under different support conditions. A reasonable value of ${\alpha }_1$ can be selected to calculate the MSD solution for different excavation types.

Figure 10 shows the relationship between excavation depth and different increments under different ${\alpha }_1$. From Figure 10, ${\alpha }_1$ was positively correlated with the strength mobilization coefficient and maximum displacement increment. At the same excavation depth, the potential energy increment of the ABCF zone decreased and the maximum displacement increment gradually increased with increasing ${\alpha }_1$. This result was because ${\alpha }_1$ stretched the wavelength of the deformation function by ${\alpha }_{1}h$ such that the actual integral area of the deformation function decreased. As ${\alpha }_1$ increased, the more integral area was lost, and the internal energy incremental decreased more than the potential energy incremental. When the excavation depth was certain, compared to ${\alpha }_1=0$, the potential energy increment was reduced by 3%, 6.3%, and 10.6% when ${\alpha }_1=0.2,0.4,0.6$, while the internal energy increment was reduced by 8.9%, 16.7%, and 19.2%, respectively.

5.2. Effect of the Influence Range Coefficient of Excavation Depth

Four sets of calculations with ${\alpha }_2$ equal to 1, k (=1–1.5), 1.5, and 2, respectively, were taken to solve for the pile deformation by the IMSD method, where ${\alpha }_2=k$ referred to the influence range determined by the method of calculating embedment depth. Similarly, the effect of the retaining structure and internal support on the deformation was not considered during the calculation of the MSD solution.

Figure 11 shows pile deformation under different ${\alpha }_2$. From Figure 11, the increase of ${\alpha }_2$ not only made the pile deformation pattern change but also significantly affected the maximum pile deformation. Contrary to the effect of ${\alpha }_1$, ${\alpha }_2$ caused a reduction in the depth of the maximum pile deformation. The pile bottom deformation, which was positively correlated with ${\alpha }_2$, was zero for the pile with good embedment depth (${\alpha }_2 < k$). When the excavation influence area exceeds the pile depth, the pile bottom may undergo a “kicking” deformation pattern with ${\alpha }_2=2$. When ${\alpha }_2$ doubled, the maximum pile deformation was about a quarter of the original pile deformation. The depth of maximum pile deformation decreased by about $0.5h$, and the deformation influence range increased dramatically. In addition, ${\alpha }_2$ can reflect the excavation deformation under different soil conditions. When ${\alpha }_2$ was small, the deformation pattern was consistent with the excavation under better soil conditions. As ${\alpha }_2$ increased, the deformation pattern was more similar to the excavation deformation characteristics in soft soil areas [5,19].

Figure 12 shows the variation in energy increments with ${\alpha }_2$ in different regions. From Figure 12, as ${\alpha }_2$ kept increasing, the potential energy increment and internal energy increment in each zone increased, while the maximum displacement increment gradually became smaller. This change was because ${\alpha }_2$ significantly increased the negative potential energy increment in the EFI zone. The increase in the sum of the potential energy increments in each zone was less than that of the internal energy increments. According to the formula ${\displaystyle {\int }_A{\gamma }_{sat}}△{\delta }_{\mathrm{v}}dA={\displaystyle {\int }_A\beta s_u\delta \gamma }dA$, the mobilized soil strength coefficient, $\beta$, decreased with the increase of ${\alpha }_2$, which resulted in the consequent decrease in the maximum pile displacement. The internal energy increment in the ABCD, CEF, FEH, and EFI regions when ${\alpha }_2=2$ was 1.3, 3, 3.1, and 3.3 times that of the ${\alpha }_2=1$, respectively.

5.3. Effect of Pile Bending Stiffness and Internal Support Stiffness

To analyze the role of ${\alpha }_1$ and ${\alpha }_2$ in the deformation mechanism, the effect of the retaining structure and internal support on the deformation energy was neglected. However, the deformation energy stored by retaining structure and internal support cannot be ignored in actual excavation engineering. The improved MSD method using an iterative method to calculate the maximum incremental deformation is a necessary improvement.

Figure 13 shows the relationship between the maximum displacement increment and excavation depth under different ${\alpha }_1$ and ${\alpha }_2$. From Figure 13a, the maximum displacement increment was inversely proportional to ${\alpha }_1$. This result was because ${\alpha }_1$ increased the pile top deformation and decreased the deformation influence wavelength. This action increased the incremental deformation energy of both the supporting structure and internal support, resulting in a decrease in the soil strength mobilization factor, which in turn suppressed the maximum incremental deformation in the MSD solution. From Figure 13b, the maximum displacement increment decreased only when ${\alpha }_2$ was small. This change was because ${\alpha }_2$ caused the deformation function value at the internal support position to decrease substantially, and the decrease degree was greater as ${\alpha }_2$ increased. Therefore, the bending deformation energy incremental of the retaining structure and the compression deformation energy incremental of the internal support were positively and negatively correlated with ${\alpha }_1$ and ${\alpha }_2$, respectively.

Figure 14 shows the effect of different stiffness of the retaining structure and internal support on excavation deformation, where $k_u$ and $k_v$ refer to the bending stiffness discount factor of the pile and support, respectively. From Figure 14, the compressive stiffness of the internal support inhibited the excavation deformation more than the pile stiffness. When the excavation depth was small, the pile bending stiffness could play a better inhibiting deformation effect more than the compressive stiffness of internal support, and the pattern was reversed for larger excavation depths. This action was because the bending deformation energy incremental depended on the pile deformation curvature, and the compression deformation energy incremental of the internal support depended on the displacement at the support.

5.4. Effect of Soil Friction Angle

In the modified MSD method, the angle of internal friction not only affects the plastic deformation mechanism and soil shear strength but also affects the deformation-influenced wavelength. Three sets of calculations with $\phi ={20}^{°},{24}^{°},{28}^{°}$ were taken to solve for pile deformation by the IMSD method.

Figure 15 shows the pile deformation under different internal friction angles. From Figure 15, pile deformation was inversely correlated with the soil friction angle. Pile deformation increased nonlinearly with decreasing internal friction angle, which was consistent with Liu et al. [41]. Relevant measured results [46] show that the maximum pile displacement was about 0.25H above the excavation surface when the soil conditions were good. While the soil condition was poor, the maximum pile displacement was near the excavation surface. When $\phi ={20}^{°}$, the deformation influence range was more than the embedment depth of the pile. Horizontal displacements were observed at the pile bottom, which was a typical deformation pattern when the embedment depth was insufficient.

5.5. Effect of Soil Shear Strength Coefficients

The soil shear strength coefficients play an important role in calculating the internal energy increment in each plastic deformation zone. In this paper, an undrained shear strength formula considering the soil anisotropy was used. Taking $K_0=0.5$, ${\phi }^{\prime }=33°$, $C_{\mathrm{c}}=0.33$, and $C_{\mathrm{s}}=0.027$, the ratio of the undrained shear strength to the overlying effective stress for different angles, ξ, can be obtained (see

Table 2. The ratio of undrained shear strength to overlying effective stress for different angles.

$\mathit{\theta }$	0°	30°	45°	60°	90°	Average Value
$S_{u\alpha }/{\sigma }_v$	0.255	0.313	0.300	0.262	0.153	0.257

). The effect of four sets of calculations with soil shear strength coefficients taken as k, 0.2, 0.25, and 0.3 on the excavation MSD solution was investigated.

Figure 16 shows the pile deformation under different soil shear strength coefficients. From Figure 16, the pile deformation was inversely proportional to ξ. The undrained soil shear strength formula used in this paper yielded pile deformation very close to ξ = 0.25. This result was because the ξ of this example was closer to the average value of ξ. The change in the soil shear strength coefficients directly affected the internal energy increment in the MSD method. A smaller ξ meant that the soil shear strength was lower; the greater the mobilization coefficient produced by exerting the same shear strength. The soil body needed to produce large deformations for the deformation mechanism to reach equilibrium.

6. Discussion

The deformation mechanisms proposed by Lam and Bolton [21] and Liu et al. [41] suffered from the problems of uncoordinated displacements and the inability to consider the effects of soil parameters. And the incremental deformation function ignored the compressive deformation energy of the first support. Therefore, based on the studies of Lam and Bolton [21] and Liu et al. [41], this paper proposes an improved plastic deformation field that can consider the effect of soil friction angle. Meanwhile, two parameters, the virtual fixed point location coefficient and the influence range coefficient of excavation depth, were introduced into the incremental deformation function, so as to reflect the influence of excavation depth and the first support on the deformation energy. Compared with the models proposed by Lam and Bolton [21] and Liu et al. [41], the accuracy of the proposed model in this paper was improved by 20–30%. This result indicated that the improved plastic deformation mechanism and incremental deformation function in this paper better reflected the actual incremental deformation of the excavation.

The model proposed in this paper is suitable for use in excavation projects with multiple internal supports combined with flexible retaining walls under stratified soil conditions. However, the proposed model does not consider the influence of complex factors such as overloading. And the selection of the soil pressure model in the deformation influence wavelength leads to a certain error in the calculated value, which needs further research.

7. Limitation

Massive field instrumentation data have already demonstrated that basal heave within an excavation enclosure results largely from the rebound of basal strata below excavation levels rather than the flow of ground behind the retaining wall into the excavation, the ground settlement behind the retaining wall derived from lateral wall movement rather than contributing to the basal heave [47,48]. For most cases in practice, the soft basal strata would be improved by ground treatment prior to soil removal, and the retaining wall would be designed with an adequate embedment length [49,50,51]. The ground movement pattern assumed in the MSD model (Figure 1) only works for a special condition, i.e., the basal strata is too soft or the length of embedded retaining wall is insufficient to constrain a kicking-out failure of the retaining wall, and then the ground behind the wall would slide into the excavation pushing the basal strata upwards. Therefore, further research and improvement of the ground movement pattern assumed in the MSD model (Figure 1) is needed in the future.

8. Conclusions

This paper mainly proposes a multi-support excavation deformation analysis method based on the improved MSD method. The main conclusions follow:

(1) Based on the traditional MSD method, an improved mobilized strength design (IMSD) method for multi-support excavation deformation analysis is proposed. The proposed model compensates for the shortcomings of the traditional MSD method; the improved incremental deformation mechanism further considers the effect of the soil friction angle, and the effect of excavation depth and the first support on deformation energy is also considered. The proposed IMSD method effectively improves the accuracy and reliability of soil deformation predictions by establishing virtual fixed points and excavation depth influence coefficients.

(2) Through an engineering case study, we found that the traditional MSD method underestimates the excavation deformation by 25% and the surface settlement by 13% compared to the measured results. The accuracy of the improved MSD solution has increased by 20–30%, which is more consistent with the measured data. This finding indicates that the IMSD method, considering the effect of soil parameters and boundary conditions, is more reasonable.

(3) The effect of internal and external parameters on the excavation deformation is discussed. The influence degree of internal and external parameters on excavation deformation is 30–40% and 10–20%, respectively. The virtual fixation point location coefficient and the soil friction angle have the greatest influence degree on the excavation deformation. The ${\alpha }_1$ mainly affects the deformation pattern, while ${\alpha }_2$ and the stiffness of the retaining structure mainly affects the maximum displacement increment.

(4) Potential applications of this method in different types of underground projects (such as subway construction and high-rise building foundations) can be explored, along with how to address challenges posed by complex environments. Furthermore, optimizing parameter selection, integrating new technologies (like artificial intelligence or machine learning) to enhance the intelligence level of the model, or developing adaptive models to respond to changes in different geological conditions will be studied in the future.