My Project
Loading...
Searching...
No Matches
NlddReporting.hpp
1/*
2 Copyright 2025, SINTEF Digital
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef OPM_NLDD_REPORTING_HEADER_INCLUDED
21#define OPM_NLDD_REPORTING_HEADER_INCLUDED
22
23#include <dune/grid/common/gridenums.hh>
24#include <dune/grid/common/partitionset.hh>
25
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>
31
32#include <fmt/format.h>
33
34#include <algorithm>
35#include <cmath>
36#include <filesystem>
37#include <fstream>
38#include <limits>
39#include <numeric>
40#include <sstream>
41#include <vector>
42
43namespace Opm {
44
53void reportNlddStatistics(const std::vector<SimulatorReport>& domain_reports,
54 const SimulatorReport& local_report,
55 const bool output_cout,
56 const Parallel::Communication& comm);
57
68template <class Grid, class Domain, class ElementMapper, class CartMapper>
70 const std::filesystem::path& odir,
71 const std::vector<Domain>& domains,
72 const std::vector<SimulatorReport>& domain_reports,
73 const Grid& grid,
74 const ElementMapper& elementMapper,
76{
77 const auto& dims = cartMapper.cartesianDimensions();
78 const auto total_size = dims[0] * dims[1] * dims[2];
79 const auto& comm = grid.comm();
80 const int rank = comm.rank();
81
82 // Create a cell-to-iterations mapping for this process
83 std::vector<int> cell_iterations(grid.size(0), 0);
84
85 // Populate the mapping with iteration counts for each domain
86 const auto ds = domains.size();
87 for (auto domain_idx = 0*ds; domain_idx < ds; ++domain_idx) {
88 const auto& domain = domains[domain_idx];
89 const auto& report = domain_reports[domain_idx];
90 const int iterations = report.success.total_newton_iterations +
91 report.failure.total_newton_iterations;
92
93 for (const int cell_idx : domain.cells) {
94 cell_iterations[cell_idx] = iterations;
95 }
96 }
97
98 // Create a full-sized vector initialized with zeros (indicating inactive cells)
99 auto full_iterations = std::vector<int>(total_size, 0);
100
101 // Convert local cell indices to cartesian indices
102 const auto& gridView = grid.leafGridView();
103 for (const auto& cell : elements(gridView, Dune::Partitions::interior)) {
104 const int cell_idx = elementMapper.index(cell);
105 const int cart_idx = cartMapper.cartesianIndex(cell_idx);
107 }
108
109 // Gather all iteration data using max operation
110 comm.max(full_iterations.data(), full_iterations.size());
111
112 // Only rank 0 writes the file
113 if (rank == 0) {
114 auto fname = odir / "ResInsight_nonlinear_iterations.txt";
115 std::ofstream resInsightFile { fname };
116 // Write header
117 resInsightFile << "NLDD_ITER" << '\n';
118
119 // Write all cells, including inactive ones
120 for (const auto& val : full_iterations) {
121 resInsightFile << val << '\n';
122 }
123 resInsightFile << "/" << '\n';
124 }
125}
126
136template <class Grid, class Domain, class ElementMapper, class CartMapper>
138 const std::filesystem::path& odir,
139 const std::vector<Domain>& domains,
140 const Grid& grid,
141 const ElementMapper& elementMapper,
142 const CartMapper& cartMapper)
143{
144 const auto& dims = cartMapper.cartesianDimensions();
145 const auto total_size = dims[0] * dims[1] * dims[2];
146 const auto& comm = grid.comm();
147 const int rank = comm.rank();
148
150
151 // Create a full-sized vector initialized with -1 (indicating inactive cells)
152 auto full_partition = std::vector<int>(total_size, -1);
153
154 // Fill active cell values for this rank
155 auto i = 0;
156 for (const auto& cell : elements(grid.leafGridView(), Dune::Partitions::interior)) {
157 full_partition[cartMapper.cartesianIndex(elementMapper.index(cell))] = partition_vector[i++];
158 }
159
160 // Gather all partitions using max operation
161 comm.max(full_partition.data(), full_partition.size());
162
163 // Only rank 0 writes the file
164 if (rank == 0) {
165 auto fname = odir / "ResInsight_compatible_partition.txt";
166 std::ofstream resInsightFile { fname };
167
168 // Write header
169 resInsightFile << "NLDD_DOM" << '\n';
170
171 // Write all cells, including inactive ones
172 for (const auto& val : full_partition) {
173 resInsightFile << val << '\n';
174 }
175 resInsightFile << "/" << '\n';
176 }
177
178 const auto nDigit = 1 + static_cast<int>(std::floor(std::log10(comm.size())));
179 auto partition_fname = odir / fmt::format("{1:0>{0}}", nDigit, rank);
180 std::ofstream pfile { partition_fname };
181
182 auto cell_index = 0;
183 for (const auto& cell : elements(grid.leafGridView(), Dune::Partitions::interior)) {
184 pfile << rank << ' '
185 << cartMapper.cartesianIndex(elementMapper.index(cell)) << ' '
186 << partition_vector[cell_index++] << '\n';
187 }
188}
189
200template <class Grid, class Domain>
202 const std::vector<int>& partition_vector,
203 const std::vector<Domain>& domains,
205 std::vector<SimulatorReport>& domain_reports_accumulated,
206 const Grid& grid,
207 int num_wells)
208{
209 const auto& gridView = grid.leafGridView();
210 const auto& comm = grid.comm();
211 const int rank = comm.rank();
212
213 const int num_domains = domains.size();
214 const int owned_cells = partition_vector.size();
215
216 // Count overlap cells using grid view iteration
217 int overlap_cells = std::count_if(elements(gridView).begin(), elements(gridView).end(),
218 [](const auto& cell) { return cell.partitionType() == Dune::OverlapEntity; });
219
220 // Store data for summary output
221 local_reports_accumulated.success.num_wells = num_wells;
222 local_reports_accumulated.success.num_domains = num_domains;
223 local_reports_accumulated.success.num_overlap_cells = overlap_cells;
224 local_reports_accumulated.success.num_owned_cells = owned_cells;
225
226 // Set statistics for each domain report
227 const auto dr_size = domain_reports_accumulated.size();
228 for (auto i = 0*dr_size; i < dr_size; ++i) {
230 domain_report.success.num_domains = 1;
231 domain_report.success.num_overlap_cells = 0;
232 domain_report.success.num_owned_cells = domains[i].cells.size();
233 }
234
235 // Gather data from all ranks
236 std::vector<int> all_owned(comm.size());
237 std::vector<int> all_overlap(comm.size());
238 std::vector<int> all_wells(comm.size());
239 std::vector<int> all_domains(comm.size());
240
241 comm.gather(&owned_cells, all_owned.data(), 1, 0);
242 comm.gather(&overlap_cells, all_overlap.data(), 1, 0);
243 comm.gather(&num_wells, all_wells.data(), 1, 0);
244 comm.gather(&num_domains, all_domains.data(), 1, 0);
245
246 if (rank == 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";
251
252 int total_owned = 0;
253 int total_overlap = 0;
254 int total_wells = 0;
255 int total_domains = 0;
256
257 for (int r = 0; r < comm.size(); ++r) {
258 ss << std::setw(6) << r
259 << std::setw(13) << all_owned[r]
260 << std::setw(15) << all_overlap[r]
261 << std::setw(14) << (all_owned[r] + all_overlap[r])
262 << std::setw(8) << all_wells[r]
263 << std::setw(9) << all_domains[r] << '\n';
264
269 }
270
271 ss << "--------------------------------------------------------------------\n"
272 << " sum"
273 << std::setw(13) << total_owned
274 << std::setw(15) << total_overlap
275 << std::setw(14) << (total_owned + total_overlap)
276 << std::setw(8) << total_wells
277 << std::setw(9) << total_domains << '\n';
278
279 OpmLog::info(ss.str());
280 }
281}
282
291template <class Grid, class Domain>
293 const std::vector<Domain>& domains,
294 const Grid& grid)
295{
296 const auto& comm = grid.comm();
297 const int rank = comm.rank();
298
299 auto numD = std::vector<int>(comm.size() + 1, 0);
300 numD[rank + 1] = static_cast<int>(domains.size());
301 comm.sum(numD.data(), numD.size());
302 std::partial_sum(numD.begin(), numD.end(), numD.begin());
303
304 auto p = std::vector<int>(grid.size(0));
305 auto maxCellIdx = std::numeric_limits<int>::min();
306
307 auto d = numD[rank];
308 for (const auto& domain : domains) {
309 for (const auto& cell : domain.cells) {
310 p[cell] = d;
311 if (cell > maxCellIdx) {
312 maxCellIdx = cell;
313 }
314 }
315
316 ++d;
317 }
318
319 p.erase(p.begin() + maxCellIdx + 1, p.end());
320 return p;
321}
322
323} // namespace Opm
324
325#endif // OPM_NLDD_REPORTING_HEADER_INCLUDED
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