CN106407678B - One kind being based on phased nonparametric anisotropy variogram construction method - Google Patents

Abstract

The invention discloses a method for constructing a non-parametric anisotropic variogram based on phase control. Different anisotropic variogram fittings are carried out based on different sedimentary facies of the formation, and the construction and solution of the exponential variogram are realized. Through Reasonable facies modeling has the following advantages: Compared with the traditional variogram construction method, the variogram fitted by the present invention is closer to the real geological situation, effectively avoiding the unreasonable construction of the variogram. The error produced in the subsequent random simulation work; the present invention has established a kind of point-to-random selection method based on the inversion data, by randomly selecting the point-to-point in each angle interval as the parameter input when the variogram is constructed, it can be more The anisotropy characteristic of good response variogram; The present invention proposes a kind of variogram parameter fitting method based on ant colony algorithm, by transforming the problem, can improve the accuracy of variogram parameter fitting, for Future work provides technical support.

Description

A Construction Method of Nonparametric Anisotropic Variation Function Based on Phase Control

technical field

The invention relates to the field of geostatistics, in particular to a method for constructing a variation function.

Background technique

The main tool for stochastic modeling is geostatistics. Geostatistics is a statistical method that uses variograms to count the spatial relationship between sampled data, and then uses this statistical rule to simulate the value of random variables at unsampled locations. On the basis of studying the regularity of spatial distribution structure characteristics of regionalized variables, geostatistics selects various appropriate Kriging methods to achieve more accurate estimation or conditional simulation of regionalized variables. It generally includes three basic components: correlation analysis of spatial functions, Kriging estimation and conditional simulation. Spatial function correlation analysis is the analysis of variogram function and covariance function, including their definition and estimated variance, which is the basis of Kriging estimation and conditional simulation.

Since the advent of geostatistics, the regionalization variable theory has been used to characterize the spatial variability of resource information in earth sciences. , bio-environment and engineering technology have had a wide range of impacts, and the research on the distribution characteristics of spatial information or the simulation of its discreteness and volatility has been formed, which can be carried out with geostatistical theory, and geostatistics has also become an important tool for evaluating various regional characteristics. Natural phenomena, natural resources, and new engineering sciences that reproduce their fluctuating processes.

As a basic tool of geostatistics, variogram is a key factor in stochastic modeling. Its main function is to describe the spatial relationship structure of spatial variables, and it can help people analyze and obtain characteristics of reservoirs in modeling. The variogram can fully describe the spatial correlation of reservoir parameters. Under different depositional conditions, the spatial distribution characteristics of reservoirs are different, especially in areas with strong heterogeneity. Anisotropic characteristics of layers. Therefore, the construction of variogram is the core technology of stochastic simulation and stochastic inversion.

Variation function

Regionalized variables are used in geostatistics to represent the object of study. The regionalization variable is a random field Z(x _u ,x _v ,x _w )=Z(x) whose independent variable is the three-dimensional Cartesian coordinates (x _u ,x _v ,x _w ) of the spatial point X, which reflects the geological variable Unique dual features, namely random features and spatial structure features.

The variogram at a point x is defined as half the variance of the difference between regionalized random variables at x with a distance of h, denoted as γ(x,h), that is

r(x,h)=0.5*D ² [Z(x)-Z(x+h)]

＝0.5*E[Z(x)-Z(x+h)] ² -0.5*{E[Z(x)-Z(x+h)]} ² (1-1)

The variogram has three basic parameters: range a, nugget value c ₀ and sill value c.

① variable range

The range a reflects the influence range of the regionalized variable, and there is correlation within the range range, and the correlation gradually weakens with the increase of the distance, and there is no correlation beyond this range.

② nugget value

The nugget value c ₀ = γ(0), theoretically the value of the variogram at the origin should be zero, but due to measurement errors and microscopic variability, random regionalization variables appear within a small distance major changes.

③Abutment value

The sill value reflects the intensity of variation of a regionalized variable within the study area.

The variogram can well describe the characteristics of geological variables. For example, the range reflects the range of influence of the variable, the sill value reflects the variation intensity of the variable, and the slope of the initial section of the variogram near zero reflects the influence of the variable in a smaller range. Within the distance, whether the change of the variable is sharp or gentle, and the variogram in different directions can reflect the anisotropy of the variable and so on ^[1] .

stochastic simulation

