Mechanical Analysis of Frozen Wall with Trapezoidal Temperature Field Distribution Based on Unified Strength Theory

Abstract

In order to study the elastic–plastic stress field distribution of a double-row-pipe frozen wall, the temperature field of the double-row-pipe frozen wall is equivalent to a trapezoidal distribution, and the frozen wall is regarded as an elastic–plastic thick-walled cylinder with functionally gradient material (FGM) characteristics in the radial direction. Considering that the elastic modulus and cohesion of the frozen wall material change linearly with the radius, the elastic–plastic analysis of the frozen wall is carried out based on unified strength theory. The analytical solutions of the elastic–plastic stress field distribution, the elastic ultimate bearing capacity, the plastic ultimate bearing capacity, and the relative radius of the plastic zone of the frozen wall are derived. The analytical solution is calculated based on the engineering case and compared with the numerical solution obtained based on COMSOL. At the same time, the influence of strength theory parameters on the mechanical properties of heterogeneous and homogeneous frozen walls is analyzed. The results show that the analytical solution and the numerical solution are in good agreement, and their accuracy is mutually verified. The external load on the frozen wall of the selected layer is greater than its elastic ultimate bearing capacity and less than its plastic ultimate bearing capacity, which indicates that the frozen wall is in a safe state of stress. The radial stress increases with the increase in the strength theoretical parameter b and the relative radius r, the tangential stress increases with the increase in the strength theoretical parameter b, and first increases and then decreases with the increase in the relative radius r. The larger the strength theoretical parameter b, the smaller the relative radius of the plastic zone of the frozen wall. The strength theoretical parameter b increases from 0 to 1, the elastic ultimate bearing capacity and plastic ultimate bearing capacity of the heterogeneous frozen wall increase by 33.3% and 40.8%, respectively, and the elastic ultimate bearing capacity and plastic ultimate bearing capacity of the homogeneous frozen wall increase by 33.3% and 41.0%, respectively. Therefore, considering the influence of intermediate principal stress, the potential of materials can be fully exerted and the ultimate bearing capacity of frozen walls can be improved. This study can provide theoretical reference for the design and construction of frozen wall.

1. Introduction

The artificial ground freezing method was first proposed by the German engineer Borzsch in 1883 [1,2]. Its principle is to circulate refrigerant in the freezing pipe driven into the stratum to freeze the water in the stratum, turn the natural rock and soil into frozen soil, enhance its strength and stability, and isolate the connection between groundwater and underground work, so as to carry out safe construction under the protection of the freezing wall [3,4,5]. The artificial ground freezing method is applicable to the construction of underground work in water rich areas. This construction method has been widely used in the construction of new coal mine shafts [6,7] and damaged shafts [8], the construction of highway tunnels in coastal cities [9,10], the reinforcement of urban subway tunnel shield tunnel entrances and exits [11], the construction of subway tunnel connecting passages [12,13,14], foundation pit works [15], underground restricted natural gas storage tanks [16], etc.

With the increasing demand for energy, China’s shallow coal resources are gradually being exhausted, and people are constantly mining coal resources in deep strata. The vertical shafts of coal mines are getting deeper and deeper, often crossing multiple aquifers. The freezing method is often used in the process of coal mine shaft construction. On the premise of meeting the requirements of isolated groundwater, the strength and stability of the frozen wall are also decisive factors for the safety of shaft excavation. Therefore, it is particularly important to study the mechanical characteristics of the frozen wall.

In elastic–plastic theoretical research of frozen walls, the frozen wall is basically regarded as a homogeneous or radial heterogeneous thick-walled cylinder. In terms of homogeneous frozen wall research, Wang et al. [17] carried out a mechanical theoretical analysis on the unlined frozen wall based on the Moore–Coulomb strength criterion. Yang et al. [18,19,20] considered the interaction between a frozen wall and the surrounding rock, analyzed the frozen wall in the elastic state, plastic state, and elastic–plastic state, respectively, in terms of mechanical theory, and gave analytical formulas for the elastic critical thickness and plastic critical thickness of frozen walls. Zhang et al. [21] conducted theoretical research on the frozen wall based on the unified strength theory and obtained analytical solutions of the stress field and displacement field. Hu [22,23] gave the mechanical model of a frozen wall considering excavation unloading and interaction between the frozen wall and surrounding soil, as well as the expression of external load of the frozen wall in this case. However, the frozen wall is essentially a functionally graded material with temperature gradient, and its mechanical properties and material parameters will also vary with temperature. Hu et al. [24,25,26,27] simplified the temperature field of a double-row-pipe and triple-row-pipe into parabola, equivalent trapezoid and trapezoid parabola superposition models, respectively. The simplified models were verified, which played a huge role in promoting the research on mechanical properties of heterogeneous frozen walls. In terms of heterogeneous frozen wall research, Hu et al. [28,29] regarded a frozen wall of double-row-pipe as a heterogeneous frozen wall with a function gradient along the radial direction, and used the Moore–Coulomb strength criterion to conduct an elastoplastic analysis of the frozen wall. Rong et al. [30] regarded the properties of the frozen wall materials as parabola changes along the radial direction, and used Drucker–Prager strength criterion to conduct elastic–plastic theoretical analysis of the frozen wall. Wang et al. [31,32] regarded the properties of the frozen wall material as parabola or linear change along the radial direction, and deduced the analytical solution of the elastoplastic stress of the frozen wall by considering the interaction between the frozen wall and the surrounding soil. Cao et al. [33] regarded the material properties of the frozen wall as a parabola along the radial direction of the frozen wall, and conducted an elastoplastic analysis of the frozen wall based on the unified strength theory.

