### Delft University of Technology

### Quantification of the impact of ensemble size on the quality of an ensemble gradient using

### principles of hypothesis testing

Fonseca, R. M.; Kahrobaei, S. S.; Van-Gastel, L. J T; Leeuwenburgh, O.; Jansen, J. D.

Publication date 2015

Document Version Final published version Published in

Society of Petroleum Engineers - SPE Reservoir Simulation Symposium 2015

Citation (APA)

Fonseca, R. M., Kahrobaei, S. S., Van-Gastel, L. J. T., Leeuwenburgh, O., & Jansen, J. D. (2015). Quantification of the impact of ensemble size on the quality of an ensemble gradient using principles of hypothesis testing. In Society of Petroleum Engineers - SPE Reservoir Simulation Symposium 2015 (Vol. 2, pp. 804-828). [SPE-173236-MS] Society of Petroleum Engineers.

Important note

To cite this publication, please use the final published version (if applicable). Please check the document version above.

Copyright

Other than for strictly personal use, it is not permitted to download, forward or distribute the text or part of it, without the consent of the author(s) and/or copyright holder(s), unless the work is under an open content license such as Creative Commons. Takedown policy

Please contact us and provide details if you believe this document breaches copyrights. We will remove access to the work immediately and investigate your claim.

This work is downloaded from Delft University of Technology.

**Green Open Access added to TU Delft Institutional Repository **

**Green Open Access added to TU Delft Institutional Repository**

**‘You share, we take care!’ – Taverne project **

**‘You share, we take care!’ – Taverne project**

**https://www.openaccess.nl/en/you-share-we-take-care **

**https://www.openaccess.nl/en/you-share-we-take-care**

### Otherwise as indicated in the copyright section: the publisher

### is the copyright holder of this work and the author uses the

### Dutch legislation to make this work public.

**SPE 173236-MS **

**Quantification of the Impact of Ensemble Size on the Quality of an Ensemble **

**Gradient Using Principles of Hypothesis Testing **

R.M. Fonseca, S.S. Kahrobaei, L.J.T.van Gastel, Delft University of Technology (TU Delft), O.Leeuwenburgh, TNO, and J.D. Jansen, TU Delft

Copyright 2015, Society of Petroleum Engineers

This paper was prepared for presentation at the SPE Reservoir Simulation Symposium held in Houston, Texas, USA, 23–25 February 2015.

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

**Abstract **

With an increase in the number of applications of ensemble optimization (EnOpt) for production optimization, the theoretical
understanding of the gradient quality has received little attention. An important factor that influences the quality of the
gradient estimate is the number of samples. In this study we use principles from statistical hypothesis testing to quantify the
number of samples needed to estimate an ensemble gradient that is comparable in quality to an accurate adjoint gradient. We
develop a methodology to estimate the necessary ensemble size to obtain an approximate gradient that is within a predefined
angle compared to the adjoint gradient, with a predefined statistical confidence. The method is first applied to the Rosenbrock
function (a standard optimization test problem), for a single realization, and subsequently for a case with uncertainty,
**represented by multiple realizations (robust optimization). The maximum allowed error applied in both experiments is a 10° **
angle between the directions of the EnOpt gradient and the exact gradient. For the single-realization case we need, depending
on the perturbation size, 900, 5 and 3 samples to estimate a “good” gradient with 95% confidence at 50 points in the
optimization space for 50 different random sequences. For the robust case, the conventional EnOpt approach is to couple one
model realization with one control sample, which leads to a computationally efficient technique to estimate a mean gradient.
However, our results show that in order to be 95% confident the original one-to-one model realization to control sample ratio
formulation is not sufficient. To achieve the required confidence requires a ratio of 1:1100, i.e. each model realization is paired
**with 1100 control samples using the original formulation. However, using a modified formulation we need a ratio of 1:10 to **
**stay within the maximum allowed error for 95% of the points in space, though a 1:1 ratio is sufficient for 85% of the points. **
We also tested our methodology on a reservoir case for deterministic and robust cases, where we observe similar trends in the
results. Our results provide insight into the necessary number of samples required for EnOpt, in particular for robust
optimization, to achieve a gradient comparable to an adjoint gradient.

**Introduction **

Multiple studies have shown the successful application of various optimization algorithms to maximize hydrocarbon recovery or net present value (NPV) over the producing life of a hydrocarbon reservoir. For such problems, gradient-based techniques, in terms of accuracy and computational efficiency, are the most successful and widely applied. The adjoint method provides the most accurate gradient and is computationally the most efficient method. However, the adjoint method requires access to a reservoir simulator's source code which for commercial simulators is practically impossible. Additionally it requires a considerable amount of time to implement and is not very flexible in adaptation to different control types. These limitations of the adjoint method have led to the development of alternative gradient-based techniques. One such alternative technique, Ensemble Optimization (EnOpt), inspired by the Ensemble Kalman Filter (EnKF) method was first introduced by Lorentzen et al. (2006) and Nwaozo (2006).

Chen (2008) proposed the now standard formulation of the EnOpt method which uses an ensemble of randomly perturbed control vectors to approximate a gradient of the objective function with respect to some specific controls. The major advantages of EnOpt are its ease of implementation, flexibility to adapt to different control types and ability to be used with any reservoir simulator. The major drawback of this method, relative to the adjoint method, is its computational inefficiency and inaccuracy of the gradient approximation. Nonetheless, recently many studies such as Chen et al. (2009), Chen and Oliver (2010), and Leeuwenburgh et al. (2010) have demostrated the applicability of EnOpt for large-scale production optimization

2 SPE 173236-MS

problems. Most of these papers have focused on deterministic optimization problems starting from a single reservoir model. However, in reality the geological and reservoir modeling process is fraught with uncertainties since a reservoir is modeled using uncertain interpretations based on uncertain data sources such as seismics, well logs etc. Incorporating these uncertainties into the optimization framework is vital to achieve results of any practical significance.

Van Essen et al. (2009) introduced a ‘robust optimization’ methodology in conjunction with the adjoint method to include the effect of uncertainties into the optimization framework. They used an ensemble of equi-probable reservoir models with differing geology and maximized the expectation of the objective function over this ensemble of models. Chen (2008) introduced this robust optimization concept within the ensemble optimization framework. They proposed the use of an ensemble of controls of equal size as the ensemble of geological models. Coupling of one member from the control ensemble with one member of the geological ensemble, a mean gradient can be approximated with the EnOpt formulation. This formulation, while computationally very attractive for robust optimization, has received scant attention with respect to its theoretical understanding. Recently Fonseca et al. (2014) demonstrated a case wherein the original formulation for ensemble-based robust optimization leads to inferior results and suggested a modified gradient formulation.

For EnOpt the two main inputs which influence the quality of the approximate gradient are the covariance matrix used to create the ensemble of perturbed controls and the number of control samples created, i.e. the ensemble size. The effect of the covariance matrix has been investigated recently in Fonseca et al. (2013) and a theoretical foundation for the use of a varying covariance matrix has been provided in Stordal et al. (2014). Sarma and Chen (2014) have investigated the applicability of different sampling techniques to improve the quality of a gradient estimate. However none of those studies have performed a detailed investigation into the effect of ensemble size on the estimated ensemble gradient quality.

In this study we aim to quantify the ensemble size required to approximate a gradient comparable to the adjoint gradient especially for robust optimization problems using principles from hypothesis testing and statistical analysis. In this paper we first provide an introduction of the hypothesis testing methodology and the different test statistics used. This will be followed by a detailed set of experiments on the widely used Rosenbrock function for cases with and without model uncertainty. Finally we test the proposed methodology on a medium-sized reservoir model, again with and without geological uncertainty

**Theory **

The two most commonly used objective functions for production optimization are ultimate recovery or an economic
*objective such as Net Present Value (NPV). In this work we chose the objective function J to be the NPV, defined in the usual *
fashion as

##

, , ,##

1 ( ) ( ) ( ) , (1 )*k*

*t*

*K*

_{o k}

_{o}

_{wp k}

_{wp}

_{wi k}

_{wi}

_{k}*t*

*k*

*q*

*r*

*q*

*r*

*q*

*r*

*t*

*b*

*J*

_{}

_{ }

_{}

_{ }

_{}

_{}

_{}

_{ }

_{}

_{}

