Transient transfer shape factor for fractured tight reservoirs: Effect of the dynamic threshold pressure gradient in unsteady flow

In fractured tight reservoirs, the seepage capacity of the matrix is poor, and fluid migration mainly depends on matrix/fracture transfer. An accurate understanding of the matrix/fracture flow is the basis of well tests and numerical simulations for tight reservoirs. In this paper, the unsteady flow equation for tight reservoirs is deduced based on the boundary layer theory, which can reflect the effect of the dynamic threshold pressure gradient, and the theoretical flow equation is verified by seepage experiments. Based on the study of the flow equation, the approximate semi‐analytical solution of the matrix/fracture unsteady transfer shape factor and the transfer function for tight reservoirs are established considering the early stage of matrix/fracture transfer (pressure does not propagate to the matrix center) and the late stage of matrix/fracture transfer (pseudo‐steady state). The results show that the shape factor of the tight reservoir is mainly affected by three factors (minimum threshold pressure gradient, static boundary layer thickness, and sensitivity coefficient of the fluid boundary layer), and the theoretical curves show that the intermediate transfer enters the pseudo‐steady state when the dimensionless time reaches approximately 0.14. The higher the minimum threshold pressure gradient is, the larger the transfer shape factor. The larger the static boundary layer thickness is, the larger the transfer shape factor; additionally, the larger the sensitivity coefficient of the fluid boundary layer is, the faster the change rate of the shape factor. Finally, the transfer shape factor is applied to a well test interpretation. Examples prove that the fitting accuracy of the new curve type is improved by 34.2% compared with the curve type for the conventional well test interpretation method.


