Method for identifying thermal response characteristic temperature zone of heavy oil aggregation structure based on molecular simulation

Abstract

The application adopts a method for identifying a thermal response characteristic temperature zone of heavy oil aggregation structure based on molecular simulation, identifies the thermal response characteristic temperature zone of heavy oil aggregation structure by processing trajectory data, and belongs to the technical field of computational material science, molecular simulation and heavy oil flowability evaluation. The method comprises the following steps: obtaining component composition and experimental characterization parameters of a heavy oil sample, constructing a representative SARA four-component heavy oil molecular model; setting different temperature conditions and performing molecular dynamics simulation; calculating the viscosity, diffusion coefficient, solvent accessible surface area, pi-pi stacking aggregate size, aggregate cohesive energy and interaction energy between the aggregate and surrounding molecules of the heavy oil from the simulation trajectory; and determining the thermal response characteristic temperature zone of the heavy oil aggregation structure according to the progressive response relationship between the above indexes. The application can quantitatively identify the thermal response transition characteristics of the heavy oil aggregation structure from the molecular scale, and provide a molecular scale basis for temperature window division and subsequent flowability regulation in the process of multi-element thermal combination flooding.

Description

Technical Field

[0001] This invention belongs to the fields of computational materials science, molecular simulation and heavy oil flowability evaluation technology, specifically involving a method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation. Background Technology

[0002] Heavy oil, due to its high content of heavy components such as gums and asphaltenes, typically exhibits high viscosity, poor fluidity, and significant development challenges. In the development of deep, high-pressure, and highly heterogeneous heavy oil reservoirs, multi-element thermal flooding, through the synergistic effect of heat, injected gas, and auxiliary control measures, improves the fluidity and sweep efficiency of heavy oil, and has become one of the important technical directions for enhancing heavy oil recovery.

[0003] In multi-component thermal flooding processes, increasing temperature can significantly reduce the viscosity of heavy oil and enhance the migration ability of heavy oil molecules. However, the macroscopic phenomenon of heavy oil viscosity decreasing with increasing temperature does not necessarily mean that the aggregate structure of heavy components within the heavy oil decreases monotonically with increasing temperature. The fluidity of heavy oil is not only directly affected by temperature but also regulated by the aggregate structure of heavy components such as gums and asphaltenes. Heavy component molecules typically contain aromatic rings and polar heteroatoms, which can form aggregate structures through π-π stacking, van der Waals interactions, and electrostatic interactions. The size, stability, surface exposure, and interaction of these aggregate structures with surrounding molecules all affect the viscous resistance of the heavy oil system.

[0004] Existing studies mostly evaluate the impact of heating on the fluidity of heavy oil from the perspective of macroscopic viscosity-temperature relationships or rheological experiments. While these studies can reveal the overall pattern of heavy oil viscosity changes with temperature, they struggle to directly track the dynamic evolution of the aggregation structure of heavy components during heating, and also find it difficult to quantitatively determine the dominant source of viscous drag in different temperature ranges. In particular, when heavy oil viscosity decreases monotonically with increasing temperature, the aggregation structure of heavy components may exhibit non-monotonic evolution characteristics, making it difficult to identify the characteristic temperature range where the aggregation structure transforms based solely on macroscopic viscosity changes.

[0005] Molecular dynamics simulations can resolve the spatial distribution, diffusion behavior, aggregation structure, and intermolecular interactions of different components in heavy oil systems at the molecular scale, providing a quantifiable means to identify the thermal response of aggregate structures of heavy components. Existing molecular simulation methods can be used to study problems such as reservoir fluid occurrence, CO2 flooding, interfacial interactions, diffusion behavior, and changes in interaction energy. For example, by constructing models, performing molecular dynamics simulations, and calculating indices such as density distribution, mean square displacement, radial distribution function, and interaction energy, microscopic mechanisms can be revealed.

[0006] However, the existing technology still lacks a molecular simulation method for the heating process of heavy oil multi-component thermal composite drive, which can comprehensively identify the characteristic temperature range of the thermal response of heavy oil heavy component aggregation structure by integrating indicators such as heavy oil viscosity, diffusion coefficient, solvent-accessible surface area of ​​heavy component components, π-π aggregate size, aggregate internal cohesion energy, and interaction energy between heavy component components and surrounding molecules.

