Refactor application of NNC / EDITNNC

- The edit manipulations from EDITNNC have already been applied to the NNC data
  from opm-common

- The NNC data structures are guaranteed to be ordered, both with cell1 <= cell2
  and the NNCs are in ascending order

- The NNC output to EGRID / INIT files is based on std::vector<NNCdata>
This commit is contained in:
Joakim Hove 2020-11-10 10:48:53 +01:00
parent 9e08e1324c
commit da2bc2affd
7 changed files with 43 additions and 227 deletions

View File

@ -23,7 +23,6 @@
# originally generated with the command:
# find opm -name '*.c*' -printf '\t%p\n' | sort
list (APPEND MAIN_SOURCE_FILES
ebos/nncsorter.cpp
opm/core/props/satfunc/RelpermDiagnostics.cpp
opm/simulators/timestepping/SimulatorReport.cpp
opm/simulators/flow/MissingFeatures.cpp
@ -83,7 +82,6 @@ list (APPEND TEST_SOURCE_FILES
tests/test_vfpproperties.cpp
tests/test_milu.cpp
tests/test_multmatrixtransposed.cpp
tests/test_nncsorter.cpp
tests/test_wellmodel.cpp
tests/test_deferredlogger.cpp
tests/test_timer.cpp

View File

@ -28,8 +28,6 @@
#ifndef EWOMS_ECL_TRANSMISSIBILITY_HH
#define EWOMS_ECL_TRANSMISSIBILITY_HH
#include <ebos/nncsorter.hpp>
#include <opm/models/utils/propertysystem.hh>
#include <opm/models/common/multiphasebaseproperties.hh>
@ -50,6 +48,7 @@
#include <dune/common/fvector.hh>
#include <dune/common/fmatrix.hh>
#include <fmt/format.h>
#include <array>
#include <vector>
#include <unordered_map>
@ -733,13 +732,11 @@ private:
// First scale NNCs with EDITNNC.
std::vector<Opm::NNCdata> unprocessedNnc;
std::vector<Opm::NNCdata> processedNnc;
const auto& nnc = vanguard_.eclState().getInputNNC();
if (!nnc.hasNNC())
const auto& nnc_input = vanguard_.eclState().getInputNNC().input();
if (nnc_input.empty())
return make_tuple(processedNnc, unprocessedNnc);
auto nncData = sortNncAndApplyEditnnc(nnc.data(), vanguard_.eclState().getInputEDITNNC().data());
for (const auto& nncEntry : nncData) {
for (const auto& nncEntry : nnc_input) {
auto c1 = nncEntry.cell1;
auto c2 = nncEntry.cell2;
auto low = cartesianToCompressed[c1];
@ -766,30 +763,50 @@ private:
if (candidate == trans_.end())
// This NNC is not resembled by the grid. Save it for later
// processing with local cell values
unprocessedNnc.push_back({c1, c2, nncEntry.trans});
unprocessedNnc.push_back(nncEntry);
else {
// NNC is represented by the grid and might be a neighboring connection
// In this case the transmissibilty is added to the value already
// set or computed.
candidate->second += nncEntry.trans;
processedNnc.push_back({c1, c2, nncEntry.trans});
processedNnc.push_back(nncEntry);
}
}
return make_tuple(processedNnc, unprocessedNnc);
}
/// \brief Multiplies the grid transmissibilities according to EDITNNC.
void applyEditNncToGridTrans_(const std::vector<int>& globalToLocal)
{
const auto& editNnc = vanguard_.eclState().getInputEDITNNC();
const auto& nnc_input = vanguard_.eclState().getInputNNC();
const auto& editNnc = nnc_input.edit();
if (editNnc.empty())
return;
const auto& cartMapper = vanguard_.cartesianIndexMapper();
const auto& cartDims = cartMapper.cartesianDimensions();
auto format_ijk = [&cartDims](std::size_t cell) -> std::string {
auto i = cell % cartDims[0]; cell /= cartDims[0];
auto j = cell % cartDims[1];
auto k = cell / cartDims[1];
return fmt::format("({},{},{})", i + 1,j + 1,k + 1);
};
auto make_warning = [&format_ijk] (const Opm::KeywordLocation& location, const Opm::NNCdata& nnc) -> std::string {
return fmt::format("Problem with EDITNNC keyword\n"
"In {} line {} \n"
"No NNC defined for connection {} -> {}", location.filename, location.lineno, format_ijk(nnc.cell1), format_ijk(nnc.cell2));
};
// editNnc is supposed to only reference non-neighboring connections and not
// neighboring connections. Use all entries for scaling if there is an NNC.
// variable nnc incremented in loop body.
auto nnc = editNnc.data().begin();
auto end = editNnc.data().end();
auto nnc = editNnc.begin();
auto end = editNnc.end();
std::size_t warning_count = 0;
while (nnc != end) {
auto c1 = nnc->cell1;
auto c2 = nnc->cell2;
@ -800,11 +817,11 @@ private:
auto candidate = trans_.find(isId_(low, high));
if (candidate == trans_.end()) {
std::ostringstream sstr;
sstr << "Cannot edit NNC from " << c1 << " to " << c2
<< " as it does not exist";
Opm::OpmLog::warning(sstr.str());
const auto& location = nnc_input.edit_location( *nnc );
auto warning = make_warning(location, *nnc);
Opm::OpmLog::warning("EDITNNC", warning);
++nnc;
warning_count++;
}
else {
// NNC exists
@ -814,6 +831,12 @@ private:
}
}
}
if (warning_count > 0) {
auto warning = fmt::format("Problems with EDITNNC keyword\n"
"A total of {} connections not defined in grid", warning_count);
Opm::OpmLog::warning(warning);
}
}
void extractPermeability_()

View File

@ -39,8 +39,6 @@
#include <opm/models/io/baseoutputwriter.hh>
#include <opm/models/parallel/tasklets.hh>
#include <ebos/nncsorter.hpp>
#include <opm/output/eclipse/EclipseIO.hpp>
#include <opm/output/eclipse/RestartValue.hpp>
@ -549,12 +547,11 @@ private:
{"TRANZ", tranz}};
}
Opm::NNC exportNncStructure_(const std::unordered_map<int,int>& cartesianToActive) const
std::vector<NNCdata> exportNncStructure_(const std::unordered_map<int,int>& cartesianToActive) const
{
std::size_t nx = eclState().getInputGrid().getNX();
std::size_t ny = eclState().getInputGrid().getNY();
auto nncData = sortNncAndApplyEditnnc(eclState().getInputNNC().data(),
eclState().getInputEDITNNC().data());
auto nncData = eclState().getInputNNC().input();
const auto& unitSystem = simulator_.vanguard().eclState().getDeckUnitSystem();
std::vector<Opm::NNCdata> outputNnc;
std::size_t index = 0;
@ -574,13 +571,6 @@ private:
++index;
}
auto nncCompare = []( const Opm::NNCdata& nnc1, const Opm::NNCdata& nnc2){
return nnc1.cell1 < nnc2.cell1 ||
( nnc1.cell1 == nnc2.cell1 && nnc1.cell2 < nnc2.cell2);};
// Sort the nncData values from the deck as they need to be
// Checked when writing NNC transmissibilities from the simulation.
std::sort(nncData.begin(), nncData.end(), nncCompare);
typedef typename EquilGrid :: LeafGridView GlobalGridView;
const GlobalGridView& globalGridView = globalGrid().leafGridView();
typedef Dune::MultipleCodimMultipleGeomTypeMapper<GlobalGridView> ElementMapper;
@ -603,7 +593,6 @@ private:
// Cartesian index mapper for the serial I/O grid
const auto& equilCartMapper = simulator_.vanguard().equilCartesianIndexMapper();
const auto& cartDims = simulator_.vanguard().cartesianIndexMapper().cartesianDimensions();
auto elemIt = globalGridView.template begin</*codim=*/0>();
const auto& elemEndIt = globalGridView.template end</*codim=*/0>();
@ -639,7 +628,7 @@ private:
// We need to check whether an NNC for this face was also specified
// via the NNC keyword in the deck (i.e. in the first origNncSize entries.
auto t = globalTrans->transmissibility(c1, c2);
auto candidate = std::lower_bound(nncData.begin(), nncData.end(), Opm::NNCdata(cc1, cc2, 0.0), nncCompare);
auto candidate = std::lower_bound(nncData.begin(), nncData.end(), Opm::NNCdata(cc1, cc2, 0.0));
while ( candidate != nncData.end() && candidate->cell1 == cc1
&& candidate->cell2 == cc2) {
@ -655,10 +644,7 @@ private:
}
}
}
Opm::NNC ret;
for(const auto& nncItem: outputNnc)
ret.addNNC(nncItem.cell1, nncItem.cell2, nncItem.trans);
return ret;
return outputNnc;
}
struct EclWriteTasklet

