CN115081210A - A method for constrained walk simulation to calculate flow tortuosity in porous media - Google Patents

Abstract

The invention discloses a method for calculating flow tortuosity in a porous medium by constrained walk simulation, which includes simplifying and correcting sample digital cores based on a watershed algorithm; including microscopic-scale continuous phase particles characterized by grayscale values of digital cores in different regions. Collision probability of moving direction; including continuous phase operator migration direction for different digital core resolutions and digital core mirror expansion mapping method; by simulating the stable diffusion coefficient under infinite mirror mapping, combined with the above operators in the real situation transport mode to quantify the tortuosity in different states. Based on the digital core, the present invention introduces the complex pore structure in the migration, compared with the existing methods of ideal mathematical models such as koponen, sousa and barrande, by simulating the tortuosity data obtained by the random walk simulation method in the porous medium. And under the influence of force, the tortuosity obtained by the new method is simulated and calculated, and the coincidence rate of the result is increased by 11% compared with the tortuosity of the traditional method.

Description

A method for constrained walk simulation to calculate flow tortuosity in porous media

technical field

The invention relates to the technical field of gas reservoir development, in particular to the field of unconventional oil and gas development, in particular to a method for calculating flow tortuosity in porous media through restricted walk simulation.

Background technique

As an important parameter to describe the seepage channel, the pore tortuosity is defined as the ratio of the actual length of the specified migration to the macroscopic length of the seepage channel during the seepage process. With the development of unconventional oil and gas reservoirs in recent years, the process of oil and gas seepage has become more and more complicated. The tortuosity obtained by the conventional tortuosity method obtained by simple calculation of pore structure has large errors in explaining the seepage process.

In recent years, experts and scholars tend to use high-precision CT inversion to make digital cores to simulate migration in complex unconventional oil and gas reservoirs, but the interpretation of tortuosity only stays on the quantification of different average particle sizes of minerals in porous media , ignoring the inter-particle force and the anisotropy of the reservoir rock during the real seepage process.

SUMMARY OF THE INVENTION

In view of the above problems, the present invention aims to provide a technique to accurately characterize the tortuosity of complex porous media, so that the results are closer to the real seepage situation.

The technical scheme of the present invention is as follows:

A method for constrained walk simulation to calculate flow tortuosity in porous media, characterized by the following steps:

Step1: Scan the rock sample in 360 degrees, and perform three-dimensional correction processing based on the multi-angle image of the rock sample;

Step 2: The processed three-dimensional grayscale data of the rock sample is segmented by image threshold to separate the pores and the skeleton of the rock sample. Compare the matrix statistical porosity with the actual test porosity adjustment threshold data to obtain the data characterizing the real rock sample skeleton and pores;

Step3: Combined with the field experiment, use the watershed algorithm to divide the Euclidean distance matrix of the void phase, convert the complex matrix into a relatively simple set of corresponding sample space connection points, and set the matrix area as the collision-free area in the random walk.

Step4: By continuously removing the boundary voxels of the rock void facies, combined with the field experiment, the collision and rebound probability of a single voxel of the corresponding rock sample at the corresponding scanning resolution will be obtained correspondingly.

Step5: Select any point in the pores of the sample digital core matrix to perform random walk simulation, record the corresponding number of steps corresponding to the walk distance during the walk, and mirror the matrix boundary during the walk. The random walk simulation coordinates take different measures to limit the number of steps or walk distance in contact with the digital core skeleton.

Step6: Characterize the tortuosity of the sample reservoir through the relationship between the number of random walk steps and the distance.

Further, the specific process of Step1 is:

Use the computer to traverse the reconstructed three-dimensional digital core matrix, and divide the grayscale data into two frequency bands with the amount of grayscale data as the standard. The reconstruction algorithm converts the gray value to the density data and replaces its gray value data.

Further, the specific process of the Step2 is:

Step 201: Determine the porosity of the sample by the liquid saturation drainage method or the helium method and other rock sample porosity test methods.

Step202: Based on the corrected 3D digital core data, divide the skeleton and the pores by setting the grayscale or density threshold, traverse the divided binarized pores, and perform the bwlabel algorithm on the coordinates marked as pores in the matrix to calculate their independence, and mark the pores. And calculate the volume occupied by the disconnected pores after adjusting the threshold. Compared with the porosity φ ₁ obtained in the experiment in Step 201, the porosity φ ₂ and the disconnected pores φ _u are calculated after adjusting the threshold. By continuously adjusting the threshold to make φ ₁ =φ ₂ -φ _u , the data that can reflect the real core pore and framework structure can be determined.

Further, the specific process of the Step3 is:

Perform Gaussian smoothing operation on the gray value in the area divided into pores in the process of Step2 in the digital core. After erasing the minimum value, adjust the grayscale threshold of the pores in the watershed algorithm, and divide the digital core into multiple similar cells according to the grayscale. And according to the size of its gray value, number, number, voxel radius R and its corresponding threshold Th are stored in the corresponding [x, y, z] position matrix of the corresponding digital core.

Based on the Boltzmann equation for the effect of the continuous phase velocity on the collision phase:

Where: f is the dimensionless external force, X _α is the corresponding position, ξ _α is the diffusion velocity of the dimensionless particle

Considering that the calculation of tortuosity in porous media ignores the external force, combined with the Boltzmann equation for the velocity field and the collision phase equation, the diffusion distribution probability of a pore point in the digital core can be equivalent to:

Considering that the pore network in the digital core can be equivalent to an oblong tube, combining the Poiseuille equation in the flow process and the independent pores separated by the watershed algorithm in the Step3 process, the distribution probability of each voxel in the rock sample flow process can be expressed by the following equation:

Simplified, the probabilities of each direction of a single voxel in the core pores at the microscale during free diffusion flow are:

Further, the specific process of the Step5 is:

Step501: Based on the pore network model obtained in Step 2, randomly generate an operator in the pore network, read the data matrix obtained in the process of Step 3, and combine the process of Step 4 to calculate the surrounding N individuals of the operator in the next iteration time step The migration probability of voxels, if the resolution of the digital core is low, that is, the size of a single voxel is large, due to the limitation of accuracy, the model cannot express the pore structure finely, so there may be a pore structure with lower accuracy, considering the surrounding 27 voxels as Migration direction; if the resolution of the digital core is high and can characterize the fine pore structure, then take 14 or 6 directions around a single voxel. The corresponding voxel transport probability is:

When the operator moves to the skeleton position in the digital core matrix during the random walk, the operator position remains unchanged, and the current migration position probability is ignored in the next iterative migration process.

Step502: Obtain the probabilities of different directions obtained in the digital core matrix through the process of Step501, call the random number library in the programming software, and make the operator perform random walks with corresponding direction probabilities in the corresponding pore positions of the digital core matrix; The core matrix performs horizontal and vertical mirror mapping to expand the boundary. The process of mirror mapping corresponding to the boundary is as follows:

Among them: x _i+1 , y _i+1 represent the coordinates in the matrix where the operator before mapping is located; REF(x, y) represent the digital core matrix coordinates after mapping; b is the size of the matrix boundary.

Step503: Repeat the iterative process of Step501-Step502 using the programming software. Simultaneously run a migration simulation in an unrestricted space, that is, all pores in the matrix; design and record the positions of the operators at different time steps according to different needs.

Further, in the described Step503 process, according to different needs, different steps are selected to obtain the specific positions of different operators. The sub-process can be explained as: if the research process is a plane two-dimensional movement, then record and record different iteration steps according to the iterative process. Corresponds to the distance r in two-dimensional coordinates; if the overall tortuosity is considered, the distance r in the corresponding three-dimensional coordinates under different iteration steps is recorded according to the iterative process. If the tortuosity under different migration scales is considered, the corresponding number of steps for the operator to reach the specified migration distance is recorded in combination with the corresponding digital core resolution.