[0007] Therefore, it is necessary to provide a method for identifying the characteristic temperature range of the thermal response of heavy oil aggregate structures based on molecular simulation, so as to identify the evolution characteristics of heavy oil aggregate structures during heating at the molecular scale. Summary of the Invention

[0008] The purpose of this invention is to provide a method for identifying the characteristic temperature range of the thermal response of heavy oil aggregate structures based on molecular simulation. This method addresses the problem in existing technologies where it is difficult to identify the thermal response transition characteristics of heavy oil aggregate structures and determine their characteristic temperature ranges based solely on macroscopic viscosity-temperature relationships. This method achieves the identification of π-π aggregates of heavy oil recombinant components and their characteristic temperature ranges through molecular dynamics trajectory analysis, clustering, structural statistics, and energy feature extraction.

[0009] To achieve the above objectives, the technical solution of the present invention is implemented as follows:

[0010] A method for identifying the characteristic temperature range of thermal response of heavy oil aggregate structures based on molecular simulation, the method includes the following steps:

[0011] S1. Obtain experimental characterization data of the heavy oil sample to be evaluated, and construct a representative SARA heavy oil molecular model containing saturated components, aromatic components, gums and asphaltenes based on the experimental characterization data, wherein gums and asphaltenes are used as heavy components.

[0012] S2. Based on the representative SARA heavy oil molecular model, a heavy oil molecular dynamics simulation system is constructed, and molecular dynamics simulations are performed at multiple target temperatures to obtain the molecular dynamics simulation trajectories at different temperatures.

[0013] S3. Analyze the molecular dynamics simulation trajectory and extract the molecular coordinates, aromatic ring structure, aromatic ring centroid coordinates, and periodic boundary information of the simulation box;

[0014] S4. Identify the π-π stacking contacts between different recombinant molecules based on the centroid coordinates of the aromatic ring, and use a clustering algorithm to obtain the π-π stacking aggregates of the recombinant molecules;

[0015] S5. Based on the molecular dynamics simulation trajectory and the π-π aggregates of the heavy components, calculate the viscosity of the heavy oil, the diffusion coefficient, the solvent-accessible surface area of ​​the heavy components, the size or size distribution of the π-π aggregates, the cohesive energy of the heavy component aggregates, the growth rate of the diffusion coefficient of the heavy components, and the interaction energy between the heavy component aggregates and the surrounding molecules.

[0016] S6. Based on the progressive response relationship between the solvent-accessible surface area of ​​the heavy component, the size or size distribution of the π-π aggregates, the cohesive energy of the heavy component aggregates, the growth rate of the diffusion coefficient of the heavy component, and the interaction energy between the heavy component aggregates and the surrounding molecules, determine the characteristic temperature range of the thermal response of the heavy oil heavy component aggregate structure.

[0017] Furthermore, step S1 specifically includes the following steps:

[0018] S1.1 The experimental characterization data includes one or more of the following: heavy oil density, mass fraction of SARA four components, molecular weight of heavy components by gel permeation chromatography, elemental analysis data of heavy components, nuclear magnetic resonance hydrogen spectrum data, and Fourier transform infrared spectrum data.

[0019] S1.2. The improved Brown-Ladner method is used to construct representative molecular structures of resins and asphaltenes. The improved Brown-Ladner method is used to determine the number of aromatic rings, cycloalkane rings, aromatic carbons connecting hydrocarbon chains, cycloalkane carbons connecting hydrocarbon chains, aliphatic heterocycles, aromatic heterocycles, and heteroatom distribution characteristics in the heavy component molecules.

[0020] S1.3. Based on the mass fractions of the four SARA components and the molecular weights of each representative molecule, determine the number of molecules of saturated fraction, aromatic fraction, resin, and asphaltenes in the representative SARA heavy oil molecular model.

[0021] Furthermore, the multiple target temperatures mentioned in step S2 include 333.15K, 353.15K, 373.15K, 393.15K, 413.15K, 433.15K, 453.15K, and 473.15K; energy minimization, isothermal and isobaric pre-equilibrium, simulated annealing, and equilibrium molecular dynamics simulation at the target temperatures are performed sequentially on the heavy oil molecular dynamics simulation system to obtain the molecular dynamics simulation trajectory for subsequent trajectory analysis and index calculation.