View File

@ -1,95 +0,0 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
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 2 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/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
#include <config.h>
#include <ebos/nncsorter.hpp>
#include <opm/common/OpmLog/OpmLog.hpp>
#include <sstream>
#include <algorithm>
#include <iostream>
namespace Opm
{
std::vector<Opm::NNCdata> sortNncAndApplyEditnnc(const std::vector<Opm::NNCdata>& nncDataIn, std::vector<Opm::NNCdata> editnncData,
bool log )
{
auto nncLess =
[](const Opm::NNCdata& d1, const Opm::NNCdata& d2) {
return
(d1.cell1 < d2.cell1)
|| (d1.cell1 == d2.cell1 && d1.cell2 < d2.cell2);
};
auto makeCell1LessCell2 =
[](const Opm::NNCdata& entry) {
if ( entry.cell2 < entry.cell1)
return Opm::NNCdata(entry.cell2, entry.cell1, entry.trans);
else
return entry;
};
// We need to make sure that for each entry cell1<=cell2 holds. Otherwise sorting
// will not make the search more accurate if the engineer chooses to define NNCs
// differently.
std::vector<Opm::NNCdata> nncData(nncDataIn);
std::transform(nncData.begin(), nncData.end(), nncData.begin(), makeCell1LessCell2);
std::transform(editnncData.begin(), editnncData.end(), editnncData.begin(), makeCell1LessCell2);
std::sort(nncData.begin(), nncData.end(), nncLess);
auto candidate = nncData.begin();
for (const auto& edit: editnncData) {
auto printNncWarning =
[](int c1, int c2) {
std::ostringstream sstr;
sstr << "Cannot edit NNC from " << c1 << " to " << c2
<< " as it does not exist";
Opm::OpmLog::warning(sstr.str());
};
if (candidate == nncData.end() && log) {
// no more NNCs left
printNncWarning(edit.cell1, edit.cell2);
continue;
}
if (candidate->cell1 != edit.cell1 || candidate->cell2 != edit.cell2) {
candidate = std::lower_bound(nncData.begin(), nncData.end(), Opm::NNCdata(edit.cell1, edit.cell2, 0), nncLess);
if (candidate == nncData.end() && log) {
// no more NNCs left
printNncWarning(edit.cell1, edit.cell2);
continue;
}
}
auto firstCandidate = candidate;
while (candidate != nncData.end()
&& candidate->cell1 == edit.cell1
&& candidate->cell2 == edit.cell2)
{
candidate->trans *= edit.trans;
++candidate;
}
// start with first match in next iteration to catch case where next
// EDITNNC is for same pair.
candidate = firstCandidate;
}
return nncData;
}
} // end namespace Opm

