Deep fracture plane reconstruction method based on microseismic signals

Abstract

The present application belongs to the technical field of deep rock mass engineering, and discloses a deep fracture surface reconstruction method based on microseismic signals, which comprises the following steps: firstly, based on two-step DBSCAN algorithm, denoising and classification are performed on microseismic events to obtain more than one fracture surface point set; then, the normalized moment magnitude of each microseismic event is taken as the weight of the microseismic event, and the positions of the microseismic events in each fracture surface point set are fitted by using the weighted least squares method; finally, based on the projection method, the fracture surface boundary obtained by fitting is determined, and the deep fracture surface is reconstructed. The present application is a new method for reconstructing the fracture surface by using the spatial distribution of the positions and moment magnitudes of microseismic events, which can not only improve the accuracy and effectiveness of the fracture surface reconstruction, but also significantly reduce the calculation cost.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of deep rock mass engineering, and relates to a deep fracture surface reconstruction method based on microseismic signals in the hydraulic fracturing and microseismic detection in deep resource mining. BACKGROUND

[0002] For deep unconventional shale gas mining, hydraulic fracturing technology greatly improves the efficiency and feasibility of unconventional reservoir mining. Due to the deep burial depth of deep reservoirs, conventional detection techniques cannot effectively obtain fracture information after fracturing. Microseismic monitoring is applied to the monitoring and analysis of the hydraulic fracturing process because it can obtain wave signals in thousands of meters deep. The microseismic monitoring method is based on the formation of shear fractures and tensile fractures by injecting fracturing fluid into the reservoir. When the formation is broken, seismic waves are generated. These seismic waves can be captured and recorded by geophones. By analyzing the seismic waves, the source location and intensity of the microseismic event can be obtained. Therefore, the extension range and extension direction of the fracture network can be determined by using the microseismic event, and the fracture morphology can be reconstructed, which can be used to evaluate the fracturing performance and optimize the design and decision-making. However, in actual monitoring, sensors are usually located on the ground and are subject to interference noise events caused by many factors. It is still a great challenge to reconstruct the deep fracture surface.

[0003] In deep shale reservoirs, most microseismic events are recorded under low signal-to-noise ratio (SNR) conditions, resulting in a large number of low-magnitude noise events in the interpretation results. According to the generation mechanism of noise events, these noise events can be divided into two categories: the first category is caused by environmental noise, which is usually random and unrelated to fracturing activities; the second category is caused by the slip of existing fractures far from the hydraulic fracture, which may reflect the stress disturbance in the reservoir. Although the first category of environmental noise can be eliminated by the signal-to-noise ratio limit, it is difficult to distinguish the second category of non-environmental noise from effective events. Traditional methods analyze all microseismic events directly to obtain the fracture structure, which leads to an overestimation of the number of fractures and misleading judgments of the fracture geometry.

[0004] In recent years, the density-based spatial clustering (DBSCAN) algorithm has been frequently proposed for removing noise events, which is also a powerful denoising technique. However, there are still two main problems in using DBSCAN combined with other methods to construct a microseismic event-based fracture configuration: (1) how to ensure that the fracture configuration is strongly related to the characteristics of the microseismic event; (2) in the analysis of microseismic events, although the parameters of the DBSCAN algorithm have a great influence on the clustering results, it is still difficult to quickly determine them.

[0005] In summary, how to realize fast and effective reconstruction of the fracture surface in microseismic detection based on the DBSCAN algorithm is still a technical problem to be solved. SUMMARY

[0006] The purpose of the present application aims at the technical problems existing in the prior art, and provides a deep fracture surface reconstruction method based on microseismic signals, which can effectively remove noise and realize effective reconstruction of the deep fracture surface.

[0007] Hydraulic fracturing greatly improves the permeability of the reservoir by creating fractures around the wellbore or establishing new hydraulic fractures connected to natural fractures, usually extending tens to hundreds of meters around the wellbore. Induced seismic activity is usually concentrated near the wellbore, but pressure disturbances generated during the fracturing process can propagate through highly permeable channels to more distant areas, possibly inducing earthquakes hundreds or even thousands of meters away from the wellbore. This increase in permeability is reflected not only in the spatial distribution of fractures but also in the spatial pattern of induced seismic activity. Therefore, locating the fracture surface by the spatial distribution of microseismic events is an effective method to determine the fracture geometry. This method is crucial for evaluating the fracture extension during hydraulic fracturing and managing and mitigating potential geological disasters.

[0008] The present application adopts normalized moment magnitude as the weight factor for constructing the fracture surface, and the moment magnitude of each microseismic event plays a key role in the final fracture surface geometry. Moment magnitude is a standard for measuring the energy released by microseismic events. Magnitude M w <3> earthquakes, usually not felt on the surface, are called "microseismic events" or "microseismic activity". Most of such events are in the range of -3.0 w ≤0. The higher the moment magnitude, the larger the scale of rock damage, and the greater the possibility of inducing major fractures. Previous studies have shown that using moment magnitude as a key parameter to determine the location of the fracture surface is a reasonable approach. It is crucial that, due to the complex mechanism of hydraulic fracturing-induced earthquakes, microseismic positioning alone cannot capture the multi-scale, complex geometry of the final fracture surface.

[0009] In view of this, the present application uses the spatial distribution of the positions and moment magnitudes of microseismic events to derive the geometry of the main fracture structure.

[0010] The present application provides a deep fracture surface reconstruction method based on microseismic signals, comprising the following steps:

[0011] S1, based on a two-step DBSCAN algorithm, denoising and classifying microseismic events to obtain more than one fracture surface point set;

[0012] S2, taking the normalized moment magnitude of each microseismic event as the weight of the microseismic event, using weighted least squares to fit the positions of the microseismic events in each fracture surface point set in step S1;

[0013] S3, based on the projection method, determining the fracture surface boundary obtained by fitting, generating the reconstructed deep fracture surface; for any fracture surface point set obtained by classification in step S1, the following steps are performed:

[0014] S31, projecting all microseismic event positions in the fracture surface point set to the coordinate plane with the largest projection area to obtain a projected point set;

[0015] S32, identifying the boundary line of the projected point set;

[0016] S33, projecting the boundary line of the projected point set to the fracture surface obtained by fitting in step S2 to obtain the fitted fracture surface boundary, generating the reconstructed deep fracture surface.

[0017] In the above step S1, the DBSCAN algorithm is used twice for noise identification and removal, and the microseismic events are classified. Based on the first step DBSCAN algorithm, the noise points and isolated points far from the microseismic events are filtered out, and the remaining effective microseismic events are more concentrated and uniformly distributed in space; then, based on the second step DBSCAN algorithm, the denoised microseismic events are classified, and the noise points in the low-density area can be further identified and removed, improving the classification accuracy, especially for microseismic events with overlapping in space or time and complex structure.

[0018] The neighborhood search radius epsilion (Eps, ) and the minimum point threshold (MinPts) are two key hyperparameters. A larger can cluster events that are far apart in space but still belong to the same geological structure or source area, thereby preventing the deletion of potential valid signals; a smaller can more accurately identify the spatial relationship between microseismic events; therefore, the present application first uses a larger to identify and remove noise that exists alone in microseismic events, and to cluster microseismic events to effectively filter out noise points and isolated points far from the main event cluster; then a smaller is used to identify and remove noise in the fracture boundary, further cluster microseismic events and obtain more than one fracture surface point set.

[0019] The present application further determines the neighborhood search radius and the minimum point threshold of the two-step DBSCAN algorithm by a parameter optimization method based on the contour coefficient.

[0020] Specifically, starting from the minimum value of the set MinPts value, the following operations are performed iteratively; select the and that maximize the contour coefficient S:

[0021] For the current MinPts value, let k = MinPts-1, calculate the distance from each point to its k nearest neighbors, and sort the resulting distances from small to large, get a k-distance curve that reflects the trend of spatial density change, that is, the k-distance distribution diagram; Then take the slope mutation point (i.e. the point on the curve where the slope changes significantly, which can be determined by difference quotient method or derivative method) in the k-distance curve as the inflection point, and the distance corresponding to the inflection point is the corresponding ;

