Commit acbaf624 authored by Praetorius, Simon's avatar Praetorius, Simon
Browse files

Add constructor with a components sub-set parameter and add a test for Vtk::Function

add AUTO RangeType

Add setFieldInfo member
parent b8b41e2a
#pragma once
#include <numeric>
#include <type_traits>
#include <dune/common/typetraits.hh>
......@@ -37,73 +38,88 @@ namespace Dune
private:
template <class T, int N>
static auto sizeOfImpl (FieldVector<T,N> const&)
-> std::integral_constant<int, N> { return {}; }
static auto sizeOfImpl (FieldVector<T,N>) -> 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 (FieldMatrix<T,N,M>) -> std::integral_constant<int, N*M> { return {}; }
static auto sizeOfImpl (...)
-> std::integral_constant<int, 1> { return {}; }
static auto sizeOfImpl (...) -> std::integral_constant<int, 1> { return {}; }
template <class T>
static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }
static std::vector<int> allComponents(int n)
{
std::vector<int> components(n);
std::iota(components.begin(), components.end(), 0);
return components;
}
public:
/// Construct from a LocalFunction directly
/// (1) Construct from a LocalFunction directly
/**
* \param localFct A local-function, providing a `bind(Element)` and an `operator()(LocalDomain)`
* \param name The name to use as 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.
* \param localFct A local-function, providing a `bind(Element)` and an `operator()(LocalDomain)`
* \param name The name to use as identification in the VTK file
* \param components A vector of component indices to extract from the range type
* \param category The type category for the range, i.e., vector, tensor or unspecified [AUTO]
* \param dataType The \ref Vtk::DataTypes used in the output. [FLOAT32]
*
* NOTE: Stores the localFunction by value.
**/
template <class LF,
class = void_t<HasBind<LF,Element>>,
class R = Range<LF,LocalDomain> >
Function (LF&& localFct, std::string name, int ncomps = sizeOf<R>(),
Vtk::DataTypes type = Vtk::dataTypeOf<R>())
Function (LF&& localFct, std::string name, std::vector<int> components,
Vtk::RangeTypes rangeType = Vtk::RangeTypes::AUTO,
Vtk::DataTypes dataType = Vtk::DataTypes::FLOAT32)
: localFct_(std::forward<LF>(localFct))
, name_(std::move(name))
{
setNumComponents(ncomps);
setDataType(type);
setComponents(std::move(components));
setRangeType(rangeType == Vtk::RangeTypes::AUTO ? rangeTypeOf(components.size()) : rangeType);
setDataType(dataType);
}
/// Construct from a GridFunction
/// (2) Construct from a LocalFunction directly
/**
* \param localFct A local-function, providing a `bind(Element)` and an `operator()(LocalDomain)`
* \param name The name to use as 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.
*
* Forwards all the other parmeters to the constructor (1)
*
* NOTE: Stores the localFunction by value.
**/
template <class LF, class... Args,
class = void_t<HasBind<LF,Element>>,
class R = Range<LF,LocalDomain> >
Function (LF&& localFct, std::string name, int ncomps = sizeOf<R>(), Args&&... args)
: Function(std::forward<LF>(localFct), std::move(name), allComponents(ncomps), std::forward<Args>(args)...)
{}
/// (3) Construct from a GridFunction
/**
* \param fct A Grid(View)-function, providing a `localFunction(fct)`
* \param name The name to use as identification in the VTK file
*
* Forwards all other arguments to the local-function constructor.
* Forwards all other arguments to the constructor (1) or (2).
*
* NOTE: Stores the localFunction(fct) by value.
*/
template <class F, class... Args,
disableCopyMove<Function,F> = 0,
class = void_t<LocalFunction<F>> >
Function (F&& fct, std::string name, Args&&... args)
: Function(localFunction(std::forward<F>(fct)), std::move(name), std::forward<Args>(args)...)
{}
/// Constructor that forwards the number of components and data type to the other constructor
/// (4) Constructor that forwards the number of components and data type to the other constructor
template <class F>
Function (F&& fct, Vtk::FieldInfo info)
: Function(std::forward<F>(fct), info.name(), info.numComponents(), info.dataType())
{}
/// Constructor that forwards the number of components and data type to the other constructor
template <class F>
Function (F&& fct, Dune::VTK::FieldInfo info)
: Function(std::forward<F>(fct), info.name(), info.size(), dataTypeOf(info.precision()))
: Function(std::forward<F>(fct), info.name(), info.size(), info.rangeType(), info.dataType())
{}
/// Construct from legacy VTKFunction
/// (5) Construct from legacy VTKFunction
/**
* \param fct The Dune::VTKFunction to wrap
**/
......@@ -111,10 +127,12 @@ namespace Dune
: localFct_(fct)
, name_(fct->name())
{
setNumComponents(fct->ncomps());
setComponents(fct->ncomps());
setDataType(dataTypeOf(fct->precision()));
setRangeType(rangeTypeOf(fct->ncomps()));
}
/// (7) Default constructor. After construction, the function is an an invalid state.
Function () = default;
/// Create a LocalFunction
......@@ -138,14 +156,22 @@ namespace Dune
/// Return the number of components of the Range as it is written to the file
int numComponents () const
{
return ncomps_ > 3 ? 9 : ncomps_ > 1 ? 3 : 1; // tensor, vector, scalar
return rangeType_ == Vtk::RangeTypes::SCALAR ? 1 :
rangeType_ == Vtk::RangeTypes::VECTOR ? 3 :
rangeType_ == Vtk::RangeTypes::TENSOR ? 9 : int(components_.size());
}
/// Set the number of components of the Range
void setNumComponents (int ncomps)
/// Set the components of the Range to visualize
void setComponents (std::vector<int> components)
{
ncomps_ = ncomps;
localFct_.setNumComponents(ncomps_);
components_ = components;
localFct_.setComponents(components_);
}
/// Set the number of components of the Range and generate component range [0...ncomps)
void setComponents (int ncomps)
{
setComponents(allComponents(ncomps));
}
/// Return the VTK Datatype associated with the functions range type
......@@ -160,11 +186,33 @@ namespace Dune
dataType_ = type;
}
/// The category of the range, SCALAR, VECTOR, TENSOR, or UNSPECIFIED
Vtk::RangeTypes rangeType () const
{
return rangeType_;
}
/// Set the category of the range, SCALAR, VECTOR, TENSOR, or UNSPECIFIED
void setRangeType (Vtk::RangeTypes type)
{
rangeType_ = type;
}
/// Set all the parameters from a FieldInfo object
void setFieldInfo (Vtk::FieldInfo info)
{
setName(info.name());
setComponents(info.size());
setRangeType(info.rangeType());
setDataType(info.dataType());
}
private:
Vtk::LocalFunction<GridView> localFct_;
std::string name_;
int ncomps_ = 1;
std::vector<int> components_;
Vtk::DataTypes dataType_ = Vtk::DataTypes::FLOAT32;
Vtk::RangeTypes rangeType_ = Vtk::RangeTypes::UNSPECIFIED;
};
} // end namespace Vtk
......
......@@ -3,7 +3,7 @@
#include <memory>
#include <type_traits>
#include <dune/common/std/type_traits.hh>
#include <dune/common/typetraits.hh>
#include "localfunctioninterface.hh"
#include "legacyvtkfunction.hh"
......@@ -31,28 +31,30 @@ namespace Dune
private:
struct RangeProxy
{
using value_type = double;
using field_type = double;
RangeProxy (LocalFunctionInterface<GridView> const& localFct,
std::size_t size, LocalCoordinate const& local)
std::vector<int> const& components,
LocalCoordinate const& local)
: localFct_(localFct)
, size_(size)
, components_(components)
, local_(local)
{}
std::size_t size () const
{
return size_;
return components_.size();
}
double operator[] (std::size_t i) const
{
return i < size_ ? localFct_.evaluate(i, local_) : 0.0;
return i < size() ? localFct_.evaluate(components_[i], local_) : 0.0;
}
private:
LocalFunctionInterface<GridView> const& localFct_;
std::size_t size_;
std::vector<int> const& components_;
LocalCoordinate local_;
};
......@@ -61,18 +63,17 @@ namespace Dune
template <class LF,
disableCopyMove<Self, LF> = 0,
class = void_t<HasBind<LF,Entity>> >
explicit LocalFunction (LF&& lf, int ncomps = 1)
explicit LocalFunction (LF&& lf)
: localFct_(std::make_shared<LocalFunctionWrapper<GridView,LF>>(std::forward<LF>(lf)))
, ncomps_(ncomps)
{}
/// Construct a Vtk::LocalFunction from a legacy VTKFunction
explicit LocalFunction (std::shared_ptr<VTKFunction<GridView> const> const& lf, int ncomps = 1)
explicit LocalFunction (std::shared_ptr<VTKFunction<GridView> const> const& lf)
: localFct_(std::make_shared<VTKLocalFunctionWrapper<GridView>>(lf))
, ncomps_(ncomps)
{}
/// Allow the default construction of a Vtk::LocalFunction
/// Allow the default construction of a Vtk::LocalFunction. After construction, the
/// LocalFunction is in an invalid state.
LocalFunction () = default;
/// Bind the function to the grid entity
......@@ -93,24 +94,24 @@ namespace Dune
RangeProxy operator() (LocalCoordinate const& xi) const
{
assert(bool(localFct_));
return {*localFct_, ncomps_, xi};
return {*localFct_, components_, xi};
}
/// Evaluate the `comp` component of the Range value at local coordinate `xi`
double evaluate (int comp, LocalCoordinate const& xi) const
/// Evaluate the `c`th component of the Range value at local coordinate `xi`
double evaluate (int c, LocalCoordinate const& xi) const
{
assert(bool(localFct_));
return comp < ncomps_ ? localFct_->evaluate(comp, xi) : 0.0;
return c < components_.size() ? localFct_->evaluate(components_[c], xi) : 0.0;
}
void setNumComponents (int ncomps)
void setComponents (std::vector<int> components)
{
ncomps_ = ncomps;
components_ = std::move(components);
}
private:
std::shared_ptr<LocalFunctionInterface<GridView>> localFct_ = nullptr;
int ncomps_ = 1;
std::vector<int> components_;
};
} // end namespace Vtk
......
dune_add_test(SOURCES test-map-datatypes.cc
LINK_LIBRARIES dunevtk)
dune_add_test(SOURCES test-function.cc
LINK_LIBRARIES dunevtk)
dune_add_test(SOURCES test-typededuction.cc
LINK_LIBRARIES dunevtk)
......
#include <config.h>
#include <optional>
#include <dune/grid/io/file/vtk/common.hh>
#include <dune/grid/utility/structuredgridfactory.hh>
#include <dune/vtk/function.hh>
#include <dune/vtk/vtkwriter.hh>
#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
// 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);
auto grid = StructuredGridFactory<GridType>::createCubeGrid({0.0,0.0}, {1.0,2.0}, {2u,4u});
using GridView = typename GridType::LeafGridView;
GridView gridView = grid->leafGridView();
VtkWriter writer{gridView, Vtk::FormatTypes::ASCII};
// 1. add (legacy) VTKFunction
std::vector<double> p1function(gridView.size(2), 1.0);
using P1Function = P1VTKFunction<GridView,std::vector<double>>;
std::shared_ptr<VTKFunction<GridView> const> f1(new P1Function(gridView, p1function, "p1"));
writer.addPointData(f1);
// 2. Add GlobalFunction
auto f2 = GlobalFunction{gridView, [](auto x) { return x[0]+x[1]; }};
writer.addPointData(f2, "global");
// 3. Add data + Vtk::FieldInfo description
auto f3 = GlobalFunction{gridView, [](auto x) { return x; }};
writer.addPointData(f3, Dune::Vtk::FieldInfo{"vector1", 3, Vtk::RangeTypes::VECTOR});
// 4. Add data + Dune::VTK::FieldInfo description
writer.addPointData(f3, Dune::VTK::FieldInfo{"vector2", Dune::VTK::FieldInfo::Type::vector, 3, Dune::VTK::Precision::float64});
// 5. Wrap function explicitly into Vtk::Function
writer.addPointData(Vtk::Function<GridView>{f3, "vector3", 3}); // calls copy-constructor of Vtk::Function
writer.addPointData(Vtk::Function<GridView>{f3, "vector3", 3}, "vector4", 3); // calls grid-function constructor
writer.addPointData(Vtk::Function<GridView>{f3, "vector5", {0,1}});
writer.addPointData(Vtk::Function<GridView>{f3, "vector5"}, "vector6", std::vector{0}, Vtk::RangeTypes::VECTOR);
// test default constructible
Vtk::Function<GridView> func0;
// test constructible by (legacy) VTKFunction
Vtk::Function<GridView> func1{f1};
writer.addPointData(func1);
// test constructible by local-function
Vtk::Function<GridView> func2{f2, "f2"};
writer.addPointData(func2);
Vtk::Function<GridView> func3{f3, "f3", 3};
writer.addPointData(func3);
// test constructible with component vector
Vtk::Function<GridView> func4a{f3, "f4a", 1};
Vtk::Function<GridView> func4b{f3, "f4b", {1}}; // == func4a
Vtk::Function<GridView> func4c{f3, "f4c", std::vector{1}};
writer.addPointData(func4a);
writer.addPointData(func4b);
writer.addPointData(func4c);
// Test copy-constructible
auto func5{func3};
auto func6 = func3;
// Test move-constructible
auto func7{std::move(func5)};
writer.addPointData(func7, "f7");
auto func8 = std::move(func6);
writer.addPointData(func8, "f8");
writer.write("test-function.vtu");
}
\ No newline at end of file
......@@ -8,48 +8,70 @@
namespace Dune {
namespace Vtk {
std::string to_string (FormatTypes type)
std::string to_string (Vtk::FormatTypes type)
{
switch (type) {
case ASCII: return "ascii";
case BINARY: return "binary";
case COMPRESSED: return "compressed";
case APPENDED: return "appended";
case Vtk::FormatTypes::ASCII: return "ascii";
case Vtk::FormatTypes::BINARY: return "binary";
case Vtk::FormatTypes::COMPRESSED: return "compressed";
case Vtk::FormatTypes::APPENDED: return "appended";
default:
DUNE_THROW(RangeError, "FormatType not found.");
std::abort();
}
}
FormatTypes formatTypeOf (Dune::VTK::OutputType o)
Vtk::FormatTypes formatTypeOf (Dune::VTK::OutputType o)
{
switch (o) {
case Dune::VTK::ascii: return Vtk::ASCII;
// case Dune::VTK::base64: return Vtk::BASE64;
case Dune::VTK::appendedraw: return Vtk::BINARY;
// case Dune::VTK::appendedbase64: return Vtk::APPENDED_BASE64;
// case Dune::VTK::binarycompressed: return Vtk::COMPRESSED;
// case Dune::VTK::compressedappended: return Vtk::APPENDED_COMPRESSED;
case Dune::VTK::ascii: return Vtk::FormatTypes::ASCII;
// case Dune::VTK::base64: return Vtk::FormatTypes::BASE64;
case Dune::VTK::appendedraw: return Vtk::FormatTypes::BINARY;
// case Dune::VTK::appendedbase64: return Vtk::FormatTypes::BASE64;
// case Dune::VTK::binarycompressed: return Vtk::FormatTypes::COMPRESSED;
// case Dune::VTK::compressedappended: return Vtk::FormatTypes::COMPRESSED;
default:
DUNE_THROW(RangeError, "OutputType not supported.");
std::abort();
}
}
// Map a dune-grid FieldInfo::Type to ValueTypes
Vtk::RangeTypes rangeTypeOf (Dune::VTK::FieldInfo::Type t)
{
switch (t) {
case Dune::VTK::FieldInfo::Type::scalar: return Vtk::RangeTypes::UNSPECIFIED;
case Dune::VTK::FieldInfo::Type::vector: return Vtk::RangeTypes::VECTOR;
case Dune::VTK::FieldInfo::Type::tensor: return Vtk::RangeTypes::TENSOR;
default:
DUNE_THROW(RangeError, "FieldInfo::Type not supported.");
std::abort();
}
}
// Map a number of components to a corresponding value type
Vtk::RangeTypes rangeTypeOf (int ncomps)
{
return ncomps > 9 ? Vtk::RangeTypes::UNSPECIFIED :
ncomps > 3 ? Vtk::RangeTypes::TENSOR :
ncomps > 1 ? Vtk::RangeTypes::VECTOR :
Vtk::RangeTypes::SCALAR;
}
std::string to_string (DataTypes type)
std::string to_string (Vtk::DataTypes type)
{
switch (type) {
case INT8: return "Int8";
case UINT8: return "UInt8";
case INT16: return "Int16";
case UINT16: return "UInt16";
case INT32: return "Int32";
case UINT32: return "UInt32";
case INT64: return "Int64";
case UINT64: return "UInt64";
case FLOAT32: return "Float32";
case FLOAT64: return "Float64";
case Vtk::DataTypes::INT8: return "Int8";
case Vtk::DataTypes::UINT8: return "UInt8";
case Vtk::DataTypes::INT16: return "Int16";
case Vtk::DataTypes::UINT16: return "UInt16";
case Vtk::DataTypes::INT32: return "Int32";
case Vtk::DataTypes::UINT32: return "UInt32";
case Vtk::DataTypes::INT64: return "Int64";
case Vtk::DataTypes::UINT64: return "UInt64";
case Vtk::DataTypes::FLOAT32: return "Float32";
case Vtk::DataTypes::FLOAT64: return "Float64";
default:
DUNE_THROW(RangeError, "DataType not found.");
std::abort();
......@@ -59,11 +81,11 @@ std::string to_string (DataTypes type)
Vtk::DataTypes dataTypeOf (Dune::VTK::Precision p)
{
switch (p) {
case Dune::VTK::Precision::int32: return Vtk::INT32;
case Dune::VTK::Precision::uint8: return Vtk::UINT8;
case Dune::VTK::Precision::uint32: return Vtk::UINT32;
case Dune::VTK::Precision::float32: return Vtk::FLOAT32;
case Dune::VTK::Precision::float64: return Vtk::FLOAT64;
case Dune::VTK::Precision::int32: return Vtk::DataTypes::INT32;
case Dune::VTK::Precision::uint8: return Vtk::DataTypes::UINT8;
case Dune::VTK::Precision::uint32: return Vtk::DataTypes::UINT32;
case Dune::VTK::Precision::float32: return Vtk::DataTypes::FLOAT32;
case Dune::VTK::Precision::float64: return Vtk::DataTypes::FLOAT64;
default:
DUNE_THROW(RangeError, "Precision not supported.");
std::abort();
......@@ -73,19 +95,19 @@ Vtk::DataTypes dataTypeOf (Dune::VTK::Precision p)
Vtk::DataTypes dataTypeOf (std::string s)
{
static const std::map<std::string, Vtk::DataTypes> to_datatype{
{"Int8", INT8},
{"UInt8", UINT8},
{"Int16", INT16},
{"UInt16", UINT16},
{"Int32", INT32},
{"UInt32", UINT32},
{"Int64", INT64},
{"UInt64", UINT64},
{"Float32", FLOAT32},
{"Float64", FLOAT64}
{"Int8", Vtk::DataTypes::INT8},
{"UInt8", Vtk::DataTypes::UINT8},
{"Int16", Vtk::DataTypes::INT16},
{"UInt16", Vtk::DataTypes::UINT16},
{"Int32", Vtk::DataTypes::INT32},
{"UInt32", Vtk::DataTypes::UINT32},
{"Int64", Vtk::DataTypes::INT64},
{"UInt64", Vtk::DataTypes::UINT64},
{"Float32", Vtk::DataTypes::FLOAT32},
{"Float64", Vtk::DataTypes::FLOAT64}
};
auto it = to_datatype.find(s);
return it != to_datatype.end() ? it->second : Vtk::UNKNOWN;
return it != to_datatype.end() ? it->second : Vtk::DataTypes::UNKNOWN;
}
......
......@@ -15,16 +15,33 @@ namespace Dune
{
namespace Vtk
{
/// Type used for representing the output format
enum FormatTypes {
ASCII = 1<<0,
BINARY = 1<<1,
COMPRESSED = 1<<2,
APPENDED = BINARY | COMPRESSED
};
std::string to_string (FormatTypes);
std::string to_string (Vtk::FormatTypes);
/// Map the dune-grid OutputType to FormatTypes
FormatTypes formatTypeOf(Dune::VTK::OutputType);
Vtk::FormatTypes formatTypeOf(Dune::VTK::OutputType);
/// Type used to determine whether to limit output components to e.g. 3 (vector), or 9 (tensor)
enum class RangeTypes {
UNSPECIFIED, //< The output components are not restricted
AUTO, //< Detect the category automatically from number of components
SCALAR, //< Use exactly 1 component
VECTOR, //< Use exactly 3 components
TENSOR //< Use exactly 9 components
};
// Map a dune-grid FieldInfo::Type to ValueTypes
Vtk::RangeTypes rangeTypeOf (Dune::VTK::FieldInfo::Type);
// Map a number of components to a corresponding value type
Vtk::RangeTypes rangeTypeOf (int ncomps);
enum DataTypes {
......@@ -36,7 +53,7 @@ namespace Dune
FLOAT32 = 32,
FLOAT64 = 64
};