21#include <dune/grid/common/rangegenerators.hh>
22#include <dune/grid/io/file/vtk/common.hh>
24#include <opm/common/ErrorMacros.hpp>
31namespace Opm::GridDataOutput {
33template <
class Gr
idView,
unsigned int partitions>
39 dimw_ = GridView::dimension;
44template <
class Gr
idView,
unsigned int partitions>
47 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
49 if (Dune::VTK::geometryType(
corner_geom.type()) == Dune::VTK::polyhedron) {
56template <
class Gr
idView,
unsigned int partitions>
73template <
class Gr
idView,
unsigned int partitions>
78 if (max_size < nvertices_) {
80 "Opm::GridDataOutput::writeGridPoints( T& x_inout, T& "
81 "y_inout, T& z_inout ) "
82 +
" Input objects max_size (" + std::to_string(max_size)
83 +
") is not sufficient to fit the nvertices_ values (" + std::to_string(nvertices_) +
")");
88 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
95 }
else if (dimw_ == 2) {
96 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
100 z_inout[i] =
static_cast<T
>(0.0);
109template <
class Gr
idView,
unsigned int partitions>
110template <
typename VectType>
120 if ((
check_size_x <
static_cast<std::size_t
>(nvertices_)) ||
122 (
check_size_z <
static_cast<std::size_t
>(nvertices_))) {
125 "Opm::GridDataOutput::writeGridPoints( VectType& x_inout, VectType& "
126 "y_inout, VectType& z_inout ) At least one of the inputs"
127 +
" object x size " + std::to_string(
check_size_x) +
" object y size "
129 +
" is not sufficient to fit the nvertices_ values( " + std::to_string(nvertices_) +
" )");
134 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
141 }
else if (dimw_ == 2) {
143 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
156template <
class Gr
idView,
unsigned int partitions>
161 if (max_size < nvertices_ * 3) {
162 assert(max_size >= nvertices_ * 3);
164 "Opm::GridDataOutput::writeGridPoints_AOS( T* xyz_inout ) " +
" Input objects max_size ("
165 + std::to_string(max_size) +
") is not sufficient to fit the nvertices_ * 3 values ("
166 + std::to_string(nvertices_ * 3) +
")");
171 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
177 }
else if (dimw_ == 2) {
178 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
188template <
class Gr
idView,
unsigned int partitions>
189template <
typename VectType>
197 if (
check_size <
static_cast<std::size_t
>(nvertices_ * 3)) {
200 "Opm::GridDataOutput::writeGridPoints_AOS( VectType& xyz_inout ) "
201 +
" Input objects check_size (" + std::to_string(
check_size)
202 +
") is not sufficient to fit the nvertices_ * 3 values (" + std::to_string(nvertices_ * 3)
208 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
214 }
else if (dimw_ == 2) {
216 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
226template <
class Gr
idView,
unsigned int partitions>
231 if (max_size < nvertices_ * 3) {
234 "Opm::GridDataOutput::writeGridPoints_SOA( T& xyz_inout ) " +
" Input objects max_size ("
235 + std::to_string(max_size) +
") is not sufficient to fit the nvertices_ * 3 values ("
236 + std::to_string(nvertices_ * 3) +
")");
245 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
252 }
else if (dimw_ == 2) {
253 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
264template <
class Gr
idView,
unsigned int partitions>
265template <
typename VectType>
271 if (
check_size <
static_cast<std::size_t
>(nvertices_ * 3)) {
274 "Opm::GridDataOutput::writeGridPoints_SOA( VectType& xyz_inout ) "
275 +
" Input objects check_size (" + std::to_string(
check_size)
276 +
") is not sufficient to fit the nvertices_ * 3 values (" + std::to_string(nvertices_ * 3)
288 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
295 }
else if (dimw_ == 2) {
297 for (
const auto&
vit :
vertices(gridView_, Dune::Partitions::all)) {
308template <
class Gr
idView,
unsigned int partitions>
309template <
typename Integer>
314 if (max_size < ncorners_) {
317 "Opm::GridDataOutput::writeConnectivity( T* connectivity_inout,... ) "
318 +
" Input max_size value (" + std::to_string(max_size)
319 +
") is not sufficient to fit the ncorners_ values (" + std::to_string(ncorners_) +
")");
325 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
328 const int vxIdx = gridView_.indexSet().subIndex(
cit,
vx, 3);
335 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
338 const int vxIdx = gridView_.indexSet().subIndex(
cit,
vx, 3);
349template <
class Gr
idView,
unsigned int partitions>
350template <
typename VectType>
357 if (
check_size <
static_cast<std::size_t
>(ncorners_)) {
360 "Opm::GridDataOutput::writeConnectivity( VectType& "
361 "connectivity_inout ) "
362 +
" Input objects size (" + std::to_string(
check_size)
363 +
") is not sufficient to fit the ncorners_ values (" + std::to_string(ncorners_) +
")");
369 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
372 const int vxIdx = gridView_.indexSet().subIndex(
cit,
vx, 3);
379 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
382 const int vxIdx = gridView_.indexSet().subIndex(
cit,
vx, 3);
393template <
class Gr
idView,
unsigned int partitions>
394template <
typename Integer>
398 if (max_size < ncells_) {
401 "Opm::GridDataOutput::writeOffsetsCells( T* offsets_inout ) " +
" Input objects max_size ("
402 + std::to_string(max_size) +
") is not sufficient to fit the ncells_ values ("
403 + std::to_string(ncells_) +
")");
407 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
415template <
class Gr
idView,
unsigned int partitions>
416template <
typename VectType>
421 if (
check_size <
static_cast<std::size_t
>(ncells_)) {
424 "Opm::GridDataOutput::writeOffsetsCells( VectType& "
426 +
" Input objects size (" + std::to_string(
offsets_inout.size())
427 +
") is not sufficient to fit the ncells_ values (" + std::to_string(ncells_) +
")");
434 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
442template <
class Gr
idView,
unsigned int partitions>
443template <
typename Integer>
447 if (max_size < ncells_) {
450 "Opm::GridDataOutput::writeCellTypes( T* types_inout ) " +
" Input objects max_size ("
451 + std::to_string(max_size) +
") is not sufficient to fit the ncells_ values ("
452 + std::to_string(ncells_) +
")");
455 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
462template <
class Gr
idView,
unsigned int partitions>
463template <
typename VectType>
469 if (
check_size <
static_cast<std::size_t
>(ncells_)) {
471 "Opm::GridDataOutput::writeCellTypes( VectType& types_inout ) " +
" Input objects check_size ("
472 + std::to_string(
check_size) +
") is not sufficient to fit the ncells_ values ("
473 + std::to_string(ncells_) +
")");
477 for (
const auto&
cit :
elements(gridView_, dunePartition_)) {
478 int vtktype =
static_cast<int>(Dune::VTK::geometryType(
cit.type()));
484template <
class Gr
idView,
unsigned int partitions>
488 if (this->dunePartition_ == Dune::Partitions::all) {
489 return "Dune::Partitions::all";
491 if (this->dunePartition_ == Dune::Partitions::interior) {
492 return "Dune::Partitions::interior";
494 if (this->dunePartition_ == Dune::Partitions::interiorBorder) {
495 return "Dune::Partitions::interiorBorder";
497 if (this->dunePartition_ == Dune::Partitions::interiorBorderOverlap) {
498 return "Dune::Partitions::interiorBorderOverlap";
500 if (this->dunePartition_ == Dune::Partitions::front) {
501 return "Dune::Partitions::front";
503 if (this->dunePartition_ == Dune::Partitions::interiorBorderOverlapFront) {
504 return "Dune::Partitions::InteriorBorderOverlapFront";
506 if (this->dunePartition_ == Dune::Partitions::border) {
507 return "Dune::Partitions::border";
509 if (this->dunePartition_ == Dune::Partitions::ghost) {
510 return "Dune::Partitions::ghost";
513 return "Unknown Dune::PartitionSet<>";
516template <
class Gr
idView,
unsigned int partitions>
517void SimMeshDataAccessor<GridView,partitions>::
518printGridDetails(std::ostream&
outstr)
const
520 outstr <<
"Dune Partition = " << partition_value_ <<
", " << getPartitionTypeString() << std::endl;
521 outstr <<
"ncells_: " << getNCells() << std::endl;
522 outstr <<
"nvertices_: " << getNVertices() << std::endl;
523 outstr <<
"ncorners_: " << getNCorners() << std::endl;
Allows model geometry data to be passed to external code - via a copy direct to input pointers.
ConnectivityVertexOrder
Allows selection of order of vertices in writeConnectivity()
Definition GridDataOutput.hpp:93
Definition GridDataOutput.hpp:97
long writeGridPoints_AOS(T *xyz_inout, long max_size=0) const
Write the positions of vertices - directly to the pointers given in parameters as Array of Structures...
Definition GridDataOutput_impl.hpp:159
SimMeshDataAccessor(const GridView &gridView, Dune::PartitionSet< partitions > dunePartition)
Construct a SimMeshDataAccessor working on a specific GridView and specialize to a Dune::PartitionSet...
Definition GridDataOutput_impl.hpp:35
void countEntities()
Count the vertices, cells and corners.
Definition GridDataOutput_impl.hpp:57
long writeGridPoints(T *x_inout, T *y_inout, T *z_inout, long max_size=0) const
Write the positions of vertices - directly to the pointers given in parameters.
Definition GridDataOutput_impl.hpp:76
long writeOffsetsCells(Integer *offsets_inout, long max_size=0) const
Write the offsets values - directly to the pointer given in parameter 1.
Definition GridDataOutput_impl.hpp:396
long writeConnectivity(Integer *connectivity_inout, ConnectivityVertexOrder whichOrder, long max_size=0) const
Write the connectivity array - directly to the pointer given in parameter 1 Reorders the indices as s...
Definition GridDataOutput_impl.hpp:311
bool polyhedralCellPresent() const
Checks for cells that have polyhedral type within the current partition of cells.
Definition GridDataOutput_impl.hpp:45
long writeCellTypes(Integer *types_inout, long max_size=0) const
Write the cell types values - directly to the pointer given in parameter 1.
Definition GridDataOutput_impl.hpp:445
long writeGridPoints_SOA(T *xyz_inout, long max_size=0) const
Write the positions of vertices - directly to the pointers given in parameters as Structure of Arrays...
Definition GridDataOutput_impl.hpp:229
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242