function.hh 8.78 KB
Newer Older
Stenger, Florian's avatar
Stenger, Florian committed
1
2
#pragma once

3
#include <numeric>
Stenger, Florian's avatar
Stenger, Florian committed
4
5
#include <type_traits>

6
#include <dune/common/typetraits.hh>
7
#include <dune/common/version.hh>
Stenger, Florian's avatar
Stenger, Florian committed
8

9
10
11
#include <dune/vtk/localfunction.hh>
#include <dune/vtk/types.hh>
#include <dune/vtk/utility/arguments.hh>
Stenger, Florian's avatar
Stenger, Florian committed
12
13
14

namespace Dune
{
15
  // forward declarations
Stenger, Florian's avatar
Stenger, Florian committed
16
17
18
19
20
21
22
23
24
25
26
27
  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
    {
28
29
30
      using Element = typename GridView::template Codim<0>::Entity;
      using LocalDomain = typename Element::Geometry::LocalCoordinate;

31
      template <class GF>
32
      using IsGridFunction = decltype(localFunction(std::declval<GF>()));
33
34

      template <class LocalFunction, class LF = std::decay_t<LocalFunction>>
35
      using IsLocalFunction = decltype((
36
37
38
        std::declval<LF&>().bind(std::declval<Element>()),
        std::declval<LF&>().unbind(),
        std::declval<LF>()(std::declval<LocalDomain>()),
39
      0));
40

41
42
      template <class F, class D>
      using Range = std::decay_t<std::result_of_t<F(D)>>;
Stenger, Florian's avatar
Stenger, Florian committed
43
44
45
46

    private:

      template <class T, int N>
47
      static auto sizeOfImpl (FieldVector<T,N>) -> std::integral_constant<int, N> { return {}; }
Stenger, Florian's avatar
Stenger, Florian committed
48
49

      template <class T, int N, int M>
50
      static auto sizeOfImpl (FieldMatrix<T,N,M>) -> std::integral_constant<int, N*M> { return {}; }
Stenger, Florian's avatar
Stenger, Florian committed
51

52
      static auto sizeOfImpl (...) -> std::integral_constant<int, 1> { return {}; }
Stenger, Florian's avatar
Stenger, Florian committed
53
54
55
56

      template <class T>
      static constexpr int sizeOf () { return decltype(sizeOfImpl(std::declval<T>()))::value; }

57
58
59
60
61
62
63
      static std::vector<int> allComponents(int n)
      {
        std::vector<int> components(n);
        std::iota(components.begin(), components.end(), 0);
        return components;
      }

Stenger, Florian's avatar
Stenger, Florian committed
64
    public:
65
      /// (1) Construct from a LocalFunction directly
Stenger, Florian's avatar
Stenger, Florian committed
66
      /**
67
68
69
      * \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
70
71
72
73
      * \param category    The \ref Vtk::RangeTypes category for the range. [Vtk::RangeTypes::AUTO]
      * \param dataType    The \ref Vtk::DataTypes used in the output. [Vtk::DataTypes::FLOAT32]
      *
      * The arguments `category` and `dataType` can be passed in any order.
74
75
      *
      * NOTE: Stores the localFunction by value.
Stenger, Florian's avatar
Stenger, Florian committed
76
      **/
77
      template <class LF, class... Args,
78
79
        class = IsLocalFunction<LF>>
      Function (LF&& localFct, std::string name, std::vector<int> components, Args&&... args)
80
        : localFct_(std::forward<LF>(localFct))
Stenger, Florian's avatar
Stenger, Florian committed
81
82
        , name_(std::move(name))
      {
83
        setComponents(std::move(components));
84
        setRangeType(getArg<Vtk::RangeTypes>(args..., Vtk::RangeTypes::UNSPECIFIED), components_.size());
85
        setDataType(getArg<Vtk::DataTypes>(args..., Vtk::DataTypes::FLOAT64));
Stenger, Florian's avatar
Stenger, Florian committed
86
87
      }

88
89
90
91
92
93
94
95
96
97
98
      /// (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.
      **/
99
100
101
102
103
104
105
106
107
108
109
110
111
      template <class LF, class... Args,
        class = IsLocalFunction<LF>>
      Function (LF&& localFct, std::string name, int ncomps, Args&&... args)
        : Function(std::forward<LF>(localFct), std::move(name), allComponents(ncomps),
                   std::forward<Args>(args)...)
      {}

      /// (3) Construct from a LocalFunction directly.
      /**
       * Same as Constructor (1) or (2) but deduces the number of components from
       * the static range type of the local-function. This defaults to 1 of no
       * static size information could be extracted.
       **/
112
      template <class LF, class... Args,
113
        class = IsLocalFunction<LF>,
114
        class R = Range<LF,LocalDomain> >
115
116
117
      Function (LF&& localFct, std::string name, Args&&... args)
        : Function(std::forward<LF>(localFct), std::move(name), sizeOf<R>(),
                   std::forward<Args>(args)...)
118
119
      {}

120
      /// (4) Construct from a Vtk::Function
121
      template <class... Args>
122
123
124
      explicit Function (Function<GridView> const& fct, Args&&... args)
        : Function(fct.localFct_,
                   getArg<std::string, char const*>(args..., fct.name_),
125
126
127
                   getArg<int,unsigned int,long,unsigned long,std::vector<int>>(args..., fct.components_),
                   getArg<Vtk::RangeTypes>(args..., fct.rangeType_),
                   getArg<Vtk::DataTypes>(args..., fct.dataType_))
128
129
      {}

130
      /// (5) Construct from a GridFunction
131
132
133
134
      /**
      * \param fct   A Grid(View)-function, providing a `localFunction(fct)`
      * \param name  The name to use as identification in the VTK file
      *
135
      * Forwards all other arguments to the constructor (1) or (2).
136
137
138
      *
      * NOTE: Stores the localFunction(fct) by value.
      */
139
140
141
142
143
      template <class GF, class... Args,
        disableCopyMove<Function, GF> = 0,
        class = IsGridFunction<GF> >
      Function (GF&& fct, std::string name, Args&&... args)
        : Function(localFunction(std::forward<GF>(fct)), std::move(name), std::forward<Args>(args)...)
144
145
      {}

146
      /// (6) Constructor that forwards the number of components and data type to the other constructor
147
      template <class F>
148
      Function (F&& fct, Vtk::FieldInfo info, ...)
149
        : Function(std::forward<F>(fct), info.name(), info.size(), info.rangeType(), info.dataType())
Stenger, Florian's avatar
Stenger, Florian committed
150
151
      {}

152
      /// (7) Construct from legacy VTKFunction
153
154
155
      /**
      * \param fct  The Dune::VTKFunction to wrap
      **/
156
      explicit Function (std::shared_ptr<VTKFunction<GridView> const> const& fct, ...)
157
158
159
        : localFct_(fct)
        , name_(fct->name())
      {
160
        setComponents(fct->ncomps());
161
#if DUNE_VERSION_LT(DUNE_GRID,2,7)
162
        setDataType(Vtk::DataTypes::FLOAT64);
163
#else
164
        setDataType(dataTypeOf(fct->precision()));
165
#endif
166
        setRangeType(rangeTypeOf(fct->ncomps()));
167
168
      }

169
      /// (8) Default constructor. After construction, the function is an an invalid state.
Stenger, Florian's avatar
Stenger, Florian committed
170
171
172
173
174
175
176
177
178
179
180
181
182
183
      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_;
      }

184
185
186
187
188
189
190
191
      /// Set the function name
      void setName (std::string name)
      {
        name_ = std::move(name);
      }

