Evaluation of Gas Migration and Rock Damage Characteristics for Underground Nuclear Waste Storage Based on a Coupled Model

1State Key Laboratory of Eco-Hydraulics in Northwest Arid Region, Xi’an University of Technology, Xi’an 710048, China 2Institute of Geotechnical Engineering, Shaanxi Provincial Key Laboratory of Loess Mechanics and Engineering, Xi’an University of Technology, Xi’an 710048, China 3Jiangsu Key Laboratory of Advanced Manufacturing Technology, Huaiyin Institute of Technology, Huai’an, Jiangsu 223003, China 4School of Civil Engineering, Henan Polytechnic University, Jiaozuo 454003, China


Introduction
Nuclear waste contains strong radioactivity, large calorific value, high toxicity, and long half-life nuclides, which should be isolated from the environment for human residence in a reliable way.In a large number of disposal schemes, the deep geological disposal of nuclear waste is a widely accepted and feasible disposal scheme at present [1][2][3][4].The geological disposal of nuclear waste is to lay waste in the geological body 500∼1000 m away from the surface, so that it can be isolated from the living environment of mankind permanently.The underground engineering of high level radioactive waste is called nuclear waste repository.The nuclear waste disposal repository adopts the design idea of "multiple barrier system", that is, to store wastes in waste cans, outside the buffer materials, and then outward to the surrounding rock.Nuclear waste repository is a special deep underground engineering.Therefore, it is necessary to solve rock mechanics problems in the location, design, construction, and stability evaluation of repository [5][6][7][8][9][10].The disposal of nuclear waste in deep sedimentary rocks has been studied by some countries.For example, in Canada, low radioactive waste and medium radioactive waste are being put forward to be disposed of in the sedimentary rocks of Ontario province.A large amount of radioactive waste can be produced in the underground storage.The gas will be generated in underground nuclear waste repository due to the metal corrosion.
The generated gas tends to migrate through underground barrier.The increase in pressure caused by the generated gas induces the formation of factures and the expansion of the pore of the barrier.This affects the long-term function of waste repository and causes the risk to ecological environment, as shown in Figure 1.Therefore, the evaluation of the safety of the waste reservoir in the sedimentary rocks should be researched systematically.

Science and Technology of Nuclear Installations
Mechanical damage due to pressure rise: Disruption to near field barriers; Fracture creation in host rock.
Gas release to biosphere: Transport of gaseous radionuclidesto biosphere; Flammability hazards at surface; e release of toxic gases; Human intrusion leading to rapid gas releases.

Modification to groundwater flow and transport:
Forcing of contaminated water from repository or canisters; Induced groundwater flows; Transport associated with bubbles; Providing a barrier to water flow.
Chemical effects of gas: Alteration of chemical environment in repository, potentially affecting some repository processes, e.g.spent fuel modification

GAS GENERATION
Figure 1: Effects of gas generation in a repository [9].
For all major international waste separation schemes, the study of gas migration and its effect on the rock has become a key goal.Horseman et al. 1999 carried out the gas injection test of the precompacted bentonite and suggested that the gas entry and breakthrough were accompanied by the development of the propagation path of clay.The initial saturated bentonite is impermeable without the absence of pressureinduced gas channels.Davy et al. [11] experimentally found that the breakthrough pressure is roughly equivalent to the sum of the expansion pressure of the clay and the external water pressure, and it is considered that gas is the key problem of the clay engineering barrier performance.Alonso et al. [1] established the THM model to describe the discontinuity of gas transportation and realized the numerical solution of the model by compiling the finite element algorithm code.They suggested that the change of the properties of buffer material led to the development of transportation channel.Du et al. [12] applied flow model for the heterogeneous unsaturated zone to DECOVALEX Task D and analyzed the sensitivity of parameters.Gerard et al. [13] studied the numerical model of gas migration under isothermal condition and analyzed the migration process of hydrogen caused by the corrosion of steel container in engineering of long-term storage of radioactive waste.These studies show that the migration of gases, such as sedimentary rocks, is a complex process.However, these modeling studies do not take into account the uncoupled processes of hydraulic and mechanical processes or the major interactions involved.
The objective of this paper is to establish a hydromechanical-damage model to predict the gas migration in the sedimentary rock.The model considering the characteristics of gas migration, especially the microcracks caused by the degradation of elastic modulus and damage, the effect of damage on permeability, and the relationship between the damage and the rock deformation and failure of fractures.Based on this model, we can simulate the gas migration process and study the mechanism of gas flow and its influence on rock stability in sedimentary rocks.

