P90 Calculation 4pa2t

SPE 162875 Analytical Probabilistic Reserve Estimation by Volumetric Method and Aggregation of Resources C. Karacaer, SPE, Colorado School of Mines, and M. Onur, SPE, Universiti Teknologi Petronas

†

Copyright 2012, Society of Petroleum Engineers This paper was prepared for presentation at the 2012 SPE Hydrocarbon, Economics, and Evaluation Symposium held in Calgary, Alberta, Calgary, 24–25 September 2012. This paper was selected for presentation by an SPE program committee following review of information contained in an abstract submitted by the author(s). Contents of the paper have not been reviewed by the Society of Petroleum Engineers and are subject to correction by the author(s). The material does not necessar ily reflect any position of the Society of Petroleum Engineers, its officers, or . Electronic reproduction, distribution, or storage of any part of this paper without the written consent of the Society of Petroleum Engineers is prohibited. Permission to reproduce in print is restricted to an abstract of not more than 300 words; illustrations may not be copied. The abstract must contain conspicuous acknowledgment of SPE copyright.

Abstract In oil and natural gas production projects, many investment and development plans are based on oil and gas reserve estimates. There is a large uncertainty in the calculation of hydrocarbon reserves because the input variables always contain uncertainties to some degree that propagate into reserve estimates. From the view point of a field investment, an accurate assessment of uncertainty in reserves is crucial for making decisions that will create value and/or mitigate loss in value. Therefore, to make good decisions, one must be able to accurately assess and manage the uncertainties and risks. In this study, we present an analytical uncertainty propagation method (AUPM) for modeling of uncertainties on volumetric reserve estimations. Analytical uncertainty propagation equations (AUPEs) are derived based on a Taylor-series expansion around the mean values of the input variables. The AUPEs are general in that correlation among the input variables, if it exists, can also be ed for on the resulting uncertainty. Comparative studies that we have conducted show that the AUPM is as accurate as the Monte Carlo method (MCM). The AUPM provides a fast alternative to Monte Carlo simulation for accurately characterizing uncertainty markers such as variance, P90, P50, and P10. In addition, we present uncertainty percentage coefficient for simulating uncertainty contribution of each parameter and correlated parameter pairs to the total uncertainty in volumetric calculations. We also discuss the problem of probabilistic aggregation of reserves for projects involving more than one reservoir or field. We provide a general analytical formulation for estimating the values of mean, variance, P10, P50 and P90 for aggregated estimates. Probabilistic aggregation requires the knowledge of pair-wise correlation of the fields. In this study, we propose uncertainty sorting method (USM) to determine pair-wise correlation coefficients for multiple resources. The method provides a simple and fast analytical approach based on uncertainty percentage coefficient of individual field parameters. Proposed analytical models can be used as a fast tool eliminating the need for MCM. 1. Introduction Uncertainty is inherent in estimation of any type of resources from underground energy systems. Unfortunately, this is also true regardless of any method used for estimation, e.g., volumetric, decline curve, or reservoir simulation methods because the input variables required for the estimation problem always contain uncertainties to some degree that propagate into estimates. Therefore, to make good decisions, one must be able to accurately assess and manage the uncertainties and risks. In this work, we limit our study to the assessment of uncertainty in oil and gas reserves estimated and probabilistic resource aggregation by the volumetric methods. Volumetric methods are usually used in the early life of oil and gas reservoirs. Estimation of the resources/reserves requires geological, geophysical, and petro-physical data including reservoir area, thickness, porosity, saturations, etc. The values of these input variables have usually large uncertainties associated with them; hence it is very important to propagate these uncertainties on to the estimation of hydrocarbon reserves or resources. Although we use the reserve and resources interchangeably here and throughout this paper, their definitions differ indeed; normally reserve is defined as the economically recoverable part of a resource. From the view point of a field investment, an accurate assessment of uncertainty in recoverable and in place reserves is a crucial task from which to make decisions that will create value and/or mitigate loss in value (risk). In the past, various authors have considered assessing uncertainty in estimation oil and gas resources from the volumetric method by using the Monte Carlo method (MCM). It is no doubt that the MCM is used as a general approach for assessing the uncertainty. In this work, however, we show that there is a simple and fast alternative method – which we refer to it as the †

on sabbatical leave from Istanbul Technical University

2

SPE 162875

analytic uncertainty propagation method (AUPM) – to the MC method for characterizing uncertainty. The validity of the AUPM for accurately characterizing uncertainty results from the fact that the distributions of recoverable or in-place hydrocarbons for a zone, well, or field to be computed from the volumetric method always tend to be log-normal. This result simply follows from the fundamental theorem of statistics and probability  the central limit theorem (CLT) (e.g., see Parzen, 1962). The CLT states that the sum of a sufficiently large number of independent random variables each with finite mean and variance will be approximately normally distributed. As a consequence of this theorem and the functional relationship of the volumetric method which involves a product/quotient of several random variables for computing oil and gas reserves, one should expect that the resulting distribution of the reserves is to be nearly log-normal as the number of input random variables increases (Capen 1996 and 2001). This result is in fact valid no matter what form of uncertainty the input variables assume. Furthermore, we discuss the problem of aggregating individual fields’ resources for oil and gas projects involving many diverse fields, zones, wells or a country’s total resources. It is shown that the arithmetic addition underestimates the P90 that would be obtained by the probabilistic addition, which is the statistically proper method for adding resources or reserves unless all fields considered in aggregation are perfectly correlated with each other. We present a method referred to as the uncertainty sorting method (USM) to determine dependency (correlation coefficient) between two fields-zones which is indispensable component of probabilistic aggregation process. 2. Uncertainty on Volumetric Method Based on the volumetric method, we can estimate recoverable gas in-place (RGIP) or oil in-place (ROIP). For instance, for a single-phase under-saturated gas reservoir, recoverable gas in-place (RGIP) is computed from the well-known following volumetric equation: RGIP 

Ah  n g 1  S wi  RF Bgi

(1)

If we replace Bgi by Boi in Eq. 1, then we obtain the volumetric equation for the recoverable oil in place ROIP. As the equations for calculations of RGIP and ROIP are of the same form, we will demonstrate the methods proposed in this work by considering only the RGIP given by Eq. 1. Although it is important how to choose the variables and their ranges for accurately assessing uncertainty in RGIP (Eq. 1), our purpose here is not to get into a detailed discussion of appropriate sources of data and how one could determine the appropriate values of the variables and the associated uncertainties. Assessment of uncertainties in the input variables itself is a notoriously difficult problem, because, to our knowledge, there is no standard rule for characterizing uncertainty in the input variables. Uncertainty in the recoverable reserves results from our lack of knowledge in most of the input variables in Eq. 1. Quantification of uncertainty is inevitably subjective because knowledge about the input variables is dependent on available data and personal experience of the interpreter (Capen 1976, Welsh et al. 2007). As well stated by Welsh et al. (2007), it is quite possible for two people to have different probability estimates for the same input variable, based on their differing knowledge. Thus, there is no single “correct” probability distribution of input variables, unless all people have identical experience and information, and process it in the same way (Welsh et al. 2007). As we will show later based on the CLT, the resulting distributions of RGIP tend to be log-normal regardless of the types of probability distributions chosen for the input variables. It is very important that some input variables involved in Eq. 1 can be statistically correlated. For example, A , the reservoir area, is expected to be dependent on the value of h, reservoir gross thickness. We would expect that increasing A increases h , which indicates that these two variables are positively correlated. Hence, this indicates that A and h may not be treated as two independent variables in Eq. 1. A further discussion of the correlation between the input parameters used in the volumetric method is provided in Holtz (1993) and Hawkins et al. (2002). Our point is that ignoring existing correlations between input variable pairs may lead to an incorrect characterization of uncertainty. If data and available information permit, one should make scatter plots of input variable pairs to identify the correlation between them, if any, and then include these correlations into the uncertainty assessment procedure. In the absence of data analogies can be used. 3. Quantification of Uncertainty In this section, we briefly review some basic equations and methods used for quantification of uncertainty in the estimation of reserves by the volumetric method (Eq. 1). We first consider the Monte Carlo method (MCM) and then the analytical uncertainty propagation method (AUPM). However, before getting into details of these methods, we first note that Eq. 1 for estimating RGIP involves a product of seven random variables; A , h ,  ,  n g  , 1  Swi  , RF and a quotient of one variable; B . Taking natural logarithm of RGIP given by Eq. 1, then we obtain the following equations, respectively:

SPE 162875

3

ln  RGIP   ln A  ln h  ln   ln  n g   ln 1  Swi   ln RF  ln Bgi

(2)

Eq. 2 clearly indicate that ln(RGIP), can be written as a sum of the natural logarithms of all independent random variables. If all these random variables are treated as independent, then it follows from the CLT, discussed previously, that the resulting distribution of ln(RGIP), will tend to be normal. It is important to note that this is true no matter what type of distribution for the input random variables assume. Hence, if ln(RGIP) will tend to be normal, then the corresponding distributions of RGIP will tend to be log-normal. Note that the CLT promises that ln(RGIP) will be normal if all the random variables are independent. However, our simulations shows that the resulting distributions of ln(RGIP) still tend to be normal even if the input variables are treated as dependent. Similar observations have also been reported by Capen (2001). Finally, we would like to state our definitions to be used for characterizing uncertainty in RGIP. For this characterization, we will refer to P90 as “proved or 1P”, P50 as “proved +probable or 2P”, and P10 as “proved +probable possible or 3P”, where P90, P50 and P10 correspond to 10th, 50th and 90th percentiles of the cumulative distribution function, respectively, for RGIP. Here, we comply with the definition used by SPE as most of the SPE papers use these definitions, though these definitions are not strictly consistent with the definitions used in the statistics and probability (Capen 1996 and 2001). Based on the CLT theorem, the uncertainty on the resource RGIP of a single field-reservoir or zone will tend to be a lognormal distribution. So, lnRGIP follows a normal distribution characterized by its mean ( ln RGIP ) and variance (  ln2 RGIP ), 2 while RGIP is a log-normal by its mean ( RGIP ) and variance (  RGIP ). The two sets of parameters are related by the equations (Capen 2001):   2 (3) RGIP  exp  ln RGIP  ln RGIP  2  

and 2 2  RGIP  RGIP  ln exp  ln2 RGIP   1 2

1

(4)

We can also derive the following equations from Eqs. 3 and 4: 1 2



ln RGIP  ln RGIP  ln 1  

2   RGIP  2 RGIP 

(5)

and 

 ln2 RGIP  ln 1  

2   RGIP  2 RGIP 

(6)

The 10th, 50th and 90th percentiles of RGIP are computed from the following equations (Capen 2001):



P90  exp ln RGIP  1.28  ln2 RGIP



P50  exp  ln RGIP 

and



P10  exp ln RGIP  1.28  ln2 RGIP

(7) (8)



(9)

In following subsections, we review some basic equations and methods used for quantification of uncertainty (i.e., computing mean, variance, P90, P50, and P10 given by Eqs. 3-9) in RGIP (or lnRGIP) by the volumetric method (Eq. 1). We first consider the Monte Carlo method (MCM) and then the analytical uncertainty propagation method (AUPM). 3.1. Monte Carlo Method (MCM) The MCM relies on a specified probability distribution of each of the input variables and generates an estimate of the overall uncertainty in the prediction due to all uncertainties in the variables (Kalos and Withlock 2008). As it does not require a linearization of the function and a continuity of the random variables, it is a more general approach for characterizing the uncertainty for any given nonlinear random function f. In our case, f represents RGIP given by Eq. 1. In the applications to be given, we perform Monte Carlo simulations by using @RISKTM spreadsheet-based software (@RISKTM 2004).

4

SPE 162875

3.2. Analytical Uncertainty Propagation Method (AUPM) The AUPM provides a simple approach for estimating the variance of a function defined by several random variables particularly so, of a function defined by products and quotients of random variables, whether they are independent or correlated. The analytical uncertainty propagation equation is based on a Taylor series approximation of the function around the mean values of the variables up to its first derivatives with respect to each of the input variables. As a consequence of this approximation, the uncertainty propagation equation provides a linearization of the function in of its input random variables (Barlow 1989; Coleman and Steele 1999; Onur et al. 2010) The method does not assume a specific type of distribution for the input variables and all needed to use the AUPM are the statistical properties of the distribution of each random variable; specifically the mean, variance (or std. dev.), and the covariance (or correlation coefficient) among variable pairs if the random variables are correlated. Before we present the derivation of the AUPM, it is worth noting that the AUPM provides an exact result for the mean and variance of a random function f if f is linear with respect to the input random variables. Otherwise, i.e., if f is nonlinear, then the AUPM provides only approximate estimates of the mean and variance of f. The approximation gets better if the nonlinear f can be well approximated by a linear function near the means of the input random variables. A few different approaches may be considered when the AUPM is used to derive approximations for the mean and variance of RGIP. For example, we may directly apply the AUPM to the RGIP given by Eq. 1 as a function of the input variables (Approach 1) or may directly apply the AUPM to the ln(RGIP) function given by Eq. 2 by treating it as a function of natural logarithms of the input variables (Approach 2) or treat RGIP given by Eq. 1 as a function of the natural logarithms of the input variables (Approach 3). We have found that Approaches 1 and 2 provide more accurate estimates of the mean, variance, P90, P50, and P10 of the lnRGIP than does Approach 3 In this paper, we present the equations for the AUPM based on only Approach 1 and Approach 2. 3.2.1. AUPM for RGIP (Approach 1). Let’s consider a random function f of M variables, Xi, i=1,2,…,M, i.e., f = f ( X1, X2 ,..., XM ) . Then, expanding f around the mean (or true) values of Xis (denoted by  X , i=1,2,...,M) by using Taylor series up i

to first derivatives, we obtain: M  f  f  X 1 , X 2 ,..., X M   f  X1 ,  X 2 ,...,  X M   X i   X i   i 1  X i  X i  Xi ,i 1,..., M









(10)

It can be shown that the mean  f and variance  2f of f are approximated, respectively, by:  f  f   X ,  X ,...,  X 1

2

M



(11)

and M

M 1 M

i 1

i 1 j i 1

 2f  i2 X2 i  2  i j cov  X i , X j 

(12)

where cov(Xi,Xj) represents the covariance between the variable pairs Xi and Xj, and if we use the relation between covariance and correlation coefficient, Xi ,X j  X i , X j , then we can express Eq. 12 in of the correlation coefficient as: M 1 M

M

 2f  i2 X2  2  i j  X , X  X2  X2 i 1

i

i 1 j i 1

i

j

i

j

(13)

In Eqs. 12 and 13,  i is the derivative of f with respect to the variable Xi, i.e.,  f    X i  X i  Xi ,i 1,...,M

i  

(14)

Note that θi represents the sensitivity of f to the variable Xi evaluated at the mean values of all the variables. It can be noticed (from Eqs. 12 and 13) that the uncertainty propagation on to f is determined not only by the variances of the variables and correlation among them, but also the sensitivity of f to each variable in the volumetric method for RGIP. As mentioned previously, for the problem of interest, f in Eqs. 10-14 represents RGIP given by Eq. 1. The sensitivities of RGIP (i.e., θis ) required in Eqs. 12 and 13 can be obtained by analytical differentiation of Eq. 1 with respect to the input variables in Eq. 1. These sensitivities are presented in Table 1.

SPE 162875

5

TABLE 1—SENSITIVITY OF RGIP WITH RESPECT TO A GIVEN INPUT VARIABLE Xi in EQ. 1 Sensitivity

Variable Xi

i  RGIP X i

h   n g  1 S   RF wi

A

B

gi

 A   n g  1 S   RF wi

h

B

gi

 A h  n g  1 S   RF



wi

B

gi

 A h  1 S   RF wi

(n/g)

B

gi

 A h   n g   RF

(1-Swi)

B

gi

 A h   n g  1 S

RF

wi



B

gi

h   n g  1 S   RF wi

Bgi

B

gi

2 Once the mean  RGIP  and variance  RGIP  of RGIP are computed by the use of Eqs. 11 and 13, we then use these

values in Eqs. 5 and 6 to compute the mean  ln RGIP  and variance  ln2 RGIP  of lnRGIP, which is normal based on the CLT theorem. The other uncertainty markers such as P90, P50, and P10 can be computed from Eqs. 7-9.

3.2.2. AUPM for lnRGIP (Approach 2). The second approach is based on the Taylor series expansion of lnRGIP around the mean values of natural log of the input variables; i.e. ln X i s . For this case, the AUP equations are given by Eqs. 10-14 with f replaced by lnf, Xis by lnXis, and ln X i s . So, in this approach, the mean and variance of lnRGIP are computed from: ln RGIP  ln RGIP  ln X , ln X ,..., ln X 1

2

M



(15)

and M 1 M

M

 ln2 RGIP  i2 ln2 X  2  i j ln X ,ln X  ln2 X  ln2 X i 1

i

i 1 j i 1

i

j

i

j

(16)

where the sensitivities θis in Eq. 16 are given by:   ln RGIP     ln X i  X 

i  

i

(17) X i ,i 1,..., M

and are tabulated in Table 2. Once the mean  ln RGIP  and variance  ln2 RGIP  of RGIP are computed by the use of Eqs. 15 and 16, we then use these 2 values in Eqs. 3 and 4 to compute the mean  RGIP  and variance  RGIP  of RGIP. The other uncertainty markers such as

P10, P50, and P90 can be computed from Eqs. 7-9 using the values ln RGIP and  ln2 RGIP .

6

SPE 162875

TABLE 2—SENSITIVITY OF lnRGIP (EQ. 2) WITH RESPECT TO NATURAL LOGARITHM OF A GIVEN INPUT VARIABLE Xi in EQ. 2 Sensitivity

Variable Xi

i   ln RGIP  ln X i

A

1

h

1



1

(n/g)

1

(1-Swi)

1

RF

1

Bgi

-1

Finally, a few remarks are in order for the AUP equations based on the Approach 1 or 2. Our numerical results indicate that Approach 2 provides a slightly better estimate of the variance than does Approach 1. However, unlike the Approach 1, the Approach 2 requires us to work with the means and variances of the natural-log of the input model variables, i.e., ln X i and  ln2 X i . In the correlated case, we will also need to convert to correlation coefficient between two pairs, say ln X i ,ln X j to

 X , X , but our results show that ln X ,ln X   X , X for all practical purposes. If the distribution of Xi is chosen as a log-normal i

j

i

j

i

j

with mean  X i and variance  X2 i , then lnXi is normal with the mean ln X i and variance  ln2 X i which can be simply computed from: 2 1  X ln X  ln  X  ln 1  2 2   X i

i

i

i

   

(18)

and



2 ln X i

  X2 i  ln 1  2  X i 

   

(19)

If the chosen distribution for the input variable Xi is not log-normal, then we can use descriptive statistics on the available data to compute ln X i and  ln2 X i . If such exhaustive data are not available, then we may generate samples from a known distribution and use descriptive statistics on these samples to compute ln X i and  ln2 X i . 3.2.3 Sensitivities of Parameters. Here we present uncertainty percentage coefficient, UPC, for uncertainty contribution of each parameter and/or correlated parameter pairs to the total uncertainty in RGIP. UPCi is a measure of the uncertainty contribution of individual parameter X i to the total uncertainty:   ln X i  UPCi  i2     ln RGIP 

2

(20)

Also UPCi,j is a measure of the uncertainty contribution of correlated parameters X i and X j :   ln X i    ln X j UPCi , j   i    j   ln RGIP    ln RGIP 2

2

  ln X i ,ln X j 

(21)

The sensitivities  e.g. i  are obtained by analytical differentiation of Eq. 2 with respect to the natural logarithm of the input variables in Eq. 2. These sensitivities are presented in Table 2. Based on Eq. 20 and 21, it can be written as following:

SPE 162875

7 M 1 M

M

UPC    UPC i

i 1

1

i, j

i 1 j i 1

(22)

As a matter of fact, UPC gives almost the same values of squared correlation coefficient obtained from MC simulations. 3.3. Example Applications In this section, we consider some example applications comparing the results obtained from MC and AUP approaches for predicting the uncertainty in reserve estimates for a single field. 3.3.1. Example Application 1. The first example application pertains to a case where all input variables are independent in RGIP given by Eq. 1. For the purpose of this example, we follow two steps: Step 1 assumes that the distribution of each input variable in Eq. 1 can be characterized by a triangular distribution, whereas Step 2 assumes that the distribution of each input variable in Eq. 1 can be characterized by a log-normal distribution. For Step 1, the minimum, most likely (mode), and maximum values of the input variables are given in Table 3. The values of mean and variance given in Table 3 were computed from the well-known formulas for a triangular distribution:

X 

MinX i  MaxX i  ModeX i

and

 Min    Max    Mode   2



2 Xi

(23)

3

i

2

Xi

Xi

18

Xi

2



MinX i  MaxX i  MinX i  ModeX i  MaxX i  ModeX i 18

(24)

For Step 2, we assume that the distribution of each input variable X i is log-normal with the mean and variances computed from Eqs. 23 and 24. Table 4 presents the values of the mean and variances for ln X i computed by using the values of mean

   and variance   for 2 Xi

Xi

X i given in Eqs. 18 and 19.

Table 5 compares the estimates of mean, variance, P90, P50, and P10 computed for RGIP from the MCM and AUPM based on Approaches 1 and 2 as described earlier.. As can be seen from Table 5, the computed values of means and uncertainty markers from the MCM and AUPM for RGIP agree well. We also notice that the values of means and variances computed from the AUPM based on Approach 2 better agrees with corresponding ones from the MCM. This is not surprising though, and is an expected result because as mentioned previously; the AUPM provides a linear approximation to a nonlinear random function around the mean values of the input variables and can provide exact results for the mean and variance for a function f if that function is linear with the input random variables. In our case, RGIP (Eq. 1) are in fact nonlinear functions of their input variables. On the other hand, lnRGIP (Eq. 2) are almost linear (or weakly nonlinear) functions of the input variables. Nonetheless, the AUPM (based on either Approach 1 or 2) provides estimate of mean and variance for RGIP agrees very well with those computed from the MCM.

TABLE 3—DISTRIBUTION PARAMETERS OF THE INPUT VARIABLES; TRIANGULAR DISTRIBUTIONS; STEP 1 Variable X i

Min

Max 6

0.5

0.5

6.67×10-3

0.95

0.9

0.717

2.18×10-2

0.5

0.9

0.85

0.750

7.92×10-3

7.317×10-3

7.386×10-3

7.318×10-3

7.340×10-3

2.58×10-10



5.50×10 122.50 0.15

n g 

0.3

0.7

1  Swi 

0.3

h

RF Bg

6

Variance,  X2 6.83×10-11 447.37 1.67×10-3

6.50×10 194 0.25

6

i

4.68×10 136.60 0.150

2.50×10 93.3 0.05

A

Mean,  X

Mode 6

i

8

SPE 162875

TABLE 4—DISTRIBUTION PARAMETERS OF THE INPUT VARIABLES; LOG-NORMAL DISTRIBUTIONS; STEP 2 Variable X i

Mean,  X

A

6.8310-11 447.37



4.68106 136.60 0.150

n g 

0.5

h

1  Swi  RF Bg

Variance,  X2

i

i

Mean, ln X

Variance,  ln2 X

i

1.6710-3

15.344 4.905 -1.933

0.031 0.024 0.071

6.6710-3

-0.706

0.026

0.717

-2

2.1810

-0.354

0.042

0.750

7.9210-3

-0.295

0.014

-10

-4.914

0.000

-3

7.34010

2.5810

i

TABLE 5—COMPARISON OF THE STATISTICAL PARAMETERS COMPUTED FROM THE AUP AND MC METHODS AUPM

Values

MCM

Approach 1

Approach 2

P90

1.82109

1.77109

1.78109

P50

9

9

3.15109

9

5.71109

9

3.52109

18

2.921018

3.1910

3.1710

9

P10

5.5910

5.6710

9

Mean

3.5110

3.5110

18

Variance

2.6210

2.8510

As stated earlier, the AUPM also provides sensitivity analysis for individual parameters. This sensitivity analysis is useful to indentify which of the parameters contribute more to the total uncertainty in the estimated reserves from the volumetric method. For instance, Table 6 presents the results of such an sensitivity analysis of Example 1 obtained from both the AUPM and MCM methods in comparison. As can be seen from Table 6, the sensitivity values (UPC for AUPM and correlation coefficient) computed from the AUPM and MCM methods agree quite well. A comparison of the values of UPC computed for the parameters indicate that the porosity is the main source of uncertainty, contributing 34% to the variance of RGIP or the total uncertainty in RGIP. TABLE 6—SENSITIVITY ANALYSIS OF INDIVIDUAL PARAMETERS; UNCORRELATED CASE Variable X i

AUPM, UPCi (fraction)

A 

0.15 0.11 0.34

MCM, Squared Correlation Coefficient 0.13 0.10 0.31

n g 

0.13

0.12

1  Swi 

0.20

0.18

0.07

0.06

0.00

0.00

h

RF Bg

3.3.2 Example Application 2. Our next example application pertains to a case where some of the input variables in Eq. 1 are correlated. As mentioned previously, it is possible that various input variables can be correlated with each other. We refer the readers to the work of Holtz (1993) and Hawkins et al. (2002) for further details regarding correlation relationships of input parameters in Eq. 1. For this investigation, we use the same input data given in Table 3, but assume correlation between the five correlated pairs and the correlation coefficients are given in Table 7. Table 8 compares the value of mean, variance and uncertainty markers computed from the MCM and AUPM (Approach 2) methods for the correlated case. As noted before, this approach requires that we work with the correlation coefficient between the pairs in of the natural-log of input random variables; i.e., ln X i ,ln X j Our results indicate ln Xi ,ln X j   X i , X j . That is; one can use the correlation coefficients based  X i , X j when using the AUPM method based on Approach 2. Table 8 compares the means and uncertainty markers obtained from MCM and AUPM based on Approach 2 for RGIP. Again, there is a very good agreement in the uncertainty markers computed from both methods.

SPE 162875

9

TABLE 7— CORRELATED VARIABLE PAIRS AND THEIR CORRELATION COEFFICIENTS Correlated Variable Pairs

Correlation coefficient X ,X

X ,X  i

j

i

 h, RF   A,    A, h 

j

0.3 0.2 0.3

 ,1  Swi   h, n g 

-0.1 0.2

TABLE 8—COMPARISON OF THE STATISTICAL PARAMETERS FROM THE MC AND AUP METHODS FOR THE CORRELATED CASE Values P90 P50 P10 Mean

AUPM 1.66109 3.17109 6.03109 3.59109

MC 1.66109 3.19109 5.96109 3.58109

Variance

3.711018

3.601018

In the correlated case the variance of RGIP is increased significantly (about 25%) and also P90, P50 and P10 is changed slightly (about 5-10%) compared to the corresponding results for the uncorrelated case (Table 8). In this case, AUPM has an advantage that it can give sensitivities both for individual parameters and correlated parameter pairs. Table 9 shows the individual sensitivity parameters for Example Application 2 obtained by the AUPM and MCM methods. Table 10 presents the values of UPC computed for correlated parameter pairs (see the second term in the left-hand side of Eq. 22). TABLE 9—SENSITIVITY ANALYSIS OF INDIVIDUAL PARAMETERS; CORRELATED CASE

A

0.121

MCM, Squared Correlation Coefficient 0.116

h

0.094

0.086



0.283

0.261

n g 

0.104

0.095

1  Swi 

0.165

0.154

RF

0.055

0.050

Bg

0.000

0.000

Variable X i

AUPM, UPCi (fraction)

TABLE 10—SENSITIVITY ANALYSIS OF CORRELATED VARIABLE PAIRS; CORRELATED CASE Correlated Variable Pairs  X i , X j 

 h, RF   A,    A, h   ,1  Swi   h, n g 

AUPM, UPCi,j (fraction) 0.0432 0.0741 0.0640 -0.0432 0.0395

As shown in Table 10, the correlated pairs may have a significant impact (about 18% in our case) on the total uncertainty, i.e., variance of the RGIP. Our results show that correlation among variables if they exist and that data available permits one to identify correlation among variables, should be ed for accurate characterization of uncertainty in the RGIP. Sensitivity analysis of individual input parameters and parameter pairs actually guide us to which parameter we should focus on to characterize much better to reduce uncertainty in the RGIP.

10

SPE 162875

3.4 Dependence on Distributions The MC simulations given in the previous sections are based on an assumption that distributions of all the input variables in Eq. 1 are log-normal. In this section, we will investigate the impact of different distribution types of input variables to the total uncertainty. Therefore several MC simulations are conducted using the same input data given in Table 3 but with different distribution combinations. Table 11 presents six scenarios considered for different distribution combination of individual input parameters. It should be noted here that the mean and/or variance of each input parameter may change from scenario to scenario depending the on distribution assumed for the parameter as we keep the min, max, and mode values the same as given in Table 3. For example, the mean and variance of  (which is assumed to be uniform) in Scenario 4 are 0.15 and 3.3310-3, respectively, whereas the corresponding values for a triangular distribution of  in Scenario 3 are 0.15 and 1.6710-3. Table 12 compares mean, variance, and uncertainty markers of six different scenarios, computed from both the MCM and AUPM (Approach 2) methods. It is should be noted that the values of AUPM given in the Table 12 are computed by using the mean and variances given in Table 4 and it is included in Table 12 as a reference case for comparison. Of course, the mean and/or variance can change for the input parameter as we change the distribution for a given parameter for each scenario considered in Table 11. Although not shown here, the results of our MC simulations show that the RGIP is well approximated by log-normal distribution for each scenario in Table 11, and the mean and the uncertainty markers computed from the AUPM for each scenario in Table 11 agree well with the corresponding ones predicted from the MCM for each scenario in Table 11. Here, our objective is to show that the type of distribution for input parameter (unless its mean and variance are not kept same from distribution to distribution) has an effect on the total uncertainty, as clearly see in Table 12. So,the shapes of distributions for individual input parameters have a remarkable effect particularly on the value of P90. As can be seen from Table 12, distribution shapes may have an impact on calculated uncertainty markers. AUPM almost identical values for uncertainty markers in the case that distributions of all the input variables are log-normal and it gives approximate results when all input variables have normal distribution. As stated in PRMS Guide (2011) normal and lognormal distributions are the most general distributions to represent uncertainty in reserve calculations. It means that AUPM is a powerful tool when decision maker has not enough data for individual parameter distributions.

TABLE 11—SCENARIOS CONSIDERED FOR DIFFERENT DISTRIBUTION TYPES FOR INPUT PARAMETERS Variable X i

Scenario 1

Scenario 2

Scenario 3

Scenario 4

Scenario 5

Scenario 6

A

Log-normal

Normal

Triangle

Uniform

Log-normal

Uniform

h

Log-normal

Normal

Triangle

Uniform

Log-normal

Normal



Log-normal

Normal

Triangle

Uniform

Normal

Normal

n g 

Log-normal

Normal

Triangle

Uniform

Log-normal

Normal

1  Swi 

Log-normal

Normal

Triangle

Uniform

Log-normal

Log-normal

RF

Log-normal

Normal

Triangle

Uniform

Normal

Log-normal

Bg

Log-normal

Normal

Triangle

Uniform

Uniform

Uniform

TABLE 12—MEAN, VARIANCE, AND UNCERTAINTY MARKERS FOR RGIP FOR SIX DIFFERENT SCENARIOS OF TABLE 11 MCM Values

AUPM

P90

Scenario 1

Scenario 2

Scenario 3

Scenario 4

Scenario 5

Scenario 6

1.77×109

1.78×109

1.66×109

1.62×109

9.14×109

1.69×109

1.35×109

P50

3.17×109

3.15×109

3.22×109

3.22×109

2.34×109

3.18×109

2.80×109

P10

5.67×109

5.71×109

5.67×109

5.80×109

5.56×109

5.72×109

5.41×109

 RGIP

3.51×109

3.52×109

3.50×109

3.52×109

2.88×109

3.51×109

3.16×109

2  RGIP

2.85×109

2.92×109

2.75×109

2.88×109

4.21×109

2.89×109

2.91×109

SPE 162875

11

4. Aggregation of Resources Here, we consider the problem of aggregating (or adding) resources of diverse fields, zones or wells. It was shown that the simple arithmetic sum may significantly underestimate the P90 and significantly overestimates the P10, relative to the corresponding ones estimated by probabilistic sum. As is well known [for example see, Capen (1996, 2001), Carter and Morales (1998), and Delfiner and Barrier (2008)], the correct procedure when aggregating the resources of many diverse fields is probabilistic addition whether the fields to be aggregated are independent or not. As to be shown mathematically, the arithmetic sum assumes that all fields considered in aggregation are fully correlated (i.e., pair-wise correlation coefficients for all fields is equal to unity). In other words, probabilistic sum will be equal to the arithmetic sum if all fields are fully correlated. On the other hand, if we assume that all fields are independent and apply the probabilistic addition, then the probabilistic sum overestimates the P10 and underestimates P90 of the all fields used in aggregation. As well stated by Delfiner and Barrier (2008), in reality fields are neither perfectly dependent nor perfectly independent, but instead are correlated. So, although the probabilistic sum is the most general approach that one should use whether the fields are independent or not, however it requires the knowledge of pair-wise correlation coefficients for the fields. 4.1. Aggregation of Means, Variances, P90, P50 and P10 In the following, we provide a general formulation for estimating the values of mean, variance, P90, P50 and P10 for a total of n resources, each following a log-normal distribution characterized by its mean



2 RGIPj





RGIPJ



, j  1, 2,..., n and variance

, j  1, 2,..., n . Note that each RGIP is a lognormal distribution, then lnRGIP is a normal distribution with a mean

2 , j  1, 2,..., n . Now suppose that we are interested in the uncertainty of the equal to RGIPJ , j  1, 2,..., n and a variance  RGIP j

total resources (denoted by RGIPs). It is not difficult to show that the mean (denoted by RGIPs ) and the variance of RGIPs 2 (denoted by  RGIPs ) are given by the following equations, respectively: n

 RGIPs    RGIP j 1

(25)

j

and n1

n

n

2 2 2 2  RGIPs    RGIP  2   RGIP , RGIP  RGIP  RGIP j 1

j

i 1 j i 1

i

j

i

j

(26)

where  RGIPi , RGIPj represents the pair-wise correlation coefficient between resources i and j. A few remarks are in order for Eq. 25 and 26: Eq. 25 indicates that the mean of the sum of the resources is equal to the arithmetic sum of the mean of each resource RGIPj, whether field resources used in aggregation are correlated or not. So, this result indicates that we can add the mean of each field resources to find the mean of the sum of all fields’ resources. Eq. 26 indicates that the variance of the sum of the resources will not be equal to the sum of the variances of individual resources unless all field resources are independent, i.e., RGIPi , RGIPj  0 for all field resource pairs i and j. So if all fields are independent Eq. 26 reduces to: n

2 2  RGIPs    RGIP j 1

(27)

j

On the other hand, if we assume that all pair-wise correlations in Eq. 28 are equal to unity, i.e., RGIPi , RGIPj  1 for all i and j such that i ≠ j, then it is not difficult to show that Eq. 26 reduces to:



2 RGIPs

 n      RGIPj   j 1 

2

(28)

As stated earlier, the sum of a sufficiently large number of independent random variables each with finite mean and variance (e.g. aggregation of reserves) will be approximately normally distributed. In other words, distribution of sum of individual fields converges to normal distribution. Capen (2001) and Delfiner and Barrier (2008) state that convergence is very slow so that log-normal approximation gives better results for aggregation purposes. To investigate the validity of this approximation we conduct a numerical experiment. Our experiment is based on aggregation of different number of random synthetic fields with log-normal RGIP distributions and varying skewness coefficient between 1.0 and 1.8. Monte Carlo simulations are accepted as correct one and relative errors are calculated from values obtained by analytical aggregation

12

SPE 162875

formulas (Eqs. 29- 36) for two cases that the sum of lognormal distributions is log-normal or normal. Figure 1 shows P90 relative errors for both assumptions. As can be seen from the Fig. 1, convergence to normal distribution is very slow because of the skewed shape of the log-normal distribution. Log-normal distribution model is a better approximation for the aggregation of field reserves than normal distribution model for a finite number of fields. This is also valid for P10.

Figure 1. Aggregation of different number of random synthetic fields.

If we will assume the sum of log-normal distribution of reserves is normal then, the P90, P50, and P10 of the RGIPs are given by: 2 P90s  RGIPs  1.28  RGIPs

(29)

P50s  RGIPs

(30)

2 P10s  RGIPs  1.28  RGIPs

(31)

and

On the other hand, if we assume that the sum of n field log-normal reserves is still log-normal, i.e., RGIPs will have a lognormal distribution (as discussed before, in practice with a finite number of fields, this assumption provides better results for P90s, P50S, and P10S), then the RGIPs follows a log-normal distribution with mean equal to  RGIPs (Eq. 25) and variance 2 equal to  RGIPs (Eq. 26). These parameters are related to the mean ln RGIPs and variance  ln2 RGIPs of lnRGIPs by:

1 2



ln RGIPs  ln RGIPs  ln 1  

2   RGIPs  2 RGIPs 

(32)

and 

2  lnRGIPs  ln 1 



2   RGIPs  2 RGIPs 

(33)

SPE 162875

13

The P10, P50, and P90 of the RGIPs, based on the assumption that RGIPs is log-normal, are given by



P90s  exp ln RGIPs  1.28  ln2 RGIPs



(34)

P50s  exp  ln RGIPs 

and

(35)



P10s  exp ln RGIPs  1.28  ln2 RGIPs



(36)

Now, we consider a synthetic example application to our theoretical findings given for aggregation of field resources of RGIP. We consider a project consisting of 2 fields, each having its log-normal distribution of RGIP with its values of mean, variance, P90, P50, and P10 as given in Table 13. TABLE 13—MEAN AND UNCERTAINTY MARKERS FOR TWO FIELDS Values

Field A

 RGIP

Field B

9

8.28×109

18

3.58×10

5.05×1018

P90

3.35×109

5.68×109

P50

5.15×109

8.00×109

P10

9

1.12×1010



5.46×10

2 RGIP

7.93×10

In Tables 14, 15, and 16, we compare the values of mean, variance, P90, P50 ,and P10 for the sum of field’s RGIP for the two fields, obtained from arithmetic and probabilistic sums (based on MC simulations using @RISK, Eqs. 27, 28, 31-33, and Eqs. 37-41). The results given in Table 14 are for the case assuming that the two fields are independent. Table 15 presents the results for the case assuming that the two fields are perfectly correlated. Table 16 presents the results for the case assuming that the two fields are correlated with a correlation coefficient equal to 0.5. The results given in Tables 14-15 and 16 indicate that (i) the arithmetic sum provides almost identical values of P90, P50, and P10 for the sum of 2 fields’ resources if all fields are fully or perfectly correlated with each other, (ii) if all fields are independent or correlated with pair-wise positive correlation coefficients different from unity, then the arithmetic sum will underestimate the value of P90 and P50, but overestimate the value of P10, and (iii) the probabilistic addition assuming normality for the sum of all fields’ resources does not provide as accurate estimates of P10, P50, and P90 as the probabilistic addition based on the assumption of lognormality for RGIPs. TABLE 14—SUM OF THE FIELDS RESOURCES FOR THE CASE WHERE FIELDS RESOURCES ARE INDEPENDENT Methods

P90

P50 9

P10 10

1.31×10

1.92×1010

Arithmetic Sum

9.03×10

Probabilistic Sum (Normal Dist.)

9.98×109

1.37×1010

1.75×1010

Probabilistic Sum (Lognormal Dist.)

1.03×1010

1.34×1010

1.76×1010

Probabilistic Sum (Monte Carlo)

1.03×1010

1.34×1010

1.76×1010

TABLE 15—SUM OF THE FIELDS RESOURCES FOR THE CASE WHERE FIELDS RESOURCES ARE PERFECTLY CORRELATED Methods

P90

P50 9

P10 10

Arithmetic Sum

9.03×10

1.315×10

1.92×1010

Probabilistic Sum (Normal Dist.)

8.44×109

1.37×1010

1.90×1010

Probabilistic Sum (Lognormal Dist.)

9.02×109

1.32×1010

1.92×1010

Probabilistic Sum (Monte Carlo)

9.03×109

1.31×1010

1.92×1010

14

SPE 162875

TABLE 16—SUM OF THE FIELDS RESOURCES FOR THE CASE WHERE FIELDS RESOURCES ARE CORRELATED WITH A CORRELATION COEFFICIENT OF 0.5 Methods

P90

P50

P10

9.03×10

9

Probabilistic Sum (Normal Dist.)

9.15×10

9

Probabilistic Sum (Lognormal Dist.)

9.57×109

1.33×1010

1.85×1010

Probabilistic Sum (Monte Carlo)

9.56×109

1.33×1010

1.85×1010

Arithmetic Sum

10

1.92×1010

10

1.37×10

1.83×1010

1.31×10

In summary, we can state that the arithmetic addition assumes that all fields resources considered in aggregation are perfectly correlated and provides a “pessimistic” estimate of P90. On the other hand, the probabilistic addition based on the assumption that all fields considered in aggregation are independent provides an “optimistic” estimate of P90. In reality, the correct P90 value should be between the P90 values estimated from the arithmetic sum which always assumes that all fields are perfectly correlated and the probabilistic addition based on the assumption that all fields are independent because fields are neither perfectly dependent nor perfectly independent, but instead are correlated. 4.2. Estimation of Pair-wise Correlation Coefficient Estimation of pair-wise correlation coefficients of the resources involved in aggregation process is very critical. Qualified persons-experts interested in the subject try to make the best prediction based on personal experience and judgment. This is the reason why different people working on the same project may give different dependency estimations. Some authors have proposed simplified and pragmatic approaches for aggregating field resources in a more subjective and systematic way to overcome this problem [Carter and Morales (1998), van Elk et al. (2000), and Delfiner and Barrier (2008)]. The method proposed by van Elk et al. (2000) is based on single factor sensitivity analysis of reservoir parameters and tornado diagrams. They use experimental design to determine uncertain reservoir parameters and then define correlation coefficients between those parameters for different resources to estimate overall dependency between two fields. Our approach – uncertainty sorting method (USM) - is based on the similar idea yet we only consider individual parameters of Eq.1 and whole process is completely analytical and straight forward. Instead of tornado diagrams, we consider UPC diagrams obtained from Eq. 20 to rank uncertainty contribution of individual parameters to total uncertainty. Here we present an equation based on Taylor expansion of Eq.1 to calculate correlation coefficient between two resources: n

 A, B 

n

  i 1 j 1

  n

i 1

A Xi



A Xi



2

n 1

n

 2     i 1 j i 1

A Xi

A Xj

A Xi ,X j

A Xi

 XB  X

  A Xi

j

A B i ,Xj

 XA  XB

 

i

n

A Xj

i 1

B Xi



B Xi

j



2

n 1

(37)

n

 2     i 1 j i 1

B Xi

B Xj

B Xi ,X j

  B Xi

B Xj

Here, A and B represents two different resources let’s say two different fields and sensitivity values of individual parameters,  , can be calculated from Table 1. Xi and Xj represent individual parameters of Eq.1 (e.g. Xi=2 represents h etc., n=7 for Eq. 1.) Note that Eq. 37 can be applied only for the equations which involve product/quotient of several random variables such as Eq. 1. USM allows including interior (different reservoir parameters for the same field) and exterior (same parameters of different fields) correlations. Different reservoir parameters for different fields (e.g. h for Field A and n/g for Field B) can also be included but generally those correlations are not taken into . Interior and exterior correlations are required to be determined by experts. Our experience shows that parameters that have UPC under the 5% for both fields do not have a significant contribution to correlation coefficient between resources. Thus one can eliminate the parameters from correlation coefficient estimation process (USM) that have UPC under the 5% for both fields. 4.2.1. Example Application. Field A and Field B are two synthetic fields and their UPC values can be seen from Figure 2. Here Field A and Field B have exterior correlations between the parameters A, h,  and (1-Swi). Besides, Field A has an interior correlation between the parameters: A and h; Field B also has an interior correlation between parameters: A and RF. As can be seen from Figure 2, (1-Swi) has very low UPC values for both fields. Thus we will ignore the correlation of (1-Swi) between the fields. Assigned correlation coefficients and required parameters are given in Table 17 and Table 18, respectively. From Eq. 37, correlation coefficient between the fields can be calculated as  A,B  0.386 . It is worth to note that the parameter with higher UPC has higher effect on correlation coefficient (A for this example) and also positive interior correlation coefficients between the field’s parameters itself reduces the correlation between the fields and negative correlations increases correlation coefficient between the fields.

SPE 162875

15

Figure 2. UPC values of Fields and Correlated Parameter Pairs.

TABLE 17—CORRELATION COEFFICIENTS OF TWO FIELDS Correlated Parameters AA  AB

Correlation Coefficient

hA  hB

0.4

A  B

0.6

AA  hA

0.3

AB  RFB

0.2

0.5

TABLE 18- VALUES REQUIRED FOR EQ.37 - EXAMPLE APPLICATION Field A Xi

A h



Field B



A



A

B

B

1.22×106

9.07×102

1.16×106

1.21×103

1

7

1

5.29×107

-2

1.65×10

5.29×1010

2.12×10

3.99×10

1.42×10



2.04×10

3.64×10

1  Swi 

5.40×10-2

6.42×109

5.40×10-2

9.75×109

-2

9

-2

-2

10

RF Bg

8.25×10

7.61×10

8.69×10

1.17×1010

4.41×10-5

-7.43×1011

4.41×10-2

-1.13×1012

n g 

6.24×10-2

1.02×1010

4.08×10-2

1.38×1010

Conclusions Analytical uncertainty propagation equations (AUPEs) were presented for computing the mean and variance of the recoverable gas or oil-in place. The AUPM method, when combined with the assumption of log-normality for the recoverable gas/oil in-place, provides a fast alternative to the Monte Carlo simulation for accurately characterizing uncertainty markers such as variance, P90, P50, and P10. The derived AUPEs are quite general in that it can for correlation among the input variables used in the volumetric equation. We discussed the problem of probabilistic aggregation of resources or reserves. We provide a general analytical formulation for estimating the values of mean, variance, P90, P50 and P10 for aggregated estimates. We proposed a method so-called the uncertainty sorting method to determine pair-wise correlation coefficients for multiple resources. This method

16

SPE 162875

provides a simple and fast analytical approach based on uncertainty percentage coefficient of individual field parameters it is general for the equations which involve only product/quotient of several random variables. Both AUPM and USM can be used as a fast tool eliminating the need for MCM. Nomenclature A

= Area, m2

AUPM

= Analytical uncertainty propagation method

Bgi

= Initial formation volume factor, m3/ sm3

cov(X, Y)

= Covariance of X and Y

f

= A random function

h

= Reservoir gross thickness

M

= Number of variables

MCM

= Monte Carlo method

n

= Number of total resources

(n/g)

= Net to gross thickness, fraction

PXX

= Value which is likely to be exceeded XX% probability

PXXs

= PXX value of total resources

RGIP

= Recoverable gas in place, m3

RGIPs

= Recoverable gas in place of total resources, m3

RF

= Recovery factor, fraction

Swi

= Initial water saturation, fraction

UPC

= Uncertainty percentage coefficient

X

= Variable

i

= Sensitivity of variable Xi

 x,y

= Correlation coefficient between x and y

x

= Standard deviation of x

x

= Mean of x



= Porosity, fraction

References Barlow, R.J.1989. Statistics: A Guide to the Use of Statistical Methods in the Physical Sciences, John Wiley & Sons. Capen, E.C. 1996. A Consistent Probabilistic Definition of Reserves. SPE Reservoir Engineering 11(1): 23-28.SPE-25830-PA. http://dx.doi.org/10.2118/25830-PA. Capen, E.C. 2001. Probabilistic Reserves! Here at Last? SPE Reservoir Evaluation & Engineering 4(5): 387-394. SPE-73828-PA. http://dx.doi.org/10.2118/73828-PA Capen, E.C. 1976. The Difficulty of Assessing Uncertainty. JPT 28(8): 843-850. SPE-5579-PA. http://dx.doi.org/10.2118/5579-PA Carter, P.J. and Morales, E. 1998. Probabilistic Addition of Gas Reserves Within a Major Gas Project. Paper SPE 50113, presented SPE Asia Pacific Oil & Gas Conference and Exhibition, Perth, Australia, 12-14 October. http://dx.doi.org/10.2118/50113-MS Coleman, H.W. and Steele, W.G. 1999. Experimentation and Uncertainty Analysis for Engineers, second edition, Wiley& Sons. Delfiner, P. and Barrier, R. 2008 . Partial Probabilistic Addition: A Practical Approach for Aggregating Gas Resources. SPE Reservoir Evaluation & Engineering 11(2): 379-385. SPE-90129-PA. http://dx.doi.org/10.2118/90129-PA Hawkins, J.T., Coopersmith E. M. and Cunningham P.C. 2002. Improving Stochastic Evaluations Using Objective Data Analysis and Expert Interviewing Techniques. Paper SPE 77421 presented at SPE Annual Technical Conference and Exhibition, San Antonio, Texas, U.S.A., 29 September-2 October. http://dx.doi.org/10.2118/77421-MS Holtz M.H. 1993, Estimating Oil Reserve Variability By Combining Geologic and Engineering Parameters, Paper SPE 25827 presented at the SPE Hydrocarbon Economics and Evaluation Symposium, Dallas, Texas, U.S.A., 29-30 March. http://dx.doi.org/10.2118/25827MS Kalos, M.V. and Whitlock P.A. 2008. Monte Carlo Methods, Willey-Blackwell.

SPE 162875

17

Onur, M., Sarak, H. and Türeyen, Ö.Đ. 2010. Probabilistic Resource Estimation Of Stored and Recoverable Thermal Energy For Geothermal Systems By Volumetrics Methods, proceedings World Geothermal Congress, Bali, Indonesia, 25-29 April Parzen, E. 1962. Modern Probability Theory and Its Applications, John Wiley & Sons. @RISK, version 4.5.5. 2004. NY, US: Palisade Corporation. SPE/AAPG/WPC/SPEE/SEG Guidelines for Application of Petroleum Resources Management System, 2011 http://www.spe.org/industry/docs/PRMS_Guidelines_Nov2011.pdf (ed 10 October 2011) van Elk, J. F., Vijayan, K., and Gupta, R. 2000. Probabilistic Addition of Reserves. Paper SPE 64454 presented at the SPE Asia Pacific Oil & Gas Conference and Exhibition, Brisbane, Australia, 16-18 October. http://dx.doi.org/10.2118/64454-MS van Elk J.F., Gupta, R. and Wann D. 2008. Probabilistic Aggregation of Oil and Gas Field Resource Estimates and Project Portfolio Analysis. Paper SPE 116395 presented at the SPE Asia Pacific Oil & Gas Conference and Exhibition, Perth, Australia, 20-22 October. http://dx.doi.org/10.2118/116395-MS Welsh, M.B., Begg, S.H., and Bratvold, R.B. 2007. Modeling the Economic Impact of Cognitive Biases on Oil and Gas Decisions. Paper SPE 110765 presented at the SPE Annual Technical Conference and Exhibition, Anaheim, CA, USA, 11-14 November. http://dx.doi.org/10.2118/110765-MS