Analysis of the creep and dilatant behavior of a salt cavern in long-term using Brazilian geotechnical properties

Abstract The geomechanical behavior of a salt rock cavern was studied using Brazilian geotechnical properties. To design the finite element model, it was necessary to implement, using Fortran language, a constitutive model representing the creep behavior, since the model used is not native to the program. The constitutive model implemented was the Multiple Deformation Mechanism Model. This model was chosen for being a robust model that represents the primary and secondary phases of creep and presents good agreement with the Brazilian salt rocks. The analyzes considered a period of 50 years after the mine closure and five internal pressures acting in the analyzed cave. The pressures considered correspond to 40%, 50%, 60%, 70% and 80% of the vertical stress at the top of the mine. The creep and dilation behaviors were analyzed, and the creep deformations obtained in the simulations was acceptable in relation to the failure criterion adopted in this paper. However, only the design of experiment that considered the two biggest internal pressure resulted in a permissible micro-crack ratio value.


Introduction
Due to the products in its chain, chlorine and caustic soda are very relevant for society.Examples include polyvinyl chloride (PVC), agrochemicals, products for water and sewage treatment, pharmaceutical inputs, and lubricants (Jörissen, 2014).
Chlorine-soda is obtained by the electrolysis of a solution of rock salt (halite) in water.One of the rock salt mining techniques is dissolution extraction using recovery wells.The wells must reach the depth of saline material that can reach thousands of meters below the earth's surface.Afterwards, water is injected by these wells and when it comes into contact with rock salt, it dissolves the rock.Thus, it is possible to extract the brine used to obtain chlorine-soda.Curi (2017) presents a more detailed explanation about the extraction process.
Inevitably, at some point, the chlorine-soda extraction ceases.The most common technique of abandonment is to fill the cavern with brine and seal the well permanently.Thus, the pressure applied by the brine maintains the structural integrity of the cavern (Crotogino & Kepplinger, 2006).However, failures in the well seal, leakage through micro-cracks in the walls of the cave, tectonic movements, and other factors, can cause loss of internal pressure in the cavern.Recent researches have been studied the possibility of using abandoned caverns for the storage of hydrocarbons (Thoraval et al., 2015;Wei et al., 2016a, b;Zhang et al., 2020a, b).
The study of cave abandonment involves three important aspects.The first aspect is the experimental characterization of the saline material.The second aspect is in-situ tests to investigate mine geometry.The third aspect is the applying constitutive models to simulate the material behavior.These models can predict the time-dependent behavior of the material (Thoraval et al., 2015).Some evaporitic rocks, especially halite, have a very particular behavior.The material suffers permanent deformations under constant loads and/or high temperatures.These deformations change throughout time.This timedependent behavior is known as creep (Costa, 1984).In the case of caverns, built on salt rocks, the material that surrounds the cavern flows into the cavern.This causes mass movement that reduces the cavern volume.As consequence, the mass movement can cause lowering of the terrain surface (subsidence).Subsidences and instabilities in salt-caverns regions are old geotechnical problems and found in several countries (Bell et al., 2000;Autin, 2002;Whyatt & Varley, 2008;Yerro et al., 2014;Wang et al., 2018;Vassileva et al., 2021).
Creep occurs in three stages.The primary creep is where the strain rate decreases.The secondary creep follows the primary phase, with a constant strain rate.Finally, the tertiary phase.Literature does not characterize this phase well.But, in this phase, there may be an acceleration of the strain rate followed by the material rupture (Poiate Junior, 2012).
The dilatancy behavior of the saline rock is very important.A high level of stress can cause the generation of micro-cracks (Firme et al., 2019).Two domains define the salt rocks behavior.The first domain is compression.
Here, the micro-cracks heal and guarantee the material impermeability.The second domain is dilatancy.The salt rock enters this domain when the active deviatoric stress exceeds the limit deviatoric stress.In this case, the expansion of the material can increase the permeability of the material (due to the occurrence of micro-cracks) and even cause damage (Van Sambeek et al., 1993).
Both dilatancy and creep behaviors can cause instabilities in abandoned caverns.This can lead to instabilities on the ground surface, such as subsidence.Although it is difficult to control the variables of a mine-abandonment project, due to its depth (it can be thousands of meters), maintaining an adequate internal pressure can help to avoid geotechnical instability.Therefore, the objective of this paper is to analyze the influence of internal pressure on the creep and dilation behavior of a rock salt mine with a geometry similar to those existing in the chlorine-soda mining fields in Brazil.