[0022] According to the determined and the current MinPts value, cluster all microseismic events, and calculate the silhouette coefficient of each cluster category using the following formula :

[0023] ;

[0024] In the formula, is-1~1, , represents the average distance from point j to other points in the same cluster; , represents the average distance from point j to all points in the nearest cluster, and N is the number of all microseismic events.

[0025] The purpose of the above step S2 is to fit the positions of microseismic events in each fracture surface point set to obtain the best fitting surface.

[0026] According to the coordinates of each microseismic event in each fracture surface point set determined in step S1, the positions of microseismic events on each fracture surface are fitted by weighted least squares method, and the fitting surface control equation is as follows:

[0027] ;

[0028] In the formula, x and y represent the coordinates of microseismic events in XOY, and z predicted represents the fitting value of Z coordinate, , , represents the fitting parameter.

[0029] Further, the following fitting error is defined as the objective function J:

[0030] ;

[0031] In the formula, n is the number of microseismic events in the fracture surface, M wl,i is the weight of the i th microseismic event, z predicted,i represents the Z coordinate fitting value of the i th microseismic event, z actual,iZ coordinate actual value of the i th microseismic event. By minimizing the objective function, the best fitting fracture surface is obtained. Specifically, different values of a, b and c are given to obtain different surface predicted values z predicted When the fitting error between the actual value and the predicted value of the surface is minimized, that is, when the J value is minimized, the best fitting fracture surface is obtained.

[0032] The normalized moment magnitude is used As a weight assigned to the i th microseismic event,

[0033] ;

[0034] In the formula, The i th moment magnitude, , The seismic moment of the i th microseismic event, ), , 、 And The thermal modulus, fault slip area and average fault slip displacement of the i th microseismic event corresponding to the fault medium are represented by The range value is-3~3.

[0035] The purpose of the above step S3 is to generate the boundary of each fracture surface.

[0036] In the above step S31, the fracture surface point set is projected on all microseismic event position points to obtain the maximum projection plane, which is realized by setting the maximum projection plane normal axis coordinate in the three-dimensional coordinates of all microseismic event positions to zero.

[0037] In the above step S32, the 2D Alpha-shape algorithm is used to identify the boundary of the projection point set.

[0038] In the above step S33, the boundary line coordinates of the projection point set are substituted into the fitting surface control equation to obtain the final boundary of the fracture surface, and then the reconstructed deep fracture surface is generated by combining the fracture surface fitting surface.

[0039] Compared with the prior art, the deep fracture surface reconstruction method based on microseismic signals provided by the present application has the following beneficial effects:

[0040] (1) The present application first removes noise and classifies microseismic events based on a two-step DBSCAN algorithm; then fits each fracture surface obtained by clustering; finally, the boundary line is identified by two-step projection combined with 2D Alpha-shape algorithm, and the reconstructed deep fracture surface is obtained by projecting the boundary line back to the fitted fracture surface; not only can improve the accuracy and effectiveness of fracture surface reconstruction, but also can significantly reduce the calculation cost;

[0041] (2) The application can obtain an accurate crack surface area through application of a two-step DBSCAN algorithm;

[0042] (3) The application introduces a magnitude as a weight factor in the process of fitting the crack surface, so that the contribution of microseismic events of different magnitudes in the crack surface reconstruction can reflect the physical meaning, thereby improving the accuracy and robustness of fitting;

[0043] (4) The application can effectively optimize the fitting accuracy of the crack surface geometry and data points through minimization of the objective function in the fitting process;

[0044] (5) The application solves the deficiency of the traditional convex hull method in processing irregular crack geometry by adaptively constructing the boundary of the data point set through the 2D Alpha-shape algorithm, thereby improving the accuracy and adaptability of the boundary extraction. BRIEF DESCRIPTION OF DRAWINGS

[0045] Figure 1 The flowchart of the deep crack surface reconstruction method based on microseismic signals provided for the embodiments of the application;

[0046] Figure 2 The flowchart for determining the crack surface boundary;

[0047] Figure 3 The k-distance distribution diagram when the MinPts value is 44;

