31#ifndef OPM_DAMARIS_WRITER_HPP
32#define OPM_DAMARIS_WRITER_HPP
34#include <dune/grid/common/partitionset.hh>
37#include <opm/common/OpmLog/OpmLog.hpp>
49#include <fmt/format.h>
59namespace DamarisOutput {
64int write(
const char* field,
int rank,
const void* data);
66 const int n_elements_local_grid,
67 std::vector<unsigned long long>& elements_rank_offsets);
82template <
class TypeTag>
84 GetPropType<TypeTag, Properties::EquilGrid>,
85 GetPropType<TypeTag, Properties::GridView>,
86 GetPropType<TypeTag, Properties::ElementMapper>,
87 GetPropType<TypeTag, Properties::Scalar>>
89 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
90 using GridView = GetPropType<TypeTag, Properties::GridView>;
91 using Grid = GetPropType<TypeTag, Properties::Grid>;
92 using EquilGrid = GetPropType<TypeTag, Properties::EquilGrid>;
93 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
94 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
95 using Element =
typename GridView::template Codim<0>::Entity;
96 using ElementMapper = GetPropType<TypeTag, Properties::ElementMapper>;
106 Parameters::registerParam<TypeTag, Properties::DamarisOutputHdfCollective>
107 (
"Write output via Damaris using parallel HDF5 to "
108 "get single file and dataset per timestep instead "
109 "of one per Damaris core with multiple datasets.");
110 Parameters::registerParam<TypeTag, Properties::DamarisSaveToHdf>
111 (
"Set to false to prevent output to HDF5. "
112 "Uses collective output by default or "
113 "set --enable-damaris-collective=false to"
114 "use file per core (file per Damaris server).");
115 Parameters::registerParam<TypeTag, Properties::DamarisSaveMeshToHdf>
116 (
"Saves the mesh data to the HDF5 file (1st iteration only). "
117 "Will set --damaris-output-hdf-collective to false "
118 "so will use file per core (file per Damaris server) output "
119 "(global sizes and offset values of mesh variables are not being provided as yet).");
120 Parameters::registerParam<TypeTag, Properties::DamarisPythonScript>
121 (
"Set to the path and filename of a Python script to run on "
122 "Damaris server resources with access to OPM flow data.");
123 Parameters::registerParam<TypeTag, Properties::DamarisPythonParaviewScript>
124 (
"Set to the path and filename of a Paraview Python script "
125 "to run on Paraview Catalyst (1 or 2) on Damaris server "
126 "resources with access to OPM flow data.");
127 Parameters::registerParam<TypeTag, Properties::DamarisSimName>
128 (
"The name of the simulation to be used by Damaris. "
129 "If empty (the default) then Damaris uses \"opm-sim-<random-number>\". "
130 "This name is used for the Damaris HDF5 file name prefix. "
131 "Make unique if writing to the same output directory.");
132 Parameters::registerParam<TypeTag, Properties::DamarisLogLevel>
133 (
"The log level for the Damaris logging system (boost log based). "
134 "Levels are: [trace, debug, info, warning, error, fatal]. "
135 "Currently debug and info are useful. ");
136 Parameters::registerParam<TypeTag, Properties::DamarisDaskFile>
137 (
"The name of a Dask json configuration file (if using Dask for processing).");
138 Parameters::registerParam<TypeTag, Properties::DamarisDedicatedCores>
139 (
"Set the number of dedicated cores (MPI processes) "
140 "that should be used for Damaris processing (per node). "
141 "Must divide evenly into the number of simulation ranks (client ranks).");
142 Parameters::registerParam<TypeTag, Properties::DamarisDedicatedNodes>
143 (
"Set the number of dedicated nodes (full nodes) "
144 "that should be used for Damaris processing (per simulation). "
145 "Must divide evenly into the number of simulation nodes.");
146 Parameters::registerParam<TypeTag, Properties::DamarisSharedMemorySizeBytes>
147 (
"Set the size of the shared memory buffer used for IPC "
148 "between the simulation and the Damaris resources. "
149 "Needs to hold all the variables published, possibly over "
150 "multiple simulation iterations.");
151 Parameters::registerParam<TypeTag, Properties::DamarisSharedMemoryName>
152 (
"The name of the shared memory area to be used by Damaris for the current. "
153 "If empty (the default) then Damaris uses \"opm-damaris-<random-string>\". "
154 "This name should be unique if multiple simulations are running on "
155 "the same node/server as it is used for the Damaris shmem name and by "
156 "the Python Dask library to locate sections of variables.");
163 :
BaseType(simulator.vanguard().schedule(),
164 simulator.vanguard().eclState(),
165 simulator.vanguard().summaryConfig(),
166 simulator.vanguard().grid(),
167 ((simulator.vanguard().grid().comm().rank() == 0)
168 ? &simulator.vanguard().equilGrid()
170 simulator.vanguard().gridView(),
171 simulator.vanguard().cartesianIndexMapper(),
172 ((simulator.vanguard().grid().comm().rank() == 0)
173 ? &simulator.vanguard().equilCartesianIndexMapper()
176 , simulator_(simulator)
178 this->damarisUpdate_ = true ;
180 this->rank_ = this->simulator_.vanguard().grid().comm().rank() ;
181 this->nranks_ = this->simulator_.vanguard().grid().comm().size();
183 this->elements_rank_offsets_.resize(this->nranks_);
193 const auto& gridView = this->simulator_.gridView();
194 const auto& interior_elements = elements(gridView, Dune::Partitions::interior);
196 this->numElements_ = std::distance(interior_elements.begin(), interior_elements.end());
199 if (this->nranks_ > 1) {
200 auto smryCfg = (this->rank_ == 0)
201 ? this->
eclIO_->finalSummaryConfig()
204 eclBroadcast(this->simulator_.vanguard().grid().comm(), smryCfg);
206 this->damarisOutputModule_ = std::make_unique<OutputBlackOilModule<TypeTag>>
210 this->damarisOutputModule_ = std::make_unique<OutputBlackOilModule<TypeTag>>
221 const int reportStepNum = simulator_.episodeIndex() + 1;
225 this->damarisOutputModule_->invalidateLocalData() ;
226 this->prepareLocalCellData(isSubStep, reportStepNum);
227 this->damarisOutputModule_->outputErrorLog(simulator_.gridView().comm());
230 auto localWellData = simulator_.problem().wellModel().wellData();
234 if (localCellData.size() == 0) {
235 this->damarisOutputModule_->assignToSolution(localCellData);
239 this->damarisOutputModule_->addRftDataToWells(localWellData, reportStepNum);
242 if (damarisUpdate_ ==
true) {
248 numElements_, elements_rank_offsets_);
251 this->setGlobalIndexForDamaris() ;
255 this->writeDamarisGridOutput() ;
259 this->damarisUpdate_ =
false;
262 if (this->damarisOutputModule_->getPRESSURE_ptr() !=
nullptr)
265 this->elements_rank_offsets_[rank_]);
267 this->damarisOutputModule_->getPRESSURE_ptr());
280 Simulator& simulator_;
281 std::unique_ptr<OutputBlackOilModule<TypeTag>> damarisOutputModule_;
282 std::vector<unsigned long long> elements_rank_offsets_ ;
283 bool damarisUpdate_ =
false;
285 static bool enableDamarisOutput_()
287 return Parameters::get<TypeTag, Properties::EnableDamarisOutput>();
290 void setGlobalIndexForDamaris ()
295 const std::vector<int>& local_to_global =
299 std::vector<int> local_to_global_filled ;
300 local_to_global_filled.resize(this->numElements_) ;
301 std::iota(local_to_global_filled.begin(), local_to_global_filled.end(), 0);
310 DamarisVarInt mpi_rank_var_test(1, {
"n_elements_local"},
"MPI_RANK", rank_);
311 mpi_rank_var_test.setDamarisParameterAndShmem( {this->numElements_ } ) ;
313 std::fill(mpi_rank_var_test.data(), mpi_rank_var_test.data() + numElements_, rank_);
316 void writeDamarisGridOutput()
318 const auto& gridView = simulator_.gridView();
319 GridDataOutput::SimMeshDataAccessor geomData(gridView, Dune::Partitions::interior) ;
322 const bool hasPolyCells = geomData.polyhedralCellPresent() ;
323 if ( hasPolyCells ) {
324 OpmLog::error(fmt::format(
"ERORR: rank {} The DUNE geometry grid has polyhedral elements - These elements are currently not supported.", rank_ ));
338 DamarisVarDbl var_x(1, {
"n_coords_local"},
"coordset/coords/values/x", rank_) ;
342 var_x.setDamarisParameterAndShmem( { geomData.getNVertices() } ) ;
344 DamarisVarDbl var_y(1, {
"n_coords_local"},
"coordset/coords/values/y", rank_) ;
345 var_y.setDamarisParameterAndShmem( { geomData.getNVertices() } ) ;
347 DamarisVarDbl var_z(1, {
"n_coords_local"},
"coordset/coords/values/z", rank_) ;
348 var_z.setDamarisParameterAndShmem( { geomData.getNVertices() } ) ;
351 if ( geomData.writeGridPoints(var_x, var_y, var_z) < 0)
352 DUNE_THROW(Dune::IOError, geomData.getError() );
367 DamarisVarInt var_connectivity(1, {
"n_connectivity_ph"},
368 "topologies/topo/elements/connectivity", rank_) ;
369 var_connectivity.setDamarisParameterAndShmem({ geomData.getNCorners()}) ;
370 DamarisVarInt var_offsets(1, {
"n_offsets_types_ph"},
371 "topologies/topo/elements/offsets", rank_) ;
372 var_offsets.setDamarisParameterAndShmem({ geomData.getNCells()+1}) ;
373 DamarisVarChar var_types(1, {
"n_offsets_types_ph"},
374 "topologies/topo/elements/types", rank_) ;
375 var_types.setDamarisParameterAndShmem({ geomData.getNCells()}) ;
381 i = geomData.writeConnectivity(var_connectivity, vtkorder) ;
382 if ( i != geomData.getNCorners())
383 DUNE_THROW(Dune::IOError, geomData.getError());
385 i = geomData.writeOffsetsCells(var_offsets);
386 if ( i != geomData.getNCells()+1)
387 DUNE_THROW(Dune::IOError,geomData.getError());
389 i = geomData.writeCellTypes(var_types) ;
390 if ( i != geomData.getNCells())
391 DUNE_THROW(Dune::IOError,geomData.getError());
393 catch (std::exception& e)
395 OpmLog::error(e.what());
399 void prepareLocalCellData(
const bool isSubStep,
400 const int reportStepNum)
402 OPM_TIMEBLOCK(prepareLocalCellData);
403 if (damarisOutputModule_->localDataValid()) {
407 const auto& gridView = simulator_.vanguard().gridView();
412 damarisOutputModule_->allocBuffers(num_interior, reportStepNum,
413 isSubStep, log,
false);
415 ElementContext elemCtx(simulator_);
418 OPM_TIMEBLOCK(prepareCellBasedData);
419 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
420 elemCtx.updatePrimaryStencil(elem);
421 elemCtx.updatePrimaryIntensiveQuantities(0);
423 damarisOutputModule_->processElement(elemCtx);
426 if(!simulator_.model().linearizer().getFlowsInfo().empty()){
427 OPM_TIMEBLOCK(prepareFlowsData);
428 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
429 elemCtx.updatePrimaryStencil(elem);
430 elemCtx.updatePrimaryIntensiveQuantities(0);
431 damarisOutputModule_->processElementFlows(elemCtx);
435 OPM_TIMEBLOCK(prepareBlockData);
436 for (
const auto& elem : elements(gridView, Dune::Partitions::interior)) {
437 elemCtx.updatePrimaryStencil(elem);
438 elemCtx.updatePrimaryIntensiveQuantities(0);
439 damarisOutputModule_->processElementBlockData(elemCtx);
443 OPM_TIMEBLOCK(prepareFluidInPlace);
445#pragma omp parallel for
447 for (
int dofIdx=0; dofIdx < num_interior; ++dofIdx){
448 const auto& intQuants = *(simulator_.model().cachedIntensiveQuantities(dofIdx, 0));
449 const auto totVolume = simulator_.model().dofTotalVolume(dofIdx);
450 damarisOutputModule_->updateFluidInPlace(dofIdx, intQuants, totVolume);
453 damarisOutputModule_->validateLocalData();
#define OPM_END_PARALLEL_TRY_CATCH(prefix, comm)
Catch exception and throw in a parallel try-catch clause.
Definition: DeferredLoggingErrorHelpers.hpp:172
#define OPM_BEGIN_PARALLEL_TRY_CATCH()
Macro to setup the try of a parallel try-catch.
Definition: DeferredLoggingErrorHelpers.hpp:138
Allows model geometry data to be passed to external code - via a copy direct to input pointers.
bool isParallel() const
Definition: CollectDataOnIORank.hpp:128
bool isIORank() const
Definition: CollectDataOnIORank.hpp:125
const std::vector< int > & localIdxToGlobalIdxMapping() const
Definition: CollectDataOnIORank.hpp:133
Definition: DamarisVar.hpp:101
Collects necessary output values and pass them to Damaris server processes.
Definition: DamarisWriter.hpp:88
void writeOutput(data::Solution &localCellData, bool isSubStep)
Writes localCellData through to Damaris servers. Sets up the unstructured mesh which is passed to Dam...
Definition: DamarisWriter.hpp:218
static void registerParameters()
Definition: DamarisWriter.hpp:104
DamarisWriter(Simulator &simulator)
Definition: DamarisWriter.hpp:162
Definition: EclGenericWriter.hpp:65
CollectDataOnIORankType collectOnIORank_
Definition: EclGenericWriter.hpp:152
std::unique_ptr< EclipseIO > eclIO_
Definition: EclGenericWriter.hpp:157
int setParameter(const char *field, int rank, int value)
int setPosition(const char *field, int rank, int64_t pos)
int endIteration(int rank)
int write(const char *field, int rank, const void *data)
int setupWritingPars(Parallel::Communication comm, const int n_elements_local_grid, std::vector< unsigned long long > &elements_rank_offsets)
ConnectivityVertexOrder
Definition: GridDataOutput.hpp:93
@ VTK
Definition: GridDataOutput.hpp:93
Dune::Communication< MPIComm > Communication
Definition: ParallelCommunication.hpp:30
std::size_t countLocalInteriorCellsGridView(const GridView &gridView)
Get the number of local interior cells in a grid view.
Definition: countGlobalCells.hpp:47
Definition: BlackoilPhases.hpp:27
void eclBroadcast(Parallel::Communication, T &)