CN115114834A - A Well Testing Simulation Method for Fracturing Wells in Complex Conditions - Google Patents

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

A Well Testing Simulation Method for Fracturing Wells in Complex Conditions

technical field

The invention relates to the technical field of a fracturing well bottom hole pressure simulation method, in particular to the technical field of a fracturing well testing model building method based on the boundary element method (BEM).

Background technique

Oil and gas well testing technology is a very important technology for understanding oil and gas reservoirs, and is known as "the eyes of oil and gas reservoir development" by the oil industry. Through well testing, some important dynamic parameters of deeply buried oil and gas reservoirs can be obtained, such as formation permeability, wellbore pollution coefficient, etc., and a reasonable unstable well testing model can be established (unstable well testing model is one of the unstable seepage models). species) etc.

Accurate simulation of the bottom hole pressure transient of oil and gas wells is the theoretical basis and necessary premise of well testing analysis technology. 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 regular-shaped reservoirs, while numerical methods can solve more complex problems.

Numerical methods can be divided into two categories: regional methods and boundary methods. Both the finite difference method and the finite element method belong to the first category, while the boundary element method belongs to the second category. Among them, the boundary element method is superior to the regional method in some aspects, such as: (1) the boundary element method retains the analytical characteristics of the understanding to a large extent, so it has higher accuracy than the finite difference, finite element and other regional methods, It can meet the accuracy requirements when solving the well test model; (2) the finite element method and the finite difference method need to divide the elements in the whole area, while the boundary element method only needs to divide the elements on the boundary of the area, so the boundary element method has the advantages of reducing Due to the characteristics of dimension, the order of the formed matrix is also small.

In the existing technology research, the research on using the boundary element method to solve the seepage problem of oil and gas reservoirs mainly focuses on the uncracked wells. In recent years, some scholars have further studied the seepage problem of fracturing wells by using the boundary element method and have made some progress. However, up to now, there is still a lack of accurate construction and simulation of well testing models with arbitrary shapes, fracture structures or flow asymmetry in gas reservoirs, and reservoirs with multi-medium characteristics. Interpretation of fracturing well testing in the above complex situations.

SUMMARY OF THE INVENTION

The purpose of the present invention is to construct a well test model that takes into account the arbitrary shape of the gas reservoir, the physical properties of natural gas at high pressure, the dual medium characteristics of the reservoir, the conductivity of the fracturing fracture, the asymmetry of the two flanks of the fracturing fracture, and the asymmetry of the flow rate of the fracturing fracture. The constructed unstable well test model obtains a simulation method that can accurately calculate the bottom hole pressure transient of the gas well in this complex situation. The simulation method is in Laplace space, the boundary element method (BEM) and Fredholm integral equation Effectively combined, the calculation is efficient and accurate.

The technical scheme of the present invention is as follows:

A well testing simulation method for fracturing wells under complex conditions, comprising:

The construction of S1 considers the physical model of the fracturing well in which the outer boundary of the gas reservoir is of any shape, the reservoir is a dual pore system with dual pore structures, the left and right flanks of the fracturing fracture can be symmetrical or asymmetric, and the fracturing fracture has limited conductivity. , the dual pore structure means that the pore structure of the reservoir includes two kinds of natural fractures and matrix pores, and the entire reservoir contains two pore systems, a natural fracture system formed by natural fractures and a matrix pore system formed by matrix pores;

S2 builds the main control model of gas reservoir seepage of the dual pore system, including:

S21 utilizes Dirac generalized functions and integral equations to couple the inner boundary pressure conditions of the pressure fractures of the strip-shaped sinks with the mass conservation equation of the natural fracture system, and couples it with the motion equations of natural gas in the natural fracture system, The equation of state and the channeling equation of natural gas in the matrix pore system and the natural fracture system are combined, and the seepage master differential equation of the natural fracture system is derived;

S22 uses the gas seepage differential equation under the pore structure to represent the seepage master differential equation of the matrix pore system, which is composed of the seepage master differential equation of the natural fracture system and the seepage master differential equation of the matrix pore system. The main control differential equation of the gas reservoir seepage of the pore system, that is, the main control model of the gas reservoir seepage;

S3 constructs a dimensioned formation seepage model of the dual pore system, including:

S31 sets the initial pressure condition equation of the gas reservoir and the outer boundary pressure condition equation under different outer boundary conditions;

S32 combines the initial pressure condition equation and the outer boundary pressure condition equation with the main control differential equation of the gas reservoir seepage to obtain the dimensioned formation seepage model;

S4 introduces dimensionless quantities, performs dimensionless transformation on the dimensional formation seepage model, and obtains a dimensionless formation seepage model;

S5 obtains the gas reservoir outer boundary seepage model of the dual pore system, including:

S51 performs Laplace transformation on the dimensionless formation seepage model, and substitutes the seepage master differential equation of the dimensionless matrix pore system into the seepage master differential equation of the dimensionless natural fracture system, and eliminates the matrix pore system The pressure parameter of the transformed natural fracture system is obtained to obtain the main control differential equation of seepage;

S52, based on the boundary element solution method, transform the main control differential equation of seepage of the transformed natural fracture system into the integral equation of seepage at the outer boundary of the gas reservoir, and obtain the seepage model of the outer boundary of the gas reservoir;

S6 performs unit discrete processing on the seepage model at the outer boundary of the gas reservoir to obtain a first linear equation system;

S7 builds a gas reservoir seepage model for the two flanks of the fracturing fracture, which takes into account the effects of the limited conductivity of the fracturing fracture, the unequal lengths of the two flanks of the fracturing fracture, and the asymmetric flow distribution of the two flanks, and the seepage model of the fracturing fracture is obtained;

S8 utilizes Laplace transformation and double integration to convert the fracturing seepage model into a Fredholm integral equation, and performs unit discretization processing to obtain a second linear equation system;

S9 Simultaneously combines the first linear equation system and the second linear equation system to obtain a closed matrix;

S10 solves the closed matrix, and uses numerical inversion to obtain a simulation result.

According to some specific embodiments of the present invention, in S2, the coupling of the inner boundary pressure conditions of the pressure fractures of the strip-shaped sinks and the mass conservation equation of the natural fracture system, the following coupling is obtained. Mass conservation equation:

Among them, ρ represents the density of natural gas in the natural fracture system; v _x represents the seepage velocity of natural gas in the natural fracture system in the x direction; v _y represents the seepage velocity of the natural gas in the natural fracture system in the y direction; x and y represent x, respectively Coordinate and y coordinate; h represents the thickness of the reservoir; F represents the line integral area corresponding to the trajectory of the fracture; ρ _sc represents the density of natural gas under ground conditions; q _fsc represents the linear density flow function of the fracture, that is, converted to the surface After the condition, the flow rate from the natural fracture system of the formation vertically flowing into the fractured fracture of unit length is a function of changing position, that is, q _fsc = q _fsc (W); dl is the micro-element length; W is any position on the fractured fracture point, W=W(x _w , y _w ); δ(xx _w , yy _w ) is a two-dimensional Dirac function; x _w , y _w represent the x and y coordinates of any point on the hydraulic fracture, t represents time, φ represents the natural fracture system porosity, and q ^* represents the channeling flow.

According to some specific embodiments of the present invention, in S2, the main control model of gas reservoir seepage of the dual pore system includes:

The seepage control differential equation of the natural fracture system:

The seepage master differential equation of the matrix pore system:

表示天然裂缝系统的拟压力；K表示天然裂缝系统的气藏渗透率；h表示储层厚度；T_sc表示地面条件下的温度；T表示气藏温度；W表示压裂缝上任一位置点，W＝W(x_w,y_w)；p_sc表示地面标况压力；q_fsc表示压裂缝的线密度流量函数；δ(x-x_w,y-y_w)为二维狄拉克函数；t表示时间；μ_i表示天然气在原始条件下的粘度；C_gi表示天然气在原始条件下的压缩系数；α表示形状因子；k_m表示基质孔隙系统渗透率；ψ_m表示基质孔隙系统的拟压力；φ_m是基质孔隙系统孔隙度；p表示压力；p₀表示参考压力；Z是天然气偏差因子；μ是天然气粘度；x、y分别为x坐标和y坐标；x_w、y_w分别表示水力压裂缝上任一点的x及y坐标；F为表示水力压裂缝轨迹的线积分区域；dl是微元长度。in:

represents the pseudo pressure of the natural fracture system; K represents the gas reservoir permeability of the natural fracture system; h represents the thickness of the reservoir; T _sc represents the temperature under ground conditions; T represents the temperature of the gas reservoir; W represents any point on the fracture, W =W(x _w , y _w ); p _sc represents the surface standard pressure; q _fsc represents the linear density flow function of the fracture; δ(xx _w , yy _w ) is a two-dimensional Dirac function; t represents time; μ _i is the viscosity of natural gas under original conditions; C _gi is the compressibility of natural gas under original conditions; α is the shape factor; _{km is the permeability of the matrix pore system; ψ m is the pseudo-pressure of the matrix pore system; φ m is the} matrix pore system System porosity; p represents pressure; p ₀ represents reference pressure; Z is natural gas deviation factor; μ is natural gas viscosity; x and _y are x-coordinate and _y -coordinate respectively; and y-coordinate; F is the line integral area representing the hydraulic fracture trajectory; dl is the element length.

According to some specific embodiments of the present invention, in S3, the outer boundary pressure condition equations under different outer boundary conditions are as follows:

If it is a closed boundary, the outer boundary pressure condition equation is:

是外边界上的向外法向向量；Among them, Γ represents the outer boundary of the gas reservoir, p represents the pressure,

is the outward normal vector on the outer boundary;

If it is a constant pressure boundary, the outer boundary pressure condition equation is:

p| _Γ = p _i (13)

Among them, _pi represents the original pressure evenly distributed in the gas reservoir;

If the outer boundary is a mixed boundary, the outer boundary pressure condition equation is:

Among them, γ ₁ , γ ₂ and γ ₂ are combination constants.

According to some specific embodiments of the present invention, in S3, the outer boundary pressure condition equation is the following pseudo-pressure form of the outer boundary pressure condition equation under the closed boundary:

The initial pressure condition equation is the following pseudo-pressure form:

ψ| _t=0 = ψ _m | _t=0 = ψ _i (11)

表示外边界上的向外法向向量

表示天然裂缝系统拟压力，

表示气藏原始拟压力，t表示时间，ψ_m表示基质孔隙系统的拟压力，p₀表示参考压力，p_i表示气藏中均匀分布的原始压力，Z表示天然气偏差因子，μ表示天然气粘度。Among them, Γ represents the outer boundary of the gas reservoir, p represents the pressure,

represents the outward normal vector on the outer boundary

represents the pseudo-pressure of the natural fracture system,

is the original pseudo pressure of the gas reservoir, t is the time, ψ _m is the pseudo pressure of the matrix pore system, p ₀ is the reference pressure, _pi is the original pressure evenly distributed in the gas reservoir, Z is the natural gas deviation factor, and μ is the natural gas viscosity.

According to some specific embodiments of the present invention, in S4, the dimensionless formation seepage model includes:

in:

是外边界上的无因次向外法向向量，t_D表示无因次时间，T_sc是地面条件下的温度，q_fsc表示压裂缝的线密度流量函数，T是气藏温度，C_gi是天然气在原始条件下的压缩系数，φ_m是基质孔隙系统孔隙度，φ表示天然裂缝系统孔隙度，α是形状因子，k_m是基质孔隙系统渗透率，K是天然裂缝系统渗透率，L_ref是压裂缝参考长度，x_fL是压裂缝左翼长度，x_fR是压裂缝右翼长度，μ_i是天然气在原始条件下的粘度，dl是微元长度。Among them, ψ _D is the dimensionless pseudo pressure in the natural fracture system of the gas reservoir, ψ _Dm is the dimensionless pseudo pressure in the gas reservoir matrix pore system, ψ _i is the original pseudo pressure of the gas reservoir, ω is the elastic storage capacity ratio, λ is the channeling coefficient, x _D is the dimensionless x coordinate, y _D is the dimensionless y coordinate, q _fD is the dimensionless linear density flow of the fracture, dl _D is the dimensionless element length,

is the dimensionless outward normal vector on the outer boundary, t _D is the dimensionless time, T _sc is the temperature under surface conditions, q _fsc is the linear density flow function of the fracture, T is the gas reservoir temperature, C _gi is the compressibility of natural gas under original conditions, φ _m is the matrix pore system porosity, φ is the natural fracture system porosity, α is the shape factor, _km is the matrix pore system permeability, K is the natural fracture system permeability, L _ref is the reference length of the fracturing fracture, x _fL is the length of the left flank of the fracturing fracture, x _fR is the length of the right flank of the fracturing fracture, μ _i is the viscosity of the natural gas under the original condition, and dl is the microelement length.

According to some specific embodiments of the present invention, in S5, the seepage master differential equation of the transformed natural fracture system is as follows:

f(s)表示与参数ω、λ和Laplace变量s有关的变换函数；s表示Laplace变量；

表示经Laplace变换后的气藏天然裂缝系统中的无因次拟压力；

表示经Laplace变换后的W点处的压裂缝无因次线密度流量；F_D表示与压裂缝的轨迹对应的无因次线积分区域，x_wD表示W点的无因次横坐标；y_wD表示W点的无因次纵坐标，即W(x_w,y_w)。in,

f(s) represents the transformation function related to the parameters ω, λ and the Laplace variable s; s represents the Laplace variable;

represents the dimensionless pseudo-pressure in the natural fracture system of the gas reservoir after Laplace transformation;

represents the dimensionless linear density flow of the fracture at point W after Laplace transformation; F _D represents the dimensionless line integral area corresponding to the trajectory of the fracture, x _wD represents the dimensionless abscissa of point W; y _wD Represents the dimensionless ordinate of the W point, namely W(x _w , y _w ).

The seepage model at the outer boundary of the gas reservoir is as follows:

中P点选择外边界Γ上的任意一点P′时的解，其中K₀表示零阶变形贝塞尔函数，r_D表示无因次径向距离；

表示经Laplace变换后的气藏天然裂缝系统中Q点的无因次拟压力；

表示经Laplace变换后的气藏天然裂缝系统中P′点的无因次拟压力；W为压裂缝上任一位置点；θ为与Q点处几何形状有关的常数；β是外边界在Q点处的左右切线的内角。Among them: P' is any point on the outer boundary of the formation, Q is any point in the formation (including the inner and outer boundaries of the region), and G(P', Q, s) represents the basic solution of the boundary element

The solution when the middle point P selects any point P' on the outer boundary Γ, where K ₀ represents the zero-order deformed Bessel function, and r _D represents the dimensionless radial distance;

represents the dimensionless pseudo-pressure at point Q in the natural fracture system of the gas reservoir after Laplace transformation;

represents the dimensionless quasi-pressure at point P' in the natural fracture system of the gas reservoir after Laplace transformation; W is any point on the fracture; θ is a constant related to the geometry at point Q; β is the outer boundary at point Q The inner corners of the left and right tangents at .

均为空间位置点和Laplace变量s的函数，其中，

中空间位置点未定，

中空间位置点为Q点，Q点既可为区域Ω内任一点，也可以为外边界Γ上任一点；

中空间位置点为P′点，P′为外边界Γ上任一点。In the above embodiment,

are functions of the spatial location point and the Laplace variable s, where,

The mid-space location is undetermined,

The mid-space position point is the Q point, and the Q point can be either any point in the region Ω or any point on the outer boundary Γ;

The mid-space position point is P', and P' is any point on the outer boundary Γ.

According to some specific embodiments of the present invention, in the S7, the fracturing seepage model is as follows:

in,

where ψ _fD is the dimensionless pseudo pressure of the fracture, ψ _f is the pseudo pressure of the fracture, x _fLD is the dimensionless length of the left flank of the fracture, x _fRD is the dimensionless length of the right flank of the fracture, and C _fD is the dimensionless length of the fracture Dimensionless conductivity coefficient, K _f is the fracture permeability, q _LD is the dimensionless production of the left flank of the fracture, q _L is the production of the left flank of the fracture, q _RD is the dimensionless production of the right flank of the fracture, q _R is the production of the right flank of the fracture, W _fD is the dimensionless fracture width, and W _f is the fracture width.

According to some specific embodiments of the present invention, in the S9, the closed matrix includes:

The first system of linear equations:

Second system of linear equations:

和

为第一线性方程组系数，N_FL、N_FR分别表示压裂缝左右两翼的离散单元数量，Δx_DL＝x_fLD/N_FL表示左翼离散单元长度为，Δx_DR＝x_fRD/N_FR表示右翼离散单元长度，i＝1,2,3……N_b表示i取遍N_b个外边界离散单元，i'＝1,2,3……N_FL+N_FR表示i'取遍压裂缝的N_FL+N_FR个离散单元，k＝1,2,3……N_FL+N_FR表示k取遍压裂缝的N_FL+N_FR个离散单元，Q_k表示压裂缝上第k个离散单元的中点。where N _b represents the number of discrete units at the outer boundary Γ of the gas reservoir,

and

are the coefficients of the first linear equation system, N _FL and N _FR respectively represent the number of discrete units on the left and right flanks of the fracture, Δx _DL =x _fLD /N _FL represents the length of the left wing discrete unit, Δx _DR =x _fRD /N _FR represents the right wing discrete unit Element length, i=1, 2, 3...N _b means i takes N _b outer boundary discrete elements, i'=1,2,3...N _FL +N _FR means i' takes N of pressure fractures _FL +N _FR discrete units, k=1, 2, 3...N _FL +N _FR means k is the N _FL +N _FR discrete units of the fracture, and Q _k represents the kth discrete unit on the fracture. midpoint.

According to the above simulation methods, a fracturing well testing simulation device under complex conditions can be obtained, including a storage medium storing programs, algorithms and/or models that can implement any of the above simulation methods.

The simulation method of the invention combines the boundary element method (BEM) and the Fredholm integral equation in the Laplace space, and establishes a simulation method for the fracturing gas well in the gas reservoir of any shape, which can comprehensively consider the arbitrary shape of the gas reservoir, the high-pressure physical properties of natural gas, and the reservoir. An unstable well test model influenced by dual medium characteristics, fracture conductivity, asymmetry between the two flanks of the fracture, and asymmetry in the flow distribution of the fracture, and the model was successfully solved. According to the simulation results, high-quality dual The logarithmic curve accurately describes the bottom hole pressure dynamics, which provides important application support for the interpretation of fracturing well testing in this complex situation.

Description of drawings

FIG. 1 is a physical model of a gas reservoir with a dual pore system of arbitrary shape and a gas well containing limited conductivity pressure fractures involved in the specific embodiment.

FIG. 2 is a schematic diagram of the discrete processing of the outer boundary unit of a gas reservoir with a dual pore system of any shape involved in the specific embodiment.

Fig. 3 is a physical model of a finite conductivity pressure fracture involved in the specific embodiment.

FIG. 4 is a double logarithmic curve diagram of a complex-shaped dual-porosity system gas reservoir obtained in a specific embodiment.

FIG. 5 is a double logarithmic curve obtained when the vertical fractured wells are located at different positions in the specific embodiment.

Fig. 6 is a graph showing the influence of the obtained pressure fracture conductivity C _fD on the double logarithmic curve of a vertical fracture well with limited conductivity in a gas reservoir with a complex shape dual pore system in a specific embodiment.

FIG. 7 is a graph showing the influence of the asymmetry of the left and right wings of the pressure fracture obtained in the specific embodiment on the double logarithmic curve of a vertical fracture well with limited conductivity in a gas reservoir with a complex shape dual pore system.

FIG. 8 is a physical model of a fracture with equal lengths of the two wings involved in the specific embodiment.

FIG. 9 is a flow distribution diagram of a fracture with equal lengths of two wings obtained in a specific embodiment.

Detailed ways

The present invention will be described in detail below with reference to the embodiments and drawings, but it should be understood that the embodiments and drawings are only used to describe the present invention by way of example, but do not limit the protection scope of the present invention. All reasonable transformations and combinations included within the scope of the inventive concept of the present invention fall into the protection scope of the present invention.

According to the technical solution of the present invention, some specific embodiments include the following processes:

S1 builds the physical model of the fracturing well considering that the outer boundary of the gas reservoir is of any shape, the reservoir is a dual pore system with dual pore structures, the left and right flanks of the fracturing fracture can be symmetrical or asymmetric, and the fracturing fracture has limited conductivity;

Further, the physical model of the fracturing well can also be set as follows: the gas reservoir undergoes isothermal seepage, the original pressure distribution in the gas reservoir is uniform, and the gas reservoir flow satisfies Darcy's law.

Wherein, the dual pore structure means that the pore structure of the reservoir includes natural fractures and matrix pores, and the entire reservoir contains two pore systems, a natural fracture system formed by natural fractures and a matrix pore system formed by matrix pores. (Dual media system shown in Figure 1).

Some more specific embodiments are, for example, constructing a physical model of a fracturing well as shown in FIG. 1 , in this model, the fracturing well is a vertical gas well formed by hydraulic fracturing, and a pair of artificial fracturing fractures is formed in the reservoir, Its gas reservoir has an arbitrarily shaped outer boundary, the gas well produces at a constant surface production q _sc , and:

The reservoir has a dual pore structure, including natural fracture system and matrix pore system, in which natural fractures are the main flow channels, matrix pores are the main storage space, and matrix pores form the supply source of natural fractures;

The left and right wings of the pressure crack can be symmetrical or asymmetrical, assuming that the length of the left wing is x _fL and the length of the right wing is x _fR ;

The pressure fracture is a vertical fracture with limited flow conductivity with a width W _f ;

The gas reservoir temperature is T, and isothermal seepage is used;

The original pressure distribution in the gas reservoir is uniform, which is _pi ;

The flow of gas reservoirs satisfies Darcy's law.

Based on the above physical model, the simulation method of the present invention further utilizes the boundary element method (BEM), Laplace transform, double integral, and unit discretization of the Fredholm integral equation to establish a fracturing gas well for a gas reservoir with arbitrary shape and double pores. A more comprehensive and stricter unstable well testing model that can comprehensively consider the effects of any shape of the gas reservoir, physical properties of natural gas at high pressure, conductivity of fracturing fractures, asymmetry between the flanks of fracturing fractures, and asymmetry of flow distribution in fracturing fractures, and implement the model. It can solve the problem and draw a high-quality pressure dynamic curve, which provides an important realization method for the well test interpretation of fracturing wells in this complex situation.

Based on the settings of the above fracturing well physical model, proceed as follows:

S2 constructs the main control model of gas reservoir seepage of the dual pore system, which specifically includes:

S21 considers the compressibility of natural gas in the gas reservoir, and uses Dirac generalized functions and integral equations to couple the internal boundary pressure conditions of the pressure fractures of the strip-shaped sinks with the mass conservation equation of the natural fracture system, and combine it with the natural gas in natural The equation of motion, the equation of state in the fracture system and the channeling equation of natural gas in the matrix pore system and the natural fracture system are combined, and the main control differential equation of the natural fracture system is derived;

S22 For the matrix pore system, the traditional gas reservoir seepage differential equation is used as the seepage master differential equation of the matrix pore system, and the seepage master differential equation of the natural fracture system and the seepage master differential equation of the matrix pore system are composed. The main control differential equation of gas reservoir seepage of the system is the main control model of gas reservoir seepage.

S3 constructs a dimensioned formation seepage model of the dual pore system, which specifically includes:

S31 sets the initial pressure condition equation of the gas reservoir and the 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 main control differential equation of the gas reservoir seepage, the formation seepage model is obtained. A Dimensional Formation Seepage Model for Gas Reservoirs of Arbitrary Shape, Dual Pore System.

S4 introduces dimensionless quantities, performs dimensionless transformation on the dimensional formation seepage model, and obtains a dimensionless formation seepage model, which includes: the seepage control differential equation of the dimensionless natural fracture system , the main control differential equation of seepage of the dimensionless matrix pore system, the dimensionless initial pressure condition equation and the dimensionless outer boundary pressure condition equation.

S5 obtains the gas reservoir outer boundary seepage model of the dual pore system, which specifically includes:

S51 Laplace transform is performed on the obtained dimensionless formation seepage model, and the seepage control differential equation of the dimensionless matrix pore system is substituted into the seepage control differential equation of the dimensionless natural fracture system, and the matrix pore system is eliminated. The pressure parameter is used to obtain the seepage master differential equation of the natural fracture system in the region, which eliminates the matrix pressure and couples the boundary condition of the strip-shaped sink such as the fracturing fracture, that is, the seepage master differential equation of the transformed natural fracture system. equation;

S52 , based on the boundary element solution method, transform the seepage master differential equation of the natural fracture system into the outer boundary seepage integral equation of the gas reservoir, and obtain the gas reservoir outer boundary seepage model.

S6 performs unit discretization processing on the obtained gas reservoir outer boundary seepage model to obtain a first linear equation system.

S7 builds a gas reservoir seepage model for the two flanks of the fracturing fracture, which takes into account the limited conductivity of the fracturing fracture, the unequal lengths of the two flanks of the fracturing fracture, and the asymmetry of the flow distribution between the two flanks, and the seepage model of the fracturing fracture is obtained.

S8 uses Laplace transform and double integral to transform the obtained fracture flow model into Fredholm integral equation, and performs element discretization to obtain the second linear equation system.

S9 combines the first linear equation system and the second linear equation system to obtain a closed matrix in which the number of equations is equal to the number of unknowns.

S10 solves the obtained sealing matrix, and uses numerical inversion to draw the double logarithmic curve of the dimensionless bottomhole pressure and the derivative of the dimensionless bottomhole pressure of the fracturing well. According to the obtained double logarithmic curve and the seepage mechanism, the gas The reservoir flow stage is divided, and the position of the fracturing well in the complex-shaped gas reservoir, the dimensionless conductivity coefficient of the fracturing fracture, and the influence of the asymmetry of the two flanks of the fracturing fracture on the double logarithmic curve are analyzed to obtain the simulation results.

Further, in some specific embodiments, the seepage control differential equation of the natural fracture system is constructed as follows:

According to the principle of mass conservation, combined with the properties of the δ generalized function and the connotation of the integral, the mass conservation equation in the natural fracture system coupled with the boundary conditions of the vertical fracture well is obtained, as follows:

In order to distinguish the natural fracture system from the matrix pore system and artificial fractures, the pressure and other parameters of the natural fracture system are not subscripted, the pressure and other parameters of the matrix pore system are all subscripted "m", and the artificial fractures are Parameters such as pressure are all subscripted "f".

Among them, the mass conservation equation of the matrix pore system still adopts the traditional form:

The channeling equation of natural gas in matrix pore system and natural fracture system is as follows:

Among them, ρ and ρ _m are the natural gas density in the natural fracture system and the matrix pore system, respectively, kg/(m ³ ); v _x is the seepage velocity of natural gas in the natural fracture system in the x direction, m/s; v _y is The seepage velocity of natural gas in the natural fracture system in the y direction, m/s; x is the x coordinate, m; y is the y coordinate, m; h is the thickness of the reservoir, m; F represents the line integral corresponding to the trajectory of the pressure fracture area; ρ _sc is the density of natural gas under ground conditions, kg/(m ³ ); q _fsc is the linear density flow rate of the fracturing, that is, the flow rate from the natural fracture system of the formation vertically flowing into the hydraulic fracture per unit length (converted to ground conditions ), is a function that changes with position, that is, q _fsc = q _fsc (W), the unit is m ³ /(sm); dl is the micro-element length, m; W is any position point on the fracture, W=W(x _w , y _w ), m ³ /(sm); δ(xx _w , yy _w ) is a two-dimensional Dirac function; φ and φ _m are the porosity of natural fracture system and matrix pore system, respectively, fraction; t is time, s; q ^* is the channeling flow, kg/(m ³ ·s); ρ _r is the natural gas density under the average pressure of the matrix pore system and the natural fracture system, kg/(m ³ ); α is the shape factor, m ^{− 2} ; km is the permeability of the matrix system, _m ² ; μ _r is the natural gas viscosity under the average pressure of the matrix pore system and the natural fracture system, mPa·s.

The equation of motion of gas in natural fracture system adopts Darcy's law, as follows:

Where: K is the permeability of the natural fracture system, m ² ; μ is the natural gas viscosity, Pa·s.

The equation of state of the gas in the natural fracture system is:

pV=nZRT (6)

Where: V is the gas volume, m ³ ; n is the number of moles, kmol; Z is the natural gas deviation factor, dimensionless; R is the universal gas constant, J/(kmol·K); T is the gas reservoir temperature, K.

The above formula can be written as:

Where: M _g is the molecular weight of natural gas, kg/kmol.

经过线性化处理后，可得耦合了压裂缝这种条带状汇内边界条件的任意形状、双重孔隙系统的气藏天然裂缝系统的渗流主控微分方程，如下：Substitute the equation of motion of gas in the natural fracture system, the equation of state of gas in the natural fracture system, and the gas channeling equation (3) between the matrix pore system and the natural fracture system into the above-mentioned natural fracture coupled with the internal boundary conditions of the fracturing fracture System mass conservation equation (1), and considering the influence of high pressure physical properties of natural gas, the pseudo pressure is introduced

After linearization, the main seepage control differential equation of the natural fracture system of the gas reservoir with arbitrary shape and dual pore system coupled with the internal boundary condition of the strip-shaped sink such as fracturing can be obtained, as follows:

Where: K is the permeability of the gas reservoir, m ² ; T _sc is the temperature under the surface condition, K; T is the temperature under the gas reservoir condition, K; μ _i is the viscosity of the natural gas under the original condition, Pa·s; C _gi is the compressibility of natural gas at original conditions, Pa ^-1 .

Further, in some specific embodiments, the seepage control differential equation of the matrix pore system adopts the following traditional form seepage differential equation:

Further, in some specific embodiments, the initial pressure condition equation of the gas reservoir is set as follows:

Assuming that the formation pressure distribution is uniform at the initial moment, which is _pi , the initial pressure condition equation of the gas reservoir is:

p| _t=0 = p _m | _t=0 = p _i (10)

Using quasi-pressure, the above formula can be transformed into:

ψ| _t=0 = ψ _m | _t=0 = ψ _i (11)

Further, in some specific embodiments, the outer boundary pressure condition equation of the gas reservoir is set as follows:

If it is a closed boundary, the outer boundary pressure condition equation is:

是外边界上的向外法向向量。in:

is the outward normal vector on the outer boundary.

If it is a constant pressure boundary, the outer boundary pressure condition equation is:

p| _Γ = p _i (13)

If the outer boundary is a mixed boundary, the outer boundary pressure condition equation is:

Among them: γ ₁ , γ ₂ and γ ₂ are combination constants.

In the above specific embodiments, for a gas reservoir, during the test time, the outer boundary mostly reflects the characteristics of a closed boundary or an infinite boundary, and the situation of a constant pressure boundary rarely occurs. In the mathematical model, when the outer boundary radius is sufficiently large, the closed boundary can be used to describe the pressure characteristics under the infinite boundary. Therefore, in the specific implementation, the outer boundary pressure condition equation ( 12).

Using quasi-pressure, equation (12) can be transformed into:

The above formula is the outer boundary pressure condition equation in the form of quasi-pressure.

The dimensioned stratigraphic mathematical model of the gas reservoir with dual pore system of arbitrary shape coupled with the inner boundary condition of the strip-like sink such as fracturing is composed of equations (8), (9), (11) and (15) .

Further, in some specific embodiments, the dimensionless formation seepage model is constructed as follows:

Define the following quantities:

是外边界上的无因次向外法向向量，无量纲；L_ref是参考长度，可取如裂缝总长度之半，即L_ref＝(x_fL+x_fR)/2，m；x_fL是压裂缝左翼长度，m；x_fR是压裂缝右翼长度，m；q_fD为压裂缝无因次线密度流量，无量纲；ω是弹性储容比，无量纲；λ是窜流系数，无量纲。Among them, ψ _D is the dimensionless pseudo pressure in the natural fracture system of the gas reservoir, dimensionless; ψ _Dm is the dimensionless pseudo pressure in the gas reservoir matrix pore system, dimensionless; ψ _i is the original pseudo pressure of the gas reservoir, Pa /s; p ₀ is the reference pressure, Pa; p _f is the pressure in the fracture, Pa; x _D is the dimensionless x coordinate, dimensionless; y _D is the dimensionless y coordinate, dimensionless; dl _D is the dimensionless Dimensional element length, dimensionless;

is the dimensionless outward normal vector on the outer boundary, dimensionless; L _ref is the reference length, which can be taken as half of the total length of the crack, that is, L _ref =(x _fL +x _fR )/2, m; x _fL is Length of the left flank of the fracture, m; x _fR is the length of the right flank of the fracture, m; q _fD is the dimensionless linear density flow of the fracture, dimensionless; ω is the elastic storage capacity ratio, dimensionless; λ is the channeling coefficient, dimensionless .

Then equations (8), (9), (11) and (15) can be transformed into the following dimensionless models:

The equations (16) to (19) are the dimensionless formation seepage model of the gas reservoir with arbitrary shape dual pore system coupled with the inner boundary condition of the strip-like sink such as fracturing. Among them, Equation (16) is the seepage master differential equation of the dimensionless natural fracture system, Equation (17) is the seepage master differential equation of the dimensionless matrix pore system, and Equation (18) is the The dimensionless initial pressure condition equation, Equation (19) is the dimensionless outer boundary pressure condition equation.

Further, in some specific embodiments, the gas reservoir outer boundary seepage model is constructed as follows:

(1) Introducing the Laplace transform based on t _D , the dimensionless formation seepage model can be transformed into:

Among them: ▽ ² is the Laplace operator, there are:

From (21) we get:

Further, based on the idea of solving problems with boundary elements, the differential equations in the region are turned into integral equations on the boundary, and then the boundary can be divided into boundary elements, and the integral equations on the boundary can be discretized into linear algebraic equations, which will be solved. The problem of differential equations is transformed into the problem of solving algebraic equations, specifically, including:

Substituting equation (21) into equation (20), we get:

in:

Then, Equation (25) is the master differential equation governing the seepage of the natural fracture system in the region where the pressure parameters of the matrix pore system are eliminated and the boundary conditions within the strip-like sink such as fracturing fractures are coupled.

(2) Transform equation (25) into the seepage integral equation on the outer boundary of a gas reservoir with a dual pore system of any shape, as follows:

Suppose the basic solution G(P,Q,s) when applying the boundary element method satisfies the following equation:

▽ ² G(P,Q,s)-f(s)G(P,Q,s)+πδ(P,Q)＝0 (26)

Among them: P and Q are any two points in the formation.

The basic solution of the above equation is as follows:

After a series of derivations, the integral equation of formula (25) can be obtained as:

In the above formula, P' is any point on the boundary Γ, W is the well point, and Q point is any point in the region Ω.

If the Q point can be taken not only in the region Ω, but also on the boundary Γ, the above formula can be replaced by the following formula to obtain the gas reservoir outer boundary seepage model:

where: θ is a constant related to the geometry at point Q,

where: β is the interior angle of the left and right tangents of the outer boundary at point Q.

Further, in some specific embodiments, referring to FIG. 2 , the process of obtaining the first linear equation system by performing cell discretization processing on the gas reservoir outer boundary seepage model includes:

The outer boundary Γ of the gas reservoir is divided into N _b discrete units, the endpoints of each discrete unit are nodes, and the pressure distribution on the outer boundary is represented by a linear discrete unit, that is, the pressure value on any discrete unit on the outer boundary is represented by the The pressure values at the nodes at both ends of the element are obtained by interpolation, and the outer boundary can be specified to take the counterclockwise direction as the positive direction in the calculation.

The left and right flanks of the fracture are divided into _NFL and _NFR discrete units, respectively, and the flow distribution on the fracture is represented by a constant discrete unit, that is, the linear density flow of each part on the same discrete unit of the fracture is equal, then the formula ( 29) can be transformed into the following first system of linear equations:

和

为线性方程组系数。in:

and

are the coefficients of the system of linear equations.

共(N_b+N_FL+N_FR)个，方程个数小于未知数个数，还不能求解，还需要其他方程。Taking k in equation (31) through the nodes of the outer boundary (k=1, 2,...N _b ), N _b linear algebraic equations can be obtained, and the current unknowns include

There are (N _b +N _FL +N _FR ) in total, and the number of equations is less than the number of unknowns, so it cannot be solved yet, and other equations are needed.

Using the quasi-pressure, if the point Q _k is taken in the research domain instead of the outer boundary, the dimensionless quasi-pressure solution of any point Q in the research domain is:

显然，此处可将

简写为

其中k＝N_b+1,N_b+2,……，N_b+N_FL+N_FR。式(31)和式(32)中的未知数个数现在共有[N_b+2(N_FL+2N_FR)]个，分别为

和

由于方程个数小于未知数个数，故还不能求解，还需要联立其它方程，因此，其后引入有限导流的压裂缝渗流模型。When the Q _k point in equation (32) is taken at the midpoint of the fracture, (N _FL +N _FR ) linear algebraic equations can be obtained. Therefore, equations (31) and (32) together represent N _b +N _FL +N _FR linear algebraic equations. However, at this time, a new unknown number of N _FL +N _FR is added, that is, the dimensionless pseudo-pressure at the midpoint of the fracture

Obviously, here

abbreviated as

where k=N _b +1, N _b +2, . . . , N _b +N _FL +N _FR . The number of unknowns in equations (31) and (32) is now [N _b +2(N _FL +2N _FR )], which are respectively

and

Since the number of equations is less than the number of unknowns, it cannot be solved yet, and other equations need to be established simultaneously. Therefore, a fractured seepage model with limited conductivity is introduced later.

Further, in some specific embodiments, referring to FIG. 3 , the construction process of the fracturing seepage model includes:

Define the following quantities:

where ψ _fD is the dimensionless quasi-pressure of the fracture, dimensionless; ψ _f is the quasi-pressure of the fracture, Pa/s; x _fLD is the dimensionless length of the left flank of the fracture, dimensionless; x _fRD is the dimensionless length of the right flank of the fracture Dimensional length, dimensionless; x _fL is the length of the left flank of the fracturing fracture, m; x _fR is the length of the right flank of the fracturing fracture, m; C _fD is the dimensionless conductivity coefficient of the fracturing fracture; K _f is the permeability of the fracturing fracture, m ² ; q _LD is the dimensionless production of the left flank of the fracture, dimensionless; q _L is the production of the left flank of the fracture, m ³ /s; q _RD is the dimensionless production of the right flank of the fracture, dimensionless; q _R is the right flank of the fracture The output of , m ³ /s; W _fD is the dimensionless fracture width, dimensionless; W _f is the fracture width, m.

According to the above-defined quantities, the fracture flow model considering the limited conductivity of the fracture, the unequal lengths of the two wings of the fracture, and the asymmetric flow distribution of the two wings is obtained, as follows:

Further, in some specific embodiments, the obtaining of the second system of linear equations includes:

Laplace transform is performed on the fracture seepage model composed of equations (33) to (39), and double integration is used to obtain the Fredholm integral equations describing the pressure distribution of the left and right flanks of the fracture as follows:

为Laplace空间中的无因次井底拟压力，无量纲；其它上部带横线的量都是Laplace空间中的相应量。In the formula,

is the dimensionless bottom hole pseudo pressure in Laplace space, dimensionless; other quantities with horizontal lines in the upper part are corresponding quantities in Laplace space.

为第j个裂缝离散单元的中点，x_Dj为第j个端点(节点)，在同一离散单元上的线密度流量是均匀的，则对于裂缝左翼上的第k个离散单元(1≤k≤N_FL)，式(40)变为：As mentioned above, the left and right flanks of the fracture are divided into N _FL and N _FR discrete units respectively, and the entire fracture is divided into N _FL +N _FR discrete units, then the length of the left-wing discrete unit is Δx _DL = x _fLD / N _FL , the right-wing discrete unit length is Δx _DR =x _fRD /N _FR . Assumption

is the midpoint of the jth fracture discrete unit, x _Dj is the jth endpoint (node), and the linear density flow on the same discrete unit is uniform, then for the kth discrete unit on the left flank of the fracture (1≤k ≤N _FL ), equation (40) becomes:

Similarly, for the kth discrete unit on the right flank of the crack (N _FL +1≤k≤N _FL +N _FR ), Equation (41) becomes:

为第k个单元的中点处的无因次拟压力，无量纲。in,

is the dimensionless quasi-pressure at the midpoint of the kth element, dimensionless.

与式(32)中的

的对应关系如下：Taking k in Equation (42) and Equation (43) as the midpoint of the discrete unit of pervasive pressure fractures, (N _FL +N _FR ) linear algebraic equations can be obtained. here

and in formula (32)

The corresponding relationship is as follows:

而式(42)和(43)却代表了(N_FL+N_FR)个线性代数方程。Therefore, in fact, the unknowns in equations (42) and (43) are only three more based on the previous equations (31) and (32):

However, equations (42) and (43) represent (N _FL +N _FR ) linear algebraic equations.

In addition, there are flow conditions:

Equations (44) to (47) represent three linear algebraic equations, that is, the second linear equation system.

Therefore, equations (42), (43) and (45)-(47) represent (N _FL +N _FR +3) linear algebraic equations in total.

Further, in some specific embodiments, the obtaining of the closed matrix includes:

方程个数与未知数个数相等，且均为线性代数方程，故可以求解。其联立获得的矩阵为小型稠密矩阵，可采用高斯消元法求解。It can be seen from the above that equations (31), (32), (42), (43), (45) to (47) represent [N _b +2(N _FL +N _FR )+3] linear algebras in total equation, and the number of unknowns is also {N _b +2(N _FL +N _FR )+3}, which are

The number of equations is equal to the number of unknowns, and they are both linear algebraic equations, so they can be solved. The matrices obtained simultaneously are small dense matrices, which can be solved by the Gaussian elimination method.

The above-mentioned embodiment of the present invention avoids the dependence on the finite element method or the finite difference method, and the obtained model comprehensively considers the arbitrary shape of the gas reservoir, the physical properties of natural gas at high pressure, the dual medium characteristics of the reservoir, the conductivity of the fracturing fracture, and the length of the two flanks of the fracturing fracture. Inequality and the influence of asymmetric flow distribution on both wings.

未考虑井储效应和表皮效应的影响。可根据Van Everdingen等人(Van-Everdingen,A.F.,1953.The skin effect and its influenceon the productive capacity of a well.J.Pet.Technol.5(6),171–176.；Kucuk,F.,Ayestaran,L.,1985.Analysis of simultaneously measured pressure and sandfaceflow rate in transient well testing(includes associated papers 13937 and14693).J.Pet.Technol.37(2),323–334)的研究，对无因次井储系数C_D和表皮系数S通过下式考虑:Further, by the above method, the

The well-reservoir effect and the skin effect are not considered. According to Van Everdingen et al. (Van-Everdingen, AF, 1953. The skin effect and its influence on the productive capacity of a well. J. Pet. Technol. 5(6), 171–176.; Kucuk, F., Ayestaran , L., 1985. Analysis of simultaneously measured pressure and sandfaceflow rate in transient well testing (includes associated papers 13937 and 14693). J.Pet.Technol.37(2), 323–334) Research on dimensionless well storage The coefficient C _D and the skin coefficient S are considered by the following equations:

From the above, the dimensionless bottomhole pseudo-pressure considering CD and _S in Laplace space can be obtained

In some specific embodiments, according to the above-mentioned embodiments, the present invention can obtain the logarithmic curve of the dimensionless bottom hole pressure (hereinafter referred to as the pressure curve) and the logarithmic curve of the non-dimensional bottom hole pressure derivative as shown in FIG. 4 . (hereinafter referred to as the derivative curve) double logarithmic curve, the gas well corresponding to the double logarithmic curve is located in the well position 1 (see Figure 1 or Figure 2) in the gas reservoir with complex shape double pore system, the upper part is the pressure curve, and the lower part is the Derivative curve. Obviously, this problem cannot be solved by traditional analytical or semi-analytical methods.

It can be seen from Figure 4 that the seepage process can be divided into 9 seepage stages: the first stage is the pure well storage section, the pressure curve and the derivative curve of this section are straight lines with a slope of 1; the second stage is the transition section, the derivative The curve is an upward "hump"; the third stage is a bilinear flow section, and the pressure derivative curve of this section is a straight line with a "1/4" slope; the fourth stage is a linear flow section, and the pressure derivative curve of this section is "1/2". "Slope straight line; Stage V is the channeling section of matrix to natural fractures, there is a downward depression on the derivative curve of this section, and the depth and morning and evening of the depression are related to the value of storage capacity ratio ω and channeling coefficient λ respectively; Stage VI is the quasi-radial flow section, and the pressure derivative curve of this section is a horizontal line with a height of 0.5; Stage VII is the right straight boundary reflection section closest to the well, and the pressure derivative curve of this section rises; Section VIII is The pressure wave propagates not only to the nearest boundary on the right, but also to the reflection segments when the upper and lower boundaries are present. The pressure derivative curve of this section continues to rise; the IX section is the reflection section of the entire outer boundary, the pressure wave in this section has propagated to the farthest boundary, and the pressure curve and the pressure derivative curve are upturned, showing a straight line with a slope of "1".

In some specific embodiments, according to the above-mentioned embodiments, the present invention can obtain a double logarithmic graph including a pressure curve and a derivative curve as shown in FIG. Well positions 1, 2, and 3 in the gas reservoir.

It can be seen from Figure 5 that when the gas well is located in different positions of the gas reservoir, the characteristics of its boundary reflection section are also different, as follows:

When the gas well is located at well position 1, the pressure derivative curve rises the earliest, because the distance from well position 1 to its nearest boundary is smaller than the distances from well positions 2 and 3 to their respective nearest boundaries, and the pressure wave first reaches the right of well position 1 As the pressure wave continues to propagate, and then reaches the upper and lower boundaries, the derivative curve further rises, and finally the pressure wave reaches the farthest left boundary, and its slope gradually reaches 1, which reflects is the outer boundary of the entire gas reservoir.

Well position 2 is closer to the upper boundary, and the distance to other boundaries is basically the same, so the influence of the upper boundary is reflected first, and then it transitions to the outer boundary of the entire gas reservoir.

Well location 3 is farther from its nearest boundary than well location 1 and well location 2, so it reflects the influence of the boundary at the latest.

Although the magnitude of the upward warping of the derivative curves at the time of reflecting the influence of the boundary is different, when the pressure wave propagates to the entire gas reservoir boundary, the derivative curves in all three cases become a straight line with a slope of 1.

The pressure curve has similar characteristics as above, but it is not as obvious as the derivative curve.

In some specific embodiments, according to the above-mentioned embodiments, the present invention can obtain the influence diagram of the fracture conductivity C _fD on the double logarithmic curve as shown in FIG. Well location 1 in the gas reservoir of the system.

It can be seen from Figure 6 that the larger the C _fD , the shorter the duration of the 1/4 slope bilinear segment on the derivative curve, and the lower the position.

In some specific embodiments, according to the above embodiments, the present invention can obtain the influence diagram of the asymmetry of the left and right wings of the pressure fracture on the double logarithmic curve of the vertical fracture well with limited conductivity in the complex shape gas reservoir as shown in FIG. 7 , The gas well corresponding to this figure is located at well position 1 in a gas reservoir with a complex shape dual pore system.

It can be seen from Figure 7 that when the two ends of the fracture are asymmetrical, the double logarithmic curve is also different. Under the condition of a certain total length, if the asymmetry of the left and right flanks of the fracture is stronger, the 1/4 slope bilinear flow characteristics on the derivative curve will be less obvious, while the 1/2 slope linear flow characteristics will be more obvious. In addition, it can be seen that the derivative curve when the left and right flanks of the fracture are asymmetric is above the derivative curve when they are symmetrical.

In addition to the above examples, the inventor further provides the necessity of considering the unequal fracturing lengths and the asymmetric flow distribution between the two wings in the above-mentioned embodiments. In the model, the obtained fracture flow distribution is shown in Figure 9.

It can be seen from Figure 9 that even if the lengths of the two flanks of the fracture are equal, the flow distribution of the two flanks (reflected by the linear density flow of each discrete unit of the fracture) is only symmetrical in the early stage, and the flow distribution in the early stage shows a lower flow rate at both ends. , the characteristics of the middle high. With the outward propagation of the pressure wave, the flow density at both ends of the fracture will gradually be greater than the flow density in the middle of the fracture; when the pressure wave reaches the outer boundary, the outer boundary will have an impact on the flow distribution on both sides of the fracture, specifically as follows: The flow density of the pressure fracture on the side of the closed boundary is greater than the flow density of the pressure fracture on the side far from the closed boundary. Assuming that A and B are used to represent the two wings of the pressure fracture, then q _B > q _A .

It can be seen from the above that even if the lengths of the two flanks of the fracture are equal, due to the influence of the outer boundary, the flow distribution of the two flanks in the later stage is asymmetric, not to mention the situation that the lengths of the two flanks of the fracture are not equal.

The above embodiments are only preferred embodiments of the present invention, and the protection scope of the present invention is not limited to the above embodiments. All the technical solutions under the idea of the present invention belong to the protection scope of the present invention. It should be pointed out that for those skilled in the art, improvements and modifications without departing from the principles of the present invention should also be regarded as the protection scope of the present 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 ₀ 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, γ ₁ 、γ ₂ And gamma ₂ 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 ₀ 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 ₀ 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.

Images