Proxy Modeling of the Production Profiles of SAGD Reservoirs Based

Jul 15, 2015 - ... S. A. ; Nikolaou , M. MPCI: A New Approach to Closed-loop Identification and Adaptive Control. Presented at FOCAPO, Snowbird, UT, 1...
0 downloads 0 Views 2MB Size
Subscriber access provided by University of Victoria Libraries

Article

Proxy Modeling of the Production Profiles of SAGD Reservoirs Based on System Identification Song Yao, Japan J Trivedi, and Vinay Prasad Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/ie502258z • Publication Date (Web): 15 Jul 2015 Downloaded from http://pubs.acs.org on August 14, 2015

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Industrial & Engineering Chemistry Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Proxy Modeling of the Production Profiles of SAGD Reservoirs Based on System Identification Song Yao1, Japan J Trivedi1* and Vinay Prasad2 School of Mining and Petroleum, Department of Civil and Environmental Engineering1 Department of Chemical and Materials Engineering2 University of Alberta, Edmonton, Alberta: T6G 2G6, Canada.

KEYWORDS: Steam assisted gravity drainage, reservoir simulation, proxy model, system identification, prediction error method, input design, subspace identification. ABSTRACT Large scale physics-based reservoir models are employed routinely in the prediction of the behaviour of steam assisted gravity drainage (SAGD) processes under different operational situations. However, parametric uncertainty persists in these models even after history matching with production data. This uncertainty, and the computational cost associated with the full-scale reservoir simulations, makes it challenging to use reservoir simulators in closed-loop control of reservoirs. As an alternative strategy, we present in this work a dynamic proxy model for the *

Corresponding Author: Japan J Trivedi

E-mail: [email protected]

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

reservoirs based on system identification and the prediction error method using only injection and production data. These proxy models are validated against field data from a SAGD reservoir and simulated synthetic reservoir data and shown to be appropriate for use in model predictive control. We also provide evidence that the predictive power of these models can be improved by the appropriate design of input signals (injection rates and pressures).

1. INTRODUCTION The recent past has seen large improvements in the development of commercial simulators for modeling the behaviour of petroleum reservoirs. These simulators have been used for history matching1 and long-term optimization of reservoir performance, 2 among other tasks. A major limiting factor in expanding their use is the computational cost of performing a high-fidelity reservoir simulation. Furthermore, as reported by Van Essen et al.,3 reservoir models have some shortcomings, including difficulty in representing near-well bore reservoir dynamics such as gas or water coning and accounting for operational activities such as breakdown maintenance or well interventions. If the primary objective is closed-loop control of the reservoir on the time scale of days, weeks or months, an alternative to using reservoir simulations is to develop proxy models based purely on dynamic injection and production data obtained from the reservoir itself. That is the approach explored in this investigation, with purely data driven models being developed that replicate reservoir dynamics with respect to selected input variables. Proxy models are widely used in different areas of reservoir modeling such as sensitivity analysis, uncertainty propagation, probabilistic forecasting and reservoir management and optimization. For example, Jalali et al.4 proposed a technique that used artificial neural networks (ANN) to build a surrogate reservoir model (SRM) to perform uncertainty analysis on a coal bed

ACS Paragon Plus Environment

Page 2 of 37

Page 3 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

methane reservoir. A small number of realizations of the reservoir were used to develop the SRM, whose use reduced the simulation time significantly. Yeten et al.5 used experiment design and response surface methodology along with linear regression and ANNs to build statistical proxy models. Slotte et al.6 used multi-dimensional polynomials and kriging to construct a proxy model for assisted history matching and uncertainty assessment. Vanegas et al.7 built a proxy model based on Butler’s theory to perform uncertainty assessment of SAGD performance. Fedutenko et al.8 proposed a method for predicting production performance during SAGD process using Green’s functions and an interpolation-based proxy model. Azad et al.9 recommended an analytical proxy model to mimic the essential features of a SAGD field for history matching. Ghasemi et al.10 described an alternative proxy approach to model SAGD with an isothermal black-oil reservoir simulator by using the solution gas to oil ratio as a surrogate for temperature. Mohajer et al.11 used a self-organizing map based proxy reservoir model in a realtime optimization framework. Proxy models have also been used in probabilistic forecasting of reservoir performance with Monte Carlo sampling,

12

and in reservoir management and

optimization. The most widely used advanced control strategy in the chemical process industry is model predictive control (MPC),13 which has become popular on the basis of theoretical results provided by the academic community and successful installations in the industry.14 Model predictive control works in an optimization framework and thus allows for the incorporation of constraints into the controller design. There have also been a few studies, mainly in simulation, applying model predictive control to reservoirs. Proxy model development is a necessary requirement for the implementation of MPC in reservoir control in order to reduce the computational costs associated with objective function evaluation in the optimization. Nikolaou

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

