Merge pull request #1932 from blattms/tran-calculator-2

Adds a calculator for modifying transmissibilities.
This commit is contained in:
Markus Blatt
2020-09-22 19:56:12 +02:00
committed by GitHub
9 changed files with 815 additions and 263 deletions

View File

@@ -572,7 +572,10 @@ if(ENABLE_ECL_INPUT)
opm/parser/eclipse/EclipseState/Util/OrderedMap.hpp
opm/parser/eclipse/EclipseState/SummaryConfig/SummaryConfig.hpp
opm/parser/eclipse/EclipseState/Edit/EDITNNC.hpp
opm/parser/eclipse/EclipseState/Grid/FieldData.hpp
opm/parser/eclipse/EclipseState/Grid/Keywords.hpp
opm/parser/eclipse/EclipseState/Grid/GridDims.hpp
opm/parser/eclipse/EclipseState/Grid/TranCalculator.hpp
opm/parser/eclipse/EclipseState/Grid/TransMult.hpp
opm/parser/eclipse/EclipseState/Grid/PinchMode.hpp
opm/parser/eclipse/EclipseState/Grid/MULTREGTScanner.hpp

View File

@@ -0,0 +1,146 @@
/*
Copyright 2020 Equinor AS.
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/>.
*/
#ifndef FIELD_DATA_HPP
#define FIELD_DATA_HPP
#include <opm/parser/eclipse/EclipseState/Grid/Box.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/Keywords.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FieldData.hpp>
#include <opm/parser/eclipse/Deck/value_status.hpp>
#include<string>
#include<vector>
#include<optional>
#include<array>
#include<algorithm>
namespace Opm
{
namespace Fieldprops
{
template<typename T>
static void compress(std::vector<T>& data, const std::vector<bool>& active_map) {
std::size_t shift = 0;
for (std::size_t g = 0; g < active_map.size(); g++) {
if (active_map[g] && shift > 0) {
data[g - shift] = data[g];
continue;
}
if (!active_map[g])
shift += 1;
}
data.resize(data.size() - shift);
}
template<typename T>
struct FieldData {
std::vector<T> data;
std::vector<value::status> value_status;
keywords::keyword_info<T> kw_info;
std::optional<std::vector<T>> global_data;
std::optional<std::vector<value::status>> global_value_status;
mutable bool all_set;
FieldData() = default;
FieldData(const keywords::keyword_info<T>& info, std::size_t active_size, std::size_t global_size) :
data(std::vector<T>(active_size)),
value_status(active_size, value::status::uninitialized),
kw_info(info),
all_set(false)
{
if (global_size != 0) {
this->global_data = std::vector<T>(global_size);
this->global_value_status = std::vector<value::status>(global_size, value::status::uninitialized);
}
if (info.scalar_init)
this->default_assign( *info.scalar_init );
}
std::size_t size() const {
return this->data.size();
}
bool valid() const {
if (this->all_set)
return true;
static const std::array<value::status,2> invalid_value = {value::status::uninitialized, value::status::empty_default};
const auto& it = std::find_first_of(this->value_status.begin(), this->value_status.end(), invalid_value.begin(), invalid_value.end());
this->all_set = (it == this->value_status.end());
return this->all_set;
}
void compress(const std::vector<bool>& active_map) {
Fieldprops::compress(this->data, active_map);
Fieldprops::compress(this->value_status, active_map);
}
void copy(const FieldData<T>& src, const std::vector<Box::cell_index>& index_list) {
for (const auto& ci : index_list) {
this->data[ci.active_index] = src.data[ci.active_index];
this->value_status[ci.active_index] = src.value_status[ci.active_index];
}
}
void default_assign(T value) {
std::fill(this->data.begin(), this->data.end(), value);
std::fill(this->value_status.begin(), this->value_status.end(), value::status::valid_default);
if (this->global_data) {
std::fill(this->global_data->begin(), this->global_data->end(), value);
std::fill(this->global_value_status->begin(), this->global_value_status->end(), value::status::valid_default);
}
}
void default_assign(const std::vector<T>& src) {
if (src.size() != this->size())
throw std::invalid_argument("Size mismatch got: " + std::to_string(src.size()) + " expected: " + std::to_string(this->size()));
std::copy(src.begin(), src.end(), this->data.begin());
std::fill(this->value_status.begin(), this->value_status.end(), value::status::valid_default);
}
void default_update(const std::vector<T>& src) {
if (src.size() != this->size())
throw std::invalid_argument("Size mismatch got: " + std::to_string(src.size()) + " expected: " + std::to_string(this->size()));
for (std::size_t i = 0; i < src.size(); i++) {
if (!value::has_value(this->value_status[i])) {
this->value_status[i] = value::status::valid_default;
this->data[i] = src[i];
}
}
}
void update(std::size_t index, T value, value::status status) {
this->data[index] = value;
this->value_status[index] = status;
}
};
} // end namespace Fieldprops
} // end namespace Opm
#endif // FIELD_DATA_HPP

View File

