Method for predicting diffusion behavior of underground hydrogen in pore scale of water-containing clay mineral and application

Abstract

The invention discloses a method for predicting the diffusion behavior of underground hydrogen in the pore scale of water-containing clay minerals and application, and belongs to the technical field of underground hydrogen energy storage and unconventional gas migration prediction. Performing energy minimization and pre-balance on the initial pore model of the water-containing clay mineral through a molecular dynamics method; filling hydrogen molecules under given temperature and pressure conditions through a giant regular Monte Carlo method, and correcting the hydrogen pressure through an SRK-EOS state equation; according to the method, the temperature and the fugacity are used as input parameters, hydrogen diffusion simulation is carried out through a molecular dynamics method, the mean square displacement of hydrogen in the pores is calculated, and the self-diffusion coefficient is obtained, and the diffusion rate change rule of the hydrogen in the clay pores under the conditions of different water film thicknesses, salinity and temperature pressure can be quantitatively represented; diffusion inflection points, diffusion suppression stages and migration limitation mechanisms caused by water film or water bridge structures are disclosed.

Description

Technical Field

[0001] This invention belongs to the field of underground hydrogen energy storage and unconventional gas reservoir behavior prediction technology. Specifically, it discloses a molecular-scale simulation method for simulating and predicting the diffusion behavior of hydrogen in the pores of hydrous clay minerals. This method is applicable to the study of hydrogen diffusion behavior in natural hydrogen reservoirs, underground hydrogen storage facilities, shale-clay capping systems, and geological bodies with clay minerals as the main porous medium. It provides a theoretical basis and microscopic mechanism support for the safety assessment of underground hydrogen storage and the prediction of hydrogen retention behavior. Background Technology

[0002] The migration and redistribution of hydrogen in underground formations are key control factors affecting hydrogen resource extraction efficiency, reservoir injectability, and hydrogen storage safety. During this process, hydrogen is influenced by numerous factors, including pore structure, water state, interfacial potential energy distribution, mineral surface charge, and dissolved ions, resulting in complex nonlinear response characteristics in its diffusion process. Underground reservoirs and caprocks are generally rich in clay minerals, whose layered structure, strong hydrophilicity, and surface charge characteristics determine their important micro-regulatory role in hydrogen migration. Especially in the water-bearing state, water films, water bridges, and interfacial water-bound structures easily form in clay pores. The presence of these water-mineral-gas multiphase interfaces significantly affects the hydrogen diffusion path and rate. Therefore, in-depth research on the diffusion behavior of hydrogen in the pores of water-bearing clay is of significant scientific importance for understanding underground hydrogen loss, predicting hydrogen migration range, and assessing reservoir sealing performance.

[0003] Existing research primarily relies on experimental diffusion tests, pore-scale theoretical models, or macroscopic transport equations (such as the Fick diffusion model) for analysis. However, these methods are usually based on idealized assumptions about pores, making it difficult to accurately reflect the water-gas-mineral interface effects in nanoscale pores, and also unable to effectively characterize microscopic mechanisms such as water film thickness variations, interfacial potential barriers, and diffusion channel contraction. Furthermore, actual underground conditions are characterized by large variations in water saturation, wide temperature and pressure ranges, and complex ionic compositions, resulting in significant nonlinear and phased changes in hydrogen diffusion behavior under different pore states. Existing studies have largely focused on dry or fully water-saturated systems, neglecting the gas-water coexistence state commonly found in underground strata. In this unsaturated pore environment, hydrogen diffusion is not only controlled by free gas diffusion but also influenced by interfacial trapping, localized retention, and pore water rearrangement, making its migration characteristics even more complex. Molecular simulation methods (including GCMC and MD) provide important technical means for studying the diffusion of various gases such as CO2 and shale gas in clay pores. However, a systematic model and simulation framework for the key issue of hydrogen diffusion behavior in water-bearing clay pores is currently lacking. Therefore, there is an urgent need for a methodological system that can construct realistic water-bearing clay pore structures at the nanoscale, simulate hydrogen diffusion trajectories, and quantify interface effects. This would reveal the microscopic control mechanisms of hydrogen migration underground, providing a scientific basis for underground hydrogen storage safety analysis, diffusion coefficient calibration, and long-term reservoir stability evaluation. Summary of the Invention

