Investigation of the dissipative structures following microseismic diffusion during hydraulic fracturing of methane-hydrate-bearing sand

Hydraulic fracturing is a prospective technology for methane hydrate deposit exploitation. The evolution of hydraulically stimulated fractures around the point of liquid injection is simulated. For this purpose, the FLAC3D computer model is used because of its explicit calculation cycle that imitates real physics, prevents numerical instability, and reproduces a realistic path during simulation of the nonlinear rock massif behavior. The results of the simulation provide for new findings, namely, the spatial asymmetry and synchronism violation, spatial deviation, discontinuity, and recurrence during microseismic diffusion, which follow the process of hydraulic fracturing. In addition, dissipative structures were developed due to entropy production, since gas hydrate strata are an open thermodynamic system, which transforms and dissipates the energy of the injected liquid. The process of dissipative structure evolution should be controlled to enhance the gas yield from the hydrates.


Introduction
Methane hydrate is deposited in permafrost regions and deep oceanic environments. The global resource of methane in the gas hydrate exceeds the other hydrocarbon reserves. Therefore, methane hydrate can play a dominant role in the global carbon cycle in the future [1]. Na et al. [2] reviewed numerous natural gas hydrate (NGH) production technologies. Thermal stimulation, depressurization, implementation of a chemical inhibitor, and CO 2 -CH 4 exchange are compared. These conventional technologies have certain advantages and disadvantages [3,4]; however, they do not satisfy economic requirements and environmental conditions. That is why specialists are searching for new innovative technologies to exploit gas hydrate deposits. Sasaki et al. [5] presented a system of hot water injection with a pair of horizontal wells to exploit the NGH. The technology of hot water cyclic steam stimulation with a single well that has been traditionally used to improve the recovery effect in a superheavy oil reservoir [6] proved to be effective for NGH production augmentation. Other conventional in the oil industry technologies such as partial oxidation, electrical heating assisted by depressurization, and CO 2 swapping depressurization promised notable progress in the development of NGH deposits [7][8][9][10].
The hydraulic fracturing process (HF) has become very popular among practitioners who extract resources from hydrocarbon and hydrothermal deposits [11][12][13]. This technology involves a high-pressure fluid injection into geologic deep strata that creates cracks in the deep-rock formations to stimulate production from the wells. Recently, HF has been successfully introduced for the NGH production stimulation [14]. Konno et. al [15] investigated HF as a well stimulation method during gas recovery from gas hydrate reservoirs. They injected distillate water into methane-hydrate-bearing sand, which was in a three-dimensional confining stress state. Pore pressure has increased rapidly, but suddenly dropped due to the delamination of the tested volume. Hydraulically stimulated fractures are developed under the action of tensile stress. The permeability of the stimulated volume increased after fracturing and was maintained even after re-confining and closing the fractures. This indicated that HF is a promising method for the well stimulation of the lowpermeable gas hydrate reservoirs.
However, requirements claimed for the stimulation of a NGH deposit are much stricter than for the activation of other unconventional gases. The methane-hydrate-bearing strata should be stimulated more uniformly to activate the process of dissociation. On the other hand, discrete fracture networks (DFN) evolve mostly stochastically and are almost not amenable to control. Is it possible to produce uniform fracturing of a body, which was stimulated by HF? Can the HF process be effectively controlled? We try to answer these questions in this paper.

