Molecular simulation calculation method for mechanical properties of natural gas hydrate and ice nanoindentation

Abstract

The invention belongs to the technical field of natural gas hydrate and ice nano-scale test methods, and provides a natural gas hydrate and ice nano-indentation mechanical property molecular simulation test method, which comprises: S1, determining the size of a molecular dynamics simulation system, the geometric shape and size of an indenter, and simulation parameters; s2, constructing a nanoindentation molecular simulation model; s3, carrying out nanoindentation molecular dynamics simulation; and S4, analyzing and calculating the indentation modulus and the indentation hardness. The method disclosed by the invention can be used for analyzing the nanoindentation modulus, hardness and damage mechanism of the ice and the hydrate at the molecular level, and the experiment cost is greatly reduced.

Description

Technical Field

[0001] The invention belongs to the technical field of natural gas hydrate and ice nanoscale testing methods, and in particular relates to a method for calculating the nanoindentation modulus and hardness of a phase change material based on molecular dynamics simulation. Background Art

[0002] Both gas hydrates and ice are crystalline substances composed of water molecules. The difference is that the holes composed of water molecules in the gas hydrate crystal structure are occupied by natural gas molecules, such as methane and ethane, while there are only water molecules in the ice crystal. Gas hydrates and ice are widely present in the interior and surface of the earth (such as seafloor sediments at the edge of the continents, permafrost zones and polar regions) and extraterrestrial planets (such as Titan). Gas hydrates are not only potential energy materials, but also related to planetary evolution, climate change, geological disasters, etc. Gas hydrate and ice mechanical property testing plays an important role in the exploration, development and utilization of gas hydrates, and is also crucial to the study of the evolution of the earth and planets.

[0003] The mechanical properties of natural gas hydrates and ice can be tested at the macroscopic level using a high-pressure, low-temperature triaxial device; at the micrometer scale, they can be tested using a low-temperature nanoindentation device; and at the molecular scale, they can be obtained using first-principles computational simulations.

[0004] The invention with application number 202310978359.5 discloses a method for testing the nanomechanical properties of natural gas hydrate samples. First, a pure nitrogen dry environment and local liquid nitrogen ultra-low temperature experimental conditions are constructed to keep the inside of the glove box in a slightly positive pressure state, and the temperature of the milling cutter and the positioning table are consistent; then the hydrate sample is transferred and installed in a dry nitrogen environment, and the milling cutter is driven to cut the surface of the natural gas hydrate sample to form a relatively flat surface; then the nanoindentation parameter reliability test, the indentation test with changing the indentation load value and the indentation point, the temperature sensitivity test and the indentation morphology observation are carried out, and the indentation test of the natural gas hydrate sample in an ultra-low temperature and dry environment is realized by relying on a nanoindenter that can be completely placed in an ultra-low temperature environment. On the one hand, the experimental test means of this method is controlled by the size of the indenter, and the laboratory is currently unable to achieve a full coverage test of the indenter size of 1-500nm; on the other hand, it is impossible to reveal the microscopic mechanism of ice and hydrate destruction during nanoindentation at the molecular level. At the same time, the nanoscale mechanical property test that connects the micrometer scale and the molecular scale is currently unattainable by both experiments and first-principles calculations. Moreover, monitoring of polycrystalline hydrates and ice mechanical failure modes in nature is also beyond the reach of current experiments and first-principles computational simulations. Summary of the invention

[0005] In order to solve the problem of difficulty in testing the mechanical properties of nanoscale indentation in the prior art, the present invention proposes a molecular simulation test method for the mechanical properties of natural gas hydrate and ice nanoindentation, which adopts the following technical solutions:

[0006] A molecular simulation test method for the mechanical properties of natural gas hydrate and ice nanoindentation, comprising:

[0007] S1. Determine the system size, indenter geometry and size, and simulation parameters for molecular dynamics simulation

[0008] S2, construction of nanoindentation molecular simulation model;

[0009] S21. Using the single crystal unit cell of methane hydrate sI structure or ice Ih structure as the seed, create a box of size L x ×L y ×L x nm 3 Single crystal or polycrystalline supercell on the left and right;

[0010] S22, establishing a diamond supercell and cutting it into a sphere with a radius of d;

