3.1. Assumptions and Computational Model
The Gongbei Tunnel is a typical engineering case of the FSPR method, and is used as a case study in this paper. To establish an efficient computational model, the typical parts of the FSPR structure are taken as the computational model, as shown in
Figure 3. The size of the model is 11 m × 1.977 m. For the unsteady-state conjugate heat conduction model, the left and right sides and the lower part of the model are the soil, which can be regarded as the second temperature boundary condition. The upper part of the model is the ground surface, which can be regarded as the third boundary condition [
14,
15]. Inside the model, the freezing tube as a cold source can be regarded as the first temperature boundary condition, and the temperature of the brine flowing into the main tube is used as the boundary temperature [
16]. The steel pipe is a conjugate heat transfer surface [
17]. In this computational model, the strong coupling integral computational method is selected, and the general control equation is used to find a global solution. Therefore, the initial values of the temperature field, velocity field, and pressure field need to be given for the whole region. For the initial value of the velocity field, both the solid domain and fluid domain are recorded as 0 m/s; for the initial pressure distribution, the fluid and solid domains are denoted as 1 atm. The upper boundary that represents the ground surface is regarded as the convective heat transfer boundary. The heat flux
q in the lower boundary is 0.4 W/m
2, as expressed by Equation (1). The left and right boundaries are the adiabatic boundaries.
where
is the thermal conductivity of the soil, and
is the geothermal gradient, with a value of 0.03 °C/m.
3.2. Governing Equation of Unsteady-State Conjugate Heat Transfer Model
Ground freezing is an unsteady-state heat transfer process with complex phase transition. In the solidification process of pure substances such as water, solidification occurs at a single temperature, and the solid phase and liquid phase are separated by a clear moving interface. However, the soil freezing occurs in a larger temperature range, and there is a separation of the solid and liquid phases by moving regions of two phases in the process [
18]. The phase transition problem is mathematically strongly nonlinear, meaning that the governing equation is linear, but the position of the two-phase interface must always be determined, and the energy conservation condition of the interface is nonlinear. It is not possible to use the superposition principle of solutions. Therefore, most of these problems are treated by numerical simulation methods [
19]. When using numerical methods to solve phase transition problems, there are generally two methods to deal with the moving boundary in the process of phase transition: The first focuses on the solution of the phase transition interface. After determining the interface position, the temperature distributions in the solid and liquid regions are solved. The second method is to assume the problem as a single-phase nonlinear heat conduction problem, determine the temperature or enthalpy distribution in the whole solution region, and then determine the position to reach the phase transition temperature as the phase transition interface [
20,
21,
22]. The second method is convenient and practical, and is more suitable for the soil phase transformation process, which has no clear interface. For the second method, the sensible heat capacity method is used. Assuming that the physical properties of the solid and liquid phases are spatially invariant, ignoring the possible natural convection in the liquid phase, the conjugate heat transfer interface between soil and air is the steel pipe. The sensible heat capacity method takes temperature as the function to be solved, without introducing the concept of enthalpy, and establishes a unified energy equation for the whole region. For the treatment of phase transition, the specific heat is expressed in the form of equivalent specific heat [
23,
24]. For convenience of explanation and comparison with the enthalpy method, the equivalent specific heat of phase transition that occurs at a given temperature Tm is as expressed in Equation (2):
where
is a Dirac function and, thus, has a heat capacity model as shown in Equation (3):
References [
25,
26] proved the equivalence between Equation (3) and the commonly used equations describing the phase transition problem.
For the phase transition that occurs in the temperature range near
, the influence of
should be taken into account when constructing the equivalent specific heat. The expression
should be expressed as shown in Equation (4):
When the specific heat and coefficient of thermal conductivity of the solid phase and liquid phase are constant, Equation (5) can be obtained:
References [
18,
27,
28] note that the phase change of water in frozen soil can be divided into three regions:
- (1)
Severe phase transition zone: when the temperature in this zone changes by 1 °C, the variation in unfrozen water content is greater than or equal to 1%;
- (2)
Transition zone: when the temperature in this zone changes by 1 °C, the variation in unfrozen water content is between 0.1% and 1%;
- (2)
Frozen zone: when the temperature in this zone decreases by 1 °C, the amount of the water phase becoming ice is less than 0.1%.
Accordingly, when the sensible heat capacity method is used to deal with the phase change problem of the soil freezing process, the change in the unfrozen water in the soil should be divided into at least three sections according to the experimental data—a violent phase change zone, transition zone, and frozen solid zone—and then the phase change should be treated with the equivalent specific heat in each section.
Similar to references [
29,
30], the specific heat of the soil region is as shown in Equation (6):
where
represent the specific heat of unfrozen soil and frozen soil, respectively (unit: J/(kg·K)); L is the latent heat of the phase change of water;
and
are the total water content and ice content of the frozen soil, respectively; and
and
are the upper and lower boundary temperatures of the frozen soil’s phase transition zone, respectively.
In the hollow pipe, when the air pressure is low and the temperature is high, the air can be treated as an ideal gas [
31]. The air in the fluid domain can also be regarded as a compressed fluid and a viscous fluid. At this time, the continuity equation and motion equation (Navier–Stokes equation) are changed. Considering the causes of fluid movement, the flow state of air can be assumed to be laminar flow. The continuous condition of heat flux can also be treated according to laminar flow in the conjugate heat transfer interface [
32].
According to the above assumptions, the general strong coupling control equations of the computational model can be obtained, containing the energy equation, ideal gas state equation, continuity equation, and motion equation.
The energy equation can be expressed with Equation (7):
The ideal gas state equation can be expressed with Equation (8):
The continuity equation can be expressed with Equation (9):
The motion equation can be expressed with Equation (10):
where
represents the fluid density,
c represents the specific heat shown in Equation (6),
represents the velocity vector,
λ represents the thermal conductivity,
T represents the temperature,
t represents time,
represents the coefficient of viscosity,
represents the intensity of the internal heat source, and
R represents the gas constant,
.
3.4. The Shape of the Enhancing Freezing Tubes
Research on the shape of enhancing freezing tubes, which has a significant effect on the freezing effect, is lacking. This is the basis of the follow-up work in this paper. Since the air flow in the hollow pipes has a great impact on the enhancing freezing tubes [
33], enhancing freezing tubes with different cross-sectional shapes were compared to improve the heat transfer efficiency. The enhancing freezing tubes of different shapes—such as circular, crescent, groove, and semicircular—are shown in
Figure 4. To amplify and compare the computational results of the freezing effect, these four types of enhancing freezing tubes were given the same cross-sectional area and tripled to 0.047 m
2.
The temperature cloud graphs of the enhancing tubes with different shapes over the course of 30 days during the active freezing phase are shown in
Figure 5. After 30 days of active freezing, the frozen soil curtain formed between the two pipes based on the different shapes of the enhancing tubes. The frozen soil curtain of Shape A and Shape D was more uniform, while the non-uniformity of the frozen soil curtain of Shape B and Shape C was greater than that of Shape A and Shape D, which may lead to frost heaving.
As shown in
Figure 6, the average temperature of the hollow pipe with a circular tube was the highest, while that of the hollow pipe with a crescent tube was the lowest. If the contact area between the enhancing freezing tube and the inner surface of the hollow pipe is increased as much as possible, the cooling capacity of the enhancing freezing tubes can be effectively and quickly transferred to the soil. The average temperature of the circular enhancing freezing tubes of Shape A was the highest; the freezing effect of these tubes was weaker than that of Shape D due to the small contact surface between the enhancing freezing tubes and the steel pipe. Therefore, the scheme of Shape A should be abandoned.
Although the average temperature of the enhancing freezing tubes of Shape B and Shape C was lower, the non-uniformity of their frozen soil curtain thickness was relatively greater. Taking the area below −10 °C as the strength of the frozen soil curtain [
34,
35], it can be concluded from
Table 6 that the thickness of the frozen soil curtain with Shape B was the thickest, and had the greatest degree of unevenness. Uneven frost heaving is detrimental to the pipeline and the surrounding environment [
36,
37]. To avoid uneven frost heaving, Shape B should be abandoned. The thickness of the frozen soil curtain with Shape C and Shape D was similar, but the degree of unevenness in the thickness of Shape D was far less than that of Shape C. Therefore, the semicircular enhancing freezing tube is the best scheme in terms of the freezing effect and the degree of unevenness in the thickness of the frozen soil curtain.
3.5. Calculation Scheme
Research on the layout of enhancing freezing tubes is the most important factor in the FSPR freezing scheme. The layout of the enhancing freezing tubes has a great influence on the freezing effect and time, which is the most concerned part in construction sites. To understand the influence of the layout, operating duration, and heat preservation measures of the semicircular enhancing freezing tubes on the freezing effect, three different simulation schemes were set up, as shown in
Figure 7. The initial temperature of the model was set at 20 °C, and the surface boundary was set as the third boundary condition, with a surface heat transfer coefficient of 15 W/(kg·K).
The layout of the enhancing freezing tubes has great influence on the freezing effect [
38]. A freezing scheme of enhancing freezing tubes in two different positions was considered, as shown in
Figure 5a. The enhancing freezing tubes were arranged far from or near to the excavation side with an angle of 15°, recorded as Scheme A and Scheme B, respectively. The master freezing tubes and enhancing freezing tubes were operated from the beginning to the end of the freezing process.
The operating duration of the enhancing freezing tubes was also set as a simulation scheme involving weather, to make the enhancing freezing tubes work and maintain an active freezing duration of 60 days, recorded as Scheme C and Scheme D. The master freezing tubes continued working during this period.
Most of the engineering literature on the use of the artificial ground freezing method posits that air convection has a great influence on the freezing temperature field [
39]. The air convection interface should be insulated to limit and reduce the loss of cooling capacity caused by air convection, so as to concentrate the cooling capacity on the formation and development of frozen soil. Fortunately, each section of the hollow jacking pipes has thermal insulation treatment in the longitudinal direction, and the air in the jacking pipe did not undergo directly convect with the atmosphere during the freezing construction of the Gongbei Tunnel. The convective heat transfer in the limited space has limited heat dissipation; even this kind of air convection can homogenize the cold capacity of the enhancing freezing tubes, and make the hollow pipes play the role of “freezing pipe”, which is beneficial to the freezing effect. Therefore, whether the strengthening of heat preservation measures is beneficial to the freezing effect must be discussed through simulation comparison. The simulation scheme is shown in
Figure 7b, where the models with and without heat preservation measures are recorded as Scheme E and Scheme F, respectively.
According to the engineering scheme [
40], polyurethane foam was selected as the heat preservation material. The thickness of the heat preservation material was set to 0.03 m based on construction experience [
41,
42]. The thermal conductivity of the heat preservation material was set to 0.04 W/(kg·K), the density was 34 kg/m
3, and the specific heat was 2016 J/(kg·K).