@@ -21,6 +21,9 @@
#include <memory>
#include <vector>
#include <unordered_map>
#include <opm/parser/eclipse/EclipseState/Grid/TranCalculator.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FieldData.hpp>
namespace Opm {
@@ -178,6 +181,14 @@ public:
template <typename T>
std::vector<std::string> keys() const;
const Fieldprops::FieldData<int>&
get_int_field_data(const std::string& keyword) const;
/// \brief Get double field data associated with a keyword
/// \param allow_unsupported If true we deactivate some checks on the
/// keyword and thus allow getting FieldData used by the TranCalculator.
const Fieldprops::FieldData<double>&
get_double_field_data(const std::string& keyword, bool allow_unsupported=false) const;
virtual const std::vector<int>& get_int(const std::string& keyword) const { return this->get<int>(keyword); }
virtual std::vector<int> get_global_int(const std::string& keyword) const { return this->get_global<int>(keyword); }
@@ -187,6 +198,59 @@ public:
virtual bool has_int(const std::string& keyword) const { return this->has<int>(keyword); }
virtual bool has_double(const std::string& keyword) const { return this->has<double>(keyword); }
/*
The transmissibility keywords TRANX, TRANY and TRANZ do not really fit
well in the FieldProps system. The opm codebase is based on a full
internalization in the parse phase, and then passing fully assembled
objects to the simulator. When it comes to the transmissibilities this
model breaks down because the input code in opm-common is not capable of
calculating the transmissibility, that is performed in the simulator.
The EDIT section can have modifiers on TRAN, these must be applied *after*
the initial transmissibilities are calculated. To support this all the
modifiers to the TRAN{XYZ} fields are assembled in "transmissibility
calculators", and then these modifiers can be applied to a TRAN vector
after it has been calculated in the simulator. Usage from the simulator
could look like:
const auto& fp = eclState.fieldProps();
// Calculate transmissibilities using grid and permeability
std::vector<double> tranx = ....
// Check if there are any active TRANX modifiers and apply them
if (fp.tran_active("TRANX"))
fp.apply_tran("TRANX", tranx);
*/
/*
Will check if there are any TRAN{XYZ} modifiers active in the deck.
*/
virtual bool tran_active(const std::string& keyword) const;
/*
Will apply all the TRAN modifiers which are present in the deck on the
already initialized vector tran_data. The vector tran_data should be
organised as the data vectors in the fieldpropsmanager - i.e. one element
for each active cell - in lexicographical order. The operations which are
supported by the transmissibility calculator are those given by the enum
ScalarOperation in FieldProps.hpp.
*/
virtual void apply_tran(const std::string& keyword, std::vector<double>& tran_data) const;
/*
When using MPI the FieldPropsManager is typically only assembled on the
root node and then distributed to the other nodes afterwards. These
methods are support methods for that, the real data used by the
transmissibility calculators is in the form of custom 3D fields, they are
distributed the same way the rest of the 3D fields are distributed.
*/
virtual std::vector<char> serialize_tran() const;
virtual void deserialize_tran(const std::vector<char>& buffer);
private:
/*
Return the keyword values as a std::vector<>. All elements in the return
@@ -227,6 +291,16 @@ private:
std::shared_ptr<FieldProps> fp;
};
void deserialize_tran(std::unordered_map<std::string, Fieldprops::TranCalculator>& tran,
const std::vector<char>& buffer);
template<class MapType>
void apply_tran(const std::unordered_map<std::string, Fieldprops::TranCalculator>& tran,
const MapType& double_data,
std::size_t active_size,
const std::string& keyword, std::vector<double>& data);
}
#endif

View File

@@ -0,0 +1,67 @@
/*
Copyright 2020 Equinor AS.
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/>.
*/
#ifndef OPM_KEYWORDS_HPP
#define OPM_KEYWORDS_HPP
namespace Opm
{
namespace Fieldprops
{
namespace keywords {
template <typename T>
struct keyword_info {
std::optional<std::string> unit = std::nullopt;
std::optional<T> scalar_init = std::nullopt;
bool multiplier = false;
bool top = false;
bool global = false;
keyword_info<T>& init(T init_value) {
this->scalar_init = init_value;
return *this;
}
keyword_info<T>& unit_string(const std::string& unit_string) {
this->unit = unit_string;
return *this;
}
keyword_info<T>& distribute_top(bool dtop) {
this->top = dtop;
return *this;
}
keyword_info<T>& mult(bool m) {
this->multiplier = m;
return *this;
}
keyword_info<T>& global_kw(bool g) {
this->global = g;
return *this;
}
};
} // end namespace Keywords
} // end namespace Fieldprops
} //end namespace Opm
#endif //OPM_KEYWORDS_HPP

View File

@@ -0,0 +1,107 @@
/*
Copyright 2020 Equinor AS.
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/>.
*/
#ifndef TRAN_CALCULATOR_HPP
#define TRAN_CALCULATOR_HPP
#include<string>
#include <opm/parser/eclipse/EclipseState/Grid/Keywords.hpp>
namespace Opm
{
namespace Fieldprops
{
enum class ScalarOperation {
ADD = 1,
EQUAL = 2,
MUL = 3,
MIN = 4,
MAX = 5
};
class TranCalculator {
public:
struct TranAction {
ScalarOperation op;
std::string field;
};
TranCalculator(const std::string& name_arg) :
m_name(name_arg)
{}
std::string next_name() const {
return this->m_name + std::to_string( this->actions.size() );
}
std::vector<TranAction>::const_iterator begin() const {
return this->actions.begin();
}
std::vector<TranAction>::const_iterator end() const {
return this->actions.end();
}
void add_action(ScalarOperation op, const std::string& field) {
this->actions.push_back(TranAction{op, field});
}
std::size_t size() const {
return this->actions.size();
}
const std::string& name() const {
return this->m_name;
}
keywords::keyword_info<double> make_kw_info(ScalarOperation op) {
keywords::keyword_info<double> kw_info;
switch (op) {
case ScalarOperation::MUL:
kw_info.init(1);
break;
case ScalarOperation::ADD:
kw_info.init(0);
break;
case ScalarOperation::MAX:
kw_info.init(std::numeric_limits<double>::max());
break;
case ScalarOperation::MIN:
kw_info.init(std::numeric_limits<double>::lowest());
break;
default:
break;
}
return kw_info;
}
private:
std::string m_name;
std::vector<TranAction> actions;
};
} // namespace Fieldprops
} // end namespace Opm
#endif // TRAN_CALCULATOR_HPP

View File

@@ -32,8 +32,10 @@
#include <opm/parser/eclipse/EclipseState/Tables/RtempvdTable.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/EclipseGrid.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/Box.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp>
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
#include <opm/common/utility/Serializer.hpp>
#include "FieldProps.hpp"
#include "Operate.hpp"
@@ -41,6 +43,9 @@
namespace Opm {
namespace Fieldprops
{
namespace keywords {
static const std::set<std::string> oper_keywords = {"ADD", "EQUALS", "MAXVALUE", "MINVALUE", "MULTIPLY", "OPERATE"};
@@ -48,7 +53,8 @@ static const std::set<std::string> region_oper_keywords = {"ADDREG", "EQUALREG",
static const std::set<std::string> box_keywords = {"BOX", "ENDBOX"};
template <>
keyword_info<double> global_kw_info(const std::string& name) {
keyword_info<double> global_kw_info(const std::string& name,
bool allow_unsupported) {
if (GRID::double_keywords.count(name))
return GRID::double_keywords.at(name);
@@ -67,12 +73,15 @@ keyword_info<double> global_kw_info(const std::string& name) {
if (SCHEDULE::double_keywords.count(name))
return SCHEDULE::double_keywords.at(name);
if (allow_unsupported)
return keyword_info<double>{};
throw std::out_of_range("INFO: No such keyword: " + name);
}
template <>
keyword_info<int> global_kw_info(const std::string& name) {
keyword_info<int> global_kw_info(const std::string& name, bool) {
if (GRID::int_keywords.count(name))
return GRID::int_keywords.at(name);
@@ -91,8 +100,9 @@ keyword_info<int> global_kw_info(const std::string& name) {
throw std::out_of_range("No such keyword: " + name);
}
} // end namespace keywords
}
} // end namespace Fieldprops
namespace {
@@ -151,7 +161,7 @@ void verify_deck_data(const DeckKeyword& keyword, const std::vector<T>& deck_dat
template <typename T>
void assign_deck(const keywords::keyword_info<T>& kw_info, const DeckKeyword& keyword, FieldProps::FieldData<T>& field_data, const std::vector<T>& deck_data, const std::vector<value::status>& deck_status, const Box& box) {
void assign_deck(const Fieldprops::keywords::keyword_info<T>& kw_info, const DeckKeyword& keyword, Fieldprops::FieldData<T>& field_data, const std::vector<T>& deck_data, const std::vector<value::status>& deck_status, const Box& box) {
verify_deck_data(keyword, deck_data, box);
for (const auto& cell_index : box.index_list()) {
auto active_index = cell_index.active_index;
@@ -181,7 +191,7 @@ void assign_deck(const keywords::keyword_info<T>& kw_info, const DeckKeyword& ke
template <typename T>
void multiply_deck(const keywords::keyword_info<T>& kw_info, const DeckKeyword& keyword, FieldProps::FieldData<T>& field_data, const std::vector<T>& deck_data, const std::vector<value::status>& deck_status, const Box& box) {
void multiply_deck(const Fieldprops::keywords::keyword_info<T>& kw_info, const DeckKeyword& keyword, Fieldprops::FieldData<T>& field_data, const std::vector<T>& deck_data, const std::vector<value::status>& deck_status, const Box& box) {
verify_deck_data(keyword, deck_data, box);
for (const auto& cell_index : box.index_list()) {
auto active_index = cell_index.active_index;
@@ -265,21 +275,21 @@ std::string make_region_name(const std::string& deck_value) {
throw std::invalid_argument("The input string: " + deck_value + " was invalid. Expected: O/F/M");
}
FieldProps::ScalarOperation fromString(const std::string& keyword) {
Fieldprops::ScalarOperation fromString(const std::string& keyword) {
if (keyword == ParserKeywords::ADD::keywordName || keyword == ParserKeywords::ADDREG::keywordName)
return FieldProps::ScalarOperation::ADD;
return Fieldprops::ScalarOperation::ADD;
if (keyword == ParserKeywords::EQUALS::keywordName || keyword == ParserKeywords::EQUALREG::keywordName)
return FieldProps::ScalarOperation::EQUAL;
return Fieldprops::ScalarOperation::EQUAL;
if (keyword == ParserKeywords::MULTIPLY::keywordName || keyword == ParserKeywords::MULTIREG::keywordName)
return FieldProps::ScalarOperation::MUL;
return Fieldprops::ScalarOperation::MUL;
if (keyword == ParserKeywords::MINVALUE::keywordName)
return FieldProps::ScalarOperation::MIN;
return Fieldprops::ScalarOperation::MIN;
if (keyword == ParserKeywords::MAXVALUE::keywordName)
return FieldProps::ScalarOperation::MAX;
return Fieldprops::ScalarOperation::MAX;
throw std::invalid_argument("Keyword operation not recognized");
}
@@ -324,6 +334,10 @@ FieldProps::FieldProps(const Deck& deck, const Phases& phases, const EclipseGrid
grid_ptr(&grid),
tables(tables_arg)
{
this->tran.emplace( "TRANX", Fieldprops::TranCalculator("TRANX") );
this->tran.emplace( "TRANY", Fieldprops::TranCalculator("TRANY") );
this->tran.emplace( "TRANZ", Fieldprops::TranCalculator("TRANZ") );
if (deck.hasKeyword<ParserKeywords::MULTREGP>()) {
const DeckKeyword& multregpKeyword = deck.getKeyword("MULTREGP");
for (const auto& record : multregpKeyword) {
@@ -395,17 +409,17 @@ void FieldProps::reset_actnum(const std::vector<int>& new_actnum) {
for (auto& data : this->int_data)
data.second.compress(active_map);
FieldProps::compress(this->cell_volume, active_map);
FieldProps::compress(this->cell_depth, active_map);
Fieldprops::compress(this->cell_volume, active_map);
Fieldprops::compress(this->cell_depth, active_map);
this->m_actnum = std::move(new_actnum);
this->active_size = new_active_size;
}
void FieldProps::distribute_toplayer(FieldProps::FieldData<double>& field_data, const std::vector<double>& deck_data, const Box& box) {
void FieldProps::distribute_toplayer(Fieldprops::FieldData<double>& field_data, const std::vector<double>& deck_data, const Box& box) {
const std::size_t layer_size = this->nx * this->ny;
FieldProps::FieldData<double> toplayer(field_data.kw_info, layer_size, 0);
Fieldprops::FieldData<double> toplayer(field_data.kw_info, layer_size, 0);
for (const auto& cell_index : box.index_list()) {
if (cell_index.global_index < layer_size) {
toplayer.data[cell_index.global_index] = deck_data[cell_index.data_index];
@@ -436,19 +450,19 @@ void FieldProps::distribute_toplayer(FieldProps::FieldData<double>& field_data,
template <>
bool FieldProps::supported<double>(const std::string& keyword) {
if (keywords::GRID::double_keywords.count(keyword) != 0)
if (Fieldprops::keywords::GRID::double_keywords.count(keyword) != 0)
return true;
if (keywords::EDIT::double_keywords.count(keyword) != 0)
if (Fieldprops::keywords::EDIT::double_keywords.count(keyword) != 0)
return true;
if (keywords::PROPS::double_keywords.count(keyword) != 0)
if (Fieldprops::keywords::PROPS::double_keywords.count(keyword) != 0)
return true;
if (keywords::PROPS::satfunc.count(keyword) != 0)
if (Fieldprops::keywords::PROPS::satfunc.count(keyword) != 0)
return true;
if (keywords::SOLUTION::double_keywords.count(keyword) != 0)
if (Fieldprops::keywords::SOLUTION::double_keywords.count(keyword) != 0)
return true;
return false;
@@ -456,27 +470,26 @@ bool FieldProps::supported<double>(const std::string& keyword) {
template <>
bool FieldProps::supported<int>(const std::string& keyword) {
if (keywords::REGIONS::int_keywords.count(keyword) != 0)
if (Fieldprops::keywords::REGIONS::int_keywords.count(keyword) != 0)
return true;
if (keywords::GRID::int_keywords.count(keyword) != 0)
if (Fieldprops::keywords::GRID::int_keywords.count(keyword) != 0)
return true;
if (keywords::SCHEDULE::int_keywords.count(keyword) != 0)
if (Fieldprops::keywords::SCHEDULE::int_keywords.count(keyword) != 0)
return true;
return keywords::isFipxxx(keyword);
return Fieldprops::keywords::isFipxxx(keyword);
}
template <>
FieldProps::FieldData<double>& FieldProps::init_get(const std::string& keyword) {
Fieldprops::FieldData<double>& FieldProps::init_get(const std::string& keyword, const Fieldprops::keywords::keyword_info<double>& kw_info) {
auto iter = this->double_data.find(keyword);
if (iter != this->double_data.end())
return iter->second;
const keywords::keyword_info<double>& kw_info = keywords::global_kw_info<double>(keyword);
this->double_data[keyword] = FieldData<double>(kw_info, this->active_size, kw_info.global ? this->global_size : 0);
this->double_data[keyword] = Fieldprops::FieldData<double>(kw_info, this->active_size, kw_info.global ? this->global_size : 0);
if (keyword == ParserKeywords::PORV::keywordName)
this->init_porv(this->double_data[keyword]);
@@ -484,34 +497,42 @@ FieldProps::FieldData<double>& FieldProps::init_get(const std::string& keyword)
if (keyword == ParserKeywords::TEMPI::keywordName)
this->init_tempi(this->double_data[keyword]);
if (keywords::PROPS::satfunc.count(keyword) == 1)
if (Fieldprops::keywords::PROPS::satfunc.count(keyword) == 1)
this->init_satfunc(keyword, this->double_data[keyword]);
return this->double_data[keyword];
}
template <>
Fieldprops::FieldData<double>& FieldProps::init_get(const std::string& keyword,
bool allow_unsupported) {
Fieldprops::keywords::keyword_info<double> kw_info = Fieldprops::keywords::global_kw_info<double>(keyword, allow_unsupported);
return this->init_get(keyword, kw_info);
}
template <>
FieldProps::FieldData<int>& FieldProps::init_get(const std::string& keyword) {
Fieldprops::FieldData<int>& FieldProps::init_get(const std::string& keyword, const Fieldprops::keywords::keyword_info<int>& kw_info) {
auto iter = this->int_data.find(keyword);
if (iter != this->int_data.end())
return iter->second;
if (keywords::isFipxxx(keyword)) {
auto kw_info = keywords::keyword_info<int>{};
kw_info.init(1);
this->int_data[keyword] = FieldData<int>(kw_info, this->active_size, 0);
} else {
const keywords::keyword_info<int>& kw_info = keywords::global_kw_info<int>(keyword);
this->int_data[keyword] = FieldData<int>(kw_info, this->active_size, kw_info.global ? this->global_size : 0);
}
this->int_data[keyword] = Fieldprops::FieldData<int>(kw_info, this->active_size, kw_info.global ? this->global_size : 0);
return this->int_data[keyword];
}
template <>
Fieldprops::FieldData<int>& FieldProps::init_get(const std::string& keyword, bool) {
if (Fieldprops::keywords::isFipxxx(keyword)) {
auto kw_info = Fieldprops::keywords::keyword_info<int>{};
kw_info.init(1);
return this->init_get(keyword, kw_info);
} else {
const Fieldprops::keywords::keyword_info<int>& kw_info = Fieldprops::keywords::global_kw_info<int>(keyword);
return this->init_get(keyword, kw_info);
}
}
std::vector<Box::cell_index> FieldProps::region_index( const std::string& region_name, int region_value ) {
const auto& region = this->init_get<int>(region_name);
@@ -609,7 +630,7 @@ std::vector<double> FieldProps::extract<double>(const std::string& keyword) {
double FieldProps::getSIValue(const std::string& keyword, double raw_value) const {
const auto& kw_info = keywords::global_kw_info<double>(keyword);
const auto& kw_info = Fieldprops::keywords::global_kw_info<double>(keyword);
if (kw_info.unit) {
const auto& dim = this->unit_system.parse( *kw_info.unit );
return dim.convertRawToSi(raw_value);
@@ -621,7 +642,7 @@ double FieldProps::getSIValue(const std::string& keyword, double raw_value) cons
void FieldProps::handle_int_keyword(const keywords::keyword_info<int>& kw_info, const DeckKeyword& keyword, const Box& box) {
void FieldProps::handle_int_keyword(const Fieldprops::keywords::keyword_info<int>& kw_info, const DeckKeyword& keyword, const Box& box) {
auto& field_data = this->init_get<int>(keyword.name());
const auto& deck_data = keyword.getIntData();
const auto& deck_status = keyword.getValueStatus();
@@ -629,8 +650,8 @@ void FieldProps::handle_int_keyword(const keywords::keyword_info<int>& kw_info,
}
void FieldProps::handle_double_keyword(Section section, const keywords::keyword_info<double>& kw_info, const DeckKeyword& keyword, const Box& box) {
auto& field_data = this->init_get<double>(keyword.name());
void FieldProps::handle_double_keyword(Section section, const Fieldprops::keywords::keyword_info<double>& kw_info, const DeckKeyword& keyword, const std::string& keyword_name, const Box& box) {
auto& field_data = this->init_get<double>(keyword_name, kw_info);
const auto& deck_data = keyword.getSIDoubleData();
const auto& deck_status = keyword.getValueStatus();
@@ -649,24 +670,27 @@ void FieldProps::handle_double_keyword(Section section, const keywords::keyword_
}
}
void FieldProps::handle_double_keyword(Section section, const Fieldprops::keywords::keyword_info<double>& kw_info, const DeckKeyword& keyword, const Box& box) {
this->handle_double_keyword(section, kw_info, keyword, keyword.name(), box );
}
template <typename T>
void FieldProps::apply(ScalarOperation op, std::vector<T>& data, std::vector<value::status>& value_status, T scalar_value, const std::vector<Box::cell_index>& index_list) {
if (op == ScalarOperation::EQUAL)
void FieldProps::apply(Fieldprops::ScalarOperation op, std::vector<T>& data, std::vector<value::status>& value_status, T scalar_value, const std::vector<Box::cell_index>& index_list) {
if (op == Fieldprops::ScalarOperation::EQUAL)
assign_scalar(data, value_status, scalar_value, index_list);
else if (op == ScalarOperation::MUL)
else if (op == Fieldprops::ScalarOperation::MUL)
multiply_scalar(data, value_status, scalar_value, index_list);
else if (op == ScalarOperation::ADD)
else if (op == Fieldprops::ScalarOperation::ADD)
add_scalar(data, value_status, scalar_value, index_list);
else if (op == ScalarOperation::MIN)
else if (op == Fieldprops::ScalarOperation::MIN)
min_value(data, value_status, scalar_value, index_list);
else if (op == ScalarOperation::MAX)
else if (op == Fieldprops::ScalarOperation::MAX)
max_value(data, value_status, scalar_value, index_list);
}
@@ -685,7 +709,7 @@ double FieldProps::get_beta(const std::string& func_name, const std::string& tar
}
template <typename T>
void FieldProps::operate(const DeckRecord& record, FieldData<T>& target_data, const FieldData<T>& src_data, const std::vector<Box::cell_index>& index_list) {
void FieldProps::operate(const DeckRecord& record, Fieldprops::FieldData<T>& target_data, const Fieldprops::FieldData<T>& src_data, const std::vector<Box::cell_index>& index_list) {
const std::string& func_name = record.getItem("OPERATION").get< std::string >(0);
const std::string& target_array = record.getItem("TARGET_ARRAY").get<std::string>(0);
const double alpha = this->get_alpha(func_name, target_array, record.getItem("PARAM1").get< double >(0));
@@ -762,20 +786,37 @@ void FieldProps::handle_operation(const DeckKeyword& keyword, Box box) {
const std::string& target_kw = record.getItem(0).get<std::string>(0);
box.update(record);
if (FieldProps::supported<double>(target_kw)) {
auto& field_data = this->init_get<double>(target_kw);
if (FieldProps::supported<double>(target_kw) || this->tran.count(target_kw) > 0) {
if (keyword.name() == ParserKeywords::OPERATE::keywordName) {
auto& field_data = this->init_get<double>(target_kw);
const std::string& src_kw = record.getItem("ARRAY").get<std::string>(0);
const auto& src_data = this->init_get<double>(src_kw);
FieldProps::operate(record, field_data, src_data, box.index_list());
} else {
std::string unique_name = target_kw;
auto operation = fromString(keyword.name());
double scalar_value = record.getItem(1).get<double>(0);
Fieldprops::keywords::keyword_info<double> kw_info;
auto tran_iter = this->tran.find(target_kw);
if (tran_iter != this->tran.end()) {
kw_info = tran_iter->second.make_kw_info(operation);
unique_name = tran_iter->second.next_name();
tran_iter->second.add_action(operation, unique_name);
} else
kw_info = Fieldprops::keywords::global_kw_info<double>(target_kw);
auto& field_data = this->init_get<double>(unique_name, kw_info);
if (keyword.name() != ParserKeywords::MULTIPLY::keywordName)
scalar_value = this->getSIValue(target_kw, scalar_value);
FieldProps::apply(fromString(keyword.name()), field_data.data, field_data.value_status, scalar_value, box.index_list());
if (field_data.global_data)
FieldProps::apply(fromString(keyword.name()), *field_data.global_data, *field_data.global_value_status, scalar_value, box.global_index_list());
if (tran_iter != this->tran.end()) {
assign_scalar(field_data.data, field_data.value_status, scalar_value, box.index_list());
} else {
FieldProps::apply(operation, field_data.data, field_data.value_status, scalar_value, box.index_list());
if (field_data.global_data)
FieldProps::apply(operation, *field_data.global_data, *field_data.global_value_status, scalar_value, box.global_index_list());
}
}
continue;
@@ -834,13 +875,13 @@ void FieldProps::handle_COPY(const DeckKeyword& keyword, Box box, bool region) {
void FieldProps::handle_keyword(const DeckKeyword& keyword, Box& box) {
const std::string& name = keyword.name();
if (keywords::oper_keywords.count(name) == 1)
if (Fieldprops::keywords::oper_keywords.count(name) == 1)
this->handle_operation(keyword, box);
else if (keywords::region_oper_keywords.count(name) == 1)
else if (Fieldprops::keywords::region_oper_keywords.count(name) == 1)
this->handle_region_operation(keyword);
else if (keywords::box_keywords.count(name) == 1)
else if (Fieldprops::keywords::box_keywords.count(name) == 1)
handle_box_keyword(keyword, box);
else if (name == ParserKeywords::COPY::keywordName)
@@ -853,7 +894,7 @@ void FieldProps::handle_keyword(const DeckKeyword& keyword, Box& box) {
/**********************************************************************/
void FieldProps::init_tempi(FieldData<double>& tempi) {
void FieldProps::init_tempi(Fieldprops::FieldData<double>& tempi) {
if (this->tables.hasTables("RTEMPVD")) {
const auto& eqlnum = this->get<int>("EQLNUM");
const auto& rtempvd = this->tables.getRtempvdTables();
@@ -870,7 +911,7 @@ void FieldProps::init_tempi(FieldData<double>& tempi) {
tempi.default_assign(this->tables.rtemp());
}
void FieldProps::init_porv(FieldData<double>& porv) {
void FieldProps::init_porv(Fieldprops::FieldData<double>& porv) {
auto& porv_data = porv.data;
auto& porv_status = porv.value_status;
@@ -956,13 +997,13 @@ void FieldProps::scanGRIDSection(const GRIDSection& grid_section) {
for (const auto& keyword : grid_section) {
const std::string& name = keyword.name();
if (keywords::GRID::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::GRID, keywords::GRID::double_keywords.at(name), keyword, box);
if (Fieldprops::keywords::GRID::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::GRID, Fieldprops::keywords::GRID::double_keywords.at(name), keyword, box);
continue;
}
if (keywords::GRID::int_keywords.count(name) == 1) {
this->handle_int_keyword(keywords::GRID::int_keywords.at(name), keyword, box);
if (Fieldprops::keywords::GRID::int_keywords.count(name) == 1) {
this->handle_int_keyword(Fieldprops::keywords::GRID::int_keywords.at(name), keyword, box);
continue;
}
@@ -974,13 +1015,24 @@ void FieldProps::scanEDITSection(const EDITSection& edit_section) {
Box box(*this->grid_ptr);
for (const auto& keyword : edit_section) {
const std::string& name = keyword.name();
if (keywords::EDIT::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::EDIT, keywords::EDIT::double_keywords.at(name), keyword, box);
auto tran_iter = this->tran.find(name);
if (tran_iter!= this->tran.end()) {
auto& tran_calc = tran_iter->second;
auto unique_name = tran_calc.next_name();
Fieldprops::keywords::keyword_info<double> kw_info;
this->handle_double_keyword(Section::EDIT, kw_info, keyword, unique_name, box);
tran_calc.add_action( Fieldprops::ScalarOperation::EQUAL, unique_name );
continue;
}
if (keywords::EDIT::int_keywords.count(name) == 1) {
this->handle_int_keyword(keywords::GRID::int_keywords.at(name), keyword, box);
if (Fieldprops::keywords::EDIT::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::EDIT, Fieldprops::keywords::EDIT::double_keywords.at(name), keyword, box);
continue;
}
if (Fieldprops::keywords::EDIT::int_keywords.count(name) == 1) {
this->handle_int_keyword(Fieldprops::keywords::GRID::int_keywords.at(name), keyword, box);
continue;
}
@@ -989,7 +1041,7 @@ void FieldProps::scanEDITSection(const EDITSection& edit_section) {
}
void FieldProps::init_satfunc(const std::string& keyword, FieldData<double>& satfunc) {
void FieldProps::init_satfunc(const std::string& keyword, Fieldprops::FieldData<double>& satfunc) {
if (this->m_rtep == nullptr)
this->m_rtep = satfunc::getRawTableEndpoints(this->tables, this->m_phases,
this->m_satfuncctrl.minimumRelpermMobilityThreshold());
@@ -1008,19 +1060,19 @@ void FieldProps::scanPROPSSection(const PROPSSection& props_section) {
for (const auto& keyword : props_section) {
const std::string& name = keyword.name();
if (keywords::PROPS::satfunc.count(name) == 1) {
keywords::keyword_info<double> sat_info{};
if (Fieldprops::keywords::PROPS::satfunc.count(name) == 1) {
Fieldprops::keywords::keyword_info<double> sat_info{};
this->handle_double_keyword(Section::PROPS, sat_info, keyword, box);
continue;
}
if (keywords::PROPS::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::PROPS, keywords::PROPS::double_keywords.at(name), keyword, box);
if (Fieldprops::keywords::PROPS::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::PROPS, Fieldprops::keywords::PROPS::double_keywords.at(name), keyword, box);
continue;
}
if (keywords::PROPS::int_keywords.count(name) == 1) {
this->handle_int_keyword(keywords::PROPS::int_keywords.at(name), keyword, box);
if (Fieldprops::keywords::PROPS::int_keywords.count(name) == 1) {
this->handle_int_keyword(Fieldprops::keywords::PROPS::int_keywords.at(name), keyword, box);
continue;
}
@@ -1034,13 +1086,13 @@ void FieldProps::scanREGIONSSection(const REGIONSSection& regions_section) {
for (const auto& keyword : regions_section) {
const std::string& name = keyword.name();
if (keywords::REGIONS::int_keywords.count(name)) {
this->handle_int_keyword(keywords::REGIONS::int_keywords.at(name), keyword, box);
if (Fieldprops::keywords::REGIONS::int_keywords.count(name)) {
this->handle_int_keyword(Fieldprops::keywords::REGIONS::int_keywords.at(name), keyword, box);
continue;
}
if (keywords::isFipxxx(name)) {
auto kw_info = keywords::keyword_info<int>{};
if (Fieldprops::keywords::isFipxxx(name)) {
auto kw_info = Fieldprops::keywords::keyword_info<int>{};
kw_info.init(1);
this->handle_int_keyword(kw_info, keyword, box);
continue;
@@ -1055,8 +1107,8 @@ void FieldProps::scanSOLUTIONSection(const SOLUTIONSection& solution_section) {
Box box(*this->grid_ptr);
for (const auto& keyword : solution_section) {
const std::string& name = keyword.name();
if (keywords::SOLUTION::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::SOLUTION, keywords::SOLUTION::double_keywords.at(name), keyword, box);
if (Fieldprops::keywords::SOLUTION::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::SOLUTION, Fieldprops::keywords::SOLUTION::double_keywords.at(name), keyword, box);
continue;
}
@@ -1068,13 +1120,13 @@ void FieldProps::scanSCHEDULESection(const SCHEDULESection& schedule_section) {
Box box(*this->grid_ptr);
for (const auto& keyword : schedule_section) {
const std::string& name = keyword.name();
if (keywords::SCHEDULE::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::SCHEDULE, keywords::SCHEDULE::double_keywords.at(name), keyword, box);
if (Fieldprops::keywords::SCHEDULE::double_keywords.count(name) == 1) {
this->handle_double_keyword(Section::SCHEDULE, Fieldprops::keywords::SCHEDULE::double_keywords.at(name), keyword, box);
continue;
}
if (keywords::SCHEDULE::int_keywords.count(name) == 1) {
this->handle_int_keyword(keywords::SCHEDULE::int_keywords.at(name), keyword, box);
if (Fieldprops::keywords::SCHEDULE::int_keywords.count(name) == 1) {
this->handle_int_keyword(Fieldprops::keywords::SCHEDULE::int_keywords.at(name), keyword, box);
continue;
}
@@ -1086,6 +1138,36 @@ const std::string& FieldProps::default_region() const {
return this->m_default_region;
}
void FieldProps::apply_tran(const std::string& keyword, std::vector<double>& data) {
Opm::apply_tran(this->tran, this->double_data, this->active_size, keyword, data);
}
std::vector<char> FieldProps::serialize_tran() const {
Serializer ser;
ser.put(this->tran.size());
for (const auto& tran_pair : this->tran) {
const auto& calc = tran_pair.second;
ser.put(calc.name());
ser.put(calc.size());
for (const auto& action : calc) {
ser.put(static_cast<int>(action.op));
ser.put(action.field);
}
}
return std::move(ser.buffer);
}
void FieldProps::deserialize_tran(const std::vector<char>& buffer) {
Opm::deserialize_tran(this->tran, buffer);
}
bool FieldProps::tran_active(const std::string& keyword) const {
auto calculator = this->tran.find(keyword);
return calculator != this->tran.end() && calculator->second.size() > 0;
}
template std::vector<bool> FieldProps::defaulted<int>(const std::string& keyword);
template std::vector<bool> FieldProps::defaulted<double>(const std::string& keyword);

View File

@@ -20,6 +20,7 @@
#define FIELDPROPS_HPP
#include <memory>
#include <limits>
#include <optional>
#include <string>
#include <unordered_set>
@@ -31,6 +32,9 @@
#include <opm/parser/eclipse/EclipseState/Grid/Box.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp>
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/Keywords.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/TranCalculator.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FieldData.hpp>
namespace Opm {
@@ -38,6 +42,9 @@ class Deck;
class EclipseGrid;
class TableManager;
namespace Fieldprops
{
namespace keywords {
/*
@@ -82,42 +89,6 @@ namespace keywords {
template <typename T>
struct keyword_info {
std::optional<std::string> unit = std::nullopt;
std::optional<T> scalar_init = std::nullopt;
bool multiplier = false;
bool top = false;
bool global = false;
keyword_info<T>& init(T init_value) {
this->scalar_init = init_value;
return *this;
}
keyword_info<T>& unit_string(const std::string& unit_string) {
this->unit = unit_string;
return *this;
}
keyword_info<T>& distribute_top(bool dtop) {
this->top = dtop;
return *this;
}
keyword_info<T>& mult(bool m) {
this->multiplier = m;
return *this;
}
keyword_info<T>& global_kw(bool g) {
this->global = g;
return *this;
}
};
inline bool isFipxxx(const std::string& keyword) {
// FIPxxxx can be any keyword, e.g. FIPREG or FIPXYZ that has the pattern "FIP.+"
// However, it can not be FIPOWG as that is an actual keyword.
@@ -162,9 +133,9 @@ static const std::unordered_map<std::string, keyword_info<int>> int_keywords = {
namespace EDIT {
static const std::unordered_map<std::string, keyword_info<double>> double_keywords = {{"MULTPV", keyword_info<double>{}.init(1.0)},
{"PORV", keyword_info<double>{}.unit_string("ReservoirVolume")},
{"TRANX", keyword_info<double>{}.unit_string("Transmissibility").init(1.0)},
{"TRANY", keyword_info<double>{}.unit_string("Transmissibility").init(1.0)},
{"TRANZ", keyword_info<double>{}.unit_string("Transmissibility").init(1.0)},
{"TRANX", keyword_info<double>{}.unit_string("Transmissibility")},
{"TRANY", keyword_info<double>{}.unit_string("Transmissibility")},
{"TRANZ", keyword_info<double>{}.unit_string("Transmissibility")},
{"MULTX", keyword_info<double>{}.init(1.0).mult(true)},
{"MULTX-", keyword_info<double>{}.init(1.0).mult(true)},
{"MULTY", keyword_info<double>{}.init(1.0).mult(true)},
@@ -271,14 +242,17 @@ static const std::unordered_map<std::string, keyword_info<int>> int_keywords = {
}
template <typename T>
keyword_info<T> global_kw_info(const std::string& name);
keyword_info<T> global_kw_info(const std::string& name, bool allow_unsupported = false);
}
} // end namespace keywords
} // end namespace FieldProps
class FieldProps {
public:
using ScalarOperation = Fieldprops::ScalarOperation;
struct MultregpRecord {
int region_value;
double multiplier;
@@ -293,29 +267,7 @@ public:
};
enum class ScalarOperation {
ADD = 1,
EQUAL = 2,
MUL = 3,
MIN = 4,
MAX = 5
};
template<typename T>
static void compress(std::vector<T>& data, const std::vector<bool>& active_map) {
std::size_t shift = 0;
for (std::size_t g = 0; g < active_map.size(); g++) {
if (active_map[g] && shift > 0) {
data[g - shift] = data[g];
continue;
}
if (!active_map[g])
shift += 1;
}
data.resize(data.size() - shift);
}
enum class GetStatus {
OK = 1,
@@ -326,106 +278,14 @@ public:
template<typename T>
struct FieldData {
std::vector<T> data;
std::vector<value::status> value_status;
keywords::keyword_info<T> kw_info;
std::optional<std::vector<T>> global_data;
std::optional<std::vector<value::status>> global_value_status;
mutable bool all_set;
FieldData() = default;
FieldData(const keywords::keyword_info<T>& info, std::size_t active_size, std::size_t global_size) :
data(std::vector<T>(active_size)),
value_status(active_size, value::status::uninitialized),
kw_info(info),
all_set(false)
{
if (global_size != 0) {
this->global_data = std::vector<T>(global_size);
this->global_value_status = std::vector<value::status>(global_size, value::status::uninitialized);
}
if (info.scalar_init)
this->default_assign( *info.scalar_init );
}
std::size_t size() const {
return this->data.size();
}
bool valid() const {
if (this->all_set)
return true;
static const std::array<value::status,2> invalid_value = {value::status::uninitialized, value::status::empty_default};
const auto& it = std::find_first_of(this->value_status.begin(), this->value_status.end(), invalid_value.begin(), invalid_value.end());
this->all_set = (it == this->value_status.end());
return this->all_set;
}
void compress(const std::vector<bool>& active_map) {
FieldProps::compress(this->data, active_map);
FieldProps::compress(this->value_status, active_map);
}
void copy(const FieldData<T>& src, const std::vector<Box::cell_index>& index_list) {
for (const auto& ci : index_list) {
this->data[ci.active_index] = src.data[ci.active_index];
this->value_status[ci.active_index] = src.value_status[ci.active_index];
}
}
void default_assign(T value) {
std::fill(this->data.begin(), this->data.end(), value);
std::fill(this->value_status.begin(), this->value_status.end(), value::status::valid_default);
if (this->global_data) {
std::fill(this->global_data->begin(), this->global_data->end(), value);
std::fill(this->global_value_status->begin(), this->global_value_status->end(), value::status::valid_default);
}
}
void default_assign(const std::vector<T>& src) {
if (src.size() != this->size())
throw std::invalid_argument("Size mismatch got: " + std::to_string(src.size()) + " expected: " + std::to_string(this->size()));
std::copy(src.begin(), src.end(), this->data.begin());
std::fill(this->value_status.begin(), this->value_status.end(), value::status::valid_default);
}
void default_update(const std::vector<T>& src) {
if (src.size() != this->size())
throw std::invalid_argument("Size mismatch got: " + std::to_string(src.size()) + " expected: " + std::to_string(this->size()));
for (std::size_t i = 0; i < src.size(); i++) {
if (!value::has_value(this->value_status[i])) {
this->value_status[i] = value::status::valid_default;
this->data[i] = src[i];
}
}
}
void update(std::size_t index, T value, value::status status) {
this->data[index] = value;
this->value_status[index] = status;
}
};
template<typename T>
struct FieldDataManager {
const std::string& keyword;
GetStatus status;
const FieldData<T> * data_ptr;
using Data = Fieldprops::FieldData<T>;
const Data * data_ptr;
FieldDataManager(const std::string& k, GetStatus s, const FieldData<T> * d) :
FieldDataManager(const std::string& k, GetStatus s, const Data * d) :
keyword(k),
status(s),
data_ptr(d)
@@ -457,7 +317,7 @@ public:
return this->data_ptr->data;
}
const FieldData<T>& field_data() const {
const Data& field_data() const {
this->verify_status();
return *this->data_ptr;
}
@@ -489,14 +349,16 @@ public:
template <typename T>
FieldDataManager<T> try_get(const std::string& keyword) {
if (!FieldProps::supported<T>(keyword))
FieldDataManager<T> try_get(const std::string& keyword,
bool allow_unsupported=false) {
if (!allow_unsupported && !FieldProps::supported<T>(keyword))
return FieldDataManager<T>(keyword, GetStatus::NOT_SUPPPORTED_KEYWORD, nullptr);
const FieldData<T> * field_data;
const Fieldprops::FieldData<T> * field_data;
bool has0 = this->has<T>(keyword);
field_data = std::addressof(this->init_get<T>(keyword));
field_data = std::addressof(this->init_get<T>(keyword,
std::is_same<T,double>::value && allow_unsupported));
if (field_data->valid())
return FieldDataManager<T>(keyword, GetStatus::OK, field_data);
@@ -519,7 +381,7 @@ public:
std::vector<T> get_global(const std::string& keyword) {
const auto& managed_field_data = this->try_get<T>(keyword);
const auto& field_data = managed_field_data.field_data();
const auto& kw_info = keywords::global_kw_info<T>(keyword);
const auto& kw_info = Fieldprops::keywords::global_kw_info<T>(keyword);
if (kw_info.global)
return *field_data.global_data;
else
@@ -539,7 +401,7 @@ public:
return field_data.data;
} else {
if (global) {
const auto& kw_info = keywords::global_kw_info<T>(keyword);
const auto& kw_info = Fieldprops::keywords::global_kw_info<T>(keyword);
return this->global_copy(this->extract<T>(keyword), kw_info.scalar_init);
} else
return this->extract<T>(keyword);
@@ -584,6 +446,10 @@ public:
return this->double_data.size();
}
bool tran_active(const std::string& keyword) const;
void apply_tran(const std::string& keyword, std::vector<double>& data);
std::vector<char> serialize_tran() const;
void deserialize_tran(const std::vector<char>& buffer);
private:
void scanGRIDSection(const GRIDSection& grid_section);
void scanEDITSection(const EDITSection& edit_section);
@@ -599,29 +465,33 @@ private:
std::vector<T> extract(const std::string& keyword);
template <typename T>
void operate(const DeckRecord& record, FieldData<T>& target_data, const FieldData<T>& src_data, const std::vector<Box::cell_index>& index_list);
void operate(const DeckRecord& record, Fieldprops::FieldData<T>& target_data, const Fieldprops::FieldData<T>& src_data, const std::vector<Box::cell_index>& index_list);
template <typename T>
static void apply(ScalarOperation op, std::vector<T>& data, std::vector<value::status>& value_status, T scalar_value, const std::vector<Box::cell_index>& index_list);
template <typename T>
FieldData<T>& init_get(const std::string& keyword);
Fieldprops::FieldData<T>& init_get(const std::string& keyword, bool allow_unsupported = false);
template <typename T>
Fieldprops::FieldData<T>& init_get(const std::string& keyword, const Fieldprops::keywords::keyword_info<T>& kw_info);
std::string region_name(const DeckItem& region_item);
std::vector<Box::cell_index> region_index( const std::string& region_name, int region_value );
void handle_operation(const DeckKeyword& keyword, Box box);
void handle_region_operation(const DeckKeyword& keyword);
void handle_COPY(const DeckKeyword& keyword, Box box, bool region);
void distribute_toplayer(FieldProps::FieldData<double>& field_data, const std::vector<double>& deck_data, const Box& box);
void distribute_toplayer(Fieldprops::FieldData<double>& field_data, const std::vector<double>& deck_data, const Box& box);
double get_beta(const std::string& func_name, const std::string& target_array, double raw_beta);
double get_alpha(const std::string& func_name, const std::string& target_array, double raw_alpha);
void handle_keyword(const DeckKeyword& keyword, Box& box);
void handle_double_keyword(Section section, const keywords::keyword_info<double>& kw_info, const DeckKeyword& keyword, const Box& box);
void handle_int_keyword(const keywords::keyword_info<int>& kw_info, const DeckKeyword& keyword, const Box& box);
void init_satfunc(const std::string& keyword, FieldData<double>& satfunc);
void init_porv(FieldData<double>& porv);
void init_tempi(FieldData<double>& tempi);
void handle_double_keyword(Section section, const Fieldprops::keywords::keyword_info<double>& kw_info, const DeckKeyword& keyword, const std::string& keyword_name, const Box& box);
void handle_double_keyword(Section section, const Fieldprops::keywords::keyword_info<double>& kw_info, const DeckKeyword& keyword, const Box& box);
void handle_int_keyword(const Fieldprops::keywords::keyword_info<int>& kw_info, const DeckKeyword& keyword, const Box& box);
void init_satfunc(const std::string& keyword, Fieldprops::FieldData<double>& satfunc);
void init_porv(Fieldprops::FieldData<double>& porv);
void init_tempi(Fieldprops::FieldData<double>& tempi);
const UnitSystem unit_system;
std::size_t nx,ny,nz;
@@ -635,8 +505,10 @@ private:
const TableManager& tables;
std::shared_ptr<satfunc::RawTableEndPoints> m_rtep;
std::vector<MultregpRecord> multregp;
std::unordered_map<std::string, FieldData<int>> int_data;
std::unordered_map<std::string, FieldData<double>> double_data;
std::unordered_map<std::string, Fieldprops::FieldData<int>> int_data;
std::unordered_map<std::string, Fieldprops::FieldData<double>> double_data;
std::unordered_map<std::string, Fieldprops::TranCalculator> tran;
};
}

View File

@@ -21,6 +21,7 @@
#include <opm/parser/eclipse/EclipseState/Grid/EclipseGrid.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FieldPropsManager.hpp>
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
#include <opm/common/utility/Serializer.hpp>
#include "FieldProps.hpp"
@@ -54,6 +55,24 @@ const std::vector<T>* FieldPropsManager::try_get(const std::string& keyword) con
return nullptr;
}
const Fieldprops::FieldData<int>&
FieldPropsManager::get_int_field_data(const std::string& keyword) const
{
const auto& data = this->fp->try_get<int>(keyword);
if (!data.valid())
throw std::out_of_range("Invalid field data requested.");
return data.field_data();
}
const Fieldprops::FieldData<double>&
FieldPropsManager::get_double_field_data(const std::string& keyword,
bool allow_unsupported) const
{
const auto& data = this->fp->try_get<double>(keyword, allow_unsupported);
if (!data.valid())
throw std::out_of_range("Invalid field data requested.");
return data.field_data();
}
template <typename T>
std::vector<T> FieldPropsManager::get_global(const std::string& keyword) const {
@@ -109,6 +128,95 @@ std::size_t FieldPropsManager::active_size() const {
return this->fp->active_size;
}
void FieldPropsManager::apply_tran(const std::string& keyword, std::vector<double>& data) const {
this->fp->apply_tran(keyword, data);
}
std::vector<char> FieldPropsManager::serialize_tran() const {
return this->fp->serialize_tran();
}
void FieldPropsManager::deserialize_tran(const std::vector<char>& buffer) {
this->fp->deserialize_tran(buffer);
}
bool FieldPropsManager::tran_active(const std::string& keyword) const {
return this->fp->tran_active(keyword);
}
template<class MapType>
void apply_tran(const std::unordered_map<std::string, Fieldprops::TranCalculator>& tran,
const MapType& double_data,
std::size_t active_size,
const std::string& keyword, std::vector<double>& data)
{
const auto& calculator = tran.at(keyword);
for (const auto& action : calculator) {
const auto& action_data = double_data.at(action.field);
for (std::size_t index = 0; index < active_size; index++) {
if (action_data.value_status[index] != value::status::deck_value)
continue;
switch (action.op) {
case Fieldprops::ScalarOperation::EQUAL:
data[index] = action_data.data[index];
break;
case Fieldprops::ScalarOperation::MUL:
data[index] *= action_data.data[index];
break;
case Fieldprops::ScalarOperation::ADD:
data[index] += action_data.data[index];
break;
case Fieldprops::ScalarOperation::MAX:
data[index] = std::min(action_data.data[index], data[index]);
break;
case Fieldprops::ScalarOperation::MIN:
data[index] = std::max(action_data.data[index], data[index]);
break;
default:
throw std::logic_error("Unhandled value in switch");
}
}
}
}
void deserialize_tran(std::unordered_map<std::string, Fieldprops::TranCalculator>& tran, const std::vector<char>& buffer) {
tran.clear();
Serializer ser(buffer);
std::size_t size = ser.get<std::size_t>();
for (std::size_t calc_index = 0; calc_index < size; calc_index++) {
std::string calc_name = ser.get<std::string>();
Fieldprops::TranCalculator calc(calc_name);
std::size_t calc_size = ser.get<std::size_t>();
for (std::size_t action_index = 0; action_index < calc_size; action_index++) {
auto op = static_cast<Fieldprops::ScalarOperation>(ser.get<int>());
auto field = ser.get<std::string>();
calc.add_action(op, field);
}
tran.emplace(calc_name, std::move(calc));
}
}
template
void apply_tran(const std::unordered_map<std::string, Fieldprops::TranCalculator>&,
const std::unordered_map<std::string, Fieldprops::FieldData<double>>&,
std::size_t, const std::string&, std::vector<double>&);
template
void apply_tran(const std::unordered_map<std::string, Fieldprops::TranCalculator>&,
const std::map<std::string, Fieldprops::FieldData<double>>&,
std::size_t, const std::string&, std::vector<double>&);
template bool FieldPropsManager::supported<int>(const std::string&);
template bool FieldPropsManager::supported<double>(const std::string&);

View File

@@ -209,7 +209,7 @@ ADDREG
BOOST_AUTO_TEST_CASE(ASSIGN) {
FieldProps::FieldData<int> data({}, 100, 0);
Fieldprops::FieldData<int> data({}, 100, 0);
std::vector<int> wrong_size(50);
BOOST_CHECK_THROW( data.default_assign( wrong_size ), std::invalid_argument );
@@ -1940,3 +1940,96 @@ OPERATE
BOOST_CHECK_THROW(make_fp(invalid_region), std::logic_error);
BOOST_CHECK_THROW(make_fp(invalid_operate), std::logic_error);
}
BOOST_AUTO_TEST_CASE(TRAN_Calculator) {
std::string deck_string = R"(
GRID
PORO
1000*0.10 /
EDIT
TRANX
1000*0.10 /
BOX
1 10 1 10 1 1 /
TRANX
100*0.20 /
ENDBOX
MULTIPLY
TRANY 2.0 /
/
ADD
TRANY 3 1 10 1 10 2 2 /
/
MAXVALUE
TRANZ 0 1 10 1 10 1 1 /
/
MINVALUE
TRANZ 3 1 10 1 10 2 2 /
/
)";
UnitSystem unit_system(UnitSystem::UnitType::UNIT_TYPE_METRIC);
auto to_si = [&unit_system](double raw_value) { return unit_system.to_si(UnitSystem::measure::transmissibility, raw_value); };
std::vector<int> actnum(1000, 1);
for (std::size_t i=0; i< 1000; i += 2)
actnum[i] = 0;
EclipseGrid grid(EclipseGrid(10,10,10), actnum);
Deck deck = Parser{}.parseString(deck_string);
FieldPropsManager fpm(deck, Phases{true, true, true}, grid, TableManager());
std::vector<double> tranx( grid.getNumActive() );
std::vector<double> trany( grid.getNumActive(), to_si(1.0) );
std::vector<double> tranz( grid.getNumActive(), to_si(1.0) );
BOOST_CHECK(!fpm.has_double("TRANX"));
BOOST_CHECK_THROW(fpm.apply_tran("TRANA", tranx), std::out_of_range);
BOOST_CHECK(!fpm.tran_active("TRANA"));
fpm.apply_tran("TRANX", tranx);
fpm.apply_tran("TRANY", trany);
fpm.apply_tran("TRANZ", tranz);
for (std::size_t i=0; i < 50; i++)
BOOST_CHECK_EQUAL(tranx[i], to_si(0.20));
for (std::size_t i=50; i < tranx.size(); i++)
BOOST_CHECK_EQUAL(tranx[i], to_si(0.10));
BOOST_CHECK_EQUAL(trany[0], to_si(2.0));
for (std::size_t i=0; i < 50; i++)
BOOST_CHECK_EQUAL(trany[i], to_si(2.0));
for (std::size_t i=50; i < 100; i++)
BOOST_CHECK_EQUAL(trany[i], to_si(5.0));
for (std::size_t i=100; i < trany.size(); i++)
BOOST_CHECK_EQUAL(trany[i], to_si(2.0));
for (std::size_t i=0; i < 50; i++)
BOOST_CHECK(tranz[i]<=to_si(0));
for (std::size_t i=50; i < 10; i++)
BOOST_CHECK(tranz[i]>=to_si(3));
for (std::size_t i=100; i < tranz.size(); i++)
BOOST_CHECK_EQUAL(tranz[i], to_si(1.0));
auto buffer = fpm.serialize_tran();
fpm.deserialize_tran(buffer);
BOOST_CHECK(fpm.tran_active("TRANX"));
BOOST_CHECK(fpm.tran_active("TRANY"));
BOOST_CHECK(fpm.tran_active("TRANZ"));
}