[0011] S23, read in the supercell of natural gas hydrate or ice and the diamond sphere, combine them in the z direction, and reserve a 2 nm space between the diamond and the crystal;

[0012] S3, nanoindentation molecular dynamics simulation;

[0013] S31, performing energy minimization and pre-balancing on the above-obtained model under simulated temperature conditions;

[0014] S32, set the fixed layer and constant temperature layer in the z direction of the simulation system, select a series of pressing speeds V and pressing depths h max , perform nanoindentation molecular dynamics simulation test; delete the molecules entering the top layer in the z direction, the simulation step length is dt = 0.001ps, and the number of simulation steps is

[0015]

[0016] S4. Analysis and calculation of indentation modulus and indentation hardness.

[0017] Furthermore, the step S3 also includes a Poisson's ratio test:

[0018] (1) For the crystal obtained in step S21, NPT simulation and NVT simulation are performed at the simulation temperature to obtain the equilibrium crystal structure;

[0019] (2) In NP xy P zUnder the T ensemble, a compression simulation in the z direction is performed to make the strain ε in the load direction reach 10%, and the strain ε in the xy direction is obtained. l , using the formula v = -ε1 / ε, the Poisson's ratio v is calculated.

[0020] Furthermore, the step S4 comprises:

[0021] S41, drawing the curves of indenter displacement-time, indenter force-time and indenter displacement-indenter force, and distinguishing the elastic compression stage and the plastic deformation stage by the curve shapes;

[0022] S42. Perform curve fitting on the indentation depth-indenter force condition in the elastic compression stage. The fitting formula is:

[0023]

[0024] Where F is the force on the indenter, E* is the simplified indentation modulus, h is the indentation depth, and the elastic modulus E* is obtained by fitting. The relationship between the indentation modulus E and E* is:

[0025] E=(1-v 2 )E * #(3)

[0026] Read F and h from the curve and calculate the indentation hardness H as:

[0027]

[0028] Where A is the equivalent contact area between the indenter and the crystal.

[0029] Furthermore, it also includes: S5. Analysis of nanoindentation fracturing mechanism

[0030] S51, loading the simulation trajectory, analyzing the coordinates of oxygen atoms in the crystal, distinguishing water molecules in different phases of hydrate, ice, and liquid water, analyzing the entire simulation trajectory data, and quantitatively analyzing the indentation fracture process;

[0031] S52. Visualize and analyze the atomic potential energy, shear strain, and Von Mises stress in the output trajectory to explain the indentation rupture mechanism.

[0032] Compared with the prior art, the advantages and positive effects of the present invention are as follows:

[0033] 1. In terms of test size and resolution, compared with the nanoindentation tip size of 500nm and above in the experiment, the molecular simulation of the present invention can test the tip size of 1-500nm, and can analyze the nanoindentation modulus, hardness and failure mechanism of ice and hydrate at the molecular level.

[0034] 2. In terms of test conditions, the molecular simulation method of the present invention can be used to analyze ice and natural gas hydrates under stable temperature and pressure conditions.

[0035] 3. In terms of simulation test cost, computational simulation method is cheaper. BRIEF DESCRIPTION OF THE DRAWINGS

[0036] Figure 1 This is a model diagram of nanoindentation molecular dynamics simulation of an embodiment of the present invention;

[0037] In the figure: 1 is the indenter, 2 is the ice single crystal;

[0038] Figure 2 It is a schematic diagram of the pressure head displacement-pressure head force curve of an embodiment of the present invention. DETAILED DESCRIPTION

[0039] The embodiments of the present application are described in detail below, and examples of the embodiments are shown in the accompanying drawings, wherein the same or similar reference numerals throughout represent the same or similar elements or elements having the same or similar functions. The embodiments described below with reference to the accompanying drawings are exemplary and are only used to explain the present application, and cannot be understood as limiting the present application.

[0040] refer to Figure 1-Figure 2 This embodiment proposes a method for calculating the nanoindentation modulus and hardness of natural gas hydrate and ice based on molecular dynamics simulation, and revealing its fracturing mechanism. The specific steps are as follows:

[0041] 1. Determine the system size, indenter geometry and size, and simulation parameters of molecular dynamics simulation based on the scale of computational simulation and size effect; establish molecular models of single crystals and polycrystalline of natural gas hydrate and ice, as well as molecular models of indenters, and combine them to form the initial structure of simulation; perform Poisson's ratio test, as well as loading, holding and unloading molecular dynamics simulation; divide the loading stage into elastic indentation stage and plastic deformation stage according to the displacement-force curve morphology during loading; read the indentation depth and loading force, calculate the indenter contact area, and fit the above elastic indentation stage curve to obtain the elastic modulus and indentation hardness; analyze the plastic deformation stage and explain the fracturing mechanism. This calculation method can greatly reduce the experimental cost, obtain the mechanical properties of natural gas hydrate and ice under extreme conditions, and dynamically display the fracturing mechanism.

[0042] 1. Determine the system size, indenter geometry and size, and simulation parameters for molecular dynamics simulation:

[0043] 1. Determine the shape and size d of the indenter required for the simulation, for example, a spherical indenter with a radius of 4 nm, d = 4 nm;

[0044] 2. Determine the minimum size of the simulation box: x and y dimensions L x =L y =2d+2×4nm, where 4nm is reserved for each side to avoid the influence of molecular simulation size effect; the box size in the pressing direction (z direction) is L z =L x +2d+9nm, of which 7nm is the boundary layer and heat bath layer of the crystal phase (3nm), the total thickness of the reserved space between the indenter and the crystal (2nm on each side), and the evaporation layer is 2nm. Hydrates and ice will undergo pressure-induced melting when compressed, and the pressure melting range varies with temperature. By setting the reserved layer in the x and y directions, it can be ensured that the melting process does not affect the calculation simulation process, and more accurate results can be obtained; by setting the evaporation layer, it can be ensured that the melting during the simulation process will not cause the simulation box to increase infinitely, ensuring that the simulation time and scale are controllable and reducing the simulation burden.

[0045] The force field parameters for water molecules, methane molecules, and carbon atoms in diamond were determined as mW / Tip4p-Ice, OPLS-UA, and Tersoff, respectively, and the cross interactions between atoms were obtained using literature parameters and the Lorentz–Berthelot rule.

[0046] 2. Establish molecular models of natural gas hydrate and ice single crystals and polycrystals, as well as the molecular model of the pressure head, and combine them to form the simulated initial structure:

[0047] 1. Using ATOMSK software, take the single crystal unit cell of methane hydrate sI structure or ice Ih structure as the seed, and create a box with a size of L x ×L y ×L x nm 3 Single crystal or polycrystalline supercell on the left and right;

[0048] 2. Create a diamond supercell using ATOMSK software and cut it into a sphere with a radius of d using a python program;

[0049] 3. Use LAMMPS to read in the supercell of natural gas hydrate or ice and the diamond sphere, combine them in the z direction, and leave a 2nm space between the diamond and the crystal. The model established is as follows: Figure 1 shown.

[0050] 3. Perform Poisson’s ratio test

[0051] 1. Using the crystal obtained in the above steps, perform 5ns of NPT simulation and NVT simulation at the simulation temperature T to obtain the equilibrium crystal structure;

[0052] 2. In NP xyP z Under the T ensemble, the LAMMPS fix deform command is used to perform compression simulation in the z direction, so that the strain ε in the load direction reaches 10%, and the strain εl in the xy direction is obtained. The Poisson's ratio v is obtained by the formula v=-ε1 / ε, which provides parameter support for the subsequent indentation simulation. By using the same force field parameters as the indentation simulation to calculate the Poisson's ratio v, the entire simulation process can be guaranteed to be self-consistent and the parameters can be controlled, greatly improving the accuracy of the computational simulation.

[0053] 4. Loading, holding and unloading nanoindentation molecular dynamics simulation:

[0054] 1. The model obtained above was energy minimized at the simulation temperature T and 0 bar, and pre-balanced for 5 ns in the NPT ensemble and balanced for 5 ns in the NVT ensemble;

[0055] 2. Set the 0-1nm layer in the z direction of the simulation system as a fixed layer, and the 1-3nm layer as a constant temperature layer, and perform temperature control. Select a series of pressing speeds V, such as And the penetration depth h max , such as 1-dnm. Then use the LAMMPSfix move command to perform nanoindentation molecular dynamics simulation tests. Use the fix evaporate command to delete the molecules that enter the top 2nm in the z direction due to hydrate decomposition, etc., to keep the simulation conditions consistent with the experimental test conditions and prevent the simulation from being interrupted due to the expansion of the simulation system. The simulation step size is dt=0.001ps and the number of simulation steps is