The construction of variogram is the premise of stochastic simulation. Stochastic modeling refers to the method of generating optional, equal-probability, and high-precision reservoir models based on the collected information, using random function theory, and adopting stochastic simulation methods. This method admits that the reservoir parameters outside the control points have a certain degree of uncertainty, that is, a certain degree of randomness. Its process is shown in Figure 1.

Contents of the invention

In order to solve the above-mentioned technical problems, the present invention proposes a construction method based on phase-controlled non-parametric anisotropic variogram, combining different sedimentary facies phases in the work area, and constructing non-parametric anisotropic variograms in different sedimentary facies respectively. When combining parameters, a method based on ant colony algorithm is proposed, which can better estimate the corresponding parameters of the variogram. Provide good technical support for further simulation work.

The technical solution adopted in the present invention is: a method for constructing a non-parametric anisotropic variogram based on phase control, comprising:

S1. Sedimentary facies division of the work area, using indicator simulation to divide the sedimentary facies of the work area;

S2. Construct the variogram of different sedimentary facies, specifically including the following sub-steps:

S21. Initialization: import the inverted wave impedance data and logging data;

S22. Phase-controlled modeling: use PETREL software to conduct phase-controlled modeling based on well logging data, and obtain the label phase_label corresponding to each point;

S23. Divide the angle equally, set the number of angle intervals as M equal divisions, mark the equal division M=degree_num, and record the interval as θ, θ=360/M, then the slope corresponding to the nth angle interval is [arctan((n- 1) θ), arctan(n θ));

S24, randomly select point pairs, classify points under different labels, randomly combine any two points into point pairs, and label each group of point pairs;

S25. Calculate the slope of the selected point pair, and put the point pair into the corresponding angle interval according to its value;

S26. Repeat steps S23 to S25 for each sedimentary facies, until there are enough sampling points in each angle interval of each sedimentary facies;

S27. Calculate the variogram respectively according to the sampling data under different sedimentary facies;

S3. Solving the fitting variogram parameters, including:

S31. Encode each variable x _i of the candidate solution {x ₁ , x ₂ ,...} with a binary code string {b _N b _N-1 ...b ₁ b ₀ } with a word length of N, and perform the following formula decoding:

Among them, b∈{0,1}, j=1,2...N, b _N-1 is the highest bit, b ₀ is the lowest bit, the left boundary of the variable x _i is the real value x _imin , and the right boundary is the real value x _imax , z indicates that the left boundary of the decimal integer value corresponding to the binary code string is a real value;

S32. Transform the parameters to be fitted into a form of a directed graph;

S33. Using an exponential model to perform parameter fitting.

Further, before the step S1, it also includes: indicating transformation, specifically: discretizing the continuously determined original data into Boolean quantities 0 or 1 according to the threshold value.

Further, the step S32 is specifically:

Define a directed graph G=(C,L), where the vertex set C is

Among them, v _s is the starting point, the vertex and Respectively represent the state of bit b _j in the binary code string with values of 0 and 1, j=1,2...N, c ₀ (v _s ), Both represent the elements in the vertex set C, e=1,2,...,2N.

Further, the step S33 specifically includes the following sub-steps:

S331, nc=0, initialization of each τ _ij and Δτ _ij , placing m ants on n vertices;

Among them, nc is the number of iteration steps or search times, τ _ij is the amount of information remaining on the line i, j at time t, is the number of pheromones per unit length track left by ant k on the edge arc (i, j);

S332. Put the initial starting point of each ant in the current solution set, and for each ant k (k=1,...,m), according to the probability Move to the next vertex j; put vertex j in the current solution set;

S333. Calculate the objective function value Z _k of each ant, and record the solution that minimizes the objective function value Z _k ; k=1,...m;

S334. Let Path ^* (t) be the best path in the t search period, the objective function value corresponding to the best path is f ^* (t), and the vertex i in the edge arc (i, j) corresponds to the first candidate solution K position, then the pheromone searched by the ant colony is updated according to the following formula

Among them, f ^* (t+1) represents the objective function value corresponding to the best path in the t+1th search period, τ _ij (t,k) represents the pheromone at time t, τ _ij (t+1,k) Represents the pheromone at time t+1, L is the code length of the binary code of the candidate solution, which is a positive integer, and k represents the kth ant;