Hydro-Mechanical-Damage Model
In order to establish the governing equation, we have made the following main assumptions: (1) The rock stratum is continuous and uniform.
(2) The gas flow process in rock stratum is pseudostatic.
(4) In terms of mechanical response, the deformation is small and the strain is infinitely small.

Stress Equilibrium Equation.
According to the theory of porous elasticity, the unit of rock satisfies the following equilibrium equation: where   is the total stress component of the rock body unit;   is the body stress in the direction ; and   is the direction coordinates in the direction .
The rock is regarded as a porous medium, and the rock element satisfies the constitutive equation.It can be expressed by stress, strain, and pore pressure as follows: where   is the symbol of Kronecker;  is the Biot coefficient of rock,  = 1 − /  ;   is the bulk modulus of rock matrix;  is the volume modulus of rock;   is the component of strain tensor; and   is the component of stress tensor.According to the continuous deformation condition, the following geometric equations are obtained: After the rock adsorbs gas, the adsorption expansion strain can be expressed where   and   are the Langmuir strain constant and pressure constant, respectively.The stress equilibrium equation can be expressed by displacement, pore pressure, and adsorption expansion.

Gas Seepage Equation.
Darcy flow is widely employed in the gas migration process.The Darcy velocity of gas is expressed as where   is the coefficient of dynamic viscosity and k is the permeability of the gas.The seepage of gas follows the law of conservation of mass.
where  is the unit volume for the gas in the rock and t is the time variable.The mass of the gas  is composed of free term and adsorption term, which can be expressed as where   is the gas density at standard condition;   is the density of the rock;   is the Langmuir volume constant; and  is the porosity of the rock.Because of the compressibility of the gas, the relation between the gas density and the pressure is The continuity equation of gas seepage can be obtained as

Permeability Model of Rock.
The basic skeleton of rock is deformed under the affection of gas pressure, which changes the porosity of rock and affects the seepage of gas in rock.The rock is subjected to the double action between external stress and pore pressure.The following equation can be obtained: where    is the effective stress,  is the Biot coefficient of the rock,  is the pore pressure, and   is the Kronecker symbol.
The rock body is regarded as the porous medium.The volume of the rock  is the summation of pore   and matrix   , so the porosity is defined as  =   /.The volume deformation of the rock matrix consists of two parts: the effective stress and the matrix adsorption.Therefore, the porosity of the rock can be expressed as The permeability of the rock body is related to the porosity, which can be expressed as cubic law.
The permeability of rock can be expressed as The following equation can be obtained after finishing transposition: where  0 is the initial permeability of rock;   is the compression coefficient of fracture,   = 1/  − 1/;  is the volume modulus of rock;   is the volume modulus of rock fractures; and   is the effective stress.

Parameter Heterogeneity of Rock.
In order to describe the heterogeneity of rock materials, it is assumed that rock consists of a large number of microscopic elements.Assuming that the mechanical properties of these units obey Weibull distribution, the distribution can be defined according to the following density distribution function: where  satisfies the numerical value of the Weibull distribution function,  0 is a parameter related to the average value, and  is the shape coefficient.
The greater the parameter , the better the uniformity of the material unit, and vice versa.Therefore,  0 and  are called the distribution parameters of materials.Using (16), the inhomogeneous parameters of the rock materials can be generated in the numerical calculation.These parameters are closer to the true sample parameters in the laboratory test.

Analysis of Damage Theory.
The maximum tensile stress criterion is used to determine the tensile damage of rock, and the Mohr-Coulomb criterion is adopted to determine the shear damage of rock [16][17][18], as shown in Figure 2. where  1 and  3 are the maximum principal stress and the minimum principal stress of rock units and   and   are uniaxial tensile strength and uniaxial compressive strength of rock unit, respectively.Based on the strain, the damage variable of rock units can be expressed using the following expression: where   and   are the maximum principal strain and the minimum principal strain and  0 and  0 are the ultimate strain corresponding to tensile damage and shear damage, respectively.
The elastic modulus of rock under damage state can be expressed as follows: where  0 is the elastic modulus of undamaged state and  is the elastic modulus of the unit in the damaged state.When the rock is damaged, the effect of the rock damage on the permeability can be described as where  0 is the initial permeability,   is the compression coefficient of the rock fracture,  is the influence coefficient of damage to permeability, and   is the effective stress.
The coupled hydro-mechanical-damage model for rock is proposed and COMSOL Multiphysics and MATLAB are employed to achieve the coupled solution of solid field, fluid field, and damage field.

