Methods for estimating fluid saturation of rocks

Abstract

The present invention provides a method for estimating fluid saturation of hydrocarbon-bearing rock based on an image of the rock. The image is segmented to represent pore space or solid material in the rock. Image pore volume is estimated based on the segmented image, and a corrected pore volume is determined to account for sub-resolution pore volume missing in the image of the rock. Direct flow simulation of the rock image is used to estimate the image-derived wetting fluid saturation of the rock and corrected for the corrected pore volume. A trained model that supports backpropagation can be used to segment the image. A method that supports backpropagation can be used to estimate the fluid saturation using an image selected from a series of 2D projection images, 3D reconstructed images, and combinations thereof.

Description

Technical Field

[0001] The present invention relates to a method for estimating fluid saturation of a rock, and in particular, to a method for estimating fluid saturation of a rock based on a 3D image of the rock. Background Art

[0002] Accurately determining the fluid saturation within a hydrocarbon-bearing reservoir is an important factor in determining whether to select a hydrocarbon-bearing reservoir to be developed and in developing and managing a hydrocarbon-bearing reservoir. Based on an understanding of the pore-scale fluid displacement dynamics, pore-scale flow simulation can facilitate the calculation of Darcy-scale flow parameters, such as used as input in reservoir simulation. However, many factors affect fluid flow, including, for example, fluid density and viscosity, interfacial tension, fluid flow rate, surface wettability, and pore geometry.

[0003] However, laboratory testing for these factors is time-consuming and costly. In addition, the number of samples that can be processed is relatively limited due to the time and expense required to perform each test. Information needs to be provided more quickly so that more timely decisions can be made. Many in the industry are turning to simulation of fluid flow in reservoirs to make informed decisions in a timely manner.

[0004] Digital rock physics is a technology that has been developed to provide faster, more extensive and cheaper analysis of hydrocarbon-bearing formation rocks to determine key petrophysical characteristics of the rocks. Digital rock physics uses digital images of formation rocks to simulate rock multiphysics at the pore scale and predict the properties of complex rocks.

[0005] Direct numerical simulations can therefore be solved directly on complex pore spaces, such as those extracted from micrometer-scale X-ray computed tomography images. Alpak et al. (“Prediction of fluid topology and relative permeability in imbibition in sandstone rock by direct numerical simulation”) Advances in Water Resources122:49-59; 2018 and “Direct simulation of pore-scale two-phase visco-capillary flow on large digital rock images using a phase-field lattice Boltzmann method on general-purpose graphics processing units” Computational Geosciences 23:849-880; 2019) described an energy-based LBM (eLBM), which is a two-phase flow simulation system consisting of three modules: a forced drainage simulation module, a forced imbibition simulation module, and a steady-state relative permeability calculation module.

[0006] Leon Carrera et al. (US2018 / 0321127A1) discloses a method for providing a numerical model of a rock sample that, when used for flow simulations, reproduces porosity and permeability based on measurements obtained in the rock sample. Data recovered from a CT scan provides a statistical density function of attenuated X-ray radiation in a volume of the rock sample. Porosity is filled into a subvolume of a 3D model, where the porosity is spatially distributed between cells of the subvolume by a Gaussian simulation algorithm that responds to a porosity statistical distribution function that is an approximation of a Gaussian density function. Permeability is filled between cells as defined by a scalar function that responds to the porosity of the cells. The overall permeability is obtained by performing a numerical simulation of the 3D model. These steps are repeated until the difference between the overall permeability measured on a plug sample of the rock and the overall permeability calculated from the 3D model is less than a threshold.

[0007] Fredrich et al. (US9,070,049B2 and WO2014 / 142976A1) are directed to systems and methods for improving direct numerical simulation of material properties from rock samples. High-resolution images are used to resolve the pore space thereof. Fredrich et al. recognize that pore system heterogeneity may not always be well represented within a small imaged portion of the rock, resulting in errors in the calculated material properties due to lack of pore system representation. Therefore, Fredrich et al. identify a characterization unit volume ("REV") determined by selecting a test volume and an adjacent volume and calculating the difference in rock physical property values ​​of the adjacent volumes.

[0008] Digital rock physics modeling conventionally relies on the assumption that the pore volume of hydrocarbon-bearing rocks can be accurately determined from micron-scale images. However, the resolution of the images is limited, and therefore a significant portion of the rock's pore volume may remain unresolved.

[0009] For example, as in Saxena et al., 2019 (“Estimating Pore Volume of Rocks from Pore-Scale Imaging” Transport in Porous Media 129:403-412; 2019), due to the resolution limitations of micron-scale images currently generated by micro-CT detectors, up to half of the total pore volume may be missing in micron-scale X-ray computed tomography images of reservoir rocks. Direct numerical simulations cannot obtain the missing pore volume, which limits the applicability of digital rock to infer the true residual saturation of reservoir fluids. In order to derive meaningful results from direct simulations, at least the original fluid saturation inferred from the micron-scale image simulations must be corrected for the missing pore volume.

