A method for identifying different types of crack extension in laminated shale based on digital images

Abstract

The present invention discloses a method for identifying different types of crack extensions in laminated shale based on digital images, and relates to the field of laminated shale oil and gas exploration and development. The present invention first takes samples on site, and then performs explicit processing on sample features, uses a two-dimensional convolution model to extract lamination type features, selects laminated shale that meets the conditions, and performs speckle production on the laminated shale, and then performs mechanical experiments, and simultaneously performs digital image acquisition, uses a high-speed camera to capture synchronous images, analyzes the acquired image data and a series of experimental data obtained from mechanical test experiments, and uses a reliability-oriented initial value transfer strategy algorithm in the area of ​​interest to obtain visualized displacement data, and finally identifies the crack type. The present invention realizes the identification of laminated shale features, obtains a laminate grayscale map, selects appropriate laminated shale for digital image analysis, and obtains a digital image analysis map. Finally, the crack type is identified by calculation.

Description

Technical Field

[0001] The invention belongs to the research field of crack propagation characteristics of laminated shale, and specifically is a method for identifying different types of crack propagation of laminated shale based on digital images. Background Art

[0002] With the rapid development of the national economy, conventional oil and gas resources are difficult to meet my country's energy needs, while my country is rich in shale oil and gas resources. Shale is a porous medium with a layered structure, and my country's continental shales generally develop multiple types of laminae such as silt, tuff or dolomite, which makes the rock more anisotropic and heterogeneous. Laminated shales are rich in natural microcracks, bedding, laminations and other types of macroscopic discontinuous weak surface structures. These defects greatly affect the mechanical properties of the rock mass, resulting in complex crack morphology when the rock mass is destroyed. Different cracks require different stresses for initiation and expansion, and the paths after expansion are also very different, which has an important impact on the stability of the wellbore in actual drilling projects and the expansion of artificial cracks in fracturing projects.

[0003] At present, the research on the law of shale crack expansion mainly focuses on layered shale without considering the laminar structure, and the description of the crack expansion path remains at the macroscopic scale and qualitative description stage, which has limitations in guiding the exploration and development of laminar shale oil and gas. Based on the current situation that multiple types of laminae are widely distributed in my country's continental shale reservoirs, the present invention proposes a method for identifying the expansion of different types of cracks in laminar shale based on digital images, clarifies the dynamic characteristics of laminar shale crack expansion from a microscopic scale, and accurately classifies it according to the crack expansion trend. The classification results reflect the differences in crack expansion within the same medium and between different media of the rock mass, which is of great significance to practical engineering problems such as the prediction of the degree of wellbore collapse during drilling and the prediction of the scale of artificial crack expansion in fracturing engineering. Summary of the invention

[0004] The purpose of the present invention is to identify the crack extension conditions of different types of laminar shale. By identifying the different types of cracks in laminar shale, digital image processing is performed to obtain the full-field stress and strain analysis results. The general flow chart of the digital image analysis of the present invention is as follows: Figure 1 shown.

[0005] The technical solution of the present invention is as follows:

[0006] A method for identifying different types of crack extensions in laminated shale based on digital images, the implementation steps of which are as follows:

[0007] Step 1, laminae identification: Cut the core with the size of φ56.6mm×80mm into a cuboid rock sample of 40mm×40mm×80mm, grind and polish each surface, then scan the end face, observe the end face of the core, and test the core parameters. The laminae type identification result and grayscale image are obtained through the laminae identification algorithm, the laminae distribution characteristics of the sample are determined, and the rock samples that meet the conditions are selected for the next experiment.

[0008] Step 2, speckle production: polish and grind the surface of the rectangular rock sample, and produce speckle on the surface of the rock sample, then conduct mechanical tests and collect images.

[0009] Step 3: Image acquisition: While conducting the mechanical test experiment, a high-speed camera is used to capture synchronous images to obtain a set of images containing displacements. The acquired image data and a series of experimental data obtained from the mechanical test experiment are analyzed.

[0010] Step 4, calculate crack displacement data: In the displacement image group, select the region of interest and determine the displacement data in the grid within the region. The reliability-oriented initial value transfer strategy algorithm is used to obtain the crack displacement data.

[0011] Step 5: According to the displacement trend and displacement of the pixels on both sides of the crack data, three crack type standards are constructed, a new method for crack type discrimination is established, and different crack types are identified.