(1)
*where qo,k is the oil production rate in bbl/day, qwp,k is the water production rate in bbl/day, qwi,k* is the water injection rate in
*bbl/day, ro is the price of oil produced in $/bbl, rwp is the cost of water produced in $/bbl, rwi* is the cost of water injected in
$/bbl, t*k is the difference between consecutive time steps in days, b is the discount factor expressed as a fraction per year, tk* is
*the cumulative time in days corresponding to time step k, and t* is the reference time period for discounting, typically one year.
Gradient-based optimization requires the gradient _{}(* _{dJ d}* )

*T*

**g** **u** which is used within an optimization algorithm to iteratively

optimize the objective function. For a detailed description of various available optimization algorithms see, e.g., Nocedal and Wright (2006). Usually the elements of the control vector are required to stay within upper and lower bounds, and different approaches for such bound control problems are available. Moreover, in addition to these constraints on the inputs, there may be constraints on the outputs of the simulator, which are much more difficult to handle. However, in this paper we are only considering the quality of the gradients under the presence of simple bound constraints.

**Ensemble Optimization (EnOpt) **

**Ensemble-Based Deterministic Formulation **

**In this section we outline the standard formulation of the EnOpt algorithm as proposed by Chen et al. (2009). We take u to **
*be a single control vector containing all the control variables to be optimized. This vector has length N equal to the product of *
the controllable well parameters (number of well settings like bottom hole pressures, rates or valve settings) and the number of
control time steps. Chen et al. (2009) sample the initial mean control vector from a Gaussian distribution while, at later
iteration steps the final control vector of the previous iteration is taken as the mean control. However the initial controls can
also be chosen by the user, as will be done in our experiments.

1 2 .
*T*
*N*
*u* *u* *u*
_{} _{}
**u** (2)

**To estimate the EnOpt gradient, a multivariate, Gaussian distributed ensemble {u**1**, u**2**, …, u***M*} is generated with a distribution
**mean u and a predefined distribution covariance matrix C where M is the ensemble size. During the iterative optimization ****process, u is updated until convergence, whereas C is, traditionally, kept constant. [An alternative procedure, in which C is **

SPE 173236-MS 3

updated during the optimization process, is decribed in Fonseca et al. (2013)]. In our implementation of EnOpt the ensemble
**members u***i, i = 1, 2, …, M, are created using *

1 2 _{,}
*i* *i*
**u** **u** **C z** (3)
with
1
1 *M*
*i*
*i*
*M*

###

**u**

**u**. (4)

We use a Cholesky decomposition to calculate **C , and draw z**1 2

*i* from a univariate Gaussian distribution. To estimate the
gradient, a mean-shifted ensemble matrix is defined as

###

1 2###

. * _{M}*

**U** **u** **u** **u** **u** **u** **u** (5)

**[Note that in earlier publications we used the transposed version of U. We modified our notation to bring it in line with that of **
textbooks such as Conn et al. (2009).] A mean-shifted objective function vector is defined as

1 2 ,
*T*
*M*
*J* *J* *J* *J* *J* *J*
_{} _{}
**j** (6)

where the expectation of the objective function is given by

1
1 *M* _{.}
*i*
*i*
*J* *J*
*M*

###

(7)The approximate gradient as proposed by Chen (2008) and Chen et al. (2009) is given by

1 _{,}
*uu uJ*
**g****C c** (8)
where
1
( )
1
*T*
*uu*
*M*
**C** **UU** (9)
and
1
( )
1
*uJ*
*M*
**c** **Uj** (10)

**are ensemble (sample) covariance and cross-covariance matrices respectively. (Note that c***uJ* is a one-dimensional matrix, i.e. a
* vector.) For the usual case where M < N, matrix Cuu* is rank-deficient, and Chen (2008) and Chen et al. (2009) therefore
propose not to use expression (8) but, instead, to use

1 _{= ,}
*uu* *uu uJ* *uJ*
**g** **C C c** **c** (11)
or
.
*uu uJ*
**g** **C c** (12)

**Alternatively, the pre-multiplication in equation (12) can be performed with C , leading to **
.

*uJ*

**g** **Cc** (13)

All three expressions (11), (12) and (13) can be interpreted as modified (regularized or smoothed) approximate gradients. In the present paper we use a straight gradient, i.e. expression (8), computed as the underdetermined least squares solution

†
( *T*)

**g** **UU** **Uj** , (14)

where the superscript † indicates the Moore-Penrose pseudo inverse, which is conveniently computed using a singular value decomposition (SVD); see, e.g., Strang (2006). Moreover, we use smoothed and double-smoothed versions of equation (14):

†
( )
*T*
**g** **C UU** **Uj** , (15)
†
( )
*T*
**g** **CC UU** **Uj** , (16)

Equation (14) was also described in Dehdari and Oliver (2012), while Do and Reynolds (2013) recently demonstrated that it is akin to what is known as a ‘Simplex gradient’ in, e.g., Conn et al. (2009). Do and Reynolds (2013) also provided theoretical connections between various ensemble methods such as simulataneous perturbation stochastic approximation (SPSA),

4 SPE 173236-MS

Simplex gradient, EnOpt etc. Moreover, they proposed a modification to the gradient formulation which uses the current
control vector ** _{u}**

_{ and the corresponding objective function value }

**

_{J}_{ to calculate the control and objective function anomalies }

**U and j: **
1 2 ,
_{} _{} _{}
*M*
**U** **u** **u** **u** **u** **u** **u** (17)
1 2 ,
_{} _{}
_{M}*T*
*J* *J* *J* *J* *J* *J*
**j** (18)

*where the superscript ℓ is the optimization iteration counter. *

Equations (11-16) can all be used to estimate a gradient-based on either the original [equations (5) and (6)] or the modified [equations (17) and (18)] formulations. Thus we can estimate as many as twelve different gradient formulations for deterministic cases. Further varieties will emege when considering robust optimization.

**Ensemble-Based Robust Formulation **

Chen (2008) proposed a computationally efficient robust ensemble algorithm in which she used an ensemble of controls of
equal size as the ensemble of geological models. Coupling one member from the control ensemble with one member of the
*geological ensemble a mean gradient is approximated; see Chen (2008) for details. Therefore the ensemble size M of the *
*controls is equal to the number of geological models, i.e. a 1:1 ratio. Hence only M simulation runs are needed to approximate *
the ‘robust’ gradient of the objective function. Recently Stordal et al. (2014) reached a similar conclusion starting from a
different mathematical viewpoint. However the theoretical understanding of using this 1:1 ratio is still incomplete. As an
alternative to this formulation, Fonseca et al. (2014) propose a modified formulation for the robust gradient which no longer
uses the mean-shifted control samples and objective values, equations (5) and (6). Instead, in equation (5) the control sample
mean **u is replaced by the control vector of the current iteration step, u***ℓ*:

1 2 ,

_{} _{} _{}

*M*

**U** **u** **u** **u** **u** **u** **u** (19)

The new formulation replacing equation (6) is

1 1 2 2 ,
_{} _{} _{}*T*
*M* *M*
*J* *J* *J* *J* *J* *J*
**j** (20)

Note that equation (19) is identical to equation (17) as used in the deterministic modified expression of Do and Reynolds (2013), but that equation (20) is different from equation (18). This modified gradient formulation [based on equations (19) and (20)] will also be tested in our set of experiments. It behaves distinctly different compared to the original robust formulation [based on equations (5) and (6)]. First, because the subtractions in the objective function values in equation (20) are with respect to the individual objective function values

*i*

*J* and not with respect to the mean. Second, because for
bound-constrained control problems, ** _{u and }**

_{u}_{ may be shifted with respect to each other. Thirdly, because the effect of outliers, }

which may strongly influence the mean value our least-squares approach to estimate the gradient, is reduced in the modified formulation.

Note: all the different gradient formulations (11-16) for deterministic optimization are also applicable to the robust case. Together with the robust modified formulation [equations (19) and (20)] this leads to a total of 18 potential robust gradient formulations for the 1:1 ratio (i.e. one control perturbation for each geological realization) approach. However, another distinction can be made if we use other ratios. E.g., Raniolo et al. (2013) suggest the use of 20 control perturbations for every model realization. For every model realization, using the 1:20 ratio, they estimate an individual gradient, whereafter they take the mean of the individual gradients to obtain the robust gradient. This formulation will hereafter be referred to as the ‘Mean of Individual Gradients’ (MIG). Alternatively, one can combine all the controls and objective function anomalies to estimate a single robust gradient, i.e. not estimate individual gradients for every model realization. This approach will hereafter be referred to as the ‘Hotch-Potch Gradient’ (HPG). This additional disctinction leads to a total of 30 potential formulations [2 times 18 minus 6 because for the MIG approach there is no difference between using equations (18) and (20)].