[0022] Furthermore, the aromatic ring structure described in step S3 includes five-membered aromatic rings, six-membered aromatic rings, and aromatic heterocycles containing nitrogen, oxygen, or sulfur heteroatoms in the recombinant molecule.

[0023] Furthermore, in step S4, π-π stacking contacts between different recombinant molecules are identified. When the distance between the centroids of aromatic rings in different recombinant molecules is not greater than a preset distance threshold, it is determined that π-π stacking contacts are formed between the corresponding aromatic rings. The preset distance threshold is 0.45 nm. The clustering algorithm is the DBSCAN clustering algorithm, which takes the centroid coordinates of the aromatic rings as input, and uses the neighborhood radius parameter and the minimum number of samples parameter as clustering conditions to output the cluster label to which the centroid of the aromatic rings belongs in each frame.

[0024] Furthermore, the solvent-accessible surface area of ​​the heavy component mentioned in step S5 is used to characterize the surface exposure degree of the resin and asphaltenes molecules. In the same evaluation system, the smaller the solvent-accessible surface area of ​​the heavy component, the higher the surface burial degree of the heavy component. The size or size distribution of the π-π aggregate includes one or more of the following: the average size of the π-π aggregate, the proportion of higher-order aggregates, and the aggregate size distribution. The cohesive energy of the heavy component aggregate is obtained by statistically analyzing the van der Waals interaction energy and Coulomb interaction energy between different heavy component molecules within the aggregate. The interaction energy between the heavy component aggregate and surrounding molecules includes the non-bonded interaction energy between the heavy component aggregate and one or more of the following: saturated components, aromatic components, resin molecules that did not participate in the formation of aggregates, and asphaltenes molecules that did not participate in the formation of aggregates. The diffusion coefficient of the heavy component is obtained by linearly fitting the mean square displacement curve of the centroid of the heavy component molecules within a long-term diffusion range. The growth rate of the diffusion coefficient of the heavy component is calculated based on the change of the diffusion coefficient of the heavy component at adjacent target temperatures.

[0025] Furthermore, in step S6, firstly, candidate characteristic temperature zones are determined based on the change in the trend of decreasing solvent accessible surface area of ​​heavy components to increasing trend with increasing temperature; secondly, within the candidate characteristic temperature zones, the aggregation enhancement characteristics of heavy components are determined based on the increase in the size of π-π aggregates, the increase in the proportion of higher-order aggregates, or the change in aggregate size distribution towards larger sizes; thirdly, the change in aggregate stability corresponding to the aggregation enhancement characteristics of heavy components is confirmed based on the increase in the absolute value of the internal cohesion energy of heavy component aggregates; finally, the correlation between the changes in the aggregation structure of heavy components and the flow response of heavy oil within the candidate characteristic temperature zones is characterized by combining the growth rate of the diffusion coefficient of heavy components and the interaction energy between heavy component aggregates and surrounding molecules.

[0026] Furthermore, the method also includes outputting the solvent-accessible surface area of ​​the heavy component, the size or size distribution of π-π aggregates, the infill energy of the heavy component aggregates, the growth rate of the diffusion coefficient of the heavy component, the interaction energy between the heavy component aggregates and surrounding molecules, and the range of the characteristic temperature range of the thermal response at different target temperatures.

[0027] Compared with existing technologies, the method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation described in this invention has the following advantages:

[0028] This invention identifies the characteristic temperature range of the thermal response of heavy oil heavy component aggregation structures at the molecular scale through molecular dynamics simulations, avoiding reliance solely on macroscopic viscosity-temperature relationships to determine the thermal response behavior of heavy oil. By coupling and analyzing indicators such as heavy oil viscosity, diffusion coefficient, solvent-accessible surface area of ​​heavy components, π-π aggregate size, aggregate infill energy, and interaction energy between heavy components and surrounding molecules, a more comprehensive characterization of the changes in the aggregation structure of heavy oil heavy components during heating can be achieved. This invention can construct corresponding representative molecular models based on the SARA composition and experimental characterization data of different heavy oil samples, and identify the characteristic temperature range of the thermal response of the heavy component aggregation structures of the corresponding oil samples. It has certain applicability and provides a molecular-scale reference for temperature window segmentation and subsequent flowability control strategies in multi-component thermal composite flooding processes. Attached Figure Description

[0029] The invention will now be further described with reference to the accompanying drawings, in which:

