CO2-H2O-shale organic matter three-phase antenna determination method based on molecular simulation

Abstract

The invention discloses a molecular simulation-based CO2-H2O-shale organic matter three-phase antenna determination method, which comprises the following steps of: based on geological deposition characteristics of a research area, preferably selecting a real kerogen molecular model to construct a characterization unit, and constructing a shale organic matter flat substrate model through annealing dynamics simulation; a CO2-H2O-organic matter substrate three-phase molecular model is established, and accurate control over CO2 pressure is achieved in combination with a Voronoi algorithm; splitting an adjacent H2O molecule phase state in an H2O molecule truncation radius range by adopting number density, and dividing the phase state into an evaporation state and a pure liquid state; the height of the spherical-crown-shaped water cluster is represented by adopting a related relational expression constructed by a plane equivalent projection area and a standard liquid volume, and the contact angle is measured through a geometrical relationship rule. According to the method, the contact angle of the CO2 on the surface of the shale organic matter under the reservoir in-situ condition can be rapidly and accurately measured, and technical support is provided for injection parameter optimization and the like in CO2 geological burying.

Description

Technical Field

[0001] The invention relates to the field of oil and gas field development, and in particular to a method for measuring the contact angle of CO2-H2O-shale organic matter three-phase based on molecular simulation. Background Art

[0002] CO2 geological storage is an effective way to mitigate the greenhouse effect. Injecting CO2 into shale reservoirs can not only enhance shale gas extraction and improve shale gas recovery efficiency through adsorption and replacement mechanism, but also achieve effective CO2 storage.

[0003] The surface properties of rocks are key factors affecting the structural and residual storage of CO2. The contact angle of the three-phase CO2-H2O-shale organic matter is a macroscopic manifestation of complex microscopic mechanisms such as van der Waals forces and electrostatic forces, and is crucial for the simulation calculation of CO2 storage efficiency. The wettability of the rock surface is the result of the combined effects of multiple factors such as temperature, pressure, mineral type, and mineralization. At present, due to the limitations of indoor experimental testing methods, it is difficult to carry out contact angle measurements under in-situ reservoir conditions indoors, especially in high temperature and high pressure environments. Even if relevant indoor tests are carried out, it is difficult to calibrate and analyze the contact angle results from a mechanism perspective due to scale limitations. In addition, indoor experiments are limited by experimental samples and time, and it is difficult to quickly carry out experimental calibration for multiple factors on the same experimental sample at the same time. At this point, if we can characterize the CO2-H2O-organic matter three-phase interaction in shale reservoirs at the atomic scale based on the method of intermolecular interaction forces, track the kinetic behavior of fluid molecules in three-dimensional space, and characterize the fluid storage configuration and wettability changes under multi-factor conditions, it will be of great significance for the accurate evaluation of the CO2 storage capacity of shale reservoirs. Summary of the invention

[0004] In order to accurately characterize the three-phase contact angle after supercritical CO2 is injected into the shale reservoir, the present invention proposes a method for measuring the three-phase contact angle of CO2-H2O-shale organic matter based on molecular simulation, thereby realizing the study of the wettability evolution law under the coupled gas-liquid-solid effect and multiple factors at the nanoscale, laying a theoretical foundation for the storage of CO2.

[0005] In order to achieve the above technical objectives, the present invention adopts the following technical solutions:

[0006] A method for determining the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation comprises the following steps:

[0007] S1: According to the element analysis results of shale samples in the study area, the reservoir sedimentary type, kerogen type and maturity are determined, and a three-dimensional kerogen molecular characterization unit is constructed;

[0008] S2: calibrating the atomic coordinates and force field interaction information of the kerogen molecular monomer model established in step S1 to obtain the kerogen monomer topological structure;

[0009] S3: Using the canonical ensemble (NVT) to obtain the basement of the shale organic wall through the annealing kinetic simulation method;

[0010] S4: According to the geometric spatial distribution characteristics, sufficient fluid molecular models, including H2O molecules and CO2 molecules, are added to the wall surface, and the stable configuration of the fluid distribution is obtained through equilibrium molecular dynamics (EMD);