and coworkers15-21 have used model predictive control in conjunction with proxy models based on subspace identification to develop constrained model predictive control of production in conventional reservoirs. Jansen and coworkers3 have also investigated the use of subspace identification based proxy models in MPC and optimization. In this work, we present a method of dynamic system identification based on the prediction error method to model the dynamic behaviour of SAGD reservoirs. In contrast to the subspace identification procedure, this allows us to extend the models easily to recursive and nonlinear structures, which will be shown to provide improved predictive capabilities. Also, the identified models retain the ability to be converted to the state space form suitable for controller design. To the best of our knowledge, this method has not been explored thoroughly in the literature on reservoir proxy modeling except in the work of Renard et al.,22 who applied it for water cut analysis in water flooded layered reservoirs. We analyze the prospects for the use of such models in characterizing SAGD operation, including investigations on the selection of input variables and their excitation for improved identification, the structure of the identified models, and the performance of the models in short and long term prediction at different parts of the reservoir’s lifetime. Thus, our main contribution lies in a careful assessment of the suitability of system identification for model development for SAGD process control. 2. METHODOLOGY In this work, we address the issue of proxy model development for reservoirs based on the identification of empirical low order models from injection and production data. This section describes the methods for model development and the metrics for evaluating the performance of models, the type of data used in the identification, and methods of input design with the potential to generate data that significantly improves the accuracy and predictive capabilities of the proxy

ACS Paragon Plus Environment

Page 4 of 37

Page 5 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

models. The methods of input design are demonstrated using simulations with a synthetic reservoir model, and the details of the reservoir model are also presented in this section.

2.1. Model development In the tradition of classical system identification, the prediction error methodology23,24 is used in this work. In this approach, the model is described as a predictor of the next output, and suitable parameterizations are chosen in order to construct the predictor. Given a model parameterization and an observed data set, the estimates of the parameters of the model are chosen such that the error between the observed data and the model predictions is minimized in a suitable norm. The model can be described as

y (k ) = G ( z −1 , θ )u (k ) + H ( z −1 , θ )e(k ) ............................................................................................(1) where G represents the process transfer function, characterizing the relation between the inputs u and the outputs y, and H is the disturbance transfer function, which describes the relation between the noise/disturbances and the outputs. The disturbance/noise term is assumed to be driven by a Gaussian white noise process e(k). The index ‘k’ indicates the current time step, and z −1 is the backward shift operator. By expressing G and H as functions of the backward shift operator, the process and transfer functions can include the effect of inputs and disturbances from previous times on the current output. A schematic of this approach to specify the structure of a model for dynamic systems is provided in Figure 1. As mentioned above, we choose the prediction error method (PEM) in our work because it provides the option to be used with a wide variety of model parameterizations, including linear and nonlinear choices. The particular choices explored in our work are described below.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 37

2.2. Linear model structure In the linear models considered, the process and disturbance transfer functions, G ( z −1 ) and

H ( z −1 ) respectively, are specified as rational functions of the backward shift operator of the form

B( z −1 ) G( z ) = ..................................................................................................................(2) A( z −1 ) F ( z −1 ) −1

H ( z −1 ) =

C ( z −1 ) ...................................................................................................................................... (3) A( z −1 ) D( z −1 )

where

A( z −1 ) = 1 + a1 z −1 + a2 z −2 + L + anA z − nA , B( z −1 ) = b0 + b1 z −1 + b2 z −2 + L + bnB z − nB , C ( z −1 ) = 1 + c1 z −1 + c2 z −2 + L + cnC z − nC , .......................................................................................(4) D( z −1 ) = 1 + d1 z −1 + d 2 z −2 + L + d nD z − nD , F ( z −1 ) = 1 + f1 z −1 + f 2 z −2 + L + f nF z − nF are polynomials in  , the backward time-shift operator. Commonly used variants of the linear model structure are the ARX (autoregressive with exogenous input) model, the ARMAX (autoregressive moving average with exogenous input), the OE (output error) model and the BJ (Box-Jenkins) model. The choices of the polynomials that lead to these model structures are provided below:

ACS Paragon Plus Environment

Page 7 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

ARX model: G ( z −1 ) =

B( z −1 ) 1 , H ( z −1 ) = ...............................................................................(5) −1 A( z ) A( z −1 )

ARMAX model: G( z −1 ) =

B( z −1 ) C ( z −1 ) −1 , H ( z ) = .........................................................................(6) A( z −1 ) A( z −1 )

B( z −1 ) , H ( z −1 ) = 1 ............................................................................................(7) OE model: G( z ) = −1 A( z ) −1

BJ model: G( z −1 ) =

B( z −1 ) C ( z −1 ) −1 H z = , ( ) ...................................................................................(8) A( z −1 ) D( z −1 )

When multiple inputs are used in the model, a separate process model G ( z −1 ) is specified with respect to each input. If there is a delay of  samples before the inputs first affect the outputs, k leading terms of the corresponding polynomials will be absent. The p-step ahead prediction is of the form

yˆ (k | k − p) = W p ( z −1 )G( z −1 )u (k ) + [1 − Wp ( z −1 )] y (k ) p −1

Wp ( z −1 ) = H p ( z −1 ) H −1 ( z −1 ), H p ( z −1 ) = ∑ h(i ) z −i ................................................................................... (9) i =0

The one-step ahead predictor is derived as a special case of the p-step ahead predictor. In the prediction error method, the model described by equation (9) is used to generate a one-step ahead predictor for the output based on data available up to the current sampling instant, and equation (9) is used to generate the p-step ahead predictions (Ljung, 1987). The p-step ahead predictions of the model are compared with process data in a batchwise fashion, i.e., predictions over a specified time interval are compared with the data.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 37