[0030] Figure 1 This invention describes the characteristic of the degree of aggregation of heavy components in the heavy oil system as temperature increases.

[0031] Figure 2 A schematic diagram illustrating the meaning of different parameters in the improved Brown-Ladner method of this invention;

[0032] Figure 3 This invention presents the simulated configuration of heavy oil after equilibrium and the molecular topology of SARA components.

[0033] Figure 4 The shear viscosity and diffusion coefficient of the heavy oil simulation system of this invention as temperature increases;

[0034] Figure 5 The following figures represent the solvent-accessible surface areas of the heavy oil system of the present invention at different temperatures: (a) is a schematic diagram of the solvent-accessible surface areas of the heavy oil system of the present invention at different temperatures; (b) is a curve of the solvent-accessible surface areas of asphaltenes in the heavy oil system of the present invention at different temperatures; and (c) is a curve of the solvent-accessible surface areas of gum molecules in the heavy oil system of the present invention at different temperatures.

[0035] Figure 6 The dimensions of the heavy component aggregates in the heavy oil system of the present invention at different temperatures are: (a) the size distribution of the heavy component aggregates in the heavy oil system of the present invention at different temperatures; and (b) the average size of the heavy component aggregates in the heavy oil system of the present invention at different temperatures.

[0036] Figure 7 This invention provides the energy accumulation of heavy component aggregates in the heavy oil system at different temperatures.

[0037] Figure 8 This refers to the interaction energy between the aggregates of heavy components and surrounding molecules in the heavy oil system of this invention at different temperatures;

[0038] Figure 9 The figures show the simulated average size of the heavy oil system of the present invention at different temperatures, ranging from 20 to 50 ns. Specifically: (a) is the simulated average size of the heavy oil system of the present invention at temperatures ranging from 333.15 K to 393.15 K, ranging from 20 to 50 ns; (b) is the simulated average size of the heavy oil system of the present invention at temperatures ranging from 413.15 K to 473.15 K, ranging from 20 to 50 ns.

[0039] Figure 10 The diffusion capacity of heavy oil as a whole and heavy components under different temperatures is calculated for the present invention, wherein: (a) is the diffusion coefficient of heavy oil as a whole and heavy components under different temperatures; and (b) is the growth rate of the diffusion coefficient of heavy oil as a whole and heavy components under different temperatures. Detailed Implementation

[0040] The present invention will be further described below with reference to specific embodiments. It should be understood that the following embodiments are only for illustrating the technical solutions of the present invention and are not intended to limit the scope of protection of the present invention. All equivalent substitutions, modifications, or improvements made based on the technical solutions of the present invention should be included within the scope of protection of the present invention.

[0041] This embodiment uses heavy oil samples from the Zheng 364 block of the Wangzhuang Oilfield as an example to illustrate a method for determining the characteristic temperature range of the thermal response of heavy oil heavy component aggregate structures based on molecular simulation. The heavy components include resins and asphaltenes, and the characteristic temperature range of the thermal response is determined through a progressive chain of indicators, including the solvent-accessible surface area of ​​the heavy components, the size of π-π aggregates, the internal cohesion energy of the aggregates, the diffusion response of the heavy components, and the interaction energy between the aggregates and surrounding molecules.

[0042] Example 1: Construction and optimization of a representative SARA heavy oil molecular model

[0043] The molecular structure of the four components of heavy oil was determined based on the experimental characterization results of oil samples from the Zheng 364 block of the Wangzhuang oilfield. The density of the oil sample at 20℃ was 0.9546 g / mL. The mass percentages of n-heptane asphaltenes and resins were 9.09% and 33.71%, respectively, with aromatics and saturated components totaling 57.2%. The molecular weights of asphaltenes and resins were determined by gel permeation chromatography (GPC). The elemental composition of asphaltenes and resins is shown in the table below. The oxygen content was calculated by the difference method. The content of the four types of hydrogen atoms was quantitatively analyzed by proton nuclear magnetic resonance spectroscopy, and the results are shown in Table 1.

[0044] Table 1. Elemental composition, molecular weight, and content of four types of hydrogen atoms in heavy oil.