[0010] Direct simulations of fluid flow need to be correctly parameterized to include the effects of sub-resolution pore volumes in order to correct the inferred fluid saturations from direct simulations of fluid flow. Summary of the invention

[0011] According to one aspect of the present invention, a method for estimating the fluid saturation of hydrocarbon-bearing rocks based on rock images is provided, which includes: obtaining a 3D image of a rock from a hydrocarbon-bearing formation in an oil field, wherein the 3D image is composed of multiple voxels and the 3D image has a resolution; processing the 3D image by selecting each voxel of the 3D image to segment the 3D image to represent pore space in the rock or solid matter in the rock; estimating the image pore volume of the rock based on the segmented 3D image, wherein the image pore volume lacks the sub-resolution pore volume of the rock; determining a corrected pore volume suitable for taking into account the sub-resolution pore volume of the rock; using direct flow simulation of the 3D image to estimate the image-derived wetting fluid saturation of the rock; and determining the corrected wetting fluid saturation of the rock based on the image-derived wetting fluid saturation of the rock, the image pore volume, and the corrected pore volume.

[0012] According to another aspect of the present invention, a back-propagation-enabled method for estimating fluid saturation of hydrocarbon-bearing rock is provided, which includes the following steps: obtaining a 3D image of a rock from a hydrocarbon-bearing formation in an oil field, the 3D image having a resolution; applying a trained model that supports back-propagation to segment the 3D image; estimating an image pore volume of the rock based on the segmented 3D image, the image pore volume lacking a sub-resolution pore volume of the rock; determining a corrected pore volume suitable for accounting for the sub-resolution pore volume of the rock; estimating an image-derived wetting fluid saturation of the rock using direct flow simulation of the 3D image; and determining the corrected wetting fluid saturation of the rock based on the image-derived wetting fluid saturation of the rock, the image pore volume, and the corrected pore volume.

[0013] According to another aspect of the present invention, a back-propagation-enabled method for estimating the fluid saturation of a rock based on an image of the rock is provided, comprising the following steps: obtaining an image of a rock from a hydrocarbon-bearing formation in an oil field, the image being selected from a group consisting of a series of 2D projection images, 3D reconstructed images, and a combination thereof; and applying a back-propagation-enabled trained model to obtain the fluid saturation of the rock. BRIEF DESCRIPTION OF THE DRAWINGS

[0014] The present invention will be better understood by reference to the following detailed description of the preferred embodiments and the accompanying drawings referred to therein, in which:

[0015] The figure is a graphical representation of relative permeability in the oil and water phases versus fluid saturation, which illustrates the imbibition of the test rock volume. DETAILED DESCRIPTION

[0016] According to the method of the present invention, the fluid saturation of the rock can be estimated from the direct simulation of the fluid flow in the 3D image faster than making laboratory measurements. In addition, by taking into account the missing pore volume in the segmented 3D image, the directional simulation of the fluid flow is more accurate.

[0017] Specifically, the method of the present invention estimates fluid saturation based on the image pore volume of the rock, the corrected pore volume, and the wetting fluid saturation of the rock derived from direct simulation of the 3D rock image. More specifically, a meaningful estimate of fluid saturation is determined by correcting the image-derived wetting fluid saturation determined from the micron-scale image. In a preferred embodiment, the image-derived wetting fluid saturation is corrected for sub-resolution pore volume missing in the lower resolution image.

[0018] The inventors have unexpectedly discovered that capillary physics in rocks can be used to quantify the effect of image resolution on pore volume. This discovery enables digital rock or digital rock physics to provide the ability to determine wetting fluid saturation from 3D images of rocks generated by micro-CT techniques. Advantageously, the corrected wetting fluid saturation, and therefore fluid saturation, compensates for limited image resolution without the need for higher resolution imaging that is only possible at the expense of image field of view or physical laboratory measurements.

[0019] The present invention provides a method for more accurately estimating the fluid saturation of a rock based on a raw 3D pore-scale image of the rock that has a limited resolution relative to the actual pore structure of the rock. Contrary to current assumptions in digital rock physics modeling, the present inventors recognize that a significant portion of the pore volume of hydrocarbon-bearing rocks is contained in pores of a size below the image resolution provided by the 3D pore-scale imaging techniques typically used to provide images of such rocks. As a result, conventional digital rock physics modeling substantially underestimates the effective wetting fluid saturation of the rock, and therefore the fluid saturation prediction, by failing to account for pores smaller than the image resolution of the pore-scale imaging techniques.

[0020] The present invention also provides a back-propagation enabled method for estimating fluid saturation of a rock from a 3D image of the rock. A back-propagation enabled trained model is applied to the 3D image to segment the 3D image of the rock.