S335. For each arc edge (i, j), set Δτ _ij =0; nc=nc+1;

S336. If nc<predetermined number of iterations, go to step S332.

Beneficial effects of the present invention: the present invention carries out different anisotropic variogram fittings based on different sedimentary facies of the formation, realizes the construction and solution of the exponential variogram, and has the following advantages through reasonable facies modeling:

(1) The variogram constructed by the present invention is based on the anisotropic variogram of different sedimentary facies. Compared with the traditional variogram construction method, the variogram fitted by the present invention is closer to the real geological situation, effectively avoiding errors in the subsequent random simulation work due to the unreasonable construction of the variogram;

(2) The present invention establishes a method for random selection of point pairs based on inversion data, through the random selection of point pairs in each angle interval as the parameter input when the variogram is constructed, which can better reflect the various parameters of the variogram. Anisotropic characteristics;

(3) The present invention proposes a variogram parameter fitting method based on the ant colony algorithm. By transforming the problem, the accuracy of variogram parameter fitting can be improved, providing technical support for future work.

Description of drawings

Figure 1 is a flowchart of the stochastic simulation.

Fig. 2 is a flow chart of the construction of the variogram provided by the present invention.

Detailed ways

In order to facilitate those skilled in the art to understand the technical content of the present invention, the content of the present invention will be further explained below in conjunction with the accompanying drawings.

As shown in Fig. 1, it is a scheme flowchart of the present invention, and the technical scheme of the present invention is: a kind of construction method based on phase-controlled non-parametric anisotropic variogram, comprising:

S1. Sedimentary facies division in the work area, use indicator simulation to divide the sedimentary facies in the work area; indicator simulation is an effective method to divide the sedimentary facies in the work area. Before performing the simulation, the application needs to carry out instruction transformation, that is, the process of discretizing the continuously determined original data into a Boolean quantity 0 or 1 according to a series of threshold values. Let {Z(x _α ),α=1,…,N} be a set of observation data on the studied reservoir, N means that N observations have been carried out and N different values have been obtained, and a threshold value is given Z, can convert all observed data into indicator value I(x _α ; Z);

Then the probability estimate of the regionalized variable I ^* (x, Z) of the area to be estimated centered on x ₀ for

The weight coefficient λ _α in the formula can be obtained by indicating the Kriging equations, and finally the estimated value Z ^* (X) of the area to be estimated can be obtained;

The weight coefficient λ _b in the formula can be obtained by indicating Kriging equations, and Z(X) represents the regional variation, which is a random variable.

In the PETREL software, you can choose to use the method of indicator simulation to divide the sedimentary facies. After importing the logging data in the work area, the PETREL software will perform clustering according to the logging data. Generally, three different sedimentary facies can meet the needs. After modeling, the phase label corresponding to each point in the work area can be obtained, and the label of each point is stored in the matrix of phase_label, and the label value is 1, 2, 3.

S2. Construct the variograms of different sedimentary facies. The traditional variograms are only constructed using well logging data and cannot fit the variogram curves well. The present invention adopts the data that has been inverted, utilizes relatively abundant data, proposes a method of random selection of point pairs, and constructs the variogram respectively for different sedimentary facies labels, as shown in Figure 2, specifically includes the following Step by step:

S21. Initialization: import the inverted wave impedance data and logging data;

S22. Phase-controlled modeling: use PETREL software to conduct phase-controlled modeling based on well logging data, and obtain the label phase_label corresponding to each point;

S23. Divide the angle equally, set the number of angle intervals as M equal divisions, M=degree_num, and the interval is recorded as θ, θ=360/M, then the slope corresponding to the nth angle interval is [arctan((n-1)· θ), arctan(n θ));

S24. Randomly select point pairs, classify the points under different labels, randomly combine any two points into point pairs, and label each group of point pairs; this application uses the rsenne Twister random number generator to generate one for each group of point pairs random number, the recursive formula of the rsenne Twister random number generator is

Among them, the initial value is (x ₀ ,x ₁ ,...,x _n-1 ), which is an n-dimensional row vector, each x=x ^(w) x ^(w-1) ...x ⁽⁰⁾ It is the binary representation of w bits, and the left side is the high bit; Represents the first wr bit of x _s (u is upper), Represents the last r bits of x _s+1 (l is lower), It is a new w-bit data consisting of the first wr bits of x _s and the last r bits of x _s+1 ; 1≤m<n, r is the character segmentation point, 0≤r≤w-1; "|" is or operator, "⊕" is an XOR operator, ">>" is a right shift, and "<<" is a left shift; A is the convolution transformation matrix of the MT algorithm, and its size is w×w; A representation that means the definition is equal to.