| INTRODUCTION
The physical properties of fractured tight reservoir matrixes are poor, the seepage is difficult, and the seepage law is very complicated. Through experiments, scholars have confirmed that fluid seepage in tight reservoirs does not follow Darcy's law. [1][2][3] After new oil wells are put into operation, the fluid in the fracture first enters the wellbore, and the fracture pressure gradually decreases. Then, there is a pressure difference between the matrix and the fracture, and a fluid transfer flow occurs between the matrix and the fracture. 4 In a tight reservoir, fluid migration is mainly achieved by the transfer flow, and it can be seen in Figure 1.
Based on the theory of fluid bridge continuity, Dejam et al and Zhou et al [5][6][7][8][9][10] studied the capillary forces between matrix blocks and matrix blocks in fractured reservoirs, and a mathematical model of fractured reservoirs was established. Warren and Root proposed a dual-medium model to approximate the fractured reservoir, which is shown in Figure 1, and a corresponding transfer function was presented: According to Equation (1), transfer flow is mainly controlled by the shape factor (α). Therefore, the accurate solution of the shape factor is the basis for the quantitative calculation of the fluid exchange between the matrix and fractures; moreover, it is also the basis of the dual-medium mathematical model.
To date, scholars have carried out much research on the matrix/fracture transfer flow law by using theoretical methods and experimental methods. For the experimental studies of the transfer shape factor, Chang and Chaobal et al [11][12][13] used experimental methods to establish the empirical formula of the transfer shape factor considering the influence of time and determined that the transfer flow enters the quasi-stable period when the dimensionless time reaches approximately 0.1. A theoretical expression of the transfer shape factor is deduced by Civan et al. They simulated the matrix and fracture transfer flow with the condition of two-phase flow and radial flow, respectively, by taking laboratory experiments. 14 The fluid state in the fractures during different stages is studied visually from the microscopic point of view by Rangel-German. They simulated the physical model of the matrix/ fracture transfer flow using a square core gripper, and it is determined that the transfer flow enters the quasi-stable period when the dimensionless time reaches approximately 0.1. 15,16 The mass transfer law between the matrix and fracture for the flow of two phases of gas-liquid is studied by Rangel-German 13 based on a two-dimensional pore micro model. The transfer flow shape factor expressions are established by Rangel-German 17 using the theoretical methods and the experimental methods in his doctoral thesis, and the distribution states of fluid in different transfer flow stages are verified by the large number of CT structure scanning experiments.
For the theoretical studies of the transfer shape factor, scholars mainly adopt analytical methods and numerical simulation methods. A dual-porous medium model for describing the matrix/fracture transfer law was proposed, and a constant shape factor expression was given by Warren and Root. 18 The finite difference method was used by Kazemi 19 to solve the 3D shape factor expression. Combined with analyzed well test data, Coats et al solved a constant expression of the shape factor by establishing a mathematical model, and its value is twice that of Kazemi's shape factor. 20 On the basis of Coats' research results, Lim and Aziz 21 obtained an anisotropic 3D shape factor expression. The transfer shape factor for twophase gas fluids is given by Peaceman, 22 and it is demonstrated that the shape factor is independent of the number of fractures. The constant transfer flow shape factors are also solved by Thomas [23][24][25][26] The influences of different block sizes and shapes on the exchange flow between the matrix and fractures are studied by Zimmerman et al, 27 and asymptotic expressions are derived for the cumulative uptake as a function of the mean and the standard deviation of the radius distribution function for both normal and lognormal radius distributions. A general shape factor expression for irregular geometric blocks was established by Wuthicharn et al 28 based on the dual-pore medium model. A new linear transfer function for improving the convergence rate of the mathematical model of a matrix/fracture dual medium was proposed by Huang et al. 29 The influence of different constant shape factors on (1) q = k m (p m − p f ) F I G U R E 1 Schematic diagram of the fractured reservoir the permeability of the matrix/fracture system was studied by Stephan. 30 The shape factor of the matrix/fracture system is determined by Ranjbar et al 31 using a semi-analytical method, and the influence of block distributions with different sizes (uniform distribution, linear distribution, exponential distribution, and Gaussian distribution) on the transfer flow of the whole matrix/fracture system is studied. A transient transfer shape factor expression for the two-phase flow of oil-water in a fractured reservoir is deduced by Saboorian-Jooybari et al 32 considering the effects of gravity and capillary force, and the shape factor solution is simplified by the convergence theorem. A transient transfer shape factor expression of the nonlinear channeling model is solved by Ranjbar et al 33,34 combining the heat integral method, the method of moments, and Duhamel's theorem. Transient shape factor expressions for the conditions of constant pressure and constant flow boundary are solved by Hassanzadeh et al 35 using the Laplace transform. The transfer shape factor considering gravity and capillary force is deduced by Abbasi et al. 36 The transient shape factor expression is solved by Wang et al 37 using the Boltzmann transform function considering the stress-sensitive effect of the tight reservoir. Gong et al proposed that the shape factor is influenced by the network behaviors of fractures and pore size, and it should be different in the primary oil recovery stage and the secondary oil recovery stage. They modified the matrix/fracture transfer model based on the complex two-dimensional network behaviors of fractures. 38 Considering the influence of the threshold pressure gradient of a low permeability fractured reservoir, the matrix/ fracture transfer model was established by He et al, 39 and the analytical solution of the time-varying transfer shape factor was presented by the integral method. A quadratic flow equation is determined by Shu et al 40 by fitting the experimental data of differential pressure and flow, and then, an unsteady time-varying shape factor is solved by introducing the flow equation into the transfer model. The semi-analytical solutions of the concentration in the tight porous medium, the average concentration in the fracture, and the average diffusion flux at the interface are derived by Dejam et al 10 46,47 and it was determined that the time-dependent transfer shape factor eventually converges to values derived by Lim and Aziz. 21 The single-phase flow transfer shape factors under different discrete fracture conditions are solved by Sarda et al 48 using a numerical simulation method. Mora et al 49 confirmed that there is a large error between the calculated results of the commonly used transfer shape factor model and the numerical simulation results, and they proposed a modified formula for calculating the transfer shape factor. By solving the capillary diffusion equation under different boundary conditions, Saboorian-Jooybari et al obtained the analytical solutions of the fluid saturation distribution in the matrix block and verified their solutions in combination with fine finite element analysis and experimental data. They proposed that the shape factor of the quasi-stationary period is stable to a constant value; however, the nonstationary period is completely phase-sensitive, and they studied the dependence of the shape factor on boundary conditions. 50 After giving the pressure distribution formulas for constant fracture pressure and fracture pressure varying with time, Ranjbar et al 51 studied the influences of different bedrock shapes, sizes, and decision regimes on shape factors. Based on the streamline model, Donato et al established a reactive transport model for two phases of oil-water in a dual medium. They solved the matrix/fracture transfer function using the difference method and compared the results with Kazemi's model. 52 The previous studies on constant shape factors and the corresponding solutions as well as the applicable reservoir types are summarized in Table 1. It can be seen from the summarized results of Table 1 that different scholars have studied the transfer flow law and the impact factors on dual-medium matrix/fracture from different perspectives at present, but most of the studies are based on linear flow equations, which causes the obtained transfer shape factors to only be suitable for medium and high permeability reservoirs and cannot be directly applied to tight reservoirs. Through investigation, only the shape factors pre- (influenced by the exponential stress-sensitive model). Shu et al determined the quadratic flow equation by fitting the experimental flow data. After being introduced into the transfer model, the transient transfer shape factor for the unsteady state is solved. 40 However, the quadratic flow equation fits the low-speed non-Darcy flow well and fits the high-speed Darcy flow poorly, and there is no theoretical basis for fitting the experimental data directly, which is not conducive to commercialization and application.
This study aims to determine a transfer function and a type of shape factor suitable for tight reservoirs and to establish a foundation for future well tests and numerical simulations of tight reservoirs. The sequence of our research is as follows: First, the derivation of the unsteady flow equation suitable for tight reservoirs is presented based on boundary layer theory. Then, the matrix pressure distributions during the process of transfer are obtained by using an analytical method after introducing the unsteady flow equation into the governing equation, and a new transfer function and shape factor are determined. Finally, the results are applied to well tests for tight reservoirs, and these results are beneficial to improve the fitting accuracy of the sample curve. The novelties of this study are as follows: The dynamic threshold pressure gradient is considered when establishing the flow equation, and the determinations of pressure distributions, transfer function, and shape factor are all derived based on the new nonlinear flow equation. The impact factors of the dynamic threshold pressure gradient (minimum threshold pressure gradient, boundary layer thickness, and sensitivity coefficient of the fluid boundary layer) on the transfer shape factor are analyzed. The new transfer flow function and shape factor are applied to the well test for tight reservoirs, and unsteady curve types of the matrix/fracture transfer stage are established.