Constitutive modeling of rock salt
The constitutive relation applied to describe the rock salt behavior follows the Multiple Deformation Mechanism developed by Munson & Dawson (1979).The Multiple Deformation Mechanism allows simulating the creep behavior of rock salt under deviatoric stress.This model simulates transient (primary) and steady-state (secondary) creep phases.Different factors such as stress state, temperature, and material chemical composition, influence dislocations in salt rock (Firme et al., 2016).The Multiple Deformation Mechanism model uses three mechanisms to determine the steady-state phase creep: dislocation climb, undefined mechanism, and dislocation glid.The sum of them quantifies the secondary creep rate.
The dislocation climb can be described as a mechanism of thermal activation that depends on the stress intensity.This mechanism has a meaningful significance at high temperatures (approximately 50 °C).The high deformation energy generated during the dislocations leads to plastic deformation and hardening.The climb mechanism allows the movement of the dislocations (Firme et al., 2016).
An Arrhenius expression determines the contribution of dislocation climb on steady-state creep rate ( ) In rocks, such as halite, the contribution of the dislocation climb is considerable when the temperature is above 40% of the melting temperature of the material.The second mechanism is known as undefined mechanism and its rate ( ) UM ε is determined using an Arrhenius expression as the one used to determine the dislocation climb, defined as: The third mechanism used is the dislocation glide ( ) DGL .It is a slow creep mechanism and arises from displacements in the material crystal lattice.This mechanism occurs at any stress level unbalanced by deviatoric stress.Shear stress also can influence the occurrence of this mechanism.Consequently, there is plastic strain and material hardening.The rate of this mechanism ( ) In this equation, ( ) − is a Heaviside step function with argument ( ) . This argument limits the occurrence of DGL for a deviatoric (or equivalent) stress higher than the reference stress value of the mechanism ( ) The primary creep was determined by performing a retro analysis with the secondary creep rate.Thus, it is possible to estimate a limit strain during the transient creep ( ) ε it is possible to determine the function F as: In this equation, ς is an internal isotropic hardening variable subject to an evolutionary rate ( ) ς defined as: The value of SS ε is obtained by the summetion of the three micromechanical creep mechanisms presented in Equations 1, 2, and 3.
Finally, the total creep rate ( ) ε is determine, including the primary and secondary creeps, by using: After integrating ε, it is possible determine the total creep strain ( ) ε as:

Validation
The Multiple Deformation Mechanism model is very effective to simulate the salt rock behavior.Munson (1979) presents accurate results for simulations involving experimental and field situations.Thus, this paper does not prove the model efficiency, which is very used and established.The purpose of validation is to verify that the model implementation using the Fortran language was performed correctly.This implementation is necessary as this model is not native to the software.
Three triaxial creep tests was simulated to validate the implementation.During the validation simulation the samples was kept under confining stress of 10 MPa and deviatoric stress of 10, 14, and 17 MPa.Poiate Junior et al. (2006) presented the test with deviatoric stress of 10 MPa, while Costa et al. (2005) presented the other tests.
All tests simulated was performed using cylindrical body tests of 88.9 mm (3.5 inches) diameter and 177.8 mm (7.0 inches) length.The aspect ratio of 2 for 1 follows the recommendation of the International Society for Rock Mechanics.During the test, the temperature was controlled and kept at 86 °C.
Figure 1 presents both experimental and numerical results obtained for the samples.Table 1 presents the parameters used to perform the numerical simulations.Firme et al. (2016).calibrated the parameters using tests of Sergipe-Alagoas-Brazil basin.Adequated agreement was observed between numerical and experimental results.
To characterize the elastic phase of the material was adopted the values presented by Poiate Junior (2012).
The author used 20.4 GPa and 21.29 3 / kN m to describe the Halite's dynamic elastic modulus, and unit specific weight.Finally, we used the dynamic Poisson's ratio equal to 0.36 as presented by Costa et al. (2005).