[0045] Structural information was obtained through proton nuclear magnetic resonance (HMR) spectroscopy. Combined with elemental composition analysis, considering the types and amounts of heteroatoms and conventional heterocyclic structures, the Brown-Ladner algorithm was improved to obtain more accurate fused ring and fused heterocyclic structures, including: the number of aromatic rings in the fused ring unit, the number of cycloalkane rings in the fused ring unit, the number of aromatic carbons connecting the hydrocarbon chain, the number of cycloalkane carbons connecting the hydrocarbon chain, the number of aliphatic heterocycles, and the number of aromatic heterocycles. This information was corroborated by Fourier transform infrared (FT-IR) spectroscopy results, and characteristic molecular structures of heavy oil that matched experimental characterization were selected. The positions of different parameters in the molecule in the BL method are shown below. Figure 2 As shown, the calculation formulas involved are shown in Table 2.

[0046] Table 2 Structural parameters of heavy oil components and their significance

[0047] Structural parameters of heavy oil heavy components (asphaltene and resins) were constructed using an improved BL method, and equivalent molecular models of heavy oil were established. All molecular models were constructed using Materials Studio 19.1 software, such as... Figure 3 As shown in Table 3, the number of molecules and molecular formulas of the four components of the heavy oil system are listed. The two colloidal molecules are isomers. The molecular structures of the light components in the system are derived from previous work, with the aromatic and saturated components constructing three types of molecular topologies, respectively.

[0048] Table 3. Molecular formulas, molecular weights, and number of molecules of SARA components in the heavy oil simulation system obtained by the improved BL method.

[0049] The initial configuration of heavy oil molecules was geometrically optimized at the B3LYP / 6-311G(d,p) level using density functional theory (DFT), and Grimme's D3 dispersion correction (Becke-Johnson damping, GD3BJ) was introduced. The calculations were performed using the GaussianE.01 package. Furthermore, the RESP charge was calculated from the optimized configuration and wavefunction information using the Multiwfn program using the restricted electrostatic potential (RESP) method to correct the atomic charges under the initial force field, thereby achieving an accurate description of the intermolecular Coulomb interactions.

[0050] Example 2: Equilibrium Molecular Dynamics (EMD) Simulation

[0051] Representative molecules of the four components of heavy oil constructed in Example 1 were randomly filled into a cubic simulation box with a side length of 10 nm using Packmol software according to a preset number of molecules, thus constructing an initial molecular dynamics model for heavy oil. This initial model was then imported into the GROMACS 2025.2 software package for equilibrium molecular dynamics simulation, with a time step of 1 fs. The CHARMM universal force field was preferably used as the molecular force field parameter. First, the energy of the system was minimized using the steepest descent algorithm and the conjugate gradient algorithm, causing the maximum force on the atoms in the system to converge to 10 kJ·mol⁻¹. -1 ·nm -1 The following steps are performed. A 20 ns pre-equilibrium simulation is then conducted under the NPT ensemble to bring the system volume and density to a stable state. After pre-equilibrium, the system undergoes a 100 ns simulated annealing process, followed by further equilibrium simulations using a Berendsen pressure coupler and a C-rescale pressure coupler to obtain the equilibrium heavy oil system for subsequent thermal response characteristic analysis. During the simulation, three-dimensional periodic boundary conditions are used; chemical bonds connected to hydrogen atoms are constrained using the LINCS algorithm; long-range electrostatic interactions are calculated using the particle mesh Ewald method; and the cutoff radius for short-range non-bonded interactions is set to 1.2 nm. The software, force fields, and parameters described above are merely preferred implementations of this embodiment and do not constitute a limitation on the scope of protection of this invention.

[0052] Example 3: Calculation of Temperature Range Indicators for Identifying the Aggregation Structure of Heavy Oil Components

[0053] Heavy oil is a non-Newtonian fluid exhibiting shear-thinning properties. This patent employs the Couette shear method to calculate the system viscosity, thereby characterizing the rheological response of the heavy oil model. Shear flow in the x-direction is applied to the equilibrium heavy oil model system, creating a stable velocity gradient in the y-direction. The shear rate is defined as γ = ∂μ x / ∂y。 After reaching a non-equilibrium steady state, according to the principles of statistical mechanics, the shear viscosity η can be expressed by the pressure tensor component P. xy The time average is obtained as follows:

[0054]

