CN113468639A - Method for establishing fracture network three-dimensional visualization model based on occurrence state - Google Patents

Method for establishing fracture network three-dimensional visualization model based on occurrence state Download PDF

Info

Publication number
CN113468639A
CN113468639A CN202110699601.6A CN202110699601A CN113468639A CN 113468639 A CN113468639 A CN 113468639A CN 202110699601 A CN202110699601 A CN 202110699601A CN 113468639 A CN113468639 A CN 113468639A
Authority
CN
China
Prior art keywords
data
disc
fracture
attitude
model
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202110699601.6A
Other languages
Chinese (zh)
Inventor
张云辉
姚荣文
李亚博
王语雁
覃礼貌
王鹰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Southwest Jiaotong University
Original Assignee
Southwest Jiaotong University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Southwest Jiaotong University filed Critical Southwest Jiaotong University
Priority to CN202110699601.6A priority Critical patent/CN113468639A/en
Publication of CN113468639A publication Critical patent/CN113468639A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/10—Geometric CAD
    • G06F30/13—Architectural design, e.g. computer-aided architectural design [CAAD] related to design of buildings, bridges, landscapes, production plants or roads
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G06F18/232—Non-hierarchical techniques
    • G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G06F18/232—Non-hierarchical techniques
    • G06F18/2321—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions
    • G06F18/23213—Non-hierarchical techniques using statistics or function optimisation, e.g. modelling of probability density functions with fixed number of clusters, e.g. K-means clustering
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/27—Design optimisation, verification or simulation using machine learning, e.g. artificial intelligence, neural networks, support vector machines [SVM] or training a model
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/08—Probabilistic or stochastic CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling

Abstract

The invention relates to the technical field of fracture networks, in particular to a method for establishing a fracture network three-dimensional visualization model based on occurrence, which comprises the following steps: 1) carrying out digital processing on fracture data; converting the attitude data into a spherical coordinate representation; expressing the attitude of the structural surface by using a normal vector; projecting a normal vector of the structural plane onto a red plane by adopting a dimension reduction method and utilizing a red plane projection to obtain a polar point diagram; realizing automatic grouping by using a clustering algorithm and an elbow method, and obtaining pole centers; 2) simulating fracture data; 3) visual modeling; obtaining a disc coordinate by using a coordinate transformation method, and constructing a disc by using a triangular surface; the three-dimensional visualization of the crack is carried out by utilizing the disc, the direction of the disc is represented by the occurrence shape, the size of the disc is represented by the trace length, and the spatial position of the position representation disc is measured. The invention can preferably realize three-dimensional visual modeling of the fracture network by the disc model.

Description