[0004] The purpose of this invention is to propose a method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals. By constructing a pore model of hydrous clay minerals and combining large canonical Monte Carlo and molecular dynamics simulation methods, the diffusion behavior of hydrogen under different water saturation, temperature, pressure and salinity conditions is quantified, revealing the interaction mechanism of the hydrogen-water-clay interface, thereby providing microscopic mechanism support for the safety assessment of underground hydrogen storage and the prediction of reservoir diffusion capacity.

[0005] The method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals according to some embodiments of this application includes the following steps:

[0006] Based on the given surface hydroxyl distribution and force field charge distribution of clay mineral crystal structure parameters, and given the pore interlayer spacing, water molecules and salt ions are filled into the pores according to the target water saturation and salinity to obtain the initial pore model of hydrous clay minerals under different water saturation and salinity conditions.

[0007] The initial pore model of hydrous clay minerals was minimized and pre-equilibriumed using molecular dynamics methods to obtain a stable pore model of hydrous clay minerals.

[0008] Based on the pore stability model of hydrous clay minerals, hydrogen molecules were filled under given temperature and pressure conditions using the giant canonical Monte Carlo method, and the hydrogen pressure was corrected using the SRK-EOS equation of state to obtain the initial model for diffusion simulation.

[0009] Based on the initial diffusion simulation model, with temperature and fugacity as input parameters, hydrogen diffusion simulation was carried out using molecular dynamics methods to calculate the mean square displacement of hydrogen in the pores and obtain the self-diffusion coefficient.

[0010] The method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals according to some embodiments of this application also includes the density distribution of the output hydrogen, and / or

[0011] The self-diffusion coefficient is used to characterize the diffusion behavior of hydrogen in the pores of water-bearing clay, by varying the relationship between water saturation, temperature, pressure and salinity.

[0012] According to the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals in some embodiments of this application, the clay minerals include kaolinite.

[0013] The unit cell parameters of the kaolinite model are a = 5.149 Å, b = 8.934 Å, c = 7.384 Å, α = 91.93°, β = 105.042°, and γ = 89.7910°; where a represents the unit vector length along the x-axis, b represents the unit vector length along the y-axis, and c represents the unit vector length along the z-axis. This represents the angle between the b-axis and the c-axis. This represents the angle between the a-axis and the c-axis. This represents the angle between the a-axis and the b-axis;

[0014] In MaterialsStudio software, the kaolinite model was hydrogenated, with three hydroxyl groups facing the aluminum oxide outer surface and the fourth hydroxyl group facing the siloxane surface, forming a surface hydroxyl distribution. Force field parameters and charge parameters were assigned, and nano-slit pores were constructed by cutting along the (001) crystal plane, where the (001) crystal plane is a set of planes parallel to the clay mineral structure layer plane. A vacuum layer was set and 8×4×2 cell expansion was performed to create two clay mineral substrates, with a given pore size of 50 Å, forming an upper and lower symmetrical wall model with x, y, and z direction dimensions of 41.192 Å, 35.736 Å, and 81.536 Å, respectively.

[0015] Based on the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals in some embodiments of this application, the number of water molecules and salt ions is calculated based on the target water saturation and salinity and added to the pores to obtain an initial pore model of hydrous clay minerals under different water saturation and salinity conditions. The calculation formula is as follows:

[0016]

[0017] in, The average density of water molecules; This represents the molar mass of a water molecule. The pore surface area; The aperture size; is Avogadro's constant.

[0018] The salt ion number is obtained based on the following calculation: the volume of 1 mol of water is calculated to be 18 ml / mol, and the volume of 1 water molecule is 3 × 10⁻⁶. -26 L is determined by the number of water molecules and salinity to determine the number of NaCl molecules.

[0019] According to the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals in some embodiments of this application, the potential energy function of the initial model of hydrous clay minerals consists of bond stretching energy, bond angle energy, van der Waals potential energy, and Coulomb potential energy, and the total potential energy U is shown in the following formula:

[0020]

[0021] In the formula, Represents total energy. Indicates the bonding energy. Indicates bond angle energy. It means that van der Waals can, Represents Coulomb potential energy;

[0022] Among them, bonds and energy It is expressed as follows:

[0023]

[0024] In the formula, The bond stretching force constant; This is the actual bond length; Standard bond length;