[0055] Heavy oil exhibits shear-thinning physical characteristics consistent with non-Newtonian fluids under Couette conditions. To simultaneously cover both the Newtonian plateau region at low shear rates and the thinning region at high shear rates, the shear rate was set to 2.5 × 10⁻⁶. -5 Up to 0.1ps -1 The obtained nonlinear curves are quantitatively described using the Carreau model equations:

[0056]

[0057] In the formula, Indicates zero shear viscosity. Viscosity at the high shear rate limit As a characteristic time scale, The apparent viscosity at the zero shear limit can be determined relatively accurately by fitting viscosity data obtained at different shear rates. .

[0058] Statistical analysis and linear fitting were performed on the mean square displacement (MSD) of the molecular centroid of the heavy oil system in equilibrium MD simulations to obtain the diffusion coefficient (D) to examine the flow diffusion capability of the entire system. The formula for calculating the diffusion coefficient is shown below:

[0059]

[0060]

[0061] Where N is the number of particles calculated in the simulation system, r i (0) is the reference position of particle i, r i (t) is the position of particle i after a time interval t from the reference point.

[0062] To quantitatively characterize the π-π packing among heavy components (asphaltene and resins), DBSCAN clustering was used for identification. Atomic coordinates, connectivity data, and periodic boundary conditions were extracted. The neighborhood radius between aromatic rings was 0.45 nm, and the minimum number of points was 2 to ensure that at least two adjacent rings could form effective aggregates. Intramolecular aggregates formed by aromatic rings within the same molecule were excluded during the clustering process to ensure that the identified π-π packing aggregates reflected the aggregation behavior between different heavy component molecules, thus guaranteeing the rationality of the clustering results. The identified nanoaggregate structures were extracted from molecular dynamic trajectories, and their aggregate size distribution, aggregate quantity, and microstructure were analyzed.

[0063] Example 4: Determining the characteristic temperature range of thermal response of heavy oil heavy component aggregation structure

[0064] Based on the simulated trajectories of the heavy oil system at different temperatures obtained in Example 2, this example sequentially calculates the flowability index and the response index of the heavy component aggregation structure, and determines the characteristic temperature range of the thermal response of the heavy component aggregation structure according to the progressive response relationship between the indices. First, the viscosity and diffusion coefficient of the heavy oil system at different temperatures are calculated, and the results are as follows: Figure 4As shown, the viscosity of heavy oil decreases overall and the diffusion coefficient increases overall with increasing temperature, indicating that heating enhances the overall flowability of the heavy oil system. Furthermore, the solvent-accessible surface areas of the resins and asphaltenes at different temperatures were calculated, and the results are as follows... Figure 5 As shown, the solvent-accessible surface area of ​​the heavy component exhibits a characteristic of first decreasing and then increasing with increasing temperature. Specifically, it gradually decreases in the lower temperature range and then gradually increases with further heating, indicating a significant change in the surface exposure of the heavy component within this temperature range. Therefore, the area around 393.15 K to 413.15 K was initially identified as a candidate temperature range for thermal response characteristics. Subsequently, the π-π packing contacts between heavy component molecules at different temperatures were identified, and the size of the π-π aggregates was statistically analyzed. The results are shown below. Figure 6 As shown, near the candidate characteristic temperature range, the average size of the π-π aggregates of heavy components increases, and the proportion of higher-order aggregates increases, indicating that the local aggregation degree between heavy component molecules is enhanced in this temperature range. Furthermore, based on the identified π-π aggregates of heavy components, the integrative energy of the heavy component aggregates at different temperatures is calculated, and the results are as follows: Figure 7 As shown, within the candidate characteristic temperature range, the absolute value of the integrative energy of the recombinant component aggregates increases, indicating that the increase in aggregate size is not simply caused by instantaneous contact, but is accompanied by enhanced internal interactions and improved structural stability. Subsequently, the interaction energies between the recombinant component aggregates and surrounding molecules at different temperatures were calculated, and the results are as follows: Figure 8 As shown, within the candidate characteristic temperature range, the interaction energy between the heavy component aggregates and surrounding molecules weakens, indicating that while the internal interactions of the heavy component aggregates increase within this temperature range, their interaction response with surrounding molecules also changes accordingly. Furthermore, the variation curves of the average size of the heavy component aggregates in the heavy oil system over a simulation time of 20-50 ns at different temperatures were calculated, and the results are shown below. Figure 9 As shown, the average size of the heavy component aggregates near the candidate characteristic temperature range exhibits significant dynamic fluctuations over time, indicating that the heavy component aggregate structures within this temperature range undergo relatively active formation, rearrangement, and adjustment processes. Finally, the diffusion coefficients and growth rates of the heavy components at different temperatures are calculated, and the results are as follows: Figure 10 As shown, the diffusion coefficient of heavy components generally increases with increasing temperature, but its growth rate decreases or becomes passive within the candidate characteristic temperature range, indicating a correlation between the changes in the aggregation structure of heavy components and their kinetic response within this temperature range. (Summary) Figures 5 to 10Based on the analysis results, this embodiment defines the region of 393.15K-413.15K as the characteristic temperature range of the heavy component aggregation structure thermal response for this heavy oil sample. Within this characteristic temperature range, the heavy components exhibit a chain response characteristic characterized by a decrease in solvent-accessible surface area, an increase in the size of π-π aggregates, an increase in the proportion of higher-order aggregates, an increase in the absolute value of the aggregation energy within the aggregates, and a weakening of the interaction energy between the aggregates and surrounding molecules. Therefore, the characteristic temperature range is not determined based on a single viscosity change, but rather through a progressive analysis of the degree of surface exposure of the heavy components, the π-π aggregation structure, aggregate stability, dynamic fluctuation characteristics, diffusion response, and the interaction between the aggregates and surrounding molecules. It should be noted that 393.15-413.15K is the characteristic temperature range of the thermal response obtained under the simulated conditions for the Zheng 364 block heavy oil sample in this embodiment; for other heavy oil samples, the corresponding characteristic temperature range of the heavy component aggregation structure thermal response can be re-determined according to the same analytical procedure.

