Method for simulating seepage pressure field of crack type porous medium and related device

Abstract

The present application belongs to the technical field of seepage simulation in fractured porous media, and discloses a seepage pressure field simulation method for fractured porous media and related devices, which comprises the following steps: constructing a seepage control equation for fractured porous media; integrating the seepage control equation in a control volume to obtain an integral form control equation; dividing embedded discrete fracture grids in the fractured porous media, identifying the interaction areas between the embedded discrete fracture grids, and establishing a pressure distribution function inside the interaction areas; according to the pressure distribution function, solving the discrete expressions of each integral term in the integral form control equation in the interaction areas to obtain the discrete expressions of each integral term; substituting the discrete expressions of each integral term into the integral form control equation to obtain a discrete equation set and solve it to obtain the seepage pressure field simulation value of the fractured porous media. The present application can simulate the situation when the fractures are blocked, and can improve the calculation accuracy of the bedrock and fracture mass exchange term.

Description

Technical Field

[0001] This invention belongs to the field of seepage simulation technology in fractured porous media, and relates to a method and related apparatus for simulating the seepage pressure field in fractured porous media. Background Technology

[0002] Accurate and efficient simulation of seepage flow in fractured porous media is crucial for underground engineering projects such as oil and gas recovery and geothermal extraction. Fractures present two challenges to seepage simulation in porous media: firstly, the pore size of fractures is very small, but the length of fractures can span the entire computational domain, making it difficult to characterize fractures during mesh generation; secondly, the permeability of fractures differs greatly from that of bedrock, spanning several orders of magnitude. When the permeability of fractures is greater than that of bedrock, they act as conductive channels, while when the permeability of fractures is less than that of bedrock, they act as sealing barriers.

[0003] To address Challenge 1, scholars have recently proposed embedded discrete fracture models, which use embedded discrete fracture meshes to divide the computational domain. This means that the bedrock and fracture regions use independent meshes, significantly reducing the difficulty of mesh generation. In the classic embedded discrete fracture model, the mass exchange between bedrock and fractures is introduced into the discrete equations through a source term. This source term is derived from the assumption that the pressure in the bedrock surrounding the fracture is linearly related to the distance to the fracture. This assumption only holds true when the fracture permeability is much greater than that of the bedrock; otherwise, it deviates significantly from reality. Therefore, the classic embedded discrete fracture model has low accuracy in calculating the mass exchange between bedrock and fractures and is unsuitable for the second challenge, where the fracture acts as a sealing barrier. Therefore, there is an urgent need to develop a new seepage field simulation method based on embedded discrete fracture meshes, applicable to cases where the fracture is sealed, and to improve the computational accuracy of the bedrock-fracture mass exchange term. Summary of the Invention

[0004] The purpose of this invention is to overcome the shortcomings of the prior art and provide a method and related apparatus for simulating the seepage pressure field of fractured porous media.

[0005] To achieve the above objectives, the present invention employs the following technical solution:

[0006] In a first aspect, the present invention provides a method for simulating the seepage pressure field of fractured porous media, comprising:

[0007] The pressure and flux on both sides of the fracture in the fractured porous medium are modeled separately, and the seepage control equation of the fractured porous medium is constructed.

[0008] Integrating the seepage control equations within the control volume yields the integral form of the control equations.

[0009] Embedded discrete crack meshes are generated within fractured porous media, and the interaction regions between the embedded discrete crack meshes are identified. Based on the characteristics that the introduced variable values ​​and pressure gradients are discontinuous, a pressure distribution function is established within the interaction region.

[0010] Based on the pressure distribution function, the discrete expressions of each integral term in the integral form of the control equation in the interaction region are solved to obtain the discrete expressions of each integral term.

[0011] Substituting the discrete expressions of each integral term into the integral form of the governing equations, we obtain a set of discrete equations and solve them to obtain the simulated values ​​of the seepage pressure field of the fractured porous medium.

[0012] Optionally, the flow control equation for constructing the fractured porous medium includes:

[0013] Establish the mass conservation equation for bedrock in fractured porous media:

[0014]

[0015] Where the superscript m indicates bedrock; K m p is the bedrock permeability tensor. m For bedrock pressure; Ω m For bedrock region; q f-m This refers to the mass exchange term from the fracture to the bedrock.

[0016] Establish the mass conservation equation in the fractures of fractured porous media:

[0017]

[0018] Wherein, the subscript τ represents the axial direction of the crack; The gradient is along the crack axis; the superscript f indicates the crack. d represents the permeability along the fracture axis. f P is the crack diameter; f The average pressure across the crack section, i.e. γ represents the crack wall surface, and γ1 and γ2 represent the two side walls of the crack; u m denoted as Darcy velocity within the bedrock; n is the unit vector of the fracture normal, its direction being from the fracture towards the bedrock relative to the γ2 surface;

[0019] Construct the relationship between the sum of fluxes through the crack walls and the wall pressure:

[0020]

[0021] in, The permeability is denoted by the fracture normal direction.

[0022] Construct the relationship between the difference in flux through the fracture wall and the wall pressure and the pressure inside the fracture:

[0023]

[0024] Set the boundary conditions as follows:

[0025]

[0026] Among them, Γ m and Γ f These represent the boundaries of the bedrock and the fractured area, respectively. For boundary pressure; For the boundary Darcy velocity;

[0027] Equations (1) to (5) are used as the seepage control equations for fractured porous media.

[0028] Optionally, integrating the seepage control equations within the control volume yields integral form control equations including:

[0029] Integrating equation (1) over the control volume, the left-hand side of equation (1) can be obtained using the divergence theorem:

[0030]

[0031] Where Ω′ represents the control volume region; n A dA is the normal vector of the control volume interface; k is the differential control volume interface; xx k xy and k yy These are components of the bedrock permeability tensor; and denoted as the partial derivatives of the bedrock pressure in the x and y directions; the subscripts e′, w′, n′, and s′ denote the control volume interfaces;