**Adjoint method **

The adjoint method has been investigated extensively for use in data assimilation and production optimization. Detailed derivations for the production optimization case can be found in, e.g., Brouwer and Jansen (2004), Sarma et al. (2005), Kraaijevanger et al. (2007) and Jansen (2011). The adjoint method is the most accurate and computationally efficient method for computing a gradient. Computation of the gradient only requires one forward simulation and one fast backward computation. Therefore the number of simulation runs is independent of the number of controls. However for robust optimization using the adjoint requires running the forward and backward simulation for every geological realization, thus to compute the robust gradient, the same number of simulation runs will be performed as required for the robust EnOpt gradient using the 1:1 ratio. For our experiments we assume that the adjoint gradient is the exact gradient, which the EnOpt method

SPE 173236-MS 5

tries to approximate. In this study the adjoint module available in the Shell in-house simulator was used (Kraaijevanger et al. 2007).

**Hypothesis testing **

We use principles from hypothesis testing to validate the research goal of this paper, namely to test if the approximate EnOpt gradient is comparable in quality to the adjoint gradient. To be able to determine the difference in gradients we compute the angle between them by using the dot product:

.
cos( ) *adj* *ens* .

*adj* *ens*

**g** **g**

**g** **g** (21)

Another measure that describes the difference in direction of two vectors is the length of the difference between the normalized adjoint and EnOpt gradients, defined as

.
*adj* *ens*
*ens*
*adj*
**g** **g**
**g**
**g** (22)

To eliminate the effect of a difference in gradient magnitude between both methods the gradient vectors are made into unit
vectors by dividing them by their norm. When goes to zero, the two gradients will point in the same direction, just like if the
angle goes to zero. These two test parameters can be used to cross-validate each other, because a difference of 10◦_{ corresponds }
to a dimensionless length difference of 0.175. Thus the two equivalent null hypotheses used are

0: 10 ,
: 0.175 .
*H*
(23)
The statistical inference method used is based on pre-defined confidence intervals for the testing parameters defined above.

**Confidence intervals **

Creating a confidence interval is a method to define a range at and the certainty that the true value of an estimated parameter lies within it, based on the knowledge of the sampling distribution (Dekking et al., 2005). In our numerical experiments we create a dataset of our parameters, given in equations (21) and (22). The parameter of interest is the maximum allowable deviation of the EnOpt gradient with regard to the exact or adjoint gradient. As it is virtually impossible to achieve a 100% confidence, we apply a confidence level of = 0.95. The general definition of a confidence interval assumes a two-sided interval, i.e. an upper and a lower limit. However it is also possible to have a one-sided interval. As the test parameters used for the numerical experiments are absolute values of deviations we only want to find the confidence interval of the maximum deviation, thus the upper limit. Using a one-sided interval, the confidence interval is just the integral of the probability distribution, i.e. the cumulative density function (CDF).

**Beta distribution **

The EnOpt method samples random points from a normal distribution with a user-defined standard deviation. The
distribution of the test parameters ( and ) are, however, not normally distributed due to the non-linearity of the objective
function and the function for the test statistics. In order to determine a confidence interval the distribution of the underlying
parameters needs to be known. The two parameters and are , by definition, ratios of two different gradients, and the beta
distribution is suitable to fit data sets that are ratios. The beta distribution forms a class of continuous probability distributions,
*parameterized by two shape parameters a and b that define the shape of the distribution. AbouRizk et al. (1994) outline several *
methods to determine these shape parameters, of which we have chosen to use the maximum likelihood estimator. AbouRizk
et al. (1994) also demonstrate that there is virtually no difference in the results when using different fitting methods. Note that
the beta distribution is always bounded in the interval [0,1]. However, in this study, because we use the cosine function, our
data is bounded between [-1 1]. Thus in order to use the beta distribution for angles greater than 90 degrees (i.e. cos() < 0) we
simply reset the angle to 90 degrees (i.e. cos() =0).

**Methodology using Beta Distributions **

The methodology proposed in this work to either accept or falsify our null hypothesis and estimate the necessary ensemble size to accept our null hypothesis is as follows:

Sample points in the control space and compute the adjoint and EnOpt gradients at each point.
* Estimate the fitting parameters a and b of the beta distribution using the data. *

Plot the CDF for the beta distribution and compute the confidence interval for the predefined error. Repeat for varying ensemble sizes.

6 SPE 173236-MS

We have chosen the 95% confidence interval to test our methodology, however any different confidence interval can be chosen within the same workflow. Varying the desired confidence interval will automatically vary the ensemble size needed to accept our hypothesis. In essence an ensemble size is found that gives an accurate gradient approximation in a number of times equal to the confidence interval, i.e. if the numerical experiment would be repated many times, 95% of the experiments would result in an approximate gradient within the error margin compared to the true (adjoint) gradient. For a more detailed explanation of confidence intervals see, e.g. Dekking et al. (2005).

**Methodology with Traditional Hypothesis Testing Principles **

In addition to using the beta distribution to test our null hypothesis we can also use a traditional hypothesis testing approach to either accept or falsify our hypothesis and estimate the necessary ensemble size as follows:

Sample points in the control space and compute the adjoint and EnOpt gradients at each point. Count the number of points that satisfy the null hypothesis

Compute the confidence interval for the predefined maximum allowable error. Repeat for varying ensemble sizes until the confidence interval is achieved.

In this paper we compare the results of the two methodologies for the different test statistics in the following section.

**Numerical Example **
**Rosenbrock Function **

The methodology is first applied to the non-linear Rosenbrock function named after the mathematician who first used it to demonstrate his optimization algorithm, see Rosenbrock (1960). This analytical function has since been used as a standard test case in mathematical optimization. The Rosenbrock function consists of a curved narrow valley which most algorithms have little difficulty finding. However once found, the difficulty lies in finding the global optimum which is situated inside the valley. Equation (24) is a slightly altered version, as it is multiplied by -1, making it a maximization problem opposed to a minimization problem:

2 2 2

1 2 2 1 1

( , ) 100( ) (1 )

*J u u* *u* *u* *u* . (24)

*The Rosenbrock function has an optimal solution J=0 at (u*1*,u*2**) = [1 1]. This point lies on a long curved ridge; see Fig. 1 for a **
contour plot. Since it is an analytical function it is possible to compute the exact gradient. The red dots in Fig. 1 are 50
randomly distributed points in space which will be used to test our methodology. All the numerical experiments will be carried
out over the same set of points for different scenarios. Since we are working with approximate gradient techniques the effect of
different random sequences also needs to be accounted for. Therefore, all the experiments are repeated for 50 different random
sequences. The results presented below are the mean values for the 50 points in space over the 50 different random sequences.

**Fig. 1: Contour plot of Rosenbrock function given by equation (24). Red dots are 50 points randomly distributed in space. **

**Deterministic Case **

We first test the two hypothesis testing methodologies which consist of calculating the angle between the gradient estimated by EnOpt and the exact gradient for different ensemble sizes. The points are uniformly distributed so as to capture the effect of the spatial variability in the objective space on the gradient quality, with many points that are on the ridge or on

*u*
*1*
*u2*
-2 -1.5 -1 -0.5 0 0.5 1 1.5 2
-1
-0.5
0
0.5
1
1.5
2
2.5
3

SP
the
be
ens
me
ens
sm
is a
est
rec
gra
at a
**Fig**
**and**
**(blu**
ens
var
the
inc
allo
the
con
in
38
ens
giv
(le
E 173236-MS
e edge of the r
comparable
semble gradie
ean angle ove
semble size fo
maller ensembl
a two-control
timate a ‘good
confirmation t
adient. For the
all the 50 poin

**g. 2: Illustration**
**d varying pertu**
**ue dots). **
One of the m
**semble size. F**
rying perturba
e spread in th
creasing ensem
owable differ
e different en
nfidence inter
this paper. As
of the 50 poi
semble size o
ve almost the s
ftmost blue cu
ridge. This wi
to the true gr
ent is the pert
er the 50 spati
or varying per
le sizes while
problem, thus
d’ gradient w
that EnOpt, f
e largest value
nts in space.
**n of the trend **
**urbation sizes **
methodologies

**Fig. 3 and Fig**

ation sizes. W
he angles at th
mble sizes the
ence between
nsemble size.
rval of 78% w
s a compariso
ints satisfy ou
f 300 and a p
same results w
urve in Fig. 3)
ill to some de
radient almos
turbation size
ial points and
rturbation size
for an increa
s we would ne
with the finite
for a small pe
e of used we
**in the mean a**
**(****). Lower va**
s we propose i
**g. 4 are the cu**
While Fig. 2 re
he different p
e curves shift
n the two grad
We see that
which is much
n, if we use th
ur null hypoth
perturbation si
which is a cro
) we cannot ac
gree ensure th
st everywhere
used to gene
d 50 random s
es . We obse
asing ensemble
eed a maximu
difference m
erturbation, al
e observe that

