CN111832164B

CN111832164B - Shale gas yield prediction method, device and equipment

Info

Publication number: CN111832164B
Application number: CN202010631132.XA
Authority: CN
Inventors: 黄世军; 程林松; 吴永辉; 王鹏; 杨彦伟
Original assignee: China University of Petroleum Beijing
Current assignee: China University of Petroleum Beijing
Priority date: 2020-07-03
Filing date: 2020-07-03
Publication date: 2023-05-12
Anticipated expiration: 2040-07-03
Also published as: CN111832164A

Abstract

The embodiment of the specification discloses a shale gas yield prediction method, a shale gas yield prediction device and shale gas yield prediction equipment, wherein the method comprises the steps of obtaining a simulated seepage model corresponding to a target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit comprises a matrix; calculating pressure values on a region unit and a fracture unit which are connected with a designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function; determining the flow rate of a regional unit and a fracture unit connected with a designated shale gas well of the target reservoir to the designated shale gas well according to the calculated pressure value; and predicting the yield of the designated shale gas well according to the flow rate. Yield prediction efficiency and accuracy can be improved using the various embodiments of the present disclosure.

Description

Shale gas yield prediction method, device and equipment

Technical Field

The specification relates to the technical field of petroleum and natural gas exploration, in particular to a shale gas yield prediction method, device and equipment.

Background

Shale gas wells generally have no natural productivity due to the low porosity of the shale matrix, low permeability and poor reservoir physical properties. Currently, the industrial exploitation of shale gas generally utilizes long horizontal well drilling technology and large-scale hydraulic fracturing technology. After large-scale hydraulic fracturing, shale reservoirs form a very complex fracture network, so accurate seepage simulation of shale gas in a dense matrix and the complex fracture network is very critical to efficient development of shale gas.

At present, a discrete fracture numerical simulation method is mostly adopted to simulate seepage of a fractured hydrocarbon reservoir. The basic idea of the discrete fracture numerical simulation method is to simulate the fracture more realistically by encrypting the grid, using a narrower grid width to simulate the fracture explicitly, or subjecting the fracture to a dimension reduction process. The existing discrete fracture numerical simulation method is mainly divided into a discrete fracture model and an embedded discrete fracture model.

The discrete fracture model is mainly characterized by encrypting grids, so that the problem of larger pressure and saturation gradient near the fracture can be simulated. However, the discrete fracture model has a large number of computational grids and is very computationally inefficient. The embedded discrete fracture model adopts a background grid to mesh the reservoir, then the dimension-reduced fracture is directly embedded into the background grid, fracture units are split according to the intersection condition of the background grid and the fracture, and then the conduction coefficient between the fracture units and the background grid is calculated to simulate the complex fracture. However, in the numerical simulation of the embedded discrete fracture model, finer computational grids are also required, especially early solving accuracy, so that the embedded discrete fracture model is often required to be weighted in accuracy and efficiency, the computational accuracy of the sparse grid is lower, and the computational efficiency of the fine grid is lower. Therefore, there is a need for an efficient development scheme for determining shale gas more accurately and efficiently by combining computational accuracy and efficient seepage simulation method.

Disclosure of Invention

The embodiment of the specification aims to provide a shale gas yield prediction method, device and equipment, which can improve the efficiency and accuracy of shale gas yield prediction.

The specification provides a shale gas yield prediction method, device and equipment, which are realized in the following modes:

a shale gas production prediction method, comprising: obtaining a simulated seepage model corresponding to a target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit includes a substrate. Calculating pressure values of a region unit and a fracture unit connected with a designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function. And determining the flow rate of the regional unit and the fracture unit connected with the appointed shale gas well of the target reservoir to the appointed shale gas well according to the calculated pressure value. And predicting the yield of the designated shale gas well according to the flow rate.

In other embodiments of the methods provided herein, the zone unit further comprises an equivalent continuous medium, and the simulated seepage model is constructed accordingly by: and carrying out continuous medium equivalent treatment on the flow characteristic parameters of the first cracks with the crack dimensions smaller than a preset threshold value in the target reservoir to obtain an equivalent continuous medium. And carrying out fusion treatment on the flow characteristic parameters of the equivalent continuous medium and the matrix of the target reservoir, and taking a medium model represented by the fused flow characteristic parameters as a reservoir background model. And carrying out discrete medium equivalent treatment on the flow characteristic parameters of the second cracks with the crack dimensions larger than or equal to a preset threshold value in the target reservoir to obtain discrete cracks. And carrying out subdivision treatment on the reservoir background model and the discrete cracks to obtain the simulated seepage model.

In other embodiments of the methods provided herein, the dissecting the reservoir background model and discrete fractures comprises: and splitting the reservoir background model based on a preset grid step length to obtain a plurality of area units. And embedding the discrete cracks into the area units, and carrying out unit division on the discrete cracks according to the intersection conditions between the discrete cracks and the boundary lines between the discrete cracks and the area units to obtain a plurality of crack units.

In other embodiments of the methods provided herein, the method further comprises: and carrying out encryption treatment on the crack unit to obtain the encrypted crack unit. Correspondingly, calculating the pressure values of a region unit connected with the designated shale gas well of the target reservoir and the encrypted fracture unit by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating the pressure values of the central point of the boundary line of the area unit and the point on the encrypted crack unit based on the unsteady-state green function.

In other embodiments of the methods provided herein, the method further comprises: encrypting the area unit where the crack unit is located based on a preset encryption method to obtain an encrypted area unit; the preset encryption method comprises the step of equally dividing the area unit where the crack unit is located or the step of equally dividing the boundary line of the area unit where the crack unit is located. And re-dividing the discrete cracks according to the intersecting conditions between the discrete cracks and the boundary lines of the encrypted area units and obtaining encrypted crack units. Correspondingly, calculating the pressure values of a region unit, an encrypted region unit and an encrypted fracture unit which are connected with a designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating the pressure values of the central points of the boundary lines of the area unit and the encrypted area unit and the points on the encrypted crack unit based on the unsteady-state green function.

In other embodiments of the methods provided herein, the preset pressure calculation model includes:

wherein ,

for solving the pressure value at point i in zone element e, θ _i For the angle value of point i, n represents the normal vector of point i, +>

Pressure value representing the center point of the j-th boundary line of the area unit e, < >>

The kth crack embedded in zone unit eFlow source and sink item per unit length of seam unit, < ->

Flow source sink item per unit length for jth boundary line of zone element e, +.>

For the flow source sink at shale gas well, +.>

For the permeability of the zone unit e +.>

For the total number of boundary lines of the area unit e, +.>

Represents the total number of crack units embedded in the area unit e, < >>

Representing the number of shale gas wells in zone unit e, m being the mth shale gas well in zone unit e, Γ ^bej The j-th boundary line Γ representing the area element e ^Fek Represents the kth crack element in the e-region element,/->

The green's basic solution is the pressure value at any point in the area element e.

In another aspect, embodiments of the present disclosure further provide a shale gas yield prediction apparatus, including: the seepage model acquisition module is used for acquiring a simulated seepage model corresponding to the target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit includes a substrate. The calculation module is used for calculating pressure values of a region unit and a fracture unit connected with the designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function. And the flow determining module is used for determining the flow of the area unit and the fracture unit connected with the designated shale gas well of the target reservoir to the designated shale gas well according to the calculated pressure value. And the yield prediction module is used for predicting the yield of the specified shale gas well according to the flow.

In other embodiments of the apparatus provided herein, the apparatus further comprises a model building module, wherein the model building module comprises: and the first equivalent unit is used for carrying out continuous medium equivalent treatment on the flow characteristic parameters of the first cracks with the crack dimensions smaller than a preset threshold value in the target reservoir to obtain an equivalent continuous medium. And the fusion unit is used for carrying out fusion treatment on the equivalent continuous medium and the flow characteristic parameters of the matrix of the target reservoir, and taking a medium model represented by the fused flow characteristic parameters as a reservoir background model. And the second equivalent unit is used for carrying out discrete medium equivalent treatment on the flow characteristic parameters of the second cracks with the crack dimensions larger than or equal to a preset threshold value in the target reservoir to obtain discrete cracks. And the model construction unit is used for carrying out subdivision treatment on the reservoir background model and the discrete cracks to obtain a simulated seepage model.

In other embodiments of the apparatus provided in the present specification, the model building unit includes: and the subdivision subunit is used for subdividing the reservoir background model based on a preset grid step length to obtain a plurality of area units. The crack dividing unit is used for embedding the discrete cracks into the area unit, and dividing the discrete cracks into units according to the intersecting condition between the discrete cracks and the boundary line between the discrete cracks and the area unit to obtain a plurality of crack units.

In another aspect, embodiments of the present disclosure further provide a shale gas production prediction apparatus, including a processor and a memory for storing processor-executable instructions, which when executed by the processor implement the steps of the method of any one or more of the embodiments described above.

Shale gas production prediction methods, apparatus, and devices provided by one or more embodiments of the present disclosure may construct a simulated seepage model characterizing heterogeneous flow characteristics of a reservoir by separately characterizing a matrix and a fracture. Then, a solution point of the simulated seepage model can be set at a central point of the unit, and the pressure distribution in the simulated seepage model can be solved based on an unsteady state green function method. By the aid of the solving mode, data processing complexity can be reduced, and data processing efficiency is greatly improved on the basis of guaranteeing solving accuracy. And then, determining flow distribution based on pressure distribution, and further accurately and efficiently predicting the yield of the shale gas well.

Drawings

In order to more clearly illustrate the embodiments of the present description or the technical solutions in the prior art, the drawings that are required in the embodiments or the description of the prior art will be briefly described below, it being obvious that the drawings in the following description are only some of the embodiments described in the present description, and that other drawings may be obtained according to these drawings without inventive effort for a person skilled in the art. In the drawings:

Fig. 1 is a schematic flow chart of a shale gas yield prediction method according to an embodiment of the present disclosure;

FIG. 2 is a schematic diagram of a simulated seepage model in one embodiment provided herein;

FIG. 3 is a schematic diagram of a crack cell distribution in a simulated seepage model in another embodiment provided herein;

FIG. 4 is a schematic illustration of flow exchange in a simulated seepage model in another embodiment provided herein;

FIG. 5 is a schematic diagram of coupling relationships between area units in another embodiment provided in the present specification;

FIG. 6 is a schematic diagram of a yield prediction flow in another embodiment provided herein;

FIG. 7 is a schematic illustration of a plate joint of a frac horizontal well model in another embodiment provided herein;

FIG. 8 is a complex, slotted network schematic of a frac horizontal well model in another embodiment provided herein;

FIG. 9 is a diagram showing comparison of yield predictions in another embodiment provided herein;

FIG. 10 is a schematic diagram showing yield prediction results in another embodiment provided in the present specification;

FIG. 11 is a schematic diagram illustrating encryption of a crack unit and a region unit according to another embodiment provided in the present disclosure;

Fig. 12 is a schematic block diagram of a shale gas yield prediction apparatus according to an embodiment of the present disclosure.

Detailed Description

In order that those skilled in the art will better understand the technical solutions in this specification, a clear and complete description of the technical solutions in one or more embodiments of this specification will be provided below with reference to the accompanying drawings in one or more embodiments of this specification, and it is apparent that the described embodiments are only some embodiments of the specification and not all embodiments. All other embodiments, which may be made by one or more embodiments of the disclosure without undue effort by one of ordinary skill in the art, are intended to be within the scope of the embodiments of the disclosure.

Fig. 1 illustrates a shale gas production prediction method provided in some embodiments of the present disclosure. The method may be applied to an apparatus, such as a server, that performs shale gas production predictions. As shown in fig. 1, the method may include the following steps.

S20: obtaining a simulated seepage model corresponding to a target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit includes a substrate.

The server may obtain a simulated seepage model corresponding to the target reservoir. The target reservoir may be any work area for which production can be predicted or reservoir characterization. The simulated seepage model may refer to a mathematical model for simulating the seepage state of a target reservoir so as to simulate and analyze the fluid flow characteristics, pressure distribution characteristics, producibility and the like of the target reservoir.