[0012] Further, the specific steps of step 1 are as follows:

[0013] Step 1.1: Scan the core end face, and the following conditions must be met: Figure 6 As shown, D represents the thickness of the shale matrix between two adjacent laminae, and d represents the thickness of a single laminae. The experimental samples are selected to ensure that: Max(d)<Max(D / 2), ∑d is less than 0.3×∑D. When the actual laminae are not parallel to the rock samples, the Monte Carlo algorithm is used to calculate the irregular laminae area, as shown in Figure 2 As shown, the following conditions need to be met: S represents the area of ​​the shale matrix between two adjacent laminae, and s represents the area of ​​a single laminae. By Max(s)<Max(S / 2), ∑s is less than 0.3×∑S.

[0014] Step 1.2: After screening and obtaining rock samples that meet the conditions, first obtain the grayscale image of the rock core, and its conversion formula is:

[0015] Gray=0.299*R+0.587*G+0.114*B (1)

[0016] Where Gray is the gray value of the core gray image; R is the red component in the core gray image; G is the green component in the core gray image; and B is the blue component in the core gray image.

[0017] Then, the grayscale image is enhanced and the laminae features are extracted using a two-dimensional convolution model. The grayscale curve is extracted at the quartile of the core image. After the average grayscale value of the grayscale curve is used to construct the average grayscale curve, the average grayscale curve is subjected to wavelet transform to remove the noise signal and obtain the average grayscale curve after wavelet transform. The grayscale difference between the average grayscale curves before and after the wavelet transform is calculated to obtain the grayscale difference curve. The average grayscale value and grayscale difference before and after the wavelet transform are used as feature values ​​for cluster analysis. The laminae type recognition results are obtained by identifying dark laminae in the first round and light laminae and dark laminae in the second round.

[0018] Step 1.3, use Python database to calculate the gray value of each row of rock image at the quartile and convert it to integer type to get the quartile, store the depth and quartile of the rock image in Python database, generate three gray curves, normalize the average gray curve and gray difference curve, set the number of clusters of K-Means algorithm, use K-Means algorithm to perform cluster analysis on the standardized eigenvalues, calculate the adjacent depth difference to represent the thickness of the laminae, and save the clustering results, depth difference and eigenvalue obtained by cluster analysis into Python database, and finally draw the difference frequency histogram according to the clustering results. The laminae type recognition result and grayscale image are obtained through the algorithm.

[0019] Further, the specific steps of step 2 are as follows:

[0020] Step 2.1, speckle processing: speckle production is performed on the surface of the laminated shale, using black speckle and a speckle density ratio of 50%, that is, the areas occupied by bright and dark pixels in any sub-region are the same.

[0021] Step 2.2, experimental image acquisition: When conducting uniaxial compression experiments on laminated shale, two white light sources are used to continuously illuminate the opposite side of the rock at an interval of 60°. The sample surface should be perpendicular to the camera optical axis, the ambient temperature should be controlled at 15℃~25℃, and the shooting environment should have no windows or observation windows.

[0022] Further, the specific steps of step 3 are as follows:

[0023] Step 3.1, Image Analysis requires using optical images (excluding images taken by scanning electron microscopy, atomic force microscopy, or X-rays), adjusting high-speed camera settings (camera resolution needs to be greater than 1080p HD / 30fps) and capturing synchronized images.

[0024] Step 3.2: When conducting mechanical experiments, the loading rate must be controlled within 0.5-1.0 MPa / s.

[0025] Step 3.3, frame rate selection in image acquisition: To prevent excessive inter-frame displacement, the following conditions need to be met: the sub-area size is 25 pixels and the image size is 20 pixels·mm -1 The maximum displacement between two adjacent images should be less than 1.25 mm. The sample movement speed should be less than 1 mm·s -1 , the minimum frame rate is set to 0.8Hz.

[0026] Step 3.4, establish a measure of similarity between the reference subset and the current subset, which is achieved by comparing the grayscale value of the final reference subset point with the grayscale value of the final current subset point. The following two indicators are two indicators in digital image analysis, and the calculation formula is as follows:

[0027]

[0028]

[0029] Where f and g are the reference image function and the current image function, respectively, and return the grayscale value corresponding to the specified (x, y) point. FM and GM are the average grayscale values ​​of the final reference and current subset, respectively, and are defined as follows:

[0030]

