VectorFacade.hpp 12.7 KB
Newer Older
Praetorius, Simon's avatar
Praetorius, Simon committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
#pragma once

#include <string>
#include <type_traits>
#include <utility>

#include <dune/common/classname.hh>
#include <dune/common/deprecated.hh>
#include <dune/common/hybridutilities.hh>
#include <dune/common/rangeutilities.hh>
#include <dune/common/shared_ptr.hh>
#include <dune/common/std/type_traits.hh>
#include <dune/functions/functionspacebases/concepts.hh>
#include <dune/functions/functionspacebases/sizeinfo.hh>

#include <amdis/Output.hpp>
#include <amdis/common/Concepts.hpp>
18
#include <amdis/common/ConceptsBase.hpp>
Praetorius, Simon's avatar
Praetorius, Simon committed
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
#include <amdis/common/FakeContainer.hpp>
#include <amdis/common/TypeTraits.hpp>
#include <amdis/functions/NodeIndices.hpp>
#include <amdis/operations/Assigner.hpp>
#include <amdis/typetree/MultiIndex.hpp>

namespace AMDiS
{
  /// \brief The basic container that stores a base vector and a corresponding basis
  /**
   * A vector storing all the assembled Operators indexed with DOF indices. The
   * vector data is associated to a global basis.
   *
   * \tparam T            The coefficient type of the vector
   * \tparam VectorImpl   A linear-algebra backend for the vector storage
   **/
  template <class T, template <class> class VectorImpl>
  class VectorFacade
  {
    using Self = VectorFacade;
    using Impl = VectorImpl<T>;

  private:
    // States the vector can be in. Is changed on various insertion or gathering methods.
    enum class VectorState
    {
      unknown = 0,
      synchronized = 1,
      insert_values = 2,
      add_values = 3
    };

    template <class A>
    using VectorStateOf_t = std::conditional_t<std::is_same<A,Assigner::plus_assign>::value,
      std::integral_constant<VectorState, VectorState::add_values>,
      std::integral_constant<VectorState, VectorState::insert_values>>;

  public:
    using size_type = typename Impl::size_type;
    using value_type = typename Impl::value_type;

60
    /// Constructor. Forwards the basis to the implementation class and
Praetorius, Simon's avatar
Praetorius, Simon committed
61
62
63
64
65
66
    /// constructs a (type-erased) size-info object.
    template <class GlobalBasis,
      class = void_t<decltype(std::declval<GlobalBasis>().dimension())> >
    VectorFacade(GlobalBasis const& basis)
      : impl_(basis)
    {
67
      resizeZero(sizeInfo(basis));
Praetorius, Simon's avatar
Praetorius, Simon committed
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
    }

    /// Return the underlying linear algebra backend
    Impl const& impl() const { return impl_; }
    Impl&       impl()       { return impl_; }


    template <class V>
    using HasLocalSize = decltype(std::declval<V>().localSize());

    template <class V>
    using HasGlobalSize = decltype(std::declval<V>().globalSize());

    /// Return the number of entries in the local part of the vector
    std::size_t localSize() const
    {
      return Dune::Hybrid::ifElse(Dune::Std::is_detected<HasLocalSize,Impl>{},
        [&](auto id) { return id(impl_).localSize(); },
        [&](auto id) { return id(impl_).size(); });
    }

    /// Return the number of entries in the global vector
    std::size_t globalSize() const
    {
      return Dune::Hybrid::ifElse(Dune::Std::is_detected<HasGlobalSize,Impl>{},
        [&](auto id) { return id(impl_).globalSize(); },
        [&](auto id) { return id(impl_).size(); });
    }

    /// Resize the \ref vector to the size of the \ref basis
98
99
    template <class SizeInfo>
    void resize(SizeInfo const& sizeInfo)
Praetorius, Simon's avatar
Praetorius, Simon committed
100
    {
101
      init(sizeInfo, false);
Praetorius, Simon's avatar
Praetorius, Simon committed
102
103
104
    }

    /// Resize the \ref vector to the size of the \ref basis and set to zero
105
106
    template <class SizeInfo>
    void resizeZero(SizeInfo const& sizeInfo)
Praetorius, Simon's avatar
Praetorius, Simon committed
107
    {
108
      init(sizeInfo, true);
Praetorius, Simon's avatar
Praetorius, Simon committed
109
110
111
112
    }