Method for establishing fracture network three-dimensional visualization model based on occurrence state
Technical Field
The invention relates to the technical field of fracture networks, in particular to a method for establishing a fracture network three-dimensional visualization model based on occurrence.
Background
The crack is formed by the fracture and deformation of hard rock under the action of stress, has important significance in the research of the fields of engineering geology, tectonic geology, hydrogeology, petroleum geology and the like, and has relation to slope stability, tectonic stress field distribution, underground water migration, petroleum migration and the like. Zhang (2013) visualizes the structural features of fractured rock mass and evaluates the rock mass structure of the Wudongde hydropower station. Field geological surveys usually reveal only a limited structural surface, so it is necessary to build fracture network models. In field geological survey, a window measurement method is often adopted to measure shallow fracture data, wherein occurrence is important data, but at present, rosettes, Wu's network and the like are mostly adopted to characterize fracture occurrence, so that the limitation is large, and the development characteristics of a regional fracture network are difficult to visually reflect. Therefore, it is necessary to write a computer program to analyze the occurrence data of the fractures, construct a numerical simulation model and realize visualization of the fracture network, and theoretical basis and corresponding reference can be provided for engineering construction, railway and highway line selection, underground water research and the like.
Disclosure of Invention
It is an object of the present invention to provide a method for building a three-dimensional visualization model of a fracture network based on attitude that overcomes some or some of the deficiencies of the prior art.
The invention relates to a method for establishing a fracture network three-dimensional visualization model based on occurrence, which comprises the following steps:
(1) carrying out digital processing on fracture data;
(1.1) converting the attitude data into a spherical coordinate representation;
(1.2) expressing the attitude of the structural plane by using a normal vector; projecting a normal vector of the structural plane onto a red plane by adopting a dimension reduction method and utilizing a red plane projection to obtain a polar point diagram;
(1.3) realizing automatic grouping by utilizing a clustering algorithm and an elbow method, and obtaining a pole center;
(2) simulating fracture data;
simulating the occurrence data by using a Fisher model, obtaining the average trace length of known data, obtaining simulation data of the radius of the disc by using negative exponential distribution, and simulating the position data of the fracture by using uniform distribution;
(3) visual modeling;
obtaining a disc coordinate by using a coordinate transformation method, and constructing a disc by using a triangular surface; the three-dimensional visualization of the crack is carried out by utilizing the disc, the direction of the disc is represented by the occurrence shape, the size of the disc is represented by the trace length, and the spatial position of the position representation disc is measured.
Preferably, the step (1) specifically comprises the following steps:
(1.1) the occurrence is converted into a spherical coordinate representation with a radius r of 1 according to the following formula:
Figure BDA0003129741500000021
theta is a horizontal angle of the horizontal,
Figure BDA0003129741500000022
is a vertical angle;
(1.2) the attitude of the structural plane is expressed by a normal vector, and the direction cosine of the normal vector can be calculated by the following formula:
Figure BDA0003129741500000023
projecting a normal vector of the structural plane onto a red plane by adopting a dimension reduction method and utilizing a red plane projection to obtain a polar point diagram; and (3) performing the erythro-plano projection by using the following formula to obtain the plane rectangular coordinates of each pole:
Figure BDA0003129741500000024
(1.3) realizing automatic grouping by using a clustering algorithm and an elbow method, obtaining pole centers, regarding the centers obtained by clustering as the average vector direction of the occurrence, and marking the pole center positions obtained by clustering as
Figure BDA0003129741500000025
Preferably, in the step (2), the fracture data simulation comprises the following specific steps:
according to the central coordinates of each grouping pole
Figure BDA0003129741500000026
Center attitude of the available fractures
Figure BDA0003129741500000027
Taking one set of data, recording the center as
Figure BDA0003129741500000031
Having a direction cosine of
Figure BDA0003129741500000032
Let the group have n data, each data having a direction cosine of (l)j,mj,nj) (ii) a The degree of dispersion can be expressed as κ, with the larger the κ value, the more centered the data is, when κ>At 6 th hour[9]This can be calculated using the following formula:
Figure BDA0003129741500000033
the attitude data obeys a Fisher distribution, and the probability density function thereof is shown as the following formula, wherein the center of the formula is the origin or z-axis direction of the declination projection:
Figure BDA0003129741500000034
the Fisher probability density function is a integrable function, and the integral of the above formula is obtained
Figure BDA0003129741500000035
θ*Cumulative density function of (d):
Figure BDA0003129741500000036
since the cumulative density function follows the uniform distribution of the (0,1) interval, the uniform distribution U-U (0,1) is introduced,
substituting u into the above formula and negating the function to obtain a sampling formula:
Figure BDA0003129741500000037
the sampling formula is utilized to generate a large amount of occurrence data which obey Fisher distribution
Figure BDA0003129741500000038
Can calculate the direction cosine
Figure BDA0003129741500000039
Coordinate transformation is performed by the following formula to obtain required data (l)i,mi,ni)：
Figure BDA00031297415000000310
Then the calculation is carried out reversely to obtain
Figure BDA00031297415000000311
And (alpha)j,βj) (ii) a Therefore, simulation of fracture occurrence data is completed;
the fracture path length L and the disc radius r have the following relationship:
Figure BDA00031297415000000312
from the negative exponential distribution of the fissure vestige, the relationship of the following formula is obtained:
Figure BDA00031297415000000313
integrating and solving a reverse function to obtain a sampling formula of the disc radius:
Figure BDA0003129741500000041
the fracture position data is then simulated with a uniform distribution.
As a preference, the first and second liquid crystal compositions are,
Figure BDA0003129741500000042
in addition to being available from clustering algorithms, the following can be calculated:
Figure BDA0003129741500000043
preferably, the degree of dispersion can also be iterated by a dichotomy, resulting in a more accurate value from the maximum likelihood estimation.
Preferably, before generating the data using the sampling formula, χ of the data is required2Fitting and testing, wherein points with large individual errors can be removed in the testing process, and if the tested data obeys Fisher distribution, the data can be generated by using the model; non-compliance requires the use of other models.
Preferably, in the case of a mixture of a plurality of Fisher models, we consider the mixture of Fisher models with different average occurrence and different dispersion degrees, that is, the average occurrence and the dispersion degree are parameters, and then we can obtain a Fisher mixture model shown in the following formula; the K-Means algorithm is added, and the solution can be carried out continuously by iteration;
Figure BDA0003129741500000044
in addition, the maximum likelihood estimation value can also be solved by utilizing the EM algorithm of the Gaussian mixture model in an iterative mode.
Preferably, in the step (3), the specific steps of visual modeling are as follows:
firstly, analyzing a single crack, setting a structural plane as a disc, establishing a disc with the center of circle positioned at an origin and the size of r on a two-dimensional plane, and obtaining a fracture model with horizontal attitude according to the analytical formula:
Figure BDA0003129741500000045
on the basis of the model, the circle is wound around the y-axis by the rotation angle according to the following formula
Figure BDA0003129741500000046
Then obtain the inclination angle of
Figure BDA0003129741500000047
The structural surface disc of (1);
Figure BDA0003129741500000051
the disc is rotated around the z-axis by an angle theta to obtain an inclination angle theta
Figure BDA0003129741500000053
The disk of (2):
Figure BDA0003129741500000052
then translating the disc to the collected position coordinates (a, b, h) to obtain the disc at any spatial position, thereby completing the three-dimensional visual modeling of a single fracture;
after the disks of a plurality of fractures are superposed, a three-dimensional model of the fracture network can be obtained.
The method realizes grouping and simulation of fracture occurrence data based on a K-Means clustering algorithm and a Fisher distribution model, and realizes three-dimensional visual modeling of a fracture network by a disc model. The invention has the advantages that:
the application of the spherical coordinate and the plano-declination projection realizes the dimension reduction and the digital processing of the occurrence data, and is convenient for the computer processing.
The invention can improve the convergence speed and the operation speed and better realize the grouping of the occurrence data.
The method is based on statistics and Monte Carlo principle, successfully utilizes Fisher model to generate a large amount of data, and provides data dependence for constructing the fracture network.
According to the invention, the three-dimensional model of the fracture network is successfully constructed by utilizing the triangular surface in the PyQtGraph, and a certain reference is provided for the future research of the fracture network.
Drawings
FIG. 1 is a flowchart of a method for building a three-dimensional visualization model of a fracture network based on occurrence status in example 1;
FIG. 2 is a schematic view of a structural plane in spherical coordinates in example 1;
FIG. 3 is a polar diagram in example 1;
FIG. 4 is a schematic view of the elbow rule in example 1;
FIG. 5 is a diagram showing grouping of extreme points in example 1;
FIG. 6 is a schematic view of a simulation of the occurrence of example 1;
FIG. 7 is a schematic view of a fracture model of the level of production in example 1;
FIG. 8 shows the inclination angle of example 1
Figure BDA0003129741500000061
Structural surface disc schematic of (1);
FIG. 9 shows the inclination angle θ of example 1
Figure BDA0003129741500000062
A disc schematic of (a);
FIG. 10 is a schematic diagram of a three-dimensional model of a fracture network in example 1;
FIG. 11 is a polar diagram in example 2;
FIG. 12 is a schematic view of an elbow rule in embodiment 2;
FIG. 13 is a diagram showing grouping of extreme points in example 2;
FIG. 14 is a schematic view of a simulation of the occurrence of example 2;
FIG. 15 is a schematic view of an actual fracture in example 2;
FIG. 16 is a schematic view of a simulated fracture in example 2;
FIG. 17 is a schematic diagram of an actual fracture of the actual fractures of PD4 and PD6 in example 2;
FIG. 18 is a schematic view of simulated fractures of PD4 and PD6 in example 2;
FIG. 19 is a polar view of all the footrill data in example 2;
figure 20 is a schematic view of all the pilot-modeled fractures of example 2.
Detailed Description
For a further understanding of the invention, reference should be made to the following detailed description taken in conjunction with the accompanying drawings and examples. It is to be understood that the examples are illustrative of the invention and not limiting.
Example 1
As shown in fig. 1, the present embodiment provides a method for building a three-dimensional visualization model of a fracture network based on a parturition, which includes the following steps:
(1) carrying out digital processing on fracture data;
(1.1) converting the attitude data into a spherical coordinate representation;
(1.2) expressing the attitude of the structural plane by using a normal vector; projecting a normal vector of the structural plane onto a red plane by adopting a dimension reduction method and utilizing a red plane projection to obtain a polar point diagram;
(1.3) realizing automatic grouping by utilizing a clustering algorithm and an elbow method, and obtaining a pole center;
(2) simulating fracture data;
simulating the occurrence data by using a Fisher model, obtaining the average trace length of known data, obtaining simulation data of the radius of the disc by using negative exponential distribution, and simulating the position data of the fracture by using uniform distribution;
(3) visual modeling;
obtaining a disc coordinate by using a coordinate transformation method, and constructing a disc by using a triangular surface; the three-dimensional visualization of the crack is carried out by utilizing the disc, the direction of the disc is represented by the occurrence shape, the size of the disc is represented by the trace length, and the spatial position of the position representation disc is measured.
The following is a detailed description:
data processing and simulation
Fracture data digital processing
The data used in the visual modeling are the attitude, trajectory length and position of the fracture, with dip and dip in the fracture attitude represented by α and β, trajectory length represented by L, and position represented by (x, y). We calculated and analyzed the data obtained using the windowing method, as shown in table 1:
TABLE 1 fracture data
Figure BDA0003129741500000071
The first step of modeling is to carry out digital processing on the fracture data, and according to the expression method of data on the Fangketai sphere, the occurrence is converted into a spherical coordinate expression (r radius, theta horizontal angle,
Figure BDA0003129741500000073
vertical angle), for ease of calculation, we specify that r is 1, taking the unit circle, as shown in fig. 2.
Figure BDA0003129741500000072
One plane can be characterized by a normal vector, so that the normal vector is used for representing the attitude of the structural plane, and the direction cosine of the normal vector can be calculated by the formula (2-2). Because the occurrence states are not convenient to carry out statistical analysis in a three-dimensional space, a dimension reduction method is adopted, and the normal vector of the structural plane is projected onto a red plane by using a red flat projection to obtain a polar point diagram, so that statistics and calculation can be carried out on a two-dimensional plane. The upper hemisphere projection is adopted, the formula (2-3) is used for the erythro-plano projection, the plane rectangular coordinate of each pole is obtained, and the machine learning algorithm is convenient for data processing.
Figure BDA0003129741500000081
Figure BDA0003129741500000082
The polar diagram shown in fig. 3 was obtained by subjecting the data of table 1 to the declination projection according to the formula (2-3), and it was clearly found that the data in the polar diagram can be divided into two groups. Usually, we will draw a rose diagram to find the dominant occurrence and then group, and chenjianpin (2007) and the like divide the dominant occurrence based on the schmidt network, but in order to realize the automatic grouping of the computers, we adopt a clustering algorithm in machine learning to realize.
At present, five major clustering algorithms are respectively based on division, grids, density, layers and models, a K-Means clustering algorithm (K-Means clustering algorithm) is selected for grouping poles, and the clustering algorithm is based on division and has the advantages of simple algorithm, extremely high convergence rate and suitability for clustering of a large amount of data. The idea of K-Means is to select K initial clustering centers from scattered points, to mark the point with the minimum distance to a certain clustering center to the cluster of the clustering center, and to find out the proper clustering center and the corresponding cluster through continuous iteration, thereby realizing grouping. The K-Means clustering algorithm has the defects that the clustering number K and the clustering center thereof need to be selected in advance, and the clustering result is easy to fall into a local optimal solution.
The K-Means algorithm is adopted to solve the problem of trapping a local optimal solution, and the operation speed is improved. The selection of the clustering number can be solved by utilizing evaluation indexes, the conventional contour coefficient and the elbow rule at present. The elbow rule is adopted, the principle is that the cost function is defined as the sum of distances from all data points to the center of a cluster of the data points, and when the cost function is lost to the minimum, the cost function is a proper k value. As shown in fig. 4, a graph is plotted of the loss function and the different values of k, and when k <2, the loss function drops quickly; when k is greater than 2, the loss function is not changed greatly, and the meaning of continuously dividing a plurality of groups is not great; when k is 2, the elbow (turning point) is located, so the optimal classification number is 2, which is consistent with the fact that the pole point diagram judged by human is divided into two groups, and the correctness can be shown. In order to allow the computer to automatically identify the turning point of the elbow rule, the relative decrease (increase) amount can be calculated by taking the ratio of the slope of the previous segment to the slope of the next segment of the polygonal line, wherein the turning point is the largest amount of change.
The pole is easily classified into 2 groups by introducing a three-party library scimit-lean for python machine learning, and using a K-Means clustering algorithm in the library and adding an elbow rule, and the center of the pole is obtained, as shown in FIG. 5. The classification result in the figure obviously accords with the fact, and the accuracy is higher. The center obtained by clustering can be regarded as the average vector direction of the attitude, and the center obtained by clustering is recorded as
Figure BDA0003129741500000091
Fracture data simulation
Through the foregoing data processing, the center coordinates of the respective grouping poles are obtained
Figure BDA0003129741500000092
Reversely solving the central attitude of the fracture by using the formula (2-3)
Figure BDA0003129741500000093
For convenient calculation, one group of data is taken and its center is recorded as
Figure BDA0003129741500000094
Having a direction cosine of
Figure BDA0003129741500000095
Let the group have n data, each data having a direction cosine of(lj,mj,nj) Center of which
Figure BDA0003129741500000096
In addition to being obtained by a clustering algorithm, the method can also be calculated according to the formulas (2-4) and (2-5), wherein | R | is the length of vector sum and is used for averaging the direction cosines. The degree of dispersion can be represented by k, the larger the value of k, the more concentrated the data is in the center, when k is>When 6, it can be calculated by the formula (2-6). In addition, the dichotomy iteration can be used, and a more accurate value is obtained by the maximum likelihood estimation.
Figure BDA0003129741500000097
Figure BDA0003129741500000098
Figure BDA0003129741500000099
Fisher (1953) indicates that the attitude data follows a Fisher distribution with a probability density function as shown in equations (2-7), where the center of the equation is the origin or z-axis direction of the equatorial projection.
Figure BDA0003129741500000101
The Fisher probability density function is a integrable function, and integration of the formula (2-7) can be obtained
Figure BDA0003129741500000102
θ*Cumulative density function of (d):
Figure BDA0003129741500000103
obeying a uniform distribution of (0,1) intervals due to the cumulative density function[9]Therefore, a uniform distribution U-U (0,1) is introduced, U is taken into formula (2-8) and the inverse function is solved to obtain a sampling formula:
Figure BDA0003129741500000104
thus, according to the Monte Carlo simulation method, a large amount of Fisher-distributed attitude data can be generated by using sampling formula (2-9)
Figure BDA0003129741500000105
The direction cosine is obtained according to the formula (2-2)
Figure BDA0003129741500000106
However, before generating the data, χ must be performed on the data2Fitting and testing, wherein extremely individual points with extremely large errors can be removed in the testing process, and if the tested data obeys Fisher distribution, the data can be generated by using the above model; non-compliance requires the use of other models.
The data generated by the Fisher model obeys standard Fisher distribution, the center of the data is the z-axis, and the actual data is hardly concentrated on the z-axis, so that the data generated by simulation needs to be subjected to coordinate conversion and converted to the center of the actual data. Firstly, the coordinate transformation is carried out by the formula (2-10) to obtain the required data (l)i，mi,ni) Then, the reaction is performed according to the formula (2-2) and the formula (2-1) to obtain
Figure BDA0003129741500000107
And (alpha)j,βj). Thus, simulation of fracture occurrence data is completed. Obtained from the preceding
Figure BDA0003129741500000108
The results shown in fig. 6 were obtained by performing numerical simulations, and the simulation results were similar and satisfactory when compared with fig. 2-2 and fig. 2-4.
Figure BDA0003129741500000109
The crack trace length L and the disc radius r have a relation of a formula (2-11), the crack trace length is distributed according to a negative exponential to obtain a relation of a formula (2-12), and a sampling formula (2-13) of the disc radius is obtained by integrating and solving a negative function. The fracture position data is then simulated with a uniform distribution.
Figure BDA0003129741500000111
Figure BDA0003129741500000112
Figure BDA0003129741500000113
So far, the data required for fracture visualization have all been simulated.
For the case of mixing a plurality of Fisher models, the Fisher models are regarded as being mixed by Fisher models with different average occurrence and different dispersion degrees, namely the average occurrence and the dispersion degrees are parameters, so that the Fisher mixed model shown as the formula (2-14) can be obtained. And the K-Means algorithm is added, and the solution can be carried out continuously by iteration. In addition, the maximum likelihood estimation value can also be solved by utilizing the EM algorithm of the Gaussian mixture model in an iterative mode.
Figure BDA0003129741500000114
Visual modeling
The visual modeling can vividly and vividly represent the development condition of the fractures in the rock body, visually display the distribution characteristics of the regional fracture network, and provide corresponding reference for engineering construction and the like.
According to the disc model, the three-dimensional visualization of the crack is carried out by using the disc, the direction of the disc is represented by the attitude, the size of the disc is represented by the trace length, and the spatial position of the position representation disc is measured. The disc coordinate is obtained by a coordinate transformation method, and then the disc is constructed by a triangular surface.
Firstly, analyzing a single crack, assuming that the structural plane of the crack is a disc, establishing a disc with the center of the circle positioned at the origin and the size of r on a two-dimensional plane, wherein the analytic expression of the disc is shown as (3-1), and thus obtaining the fracture model with the attitude level shown in fig. 7.
Figure BDA0003129741500000121
Based on the model shown in FIG. 7, the circle is rotated around the y-axis by the rotation angle according to the formula (3-2)
Figure BDA0003129741500000125
(rotate in the positive direction as determined by the right-hand rule, the same applies below) to obtain a tilt angle of
Figure BDA0003129741500000126
The structural face disk of (1), as shown in fig. 8.
On the basis of the model shown in FIG. 8, the disc is rotated by an angle theta around the z-axis according to equation (3-3) to obtain an inclination angle theta
Figure BDA0003129741500000127
The disc is translated to the acquired position coordinates (a, b, h) as shown in fig. 9, so that a disc at any position in space is obtained, and therefore, three-dimensional visual modeling of a single fracture is completed.
Figure BDA0003129741500000122
Figure BDA0003129741500000123
For a plurality of fractures, a single disc needs to be established first, and then superposition is carried out, so that a three-dimensional model of the fracture network can be obtained. For some disks which are arranged alternately and have the same shape, position coordinates are calculated according to the same group of intervals, and modeling is carried out according to the steps.
Zhang indicates that data simulation of the fracture network needs to be carried out in a certain space range, so when modeling is carried out by using a computer program, firstly, a simulated space range and data volume are determined by a square box, then, simulated data (attitude, trace length, position or distance) are generated by measured data, then, a circle of a formula (3-1) (figure 7) is cut into n nodes, then, coordinate transformation is carried out on the data of the n nodes according to the formulas (3-2) and (3-3), and finally, a circular disc is spliced by using a triangular curved surface. The previous steps are repeated, and a plurality of discs are drawn, so that a three-dimensional model of the fracture network shown in the figure 10 can be established.
Example 2
This example selects a data (as in table 2) for visualization in a adit (No. PD 4).
TABLE 2 Pilotte PD4 data
Figure BDA0003129741500000124
Figure BDA0003129741500000131
First, the attitude data is converted into spherical coordinates by using the formula (2-1), and the polar diagram shown in fig. 11 is obtained by performing the erythroplanic projection by using the formula (2-3). Then, the best classification number is determined to be 3 by elbow rule (as shown in fig. 12) by using K-Means algorithm for automatic grouping, so as to obtain the center coordinates of the clusters as (0.13,0.01), (-0.76, -0.13), (-0.14, -0.84), and the final classification results are shown in fig. 13, with the discrete degrees of 7.75, 48.35, and 36.63. Dividing each group of data into 8 intervals, and then dividing the Pearson statistic chi of the three groups of data2Respectively as follows: 13.686, 4.158 and 9.333, selecting a significant level alpha which is 0.05, and looking up the table to obtain:
Figure BDA0003129741500000132
therefore, it is
Figure BDA0003129741500000133
All three sets of data obey Fisher's distributions at a significance level α -0.05.
The Fisher model was used to obtain the attitude simulation data shown in FIG. 14, according to the Monte Carlo principle. The average track length of the known data was determined to be 0.485, 1.327, 1.425, and the negative exponential distribution was used to obtain the simulated data of the disc radius. Assuming that the positions of the fractures are subject to uniform distribution, the corresponding data can be obtained from the uniform distribution and are not repeated here.
And according to the step of computer modeling, converting the actual data into a disc so as to realize visualization. In order to facilitate software packaging, a three-dimensional model is constructed by adopting PyQtGraph based on OpenGL and PyQt5 under a Python platform to obtain the model shown in FIG. 15, and only partial discs are provided due to data segmented acquisition. The simulation data was then limited to a 10m × 5m × 5m rectangular parallelepiped, and the model shown in fig. 16 was obtained in the same manner.
The data simulation and visualization are only aimed at one adit, and the data of the actual rock mass is much more complex than that in the model. In a project, a plurality of adits are usually drilled in a region to acquire more detailed fracture data, so that the data of the plurality of adits are required to be simulated to really reflect the development condition of the fracture in the region. The realization idea is consistent with the idea of establishing the disc network model, and the data of a single adit is generated firstly by adopting a superposition method and then superposed to obtain the fracture network model combined by a plurality of adits. Because the data are collected in the same area, the simulation data of a plurality of adits are put together, and the actual situation of the area can be reflected better. Two adits numbered PD4 and PD6 were selected as follows, wherein the PD4 had an azimuth of 80 ° and a depth of 69.5 m; the PD6 has an azimuth angle of 335 deg. and a depth of 94.4 m. Modeling is then performed according to the above grouping re-simulation procedure, and an actual model shown in fig. 17 and a simulation model shown in fig. 18 are obtained.
Because the crack data volume of some footriles is less and cannot meet the requirement of statistics, a Bootstrap method is adopted to solve the problem in the calculation process. The principle is that a certain number of new samples are obtained by resampling in original samples, then statistics needing to be estimated is calculated based on the new samples, and a large amount of data can be obtained by continuous repetition, so that the problem of insufficient data volume is solved.
The data of 14 footrills were studied this time, and as shown in table 3, fracture data in the 14 footrills were integrated and processed to obtain a polar diagram as shown in fig. 19.
TABLE 3 Pionee data
Figure BDA0003129741500000141
On the basis of the polar diagram, basic data for modeling is obtained through data grouping and simulation, and finally, a model as shown in fig. 20 is established by using a disk model. The model integrates data of 14 galleries, in the figure, the combination relationship between the disks can be clearly seen, when the two disks are intersected, the two fractures are communicated, and the structure has the water guiding capacity; when a plurality of disks intersect to form a wedge, it is considered that there is a dangerous rock mass in the area, and the possibility of occurrence of an unfavorable geological disaster such as collapse is high. Therefore, on the basis of the model, through further analysis and research, the defects of the traditional method can be avoided, the crack development condition in one area can be clarified, and the engineering problems of crack connectivity, rock mass structure surface density, underground water migration, slope stability and the like can be analyzed and researched, so that the method has obvious advantages in cost and visualization.
The present invention and its embodiments have been described above schematically, without limitation, and what is shown in the drawings is only one of the embodiments of the present invention, and the actual structure is not limited thereto. Therefore, if the person skilled in the art receives the teaching, without departing from the spirit of the invention, the person skilled in the art shall not inventively design the similar structural modes and embodiments to the technical solution, but shall fall within the scope of the invention.