      /// Return the number of components of the Range as it is written to the file
      int numComponents () const
Stenger, Florian's avatar
Stenger, Florian committed
192
      {
193
194
195
        return rangeType_ == Vtk::RangeTypes::SCALAR ? 1 :
               rangeType_ == Vtk::RangeTypes::VECTOR ? 3 :
               rangeType_ == Vtk::RangeTypes::TENSOR ? 9 : int(components_.size());
Stenger, Florian's avatar
Stenger, Florian committed
196
197
      }

198
199
      /// Set the components of the Range to visualize
      void setComponents (std::vector<int> components)
200
      {
201
202
203
204
205
206
207
208
        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));
209
210
      }

Stenger, Florian's avatar
Stenger, Florian committed
211
      /// Return the VTK Datatype associated with the functions range type
212
213
214
215
216
217
218
      Vtk::DataTypes dataType () const
      {
        return dataType_;
      }

      /// Set the data-type for the components
      void setDataType (Vtk::DataTypes type)
Stenger, Florian's avatar
Stenger, Florian committed
219
      {
220
        dataType_ = type;
Stenger, Florian's avatar
Stenger, Florian committed
221
222
      }

223
224
225
226
227
228
229
      /// 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
230
      void setRangeType (Vtk::RangeTypes type, std::size_t ncomp = 1)
231
232
      {
        rangeType_ = type;
233
234
        if (type == Vtk::RangeTypes::AUTO)
          rangeType_ = rangeTypeOf(ncomp);
235
236
237
238
239
240
241
242
243
244
245
      }

      /// 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());
      }

Stenger, Florian's avatar
Stenger, Florian committed
246
247
248
    private:
      Vtk::LocalFunction<GridView> localFct_;
      std::string name_;
249
      std::vector<int> components_;
250
      Vtk::DataTypes dataType_ = Vtk::DataTypes::FLOAT64;
251
      Vtk::RangeTypes rangeType_ = Vtk::RangeTypes::UNSPECIFIED;
Stenger, Florian's avatar
Stenger, Florian committed
252
253
254
255
    };

  } // end namespace Vtk
} // end namespace Dune