Generalized Automatic Global Predictive Optimization
and Forecasting*
October 2024
Also See:
Our new optimization page
describing our new optimizer SensorOpt is coming soon. We're
finishing up our probabilistic workflows to quantify uncertainty in
optimizations, predictions, and forecasts. For now, see this SensorOpt
log file, Spe5s.txt, for the solution to
the deterministic example problem given on our old optimization page (below)
to determine the optimal production strategy for SPE5^{1},
obtained in a few minutes as opposed to about an hour to make several
optimization runs for the previous solution: It maximizes NPV to US
$1.9438 billion at 11081 days and recovers 100% of the oil in place.
Note that the solution gives an optimal gas flood producing bhp of 3755 psia,
the same as the previous solution. This is not a coincidence. It
places the reservoir almost entirely in the firstcontact miscible region
during the optimal gas flood (firstcontact fluid saturation pressure = 3874
psia, drawdown is about 80 psia).
Also see this SensorOpt log file for the same problem,
but with the added constraint that total runtime is less than or equal to
9000 days: Spe5c.txt. The unconstrained
problem gives the solution that absolutely maximizes NPV and oil recovery
regardless of the time taken to achieve it, and is mostly of academic
interest. The constrained optimized solution maximizes NPV to $1.9148
billion. at 8983.6 days, recovers 98.1% of the oil in place, and gives a
12.56% average annual net return on investment (ANROI is the yield on top of
inflation, assuming everytthing inflates at the same rate). The
unconstrained solution gives an ANROI of 8.5%.
We can easily and accurately determine the optimal
production strategy for our reservoirs, including the optimal well and
fracture spacings that are part of it, using simple models incorporating our
reservoir and fluid properties. Numerical error (dispersion) is
eliminated through the use of fine gridding or proper upscaling. See
our SPE3 and SPE8 pages demonstrating the critical need for proper upscaling
in our coarse grid reservoir models. In these spe5 optimization
problems, we have thus far assumed that the originally given grid and
properties are properly upscaled to reproduce the correct finescale
solution.
Generalized Automatic Global Predictive Optimization and Forecasting*
(2013) Our associate Dr. Curtis Whitson
and his team developed workflow integration, automation, and
iterative optimization software called ECIAM^{®}
(Tieto). ECIAM^{® }can
integrate any number of linked batch applications to create automated
workflows, to which various available iterative optimization methods can be
applied. Any program inputs can be automatically and iteratively
adjusted to maximize or minimize any benefit function that can be defined
from the workflow outputs. They have solved the generalized automatic
optimization problem, which naturally leads to the greatest advances in
reservoir modeling in decades, when applied to reservoir simulation.
Using our reservoir simulator Sensor, we have created automatic workflows in
ECIAM^{®}
for:

deterministic predictive
optimization and forecasting (SensorCast)

probabilistic predictive
optimization and forecasting (SensorPCast)

deterministic upscaling and history
matching (SensorMatch)