[0032] The right-hand side of equation (1) equals the flux flowing into the bedrock region at the fracture interface minus the flux flowing out of the bedrock region:

[0033]

[0034] Where γ′ represents the crack segment within the control volume; γ′1 and γ′2 are the segments on both sides of the crack wall; dl is the differential crack segment; n x and n y These are the components of vector n;

[0035] Integrating equation (2) along the crack segment within the control volume, the divergence theorem can be applied to the left-hand side of equation (2) to obtain:

[0036]

[0037] The integral of the right-hand side of equation (2) is expressed as:

[0038]

[0039] Integrating equation (3) along the crack segment within the control volume, the left-hand side of equation (3) is expressed as:

[0040]

[0041] The integral of the right-hand side of equation (3) is expressed as:

[0042]

[0043] Integrating equation (4) along the crack segment within the control volume, it can be expressed as:

[0044]

[0045] The left-hand side of equation (12) can be represented by equation (7);

[0046] By rearranging equations (6) to (12), we can obtain the integral form of the governing equations as follows:

[0047]

[0048]

[0049]

[0050]

[0051] Optionally, the step of dividing the fractured porous medium into embedded discrete fracture grids and identifying the interaction regions between the embedded discrete fracture grids includes:

[0052] The bedrock region of the fractured porous medium is discretized using a structured rectangular grid. The location of the control volume interface is determined, with the element node being the grid center. Interaction regions are delineated, offset from the control volume by a quarter of the grid. Each interaction region is rectangular, and its vertices are composed of element nodes. Each control volume is covered by four interaction regions, and the interface of the control volume is divided into eight half-interfaces, with each interaction region covering two half-interfaces. The fractures are divided using a one-dimensional grid. A fracture segment covered by a single interaction region is a fracture grid, and the interface of the fracture grid is located at the intersection of the fracture segment and the boundary of the interaction region, with the grid center being an element node. An additional degree of freedom is added to each fracture grid on both sides.

[0053] Optionally, the establishment of the pressure distribution function within the interaction region includes:

[0054] Integral terms in the inductive integral form of the governing equations:

[0055]

[0056] Where i∈{x,y}, k∈{e,w,n,s} and j∈{1,2}, e, w, n, and s represent the four half-interfaces in the interactive area; D i,k L is the integral of the partial derivative of bedrock pressure along the control volume interface; i,j M is the integral of the pressure partial derivative along one side of the crack; j This is the integral of the pressure value along one side of the crack;

[0057] The pressure distribution function within the interaction region is defined as follows:

[0058]

[0059] in, For the pressure at the four bedrock unit nodes; For the crack element, these are two additional degrees of freedom; their value is always 0 if no crack passes through the interaction region. h (x) is the interpolation shape function of the rectangular element, and its expression in the local coordinate system with the lower left corner element node of the interaction region as the origin is:

[0060]

[0061] Where δx and δy are the width and height of the interaction area, respectively;

[0062] In equation (18), B h (x), C h (x) and These are three step functions used to characterize the discontinuity of variable values ​​on both sides of the crack, and their expressions are as follows:

[0063]

[0064] Here, sgn(x) is the signed distance function for the crack within the interaction region, with a value of 1 on one side of the crack and a value of -1 on the other side; x h These are the coordinates of the bedrock unit nodes; For in crack γ j The coordinates of any point on the side.

[0065] Optionally, the step of solving for the discrete expressions of each integral term in the integral form of the control equation within the interaction region based on the pressure distribution function includes:

[0066] Solve for the nodal pressure of the bedrock element with respect to the interface integral D i,k The contributions are as follows:

[0067]

[0068]

[0069] in, It is the pressure vector of the bedrock unit node in the interaction region φ;

[0070] Solve for the integral L of the nodal pressure of the bedrock element along the fracture segment. i,j and M j Contributions:

[0071]

[0072]

[0073]

[0074] Where (fx, fy) are the coordinates of the crack segment in the local coordinate system;

[0075] Solve for the interface integral D of the additional degrees of freedom of the crack. i,k Contributions:

[0076]

[0077]

[0078]

[0079] in, , which is the additional degree of freedom vector for the crack segment in the interaction region φ;

[0080] Solve for the additional degree of freedom of the crack by integrating along the crack segment L. i,j and M j Contributions:

[0081]

[0082]

[0083]

[0084] The discrete expressions for each integral term are obtained using the following formula:

[0085] and

[0086] Optionally, the step of substituting the discrete expressions of each integral term into the integral form of the governing equations to obtain a discrete set of equations and solving it to obtain the simulated values ​​of the seepage pressure field of the fractured porous medium includes:

[0087] For equation (13), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0088]

[0089] In this matrix, rows 1 to 4 correspond to the control volumes of unit nodes x1 to x4, respectively. and The superscript number indicates the value of the integral along the crack segment in the relevant control volume, satisfying... and

[0090] The coefficient matrix in equation (32) can be simplified to the following form:

[0091]

[0092] Rearranging equation (33) further, we get:

[0093]

[0094] By superimposing all interactive regions, the matrix of equation (34) is integrated to obtain:

[0095]

[0096] Where, φ tot p represents the number of interactive areas. m p represents the pressure vector at all bedrock unit nodes. fa Add a degree-of-freedom vector to all cracks;

[0097] For equation (14), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0098]

[0099] in, φ represents the pressure vector of the crack element and its adjacent crack element nodes in the interactive region φ. The discrete matrix for axial flow in the crack is discretized using central difference.

[0100] Equation (36) can be simplified as follows:

[0101]

[0102] Rearranging equation (37) further, we get:

[0103]

[0104] By superimposing all interactive regions, the matrix of equation (37) is integrated to obtain:

[0105]

[0106] Where, p fFor the pressure vector of all crack element nodes;

[0107] For equation (15), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0108]

[0109] Equation (40) can be simplified as follows:

[0110]

[0111] Rearranging equation (41) further, we get:

[0112]

[0113] By superimposing all interactive regions, the matrix of equation (41) is integrated to obtain:

[0114]

[0115] For equation (16), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0116]