[0011] S5: Extract the coordinates of fluid molecular trajectories, take the O atom in the H2O molecule as the center, capture the density distribution relationship of water clusters on the longitudinal plane, and divide the liquid-solid and gas-liquid phase interfaces;

[0012] S6: fitting and solving the water cluster density distribution depicted in step S5 to quantitatively characterize the contact angle.

[0013] Furthermore, the kerogen monomer molecule in step S2 is geometrically optimized by the conjugate gradient method until the energy reaches 1.0×10 -4 The convergence value of kcal / mol results in stable kerogen monomer molecules.

[0014] Furthermore, in step S3, the van der Waals term is described by the Lennard-Jones 12-6 potential energy function model, the unpaired atomic interactions are approximated by the Lorentz-Berthelot (LB) mixing rule, and the temperature field is controlled and constrained by the Nosé-Hoover algorithm.

[0015] Furthermore, in step S4, periodic boundary conditions are applied to the system in three directions, with the aim of overcoming the limitation of the small scale of the simulation system and characterizing the infinite system with a finite system; in the vertical direction, in order to avoid the influence of mirror mapping, a virtual reflection boundary controlled by the Lennard-Jones 9-6 model is introduced at the upper and lower boundaries.

[0016] Furthermore, in step S4, the gas pressure is obtained by respectively calculating the atomic stress tensor and the gas molecular volume; the volume of the CO2 gas molecule is calculated using the Voronoi algorithm, and a continuous polygonal structure is constructed by dividing the Thiessen polygons to construct the perpendicular bisectors of the straight lines connecting two adjacent points, so as to obtain the volume of the supercritical CO2 system and the local pressure.

[0017] Furthermore, the method for dividing the liquid state and the evaporated state of H2O molecules in step S5 is:

[0018] A neighbor list is established with a single H2O molecule as the origin to retrieve and count the H2O molecules within the cutoff radius, and the number density distribution of the H2O molecules in each calculation unit is calculated;

[0019] Based on the standard number density of H2O molecules measured experimentally under corresponding temperature and pressure conditions, the equilibrium distance between liquid H2O molecules can be calculated:

[0020]

[0021] In the above formula, d equ represents the equilibrium distance between molecules under liquid H2O conditions, ρ liquid It indicates the experimentally measured number density of liquid H2O. The number density of H2O under different temperature and pressure conditions can be found in the database of the National Institute of Standards and Technology (NIST).

[0022] Considering the error caused by atomic offset due to transient imbalance of intermolecular interactions, the minimum number of H2O molecules per unit volume can be calculated using 2 times the equilibrium distance as a constraint:

[0023]

[0024] In the above formula, N liquid represents the number of H2O molecules per cubic angstrom volume after taking into account the local energy error, .

[0025] If we determine the number of H2O molecules N in a spatial unit <N liquid , the H2O molecules in the unit are classified as evaporated and deleted from the statistical list to save computing resources. Otherwise, they are marked as liquid H2O molecules and retained for contact angle calculation.

[0026] Further, for the contact angle measurement in step S6, an indirect solution is performed by the following method:

[0027] For the retained liquid H2O molecules, the system is discretized in the plane space. The area occupied by liquid H2O molecules is characterized by a solid structure, and the area without liquid H2O molecules is a blank structure. In this way, a nearly circular solid area can be obtained on the two-dimensional projection surface. The area corresponding to this area is S, and the radius is r.

[0028] The liquid H2O molecule region is retained, and its volume V can be approximately calculated by the number density. liduid :

[0029]

[0030] For the volume of the regional water clusters in the region being solved, the volume model of the local spherical cap is constructed through the idea of ​​calculus, and the correlation between the height h of the water cluster and the volume V is established to determine the height h;

[0031]

[0032] The contact angle θ can be solved through simple geometric relationships.

[0033]