2.3. Parameter estimation The set of parameter values that minimizes the error between the one-step ahead model predictions (except in the case of the OE model) and the measured outputs in a least squares sense is chosen as the optimal parameter set. For the OE model, the minimization of the errors between measurement data and model simulation (i.e., the so-called infinite horizon prediction) is sought. The one-step ahead prediction error is given below: ) e(k ) = y ( k ) − y (k | k − 1) ...............................................................................................................(10) Usually, the sum of squared prediction errors is used as the loss function:

L = ∑i =1 e 2 (i ) ................................................................................................................................(11) N

where N is the number of data samples.

2.4. Model validation using p-step ahead predictions While the parameter estimation is based on one-step ahead errors, model validation is based on the goodness of fit for the p-step ahead model predictions, which is characterized by the following measure:

 ∑N ( y ( j ) − yˆ ( j | j − p) 2 )   .............................................................................. (12) fit % = 100 × 1 − j =1 N 2  ∑ j =1 ( y( j) − y ( j ))  

where $y denotes the estimated output, y is the measured output, y is the mean of the measured output data and N represents the number of data samples. In addition, the autocorrelation of the residuals (errors between the model predictions and the measured values), and their cross-

ACS Paragon Plus Environment

Page 9 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

correlation against the inputs are investigated to ensure that there is no systematic dynamic behaviour left unmodeled. The Akaike information criterion (AIC) is used as the measure to compare models with different structures, with the model with the lowest AIC being chosen.

AIC = log V +

2d .........................................................................................................................(13) N

where V is the loss function (a quadratic measure of the residuals of the model with respect to the entire data set), d is the number of parameters in the model and N is the number of data points in the estimation data set.

2.5. Recursive ARX models When linear time-invariant models are inadequate to represent the behaviour of a process, the alternatives that are available are nonlinear models. This can be achieved using pseudo-linear parameter varying models in which the parameters of linear models such as the ARX model are updated recursively with each new data point that is acquired.25 This is applied with the assumption that the unknown parameters are time-varying. Even if the true parameters are constant (but unknown), allowing the estimated parameters to vary with time provides a recursive update so that the parameter estimates can converge to the true values as additional data is evaluated. In this work, we explore the recursive ARX model structure:

y (k ) = −∑i=a1 ai (k ) y (k − i ) + ∑i=b1 bi (k )u (k − i ) ..............................................................................(14) n

n

θ (k ) = [a1 (k ) L an (k ) b1 (k ) L bn (k )] .......................................................................................(15) a

b

) and θ (k | k − 1) is a complex function of the past inputs and outputs. The general recursive identification algorithm is given by the following equation:

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

θˆ(k ) = θˆ(k − 1) +

Page 10 of 37

P(k − 1)Ψ (k ) [ y(k ) − y)(k | k − 1)] ......................................................(16) T λ + Ψ (k ) P(k − 1)Ψ (k )

...................................................................................................................................................... where θ$ is the parameter estimate at time k. The model prediction is given by

yˆ (k | k − 1) = Ψ T (k )θˆ(k − 1) ...........................................................................................................(17) ...................................................................................................................................................... The estimation algorithm minimizes the prediction error, and P(k) is given by

1 P(k − 1)Ψ(k )ΨT (k ) P(k − 1)  P ( k ) =  P ( k − 1) −  ......................................................................(18) λ λ + ΨT (k ) P(k − 1)Ψ(k )  where λ is a forgetting factor ( 0 ≤ λ ≤ 1 ), which results in older prediction errors having lesser weight in the regression than newer prediction errors. Ψ (k ) is the regression vector that is computed based on previous values of measured inputs and outputs, whose specific form depends on the structure of the polynomial model:

Ψ(k ) = [u (k − 1) L u (k − nb ) − y(k − 1) L − y (k − na )] .................................................(19) T

The results described in this work have been obtained using a forgetting factor of λ = 0.98 .

2.6. Reservoir data and models The data used in this study was obtained from 15-year operation with sampling time of 1 day of a SAGD reservoir with three horizontal well pairs located in northern Alberta. The depth from the surface was approximately 125-175 meters, with net pay thickness approximately 15 meters. The

ACS Paragon Plus Environment

Page 11 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

sand porosity and permeability were approximately 35% and 5-12 Darcy, respectively, and average bitumen saturation was approximately 85%. All the results are presented in terms of normalized variables for the reason of confidentiality of the raw data. The production rates presented throughout are subtracted by the linear trends. Reservoir simulations were also conducted to demonstrate the potential effect of input design on the predictive capabilities of the proxy models. The Advanced Process and Thermal Reservoir Simulator STARS 2012.1027 was used to model the SAGD process. These simulations were carried out on a synthetic reservoir that closely matched the characteristics of the reservoir from which the production data was obtained; however, the reservoir models were not rigorously history-matched. A 3-D SAGD reservoir model, with 11 × 10 × 10 grid blocks with model dimensions of 44 meters × 1400 meters ×30 meters (I×J×K), has been used in this case study. A pair of injector and producer wells is used in this simulation model. As we discussed before, steam injection rate and production pressure will be used as in the input signals, and the output signals are water rate and bitumen rate in the producer. A schematic of the SAGD process and variables used for the proxy modeling are shown in Figure 2.

2.7. Input design Proper design of the input sequence for which output data is collected is essential to ensure good estimation performance. The input sequence must be persistently exciting26 to the order of the linear model being identified. A signal is said to be persistently exciting of order n if the following limit exists: 1 N →∞ N