[0117] Equation (44) can be simplified as follows:

[0118]

[0119] Rearranging equation (45) further, we get:

[0120]

[0121] By superimposing all interactive regions, the matrix of equation (46) is integrated to obtain:

[0122]

[0123] By assembling equations (35), (39), (43), and (47), we obtain the overall discrete equation system:

[0124]

[0125] Solving equation (48) yields the simulated values ​​of the seepage pressure field in the fractured porous medium.

[0126] In a second aspect, the present invention provides a system for simulating the seepage pressure field of fractured porous media, comprising:

[0127] The module is used to model the pressure and flux on both sides of the fracture in the fractured porous medium, and to construct the seepage control equation for the fractured porous medium.

[0128] The integral module is used to integrate the seepage control equations within the control volume to obtain the integral form of the control equations.

[0129] The mesh embedding module is used to divide embedded discrete crack meshes in cracked porous media, identify the interaction regions between embedded discrete crack meshes, and establish the pressure distribution function inside the interaction region based on the characteristics that the introduced variable values ​​and pressure gradients are discontinuous.

[0130] The discrete module is used to solve the discrete expressions of each integral term in the integral form of the control equation within the interaction region based on the pressure distribution function, and to obtain the discrete expressions of each integral term.

[0131] The solver module is used to substitute the discrete expressions of each integral term into the integral form of the governing equations to obtain a discrete set of equations and solve them to obtain the simulated values ​​of the seepage pressure field of the fractured porous medium.

[0132] In a third aspect, the present invention provides a computer device, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the steps of the above-described method for simulating the seepage pressure field of fractured porous media.

[0133] In a fourth aspect, the present invention provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the steps of the above-described method for simulating the seepage pressure field of a fractured porous medium.

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

[0135] This invention presents a method for simulating the seepage pressure field in fractured porous media. It involves constructing a seepage control equation for the fractured porous media, integrating this equation within a control volume to obtain an integral form control equation, further dividing the fractured porous media into embedded discrete fracture grids, identifying interaction regions between these grids, and establishing a pressure distribution function within each interaction region. Based on this pressure distribution function, the method solves for the discrete expressions of each integral term in the integral form control equation within the interaction region. Substituting these discrete expressions into the integral form control equation yields a set of discrete equations, which are then solved to obtain the simulated seepage pressure field values ​​for the fractured porous media. By modeling the pressure and flux on both sides of the fracture separately in the control equations, and introducing discontinuous variable values ​​and gradients into the pressure distribution function, the method accurately reflects the reality of discontinuous pressure and flux on both sides of a sealed fracture, thus simulating a sealed fracture. Furthermore, by supplementing the equations with two relationships between pressure and flux on both sides and inside the fracture, the method can more accurately describe the interaction between the bedrock and the fracture, thereby improving the calculation accuracy of the bedrock-fracture mass exchange term. Meanwhile, based on the embedded discrete fracture mesh discretization method, it has a small memory footprint and low computation time, which can provide more accurate prediction and guidance for the seepage state of underground reservoirs in the oil and gas extraction and geothermal development industries. Attached Figure Description

[0136] Figure 1 This is a flowchart illustrating the seepage pressure field simulation method for fractured porous media according to an embodiment of the present invention.

[0137] Figure 2 This is a schematic diagram illustrating the parameter identification of the control equations in an embodiment of the present invention.

[0138] Figure 3 Figure (a) is a schematic diagram of the embedded discrete fracture mesh and interactive region division in an embodiment of the present invention, and Figure (b) is a schematic diagram of the mesh division of the bedrock and fracture regions, and Figure (b) is a schematic diagram of the additional degrees of freedom of the fracture mesh in the interactive region.

[0139] Figure 4 This is a schematic diagram of the computational region information according to an embodiment of the present invention.

[0140] Figure 5 This is a schematic diagram of the computational grid according to an embodiment of the present invention.

[0141] Figure 6 Figure 1 is a schematic diagram of the seepage pressure field distribution when the crack is conductive according to an embodiment of the present invention; wherein, Figure (a) is the pressure field obtained by the method of the present invention, Figure (b) is the pressure field obtained by the classical EDFM method, and Figure (c) is the reference solution.

[0142] Figure 7This is a schematic diagram of the relative error distribution in the bedrock region when the fracture is conductive according to an embodiment of the present invention; wherein, Figure (a) is the error distribution of the method of the present invention, and Figure (b) is the error distribution of the classical EDFM method.

[0143] Figure 8 Figure (a) is a schematic diagram of the pressure distribution inside the crack when the crack is conducting, and Figure (b) is a schematic diagram of the pressure distribution inside the crack when crack f1 is conducting.

[0144] Figure 9 This is a schematic diagram showing the variation of the average relative error between the bedrock and the fracture region with the number of grids when the fracture is connected according to an embodiment of the present invention; wherein, Figure (a) shows the error in the bedrock region and Figure (b) shows the error in the fracture region.

[0145] Figure 10 Figure 1 is a schematic diagram of the seepage pressure field distribution during crack sealing according to an embodiment of the present invention; wherein, Figure (a) is the pressure field obtained by the method of the present invention, Figure (b) is the pressure field obtained by the classical EDFM method, and Figure (c) is the reference solution.

[0146] Figure 11 This is a schematic diagram showing the variation of the average relative error between the bedrock and the fracture region with the number of grids during fracture sealing according to an embodiment of the present invention; wherein, Figure (a) shows the error in the bedrock region and Figure (b) shows the error in the fracture region. Detailed Implementation

[0147] To enable those skilled in the art to better understand the present invention, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of the present invention.

[0148] It should be noted that the terms "first," "second," etc., in the specification, claims, and accompanying drawings of this invention are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that such data can be interchanged where appropriate so that the embodiments of the invention described herein can be implemented in orders other than those illustrated or described herein. Furthermore, the terms "comprising" and "having," and any variations thereof, are intended to cover a non-exclusive inclusion; for example, a process, method, system, product, or apparatus that comprises a series of steps or units is not necessarily limited to those steps or units explicitly listed, but may include other steps or units not explicitly listed or inherent to such processes, methods, products, or apparatus.

