21#include <dune/grid/common/rangegenerators.hh>
22#include <dune/grid/io/file/vtk/common.hh>
24#include <opm/common/ErrorMacros.hpp>
33template <
class Gr
idView,
unsigned int partitions>
37 , dunePartition_(dunePartition)
39 dimw_ = GridView::dimension;
44template <
class Gr
idView,
unsigned int partitions>
47 for (
const auto& cit : elements(gridView_, dunePartition_)) {
48 auto corner_geom = cit.geometry();
49 if (Dune::VTK::geometryType(corner_geom.type()) == Dune::VTK::polyhedron) {
56template <
class Gr
idView,
unsigned int partitions>
60 const auto& vert_partition = vertices(gridView_, Dune::Partitions::all);
61 nvertices_ = std::distance(vert_partition.begin(), vert_partition.end());
63 const auto& cell_partition = elements(gridView_, dunePartition_);
66 for (
const auto& cit : cell_partition) {
67 auto corner_geom = cit.geometry();
68 ncorners_ += corner_geom.corners();
73template <
class Gr
idView,
unsigned int partitions>
78 if (max_size < nvertices_) {
79 OPM_THROW(std::runtime_error,
80 "Opm::GridDataOutput::writeGridPoints( T& x_inout, T& "
81 "y_inout, T& z_inout ) "
83 +
") is not sufficient to fit the nvertices_ values (" +
std::to_string(nvertices_) +
")");
88 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
89 auto xyz_local = vit.geometry().corner(0);
90 x_inout[i] =
static_cast<T
>(xyz_local[0]);
91 y_inout[i] =
static_cast<T
>(xyz_local[1]);
92 z_inout[i] =
static_cast<T
>(xyz_local[2]);
95 }
else if (dimw_ == 2) {
96 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
97 auto xyz_local = vit.geometry().corner(0);
98 x_inout[i] =
static_cast<T
>(xyz_local[0]);
99 y_inout[i] =
static_cast<T
>(xyz_local[1]);
100 z_inout[i] =
static_cast<T
>(0.0);
104 assert(i == nvertices_);
109template <
class Gr
idView,
unsigned int partitions>
110template <
typename VectType>
112writeGridPoints(VectType& x_inout, VectType& y_inout, VectType& z_inout)
const
114 const std::size_t check_size_x = x_inout.size();
115 const std::size_t check_size_y = y_inout.size();
116 const std::size_t check_size_z = z_inout.size();
118 using VT =
decltype(x_inout.data()[0]);
120 if ((check_size_x <
static_cast<std::size_t
>(nvertices_)) ||
121 (check_size_y < static_cast<std::size_t>(nvertices_)) ||
122 (check_size_z <
static_cast<std::size_t
>(nvertices_))) {
124 OPM_THROW(std::runtime_error,
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)) {
135 auto xyz_local = vit.geometry().corner(0);
136 x_inout.data()[i] =
static_cast<VT
>(xyz_local[0]);
137 y_inout.data()[i] =
static_cast<VT
>(xyz_local[1]);
138 z_inout.data()[i] =
static_cast<VT
>(xyz_local[2]);
141 }
else if (dimw_ == 2) {
143 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
144 auto xyz_local = vit.geometry().corner(0);
145 x_inout.data()[i] =
static_cast<VT
>(xyz_local[0]);
146 y_inout.data()[i] =
static_cast<VT
>(xyz_local[1]);
147 z_inout.data()[i] =
static_cast<VT
>(td);
151 assert(i == nvertices_);
156template <
class Gr
idView,
unsigned int partitions>
161 if (max_size < nvertices_ * 3) {
162 assert(max_size >= nvertices_ * 3);
163 OPM_THROW(std::runtime_error,
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 ("
171 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
172 auto xyz_local = vit.geometry().corner(0);
173 xyz_inout[i++] =
static_cast<T
>(xyz_local[0]);
174 xyz_inout[i++] =
static_cast<T
>(xyz_local[1]);
175 xyz_inout[i++] =
static_cast<T
>(xyz_local[2]);
177 }
else if (dimw_ == 2) {
178 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
179 auto xyz_local = vit.geometry().corner(0);
180 xyz_inout[i++] =
static_cast<T
>(xyz_local[0]);
181 xyz_inout[i++] =
static_cast<T
>(xyz_local[1]);
182 xyz_inout[i++] =
static_cast<T
>(0.0);
188template <
class Gr
idView,
unsigned int partitions>
189template <
typename VectType>
193 const std::size_t check_size = xyz_inout.size();
195 using VT =
decltype(xyz_inout.data()[0]);
197 if (check_size <
static_cast<std::size_t
>(nvertices_ * 3)) {
198 assert(check_size >= nvertices_ * 3);
199 OPM_THROW(std::runtime_error,
200 "Opm::GridDataOutput::writeGridPoints_AOS( VectType& xyz_inout ) "
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)) {
209 auto xyz_local = vit.geometry().corner(0);
210 xyz_inout.data()[i++] =
static_cast<VT
>(xyz_local[0]);
211 xyz_inout[i++] =
static_cast<VT
>(xyz_local[1]);
212 xyz_inout[i++] =
static_cast<VT
>(xyz_local[2]);
214 }
else if (dimw_ == 2) {
216 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
217 auto xyz_local = vit.geometry().corner(0);
218 xyz_inout[i++] =
static_cast<VT
>(xyz_local[0]);
219 xyz_inout[i++] =
static_cast<VT
>(xyz_local[1]);
220 xyz_inout[i++] =
static_cast<VT
>(td);
226template <
class Gr
idView,
unsigned int partitions>
231 if (max_size < nvertices_ * 3) {
233 OPM_THROW(std::runtime_error,
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 ("
241 T* xyz_inout_y = xyz_inout + nvertices_;
242 T* xyz_inout_z = xyz_inout + (2 * nvertices_);
245 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
246 auto xyz_local = vit.geometry().corner(0);
247 xyz_inout[i] =
static_cast<T
>(xyz_local[0]);
248 xyz_inout_y[i] =
static_cast<T
>(xyz_local[1]);
249 xyz_inout_z[i] =
static_cast<T
>(xyz_local[2]);
252 }
else if (dimw_ == 2) {
253 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
254 auto xyz_local = vit.geometry().corner(0);
255 xyz_inout[i] =
static_cast<T
>(xyz_local[0]);
256 xyz_inout_y[i] =
static_cast<T
>(xyz_local[1]);
257 xyz_inout_z[i] =
static_cast<T
>(0.0);
264template <
class Gr
idView,
unsigned int partitions>
265template <
typename VectType>
269 const std::size_t check_size = xyz_inout.size();
271 if (check_size <
static_cast<std::size_t
>(nvertices_ * 3)) {
273 OPM_THROW(std::runtime_error,
274 "Opm::GridDataOutput::writeGridPoints_SOA( VectType& xyz_inout ) "
276 +
") is not sufficient to fit the nvertices_ * 3 values (" +
std::to_string(nvertices_ * 3)
280 using VT =
decltype(xyz_inout.data()[0]);
284 VT* xyz_inout_y = xyz_inout.data() + nvertices_;
285 VT* xyz_inout_z = xyz_inout.data() + (2 * nvertices_);
288 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
289 auto xyz_local = vit.geometry().corner(0);
290 xyz_inout.data()[i] =
static_cast<VT
>(xyz_local[0]);
291 xyz_inout_y[i] =
static_cast<VT
>(xyz_local[1]);
292 xyz_inout_z[i] =
static_cast<VT
>(xyz_local[2]);
295 }
else if (dimw_ == 2) {
297 for (
const auto& vit : vertices(gridView_, Dune::Partitions::all)) {
298 auto xyz_local = vit.geometry().corner(0);
299 xyz_inout.data()[i] =
static_cast<VT
>(xyz_local[0]);
300 xyz_inout_y[i] =
static_cast<VT
>(xyz_local[1]);
301 xyz_inout_z[i] =
static_cast<VT
>(td);
308template <
class Gr
idView,
unsigned int partitions>
309template <
typename Integer>
314 if (max_size < ncorners_) {
316 OPM_THROW(std::runtime_error,
317 "Opm::GridDataOutput::writeConnectivity( T* connectivity_inout,... ) "
319 +
") is not sufficient to fit the ncorners_ values (" +
std::to_string(ncorners_) +
")");
323 if (whichOrder ==
DUNE) {
325 for (
const auto& cit : elements(gridView_, dunePartition_)) {
326 auto cell_corners = cit.geometry().corners();
327 for (
auto vx = 0; vx < cell_corners; ++vx) {
328 const int vxIdx = gridView_.indexSet().subIndex(cit, vx, 3);
329 connectivity_inout[i + vx] = vxIdx;
335 for (
const auto& cit : elements(gridView_, dunePartition_)) {
336 auto cell_corners = cit.geometry().corners();
337 for (
auto vx = 0; vx < cell_corners; ++vx) {
338 const int vxIdx = gridView_.indexSet().subIndex(cit, vx, 3);
340 vtkOrder = Dune::VTK::renumber(cit.type(), vx);
341 connectivity_inout[i + vtkOrder] = vxIdx;
349template <
class Gr
idView,
unsigned int partitions>
350template <
typename VectType>
355 const std::size_t check_size = connectivity_inout.size();
357 if (check_size <
static_cast<std::size_t
>(ncorners_)) {
359 OPM_THROW(std::runtime_error,
360 "Opm::GridDataOutput::writeConnectivity( VectType& "
361 "connectivity_inout ) "
363 +
") is not sufficient to fit the ncorners_ values (" +
std::to_string(ncorners_) +
")");
367 if (whichOrder ==
DUNE) {
369 for (
const auto& cit : elements(gridView_, dunePartition_)) {
370 auto cell_corners = cit.geometry().corners();
371 for (
auto vx = 0; vx < cell_corners; ++vx) {
372 const int vxIdx = gridView_.indexSet().subIndex(cit, vx, 3);
373 connectivity_inout.data()[i + vx] = vxIdx;
379 for (
const auto& cit : elements(gridView_, dunePartition_)) {
380 auto cell_corners = cit.geometry().corners();
381 for (
auto vx = 0; vx < cell_corners; ++vx) {
382 const int vxIdx = gridView_.indexSet().subIndex(cit, vx, 3);
384 vtkOrder = Dune::VTK::renumber(cit.type(), vx);
385 connectivity_inout.data()[i + vtkOrder] = vxIdx;
393template <
class Gr
idView,
unsigned int partitions>
394template <
typename Integer>
398 if (max_size < ncells_) {
400 OPM_THROW(std::runtime_error,
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 ("
406 offsets_inout[0] = 0;
407 for (
const auto& cit : elements(gridView_, dunePartition_)) {
408 auto cell_corners = cit.geometry().corners();
409 offsets_inout[i] = offsets_inout[i - 1] + cell_corners;
415template <
class Gr
idView,
unsigned int partitions>
416template <
typename VectType>
420 const std::size_t check_size = offsets_inout.size();
421 if (check_size <
static_cast<std::size_t
>(ncells_)) {
423 OPM_THROW(std::runtime_error,
424 "Opm::GridDataOutput::writeOffsetsCells( VectType& "
427 +
") is not sufficient to fit the ncells_ values (" +
std::to_string(ncells_) +
")");
433 offsets_inout.data()[0] = 0;
434 for (
const auto& cit : elements(gridView_, dunePartition_)) {
435 auto cell_corners = cit.geometry().corners();
436 offsets_inout.data()[i] = offsets_inout.data()[i - 1] + cell_corners;
442template <
class Gr
idView,
unsigned int partitions>
443template <
typename Integer>
447 if (max_size < ncells_) {
449 OPM_THROW(std::runtime_error,
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 ("
455 for (
const auto& cit : elements(gridView_, dunePartition_)) {
456 Integer vtktype =
static_cast<Integer
>(Dune::VTK::geometryType(cit.type()));
457 types_inout[i++] = vtktype;
462template <
class Gr
idView,
unsigned int partitions>
463template <
typename VectType>
467 const std::size_t check_size = types_inout.size();
469 if (check_size <
static_cast<std::size_t
>(ncells_)) {
470 OPM_THROW(std::runtime_error,
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 ("
477 for (
const auto& cit : elements(gridView_, dunePartition_)) {
478 int vtktype =
static_cast<int>(Dune::VTK::geometryType(cit.type()));
479 types_inout.data()[i++] = vtktype;
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>
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.
long writeGridPoints_AOS(T *xyz_inout, long max_size=0) const
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()
Definition: GridDataOutput_impl.hpp:57
long writeGridPoints(T *x_inout, T *y_inout, T *z_inout, long max_size=0) const
Definition: GridDataOutput_impl.hpp:76
long writeOffsetsCells(Integer *offsets_inout, long max_size=0) const
Definition: GridDataOutput_impl.hpp:396
long writeConnectivity(Integer *connectivity_inout, ConnectivityVertexOrder whichOrder, long max_size=0) const
Definition: GridDataOutput_impl.hpp:311
unsigned int partition_value_
Definition: GridDataOutput.hpp:359
void printGridDetails(std::ostream &outstr) const
Definition: GridDataOutput_impl.hpp:518
bool polyhedralCellPresent() const
Definition: GridDataOutput_impl.hpp:45
long writeCellTypes(Integer *types_inout, long max_size=0) const
Definition: GridDataOutput_impl.hpp:445
int dimw_
Definition: GridDataOutput.hpp:374
long writeGridPoints_SOA(T *xyz_inout, long max_size=0) const
Definition: GridDataOutput_impl.hpp:229
std::string getPartitionTypeString() const
Definition: GridDataOutput_impl.hpp:486
Definition: GridDataOutput.hpp:89
ConnectivityVertexOrder
Definition: GridDataOutput.hpp:93
@ DUNE
Definition: GridDataOutput.hpp:93
std::string to_string(const ConvergenceReport::ReservoirFailure::Type t)