Improving gas capacity forecast model by mapping real time data to curves

(1)

Improving gas capacity forecast model

by mapping real time data to curves

Marina de Wolff

Groningen, December 8, 2011 Master thesis EORAS

Specialization: Operations Research University of Groningen

Supervisors:

prof. dr. R.H. Teunter, RUG G. de Boer, ORTEC

(2)

Abstract

This thesis is written in the field of Econometrics, Operations Research and Actuarial Studies and will be about well behaviour of a gasfield in the northern part of the Netherlands. The modeling of well behaviour is based on data measured at the surface of a well. With this data a model is formed and used by the NAM to make forecasts about the gasfield. This well modeling is of great importance, since the NAM makes agreements with gas delivery companies. This model is calibrated once a year and is mainly based on well configuration and geological measurements. Lately it is noticed that this model is not a good reflection of the well behaviour, and therefore not reliable enough. At the moment the model is adjusted manually in order to come to the right solution. This is not only time consuming, but might also lead to the wrong corrections. This thesis will focus on two aspects; filtering of the incoming data (not all measured data can be used) and finding a good well behaviour model by using filtered real time data.

For filtering the data two methods are used: the Signal to Noise Ratio and the Mean Absolute Deviation. Both filter techniques lead to similar results. These filters are easy to use and easy to implement, which is a great advantage when comparing it with the currently used method. Next the well modeling method used by the NAM is analyzed. We show that by regressing on real-time filtered data, a much better model fit is attained compared to the current practice of modeling well behaviour, at least in the short term. Moreover, even if we restrict the regres-sion such that some aspects of well modeling have to remain unchanged (in line with geological knowledge), there is still a considerable improvement in fit.

(3)

Introduction

There are different types of fuels which can generate energy. The three main types are

• Biofuels; such as bioethanol and biodiesel.

• Fossil fuels; such as coal, petroleum and natural gas.

• Nuclear fuels; such as uranium-235 and plutonium-239.