Claims (8)

1. The method for establishing the three-dimensional visualization model of the fracture network based on the occurrence state is characterized by comprising the following steps of: the method comprises the following steps:
(1) carrying out digital processing on fracture data;
(1.1) converting the attitude data into a spherical coordinate representation;
(1.2) expressing the attitude of the structural plane by using a normal vector; projecting a normal vector of the structural plane onto a red plane by adopting a dimension reduction method and utilizing a red plane projection to obtain a polar point diagram;
(1.3) realizing automatic grouping by utilizing a clustering algorithm and an elbow method, and obtaining a pole center;
(2) simulating fracture data;
simulating the occurrence data by using a Fisher model, obtaining the average trace length of known data, obtaining simulation data of the radius of the disc by using negative exponential distribution, and simulating the position data of the fracture by using uniform distribution;
(3) visual modeling;
obtaining a disc coordinate by using a coordinate transformation method, and constructing a disc by using a triangular surface; the three-dimensional visualization of the crack is carried out by utilizing the disc, the direction of the disc is represented by the occurrence shape, the size of the disc is represented by the trace length, and the spatial position of the position representation disc is measured.
2. The method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 1, wherein: the step (1) comprises the following steps:
(1.1) the occurrence is converted into a spherical coordinate representation with a radius r of 1 according to the following formula:
Figure FDA0003129741490000011
theta is a horizontal angle of the horizontal,
Figure FDA0003129741490000012
is a vertical angle;
(1.2) the attitude of the structural plane is expressed by a normal vector, and the direction cosine of the normal vector can be calculated by the following formula:
Figure FDA0003129741490000013
projecting a normal vector of the structural plane onto a red plane by adopting a dimension reduction method and utilizing a red plane projection to obtain a polar point diagram; and (3) performing the erythro-plano projection by using the following formula to obtain the plane rectangular coordinates of each pole:
Figure FDA0003129741490000021
(1.3) realizing automatic grouping by using a clustering algorithm and an elbow method, obtaining pole centers, regarding the centers obtained by clustering as the average vector direction of the occurrence, and marking the pole center positions obtained by clustering as
Figure FDA0003129741490000022
3. The method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 2, wherein: in the step (2), the fracture data simulation comprises the following specific steps:
according to the central coordinates of each grouping pole
Figure FDA0003129741490000023
Center attitude of the available fractures
Figure FDA0003129741490000024
Taking one set of data, recording the center as
Figure FDA0003129741490000025
Having a direction cosine of
Figure FDA0003129741490000026
Let the group have n data, each data having a direction cosine of (l)j,mj,nj) (ii) a The degree of dispersion can be expressed as κ, with the larger the κ value, the more centered the data is, when κ>At 6 th hour[9]This can be calculated using the following formula:
Figure FDA0003129741490000027
the attitude data obeys a Fisher distribution, and the probability density function thereof is shown as the following formula, wherein the center of the formula is the origin or z-axis direction of the declination projection:
Figure FDA0003129741490000028
the Fisher probability density function is a integrable function, and the integral of the above formula is obtained
Figure FDA0003129741490000029
θ*Cumulative density function of (d):
Figure FDA00031297414900000210
since the cumulative density function follows the uniform distribution of the (0,1) interval, the uniform distribution U-U (0,1) is introduced,
substituting u into the above formula and negating the function to obtain a sampling formula:
Figure FDA00031297414900000211
the sampling formula is utilized to generate a large amount of occurrence data which obey Fisher distribution
Figure FDA00031297414900000212
Can calculate the direction cosine
Figure FDA00031297414900000213
Coordinate transformation is performed by the following formula to obtain required data (l)i,mi,ni)：
Figure FDA0003129741490000031
Then the calculation is carried out reversely to obtain
Figure FDA0003129741490000032
And (alpha)j,βj) (ii) a Therefore, simulation of fracture occurrence data is completed;
the fracture path length L and the disc radius r have the following relationship:
Figure FDA0003129741490000033
from the negative exponential distribution of the fissure vestige, the relationship of the following formula is obtained:
Figure FDA0003129741490000034
integrating and solving a reverse function to obtain a sampling formula of the disc radius:
Figure FDA0003129741490000035
the fracture position data is then simulated with a uniform distribution.
4. The method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 3, wherein:
Figure FDA0003129741490000036
in addition to being available from clustering algorithms, the following can be calculated:
Figure FDA0003129741490000037
5. the method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 4, wherein: the degree of dispersion can also be iterated by dichotomy, resulting in more accurate values from maximum likelihood estimation.
6. The method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 5, wherein: it is necessary to χ' data before generating data using the sampling formula2Fitting and testing, wherein points with large individual errors can be removed in the testing process, and if the tested data obeys Fisher distribution, the data can be generated by using the model; non-compliance requires the use of other models.
7. The method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 6, wherein: for the case of mixing a plurality of Fisher models, the Fisher models are regarded as being mixed by Fisher models with different average occurrence and different dispersion degrees, namely the average occurrence and the dispersion degrees are parameters, so that the Fisher mixed model shown as the following formula can be obtained; the K-Means algorithm is added, and the solution can be carried out continuously by iteration;
Figure FDA0003129741490000041
in addition, the maximum likelihood estimation value can also be solved by utilizing the EM algorithm of the Gaussian mixture model in an iterative mode.
8. The method for building the three-dimensional visualization model of the fracture network based on the attitude as claimed in claim 7, wherein: in the step (3), the specific steps of visual modeling are as follows:
firstly, analyzing a single crack, setting a structural plane as a disc, establishing a disc with the center of circle positioned at an origin and the size of r on a two-dimensional plane, and obtaining a fracture model with horizontal attitude according to the analytical formula:
Figure FDA0003129741490000042
on the basis of the model, the circle is wound around the y-axis by the rotation angle according to the following formula
Figure FDA0003129741490000043
Then obtain the inclination angle of
Figure FDA0003129741490000044
The structural surface disc of (1);
Figure FDA0003129741490000045
the disc is rotated around the z-axis by an angle theta to obtain an inclination angle theta
Figure FDA0003129741490000046
The disk of (2):
Figure FDA0003129741490000047
then translating the disc to the collected position coordinates (a, b, h) to obtain the disc at any spatial position, thereby completing the three-dimensional visual modeling of a single fracture;
after the disks of a plurality of fractures are superposed, a three-dimensional model of the fracture network can be obtained.
CN202110699601.6A 2021-06-23 2021-06-23 Method for establishing fracture network three-dimensional visualization model based on occurrence state Pending CN113468639A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110699601.6A CN113468639A (en) 2021-06-23 2021-06-23 Method for establishing fracture network three-dimensional visualization model based on occurrence state

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110699601.6A CN113468639A (en) 2021-06-23 2021-06-23 Method for establishing fracture network three-dimensional visualization model based on occurrence state