[0021] In a preferred embodiment, a trained model is generated by providing a training set of images of rock, segmenting the images into a plurality of labeled voxels representing pore space and solid material in the rock, and training the model using the labeled voxels via back propagation.

[0022] The image training set of the rock may include, for example, 2D projection images obtained from pore-scale imaging techniques, 3D images reconstructed from 2D projection images, synthetic 2D images, synthetic 3D images, and combinations thereof. In a preferred embodiment, the training set of images is obtained from a cloud-based tool adapted to store 2D projection images from pore space imaging techniques, in particular from micro-CT and thin sections. The tool is adapted to process the 2D projection images to produce reconstructed 3D images. The tool is also adapted to store the resulting 3D images.

[0023] Examples of processes that support back propagation include, but are not limited to, artificial intelligence, machine learning, and deep learning. Those skilled in the art will appreciate that advances in processes that support back propagation continue rapidly. Even under different names, the methods of the present invention are expected to be applicable to those advances. Therefore, even if not explicitly cited herein, the methods of the present invention are also applicable to other advances in processes that support back propagation.

[0024] A preferred embodiment of a process supporting back-propagation is a deep learning process, including but not limited to a convolutional neural network.

[0025] The process supporting back-propagation can be supervised, semi-supervised, unsupervised, or a combination thereof. In one embodiment, a supervised method is made semi-supervised by adding unsupervised techniques.

[0026] In a supervised back-propagation process, a training set of images is labeled to provide examples of pore space and solid matter of interest. In an unsupervised back-propagation process, pore space and / or solid matter of interest can be identified by, for example, drawing a polygon around the image of interest in the image. The trained process will then identify regions of interest with similar latent space characteristics. When the training set is labeled as images, the dimensions of the labels can be 1D-3D.

[0027] In one embodiment, the supervised back-propagation-enabled process is a classification process. The classification process can be performed voxel-wise, slice-wise, and / or volume-wise.

[0028] In another embodiment, the unsupervised back-propagation-enabled process is a clustering process. The clustering process can be performed voxel-wise, slice-wise, and / or volume-wise.

[0029] In another embodiment, the unsupervised back-propagation-enabled process is a generative process. The generative process can be performed voxel-wise, slice-wise, and / or volume-wise.

[0030] Preferably, the process supporting back-propagation is a segmentation process.

[0031] In a preferred embodiment, the training step includes validation and testing.

[0032] In the method of the present invention, the petrophysical characteristics of the rock, in particular the fluid saturation of the rock, can be estimated from the image of the rock. The rock image is obtained from a rock from a hydrocarbon-bearing formation, and the petrophysical characteristics of the formation or a portion thereof are of interest. Preferably, the rock can be a sandstone, a carbonate rock, a shale, and a combination thereof from a hydrocarbon-bearing formation. The rock can be obtained by conventional means for obtaining rock samples from hydrocarbon formations. In a preferred embodiment, a core sample of the rock is obtained by coring a portion of the formation from a well in the formation. Alternatively, the rock sample can be obtained from drill cuttings produced in a well being drilled in the formation. The rock can be obtained from the same borehole as the resistivity log. Alternatively, the rock can be obtained from another borehole in the same oil field as the borehole in which the resistivity log was produced.

[0033] The rock sample should be of sufficient size to obtain a 3D image with sufficient volume at the scale at which the image is to be generated. Specifically, the rock sample should be of sufficient size so that features of the bulk of the sample dominate over features of the edges of the sample at the scale or field of view of the image to be generated.

[0034] A 3D image consisting of a plurality of voxels is obtained from a rock sample. A 3D image of a rock may be obtained using pore-scale imaging techniques. A 3D image of a rock may be obtained by X-ray computed tomography, including but not limited to X-ray micro-computed tomography (micro-CT) and X-ray nano-computed tomography (nano-CT), acoustic microscopy, or magnetic resonance imaging. Most preferably, a 3D image of a rock is obtained by micro-CT to provide a sufficient field of view of the rock to avoid edge pores distorting the overall pore volume of the resulting image, as well as to reduce the scanning time and computational requirements required for higher resolution tomography (e.g., nano-CT).

[0035] In a preferred embodiment, the 3D image is obtained from a cloud-based tool adapted to store 2D projection images from pore space imaging techniques, in particular from micro-CT and thin sections. The tool is adapted to process the 2D projection images to produce a reconstructed 3D image. The tool is also adapted to store the resulting 3D image.

[0036] A 3D image of a rock obtained by pore-scale imaging techniques has a resolution. The voxels of a 3D image define the resolution of the image. An image is composed of a plurality of voxels, wherein the volume defined by each voxel represents the maximum resolution of the image. The resolution of the image should be selected to provide a voxel size at which the major pore throats for fluid flow in the rock are adequately resolved and to provide a sufficient field of view to represent the entire rock to be analyzed for a given rock physical property, such as pore volume and fluid saturation. For the purposes of this paper, the major pore throat size (D d ) is the non-wetting liquid pressure at the pore entrance (P d ), where the pore entrance pressure is the minimum pressure required before a non-wetting liquid can begin to invade the pore structure of the rock.

