Using experimental design and response surface methodology to model induced fracture geometry in Shublik shale

Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

Contents lists available at ScienceDirect

Journal of Unconventional Oil and Gas Resources journal homepage: www.elsevier.com/locate/juogr

Using experimental design and response surface methodology to model induced fracture geometry in Shublik shale S. Poludasu ⇑, O. Awoleke, M. Ahmadi, C. Hanks University of Alaska, Fairbanks, United States

a r t i c l e

i n f o

Article history: Received 10 March 2015 Revised 26 April 2016 Accepted 26 April 2016 Available online 17 May 2016

a b s t r a c t In this study, we developed a methodology for identifying the critical variables needed for accurate planning of a hydraulic fracturing treatment in a shale resource play where much of the properties required for hydraulic fracture modeling remain unknown. The critical variables identified can thereafter be used to develop a proxy model that can be used in lieu of a numerical simulator. This study was conducted in two stages. In the first stage, we used 2-level fractional factorial designs and a pseudo-3D simulator to identify the most important variables affecting the simulated fracture geometry. The variables investigated included geologic, mechanical and treatment design parameters. Using the three most significant variables for each fracture geometry component identified from the first stage, the second stage of this study applied Box-Behnken experimental design and response surface methodology to quantify functional relationships between input variables and the fracture geometry. These proxy models, typically polynomial equations, can be used to predict the fracture geometry with very little computational time. The use of experimental design drastically reduces the number of simulations required to evaluate large number of variables. With only 137 simulations, 26 variables were ranked based on their statistical significance and non-linear proxy models were developed for the nine fracture geometry variables. Predicted values of the fracture geometry using the proxy models were in good agreement with the simulated values (R2 value of 0.99 for fracture length and fracture height and R2 value of 0.96 for fracture width). These linear and non-linear proxy models were validated by comparing the results from the proxies and the actual simulator using a random value dataset within the design space. The results indicate a good match for the width at the top and bottom of the fracture and propped fracture height/length. Engineers can use the results described here for quick estimates of fracture dimensions and the methodology outlined here can be used with more complicated fracturing models. Ó 2016 Elsevier Ltd. All rights reserved.

Introduction The ultra-low permeability of unconventional shale reservoirs demands large-scale stimulation treatments (multi-stage hydraulic fracturing of horizontal wells) in order to produce economically. Even with the significant technological advances in modern day multi-stage fracturing, development risk has not been eliminated for developing these shale reservoirs. The risk is largely associated with the limited knowledge of reservoir geology, presence, conductivity and connectivity of natural fractures and their influence on the fracture geometry. Typically, a hydraulic fracture design engineer uses a numerical simulator to predict the fracture geometry for a given reservoir. The choice of the simulator can vary from relatively simple and compu⇑ Corresponding author at: 6607 Lake Woodlands Dr APT 414, The Woodlands, TX 77382, United States. E-mail address: [email protected] (S. Poludasu). http://dx.doi.org/10.1016/j.juogr.2016.04.002 2213-3976/Ó 2016 Elsevier Ltd. All rights reserved.

tationally inexpensive 2D models to more complex full 3D models. 2D fracture propagation models assume that an induced fracture will extend vertically to the entire height of the pay zone, and remain within the pay zone while propagating laterally (Zeng, 2002). Pseudo 3D fracture propagation models are similar to 2D models, except that the fracture height is not constrained to the payzone thickness. These models also assume that the fracturing fluid flows in one dimension (from perforations to fracture tips) to induce an elliptical fracture (Zeng, 2002). Lastly, 3D fracture propagation models have no assumptions about the orientation of the fracture. They use the local stress field and fracture mechanics criteria to estimate the fracture propagation direction (Zeng, 2002). The choice of hydraulic fracturing model used entirely depends on the complexity of the reservoir geology and the availability of pertinent data. For simple systems, 2D equations can be used to estimate the fracture geometry. For more complicated systems, the use of either pseudo 3D or full 3D models is common or required. There is however a positive correlation between model complexity and

44

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

data requirements—as such, it is best to go with the simplest model that can fulfill the technical objectives. For a detailed treatment of existing hydraulic fracture geometry simulators, the reader is referred to the SPE Monograph on Hydraulic Fracturing. This work is focused on the development of a methodology or workflow that identifies critical variables needed for planning hydraulic fracture treatments in shale where little geologic information is available. We also developed linear and non-linear proxies for a pseudo-3D hydraulic fracture simulator. Given the dependence of fracture geometry on a large number of reservoir and treatment variables, using a typical numerical simulator to evaluate all the possible development scenarios would be time-intensive. Therefore, using a proxy as described above can greatly reduce the computational time and can be used to screen various scenarios. Experimental design coupled with response surface methodology is very efficient in extracting the maximum amount of information from relatively small subset of the simulation space. Experimental designs have been widely used in petroleum engineering studies (Awoleke et al., 2012; Segnini et al., 2014; Ambastha, 2014; Yu and Sepehrnoori, 2014). In order to fully investigate the space of factorial experiments with ‘n’ variables, we would require 2n simulation experiments (full factorial design). However, as the value of ‘n’ increases, the number of simulations required would also increase exponentially. If we limit our inquiry to being able to uniquely characterize the effects of each investigated variable (fractional factorial design), we can drastically reduce the number of simulations required by several orders of magnitude. The endpoint for this part of the work is to identify the statistically significant variables from our input variable set. In essence, we are using fractional factorial designs to identify and eliminate the non-statistically significant variables. We also developed a linear proxy based on the results of our simulations. However, because factorial and fractional factorial designs assume linearity, we ran another set of simulations using Box-Behnken designs and the three most significant variables (for each response variable) from the initial set of results. Using this second set of simulations, we developed a non-linear proxy. We concede that selecting only the top three significant variables would mean sacrificing some of the accuracy of the proxy for some of the response variables, as we will discuss later. Thus, in this study, we develop some functional relationships between the fracture design variables and the predicted fracture geometry by using experimental design, a pseudo 3D simulator and data from the Shublik shale (with the Eagle Ford of Texas as an analog whenever Shublik data is unavailable) of the Alaskan North Slope. The relationships developed can be used in lieu of the numerical simulator to quickly evaluate and rank various development scenarios.