[0025] Among them, bond angle energy It is expressed as follows:

[0026]

[0027] In the formula, The constant of the key angle bending force; This refers to the actual bond angle; Standard bond angle;

[0028] Among them, van der Waals can It is expressed as follows:

[0029]

[0030] In the formula, and For different particles, For particles and particles Between potential well depths, For the same particles Between potential well depths, For the same particles Between potential well depths, For particles and particles Collision distance between them For the same particles Collision distance between them For the same particles The collision distance between particles was calculated using the Lorentz-Berthelot mixing rule, where the LJ parameters between different particles were... , , For particles and particles The distance between them;

[0031] Among them, Coulomb potential energy It is expressed as follows:

[0032]

[0033] In the formula, Where is the dielectric constant. and It is a particle and particles The charge, For particles and particles The distance between them.

[0034] According to the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals in some embodiments of this application, an energy minimization and pre-equilibrium method is used to perform energy minimization and pre-equilibrium on the initial pore model of hydrous clay minerals to obtain a stable pore model of hydrous clay minerals, including:

[0035] Clay mineral force field setting Force field, water molecule model set The model, the hydrogen model is set to The unit point model, with the unit system set as follows: , Periodic boundary conditions are set in all three directions, and the atomic style is set to... The time step is set to 1.0. Runtime set to 1 The mineral substrate is set as a rigid body, and the bonds and bond angles of water molecules are set as rigid, while maintaining the number of particles in the system. ,volume and temperature Under certain conditions, in large-scale atomic and molecular parallel simulators A simulation was conducted to obtain a pore stability model for water-bearing clay minerals.

[0036] According to the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals in some embodiments of this application, the giant canonical Monte Carlo method is used to fill hydrogen molecules under given temperature and pressure conditions, including setting the temperature to 300–420°C. Pressure 50–250 The time step is set to 1.0. The runtime is 1 Execute once every 100 steps The cycle consists of 100 attempts to insert or remove hydrogen molecules, and is based on... The criteria determine whether to accept or reject.

[0037] The SRK-EOS equation of state is corrected for hydrogen pressure, and the calculation is as follows:

[0038]

[0039]

[0040]

[0041]

[0042]

[0043]

[0044]

[0045] In the formula, For the fugacity, For pressure, As the compression factor, The fugacity coefficient. The critical pressure. Thermodynamic temperature The critical temperature. The eccentricity factor for hydrogen;

[0046] Based on the initial diffusion simulation model, using temperature and fugacity as input parameters, hydrogen diffusion simulation was conducted using molecular dynamics methods. This included: at a given temperature and fugacity, performing simulations using only molecular dynamics methods to obtain the diffusion behavior and distribution characteristics of hydrogen in mineral pores, with a time step of 1.0. ,exist Running under the ensemble 4 and adopt The thermostat maintains a constant temperature, the first 2 Used for balance, the last 2 Used for trajectory analysis.

[0047] 8. The method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals according to claim 7, characterized in that the self-diffusion coefficient is expressed as follows:

[0048]

[0049] In the formula, For particles diffusion coefficient, For the relevant time, For particles The number of molecules, For particles No. molecule Location at any given moment For particles No. The initial position of each molecule;

[0050] in: This represents the mean square displacement of a hydrogen molecule.

[0051] According to the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals in some embodiments of this application, the hydrogen self-diffusion coefficient is as follows under the external force caused by the concentration difference:

[0052]

[0053] In the formula, For particles diffusion coefficient, The gas constant is... Thermodynamic temperature This refers to the external force applied to hydrogen molecules, with a numerical range of 0.2-0.45. Let t be the position of the hydrogen molecule at time t.

[0054] The method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals according to some embodiments of this application is used to identify the nonlinear response of hydrogen diffusion caused by changes in water saturation, temperature, pressure and salinity.