[0065] The embodiments described above are only used to illustrate the technical solutions of the present invention and are not intended to limit the scope of protection of the present invention. For those skilled in the art, appropriate adjustments, substitutions, combinations or improvements can be made to the specific embodiments based on the content disclosed in this specification without departing from the concept and technical essence of the present invention; all equivalent changes, equivalent substitutions or conventional improvements made based on the technical solutions of the present invention should fall within the scope of protection of the present invention.

Claims

1. A method for identifying the characteristic temperature range of thermal response of heavy oil aggregate structures based on molecular simulation, characterized in that, The method includes the following steps: S1. Obtain experimental characterization data of the heavy oil sample to be evaluated, and construct a representative SARA heavy oil molecular model containing saturated components, aromatic components, gums and asphaltenes based on the experimental characterization data, wherein gums and asphaltenes are used as heavy components. S2. Based on the representative SARA heavy oil molecular model, a heavy oil molecular dynamics simulation system is constructed, and molecular dynamics simulations are performed at multiple target temperatures to obtain the molecular dynamics simulation trajectories at different temperatures. S3. Analyze the molecular dynamics simulation trajectory and extract the molecular coordinates, aromatic ring structure, aromatic ring centroid coordinates, and periodic boundary information of the simulation box; S4. Identify the π-π stacking contacts between different recombinant molecules based on the centroid coordinates of the aromatic ring, and use a clustering algorithm to obtain the π-π stacking aggregates of the recombinant molecules; S5. Based on the molecular dynamics simulation trajectory and the π-π stacked aggregates of the heavy components, calculate the solvent-accessible surface area of ​​the heavy component, the size or size distribution of the π-π stacked aggregates, the infill energy of the heavy component aggregates, the growth rate of the diffusion coefficient of the heavy component, and the interaction energy between the heavy component aggregates and surrounding molecules. S6. Based on the progressive response relationship between the solvent-accessible surface area of ​​the heavy component, the size or size distribution of the π-π aggregates, the cohesive energy of the heavy component aggregates, the growth rate of the diffusion coefficient of the heavy component, and the interaction energy between the heavy component aggregates and the surrounding molecules, determine the characteristic temperature range of the thermal response of the heavy oil heavy component aggregate structure.

2. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, Step S1 specifically includes the following steps: S1.1 The experimental characterization data includes one or more of the following: heavy oil density, mass fraction of SARA four components, molecular weight of heavy components by gel permeation chromatography, elemental analysis data of heavy components, nuclear magnetic resonance hydrogen spectrum data, and Fourier transform infrared spectrum data. S1.