On the basis of ensuring easy implementation, the present invention is more practical in screening the main control factors affecting the production of tight gas than the existing method, and has extremely far-reaching significance for the subsequent research on the production of tight gas.

Description of drawings

In order to explain the embodiments of the present invention or the technical solutions in the prior art more clearly, the following briefly introduces the accompanying drawings that need to be used in the description of the embodiments or the prior art. Obviously, the accompanying drawings in the following description are only These are some embodiments of the present invention, and for those of ordinary skill in the art, other drawings can also be obtained from these drawings without any creative effort.

1 is a schematic flowchart of a method for calculating the flow tortuosity in a porous medium by simulating a restricted walk of the present invention;

Fig. 2 is the digital core matrix after Avizo visualization obtained by embodiment 1 calculation;

Fig. 3 is the matrix data (Avizo visualization) after calculating and processing based on the watershed algorithm that embodiment 1 calculates;

Fig. 4 calculates and obtains the diffusion coefficient under different time steps in embodiment 1;

Detailed ways

The present invention will be further described below in conjunction with the accompanying drawings and embodiments. It should be noted that the embodiments in the present application and the technical features in the embodiments may be combined with each other under the condition of no conflict. It should be noted that, unless otherwise specified, all technical and scientific terms used in this application have the same meaning as commonly understood by one of ordinary skill in the art to which this application belongs. The use of "comprising" or "comprising" and similar words in the present disclosure means that the elements or items appearing before the word encompass the elements or items listed after the word and their equivalents, but do not exclude other elements or items.

Example 1

As shown in FIG. 1, a method for accurately measuring the dynamic tortuosity of seepage flow in porous media based on random walk, including the following steps:

Step1: Perform 360-degree scanning on the rock sample, and perform three-dimensional correction processing based on the multi-angle image of the rock sample; the transformed three-dimensional digital core matrix is shown in Table 1, and the data read visually by Avizo is shown in Figure 2 of the manual.

Table 1 Converted 3D digital core data matrix

Step 2: The processed three-dimensional grayscale data of the rock sample is segmented by image threshold to separate the pores and the skeleton of the rock sample. Compare the matrix statistical porosity with the actual test porosity adjustment threshold data to obtain the data characterizing the real rock sample skeleton and pores;

After adjusting the threshold, the porosity of 0.143 obtained by setting the threshold value 4817 is the closest to the test result, so the porosity is divided by this threshold value. Figure 2 shows the result of the processed matrix data visualized by Avizo.

Step3: Combined with the field experiment, use the watershed algorithm to divide the Euclidean distance matrix of the void phase, convert the complex matrix into a relatively simple set of corresponding sample space connection points, and set the matrix area as the collision-free area in the random walk.

The collision-free area of random walk after the watershed algorithm and correction is shown in the three blue positions in Figure 3.

Step4: By continuously removing the boundary voxels of the rock void facies, combined with the field experiment, the collision and rebound probability of a single voxel of the corresponding rock sample at the corresponding scanning resolution will be obtained correspondingly.

Taking the digital core matrix coordinates [77, 77, 77] as an example, based on its gray value, the migration probabilities in 9 directions are calculated as 163/814, 5/906, 16/103, 49/370, 5/56 , 25/222, 24/299, 50/841, 35/213.

Step5: Select any point in the pores of the sample digital core matrix to perform random walk simulation, record the corresponding number of steps corresponding to the walk distance during the walk, and mirror the matrix boundary during the walk. The random walk simulation coordinates take different measures to limit the number of steps or walk distance in contact with the digital core skeleton.

There are different particle sizes in the digital core samples, and the digital core samples are divided into large particle size, small particle size overall, and 10,000 times*100,000 steps in an unconstrained space, respectively, and the migration distance under different conditions is recorded.

Step6: Characterize the tortuosity of the sample reservoir through the relationship between the number of random walk steps and the distance.