[0037] The resolution of the micro-CT images may be selected based on the size of the rock sample, the relative average pore size of the rock type, the time required for imaging, and the computational power required to store the image data and perform additional computational activities on the image data. The image resolution should be selected to be detailed enough so that a non-wetting liquid capillary injection curve can be plotted based on a segmented image generated from the image, while maintaining a sufficient field of view to avoid edge pores distorting the overall pore volume of the resulting image. In a preferred embodiment, the image resolution is selected to require little computational power to store the image and perform additional computational activities on the image, while providing sufficient detail to construct a capillary injection curve based on the segmented image. The image resolution may be selected based on the type of rock, with sandstones typically having a larger pore structure than carbonates and requiring a lower image resolution than carbonates, and carbonates having a larger pore structure than shale and requiring a lower image resolution than shale. Micro-CT image resolution may be as low as 0.1 μm per voxel. 3 Up to 30μm 3 For sandstone, micro-CT images are preferably acquired at 1 μm per voxel. 3 Up to 25 μm 3 or 2.5 μm per voxel 3 Up to 15μm 3 For carbonate rocks, the resolution of micro-CT images is preferably 0.5 μm. 3 Up to 20μm 3 or 1μm 3 Up to 10μm 3 and for shale, the resolution of micro-CT (or nano-CT) images is preferably within 0.1 μm 3 Up to 10μm 3 or 0.5μm 3 Up to 5μm 3 within the range.

[0038] In a preferred embodiment, the acquired image may be processed to reduce noise and image artifacts. Noise may be filtered out of the acquired image by filtering using a local mean filter to reduce noise. Imaging artifacts primarily at the outer edges of the acquired image may be reduced by processing the image while excluding the outer edges of the image.

[0039] Processing a 3D image obtained for a rock to segment voxels of the image into voxels representing pore space in the rock or solid matter in the rock, thereby producing a binary image in which pore voxels have a value of 0 and solid matter voxels have a value of 1 (or vice versa). The image can be a grayscale image, and processing the voxels of the image to segment the image into voxels representing pore space or solid matter can be achieved by assigning the voxels a name of pore space or solid matter based on a threshold, wherein voxels with image intensity above the threshold can be assigned a value representing pores (or solid matter) and voxels with image intensity below the threshold can be assigned a value representing solid matter (or pores). The Otsu (“A Threshold Selection Method from Gray-level Histogram”) method can be used. IEEE Trans.SMC The threshold is calculated using the Otsu method described in 9:62-66; 1979), or other threshold calculation algorithms known in the art.

[0040] The 3D image of the rock can be processed using a segmentation algorithm known in the art to segment the voxels into pore space voxels and solid material voxels. Preferably, the segmentation method is selected to distinguish between conductive pores and non-conductive minerals. Examples of segmentation methods are described in Otsu (1979), Andra et al. ("Digital Rock Physics Benchmarks-Part II: Computing Effective Properties"), "count Computer and Earth Sciences ( Computers and Geosciences 50:33-43;2013), Saxena et al. (“Effect of Image Segmentation & Voxel Size on Micro-CT Computed Effective Transport & Elastic Properties” Ship Marine and Petroleum Geology 86:972-990; 2019), and Chuang et al. (“Fuzzy C-Means Clustering with SpatialInformation for Image Segmentation,” Computerized Medical Imaging Graphics Comput.Med.Imaging Graph.)》30:9-15; 2006). Those skilled in the art will understand the desired segmentation selection. Preferably, segmentation using a segmentation algorithm is performed automatically using a data processing system.

[0041] After the image has been segmented, the image-derived pore volume, φ, is estimated from the segmented 3D image of the rock. I . The image pore volume of the rock can be estimated by summing the number of voxels representing pore space in the segmented image, summing the total number of voxels in the segmented image (or obtaining the total number of voxels from the imaging parameters), and then dividing the sum of the number of voxels representing pore space in the segmented image by the total number of voxels in the segmented image. The sum of the number of voxels representing pore space in the segmented image can be determined by adding the number of voxels assigned a binary value (e.g., 1 or 0) representing pore space. The sum of the total number of voxels in the segmented image can be determined by adding the total number of voxels assigned a binary value (pore space voxels and solid matter voxels). The image pore volume lacks the sub-resolution pore volume of the rock.

[0042] According to the present invention, a corrected pore volume suitable for taking into account the sub-resolution pore volume of the rock is determined. In one embodiment of the present invention, the corrected pore volume is measured using conventional laboratory measurements. In another embodiment of the present invention, the corrected pore volume is determined by applying a correction factor to the image pore volume, for example by using the image-derived pore volume (φ I ) and image-derived correction factors (α R ) to estimate the actual pore volume (φ ∞ ) to determine:

[0043] φ ∞ =φ1 / α R (1)

[0044] The correction factor can be estimated by a transformation. In a preferred embodiment, the correction factor is determined by determining a non-wetting liquid capillary pressure curve from a segmented 3D image of the rock for pores identifiable in the segmented image at pressures up to the image limit pressure. Preferably, mercury or Wood's metal is selected as the non-wetting liquid. The non-wetting liquid capillary pressure curve can be determined from the segmented image by plotting the porosity of the rock occupied by the non-wetting liquid at selected pressures up to the image limit pressure based on a simulation of the non-wetting liquid filling the pore space of the image.

[0045] Correction factor (α R ) depends on the choice of mercury intrusion capillary pressure (MICP) model. For the Thomeer model of the capillary pressure curve, α R Given by:

[0046]

[0047] where G is the pore geometry factor that captures the shape of the MICP curve, and N is the pore throat resolution parameter.

[0048] The pore resolution parameter N is determined based on the non-wetting liquid capillary pressure curve derived from the image and the image resolution determined based on the size of the voxel. The pore resolution parameter Dd is the ratio of the pore throat size (Dd) into which the non-wetting liquid enters at the inlet pressure (Pd) to the voxel size (Δx), N = (Dd / Δx). The size of the voxel can be determined based on the parameters of the 3D imaging, i.e., the resolution of the image. The pore throat size (Dd) of the pore into which the non-wetting liquid enters at the inlet pressure (Pd) can be determined based on the non-wetting liquid capillary pressure curve, which is derived from the segmented 3D image of the rock according to the following equation:

[0049] Dd=4σcosθ / Pd (3)

[0050] The pore geometry factor G is determined from the non-wetting liquid capillary pressure curve derived from the image. The pore geometry factor G can be determined by plotting a best fit curve onto the non-wetting liquid capillary pressure curve simulated from the segmented image and determining the pore geometry factor based on the shape of the curve. The best fit curve can be plotted by the least squares method or by any conventional curve fitting method.

[0051] The parameter G, which increases with pore volume complexity, is exactly 0 for pipes, about 0.1 to 0.3 for sandstones, and about 0.3 to 0.5 for carbonates. The parameters N, G, and α can also be estimated directly from the limited-resolution 3D images by analyzing the response of numerically simulated mercury injection. R , without the need for laboratory measurements.

[0052] In order to derive meaningful results from the direct simulations, and in accordance with the present invention, the original fluid saturation inferred from the microscale image simulations is corrected for the missing pore volume.

[0053] The following expression relates the fraction of the rock pore volume occupied by mercury to the number of voxels (ξ) available to resolve the last pore throat penetrated during the MICP simulation at a given capillary pressure P:

[0054] φ(ξ)=α(ξ)φ ∞ (4)

[0055] in

[0056]

[0057] In equation (5), Δx is the image voxel size, D(P) is the size of the pore throat penetrated at a given pressure P in the MICP simulation, σ is the mercury-air surface tension, and θ is the contact angle. If all pore throats are invaded, then ξ = 1, and therefore the pore volume entered by mercury will be equal to the pore volume of the image, i.e., φ(1) = φ I and α(1) = α R The expression for α(ξ) depends on the type of MICP model that best describes the behavior of the simulation curve. For Thomeer's MICP, the functional form α(ξ) is given by:

[0058]

[0059] The corresponding expression for the Brooks-Corey fit is as follows:

[0060]

[0061] where λ is a fitting parameter similar to parameter G.

[0062] Using equation (4), the air saturation in the MICP simulation is related to the ξ value:

[0063]

[0064] Where S wI and S nwI are the image-derived saturations of the wetting phase (air; denoted by subscript E) and the non-wetting phase (mercury; denoted by subscript FE), respectively. Symbol V nwI represents the volume fraction occupied by mercury in the image. Fluid saturation is relative to the actual rock pore volume (φ ∞ ) is given by:

[0065]

[0066] Equations (8) and (9) assume that the mercury volumes estimated from direct simulations are accurate for pores and pore throats corresponding to sufficiently large values ​​of ξ. Using these equations, the image-derived wetting fluid saturation is transformed to the corrected true wetting fluid saturation:

[0067] S w∞ (ξ)=1-(1-S wI (ξ))α(1)(drainage) (10)

[0068] Equations (8) to (10) assume that all the unresolved porosity is filled with a wetting phase (corresponding to either a water-wet or oil-wet condition). For drainage processes, we can assume that the rock is essentially water-wet, so the previous assumption is valid and equation (10) can be applied to drainage processes. However, for imbibition processes, mixed wettability can be formed in the rock due to the presence of hydrocarbons. Under this condition, it cannot be assumed that the unresolved porosity is filled with only one type of fluid. It is not clear what percentage of fluid is present in the unresolved pores. As an initial estimate, the average value between α(1) and 1 can be used:

[0069] S w∞ (ξ)=1-(1-S wI (ξ))α imb (Imbibition) (11)

[0070] where α imb is the average value between α(1) and 1.

[0071] Equations (10) and (11) represent the inventors' recognition that image-derived properties, especially for lower resolution images, do not match those measured in the laboratory. The inventors recognized that without correcting for the missing sub-resolution pore volume, due to φ I <φ ∞ , thus underestimating the fluid flow predictions derived from direct flow simulations using segmented micro-CT images.

[0072] If the resolution of the sandstone reservoir (e.g., G = 0.2) is moderate (e.g., N = 5) and the air saturation inferred from the image is low (e.g., S wI =0.4), then α(1)=0.75, and therefore the corrected air saturation will be S w∞ = 0.55. This shows that correction of the image-derived fluid saturation is crucial before comparing any simulation results with laboratory measurements. In other words, the saturation S w∞ The saturation is comparable to the lab measured saturation, while the original image saturation S wI Not comparable to the saturation measured in the laboratory.

[0073] The parameters (N, G, and λ) can be directly obtained from Saxena et al. (“Estimating Pore Volume of Rocks from Pore-Scale Imaging”) Transport in Porous MediaComparing raw MICP simulation results (fluid volumes derived from resolved porosity-normalized images in 3D images of rock) with those from laboratory measurements can lead to false negative comparisons.

[0074] From an imaging perspective, pore throats corresponding to larger values ​​of ξ are better resolved and can therefore also be segmented with higher confidence. For any direct simulation, the fluid saturations corresponding to larger values ​​of ξ will have higher fidelity than those estimated for lower values ​​of ξ. If ξ = 2, then only 2 voxels (in 1D) are available for non-wetting fluids to penetrate. Therefore, the associated fluid volume within such low-resolution pores remains uncertain, and therefore the associated fluid distribution within these narrow pore throats cannot be accurately captured by direct simulation. To ensure that ξ ≥ 2, the following inequality needs to be satisfied:

[0075]

[0076] as well as

[0077] S w∞ ≥1-α(2) or S nw∞ ≤α(2) (13)

[0078] The inequalities in equations (12) and (13) follow directly from the results in equations (8) and (9). Before comparing the image simulation-derived wetting / non-wetting fluid (e.g., air / mercury) saturations with laboratory-measured saturations, another aspect that needs to be considered is the pressure at which the non-wetting fluid begins to enter the digital rock. In the numerical simulation, the non-wetting fluid first enters the digital rock sample via the "external" pores rather than the pore throats, because the pores can enter at much lower pressures than the pore throats - compensating for the pressure and volume biases introduced by the external pores to the flow characteristics is referred to here as closure correction. Although this phenomenon also occurs in physical MICP measurements, the volume contribution of the external pores is negligible compared to the total volume of non-wetting fluid injected into the physical sample, simply because the sample size is significantly larger. Therefore, in the initial drainage simulation, the pressure and flow data are only meaningful after the external pores are filled with non-wetting fluid (e.g., mercury). This condition is guaranteed if the following equation is satisfied:

[0079]

[0080] as well as

[0081]

[0082] The saturation cutoff in equations (12) and (13) is related to the image voxel size, i.e., the finer the voxel size, the larger the parameter N, and therefore the less restrictive the cutoff value. At the same time, the saturation cutoff in equations (14) and (15) is related to the field of view, where the larger the field of view, the smaller the value V closure The smaller the value, the less restrictive the cutoff value.

[0083] In another embodiment, the saturation cutoffs associated with image resolution in equations (12) and (13) can be relaxed by re-imaging the rock sample at a finer voxel size. However, this will further limit the saturation cutoff associated with the field of view, since state-of-the-art micro-CT detectors can only capture a limited number of voxels (typically less than 5000x5000x5000 voxels) and can only image at a finer voxel size (thus increasing the value of parameter N) at the expense of field of view (thus producing larger values ​​of V closure ).

[0084] Multiphase flow is a complex phenomenon in which viscous and capillary forces interact on a micrometer scale. The degree of resolution of the pore space further complicates the numerical simulation of this interaction between these two forces and can therefore negatively impact the interpretation of the simulation results. Before comparing the results of direct numerical simulations with laboratory measured data, any limitations due to the missing pore volume must be taken into account. According to the method of the present invention, the fluid saturation inferred using direct flow simulations is corrected for the sub-resolution pore volume to recover the true saturation.

[0085] In one embodiment of the invention, the method of the invention is used to correct the brine saturation inferred from direct numerical simulation of multiphase flow using segmented micro-CT images.

[0086] For example, the simulation can be limited to the initial discharge of a non-wetting oil phase into a completely water-wet rock. It should be noted that the results in equations (10) to (14) are independent of surface tension and wetting angle, and are therefore applicable to initial discharge of other fluid systems (e.g., water-oil). However, the correction assumes that the non-wetting fluid (oil) cannot enter any sub-resolution pores, and therefore does not apply if the simulation is performed only on the same field of view imaged at infinite resolution, so that the numerical simulator can enter all pores, and if the non-wetting phase otherwise invades the sub-resolution pore volume.