[0056]

[0057] Here, 2 nm is the distance between the indenter and the crystal, and 3 means that the number of steps for pressing, holding, and unloading is the same.

[0058] Conventional nanoindentation molecular simulation uses no boundary conditions in the z direction, which will cause the simulation box to expand infinitely when used for phase change materials, greatly increasing the burden of simulation and making it impossible to use for large system simulation. At the same time, it will cause anomalies when analyzing a series of parameters such as molecular displacement vectors. Using the fix zwalls command will increase the partial pressure of methane molecules during the simulation, which does not meet the experimental test conditions. This simulation method can save simulation time and improve efficiency, simulate the actual experimental process, and ensure the smooth progress of subsequent analysis.

[0059] The optimized V and h are obtained through the above steps. max, repeat the above nanoindentation molecular dynamics simulation at different temperatures, output the simulation trajectory, indentation depth, interaction force between the indenter and the crystal, as well as the potential energy, kinetic energy, shear deformation, Von Mises stress and other results of the atoms in the crystal, ensure the constancy of the system pressure and size during the simulation process, eliminate the influence of phase change on the pressure and size of the simulated system, reproduce the experimental process, and cross-validate with the experimental results.

[0060] 5. Analysis and calculation of indentation modulus and indentation hardness. According to the morphology of the displacement-force curve during loading, the loading stage is divided into an elastic indentation stage and a plastic deformation stage. The indentation depth and loading force are read, the indenter contact area is calculated, and the elastic indentation stage curve is fitted to obtain the elastic modulus and indentation hardness.

[0061] 1. Draw the curves of pressure head displacement-time, pressure head force-time and pressure head displacement-pressure head force (such as Figure 2 shown), distinguished by the shape of the curve Figure 2 The left part of the middle dashed line is the elastic compression stage, and the right part is the plastic deformation stage;

[0062] 2. For the indentation depth-indenter force in the elastic compression stage, the Hertz theory is used to perform curve fitting, and the fitting formula is:

[0063]

[0064] Where F is the force on the indenter, E* is the simplified indentation modulus, and h is the indentation depth. The elastic modulus E* is obtained by fitting, and the fitting curve is as follows: Figure 2 The relationship between the indentation modulus E and E* is:

[0065] E=(1-v 2 )E * #(3)

[0066] Where v is Poisson's ratio, thus obtaining the indentation modulus;

[0067] from Figure 2 The middle curve reads F and h, which are the indenter force and indentation depth respectively. The indentation hardness H is:

[0068]

[0069] Where A is the equivalent contact area between the indenter and the crystal. The above data are plotted against H and h. The curve first rises and then reaches equilibrium. The equilibrium value is the indentation hardness H.

[0070] 6. Analysis of Nanoindentation Fracturing Mechanism

[0071] Load the simulation trajectory and use the CHILL+ algorithm in the Ovito software to analyze the coordinates of oxygen atoms in the crystal and distinguish water molecules in different phases such as hydrates, ice, and liquid water. Analyze the entire simulation trajectory data and quantitatively analyze the indentation fracture process;

[0072] Ovito software is used to visualize the atomic potential energy, shear strain, and Von Mises stress in the output trajectory to explain the indentation rupture mechanism.

[0073] VII. Verification of simulation results

[0074] The molecular dynamics simulation method was used to establish a molecular simulation of nanoindentation with a tip size of 500nm, and the calculated results were compared and verified with the hardness and modulus results of hydrates and ice reported in the existing literature. The water molecule models such as mW, Tip4p / Ice, Tip4p / 2005, and DeepMD were selected for testing to obtain a molecular model that can reproduce the experimental results and be used for subsequent calculation simulations.

[0075] The method of the present invention can reveal the microscopic mechanism of ice and hydrate destruction during nanoindentation at the molecular level, which is impossible to achieve with existing experimental tests. At the same time, it overcomes the experimental errors caused by experimental equipment, such as indentation equipment, temperature sensors, etc. The computational simulation process of the present invention has good repeatability, and better test accuracy can be obtained through multiple parallel simulations.