Mathematical representation and description of using the airflow energy
Computer models are the best tool to investigate new geomechanical processes greatly augmenting the development of modern rock mechanics. Jing [16] classified existent geomechanical models to finite element, boundary element, finite difference, and discrete element methods. Lei et al. [17] discussed the state-of-the-art on the use of discrete fracture networks (DFN), which are relevant for modeling structural characteristics, geomechanical evolution, and hydro-mechanical behavior of stochastically generated fracture networks in a rock mass. However, the stochastic nature of DFN introduces uncertainty to fracture distribution and evolution. Consequently, DFN models can hide important fracture features. To achieve the goal, the FLAC3D model was chosen, which can simulate an irreversible behavior of the rock mass stimulated by HF [18,19]. FLAC3D simulates flow in parallel with the mechanical modeling, in order to capture the effects of fluid/solid interaction. HF generates microseismicity characterized by a seismic moment [20], which is defined by the equation: where G -the shear modulus of the rocks involved in the micro-earthquake, Pa; A -the area of the rupture along the fracture, produced by HF, m 2 ; D -the average slip or displacement offset between the two sides of the fracture, m.
The seismic moment has an energy dimension, and transforming the formula, we use the incremental deformation of the volume of the rock surrounding the crack: FLAC3D cannot simulate a set of fractures explicitly. Thus, we may simplify Formula (2) by eliminating the factor V. Let us consider this expression as a specific indicator of seismic activity or specific seismic moment, namely total seismic moment divided by volume of the zone where the incremental strain was calculated by FLAC3D.
The process of fracture initiation is simulated using Hubbert and Willis formula [21]: where P b -breakdown pressure; σ h and σ H -minimum and maximum horizontal in-situ stresses, respectively; T -the tensile strength of a rock layer; p o -the pore pressure; η -the poroelastic parameter that varies in the range of 0 to 0.5. Recalculated permeability of the rock volume after fracture generation according to recommendations of Min et al. [22] reads: where k 1 and k 1 ′ -permeability along direction 1 before and after HF; ΔS 1 -the normal component of strain increment; ΔS 12 and ΔS 13 -shear components of strain increments in planes 12 and 13; f norm and f tan -empirical coefficients 0.83 and 0.095, respectively.
The other components of orthotropic permeability k 2 and k 3 were recalculated replacing corresponding indexes.