View File

@ -1,42 +0,0 @@
// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
// vi: set et ts=4 sw=4 sts=4:
/*
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 2 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/>.
Consult the COPYING file in the top-level source directory of this
module for the precise wording of the license and the list of
copyright holders.
*/
#ifndef EWOMS_EBOS_NNCSORTER_HPP
#define EWOMS_EBOS_NNCSORTER_HPP
#include <opm/parser/eclipse/EclipseState/Grid/NNC.hpp>
#include <vector>
namespace Opm
{
/// \brief Scale NNC data wit informtion form EDITNNC and sort it.
/// \param nncData The NNC data as provided by the deck.
/// \param editnncData The EDITNNC data as provided by the deck.
/// \return A lexicographically sorted vector of the scaled NNC data.
/// For each entry entry.cell1<entry.cell2 will hold for convenience.
std::vector<Opm::NNCdata> sortNncAndApplyEditnnc(const std::vector<Opm::NNCdata>& nncData, std::vector<Opm::NNCdata> editnncData,
bool log = true);
}
#endif

View File

@ -37,7 +37,6 @@
#include <opm/parser/eclipse/EclipseState/EclipseConfig.hpp>
#include <opm/parser/eclipse/EclipseState/Runspec.hpp>
#include <opm/parser/eclipse/EclipseState/TracerConfig.hpp>
#include <opm/parser/eclipse/EclipseState/Edit/EDITNNC.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FaceDir.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/Fault.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/FaultCollection.hpp>
@ -483,7 +482,6 @@ TEST_FOR_TYPE(DenT)
TEST_FOR_TYPE(Dimension)
TEST_FOR_TYPE(EclHysterConfig)
TEST_FOR_TYPE(EclipseConfig)
TEST_FOR_TYPE(EDITNNC)
TEST_FOR_TYPE(EndpointScaling)
TEST_FOR_TYPE(Eqldims)
TEST_FOR_TYPE(Equil)