complex, nonlinear systems described by computer models (Parikh, 2003). Experimental designs have also been used in petroleum engineering studies (Awoleke et al., 2012; Segnini et al., 2014; Ambastha, 2014; Yu and Sepehrnoori, 2014). Factorial design Consider a simulation study with ‘k’ input variables. Each input parameter is assigned a maximum or minimum value based on our engineering judgment. In other words, we have ‘k’ input variables in two levels (a higher value denoted by ‘+1’, and a lower value denoted by ‘
1’). The factorial design considers all the possible combinations of the input variables on both levels. This implies that the total number of simulations required in a factorial design with k factors is 2k. This design considers all the main effects and interaction effects of all the input variables. Main effect of an input parameter is the quantification of the variation in response with change in that input parameter alone. An interaction effect signifies the relative dependence of two or more input variables among themselves based on their shared effect on the response (Parikh, 2003). Fractional factorial design As the number of input variables increase, the number of simulation runs required using factorial design also increases exponentially. For such cases, fractional factorial designs are utilized. This design assumes that only main effects and few of the two/threefactor interactions of input variables have significant effect on the responses. By considering, only a subset of the factorial design, fractional factorial designs drastically reduces the number of simulations required to uniquely estimate the significance all the input variables on the responses (Parikh, 2003). A disadvantage of fractional factorial designs is that it assumes linearity between the input and response variables. Box-Behnken design Box-Behnken design is a rotatable quadratic design based on 3level fractional factorial design (Aslan and Cebeci, 2007). Each input factor is placed at one of the three equally spaced values, generally coded as
1, 0, +1 (lower, middle, and higher values of the input parameter range) as seen in Fig. 1. At least three levels are required for these designs as this design fits the data into a quadratic model. Since the design is quadratic, it does not assume linearity between the input and response variables.

Experimental design concepts Numerical models are widely used in engineering and scientific studies with the help of high performance computers. As a result, researchers have shifted to intricate mathematical models to simulate complex systems. The computer models often have multidimensional inputs, like scalars or functions. The output may also be multidimensional. Making a number of simulation runs at various input conditions is what is called a simulation experiment. ’Experimental design’ (ED) builds a response surface which is an empirical fit of computed responses as a function of input variables. ED is an efficient way to choose the input conditions that minimize the number of computer simulation runs required for data analysis, inversion problems and input uncertainty assessment and has been used in diverse areas such as aerospace, civil engineering and electronics for analysis and optimization of

Fig. 1. Box-Behnken design.

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

45

We used regression analysis to develop the proxies for both the simulation runs suing the fractional factorial design and the Box Behnken design. The typical regression equation is shown below.

^ y ¼ Xb

2

3 y1 6 y2 7 6 7 7 where y ¼ 6 6 : 7 is the response matrix computed from the simula4 : 5 yN 2 3 X11 X12 X13 ::::X1M 6 X21 X22 X23 ::::X2M 7 6 7 7 tor, X ¼ 6 6 :::::::::::::::::::::::::::: 7 is the matrix of simulation conditions 4 :::::::::::::::::::::::::::: 5 XN1 XN2 XN3 :::XNM determined using fractional factorial or Box-Behnken design. N is the number of simulation runs and M is the number of independent 2 3 ^1 b 6b ^ 7 6 27 ^¼6 : 7 variables and their interactions that are investigated, and b 6 7 4 : 5 ^N b
1 T T ^ ¼ ðX XÞ X y. is the regression coefficient matrix calculated using b ^ represents the statistical The regression coefficient matrix, b, effect of the input parameter on the outcomes in y. The significance of individual input variable on the outcome can be estimated through the absolute magnitude of its regression coefficient and by use of the p-value, as we will show later. A more detailed explanation about experimental design concepts/methodologies can be seen in the work of Poludasu (2014). Methodology The workflow (Fig. 2) followed in this study is divided into two stages. The first stage focuses on identifying the significant variables affecting the hydraulic fracture geometry. The second stage of the study uses the three most significant variables identified in the first stage to develop a functional relationship between these variables and the predicted fracture geometry using Box-Behnken experimental design methodology. The workflow can be summarized as follows:  Stage 1: Significant variable identification o Determine the input variables o Determine the response variables o Choose a range (minimum, maximum) for the input variables, based on a detailed literature search o Use fractional factorial design to plan the number of simulations o Estimate the fracture geometry variables (response variables) using the pseudo-3D numerical simulator o Identify the significant parameters for each response variable o Develop a proxy linear model.  Stage 2: Non-linear proxy model development o Based on the linear model, determine the statistical significance of the input variables evaluated in Stage 1. o Select the top three statistically significant variables for all the nine fracture geometry variables. o Using these significant variables, plan a 3-level BoxBehnken design. o Perform the planned simulations using pseudo 3D numerical simulator and generate the non-linear response surface or quadratic proxy model o Estimate the prediction accuracy of the developed proxy models using the validation dataset.

Fig. 2. Two-stage workflow used in this study using experimental design and response surface methodology.

Geology The Triassic Shublik Formation is calcareous shale and phosphatic limestone interval and a major source rock for hydrocarbon accumulations on the North Slope. At Prudhoe Bay, the Shublik Formation is thin (varying between 0 and 585 ft) and is bounded by the Eileen and the Sag River Sandstones (Parrish et al., 2001; Kelly et al., 2007). The Shublik is subdivided into four distinct zones in the subsurface (Fig. 3; Kupecz, 1995; Hulm, 1999). These zones are labelled A through D, from the top to base of the section. Zones A and C are organic rich, consist of black shale and dark grey limestone and are the target zones for stimulation. The thicknesses of zones A and C range from 0 to 83 ft and 0 to 46 ft, respectively. Zone B varies from 0 to 28 ft in thickness and is mainly composed of phosphorite, phosphatic carbonates and siliciclastic rocks. Lastly, Zone D is fine-to-medium-grained phosphatic sandstone with thicknesses ranging from 0 to 24 ft (Hulm, 1999). Development of the Shublik as a shale resource play has been hindered partly because of the unavailability of accurate estimates of reservoir and rock mechanical properties. For this reason, when the reservoir and mechanical properties of Shublik needed for simulation were unavailable, the properties of Eagle Ford shale, a geologic analog, were used. The Eagle Ford shale is a Cretaceous-age heterogeneous calcareous shale formation that is the source rock for the Austin Chalk Formation and the East Texas oil fields (Jiang, 1989). As seen in Table 1, Eagle Ford and Shublik appear to have similar Total

46

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

Organic Carbon (TOC) values and kerogen types. Both Eagle Ford and Shublik shales are calcareous shales, are brittle and have natural fractures. These similarities make the Eagle Ford a reasonable geologic analog for the Shublik, suggesting that hydraulic fracture treatments effective in the Eagle Ford could also be effective in the Shublik (Decker, 2011). These similarities also suggest that well and production data from the Eagle Ford can be used to investigate the Shublik’s response to a simulated fracture treatments (Hutton et al., 2012).