Publications (1)

Publication Number Publication Date
CN113468639A true CN113468639A (en) 2021-10-01

Family

ID=77872505

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110699601.6A Pending CN113468639A (en) 2021-06-23 2021-06-23 Method for establishing fracture network three-dimensional visualization model based on occurrence state

Country Status (1)

Country Link
CN (1) CN113468639A (en)

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2014208886A1 (en) * 2013-06-26 2014-12-31 부경대학교 산학협력단 Three-dimensional imaging method and three-dimensional imaging system of discontinuity network structure in fractured rock
CN107145633A (en) * 2017-04-07 2017-09-08 中国地质大学（武汉） A kind of Forecasting Methodology of the three-dimensional statistical distribution of rock fracture network occurrence

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2014208886A1 (en) * 2013-06-26 2014-12-31 부경대학교 산학협력단 Three-dimensional imaging method and three-dimensional imaging system of discontinuity network structure in fractured rock
CN107145633A (en) * 2017-04-07 2017-09-08 中国地质大学（武汉） A kind of Forecasting Methodology of the three-dimensional statistical distribution of rock fracture network occurrence

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
YUNHUI ZHANG等: "《Nature and evolution of the lithospheric mantle revealed by water contents and He-Ar isotopes of peridotite xenoliths from Changbaishan and Longgang basalts in Northeast China》", 《SCIENCE BULLETIN》 *
占洁伟: "《复杂岩体结构的几何特征精细描述方法研究》", 《中国优秀博硕士学位论文全文数据库(博士)》 *
张亚等: "《基于Monte-Carlo的岩体三维裂隙网络模型研究》", 《甘肃科学学报》 *
李志浩等: "《基于纵向裂隙修正的岩体三维裂隙模拟及连通性分析》", 《隧道建设》 *
郑健: "《岩体结构面几何参数获取的数字化方法及三维网络模型研究》", 《中国优秀博硕士学位论文全文数据库(硕士)》 *
陈剑平，肖树芳，王清编著: "《随机不连续面三维网络计算机模拟原理》", 31 December 1995, 长春：东北师范大学出版社 *

