2016-04-20 09:55:17 +02:00
|
|
|
/*
|
2019-10-09 00:25:45 -05:00
|
|
|
Copyright 2019 Equinor ASA.
|
2016-04-20 09:55:17 +02:00
|
|
|
Copyright 2016 Statoil ASA.
|
|
|
|
|
|
|
|
|
|
This file is part of the Open Porous Media project (OPM).
|
|
|
|
|
|
|
|
|
|
OPM is free software: you can redistribute it and/or modify
|
|
|
|
|
it under the terms of the GNU General Public License as published by
|
|
|
|
|
the Free Software Foundation, either version 3 of the License, or
|
|
|
|
|
(at your option) any later version.
|
|
|
|
|
|
|
|
|
|
OPM is distributed in the hope that it will be useful,
|
|
|
|
|
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
|
|
|
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
|
|
|
GNU General Public License for more details.
|
|
|
|
|
|
|
|
|
|
You should have received a copy of the GNU General Public License
|
|
|
|
|
along with OPM. If not, see <http://www.gnu.org/licenses/>.
|
|
|
|
|
*/
|
2016-10-28 14:42:07 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
#include <opm/output/eclipse/Summary.hpp>
|
2019-03-08 07:12:07 +01:00
|
|
|
|
2017-06-22 22:07:59 +02:00
|
|
|
#include <opm/common/OpmLog/OpmLog.hpp>
|
2019-11-07 12:03:38 +01:00
|
|
|
#include <opm/common/OpmLog/Location.hpp>
|
2017-06-22 22:07:59 +02:00
|
|
|
|
2016-04-20 09:55:17 +02:00
|
|
|
#include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
|
|
|
|
|
#include <opm/parser/eclipse/EclipseState/IOConfig/IOConfig.hpp>
|
2016-11-01 11:19:59 +01:00
|
|
|
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
|
2019-11-12 08:29:28 +01:00
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/Group/Group.hpp>
|
2016-11-01 11:19:59 +01:00
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/Schedule.hpp>
|
2019-10-09 00:25:45 -05:00
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/SummaryState.hpp>
|
|
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/UDQ/UDQConfig.hpp>
|
|
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/UDQ/UDQContext.hpp>
|
|
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/Well/WellProductionProperties.hpp>
|
|
|
|
|
#include <opm/parser/eclipse/EclipseState/Schedule/Well/WellInjectionProperties.hpp>
|
2016-04-20 09:55:17 +02:00
|
|
|
#include <opm/parser/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp>
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2016-04-20 09:55:17 +02:00
|
|
|
#include <opm/parser/eclipse/Units/UnitSystem.hpp>
|
2019-10-09 00:25:45 -05:00
|
|
|
#include <opm/parser/eclipse/Units/Units.hpp>
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
#include <opm/io/eclipse/EclOutput.hpp>
|
|
|
|
|
#include <opm/io/eclipse/OutputStream.hpp>
|
|
|
|
|
|
|
|
|
|
#include <opm/output/data/Wells.hpp>
|
2020-03-16 13:51:05 +01:00
|
|
|
#include <opm/output/data/Groups.hpp>
|
2016-10-11 14:39:15 +02:00
|
|
|
#include <opm/output/eclipse/RegionCache.hpp>
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
#include <algorithm>
|
|
|
|
|
#include <array>
|
|
|
|
|
#include <cassert>
|
|
|
|
|
#include <chrono>
|
|
|
|
|
#include <cstddef>
|
|
|
|
|
#include <cctype>
|
|
|
|
|
#include <ctime>
|
|
|
|
|
#include <exception>
|
|
|
|
|
#include <functional>
|
|
|
|
|
#include <initializer_list>
|
|
|
|
|
#include <iterator>
|
|
|
|
|
#include <memory>
|
|
|
|
|
#include <numeric>
|
2020-03-19 13:28:17 +01:00
|
|
|
#include <regex>
|
2019-10-09 00:25:45 -05:00
|
|
|
#include <stdexcept>
|
|
|
|
|
#include <string>
|
|
|
|
|
#include <unordered_map>
|
|
|
|
|
#include <utility>
|
|
|
|
|
#include <vector>
|
2019-05-19 15:39:03 +02:00
|
|
|
|
2018-06-25 13:47:34 +02:00
|
|
|
namespace {
|
2020-03-19 10:53:46 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
struct ParamCTorArgs
|
2018-10-09 13:59:13 +02:00
|
|
|
{
|
2019-10-09 00:25:45 -05:00
|
|
|
std::string kw;
|
2020-03-13 14:24:29 +01:00
|
|
|
Opm::SummaryConfigNode::Type type;
|
2018-10-09 13:59:13 +02:00
|
|
|
};
|
|
|
|
|
|
2020-03-16 13:51:05 +01:00
|
|
|
using p_cmode = Opm::Group::ProductionCMode;
|
|
|
|
|
const std::map<p_cmode, int> pCModeToPCntlMode = {
|
|
|
|
|
{p_cmode::NONE, 0},
|
|
|
|
|
{p_cmode::ORAT, 1},
|
|
|
|
|
{p_cmode::WRAT, 2},
|
|
|
|
|
{p_cmode::GRAT, 3},
|
|
|
|
|
{p_cmode::LRAT, 4},
|
|
|
|
|
{p_cmode::CRAT, 9},
|
|
|
|
|
{p_cmode::RESV, 5},
|
|
|
|
|
{p_cmode::PRBL, 6},
|
|
|
|
|
{p_cmode::FLD, 0}, // same as NONE
|
|
|
|
|
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
using i_cmode = Opm::Group::InjectionCMode;
|
|
|
|
|
const std::map<i_cmode, int> iCModeToICntlMode = {
|
|
|
|
|
{i_cmode::NONE, 0},
|
|
|
|
|
{i_cmode::RATE, 1},
|
|
|
|
|
{i_cmode::RESV, 2},
|
|
|
|
|
{i_cmode::REIN, 3},
|
|
|
|
|
{i_cmode::VREP, 4},
|
|
|
|
|
{i_cmode::FLD, 0}, // same as NONE
|
|
|
|
|
{i_cmode::SALE, 0}, // not used in E100
|
|
|
|
|
};
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
std::vector<ParamCTorArgs> requiredRestartVectors()
|
2018-06-25 13:47:34 +02:00
|
|
|
{
|
2020-03-13 14:24:29 +01:00
|
|
|
using Type = ::Opm::SummaryConfigNode::Type;
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2018-06-25 13:47:34 +02:00
|
|
|
return {
|
2019-10-09 00:25:45 -05:00
|
|
|
// Production
|
2020-03-19 14:01:29 +01:00
|
|
|
{ "OPR" , Type::Rate },
|
|
|
|
|
{ "WPR" , Type::Rate },
|
|
|
|
|
{ "GPR" , Type::Rate },
|
|
|
|
|
{ "VPR" , Type::Rate },
|
|
|
|
|
{ "OPP" , Type::Rate },
|
|
|
|
|
{ "WPP" , Type::Rate },
|
|
|
|
|
{ "GPP" , Type::Rate },
|
2020-03-19 10:49:41 +01:00
|
|
|
{ "OPT" , Type::Total },
|
|
|
|
|
{ "WPT" , Type::Total },
|
|
|
|
|
{ "GPT" , Type::Total },
|
|
|
|
|
{ "VPT" , Type::Total },
|
|
|
|
|
{ "OPTH", Type::Total },
|
|
|
|
|
{ "WPTH", Type::Total },
|
|
|
|
|
{ "GPTH", Type::Total },
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
// Flow rate ratios (production)
|
2020-03-19 10:49:41 +01:00
|
|
|
{ "WCT" , Type::Ratio },
|
|
|
|
|
{ "GOR" , Type::Ratio },
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
// injection
|
2020-03-19 10:49:41 +01:00
|
|
|
|
2020-03-19 14:01:29 +01:00
|
|
|
{ "WIR" , Type::Rate },
|
|
|
|
|
{ "GIR" , Type::Rate },
|
|
|
|
|
{ "OPI" , Type::Rate },
|
|
|
|
|
{ "WPI" , Type::Rate },
|
|
|
|
|
{ "GPI" , Type::Rate },
|
2020-03-19 10:49:41 +01:00
|
|
|
{ "WIT" , Type::Total },
|
|
|
|
|
{ "GIT" , Type::Total },
|
|
|
|
|
{ "WITH", Type::Total },
|
|
|
|
|
{ "GITH", Type::Total },
|
2018-06-25 13:47:34 +02:00
|
|
|
};
|
|
|
|
|
}
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
std::vector<Opm::SummaryConfigNode>
|
2018-06-25 13:47:34 +02:00
|
|
|
requiredRestartVectors(const ::Opm::Schedule& sched)
|
|
|
|
|
{
|
2020-03-19 11:12:59 +01:00
|
|
|
std::vector<Opm::SummaryConfigNode> entities {} ;
|
2020-03-19 11:52:46 +01:00
|
|
|
const auto& restartVectors { requiredRestartVectors() } ;
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
using SN = ::Opm::SummaryConfigNode;
|
2018-06-25 13:47:34 +02:00
|
|
|
|
2020-03-19 11:52:46 +01:00
|
|
|
auto makeEntities = [&restartVectors, &entities]
|
2019-10-09 00:25:45 -05:00
|
|
|
(const char kwpref,
|
|
|
|
|
const SN::Category cat,
|
|
|
|
|
const std::string& name) -> void
|
2018-06-25 13:47:34 +02:00
|
|
|
{
|
2020-03-19 11:52:46 +01:00
|
|
|
for (const auto& restartVector : restartVectors) {
|
|
|
|
|
entities.emplace_back(kwpref + restartVector.first, cat, ::Opm::Location());
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
entities.back().namedEntity(name)
|
2020-03-19 11:52:46 +01:00
|
|
|
.parameterType(restartVector.second);
|
2018-06-25 13:47:34 +02:00
|
|
|
}
|
|
|
|
|
};
|
|
|
|
|
|
2019-04-19 07:24:10 +02:00
|
|
|
for (const auto& well_name : sched.wellNames()) {
|
2019-10-09 00:25:45 -05:00
|
|
|
makeEntities('W', SN::Category::Well, well_name);
|
2018-06-25 13:47:34 +02:00
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
entities.emplace_back("WBHP", SN::Category::Well, ::Opm::Location());
|
2019-10-09 00:25:45 -05:00
|
|
|
entities.back().namedEntity(well_name)
|
|
|
|
|
.parameterType(SN::Type::Pressure);
|
|
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
entities.emplace_back("WGVIR", SN::Category::Well, ::Opm::Location());
|
2019-10-09 00:25:45 -05:00
|
|
|
entities.back().namedEntity(well_name)
|
|
|
|
|
.parameterType(SN::Type::Rate);
|
|
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
entities.emplace_back("WWVIR", SN::Category::Well, ::Opm::Location());
|
2019-10-09 00:25:45 -05:00
|
|
|
entities.back().namedEntity(well_name)
|
|
|
|
|
.parameterType(SN::Type::Rate);
|
2018-06-25 13:47:34 +02:00
|
|
|
}
|
|
|
|
|
|
2019-07-07 20:14:41 +02:00
|
|
|
for (const auto& grp_name : sched.groupNames()) {
|
2020-03-16 13:51:05 +01:00
|
|
|
if (grp_name != "FIELD") {
|
2019-10-09 00:25:45 -05:00
|
|
|
makeEntities('G', SN::Category::Group, grp_name);
|
2020-03-16 13:51:05 +01:00
|
|
|
|
|
|
|
|
entities.emplace_back("GMCTP", SN::Category::Group, ::Opm::Location());
|
|
|
|
|
entities.back().namedEntity(grp_name)
|
|
|
|
|
.parameterType(SN::Type::Mode);
|
|
|
|
|
|
|
|
|
|
entities.emplace_back("GMCTW", SN::Category::Group, ::Opm::Location());
|
|
|
|
|
entities.back().namedEntity(grp_name)
|
|
|
|
|
.parameterType(SN::Type::Mode);
|
|
|
|
|
|
|
|
|
|
entities.emplace_back("GMCTG", SN::Category::Group, ::Opm::Location());
|
|
|
|
|
entities.back().namedEntity(grp_name)
|
|
|
|
|
.parameterType(SN::Type::Mode);
|
|
|
|
|
}
|
2018-06-25 13:47:34 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
makeEntities('F', SN::Category::Field, "FIELD");
|
2018-06-25 13:47:34 +02:00
|
|
|
|
2020-03-16 13:51:05 +01:00
|
|
|
entities.emplace_back("FMCTP", SN::Category::Field, ::Opm::Location());
|
|
|
|
|
entities.back().namedEntity("FIELD")
|
|
|
|
|
.parameterType(SN::Type::Mode);
|
|
|
|
|
|
|
|
|
|
entities.emplace_back("FMCTW", SN::Category::Field, ::Opm::Location());
|
|
|
|
|
entities.back().namedEntity("FIELD")
|
|
|
|
|
.parameterType(SN::Type::Mode);
|
|
|
|
|
|
|
|
|
|
entities.emplace_back("FMCTG", SN::Category::Field, ::Opm::Location());
|
|
|
|
|
entities.back().namedEntity("FIELD")
|
|
|
|
|
.parameterType(SN::Type::Mode);
|
|
|
|
|
|
2018-06-25 13:47:34 +02:00
|
|
|
return entities;
|
|
|
|
|
}
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
std::vector<Opm::SummaryConfigNode>
|
2018-10-09 13:59:13 +02:00
|
|
|
requiredSegmentVectors(const ::Opm::Schedule& sched)
|
|
|
|
|
{
|
2020-03-13 14:24:29 +01:00
|
|
|
using SN = Opm::SummaryConfigNode;
|
2019-10-09 00:25:45 -05:00
|
|
|
auto ret = std::vector<SN>{};
|
|
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
auto sofr = SN{ "SOFR", SN::Category::Segment, ::Opm::Location() }
|
2019-10-09 00:25:45 -05:00
|
|
|
.parameterType(SN::Type::Rate);
|
|
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
auto sgfr = SN{ "SGFR", SN::Category::Segment, ::Opm::Location() }
|
2019-10-09 00:25:45 -05:00
|
|
|
.parameterType(SN::Type::Rate);
|
|
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
auto swfr = SN{ "SWFR", SN::Category::Segment, ::Opm::Location() }
|
2019-10-09 00:25:45 -05:00
|
|
|
.parameterType(SN::Type::Rate);
|
|
|
|
|
|
2020-01-03 08:20:20 +01:00
|
|
|
auto spr = SN{ "SPR", SN::Category::Segment, ::Opm::Location() }
|
2019-10-09 00:25:45 -05:00
|
|
|
.parameterType(SN::Type::Pressure);
|
2018-10-09 13:59:13 +02:00
|
|
|
|
|
|
|
|
auto makeVectors =
|
2019-10-09 00:25:45 -05:00
|
|
|
[&](const std::string& well,
|
|
|
|
|
const int segNumber) -> void
|
2018-10-09 13:59:13 +02:00
|
|
|
{
|
2019-10-09 00:25:45 -05:00
|
|
|
ret.push_back(sofr.namedEntity(well).number(segNumber));
|
|
|
|
|
ret.push_back(sgfr.namedEntity(well).number(segNumber));
|
|
|
|
|
ret.push_back(swfr.namedEntity(well).number(segNumber));
|
|
|
|
|
ret.push_back(spr .namedEntity(well).number(segNumber));
|
2018-10-09 13:59:13 +02:00
|
|
|
};
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
for (const auto& wname : sched.wellNames()) {
|
2019-11-12 08:29:28 +01:00
|
|
|
const auto& well = sched.getWellatEnd(wname);
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
if (! well.isMultiSegment()) {
|
2018-10-09 13:59:13 +02:00
|
|
|
// Don't allocate MS summary vectors for non-MS wells.
|
|
|
|
|
continue;
|
|
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
const auto nSeg = well.getSegments().size();
|
2018-10-09 13:59:13 +02:00
|
|
|
|
|
|
|
|
for (auto segID = 0*nSeg; segID < nSeg; ++segID) {
|
|
|
|
|
makeVectors(wname, segID + 1); // One-based
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
return ret;
|
|
|
|
|
}
|
|
|
|
|
|
2016-04-20 09:55:17 +02:00
|
|
|
/*
|
2016-06-15 10:45:07 +02:00
|
|
|
* This class takes simulator state and parser-provided information and
|
|
|
|
|
* orchestrates ert to write simulation results as requested by the SUMMARY
|
|
|
|
|
* section in eclipse. The implementation is somewhat compact as a lot of the
|
|
|
|
|
* requested output may be similar-but-not-quite. Through various techniques
|
|
|
|
|
* the compiler writes a lot of this code for us.
|
2016-04-20 09:55:17 +02:00
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
using rt = Opm::data::Rates::opt;
|
|
|
|
|
using measure = Opm::UnitSystem::measure;
|
2016-11-02 12:35:11 +01:00
|
|
|
constexpr const bool injector = true;
|
|
|
|
|
constexpr const bool producer = false;
|
2018-06-22 14:09:08 +02:00
|
|
|
constexpr const bool polymer = true;
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
/* Some numerical value with its unit tag embedded to enable caller to apply
|
2016-07-19 13:18:06 +02:00
|
|
|
* unit conversion. This removes a lot of boilerplate. ad-hoc solution to poor
|
|
|
|
|
* unit support in general.
|
2016-06-15 10:45:07 +02:00
|
|
|
*/
|
2016-07-19 13:18:06 +02:00
|
|
|
measure div_unit( measure denom, measure div ) {
|
|
|
|
|
if( denom == measure::gas_surface_rate &&
|
|
|
|
|
div == measure::liquid_surface_rate )
|
|
|
|
|
return measure::gas_oil_ratio;
|
|
|
|
|
|
|
|
|
|
if( denom == measure::liquid_surface_rate &&
|
|
|
|
|
div == measure::gas_surface_rate )
|
|
|
|
|
return measure::oil_gas_ratio;
|
|
|
|
|
|
|
|
|
|
if( denom == measure::liquid_surface_rate &&
|
|
|
|
|
div == measure::liquid_surface_rate )
|
|
|
|
|
return measure::water_cut;
|
|
|
|
|
|
|
|
|
|
if( denom == measure::liquid_surface_rate &&
|
|
|
|
|
div == measure::time )
|
|
|
|
|
return measure::liquid_surface_volume;
|
|
|
|
|
|
|
|
|
|
if( denom == measure::gas_surface_rate &&
|
|
|
|
|
div == measure::time )
|
|
|
|
|
return measure::gas_surface_volume;
|
|
|
|
|
|
2018-06-22 14:09:08 +02:00
|
|
|
if( denom == measure::mass_rate &&
|
|
|
|
|
div == measure::time )
|
|
|
|
|
return measure::mass;
|
|
|
|
|
|
2016-07-19 13:18:06 +02:00
|
|
|
return measure::identity;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
measure mul_unit( measure lhs, measure rhs ) {
|
|
|
|
|
if( lhs == rhs ) return lhs;
|
|
|
|
|
|
|
|
|
|
if( ( lhs == measure::liquid_surface_rate && rhs == measure::time ) ||
|
|
|
|
|
( rhs == measure::liquid_surface_rate && lhs == measure::time ) )
|
|
|
|
|
return measure::liquid_surface_volume;
|
|
|
|
|
|
|
|
|
|
if( ( lhs == measure::gas_surface_rate && rhs == measure::time ) ||
|
|
|
|
|
( rhs == measure::gas_surface_rate && lhs == measure::time ) )
|
|
|
|
|
return measure::gas_surface_volume;
|
|
|
|
|
|
2017-12-14 10:16:03 +01:00
|
|
|
if( ( lhs == measure::rate && rhs == measure::time ) ||
|
|
|
|
|
( rhs == measure::rate && lhs == measure::time ) )
|
|
|
|
|
return measure::volume;
|
|
|
|
|
|
2018-06-22 14:09:08 +02:00
|
|
|
if( lhs == measure::mass_rate && rhs == measure::time)
|
|
|
|
|
return measure::mass;
|
|
|
|
|
|
2016-07-19 13:18:06 +02:00
|
|
|
return lhs;
|
|
|
|
|
}
|
|
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
struct quantity {
|
|
|
|
|
double value;
|
2019-10-07 13:32:31 +02:00
|
|
|
Opm::UnitSystem::measure unit;
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
quantity operator+( const quantity& rhs ) const {
|
2016-07-19 13:18:06 +02:00
|
|
|
assert( this->unit == rhs.unit );
|
2016-06-15 10:45:07 +02:00
|
|
|
return { this->value + rhs.value, this->unit };
|
|
|
|
|
}
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
quantity operator*( const quantity& rhs ) const {
|
2016-07-19 13:18:06 +02:00
|
|
|
return { this->value * rhs.value, mul_unit( this->unit, rhs.unit ) };
|
2016-06-15 10:45:07 +02:00
|
|
|
}
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
quantity operator/( const quantity& rhs ) const {
|
2016-07-19 13:18:06 +02:00
|
|
|
const auto res_unit = div_unit( this->unit, rhs.unit );
|
|
|
|
|
|
|
|
|
|
if( rhs.value == 0 ) return { 0.0, res_unit };
|
|
|
|
|
return { this->value / rhs.value, res_unit };
|
2016-06-15 10:45:07 +02:00
|
|
|
}
|
2016-10-24 12:44:56 +02:00
|
|
|
|
|
|
|
|
quantity operator/( double divisor ) const {
|
|
|
|
|
if( divisor == 0 ) return { 0.0, this->unit };
|
|
|
|
|
return { this->value / divisor , this->unit };
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
quantity& operator/=( double divisor ) {
|
|
|
|
|
if( divisor == 0 )
|
|
|
|
|
this->value = 0;
|
|
|
|
|
else
|
|
|
|
|
this->value /= divisor;
|
|
|
|
|
|
|
|
|
|
return *this;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
2016-09-23 15:57:10 +02:00
|
|
|
quantity operator-( const quantity& rhs) const {
|
|
|
|
|
return { this->value - rhs.value, this->unit };
|
|
|
|
|
}
|
2016-06-15 10:45:07 +02:00
|
|
|
};
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2020-03-16 13:51:05 +01:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
/*
|
|
|
|
|
* All functions must have the same parameters, so they're gathered in a struct
|
|
|
|
|
* and functions use whatever information they care about.
|
|
|
|
|
*
|
2016-07-19 08:46:09 +02:00
|
|
|
* schedule_wells are wells from the deck, provided by opm-parser. active_index
|
|
|
|
|
* is the index of the block in question. wells is simulation data.
|
2016-06-15 10:45:07 +02:00
|
|
|
*/
|
|
|
|
|
struct fn_args {
|
2019-11-12 08:29:28 +01:00
|
|
|
const std::vector<Opm::Well>& schedule_wells;
|
2020-03-16 13:51:05 +01:00
|
|
|
const std::string group_name;
|
2016-06-15 10:45:07 +02:00
|
|
|
double duration;
|
2018-03-07 08:57:58 +01:00
|
|
|
const int sim_step;
|
2016-08-23 14:09:44 +02:00
|
|
|
int num;
|
2019-10-07 13:32:31 +02:00
|
|
|
const Opm::SummaryState& st;
|
|
|
|
|
const Opm::data::Wells& wells;
|
2020-03-16 13:51:05 +01:00
|
|
|
const Opm::data::Group& group;
|
2019-10-07 13:32:31 +02:00
|
|
|
const Opm::out::RegionCache& regionCache;
|
|
|
|
|
const Opm::EclipseGrid& grid;
|
2018-01-31 08:16:31 +01:00
|
|
|
const std::vector< std::pair< std::string, double > > eff_factors;
|
2016-06-15 10:45:07 +02:00
|
|
|
};
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
/* Since there are several enums in opm scattered about more-or-less
|
|
|
|
|
* representing the same thing. Since functions use template parameters to
|
|
|
|
|
* expand into the actual implementations we need a static way to determine
|
|
|
|
|
* what unit to tag the return value with.
|
|
|
|
|
*/
|
|
|
|
|
template< rt > constexpr
|
|
|
|
|
measure rate_unit() { return measure::liquid_surface_rate; }
|
2019-10-07 13:32:31 +02:00
|
|
|
template< Opm::Phase > constexpr
|
2016-06-15 10:45:07 +02:00
|
|
|
measure rate_unit() { return measure::liquid_surface_rate; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::gas >() { return measure::gas_surface_rate; }
|
|
|
|
|
template<> constexpr
|
2019-10-07 13:32:31 +02:00
|
|
|
measure rate_unit< Opm::Phase::GAS >() { return measure::gas_surface_rate; }
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2016-09-23 15:57:10 +02:00
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::solvent >() { return measure::gas_surface_rate; }
|
|
|
|
|
|
2017-12-14 10:16:03 +01:00
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::reservoir_water >() { return measure::rate; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::reservoir_oil >() { return measure::rate; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::reservoir_gas >() { return measure::rate; }
|
|
|
|
|
|
2018-11-02 15:21:28 +01:00
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit < rt::productivity_index_water > () { return measure::liquid_productivity_index; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit < rt::productivity_index_oil > () { return measure::liquid_productivity_index; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit < rt::productivity_index_gas > () { return measure::gas_productivity_index; }
|
|
|
|
|
|
2018-11-07 09:14:06 +01:00
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::well_potential_water >() { return measure::liquid_surface_rate; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::well_potential_oil >() { return measure::liquid_surface_rate; }
|
|
|
|
|
|
|
|
|
|
template<> constexpr
|
|
|
|
|
measure rate_unit< rt::well_potential_gas >() { return measure::gas_surface_rate; }
|
|
|
|
|
|
2018-01-31 08:16:31 +01:00
|
|
|
double efac( const std::vector<std::pair<std::string,double>>& eff_factors, const std::string& name ) {
|
|
|
|
|
auto it = std::find_if( eff_factors.begin(), eff_factors.end(),
|
|
|
|
|
[&] ( const std::pair< std::string, double > elem )
|
|
|
|
|
{ return elem.first == name; }
|
|
|
|
|
);
|
|
|
|
|
|
|
|
|
|
return (it != eff_factors.end()) ? it->second : 1;
|
|
|
|
|
}
|
|
|
|
|
|
2018-06-22 14:09:08 +02:00
|
|
|
template< rt phase, bool injection = true, bool polymer = false >
|
2016-06-15 10:45:07 +02:00
|
|
|
inline quantity rate( const fn_args& args ) {
|
|
|
|
|
double sum = 0.0;
|
|
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
for( const auto& sched_well : args.schedule_wells ) {
|
|
|
|
|
const auto& name = sched_well.name();
|
2016-09-25 18:31:34 +02:00
|
|
|
if( args.wells.count( name ) == 0 ) continue;
|
2018-01-31 08:16:31 +01:00
|
|
|
|
|
|
|
|
double eff_fac = efac( args.eff_factors, name );
|
|
|
|
|
|
2018-06-22 14:09:08 +02:00
|
|
|
double concentration = polymer
|
2019-05-04 12:00:32 +02:00
|
|
|
? sched_well.getPolymerProperties().m_polymerConcentration
|
2018-06-22 14:09:08 +02:00
|
|
|
: 1;
|
|
|
|
|
|
|
|
|
|
const auto v = args.wells.at(name).rates.get(phase, 0.0) * eff_fac * concentration;
|
2018-01-31 08:16:31 +01:00
|
|
|
|
2016-07-18 12:02:29 +02:00
|
|
|
if( ( v > 0 ) == injection )
|
|
|
|
|
sum += v;
|
2016-06-15 10:45:07 +02:00
|
|
|
}
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-07-18 12:02:29 +02:00
|
|
|
if( !injection ) sum *= -1;
|
2018-06-22 14:09:08 +02:00
|
|
|
|
|
|
|
|
if( polymer ) return { sum, measure::mass_rate };
|
2016-06-15 10:45:07 +02:00
|
|
|
return { sum, rate_unit< phase >() };
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
template< bool injection >
|
2016-11-02 12:18:39 +01:00
|
|
|
inline quantity flowing( const fn_args& args ) {
|
|
|
|
|
const auto& wells = args.wells;
|
2019-11-12 08:29:28 +01:00
|
|
|
auto pred = [&wells]( const Opm::Well& w ) {
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto& name = w.name();
|
|
|
|
|
return w.isInjector( ) == injection
|
2016-11-02 12:18:39 +01:00
|
|
|
&& wells.count( name ) > 0
|
|
|
|
|
&& wells.at( name ).flowing();
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
return { double( std::count_if( args.schedule_wells.begin(),
|
|
|
|
|
args.schedule_wells.end(),
|
|
|
|
|
pred ) ),
|
|
|
|
|
measure::identity };
|
|
|
|
|
}
|
|
|
|
|
|
2018-06-26 14:24:17 +02:00
|
|
|
template< rt phase, bool injection = true, bool polymer = false >
|
2016-07-15 15:20:56 +02:00
|
|
|
inline quantity crate( const fn_args& args ) {
|
2016-07-19 09:45:50 +02:00
|
|
|
const quantity zero = { 0, rate_unit< phase >() };
|
2016-10-25 16:25:46 +02:00
|
|
|
// The args.num value is the literal value which will go to the
|
|
|
|
|
// NUMS array in the eclispe SMSPEC file; the values in this array
|
|
|
|
|
// are offset 1 - whereas we need to use this index here to look
|
|
|
|
|
// up a completion with offset 0.
|
2018-02-08 10:40:13 +01:00
|
|
|
const size_t global_index = args.num - 1;
|
2016-07-19 09:45:50 +02:00
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
2016-07-15 15:20:56 +02:00
|
|
|
|
2018-06-26 14:24:17 +02:00
|
|
|
const auto& well = args.schedule_wells.front();
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto& name = well.name();
|
2016-09-25 18:31:34 +02:00
|
|
|
if( args.wells.count( name ) == 0 ) return zero;
|
2016-07-18 12:02:29 +02:00
|
|
|
|
2018-06-26 14:24:17 +02:00
|
|
|
const auto& well_data = args.wells.at( name );
|
|
|
|
|
const auto& completion = std::find_if( well_data.connections.begin(),
|
|
|
|
|
well_data.connections.end(),
|
2019-10-07 13:32:31 +02:00
|
|
|
[=]( const Opm::data::Connection& c ) {
|
2018-02-08 10:40:13 +01:00
|
|
|
return c.index == global_index;
|
2016-10-24 16:52:26 +02:00
|
|
|
} );
|
|
|
|
|
|
2018-06-26 14:24:17 +02:00
|
|
|
if( completion == well_data.connections.end() ) return zero;
|
2018-01-31 08:16:31 +01:00
|
|
|
|
|
|
|
|
double eff_fac = efac( args.eff_factors, name );
|
2018-06-26 14:24:17 +02:00
|
|
|
double concentration = polymer
|
2019-05-04 12:00:32 +02:00
|
|
|
? well.getPolymerProperties().m_polymerConcentration
|
2018-06-26 14:24:17 +02:00
|
|
|
: 1;
|
2018-01-31 08:16:31 +01:00
|
|
|
|
2018-06-26 14:24:17 +02:00
|
|
|
auto v = completion->rates.get( phase, 0.0 ) * eff_fac * concentration;
|
2016-07-19 09:45:50 +02:00
|
|
|
if( ( v > 0 ) != injection ) return zero;
|
2018-06-26 14:24:17 +02:00
|
|
|
if( !injection ) v *= -1;
|
2016-07-19 09:45:50 +02:00
|
|
|
|
2018-06-26 14:24:17 +02:00
|
|
|
if( polymer ) return { v, measure::mass_rate };
|
2016-07-19 09:45:50 +02:00
|
|
|
return { v, rate_unit< phase >() };
|
2016-07-15 15:20:56 +02:00
|
|
|
}
|
|
|
|
|
|
2018-10-05 12:59:42 +02:00
|
|
|
template< rt phase, bool polymer = false >
|
|
|
|
|
inline quantity srate( const fn_args& args ) {
|
|
|
|
|
const quantity zero = { 0, rate_unit< phase >() };
|
|
|
|
|
// The args.num value is the literal value which will go to the
|
|
|
|
|
// NUMS array in the eclispe SMSPEC file; the values in this array
|
|
|
|
|
// are offset 1 - whereas we need to use this index here to look
|
|
|
|
|
// up a completion with offset 0.
|
|
|
|
|
const size_t segNumber = args.num;
|
|
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
|
|
|
|
|
|
|
|
|
const auto& well = args.schedule_wells.front();
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto& name = well.name();
|
2018-10-05 12:59:42 +02:00
|
|
|
if( args.wells.count( name ) == 0 ) return zero;
|
|
|
|
|
|
2018-10-29 09:12:50 +01:00
|
|
|
const auto& well_data = args.wells.at( name );
|
|
|
|
|
|
2018-10-05 12:59:42 +02:00
|
|
|
const auto& segment = well_data.segments.find(segNumber);
|
|
|
|
|
|
|
|
|
|
if( segment == well_data.segments.end() ) return zero;
|
|
|
|
|
|
|
|
|
|
double eff_fac = efac( args.eff_factors, name );
|
|
|
|
|
double concentration = polymer
|
2019-05-04 12:00:32 +02:00
|
|
|
? well.getPolymerProperties().m_polymerConcentration
|
2018-10-05 12:59:42 +02:00
|
|
|
: 1;
|
|
|
|
|
|
|
|
|
|
auto v = segment->second.rates.get( phase, 0.0 ) * eff_fac * concentration;
|
|
|
|
|
//switch sign of rate - opposite convention in flow vs eclipse
|
|
|
|
|
v *= -1;
|
|
|
|
|
|
|
|
|
|
if( polymer ) return { v, measure::mass_rate };
|
|
|
|
|
return { v, rate_unit< phase >() };
|
|
|
|
|
}
|
|
|
|
|
|
2018-09-13 10:21:57 +02:00
|
|
|
inline quantity trans_factors ( const fn_args& args ) {
|
|
|
|
|
const quantity zero = { 0, measure::transmissibility };
|
|
|
|
|
|
|
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
|
|
|
|
// Like completion rate we need to look
|
|
|
|
|
// up a connection with offset 0.
|
|
|
|
|
const size_t global_index = args.num - 1;
|
|
|
|
|
|
|
|
|
|
const auto& well = args.schedule_wells.front();
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto& name = well.name();
|
2018-09-13 10:21:57 +02:00
|
|
|
if( args.wells.count( name ) == 0 ) return zero;
|
|
|
|
|
|
|
|
|
|
const auto& grid = args.grid;
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto& connections = well.getConnections();
|
2018-09-13 10:21:57 +02:00
|
|
|
|
|
|
|
|
const auto& connection = std::find_if(
|
|
|
|
|
connections.begin(),
|
|
|
|
|
connections.end(),
|
|
|
|
|
[=]( const Opm::Connection& c ) {
|
|
|
|
|
return grid.getGlobalIndex(c.getI(), c.getJ(), c.getK()) == global_index;
|
|
|
|
|
} );
|
|
|
|
|
|
|
|
|
|
if( connection == connections.end() ) return zero;
|
|
|
|
|
|
2018-10-29 09:12:50 +01:00
|
|
|
const auto& v = connection->CF() * connection->wellPi();
|
2018-09-13 10:21:57 +02:00
|
|
|
return { v, measure::transmissibility };
|
|
|
|
|
}
|
|
|
|
|
|
2018-10-05 12:59:42 +02:00
|
|
|
inline quantity spr ( const fn_args& args ) {
|
|
|
|
|
const quantity zero = { 0, measure::pressure };
|
|
|
|
|
|
|
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
|
|
|
|
// Like completion rate we need to look
|
|
|
|
|
// up a connection with offset 0.
|
|
|
|
|
const size_t segNumber = args.num;
|
|
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
|
|
|
|
|
|
|
|
|
const auto& well = args.schedule_wells.front();
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto& name = well.name();
|
2018-10-05 12:59:42 +02:00
|
|
|
if( args.wells.count( name ) == 0 ) return zero;
|
|
|
|
|
|
|
|
|
|
const auto& well_data = args.wells.at( name );
|
|
|
|
|
|
|
|
|
|
const auto& segment = well_data.segments.find(segNumber);
|
|
|
|
|
|
|
|
|
|
if( segment == well_data.segments.end() ) return zero;
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
const auto& v = segment->second.pressure;
|
|
|
|
|
return { v, measure::pressure };
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
inline quantity bhp( const fn_args& args ) {
|
2016-07-19 13:00:25 +02:00
|
|
|
const quantity zero = { 0, measure::pressure };
|
|
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto p = args.wells.find( args.schedule_wells.front().name() );
|
2016-09-25 18:31:34 +02:00
|
|
|
if( p == args.wells.end() ) return zero;
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
return { p->second.bhp, measure::pressure };
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
inline quantity thp( const fn_args& args ) {
|
2016-07-19 13:00:25 +02:00
|
|
|
const quantity zero = { 0, measure::pressure };
|
|
|
|
|
if( args.schedule_wells.empty() ) return zero;
|
2016-04-28 09:46:01 +02:00
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
const auto p = args.wells.find( args.schedule_wells.front().name() );
|
2016-09-25 18:31:34 +02:00
|
|
|
if( p == args.wells.end() ) return zero;
|
2016-04-28 09:46:01 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
return { p->second.thp, measure::pressure };
|
2016-04-28 09:46:01 +02:00
|
|
|
}
|
|
|
|
|
|
2017-11-27 14:33:13 +01:00
|
|
|
inline quantity bhp_history( const fn_args& args ) {
|
2018-03-07 08:57:58 +01:00
|
|
|
if( args.schedule_wells.empty() ) return { 0.0, measure::pressure };
|
2017-11-27 14:33:13 +01:00
|
|
|
|
2019-11-12 08:29:28 +01:00
|
|
|
const Opm::Well& sched_well = args.schedule_wells.front();
|
2017-11-27 14:33:13 +01:00
|
|
|
|
|
|
|
|
double bhp_hist;
|
2019-05-04 12:00:32 +02:00
|
|
|
if ( sched_well.isProducer( ) )
|
|
|
|
|
bhp_hist = sched_well.getProductionProperties().BHPH;
|
2017-11-27 14:33:13 +01:00
|
|
|
else
|
2019-05-04 12:00:32 +02:00
|
|
|
bhp_hist = sched_well.getInjectionProperties().BHPH;
|
2017-11-27 14:33:13 +01:00
|
|
|
|
|
|
|
|
return { bhp_hist, measure::pressure };
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
inline quantity thp_history( const fn_args& args ) {
|
2018-03-07 08:57:58 +01:00
|
|
|
if( args.schedule_wells.empty() ) return { 0.0, measure::pressure };
|
2017-11-27 14:33:13 +01:00
|
|
|
|
2019-11-12 08:29:28 +01:00
|
|
|
const Opm::Well& sched_well = args.schedule_wells.front();
|
2017-11-27 14:33:13 +01:00
|
|
|
|
|
|
|
|
double thp_hist;
|
2019-05-04 12:00:32 +02:00
|
|
|
if ( sched_well.isProducer() )
|
|
|
|
|
thp_hist = sched_well.getProductionProperties().THPH;
|
2017-11-27 14:33:13 +01:00
|
|
|
else
|
2019-05-04 12:00:32 +02:00
|
|
|
thp_hist = sched_well.getInjectionProperties().THPH;
|
2017-11-27 14:33:13 +01:00
|
|
|
|
|
|
|
|
return { thp_hist, measure::pressure };
|
|
|
|
|
}
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
template< Opm::Phase phase >
|
2016-06-15 10:45:07 +02:00
|
|
|
inline quantity production_history( const fn_args& args ) {
|
2016-06-10 09:54:44 +02:00
|
|
|
/*
|
|
|
|
|
* For well data, looking up historical rates (both for production and
|
|
|
|
|
* injection) before simulation actually starts is impossible and
|
|
|
|
|
* nonsensical. We therefore default to writing zero (which is what eclipse
|
2018-05-14 09:15:53 +02:00
|
|
|
* seems to do as well).
|
2016-06-10 09:54:44 +02:00
|
|
|
*/
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
double sum = 0.0;
|
2019-05-04 12:00:32 +02:00
|
|
|
for( const auto& sched_well : args.schedule_wells ){
|
2018-03-22 15:58:55 +01:00
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
double eff_fac = efac( args.eff_factors, sched_well.name() );
|
2019-06-21 21:18:46 +02:00
|
|
|
sum += sched_well.production_rate( args.st, phase ) * eff_fac;
|
2018-03-22 15:58:55 +01:00
|
|
|
}
|
|
|
|
|
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
return { sum, rate_unit< phase >() };
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
template< Opm::Phase phase >
|
2016-06-15 10:45:07 +02:00
|
|
|
inline quantity injection_history( const fn_args& args ) {
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
double sum = 0.0;
|
2019-05-04 12:00:32 +02:00
|
|
|
for( const auto& sched_well : args.schedule_wells ){
|
|
|
|
|
double eff_fac = efac( args.eff_factors, sched_well.name() );
|
2019-06-21 21:18:46 +02:00
|
|
|
sum += sched_well.injection_rate( args.st, phase ) * eff_fac;
|
2018-03-22 15:58:55 +01:00
|
|
|
}
|
|
|
|
|
|
2016-05-02 13:04:43 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
return { sum, rate_unit< phase >() };
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2017-12-20 10:37:24 +01:00
|
|
|
inline quantity res_vol_production_target( const fn_args& args ) {
|
|
|
|
|
|
|
|
|
|
double sum = 0.0;
|
2019-11-12 08:29:28 +01:00
|
|
|
for( const Opm::Well& sched_well : args.schedule_wells )
|
2019-05-04 12:00:32 +02:00
|
|
|
if (sched_well.getProductionProperties().predictionMode)
|
2020-01-03 12:58:55 +01:00
|
|
|
sum += sched_well.getProductionProperties().ResVRate.getSI();
|
2017-12-20 10:37:24 +01:00
|
|
|
|
|
|
|
|
return { sum, measure::rate };
|
|
|
|
|
}
|
|
|
|
|
|
2016-07-19 13:18:06 +02:00
|
|
|
inline quantity duration( const fn_args& args ) {
|
|
|
|
|
return { args.duration, measure::time };
|
|
|
|
|
}
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2016-10-25 16:26:14 +02:00
|
|
|
template<rt phase , bool injection>
|
|
|
|
|
quantity region_rate( const fn_args& args ) {
|
|
|
|
|
double sum = 0;
|
2018-06-10 20:54:34 +02:00
|
|
|
const auto& well_connections = args.regionCache.connections( args.num );
|
2018-01-31 08:16:31 +01:00
|
|
|
|
2018-06-10 20:54:34 +02:00
|
|
|
for (const auto& pair : well_connections) {
|
2018-01-31 08:16:31 +01:00
|
|
|
|
|
|
|
|
double eff_fac = efac( args.eff_factors, pair.first );
|
|
|
|
|
|
|
|
|
|
double rate = args.wells.get( pair.first , pair.second , phase ) * eff_fac;
|
2016-10-25 16:26:14 +02:00
|
|
|
|
|
|
|
|
// We are asking for the production rate in an injector - or
|
|
|
|
|
// opposite. We just clamp to zero.
|
|
|
|
|
if ((rate > 0) != injection)
|
|
|
|
|
rate = 0;
|
|
|
|
|
|
|
|
|
|
sum += rate;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
if( injection )
|
|
|
|
|
return { sum, rate_unit< phase >() };
|
|
|
|
|
else
|
|
|
|
|
return { -sum, rate_unit< phase >() };
|
|
|
|
|
}
|
|
|
|
|
|
2018-12-07 09:23:01 +01:00
|
|
|
template < rt phase, bool outputProducer = true, bool outputInjector = true>
|
2019-01-11 09:51:36 +01:00
|
|
|
inline quantity potential_rate( const fn_args& args ) {
|
2018-12-07 09:23:01 +01:00
|
|
|
double sum = 0.0;
|
|
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
for( const auto& sched_well : args.schedule_wells ) {
|
|
|
|
|
const auto& name = sched_well.name();
|
2018-12-07 09:23:01 +01:00
|
|
|
if( args.wells.count( name ) == 0 ) continue;
|
|
|
|
|
|
2019-05-04 12:00:32 +02:00
|
|
|
if (sched_well.isInjector() && outputInjector) {
|
2018-12-07 09:23:01 +01:00
|
|
|
const auto v = args.wells.at(name).rates.get(phase, 0.0);
|
|
|
|
|
sum += v;
|
|
|
|
|
}
|
2019-05-04 12:00:32 +02:00
|
|
|
else if (sched_well.isProducer() && outputProducer) {
|
2018-12-07 09:23:01 +01:00
|
|
|
const auto v = args.wells.at(name).rates.get(phase, 0.0);
|
|
|
|
|
sum += v;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
return { sum, rate_unit< phase >() };
|
2018-11-02 15:21:28 +01:00
|
|
|
}
|
2018-11-07 09:14:06 +01:00
|
|
|
|
2020-03-16 15:29:08 +01:00
|
|
|
template < bool isGroup, bool Producer, bool waterInjector, bool gasInjector>
|
2020-03-16 13:51:05 +01:00
|
|
|
inline quantity group_control( const fn_args& args ) {
|
|
|
|
|
|
|
|
|
|
std::string g_name = "";
|
|
|
|
|
if (isGroup) {
|
|
|
|
|
const quantity zero = { static_cast<double>(0), Opm::UnitSystem::measure::identity};
|
|
|
|
|
if( args.group_name.empty() ) return zero;
|
|
|
|
|
|
|
|
|
|
g_name = args.group_name;
|
|
|
|
|
}
|
|
|
|
|
else {
|
|
|
|
|
g_name = "FIELD";
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
int cntl_mode = 0;
|
|
|
|
|
|
|
|
|
|
// production control
|
|
|
|
|
if (Producer) {
|
|
|
|
|
const auto it_g = args.group.find(g_name);
|
|
|
|
|
if (it_g != args.group.end()) {
|
|
|
|
|
const auto& value = it_g->second.currentProdConstraint;
|
|
|
|
|
auto it_c = pCModeToPCntlMode.find(value);
|
|
|
|
|
if (it_c == pCModeToPCntlMode.end()) {
|
|
|
|
|
std::stringstream str;
|
|
|
|
|
str << "unknown control CMode: " << static_cast<int>(value);
|
|
|
|
|
throw std::invalid_argument(str.str());
|
|
|
|
|
}
|
|
|
|
|
cntl_mode = it_c->second;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
// water injection control
|
|
|
|
|
else if (waterInjector){
|
|
|
|
|
const auto it_g = args.group.find(g_name);
|
|
|
|
|
if (it_g != args.group.end()) {
|
|
|
|
|
const auto& value = it_g->second.currentWaterInjectionConstraint;
|
|
|
|
|
auto it_c = iCModeToICntlMode.find(value);
|
|
|
|
|
if (it_c == iCModeToICntlMode.end()) {
|
|
|
|
|
std::stringstream str;
|
|
|
|
|
str << "unknown control CMode: " << static_cast<int>(value);
|
|
|
|
|
throw std::invalid_argument(str.str());
|
|
|
|
|
}
|
|
|
|
|
cntl_mode = it_c->second;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
// gas injection control
|
|
|
|
|
else if (gasInjector){
|
|
|
|
|
const auto it_g = args.group.find(g_name);
|
|
|
|
|
if (it_g != args.group.end()) {
|
|
|
|
|
const auto& value = it_g->second.currentGasInjectionConstraint;
|
|
|
|
|
auto it_c = iCModeToICntlMode.find(value);
|
|
|
|
|
if (it_c == iCModeToICntlMode.end()) {
|
|
|
|
|
std::stringstream str;
|
|
|
|
|
str << "unknown control CMode: " << static_cast<int>(value);
|
|
|
|
|
throw std::invalid_argument(str.str());
|
|
|
|
|
}
|
|
|
|
|
cntl_mode = it_c->second;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
return {static_cast<double>(cntl_mode), Opm::UnitSystem::measure::identity};
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
/*
|
|
|
|
|
* A small DSL, really poor man's function composition, to avoid massive
|
|
|
|
|
* repetition when declaring the handlers for each individual keyword. bin_op
|
|
|
|
|
* and its aliases will apply the pair of functions F and G that all take const
|
|
|
|
|
* fn_args& and return quantity, making them composable.
|
|
|
|
|
*/
|
|
|
|
|
template< typename F, typename G, typename Op >
|
|
|
|
|
struct bin_op {
|
|
|
|
|
bin_op( F fn, G gn = {} ) : f( fn ), g( gn ) {}
|
|
|
|
|
quantity operator()( const fn_args& args ) const {
|
|
|
|
|
return Op()( f( args ), g( args ) );
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
|
|
|
|
F f;
|
|
|
|
|
G g;
|
|
|
|
|
};
|
|
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
template< typename F, typename G >
|
2016-07-19 13:18:06 +02:00
|
|
|
auto mul( F f, G g ) -> bin_op< F, G, std::multiplies< quantity > >
|
|
|
|
|
{ return { f, g }; }
|
2016-06-15 10:45:07 +02:00
|
|
|
|
|
|
|
|
template< typename F, typename G >
|
2016-07-19 13:18:06 +02:00
|
|
|
auto sum( F f, G g ) -> bin_op< F, G, std::plus< quantity > >
|
|
|
|
|
{ return { f, g }; }
|
2016-06-15 10:45:07 +02:00
|
|
|
|
|
|
|
|
template< typename F, typename G >
|
2016-07-19 13:18:06 +02:00
|
|
|
auto div( F f, G g ) -> bin_op< F, G, std::divides< quantity > >
|
|
|
|
|
{ return { f, g }; }
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2016-09-23 15:57:10 +02:00
|
|
|
template< typename F, typename G >
|
|
|
|
|
auto sub( F f, G g ) -> bin_op< F, G, std::minus< quantity > >
|
|
|
|
|
{ return { f, g }; }
|
|
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
using ofun = std::function< quantity( const fn_args& ) >;
|
|
|
|
|
|
|
|
|
|
static const std::unordered_map< std::string, ofun > funs = {
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "WWIR", rate< rt::wat, injector > },
|
|
|
|
|
{ "WOIR", rate< rt::oil, injector > },
|
|
|
|
|
{ "WGIR", rate< rt::gas, injector > },
|
|
|
|
|
{ "WNIR", rate< rt::solvent, injector > },
|
2018-06-22 14:09:08 +02:00
|
|
|
{ "WCIR", rate< rt::wat, injector, polymer > },
|
2019-11-06 14:58:35 +01:00
|
|
|
{ "WVIR", sum( sum( rate< rt::reservoir_water, injector >, rate< rt::reservoir_oil, injector > ),
|
|
|
|
|
rate< rt::reservoir_gas, injector > ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "WWIT", mul( rate< rt::wat, injector >, duration ) },
|
|
|
|
|
{ "WOIT", mul( rate< rt::oil, injector >, duration ) },
|
|
|
|
|
{ "WGIT", mul( rate< rt::gas, injector >, duration ) },
|
|
|
|
|
{ "WNIT", mul( rate< rt::solvent, injector >, duration ) },
|
2018-06-22 14:09:08 +02:00
|
|
|
{ "WCIT", mul( rate< rt::wat, injector, polymer >, duration ) },
|
2018-02-02 10:41:35 +01:00
|
|
|
{ "WVIT", mul( sum( sum( rate< rt::reservoir_water, injector >, rate< rt::reservoir_oil, injector > ),
|
|
|
|
|
rate< rt::reservoir_gas, injector > ), duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "WWPR", rate< rt::wat, producer > },
|
|
|
|
|
{ "WOPR", rate< rt::oil, producer > },
|
|
|
|
|
{ "WGPR", rate< rt::gas, producer > },
|
|
|
|
|
{ "WNPR", rate< rt::solvent, producer > },
|
|
|
|
|
|
2017-11-15 09:35:36 +01:00
|
|
|
{ "WGPRS", rate< rt::dissolved_gas, producer > },
|
|
|
|
|
{ "WGPRF", sub( rate< rt::gas, producer >, rate< rt::dissolved_gas, producer > ) },
|
2018-01-23 08:16:15 +01:00
|
|
|
{ "WOPRS", rate< rt::vaporized_oil, producer > },
|
|
|
|
|
{ "WOPRF", sub (rate < rt::oil, producer >, rate< rt::vaporized_oil, producer > ) },
|
2018-02-02 12:11:23 +01:00
|
|
|
{ "WVPR", sum( sum( rate< rt::reservoir_water, producer >, rate< rt::reservoir_oil, producer > ),
|
|
|
|
|
rate< rt::reservoir_gas, producer > ) },
|
2018-02-13 16:32:06 +01:00
|
|
|
{ "WGVPR", rate< rt::reservoir_gas, producer > },
|
2017-11-15 09:35:36 +01:00
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "WLPR", sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) },
|
|
|
|
|
{ "WWPT", mul( rate< rt::wat, producer >, duration ) },
|
|
|
|
|
{ "WOPT", mul( rate< rt::oil, producer >, duration ) },
|
|
|
|
|
{ "WGPT", mul( rate< rt::gas, producer >, duration ) },
|
|
|
|
|
{ "WNPT", mul( rate< rt::solvent, producer >, duration ) },
|
|
|
|
|
{ "WLPT", mul( sum( rate< rt::wat, producer >, rate< rt::oil, producer > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2016-07-18 12:02:29 +02:00
|
|
|
|
2017-11-15 15:01:13 +01:00
|
|
|
{ "WGPTS", mul( rate< rt::dissolved_gas, producer >, duration )},
|
|
|
|
|
{ "WGPTF", sub( mul( rate< rt::gas, producer >, duration ),
|
|
|
|
|
mul( rate< rt::dissolved_gas, producer >, duration ))},
|
2017-11-16 14:08:59 +01:00
|
|
|
{ "WOPTS", mul( rate< rt::vaporized_oil, producer >, duration )},
|
|
|
|
|
{ "WOPTF", sub( mul( rate< rt::oil, producer >, duration ),
|
|
|
|
|
mul( rate< rt::vaporized_oil, producer >, duration ))},
|
2018-02-02 12:11:23 +01:00
|
|
|
{ "WVPT", mul( sum( sum( rate< rt::reservoir_water, producer >, rate< rt::reservoir_oil, producer > ),
|
|
|
|
|
rate< rt::reservoir_gas, producer > ), duration ) },
|
2017-11-15 15:01:13 +01:00
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "WWCT", div( rate< rt::wat, producer >,
|
|
|
|
|
sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) ) },
|
|
|
|
|
{ "GWCT", div( rate< rt::wat, producer >,
|
|
|
|
|
sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) ) },
|
|
|
|
|
{ "WGOR", div( rate< rt::gas, producer >, rate< rt::oil, producer > ) },
|
|
|
|
|
{ "GGOR", div( rate< rt::gas, producer >, rate< rt::oil, producer > ) },
|
|
|
|
|
{ "WGLR", div( rate< rt::gas, producer >,
|
|
|
|
|
sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) ) },
|
2016-06-15 10:45:07 +02:00
|
|
|
|
|
|
|
|
{ "WBHP", bhp },
|
|
|
|
|
{ "WTHP", thp },
|
2018-02-02 12:11:23 +01:00
|
|
|
{ "WVPRT", res_vol_production_target },
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "GWIR", rate< rt::wat, injector > },
|
2018-02-13 16:32:06 +01:00
|
|
|
{ "WGVIR", rate< rt::reservoir_gas, injector >},
|
|
|
|
|
{ "WWVIR", rate< rt::reservoir_water, injector >},
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "GOIR", rate< rt::oil, injector > },
|
|
|
|
|
{ "GGIR", rate< rt::gas, injector > },
|
|
|
|
|
{ "GNIR", rate< rt::solvent, injector > },
|
2018-06-26 11:40:41 +02:00
|
|
|
{ "GCIR", rate< rt::wat, injector, polymer > },
|
2018-02-02 10:41:35 +01:00
|
|
|
{ "GVIR", sum( sum( rate< rt::reservoir_water, injector >, rate< rt::reservoir_oil, injector > ),
|
|
|
|
|
rate< rt::reservoir_gas, injector > ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "GWIT", mul( rate< rt::wat, injector >, duration ) },
|
|
|
|
|
{ "GOIT", mul( rate< rt::oil, injector >, duration ) },
|
|
|
|
|
{ "GGIT", mul( rate< rt::gas, injector >, duration ) },
|
|
|
|
|
{ "GNIT", mul( rate< rt::solvent, injector >, duration ) },
|
2018-06-26 11:40:41 +02:00
|
|
|
{ "GCIT", mul( rate< rt::wat, injector, polymer >, duration ) },
|
2018-02-02 10:41:35 +01:00
|
|
|
{ "GVIT", mul( sum( sum( rate< rt::reservoir_water, injector >, rate< rt::reservoir_oil, injector > ),
|
|
|
|
|
rate< rt::reservoir_gas, injector > ), duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "GWPR", rate< rt::wat, producer > },
|
|
|
|
|
{ "GOPR", rate< rt::oil, producer > },
|
|
|
|
|
{ "GGPR", rate< rt::gas, producer > },
|
|
|
|
|
{ "GNPR", rate< rt::solvent, producer > },
|
2017-11-17 14:56:57 +01:00
|
|
|
{ "GOPRS", rate< rt::vaporized_oil, producer > },
|
|
|
|
|
{ "GOPRF", sub (rate < rt::oil, producer >, rate< rt::vaporized_oil, producer > ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "GLPR", sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) },
|
2017-12-14 10:16:03 +01:00
|
|
|
{ "GVPR", sum( sum( rate< rt::reservoir_water, producer >, rate< rt::reservoir_oil, producer > ),
|
|
|
|
|
rate< rt::reservoir_gas, producer > ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "GWPT", mul( rate< rt::wat, producer >, duration ) },
|
|
|
|
|
{ "GOPT", mul( rate< rt::oil, producer >, duration ) },
|
|
|
|
|
{ "GGPT", mul( rate< rt::gas, producer >, duration ) },
|
|
|
|
|
{ "GNPT", mul( rate< rt::solvent, producer >, duration ) },
|
2017-11-17 14:56:57 +01:00
|
|
|
{ "GOPTS", mul( rate< rt::vaporized_oil, producer >, duration ) },
|
|
|
|
|
{ "GOPTF", mul( sub (rate < rt::oil, producer >,
|
|
|
|
|
rate< rt::vaporized_oil, producer > ),
|
|
|
|
|
duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "GLPT", mul( sum( rate< rt::wat, producer >, rate< rt::oil, producer > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2017-12-15 12:29:49 +01:00
|
|
|
{ "GVPT", mul( sum( sum( rate< rt::reservoir_water, producer >, rate< rt::reservoir_oil, producer > ),
|
|
|
|
|
rate< rt::reservoir_gas, producer > ), duration ) },
|
2018-12-05 14:29:32 +01:00
|
|
|
// Group potential
|
2019-01-11 09:51:36 +01:00
|
|
|
{ "GWPP", potential_rate< rt::well_potential_water , true, false>},
|
|
|
|
|
{ "GOPP", potential_rate< rt::well_potential_oil , true, false>},
|
|
|
|
|
{ "GGPP", potential_rate< rt::well_potential_gas , true, false>},
|
|
|
|
|
{ "GWPI", potential_rate< rt::well_potential_water , false, true>},
|
|
|
|
|
{ "GOPI", potential_rate< rt::well_potential_oil , false, true>},
|
|
|
|
|
{ "GGPI", potential_rate< rt::well_potential_gas , false, true>},
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2020-03-16 13:51:05 +01:00
|
|
|
//Group control mode
|
|
|
|
|
{ "GMCTP", group_control< true, true, false, false >},
|
|
|
|
|
{ "GMCTW", group_control< true, false, true, false >},
|
|
|
|
|
{ "GMCTG", group_control< true, false, false, true >},
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "WWPRH", production_history< Opm::Phase::WATER > },
|
|
|
|
|
{ "WOPRH", production_history< Opm::Phase::OIL > },
|
|
|
|
|
{ "WGPRH", production_history< Opm::Phase::GAS > },
|
|
|
|
|
{ "WLPRH", sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) },
|
|
|
|
|
|
|
|
|
|
{ "WWPTH", mul( production_history< Opm::Phase::WATER >, duration ) },
|
|
|
|
|
{ "WOPTH", mul( production_history< Opm::Phase::OIL >, duration ) },
|
|
|
|
|
{ "WGPTH", mul( production_history< Opm::Phase::GAS >, duration ) },
|
|
|
|
|
{ "WLPTH", mul( sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "WWIRH", injection_history< Opm::Phase::WATER > },
|
|
|
|
|
{ "WOIRH", injection_history< Opm::Phase::OIL > },
|
|
|
|
|
{ "WGIRH", injection_history< Opm::Phase::GAS > },
|
|
|
|
|
{ "WWITH", mul( injection_history< Opm::Phase::WATER >, duration ) },
|
|
|
|
|
{ "WOITH", mul( injection_history< Opm::Phase::OIL >, duration ) },
|
|
|
|
|
{ "WGITH", mul( injection_history< Opm::Phase::GAS >, duration ) },
|
2016-06-15 10:45:07 +02:00
|
|
|
|
|
|
|
|
/* From our point of view, injectors don't have water cuts and div/sum will return 0.0 */
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "WWCTH", div( production_history< Opm::Phase::WATER >,
|
|
|
|
|
sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) ) },
|
2016-04-26 12:26:38 +02:00
|
|
|
|
2016-06-15 10:45:07 +02:00
|
|
|
/* We do not support mixed injections, and gas/oil is undefined when oil is
|
|
|
|
|
* zero (i.e. pure gas injector), so always output 0 if this is an injector
|
|
|
|
|
*/
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "WGORH", div( production_history< Opm::Phase::GAS >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) },
|
|
|
|
|
{ "WGLRH", div( production_history< Opm::Phase::GAS >,
|
|
|
|
|
sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) ) },
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2017-11-27 14:33:13 +01:00
|
|
|
{ "WTHPH", thp_history },
|
|
|
|
|
{ "WBHPH", bhp_history },
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "GWPRH", production_history< Opm::Phase::WATER > },
|
|
|
|
|
{ "GOPRH", production_history< Opm::Phase::OIL > },
|
|
|
|
|
{ "GGPRH", production_history< Opm::Phase::GAS > },
|
|
|
|
|
{ "GLPRH", sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) },
|
|
|
|
|
{ "GWIRH", injection_history< Opm::Phase::WATER > },
|
|
|
|
|
{ "GOIRH", injection_history< Opm::Phase::OIL > },
|
|
|
|
|
{ "GGIRH", injection_history< Opm::Phase::GAS > },
|
|
|
|
|
{ "GGORH", div( production_history< Opm::Phase::GAS >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) },
|
|
|
|
|
{ "GWCTH", div( production_history< Opm::Phase::WATER >,
|
|
|
|
|
sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) ) },
|
|
|
|
|
|
|
|
|
|
{ "GWPTH", mul( production_history< Opm::Phase::WATER >, duration ) },
|
|
|
|
|
{ "GOPTH", mul( production_history< Opm::Phase::OIL >, duration ) },
|
|
|
|
|
{ "GGPTH", mul( production_history< Opm::Phase::GAS >, duration ) },
|
2017-11-14 12:52:10 +01:00
|
|
|
{ "GGPRF", sub( rate < rt::gas, producer >, rate< rt::dissolved_gas, producer > )},
|
2017-11-15 14:09:55 +01:00
|
|
|
{ "GGPRS", rate< rt::dissolved_gas, producer> },
|
|
|
|
|
{ "GGPTF", mul( sub( rate < rt::gas, producer >, rate< rt::dissolved_gas, producer > ),
|
|
|
|
|
duration ) },
|
|
|
|
|
{ "GGPTS", mul( rate< rt::dissolved_gas, producer>, duration ) },
|
2016-11-08 14:21:27 +01:00
|
|
|
{ "GGLR", div( rate< rt::gas, producer >,
|
|
|
|
|
sum( rate< rt::wat, producer >,
|
|
|
|
|
rate< rt::oil, producer > ) ) },
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "GGLRH", div( production_history< Opm::Phase::GAS >,
|
|
|
|
|
sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) ) },
|
|
|
|
|
{ "GLPTH", mul( sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "GWITH", mul( injection_history< Opm::Phase::WATER >, duration ) },
|
|
|
|
|
{ "GGITH", mul( injection_history< Opm::Phase::GAS >, duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "GMWIN", flowing< injector > },
|
|
|
|
|
{ "GMWPR", flowing< producer > },
|
|
|
|
|
|
2017-12-20 10:37:24 +01:00
|
|
|
{ "GVPRT", res_vol_production_target },
|
|
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "CWIR", crate< rt::wat, injector > },
|
|
|
|
|
{ "CGIR", crate< rt::gas, injector > },
|
2018-06-26 14:24:17 +02:00
|
|
|
{ "CCIR", crate< rt::wat, injector, polymer > },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "CWIT", mul( crate< rt::wat, injector >, duration ) },
|
|
|
|
|
{ "CGIT", mul( crate< rt::gas, injector >, duration ) },
|
|
|
|
|
{ "CNIT", mul( crate< rt::solvent, injector >, duration ) },
|
|
|
|
|
|
|
|
|
|
{ "CWPR", crate< rt::wat, producer > },
|
|
|
|
|
{ "COPR", crate< rt::oil, producer > },
|
|
|
|
|
{ "CGPR", crate< rt::gas, producer > },
|
2020-03-17 22:24:02 +01:00
|
|
|
{ "CGFR", sub(crate<rt::gas, producer>, crate<rt::gas, injector>) },
|
|
|
|
|
{ "COFR", sub(crate<rt::oil, producer>, crate<rt::oil, injector>) },
|
|
|
|
|
{ "CWFR", sub(crate<rt::wat, producer>, crate<rt::wat, injector>) },
|
|
|
|
|
{ "CWCT", div( crate< rt::wat, producer >,
|
|
|
|
|
sum( crate< rt::wat, producer >, crate< rt::oil, producer > ) ) },
|
|
|
|
|
{ "CGOR", div( crate< rt::gas, producer >, crate< rt::oil, producer > ) },
|
2016-09-23 15:57:10 +02:00
|
|
|
// Minus for injection rates and pluss for production rate
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "CNFR", sub( crate< rt::solvent, producer >, crate<rt::solvent, injector >) },
|
|
|
|
|
{ "CWPT", mul( crate< rt::wat, producer >, duration ) },
|
|
|
|
|
{ "COPT", mul( crate< rt::oil, producer >, duration ) },
|
|
|
|
|
{ "CGPT", mul( crate< rt::gas, producer >, duration ) },
|
|
|
|
|
{ "CNPT", mul( crate< rt::solvent, producer >, duration ) },
|
2018-06-26 14:24:17 +02:00
|
|
|
{ "CCIT", mul( crate< rt::wat, injector, polymer >, duration ) },
|
2019-11-11 09:26:36 +01:00
|
|
|
{ "CCPT", mul( crate< rt::wat, producer, polymer >, duration ) },
|
2018-09-13 10:21:57 +02:00
|
|
|
{ "CTFAC", trans_factors },
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "FWPR", rate< rt::wat, producer > },
|
|
|
|
|
{ "FOPR", rate< rt::oil, producer > },
|
|
|
|
|
{ "FGPR", rate< rt::gas, producer > },
|
|
|
|
|
{ "FNPR", rate< rt::solvent, producer > },
|
2017-12-18 13:26:54 +01:00
|
|
|
{ "FVPR", sum( sum( rate< rt::reservoir_water, producer>, rate< rt::reservoir_oil, producer >),
|
|
|
|
|
rate< rt::reservoir_gas, producer>)},
|
2017-12-19 08:48:12 +01:00
|
|
|
{ "FGPRS", rate< rt::dissolved_gas, producer > },
|
|
|
|
|
{ "FGPRF", sub( rate< rt::gas, producer >, rate< rt::dissolved_gas, producer > ) },
|
2018-01-23 14:36:56 +01:00
|
|
|
{ "FOPRS", rate< rt::vaporized_oil, producer > },
|
|
|
|
|
{ "FOPRF", sub (rate < rt::oil, producer >, rate< rt::vaporized_oil, producer > ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "FLPR", sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) },
|
|
|
|
|
{ "FWPT", mul( rate< rt::wat, producer >, duration ) },
|
|
|
|
|
{ "FOPT", mul( rate< rt::oil, producer >, duration ) },
|
|
|
|
|
{ "FGPT", mul( rate< rt::gas, producer >, duration ) },
|
2018-05-08 11:25:47 +02:00
|
|
|
{ "FNPT", mul( rate< rt::solvent, producer >, duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "FLPT", mul( sum( rate< rt::wat, producer >, rate< rt::oil, producer > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2017-12-18 13:26:54 +01:00
|
|
|
{ "FVPT", mul(sum (sum( rate< rt::reservoir_water, producer>, rate< rt::reservoir_oil, producer >),
|
|
|
|
|
rate< rt::reservoir_gas, producer>), duration)},
|
2017-12-19 08:54:55 +01:00
|
|
|
{ "FGPTS", mul( rate< rt::dissolved_gas, producer > , duration )},
|
|
|
|
|
{ "FGPTF", mul( sub( rate< rt::gas, producer >, rate< rt::dissolved_gas, producer > ), duration )},
|
2018-01-23 14:36:56 +01:00
|
|
|
{ "FOPTS", mul( rate< rt::vaporized_oil, producer >, duration ) },
|
|
|
|
|
{ "FOPTF", mul( sub (rate < rt::oil, producer >,
|
|
|
|
|
rate< rt::vaporized_oil, producer > ),
|
|
|
|
|
duration ) },
|
2016-07-18 12:57:06 +02:00
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "FWIR", rate< rt::wat, injector > },
|
|
|
|
|
{ "FOIR", rate< rt::oil, injector > },
|
|
|
|
|
{ "FGIR", rate< rt::gas, injector > },
|
|
|
|
|
{ "FNIR", rate< rt::solvent, injector > },
|
2018-06-26 11:41:10 +02:00
|
|
|
{ "FCIR", rate< rt::wat, injector, polymer > },
|
2019-11-11 09:26:36 +01:00
|
|
|
{ "FCPR", rate< rt::wat, producer, polymer > },
|
2017-12-18 13:26:54 +01:00
|
|
|
{ "FVIR", sum( sum( rate< rt::reservoir_water, injector>, rate< rt::reservoir_oil, injector >),
|
|
|
|
|
rate< rt::reservoir_gas, injector>)},
|
2016-11-02 12:35:11 +01:00
|
|
|
|
|
|
|
|
{ "FLIR", sum( rate< rt::wat, injector >, rate< rt::oil, injector > ) },
|
|
|
|
|
{ "FWIT", mul( rate< rt::wat, injector >, duration ) },
|
|
|
|
|
{ "FOIT", mul( rate< rt::oil, injector >, duration ) },
|
|
|
|
|
{ "FGIT", mul( rate< rt::gas, injector >, duration ) },
|
|
|
|
|
{ "FNIT", mul( rate< rt::solvent, injector >, duration ) },
|
2018-06-26 11:41:10 +02:00
|
|
|
{ "FCIT", mul( rate< rt::wat, injector, polymer >, duration ) },
|
2019-11-11 09:26:36 +01:00
|
|
|
{ "FCPT", mul( rate< rt::wat, producer, polymer >, duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "FLIT", mul( sum( rate< rt::wat, injector >, rate< rt::oil, injector > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2017-12-18 13:26:54 +01:00
|
|
|
{ "FVIT", mul( sum( sum( rate< rt::reservoir_water, injector>, rate< rt::reservoir_oil, injector >),
|
|
|
|
|
rate< rt::reservoir_gas, injector>), duration)},
|
2018-12-05 14:29:32 +01:00
|
|
|
// Field potential
|
2019-01-11 09:51:36 +01:00
|
|
|
{ "FWPP", potential_rate< rt::well_potential_water , true, false>},
|
|
|
|
|
{ "FOPP", potential_rate< rt::well_potential_oil , true, false>},
|
|
|
|
|
{ "FGPP", potential_rate< rt::well_potential_gas , true, false>},
|
|
|
|
|
{ "FWPI", potential_rate< rt::well_potential_water , false, true>},
|
|
|
|
|
{ "FOPI", potential_rate< rt::well_potential_oil , false, true>},
|
|
|
|
|
{ "FGPI", potential_rate< rt::well_potential_gas , false, true>},
|
2018-12-05 14:29:32 +01:00
|
|
|
|
2016-07-18 12:57:06 +02:00
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "FWPRH", production_history< Opm::Phase::WATER > },
|
|
|
|
|
{ "FOPRH", production_history< Opm::Phase::OIL > },
|
|
|
|
|
{ "FGPRH", production_history< Opm::Phase::GAS > },
|
|
|
|
|
{ "FLPRH", sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) },
|
|
|
|
|
{ "FWPTH", mul( production_history< Opm::Phase::WATER >, duration ) },
|
|
|
|
|
{ "FOPTH", mul( production_history< Opm::Phase::OIL >, duration ) },
|
|
|
|
|
{ "FGPTH", mul( production_history< Opm::Phase::GAS >, duration ) },
|
|
|
|
|
{ "FLPTH", mul( sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ),
|
2016-07-19 13:18:06 +02:00
|
|
|
duration ) },
|
2016-07-18 12:57:06 +02:00
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "FWIRH", injection_history< Opm::Phase::WATER > },
|
|
|
|
|
{ "FOIRH", injection_history< Opm::Phase::OIL > },
|
|
|
|
|
{ "FGIRH", injection_history< Opm::Phase::GAS > },
|
|
|
|
|
{ "FWITH", mul( injection_history< Opm::Phase::WATER >, duration ) },
|
|
|
|
|
{ "FOITH", mul( injection_history< Opm::Phase::OIL >, duration ) },
|
|
|
|
|
{ "FGITH", mul( injection_history< Opm::Phase::GAS >, duration ) },
|
2016-07-18 12:57:06 +02:00
|
|
|
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "FWCT", div( rate< rt::wat, producer >,
|
|
|
|
|
sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) ) },
|
|
|
|
|
{ "FGOR", div( rate< rt::gas, producer >, rate< rt::oil, producer > ) },
|
|
|
|
|
{ "FGLR", div( rate< rt::gas, producer >,
|
|
|
|
|
sum( rate< rt::wat, producer >, rate< rt::oil, producer > ) ) },
|
2019-10-07 13:32:31 +02:00
|
|
|
{ "FWCTH", div( production_history< Opm::Phase::WATER >,
|
|
|
|
|
sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) ) },
|
|
|
|
|
{ "FGORH", div( production_history< Opm::Phase::GAS >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) },
|
|
|
|
|
{ "FGLRH", div( production_history< Opm::Phase::GAS >,
|
|
|
|
|
sum( production_history< Opm::Phase::WATER >,
|
|
|
|
|
production_history< Opm::Phase::OIL > ) ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "FMWIN", flowing< injector > },
|
|
|
|
|
{ "FMWPR", flowing< producer > },
|
2018-02-02 12:11:23 +01:00
|
|
|
{ "FVPRT", res_vol_production_target },
|
2016-08-23 14:09:44 +02:00
|
|
|
|
2020-03-16 13:51:05 +01:00
|
|
|
//Field control mode
|
|
|
|
|
{ "FMCTP", group_control< false, true, false, false >},
|
|
|
|
|
{ "FMCTW", group_control< false, false, true, false >},
|
|
|
|
|
{ "FMCTG", group_control< false, false, false, true >},
|
|
|
|
|
|
2016-08-23 14:09:44 +02:00
|
|
|
/* Region properties */
|
2016-11-08 14:37:25 +01:00
|
|
|
{ "ROIR" , region_rate< rt::oil, injector > },
|
|
|
|
|
{ "RGIR" , region_rate< rt::gas, injector > },
|
|
|
|
|
{ "RWIR" , region_rate< rt::wat, injector > },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "ROPR" , region_rate< rt::oil, producer > },
|
|
|
|
|
{ "RGPR" , region_rate< rt::gas, producer > },
|
|
|
|
|
{ "RWPR" , region_rate< rt::wat, producer > },
|
2016-11-08 14:37:25 +01:00
|
|
|
{ "ROIT" , mul( region_rate< rt::oil, injector >, duration ) },
|
|
|
|
|
{ "RGIT" , mul( region_rate< rt::gas, injector >, duration ) },
|
|
|
|
|
{ "RWIT" , mul( region_rate< rt::wat, injector >, duration ) },
|
2016-11-02 12:35:11 +01:00
|
|
|
{ "ROPT" , mul( region_rate< rt::oil, producer >, duration ) },
|
|
|
|
|
{ "RGPT" , mul( region_rate< rt::gas, producer >, duration ) },
|
|
|
|
|
{ "RWPT" , mul( region_rate< rt::wat, producer >, duration ) },
|
2018-10-05 12:59:42 +02:00
|
|
|
//Multisegment well segment data
|
|
|
|
|
{ "SOFR", srate< rt::oil > },
|
|
|
|
|
{ "SWFR", srate< rt::wat > },
|
|
|
|
|
{ "SGFR", srate< rt::gas > },
|
|
|
|
|
{ "SPR", spr },
|
2018-11-07 09:14:06 +01:00
|
|
|
// Well productivity index
|
2019-01-11 09:51:36 +01:00
|
|
|
{ "WPIW", potential_rate< rt::productivity_index_water >},
|
|
|
|
|
{ "WPIO", potential_rate< rt::productivity_index_oil >},
|
|
|
|
|
{ "WPIG", potential_rate< rt::productivity_index_gas >},
|
|
|
|
|
{ "WPIL", sum( potential_rate< rt::productivity_index_water >, potential_rate< rt::productivity_index_oil>)},
|
2018-11-07 09:14:06 +01:00
|
|
|
// Well potential
|
2019-01-11 09:51:36 +01:00
|
|
|
{ "WWPP", potential_rate< rt::well_potential_water , true, false>},
|
|
|
|
|
{ "WOPP", potential_rate< rt::well_potential_oil , true, false>},
|
|
|
|
|
{ "WGPP", potential_rate< rt::well_potential_gas , true, false>},
|
|
|
|
|
{ "WWPI", potential_rate< rt::well_potential_water , false, true>},
|
|
|
|
|
{ "WOPI", potential_rate< rt::well_potential_oil , false, true>},
|
|
|
|
|
{ "WGPI", potential_rate< rt::well_potential_gas , false, true>},
|
2016-06-15 10:45:07 +02:00
|
|
|
};
|
2016-04-26 12:26:38 +02:00
|
|
|
|
2016-10-24 12:44:56 +02:00
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
static const std::unordered_map< std::string, Opm::UnitSystem::measure> single_values_units = {
|
|
|
|
|
{"TCPU" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"ELAPSED" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"NEWTON" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"NLINERS" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"NLINSMIN" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"NLINSMAX" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"MLINEARS" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"MSUMLINS" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"MSUMNEWT" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"TCPUTS" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"TIMESTEP" , Opm::UnitSystem::measure::time },
|
|
|
|
|
{"TCPUDAY" , Opm::UnitSystem::measure::time },
|
|
|
|
|
{"STEPTYPE" , Opm::UnitSystem::measure::identity },
|
|
|
|
|
{"TELAPLIN" , Opm::UnitSystem::measure::time },
|
|
|
|
|
{"FWIP" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"FOIP" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"FGIP" , Opm::UnitSystem::measure::gas_surface_volume },
|
|
|
|
|
{"FOIPL" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"FOIPG" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"FGIPL" , Opm::UnitSystem::measure::gas_surface_volume },
|
|
|
|
|
{"FGIPG" , Opm::UnitSystem::measure::gas_surface_volume },
|
|
|
|
|
{"FPR" , Opm::UnitSystem::measure::pressure },
|
2018-01-17 15:23:23 +01:00
|
|
|
|
|
|
|
|
};
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
static const std::unordered_map< std::string, Opm::UnitSystem::measure> region_units = {
|
|
|
|
|
{"RPR" , Opm::UnitSystem::measure::pressure},
|
|
|
|
|
{"ROIP" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"ROIPL" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"ROIPG" , Opm::UnitSystem::measure::liquid_surface_volume },
|
|
|
|
|
{"RGIP" , Opm::UnitSystem::measure::gas_surface_volume },
|
|
|
|
|
{"RGIPL" , Opm::UnitSystem::measure::gas_surface_volume },
|
|
|
|
|
{"RGIPG" , Opm::UnitSystem::measure::gas_surface_volume },
|
|
|
|
|
{"RWIP" , Opm::UnitSystem::measure::liquid_surface_volume }
|
2017-06-06 13:02:46 +02:00
|
|
|
};
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
static const std::unordered_map< std::string, Opm::UnitSystem::measure> block_units = {
|
|
|
|
|
{"BPR" , Opm::UnitSystem::measure::pressure},
|
|
|
|
|
{"BPRESSUR" , Opm::UnitSystem::measure::pressure},
|
|
|
|
|
{"BSWAT" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BWSAT" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BSGAS" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BGSAT" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BOSAT" , Opm::UnitSystem::measure::identity},
|
2020-03-19 14:18:00 +01:00
|
|
|
{"BWKR" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BOKR" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BKRO" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BGKR" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BKRG" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BKRW" , Opm::UnitSystem::measure::identity},
|
|
|
|
|
{"BWPC" , Opm::UnitSystem::measure::pressure},
|
|
|
|
|
{"BGPC" , Opm::UnitSystem::measure::pressure},
|
2019-10-07 13:32:31 +02:00
|
|
|
{"BVWAT" , Opm::UnitSystem::measure::viscosity},
|
|
|
|
|
{"BWVIS" , Opm::UnitSystem::measure::viscosity},
|
|
|
|
|
{"BVGAS" , Opm::UnitSystem::measure::viscosity},
|
|
|
|
|
{"BGVIS" , Opm::UnitSystem::measure::viscosity},
|
|
|
|
|
{"BVOIL" , Opm::UnitSystem::measure::viscosity},
|
|
|
|
|
{"BOVIS" , Opm::UnitSystem::measure::viscosity},
|
2018-01-19 15:27:28 +01:00
|
|
|
};
|
|
|
|
|
|
2019-11-12 08:29:28 +01:00
|
|
|
inline std::vector<Opm::Well> find_wells( const Opm::Schedule& schedule,
|
2020-03-19 14:17:16 +01:00
|
|
|
const Opm::EclIO::SummaryNode& node,
|
|
|
|
|
const int sim_step,
|
|
|
|
|
const Opm::out::RegionCache& regionCache ) {
|
|
|
|
|
switch (node.category) {
|
|
|
|
|
case Opm::SummaryConfigNode::Category::Well: [[fallthrough]];
|
|
|
|
|
case Opm::SummaryConfigNode::Category::Connection: [[fallthrough]];
|
|
|
|
|
case Opm::SummaryConfigNode::Category::Segment: {
|
|
|
|
|
const auto& name = node.wgname;
|
|
|
|
|
|
|
|
|
|
if (schedule.hasWell(node.wgname, sim_step)) {
|
|
|
|
|
return { schedule.getWell( name, sim_step ) };
|
|
|
|
|
} else {
|
2019-05-04 12:00:32 +02:00
|
|
|
return {};
|
2020-03-19 14:17:16 +01:00
|
|
|
}
|
2016-04-26 12:26:38 +02:00
|
|
|
}
|
2016-04-27 09:57:39 +02:00
|
|
|
|
2020-03-19 14:17:16 +01:00
|
|
|
case Opm::SummaryConfigNode::Category::Group: {
|
|
|
|
|
const auto& name = node.wgname;
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2016-10-05 14:52:25 +02:00
|
|
|
if( !schedule.hasGroup( name ) ) return {};
|
|
|
|
|
|
2019-09-01 07:57:48 +02:00
|
|
|
return schedule.getChildWells2( name, sim_step);
|
2016-04-26 12:26:38 +02:00
|
|
|
}
|
2016-04-27 09:57:39 +02:00
|
|
|
|
2020-03-19 14:17:16 +01:00
|
|
|
case Opm::SummaryConfigNode::Category::Field:
|
2019-11-12 08:29:28 +01:00
|
|
|
return schedule.getWells(sim_step);
|
2016-07-18 12:57:06 +02:00
|
|
|
|
2020-03-19 14:17:16 +01:00
|
|
|
case Opm::SummaryConfigNode::Category::Region: {
|
2019-11-12 08:29:28 +01:00
|
|
|
std::vector<Opm::Well> wells;
|
2018-02-02 10:10:35 +01:00
|
|
|
|
2020-03-19 14:17:16 +01:00
|
|
|
const auto region = node.number;
|
2018-02-02 10:10:35 +01:00
|
|
|
|
2018-06-10 20:54:34 +02:00
|
|
|
for ( const auto& connection : regionCache.connections( region ) ){
|
|
|
|
|
const auto& w_name = connection.first;
|
2019-05-04 12:00:32 +02:00
|
|
|
if (schedule.hasWell(w_name, sim_step)) {
|
2019-11-12 08:29:28 +01:00
|
|
|
const auto& well = schedule.getWell( w_name, sim_step );
|
2019-05-04 12:00:32 +02:00
|
|
|
|
|
|
|
|
const auto& it = std::find_if( wells.begin(), wells.end(),
|
2019-11-12 08:29:28 +01:00
|
|
|
[&] ( const Opm::Well& elem )
|
2019-05-04 12:00:32 +02:00
|
|
|
{ return elem.name() == well.name(); });
|
|
|
|
|
if ( it == wells.end() )
|
2019-10-09 00:25:45 -05:00
|
|
|
wells.push_back( well );
|
2019-05-04 12:00:32 +02:00
|
|
|
}
|
2018-02-02 10:10:35 +01:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
return wells;
|
2020-03-19 14:17:16 +01:00
|
|
|
|
2018-02-02 10:10:35 +01:00
|
|
|
}
|
|
|
|
|
|
2020-03-19 14:17:16 +01:00
|
|
|
default:
|
|
|
|
|
return {};
|
|
|
|
|
}
|
2016-04-26 12:26:38 +02:00
|
|
|
}
|
|
|
|
|
|
2020-03-19 13:28:17 +01:00
|
|
|
bool need_wells(const Opm::EclIO::SummaryNode& node) {
|
|
|
|
|
static const std::regex region_keyword_regex { "R[OGW][IP][RT]" };
|
|
|
|
|
|
|
|
|
|
switch (node.category) {
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Connection:
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Field:
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Group:
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Segment:
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Well:
|
|
|
|
|
return true;
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Region:
|
|
|
|
|
return std::regex_match(node.keyword, region_keyword_regex);
|
|
|
|
|
default:
|
|
|
|
|
return false;
|
2019-05-02 10:14:15 +02:00
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
2019-10-07 13:32:31 +02:00
|
|
|
void eval_udq(const Opm::Schedule& schedule, std::size_t sim_step, Opm::SummaryState& st)
|
2019-01-31 19:33:45 +01:00
|
|
|
{
|
2019-10-07 13:32:31 +02:00
|
|
|
using namespace Opm;
|
|
|
|
|
|
2019-07-26 17:17:06 +02:00
|
|
|
const UDQConfig& udq = schedule.getUDQConfig(sim_step);
|
2019-03-07 10:32:35 +01:00
|
|
|
const auto& func_table = udq.function_table();
|
2019-09-18 12:58:55 +02:00
|
|
|
UDQContext context(func_table, st);
|
2019-12-17 16:45:53 +01:00
|
|
|
{
|
|
|
|
|
const std::vector<std::string> wells = st.wells();
|
|
|
|
|
|
|
|
|
|
for (const auto& assign : udq.assignments(UDQVarType::WELL_VAR)) {
|
|
|
|
|
auto ws = assign.eval(wells);
|
|
|
|
|
for (const auto& well : wells) {
|
|
|
|
|
const auto& udq_value = ws[well];
|
|
|
|
|
if (udq_value)
|
|
|
|
|
st.update_well_var(well, ws.name(), udq_value.value());
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
for (const auto& def : udq.definitions(UDQVarType::WELL_VAR)) {
|
|
|
|
|
auto ws = def.eval(context);
|
|
|
|
|
for (const auto& well : wells) {
|
|
|
|
|
const auto& udq_value = ws[well];
|
|
|
|
|
if (udq_value)
|
|
|
|
|
st.update_well_var(well, def.keyword(), udq_value.value());
|
|
|
|
|
}
|
2019-02-21 10:44:47 +01:00
|
|
|
}
|
|
|
|
|
}
|
2019-01-31 19:33:45 +01:00
|
|
|
|
2019-12-17 16:45:53 +01:00
|
|
|
{
|
|
|
|
|
const std::vector<std::string> groups = st.groups();
|
|
|
|
|
|
|
|
|
|
for (const auto& assign : udq.assignments(UDQVarType::GROUP_VAR)) {
|
|
|
|
|
auto ws = assign.eval(groups);
|
|
|
|
|
for (const auto& group : groups) {
|
|
|
|
|
const auto& udq_value = ws[group];
|
|
|
|
|
if (udq_value)
|
|
|
|
|
st.update_group_var(group, ws.name(), udq_value.value());
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
for (const auto& def : udq.definitions(UDQVarType::GROUP_VAR)) {
|
|
|
|
|
auto ws = def.eval(context);
|
|
|
|
|
for (const auto& group : groups) {
|
|
|
|
|
const auto& udq_value = ws[group];
|
|
|
|
|
if (udq_value)
|
|
|
|
|
st.update_group_var(group, def.keyword(), udq_value.value());
|
|
|
|
|
}
|
2019-03-07 10:32:35 +01:00
|
|
|
}
|
|
|
|
|
}
|
2019-05-24 15:50:40 +02:00
|
|
|
|
|
|
|
|
for (const auto& def : udq.definitions(UDQVarType::FIELD_VAR)) {
|
|
|
|
|
auto field_udq = def.eval(context);
|
|
|
|
|
if (field_udq[0])
|
|
|
|
|
st.update(def.keyword(), field_udq[0].value());
|
|
|
|
|
}
|
2019-03-07 10:32:35 +01:00
|
|
|
}
|
2019-10-07 13:44:51 +02:00
|
|
|
|
2020-03-19 13:34:56 +01:00
|
|
|
void updateValue(const Opm::EclIO::SummaryNode& node, const double value, Opm::SummaryState& st) {
|
|
|
|
|
switch (node.category) {
|
|
|
|
|
case Opm::EclIO::SummaryNode::Category::Well:
|
|
|
|
|
return st.update_well_var(node.wgname, node.keyword, value);
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2020-03-19 13:34:56 +01:00
|
|
|
case Opm::SummaryConfigNode::Category::Group:
|
|
|
|
|
return st.update_group_var(node.wgname, node.keyword, value);
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2020-03-19 13:34:56 +01:00
|
|
|
default:
|
|
|
|
|
return st.update(node.unique_key(), value);
|
|
|
|
|
}
|
2019-10-09 00:25:45 -05:00
|
|
|
}
|
|
|
|
|
|
2019-10-07 13:44:51 +02:00
|
|
|
/*
|
|
|
|
|
* The well efficiency factor will not impact the well rate itself, but is
|
|
|
|
|
* rather applied for accumulated values.The WEFAC can be considered to shut
|
|
|
|
|
* and open the well for short intervals within the same timestep, and the well
|
|
|
|
|
* is therefore solved at full speed.
|
|
|
|
|
*
|
|
|
|
|
* Groups are treated similarly as wells. The group's GEFAC is not applied for
|
|
|
|
|
* rates, only for accumulated volumes. When GEFAC is set for a group, it is
|
|
|
|
|
* considered that all wells are taken down simultaneously, and GEFAC is
|
|
|
|
|
* therefore not applied to the group's rate. However, any efficiency factors
|
|
|
|
|
* applied to the group's wells or sub-groups must be included.
|
|
|
|
|
*
|
|
|
|
|
* Regions and fields will have the well and group efficiency applied for both
|
|
|
|
|
* rates and accumulated values.
|
|
|
|
|
*
|
|
|
|
|
*/
|
2019-10-09 00:25:45 -05:00
|
|
|
struct EfficiencyFactor
|
|
|
|
|
{
|
|
|
|
|
using Factor = std::pair<std::string, double>;
|
|
|
|
|
using FacColl = std::vector<Factor>;
|
|
|
|
|
|
|
|
|
|
FacColl factors{};
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
void setFactors(const Opm::SummaryConfigNode& node,
|
2019-10-09 00:25:45 -05:00
|
|
|
const Opm::Schedule& schedule,
|
2019-11-12 08:29:28 +01:00
|
|
|
const std::vector<Opm::Well>& schedule_wells,
|
2019-10-09 00:25:45 -05:00
|
|
|
const int sim_step);
|
|
|
|
|
};
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
void EfficiencyFactor::setFactors(const Opm::SummaryConfigNode& node,
|
2019-10-09 00:25:45 -05:00
|
|
|
const Opm::Schedule& schedule,
|
2019-11-12 08:29:28 +01:00
|
|
|
const std::vector<Opm::Well>& schedule_wells,
|
2019-10-09 00:25:45 -05:00
|
|
|
const int sim_step)
|
|
|
|
|
{
|
|
|
|
|
this->factors.clear();
|
|
|
|
|
|
|
|
|
|
if (schedule_wells.empty()) { return; }
|
2019-10-07 13:44:51 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
const auto cat = node.category();
|
2020-03-13 14:24:29 +01:00
|
|
|
if( cat != Opm::SummaryConfigNode::Category::Group
|
|
|
|
|
&& cat != Opm::SummaryConfigNode::Category::Field
|
|
|
|
|
&& cat != Opm::SummaryConfigNode::Category::Region
|
|
|
|
|
&& (node.type() != Opm::SummaryConfigNode::Type::Total))
|
2019-10-09 00:25:45 -05:00
|
|
|
return;
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
const bool is_group = (cat == Opm::SummaryConfigNode::Category::Group);
|
|
|
|
|
const bool is_rate = (node.type() != Opm::SummaryConfigNode::Type::Total);
|
2019-10-07 13:44:51 +02:00
|
|
|
|
|
|
|
|
for( const auto& well : schedule_wells ) {
|
|
|
|
|
if (!well.hasBeenDefined(sim_step))
|
|
|
|
|
continue;
|
|
|
|
|
|
|
|
|
|
double eff_factor = well.getEfficiencyFactor();
|
2019-11-12 08:29:28 +01:00
|
|
|
const auto* group_ptr = std::addressof(schedule.getGroup(well.groupName(), sim_step));
|
2019-10-07 13:44:51 +02:00
|
|
|
|
|
|
|
|
while(true){
|
|
|
|
|
if(( is_group
|
|
|
|
|
&& is_rate
|
2019-10-09 00:25:45 -05:00
|
|
|
&& group_ptr->name() == node.namedEntity() ))
|
2019-10-07 13:44:51 +02:00
|
|
|
break;
|
|
|
|
|
eff_factor *= group_ptr->getGroupEfficiencyFactor();
|
|
|
|
|
|
|
|
|
|
if (group_ptr->name() == "FIELD")
|
|
|
|
|
break;
|
2019-11-12 08:29:28 +01:00
|
|
|
group_ptr = std::addressof( schedule.getGroup( group_ptr->parent(), sim_step ) );
|
2019-10-07 13:44:51 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
this->factors.emplace_back( well.name(), eff_factor );
|
|
|
|
|
}
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
namespace Evaluator {
|
|
|
|
|
struct InputData
|
|
|
|
|
{
|
|
|
|
|
const Opm::EclipseState& es;
|
|
|
|
|
const Opm::Schedule& sched;
|
|
|
|
|
const Opm::EclipseGrid& grid;
|
|
|
|
|
const Opm::out::RegionCache& reg;
|
|
|
|
|
};
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
struct SimulatorResults
|
|
|
|
|
{
|
|
|
|
|
const Opm::data::WellRates& wellSol;
|
2020-03-16 13:51:05 +01:00
|
|
|
const Opm::data::Group& groupSol;
|
2019-10-09 00:25:45 -05:00
|
|
|
const std::map<std::string, double>& single;
|
|
|
|
|
const std::map<std::string, std::vector<double>>& region;
|
|
|
|
|
const std::map<std::pair<std::string, int>, double>& block;
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class Base
|
|
|
|
|
{
|
2016-06-15 10:45:07 +02:00
|
|
|
public:
|
2019-10-09 00:25:45 -05:00
|
|
|
virtual ~Base() {}
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
virtual void update(const std::size_t sim_step,
|
|
|
|
|
const double stepSize,
|
|
|
|
|
const InputData& input,
|
|
|
|
|
const SimulatorResults& simRes,
|
|
|
|
|
Opm::SummaryState& st) const = 0;
|
|
|
|
|
};
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
class FunctionRelation : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
2020-03-13 14:24:29 +01:00
|
|
|
explicit FunctionRelation(Opm::SummaryConfigNode node, ofun fcn)
|
2019-10-09 00:25:45 -05:00
|
|
|
: node_(std::move(node))
|
|
|
|
|
, fcn_ (std::move(fcn))
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void update(const std::size_t sim_step,
|
|
|
|
|
const double stepSize,
|
|
|
|
|
const InputData& input,
|
|
|
|
|
const SimulatorResults& simRes,
|
|
|
|
|
Opm::SummaryState& st) const override
|
|
|
|
|
{
|
|
|
|
|
const auto get_wells =
|
2020-03-19 13:28:17 +01:00
|
|
|
need_wells(node_);
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
const auto wells = get_wells
|
|
|
|
|
? find_wells(input.sched, this->node_,
|
|
|
|
|
static_cast<int>(sim_step), input.reg)
|
2019-11-12 08:29:28 +01:00
|
|
|
: std::vector<Opm::Well>{};
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
if (get_wells && wells.empty())
|
|
|
|
|
// Parameter depends on well information, but no active
|
|
|
|
|
// wells apply at this sim_step. Nothing to do.
|
|
|
|
|
return;
|
|
|
|
|
|
2020-03-16 13:51:05 +01:00
|
|
|
std::string group_name = this->node_.category() == Opm::SummaryConfigNode::Category::Group ? this->node_.namedEntity() : "";
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
EfficiencyFactor efac{};
|
|
|
|
|
efac.setFactors(this->node_, input.sched, wells, sim_step);
|
|
|
|
|
|
|
|
|
|
const fn_args args {
|
2020-03-16 13:51:05 +01:00
|
|
|
wells, group_name, stepSize, static_cast<int>(sim_step),
|
2019-10-09 00:25:45 -05:00
|
|
|
std::max(0, this->node_.number()),
|
2020-03-16 13:51:05 +01:00
|
|
|
st, simRes.wellSol, simRes.groupSol, input.reg, input.grid,
|
2019-10-09 00:25:45 -05:00
|
|
|
std::move(efac.factors)
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
const auto& usys = input.es.getUnits();
|
|
|
|
|
const auto prm = this->fcn_(args);
|
|
|
|
|
|
|
|
|
|
updateValue(this->node_, usys.from_si(prm.unit, prm.value), st);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
2020-03-13 14:24:29 +01:00
|
|
|
Opm::SummaryConfigNode node_;
|
2019-10-09 00:25:45 -05:00
|
|
|
ofun fcn_;
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class BlockValue : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
2020-03-13 14:24:29 +01:00
|
|
|
explicit BlockValue(Opm::SummaryConfigNode node,
|
2019-10-09 00:25:45 -05:00
|
|
|
const Opm::UnitSystem::measure m)
|
|
|
|
|
: node_(std::move(node))
|
|
|
|
|
, m_ (m)
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void update(const std::size_t /* sim_step */,
|
|
|
|
|
const double /* stepSize */,
|
|
|
|
|
const InputData& input,
|
|
|
|
|
const SimulatorResults& simRes,
|
|
|
|
|
Opm::SummaryState& st) const override
|
|
|
|
|
{
|
|
|
|
|
auto xPos = simRes.block.find(this->lookupKey());
|
|
|
|
|
if (xPos == simRes.block.end()) {
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
const auto& usys = input.es.getUnits();
|
|
|
|
|
updateValue(this->node_, usys.from_si(this->m_, xPos->second), st);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
2020-03-13 14:24:29 +01:00
|
|
|
Opm::SummaryConfigNode node_;
|
2019-10-09 00:25:45 -05:00
|
|
|
Opm::UnitSystem::measure m_;
|
|
|
|
|
|
|
|
|
|
Opm::out::Summary::BlockValues::key_type lookupKey() const
|
|
|
|
|
{
|
|
|
|
|
return { this->node_.keyword(), this->node_.number() };
|
|
|
|
|
}
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class RegionValue : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
2020-03-13 14:24:29 +01:00
|
|
|
explicit RegionValue(Opm::SummaryConfigNode node,
|
2019-10-09 00:25:45 -05:00
|
|
|
const Opm::UnitSystem::measure m)
|
|
|
|
|
: node_(std::move(node))
|
|
|
|
|
, m_ (m)
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void update(const std::size_t /* sim_step */,
|
|
|
|
|
const double /* stepSize */,
|
|
|
|
|
const InputData& input,
|
|
|
|
|
const SimulatorResults& simRes,
|
|
|
|
|
Opm::SummaryState& st) const override
|
|
|
|
|
{
|
|
|
|
|
if (this->node_.number() < 0)
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
auto xPos = simRes.region.find(this->node_.keyword());
|
|
|
|
|
if (xPos == simRes.region.end())
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
const auto ix = this->index();
|
|
|
|
|
if (ix >= xPos->second.size())
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
const auto val = xPos->second[ix];
|
|
|
|
|
const auto& usys = input.es.getUnits();
|
|
|
|
|
|
|
|
|
|
updateValue(this->node_, usys.from_si(this->m_, val), st);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
2020-03-13 14:24:29 +01:00
|
|
|
Opm::SummaryConfigNode node_;
|
2019-10-09 00:25:45 -05:00
|
|
|
Opm::UnitSystem::measure m_;
|
|
|
|
|
|
|
|
|
|
std::vector<double>::size_type index() const
|
|
|
|
|
{
|
|
|
|
|
return this->node_.number() - 1;
|
|
|
|
|
}
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class GlobalProcessValue : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
2020-03-13 14:24:29 +01:00
|
|
|
explicit GlobalProcessValue(Opm::SummaryConfigNode node,
|
2019-10-09 00:25:45 -05:00
|
|
|
const Opm::UnitSystem::measure m)
|
|
|
|
|
: node_(std::move(node))
|
|
|
|
|
, m_ (m)
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void update(const std::size_t /* sim_step */,
|
|
|
|
|
const double /* stepSize */,
|
|
|
|
|
const InputData& input,
|
|
|
|
|
const SimulatorResults& simRes,
|
|
|
|
|
Opm::SummaryState& st) const override
|
|
|
|
|
{
|
|
|
|
|
auto xPos = simRes.single.find(this->node_.keyword());
|
|
|
|
|
if (xPos == simRes.single.end())
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
const auto val = xPos->second;
|
|
|
|
|
const auto& usys = input.es.getUnits();
|
|
|
|
|
|
|
|
|
|
updateValue(this->node_, usys.from_si(this->m_, val), st);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
2020-03-13 14:24:29 +01:00
|
|
|
Opm::SummaryConfigNode node_;
|
2019-10-09 00:25:45 -05:00
|
|
|
Opm::UnitSystem::measure m_;
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class UserDefinedValue : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
void update(const std::size_t /* sim_step */,
|
|
|
|
|
const double /* stepSize */,
|
|
|
|
|
const InputData& /* input */,
|
|
|
|
|
const SimulatorResults& /* simRes */,
|
|
|
|
|
Opm::SummaryState& /* st */) const override
|
|
|
|
|
{
|
|
|
|
|
// No-op
|
|
|
|
|
}
|
|
|
|
|
};
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
class Time : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
explicit Time(std::string saveKey)
|
|
|
|
|
: saveKey_(std::move(saveKey))
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void update(const std::size_t /* sim_step */,
|
|
|
|
|
const double stepSize,
|
|
|
|
|
const InputData& input,
|
|
|
|
|
const SimulatorResults& /* simRes */,
|
|
|
|
|
Opm::SummaryState& st) const override
|
|
|
|
|
{
|
|
|
|
|
const auto& usys = input.es.getUnits();
|
|
|
|
|
|
|
|
|
|
const auto m = ::Opm::UnitSystem::measure::time;
|
|
|
|
|
const auto val = st.get_elapsed() + stepSize;
|
|
|
|
|
|
|
|
|
|
st.update(this->saveKey_, usys.from_si(m, val));
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
|
|
|
|
std::string saveKey_;
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class Years : public Base
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
explicit Years(std::string saveKey)
|
|
|
|
|
: saveKey_(std::move(saveKey))
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void update(const std::size_t /* sim_step */,
|
|
|
|
|
const double stepSize,
|
|
|
|
|
const InputData& /* input */,
|
|
|
|
|
const SimulatorResults& /* simRes */,
|
|
|
|
|
Opm::SummaryState& st) const override
|
|
|
|
|
{
|
|
|
|
|
using namespace ::Opm::unit;
|
|
|
|
|
|
|
|
|
|
const auto val = st.get_elapsed() + stepSize;
|
|
|
|
|
|
|
|
|
|
st.update(this->saveKey_, convert::to(val, year));
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
private:
|
|
|
|
|
std::string saveKey_;
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
class Factory
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
struct Descriptor
|
|
|
|
|
{
|
|
|
|
|
std::string uniquekey{};
|
|
|
|
|
std::string unit{};
|
|
|
|
|
std::unique_ptr<Base> evaluator{};
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
explicit Factory(const Opm::EclipseState& es,
|
|
|
|
|
const Opm::EclipseGrid& grid,
|
|
|
|
|
const Opm::SummaryState& st,
|
|
|
|
|
const Opm::UDQConfig& udq)
|
|
|
|
|
: es_(es), grid_(grid), st_(st), udq_(udq)
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
~Factory() = default;
|
|
|
|
|
|
|
|
|
|
Factory(const Factory&) = delete;
|
|
|
|
|
Factory(Factory&&) = delete;
|
|
|
|
|
Factory& operator=(const Factory&) = delete;
|
|
|
|
|
Factory& operator=(Factory&&) = delete;
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
Descriptor create(const Opm::SummaryConfigNode&);
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
private:
|
|
|
|
|
const Opm::EclipseState& es_;
|
|
|
|
|
const Opm::EclipseGrid& grid_;
|
|
|
|
|
const Opm::SummaryState& st_;
|
|
|
|
|
const Opm::UDQConfig& udq_;
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
const Opm::SummaryConfigNode* node_;
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
Opm::UnitSystem::measure paramUnit_;
|
|
|
|
|
ofun paramFunction_;
|
|
|
|
|
|
|
|
|
|
Descriptor functionRelation();
|
|
|
|
|
Descriptor blockValue();
|
|
|
|
|
Descriptor regionValue();
|
|
|
|
|
Descriptor globalProcessValue();
|
|
|
|
|
Descriptor userDefinedValue();
|
|
|
|
|
Descriptor unknownParameter();
|
|
|
|
|
|
|
|
|
|
bool isBlockValue();
|
|
|
|
|
bool isRegionValue();
|
|
|
|
|
bool isGlobalProcessValue();
|
|
|
|
|
bool isFunctionRelation();
|
|
|
|
|
bool isUserDefined();
|
|
|
|
|
|
|
|
|
|
std::string functionUnitString() const;
|
|
|
|
|
std::string directUnitString() const;
|
|
|
|
|
std::string userDefinedUnit() const;
|
|
|
|
|
};
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
Factory::Descriptor Factory::create(const Opm::SummaryConfigNode& node)
|
2019-10-09 00:25:45 -05:00
|
|
|
{
|
|
|
|
|
this->node_ = &node;
|
|
|
|
|
|
|
|
|
|
if (this->isUserDefined())
|
|
|
|
|
return this->userDefinedValue();
|
|
|
|
|
|
|
|
|
|
if (this->isBlockValue())
|
|
|
|
|
return this->blockValue();
|
|
|
|
|
|
|
|
|
|
if (this->isRegionValue())
|
|
|
|
|
return this->regionValue();
|
|
|
|
|
|
|
|
|
|
if (this->isGlobalProcessValue())
|
|
|
|
|
return this->globalProcessValue();
|
|
|
|
|
|
|
|
|
|
if (this->isFunctionRelation())
|
|
|
|
|
return this->functionRelation();
|
|
|
|
|
|
|
|
|
|
return this->unknownParameter();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Factory::Descriptor Factory::functionRelation()
|
|
|
|
|
{
|
|
|
|
|
auto desc = this->unknownParameter();
|
|
|
|
|
|
|
|
|
|
desc.unit = this->functionUnitString();
|
|
|
|
|
desc.evaluator.reset(new FunctionRelation {
|
|
|
|
|
*this->node_, std::move(this->paramFunction_)
|
|
|
|
|
});
|
|
|
|
|
|
|
|
|
|
return desc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Factory::Descriptor Factory::blockValue()
|
|
|
|
|
{
|
|
|
|
|
auto desc = this->unknownParameter();
|
|
|
|
|
|
|
|
|
|
desc.unit = this->directUnitString();
|
|
|
|
|
desc.evaluator.reset(new BlockValue {
|
|
|
|
|
*this->node_, this->paramUnit_
|
|
|
|
|
});
|
|
|
|
|
|
|
|
|
|
return desc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Factory::Descriptor Factory::regionValue()
|
|
|
|
|
{
|
|
|
|
|
auto desc = this->unknownParameter();
|
|
|
|
|
|
|
|
|
|
desc.unit = this->directUnitString();
|
|
|
|
|
desc.evaluator.reset(new RegionValue {
|
|
|
|
|
*this->node_, this->paramUnit_
|
|
|
|
|
});
|
|
|
|
|
|
|
|
|
|
return desc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Factory::Descriptor Factory::globalProcessValue()
|
|
|
|
|
{
|
|
|
|
|
auto desc = this->unknownParameter();
|
|
|
|
|
|
|
|
|
|
desc.unit = this->directUnitString();
|
|
|
|
|
desc.evaluator.reset(new GlobalProcessValue {
|
|
|
|
|
*this->node_, this->paramUnit_
|
|
|
|
|
});
|
|
|
|
|
|
|
|
|
|
return desc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Factory::Descriptor Factory::userDefinedValue()
|
|
|
|
|
{
|
|
|
|
|
auto desc = this->unknownParameter();
|
|
|
|
|
|
|
|
|
|
desc.unit = this->userDefinedUnit();
|
|
|
|
|
desc.evaluator.reset(new UserDefinedValue {});
|
|
|
|
|
|
|
|
|
|
return desc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Factory::Descriptor Factory::unknownParameter()
|
|
|
|
|
{
|
|
|
|
|
auto desc = Descriptor{};
|
|
|
|
|
|
|
|
|
|
desc.uniquekey = this->node_->uniqueNodeKey();
|
|
|
|
|
|
|
|
|
|
return desc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
bool Factory::isBlockValue()
|
|
|
|
|
{
|
|
|
|
|
auto pos = block_units.find(this->node_->keyword());
|
|
|
|
|
if (pos == block_units.end())
|
|
|
|
|
return false;
|
|
|
|
|
|
|
|
|
|
if (! this->grid_.cellActive(this->node_->number() - 1))
|
|
|
|
|
// 'node_' is a block value, but it is configured in a
|
|
|
|
|
// deactivated cell. Don't create an evaluation function.
|
|
|
|
|
return false;
|
|
|
|
|
|
|
|
|
|
// 'node_' represents a block value in an active cell.
|
|
|
|
|
// Capture unit of measure and return true.
|
|
|
|
|
this->paramUnit_ = pos->second;
|
|
|
|
|
return true;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
bool Factory::isRegionValue()
|
|
|
|
|
{
|
|
|
|
|
auto pos = region_units.find(this->node_->keyword());
|
|
|
|
|
if (pos == region_units.end())
|
|
|
|
|
return false;
|
|
|
|
|
|
|
|
|
|
// 'node_' represents a region value. Capture unit
|
|
|
|
|
// of measure and return true.
|
|
|
|
|
this->paramUnit_ = pos->second;
|
|
|
|
|
return true;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
bool Factory::isGlobalProcessValue()
|
|
|
|
|
{
|
|
|
|
|
auto pos = single_values_units.find(this->node_->keyword());
|
|
|
|
|
if (pos == single_values_units.end())
|
|
|
|
|
return false;
|
|
|
|
|
|
|
|
|
|
// 'node_' represents a single value (i.e., global process)
|
|
|
|
|
// value. Capture unit of measure and return true.
|
|
|
|
|
this->paramUnit_ = pos->second;
|
|
|
|
|
return true;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
bool Factory::isFunctionRelation()
|
|
|
|
|
{
|
|
|
|
|
auto pos = funs.find(this->node_->keyword());
|
|
|
|
|
if (pos == funs.end())
|
|
|
|
|
return false;
|
|
|
|
|
|
|
|
|
|
// 'node_' represents a functional relation.
|
|
|
|
|
// Capture evaluation function and return true.
|
|
|
|
|
this->paramFunction_ = pos->second;
|
|
|
|
|
return true;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
bool Factory::isUserDefined()
|
|
|
|
|
{
|
|
|
|
|
return this->node_->isUserDefined();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
std::string Factory::functionUnitString() const
|
|
|
|
|
{
|
|
|
|
|
const auto reg = Opm::out::RegionCache{};
|
|
|
|
|
|
|
|
|
|
const fn_args args {
|
2020-03-16 13:51:05 +01:00
|
|
|
{}, "", 0.0, 0, std::max(0, this->node_->number()),
|
|
|
|
|
this->st_, {}, {}, reg, this->grid_,
|
2019-10-09 00:25:45 -05:00
|
|
|
{}
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
const auto prm = this->paramFunction_(args);
|
|
|
|
|
|
|
|
|
|
return this->es_.getUnits().name(prm.unit);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
std::string Factory::directUnitString() const
|
|
|
|
|
{
|
|
|
|
|
return this->es_.getUnits().name(this->paramUnit_);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
std::string Factory::userDefinedUnit() const
|
|
|
|
|
{
|
|
|
|
|
const auto& kw = this->node_->keyword();
|
|
|
|
|
|
|
|
|
|
return this->udq_.has_unit(kw)
|
|
|
|
|
? this->udq_.unit(kw) : "?????";
|
|
|
|
|
}
|
|
|
|
|
} // namespace Evaluator
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
void reportUnsupportedKeywords(std::vector<Opm::SummaryConfigNode> keywords)
|
2019-10-09 00:25:45 -05:00
|
|
|
{
|
2020-01-24 11:13:05 +01:00
|
|
|
// Sort by location first, then keyword.
|
2020-03-13 14:24:29 +01:00
|
|
|
auto loc_kw_ordering = [](const Opm::SummaryConfigNode& n1, const Opm::SummaryConfigNode& n2) {
|
2020-01-24 11:13:05 +01:00
|
|
|
if (n1.location() == n2.location()) {
|
|
|
|
|
return n1.keyword() < n2.keyword();
|
|
|
|
|
}
|
|
|
|
|
if (n1.location().filename == n2.location().filename) {
|
|
|
|
|
return n1.location().lineno < n2.location().lineno;
|
|
|
|
|
}
|
|
|
|
|
return n1.location().filename < n2.location().filename;
|
|
|
|
|
};
|
|
|
|
|
std::sort(keywords.begin(), keywords.end(), loc_kw_ordering);
|
|
|
|
|
|
|
|
|
|
// Reorder to remove duplicate { keyword, location } pairs, since
|
|
|
|
|
// that will give duplicate and therefore useless warnings.
|
2020-03-13 14:24:29 +01:00
|
|
|
auto same_kw_and_loc = [](const Opm::SummaryConfigNode& n1, const Opm::SummaryConfigNode& n2) {
|
2020-01-24 11:13:05 +01:00
|
|
|
return (n1.keyword() == n2.keyword()) && (n1.location() == n2.location());
|
|
|
|
|
};
|
|
|
|
|
auto uend = std::unique(keywords.begin(), keywords.end(), same_kw_and_loc);
|
2019-10-09 00:25:45 -05:00
|
|
|
|
2019-11-07 12:03:38 +01:00
|
|
|
for (auto node = keywords.begin(); node != uend; ++node) {
|
|
|
|
|
const auto& location = node->location();
|
|
|
|
|
::Opm::OpmLog::warning("Unhandled summary keyword '" + node->keyword() + "' at " + location.filename + ", line " + std::to_string(location.lineno));
|
|
|
|
|
}
|
2019-10-09 00:25:45 -05:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
std::string makeWGName(std::string name)
|
2016-04-20 09:55:17 +02:00
|
|
|
{
|
2019-10-09 00:25:45 -05:00
|
|
|
// Use default WGNAME if 'name' is empty or consists
|
|
|
|
|
// exlusively of white-space (space and tab) characters.
|
|
|
|
|
//
|
|
|
|
|
// Use 'name' itself otherwise.
|
2017-06-22 22:07:59 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
const auto use_dflt = name.empty() ||
|
|
|
|
|
(name.find_first_not_of(" \t") == std::string::npos);
|
Decouple SummaryConfig From LibECL
This commit reimplements the SummaryNode class in order not to use
the ecl::smspec_node class from LibECL. Consequently, we remove
class SummaryConfig's binding to LibECL.
Class SummaryNode maintains the same information as before. We also
implement a "named constructor" strategy to assign data members that
only make sense for a subset of the node categories. The previous
member function 'type' is renamed to 'category' to identify the
attachment category (e.g., Well, Group, Field, Block, Region). In
turn, we introduce a new 'type' member function to identify the
parameter kind (e.g, pressure, rate, cumulative total, well count)
represented by a given node. We furthermore capture whether or not
the node is a user defined quantity (i.e., a UDQ).
We reimplement the keyword classifier operations that are currently
needed as free functions named 'is_*()' in SummaryConfig.cpp.
Note that in addition to the renamed member functions of class
SummaryNode, this commit also switches the summary key strategy for
block parameters. Rather than capturing the 'ijk' values
individually as "BOSAT:3,3,6", we now store the equivalent global
(Cartesian) index (i.e., as "BOSAT:523"). Code that directly
constructs block parameter keys must be updated accordingly.
Chase the API change in the 'Summary' constructor and update unit
tests as needed.
2019-09-27 16:50:42 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
return use_dflt ? std::string(":+:+:+:+") : std::move(name);
|
|
|
|
|
}
|
2017-06-22 22:07:59 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
class SummaryOutputParameters
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
using EvalPtr = std::unique_ptr<Evaluator::Base>;
|
|
|
|
|
using SMSpecPrm = Opm::EclIO::OutputStream::
|
|
|
|
|
SummarySpecification::Parameters;
|
|
|
|
|
|
|
|
|
|
SummaryOutputParameters() = default;
|
|
|
|
|
~SummaryOutputParameters() = default;
|
|
|
|
|
|
|
|
|
|
SummaryOutputParameters(const SummaryOutputParameters& rhs) = delete;
|
|
|
|
|
SummaryOutputParameters(SummaryOutputParameters&& rhs) = default;
|
|
|
|
|
|
|
|
|
|
SummaryOutputParameters&
|
|
|
|
|
operator=(const SummaryOutputParameters& rhs) = delete;
|
|
|
|
|
|
|
|
|
|
SummaryOutputParameters&
|
|
|
|
|
operator=(SummaryOutputParameters&& rhs) = default;
|
|
|
|
|
|
|
|
|
|
void makeParameter(std::string keyword,
|
|
|
|
|
std::string name,
|
|
|
|
|
const int num,
|
|
|
|
|
std::string unit,
|
|
|
|
|
EvalPtr evaluator)
|
|
|
|
|
{
|
|
|
|
|
this->smspec_.add(std::move(keyword), std::move(name),
|
|
|
|
|
std::max (num, 0), std::move(unit));
|
|
|
|
|
|
|
|
|
|
this->evaluators_.push_back(std::move(evaluator));
|
2017-06-22 22:07:59 +02:00
|
|
|
}
|
2018-01-19 15:27:28 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
const SMSpecPrm& summarySpecification() const
|
|
|
|
|
{
|
|
|
|
|
return this->smspec_;
|
|
|
|
|
}
|
2018-01-19 15:27:28 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
const std::vector<EvalPtr>& getEvaluators() const
|
|
|
|
|
{
|
|
|
|
|
return this->evaluators_;
|
|
|
|
|
}
|
2017-06-06 13:02:46 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
private:
|
|
|
|
|
SMSpecPrm smspec_{};
|
|
|
|
|
std::vector<EvalPtr> evaluators_{};
|
|
|
|
|
};
|
2017-06-06 13:02:46 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
class SMSpecStreamDeferredCreation
|
|
|
|
|
{
|
|
|
|
|
private:
|
|
|
|
|
using Spec = ::Opm::EclIO::OutputStream::SummarySpecification;
|
|
|
|
|
|
|
|
|
|
public:
|
|
|
|
|
using ResultSet = ::Opm::EclIO::OutputStream::ResultSet;
|
|
|
|
|
using Formatted = ::Opm::EclIO::OutputStream::Formatted;
|
2017-06-06 13:02:46 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
explicit SMSpecStreamDeferredCreation(const Opm::InitConfig& initcfg,
|
|
|
|
|
const Opm::EclipseGrid& grid,
|
|
|
|
|
const std::time_t start,
|
|
|
|
|
const Opm::UnitSystem::UnitType utype);
|
2017-06-06 13:02:46 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
std::unique_ptr<Spec>
|
|
|
|
|
createStream(const ResultSet& rset, const Formatted& fmt) const
|
|
|
|
|
{
|
|
|
|
|
return std::make_unique<Spec>(rset, fmt, this->uconv(),
|
|
|
|
|
this->cartDims_, this->restart_,
|
|
|
|
|
this->start_);
|
|
|
|
|
}
|
2017-06-06 13:02:46 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
private:
|
|
|
|
|
Opm::UnitSystem::UnitType utype_;
|
|
|
|
|
std::array<int,3> cartDims_;
|
|
|
|
|
Spec::StartTime start_;
|
|
|
|
|
Spec::RestartSpecification restart_{};
|
2019-03-09 10:19:27 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
Spec::UnitConvention uconv() const;
|
|
|
|
|
};
|
2019-03-09 10:19:27 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
SMSpecStreamDeferredCreation::
|
|
|
|
|
SMSpecStreamDeferredCreation(const Opm::InitConfig& initcfg,
|
|
|
|
|
const Opm::EclipseGrid& grid,
|
|
|
|
|
const std::time_t start,
|
|
|
|
|
const Opm::UnitSystem::UnitType utype)
|
|
|
|
|
: utype_ (utype)
|
|
|
|
|
, cartDims_(grid.getNXYZ())
|
|
|
|
|
, start_ (std::chrono::system_clock::from_time_t(start))
|
|
|
|
|
{
|
|
|
|
|
if (initcfg.restartRequested()) {
|
|
|
|
|
this->restart_.root = initcfg.getRestartRootName();
|
|
|
|
|
this->restart_.step = initcfg.getRestartStep();
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
2019-10-09 00:25:45 -05:00
|
|
|
}
|
|
|
|
|
|
|
|
|
|
SMSpecStreamDeferredCreation::Spec::UnitConvention
|
|
|
|
|
SMSpecStreamDeferredCreation::uconv() const
|
|
|
|
|
{
|
|
|
|
|
using UType = ::Opm::UnitSystem::UnitType;
|
|
|
|
|
|
|
|
|
|
if (this->utype_ == UType::UNIT_TYPE_METRIC)
|
|
|
|
|
return Spec::UnitConvention::Metric;
|
|
|
|
|
|
|
|
|
|
if (this->utype_ == UType::UNIT_TYPE_FIELD)
|
|
|
|
|
return Spec::UnitConvention::Field;
|
|
|
|
|
|
|
|
|
|
if (this->utype_ == UType::UNIT_TYPE_LAB)
|
|
|
|
|
return Spec::UnitConvention::Lab;
|
|
|
|
|
|
|
|
|
|
if (this->utype_ == UType::UNIT_TYPE_PVT_M)
|
|
|
|
|
return Spec::UnitConvention::Pvt_M;
|
|
|
|
|
|
|
|
|
|
throw std::invalid_argument {
|
|
|
|
|
"Unsupported Unit Convention (" +
|
|
|
|
|
std::to_string(static_cast<int>(this->utype_)) + ')'
|
|
|
|
|
};
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
std::unique_ptr<SMSpecStreamDeferredCreation>
|
|
|
|
|
makeDeferredSMSpecCreation(const Opm::EclipseState& es,
|
|
|
|
|
const Opm::EclipseGrid& grid,
|
|
|
|
|
const Opm::Schedule& sched)
|
|
|
|
|
{
|
|
|
|
|
return std::make_unique<SMSpecStreamDeferredCreation>
|
|
|
|
|
(es.cfg().init(), grid, sched.posixStartTime(),
|
|
|
|
|
es.getUnits().getType());
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
std::string makeUpperCase(std::string input)
|
|
|
|
|
{
|
|
|
|
|
for (auto& c : input) {
|
|
|
|
|
const auto u = std::toupper(static_cast<unsigned char>(c));
|
|
|
|
|
c = static_cast<std::string::value_type>(u);
|
2017-11-30 14:52:22 +01:00
|
|
|
}
|
2018-05-20 14:50:31 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
return input;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Opm::EclIO::OutputStream::ResultSet
|
|
|
|
|
makeResultSet(const Opm::IOConfig& iocfg, const std::string& basenm)
|
|
|
|
|
{
|
|
|
|
|
const auto& base = basenm.empty()
|
|
|
|
|
? makeUpperCase(iocfg.getBaseName())
|
|
|
|
|
: basenm;
|
|
|
|
|
|
|
|
|
|
return { iocfg.getOutputDir(), base };
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
} // Anonymous namespace
|
|
|
|
|
|
|
|
|
|
class Opm::out::Summary::SummaryImplementation
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
explicit SummaryImplementation(const EclipseState& es,
|
|
|
|
|
const SummaryConfig& sumcfg,
|
|
|
|
|
const EclipseGrid& grid,
|
|
|
|
|
const Schedule& sched,
|
|
|
|
|
const std::string& basename);
|
|
|
|
|
|
|
|
|
|
SummaryImplementation(const SummaryImplementation& rhs) = delete;
|
|
|
|
|
SummaryImplementation(SummaryImplementation&& rhs) = default;
|
|
|
|
|
SummaryImplementation& operator=(const SummaryImplementation& rhs) = delete;
|
|
|
|
|
SummaryImplementation& operator=(SummaryImplementation&& rhs) = default;
|
|
|
|
|
|
|
|
|
|
void eval(const EclipseState& es,
|
|
|
|
|
const Schedule& sched,
|
|
|
|
|
const int sim_step,
|
|
|
|
|
const double duration,
|
|
|
|
|
const data::WellRates& well_solution,
|
2020-03-16 13:51:05 +01:00
|
|
|
const data::Group& group_solution,
|
2019-10-09 00:25:45 -05:00
|
|
|
const GlobalProcessParameters& single_values,
|
|
|
|
|
const RegionParameters& region_values,
|
|
|
|
|
const BlockValues& block_values,
|
|
|
|
|
SummaryState& st) const;
|
|
|
|
|
|
|
|
|
|
void internal_store(const SummaryState& st, const int report_step);
|
|
|
|
|
void write();
|
|
|
|
|
|
|
|
|
|
private:
|
|
|
|
|
struct MiniStep
|
2018-06-25 13:47:34 +02:00
|
|
|
{
|
2019-10-09 00:25:45 -05:00
|
|
|
int id{0};
|
|
|
|
|
int seq{-1};
|
|
|
|
|
std::vector<float> params{};
|
|
|
|
|
};
|
2018-06-25 13:47:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
using EvalPtr = SummaryOutputParameters::EvalPtr;
|
2018-06-25 13:47:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
std::reference_wrapper<const Opm::EclipseGrid> grid_;
|
|
|
|
|
Opm::out::RegionCache regCache_;
|
2018-10-09 13:59:13 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
std::unique_ptr<SMSpecStreamDeferredCreation> deferredSMSpec_;
|
2018-10-09 13:59:13 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
Opm::EclIO::OutputStream::ResultSet rset_;
|
|
|
|
|
Opm::EclIO::OutputStream::Formatted fmt_;
|
|
|
|
|
Opm::EclIO::OutputStream::Unified unif_;
|
2018-10-09 13:59:13 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
int miniStepID_{0};
|
|
|
|
|
int prevCreate_{-1};
|
|
|
|
|
int prevReportStepID_{-1};
|
|
|
|
|
std::vector<MiniStep>::size_type numUnwritten_{0};
|
|
|
|
|
|
|
|
|
|
SummaryOutputParameters outputParameters_{};
|
|
|
|
|
std::vector<EvalPtr> requiredRestartParameters_{};
|
|
|
|
|
std::vector<std::string> valueKeys_{};
|
|
|
|
|
std::vector<MiniStep> unwritten_{};
|
|
|
|
|
|
|
|
|
|
std::unique_ptr<Opm::EclIO::OutputStream::SummarySpecification> smspec_{};
|
|
|
|
|
std::unique_ptr<Opm::EclIO::EclOutput> stream_{};
|
|
|
|
|
|
|
|
|
|
void configureTimeVectors(const EclipseState& es);
|
|
|
|
|
|
|
|
|
|
void configureSummaryInput(const EclipseState& es,
|
|
|
|
|
const SummaryConfig& sumcfg,
|
|
|
|
|
const EclipseGrid& grid,
|
|
|
|
|
const Schedule& sched);
|
|
|
|
|
|
|
|
|
|
void configureRequiredRestartParameters(const SummaryConfig& sumcfg,
|
|
|
|
|
const Schedule& sched);
|
|
|
|
|
|
|
|
|
|
MiniStep& getNextMiniStep(const int report_step);
|
|
|
|
|
const MiniStep& lastUnwritten() const;
|
|
|
|
|
|
|
|
|
|
void write(const MiniStep& ms);
|
|
|
|
|
|
|
|
|
|
void createSMSpecIfNecessary();
|
|
|
|
|
void createSmryStreamIfNecessary(const int report_step);
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
SummaryImplementation(const EclipseState& es,
|
|
|
|
|
const SummaryConfig& sumcfg,
|
|
|
|
|
const EclipseGrid& grid,
|
|
|
|
|
const Schedule& sched,
|
|
|
|
|
const std::string& basename)
|
|
|
|
|
: grid_ (std::cref(grid))
|
2020-01-29 15:46:58 +01:00
|
|
|
, regCache_ (es.globalFieldProps().get_int("FIPNUM"), grid, sched)
|
2019-10-09 00:25:45 -05:00
|
|
|
, deferredSMSpec_(makeDeferredSMSpecCreation(es, grid, sched))
|
|
|
|
|
, rset_ (makeResultSet(es.cfg().io(), basename))
|
|
|
|
|
, fmt_ { es.cfg().io().getFMTOUT() }
|
|
|
|
|
, unif_ { es.cfg().io().getUNIFOUT() }
|
|
|
|
|
{
|
|
|
|
|
this->configureTimeVectors(es);
|
|
|
|
|
this->configureSummaryInput(es, sumcfg, grid, sched);
|
|
|
|
|
this->configureRequiredRestartParameters(sumcfg, sched);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
internal_store(const SummaryState& st, const int report_step)
|
|
|
|
|
{
|
|
|
|
|
auto& ms = this->getNextMiniStep(report_step);
|
|
|
|
|
|
|
|
|
|
const auto nParam = this->valueKeys_.size();
|
|
|
|
|
|
|
|
|
|
for (auto i = decltype(nParam){0}; i < nParam; ++i) {
|
|
|
|
|
if (! st.has(this->valueKeys_[i]))
|
|
|
|
|
// Parameter not yet evaluated (e.g., well/group not
|
|
|
|
|
// yet active). Nothing to do here.
|
|
|
|
|
continue;
|
|
|
|
|
|
|
|
|
|
ms.params[i] = st.get(this->valueKeys_[i]);
|
2018-06-25 13:47:34 +02:00
|
|
|
}
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
void
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
eval(const EclipseState& es,
|
|
|
|
|
const Schedule& sched,
|
|
|
|
|
const int sim_step,
|
|
|
|
|
const double duration,
|
|
|
|
|
const data::WellRates& well_solution,
|
2020-03-16 13:51:05 +01:00
|
|
|
const data::Group& group_solution,
|
2019-10-09 00:25:45 -05:00
|
|
|
const GlobalProcessParameters& single_values,
|
|
|
|
|
const RegionParameters& region_values,
|
|
|
|
|
const BlockValues& block_values,
|
|
|
|
|
Opm::SummaryState& st) const
|
|
|
|
|
{
|
|
|
|
|
const Evaluator::InputData input {
|
|
|
|
|
es, sched, this->grid_, this->regCache_
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
const Evaluator::SimulatorResults simRes {
|
2020-03-16 13:51:05 +01:00
|
|
|
well_solution, group_solution, single_values, region_values, block_values
|
2019-10-09 00:25:45 -05:00
|
|
|
};
|
|
|
|
|
|
|
|
|
|
for (auto& evalPtr : this->outputParameters_.getEvaluators()) {
|
|
|
|
|
evalPtr->update(sim_step, duration, input, simRes, st);
|
2018-11-02 13:15:43 +01:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
for (auto& evalPtr : this->requiredRestartParameters_) {
|
|
|
|
|
evalPtr->update(sim_step, duration, input, simRes, st);
|
|
|
|
|
}
|
|
|
|
|
}
|
2018-03-07 08:57:58 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
void Opm::out::Summary::SummaryImplementation::write()
|
|
|
|
|
{
|
|
|
|
|
const auto zero = std::vector<MiniStep>::size_type{0};
|
|
|
|
|
if (this->numUnwritten_ == zero)
|
|
|
|
|
// No unwritten data. Nothing to do so return early.
|
|
|
|
|
return;
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
this->createSMSpecIfNecessary();
|
2016-06-15 10:45:07 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
if (this->prevReportStepID_ < this->lastUnwritten().seq) {
|
|
|
|
|
this->smspec_->write(this->outputParameters_.summarySpecification());
|
2017-06-06 13:02:46 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
for (auto i = 0*this->numUnwritten_; i < this->numUnwritten_; ++i)
|
|
|
|
|
this->write(this->unwritten_[i]);
|
|
|
|
|
|
|
|
|
|
// Eagerly output last set of parameters to permanent storage.
|
|
|
|
|
this->stream_->flushStream();
|
|
|
|
|
|
|
|
|
|
// Reset "unwritten" counter to reflect the fact that we've
|
|
|
|
|
// output all stored ministeps.
|
|
|
|
|
this->numUnwritten_ = zero;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void Opm::out::Summary::SummaryImplementation::write(const MiniStep& ms)
|
|
|
|
|
{
|
|
|
|
|
this->createSmryStreamIfNecessary(ms.seq);
|
|
|
|
|
|
|
|
|
|
if (this->prevReportStepID_ < ms.seq) {
|
|
|
|
|
// XXX: Should probably write SEQHDR = 0 here since
|
|
|
|
|
/// we do not know the actual encoding needed.
|
|
|
|
|
this->stream_->write("SEQHDR", std::vector<int>{ ms.seq });
|
|
|
|
|
this->prevReportStepID_ = ms.seq;
|
2018-01-17 15:23:23 +01:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
this->stream_->write("MINISTEP", std::vector<int>{ ms.id });
|
|
|
|
|
this->stream_->write("PARAMS" , ms.params);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
configureTimeVectors(const EclipseState& es)
|
|
|
|
|
{
|
|
|
|
|
const auto dfltwgname = std::string(":+:+:+:+");
|
|
|
|
|
const auto dfltnum = 0;
|
|
|
|
|
|
|
|
|
|
// XXX: Save keys might need/want to include a random component too.
|
|
|
|
|
auto makeKey = [this](const std::string& keyword) -> void
|
|
|
|
|
{
|
|
|
|
|
this->valueKeys_.push_back(
|
|
|
|
|
"SMSPEC.Internal." + keyword + ".Value.SAVE"
|
|
|
|
|
);
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
// TIME
|
|
|
|
|
{
|
|
|
|
|
const auto& kw = std::string("TIME");
|
|
|
|
|
makeKey(kw);
|
|
|
|
|
|
|
|
|
|
const auto* utime = es.getUnits().name(UnitSystem::measure::time);
|
|
|
|
|
auto eval = std::make_unique<Evaluator::Time>(this->valueKeys_.back());
|
|
|
|
|
|
|
|
|
|
this->outputParameters_
|
|
|
|
|
.makeParameter(kw, dfltwgname, dfltnum, utime, std::move(eval));
|
2018-01-19 15:27:28 +01:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
#if NOT_YET
|
|
|
|
|
// YEARS
|
|
|
|
|
{
|
|
|
|
|
const auto& kw = std::string("YEARS");
|
|
|
|
|
makeKey(kw);
|
|
|
|
|
|
|
|
|
|
auto eval = std::make_unique<Evaluator::Years>(this->valueKeys_.back());
|
|
|
|
|
|
|
|
|
|
this->outputParameters_
|
|
|
|
|
.makeParameter(kw, dfltwgname, dfltnum, kw, std::move(eval));
|
2016-04-26 12:26:38 +02:00
|
|
|
}
|
2019-10-09 00:25:45 -05:00
|
|
|
#endif // NOT_YET
|
2019-05-11 11:26:34 +02:00
|
|
|
}
|
2018-06-25 13:47:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
void
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
configureSummaryInput(const EclipseState& es,
|
|
|
|
|
const SummaryConfig& sumcfg,
|
|
|
|
|
const EclipseGrid& grid,
|
|
|
|
|
const Schedule& sched)
|
|
|
|
|
{
|
|
|
|
|
const auto st = SummaryState {
|
|
|
|
|
std::chrono::system_clock::from_time_t(sched.getStartTime())
|
|
|
|
|
};
|
2019-05-11 11:26:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
Evaluator::Factory fact {
|
|
|
|
|
es, grid, st, sched.getUDQConfig(sched.size() - 1)
|
|
|
|
|
};
|
|
|
|
|
|
2020-03-13 14:24:29 +01:00
|
|
|
auto unsuppkw = std::vector<SummaryConfigNode>{};
|
2019-10-09 00:25:45 -05:00
|
|
|
for (const auto& node : sumcfg) {
|
|
|
|
|
auto prmDescr = fact.create(node);
|
|
|
|
|
|
|
|
|
|
if (! prmDescr.evaluator) {
|
|
|
|
|
// No known evaluation function/type for this keyword
|
2019-11-07 12:03:38 +01:00
|
|
|
unsuppkw.push_back(node);
|
2019-05-11 11:26:34 +02:00
|
|
|
continue;
|
2019-10-09 00:25:45 -05:00
|
|
|
}
|
2019-05-11 11:26:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
// This keyword has a known evaluation method.
|
2019-05-11 11:26:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
this->valueKeys_.push_back(std::move(prmDescr.uniquekey));
|
|
|
|
|
|
|
|
|
|
this->outputParameters_
|
|
|
|
|
.makeParameter(node.keyword(),
|
|
|
|
|
makeWGName(node.namedEntity()),
|
|
|
|
|
node.number(),
|
|
|
|
|
std::move(prmDescr.unit),
|
|
|
|
|
std::move(prmDescr.evaluator));
|
2018-06-25 13:47:34 +02:00
|
|
|
}
|
2019-10-09 00:25:45 -05:00
|
|
|
|
|
|
|
|
if (! unsuppkw.empty())
|
|
|
|
|
reportUnsupportedKeywords(std::move(unsuppkw));
|
2019-05-11 11:26:34 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
void
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
configureRequiredRestartParameters(const SummaryConfig& sumcfg,
|
|
|
|
|
const Schedule& sched)
|
|
|
|
|
{
|
2020-03-13 14:24:29 +01:00
|
|
|
auto makeEvaluator = [&sumcfg, this](const SummaryConfigNode& node) -> void
|
2019-10-09 00:25:45 -05:00
|
|
|
{
|
|
|
|
|
if (sumcfg.hasSummaryKey(node.uniqueNodeKey()))
|
|
|
|
|
// Handler already exists. Don't add second evaluation.
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
auto fcnPos = funs.find(node.keyword());
|
|
|
|
|
assert ((fcnPos != funs.end()) &&
|
|
|
|
|
"Internal error creating required restart vectors");
|
|
|
|
|
|
|
|
|
|
auto eval = std::make_unique<
|
|
|
|
|
Evaluator::FunctionRelation>(node, fcnPos->second);
|
|
|
|
|
|
|
|
|
|
this->requiredRestartParameters_.push_back(std::move(eval));
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
for (const auto& node : requiredRestartVectors(sched))
|
|
|
|
|
makeEvaluator(node);
|
2019-05-11 11:26:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
for (const auto& node : requiredSegmentVectors(sched))
|
|
|
|
|
makeEvaluator(node);
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
Opm::out::Summary::SummaryImplementation::MiniStep&
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::getNextMiniStep(const int report_step)
|
|
|
|
|
{
|
|
|
|
|
if (this->numUnwritten_ == this->unwritten_.size())
|
|
|
|
|
this->unwritten_.emplace_back();
|
|
|
|
|
|
|
|
|
|
assert ((this->numUnwritten_ < this->unwritten_.size()) &&
|
|
|
|
|
"Internal inconsistency in 'unwritten' counter");
|
2019-05-11 11:26:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
auto& ms = this->unwritten_[this->numUnwritten_++];
|
|
|
|
|
|
|
|
|
|
ms.id = this->miniStepID_++; // MINSTEP IDs start at zero.
|
|
|
|
|
ms.seq = report_step;
|
|
|
|
|
|
|
|
|
|
ms.params.resize(this->valueKeys_.size(), 0.0f);
|
|
|
|
|
|
|
|
|
|
std::fill(ms.params.begin(), ms.params.end(), 0.0f);
|
|
|
|
|
|
|
|
|
|
return ms;
|
2016-04-20 09:55:17 +02:00
|
|
|
}
|
|
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
const Opm::out::Summary::SummaryImplementation::MiniStep&
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::lastUnwritten() const
|
|
|
|
|
{
|
|
|
|
|
assert (this->numUnwritten_ <= this->unwritten_.size());
|
|
|
|
|
assert (this->numUnwritten_ > decltype(this->numUnwritten_){0});
|
2019-05-11 11:26:34 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
return this->unwritten_[this->numUnwritten_ - 1];
|
|
|
|
|
}
|
2016-04-20 09:55:17 +02:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
void Opm::out::Summary::SummaryImplementation::createSMSpecIfNecessary()
|
|
|
|
|
{
|
|
|
|
|
if (this->deferredSMSpec_) {
|
|
|
|
|
// We need an SMSPEC file and none exists. Create it and release
|
|
|
|
|
// the resources captured to make the deferred creation call.
|
|
|
|
|
this->smspec_ = this->deferredSMSpec_
|
|
|
|
|
->createStream(this->rset_, this->fmt_);
|
|
|
|
|
|
|
|
|
|
this->deferredSMSpec_.reset();
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
Opm::out::Summary::SummaryImplementation::
|
|
|
|
|
createSmryStreamIfNecessary(const int report_step)
|
|
|
|
|
{
|
|
|
|
|
// Create stream if unset or if non-unified (separate) and new step.
|
|
|
|
|
|
|
|
|
|
assert ((this->prevCreate_ <= report_step) &&
|
|
|
|
|
"Inconsistent Report Step Sequence Detected");
|
|
|
|
|
|
|
|
|
|
const auto do_create = ! this->stream_
|
|
|
|
|
|| (! this->unif_.set && (this->prevCreate_ < report_step));
|
2018-11-02 15:45:33 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
if (do_create) {
|
|
|
|
|
this->stream_ = Opm::EclIO::OutputStream::
|
|
|
|
|
createSummaryFile(this->rset_, report_step,
|
|
|
|
|
this->fmt_, this->unif_);
|
2018-11-02 15:45:33 +01:00
|
|
|
|
2019-10-09 00:25:45 -05:00
|
|
|
this->prevCreate_ = report_step;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
namespace {
|
|
|
|
|
|
|
|
|
|
void validateElapsedTime(const double secs_elapsed,
|
|
|
|
|
const Opm::EclipseState& es,
|
|
|
|
|
const Opm::SummaryState& st)
|
|
|
|
|
{
|
|
|
|
|
if (! (secs_elapsed < st.get_elapsed()))
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
const auto& usys = es.getUnits();
|
|
|
|
|
const auto elapsed = usys.from_si(measure::time, secs_elapsed);
|
|
|
|
|
const auto prev_el = usys.from_si(measure::time, st.get_elapsed());
|
|
|
|
|
const auto unt = '[' + std::string{ usys.name(measure::time) } + ']';
|
|
|
|
|
|
|
|
|
|
throw std::invalid_argument {
|
|
|
|
|
"Elapsed time ("
|
|
|
|
|
+ std::to_string(elapsed) + ' ' + unt
|
|
|
|
|
+ ") must not precede previous elapsed time ("
|
|
|
|
|
+ std::to_string(prev_el) + ' ' + unt
|
|
|
|
|
+ "). Incorrect restart time?"
|
|
|
|
|
};
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
} // Anonymous namespace
|
|
|
|
|
|
|
|
|
|
namespace Opm { namespace out {
|
|
|
|
|
|
|
|
|
|
Summary::Summary(const EclipseState& es,
|
|
|
|
|
const SummaryConfig& sumcfg,
|
|
|
|
|
const EclipseGrid& grid,
|
|
|
|
|
const Schedule& sched,
|
|
|
|
|
const std::string& basename)
|
|
|
|
|
: pImpl_(new SummaryImplementation(es, sumcfg, grid, sched, basename))
|
|
|
|
|
{}
|
|
|
|
|
|
|
|
|
|
void Summary::eval(SummaryState& st,
|
|
|
|
|
const int report_step,
|
|
|
|
|
const double secs_elapsed,
|
|
|
|
|
const EclipseState& es,
|
|
|
|
|
const Schedule& schedule,
|
|
|
|
|
const data::WellRates& well_solution,
|
2020-03-16 13:51:05 +01:00
|
|
|
const data::Group& group_solution,
|
2019-10-09 00:25:45 -05:00
|
|
|
const GlobalProcessParameters& single_values,
|
|
|
|
|
const RegionParameters& region_values,
|
|
|
|
|
const BlockValues& block_values) const
|
|
|
|
|
{
|
|
|
|
|
validateElapsedTime(secs_elapsed, es, st);
|
|
|
|
|
|
|
|
|
|
const double duration = secs_elapsed - st.get_elapsed();
|
|
|
|
|
|
|
|
|
|
/* report_step is the number of the file we are about to write - i.e. for instance CASE.S$report_step
|
|
|
|
|
* for the data in a non-unified summary file.
|
|
|
|
|
* sim_step is the timestep which has been effective in the simulator, and as such is the value
|
|
|
|
|
* necessary to use when consulting the Schedule object. */
|
|
|
|
|
const auto sim_step = std::max( 0, report_step - 1 );
|
|
|
|
|
|
|
|
|
|
this->pImpl_->eval(es, schedule, sim_step, duration,
|
2020-03-16 13:51:05 +01:00
|
|
|
well_solution, group_solution, single_values,
|
2019-10-09 00:25:45 -05:00
|
|
|
region_values, block_values, st);
|
|
|
|
|
|
|
|
|
|
eval_udq(schedule, sim_step, st);
|
|
|
|
|
|
|
|
|
|
st.update_elapsed(duration);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void Summary::add_timestep(const SummaryState& st, const int report_step)
|
|
|
|
|
{
|
|
|
|
|
this->pImpl_->internal_store(st, report_step);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void Summary::write() const
|
|
|
|
|
{
|
|
|
|
|
this->pImpl_->write();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
Summary::~Summary() {}
|
2019-05-12 10:01:26 +02:00
|
|
|
|
2018-06-25 13:47:34 +02:00
|
|
|
}} // namespace Opm::out
|