    /// Prepare the vector for insertion of values, finish the insertion with
    /// \ref finish().
113
114
    template <class SizeInfo>
    void init(SizeInfo const& sizeInfo, bool clear)
Praetorius, Simon's avatar
Praetorius, Simon committed
115
    {
116
      impl_.init(sizeInfo, clear);
Praetorius, Simon's avatar
Praetorius, Simon committed
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
      state_ = clear ? VectorState::synchronized : VectorState::unknown;
    }

    /// Finish the insertion of values, see \ref init()
    void finish()
    {
      impl_.finish();
      state_ = VectorState::unknown;
    }

    /// Return the value of the vector at the local index \ref idx
    template <class Index,
      REQUIRES(Concepts::MultiIndex<Index>)>
    typename Impl::value_type at(Index const& idx) const
    {
      const_cast<Self*>(this)->synchronize();
      return impl_.at(idx);
    }

    /// \brief Insert a single value into the matrix (add or overwrite to existing value)
    /**
     * Inserts or adds a value into a certain location \p dof (given as dof multi-index)
     * of a vector. The insertion mode is determined by the \p assign functor. Use
     * \ref Assigner::plus_assign for adding values (default) or \ref Assigner::assign
     * for overwriting (setting) values. Different insertion modes can not be mixed!
     *
     * Insertion must be closed with a call to \ref finish().
     *
     * [[not collective]]
     */
    template <class Index, class Assign = Assigner::plus_assign,
      REQUIRES(Concepts::MultiIndex<Index>)>
    void insert(Index const& idx, typename Impl::value_type const& value, Assign assign = {})
    {
      test_exit_dbg(state_ == VectorStateOf_t<Assign>::value ||
                    state_ == VectorState::unknown ||
                    state_ == VectorState::synchronized,
        "Vector in invalid state {} for insertion by {}.", to_string(state_), Dune::className<Assign>());

      impl_.insert(idx, value, assign);

      // set the state to insert_values or add_values
      state_ = VectorStateOf_t<Assign>::value;
    }

    /// See \ref insert for assignment operation \ref Assigner::assign
    template <class Index,
      REQUIRES(Concepts::MultiIndex<Index>)>
    void set(Index const& idx, typename Impl::value_type const& value)
    {
      insert(idx, value, Assigner::assign{});
    }

    /// See \ref insert for assignment operation \ref Assigner::plus_assign
    template <class Index,
      REQUIRES(Concepts::MultiIndex<Index>)>
    void add(Index const& idx, typename Impl::value_type const& value)
    {
      insert(idx, value, Assigner::plus_assign{});
    }