[0031] Where n(S) is the number of elements in S.

[0032] Step 3.5, nonlinear optimization, is mainly divided into three parts: providing an initial guess, Gauss-Newton (GN) iterative optimization scheme, and Lagrangian interpolation. Initial guess, because the iterative optimization scheme converges to a local maximum / minimum; therefore, an initial guess needs to be made near the global maximum / minimum. In the GN iterative optimization scheme section, two different schemes are discussed, the forward addition Gauss-Newton method (FA-GN) and the inverse combined Gauss-Newton method (IC-GN). The FA-GN method is discussed first, because this is a commonly used algorithm in digital image analysis (and optimization in general). On the other hand, the IC-GN method is a special case of nonlinear optimization (usually used in image alignment algorithms), which is faster than the standard FA-GN method and has recently been applied to digital image analysis. Finally, interpolation is discussed, because grayscale values ​​and gradients at sub-pixel locations are needed to calculate the quantities used in the FA-GN and IC-GN methods.

[0033] Step 3.6, Gauss-Newton nonlinear iterative least squares method, this scheme can be inferred to be an optimization by finding the roots of the derivatives of the function. Furthermore, it can be generalized to multivariate optimization by using gradients instead of derivatives and then determining where the gradient norm converges to zero.

[0034] Step 3.7, Obtaining Full-Field Displacement, enables the calculation of displacement data for a single material point located at the center of the initial reference subset. The next step is to obtain a complete array of displacement values. Typically, this is accomplished by first selecting a region of interest and then determining the displacement data in a grid within the region of interest. The displacements are then reduced or interpolated to form a "continuous" displacement field. Figure 7 A sample with a complex region of interest is shown, in this case a sample with a crack. The region of interest is the white highlighted area. The sample subset is superimposed on the region of interest. The subset near the crack is a special case that shows that the subset needs to be continuous within the region of interest.

[0035] Step 3.8. Obtain the Green-Lagrangian strain by using the four displacement gradients as follows:

[0036]

[0037] The displacement gradients are actually obtained directly by the IC-GN scheme. However, these values ​​are noisy, so they must be "smoothed" in some way before calculating the strain field. One possible approach is to use a low-pass filter to smooth the displacement gradient field before using the above equation. The method used, which uses a least squares plane fit of a subset of the displacement data (u and v) to find the plane parameters in the following equation.

[0038]

[0039] The idea here is to form an overconstrained system of equations of the form:

[0040]

[0041] Further, the specific steps of step 4 are as follows:

[0042] Step 4.1: Reliability-oriented initial value transfer strategy algorithm (such as Fig. 9 ) is accomplished by first selecting a seed point, calculating the corresponding deformation parameters of the surrounding subset centered at this point, and then marking each point that has been traversed and deactivating it.

[0043] Step 4.2, the subsets form a queue, and then the deformation parameters and the position of the subset center are stored in the queue. The program will enter a while loop, and at each iteration, the top element in the queue will be copied and deleted from the queue, and then the data in the loaded queue will be added to the data graph.

[0044] Step 4.3. Next, use the displacement data from the loaded queue point as an initial guess for the nonlinear optimization scheme and analyze the four surrounding pixels (if they are not already marked or within the range of interest). Repeat this procedure until the queue is empty, at which point all consecutive points of interest will be analyzed.

[0045] Step 4.4, in the adjacent calculation stage, the white cross in the middle is the material point loaded from the queue. The black cross is the unanalyzed point added to the queue. The gray cross has been analyzed and is skipped. This process has high accuracy. The poor data points with low CLS are processed first, where CLS is a deformation parameter, thereby preventing the data from being used as an initial guess for adjacent points. In addition, the smaller the subset spacing, the better the initial guess provided by the adjacent subset because the displacement field is smooth, the region of interest is partitioned, and multiple seed points are selected for parallel processing. Fig. 9 It is a reliability-oriented initial value transfer strategy algorithm flow chart.

[0046] Further, the specific steps of step 5 are as follows:

