Merge branch 'feature/test_reader_gridfunction' into 'master'

Add test for compatibility of grid-functions with Vtk::Function

See merge request extensions/dune-vtk!15

dune/vtk/function.hh

@@ -149,6 +149,19 @@ namespace Dune
         : Function(std::forward<F>(fct), info.name(), info.size(), info.rangeType(), info.dataType())
       {}
 
+      /// (6) Automatically extract name and num components from GridFunction if available
+      template <class F, class... Args,
+        disableCopyMove<Function, F> = 0,
+        class = decltype(localFunction(std::declval<F>())),
+        class = decltype(std::declval<F>().name()),
+        class = decltype(std::declval<F>().numComponents()),
+        class = decltype(std::declval<F>().dataType()) >
+      explicit Function (F&& fct, ...)
+        : Function(localFunction(std::forward<F>(fct)), fct.name(), fct.numComponents(),
+                   Vtk::RangeTypes::UNSPECIFIED, fct.dataType())
+      {}
+
+
       /// (7) Construct from legacy VTKFunction
       /**
       * \param fct  The Dune::VTKFunction to wrap

dune/vtk/gridcreators/lagrangegridcreator.hh

@@ -276,6 +276,21 @@ namespace Dune
         return LocalFunction{gridCreator};
       }
 
+      std::string name () const
+      {
+        return "LagrangeParametrization";
+      }
+
+      int numComponents () const
+      {
+        return GlobalCoordinate::size();
+      }
+
+      Vtk::DataTypes dataType () const
+      {
+        return dataTypeOf<typename Element::Geometry::ctype>();
+      }
+
       struct EntitySet
       {
         using Grid = GridType;

dune/vtk/gridfunctions/CMakeLists.txt

@@ -4,3 +4,5 @@ install(FILES
   continuousgridfunction.hh
   lagrangegridfunction.hh
   DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/dune/vtk/gridfunctions)
+
+add_subdirectory(test)
\ No newline at end of file

dune/vtk/gridfunctions/continuousgridfunction.hh

@@ -159,13 +159,16 @@ namespace Dune
 
     public:
       template <class GridCreator>
-      ContinuousGridFunction (GridCreator const& creator, std::vector<Field> const& values, unsigned int comp,
+      ContinuousGridFunction (GridCreator const& creator, std::vector<Field> const& values,
+                              std::string name, unsigned int ncomps, Vtk::DataTypes dataType,
                               std::vector<std::uint8_t> const& /*types*/,
                               std::vector<std::int64_t> const& /*offsets*/,
                               std::vector<std::int64_t> const& /*connectivity*/)
         : factory_(&creator.factory())
         , values_(&values)
-        , comp_(comp)
+        , name_(std::move(name))
+        , ncomps_(ncomps)
+        , dataType_(dataType)
       {}
 
       ContinuousGridFunction () = default;
@@ -173,7 +176,7 @@ namespace Dune
       Range operator() (Domain const& global) const
       {
         DUNE_THROW(Dune::NotImplemented, "Evaluation in global coordinates not implemented.");
-        return Range(comp_, 0);
+        return Range(ncomps_, 0);
       }
 
       EntitySet const& entitySet () const
@@ -181,15 +184,32 @@ namespace Dune
         return entitySet_;
       }
 
+      std::string const& name () const
+      {
+        return name_;
+      }
+
+      int numComponents () const
+      {
+        return ncomps_;
+      }
+
+      Vtk::DataTypes dataType () const
+      {
+        return dataType_;
+      }
+
       friend LocalFunction<typename EntitySet::Element> localFunction (ContinuousGridFunction const& gf)
       {
-        return {gf.factory_, gf.values_, gf.comp_};
+        return {gf.factory_, gf.values_, gf.ncomps_};
       }
 
     private:
       Factory const* factory_;
       std::vector<Field> const* values_ = nullptr;
-      unsigned int comp_ = 0;
+      std::string name_ = "GridFunction";
+      unsigned int ncomps_ = 0;
+      Vtk::DataTypes dataType_ = Vtk::DataTypes::UNKNOWN;
 
       EntitySet entitySet_;
     };

dune/vtk/gridfunctions/lagrangegridfunction.hh

@@ -212,13 +212,16 @@ namespace Dune
     public:
       /// Construct a grid-function. Passed in data is stroed by reference, thus must have
       /// a life-time greater than that of the grid-function and corresponding local-function.