[0087] A 3D micro-CT image is obtained having an image and a voxel size. The 3D image is processed to segment the image by selecting each voxel of the 3D image to represent pore space in the rock or solid material in the rock.

[0088] Image porosity is determined from the segmented three-dimensional image of the rock by summing the number of voxels representing pore space in the segmented image and dividing by the total number of voxels in the image determined from the imaging parameters.

[0089] The mercury capillary injection curve is then determined from the segmented 3D images of the rock. An image-based curve is created by plotting the volume fraction occupied by each segmented rock image versus pressure according to equation (3).

[0090] The pore throat resolution parameter N is determined from the mercury capillary injection curve and the voxel size. The pore throat size into which mercury enters at the inlet pressure (Dd) is determined from the mercury capillary injection curve and equation (3), and the pore throat resolution parameter N = Dd / Δx is calculated.

[0091] The pore geometry factor G and / or the Brooks-Corey fitting parameter λ are determined from the mercury capillary injection curve by curve fitting using least squares curve fitting, where a higher weight is given to pores with better resolution in the image. Using G and / or λ, the rock sample (φ ∞ ) pore volume correction factor (α R ) and the corrected pore volume, thereby providing an indication of an underestimate of the true pore volume, for example based on the image pore volume.

[0092] A proprietary direct pore-scale eLBM flow simulator is used to numerically compute relative permeabilities. The direct current simulator takes into account drainage imbibition hysteresis. This approach represents a steady-state type displacement method for relative permeability calculations. A forced drainage process is first simulated. During forced drainage, the in-situ wetting phase that fully saturates the pore space is replaced by an injected non-wetting phase. The rock model is completely saturated with the wetting phase. The drainage process is simulated assuming hydrophilic behavior (constant contact angle of 40 degrees). A buffer layer saturated with oil is placed at the inlet, and a buffer layer saturated with water is placed at the outlet. Oil is injected from the inlet buffer into the pore space at a constant prescribed rate until the average water saturation in the pore space shows no significant change (convergence), at which point the drainage calculation is terminated.

[0093] It should be noted that the water saturation obtained in this way may not necessarily correspond to the irreducible water saturation or coexisting water saturation characterizing the unit volume, because the viscous pressure drop obtained in a rather small computational domain size (much smaller than the core) is significantly smaller than the capillary pressure established in the porous plate or centrifuge method, which corresponds to the saturation height in the reservoir.

[0094] Snapshots of the fluid configurations at different average water saturation values ​​are used for drained relative permeability calculations. The drained fluid configuration is then used as the initial state for a forced imbibition simulation, where the non-wetting phase is replaced by the wetting phase. Prior to the forced imbibition calculation, the fluids in the inlet and outlet buffers are swapped. Forced imbibition calculations are again continued until the average water saturation in the porous domain converges.

[0095] The fluid configuration snapshots at different average water saturation values ​​acquired during the forced imbibition simulations are then used for the relative permeability calculations. To ensure steady-state conditions in the relative permeability calculations, the calculations are performed with cyclic-type periodic boundary conditions applied on all faces of the domain for a given saturation snapshot. A complete mirroring of the porous domain in the mainstream direction is applied so that the periodic boundary conditions are fully satisfied on the outlet faces of the domain. The system is then driven to a steady state, which is very similar to the steady-state relative permeability experiments.

[0096] After the system converges to a steady state, the effective wetting and nonwetting phase permeabilities calculated by eLBM are normalized by the absolute permeabilities calculated using the accurate and efficient MRT-LBM code prior to the relative permeability simulation (Alpak et al. “Prediction of fluid topology and relative permeability in imbibition in sandstone rock by direct numerical simulation”). Advances in Water Resources 122:49-59; 2018). The fluxes used to calculate the effective wetting and nonwetting phase permeabilities are derived by monitoring the average velocities of the wetting and nonwetting phases, respectively, over a portion of the domain corresponding to the porous rock. In other words, the buffer is excluded from the flux calculation. Relatively long buffers are used for forced drainage, and in particular for imbibition process simulations, where the goal is to use the resulting fluid configuration for subsequent relative permeability calculations to minimize boundary effects. Furthermore, additional rings of the interior (rock) domain adjacent to the buffer can be excluded from the effective property calculations to further minimize boundary effects. This approach is typically used for forced imbibition simulations.

[0097] The process is then repeated for all saturation snapshots covering the range of saturations observed during forced drainage or forced imbibition. The relative permeability calculation process is accurate, physically robust, and computationally demanding. The computational demand is directly proportional to the saturation resolution (i.e., the number of saturation snapshots) required to accurately capture the relative permeability function.