FLOW EQUATION FOR TIGHT RESERVOIRS
The unsteady flow equation for tight reservoirs is established based on the boundary layer theory, 55,56 and the accuracy of the flow equation is verified by conducting the core seepage experiment.

| Derivation of flow equation
In the porous medium of a tight reservoir, the pore radius and fluid boundary thickness are on the same order of magnitude. The interaction between fluids and porous medium is strong, and the macroscopic performance is the unsteady flow characteristics. 57 According to the boundary layer theory, for the same capillary, the displacement pressure gradient is inversely proportional to the thickness of the boundary layer. The larger the differential pressure is, the thinner the viscous layer on the capillary wall, and there is an exponential decrease in the relation compliance. 58,59 The relation between the thickness of the fluid boundary layer and the displacement pressure gradient is as follows: There are many fine throats in tight reservoirs, and fluids passing through throats behave as non-Newtonian Bingham fluids because of the influence of the boundary layer. That is, the fluids begin to flow when the displacement pressure gradient is greater than the yield stress τ 0 . The constitutive equation of Bingham is as follows: where du dr is the velocity gradient, that is, the shear rate, s −1 . This equation is usually used in describing the flow characteristics of fluids in porous media in tight reservoirs. 60,61 According to capillary theory, the radial flow velocity of a Bingham fluid in capillaries is: Considering the influence of the boundary layer, the fluid in the viscous layer does not flow, and then, the equivalent flow radius is: The distribution equation of flow velocity can be obtained by integrating the whole flow interval: The flow rate in a single capillary can be obtained by integrating the flow velocity on the cross-section.
That is, According to Equation (9), the fluid passes a capillary of radius r 0 due to the action of the displacement pressure gradient ∇p; its flow rate consists of two parts: The first part is the flow caused by the displacement pressure gradient, and the second part is the flow loss caused by the yield stress of the boundary layer. According to the definition of the threshold pressure gradient, the pressure gradient when the flow rate in the capillary is greater than the critical value zero is: The expression of the dynamic threshold pressure gradient can be obtained by substituting Equation (2) into Equation (10).
According to Equation (11), the dynamic threshold pressure gradient is related to the static parameters (yield stress τ 0 , throat radius r 0 , and static boundary layer thickness δ 0 ) and dynamic parameters (displacement pressure gradient ▽p). Substituting Equation (2) into Equation (9) and rearranging gives: Based on the capillary model, the Poiseuille equation and the Darcy equation are simultaneously: The expression of permeability is obtained.
A new unsteady flow equation considering the influence of the throat boundary layer for tight reservoirs can be obtained by substituting Equation (11) and Equation (15) into the capillary flow equation Equation (12).
According to Equation (16), the throat boundary layer thickness δ 0 and fluid boundary layer sensitivity coefficient c in Equation (16) reflect the change in the dynamic threshold pressure gradient with different pressures. The parameters can be fitted by flow experiments of tight reservoir cores.