[0055] Beneficial Effects: The proposed method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals can accurately reveal the migration mechanism of hydrogen molecules within the pores of hydrous clay based on hydrogen diffusion characteristics under different water saturation, temperature, pressure, and salinity conditions, and achieve refined prediction of the diffusion coefficient at the pore scale. Compared with existing studies that rely on macroscopic diffusion models or idealized pore assumptions, this invention constructs a molecular-scale simulation system that can realistically reflect the microenvironment of underground clay pores. It comprehensively considers the combined influence of water film structure, pore geometric constraints, and interfacial potential barriers on the hydrogen diffusion process, thereby more comprehensively capturing key microscopic phenomena such as interfacial trapping, diffusion inhibition zones, and diffusion path rearrangement during hydrogen migration. The method of this invention can identify the nonlinear stages and diffusion inflection points in the diffusion behavior, significantly improving the reliability of underground hydrogen diffusion prediction under complex water-gas coexistence conditions. Furthermore, this invention can be flexibly extended to different water saturation levels, pore structure types, and underground working conditions, and has wide applicability. It can provide important microscopic mechanism foundations and theoretical support for key engineering problems such as hydrogen loss assessment, effective diffusion coefficient calibration, and caprock sealing analysis during underground hydrogen storage. Attached Figure Description

[0056] Figure 1 This is a flowchart illustrating the method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals according to an embodiment of the present invention.

[0057] Figure 2 The following are simulation diagrams: (a) is a simulation diagram under dry conditions in an embodiment of the present invention; (b) is a simulation diagram under 0.3 water saturation conditions in an embodiment of the present invention; and (c) is a simulation diagram under 0.7 water saturation conditions in an embodiment of the present invention.

[0058] Figure 3 This is a schematic diagram of the hydrogen self-diffusion coefficient calculated under different water saturation levels in an embodiment of the present invention.

[0059] Figure 4 This is a schematic diagram of the hydrogen self-diffusion coefficient calculated under different salinities in an embodiment of the present invention.

[0060] Figure 5 This is a schematic diagram of the hydrogen self-diffusion coefficient calculated at different temperatures in an embodiment of the present invention.

[0061] Figure 6This is a schematic diagram of the hydrogen self-diffusion coefficient calculated under different pressures in an embodiment of the present invention. Detailed Implementation

[0062] The embodiments of the present invention will be described in further detail below with reference to the accompanying drawings and examples. The following examples are for illustrative purposes only and should not be construed as limiting the scope of the invention.

[0063] Example: This example provides a method for predicting the diffusion behavior of underground hydrogen at the pore scale of hydrous clay minerals, including the following steps:

[0064] S1. Constructing an initial pore model for hydrous clay minerals: Based on the given surface hydroxyl distribution and force field charge distribution of the clay mineral crystal structure parameters, and the given pore interlayer spacing, water molecules and salt ions are filled into the pores according to the target water saturation and salinity to obtain an initial pore model for hydrous clay minerals under different water saturation and salinity conditions. Specifically, step S1 includes the following steps:

[0065] S101. Kaolinite was selected as the representative mineral. The unit cell parameters were derived from an inorganic crystal structure database to ensure high accuracy in the model construction. The individual unit cell parameters are a = 5.149 Å, b = 8.934 Å, c = 7.384 Å, α = 91.93°, β = 105.042°, and γ = 89.7910°; where a represents the unit vector length along the x-axis, b represents the unit vector length along the y-axis, and c represents the unit vector length along the z-axis. This represents the angle between the b-axis and the c-axis. This represents the angle between the a-axis and the c-axis. This represents the angle between the a-axis and the b-axis.

[0066] S102. In MaterialsStudio software, the kaolinite model is hydrogenated, with three hydroxyl groups facing the aluminum oxide outer surface and the fourth hydroxyl group facing the siloxane surface to form a surface hydroxyl distribution. Force field parameters and charge parameters are assigned, and nano-slit pores are constructed by cutting along the (001) crystal plane. The (001) crystal plane is a set of planes parallel to the clay mineral structure layer plane. A vacuum layer is set and 8×4×2 cell expansion is performed. Then, two clay mineral substrates are created, with a given pore size of 50 Å, thus forming an upper and lower symmetrical wall model with x, y, and z direction dimensions of 41.192 Å, 35.736 Å, and 81.536 Å, respectively.

[0067] S103. Based on the water molecule density and target water saturation and salinity, calculate the number of water molecules and salt ions, and add them to the pores to obtain the initial pore model of water-bearing clay minerals under different water saturation and salinity conditions. The water molecule calculation formula is as follows:

[0068]