[0149] The present invention will now be described in further detail with reference to the accompanying drawings:

[0150] See Figure 1 In one embodiment of the present invention, a method for simulating the seepage pressure field of fractured porous media is provided. This method can simulate the situation where the fractures are blocked and can improve the calculation accuracy of the mass exchange term between bedrock and fractures. Specifically, the method for simulating the seepage pressure field of fractured porous media includes the following steps:

[0151] S1: Model the pressure and flux on both sides of the fracture in the fractured porous medium and construct the seepage control equation for the fractured porous medium.

[0152] S2: Integrate the seepage control equation within the control volume to obtain the integral form of the control equation.

[0153] S3: In the fractured porous medium, an embedded discrete crack mesh is defined, and the interaction region between the embedded discrete crack meshes is identified. Based on the characteristics that the introduced variable values ​​and pressure gradients are discontinuous, a pressure distribution function is established within the interaction region.

[0154] S4: Based on the pressure distribution function, solve for the discrete expressions of each integral term in the integral form of the control equation within the interaction region to obtain the discrete expressions of each integral term.

[0155] S5: Substitute the discrete expressions of each integral term into the integral form of the governing equations to obtain a discrete set of equations and solve them to obtain the simulated values ​​of the seepage pressure field of the fractured porous medium.

[0156] In summary, the method for simulating the seepage pressure field of fractured porous media in this invention constructs the seepage control equations for the fractured porous media, integrates these equations within a control volume to obtain integral control equations, further divides the fractured porous media into embedded discrete fracture grids, identifies the interaction regions between these grids, establishes pressure distribution functions within these interaction regions, and solves for the discrete expressions of each integral term in the integral control equations within the interaction regions based on these pressure distribution functions. Substituting these discrete expressions into the integral control equations yields a set of discrete equations, which are then solved to obtain the simulated seepage pressure field values ​​for the fractured porous media. By modeling the pressure and flux on both sides of the fracture separately in the control equations, and introducing the discontinuous nature of variable values ​​and gradients in the pressure distribution function, the method accurately reflects the reality of discontinuous pressure and flux on both sides of a sealed fracture, thus simulating a sealed fracture. Furthermore, by supplementing the equations with two additional equations relating pressure and flux on both sides and inside the fracture, the method can more accurately describe the interaction between the bedrock and the fracture, thereby improving the calculation accuracy of the bedrock-fracture mass exchange term. Meanwhile, based on the embedded discrete fracture mesh discretization method, it has a small memory footprint and low computation time, which can provide more accurate prediction and guidance for the seepage state of underground reservoirs in the oil and gas extraction and geothermal development industries.

[0157] In one possible implementation, the flow control equation for constructing the fractured porous medium includes:

[0158] See Figure 2 Based on the parameter identifiers, establish the mass conservation equation for bedrock in fractured porous media:

[0159]

[0160] Where the superscript m indicates bedrock; K m p is the bedrock permeability tensor. m For bedrock pressure; Ω m For bedrock region; q f-m This is the mass exchange term from the fracture to the bedrock.

[0161] Establish the mass conservation equation in the fractures of fractured porous media:

[0162]

[0163] Wherein, the subscript τ represents the axial direction of the crack; The gradient is along the crack axis; the superscript f indicates the crack. d represents the permeability along the fracture axis. f P is the crack diameter; f The average pressure across the crack section, i.e. γ represents the crack wall surface, and γ1 and γ2 represent the two side walls of the crack; u m denoted as Darcy velocity within the bedrock; n is the unit vector of the fracture normal, its direction being relative to the γ2 surface from the fracture towards the bedrock.

[0164] Construct the relationship between the sum of fluxes through the crack walls and the wall pressure:

[0165]

[0166] in, The permeability is denoted by the normal direction of the crack.

[0167] Construct the relationship between the difference in flux through the fracture wall and the wall pressure and the pressure inside the fracture:

[0168]

[0169] Set the boundary conditions as follows:

[0170]

[0171] Among them, Γ m and Γ f These represent the boundaries of the bedrock and the fractured area, respectively. For boundary pressure; For the boundary Darcy velocity.

[0172] Equations (1) to (5) are used as the seepage control equations for fractured porous media.

[0173] Integrating the above governing equations within the control volume yields integral-form governing equations. Specifically, integrating the seepage governing equations within the control volume yields integral-form governing equations including:

[0174] Integrating equation (1) over the control volume, the left-hand side of equation (1) can be obtained using the divergence theorem:

[0175]

[0176] Where Ω′ represents the control volume region; n A dA is the normal vector of the control volume interface; k is the differential control volume interface; xx k xy and k yy These are components of the bedrock permeability tensor; and denoted as the partial derivatives of the bedrock pressure in the x and y directions; the subscripts e′, w′, n′, and s′ denote the control volume interfaces.

[0177] The right-hand side of equation (1) equals the flux flowing into the bedrock region at the fracture interface minus the flux flowing out of the bedrock region:

[0178]

[0179] Where γ′ represents the crack segment within the control volume; γ′1 and γ′2 are the segments on both sides of the crack wall; dl is the differential crack segment; n x and n y These are the components of vector n.

[0180] Integrating equation (2) along the crack segment within the control volume, the divergence theorem can be applied to the left-hand side of equation (2) to obtain:

[0181]

[0182] The integral of the right-hand side of equation (2) is expressed as:

[0183]

[0184] Integrating equation (3) along the crack segment within the control volume, the left-hand side of equation (3) is expressed as:

[0185]

[0186] The integral of the right-hand side of equation (3) is expressed as:

[0187]

[0188] Integrating equation (4) along the crack segment within the control volume, it can be expressed as:

[0189]

[0190] The left-hand side of equation (12) can be represented by equation (7).

[0191] By rearranging equations (6) to (12), we can obtain the integral form of the governing equations as follows:

[0192]

[0193]

[0194]

[0195]

