Fractured shale oil seepage-stress field coupling simulation method based on virtual element
By using a virtual element-based coupled simulation method for seepage-stress field in fractured shale oil, the problem of accurately characterizing the discontinuous displacement characteristics of complex fracture networks after fracturing in existing technologies is solved, enabling efficient simulation and production capacity prediction of shale oil reservoirs.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-12-30
- Publication Date
- 2026-04-10
AI Technical Summary
Existing numerical simulation technologies for shale oil fracturing development are insufficient to accurately characterize the discontinuous displacement features of complex fracture networks after fracturing, and cannot effectively capture the dynamic coupling mechanism of seepage-stress field, resulting in insufficient accuracy in production capacity prediction.
A virtual element-based simulation method for coupled flow-stress field in fractured shale oil is adopted. By constructing a solid deformation model and a fluid flow model of the fractured shale oil reservoir, and performing iterative calculations using the virtual element method, a coupled flow-stress field model is established. The permeability is updated using the cubic relationship between porosity and permeability, integrating the elastic deformation and stress-sensitive effects of porous media.
It enables efficient characterization of complex multi-scale fracture networks, enhances the ability to characterize highly heterogeneous shale reservoirs, improves the reliability of production capacity and reservoir dynamic prediction, and accurately reflects the time-varying laws of reservoir physical parameters during development.
Smart Images

