diff --git a/dune/vtk/CMakeLists.txt b/dune/vtk/CMakeLists.txt
index 4779d6df3f8bbc8c4e17129bb24eeec1dd8b8f95..49f14785ec3b546a35e00ed0f040094c2b42e501 100644
--- a/dune/vtk/CMakeLists.txt
+++ b/dune/vtk/CMakeLists.txt
@@ -1,5 +1,5 @@
dune_add_library("vtktypes" OBJECT
- vtktypes.cc)
+ types.cc)
#install headers
install(FILES
@@ -8,21 +8,23 @@ install(FILES
defaultvtkfunction.hh
filereader.hh
filewriter.hh
+ forward.hh
+ gridcreatorinterface.hh
legacyvtkfunction.hh
pvdwriter.hh
pvdwriter.impl.hh
- vtkfunction.hh
- vtklocalfunction.hh
- vtklocalfunctioninterface.hh
+ function.hh
+ localfunction.hh
+ localfunctioninterface.hh
vtkreader.hh
vtkreader.impl.hh
vtktimeserieswriter.hh
vtktimeserieswriter.impl.hh
- vtktypes.hh
+ types.hh
vtkwriter.hh
vtkwriterinterface.hh
vtkwriterinterface.impl.hh
- DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter)
+ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtk)
add_subdirectory(datacollectors)
add_subdirectory(gridcreators)
diff --git a/dune/vtk/datacollectorinterface.hh b/dune/vtk/datacollectorinterface.hh
index d37ca636ba8cb36c6e24b4c29d562df54b5ded54..a2153bc5b8aa2afaf2a3078f64ecf0d0b2922f20 100644
--- a/dune/vtk/datacollectorinterface.hh
+++ b/dune/vtk/datacollectorinterface.hh
@@ -5,130 +5,133 @@
namespace Dune
{
- /// Base class for data collectors in a CRTP style.
- /**
- * \tparam GridViewType Model of Dune::GridView
- * \tparam Derived Implementation of a concrete DataCollector.
- * \tparam Partition Dune::PartitionType [Partitions::InteriorBorder]
- **/
- template <class GridViewType, class Derived, class Partition>
- class DataCollectorInterface
+ namespace Vtk
{
- public:
- /// The partitionset to collect data from
- static constexpr auto partition = Partition{};
-
- using GridView = GridViewType;
-
- /// The dimension of the grid
- enum { dim = GridView::dimension };
-
- /// The dimension of the world
- enum { dow = GridView::dimensionworld };
-
- public:
- /// Store a copy of the GridView
- DataCollectorInterface (GridView const& gridView)
- : gridView_(gridView)
- {}
-
- /// Update the DataCollector on the current GridView
- void update ()
- {
- asDerived().updateImpl();
- }
-
- /// Return the number of ghost elements
- int ghostLevel () const
- {
- return asDerived().ghostLevelImpl();
- }
-
- /// Return the number of cells in (this partition of the) grid
- std::uint64_t numCells () const
- {
- return asDerived().numCellsImpl();
- }
-
- /// Return the number of points in (this partition of the) grid
- std::uint64_t numPoints () const
- {
- return asDerived().numPointsImpl();
- }
-
- /// \brief Return a flat vector of point coordinates
- /**
- * All coordinates are extended to 3 components and concatenated.
- * [p0_x, p0_y, p0_z, p1_x, p1_y, p1_z, ...]
- * If the GridView::dimensionworld < 3, the remaining components are
- * set to 0
- **/
- template <class T>
- std::vector<T> points () const
- {
- return asDerived().template pointsImpl<T>();
- }
-
- /// \brief Return a flat vector of function values evaluated at the points.
+ /// Base class for data collectors in a CRTP style.
/**
- * In case of a vector valued function, flat the vector entries:
- * [fct(p0)_0, fct(p0)_1, fct(p0)_2, fct(p1)_0, ...]
- * where the vector dimension must be 3 (possible extended by 0s)
- * In case of tensor valued function, flat the tensor row-wise:
- * [fct(p0)_00, fct(p0)_01, fct(p0)_02, fct(p0)_10, fct(p0)_11, fct(p0)_12, fct(p0)_20...]
- * where the tensor dimension must be 3x3 (possible extended by 0s)
+ * \tparam GridViewType Model of Dune::GridView
+ * \tparam Derived Implementation of a concrete DataCollector.
+ * \tparam Partition Dune::PartitionType [Partitions::InteriorBorder]
**/
- template <class T, class VtkFunction>
- std::vector<T> pointData (VtkFunction const& fct) const
+ template <class GridViewType, class Derived, class Partition>
+ class DataCollectorInterface
{
- return asDerived().template pointDataImpl<T>(fct);
- }
-
- /// \brief Return a flat vector of function values evaluated at the cells in the order of traversal.
- /**
- * \see pointData.
- * Note: Cells might be described explicitly by connectivity, offsets, and types, e.g. in
- * an UnstructuredGrid, or might be described implicitly by the grid type, e.g. in
- * StructuredGrid.
- */
- template <class T, class VtkFunction>
- std::vector<T> cellData (VtkFunction const& fct) const
- {
- return asDerived().template cellDataImpl<T>(fct);
- }
-
- protected: // cast to derived type
-
- Derived& asDerived ()
- {
- return static_cast<Derived&>(*this);
- }
-
- const Derived& asDerived () const
- {
- return static_cast<const Derived&>(*this);
- }
-
- public: // default implementations
-
- void updateImpl ()
- {
- /* do nothing */
- }
-
- int ghostLevelImpl () const
- {
- return gridView_.overlapSize(0);
- }
-
- // Evaluate `fct` in center of cell.
- template <class T, class VtkFunction>
- std::vector<T> cellDataImpl (VtkFunction const& fct) const;
-
- protected:
- GridView gridView_;
- };
-
+ public:
+ /// The partitionset to collect data from
+ static constexpr auto partition = Partition{};
+
+ using GridView = GridViewType;
+
+ /// The dimension of the grid
+ enum { dim = GridView::dimension };
+
+ /// The dimension of the world
+ enum { dow = GridView::dimensionworld };
+
+ public:
+ /// Store a copy of the GridView
+ DataCollectorInterface (GridView const& gridView)
+ : gridView_(gridView)
+ {}
+
+ /// Update the DataCollector on the current GridView
+ void update ()
+ {
+ asDerived().updateImpl();
+ }
+
+ /// Return the number of ghost elements
+ int ghostLevel () const
+ {
+ return asDerived().ghostLevelImpl();
+ }
+
+ /// Return the number of cells in (this partition of the) grid
+ std::uint64_t numCells () const
+ {
+ return asDerived().numCellsImpl();
+ }
+
+ /// Return the number of points in (this partition of the) grid
+ std::uint64_t numPoints () const
+ {
+ return asDerived().numPointsImpl();
+ }
+
+ /// \brief Return a flat vector of point coordinates
+ /**
+ * All coordinates are extended to 3 components and concatenated.
+ * [p0_x, p0_y, p0_z, p1_x, p1_y, p1_z, ...]
+ * If the GridView::dimensionworld < 3, the remaining components are
+ * set to 0
+ **/
+ template <class T>
+ std::vector<T> points () const
+ {
+ return asDerived().template pointsImpl<T>();
+ }
+
+ /// \brief Return a flat vector of function values evaluated at the points.
+ /**
+ * In case of a vector valued function, flat the vector entries:
+ * [fct(p0)_0, fct(p0)_1, fct(p0)_2, fct(p1)_0, ...]
+ * where the vector dimension must be 3 (possible extended by 0s)
+ * In case of tensor valued function, flat the tensor row-wise:
+ * [fct(p0)_00, fct(p0)_01, fct(p0)_02, fct(p0)_10, fct(p0)_11, fct(p0)_12, fct(p0)_20...]
+ * where the tensor dimension must be 3x3 (possible extended by 0s)
+ **/
+ template <class T, class VtkFunction>
+ std::vector<T> pointData (VtkFunction const& fct) const
+ {
+ return asDerived().template pointDataImpl<T>(fct);
+ }
+
+ /// \brief Return a flat vector of function values evaluated at the cells in the order of traversal.
+ /**
+ * \see pointData.
+ * Note: Cells might be described explicitly by connectivity, offsets, and types, e.g. in
+ * an UnstructuredGrid, or might be described implicitly by the grid type, e.g. in
+ * StructuredGrid.
+ */
+ template <class T, class VtkFunction>
+ std::vector<T> cellData (VtkFunction const& fct) const
+ {
+ return asDerived().template cellDataImpl<T>(fct);
+ }
+
+ protected: // cast to derived type
+
+ Derived& asDerived ()
+ {
+ return static_cast<Derived&>(*this);
+ }
+
+ const Derived& asDerived () const
+ {
+ return static_cast<const Derived&>(*this);
+ }
+
+ public: // default implementations
+
+ void updateImpl ()
+ {
+ /* do nothing */
+ }
+
+ int ghostLevelImpl () const
+ {
+ return gridView_.overlapSize(0);
+ }
+
+ // Evaluate `fct` in center of cell.
+ template <class T, class VtkFunction>
+ std::vector<T> cellDataImpl (VtkFunction const& fct) const;
+
+ protected:
+ GridView gridView_;
+ };
+
+ } // end namespace Vtk
} // end namespace Dune
#include "datacollectorinterface.impl.hh"
diff --git a/dune/vtk/datacollectorinterface.impl.hh b/dune/vtk/datacollectorinterface.impl.hh
index 957a6c0d4714e752edafa7885ebc6986de5ee05a..ada22528487829a922dd787b76cbb3e2ee6e985e 100644
--- a/dune/vtk/datacollectorinterface.impl.hh
+++ b/dune/vtk/datacollectorinterface.impl.hh
@@ -3,6 +3,7 @@
#include <dune/geometry/referenceelements.hh>
namespace Dune {
+namespace Vtk {
template <class GV, class D, class P>
template <class T, class VtkFunction>
@@ -23,4 +24,4 @@ std::vector<T> DataCollectorInterface<GV,D,P>
return data;
}
-} // end namespace Dune
+} } // end namespace Dune::Vtk
diff --git a/dune/vtk/datacollectors/CMakeLists.txt b/dune/vtk/datacollectors/CMakeLists.txt
index 0b6e52ac5ee3a43ee6d65cb8a59427310408748f..e7af67657c495009bece8a542d20efa91427ac4d 100644
--- a/dune/vtk/datacollectors/CMakeLists.txt
+++ b/dune/vtk/datacollectors/CMakeLists.txt
@@ -7,6 +7,6 @@ install(FILES
structureddatacollector.hh
unstructureddatacollector.hh
yaspdatacollector.hh
- DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter/datacollectors)
+ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtk/datacollectors)
add_subdirectory(test)
\ No newline at end of file
diff --git a/dune/vtk/datacollectors/continuousdatacollector.hh b/dune/vtk/datacollectors/continuousdatacollector.hh
index e6460d1accd6b41d2a44198954c0627a669f6d41..4dc1de48142656ea1da75dc46fbda8277284167f 100644
--- a/dune/vtk/datacollectors/continuousdatacollector.hh
+++ b/dune/vtk/datacollectors/continuousdatacollector.hh
@@ -7,140 +7,142 @@
#include <dune/grid/utility/globalindexset.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include "unstructureddatacollector.hh"
namespace Dune
{
+ namespace Vtk
+ {
+ /// Implementation of \ref DataCollector for linear cells, with continuous data.
+ template <class GridView, class Partition>
+ class ContinuousDataCollector
+ : public UnstructuredDataCollectorInterface<GridView, ContinuousDataCollector<GridView,Partition>, Partition>
+ {
+ using Self = ContinuousDataCollector;
+ using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
+
+ public:
+ using Super::dim;
+ using Super::partition;
+
+ public:
+ ContinuousDataCollector (GridView const& gridView)
+ : Super(gridView)
+ {}
+
+ /// Collect the vertex indices
+ void updateImpl ()
+ {
+ numPoints_ = 0;
+ indexMap_.resize(gridView_.size(dim));
+ auto const& indexSet = gridView_.indexSet();
+ for (auto const& vertex : vertices(gridView_, partition))
+ indexMap_[indexSet.index(vertex)] = std::int64_t(numPoints_++);
+
+ if (gridView_.comm().size() > 1) {
+ auto&& e = elements(gridView_, partition);
+ numCells_ = std::distance(std::begin(e), std::end(e));
+ } else {
+ numCells_ = gridView_.size(0);
+ }
+ }
-/// Implementation of \ref DataCollector for linear cells, with continuous data.
-template <class GridView, class Partition>
-class ContinuousDataCollector
- : public UnstructuredDataCollectorInterface<GridView, ContinuousDataCollector<GridView,Partition>, Partition>
-{
- using Self = ContinuousDataCollector;
- using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
+ /// Return number of grid vertices
+ std::uint64_t numPointsImpl () const
+ {
+ return numPoints_;
+ }
-public:
- using Super::dim;
- using Super::partition;
+ /// Return the coordinates of all grid vertices in the order given by the indexSet
+ template <class T>
+ std::vector<T> pointsImpl () const
+ {
+ std::vector<T> data;
+ data.reserve(numPoints_ * 3);
+ for (auto const& vertex : vertices(gridView_, partition)) {
+ auto v = vertex.geometry().center();
+ for (std::size_t j = 0; j < v.size(); ++j)
+ data.emplace_back(v[j]);
+ for (std::size_t j = v.size(); j < 3u; ++j)
+ data.emplace_back(0);
+ }
+ return data;
+ }
-public:
- ContinuousDataCollector (GridView const& gridView)
- : Super(gridView)
- {}
+ /// Return a vector of global unique ids of the points
+ std::vector<std::uint64_t> pointIdsImpl () const
+ {
+ std::vector<std::uint64_t> data;
+ data.reserve(numPoints_);
+ GlobalIndexSet<GridView> globalIndexSet(gridView_, dim);
+ for (auto const& vertex : vertices(gridView_, partition)) {
+ data.emplace_back(globalIndexSet.index(vertex));
+ }
+ return data;
+ }
- /// Collect the vertex indices
- void updateImpl ()
- {
- numPoints_ = 0;
- indexMap_.resize(gridView_.size(dim));
- auto const& indexSet = gridView_.indexSet();
- for (auto const& vertex : vertices(gridView_, partition))
- indexMap_[indexSet.index(vertex)] = std::int64_t(numPoints_++);
-
- if (gridView_.comm().size() > 1) {
- auto&& e = elements(gridView_, partition);
- numCells_ = std::distance(std::begin(e), std::end(e));
- } else {
- numCells_ = gridView_.size(0);
- }
- }
-
- /// Return number of grid vertices
- std::uint64_t numPointsImpl () const
- {
- return numPoints_;
- }
+ /// Return number of grid cells
+ std::uint64_t numCellsImpl () const
+ {
+ return numCells_;
+ }
- /// Return the coordinates of all grid vertices in the order given by the indexSet
- template <class T>
- std::vector<T> pointsImpl () const
- {
- std::vector<T> data;
- data.reserve(numPoints_ * 3);
- for (auto const& vertex : vertices(gridView_, partition)) {
- auto v = vertex.geometry().center();
- for (std::size_t j = 0; j < v.size(); ++j)
- data.emplace_back(v[j]);
- for (std::size_t j = v.size(); j < 3u; ++j)
- data.emplace_back(0);
- }
- return data;
- }
-
- /// Return a vector of global unique ids of the points
- std::vector<std::uint64_t> pointIdsImpl () const
- {
- std::vector<std::uint64_t> data;
- data.reserve(numPoints_);
- GlobalIndexSet<GridView> globalIndexSet(gridView_, dim);
- for (auto const& vertex : vertices(gridView_, partition)) {
- data.emplace_back(globalIndexSet.index(vertex));
- }
- return data;
- }
-
- /// Return number of grid cells
- std::uint64_t numCellsImpl () const
- {
- return numCells_;
- }
+ /// Return the types, offsets and connectivity of the cells, using the same connectivity as
+ /// given by the grid.
+ Cells cellsImpl () const
+ {
+ auto const& indexSet = gridView_.indexSet();
+ auto types = indexSet.types(0);
+ int maxVertices = std::accumulate(types.begin(), types.end(), 1, [](int m, GeometryType t) {
+ auto refElem = referenceElement<double,dim>(t);
+ return std::max(m, refElem.size(dim));
+ });
+
+ Cells cells;
+ cells.connectivity.reserve(numCells_ * maxVertices);
+ cells.offsets.reserve(numCells_);
+ cells.types.reserve(numCells_);
+
+ std::int64_t old_o = 0;
+ for (auto const& c : elements(gridView_, partition)) {
+ Vtk::CellType cellType(c.type());
+ for (unsigned int j = 0; j < c.subEntities(dim); ++j)
+ cells.connectivity.emplace_back(indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)]);
+ cells.offsets.push_back(old_o += c.subEntities(dim));
+ cells.types.push_back(cellType.type());
+ }
+ return cells;
+ }
- /// Return the types, offsets and connectivity of the cells, using the same connectivity as
- /// given by the grid.
- Cells cellsImpl () const
- {
- auto const& indexSet = gridView_.indexSet();
- auto types = indexSet.types(0);
- int maxVertices = std::accumulate(types.begin(), types.end(), 1, [](int m, GeometryType t) {
- auto refElem = referenceElement<double,dim>(t);
- return std::max(m, refElem.size(dim));
- });
-
- Cells cells;
- cells.connectivity.reserve(numCells_ * maxVertices);
- cells.offsets.reserve(numCells_);
- cells.types.reserve(numCells_);
-
- std::int64_t old_o = 0;
- for (auto const& c : elements(gridView_, partition)) {
- Vtk::CellType cellType(c.type());
- for (unsigned int j = 0; j < c.subEntities(dim); ++j)
- cells.connectivity.emplace_back(indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)]);
- cells.offsets.push_back(old_o += c.subEntities(dim));
- cells.types.push_back(cellType.type());
- }
- return cells;
- }
-
- /// Evaluate the `fct` at the corners of the elements
- template <class T, class GlobalFunction>
- std::vector<T> pointDataImpl (GlobalFunction const& fct) const
- {
- std::vector<T> data(numPoints_ * fct.ncomps());
- auto const& indexSet = gridView_.indexSet();
- auto localFct = localFunction(fct);
- for (auto const& e : elements(gridView_, partition)) {
- localFct.bind(e);
- Vtk::CellType cellType{e.type()};
- auto refElem = referenceElement(e.geometry());
- for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
- std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e,cellType.permutation(j),dim)];
- for (int comp = 0; comp < fct.ncomps(); ++comp)
- data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+ /// Evaluate the `fct` at the corners of the elements
+ template <class T, class GlobalFunction>
+ std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+ {
+ std::vector<T> data(numPoints_ * fct.ncomps());
+ auto const& indexSet = gridView_.indexSet();
+ auto localFct = localFunction(fct);
+ for (auto const& e : elements(gridView_, partition)) {
+ localFct.bind(e);
+ Vtk::CellType cellType{e.type()};
+ auto refElem = referenceElement(e.geometry());
+ for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
+ std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e,cellType.permutation(j),dim)];
+ for (int comp = 0; comp < fct.ncomps(); ++comp)
+ data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+ }
+ localFct.unbind();
+ }
+ return data;
}
- localFct.unbind();
- }
- return data;
- }
-
-protected:
- using Super::gridView_;
- std::uint64_t numPoints_ = 0;
- std::uint64_t numCells_ = 0;
- std::vector<std::int64_t> indexMap_;
-};
+ protected:
+ using Super::gridView_;
+ std::uint64_t numPoints_ = 0;
+ std::uint64_t numCells_ = 0;
+ std::vector<std::int64_t> indexMap_;
+ };
+
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/discontinuousdatacollector.hh b/dune/vtk/datacollectors/discontinuousdatacollector.hh
index eab3952372315889a5a226ad892a65c3a9c96c8a..8d75ae890bdfd56cdfcc21d5f2b9fb1bc4d4bc6a 100644
--- a/dune/vtk/datacollectors/discontinuousdatacollector.hh
+++ b/dune/vtk/datacollectors/discontinuousdatacollector.hh
@@ -3,126 +3,128 @@
#include <vector>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include "unstructureddatacollector.hh"
namespace Dune
{
-
-/// Implementation of \ref DataCollector for linear cells, with discontinuous data.
-template <class GridView, class Partition>
-class DiscontinuousDataCollector
- : public UnstructuredDataCollectorInterface<GridView, DiscontinuousDataCollector<GridView,Partition>, Partition>
-{
- using Self = DiscontinuousDataCollector;
- using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
-
-public:
- using Super::dim;
- using Super::partition;
-
-public:
- DiscontinuousDataCollector (GridView const& gridView)
- : Super(gridView)
- {}
-
- /// Create an index map the uniquely assigns an index to each pair (element,corner)
- void updateImpl ()
- {
- numPoints_ = 0;
- numCells_ = 0;
- indexMap_.resize(gridView_.size(dim));
- std::int64_t vertex_idx = 0;
- auto const& indexSet = gridView_.indexSet();
- for (auto const& c : elements(gridView_, partition)) {
- numCells_++;
- numPoints_ += c.subEntities(dim);
- for (unsigned int i = 0; i < c.subEntities(dim); ++i)
- indexMap_[indexSet.subIndex(c, i, dim)] = vertex_idx++;
- }
- }
-
- /// The number of points approx. #cell * #corners-per-cell
- std::uint64_t numPointsImpl () const
+ namespace Vtk
{
- return numPoints_;
- }
+ /// Implementation of \ref DataCollector for linear cells, with discontinuous data.
+ template <class GridView, class Partition>
+ class DiscontinuousDataCollector
+ : public UnstructuredDataCollectorInterface<GridView, DiscontinuousDataCollector<GridView,Partition>, Partition>
+ {
+ using Self = DiscontinuousDataCollector;
+ using Super = UnstructuredDataCollectorInterface<GridView, Self, Partition>;
+
+ public:
+ using Super::dim;
+ using Super::partition;
+
+ public:
+ DiscontinuousDataCollector (GridView const& gridView)
+ : Super(gridView)
+ {}
+
+ /// Create an index map the uniquely assigns an index to each pair (element,corner)
+ void updateImpl ()
+ {
+ numPoints_ = 0;
+ numCells_ = 0;
+ indexMap_.resize(gridView_.size(dim));
+ std::int64_t vertex_idx = 0;
+ auto const& indexSet = gridView_.indexSet();
+ for (auto const& c : elements(gridView_, partition)) {
+ numCells_++;
+ numPoints_ += c.subEntities(dim);
+ for (unsigned int i = 0; i < c.subEntities(dim); ++i)
+ indexMap_[indexSet.subIndex(c, i, dim)] = vertex_idx++;
+ }
+ }
- /// Return the coordinates of the corners of all cells
- template <class T>
- std::vector<T> pointsImpl () const
- {
- std::vector<T> data(numPoints_ * 3);
- auto const& indexSet = gridView_.indexSet();
- for (auto const& element : elements(gridView_, partition)) {
- for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
- std::size_t idx = 3 * indexMap_[indexSet.subIndex(element, i, dim)];
- auto v = element.geometry().corner(i);
- for (std::size_t j = 0; j < v.size(); ++j)
- data[idx + j] = T(v[j]);
- for (std::size_t j = v.size(); j < 3u; ++j)
- data[idx + j] = T(0);
+ /// The number of points approx. #cell * #corners-per-cell
+ std::uint64_t numPointsImpl () const
+ {
+ return numPoints_;
}
- }
- return data;
- }
- /// Return number of grid cells
- std::uint64_t numCellsImpl () const
- {
- return numCells_;
- }
+ /// Return the coordinates of the corners of all cells
+ template <class T>
+ std::vector<T> pointsImpl () const
+ {
+ std::vector<T> data(numPoints_ * 3);
+ auto const& indexSet = gridView_.indexSet();
+ for (auto const& element : elements(gridView_, partition)) {
+ for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
+ std::size_t idx = 3 * indexMap_[indexSet.subIndex(element, i, dim)];
+ auto v = element.geometry().corner(i);
+ for (std::size_t j = 0; j < v.size(); ++j)
+ data[idx + j] = T(v[j]);
+ for (std::size_t j = v.size(); j < 3u; ++j)
+ data[idx + j] = T(0);
+ }
+ }
+ return data;
+ }
- /// Connect the corners of each cell. The leads to a global discontinuous grid
- Cells cellsImpl () const
- {
- Cells cells;
- cells.connectivity.reserve(numPoints_);
- cells.offsets.reserve(numCells_);
- cells.types.reserve(numCells_);
-
- std::int64_t old_o = 0;
- auto const& indexSet = gridView_.indexSet();
- for (auto const& c : elements(gridView_, partition)) {
- Vtk::CellType cellType(c.type());
- for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
- std::int64_t vertex_idx = indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)];
- cells.connectivity.push_back(vertex_idx);
+ /// Return number of grid cells
+ std::uint64_t numCellsImpl () const
+ {
+ return numCells_;
}
- cells.offsets.push_back(old_o += c.subEntities(dim));
- cells.types.push_back(cellType.type());
- }
- return cells;
- }
+ /// Connect the corners of each cell. The leads to a global discontinuous grid
+ Cells cellsImpl () const
+ {
+ Cells cells;
+ cells.connectivity.reserve(numPoints_);
+ cells.offsets.reserve(numCells_);
+ cells.types.reserve(numCells_);
+
+ std::int64_t old_o = 0;
+ auto const& indexSet = gridView_.indexSet();
+ for (auto const& c : elements(gridView_, partition)) {
+ Vtk::CellType cellType(c.type());
+ for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
+ std::int64_t vertex_idx = indexMap_[indexSet.subIndex(c,cellType.permutation(j),dim)];
+ cells.connectivity.push_back(vertex_idx);
+ }
+ cells.offsets.push_back(old_o += c.subEntities(dim));
+ cells.types.push_back(cellType.type());
+ }
+
+ return cells;
+ }
- /// Evaluate the `fct` in the corners of each cell
- template <class T, class GlobalFunction>
- std::vector<T> pointDataImpl (GlobalFunction const& fct) const
- {
- std::vector<T> data(numPoints_ * fct.ncomps());
- auto const& indexSet = gridView_.indexSet();
- auto localFct = localFunction(fct);
- for (auto const& e : elements(gridView_, partition)) {
- localFct.bind(e);
- Vtk::CellType cellType{e.type()};
- auto refElem = referenceElement(e.geometry());
- for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
- std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e, cellType.permutation(j), dim)];
- for (int comp = 0; comp < fct.ncomps(); ++comp)
- data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+ /// Evaluate the `fct` in the corners of each cell
+ template <class T, class GlobalFunction>
+ std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+ {
+ std::vector<T> data(numPoints_ * fct.ncomps());
+ auto const& indexSet = gridView_.indexSet();
+ auto localFct = localFunction(fct);
+ for (auto const& e : elements(gridView_, partition)) {
+ localFct.bind(e);
+ Vtk::CellType cellType{e.type()};
+ auto refElem = referenceElement(e.geometry());
+ for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
+ std::size_t idx = fct.ncomps() * indexMap_[indexSet.subIndex(e, cellType.permutation(j), dim)];
+ for (int comp = 0; comp < fct.ncomps(); ++comp)
+ data[idx + comp] = T(localFct.evaluate(comp, refElem.position(cellType.permutation(j),dim)));
+ }
+ localFct.unbind();
+ }
+ return data;
}
- localFct.unbind();
- }
- return data;
- }
-
-protected:
- using Super::gridView_;
- std::uint64_t numCells_ = 0;
- std::uint64_t numPoints_ = 0;
- std::vector<std::int64_t> indexMap_;
-};
+ protected:
+ using Super::gridView_;
+ std::uint64_t numCells_ = 0;
+ std::uint64_t numPoints_ = 0;
+ std::vector<std::int64_t> indexMap_;
+ };
+
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/lagrangedatacollector.hh b/dune/vtk/datacollectors/lagrangedatacollector.hh
index d38f3411026548a241761bf6924c7cda87154a57..684c3d9890c95fca61c8088ddc8949e07b1a380f 100644
--- a/dune/vtk/datacollectors/lagrangedatacollector.hh
+++ b/dune/vtk/datacollectors/lagrangedatacollector.hh
@@ -8,189 +8,192 @@
#include <dune/grid/common/partitionset.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/utility/lagrangepoints.hh>
#include "unstructureddatacollector.hh"
-namespace Dune {
-
-/// Implementation of \ref DataCollector for lagrange cells
-template <class GridView, int ORDER = -1>
-class LagrangeDataCollector
- : public UnstructuredDataCollectorInterface<GridView, LagrangeDataCollector<GridView,ORDER>, Partitions::All>
+namespace Dune
{
- using Self = LagrangeDataCollector;
- using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
-
-public:
- static_assert(ORDER != 0, "Order 0 not supported");
- using Super::dim;
- using Super::partition; // NOTE: lagrange data-collector currently implemented for the All partition only
-
-public:
- LagrangeDataCollector (GridView const& gridView, int order = (ORDER < 0 ? 2 : ORDER))
- : Super(gridView)
- , order_(order)
+ namespace Vtk
{
- assert(order > 0 && "Order 0 not supported");
- assert(ORDER < 0 || order == ORDER);
- }
+ /// Implementation of \ref DataCollector for lagrange cells
+ template <class GridView, int ORDER = -1>
+ class LagrangeDataCollector
+ : public UnstructuredDataCollectorInterface<GridView, LagrangeDataCollector<GridView,ORDER>, Partitions::All>
+ {
+ using Self = LagrangeDataCollector;
+ using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
+
+ public:
+ static_assert(ORDER != 0, "Order 0 not supported");
+ using Super::dim;
+ using Super::partition; // NOTE: lagrange data-collector currently implemented for the All partition only
+
+ public:
+ LagrangeDataCollector (GridView const& gridView, int order = (ORDER < 0 ? 2 : ORDER))
+ : Super(gridView)
+ , order_(order)
+ {
+ assert(order > 0 && "Order 0 not supported");
+ assert(ORDER < 0 || order == ORDER);
+ }
- /// Construct the point sets
- void updateImpl ()
- {
- auto const& indexSet = gridView_.indexSet();
+ /// Construct the point sets
+ void updateImpl ()
+ {
+ auto const& indexSet = gridView_.indexSet();
- pointSets_.clear();
- for (auto gt : indexSet.types(0))
- pointSets_.emplace(gt, order_);
+ pointSets_.clear();
+ for (auto gt : indexSet.types(0))
+ pointSets_.emplace(gt, order_);
- for (auto& pointSet : pointSets_)
- pointSet.second.build(pointSet.first);
+ for (auto& pointSet : pointSets_)
+ pointSet.second.build(pointSet.first);
- numPoints_ = indexSet.size(dim);
- for (auto const& pointSet : pointSets_) {
- auto gt = pointSet.first;
- auto refElem = referenceElement<double,dim>(gt);
- numPoints_ += (pointSet.second.size() - refElem.size(dim)) * indexSet.size(gt);
- }
- }
+ numPoints_ = indexSet.size(dim);
+ for (auto const& pointSet : pointSets_) {
+ auto gt = pointSet.first;
+ auto refElem = referenceElement<double,dim>(gt);
+ numPoints_ += (pointSet.second.size() - refElem.size(dim)) * indexSet.size(gt);
+ }
+ }
- /// Return number of lagrange nodes
- std::uint64_t numPointsImpl () const
- {
- return numPoints_;
- }
-
- /// Return a vector of point coordinates.
- /**
- * The vector of point coordinates is composed of vertex coordinates first and second
- * edge center coordinates.
- **/
- template <class T>
- std::vector<T> pointsImpl () const
- {
- std::vector<T> data(this->numPoints() * 3);
- auto const& indexSet = gridView_.indexSet();
-
- std::size_t shift = indexSet.size(dim);
-
- for (auto const& element : elements(gridView_, partition)) {
- auto geometry = element.geometry();
- auto refElem = referenceElement<T,dim>(element.type());
-
- auto const& pointSet = pointSets_.at(element.type());
- unsigned int vertexDOFs = refElem.size(dim);
- unsigned int innerDOFs = pointSet.size() - vertexDOFs;
-
- for (std::size_t i = 0; i < pointSet.size(); ++i) {
- auto const& p = pointSet[i];
- if (i < vertexDOFs)
- assert(p.localKey().codim() == dim);
-
- auto const& localKey = p.localKey();
- std::size_t idx = 3 * (localKey.codim() == dim
- ? indexSet.subIndex(element, localKey.subEntity(), dim)
- : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
-
- auto v = geometry.global(p.point());
- for (std::size_t j = 0; j < v.size(); ++j)
- data[idx + j] = T(v[j]);
- for (std::size_t j = v.size(); j < 3u; ++j)
- data[idx + j] = T(0);
+ /// Return number of lagrange nodes
+ std::uint64_t numPointsImpl () const
+ {
+ return numPoints_;
}
- }
- return data;
- }
- /// Return number of grid cells
- std::uint64_t numCellsImpl () const
- {
- return gridView_.size(0);
- }
-
- /// \brief Return cell types, offsets, and connectivity. \see Cells
- /**
- * The cell connectivity is composed of cell vertices first and second cell edges,
- * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
- **/
- Cells cellsImpl () const
- {
- Cells cells;
- cells.connectivity.reserve(this->numPoints());
- cells.offsets.reserve(this->numCells());
- cells.types.reserve(this->numCells());
-
- auto const& indexSet = gridView_.indexSet();
- std::size_t shift = indexSet.size(dim);
-
- std::int64_t old_o = 0;
- for (auto const& element : elements(gridView_, partition)) {
- auto refElem = referenceElement<double,dim>(element.type());
- Vtk::CellType cellType(element.type(), Vtk::LAGRANGE);
-
- auto const& pointSet = pointSets_.at(element.type());
- unsigned int vertexDOFs = refElem.size(dim);
- unsigned int innerDOFs = pointSet.size() - vertexDOFs;
-
- for (std::size_t i = 0; i < pointSet.size(); ++i) {
- auto const& p = pointSet[i];
- auto const& localKey = p.localKey();
- std::size_t idx = (localKey.codim() == dim
- ? indexSet.subIndex(element, localKey.subEntity(), dim)
- : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
- cells.connectivity.push_back(idx);
+ /// Return a vector of point coordinates.
+ /**
+ * The vector of point coordinates is composed of vertex coordinates first and second
+ * edge center coordinates.
+ **/
+ template <class T>
+ std::vector<T> pointsImpl () const
+ {
+ std::vector<T> data(this->numPoints() * 3);
+ auto const& indexSet = gridView_.indexSet();
+
+ std::size_t shift = indexSet.size(dim);
+
+ for (auto const& element : elements(gridView_, partition)) {
+ auto geometry = element.geometry();
+ auto refElem = referenceElement<T,dim>(element.type());
+
+ auto const& pointSet = pointSets_.at(element.type());
+ unsigned int vertexDOFs = refElem.size(dim);
+ unsigned int innerDOFs = pointSet.size() - vertexDOFs;
+
+ for (std::size_t i = 0; i < pointSet.size(); ++i) {
+ auto const& p = pointSet[i];
+ if (i < vertexDOFs)
+ assert(p.localKey().codim() == dim);
+
+ auto const& localKey = p.localKey();
+ std::size_t idx = 3 * (localKey.codim() == dim
+ ? indexSet.subIndex(element, localKey.subEntity(), dim)
+ : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
+
+ auto v = geometry.global(p.point());
+ for (std::size_t j = 0; j < v.size(); ++j)
+ data[idx + j] = T(v[j]);
+ for (std::size_t j = v.size(); j < 3u; ++j)
+ data[idx + j] = T(0);
+ }
+ }
+ return data;
+ }
+
+ /// Return number of grid cells
+ std::uint64_t numCellsImpl () const
+ {
+ return gridView_.size(0);
}
- cells.offsets.push_back(old_o += pointSet.size());
- cells.types.push_back(cellType.type());
- }
- return cells;
- }
+ /// \brief Return cell types, offsets, and connectivity. \see Cells
+ /**
+ * The cell connectivity is composed of cell vertices first and second cell edges,
+ * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
+ **/
+ Cells cellsImpl () const
+ {
+ Cells cells;
+ cells.connectivity.reserve(this->numPoints());
+ cells.offsets.reserve(this->numCells());
+ cells.types.reserve(this->numCells());
+
+ auto const& indexSet = gridView_.indexSet();
+ std::size_t shift = indexSet.size(dim);
+
+ std::int64_t old_o = 0;
+ for (auto const& element : elements(gridView_, partition)) {
+ auto refElem = referenceElement<double,dim>(element.type());
+ Vtk::CellType cellType(element.type(), Vtk::LAGRANGE);
+
+ auto const& pointSet = pointSets_.at(element.type());
+ unsigned int vertexDOFs = refElem.size(dim);
+ unsigned int innerDOFs = pointSet.size() - vertexDOFs;
+
+ for (std::size_t i = 0; i < pointSet.size(); ++i) {
+ auto const& p = pointSet[i];
+ auto const& localKey = p.localKey();
+ std::size_t idx = (localKey.codim() == dim
+ ? indexSet.subIndex(element, localKey.subEntity(), dim)
+ : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
+ cells.connectivity.push_back(idx);
+ }
+
+ cells.offsets.push_back(old_o += pointSet.size());
+ cells.types.push_back(cellType.type());
+ }
+ return cells;
+ }
- /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
- template <class T, class GlobalFunction>
- std::vector<T> pointDataImpl (GlobalFunction const& fct) const
- {
- int nComps = fct.ncomps();
- std::vector<T> data(this->numPoints() * nComps);
- auto const& indexSet = gridView_.indexSet();
-
- std::size_t shift = indexSet.size(dim);
-
- auto localFct = localFunction(fct);
- for (auto const& element : elements(gridView_, partition)) {
- localFct.bind(element);
- auto refElem = referenceElement<T,dim>(element.type());
-
- auto const& pointSet = pointSets_.at(element.type());
- unsigned int vertexDOFs = refElem.size(dim);
- unsigned int innerDOFs = pointSet.size() - vertexDOFs;
-
- for (std::size_t i = 0; i < pointSet.size(); ++i) {
- auto const& p = pointSet[i];
- auto const& localKey = p.localKey();
- std::size_t idx = nComps * (localKey.codim() == dim
- ? indexSet.subIndex(element, localKey.subEntity(), dim)
- : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
-
- for (int comp = 0; comp < nComps; ++comp)
- data[idx + comp] = T(localFct.evaluate(comp, p.point()));
+ /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
+ template <class T, class GlobalFunction>
+ std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+ {
+ int nComps = fct.ncomps();
+ std::vector<T> data(this->numPoints() * nComps);
+ auto const& indexSet = gridView_.indexSet();
+
+ std::size_t shift = indexSet.size(dim);
+
+ auto localFct = localFunction(fct);
+ for (auto const& element : elements(gridView_, partition)) {
+ localFct.bind(element);
+ auto refElem = referenceElement<T,dim>(element.type());
+
+ auto const& pointSet = pointSets_.at(element.type());
+ unsigned int vertexDOFs = refElem.size(dim);
+ unsigned int innerDOFs = pointSet.size() - vertexDOFs;
+
+ for (std::size_t i = 0; i < pointSet.size(); ++i) {
+ auto const& p = pointSet[i];
+ auto const& localKey = p.localKey();
+ std::size_t idx = nComps * (localKey.codim() == dim
+ ? indexSet.subIndex(element, localKey.subEntity(), dim)
+ : innerDOFs*indexSet.index(element) + (i - vertexDOFs) + shift);
+
+ for (int comp = 0; comp < nComps; ++comp)
+ data[idx + comp] = T(localFct.evaluate(comp, p.point()));
+ }
+ localFct.unbind();
+ }
+ return data;
}
- localFct.unbind();
- }
- return data;
- }
-protected:
- using Super::gridView_;
+ protected:
+ using Super::gridView_;
- unsigned int order_;
- std::uint64_t numPoints_ = 0;
+ unsigned int order_;
+ std::uint64_t numPoints_ = 0;
- using PointSet = VtkLagrangePointSet<typename GridView::ctype, GridView::dimension>;
- std::map<GeometryType, PointSet> pointSets_;
-};
+ using PointSet = LagrangePointSet<typename GridView::ctype, GridView::dimension>;
+ std::map<GeometryType, PointSet> pointSets_;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/quadraticdatacollector.hh b/dune/vtk/datacollectors/quadraticdatacollector.hh
index 0ddee00b9028c156d017ca9f430681439b4bbb8b..f4fdaeb040d47972f7ca44996892f0812dc011e3 100644
--- a/dune/vtk/datacollectors/quadraticdatacollector.hh
+++ b/dune/vtk/datacollectors/quadraticdatacollector.hh
@@ -6,140 +6,142 @@
#include <dune/grid/common/partitionset.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include "unstructureddatacollector.hh"
namespace Dune
{
-
-/// Implementation of \ref DataCollector for quadratic cells, with continuous data.
-template <class GridView>
-class QuadraticDataCollector
- : public UnstructuredDataCollectorInterface<GridView, QuadraticDataCollector<GridView>, Partitions::All>
-{
- using Self = QuadraticDataCollector;
- using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
-
-public:
- using Super::dim;
- using Super::partition; // NOTE: quadratic data-collector currently implemented for the All partition only
-
-public:
- QuadraticDataCollector (GridView const& gridView)
- : Super(gridView)
- {}
-
- /// Return number of vertices + number of edge
- std::uint64_t numPointsImpl () const
- {
- return gridView_.size(dim) + gridView_.size(dim-1);
- }
-
- /// Return a vector of point coordinates.
- /**
- * The vector of point coordinates is composed of vertex coordinates first and second
- * edge center coordinates.
- **/
- template <class T>
- std::vector<T> pointsImpl () const
+ namespace Vtk
{
- std::vector<T> data(this->numPoints() * 3);
- auto const& indexSet = gridView_.indexSet();
- for (auto const& element : elements(gridView_, partition)) {
- auto geometry = element.geometry();
- auto refElem = referenceElement<T,dim>(element.type());
-
- // vertices
- for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
- std::size_t idx = 3 * indexSet.subIndex(element, i, dim);
- auto v = geometry.global(refElem.position(i,dim));
- for (std::size_t j = 0; j < v.size(); ++j)
- data[idx + j] = T(v[j]);
- for (std::size_t j = v.size(); j < 3u; ++j)
- data[idx + j] = T(0);
- }
- // edge centers
- for (unsigned int i = 0; i < element.subEntities(dim-1); ++i) {
- std::size_t idx = 3 * (indexSet.subIndex(element, i, dim-1) + gridView_.size(dim));
- auto v = geometry.global(refElem.position(i,dim-1));
- for (std::size_t j = 0; j < v.size(); ++j)
- data[idx + j] = T(v[j]);
- for (std::size_t j = v.size(); j < 3u; ++j)
- data[idx + j] = T(0);
+ /// Implementation of \ref DataCollector for quadratic cells, with continuous data.
+ template <class GridView>
+ class QuadraticDataCollector
+ : public UnstructuredDataCollectorInterface<GridView, QuadraticDataCollector<GridView>, Partitions::All>
+ {
+ using Self = QuadraticDataCollector;
+ using Super = UnstructuredDataCollectorInterface<GridView, Self, Partitions::All>;
+
+ public:
+ using Super::dim;
+ using Super::partition; // NOTE: quadratic data-collector currently implemented for the All partition only
+
+ public:
+ QuadraticDataCollector (GridView const& gridView)
+ : Super(gridView)
+ {}
+
+ /// Return number of vertices + number of edge
+ std::uint64_t numPointsImpl () const
+ {
+ return gridView_.size(dim) + gridView_.size(dim-1);
}
- }
- return data;
- }
- /// Return number of grid cells
- std::uint64_t numCellsImpl () const
- {
- return gridView_.size(0);
- }
-
- /// \brief Return cell types, offsets, and connectivity. \see Cells
- /**
- * The cell connectivity is composed of cell vertices first and second cell edges,
- * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
- **/
- Cells cellsImpl () const
- {
- Cells cells;
- cells.connectivity.reserve(this->numPoints());
- cells.offsets.reserve(this->numCells());
- cells.types.reserve(this->numCells());
-
- std::int64_t old_o = 0;
- auto const& indexSet = gridView_.indexSet();
- for (auto const& c : elements(gridView_, partition)) {
- Vtk::CellType cellType(c.type(), Vtk::QUADRATIC);
- for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
- int k = cellType.permutation(j);
- std::int64_t point_idx = indexSet.subIndex(c,k,dim);
- cells.connectivity.push_back(point_idx);
+ /// Return a vector of point coordinates.
+ /**
+ * The vector of point coordinates is composed of vertex coordinates first and second
+ * edge center coordinates.
+ **/
+ template <class T>
+ std::vector<T> pointsImpl () const
+ {
+ std::vector<T> data(this->numPoints() * 3);
+ auto const& indexSet = gridView_.indexSet();
+ for (auto const& element : elements(gridView_, partition)) {
+ auto geometry = element.geometry();
+ auto refElem = referenceElement<T,dim>(element.type());
+
+ // vertices
+ for (unsigned int i = 0; i < element.subEntities(dim); ++i) {
+ std::size_t idx = 3 * indexSet.subIndex(element, i, dim);
+ auto v = geometry.global(refElem.position(i,dim));
+ for (std::size_t j = 0; j < v.size(); ++j)
+ data[idx + j] = T(v[j]);
+ for (std::size_t j = v.size(); j < 3u; ++j)
+ data[idx + j] = T(0);
+ }
+ // edge centers
+ for (unsigned int i = 0; i < element.subEntities(dim-1); ++i) {
+ std::size_t idx = 3 * (indexSet.subIndex(element, i, dim-1) + gridView_.size(dim));
+ auto v = geometry.global(refElem.position(i,dim-1));
+ for (std::size_t j = 0; j < v.size(); ++j)
+ data[idx + j] = T(v[j]);
+ for (std::size_t j = v.size(); j < 3u; ++j)
+ data[idx + j] = T(0);
+ }
+ }
+ return data;
}
- for (unsigned int j = 0; j < c.subEntities(dim-1); ++j) {
- int k = cellType.permutation(c.subEntities(dim) + j);
- std::int64_t point_idx = (indexSet.subIndex(c,k,dim-1) + gridView_.size(dim));
- cells.connectivity.push_back(point_idx);
+
+ /// Return number of grid cells
+ std::uint64_t numCellsImpl () const
+ {
+ return gridView_.size(0);
}
- cells.offsets.push_back(old_o += c.subEntities(dim)+c.subEntities(dim-1));
- cells.types.push_back(cellType.type());
- }
- return cells;
- }
-
- /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
- template <class T, class GlobalFunction>
- std::vector<T> pointDataImpl (GlobalFunction const& fct) const
- {
- std::vector<T> data(this->numPoints() * fct.ncomps());
- auto const& indexSet = gridView_.indexSet();
- auto localFct = localFunction(fct);
- for (auto const& e : elements(gridView_, partition)) {
- localFct.bind(e);
- Vtk::CellType cellType{e.type(), Vtk::QUADRATIC};
- auto refElem = referenceElement(e.geometry());
- for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
- int k = cellType.permutation(j);
- std::size_t idx = fct.ncomps() * indexSet.subIndex(e, k, dim);
- for (int comp = 0; comp < fct.ncomps(); ++comp)
- data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim)));
+
+ /// \brief Return cell types, offsets, and connectivity. \see Cells
+ /**
+ * The cell connectivity is composed of cell vertices first and second cell edges,
+ * where the indices are grouped [vertex-indices..., (#vertices)+edge-indices...]
+ **/
+ Cells cellsImpl () const
+ {
+ Cells cells;
+ cells.connectivity.reserve(this->numPoints());
+ cells.offsets.reserve(this->numCells());
+ cells.types.reserve(this->numCells());
+
+ std::int64_t old_o = 0;
+ auto const& indexSet = gridView_.indexSet();
+ for (auto const& c : elements(gridView_, partition)) {
+ Vtk::CellType cellType(c.type(), Vtk::QUADRATIC);
+ for (unsigned int j = 0; j < c.subEntities(dim); ++j) {
+ int k = cellType.permutation(j);
+ std::int64_t point_idx = indexSet.subIndex(c,k,dim);
+ cells.connectivity.push_back(point_idx);
+ }
+ for (unsigned int j = 0; j < c.subEntities(dim-1); ++j) {
+ int k = cellType.permutation(c.subEntities(dim) + j);
+ std::int64_t point_idx = (indexSet.subIndex(c,k,dim-1) + gridView_.size(dim));
+ cells.connectivity.push_back(point_idx);
+ }
+ cells.offsets.push_back(old_o += c.subEntities(dim)+c.subEntities(dim-1));
+ cells.types.push_back(cellType.type());
+ }
+ return cells;
}
- for (unsigned int j = 0; j < e.subEntities(dim-1); ++j) {
- int k = cellType.permutation(e.subEntities(dim) + j);
- std::size_t idx = fct.ncomps() * (indexSet.subIndex(e, k, dim-1) + gridView_.size(dim));
- for (int comp = 0; comp < fct.ncomps(); ++comp)
- data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim-1)));
+
+ /// Evaluate the `fct` at element vertices and edge centers in the same order as the point coords.
+ template <class T, class GlobalFunction>
+ std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+ {
+ std::vector<T> data(this->numPoints() * fct.ncomps());
+ auto const& indexSet = gridView_.indexSet();
+ auto localFct = localFunction(fct);
+ for (auto const& e : elements(gridView_, partition)) {
+ localFct.bind(e);
+ Vtk::CellType cellType{e.type(), Vtk::QUADRATIC};
+ auto refElem = referenceElement(e.geometry());
+ for (unsigned int j = 0; j < e.subEntities(dim); ++j) {
+ int k = cellType.permutation(j);
+ std::size_t idx = fct.ncomps() * indexSet.subIndex(e, k, dim);
+ for (int comp = 0; comp < fct.ncomps(); ++comp)
+ data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim)));
+ }
+ for (unsigned int j = 0; j < e.subEntities(dim-1); ++j) {
+ int k = cellType.permutation(e.subEntities(dim) + j);
+ std::size_t idx = fct.ncomps() * (indexSet.subIndex(e, k, dim-1) + gridView_.size(dim));
+ for (int comp = 0; comp < fct.ncomps(); ++comp)
+ data[idx + comp] = T(localFct.evaluate(comp, refElem.position(k, dim-1)));
+ }
+ localFct.unbind();
+ }
+ return data;
}
- localFct.unbind();
- }
- return data;
- }
-protected:
- using Super::gridView_;
-};
+ protected:
+ using Super::gridView_;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/spdatacollector.hh b/dune/vtk/datacollectors/spdatacollector.hh
index 25ab8be2ab1d222f8db359ac122e819ff712a010..0628d342950d0abbd419c3b175ca81c220fe8fd9 100644
--- a/dune/vtk/datacollectors/spdatacollector.hh
+++ b/dune/vtk/datacollectors/spdatacollector.hh
@@ -16,86 +16,90 @@
namespace Dune
{
-#if HAVE_DUNE_SPGRID
-
-// Specialization for SPGrid
-template <class GridView>
-class SPDataCollector
- : public StructuredDataCollectorInterface<GridView, SPDataCollector<GridView>>
-{
- using Self = SPDataCollector;
- using Super = StructuredDataCollectorInterface<GridView, Self>;
- using ctype = typename GridView::ctype;
-
-public:
- using Super::dim;
- using Super::partition;
-
-public:
- SPDataCollector (GridView const& gridView)
- : Super(gridView)
- , wholeExtent_(filledArray<6,int>(0))
- , extent_(filledArray<6,int>(0))
- , origin_(0.0)
- , spacing_(0.0)
- {}
-
- std::array<int, 6> const& wholeExtentImpl () const
- {
- return wholeExtent_;
- }
-
- std::array<int, 6> const& extentImpl () const
- {
- return extent_;
- }
-
- auto const& originImpl () const
+ namespace Vtk
{
- return origin_;
- }
-
- auto const& spacingImpl () const
- {
- return spacing_;
- }
-
- void updateImpl ()
- {
- Super::updateImpl();
-
- auto const& localMesh = gridView_.impl().gridLevel().localMesh();
- auto const& begin = gridView_.impl().gridLevel().globalMesh().begin();
- auto const& end = gridView_.impl().gridLevel().globalMesh().end();
- auto const& cube = gridView_.grid().domain().cube();
-
- for (int i = 0; i < dim; ++i) {
- wholeExtent_[2*i] = begin[i];
- wholeExtent_[2*i+1] = end[i];
- extent_[2*i] = localMesh.begin()[i];
- extent_[2*i+1] = localMesh.end()[i];
- spacing_[i] = cube.width()[i] / (end[i] - begin[i]);
- origin_[i] = cube.origin()[i];
+ #if HAVE_DUNE_SPGRID
+
+ // Specialization for SPGrid
+ template <class GridView>
+ class SPDataCollector
+ : public StructuredDataCollectorInterface<GridView, SPDataCollector<GridView>>
+ {
+ using Self = SPDataCollector;
+ using Super = StructuredDataCollectorInterface<GridView, Self>;
+ using ctype = typename GridView::ctype;
+
+ public:
+ using Super::dim;
+ using Super::partition;
+
+ public:
+ SPDataCollector (GridView const& gridView)
+ : Super(gridView)
+ , wholeExtent_(filledArray<6,int>(0))
+ , extent_(filledArray<6,int>(0))
+ , origin_(0.0)
+ , spacing_(0.0)
+ {}
+
+ std::array<int, 6> const& wholeExtentImpl () const
+ {
+ return wholeExtent_;
+ }
+
+ std::array<int, 6> const& extentImpl () const
+ {
+ return extent_;
+ }
+
+ auto const& originImpl () const
+ {
+ return origin_;
+ }
+
+ auto const& spacingImpl () const
+ {
+ return spacing_;
+ }
+
+ void updateImpl ()
+ {
+ Super::updateImpl();
+
+ auto const& localMesh = gridView_.impl().gridLevel().localMesh();
+ auto const& begin = gridView_.impl().gridLevel().globalMesh().begin();
+ auto const& end = gridView_.impl().gridLevel().globalMesh().end();
+ auto const& cube = gridView_.grid().domain().cube();
+
+ for (int i = 0; i < dim; ++i) {
+ wholeExtent_[2*i] = begin[i];
+ wholeExtent_[2*i+1] = end[i];
+ extent_[2*i] = localMesh.begin()[i];
+ extent_[2*i+1] = localMesh.end()[i];
+ spacing_[i] = cube.width()[i] / (end[i] - begin[i]);
+ origin_[i] = cube.origin()[i];
+ }
+ }
+
+ protected:
+ using Super::gridView_;
+ std::array<int, 6> wholeExtent_;
+ std::array<int, 6> extent_;
+ FieldVector<ctype,3> origin_;
+ FieldVector<ctype,3> spacing_;
+ std::vector<std::size_t> indexMap_;
+ };
+
+ namespace Impl
+ {
+ template <class GridView, class ct, int dim, template< int > class Ref, class Comm>
+ struct StructuredDataCollectorImpl<GridView, SPGrid<ct,dim,Ref,Comm>>
+ {
+ using type = SPDataCollector<GridView>;
+ };
}
- }
-
-protected:
- using Super::gridView_;
- std::array<int, 6> wholeExtent_;
- std::array<int, 6> extent_;
- FieldVector<ctype,3> origin_;
- FieldVector<ctype,3> spacing_;
- std::vector<std::size_t> indexMap_;
-};
-
-namespace Impl
-{
- template <class GridView, class ct, int dim, template< int > class Ref, class Comm>
- struct StructuredDataCollectorImpl<GridView, SPGrid<ct,dim,Ref,Comm>>
- {
- using type = SPDataCollector<GridView>;
- };
-}
-#endif // HAVE_DUNE_SPGRID
+ #endif // HAVE_DUNE_SPGRID
+
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/structureddatacollector.hh b/dune/vtk/datacollectors/structureddatacollector.hh
index 0b1d5ab5d56f42f7a5aaa932a192c1eb0cb44e6d..eb69931f439e00d357114b9da06f4ffc339ba0e5 100644
--- a/dune/vtk/datacollectors/structureddatacollector.hh
+++ b/dune/vtk/datacollectors/structureddatacollector.hh
@@ -12,225 +12,228 @@
namespace Dune
{
-/// The Interface for structured data-collectors
-template <class GridView, class Derived>
-class StructuredDataCollectorInterface
- : public DataCollectorInterface<GridView, Derived>
-{
-protected:
- using Super = DataCollectorInterface<GridView, Derived>;
- using SubDataCollector = ContinuousDataCollector<GridView>;
- using ctype = typename GridView::ctype;
-
-public:
- using Super::dim;
- using Super::partition;
-
-public:
- StructuredDataCollectorInterface (GridView const& gridView)
- : Super(gridView)
- , subDataCollector_(gridView)
- {}
-
- /// Sequence of Index pairs [begin, end) for the cells in each direction
- std::array<int, 6> wholeExtent () const
- {
- return this->asDerived().wholeExtentImpl();
- }
+ namespace Vtk
+ {
+ /// The Interface for structured data-collectors
+ template <class GridView, class Derived>
+ class StructuredDataCollectorInterface
+ : public DataCollectorInterface<GridView, Derived>
+ {
+ protected:
+ using Super = DataCollectorInterface<GridView, Derived>;
+ using SubDataCollector = ContinuousDataCollector<GridView>;
+ using ctype = typename GridView::ctype;
+
+ public:
+ using Super::dim;
+ using Super::partition;
+
+ public:
+ StructuredDataCollectorInterface (GridView const& gridView)
+ : Super(gridView)
+ , subDataCollector_(gridView)
+ {}
+
+ /// Sequence of Index pairs [begin, end) for the cells in each direction
+ std::array<int, 6> wholeExtent () const
+ {
+ return this->asDerived().wholeExtentImpl();
+ }
- /// Sequence of Index pairs [begin, end) for the cells in each direction of the local partition
- std::array<int, 6> extent () const
- {
- return this->asDerived().extentImpl();
- }
+ /// Sequence of Index pairs [begin, end) for the cells in each direction of the local partition
+ std::array<int, 6> extent () const
+ {
+ return this->asDerived().extentImpl();
+ }
- /// Call the `writer` with extent
- template <class Writer>
- void writeLocalPiece (Writer const& writer) const
- {
- this->asDerived().writeLocalPieceImpl(writer);
- }
+ /// Call the `writer` with extent
+ template <class Writer>
+ void writeLocalPiece (Writer const& writer) const
+ {
+ this->asDerived().writeLocalPieceImpl(writer);
+ }
- /// Call the `writer` with piece number and piece extent
- template <class Writer>
- void writePieces (Writer const& writer) const
- {
- this->asDerived().writePiecesImpl(writer);
- }
+ /// Call the `writer` with piece number and piece extent
+ template <class Writer>
+ void writePieces (Writer const& writer) const
+ {
+ this->asDerived().writePiecesImpl(writer);
+ }
- /// Interface for ImageData:
- /// @{
+ /// Interface for ImageData:
+ /// @{
- /// Lower left corner of the grid
- FieldVector<ctype, 3> origin () const
- {
- return this->asDerived().originImpl();
- }
+ /// Lower left corner of the grid
+ FieldVector<ctype, 3> origin () const
+ {
+ return this->asDerived().originImpl();
+ }
- /// Constant grid spacing in each coordinate direction
- FieldVector<ctype, 3> spacing () const
- {
- return this->asDerived().spacingImpl();
- }
+ /// Constant grid spacing in each coordinate direction
+ FieldVector<ctype, 3> spacing () const
+ {
+ return this->asDerived().spacingImpl();
+ }
- /// @}
+ /// @}
- /// Interface for RectilinearGrid
- /// @{
+ /// Interface for RectilinearGrid
+ /// @{
- /// The coordinates defines point coordinates for an extent by specifying the ordinate along each axis.
- template <class T>
- std::array<std::vector<T>, 3> coordinates () const
- {
- return this->asDerived().template coordinatesImpl<T>();
- }
+ /// The coordinates defines point coordinates for an extent by specifying the ordinate along each axis.
+ template <class T>
+ std::array<std::vector<T>, 3> coordinates () const
+ {
+ return this->asDerived().template coordinatesImpl<T>();
+ }
- /// @}
+ /// @}
-public: // default implementation:
+ public: // default implementation:
- /// \copyref DefaultDataCollector::update.
- void updateImpl ()
- {
- subDataCollector_.update();
+ /// \copyref DefaultDataCollector::update.
+ void updateImpl ()
+ {
+ subDataCollector_.update();
-#if HAVE_MPI
- int rank = gridView_.comm().rank();
- int numRanks = gridView_.comm().size();
+ #if HAVE_MPI
+ int rank = gridView_.comm().rank();
+ int numRanks = gridView_.comm().size();
- if (rank == 0) {
- extents_.resize(numRanks);
- requests_.resize(numRanks, MPI_REQUEST_NULL);
- for (int i = 1; i < numRanks; ++i)
- MPI_Irecv(extents_[i].data(), extents_[i].size(), MPI_INT, i, /*tag=*/6, gridView_.comm(), &requests_[i]);
- }
+ if (rank == 0) {
+ extents_.resize(numRanks);
+ requests_.resize(numRanks, MPI_REQUEST_NULL);
+ for (int i = 1; i < numRanks; ++i)
+ MPI_Irecv(extents_[i].data(), extents_[i].size(), MPI_INT, i, /*tag=*/6, gridView_.comm(), &requests_[i]);
+ }
- sendRequest_ = MPI_REQUEST_NULL;
-#endif
- }
+ sendRequest_ = MPI_REQUEST_NULL;
+ #endif
+ }
- /// Return number of grid cells
- std::uint64_t numCellsImpl () const
- {
- auto extent = this->extent();
- std::uint64_t num = 1;
- for (int d = 0; d < dim; ++d)
- num *= extent[2*d+1] - extent[2*d];
- return num;
- }
-
- /// Return number of grid vertices
- std::uint64_t numPointsImpl () const
- {
- return subDataCollector_.numPoints();
- }
+ /// Return number of grid cells
+ std::uint64_t numCellsImpl () const
+ {
+ auto extent = this->extent();
+ std::uint64_t num = 1;
+ for (int d = 0; d < dim; ++d)
+ num *= extent[2*d+1] - extent[2*d];
+ return num;
+ }
- /// \copydoc DefaultDataCollector::points.
- template <class T>
- std::vector<T> pointsImpl () const
- {
- return subDataCollector_.template points<T>();
- }
+ /// Return number of grid vertices
+ std::uint64_t numPointsImpl () const
+ {
+ return subDataCollector_.numPoints();
+ }
- /// \copydoc DefaultDataCollector::pointData
- template <class T, class GlobalFunction>
- std::vector<T> pointDataImpl (GlobalFunction const& fct) const
- {
- return subDataCollector_.template pointData<T>(fct);
- }
+ /// \copydoc DefaultDataCollector::points.
+ template <class T>
+ std::vector<T> pointsImpl () const
+ {
+ return subDataCollector_.template points<T>();
+ }
- // Calculates the extent and communicates it to rank 0.
- template <class Writer>
- void writeLocalPieceImpl (Writer const& writer) const
- {
- auto&& extent = this->extent();
+ /// \copydoc DefaultDataCollector::pointData
+ template <class T, class GlobalFunction>
+ std::vector<T> pointDataImpl (GlobalFunction const& fct) const
+ {
+ return subDataCollector_.template pointData<T>(fct);
+ }
-#if HAVE_MPI
- int sendFlag = 0;
- MPI_Status sendStatus;
- MPI_Test(&sendRequest_, &sendFlag, &sendStatus);
+ // Calculates the extent and communicates it to rank 0.
+ template <class Writer>
+ void writeLocalPieceImpl (Writer const& writer) const
+ {
+ auto&& extent = this->extent();
+
+ #if HAVE_MPI
+ int sendFlag = 0;
+ MPI_Status sendStatus;
+ MPI_Test(&sendRequest_, &sendFlag, &sendStatus);
+
+ if (sendFlag) {
+ int rank = gridView_.comm().rank();
+ if (rank != 0) {
+ MPI_Isend(extent.data(), extent.size(), MPI_INT, 0, /*tag=*/6, gridView_.comm(), &sendRequest_);
+ } else {
+ extents_[0] = extent;
+ }
+ }
+ #endif
+
+ writer(extent);
+ }
- if (sendFlag) {
- int rank = gridView_.comm().rank();
- if (rank != 0) {
- MPI_Isend(extent.data(), extent.size(), MPI_INT, 0, /*tag=*/6, gridView_.comm(), &sendRequest_);
- } else {
- extents_[0] = extent;
+ // Receive extent from all ranks and call the `writer` with the rank's extent vector
+ template <class Writer>
+ void writePiecesImpl (Writer const& writer) const
+ {
+ #if HAVE_MPI
+ writer(0, extents_[0], true);
+
+ int numRanks = gridView_.comm().size();
+ for (int p = 1; p < numRanks; ++p) {
+ int idx = -1;
+ MPI_Status status;
+ MPI_Waitany(numRanks, requests_.data(), &idx, &status);
+ if (idx != MPI_UNDEFINED) {
+ assert(idx == status.MPI_SOURCE && status.MPI_TAG == 6);
+ writer(idx, extents_[idx], true);
+ }
+ }
+ #else
+ writer(0, this->extent(), true);
+ #endif
}
- }
-#endif
- writer(extent);
- }
+ // Origin (0,0,0)
+ FieldVector<ctype, 3> originImpl () const
+ {
+ FieldVector<ctype, 3> vec; vec = ctype(0);
+ return vec;
+ }
- // Receive extent from all ranks and call the `writer` with the rank's extent vector
- template <class Writer>
- void writePiecesImpl (Writer const& writer) const
- {
-#if HAVE_MPI
- writer(0, extents_[0], true);
-
- int numRanks = gridView_.comm().size();
- for (int p = 1; p < numRanks; ++p) {
- int idx = -1;
- MPI_Status status;
- MPI_Waitany(numRanks, requests_.data(), &idx, &status);
- if (idx != MPI_UNDEFINED) {
- assert(idx == status.MPI_SOURCE && status.MPI_TAG == 6);
- writer(idx, extents_[idx], true);
- }
- }
-#else
- writer(0, this->extent(), true);
-#endif
- }
-
- // Origin (0,0,0)
- FieldVector<ctype, 3> originImpl () const
- {
- FieldVector<ctype, 3> vec; vec = ctype(0);
- return vec;
- }
+ // Grid spacing (0,0,0)
+ FieldVector<ctype, 3> spacingImpl () const
+ {
+ FieldVector<ctype, 3> vec; vec = ctype(0);
+ return vec;
+ }
- // Grid spacing (0,0,0)
- FieldVector<ctype, 3> spacingImpl () const
- {
- FieldVector<ctype, 3> vec; vec = ctype(0);
- return vec;
- }
+ // Ordinate along each axis with constant \ref spacing from the \ref origin
+ template <class T>
+ std::array<std::vector<T>, 3> coordinatesImpl () const
+ {
+ auto origin = this->origin();
+ auto spacing = this->spacing();
+ auto extent = this->extent();
+
+ std::array<std::vector<T>, 3> ordinates{};
+ for (int d = 0; d < dim; ++d) {
+ auto s = extent[2*d+1] - extent[2*d] + 1;
+ ordinates[d].resize(s);
+ for (int i = 0; i < s; ++i)
+ ordinates[d][i] = origin[d] + (extent[2*d] + i)*spacing[d];
+ }
+
+ for (int d = dim; d < 3; ++d)
+ ordinates[d].resize(1, T(0));
+
+ return ordinates;
+ }
- // Ordinate along each axis with constant \ref spacing from the \ref origin
- template <class T>
- std::array<std::vector<T>, 3> coordinatesImpl () const
- {
- auto origin = this->origin();
- auto spacing = this->spacing();
- auto extent = this->extent();
-
- std::array<std::vector<T>, 3> ordinates{};
- for (int d = 0; d < dim; ++d) {
- auto s = extent[2*d+1] - extent[2*d] + 1;
- ordinates[d].resize(s);
- for (int i = 0; i < s; ++i)
- ordinates[d][i] = origin[d] + (extent[2*d] + i)*spacing[d];
- }
-
- for (int d = dim; d < 3; ++d)
- ordinates[d].resize(1, T(0));
-
- return ordinates;
- }
-
-protected:
- using Super::gridView_;
- SubDataCollector subDataCollector_;
-
-#if HAVE_MPI
- mutable std::vector<std::array<int,6>> extents_;
- mutable std::vector<MPI_Request> requests_;
- mutable MPI_Request sendRequest_ = MPI_REQUEST_NULL;
-#endif
-};
+ protected:
+ using Super::gridView_;
+ SubDataCollector subDataCollector_;
+
+ #if HAVE_MPI
+ mutable std::vector<std::array<int,6>> extents_;
+ mutable std::vector<MPI_Request> requests_;
+ mutable MPI_Request sendRequest_ = MPI_REQUEST_NULL;
+ #endif
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc b/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc
index 99d22c5fca646eebb7005f5d2f99bb847790ee6b..32d6870db952f757fd1aecd6cee108ba1cc3e625 100644
--- a/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc
+++ b/dune/vtk/datacollectors/test/test-lagrangedatacollector.cc
@@ -18,7 +18,7 @@
#include <dune/vtk/datacollectors/lagrangedatacollector.hh>
/**
- * This test compares LagrangeDataCollector<1> against the ContinuousDataCollector, thus
+ * This test compares LagrangeDataCollector<1> against the ContinuousDataCollector, thus
* just compares the linear order implementations
**/
@@ -34,7 +34,7 @@ void testDataCollector (std::string prefix, Dune::TestSuite& testSuite, DataColl
// check that number of points and cells are equal
testSuite.check(dataCollector1.numPoints() == dataCollector2.numPoints(), prefix + "_numPoints");
testSuite.check(dataCollector1.numCells() == dataCollector2.numCells(), prefix + "_numCells");
-
+
auto points1 = dataCollector1.template points<double>();
auto points2 = dataCollector2.template points<double>();
@@ -71,14 +71,14 @@ void testGridView (std::string prefix, Dune::TestSuite& testSuite, GridView cons
{
// 1. test linear order lagrange data-collector
{
- Dune::ContinuousDataCollector<GridView> linearDataCollector(gridView);
-
+ Dune::Vtk::ContinuousDataCollector<GridView> linearDataCollector(gridView);
+
// data collector with template order
- Dune::LagrangeDataCollector<GridView, 1> lagrangeDataCollector1a(gridView);
+ Dune::Vtk::LagrangeDataCollector<GridView, 1> lagrangeDataCollector1a(gridView);
testDataCollector(prefix + "_linear_template", testSuite, lagrangeDataCollector1a, linearDataCollector);
// data collector with runtime order
- Dune::LagrangeDataCollector<GridView> lagrangeDataCollector1b(gridView, 1);
+ Dune::Vtk::LagrangeDataCollector<GridView> lagrangeDataCollector1b(gridView, 1);
testDataCollector(prefix + "_linear_runtime", testSuite, lagrangeDataCollector1b, linearDataCollector);
}
}
diff --git a/dune/vtk/datacollectors/unstructureddatacollector.hh b/dune/vtk/datacollectors/unstructureddatacollector.hh
index b5ff27a5dc786649ecb719887e3d08025b369405..070177a8a207f5a67e7c92cee03fc5d4960383b3 100644
--- a/dune/vtk/datacollectors/unstructureddatacollector.hh
+++ b/dune/vtk/datacollectors/unstructureddatacollector.hh
@@ -6,50 +6,53 @@
#include <dune/vtk/forward.hh>
#include <dune/vtk/datacollectorinterface.hh>
-namespace Dune {
-
-struct Cells
-{
- std::vector<std::uint8_t> types;
- std::vector<std::int64_t> offsets;
- std::vector<std::int64_t> connectivity;
-};
-
-template <class GridView, class Derived, class Partition>
-class UnstructuredDataCollectorInterface
- : public DataCollectorInterface<GridView, Derived, Partition>
+namespace Dune
{
- using Super = DataCollectorInterface<GridView, Derived, Partition>;
-
-public:
- using Super::dim;
- using Super::partition;
-
-public:
- UnstructuredDataCollectorInterface (GridView const& gridView)
- : Super(gridView)
- {}
-
- /// \brief Return cell types, offsets, and connectivity. \see Cells
- Cells cells () const
+ namespace Vtk
{
- return this->asDerived().cellsImpl();
- }
-
- std::vector<std::uint64_t> pointIds () const
- {
- return this->asDerived().pointIdsImpl();
- }
-
-protected:
- // default implementation
- std::vector<std::uint64_t> pointIdsImpl () const
- {
- return {};
- }
-
-protected:
- using Super::gridView_;
-};
-
+ struct Cells
+ {
+ std::vector<std::uint8_t> types;
+ std::vector<std::int64_t> offsets;
+ std::vector<std::int64_t> connectivity;
+ };
+
+ template <class GridView, class Derived, class Partition>
+ class UnstructuredDataCollectorInterface
+ : public DataCollectorInterface<GridView, Derived, Partition>
+ {
+ using Super = DataCollectorInterface<GridView, Derived, Partition>;
+
+ public:
+ using Super::dim;
+ using Super::partition;
+
+ public:
+ UnstructuredDataCollectorInterface (GridView const& gridView)
+ : Super(gridView)
+ {}
+
+ /// \brief Return cell types, offsets, and connectivity. \see Cells
+ Cells cells () const
+ {
+ return this->asDerived().cellsImpl();
+ }
+
+ std::vector<std::uint64_t> pointIds () const
+ {
+ return this->asDerived().pointIdsImpl();
+ }
+
+ protected:
+ // default implementation
+ std::vector<std::uint64_t> pointIdsImpl () const
+ {
+ return {};
+ }
+
+ protected:
+ using Super::gridView_;
+ };
+
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/datacollectors/yaspdatacollector.hh b/dune/vtk/datacollectors/yaspdatacollector.hh
index 818f200f1a5410295ac198d8ae941acee681a203..196c921630f9bc34060f891a43ac5d06d8a0d1b7 100644
--- a/dune/vtk/datacollectors/yaspdatacollector.hh
+++ b/dune/vtk/datacollectors/yaspdatacollector.hh
@@ -17,153 +17,156 @@
namespace Dune
{
-// Specialization for YaspGrid
-template <class GridView>
-class YaspDataCollector
- : public StructuredDataCollectorInterface<GridView, YaspDataCollector<GridView>>
-{
- using Self = YaspDataCollector;
- using Super = StructuredDataCollectorInterface<GridView, Self>;
- using ctype = typename GridView::ctype;
-
-public:
- using Super::dim;
- using Super::partition;
-
-public:
- YaspDataCollector (GridView const& gridView)
- : Super(gridView)
- , wholeExtent_(filledArray<6,int>(0))
- , extent_(filledArray<6,int>(0))
- , origin_(0.0)
- , spacing_(0.0)
- , level_(0)
- {}
-
- std::array<int, 6> const& wholeExtentImpl () const
- {
- return wholeExtent_;
- }
-
- std::array<int, 6> const& extentImpl () const
- {
- return extent_;
- }
-
- auto const& originImpl () const
- {
- return origin_;
- }
-
- auto const& spacingImpl () const
- {
- return spacing_;
- }
-
- void updateImpl ()
- {
- Super::updateImpl();
-
- level_ = gridView_.template begin<0,All_Partition>()->level();
- for (int i = 0; i < dim; ++i) {
- wholeExtent_[2*i] = 0;
- wholeExtent_[2*i+1] = grid(gridView_.grid()).levelSize(level_,i);
- }
-
- auto const& gl = *grid(gridView_.grid()).begin(level_);
- auto const& g = gl.interior[0];
- auto const& gc = *g.dataBegin();
- for (int i = 0; i < dim; ++i) {
- extent_[2*i] = gc.min(i);
- extent_[2*i+1] = gc.max(i)+1;
- }
-
- auto it = grid(gridView_.grid()).begin(level_);
- initGeometry(it->coords);
- }
-
- void initGeometry (EquidistantCoordinates<ctype,dim> const& coords)
- {
- for (int i = 0; i < dim; ++i) {
- spacing_[i] = coords.meshsize(i,0);
- origin_[i] = 0;
- }
- }
-
- void initGeometry (EquidistantOffsetCoordinates<ctype,dim> const& coords)
+ namespace Vtk
{
- for (int i = 0; i < dim; ++i) {
- spacing_[i] = coords.meshsize(i,0);
- origin_[i] = coords.origin(i);
+ // Specialization for YaspGrid
+ template <class GridView>
+ class YaspDataCollector
+ : public StructuredDataCollectorInterface<GridView, YaspDataCollector<GridView>>
+ {
+ using Self = YaspDataCollector;
+ using Super = StructuredDataCollectorInterface<GridView, Self>;
+ using ctype = typename GridView::ctype;
+
+ public:
+ using Super::dim;
+ using Super::partition;
+
+ public:
+ YaspDataCollector (GridView const& gridView)
+ : Super(gridView)
+ , wholeExtent_(filledArray<6,int>(0))
+ , extent_(filledArray<6,int>(0))
+ , origin_(0.0)
+ , spacing_(0.0)
+ , level_(0)
+ {}
+
+ std::array<int, 6> const& wholeExtentImpl () const
+ {
+ return wholeExtent_;
+ }
+
+ std::array<int, 6> const& extentImpl () const
+ {
+ return extent_;
+ }
+
+ auto const& originImpl () const
+ {
+ return origin_;
+ }
+
+ auto const& spacingImpl () const
+ {
+ return spacing_;
+ }
+
+ void updateImpl ()
+ {
+ Super::updateImpl();
+
+ level_ = gridView_.template begin<0,All_Partition>()->level();
+ for (int i = 0; i < dim; ++i) {
+ wholeExtent_[2*i] = 0;
+ wholeExtent_[2*i+1] = grid(gridView_.grid()).levelSize(level_,i);
+ }
+
+ auto const& gl = *grid(gridView_.grid()).begin(level_);
+ auto const& g = gl.interior[0];
+ auto const& gc = *g.dataBegin();
+ for (int i = 0; i < dim; ++i) {
+ extent_[2*i] = gc.min(i);
+ extent_[2*i+1] = gc.max(i)+1;
+ }
+
+ auto it = grid(gridView_.grid()).begin(level_);
+ initGeometry(it->coords);
+ }
+
+ void initGeometry (EquidistantCoordinates<ctype,dim> const& coords)
+ {
+ for (int i = 0; i < dim; ++i) {
+ spacing_[i] = coords.meshsize(i,0);
+ origin_[i] = 0;
+ }
+ }
+
+ void initGeometry (EquidistantOffsetCoordinates<ctype,dim> const& coords)
+ {
+ for (int i = 0; i < dim; ++i) {
+ spacing_[i] = coords.meshsize(i,0);
+ origin_[i] = coords.origin(i);
+ }
+ }
+
+ void initGeometry (TensorProductCoordinates<ctype,dim> const& coords)
+ {
+ for (int i = 0; i < dim; ++i) {
+ spacing_[i] = coords.meshsize(i,0); // is not constant, but also not used.
+ origin_[i] = coords.coordinate(i,0);
+ }
+ }
+
+
+ /// Extract the ordinates from the coordinates object of the current level
+ template <class T>
+ std::array<std::vector<T>, 3> coordinatesImpl () const
+ {
+ auto it = grid(gridView_.grid()).begin(level_);
+ auto const& coords = it->coords;
+
+ std::array<std::vector<T>, 3> ordinates{};
+ for (int d = 0; d < dim; ++d) {
+ auto s = extent_[2*d+1] - extent_[2*d] + 1;
+ ordinates[d].resize(s);
+ for (int i = 0; i < s; ++i)
+ ordinates[d][i] = coords.coordinate(d, extent_[2*d] + i);
+ }
+
+ for (int d = dim; d < 3; ++d)
+ ordinates[d].resize(1, T(0));
+
+ return ordinates;
+ }
+
+
+ private:
+
+ template <class G>
+ using HostGrid = decltype(std::declval<G>().hostGrid());
+
+ template <class G,
+ std::enable_if_t<not Std::is_detected<HostGrid, G>::value, int> = 0>
+ auto const& grid (G const& g) const
+ {
+ return g;
+ }
+
+ template <class G,
+ std::enable_if_t<Std::is_detected<HostGrid, G>::value, int> = 0>
+ auto const& grid (G const& g) const
+ {
+ return grid(g.hostGrid());
+ }
+
+ protected:
+ using Super::gridView_;
+ std::array<int, 6> wholeExtent_;
+ std::array<int, 6> extent_;
+ FieldVector<ctype,3> origin_;
+ FieldVector<ctype,3> spacing_;
+ int level_;
+ };
+
+ namespace Impl
+ {
+ template <class GridView, int dim, class Coordinates>
+ struct StructuredDataCollectorImpl<GridView, YaspGrid<dim,Coordinates>>
+ {
+ using type = YaspDataCollector<GridView>;
+ };
}
- }
-
- void initGeometry (TensorProductCoordinates<ctype,dim> const& coords)
- {
- for (int i = 0; i < dim; ++i) {
- spacing_[i] = coords.meshsize(i,0); // is not constant, but also not used.
- origin_[i] = coords.coordinate(i,0);
- }
- }
-
-
- /// Extract the ordinates from the coordinates object of the current level
- template <class T>
- std::array<std::vector<T>, 3> coordinatesImpl () const
- {
- auto it = grid(gridView_.grid()).begin(level_);
- auto const& coords = it->coords;
-
- std::array<std::vector<T>, 3> ordinates{};
- for (int d = 0; d < dim; ++d) {
- auto s = extent_[2*d+1] - extent_[2*d] + 1;
- ordinates[d].resize(s);
- for (int i = 0; i < s; ++i)
- ordinates[d][i] = coords.coordinate(d, extent_[2*d] + i);
- }
-
- for (int d = dim; d < 3; ++d)
- ordinates[d].resize(1, T(0));
-
- return ordinates;
- }
-
-
-private:
-
- template <class G>
- using HostGrid = decltype(std::declval<G>().hostGrid());
-
- template <class G,
- std::enable_if_t<not Std::is_detected<HostGrid, G>::value, int> = 0>
- auto const& grid (G const& g) const
- {
- return g;
- }
-
- template <class G,
- std::enable_if_t<Std::is_detected<HostGrid, G>::value, int> = 0>
- auto const& grid (G const& g) const
- {
- return grid(g.hostGrid());
- }
-
-protected:
- using Super::gridView_;
- std::array<int, 6> wholeExtent_;
- std::array<int, 6> extent_;
- FieldVector<ctype,3> origin_;
- FieldVector<ctype,3> spacing_;
- int level_;
-};
-
-namespace Impl
-{
- template <class GridView, int dim, class Coordinates>
- struct StructuredDataCollectorImpl<GridView, YaspGrid<dim,Coordinates>>
- {
- using type = YaspDataCollector<GridView>;
- };
-}
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/defaultvtkfunction.hh b/dune/vtk/defaultvtkfunction.hh
index 32a02225e906ef974bdcd37789d4be0ed2fdb5ed..2dd750bf9167748c21be6a25ad0d7c28e0c31403 100644
--- a/dune/vtk/defaultvtkfunction.hh
+++ b/dune/vtk/defaultvtkfunction.hh
@@ -5,76 +5,79 @@
#include <dune/common/fmatrix.hh>
#include <dune/common/fvector.hh>
-#include "vtklocalfunctioninterface.hh"
+#include "localfunctioninterface.hh"
namespace Dune
{
- /// Type erasure for dune-functions LocalFunction interface
- template <class GridView, class LocalFunction>
- class LocalFunctionWrapper final
- : public VtkLocalFunctionInterface<GridView>
+ namespace Vtk
{
- using Self = LocalFunctionWrapper;
- using Interface = VtkLocalFunctionInterface<GridView>;
- using Entity = typename Interface::Entity;
- using LocalCoordinate = typename Interface::LocalCoordinate;
+ /// Type erasure for dune-functions LocalFunction interface
+ template <class GridView, class LocalFunction>
+ class LocalFunctionWrapper final
+ : public LocalFunctionInterface<GridView>
+ {
+ using Self = LocalFunctionWrapper;
+ using Interface = LocalFunctionInterface<GridView>;
+ using Entity = typename Interface::Entity;
+ using LocalCoordinate = typename Interface::LocalCoordinate;
- template <class F, class D>
- using Range = std::decay_t<decltype(std::declval<F>()(std::declval<D>()))>;
+ template <class F, class D>
+ using Range = std::decay_t<decltype(std::declval<F>()(std::declval<D>()))>;
- public:
- /// Constructor. Stores a copy of the passed `localFct` in a local variable.
- template <class LocalFct,
- disableCopyMove<Self, LocalFct> = 0>
- LocalFunctionWrapper (LocalFct&& localFct)
- : localFct_(std::forward<LocalFct>(localFct))
- {}
+ public:
+ /// Constructor. Stores a copy of the passed `localFct` in a local variable.
+ template <class LocalFct,
+ disableCopyMove<Self, LocalFct> = 0>
+ LocalFunctionWrapper (LocalFct&& localFct)
+ : localFct_(std::forward<LocalFct>(localFct))
+ {}
- /// Bind the LocalFunction to the Entity
- virtual void bind (Entity const& entity) override
- {
- localFct_.bind(entity);
- }
+ /// Bind the LocalFunction to the Entity
+ virtual void bind (Entity const& entity) override
+ {
+ localFct_.bind(entity);
+ }
- /// Unbind the LocalFunction from the Entity
- virtual void unbind () override
- {
- localFct_.unbind();
- }
+ /// Unbind the LocalFunction from the Entity
+ virtual void unbind () override
+ {
+ localFct_.unbind();
+ }
- /// Evaluate the LocalFunction in LocalCoordinates
- virtual double evaluate (int comp, LocalCoordinate const& xi) const override
- {
- return evaluateImpl(comp, localFct_(xi));
- }
+ /// Evaluate the LocalFunction in LocalCoordinates
+ virtual double evaluate (int comp, LocalCoordinate const& xi) const override
+ {
+ return evaluateImpl(comp, localFct_(xi));
+ }
- private:
- // Evaluate a component of a vector valued data
- template <class T, int N, int M>
- double evaluateImpl (int comp, FieldMatrix<T,N,M> const& mat) const
- {
- int r = comp / 3;
- int c = comp % 3;
- return r < N && c < M ? mat[r][c] : 0.0;
- }
+ private:
+ // Evaluate a component of a vector valued data
+ template <class T, int N, int M>
+ double evaluateImpl (int comp, FieldMatrix<T,N,M> const& mat) const
+ {
+ int r = comp / 3;
+ int c = comp % 3;
+ return r < N && c < M ? mat[r][c] : 0.0;
+ }
- // Evaluate a component of a vector valued data
- template <class T, int N>
- double evaluateImpl (int comp, FieldVector<T,N> const& vec) const
- {
- return comp < N ? vec[comp] : 0.0;
- }
+ // Evaluate a component of a vector valued data
+ template <class T, int N>
+ double evaluateImpl (int comp, FieldVector<T,N> const& vec) const
+ {
+ return comp < N ? vec[comp] : 0.0;
+ }
- // Return the scalar values
- template <class T>
- double evaluateImpl (int comp, T const& value) const
- {
- assert(comp == 0);
- return value;
- }
+ // Return the scalar values
+ template <class T>
+ double evaluateImpl (int comp, T const& value) const
+ {
+ assert(comp == 0);
+ return value;
+ }
- private:
- LocalFunction localFct_;
- };
+ private:
+ LocalFunction localFct_;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/filereader.hh b/dune/vtk/filereader.hh
index 77c40b978f15d873dcf6896073c0b8cc24c07596..8f021d7ae08c9c021f0c6b0a1bc52b7c364cfce8 100644
--- a/dune/vtk/filereader.hh
+++ b/dune/vtk/filereader.hh
@@ -10,66 +10,73 @@
namespace Dune
{
- template <class Grid, class FilerReaderImp>
- class FileReader
+ namespace Vtk
{
- private:
- // type of underlying implementation, for internal use only
- using Implementation = FilerReaderImp;
-
- /// \brief An accessor class to call protected members of reader implementations.
- struct Accessor : public Implementation
+ template <class Grid, class FilerReaderImp>
+ class FileReader
{
+ private:
+ // type of underlying implementation, for internal use only
+ using Implementation = FilerReaderImp;
+
+ /// \brief An accessor class to call protected members of reader implementations.
+ struct Accessor : public Implementation
+ {
+ template <class... Args>
+ static std::unique_ptr<Grid> createGridFromFileImpl (Args&&... args)
+ {
+ return Implementation::createGridFromFileImpl(std::forward<Args>(args)...);
+ }
+
+ template <class... Args>
+ static void fillFactoryImpl (Args&&... args)
+ {
+ return Implementation::fillFactoryImpl(std::forward<Args>(args)...);
+ }
+ };
+
+ public:
+ /// Reads the grid from a file with filename and returns a unique_ptr to the created grid.
+ /// Redirects to concrete implementation of derivated class.
template <class... Args>
- static std::unique_ptr<Grid> readImp (Args&&... args)
+ static std::unique_ptr<Grid> createGridFromFile (const std::string &filename, Args&&... args)
{
- return Implementation::readImpl(std::forward<Args>(args)...);
+ return Accessor::createGridFromFileImpl(filename, std::forward<Args>(args)...);
}
+ /// Reads the grid from a file with filename into a grid-factory.
+ /// Redirects to concrete implementation of derivated class.
template <class... Args>
- static void readFactoryImpl (Args&&... args)
+ static void fillFactory (GridFactory<Grid> &factory,
+ const std::string &filename,
+ Args&&... args)
{
- return Implementation::readFactoryImpl(std::forward<Args>(args)...);
+ Accessor::fillFactoryImpl(factory, filename, std::forward<Args>(args)...);
}
- };
-
- public:
- /// Reads the grid from a file with filename and returns a unique_ptr to the created grid.
- /// Redirects to concrete implementation of derivated class.
- template <class... Args>
- static std::unique_ptr<Grid> read (const std::string &filename, Args&&... args)
- {
- return Accessor::readImpl(filename, std::forward<Args>(args)...);
- }
-
- /// Reads the grid from a file with filename into a grid-factory.
- /// Redirects to concrete implementation of derivated class.
- template <class... Args>
- static void read (GridFactory<Grid> &factory, const std::string &filename, Args&&... args)
- {
- Accessor::readFactoryImpl(factory, filename, std::forward<Args>(args)...);
- }
- protected: // default implementations
+ protected: // default implementations
- // Default implementation, redirects to factory read implementation.
- template <class... Args>
- static std::unique_ptr<Grid> readImpl (const std::string &filename, Args&&... args)
- {
- GridFactory<Grid> factory;
- read(factory, filename, std::forward<Args>(args)...);
+ // Default implementation, redirects to fillFactory implementation.
+ template <class... Args>
+ static std::unique_ptr<Grid> createGridFromFileImpl (const std::string &filename,
+ Args&&... args)
+ {
+ GridFactory<Grid> factory;
+ fillFactory(factory, filename, std::forward<Args>(args)...);
- return std::unique_ptr<Grid>{ factory.createGrid() };
- }
+ return std::unique_ptr<Grid>{ factory.createGrid() };
+ }
- // Default implementation for reading into grid-factory: produces a runtime-error.
- template <class... Args>
- static void readFactoryImpl (GridFactory<Grid> &/*factory*/, const std::string &/*filename*/,
- Args&&... /*args*/)
- {
- DUNE_THROW(NotImplemented,
- "GridReader using a factory argument not implemented for concrete reader implementation.");
- }
- };
+ // Default implementation for reading into grid-factory: produces a runtime-error.
+ template <class... Args>
+ static void fillFactoryImpl (GridFactory<Grid> &/*factory*/,
+ const std::string &/*filename*/,
+ Args&&... /*args*/)
+ {
+ DUNE_THROW(NotImplemented,
+ "GridReader using a factory argument not implemented for concrete reader implementation.");
+ }
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/filewriter.hh b/dune/vtk/filewriter.hh
index 803fe2566dc5a2bc2d93183e7145508497854ef0..840df5dc53f4866f496e6cacfded5cad0ac04cf0 100644
--- a/dune/vtk/filewriter.hh
+++ b/dune/vtk/filewriter.hh
@@ -7,14 +7,17 @@
namespace Dune
{
- class FileWriter
+ namespace Vtk
{
- public:
- /// Virtual destructor
- virtual ~FileWriter () = default;
+ class FileWriter
+ {
+ public:
+ /// Virtual destructor
+ virtual ~FileWriter () = default;
- /// Write to file given by `filename` and (optionally) store additional data in `dataDir`
- virtual std::string write (std::string const& filename, std::optional<std::string> dataDir = {}) const = 0;
- };
+ /// Write to file given by `filename` and (optionally) store additional data in `dataDir`
+ virtual std::string write (std::string const& filename, std::optional<std::string> dataDir = {}) const = 0;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/forward.hh b/dune/vtk/forward.hh
index f3946eab572f9d9362f38a34d9e480788c09e838..eb90b64109807ce38f20c6b05d55632f17e14a4c 100644
--- a/dune/vtk/forward.hh
+++ b/dune/vtk/forward.hh
@@ -4,81 +4,81 @@
namespace Dune
{
- // forward declaration of all classes in dune-vtk
-
- template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
- class DataCollectorInterface;
+ namespace Vtk
+ {
+ // forward declaration of all classes in dune-vtk
- // @{ datacollectors
- template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
- class UnstructuredDataCollectorInterface;
+ template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
+ class DataCollectorInterface;
- // @{ unstructured-datacollectors
- template <class GridView, class Partition = Partitions::InteriorBorder>
- class ContinuousDataCollector;
+ // @{ datacollectors
+ template <class GridView, class Derived, class Partition = Partitions::InteriorBorder>
+ class UnstructuredDataCollectorInterface;
- template <class GridView, class Partition = Partitions::InteriorBorder>
- class DiscontinuousDataCollector;
+ // @{ unstructured-datacollectors
+ template <class GridView, class Partition = Partitions::InteriorBorder>
+ class ContinuousDataCollector;
- template <class GridView>
- class QuadraticDataCollector;
- // @} unstructured-datacollectors
+ template <class GridView, class Partition = Partitions::InteriorBorder>
+ class DiscontinuousDataCollector;
- template <class GridView, class Derived>
- class StructuredDataCollectorInterface;
+ template <class GridView>
+ class QuadraticDataCollector;
+ // @} unstructured-datacollectors
- namespace Impl
- {
- // Should be specialized for concrete structured grid
- template <class GridView, class Grid>
- struct StructuredDataCollectorImpl;
- }
+ template <class GridView, class Derived>
+ class StructuredDataCollectorInterface;
- template <class GridView>
- using StructuredDataCollector = typename Impl::StructuredDataCollectorImpl<GridView, typename GridView::Grid>::type;
+ namespace Impl
+ {
+ // Should be specialized for concrete structured grid
+ template <class GridView, class Grid>
+ struct StructuredDataCollectorImpl;
+ }
- // @{ structured-datacollectors
- template <class GridView>
- class SPDataCollector;
+ template <class GridView>
+ using StructuredDataCollector = typename Impl::StructuredDataCollectorImpl<GridView, typename GridView::Grid>::type;
- template <class GridView>
- class YaspDataCollector;
- // @} structured-datacollectors
+ // @{ structured-datacollectors
+ template <class GridView>
+ class SPDataCollector;
- // @} datacollectors
+ template <class GridView>
+ class YaspDataCollector;
+ // @} structured-datacollectors
- template <class Grid, class Derived>
- class GridCreatorInterface;
+ // @} datacollectors
- template <class GridCreator, class Derived>
- struct DerivedGridCreator;
+ template <class Grid, class Derived>
+ class GridCreatorInterface;
- // @{ gridcreators
- template <class Grid>
- struct ContinuousGridCreator;
+ template <class GridCreator, class Derived>
+ struct DerivedGridCreator;
- template <class Grid>
- struct DiscontinuousGridCreator;
+ // @{ gridcreators
+ template <class Grid>
+ struct ContinuousGridCreator;
- template <class Grid>
- struct ParallelGridCreator;
+ template <class Grid>
+ struct DiscontinuousGridCreator;
- template <class Grid>
- struct SerialGridCreator;
+ template <class Grid>
+ struct ParallelGridCreator;
- template <class Grid>
- class LagrangeGridCreator;
- // @} gridcreators
+ template <class Grid>
+ struct SerialGridCreator;
+ template <class Grid>
+ class LagrangeGridCreator;
+ // @} gridcreators
- template <class Grid, class FilerReaderImp>
- class FileReader;
- template <class Grid, class GridCreator = ContinuousGridCreator<Grid>>
- class VtkReader;
+ template <class Grid, class FilerReaderImp>
+ class FileReader;
+ class FileWriter;
- class FileWriter;
+ } //end namespace Vtk
// @{ filewriters
template <class VtkWriter>
@@ -91,19 +91,22 @@ namespace Dune
class VtkWriterInterface;
// @{ vtkwriters
- template <class GridView, class DataCollector = StructuredDataCollector<GridView>>
+ template <class GridView, class DataCollector = Vtk::StructuredDataCollector<GridView>>
class VtkImageDataWriter;
- template <class GridView, class DataCollector = StructuredDataCollector<GridView>>
+ template <class GridView, class DataCollector = Vtk::StructuredDataCollector<GridView>>
class VtkRectilinearGridWriter;
- template <class GridView, class DataCollector = StructuredDataCollector<GridView>>
+ template <class GridView, class DataCollector = Vtk::StructuredDataCollector<GridView>>
class VtkStructuredGridWriter;
- template <class GridView, class DataCollector = ContinuousDataCollector<GridView>>
+ template <class GridView, class DataCollector = Vtk::ContinuousDataCollector<GridView>>
class VtkUnstructuredGridWriter;
// @} vtkwriters
// @} filewriters
+ template <class Grid, class GridCreator = Vtk::ContinuousGridCreator<Grid>>
+ class VtkReader;
+
} // end namespace Dune
diff --git a/dune/vtk/function.hh b/dune/vtk/function.hh
new file mode 100644
index 0000000000000000000000000000000000000000..0e99387d777b65e839135a2c474d5614ffc4457f
--- /dev/null
+++ b/dune/vtk/function.hh
@@ -0,0 +1,131 @@
+#pragma once
+
+#include <optional>
+#include <type_traits>
+
+#include <dune/common/std/type_traits.hh>
+
+#include "localfunction.hh"
+#include "types.hh"
+
+namespace Dune
+{
+ // forward declarations
+ template <class T, int N>
+ class FieldVector;
+
+ template <class T, int N, int M>
+ class FieldMatrix;
+
+ namespace Vtk
+ {
+ /// Wrapper class for functions allowing local evaluations.
+ template <class GridView>
+ class Function
+ {
+ template <class F>
+ using LocalFunction = decltype(localFunction(std::declval<F>()));
+
+ using Domain = typename GridView::template Codim<0>::Entity::Geometry::LocalCoordinate;
+
+ template <class F>
+ using Range = std::decay_t<std::result_of_t<F(Domain)>>;
+
+ private:
+
+ template <class T, int N>
+ static auto sizeOfImpl (FieldVector<T,N> const&)
+ -> std::integral_constant<int, N> { return {}; }
+
+ template <class T, int N, int M>
+ static auto sizeOfImpl (FieldMatrix<T,N,M> const&)
+ -> std::integral_constant<int, N*M> { return {}; }
+
+ static auto sizeOfImpl (...)
+ -> std::integral_constant<int, 1> { return {}; }
+
+ template <class T>
+ static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }
+
+ public:
+ /// Constructor VtkFunction from legacy VTKFunction
+ /**
+ * \param fct The VTKFunction to wrap
+ * \param type The VTK datatype how to write the function values to the output [Vtk::FLOAT64]
+ **/
+ Function (std::shared_ptr<VTKFunction<GridView> const> const& fct,
+ std::optional<Vtk::DataTypes> type = {})
+ : localFct_(fct)
+ , name_(fct->name())
+ , ncomps_(fct->ncomps())
+ , type_(type ? *type : Vtk::FLOAT64)
+ {}
+
+ /// Construct VtkFunction from dune-functions GridFunction with Signature
+ // NOTE: Stores the localFunction(fct) by value.
+ /**
+ * \param fct A Grid(View)-function, providing a `localFunction(fct)`
+ * \param name The name to use component identification in the VTK file
+ * \param ncomps Number of components of the pointwise data. Is extracted
+ * from the range type of the GridFunction if not given.
+ * \param type The \ref Vtk::DataTypes used in the output. E.g. FLOAT32,
+ * or FLOAT64. Is extracted from the range type of the
+ * GridFunction if not given.
+ **/
+ template <class F,
+ class = void_t<LocalFunction<F>> >
+ Function (F&& fct, std::string name,
+ std::optional<int> ncomps = {},
+ std::optional<Vtk::DataTypes> type = {})
+ : localFct_(localFunction(std::forward<F>(fct)))
+ , name_(std::move(name))
+ {
+ using R = Range<LocalFunction<F>>;
+
+ ncomps_ = ncomps ? *ncomps : sizeOf<R>();
+ type_ = type ? *type : Vtk::Map::type<R>();
+ }
+
+ /// Constructor that forward the number of components and data type to the other constructor
+ template <class F,
+ class = void_t<LocalFunction<F>> >
+ Function (F&& fct, Vtk::FieldInfo fieldInfo,
+ std::optional<Vtk::DataTypes> type = {})
+ : Function(std::forward<F>(fct), fieldInfo.name(), fieldInfo.ncomps(), type)
+ {}
+
+ Function () = default;
+
+ /// Create a LocalFunction
+ friend Vtk::LocalFunction<GridView> localFunction (Function const& self)
+ {
+ return self.localFct_;
+ }
+
+ /// Return a name associated with the function
+ std::string const& name () const
+ {
+ return name_;
+ }
+
+ /// Return the number of components of the Range
+ int ncomps () const
+ {
+ return ncomps_ > 3 ? 9 : ncomps_ > 1 ? 3 : 1; // tensor, vector, scalar
+ }
+
+ /// Return the VTK Datatype associated with the functions range type
+ Vtk::DataTypes type () const
+ {
+ return type_;
+ }
+
+ private:
+ Vtk::LocalFunction<GridView> localFct_;
+ std::string name_;
+ int ncomps_ = 1;
+ Vtk::DataTypes type_ = Vtk::FLOAT32;
+ };
+
+ } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/gridcreatorinterface.hh b/dune/vtk/gridcreatorinterface.hh
index 02def76089d8299890dbea7c7fac6024122d9647..93c68245a118b2baf906b9a9cc5ef0092ee7a5ce 100644
--- a/dune/vtk/gridcreatorinterface.hh
+++ b/dune/vtk/gridcreatorinterface.hh
@@ -12,102 +12,110 @@
namespace Dune
{
-
- /// Base class for grid creators in a CRTP style.
- /**
- * Construct a grid from data read from VTK files.
- *
- * \tparam GridType Model of Dune::Grid
- * \tparam GlobalCoordType Type of the global coordinates.
- * \tparam DerivedType Implementation of a concrete GridCreator.
- **/
- template <class GridType, class DerivedType>
- class GridCreatorInterface
+ namespace Vtk
{
- public:
- using Grid = GridType;
- using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
- using Derived = DerivedType;
-
- public:
- /// Constructor. Stores a reference to the passed GridFactory
- GridCreatorInterface (GridFactory<Grid>& factory)
- : factory_(&factory)
- {}
-
- /// Insert all points as vertices into the factory
- void insertVertices (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& point_ids)
- {
- asDerived().insertVerticesImpl(points, point_ids);
- }
-
- /// Create elements based on type and connectivity description
- void insertElements (std::vector<std::uint8_t> const& types,
- std::vector<std::int64_t> const& offsets,
- std::vector<std::int64_t> const& connectivity)
- {
- asDerived().insertElementsImpl(types, offsets, connectivity);
- }
-
- /// Insert part of a grid stored in file into factory
- void insertPieces (std::vector<std::string> const& pieces)
- {
- asDerived().insertPiecesImpl(pieces);
- }
-
- /// Return the associated GridFactory
- GridFactory<Grid>& factory ()
+ /// Base class for grid creators in a CRTP style.
+ /**
+ * Construct a grid from data read from VTK files.
+ *
+ * \tparam GridType Model of Dune::Grid
+ * \tparam GlobalCoordType Type of the global coordinates.
+ * \tparam DerivedType Implementation of a concrete GridCreator.
+ **/
+ template <class GridType, class DerivedType>
+ class GridCreatorInterface
{
- return *factory_;
- }
-
- /// Return the associated (const) GridFactory
- GridFactory<Grid> const& factory () const
- {
- return *factory_;
- }
-
- /// Return the mpi collective communicator
- auto comm () const
- {
- return MPIHelper::getCollectiveCommunication();
- }
-
- protected: // cast to derived type
-
- Derived& asDerived ()
- {
- return static_cast<Derived&>(*this);
- }
-
- const Derived& asDerived () const
- {
- return static_cast<const Derived&>(*this);
- }
-
- public: // default implementations
-
- void insertVerticesImpl (std::vector<GlobalCoordinate> const&,
- std::vector<std::uint64_t> const&)
- {
- /* do nothing */
- }
-
- void insertElementsImpl (std::vector<std::uint8_t> const&,
- std::vector<std::int64_t> const&,
- std::vector<std::int64_t> const&)
- {
- /* do nothing */
- }
-
- void insertPiecesImpl (std::vector<std::string> const&)
- {
- /* do nothing */;
- }
-
- protected:
- GridFactory<Grid>* factory_;
- };
-
+ public:
+ using Grid = GridType;
+ using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
+ using Derived = DerivedType;
+
+ public:
+ /// Constructor. Stores a reference to the passed GridFactory
+ GridCreatorInterface (GridFactory<Grid>& factory)
+ : factory_(&factory)
+ {}
+
+ /// Insert all points as vertices into the factory
+ void insertVertices (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& point_ids)
+ {
+ asDerived().insertVerticesImpl(points, point_ids);
+ }
+
+ /// Create elements based on type and connectivity description
+ void insertElements (std::vector<std::uint8_t> const& types,
+ std::vector<std::int64_t> const& offsets,
+ std::vector<std::int64_t> const& connectivity)
+ {
+ asDerived().insertElementsImpl(types, offsets, connectivity);
+ }
+
+ /// Insert part of a grid stored in file into factory
+ void insertPieces (std::vector<std::string> const& pieces)
+ {
+ asDerived().insertPiecesImpl(pieces);
+ }
+
+ /// Construct the actual grid using the GridFactory
+ std::unique_ptr<Grid> createGrid () const
+ {
+ return factory_->createGrid();
+ }
+
+ /// Return the associated GridFactory
+ GridFactory<Grid>& factory ()
+ {
+ return *factory_;
+ }
+
+ /// Return the associated (const) GridFactory
+ GridFactory<Grid> const& factory () const
+ {
+ return *factory_;
+ }
+
+ /// Return the mpi collective communicator
+ auto comm () const
+ {
+ return MPIHelper::getCollectiveCommunication();
+ }
+
+ protected: // cast to derived type
+
+ Derived& asDerived ()
+ {
+ return static_cast<Derived&>(*this);
+ }
+
+ const Derived& asDerived () const
+ {
+ return static_cast<const Derived&>(*this);
+ }
+
+ public: // default implementations
+
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const&,
+ std::vector<std::uint64_t> const&)
+ {
+ /* do nothing */
+ }
+
+ void insertElementsImpl (std::vector<std::uint8_t> const&,
+ std::vector<std::int64_t> const&,
+ std::vector<std::int64_t> const&)
+ {
+ /* do nothing */
+ }
+
+ void insertPiecesImpl (std::vector<std::string> const&)
+ {
+ /* do nothing */;
+ }
+
+ protected:
+ GridFactory<Grid>* factory_;
+ };
+
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/CMakeLists.txt b/dune/vtk/gridcreators/CMakeLists.txt
index b605ec6eaf1b42d5fafed1a0d5fc8ddbfb9c9d54..4eed8729fcbbb4670165564878a0f0dfe7c4b18c 100644
--- a/dune/vtk/gridcreators/CMakeLists.txt
+++ b/dune/vtk/gridcreators/CMakeLists.txt
@@ -7,4 +7,4 @@ install(FILES
lagrangegridcreator.hh
parallelgridcreator.hh
serialgridcreator.hh
- DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter/gridcreators)
+ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtk/gridcreators)
diff --git a/dune/vtk/gridcreators/common.hh b/dune/vtk/gridcreators/common.hh
index e47eb4dcd50ee57d37b7bd746972ab52c694b00f..6031cec47d2c8ee006fb38791408b39007970b63 100644
--- a/dune/vtk/gridcreators/common.hh
+++ b/dune/vtk/gridcreators/common.hh
@@ -4,19 +4,23 @@
namespace Dune
{
- template <class Factory, class... Args>
- using HasInsertVertex = decltype( std::declval<Factory>().insertVertex(std::declval<Args>()...) );
-
- namespace Impl
+ namespace Vtk
{
- template <class GF, class = void>
- struct VertexIdType { using type = unsigned int; };
- template <class GF>
- struct VertexIdType<GF, typename GF::VertexId> { using type = typename GF::VertexId; };
- }
+ template <class Factory, class... Args>
+ using HasInsertVertex = decltype( std::declval<Factory>().insertVertex(std::declval<Args>()...) );
+
+ namespace Impl
+ {
+ template <class GF, class = void>
+ struct VertexIdType { using type = unsigned int; };
- template <class GF>
- using VertexId_t = typename Impl::VertexIdType<GF>::type;
+ template <class GF>
+ struct VertexIdType<GF, typename GF::VertexId> { using type = typename GF::VertexId; };
+ }
+
+ template <class GF>
+ using VertexId_t = typename Impl::VertexIdType<GF>::type;
+ } //end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/continuousgridcreator.hh b/dune/vtk/gridcreators/continuousgridcreator.hh
index 3db1d0cabcf7ce53d23f2506309253b0517e7821..6e77e60c3bef7e8513d532f57b1a1f1a02c4ae28 100644
--- a/dune/vtk/gridcreators/continuousgridcreator.hh
+++ b/dune/vtk/gridcreators/continuousgridcreator.hh
@@ -9,62 +9,65 @@
#include <dune/common/hybridutilities.hh>
#include <dune/grid/common/gridfactory.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/gridcreatorinterface.hh>
namespace Dune
{
- // Create a grid where the input points and connectivity is already
- // connected correctly.
- template <class Grid>
- struct ContinuousGridCreator
- : public GridCreatorInterface<Grid, ContinuousGridCreator<Grid>>
+ namespace Vtk
{
- using Self = ContinuousGridCreator;
- using Super = GridCreatorInterface<Grid, ContinuousGridCreator>;
- using GlobalCoordinate = typename Super::GlobalCoordinate;
-
- ContinuousGridCreator (GridFactory<Grid>& factory)
- : Super(factory)
- {}
+ // Create a grid where the input points and connectivity is already
+ // connected correctly.
+ template <class Grid>
+ struct ContinuousGridCreator
+ : public GridCreatorInterface<Grid, ContinuousGridCreator<Grid>>
+ {
+ using Self = ContinuousGridCreator;
+ using Super = GridCreatorInterface<Grid, ContinuousGridCreator>;
+ using GlobalCoordinate = typename Super::GlobalCoordinate;
- using Super::factory;
+ ContinuousGridCreator (GridFactory<Grid>& factory)
+ : Super(factory)
+ {}
- void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& /*point_ids*/)
- {
- for (auto const& p : points)
- factory().insertVertex(p);
- }
+ using Super::factory;
- void insertElementsImpl (std::vector<std::uint8_t> const& types,
- std::vector<std::int64_t> const& offsets,
- std::vector<std::int64_t> const& connectivity)
- {
- std::size_t idx = 0;
- for (std::size_t i = 0; i < types.size(); ++i) {
- auto type = Vtk::to_geometry(types[i]);
- Vtk::CellType cellType{type};
- auto refElem = referenceElement<double,Grid::dimension>(type);
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& /*point_ids*/)
+ {
+ for (auto const& p : points)
+ factory().insertVertex(p);
+ }
- int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
- assert(nNodes == refElem.size(Grid::dimension));
- std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
- for (int j = 0; j < nNodes; ++j)
- vtk_cell.push_back( connectivity[idx++] );
+ void insertElementsImpl (std::vector<std::uint8_t> const& types,
+ std::vector<std::int64_t> const& offsets,
+ std::vector<std::int64_t> const& connectivity)
+ {
+ std::size_t idx = 0;
+ for (std::size_t i = 0; i < types.size(); ++i) {
+ auto type = Vtk::to_geometry(types[i]);
+ Vtk::CellType cellType{type};
+ auto refElem = referenceElement<double,Grid::dimension>(type);
- if (cellType.noPermutation())
- factory().insertElement(type,vtk_cell);
- else {
- // apply index permutation
- std::vector<unsigned int> cell(nNodes);
+ int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
+ assert(nNodes == refElem.size(Grid::dimension));
+ std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
for (int j = 0; j < nNodes; ++j)
- cell[j] = vtk_cell[cellType.permutation(j)];
+ vtk_cell.push_back( connectivity[idx++] );
+
+ if (cellType.noPermutation())
+ factory().insertElement(type,vtk_cell);
+ else {
+ // apply index permutation
+ std::vector<unsigned int> cell(nNodes);
+ for (int j = 0; j < nNodes; ++j)
+ cell[j] = vtk_cell[cellType.permutation(j)];
- factory().insertElement(type,cell);
+ factory().insertElement(type,cell);
+ }
}
}
- }
- };
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/derivedgridcreator.hh b/dune/vtk/gridcreators/derivedgridcreator.hh
index 91fa8463f028872f9160a8cce9638945615399c3..38b2ff97a28d50ce2947008ca27475d9ca785c94 100644
--- a/dune/vtk/gridcreators/derivedgridcreator.hh
+++ b/dune/vtk/gridcreators/derivedgridcreator.hh
@@ -12,40 +12,43 @@
namespace Dune
{
- template <class GridCreator, class Derived>
- struct DerivedGridCreator
- : public GridCreatorInterface<typename GridCreator::Grid, Derived>
+ namespace Vtk
{
- using Self = DerivedGridCreator;
- using Super = GridCreatorInterface<typename GridCreator::Grid, Derived>;
- using Grid = typename GridCreator::Grid;
- using GlobalCoordinate = typename Super::GlobalCoordinate;
-
- DerivedGridCreator (GridFactory<Grid>& factory)
- : Super(factory)
- , gridCreator_(factory)
- {}
-
- void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& point_ids)
+ template <class GridCreator, class Derived>
+ struct DerivedGridCreator
+ : public GridCreatorInterface<typename GridCreator::Grid, Derived>
{
- gridCreator_.insertVertices(points, point_ids);
- }
-
- void insertElementsImpl (std::vector<std::uint8_t> const& types,
- std::vector<std::int64_t> const& offsets,
- std::vector<std::int64_t> const& connectivity)
- {
- gridCreator_.insertElements(types, offsets, connectivity);
- }
-
- void insertPiecesImpl (std::vector<std::string> const& pieces)
- {
- gridCreator_.insertPieces(pieces);
- }
-
- private:
- GridCreator gridCreator_;
- };
-
+ using Self = DerivedGridCreator;
+ using Super = GridCreatorInterface<typename GridCreator::Grid, Derived>;
+ using Grid = typename GridCreator::Grid;
+ using GlobalCoordinate = typename Super::GlobalCoordinate;
+
+ DerivedGridCreator (GridFactory<Grid>& factory)
+ : Super(factory)
+ , gridCreator_(factory)
+ {}
+
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& point_ids)
+ {
+ gridCreator_.insertVertices(points, point_ids);
+ }
+
+ void insertElementsImpl (std::vector<std::uint8_t> const& types,
+ std::vector<std::int64_t> const& offsets,
+ std::vector<std::int64_t> const& connectivity)
+ {
+ gridCreator_.insertElements(types, offsets, connectivity);
+ }
+
+ void insertPiecesImpl (std::vector<std::string> const& pieces)
+ {
+ gridCreator_.insertPieces(pieces);
+ }
+
+ private:
+ GridCreator gridCreator_;
+ };
+
+ } // end namespace Vtk;
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/discontinuousgridcreator.hh b/dune/vtk/gridcreators/discontinuousgridcreator.hh
index c285ba51a44036d26cb8adb99372a4479e7c88f5..086a9fb04865b112da71ca86925b4883532f99f5 100644
--- a/dune/vtk/gridcreators/discontinuousgridcreator.hh
+++ b/dune/vtk/gridcreators/discontinuousgridcreator.hh
@@ -9,89 +9,92 @@
#include <dune/common/hybridutilities.hh>
#include <dune/grid/common/gridfactory.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/gridcreatorinterface.hh>
namespace Dune
{
- // Create a grid where the input points are not connected and the connectivity
- // describes separated elements.
- template <class Grid>
- struct DiscontinuousGridCreator
- : public GridCreatorInterface<Grid, DiscontinuousGridCreator<Grid>>
+ namespace Vtk
{
- using Self = DiscontinuousGridCreator;
- using Super = GridCreatorInterface<Grid, DiscontinuousGridCreator>;
- using GlobalCoordinate = typename Super::GlobalCoordinate;
-
- struct CoordLess
+ // Create a grid where the input points are not connected and the connectivity
+ // describes separated elements.
+ template <class Grid>
+ struct DiscontinuousGridCreator
+ : public GridCreatorInterface<Grid, DiscontinuousGridCreator<Grid>>
{
- template <class T, int N>
- bool operator() (FieldVector<T,N> const& lhs, FieldVector<T,N> const& rhs) const
+ using Self = DiscontinuousGridCreator;
+ using Super = GridCreatorInterface<Grid, DiscontinuousGridCreator>;
+ using GlobalCoordinate = typename Super::GlobalCoordinate;
+
+ struct CoordLess
{
- for (int i = 0; i < N; ++i) {
- if (std::abs(lhs[i] - rhs[i]) < std::numeric_limits<T>::epsilon())
- continue;
- return lhs[i] < rhs[i];
+ template <class T, int N>
+ bool operator() (FieldVector<T,N> const& lhs, FieldVector<T,N> const& rhs) const
+ {
+ for (int i = 0; i < N; ++i) {
+ if (std::abs(lhs[i] - rhs[i]) < std::numeric_limits<T>::epsilon())
+ continue;
+ return lhs[i] < rhs[i];
+ }
+ return false;
}
- return false;
- }
- };
+ };
- DiscontinuousGridCreator (GridFactory<Grid>& factory)
- : Super(factory)
- {}
+ DiscontinuousGridCreator (GridFactory<Grid>& factory)
+ : Super(factory)
+ {}
- using Super::factory;
- void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& /*point_ids*/)
- {
- points_ = &points;
- uniquePoints_.clear();
- std::size_t idx = 0;
+ using Super::factory;
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& /*point_ids*/)
+ {
+ points_ = &points;
+ uniquePoints_.clear();
+ std::size_t idx = 0;
- for (auto const& p : points) {
- auto b = uniquePoints_.emplace(std::make_pair(p,idx));
- if (b.second) {
- factory().insertVertex(p);
- ++idx;
+ for (auto const& p : points) {
+ auto b = uniquePoints_.emplace(std::make_pair(p,idx));
+ if (b.second) {
+ factory().insertVertex(p);
+ ++idx;
+ }
}
}
- }
- void insertElementsImpl (std::vector<std::uint8_t> const& types,
- std::vector<std::int64_t> const& offsets,
- std::vector<std::int64_t> const& connectivity)
- {
- assert(points_ != nullptr);
- std::size_t idx = 0;
- for (std::size_t i = 0; i < types.size(); ++i) {
- auto type = Vtk::to_geometry(types[i]);
- Vtk::CellType cellType{type};
+ void insertElementsImpl (std::vector<std::uint8_t> const& types,
+ std::vector<std::int64_t> const& offsets,
+ std::vector<std::int64_t> const& connectivity)
+ {
+ assert(points_ != nullptr);
+ std::size_t idx = 0;
+ for (std::size_t i = 0; i < types.size(); ++i) {
+ auto type = Vtk::to_geometry(types[i]);
+ Vtk::CellType cellType{type};
- int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
- assert(nNodes > 0);
- std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
- for (int j = 0; j < nNodes; ++j) {
- std::size_t v_j = connectivity[idx++];
- std::size_t new_idx = uniquePoints_[(*points_)[v_j]];
- vtk_cell.push_back(new_idx);
- }
+ int nNodes = offsets[i] - (i == 0 ? 0 : offsets[i-1]);
+ assert(nNodes > 0);
+ std::vector<unsigned int> vtk_cell; vtk_cell.reserve(nNodes);
+ for (int j = 0; j < nNodes; ++j) {
+ std::size_t v_j = connectivity[idx++];
+ std::size_t new_idx = uniquePoints_[(*points_)[v_j]];
+ vtk_cell.push_back(new_idx);
+ }
- if (cellType.noPermutation()) {
- factory().insertElement(type,vtk_cell);
- } else {
- // apply index permutation
- std::vector<unsigned int> cell(nNodes);
- for (int j = 0; j < nNodes; ++j)
- cell[j] = vtk_cell[cellType.permutation(j)];
- factory().insertElement(type,cell);
+ if (cellType.noPermutation()) {
+ factory().insertElement(type,vtk_cell);
+ } else {
+ // apply index permutation
+ std::vector<unsigned int> cell(nNodes);
+ for (int j = 0; j < nNodes; ++j)
+ cell[j] = vtk_cell[cellType.permutation(j)];
+ factory().insertElement(type,cell);
+ }
}
}
- }
- private:
- std::vector<GlobalCoordinate> const* points_ = nullptr;
- std::map<GlobalCoordinate, std::size_t, CoordLess> uniquePoints_;
- };
+ private:
+ std::vector<GlobalCoordinate> const* points_ = nullptr;
+ std::map<GlobalCoordinate, std::size_t, CoordLess> uniquePoints_;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/lagrangegridcreator.hh b/dune/vtk/gridcreators/lagrangegridcreator.hh
index 8901e897b3bfe026b1af568b66331b7c8f7a8e56..4f183784d6a8db695a1664a4ade9fa6f73778b96 100644
--- a/dune/vtk/gridcreators/lagrangegridcreator.hh
+++ b/dune/vtk/gridcreators/lagrangegridcreator.hh
@@ -9,362 +9,381 @@
#include <dune/common/exceptions.hh>
#include <dune/common/hybridutilities.hh>
#include <dune/geometry/utility/typefromvertexcount.hh>
+#include <dune/geometry/multilineargeometry.hh>
#include <dune/localfunctions/lagrange.hh>
#include <dune/grid/common/gridfactory.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/gridcreatorinterface.hh>
#include <dune/vtk/utility/lagrangepoints.hh>
namespace Dune
{
- // \brief Create a grid from data that represents higher (lagrange) cells.
- /**
- * The grid is created from the first nodes of a cell parametrization, representing
- * the corner vertices. Thus a piecewise "flat" grid is constructed. The
- * parametrization is 1. passed as a local element parametrization to the
- * `insertElement()` function of a gridFactory to allow the grid itself to handle the
- * parametrization and 2. is stored internally that can be accessed by using this
- * GridCreator object as a grid function, or by extracting locally the parametrization
- * on each existing grid element after creation of the grid.
- *
- * So, the LagrangeGridCreator models both, a `GridCreator` and a `GridFunction`.
- **/
- template <class GridType>
- struct LagrangeGridCreator
- : public GridCreatorInterface<GridType, LagrangeGridCreator<GridType>>
+ namespace Vtk
{
- using Self = LagrangeGridCreator;
- using Super = GridCreatorInterface<GridType, Self>;
- using GlobalCoordinate = typename Super::GlobalCoordinate;
-
- using Nodes = std::vector<GlobalCoordinate>;
-
- struct ElementParametrization
+ // \brief Create a grid from data that represents higher (lagrange) cells.
+ /**
+ * The grid is created from the first nodes of a cell parametrization, representing
+ * the corner vertices. Thus a piecewise "flat" grid is constructed. The
+ * parametrization is 1. passed as a local element parametrization to the
+ * `insertElement()` function of a gridFactory to allow the grid itself to handle the
+ * parametrization and 2. is stored internally that can be accessed by using this
+ * GridCreator object as a grid function, or by extracting locally the parametrization
+ * on each existing grid element after creation of the grid.
+ *
+ * So, the LagrangeGridCreator models both, a `GridCreator` and a `GridFunction`.
+ **/
+ template <class GridType>
+ struct LagrangeGridCreator
+ : public GridCreatorInterface<GridType, LagrangeGridCreator<GridType>>
{
- GeometryType type; //< Geometry type of the element
- std::vector<std::int64_t> nodes; //< Indices of the w.r.t. `nodes_` vector
- std::vector<unsigned int> corners; //< Insertion-indices of the element corner nodes
- };
-
- using Parametrization = std::vector<ElementParametrization>;
- using Element = typename GridType::template Codim<0>::Entity;
- using LocalCoordinate = typename Element::Geometry::LocalCoordinate;
-
- class LocalParametrization;
- class LocalFunction;
-
- public:
- using Super::factory;
-
- LagrangeGridCreator (GridFactory<GridType>& factory)
- : Super(factory)
- {}
+ using Self = LagrangeGridCreator;
+ using Super = GridCreatorInterface<GridType, Self>;
+ using GlobalCoordinate = typename Super::GlobalCoordinate;
+
+ using Nodes = std::vector<GlobalCoordinate>;
+
+ struct ElementParametrization
+ {
+ GeometryType type; //< Geometry type of the element
+ std::vector<std::int64_t> nodes; //< Indices of the w.r.t. `nodes_` vector
+ std::vector<unsigned int> corners; //< Insertion-indices of the element corner nodes
+ };
+
+ using Parametrization = std::vector<ElementParametrization>;
+ using Element = typename GridType::template Codim<0>::Entity;
+ using LocalCoordinate = typename Element::Geometry::LocalCoordinate;
+
+ class LocalParametrization;
+ class LocalFunction;
+
+ public:
+ using Super::factory;
+
+ LagrangeGridCreator (GridFactory<GridType>& factory)
+ : Super(factory)
+ {}
+
+ /// Implementation of the interface function `insertVertices()`
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& /*point_ids*/)
+ {
+ // store point coordinates in member variable
+ nodes_ = points;
+ }
- /// Implementation of the interface function `insertVertices()`
- void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& /*point_ids*/)
- {
- // store point coordinates in member variable
- nodes_ = points;
- }
-
- template <class F>
- using HasParametrizedElements = decltype(std::declval<F>().insertElement(std::declval<GeometryType>(),
- std::declval<std::vector<unsigned int> const&>(), std::declval<std::function<GlobalCoordinate(LocalCoordinate)>>()));
-
- /// Implementation of the interface function `insertElements()`
- void insertElementsImpl (std::vector<std::uint8_t> const& types,
- std::vector<std::int64_t> const& offsets,
- std::vector<std::int64_t> const& connectivity)
- {
- assert(nodes_.size() > 0);
-
- // mapping of node index to element-vertex index
- std::vector<std::int64_t> elementVertices(nodes_.size(), -1);
- parametrization_.reserve(types.size());
-
- std::int64_t vertexIndex = 0;
- for (std::size_t i = 0; i < types.size(); ++i) {
- auto type = Vtk::to_geometry(types[i]);
- if (type.dim() != GridType::dimension)
- continue;
-
- Vtk::CellType cellType{type};
- auto refElem = referenceElement<double,GridType::dimension>(type);
-
- std::int64_t shift = (i == 0 ? 0 : offsets[i-1]);
- int nNodes = offsets[i] - shift;
- int nVertices = refElem.size(GridType::dimension);
-
- // insert vertices into grid and construct element vertices
- std::vector<unsigned int> element(nVertices);
- for (int j = 0; j < nVertices; ++j) {
- auto index = connectivity.at(shift + j);
- auto& vertex = elementVertices.at(index);
- if (vertex < 0) {
- factory().insertVertex(nodes_.at(index));
- vertex = vertexIndex++;
+ template <class F>
+ using HasParametrizedElements = decltype(std::declval<F>().insertElement(std::declval<GeometryType>(),
+ std::declval<std::vector<unsigned int> const&>(), std::declval<std::function<GlobalCoordinate(LocalCoordinate)>>()));
+
+ /// Implementation of the interface function `insertElements()`
+ void insertElementsImpl (std::vector<std::uint8_t> const& types,
+ std::vector<std::int64_t> const& offsets,
+ std::vector<std::int64_t> const& connectivity)
+ {
+ assert(nodes_.size() > 0);
+
+ // mapping of node index to element-vertex index
+ std::vector<std::int64_t> elementVertices(nodes_.size(), -1);
+ parametrization_.reserve(types.size());
+
+ std::int64_t vertexIndex = 0;
+ for (std::size_t i = 0; i < types.size(); ++i) {
+ auto type = Vtk::to_geometry(types[i]);
+ if (type.dim() != GridType::dimension)
+ continue;
+
+ Vtk::CellType cellType{type};
+ auto refElem = referenceElement<double,GridType::dimension>(type);
+
+ std::int64_t shift = (i == 0 ? 0 : offsets[i-1]);
+ int nNodes = offsets[i] - shift;
+ int nVertices = refElem.size(GridType::dimension);
+
+ // insert vertices into grid and construct element vertices
+ std::vector<unsigned int> element(nVertices);
+ for (int j = 0; j < nVertices; ++j) {
+ auto index = connectivity.at(shift + j);
+ auto& vertex = elementVertices.at(index);
+ if (vertex < 0) {
+ factory().insertVertex(nodes_.at(index));
+ vertex = vertexIndex++;
+ }
+ element[j] = vertex;
}
- element[j] = vertex;
- }
-
- // permute element indices
- if (!cellType.noPermutation()) {
- // apply index permutation
- std::vector<unsigned int> cell(element.size());
- for (int j = 0; j < element.size(); ++j)
- cell[j] = element[cellType.permutation(j)];
- std::swap(element, cell);
- }
- // fill vector of element parametrizations
- parametrization_.push_back(ElementParametrization{type});
- auto& param = parametrization_.back();
-
- param.nodes.resize(nNodes);
- for (int j = 0; j < nNodes; ++j)
- param.nodes[j] = connectivity.at(shift + j);
- param.corners = element;
+ // permute element indices
+ if (!cellType.noPermutation()) {
+ // apply index permutation
+ std::vector<unsigned int> cell(element.size());
+ for (int j = 0; j < element.size(); ++j)
+ cell[j] = element[cellType.permutation(j)];
+ std::swap(element, cell);
+ }
- // try to create element with parametrization
- if constexpr (Std::is_detected_v<HasParametrizedElements, GridFactory<GridType>>) {
- try {
- factory().insertElement(type, element, localParametrization(parametrization_.size()-1));
- } catch (Dune::GridError const& /* notImplemented */) {
+ // fill vector of element parametrizations
+ parametrization_.push_back(ElementParametrization{type});
+ auto& param = parametrization_.back();
+
+ param.nodes.resize(nNodes);
+ for (int j = 0; j < nNodes; ++j)
+ param.nodes[j] = connectivity.at(shift + j);
+ param.corners = element;
+
+ // try to create element with parametrization
+ if constexpr (Std::is_detected_v<HasParametrizedElements, GridFactory<GridType>>) {
+ try {
+ factory().insertElement(type, element,
+ localParametrization(parametrization_.size()-1));
+ } catch (Dune::GridError const& /* notImplemented */) {
+ factory().insertElement(type, element);
+ }
+ } else {
factory().insertElement(type, element);
}
- } else {
- factory().insertElement(type, element);
}
}
- }
- /// \brief Construct an element parametrization
- /**
- * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
- * where `LocalCoordinate is w.r.t. the local coordinate system in an element with
- * given `insertionIndex` (defined by the inserted corner vertices) and
- * `GlobalCoordinate` a world coordinate in the parametrized grid.
- **/
- LocalParametrization localParametrization (unsigned int insertionIndex) const
- {
- assert(!nodes_.empty() && !parametrization_.empty());
- auto const& localParam = parametrization_.at(insertionIndex);
- return LocalParametrization{nodes_, localParam, order(localParam)};
- }
-
- /// \brief Construct an element parametrization
- /**
- * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
- * where `LocalCoordinate is w.r.t. the local coordinate system in the passed element
- * and `GlobalCoordinate` a world coordinate in the parametrized grid.
- *
- * Note, when an element is passed, it might have a different local coordinate system
- * than the coordinate system used to defined the element parametrization. Thus
- * coordinate transform is internally chained to the evaluation of the local
- * parametrization. This local geometry transform is obtained by figuring out the
- * permutation of corners in the element corresponding to the inserted corner
- * vertices.
- **/
- LocalParametrization localParametrization (Element const& element) const
- {
- assert(!nodes_.empty() && !parametrization_.empty());
-
- unsigned int insertionIndex = factory().insertionIndex(element);
- auto const& localParam = parametrization_.at(insertionIndex);
- assert(element.type() == localParam.type);
-
- // collect indices of vertices
- std::vector<unsigned int> indices(element.subEntities(GridType::dimension));
- for (unsigned int i = 0; i < element.subEntities(GridType::dimension); ++i)
- indices[i] = factory().insertionIndex(element.template subEntity<GridType::dimension>(i));
-
- // calculate permutation vector
- std::vector<unsigned int> permutation(indices.size());
- for (std::size_t i = 0; i < indices.size(); ++i) {
- auto it = std::find(localParam.corners.begin(), localParam.corners.end(), indices[i]);
- assert(it != localParam.corners.end());
- permutation[i] = std::distance(localParam.corners.begin(), it);
+ /// \brief Construct an element parametrization
+ /**
+ * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
+ * where `LocalCoordinate is w.r.t. the local coordinate system in an element with
+ * given `insertionIndex` (defined by the inserted corner vertices) and
+ * `GlobalCoordinate` a world coordinate in the parametrized grid.
+ **/
+ LocalParametrization localParametrization (unsigned int insertionIndex) const
+ {
+ assert(!nodes_.empty() && !parametrization_.empty());
+ auto const& localParam = parametrization_.at(insertionIndex);
+ return LocalParametrization{nodes_, localParam, order(localParam)};
}
- return LocalParametrization{nodes_, localParam, order(localParam), permutation};
- }
+ /// \brief Construct an element parametrization
+ /**
+ * The returned LocalParametrization is a mapping `GlobalCoordinate(LocalCoordinate)`
+ * where `LocalCoordinate is w.r.t. the local coordinate system in the passed element
+ * and `GlobalCoordinate` a world coordinate in the parametrized grid.
+ *
+ * Note, when an element is passed, it might have a different local coordinate system
+ * than the coordinate system used to defined the element parametrization. Thus
+ * coordinate transform is internally chained to the evaluation of the local
+ * parametrization. This local geometry transform is obtained by figuring out the
+ * permutation of corners in the element corresponding to the inserted corner
+ * vertices.
+ **/
+ LocalParametrization localParametrization (Element const& element) const
+ {
+ assert(!nodes_.empty() && !parametrization_.empty());
+
+ unsigned int insertionIndex = factory().insertionIndex(element);
+ auto const& localParam = parametrization_.at(insertionIndex);
+ assert(element.type() == localParam.type);
+
+ // collect indices of vertices
+ std::vector<unsigned int> indices(element.subEntities(GridType::dimension));
+ for (unsigned int i = 0; i < element.subEntities(GridType::dimension); ++i)
+ indices[i] = factory().insertionIndex(element.template subEntity<GridType::dimension>(i));
+
+ // calculate permutation vector
+ std::vector<unsigned int> permutation(indices.size());
+ for (std::size_t i = 0; i < indices.size(); ++i) {
+ auto it = std::find(localParam.corners.begin(), localParam.corners.end(), indices[i]);
+ assert(it != localParam.corners.end());
+ permutation[i] = std::distance(localParam.corners.begin(), it);
+ }
- /// Determine lagrange order from number of points
- template <class LocalParam>
- unsigned int order (LocalParam const& localParam) const
- {
- GeometryType type = localParam.type;
- std::size_t nNodes = localParam.nodes.size();
- for (unsigned int o = 1; o <= nNodes; ++o)
- if (numLagrangePoints(type.id(), type.dim(), o) == nNodes)
- return o;
+ return LocalParametrization{nodes_, localParam, order(localParam), permutation};
+ }
- return 1;
- }
+ /// Determine lagrange order from number of points
+ template <class LocalParam>
+ unsigned int order (LocalParam const& localParam) const
+ {
+ GeometryType type = localParam.type;
+ std::size_t nNodes = localParam.nodes.size();
+ for (unsigned int o = 1; o <= nNodes; ++o)
+ if (numLagrangePoints(type.id(), type.dim(), o) == nNodes)
+ return o;
+
+ return 1;
+ }
- /// Determine lagrange order from number of points from the first element parametrization
- unsigned int order () const
- {
- assert(!parametrization_.empty());
- return order(parametrization_.front());
- }
+ /// Determine lagrange order from number of points from the first element parametrization
+ unsigned int order () const
+ {
+ assert(!parametrization_.empty());
+ return order(parametrization_.front());
+ }
- public:
- /// \brief Local function representing the parametrization of the grid.
- /**
- * The returned object models Functions::Concept::LocalFunction
- * and can thus be bound to an element of the created grid and evaluated in
- * the local coordinates of the bound element.
- *
- * It is implemented in terms of the \ref LocalParametrization function
- * returned by the method \ref localParametrization(element). See comments
- * there for further details.
- *
- * Note, this methods requires the GridCreator to be based by
- * lvalue-reference. This is necessary, since we want to guarantee that all
- * internal storage is preserved while evaluating the local function.
- **/
- friend LocalFunction localFunction (LagrangeGridCreator& gridCreator)
- {
- return LocalFunction{gridCreator};
- }
+ public:
+ /// \brief Local function representing the parametrization of the grid.
+ /**
+ * The returned object models Functions::Concept::LocalFunction
+ * and can thus be bound to an element of the created grid and evaluated in
+ * the local coordinates of the bound element.
+ *
+ * It is implemented in terms of the \ref LocalParametrization function
+ * returned by the method \ref localParametrization(element). See comments
+ * there for further details.
+ *
+ * Note, this methods requires the GridCreator to be based by
+ * lvalue-reference. This is necessary, since we want to guarantee that all
+ * internal storage is preserved while evaluating the local function.
+ **/
+ friend LocalFunction localFunction (LagrangeGridCreator& gridCreator)
+ {
+ return LocalFunction{gridCreator};
+ }
- struct EntitySet
- {
- using Grid = GridType;
- };
+ friend LocalFunction localFunction (LagrangeGridCreator const& gridCreator)
+ {
+ return LocalFunction{gridCreator};
+ }
- /// Dummy function returning a placeholder entityset
- EntitySet entitySet () const
- {
- assert(false && "Should not be used!");
- return EntitySet{};
- }
+ friend LocalFunction localFunction (LagrangeGridCreator&& gridCreator)
+ {
+ DUNE_THROW(Dune::Exception, "Cannot pass temporary LagrangeGridCreator to localFunction(). Pass an lvalue-reference instead.");
+ return LocalFunction{gridCreator};
+ }
- /// Dummy function returning a placeholder entityset
- GlobalCoordinate operator() (GlobalCoordinate const&) const
- {
- assert(false && "Should not be used!");
- return GlobalCoordinate{};
- }
+ struct EntitySet
+ {
+ using Grid = GridType;
+ using GlobalCoordinate = typename Self::GlobalCoordinate;
+ };
+
+ /// Dummy function returning a placeholder entityset
+ EntitySet entitySet () const
+ {
+ assert(false && "Should not be used!");
+ return EntitySet{};
+ }
- private:
- /// All point coordinates inclusing the higher-order lagrange points
- Nodes nodes_;
+ /// Dummy function returning a placeholder entityset
+ GlobalCoordinate operator() (GlobalCoordinate const&) const
+ {
+ assert(false && "Should not be used!");
+ return GlobalCoordinate{};
+ }
- /// Parametrization for all elements
- Parametrization parametrization_;
- };
+ private:
+ /// All point coordinates inclusing the higher-order lagrange points
+ Nodes nodes_;
+ /// Parametrization for all elements
+ Parametrization parametrization_;
+ };
- template <class Grid>
- class LagrangeGridCreator<Grid>::LocalParametrization
- {
- using ctype = typename Grid::ctype;
-
- using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
- using LocalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::LocalCoordinate;
- using LocalGeometry = MultiLinearGeometry<ctype,Grid::dimension,Grid::dimension>;
-
- using LocalFE = LagrangeLocalFiniteElement<VtkLagrangePointSet, Grid::dimension, ctype, ctype>;
- using LocalBasis = typename LocalFE::Traits::LocalBasisType;
- using LocalBasisTraits = typename LocalBasis::Traits;
-
- public:
- /// Construct a local element parametrization
- template <class Nodes, class LocalParam>
- LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order)
- : localFE_(param.type, order)
- , localNodes_(param.nodes.size())
- {
- for (std::size_t i = 0; i < localNodes_.size(); ++i)
- localNodes_[i] = nodes[param.nodes[i]];
- }
-
- /// Construct a local element parametrization for elements with permuted corners
- template <class Nodes, class LocalParam, class Permutation>
- LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order, Permutation const& permutation)
- : LocalParametrization(nodes, param, order)
+
+ template <class Grid>
+ class LagrangeGridCreator<Grid>::LocalParametrization
{
- auto refElem = referenceElement<ctype,Grid::dimension>(param.type);
- std::vector<LocalCoordinate> corners(permutation.size());
- for (std::size_t i = 0; i < permutation.size(); ++i)
- corners[i] = refElem.position(permutation[i], Grid::dimension);
+ using ctype = typename Grid::ctype;
+
+ using GlobalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::GlobalCoordinate;
+ using LocalCoordinate = typename Grid::template Codim<0>::Entity::Geometry::LocalCoordinate;
+ using LocalGeometry = MultiLinearGeometry<ctype,Grid::dimension,Grid::dimension>;
+
+ using LocalFE = LagrangeLocalFiniteElement<Vtk::LagrangePointSet, Grid::dimension, ctype, ctype>;
+ using LocalBasis = typename LocalFE::Traits::LocalBasisType;
+ using LocalBasisTraits = typename LocalBasis::Traits;
+
+ public:
+ /// Construct a local element parametrization
+ template <class Nodes, class LocalParam>
+ LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order)
+ : localFE_(param.type, order)
+ , localNodes_(param.nodes.size())
+ {
+ for (std::size_t i = 0; i < localNodes_.size(); ++i)
+ localNodes_[i] = nodes[param.nodes[i]];
+ }
- localGeometry_.emplace(param.type, corners);
- }
+ /// Construct a local element parametrization for elements with permuted corners
+ template <class Nodes, class LocalParam, class Permutation>
+ LocalParametrization (Nodes const& nodes, LocalParam const& param, unsigned int order, Permutation const& permutation)
+ : LocalParametrization(nodes, param, order)
+ {
+ auto refElem = referenceElement<ctype,Grid::dimension>(param.type);
+ std::vector<LocalCoordinate> corners(permutation.size());
+ for (std::size_t i = 0; i < permutation.size(); ++i)
+ corners[i] = refElem.position(permutation[i], Grid::dimension);
+
+ localGeometry_.emplace(param.type, corners);
+ }
- /// Evaluate the local parametrization in local coordinates
- template <class LocalCoordinate>
- GlobalCoordinate operator() (LocalCoordinate const& local) const
- {
- // map coordinates if element corners are permuted
- LocalCoordinate x = localGeometry_ ? localGeometry_->global(local) : local;
+ /// Evaluate the local parametrization in local coordinates
+ template <class LocalCoordinate>
+ GlobalCoordinate operator() (LocalCoordinate const& local) const
+ {
+ // map coordinates if element corners are permuted
+ LocalCoordinate x = localGeometry_ ? localGeometry_->global(local) : local;
- LocalBasis const& localBasis = localFE_.localBasis();
- localBasis.evaluateFunction(x, shapeValues_);
- assert(shapeValues_.size() == localNodes_.size());
+ LocalBasis const& localBasis = localFE_.localBasis();
+ localBasis.evaluateFunction(x, shapeValues_);
+ assert(shapeValues_.size() == localNodes_.size());
- GlobalCoordinate out(0);
- for (std::size_t i = 0; i < shapeValues_.size(); ++i)
- out.axpy(shapeValues_[i], localNodes_[i]);
+ GlobalCoordinate out(0);
+ for (std::size_t i = 0; i < shapeValues_.size(); ++i)
+ out.axpy(shapeValues_[i], localNodes_[i]);
- return out;
- }
+ return out;
+ }
- private:
- LocalFE localFE_;
- std::vector<GlobalCoordinate> localNodes_;
- std::optional<LocalGeometry> localGeometry_;
+ private:
+ LocalFE localFE_;
+ std::vector<GlobalCoordinate> localNodes_;
+ std::optional<LocalGeometry> localGeometry_;
- mutable std::vector<typename LocalBasisTraits::RangeType> shapeValues_;
- };
+ mutable std::vector<typename LocalBasisTraits::RangeType> shapeValues_;
+ };
- template <class Grid>
- class LagrangeGridCreator<Grid>::LocalFunction
- {
- using ctype = typename Grid::ctype;
- using LocalContext = typename Grid::template Codim<0>::Entity;
- using GlobalCoordinate = typename LocalContext::Geometry::GlobalCoordinate;
- using LocalCoordinate = typename LocalContext::Geometry::LocalCoordinate;
- using LocalParametrization = typename LagrangeGridCreator::LocalParametrization;
-
- public:
- explicit LocalFunction (LagrangeGridCreator& gridCreator)
- : gridCreator_(&gridCreator)
- {}
-
- /// Collect a local parametrization on the element
- void bind (LocalContext const& element)
+ template <class Grid>
+ class LagrangeGridCreator<Grid>::LocalFunction
{
- localContext_ = element;
- localParametrization_.emplace(gridCreator_->localParametrization(element));
- }
+ using ctype = typename Grid::ctype;
+ using LocalContext = typename Grid::template Codim<0>::Entity;
+ using GlobalCoordinate = typename LocalContext::Geometry::GlobalCoordinate;
+ using LocalCoordinate = typename LocalContext::Geometry::LocalCoordinate;
+ using LocalParametrization = typename LagrangeGridCreator::LocalParametrization;
+
+ public:
+ explicit LocalFunction (LagrangeGridCreator const& gridCreator)
+ : gridCreator_(&gridCreator)
+ {}
+
+ explicit LocalFunction (LagrangeGridCreator&& gridCreator) = delete;
+
+ /// Collect a local parametrization on the element
+ void bind (LocalContext const& element)
+ {
+ localContext_ = element;
+ localParametrization_.emplace(gridCreator_->localParametrization(element));
+ }
- void unbind () { /* do nothing */ }
+ void unbind () { /* do nothing */ }
- /// Evaluate the local parametrization in local coordinates
- GlobalCoordinate operator() (LocalCoordinate const& local) const
- {
- assert(!!localParametrization_);
- return (*localParametrization_)(local);
- }
+ /// Evaluate the local parametrization in local coordinates
+ GlobalCoordinate operator() (LocalCoordinate const& local) const
+ {
+ assert(!!localParametrization_);
+ return (*localParametrization_)(local);
+ }
- /// Return the bound element
- LocalContext const& localContext () const
- {
- return localContext_;
- }
+ /// Return the bound element
+ LocalContext const& localContext () const
+ {
+ return localContext_;
+ }
- private:
- LagrangeGridCreator const* gridCreator_;
+ private:
+ LagrangeGridCreator const* gridCreator_;
- LocalContext localContext_;
- std::optional<LocalParametrization> localParametrization_;
- };
+ LocalContext localContext_;
+ std::optional<LocalParametrization> localParametrization_;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/parallelgridcreator.hh b/dune/vtk/gridcreators/parallelgridcreator.hh
index 73cb46058d79f26b2b503214a5073f89604bdf9e..89c4584b2815b6ab8ee77ef3a95b57f604857576 100644
--- a/dune/vtk/gridcreators/parallelgridcreator.hh
+++ b/dune/vtk/gridcreators/parallelgridcreator.hh
@@ -13,38 +13,41 @@
namespace Dune
{
- // create a distributed grid in parallel. Currently only supported by ALUGrid
- template <class Grid>
- struct ParallelGridCreator
- : public DerivedGridCreator<ContinuousGridCreator<Grid>, ParallelGridCreator<Grid>>
+ namespace Vtk
{
- using Self = ParallelGridCreator;
- using Super = DerivedGridCreator<ContinuousGridCreator<Grid>, Self>;
- using GlobalCoordinate = typename Super::GlobalCoordinate;
- using VertexId = VertexId_t<GridFactory<Grid>>;
-
- // The GridFactory must support insertion of global vertex IDs
- static_assert(Std::is_detected<HasInsertVertex, GridFactory<Grid>, GlobalCoordinate, VertexId>{}, "");
-
- ParallelGridCreator (GridFactory<Grid>& factory)
- : Super(factory)
- {}
-
- void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& point_ids)
+ // create a distributed grid in parallel. Currently only supported by ALUGrid
+ template <class Grid>
+ struct ParallelGridCreator
+ : public DerivedGridCreator<ContinuousGridCreator<Grid>, ParallelGridCreator<Grid>>
{
- assert(point_ids.size() == points.size());
- for (std::size_t i = 0; i < points.size(); ++i)
- this->factory().insertVertex(points[i], VertexId(point_ids[i]));
- }
+ using Self = ParallelGridCreator;
+ using Super = DerivedGridCreator<ContinuousGridCreator<Grid>, Self>;
+ using GlobalCoordinate = typename Super::GlobalCoordinate;
+ using VertexId = VertexId_t<GridFactory<Grid>>;
+
+ // The GridFactory must support insertion of global vertex IDs
+ static_assert(Std::is_detected<HasInsertVertex, GridFactory<Grid>, GlobalCoordinate, VertexId>{}, "");
+
+ ParallelGridCreator (GridFactory<Grid>& factory)
+ : Super(factory)
+ {}
+
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& point_ids)
+ {
+ assert(point_ids.size() == points.size());
+ for (std::size_t i = 0; i < points.size(); ++i)
+ this->factory().insertVertex(points[i], VertexId(point_ids[i]));
+ }
- void insertPiecesImpl (std::vector<std::string> const& pieces)
- {
- if (int(pieces.size()) == this->comm().size()) {
- VtkReader<Grid, Self> pieceReader(this->factory());
- pieceReader.readFromFile(pieces[this->comm().rank()], true);
+ void insertPiecesImpl (std::vector<std::string> const& pieces)
+ {
+ if (int(pieces.size()) == this->comm().size()) {
+ VtkReader<Grid, Self> pieceReader(this->factory());
+ pieceReader.read(pieces[this->comm().rank()], true);
+ }
}
- }
- };
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/gridcreators/serialgridcreator.hh b/dune/vtk/gridcreators/serialgridcreator.hh
index d729149529ef8c64a7c6a30adf2dc0973ed243d2..f438fdad0215afebd6ebf31e7a03cfc494fd8d23 100644
--- a/dune/vtk/gridcreators/serialgridcreator.hh
+++ b/dune/vtk/gridcreators/serialgridcreator.hh
@@ -11,64 +11,67 @@
namespace Dune
{
- // create a distributed grid on rank 0. Needs to be load balanced afterwards.
- template <class Grid>
- struct SerialGridCreator
- : public GridCreatorInterface<Grid, SerialGridCreator<Grid>>
+ namespace Vtk
{
- using Self = SerialGridCreator;
- using Super = GridCreatorInterface<Grid, Self>;
- using GlobalCoordinate = typename Super::GlobalCoordinate;
+ // create a distributed grid on rank 0. Needs to be load balanced afterwards.
+ template <class Grid>
+ struct SerialGridCreator
+ : public GridCreatorInterface<Grid, SerialGridCreator<Grid>>
+ {
+ using Self = SerialGridCreator;
+ using Super = GridCreatorInterface<Grid, Self>;
+ using GlobalCoordinate = typename Super::GlobalCoordinate;
- SerialGridCreator (GridFactory<Grid>& factory)
- : Super(factory)
- {}
+ SerialGridCreator (GridFactory<Grid>& factory)
+ : Super(factory)
+ {}
- void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
- std::vector<std::uint64_t> const& /*point_ids*/)
- {
- shift_.push_back(points_.size());
- points_.reserve(points_.size() + points.size());
- points_.insert(points_.end(), points.begin(), points.end());
- }
+ void insertVerticesImpl (std::vector<GlobalCoordinate> const& points,
+ std::vector<std::uint64_t> const& /*point_ids*/)
+ {
+ shift_.push_back(points_.size());
+ points_.reserve(points_.size() + points.size());
+ points_.insert(points_.end(), points.begin(), points.end());
+ }
- void insertElementsImpl (std::vector<std::uint8_t> const& types,
- std::vector<std::int64_t> const& offsets,
- std::vector<std::int64_t> const& connectivity)
- {
- types_.reserve(types_.size() + types.size());
- types_.insert(types_.end(), types.begin(), types.end());
+ void insertElementsImpl (std::vector<std::uint8_t> const& types,
+ std::vector<std::int64_t> const& offsets,
+ std::vector<std::int64_t> const& connectivity)
+ {
+ types_.reserve(types_.size() + types.size());
+ types_.insert(types_.end(), types.begin(), types.end());
- offsets_.reserve(offsets_.size() + offsets.size());
- std::transform(offsets.begin(), offsets.end(), std::back_inserter(offsets_),
- [shift=offsets_.empty() ? 0 : offsets_.back()](std::int64_t o) { return o + shift; });
+ offsets_.reserve(offsets_.size() + offsets.size());
+ std::transform(offsets.begin(), offsets.end(), std::back_inserter(offsets_),
+ [shift=offsets_.empty() ? 0 : offsets_.back()](std::int64_t o) { return o + shift; });
- connectivity_.reserve(connectivity_.size() + connectivity.size());
- std::transform(connectivity.begin(), connectivity.end(), std::back_inserter(connectivity_),
- [shift=shift_.back()](std::int64_t idx) { return idx + shift; });
- }
+ connectivity_.reserve(connectivity_.size() + connectivity.size());
+ std::transform(connectivity.begin(), connectivity.end(), std::back_inserter(connectivity_),
+ [shift=shift_.back()](std::int64_t idx) { return idx + shift; });
+ }
- void insertPiecesImpl (std::vector<std::string> const& pieces)
- {
- if (this->comm().rank() == 0) {
- VtkReader<Grid, Self> pieceReader(*this);
- for (std::string const& piece : pieces) {
- pieceReader.readFromFile(piece, false);
- pieceReader.createGrid(false);
- }
+ void insertPiecesImpl (std::vector<std::string> const& pieces)
+ {
+ if (this->comm().rank() == 0) {
+ VtkReader<Grid, Self> pieceReader(*this);
+ for (std::string const& piece : pieces) {
+ pieceReader.read(piece, false);
+ pieceReader.fillGridCreator(false);
+ }
- DiscontinuousGridCreator<Grid> creator(this->factory());
- creator.insertVertices(points_, {});
- creator.insertElements(types_, offsets_, connectivity_);
+ DiscontinuousGridCreator<Grid> creator(this->factory());
+ creator.insertVertices(points_, {});
+ creator.insertElements(types_, offsets_, connectivity_);
+ }
}
- }
- private:
- std::vector<GlobalCoordinate> points_;
- std::vector<std::uint8_t> types_;
- std::vector<std::int64_t> offsets_;
- std::vector<std::int64_t> connectivity_;
- std::vector<std::int64_t> shift_;
- };
+ private:
+ std::vector<GlobalCoordinate> points_;
+ std::vector<std::uint8_t> types_;
+ std::vector<std::int64_t> offsets_;
+ std::vector<std::int64_t> connectivity_;
+ std::vector<std::int64_t> shift_;
+ };
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/legacyvtkfunction.hh b/dune/vtk/legacyvtkfunction.hh
index c0339e8496b11f1cce61be9dc3666ca45a6e650f..47faf6fe9099911037f2ce235ee56545a868e289 100644
--- a/dune/vtk/legacyvtkfunction.hh
+++ b/dune/vtk/legacyvtkfunction.hh
@@ -4,46 +4,49 @@
#include <dune/grid/io/file/vtk/function.hh>
-#include "vtklocalfunctioninterface.hh"
+#include "localfunctioninterface.hh"
namespace Dune
{
- /// Type erasure for Legacy VTKFunction
- template <class GridView>
- class VTKLocalFunctionWrapper final
- : public VtkLocalFunctionInterface<GridView>
+ namespace Vtk
{
- using Interface = VtkLocalFunctionInterface<GridView>;
- using Entity = typename Interface::Entity;
- using LocalCoordinate = typename Interface::LocalCoordinate;
-
- public:
- /// Constructor. Stores a shared pointer to the passed Dune::VTKFunction
- VTKLocalFunctionWrapper (std::shared_ptr<VTKFunction<GridView> const> const& fct)
- : fct_(fct)
- {}
-
- /// Stores a pointer to the passed entity
- virtual void bind (Entity const& entity) override
+ /// Type erasure for Legacy VTKFunction
+ template <class GridView>
+ class VTKLocalFunctionWrapper final
+ : public LocalFunctionInterface<GridView>
{
- entity_ = &entity;
- }
-
- /// Unsets the stored entity pointer
- virtual void unbind () override
- {
- entity_ = nullptr;
- }
-
- /// Evaluate the Dune::VTKFunction in LocalCoordinates on the stored Entity
- virtual double evaluate (int comp, LocalCoordinate const& xi) const override
- {
- return fct_->evaluate(comp, *entity_, xi);
- }
-
- private:
- std::shared_ptr<VTKFunction<GridView> const> fct_;
- Entity const* entity_;
- };
-
+ using Interface = LocalFunctionInterface<GridView>;
+ using Entity = typename Interface::Entity;
+ using LocalCoordinate = typename Interface::LocalCoordinate;
+
+ public:
+ /// Constructor. Stores a shared pointer to the passed Dune::VTKFunction
+ VTKLocalFunctionWrapper (std::shared_ptr<VTKFunction<GridView> const> const& fct)
+ : fct_(fct)
+ {}
+
+ /// Stores a pointer to the passed entity
+ virtual void bind (Entity const& entity) override
+ {
+ entity_ = &entity;
+ }
+
+ /// Unsets the stored entity pointer
+ virtual void unbind () override
+ {
+ entity_ = nullptr;
+ }
+
+ /// Evaluate the Dune::VTKFunction in LocalCoordinates on the stored Entity
+ virtual double evaluate (int comp, LocalCoordinate const& xi) const override
+ {
+ return fct_->evaluate(comp, *entity_, xi);
+ }
+
+ private:
+ std::shared_ptr<VTKFunction<GridView> const> fct_;
+ Entity const* entity_;
+ };
+
+ } //end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/localfunction.hh b/dune/vtk/localfunction.hh
new file mode 100644
index 0000000000000000000000000000000000000000..e8f8343f4483a47dfeef9067d9f24b1cf636ae4d
--- /dev/null
+++ b/dune/vtk/localfunction.hh
@@ -0,0 +1,74 @@
+#pragma once
+
+#include <memory>
+#include <type_traits>
+
+#include <dune/common/std/type_traits.hh>
+
+#include "localfunctioninterface.hh"
+#include "legacyvtkfunction.hh"
+#include "defaultvtkfunction.hh"
+
+namespace Dune
+{
+ namespace Vtk
+ {
+ /// \brief A Vtk::LocalFunction is a function-like object that can be bound to a grid element
+ /// an that provides an evaluate method with a component argument.
+ /**
+ * Stores internally a Vtk::LocalFunctionInterface object for the concrete evaluation.
+ **/
+ template <class GridView>
+ class LocalFunction
+ {
+ using Self = LocalFunction;
+ using Entity = typename GridView::template Codim<0>::Entity;
+ using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
+
+ template <class LF, class E>
+ using HasBind = decltype(std::declval<LF>().bind(std::declval<E>()));
+
+ public:
+ /// Construct the Vtk::LocalFunction from any function object that has a bind(element) method.
+ template <class LF,
+ disableCopyMove<Self, LF> = 0,
+ class = void_t<HasBind<LF,Entity>> >
+ LocalFunction (LF&& lf)
+ : localFct_(std::make_shared<LocalFunctionWrapper<GridView,LF>>(std::forward<LF>(lf)))
+ {}
+
+ /// Construct a Vtk::LocalFunction from a legacy VTKFunction
+ LocalFunction (std::shared_ptr<VTKFunction<GridView> const> const& lf)
+ : localFct_(std::make_shared<VTKLocalFunctionWrapper<GridView>>(lf))
+ {}
+
+ /// Allow the default construction of a Vtk::LocalFunction
+ LocalFunction () = default;
+
+ /// Bind the function to the grid entity
+ void bind (Entity const& entity)
+ {
+ assert(bool(localFct_));
+ localFct_->bind(entity);
+ }
+
+ /// Unbind from the currently bound entity
+ void unbind ()
+ {
+ assert(bool(localFct_));
+ localFct_->unbind();
+ }
+
+ /// Evaluate the `comp` component of the Range value at local coordinate `xi`
+ double evaluate (int comp, LocalCoordinate const& xi) const
+ {
+ assert(bool(localFct_));
+ return localFct_->evaluate(comp, xi);
+ }
+
+ private:
+ std::shared_ptr<LocalFunctionInterface<GridView>> localFct_ = nullptr;
+ };
+
+ } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/localfunctioninterface.hh b/dune/vtk/localfunctioninterface.hh
new file mode 100644
index 0000000000000000000000000000000000000000..69682e62d17b874c13adba2a38f7997c2a4b4c06
--- /dev/null
+++ b/dune/vtk/localfunctioninterface.hh
@@ -0,0 +1,30 @@
+#pragma once
+
+namespace Dune
+{
+ namespace Vtk
+ {
+ /// \brief An abstract base class for LocalFunctions that can be bound to an element and
+ /// evaluated in local coordinates w.r.t. to a component of its value.
+ template <class GridView>
+ class LocalFunctionInterface
+ {
+ public:
+ using Entity = typename GridView::template Codim<0>::Entity;
+ using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
+
+ /// Bind the function to the grid entity
+ virtual void bind (Entity const& entity) = 0;
+
+ /// Unbind from the currently bound entity
+ virtual void unbind () = 0;
+
+ /// Evaluate single component comp in the entity at local coordinates xi
+ virtual double evaluate (int comp, LocalCoordinate const& xi) const = 0;
+
+ /// Virtual destructor
+ virtual ~LocalFunctionInterface () = default;
+ };
+
+ } // end namespace Vtk
+} // end namespace Dune
diff --git a/dune/vtk/pvdwriter.hh b/dune/vtk/pvdwriter.hh
index 7b256cc802e8f292219e46177a555a1f11ffa1cf..a3d78ecb7250208247be225154fb9982112a86d0 100644
--- a/dune/vtk/pvdwriter.hh
+++ b/dune/vtk/pvdwriter.hh
@@ -6,7 +6,7 @@
#include <vector>
#include <tuple>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
@@ -15,7 +15,7 @@ namespace Dune
/// File-Writer for ParaView .pvd files
template <class VtkWriter>
class PvdWriter
- : public FileWriter
+ : public Vtk::FileWriter
{
using Self = PvdWriter;
@@ -53,7 +53,7 @@ namespace Dune
**/
virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
- /// Attach point data to the writer, \see VtkFunction for possible arguments
+ /// Attach point data to the writer, \see Vtk::Function for possible arguments
template <class Function, class... Args>
PvdWriter& addPointData (Function const& fct, Args&&... args)
{
@@ -61,7 +61,7 @@ namespace Dune
return *this;
}
- /// Attach cell data to the writer, \see VtkFunction for possible arguments
+ /// Attach cell data to the writer, \see Vtk::Function for possible arguments
template <class Function, class... Args>
PvdWriter& addCellData (Function const& fct, Args&&... args)
{
diff --git a/dune/vtk/pvdwriter.impl.hh b/dune/vtk/pvdwriter.impl.hh
index 7f84bd1ffdbd945fb99941010355c4eb3f410084..b26672083559d1f450ba896b8001406dbb13187b 100644
--- a/dune/vtk/pvdwriter.impl.hh
+++ b/dune/vtk/pvdwriter.impl.hh
@@ -12,13 +12,13 @@ template <class W>
void PvdWriter<W>
::writeTimestep (double time, std::string const& fn, std::optional<std::string> dir, bool writeCollection) const
{
- auto p = filesystem::path(fn);
+ auto p = Vtk::Path(fn);
auto name = p.stem();
- p.remove_filename();
+ p.removeFilename();
- filesystem::path fn_dir = p;
- filesystem::path data_dir = dir ? filesystem::path(*dir) : fn_dir;
- filesystem::path rel_dir = filesystem::relative(data_dir, fn_dir);
+ Vtk::Path fn_dir = p;
+ Vtk::Path data_dir = dir ? Vtk::Path(*dir) : fn_dir;
+ Vtk::Path rel_dir = Vtk::relative(data_dir, fn_dir);
std::string pvd_fn = fn_dir.string() + '/' + name.string();
std::string seq_fn = data_dir.string() + '/' + name.string() + "_t" + std::to_string(timesteps_.size());
@@ -52,9 +52,9 @@ template <class W>
std::string PvdWriter<W>
::write (std::string const& fn, std::optional<std::string> /*dir*/) const
{
- auto p = filesystem::path(fn);
+ auto p = Vtk::Path(fn);
auto name = p.stem();
- p.remove_filename();
+ p.removeFilename();
p /= name.string();
std::string outputFilename;
diff --git a/dune/vtk/test/test-typededuction.cc b/dune/vtk/test/test-typededuction.cc
index 0cc21c23430ac1bcdce0bf343f38209da0324ffe..23baf7ec93ad85bd0da53a70630d023bb3310393 100644
--- a/dune/vtk/test/test-typededuction.cc
+++ b/dune/vtk/test/test-typededuction.cc
@@ -28,7 +28,7 @@ int main (int argc, char** argv)
VtkUnstructuredGridWriter writer1(grid->leafGridView());
// 2. construct writer from datacollector
- ContinuousDataCollector dataCollector1(grid->leafGridView());
+ Vtk::ContinuousDataCollector dataCollector1(grid->leafGridView());
VtkUnstructuredGridWriter writer2(dataCollector1);
VtkUnstructuredGridWriter writer3(stackobject_to_shared_ptr(dataCollector1));
@@ -40,6 +40,6 @@ int main (int argc, char** argv)
VtkReader reader1(factory);
// 5. construct reader from grid-creator
- ContinuousGridCreator creator(factory);
+ Vtk::ContinuousGridCreator creator(factory);
VtkReader reader2(creator);
}
\ No newline at end of file
diff --git a/dune/vtk/vtktypes.cc b/dune/vtk/types.cc
similarity index 98%
rename from dune/vtk/vtktypes.cc
rename to dune/vtk/types.cc
index 90f993ec79f5442974039d02f6bac52435e850db..768408cac416e7beb209f5165940d27ef4f853bb 100644
--- a/dune/vtk/vtktypes.cc
+++ b/dune/vtk/types.cc
@@ -1,4 +1,4 @@
-#include "vtktypes.hh"
+#include "types.hh"
#include <iostream>
@@ -57,7 +57,7 @@ GeometryType to_geometry (std::uint8_t cell)
case QUADRATIC_QUAD: return GeometryTypes::quadrilateral;
case QUADRATIC_TETRA: return GeometryTypes::tetrahedron;
case QUADRATIC_HEXAHEDRON: return GeometryTypes::hexahedron;
-
+
// Arbitrary order Lagrange elements
case LAGRANGE_CURVE: return GeometryTypes::line;
case LAGRANGE_TRIANGLE: return GeometryTypes::triangle;
@@ -197,4 +197,4 @@ CellType::CellType (GeometryType const& t, CellParametrization parametrization)
}
}
-}} // end namespace Dune::Vtk
+} } // end namespace Dune::Vtk
diff --git a/dune/vtk/vtktypes.hh b/dune/vtk/types.hh
similarity index 100%
rename from dune/vtk/vtktypes.hh
rename to dune/vtk/types.hh
diff --git a/dune/vtk/utility/CMakeLists.txt b/dune/vtk/utility/CMakeLists.txt
index 048e6f91ecc2a129debba3daa4503fb5c9de6b0b..efed088eecbd4209a6ed113647ea3500efe5b064 100644
--- a/dune/vtk/utility/CMakeLists.txt
+++ b/dune/vtk/utility/CMakeLists.txt
@@ -6,8 +6,9 @@ install(FILES
enum.hh
filesystem.hh
lagrangepoints.hh
+ lagrangepoints.impl.hh
string.hh
uid.hh
- DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter/utility)
+ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtk/utility)
add_subdirectory(test)
\ No newline at end of file
diff --git a/dune/vtk/utility/enum.hh b/dune/vtk/utility/enum.hh
index 6c0b3fc16f71845d9215c2b25069ea336921b00d..025719e97b22448052c3de059483d4a474000d1c 100644
--- a/dune/vtk/utility/enum.hh
+++ b/dune/vtk/utility/enum.hh
@@ -4,11 +4,15 @@
namespace Dune
{
- template <class E, class Integer,
- std::enable_if_t<std::is_enum<E>::value, int> = 0>
- constexpr bool is_a(E a, Integer b)
+ namespace Vtk
{
- return (int(a) & int(b)) != 0;
- }
+ template <class E, class Integer,
+ std::enable_if_t<std::is_enum<E>::value, int> = 0>
+ constexpr bool is_a(E a, Integer b)
+ {
+ return (int(a) & int(b)) != 0;
+ }
+
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/utility/filesystem.cc b/dune/vtk/utility/filesystem.cc
index e8c4c06edf5e02dbd14d11ebb93f4fa694fd4fdc..42eb07cb48d16efb9396e42de4ff2ae9921ee9db 100644
--- a/dune/vtk/utility/filesystem.cc
+++ b/dune/vtk/utility/filesystem.cc
@@ -21,9 +21,10 @@
template <class... Args>
void inline _ignore_(Args&&...) {}
-namespace Dune { namespace filesystem {
+namespace Dune {
+namespace Vtk {
-std::string path::string() const
+std::string Path::string() const
{
if (empty())
return ".";
@@ -31,18 +32,18 @@ std::string path::string() const
auto it = begin();
auto result = *it;
for (++it; it != end(); ++it)
- result += preferred_separator + *it;
+ result += preferredSeparator + *it;
return result;
}
-void path::split(std::string p)
+void Path::split(std::string p)
{
std::string separators = "/\\";
bool relative = true;
trim(p);
- Dune::split(p.begin(), p.end(), separators.begin(), separators.end(),
+ Dune::Vtk::split(p.begin(), p.end(), separators.begin(), separators.end(),
[this,&relative](auto first, auto end)
{
auto token = std::string(first, end);
@@ -62,7 +63,7 @@ void path::split(std::string p)
}
-path path::stem() const
+Path Path::stem() const
{
auto f = filename().string();
auto pos = f.find_last_of('.');
@@ -73,7 +74,7 @@ path path::stem() const
}
-path path::extension() const
+Path Path::extension() const
{
auto f = filename().string();
auto pos = f.find_last_of('.');
@@ -84,7 +85,7 @@ path path::extension() const
}
-bool path::is_absolute(std::string p)
+bool Path::isAbsolute(std::string p)
{
if (p[0] == '/')
return true;
@@ -97,15 +98,15 @@ bool path::is_absolute(std::string p)
}
-path& path::operator/=(path const& p)
+Path& Path::operator/=(Path const& p)
{
insert(end(), p.begin(), p.end());
- original += preferred_separator + p.original;
+ original += preferredSeparator + p.original;
return *this;
}
-bool path::is_file() const
+bool Path::isFile() const
{
std::string p = this->string();
struct stat info;
@@ -113,7 +114,7 @@ bool path::is_file() const
}
-bool path::is_directory() const
+bool Path::isDirectory() const
{
std::string p = this->string();
struct stat info;
@@ -121,7 +122,7 @@ bool path::is_directory() const
}
-path current_path()
+Path currentPath()
{
char cwd_[FILENAME_MAX];
_ignore_(GET_CURRENT_DIR(cwd_, sizeof(cwd_)));
@@ -130,20 +131,20 @@ path current_path()
}
-bool exists(path const& p)
+bool exists(Path const& p)
{
- return p.is_file() || p.is_directory();
+ return p.isFile() || p.isDirectory();
}
-bool create_directories(path const& p)
+bool createDirectories(Path const& p)
{
- if (p.is_directory())
+ if (p.isDirectory())
return true;
- auto parent = p.parent_path();
- if (!parent.empty() && !parent.is_directory())
- create_directories(parent);
+ auto parent = p.parentPath();
+ if (!parent.empty() && !parent.isDirectory())
+ createDirectories(parent);
#ifdef _WIN32
int ret = _mkdir(p.string().c_str());
@@ -169,7 +170,7 @@ bool create_directories(path const& p)
}
}
-path relative(path const& a, path const& b)
+Path relative(Path const& a, Path const& b)
{
// find common base path
auto a_it = a.begin();
@@ -180,11 +181,11 @@ path relative(path const& a, path const& b)
}
// combine remaining parts of a to result path
- path rel(".");
+ Path rel(".");
for (; a_it != a.end(); ++a_it)
rel /= *a_it;
return rel;
}
-} } // end namespace Dune::filesystem
+} } // end namespace Dune::Vtk
diff --git a/dune/vtk/utility/filesystem.hh b/dune/vtk/utility/filesystem.hh
index 12aacfd31fe6a614b079e3ba149f51fd0cf1c9c4..827273a2ceef956bb1bbcf37fe2cf27a1d939fb7 100644
--- a/dune/vtk/utility/filesystem.hh
+++ b/dune/vtk/utility/filesystem.hh
@@ -7,10 +7,10 @@
namespace Dune
{
- namespace filesystem
+ namespace Vtk
{
// A minimalistic filesystem class
- class path
+ class Path
: public std::vector<std::string>
{
using Super = std::vector<std::string>;
@@ -19,58 +19,58 @@ namespace Dune
public:
#ifdef _WIN32
- static constexpr char preferred_separator = '\\';
+ static constexpr char preferredSeparator = '\\';
#else
- static constexpr char preferred_separator = '/';
+ static constexpr char preferredSeparator = '/';
#endif
public:
- path() = default;
+ Path() = default;
// NOTE: implicit conversion is allowed here
template <class String>
- path(String const& p)
+ Path(String const& p)
: original(p)
{
split(p);
}
template <class InputIt>
- path(InputIt it, InputIt end_it)
+ Path(InputIt it, InputIt end_it)
: Super(it, end_it)
{
original = this->string();
}
template <class String>
- path(std::initializer_list<String> const& list)
- : path(list.begin(), list.end())
+ Path(std::initializer_list<String> const& list)
+ : Path(list.begin(), list.end())
{}
/// Removes filename path component
- path& remove_filename()
+ Path& removeFilename()
{
this->pop_back();
return *this;
}
/// Returns the path of the parent path
- path parent_path() const
+ Path parentPath() const
{
- return empty() ? path() : path(begin(), --end());
+ return empty() ? Path() : Path(begin(), --end());
}
/// Returns filename path component
- path filename() const
+ Path filename() const
{
- return empty() ? path() : path(back());
+ return empty() ? Path() : Path(back());
}
/// Returns the stem path component
- path stem() const;
+ Path stem() const;
/// Returns the file extension path component
- path extension() const;
+ Path extension() const;
/// Return the path as string
std::string string() const;
@@ -79,30 +79,30 @@ namespace Dune
/** In Linux, test whether the path starts with `/`, in Windows whether it starts
* with `[a-z]:\\`.
**/
- static bool is_absolute(std::string p);
+ static bool isAbsolute(std::string p);
- bool is_absolute() const { return is_absolute(original); }
+ bool isAbsolute() const { return isAbsolute(original); }
- bool is_relative() const { return !is_absolute(); }
+ bool isRelative() const { return !isAbsolute(); }
/// Check whether path is a regular file
- bool is_file() const;
+ bool isFile() const;
/// Check whether path is a regular file
- bool is_directory() const;
+ bool isDirectory() const;
/// Lexicographically compares two paths
- bool operator==(path const& p)
+ bool operator==(Path const& p)
{
return this->string() == p.string();
}
/// Appends elements to the path
- path& operator/=(path const& p);
+ Path& operator/=(Path const& p);
/// output of the path
template <class CharT, class Traits>
- friend std::basic_ostream<CharT, Traits>& operator<<(std::basic_ostream<CharT, Traits>& out, path const& p)
+ friend std::basic_ostream<CharT, Traits>& operator<<(std::basic_ostream<CharT, Traits>& out, Path const& p)
{
out << '"' << p.string() << '"';
return out;
@@ -119,16 +119,16 @@ namespace Dune
};
/// Test whether the path is a valid (existing and accessible) file / directory
- bool exists(path const&);
+ bool exists(Path const&);
/// Create directory and non existing parent directories.
- bool create_directories(path const&);
+ bool createDirectories(Path const&);
/// Returns the current path
- path current_path();
+ Path currentPath();
/// Find the path of `a` relative to directory of `b`
- path relative(path const& a, path const& b);
+ Path relative(Path const& a, Path const& b);
- } // end namespace filesystem
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/utility/lagrangepoints.hh b/dune/vtk/utility/lagrangepoints.hh
index b598c89527500007b91af4a8a4cfc46316678826..ebda327d43c92e1cf3f92ae7305e0dcb5a9e1807 100644
--- a/dune/vtk/utility/lagrangepoints.hh
+++ b/dune/vtk/utility/lagrangepoints.hh
@@ -7,599 +7,193 @@
#include <dune/geometry/type.hh>
#include <dune/localfunctions/lagrange/equidistantpoints.hh>
-namespace Dune {
-
-namespace Impl {
- // forward declaration
- template <class K, unsigned int dim>
- class VtkLagrangePointSetBuilder;
-}
-
-
-/// \brief A set of lagrange points compatible with the numbering of VTK
-/**
- * \tparam K Field-type for the coordinates
- * \tparam dim Dimension of the coordinates
- **/
-template <class K, unsigned int dim>
-class VtkLagrangePointSet
- : public EmptyPointSet<K, dim>
+namespace Dune
{
- using Super = EmptyPointSet<K, dim>;
-
-public:
- static const unsigned int dimension = dim;
-
- VtkLagrangePointSet (std::size_t order)
- : Super(order)
- {
- assert(order > 0);
- }
-
- /// Fill the lagrange points for the given geometry type
- void build (GeometryType gt)
- {
- assert(gt.dim() == dimension);
- builder_(gt, order(), points_);
- }
-
- /// Fill the lagrange points for the given topology type `Topology`
- template <class Topology>
- bool build ()
- {
- build(GeometryType(Topology{}));
- return true;
- }
-
- /// Returns whether the point set support the given topology type `Topology` and can
- /// generate point for the given order.
- template <class Topology>
- static bool supports (std::size_t order)
- {
- return true;
- }
-
- using Super::order;
-
-private:
- using Super::points_;
- Impl::VtkLagrangePointSetBuilder<K,dim> builder_;
-};
-
-
-namespace Impl {
-
-// Build for lagrange point sets in different dimensions
-// Specialized for dim=1,2,3
-template <class K, unsigned int dim>
-class VtkLagrangePointSetBuilder
-{
-public:
- template <class Points>
- void operator()(GeometryType, unsigned int, Points& points) const
- {
- DUNE_THROW(Dune::NotImplemented,
- "Lagrange points not yet implemented for this GeometryType.");
- }
-};
-
-/**
- * The implementation of the point set builder is directly derived from VTK.
- * Modification are a change in data-types and naming scheme. Additionally
- * a LocalKey is created to reflect the concept of a Dune PointSet.
- *
- * Included is the license of the BSD 3-clause License included in the VTK
- * source code from 2020/04/13 in commit b90dad558ce28f6d321420e4a6b17e23f5227a1c
- * of git repository https://gitlab.kitware.com/vtk/vtk.
- *
- Program: Visualization Toolkit
- Module: Copyright.txt
-
- Copyright (c) 1993-2015 Ken Martin, Will Schroeder, Bill Lorensen
- All rights reserved.
-
- Redistribution and use in source and binary forms, with or without
- modification, are permitted provided that the following conditions are met:
-
- * Redistributions of source code must retain the above copyright notice,
- this list of conditions and the following disclaimer.
-
- * Redistributions in binary form must reproduce the above copyright notice,
- this list of conditions and the following disclaimer in the documentation
- and/or other materials provided with the distribution.
-
- * Neither name of Ken Martin, Will Schroeder, or Bill Lorensen nor the names
- of any contributors may be used to endorse or promote products derived
- from this software without specific prior written permission.
-
- THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS ``AS IS''
- AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS OR CONTRIBUTORS BE LIABLE FOR
- ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
- SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
- CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
- OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
- OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
- *
- **/
-
-// Lagrange points on point geometries
-template <class K>
-class VtkLagrangePointSetBuilder<K,0>
-{
- static constexpr int dim = 0;
- using LP = LagrangePoint<K,dim>;
- using Vec = typename LP::Vector;
- using Key = LocalKey;
-
-public:
- template <class Points>
- void operator()(GeometryType gt, int /*order*/, Points& points) const
- {
- assert(gt.isVertex());
- points.push_back(LP{Vec{},Key{0,0,0}});
- }
-};
-
-
-// Lagrange points on line geometries
-template <class K>
-class VtkLagrangePointSetBuilder<K,1>
-{
- static constexpr int dim = 1;
- using LP = LagrangePoint<K,dim>;
- using Vec = typename LP::Vector;
- using Key = LocalKey;
-
-public:
- template <class Points>
- void operator()(GeometryType gt, int order, Points& points) const
- {
- assert(gt.isLine());
-
- // Vertex nodes
- points.push_back(LP{Vec{0.0}, Key{0,dim,0}});
- points.push_back(LP{Vec{1.0}, Key{1,dim,0}});
-
- if (order > 1) {
- // Inner nodes
- Vec p{0.0};
- for (unsigned int i = 0; i < order-1; i++)
- {
- p[0] += 1.0 / order;
- points.push_back(LP{p,Key{0,dim-1,i}});
- }
- }
- }
-};
-
-
-// Lagrange points on 2d geometries
-template <class K>
-class VtkLagrangePointSetBuilder<K,2>
-{
- static constexpr int dim = 2;
- using LP = LagrangePoint<K,dim>;
- using Vec = typename LP::Vector;
- using Key = LocalKey;
-
- friend class VtkLagrangePointSetBuilder<K,3>;
-
-public:
- template <class Points>
- void operator()(GeometryType gt, int order, Points& points) const
- {
- std::size_t nPoints = numLagrangePoints(gt.id(), dim, order);
-
- if (gt.isTriangle())
- buildTriangle(nPoints, order, points);
- else if (gt.isQuadrilateral())
- buildQuad(nPoints, order, points);
- else {
- DUNE_THROW(Dune::NotImplemented,
- "Lagrange points not yet implemented for this GeometryType.");
- }
-
- assert(points.size() == nPoints);
- }
-
-private:
- // Construct the point set in a triangle element.
- // Loop from the outside to the inside
- template <class Points>
- void buildTriangle (std::size_t nPoints, int order, Points& points) const
+ namespace Vtk
{
- points.reserve(nPoints);
-
- const int nVertexDOFs = 3;
- const int nEdgeDOFs = 3 * (order-1);
-
- static const unsigned int vertexPerm[3] = {0,1,2};
- static const unsigned int edgePerm[3] = {0,2,1};
-
- auto calcKey = [&](int p) -> Key
+ namespace Impl
{
- if (p < nVertexDOFs) {
- return Key{vertexPerm[p], dim, 0};
- }
- else if (p < nVertexDOFs+nEdgeDOFs) {
- unsigned int entityIndex = (p - nVertexDOFs) / (order-1);
- unsigned int index = (p - nVertexDOFs) % (order-1);
- return Key{edgePerm[entityIndex], dim-1, index};
- }
- else {
- unsigned int index = p - (nVertexDOFs + nEdgeDOFs);
- return Key{0, dim-2, index};
- }
- };
-
- std::array<int,3> bindex;
-
- double order_d = double(order);
- for (std::size_t p = 0; p < nPoints; ++p) {
- barycentricIndex(p, bindex, order);
- Vec point{bindex[0] / order_d, bindex[1] / order_d};
- points.push_back(LP{point, calcKey(p)});
+ // forward declaration
+ template <class K, unsigned int dim>
+ class LagrangePointSetBuilder;
}
- }
-
- // "Barycentric index" is a triplet of integers, each running from 0 to
- // <Order>. It is the index of a point on the triangle in barycentric
- // coordinates.
- static void barycentricIndex (int index, std::array<int,3>& bindex, int order)
- {
- int max = order;
- int min = 0;
- // scope into the correct triangle
- while (index != 0 && index >= 3 * order)
- {
- index -= 3 * order;
- max -= 2;
- min++;
- order -= 3;
- }
- // vertex DOFs
- if (index < 3)
- {
- bindex[index] = bindex[(index + 1) % 3] = min;
- bindex[(index + 2) % 3] = max;
- }
- // edge DOFs
- else
+ /// \brief A set of lagrange points compatible with the numbering of VTK and Gmsh
+ /**
+ * \tparam K Field-type for the coordinates
+ * \tparam dim Dimension of the coordinates
+ **/
+ template <class K, unsigned int dim>
+ class LagrangePointSet
+ : public EmptyPointSet<K, dim>
{
- index -= 3;
- int dim = index / (order - 1);
- int offset = (index - dim * (order - 1));
- bindex[(dim + 1) % 3] = min;
- bindex[(dim + 2) % 3] = (max - 1) - offset;
- bindex[dim] = (min + 1) + offset;
- }
- }
-
-
- // Construct the point set in the quad element
- // 1. build equispaced points with index tuple (i,j)
- // 2. map index tuple to DOF index and LocalKey
- template <class Points>
- void buildQuad(std::size_t nPoints, int order, Points& points) const
- {
- points.resize(nPoints);
+ using Super = EmptyPointSet<K, dim>;
- std::array<int,2> orders{order, order};
- std::array<Vec,4> nodes{Vec{0., 0.}, Vec{1., 0.}, Vec{1., 1.}, Vec{0., 1.}};
+ public:
+ static const unsigned int dimension = dim;
- for (int n = 0; n <= orders[1]; ++n) {
- for (int m = 0; m <= orders[0]; ++m) {
- // int idx = pointIndexFromIJ(m,n,orders);
-
- const double r = double(m) / orders[0];
- const double s = double(n) / orders[1];
- Vec p = (1.0 - r) * (nodes[3] * s + nodes[0] * (1.0 - s)) +
- r * (nodes[2] * s + nodes[1] * (1.0 - s));
-
- auto [idx,key] = calcQuadKey(m,n,orders);
- points[idx] = LP{p, key};
- // points[idx] = LP{p, calcQuadKey(n,m,orders)};
+ LagrangePointSet (std::size_t order)
+ : Super(order)
+ {
+ assert(order > 0);
}
- }
- }
- // Obtain the VTK DOF index of the node (i,j) in the quad element
- // and construct a LocalKey
- static std::pair<int,Key> calcQuadKey (int i, int j, std::array<int,2> order)
- {
- const bool ibdy = (i == 0 || i == order[0]);
- const bool jbdy = (j == 0 || j == order[1]);
- const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0); // How many boundaries do we lie on at once?
-
- int dof = 0;
- unsigned int entityIndex = 0;
- unsigned int index = 0;
-
- if (nbdy == 2) // Vertex DOF
- {
- dof = (i ? (j ? 2 : 1) : (j ? 3 : 0));
- entityIndex = (j ? (i ? 3 : 2) : (i ? 1 : 0));
- return std::make_pair(dof,Key{entityIndex, dim, 0});
- }
+ /// Fill the lagrange points for the given geometry type
+ void build (GeometryType gt)
+ {
+ assert(gt.dim() == dimension);
+ builder_(gt, order(), points_);
+ }
- int offset = 4;
- if (nbdy == 1) // Edge DOF
- {
- if (!ibdy) {
- dof = (i - 1) + (j ? order[0]-1 + order[1]-1 : 0) + offset;
- entityIndex = j ? 3 : 2;
- index = i-1;
+ /// Fill the lagrange points for the given topology type `Topology`
+ template <class Topology>
+ bool build ()
+ {
+ build(GeometryType(Topology{}));
+ return true;
}
- else if (!jbdy) {
- dof = (j - 1) + (i ? order[0]-1 : 2 * (order[0]-1) + order[1]-1) + offset;
- entityIndex = i ? 1 : 0;
- index = j-1;
+
+ /// Returns whether the point set support the given topology type `Topology` and can
+ /// generate point for the given order.
+ template <class Topology>
+ static bool supports (std::size_t order)
+ {
+ return true;
}
- return std::make_pair(dof, Key{entityIndex, dim-1, index});
- }
- offset += 2 * (order[0]-1 + order[1]-1);
+ using Super::order;
- // nbdy == 0: Face DOF
- dof = offset + (i - 1) + (order[0]-1) * ((j - 1));
- Key innerKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
- return std::make_pair(dof, Key{0, dim-2, innerKey.index()});
- }
-};
+ private:
+ using Super::points_;
+ Impl::LagrangePointSetBuilder<K,dim> builder_;
+ };
-// Lagrange points on 3d geometries
-template <class K>
-class VtkLagrangePointSetBuilder<K,3>
-{
- static constexpr int dim = 3;
- using LP = LagrangePoint<K,dim>;
- using Vec = typename LP::Vector;
- using Key = LocalKey;
-
-public:
- template <class Points>
- void operator() (GeometryType gt, unsigned int order, Points& points) const
- {
- std::size_t nPoints = numLagrangePoints(gt.id(), dim, order);
-
- if (gt.isTetrahedron())
- buildTetra(nPoints, order, points);
- else if (gt.isHexahedron())
- buildHex(nPoints, order, points);
- else {
- DUNE_THROW(Dune::NotImplemented,
- "Lagrange points not yet implemented for this GeometryType.");
- }
+ namespace Impl
+ {
+ // Build for lagrange point sets in different dimensions
+ // Specialized for dim=1,2,3
+ template <class K, unsigned int dim>
+ class LagrangePointSetBuilder
+ {
+ public:
+ template <class Points>
+ void operator()(GeometryType, unsigned int, Points& points) const
+ {
+ DUNE_THROW(Dune::NotImplemented,
+ "Lagrange points not yet implemented for this GeometryType.");
+ }
+ };
- assert(points.size() == nPoints);
- }
-private:
- // Construct the point set in the tetrahedron element
- // 1. construct barycentric (index) coordinates
- // 2. obtains the DOF index, LocalKey and actual coordinate from barycentric index
- template <class Points>
- void buildTetra (std::size_t nPoints, int order, Points& points) const
- {
- points.reserve(nPoints);
+ // Lagrange points on point geometries
+ template <class K>
+ class LagrangePointSetBuilder<K,0>
+ {
+ static constexpr int dim = 0;
+ using LP = LagrangePoint<K,dim>;
+ using Vec = typename LP::Vector;
+ using Key = LocalKey;
- const int nVertexDOFs = 4;
- const int nEdgeDOFs = 6 * (order-1);
- const int nFaceDOFs = 4 * (order-1)*(order-2)/2;
+ public:
+ template <class Points>
+ void operator()(GeometryType gt, int /*order*/, Points& points) const;
+ };
- static const unsigned int vertexPerm[4] = {0,1,2,3};
- static const unsigned int edgePerm[6] = {0,2,1,3,4,5};
- static const unsigned int facePerm[4] = {1,2,0,3};
- auto calcKey = [&](int p) -> Key
- {
- if (p < nVertexDOFs) {
- return Key{vertexPerm[p], dim, 0};
- }
- else if (p < nVertexDOFs+nEdgeDOFs) {
- unsigned int entityIndex = (p - nVertexDOFs) / (order-1);
- unsigned int index = (p - nVertexDOFs) % (order-1);
- return Key{edgePerm[entityIndex], dim-1, index};
- }
- else if (p < nVertexDOFs+nEdgeDOFs+nFaceDOFs) {
- unsigned int index = (p - (nVertexDOFs + nEdgeDOFs)) % ((order-1)*(order-2)/2);
- unsigned int entityIndex = (p - (nVertexDOFs + nEdgeDOFs)) / ((order-1)*(order-2)/2);
- return Key{facePerm[entityIndex], dim-2, index};
- }
- else {
- unsigned int index = p - (nVertexDOFs + nEdgeDOFs + nFaceDOFs);
- return Key{0, dim-3, index};
- }
- };
+ // Lagrange points on line geometries
+ template <class K>
+ class LagrangePointSetBuilder<K,1>
+ {
+ static constexpr int dim = 1;
+ using LP = LagrangePoint<K,dim>;
+ using Vec = typename LP::Vector;
+ using Key = LocalKey;
- std::array<int,4> bindex;
-
- double order_d = double(order);
- for (std::size_t p = 0; p < nPoints; ++p) {
- barycentricIndex(p, bindex, order);
- Vec point{bindex[0] / order_d, bindex[1] / order_d, bindex[2] / order_d};
- points.push_back(LP{point, calcKey(p)});
- }
- }
+ public:
+ template <class Points>
+ void operator()(GeometryType gt, int order, Points& points) const;
+ };
- // "Barycentric index" is a set of 4 integers, each running from 0 to
- // <Order>. It is the index of a point in the tetrahedron in barycentric
- // coordinates.
- static void barycentricIndex (std::size_t p, std::array<int,4>& bindex, int order)
- {
- const int nVertexDOFs = 4;
- const int nEdgeDOFs = 6 * (order-1);
- static const int edgeVertices[6][2] = {{0,1}, {1,2}, {2,0}, {0,3}, {1,3}, {2,3}};
- static const int linearVertices[4][4] = {{0,0,0,1}, {1,0,0,0}, {0,1,0,0}, {0,0,1,0}};
- static const int vertexMaxCoords[4] = {3,0,1,2};
- static const int faceBCoords[4][3] = {{0,2,3}, {2,0,1}, {2,1,3}, {1,0,3}};
- static const int faceMinCoord[4] = {1,3,0,2};
+ // Lagrange points on 2d geometries
+ template <class K>
+ class LagrangePointSetBuilder<K,2>
+ {
+ static constexpr int dim = 2;
+ using LP = LagrangePoint<K,dim>;
+ using Vec = typename LP::Vector;
+ using Key = LocalKey;
- int max = order;
- int min = 0;
+ friend class LagrangePointSetBuilder<K,3>;
- // scope into the correct tetra
- while (p >= 2 * (order * order + 1) && p != 0 && order > 3)
- {
- p -= 2 * (order * order + 1);
- max -= 3;
- min++;
- order -= 4;
- }
+ public:
+ template <class Points>
+ void operator()(GeometryType gt, int order, Points& points) const;
- // vertex DOFs
- if (p < nVertexDOFs)
- {
- for (int coord = 0; coord < 4; ++coord)
- bindex[coord] = (coord == vertexMaxCoords[p] ? max : min);
- }
- // edge DOFs
- else if (p < nVertexDOFs+nEdgeDOFs)
- {
- int edgeId = (p - nVertexDOFs) / (order-1);
- int vertexId = (p - nVertexDOFs) % (order-1);
- for (int coord = 0; coord < 4; ++coord)
- {
- bindex[coord] = min +
- (linearVertices[edgeVertices[edgeId][0]][coord] * (max - min - 1 - vertexId) +
- linearVertices[edgeVertices[edgeId][1]][coord] * (1 + vertexId));
- }
- }
- // face DOFs
- else
- {
- int faceId = (p - (nVertexDOFs+nEdgeDOFs)) / ((order-2)*(order-1)/2);
- int vertexId = (p - (nVertexDOFs+nEdgeDOFs)) % ((order-2)*(order-1)/2);
-
- std::array<int,3> projectedBIndex;
- if (order == 3)
- projectedBIndex[0] = projectedBIndex[1] = projectedBIndex[2] = 0;
- else
- VtkLagrangePointSetBuilder<K,2>::barycentricIndex(vertexId, projectedBIndex, order-3);
-
- for (int i = 0; i < 3; i++)
- bindex[faceBCoords[faceId][i]] = (min + 1 + projectedBIndex[i]);
-
- bindex[faceMinCoord[faceId]] = min;
- }
- }
-
-private:
- // Construct the point set in the heyhedral element
- // 1. build equispaced points with index tuple (i,j,k)
- // 2. map index tuple to DOF index and LocalKey
- template <class Points>
- void buildHex (std::size_t nPoints, int order, Points& points) const
- {
- points.resize(nPoints);
-
- std::array<int,3> orders{order, order, order};
- std::array<Vec,8> nodes{Vec{0., 0., 0.}, Vec{1., 0., 0.}, Vec{1., 1., 0.}, Vec{0., 1., 0.},
- Vec{0., 0., 1.}, Vec{1., 0., 1.}, Vec{1., 1., 1.}, Vec{0., 1., 1.}};
-
- for (int p = 0; p <= orders[2]; ++p) {
- for (int n = 0; n <= orders[1]; ++n) {
- for (int m = 0; m <= orders[0]; ++m) {
- const double r = double(m) / orders[0];
- const double s = double(n) / orders[1];
- const double t = double(p) / orders[2];
- Vec point = (1.0-r) * ((nodes[3] * (1.0-t) + nodes[7] * t) * s + (nodes[0] * (1.0-t) + nodes[4] * t) * (1.0-s)) +
- r * ((nodes[2] * (1.0-t) + nodes[6] * t) * s + (nodes[1] * (1.0-t) + nodes[5] * t) * (1.0-s));
-
- auto [idx,key] = calcHexKey(m,n,p,orders);
- points[idx] = LP{point, key};
- }
- }
- }
- }
+ private: // implementation details
- // Obtain the VTK DOF index of the node (i,j,k) in the hexahedral element
- static std::pair<int,Key> calcHexKey (int i, int j, int k, std::array<int,3> order)
- {
- const bool ibdy = (i == 0 || i == order[0]);
- const bool jbdy = (j == 0 || j == order[1]);
- const bool kbdy = (k == 0 || k == order[2]);
- const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0) + (kbdy ? 1 : 0); // How many boundaries do we lie on at once?
+ // Construct the point set in a triangle element.
+ // Loop from the outside to the inside
+ template <class Points>
+ void buildTriangle (std::size_t nPoints, int order, Points& points) const;
- int dof = 0;
- unsigned int entityIndex = 0;
- unsigned int index = 0;
+ // "Barycentric index" is a triplet of integers, each running from 0 to
+ // <Order>. It is the index of a point on the triangle in barycentric
+ // coordinates.
+ static void barycentricIndex (int index, std::array<int,3>& bindex, int order);
- if (nbdy == 3) // Vertex DOF
- {
- dof = (i ? (j ? 2 : 1) : (j ? 3 : 0)) + (k ? 4 : 0);
- entityIndex = (i ? 1 : 0) + (j ? 2 : 0) + (k ? 4 : 0);
- return std::make_pair(dof, Key{entityIndex, dim, 0});
- }
+ // Construct the point set in the quad element
+ // 1. build equispaced points with index tuple (i,j)
+ // 2. map index tuple to DOF index and LocalKey
+ template <class Points>
+ void buildQuad(std::size_t nPoints, int order, Points& points) const;
- int offset = 8;
- if (nbdy == 2) // Edge DOF
- {
- entityIndex = (k ? 8 : 4);
- if (!ibdy)
- { // On i axis
- dof = (i - 1) + (j ? order[0]-1 + order[1]-1 : 0) + (k ? 2 * (order[0]-1 + order[1]-1) : 0) + offset;
- index = i;
- entityIndex += (i ? 1 : 0);
- }
- else if (!jbdy)
- { // On j axis
- dof = (j - 1) + (i ? order[0]-1 : 2 * (order[0]-1) + order[1]-1) + (k ? 2 * (order[0]-1 + order[1]-1) : 0) + offset;
- index = j;
- entityIndex += (j ? 3 : 2);
- }
- else
- { // !kbdy, On k axis
- offset += 4 * (order[0]-1) + 4 * (order[1]-1);
- dof = (k - 1) + (order[2]-1) * (i ? (j ? 3 : 1) : (j ? 2 : 0)) + offset;
- index = k;
- entityIndex = (i ? 1 : 0) + (j ? 2 : 0);
- }
- return std::make_pair(dof, Key{entityIndex, dim-1, index});
- }
-
- offset += 4 * (order[0]-1 + order[1]-1 + order[2]-1);
- if (nbdy == 1) // Face DOF
- {
- Key faceKey;
- if (ibdy) // On i-normal face
- {
- dof = (j - 1) + ((order[1]-1) * (k - 1)) + (i ? (order[1]-1) * (order[2]-1) : 0) + offset;
- entityIndex = (i ? 1 : 0);
- faceKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(j-1,k-1,{order[1]-2, order[2]-2}).second;
- }
- else {
- offset += 2 * (order[1] - 1) * (order[2] - 1);
- if (jbdy) // On j-normal face
- {
- dof = (i - 1) + ((order[0]-1) * (k - 1)) + (j ? (order[2]-1) * (order[0]-1) : 0) + offset;
- entityIndex = (j ? 3 : 2);
- faceKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(i-1,k-1,{order[0]-2, order[2]-2}).second;
- }
- else
- { // kbdy, On k-normal face
- offset += 2 * (order[2]-1) * (order[0]-1);
- dof = (i - 1) + ((order[0]-1) * (j - 1)) + (k ? (order[0]-1) * (order[1]-1) : 0) + offset;
- entityIndex = (k ? 5 : 4);
- faceKey = VtkLagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
- }
- }
- return std::make_pair(dof, Key{entityIndex, dim-2, faceKey.index()});
- }
+ // Obtain the VTK DOF index of the node (i,j) in the quad element
+ // and construct a LocalKey
+ static std::pair<int,Key> calcQuadKey (int i, int j, std::array<int,2> order);
+ };
- // nbdy == 0: Body DOF
- offset += 2 * ((order[1]-1) * (order[2]-1) + (order[2]-1) * (order[0]-1) + (order[0]-1) * (order[1]-1));
- dof = offset + (i - 1) + (order[0]-1) * ((j - 1) + (order[1]-1) * ((k - 1)));
- Key innerKey = VtkLagrangePointSetBuilder<K,3>::calcHexKey(i-1,j-1,k-1,{order[0]-2, order[1]-2, order[2]-2}).second;
- return std::make_pair(dof, Key{0, dim-3, innerKey.index()});
- }
-};
-}} // end namespace Dune::Impl
+ // Lagrange points on 3d geometries
+ template <class K>
+ class LagrangePointSetBuilder<K,3>
+ {
+ static constexpr int dim = 3;
+ using LP = LagrangePoint<K,dim>;
+ using Vec = typename LP::Vector;
+ using Key = LocalKey;
+
+ public:
+ template <class Points>
+ void operator() (GeometryType gt, unsigned int order, Points& points) const;
+
+ private: // implementation details
+
+ // Construct the point set in the tetrahedron element
+ // 1. construct barycentric (index) coordinates
+ // 2. obtains the DOF index, LocalKey and actual coordinate from barycentric index
+ template <class Points>
+ void buildTetra (std::size_t nPoints, int order, Points& points) const;
+
+ // "Barycentric index" is a set of 4 integers, each running from 0 to
+ // <Order>. It is the index of a point in the tetrahedron in barycentric
+ // coordinates.
+ static void barycentricIndex (std::size_t p, std::array<int,4>& bindex, int order);
+
+ // Construct the point set in the heyhedral element
+ // 1. build equispaced points with index tuple (i,j,k)
+ // 2. map index tuple to DOF index and LocalKey
+ template <class Points>
+ void buildHex (std::size_t nPoints, int order, Points& points) const;
+
+ // Obtain the VTK DOF index of the node (i,j,k) in the hexahedral element
+ static std::pair<int,Key> calcHexKey (int i, int j, int k, std::array<int,3> order);
+ };
+
+ } // end namespace Impl
+ } // end namespace Vtk
+} // end namespace Dune
+
+#include <dune/vtk/utility/lagrangepoints.impl.hh>
diff --git a/dune/vtk/utility/lagrangepoints.impl.hh b/dune/vtk/utility/lagrangepoints.impl.hh
new file mode 100644
index 0000000000000000000000000000000000000000..e49934513c16e2161bdebad4141c597e9e7b341b
--- /dev/null
+++ b/dune/vtk/utility/lagrangepoints.impl.hh
@@ -0,0 +1,507 @@
+#pragma once
+
+#include <cassert>
+#include <array>
+
+#include <dune/common/exceptions.hh>
+#include <dune/geometry/type.hh>
+#include <dune/localfunctions/lagrange/equidistantpoints.hh>
+
+namespace Dune {
+namespace Vtk {
+namespace Impl {
+
+/**
+ * The implementation of the point set builder is directly derived from VTK.
+ * Modification are a change in data-types and naming scheme. Additionally
+ * a LocalKey is created to reflect the concept of a Dune PointSet.
+ *
+ * Included is the license of the BSD 3-clause License included in the VTK
+ * source code from 2020/04/13 in commit b90dad558ce28f6d321420e4a6b17e23f5227a1c
+ * of git repository https://gitlab.kitware.com/vtk/vtk.
+ *
+ Program: Visualization Toolkit
+ Module: Copyright.txt
+
+ Copyright (c) 1993-2015 Ken Martin, Will Schroeder, Bill Lorensen
+ All rights reserved.
+
+ Redistribution and use in source and binary forms, with or without
+ modification, are permitted provided that the following conditions are met:
+
+ * Redistributions of source code must retain the above copyright notice,
+ this list of conditions and the following disclaimer.
+
+ * Redistributions in binary form must reproduce the above copyright notice,
+ this list of conditions and the following disclaimer in the documentation
+ and/or other materials provided with the distribution.
+
+ * Neither name of Ken Martin, Will Schroeder, or Bill Lorensen nor the names
+ of any contributors may be used to endorse or promote products derived
+ from this software without specific prior written permission.
+
+ THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS ``AS IS''
+ AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+ IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
+ ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHORS OR CONTRIBUTORS BE LIABLE FOR
+ ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+ DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+ SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+ CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+ OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+ *
+ **/
+
+
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,0>::operator() (GeometryType gt, int /*order*/, Points& points) const
+{
+ assert(gt.isVertex());
+ points.push_back(LP{Vec{},Key{0,0,0}});
+}
+
+
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,1>::operator() (GeometryType gt, int order, Points& points) const
+{
+ assert(gt.isLine());
+
+ // Vertex nodes
+ points.push_back(LP{Vec{0.0}, Key{0,dim,0}});
+ points.push_back(LP{Vec{1.0}, Key{1,dim,0}});
+
+ if (order > 1) {
+ // Inner nodes
+ Vec p{0.0};
+ for (unsigned int i = 0; i < order-1; i++)
+ {
+ p[0] += 1.0 / order;
+ points.push_back(LP{p,Key{0,dim-1,i}});
+ }
+ }
+}
+
+
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,2>::operator() (GeometryType gt, int order, Points& points) const
+{
+ std::size_t nPoints = numLagrangePoints(gt.id(), dim, order);
+
+ if (gt.isTriangle())
+ buildTriangle(nPoints, order, points);
+ else if (gt.isQuadrilateral())
+ buildQuad(nPoints, order, points);
+ else {
+ DUNE_THROW(Dune::NotImplemented,
+ "Lagrange points not yet implemented for this GeometryType.");
+ }
+
+ assert(points.size() == nPoints);
+}
+
+
+// Construct the point set in a triangle element.
+// Loop from the outside to the inside
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,2>::buildTriangle (std::size_t nPoints, int order, Points& points) const
+{
+ points.reserve(nPoints);
+
+ const int nVertexDOFs = 3;
+ const int nEdgeDOFs = 3 * (order-1);
+
+ static const unsigned int vertexPerm[3] = {0,1,2};
+ static const unsigned int edgePerm[3] = {0,2,1};
+
+ auto calcKey = [&](int p) -> Key
+ {
+ if (p < nVertexDOFs) {
+ return Key{vertexPerm[p], dim, 0};
+ }
+ else if (p < nVertexDOFs+nEdgeDOFs) {
+ unsigned int entityIndex = (p - nVertexDOFs) / (order-1);
+ unsigned int index = (p - nVertexDOFs) % (order-1);
+ return Key{edgePerm[entityIndex], dim-1, index};
+ }
+ else {
+ unsigned int index = p - (nVertexDOFs + nEdgeDOFs);
+ return Key{0, dim-2, index};
+ }
+ };
+
+ std::array<int,3> bindex;
+
+ double order_d = double(order);
+ for (std::size_t p = 0; p < nPoints; ++p) {
+ barycentricIndex(p, bindex, order);
+ Vec point{bindex[0] / order_d, bindex[1] / order_d};
+ points.push_back(LP{point, calcKey(p)});
+ }
+}
+
+
+// "Barycentric index" is a triplet of integers, each running from 0 to
+// <Order>. It is the index of a point on the triangle in barycentric
+// coordinates.
+template <class K>
+void LagrangePointSetBuilder<K,2>::barycentricIndex (int index, std::array<int,3>& bindex, int order)
+{
+ int max = order;
+ int min = 0;
+
+ // scope into the correct triangle
+ while (index != 0 && index >= 3 * order)
+ {
+ index -= 3 * order;
+ max -= 2;
+ min++;
+ order -= 3;
+ }
+
+ // vertex DOFs
+ if (index < 3)
+ {
+ bindex[index] = bindex[(index + 1) % 3] = min;
+ bindex[(index + 2) % 3] = max;
+ }
+ // edge DOFs
+ else
+ {
+ index -= 3;
+ int dim = index / (order - 1);
+ int offset = (index - dim * (order - 1));
+ bindex[(dim + 1) % 3] = min;
+ bindex[(dim + 2) % 3] = (max - 1) - offset;
+ bindex[dim] = (min + 1) + offset;
+ }
+}
+
+
+// Construct the point set in the quad element
+// 1. build equispaced points with index tuple (i,j)
+// 2. map index tuple to DOF index and LocalKey
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,2>::buildQuad(std::size_t nPoints, int order, Points& points) const
+{
+ points.resize(nPoints);
+
+ std::array<int,2> orders{order, order};
+ std::array<Vec,4> nodes{Vec{0., 0.}, Vec{1., 0.}, Vec{1., 1.}, Vec{0., 1.}};
+
+ for (int n = 0; n <= orders[1]; ++n) {
+ for (int m = 0; m <= orders[0]; ++m) {
+ // int idx = pointIndexFromIJ(m,n,orders);
+
+ const double r = double(m) / orders[0];
+ const double s = double(n) / orders[1];
+ Vec p = (1.0 - r) * (nodes[3] * s + nodes[0] * (1.0 - s)) +
+ r * (nodes[2] * s + nodes[1] * (1.0 - s));
+
+ auto [idx,key] = calcQuadKey(m,n,orders);
+ points[idx] = LP{p, key};
+ // points[idx] = LP{p, calcQuadKey(n,m,orders)};
+ }
+ }
+}
+
+
+// Obtain the VTK DOF index of the node (i,j) in the quad element
+// and construct a LocalKey
+template <class K>
+std::pair<int,typename LagrangePointSetBuilder<K,2>::Key>
+LagrangePointSetBuilder<K,2>::calcQuadKey (int i, int j, std::array<int,2> order)
+{
+ const bool ibdy = (i == 0 || i == order[0]);
+ const bool jbdy = (j == 0 || j == order[1]);
+ const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0); // How many boundaries do we lie on at once?
+
+ int dof = 0;
+ unsigned int entityIndex = 0;
+ unsigned int index = 0;
+
+ if (nbdy == 2) // Vertex DOF
+ {
+ dof = (i ? (j ? 2 : 1) : (j ? 3 : 0));
+ entityIndex = (j ? (i ? 3 : 2) : (i ? 1 : 0));
+ return std::make_pair(dof,Key{entityIndex, dim, 0});
+ }
+
+ int offset = 4;
+ if (nbdy == 1) // Edge DOF
+ {
+ if (!ibdy) {
+ dof = (i - 1) + (j ? order[0]-1 + order[1]-1 : 0) + offset;
+ entityIndex = j ? 3 : 2;
+ index = i-1;
+ }
+ else if (!jbdy) {
+ dof = (j - 1) + (i ? order[0]-1 : 2 * (order[0]-1) + order[1]-1) + offset;
+ entityIndex = i ? 1 : 0;
+ index = j-1;
+ }
+ return std::make_pair(dof, Key{entityIndex, dim-1, index});
+ }
+
+ offset += 2 * (order[0]-1 + order[1]-1);
+
+ // nbdy == 0: Face DOF
+ dof = offset + (i - 1) + (order[0]-1) * ((j - 1));
+ Key innerKey = LagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
+ return std::make_pair(dof, Key{0, dim-2, innerKey.index()});
+}
+
+
+// Lagrange points on 3d geometries
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,3>::operator() (GeometryType gt, unsigned int order, Points& points) const
+{
+ std::size_t nPoints = numLagrangePoints(gt.id(), dim, order);
+
+ if (gt.isTetrahedron())
+ buildTetra(nPoints, order, points);
+ else if (gt.isHexahedron())
+ buildHex(nPoints, order, points);
+ else {
+ DUNE_THROW(Dune::NotImplemented,
+ "Lagrange points not yet implemented for this GeometryType.");
+ }
+
+ assert(points.size() == nPoints);
+}
+
+
+// Construct the point set in the tetrahedron element
+// 1. construct barycentric (index) coordinates
+// 2. obtains the DOF index, LocalKey and actual coordinate from barycentric index
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,3>::buildTetra (std::size_t nPoints, int order, Points& points) const
+{
+ points.reserve(nPoints);
+
+ const int nVertexDOFs = 4;
+ const int nEdgeDOFs = 6 * (order-1);
+ const int nFaceDOFs = 4 * (order-1)*(order-2)/2;
+
+ static const unsigned int vertexPerm[4] = {0,1,2,3};
+ static const unsigned int edgePerm[6] = {0,2,1,3,4,5};
+ static const unsigned int facePerm[4] = {1,2,0,3};
+
+ auto calcKey = [&](int p) -> Key
+ {
+ if (p < nVertexDOFs) {
+ return Key{vertexPerm[p], dim, 0};
+ }
+ else if (p < nVertexDOFs+nEdgeDOFs) {
+ unsigned int entityIndex = (p - nVertexDOFs) / (order-1);
+ unsigned int index = (p - nVertexDOFs) % (order-1);
+ return Key{edgePerm[entityIndex], dim-1, index};
+ }
+ else if (p < nVertexDOFs+nEdgeDOFs+nFaceDOFs) {
+ unsigned int index = (p - (nVertexDOFs + nEdgeDOFs)) % ((order-1)*(order-2)/2);
+ unsigned int entityIndex = (p - (nVertexDOFs + nEdgeDOFs)) / ((order-1)*(order-2)/2);
+ return Key{facePerm[entityIndex], dim-2, index};
+ }
+ else {
+ unsigned int index = p - (nVertexDOFs + nEdgeDOFs + nFaceDOFs);
+ return Key{0, dim-3, index};
+ }
+ };
+
+ std::array<int,4> bindex;
+
+ double order_d = double(order);
+ for (std::size_t p = 0; p < nPoints; ++p) {
+ barycentricIndex(p, bindex, order);
+ Vec point{bindex[0] / order_d, bindex[1] / order_d, bindex[2] / order_d};
+ points.push_back(LP{point, calcKey(p)});
+ }
+}
+
+
+// "Barycentric index" is a set of 4 integers, each running from 0 to
+// <Order>. It is the index of a point in the tetrahedron in barycentric
+// coordinates.
+template <class K>
+void LagrangePointSetBuilder<K,3>::barycentricIndex (std::size_t p, std::array<int,4>& bindex, int order)
+{
+ const int nVertexDOFs = 4;
+ const int nEdgeDOFs = 6 * (order-1);
+
+ static const int edgeVertices[6][2] = {{0,1}, {1,2}, {2,0}, {0,3}, {1,3}, {2,3}};
+ static const int linearVertices[4][4] = {{0,0,0,1}, {1,0,0,0}, {0,1,0,0}, {0,0,1,0}};
+ static const int vertexMaxCoords[4] = {3,0,1,2};
+ static const int faceBCoords[4][3] = {{0,2,3}, {2,0,1}, {2,1,3}, {1,0,3}};
+ static const int faceMinCoord[4] = {1,3,0,2};
+
+ int max = order;
+ int min = 0;
+
+ // scope into the correct tetra
+ while (p >= 2 * (order * order + 1) && p != 0 && order > 3)
+ {
+ p -= 2 * (order * order + 1);
+ max -= 3;
+ min++;
+ order -= 4;
+ }
+
+ // vertex DOFs
+ if (p < nVertexDOFs)
+ {
+ for (int coord = 0; coord < 4; ++coord)
+ bindex[coord] = (coord == vertexMaxCoords[p] ? max : min);
+ }
+ // edge DOFs
+ else if (p < nVertexDOFs+nEdgeDOFs)
+ {
+ int edgeId = (p - nVertexDOFs) / (order-1);
+ int vertexId = (p - nVertexDOFs) % (order-1);
+ for (int coord = 0; coord < 4; ++coord)
+ {
+ bindex[coord] = min +
+ (linearVertices[edgeVertices[edgeId][0]][coord] * (max - min - 1 - vertexId) +
+ linearVertices[edgeVertices[edgeId][1]][coord] * (1 + vertexId));
+ }
+ }
+ // face DOFs
+ else
+ {
+ int faceId = (p - (nVertexDOFs+nEdgeDOFs)) / ((order-2)*(order-1)/2);
+ int vertexId = (p - (nVertexDOFs+nEdgeDOFs)) % ((order-2)*(order-1)/2);
+
+ std::array<int,3> projectedBIndex;
+ if (order == 3)
+ projectedBIndex[0] = projectedBIndex[1] = projectedBIndex[2] = 0;
+ else
+ LagrangePointSetBuilder<K,2>::barycentricIndex(vertexId, projectedBIndex, order-3);
+
+ for (int i = 0; i < 3; i++)
+ bindex[faceBCoords[faceId][i]] = (min + 1 + projectedBIndex[i]);
+
+ bindex[faceMinCoord[faceId]] = min;
+ }
+}
+
+
+// Construct the point set in the hexahedral element
+// 1. build equispaced points with index tuple (i,j,k)
+// 2. map index tuple to DOF index and LocalKey
+template <class K>
+ template <class Points>
+void LagrangePointSetBuilder<K,3>::buildHex (std::size_t nPoints, int order, Points& points) const
+{
+ points.resize(nPoints);
+
+ std::array<int,3> orders{order, order, order};
+ std::array<Vec,8> nodes{Vec{0., 0., 0.}, Vec{1., 0., 0.}, Vec{1., 1., 0.}, Vec{0., 1., 0.},
+ Vec{0., 0., 1.}, Vec{1., 0., 1.}, Vec{1., 1., 1.}, Vec{0., 1., 1.}};
+
+ for (int p = 0; p <= orders[2]; ++p) {
+ for (int n = 0; n <= orders[1]; ++n) {
+ for (int m = 0; m <= orders[0]; ++m) {
+ const double r = double(m) / orders[0];
+ const double s = double(n) / orders[1];
+ const double t = double(p) / orders[2];
+ Vec point = (1.0-r) * ((nodes[3] * (1.0-t) + nodes[7] * t) * s + (nodes[0] * (1.0-t) + nodes[4] * t) * (1.0-s)) +
+ r * ((nodes[2] * (1.0-t) + nodes[6] * t) * s + (nodes[1] * (1.0-t) + nodes[5] * t) * (1.0-s));
+
+ auto [idx,key] = calcHexKey(m,n,p,orders);
+ points[idx] = LP{point, key};
+ }
+ }
+ }
+}
+
+
+// Obtain the VTK DOF index of the node (i,j,k) in the hexahedral element
+template <class K>
+std::pair<int,typename LagrangePointSetBuilder<K,3>::Key>
+LagrangePointSetBuilder<K,3>::calcHexKey (int i, int j, int k, std::array<int,3> order)
+{
+ const bool ibdy = (i == 0 || i == order[0]);
+ const bool jbdy = (j == 0 || j == order[1]);
+ const bool kbdy = (k == 0 || k == order[2]);
+ const int nbdy = (ibdy ? 1 : 0) + (jbdy ? 1 : 0) + (kbdy ? 1 : 0); // How many boundaries do we lie on at once?
+
+ int dof = 0;
+ unsigned int entityIndex = 0;
+ unsigned int index = 0;
+
+ if (nbdy == 3) // Vertex DOF
+ {
+ dof = (i ? (j ? 2 : 1) : (j ? 3 : 0)) + (k ? 4 : 0);
+ entityIndex = (i ? 1 : 0) + (j ? 2 : 0) + (k ? 4 : 0);
+ return std::make_pair(dof, Key{entityIndex, dim, 0});
+ }
+
+ int offset = 8;
+ if (nbdy == 2) // Edge DOF
+ {
+ entityIndex = (k ? 8 : 4);
+ if (!ibdy)
+ { // On i axis
+ dof = (i - 1) + (j ? order[0]-1 + order[1]-1 : 0) + (k ? 2 * (order[0]-1 + order[1]-1) : 0) + offset;
+ index = i;
+ entityIndex += (i ? 1 : 0);
+ }
+ else if (!jbdy)
+ { // On j axis
+ dof = (j - 1) + (i ? order[0]-1 : 2 * (order[0]-1) + order[1]-1) + (k ? 2 * (order[0]-1 + order[1]-1) : 0) + offset;
+ index = j;
+ entityIndex += (j ? 3 : 2);
+ }
+ else
+ { // !kbdy, On k axis
+ offset += 4 * (order[0]-1) + 4 * (order[1]-1);
+ dof = (k - 1) + (order[2]-1) * (i ? (j ? 3 : 1) : (j ? 2 : 0)) + offset;
+ index = k;
+ entityIndex = (i ? 1 : 0) + (j ? 2 : 0);
+ }
+ return std::make_pair(dof, Key{entityIndex, dim-1, index});
+ }
+
+ offset += 4 * (order[0]-1 + order[1]-1 + order[2]-1);
+ if (nbdy == 1) // Face DOF
+ {
+ Key faceKey;
+ if (ibdy) // On i-normal face
+ {
+ dof = (j - 1) + ((order[1]-1) * (k - 1)) + (i ? (order[1]-1) * (order[2]-1) : 0) + offset;
+ entityIndex = (i ? 1 : 0);
+ faceKey = LagrangePointSetBuilder<K,2>::calcQuadKey(j-1,k-1,{order[1]-2, order[2]-2}).second;
+ }
+ else {
+ offset += 2 * (order[1] - 1) * (order[2] - 1);
+ if (jbdy) // On j-normal face
+ {
+ dof = (i - 1) + ((order[0]-1) * (k - 1)) + (j ? (order[2]-1) * (order[0]-1) : 0) + offset;
+ entityIndex = (j ? 3 : 2);
+ faceKey = LagrangePointSetBuilder<K,2>::calcQuadKey(i-1,k-1,{order[0]-2, order[2]-2}).second;
+ }
+ else
+ { // kbdy, On k-normal face
+ offset += 2 * (order[2]-1) * (order[0]-1);
+ dof = (i - 1) + ((order[0]-1) * (j - 1)) + (k ? (order[0]-1) * (order[1]-1) : 0) + offset;
+ entityIndex = (k ? 5 : 4);
+ faceKey = LagrangePointSetBuilder<K,2>::calcQuadKey(i-1,j-1,{order[0]-2, order[1]-2}).second;
+ }
+ }
+ return std::make_pair(dof, Key{entityIndex, dim-2, faceKey.index()});
+ }
+
+ // nbdy == 0: Body DOF
+ offset += 2 * ((order[1]-1) * (order[2]-1) + (order[2]-1) * (order[0]-1) + (order[0]-1) * (order[1]-1));
+ dof = offset + (i - 1) + (order[0]-1) * ((j - 1) + (order[1]-1) * ((k - 1)));
+ Key innerKey = LagrangePointSetBuilder<K,3>::calcHexKey(i-1,j-1,k-1,{order[0]-2, order[1]-2, order[2]-2}).second;
+ return std::make_pair(dof, Key{0, dim-3, innerKey.index()});
+}
+
+}}} // end namespace Dune::Vtk::Impl
diff --git a/dune/vtk/utility/string.hh b/dune/vtk/utility/string.hh
index 7724911a9df42db7474b297b0f52424913242695..be35b43f5249ea496548cf2d0c289ca94d937dbb 100644
--- a/dune/vtk/utility/string.hh
+++ b/dune/vtk/utility/string.hh
@@ -8,113 +8,116 @@
namespace Dune
{
- /// convert all characters in a string to upper case
- inline std::string to_upper(std::string input)
+ namespace Vtk
{
- for (auto& c : input)
- c = toupper(c);
- return input;
- }
-
- /// convert all characters in a string to upper case
- inline std::string to_lower(std::string input)
- {
- for (auto& c : input)
- c = tolower(c);
- return input;
- }
+ /// convert all characters in a string to upper case
+ inline std::string to_upper(std::string input)
+ {
+ for (auto& c : input)
+ c = toupper(c);
+ return input;
+ }
- /// trim a string from the left
- inline std::string& ltrim(std::string& str)
- {
- auto it = std::find_if(str.begin(), str.end(), [](char ch)
+ /// convert all characters in a string to upper case
+ inline std::string to_lower(std::string input)
{
- return !std::isspace<char>(ch, std::locale::classic());
- });
- str.erase(str.begin() , it);
- return str;
- }
-
- /// trim a string from the right
- inline std::string& rtrim(std::string& str)
- {
- auto it = std::find_if(str.rbegin(), str.rend(), [](char ch)
+ for (auto& c : input)
+ c = tolower(c);
+ return input;
+ }
+
+ /// trim a string from the left
+ inline std::string& ltrim(std::string& str)
{
- return !std::isspace<char>(ch, std::locale::classic());
- });
- str.erase(it.base(), str.end());
- return str;
- }
-
- /// trim a string from both sides
- inline std::string& trim(std::string& str)
- {
- return ltrim(rtrim(str));
- }
+ auto it = std::find_if(str.begin(), str.end(), [](char ch)
+ {
+ return !std::isspace<char>(ch, std::locale::classic());
+ });
+ str.erase(str.begin() , it);
+ return str;
+ }
- /// trim a (copy of the) string from both sides
- inline std::string trim_copy(std::string const& str)
- {
- auto s = str;
- return trim(s);
- }
+ /// trim a string from the right
+ inline std::string& rtrim(std::string& str)
+ {
+ auto it = std::find_if(str.rbegin(), str.rend(), [](char ch)
+ {
+ return !std::isspace<char>(ch, std::locale::classic());
+ });
+ str.erase(it.base(), str.end());
+ return str;
+ }
+ /// trim a string from both sides
+ inline std::string& trim(std::string& str)
+ {
+ return ltrim(rtrim(str));
+ }
- template <class InputIter, class T, class Func>
- void split(InputIter first, InputIter end, T const& t, Func f)
- {
- if (first == end)
- return;
-
- while (true) {
- InputIter found = std::find(first, end, t);
- f(first, found);
- if (found == end)
- break;
- first = ++found;
+ /// trim a (copy of the) string from both sides
+ inline std::string trim_copy(std::string const& str)
+ {
+ auto s = str;
+ return trim(s);
}
- }
- template <class InputIter, class SeparatorIter, class Func>
- void split(InputIter first, InputIter end, SeparatorIter s_first, SeparatorIter s_end, Func f)
- {
- if (first == end)
- return;
-
- while (true) {
- InputIter found = std::find_first_of(first, end, s_first, s_end);
- f(first, found);
- if (found == end)
- break;
- first = ++found;
+
+ template <class InputIter, class T, class Func>
+ void split(InputIter first, InputIter end, T const& t, Func f)
+ {
+ if (first == end)
+ return;
+
+ while (true) {
+ InputIter found = std::find(first, end, t);
+ f(first, found);
+ if (found == end)
+ break;
+ first = ++found;
+ }
}
- }
- /// Replace all occurences of substring `from` with `to` in source `str`.
- inline void replaceAll(std::string& str, std::string const& from, std::string const& to)
- {
- if (from.empty())
- return;
- std::size_t start_pos = 0;
- while ((start_pos = str.find(from, start_pos)) != std::string::npos)
+ template <class InputIter, class SeparatorIter, class Func>
+ void split(InputIter first, InputIter end, SeparatorIter s_first, SeparatorIter s_end, Func f)
{
- str.replace(start_pos, from.length(), to);
- start_pos += to.length();
+ if (first == end)
+ return;
+
+ while (true) {
+ InputIter found = std::find_first_of(first, end, s_first, s_end);
+ f(first, found);
+ if (found == end)
+ break;
+ first = ++found;
+ }
}
- }
+ /// Replace all occurences of substring `from` with `to` in source `str`.
+ inline void replaceAll(std::string& str, std::string const& from, std::string const& to)
+ {
+ if (from.empty())
+ return;
+ std::size_t start_pos = 0;
+ while ((start_pos = str.find(from, start_pos)) != std::string::npos)
+ {
+ str.replace(start_pos, from.length(), to);
+ start_pos += to.length();
+ }
+ }
- template <class InputIter>
- std::string join (InputIter first, InputIter end, std::string sep = " ")
- {
- if (first == end)
- return "";
- std::ostringstream os;
- os << *first++;
- while (first != end)
- os << sep << *first++;
- return os.str();
- }
+ template <class InputIter>
+ std::string join (InputIter first, InputIter end, std::string sep = " ")
+ {
+ if (first == end)
+ return "";
+
+ std::ostringstream os;
+ os << *first++;
+ while (first != end)
+ os << sep << *first++;
+ return os.str();
+ }
+ } // end namespace Vtk
} // end namspace Dune
diff --git a/dune/vtk/utility/test/test-lagrange.cc b/dune/vtk/utility/test/test-lagrange.cc
index 62ed50a11e6a44cb4bc73ace86ccffc31a2b2fdc..f25c2cb55826ebfc2c5435e6d14c6904bbf69dd3 100644
--- a/dune/vtk/utility/test/test-lagrange.cc
+++ b/dune/vtk/utility/test/test-lagrange.cc
@@ -75,8 +75,8 @@ bool test(const Basis &basis, const Points &points, bool verbose)
template <class Topology>
bool test(unsigned int order, bool verbose = false)
{
- typedef Dune::LagrangeBasisFactory<Dune::VtkLagrangePointSet,Topology::dimension,StorageField,ComputeField> BasisFactory;
- typedef Dune::LagrangeCoefficientsFactory< Dune::VtkLagrangePointSet, Topology::dimension,double > LagrangeCoefficientsFactory;
+ typedef Dune::LagrangeBasisFactory<Dune::Vtk::LagrangePointSet,Topology::dimension,StorageField,ComputeField> BasisFactory;
+ typedef Dune::LagrangeCoefficientsFactory< Dune::Vtk::LagrangePointSet, Topology::dimension,double > LagrangeCoefficientsFactory;
bool ret = true;
@@ -154,14 +154,14 @@ int main ( int argc, char **argv )
#endif
#ifdef CHECKDIM2
- tests &= test<Prism<Prism<Point> > > (order); // quad
- tests &= test<Pyramid<Pyramid<Point> > >(order); // triangle
+ //tests &= test<Prism<Prism<Point> > > (order); // quad
+ //tests &= test<Pyramid<Pyramid<Point> > >(order); // triangle
#endif
#ifdef CHECKDIM3
- tests &= test<Prism<Prism<Prism<Point> > > >(order); // hexahedron
+ // tests &= test<Prism<Prism<Prism<Point> > > >(order); // hexahedron
// tests &= test<Prism<Pyramid<Pyramid<Point> > > >(order);
- tests &= test<Pyramid<Pyramid<Pyramid<Point> > > >(order); // tetrahedron
+ // tests &= test<Pyramid<Pyramid<Pyramid<Point> > > >(order); // tetrahedron
#endif
return (tests ? 0 : 1);
diff --git a/dune/vtk/utility/uid.hh b/dune/vtk/utility/uid.hh
index 22c2bb7e7bf4eb72c4aa2062cea325b39c428825..8fbe700e92f2d725610ab48965e91fd554775ba8 100644
--- a/dune/vtk/utility/uid.hh
+++ b/dune/vtk/utility/uid.hh
@@ -7,16 +7,19 @@
namespace Dune
{
- inline std::string uid (std::size_t len = 8)
+ namespace Vtk
{
- static const auto digits = "0123456789abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ";
- static const int N = std::strlen(digits);
+ inline std::string uid (std::size_t len = 8)
+ {
+ static const auto digits = "0123456789abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ";
+ static const int N = std::strlen(digits);
- std::string id(len,' ');
- for (std::size_t i = 0; i < len; ++i)
- id[i] = digits[std::rand()%N];
+ std::string id(len,' ');
+ for (std::size_t i = 0; i < len; ++i)
+ id[i] = digits[std::rand()%N];
- return id;
- }
+ return id;
+ }
+ } // end namespace Vtk
} // end namespace Dune
diff --git a/dune/vtk/vtkfunction.hh b/dune/vtk/vtkfunction.hh
deleted file mode 100644
index d51f64c9af1a69c6affc6ff4481031a289fdeba9..0000000000000000000000000000000000000000
--- a/dune/vtk/vtkfunction.hh
+++ /dev/null
@@ -1,128 +0,0 @@
-#pragma once
-
-#include <optional>
-#include <type_traits>
-
-#include <dune/common/std/type_traits.hh>
-
-#include "vtklocalfunction.hh"
-#include "vtktypes.hh"
-
-namespace Dune
-{
- template <class T, int N>
- class FieldVector;
-
- template <class T, int N, int M>
- class FieldMatrix;
-
-
- /// Wrapper class for functions allowing local evaluations.
- template <class GridView>
- class VtkFunction
- {
- template <class F>
- using LocalFunction = decltype(localFunction(std::declval<F>()));
-
- using Domain = typename GridView::template Codim<0>::Entity::Geometry::LocalCoordinate;
-
- template <class F>
- using Range = std::decay_t<std::result_of_t<F(Domain)>>;
-
- private:
-
- template <class T, int N>
- static auto sizeOfImpl (FieldVector<T,N> const&)
- -> std::integral_constant<int, N> { return {}; }
-
- template <class T, int N, int M>
- static auto sizeOfImpl (FieldMatrix<T,N,M> const&)
- -> std::integral_constant<int, N*M> { return {}; }
-
- static auto sizeOfImpl (...)
- -> std::integral_constant<int, 1> { return {}; }
-
- template <class T>
- static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }
-
- public:
- /// Constructor VtkFunction from legacy VTKFunction
- /**
- * \param fct The VTKFunction to wrap
- * \param type The VTK datatype how to write the function values to the output [Vtk::FLOAT64]
- **/
- VtkFunction (std::shared_ptr<VTKFunction<GridView> const> const& fct,
- std::optional<Vtk::DataTypes> type = {})
- : localFct_(fct)
- , name_(fct->name())
- , ncomps_(fct->ncomps())
- , type_(type ? *type : Vtk::FLOAT64)
- {}
-
- /// Construct VtkFunction from dune-functions GridFunction with Signature
- // NOTE: Stores the localFunction(fct) by value.
- /**
- * \param fct A Grid(View)-function, providing a `localFunction(fct)`
- * \param name The name to use component identification in the VTK file
- * \param ncomps Number of components of the pointwise data. Is extracted
- * from the range type of the GridFunction if not given.
- * \param type The \ref Vtk::DataTypes used in the output. E.g. FLOAT32,
- * or FLOAT64. Is extracted from the range type of the
- * GridFunction if not given.
- **/
- template <class F,
- class = void_t<LocalFunction<F>> >
- VtkFunction (F&& fct, std::string name,
- std::optional<int> ncomps = {},
- std::optional<Vtk::DataTypes> type = {})
- : localFct_(localFunction(std::forward<F>(fct)))
- , name_(std::move(name))
- {
- using R = Range<LocalFunction<F>>;
-
- ncomps_ = ncomps ? *ncomps : sizeOf<R>();
- type_ = type ? *type : Vtk::Map::type<R>();
- }
-
- /// Constructor that forward the number of components and data type to the other constructor
- template <class F,
- class = void_t<LocalFunction<F>> >
- VtkFunction (F&& fct, Vtk::FieldInfo fieldInfo,
- std::optional<Vtk::DataTypes> type = {})
- : VtkFunction(std::forward<F>(fct), fieldInfo.name(), fieldInfo.ncomps(), type)
- {}
-
- VtkFunction () = default;
-
- /// Create a LocalFunction
- friend VtkLocalFunction<GridView> localFunction (VtkFunction const& self)
- {
- return self.localFct_;
- }
-
- /// Return a name associated with the function
- std::string const& name () const
- {
- return name_;
- }
-
- /// Return the number of components of the Range
- int ncomps () const
- {
- return ncomps_ > 3 ? 9 : ncomps_ > 1 ? 3 : 1; // tensor, vector, scalar
- }
-
- /// Return the VTK Datatype associated with the functions range type
- Vtk::DataTypes type () const
- {
- return type_;
- }
-
- private:
- VtkLocalFunction<GridView> localFct_;
- std::string name_;
- int ncomps_ = 1;
- Vtk::DataTypes type_ = Vtk::FLOAT32;
- };
-
-} // end namespace Dune
diff --git a/dune/vtk/vtklocalfunction.hh b/dune/vtk/vtklocalfunction.hh
deleted file mode 100644
index e636b0b4396adb03db94ab4366227ed35cce7c3d..0000000000000000000000000000000000000000
--- a/dune/vtk/vtklocalfunction.hh
+++ /dev/null
@@ -1,71 +0,0 @@
-#pragma once
-
-#include <memory>
-#include <type_traits>
-
-#include <dune/common/std/type_traits.hh>
-
-#include "vtklocalfunctioninterface.hh"
-#include "legacyvtkfunction.hh"
-#include "defaultvtkfunction.hh"
-
-namespace Dune
-{
- /// \brief A VtkLocalFunction is a function-like object that can be bound to a grid element
- /// an that provides an evaluate method with a component argument.
- /**
- * Stores internally a VtkLocalFunctionInterface object for the concrete evaluation.
- **/
- template <class GridView>
- class VtkLocalFunction
- {
- using Self = VtkLocalFunction;
- using Entity = typename GridView::template Codim<0>::Entity;
- using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
-
- template <class LF, class E>
- using HasBind = decltype(std::declval<LF>().bind(std::declval<E>()));
-
- public:
- /// Construct the VtkLocalFunction from any function object that has a bind(element) method.
- template <class LF,
- disableCopyMove<Self, LF> = 0,
- class = void_t<HasBind<LF,Entity>> >
- VtkLocalFunction (LF&& lf)
- : localFct_(std::make_shared<LocalFunctionWrapper<GridView,LF>>(std::forward<LF>(lf)))
- {}
-
- /// Construct a VtkLocalFunction from a legacy VTKFunction
- VtkLocalFunction (std::shared_ptr<VTKFunction<GridView> const> const& lf)
- : localFct_(std::make_shared<VTKLocalFunctionWrapper<GridView>>(lf))
- {}
-
- /// Allow the default construction of a VtkLocalFunction
- VtkLocalFunction () = default;
-
- /// Bind the function to the grid entity
- void bind (Entity const& entity)
- {
- assert(bool(localFct_));
- localFct_->bind(entity);
- }
-
- /// Unbind from the currently bound entity
- void unbind ()
- {
- assert(bool(localFct_));
- localFct_->unbind();
- }
-
- /// Evaluate the `comp` component of the Range value at local coordinate `xi`
- double evaluate (int comp, LocalCoordinate const& xi) const
- {
- assert(bool(localFct_));
- return localFct_->evaluate(comp, xi);
- }
-
- private:
- std::shared_ptr<VtkLocalFunctionInterface<GridView>> localFct_ = nullptr;
- };
-
-} // end namespace Dune
diff --git a/dune/vtk/vtklocalfunctioninterface.hh b/dune/vtk/vtklocalfunctioninterface.hh
deleted file mode 100644
index 9ba2afeccea47eb143edfab28d2d84e35cd0bc1e..0000000000000000000000000000000000000000
--- a/dune/vtk/vtklocalfunctioninterface.hh
+++ /dev/null
@@ -1,27 +0,0 @@
-#pragma once
-
-namespace Dune
-{
- /// \brief An abstract base class for LocalFunctions that can be bound to an element and
- /// evaluated in local coordinates w.r.t. to a component of its value.
- template <class GridView>
- class VtkLocalFunctionInterface
- {
- public:
- using Entity = typename GridView::template Codim<0>::Entity;
- using LocalCoordinate = typename Entity::Geometry::LocalCoordinate;
-
- /// Bind the function to the grid entity
- virtual void bind (Entity const& entity) = 0;
-
- /// Unbind from the currently bound entity
- virtual void unbind () = 0;
-
- /// Evaluate single component comp in the entity at local coordinates xi
- virtual double evaluate (int comp, LocalCoordinate const& xi) const = 0;
-
- /// Virtual destructor
- virtual ~VtkLocalFunctionInterface () = default;
- };
-
-} // end namespace Dune
diff --git a/dune/vtk/vtkreader.hh b/dune/vtk/vtkreader.hh
index e3b9c74fdc1d2208d4140f8405782ce80494ba0a..ea4ded7e57bb04bbe5304a1fdd4e95c09f4fe30a 100644
--- a/dune/vtk/vtkreader.hh
+++ b/dune/vtk/vtkreader.hh
@@ -10,7 +10,7 @@
#include <dune/vtk/filereader.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
// default GridCreator
#include <dune/vtk/gridcreators/continuousgridcreator.hh>
@@ -26,7 +26,7 @@ namespace Dune
**/
template <class Grid, class GridCreator>
class VtkReader
- : public FileReader<Grid, VtkReader<Grid, GridCreator>>
+ : public Vtk::FileReader<Grid, VtkReader<Grid, GridCreator>>
{
// Sections visited during the xml parsing
enum Sections {
@@ -61,18 +61,18 @@ namespace Dune
: creator_(std::move(creator))
{}
- /// Read the grid from file with `filename` into the GridFactory \ref factory_
+ /// Read the grid from file with `filename` into the GridCreator
/**
* \param filename The name of the input file
- * \param create If `false`, only fill internal data structures, if `true`, also create the grid. [true]
+ * \param fillCreator If `false`, only fill internal data structures, if `true`, pass the internal data to the GridCreator. [true]
**/
- void readFromFile (std::string const& filename, bool create = true);
+ void read (std::string const& filename, bool fillCreator = true);
/// Implementation of \ref FileReader interface
- static void readFactoryImpl (GridFactory<Grid>& factory, std::string const& filename)
+ static void fillFactoryImpl (GridFactory<Grid>& factory, std::string const& filename)
{
VtkReader reader{factory};
- reader.readFromFile(filename);
+ reader.read(filename);
}
/// Obtains the creator of the grid
@@ -81,6 +81,12 @@ namespace Dune
return *creator_;
}
+ /// Construct the actual grid using the GridCreator
+ std::unique_ptr<Grid> createGrid () const
+ {
+ return creator_->createGrid();
+ }
+
/// Advanced read methods
/// @{
@@ -98,9 +104,9 @@ namespace Dune
**/
void readParallelFileFromStream (std::ifstream& input, int rank, int size, bool create = true);
- /// Construct a grid using the GridCreator
- // NOTE: requires the internal data structures to be filled by an aforgoing call to readFromFile
- void createGrid (bool insertPieces = true);
+ /// Insert all internal data to the GridCreator
+ // NOTE: requires the internal data structures to be filled by an aforgoing call to read()
+ void fillGridCreator (bool insertPieces = true);
/// @}
@@ -197,7 +203,7 @@ namespace Dune
// deduction guides
template <class Grid>
VtkReader (GridFactory<Grid>&)
- -> VtkReader<Grid, ContinuousGridCreator<Grid>>;
+ -> VtkReader<Grid, Vtk::ContinuousGridCreator<Grid>>;
template <class GridCreator,
class = std::void_t<typename GridCreator::Grid>>
diff --git a/dune/vtk/vtkreader.impl.hh b/dune/vtk/vtkreader.impl.hh
index 604dd0caf6d330059d1062d256e8372013acc54f..c3a0bf01d1ff337891f3bc679c43d98ea203927c 100644
--- a/dune/vtk/vtkreader.impl.hh
+++ b/dune/vtk/vtkreader.impl.hh
@@ -16,21 +16,21 @@
namespace Dune {
template <class Grid, class Creator>
-void VtkReader<Grid,Creator>::readFromFile (std::string const& filename, bool create)
+void VtkReader<Grid,Creator>::read (std::string const& filename, bool fillCreator)
{
// check whether file exists!
- if (!filesystem::exists(filename))
+ if (!Vtk::exists(filename))
DUNE_THROW(IOError, "File " << filename << " does not exist!");
std::ifstream input(filename, std::ios_base::in | std::ios_base::binary);
assert(input.is_open());
- std::string ext = filesystem::path(filename).extension().string();
+ std::string ext = Vtk::Path(filename).extension().string();
if (ext == ".vtu") {
- readSerialFileFromStream(input, create);
+ readSerialFileFromStream(input, fillCreator);
pieces_.push_back(filename);
} else if (ext == ".pvtu") {
- readParallelFileFromStream(input, comm().rank(), comm().size(), create);
+ readParallelFileFromStream(input, comm().rank(), comm().size(), fillCreator);
} else {
DUNE_THROW(IOError, "File has unknown file-extension '" << ext << "'. Allowed are only '.vtu' and '.pvtu'.");
}
@@ -38,7 +38,7 @@ void VtkReader<Grid,Creator>::readFromFile (std::string const& filename, bool cr
template <class Grid, class Creator>
-void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bool create)
+void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bool fillCreator)
{
clear();
std::string compressor = "";
@@ -49,7 +49,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
Sections section = NO_SECTION;
for (std::string line; std::getline(input, line); ) {
- ltrim(line);
+ Vtk::ltrim(line);
if (isSection(line, "VTKFile", section)) {
bool closed = false;
@@ -109,7 +109,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
data_type = Vtk::Map::to_datatype[attr["type"]];
if (!attr["Name"].empty())
- data_name = to_lower(attr["Name"]);
+ data_name = Vtk::to_lower(attr["Name"]);
else if (section == POINTS)
data_name = "points";
@@ -118,7 +118,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
data_components = std::stoul(attr["NumberOfComponents"]);
// determine FormatType
- data_format = to_lower(attr["format"]);
+ data_format = Vtk::to_lower(attr["format"]);
if (data_format == "appended") {
format_ = !compressor.empty() ? Vtk::COMPRESSED : Vtk::BINARY;
} else {
@@ -139,7 +139,7 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
if (data_format == "appended") {
if (!closed) {
while (std::getline(input, line)) {
- ltrim(line);
+ Vtk::ltrim(line);
if (line.substr(1,10) == "/DataArray")
break;
}
@@ -225,19 +225,19 @@ void VtkReader<Grid,Creator>::readSerialFileFromStream (std::ifstream& input, bo
if (section != NO_SECTION)
DUNE_THROW(IOError, "VTK-File is incomplete. It must end with </VTKFile>!");
- if (create)
- createGrid();
+ if (fillCreator)
+ fillGridCreator();
}
template <class Grid, class Creator>
-void VtkReader<Grid,Creator>::readParallelFileFromStream (std::ifstream& input, int commRank, int commSize, bool create)
+void VtkReader<Grid,Creator>::readParallelFileFromStream (std::ifstream& input, int commRank, int commSize, bool fillCreator)
{
clear();
Sections section = NO_SECTION;
for (std::string line; std::getline(input, line); ) {
- ltrim(line);
+ Vtk::ltrim(line);
if (isSection(line, "VTKFile", section)) {
bool closed = false;
@@ -276,8 +276,8 @@ void VtkReader<Grid,Creator>::readParallelFileFromStream (std::ifstream& input,
if (section != NO_SECTION)
DUNE_THROW(IOError, "VTK-File is incomplete. It must end with </VTKFile>!");
- if (create)
- createGrid();
+ if (fillCreator)
+ fillGridCreator();
}
@@ -297,7 +297,7 @@ Sections readDataArray (IStream& input, std::vector<T>& values, std::size_t max_
std::size_t idx = 0;
for (std::string line; std::getline(input, line);) {
- trim(line);
+ Vtk::trim(line);
if (line.substr(1,10) == "/DataArray")
return parent_section;
if (line[0] == '<')
@@ -318,7 +318,7 @@ template <class IStream, class Sections>
Sections skipRestOfDataArray (IStream& input, Sections section, Sections parent_section)
{
for (std::string line; std::getline(input, line);) {
- ltrim(line);
+ Vtk::ltrim(line);
if (line.substr(1,10) == "/DataArray")
return parent_section;
}
@@ -543,7 +543,7 @@ void VtkReader<Grid,Creator>::readAppended (std::ifstream& input, std::vector<T>
template <class Grid, class Creator>
-void VtkReader<Grid,Creator>::createGrid (bool insertPieces)
+void VtkReader<Grid,Creator>::fillGridCreator (bool insertPieces)
{
assert(vec_points.size() == numberOfPoints_);
assert(vec_types.size() == numberOfCells_);
diff --git a/dune/vtk/vtktimeserieswriter.hh b/dune/vtk/vtktimeserieswriter.hh
index bbac2b301ff0b4f167c0b7814b863ddee83f3a9e..039244508d0b52fbf944944083c7815f8aa4967b 100644
--- a/dune/vtk/vtktimeserieswriter.hh
+++ b/dune/vtk/vtktimeserieswriter.hh
@@ -8,7 +8,7 @@
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/utility/filesystem.hh>
#include <dune/vtk/utility/uid.hh>
@@ -22,7 +22,7 @@ namespace Dune
**/
template <class VtkWriter>
class VtkTimeseriesWriter
- : public FileWriter
+ : public Vtk::FileWriter
{
protected:
using Self = VtkTimeseriesWriter;
@@ -47,7 +47,7 @@ namespace Dune
std::srand(std::time(nullptr));
// put temporary file to /tmp directory
- tmpDir_ = filesystem::path("/tmp/vtktimeserieswriter_" + uid(10) + "/");
+ tmpDir_ = Vtk::Path("/tmp/vtktimeserieswriter_" + Vtk::uid(10) + "/");
}
/// Remove all written intermediate files and remove temporary directory
@@ -75,7 +75,7 @@ namespace Dune
*
* \param fn Filename of the Timeseries file. May contain a directory and any file extension.
* \param dir The optional parameter specifies the directory of the partition files.
- *
+ *
* \returns File name of the actual written file
**/
virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
@@ -88,7 +88,7 @@ namespace Dune
return *this;
}
- /// Attach cell data to the writer, \see VtkFunction for possible arguments
+ /// Attach cell data to the writer, \see Vtk::Function for possible arguments
template <class Function, class... Args>
VtkTimeseriesWriter& addCellData (Function const& fct, Args&&... args)
{
@@ -98,7 +98,7 @@ namespace Dune
protected:
VtkWriter vtkWriter_;
- filesystem::path tmpDir_;
+ Vtk::Path tmpDir_;
mutable bool initialized_ = false;
diff --git a/dune/vtk/vtktimeserieswriter.impl.hh b/dune/vtk/vtktimeserieswriter.impl.hh
index f3822daa7fc832cbe157a3e91fe71b06a55146fe..3f2c7737525b7162545cbaefeb16c8301cc47614 100644
--- a/dune/vtk/vtktimeserieswriter.impl.hh
+++ b/dune/vtk/vtktimeserieswriter.impl.hh
@@ -41,8 +41,8 @@ template <class W>
void VtkTimeseriesWriter<W>
::writeTimestep (double time, std::string const& fn, std::optional<std::string> tmpDir, bool writeCollection) const
{
- auto name = filesystem::path(fn).stem();
- auto tmpBase = tmpDir ? filesystem::path(*tmpDir) : tmpDir_;
+ auto name = Vtk::Path(fn).stem();
+ auto tmpBase = tmpDir ? Vtk::Path(*tmpDir) : tmpDir_;
auto tmp = tmpBase;
tmp /= name.string();
@@ -54,7 +54,7 @@ void VtkTimeseriesWriter<W>
filenameBase = tmp.string() + "_p" + std::to_string(vtkWriter_.comm().rank());
if (!initialized_) {
- filesystem::create_directories(tmpBase);
+ Vtk::createDirectories(tmpBase);
// write points and cells only once
filenameMesh_ = filenameBase + ".mesh.vtkdata";
@@ -79,13 +79,13 @@ std::string VtkTimeseriesWriter<W>
{
assert( initialized_ );
- auto p = filesystem::path(fn);
+ auto p = Vtk::Path(fn);
auto name = p.stem();
- p.remove_filename();
+ p.removeFilename();
- filesystem::path fn_dir = p;
- filesystem::path data_dir = dir ? filesystem::path(*dir) : fn_dir;
- filesystem::path rel_dir = filesystem::relative(data_dir, fn_dir);
+ Vtk::Path fn_dir = p;
+ Vtk::Path data_dir = dir ? Vtk::Path(*dir) : fn_dir;
+ Vtk::Path rel_dir = Vtk::relative(data_dir, fn_dir);
std::string serial_fn = fn_dir.string() + '/' + name.string() + "_ts";
std::string parallel_fn = data_dir.string() + '/' + name.string() + "_ts";
diff --git a/dune/vtk/vtkwriter.hh b/dune/vtk/vtkwriter.hh
index b360f86ea7e657a8da9ea1f290842c4d5a694ecb..05b437c05fb860abfbee6242ff64634d24c841dc 100644
--- a/dune/vtk/vtkwriter.hh
+++ b/dune/vtk/vtkwriter.hh
@@ -40,7 +40,7 @@ namespace Dune
: public Impl::VtkWriterImpl<GridView, typename GridView::Grid>::type
{
using Super = typename Impl::VtkWriterImpl<GridView, typename GridView::Grid>::type;
-
+
public:
using Super::Super;
};
@@ -58,7 +58,7 @@ namespace Dune
template <class GridView, int dim, class Coordinates>
struct VtkWriterImpl<GridView, YaspGrid<dim,Coordinates>>
{
- using type = VtkRectilinearGridWriter<GridView, YaspDataCollector<GridView>>;
+ using type = VtkRectilinearGridWriter<GridView, Vtk::YaspDataCollector<GridView>>;
};
#if HAVE_DUNE_SPGRID
@@ -66,7 +66,7 @@ namespace Dune
template <class GridView, class ct, int dim, template <int> class Ref, class Comm>
struct VtkWriterImpl<GridView, SPGrid<ct,dim,Ref,Comm>>
{
- using type = VtkRectilinearGridWriter<GridView, SPDataCollector<GridView>>;
+ using type = VtkRectilinearGridWriter<GridView, Vtk::SPDataCollector<GridView>>;
};
#endif
@@ -74,7 +74,7 @@ namespace Dune
template <class GridView, int dim, class ct>
struct VtkWriterImpl<GridView, YaspGrid<dim,TensorProductCoordinates<ct,dim>>>
{
- using type = VtkRectilinearGridWriter<GridView, YaspDataCollector<GridView>>;
+ using type = VtkRectilinearGridWriter<GridView, Vtk::YaspDataCollector<GridView>>;
};
// A transformed structured grid has structured connectivity but unstructured point
@@ -82,7 +82,7 @@ namespace Dune
template <class GridView, int dim, class Coordinates, class CoordFunction, class Allocator>
struct VtkWriterImpl<GridView, GeometryGrid<YaspGrid<dim,Coordinates>, CoordFunction, Allocator>>
{
- using type = VtkStructuredGridWriter<GridView, YaspDataCollector<GridView>>;
+ using type = VtkStructuredGridWriter<GridView, Vtk::YaspDataCollector<GridView>>;
};
} // end namespace Impl
diff --git a/dune/vtk/vtkwriterinterface.hh b/dune/vtk/vtkwriterinterface.hh
index 8dadb0a5ed609c1c8d9af6bb7b066935fee4c12c..7ef3e1e6894d01e62143eca3ee486e3b413bfb1a 100644
--- a/dune/vtk/vtkwriterinterface.hh
+++ b/dune/vtk/vtkwriterinterface.hh
@@ -10,8 +10,8 @@
#include <dune/common/parallel/mpihelper.hh>
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
namespace Dune
{
@@ -22,15 +22,15 @@ namespace Dune
**/
template <class GridView, class DataCollector>
class VtkWriterInterface
- : public FileWriter
+ : public Vtk::FileWriter
{
- template <class> friend class VtkTimeseriesWriter;
+ template <class> friend class TimeseriesWriter;
template <class> friend class PvdWriter;
protected:
static constexpr int dimension = GridView::dimension;
- using VtkFunction = Dune::VtkFunction<GridView>;
+ using VtkFunction = Dune::Vtk::Function<GridView>;
using Communicator = CollectiveCommunication<typename MPIHelper::MPICommunicator>;
using pos_type = typename std::ostream::pos_type;
@@ -81,7 +81,7 @@ namespace Dune
/**
* \param fn Filename of the VTK file. May contain a directory and any file extension.
* \param dir The optional parameter specifies the directory of the partition files for parallel writes.
- *
+ *
* \returns File name that is actually written.
**/
virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
@@ -90,8 +90,8 @@ namespace Dune
/**
* Attach a global function to the writer that will be evaluated at grid points
* (vertices and higher order points). The global function must be
- * assignable to the function wrapper \ref VtkFunction. Additional argument
- * for output datatype and number of components can be passed. See \ref VtkFunction
+ * assignable to the function wrapper \ref Vtk::Function. Additional argument
+ * for output datatype and number of components can be passed. See \ref Vtk::Function
* Constructor for possible arguments.
**/
template <class Function, class... Args>
@@ -104,9 +104,9 @@ namespace Dune
/// \brief Attach cell data to the writer
/**
* Attach a global function to the writer that will be evaluated at cell centers.
- * The global function must be assignable to the function wrapper \ref VtkFunction.
+ * The global function must be assignable to the function wrapper \ref Vtk::Function.
* Additional argument for output datatype and number of components can be passed.
- * See \ref VtkFunction Constructor for possible arguments.
+ * See \ref Vtk::Function Constructor for possible arguments.
**/
template <class Function, class... Args>
VtkWriterInterface& addCellData (Function&& fct, Args&&... args)
@@ -120,7 +120,7 @@ namespace Dune
void setFormat (Vtk::FormatTypes format)
{
format_ = format;
-
+
#if !HAVE_VTK_ZLIB
if (format_ == Vtk::COMPRESSED) {
std::cout << "Dune is compiled without compression. Falling back to BINARY VTK output!\n";
diff --git a/dune/vtk/vtkwriterinterface.impl.hh b/dune/vtk/vtkwriterinterface.impl.hh
index 597774c9cbacb45aedd3269e8f6d28b7777a0e39..49bd0ea5a61f922e070ccb589fbb1a480eefdc6b 100644
--- a/dune/vtk/vtkwriterinterface.impl.hh
+++ b/dune/vtk/vtkwriterinterface.impl.hh
@@ -28,13 +28,13 @@ std::string VtkWriterInterface<GV,DC>
{
dataCollector_->update();
- auto p = filesystem::path(fn);
+ auto p = Vtk::Path(fn);
auto name = p.stem();
- p.remove_filename();
+ p.removeFilename();
- filesystem::path fn_dir = p;
- filesystem::path data_dir = dir ? filesystem::path(*dir) : fn_dir;
- filesystem::path rel_dir = filesystem::relative(data_dir, fn_dir);
+ Vtk::Path fn_dir = p;
+ Vtk::Path data_dir = dir ? Vtk::Path(*dir) : fn_dir;
+ Vtk::Path rel_dir = Vtk::relative(data_dir, fn_dir);
std::string serial_fn = data_dir.string() + '/' + name.string();
std::string parallel_fn = fn_dir.string() + '/' + name.string();
@@ -44,7 +44,7 @@ std::string VtkWriterInterface<GV,DC>
serial_fn += "_p" + std::to_string(comm().rank());
std::string outputFilename;
-
+
{ // write serial file
outputFilename = serial_fn + "." + fileExtension();
std::ofstream serial_out(outputFilename, std::ios_base::ate | std::ios::binary);
diff --git a/dune/vtk/writers/CMakeLists.txt b/dune/vtk/writers/CMakeLists.txt
index 779237467679ed664f32981dfa76537d2d2cb4fb..dc57461be4b7b0956258cdb06be527610df9dfa4 100644
--- a/dune/vtk/writers/CMakeLists.txt
+++ b/dune/vtk/writers/CMakeLists.txt
@@ -8,4 +8,4 @@ install(FILES
vtkstructuredgridwriter.impl.hh
vtkunstructuredgridwriter.hh
vtkunstructuredgridwriter.impl.hh
- DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtkwriter/writers)
+ DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtk/writers)
diff --git a/dune/vtk/writers/vtkimagedatawriter.hh b/dune/vtk/writers/vtkimagedatawriter.hh
index ed7cf1fc2027909dffb6e0ac6a09241c94151734..c0d4e45e21f76b7d916e85237f5a5678c4c427bc 100644
--- a/dune/vtk/writers/vtkimagedatawriter.hh
+++ b/dune/vtk/writers/vtkimagedatawriter.hh
@@ -6,8 +6,8 @@
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/datacollectors/structureddatacollector.hh>
#include <dune/vtk/vtkwriterinterface.hh>
@@ -62,7 +62,7 @@ namespace Dune
template <class GridView,
class = std::void_t<typename GridView::IndexSet>>
VtkImageDataWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
- -> VtkImageDataWriter<GridView, StructuredDataCollector<GridView>>;
+ -> VtkImageDataWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
template <class DataCollector,
class = std::void_t<typename DataCollector::GridView>>
diff --git a/dune/vtk/writers/vtkimagedatawriter.impl.hh b/dune/vtk/writers/vtkimagedatawriter.impl.hh
index 898b076ca5fdfc0d43988cf5057b9b854a111cc8..03c1ea196b580874b2b7357320cdcd1109fa67bf 100644
--- a/dune/vtk/writers/vtkimagedatawriter.impl.hh
+++ b/dune/vtk/writers/vtkimagedatawriter.impl.hh
@@ -27,13 +27,13 @@ void VtkImageDataWriter<GV,DC>
auto const& origin = dataCollector_->origin();
auto const& spacing = dataCollector_->spacing();
out << "<ImageData"
- << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
- << " Origin=\"" << join(origin.begin(), origin.end()) << "\""
- << " Spacing=\"" << join(spacing.begin(), spacing.end()) << "\""
+ << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
+ << " Origin=\"" << Vtk::join(origin.begin(), origin.end()) << "\""
+ << " Spacing=\"" << Vtk::join(spacing.begin(), spacing.end()) << "\""
<< ">\n";
dataCollector_->writeLocalPiece([&out](auto const& extent) {
- out << "<Piece Extent=\"" << join(extent.begin(), extent.end()) << "\">\n";
+ out << "<Piece Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\">\n";
});
// Write data associated with grid points
@@ -67,9 +67,9 @@ void VtkImageDataWriter<GV,DC>
auto const& spacing = dataCollector_->spacing();
out << "<PImageData"
<< " GhostLevel=\"" << dataCollector_->ghostLevel() << "\""
- << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
- << " Origin=\"" << join(origin.begin(), origin.end()) << "\""
- << " Spacing=\"" << join(spacing.begin(), spacing.end()) << "\""
+ << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
+ << " Origin=\"" << Vtk::join(origin.begin(), origin.end()) << "\""
+ << " Spacing=\"" << Vtk::join(spacing.begin(), spacing.end()) << "\""
<< ">\n";
// Write data associated with grid points
@@ -100,7 +100,7 @@ void VtkImageDataWriter<GV,DC>
std::string piece_source = pfilename + "_p" + std::to_string(p) + "." + ext;
out << "<Piece Source=\"" << piece_source << "\"";
if (write_extent)
- out << " Extent=\"" << join(extent.begin(), extent.end()) << "\"";
+ out << " Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\"";
out << " />\n";
});
diff --git a/dune/vtk/writers/vtkrectilineargridwriter.hh b/dune/vtk/writers/vtkrectilineargridwriter.hh
index f119a0501f54adff0987b4e7112fb4d386f8300a..af88be7928da7148c14a683a5033b950b72f6210 100644
--- a/dune/vtk/writers/vtkrectilineargridwriter.hh
+++ b/dune/vtk/writers/vtkrectilineargridwriter.hh
@@ -6,8 +6,8 @@
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/datacollectors/structureddatacollector.hh>
#include <dune/vtk/vtkwriterinterface.hh>
@@ -68,7 +68,7 @@ namespace Dune
template <class GridView,
class = std::void_t<typename GridView::IndexSet>>
VtkRectilinearGridWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
- -> VtkRectilinearGridWriter<GridView, StructuredDataCollector<GridView>>;
+ -> VtkRectilinearGridWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
template <class DataCollector,
class = std::void_t<typename DataCollector::GridView>>
diff --git a/dune/vtk/writers/vtkrectilineargridwriter.impl.hh b/dune/vtk/writers/vtkrectilineargridwriter.impl.hh
index f2c822382e96f4c5d4f0c7e4d991fbf44146e483..d26428dc8724a8731d826f298626d480ccc3cc1f 100644
--- a/dune/vtk/writers/vtkrectilineargridwriter.impl.hh
+++ b/dune/vtk/writers/vtkrectilineargridwriter.impl.hh
@@ -25,11 +25,11 @@ void VtkRectilinearGridWriter<GV,DC>
auto const& wholeExtent = dataCollector_->wholeExtent();
out << "<RectilinearGrid"
- << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
+ << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
<< ">\n";
dataCollector_->writeLocalPiece([&out](auto const& extent) {
- out << "<Piece Extent=\"" << join(extent.begin(), extent.end()) << "\">\n";
+ out << "<Piece Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\">\n";
});
// Write point coordinates for x, y, and z ordinate
@@ -66,7 +66,7 @@ void VtkRectilinearGridWriter<GV,DC>
auto const& wholeExtent = dataCollector_->wholeExtent();
out << "<PRectilinearGrid"
<< " GhostLevel=\"" << dataCollector_->ghostLevel() << "\""
- << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
+ << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
<< ">\n";
// Write point coordinates for x, y, and z ordinate
@@ -104,7 +104,7 @@ void VtkRectilinearGridWriter<GV,DC>
std::string piece_source = pfilename + "_p" + std::to_string(p) + "." + ext;
out << "<Piece Source=\"" << piece_source << "\"";
if (write_extent)
- out << " Extent=\"" << join(extent.begin(), extent.end()) << "\"";
+ out << " Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\"";
out << " />\n";
});
diff --git a/dune/vtk/writers/vtkstructuredgridwriter.hh b/dune/vtk/writers/vtkstructuredgridwriter.hh
index 21d23134e8fb166a15ac0e6afdcd93261118cb99..c0b37d366935b57bd14ba16d984fe9530a5b295b 100644
--- a/dune/vtk/writers/vtkstructuredgridwriter.hh
+++ b/dune/vtk/writers/vtkstructuredgridwriter.hh
@@ -6,8 +6,8 @@
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/datacollectors/structureddatacollector.hh>
#include <dune/vtk/vtkwriterinterface.hh>
@@ -62,7 +62,7 @@ namespace Dune
template <class GridView,
class = std::void_t<typename GridView::IndexSet>>
VtkStructuredGridWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
- -> VtkStructuredGridWriter<GridView, StructuredDataCollector<GridView>>;
+ -> VtkStructuredGridWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
template <class DataCollector,
class = std::void_t<typename DataCollector::GridView>>
diff --git a/dune/vtk/writers/vtkstructuredgridwriter.impl.hh b/dune/vtk/writers/vtkstructuredgridwriter.impl.hh
index ec8707990614f24175c6248950d65d0d7148686e..0d09b88d5cead0fe046c7dd0134b61443af5b506 100644
--- a/dune/vtk/writers/vtkstructuredgridwriter.impl.hh
+++ b/dune/vtk/writers/vtkstructuredgridwriter.impl.hh
@@ -24,10 +24,10 @@ void VtkStructuredGridWriter<GV,DC>
this->writeHeader(out, "StructuredGrid");
auto const& wholeExtent = dataCollector_->wholeExtent();
- out << "<StructuredGrid WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\">\n";
+ out << "<StructuredGrid WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\">\n";
dataCollector_->writeLocalPiece([&out](auto const& extent) {
- out << "<Piece Extent=\"" << join(extent.begin(), extent.end()) << "\">\n";
+ out << "<Piece Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\">\n";
});
// Write point coordinates
@@ -64,7 +64,7 @@ void VtkStructuredGridWriter<GV,DC>
auto const& wholeExtent = dataCollector_->wholeExtent();
out << "<PStructuredGrid"
<< " GhostLevel=\"" << dataCollector_->ghostLevel() << "\""
- << " WholeExtent=\"" << join(wholeExtent.begin(), wholeExtent.end()) << "\""
+ << " WholeExtent=\"" << Vtk::join(wholeExtent.begin(), wholeExtent.end()) << "\""
<< ">\n";
// Write points
@@ -103,7 +103,7 @@ void VtkStructuredGridWriter<GV,DC>
std::string piece_source = pfilename + "_p" + std::to_string(p) + "." + ext;
out << "<Piece Source=\"" << piece_source << "\"";
if (write_extent)
- out << " Extent=\"" << join(extent.begin(), extent.end()) << "\"";
+ out << " Extent=\"" << Vtk::join(extent.begin(), extent.end()) << "\"";
out << " />\n";
});
diff --git a/dune/vtk/writers/vtkunstructuredgridwriter.hh b/dune/vtk/writers/vtkunstructuredgridwriter.hh
index a1d557ba9a35e04c43369c123c9e257d0a598482..8bd16aa190f9062205dac9b29c9d9f37520ef043 100644
--- a/dune/vtk/writers/vtkunstructuredgridwriter.hh
+++ b/dune/vtk/writers/vtkunstructuredgridwriter.hh
@@ -6,8 +6,8 @@
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/forward.hh>
-#include <dune/vtk/vtkfunction.hh>
-#include <dune/vtk/vtktypes.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/types.hh>
#include <dune/vtk/datacollectors/continuousdatacollector.hh>
#include <dune/vtk/vtkwriterinterface.hh>
@@ -94,7 +94,7 @@ namespace Dune
template <class GridView,
class = std::void_t<typename GridView::IndexSet>>
VtkUnstructuredGridWriter(GridView const&, Vtk::FormatTypes = Vtk::BINARY, Vtk::DataTypes = Vtk::FLOAT32)
- -> VtkUnstructuredGridWriter<GridView, ContinuousDataCollector<GridView>>;
+ -> VtkUnstructuredGridWriter<GridView, Vtk::ContinuousDataCollector<GridView>>;
template <class DataCollector,
class = std::void_t<typename DataCollector::GridView>>
diff --git a/src/datacollector.cc b/src/datacollector.cc
index 9f126d2c54c19e59ecc9050ff6f8e9c606cedd30..c342a24a2fe1c067459cb0aa5b750a737c40e334 100644
--- a/src/datacollector.cc
+++ b/src/datacollector.cc
@@ -21,7 +21,7 @@
#if HAVE_UG
#include <dune/grid/uggrid.hh>
-#endif
+#endif
#include <dune/grid/yaspgrid.hh>
#include <dune/grid/utility/structuredgridfactory.hh>
@@ -73,12 +73,12 @@ void write (std::string prefix, GridView const& gridView)
// write analytic function
auto p1Analytic = makeAnalyticGridViewFunction([&c](auto const& x) { return c.dot(x); }, gridView);
- write_dc<ContinuousDataCollector<GridView>>(prefix + "_continuous", gridView, p1Interpol, p1Analytic);
- write_dc<DiscontinuousDataCollector<GridView>>(prefix + "_discontinuous", gridView, p1Interpol, p1Analytic);
- write_dc<QuadraticDataCollector<GridView>>(prefix + "_quadratic", gridView, p1Interpol, p1Analytic);
+ write_dc<Vtk::ContinuousDataCollector<GridView>>(prefix + "_continuous", gridView, p1Interpol, p1Analytic);
+ write_dc<Vtk::DiscontinuousDataCollector<GridView>>(prefix + "_discontinuous", gridView, p1Interpol, p1Analytic);
+ write_dc<Vtk::QuadraticDataCollector<GridView>>(prefix + "_quadratic", gridView, p1Interpol, p1Analytic);
Hybrid::forEach(StaticIntegralRange<int,7,1>{}, [&](auto p) {
- write_dc<LagrangeDataCollector<GridView,p>>(prefix + "_lagrange_p" + std::to_string(p), gridView, p1Interpol, p1Analytic);
+ write_dc<Vtk::LagrangeDataCollector<GridView,p>>(prefix + "_lagrange_p" + std::to_string(p), gridView, p1Interpol, p1Analytic);
});
}
diff --git a/src/lagrangepoints.cc b/src/lagrangepoints.cc
index a40f83e80c6e84e0f2f324d2d6becb7032d2b39b..9839ebb40dac1f51a53422d8ba049da5a62e6d3c 100644
--- a/src/lagrangepoints.cc
+++ b/src/lagrangepoints.cc
@@ -25,7 +25,7 @@ void write (std::string prefix, index_constant<dim>)
for (int order = 1; order < 6; ++order) {
std::cout << "order: " << order << std::endl;
{
- VtkLagrangePointSet<double, dim> pointSet(order);
+ Vtk::LagrangePointSet<double, dim> pointSet(order);
pointSet.build(GeometryTypes::cube(dim));
std::size_t i = 0;
@@ -34,9 +34,9 @@ void write (std::string prefix, index_constant<dim>)
std::cout << i++ << ") p = " << p.point() << ", key = " << p.localKey() << std::endl;
}
- using BasisF = LagrangeBasisFactory<VtkLagrangePointSet, dim, double, double>;
- using CoefficientF = LagrangeCoefficientsFactory<VtkLagrangePointSet, dim, double>;
- using InterpolationF = LagrangeInterpolationFactory<VtkLagrangePointSet, dim, double>;
+ using BasisF = LagrangeBasisFactory<Vtk::LagrangePointSet, dim, double, double>;
+ using CoefficientF = LagrangeCoefficientsFactory<Vtk::LagrangePointSet, dim, double>;
+ using InterpolationF = LagrangeInterpolationFactory<Vtk::LagrangePointSet, dim, double>;
GenericLocalFiniteElement<BasisF, CoefficientF, InterpolationF> localFE(GeometryTypes::cube(dim), order);
auto const& localBasis = localFE.localBasis();
@@ -45,7 +45,7 @@ void write (std::string prefix, index_constant<dim>)
}
{
- VtkLagrangePointSet<double, dim> pointSet(order);
+ Vtk::LagrangePointSet<double, dim> pointSet(order);
pointSet.build(GeometryTypes::simplex(dim));
std::size_t i = 0;
diff --git a/src/lagrangereader.cc b/src/lagrangereader.cc
index 8a3543cda24742ac3436e296ddfc1dd23238f8f1..fb1af5e2486a737ccee6e66ee8eaba13ee123adb 100644
--- a/src/lagrangereader.cc
+++ b/src/lagrangereader.cc
@@ -42,8 +42,8 @@ int main(int argc, char** argv)
// using GridType = Dune::ALUGrid<dim,dow,Dune::simplex,Dune::conforming>;
using GridType = FoamGrid<dim,dow>;
using GridView = typename GridType::LeafGridView;
- using DataCollector = LagrangeDataCollector<GridView, order>;
- using GridCreator = LagrangeGridCreator<GridType>;
+ using DataCollector = Vtk::LagrangeDataCollector<GridView, order>;
+ using GridCreator = Vtk::LagrangeGridCreator<GridType>;
std::string filename = "triangles_" + std::to_string(dow) + "d_order" + std::to_string(order);
@@ -51,7 +51,7 @@ int main(int argc, char** argv)
{ // Test using the (static) file-reader interface
std::cout << "Test 1..." << std::endl;
- auto gridPtr = VtkReader<GridType,GridCreator>::read(GRID_PATH "/" + filename + ".vtu");
+ auto gridPtr = VtkReader<GridType,GridCreator>::createGridFromFile(GRID_PATH "/" + filename + ".vtu");
auto& grid = *gridPtr;
grid.globalRefine(2);
@@ -63,7 +63,7 @@ int main(int argc, char** argv)
std::cout << "Test 2..." << std::endl;
GridFactory<GridType> factory;
VtkReader<GridType,GridCreator> reader(factory);
- reader.readFromFile(GRID_PATH "/" + filename + ".vtu");
+ reader.read(GRID_PATH "/" + filename + ".vtu");
auto gridPtr = factory.createGrid();
auto& grid = *gridPtr;
@@ -80,7 +80,7 @@ int main(int argc, char** argv)
GridFactory<GridType> factory;
GridCreator creator(factory);
VtkReader<GridType,GridCreator> reader(creator);
- reader.readFromFile(GRID_PATH "/" + filename + ".vtu");
+ reader.read(GRID_PATH "/" + filename + ".vtu");
auto gridPtr = factory.createGrid();
auto& grid = *gridPtr;
diff --git a/src/test/mixed_element_test.cc b/src/test/mixed_element_test.cc
index 0a5b446cb5d5846c3542d26dd3325dc0bc5eb65c..6dbd9e104eb2ad5b7d7d99c0e6d07a6e00adf77e 100644
--- a/src/test/mixed_element_test.cc
+++ b/src/test/mixed_element_test.cc
@@ -82,7 +82,7 @@ template <class Grid, class Test>
void reader_test (Test& test)
{
for (auto const& test_case : test_cases) {
- auto grid = VtkReader<Grid>::read("mixed_element_test_" + std::get<0>(test_case) + ".vtu");
+ auto grid = VtkReader<Grid>::createGridFromFile("mixed_element_test_" + std::get<0>(test_case) + ".vtu");
VtkUnstructuredGridWriter<typename Grid::LeafGridView> vtkWriter(grid->leafGridView(),
std::get<1>(test_case), std::get<2>(test_case));
vtkWriter.write("mixed_element_test_" + std::get<0>(test_case) + "_2.vtu");
diff --git a/src/test/parallel_reader_writer_test.cc b/src/test/parallel_reader_writer_test.cc
index 72dd52bbf4aa832381fa056a4097822424a7add3..d46d1371595ffd3dfded11c20ff611791db996bb 100644
--- a/src/test/parallel_reader_writer_test.cc
+++ b/src/test/parallel_reader_writer_test.cc
@@ -79,7 +79,7 @@ using HasParallelGridFactory = Std::is_detected<HasParallelGridFactoryImpl, Grid
template <class Test>
-void compare (Test& test, filesystem::path const& dir, filesystem::path const& name)
+void compare (Test& test, Vtk::Path const& dir, Vtk::Path const& name)
{
test.check(compare_files(dir.string() + '/' + name.string() + ".vtu",
dir.string() + '/' + name.string() + "_2.vtu"));
@@ -119,7 +119,7 @@ void reader_writer_test(MPIHelper& mpi, TestSuite& test, std::string const& test
// Step 2: read the grid from file1 and write it back to file2
{ GridFactory<Grid> factory;
VtkReader<Grid, Creator> reader{factory};
- reader.readFromFile(base_name + ext);
+ reader.read(base_name + ext);
std::unique_ptr<Grid> grid{ Hybrid::ifElse(HasParallelGridFactory<Grid>{},
[&](auto id) { return id(factory).createGrid(std::true_type{}); },
@@ -135,7 +135,7 @@ void reader_writer_test(MPIHelper& mpi, TestSuite& test, std::string const& test
// Step 3: read the (parallel) file1 to get the piece filenames
{ GridFactory<Grid> factory;
VtkReader<Grid, Creator> reader{factory};
- reader.readFromFile(base_name + "_1" + ext);
+ reader.read(base_name + "_1" + ext);
std::unique_ptr<Grid> grid{ Hybrid::ifElse(HasParallelGridFactory<Grid>{},
[&](auto id) { return id(factory).createGrid(std::true_type{}); },
@@ -152,7 +152,7 @@ void reader_writer_test(MPIHelper& mpi, TestSuite& test, std::string const& test
// Step 4: read the (parallel) file2 to get the piece filenames
{ GridFactory<Grid> factory;
VtkReader<Grid, Creator> reader{factory};
- reader.readFromFile(base_name + "_2" + ext, false);
+ reader.read(base_name + "_2" + ext, false);
pieces2 = reader.pieces();
}
@@ -182,19 +182,19 @@ int main (int argc, char** argv)
TestSuite test{};
#if HAVE_UG
- reader_writer_test<UGGrid<2>, SerialGridCreator<UGGrid<2>>>(mpi, test, "UGGrid<2>");
- reader_writer_test<UGGrid<3>, SerialGridCreator<UGGrid<3>>>(mpi, test, "UGGrid<3>");
+ reader_writer_test<UGGrid<2>, Vtk::SerialGridCreator<UGGrid<2>>>(mpi, test, "UGGrid<2>");
+ reader_writer_test<UGGrid<3>, Vtk::SerialGridCreator<UGGrid<3>>>(mpi, test, "UGGrid<3>");
#endif
#if HAVE_DUNE_ALUGRID
// Test VtkWriter for ALUGrid.
- reader_writer_test<ALUGridType<2>, SerialGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2>");
- reader_writer_test<ALUGridType<2>, ParallelGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2, Parallel>", false);
+ reader_writer_test<ALUGridType<2>, Vtk::SerialGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2>");
+ reader_writer_test<ALUGridType<2>, Vtk::ParallelGridCreator<ALUGridType<2>>>(mpi, test, "ALUGridType<2, Parallel>", false);
- reader_writer_test<ALUGridType<3>, SerialGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3>");
+ reader_writer_test<ALUGridType<3>, Vtk::SerialGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3>");
#if DUNE_VERSION_LT(DUNE_GRID,2,7)
// Currently the 2.7 branch is not working, due to a new bisection compatibility check in 3d
- reader_writer_test<ALUGridType<3>, ParallelGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3, Parallel>", false);
+ reader_writer_test<ALUGridType<3>, Vtk::ParallelGridCreator<ALUGridType<3>>>(mpi, test, "ALUGridType<3, Parallel>", false);
#endif
#endif
diff --git a/src/test/reader_writer_test.cc b/src/test/reader_writer_test.cc
index c7edefb7c7baee5c3e46b8fac5b2948590223e49..927a660073ca8ea84d25c1cd9b3243258dca4582 100644
--- a/src/test/reader_writer_test.cc
+++ b/src/test/reader_writer_test.cc
@@ -80,7 +80,7 @@ static TestCases test_cases = {
};
template <class Test>
-void compare (Test& test, filesystem::path const& dir, filesystem::path const& name)
+void compare (Test& test, Vtk::Path const& dir, Vtk::Path const& name)
{
test.check(compare_files(dir.string() + '/' + name.string() + ".vtu",
dir.string() + '/' + name.string() + "_2.vtu"));
@@ -105,7 +105,7 @@ void reader_test (MPIHelper& mpi, Test& test)
{
GridFactory<Grid> factory;
VtkReader<Grid> reader{factory};
- reader.readFromFile("reader_writer_test_" + std::get<0>(test_case) + ext);
+ reader.read("reader_writer_test_" + std::get<0>(test_case) + ext);
std::unique_ptr<Grid> grid{factory.createGrid()};
pieces1 = reader.pieces();
@@ -118,7 +118,7 @@ void reader_test (MPIHelper& mpi, Test& test)
{
GridFactory<Grid> factory2;
VtkReader<Grid> reader2{factory2};
- reader2.readFromFile("reader_writer_test_" + std::get<0>(test_case) + "_2" + ext, false);
+ reader2.read("reader_writer_test_" + std::get<0>(test_case) + "_2" + ext, false);
pieces2 = reader2.pieces();
}
diff --git a/src/vtkreader.cc b/src/vtkreader.cc
index 5ac2fd7b0ca3e47fce50393bad1ea80887fc7959..db4dd25e48472203fbf418e429200868312633e3 100644
--- a/src/vtkreader.cc
+++ b/src/vtkreader.cc
@@ -51,7 +51,7 @@ int main(int argc, char** argv)
{
std::cout << "read 'test_ascii_float32.vtu'..." << std::endl;
- auto gridPtr = VtkReader<GridType>::read("test_ascii_float32.vtu");
+ auto gridPtr = VtkReader<GridType>::createGridFromFile("test_ascii_float32.vtu");
auto& grid = *gridPtr;
VtkUnstructuredGridWriter<GridView> vtkWriter(grid.leafGridView(), Vtk::ASCII);
@@ -60,7 +60,7 @@ int main(int argc, char** argv)
{
std::cout << "read 'test_binary_float32.vtu'..." << std::endl;
- auto gridPtr = VtkReader<GridType>::read("test_binary_float32.vtu");
+ auto gridPtr = VtkReader<GridType>::createGridFromFile("test_binary_float32.vtu");
auto& grid = *gridPtr;
VtkUnstructuredGridWriter<GridView> vtkWriter(grid.leafGridView(), Vtk::ASCII);
@@ -69,7 +69,7 @@ int main(int argc, char** argv)
{
std::cout << "read 'test_compressed_float64.vtu'..." << std::endl;
- auto gridPtr = VtkReader<GridType>::read("test_compressed_float64.vtu");
+ auto gridPtr = VtkReader<GridType>::createGridFromFile("test_compressed_float64.vtu");
auto& grid = *gridPtr;
VtkUnstructuredGridWriter<GridView> vtkWriter(grid.leafGridView(), Vtk::ASCII, Vtk::FLOAT64);
@@ -78,31 +78,31 @@ int main(int argc, char** argv)
{
std::cout << "read 'test_ascii_float32.vtu'..." << std::endl;
- auto gridPtr = VtkReader<GridType,DiscontinuousGridCreator<GridType>>::read("test_ascii_float32.vtu");
+ auto gridPtr = VtkReader<GridType,Vtk::DiscontinuousGridCreator<GridType>>::createGridFromFile("test_ascii_float32.vtu");
auto& grid = *gridPtr;
VtkUnstructuredGridWriter<GridView> vtkWriter(grid.leafGridView(), Vtk::ASCII);
vtkWriter.write("test_ascii_float32_4.vtu");
}
- if (filesystem::exists("paraview_3d.vtu")) {
+ if (Vtk::exists("paraview_3d.vtu")) {
using GridType3d = UGGrid<3>;
using GridView3d = typename GridType3d::LeafGridView;
std::cout << "read 'paraview_3d.vtu'..." << std::endl;
- auto gridPtr = VtkReader<GridType3d>::read("paraview_3d.vtu");
+ auto gridPtr = VtkReader<GridType3d>::createGridFromFile("paraview_3d.vtu");
auto& grid = *gridPtr;
VtkUnstructuredGridWriter<GridView3d> vtkWriter(grid.leafGridView(), Vtk::ASCII, Vtk::FLOAT64);
vtkWriter.write("paraview_3d_ascii.vtu");
}
- if (filesystem::exists("paraview_2d.vtu")) {
+ if (Vtk::exists("paraview_2d.vtu")) {
std::cout << "paraview_2d_ascii...\n";
using GridType2d = UGGrid<2>;
using GridView2d = typename GridType2d::LeafGridView;
std::cout << "read 'paraview_2d.vtu'..." << std::endl;
- auto gridPtr = VtkReader<GridType2d>::read("paraview_2d.vtu");
+ auto gridPtr = VtkReader<GridType2d>::createGridFromFile("paraview_2d.vtu");
auto& grid = *gridPtr;
VtkUnstructuredGridWriter<GridView2d> vtkWriter(grid.leafGridView(), Vtk::ASCII, Vtk::FLOAT64);