[0034] The present invention provides a method for determining the contact angle of three phases of CO2-H2O-shale organic matter based on molecular simulation. On the basis of the study of the geological sedimentary characteristics of the study area, a real kerogen molecular model is selected to construct a three-dimensional characterization unit, and an annealing kinetic simulation method is used to establish a shale organic matter flat plate base model. Fluid molecules such as CO2 and H2O are added on top of it, and the thermodynamically stable configuration is obtained by the equilibrium molecular dynamics simulation method. The pressure is calibrated by the Thiessen polygons combined with the Voronoi algorithm. This method identifies and calibrates the evaporated H2O molecules and liquid H2O molecules through number density division, extracts the coordinates and thermodynamic information of the pure liquid H2O molecules, and reconstructs the data volume file. The plane equivalent projected area and the standard liquid volume are combined to establish a related correlation formula, the liquid water cluster height is characterized, and the contact angle is solved by the geometric relationship law. This method can realize the rapid measurement of the contact angle of the shale surface. At the same time, through the number density division, it can get rid of the problem of inaccurate contact angle measurement caused by the density profile fitting method due to the unstable and irregular edge shape. This method can be extended to different organic matter types and maturity, different inorganic mineral crystal surfaces, etc. It has wide application value and can provide theoretical support for CO2 geological storage site selection, injection parameter optimization, etc.

[0035] Beneficial effects:

[0036] Compared with the prior art, the present invention has the following beneficial effects:

[0037] 1. Based on the real kerogen model, a shale organic matter flat plate model was established, which can realize the characterization of the occurrence mechanism under the complex interaction of gas-liquid-solid three phases under in situ conditions.

[0038] 2. Combined with the Voronoi algorithm, the accurate calibration of the local pressure of the supercritical CO2 system is achieved through the Thiessen polygons.

[0039] 3. Coarse-graining H2O molecules into particles can effectively improve the efficiency and accuracy of the number density distribution meter of water clusters in the longitudinal plane, reduce the influence of geometric unevenness due to local imbalance of the system, and have higher accuracy than the traditional density distribution contour fitting method.

[0040] 4. The contact angle measurement workflow at the atomic scale of the shale organic matter wall proposed in the present invention can realize the contact angle measurement on various mineral surfaces and in multiple temperature and pressure ranges, making up for the dilemma of limited indoor experimental conditions. BRIEF DESCRIPTION OF THE DRAWINGS

[0041] Figure 1 Flow chart of the method for calibrating the contact angle of the real rock substrate surface described in the embodiment of the present invention

[0042] Figure 2 Initialize the configuration of the fluid model described in the embodiment of the present invention

[0043] Figure 3 The coarse-grained configuration of H2O molecules in equilibrium described in the embodiment of the present invention is

[0044] Figure 4 Schematic diagram of using thresholds to divide vaporized and liquid H2O molecules in an embodiment of the present invention

[0045] Figure 5 Schematic diagram of the contact angle calibration method in an embodiment of the present invention DETAILED DESCRIPTION

[0046] In order to make the purpose, technical solution and advantages of the present invention more clearly understood, the present invention is further described in detail in conjunction with the following embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention and are not used to limit the present invention.

[0047] Example:

[0048] A method for measuring the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation was used to calculate the contact angle of CO2-H2O-shale organic matter three phases in a shale in the Sichuan Basin. Figure 1 The process shown includes the following steps:

[0049] According to the element analysis results of shale samples in the Sichuan Basin, which is a typical marine sedimentary environment, the IID type kerogen molecular unit (C 175 H 102 N4O9S2) to construct a three-dimensional characterization unit.

[0050] The atomic coordinates of the kerogen molecular monomer model established in step S1 are calibrated, and its atomic types are calibrated using a consistent force field to divide atoms into types.

[0051] Table 1 Basic parameters of organic matter model

[0052]

[0053] The van der Waals term is described by the Lennard-Jones 12-6 potential energy function model, and the Coulomb long-range force is summed using the Ewald algorithm, which is implemented through real space and Fourier space transformation.

[0054]