probabilistic discrete alternative
optimization (SensorPopt)
The iterative optimization problems
are very similar and differ only in the definition of the benefit function
to be minimized or maximized, and identity of the adjusted variables, and
the (probabilistic) consideration of uncertainties. Workflows for
characterization, upscaling, and history matching are essentially the
same. They minimize uncertainty by minimizing a mismatch function.
Global predictive optimization adjusts operational variables to maximize an
overall value function such as cumulative BOE or NPV, making most
instantaneous optimizations obsolete, other than those designed to achieve
the globally optimal solution.
The probabilistic workflows require
as input an arbitrary number of equallyprobable realizations of the
uncertainties (the number of realizations required must be determined by the
user as those that give reproducible probabilistic results). Sensor
can represent some of those uncertainties through its Uncertain Inputs
features, in which each execution of the program automatically creates an
equallyprobable realization of the specified uncertainties in the input
data. Our programs MakeSpx and SensorPx are used to make an arbitrary
number of runs representing the uncertainties and compute the probabilistic
results from them, within each iteration of the optimization.
Our first example illustrates the
use of SensorCast for automated deterministic predictive optimization
(considering no uncertainties in the input data other than the optimized
variables). It is based on the reservoir and fluid description from
the Fifth SPE Comparative Solution Project^{1}. Results from a
number of miscible flood (extended black oil) and fully compositional
simulators were compared for three given production/injection scenarios.
The Scenarios are described in detail in Appendix 1 of the reference^{ibid.}
The first Scenario is 2 years of
depletion with producing BHP=1000, followed by a WAG flood with 1 year gas
and water injection cycles. The second Scenario is a primary WAG flood
with 3 month gas and water injection cycles with higher maximum rates.
The third Scenario is 1 year of depletion using BHP=1000, followed by a 1
year waterflood, followed by a WAG flood using the same 3 month cycle and
WAG rates as Scenario 2 except that maximum gas injection rate is increased
from 20000 mcf/d to 30000 mcf/d. Scenarios 1 and 3 use a minimum
producing BHP of 1000 psia during all production periods, but Scenario 2
limits producing BHP to 3000 psia during WAG cycling, to better operate in
the miscible region. The Sensor compositional model datafiles for
these three scenarios are spe5a.dat, spe5b.dat, and spe5c.dat. These
datafiles do not include Sensor's options for treatment of oil and gas
relative permeability and Pcgo as functions of gas/oil interfacial tension,
as miscibility is approached, which may make a significant difference in
results. We used the optional maximum trapped gas saturation of .20 for Krg
hysteresis. In Scenarios 1 and 3, we used both injector and producer
for depletion during the depletion periods (for a better determination of
optimal strategy).
The problem specification includes
a labmeasured minimum miscibility pressure (MMP) of 3000 psia, which is
used in miscible flood simulators as a condition of miscibility. Miscibility
in reservoirs is generally not achieved at MMP measured in a laboratory
slimtube experiment, because of the compositional mixing in the reservoir
due to heterogeneity and multiphase flow and capillary and gravity forces
and thermodynamics, which prevents the very gradual 1D enrichment that is
required to achieve MMP well below the pressure required for firstcontact
miscibility (the FCM pressure) in the lab experiment (see our
MMP page). Optimal
minimum (producer) miscible flood pressure is generally slightly lower than
the firstcontact miscible (FCM) pressure. An FCM pressure exists for
a given fluid/solvent if there exists a closed Px phase envelope (a hubble
point, dew point, or critical point exists for all mixtures) as the maximum
saturation pressure of any mixture. Sensor prints the phase envelope
for the given oil / injected gas mixture in all output files, showing that
the FCM pressure for the given oil and injected gas compositions is 3874.8
psia at a composition of (.2 oil, .8 solvent).
Scenario 1 (spe5a.dat)
terminates at 6317.9 days when the field watercut limit is exceeded, with
oil and gas recoveries of 70.2% and 34.8%, cumulative produced oil, gas, and
water are 29398 mstb, 40117 mmscf, and 4411 mstb, respectively.
Cumulative water and gas injection are 35040 mstb and 32015 mmscf,
respectively. HC average pressure at end of run is 3286 psia.
Scenario 2 (spe5b.dat)
terminates at 6186 days when the field GOR limit is exceeded, with oil and
gas recoveries of 78.4% and 1.4%, cumulative produced oil, gas, and water
are 32824 mstb, 59655 mmscf, and 73 mstb, respectively. Cumulative water and
gas injection are 29738 mstb and 59976 mmscf, respectively. HC average
pressure at end of run is 3811 psia.
Scenario 3 (spe5c.dat)
terminates at 3531 days when the field GOR limit is exceeded, with oil and
gas recoveries of 72.0 and 17.6%, cumulative produced oil, gas, and water
are 30166 mstb, 44469 mmscf, and 2585 mstb, respectively. Cumulative
water and gas injection are 31845 mstb and 40357 mmscf, respectively. HC
average pressure at end of run is 3501 psia.
To improve accuracy of the results,
we added Sensor's options for treatment of oil and gas relperms and
capillary pressure as functions of gas/oil interfacial tension (KRIFT
and TENSION). If we assume a lab reference IFT for the air/water
system used to determine Pcgo, then the Pcgo scaling performed by our
TENSION option scales the curve down to negligible values at immiscible and
miscible reservoir conditions. So in the following datafiles we use
the KRIFT option to treat oil and gas relperms as functions of IFT, but omit
the TENSION option for Pcgo scaling. The KRIFT option without the
following PC label automatically zeroes Pcgo. The Sensor datafiles for
the 3 Scenarios using KRIFT are spe5akrift.dat, spe5bkrift.dat, and
spe5ckrift.dat. Results are:
Scenario 1 (spe5akrift.dat)
terminates at 6357.4 days when the field watercut limit is exceeded, with
oil and gas recoveries of 70.3% and 34.7%, cumulative produced oil, gas, and
water are 29430 mstb, 40574 mmscf, and 4553 mstb, respectively.
Cumulative water and gas injection are 35040 mstb and 32489 mmscf,
respectively. HC average pressure at end of run is 3218 psia.
Scenario 2 (spe5bkrift.dat)
terminates at 6478.8 days when the field GOR limit is exceeded, with oil and
gas recoveries of 79.3% and 11.0%, cumulative produced oil, gas, and water
are 33191 mstb, 64812 mmscf, and 192 mstb, respectively. Cumulative water
and gas injection are 31558 mstb and 62252 mmscf, respectively. HC average
pressure at end of run is 3664 psia.
Scenario 3 (spe5ckrift.dat)
terminates at 3376.2 days when the field GOR limit is exceeded, with oil and
gas recoveries of 71.5 and 22.0%, cumulative produced oil, gas, and water
are 29946 mstb, 42947 mmscf, and 21918 mstb, respectively. Cumulative
water and gas injection are 31280 mstb and 37822 mmscf, respectively. HC
average pressure at end of run is 3501 psia.
Adding the KRIFT option made the
most difference in gas production in Scenarios 2 and 3 (8.6% more in
Scenario 2).
The economic model parameters
specified in Example 1 below were used in the above datafiles, while
retaining the originally specified economic limits.
Reported NPV at end of run for the last 3 KRIFT runs are $1.40575E9 at
6357.4 days, $1.57631E9 at 6478.8 days, and $1.43252E9 at 3376.2 days,
respectively.
Example 1  Extension of the
SPE5 problem to determine the optimal production/injection strategy
Determine the optimal
production/injection scenario for the SPE5 problem. Allow for a
possible primary depletion period, possibly followed by a waterflood,
possibly followed by a WAG flood, possibly followed by blowdown. Use the
same maximum production rate for the producer (12000 stb oil/d) along with
the same minimum producing bhp of 1000 psia for depletion, waterflooding,
and blowdown. Assume that the maximum well injection rates given in any
period of the original Scenarios represent fixed well or facility limits.
Instead of the specified limiting GOR and WCUT values for shutting in the
field, use the following economic model variables to determine the optimal
scenario by maximizing Net Present Value:
Discount factor and all
interest and inflation rates are zero. Capex for lease and
facilities excluding wells is $10,000,000 at time 0. Capex per
well is $2,000,000 at time 0. Facilities Opex is $400/day and well
Opex is $200/day/well. The values/costs of fluid production and
injection are:
Vo = $50/stb
Vg = $2.0/mscf
Cw = $1.0/stb (cost
of water production)
Cwi = $1.2/stb (cost of
water injection)
Cgi = $2.5/mscf (cost of
gas injection)
Solution
We first define the optimization
variables and their initial values and minimum and maximum limiting values
for optimization in the ECIAM Optimizer. A link is created for each
variable to its specification in one of the data or output files. The
initial values are based on original specifications or on our best guesses
of what we think may be optimal:
Variable
Description
Min Max
Initial
TD
depletion end time, d 10
7300 365
DTWF
waterflood time, d
10 7300
365
DTWAG WAG flood time, d
10 10000
6500
DTWAGW water inj. cycle time, d
1.0 500
91.25
DTWAGG gas inj. cycle time, d
1.0 500
91.25
DTBD
blowdown time 365 7300 1370
BHPWFI wf injector max bhp
3000 6500
6000 (assume frac P > 6500)
BHPWAGP wag producer min bhp
1000 4000
3000
BHPWAGI wag injector max bhp
4000 6000 4500
The set limits allow any of the
production phases to be eliminated or to be included down to their minimums
and up to their maximum
time limits. The optimized values of these variables represent the
optimal operating strategy.
The Sensor data file prior to
optimization using the initial values of the defined variables is
spe5d0.dat. spe5d0.out reports
that the NPV at end of run for the initial guess is $1.5766E9 at 8600 days
(about the same as Scenario 2 but 3 years later), and the fractional return on investment ROI
is 5.485.
NPV = cumulative npv of
production  cumulative npv of investment, $ at t=0
ROI = cumulative npv of
production / cumulative npv of investment, fraction
Injected fluid costs are accounted
for as investment expenses.
We note that miscible recovery is
not achieved in much of the bottom 2 layers, due to gravity override of the
oil by the injected gas. So we make a test run changing the injection
well perforation from layer 1 to layer 3. NPV increases to
$1.61823E9 at 8600 days (about $42 million higher), and ROI increases
by 27% to 6.916.
We then copy this modified datafile (spe5di.dat) to Sensor.dat in our
ECIAM project directory data folder. This is an example of manual
predictive optimization, using engineering judgment. An alternative is
to define as additional optimization variables the injector (and producer) completion layer
ranges (as K1 K2 in the well data) for all or for each of the production
phases, but for this example we will consider the completed layers as fixed
in the bottom layer for both wells in all production phases.
The above variables are linked
directly to their specifications in the Sensor.dat file (the DT* variables
are defined in comments). Other "auxilliary" variables
representing endtime specifications in the Sensor datafile are defined in
the optimizer as specified functions of the optimization variables:
TW end of waterfloodng.
= TD + DTWF
TWAG end of WAG flooding = TD
+ DTWF + DTWAG
TEND end of blowdown =
TD + DTWF + DTWAG + DTBD
The auxilliary variables are automatically calculated and set in the Sensor.dat datafile by the optimizer at the beginning of each optimization
iteration.
The recurrent well data in
Sensor.dat (spe5di.dat) containing all of the defined variables is:
WELL
I J K
INJR
1 1 3
PRODR
7 7 3
WELLTYPE
INJR STBOIL ! first define the injector as a producer, for depletion using
both wells (changed May 2018), previously only PRODR was used for depletion
PRODR STBOIL
BHP
INJR 1000
PRODR 1000
INJGAS
INJR
.77 .2 .03 3*0
RATE
INJR 12000
PRODR 12000
C ** if IMPLICIT, limit first time step at injection switchover shocks
C DTSTART 4
TIME 365 ! TD, end of depletion
WELLTYPE
INJR STBWATINJ
BHP
INJR 4500 ! BHPWFI
RATE
INJR 45000
TIME 730 ! TW = TD + DTWF, end of waterflood
C DTWF 365
BHP
INJR 4500 ! BHPWAGI
PRODR 3000 ! BHPWAGP
TIMEWAG
7230 ! TWAG = TD + DTWF + DTWAG, end of wag cycling
INJR ! using well INJR as injector for repeated wag cycles
91.25 45000 3 ! DTWAGW qwwag
3=all injectors,
! max injection of 45000 stb/d water for dtwagw days, followed
by
91.25 30000 1 ! DTWAGG qgwag
1=all producers,
! max injection of 30000 mcf/d gas for dtwagg days
C DTWAG 6500
WELLTYPE
INJR STBOIL
BHP
INJR 1000
PRODR 1000
RATE
INJR 12000
MAPSFREQ 1
TIME 8600 ! TEND = TD + DTWF + DTWAG + DTBD, end of blowdown, and end of run
C DTBD 1370
END
We define the objective function in
the ECIAM Optimization window by creating a link to the
value of NPV reported at end of run in the Run Summary of the Sensor.out
file and assign it as the objective function to be maximized. Also
reported are the values of the maximum NPV achieved during the run and the
time it is achieved. That time is the optimal end time as reported by
Sensor, if it occurs before TEND. If not, the (exactly) optimal end
time for any given scenario is easily determined by increasing TEND until
NPV begins to decrease.
After setting up the above (9
optimization, 3 auxiliary, and 1 objective function) variables in the
ECIAM Optimization and Linkz windows (taking about 20 minutes), we are ready
to optimize. The general solution method to find a global maximum,
which is demonstrated here as a manual process but can be automated, is to
keep executing the same or different 'solvers' available in ECIAM in series,
and by trying different starting points until we're unable to significantly
improve NPV. This is simply automation and improvement of the method
that engineers have always used in reservoir modeling to make optimizations
between discrete alternatives, for example as examined in the original SPE5
paper for the three given production scenarios using various reservoir
simulators. We select the Reflection optimizer (or "solver", a
proprietary Petrostreamz/Tieto
implementation of the Simplex method) in the ECIAM Optimization window, and
allow for 500 maximum iterations. The number of required iterations
and convergence is highly dependent on the number of optimization variables
(and the optimization method used). We then select "Run" in the Optimizer window and the optimizer begins to
iterate.
The optimization terminates after
500 iterations with an optimized NPV of $1.89303E+9 at 10992 days, an
increase of $172 million or 17%, 6 years later. Total workflow
execution time is about 12 minutes. Spe5_opt1.xlsx gives the iterative
history of the optimization. The best found solution is repeated at
the end as iteration no. 501.. The optimized datafile Sensor.dat is copied to
spe5_opt1.dat. The values of the
optimized variables are:
TD
170.3299861
DTWF 238.0863996
DTWAG 9988.348595
DTWAGW 35.82821765
DTWAGG 495.3740961
DTBD 595.5213266
BHPWFI 4671.993396
BHPWAGP 3349.564431
BHPWAGI 4483.139139
TW 408.4163857
TWAG 10396.76498
TEND 10992.28631
NPV 1893030000
Since DTWAG has reached its maximum limit, we reset
its maximum value in the optimizer to 12000 days. We also reset minimum DTWAG
to 6000, and continue optimization allowing another 500 iterations.
NPV is increased by over $500 million to 1.94581E+9 at 12052.2 days
(spe5d_opt2.dat, spe5d_opt2.out). This solution is nearoptimal (near
a local or global optimum). NPV for this solution is within 0.1% of
the best solution that we have found for this problem upon further
optimization, as described below.
TD 136.50915184286688
DTWF 599.9401793463742
DTWAG 10659.776837735773
DTWAGW 11.928941344948099
DTWAGG 477.7340847527599
DTBD 655.9837457117383
BHPWFI 4767.714517037251
BHPWAGP 3755.2291830134295
BHPWAGI 4558.572982510338
TW 736.4493311892411
TWAG 11396.226168925015
TEND 12052.209914636755
NPV 1945810000
To check for possible further improvement, we switch to
the NOMAD optimizer in ECIAM, and allow another 500 iterations, which take
another 14 minutes to complete. NOMAD works well for improving the
accuracy of near optimal solutions, and also by itself on most cases.
NPV improves by almost a million dollars to 1.946674E9 at 11736.1 days
(spe5d_opt3.dat, spe5d_opt3.out):
TD 136.51121714380835
DTWF 599.9427757247006
DTWAG 10294.159347921015
DTWAGW 11.929684912891672
DTWAGG 477.7286329125976
DTBD 705.4821524795834
BHPWFI 4767.722027222493
BHPWAGP 3755.222388083925
BHPWAGI 4558.572267254601
TW 736.453992868509
TWAG 11030.613340789525
TEND 11736.09549326911
NPV 1946673753.4
An optimal solution is indicated mathematically by dNPV/dx
= 0 for all optimization variables x. No solution found by the
optimizer will be exactly optimal. When the optimizer converges in
less than the maximum allowable iterations, its solution is very close to a
local or global optimum When it does not converge but enough
iterations are allowed, the optimizer will usually find and 'bounce around'
a local or global optimum but without satisfying convergence tolerances.
This tends to increasingly happen as the number of optimization variables
increases, and as their codependencies become stronger. Although the
optimization variables are completely independent, their optimal values are
very highly dependent on each other (with the exception of optimal water and
gas flood operating pressures, which are mainly dependent on phase behavior).
This makes the NPV response surface very complex with multiple local optima
and the optimization extremely difficult or impossible to perform manually
when there are more than just a few variables. It also makes it
virtually impossible to replace or sufficiently estimate the
simulatorpredicted NPV function or its optimization with any proxy model,
neural network, or any 'artificial intelligence' or 'machine learning'
method. To demonstrate this, NPV results for an unlimited number of
equallyprobable realizations of this example problem were offered to those
claiming proxy or neuralnetwork or any ai/ml model applicability to
reservoir engineering. No solutions were offered based on the 100
realizations given nor were any more requested. See the SPE Reservoir
technical group discussion of "Artificial
Intelligence / Machine Learning Problem Solving Methods"
(SPE and
Reservoir group memberships are required). Also see our page on
Artificial
Intelligence.
The optimizer may find a local or a global optimum in any
given optimization run. A locally optimal solution that may be found
can depend on the iterative solution path taken by the optimizer from the
assumed initial variable values, which depends on those values and the
specified variable limits and the optimization algorithm. So to find a
global optimum, we use different starting points for optimization while
trying to allow for any reasonable values of the variables. Continuing
an optimization with the same or a different optimizer sometimes results in
an improved solution. Some optimizers are much better than others. Of
the available optimizers in ECIAM, we have found that both the Reflection
and the NOMAD solvers work well on this problem and others like it.
From the ECIAM Pipeit Manual,
The Nelder and Mead (1965) reflection simplex
algorithm is the default "real" solver provided with PipeIt. This
algorithm does not require derivatives and is very robust. It utilizes a
direct search method of optimization that works moderately well for
stochastic problems. It is based on evaluating a function at the
vertices of a simplex, then iteratively shrinking the simplex as better
points are found until some desired bound is obtained. Lagarias et al.
(1998) studied the convergience properties of the NelderMead reflection
simplex, and concludes that they are at best linear."
NOMAD (Nonlinear Optimization by Mesh Adaptive Direct
Search) is a C++ implementation of the Mesh Adaptive Direct Search (Mads)
algorithm.
http://www.gerad.ca/nomad/Project/Home.html
We can check that our optimized solutions are close to a
local or global maxima by perturbing each optimization variable in turn,
using perturbations that are large enough to have a significant effect on
the objective function. At an optimal solution, all perturbations of
all variables plus and minus should result in a decrease in NPV. Test runs
were made by manually perturbing all variables by small amounts in the
optimizer interface to confirm that spe5_opt3 is very close to an optimal
solution, and all perturbations made result in decreasing NPV. For
example, changing time of depletion TD to either 130 or 140 days results in
significant decrease in NPV. Even the very small water injection cycle
time of 11.93 for the WAG flood is found to be very accurately optimal,
which is a little surprising given the very long optimal gas injection cycle
time of 477.73 days. These tests take only a few minutes in this case.
However, since the last NOMAD optimization increased NPV by about a million
dollars, we continue with another 500 NOMAD iterations (spe5d_opt4.dat,
spe5d_opt4.out). NPV increases by about half a million dollars to
1.947266E9 at 11736.0 days. Another 1000 iterations of NOMAD
(spe5d_opt5.dat, spe5d_opt5.out) increases NPV by about $330 thousand
dollars to $1.947597E9 at 11735.92 days. The NOMAD optimizations
are slowly increasing the accuracy and NPV of the same optimal strategy that
was found in the second optimization using the Reflection solver, spe5_opt2.
Our spe5d_opt5 solution, taking a total time of about an hour to perform the
5 optimizations, is
TD 136.3965841775175
DTWF 600.0111131013296
DTWAG 10294.03529942354
DTWAGW 11.946483704052952
DTWAGG 477.67823466207545
DTBD 705.4771862199365
BHPWFI 4768.00586143727
BHPWAGP 3755.3412987055735
BHPWAGI 4558.521919736671
TWF 736.4076972788471
TWAG 11030.442996702386
TEND 11735.920182922322
NPV 1947597546.3
To see it could possibly be better
to completely deplete the reservoir first, we run case spe5d_dep.dat
simulating only depletion to 10,000 days. The only optimization
variable is the optimal end time, which Sensor reports at the end of
spe5d_dep.out as the time of maximum NPV:
MAX NPV (PRIOR TO
ABEX) = 0.6271021E+09, AT TIME = 0.7197584E+04 DAYS
A very important observation is
that if production strategy is not optimized and only depletion is
considered, NPV of the reservoir is reduced by 67.8%.
We then specify 7197.58 days as end time of depletion TD
and deactivate it as an optimization variable, and reoptimize using 500
iterations of the Reflection solver, starting with the other optimized
variable values from spe5d_opt5. The best found solution
(spe5d_dopt1.dat, spe5d_dopt1.out) is NPV = 1930191954.0 at time 18982.5
days. We then apply 500 iterations of NOMAD and NPV improves by about
$500,000 to 1.930689E9 at 19031.9 days giving our optimal solution for the
case of mandatory depletion to 7197.59 days (spe5d_dopt2.dat,
spe5d_dopt2.out).
TD 7197.58 (specified)
DTWF 512.9597784073119
DTWAG 10602.090610824198
DTWAGW 3.2521265727310475
DTWAGG 483.77406322367017
DTBD 719.2467111404744
BHPWFI 4391.326519799346
BHPWAGP 3689.1731626916003
BHPWAGI 4517.1076610443215
TWF 7710.539778407312
TWAG 18312.63038923151
TEND 19031.877100371985
NPV 1930689086.2
So requiring total depletion
first only costs 17 million dollars in final NPV, but at a very much longer
final time (19032 vs. 11736 days, 20 years longer!). If discount
factor were significant here, the lost value due to requiring initial
depletion of the reservoir would be huge.
To compare optimizers and to possibly find a better
solution, we try 4 successive applications of 500 maximum iterations of the
NOMAD optimizer, with the same spe5di.dat starting point, with maximum DTWAG
set to 12000. A slightly different locallyoptimal solution is found
(in less than an hour) with an NPV of 1944930000, with the main difference being a
much shorter depletion period of only 10 days, which is the specified
minimum (spe5d_optn4.dat, spe5d_optn4.out):
TD 10
DTWF 733.2169148921967
DTWAG 9947.705966472626
DTWAGW 1.0603598713874816
DTWAGG 125.81382192373273
DTBD 661.9549714326859
BHPWFI 4706.688094139099
BHPWAGP 3737.332844734192
BHPWAGI 5689.139604568481
TWF 743.2169148921967
TWAG 10690.922881364822
TEND 11352.877852797508
NPV 1944930000
Both spe5d_opt5 and spe5_dopt2 solutions might be
improved slightly but are tightly converged to optimal solutions, which is
confirmed with declining NPV upon small perturbations of the variables.
In spe5d_dopt2, the extremely small optimal water injection phase time of
3.2521 days (with gas injection time of 483.77 days) for the WAG flood and
its strong effect on NPV is very surprising, but accurate. Values of
3.2 and 3.3 days both give reduced NPV values of 1.929770E9 and 1.928870E9,
respectively (reductions of about 0.8%).
We believe that spe5d_opt5 represents the globally
optimal solution.
The firstcontact miscible pressure
for the oil and the specified injected solvent is given by the maximum
saturation pressure on the solvent/oil phase diagram printed n the Sensor
output file, and is equal to 3874.8 psia. It is interesting that the
optimal minimum bottomhole pressure for the producer (3755 psia, first
determined accurately in spe5d_opt2) is only just below this fcm pressure,
and not even close to the given labmeasured slimtube MMP of 3000 psia.
Any (miscible) simulator assuming a reservoir MMP applied to Scenario 2
of the original spe5 case is incorrect. Incorrect application of
labmeasured slimtube MMP to field operations is likely the main reason for
their failures to achieve high miscible recovery. The final oil saturation map
shown below for
the optimized case shows that
miscible recovery was achieved in most of the reservoir. This is not
true in any of the original Spe5 scenarios (counter to statements in the
paper). Final oil and gas recoveries in spe5d_opt5 are 100.9% and
17.8%, respectively, and final average reservoir pressure is 1002 psia.
Very high miscible recovery is achieved in the confined model, just as might
be achieved from interior wells in a pattern miscible flood. This case
might be considered to be a model of a quarter5 spot pattern, with the well
rates representing partial well values. To do so rigorously for this
example problem, only the well PI values and block sizes on the boundary
need adjustment for the partial wells placed along the outer lateral
boundaries of half and quarter blocks that are created using our EDGE
option. The constant value of the lateral block size can then be added
as an optimization variable in order to additionally optimize well spacing
for the simultaneous pattern flood. We have long known that simultaneous
pattern operation is generally required to optimize recovery from our
petroleum reservoirs.
SO
TIME = 11735.9 DAYS DATE: 0 0 0 MAP WINDOW 1