View File

@ -1,52 +0,0 @@
/*
Copyright 2019 Equinor.
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/>.
*/
#include "config.h"
#define BOOST_TEST_MODULE NNCSortTest
#include <boost/test/unit_test.hpp>
#include <ebos/nncsorter.hpp>
#include <opm/parser/eclipse/EclipseState/Grid/NNC.hpp>
#include <cmath>
BOOST_AUTO_TEST_CASE(Test1) {
std::vector<Opm::NNCdata> nncDataIn =
{ {9, 8, 10.0 }, { 1, 2, 3.0 }, { 3, 4, 2.0 }, { 2, 1, 5.0 } };
std::vector<Opm::NNCdata> editnncData =
{ {20, 5, .001}, { 2, 1, .1}, {3, 4, .001}, {0, 0, 0.0}, {4, 3, 2.0} };
std::vector<Opm::NNCdata> nncDataOut1 =
{ { 1, 2, 0.3 }, { 1, 2, 0.5 }, { 3, 4, 0.004 }, { 8, 9, 10.0 } };
std::vector<Opm::NNCdata> nncDataOut2 =
{ { 1, 2, 0.5 }, { 1, 2, 0.3 }, { 3, 4, 0.4 }, { 8, 9, 10.0 } };
auto nncDataProcessed = Opm::sortNncAndApplyEditnnc(nncDataIn, editnncData);
BOOST_CHECK(nncDataProcessed.size() == nncDataOut1.size());
auto expectedNnc1 = nncDataOut1.begin();
auto expectedNnc2 = nncDataOut2.begin();
for(const auto& entry: nncDataProcessed) {
BOOST_CHECK( entry.cell1 == expectedNnc1->cell1);
BOOST_CHECK( entry.cell2 == expectedNnc1->cell2);
BOOST_CHECK( std::abs(entry.trans - expectedNnc1->trans) / std::abs(entry.trans) < 1e-5 ||
std::abs(entry.trans - expectedNnc2->trans) / std::abs(entry.trans) < 1e-5);
++expectedNnc1;
++expectedNnc2;
}
}