**ngle over the 5**
**lues (green & **

in this paper u umulative den epresents a sin points for vary

towards an an dients i.e. 10 d

t even for an less than the he second met hesis. Thus 38

ize () equal t ss validation o chieve good re

hat once the r e in space. Th erate the ense

seeds versus erve that, as ex

e size the gra um of three (i. method provide

lso would req t even for an e

**50 points in sp**
**red dots) of **

uses beta distr nsity function ngle mean val ying perturba ngle equal to degrees. We c n ensemble si 95% confiden thodology, i.e 8/50, i.e. 76% to 0.1, we hav of the beta dis esults and this

ight ensemble
he other impo
emble of contr
an increasing
xpected, as
dient quality
*e. n+1, where*
ed the perturb
quire an ense
ensemble size
**pace and 50 ra**
** give significa**
ributions to re
of the beta di
ue for a partic
tion and ense

zero degrees
an read-off th
ize of 300 (r
nce interval w
e. a traditional
of the points
ve a confiden
stribution appr
s is also the ca
e size is determ
ortant factor
rols for the g
ensemble siz
decreases the
always impro
*e n is the prob*
bation size is
emble size equ

of 300 we do

**andom seeds f**
**ntly better res**

epresent the ca istributions fo

cular ensembl emble sizes. W

and the black he achieved co ightmost ligh which would sa l hypothesis te , satisfy our h ce interval of roach. Note: f ase for the diff

mined the En in determinin gradient estim

ze which show
e gradient qua
oves. The Ros
blem dimensio
sufficiently s
qual to 3 to e
o not achieve
**for an increasi**
**sults compared**
alculated angl
or different en
le size, Fig. 3
We observe in
k line is the d
onfidence inte
ht green curv

atisfy the null
esting approac
hypothesis, an
f 76%. The tw
for an ensemb
fferent perturb
Opt gradient w
ng a high-qua
**ate. Fig. 2 is **
ws the impac
ality increases
enbrock funct
on) simulation
small. Fig. 2
estimate a ‘go
a ‘good’ grad
**ing ensemble s**
**d to higher val**
les to estimate
nsemble sizes
and Fig. 4 sh
n Fig. 3, that
efined maxim
*erval as 1-Pcdf*
e) we achiev
l hypothesis u
ch, we count t
nd hence, with
wo methodolog
ble size equal t
ation sizes.
7
will
ality
the
t of
for
tion
ns to
is a
ood’
ient
**size **
**lues **
e an
and
how
for
mum
*f* for
ve a
used
that
h an
gies
to 2

8
**Fig**
sim
No
sam
=
in
sid
wit
ord
Th
**Fig**
**poi**
ens
lar
ens
fin
int
Th
unc
bet
**g. 3: Illustration**
The same ex
milar results, i
ote: this result
me results. To
= 0.1. Fig. 4 c
Fig. 3, the dis
de plot) shows
th an ensembl
der of magnitu
hese results are

**g. 4: (a) Cumm**
**ints. (b) Cumul**
The above fi
semble size n
ger ensemble
semble size of
nite difference
* Rosenbrock *
The quality
roduced throu
he two constan
certainty intro
tter understan

**n of the cumula**xercise is rep i.e. similar con t is not surpris o achieve the c consists of two stribution of th s the distributi le size of 5 sam ude smaller, w e also cross-va

**mulative densit**

**lative density f**igures show th needed to achi sizes to achi f 3, i.e. for an gradient.