Stratigraphy and typical geometry
The geometries of chlorine-soda mines are very variable due to their construction techniques and the material characteristics.The cavern studied in this paper has geometry and dimensions similar to existing mines in northeastern Brazil (CPRM, 2019).Figure 2 shows the geometry of the mine studied and its dimensions obtained by sonar measurement.
The stratigraphy used is based on the characteristics of the Sergipe-Alagoas sedimentary basin, also located in northeastern Brazil.In this basement, the rock salt stratum starts at approximately 900 m of depth and the overburden is composed mainly of shale and sandstone, and there may be the presence of other sedimentary rocks.
In this paper, it was considered the most recent geometry obtained by sonar measurement performed in 2018.However, the geometry was simplified by replacing the original format by an ellipse, with height and radius similar to the original geometry.Figure 3 shows the simplified geometry.
The simulation, despite being based on data of parameters and geometry of Brazilian mines, is not representative of any specific case, it is a hypothetical situation that aims to analyze the behavior for conditions similar to that simulated.
Simplifications are common in studies of underground salt caverns (Goulart et al., 2020;Li et al., 2021).They bring the advantage of avoiding the action of stress concentrators due to the corners.Another important point is the quality of the mesh used to discretize the model.A simpler geometry avoids the use of a greater number of elements in regions with sharp corners, for example.Furthermore, discretization can be performed with more robust elements.This provides an improvement in determining the variables involved in the simulation.

Geomechanical model, boundary, and initial conditions
Figure 3 presents the model used for performing the simulations.At the top of the model was applied a 19.03 MPa pressure, which characterizes the overburden.This pressure represents the sediments above the shale layer.The height of the mine is 96 m and the approximate radius is 55 m. Figure 3 also presents the others measurements of the model.
The finite element model was built using the commercial software Abaqus.An axisymmetric model represented the saline environment.It was used CAX4R elements (4-node bilinear axisymmetric quadrilateral, reduced integration) to discretize the model.Two materials composed the model: halite and shale (Figure 3).We applied the elastoplastic model proposed by Mohr-Coulomb to the shale's behavior.It was considered the elastic modulus, Poisson's ratio, friction angle, cohesion and the specific weight equal to 18.97 GPa, 0.15, 22°, 4.8 MPa and 24.00 kN/m³ (Costa et al., 2015).The geothermal gradients, i.e., the temperature rate of geotechnical materials, suggested by Costa et al. (2012) were used.The temperature is higher the greater the depth, these authors suggest 10ºC/km for saline material and 30ºC/km for the sedimentary rocks that constitute the overburden.
Figure 4a presents the mesh discretization around the cavern and the five points of interest (A, B, C, D, e E).These points were used to perform some analysis presented in section 4. Figure 4b presents the mesh and boundary conditions.The thick black line indicates the axis of revolution to build the model.Figure 4b represents a 90-degree revolution of Figure 3.However, the software understands the model as a complete 360-degree revolution.The blue and gray portions represent shale and halite.
On the right side, the horizontal displacement was constrained.This restriction represents the continuity of the layer and allows the vertical movement.Vertical displacements were constrained at the model bottom.
The simulation was performed in three steps.In the first step, the geostatic stresses are generated, and the excavation of the mine is carried out.The coefficient of earth pressure at rest of overburden was considered 0.25 and for salt rock was considered 1, cause your creep behavior.The excavation is carried out instantly with the complete removal of the material.The second step is the instantaneous elastic behavior of the material.Finally, the third step if visco-elasto-plastic behavior using the Multiple Deformation Mechanism model concerning a period of 50 years.
Given the fact that some units of the constitutive model use seconds as a time unit, the simulations were performed considering seconds as the time unit.The time for data processing of each simulation is directly linked to the internal pressure considered.Higher internal pressures are more stable conditions.So, they are more agile, while smaller internal pressures need more time to be completed.The simulations performed in this research took between 10 and 30 minutes using an Intel Core I7 computer with 8GB of RAM.