[0098] For a predetermined voxel size, the pressure and volume cutoffs are calculated from equations (12) to (15) to provide the non-wetting phase Vclosure and non-wetting phase V nw .

[0099] Examples

[0100] The following non-limiting examples of embodiments of the method of the present invention as claimed herein are provided for illustrative purposes only. The examples illustrate the application of saturation transformations. In the figure, relative permeabilities are plotted for fluid saturations of brine (wetting phase) and hydrocarbon (non-wetting phase). The original results of the direct simulation of the imbibition of the test rock volume are shown as dashed lines 12 and 14 for the oil curve and the water curve, respectively. The original simulation data are then transformed into fluid saturations inferred from the simulation results between saturation cutoffs by correcting the saturation using equation (11). After saturation correction, the transformed simulation results (lines 16, 18 for oil and water, respectively) are consistent with the steady-state relative permeabilities (lines 22, 24 for oil and water, respectively) of laboratory measurements of plugs from the same well.

[0101] The method of the present invention provides a more accurate estimate of fluid saturation and a continuous saturation characteristic curve derived therefrom, which is particularly suitable for exploration and production decision-making and oil field evaluation. The method of the present invention enables rapid decision-making. For example, the process of estimating hydrocarbon saturation using physical measurements of rock samples may take about 7 to 8 months, while the method of the present invention can provide substantially the same or better accuracy in about a few days.

[0102] The present invention is well suited to obtain the objects and advantages mentioned and the objects and advantages inherent therein. The specific embodiments disclosed above are illustrative only, because it is obvious to those skilled in the art who benefit from the teachings herein that the present invention can be modified and practiced in different but equivalent ways. In addition, except as described in the appended claims, there is no intention to limit the details of the construction or design shown herein. Therefore, it is obvious that the specific illustrative embodiments disclosed above can be changed, combined or modified, and all such changes are considered to be within the scope of the present invention. The invention disclosed illustratively herein can be appropriately practiced in the absence of any element not specifically disclosed herein and / or any optional element disclosed herein. Although compositions and methods are described according to "including", "containing" or "comprising" various components or steps, the compositions and methods can also be "essentially composed of various components and steps" or "composed thereof". All numbers and ranges disclosed above can vary by a certain amount. Whenever a numerical range with a lower limit and an upper limit is disclosed, any number and any included range falling within the range are specifically disclosed. In particular, each range of values ​​disclosed herein (in the form of "about a to about b", or equivalently "about a to b", or equivalently "about ab") should be understood to list each number and range covered within the broader range of values. In addition, the terms in the claims have their ordinary general meaning unless otherwise explicitly and clearly defined by the patentee. In addition, the indefinite articles "a" or "an" used in the claims are defined herein to mean one or more of the elements to which they refer. If there is any conflict in the use of words or terms in this specification and one or more patents or other documents that may be incorporated herein by reference, the definition consistent with this specification should be adopted.

Claims

1. A method for estimating fluid saturation of hydrocarbon-bearing rock based on a rock image, comprising: - obtaining a 3D image of rock from a hydrocarbon-bearing formation in an oil field, wherein the 3D image is composed of a plurality of voxels and the 3D image has a resolution; - processing the 3D image to segment the 3D image by selecting each voxel of the 3D image to represent pore spaces in the rock or solid material in the rock; - estimating an image pore volume of the rock based on the segmented 3D image, the image pore volume lacking a sub-resolution pore volume of the rock; - determining a corrected pore volume suitable for taking into account said sub-resolution pore volume of said rock; - estimating an image-derived wetting fluid saturation of the rock using direct flow simulation on the 3D image; as well as - determining a corrected wetting fluid saturation of the rock based on the image-derived wetting fluid saturation of the rock, the image pore volume and the corrected pore volume.

2. The method according to claim 1, wherein the 3D image of the rock is obtained by X-ray computed tomography.

3. The method according to claim 1, wherein the rock is obtained from a hydrocarbon-bearing formation consisting of sandstone, carbonate rock, shale and combinations thereof.

4. The method of claim 1, wherein the corrected pore volume of the rock is determined by laboratory measurements.

5. The method according to claim 1, wherein the corrected pore volume of the rock is determined by: deriving a non-wetting liquid capillary pressure curve from the segmented 3D image of the rock; and determining the corrected pore volume based on the segmented 3D image and the non-wetting liquid capillary pressure curve.

6. A method according to claim 5, wherein the non-wetting liquid capillary pressure curve is derived from the segmented 3D image of the rock at pressures up to an image limit pressure, wherein the image limit pressure is the minimum pressure that can be applied to the non-wetting liquid to overcome the capillary pressure of the narrowest pore throat identifiable from the segmented 3D image of the rock.

7. The method of claim 1, wherein the 3D image is obtained from a cloud-based tool adapted to store and process 2D projection images from pore-scale imaging techniques.

Images