OilfieldToolsNet/opm-common
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Merge pull request #2119 from jalvestad/glopt

Browse Source 
Enable eclipse compatible restart for wells with gas lift optimization
This commit is contained in:
Joakim Hove
2020-12-17 15:56:48 +01:00
committed by GitHub
parent 2c5403409e 229950a591
commit 93913d64a3
15 changed files with 51305 additions and 11 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
CMakeLists_files.cmake
View File

@@ -407,6 +407,7 @@ if(ENABLE_ECL_OUTPUT)
tests/test_InteHEAD.cpp
tests/test_LinearisedOutputTable.cpp
tests/test_LogiHEAD.cpp
tests/test_LGOData.cpp
tests/test_OutputStream.cpp
tests/test_regionCache.cpp
tests/test_PaddedOutputString.cpp
@@ -460,6 +461,7 @@ if(ENABLE_ECL_OUTPUT)
tests/SPE1CASE2.X0060
tests/PYACTION.DATA
tests/0A4_GRCTRL_LRAT_LRAT_GGR_BASE_MODEL2_MSW_ALL.DATA
tests/2_WLIFT_MODEL5_NOINC.DATA
tests/act1.py
tests/MSW.DATA
tests/EXIT_TEST.DATA

7
opm/output/eclipse/DoubHEAD.hpp
View File