Numerical experiments and stability criteria
Five simulations were performed by considering different internal pressures applied on the cavern (Table 2).The values of the internal pressures are a percentage of the effective vertical tension at the top of the cavern.Those values ranging from 40% to 80% of the value of the effective vertical tension.
The analyzed internal pressures were based on literature to choose, which indicates values between 30% and 85% of the effective vertical stress at the top of the cavern.Pressures less than 30% cause instability in the caverns.While pressures greater than 85% can generate tensile stress.Salt rocks do not resist this type of effort (Bérest et al., 2020;Costa et al., 2015;Yuan et al., 2021).Stability criteria considered both the creep and dilatation behaviors.

Creep criterion
The creep behavior causes a decrease in the volume of salt caverns.Thus, the creep criterion relates to the cavern volume loss.The largest shrinkage volume limit was defined as 10% of the initial cavern volume throughout the analysis (50 years).This is a value often used in salt cavern stability analyses (Yuan et al., 2021).

Dilatancy criterion
In this research, the tertiary (accelerated) stage of creep was not considered.However, an expansion criterion monitored the damage in the salt rock.Salt rocks in situ are usually impermeable.But, when an increase in deformation happens, there is a rise in permeability of the rock due to damage (Van Sambeek et al., 1993;Firme et al., 2016).Ratigan's criterion (Van Sambeek et al., 1993) was used, since this criterion is widely used to model the Brazilian halite (Costa et al., 2015;Firme et al., 2016;Firme et al., 2019).This criterion considers the mean stress ( )   The deviatoric effective stress ( ) d σ considered in this paper is the equivalent von Mises stress ( ) that can be described from the second invariant of the deviatoric stress tensors ( ) 2 J .

Salt cavern convergence
As highlighted above, the creep behavior is responsible for the convergence of salt caverns.That is, the volume of shrinkage of the cavern.The strain at each point of the model depends on the acting deviatoric stress (von Mises equivalent stress), as shown in Equations 1, 2, and 3. Thus, the shrinkage volume is dependent on the von Mises equivalent stress.
The Figure 5 shows the von Mises equivalent stress for IP40.One can observe that the cap rock causes a stress arching effect.This relieves the stresses in the saline material.Thus, cap rock is very important for cavern simulations.Zhang et al. (2020a) also verified the importance of the surrounding rocks.
Figure 6 shows the equivalent creep strains for all internal pressures presented in Table 2.The point C was presented the highest level of strain among five points.
Figure 7 shows the wall equivalent strain rate at point C for all analyzed internal pressures.It was analyzed this point because it is at the region with the highest strain.Costa et al. (2015) indicates an admissible creep rate in steady-state regime equal to 0.5x10 -06 per hour.One can note that the strain rates (per hour) are between 10 -07 and 10 -10 at the end of analysis.Instantly the creep strain rate has values greater than the admissible value for all internal pressures.However, in a small portion of time, the creep rate becomes higher than the admissible value.Afterwards, one can observe the deceleration phase of the strain rate until reaching the steady-state.There is an unstable condition for the IP80.However, the creep rate presents the smallest values at the end of simulation.
Figure 8 shows the wall equivalent strain at point C. In the numerical analysis performed, the maximum strain was 0.038 for IP40.The simulations did not reach the maximum value allowable considered equal to 0.100 (Costa et al., 2015).This indicates that shrinkage caused by salt creep is not meaningful even for IP40, which presents an increasing behavior during the simulation.
Figure 9 shows the percentage of volume shrinkage over time.None of the numerical experiments reached the maximum shrinkage volume of 10% of the initial volume.
Figure 10 presents the numerical results of displacements norm resulting from creep strain.Clearly, displacements decreased when the internal pressure increases.
A large displacement increase was noted when the internal pressure was reduced from 80% to 70% of the in situ vertical stress on the mine roof.Bérest et al. (2020) performed several studies on the internal pressure of salt caverns.They indicate that the adequate pressure is between 80 and 85% of the in situ vertical stress on the mine roof.However, it is not always possible to ensure this pressurization, due to several factors such as reactivation of geological faults, cracks due to the low tensile strength of the saline material or, even, due to leakage due to the filling of the well.