*of the ‘robust ugh equation (*

**Function Wi***nts c*1

*and c*2 a oduced. This nd these formu

**ative density fu**eated using th nfidence inter sing because e confidence in o plots for dif he angles calc ion of the calc mples (Fig. 4a we achieve a

alidated using

**y functions fo**
**functions for **

hat the gradie ieve good qua

eve the 95%
*n n dimensiona*
* ith Uncertaint*
t’ EnOpt grad
(25). The unce
(

*J u*are drawn from

provides us a
ulations. The
**unctions of the**
he norm in
rvals, compare
essentially cal
nterval define
fferent perturb
culated for the
culated angles
a, purple curv
100% confide
g the traditiona
**or **** = 0.01 whe**
** = 0.001, where**
ent quality is s
ality gradients
confidence in
al problem we
**ty **

dient has not b ertainty is intr

1, ) = 100(2

*u u*

m a Normal d a fast and accu

*value for c*1 m
**e beta distribut**
nstead of the
ed to the resu
lculating an an
ed in this work
bation sizes. F
e 50 points usi
s for the same
ve) we achieve
ence interval
al hypothesis

**ere for an ens**
**e for an ensem**
strongly affec
s. We observe
nterval. For si
e need an ense
been investiga
roduced to mim
2 2
1 2 1
(*c u* *u* )
distribution, w
urate way to
mainly affects

**tion for varying**

angle as the ults shown abo

ngle or the
k we would ne
ig. 4a (left-sid
ing different e
e 50 points wi
e a 100% conf
with an ensem
testing metho
**semble size of**
**mble size of 3 a**

ted by the per e that, for larg ignificantly sm emble size equ

ated before, so mic geologica

2

sin( ) (1*c*

where the stand test the vario s the steepness

**g ensemble siz**

data. We obs
ove with exac
norm with the
eed approxim
de plot) depict
ensemble sizes
ith = 0.001.
fidence interva
mble size of 3
d as well as u
**5 or higher th**
**nd higher the h**
rturbation size
ge perturbatio
maller perturb
*ual to n+1, jus*
o uncertainty
al uncertainty.
2
1) .
*u*
dard deviation
us ‘robust’ En
s of the object
**zes for **** = 0.1. **
serve (results
ctly the same
he same data s
mately 900 sam
ts, similar to t
s with = 0.0
. We observe
al while for
3 samples (Fig
using the nor

**he hypothesis**
**hypothesis is s**
e which in tur
on sizes, we n
bation sizes w
st like if we w
in the Rosenb
.
n reflects the
EnOpt gradien
tive space, wh
SPE
173236-not shown) v
trend in the d
should lead to
mples when us
the curves sho
01. Fig. 4b (rig
that for = 0
= 0.001, i.e.
g. 4b, red curv
rm.
** is satisfied at**
**satisfied. **
rn determines
need significan
we would need
were to estima
brock function
(2
magnitude of
t formulation
*hile c*2 rotates
-MS
very
data.
the
sing
own
ght-0.01
one
ve).
**t all **
the
ntly
d an
ate a
n is
25)
f the
s to
the

SP
spa
Th
**Fig**
**unc**
**Fig**
**unc**
In
wit
Ch
ma
ens
exp
we
ens
me
sat
on
sig
hav
**Fig**
**(blu**
Th
we
me
wh
E 173236-MS
**ace. Fig. 5 an**
his was done to

**g. 5: Contour p**
**certainty. **
**g. 6: Contour p**
**certainty. **
view of the r
th the tradition
hen (2008) fir
athematical re
sembles for bo
periments the
ere reasoned t
semble size o
ean angle sign
tisfy our null h
equations (19
gnificantly bet

ve satisfied ou

**g. 7: Comparis**
**ue) and modifi**

hese results ar e usually use a ean angle of 4 hen using the 1

**nd Fig. 6 show**
o investigate t

**plot of five rea**

**plot of five rea**

results from th nal hypothesis rst proposed

asoning for th oth the model e original 1:1 to be the lack f model realiz nificantly redu hypothesis wh 9) and (20) (re tter. We obser ut null hypoth

**son of the tren**
**ed form (red) u**
e based on th
an ensemble
40 degrees wh
1:1 ratio for ro
w five realizat
the impact of u
**alizations of th**
**alizations of th**
he determinis
s approach to
the idea of th
he applicabilit
ls and consequ
formulation i
k of a good qu
zations and, c
uces for an in
hen using the o
ed curve) whi
rve that for an
hesis and achie

**nd in the mean**
**using the comp**

he case represe size of 100 m hile the modifi obust optimiza

tions each for uncertainty on

**he Rosenbroc**

**e Rosenbrock**

tic case, for t
determine the
he 1:1 ratio e
ty of this 1:1 r
uently the con
in many insta
uality gradien
consequently,
ncreasing ens
original formu
ch still retains
n ensemble siz
eved the 95%
**n angle over th**
**putationally att**

enting the hig models to capt fied formulatio ation, the mod

cases which n the gradient

**k function giv**

**function give**

the robust cas
e ensemble siz
explained in
ratio provided
ntrols. In realit
ances indeed p
**nt. Fig. 7 illus**
also control r
semble size. H
ulation (blue c
s the computa
ze of 1000 w
confidence in
**he 50 spatial p**
**tractive 1:1 rat**
ghest uncertain
ture the uncer
on achieves a
dified formula
are representa
estimate.
**ven by equatio**
**n by equation **
se we only us
ze required to
the theory se
d by Chen (20
ty, however, w
produces resu
strates the eff
realizations. W
However, even
curve). On the
ational attracti
we have a mea
nterval.

**points with inc**
**io. **
nty in the mo
rtainty. For th
mean angle o
ation based on
ative of low a
**on (25) for a ca**
**(25) for a cas**
e the angles a
estimate a go
ection to estim
08) is only ap
we usually wo
ults of insuffic
fect of using t
We observe, in
n for an ense
e other hand th
iveness of the
n angle of app
**creasing ensem**
dels. In our re
hat size, the or
of 7 degrees. W
n equations (19

and high unce

**case mimicking**

**se mimicking a**

as our testing ood quality ‘ro

mate a ‘robus pplicable if we ork with 100 m cient quality. the 1:1 ratio f n line with th emble size of he modified fo original form proximately 4

**mble sizes for**

reservoir optim
riginal formu
We recommen
9) and (20) be
rtainty scenar
**g a low degree**
**a higher degree**
parameter al
obust’ gradien
st’ gradient. T
e have very la
models and in
The poor res
for an increas
e theory, that
1000 we do
formulation ba
mulation perfor
4 degrees and
**r the original fo**
mization prob
ulation achieve
nd therefore t
e applied.
9
rios.
**e of **
**e of **
ong
nt.
The
arge
our
sults
sing
the
not
ased
rms
d we
**orm **
lem
es a
that,

10
stu
the
mo
for
aga
ori
sca
rea
bet
com
.
**Fig**
**for**
des
Fig
the
all
the
rea
sho
mo
of
me
the
* Effect of Hig*
The results i
udies, e.g. Ran
e use of higher
odel and cont
rmulation usin
ainst the ense
iginal formula
ale). In both p
alizations. Onc
tter gradient e
mputational c

**g. 8: (a). Illustr**

**rmulation. (b) Il**

*The effect of sired confiden g. 9(a) (left-si e lowest uncer the results thu e modified fo alizations was*

**Effect of Un**own) that for
odified formul
uncertainty, t
ean angle less
e necessary co
* gher Ratios *
in Fig. 7 were
niolo et al. (20
r ratios (1:20
trol perturbat
ng a perturbat
emble size of
ation, and Fig.
plots we observ

ce again, the m estimates can

osts of using t

**rates the effec**
**llustrates the e**

**ncertainty **

f uncertainty i
**nce interval. F**
ide plot) depic
rtainty scenar
us far, that the
ormulation gr
s kept constan
r the original

lation the grad the original fo than 10 degr onfidence inter

e for the case
013) and Li et
etc.) to find b
tions on the
tion size equa
model realiza
8(b) which d
ve that an incr
modified form
be achieved w
these ratios m
**ct of using hig**
**effect when usi**

is investigated
**Fig. 9 is an ill**
cts the results
io for both th
e modified for
radient estima
nt at 100 for t
formulation
dient estimate
ormulation w
rees. Note: if t
rval.
representativ
al. (2013), did
etter gradient
quality of th
**al to 0.01. Fi**
ation with the
displays the res

reased ratio g mulation outpe while using hi must also be ac

**gher ratios on**
**ing the modifie**

d in conjunctio lustration of t

for the highe e original and rmulation resu ate is less se his exercise a using a large e improved for ould need a r the mean ang

ve of the highe
d not achieve
estimates. W
he gradient e
**ig. 8 consists **
e curves repre
sults using the
gives better gra
erforms the o
igher ratios. N
ccounted for, e
**n the gradient **
**ed formulation **
on with the ra
the impact tha
est uncertainty
d modified for
ults in signific
ensitive to th
and a perturba
er perturbatio
r smaller pertu
ratio of 1:20
gle is less than

est model unc
results of pra
e investigate h
estimate. Note
of two plots,
esenting the ra
e modified for
adient estimat
riginal formul
Note though th
especially for
**quality when **
**with the HPG f**
atios needed to
at uncertainty
y case while F
rmulations. In
cantly better g
he effect of u
ation size of 0
on size resulte
urbation sizes
or 1:50 i.e. 2
n 10 degrees i
certainty using
ctical value w
here the impac
e: the results
Fig. 8(a) wh
atio used for t
rmulation (not
tes for the diff
lation also for
hat while it is
large-scale hi

**using the orig**
**formulation. N**
o satisfy the n
has on the qu
**Fig. 9(b) (righ**
n both plots w
gradients than
uncertainty. T
0.01 was used
ed in better m
s. We observe
000 or 5000
it does not gua

g the 1:1 ratio
with the 1:1 rat
ct of varying t
are obtained
hich is the me
the gradient e
te the differen
ferent ensemb
r larger ratios
better to use
igh-dimension
**ginal formulat**
**ote the differe**

null hypothesi uality of the g ht-side plot) d we observe, in n the original f The ensemble d. We also no mean angles, e that, dependi function eval uarantee that w SPE 173236-o. Some previ tio and sugges the ratio betw d with the H ean angle plot estimate with nce in the vert ble sizes of mo

. Thus in gen higher ratios, nal problems.

**ion with the H**
**nt vertical scal**
s and achieve
gradient estim
depicts results
accordance w
formulation. A
e size of mo
ticed (results
whereas for
ing on the deg
luations to fin
we have achie
-MS
ious
sted
ween
HPG
tted
the
tical
odel
eral
the
**HPG **
**les. **
the
mate.
for
with
Also
odel
not
the
gree
nd a
eved

SP
**Fig**
**Illu**
**ang**
con
tha
ori
als
for
for
**Fig**
**unc**
**hyp**
we
dim
‘ro
the
pro
Fig
rat
dif
the
the
E 173236-MS
**g. 9: (a) shows**
**ustrates the res**

**gles for the hig**

**Hypothesis T****Fig. 10 illus**
nfidence inter
at for a smalle
iginal formula
so observe tha
r the highest u
rmulation for a
**g. 10: Hypothes**
**certainty case **
**pothesis and a**
* Mean of Ind*
The results p
e would requi
mensional con
obust’ gradien
e adjoint form
oblem we shou
g. 2. Using th
io of 1:5 to e
fference betwe
e ratio needed
e same trend is

**s the effect of**

**sults for the lo**

**gh uncertainty**

**Test Results **

strates the ra rval for the di er ensemble s ation while fo at the required uncertainty sce

an ensemble s

**sis testing resu**
**with the HPG**
**achieve the des**

* dividual Grad*
presented in F
re ratios as h
ntrol problem
nt when using
mulation and al
uld need 3 fun
he MIG formu
estimate a me
een the origin
compared to
s observed for

**a high degree**
**owest uncertain**

**case are highe**

atio necessary fferent gradie size of model or the same en d ratio decreas enario. For th size of 10 mod

**ults using the **
**G formulation. **
**sired confidenc**
* dients (MIG) *
Figs. 8-10, wh
high as 1:1100
these results
a ratio other t
lso by Raniol
nction evaluat
ulation we obs
ean angle less

al and modifie the results wi r the low unce

**e of uncertaint**
**nty case. Note**
**er than the low**

y to satisfy th
ent formulation
l realizations,
nsemble size,
ses with an in
e lowest unce
del realization
**traditional met**
**A lower ratio **
**ce interval. **

ich are obtain
0 to satisfy th
are complete
than 1:1. Note
o et al. (2013
tions to estima
**serve in Fig. **
s than 10 degr
ed formulation
ith the HPG fo
ertainty case.
**ty on the mea**
**e: for both thes**
**w uncertainty ca**

he null hypot
ns and the dif
e.g. 10, we n
we would re
ncrease in the
ertainty scenar
ns, and a 1:5 ra
**thodology for **
**is needed for **

ned using the H
he hypothesis
ely counter-in
e: This approa
) using EnOp
ate a gradient,
**11 that for a **
rees. Increasi
n is minor at b
ormulation. T
**n angle for in**
**se plots an ens**

**ase. (Note the **

thesis used in fferent ensemb need significa quire a ratio ensemble size rio we would

atio for the mo

**the original (b**
**higher ensem**

HPG formula and achieve ntuitive. Thus

ach was follow
pt gradients. S
, for a sufficie
perturbation
ng the ratio w
best. These re
he results are
**creasing ratios**
**semble of 100 **
**different vertic**
n this paper a
ble sizes of m
antly higher ra
of 1:10 with t
e of model rea
require lower
odified formu

**lue) and modif**
**mble sizes of m**
tion, illustrate
the desired co
we test the M
wed by Van E
ince we are d
ently small per
size equal to
will improve t
esults are a ma
based on the
**s with the HPG**
**model realizat**
**cal scales). **
and achieve
model realizati
atios (1:1100)
the modified
alizations. All
r ratios: 1:800
ulation.

**fied (red) form **
**model realizatio**

e that for sma
onfidence inte
MIG formulat
Essen et al. (2
dealing with a
rturbation size
0.01 we wou
the gradient e
arked improve
highest uncer
**G formulation.**
**tions is used. **
the desired 9
ions. We obse
) when using
formulation.
l these results

with the orig

**using the high**
**ons to satisfy**
ll ensemble si
erval. For a tw
ion to estimat
2009) albeit us
two dimensio
e, as illustrated
uld at best nee

estimate, and
ement in term
rtainty case wh
11
**. (b) **
**The **
95%
erve
the
We
are
inal
**hest **
**the **
izes
wo-te a
sing
onal
d in
ed a
the
s of
hile

12
**Fig**
**req**
det
suf
con
the
**Fig**
**wo**
In
sm
**3D**
illu
we
sys
act
flu
is a
the
is
inj
En
**g. 11: Illustrati**
**quired even for**

The effect o
terministic ca
fficiently sma
nfidence inter
e ratio increase
**g. 12: Hypothes**
**ould require sig**

summary, for mall ensembles

**D Reservoir M**

The ‘Egg M ustration of th ells (blue) and

stem in the fo
tive. A detaile
uid properties u
assumed to be
e injectors are
simulated for
ection rates as
We use a com
nOpt gradient
**on of the mea**
**r small ensemb**
of perturbatio
ase, that a dec
all perturbation
rval even for

es, although, i

**sis testing res**
**gnificantly sma**

r this case usin s of model rea

**Model: “Egg M**

odel’, first int
he permeabilit
d four produc
orm of discre
ed description
used for all en
e incompressi
rate-controlle
a period of 3
s controls we
mmercial fully
*of J with re*
**an angle for in**
**ble sizes comp**

n size estima creasing pertu n size ( = 0 small ensemb it is significan

**sults using the**
**aller ratio comp**

ng the MIG fo
alizations.
**Model” **
troduced by v
ty field of a si
ction wells (r
ete permeabili
n of a standard
nsemble mem
ible. The botto
ed between 1
3600 days or
have of 40x8
y implicit blac
espect to the
**ncreasing ens**
**pared to Fig (10**

ates from the urbation size 0.001) a ratio ble sizes of m ntly lower than

**e traditional me**
**pared to the HP**
ormulation giv
van Essen (200
ingle model r
ed) completed
ity fields mod
dized version
mbers are given
om hole press
and 79.5 m3_{/d}
slightly less t
= 320 control
ck oil simulat
**controls u. A**
**semble sizes a**
**0) for the highe**

MIG formul results in an of 1:3 is suffi model realizati n the ratio req

**ethodology for**
**PG formulation**

ves much bett

09), is a chan
realization wit
d in all the l
deled with 60
of this Egg M
n in Table 1. N
sures of the pr
day. The initia
than 10 years
ls, i.e. a 320 d
tor (Eclipse, 2
An in-house
**and ratios. We**
**est uncertainty **
lation is show
increase in th
ficient to satisf
ons. We obse
quired when us
**r different pert**
**n. **
er angles for s
nelized reserv
th the location
ayers. The m
0 × 60 × 7 = 2
Model is given
No capillary p
roducers are c
al reservoir pr
s. There are 4
dimensional pr
011) for the r
simulator is u
**e observe a si**
**case. **
**wn in Fig. 12**
he quality of
fy our null hy
erve that for a
sing the HPG

**turbation sizes**

smaller ratios

voir with seve
ns of the eight
model represen
5.200 grid ce
n in Jansen et
pressures are in
constrained be
essure is at 40
0 control time
roblem.
eservoir simu
used to calcu
**ignificant redu**
**2. We observ**
the gradient
ypothesis and
an increasing
formulation.
**s (****) and ratios**
especially wh
en vertical lay
ht (mainly peri
nts a channeli
ells of which
t al. (2013). T
included and t
etween 385 an
00 bars. Produ
me steps of 90
ulations to esti

ulate the adjo

SPE
**173236-uction in the r**
ve, ‘akin’ to
estimated. Fo
achieve the 9
perturbation s
**s. Lower **** val**
hen working w
**yers. Fig. 13 is**
ipheral) inject
ized depositio
18.553 cells
The reservoir
the reservoir r
nd 400 bar, wh
uction of the f
days, thus us
mate the (robu
oint gradient,
-MS
**ratio **
the
or a
95%
size
**lues **
with
s an
tion
onal
are
and
rock
hile
field
sing
ust)
see

SP
Kr
can
**Fig**
det
spa
cov
‘re
dif
to
con
con
sig
pla
dif
fun
**Fig**
**(ind**
ang
wit
ob
E 173236-MS
raaijevanger et
n be found in
**g. 13: Reservoi**
* Deterministi*
The methodo
terministic ca
ace. Randoml
vered all featu
elevant’ points
fferent initial s
79.5 m3

_{/day. }ntrol points, a ntrol trajector gnificant scope ateau in the ob fferent parame nction value a

**g. 14: Illustratio**

**dicative of the**

*The Steep Re*

**Fig. 15 disp**

gle for an inc
th different ra
serve that we
t al. (2007) fo
Jansen et al. (
**r model displa**
* ic case *
ology present
ase. When wo
ly creating po
ures of the obj
s which captur

**strategies. Fig**Thus we have as was done f ry. Fig. 14 als e for optimiza bjective funct eters are requ at iteration num

**on of the optim**
**steeper region**

*egion *
**lays the resul**
creasing ensem
andom seeds.

need an ensem

or details. Note (2013).

**ying the positi**

ted above is fi
rking with a
oints to evalu
jective functio
re the nature o
**g. 14 shows th**
e 35 points to
for the Rosen
so shows two
ation, and a fl
tion space. W
uired to achiev
mber 7 is an ar
**mization proces**
**n of the curve) **

lts for the poi mble size with

The left subpl mble size of a

e that the mod

**on of the injec**

first tested on 320-control p ate our metho on space, due

of the objectiv
he optimization
test the meth
nbrock functio
dashed regio
atter region, i
While using the
ve the necessa
rtifact of the p
**ss (blue curve)**
**and a black el**
ints encapsula
h a constant p
lot of Fig. 15
approximately

del has been b

**ctors (blue) and**

a single reali
problem it is,
odology is al
to the ‘emptin
ve function sp
n process for
hodology, i.e.,
on, we now te
ons, a steep re
indicated with
e same hypoth
ary confidenc
plotting algori
**) using adjoint**
**lipse (indicativ**
ated by the re
perturbation s
illustrates the
y 150 to satisfy
benchmarked f
**d producers (re**
ization of the
of course, im
so risky since
ness’ of a 320
pace we perfo
35 iterations f
, instead of us
est the gradien
egion, indicat
h the black ell
hesis, we inve
ce intervals. N
ithm.
**t gradients. Th**
**ve of the flatter**
ed ellipse in F
ize, = 0.01.
e mean angle o

y the null hyp

for the differe

**ed) **
Egg Model a
mpossible to vi
e we will not
-dimensional
rm an adjoint
for an initial s
sing a large nu
nt quality for
ed with a red
ipse, which is
estigate the tw
Note: the appa

**e curve is divi**
**r region of the **
Fig. 14. It dep
The results a
of the differen
othesis. We a
ent simulators
and thus we ar
isualize the o
t be sure that
space. Thus, i
t- based optim
strategy of con
umber of rand
r 35 poins alo
d ellipse, in w
s indicative of
wo regions sep
arent decrease

**ided into two p**
**curve) **

picts the decre are the mean

nt points in th also observe th , details of wh re dealing wit bjective funct t our points h in order to obt mization from t nstant rates eq domly distribu ong a pre-defi which there ex f a peak, ridge parately, beca e in the object

**parts, a red elli**

ease in the m
of 5 experime
he red ellipse.
hat increasing
13
hich
th a
tion
have
tain
two
qual
uted
ined
xists
e or
ause
tive
**ipse **
mean
ents
We
the

14
ens
est
gra
num
400
res
Th
Re
[eq
gra
for
des
We
neg
stra
equ
ens
cov
**Fig**
**(14**
the
Sm
reg
fla
ach
cor
hyp
her
a s
equ
fur
bet
imp
not
gre
bec
dir
semble size le
timate a
high-adient estimat
merical round
0 and 300 sam
sults we have
he results pres
eynolds (2013
quations (5) a
adients and the
r different per
sired confiden
e observe that
gative impact
aight gradient
uation (11). W
semble sizes
variance vecto
**g. 15: (a) illustr**
**4), while (b) dis**
*The Flatter R*
The general
e straight grad
moothing of th
gions the mod

tter region is t
hieve the 95%
rrelation lengt
pothesis at mo
re were obtain
smoothed grad
**uation (11). F**
rther single sm
tter quality gr
proves with a
t observe a si
eater than 90
cause it impli
rection.
eads to higher
-quality gradi
te decreases, w
d-off errors. E
mples respect
used a time-c
sented in Fig
3) [equations
and (6)]. Irres
e modified an
rturbation size
nce intervals w
t smoothing o
t on the gradi
t obtained fro
We also note th
and inferior
or, i.e. equatio

**rates the perfo**
**splays the effec**

*Region *
trend in the re
dient which g
he gradients,
dified formula
the sensitivity
% confidence
ths varying fr
ost of the poin
ned through an
dient given by
**Fig. 16 shows**
moothing, i.e.
radients for al
an increase in
ignificant imp
degrees, i.e.
ies that even

r-quality gradi
ent. We also
while for sma
E.g., with =
tively to satisf
correlated cov
g. 15(a) are o
(17) and (18
spective of th
nd original form
es and correla
will then be di
of the gradien
ent quality fo
om equation (
hat a double s
results for lar
on (13), achiev
**ormance of the**
**ct of using diffe**
esults observe
gives good es
on the other
ation of Do an
y of the gradie
interval for
rom 8 to 12 c
nts. It is virtu
n exhaustive t
y equation (15
s the quality o
equation (13)
ll ensemble si
ensemble siz
provement in g
within the f
for smaller e
ients. Thus fo
tested (result
ller perturbati
0.1 (i.e. 10 t
fy the null hy
variance matri
obtained with
8)]. This perf
he gradient fo
mulation grad
ation lengths (
ifferent. Fig. 1
nt, i.e. pre-mul
or this set of p
14), i.e. equa
smoothing of
rger ensemble
ves better resu

**e two different **
**erent versions **

ed for the stee timates for th hand, achiev nd Reynolds (2 ent estimate qu this region w control time s ally impossib trial and error

) achieves a b of the differen ), which is aki zes. Note: irre ze. We also ob gradient quali first quadrant ensemble size or a 320-dimen ts not shown) ion sizes we a times larger), ypothesis and ix with a corr h equation (14 forms better t ormulation us dients converg (results not sh 15(b) shows th ltiplication of points and a ation (15), per the straight g e sizes. Note ults than those

**formulations **
** of the ‘smooth**
eper region is
he steeper reg
es much-high
2013) perform
quality to the c
we observe th
steps in conjun
le to know
a-approach. Sim
better angle co
nt smoothed g
kin to a double
espective of t
bserve that fo
ity. Also note
with respect
s we are able

nsional proble that for incr also obtain inf

and =0.000
achieve the d
relation length
4) and the m
than using eq
ed, an increa
ge for larger en
hown), althoug
he effect of sm
f the gradient
= 0.01. We
rforms better
radient, i.e. eq
: pre-multipli
e obtained wit
**used to estima**
**hed’ gradient a**
also observed
gion does not
her-quality gra
ms best. A maj
choice of the c
hat we need s
nction with an
priori the corr
milar to the re
ompared to di
gradient formu
e smoothing o
the different g
or this region w
e that, irrespec
to the adjoin
e to approxim
em with appro
easing perturb
ferior gradien
01, (i.e. 100 ti
desired confid
h equal to 20
modified formu
quation (14) w
sing ensembl
nsemble sizes
gh the ensemb
moothing on th
by the covari
e also observe
than using th
quation (16), g
ication of the
h equation (11

**ate the ‘straigh**
**and their relativ**

d for the flatte achieve good adients in this ajor difference correlation len smaller pertur n ensemble si rect correlatio esults for the s rectly estimat ulations. Whe of equation (14 gradient estim with ensemble ctive of ensem nt gradient. T mate, in gener oximately 150 bation sizes t nt estimates, m imes smaller) dence interval which was ar mulation sugge with the orig le size leads s. This is trend

ble size neede
he different gr
iance matrix,
e that single s
he smoothed g
gives better re
e gradient giv
1).
**ht gradient’ g g**
**ve behavior. **
er region. How
d gradients in
s flatter regio
e in gradient e
ngth used. Alt
rbation sizes,
ize equal to 3
on length, and
steeper region
ting a smoothe
en using equa
4), i.e. equatio
mates used, the

le sizes higher mble size, the This is particu ral, a roughly SPE 173236-0 samples we he quality of most likely du ), we would n . To obtain th rbitrarily chos ested by Do ginal formulat to higher-qua d is also obser ed to achieve radient estima has a margin smoothing of gradient given esults for sma ven by the cro

**given by equa**
wever, estimat
n the flatter p
on. Note, in b
estimation for
though we do
= 0.001,
300 to satisfy
d the lengths u
n, we observe t
ed gradient us
ation (11), the
on (16), achie
e gradient qua
r than 300 we
angles are ne
ularly import
‘correct’
up--MS
can
the
e to
need
hese
sen.
and
tion
ality
rved
the
ates.
ally
f the
n by
aller
**oss-tion **
ting
part.
both
the
not
and
the
used
that
sing
en a
eves
ality
e do
ever
tant,
-hill

SP
**Fig**
req
ens
the
det
mo
gen
opt
thi
25
two
**Fig**
**Fig**
**usi**
**elli**
E 173236-MS
**g. 16: Comparis**
* Robust case *
The main pu
quired to estim
semble of 100
e different pe

tails. The hyp
odel. The cont
nerate the poi
timization usi
s control traje
iterations due
o regions usin
**g. 17: Six rando**
.
**g. 18: Illustratio**
**ing adjoint gra**
**pse (indicative**

**son of the effec**

urpose of this
mate a
high-0 equi-probab
rmeability fie
pothesis tests
trols and the f
ints to test ou
ing the adjoint
ectory. The re
e to the comp
ng the same re
**omly chosen re**
**on of the expe**
**adients. The cu**
**e of the flatter r**
**ct of smoothin**
study is to qu
quality gradie
ble geological

elds and the d
are performe
fluid model a
ur methodolog
t method as de
esult of the rob
putational com
easoning as ex
**ealizations disp**
**ected objective**
**urve is divided**
**region of the c**
**g and different**

uantify the qua ent for a rea

models, six o different direc ed to estimate s well as the gy, the same a

escribed in Va bust optimiza mplexity invol xplained above

**playing the unc**

**e function valu**
**d into two parts**
**curve). **

**t gradient estim**

ality of the ‘ro
alistic reservo
of which are d
ctions and or
e a ‘robust’ e
other properti
approach is u
an Essen et al
ation is display
lved. Again, f
e.
**certainty in the**
**ue over 100 rea**
**ts, a red ellipse**

**mates for the f**

obust’ ensemb
ir test case. T
**displayed in F**
rientations of
ensemble grad
ies are exactly
used as describ
l. (2009), and

**yed in Fig. 18**
for the analys

**e geological m**
**alizations with**
**e (indicative of**
**latter region (b**
ble gradient an
Thus, to test
**Fig. 17. The u**
the channels
dient using al
y the same as
bed above. Th
then assess th
**8. Note that th**
sis, we divide

**odels; see Jan**

**the robust op**
**f the steeper r**
**black ellipse) o**
nd determine
our methodo
uncertainty is
s, see Jansen
ll hundred rea

for the determ
hat is, we firs
he quality of th
he optimizatio
e the optimiza
**nsen et al. (201**
**ptimization pro**
**region of the c**
**of Fig. (14). **
the optimal r
ology we use
captured throu
et al. (2013)
alizations of
ministic case.
st perform rob
he gradient al
on was limited
ation process i
**3) for details.**
**ocess (blue cu**
**curve) and a bl**
15
atio
e an
ugh
for
this
To
bust
ong
d to
into
**rve) **
**lack **

16
enc
the
all
rob
equ
bet
of
13
of
in
tha
adv
obj
thi
com
qua
**Fig**
**mo**
Fo
for
gra
wh
We
usi
for
sat
*The Steep Re*
Following th
capsulated by
e original and

the five poin
bust formulati
uations (5) an
tween 80 and
5 experiment
,15 and 16), w
the gradient l
that the smoo
at using a rela
vantage of us
jective functio
s region. Note
mpared to usi
ality gradients
**g. 19: Illustrate**
**odel realization**
r the Rosenb
rmulations. Fo
adient quality
here, similar to
e also observ
ing the modif
rmulations, al
tisfies the hyp