ru (τ ) = lim

N

∑u

t +τ

utT

(20)

t =1

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 37

and the matrix given below is positive definite:

ru (1) … ru (n − 1)   ru (0)  r (−1) r (0) … r (n − 2)  u u  Ru ( n ) =  u  M  M M M   ru (0)   ru (1 − n) … …

(21)

Input sequences such as pseudorandom binary sequences (PRBS) are known to be persistently exciting for finite order models.23

3. RESULTS AND DISCUSSION The workflow for proxy modeling used in this work is shown in Figure 3. First, the input variables that are most useful in estimating the chosen output data are identified. Next, the model structure (OE, ARX, ARMAX or BJ) that is most effective at describing the observed output behaviour is identified. Input variable selection is performed using one-step ahead model predictions, since it is used to screen the most promising approaches to system identification from reservoir production data. Once the set of input variables to be used in model development is identified, the determination of the structure of the model and the estimation of its parameters is based on multi-step ahead predictions. In the case of OE models, the estimation is based on infinite horizon prediction. Next, the quality of input datasets for system identification is investigated in two respects: first, the period in the reservoir’s lifetime that is best suited for proxy model development using system identification is identified; second, the time periods over which identified linear models are valid without re-identification using new data are investigated. These investigations are also conducted using one-step ahead model predictions, since they are also screening studies. Finally, proxy model estimation and validation based on the chosen input variables and model structures are carried out using multi-step ahead predictions. The results of

ACS Paragon Plus Environment

Page 13 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

7-step ahead predictions are reported, although we have also investigated the effect of an increase in the prediction horizon, since that is a key variable in the design of model predictive controllers that would use these identified models. Further data preprocessing to improve proxy model identification is then investigated, along with input designs known to provide improved identification. In the proxy model estimation and validation using linear models without recursive updating, the testing data set is 2 months, 6 months or 1 year long, and the validation data set is of that same length or a longer time period. As noted earlier, all the data sets have a sampling period of 1 day. We have employed the System Identification Toolbox of MATLAB to develop the models, except in the case of recursive estimation, which we have coded in MATLAB.

3.1. Selection of input variables The first step towards developing a model for a SAGD reservoir is to assign input and output variables. The water and bitumen production rates are natural choices for the output variables of the proxy model. System identification is performed for each of these outputs separately considering multiple inputs, i.e., identification is performed in a multiple input-single output (MISO) framework. It has been established that real-time field data can be used in the development of reservoir proxy models.,15. Van Essen et al.3 suggest the use of steam injection rate in the injector and bottomhole pressure in the producer well as input variables. In addition, we also considered the injection pressure as a potential input. MISO models were constructed considering paired combinations of the input variables, and also considering all three of the input variables in the model. An ARX model structure was chosen for this exercise, since the objective was just to

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 37

screen for the most appropriate combination of input variables. Figures 4 and 5 show the onestep ahead predictive ability (over 70 days) of the best ARX model that was obtained, which included steam injection flow rate and the bottom-hole pressure in the producer well as the inputs. Figure 4 shows the performance of the ARX model on the training data, and Figure 5 is performance on the validation data set. The input data is depicted in Figure 6. The input data of the first 70 days is used for model training and the data from day 71 to day 140 is for model validation. The input sequence was confirmed to have persistent excitation for the models tested based on the conditions described in section 2.5. The identification performance with other input combinations on the training data can be found in Table 1, and the combination with the lowest value of AIC was chosen. The orders of the different models investigated for input variable selection are shown in Table 2. The performance of the ARX model was satisfactory, indicating that the input variables used here are appropriate to capture the reservoir dynamics. Note that the production data is mean-centered and has linear trends subtracted, and the plots display the deviation variables thus obtained, which is why negative values are possible. Another test that was performed to evaluate the suitability of the inputs for inclusion in the proxy model was to evaluate the prediction performance of an autoregressive (AR) model without an exogenous input. If the AR models were to provide adequate prediction performance, this would imply that the dynamic trends in the outputs can be explained purely in terms of the values of the output variables in the recent past. However, the AR models were extremely poor in prediction of both water and bitumen production rates, indicating that using the specified input variables is critical for explaining the dynamic variation in the outputs.

3.2. Model structure selection

ACS Paragon Plus Environment

