Applicability of Vyalov’s equations to ice wall strength estimation

A simple analytical relations are commonly used in engineering practice to calculate ice wall thickness. One of them is Vyalov’s relation that takes into account the features of a real technological process of tubing lining and the inelastic deformation associated with frozen soil creep. An estimation of applicability and margin of safety of this equation is an issue of engineering mechanics. In this paper, we propose a mathematical model for description of ice wall deformation under natural external loading and present the results of the computational experiments in which an optimal thickness for the ice wall is determined. Based on this simulation, we modify the existing analytical relation, which makes it possible to calculate the thickness of an ice wall of unlimited height.


INTRODUCTION
he ice wall is an important engineering structure the integrity of which determines the safety of sinking operations in vertical mine shafts under construction. The effectiveness of the ice wall as a temporary protective lining depends on the adequate evaluation of its thickness. Optimal assessment of the wall thickness comes from the need to reduce the cost of mining operations at the design stage. This will help to exclude the failure of an ice wall and the breakthrough of groundwater into the mine. Furthermore, the savings in artificial freezing costs may be significant. Engineering calculations for ice walls are currently performed for two limit states -strain condition and strength condition. In strain calculations, it is necessary to determine the minimum thickness of the ice wall at which its deformation in the design phase does not exceed the value acceptable for shaft construction and promotes no dangerous movements of freezing columns. The objective of the limit strength design is to calculate the optimal minimum thickness of the ice wall which T prevents its failure under external loads, i.e. no cracks occur in the wall. At this thickness value, the ice wall is in the ultimate stress-strain state, i.e., in the state where the stresses at a given instant of time do not exceed the frozen soil strength at the same instant. The most common methods for calculating of ice wall thickness are those proposed by Domke [1] and Lame-Gadolin [1]. The Lame-Godolin formula was obtained in a linear elastic approximation and did not prove effectiveness at soil depths of more than 50 m. The Domke equation takes into account the plastic deformations of frozen soil, the onset of which is determined according to the criterion for the difference between the largest and smallest principal normal stresses. Because the strength of saturated soils increases with an increase in the mechanical pressure up to a certain value, the application of the Mohr-Coulomb and Druker-Prager criteria makes it possible to calculate plastic strain with higher accuracy and to take into account the fact that the frozen soil strength increases with decreasing temperature [2]. The triggering of high hydrostatic pressure leads to melting and failure of the ice contained in the pores. This has given impetus to different modifications of these criteria [3]. The thermo-hydro-dynamical models of soil freezing which take into account plastic strains are presented in [4][5][6][7]. In [4][5][6], the so-called Barcelona Basic Model [8] is used to describe the plastic deformations of saturated frozen soils. According to this model, the stress-strain state of soils is determined by analyzing two state variables -effective stress and a parameter that describes the moisture migration effect. For the non-frozen soil, the Barcelona Basic Model is simplified and reduced to the modified Cam-Clay model. The thermo-hydro-dynamical models based on this approach can be employed to consider the influence of moisture migration and mechanical pressure on the plastic deformation of frozen soils. In [6], the Barcelona Basic Model is generalized to the case of large deformations. Another way [7] has been proposed to derive the constitutive relations for describing the elastoplastic behavior of frozen soils in the framework of a thermodynamic approach. It has been suggested in [9] that the yield surface defined in terms of a double clay hardening model can be used to assess plastic strains. Elastic and plastic strains characterize the instantaneous response of the material to the applied loads. However, due to the pronounced rheological behavior of saturated soils under long-term loads, their deformation increases and strength reduces. This feature of soils is of prime importance for shaft sinking where the loads from the surrounding rocks and ground water are taken over by the lateral surface of a mine during the entire time required for lining construction. In order to calculate the optimal thickness of an ice wall () and to take into account the rheological behavior of frozen soils, S.S.Vyalov suggested the following relations [10,11]: where  is the optimal ice wall thickness, a is the cylinder inner radius, k' is the coefficient dependent on soil compaction conditions, m, A(t p , ) are the parameters characterizing the rheological behavior of the soil, t p is the loading period, p is the load applied to the ice wall, h is the height of unfixed part of mine shaft, is the maximum radial movement of the frozen soil, s(tp) is the ultimate strength of frozen soil to uniaxial compression, c(tp) is the cohesion, and  is the angle of internal friction.
Eqns. (1)-(2) estimate for the cylinder wall thickness according to the criterion for the maximum radial movement of its internal wall, and Eqns. (3)-(4) -the ultimate stress criterion. We notice that Eqns. (1) and (3) are obtained on the assumption that the ice wall is a cylinder of finite height, and Eqns. (2) and (4) -on the assumption that the ice wall is a cylinder of unlimited height. It is also worth noting that formula (3) is obtained based on relationship (1) assuming that the frozen soil exhibits ideal plastic behavior described by the von Mises model. In [12], the applicability of formula (1) for calculating the frozen wall thickness was evaluated. The aim of our work is to analyze the correctness of analytical relationships (3) and (4) for calculating the frozen wall thickness according to the ultimate stress criterion. For this purpose, the ice wall thickness calculated by these formulas was compared with the numerical simulation results obtained by the formulation whose geometry is close to the real shaft sinking conditions. Using the obtained data, we propose a modification of formula (4) because it shows better agreement with the numerical simulation results. In addition, we perform a comparative analysis of the values of ice wall thickness calculated by formulas (1)-(4) and by their modifications, which provides guidelines for using each of these formulas.