*egion *
he analysis pr
y the red ellips
modified form
nts we are abl
ion [based on
nd (6)] we nev
90 degrees fo
ts with differe
we observe mu
leads to better
othed version
atively larger
sing the mod
on) is from wi
e : when usin
ing smoothed
s for both grad

**s the differenc**
**n using the HPG**

rock function
or this case we
for the two
o the Rosenbr
e, as expected
fied formulat
lthough the M
othesis at all t
rovided for th
se as shown in
mulation usin
le to satisfy o
n equations (1
ver satisfy ou
or all the poin
ent random se
uch better ang
r estimates of
of equation (
r perturbation
dified formula
ithin this regio
ng non-smooth
gradients. No
dient formulat
**ce in gradient q**
**G ratio formula**
n we observed
e observe a si
ratio formula
rock results, w
d, that as the
ion. The resu
MIG formulati
the points, eve

he determinis n Fig. 18. We ng the computa our null hypot

9) and (20)]. ur hypothesis a nts when using eeds.) If we u gles as shown the ensemble (14), i.e. equa size ( = 1) ation with th on so achievin hed gradients, ote: while the tions with the

**quality betwee**
**ation. **
d significant
milar trend in
tions. Fig. 20
we always ach
ratio increas
ults indicate t
ion still perfo
en for the 1:1