Page 15 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Once the input and output variables have been identified, the next step is to identify the most suitable model structure that best describes the reservoir dynamics. We evaluate the four linear model structures described earlier (OE, ARX, ARMAX and BJ) using infinite horizon predictions, and the results for 70 days of prediction for the validation data set are presented in Figure 7. And because each model has different number of parameters, it is necessary to employ the Akaike Information Criterion (AIC) to take both the number of model parameters and prediction errors into consideration. The AIC values are shown in Table 3, and the most accurate model has the smallest AIC. The Box-Jenkins model has the best performance, though the ARX and ARMAX models also provide reasonable prediction performance. The OE model proved to be completely unsuitable for use as a proxy model. Based on these results, the Box-Jenkins model is the linear non-recursive model of choice for further investigations. The orders of all models are reported in Table 4. Residuals of the models are also investigated to ensure that there is no systematic dynamic behaviour left unmodeled, as mentioned in section 2.2. The residuals for the ARX model on the identification data set are depicted in Figure 8. Also, the gains and time scales of the identified models are consistent with the known physics of the system. For example, the gain for the water production rate with respect to the steam injection rate is positive, since more water is produced at higher injection rates, as is shown in Table 5. Also, when the other input (the production pressure) is increased, due to the fact that the injection pressure is also increased correspondingly, the water production rate increases, and these gains are positive, too. However, most of the models describing the dependence of the bitumen production rate on the inputs have negative gains as a result of data detrending. The untreated raw data have positive gains (not shown here because of confidentiality). The time scales of most of the identified models (shown

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 37

in Table 6) range from 1-16 days, which is reasonable for the response of such reservoirs to changes in the injection rate and the injection/production pressure. The first term in nB represents the order of the process models with respect to the steam injection rate and the second its order with respect to the bottom-hole pressure. A delay of 1 time unit with respect to each input was identified for the water production rate model, while a delay of 3 time units was identified with respect to each input for the bitumen production rate model. Subspace identification methods were also found to identify models of comparable accuracy in one-step ahead prediction for this data set (see Supplementary Information, Figure S1); however, we chose to use the models described above to retain the possibility for recursive updating of the model’s parameters. The best subspace models identified were a second order subspace model for the water rate and a third order model for the bitumen rate.

3.3. Time period for identification The quality of identification obtained using linear non-recursive models was investigated next. First, the period in the reservoir’s lifetime over which the best identification performance can be achieved was investigated. Production data was available over 15 years for the SAGD reservoir, and BJ models were identified over two month periods in the first, third, eighth and thirteenth years of the reservoir’s lifetime. Figures 9 and 10 describe the one-step ahead prediction performance of BJ models on validation data sets in the first and third year of the reservoir’s lifetime, respectively. The order of water rate model is nB = [2 2], nC = 2, nD = 2, nF = [2 2] and the order of bitumen rate model is nB = [2 2], nC = 2, nD = 2, nF = [2 2] in Figure 9. The order of water rate model is nB = [1 1], nC = 1, nD = 1, nF = [2 2] and the order of bitumen rate model is nB = [6 1], nC = 1, nD = 1, nF = [1 1] in Figure 10. While the model is reasonably

ACS Paragon Plus Environment

Page 17 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

accurate in describing the production in the third year, it is much less successful in describing the production data in the first year. While the results are not shown, the model performed adequately in the eighth year of production, and the prediction performance in the thirteenth year was comparable to that in the first year. This indicates that linear non-recursive models are inadequate in the early part of the reservoir’s lifetime, when steam chambers and production profiles are yet to fully stabilize, and in the final part of its lifetime, when saturation and production profiles reduce significantly. For those periods, nonlinear or recursive linear models may be required even for single-step ahead predictions. The other aspect investigated was the time period over which non-recursive linear models are valid without re-identification using new data. Figures 11 and 12 show the prediction performance of BJ models over a six month and a one year period, respectively, and indicate that these models are adequate for one-step ahead prediction over these time periods. Note that this conclusion is drawn using data taken from the intermediate portion of the reservoir’s lifetime (i.e., between years 3 and 10).

3.4. Multi-step ahead prediction While one-step ahead predictions are useful in screening studies to identify appropriate inputs and model structures, the primary objective in building a dynamic proxy model is to obtain good multi-step predictions. Model predictive controllers use models to predict and optimize the effect of future manipulated variable changes on the outputs over a prediction horizon into the future. A prediction horizon of one week is typical in MPC for reservoir control;3 thus, with a sample time of 1 day, we choose 7-step ahead predictions as our basis for evaluating models. The linear non-

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 37

recursive models described above were incapable of providing meaningful predictions in 7-step ahead prediction, which led us to investigate recursive ARX models. In contrast, the recursive ARX models showed considerable promise in multi-step ahead prediction. Before providing results on multi-step prediction, we first provide results for one-step ahead prediction, noting that the recursive models are developed based on the minimization of one-step ahead predictions. Figure 13 shows that the recursive ARX model was able to provide highly accurate one-step ahead prediction for production data over one year of the reservoir’s lifetime, and Figure 14 shows that the parameters of the model are adapted significantly over that time period to allow the model to track the production data accurately. The recursive ARX models also provide accurate one-step ahead prediction over longer periods of two and five years (see Supplementary Information, Figures S2-S5). The order of water rate model is nA = 2, nB = [1 1] and the order of bitumen rate is nA = 1, nB = [1 1]. Next, we present the results of 7-step ahead prediction using recursive ARX models over 5 years of the reservoir’s production in Figure 15. The degree of fit was 59.1% for prediction of the water production rate, and 47.6% for prediction of the bitumen production rate. While a better fit would obviously be advantageous when using the models in closed-loop control, the use of state estimation with MPC is expected to mitigate the effects of model-plant mismatch on the control performance. Also, we investigated 10-step ahead prediction using these recursive ARX models; while the results are not presented here, the degrees of fit dropped to 45.6% for the water production rate and 33.7% for the bitumen production rate, which may still be adequate for use in MPC28. Normal operation of the reservoir is often interrupted by production shutdowns, which may occur for many reasons including maintenance, repair and safety considerations. The large downward steps (of the order of 100-200 T/d) in the water and bitumen production rates

ACS Paragon Plus Environment

Page 19 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