Similar Documents

Publication Publication Date Title
CN104573198A (en) Method for reconstructing digital rock core and pore network model based on random fractal theory
US20140058713A1 (en) Seismic modeling system and method
Li et al. An improved computing method for 3D mechanical connectivity rates based on a polyhedral simulation model of discrete fracture network in rock masses
CN113839725B (en) Method and device for predicting wireless signal propagation
CN114088015A (en) Rapid intelligent generation and sectioning method for rock three-dimensional fracture network model
CN112464513A (en) RQD based on photogrammetry and RQD inversiontOptimal threshold t solving method
CN113468639A (en) Method for establishing fracture network three-dimensional visualization model based on occurrence state
CN112132407A (en) Space RQD based on BQ inversion optimal threshold ttSolving method
CN114088014A (en) Intelligent drawing and sectioning system for three-dimensional fracture network model of rock mass
CN112132411A (en) Based on laser scanning, BQ and RQDtMethod for solving Q anisotropy of anisotropy
Riley et al. Void space and scaffold analysis of packed particles: applications in granular biomaterials
CN109933921B (en) Rolling rock disaster risk assessment method, device and system and storage medium
CN112464514A (en) Based on photogrammetry, RQD and RQDtMethod for solving unfavorable position of roadway excavation
CN112464516A (en) Space RQD based on laser scanning and RQD inversion optimal threshold ttSolving method
CN112132408A (en) Space RQD based on laser scanning and BQ inversion optimal threshold ttSolving method
CN112200431A (en) Dynamic evaluation method for stability of empty zone based on laser scanning, BQ and numerical simulation
CN112132410A (en) Space RQD based on photogrammetry and BQ inversion optimal threshold ttSolving method
CN114088013A (en) Automatic drawing and sectioning system for three-dimensional fracture network model of rock mass
CN115469361B (en) Clastic rock stratum three-dimensional geological modeling method
CN114092628A (en) Intelligent solving and drawing system for generalized RQD of rock mass
Wood Joint spacing in the Caples Lake granodiorite of the Sierra Nevada Batholith in Eldorado National Forest, California: A comparative analysis of joint sets and data resolution
CN113607920B (en) Reconstruction analysis method, experimental device and medium for sedimentary basin by magma bottom wall
CN114088016A (en) Rapid and automatic generation and sectioning method for rock three-dimensional fracture network model
CN114092585A (en) Rapid intelligent solution drawing method for generalized RQD of rock mass
CN114092581A (en) Automatic solution drawing system for rock RQDt

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
RJ01 Rejection of invention patent application after publication
RJ01 Rejection of invention patent application after publication

Application publication date: 20211001