Figure CN121835494A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the technical field of oilfield development, and in particular to a virtual element-based method for simulating the flow-stress field coupling of fractured shale oil. Background Technology
[0002] As a key unconventional resource in the global energy structure transformation, the commercial development of shale oil heavily relies on large-scale hydraulic fracturing technology to modify reservoirs. After fracturing, shale reservoirs form a complex multi-scale fracture network consisting of primary natural fractures, secondary hydraulic fractures, and nanopores. This network exhibits strong heterogeneity, anisotropy, and significant stress sensitivity, leading to a strong flow field-stress field coupling effect between fluid flow and reservoir deformation. This significantly increases the uncertainty in predicting the dynamics of shale oil production.
[0003] Existing numerical simulation technologies for shale oil fracturing development suffer from several key shortcomings: traditional fluid-structure interaction models, limited by structured mesh generation, struggle to accurately characterize the discontinuous displacement characteristics of complex fracture networks after fracturing, and fail to effectively capture the dynamic coupling mechanism of seepage-stress fields; furthermore, most models treat shale reservoirs as rigid porous media, neglecting the interaction between porous media elastic deformation, stress-sensitive effects, and fluid seepage, leading to a systematic overestimation of production capacity. While existing coupling models (such as equivalent continuous medium models and dual-medium models) have been applied, they fail to fully quantify the synergistic seepage mechanism of multi-scale pore structures and lack dynamic characterization of the stress-sensitive effect of fracture permeability, failing to accurately reflect the time-varying patterns of reservoir physical parameters during development, resulting in insufficient accuracy in predicting the dynamics of shale oil depletion development. Summary of the Invention
[0004] The purpose of this invention is to provide a virtual element-based simulation method for the coupled flow-stress field of fractured shale oil, thereby solving the above-mentioned problems.
[0005] To achieve the above objectives, this invention provides a virtual element-based method for coupled simulation of seepage and stress fields in fractured shale oil, comprising the following steps: S1. Construct a solid deformation model for fractured shale oil reservoirs; S2. Construct a fluid flow model for fractured shale oil reservoirs; S3. By correlating the change in porosity with volumetric strain, the permeability is updated using the cubic relationship between porosity and permeability. Combined with the solid deformation model and the fluid seepage model, a seepage-stress field coupling model is established. S4. Use the virtual element method to solve the seepage-stress field coupling model, perform iterative calculations, and output the simulation results.
[0006] Preferably, step S1 specifically includes the following steps: S11. Based on the elastoplastic deformation characteristics of shale, an incremental constitutive equation is established using Hook's law to describe the relationship between effective stress and strain. The constitutive equation is as follows: ; in, Represents effective stress. Represents response, The tensor representing the elastic-plastic coefficient matrix; S12. Establish the stress balance equation for the spatial eight-node element, ensuring that the sum of the force vectors in all directions of the micro-element is zero; S13. Based on the small deformation assumption of continuum mechanics, establish the geometric equations: ; ; in, Represents the Laplace operator, Represents the stress tensor. Represents the displacement field components, Represents effective stress. Represents the Lamé coefficient, This represents the shear modulus.
[0007] Preferably, step S2 specifically includes the following steps: S21. Construct the control equation for the two-phase flow of oil and water in the fracture. The control equation is as follows: ; in, Represents the divergence operator. Represents the porosity of shale. Represents the pressure gradient. Represents the absolute permeability of the reservoir. Represents the saturation of oil-phase fluids in the reservoir. The density of the oil phase fluid. Represents the relative permeability of the oil phase. The viscosity represents the viscosity of the oil phase fluid. Representing the oil source and sink phases in the fracture system. Represents the saturation of aqueous fluids in the reservoir. The density of an aqueous fluid. Represents the relative permeability of the aqueous phase. Viscosity represents the viscosity of an aqueous fluid. The water source and sink phases represent the fracture system; S22. Considering stress sensitivity, an exponential equation representing the permeability of shale reservoirs is established, and its formula is: ; in, This represents the effective permeability after considering the stress-sensitive effect of shale reservoirs. Represents the initial penetration rate. This represents the stress sensitivity coefficient.
[0008] Preferably, in step S3, the formula for relating porosity changes to volumetric strain is: ; in, Represents the initial porosity of shale. Represents volumetric strain; The formula for updating permeability using the cubic relationship between porosity and permeability is: ; in, Represents shale permeability, This represents the initial permeability of shale.
[0009] Preferably, in step S3, by simultaneously solving the seepage equation and the solid mechanics equation, the linear elastic deformation equation for the fluid-structure interaction problem of shale oil is obtained, and its formula is as follows: ; in, Represents shale permeability, Represents boundary forces.
[0010] Preferably, step S4 specifically includes the following steps: S41. Divide fractured shale oil reservoirs into unstructured polyhedral grid units; S42. Decompose the motion of each mesh element into rigid body motion, constant strain motion and higher-order motion forms, where rigid body motion includes translation and rotation, and constant strain motion includes axial strain and shear strain. S43. Constructing Rigid Body Motion Projection Operators And constant strain motion projection operator The decomposition relationship of the unit displacement field is characterized by the projection matrix; S44. Introduce stabilizing terms to suppress oscillations in higher-order motion forms, construct the element stiffness matrix based on the principle of minimum energy, and then assemble the global stiffness matrix. S45. Combining stress boundary conditions and displacement boundary conditions, the displacement field descriptions for all elements are obtained: ; in, Represents the global stiffness matrix. The displacement matrix representing all elements. This represents the load matrix acting on the element; S46. Calculate the stress and strain fields using the geometric and constitutive equations of the solid deformation model; S47. Update the porosity and permeability of the shale reservoir based on the strain field, substitute the updated porosity and permeability into the flow control equation, and repeat steps S41-S47 for iterative calculation until the preset maximum time step is reached. Output the evolution law of shale porosity, permeability, fluid pressure and dynamic shale oil production capacity.
[0011] Preferably, in step S44, the element stiffness matrix is: ; in, Represents the volume of a mesh cell. The weight mapping matrix representing the constant strain mode, Represents the elastic tensor matrix. Represents a unit array, Represents the stable term. This represents the decomposition of the displacement field into rigid body displacement and constant strain motion, i.e. .
[0012] Therefore, the present invention employs the above-mentioned virtual element-based coupled simulation method for seepage-stress field in fractured shale oil, which has the following advantages: (1) In this invention, relying on the natural adaptability of the Virtual Element Method (VEM) to unstructured grids, it supports arbitrary polyhedral element division. It can efficiently characterize the complex multi-scale fracture network composed of primary natural fractures and secondary hydraulic fractures after shale fracturing without relying on structured grids. It can effectively capture the discontinuous displacement characteristics at the fractures, solve the adaptation problem of traditional numerical methods in the simulation of complex reservoir geometry and discrete fracture elements, and greatly improve the ability to characterize highly heterogeneous shale reservoirs.
[0013] (2) In this invention, the elastic deformation of porous media, stress-sensitive effect and seepage dynamic equation are integrated by equivalent continuous medium theory to construct a two-way coupled control system with high simulation accuracy and strong reliability of production capacity and reservoir dynamic prediction.
[0014] The technical solution of the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. Attached Figure Description
[0015] Figure 1 This is a schematic diagram of the process of the virtual element-based simulation method for fractured shale oil seepage-stress field coupling in this invention. Figure 2 This is a comparison and verification diagram between the simulation results of the present invention and the research results of Gudala in the embodiments; Figure 3The figure shows the simulation results of the vertical displacement of the model of the present invention and the study of Gudala in the embodiment; Figure 4 This is a comparison chart of the yields of the fluid-structure interaction model, the black oil model, and the weakened model in the embodiments. Detailed Implementation
[0016] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. The components of the embodiments of the present invention described and shown in the accompanying drawings can generally be arranged and designed in various different configurations. Specific model specifications need to be selected and determined according to the actual specifications of the device, etc. The specific selection calculation method adopts existing technology in the art, and therefore will not be described in detail.
[0017] Example like Figure 1 As shown, this invention provides a virtual element-based method for simulating the coupled flow-stress field of fractured shale oil, comprising the following steps: S1. Construct a solid deformation model of fractured shale oil reservoirs, specifically including the following steps: S11. Based on the elastoplastic deformation characteristics of shale, an incremental constitutive equation is established using Hook's law to describe the relationship between effective stress and strain. The constitutive equation is as follows: ; in, Represents effective stress. Represents response, The tensor representing the elastic-plastic coefficient matrix; S12. Establish the stress balance equation for the spatial eight-node element, ensuring that the sum of the force vectors in all directions of the micro-element is zero; S13. Based on the small deformation assumption of continuum mechanics, establish the geometric equations: ; ; in, Represents the Laplace operator, Represents the stress tensor. Represents the displacement field components, Represents effective stress. Represents the Lamé coefficient, This represents the shear modulus.
[0018] S2. Construct a fluid flow model for fractured shale oil reservoirs, specifically including the following steps: S21. Construct the governing equations for the two-phase oil-water flow in the fracture. The governing equations are: ; in, Represents the divergence operator. Represents the porosity of shale. Represents the pressure gradient. Represents the absolute permeability of the reservoir. Represents the saturation of oil-phase fluids in the reservoir. The density of the oil phase fluid. Represents the relative permeability of the oil phase. The viscosity represents the viscosity of the oil phase fluid. Representing the oil source and sink phases in the fracture system. Represents the saturation of aqueous fluids in the reservoir. The density of an aqueous fluid. Represents the relative permeability of the aqueous phase. Viscosity represents the viscosity of an aqueous fluid. The water source and sink phases represent the fracture system; S22. Considering stress sensitivity, an exponential equation representing the permeability of shale reservoirs is established, and its formula is: ; in, This represents the effective permeability after considering the stress-sensitive effect of shale reservoirs. Represents the initial penetration rate. This represents the stress sensitivity coefficient.
[0019] S3. By correlating the change in porosity with volumetric strain, the permeability is updated using the cubic relationship between porosity and permeability. Combined with the solid deformation model and the fluid seepage model, a seepage-stress field coupling model is established. The formula relating volumetric strain to porosity change is: ; in, Represents the initial porosity of shale. Represents volumetric strain; The formula for updating permeability using the cubic relationship between porosity and permeability is: ; in, Represents shale permeability, This represents the initial permeability of shale.
[0020] By simultaneously applying the seepage equation and the solid mechanics equation, the linear elastic deformation equation for the fluid-structure interaction problem of shale oil is obtained, and its formula is as follows: ; in, Represents shale permeability, Represents boundary forces.
[0021] S4. Solve the seepage-stress field coupling model using the virtual element method, perform iterative calculations, and output the simulation results. This includes the following steps: S41. Divide fractured shale oil reservoirs into unstructured polyhedral grid units; S42. Decompose the motion of each mesh element into rigid body motion, constant strain motion and higher-order motion forms, where rigid body motion includes translation and rotation, and constant strain motion includes axial strain and shear strain. For rigid body translational motion, take the unit vectors in the x, y, and z directions. , and Let represent the rigid body translation along the x, y, and z directions. Then, the rigid body translation is expressed as: ; For rigid body rotational motion, the displacement field of the rigid body rotation is expressed as: ; in, Represents the angular velocity vector. Represents the coordinate vector of the grid nodes. Represents the centroid coordinates of the grid cells; Rotations about the x-axis, y-axis, and z-axis can be expressed as follows: ; ; ; Where T represents the transpose matrix, and the subscripts (1), (2) and (3) represent the first, second and third components of the vector; S43. Constructing Rigid Body Motion Projection Operators And constant strain motion projection operator The decomposition relationship of the unit displacement field is characterized by the projection matrix; Rigid body motion can be achieved through the projection operator. Represented as: ; Constant strain motion can be achieved through the projection operator Represented as: ; S44. Introducing stabilizing terms to suppress oscillations in higher-order motion forms, constructing the element stiffness matrix based on the principle of minimum energy, and then assembling the global stiffness matrix, the element stiffness matrix is: ; in, Represents the volume of a mesh cell. The weight mapping matrix representing the constant strain mode, Represents the elastic tensor matrix. Represents a unit array, Represents the stable term. This represents the decomposition of the displacement field into rigid body displacement and constant strain motion, i.e. ; S45. Combining stress boundary conditions and displacement boundary conditions, the displacement field descriptions for all elements are obtained: ; in, Represents the global stiffness matrix. The displacement matrix representing all elements. This represents the load matrix acting on the element; S46. Calculate the stress and strain fields using the geometric and constitutive equations of the solid deformation model; S47. Update the porosity and permeability of the shale reservoir based on the strain field, substitute the updated porosity and permeability into the flow control equation, and repeat steps S41-S47 for iterative calculation until the preset maximum time step is reached. Output the evolution law of shale porosity, permeability, fluid pressure and dynamic shale oil production capacity.
[0022] In the model validation phase, this embodiment compares and analyzes the solution of the fluid-structure interaction stress field with existing research by Gudala et al. using a numerical example. The parameters used in the comparison validation are shown in Table 1. The boundary conditions of the model are top drainage, with other boundaries fixed. The model validation results are as follows: Figure 2 As shown, it can be observed that the solution results of this invention agree well with the research on Gudala; the simulation results of the vertical displacement of the model of this invention and Gudala are as follows. Figure 3 As shown, the vertical displacements of the two are well matched, with an average error of 1.8%.
[0023] Table 1. List of parameters used in the verification process comparing this embodiment with the Gudala study.
[0024] In terms of verifying seepage problems, by constructing benchmark cases with the same porosity, permeability, and fluid properties, the differences in daily oil and water production were compared and analyzed between the fluid-structure interaction model established in this invention, the black oil model, and the weakened model of this invention (which weakens the fluid-structure interaction effect). This verified the correctness of the model of this invention. Figure 4As shown, the difference in oil production between the weakened model and the black oil model is less than 2%, indicating that the model of this invention is comparable to the conventional black oil model after weakening the fluid-structure interaction effect. Because the model of this invention quantifies the stress-sensitive effect of reservoir properties under the coupling effect of rock deformation and fluid flow, it exhibits a continuous production suppression phenomenon compared to the black oil model, with a maximum production difference of up to 18.7%. This reflects the regulation of reservoir fluid seepage by pore pressure redistribution. The low production capacity of the fluid-structure interaction model stems from its detailed characterization of elastoplastic mechanisms such as stress redistribution and stress sensitivity, confirming the necessity of the fluid-structure interaction model for shale oil depletion development.
[0025] Therefore, this invention adopts the above-mentioned virtual element-based fractured shale oil seepage-stress field coupling simulation method. By modifying the seepage equation and the effective stress principle, a shale oil fluid-solid coupling model considering multi-scale pore structure and stress-sensitive effects is established. The simulation accuracy is high, and the reliability of production capacity and reservoir dynamic prediction is strong.
[0026] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit them. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the technical solutions of the present invention, and these modifications or equivalent substitutions cannot cause the modified technical solutions to deviate from the spirit and scope of the technical solutions of the present invention.
Claims
1. A virtual element-based simulation method for coupled flow-stress field in fractured shale oil, characterized in that: Includes the following steps: S1. Construct a solid deformation model for fractured shale oil reservoirs; S2. Construct a fluid flow model for fractured shale oil reservoirs; S3. By correlating the change in porosity with volumetric strain, the permeability is updated using the cubic relationship between porosity and permeability. Combined with the solid deformation model and the fluid seepage model, a seepage-stress field coupling model is established. S4. Use the virtual element method to solve the seepage-stress field coupling model, perform iterative calculations, and output the simulation results.
2. The virtual element-based simulation method for coupled flow-stress field in fractured shale oil as described in claim 1, characterized in that: Step S1 specifically includes the following steps: S11. Based on the elastoplastic deformation characteristics of shale, an incremental constitutive equation is established using Hook's law to describe the relationship between effective stress and strain. The constitutive equation is as follows: ; in, Represents effective stress. Represents response, The tensor representing the elastic-plastic coefficient matrix; S12. Establish the stress balance equation for the spatial eight-node element, ensuring that the sum of the force vectors in all directions of the micro-element is zero; S13. Based on the small deformation assumption of continuum mechanics, establish the geometric equations: ; ; in, Represents the Laplace operator, Represents the stress tensor. Represents the displacement field components, Represents effective stress. Represents the Lamé coefficient, This represents the shear modulus.
3. The virtual element-based coupled simulation method for seepage-stress field in fractured shale oil according to claim 2, characterized in that: Step S2 specifically includes the following steps: S21. Construct the control equation for the two-phase flow of oil and water in the fracture. The control equation is as follows: ; in, Represents the divergence operator. Represents the porosity of shale. Represents the pressure gradient. Represents the absolute permeability of the reservoir. Represents the saturation of oil-phase fluids in the reservoir. The density of the oil phase fluid. Represents the relative permeability of the oil phase. The viscosity represents the viscosity of the oil phase fluid. Representing the oil source and sink phases in the fracture system. Represents the saturation of aqueous fluids in the reservoir. The density of an aqueous fluid. Represents the relative permeability of the aqueous phase. Viscosity represents the viscosity of an aqueous fluid. The water source and sink phases represent the fracture system; S22. Considering stress sensitivity, an exponential equation representing the permeability of shale reservoirs is established, and its formula is: ; in, This represents the effective permeability after considering the stress-sensitive effect of shale reservoirs. Represents the initial penetration rate. This represents the stress sensitivity coefficient.
4. The virtual element-based coupled simulation method for seepage-stress field in fractured shale oil according to claim 3, characterized in that: In step S3, the formula for relating volumetric strain to porosity change is: ; in, Represents the initial porosity of shale. Represents volumetric strain; The formula for updating permeability using the cubic relationship between porosity and permeability is: ; in, Represents shale permeability, This represents the initial permeability of shale.
5. The virtual element-based simulation method for coupled flow-stress field in fractured shale oil as described in claim 4, characterized in that: In step S3, by simultaneously solving the seepage equation and the solid mechanics equation, the linear elastic deformation equation for the fluid-structure interaction problem of shale oil is obtained, and its formula is as follows: ; in, Represents shale permeability, Represents boundary forces.
6. The virtual element-based coupled simulation method for seepage-stress field in fractured shale oil according to claim 5, characterized in that: Step S4 specifically includes the following steps: S41. Divide fractured shale oil reservoirs into unstructured polyhedral grid units; S42. Decompose the motion of each mesh element into rigid body motion, constant strain motion and higher-order motion forms, where rigid body motion includes translation and rotation, and constant strain motion includes axial strain and shear strain. S43. Constructing the rigid body motion projection operator And constant strain motion projection operator The decomposition relationship of the unit displacement field is characterized by the projection matrix; S44. Introduce stabilizing terms to suppress oscillations in higher-order motion forms, construct the element stiffness matrix based on the principle of minimum energy, and then assemble the global stiffness matrix. S45. Combining stress boundary conditions and displacement boundary conditions, the displacement field descriptions for all elements are obtained: ; in, Represents the global stiffness matrix. The displacement matrix representing all elements. This represents the load matrix acting on the element; S46. Calculate the stress and strain fields using the geometric and constitutive equations of the solid deformation model; S47. Update the porosity and permeability of the shale reservoir based on the strain field, substitute the updated porosity and permeability into the flow control equation, and repeat steps S41-S47 for iterative calculation until the preset maximum time step is reached. Output the evolution law of shale porosity, permeability, fluid pressure and dynamic shale oil production capacity.
7. The virtual element-based coupled simulation method for seepage-stress field in fractured shale oil according to claim 6, characterized in that: In step S44, the element stiffness matrix is: ; in, Represents the volume of a mesh cell. The weight mapping matrix representing the constant strain mode, Represents the elastic tensor matrix. Represents a unit array, Represents the stable term. This represents the decomposition of the displacement field into rigid body displacement and constant strain motion, i.e. .