observed in Figure 15 are due to these production shutdowns. In the results shown thus far, the proxy models have been evaluated and recursively identified (where applicable) even during these shutdowns. However, it is reasonable to assume that knowledge of the shutdowns would be available to the modeler in the field, which implies that the model identification and evaluation can be conducted by freezing the model’s parameters over the period of the shutdown and excluding that time period from the identification. The results of 7-step ahead model prediction with recursive ARX models after this additional data preprocessing are presented in Figure 16. It is clear that the accuracy of model predictions is enhanced considerably with the exclusion of the shutdown periods. The same recursive ARX model structure was used to identify proxy models for the other two well pairs in the reservoir, and similar 7-step ahead prediction performance was achieved (results not shown). While the parameter values that were identified were different for the other well pairs, the same model structure was shown to be general enough to capture the dynamic production behaviour of other parts of the reservoir.

3.5. Input design for improved prediction The rank condition on the input data for persistent excitation was described in section 2.5. The input data obtained from the reservoir meets the condition for persistent excitation for models up to order 8, considering data over two months or longer time periods, which implies that the identification of lower order models is possible. This is, of course, confirmed based on the results described above. However, if carefully designed inputs are used to generate the dynamic production data, it is expected that much more accurate models can be developed. Since we do not have access to a real reservoir, we demonstrated this using synthetic production data from

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 37

SAGD reservoir simulations using a flex-well model developed in CMG STARS.27 A pseudorandom binary sequence (PRBS) was applied to the inputs (steam injection rate and bottom-hole pressure) of the reservoir simulation, and the water and bitumen production rate were modeled using non-recursive linear models. The average amplitude of the PRBS signals were less than 20% of the nominal values of the inputs, and its frequency bounds were calculated based on the Nyquist frequency and estimates of the corner frequencies using step responses of the reservoir model (Ljung, 1987). Box-Jenkins models of the same structure as those used for identification provided the best fits for the synthetic production data. The model performance was evaluated in open loop simulation, also referred to as infinite horizon prediction. This means that no data is presented to the model for use in predictions other than the initial condition. Figure 17 shows the prediction performance of the BJ model in open loop simulation for a period of two months, with 72.8% and 64.4% fits for the water and bitumen production rates, respectively. The model order of water rate is nB = [1 1], nC = 1, nD = 1, nF = [1 1] and the model order of bitumen is nB = [1 1], nC = 1, nD = 1, nF = [1 1]. The figure shows that the dynamics of the production data are captured very well, and only the process gain of the BJ model is slightly different from that for the reservoir simulation. In contrast, the model fits for BJ models identified using reservoir simulations with input profiles similar to that used in the real reservoir were much poorer (fits of 26.0% and 19.6% for water and bitumen production rates, respectively), indicating that the use of PRBS inputs to generate data for system identification is necessary to improve the accuracy of the identified models. A multi-level PRBS input was also investigated (see Supplementary Information, Figure S6); however, the single level PRBS input resulted in more accurate proxy model identification. The good performance of the model in infinite horizon prediction indicates, of course, that it can be used with very long prediction

ACS Paragon Plus Environment

Page 21 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

horizons in model predictive control. We expect the gains in model accuracy and length of prediction horizon to hold in practice, but that can only be proven using field tests.

4. CONCLUSIONS In this work, we have investigated the possibility of using system identification methods for proxy model development to describe the dynamics of SAGD operation using production data from an operational SAGD reservoir. The main findings are: it is possible to develop proxy models with multi-step predictive capability of sufficient accuracy to be used for model predictive control. Box-Jenkins models provide the most accurate single-step predictions among linear non-recursive models, but recursive ARX (autoregressive with exogenous input) models provide the best multi-step ahead predictions. The recursive ARX models provided accurate 7day ahead predictions of production data over periods as long as 5 years. Steam injection rates and bottom-hole production pressures were found to be the most appropriate inputs to build proxy models to describe water and bitumen production rates. Most of the proxy models were slightly more accurate in the prediction of water production rates than in the prediction of bitumen production rates. The developed proxy models are least accurate at the start of the reservoir’s lifetime (in the first year of production) and at its end (in the last two years); however, effective proxy models were obtained over the majority of the reservoir’s lifetime. For the reservoir we studied, the same model structure provided accurate predictions for production data from all the well pairs in the reservoir, indicating similarity in the qualitative dynamic behaviour of these well pairs. Finally, through proxy modeling of synthetic production data, we have shown that using a pseudo-random binary sequence in the input variables provides proxy models of much higher accuracy over very long prediction horizons, and recommend that these sequences be used for proxy model identification.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 37

ACKNOWLEDGEMENT The authors acknowledge financial support from the Natural Sciences and Engineering Research Council, Canada through Discovery Grants. The thermal reservoir simulator STARS was provided by Computer Modelling Group.

ASSOCIATED CONTENT Supporting Information Predictions using subspace identification, nonlinear ARX and recursive ARX models. The Supporting Information is available free of charge on the ACS Publications website at http://pubs.acs.org/.