The simulated seepage model may include a plurality of zone units and one or more fracture units embedded in the zone units. The area unit may include a matrix. Typically, reservoirs in shale fracking sites are composed primarily of matrix, a large number of microcracks, and major fractures. The main cracks are cracks with larger crack dimensions, and the micro cracks are cracks with smaller crack dimensions. The matrix and fracture in the reservoir may be split to obtain a plurality of zone units and fracture units. In some embodiments, different flow characteristic parameter values can be respectively given to each area unit or each fracture unit to simulate the seepage distribution characteristics of the target reservoir, so that a simulated seepage model is constructed. The flow characteristic parameter may be, for example, a shape factor, porosity, permeability, etc.

In some embodiments, the zone unit may also include an equivalent continuous medium. Accordingly, the simulated seepage model can be constructed in the following manner. And carrying out continuous medium equivalent treatment on the flow characteristic parameters of the first cracks with the crack dimensions smaller than a preset threshold value in the target reservoir to obtain an equivalent continuous medium. And carrying out fusion treatment on the flow characteristic parameters of the equivalent continuous medium and the matrix of the target reservoir, and taking a medium model represented by the fused flow characteristic parameters as a reservoir background model. And carrying out discrete medium equivalent treatment on the flow characteristic parameters of the second cracks with the crack dimensions larger than or equal to a preset threshold value in the target reservoir to obtain discrete cracks. And carrying out subdivision treatment on the reservoir background model and the discrete cracks to obtain the simulated seepage model.

In some embodiments, for ease of description, a fracture with a fracture dimension less than a preset threshold may be described as a first fracture, i.e., a micro-fracture is described as a first fracture. A crack with a crack size greater than a preset threshold is described as a second crack, i.e., a main crack is described as a second crack. The fracture dimension may be a lateral dimension of the fracture. The magnitude of the preset threshold can be set according to the actual application scene.

As shown in fig. 2, fig. 2 shows a schematic diagram of a simulated seepage model of a hydraulic fracturing reformation work area. Fig. 2 (a) is a diagram showing microseismic monitoring results of a work area. Fig. 2 (b) shows a schematic diagram of a mesh model quantized based on microseismic monitoring results. The large-scale crack (second crack) can be extracted, and discrete medium treatment is adopted to obtain the discrete crack. Such as large scale fractures can be extracted from the target reservoir, and individually reduced to smooth parallel plates. The small-scale fracture (first fracture) may also be roughened into a continuous medium. Then, the equivalent flow characteristic parameters of the continuous medium and the flow characteristic parameters of the matrix are fused. For example, the flow characteristic parameter data representation model of the target reservoir at different positions can be constructed by combining the equivalent data of the flow characteristic parameters of the first fracture and the flow characteristic parameters of the matrix at different positions of the target reservoir and the seepage relation between the pores of the matrix and the first fracture, so as to form a dual medium model. The continuous medium equivalent treatment and the specific construction process of the dual medium model can be performed with reference to the prior art, and are not described herein. Thus obtaining a simulated seepage conceptual model of the reservoir formation in the graph (c) of fig. 2.

Correspondingly, the graph (c) in fig. 2 includes discrete fractures, heterogeneous reservoirs, and homogeneous reservoirs. The heterogeneous reservoir background model is a double medium consisting of small cracks equivalent to continuous medium and a matrix, and the homogeneous reservoir background model is a region formed by the matrix.

The heterogeneous reservoir, homogeneous reservoir formed by the matrix and the equivalent continuous medium may be collectively referred to as a reservoir background model. Accordingly, the target reservoir after the treatment is simplified into a reservoir background model and discrete fractures. Different reservoir and flow parameter values can be respectively assigned to different positions of a reservoir background model and a discrete fracture to simulate the seepage distribution characteristics of a target reservoir.

The server can divide a reservoir background model and discrete cracks of the target reservoir based on a preset grid step length to obtain a plurality of area units. The grid form of the area unit may be a rectangular grid, a triangular grid, an arbitrary quadrilateral grid, or the like. For example, the COMSOL can be utilized to perform the subdivision processing, so that the subdivision processing can be performed more accurately and efficiently. In some embodiments, it may be preferable to split the reservoir background model using a rectangular grid. The rectangular grid is used, the solving precision is not affected, the subdivision processing is simpler, and the assistance of COMSOL and other tools is not needed. Meanwhile, intersection points and crack dispersion are also simpler after the discrete cracks are embedded. And the rectangular grid system can be better combined with commercial geological modeling software, such as Petrel and the like.

The preset grid step length can be determined according to the size of the work area to be analyzed, the application scene and the like. The reservoir background model may be partitioned based on a preset mesh step size to obtain a plurality of region units. Geometric information of the respective region units may then be determined. Such as vertex coordinates, length and width of the grid cells, etc.

In some embodiments, the server may embed discrete cracks in the area units. Then, the discrete cracks can be subjected to unit division according to the intersecting conditions between the discrete cracks and the boundary lines of the discrete cracks and the area grid, so as to obtain a plurality of crack units.

For example, the server may embed each discrete crack into the divided area unit according to the position information of the discrete crack. Then, the server may divide the discrete cracks into units according to the intersections between the discrete cracks and the boundary lines between the discrete cracks and the area units, to obtain a plurality of crack units. As shown in fig. 3, fig. 3 (a) shows a schematic diagram of the plurality of fracture cells obtained. 1, 2, and 3 in the (a) diagram in fig. 3 represent the obtained fracture unit, respectively.

S22: calculating pressure values of a region unit and a fracture unit connected with a designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function.

And the server can calculate the pressure values of the area units and the fracture units connected with the designated shale gas well of the target reservoir by using a preset pressure calculation model. The preset pressure calculation model may include a calculation model that calculates pressure values of a center point of a boundary line of the area unit and a point on the fracture unit based on an unsteady green function.

As shown in fig. 4, in the simulated seepage model given in the above embodiment of the present specification, if the area unit is a dual medium, only flow exchange occurs between the matrix m and the corresponding first fracture f. The morphology and seepage characteristics of the first cracks in each area unit may have large differences, so that the distribution characteristics of the flow characteristic parameter data endowed by each area unit also have large differences. The flow characteristic characterization models of the different region units can be characterized in Laplace space with f(s). Where m is a Laplace space parameter corresponding to time, which is a constant when the time point is fixed.

