CN115114834A - Fracturing well testing simulation method under complex conditions - Google Patents
Fracturing well testing simulation method under complex conditions Download PDFInfo
- Publication number
- CN115114834A CN115114834A CN202210902933.4A CN202210902933A CN115114834A CN 115114834 A CN115114834 A CN 115114834A CN 202210902933 A CN202210902933 A CN 202210902933A CN 115114834 A CN115114834 A CN 115114834A
- Authority
- CN
- China
- Prior art keywords
- fracture
- pressure
- seepage
- dimensionless
- equation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 50
- 238000012360 testing method Methods 0.000 title claims abstract description 33
- 238000004088 simulation Methods 0.000 title claims abstract description 29
- 239000007789 gas Substances 0.000 claims abstract description 144
- VNWKTOKETHGBQD-UHFFFAOYSA-N methane Chemical compound C VNWKTOKETHGBQD-UHFFFAOYSA-N 0.000 claims abstract description 96
- 239000003345 natural gas Substances 0.000 claims abstract description 49
- 230000009977 dual effect Effects 0.000 claims abstract description 30
- 238000009826 distribution Methods 0.000 claims abstract description 19
- 239000011148 porous material Substances 0.000 claims description 103
- 239000011159 matrix material Substances 0.000 claims description 69
- 230000006870 function Effects 0.000 claims description 22
- 230000009466 transformation Effects 0.000 claims description 19
- 230000015572 biosynthetic process Effects 0.000 claims description 16
- 230000035699 permeability Effects 0.000 claims description 16
- 238000012545 processing Methods 0.000 claims description 8
- 239000000243 solution Substances 0.000 claims description 8
- 230000005465 channeling Effects 0.000 claims description 7
- 239000013598 vector Substances 0.000 claims description 7
- 230000008878 coupling Effects 0.000 claims description 6
- 238000010168 coupling process Methods 0.000 claims description 6
- 238000005859 coupling reaction Methods 0.000 claims description 6
- 239000003637 basic solution Substances 0.000 claims description 4
- 230000010354 integration Effects 0.000 claims description 4
- 230000004907 flux Effects 0.000 claims description 2
- 230000000704 physical effect Effects 0.000 abstract description 6
- 238000004458 analytical method Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 4
- 238000005325 percolation Methods 0.000 description 4
- 230000014509 gene expression Effects 0.000 description 3
- 230000008569 process Effects 0.000 description 3
- 238000003860 storage Methods 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 2
- 238000010276 construction Methods 0.000 description 2
- 230000002500 effect on skin Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 238000013178 mathematical model Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 230000000644 propagated effect Effects 0.000 description 2
- 239000000758 substrate Substances 0.000 description 2
- 241000282465 Canis Species 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 239000006185 dispersion Substances 0.000 description 1
- 230000008030 elimination Effects 0.000 description 1
- 238000003379 elimination reaction Methods 0.000 description 1
- 230000005251 gamma ray Effects 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 239000003208 petroleum Substances 0.000 description 1
- 238000003825 pressing Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 230000002195 synergetic effect Effects 0.000 description 1
- 230000001052 transient effect Effects 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/13—Architectural design, e.g. computer-aided architectural design [CAAD] related to design of buildings, bridges, landscapes, production plants or roads
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- Mathematical Physics (AREA)
- Mathematical Analysis (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Computational Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Evolutionary Computation (AREA)
- Structural Engineering (AREA)
- Civil Engineering (AREA)
- Architecture (AREA)
- Operations Research (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
Abstract
The invention discloses a fracturing well testing simulation method under a complex condition, which can combine a Boundary Element Method (BEM) and a Fredholm integral equation in a Laplace space, establish an unstable well testing model for a fracturing well in a gas reservoir in any shape, which can comprehensively consider the influences of any shape of the gas reservoir, natural gas high-pressure physical properties, reservoir dual medium characteristics, fracturing flow conductivity, the asymmetry of two wings of a fracturing and the asymmetry of fracturing flow distribution, and obtain a simulation result according to the solution of the model, wherein the simulation result can comprise a double logarithmic curve and a fracturing flow distribution curve of non-dimensional bottom hole pressure and non-dimensional bottom hole pressure derivative. The invention can provide important support for well testing interpretation of the fractured well under complex conditions.
Description
Technical Field
The invention relates to the technical field of a fracturing well bottom pressure simulation method, in particular to the technical field of a fracturing well testing model construction method based on a Boundary Element Method (BEM).
Background
The well testing technology of oil and gas wells is a very important technology for recognizing oil and gas reservoirs, and is known as the eye of oil and gas reservoir development by the petroleum industry. Through well testing, some important dynamic parameters of the oil and gas reservoir deeply buried underground such as stratum permeability, shaft pollution coefficient and the like can be obtained, and a reasonable unstable well testing model (the unstable well testing model is one of the unstable seepage models) and the like are established.
Accurate simulation of bottom hole pressure transients in oil and gas wells is the theoretical basis and necessary premise for well testing analysis techniques. Conventional analytical or semi-analytical methods can only be used to solve and characterize the bottom hole pressure transients of oil and gas wells in regularly shaped reservoirs, while numerical methods can solve more complex problems.
Numerical methods can be divided into two broad categories: a region method and a boundary method. Both finite difference and finite element methods belong to the first category, while boundary element methods belong to the second category. Among these, the boundary element method is superior to the area method in some aspects, such as: (1) the boundary element method retains the understood analytic characteristics to a greater extent, so that the boundary element method has higher precision than regional methods such as finite difference, finite element and the like, and can meet the requirement on precision when solving a well test model; (2) the finite element method and the finite difference method need to divide units in the whole area, and the boundary element method only needs to divide the units on the boundary of the area, so the boundary element method has the characteristic of dimension reduction, and the formed matrix order is smaller.
In the prior art, the research for solving the seepage problem of the oil and gas reservoir by using the boundary element method mainly focuses on the aspect of non-fractured wells. In recent years, some researchers have further studied the seepage problem of fractured wells by using the boundary element method and made some progress. However, at present, the prior art still lacks accurate construction and simulation of a well testing model with any shape of a gas reservoir, asymmetric fracturing structure or flow and multiple medium characteristics of a reservoir, and cannot obtain fracturing well testing explanations which often exist in practice and under the complex conditions.
Disclosure of Invention
The invention aims to construct a well testing model considering the influences of any shape of a gas reservoir, natural gas high-pressure physical properties, reservoir dual medium characteristics, fracture fracturing flow conductivity, fracture two-wing asymmetry and fracture flow asymmetry, and obtain a simulation method capable of accurately calculating the bottom hole pressure transient state of a gas well under the complex condition based on the constructed unstable well testing model.
The technical scheme of the invention is as follows:
a well testing simulation method for a fractured well under complex conditions comprises the following steps:
s1, constructing a physical model of a fracturing well, wherein the physical model considers that the outer boundary of a gas reservoir is in any shape, the reservoir is a double-pore system containing a double-pore structure, the left wing and the right wing of a fracture can be symmetrical or asymmetrical, and the fracture has limited flow conductivity; wherein, the dual pore structure means that the pore structure of the reservoir comprises two pore systems, namely a natural fracture system formed by the natural fracture and a matrix pore system formed by the matrix pore;
s2, constructing a gas reservoir seepage master control model of the dual pore system, wherein the model comprises the following steps:
s21, coupling the inner boundary pressure condition of the zonally converged fracturing with the mass conservation equation of the natural fracture system by utilizing a Dirac generalized function and an integral equation, and simultaneously establishing a motion equation and a state equation of natural gas in the natural fracture system and a channeling equation of the natural gas in the matrix pore system and the natural fracture system to derive a seepage master control differential equation of the natural fracture system;
s22, a gas seepage master control differential equation of the matrix pore system is expressed by using a gas seepage differential equation under a pore structure, and a gas reservoir seepage master control differential equation of the dual pore system, namely the gas reservoir seepage master control model, is formed by the seepage master control differential equation of the natural fracture system and the seepage master control differential equation of the matrix pore system;
s3, constructing a factorial formation seepage model of the dual pore system, comprising:
s31 setting an initial pressure condition equation of the gas reservoir and an outer boundary pressure condition equation under different outer boundary conditions;
s32, combining the initial pressure condition equation and the outer boundary pressure condition equation with the gas reservoir seepage master control differential equation to obtain the factorial formation seepage model;
s4, introducing a dimensionless quantity, and performing dimensionless transformation on the dimensionless stratum seepage model to obtain a dimensionless stratum seepage model;
s5 obtaining a gas reservoir outer boundary seepage model of the dual pore system, comprising:
s51, Laplace transformation is carried out on the dimensionless stratum seepage model, the seepage master control differential equation of the dimensionless matrix pore system is substituted into the seepage master control differential equation of the dimensionless natural fracture system, the pressure parameter of the matrix pore system is eliminated, and the seepage master control differential equation of the transformed natural fracture system is obtained;
s52, converting the transformed seepage master control differential equation of the natural fracture system into an outer boundary seepage integral equation of the gas reservoir based on a boundary element solution to obtain an outer boundary seepage model of the gas reservoir;
s6, carrying out unit discrete processing on the gas reservoir outer boundary seepage model to obtain a first linear equation set;
s7, respectively constructing gas reservoir seepage models for two wings of the fracture, wherein the gas reservoir seepage models take the influences of limited flow conductivity of the fracture, unequal lengths of the two wings of the fracture and asymmetric flow distribution of the two wings into consideration, and obtaining fracture seepage models;
s8, converting the fracturing fracture seepage model into a Fredholm integral equation by using Laplace transformation and double integration, and performing unit discrete processing to obtain a second linear equation set;
s9, the first linear equation set and the second linear equation set are combined to obtain a closed matrix;
s10, solving the closed matrix, and obtaining a simulation result by numerical inversion.
According to some embodiments of the invention, in S2, the coupling of the inner boundary pressure condition of the zonal confluent fractures and the mass conservation equation of the natural fracture system results in the following coupled mass conservation equation:
where ρ represents the natural gas density in the natural fracture system; v. of x Representing the seepage velocity of natural gas in the x direction in a natural fracture system; v. of y The seepage velocity of natural gas in the y direction in the natural fracture system is represented; x and y respectively represent an x coordinate and a y coordinate; h represents reservoir thickness; f represents a line integral region corresponding to the trajectory of the fracturing fracture; rho sc Represents the density of natural gas at surface conditions; q. q.s fsc The linear density flow rate function of the fracture, i.e. the flow rate of the fracture per unit length vertically flowing from the natural fracture system of the formation, converted to surface conditions, is a function of the location, i.e. q fsc ＝q fsc (W); dl is the infinitesimal length; w is any point on the fracture, W ═ W (x) w ,y w )；δ(x-x w ,y-y w ) Is a two-dimensional dirac function; x is the number of w 、y w X and y coordinates of any point on the hydraulic fracture are respectively expressed, t represents time, phi represents the porosity of a natural fracture system, and q represents the porosity of the natural fracture system * Indicating the amount of cross flow.
According to some embodiments of the present invention, in S2, the gas reservoir seepage control model of the dual pore system includes:
the seepage master differential equation of the natural fracture system:
seepage master control differential equation for the matrix pore system:
wherein:
representing the pseudo-pressure of the natural fracture system; k represents the gas reservoir permeability of the natural fracture system; h represents reservoir thickness; t is sc Represents the temperature at ground conditions; t represents the gas reservoir temperature; w represents any point on the fracture, W ═ W (x) w ,y w )；p sc Representing the ground standard condition pressure; q. q.s fsc A linear density flow function representing the fracture; delta (x-x) w ,y-y w ) Is a two-dimensional dirac function; t represents time; mu.s i Represents the viscosity of natural gas under the original condition; c gi Representing the compressibility of natural gas under original conditions; α represents a shape factor; k is a radical of m Represents the permeability of the matrix pore system; psi m Representing the pseudo-pressure of the matrix pore system; phi is a m Is the porosity of the matrix pore system; p represents pressure; p is a radical of 0 Represents a reference pressure; z is a natural gas deviation factor; μ is the natural gas viscosity; x and y are x coordinate and y coordinate respectively; x is the number of w 、y w Respectively representing x and y coordinates of any point on the hydraulic fracture; f is a line integral area representing a hydraulic fracture trajectory; dl is the infinitesimal length.
According to some embodiments of the invention, in S3, the external boundary pressure condition equations for the different external boundary situations are as follows:
if the boundary is a closed boundary, the external boundary pressure condition equation is as follows:
where Γ represents the gas reservoir outer boundary, p represents pressure,
is an outward normal vector on the outer boundary;
if the pressure is the constant pressure boundary, the conditional equation of the pressure of the outer boundary is as follows:
p| Γ ＝p i (13)
wherein p is i Representing an evenly distributed original pressure in the gas reservoir;
if the outer boundary is a mixed boundary, the outer boundary pressure condition equation is as follows:
wherein, γ 1 、γ 2 And gamma 2 Is a combination constant.
According to some embodiments of the present invention, in S3, the outer boundary pressure condition equation is a pseudo-pressure form of the outer boundary pressure condition equation under the closed boundary as follows:
the initial pressure condition equation is in the form of a pseudo pressure as follows:
ψ| t＝0 ＝ψ m | t＝0 ＝ψ i (11)
where Γ represents the gas reservoir outer boundary, p represents pressure,
representing outward normal vectors on an outer boundary
The pseudo-pressure of the natural fracture system is shown,
representing the original pseudo-pressure of the gas reservoir, t representing time,. phi m Denotes the pseudo-pressure, p, of the pore system of the substrate 0 Denotes a reference pressure, p i Representing the original pressure of the gas reservoir distributed uniformly, Z represents the natural gas deviation factor, and μ represents the natural gas viscosity.
According to some embodiments of the invention, in S4, the dimensionless formation seepage model includes:
wherein:
wherein psi D Is dimensionless pseudo-pressure, # in a gas reservoir natural fracture system Dm Is dimensionless pseudo-pressure,. psi.in the pore system of the gas reservoir matrix i Is the original pseudo pressure of the gas reservoir, omega is the elastic storage-volume ratio, lambda is the cross-flow coefficient, x D Is a dimensionless x coordinate, y D Is a dimensionless y coordinate, q fD Dimensionless linear density flux, dl, for fracturing D Is a dimensionless infinitesimal length,
is a dimensionless outward normal on the outer boundaryVector, t D Denotes dimensionless time, T sc Is the temperature under ground conditions, q fsc Linear density flow function representing fracture, T is gas reservoir temperature, C gi Is the compressibility factor, phi, of natural gas under the original conditions m Is the porosity of the matrix pore system, phi denotes the porosity of the natural fracture system, alpha is the shape factor, k m Is the permeability of the matrix pore system, K is the permeability of the natural fracture system, L ref Is the reference length of the fracture, x fL Is the left wing length of the crack, x fR Is the right wing length of the crack, μ i Is the viscosity of natural gas at the original condition and dl is the infinitesimal length.
According to some embodiments of the invention, in S5, the seepage-dominated differential equation of the transformed natural fracture system is as follows:
wherein,
f(s) represents a transformation function related to the parameters ω, λ and the Laplace variable s; s represents a Laplace variable;
expressing dimensionless simulated pressure in the gas reservoir natural fracture system after Laplace transformation;
showing the dimensionless linear density flow of the fracture at the W point after Laplace transformation; f D Representing the dimensionless line integral region, x, corresponding to the trajectory of the fracture wD A dimensionless abscissa representing the W point; y is wD Expressing the dimensionless ordinate of the W point, i.e. W (x) w ,y w )。
The gas reservoir outer boundary seepage model is as follows:
wherein: p 'is any point on the outer boundary of the stratum, Q is any point in the stratum (including in the region and on the outer boundary), and G (P', Q, s) represents the basic solution of the boundary element
The middle P point selects the solution at any point P' on the outer boundary gamma, where K is 0 Representing a zero order deformed Bessel function, r D Represents a dimensionless radial distance;
expressing dimensionless simulated pressure of a Q point in the gas reservoir natural fracture system after Laplace transformation;
representing dimensionless simulated pressure of a P' point in the gas reservoir natural fracture system after Laplace transformation; w is any position point on the fracturing crack; θ is a constant related to the geometry at point Q; beta is the inner angle of the left and right tangents to the outer boundary at point Q.
In the above-described embodiments, the first and second embodiments,
are both functions of spatial location points and Laplace variables s, where,
the mid-space location point is not determined,
the middle space position point is a point Q, and the point Q can be any point in the region omega or any point on the outer boundary gamma;
the middle space position point is P 'point, and P' is any point on the outer boundary gamma.
According to some embodiments of the invention, in S7, the fracture-seepage model is as follows:
wherein,
wherein psi fD Dimensionless pseudo-pressure for fracturing f Pseudo-pressure, x, for fracturing fLD Is the dimensionless length, x, of the left wing of the crack fRD Dimensionless length of right wing of crack press C fD Is a fracture dimensionless diversion coefficient, K f Is the fracture permeability, q LD Is the dimensionless yield of the left wing of the fracture, q L Is the yield of the left wing of the fracturing fracture, q RD Is the dimensionless yield of the right flank of the fracture, q R Is the yield of the right wing of the fracturing, W fD Is the dimensionless crack width, W f Is the crack width.
According to some embodiments of the invention, in S9, the closed matrix includes:
a first set of linear equations:
a second linear system of equations:
wherein N is b Representing the number of discrete cells of the outer boundary Γ of the gas reservoir,
and
is a coefficient of a first linear system of equations, N FL 、N FR Respectively representing the number of discrete elements, Deltax, of the left and right wings of the fracturing fracture DL ＝x fLD /N FL Denotes the left wing discrete element length, Δ x DR ＝x fRD /N FR Indicating the length of the right wing discrete element, i ═ 1,2,3 … … N b Denotes i taken over N b Discrete units of outer boundary, i ═ 1,2,3 … … N FL +N FR N representing i' Take through crack FL +N FR A discrete unit, k is 1,2,3 … … N FL +N FR N representing k-cut through fracture FL +N FR A discrete unit, Q k The midpoint of the kth discrete cell on the fracture is shown.
According to the simulation method, the fracturing well testing simulation device under the complex condition can be obtained, and comprises a storage medium which stores programs, algorithms and/or models capable of realizing any simulation method.
The simulation method combines a Boundary Element Method (BEM) and a Fredholm integral equation in a Laplace space, establishes an unstable well testing model for fracturing gas wells in gas reservoirs with any shapes, which can comprehensively consider the influences of any shapes of the gas reservoirs, high-pressure physical properties of natural gas, double medium characteristics of reservoirs, fracturing flow conductivity, asymmetry of two wings of fracturing and flow distribution asymmetry of fracturing, successfully solves the model, can draw a high-quality log-log curve according to simulation results, accurately describes the bottom hole pressure dynamics, and provides important application support for the interpretation of the fracturing well testing under complex conditions.
Drawings
FIG. 1 is a physical model of a gas reservoir of an arbitrarily shaped dual pore system and a gas well containing a confined flow fracture involving an embodiment.
FIG. 2 is a schematic illustration of discrete processing of the outer boundary elements of a gas reservoir of an arbitrarily shaped dual pore system as contemplated in the embodiments.
Fig. 3 is a physical model of a limited fracture conductivity fracturing system according to an embodiment.
Fig. 4 is a log-log plot of a complex-shaped dual pore system gas reservoir obtained in an embodiment.
FIG. 5 is a log-log plot of vertical fracture wells at different locations in an embodiment.
FIG. 6 shows the fracture conductivity C obtained in the present embodiment fD And (3) an influence graph of the log curves of the limited diversion vertical fracture well of the gas reservoir of the complex-shaped dual-pore system.
Fig. 7 is a graph showing the effect of asymmetry of left and right wings of a fracture on a log-log curve of a finite flow guide vertical fracture well of a complex-shaped dual-pore system gas reservoir in a specific embodiment.
Fig. 8 is a physical model of a fracture with equal wing length according to an embodiment.
FIG. 9 is a graph of the fracture flow profile with equal wing length according to an embodiment.
Detailed Description
The present invention is described in detail below with reference to the following embodiments and the attached drawings, but it should be understood that the embodiments and the attached drawings are only used for the illustrative description of the present invention and do not limit the protection scope of the present invention in any way. All reasonable variations and combinations that fall within the spirit of the invention are intended to be within the scope of the invention.
According to the technical scheme of the invention, some specific embodiments comprise the following processes:
s1, constructing a physical model of a fracturing well, wherein the physical model considers that the outer boundary of a gas reservoir is in any shape, the reservoir is a double-pore system containing a double-pore structure, the left wing and the right wing of a fracture can be symmetrical or asymmetrical, and the fracture has limited flow conductivity;
further, the physical model of the fractured well can also set: isothermal seepage is carried out on the gas reservoir, the original pressure in the gas reservoir is uniformly distributed, and the flow of the gas reservoir meets the Darcy law.
Wherein, the dual pore structure means that the pore structure of the reservoir comprises both natural fractures and matrix pores, and the whole reservoir comprises both pore systems of natural fractures formed by the natural fractures and matrix pores formed by the matrix pores (the dual medium system shown in fig. 1).
More particularly, some embodiments are directed to constructing a physical model of a fractured well as shown in fig. 1, wherein the fractured well is formed by hydraulic fracturing of a vertical gas well, a double-wing artificial fracture is formed in a reservoir, the gas reservoir has an arbitrarily shaped outer boundary, and the gas well has a constant surface production q sc The production is carried out, and:
the reservoir is of a dual-pore structure and comprises a natural fracture system and a matrix pore system, wherein the natural fracture is a main flow channel, the matrix pore is a main reservoir space, and the matrix pore forms a supply source of the natural fracture;
the left wing and the right wing of the fracturing crack can be symmetrical or asymmetrical, and the length of the left wing is assumed to be x fL The length of the right wing is x fR ；
The fracture is of width W f Vertical fractures with limited conductivity;
the gas reservoir temperature is T, and isothermal seepage is adopted;
the original pressure in the reservoir is uniformly distributed and is p i ；
The reservoir flow satisfies darcy's law.
Based on the physical model, the simulation method further utilizes the unit dispersion of Boundary Element Method (BEM), Laplace transformation, double integral and Fredholm integral equations to establish a more comprehensive and stricter unstable well testing model for the fractured gas well of the gas reservoir with any shape and double pores, which can comprehensively consider the influences of any shape of the gas reservoir, natural gas high-pressure physical properties, fracture flow conductivity, fracture two-wing asymmetry and fracture flow distribution asymmetry, solve the model and draw a high-quality pressure dynamic curve, and provide an important implementation mode for the fractured well testing explanation under the complex condition.
And setting a physical model of the fracturing well based on the steps of:
s2, constructing a gas reservoir seepage master control model of the dual pore system, which specifically comprises the following steps:
s21, considering compressibility of natural gas in the gas reservoir, coupling the inner boundary pressure condition of the zonally converged fracturing with a mass conservation equation of the natural fracture system by using a Dirac generalized function and an integral equation, and combining the pressure condition with a motion equation and a state equation of natural gas in the natural fracture system and a channeling equation of the natural gas in a matrix pore system and the natural fracture system to derive a seepage master control differential equation of the natural fracture system;
s22, aiming at the matrix pore system, the traditional gas reservoir seepage differential equation is used as the seepage master control differential equation of the matrix pore system, and the gas reservoir seepage master control differential equation of the dual pore system, namely the gas reservoir seepage master control model, is composed of the seepage master control differential equation of the natural fracture system and the seepage master control differential equation of the matrix pore system.
S3, constructing a factorial formation seepage model of the dual pore system, wherein the factorial formation seepage model specifically comprises the following steps:
s31 setting an initial pressure condition equation of the gas reservoir and an outer boundary pressure condition equation under different outer boundary conditions;
s32, combining the initial pressure condition equation and the outer boundary pressure condition equation with the gas reservoir seepage master control differential equation to obtain the stratum seepage model, wherein the model is a factorial stratum seepage model of the gas reservoir of any shape and double pore system coupled with the strip-shaped confluent inner boundary condition of the fracture.
S4, introducing a dimensionless quantity, and performing dimensionless transformation on the dimensionless stratum seepage model to obtain a dimensionless stratum seepage model, wherein the dimensionless stratum seepage model comprises the following steps: a dimensionless natural fracture system seepage master differential equation, a dimensionless matrix pore system seepage master differential equation, a dimensionless initial pressure condition equation, and a dimensionless outer boundary pressure condition equation.
S5, obtaining a gas reservoir outer boundary seepage model of the dual pore system, which specifically includes:
s51, Laplace transformation is carried out on the obtained dimensionless stratum seepage model, the seepage master control differential equation of the dimensionless matrix pore system is substituted into the seepage master control differential equation of the dimensionless natural fracture system, the pressure parameter of the matrix pore system is eliminated, and the seepage master control differential equation of the natural fracture system in the region, which is obtained by eliminating the matrix pressure and is coupled with the strip convergence boundary condition of the fracturing, is obtained, namely the seepage master control differential equation of the transformed natural fracture system;
s52, converting the seepage master control differential equation of the natural fracture system into an outer boundary seepage integral equation of the gas reservoir based on a boundary element solution to obtain the gas reservoir outer boundary seepage model.
S6, unit discrete processing is carried out on the obtained gas reservoir outer boundary seepage model to obtain a first linear equation set.
S7, respectively constructing gas reservoir seepage models considering limited flow conductivity of the fracturing, unequal lengths of two wings of the fracturing and asymmetric influence of flow distribution of the two wings of the fracturing to obtain the fracturing seepage models.
S8, converting the obtained fracture seepage model into a Fredholm integral equation by Laplace transformation and double integration, and performing unit discrete processing to obtain a second linear equation set.
S9, the first linear equation set and the second linear equation set are combined to obtain a closed matrix with the equation number equal to the unknown number.
S10, solving the obtained closed matrix, utilizing numerical inversion to draw a double logarithmic curve of the dimensionless bottom hole pressure derivative and the dimensionless bottom hole pressure derivative of the obtained fractured well, dividing the flowing stage of the gas reservoir according to the obtained double logarithmic curve and a seepage mechanism, analyzing and obtaining the position of the fractured well in the gas reservoir with a complex shape, the dimensionless diversion coefficient of the fracturing well, the influence of the asymmetry of two wings of the fracturing well on the double logarithmic curve and the like, and obtaining a simulation result.
Further, in some embodiments, the seepage-dominated differential equation for the natural fracture system is constructed as follows:
according to the mass conservation principle, the property of a delta generalized function and the connotation of integral are combined to obtain a mass conservation equation coupling the boundary conditions in the vertical fracture well in the natural fracture system, and the mass conservation equation comprises the following steps:
in order to distinguish the natural fracture system from the matrix pore system and the artificial fracture, the pressure and other parameters of the natural fracture system are not provided with subscripts, the pressure and other parameters of the matrix pore system are provided with subscripts "m", and the pressure and other parameters of the artificial fracture are provided with subscripts "f".
Wherein the mass conservation equation of the matrix pore system still takes the traditional form:
the channeling equation of natural gas in a matrix pore system and a natural fracture system is as follows:
where ρ, ρ m The natural gas density in the natural fracture system and the matrix pore system, respectively, is kg/(m) 3 )；v x The seepage velocity of natural gas in a natural fracture system in the x direction is m/s; v. of y The seepage velocity of natural gas in a natural fracture system in the y direction is m/s; x is the x coordinate, m; y is the y coordinate, m; h is reservoir thickness, m; f represents a line integral region corresponding to the trajectory of the fracturing fracture; rho sc Is the density of natural gas under ground conditions, kg/(m) 3 )；q fsc Is the linear density flow of the fracture, i.e., the flow (converted to surface conditions) from the natural fracture system of the formation into hydraulic fractures per unit length, as a function of location, i.e., q fsc ＝q fsc (W) in m 3 V (s.m); dl is the infinitesimal length, m; w is any point on the fracture, W is W (x) w ,y w )，m 3 /(s.m)；δ(x-x w ,y-y w ) Is a two-dimensional dirac function; phi, phi m Respectively, natural fracture system porosity and matrix pore system porosity, fraction; t is time, s; q. q.s * Is the flow rate, kg/(m) 3 ·s)；ρ r Is the natural gas density under the average pressure of a matrix pore system and a natural fracture system, and is kg/(m) 3 ) (ii) a α is the shape factor, m -2 ；k m Is the permeability of the matrix system, m 2 ；μ r Is the natural gas viscosity, mPa · s, at the average pressure of the matrix pore system and the natural fracture system.
The equation of motion of a gas in a natural fracture system uses darcy's law as follows:
wherein: k is the natural fracture system permeability, m 2 (ii) a μ is the natural gas viscosity, pas.
The equation of state of the gas in a natural fracture system is:
pV＝nZRT (6)
wherein: v is the volume of gas, m 3 (ii) a n is the number of moles, kmol; z is a natural gas deviation factor, and is dimensionless; r is the general gas constant, J/(kmol. K); t is the gas reservoir temperature, K.
The above formula can be written as:
wherein: m g Is the molecular weight of natural gas, kg/kmol.
Substituting the equation of motion of gas in a natural fracture system, the equation of state of gas in the natural fracture system and the channeling equation (3) of gas between a matrix pore system and the natural fracture system into the natural fracture system mass conservation equation (1) coupled with the inner boundary condition of the fracturing fracture, and introducing pseudo-pressure into the natural fracture system mass conservation equation (1) in consideration of the influence of the high-pressure physical property of natural gas
After linearization treatment, the seepage master control differential equation of the gas reservoir natural fracture system of the dual-pore system with any shape and any shape coupled with the strip-shaped confluent boundary condition of the fracture can be obtained as follows:
wherein: k is the gas reservoir permeability, m 2 ；T sc Is the temperature under ground conditions, K; t is the temperature under gas reservoir conditions, K; mu.s i Is the viscosity of natural gas under the original conditions, Pa · s; c gi Is the compressibility factor of natural gas under original conditions, Pa -1 。
Further, in some embodiments, the percolation primary differential equation for the matrix pore system is in the following conventional form percolation differential equation:
further, in some embodiments, the initial pressure condition equation for the gas reservoir is set as follows:
setting the formation pressure distribution at the initial moment to be even and p i Then the initial pressure condition equation of the gas reservoir is:
p| t＝0 ＝p m | t＝0 ＝p i (10)
with pseudo-pressure, the above equation can be:
ψ| t＝0 ＝ψ m | t＝0 ＝ψ i (11)
further, in some embodiments, the equation for the pressure condition of the outer boundary of the gas reservoir is set as follows:
if the boundary is a closed boundary, the external boundary pressure condition equation is as follows:
wherein:
is the outward normal vector on the outer boundary.
If the pressure is the constant pressure boundary, the conditional equation of the pressure of the outer boundary is as follows:
p| Γ ＝p i (13)
if the outer boundary is a mixed boundary, the outer boundary pressure condition equation is as follows:
wherein: gamma ray 1 、γ 2 And gamma 2 Is a combination constant.
In the above embodiments, for the gas reservoir, the outer boundary mostly reflects the characteristics of the closed boundary or the infinite boundary within the test time, and the situation of the constant pressure boundary rarely occurs. In the mathematical model, when the outer boundary radius is large enough, the pressure characteristic under the infinite boundary can be characterized by the closed boundary, so that, in the specific implementation, the outer boundary pressure condition equation (12) under the condition of the closed boundary can be directly and preferably used.
Using pseudo-pressure, equation (12) can be transformed:
the above formula is an outer boundary pressure condition equation in the form of pseudo pressure.
The factorial formation mathematical model of an arbitrarily shaped dual pore system reservoir that couples the band-like in-sink boundary conditions of a fracture is composed of equations (8), (9), (11) and (15).
Further, in some embodiments, the dimensionless formation permeability model is constructed as follows:
the following quantities are defined:
wherein psi D Is dimensionless pseudo-pressure in a gas reservoir natural fracture system and dimensionless; psi Dm Is dimensionless pseudo-pressure in the pore system of the gas reservoir matrix; psi i Is the original pseudo-pressure of the gas reservoir, Pa/s; p is a radical of 0 Is a reference pressure, Pa; p is a radical of formula f Is the pressure in the pack, Pa; x is the number of D Is a dimensionless x coordinate, dimensionless; y is D Is a dimensionless y coordinate, dimensionless; dl (dl) D The length of the dimensionless infinitesimal is dimensionless;
dimensionless outward normal vectors on the outer boundary, dimensionless; l is ref Is a reference length, and may be, for example, half the total length of the crack, i.e., L ref ＝(x fL +x fR )/2，m；x fL Is the crack left wing length, m; x is the number of fR Is the crack right wing length, m; q. q.s fD The fracture is dimensionless and has no linear density flow and dimension; omega is elastic storage-capacity ratio and is dimensionless; λ is the cross-flow coefficient, dimensionless.
Equations (8), (9), (11) and (15) can be transformed into the following dimensionless model:
equations (16) - (19) are the dimensionless formation seepage model for an arbitrarily shaped dual pore system reservoir coupled with the band-like in-plane boundary condition of a fracture. Wherein, formula (16) is a seepage master differential equation of the dimensionless natural fracture system, formula (17) is a seepage master differential equation of the dimensionless matrix pore system, formula (18) is an initial pressure condition equation of the dimensionless, and formula (19) is an outer boundary pressure condition equation of the dimensionless.
Further, in some embodiments, the model of gas reservoir outer boundary seepage is constructed as follows:
(1) introduction based on t D The dimensionless formation seepage model can be transformed into:
wherein: v 2 For the Laplace operator, there are:
obtained from (21):
further, based on the idea of solving the problem with the boundary element, the differential equation in the region is changed into the integral equation on the boundary, and then the boundary can be divided into boundary units, the integral equation on the boundary is dispersed into a linear algebraic equation, and the problem of solving the differential equation is converted into a problem of solving the algebraic equation, which specifically comprises:
substituting formula (21) for formula (20) yields:
wherein:
equation (25) is then the natural fracture system seepage dominated differential equation in the region coupled with the band-like confluent boundary condition of the fracture, with the matrix pore system pressure parameters removed.
(2) Equation (25) is converted to the integral equation of the seepage over the outer boundary of an arbitrarily shaped dual pore system gas reservoir as follows:
let the basic solution G (P, Q, s) when applying the boundary element method satisfy the following equation:
▽ 2 G(P,Q,s)-f(s)G(P,Q,s)+πδ(P,Q)＝0 (26)
wherein: p, Q are any two points in the formation.
The basic solution of the above formula is as follows:
through a series of derivation, the integral equation of equation (25) can be obtained as:
in the above formula, P' is an arbitrary point on the boundary Γ, W is a well point, and Q is an arbitrary point within the region Ω.
If the Q point is not only within the region Ω but also on the boundary Γ, the above formula is replaced by the following formula to obtain the gas reservoir outer boundary seepage model:
wherein: theta is a constant related to the geometry at point Q,
wherein: beta is the inner angle of the left and right tangents to the outer boundary at point Q.
Further, in some embodiments, with reference to fig. 2, the obtaining the first linear equation set by performing a unit discrete processing on the gas reservoir outer boundary seepage model includes:
splitting the gas reservoir outer boundary Γ into N b The pressure distribution on the outer boundary is expressed by adopting linear discrete units, namely, the pressure value on any discrete unit of the outer boundary is obtained by interpolating the pressure values on the nodes at the two ends of the unit, and the outer boundary can be specified to take the anticlockwise direction as the positive direction in the calculation.
The left wing and the right wing of the crack are respectively divided into N FL 、N FR Discrete units, whereby the flow distribution over the fracture is represented by the normally discrete units, i.e. fracturingWhen the linear density flow rates at each location on the same discrete element are equal, equation (29) can be converted to a first linear equation set as follows:
wherein:
and
are coefficients of a linear system of equations.
The node where k in equation (31) passes through the outer boundary (k is 1,2, … N) b ) Obtaining N b A linear algebraic equation, and the current unknowns include
In total (N) b +N FL +N FR ) And the number of equations is less than the number of unknowns, so that the solution cannot be carried out, and other equations are needed.
Using pseudo-pressure, if Q k Points are taken within the study domain rather than on the outer boundary, the dimensionless pseudo-pressure solution for any point Q within the study domain is:
q in formula (32) k When the point is located at the midpoint of the fracture, (N) is obtained FL +N FR ) A linear algebraic equation, such that equation (31) and equation (32) together represent N b +N FL +N FR A linear algebraic equation. However, N is increased again at this time FL +N FR New unknowns, dimensionless pseudo-pressure at mid-point of fracture
It is obvious that here one can
Is abbreviated as
Wherein k is N b +1,N b +2,……，N b +N FL +N FR . The unknown numbers in the formulae (31) and (32) now share [ N b +2(N FL +2N FR )]Each is respectively
And
because the number of the equations is less than the number of the unknown numbers, the equations cannot be solved, and other equations need to be simultaneously established, so that a fracture seepage model with limited diversion is introduced afterwards.
Further, in some embodiments, referring to fig. 3, the fracture-seepage model building process includes:
the following quantities are defined:
wherein psi fD The method is dimensionless pressure simulation and dimensionless for pressing cracks; psi f Is the pseudo-pressure of the fracturing fracture, Pa/s; x is the number of fLD The left wing of the crack is dimensionless in length and dimension; x is the number of fRD The right wing of the fracturing crack has dimensionless length and dimension; x is the number of fL Is the crack left wing length, m; x is the number of fR Is the crack right wing length, m; c fD The fracture non-dimensional flow conductivity coefficient is obtained; k f Is the fracture permeability, m 2 ；q LD Is dimensionless yield and dimension of the left wing of the fracturing fracture; q. q.s L Is the yield of the left wing of the fracturing fracture, m 3 /s；q RD Is the dimensionless yield of the right flank of the fracture,dimensionless; q. q.s R Yield of right wing of fracturing, m 3 /s；W fD The width of the dimensionless crack is zero, and the dimension is zero; w f Is the crimp gap width, m.
According to the definition, the fracture seepage model considering the influences of limited flow conductivity of the fracture, unequal lengths of two wings of the fracture and asymmetric flow distribution of the two wings is obtained as follows:
further, in some embodiments, obtaining the second system of linear equations includes:
laplace transformation is carried out on the fracture seepage model formed by the formulas (33) to (39), and by means of double integration, Fredholm integral equations respectively describing the pressure distribution of the left wing and the right wing of the fracture are obtained as follows:
in the formula,
dimensionless bottom hole pseudo pressure in Laplace space is dimensionless; the other upper band transverse lines are all corresponding quantities in the Laplace space.
As described above, the left and right wings of the fracturing crack are divided into N FL 、N FR A discrete unit, the whole crack is divided into N FL +N FR A discrete unit, the left wing discrete unit length is Deltax DL ＝x fLD /N FL The length of the right wing discrete unit is Deltax DR ＝x fRD /N FR . Suppose that
Is the midpoint of the j-th slit discrete cell, x Dj For the jth endpoint (node), the linear density flow on the same discrete unit is uniform, and for the kth discrete unit (k is more than or equal to 1 and less than or equal to N) on the left wing of the crack FL ) Equation (40) becomes:
similarly, for the kth discrete cell (N) on the right wing of the crack FL +1≤k≤N FL +N FR ) Equation (41) becomes:
wherein,
dimensionless pseudo pressure, dimensionless, at the midpoint of the kth cell.
(N) is obtained by dividing k in the formulae (42) and (43) by the midpoint of the discrete fracture cell FL +N FR ) A linear algebraic equation. Herein, the
And in formula (32)
The correspondence of (a) is as follows:
therefore, in practice, the unknowns in the expressions (42) and (43) are merely 3 more than the unknowns in the expressions (31) and (32):
while the formulae (42) and (43) represent (N) FL +N FR ) A linear algebraic equation.
Further, there are flow conditions:
equations (44) to (47) represent 3 linear algebraic equations, i.e. the second linear system of equations.
Accordingly, the expressions (42), (43) and (45) to (47) collectively represent (N) FL +N FR +3) linear algebraic equations.
Further, in some embodiments, obtaining the closed matrix comprises:
as is clear from the above, the formulae (31), (32), (42), (43), (45) to (47) collectively represent [ N b +2(N FL +N FR )+3]A linear algebraic equation, and the number of unknowns is also { N } b +2(N FL +N FR ) +3 are respectively
The number of the equations is equal to the number of the unknown numbers, and the equations are linear algebraic equations, so that the solution can be realized. The matrix obtained by simultaneous calculation is a small dense matrix which can be solved by adopting a Gaussian elimination method.
The embodiment of the invention avoids the dependence on a finite element method or a finite difference method, and the obtained model comprehensively considers the influences of any shape of a gas reservoir, high-pressure physical properties of natural gas, dual medium characteristics of a reservoir, fracture conductivity, unequal fracture wing lengths and asymmetric fracture wing flow distribution.
Further, obtained by the above method
The effects of the reservoir effect and the skin effect are not considered. The dimensionless well storage coefficient C.sub.F.can be studied according to Van Evandingen et al (Van-Evandingen, A.F.,1953.The skin effect and its induced property on The productive capacity of a well.J.Pet.Technol.5(6), 171-176.; Kucuk, F., Ayestaran, L.,1985.Analysis of a synergistic measured pressure and yield flow in transformed well testing (included associated papers 13937 and 14693). J.Pet.Technol.37(2), 323-334) D Coefficient of and epidermisConsider the following equation:
from top to bottomObtaining C consideration in Laplace space D Dimensionless bottom hole pseudo pressure of S
In some embodiments, according to the above embodiments, the present invention can obtain a log-log curve including a dimensionless log curve of bottom hole pressure (hereinafter referred to as a pressure curve) and a dimensionless log curve of bottom hole pressure derivative (hereinafter referred to as a derivative curve), as shown in fig. 4, corresponding to a well location 1 (see fig. 1 or fig. 2) of a gas well in a gas reservoir of a complex-shaped dual pore system, with the pressure curve on the top and the derivative curve on the bottom. Obviously, the problem cannot be solved by using the conventional analytic method or the semi-analytic method.
As can be seen from fig. 4, the percolation process can be divided into 9 percolation stages: the I stage is a pure well reservoir section, and the pressure curve and the derivative curve of the section both form a straight line with the slope of 1; stage II is the transition section, and the derivative curve is an upward "hump"; stage III is a bi-linear flow segment, the pressure derivative curve of which is a '1/4' slope line; stage IV is a linear flow segment with a pressure derivative curve that follows a "1/2" slope line; the V stage is a substrate channeling section towards a natural fracture, a downward concave seed is arranged on a derivative curve of the section, and the depth and the morning and evening of the concave seed are respectively related to the values of a storage capacity ratio omega and a channeling coefficient lambda; stage VI is a pseudo-radial flow section, the pressure derivative curve of which is a horizontal line of 0.5 height; the VII section is a right linear boundary reflecting section closest to the well, and the pressure derivative curve of the section rises; the VIII-th segment is a reflection segment when the pressure wave has propagated not only to the right nearest boundary but also to the upper and lower boundaries. The pressure derivative curve of the segment continues to rise; the section IX is the entire outer boundary reflection section, which has propagated the pressure wave to the farthest boundary, and the pressure curve and the pressure derivative curve of which are upwarped and form a straight line with a slope of "1".
In some embodiments, the present invention may be used to obtain a log-log plot of the pressure curve and derivative curve for wells 1,2,3 in a complex-shaped dual-pore system reservoir, as shown in FIG. 5, according to the embodiments described above.
As can be seen from fig. 5, when a gas well is located in different locations of a gas reservoir, the characteristics of its boundary reflection section are different, as follows:
when the gas well is at well location 1, the pressure derivative curve rises earliest because well location 1 is closer to its nearest boundary than well locations 2 and 3 are to their respective nearest boundaries, and the pressure wave first reaches the nearest boundary to the right of well location 1, so its pressure derivative curve rises earliest; as the pressure wave continues to propagate, then reaches the upper and lower boundaries, the derivative curve rises further, and finally the pressure wave reaches the furthest left boundary, with its slope gradually reaching 1, which reflects the entire gas reservoir outer boundary.
Well site 2 is closer to the upper boundary and substantially equidistant from the other boundaries, thus reflecting the influence of the upper boundary first and then transitioning to the entire gas reservoir outer boundary.
Well site 3 is a greater distance from its nearest boundary than well site 1 and well site 2, and therefore reflects the influence of the boundary at the latest.
Although the derivative curves just reflecting the boundary effect are different in magnitude of the upwarp, eventually the derivative curves in all three cases become straight lines with a slope of 1 as the pressure wave propagates across the boundary of the gas reservoir.
The pressure curve also has similar characteristics as described above, except that no derivative curve is apparent.
In some embodiments, the fracture conductivity C shown in FIG. 6 can be obtained according to the above embodiments fD Graph of influence on a log-log curve, the graph pairThe gas well is located at well site 1 in a reservoir of complex-shaped dual pore system.
As can be seen from FIG. 6, C fD The larger the 1/4 slope bilinear segment on the derivative curve the shorter the duration and the lower the position.
In some embodiments, according to the above embodiments, the present invention can obtain the influence graph of the asymmetry of the left and right wings of a fracture on the log curve of a limited diversion vertical fractured well of a complex-shaped gas reservoir as shown in fig. 7, where the gas well is located at well position 1 in the gas reservoir of the complex-shaped dual pore system.
As can be seen from fig. 7, when the two ends of the fracture are asymmetric, the log-log curves are different. For a given total length, the greater the asymmetry between the left and right wings of the fracture, the less pronounced the 1/4 slope versus 1/2 slope linear flow characteristic on the derivative curve. Furthermore, it can be seen that the derivative curve for the left and right wings of the fracture when asymmetric is above the derivative curve for the symmetric.
In addition to the above examples, the inventors further provide the necessity of considering the unequal fracture lengths and the asymmetry of the two-wing flow distribution in the above embodiments, for example, when a fracture model with equal wing lengths as shown in fig. 8 is used, the resulting fracture flow distribution is shown in fig. 9.
As can be seen from fig. 9, even if the length of both wings of the fracturing fracture is equal, the flow distribution of both wings (reflected by the linear density flow of each discrete unit of the fracturing fracture) is symmetrical only at the early stage, and the flow distribution at the early stage has the characteristics of low both ends and high middle. As the pressure wave propagates outward, the flow density at both ends of the fracture will gradually be greater than the flow density in the middle of the fracture; after the pressure wave reaches the outer boundary, the outer boundary will have an influence on the flow distribution of the two wings of the crack, which is specifically shown as follows: the flow density of the fracturing seam on the side close to the closed boundary is greater than that of the fracturing seam on the side far from the closed boundary, and if A and B respectively represent two wings of the fracturing seam, q is B ＞q A 。
From the above, even if the lengths of the two wings of the fracturing fracture are equal, the flow distribution of the two wings at the later stage is also asymmetric due to the influence of the outer boundary, not to mention the situation that the lengths of the two wings of the fracturing fracture are unequal.
The above examples are merely preferred embodiments of the present invention, and the scope of the present invention is not limited to the above examples. All technical schemes belonging to the idea of the invention belong to the protection scope of the invention. It should be noted that modifications and embellishments within the scope of the invention may be made by those skilled in the art without departing from the principle of the invention, and such modifications and embellishments should also be considered as within the scope of the invention.
Claims (9)
1. A fracturing well testing simulation method under complex conditions is characterized by comprising the following steps:
s1, constructing a physical model of a fracturing well, wherein the physical model considers that the outer boundary of a gas reservoir is in any shape, the reservoir is a double-pore system containing a double-pore structure, the left wing and the right wing of a fracture can be symmetrical or asymmetrical, and the fracture has limited flow conductivity; wherein, the dual pore structure means that the pore structure of the reservoir comprises two pore systems, namely a natural fracture system formed by the natural fracture and a matrix pore system formed by the matrix pore;
s2, constructing a gas reservoir seepage master control model of the dual pore system, including:
s21, coupling the inner boundary pressure condition of the zonally converged fracturing with the mass conservation equation of the natural fracture system by utilizing a Dirac generalized function and an integral equation, and simultaneously establishing a motion equation and a state equation of natural gas in the natural fracture system and a channeling equation of the natural gas in the matrix pore system and the natural fracture system to derive a seepage master control differential equation of the natural fracture system;
s22, a gas seepage master control differential equation of the matrix pore system is expressed by using a gas seepage differential equation under a pore structure, and a gas reservoir seepage master control differential equation of the dual pore system, namely the gas reservoir seepage master control model, is formed by the seepage master control differential equation of the natural fracture system and the seepage master control differential equation of the matrix pore system;
s3, constructing a factorial formation seepage model of the dual pore system, comprising:
s31 setting an initial pressure condition equation of the gas reservoir and an outer boundary pressure condition equation under different outer boundary conditions;
s32, combining the initial pressure condition equation and the outer boundary pressure condition equation with the gas reservoir seepage master control differential equation to obtain the factorial formation seepage model;
s4, introducing a dimensionless quantity, and performing dimensionless transformation on the dimensionless stratum seepage model to obtain a dimensionless stratum seepage model;
s5 obtaining a gas reservoir outer boundary seepage model of the dual pore system, comprising:
s51, Laplace transformation is carried out on the dimensionless stratum seepage model, the seepage master control differential equation of the dimensionless matrix pore system is substituted into the seepage master control differential equation of the dimensionless natural fracture system, the pressure parameter of the matrix pore system is eliminated, and the seepage master control differential equation of the transformed natural fracture system is obtained;
s52, converting the transformed seepage master control differential equation of the natural fracture system into an outer boundary seepage integral equation of the gas reservoir based on a boundary element solution to obtain an outer boundary seepage model of the gas reservoir;
s6, carrying out unit discrete processing on the gas reservoir outer boundary seepage model to obtain a first linear equation set;
s7, respectively constructing gas reservoir seepage models considering the influences of limited flow conductivity of the fracturing, unequal lengths of two wings of the fracturing and asymmetric flow distribution of the two wings of the fracturing to obtain fracturing seepage models;
s8, converting the fracturing fracture seepage model into a Fredholm integral equation by using Laplace transformation and double integration, and performing unit discrete processing to obtain a second linear equation set;
s9, the first linear equation set and the second linear equation set are combined to obtain a closed matrix;
s10, solving the closed matrix, and obtaining a simulation result by numerical inversion.
2. The well testing simulation method of claim 1, wherein the coupling of the inner boundary pressure condition of the zonal confluent fractures and the mass conservation equation of the natural fracture system in S2 results in the following coupled mass conservation equation:
where ρ represents the natural gas density in the natural fracture system; v. of x Representing the seepage velocity of natural gas in the x direction in a natural fracture system; v. of y Representing the seepage velocity of natural gas in the y direction in the natural fracture system; x and y respectively represent an x coordinate and a y coordinate; h is the reservoir thickness; f represents a line integral region corresponding to the trajectory of the fracturing fracture; rho sc Represents the density of natural gas at surface conditions; q. q.s fsc A linear density flow function representing the fracture; dl represents the infinitesimal length; w represents any point on the hydraulic fracture, and W is W (x) w ,y w )；δ(x-x w ,y-y w ) Is a two-dimensional Dirac function, where x w 、y w X and y coordinates representing any point on the fracture, respectively, t represents time, phi represents natural fracture system porosity, q * Indicating the amount of cross flow.
3.The well testing simulation method of claim 1, wherein in the S2, the gas reservoir seepage master control model of the dual pore system comprises:
the seepage master control differential equation of the natural fracture system:
seepage master control differential equation for the matrix pore system:
wherein:
representing the pseudo-pressure of the natural fracture system; k represents the gas reservoir permeability of the natural fracture system; h is the reservoir thickness; t is sc Is the temperature at ground conditions; t represents the gas reservoir temperature; w represents any point on the fracture, W ═ W (x) w ,y w )；p sc Representing the ground standard condition pressure; q. q.s fsc A linear density flow function representing the fracture; delta (x-x) w ,y-y w ) Is a two-dimensional dirac function; t represents time; mu.s i Represents the viscosity of natural gas under the original conditions; c gi Representing the compressibility of natural gas under original conditions; α represents a shape factor; k is a radical of m Represents the permeability of the matrix pore system; psi m Representing the pseudo-pressure of the matrix pore system; phi is a m Represents the porosity of the matrix pore system; p represents pressure; p is a radical of 0 Represents a reference pressure; z represents a natural gas deviation factor; μ represents the natural gas viscosity; x and y are x coordinate and y coordinate respectively; x is the number of w 、y w Respectively representing the x and y coordinates of any point on the fracturing fracture; f is a line integral area representing a fracture trajectory; dl is the infinitesimal length.
4. The well testing simulation method of claim 1, wherein in the S3, the external boundary pressure condition equations for the different external boundary situations are as follows:
if the boundary is a closed boundary, the external boundary pressure condition equation is as follows:
wherein Γ represents the gas reservoir outer boundary, p represents the pressure,
is an outward normal vector on the outer boundary;
if the pressure boundary is the constant pressure boundary, the conditional equation of the pressure of the outer boundary is as follows:
p| Γ ＝p i (13)
wherein p is i Representing an evenly distributed original pressure in the gas reservoir;
if the outer boundary is a mixed boundary, the outer boundary pressure condition equation is as follows:
wherein, γ 1 、γ 2 And gamma 2 Are the combination constants.
5. The well testing simulation method of claim 1, wherein in the S3, the outer boundary pressure condition equation is a pseudo-pressure form of the outer boundary pressure condition equation under a closed boundary as follows:
the initial pressure condition equation is in the form of a pseudo pressure as follows:
ψ| t＝0 ＝ψ m | t＝0 ＝ψ i (11)
wherein Γ represents the gas reservoir outer boundary, p represents the pressure,
indicating an outward normal vector on the outer boundary,
the pseudo-pressure of the natural fracture system is shown,
representing the original pseudo-pressure of the gas reservoir, t representing time,. phi m Denotes the pseudo-pressure of the matrix pore system, p 0 Denotes a reference pressure, p i Representing the original pressure of the gas reservoir distributed uniformly, Z represents the natural gas deviation factor, and μ represents the natural gas viscosity.
6. The well testing simulation method of claim 5, wherein in the S4, the dimensionless formation seepage model comprises:
wherein:
wherein psi D Is dimensionless in gas reservoir natural fracture systemsPseudo pressure,. psi Dm Is dimensionless pseudo-pressure,. psi.in the pore system of the gas reservoir matrix i Is the original pseudo pressure of the gas reservoir, omega is the elastic storage-volume ratio, lambda is the cross-flow coefficient, x D Is a dimensionless x coordinate, y D Is a dimensionless y coordinate, q fD Dimensionless linear density flux, dl, for fracturing D Is a dimensionless infinitesimal length,
is a dimensionless outward normal vector, t, on the outer boundary D Denotes dimensionless time, T sc Is the temperature under ground conditions, q fsc Linear density flow function representing fracture, T is gas reservoir temperature, C gi Is the compressibility factor, phi, of natural gas under the original conditions m Is the porosity of the matrix pore system, phi denotes the porosity of the natural fracture system, alpha is the shape factor, k m Is the permeability of the matrix pore system, K is the permeability of the natural fracture system, L ref Is the reference length of the fracture, x fL Is the left wing length of the crack, x fR Is the right wing length of the crack, μ i Is the viscosity of natural gas at the original condition and dl is the infinitesimal length.
7. The well testing simulation method of claim 6, wherein in the S5, the seepage master differential equation of the transformed natural fracture system is as follows:
wherein,
i.e.,(s) represents a function related to the parameters ω, λ and the Laplace variable s; s represents a Laplace variable;
in a gas reservoir natural fracture system after Laplace transformationDimensionless pseudo-pressure of (a);
showing the dimensionless linear density flow of the fracture at the W point after Laplace transformation; f D Representing the dimensionless line integral region, x, corresponding to the trajectory of the fracture wD A dimensionless abscissa representing the W point; y is wD Expressing the dimensionless ordinate of the W point, i.e. W (x) w ,y w )；
The gas reservoir outer boundary seepage model is as follows:
wherein: p 'is any point on the outer boundary of the stratum, Q is any point in the stratum, and G (P', Q, s) represents the basic solution of the boundary element
The middle P point selects the solution at any point P' on the outer boundary gamma, where K is 0 Representing a zero order deformed Bessel function, r D Represents a dimensionless radial distance;
expressing dimensionless simulated pressure of a Q point in the gas reservoir natural fracture system after Laplace transformation;
representing dimensionless simulated pressure of a P' point in the gas reservoir natural fracture system after Laplace transformation; w is any position point on the fracturing crack; θ is a constant related to the geometry at point Q; beta is the inner angle of the left and right tangents to the outer boundary at point Q.
8. The well testing simulation method of claim 6, wherein in the step S7, the fracture seepage model is as follows:
wherein,
wherein psi fD Dimensionless pseudo-pressure for fracturing f Pseudo-pressure, x, for fracturing fLD Is the dimensionless length, x, of the left wing of the crack fRD Dimensionless length of right wing of crack press C fD Is a fracture dimensionless diversion coefficient, K f Is the fracture permeability, q LD Is the dimensionless yield of the left wing of the fracture, q L Is the yield of the left wing of the fracturing fracture, q RD Is dimensionless yield of right flank of fracture, q R Is the yield of the right wing of the fracturing, W fD Is the dimensionless crack width, W f Is the crimp gap width.
9. The well testing simulation method of claim 8, wherein in S9, the closed matrix comprises:
a first set of linear equations:
a second linear system of equations:
wherein N is b Representing the number of discrete cells of the outer boundary Γ of the gas reservoir,
and
is a coefficient of a first linear system of equations, N FL 、N FR Respectively representing the number of discrete elements, Deltax, of the left and right wings of the fracturing fracture DL ＝x fLD /N FL Denotes the left wing discrete element length, Δ x DR ＝x fRD /N FR Indicating the length of the right wing discrete element, i ═ 1,2,3 … … N b Denotes that i takes over N b Discrete units of outer boundary, i ═ 1,2,3 … … N FL +N FR N representing i' Take through crack FL +N FR A discrete unit, k is 1,2,3 … … N FL +N FR N representing k-cut through fracture FL +N FR A discrete unit, Q k The midpoint of the kth discrete cell on the fracture is shown.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210902933.4A CN115114834B (en) | 2022-07-29 | 2022-07-29 | Fracturing well test simulation method under complex condition |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210902933.4A CN115114834B (en) | 2022-07-29 | 2022-07-29 | Fracturing well test simulation method under complex condition |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN115114834A true CN115114834A (en) | 2022-09-27 |
| CN115114834B CN115114834B (en) | 2024-02-23 |
Family
ID=83334254
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210902933.4A Active CN115114834B (en) | 2022-07-29 | 2022-07-29 | Fracturing well test simulation method under complex condition |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN115114834B (en) |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116629154A (en) * | 2023-05-24 | 2023-08-22 | 西南石油大学 | Fractal composite gas reservoir fracturing well transient pressure calculation method, system and equipment |
| CN116861818A (en) * | 2023-07-21 | 2023-10-10 | 西南石油大学 | Multilayer gas reservoir pressure-dividing and pressure-mixing test well simulation method under complex condition |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111927420A (en) * | 2020-08-15 | 2020-11-13 | 西南石油大学 | Method for simulating pressure of asymmetric fractured well with limited diversion for gas reservoir in any shape |
| CN113743037A (en) * | 2021-09-15 | 2021-12-03 | 陕西延长石油（集团）有限责任公司 | Low-permeability reservoir water injection induced dynamic fracture variable flow conductivity calculation method |
| CN114201932A (en) * | 2021-12-10 | 2022-03-18 | 西南石油大学 | Well testing simulation method for tight reservoir fracturing well under complex condition |
-
2022
- 2022-07-29 CN CN202210902933.4A patent/CN115114834B/en active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111927420A (en) * | 2020-08-15 | 2020-11-13 | 西南石油大学 | Method for simulating pressure of asymmetric fractured well with limited diversion for gas reservoir in any shape |
| CN113743037A (en) * | 2021-09-15 | 2021-12-03 | 陕西延长石油（集团）有限责任公司 | Low-permeability reservoir water injection induced dynamic fracture variable flow conductivity calculation method |
| CN114201932A (en) * | 2021-12-10 | 2022-03-18 | 西南石油大学 | Well testing simulation method for tight reservoir fracturing well under complex condition |
Non-Patent Citations (4)
| Title |
|---|
| 何彦锋;孙伟家;符力耘;: "复杂介质地震波传播模拟中边界元法与有限差分法的比较研究", 地球物理学进展, no. 02, 15 April 2013 (2013-04-15), pages 140 - 154 * |
| 王海涛;彭倩;张烈辉;郭晶晶;聂权;: "考虑非达西渗流的复合页岩气藏试井模型", 东北石油大学学报, no. 01, 15 February 2018 (2018-02-15), pages 100 - 106 * |
| 郭大立, 曾晓慧, 赵金洲, 刘慈群: "垂直裂缝井试井分析模型和方法", 应用数学和力学, no. 05, 15 May 2005 (2005-05-15), pages 26 - 32 * |
| 陈军;杨峻懿;刘启国;周丽莎;杨学锋;: "两区复合气藏有限导流垂直裂缝井压力动态分析", 钻采工艺, no. 04, 25 July 2016 (2016-07-25), pages 44 - 46 * |
Cited By (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116629154A (en) * | 2023-05-24 | 2023-08-22 | 西南石油大学 | Fractal composite gas reservoir fracturing well transient pressure calculation method, system and equipment |
| CN116629154B (en) * | 2023-05-24 | 2024-01-09 | 西南石油大学 | Fractal composite gas reservoir fracturing well transient pressure calculation method, system and equipment |
| CN116861818A (en) * | 2023-07-21 | 2023-10-10 | 西南石油大学 | Multilayer gas reservoir pressure-dividing and pressure-mixing test well simulation method under complex condition |
| CN116861818B (en) * | 2023-07-21 | 2024-04-30 | 西南石油大学 | Multilayer gas reservoir pressure-dividing and pressure-mixing test well simulation method under complex condition |
Also Published As
| Publication number | Publication date |
|---|---|
| CN115114834B (en) | 2024-02-23 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| WO2021180189A1 (en) | Multi-element thermal fluid thermal recovery oil reservoir numerical simulation method | |
| CN115114834B (en) | Fracturing well test simulation method under complex condition | |
| CN109441422B (en) | Shale gas well spacing optimization mining method | |
| WO2020224539A1 (en) | Flow simulation and transient well analysis method based on generalized pipe flow seepage coupling | |
| CN107506948B (en) | Shale oil gas comprehensive yield analysis method based on dynamic drainage volume | |
| CN111980654B (en) | Method for calculating capacity of staged fracturing horizontal well of heterogeneous shale oil reservoir | |
| CN109522634B (en) | Numerical analysis method for compact gas multistage volume fracturing horizontal well | |
| Thomas et al. | Fractured reservoir simulation | |
| CN113076676B (en) | Unconventional oil and gas reservoir horizontal well fracture network expansion and production dynamic coupling method | |
| CN113836695B (en) | Oil reservoir numerical simulation method based on gridless connecting element | |
| CN107313759B (en) | Hypotonic heavy crude reservoir straight well thermal recovery pressure distribution forecasting method and system | |
| CN103399970B (en) | The method of digital-to-analogue measuring and calculating oil reservoir flow condition is carried out with the process of discrete fractures line | |
| CN103590824A (en) | Capacity calculation method for compact gas reservoir horizontal well after multi-section fracturing modification | |
| CN114201932B (en) | Compact oil reservoir fracturing well test simulation method under complex condition | |
| CN111734394B (en) | Method for determining unsteady bottom-hole pressure of tight oil reservoir fracturing well | |
| CN114580100B (en) | Method and device for calculating full wellbore pressure of fractured horizontal well and computer readable storage medium | |
| CN105089612A (en) | Determining method for distance of well-drain and length of pressure break of low penetration oil reservoir artificial fracture | |
| CN107462936A (en) | Utilize the method for pressure monitoring Data Inversion low permeability reservoir non-Darcy percolation law | |
| CN111353205A (en) | Method for calculating stratum pressure and dynamic capacity of water-producing gas well of tight gas reservoir | |
| CN112541287A (en) | Loose sandstone fracturing filling sand control production increase and profile control integrated design method | |
| CN111502652A (en) | Yield decreasing and production dynamic prediction method for three-hole medium gas reservoir horizontal well | |
| CN104989385B (en) | The HTHP oil gas straight well perforating parameter optimization method calculated based on skin factor | |
| CN111927420A (en) | Method for simulating pressure of asymmetric fractured well with limited diversion for gas reservoir in any shape | |
| CN115935857A (en) | EDFM-based unconventional oil and gas reservoir capacity rapid simulation method | |
| CN115345090A (en) | Calculation method for dynamic propagation of undersaturated coalbed methane reservoir pressure drop funnel |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |