1 School of Petroleum Engineering, China University of Petroleum, East China, Qingdao 266580, China
2 Department of Geological Sciences and Engineering, Missouri University of Science and Technology, Rolla, MO 65409, USA
3 State Key Laboratory of Oil and Gas Reservoir Geology and Exploitation, Southwest Petroleum University, Chengdu 610500, China
A mathematical model of multistage fractured horizontal well (MsFHW) considering stimulated reservoir volume (SRV) was presented for tight oil reservoirs. Both inner and outer regions were assumed as single porosity media but had different formation parameters. Laplace transformation method, point source function integration method, superposition principle, Stehfest numerical algorithm, and Duhamel’s theorem were used comprehensively to obtain the semianalytical solution. Different flow regimes were divided based on pressure transient analysis (PTA) curves. According to rate transient analysis (RTA), the effects of related parameters such as SRV radius, storativity ratio, mobility ratio, fracture number, fracture half-length, and fracture spacing were analyzed.
The presented model and obtained results in this paper enrich the performance analysis models of MsFHW considering SRV.
Development of unconventional gas and oil has accelerated in recent years as conventional reservoirs have become increasingly difficult to find and produce [1–3]. Unconventional resources mainly include shale gas and oil, coal bed methane, tight gas and oil, and heavy oil. As one part of “unconventionals,” tight oil has garnered a lot of attention both in North American and Asia. Some of the noteworthy tight oil plays in North America include the Barnett, Haynesville, Marcellus, Eagle Ford, and Bakken. In China, tight oil also distributes widely such as in Ordos, Junggar, Songliao, Sichuan, and Qaidam Basin. The amount of reserves is about (80∼100) ×108 t [4–7]. To economically produce these hydrocarbons, unconventional methods are required. As an effective technique, multistage fractured horizontal well (MsFHW) has been widely used due to its advantages, such as creating high flow channels for liquids to flow into the well and increasing drainage area.
The creation of large complex fracture networks by hydraulic fracturing is imperative in many unconventional reservoirs by microseismic mapping observation. These networks are defined as stimulated reservoir volume (SRV) [8,9]. The fracture networks can happen when the differences between the principle stresses are small. SRV benefits oil or gas production because of its high conductivity [10, 11]. Some tight oil or gas reservoirs such as Eagle ford, Sulige, and Barnett have obtained high production by applying MsFHW with SRV [12–15].
To study flow mechanisms of MsFHW with SRV, both analytical and numerical methods are used. (1) Analytical method: Ozkan’ research group [16–19] proposed the “trilinear” model to study MsFHW performance. Pressure transient analysis was obtained in unconventional gas reservoirs. Stalgorova and Mattar [20, 21] improved the trilinear flow model. Five regions were defined to simulate the simulated reservoir volume. In their model, SRV was simplified into a simulated region of limited width. Both “trilinear” model and “five regions” models are based on the assumption that flow obeys “linear flow” in different regions. This assumption may ignore some flow regimes for MsFHW. Then, Ketineni and Ertekin  used equivalent flow model to describe SRV. In their model, reservoir was approximated as composite naturally fractured. The mathieu modified functions were used to solve the elliptical flow problem and some key factors were analyzed, such as mobility ratio, diffusivity ratio, storativity ratio, and interporosity flow coefficient ratios. Similar to this model, Zhao et al.  used a circular region to characterize the SRV in tight gas reservoir and the pressure transient response was given considering the effect of SRV. (2) Numerical simulation method: Mayerhofer et al. [8, 15] used numerical simulator to characterize SRV explicitly. Impacts of fracture network properties including fracture network size, network density, fracture conductivity, matrix permeability, and gaps in the network on well performance were studies. Meyer and Bazan  provided the foundation for predicting the behavior of discrete fracture networks. Wang et al.  used numerical simulator to study flow regimes considering SRV. Five regimes were divided and the characteristics in each regime were given.
All above work is meaningful to understand MsFHW performance with SRV. To our best knowledge, there are few models which can successfully calculate the performance behavior including pressure transient analysis (PTA) and rate transient analysis (RTA) for MsFHW with SRV in tight oil reservoirs. We try to solve this problem using the composite reservoir model. Lots of work has been done on RTA or PTA for composite reservoirs [26–29]. These studies mainly focused on vertical well or horizontal well and few are related to multistage fractured horizontal well. In our paper, we extended the composite model to MsFHW with SRV in tight oil reservoirs. We assumed the inner region was a single porosity medium to characterize SRV by equivalent continuum model [30, 31]. Another contribution in our work is that we try to obtain the solution with point function. Line source function was obtained by point function integration. Compared with the line source solution proposed by Zhao et al.  which can only be used to fully penetrating fracture, point function is more practical for fully/partially penetrating fracture or planar/bending fracture. Besides, we try to give the flow regimes for MsFHW and parameters effect analysis which is meaningful for well performance diagnosis. This paper is organized as follows: in Section 2, we established the composite model to describe the SRV for tight oil reservoirs; in Section 3, the point source function was obtained and used for semi analytical solution o fthe proposed model; in Section 4, the pressure transient analysis (PTA) and rate transient analysis (RTA) curves for MsFHW were discussed and the effects of related parameters were analyzed as well.
The schematic diagram for MsFHW with SRV is shown in Figure 1 . The reservoir has two regions: inner and outer region, which have different reservoir properties. The inner region is a medium including matrix and induced fracture. The equivalent continuum model concept is used to describe the inner region as a single porosity medium which has high flow capacity. The outer region is another single porosity medium which is not influenced by hydraulic fracturing. The model assumptions are as follows: (1) the outer region of a circular reservoir is infinite and the inner region radius is r1 ; (2) the reservoir is horizontal with uniform thickness of h and original pressure p𝑖; (3) for the inner region, the horizontal permeability is Kh1, the vertical permeability is KV1 , the compressibility is C𝑡1 , and the porosity is 𝜑1 ; while for the outer region, they are Kh2 , KV2 , C𝑡2 , and 𝜑2 . (4) The influence of gravity and capillary forces are ignored; (5) Wellbore storage effect and formation damage are taken into account.
The purpose of hydraulic fracturing is to create high conductivity around the wellbore that means inner region conductivity is higher than outer region conductivity. The left picture of Figure 1 shows the geology schematic after hydraulic fracturing. After hydraulic fracturing, permeability near the fractures will increase because the high flow channels (region in the green circle) form near the multistage fractured horizontal well. The outer region (region out of the green circle) is not affected by hydraulic fracturing. Thus the flow capacity is lower than the inner region. This analysis suggests the scenario that the MsFHW completely in the inner region should be considered (Figure 1). The proposed model is effective when conductivity capacity around the well improves greatly after hydraulic fracturing. It is worth noting that this paper assumes the inner region is a circle one which is used to characterize the SRV. In fact, the SRV is rather complex and has no regular shape in tight oil reservoirs. On the other hand, the inner region is assumed as an equivalent continuum single porosity model which has high conductivity, so the induced fracture cannot be characterized explicitly. It is necessary to conduct some further researches considering the above issues. These are some limitations of the model.
With radial cylindrical coordinate system, the flow equation can be expressed as follows. Flow equation in inner region is
Flow equation in outer region is
Inner boundary condition is
point source in inner region. Outer boundary is
Interface condition is
Dimensionless variables are defined as follows:
Using the dimensionless variables and Laplace transformation to (1)–(5), the following equations can be derived:
Inner boundary is:
point source in inner region. Outer boundary is
Interface condition is
where q is the point source, m3/s.
From the characteristics of the Bessel functions and the point source functions in pure medium described previously [32, 33], pressure distribution in the inner region can be expressed as follows when the point source is in the inner region:
where (rD, θ,z1D) is the pressure observation point and (r'D, θ',z'1D) point source position.
By integrating of the point source function with respect to z1D over the interval 0 to h1D, we can obtain the pressure distribution for a line source function:
where qL is the line source, m3/s and RD is the distance between the line source and pressure observation point in horizontal direction.
To obtain the pressure or rate behaviors for MsFHW, the assumptions are made as follows [34–36]: the wellbore is intersected by N fractures and all fractures are transverse to the wellbore as shown in Figure 2; flow from the reservoir to the wellbore is negligible compared to the flow from the hydraulic fractures; the well is assumed to produce at a constant wellbore pressure or at a constant production rate and the horizontal well is treated as infinite conductivity one.
The fracture is discretized to make the solution possible as shown in Figure 3. Each fracture includes 2n unites. Each unit can be taken as a line source. For the jth discrete unit of the ith fracture, the pressure distribution caused by the line source is
Equation (14) becomes
By applying the principle of superposition, the pressure response for N ×2n discrete segments can be obtained as follows:
Thus, the pressure response at the discrete segment (K, v) can be obtained:
where (xDk,v, yDk,v) is the center coordinate of the jth discrete unit in the ith fracture.The assumptions of infinite conductivity wellbore and fractures result in that the pressure at each point within the fractures and the horizontal wellbore is identical to bottomhole pressure, pωD.Thus, the wellbore pressure drop may be expressed as
The assumption of constant flow rate gives the following condition that must hold at anytime:
There are 2n∗(N+1)equations which can solve the 2n∗(N+1) unknowns of pwD, qD1, qD2,...,qD2n*N. A system of equations which is formed by (19) and (20) is now obtained and solution of such a system produces values for bottom hole pressure distribution as well as flux distribution for each fracture in Laplace domain.
The solution can be inverted back into real time domain using Stehfest algorithm. The Gauss elimination method was used to solve the system of equations. Then the Stehfest numerical inversion algorithm  was chosen to calculate the dimensionless bottom hole pressure as well as the dimensionless production rate distribution in real time space. For the numerical integration, Gauss-Legendre method was used [38 ].
Using Duhamel’s principle , we can obtain the solution considering well storage effect and skin effect:
When the well is producing at a constant bottom-hole pressure, the dimensionless well production rate can be defined as follows :
The code of the framework was programed by MATLAB2013a. Then type curves can be analyzed.
In this section, the dimensionless pressure and derivative responses for a multistage fractured horizontal well in tight oil reservoirs are calculated with the model proposed above. The effects of relevant parameters on rate transient responses are studied.
As shown in Figure 4, the type curve of MsFHW considering SRV can be divided into the following nine regimes: (1) the early wellbore storage characterized by a slope of 1 in the pressure curve and derivative curve;(2) the first transition flow period between wellbore storage and the early linear flow; (3) the early linear flow characterized by a slope of 1/2 in the pressure derivative curve (during the first linear flow, each fracture produces independently with the oil flowing perpendicular to the fracture as shown in Figure 5(a) ); (4) the second transition flow between the early linear flow and the early radial flow; (5) the early radial period around individual fractures marked by a constant pressure derivative value which is 1/(2*N*M) (this regime can occur if the fracture spacing is large enough compared to the fracture half-length as shown in Figure 5(b); this flow regime may not be observed if the fracture length is relatively long compared to the fracture spacing); (6) the second linear flow period characterized by a slope of 1/2 in the pressure derivative curve (this flow regime is shown in Figure 5(c)); (7) the second radial flow characterized by a horizontal straight line in the derivative curve and the value is 1/(2M) (this flow regime is shown in Figure 5(d); this regime can occur if the SRV radius is large enough; otherwise, it may be covered up such as when r1D=8000 as shown in Figure 6); (8) the fourth transition flow period between the second radial flow and the third radial flow; (9) the third radial flow characterized by a 1/2 horizontal straight line. This flow regime is shown in Figure 5(e). It can be seen that compared with the type curves (the dashed blue lines in Figure 4) of homogenous reservoir model, two more flow regimes appear—regimes (8) and (9) for MsFHW considering SRV. It should be point out that all these regimes are the reflection of the MsFHW and the formation parameters. If conditions change, these regimes will not be complete. An example is shown in Figure 6,when r1D decreases from 15000 to 8000, the second radial flow will disappear, and when decreasing to 3500, the second linear flow and second radial flow will disappear. Figure 7 is another example which presents the effect of fracture half-length on pressure transient responses.
When fracture half-length increases from 300 to 600, the first radial flow regime will disappear. From the two examples we can see that if parameters (fracture parameters and formation parameters) change, not all nine regimes will appear but depend on the combination of those related parameters.
Figure 8 shows effect of SRV radius on rate transient curves when the well produces at a constant wellbore pressure. The values of relevant parameters are listed as follows: CD=1000, S=0,02, Lf=100m, N=3M, M=15, L=0,1m, DL= 300m, η=15, r1=350m, 550m, 750m, 950m. It can be seen the SRV radius has significant effect on the production rate during the practical life of the well when all other parameters remain constant. As the SRV radius becomes larger, there will be a greater production rate. For the tight oil reservoirs, because of the extremely low permeability, flow resistance is high for oil. SRV is where the induced fractures exist and these fractures make the oil flow capacity increase to obtain economic development. Large SRV radius leads to large “high permeable” area and thus decreases the flow resistance. Therefore, in the development of tight oil reservoir it is desirable to obtain a large SRV radius to obtain a high production rate.
Figure 9 shows effect of storativity ratio on rate transient curves when the well produces at a constant wellbore pressure. The values of relevant parameters are listed as follows: 𝐶𝐷=10000, 𝑆=0.02, 𝐿𝑓=100m, N=3, M=15, L=0.1m, D𝐿=300m, r1=350m, 𝜔=0.2, 0.5, 1, 2, 5. Large storativity ratio leads to large production rate. Large storativity ratio means that the compressibility of the inner region is large. The production rate will become high under the same pressure difference.
Figure 10 shows effect of mobility ratio on rate transient curves when the well produces at a constant wellbore pressure. The values of relevant parameters are listed as follows: C𝐷=10000, S=0.02, L𝑓= 100 m, N=3, L=0.1m, D𝐿 =300 m, r1=350m, 𝜔=1, M=1, 5, 10, 15. As can be seen, when the mobility ratio is bigger than one, the production rate increases greatly. The blue line is rate transient curve when the mobility ratio is equal to one which means the reservoir is homogenous and there is no SRV near the well. The production rate is very low for the homogenous reservoir which indicates the necessity of the SRV. Another conclusion can be drawn from Figure 10 is that larger mobility ratio leads to larger production rate. Large mobility ratio represents large inner permeability and thus decreases flow resistance. So production rate increases effectively. When developing tight oil reservoirs, we hope to improve the inner region permeability as possible as we can to obtain high production rate.
Figure 11 shows effect of fracture number on rate transient curves when the well produces at a constant wellbore pressure. The values of relevant parameters are listed as follows: C𝐷=10000, S= 0.02, L𝑓=100m, L=0.1 m, D𝐿=80m, r1=350m, ω=1, N=3,5,7, 9. The fractures are assumed to be equally spaced and the properties of the fractures are identical in this case. As shown in Figure 11, increasing the number of fractures mainly influences the early production rate. Increasing the number of hydraulic fractures will improve the permeability around the wellbore. Flow resistance in the vicinity of the wellbore will be small. At the same time, more fractures increase the contact area with the formation for well and thus high production rate can be obtained.
Figure 12 shows effect of fracture half-length on rate transient curves when the well produces at a constant wellbore pressure. The values of relevant parameters are listed as follows: C𝐷=10000, S=0.02, N=3, M=15, L=0.1m, D𝐿=300m, r1=350m, ω=1, L𝑓=20m,50m, 80m,120m. The fractures are assumed to be equally spaced and the properties of the fractures are identical in this case. As shown in Figure 12, fracture half-length significantly influences the early production rate. Increasing fracture half-length can increase the contact area between the MsFHW and formation which is beneficial for the oil production.
Figure 13 shows effect of fracture spacing on rate transient curves when the well produces at a constant wellbore pressure. The values of relevant parameters are listed as follows: C𝐷=10000, S= 0.02, L𝑓=100 m, N=3, M=15, L= 0.1m, r1=350m, ω=1, D𝐿=50 m, 100m, 200m, 300m. The properties of the fractures are identical in this case. The well production performance indicates that large fracture spacing increases the production rate in early time. The reason is that when fracture number remains unchanged, larger fracture spacing leads to larger drainage area, and oil production rate will increase.
In this paper, a new model was presented for multistage fractured horizontal well considering simulated reservoir volume in tight oil reservoirs. The solution was obtained with point function. Pressure transient responses and rate transient responses were discussed. The following conclusions can be drawn.
(1) Oil flow in tight reservoir is a comprehensive result coupling the MsFHW and formation. The mathematical model is verified to describe the flow in both MsFHW and formation. Specifically, the SRV is taken into consideration compared to the existing PTA and RTA methods in tight oil reservoirs. Nine flow regimes are identified from the transient pressure curves. These regimes may not be complete when the formation properties and hydraulic fracture parameters change.
(2) The SRV has significant effect on rate transient curves. Large SRV radius leads to high production rate. That means a large SRV should be created for hydraulic fracture treatment in development of tight oil reservoirs. Permeability in the SRV region affects not only the early production rate, but also the production rate in later flow period. If mobility ratio is larger than 1, the rate is much bigger than that of the homogenous model. Large storativity ratio leads to large compressibility of inner region and high production rate can be obtained.
(3) The hydraulic fracture properties have significant effect on well early production rate. But the late-time behaviors are not affected. Increasing fracture number, fracture half-length, fracture spacing will increase the drainage area and improve the permeability around the wellbore. So if high early production rate is needed, the hydraulic fracture properties should be considered.
Ruizhong Jiang, Jianchun Xu, Zhaobo Sun, Chaohua Guo and Yulong Zhao declare that there is no conflict of interests regarding the publication of this paper.
This work is supported by the National Natural Science Foundation of China (Grant no. 51174223, no. 51374181, and E0403), the Fundamental Research Funds for the Central University under 14CX06087A, the graduate innovation fund of China University of Petroleum (East China) under CX-1210 and YCX2014017, the National Science and Technology Major Project under 2011ZX05013-006 and 2011ZX05051. The authors would also like to express their gratitude for reviewers for their careful review of this paper.
Correspondence should be addressed to Zhaobo Sun; firstname.lastname@example.org
Received 6 August 2014; Revised 23 September 2014; Accepted 7 October 2014; Published 10 November 2014
Academic Editor: Yang Tang
Copyright © 2014 Ruizhong Jiang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
 C. Jia, M. Zheng, and Y. Zhang, “Unconventional hydrocarbon resources in China and the prospect of exploration and development,” Petroleum Exploration and Development, vol. 39, no. 2, pp. 139–146, 2012.
 P. H. Nelson, “Pore-throat sizes in sandstones, tight sandstones, and shales,” AAPG Bulletin, vol. 93, no. 3, pp. 329–340, 2009.
 L. Xu, Y. Shi, C. Xu, Y. Yang, H. Li, and Z. Chai, “Influences of feldspars on the storage and permeability conditions in tight oil reservoirs: A case study of Chang-6 group, Ordos Basin,” Petroleum Exploration and Development, vol. 40, no. 4, pp. 481–487, 2013.
 L. Kuang, Y. Tang, D. Lei et al., “Formation conditions and exploration potential of tight oil in the Permian saline lacustrine dolomitic rock, Junggar Basin, NW China,” Petroleum Exploration and Development, vol. 39, no. 6, pp. 700–711, 2012.
 National Energy Administration, SY/T 6832-2011 Oil and Gas Industry Standard of the People’s Republic of China, Petroleum Industry Press, Beijing, China, 2011.
 W. Mingjian, H. Dengfa, B. Hongping, L. Renqi, and G. Baoling, “Upper Palaeozoic gas accumulations of the Yimeng Uplift, Ordos Basin,” Petroleum Exploration and Development, vol. 38, no. 1, pp. 30–39, 2011.
 W. Hu, Theory of Low-Permeability Reservoir, Petroleum Indus- try Press, Beijing, China, 2009.
 M. J. Mayerhofer, E. P. Lolon, C. Rightmire, D. Walser, C. L. Cipolla, and N. R. Warplnskl, “What is stimulated reservoir volume?” SPE Production and Operations, vol. 25, no. 1, pp. 89–98, 2010.
 H. Wang, X. Liao, N. Lu et al., “A study on development effect of horizontal well with SRV in unconventional tight oil reservoir,” Journal of the Energy Institute, vol. 87, no. 2, pp. 114–120, 2014.
 E. Stalgorova and L. Mattar, “Analytical model for history matching and forecasting production in multifrac composite systems,” in Proceedings of the SPE Canadian Unconventional Resources Conference (CURC ’12), pp. 450–466, November 2012.
 C. R. Clarkson, “Production data analysis of unconventional gas wells: review of theory and best practices,” International Journal of Coal Geology, vol. 109-110, pp. 101–146, 2013.
 R. Hull, H. Bello, L. P. Richmond et al., “Variable Stimulated Reservoir Volume (SRV) simulation: eagle ford shale case study,” in SPE Unconventional Resources Conference-USA, Society of Petroleum Engineers, The Woodlands, Tex, USA, 2013.
 D. K. Agboada and M. Ahmadi, “Production decline and numerical simulation model analysis of the eagle ford shale play,” in Proceedings of the SPE Western Regional & AAPG Pacific Section Meeting Joint Technical Conference, Society of Petroleum Engineers, 2013.
 H. Wang, C. Liao, H. Ye et al., “The performance evaluation of old well after SRV in Ordos basin tight oil reservoir,” in Proceedings of the SPE Energy Resources Conference, SPE-169968-MS, Society of Petroleum Engineers, Port of Spain, Trinidad and Tobago, June 2014.
 M. J. Mayerhofer, E. P. Lolon, J. E. Youngblood, and J. R. Heinze, “Integration of microseismic fracture mapping results with numerical fracture network production modeling in the Barnett shale,” in Proceedings of the SPE Annual Technical Conference and Exhibition (ATCE ’06), pp. 976–983, San Antonio, Tex, USA, September 2006.
 E. Ozkan, M. Brown, R. Raghavan, and H. Kazemi, “Comparison of fractured horizontal-well performance in conventional and unconventional reservoirs,” in Proceedings of the SPE Western Regional Meeting, pp. 345–360, Society of Petroleum Engineers, March 2009.
 E. Ozkan, L. M. Brown, R. Raghavan et al., “Comparison of fractured-horizontal-well performance in tight sand and shale reservoirs,” SPE Reservoir Evaluation & Engineering, vol. 14, no.02, pp. 248–259, 2011.
 M. Brown, E. Ozkan, R. Raghavan, and H. Kazemi, “Practical solutions for pressure-transient responses of fractured hori- zontal wells in unconventional shale reservoirs,” SPE Reservoir Evaluation and Engineering, vol. 14, no. 6, pp. 663–676, 2011.
 O. G. Apaydin, E. Ozkan, and R. Raghavan, “Effect of discontinuous microfractures on ultratight matrix permeability of a dual- porosity medium,” SPE Reservoir Evaluation & Engineering, vol.
15, no. 4, pp. 473–485, 2012.
 E. Stalgorova and L. Mattar, “Practical analytical model to simulate production of horizontal wells with branch fractures,” in Proceedings of SPE Canadian Unconventional Resources Conference, Society of Petroleum Engineers, Alberta, Canada, October-November 2012.
 E. Stalgorova and L. Mattar, “Analytical model for history matching and forecasting production in multifrac composite systems,” in Proceedings of the SPE Canadian Unconventional Resources Conference (CURC ’12), pp. 450–466, Society of Petroleum Engineers, November 2012.
 S. P. Ketineni and T. Ertekin, “Analysis of production decline characteristics of a multistage hydraulically fractured horizontal well in a naturally fractured reservoir,” in Proceedings of the SPE Eastern Regional Meeting, Society of Petroleum Engineers, Lexington, Ky, USA, October 2012.
 Y.-L. Zhao, L.-H. Zhang, J.-X. Luo, and B.-N. Zhang, “Performance of fractured horizontal well with stimulated reservoir volume in unconventional gas reservoir,” Journal of Hydrology, vol. 512, pp. 447–456, 2014.
 B. R. Meyer and L. W. Bazan, “A discrete fracture network model for hydraulically induced fractures: theory, parametric and case studies,” in Proceedings of the SPE Hydraulic Fracturing Technology Conference, pp. 571–606, January 2011.
 H. Wang, X. Liao, and X. Zhao, “Study of tight oil reservoir flow regimes in different treated horizontal,” Journal of the Energy Institute, 2014.
 J. S. Olarewaju and J. W. Lee, “A comprehensive application of a composite reservoir model to pressure-transient analysis,” SPE Reservoir Engineering, vol. 4, no. 3, pp. 325–331, 1989.
 O. Ezulike and A. Igbokoyi, “Horizontal well pressure transient analysis in anisotropic composite reservoirs—a three dimen- sional semi-analytical approach,” Journal of Petroleum Science and Engineering, vol. 96-97, pp. 120–139, 2012.
 A. Satman, “Pressure-transient analysis of a composite naturally fractured reservoir,” SPE Formation Evaluation, vol. 6, no. 2, pp. 169–175, 1991.
 L.-H. Zhang, J.-J. Guo, and Q.-G. Liu, “A new well test model for a two-zone linear composite reservoir with varied thicknesses,” Journal of Hydrodynamics, vol. 22, no. 6, pp. 804–809, 2010.
 D. T. Snow, “Rock fracture spacings, openings, and porosities,” Journal of the Soil Mechanics and Foundations Division, vol. 94, no. 1, pp. 73–92, 1968.
 D. T. Snow, “Anisotropie permeability of fractured media,” Water Resources Research, vol. 5, no. 6, pp. 1273–1289, 1969.
 E. Ozkan and R. Raghavan, “New solutions for well-test- analysis problems. Part 1. Analytical considerations,” SPE Formation Evaluation, vol. 6, no. 3, pp. 359–368, 1991.
 E. Ozkan and R. Raghavan, “New solutions for well-test- analysis problems: part III-additional algorithms,” in Proceedings of the SPE Annual Technical Conference and Exhibition, Society of Petroleum Engineers, 1994.
 Y.-L. Zhao, L.-H. Zhang, J.-Z. Zhao, J.-X. Luo, and B.-N. Zhang, ““Triple porosity” modeling of transient well test and rate decline analysis for multi-fractured horizontal well in shale gas reservoirs,” Journal of Petroleum Science and Engineering, vol.110, pp. 253–262, 2013.
 H.-T. Wang, “Performance of multiple fractured horizontal wells in shale gas reservoirs with consideration of multiple mechanisms,” Journal of Hydrology, vol. 510, pp. 299–312, 2014.
 J. Guo, L. Zhang, H. Wang, and G. Feng, “Pressure transient analysis for multi-stage fractured horizontal wells in shale gas reservoirs,” Transport in Porous Media, vol. 93, no. 3, pp. 635–653, 2012.
 H. Stehfest, “Algorithm 368: numerical inversion of Laplace transforms,” Communications of the ACM, vol. 13, no. 1, pp. 47–49, 1970.
 A. Iserles, A First Course in the Numerical Analysis of Differential Equations, Cambridge University Press, Cambridge, UK, 1996.
 A. F. van Everdingen and W. Hurst, “The application of the Laplace transformation to flow problems in reservoirs,” Journal of Petroleum Technology, vol. 1, no. 12, pp. 305–324, 1949.