{'fraction_non_alphanumeric': 0.02219565443393284, 'fraction_numerical': 0.006653498284644973, 'mean_word_length': 4.99159140454687, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Characterization of Fault-slip Taking Place in a Fractured Fault Damage Zone

Significant efforts have been made in numerically simulating anthropogenic seismicity in fault zones. However, such studies do not adequately consider the complex geological structure and stress state within a fault damage zone. The present study aims at investigating the effect of the near-fault complex stress state on the intensity of fault-slip. To this end, a three-dimensional numerical model including a fault core and its damage zone is constructed with an anisotropic compliance tensor equivalent to the elasticity of fractures within the damage zone as well as interface elements to model fault-slip occurring along the core. The analysis result shows that the model yields a fault-slip with the maximum shear displacement as much as 0.14 m (M w = 2.55), whilst the maximum shear displacement produced in a homogeneous, isotropic model is merely 0.038 m (M w = 2.03), implying the significant influence of stress concentrations generated in the fault damage zone on the shear behaviour of the fault. This study emphasizes the importance of considering the complex mechanics of fault damage zone in estimating the intensity of induced seismicity and provides an insight into the development of methodologies to predict source parameters of induced seismicity more accurately.


Introduction
Fault-slip is an important physical phenomenon pertinent to various engineering projects that could induce significant stress change at great depths, thereby re-activating pre-existing geological structures, e.g., deep underground mining, geothermal energy development, oil and gas production, and tunnelling. Importantly, seismic waves arising from fault-slip as well as the aseismic deformation of the fault undergoing slip can inflict devastating damage to underground openings and/or surface facilities. This could pose a life-threatening risk. Furthermore, the shear displacement induced by fault-slip can change the hydraulic properties of the fault, thus exerting a direct influence on fluid and gas flow in the fault zone and thus affecting the productivity of fluid and gas and their sustainability. Hence, an accurate estimation of seismic source parameters of fault-slip, such as location, magnitude, and dynamic parameters, is crucial to accomplish the sustainable development of underground resources whilst mitigating the risk for fatal accidents.
Over the past decades, a number of studies have been undertaken with numerical simulations in an effort to evaluate the potential for fault-slip induced by anthropogenic stress changes, reproduce the dynamic and quasi-static shear movements, and estimate the severity of the damage inflicted. In such simulations, attentions have been predominantly paid to the mechanical properties of causative faults and their constitutive frictional behaviour during aseismic and seismic periods. Although the studies have contributed to a better understanding of fault-slip mechanism and its influence on the surrounding rock mass, it is still challenging to identify where an anthropogenic seismic event is triggered and estimate its source parameters. This is partially attributed to the fact that most of the numerical simulations assume that the initial stress regime varies only with depth and is homogeneous at the same depth. In fault-slip simulations, this assumption could lead to an inaccurate result. In fact, seismic monitoring in deep underground mines indicates that severe seismic activities can take place in a region away from an active mining area, despite the infinitesimal magnitude of mining-induced stress change, implying the presence of stress anomalies in the fractured rock mass, which cannot be considered in conventional numerical simulations.
In recent years, several studies have numerically simulated the heterogeneous pre-mining stress state by taking into account the difference in stiffness amongst large-scale geological structures. Based on the studies, Sainoki et al. [1] proposed a methodology to reproduce a complex stress state with metre-scale stress anomalies resulting from the difference in fracture density in a large-scale model with dimensions of several hundreds of meters in the framework of continuum analysis. The method has made it possible to simulate a complex, heterogeneous stress state in a large-scale fault damage zone (fractured rock mass) composed of millions of fractures whilst reducing the amount of PC memory required compared to the discrete element method.
In light of the recent advancement in the simulation method of the in-situ stress state in deep underground, the present study investigates the influence of the complex, heterogeneous stress state in a fault damage zone on the intensity of fault-slip taking place along the fault core by applying the method to fault-slip simulation. It is then aimed to quantify the difference in seismic source parameters between fault-slips taking place in heterogeneous and homogeneous stress states. This lays the foundation for a more accurate estimation of seismic source parameters of anthropogenic seismicity.

Procedure to construct an equivalent continuum model Examples
The presents study employs the method developed by Sainoki et al. [1] to construct an equivalent continuum model and simulate fracture-induced heterogeneous stress state within a fault damage zone. The method employs the crack tensor theory [2] for the model construction. As the first step, a discrete fracture network (DFN) composed of millions of fractures is statistically generated. Simultaneously, a continuum model composed of parallelepiped zones is generated in the framework of the finite difference method (FDM) with FLAC3D [3]. Then, for each zone in the continuum model, a full anisotropic compliance tensor is calculated based on the characteristic of the DFN located in a parallelepiped region corresponding to each zone in the FDM model. Eventually, the heterogeneity and anisotropy of the rock mass stiffness are realized with the anisotropic compliance tensor equivalent to the elasticity of fractures in the fault damage zone.
The anisotropic compliance tensor is calculated as follows according to the crack tensor theory.
where F and n are the crack tensor and the unit normal vector of fractures intersecting with a zone in the DFM model, respectively; C, M, and T represent the anisotropic compliance tensor equivalent to the elasticity of the fractures, the isotropic compliance tensor of the rock mass, and the anisotropic compliance tensor of fault damage zone, respectively; NCR is the total number of fractures intersecting with the zone; Kn and Ks are shear and normal stiffnesses of the fracture, respectively; δ is the Kronecker's delta; and D denotes the diameter of the fracture.

Fracture characteristics in a fault zone
To construct a reliable equivalent continuum model, it is of paramount importance to reproduce a DFN with realistic fracture distribution in a fault damage zone. Figure 1(a) illustrates a typical fault zone, which is generally divided into a fault core and the surrounding fault damage zone depending on the physical properties of the rock mass. The fault damage zone is further divided into inner and outer damage zones based on the characteristic of fracture density. Previous studies report that fracture density within a fault damage zone decreases with the distance from the fault core , which can be approximated with a power law as shown in figure 1(b). In the presents study, only metre-scale macrofractures are statistically generated for the DFN whilst applying the decay law, and the fractures are considered in the calculation of the crack tensor and [C] in equations (1), (2), and (3). This is based on an assumption that the influence of macro-fractures on the in-situ stress state is more significant than that of meso-and micro-fractures in the order of cm or less.
As for meso-fractures in the host rock, their characteristics are assumed to be homogeneous and are implicitly considered by applying the concept of geological strength index (GSI) , i.e., the property of the intact rock is reduced according to a GSI assumed, which gives the compliance matrix of the rock mass [M] in equation (4). Note that the compliance matrix for the rock mass does not vary in the model. Specific input parameters are shown in a later section.  Table 1 lists parameters used to statistically generate fractures in a fault damage zone. As shown, the maximum and minimum fracture lengths are set at 10 m and 2 m, respectively, because only metrescale macro-fractures are considered. Regarding the statistical variation of the fracture length, a power law function is applied with a constant of 2.5, based on previous study [6] . The statistical variation of the dip angle and dip direction of fractures is assumed to follow the Fisher distribution. Regarding fracture density in the fault damage zone that varies with the distance from the fault core, P10 is employed as the index is commonly measured in the field, thus frequently reported in previous studies. The coefficient of the power law function for the fracture density represents the maximum fracture density near the fault core, and d denotes the distance from the fault core. The value falls within a reasonable range of near-fault fracture density that varies from 10 /m or less to 100 /m [7] .    Table 2 summarizes model input parameters. The elastic modulus is determined by reducing the Young's modulus of intact granite (60 GPa) [8] while assuming a range of GSI from 60 to 70 for the host rock excluding the effect of macro-fractures in the faulted region. For granite, Poisson's ratio ranges from 0.17 to 0.3 [8,9], according to which, 0.2 is assumed. The joint normal and shear stiffnesses of 30 GPa/m and 10 GPa/m for the crack tensor model are the ones substituted to equation (3) to calculate the equivalent compliance tensor. It is to be noted that the crack tensor theory embeds the scale-dependency of joint stiffness in the equation in the form of inverse relation. Such that, the stiffnesses in the table correspond to those of a joint with a unit length. As for the interface elements representing the fault core, 1 GPa/m and 0.33 GPa/m are applied as joint normal and shear stiffnesses. Considering the difference in scale between fractures in the fault damage zone and the fault core [10], the normal and shear stiffnesses of the fault core are reduced compared to those for the crack tensor model. As for the friction angle of the fault core, 30°is assumed, which falls within a reasonable range of frictional coefficients of rock joints [11].  Figure 3 shows the analysis and boundary conditions. As shown in figure 3(a), the extremely high friction angle and cohesive strength are applied to the interfaces coloured in red with the aim of achieving mechanical convergence by preventing the entire fault core from undergoing slip. Figure 3(b) shows stresses applied on the model boundaries to simulate an initial stress condition. This method is called boundary traction method [12] and have been used in previous studies to simulate burst-prone stress conditions induced by stiffness difference between macro-scale geological structures. Sainoki et al. [1] combined the method with the crack tensor theory to reproduce the heterogeneous stress state in a fault damage zone with the purpose of elucidating the mechanism of severe seismicity taking place in a region away from active mining areas.

Model input parameters and boundary condition
As can be seen from figure 3(b), theoretically no slip takes place under the frictional strength and the initial stress regime. To initiate slip on the interface, after simulating the initial stress state, the shear stress (σxy) acting on the model lateral boundaries is increased from 5 MPa to 6.5 MPa by assuming anthropogenic activities in deep underground, such as large-scale ore extraction and underground cavern excavation for a hydraulic power plant. As the friction angle of the fault is 30°, the shear stress is large enough to cause slip along the fault under the confining stress of 10 MPa. This simple analysis procedure is intended to focus on the influence of the heterogeneous stress state on the characteristic of slip on the fault by excluding other influential factors, e.g., stress re-distribution induced by mining activity that varies from place to place, thus making the evaluation more complicated.   Figure 4 shows the stress state on cross-sections in the initial stress state. It is found from the figure that the heterogeneous stress state is successfully simulated in the fault damage zone, showing metrescale stress anomalies exceeding the maximum compressive stress applied on the model boundaries. More specifically, the red-coloured region in figure 4(a) corresponds to zones with the maximum compressive stress of approximately 20 MPa, which is twice as large as the stress applied on the model boundaries. Likewise, there are regions with abnormal minimum stress states that are significantly lower than the applied stress level. These stress anomalies are caused by the anisotropic behaviour of zones with the crack tensor model and the stiffness heterogeneity. As implied from the previous study [1], it is reasonable to conceive that these stress anomalies contribute to the occurrence of intense seismic activity in a fault damage zone. Under this initial stress state with a significant degree of heterogeneity, fault-slip simulation is performed.  Figure 5 depicts slip distribution on the fault simulated with the shear stress increased to 6.5 MPa. For comparison, the same numerical simulation was performed for a homogeneous model, to which only the elastic model with an elastic modulus of 10 GPa and Poisson's ratio of 0.2 was applied. The elastic modulus of the rock mass (40 GPa) is uniformly reduced to 10 GPa to take into account the influence of macro-fractures in the fault damage zone on the rock mass stiffness. Figure 5 clearly indicates the impact of the heterogeneous stress state on the magnitude of fault-slip. Specifically, the maximum slip displacement is approximately 0.15 m in the heterogeneous model, whilst the maximum magnitude of shear displacement is as large as 0.03 m in the homogeneous model. More interestingly, the slip distribution in the heterogeneous model is non-uniform, resulting in the generation of regions with intense slips and irregular slip propagation from the hypocentre with the largest slip magnitude. This implies the significant impact of the heterogeneity of the initial stress state on the characteristic of fault-slip. In other words, this result suggests that an accurate estimation of the initial stress state in a fault zone be achieved in order to evaluate the severity of fault-slip and its risk, especially when the fault core is surrounded by a developed damage zone composed of fractured rock masses.

Comparison of shear displacement between the anisotropic model and a homogeneous model
For further comparison, seismic moment is calculated for the two models. It was then shown that the seismic moment of the fault-slip in the heterogeneous model is 7.45 × 10 12 N·m, whilst that of the homogeneous model is 1.26 × 10 12 N·m. The difference between the two models becomes more pronounced in terms of seismic moment, indicating the considerable influence of pre-existing stress anomalies in the initial stress state on the magnitude of resultant fault-slip.

Discussion
To delve into the impact of the stress anomalies on the slip distribution and the magnitude of fault-slip, the stress state on the fault is examined in the initial stress regime. Figures 6(a) and (b) illustrate the shear and normal stresses acting on the fault in the initial stress state. As can be seen, metre-scale stress anomalies with higher shear and normal stresses than the values applied on the model boundaries are present throughout the fault plane. Note that the high stress region near the edge is related to the boundary effect and does not exert any influence on the characteristic of fault-slip simulated since the region is included in the no-slip region. Figures 6(a) and (b) indicate that although there are highly stressed regions, such regions do not coincide with areas undergoing intense slip shown in figure 5(a). This implies that the degree of stress concentration is not directly related to the intensity of fault-slip. In light of this result, fault-slip potential defined as the ratio of shear stress to shear strength determined by the friction law is computed. The result is depicted in figure 6(c). If fault-slip potential is equal or close to 1, it indicates the occurrence of slip. It is remarked that an extensive area undergoes slip in the initial stress state, albeit the stress state applied on the model boundaries that does not induce slip theoretically. Unlike the shear and normal stress distribution, the correlation between fault-slip potential and the magnitude of fault-slip is obvious, that is, the area with high fault-slip potential in the initial stress state is subjected to relatively large magnitudes of shear displacement when the shear stress is increased, as shown in figure 5(a). This result provides insight into the occurrence of intense seismic activity in a region not significantly affected by anthropogenic activity in deep underground and sheds light on the development of more advanced fault-slip simulation methodology considering stress anomalies on causative faults.

Conclusion
The present study examines the influence of the heterogeneous stress state generated in a fault damage zone on the intensity of fault-slip induced by anthropogenic activity. An equivalent continuum model representing a fault zone is constructed with anisotropic compliance tensor derived from a discrete fracture network, whilst interface elements are used to model the fault core. Initial stress analysis is then performed whilst applying stresses to the model boundaries. Thereafter, fault-slip is initiated by increasing the shear stress acting on the model boundaries. For comparison, the same simulation is performed for a homogeneous model. The result shows that there is a fivefold difference in the maximum shear displacement between the two models. In terms of seismic moment, the difference becomes more pronounced. These results indicate the significant influence of the heterogeneous stress state on the intensity of fault-slip. To verify this, fault-slip potential defined as the ratio of shear stress to shear strength is calculated for the fault core in the initial stress state, which demonstrates that high fault-slip potential regions coincide with the area undergoing intense slip. This result suggests that an accurate estimation of the in-situ stress state in a fault damage zone is indispensable in fault-slip simulation to estimate seismic source parameters accurately. Furthermore, the analysis result gives insight into the mechanism of intense seismic activity taking place in a region not significantly affected by stress change resulting from anthropogenic activity in deep underground.