K = 1

J I= 1 2 3 4 5 6 7
1 0.0000 0.0000 0.0000 0.0000 0.0001 0.0014 0.0076
2 0.0000 0.0000 0.0000 0.0000 0.0001 0.0008 0.0067
3 0.0000 0.0000 0.0000 0.0000 0.0000 0.0005 0.0059
4 0.0000 0.0000 0.0000 0.0000 0.0000 0.0001 0.0056
5 0.0000 0.0001 0.0000 0.0000 0.0000 0.0001 0.0033
6 0.0015 0.0009 0.0002 0.0003 0.0001 0.0000 0.0005
7 0.0079 0.0067 0.0057 0.0055 0.0033 0.0005 0.0075
K = 2

J I= 1 2 3 4 5 6 7
1 0.0000 0.0000 0.0000 0.0000 0.0001 0.0025 0.0196
2 0.0000 0.0000 0.0000 0.0000 0.0001 0.0018 0.0220
3 0.0000 0.0000 0.0000 0.0000 0.0001 0.0007 0.0209
4 0.0000 0.0000 0.0000 0.0000 0.0000 0.0003 0.0145
5 0.0001 0.0001 0.0001 0.0000 0.0000 0.0001 0.0113
6 0.0026 0.0018 0.0008 0.0004 0.0002 0.0001 0.0034
7 0.0201 0.0195 0.0222 0.0155 0.0119 0.0034 0.0108
K = 3