these layers was also included to evaluate the effects of fluid leak-off on fracture geometry. Treatment parameters like fracturing fluid type, proppant size/type, injection rate, and position of horizontal leg within the formation were also considered. Table 2 lists the 26 parameters evaluated in this study to model the nine fracture geometry variables (Fig. 5; Table 3). Each parameter is sampled at two levels and encoded as
1 (lower level) or +1 (higher level).

Input and output parameter selection

For the chosen 26 input variables as seen in Table 2, the total number of simulations required to evaluate the full design space (full factorial method) including the main factor effects and the interaction effects are 67108864 (226). By using fractional factorial design, we can uniquely compute the effect of each of the investigated factors (main effect) with just 32 simulations (226-21 fractional factorial design). Table 4 shows the 25 planning table for the 26 input variables. Each row of the planning table represents a simulation run. These planned simulations were implemented using a pseudo-3D fracturing simulator. For each simulation run, nine response variables (Table 3) were recorded. Using linear regression, a linear response surface model was developed for the 32 simulations. The absolute regression coefficients from this exercise can be seen in Table 5 presented in a Pareto chart format. The rows and columns of Table 5 are the input (independent) and the response (dependent) variables respectively. Thus, for example, the most important variable affecting the width at the top of the fracture is hu, the upper boundary layer thickness while the most important property affecting the width at the middle of the fracture is hr, the reservoir thickness. The rows of Table 6 represent the p-values from this analysis. The p-values are indicative of the relative significance of the parameters. The significant variables are also highlighted (at a significance level of 0.05) for the 26 parameters. Table 7 summarizes the linear proxy models with their respective R2 values. It is entirely possible to use this linear model to predict fracture dimensions and thus mimic the numerical simulator. But, this is not advisable as the 2-level design assumes a linear relationship between input and output parameters which is not true. However, using results from the 2-level linear model, we can rank the investigated parameters in order of the influence on the response. This step will greatly reduce the number of parameters considered for the next step, which is 3-level non-linear analysis. The top three significant parameters selected among the 26 parameters based on magnitude of the p-value (Table 6) from the linear model as seen in Table 8. These top 3 significant parameters screened are further used in 3-level Box-Behnken design to model the nine fracture geometry variables.

For this study, we selected 26 parameters including geologic, mechanical and treatment design parameters. The Shublik shale was assumed to be subdivided into four distinct zones of equal thickness (labelled A through D in Fig. 4) following (Kupecz, 1995, and Hulm, 1999). Zones A and C are organic rich, consist of organic-rich black shale and dark grey limestone and are considered as the target zones for stimulation (Hulm, 1999). It is also assumed that sandstones (Parrish et al., 2001) overlie and underlie the Shublik Formation. The mechanical properties considered for this study include Young’s modulus and Poisson’s ratio of the four Shublik sub-layers and the boundary layers. The permeability of

Stage 1—Significant parameter identification

Stage 2—Proxy model development Out of the 26 variables evaluated in Stage 1, based on the magnitude of the p-values from the 2-level linear model (Table 6), the Table 1 Similar characteristics of eagle ford shale and shublik shale (Decker, 2011).

Total organic carbon Main Kerogen types Oil gravity, API Thickness Thermal maturity Lithology and variability Natural fractures Overpressure Fig. 3. Lithostratigraphy and corresponding gamma ray log response for Shublik shale based on one of the wells in Prudhoe Bay (Hulm, 1999).

*

At a given outcrop.

Eagle Ford

Shublik

2–7% I/II (oil) 30–50 API 50–250 ft Imm-Oil–Gas Sh-Slts-Sh Yes Yes

*

2.40% (at a given outcrop) I/II-S (oil) 24 API 0–600 ft Imm-Oil–Gas Sh-Slts-Ls Yes Locally

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

47