MATHEMATICAL FORMULATION
n order to evaluate the correctness of relations (3) and (4), we have computed the problem of deformation of ice wall exposed to the external radial load p applied to the lateral surface of the cylinder and the vertical load p' acting on the upper end of the cylinder. The calculation scheme of shaft sinking is given in Fig.1. An ice wall is modeled as a twopiece cylinder. One part is a hollow cylinder of a height equal to that of the ice wall section without lining, and another part is a solid section. The soil inside the ice wall is assumed to be a thawed soil that experiences only elastic deformations. The load applied to the lateral surface of the cylinder is calculated by the formula for mining and hydrostatic pressure, which takes into account the influence of cohesive forces and the results of hydrogeological observations. The mathematical model is based on the following hypotheses and assumptions: 1) The computational domain has radial symmetry, and therefore calculations are made in axisymmetric formulation.
2) The rocks under consideration are isotropic at macrolevel. 3) Deformations as well as strain increments at each step are small. 4) The time and temperature effects on the strength and elastic characteristics of rocks are taken into account parametrically. 5) The surrounding thawed and burden formations generate constant rock pressure. 6) A linear relation exists between the stress and elastic strain tensors. 7) The minimum value for the ice wall thickness is determined on the assumption that no plastic deformations are present on its internal wall. 8) The Mohr-Coulomb criterion is used as a yield criterion that is commonly applied to describe soil failure (especially, that of cohesionless soil). By this criterion, the rock failure begins when the maximum tangential stress n reaches the critical value which is dependent on the normal stress  n in the same area.
The mathematical statement of the problem includes equilibrium Eqn. (5), Hooke's law (6), geometric relation for small strain tensor (7), yield condition (8), and associated plastic flow rule (9): where ={ r ,   , z ,  rz } is the stress tensor, ={ r ,   , z ,  rz } is the strain tensor, С is the tensor of elastic constants which reduces in the isotropic case to two elastic constants (E′ -Young's modulus,  -Poisson's ratio), u={u r , u z } is the displacement vector,  p ={ p r,  p , p z,  p rz} is the plastic strain tensor, F is the Mohr-Coulomb yield criterion, n is the maximum tangential stress in the area with a normal n,  n is the normal stress operating in the same area,  is the indefinite multiplier determined using Prager's compatibility conditions, and the dots over the symbols denote time derivatives. The Prager compatibility conditions are written as The system of Eqns. (5)-(11) is supplemented with boundary conditions that correspond to the calculation scheme from Fig.1.
where the vectors P={p,0} and P′={p′,0} correspond to the calculated loads obtained by the formulas: In formulas (15) -(18), the following notation is used:  m =2000 kg/m 3 -average density of the material, h′ -bed rock density, w=1000 kg/m 3 -water density, g=10 m/s 2 -gravitational constant, hgw=1.5 m -groundwater depth. The load ph is considered only for saturated rocks.

MATERIALS AND METHOD
he allowable ice wall thickness was determined for three typical materials: sand, chalk and clay. We considered the bed rock depth values of 100m, 200m, 300m, and 500m. Tabs. 1 and 2 give the mechanical parameters of the soil in thawed and frozen states which were provided by the Institute of Nature Management of the NAS of Belarus. The parameters of the frozen soil were obtained at T=-8 0 C and for load duration of 12 h. These data were used for calculating the load acting on the ice wall according to formulas (15) -(18). The obtained data are given in Tabs. 3     where i[0.02; 1.28] indicates the size of the element used to determine the von Mises stress value m on the internal wall of the ice cylinder, and  min m is the von Mises stress value on the internal wall of the ice cylinder obtained in the solution of the problem with minimum element size. Fig. 2 shows the effect of size of the finite element of the grid on relative error. The results demonstrate that the optimal size of the element is 0.04m. It is interesting that if the element size is maximum (1.28m), then the relative error does not exceed 8.3%, which is well within the permissible limits of an engineering error.