J I= 1 2 3 4 5 6 7
1 0.0000* 0.0000 0.0000 0.0000 0.0001 0.0044 0.0715
2 0.0000 0.0000 0.0000 0.0000 0.0002 0.0040 0.0766
3 0.0000 0.0000 0.0000 0.0000 0.0002 0.0034 0.0805
4 0.0000 0.0000 0.0000 0.0000 0.0002 0.0023 0.0683
5 0.0001 0.0002 0.0002 0.0002 0.0002 0.0014 0.0606
6 0.0046 0.0041 0.0036 0.0025 0.0015 0.0013 0.0153
7 0.0731 0.0764 0.0820 0.0710 0.0633 0.0154 0.0568*
The SensorCast workflow is very simple and consists only
of the Sensor datafile being executed by Sensor to produce its ascii text
output file. The optimization variables are graphically linked to
their values in Sensor.dat and the NPV objective function is graphically
linked to its reported value at the end of Sensor.out:
Example 1b
What is the value of the final
blowdown period in the optimized solution of Example 1?
Answer: Fort.75 from run spe5d_opt5 gives NPV vs. time
output from the economic model, showing that NPV increased by $27.03 million
in the final 705.4 days of the run during blowdown.
1. Killough, J. E., "Fifth
Comparative Solution Project: Evaluation of Miscible Flood Simulators", SPE
16000 presented at the Ninth SPE Symposium on Reservoir Simulation held in
San Antonio, Texas, February 14, 1987
* This example problem was first
presented in December 2020 in the SPE Connect Reservoir Group Discussion of
"Artificial Intelligence / Machine Learning Problem Solving Methods" in
which it was proven that no such useful methods exist in reservoir
engineering, invalidating very many published claims, and SPE technical
sections and courses. A spreadsheet giving NPV results for some number
of realizations of the unknowns was given with any number offered, but
nobody was able to predict an optimal solution.