[0048] Figure 4 The process diagram for determining the crack surface boundary based on the 2D Alpha-shape algorithm; wherein (a) is to generate the maximum projection plane, (b) is to generate the boundary of the identified projection point set, and (c) is to generate the crack surface boundary;

[0049] Figure 5 The spatial distribution diagram of the microseismic events recorded during hydraulic fracturing (scaled according to the moment magnitude): (a) is the map view of the microseismic events, (b) is the east-west direction depth view, and (c) is the north-south direction depth view;

[0050] Figure 6 The moment magnitude distribution diagram of the microseismic events, wherein (a) is the microseismic event distribution before and after denoising by using the first-step DBSCAN; (b) is the microseismic event distribution diagram corresponding to the three crack surfaces obtained by classification;

[0051] Figure 7 The spatial distribution diagram of the microseismic events corresponding to the three reconstructed crack surfaces;

[0052] ​​Figure 8 Figure 16 is a plot of the Euclidean distance between the three fracture planes and the corresponding microseismic events;

[0053] Figure 9 Figure 17 is a plot of the strike and dip of the focal mechanism solution of the microseismic events and fracture data; (a) and (d) correspond to the strike and dip of fracture plane 1; (b) and (e) correspond to the strike and dip of fracture plane 2; (c) and (f) correspond to the strike and dip of fracture plane 3. DETAILED DESCRIPTION

[0054] The technical solutions of the embodiments of the present application will be described clearly and completely below with reference to the drawings. Obviously, the described embodiments are only some of the embodiments of the present application, but not all the embodiments. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative work are within the present application.

[0055] Embodiment 1

[0056] The present embodiment provides a deep fracture plane reconstruction method based on microseismic signals, as shown in Figure 1, which comprises the following steps: Figures 1-2

[0057] S1, based on the two-step DBSCAN algorithm, denoising and classifying the microseismic events to obtain more than one fracture plane point set.

[0058] In this step, the microseismic events are formatted into strings, and then the DBSCAN algorithm (see Ester, M., Kriegel, H.P., Sander, J., Xiaowei, X., 1996. A density-based algorithm for discovering clusters in large spatial databases with noise. KDD-96 Proceedings. Second International Conference on Knowledge Discovery and Data Mining, 226-231) is used twice to identify and remove noise and cluster the microseismic events. Based on the first-step DBSCAN algorithm, the microseismic events are denoised to filter out noise points and isolated points far from the microseismic events, and the remaining effective microseismic events are more concentrated and uniformly distributed in space; then, based on the second-step DBSCAN algorithm, the denoised microseismic events are classified, and noise points in low-density areas can be further identified and removed to improve clustering accuracy, especially for microseismic events with overlapping in space or time and complex structure. ​