-      LagrangeGridFunction (GridCreator const& creator, std::vector<Field> const& values, unsigned int comp,
+      LagrangeGridFunction (GridCreator const& creator, std::vector<Field> const& values,
+                            std::string name ,unsigned int ncomps, Vtk::DataTypes dataType,
                             std::vector<std::uint8_t> const& types,
                             std::vector<std::int64_t> const& offsets,
                             std::vector<std::int64_t> const& connectivity)
         : creator_(&creator)
         , values_(&values)
-        , comp_(comp)
+        , name_(std::move(name))
+        , ncomps_(ncomps)
+        , dataType_(dataType)
         , types_(&types)
         , offsets_(&offsets)
         , connectivity_(&connectivity)
@@ -230,7 +233,7 @@ namespace Dune
       Range operator() (Domain const& global) const
       {
         DUNE_THROW(Dune::NotImplemented, "Evaluation in global coordinates not implemented.");
-        return Range(comp_, 0);
+        return Range(ncomps_, 0);
       }
 
       /// Return a type that defines the element that can be iterated.
@@ -239,17 +242,34 @@ namespace Dune
         return entitySet_;
       }
 
+      std::string const& name () const
+      {
+        return name_;
+      }
+
+      int numComponents () const
+      {
+        return ncomps_;
+      }
+
+      Vtk::DataTypes dataType () const
+      {
+        return dataType_;
+      }
+
       /// Construct a local-function depending on the Context type either PointDataLocalFunction
       /// or CellDataLocalFunction
       friend LocalFunction<typename EntitySet::Element> localFunction (LagrangeGridFunction const& gf)
       {
-        return {gf.creator_, gf.values_, gf.comp_, gf.types_, gf.offsets_, gf.connectivity_};
+        return {gf.creator_, gf.values_, gf.ncomps_, gf.types_, gf.offsets_, gf.connectivity_};
       }
 
     private:
       GridCreator const* creator_ = nullptr;
       std::vector<Field> const* values_ = nullptr;
-      unsigned int comp_ = 0;
+      std::string name_ = "GridFunction";
+      unsigned int ncomps_ = 0;
+      Vtk::DataTypes dataType_ = Vtk::DataTypes::UNKNOWN;
       std::vector<std::uint8_t> const* types_ = nullptr;
       std::vector<std::int64_t> const* offsets_ = nullptr;
       std::vector<std::int64_t> const* connectivity_ = nullptr;

dune/vtk/gridfunctions/test/CMakeLists.txt

+dune_add_test(SOURCES test-function-gridfunction.cc
+  LINK_LIBRARIES dunevtk)

dune/vtk/gridfunctions/test/test-function-gridfunction.cc

+#include <config.h>
+
+#if HAVE_DUNE_UGGRID
+  #include <dune/grid/uggrid.hh>
+  using GridType = Dune::UGGrid<2>;
+#else
+  #include <dune/grid/yaspgrid.hh>
+  using GridType = Dune::YaspGrid<2>;
+#endif
+
+#include <dune/grid/utility/structuredgridfactory.hh>
+#include <dune/vtk/vtkreader.hh>
+#include <dune/vtk/vtkwriter.hh>
+#include <dune/vtk/function.hh>
+#include <dune/vtk/gridcreators/lagrangegridcreator.hh>
+#include <dune/vtk/utility/errors.hh>
+
+
+// Wrapper for global-coordinate functions F
+template <class GridView, class F>
+class GlobalFunction
+{
+  using Element = typename GridView::template Codim<0>::Entity;
+  using Geometry = typename Element::Geometry;
+
+public:
+  GlobalFunction (GridView const& gridView, F const& f)
+    : gridView_(gridView)
+    , f_(f)
+  {}
+
+  void bind(Element const& element) { geometry_.emplace(element.geometry()); }
+  void unbind() { geometry_.reset(); }
+
+  auto operator() (typename Geometry::LocalCoordinate const& local) const
+  {
+    assert(!!geometry_);
+    return f_(geometry_->global(local));
+  }
+
+private:
+  GridView gridView_;
+  F f_;
+  std::optional<Geometry> geometry_;
+};
+
+int main(int argc, char** argv)
+{
+  using namespace Dune;
+  MPIHelper::instance(argc, argv);
+
+  { // write point and cell data to a file
+
+    auto grid = StructuredGridFactory<GridType>::createCubeGrid({0.0,0.0}, {1.0,2.0}, {2u,4u});
+    auto gridView = grid->leafGridView();
+
+    VtkUnstructuredGridWriter writer(gridView);
+
+    auto f = GlobalFunction{gridView, [](auto x) { return x[0] + x[1]; }};
+    writer.addPointData(f, "pointdata");
+    writer.addCellData(f, "celldata");
+    writer.write("test-function-gridfunction.vtu");
+  }
+
+  { // read data from file
+
+    Vtk::LagrangeGridCreator<GridType> gridCreator;
+    VtkReader reader(gridCreator);
+    reader.read("test-function-gridfunction.vtu");
+
+    auto grid = reader.createGrid();
+    auto gridView = grid->leafGridView();
+    using GridView = decltype(gridView);
+
+    VtkUnstructuredGridWriter writer(gridView);
+
+    // store point-data and cell-data gridfunction into a Vtk::Function
+    Vtk::Function<GridView> pd{ reader.getPointData("pointdata") };
+    Vtk::Function<GridView> cd{ reader.getCellData("celldata") };
+
+    // store the grid-parametrization grid-function into a Vtk::Function
+    Vtk::Function<GridView> param{ gridCreator };
+
+    writer.addPointData(reader.getPointData("pointdata"));
+    writer.addCellData(reader.getCellData("celldata"));
+    writer.addPointData(gridCreator);
+
+    auto localParam = localFunction(param);
+    auto localPd = localFunction(pd);
+    auto localCd = localFunction(cd);
+    for (auto const& element : elements(gridView)) {
+      localParam.bind(element);
+      localPd.bind(element);
+      localCd.bind(element);
+
+      auto refElem = Dune::referenceElement(element);
+
+      // evaluate grid-functions in local coordinates
+      auto x0 = localParam(refElem.position(0,0));
+      VTK_ASSERT(x0.size() == 2);
+
+      auto x1 = localPd(refElem.position(0,0));
+      VTK_ASSERT(x1.size() == 1);
+
+      auto x2 = localCd(refElem.position(0,0));
+      VTK_ASSERT(x2.size() == 1);
+    }
+  }
+}
\ No newline at end of file