[0047] Step 5.1: Select several groups of symmetrical points H about the middle crack near the crack. 0 and H 1 , transform the displacement of the symmetric point in the uv coordinate into the newly established u'-v' coordinate system (such as Fig.10 As shown), through point H 0 and H 1 The u' and v' of the cracks are divided into three types (such as Fig.11 As shown in Figure 2): When only the u' of the symmetrical point is positive and negative, it means that the points on both sides of the crack move away from the crack in the vertical direction of the crack, and the crack is a tensile crack; when the symmetrical point H 0 and H 1 When the difference of v' is greater than 0.15, it means that the points on both sides of the crack are moving in the opposite direction parallel to the crack direction, and the crack is a shear crack. 0 and H 1 When u' is one positive and one negative and the difference of v' is >0.15, the crack is a tension-shear composite crack.

[0048] u'=u·sinβ+v·cosβ (12)

[0049] v'=u·sinβ+v·sinβ (13)

[0050] In summary, due to the adoption of the above technical solution, the beneficial effects of the present invention are:

[0051] By identifying different types of laminae and different crack types, the stress and strain analysis of crack propagation can be performed more accurately. The stress and strain can be calculated according to different laminae, and the crack propagation can be predicted by identifying different crack types. Through the initial guess at the beginning, the algorithm is further optimized by nonlinear optimization. A series of methods such as Gauss-Newton nonlinear iterative least squares method, forward addition, and inverse combination method are used to finally obtain the full-field displacement and full-field strain of different media, and obtain strain images, two-dimensional cloud maps, etc. Finally, the rock laminae and crack types are identified to achieve more accurate stress and strain analysis, which is of great significance to practical engineering problems such as the degree of wellbore collapse during drilling and the prediction of the expansion scale of artificial cracks in fracturing projects. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] Figure 1 This is a digital image analysis flow chart of the present invention

[0053] Figure 2 Schematic diagram of irregular characteristic distribution of laminae in the present invention

[0054] Figure 3 The grayscale image of the laminae of the present invention

[0055] Figure 4 Flow chart of the algorithm for identifying lamination types in the present invention

[0056] Figure 5 Schematic diagram of image acquisition for mechanical experiment of the present invention

[0057] Figure 6 Schematic diagram of the regular characteristic distribution of the laminae of the present invention

[0058] Figure 7 Schematic diagram of selecting samples of interest in the region with cracks in the present invention

[0059] Figure 8 The reliability-oriented initial value transfer strategy algorithm of the present invention

[0060] Fig. 9 The reliability-oriented initial value transfer strategy algorithm flow chart of the present invention

[0061] Fig.10 Mechanical analysis of the crack sites on both sides of the present invention

[0062] Fig.11 The mechanical models of different types of cracks in this invention are

[0063] Fig.12Digital image analysis strain map of the present invention

[0064] Fig.13 This is the crack type judgment diagram of the present invention Specific implementation methods

[0065] Step 1: Cut the core of φ56.6mm×80mm into a 40mm×40mm×80mm cuboid rock sample, grind and polish each surface, and then scan the end face, observe the end face of the core, and test the core parameters. The laminae type recognition result and grayscale image are obtained through the laminae recognition algorithm, the laminae distribution characteristics of the sample are determined, and the rock samples that meet the conditions are selected for the next experiment, and the laminae type recognition result and grayscale image are obtained through the algorithm.

[0066] Specifically, the details of step 1 are as follows:

[0067] Step 1.1: Determine the laminae distribution characteristics of the sample. The following conditions must be met: Figure 6 As shown, D represents the thickness of the shale matrix between two adjacent laminae, and d represents the thickness of a single laminae. The experimental samples are selected to ensure that: Max(d)<Max(D / 2), ∑d is less than 0.3*∑D. When the actual laminae are not parallel to the rock samples, the Monte Carlo algorithm is used to calculate the irregular laminae area, as shown in Figure 2 As shown, the following conditions need to be met: S represents the area of ​​the shale matrix between two adjacent laminae, and s represents the area of ​​a single laminae. By Max(s)<Max(S / 2), ∑s is less than 0.3*∑S.

[0068] After screening and obtaining rock samples that meet the conditions, the core grayscale image is first obtained, and its conversion formula is:

[0069] Gray=0.299*R+0.587*G+0.114*B (1)

[0070] Where Gray is the gray value of the core gray image; R is the red component in the core gray image; G is the green component in the core gray image; and B is the blue component in the core gray image.