REFERENCES (1) Oliver, D.S.; Chen, Y. Recent Progress on Reservoir History Matching: A Review. Comput. Geosci. 2011, 15(1), 185. (2) Sarma, P.; Durlofsky, L.J., Aziz, K.; Chen, W.H. Efficient Real-time Reservoir Management Using Adjoint-based Optimal Control and Model Updating. Comput. Geosci. 2006, 10(1), 3. (3) Van Essen, G.M.; Van den Hof, P.M.J.; Jansen, J.D. A Two-level Strategy to Realize Lifecycle Production Optimization in an Operational Setting. Presented at the SPE Intelligent Energy International Conference, Utrecht, Netherlands, March 27-29, 2012. (4) Jalali, J.; Mohaghegh, S.D.; Gaskari, R. Coalbed Methane Reservoir Simulation and Uncertainty Analysis with Artificial Neural Networks. Sci. Iran. Trans. C: Chem. Chem. Eng.

2010, 17(1), 65. (5) Yeten, B.; Castellini, A.; Guyaguler, B.; Chen, W.H. A Comparison Study on Experimental Design and Response Surface Methodologies. Presented at the SPE Reservoir Simulation Symposium, Houston, TX, January 31–February 2, 2005.

ACS Paragon Plus Environment

Page 23 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(6) Slotte, P.A.; Smorgrav, E. Response Surface Methodology Approach for History Matching and Uncertainty Assessment of Reservoir Simulation Models. Presented at the SPE Europec/EAGE Annual Conference and Exhibition, Rome, Italy, June 9-12, 2008. (7) Vanegas, J.W.; Deutsch, C.V.; Cunha, L.B. Uncertainty Assessment of SAGD Performance Using a Proxy Model based on Butler’s Theory. Presented at the SPE Annual Technical Conference and Exhibition, Denver, CO, September 21-24, 2008. (8) Fedutenko, E.; Yang, C.; Nghiem, L. Forecasting SAGD Process Under Geological Uncertainties Using Data-driven Proxy Model. Presented at the SPE Heavy Oil Conference, Calgary, Canada, June 12-14, 2012. (9) Azad, A.; Chalaturnyk, R.J. Application of Analytical Proxy Models in Reservoir Estimation for SAGD Process: UTF-project Case Study. J. Can. Petrol. Tech. 2013, 52(3), 219. (10) Ghasemi, M.; Whitson, C.H. Modeling Steam-assisted Gravity Drainage with a Black-oil Proxy. Presented at the SPE Annual Technical Conference and Exhibition, Denver, CO, October 30-November 2, 2011. (11) Mohajer, M.; Damas, C.; Silva, A.J.B.; Al-Kinani, A. An Integrated Framework for SAGD Real-time Optimization. Presentation at the SPE Intelligent Energy Conference and Exhibition, Utrecht, Netherlands, March 23-25, 2010. (12) Fishman, G.S. Monte Carlo: Concepts, Algorithms and Applications. Springer-Verlag, NY, 1996. (13) Bemporad, A.; Morari, M.; Ricker, N.L. Model Predictive Control Toolbox: User’s Guide, MathWorks, Natick, MA, 2014.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 37

(14) Borgesen, J.F. Efficient Optimization for Model Predictive Control in Reservoir Models. MSc. Dissertation, Norwegian University of Science and Technology, Trondheim, 2009. (15) Genceli, H.; Nikolaou, M. New Approach to Constrained Predictive Control with Simultaneous Model Identification. AIChE J. 1996, 42(10), 2857. (16) Saputelli, L.; Nikolaou, M.; Economides, M.J. Self-learning Reservoir Management. SPE Reservoir Eval. Eng. 2005, 8(6), 534. (17) Awasthi A.; Sankaran, S.; Nikolaou, M.; Saputelli, L.; Mijares, G. Closing the Gap Between Reservoir Simulation and Production Optimization. Presentation at the SPE Digital Energy Conference and Exhibition, Houston, USA, April 11-12, 2007. (18) Eker, S.A.; Nikolaou, M. Adaptive Control through On-line Optimization: The MPCI Paradigm. Appl. Math. Comp. Sci. 1999, 9(1), 101-128. (19) Nikolaou, M. Model Predictive Controllers: A Critical Synthesis of Recent Developments and Industrial Needs. In Advances in Chemical Engineering Series, Academic Press, Burlington, MA, 2001. (20) Schwarm, A.; Eker, S.A.; Nikolaou, M. MPCI: A New Approach to Closed-loop Identification and Adaptive Control. FOCAPO, Snowbird, Utah, 1998. (21) Schwarm, A.; Nikolaou, M. Chance-constrained Model Predictive Control, AIChE J. 1999, 45(8), 1743. (22) Renard, G.; Dembele, D.; Lessi, J; Mari, J. System Identification Approach to Watercut Analysis in Waterflooded Layered Reservoirs. Presented at the SPE/DOE Improved Oil Recovery Symposium, Tulsa, OK, April 19-22, 1998.

ACS Paragon Plus Environment

Page 25 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(23) Ljung, L. System Identification: Theory for the User; Prentice Hall: Englewood Cliffs, NJ, 1987. (24) Nells, O. Nonlinear System Identification: From Classical Approaches to Neural Networks and Fuzzy Models; Springer: Berlin, 2001. (25) Wellstead, P.E.; Zarrop, M.B. Self-Tuning Systems: Control and Signal Processing; Wiley: Chichester, New York, 1996. (26) Soderstrom, T.; Stoica, P. System Identification; Prentice Hall: Englewood Cliffs, NJ (1994). (27) Computer Modelling Group STARS Version 2012.10 User's Guide. Computer Modelling Group, Calgary, Alberta, Canada, 2012. (28) Vembadi, S.; Trivedi, J.; Prasad, V. Real-Time Feedback Control of SAGD Wells using Model Predictive Control to Optimize Steam Chamber Development under Uncertainty. Presented at the World Heavy Oil Congress, Edmonton, AB, March 24-26, 2015.

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 37