[0055] Among them, r ij is the distance between the centers of two atoms, σ ij is the potential parameter related to the atomic size, ε ij is the depth of the potential well of the interaction potential between two atoms, kcal / mol; q i and q j is the charge carried by the atom; ε0 is the dielectric constant of free space.

[0056] The interactions between different atoms are approximately characterized by the Lorentz-Berthelot mixing law.

[0057]

[0058] For the organic matter system, the initial velocity field of the atoms in the system is created according to the initial temperature of the model, and its distribution characteristics obey the Gaussian distribution. The canonical ensemble (NVT) is applied to the kerogen molecules in the entire system to carry out annealing dynamics simulation. The temperature of the system needs to be controlled and constrained by the Nosé-Hoover algorithm during the simulation process.

[0059] like Figure 2 As shown in the figure, 960 H2O molecules are first added in the middle position above the organic substrate, and the extended single point charge (SPC / E) model is used to describe its properties; then a certain number of CO2 molecules are arrayed around the H2O molecules, which are characterized by the EPM2 model. Periodic boundary conditions are imposed on the three directions of the simulation system in order to overcome the limitation of the simulation system being too small and to characterize the infinite with the finite. In the vertical direction, in order to avoid the influence of mirror mapping, a "virtual reflection boundary" controlled by the LJ 9-6 model is introduced at the upper and lower boundaries.

[0060]

[0061] Under the conditions of a canonical system, equilibrium molecular dynamics (EMD) was performed on H2O and CO2 molecules until the system reached thermodynamic equilibrium, and key information such as atomic coordinates of the system was recorded at each time step. During the simulation, CO2 pressure was obtained by calculating the atomic stress tensor and gas molecular volume respectively; the volume of CO2 gas molecules was obtained using the Voronoi algorithm, and the perpendicular bisector of the straight line connecting two adjacent points was constructed into a continuous polygonal structure by dividing the Thiessen polygon to obtain the volume of the supercritical CO2 system. In order to obtain the local pressure, the pressure was adjusted to the research range by adjusting the number of CO2 molecules. During the simulation calculation, under the condition of a 1fs time step, the first 2.4ns was used for the system to reach thermodynamic equilibrium, and the last 4.6ns was used for data collection, with a total simulation time of 6ns. In the data collection stage, the pressure was averaged over time, and the collection frequency was 10ps / point within every 400ps, and a total of 400 data points were collected for ensemble averaging. Through simulation, it was found that under the conditions of 330K and 20MPa respectively, the number of CO2 was 3100, and the simulated pressure was basically consistent with the set pressure, with high accuracy.

[0062] Table 2 CO2 system pressure control results

[0063] Simulation steps Pressure(MPa) 2400000 19.94 2800000 20.18 3200000 20.36 3600000 20.84 4000000 20.24 4400000 20.32 4800000 20.35 5200000 20.10 5600000 20.29 6000000 20.28

[0064] Extract the molecular dynamic evolution trajectory file under the simulation system. In order to save computing resources, first delete the organic base wall and CO2 molecules, and only keep H2O molecules. Figure 3 As shown, H2O is coarse-grained into a particle with geometric configuration as constraint, and a neighbor list is established with a single H2O molecule as the origin to perform frame-by-frame search and statistics on the H2O molecules within the truncation radius. The number density distribution of H2O molecules in each calculation unit is calculated, and the number density distribution law is stored and calibrated.

[0065] Based on the standard number density of H2O molecules measured experimentally under the research temperature and pressure conditions, the equilibrium distance between liquid H2O molecules can be calculated:

[0066]

[0067] In the above formula, d equ represents the equilibrium distance between molecules under liquid H2O conditions, ρ liquid It indicates the experimentally measured number density of liquid H2O. The number density of H2O under different temperature and pressure conditions can be found in the database of the National Institute of Standards and Technology (NIST). Under 330K and 20MPa, the number density of H2O molecules is about Therefore, the equilibrium distance of H2O under this condition is approximately