| Simplification of the flow equation
Determination of the flow equation requires that the theoretical calculation results coincide with the "differential pressure-flow velocity" curve of the flow experiment, that is, the accuracy of the equation; on the other hand, the determination of the flow equation requires that the equation can be solved analytically and applied easily after substituting it into the matrix/fracture transfer model.
The flow in the tight reservoir matrix includes both porous flow and throat flow. The sponsoring agency of this project provided some interpretation results of rate-controlled mercury penetration experiments of pore throats in a tight reservoir block, as shown in Table 2.
Through analyzing the interpretation results of pore throats in tight reservoirs, it can be seen that the average pore radius is 7.529 μm, while the main throat radius is only 0.8948 μm, which is approximately one order of magnitude different from each other. Therefore, the flow law of fluid in pores is quite different than that in throats.
The flow equation Equation (16) we established is based on the boundary layer theory. In other words, the flow equation is an idealized description of fluid flow in capillaries. According to the sketch map (Figure 2), the flow equation can relatively accurately describe fluid flow in the throat (the capillary radius is very small, and the thickness of the boundary layer accounts for a high proportion); however, when the flow equation is used to describe the flow of fluid in a porous medium (pores and throats), the initial curvature of the flow curve with an exponent of 4 is too large, especially when the pressure gradient is very small. Therefore, we proposed a mathematical method to simplify the flow equation and expand it by the Taylor formula, and then, the simplified equation becomes the following form: The theoretical curves of the flow equation before and after simplification are proposed, as shown in Figure 2. After simplification, the curvature of the flow curve decreases when the pressure gradient is small. In terms of the phenomenon, this mathematical treatment has achieved a certain effect.
Moreover, the purpose of simplification is to ensure the solvability of the subsequent pressure governing equation. The initial form of the flow equation Equation (16) } includes both an exponent of 4 and the exponential term. Without simplification, the highest exponent of the basic differential equation is 7. For analytical solutions, such an equation cannot be solved; for numerical solutions, the difference in equations after discretization are very complex, there are also high-order problems, and the computer solving speed is very slow. Therefore, we adhere to the principle of simplifying reasonably, and the simplified flow equation is proposed above. Next, the simplification is proven to be feasible through fitting flow experiment data.

| Experimental verification of flow equation
After theoretically explaining the rationality and necessity of the simplification, we further carried out flow experiments by using the experimental method of the "capillary balance method" and "differential pressure flow method" to record the flow velocity under different differential pressures, and both the original flow equation and the simplified flow equation are used to fit the experimental data.
(1) The experimental device and principle Experimental equipment: high-precision double cylinder constant speed constant pressure displacement pump produced by ISCO company (Type A260), rotary-vane vacuum pump (Type 2XZ-8B), thermostat (Type ZJ-HK), confining pressure controller by Glocom company of the USA (Type Cpc-110), DXD high-precision digital pressure sensor, microflowmeter (Type YRD-WL-003), triaxial core holder, and simulated formation oil (mixture of kerosene and white oil with viscosity of 1.5 mPa s).
The principle of testing: the steady state method.
(2) The experimental procedures Step 1: The core is dried for 48 hours, and the length, diameter, and gas permeability of the cores are measured.
Step 2: The core is vacuumed for 12 hours using the simulated formation oil (kerosene + white oil) for saturation, and then, the effective porosities of the cores are calculated by weight.
Step 3: The core is put into the core holder, the instruments are connected according to the flowchart (Figure 3), the zero point of the instruments is calibrated, and then, the confining pressure controller is used to gradually increase the confining pressure value to the effective stress range of the reservoir (25 MPa) in steps of 2 MPa and intervals of 10 min.
Step 4: The liquid in the outlet pipe is drained, and the pipe is placed under the liquid level of the collection beaker.
Step 5: The displacement pump is turned on, simulated formation oil is injected into the core at a flow rate of 0.01 mL/min, and the outlet was observed to determine whether there is bubble emerging from it.
Step 6: When the first bubble comes out of the outlet pipe, the inlet valve is closed at the front end of the gripper, the pump is stopped, and it remained standing until the data collected by the differential pressure sensor are stable for more than 48 hours; then, the minimum threshold pressure is displayed on the differential pressure sensor.
Step 7: The inlet valve is opened, the effective stress is set as 25 MPa, the displacement pump is started, the simulated formation oil is injected at a constant pressure, and the flow velocity is recorded on the microflowmeter until the flow is stable.
Step 8: The displacement pressure is increased, Step 7 is repeated until all experimental design pressure points are determined, the experiment is finished, and the instruments are arranged.  The "pressure difference-flow velocity" data of three cores with permeability values less than 1 × 10 −3 μm 2 are tested, and the basic parameters of the experimental cores are shown in Table 3. The fitting results of the theoretical equations are shown in Figure 4.
The establishment of the flow equation is mainly used to describe the seepage characteristics of fluid in tight reservoirs. According to Figure 4, the two flow equations before and after simplification can both fit the experimental data with high accuracy, and the simplified seepage equation has a higher accuracy when describing the unsteady flow period with a low pressure gradient. Therefore, this simplification is both reasonable and necessary.

| MATRIX/FRACTURE TRANSFER MODEL AND ITS SOLUTION
The unsteady pressure diffusion equation for tight reservoirs is deduced by simplifying the matrix/fracture dual-porous medium model. The pressure distribution in the early and late stages of transfer is solved at the constant pressure boundary condition. The semi-analytical solution of the transfer shape factor varying with time and the transfer function are given.