[0069] in, This is the average density of water molecules; detailed values ​​can be found in the NIST database. This represents the molar mass of a water molecule. The pore surface area; The aperture size; This is Avogadro's constant. For example, there are 744 water molecules at a water saturation level of 0.3 and 1730 water molecules at a water saturation level of 0.7.

[0070] The determination of the number of salt ions is based on the following assumptions: the volume of the solution is approximately the volume of water; the density of water is 1 g / cm³. 3 The molar mass of a water molecule is 18 g / mol; Avogadro's constant N A ≈6×10 23 / mol; the volume of 1 mol of water is calculated to be 18 ml / mol, and the volume of 1 water molecule is 3 × 10⁻⁶. -26 L is then used to determine the number of NaCl molecules based on the number of water molecules and the salinity. For example, at a water saturation of 0.3 and a salinity of 5M, the number of NaCl molecules is 67, and at a water saturation of 0.7 and a salinity of 5M, the number of NaCl molecules is 156.

[0071] S2. Obtaining a stable porosity model for hydrous clay minerals: Molecular dynamics methods were used to minimize and pre-equilibrium the initial porosity model of hydrous clay minerals under different water saturation and salinity conditions in step S1. This included setting the force field, unit system, boundary conditions, atomic pattern, time step, water-mineral treatment method, and NVT relaxation step to establish a stable porosity model for hydrous clay minerals. Specifically, the settings included: based on the initial porosity model of hydrous clay minerals obtained in step S1, the CLAYFF force field was selected for the clay mineral force field; the SPC / E model, compatible with the CLAYFF force field, was used for water molecules; Michels' unit point model was used for hydrogen; the unit system was set to real; periodic boundary conditions were set in all three directions; the atomic pattern was set to full; the time step was set to 1.0 fs; the running time was 1 ns; the mineral substrate was kept rigid; and the bonds and bond angles of water molecules were kept rigid. Simulations were conducted while keeping the number of particles N, volume V, and temperature T constant.

[0072] The potential energy function of the initial model of the pore size of hydrous clay minerals consists of bond stretching energy, bond angle energy, van der Waals potential energy, and Coulomb potential energy. The total potential energy U is shown in the following equation:

[0073]

[0074] Represents total energy. Indicates the bonding energy. Indicates bond angle energy. It means that van der Waals can, This represents the Coulomb potential energy.

[0075] Among them, bonds and energy It is expressed as follows:

[0076]

[0077] in, The bond stretching force constant is expressed in units of 1 / 3. ; For actual bond length, in units ; Standard bond length, unit .

[0078] Among them, bond angle energy It is expressed as follows:

[0079]

[0080] in, The key angle bending force constant, in units of ; Actual bond angle, unit ; Standard bond angle, unit .

[0081] Among them, van der Waals can It is expressed as follows:

[0082]

[0083] in, and For different particles, For particles and particles Between potential well depths, For the same particles Between potential well depths, For the same particles Between potential well depths, For particles and particles Collision distance between them For the same particles Collision distance between them For the same particles The collision distance between particles is calculated using the Lorentz-Berthelot mixing rule, which determines the LJ parameters between different particles. , , For particles and particles The distance between them.

[0084] Among them, Coulomb potential energy It is expressed as follows:

[0085]

[0086] in, Where is the dielectric constant. and It is a particle and particles The charge, For particles and particles The distance between them.

[0087] S3. Generate the initial model for hydrogen diffusion simulation. Based on the pore stability model of hydrous clay minerals from step S2, use the giant canonical Monte Carlo method to fill the model with hydrogen molecules under given temperature, pressure, water saturation, and salinity conditions. Correct for pressure using the SRK-EOS equation of state to obtain the initial diffusion simulation model. Specifically, the model includes setting up the SRK equation to correlate the compressibility factor and fugacity coefficient, calculated as follows:

[0088]

[0089]

[0090]

[0091]

[0092]

[0093]

[0094]

[0095] in, For the fugacity, For pressure, As the compression factor, The fugacity coefficient. The critical pressure. Thermodynamic temperature The critical temperature. The eccentricity factor of hydrogen is denoted as , where the eccentricity factor of hydrogen is... Critical pressure and critical temperature The values ​​are -0.216, 1.28 MPa, and 33.2 K, respectively.