Sort the point pairs according to the size of the random number, and select N groups of point pairs as the input of the variogram calculation.

S25. Calculate the slope of the selected point pair, and put the point pair into the corresponding angle interval according to its value;

S26. Repeat steps S23 to S25 for each sedimentary facies, until the number of sampling points in each angle interval of each sedimentary facies is greater than or equal to N; judge according to the size of the divided angle intervals, if the angle intervals are larger, there will be more Some, if smaller, can be less. In this application, the number N of sampling points is set to be about 20.

S27. Calculate the variogram respectively according to the sampling data under different sedimentary facies;

S3. Solve the parameters of the fitting variogram. Although many fitting methods of the variogram have been proposed, these methods have many inconveniences and shortcomings. For example, the weighted polynomial fitting method largely relies on the experience of geological personnel, and the determination of its weight coefficient is also difficult; compared with the weighted polynomial regression method, the linear programming method only improves the method of calculating polynomial coefficients, without considering weighting. Problem; the goal programming rule is too complicated; the weighted linear programming method takes into account the weighting of the experimental variation function values obtained under different lag distances h, and also ensures the success of the fitting. At the same time, manual intervention is carried out. Proposed legally and linearly, the objective function is too complex to solve. According to the characteristics of these fitting methods, the present invention proposes a method of using an ant colony algorithm to fit the parameters of the variation function. include:

S31. Encode each variable x _i of the candidate solution {x ₁ , x ₂ ,...} with a binary code string {b _N b _N-1 ...b ₁ b ₀ } with a word length of N, and perform the following formula decoding:

Among them, b∈{0,1}, j=1,2...N, b _N-1 is the highest bit, b ₀ is the lowest bit, the left boundary of the variable x _i is the real value x _imin , and the right boundary is the real value x _imax , z indicates that the left boundary of the decimal integer value corresponding to the binary code string is a real value;

S32, convert the parameters to be fitted into the form of a directed graph; define the directed graph G=(C, L), wherein the vertex set C is

v _s is the starting point, the vertex and Respectively represent the state that the bit b _j in the binary code string takes the value of 0 and 1, j=0, 1, 2,..., N, c ₁ is the actual vertex, and the Indicates its generation method, the superscript 0, 1 indicates its two states, the subscript indicates starting from N, then the directed arc L is

That is, for j=2,3,...N, at all vertices and , respectively have and only point to Two directed arcs of ;

S33. Using an exponential model for parameter fitting, the formula is

The parameters to be fitted are a and c.

A∈(40,100) and c∈(8000,16000) can be obtained from the experimental variogram, which are respectively coded and converted into the form of a directed graph G through step S32. then for the problem

Its minimum value problem can be transformed into an ant colony algorithm minimum path finding problem with parameters a and c as a graph, where γ*(h _i ) is the function value corresponding to the lag distance of the experimental variogram h _i . The solution steps are as follows:

S331, nc=0, initialization of each τ _ij and Δτ _ij , placing m ants on n vertices;

Among them, nc is the number of iteration steps or search times, τ _ij is the amount of information remaining on the line i, j at time t, is the number of pheromones per unit length track left by ant k on the arc edge (i, j);

S332. Put the initial starting point of each ant in the current solution set, and for each ant k, k=1,...,m, according to the probability Move to the next vertex j; put vertex j in the current solution set;

S333. Calculate the objective function value Z _k of each ant, and record the solution that minimizes the objective function value Z _k ;

S334. Let Path ^* (t) be the best path in the t search period, the objective function value corresponding to this best path is f ^* (t), and the vertex i in the arc edge (i, j) corresponds to the first candidate solution K position, then the pheromone searched by the ant colony is updated according to the following formula

Among them, f ^* (t+1) represents the objective function value corresponding to the best path in the t+1th search period, τ _ij (t,k) represents the pheromone at time t, τ _ij (t+1,k) Represents the pheromone at time t+1, L is the code length of the binary code of the candidate solution, which is a positive integer, and k represents the kth ant;

S335. For each arc edge (i, j), set Δτ _ij =0; nc=nc+1;