| Establishment of the pressure diffusion equation
The transfer between the matrix and fracture system in tight reservoirs is mainly controlled by the pressure diffusion equation of fluid in matrix blocks. The fracture system is assumed to be the pressure boundary to the matrix rock block. Considering the symmetry of the matrix rock block, the center of the matrix block is taken as the origin. The simplified model is shown in Figure 5.
The continuity equation describing fluid flow in the matrix is: The unsteady pressure diffusion equation in the matrix is obtained by substituting the unsteady flow equation Equation (17) derived from the boundary layer theory into the continuity equation Equation (18). Its 1D form is as follows:

| Solution of the diffusion equation
Pressure diffusion in a matrix/fracture system can be divided into the following two stages.
a. During the early stage of transfer, the pressure does not propagate to the matrix center, and the matrix is mainly divided into two parts: the pressure excitation zone affected by the fracture system and the nonexcitation zone at the initial state. b. During the late stage of transfer (quasi-steady state), the pressure has spread to the matrix center, and the whole matrix system is in the excited state.
A dimensionless parameter related to time is introduced in this paper for describing the pressure diffusion characteristics of the matrix/fracture system shown in Equation (20). The physical meaning of this parameter is the location of the boundary of the excitation zone, and T * is the time when the pressure propagates to the matrix center.
After a qualitative analysis of the pressure diffusion law, the diffusion equation can be solved. Considering that the partial differential equation in the form of Equation (18) cannot be solved by the conventional method, in this paper,  62 The pressure governing equations, boundary conditions, and pressure analytical solutions are given, while the specific solution process can be seen in Appendix.
(1) Dimensionless treatment Considering the symmetry of the flow in each direction of the matrix block, this paper only studies the matrix/fracture system in the left half of the simplified model ( Figure 5). Dimensionless treatments of the equation are as follows: After the dimensionless treatments, the unsteady pressure diffusion equation and initial/boundary conditions are obtained. During the early stage of transfer, the boundary of the excitation zone θ (T) is less than 1, and the nonexcitation zone is still in the initial state, so the initial and boundary conditions become: Former Soviet scholar Bachiev proposed an integral method to solve the problem, and he believed that the pressure distribution in the excitation zone affected by the pressure wave could be approximated by an exponential polynomial of X.
Coefficients a 0 , a 1 , and a 2 are determined by substituting the pressure boundary conditions Equation (25) into Equation (26), and then, the pressure distribution equation related to θ (T) is obtained: The boundary of the excitation zone is 0 ≤ (T) < 1; therefore, the dimensionless time range is: The pressure distribution in the excitation zone is obtained by substituting Equation (27) into Equation (26).
To obtain the mean pressure of the matrix system, the integral of pressure is carried out over the whole matrix, including the excitation and nonexcitation zones. The average pressure of the matrix system during the early stage of transfer is: (3) Pressure distribution during the late stage of transfer The whole matrix system enters the excitation state after the pressure propagates to the matrix center. According to Equation (26), the pressure distribution in the matrix becomes as follows: Combined with the results of the pressure distribution during the early stage of unsteady transfer, the initial and boundary conditions become as follows: The equations of parameters (a 3 , a 4 , and a 5 ) can be obtained by substituting the initial conditions Equation (32a) and boundary conditions Equation (32b) and Equation (32c) into the pressure distribution formula Equation (31). Therefore, the pressure distribution solution during the late unsteady stage is obtained.
where The average pressure of the matrix system in the late stage of transfer is:

SHAPE FACTOR AND TRANSFER FLOW FUNCTION
After solving the average pressure at different states in the matrix system, the expression of the shape factor is determined by the mass conservation equation. After substituting the shape factor into the transfer function, the new expression of the transfer flow function is determined.

| Determination of the shape factor
According to the mass conservation equation, the intermediate transfer flow function can be expressed as the mass change per unit time.
Combining Warren and Root's transfer function and the unsteady flow equation Equation (16), the transfer flow function can also be expressed as: The expression of the transient transfer shape factor is obtained by simulating Equation (36) and Equation (37).
In the formula, the parameter group 1 − e −c|∇p−G| 4 1 − G ∇p is the influence of unsteady flow on transfer. Lim et al 21 proposed Equation (11) for extending the transfer shape factor directly from 1D to 3D. The extension of the shape factor to the 3D case is: That is where

| Determination of the transfer flow function
The new transfer function is determined by substituting the derived expression of the transfer shape factor Equation (40) into the transfer function Equation (1).
According to Equation (30) and Equation (35), there are: The unsteady matrix/fracture transfer function for tight reservoirs can be obtained by substituting Equation (42a) and Equation (42b) into Equation (41).
where According to Equation (43), the transfer flow between the matrix and fracture is related to the pore structure parameters (δ 0 , c), reservoir physical properties parameters (k m , μ, ρ, c m , φ, and G), the initial pressure in the matrix system, and the pressure in the fracture system.