Table 1. Comparison of model accuracy for different input combinations. SIR: steam injection rate, PP: production (bottom-hole) pressure, SIP: steam injection pressure, AIC: Akaike information criterion. Fit percentage is based on equation (12). Input Combinations

Water Rate

Bitumen Rate

Fit %

AIC

Best Fit

AIC

SIR, PP

75.0

4.2014

73.0

2.7441

SIR, SIP

73.6

4.3385

72.7

2.9143

SIP, PP

72.5

4.3145

72.8

2.7551

SIR,SIP,PP

72.6

4.3302

72.2

2.8341

Table 2. Model orders for selection of input variables Input Combinations

Water Rate

Bitumen Rate

SIR, PP

nA = 2, nB = [1 1]

nA = 2, nB = [1 1]

SIR, SIP

nA = 3, nB = [1 1]

nA = 4, nB = [3 3]

SIP, PP

nA = 2, nB = [1 1]

nA = 2, nB = [1 1]

SIR,SIP,PP

nA = 2, nB = [1 1 1]

nA = 2, nB = [1 1 1]

Table 3. Value of the Akaike information criterion (AIC) for model selection Model

AIC (Water Rate)

AIC(Bitumen Rate)

BJ

4.1071

2.6985

ARX

4.2014

2.7022

ARMAX

4.3271

2.8101

OE

4.5233

4.8711

Table 4. Model orders for model selection Water Rate BJ ARX ARMAX

nB = [1 1], nC = 1, nD = 1, nF = [2 2] nA = 2, nB = [1 1] nA = 2, nB = [2 2], nC = 2

ACS Paragon Plus Environment

Page 27 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

OE

nB = [3 1], nF = [1 1] Bitumen Rate

BJ

nB = [6 1], nC = 1, nD = 1, nF = [1 1]

ARX

nA = 1, nB = [4 1]

ARMAX

nA = 1, nB = [4 1], nC = 1

OE

nB = [2 2], nF = [2 2]

Table 5. Model steady state gains for water and bitumen production rates with respect to the inputs: steam injection rate (SIR) and production pressure (PP). The gains are unitless since the raw data has been normalized to keep it confidential. Water Production Rate Model

Bitumen Production Rate

SIR

PP

SIR

PP

BJ

0.5228

0.2058

-0.6401

-0.79

ARX

0.8008

0.4450

-2.2555

0.0208

ARMAX

0.7328

0.1254

-1.9675

-0.0343

Table 6. The dominant time constants (in days) for the models describing the dynamic relations between the outputs (water and bitumen production rates) and the inputs (steam injection rate and production pressure) Water Rate Model

Bitumen Rate

SIR

PP

SIR

PP

BJ

40

3

9

7

ARX

7

7

16

15

ARMAX

6

1

16

11

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 1. Linear model structure for system identification

ACS Paragon Plus Environment

Page 28 of 37

Page 29 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Figure 2. SAGD process and variables used for proxy modeling

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 37

Figure 3. Proxy-modeling workflow

(a) (b) Figure 4. Input variable selection: model identification on training data. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment

Page 31 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(a) (b) Figure 5. Input variable selection: model validation using one-step ahead prediction. (a) Water production rate, (b) Bitumen production rate

(a) (b) Figure 6. Input data for system identification. (a) Steam injection rate, (b) Production pressure

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 37

(a) (b) Figure 7. Model structure selection: infinite horizon prediction on validation data set. (a) Water production rate, (b) Bitumen production rate

(a)

(b)

Figure 8. Autocorrelation of ARX model residuals on the identification data set. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment

Page 33 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(a)

(b)

Figure 9. One-step ahead predictions using Box-Jenkins (BJ) model based on two months data for early part of the reservoir’s life (year 1). (a) Water production rate, (b) Bitumen production rate

(a) (b) Figure 10. One-step ahead predictions using Box-Jenkins (BJ) model based on two months data starting from year 3 of the reservoir’s life. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 37

(a) (b) Figure 11. One-step ahead predictions using Box-Jenkins (BJ) model based on 6 months data. (a) Water production rate, (b) Bitumen production rate

(a) (b) Figure 12. One-step ahead predictions using Box-Jenkins (BJ) model based on 1 year data. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment

Page 35 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(a) (b) Figure 13. One-step ahead prediction with a recursive ARX model for 1 year data. (a) Water production rate, (b) Bitumen production rate

(a) (b) Figure 14. The adaptation of model coefficients with time for the recursive ARX model for 1 year data. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 37

(a) (b) Figure 15. 7-step ahead prediction with a recursive ARX model for 5 years data. (a) Water production rate, (b) Bitumen production rate

(a) (b) Figure 16. 7-step ahead prediction with a recursive ARX model with data preprocessing over a five year period. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment

Page 37 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

(a) (b) Figure 17. Performance of open loop simulation predictions of Box-Jenkins model with PRBS input design. (a) Water production rate, (b) Bitumen production rate

ACS Paragon Plus Environment