@@ -51,6 +51,12 @@ namespace Opm { namespace RestartIO {
double damping_fact;
};
struct liftOptPar {
double min_int;
double incr;
double min_ec_grad;
};
DoubHEAD();
~DoubHEAD() = default;
@@ -72,6 +78,7 @@ namespace Opm { namespace RestartIO {
DoubHEAD& udq_param(const UDQParams& udqPar);
DoubHEAD& guide_rate_param(const guideRate& guide_rp);
DoubHEAD& lift_opt_param(const liftOptPar& lo_par);
const std::vector<double>& data() const
{

1
opm/output/eclipse/InteHEAD.hpp
View File

@@ -146,6 +146,7 @@ namespace RestartIO {
InteHEAD& variousUDQ_ACTIONXParam();
InteHEAD& nominatedPhaseGuideRate(GuideRateNominatedPhase nphase);
InteHEAD& whistControlMode(int mode);
InteHEAD& liftOptParam(int in_enc);
const std::vector<int>& data() const
{

3
opm/output/eclipse/VectorItems/doubhead.hpp
View File

@@ -53,6 +53,9 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems
GRpar_d = 90, // Guiderate parameter D
GRpar_e = 91, // Guiderate parameter E
GRpar_f = 92, // Guiderate parameter F
LOminInt = 93, // LIFTOP - Minimum interval between gas lift optimizations
LOincrsz = 95, // LIFTOPT - Increment size for lift gas injection rate
LOminEcGrad = 96, // LIFTOPT - Minimum economic gradient
GRpar_int = 97, // Guiderate parameter delay interval
ThrUPT = 99,
XxxDPR = 100,

2
opm/output/eclipse/VectorItems/intehead.hpp
View File

@@ -80,6 +80,8 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems
NGCTRL = 51, // Index indicating if group control is used or not (1 - if group control, 0 if not)
NGRNPH = 58, // Index indicating if group control is used or not (1 - if group control, 0 if not)
EACHNCITS = 59, // Index indicating if lift gas distribution optimized each of the NUPCOL first iterations or not
// 1 - optimized only first newton iteration, 2 - optimized each of NUPCOL newton iterations
DAY = 64, // Calendar day of report step (1..31)
MONTH = 65, // Calendar month of report step (1..12)

7
opm/output/eclipse/VectorItems/well.hpp
View File

@@ -50,6 +50,8 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems
HistReqWCtrl = 49, // Well's requested control mode from
// simulation deck (WCONHIST, WCONINJH)
LiftOpt = 53, // Well's lift gas injection to be calculated by optimisation or not
MsWID = 70, // Multisegment well ID
// Value 0 for regular wells
// Value 1..#MS wells for MS wells
@@ -168,6 +170,11 @@ namespace Opm { namespace RestartIO { namespace Helpers { namespace VectorItems
HistBHPTarget = 55, // Well's historical/observed bottom
// hole pressure target/limit
LOmaxRate = 56, // Well's maximum lift gas rate
LOweightFac = 57, // Well's wighting factor for preferential allocation of lift gas
LOminRate = 67, // Well's mimimum lift gas rate
};
} // SWell

3
opm/parser/eclipse/EclipseState/Schedule/GasLiftOpt.hpp
View File

@@ -204,7 +204,7 @@ public:
std::string m_name;
std::optional<double> m_max_rate;
double m_min_rate = 0;
bool m_use_glo;
bool m_use_glo = false;
double m_weight = 1;
double m_inc_weight = 0;
bool m_alloc_extra_gas = false;
@@ -219,6 +219,7 @@ public:
void gaslift_increment(double gaslift_increment);
double min_eco_gradient() const;
void min_eco_gradient(double min_eco_gradient);
double min_wait() const;
void min_wait(double min_wait);
void all_newton(double all_newton);
bool all_newton() const;

32
src/opm/output/eclipse/AggregateWellData.cpp
View File

@@ -268,6 +268,7 @@ namespace {
template <class IWellArray>
void staticContrib(const Opm::Well& well,
const Opm::GasLiftOpt& glo,
const Opm::SummaryState& st,
const std::size_t msWellID,
const std::map <const std::string, size_t>& GroupMapNameInd,
@@ -316,6 +317,14 @@ namespace {
iWell[Ix::item32] = 7;
iWell[Ix::item48] = - 1;
// integer flag indicating lift gas optimisation to be calculated
if (glo.has_well(well.name())) {
const auto& w_glo = glo.well(well.name());
iWell[Ix::LiftOpt] = (w_glo.use_glo()) ? 1 : 0;
} else {
iWell[Ix::LiftOpt] = 0;
}
// Deliberate misrepresentation. Function 'eclipseControlMode'
// returns the target control mode requested in the simulation
// deck. This item is supposed to be the well's actual, active
@@ -499,6 +508,7 @@ namespace {
};
template <class SWellArray>
void staticContrib(const Opm::Well& well,
const Opm::GasLiftOpt& glo,
const Opm::UnitSystem& units,
const std::size_t sim_step,
const Opm::Schedule& sched,
@@ -644,6 +654,19 @@ namespace {
sWell[Ix::HistBHPTarget] = sWell[Ix::BHPTarget];
}
// assign gas lift data
if (glo.has_well(well.name())) {
const auto& w_glo = glo.well(well.name());
sWell[Ix::LOmaxRate] = swprop(M::gas_surface_rate, w_glo.max_rate().value_or(0.));
sWell[Ix::LOminRate] = swprop(M::gas_surface_rate, w_glo.min_rate());
sWell[Ix::LOweightFac] = static_cast<float>(w_glo.weight_factor());
} else {
sWell[Ix::LOmaxRate] = 0.;
sWell[Ix::LOminRate] = 0.;
sWell[Ix::LOweightFac] = 0.;
}
sWell[Ix::DatumDepth] = swprop(M::length, datumDepth(well));
sWell[Ix::DrainageRadius] = swprop(M::length, well.getDrainageRadius());
sWell[Ix::EfficiencyFactor1] = well.getEfficiencyFactor();
@@ -953,6 +976,7 @@ captureDeclaredWellData(const Schedule& sched,
const std::vector<int>& inteHead)
{
const auto& wells = sched.getWells(sim_step);
const auto& step_glo = sched.glo(sim_step);
// Static contributions to IWEL array.
{
@@ -960,23 +984,23 @@ captureDeclaredWellData(const Schedule& sched,
const auto groupMapNameIndex = IWell::currentGroupMapNameIndex(sched, sim_step, inteHead);
auto msWellID = std::size_t{0};
wellLoop(wells, [&groupMapNameIndex, &msWellID, &smry, this]
wellLoop(wells, [&groupMapNameIndex, &msWellID, &step_glo, &smry, this]
(const Well& well, const std::size_t wellID) -> void
{
msWellID += well.isMultiSegment(); // 1-based index.
auto iw = this->iWell_[wellID];
IWell::staticContrib(well, smry, msWellID, groupMapNameIndex, iw);
IWell::staticContrib(well, step_glo, smry, msWellID, groupMapNameIndex, iw);
});
}
// Static contributions to SWEL array.
wellLoop(wells, [&units, &sim_step, &sched, &smry, this]
wellLoop(wells, [&units, &step_glo, &sim_step, &sched, &smry, this]
(const Well& well, const std::size_t wellID) -> void
{
auto sw = this->sWell_[wellID];
SWell::staticContrib(well, units, sim_step, sched, smry, sw);
SWell::staticContrib(well, step_glo, units, sim_step, sched, smry, sw);
});
// Static contributions to XWEL array.

14
src/opm/output/eclipse/CreateDoubHead.cpp
View File

@@ -107,6 +107,19 @@ namespace {
damping_fact
};
}
Opm::RestartIO::DoubHEAD::liftOptPar
computeLiftOptParam(const ::Opm::Schedule& sched,
const Opm::UnitSystem& units,
const std::size_t lookup_step)
{
using M = ::Opm::UnitSystem::measure;
const auto& glo = sched.glo(lookup_step);
return {
units.from_si(M::time, glo.min_wait()), units.from_si(M::gas_surface_rate, glo.gaslift_increment()),
units.from_si(M::oil_gas_ratio, glo.min_eco_gradient())
};
}
} // Anonymous
// #####################################################################
@@ -131,6 +144,7 @@ createDoubHead(const EclipseState& es,
.drsdt (sched, lookup_step, tconv)
.udq_param(rspec.udqParams())
.guide_rate_param(computeGuideRate(sched, lookup_step))
.lift_opt_param(computeLiftOptParam(sched, usys, lookup_step))
;
if (nextTimeStep > 0.0) {

9
src/opm/output/eclipse/CreateInteHead.cpp
View File

@@ -447,6 +447,14 @@ namespace {
return mode;
}
int getLiftOptPar(const ::Opm::Schedule& sched,
const std::size_t lookup_step)
{
const auto& each_nupcol = sched.glo(lookup_step).all_newton();
int in_enc = (each_nupcol) ? 2 : 1;
return in_enc;
}
} // Anonymous
// #####################################################################
@@ -494,6 +502,7 @@ createInteHead(const EclipseState& es,
.params_NAAQZ (1, 18, 24, 10, 7, 2, 4)
.stepParam (num_solver_steps, report_step)
.tuningParam (getTuningPars(sched.getTuning(lookup_step)))
.liftOptParam (getLiftOptPar(sched, lookup_step))
.wellSegDimensions (getWellSegDims(rspec, sched, report_step, lookup_step))
.regionDimensions (getRegDims(tdim, rdim))
.ngroups ({ ngmax })

21
src/opm/output/eclipse/DoubHEAD.cpp
View File

@@ -154,10 +154,10 @@ enum Index : std::vector<double>::size_type {
grpar_d = VI::doubhead::GRpar_d,
grpar_e = VI::doubhead::GRpar_e,
grpar_f = VI::doubhead::GRpar_f,
dh_093 = 93,
dh_094 = 94,
LOmini = VI::doubhead::LOminInt,
LOincr = VI::doubhead::LOincrsz,
dh_095 = 95,
dh_096 = 96,
LOmineg = VI::doubhead::LOminEcGrad,
grpar_int = VI::doubhead::GRpar_int,
dh_098 = 98,
ThrUPT = VI::doubhead::ThrUPT,
@@ -398,8 +398,9 @@ Opm::RestartIO::DoubHEAD::DoubHEAD()
this->data_[Index::dh_081] = 1.0e+20;
this->data_[grpar_e] = 0.0;
this->data_[grpar_f] = 0.0;
this->data_[Index::dh_093] = 0.0;
this->data_[Index::dh_096] = 0.0;
this->data_[LOmini] = 0.0;
this->data_[LOincr] = 0.0;
this->data_[LOmineg] = 0.0;
this->data_[Index::dh_105] = 1.0;
this->data_[Index::dh_108] = 0.0;
this->data_[Index::dh_109] = 0.0;
@@ -643,3 +644,13 @@ Opm::RestartIO::DoubHEAD::guide_rate_param(const guideRate& guide_rp)
return *this;
}
Opm::RestartIO::DoubHEAD&
Opm::RestartIO::DoubHEAD::lift_opt_param(const liftOptPar& lo_par)
{
this->data_[LOmini] = lo_par.min_int;
this->data_[LOincr] = lo_par.incr;
this->data_[LOmineg] = lo_par.min_ec_grad;
return *this;
}

11
src/opm/output/eclipse/InteHEAD.cpp
View File

@@ -79,7 +79,7 @@ enum index : std::vector<int>::size_type {
ih_056 = 56 , // 0 0
ih_057 = 57 , // 0 0
NGRNPHASE = VI::intehead::NGRNPH, // Parameter to determine the nominated phase for the guiderate
ih_059 = 59 , // 0 0
EACHNC = VI::intehead::EACHNCITS, // Index indicating if lift gas distribution optimized each of the NUPCOL first iterations or not
ih_060 = 60 , // 0 0
ih_061 = 61 , // 0 0
ih_062 = 62 , // 0 0
@@ -738,6 +738,15 @@ whistControlMode(int mode)
return *this;
}
Opm::RestartIO::InteHEAD&
Opm::RestartIO::InteHEAD::
liftOptParam(int in_enc)
{
this -> data_[EACHNC] = in_enc;
return *this;
}
// =====================================================================
// Free functions (calendar/time utilities)
// =====================================================================

4
src/opm/parser/eclipse/EclipseState/Schedule/GasLiftOpt.cpp
View File

@@ -47,6 +47,10 @@ void GasLiftOpt::min_wait(double min_wait) {
this->m_min_wait = min_wait;
}
double GasLiftOpt::min_wait() const {
return this->m_min_wait;
}
void GasLiftOpt::all_newton(double all_newton) {
this->m_all_newton = all_newton;
}

51025
tests/2_WLIFT_MODEL5_NOINC.DATA Normal file
View File

File diff suppressed because it is too large Load Diff

175
tests/test_LGOData.cpp Normal file
View File

@@ -0,0 +1,175 @@
#define BOOST_TEST_MODULE UDQ_Data
#include <boost/test/unit_test.hpp>
#include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Action/State.hpp>
#include <opm/parser/eclipse/EclipseState/Schedule/Well/Well.hpp>
#include <opm/parser/eclipse/Parser/Parser.hpp>
#include <opm/parser/eclipse/Deck/Deck.hpp>
#include <opm/output/eclipse/AggregateWellData.hpp>
#include <opm/output/eclipse/WriteRestartHelpers.hpp>
#include <opm/output/eclipse/InteHEAD.hpp>
#include <opm/output/eclipse/VectorItems/intehead.hpp>
#include <opm/output/eclipse/VectorItems/doubhead.hpp>
#include <opm/output/eclipse/VectorItems/well.hpp>
#include <opm/output/eclipse/DoubHEAD.hpp>
#include <opm/parser/eclipse/Python/Python.hpp>
#include <opm/parser/eclipse/Units/UnitSystem.hpp>
#include <opm/parser/eclipse/Units/Units.hpp>
#include <stdexcept>
#include <utility>
#include <exception>
#include <iostream>
#include <string>
#include <vector>
namespace {
Opm::Deck first_sim(std::string fname) {
return Opm::Parser{}.parseFile(fname);
}
}
Opm::SummaryState sum_state()
{
auto state = Opm::SummaryState{std::chrono::system_clock::now()};
state.update("FULPR", 460.);
return state;
}
//int main(int argc, char* argv[])
struct SimulationCase
{
explicit SimulationCase(const Opm::Deck& deck)
: es { deck }
, grid { deck }
, python { std::make_shared<Opm::Python>()}
, sched{ deck, es, python }
{}
// Order requirement: 'es' must be declared/initialised before 'sched'.
Opm::EclipseState es;
Opm::EclipseGrid grid;
std::shared_ptr<Opm::Python> python;
Opm::Schedule sched;
};
BOOST_AUTO_TEST_SUITE(LiftGasOptimization)
// test lift gas optimisation data
BOOST_AUTO_TEST_CASE (liftGasOptimzation_data)
{
const auto simCase = SimulationCase{first_sim("2_WLIFT_MODEL5_NOINC.DATA")};
Opm::EclipseState es = simCase.es;
Opm::SummaryState st = sum_state();
Opm::Schedule sched = simCase.sched;
Opm::EclipseGrid grid = simCase.grid;
//const auto& ioConfig = es.getIOConfig();
Opm::Action::State action_state;
// Report Step 1: 2020-01-01
const auto rptStep = std::size_t{2};
const auto simStep = std::size_t{1};
double secs_elapsed = 3.1536E07;
const auto ih = Opm::RestartIO::Helpers::
createInteHead(es, grid, sched, secs_elapsed,
rptStep, rptStep, simStep);
//set dummy value for next_step_size
const double next_step_size= 0.1;
const auto dh = Opm::RestartIO::Helpers::createDoubHead(es, sched, simStep,
secs_elapsed, next_step_size);
auto wellData = Opm::RestartIO::Helpers::AggregateWellData(ih);
wellData.captureDeclaredWellData(sched, es.getUnits(), simStep, action_state, st, ih);
// intehead data
auto eachnc = Opm::RestartIO::Helpers::VectorItems::intehead::EACHNCITS;
auto niwelz = Opm::RestartIO::Helpers::VectorItems::intehead::NIWELZ;
auto nswelz = Opm::RestartIO::Helpers::VectorItems::intehead::NSWELZ;
// doubhead data
auto lomini = Opm::RestartIO::Helpers::VectorItems::doubhead::LOminInt;
auto loincr = Opm::RestartIO::Helpers::VectorItems::doubhead::LOincrsz;
auto lomineg = Opm::RestartIO::Helpers::VectorItems::doubhead::LOminEcGrad;
// iwel data
auto liftopt = static_cast<std::size_t>(Opm::RestartIO::Helpers::VectorItems::IWell::index::LiftOpt);
// swel data
auto lomaxrate = static_cast<std::size_t>(Opm::RestartIO::Helpers::VectorItems::SWell::index::LOmaxRate);
auto lominrate = static_cast<std::size_t>(Opm::RestartIO::Helpers::VectorItems::SWell::index::LOminRate);
auto loweightfac = static_cast<std::size_t>(Opm::RestartIO::Helpers::VectorItems::SWell::index::LOweightFac);
{
/*
Check of InteHEAD and DoubHEAD data for LIFTOPT variables
*/
BOOST_CHECK_EQUAL(ih[eachnc] , 2);
BOOST_CHECK_EQUAL(dh[lomini] , 37.);
BOOST_CHECK_EQUAL(dh[loincr] , 12500);
BOOST_CHECK_EQUAL(dh[lomineg] , 5E-3);
}
{
// IWEL
const auto& iWel = wellData.getIWell();
auto start = 0*static_cast<std::size_t>(ih[niwelz]);
BOOST_CHECK_EQUAL(iWel[start + liftopt] , 1);
start = 1*static_cast<std::size_t>(ih[niwelz]);
BOOST_CHECK_EQUAL(iWel[start + liftopt] , 0);
}
{
// SWEL
const auto& sWel = wellData.getSWell();
auto start = 0*static_cast<std::size_t>(ih[nswelz]);
BOOST_CHECK_CLOSE(sWel[start + lomaxrate] , 150000.f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + lominrate] , -1.0f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + loweightfac] , 1.01f, 1.0e-7f);
start = 1*static_cast<std::size_t>(ih[nswelz]);
BOOST_CHECK_CLOSE(sWel[start + lomaxrate] , 150000.f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + lominrate] , 0.0f, 1.0e-7f); // default value since item 2 for this well is 'NO'
BOOST_CHECK_CLOSE(sWel[start + loweightfac] , 1.0f, 1.0e-7f); // default value since item 2 for this well is 'NO'
start = 2*static_cast<std::size_t>(ih[nswelz]);
BOOST_CHECK_CLOSE(sWel[start + lomaxrate] , 150000.f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + lominrate] , 3.0f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + loweightfac] , 1.21f, 1.0e-7f);
start = 3*static_cast<std::size_t>(ih[nswelz]);
BOOST_CHECK_CLOSE(sWel[start + lomaxrate] , 150000.f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + lominrate] , 0.f, 1.0e-7f);
BOOST_CHECK_CLOSE(sWel[start + loweightfac] , 1.01f, 1.0e-7f);
}
}
BOOST_AUTO_TEST_SUITE_END()