20#ifndef OPM_NLDD_REPORTING_HEADER_INCLUDED
21#define OPM_NLDD_REPORTING_HEADER_INCLUDED
23#include <dune/grid/common/gridenums.hh>
24#include <dune/grid/common/partitionset.hh>
26#include <opm/common/ErrorMacros.hpp>
27#include <opm/common/OpmLog/OpmLog.hpp>
28#include <opm/simulators/timestepping/SimulatorReport.hpp>
29#include <opm/simulators/utils/DeferredLogger.hpp>
30#include <opm/simulators/utils/ParallelCommunication.hpp>
32#include <fmt/format.h>
56 const Parallel::Communication& comm);
68template <
class Gr
id,
class Domain,
class ElementMapper,
class CartMapper>
70 const std::filesystem::path&
odir,
71 const std::vector<Domain>&
domains,
74 const ElementMapper& elementMapper,
79 const auto& comm = grid.comm();
80 const int rank = comm.rank();
90 const int iterations = report.success.total_newton_iterations +
91 report.failure.total_newton_iterations;
102 const auto& gridView = grid.leafGridView();
103 for (
const auto& cell :
elements(gridView, Dune::Partitions::interior)) {
104 const int cell_idx = elementMapper.index(cell);
114 auto fname =
odir /
"ResInsight_nonlinear_iterations.txt";
136template <
class Gr
id,
class Domain,
class ElementMapper,
class CartMapper>
138 const std::filesystem::path&
odir,
139 const std::vector<Domain>&
domains,
141 const ElementMapper& elementMapper,
146 const auto& comm = grid.comm();
147 const int rank = comm.rank();
156 for (
const auto& cell :
elements(grid.leafGridView(), Dune::Partitions::interior)) {
165 auto fname =
odir /
"ResInsight_compatible_partition.txt";
178 const auto nDigit = 1 +
static_cast<int>(std::floor(std::log10(comm.size())));
183 for (
const auto& cell :
elements(grid.leafGridView(), Dune::Partitions::interior)) {
185 <<
cartMapper.cartesianIndex(elementMapper.index(cell)) <<
' '
200template <
class Gr
id,
class Domain>
203 const std::vector<Domain>&
domains,
209 const auto& gridView = grid.leafGridView();
210 const auto& comm = grid.comm();
211 const int rank = comm.rank();
213 const int num_domains =
domains.size();
218 [](
const auto& cell) {
return cell.partitionType() == Dune::OverlapEntity; });
243 comm.gather(&num_wells,
all_wells.data(), 1, 0);
244 comm.gather(&num_domains,
all_domains.data(), 1, 0);
247 std::ostringstream
ss;
248 ss <<
"\nNLDD domain distribution summary:\n"
249 <<
" rank owned cells overlap cells total cells wells domains\n"
250 <<
"--------------------------------------------------------------------\n";
257 for (
int r = 0;
r < comm.size(); ++
r) {
258 ss << std::setw(6) <<
r
271 ss <<
"--------------------------------------------------------------------\n"
279 OpmLog::info(
ss.str());
291template <
class Gr
id,
class Domain>
293 const std::vector<Domain>&
domains,
296 const auto& comm = grid.comm();
297 const int rank = comm.rank();
299 auto numD = std::vector<int>(comm.size() + 1, 0);
302 std::partial_sum(
numD.begin(),
numD.end(),
numD.begin());
304 auto p = std::vector<int>(grid.size(0));
305 auto maxCellIdx = std::numeric_limits<int>::min();
309 for (
const auto& cell :
domain.cells) {
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
void printDomainDistributionSummary(const std::vector< int > &partition_vector, const std::vector< Domain > &domains, SimulatorReport &local_reports_accumulated, std::vector< SimulatorReport > &domain_reports_accumulated, const Grid &grid, int num_wells)
Prints a summary of domain distribution across ranks.
Definition NlddReporting.hpp:201
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
void reportNlddStatistics(const std::vector< SimulatorReport > &domain_reports, const SimulatorReport &local_report, const bool output_cout, const Parallel::Communication &comm)
Reports NLDD statistics after simulation.
Definition NlddReporting.cpp:44
std::vector< int > reconstitutePartitionVector(const std::vector< Domain > &domains, const Grid &grid)
Reconstructs the partition vector that maps each grid cell to its corresponding domain ID,...
Definition NlddReporting.hpp:292
void writePartitions(const std::filesystem::path &odir, const std::vector< Domain > &domains, const Grid &grid, const ElementMapper &elementMapper, const CartMapper &cartMapper)
Writes the partition vector to a file in ResInsight compatible format and a partition file for each r...
Definition NlddReporting.hpp:137
void writeNonlinearIterationsPerCell(const std::filesystem::path &odir, const std::vector< Domain > &domains, const std::vector< SimulatorReport > &domain_reports, const Grid &grid, const ElementMapper &elementMapper, const CartMapper &cartMapper)
Writes the number of nonlinear iterations per cell to a file in ResInsight compatible format.
Definition NlddReporting.hpp:69
Definition SimulatorReport.hpp:119