In this paper, the temperature field of a double-row-pipe frozen wall is simplified as a trapezoidal distribution along the radial direction, and the linear relationship between the elastic modulus E and the cohesive force c of the freezing wall material and temperature is considered. Taking Qingdong Mine as the engineering background, the elastic–plastic analysis of the frozen wall of the coal mine is carried out using unified strength theory, and the analytical solution of the stress field is obtained. The numerical simulation of the temperature field of the frozen wall is carried out based on COMSOL software, and the numerical solution of the stress field is obtained. The results obtained by the two methods are mutually verified, and the correctness of the two methods is obtained. This study evaluates whether the strength and stability of the frozen wall of the coal mine meet the requirements. At the same time, the influence of the strength theoretical parameter b on the stress distribution and stress characteristics of the heterogeneous and homogeneous frozen wall is analyzed. This study enriches the elastic–plastic theoretical analysis of frozen walls. Compared with the numerical simulation, it greatly reduces the calculation time and cost and can provide a theoretical reference for the design and construction of frozen walls. The research flow chart is shown in Figure 1.

2. Equivalent Temperature Field

According to existing research, in order to more accurately describe the temperature field of the double-row-pipe frozen wall, the temperature field on a typical section of the frozen wall (the section at 1/4 of the pipe spacing) is selected to replace the distribution of the overall temperature field. The temperature field on the typical section is equivalent to a trapezoidal distribution along the radial direction of the freezing wall. The equivalent temperature field model is shown in Figure 2 [24]. The temperature of any point in the frozen wall can be calculated by Equation (1).