[0196] The advantage of using integral form of the control equations is that the discrete equations can be guaranteed to satisfy the conservation property. Equations (13) to (16) are applicable to the case where the control volume is traversed by a crack.

[0197] The process of dividing the fractured porous medium into embedded discrete fracture grids and identifying the interaction regions between these grids includes: discretizing the bedrock region using a structured rectangular grid; firstly, determining the location of the bedrock control volume interface, with the grid center as the element node; then defining the interaction regions, which are offset from the control volume by a quarter of the grid. Each interaction region is still rectangular, and the vertices of the interaction region are composed of element nodes, such as... Figure 3 As shown in (a), the discreteness of the integral form of the control equations is performed in the interaction regions. Each control volume is covered by four interaction regions, and the interface of the control volume is divided into eight half-interfaces, with each interaction region covering two half-interfaces. Cracks are meshed using a one-dimensional grid, with the meshing location depending on the interaction region. That is, a crack segment covered by a single interaction region constitutes a crack mesh, and the interface of the crack mesh is at the intersection of the crack segment and the boundary of the interaction region. The element nodes are at the center of the mesh. Each crack mesh has an additional degree of freedom added to each of its two sides, as shown in (a). Figure 3 As shown in (b), it is used to construct the pressure distribution shape function within the interactive region.

[0198] The pressure distribution function established within the interactive region includes:

[0199] Integral terms in the inductive integral form of the governing equations:

[0200]

[0201] Where i∈{x,y}, k∈{e,w,n,s} and j∈{1,2}, e, w, n, and s represent the four half-interfaces in the interactive area, such as Figure 3 As shown; D i,k L is the integral of the partial derivative of bedrock pressure along the control volume interface; i,j M is the integral of the pressure partial derivative along one side of the crack; j This is the integral of the pressure value along one side of the crack.

[0202] The integral terms in equation (17) are discrete in the interaction region, meaning the purpose is to replace these integral terms with a combination of the values ​​of the variables to be determined (unit nodal pressure, crack additional degrees of freedom). Therefore, the pressure distribution function in the interaction region is defined based on the inductive integral terms as follows:

[0203]

[0204] in, For the pressure at the four bedrock unit nodes; For the crack element, these are two additional degrees of freedom; their value is always 0 if no crack passes through the interaction region. h (x) is the interpolation shape function of the rectangular element, and its expression in the local coordinate system with the lower left corner element node of the interaction region as the origin is:

[0205]

[0206] Where δx and δy are the width and height of the interaction area, respectively.

[0207] In equation (18), B h (x), C h (x) and C j fa (x) are three step functions used to characterize the discontinuity of variable values ​​on both sides of the crack, and their expressions are:

[0208]

[0209] Here, sgn(x) is the signed distance function for the crack within the interaction region, with a value of 1 on one side of the crack and a value of -1 on the other side; x h These are the coordinates of the bedrock unit nodes; To in the crack γ j The coordinates of any point on the side.

[0210] The step caused by equation (20) results in the pressure distribution shape function being discontinuous in both variable values ​​and gradient on both sides of the crack.

[0211] The discrete expressions of each integral term in the integral form of the control equations in the interaction region are obtained by solving the pressure distribution function, resulting in the following discrete expressions for each integral term:

[0212] Substituting equation (18) into equation (17), we can obtain the discrete expressions for the three integral terms.

[0213] Solve for the nodal pressure of the bedrock element with respect to the interface integral D i,k The contributions are as follows:

[0214]

[0215]

[0216] in, It is the pressure vector of the bedrock unit node in the interaction region φ;

[0217] Solve for the integral L of the nodal pressure of the bedrock element along the fracture segment. i,j and M j Contributions:

[0218]

[0219]

[0220]

[0221] Where (fx, fy) are the coordinates of the crack segment in the local coordinate system;

[0222] Solve for the interface integral D of the additional degrees of freedom of the crack. i,k Contributions:

[0223]

[0224]

[0225]

[0226] in, Let L be the additional degree of freedom vector of the crack segment in the interaction region φ; solve for the additional degree of freedom of the crack by integrating L along the crack segment. i,j and M j Contributions:

[0227]

[0228]

[0229]

[0230] Finally, the discrete expression for the integral term is the superposition of the contributions from nodal pressure and additional degrees of freedom of fractures in the bedrock element. The discrete expression for each integral term is obtained through the following equation:

[0231] and

[0232] The process of substituting the discrete expressions of each integral term into the integral form of the governing equations to obtain a set of discrete equations and solving them yields the simulated values ​​of the seepage pressure field in the fractured porous medium, including:

[0233] According to equations (21) to (31), the discrete results of the integral terms in the interaction region φ can be obtained, and then they can be further allocated to the control volume to which the integral terms belong.

[0234] For equation (13), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0235]

[0236] In this matrix, rows 1 to 4 correspond to the control volumes of unit nodes x1 to x4, respectively. and The superscript number indicates the value of the integral along the crack segment in the relevant control volume, satisfying... and

[0237] The coefficient matrix in equation (32) can be simplified to the following form:

[0238]

[0239] Rearranging equation (33) further, we get:

[0240]

[0241] By superimposing all interactive regions, the matrix of equation (34) is integrated to obtain:

[0242]

[0243] Where, φ tot p represents the number of interactive areas. m p represents the pressure vector at all bedrock unit nodes. fa Add a degree-of-freedom vector to all cracks.

[0244] For equation (14), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0245]

[0246] in, φ represents the pressure vector of the crack element and its adjacent crack element nodes in the interactive region φ. The discrete matrix for axial flow in the crack is discretized using central difference.

[0247] Equation (36) can be simplified as follows:

[0248]

[0249] Rearranging equation (37) further, we get:

[0250]

[0251] By superimposing all interactive regions, the matrix of equation (37) is integrated to obtain:

[0252]

[0253] Where, p f This represents the pressure vector at all crack element nodes.

[0254] For equation (15), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0255]

[0256] Equation (40) can be simplified as follows:

[0257]

[0258] Rearranging equation (41) further, we get:

[0259]

[0260] By superimposing all interactive regions, the matrix of equation (41) is integrated to obtain:

[0261]

[0262] For equation (16), substituting the discrete expressions of each integral term into the integral form of the control equation, we get:

[0263]

[0264] Equation (44) can be simplified as follows:

[0265]

[0266] Rearranging equation (45) further, we get:

[0267]

[0268] By superimposing all interactive regions, the matrix of equation (46) is integrated to obtain:

[0269]

[0270] By assembling equations (35), (39), (43), and (47), we obtain the overall discrete equation system:

[0271]

[0272] Solving equation (48) yields the simulated values ​​of the seepage pressure field in the fractured porous medium.

[0273] In one possible implementation, a porous medium with a computational area of ​​10m × 10m is used, containing two parallel, conductive fractures with a fracture diameter of 1cm. The bedrock permeability is 1 Darcy, and the fracture permeability is 10. 4 Darcy's example. The boundary conditions are calculated as follows: Figure 4 As shown, the top and bottom are closed interfaces, and the left and right sides are constant pressure conditions. The mesh is generated as follows. Figure 5 As shown, the bedrock mesh is 50×50, and the fracture mesh consists of 111 elements. The numerical solution of the discrete fracture model under a fine mesh is used as the reference solution, and the relative error is defined as follows:

[0274]

[0275]

[0276] Among them, err and These represent the relative error of a single grid and the average relative error of the entire computational domain, respectively; the superscripts m and f represent the bedrock and fracture regions, respectively; p re f The reference solution is interpolated onto the coarse grid; p coarse For the solution on the coarse grid; N m N represents the number of bedrock units. f This represents the number of crack elements.

[0277] See Figure 6 This paper demonstrates the pressure field distribution under the method of this invention, the classical embedded discrete fracture model (EDFM), and the reference solution. It can be seen that both the method of this invention and the classical EDFM method can calculate a reasonable pressure field. See also... Figure 7 This paper illustrates the relative error distribution of the method of this invention and the classical EDFM method within the bedrock region. The pressure error distribution area of ​​the method of this invention is significantly smaller than that of the classical EDFM method, indicating that the bedrock pressure calculation accuracy of the method of this invention is higher. See also Figure 8 The pressure distribution within the fracture is shown, demonstrating that the calculation results of the method described in this invention agree better with the reference solution. Since this method and classical EDFM use the same fracture axial calculation method, the difference lies in the mass exchange between the bedrock and the fracture. Therefore, it can be concluded that the method of this invention has higher accuracy than classical EDFM in calculating the mass exchange between the bedrock and the fracture. See also... Figure 9 The paper presents a convergence analysis of the method of the present invention and the classical EDFM with mesh refinement. It can be seen that, for various mesh numbers, the average relative error of the method of the present invention is smaller than that of the classical EDFM in both the bedrock region and the fracture region, which further illustrates that the method of the present invention has higher calculation accuracy in all aspects.

[0278] In one possible implementation, the fissure serves as a sealing barrier with a permeability of 10. -4 Darcy, other calculation parameters are the same as in the above implementation method. See also Figure 10 This paper illustrates the pressure field distribution under the proposed method, the classical EDFM method, and the reference solution. It shows that the proposed method can calculate a reasonable pressure field, with a step pressure value on both sides of the sealed fracture, while the classical EDFM method fails to capture the effect of the sealed fracture. See also... Figure 11 The paper presents the convergence analysis of the method of the present invention and the classical EDFM with mesh refinement. It can be seen that the average relative error of the method of the present invention decreases in the bedrock and fracture regions with mesh refinement, indicating that the method of the present invention has a mesh-independent solution, while the average relative error of the classical EDFM remains high with mesh refinement.

[0279] Overall, the method of the present invention can reasonably describe the seepage pressure field when sealing cracks.

[0280] The following are embodiments of the apparatus of the present invention, which can be used to execute embodiments of the method of the present invention. For details not disclosed in the apparatus embodiments, please refer to the embodiments of the method of the present invention.

[0281] In another embodiment of the present invention, a seepage pressure field simulation system for fractured porous media is provided, which can be used to implement the above-mentioned seepage pressure field simulation method for fractured porous media. Specifically, the seepage pressure field simulation system for fractured porous media includes a construction module, an integration module, a mesh embedding module, a discretization module, and a solution module.

[0282] The system comprises the following modules: a construction module for constructing the seepage control equations for fractured porous media; an integration module for integrating the seepage control equations within the control volume to obtain integral form control equations; a mesh embedding module for dividing the fractured porous media into embedded discrete fracture meshes, identifying the interaction regions between the embedded discrete fracture meshes, and establishing the pressure distribution function within the interaction regions; a discretization module for solving the discrete expressions of each integral term in the integral form control equations within the interaction regions based on the pressure distribution function; and a solution module for substituting the discrete expressions of each integral term into the integral form control equations to obtain a set of discrete equations and solving them to obtain the simulated values ​​of the seepage pressure field of the fractured porous media.

[0283] All relevant content of each step involved in the aforementioned embodiments of the seepage pressure field simulation method for fractured porous media can be referenced to the functional description of the corresponding functional module of the seepage pressure field simulation system for fractured porous media in the embodiments of the present invention, and will not be repeated here.

[0284] The module division in this embodiment of the invention is illustrative and represents only one logical functional division. In actual implementation, other division methods may be used. Furthermore, the functional modules in the various embodiments of the invention can be integrated into a single processor, exist as separate physical entities, or be integrated into a single module. The integrated modules described above can be implemented in hardware or as software functional modules.

[0285] In another embodiment of the present invention, a computer device is provided, comprising a processor and a memory. The memory stores a computer program, which includes program instructions. The processor executes the program instructions stored in the computer storage medium. The processor may be a Central Processing Unit (CPU), or other general-purpose processors, digital signal processors (DSPs), application-specific integrated circuits (ASICs), field-programmable gate arrays (FPGAs), or other programmable logic devices, discrete gate or transistor logic devices, discrete hardware components, etc. It is the computing and control core of the terminal, suitable for implementing one or more instructions, specifically suitable for loading and executing one or more instructions in the computer storage medium to achieve a corresponding method flow or corresponding function. The processor described in this embodiment of the present invention can be used to simulate the operation of a method using a seepage pressure field of a cracked porous medium.