Model Establishment
The in situ gas injection test was carried out at the Mont Terri underground laboratory in the Ru La mountains, northwest of Switzerland [14,15].In the area of the laboratory, the covering layer is between 230 and 320 meters.The in situ gas experiment was conducted in the EZ-A niche.This niche is located in the sedimentary rocks and was unearthed in 2003.Figure 3 shows a schematic diagram of a field gas injection test.
Based on the engineering geological characteristics, a calculation model is established.The size of the numerical model is 10m×10m, and the diameter of the borehole is 0.1m. Figure 4 shows the model and the boundary condition.The modulus and strength of the model unit are subjected to the Weibull distribution and the parameters of the rock are shown in Table 1.The normal displacement is applied on the left and the lower sides of the model, and the external boundary stress is applied on the right and the upper sides.Based on the established calculation model, the hydromechanical-damage coupling procedure was compiled.The damage and seepage evolution characteristics of surrounding rock were studied after drilling and gas injection.

Analysis of the Damage Characteristics.
The seepage distribution and the damage distribution of the rock around the borehole can be obtained through the analytical solution in previous literature [19].The analytical solution can be used to verify the validity of this model.Before the peak stress intensity of rock is reached, there is no damage.It is the linear isotropic damage evolution when rock enters the postpeak stage.The damage evolution equation is expressed as where  0 is the borehole radius;   is the radius of damage zone; and  is the distance from rock to borehole center.Because the excavation damage zone (EDZ) is mostly concentrated near the excavation wall, we mainly focus on the area of 0.5m×0.5maround the borehole.It can be seen in Figure 5, where a circular excavation damage zone appears in the surrounding rock after drilling.When it is closer to the borehole wall, the damage degree of the rock is higher and the stress release degree is higher.When it is far from the borehole wall, the influence of the excavation is smaller, and there is a gradual transition to the no damage area.For the soft rock, because of the high degree of damage of the borehole surrounding rock, the hole collapse often occurs, and it is necessary to reinforce the borehole.Figure 6 is a comparison between numerical solution result and analytical solution result of damage distribution.The numerical solution results are the damage values of each unit of rock on the monitoring line of finite element model.Because of the inhomogeneity of the initial mechanical parameters of rock, the damage value of rock units also presents a discrete distribution state.Although the damage of rock units on the monitoring line shows discrete distribution characteristics, it can be seen that damage shows a decreasing distribution with the increase of the distance to borehole.The calculation results of the hydromechanical-damage coupling model are in good agreement with the analytical results, and the correctness of the model is verified.
(a) Geological section of the Mont Terri Underground Rock Laboratory [14] (b) Schematic presentation of field gas injection test [15] Figure 3: Location of rock laboratory and field gas injection test.The damage degree of surrounding rock is relatively high, and the elastic modulus generally decreases.Although the initial elastic modulus of rock is distributed randomly, the damage degree of rock is higher in the range of distance to the borehole wall 40mm.The distribution of elastic modulus gradually changes from random distribution to regular distribution.rock is seriously damaged, and the permeability coefficient is high.The excavation disturbance produces effective pressure relief and permeability enhancement for rock.The maximum permeability around the borehole reaches 1.07 × 10 −15 m 2 , which increases by two orders of magnitude compared with the initial permeability of rock.The permeability in the red area is high, and the damage degree is high, which is a high permeability area.After the drill is excavated, the gas can quickly pass through the high permeable area, and the high permeability area also improves the gas migration efficiency.With the increase of distance to the borehole wall, the permeability gradually decreases.Although the damage degree is small, the damage still improves the permeability of rock and accelerates the seepage of gas in the rock.Figure 9 shows the comparison between numerical solution results and analytical solution results of permeability distribution of surrounding rock.The results are the permeability values of each unit of rock on the monitoring line of finite element model.Because of the inhomogeneity of the initial mechanical parameters of rock, the permeability of rock units also presents the characteristic of discrete distribution.It can be seen from Figure 9 that the calculation results of the hydro-mechanical-damage coupling model fit well with the calculated results, and the correctness of the model is verified.

Analysis of Percolation Characteristics of Borehole
Figure 10 shows the gas pressure distribution of the rock mass.The gas in the damaged area can flow into the rock quickly, and the gas pressure increases rapidly.Subsequently, the gas pressure in the undamaged area also gradually increases, and the gas pressure in the damaged area remains at a high state.The gas pressure gradient at the boundary of the damage zone is larger, consistent with the seepage law of the elasto-plastic analytical solution.After that, the gas pressure of the farther area outside the damaged area gradually increases.