[0096] The Grand Canonical Monte Carlo (GCMC) method is used to simulate the filling of hydrogen molecules to obtain an initial diffusion model. The temperature is set to 300–420 K, the pressure to 50–250 atm, the water saturation to 0.3, the salinity to 1 M, the time step to 1.0 fs, the run time to 1 ns, and a GCMC cycle is executed every 100 steps. Each cycle contains 100 attempts to insert or remove hydrogen molecules, and the acceptance or rejection is determined according to the Metropolis criterion.

[0097] S4. Conduct hydrogen diffusion simulation. Based on the initial diffusion simulation model from step S3, using temperature and fugacity as input parameters, calculate the mean square displacement (MSD) of hydrogen within the pores using molecular dynamics methods to obtain the self-diffusion coefficient. Further output the hydrogen density distribution and the relationship between the diffusion coefficient and water saturation, temperature, pressure, and salinity to characterize the diffusion behavior of hydrogen in the pores of water-bearing clay. Specifically, the settings are as follows: at a given temperature (e.g., 300 K) and fugacity (e.g., 200 atm), use only the MD method to simulate and study the diffusion behavior and distribution characteristics of hydrogen in mineral pores. The time step remains 1.0 fs, running for 4 ns under the NVT ensemble, and using a Nosé-Hoover thermostat for isothermal control. The first 2 ns are used for equilibration, and the last 2 ns are used for trajectory analysis. Specifically, the hydrogen self-diffusion coefficient is calculated using the following formula:

[0098]

[0099] in: For particles The diffusion coefficient; For the relevant time; For particles The number of molecules; For particles No. The position of each molecule, The mean square displacement of hydrogen molecules represents the average square of the deviation of hydrogen molecules from their initial positions during continuous movement. By tracking the positional changes of particles, the mean square displacement of particles can be obtained. Then, by linearly fitting the mean square displacement curve and processing the slope, the corresponding diffusion coefficient can be obtained. By outputting the mass density distribution of hydrogen along the pore height at a certain number of frames, and averaging multiple frames, the density distribution of hydrogen under a specific condition can be obtained. Then, simulations under multiple conditions are performed to obtain the relationship between density distribution and diffusion coefficient with water saturation, temperature, pressure, and salinity.

[0100] In the presence of external forces such as concentration differences, the hydrogen self-diffusion coefficient can be expressed as:

[0101]

[0102] in, For particles diffusion coefficient, The gas constant is... Thermodynamic temperature This refers to the external force applied to hydrogen molecules, with a numerical range of 0.2-0.45. Let t represent the position of the hydrogen molecules at time t. The concentration difference is due to the solvent, such as pure water or salt water. Finally, the calculated self-diffusion coefficient under external force and the self-diffusion coefficient obtained by calculating the mean square displacement are not significantly different; therefore, only the self-diffusion coefficient obtained by the mean square displacement is included in the attached figure.

[0103] In this embodiment, Figure 3 This is a schematic diagram of the hydrogen self-diffusion coefficient calculated under different water saturation levels in an embodiment of the present invention; Figure 4 This is a schematic diagram of the hydrogen self-diffusion coefficient calculated at different salinities in an embodiment of the present invention; Figure 5 This is a schematic diagram of the hydrogen self-diffusion coefficient at different temperatures calculated in this embodiment of the invention; Figure 6 This is a schematic diagram of the hydrogen self-diffusion coefficient calculated under different pressures in an embodiment of the present invention; as shown. Figure 3-6 As shown, the self-diffusion coefficient of hydrogen under different water saturation levels exhibits a stepwise decreasing trend, indicating that hydrogen migration is hindered under water-containing conditions, which is unfavorable for hydrogen extraction and flow. The self-diffusion coefficient decreases with increasing salinity, and the effect is more pronounced under low salinity conditions, indicating that the presence of salt ions also weakens hydrogen flow. Under dry conditions, the self-diffusion coefficient shows an approximately exponential increasing or decreasing nonlinear trend with increasing temperature or pressure. The presence of water weakens this nonlinear change, and the effect increases with increasing water saturation, indicating that high temperature and low pressure are conducive to hydrogen migration, but the enhancing effect weakens with increasing water saturation. This effectively characterizes the hydrogen diffusion behavior.