The first cracks F in the area cells are in flow exchange with the corresponding matrix m, and also in flow exchange with the first cracks F of the adjacent area cells through the connected area cell boundary lines and in flow exchange with the crack cells F in the grid cells, so that the cracks are regarded as sources/sinks in the area cells. And the fracture unit F exchanges flow with other fracture units F intersecting the first fracture F in the embedded area unit, in addition to the first fracture F. But the crack cell does not exchange flow with the first crack in the otherwise unembossed area cell. Accordingly, a discrete fracture may be considered a source sink within a connected region element. The flow exchange between the matrix m and the first cracks f does not occur on the boundary line of the area unit, and only the flow exchange of the first cracks between the area units occurs on the boundary line of the area unit.

For well testing problems, the unsteady green's function base solution can be expressed in the form:

in the formula ,

r _D representing arbitrary points (x _D ′,y _D ') to the current solution point (x) _D ,y _D ) Is a distance of (2); k (K) ₀ Is a second class 0-order Bessel function; s is Laplace constant; f(s) is a parameter characterizing the type of medium.

Based on the unsteady green's function, a boundary integral equation of the form of the following equation (2) can be established in one area unit:

the above formula (2) is a boundary integral equation established by any point i in the area unit e, and the physical meaning of the equation is a representation form of the pressure value of the point i constructed after considering the influence of the boundary line of the area unit e, a crack unit and the seepage characteristics of the shale gas well on the pressure of the point i in the area unit e. The pressure and flow exchange at each point i in the area unit can be solved using the above equation. The first term represents the pressure value of the calculated point i, the second term and the third term represent the influence of the boundary line of the area unit e on the point i pressure, the fourth term represents the influence of the crack unit in the area unit e on the point i pressure, and the fifth term represents the influence of the shale gas well (if the shale gas well exists) in the area unit e on the point i pressure.

The above formula (2),

the green's basic solution is the pressure value at any point within the area element e. The complex seepage problem can be approximated by using the green basic solution, the calculation complexity is reduced, and meanwhile, for the well test problem, the calculation result can be more accurate by selecting the unsteady green basic solution. / >

The specific representation of (2) is shown in the formula (1).

Is the pressure at point i. θ _i For the angle value of the point i, if the point i is taken inside the grid, theta _i =2pi; if the point i is taken on the grid edge, θ _i =pi; if the i point is taken from the grid corner point, theta _i ＝π/2；/>

A pressure value for any point within the area unit; n represents the normal vector. />

Is a flow source sink item of a unit length on the fracture unit. />

Is a flow source sink item per unit length at the shaft of the shale gas well. k (k) _fD Is the permeability of the mesh. In the area unit e Ω ^e Is the whole area of the grid Γ ^be Is the region omega ^e Is defined by the boundary of the first pair of second pair of third pair. For orthogonal grid, there is +/per region unit>

Strip edge and gamma ^bej Represents the j-th boundary line of the area unit e. Γ -shaped structure ^Fe Is a crack unit contained in the area unit e, < >>

Representing the total number of crack elements embedded in the region element e Γ ^Fek Representing the k-th slit cell in the e-zone cell. All the parameters of the superscript band e represent the amount in the region element e.

Preferably, the solution point i may be set as a center point of a boundary line of the area unit or a point on the crack unit. The solution point is arranged on the central point of the boundary line of the area units, the two area units share one boundary line, and only three unknowns can be arranged on one solution point, including the pressure of the point and the flow of the two area units at the point.

As shown in fig. 5, the coupling conditions of the two area units are: the pressure of the two area units at the solving point is equal; conservation of mass. By setting an unknown pressure at a solution point, coupling conditions of equal pressure are readily achieved. The flow exchange capacity of two adjacent area units can be obtained by using the formula (3), and the mass conservation condition can be satisfied by the flow conservation at the same solution point, namely, the formula (4). The flow terms of the two area units on the common boundary line are equal in size and opposite in direction, and the sum of the flow terms is 0. Therefore, the unknown flow term on the edge of the area unit can be eliminated by using the formula (4), so that the solving point on the boundary line of the area unit has only one unknown pressure, the solving complexity can be greatly reduced, and the computing simplicity and efficiency are improved.

It is to be noted that,

the addition of the superscript e may indicate that the flow term at this solution point is related to the area element, i.e. the normal derivatives of the pressure of different area elements at the solution point may not be equal on the boundary line of the area element. Thus, when calculating the pressure of the solution point, it means that two unknown flow terms relating to the area unit have been considered at the same time, so that the accuracy of the calculation result can be ensured. After the pressure of the solving point is obtained through calculation, the two unknown flow items can be further determined, so that the flow of the crack intersected with the shale gas well and the contribution of the grid unit to the shale gas well can be determined, and the yield of the shale gas well can be predicted.

Based on formulas (3) and (4), formula (2) can be written in the following simplified form:

in the formula (5), the amino acid sequence of the compound,

/>

wherein ,

For the flow source and sink item per unit length of the kth crack unit embedded in the area unit e,/->

For the flow source sink at shale gas well, +.>

For the permeability of the zone unit e +.>

For the total number of boundary lines of the area unit e, +.>

Represents the total number of crack units embedded in the area unit e, < >>

The pressure coupling conditions between the area unit and the crack unit are as follows: when i is taken at the midpoint of the fracture cell, the pressure on the zone cell is equal to the pressure of the corresponding fracture cell. Accordingly, if i is at the midpoint of the fracture cell,

If i is in the vertical well position, +.>

wherein ,p_FD For the midpoint pressure value, p, at the fracture cell _wD Is the pressure value at the midpoint of the well bore of the vertical well in the zone unit.

Equations established based on all the zone units and the fracture units can then be coupled for problem solving. In some embodiments, the coupling matrix may be established using boundary lines of all zone units, fracture units, and shale gas wells. Correspondingly, the number of unknowns corresponding to the coupling matrix is N _b +2N _F +N _w And the equation is N _b +N _F +N _w And each. Therefore, an additional N is also required _F And equation, the solution of the analytical solution of the coupling matrix can be realized. In some embodiments, the solution matrix corresponding to each fracture unit may be determined according to a finite difference method of the discrete fracture percolation equation, thereby constructing another N _F And equations.

If the horizontal well is produced at fixed yield, the equation corresponding to the coupling area unit, the fracture unit and the shale gas well can be obtained _b +2N _F +N _w +1 equations, corresponding to a total unknown of N _b +2N _F +N _w +1, a unique solution of the model in Laplace space can be obtained by solving the equation.

If the horizontal well is produced at constant pressure, the pressure of the horizontal well is a constant value, and equations corresponding to the direct area unit, the crack unit and the shale gas well can be obtained _b +2N _F +N _w And total unknown quantity is N _b +2N _F +N _w And obtaining a unique solution of the model in Laplace space by solving an equation.

The solution of the Laplace space obtained can be converted into real space. For example, the solution of the Laplace space obtained can be converted to real space using the Stehfest method.

Due to the nonlinearity of the model, there is nonlinearity of both the matrix and the percolation model of the fracture system, including: PVT properties of the gas, desorption of the gas, nonlinear flow mechanisms of the gas in the matrix, stress sensitivity of the fracture system conductivity, etc. The flow shown in fig. 6 may be used to solve for the iterative manner of pressure and time dependent parameters. In FIG. 6

Representing the pressure corresponding to the initial iteration at the nth time step, p _i The initial pressure preset by the target reservoir is represented, ne represents the total number of boundary lines of each area unit in the simulated seepage model, and epsilon represents the allowable error limit. N (N) _t For time steps, corresponds to a constant s in Laplace space.

S24: and determining the flow rate of a regional unit or a fracture unit connected with the designated shale gas well of the target reservoir to the designated shale gas well according to the calculated pressure value.

S26: and predicting the yield of the designated shale gas well according to the flow rate.

After the pressure of the zone units and fracture units associated with the horizontal well bore is obtained, the flow rates of each zone unit and fracture unit to the designated shale gas well may be obtained, obtaining the component contributing to the production of the horizontal well bore. The horizontal well production may be the sum of all components.

The simulated seepage model and the solution mode established on the simulated seepage model based on the scheme provided by the embodiment can greatly simplify the solution complexity, greatly improve the solution efficiency, and further improve the reservoir seepage characteristics and the processing efficiency of shale gas well yield prediction.

The present specification also provides two examples to compare the accuracy of the results of the solution simulation provided by the above embodiments and to compare the computational efficiency of producing dynamic predictions with a commercial numerical simulator. FIG. 7 shows a schematic of a plate joint of a segmented multi-cluster fracturing horizontal well model. FIG. 8 shows a complex slotted network schematic of a staged multi-cluster fractured horizontal well model. Table 1 shows an example table of parameters used in the prediction.

TABLE 1

Parameters (parameters)	Value taking	Parameters (parameters)	Value taking
				Reservoir original pressure, MPa	30	Fluid compression coefficient, MPa ^-1	9×10 ^-3
Reservoir thickness, m	30	Viscosity of fluid, mPas	0.8
				Matrix permeability, mD	1×10 ^-4	Crack width, m	0.01
Porosity of the matrix	0.05	Crack porosity	0.2
				Rock compression coefficient, MPa ^-1	1×10 ^-3	Fracture permeability, mD	500
Bottom hole flow pressure, MPa	5

FIG. 9 is a schematic diagram comparing the predicted result of the embodiment scheme with the Eclipse predicted result of the commercial simulator corresponding to the flat seam. FIG. 10 is a schematic diagram comparing the simulation results based on the embodiment scheme with the EDFM simulation results corresponding to the complex stitch net.

For ease of description, the embodiment model of the present specification may be defined as the eGEM model. Eclipse and eGEM both employ three types of grids. In the EGEM model, the mesh steps of the three types of meshes are 105m, 52.5m and 26.25m, respectively, corresponding to 28, 112 and 448 area units and 90, 90 and 180 crack units, respectively. The Eclipse grids have maximum grid sizes of 30m, 6m and 2m respectively, and the corresponding grid numbers are 2233, 33775 and 124425 respectively.

Fig. 9 (a) shows Eclipse simulation results. As can be seen from the graph (a) in fig. 9, the calculation result error is large with 2233 meshes, and since the result change is small when the number of meshes is encrypted from 33775 to 124425, the result obtained with 124425 meshes can be used as a reference solution for the present example.

Fig. 9 (b) shows the result of the eGEM simulation. It can be seen that the error of the eGEM solution for 105m mesh step and 52.5m mesh step is small, with a maximum relative error of about 15% within 8000d of production. Since the calculation result obtained by adopting the 26.25m grid step length is basically consistent with the reference solution, eGEM can be realized by using a coarse grid (the grid size is 20-50 m) for the simulation result of Eclipse under 124425 grids.

The simulation calculation time of Eclipse and eGEM is shown in FIG. 7, the calculation time of eGEM is very similar to that of Eclipse of a commercial simulator, but the grid number used by eGEM is about 1/10-1/20 of that of Eclipse, and if parallel calculation is adopted, the calculation efficiency of eGEM provided by the application is greatly improved compared with that of Eclipse.

For the complex fracture shown in fig. 10, the co-simulation in this example gives 385 fractures, with 82 fractures in direct communication with the horizontal wellbore. Three types of grid sizes are used by eGEM and EDFM, respectively. The mesh sizes of the eGEM were 50m, 25m and 12.5m, respectively, and the EDFM used were 25m, 12.5m and 3.125m, respectively. The production prediction results are shown in fig. 8, in which the graph (a) in fig. 10 represents the simulation results of the EDFM; fig. 10 (b) shows simulation results of the eDEM. It can be seen that the accuracy of EDFM is more sensitive to grid size than eGEM. The present example uses the EDFM solution under a fine grid as a reference solution with a grid step size of 3.125m.

As can be seen from the graph (b) in fig. 10, when the mesh size is 50m or less, the relative error of 1700d and 9650d is less than 10%, and when the mesh size is 25m or less, the relative error is less than 5%. Using the same grid system, the computation time of the eGEM is longer than the EDFM. However, the simulation result shows that the grid sensitivity of the eGEM is low, and the coarse grid can be used to save the calculation amount. If a relative error of 5% is allowed, the present example may use a mesh size of 25m, requiring only 41 seconds in the simulation.

As can be seen from the comparison, the above embodiments provided in the present specification can use a wider mesh size to achieve more accurate calculation accuracy. With the same accuracy, the calculation efficiency of the above embodiment provided in the present specification is far higher than that of the existing commercial numerical simulation method. Under the condition of considering both the calculated amount and the calculated precision, the high-efficiency simulation of the complex stitch net is realized.

In other embodiments, the encryption processing may be further performed on the crack unit or the area unit where the crack unit and the crack unit are located based on a preset encryption method.

In some embodiments, only the crack unit may be encrypted, so as to obtain an encrypted crack unit. As shown in fig. 3 (b), the crack unit may be subjected to encryption processing. For example, the number of divisions corresponding to each slit unit may be determined according to the longitudinal extension length of each slit unit, and the respective slit units may be further divided equally according to the determined number of divisions, to obtain the encrypted slit units 1 to 6.

In other embodiments, the encryption processing may be performed on the area unit where the crack unit is located based on a preset encryption method, so as to obtain an encrypted area unit; the preset encryption method comprises the step of equally dividing the area unit where the crack unit is located or the step of equally dividing the boundary line of the area unit where the crack unit is located. And re-dividing the discrete cracks according to the intersecting conditions between the discrete cracks and the boundary lines of the encrypted area units and obtaining encrypted crack units.

In some embodiments, encryption processing may be performed on the area unit where the fracture unit is located, so as to implement encryption processing on the fracture unit. For example, the area units where the discrete cracks are located may be equally divided according to a predetermined division number, respectively, to obtain encrypted area units, as shown in fig. 11 (c). FIG. 11 shows a schematic diagram of a simulated seepage model. Fig. 11 (a) shows a schematic diagram of a model before the encryption processing. Then, the discrete cracks can be re-divided according to the intersecting condition between the discrete cracks and the boundary line between the discrete cracks and the encrypted area unit, so as to obtain a plurality of encrypted crack units.

In other embodiments, encryption processing may be performed on the boundary line of the area unit where the crack unit is located, and then, based on the encrypted boundary line, unit division may be performed again on the discrete crack according to the intersection condition between the discrete crack and between the discrete crack and the boundary line of the encrypted area unit, so as to obtain the encrypted crack unit. Alternatively, the boundary line of the area unit where the crack unit is located and the crack unit may be encrypted. For example, the boundary line of the area unit where the discrete fracture is located and the fracture unit may be divided equally according to a predetermined division number, respectively, to obtain an encrypted area unit and fracture unit, as shown in fig. 11 (d).

Because the scheme provided by the embodiment of the specification is based on the green function to process the flow between the area unit and the crack unit and between the area units, the grid effect of numerical simulation is weak, and the area unit with larger grid scale can be used. The flow calculation among the fracture units is based on finite difference, the calculation error is relatively larger, and the fracture units can be independently encrypted. In the same grid, encryption processing can be performed only on the crack unit under the condition that the regional unit parameters are kept unchanged. Encryption processing of the slit unit may also be achieved by encrypting the region unit containing the slit unit or encrypting the boundary line of the region unit under the same basic grid.

In view of the above embodiments of the present disclosure, the flow between the area unit and the crack unit and between the area units is processed by taking the center point of the boundary line of the area unit as the solving point, so that the same boundary integral equation can be directly used for solving whether only the crack unit is encrypted or the boundary line of the area unit is encrypted. The unknown parameters and the corresponding equation numbers related to the constructed coupling matrix directly correspond to the corresponding encrypted boundary line numbers and the corresponding crack unit numbers, and the solving cannot be influenced, so that other additional processing is not needed in the solving process. Compared with the scheme that the unencrypted grid is required to be subjected to virtual encryption in the Eclipse of the commercial numerical simulator, and the conductivity between the virtual grid and the encrypted grid is calculated, and the solution is already carried out, the scheme of the embodiment of the specification can further reduce the complexity and the processing amount of data processing and improve the processing efficiency while the encryption improves the calculation precision.

Of course, the embodiment of the present specification does not exclude the scheme of performing encryption processing on all the area units to simultaneously implement encryption on all the slit units. For example, each area unit may be equally divided according to a predetermined division number, and an encrypted area unit may be obtained as shown in fig. 11 (b). Then, the discrete cracks can be re-divided according to the intersecting condition between the discrete cracks and the boundary line between the discrete cracks and the encrypted area unit, so as to obtain a plurality of encrypted crack units.

In this specification, each embodiment is described in a progressive manner, and identical and similar parts of each embodiment are all referred to each other, and each embodiment mainly describes differences from other embodiments. Specific reference may be made to the foregoing description of related embodiments of the related process, which is not described herein in detail.

The foregoing describes specific embodiments of the present disclosure. Other embodiments are within the scope of the following claims. In some cases, the actions or steps recited in the claims can be performed in a different order than in the embodiments and still achieve desirable results. In addition, the processes depicted in the accompanying figures do not necessarily require the particular order shown, or sequential order, to achieve desirable results. In some embodiments, multitasking and parallel processing are also possible or may be advantageous.

Based on one or more of the embodiments described above, a simulated seepage model may be constructed that can accurately characterize heterogeneous flow characteristics of a reservoir by equating small-scale fractures to continuous media and large-scale fractures to discrete media. And then, setting a solving point of the simulated seepage model on a central point of the unit, and solving the pressure distribution in the simulated seepage model based on an unsteady-state green function method. By the aid of the solving mode, data processing complexity can be reduced, and data processing efficiency is greatly improved on the basis of guaranteeing solving accuracy. And then, determining flow distribution based on pressure distribution, and further accurately and efficiently predicting the yield of the shale gas well.

Based on the shale gas yield prediction method, one or more embodiments of the present disclosure further provide a shale gas yield prediction device. The apparatus may include a system, software (application), module, component, server, etc. using the methods described in the embodiments of the present specification in combination with necessary hardware implementation. Based on the same innovative concepts, the embodiments of the present description provide means in one or more embodiments as described in the following embodiments. Because the implementation scheme and the method for solving the problem by the device are similar, the implementation of the device in the embodiment of the present disclosure may refer to the implementation of the foregoing method, and the repetition is not repeated. As used below, the term "unit" or "module" may be a combination of software and/or hardware that implements the intended function. While the means described in the following embodiments are preferably implemented in software, implementation in hardware, or a combination of software and hardware, is also possible and contemplated. Specifically, fig. 12 is a schematic block diagram illustrating an embodiment of a shale gas yield prediction apparatus provided in the specification, and as shown in fig. 12, the apparatus may include:

The seepage model obtaining module 102 may be configured to obtain a simulated seepage model corresponding to the target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit comprises a matrix;

the calculating module 104 may be configured to calculate pressure values of a region unit and a fracture unit connected to a specified shale gas well of the target reservoir using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function;

the flow determination module 106 may be configured to determine, according to the calculated pressure value, a flow rate from a region unit and a fracture unit connected to a specified shale gas well of the target reservoir to the specified shale gas well;

the production prediction module 108 may be configured to predict production of the specified shale gas well based on the flow rate.

In other embodiments, the apparatus may further comprise a model building module. Wherein the model building module may include:

the first equivalent unit can be used for carrying out continuous medium equivalent processing on the flow characteristic parameters of the first cracks with the crack dimensions smaller than a preset threshold value in the target reservoir to obtain an equivalent continuous medium.

And the fusion unit can be used for carrying out fusion treatment on the flow characteristic parameters of the equivalent continuous medium and the matrix of the target reservoir, and taking a medium model represented by the fused flow characteristic parameters as a reservoir background model.

The second equivalent unit can be used for carrying out discrete medium equivalent processing on the flow characteristic parameters of the second cracks with the crack scale larger than or equal to a preset threshold value in the target reservoir to obtain discrete cracks.

The model construction unit can be used for carrying out subdivision treatment on the reservoir background model and the discrete cracks to obtain a simulated seepage model.

In other embodiments, the model building unit may include:

and the subdivision subunit can be used for subdividing the reservoir background model based on a preset grid step length to obtain a plurality of area units.

The crack dividing unit can be used for embedding the discrete cracks into the area unit, and dividing the discrete cracks into units according to the intersecting conditions between the discrete cracks and the boundary lines of the discrete cracks and the area unit to obtain a plurality of crack units.

In other embodiments, the model building module may further include:

The first encryption unit can be used for carrying out encryption processing on the crack unit to obtain an encrypted crack unit;

correspondingly, the calculation module can also be used for calculating the pressure values of the area units connected with the designated shale gas well of the target reservoir and the encrypted fracture units by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating the pressure values of the central point of the boundary line of the area unit and the point on the encrypted crack unit based on an unsteady-state green function.

In other embodiments, the model building module may further include:

the second encryption unit can be used for carrying out encryption processing on the area unit where the crack unit is located based on a preset encryption method to obtain an encrypted area unit; the preset encryption method comprises the steps of equally dividing an area unit where a crack unit is located or equally dividing a boundary line of the area unit where the crack unit is located;

the crack dividing unit can be further used for carrying out unit division on the discrete cracks again according to the intersecting conditions between the discrete cracks and the boundary lines of the discrete cracks and the encrypted area units to obtain encrypted crack units;

Correspondingly, the calculation module can also be used for calculating the pressure values of the area unit, the encrypted area unit and the encrypted fracture unit which are connected with the designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating the pressure values of the central points of the boundary lines of the area unit and the encrypted area unit and the points on the encrypted crack unit based on an unsteady state green function.

It should be noted that the above description of the apparatus according to the method embodiment may also include other implementations. Specific implementation may refer to descriptions of related method embodiments, which are not described herein in detail.

According to the shale gas yield prediction device provided by one or more embodiments of the specification, the small-scale cracks can be equivalent to continuous media, the large-scale cracks can be equivalent to discrete media, and a simulated seepage model capable of accurately representing heterogeneous flow characteristics of a reservoir can be constructed. And then, setting a solving point of the simulated seepage model on a central point of the unit, and solving the pressure distribution in the simulated seepage model based on an unsteady-state green function method. By the aid of the solving mode, data processing complexity can be reduced, and data processing efficiency is greatly improved on the basis of guaranteeing solving accuracy. And then, determining flow distribution based on pressure distribution, and further accurately and efficiently predicting the yield of the shale gas well.

The method or apparatus according to the above embodiments provided in the present specification may implement service logic by a computer program and be recorded on a storage medium, where the storage medium may be read and executed by a computer, to implement the effects of the schemes described in the embodiments of the present specification. Accordingly, the present specification also provides a shale gas production prediction apparatus comprising a processor and a memory storing processor executable instructions which when executed by the processor perform steps comprising the method of any one of the embodiments described above.

The storage medium may include physical means for storing information, typically by digitizing the information before storing it in an electronic, magnetic, or optical medium. The storage medium may include: means for storing information using electrical energy such as various memories, e.g., RAM, ROM, etc.; devices for storing information using magnetic energy such as hard disk, floppy disk, magnetic tape, magnetic core memory, bubble memory, and USB flash disk; devices for optically storing information, such as CDs or DVDs. Of course, there are other ways of readable storage medium, such as quantum memory, graphene memory, etc.

It should also be noted that the terms "comprises," "comprising," or any other variation thereof, are intended to cover a non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements does not include only those elements but may include other elements not expressly listed or inherent to such process, method, article, or apparatus. Without further limitation, an element defined by the phrase "comprising one … …" does not exclude the presence of other like elements in a process, method or apparatus comprising such elements.

In this specification, each embodiment is described in a progressive manner, and identical and similar parts of each embodiment are all referred to each other, and each embodiment mainly describes differences from other embodiments. In particular, for system embodiments, since they are substantially similar to method embodiments, the description is relatively simple, as relevant to see a section of the description of method embodiments. In the description of the present specification, a description referring to terms "one embodiment," "some embodiments," "examples," "specific examples," or "some examples," etc., means that a particular feature, structure, material, or characteristic described in connection with the embodiment or example is included in at least one embodiment or example of the present specification. In this specification, schematic representations of the above terms are not necessarily directed to the same embodiment or example. Furthermore, the particular features, structures, materials, or characteristics described may be combined in any suitable manner in any one or more embodiments or examples. Furthermore, the different embodiments or examples described in this specification and the features of the different embodiments or examples may be combined and combined by those skilled in the art without contradiction.

The foregoing is merely exemplary of the present disclosure and is not intended to limit the disclosure. Various modifications and alterations to this specification will become apparent to those skilled in the art. Any modifications, equivalent substitutions, improvements, or the like, which are within the spirit and principles of the present description, are intended to be included within the scope of the claims of the present description.

Claims

1. A shale gas production prediction method, comprising:

obtaining a simulated seepage model corresponding to a target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit comprises a matrix;

calculating pressure values of a region unit and a fracture unit connected with a designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function;

determining the flow rate of a regional unit and a fracture unit connected with a designated shale gas well of the target reservoir to the designated shale gas well according to the calculated pressure value;

predicting the yield of the designated shale gas well according to the flow rate;

The area unit also comprises an equivalent continuous medium, and correspondingly, the simulated seepage model is constructed by adopting the following mode:

performing continuous medium equivalent treatment on the flow characteristic parameters of the first cracks with the crack dimensions smaller than a preset threshold value in the target reservoir to obtain an equivalent continuous medium;

carrying out fusion treatment on the flow characteristic parameters of the equivalent continuous medium and the matrix of the target reservoir, and taking a medium model represented by the fused flow characteristic parameters as a reservoir background model;

performing discrete medium equivalent treatment on the flow characteristic parameters of the second cracks with the crack dimensions larger than or equal to a preset threshold value in the target reservoir to obtain discrete cracks;

performing subdivision treatment on the reservoir background model and the discrete cracks to obtain a simulated seepage model;

the splitting treatment of the reservoir background model and the discrete fracture comprises the following steps:

dividing the reservoir background model based on a preset grid step length to obtain a plurality of area units;

embedding the discrete cracks into the area units, and carrying out unit division on the discrete cracks according to the intersecting conditions between the discrete cracks and the boundary lines between the discrete cracks and the area units to obtain a plurality of crack units;

The preset pressure calculation model includes:

wherein ,

A flow source sink item per unit length of the j-th boundary line of the area unit e,

for the flow source sink at shale gas well, +.>

For the permeability of the zone unit e +.>

For the total number of boundary lines of the area unit e, +.>

Represents the total number of crack units embedded in the area unit e, < >>

2. The method according to claim 1, wherein the method further comprises:

encrypting the crack unit to obtain an encrypted crack unit;

correspondingly, calculating the pressure values of a region unit connected with the designated shale gas well of the target reservoir and the encrypted fracture unit by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating the pressure values of the central point of the boundary line of the area unit and the point on the encrypted crack unit based on the unsteady-state green function.

3. The method according to claim 1, wherein the method further comprises:

encrypting the area unit where the crack unit is located based on a preset encryption method to obtain an encrypted area unit; the preset encryption method comprises the steps of equally dividing an area unit where a crack unit is located or equally dividing a boundary line of the area unit where the crack unit is located;

carrying out unit division again on the discrete cracks according to the intersecting conditions between the discrete cracks and boundary lines of the discrete cracks and the encrypted area units to obtain encrypted crack units;

correspondingly, calculating the pressure values of a region unit, an encrypted region unit and an encrypted fracture unit which are connected with a designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating the pressure values of the central points of the boundary lines of the area unit and the encrypted area unit and the points on the encrypted crack unit based on the unsteady-state green function.

4. A shale gas production prediction apparatus, comprising:

the seepage model acquisition module is used for acquiring a simulated seepage model corresponding to the target reservoir; the simulated seepage model comprises a plurality of area units and a plurality of crack units; the area unit comprises a matrix;

The calculation module is used for calculating pressure values of a region unit and a fracture unit connected with the designated shale gas well of the target reservoir by using a preset pressure calculation model; the preset pressure calculation model comprises a calculation model for calculating pressure values of a central point of a boundary line of the area unit and a point on the crack unit based on an unsteady-state green function;

the flow determining module is used for determining the flow of a region unit and a crack unit connected with a designated shale gas well of the target reservoir to the designated shale gas well according to the calculated pressure value;

the yield prediction module is used for predicting the yield of the designated shale gas well according to the flow;

the apparatus further comprises a model building module, wherein the model building module comprises:

the first equivalent unit is used for carrying out continuous medium equivalent treatment on the flow characteristic parameters of the first cracks with the crack dimensions smaller than a preset threshold value in the target reservoir to obtain an equivalent continuous medium;

the fusion unit is used for carrying out fusion treatment on the equivalent continuous medium and the flow characteristic parameters of the matrix of the target reservoir, and taking a medium model characterized by the fused flow characteristic parameters as a reservoir background model;

The second equivalent unit is used for carrying out discrete medium equivalent treatment on the flow characteristic parameters of the second fracture with the fracture scale larger than or equal to a preset threshold value in the target reservoir to obtain a discrete fracture;

the model construction unit is used for carrying out subdivision treatment on the reservoir background model and the discrete cracks to obtain a simulated seepage model;

the model construction unit includes:

the subdivision subunit is used for subdividing the reservoir background model based on a preset grid step length to obtain a plurality of area units;

the crack dividing unit is used for embedding the discrete cracks into the area unit, and dividing the discrete cracks into units according to the intersecting condition between the discrete cracks and the boundary line between the discrete cracks and the area unit to obtain a plurality of crack units;

the preset pressure calculation model includes:

wherein ,

For the flow source sink at shale gas well, +.>

For the permeability of the zone unit e +.>

For the total number of boundary lines of the area unit e, +.>

Represents the total number of crack units embedded in the area unit e, < >>

5. Shale gas production prediction apparatus comprising a processor and a memory for storing processor executable instructions which when executed by the processor carry out the steps of the method of any of claims 1-3.