dune/vtk/utility/CMakeLists.txt

@@ -4,6 +4,7 @@ dune_add_library("filesystem" OBJECT
 #install headers
 install(FILES
   arguments.hh
+  concepts.hh
   enum.hh
   errors.hh
   filesystem.hh

dune/vtk/utility/concepts.hh

+#pragma once
+
+#include <type_traits>
+#include <dune/geometry/type.hh>
+
+namespace Dune
+{
+  namespace Vtk
+  {
+    template <class...> struct CheckTypes {};
+
+    template <class DataCollector, class DC = std::decay_t<DataCollector>>
+    using IsDataCollector = decltype((
+      std::declval<DC&>().update(),
+      std::declval<DC>().numPoints(),
+      std::declval<DC>().numCells(),
+      CheckTypes<typename DC::GridView>{},
+    true));
+
+    template <class GridView, class GV = std::decay_t<GridView>>
+    using IsGridView = decltype((
+      std::declval<GV>().grid(),
+      std::declval<GV>().indexSet(),
+      std::declval<GV>().size(0),
+      std::declval<GV>().size(std::declval<Dune::GeometryType>()),
+      CheckTypes<typename GV::Grid, typename GV::IndexSet>{},
+    true));
+
+  } // end namespace Vtk
+} // end namespace Dune

dune/vtk/vtkreader.hh

@@ -128,8 +128,8 @@ namespace Dune
         "The data to extract is not found in point-data. Try `getCellData()` instead!");
       VTK_ASSERT(data_it->second.section == POINT_DATA);
 
-      return {*creator_, point_it->second, data_it->second.components,
-              vec_types, vec_offsets, vec_connectivity};
+      return {*creator_, point_it->second, data_it->second.name, data_it->second.components,
+              data_it->second.type, vec_types, vec_offsets, vec_connectivity};
     }
 
     /// Return a vector of DataArrayAttributes for all POINT_DATA blocks
@@ -155,8 +155,8 @@ namespace Dune
         "The data to extract is not found in cell-data. Try `getPointData()` instead!");
       VTK_ASSERT(data_it->second.section == CELL_DATA);
 
-      return {*creator_, cell_it->second, data_it->second.components,
-              vec_types, vec_offsets, vec_connectivity};
+      return {*creator_, cell_it->second, data_it->second.name, data_it->second.components,
+              data_it->second.type, vec_types, vec_offsets, vec_connectivity};
     }
 
     /// Return a vector of DataArrayAttributes for all CELL_DATA blocks

dune/vtk/writers/vtkimagedatawriter.hh

@@ -8,6 +8,7 @@
 #include <dune/vtk/function.hh>
 #include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/structureddatacollector.hh>
+#include <dune/vtk/utility/concepts.hh>
 
 #include <dune/vtk/vtkwriterinterface.hh>
 
@@ -57,19 +58,19 @@ namespace Dune
   };
 
   // deduction guides