[0071] Then, the grayscale image is enhanced and the laminae features are extracted using a two-dimensional convolution model. The grayscale curve is extracted at the quartile of the core image. After the grayscale curve is used to construct the average grayscale curve, the average grayscale curve is subjected to wavelet transform to remove noise signals. The grayscale difference between the average grayscale curves before and after the wavelet transform is calculated to obtain the grayscale difference curve. Cluster analysis is then performed using the average grayscale value and grayscale difference before and after the wavelet transform as the eigenvalues. The laminae type recognition results are obtained by identifying dark laminae in the first round and light laminae and dark laminae in the second round.

[0072] Step 1.3, use Python database to calculate the gray value of each row of rock image at the quartile and convert it into integer type to get the quartile, store the depth and quartile of the rock image in Python database, generate three gray curves, normalize the average gray curve and gray difference curve, set the number of clusters of K-Means algorithm, use K-Means algorithm to perform cluster analysis on the standardized eigenvalues, calculate the adjacent depth difference to represent the thickness of the laminae, and save the clustering results, depth difference and eigenvalue obtained by cluster analysis into Python database, finally draw the difference frequency histogram according to the clustering results to get the gray map (such as Figure 3 The specific implementation steps are as follows: Figure 4 shown.

[0073] Step 2, speckle production: polish and grind the surface of the rectangular rock sample, and produce speckle on the surface of the rock sample, then conduct mechanical tests and collect images.

[0074] Specifically, the details of step 2 are as follows:

[0075] Step 2.1, speckle processing: speckle production is performed on the surface of the laminated shale, using black speckle, and the speckle density ratio is approximately 50%, that is, the areas occupied by bright and dark pixels in any sub-area are approximately the same.

[0076] Step 2.2, experimental image acquisition preprocessing: When conducting uniaxial compression experiments on laminated shale, two white light sources are used to continuously illuminate the opposite side of the rock at an interval of 60°. The sample surface should be perpendicular to the camera optical axis, the ambient temperature should be controlled at 15℃~25℃, and the shooting environment should have no windows or observation windows.

[0077] Step 3: Image acquisition: While conducting the mechanical test experiment, a high-speed camera is used to capture synchronous images to obtain a set of images containing displacements. The acquired image data and a series of experimental data obtained from the mechanical test experiment are analyzed.

[0078] Specifically, the details of step 3 are as follows:

[0079] Step 3.1, Image Analysis requires using optical images (excluding images taken by scanning electron microscopy, atomic force microscopy, or X-rays), adjusting high-speed camera settings (camera resolution needs to be greater than 1080p HD / 30fps) and capturing synchronized images.

[0080] Step 3.2: When conducting mechanical experiments, the loading rate must be controlled within 0.5-1.0 MPa / s.

[0081] Step 3.3, frame rate selection in image acquisition: To prevent excessive inter-frame displacement, the following conditions need to be met: the sub-area size is 25 pixels and the image size is 20 pixels·mm -1 The maximum displacement between two adjacent images should be less than 1.25 mm. The sample movement speed should be less than 1 mm·s -1 , the minimum frame rate is set to 0.8Hz.

[0082] Step 4, calculate crack displacement data: In the displacement image group, select the region of interest and determine the displacement data in the grid within the region. The reliability-oriented initial value transfer strategy algorithm is used to obtain the crack displacement data.

[0083] Specifically, the details of step 4 are as follows:

[0084] Step 4.1: Reliability-oriented initial value transfer strategy algorithm (such as Fig. 9 ) is accomplished by first selecting a seed point, calculating the corresponding deformation parameters of the surrounding subset centered at this point, and then marking each point that has been traversed and deactivating it.

[0085] Step 4.2, the subsets form a queue, and then the deformation parameters and the position of the subset center are stored in the queue. The program will enter a while loop, and at each iteration, the top element in the queue will be copied and deleted from the queue, and then the data in the loaded queue will be added to the data graph.

[0086] Step 4.3. Next, use the displacement data from the loaded queue point as an initial guess for the nonlinear optimization scheme and analyze the four surrounding pixels (if they are not already marked or within the range of interest). Repeat this procedure until the queue is empty, at which point all consecutive points of interest will be analyzed.

[0087] Step 4.4, in the adjacent calculation stage, the white cross in the middle is the material point loaded from the queue. The black cross is the unanalyzed point added to the queue. The gray cross has been analyzed and is skipped. This process has high accuracy and first processes poor data points, such as data points with low CLS (CLS is a deformation parameter), thereby preventing the data from being used as an initial guess for adjacent points. In addition, the initial guess provided by the adjacent subset usually works quite well (when the subset spacing is lower than the average spacing) because the displacement field is usually relatively smooth. Partition the region of interest and choose to process multiple seed points in parallel. Fig. 9 It is a reliability-oriented initial value transfer strategy algorithm flow chart.

[0088] Step 5: According to the displacement trend and displacement of the pixels on both sides of the crack data, three crack type standards are constructed, a new method for crack type discrimination is established, and different crack types are identified.

[0089] Specifically, the details of step 5 are as follows:

[0090] Step 5.1, reliability-oriented initial value transfer strategy algorithm (such as Fig. 9 ) is accomplished by first selecting a seed point, calculating the corresponding deformation parameters of the surrounding subset centered at this point, and then marking each point that has been traversed and deactivating it.

[0091] Step 5.2, the subsets form a queue, and then the deformation parameters and the position of the subset center are stored in the queue. The program will enter a while loop, and at each iteration, the top element in the queue will be copied and deleted from the queue, and then the data in the loaded queue will be added to the data graph.

[0092] Step 5.3. Next, use the displacement data from the loaded queue point as an initial guess for the nonlinear optimization scheme and analyze the surrounding four pixels (if they are not already marked or within the range of interest). Repeat this procedure until the queue is empty, at which point all consecutive points of interest will be analyzed.

[0093] Step 5.4, in the adjacent calculation stage, the white cross in the middle is the material point loaded from the queue. The black cross is the unanalyzed point added to the queue. The gray cross has been analyzed and is skipped. This process has a high accuracy and first processes the poor data points, such as those with low CLS (CLS is a deformation parameter), thereby preventing the data from being used as an initial guess for the adjacent points. In addition, the initial guess provided by the adjacent subset usually works quite well (as long as the subset spacing is relatively small) because the displacement field is usually relatively smooth. The only problem with this algorithm is that it must be done serially. One way to solve this problem is to partition the area of ​​interest and choose to process multiple seed points in parallel.

[0094] Step 6: According to the displacement trend and displacement of the pixels on both sides of the crack data, three crack type standards are constructed, a new method for crack type discrimination is established, and different crack types are identified.

[0095] Specifically, the details of step 6 are as follows:

[0096] Select several groups of symmetrical points H about the middle crack near the crack. 0 and H 1 , transform the displacement of the symmetric point in the uv coordinate into the newly established u'-v' coordinate system (such as Fig.10 As shown), through point H 0 and H 1 The u' and v' of the cracks are divided into three types (such as Fig.11As shown in Figure 2): When only the u' of the symmetrical point is positive and negative, it means that the points on both sides of the crack move away from the crack in the vertical direction of the crack, and the crack is a tensile crack; when the symmetrical point H 0 and H 1 When the difference of v' is greater than 0.15, it means that the points on both sides of the crack are moving in the opposite direction parallel to the crack direction, and the crack is a shear crack. 0 and H 1 When u' is one positive and one negative and the difference of v' is >0.15, the crack is a tension-shear composite crack.

[0097] u'=u·sinβ+v·cosβ (12)

[0098] v'=u·sinβ+v·sinβ (13)

[0099] Example 1

[0100] (1) This embodiment uses a laminar shale rock sample as the experimental object, uses a square core experiment, and the core size is φ56.6 mm×80 mm, which is then cut into a 40 mm×40 mm×80 mm rectangular specimen. By scanning the end face of the laminar shale core, a two-dimensional convolution model is used to extract and enhance the laminar features to obtain the laminar type recognition result and obtain a grayscale image (such as Figure 3 As shown), and select the laminated shale that meets the requirements.

[0101] (2) Speckle production: Select black speckle and the speckle density ratio is approximately 50%, that is, the areas occupied by bright and dark pixels in any sub-region are approximately the same.

[0102] (3) Experimental image acquisition: The loading rate is controlled at 0.5-1.0 MPa / s. Two white light sources are used with an interval of 60° to continuously illuminate the opposite side of the rock. The sample surface should be perpendicular to the optical axis of the camera (e.g. Figure 5 The uniaxial compression test of laminated shale was carried out with the ambient temperature controlled at 15℃~25℃ and no special environmental conditions (such as no window or observation window, no explosion or shock wave). The synchronous images were captured using a high-speed camera with a camera resolution of 1080p HD / 30fps, a sub-area size of 25pixel and an image size of 20pixel·mm -1 The maximum displacement between two adjacent images is 1.25 mm. The sample movement speed is 1 mm·s -1 , the minimum frame rate is set to 0.8Hz.

[0103] (4) Digital image analysis is performed on the obtained images to calculate the displacement data of a single material point located in the region of interest. The region of interest is first selected, and then the displacement data in the grid is determined within the region of interest. The reliability-oriented initial value transfer strategy algorithm is used (the algorithm principle is as follows Figure 8 As shown), the full-field displacement is obtained and the digital image analysis strain image is obtained, as shown Fig.12 shown.

[0104] (5) According to the displacement trend and displacement of the pixels on both sides of the crack, three types of crack types are constructed, namely tensile crack, shear crack and tensile-combined-shear crack. By identifying and calculating the u' of the symmetrical points H0 and H1 in the pixel points, one is positive and the other is negative, and the difference of v' is 0.16, which meets the difference value > 0.15, it is obtained that this crack is a tensile-combined-shear crack, such as Fig.13 shown.

[0105] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A method for identifying different types of crack extension in laminated shale based on digital images, characterized in that: The following steps are involved: Step 1, laminae identification: cut the core with the size of φ56.6mm×80mm into a cuboid rock sample of 40mm×40mm×80mm, and grind and polish each surface, then perform end face scanning, observe the end face of the core, test the core parameters, obtain the laminae type identification result and grayscale image through the laminae identification algorithm, determine the laminae distribution characteristics of the sample, and select the rock sample that meets the conditions for the next step of the experiment; Step 2, speckle production: polishing the surface of the rectangular rock sample, and producing speckles on the surface of the rock sample, followed by mechanical testing and image acquisition; Step 3, image acquisition: while conducting the mechanical test experiment, a high-speed camera is used to capture synchronous images to obtain a set of images containing displacements, and the collected image data and a series of experimental data obtained from the mechanical test experiment are analyzed; Step 4, calculating crack displacement data: in the displacement image group, the region of interest is selected to determine the displacement data in the grid within the region, and the reliability-oriented initial value transfer strategy algorithm is used to obtain the crack displacement data; Step 5: According to the displacement trend and displacement of the pixels on both sides of the crack data, three crack type standards are constructed, a new method for crack type discrimination is established, and different crack types are identified; In step 4, when calculating the crack displacement data, the reliability initial value transfer strategy algorithm used specifically includes: Step 4.1, the reliability-oriented initial value transfer strategy algorithm is to first select a seed point, calculate the corresponding deformation parameters of the surrounding subset centered on this point, and then mark each point that has been traversed and deactivate it; Step 4.2, the subsets form a queue, and then the deformation parameters and the position of the subset center are stored in the queue. The program will enter a while loop. At each iteration, the top element in the queue will be copied and deleted from the queue. Then, the data in the loaded queue will be added to the data map; Step 4.

3. Next, use the displacement data from the loading queue points as the initial guess for the nonlinear optimization scheme, and analyze the four surrounding pixels when they are not yet marked or within the range of interest. Repeat the procedure until the queue is empty, at which point all consecutive points of interest will be analyzed. Step 4.4, in the adjacent calculation stage, the middle cross is the material point loaded from the queue, and the surrounding crosses are the unanalyzed points added to the queue. The marks have been analyzed and skipped. This process has higher accuracy. The poor data points with lower CLS are processed first, where CLS is the deformation parameter, thereby preventing the data from being used as the initial guess of the adjacent points. In addition, the smaller the subset spacing, the better the initial guess provided by the adjacent subset because the displacement field is smooth, the region of interest is partitioned, and multiple seed points are selected for parallel processing.

2. The method for identifying different types of crack extensions in laminated shale based on digital images according to claim 1, characterized in that: The step 1 includes the following steps: Step 1.1: Determine the laminae distribution characteristics of the sample. When the shale matrix between two adjacent laminae is parallel, the following conditions must be met: Max(d)<Max(D / 2) (1) Σd<0.3*ΣD (2) In the formula, the thickness of the shale matrix between two adjacent laminae is D, and the thickness of a single laminae is d. When the actual laminae are not parallel to the rock sample, the Monte Carlo algorithm is needed to calculate the irregular laminae area, and the following conditions must be met: Max(s)<Max(S / 2) (3) Σs<0.3*ΣS (4) In the formula, the area of ​​shale matrix between two adjacent laminae is S, and the area of ​​a single laminae is s; Step 1.2: After screening and obtaining rock samples that meet the conditions, first obtain the grayscale image of the rock core, and its conversion formula is: Gary=0.299*R+0.587*G+0.114*B (5) Where Gray is the gray value of the core gray image; R is the red component in the core gray image; G is the green component in the core gray image; B is the blue component in the core gray image; Then, the grayscale image is enhanced and the laminae features are extracted using a two-dimensional convolution model. The grayscale curve is extracted at the quartile of the core image. After the grayscale curve is used to construct the average grayscale curve, the average grayscale curve is subjected to wavelet transform to remove the noise signal. The grayscale difference between the average grayscale curve before and after the wavelet transform is calculated to obtain the grayscale difference curve. The average grayscale value and grayscale difference before and after the wavelet transform are used as the three eigenvalues ​​for cluster analysis. The dark laminae are identified in the first round, and the light laminae and dark laminae are identified in the second round, and the laminae type recognition results are obtained. Step 1.3, use the Python database to calculate the grayscale value at the quartile of each row of the rock image and convert it to an integer type to obtain the quartile, store the depth and quartile of the rock image in the Python database, generate three grayscale curves, normalize the average grayscale curve and the grayscale difference curve, set the number of clusters of the K-Means algorithm, use the K-Means algorithm to perform cluster analysis on the standardized eigenvalues, calculate the adjacent depth difference to represent the thickness of the laminae, and save the clustering results, depth difference and eigenvalue obtained by the cluster analysis into the Python database, and finally draw a difference frequency histogram based on the clustering results to obtain a grayscale image.

3. The method for identifying different types of crack extensions in laminated shale based on digital images according to claim 1, characterized in that: The step 2 includes the following steps: Step 2.1, speckle processing: speckle production is performed on the laminar shale surface, using black speckle and a speckle density ratio of 50%, that is, the areas occupied by bright and dark pixels in any sub-area are the same; Step 2.2, experimental image acquisition preprocessing: When conducting uniaxial compression experiments on laminated shale, two white light sources are used to continuously illuminate the opposite side of the rock at an interval of 60°. The sample surface should be perpendicular to the camera optical axis, the ambient temperature should be controlled at 15℃~25℃, and the shooting environment should have no windows or observation windows.

4. The method for identifying different types of crack extensions in laminated shale based on digital images according to claim 1, characterized in that: The step 3 includes the following steps: Step 3.1: Image analysis requires the use of optical images, which do not include images from scanning electron microscopes, atomic force microscopes, or X-rays. Adjust the high-speed camera settings, where the camera resolution needs to be greater than 1080p HD / 30fps, and capture synchronized images at the same time. Step 3.2: When conducting mechanical experiments, the loading rate must be controlled within 0.5-1.0 MPa / s; Step 3.3, frame rate selection in image acquisition: To prevent excessive inter-frame displacement, the following conditions need to be met: the sub-area size is 25 pixels and the image size is 20 pixels·mm -1 The maximum displacement between two adjacent images should be less than 1.25 mm, and the sample movement speed should be less than 1 mm·s -1 , the minimum frame rate is set to 0.8Hz.

5. The method for identifying different types of crack extensions in laminated shale based on digital images according to claim 1, characterized in that: The new method for determining the crack type in step 5 includes the following features: For the identification of crack types, several groups of symmetrical points H0 and H1 about the middle crack are selected along both sides of the crack. The displacement of the symmetrical points in the uv coordinate is converted to the newly established u,-v, coordinate system through calculation by the formula. The conversion formula is as follows: u'=u·sinβ+v·cosβ (6) v'=u·sinβ+v·sinβ (7) In the formula, β represents the angle between the crack and the horizontal, u represents the displacement of the pixel point under the horizontal coordinate, and v represents the displacement of the pixel point under the vertical coordinate; The cracks are divided into three types by u, and v at points H0 and H1: Tensile crack: When the u of only symmetrical points is one positive and one negative, it means that the points on both sides of the crack move away from the crack in the vertical direction of the crack; Shear crack: When the difference of v between the symmetrical points H0 and H1 is > 0.15, it means that the points on both sides of the crack move in the opposite direction parallel to the crack direction; Tensile-shear composite crack: When the v of the two symmetrical points H0 and H1 are one positive and one negative and the difference of v is >0.15.

Images