2. The improved Brown-Ladner method is used to construct representative molecular structures of resins and asphaltenes. The improved Brown-Ladner method is used to determine the number of aromatic rings, cycloalkane rings, aromatic carbons connecting hydrocarbon chains, cycloalkane carbons connecting hydrocarbon chains, aliphatic heterocycles, aromatic heterocycles, and heteroatom distribution characteristics in the heavy component molecules. S1.

3. Based on the mass fractions of the four SARA components and the molecular weights of each representative molecule, determine the number of molecules of saturated fraction, aromatic fraction, resin, and asphaltenes in the representative SARA heavy oil molecular model.

3. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, The multiple target temperatures mentioned in step S2 include 333.15K, 353.15K, 373.15K, 393.15K, 413.15K, 433.15K, 453.15K, and 473.15K. Energy minimization, isothermal and isobaric pre-equilibrium, simulated annealing, and equilibrium molecular dynamics simulation at the target temperatures are performed sequentially on the heavy oil molecular dynamics simulation system to obtain the molecular dynamics simulation trajectory for subsequent trajectory analysis and index calculation.

4. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, The aromatic ring structure described in step S3 includes five-membered aromatic rings, six-membered aromatic rings, and aromatic heterocycles containing nitrogen, oxygen, or sulfur heteroatoms in the recombinant molecule.

5. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, In step S4, π-π stacking contacts between different recombinant molecules are identified. When the distance between the centroids of aromatic rings in different recombinant molecules is not greater than a preset distance threshold, it is determined that π-π stacking contacts are formed between the corresponding aromatic rings. The preset distance threshold is 0.45 nm. The clustering algorithm is the DBSCAN clustering algorithm, which takes the centroid coordinates of the aromatic rings as input, and uses the neighborhood radius parameter and the minimum number of samples parameter as clustering conditions to output the cluster label to which the centroid of the aromatic rings belongs in each frame.

6. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, The solvent-accessible surface area of ​​the heavy component in step S5 is used to characterize the surface exposure degree of the resin and asphaltenes molecules. In the same evaluation system, the smaller the solvent-accessible surface area of ​​the heavy component, the higher the surface burial degree of the heavy component. The size or size distribution of the π-π aggregate includes one or more of the following: the average size of the π-π aggregate, the proportion of higher-order aggregates, and the aggregate size distribution. The cohesive energy of the heavy component aggregate is obtained by statistically analyzing the van der Waals interaction energy and Coulomb interaction energy between different heavy component molecules within the aggregate. The interaction energy between the heavy component aggregate and surrounding molecules includes the non-bonded interaction energy between the heavy component aggregate and one or more of the following: saturated components, aromatic components, resin molecules that did not participate in the formation of aggregates, and asphaltenes molecules that did not participate in the formation of aggregates. The diffusion coefficient of the heavy component is obtained by linearly fitting the mean square displacement curve of the centroid of the heavy component molecules within a long-term diffusion range. The growth rate of the diffusion coefficient of the heavy component is calculated based on the change of the diffusion coefficient of the heavy component at adjacent target temperatures.

7. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, In step S6, candidate characteristic temperature zones are first determined based on the change in the trend of decreasing solvent surface area of ​​heavy components to increasing trend with increasing temperature; then, within the candidate characteristic temperature zones, the aggregation enhancement characteristics of heavy components are determined based on the increase in the size of π-π aggregates, the increase in the proportion of higher-order aggregates, or the change in the size distribution of aggregates towards larger sizes. The change in aggregate stability corresponding to the enhanced aggregation feature of the heavy component was confirmed again based on the increase in the absolute value of the aggregation energy within the heavy component aggregates. Finally, the correlation between the change in the aggregation structure of the heavy component and the flow response of heavy oil within the candidate characteristic temperature range was characterized by combining the growth rate of the diffusion coefficient of the heavy component and the interaction energy between the heavy component aggregates and surrounding molecules.

8. The method for identifying the thermal response characteristic temperature range of heavy oil aggregate structures based on molecular simulation according to claim 1, characterized in that, The method also includes outputting the solvent-accessible surface area of ​​the heavy component, the size or size distribution of π-π aggregates, the cohesive energy of the heavy component aggregates, the growth rate of the diffusion coefficient of the heavy component, the interaction energy between the heavy component aggregates and surrounding molecules, and the range of the characteristic temperature range of the thermal response at different target temperatures.

Images