[0286] In another embodiment of the present invention, a storage medium is provided, specifically a computer-readable storage medium (Memory), which is a memory device in a computer device used to store programs and data. It is understood that the computer-readable storage medium here can include both the built-in storage medium in the computer device and extended storage media supported by the computer device. The computer-readable storage medium provides storage space that stores the terminal's operating system. Furthermore, the storage space also stores one or more instructions suitable for loading and execution by a processor. These instructions can be one or more computer programs (including program code). It should be noted that the computer-readable storage medium here can be a high-speed RAM memory or a non-volatile memory, such as at least one disk storage device. The processor can load and execute one or more instructions stored in the computer-readable storage medium to implement the corresponding steps of the seepage pressure field simulation method for fractured porous media in the above embodiments.

[0287] Those skilled in the art will understand that embodiments of the present invention can be provided as methods, systems, or computer program products. Therefore, the present invention can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, the present invention can take the form of a computer program product embodied on one or more computer-usable storage media (including, but not limited to, disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0288] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, apparatus (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing apparatus to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing apparatus, generate instructions for implementing the flowchart illustrations and / or block diagrams. Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0289] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing device to function in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0290] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0291] The above embodiments are only used to illustrate the technical solutions of the present invention and not to limit it. Although the present invention has been described in detail with reference to the above embodiments, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the specific implementation of the present invention. Any modifications or equivalent substitutions that do not depart from the spirit and scope of the present invention should be covered within the protection scope of the claims of the present invention.

Claims

1. A method for simulating the seepage pressure field of fractured porous media, characterized in that, include: The pressure and flux on both sides of the fracture in the fractured porous medium are modeled separately, and the seepage control equation of the fractured porous medium is constructed. Integrating the seepage control equations within the control volume yields the integral form of the control equations. Embedded discrete crack meshes are defined within fractured porous media, and the interaction regions between the embedded discrete crack meshes are identified. Based on the characteristics that the introduced variable values ​​and pressure gradients are discontinuous, a pressure distribution function is established within the interaction region. Based on the pressure distribution function, the discrete expressions of each integral term in the integral form of the control equation in the interaction region are solved to obtain the discrete expressions of each integral term. Substituting the discrete expressions of each integral term into the integral form of the governing equations, we obtain a set of discrete equations and solve them to obtain the simulated values ​​of the seepage pressure field of the fractured porous medium. The flow control equations for constructing fractured porous media include: Establish the mass conservation equation for bedrock in fractured porous media: （1） Where the superscript m indicates bedrock; K m This is the bedrock permeability tensor; p m Bedrock pressure; Bedrock region; This refers to the mass exchange term from the fracture to the bedrock. Establish the mass conservation equation in the fractures of fractured porous media: （2） Among them, subscript τ The axial direction of the crack; The gradient is along the crack axis; the superscript f indicates the crack. The permeability is axial along the fracture axis. The diameter of the crack; The average pressure across the crack section, i.e. ; For the crack wall, and Represents the two side walls of the crack; denoted as Darcy velocity within the bedrock; 'n' is the unit vector of the fracture normal, its direction being relative. The surface points from the crack towards the bedrock; Construct the relationship between the sum of fluxes through the crack walls and the wall pressure: （3） in, The permeability is denoted by the fracture normal direction. Construct the relationship between the difference in flux through the fracture wall and the wall pressure and the pressure inside the fracture: （4） Set the boundary conditions as follows: （5） in, and These represent the boundaries of the bedrock and the fractured area, respectively. For boundary pressure; For the boundary Darcy velocity; Equations (1) to (5) are used as the seepage control equations for fractured porous media.

2. The method for simulating the seepage pressure field of fractured porous media according to claim 1, characterized in that, Integrating the seepage control equations over the control volume yields the integral form of the control equations, including: Integrating equation (1) over the control volume, the left-hand side of equation (1) can be obtained using the divergence theorem: in, Indicates the controlled volume area; To control the normal vector of the volume interface; d A For differential control of volumetric interface; , and These are components of the bedrock permeability tensor; and For bedrock pressure in x direction and y Partial derivatives in direction; subscript , , and This indicates the control volume interface; The right-hand side of equation (1) equals the flux flowing into the bedrock region at the fracture interface minus the flux flowing out of the bedrock region: （7） in, Indicates the crack segment within the control volume; and d represents the line segments on both sides of the crack wall; l This is a differential crack segment; n x and n y These are the components of vector n; Integrating equation (2) along the crack segment within the control volume, the divergence theorem can be applied to the left-hand side of equation (2) to obtain: （8） The integral of the right-hand side of equation (2) is expressed as: （9） Integrating equation (3) along the crack segment within the control volume, the left-hand side of equation (3) is expressed as: （10） The integral of the right-hand side of equation (3) is expressed as: （11） Integrating equation (4) along the crack segment within the control volume, it can be expressed as: （12） The left-hand side of equation (12) can be represented by equation (7); By rearranging equations (6) to (12), we can obtain the integral form of the governing equations as follows: （13） （14） （15） （16）。 3. The method for simulating the seepage pressure field of fractured porous media according to claim 2, characterized in that, The process of dividing an embedded discrete crack mesh within a cracked porous medium and identifying the interaction regions between the embedded discrete crack meshes includes: The bedrock region of the fractured porous medium is discretized using a structured rectangular grid. The location of the control volume interface is determined, with the element node being the grid center. Interaction regions are delineated, offset from the control volume by a quarter of the grid. Each interaction region is rectangular, and its vertices are composed of element nodes. Each control volume is covered by four interaction regions, and the interface of the control volume is divided into eight half-interfaces, with each interaction region covering two half-interfaces. The fractures are divided using a one-dimensional grid. A fracture segment covered by a single interaction region is a fracture grid, and the interface of the fracture grid is located at the intersection of the fracture segment and the boundary of the interaction region, with the grid center being an element node. An additional degree of freedom is added to each fracture grid on both sides.