stic case we f
have five poi
ationally attra
thesis, i.e. the
However if w
at any of the
g equation (14
use the‘smooth
**in Fig. 19, th**
e gradient, and
ation (15), ach
) for this regi
he 1:1 ratio b
ng good angle
, i.e. equation
e 1:1 ratio sati
HPG formula
**en the original **
differences in
**n the results. F**
0(a) depicts th
hieve significa
es the gradien
that the modi
orms better. I
ratio.
first consider
ints (control s
active 1:1 ratio
e angle is less
we use the 1:1
points. With
4). (As for the
hed’ versions
hough the hypo

d we observe
hieves better r
ion leads to th
because the m
es at acceptabl
n (14), the gra
isfies our hyp
ation as depict
**and modified f**
n gradient qu
**Fig. 20 consist**
he results obta
antly better gr
nt quality imp
ified formulat
Irrespective o
the steeper p
strategies) whi
o have been c
than 10 degr
1 ratio with th
the original f
e deterministic
of the gradie
othesis is still
a similar tren
results than eq
he best result
maximum con
le computation
adient quality
othesis, using
ted in Fig. (19
**formulations fo**
uality when u
ts of two plots
ained when u
radient estima
proves. Fig. 2
tion is less s
f the ratio fo
part of the op
ich lie in this
compared. We
rees, when us
he original fo
formulation w
c case, all resu
ent as given b
l not satisfied.
nd as for the d
quation (11).
ts. These resu
ntribution (lar
nal costs is m
is inferior (re
g higher ratios
9).
**or varying rati**
using the HPG
s which show
using the orig
ates with the M

20(b) depicts sensitive to th ormulation the SPE 173236-ptimization cu region for wh e observe that sing the modif

rmulation [ba we observe ang

ults are the m by equations ( . Thus smooth deterministic c We also obse ults highlight rgest increase more important

esults not show
s leads to high
**os of controls **
G and MIG r
the differenc
ginal formulati
MIG formulati
the results w
he different r
e modified fo
-MS
urve
hich
t for
fied
ased
gles
mean
(11,
hing
case
erve
the
e in
t for
wn)
**her-per **
atio
e in
ion,
ion.
when
atio
orm,

SP
**Fig**
**hig**
for
[ba
qua
‘sm
for
var
cor
**Fig**
**mo**
Sim
fla
Fig
mo
sat
pro
allo
for
sam
wit
com
we
E 173236-MS
**g. 20: (a) Differ**
**ghlights the gra**

*The Flatter R*
The 1:1 ratio
r the steeper p
ased equations
ality as shown
moothed’ grad
rmulations. No
ries between
rrelation lengt
**g. 21: Illustrate**
**odel realization**
milar to the re
**tter region. F**
g. 22 (a) depi
odified formul
tisfied our hyp
oblem dimens

owable error rmulation, we mples per geo thin the optim mputationally e suggest the u

**rence in gradie**
**adient quality w**

*Region *
o with the mod

part of the op
s (5) and (6)]
**n in Fig. 21. W**
dient, given b
ote: similar to
8 and 17. It
th which has a
**s the differenc**
**n using the HPG**
esults shown a
**ig. 22 consist**

icts the result lation. We ob pothesis for on sion, the differ

margin of 2 would satisfy ological realiz mization proce y challenging a

use of the mod

**ent quality for **
**when using the**

dified formula
ptimization cu
, the gradient
We also obser
by equations
o the determin
requires a su
a strong impac
**ce in gradient q**
**G formulation. **

above for the ts of two plots ts obtained w bserve that the

nly a few poin rence between 20 degrees th y our hypothe zation to estim ess we would and therefore, dified formula

**the two ratio f**
**e modified grad**
ation [based on
urve. However
quality is not
rve, slightly d
(11, 13,15 an
nistic case we
ubstantial num
ct on the gradi
**quality betwee**
steeper region
s which show
when using the

e same trend i nts for this rob n 10 and 20 de en, with the esis at most of mate a ‘good’

need to perfo , since the inc ation with a 1:

**formulations (H**
**dient formulati**

n equations (1
r, for the flatt
t as high. An
different to ou
nd 16), give
use a smaller
mber of simul
ient quality.
**en the original **
n, the differen
w the differenc
e original for
in the results
bust case. We
egrees is marg
MIG ratio f
f the points fo
’ gradient. Th
orm 5000 rese
crease in obje
1 ratio.
**HPG and MIG) **
**ion. Note the d**

19) and (20)] l ter part, thoug

increase in th
ur previous ob
very similar
r perturbation
lations to per
**and modified f**
nt ratio formul
ce in gradient
rmulation and
exists, althou
e could argue t
ginal. Thus if
formulation an
or a ratio of 1
his would imp
ervoir simulat
ctive function

**for the origina**
**ifferent vertica**
leads to very h
gh still better
he ratios signi
bservations tha
results for bo
size ( = 0.01
rform an exte
**formulations fo**
ations show a
quality for th
d 22 (b) depic
ugh on the fla
that for the ro
we were to ha
nd both the o
:50, i.e. we w
ply that, to est
ions since the
n value for the

**al gradient for**
**al scales **

high quality g
r than the orig
ificantly impro
at the differen
oth the modif
1) and a corre
ensive search
**or varying rati**
a similar trend
he different ra
cts the results
atter part of th
obust case and
ave a hypothe
original or m
would need to
timate a grad
e ensemble siz
e flatter region
**mulation while**
gradient estima
ginal formulat

oves the grad nt versions of

fied and orig elation length t

for the ‘corr

**os of controls **

d for the point
atio formulatio
s when using
he curve we h
d, considering
esis which had
modified grad
apply 50 con
dient at one po
ze is 100. Thi
n is less than
17
**e (b) **
ates
tion
ient
f the
inal
that
ect’
**per **
ts in
ons.
the
have
our
d an
ient
ntrol
oint
is is
1%