Fig. 4. Geologic setting used in this study (Modified from Hulm (1999).

Table 2 26 parameters chosen to be investigated by this study. Code

Input parameter

Symbol

Lower level (
1)

Higher level (+1)

Sources of parameter ranges

A B C D E F G H I J K L M N O P Q R

Shublik A Young’s modulus Shublik A Poisson’s ratio Shublik B Young’s modulus Shublik B Poisson’s ratio Shublik C Young’s modulus Shublik C Poisson’s ratio Shublik D Young’s modulus Shublik D Poisson’s ratio Permeability in Shublik A, nD Reservoir depth, ft Upper layer Poisson’s ratio Upper layer Young’s modulus, MMpsi Upper layer thickness, ft Upper layer permeability, mD Lower layer Poisson’s ratio Lower layer Young’s modulus, MMpsi Lower layer thickness, ft Lower layer permeability, mD

ErA

6 0.26 3.6 0.23 6 0.26 2.95 0.38 800 13,500 0.38 3.55 1000 21 0.38 3.55 750 23

For the analog Eagle Ford from Manchanda et al. (2012)

El hl kl

1.5 0.22 3 0.18 1.5 0.22 1.5 0.21 1 8000 0.21 2.38 20 1.8 0.21 2.38 100 2.9

S T U V W X Y Z

Proppant type Proppant size Injection rate, bpm Fracturing fluid type Reservoir thickness, ft Well placement Permeability in Shublik C, nD Permeability in Shublik B & D, mD

Pt Ps q t hr Wp krC krBD

Sand 40/70 45 Slickwater 120 A 1 10

Ceramic 16/30 100 Cross-linked 550 C 800 100

mrA

ErB

mrB

ErC

mrC

ErD

mrD

krA D

mu Eu hu ku

ml

top three significant variables were chosen for further evaluation using the Box-Behnken experimental design (Table 8). The choice of just selecting the top three-significant variables was made because we wanted to limit the number of simulation experiments run during this phase of the study. Box-Behnken design In this study, the Box-Behnken experimental design was chosen for finding out the relationship between the top three significant variables (sampled at 3 levels) screened from the 2-level design (Table 8) and the fracture geometry variables. The Box-Behnken

The mechanical properties of typical carbonates For the analog Eagle Ford from Manchanda et al. (2012) The mechanical properties of typical sandstones For the analog Eagle Ford from Stegent et al. (2010) For the analog Eagle Ford from Centurion (2011) The mechanical properties of typical sandstones Upper layer thickness at 2 outcrops by Kelly et al. (2007) Permeability of the Ivishak sandstone from Miller et al. (2002) The mechanical properties of typical sandstones Lower layer thickness at 2 outcrops by Kelly et al. (2007) Permeability of the Sagriver sandstone from Johnston and Christenson (1998) From that data library of FracPro PT, numerical fracture simulator used in this study

Chosen from the 2 organic-rich shales from Hulm (1999) For the analog Eagle Ford from Stegent et al. (2010) Permeability range of 10–100 mD was chosen

design is a rotatable quadratic design based on 3-level fractional factorial design (Aslan and Cebeci, 2007). Each input factor is placed at one of the three equally spaced values, coded as
1, 0, +1, representing lower, middle (arithmetic average), and higher values of the input parameter range, respectively, as seen in Fig. 1. At least three levels are required for these designs as this design tries to fit into a quadratic model (non-linear model). For the 3 variables screened from the previous step, Box-Behnken design would require a total of 15 runs as seen in Table 9. After the regression analysis of results predicted by the numerical simulator for these 15 runs, the resultant non-linear model considering three variables would look in the following form:

48

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

Fig. 5. Output parameters estimation from the simulator generated fracture profile.

Table 3 Nine fracture geometry (Response) variables modeled in this study. Dependent parameter

Symbol

Width at the top of the fracture (in) Width at the middle of the fracture (in) Width at the bottom of the fracture (in) Fracture length (ft) Propped length (ft) Fracture height (ft) Propped height (ft) Upper fracture outgrowth (ft) Lower fracture outgrowth (ft)

width_top width_mid width_bot fracture_length propped_length fracture_height propped_height upper_outgrowth lower_outgrowth

y ¼ b0 þ b1 x1 þ b2 x2 þ b3 x3 þ b11 x21 þ b22 x22 þ b33 x23 þ b12 x1 x2 þ b13 x1 x3 þ b23 x2 x3 where y is predicted response (fracture geometry variables), b0 is the intercept of the model, b1 ; b2 ; b3 are the regression coefficients of the input variables x1 ; x2 ; x3 , b12 ; b13 ; b23 are the regression coefficients of two factor interaction terms, and b11 ; b22 ; b33 are the second-order regression coefficients. The highest order model

that we can develop based on the Box-Behnken design is a second order model. This is because the columns of the design matrix for the Box-Behnken design will become linearly dependent if we increase the order of the model. This in turn will lead to the development of erroneous regression models. Non-linear proxy models developed in Stage 2 A three parameter coded 3-level Box Behnken design was used to determine the relationships between the significant variables (Table 8) and nine fracture geometry variables. The non-linear models developed using the above technique is as shown in Table 8. Validation As seen in Table 7 and Table 10, both the linear and non-linear proxy models have reasonable prediction accuracy as seen from the regression coefficients. However to provide some type of validation for the models, we generated a random sample of independent variables (within the bounds of our simulation experiments),

Table 4 Planning table for the 26 parameters.

mrA

ErB

mrB

ErC

mrC

ErD

mrD

krA

D

mu

Eu

hu

ku

ml

El

hl

kl

Pt

Ps

q

t

hr

Wp

krC

krBD

1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 1.5 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6 6

0.22 0.22 0.22 0.22 0.22 0.22 0.22 0.22 0.26 0.26 0.26 0.26 0.26 0.26 0.26 0.26 0.22 0.22 0.22 0.22 0.22 0.22 0.22 0.22 0.26 0.26 0.26 0.26 0.26 0.26 0.26 0.26

3 3 3 3 3.6 3.6 3.6 3.6 3 3 3 3 3.6 3.6 3.6 3.6 3 3 3 3 3.6 3.6 3.6 3.6 3 3 3 3 3.6 3.6 3.6 3.6

0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23 0.18 0.18 0.23 0.23

1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6 1.5 6

0.22 0.26 0.26 0.22 0.26 0.22 0.22 0.26 0.26 0.22 0.22 0.26 0.22 0.26 0.26 0.22 0.26 0.22 0.22 0.26 0.22 0.26 0.26 0.22 0.22 0.26 0.26 0.22 0.26 0.22 0.22 0.26

2.95 1.5 1.5 2.95 1.5 2.95 2.95 1.5 1.5 2.95 2.95 1.5 2.95 1.5 1.5 2.95 2.95 1.5 1.5 2.95 1.5 2.95 2.95 1.5 1.5 2.95 2.95 1.5 2.95 1.5 1.5 2.95

0.38 0.21 0.21 0.38 0.21 0.38 0.38 0.21 0.38 0.21 0.21 0.38 0.21 0.38 0.38 0.21 0.21 0.38 0.38 0.21 0.38 0.21 0.21 0.38 0.21 0.38 0.38 0.21 0.38 0.21 0.21 0.38

1 800 800 1 800 1 1 800 1 800 800 1 800 1 1 800 1 800 800 1 800 1 1 800 1 800 800 1 800 1 1 800

13,500 8000 8000 13,500 13,500 8000 8000 13,500 8000 13,500 13,500 8000 8000 13,500 13,500 8000 8000 13,500 13,500 8000 8000 13,500 13,500 8000 13,500 8000 8000 13,500 13,500 8000 8000 13,500

0.21 0.38 0.38 0.21 0.21 0.38 0.38 0.21 0.38 0.21 0.21 0.38 0.38 0.21 0.21 0.38 0.21 0.38 0.38 0.21 0.21 0.38 0.38 0.21 0.38 0.21 0.21 0.38 0.38 0.21 0.21 0.38

2.38 3.55 3.55 2.38 2.38 3.55 3.55 2.38 2.38 3.55 3.55 2.38 2.38 3.55 3.55 2.38 3.55 2.38 2.38 3.55 3.55 2.38 2.38 3.55 3.55 2.38 2.38 3.55 3.55 2.38 2.38 3.55

1000 20 20 1000 1000 20 20 1000 1000 20 20 1000 1000 20 20 1000 1000 20 20 1000 1000 20 20 1000 1000 20 20 1000 1000 20 20 1000

21 1.8 21 1.8 1.8 21 1.8 21 1.8 21 1.8 21 21 1.8 21 1.8 1.8 21 1.8 21 21 1.8 21 1.8 21 1.8 21 1.8 1.8 21 1.8 21

0.21 0.38 0.21 0.38 0.38 0.21 0.38 0.21 0.38 0.21 0.38 0.21 0.21 0.38 0.21 0.38 0.21 0.38 0.21 0.38 0.38 0.21 0.38 0.21 0.38 0.21 0.38 0.21 0.21 0.38 0.21 0.38

2.38 3.55 2.38 3.55 3.55 2.38 3.55 2.38 2.38 3.55 2.38 3.55 3.55 2.38 3.55 2.38 3.55 2.38 3.55 2.38 2.38 3.55 2.38 3.55 3.55 2.38 3.55 2.38 2.38 3.55 2.38 3.55

750 100 750 100 100 750 100 750 750 100 750 100 100 750 100 750 750 100 750 100 100 750 100 750 750 100 750 100 100 750 100 750

2.9 23 2.9 23 2.9 23 2.9 23 23 2.9 23 2.9 23 2.9 23 2.9 23 2.9 23 2.9 23 2.9 23 2.9 2.9 23 2.9 23 2.9 23 2.9 23

Ceramic Sand Ceramic Sand Ceramic Sand Ceramic Sand Sand Ceramic Sand Ceramic Sand Ceramic Sand Ceramic Ceramic Sand Ceramic Sand Ceramic Sand Ceramic Sand Sand Ceramic Sand Ceramic Sand Ceramic Sand Ceramic

16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30 40/70 16/30

100 100 45 45 45 45 100 100 45 45 100 100 100 100 45 45 45 45 100 100 100 100 45 45 100 100 45 45 45 45 100 100

Slickwater Slickwater X-linked X-linked X-linked X-linked Slickwater Slickwater X-linked X-linked Slickwater Slickwater Slickwater Slickwater X-linked X-linked Slickwater Slickwater X-linked X-linked X-linked X-linked Slickwater Slickwater X-linked X-linked Slickwater Slickwater Slickwater Slickwater X-linked X-linked

120 120 550 550 550 550 120 120 120 120 550 550 550 550 120 120 550 550 120 120 120 120 550 550 550 550 120 120 120 120 550 550

C C A A A A C C C C A A A A C C C C A A A A C C C C A A A A C C

1 1 800 800 1 1 800 800 800 800 1 1 800 800 1 1 800 800 1 1 800 800 1 1 1 1 800 800 1 1 800 800

100 100 10 10 100 100 10 10 10 10 100 100 10 10 100 100 100 100 10 10 100 100 10 10 10 10 100 100 10 10 100 100

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

ErA

49

50

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

Table 5 Absolute regression coefficient values evaluated using 2-level fractional factorial design (also shows the parameter significance pareto chart).

Symbol

width_top

width_mid

width_bot

Fracture length

Propped length

Fracture height

Propped height

E rA ν rA E rB ν rB E rC ν rC E rD ν rD k rA D νu Eu hu ku νl El hl kl Pt Ps q t hr Wp k rC k rBD

0.016 0.060 0.078 0.050 0.007 0.099 0.002 0.051 0.035 0.087 0.089 0.061 0.304 0.005 0.084 0.105 0.161 0.008 0.120 0.066 0.031 0.045 0.115 0.027 0.139 0.044

0.137 0.029 0.086 0.012 0.156 0.048 0.003 0.003 0.039 0.067 0.029 0.162 0.087 0.006 0.102 0.029 0.104 0.010 0.148 0.043 0.076 0.127 0.277 0.007 0.189 0.063

0.174 0.083 0.020 0.035 0.147 0.015 0.062 0.022 0.030 0.103 0.061 0.100 0.056 0.081 0.159 0.026 0.009 0.007 0.228 0.024 0.061 0.091 0.041 0.145 0.127 0.053

0.000 0.041 0.024 0.059 0.041 0.020 0.005 0.006 0.015 0.015 0.015 0.009 0.063 0.002 0.015 0.054 0.012 0.025 0.011 0.033 0.043 0.029 0.145 0.031 0.080 0.007

0.008 0.039 0.028 0.052 0.050 0.037 0.005 0.010 0.021 0.013 0.022 0.008 0.068 0.009 0.013 0.062 0.003 0.023 0.028 0.040 0.045 0.036 0.184 0.029 0.105 0.006

0.138 0.014 0.018 0.028 0.001 0.014 0.060 0.010 0.033 0.105 0.025 0.039 0.066 0.045 0.079 0.063 0.058 0.031 0.023 0.043 0.036 0.037 0.300 0.006 0.077 0.003

0.141 0.017 0.024 0.012 0.034 0.056 0.033 0.017 0.006 0.037 0.011 0.052 0.001 0.036 0.021 0.011 0.067 0.010 0.039 0.003 0.033 0.049 0.407 0.009 0.076 0.050

Upper fracture growth 0.602 0.228 0.239 0.071 0.016 0.014 0.227 0.347 0.058 0.424 0.048 0.006 0.126 0.102 0.034 0.214 0.165 0.197 0.073 0.239 0.327 0.106 0.891 1.237 0.219 0.093

Lower fracture growth 0.440 0.097 0.020 0.041 0.247 0.214 0.204 0.247 0.135 0.140 0.010 0.015 0.064 0.010 0.140 0.064 0.015 0.204 0.214 0.135 0.097 0.054 1.283 0.440 0.020 0.092

Table 6 P-Values of the variables evaluated in phase III (significant variables (p-values <0.05) are highlighted).

and compared the output from the simulator and the linear and non-linear proxy models. The results are as shown in Figs. 6–14. We see that the results from the simulator and the non-linear proxy agree for most of the variables except for width at the middle of the fracture and the fracture outgrowths at the top and bottom of the target layer. We suspect it would be difficult to

use a proxy model to describe the relationship between the lower/upper outgrowth and the independent variables. This is because they have an essentially discontinuous relationship. For the width at the middle to the fracture, we think we will have more accuracy if more independent variables are included in the proxy development process.

Table 7 Linear proxy model summary. Linear proxy model

Prediction accuracy (R2)

Width at the top of the fracture

lnðwidth topÞ ¼
0:15
0:016ErA
0:06v rA
0:078ErB
0:05v rB þ 0:007ErC þ 0:099v rC
0:002ErD þ 0:051v rD þ 0:035krA þ 0:087D
0:089v u
0:061Eu þ 0:304hu þ 0:005ku
0:084v l
0:105El
0:161hl þ0:008kl
0:12P t þ 0:66P s þ 0:031q þ 0:045t þ 0:115hr
0:027W p þ 0:139krC þ 0:044krBD lnðwidth midÞ ¼ 0:306 þ 0:137ErA
0:029v rA
0:086ErB
0:012v rB
0:156ErC þ0:048v rC þ 0:003ErD
0:003v rD
0:039krA þ 0:067D
0:029v u þ 0:162Eu þ 0:087hu
0:006ku
0:102v l
0:029El
0:104hl
0:01kl
0:148P t þ 0:043P s
0:076q þ 0:127t þ0:277hr
0:007W p þ 0:189krC þ 0:063krBD lnðwidth botÞ ¼
0:177 þ 0:174ErA
0:083v rA
0:02ErB þ 0:035v rB
0:147ErC
0:015v rC
0:062ErD þ0:022v rD þ 0:03krA þ 0:103D
0:061v u þ 0:1Eu þ 0:056hu þ 0:081ku
0:159v l þ 0:026El
0:009hl þ0:007kl
0:228P t þ 0:024P s þ 0:061q þ 0:091t þ 0:041hr þ 0:145W p þ 0:127krC þ 0:053krBD lnðfracture lengthÞ ¼ 5:777 þ 0:0002ErA
0:041v rA
0:024ErB
0:059v rB
0:041ErC þ0:02v rC
0:005ErD þ 0:006v rD
0:015krA þ 0:015D
0:015v u
0:009Eu
0:063hu þ 0:002ku þ 0:015v l þ 0:054El þ0:012hl
0:025kl
0:011P t þ 0:033P s
0:043q þ 0:029t þ 0:145hr þ 0:031W p þ 0:08krC þ 0:007krBD lnðpropped lengthÞ ¼ 5:723 þ 0:008ErA
0:039v rA
0:028ErB
0:052v rB
0:05ErC þ0:037v rC
0:005ErD þ 0:01v rD
0:021krA þ 0:013D
0:022v u þ 0:008Eu
0:068hu þ 0:009ku þ0:013v l þ 0:062El
0:003hl
0:023kl
0:028P t þ 0:04P s
0:045q þ 0:036t þ 0:184hr þ 0:0291W p þ 0:105krC þ 0:006krBD lnðfracture heightÞ ¼ 5:874 þ 0:138ErA þ 0:014v rA
0:018ErB þ0:028v rB þ 0:001ErC
0:014v rC
0:06ErD
0:01v rD
0:033krA þ 0:105D þ 0:025v u þ 0:039Eu þ 0:066hu
0:045ku
0:079v l
0:063El
0:058hl þ 0:031kl þ0:023P t þ 0:043P s
0:036q
0:037t þ 0:3hr
0:006W p þ 0:077krC
0:003krBD lnðpropped heightÞ ¼ 5:74 þ 0:141ErA þ 0:017v rA þ 0:024ErB þ 0:012v rB
0:034ErC þ0:056v rC
0:033ErD þ 0:017v rD þ 0:006krA þ 0:037D þ 0:011v u þ 0:052Eu þ 0:001hu
0:036ku
0:021v l þ 0:011El
0:067hl þ 0:01kl
0:039P t þ0:003P s
0:003q þ 0:049t þ 0:407hr
0:009W p þ 0:076krC
0:05krBD lnðupper outgrowthÞ ¼ 3:243
0:602ErA
0:228v rA
0:239ErB
0:071v rB
0:016ErC
0:014v rC
0:227ErD þ0:347v rD
0:058krA þ 0:424D
0:048v u
0:006Eu
0:126hu þ 0:102ku þ 0:034v l
0:214El
0:165hl
0:197kl
0:073P t þ0:239P s
0:327q þ 0:106t
0:891hr
1:237W p þ 0:219krC þ 0:093krBD lnðlower outgrowthÞ ¼ 1:283
0:44ErA
0:097v rA þ 0:02ErB þ 0:041v rB
0:247ErC þ0:214v rC þ 0:204ErD þ 0:247v rD þ 0:135krA þ 0:14D
0:01v u
0:015Eu
0:064hu þ 0:01ku
0:14v l þ 0:064El þ0:015hl
0:204kl
0:214P t
0:135P s þ 0:097q
0:054t
1:283hr þ 0:44W p þ 0:02krC
0:092krBD

0.82

Width at the middle of the fracture Width at the bottom of the fracture Fracture length

Propped length

Fracture height

Propped height

Upper fracture outgrowth Lower fracture outgrowth

0.96

0.67

0.99

0.99

0.99

0.99

0.90

0.87

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

Modeled fracture geometry

51

Width at the bottom of the fracture

Fracture length

Propped length

Fracture height

Propped height

Upper fracture outgrowth

0.87

0.90

0.99

0.99

0.99

0.99

v

Lower fracture outgrowth

0.67

Width at the middle of the fracture

0.96

Upper layer thickness (hu) Lower layer thickness (hl) Permeability in Shublik C (krC) Reservoir thickness (hr) Permeability in Shublik C (krC) Upper layer Young’s Modulus (Eu) Proppant type (Pt) Shublik A Young’s modulus (ErA) Lower layer Poisson’s Ratio (ml) Reservoir thickness (hr) Permeability in Shublik C (krC) Upper layer thickness (hu) Reservoir thickness (hr) Permeability in Shublik C (krC) Upper layer thickness (hu) Reservoir thickness (hr) Shublik A Young’s modulus (ErA) Reservoir depth (D) Reservoir thickness (hr) Shublik A Young’s modulus (ErA) Permeability in Shublik C (krC) Well placement (Wp) Reservoir thickness (hr) Shublik A Young’s modulus (ErA) Reservoir thickness (hr) Shublik A Young’s modulus (ErA) Well placement (Wp)

v

Width at the top of the fracture

v

Most significant input variables

2 0:0086hu

Modeled fracture geometry

Prediction accuracy (R2)

Table 8 Three most important input variables affecting the modeled fracture geometry.

0.82

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

2 2 lnðwidth topÞ ¼ 0:141
0:0219hu þ 0:0798hl
0:036krC þ 0:0437ðhu Þðhl Þ þ 0:0614ðhu ÞðkrC Þ
0:0192ðhl ÞðkrC Þ

0:0506hl þ 0:0601krC 2 2 2 lnðwidth midÞ ¼ 0:5008 þ 0:6293hr þ 0:0068krC þ 0:0319Eu
0:1181ðhr ÞðkrC Þ
0:0026ðhr ÞðEu Þ þ 0:0834ðEu ÞðkrC Þ
0:2711hr þ 0:0358Eu
0:0557krC lnðwidth botÞ ¼ 0:2852
0:1947P t þ 0:1664ErA þ 0:1576 l
0:1676ðP t ÞðErA Þ
0:0148ðP t Þð l Þ
0:0323ðErA Þð l Þ
0:3983P 2t þ 0:2705E2rA þ 0:1112 2l 2 2 2 lnðfracture lengthÞ ¼ 5:507 þ 0:363hr
0:0023krC
0:112hl þ 0:0095ðhl ÞðkrC Þ þ 0:0017ðhr Þðhl Þ þ 0:0071ðkrC Þðhl Þ þ 0:2509hr þ 0:0062krC þ 0:0057hl 2 2 2 lnðpropped lengthÞ ¼ 5:507 þ 0:393hr
0:0023krC
0:146hu þ 0:0079ðhu ÞðkrC Þ þ 0:0048ðhr Þðhu Þ þ 0:0048ðkrC Þðhu Þ þ 0:2142hr þ 0:0033krC þ 0:0027hl 2 lnðfracture heightÞ ¼ 5:8093 þ 0:6365hr
0:0104ErA
0:0104D þ 0:0208ðhr ÞðErA Þ þ 0:008ðhr ÞðDÞ
0:0207ðErA ÞðDÞ
0:2228hr þ 0:0243E2rA
0:00468D2 2 2 2 lnðpropped heightÞ ¼ 5:7981 þ 0:647hr þ 0:009ErA
0:0017krCl
0:0019ðhr ÞðkrC Þ þ 0:0211ðhr ÞðErA Þ
0:004ðkrC ÞðErA Þ
0:2445hr þ 0:0047krC þ 0:0000565hl 2 2 2 lnðupper outgrowthÞ ¼ 3:96
0:999hr
0:835ErA
3:060W p
0:725ðhr ÞðErA Þ
0:127ðhr ÞðW p Þ
0:802ðErA ÞðW p Þ
1:736hr
0:904ErA
1:035W p 2 lnðlower outgrowthÞ ¼ 1:551
2:251hr þ 0:310ErA þ 1:090W p
0:101ðhr ÞðErA Þ
0:162ðhr ÞðW p Þ þ 0:519ðErA ÞðW p Þ
0:684hr
1:223E2rA
0:612W 2p

52

Lower fracture outgrowth

Upper fracture outgrowth

Propped height

Fracture height

Propped length

 These proxy models can be used to perform sensitivity analysis and understand the effect of all the input variables on the predicted fracture geometry and also to evaluate a particular treatment design by estimating fracture geometry in a time efficient manner. This quick estimation is useful in screening and ranking several available stimulation treatments.  They can be used to solve inversion problems (reverse calculating the input parameters based on the desired fracture geometry). This capability of the proxy models is very useful, especially in the case of Shublik shale, where the waterbearing boundary layers create a need to contain the fracture within the reservoir thickness. Therefore, reverse calculating the treatment properties by setting the fracture outgrowth to zero can be very helpful.

Table 10 Non-linear proxy model summary.

The proxy models developed for nine fracture geometry variables have the following possible applications.

Fracture length

Applications of proxy models

Width at the bottom of the fracture

1
1 1 0 0 1 0 1
1 0
1 0
1 0 0

Width at the top of the fracture

Input parameter 3

1 1 0 1 1 0 0
1 0 0 0 0
1
1
1

Width at the middle of the fracture

Input parameter 2

0 0 1 1
1
1 0 0
1 0 1 0 0 1
1

Non-linear proxy model

Input parameter 1

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15

Modeled fracture geometry

Run#

v

Table 9 Box-behnken design for three variables coded as
1 (minimum), 0 (median), +1 (maximum).

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

Fig. 6. Comparison of width_top results from non-linear and linear proxy models with the simulator for a random input data set.

Fig. 7. Comparison of width_mid results from non-linear and linear proxy models with the simulator for a random input data set.

Fig. 8. Comparison of width_bot results from non-linear and linear proxy models with the simulator for a random input data set.

 These proxy models can also be used to estimate the statistics of uncertainty. Monte Carlo simulations performed using the proxy model can estimate the P10, P50, and P90 scenarios.

General discussion about proxy models We have presented in the body of this work a way to rank data required for hydraulic fracture treatments especially in emerging unconventional plays like the Shublik. We also presented a methodology on how to develop proxy models for the same play. We did not concern ourselves with a critique of the various options engineers have as far as hydraulic fracture geometry models are concerned—as this was not our objective. However, we cannot

53

Fig. 9. Comparison of fracture_length results from non-linear and linear proxy models with the simulator for a random input data set.

Fig. 10. Comparison of propped_length results from non-linear and linear proxy models with the simulator for a random input data set.

Fig. 11. Comparison of fracture_height results from non-linear and linear proxy models with the simulator for a random input data set.

think of a reason why the methodology that we have presented here cannot be used with any other hydraulic fracture simulator, albeit with some modifications. As such, we do not attempt to validate the pseudo-3D hydraulic fracture model that we base our proxies on. It should suffice to say though that engineers in the industry routinely use this pseudo-3D simulator. We also emphasize that the proxy models that we have developed here are approximate and only useful for back of the envelope calculations. In Stage 2 of this study, only the top three statistically significant parameters were further evaluated using Box-Behnken design. This could be the reason why both the linear and nonlinear models have comparable results in the validation stage. This accuracy of the non-linear proxies can be improved by increasing

54

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

rock characteristics will play an important role in fracture initiation and propagation. However, we believe it is still very useful to use a pseudo-3D simulator to illustrate the power of this methodology. This methodology can be used with even greater impact with more complex models of hydraulic fractures. Conclusions

Fig. 12. Comparison of propped_height results from non-linear and linear proxy models with the simulator for a random input data set.

In stage 1 of this study, we used a pseudo 3D fracture propagation model and fractional factorial experiment design to identify the most significant reservoir variables affecting the modeled hydraulic fracture geometry for the Shublik Shale of the North Slope. Fractional factorial method and regression analysis were used to quantify the relative significance of each individual input parameter on the resultant modeled fracture geometry. An analysis of the results from this study indicates the following conclusions:  The statistical significance of the evaluated reservoir parameters to each fracture geometry variable was quantified based on the p-values.  Thicknesses of the shale reservoir and boundary layers are the critical variables affecting most of the modeled fracture geometry variables for the Shublik.  Fracturing fluid leak-off into the boundary layers has little or no effect on fracture geometry.  The mechanical properties of the 2 shale sub-layers have higher significance than the 2 limestone and sandstone sub-layers and the boundaries on the fracture geometry.  Well placement between the shale layers has a significant effect on the upper and lower fracture outgrowth.

Fig. 13. Comparison of outgrowth_top results from non-linear and linear proxy models with the simulator for a random input data set.

In stage 2 of this study, the Box-Behnken experimental design and response surface methodology were applied to model fracture geometry. With 137 simulation runs (32 simulations for screening and 15  7 = 105 simulations for non-linear model regression), 26 variables are evaluated and non-linear proxy models are developed for all the nine fracture geometry variables.  The non-linear model for fracture width (width_top, width_mid, and width_bot) predicted with R2 of 0.82, 0.97, and 0.67 respectively.  Similarly, R2 of 0.99 was achieved for the non-linear models for fracture length (fracture_length and propped_length) and fracture height (fracture_height and propped_height).  The non-linear model for fracture outgrowth (upper_outgrowth and lower_outgrowth) has R2 of 0.9 and 0.87 respectively.

Fig. 14. Comparison of outgrowth_bot results from non-linear and linear proxy models with the simulator for a random input data set.

the number of factors considered for Box-Behnken design. However, for the sake of reduced computation time, only three significant parameters were considered as a part of this study. Based on the simulation runs (during the pre-validation and the validation phases), we recommend that engineers do not use this methodology for response variables that tend to be discontinuous or choppy functions of the independent variables like the lower and upper fracture outgrowths. As stated above, the proxy model was developed using a pseudo-3D numerical fracturing simulator. This simulator neglect the effects of permeability anisotropy, stress shadowing, and natural fracture interactions. These effects depending on geology and

In the model validation stage, where a randomly generated input data was used to compare the results predicted by the proxy (linear and non-linear) to the results predicted by the simulator:  The validation results were a good match in the case of fracture_length, propped_length, fracture_height, and propped_height  The validation results were a bad match for both upper and lower fracture outgrowth. We suspect we would need to incorporate more independent variables to increase the levels of accuracy.

Acknowledgements The authors would like to acknowledge the financial support provided by the Alaska Department of Natural Resources during this study.

S. Poludasu et al. / Journal of Unconventional Oil and Gas Resources 15 (2016) 43–55

References Ambastha, A.K., 2014. Guidelines for Uncertainty Assessment using Reservoir Simulation Models for Green- and Brown-field Situations. In: SPE Nigeria Annual International Conference and Exhibition: Lagos, Nigeria, 5–7 August. SPE-172432-MS, 10.2118/172432-MS. Aslan, N., Cebeci, Y., 2007. Application of Box-Behnken Design and Response Surface Methodology for Modeling of Some Turkish Coals, Vol. 86. Fuel, Issues 1–2, January, Pp. 90–97, ISSN 0016-2361, http://dx.doi.org/10.1016/j.fuel.2006.06. 010. Awoleke, O.O., Romero, J.D., Zhu, D., et al., 2012. Experimental Investigation of Propped Fracture Conductivity in Tight Gas Reservoirs Using Factorial Design. In: SPE Hydraulic Fracturing Technology Conference, The Woodlands, Texas, 6– 8 February. SPE-151963-MS. http://dx.doi.org/10.2118/151963-MS. Centurion, S.M., 2011. Eagle Ford Shale: A Multistage Hydraulic Fracturing, Completion Trends and Production Outcome Study Using Practical Data Mining Techniques. In: SPE Eastern Regional Meeting, Columbus, Ohio, 17–19 August. SPE-149258-MS. http://dx.doi.org/10.2118/149258-MS. Decker, P.L., 2011. Source-Reservoir Oil Resources Alaskan North Slope. Alaska Department of Natural Resources, March 2001. (accessed 26.02.16). Hulm, E.J., 1999. Facies & Sequence Analysis of Eileen-Sag River Interval MS thesis. University of Alaska Fairbanks, Fairbanks, Alaska. Hutton, E., Agboada, D.K., Whalen, M.T., et al., 2012. Upper Triassic Shublik Formation of North Alaska – An Eagle Ford Type Shale Reservoir? AAPG Annual Convention and Exhibition, Long Beach, California, April 22–25 Search and Discovery Article #10437, . Jiang, M., 1989. Biostratigraphy and Geochronology of the Eagle Ford Shale, Austin Chalk, and Lower Taylor Marl in Texas Based on Calcareous Nanofossils MS thesis. Texas A&M University, College Station, Texas. Johnston, D.W., Christenson, R.N., 1998. The Application of BP Exploration (Alaska) Inc. to Present Testimony for Classification of a New Pool and the Establishment of Pool Rules for Development of the Sag River Oil Pool in the Milne Point Unit. Alaska Department of Administration, May 6 1998. (accessed 26.02.16). Kelly, L.N., Whalen, M.T., McRoberts, C.A., 2007. Sequence Stratigraphy and Geochemistry of the Upper Lower Through Upper Triassic of Northern Alaska:

55

Implications for Paleoredox History, Source Rock Accumulation, and Paleoceanography. Alaska Division of Geological & Geophysical Surveys Report of Investigation. http://dx.doi.org/10.14509/15773, 2007-1, p. 50. Kupecz, J.A., 1995. Depositional Setting, Sequence Stratigraphy, Diagenesis, and Reservoir Potential of a Mixed Lithology, Upwelling Deposit: The Upper Triassic Shublik Formation, Prudhoe Bay, Alaska. Am. Assoc. Pet. Geol. Bull. 79, 1301– 1319, http://archives.datapages.com/data/bulletns/1994-96/data/pg/0079/ 0009/1300/1301.htm?doi=10.1306%2F7834D4AE-1721-11D7-8645000102C18 65D. Manchanda, R., Roussel, N.P., Sharma, M.M., 2012. Factors Influencing Fracture Trajectories and Fracturing Pressure Data in a Horizontal Completion. In: 46th U.S. Rock Mechanics/Geomechanics Symposium, Chicago, Illinois, 24–27 June. ARMA-2012-633. Miller, E.L., Grantz, A., Klemperer, S., 2002. Tectonic Evolution of the Bering ShelfChukchi Sea-Arctic Margin and Adjacent Landmasses, first ed. Geological Society of America. Parrish, J.T., Whalen, M.T., Hulm, E.J., 2001. Shublik Formation Lithofacies, Environments, and Sequence Stratigraphy, Vol. 21. Petroleum Plays and Systems in the National Petroleum Reserve-Alaska, Arctic Alaska, U.S.A. http://dx.doi.org/10.2110/cor.01.01.0089, pp. 89–110. Poludasu, S., 2014. Using Experimental Design and Response Surface Methodology to Model Induced Fracture Geometry in Shublik Shale MS thesis. University of Alaska Fairbanks, Fairbanks, Alaska. Segnini, C.J., Rashwan, M., Hernandez, et al., 2014. Probabilistic Analysis by Applying Uncertainty Workflows: A Case Study of the Teak Field, East Coast Trinidad. In: SPE Energy Resources Conference: Port of Spain, Trinidad and Tobago, 9–11 June. SPE-169924-MS. http://dx.doi.org/10.2118/169924-MS. Stegent, N.A., Wagner, A.L., Mullen, J., et al., 2010. Engineering a Successful FractureStimulation Treatment in the Eagle Ford Shale. In: Tight Gas Completions Conference: San Antonio, Texas, 2–3 November. SPE-136183-MS. http://dx.doi. org/10.2118/136183-MS. Yu, W., Sepehrnoori, K., 2014. Sensitivity Study and History Matching and Economic Optimization for Marcellus Shale. In: SPE/AAPG/SEG Unconventional Resources Technology Conference: Denver, Colorado, 25-27 August. SPE-2014-1923491MS, http://dx.doi.org/10.15530/urtec-2014-1923491. Zeng, Z., 2002. Laboratory Imaging of Hydraulic Fractures Using Microseismicity MS thesis. University of Oklahoma, Norman, Oklahoma.