EFFECT OF DYNAMIC THRESHOLD PRESSURE GRADIENT ON SHAPE FACTOR
Taking the actual data of the tight reservoir in the S oilfield as an example, the exact numerical solution of the finite element method is applied to verify the calculation model of the shape factor, and then, the influence of the dynamic threshold pressure gradient on the shape factor is analyzed by plotting the chart.
LIU et aL.

| Model validation
The approximate semi-analytic solution of the shape factor and transfer flow function has been obtained in the previous article by establishing the matrix/fracture transfer model for tight reservoirs. Presently, the finite element method is usually used to verify the solution of this kind of model at home and abroad. 37,40,54 Finite element analysis simplifies the complex matrix/fracture system into multiple subsystems, which includes the governing differential Equation (19) and the volume transport flux Equation (43). The exact solution of the accumulative flux of fine grids can be obtained by a numerical method, which verifies the accuracy of the approximate analytical solution of the shape factor. The simplified finite element model and size of the matrix/fracture system are shown in Figure 6. The physical and fluid parameters of the tight reservoir in the S oilfield are shown in Table 4. The comparison between the calculated results of the transfer shape factor derived in this paper and those of Warren and Root's model is shown in Figure 7. The comparisons of cumulative flux among the different transfer flow functions (finite element model, the new transfer function established in this paper, and Warren and Root's transfer function) are shown in Figure 8.
The comparison between the different shape factor models shows that the shape factors during the early and late stages of transfer are in good continuity because the pressure distributions during the different stages of transfer use the same polynomial distributions Equation (26) and Equation (31) approximation instead of the real solution when solving the nonlinear differential equation (19). Matrix/fracture transfer enter the quasi-steady state when the dimensionless time reaches approximately 0.14; the shape factor during the quasi-steady state is stable at approximately 18.6. The theoretical critical point of dimensionless time obtained in this paper is slightly larger than the experimental results obtained by Chang et al and the theoretical calculation results given by Shu et al. This outcome occurs because the unsteady flow in tight reservoirs is seriously affected by the boundary layer and throat, and the pressure propagates more slowly in the matrix; that is, the time it takes to reach the quasi-steady state is longer. Compared with Warren and Root's constant shape factor, the calculated results of the two methods are more than one order of magnitude different during the early stage of transfer, but there is little difference during the late stage of transfer. The

| Effect of the dynamic threshold pressure gradient on the transfer shape factor
To explain this problem more succinctly and clearly, we have presented the conventional research sequence on the shape factor for tight reservoirs and the sequence we proposed. According to Figure 9, the flow equation affects the shape factor by affecting the pressure distribution in the matrix during the process of transfer. In previous studies, scholars only used one single parameter (minimum threshold pressure gradient) to modify the flow equation for tight reservoirs. In our study, the unsteady flow equation is affected by three factors (minimum threshold pressure gradient, thickness of fluid boundary layer, and sensitivity coefficient of fluid boundary layer), and these three parameters can be determined simultaneously by seepage experiments. For the minimum threshold pressure gradient, the equilibrium method can be used to measure it directly. For the thickness of the fluid boundary layer and sensitivity coefficient of the fluid boundary layer, the simplest way to obtain it is to fit the flow curve. According to the boundary layer theory, the fitted parameters also have clear physical significance, and it is more convenient for engineering applications.
In this paper, the effects of different minimum threshold pressure gradients (0.2, 0.4, and 0.6 MPa m −1 ) and static boundary layer thickness δ 0 (0.55, 0.50, and 0.45) and sensitivity coefficient of fluid boundary layer c (0.8, 0.5, and 0.2 MPa −1 ) on the transfer shape factor are analyzed in Figure 10.
According to the theoretical curves under the different conditions, the following conclusions can be obtained.
1. The higher the minimum threshold pressure gradient is, the larger the initial transfer shape factor. This outcome occurs because the higher the minimum threshold pressure gradient is, the greater the fluid yield stress, and the flow in the pore throat is more difficult. According to the simultaneous equations Equation (36) and Equation (37), the smaller the G min is, the larger the shape factor under a certain transfer flow. 2. The larger the boundary layer thickness is, the larger the transfer shape factor. This outcome occurs because the lower the permeability is, the smaller the average pore throat radius. When the fluid passes through, the higher the proportion of the boundary layer in the flow channel is, the lower the proportion of the flowable region, and then, the instantaneous flow in the pore throat is lower, which requires a larger shape factor. 3. The larger the sensitivity coefficient of the fluid boundary layer is, the faster the change in the transfer shape factor. This outcome occurs because the larger the sensitivity coefficient of the fluid boundary layer is, the more sensitive the thickness of the boundary layer to the change of the pressure difference, and the greater the change in the instantaneous flow rate in the pore throat, which leads to a faster reduction rate in the shape factor.