$$T(R)=\{\begin{array}{cc}{T}_n+\frac{T_k-T_n}{R_1-R_0}(R-R_0)& R_0\le R\le R_1\\ T_k& R_1\le R\le R_2\\ T_n+\frac{T_k-T_n}{R_1-R_0}(R-R_0)& R_2\le R\le R_e\end{array}$$

(1)

where T₀ is the freezing temperature of soil mass, °C; T_n is the temperature of the inner wall of the frozen wall after excavation, °C; T_k is the lowest temperature of trapezoidal equivalent temperature field, °C; R₀, R₁, R₂, R_e, and R are the cylindrical coordinates of the wellbore centerline as the origin, m; R₀ is the radius of the inner frozen wall, m; R₁ is the radius of the inner frozen pipe, m; R₂ is the radius of the outer frozen pipe, m; R_e is the radius of the outer frozen wall, m; and R is the radius of any point in the frozen wall, m.

3. Mechanical Model of a Frozen Wall

The mechanical model of a frozen wall selected for elastic–plastic analysis is shown in Figure 3, where P is the external load on the frozen wall, and $R_{\rho }$ is the radius of the plastic zone of the frozen wall. The elastic modulus E(R) and cohesion C(R) are directly proportional to the temperature within the scope of the engineering application.

$$\{\begin{array}{l}E\left(R\right)=aT\left(R\right)+E_0\\ C\left(R\right)=lT\left(R\right)+C_0\end{array}$$

(2)

where a and l are the coefficients of linear Equation (2); E₀ is the elastic modulus of the inner frozen wall; and C₀ is the cohesion of the inner frozen wall.

Substituting Equation (1) into E(R) and C(R), and making the relative radius r = R/R₀, then relative radius of inner frozen pipe r₁ = R₁/R₀, relative radius of outer frozen pipe r₂ = R₂/R₀, relative radius of outer frozen wall r_e = R_e/R₀, we get [29]:

$$E(r)=\{\begin{array}{cc}{a}_{1}r+b_1& 1\le r\le r_1\\ E_1& r_1\le r\le r_2\\ a_{2}r+b_2& r_2\le r\le r_e\end{array}$$

(3)

where a₁, b₁, a₂, and b₂ are the coefficients of Equation (3); E₁ is the elastic modulus of the inner frozen pipe:

$$\begin{gathered} a_1=a\frac{T_k-T_n}{r_1-1};b_1=a\left(T_n-\frac{T_k-T_n}{r_1-1}\right)+E_0;E_1=a{T}_k+E_0; \\ {a}_2=a\frac{T_k-T_0}{r_2-r_e}; \\ {b}_2=a\left(T_0-\frac{T_k-T_0}{r_2-r_e}{r}_e\right)+E_0.\hspace{1em}C(r)=\{\begin{array}{cc}{l}_{1}r+m_1& 1\le r\le r_1\\ C_1& r_1\le r\le r_2\\ l_{2}r+m_2& r_2\le r\le r_e\end{array} \end{gathered}$$

(4)

where l₁, m₁, l₂, and m₂ are the coefficients of Equation (4) [29]; C₁ is the cohesion of the inner frozen pipe; $l_1=l\frac{T_k-T_n}{r_1-1}$; $m_1=l\left(T_n-\frac{T_k-T_n}{r_1-1}\right)+C_0$;

$$C_1=l{T}_k+C_0;l_2=l\frac{T_k-T_0}{r_2-r_e};m_2=l\left(T_0-\frac{T_k-T_0}{r_2-r_e}{r}_e\right)+C_0$$

Assuming the same development rate on both sides of the frozen wall, then a₂ = −a₁.

4. Model Solving

4.1. Solution of Stress Field in the Elastic State

Under the elastic state, according to the theoretical solution of a thick-walled cylinder in elasticity [34], the analytical solution of the frozen wall stress is as follows:

$$\{\begin{array}{cc}{\sigma }_r=\frac{{\sigma }_{r_1}[a_1\left(1-\frac{1}{r}\right)+\frac{b_1}{2}\left(1-\frac{1}{r^2}\right)]}{\left(1-\frac{1}{r_1}\right)a_1+\left(1-\frac{1}{r^2}\right)\frac{b_1}{2}}& \\ {\sigma }_{\theta }=\frac{{\sigma }_{r_1}[a_1+\frac{b_1}{2}\left(1+\frac{1}{r^2}\right)]}{\left(1-\frac{1}{r_1}\right)a_1+\left(1-\frac{1}{r^2}\right)\frac{b_1}{2}}& \end{array}r\in [1,r_1]$$

(5)

$$\{\begin{array}{cc}{\sigma }_r=\frac{1}{r_2^2-r_1^2}\left[-\frac{r_1^2{r}_2^2({\sigma }_{r_2}-{\sigma }_{r_1})}{r^2}+\left(r_2^2{\sigma }_{r_2}-r_1^2{\sigma }_{r_1}\right)\right]& \\ {\sigma }_{\theta }=\frac{1}{r_2^2-r_1^2}\left[\frac{r_1^2{r}_2^2({\sigma }_{r_2}-{\sigma }_{r_1})}{r^2}+\left(r_2^2{\sigma }_{r_2}-r_1^2{\sigma }_{r_1}\right)\right]& \end{array}r\in [r_1,r_2]$$

(6)

$$\{\begin{array}{cc}{\sigma }_r=\frac{(P-{\sigma }_{r_2})\left(\frac{a_1}{r}-\frac{b_2}{2r^2}\right)-a_1\left(\frac{P}{r_2}-\frac{{\sigma }_{r_2}}{r_e}\right)+\frac{b_2}{2}\left(\frac{P}{r_2{}^2}-\frac{{\sigma }_{r_2}}{r_e{}^2}\right)}{\left(\frac{1}{r_e}-\frac{1}{r_2}\right)a_1+\frac{b_2}{2}\left(\frac{1}{r_2{}^2}-\frac{1}{r_e{}^2}\right)}& \\ {\sigma }_{\theta }=\frac{\frac{b_2}{2r^2}(P-{\sigma }_{r_2})-a_1\left(\frac{P}{r_2}-\frac{{\sigma }_{r_2}}{r_e}\right)+\frac{b_2}{2}\left(\frac{P}{r_2{}^2}-\frac{{\sigma }_{r_2}}{r_e{}^2}\right)}{\left(\frac{1}{r_e}-\frac{1}{r_2}\right)a_1+\frac{b_2}{2}\left(\frac{1}{r_2{}^2}-\frac{1}{r_e{}^2}\right)}& \end{array}r\in [r_2,r_e]$$

(7)

where ${\sigma }_r$ is the radial stress, MPa; ${\sigma }_{\theta }$ is the tangential stress, MPa; P is the external load; ${\sigma }_{r_1}$ is the radial stress of the inner frozen pipe, MPa; ${\sigma }_{r_2}$ is the radial stress of the outer frozen pipe, MPa; and P is the external load, MPa.

According to the displacement continuity condition, we solve ${\sigma }_{r_1}$ and ${\sigma }_{r_2}$, where $u_1^{-}$ is the displacement of the inner boundary when the relative radius r = r₁, m; $u_1^{+}$ is the displacement of the outer boundary when the relative radius r = r₁, m; $u_2^{-}$ is the displacement of the inner boundary when the relative radius r = r₂, m; and $u_2^{+}$ is the displacement of the outer boundary when the relative radius r = r₂, m.

According to $u_1^{-}=u_1^{+}$, we get

$$-\mu {\sigma }_{r_1}+(1-\mu )\frac{{\sigma }_{r_1}[a_1+\frac{b_1}{2}\left(1+\frac{1}{r^2}\right)]}{\left(1-\frac{1}{r_1}\right)a_1+\left(1-\frac{1}{r^2}\right)\frac{b_1}{2}}=\frac{(2-2\mu )r_2^2}{r_2^2-r_1^2}{\sigma }_{r_2}-\frac{(1-2\mu )r_1^2+r_2^2}{r_2^2-r_1^2}{\sigma }_{r_1}$$

(8)

where μ is Poisson’s ratio.

According to $u_2^{-}=u_2^{+}$, we get:

$$\begin{array}{l}\frac{r_1^2+(1-2\mu )r_2^2}{r_2^2-r_1^2}{\sigma }_{r_2}-\frac{(2-2\mu )r_2^2}{r_2^2-r_1^2}{\sigma }_{r_1}=\frac{-\mu \left[(P-{\sigma }_{r_2})\left(\frac{a_1}{r_2}-\frac{b_2}{2r_2^2}\right)-a_1\left(\frac{P}{r_2}-\frac{{\sigma }_{r_2}}{r_e}\right)+\frac{b_2}{2}\left(\frac{P}{r_2{}^2}-\frac{{\sigma }_{r_2}}{r_e{}^2}\right)\right]}{\left(\frac{1}{r_e}-\frac{1}{r_2}\right)a_1+\frac{b_2}{2}\left(\frac{1}{r_2{}^2}-\frac{1}{r_e{}^2}\right)}\\ +\frac{(1-\mu )\left[\frac{b_2}{2r_2{}^2}(P-{\sigma }_{r_2})-a_1\left(\frac{P}{r_2}-\frac{{\sigma }_{r_2}}{r_e}\right)+\frac{b_2}{2}\left(\frac{P}{r_2{}^2}-\frac{{\sigma }_{r_2}}{r_e{}^2}\right)\right]}{\left(\frac{1}{r_e}-\frac{1}{r_2}\right)a_1+\frac{b_2}{2}\left(\frac{1}{r_2{}^2}-\frac{1}{r_e{}^2}\right)}\end{array}$$

(9)

By simultaneous Equations (8) and (9), we get ${\sigma }_{r_1}$ and ${\sigma }_{r_2}$.

4.2. Elastic–Plastic Analysis of Frozen Wall

4.2.1. Solution of Stress Field in the Plastic Zone

The elasto-plastic analysis of a frozen wall is based on unified strength theory, which has different expressions. In geotechnical engineering, the cohesive force C and the angle of internal friction $\phi$ of materials are generally used to express the unified strength [35].

$$\{\begin{array}{c}F={\sigma }_1-\frac{1+\sin \phi }{(1+b)(1-\sin \phi )}\left(b{\sigma }_2+{\sigma }_3\right)=\frac{2C_0\cos \phi }{1-\sin \phi }\\ {\sigma }_2\le \frac{1}{2}\left({\sigma }_1+{\sigma }_3\right)-\frac{\sin \phi }{2}\left({\sigma }_1-{\sigma }_3\right)\\ F^{\prime }=\frac{1}{1+b}\left({\sigma }_1+b{\sigma }_2\right)-\frac{1+\sin \phi }{1-\sin \phi }{\sigma }_3=\frac{2C_0\cos \phi }{1-\sin \phi }\\ {\sigma }_2\ge \frac{1}{2}\left({\sigma }_1+{\sigma }_3\right)-\frac{\sin \phi }{2}\left({\sigma }_1-{\sigma }_3\right)\end{array}$$

(10)

where F and F′ are strength theoretical functions, respectively; ${\sigma }_1$ is the first principal stress, the compressive stress is positive, MPa; ${\sigma }_2$ is the second principal stress, MPa; ${\sigma }_3$ is the third principal stress, MPa; b is the strength theory parameter, which represents the influence of the intermediate principal shear stress and the normal stress on the surface on the material damage. When b takes different values, it can represent or linearly approximate different strength criteria. When b = 0, it degenerates into the molar strength theory; when b = 1, it degenerates into the twin-shear strength criterion; when 0 < b < 1, it is a series of new strength criteria.

When the external load applied to the frozen wall exceeds the elastic limit load of the frozen wall, the frozen wall enters the elasto-plastic stress state. According to elasticity, under the plane strain state,

$${\sigma }_z=\frac{1}{2}\left({\sigma }_r+{\sigma }_{\theta }\right)$$

(11)

where ${\sigma }_z$ is axial stress.

It can be seen from ${\sigma }_3\le {\sigma }_2\le {\sigma }_1$,

$${\sigma }_1={\sigma }_{\theta },{\sigma }_2={\sigma }_z,{\sigma }_3={\sigma }_r$$

(12)

There is ${\sigma }_2\ge \frac{1}{2}\left({\sigma }_1+{\sigma }_3\right)-\frac{\sin \phi }{2}\left({\sigma }_1-{\sigma }_3\right)$, substituting Equations (11) and (12) into Equation (10), we get:

$${\sigma }_{\theta }^p-{\sigma }_r^p=\frac{4(1+b)\left[C(r)\cos \phi +{\sigma }_r\sin \phi \right]}{(1-\sin \phi )(2+b)}$$

(13)

where ${\sigma }_{\theta }^p$ is the tangential stress in the plastic zone; ${\sigma }_r^p$ is the radial stress in the plastic zone.

The equilibrium equation of the plastic zone is:

$$\frac{d{\sigma }_r}{dr}+\frac{{\sigma }_r-{\sigma }_{\theta }}{r}=0$$

(14)

Substituting Equation (13) into Equation (14), we get:

$$\frac{d{\sigma }_r}{dr}-\frac{4(1+b)\sin \phi }{(1-\sin \phi )(2+b)}\frac{{\sigma }_r}{r}=\frac{4(1+b)\cos \phi }{(1-\sin \phi )(2+b)}\frac{C(r)}{r}$$

(15)

The general solution of non-homogeneous linear differential Equation (15) is as follows:

$${\sigma }_r^p=C{r}^{\frac{4(1+b)\sin \phi }{(1-\sin \phi )(b+2)}}+\frac{4k_1(1+b)\cos \phi }{(1+b)(1-5\sin \phi )+(1-\sin \phi )}r-\frac{k_2\cos \phi }{\sin \phi }$$

(16)

The relative radius of the plastic zone is $r_{\rho }$. When $r_{\rho }\in [1,r_1]$, k₁ = l₁ and k₂ = m₁; when r = 1, ${\sigma }_r^p=0$, and substituting it into Equation (16), we get:

$$C=-\frac{4l_1(1+b)\cos \phi }{(1+b)(1-5\sin \phi )+(1-\sin \phi )}+\frac{m_1\cos \phi }{\sin \phi }$$

(17)

Substituting Equation (17) into Equation (16), we get Equation (A1) in the Appendix A for the stress in the plastic zone.

When $r_{\rho }\in [r_1,r_2]$, k₁ = 0, k₂ = C₁ = l₁r₁ + m₁; when r = r₁, ${\sigma }_r^p={\sigma }_{r_1}^p$, ${\sigma }_{r_1}^p$ is the radial stress in the elastic zone when r = r₁; substituting it into Equations (16) and (A1), we get Equation (A2).

Substituting Equation (A2) into Equation (16), we get Equation (A3) for the stress in the plastic zone.

When $r_{\rho }\in [r_1,r_e]$, k₁ = l₂ and k₂ = m₂; when r = r₂, ${\sigma }_r^p={\sigma }_{r_2}^p$, ${\sigma }_{r_2}^p$ is the radial stress in the elastic zone when r = r₂; substituting it into Equations (16) and (A3), we get Equation (A4).

Substituting Equation (A4) into Equation (16), we get Equation (A5) for the stress in the plastic zone.

The unified strength theory is shown in Equation (18); the tangential stress ${\sigma }_{\theta }^p$ of each plastic zone can be obtained from the following equation:

$${\sigma }_{\theta }^p-{\sigma }_r^p=\frac{4(1+b)\left[C(r)\cos \phi +{\sigma }_r^p\sin \phi \right]}{(1-\sin \phi )(2+b)}$$

(18)

4.2.2. Solution of Stress Field in the Elastic Region

When $r_{\rho }\in [1,r_1]$, it is an elastic thick-walled cylinder model with inner frozen wall stress ${\sigma }_{r_{\rho }}$ and outer frozen wall stress ${\sigma }_{r_1}$ in the zone $[r_{\rho },r_1]$; ${\sigma }_{r_{\rho }}$ is the radial stress in the elastic zone when r = r_ρ, The elastic solution is as follows:

$$\{\begin{array}{l}{\sigma }_r=\frac{\left({\sigma }_{r_1}-{\sigma }_{r_{\rho }}\right)\left(-\frac{a_1}{r}-\frac{b_1}{2r^2}\right)+a_1\left(\frac{{\sigma }_{r_1}}{r_{\rho }}-\frac{{\sigma }_{r_{\rho }}}{r_1}\right)+\frac{b_1}{2}\left(\frac{{\sigma }_{r_1}}{r_{\rho }^2}-\frac{{\sigma }_{r_{\rho }}}{r_1^2}\right)}{a_1\left(\frac{1}{r_{\rho }}-\frac{1}{r_1}\right)+\frac{b_1}{2}\left(\frac{1}{r_{\rho }^2}-\frac{1}{r_1^2}\right)}\\ {\sigma }_{\theta }=\frac{\frac{b_1}{2r^2}\left({\sigma }_{r_1}-{\sigma }_{r_{\rho }}\right)+a_1\left(\frac{{\sigma }_{r_1}}{r_{\rho }}-\frac{{\sigma }_{r_{\rho }}}{r_1}\right)+\frac{b_1}{2}\left(\frac{{\sigma }_{r_1}}{r_{\rho }^2}-\frac{{\sigma }_{r_{\rho }}}{r_1^2}\right)}{a_1\left(\frac{1}{r_{\rho }}-\frac{1}{r_1}\right)+\frac{b_1}{2}\left(\frac{1}{r_{\rho }^2}-\frac{1}{r_1^2}\right)}\end{array}r\in [r_{\rho },r_1]$$

(19)

The elastic solution in zone II and zone III is the same as Equations (6) and (7).

When $r_{\rho }\in [r_1,r_2]$, it is an elastic thick-walled cylinder model with inner frozen wall stress ${\sigma }_{r_{\rho }}$ outer frozen wall stress ${\sigma }_{r_2}$ in the zone $[r_{\rho },r_2]$. The elastic solution is as follows:

$$\{\begin{array}{l}{\sigma }_r=\frac{1}{r_2^2-r_{\rho }^2}\left[-\frac{r_2^2{r}_{\rho }^2\left({\sigma }_{r_2}-{\sigma }_{r_{\rho }}\right)}{r^2}+r_2^2{\sigma }_{r_2}-r_{\rho }^2{\sigma }_{r_{\rho }}\right]\\ {\sigma }_{\theta }=\frac{1}{r_2^2-r_{\rho }^2}\left[\frac{r_2^2{r}_{\rho }^2\left({\sigma }_{r_2}-{\sigma }_{r_{\rho }}\right)}{r^2}+r_2^2{\sigma }_{r_2}-r_{\rho }^2{\sigma }_{r_{\rho }}\right]\end{array}r\in [r_{\rho },r_2]$$

(20)

The elastic solution in zone III is the same as Equation (7).

When $r_{\rho }\in [r_2,r_e]$, it is an elastic thick-walled cylinder model with inner frozen wall stress ${\sigma }_{r_{\rho }}$ and outer frozen wall stress P in the zone $[r_{\rho },r_e]$. The elastic solution is as follows:

$$\{\begin{array}{l}{\sigma }_r=\frac{\left(P-{\sigma }_{r_{\rho }}\right)\left(\frac{a_1}{r}-\frac{b_2}{2r^2}\right)+a_1\left(\frac{{\sigma }_{r_{\rho }}}{r_e}-\frac{P}{r_{\rho }}\right)+\frac{b_2}{2}\left(\frac{P}{r_{\rho }^2}-\frac{{\sigma }_{r_{\rho }}}{r_e^2}\right)}{a_1\left(\frac{1}{r_e}-\frac{1}{r_{\rho }}\right)+\frac{b_2}{2}\left(\frac{1}{r_{\rho }^2}-\frac{1}{r_e^2}\right)}\\ {\sigma }_{\theta }=\frac{\frac{b_2}{2r^2}\left(P-{\sigma }_{r_{\rho }}\right)+a_1\left(\frac{{\sigma }_{r_{\rho }}}{r_e}-\frac{P}{r_{\rho }}\right)+\frac{b_2}{2}\left(\frac{P}{r_{\rho }^2}-\frac{{\sigma }_{r_{\rho }}}{r_e^2}\right)}{a_1\left(\frac{1}{r_e}-\frac{1}{r_{\rho }}\right)+\frac{b_2}{2}\left(\frac{1}{r_{\rho }^2}-\frac{1}{r_e^2}\right)}\end{array}r\in [r_{\rho },r_e]$$

(21)

4.2.3. Solution of Relative Radius $r_{\rho }$ in the Plastic Zone

The frozen wall begins to yield from the inner wall. After entering the plastic zone, the stress state of the frozen wall meets the unified strength criterion.

When $r_{\rho }\in [1,r_1]$ and $r=r_{\rho }$, the radial stress and tangential stress meet the unified strength theory and stress continuity conditions, as shown in Equations (A6) and (A7).

Here, ${\sigma }_{{\theta }_{r_{\rho }}}^p$ is the tangential stress in the plastic zone when $r=r_{\rho }$; ${\sigma }_{{\theta }_{r_{\rho }}}^e$ is the tangential stress in the elastic zone when $r=r_{\rho }$; ${\sigma }_{r_{\rho }}^e$ is the radial stress in the elastic zone when $r=r_{\rho }$.

The relative radius $r_{\rho }$ of the plastic zone can be obtained using simultaneous Equations (A6) and (A7). The solution of ${\sigma }_{r_1}$ in Equation (A6) refers to the solution of ${\sigma }_{r_1}$ in the elastic state using the displacement continuity condition.

When $r_{\rho }\in [r_1,r_2]$ and $r=r_{\rho }$, the radial stress and tangential stress meet the unified strength theory and stress continuity conditions, as shown in Equations (A8) and (A9).

The relative radius $r_{\rho }$ of the plastic zone can be obtained using simultaneous Equations (A8) and (A9). The solution of ${\sigma }_{r_2}$ in Equation (A8) refers to the solution of ${\sigma }_{r_2}$ in the elastic state using the displacement continuity condition.

When $r_{\rho }\in [r_2,r_e]$ and $r=r_{\rho }$, the radial stress and tangential stress meet the unified strength theory and stress continuity conditions, as shown in Equations (A10) and (A11).

The relative radius $r_{\rho }$ of the plastic zone can be obtained using simultaneous Equations (A10) and (A11).

4.3. Elastic Limit Load

When the frozen wall enters the plastic state, the plastic zone of the frozen wall begins to yield from the inner wall of the frozen wall. Therefore, when the inner wall of the frozen wall just enters the plastic stage ($r_{\rho }=1$), the frozen wall is in the elastic limit state. The external load borne by the frozen wall is called the elastic limit load $P_e$ of the frozen wall.

When the inner wall of the frozen wall just enters the plastic stage ($r_{\rho }=1$), the stress at the inner wall (r = 1) meets the unified strength theory and stress continuity conditions:

$${\sigma }_{{\theta }_{r_{\rho }}}^p-{\sigma }_{r_{\rho }}^p={\sigma }_{{\theta }_{r_{\rho }}}^e-{\sigma }_{r_{\rho }}^e=\frac{{\sigma }_{r_1}\left(a_1+b_1\right)}{a_1\left(1-\frac{1}{r_1}\right)+\left(1-\frac{1}{r_1^2}\right)\frac{b_1}{2}}=\frac{4(1+b)\left[\left(l_1+m_1\right)\cos \phi +{\sigma }_{r_{\rho }}^p\sin \phi \right]}{(1-\sin \phi )(2+b)}$$

(22)

We know that ${\sigma }_{r_{\rho }}^p=0$, and substituting it into Equation (21), we get:

$${\sigma }_{r_1}=\frac{4(1+b)\left(l_1+m_1\right)\cos \phi \left[a_1\left(1-\frac{1}{r_1}\right)+\frac{b_1}{2}\left(1-\frac{1}{r_1^2}\right)\right]}{(1-\sin \phi )(b+2)\left(a_1+b_1\right)}$$

(23)

According to the displacement continuity condition $u_1^{-}=u_1^{+}$ at r₁, we get:

$$-\mu {\sigma }_{r_1}+(1-\mu )\frac{{\sigma }_{r_1}\left[a_1+\frac{b_1}{2}\left(1+\frac{1}{r_1^2}\right)\right]}{a_1\left(1-\frac{1}{r_1}\right)+\frac{b_1}{2}\left(1-\frac{1}{r_1^2}\right)}=\frac{(2-2\mu )r_2^2}{r_2^2-r_1^2}{\sigma }_{r_2}-\frac{(1-2\mu )r_1^2+r_2^2}{r_2^2-r_1^2}{\sigma }_{r_1}$$

(24)

According to the displacement continuity condition $u_2^{-}=u_2^{+}$ at r₂, we get Equation (A12).

The elastic limit load $P_e$ of the frozen wall can be obtained using simultaneous Equations (23), (24), and (A12).

4.4. Plastic Limit Load

When the plastic radius of the frozen wall increases to the radius of the outer frozen wall ($r_{\rho }=r_e$), the load borne by the frozen wall is the plastic limit external load of the frozen wall. According to Equation (A11), the plastic limit load of the frozen wall can be obtained as shown in Equation (A13).

5. Stress Field Distribution of Frozen Wall and Numerical Verification: A Case Study

The Dongfeng shaft of Qingdong Mine in the Huaibei mining area is constructed by the artificial ground freezing method. The clay layer with a vertical depth of 225 m is frozen by double-row-pipe. The inner radius of the frozen wall is R_B = 4.053 m, the wall temperature is −9 °C, the outer radius is R_e = 8.147 m, and the corresponding temperature is −2.11 °C; r_e = R_e/R_B = 2.01. The minimum temperature of the trapezoidal equivalent temperature field is −22.32 °C. The material is soil taken from the frozen wall. The internal friction angle of soil is measured by a triaxial shear test, and the Poisson’s ratio of the soil is measured by a triaxial compression test. The internal friction angle of the frozen soil is 10.5°, and the Poisson’s ratio is 0.29. The temperature field distribution of the frozen wall and the material parameters obtained from the thermophysical and mechanical test report of the frozen soil are shown in Equations (25)–(27).

$$T(r)=\{\begin{array}{cc}-53.27r+44.27& r\in [1,1.25]\\ -22.32& r\in [1.25,1.63]\\ 53.27r-109.15& r\in [1.63,2.01]\end{array}$$

(25)

$$E(r)=\{\begin{array}{cc}359.04r-284.26& r\in [1,1.25]\\ 164.54& r\in [1.25,1.63]\\ -359.04r+749.78& r\in [1.63,2.01]\end{array}$$

(26)

$$C(r)=\{\begin{array}{cc}5.27r-4.51& r\in [1,1.25]\\ 2.08& r\in [1.25,1.63]\\ -5.27r+10.67& r\in [1.63,2.01]\end{array}$$

(27)

The external load of the frozen wall is generally calculated by the heavy liquid equation [36]:

P = 0.013h

(28)

where h is the depth, m.

The calculation shows that the external load borne by the frozen wall at the vertical depth of 225 m is 2.925 MPa. On the premise that b = 0.5, according to the equations derived above, it can be obtained that the relative plastic radius of the frozen wall under external load is $r_{\rho }=1.32$. The radial stress ${\sigma }_r$ and tangential stress ${\sigma }_{\theta }$ distribution of the frozen wall are shown in Figure 3. Figure 4 is drawn by Origin software. The radial stress ${\sigma }_r$ of the frozen wall increases with the increase in relative radius r and reaches the maximum at the outer frozen wall. The tangential stress ${\sigma }_{\theta }$ first increases and then decreases with the increase in the relative radius r. When $r\in [1,r_{\rho }]$, in the plastic zone, the tangential stress ${\sigma }_{\theta }$ of the frozen wall increases with the increase in the relative radius r, and the rate of increase in the tangential stress between the inner wall of the frozen wall and the inner frozen pipe is greater than the rate of increase between the inner frozen pipe and the radius of the plastic zone; the tangential stress reaches the maximum at the plastic radius. When $r\in [r_{\rho },r_e]$, in the elastic zone, the tangential stress of the frozen wall decreases with the increase in the relative radius, and the reduction rate of the tangential stress between the plastic zone radius and the outer frozen pipe is less than the reduction rate between the outer frozen pipe and the outer wall of the frozen wall. According to the equation derived above, when the strength theoretical parameter b = 0.5, the elastic limit load and plastic limit load of the frozen wall are 1.44 MPa and 3.55 MPa, respectively. The external load of the frozen wall is greater than its elastic ultimate bearing capacity, and the external load of the frozen wall is less than its plastic ultimate bearing capacity, indicating that the stress state of the frozen wall is safe.

In order to verify the rationality of the analytical solution, COMSOL software is used for numerical simulation. A numerical calculation model is established based on the Drucker–Prager strength criterion. At the same time, the numerical calculation model is simplified according to the geometric shape, load distribution, and constraint conditions, and finally a quarter symmetric numerical calculation model is adopted. The parameters of the numerical analysis are consistent with those of the analytical solution. The calculation model and grid division are shown in Figure 5. Through the analysis of numerical simulation software, the quality of the model grid cells is 1, and the quality of the grid division is high. The radial stress and tangential stress of the model are derived through the COMSOL post-processing function, as shown in Figure 4. The plastic radius of the numerical solution is consistent with that of the analytical solution. The results of the numerical solution are highly consistent with the results of the analytical solution. The change trend in the plastic zone and the elastic zone is consistent, and the difference between the values is very small. The results of the numerical solution and the analytical solution are mutually verified. However, the analytical solution uses the unified strength theory, and the numerical solution uses the Drucker–Prager strength criterion, so there is a small error.

Different values of b, the stress distribution of frozen wall under elastic limit state and plastic limit state, are shown in Figure 6 and Figure 7. The radial stress of the frozen wall increases with the increase in the relative radius r, and the maximum radial stress is located at the outer frozen wall. The tangential stress first increases and then decreases with the increase in the relative radius, but the change rule is different under different stress states.

When $r\in [1,r_1]$, the tangential stress increases with the increase in the relative radius r. When $r\in [r_1,r_2]$, the tangential stress decreases with the increase in r in the elastic limit state, and increases with the increase in r in the plastic limit state. When $r\in [r_2,r_e]$, the tangential stress of the frozen wall decreases with the increase in r. Under the elastic limit state, the maximum tangential stress of frozen wall appears at r = r₁ = 1.25. Under the plastic limit state, the maximum tangential stress of frozen wall appears at r = r₂ = 1.63.

6. Discussion

The main content of this paper is based on the fact that a double-row-pipe frozen wall is equivalent to a trapezoidal heterogeneous frozen wall, and the existing research has already compared the double-row-pipe frozen wall to a homogeneous frozen wall [17]. Now we will introduce them briefly, and compare and analyze the two equivalent methods. According to unified strength theory (Equation (18)), the general solution of the elastic zone of the homogeneous frozen wall is

$$\{\begin{array}{l}{\sigma }_r^e=\frac{A}{r^2}+2D\\ {\sigma }_{\theta }^e=-\frac{A}{r^2}+2D\end{array}$$

(29)

where A and D are coefficients.

According to the boundary conditions, when $r=r_{\rho }$, ${\sigma }_{r_{\rho }}^e={\sigma }_{r_{\rho }}$; when $r=r_{\rho }$, ${\sigma }_{r_{\rho }}^e=P$. Substituting the boundary conditions into Equation (29), we get:

$$\{\begin{array}{l}A=\frac{r_{\rho }^2{r}_e^2(P-{\sigma }_{r_{\rho }})}{r_{\rho }^2-r_e^2}\\ 2D=\frac{r_{\rho }^2{\sigma }_{r_{\rho }}-r_e^{2}P}{r_{\rho }^2-r_e^2}\end{array}$$

(30)

The equation of the stress field in the elastic zone is:

$$\{\begin{array}{l}{\sigma }_r^e=\frac{r_{\rho }^2}{r^2}\frac{r_e^2(P-{\sigma }_{r_{\rho }})}{r_{\rho }^2-r_e^2}+\frac{r_{\rho }^2{\sigma }_{r_{\rho }}-r_e^{2}P}{r_{\rho }^2-r_e^2}\\ {\sigma }_{\theta }^e=-\frac{r_{\rho }^2}{r^2}\frac{r_e^2(P-{\sigma }_{r_{\rho }})}{r_{\rho }^2-r_e^2}+\frac{r_{\rho }^2{\sigma }_{r_{\rho }}-r_e^{2}P}{r_{\rho }^2-r_e^2}\end{array}$$

(31)

The equation of the stress field in the plastic zone is:

$$\{\begin{array}{l}{\sigma }_r^p=\frac{C\cos \phi }{\sin \phi }\left[r^{\frac{4(1+b)\sin \phi }{(2+b)(1-\sin \phi )}}-1\right]\\ {\sigma }_{\theta }^p=\frac{4(1+b)\left[C\cos \phi +{\sigma }_r^p\sin \phi \right]}{(2+b)(1-\sin \phi )}+{\sigma }_r^p\end{array}$$

(32)

At the interface of elastic–plastic zone, i.e., $r=r_{\rho }$; according to Equation (31), the radial stress and tangential stress at the interface are:

$$\{\begin{array}{l}{\sigma }_{r_{\rho }}^e=\frac{r_{\rho }^2}{r_{\rho }{}^2}\frac{r_e^2(P-{\sigma }_{r_{\rho }})}{r_{\rho }^2-r_e^2}+\frac{r_{\rho }^2{\sigma }_{r_{\rho }}-r_e^{2}P}{r_{\rho }^2-r_e^2}={\sigma }_{r_{\rho }}\\ {\sigma }_{{\theta }_{r_{\rho }}}^e=-\frac{r_{\rho }^2}{r_{\rho }{}^2}\frac{r_e^2(P-{\sigma }_{r_{\rho }})}{r_{\rho }^2-r_e^2}+\frac{r_{\rho }^2{\sigma }_{r_{\rho }}-r_e^{2}P}{r_{\rho }^2-r_e^2}=\frac{2r_e^{2}P-(r_e^2+r_{\rho }^2){\sigma }_{r_{\rho }}}{r_e^2-r_{\rho }^2}\end{array}$$

(33)

Equation (33) should satisfy strength theory Equation (18), namely:

$${\sigma }_{{\theta }_{r_{\rho }}}^e-{\sigma }_{r_{\rho }}^e=\frac{2r_e^2(P-{\sigma }_{r_{\rho }}^e)}{r_e^2-r_{\rho }^2}=\frac{4(1+b)\left[C\cos \phi +{\sigma }_{r_{\rho }}^e\sin \phi \right]}{(1-\sin \phi )(2+b)}$$

(34)

At the same time, according to the stress continuity conditions at the interface between the elastic zone and the plastic zone, it can be obtained that

$${\sigma }_{r_{\rho }}^e={\sigma }_{r_{\rho }}^p=\frac{C\cos \phi }{\sin \phi }\left[r_{\rho }{}^{\frac{4(1+b)\sin \phi }{(2+b)(1-\sin \phi )}}-1\right]$$

(35)

Substituting Equation (35) into Equation (34) obtains the relative radius of the plastic zone. In the zone $[1,r_{\rho })$, the radial stress and tangential stress are plotted by Equation (32), respectively. In the zone $[r_{\rho },r_e]$, the radial stress and tangential stress are plotted by Equation (31), respectively. These stress fields are compared with the stress field of a heterogeneous frozen wall. The comparison of stress fields of heterogeneous and homogeneous frozen walls is shown in Figure 8. It can be seen from Figure 8 that when the plastic zone of both the heterogeneous frozen wall and the homogeneous frozen wall extends to the middle of the inner and outer rows of frozen pipes, the corresponding elastic ultimate bearing capacity and plastic ultimate bearing capacity of the heterogeneous frozen wall are 2.93 MPa and 8.03 MPa, respectively. The elastic ultimate bearing capacity and plastic ultimate bearing capacity of the homogeneous frozen wall are 2.92 MPa and 6.79 MPa, respectively.

In order to analyze the influence of the strength theory parameter b on the mechanical properties of the frozen wall, the radius of the plastic zone $R_{\rho }$, the elastic ultimate bearing capacity $P_e$, the plastic ultimate bearing capacity $P_p$, and the safety factor K of the frozen wall are calculated when b takes different values; where K = P_p/P, P are the external loads borne by the frozen wall, as shown in Figure 9 and

Table 1. Mechanical characteristics of heterogeneous and homogeneous frozen walls with different values of strength theory parameter b.

b	P (MPa)	Heterogeneous Frozen Wall				Homogeneous Frozen Wall
		rρ	Pe (MPa)	Pp (MPa)	K	rρ	Pe (MPa)	Pp (MPa)	K
0	2.925	1.63	1.20	3.11	1.06	1.94	1.35	2.93	2.17
0.2		1.44	1.31	3.44	1.18	1.57	1.47	3.28	2.23
0.4		1.35	1.40	3.73	1.28	1.45	1.57	3.52	2.24
0.6		1.30	1.47	3.97	1.36	1.37	1.66	3.75	2.26
0.8		1.26	1.54	4.19	1.43	1.32	1.73	3.95	2.28
1		1.23	1.60	4.38	1.50	1.28	1.80	4.13	2.29

. When the strength theoretical parameter b increases, the radius of the plastic zone of the frozen wall decreases, and the elastic ultimate bearing capacity and plastic ultimate bearing capacity of the frozen wall increase. It can be seen from Table 1 that the strength theoretical parameter b increases from 0 to 1, and the elastic ultimate bearing capacity and plastic ultimate bearing capacity of the heterogeneous frozen wall increase by 33.3% and 40.8%, respectively. The elastic ultimate bearing capacity and plastic ultimate bearing capacity of the homogeneous frozen wall increase by 33.3% and 41.0%, respectively. The elastic ultimate bearing capacity of the heterogeneous frozen wall is smaller than that of the homogeneous frozen wall. The plastic ultimate bearing capacity of the heterogeneous frozen wall is greater than that of the homogeneous frozen wall. The safety factor of the heterogeneous frozen wall is greater than that of the homogeneous frozen wall. When the intermediate principal stress is considered, the potential of the material can be fully exploited, and the ultimate bearing capacity of the frozen wall can be improved. Therefore, the influence of intermediate principal stress should be properly considered in the design of a frozen wall, so that the potential of materials can be fully utilized, the waste of material resources can be avoided, and the economic benefits can be improved.

7. Conclusions

Based on unified strength theory, the elastic–plastic analysis of the double-row-pipe frozen wall is carried out, and the analytical solutions of the elastic–plastic stress field distribution, the elastic ultimate bearing capacity, the plastic ultimate bearing capacity, and the relative radius of the plastic zone of the frozen wall are derived. The stress field distribution of the frozen wall in the control layer of the surface soil layer in Qingdong Mine is analyzed using the analytical solution, and the numerical simulation is carried out based on COMSOL. This study determined whether the frozen wall of this layer meets the strength requirements and analyzed the influence of strength theoretical parameters on the stress characteristics of heterogeneous and homogeneous frozen walls. The following conclusions were obtained:

• Considering the FGM characteristics of the frozen wall, the temperature field of the double-row-pipe frozen wall is equivalent to a trapezoidal distribution along the radial direction, and the elastic modulus E and cohesion c of the frozen wall are considered to be linearly distributed with the radius. Based on unified strength theory, the analytical solutions of the elastic–plastic stress field, the elastic ultimate bearing capacity, the plastic ultimate bearing capacity, and the relative radius of the plastic zone of the double-row-pipe frozen wall are derived; the results can provide a theoretical reference for the design of frozen wall.
• We analyzed a frozen wall with a vertical depth of 225 m in Qingdong Mine. When b is 0.5, the elastic ultimate bearing capacity of the frozen wall is 1.44 MPa, and the plastic ultimate bearing capacity of the frozen wall is 3.55 MPa. The external load p = 2.925 MPa of the frozen wall is between the elastic ultimate bearing capacity and the plastic ultimate bearing capacity, indicating that the stress state of the frozen wall is safe. Under the action of external load, the relative radius $r_{\rho }$ of the plastic zone of the frozen wall is 1.32, which is located between the inner and outer frozen pipes.
• Based on the Qingdong Mine, the analytical solution was calculated and compared with the numerical solution, and it is concluded that the stress distribution of the two solutions is consistent, which mutually verifies the accuracy of the two solutions.
• The strength theory parameter b has an important influence on the stress state of heterogeneous and homogeneous frozen walls. The radial stress increases with the increase in b and the relative radius r, the tangential stress increases with the increase in b, and the tangential stress first increases and then decreases with the increase in r. The elastic ultimate bearing capacity P_e and P_p of the frozen wall increase with the increase in b. When b increases from 0 to 1, the elastic ultimate bearing capacity and plastic ultimate bearing capacity of the heterogeneous frozen wall increase by 33.3% and 40.8%, respectively, and the elastic ultimate bearing capacity and plastic ultimate bearing capacity of the homogeneous frozen wall increase by 33.3% and 41.0%, respectively. The relative radius of the plastic zone decreases with the increase in b.