The primary source of energy comes from fossil fuels (86.37 % of the total energy consumed in 2008 (www.eia.gov). Fossil fuels are of great importance as they produce a lot of energy. Natural gas is one of the ‘cleanest’ and safest fossil fuels. It produces less by-products then other fossil fuels. This is one of the reasons why it is so popular. When burning 1 m3 _{natural gas it}

will discharge 1.8 kg CO2. Coal for example produces 2.6 kg/m3 CO2.

(6)

Figure 1.1: Natural gas production over the world (in m3₎

The United States and Russia were the largest producers of natural gas (www.enerdata.net) in 2009 and 2010. The Netherlands is in the ninth place (2009). When only looking at Europe the Netherlands is in the top three producers. One of the largest natural gasfields in the world is situated in the Netherlands, the so-called Groningen gasfield. This gasfield has a surface of approximately 900 km2 _{and lies 3000 meter below surface.}

Gas extracting companies make a long-term forecast on how to produce natural gas in the future, to satisfy demand and to develop the gasfield in the best way. This forecast is also used to give a prediction on how long the gasfield can be used. These companies make long-term agreements with gas delivery companies, on where, when an how much gas to deliver. To make such long-term agreements the gas extracting company must have an idea how the gasfield will behave in the near future. There is a need for short-term forecasting as well, since the NAM wants to know what can be expected tomorrow based on the knowledge of today.

Figure 1.2: Well modeling

(7)

and taking into account that the amount of gas in the gasfield decreases. The main focus of this research will be on the behaviour of wells, and how these can be modeled.

The research questions are:

How is the current methodology specified when modeling the behaviour of a well? Most im-portantly, is this method sufficiently accurate? If not, can it be improved or is it necessary to define a new methodology?

(8)

Chapter 2

Natural gas

Natural gas is an example of a mixture of hydrocarbons, which is a combination of carbon (C) and hydrogen (H). There are different types of hydrocarbons, where the simplest form are alkanes. Natural gas mainly consists of methane (CH4), but also other longer hydrocarbons, nitrogen (N)

and carbon dioxide (CO2) are present. When natural gas is delivered to the consumer, it consists

mainly of CH4since all other materials have been removed as much as possible. Natural gas has

no color, taste or smell when it is in it’s pure form. The familiar smell of gas (a thiol, smells like rotten eggs) is added for safety reasons. Natural gas can be found in several forms before it is treated

• Wet gas, contains a considerable amount of longer hydrocarbons. • Dry gas, almost pure methane.

• Sour gas, contains hydrogen sulfide. • Sweet gas, very little hydrogen sulfide.

2.1 Formation of natural gas (in the Netherlands)

(9)

Figure 2.1: Sectional plane of a gas reservoir

Under the surface of the Netherlands several gasfields are present, where the Groningen gasfield is the largest one. There are also gasfields at the Wadden and North Sea. The map of the Dutch gasfields can be found in appendix B.

Natural gas plays a major role in our daily life, but the amount of gas available is finite. It is true that even now new gasfields are discovered, but it is very unlikely that new gas will be created. This makes it really important that optimal decisions are made when producing gas.

2.2 NAM

(10)

gas fields were discovered. Nowadays the main focus is natural gas. NAM is the largest gas producer of the Netherlands and produces 75 % of the gasproduction. About 50 % of the total gasproduction comes from the Groningen gasfield, and the rest from smaller gasfields.

2.2.1 Groningen gasfield

The first drilling of the Groningen gasfield took place in Slochteren on may 29th 1959 at 5.45 in the morning. This drilling started on the land of Kornelis P. Boon, under the management of Joop Boering, and reached a depth of 3000 meter below surface before finding gas. No one was really exited. They did not know that their discovery would lead to finding an enormous gasfield. In 1960 drilling took place near Delfzijl, where gas was found of the same composition and under the same pressure as the gas in Slochteren. These similarities were an indication that the gas came from one reservoir, and was estimated to have a volume of 60 billion m3_{. Two}

years later, in 1962, the Groningen gasfield was taken into use. The volume was extended to 470 billion m3_{, and kept rising. The discovery of such a large gasfield led to a fast increase of gas}

usage in the Netherlands, and therefore the Dutch appliances are tuned for the Groningen gas (table 2.1). Ten years later, 75% of the Dutch households used natural gas as a source of energy. Nowadays almost all Dutch households use natural gas.

Component1 _% _Component _% _Component _%

CH4 81.30 C4H10 0.14 CO2 0.89

C2H6 2.85 C5H12 0.04 O2 0.01

C3H8 0.37 C6H14 0.05 N2 14.35

Table 2.1: Composition of Groningen natural gas

2.2.2 Small fields

Besides the large Groningen gasfield several small gasfields are in use. In 1963 the Dutch gov-ernment sold a lot of natural gas to foreign countries as they thought that nuclear energy would become the biggest future energy source. It turned out this wasn’t the case, and in 1974 the NAM decided to search for other (smaller) gasfields to lower the burden of the Groningen gas-field. There are approximately 200 small fields in use, on land and in the North Sea (a part belongs to the Netherlands). Because of it’s success this policy is still in place. The smaller fields are used at a constant production rate, and the Groningen gasfield is used during peak periods.

(11)

Figure 2.2: Production of natural gas in the Netherlands

The blue area in Figure 2.2 represents the production of the smaller gasfields, whereas the yellow area represents the production of the Groningen gasfield. In some months the demand is so high that gas needs to be extracted from under ground storage areas (UGS’s), meaning that the green areas in the figure are filled up with gas extracted in low demand periods, for example during the summer.

(12)

Chapter 3

The research

The NAM uses GENREM2 _{to make (long-term) forecasts of the predicted capacity of the}

dif-ferent gasfields. This network simulation is used to utilize the gasfields in the best possible way. During the winter of 2010 the production at the Groningen gasfield and some smaller gasfields were high and it was noticed that the forecast did not correspond to the actual values. This mismatch was caused by the so-called THP/Q curve, which reflects the behaviour at the wells as will be explained in what follows. This curve did not correspond well enough with the actual values measured. Due to this mismatch the NAM could not rely on the GENREM forecast to fulfill their obligation to the Gasunie, which is the Dutch gas transportation company. A (high) fine will be received when the NAM cannot deliver the amount of gas nominated to the Gasunie, but the NAM will also miss out on a lot of money when they can deliver more than nominated. The additional gas extracted is a waste product, or is sold at a very low price.

At this moment the THP/Q curve is updated once a year, and the same curve is used throughout the year. It seems preferable to update the curve during the year with real time data, to capture changing behaviour of the well. Instead of updating the curve once a year it would adjust itself when real time data is added to give a better prediction of the relationship. Employees of the NAM already did some adjustments to improve the fit of the THP/Q curve. At the well they performed corrections manually to match the curve to their observations. This was based on ‘trial and error’ and is very time consuming.

Due to thorough research the NAM knows quite well how a well will behave. In general the behaviour is; the lower the THP, the higher the capacity. An important variable is the reser-voir pressure, which is the pressure below surface in the reserreser-voir. When extracting gas from a reservoir the THP needs to be lower than the reservoir pressure, so that gas will flow from the reservoir to the surface. The difference between these two pressures determines the flow rate of the gas, which is the capacity. Lowering the THP will result in a larger difference with the reservoir pressure (assuming that the reservoir pressure did not change), and thus in a higher flow rate (capacity).

Inside NAM there are a different programs used to generate THP/Q curves, all using the same underlying geological model, but using other methods to get the parameters. WellMon for ex-ample first filters data based on the state of the measurements (explained further on), and uses (linear) regression to find the THP/Q curve. An other approach is to use ‘raw’ data and find

(13)

the THP/Q curve based on well configuration and gas and reservoir characteristics.

When fitting a model to real time data it is very important to filter the dataset, since only steady data can be used. If this is not the case the model might give misleading results. The WellMon manual will be described, and used as a guideline to develop a working method which can be implemented in GENREM. The filter technique developed in this thesis will not use ad-ditional information, such as the state of the measurement. With this property it distinguishes itself from the WellMon filter.

First a distinction needs to be made between short-term forecasts (48 hours) and long-term forecasts (> 1 year), since they can require different working methods.

3.1 Short-term forecast

The reservoir pressure declines over time. For the short-term forecast this is not that influential. It is also not of interest how the well behaves outside the (currently) operating area, since the production of tomorrow is similar to the production of today. Therefore it is possible to determine a regression THP/Q curve of the data points (Q,THP) that does not take the reservoir pressure (Pres) into account.

3.2 Long-term forecast

When the NAM wants to make a long-term forecast based on the THP/Q curves, it makes use of reservoir characteristics. The THP/Q curve is corrected due to reservoir decline, gas and fluid characteristics. This might put a restriction on the possible curve fitting methods, as the reservoir pressure must be taken into account. For this reason it is preferred by the NAM to use the current method, but improve it to reduce the workload for the user. This could be done by inserting real time data, which would change the model parameters.

Besides improving the current THP/Q method, it is very usefull to manner in which the long-term predictions are done is correct. The current THP/Q curve is corrected for another reservoir pressure, meaning that the parameters in the model are corrected for a reservoir pressure which would occur in the future. The idea is to fit a THP/Q curve on data which is measured a few years ago, and correct for other reservoir pressures, which are required a few years later.

3.3 WellMon

WellMon is a cluster performance-monitoring program developed in 2005. It determines the well and cluster inflow and outflow parameters. Real time Q and THP can be plotted against the reference curve and gives an overview of the performance in time.

3.3.1 Data filtering

(14)

Step 1: When at least one well has a quality indicator ‘bad’ or ‘uncertain’ all data points (of all wells in that time instant) are removed.

Step 2: Determine the state for each well for each data line. This can be no production (NP), free flow (FF) or compression (CO).

Step 3: Search for time periods that meet the following input parameters; the minimum closed in period, minimum time free flowing, minimum time with compressor and maximum time of production.

This will lead to data where only steady-state situations are used.

3.3.2 AFBC model and THP/Q curve

To model the well behaviour the NAM makes use of the AFBC model, as discussed in Manual Prototype WellMon (Version 1.0.0.). This method finds parameter values for A, F, B and C. All these parameters represent a geological measure. For finding these values data is needed of (Q,THP), but also the reservoir pressure. Normally the reservoir pressure used comes from GENREM. To make sure that this is correct the reservoir is calibrated once a year. After this calibration the THP/Q curve is updated for all wells of that reservoir.

The AFBC model is based on the following equations:

(Cullender-Smith) BHP2 = B × THP2+ C × Q2 (Forcheimer) P2_{res − BHP}2 = A × Q + F × Q2

In the first equation the capacity is based on the difference between the bottom hole pressure and the tubing head pressure, and in the second equation the capacity is based on the difference between the reservoir pressure and the bottom hole pressure. The B parameter is defined by

B = CIBHP

2

CITHP2 where

A Darcy flow coefficient(bar2/(103_m3_/d))

B Static pressure correction coefficient BHP Bottom Hole Pressure (bar)

C Tubing friction coefficient(bar2/(103m3/d)2) CIBHP Closed in bottom hole pressure

CITHP Closed in tubing head pressure

(15)

The Cullender-Smith equation represents the performance of the vertical flow, which is the flow in the tube. The Forcheimer equation represents the performance of the horizontal flow, which is the flow from the reservoir to the well. It is important to note that the parameters A, F, B and C are nonnegative.

When the B parameter is determined (using measurements for CIBHP and CITHP), the other parameters in the AFBC-model can be found by combining the Cullender-Smith and Forcheimer equations to

P2_{res − B × THP}2

Q = A + (F + C) × Q

and performing linear regression (Rice (2007)). The intercept corresponds to A, and the slope to (F + C). When generating the THP/Q curve it does not matter which values F and C individually have. The THP/Q curve is given by

THP =

s

P2_{res − (A × Q) − ((F+C)Q}2) B

This curve is calibrated once a year, and used throughout the year. However, natural gas is not an ideal gas and the AFBC parameters are not constant for other reservoir pressures, and therefore changes when making long-term forecasts3_.

3.3.3 Correcting THP/Q curve manually

During the winter of 2010 the AFBC values were adjusted in such a way that the AFBC-model fitted better to the data, and these corrections were done manually. This is not only very time consuming for the user but also quite subjective as we will explain in this subsection. Let us consider the following example: A = 10, F = 0.0015, B = 1.5 and C = 0.007, and assuming a constant reservoir pressure of 80 bar.

(16)

Figure 3.1: Correcting A, F, B and C manually

(17)

Chapter 4

Data set

The data provided by the NAM consists of tubing head pressures and capacity flows of 8 wells in the north of the Netherlands (due to confidential information the exact location and the names of the wells remain unknown); well W1, up to and including well W8. The AFBC values and the reservoir pressure used that time to generate the THP/Q curve are also available. These 8 wells are part of some smaller gasfields in the Netherlands.

The timespan of the data is approximately one month, march 2011, where one data point repre-sents one minute. The size of a monthly dataset is therefore approximately 42.000 data points. This data set consists of ‘raw’ data, which means that it has not been filtered.

4.1 AFBC values

The NAM provided the Pres and AFBC values of the GENREM simulation during march 2011. The parameter values for A, F, B and C are not based on the technique used in WellMon, but derived from other interpretation techniques and software.

Well A F B C _Pres W1 14.194 0.0012 1.5479 0.00331 63.9 W2 NA NA NA NA NA W3 15.834 0.00721 1.6303 0.01281 83.4 W4 3.1708 0.001 1.5679 0.0017 64.2 W5 28.982 0 1.6208 0.01312 248.7 W6 3.6413 0.014916 1.7795 0.012272 140.0 W7 11.838 0 1.6431 0.00811 310.9 W8 13.98 0.00336 1.5361 0.0093 450.1

Table 4.1: ‘Old’ A, F, B, C and Pres, march 2011

(18)

(19)

Figure 4.1: AFBC models provided by the NAM

It is not difficult to see that this curve is not always a good reflection of the behaviour of all wells. Modifying the parameters could improve the fit. It is also obvious from Figure 4.1 that data filtering is needed. This is necessary for finding a model which reflects the well best, since using the raw data gives misleading conclusions. When the THP is adjusted at the well it needs some time to achieve the corresponding capacity. This unstable situation needs to be filtered out of the dataset. It is hard to say from Figure 4.1 which data points will be in the filtered dataset, since it is unknown which data points occurred when. How to filter the dataset will be described in the next section.

4.2 Data filtering

The data cannot be used in its raw form for the research. It needs to go through a filter first. A filter technique needs to be implemented in the software program before model fitting can be done, and in this research the software program R is used. It is preferred to use a technique in line with WellMon, thus where steady-state data will be in the filtered dataset.

Kim et al. (2008) describe a steady state detection technique, which is based on a moving window standard deviation. The data which is selected as steady state falls between ± 3σ of its average value. A suggestion made, is to use ± 1σ for very restrictive indication of steady-state. Applying this method for the data in this thesis, where the moving window is 60 minutes, almost all data was indicated as stable. Another approach is to fit linear models on subsets of the data as in Pilar Moreno (2004), again using a moving window. If the slope of the line is smaller than a given threshold, the data is indicated as stable. A slope close to zero would mean that the data is stable. However, this does not tell us if the data is volatile or not. As an addition statistical measures, such as the R2_{, could be used to address this volatility. Taiwen et al. (2003) describe}

(20)

indicated as stable by the NAM.

A statistical measure indicating the variability in the data is the coefficient of variation (CV) or signal to noise ratio (SNR), where SNR = CV−1. The main advantage of this measurement is that it is dimensionless, and therefore possible to compare the CV or SNR of datasets with different units. This is what is preferred here, since both THP and Q need to be stable, and have different units. In Aja-Fern´andez (2006) the CV was used as a component of a filter, to detect noise in imaging sytems. In this filtering technique the variation of the noise and image was used to filter images, and increase its quality. In Baldwin et al (2004) variability was examined using moving-window coefficient of variation. The CV was used to measure the change in variability after a sign. This measure of variability is the base of one of the filters in this thesis. To my knowledge the CV or SNR is not used as a filter for data types similar to the dataset in this thesis. WellMon currently uses additional information, such as the state of the well, where the fil-tering techniques developed in this thesis will not use this information. This is a big advantage. These filters will only require raw data, and can therefore be easily implemented in other software programs. Another advantage of the filtering techniques developed in this thesis is that it only takes a couple of minutes to filter the dataset (even seconds in the software program used by ORTEC), and therefore saving a lot of time. The filter will consists of two main steps.

Step 1: Delete the data points where no value or values are available. Step 2: Only include points were the data is in a stable situation. Two filter techniques will be described in the following subsections.

4.2.1 SNR filter

For this filter we use the Signal to Noise Ratio (SNR) to determine which data will be in the final dataset. The SNR, at time t, is calculated as follows

SNR = µˆ ˆ σ

where ˆµ is the mean estimate and ˆσ is the estimated standard deviation of the noise, both in the interval [t − h + 1, t]. A high SNR indicates that the data is very stable, thus gives meaningful information, whereas a low SNR indicates the data being volatile. Subsets with a high mean are allowed to be more scattered. When two subsets have the same SNR and if the first subset has a much lower mean than the second, this implies a higher variance for the second subset. When a well is closed for a longer period (longer than h), the SNR for Q gives ‘Not a Number (NaN)’, this data is not taken into account. It is of course possible to include data for closed wells. A criterion is that the SNR of THP is above the threshold. The filter works as follows;

• Decide if closed in data needs to be included. • Determine the length of the subset (h).

• Set a threshold for the SNR. Data with a SNR above the threshold will end up in the final dataset.

(21)

threshold is set at 100, which works as desired for this dataset. It is always necessary to check if the final dataset is indeed as desired, and if not it might be wise to adjust the length of the subset and/or the threshold. Both the THP and the Q are filtered with this technique. The filtered dataset contains data where both THP and Q have a SNR above the threshold.

This technique gives us the following filtered data for well W1, if the length of the subset is 60 and the threshold is 100 for Q and THP.

Figure 4.2: SNR filter for well W1

This filtered data looks very reasonable for well W1, and might be a good replacement for Well-Mon. The SNR filter technique works not only as desired, but is also time-saving and a lot easier to use.

(22)

4.2.2 MAD filter

The routine of this filter works quite the same as of the SNR filter. The median absolute deviation (MAD), at time t, is calculated as follows;

MADQ = _h1 t X i=t−h+1 |Q_i− ˆ_µQ| MADTHP = _h1 t X i=t−h+1 |THPi− ˆµTHP| ˆ

µQ and ˆµTHP are respectively the estimate of the mean of Q and THP in the interval [t−h+1, t]. The MAD is a more robust filter, where outliers have less influence. Here a low MAD indicates stable data and a high MAD indicates the data being unstable. The filter works as follows

• Determine the length of the subset (h).

• Set a threshold for the MAD. Data with a MAD below the threshold will end up in the final dataset.

The choice of both h and the threshold determine the final dataset. It is again necessary to check if the final dataset is indeed as desired, and if not it might be wise to adjust the length of the subset and/or the threshold. Both the THP and the Q are filtered with this technique. The filtered dataset contains data where both THP and Q have a low MAD.

(23)

Figure 4.3: MAD filter for well W1

This technique takes closed in values into account, since the MADQ is zero and will always be below the threshold. For this dataset these values are likely to be incorrect for the same reasoning described in section 4.2.1. This filter technique might work better for other wells, where the well is closed for a longer period.

(24)

Chapter 5

Results

In this chapter the results of the models fitted will be discussed. The dataset of well W1 will be used in this section. The results of the other wells can be found in appendix C.

In this chapter the AFBC model will be inspected in detail, and can be used for short and long-term predictions. Several linear, nonlinear and nonparametric models were fitted to the data, which can only be used for short-term prediction. These models did not make any im-provement when comparing it with the updated AFBC-model. Therefore the result of these models are not included in the main text, but can be found in appendix C.

5.1 AFBC model

Recall that the AFBC model reads as follows;

THP = s P2_{res − (A × Q) − ((F + C)Q}2) B + e THP2 = P 2 res B − A BQ − F + C B Q 2_{+ e}0

Where e represents an error term. The linear model fitted is

THP2 = β0+ β1Q + β2Q2+ e (5.1)

The model might look nonlinear rather than a linear model due to the squares in the expression, but this is just a matter of notation. The model is linear in the parameters, and is therefore a linear model.

(25)

the following result is obtained

Figure 5.1: Model (5.1)

Model β0 β1 β2 R2

(5.1) 1506.14 -2.09 -1.12E−3 _0.983

Table 5.1: Result model (5.1)

Comparing the old AFBC model with model (5.1) a major improvement has taken place. As model (5.1) is significantly different from the line based on the currently assumed values of Pres, A, F, B and C, it is obvious that one or more of these parameters have to be modified to obtain the results of model (5.1). Note that there is no unique set of values, since there are 4 variables (F+C is seen as one) and 3 equations.

Pres A B F+C

63.9 5.66 2.71 3.03E−3

78.00 8.43 4.04 4.51E−3 48.28 3.23 1.5479 1.73E−3

(26)

Fixing Pres lead to considerable changes in the AFBC parameters when comparing with the ‘old’ AFBC values. These parameter results are not reasonable in practice, since they are based on knowledge of the well and reservoir. It is for example very unlikely that the F and C parameters are incorrect. Also, a high B indicates that a lot of water is present close to the well in the reservoir, which is probably not the case here. Fixing (F+C) gives an even more implausible result. The B factor has even increased more than with a fixed Pres. This is due to the increase of the (F+C) factor, which forces to increase B as well. Fixing the B parameter gives a more plausible result with indeed a lower reservoir pressure.

Next a more restricted model will be fitted, where the B,C and F parameters are fixed and the Pres and A are the varying parameters.

THP2 = β0+ β1Q −

F + C

B Q

2_{+ e} _(5.2)

The following result is obtained

Figure 5.2: Model (5.1) and (5.2)

Model β0 β1 R2

(27)

Pres A B F+C 43.41 0.97 1.5479 4.51E−3

Table 5.4: Fixing F, B and C

Adding the restriction of fixed F, B and C gives again a nice fitting curve, where the statistically seen the fit of model (5.2) will always be equal or lower than the fit of model (5.1). Both models give back different AFBC-Pres combinations. The best AFBC curve is thus the one which is most plausible in practice. It is difficult to draw a conclusion about these curves as a non-reservoir engineer, therefore it is best to generate different AFBC-Pres curves under certain restrictions.

5.2 Coefficient corrections for long-term forecasting

The AFBC-model is used for short-term as well as long-term forecasting. When using the model for long-term forecasting the AFBC parameters are corrected for other gas properties4_{. In this}

section the ‘old’ AFBC model, model (5.1) and (5.2) are fitted on a dataset, and the fit of the corrected models in the subsequent dataset will be examined. Thus if these corrected curves do indeed match with real time data.

The data for this part provided by the NAM is data of well W4. Different years of data were available. The results of 2007-2011 will be discussed, since other datasets gave similar results. The corresponding A, F, B, C and Pres provided by the NAM are

Year A F B C _Pres

2007 3.31 0.001 1.57 0.0017 97.8 2011 3.17 0.001 1.57 0.0017 64.2 Table 5.5: ‘Old’ A, F, B, C and Pres well W4

Note that the A, F, B and C parameters used by the NAM did not change (a lot) over time. First the datasets need to be filtered as described in 4.2.1, all with a threshold of 100. This is of great importance, since the fitted models might otherwise reflect the wrong behaviour. To examine the coefficient corrections done by the NAM model (5.1) and (5.2) will be fitted to the filtered data of 2007 and 2011. As been addressed earlier, model (5.1) has no unique solution for the A, F, B, C and Pres. In this section we will only consider the solution where Pres is fixed. If we would assume that the reservoir pressure might be incorrect, it is hard to believe that the reservoir pressure in the 2011 dataset will be correct. Using model (5.2) for corrections is therefore a bit tricky, since we already assume that the reservoir pressure is incorrect. With the fitted models of 2007 coefficient correction will be done for the reservoir situation in 2011, and examine the fit of the models to this dataset.

When fitting model (5.1) and (5.2) to the filtered data of 2007 and 2011 the following results are obtained

(28)

Figure 5.3: Model (5.1) and (5.2) Year Model β0 β1 β2 R2 2007 (5.1) 5181.27 -1.63 -9.87E−4 0.970 (5.2) 4783.96 -0.52 -1.72E−3 0.968 2011 (5.1) 1852.01 -2.64 0 0.967 (5.2) 1812.46 -1.69 -1.72E−3 NA Table 5.6: Result model (5.1) and (5.2)

It can be seen that the models fitted to the data are a major improvement compared to the old model. It is obvious that one or more of the AFBC-Pres parameters have to be modified to obtain a nice fit. It is therefore not likely that the AFBC values stay constant over the years, which is the case for the parameter values provided by the NAM.

Year Model _Pres A B F+C

2007 ‘old’ 97.8 3.31 1.57 0.0027 (5.1) 97.8 3.00 1.85 1.82E−3 (5.2) 86.66 0.82 1.57 0.0027 2011 ‘old’ 64.2 3.17 1.57 0.0027 (5.1) 64.2 5.87 2.23 0 (5.2) 53.34 2.65 1.57 0.0027 Table 5.7: Fixing parameter values

For simplicity we will assume that the corrected AFBC factors will be very close to the old ones, and therefore only a different reservoir pressure is taken into account.

(29)

Figure 5.4: Model (5.1) and (5.2)

Model (5.1) follows the pattern of 2011 the best, but the curve lies above the data. The old AFBC model and model (5.2) are not closer to the data, and show the wrong shape. To draw right conclusions about what happens the results of model (5.1) are combined in one table

data curve β0 β1 β2

2007 2007 5181.27 -1.63 -9.87E−4 2011 2007 2232.70 -1.63 -9.87E−4

2011 2011 1852.01 -2.64 0

Table 5.8: Comparing parameter values

It is safe to say the the β2parameter is constant over the years, it is almost always zero, or very

close to zero. On the other hand, the β1 parameter is different is both years, so it is hard to

believe that this parameter stays constant over the years. The β0 parameters, the one which

(30)

Chapter 6

Conclusion

In this master thesis real time data was filtered and the well behaviour was modeled. A distinc-tion was made between short-term and long-term forecasting. For the short-term forecasting the fit of the model in the operating area was most important, whereas for long-term forecasting it was important that the behaviour of the model outside the operating area was also as expected. Besides that, for long-term prediction the reservoir pressure is an important factor which needs to be taken into account.

First the data, provided by the NAM, was filtered. Two different methods were applied, the SNR and MAD filter, where the SNR had an additional option; include or exclude stable closed in data. Both methods worked equally well, but due to the additional option the SNR filter was preferred over the MAD filter. The filter techniques described in this thesis are fast, easy to use and easy to implement. These are great advantages, since in the old situation different programs were needed to get to the desired results.

It was shown when fitting regression models to real time data that the models gave a much better fit compared with the currently used method, especially for short-term forecasting. Dif-ferent types of regression models were inspected, such as linear, nonlinear and nonparametric models. The AFBC model was handled as a linear regression model with restrictions (all AFBC parameters need to be nonnegative). Within this real time data modeling the AFBC model gave the best result, and no other regression model improved the fit. In addition to this the NAM prefers this model above the other regression models, since it can be used for short-term and long-term forecasting.

(31)

Figure 6.1: Model (5.1) Figure 6.2: Model (5.2) In addition for the future, when implementing this method in GENREM, different AFBC-Pres restrictions (upper and/or lower bound based on observations and petroleum operational knowl-edge) can be set when fitting the model to the data. This of course will deteriorate the fit, but might still be a major improvement compared to the model currently used, especially for short term forecasting.

Further research

The first thing that needs to be remarked is that this research is done for eight wells (where there are 200 wells in the Netherlands) based on a monthly dataset. These eight wells do not give enough information to draw a general conclusion for all other wells. It might be wise to fit the models described in this thesis on data for other wells and other months.

(32)

Bibliography

http://nl.wikipedia.org/wiki/Aardgas. http://nl.wikipedia.org/wiki/Gronings_gas. http://www.naturalgas.org. http://www.nam.nl. http://www.netlib.org/port/. http://www.eia.gov. http://www.enerdata.net.

Aja-Fern´andez, Santiago and Carlos Alberola-L´opez (2006). On the estimation of the coefficient of variation for anisotropic diffusion speckle filtering. IEEE TRANSACTIONS ON IMAGE PROCESSING 15, 2694–2701.

Baldwin, D.N., Suke B. Pillow J.J. Roiha H.L. Minocchieri S. and U. Frey (2004). Effect of sighs on breathing memory and dynamics in healthy infants. J Appl Physiol 97, 1830–1839. Faraway, Julian J. (2005, December). Extending the Linear Model with R: Generalized Linear,

Mixed Effects and Nonparametric Regression Models (1 ed.). Chapman and Hall/CRC. Gasterra (2009). Onzichbaar goud. Waanders B.V., Uitgeverij.

Kennedy, William J and James E. Gentle (1980, March). Statistical Computing (1 ed.). Marcel Dekker, Inc.

Kim, M., Ho Yoon S. Domanski P. and V. Payne (2008). Design of a steady-state detector for fault detection and diagnosis of a residential air conditioner. Elsevier 31, 790–799.

NAM. Manual Prototype WellMon (Version 1.0.0. ed.). NAM.

ORTEC (2011, March). GENeralized Reservoir Evaluation Model, user manual (Version 2011.0.1 ed.). ORTEC.

del Pilar Moreno, Rocio (2004). Steady State Detection, Data Reconciliation, and Gross Error Detection: Development for Industrial Processes. Master’s thesis, Francisco Jose de Caldas District University.

(33)

Rice, J. (2007). Mathematical Statistics and Data Analysis. Thomson.

Ritz, Christian and Jens Carl Streibig (2008, November). Nonlinear Regression with R (1 ed.). Springer.

(34)

Appendix A

Glossary

A Darcy flow coefficient (bar2_/(103_m3_/d))

B Static pressure correction coefficient BHP Bottom Hole Pressure (bar)

C Carbon

C Tubing friction coefficient (bar2_/(103_m3_/d)2₎

CH4 Methane C2H6 Ethane C3H8 Propane C4H10 Butane C5H12 Pentane C6H14 Hexane

CIBHP Closed in bottom hole pressure CITHP Closed in tubing head pressure

CO compression

CO2 Carbon Dioxide

F Non-Darcy flow coefficient (bar2_/(103_m3_/d))

(35)

GENREM GENeralized Reservoir Evaluation Model

H2 Hydrogen

MAD Mean Absolute Deviation

µ Mean of a dataset

m number of independent variables

n number of observations

N2 Nitrogen

NA Not Available

NAM Nederlandse Aardolie Maatschappij N2O Nitrous oxide

NP no production

O2 Oxygen

p order of polynomial

Pres Reservoir pressure Q Capacity (103_m3_/d)

SNR Signal to Noise Ratio

SS Sum of Squares

σ Standard deviation of a dataset THP Tubing Head Pressure (bar)

UGS Under Ground Storage

(36)

(37)

Appendix B

Map of Dutch gasfields

(38)

Appendix C

Data and models

Here the AFBC model for all other wells/datasets provided by the NAM will be discussed, as well as the other models fitted to the data of well W1.

C.1 Well W1

This section will handle the linear, nonlinear and nonparametric models fitted to the data of well W1. These models are used for making short-term predictions, where the focus will be on the (currently) operating area.

C.1.1 Linear Models

Several linear models were fitted to filtered data of well W1, the SNR filter without closed in data is used.

A short overview of the models used

THP = β0+ β1Q + e (C.1)

THP = β0+ β1Q2+ e (C.2)

THP = β0+ β1Q3+ e (C.3)

THP = β0+ β1Q + β2Q2+ β3Q3+ e (C.4)

(39)

Figure C.1: Linear models Model β0 β1 β2 β3 R2 (C.1) 50.21 -0.07 - - 0.974 (C.2) 36.08 -8.83E−5 - - 0.981 (C.3) 31.45 -1.43E−7 - - 0.982 (C.4) 219.12 -1.44 3.63E−3 -3.16E−6 0.986 (C.5) 23.44 -1.02E−216 - - 3.48E−4

Table C.1: Result linear models

Looking at the R2 5_{values in table C.1 and the models in figure C.1 almost all (except (C.5))}

give a good reflection of the well behaviour at the operating area. Comparing it with the R2 values of the AFBC model there is no (big) increase, and therefore not a big improvements of the fit in the operating area.

C.1.2 Nonlinear Models

Several nonlinear models were fitted to filtered data of well W1. A short overview of the models used

(40)

THP = β0+ β1Qβ2+ e (C.6)

THP = β0+ β1exp(Q/β2) + e (C.7)

THP = β0+ β1exp(Q/β2) + β3Q + e (C.8)

When fitting nonlinear models to a dataset the Gauss-Newton method is used. A big downside of this procedure is that starting values need to be given in order to find the parameters. It is important that these values are close to the optimal values, otherwise a local minimum or no solution is found.

For model (C.6) starting values are choses as the parameter values of model (C.2) or model (C.3), with β2 = 2 or 3. For model (C.7) starting values of model (C.2) are used, β2 is varied

with a factor of 10. The β2 value of the curve closest to the data is used as a starting value.

Model (C.8) is an expansion of model (C.7), the same starting values are used for βi(i= 0, ... ,

2) and a small negative number is chosen for β3.

Figure C.2: Starting values Figure C.3: Nonlinear models

Model β0 β1 β2 β3 R2

(C.6) 32.08 -4.90E−7 2.81 - 0.982

(C.7) 37.11 -2.37 217.74 - 0.982

(C.8) 44.08 -1.24E−5 39.11 -0.05 0.984 Table C.2: Result nonlinear models

All models have a high R2_{, meaning that the model fits nicely to the filtered data. Comparing}

it with the R2 _{values of the AFBC model there is no (big) increase, and therefore not a big}

(41)

C.1.3 Nonparametric Models

In R the function loess() is used to fit nonparametric models to the data. The default weight function, the tricube weight function, is used. In addition extrapolation is used, since the range of the filtered data is rather small (Faraway (2005)).

Figure C.4: Loess model (p=1) Figure C.5: Loess model (p=2)

(42)

C.2 Well W2

For this well data of march 2011 is available.

Figure C.6: Raw data well W2

(43)

C.3 Well W3

For this well data of march 2011 is available. This well is closed for a longer period, approximately two weeks. When the well is closed, the THP is also quite stable and it would make sense to include the closed in data for the SNR filter. If not, it will be harder to fit a right model for the well, since the model needs to go through one large cloud of points. This filter gives us the following dataset for well W3 if we set h=60 and the threshold of Q at 20, and the threshold of THP at 100.

Figure C.7: SNR filter for well W3

(44)

Figure C.8: Model (5.1) and (5.2)

(5.1) 6802.97 -28.40 0 0.967 (5.2) 6802.23 -26.53 - 0.965 Table C.3: Result model (5.1) and (5.2)

model _Pres A B F+C

‘old’ 83.4 15.834 1.6303 0.02002

(5.1) 83.4 29.04 1.02 0

(5.1) 105.31 46.30 1.6303 0

(5.2) 105.31 43.26 1.6303 0.02002 Table C.4: Fixing parameter values

When only fixing F+C it is not possible to obtain the rest of the parameters, since β2 is equal

(45)

C.4 Well W5

For this well data of march 2011 is available. This well is not closed during this period, both SNR filters will give exactly the same dataset. The SNR filter gives us the following dataset for well W5 if we set h=60 and the threshold at 120.

(46)

(5.1) 70121.01 -63.24067 0 0.963 (5.2) 63558.83 -48.61848 - 0.962 Table C.5: Result model (5.1) and (5.2)

model _Pres A B F+C

‘old’ 248.7 28.982 1.6208 0.01312

(5.1) 248.7 55.78 0.88 0

(5.1) 337.12 102.50 1.6208 0

(47)

C.5 Well W6

For this well data of march 2011 is available. This well is not closed during this period. The SNR filter gives us the following dataset for well W6 if we set h=60 and the threshold at 120.

(48)

(5.1) 14215.76 -9.942032 0 0.327

(5.2) 13673.04 0 - 0.039

Table C.7: Result model (5.1) and (5.2)

model _Pres A B F+C

‘old’ 140.04 3.6413 1.7795 0.027188

(5.1) 140.04 13.71 1.38 0

(5.1) 159.05 17.69 1.7795 0

(5.2) 155.98 0 1.7795 0.027188

(49)

C.6 Well W7

For this well data of march 2011 is available. Looking at the raw dataset it is very unlikely that the SNR filter with closed in data included will give a nice result, lowering the capacity will lead to a higher tubing head pressure, which is not the case for this dataset. Therefore the ordinary SNR filter is used with h=60 and a threshold of 100

(50)

(5.1) 52720 -3.683 -9.383E−3 0.735

(5.2) 63420.45 -18.18071 - 0.734

Table C.9: Result model (5.1) and (5.2)

(51)

C.7 Well W8

For this well data of march 2011 is available. In this dataset the well is closed for a really short period, the SNR filter without or with closed in data gives back the same results. The SNR filter gives us the following dataset for well W8 if we set h=60 and the threshold at 100.

(52)

(5.1) 189965.6 -54.2 0 0.997

(5.2) 115300.7 -4.51974 - 0.997 Table C.11: Result model (5.1) and (5.2)

model _Pres A B F+C

‘old’ 450.1 13.98 1.5361 0.01266

(5.1) 450.1 57.80 1.07 0

(5.1) 540.19 83.25 1.53611 0

(53)

Appendix D

GENREM

The GENeralized Reservoir Evaluation Model (GENREM) is a planning tool used by the NAM to make certain decisions in the gasfields. The model makes forecasts of the predicted capac-ity of the different gasfields. It models both the physical parameters of the subsurface fields and the surface network for transporting the gas to Gasunie. The well model used couples the subsurface and surface parts, and it is obvious to say that the well model is one of the key factors. The first version of GENREM was developed in 1989 by ORTEC. Nowadays ORTEC supplies support and maintenance of the software, but is also responsible for the development and imple-mentation of new GENREM functionalities.

D.1 Subsurface

In the subsurface part a reservoir model calculates the pressures in the reservoir based on the (inter)gridblock flows and well takeoffs. There are four options available to model the gasfield:

1. Use the relation between the cumulative production and the reservoir pressure. 2. No reservoir connection but only fixed well rates.

3. ResMod, a two dimensional one-phase (gas) simulator.

4. MoRes, a three dimensional three-phase (oil, water, gas) simulator.

ResMod is a fast and relatively simple model and is therefore preferred by the NAM as a reservoir simulator.

D.2 Well models

(54)

D.3 Surface network

The surface network is very important, since this is the most restraining part for the NAM. The GENREM program contains many building blocks to allow for an adequate modeling of reality. The network structure is set up by nodes and segments, where segments can be pipelines, con-densers or compressors. The nodes connect these segments. The user can specify a number of pressure and flow restrictions on every node in the surface network.

Before running the GENREM model the user must decide what time step to use, which is the duration of a simulation step. The default time step is one month, and this is used throughout this thesis. Within a time step GENREM will execute the following important steps.

First the reservoir pressure will be determined, which is the simulated or measured pressure below surface in the gasfield based on previous production volumes (subsurface). The reservoir pressure can be seen as an average pressure of the entire gasfield.

In the second step the maximum well capacities (the maximum daily production) are calcu-lated. This maximum capacity (at THP = 0) is calculated under fixed reservoir pressure, fixed delivery pressure (pressure at endpoints of the surface network) and all the restrictions in the network. With this maximum capacity the maximum offtake during a time step can be calculated as

maximum offtake = maximum capacity × timestep duration

Assuming that the demand is lower than the total of all maximum well capacities, not all max-imum capacity is needed. The partition of where to produce which amount of gas is done in the production step, where optimization takes place. This production step is based on a linear programming (LP) model, the objective is to maximize

X

i

ciqixi,

where ci stands for the priority of well i, some wells are more popular to use. qi stand for the

maximum capacity of well i, xi is the decision variable indicating the fraction to use of the

maximum capacity and is always smaller than or equal to one. The capacity, denoted by Q, for a well i in section D.2 is thus equal to qixi. To find the optimal allocation of the wells different

restrictions must be satisfied in the GENREM LP model, among others quality restrictions need to be satisfied. After solving the LP model production of the wells takes place (GENeralized Reservoir Evaluation Model, user manual (Version 2011.0.1.)).

(55)

(56)

Appendix E

Coefficient corrections for

long-term forecasts

When making a prediction on how the well will behave a few years from now, the AFBC curve needs to be adjusted for different reservoir pressure.

A = Aoriginal ZoriginalµoriginalZnewµnew F = Foriginal ZoriginalZnew B = _Boriginal C = CoriginalZ 2 new Z2_original

(57)

Appendix F

Gauss-Newton method

The Gauss-Newton method is a modification of the NewtonRaphson method, which minimizes f(x), as described in (Ratkowsky (1989)). The Gauss-Newton method estimates the parameters of the SS of a regression function, and makes use of the Jacobian. This method is used in finding the parameter values of a nonlinear regression model.

Recall that the (nonlinear) regression model is given by THP = f (Q, β) + e and the sum of squares by

SS = n X i=1 (THPi− f (Qi, β))2 or SS = (THP − f (Q, β))|(THP − f (Q, β))

Where THP and f (Q,β) are vectors of length n, this notation will be used from now on. The Jacobian matrix is then defined by

J(Q, β) =         ∂f(Q1,β) ∂β1 ∂f(Q1,β) ∂β2 · · · ∂f(Q1,β) ∂βp ∂f(Q2,β) ∂β1 ∂f(Q2,β) ∂β2 · · · ∂f(Q2,β) ∂βp .. . ... ... ... ∂f(Qm,β) ∂β1 ∂f(Qm,β) ∂β2 · · · ∂f(Qm,β) ∂βp        

Using the first two terms of a Taylor series

f (Q, β) ∼= f (Q, βk) + J(Q, βk)(β − βk)

where k indicates the kth iteration. Using this for f (Q,β) in the sum of squares results in the following expression

SS ∼= (THP − f (Q, βk) − J(Q, βk)(β − βk))|(THP − f (Q, βk) − J(Q, βk)(β − βk)) ∼

= SSk− 2(THP − f (Q, βk₎₎|_{J(Q, β}k_{)(β − β}k_{) +}

(58)

The last two terms together present the gradient vector and needs to be equal to zero to obtain the minimum for SS. Rearranging this leads to the following expression

(59)

Appendix G

Coefficient of Determination

Usually the coefficient of determination, the statistic R2, is used to give a conclusion about the fit of a model. It measures the strength of the model fitted by least squares and is calculated by as described in (Rice (2007)) R2 = 1 −SSerr SStot SSerr = n X i=1 (yi− fi)2 SStot = n X i=1 (yi− ¯y)2

Comparing R2 values of different models, it is not necessarily true that the model with the highest R2 has the best fit. Comparing different models with each other can be done by using the adjusted R2, and is defined by

¯

R2 = 1 − (1 − R2) n − 1 n − p − 1

The datasets used in this thesis are very large, around 30.000 data points. The fraction in this expression will be very close to one, and would lead to the following

¯

R2 ≈ R2

for this reason the R2 _{of the different models will be compared, since it will give to the same}

Improving gas capacity forecast model by mapping real time data to curves