[0068] At the same time, considering the error caused by atomic offset due to transient imbalance of intermolecular interactions, the minimum number of H2O molecules per unit volume can be calculated using 2 times the equilibrium distance as a constraint.

[0069]

[0070] In the above formula, N liquid It represents the number of H2O molecules per cubic angstrom volume after considering the local energy error. It can be concluded that the number of H2O molecules per cubic angstrom space is about 0.037 after calibration by the molecular equilibrium distance under 330K 20MPa conditions.

[0071] like Figure 4 As shown in the figure, the phase state of H2O molecules is divided and determined based on the minimum number of H2O molecules in a unit volume. If the number of H2O molecules N in the space unit is determined <N liquid , the H2O molecules in the unit are classified as evaporated, deleted from the statistical list, and not used for contact angle calculation; otherwise, they are marked as liquid H2O molecules and retained for contact angle calculation.

[0072] For the retained liquid H2O molecules, the system is discretized in the plane space. The area occupied by liquid H2O molecules is characterized by a solid structure, and the area without liquid H2O molecules is a blank structure. In this way, a nearly circular solid area can be obtained on the two-dimensional projection surface. The area corresponding to this area is S, and the radius is r.

[0073] The liquid H2O molecule region is retained, and its volume V can be approximately calculated by the number density. liquid :

[0074]

[0075] For the volume of the regional water clusters to be solved, the volume model of the local spherical cap is constructed through the calculus idea, and the correlation between the height h of the water cluster and the volume V is established to determine the height h:

[0076]

[0077] like Figure 5 As shown, the contact angle can be solved through simple geometric relationships:

[0078]

[0079] The contact angle measured in this example is 175°, showing a strong hydrophobic property.

[0080] The present invention provides a method for determining the contact angle of three phases of CO2-H2O-shale organic matter based on molecular simulation. On the basis of the study of the geological sedimentary characteristics of the study area, a real kerogen molecular model is selected to construct a three-dimensional characterization unit, and an annealing kinetic simulation method is used to establish a shale organic matter flat plate base model. Fluid molecules such as CO2 and H2O are added on top of it, and the thermodynamically stable configuration is obtained by the equilibrium molecular dynamics simulation method. The pressure is calibrated by the Thiessen polygons combined with the Voronoi algorithm. This method identifies and calibrates the evaporated H2O molecules and liquid H2O molecules through number density division, extracts the coordinates and thermodynamic information of the pure liquid H2O molecules, and reconstructs the data volume file. The plane equivalent projected area and the standard liquid volume are combined to establish a related correlation formula, the liquid water cluster height is characterized, and the contact angle is solved by the geometric relationship law. This method can realize the rapid measurement of the contact angle of the shale surface. At the same time, through the number density division, it can get rid of the problem of inaccurate contact angle measurement caused by the density profile fitting method due to the unstable and irregular edge shape. This method can be extended to different organic matter types and maturity, different inorganic mineral crystal surfaces, etc. It has wide application value and can provide theoretical support for CO2 geological storage site selection, injection parameter optimization, etc.

[0081] The above description does not limit the present invention in any form. Although the present invention discloses the effect of the implementation of the method as described in the above embodiment, the embodiment is not used to limit the present invention. Any technician familiar with the profession can make some changes or modify the above disclosed technical contents into equivalent embodiments without departing from the scope of the technical solution of the present invention. However, any simple modification, equivalent change and modification made to the above embodiments according to the technical essence of the present invention without departing from the content of the technical solution of the present invention are within the scope of the technical solution of the present invention.

Claims