Results and discussion
Initial data. Boundary conditions.  The model was divided into 9456 zones containing 12829 grid nodes. The geometry of the model was symmetrical relatively all axis. The displacements that are normal to lateral walls were fixed. All components of displacements were prescribed to zero on the bottom. The initial stress state in the model corresponded to the depth of 800 m at the level of the liquid injection: components of the geostatic stress σ y , σ x , and σ z were 20, 10, and 10 MPa respectively in the point of injection.
Orthotropic rock mass had three basic fracture systems, namely two sub-vertical cleavages and one sub-horizontal lamination of sedimentary rocks [23]. Axial symmetry divided the model into two mirror symmetries, relatively horizontal or X-and Z-axis. It does mean that the hydraulic fractures were expected to expand synchronously in time and symmetrically in space relatively the corresponding horizontal axes, thereat along the X-axis quicker because permeability k x was five times more than k z . Mechanic and hydraulic properties of gas hydrate-bearing strata are depicted in Tables 1 and 2. The liquid is injected under a constant rate. Results of computer simulation. Cundall and Strack [24] proposed an explicit calculation cycle (ECC). The calculation of motion (Newton's second law) was solved ahead of the constitutive equation (stress-strain relation, including nonlinear behavior of the rock). This approach imitated real physics because the velocity of a disturbing wave is always limited in solid and liquid. Such a tactic has provided success in preventing numerical instability and reproducing a realistic path during simulation of nonlinear behavior of the rock mass. Furthermore, ECC approach provided accounting for the path of loading what is a fundamental feature of the irreversible processes [25].
ECC delivered essentially new results concerning the irreversible behavior of the fractured rock mass. These results are summarized in the following.
Spatial asymmetry and synchronism violation. Let us recall, there were both physical and geometrical symmetries in the model relatively axis X and Z, and the liquid was injected at a constant rate. Therefore, expansion of the fractured body was expected synchronously and symmetrically in space but it was not so. Fig. 2 depicts a set of consecutive states of the fractured bodies around the coordinate origin where fluid was injected. Every state is marked with a number of cycles beginning from the start of the injection. Ten cycles correspond approximately one minute in situ. Generally, a symmetric development of the fractured body is evident but there is a tangible deviation from the symmetry. For example, set 1 residing at the positive half of the X-axis on the 130 th cycle of injection encompassed more fractured zones than antagonistic set 2 located at the negative part of this axis. The same situation concerns the sets 3 and 4 on the 140 th cycle: set 3, which is on the positive end of the Z-axis has two fractured zones against three The rate of the spatial asymmetry is not big and keeps in the range of 10% but it is much more than the error of calculation that was less than 1·10 -5 .
Spatial deviation, discontinuity, recurrence. Fig. 2 demonstrates that the fractured body evolution is discontinuing and discrete in space. There are blanks in the 3D fractured body where fractures are absent at the initial stage of HF process and patches of the virgin intact rock remain a long time. Geophysical monitoring of the microseismicity following HF process proved that the discrete pattern of hydraulic fracture development is a usual and natural phenomenon [26] (Fig. 2, bottom fragment). Eventually, HF process returned to the blanked spots and disintegrated them emphasizing the recurrence of HF. It does mean the gas hydrate strata can disintegrate by turn, with separate portions in space and in time.
As a contrast, the pore pressure diffused and expanded continuously what Fig. 3 demonstrates. The shape of the pore pressure cloud tends to correspond to the permeability anisotropy because the dimensions of the cloud in the horizontal plane intersecting the point of injection relate as 1:3 that is approximately reversely proportional to the ratio of the permeability along Z-and X-axis. Dissipative structures development. Specialists in microseismicity use a diagram in coordinate "distance from the injection point" -"time" to characterize the process of microseismicity diffusion (Fig. 4). Fragment (a) demonstrates this diagram with a wireframe surface, whereas fragment (c) shows the distribution of seismic moments for a section along the horizontal plane going through the point of injection. Apparently, HF process went spontaneously to the self-organization state. Despite the uniform distribution of mechanic and hydraulic parameters of the gas hydrate-bearing strata, the microseismic intensity evolved to the complex structure in time and space.
Periodic maximums of microseismicity (indicated with white arrows in fragment (a) interchanged with a less intensive manifestation of ground/liquid pressure that perfectly coincides with the results of geophysical monitoring [27] (fragment (b) in Fig. 4 -black arrows). Furthermore, both the parabolic envelope and the boundary of the back front match on the computer (a) and experimental (b) diagrams. Such a self-organization is natural and physically substantiated because the gas hydrate strata is a typical open thermodynamic system that transforms energy of the liquid pressure dissipating it in a form of microseismic energy, the surface energy of the hydraulically stimulated fractures, and finally heat. According to [28] such a system transits to the self-organization state if the transformed energy flow is sufficiently big.
Thermodynamics of irreversible processes has evolved over the past century, and several Nobel Prizes have marked its theoretical basis. It is important from the position of gas hydrate recovery that the rock mass is susceptible to form the dissipative structures during HS deployment even if it is geologically homogeneous and has a uniform distribution of mechanic properties. The law of the minimum entropy production [28]  The study showed [23] that there are natural sources of feedback that enhance the selfoscillatory process of periodic amplification and decay of HF, which is clearly shown in Fig. 4. There is a wide range of the dissipative structure patterns in a rock mass [29]: unusual irreversible ground movement sort of rotors, torrents, sinks, and sources; sequential asymmetrical expansion of disintegrated and loosening rock mass around underground opening [30]; yielding of the frame support clutches that proceeds by turn, one after another.
Computer simulation has demonstrated that HF can produce dissipative structures, which cause irregular fracturing of the gas hydrate-bearing strata. This reduces the positive effect of gas production stimulation. This disadvantage can increase because of the natural variation of the rock mass strength. Therefore, dissipative structures must be controlled and, first of all, slowed down or eliminated [23]. This is a subject for future research work.

Conclusions
The FLAC3D model has been used to simulate the deployment of HF in gas hydrate reservoirs at a depth of 800 m. It was found that a fractured zone around the point of liquid injection develops asymmetrically in space. This zone increases in parts and fragments one after another in different directions. The rate of spatial asymmetry is not big and keeps in the range of 10% but it is much more than the error of calculation that was less than 1·10 -5 . The fractured body evolution is discontinuing and discrete in space. There are blanks in the three-dimensional fractured body where fractures are absent at the initial stage of HF process. The patches of the virgin intact rock can remain a long time. Eventually, HF process returns to the blanked spots and disintegrates them emphasizing the recurrence of HF. At the last stage of injection, most of the blind spots were processed with HF but some intact patches could persist.
As a contrast, the pore pressure diffused and expanded continuously shaping a cloud, dimensions of which were reversely proportional to permeability in the orthotropic rock mass.
Gas hydrate-bearing strata are an open thermodynamic system, which transforms and dissipates the energy of injected liquid and ground pressure. Dissipative structures follow the HF process if the energy flow is sufficient. The dissipative structures are exposed on the twodimensional diagram in coordinates 'distance from the injection point' -'time' as periodic maximums of microseismicity interchanging with a less intensive manifestation of ground pressure. The inclination of the rock mass to generate the dissipative structures is governed by the minimum entropy production law. The dissipative structure development is triggered by small thermodynamic fluctuations, for example, a variation of ground pressure, instability of the temperature, or deviation of rock mass strength relatively average level.
The dissipative structures reduce the positive effect of gas production stimulation. This disadvantage can further increase because of the natural variation of the rock mass strength. Therefore, the dissipative structures should be controlled and, first of all, put a brake on or inhibited.