[0059] Neighborhood search radius epsilon(Eps, The minimum number of points (MinPts) and the minimum point threshold are two key hyperparameters. Larger values... It can cluster events that are spatially far apart but still belong to the same geological structure or seismic source region, thus preventing the deletion of potentially valid signals; smaller This allows for more precise identification of spatial relationships between microseismic events; therefore, this invention first employs larger [seismic events]. Identify and remove isolated noise in microseismic events, and cluster the microseismic events to effectively filter out noise points and isolated points far from the main event clusters; then use a smaller... Noise in the crack boundaries is identified and removed, and microseismic events are further clustered to obtain more than one crack surface point set.

[0060] This embodiment determines the neighborhood search radius of the two-step DBSCAN algorithm using a parameter optimization method based on the silhouette coefficient. and minimum number of points threshold .

[0061] Specifically, starting from setting the minimum value of MinPts, perform an iterative search as follows: select the option that maximizes the profile coefficient S. and :

[0062] For the current MinPts value, let k = MinPts-1, calculate the distance from each point to its k nearest neighbors, and sort the obtained distances in ascending order to obtain a k-distance curve that reflects the trend of spatial density change, i.e., the k-distance distribution map. Figure 3 The k-distance distribution is given when MinPts = 44; then, the point where the slope of the k-distance curve changes abruptly is taken as the inflection point (e.g., Figure 3 (3436 in the text), the distance corresponding to the inflection point is the distance corresponding to the current MinPts value. ;

[0063] Based on the determined Based on the current MinPts value, all microseismic events are clustered, and the silhouette coefficient of each cluster category is calculated using the following formula. :

[0064] ;

[0065] In the formula, The range is -1 to 1. This represents the average distance from point j to other points in the same cluster; Let J represent the average distance from point j to all points in the nearest cluster, and N be the total number of microseismic events.

[0066] The first step DBSCAN algorithm and the second step DBSCAN algorithm are respectively optimized by the above-mentioned parameter optimization method based on the profile coefficient. And The value is optimized.

[0067] S2, taking the normalized moment magnitude of each microseismic event as the weight of the microseismic event, fitting the position of the microseismic event in each crack surface point set in step S1 by using the weighted least square method.

[0068] According to step S1, the coordinates of each microseismic event in each crack surface point set are determined, and the position of the microseismic event on each crack surface is fitted by using the weighted least square method. The fitting surface control equation is as follows:

[0069] ;

[0070] Wherein, x, y represent the coordinates of the microseismic event in XOY, z predicted represents the Z coordinate fitting value, a, b, c represent the fitting parameters.

[0071] The following fitting error is defined as the objective function J:

[0072] ;

[0073] In the formula, n is the number of microseismic events in the crack surface, M wl,i is the weight of the i-th microseismic event, z predicted,i represents the Z coordinate fitting value of the i-th microseismic event, and z actual,i represents the actual value of the Z coordinate of the i-th microseismic event.

[0074] By introducing an optimization algorithm, the plane parameters are gradually adjusted, so that the value of the objective function gradually decreases. The objective function provides a clear evaluation standard for the optimization problem, which is used to measure the pros and cons of the model parameters. By optimizing the objective function, a set of parameters can be found to minimize the value of the objective function, so that the fitting effect or performance is optimal.

[0075] Specifically, different values are given to a, b, and c to obtain different surface prediction values z predicted When the fitting error between the actual value and the predicted value of the surface is minimized, that is, when the value of J is minimized, the crack surface with the best fitting effect is obtained.

[0076] The moment magnitude can accurately represent the energy release of the microseismic event, and provides important geophysical information in microseismic monitoring of deep shale gas development and mining. In this embodiment, the normalized moment magnitude is used as the weight of the i-th microseismic event,

[0077] ;

[0078] wherein, is the i-th moment magnitude, , is the seismic moment of the i-th microseismic event, ), , , and respectively represent the thermal modulus, the fault slip area and the average fault slip displacement of the i-th microseismic event corresponding to the fault medium; The range value is -3-3.

[0079] S3, based on the projection method, determines the fracture surface boundary obtained by fitting, and generates a reconstructed deep fracture surface.

[0080] The embodiment proposes a method of applying 2D Alpha-shape algorithm and two-step projection to obtain the fracture surface boundary. For any fracture surface point set classified in step S1, as shown in Figure 4 , the following steps are performed:

[0081] S31, project all microseismic event positions in the fracture surface point set to the coordinate plane with the largest projection area to obtain a projected point set.

[0082] Project all microseismic event position points in the fracture surface point set to obtain the maximum projection plane. Here, it is realized by setting the maximum projection plane normal axis coordinate in the three-dimensional coordinates of all microseismic event positions to zero. As shown in Figure 4 (a), when the projection area of the fracture surface in the XOY plane is the largest, only the Z coordinate of each microseismic event on the fracture surface needs to be set to 0.

[0083] S32, identify the boundary line of the projected point set.

[0084] This step uses the 2D Alpha-shape algorithm to identify the boundary of the projected point set. Alpha-shape is a generalized convex hull algorithm for a given point set, which can more accurately represent the shape of the point set. In this embodiment, according to the set radius a, the Alpha-shape algorithm is used to identify the microseismic event point boundary on the two-dimensional plane after rolling the edge, and the specific operation can be referred to in H. Edelsbrunner, D. Kirkpatrick and R. Seidel, "On the shape of a set of points in the plane," in IEEE Transactions on Information Theory, vol. 29, no. 4, pp. 551-559, July 1983, doi: 10.1109 / TIT.1983.1056714.

[0085] It is worth noting that although this method provides high accuracy for polygonal fracture configuration, careful consideration of the projection direction is essential. This step can ensure accurate representation of the fracture boundary, thereby minimizing distortion in the final result.

[0086] S33, project the boundary line of the projected point set onto the fracture surface fitted in step S2 to obtain the fitted fracture surface boundary and generate the reconstructed deep fracture surface.

[0087] This step substitutes the boundary line coordinates of the projected point set into the fitted surface control equation to obtain the final boundary of the fracture surface, and then generates the reconstructed deep fracture surface in combination with the fracture surface fitting surface.

[0088] Application example

[0089] This application example compares the fracture network formed after multi-stage hydraulic fracturing of the A well section in Luzhou Oilfield with the fracture surface obtained by the deep fracture surface reconstruction method based on microseismic signals provided by the present application. The foregoing will be described in detail below.

[0090] This well section is located in Luzhou Oilfield, with a depth of 3100m to 3400m and a horizontal length of about 1550m. Multi-stage hydraulic fracturing is performed to form a fracture network and maximize reservoir permeability by maintaining high flow and continuous sand injection. During the stimulation process, a radial configuration of surface geophones is used to record microseismic events. The configuration consists of 1312 geophone stations distributed on 12 measurement lines.

[0091] As Figure 5The image shows the spatial distribution of microseismic events recorded during hydraulic fracturing (scaled by moment magnitude), containing 4110 microseismic events. The moment magnitude can be seen from the image. The range was between -2.72 and 2.19, with numerous microseismic events observed occurring in clusters around well section A, while the density of microseismic events was lower further away from well section A. Furthermore, a significant portion of the microseismic events within the cluster exhibited moment magnitude. The trend of values ​​greater than 0 suggests that the fracture surface may be roughly perpendicular to the Earth's surface.

[0092] Following the deep fracture surface reconstruction method based on microseismic signals presented earlier, the hyperparameters of the two-step DBSCAN algorithm are first determined using the parameter optimization method based on profile coefficients presented earlier. and As shown in Table 1.

[0093] Then, the hyperparameters shown in Table 1 are used. and The microseismic events are denoised and classified according to the two-step DBSCAN algorithm given in step S1. In this application example, the minimum value of MinPts is 2 in both steps of the DBSCAN algorithm.

[0094] Table 1 DBSCAN Parameter Selection Results

[0095] .

[0096] Moment magnitude of microseismic events Distribution as Figure 6 As shown, from Figure 6 As can be seen from this, the moment magnitude of most microseismic events... All values ​​are less than -1, with a significant portion falling between -3 and -2. After two steps of the DBSCAN algorithm for denoising and classification, the remaining microseismic events still exhibit similar moment magnitude distributions. Figure 6 (a) shows that a significant portion of the microseismic events (more than 1000 events) in the denoised dataset were filtered out, and the moment magnitudes of these filtered microseismic events... Between -3 and -2. Figure 6 (b) Shows the moment magnitudes of the microseismic events corresponding to the three fracture surfaces. The distribution shows that the moment magnitudes of most microseismic events corresponding to the three fracture surfaces are below -1, and their distribution pattern is similar to... Figure 6 The event distribution after the first step of DDSCAN denoising is similar to that shown in (a).

[0097] Based on the identified moment magnitude Using the normalized moment magnitude given above The formula for calculating the normalized moment magnitude of each microseismic event .

[0098] Then the best fitting fracture surface is obtained by the least square method combined with the fitting control equation and the objective function in step S2.

[0099] Then the reconstructed deep fracture surface is generated according to the step S3 given above combined with the 2D Alpha-shape algorithm and two-step projection, as shown in Figure 7 . Figure 7 Three polygonal models representing the fracture geometry are given.

[0100] Figure 8 The Euclidean distance between microseismic events on the three fracture surfaces is compared in the figure. Figure 8 From the box plot and trend line in , it can be seen that most of the microseismic events on the three fracture surfaces are concentrated within 20 meters; at the same time, the median distances of fracture surfaces 1, 2 and 3 are 12.17, 15.42 and 20.35 meters, respectively. 75% of the microseismic events on the three fracture surfaces are located within the interquartile range of 21.38 meters, 23.22 meters and 37.59 meters, respectively, indicating a strong correlation between the fracture surface geometry and the microseismic events. Fracture surfaces 1 and 2 have 19 and 1 outliers, respectively, while fracture surface 3 has 77 outliers. The number of microseismic events on fracture surfaces 1, 2 and 3 is 632, 186 and 1403, respectively, and the proportion of outliers is 3.01%, 0.54% and 5.49%, respectively.

[0101] The focal mechanisms of the three groups of microseismic events around fracture surface 1, fracture surface 2 and fracture surface 3 in Figure 7 are calculated by moment tensor inversion (P-wave polarity inversion), and the specific calculation process is described in (Hardebeck, J.L., Shearer, P.M., 2002. A new method for determining first-motion focal mechanisms. BULLETIN OF THE SEISMOLOGICAL SOCIETY OF AMERICA 92(6), 2264-2276).

[0102] The strike and dip of the fracture surface structure are calculated by selecting 66, 47 and 48 microseismic events with relatively small positioning errors in fracture surface 1, fracture surface 2 and fracture surface 3, respectively, and the calculation formula is as follows:

[0103] ;

[0104] ;

[0105] wherein, is the dip or strike of the final calculated fracture plane; is the dip or strike of the fracture plane of the i-th microseismic event solved based on the focal mechanism of the i-th microseismic event, is the Euclidean distance of the i-th microseismic event to the fracture plane fitted in step S2; is used as a hyperparameter to control the impact of distance on weight, in the present application is set to 1.

[0106] The strike and dip of the fracture plane structure are shown in Table 2. The microseismic data is shown in Figure 9 to compare the results of the crack reconstruction method proposed in the present application and the moment tensor inversion method on the spatial trend of crack propagation. From Figure 9 it can be seen that the focal mechanism results of the clusters around fracture plane 1 and fracture plane 2 are in good agreement. Although the strike data of fracture plane 3 is consistent, the dip data shows greater variability and complexity. The strike of fracture plane 1 is 38.41°, and the dip is 64.13°, which is in good agreement with the dominant direction observed in microseismic. The strike of fracture plane 2 is 50.45°, and the dip is 68.07°, which is in good agreement with the dominant direction of microseismic. The strike of fracture plane 3 is 19.26°, and the dip is 76.75°, which is in good agreement with the strike result of microseismic. In summary, the crack configuration direction obtained by focal mechanism analysis is consistent with the main direction of the calculated crack plane.

[0107] Table 2. Fracture structure: event count, strike and dip

[0108] .

[0109] Therefore, in the present application, the DBSCAN method is used in two steps to denoise and classify, which avoids over-sensitivity to event data in a single processing step, thereby reducing the risk of over-fitting (i.e. over-reliance on details) or under-fitting (i.e. failure to identify key structures). After initial denoising, the amount of data is reduced, thereby providing a more streamlined data set for subsequent classification. This not only improves computational efficiency, but also reduces the risk of false classification caused by noise points in the classification process. The present application also uses an optimization method that combines the use of k-distance maps with a contour coefficient-based parameter optimization method to select the hyperparameters of DBSCAN. Compared with the traditional grid search method, this method provides a clear indicator for evaluating the quality of clustering and reduces the computational cost. In addition, for data sets with an undetermined number of clusters, this method can improve the accuracy of clustering.

[0110] A key difference between the present invention and previous studies is the use of a weighted least trimmed squares optimization algorithm based on moment magnitudes to determine the spatial extension of the fracture structure. In contrast, previous studies often relied on the random sample consensus (RANSAC) method or simpler optimization techniques based on this method. The least trimmed squares minimization approach is advantageous for more accurate fitting of significant events, especially when these events exhibit different energy levels. The addition of moment magnitude weights allows the fitting process to better capture the impact of large magnitude events on the model. Furthermore, compared to RANSAC, the least squares method generally has higher computational efficiency, especially when dealing with large datasets, as it does not require a large number of random samples and iterations. Most importantly, the performance of RANSAC is highly dependent on the selection of two key parameters, the maximum number of iterations and the tolerance range, given that DBSCAN itself requires parameter adjustment, combining the two algorithms and selecting the best parameters becomes increasingly complex.

[0111] Moreover, in previous studies, the fracture structure was often assumed to have various shapes, such as rectangular or tree-like fracture models, but in reality, the shape of the fracture is generally irregular. Traditional boundary extraction methods such as convex hull method tend to oversimplify the fracture surface topography, while the two-step projection combined with 2D Alpha-shape algorithm used in the present invention has better anti-interference ability when dealing with noisy data, ensuring that the fitting results accurately reflect the complexity of the actual geological structure.

[0112] Those skilled in the art will realize that the embodiments described herein are for the purpose of helping the reader understand the principles of the present invention and should be understood as not limiting the scope of protection of the present invention to such specific statements and embodiments. Those skilled in the art can make various other specific modifications and combinations according to the technical inspiration disclosed in the present invention without departing from the spirit of the present invention, and these modifications and combinations are still within the scope of protection of the present invention.

Claims

1. A method for deep fracture plane reconstruction based on microseismic signals, characterized in that, The method comprises the following steps: S1, denoising and classifying microseismic events to obtain more than one fracture surface point set based on a two-step DBSCAN algorithm; denoising microseismic events based on a first-step DBSCAN algorithm; classifying the denoised microseismic events based on a second-step DBSCAN algorithm; the neighborhood search radius in the first-step DBSCAN algorithm is greater than the neighborhood search radius in the second-step DBSCAN algorithm; The neighborhood search radius and the minimum point threshold MinPts of the two-step DBSCAN algorithm are determined by a parameter optimization method based on contour coefficients respectively: The iterative search is performed as follows starting with a minimum value of MinPts; select the cluster with the largest silhouette coefficient S and : For the current MinPts value, let k = MinPts-1, calculate the distance from each point to its k nearest neighbors, and sort the resulting distances in ascending order to obtain a k-distance curve; then take the slope mutation point in the k-distance curve as the inflection point, and the distance corresponding to the inflection point is the corresponding ; According to the determined and the current MinPts value, all microseismic events are clustered, and the silhouette coefficient of each cluster category is calculated by the following formula : ； wherein is -1 ~ 1, denotes the average distance of point j to other points of the same cluster; denotes the average distance of point j to the closest cluster of all points, N is the number of all microseismic events; S2, taking the normalized magnitude of each microseismic event as the weight of the microseismic event, and fitting the position of the microseismic event in each fracture surface point set in step S1 by using a weighted least square method; S3, determining the boundary of the fitted fracture surface based on a projection method to generate a reconstructed deep fracture surface; for any fracture surface point set classified in step S1, the following steps are performed: S31, projecting all microseismic event positions in the fracture surface point set to the coordinate plane with the largest projection area to obtain a projected point set; S32, identifying the boundary line of the projected point set; S33, projecting the boundary line of the projected point set to the fracture surface fitted in step S2 to obtain the boundary of the fitted fracture surface and generate a reconstructed deep fracture surface.

2. The method of claim 1, wherein, In step S2, the coordinates of each microseismic event in each fracture surface point set are determined according to step S1, and the positions of the microseismic events on each fracture surface are fitted by using a weighted least square method, and the fitting surface control equation is as follows: ； where x, y represent the coordinates of the microseismic event in XOY, z predicted represents the Z coordinate fitting value, , , represents the fitting parameter.

3. The method of claim 2, wherein, The following fitting error is defined as the objective function J: ； where n is the number of microseismic events in the fracture plane, M wl,i is the weight of the i-th microseismic event, z predicted,i represents the Z coordinate fitted value of the i-th microseismic event, z actual,i represents the Z coordinate actual value of the i-th microseismic event.

4. The method of claim 3, wherein, Using normalized moment magnitude as a weight assigned to the i-th microseismic event, ； wherein is the i-th moment magnitude, , is the seismic moment of the i-th microseismic event.

5. The method of claim 1, wherein, In step S32, a 2D alpha-shape algorithm is used to identify the boundary of the projected point set.

Images