Dilatant behavior
Figure 11 presents the micro-cracking (MC) ratio obtained for each internal pressure considered.One observes that the regions close to points B and C provide the necessary condition for initiation the micro cracks.
At the end of the 50 years analyzed by finite element modeling, the micro-crack ratio was acceptable only for IP70 and IP80. Figure 12 shows the micro-cracking ratio values, considering 50 years of analysis, for point C for all simulated internal pressures.In this sense, the knowledge of this complex behavior is fundamental.The use of Multiple Deformation Mechanism model combined with a suitable dilatation criterion allows modeling the primary creep, the steady-state creep and damage

Conclusions
The high level of stress due to the creep phenomenon in salt rocks can cause the failure of abandoned caverns. of a salt cavern.Thus, one can perform numerical studies by considering the modifications of the cavern's geometry or the internal pressure.
Simulations was carried out, using the finite element method, to analyze the behavior of caverns located in saline environments.Creep and dilation were investigated.Despite the occurrence of displacements with high magnitude (greater than 1 m), the creep deformation for all analyzed internal pressure was considered acceptable, as the shrinkage volume was smaller than the maximum allowable shrinkage volume.Another important aspect simulated was the micro-crack ratio.It was considered as an acceptable the value of 0.8.Values greater than this indicate the occurrence of micro-cracks in the cavern walls.This paper indicates that only the IP70 and IP80 condition provided acceptable results.
Several research have been carried out in order to analyze the behavior of salt-rock caverns, but there are difficulties in representing the behavior of this material over time.It is necessary to emphasize that the combination of a model that simulates primary and secondary creep with a damage criterion, as done in this research, is one of the most robust ways of analyzing rock salt behavior over time.
Finally, the integrated analysis between creep and dilatant behavior was considered an important tool for reliability analysis, as it allows an accurate analysis of the behavior of rock salts.Future research should focus on analyzing the variability of geotechnical parameters to also provide a probabilistic approach.

Figure 1 .
Figure 1.Comparison between numerical results and experimental results of triaxial tests on halite.

Figure 2 .
Figure 2. Shape and dimensions of cavern obtained by the sonar survey in 2018 (CPRM, 2019).

Figure 3 .
Figure 3. Geometry used for representing the cavern.
a receives the value of 0.81.Damage will occur if the ratio between the deviatoric ( ) d σ effective stress and the critical deviatoric stress ( ) , d cr σ is greater than the critical value indicated by the criterion.This value corresponds to 0.8, i.e., the ratio between the acting deviatoric effective stress and critical deviatoric stress is the micro-crack ratio ( ) MC defined as:

Figure 4 .
Figure 4. Axisymmetric model and finite element mesh discretization: (a) detail of the analysis points; (b) partial revolution of the analysis plane.

Figure 6 .
Figure 6.Equivalent creep strains (CEEQ) around the cavern for the internal pressures studied.

Figure 7 .
Figure 7. Wall equivalent strain rate at point C.

Figure 8 .
Figure 8.Wall equivalent strain at point C.

Figure 9 .
Figure 9. Percentage of volume shrinkage over time.

Figure 10 .
Figure 10.Magnitude of displacements at the analyzed points.

Figure 11 .
Figure 11.Micro-cracking ratio (MC) around the cavern for the internal pressures studied.

Figure 12 .
Figure 12.MC ratio at point C after 50 years of analysis.Dashed line represents the dilatancy criterion critical value.

Table 1 .
Parameters used to characterize the creep behavior of the rock salt.

Table 2 .
Design of experiments.