-  template <class GridView,
-    class = std::void_t<typename GridView::IndexSet>>
-  VtkImageDataWriter(GridView const&, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class GridView, class... Args,
+    Vtk::IsGridView<GridView> = true>
+  VtkImageDataWriter(GridView, Args...)
     -> VtkImageDataWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
 
-  template <class DataCollector,
-    class = std::void_t<typename DataCollector::GridView>>
-  VtkImageDataWriter(DataCollector&, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class DataCollector, class... Args,
+    Vtk::IsDataCollector<DataCollector> = true>
+  VtkImageDataWriter(DataCollector&, Args...)
     -> VtkImageDataWriter<typename DataCollector::GridView, DataCollector>;
 
-  template <class DataCollector,
-    class = std::void_t<typename DataCollector::GridView>>
-  VtkImageDataWriter(std::shared_ptr<DataCollector>, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class DataCollector, class... Args,
+    Vtk::IsDataCollector<DataCollector> = true>
+  VtkImageDataWriter(std::shared_ptr<DataCollector>, Args...)
     -> VtkImageDataWriter<typename DataCollector::GridView, DataCollector>;
 
 } // end namespace Dune

dune/vtk/writers/vtkrectilineargridwriter.hh

@@ -8,6 +8,7 @@
 #include <dune/vtk/function.hh>
 #include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/structureddatacollector.hh>
+#include <dune/vtk/utility/concepts.hh>
 
 #include <dune/vtk/vtkwriterinterface.hh>
 
@@ -63,19 +64,19 @@ namespace Dune
   };
 
   // deduction guides
-  template <class GridView,
-    class = std::void_t<typename GridView::IndexSet>>
-  VtkRectilinearGridWriter(GridView const&, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class GridView, class... Args,
+    Vtk::IsGridView<GridView> = true>
+  VtkRectilinearGridWriter(GridView, Args...)
     -> VtkRectilinearGridWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
 
-  template <class DataCollector,
-    class = std::void_t<typename DataCollector::GridView>>
-  VtkRectilinearGridWriter(DataCollector&, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class DataCollector, class... Args,
+    Vtk::IsDataCollector<DataCollector> = true>
+  VtkRectilinearGridWriter(DataCollector&, Args...)
     -> VtkRectilinearGridWriter<typename DataCollector::GridView, DataCollector>;
 
-  template <class DataCollector,
-    class = std::void_t<typename DataCollector::GridView>>
-  VtkRectilinearGridWriter(std::shared_ptr<DataCollector>, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class DataCollector, class... Args,
+    Vtk::IsDataCollector<DataCollector> = true>
+  VtkRectilinearGridWriter(std::shared_ptr<DataCollector>, Args...)
     -> VtkRectilinearGridWriter<typename DataCollector::GridView, DataCollector>;
 
 } // end namespace Dune

dune/vtk/writers/vtkstructuredgridwriter.hh

@@ -8,6 +8,7 @@
 #include <dune/vtk/function.hh>
 #include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/structureddatacollector.hh>
+#include <dune/vtk/utility/concepts.hh>
 
 #include <dune/vtk/vtkwriterinterface.hh>
 
@@ -57,19 +58,19 @@ namespace Dune
   };
 
   // deduction guides
-  template <class GridView,
-    class = std::void_t<typename GridView::IndexSet>>
-  VtkStructuredGridWriter(GridView const&, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class GridView, class... Args,
+    Vtk::IsGridView<GridView> = true>
+  VtkStructuredGridWriter(GridView, Args...)
     -> VtkStructuredGridWriter<GridView, Vtk::StructuredDataCollector<GridView>>;
 
-  template <class DataCollector,
-    class = std::void_t<typename DataCollector::GridView>>
-  VtkStructuredGridWriter(DataCollector&, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class DataCollector, class... Args,
+    Vtk::IsDataCollector<DataCollector> = true>
+  VtkStructuredGridWriter(DataCollector&, Args...)
     -> VtkStructuredGridWriter<typename DataCollector::GridView, DataCollector>;
 
-  template <class DataCollector,
-    class = std::void_t<typename DataCollector::GridView>>
-  VtkStructuredGridWriter(std::shared_ptr<DataCollector>, Vtk::FormatTypes = Vtk::FormatTypes::BINARY, Vtk::DataTypes = Vtk::DataTypes::FLOAT32)
+  template <class DataCollector, class... Args,
+    Vtk::IsDataCollector<DataCollector> = true>
+  VtkStructuredGridWriter(std::shared_ptr<DataCollector>, Args...)
     -> VtkStructuredGridWriter<typename DataCollector::GridView, DataCollector>;
 
 } // end namespace Dune

dune/vtk/writers/vtkunstructuredgridwriter.hh

@@ -8,6 +8,7 @@
 #include <dune/vtk/function.hh>
 #include <dune/vtk/types.hh>
 #include <dune/vtk/datacollectors/continuousdatacollector.hh>