S336. If nc<predetermined number of iterations, go to step S332.

Those skilled in the art will appreciate that the embodiments described here are to help readers understand the principles of the present invention, and it should be understood that the protection scope of the present invention is not limited to such specific statements and embodiments. Various modifications and variations of the present invention will occur to those skilled in the art. Any modifications, equivalent replacements, improvements, etc. made within the spirit and principles of the present invention shall be included within the scope of the claims of the present invention.

Claims (4)

1. A method for constructing a non-parametric anisotropic variogram based on phase control, characterized in that, comprising: S1. Sedimentary facies division of the work area, using indicator simulation to divide the sedimentary facies of the work area; S2. Construct the variogram of different sedimentary facies, specifically including the following sub-steps: S21. Initialization: import the inverted wave impedance data and logging data; S22. Phase-controlled modeling: use PETREL software to perform phase-controlled modeling on the logging data, and obtain the label phase_label corresponding to each point; S23. Divide the angle equally, set the number of angle intervals as M equal divisions, mark the equal division M=degree_num, and record the interval as θ, θ=360/M, then the slope corresponding to the nth angle interval is [arctan((n- 1) θ), arctan(n θ)]; S24, randomly select point pairs, classify points under different labels, randomly combine any two points into point pairs, and label each group of point pairs; S25. Calculate the slope of the selected point pair, and put the point pair into the corresponding angle interval according to its value; S26. Repeat step S23 to step S25 for each sedimentary facies, until there are enough sampling points in each angle interval of each sedimentary facies; S27. Calculate the variogram respectively according to the sampling data under different sedimentary facies; S3. Solving the fitting variogram parameters, including: S31. Encode each variable x _i of the candidate solution {x ₁ , x ₂ ,...} with a binary code string {b _N b _N-1 ... b ₁ } whose word length is N, and decode according to the following formula: Among them, b _j ∈ {0,1}, j=1,2...N, b _N is the highest bit, b ₁ is the lowest bit, the left boundary of the variable x _i is the real value x _imin , and the right boundary is the real value x _imax , z indicates that the left boundary of the decimal integer value corresponding to the binary code string is a real value; S32. Transform the parameters to be fitted into a form of a directed graph; S33. Using an exponential model to perform parameter fitting. 2. A method for constructing a non-parametric anisotropic variogram based on phase control according to claim 1, characterized in that, before the step S1, it also includes: indicating transformation, specifically: according to the threshold value, the continuous The determined raw data is discretized into Boolean 0 or 1. 3. A kind of variogram construction method based on phase-controlled non-parametric anisotropy according to claim 1, is characterized in that, described step S32 is specifically: Define a directed graph G=(C,L), where the vertex set C is Among them, L is a directed arc, v _s is the starting point, and the vertex and Respectively represent the state of bit b _j in the binary code string with values of 0 and 1, j=1,2...N, c ₀ (v _s ), Both represent the elements in the vertex set C, e=1,2,...,2N. 4. A kind of variogram construction method based on phase-controlled nonparametric anisotropy according to claim 1, is characterized in that, described step S33 specifically comprises the following sub-steps: S331, nc=0, each τ _ij and The initialization of , put m ants on n vertices; Among them, nc is the number of iteration steps or search times, τ _ij is the amount of information remaining on the line i, j at time t, is the number of pheromones per unit length track left by ant k on the edge arc (i, j); S332. Put the initial starting point of each ant in the current solution set, and for each ant k, k=1,...,m, according to the probability Move to the next vertex j; put vertex j in the current solution set; S333. Calculate the objective function value Z _k of each ant, and record the solution that minimizes the objective function value Z _k ; k=1,...m; S334. Let Path ^* (t) be the best path in the t search period, the objective function value corresponding to the best path is f ^* (t), and the vertex i in the edge arc (i, j) corresponds to the first candidate solution K bit, the pheromone searched by the ant colony is updated according to the following formula Among them, f ^* (t+1) represents the objective function value corresponding to the best path in the t+1th search period, τ _ij (t,k) represents the pheromone at time t, τ _ij (t+1,k) Represents the pheromone at time t+1, L is the code length of the binary code of the candidate solution, which is a positive integer, and k represents the kth ant; S335. For each arc edge (i, j), set nc=nc+1; S336. If nc<predetermined number of iterations, go to step S332.