    /// \brief Extract values from the vector referring to the given local indices
    /// and store it into a buffer
    /**
     * Collect value of indices and store them into a buffer. The buffer must be
     * a vector-like container with `buffer.resize()` and `buffer.begin()`. The
     * indices must be stored in an iterable container.
     *
     * If the vector is not in synchronized state, collects all necessary values
     * possibly from neighbouring processors.
     *
     * [[expects: localView is bound to an element]]
     * [[expects: node is in localView.tree()]]
     * [[possibly neighbor-wise collective operation]]
     */
    template <class LocalView, class Node, class Buffer,
194
      REQUIRES(Concepts::LocalView<LocalView>),
195
      REQUIRES(Concepts::BasisNode<Node>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
196
197
198
199
200
201
202
203
204
205
206
207
208
209
    void gather(LocalView const& localView, Node const& node, Buffer& buffer) const
    {
      test_exit(state_ == VectorState::unknown ||
                state_ == VectorState::synchronized,
        "Vector in invalid state {} for gather operations.", to_string(state_));

      const_cast<Self*>(this)->synchronize();

      buffer.resize(node.size());
      impl_.gather(nodeIndices(localView, node), buffer.begin());
    }

    // [[expects: localView is bound to an element]]
    template <class LocalView, class Buffer,
210
      REQUIRES(Concepts::LocalView<LocalView>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
    void gather(LocalView const& localView, Buffer& buffer) const
    {
      gather(localView, localView.tree(), buffer);
    }

    /// Insert a block of values into the vector (add or overwrite to existing values)
    /**
     * Inserts or adds values into certain locations of a vector. Insertion indices
     * are extracted from the given \p localView. The insertion mode is determined
     * by the \p assign functor. Use \ref Assigner::plus_assign for adding values
     * (default) or \ref Assigner::assign for overwriting (setting) values. Different
     * insertion modes can not be mixed! The \p localVector is assumed to be a continuous
     * memory container with a `data()` method to get a pointer to the beginning.
     *
     * The \p mask models a boolean range with at least a `begin()` method. Must
     * be forward iterable for at least `localVector.size()` elements. Does not
     * need an `end()` method. See, e.g. \ref FakeContainer.
     *
     * Insertion must be closed with a call to \ref finish(). It is not allowed
     * to switch insertion mode before calling `finish()`.
     *
     * [[expects: localView is bound to an element]]
     * [[expects: node is in localView.tree()]]
     * [[not collective]]
     */
    template <class LocalView, class Node, class NodeVector, class MaskRange, class Assign,
237
      REQUIRES(Concepts::LocalView<LocalView>),
238
      REQUIRES(Concepts::BasisNode<Node>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
    void scatter(LocalView const& localView, Node const& node, NodeVector const& localVector,
                 MaskRange const& mask, Assign assign)
    {
      test_exit(state_ == VectorStateOf_t<Assign>::value ||
                state_ == VectorState::unknown ||
                state_ == VectorState::synchronized,
        "Vector in invalid state {} for insertion by {}.", to_string(state_), Dune::className<Assign>());

      assert(localVector.size() == node.size());

      // create a range of DOF indices on the node
      impl_.scatter(nodeIndices(localView, node), localVector, mask, assign);

      // set the state to insert_values or add_values
      state_ = VectorStateOf_t<Assign>::value;
    }

    /// Call \ref scatter with `MaskRange` set to \ref FakeContainer.
    // [[expects: localView is bound to an element]]
    // [[expects: node is in localView.tree()]]
    template <class LocalView, class Node, class NodeVector, class Assign,
260
      REQUIRES(Concepts::LocalView<LocalView>),
261
      REQUIRES(Concepts::BasisNode<Node>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
262
263
264
265
266
267
268
269
    void scatter(LocalView const& localView, Node const& node, NodeVector const& localVector, Assign assign)
    {
      scatter(localView, node, localVector, FakeContainer<bool,true>{}, assign);
    }

    /// Call \ref scatter with `Node` given by the tree of the \p localView.
    // [[expects: localView is bound to an element]]
    template <class LocalView, class LocalVector, class Assign,
270
      REQUIRES(Concepts::LocalView<LocalView>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
    void scatter(LocalView const& localView, LocalVector const& localVector, Assign assign)
    {
      scatter(localView, localView.tree(), localVector, assign);
    }

    /// Apply \p func to each value at given indices \p localInd
    /**
     * First, synchronizes the values of the vector, then applies the functor to
     * each local value associated to the local indices \p localInd.
     *
     * [[mutable]]
     **/
    template <class LocalInd, class Func>
    void forEach(LocalInd const& localInd, Func&& func)
    {
      synchronize();
      impl_.forEach(localInd, FWD(func));
    }

    /// Call \ref forEach for all entries in the vector.
    template <class Func>
    void forEach(Func&& func)
    {
      synchronize();
      impl_.forEach(FWD(func));
    }

    /// Apply \p func to each value at given indices \p localInd
    /**
     * First, synchronizes the values of the vector, then applies the functor to
     * each local value associated to the local indices \p localInd.
     *
     * [[const]]
     **/
    template <class LocalInd, class Func>
    void forEach(LocalInd const& localInd, Func&& func) const
    {
      const_cast<Self*>(this)->synchronize();
      impl_.forEach(localInd, FWD(func));
    }

    /// Call \ref forEach for all entries in the vector.
    template <class Func>
    void forEach(Func&& func) const
    {
      const_cast<Self*>(this)->synchronize();
      impl_.forEach(FWD(func));
    }


  private:
    // synchronizes ghost values. Does not touch owned values
    void synchronize()
    {
      if (state_ != VectorState::synchronized)
        impl_.synchronize();

      state_ = VectorState::synchronized;
    }

    // print the vector state to string for debugging purposes
    static std::string to_string(VectorState state)
    {
      switch (state) {
        case VectorState::synchronized: return "synchronized";
        case VectorState::insert_values: return "insert_values";
        case VectorState::add_values: return "add_values";
        default: return "unknown";
      }
    }

  private:
    /// Data backend
    Impl impl_;

    /// The current state of the vector, one of {synchronized, insert_values,
    /// add_values, unknown}
    VectorState state_ = VectorState::unknown;
  };


  namespace Impl
  {
    template <class T, class Facade>
    struct VectorTypeImpl;

    template <class T, class S, template <class> class Impl>
    struct VectorTypeImpl<T,VectorFacade<S,Impl>>
    {
      using type = VectorFacade<T,Impl>;
    };

  } // end namespace Impl

  /// Type transformation changing the value type of the vector
  template <class T, class Facade>
  using VectorType_t = typename Impl::VectorTypeImpl<T,Facade>::type;

} // end namespace AMDiS