[0104] This invention belongs to the field of underground hydrogen energy storage and unconventional gas migration prediction technology, proposing a method for predicting the diffusion behavior of underground hydrogen in the pore scale of hydrous clay minerals. Based on the actual pore structure of clay minerals, the characteristics of the water-gas coupling interface, and the microscopic interaction mechanism of hydrogen on the mineral surface, this method constructs a pore-scale hydrogen diffusion dynamics model. Utilizing a strategy combining statistical mechanics and molecular simulation theory, the initial filling configuration of hydrogen under different water saturation, temperature, pressure, and salinity conditions is generated using the grand canonical Monte Carlo method (GCMC). The mean square displacement of hydrogen in the pores is calculated using molecular dynamics (MD) methods to obtain the self-diffusion coefficient. This invention can quantitatively characterize the variation law of hydrogen diffusion rate in clay pores under different water film thicknesses, salinity, and temperature and pressure conditions, revealing diffusion inflection points, diffusion inhibition stages, and migration-restricted mechanisms caused by water film or water bridge structures.

[0105] The embodiments of the present invention are given for illustrative and descriptive purposes only, and are not intended to be exhaustive or to limit the invention to the forms disclosed. Many modifications and variations will be apparent to those skilled in the art. The embodiments were chosen and described in order to better illustrate the principles and practical application of the invention, and to enable those skilled in the art to understand the invention and to design various embodiments with various modifications suitable for a particular purpose.

Claims

1. A method for predicting the diffusion behavior of subsurface hydrogen in the pore scale of water-bearing clay minerals, characterized by, The method comprises the following steps: The surface hydroxyl distribution and the force field charge distribution are given based on the clay mineral crystal structure parameters, the interlayer spacing of the pores is given, water molecules and salt ions are filled into the pores according to the target water saturation and salinity, and an initial model of the water-containing clay mineral pores under different water saturations and different salinities is obtained; The initial model of the water-containing clay mineral pores is subjected to energy minimization and pre-equilibrium by a molecular dynamics method, and a stable model of the water-containing clay mineral pores is obtained; Based on the stable model of the water-containing clay mineral pores, hydrogen molecules are filled in under given temperature and pressure conditions by a grand canonical Monte Carlo method, and the hydrogen pressure is corrected by an SRK-EOS state equation, and an initial model for diffusion simulation is obtained; Based on the initial model for diffusion simulation, hydrogen diffusion simulation is carried out by a molecular dynamics method with temperature and fugacity as input parameters, the mean square displacement of hydrogen in the pores is calculated, and the self-diffusion coefficient is obtained.

2. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 1, characterized by, The method also comprises outputting the density distribution of hydrogen, and / or The self-diffusion coefficient changes with at least one of water saturation, temperature, pressure and salinity, and is used to characterize the diffusion behavior of hydrogen in the water-containing clay pores.

3. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 1, characterized by, The clay mineral comprises kaolinite; The single cell parameters of the kaolinite model are a = 5.149 A, b = 8.934 A, c = 7.384 A, α = 91.93°, β = 105.042°, γ = 89.7910°; wherein a represents the unit vector length along the x-axis direction, b represents the unit vector length along the y-axis direction, and c represents the unit vector length along the z-axis direction, represents an angle between the b-axis and the c-axis, represents an angle between the a-axis and the c-axis, represents an angle between the a-axis and the b-axis; The kaolinite model is subjected to hydrogenation treatment in the Materials Studio software, three hydroxyl groups are directed to the aluminum-oxygen outer surface, the fourth hydroxyl group is directed to the siloxane surface, the surface hydroxyl distribution is formed, the force field parameters and the charge parameters are distributed, and a nano-slit pore is constructed along the (001) crystal surface, wherein the (001) crystal surface is a set of planes parallel to the structural layer plane of the clay mineral; a vacuum layer is set and 8*4*2 cell expansion is carried out, an upper and lower clay mineral substrate is created, the pore size is set to 50A, and an upper and lower symmetric wall model with x, y and z direction sizes of 41.192A, 35.736A and 81.536A is formed.

4. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 3, characterized by, According to the target water saturation and salinity, the number of water molecules and salt ions is calculated and the pores are added to obtain the initial model of the water-containing clay mineral pores under different water saturations and salinities, the number of water molecules The calculation formula is as follows: wherein, is the average density of water molecules; is the molar mass of water molecules; is the pore surface area; is the pore size; is the Avogadro constant; The number of salt ions was obtained by calculating the volume of 1 mol of water as 18 ml / mol, the volume of 1 water molecule as 3 x 10 -26 L, and the number of NaCl from the number of water molecules and the salinity. Example 2 5. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 1, characterized by, The potential energy function of the initial model of the water-containing clay mineral comprises bond stretching energy, bond angle energy, van der Waals potential energy and Coulomb potential energy, and the total potential energy U is shown in the following formula: wherein represents the total energy, represents the bond energy, represents the bond angle energy, represents the van der Waals energy, represents the Coulomb potential energy; where the bond and energy are represented as follows: wherein is the bond-stretching force constant; is the actual bond length; is the standard bond length; where the bond angle energy is represented as follows: wherein is the bond angle bending force constant; is the actual bond angle; is the standard bond angle; wherein the van der Waals energy is represented as follows: wherein and are different particles, is a particle and a particle the potential well depth between, are the same particles the potential well depth between, are the same particles the potential well depth between, is a particle and a particle the collision distance between, are the same particles the collision distance between, are the same particles the collision distance between, the LJ parameters between different particles are calculated using the Lorentz-Berthelot mixing rules, where, , , is the distance between a particle and a particle ; where the Coulomb potential energy is represented as follows: wherein is the dielectric constant, and is the particle and the charge of the particle is the distance between the particle and the particle ​​ 6. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 5, characterized by, The initial model of the water-containing clay mineral pores is subjected to energy minimization and pre-equilibrium by a molecular dynamics method, and a stable model of the water-containing clay mineral pores is obtained, including: The force field of clay mineral is set as The model of water molecule is set as The model of hydrogen is set as The unit of the unit point model is set as , Periodic boundary conditions are set in three directions, and the atomic pattern is set as , the time step is set as 1.0 , the running time is set as 1 , the mineral base is set as a rigid body, and the bond and bond angle of water molecules are set as rigid, under the condition of keeping the particle number , volume and temperature , the simulation is carried out in a large-scale atomic and molecular parallel simulator to obtain a pore stability model of water-containing clay mineral.

7. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 1, characterized by, The grand canonical Monte Carlo method was used to fill the hydrogen gas molecules under given temperature and pressure conditions, including setting the temperature 300-420 , the pressure 50-250 , the time step 1.0 , the running time 1 , performing once every 100 steps , the loop containing 100 insertion or removal attempts of hydrogen molecules per cycle, and deciding acceptance or rejection according to the criteria; The SRK-EOS state equation is used to correct the hydrogen pressure, and the calculation is as follows: wherein is the fugacity, is the pressure, is the compressibility factor, is the fugacity coefficient, is the critical pressure, is the thermodynamic temperature, is the critical temperature, is the eccentric factor for hydrogen; Based on the initial model of diffusion simulation, the hydrogen diffusion simulation is carried out by the molecular dynamics method with temperature and fugacity as input parameters, including: at a given temperature and fugacity, only using the molecular dynamics method to simulate, obtaining the diffusion behavior and distribution characteristics of hydrogen in the mineral pores, the time step is still 1.0 , in , the system is run for 4 , and the thermostat is used for constant temperature control, the first 2 are used for equilibrium, and the last 2 are used for trajectory analysis.

8. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 7, characterized by, The self-diffusion coefficient is represented as follows: wherein is the diffusion coefficient of the particles , is the correlation time, is the number of molecules of the particles , is the position of the particles at the time , is the position of the particles at the time , is the position of the particles at the initial time wherein: is the mean squared displacement of hydrogen molecules.

9. The method for predicting the diffusion behavior of hydrogen gas in the pore scale of water-bearing clay minerals according to claim 7, characterized by, Under the external force caused by the concentration difference, the hydrogen self-diffusion coefficient is: wherein is the diffusion coefficient of the particles , is the gas constant, is the thermodynamic temperature, is the external force applied to the hydrogen molecule, with a value ranging from 0.2 to 0.45, is the position of the hydrogen molecule at time t.

10. The method for predicting the diffusion behavior of underground hydrogen in water-containing clay mineral pores, for identifying the nonlinear response of hydrogen diffusion caused by changes in water saturation, temperature, pressure and salinity.

Images