| APPLICATION OF THE MODEL IN WELL TEST INTERPRETATIONS
Applications of the new unsteady transfer shape factor model can be applied to the numerical simulation of tight reservoirs to improve the simulation accuracy; on the other hand, it can be applied to the well test interpretation of tight reservoirs to improve the fitting accuracy of typical well test curves. In this paper, the new factor model is applied to a well test interpretation, and the new theory of the well test interpretation is validated by taking the actual pressure recovery data of a single well in the S oilfield as an example. (

1) New Model of the Well Test Interpretation Theory
After opening a well, the flow after fractured tight reservoirs can be divided into the following three stages. First, the stage of fluid in the fracture flowing to the wellbore is when the fracture pressure decreases continuously. Second, the stage of intermediate unsteady flow is when there is a pressure difference between the matrix and fracture and matrix transfer to fracture. Third, both the fluid flow from the matrix to the fracture and the fluid flow from the fracture to the wellbore exist when the matrix pressure p m decreases to the fracture pressure p f . At this stage, the decrease in p m and p f is the same.
Current well test theory [63][64][65] usually interprets this model by using two sets of curve types, one of which is a homogeneous reservoir curve type. Each curve corresponds to a C D e 2S , which is used to fit the first stage {CDe 2S =(C D e 2S ) f } and the third stage {CDe 2S =(C D e 2S ) f+m }. The other set is the unsteady flow curve type of the matrix/fracture system. Each curve corresponds to a β′, and the β′ curve type is used to fit the second stage. The value of ε is 1.8914 for the Warren and Root dual-porosity medium model, and the value of λ is the transfer coefficient.
A new well test interpretation model for intermediate unsteady transfer in fractured tight reservoirs is obtained by substituting the transfer shape factor established in this paper Equation (40) into Equation (45).
(2) Fitting results of typical well test curves The pressure build-up curve of the 10R108-68 single well in the S block of the domestic tight reservoir is taken as an example, and the actual pressure build-up curve is fitted with the classic curve type Equation (44) Figure 11, and the basic data for the well test interpretation are shown in Table 5.
The fitting results of the curve types for the unsteady intermediate transfer stage show that the fitting error of the new curve type is 8.6%, which is 34.2% lower than that of the classic curve type. This outcome occurs because the new well test interpretation model Equation (46) established in this paper is based on the unsteady flow equation of tight reservoirs and takes the effect of the dynamic threshold pressure gradient into account, which is more in line with the real situation.

| CONCLUSIONS
1. Based on the boundary layer theory, the unsteady flow equation for tight reservoirs is deduced, which can reflect the effect of the dynamic threshold pressure gradient. It is clear that the dynamic threshold pressure gradient is mainly affected by three factors (minimum threshold pressure gradient, thickness of boundary layer, and the sensitivity coefficient of the fluid boundary layer). Moreover, the theoretical "differential pressure-seepage velocity" curve is in good agreement with the experimental data. 2. Considering the early stage of transfer (pressure does not propagate to the center of the matrix) and the late stage of transfer (quasi-steady state), the approximate analytical solution of the transient shape factor is obtained by combining the integral method with the moment method. It is clear that the intermediate transfer enters the quasisteady state when the dimensionless time reaches 0.14. Moreover, the matrix/fracture transfer flow function is determined, and the calculated cumulative volume flux is in good agreement with the exact numerical solution of the finite element model. 3. The influence of the dynamic threshold pressure gradient sensitivity parameters on the shape factor is studied. That is, the larger the minimum threshold pressure gradient is, the larger the shape factor, the larger the thickness of static boundary layer, and the larger the shape factor; in addition, the larger the sensitivity coefficient of fluid boundary layer, the faster the change rate in the transfer shape factor with time. 4. The transient transfer shape factor is applied to well testing, and a new curve type for intermediate unsteady transfer in fractured reservoirs is established. Taking the pressure buildup data of a single well in the S block of a domestic tight reservoir as an example, the fitting accuracy of the new curve type is improved by 34.2% compared with the conventional well test curve type. The matrix permeability, 10 −3 μm 2 μ The fluid viscosity, mPa s α The transfer shape factor, m −2 P m The average matrix pressure, MPa p f The fracture pressure, MPa Δ The thickness of fluid boundary layer, m r 0 The capillary radius, m δ 0 The ratio of the thickness of the fluid boundary layer to the radius of the capillary when the differential pressure is 0, also called static boundary layer thickness c The fluid boundary layer sensitivity coefficient, The displacement pressure gradient, MPa m −1 τ The shear stress on fluid, MPa τ 0 The yield stress, MPa du dr The velocity gradient, that is, the shear rate, s −1 r The equivalent flow radius of fluid in capillary, m q r The flow rate in a single capillary, m 3 s −1 G The threshold pressure gradient, MPa m −1 △L The single capillary length, m v The flow velocity, m s −1 φ The porosity v x , v y , v z The flow velocities in x, y, and z directions, respectively, m s −1 t The transfer time, s T, P, X The dimensionless time, pressure, and distance to the boundary of fracture, respectively T * The time of pressure propagation to matrix center, s θ (T) The location of the boundary of the excitation zone, dimensionless p mi The matrix initial pressure, MPa L x , L y , L z The dimensions of matrix blocks in all directions, m C t The total compressibility of reservoir, MPa −1 α x , α y , α z The transfer shape factors in x, y, and z directions, respectively, m −2 ξ 1 The parameters group, ξ 1 = 2δ 0 c ξ 2 the parameters group, ξ 2 = (1 − δ 0 ) − 2δ 0 cG m The parameters group, m = 2ξ 1 (p f − p mi )/L x n The parameter, n = ξ 2 ω The dimensionless parameter R The residual error C 0 The unknown parameters of equation after indefinite integration C D The dimensionless wellbore reservoir coefficient, The skin factor ε The shape parameters of dual-porous medium model, the ε is 1.8914 when the matrix block is square λ The flow factors in porous medium, m 2 Considering the symmetry of the flow in each direction of the matrix block, this paper only studies the matrix/fracture system in the left half of the simplified model ( Figure 2). The initial and boundary conditions at constant pressure conditions are as follows.