As shown in Figure 4, when the number of walking steps is higher than 70,000 times, the diffusion coefficients in the four different cases tend to be constant. Therefore, the large particle size area, the small particle size area and the overall tortuosity in this sample are respectively , 1.7, 2.1 and 1.9. Bring this value into the Poiseuille equation for flow simulation, the simulation result differs by 8% from the real result, and the simulation error is reduced by 11% compared with the tortuosity value calculated by porosity.

The above are only preferred embodiments of the present invention, and do not limit the present invention in any form. Although the present invention has been disclosed above with preferred embodiments, it is not intended to limit the present invention. Technical personnel, within the scope of the technical solution of the present invention, can make some changes or modifications to equivalent embodiments of equivalent changes by using the technical content disclosed above, but any content that does not depart from the technical solution of the present invention, according to the present invention Any simple modifications, equivalent changes and modifications made to the above embodiments still fall within the scope of the technical solutions of the present invention.

Claims (8)

1. A method of constrained migration simulation to calculate flow tortuosity in a porous medium, comprising the steps of:

step 1: scanning the rock sample by 360 degrees, and performing three-dimensional correction processing based on the multi-angle image of the rock sample;

step 2: separating rock sample pores and rock sample frameworks from the processed three-dimensional gray data of the rock sample through image threshold segmentation, and comparing matrix statistics on the pores and actual test porosity adjustment threshold data to obtain data representing real rock sample frameworks and pores;

step 3: in combination with a field experiment, a watershed algorithm is used for segmenting a space phase Euclidean distance matrix, a complex matrix is converted into a simpler corresponding sample space connection point set, and the matrix area is set as a random walk non-collision area;

step 4: continuously removing rock gap phase boundary voxels, and correspondingly acquiring the collision and rebound probability of a single voxel of a corresponding rock sample under the corresponding scanning resolution by combining a field experiment;

step 5: selecting any point in a pore in a sample digital core matrix, performing random walk simulation, recording a corresponding step number and a corresponding walk distance in the walk process, performing mirror image mapping on a matrix boundary in the walk process, and taking different measures to limit the step number or the walk distance when a random walk simulation coordinate is in contact with a digital core framework;

step 6: and characterizing the reservoir tortuosity of the sample through the relation between the random walk step number and the distance.

2. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation is characterized in that Step1 is specifically processed by the following steps: and traversing the reconstructed three-dimensional digital core matrix by using a computer, dividing the gray data into two frequency bands by taking the gray data amount as a standard, and converting the gray value into density data on the part which is difficult to distinguish based on a reconstruction algorithm of the corresponding gray value and the density for the data after the frequency bands are divided, and replacing the gray value data.

3. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation is characterized in that Step2 is specifically processed by the following steps:

step 201: measuring the porosity of the sample by a rock sample porosity experimental measurement method such as a liquid saturation drainage method or a helium method;

step 202: based on the corrected three-dimensional digital core data, a framework and pores are divided by setting a gray level or density threshold value, divided binary pores are traversed, coordinates marked as pores in a matrix are subjected to a bwlabel algorithm to calculate the independence of the coordinates, the volumes of the unconnected pores are marked and calculated after the threshold value is adjusted, and the porosity phi obtained in the experiment in Step201 is compared ₁ After adjusting the threshold, the porosity phi is counted ₂ And a disconnected pore phi _u By continuously adjusting the threshold value to phi ₁ ＝φ ₂ -φ _u And determining data capable of reflecting the pore and skeleton structure of the real core.

4. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation is characterized in that Step3 is specifically processed by the following steps: and performing Gaussian smoothing operation on gray values in the region divided into pores in the digital core according to Step2, adjusting a pore gray threshold in a watershed algorithm after erasing a minimum value, dividing the region into a plurality of similar cells according to the gray values in the digital core, and numbering according to the gray values, wherein the voxel radius R and the corresponding threshold Th are stored in a corresponding [ x, y, z ] position matrix of the corresponding digital core.

5. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation is characterized in that Step4 is specifically processed by the following steps: based on the effect of continuous phase velocity on the collision phase in the Boltzmann equation:

wherein: f is a dimensionless external force, X _α Is corresponding to the position, xi _α To a dimensionless particle diffusion velocity

Considering the influence of the external force neglected by the tortuosity calculation in the porous medium, and combining the equation of the boltzmann equation for the velocity field in the collision phase, the diffusion distribution probability of a certain pore point in the digital core can be equivalent to:

considering that the pore network in the digital core can be equivalent to a long circular tube, and therefore, in combination with the poisson equation in the flow process and the independent pores separated by the watershed algorithm in the Step3 process, the distribution probability of each voxel in the rock sample flow process can be expressed by the following equation:

simplifying, the probability of each direction of a single voxel in the core pore in the microscopic scale during free diffusion flow is as follows:

wherein: th _N The value is the grey value in the digital rock core; r _i Calculating an equivalent radius for the watershed algorithm; theta.theta. ^* Is the contact angle; _μ is the fluid viscosity.

6. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation, which is characterized in that the Step5 comprises the following specific processes:

step501, randomly generating an operator in a pore network based on the pore network model obtained in the Step2, reading a data matrix obtained in the Step3 process, calculating the migration probability of N surrounding voxels of the random voxel in the next iteration time Step by combining the Step4 process, and if the resolution of the digital core is low, namely the size of a single voxel is large and is limited to the precision limit, the model cannot finely express the pore structure, so that a pore structure lower than the precision possibly exists, and the 27 surrounding voxels are considered as the migration direction; if the resolution of the digital core is high and the fine pore structure can be represented, taking 14 or 6 directions around a single voxel, wherein the corresponding voxel migration probability is as follows:

when the operator is transported to the skeleton position in the digital core matrix in the random migration process, the operator position is kept still, and the current migration position probability is ignored in the next iterative migration process;

step 502: acquiring the probabilities in different directions acquired in the digital core matrix through the Step501 process, and calling a random number library in programming software to enable operators to randomly walk in the positions of the corresponding holes of the digital core matrix with the probabilities in the corresponding directions; performing transverse and longitudinal mirror mapping on the digital core matrix to expand the boundary, wherein the process corresponding to the mirror mapping of the boundary is as follows:

wherein: x is the number of _i+1 ，y _i+1 Representing coordinates in a matrix where an operator before mapping is located; REF (x, y) represents the matrix coordinates of the mapped digital core; b is the matrix boundary size;

step 503: repeatedly iterating the Step501-Step502 process by using programming software, and simultaneously running migration simulation in an unlimited space, namely all pores are in the matrix; and designing and recording the positions of operators under different time step lengths according to different requirements.

7. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation is characterized in that the specific positions of different operators obtained by selecting different Step lengths according to different requirements in the Step503 process can be interpreted as follows: if the research process is planar two-dimensional migration, recording the distance r under the corresponding two-dimensional coordinates under different iteration steps according to the records of the iteration process; and if the integral tortuosity is considered, recording the distance r under the corresponding three-dimensional coordinates under different iteration steps according to the iteration process, and if the tortuosity under different migration scales is considered, combining the corresponding digital core resolution recording operator to reach the corresponding step number of the specified migration distance.

8. A method for calculating the flow tortuosity in a porous medium by using a limited walk simulation is characterized in that Step6 is specifically processed by the following steps: counting the operator average moving distance variance under the corresponding Step length acquired in the Step5 process, calculating the diffusion coefficient in the digital core matrix and the unlimited space, judging the tortuosity through the slope of the image of the average moving distance variance and the Step number, and judging the corresponding tortuosity as follows:

wherein: d _diff Diffusion parameter for free migration, d _lmtd In order to limit the diffusion coefficient of the migration,

corresponding to the image slope.

Images