From 967954fc4b7db12de2a2139e0834f86ebebea4d0 Mon Sep 17 00:00:00 2001 From: Jostein Alvestad Date: Fri, 23 Aug 2019 17:08:44 +0200 Subject: [PATCH] corrections and changes according to comments from maintainer --- opm/output/eclipse/AggregateUDQData.hpp | 2 +- src/opm/output/eclipse/AggregateUDQData.cpp | 46 +++++++++--------- src/opm/output/eclipse/CreateUdqDims.cpp | 54 ++++++++------------- 3 files changed, 43 insertions(+), 59 deletions(-) diff --git a/opm/output/eclipse/AggregateUDQData.hpp b/opm/output/eclipse/AggregateUDQData.hpp index b2987576f..c166ff129 100644 --- a/opm/output/eclipse/AggregateUDQData.hpp +++ b/opm/output/eclipse/AggregateUDQData.hpp @@ -128,7 +128,7 @@ private: /// Aggregate 'DUDW' array (Double Precision) for all UDQ data. (Dimension = max no wells * noOfUDQ's) WindowedArray dUDW_; - /// Aggregate 'DUDG' array (Double Precision) for all UDQ data. (Dimension = max no wells * noOfUDQ's) + /// Aggregate 'DUDG' array (Double Precision) for all UDQ data. (Dimension = (max no groups + 1) * noOfUDQ's) WindowedArray dUDG_; /// Aggregate 'DUDF' array (Double Precision) for all UDQ data. (Dimension = Number of FU - UDQ's, with value equal to the actual constraint) diff --git a/src/opm/output/eclipse/AggregateUDQData.cpp b/src/opm/output/eclipse/AggregateUDQData.cpp index 76dcac4e0..1f8382b61 100644 --- a/src/opm/output/eclipse/AggregateUDQData.cpp +++ b/src/opm/output/eclipse/AggregateUDQData.cpp @@ -247,18 +247,13 @@ namespace { DUDWArray& dUdw) { //initialize array to the default value for the array + for (std::size_t ind = 0; ind < nwmaxz; ind++) { + dUdw[ind] = -0.3E+21; + } for (std::size_t ind = 0; ind < wnames.size(); ind++) { if (st.has_well_var(wnames[ind], udq)) { dUdw[ind] = st.get_well_var(wnames[ind], udq); } - else { - dUdw[ind] = -0.3E+21; - } - } - if (wnames.size() < nwmaxz) { - for (std::size_t ind = wnames.size(); ind < nwmaxz; ind++) { - dUdw[ind] = -0.3E+21; - } } } } // dUdw @@ -354,18 +349,30 @@ std::pair findInVector(const std::vector & vecOfElements, const return result; } -const std::vector Opm::RestartIO::Helpers::igphData::ig_phase(const Opm::Schedule& sched, - const std::size_t simStep, - const std::vector& inteHead - ) +// Make ordered list of current groups +const std::vector currentGroups(const Opm::Schedule& sched, + const std::size_t simStep, + const std::vector& inteHead ) { std::vector curGroups(ngmaxz(inteHead), nullptr); for (const auto& group_name : sched.groupNames(simStep)) { const auto& group = sched.getGroup2(group_name, simStep); + + //The FIELD group is the first group according to the insert_index() + //In the Eclipse compatible restart file, the FILED group is put at the end of the list of groups (ngmaxz(inteHead)-1) int ind = (group.name() == "FIELD") ? ngmaxz(inteHead)-1 : group.insert_index()-1; curGroups[ind] = std::addressof(group); - } + + } + return curGroups; +} + +const std::vector Opm::RestartIO::Helpers::igphData::ig_phase(const Opm::Schedule& sched, + const std::size_t simStep, + const std::vector& inteHead ) +{ + const auto curGroups = currentGroups(sched, simStep, inteHead); std::vector inj_phase(ngmaxz(inteHead), 0); for (std::size_t ind = 0; ind < curGroups.size(); ind++) { if (curGroups[ind] != nullptr) { @@ -409,7 +416,7 @@ const std::vector iuap_data(const Opm::Schedule& sched, return wg_no; } - + Opm::RestartIO::Helpers::AggregateUDQData:: AggregateUDQData(const std::vector& udqDims) : iUDQ_ (iUdq::allocate(udqDims)), @@ -485,16 +492,9 @@ captureDeclaredUDQData(const Opm::Schedule& sched, i_wudq++; } } - - // Make ordered list of current groups - std::vector curGroups(ngmaxz(inteHead), nullptr); - for (const auto& group_name : sched.groupNames(simStep)) { - const auto& group = sched.getGroup2(group_name, simStep); - int ind = (group.name() == "FIELD") - ? ngmaxz(inteHead)-1 : group.insert_index()-1; - curGroups[ind] = std::addressof(group); - } + std::size_t i_gudq = 0; + const auto curGroups = currentGroups(sched, simStep, inteHead); const auto ngmax = ngmaxz(inteHead); for (const auto& udq_input : udqCfg.input()) { if (udq_input.var_type() == UDQVarType::GROUP_VAR) { diff --git a/src/opm/output/eclipse/CreateUdqDims.cpp b/src/opm/output/eclipse/CreateUdqDims.cpp index 966582c3a..4b46fa566 100755 --- a/src/opm/output/eclipse/CreateUdqDims.cpp +++ b/src/opm/output/eclipse/CreateUdqDims.cpp @@ -101,15 +101,10 @@ int noGroupUdqs(const Opm::Schedule& sched, const std::size_t simStep) { const auto& udqCfg = sched.getUDQConfig(simStep); - std::size_t i_gudq = 0; - for (const auto& udq_input : udqCfg.input()) { - if (udq_input.var_type() == Opm::UDQVarType::GROUP_VAR) { - i_gudq++; - } - } - return i_gudq; -} + const auto& input = udqCfg.input(); + return std::count_if(input.begin(), input.end(), [](const Opm::UDQInput inp) { return (inp.var_type() == Opm::UDQVarType::GROUP_VAR); }); +} int noFieldUdqs(const Opm::Schedule& sched, const std::size_t simStep) @@ -120,8 +115,8 @@ int noFieldUdqs(const Opm::Schedule& sched, if (udq_input.var_type() == Opm::UDQVarType::FIELD_VAR) { i_fudq++; } - return i_fudq; } + return i_fudq; } } // Anonymous @@ -139,32 +134,21 @@ createUdqDims(const Schedule& sched, const auto& udqCfg = sched.getUDQConfig(lookup_step); const auto& udqActive = sched.udqActive(lookup_step); std::vector udqDims; - // 0 - udqDims.emplace_back(udqCfg.size()); - // 1 - udqDims.emplace_back(entriesPerIUDQ()); - // 2 - udqDims.emplace_back(udqActive.IUAD_size()); - // 3 - udqDims.emplace_back(entriesPerIUAD()); - // 4 - udqDims.emplace_back(entriesPerZUDN()); - // 5 - udqDims.emplace_back(entriesPerZUDL()); - // 6 - udqDims.emplace_back(noIGphs(inteHead)); - // 7 - udqDims.emplace_back(udqActive.IUAP_size()); - // 8 - udqDims.emplace_back(nwmaxz(inteHead)); - // 9 - udqDims.emplace_back(noWellUdqs(sched, lookup_step)); - // 10 - udqDims.emplace_back(ngmaxz(inteHead)); - // 11 - udqDims.emplace_back(noGroupUdqs(sched, lookup_step)); - // 12 - udqDims.emplace_back(noFieldUdqs(sched, lookup_step)); + udqDims.resize(13,0); + + udqDims[ 0] = udqCfg.size(); + udqDims[ 1] = entriesPerIUDQ(); + udqDims[ 2] = udqActive.IUAD_size(); + udqDims[ 3] = entriesPerIUAD(); + udqDims[ 4] = entriesPerZUDN(); + udqDims[ 5] = entriesPerZUDL(); + udqDims[ 6] = noIGphs(inteHead); + udqDims[ 7] = udqActive.IUAP_size(); + udqDims[ 8] = nwmaxz(inteHead); + udqDims[ 9] = noWellUdqs(sched, lookup_step); + udqDims[10] = ngmaxz(inteHead); + udqDims[11] = noGroupUdqs(sched, lookup_step); + udqDims[12] = noFieldUdqs(sched, lookup_step); return udqDims; }