Dimensionless treatment of initial and boundary conditions is
That is, After finishing, the dimensionless unsteady pressure diffusion equation is obtained.
(2) Pressure distribution during the early stage of transfer During the early stage of transfer, the boundary of the excitation zone θ (T) is less than 1, and the nonexcitation zone is still in the initial state; thus, the initial and boundary conditions become Former Soviet scholar Bachiev proposed an integral method to solve the problem, and he believed that the pressure distribution in the excitation zone affected by the pressure wave could be approximated by an exponential polynomial of X.
Coefficients a 0 , a 1 , and a 2 are determined by substituting the pressure boundary conditions Equation (A7) into Equation (A8), and then, the pressure distribution equation related to θ (T) is obtained.

Thus,
To solve θ (T) , the integral is carried out in the whole excitation region, and the substitution of Equation (A10) into Equation (A4) occurs.
The first-order nonhomogeneous nonlinear ordinary differential equation is obtained.
The mathematical software Maple is used to solve the above equation. The analytical solution includes the Lambert W function. There are two considerations for utilizing this function: One consideration is that there are multiple solutions to solve the Lambert W function directly, which is not considered to be single mapping, and there are still multiple solutions to the function on the premise that the independent variable is a real number; the other consideration is that the convergence radius (e −1 ) needs to be considered when using the asymptotic expansion of the Taylor expansion, and the function still has two branches when using the other expansion with a very complicated form. Overall, the Lambert W function is inconvenient to use. Therefore, Equation (A12) is appropriately simplified by introducing parameter ω in this paper.
The value of parameter ω directly affects the smoothness of the inflection point at the border of the excitation zone and nonexcitation zone. We determine the parameter ω mainly from the following two aspects: first, to ensure the smoothness of the border; second, to ensure the first moment of the residual is 0 after this simplification. Then, the analytical solution of θ (T) is obtained by substituting the initial condition Equation (A7) into Equation (A13).
The pressure distribution in the excitation zone is obtained by substituting Equation (A14) into Equation (A9).
To obtain the mean pressure of the matrix system, the integral of pressure is carried out over the whole matrix, including the excitation and nonexcitation zones.
The average pressure of the matrix system during the early stage of transfer is obtained by solving Equation (A17).
(3) Pressure distribution during the late stage of transfer The whole matrix system enters the excitation state after the pressure propagates to the matrix center. According to Equation (A8), the pressure distribution in the matrix becomes as follows.
Combined with the results of the pressure distribution during the early stage of unsteady transfer, the initial and boundary conditions become as follows.
The equations of parameters (a 3 , a 4 , a 5 ) can be obtained by substituting the initial conditions Equation (A20a) and boundary conditions Equation (A20b), Equation (A20c) into the pressure distribution formula Equation (A19).
It is impossible to determine three parameters at the same time by only solving Equation (A21). There are residual errors after the model simplification.
When we make the residual expectations zero, where The relation between the derivatives of various parameters is obtained by taking the derivative of Equation (A21).
Therefore, the solution of the pressure distribution during the late unsteady stage is obtained.
where To find the average pressure of the matrix system, an integral is carried out over the whole matrix range.
The average pressure of the matrix system during the late stage of transfer is obtained by solving Equation (A33).