1. A method for determining the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation, comprising the following steps: S1: According to the element analysis results of shale samples in the study area, the reservoir sedimentary type, kerogen type and maturity are determined, and a three-dimensional kerogen molecular characterization unit is constructed; S2: calibrating the atomic coordinates and force field interaction information of the kerogen molecular monomer model established in step S1 to obtain the kerogen monomer topological structure; S3: Using the NVT canonical ensemble, the base of the shale organic wall is obtained by annealing kinetic simulation method; S4: According to the geometric spatial distribution characteristics, a sufficient amount of fluid molecular models, including H2O molecules and CO2 molecules, are added to the wall surface, and the stable configuration of the fluid distribution is obtained by equilibrium molecular dynamics; S5: Extract the coordinates of fluid molecular trajectories, take the O atom in the H2O molecule as the center, capture the density distribution relationship of water clusters on the longitudinal plane, and divide the liquid-solid and gas-liquid phase interfaces; S6: fitting and solving the water cluster density distribution depicted in step S5 to quantitatively characterize the contact angle; The kerogen monomer molecule in step S2 is geometrically optimized by the conjugate gradient method until the energy reaches 1.0×10 -4 The convergence value of kcal / mol resulted in stable kerogen monomer molecules; In step S3, the van der Waals term is described by the Lennard-Jones 12-6 potential energy function model, the unpaired atomic interactions are approximated by the Lorentz-Berthelot mixing rule, and the temperature field is controlled and constrained by the Nosé-Hoover algorithm.

2. A method for measuring the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation as claimed in claim 1, characterized in that: In step S4, periodic boundary conditions are applied to the system in three directions, in order to overcome the limitation of the small scale of the simulation system and to characterize the infinite system with a finite system; in the vertical direction, in order to avoid the influence of mirror mapping, a virtual reflection boundary controlled by the Lennard-Jones 9-6 model is introduced at the upper and lower boundaries.

3. A method for measuring the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation as claimed in claim 1, characterized in that: In step S4, the gas pressure is obtained by respectively calculating the atomic stress tensor and the gas molecular volume; the volume of the CO2 gas molecule is calculated using the Voronoi algorithm, and a continuous polygonal structure is constructed by dividing the Thiessen polygons to construct the perpendicular bisectors of the straight lines connecting two adjacent points, thereby obtaining the volume of the supercritical CO2 system and obtaining the local pressure.

4. A method for measuring the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation as claimed in claim 1, characterized in that: The method for dividing the liquid state and the evaporated state of H2O molecules in step S5 is: A neighbor list is established with a single H2O molecule as the origin to retrieve and count the H2O molecules within the cutoff radius, and the number density distribution of the H2O molecules in each calculation unit is calculated; Based on the standard number density of H2O molecules measured experimentally under corresponding temperature and pressure conditions, the equilibrium distance between liquid H2O molecules can be calculated: In the above formula, d equ represents the equilibrium distance between molecules under liquid H2O conditions, ρ liquid It indicates the experimentally measured number density of liquid H2O. The number density of H2O under different temperature and pressure conditions can be found in the database of the National Institute of Standards and Technology of the United States; Considering the error caused by atomic offset due to transient imbalance of intermolecular interactions, the minimum number of H2O molecules per unit volume can be calculated using 2 times the equilibrium distance as a constraint: In the above formula, N liquid represents the number of H2O molecules per cubic angstrom volume after considering the local energy error, ; If we determine the number of H2O molecules N in a spatial unit <N liquid , the H2O molecules in the unit are classified as evaporated and deleted from the statistical list to save computing resources. Otherwise, they are marked as liquid H2O molecules and retained for contact angle calculation.

5. A method for measuring the contact angle of CO2-H2O-shale organic matter three phases based on molecular simulation as claimed in claim 1, characterized in that: For the contact angle measurement in step S6, an indirect solution is performed by the following method: For the retained liquid H2O molecules, the system is discretized in the plane space. The area occupied by liquid H2O molecules is characterized by a solid structure, and the area without liquid H2O molecules is a blank structure. In this way, a nearly circular solid area can be obtained on the two-dimensional projection surface. The area corresponding to this area is S, and the radius is r. The liquid H2O molecule region is retained, and its volume V can be approximately calculated by the number density. liduid : For the volume of the regional water clusters to be solved, the volume model of the local spherical cap is constructed by using the calculus idea, and the correlation between the height h of the water cluster and the volume V is established to determine the height h; The contact angle θ can be solved through simple geometric relationships: Using the measured parameters, the contact angle can be calculated using the above formula.

Images