4. The method for simulating the seepage pressure field of fractured porous media according to claim 3, characterized in that, The pressure distribution function established within the interactive region includes: Integral terms in the inductive integral form of the governing equations: （17） in, i ∈{ x , y }、 k ∈{ e , w, n, s }and j ∈{1, 2}, where, e , w , n, s This represents the four and a half interfaces in the interactive area; This is the integral of the partial derivative of the bedrock pressure along the control volume interface; This is the integral of the pressure partial derivative along one side of the crack; This is the integral of the pressure value along one side of the crack; The pressure distribution function within the interaction region is defined as follows: （18） in, ( h =1, 2, 3, 4) represent the nodal pressures of four bedrock units; ( j =1,2) are two additional degrees of freedom of the crack element, and their values ​​are always 0 if no crack passes through the interaction region; The interpolation shape function for the rectangular element, expressed in the local coordinate system with the lower left corner element node of the interaction region as the origin, is as follows: （19） in, δx and δy These are the width and height of the interactive area, respectively. In formula (18) , and These are three step functions used to characterize the discontinuity of variable values ​​on both sides of the crack, and their expressions are as follows: （20） Here, sgn(x) is the signed distance function for the crack within the interaction region, with a value of 1 on one side of the crack and a value of -1 on the other side; x h These are the coordinates of the bedrock unit nodes; In the crack The coordinates of any point on the side.

5. The method for simulating the seepage pressure field of fractured porous media according to claim 4, characterized in that, The process of solving for the discrete expressions of each integral term in the integral form of the control equation within the interaction region based on the pressure distribution function includes: Solving the nodal pressure versus interface integral of the bedrock element The contributions are as follows: （21） （22） in, It is an interactive area The bedrock unit node pressure vector; Solve the integral of the nodal pressures of the bedrock element along the fracture segment. and Contributions: （23） （24） （25） in,( fx , fy () represents the coordinates of the crack segment in the local coordinate system; Solve for the interface integral of the additional degrees of freedom of the crack. Contributions: 、 （26） （27） （28） in, Interactive area The additional degree of freedom vector of the crack segment in the structure; Solve for the additional degrees of freedom of the crack by integrating along the crack segment. and Contributions: （29） （30） （31） The discrete expressions for each integral term are obtained using the following formula: , and .

6. The method for simulating the seepage pressure field of fractured porous media according to claim 5, characterized in that, The process of substituting the discrete expressions of each integral term into the integral form of the governing equations to obtain a set of discrete equations and solving them yields the simulated values ​​of the seepage pressure field in the fractured porous medium, including: For equation (13), substituting the discrete expressions of each integral term into the integral form of the control equation, we get: （32） In this matrix, rows 1 to 4 correspond to the control volumes of unit nodes x1 to x4, respectively. and The superscript number indicates the value of the integral along the crack segment in the relevant control volume, satisfying... and ; The coefficient matrix in equation (32) can be simplified to the following form: （33） Rearranging equation (33) further, we get: （34） By superimposing all interactive regions, the matrix of equation (34) is integrated to obtain: （35） in, The number of interactive areas; For all bedrock unit node pressure vectors; Add a degree-of-freedom vector to all cracks; For equation (14), substituting the discrete expressions of each integral term into the integral form of the control equation, we get: （36） in, Interactive area Pressure vectors at nodes of the intermediate fracture element and its adjacent fracture elements; The discrete matrix for axial flow in the crack is discretized using central difference. Equation (36) can be simplified as follows: （37） Rearranging equation (37) further, we get: （38） By superimposing all interactive regions, the matrix of equation (37) is integrated to obtain: （39） in, For the pressure vector of all crack element nodes; For equation (15), substituting the discrete expressions of each integral term into the integral form of the control equation, we get: （40） Equation (40) can be simplified as follows: （41） Rearranging equation (41) further, we get: （42） By superimposing all interactive regions, the matrix of equation (41) is integrated to obtain: （43） For equation (16), substituting the discrete expressions of each integral term into the integral form of the control equation, we get: （44） Equation (44) can be simplified as follows: （45） Rearranging equation (45) further, we get: （46） By superimposing all interactive regions, the matrix of equation (46) is integrated to obtain: （47） By assembling equations (35), (39), (43), and (47), we obtain the overall discrete equation system: （48） Solving equation (48) yields the simulated values ​​of the seepage pressure field in the fractured porous medium.

7. A system for simulating the seepage pressure field of fractured porous media based on the seepage pressure field simulation method of fractured porous media according to claim 1, characterized in that, include: The module is used to model the pressure and flux on both sides of the fracture in the fractured porous medium, and to construct the seepage control equation for the fractured porous medium. The integral module is used to integrate the seepage control equations within the control volume to obtain the integral form of the control equations. The mesh embedding module is used to divide embedded discrete crack meshes in cracked porous media, identify the interaction regions between embedded discrete crack meshes, and establish the pressure distribution function inside the interaction region based on the characteristics that the introduced variable values ​​and pressure gradients are discontinuous. The discrete module is used to solve the discrete expressions of each integral term in the integral form of the control equation within the interaction region based on the pressure distribution function, and to obtain the discrete expressions of each integral term. The solver module is used to substitute the discrete expressions of each integral term into the integral form of the governing equations to obtain a discrete set of equations and solve them to obtain the simulated values ​​of the seepage pressure field of the fractured porous medium.

8. A computer device comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that, When the processor executes the computer program, it implements the steps of the seepage pressure field simulation method for fractured porous media as described in any one of claims 1 to 6.

9. A computer-readable storage medium storing a computer program, characterized in that, When the computer program is executed by the processor, it implements the steps of the seepage pressure field simulation method for fractured porous media as described in any one of claims 1 to 6.

Images