[0076] The molecular simulation method of the present invention can be analyzed under stable temperature and pressure conditions of ice and natural gas hydrate (0-273K), while the test conditions of the experiment are limited by the type of gas that controls the experimental temperature, such as liquid nitrogen, which can control the temperature of 120-273K. Therefore, the test range of the computational simulation of the present invention is wider, and the experimental conditions are easy to control.

[0077] The computational simulation method of the present invention has lower cost. Taking ice nanoindentation as an example, in a 172980-atom system, 48 Intel Xeon E5-2667 CPU cores were used to perform nanoindentation simulation, which took a total of 2.5 hours. Based on the market price of 0.04 yuan / CPU core hour, the computational simulation cost a total of 4.8 yuan; the experiment took the charging standard of the solid material nanomechanics test system of East China University of Science and Technology as an example, and the single-point test cost 80 yuan / time.

[0078] The above-described embodiments of the present invention do not constitute a limitation on the protection scope of the present invention. Any modification, equivalent substitution and improvement made within the spirit and principle of the present invention shall be included in the protection scope of the claims of the present invention.

Claims

1. A molecular simulation test method for the mechanical properties of natural gas hydrate and ice nanoindentation, characterized in that: include: S1. Determine the system size, pressure head geometry and size, and simulation parameters for molecular dynamics simulation; S2, construction of nanoindentation molecular simulation model; S21. Using the single crystal unit cell of methane hydrate sI structure or ice Ih structure as the seed, create a box of size L x ×L y ×L z nm 3 Single crystal or polycrystalline supercell on the left and right; S22, establishing a diamond supercell and cutting it into a sphere with a radius of d; S23, read in the supercell of natural gas hydrate or ice and the diamond sphere, combine them in the z direction, and reserve a 2 nm space between the diamond and the crystal; S3, nanoindentation molecular dynamics simulation; S31, performing energy minimization and pre-balancing on the above-obtained model under simulated temperature conditions; S32, set the fixed layer and constant temperature layer in the z direction of the simulation system, select a series of pressing speeds V and pressing depths h max , perform nanoindentation molecular dynamics simulation test; delete the molecules entering the top layer in the z direction, the simulation step length is dt = 0.001 ps, and the number of simulation steps is calculated by formula (1) S4. Analysis and calculation of indentation modulus and indentation hardness.

2. The molecular simulation test method for the mechanical properties of natural gas hydrate and ice nanoindentation according to claim 1, characterized in that: The step S3 also includes a Poisson's ratio test: (1) For the crystal obtained in step S21, NPT simulation and NVT simulation are performed at the simulation temperature to obtain the equilibrium crystal structure; (2) In NP xy P z Under the T ensemble, a compression simulation in the z direction is performed to make the strain ε in the load direction reach 10%, and the strain ε in the xy direction is obtained. l , using the formula v = -ε1 / ε, the Poisson's ratio v is calculated.

3. The molecular simulation test method for the mechanical properties of natural gas hydrate and ice nanoindentation according to claim 2, characterized in that: The step S4 comprises: S41, drawing the curves of indenter displacement-time, indenter force-time and indenter displacement-indenter force, and distinguishing the elastic compression stage and the plastic deformation stage by the curve shapes; S42. Perform curve fitting on the indentation depth-indenter force condition in the elastic compression stage. The fitting formula is: Where F is the force on the indenter, E* is the simplified indentation modulus, h is the indentation depth, and the elastic modulus E* is obtained by fitting. The relationship between the indentation modulus E and E* is: Through the indenter force F and the indentation depth h, the indentation hardness H can be calculated as: Where A is the equivalent contact area between the indenter and the crystal.

4. The molecular simulation test method for the mechanical properties of natural gas hydrate and ice nanoindentation according to claim 1, characterized in that: Also includes: S5. Analysis of nanoindentation fracturing mechanism S51, loading the simulation trajectory, analyzing the coordinates of oxygen atoms in the crystal, distinguishing water molecules in different phases of hydrate, ice, and liquid water, analyzing the entire simulation trajectory data, and quantitatively analyzing the indentation fracture process; S52. Visualize and analyze the atomic potential energy, shear strain, and Von Mises stress in the output trajectory to explain the indentation rupture mechanism.

Images