Comparison of numerical simulation results and analytical estimates
he developed mathematical model is used for assessing the optimal ice wall thickness according to the criterion for ultimate stress-strain state. The calculation results are summarized in Tab. 6. When the optimal ice wall thickness is calculated by Eqn. (2), the ultimate uniaxial compression strength values are taken from Tab. 2. The values for the cohesion and internal friction angle used in Eqn. (4) are also given in Tab. 2. The calculation results are summarized in Tabs. 7 and 8. The obtained results demonstrate that, for all materials under study, the design thicknesses obtained by Eqn. (3) exceed those found using Eqn. (4). This can be attributed to the fact that the internal friction is ignored, which increases the strength. In general, the scatter in the estimates obtained by Eqns. (3) and (4)    A comparison of the ice wall thickness values found during the numerical solution of problem (5)-(18) with the results obtained by Eqns. (3) and (4) is illustrated in Fig.3. For all materials under study, Eqn. (3) shows the most distinction from the numerical simulation results (Fig. 3a). The greatest ice wall thickness reserve is seen for clay and sand, especially, for the layers at a depth of 500 m. This is because the ice wall is subjected to maximum loads at this depth. In this case, the deviations from the numerical results obtained for clay and sand are the greatest, 12.8 m and 13.7m, respectively. It was found that the results of numerical simulations of the ice wall thickness for chalk and sand are higher than the estimations obtained by Eqn. (4) (Fig.3 b). The most distinction from the numerical simulation results is 1.49m for sand at a depth of 300 m and 0.242 for chalk at a depth of 300 m. A comparison of the same results but obtained for clay shows that the numerical calculations give the higher ice wall thickness values for the material at depths of 100-200 m. At depths of 300-500 m, the simulation results are comparable with those obtained by formula (4); the relative error does not exceed 2.8%.

Modification of Eqn. (4)
Analysis of the obtained results has revealed that Eqn. (4) has the best correlation with the numerical simulation results. This equation was modified here in order to get a more accurate description of the relation between the designed ice wall thickness and the applied load. The modification made it possible to minimize the difference between the calculation results obtained by this formula and the numerical simulation data. Thus, we have    In addition to approximation (20), it is suggested that formula (4) can be modified as The values of parameter z for sand, chalk and clay are given in Tab. 11. The calculated results are obtained for the designed thickness by Eqn. (21) and the numerical simulation results are shown for each material in Fig.5. The approximations are nonlinear, which is caused by the nonlinear dependence of stress on the radial coordinate and, as a consequence, on the desired ice wall thickness. The nonlinear character is determined by the material parameters of the soil, e.g., for chalk and sand, the curve is convex upwards and for clay -downwards. Chalk and sand exhibit similar cohesion values (6.2 and 6.33 MPa), and clay has significantly lower cohesion value (3.78 MPa). It is interesting that almost the same nonlinear behavior is observed for Eqn. (4).   A comparison of the results of approximation by Eqn. (4) shows that in case of piecewise-linear dependence (20) they are in rather good qualitative and quantitative agreement. A deficiency of this approach is the necessity of calculating four material parameters, as well as the intervals for determining piecewise-linear correction coefficients. Therefore, it is advisable to use formula (21) in order make rapid estimation for the designed ice wall thickness. In this case, it is necessary to determine only one material parameter, whereas formula (20) can be used when there is a need for more accurate calculation.

Comparative analysis of Vyalov's formulas for ultimate displacements and stresses
where the material function g(p) is dependent on the value of the load acting on the ice wall and determined by applying the numerical simulation results. The function for the materials under study was taken from [12]. For sand (Fig.6a), the maximum discrepancy between Eqns. (4) and (20) (2) cannot be used in this case to calculate the ice wall thickness at a depth of 500 m because the base of power function becomes negative. For chalk Fig.6(b), there are two similar data blocks. The first group includes Eqns. (4) and (20), and another -Eqns. (2) (3) and (22). Analogous results are obtained for clay (Fig. 6c). Two groups of equations provide almost similar values for ice wall thickness: (4), (20) and (2)

CONCLUSION
e have performed a theoretical study to evaluate the applicability of analytical Vyalov's formulas to the calculation of the ultimate stress state of an ice-rock cylinder of unlimited and finite heights. The numerical results showed that the design values of the ice wall thickness obtained by the formula for a cylinder of limited height have a margin of safety for all considered rocks (sand, chalk, and clay). This is associated with the fact that the internal friction is ignored, which increases the strength. The application of the Coulomb-Mohr criterion made it possible to model the deformation of the ice-rock cylinder. A comparison of the thicknesses obtained numerically with the results found by Eqn. (3) showed that the greatest thickness reserve was observed for clay and sand. Using the results of numerical simulation of the stress-strain state of the ice wall, we suggested two modifications of Eqn. (4). The first variant involves the use of the piecewise-linear relation and enables describing the simulation results both qualitatively and quantitatively. The second is a simple variant and it can be used for making rapid estimation for the designed ice-wall thickness. A comprehensive comparative analysis of six formulas (four Vyalov's formulas and two modifications) has revealed that in strength calculations it is advisable to use Eqn. (20) for all soils examined here because Eqn. (3) causes the excess margin of thickness to occur, and Eqn. (4) must be used with some prescribed margin of safety. In ultimate displacement calculations for chalk, Eqn. (22) or (2) should be used at the load not exceeding 2.25 MPa, which corresponds to the soil depth of 300 m. In creep calculations for clay at the load of 2.301 MPa, which corresponds to the soil depth of 300 m, relations (2) and (22) can be used. In order to calculate the ice-wall made of frozen sand at the load less than 4.889 MPa (up to the soil depth of 300 m), it is advisable to use relation (22) and relation (2) at the load less than 3.222 MPa (up to the soil depth of 200 m). At high loads that exceed 4.899 MPa for sand and 2.301 MPa for clay there is a need for further investigations so that the existing constitutive equations describing the stress-strain state of rocks can be modified. W