Influence of Lateral Pressure Coefficient on Damage
Distribution and Seepage Characteristics.The damage distribution characteristics and gas seepage characteristics of rock under the condition of equal horizontal stress and vertical stress are analyzed in the above section; that is, the lateral pressure coefficient is equal to 1.However, under some geological conditions, the rock stratum is not in the hydrostatic pressure state and the influence of the lateral pressure coefficient should be considered.Figure 11 shows the damage distribution of the surrounding rock when the lateral pressure coefficient   = 1.5.The distribution of damage changes from circular shape to oval shape.Because the stress in horizontal direction is larger, the maximum principal stress on the upper and lower sides of the borehole increases, and it is more likely to reach the ultimate strength, resulting in the occurrence of the damage.The nonisotropic stress causes the damage distribution of surrounding rock to be elliptical and no longer presents a circular distribution.Figure 12 is the permeability distribution of the borehole surrounding rock when the lateral pressure coefficient   = 1.5.The distribution of permeability is consistent with the distribution of damage.The permeability maximum increases by up to two orders of magnitude compared to the initial value.The permeability in high damage zone increases obviously, while the permeability in low damage area increases several times.The area of permeability region is slightly less than the damage area.
Figures 13 and 14 are the damage distribution and the permeability distribution of the surrounding rock when the lateral pressure coefficient   = 3.When the lateral pressure coefficient   = 3, the damage distribution and permeability distribution of surrounding rock appear in butterfly shape.The damage develops along the 45 ∘ direction with the principal stress direction, which coincides with the theoretical research result of the plastic zone of rock mass around borehole surrounding rock.
In practical engineering of underground nuclear waste storage, the damage will appear around the excavation cavity, which will seriously affect the strength and seepage characteristics of the rock mass.The permeability of rock in excavation damage zone normally increases by 2-3 orders of magnitude compared to the initial value.The gas will quickly flow through the excavation zone, while the flow will be relatively slow in the undamaged zone.Although a lot of scholars have done research on the gas seepage process [20][21][22][23][24][25][26][27][28], including the influence of chemistry and temperature on seepage, many studies neglected the effect of damage, which may cause errors in the analysis of gas seepage characteristics.In addition, the nonuniform in situ stress and bedding of rock layers will affect the damage and seepage characteristics of the rock mass.When the geostress is nonuniform, the shape of the damage zone is no longer circular.This is also helpful to the research on gas migration in sedimentary rocks.Reinforcement of surrounding rock in damaged area can effectively prevent gas seepage.

Conclusions
The coupled hydro-mechanical-damage model is established to predict and analyze the gas migration in rocks, which   Based on the hydro-mechanical-damage coupling model of rock mass, the damage and fracture characteristics of surrounding rock and seepage characteristics in rock mass are analyzed, and the numerical solution results are compared with the analytical solution results to verify the validity of the model.In addition, the influence of the lateral pressure coefficient on the damage and the permeability of rock are  also considered.The results indicate that a circular excavation damage zone appears in the surrounding rock after drilling, and the damage degree of rock is higher within the range of 40 mm from the borehole wall.The closer to the borehole wall, the higher the damage degree of the rock.When the lateral pressure coefficient   = 1.5, the distribution of damage changes from circular shape to oval shape.When the lateral pressure coefficient   = 3, the damage distribution and permeability distribution of surrounding rock appear in butterfly shape.

Figure 2 :
Figure 2: The constitutive law of rock.

Figure 4 :
Figure 4: Numerical model of gas seepage process of borehole surrounding rock.

Figure 5 :
Figure 5: Damage distribution of borehole surrounding rock.

Figure 7
Figure 7 displays the distribution of the elastic modulus of surrounding rock.The distribution law of elastic modulus of rock is generally similar to the distribution of damage.The damage degree of surrounding rock is relatively high, and the elastic modulus generally decreases.Although the initial elastic modulus of rock is distributed randomly, the damage degree of rock is higher in the range of distance to the borehole wall 40mm.The distribution of elastic modulus gradually changes from random distribution to regular distribution.

Figure 6 :
Figure 6: Comparison between numerical results and analytical results of damage.

Figure 10 :
Figure 10: Gas pressure distribution of borehole surrounding rock.

Figure 11 :
Figure 11: Damage distribution of surrounding rock under lateral pressure coefficient   = 1.5.

Figure 12 :
Figure 12: Permeability distribution of surrounding rock under lateral pressure coefficient   = 1.5.

Table 1 :
Parameters in the calculation.