dune-vtk 2.8
Loading...
Searching...
No Matches
function.hh
Go to the documentation of this file.
1#pragma once
2
3#include <numeric>
4#include <type_traits>
5
6#include <dune/common/typetraits.hh>
7#include <dune/common/version.hh>
8
10#include <dune/vtk/types.hh>
13
14namespace Dune
15{
16 // forward declarations
17 template <class T, int N>
19
20 template <class T, int N, int M>
22
23 namespace Vtk
24 {
26 template <class GridView>
28 {
29 using Element = typename GridView::template Codim<0>::Entity;
30 using LocalDomain = typename Element::Geometry::LocalCoordinate;
31
32 template <class F, class D>
33 using Range = std::decay_t<std::result_of_t<F(D)>>;
34
35 private:
36
37 template <class T, int N>
38 static auto sizeOfImpl (FieldVector<T,N>) -> std::integral_constant<int, N> { return {}; }
39
40 template <class T, int N, int M>
41 static auto sizeOfImpl (FieldMatrix<T,N,M>) -> std::integral_constant<int, N*M> { return {}; }
42
43 static auto sizeOfImpl (...) -> std::integral_constant<int, 1> { return {}; }
44
45 template <class T>
46 static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }
47
48 static std::vector<int> allComponents(int n)
49 {
50 std::vector<int> components(n);
51 std::iota(components.begin(), components.end(), 0);
52 return components;
53 }
54
55 public:
57
68 template <class LF, class... Args,
70 Function (LF&& localFct, std::string name, std::vector<int> components, Args&&... args)
71 : localFct_(std::forward<LF>(localFct))
72 , name_(std::move(name))
73 {
74 setComponents(std::move(components));
75 setRangeType(getArg<Vtk::RangeTypes>(args..., Vtk::RangeTypes::UNSPECIFIED), components_.size());
76 setDataType(getArg<Vtk::DataTypes>(args..., Vtk::DataTypes::FLOAT64));
77 }
78
80
90 template <class LF, class... Args,
92 Function (LF&& localFct, std::string name, int ncomps, Args&&... args)
93 : Function(std::forward<LF>(localFct), std::move(name), allComponents(ncomps),
94 std::forward<Args>(args)...)
95 {}
96
98
103 template <class LF, class... Args,
105 class R = Range<LF,LocalDomain> >
106 Function (LF&& localFct, std::string name, Args&&... args)
107 : Function(std::forward<LF>(localFct), std::move(name), sizeOf<R>(),
108 std::forward<Args>(args)...)
109 {}
110
112 template <class... Args>
113 explicit Function (Function<GridView> const& fct, Args&&... args)
114 : Function(fct.localFct_,
115 getArg<std::string, char const*>(args..., fct.name_),
116 getArg<int,unsigned int,long,unsigned long,std::vector<int>>(args..., fct.components_),
117 getArg<Vtk::RangeTypes>(args..., fct.rangeType_),
118 getArg<Vtk::DataTypes>(args..., fct.dataType_))
119 {}
120
122
130 template <class GF, class... Args,
131 disableCopyMove<Function, GF> = 0,
132 IsGridFunction<GF> = true>
133 Function (GF&& fct, std::string name, Args&&... args)
134 : Function(localFunction(std::forward<GF>(fct)), std::move(name), std::forward<Args>(args)...)
135 {}
136
138 template <class F>
139 Function (F&& fct, Vtk::FieldInfo info, ...)
140 : Function(std::forward<F>(fct), info.name(), info.size(), info.rangeType(), info.dataType())
141 {}
142
144 template <class F, class... Args,
145 disableCopyMove<Function, F> = 0,
146 IsGridFunction<F> = true,
147 class = decltype(std::declval<F>().name()),
148 class = decltype(std::declval<F>().numComponents()),
149 class = decltype(std::declval<F>().dataType()) >
150 explicit Function (F&& fct, ...)
151 : Function(localFunction(std::forward<F>(fct)), fct.name(), fct.numComponents(),
152 Vtk::RangeTypes::UNSPECIFIED, fct.dataType())
153 {}
154
156
159 explicit Function (std::shared_ptr<VTKFunction<GridView> const> const& fct, ...)
160 : localFct_(fct)
161 , name_(fct->name())
162 {
163 setComponents(fct->ncomps());
164 setDataType(dataTypeOf(fct->precision()));
165 setRangeType(rangeTypeOf(fct->ncomps()));
166 }
167
169 Function () = default;
170
173 {
174 return self.localFct_;
175 }
176
178 std::string const& name () const
179 {
180 return name_;
181 }
182
184 void setName (std::string name)
185 {
186 name_ = std::move(name);
187 }
188
190 int numComponents () const
191 {
192 return rangeType_ == Vtk::RangeTypes::SCALAR ? 1 :
193 rangeType_ == Vtk::RangeTypes::VECTOR ? 3 :
194 rangeType_ == Vtk::RangeTypes::TENSOR ? 9 : int(components_.size());
195 }
196
198 void setComponents (std::vector<int> components)
199 {
200 components_ = components;
201 localFct_.setComponents(components_);
202 }
203
205 void setComponents (int ncomps)
206 {
207 setComponents(allComponents(ncomps));
208 }
209
212 {
213 return dataType_;
214 }
215
218 {
219 dataType_ = type;
220 }
221
224 {
225 return rangeType_;
226 }
227
229 void setRangeType (Vtk::RangeTypes type, std::size_t ncomp = 1)
230 {
231 rangeType_ = type;
232 if (type == Vtk::RangeTypes::AUTO)
233 rangeType_ = rangeTypeOf(ncomp);
234 }
235
238 {
239 setName(info.name());
240 setComponents(info.size());
241 setRangeType(info.rangeType());
242 setDataType(info.dataType());
243 }
244
245 private:
247 std::string name_;
248 std::vector<int> components_;
251 };
252
253 } // end namespace Vtk
254} // end namespace Dune
Definition: writer.hh:13
decltype((std::declval< LF & >().bind(std::declval< LocalContext >()), std::declval< LF & >().unbind(), std::declval< LF >()(std::declval< typename LocalContext::Geometry::LocalCoordinate >()), true)) IsLocalFunction
Definition: concepts.hh:39
Vtk::DataTypes dataTypeOf()
Definition: types.hh:79
decltype((localFunction(std::declval< GF const & >()), true)) IsGridFunction
Definition: concepts.hh:32
RangeTypes
Type used to determine whether to limit output components to e.g. 3 (vector), or 9 (tensor)
Definition: types.hh:35
DataTypes
Definition: types.hh:52
decltype(auto) getArg(Arg0 &&arg0, Args &&... args)
Definition: arguments.hh:29
Vtk::RangeTypes rangeTypeOf(Dune::VTK::FieldInfo::Type t)
Definition: types.cc:56
Definition: function.hh:18
Definition: function.hh:21
Wrapper class for functions allowing local evaluations.
Definition: function.hh:28
Function(LF &&localFct, std::string name, int ncomps, Args &&... args)
(2) Construct from a LocalFunction directly
Definition: function.hh:92
Vtk::DataTypes dataType() const
Return the VTK Datatype associated with the functions range type.
Definition: function.hh:211
void setName(std::string name)
Set the function name.
Definition: function.hh:184
friend Vtk::LocalFunction< GridView > localFunction(Function const &self)
Create a LocalFunction.
Definition: function.hh:172
std::string const & name() const
Return a name associated with the function.
Definition: function.hh:178
Function(GF &&fct, std::string name, Args &&... args)
(5) Construct from a GridFunction
Definition: function.hh:133
Function()=default
(9) Default constructor. After construction, the function is an an invalid state.
void setRangeType(Vtk::RangeTypes type, std::size_t ncomp=1)
Set the category of the range, SCALAR, VECTOR, TENSOR, or UNSPECIFIED.
Definition: function.hh:229
int numComponents() const
Return the number of components of the Range as it is written to the file.
Definition: function.hh:190
void setDataType(Vtk::DataTypes type)
Set the data-type for the components.
Definition: function.hh:217
Function(F &&fct, Vtk::FieldInfo info,...)
(6) Constructor that forwards the number of components and data type to the other constructor
Definition: function.hh:139
void setComponents(int ncomps)
Set the number of components of the Range and generate component range [0...ncomps)
Definition: function.hh:205
Function(std::shared_ptr< VTKFunction< GridView > const > const &fct,...)
(8) Construct from legacy VTKFunction
Definition: function.hh:159
Vtk::RangeTypes rangeType() const
The category of the range, SCALAR, VECTOR, TENSOR, or UNSPECIFIED.
Definition: function.hh:223
Function(Function< GridView > const &fct, Args &&... args)
(4) Construct from a Vtk::Function
Definition: function.hh:113
void setComponents(std::vector< int > components)
Set the components of the Range to visualize.
Definition: function.hh:198
Function(LF &&localFct, std::string name, Args &&... args)
(3) Construct from a LocalFunction directly.
Definition: function.hh:106
Function(LF &&localFct, std::string name, std::vector< int > components, Args &&... args)
(1) Construct from a LocalFunction directly
Definition: function.hh:70
void setFieldInfo(Vtk::FieldInfo info)
Set all the parameters from a FieldInfo object.
Definition: function.hh:237
A Vtk::LocalFunction is a function-like object that can be bound to a grid element an that provides a...
Definition: localfunction.hh:23
Definition: types.hh:228
int size() const
The number of components in the data field.
Definition: types.hh:253
std::string const & name() const
The name of the data field.
Definition: types.hh:247
Vtk::RangeTypes rangeType() const
Return the category of the stored range.
Definition: types.hh:259
Vtk::DataTypes dataType() const
Return the data tpe of the data field.
Definition: types.hh:265