AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
GlobalIdSet.hpp
1#pragma once
2
3#include <tuple>
4#include <utility>
5#include <vector>
6
7#include <dune/common/exceptions.hh>
8#include <dune/common/hybridutilities.hh>
9#include <dune/common/version.hh>
10#include <dune/functions/functionspacebases/nodes.hh>
11#include <dune/functions/functionspacebases/lagrangebasis.hh>
12#include <dune/functions/functionspacebases/lagrangedgbasis.hh>
13#include <dune/grid/common/gridenums.hh>
14#include <dune/localfunctions/common/localkey.hh>
15
16#include <amdis/Output.hpp>
17#include <amdis/common/Apply.hpp>
18#include <amdis/common/ConceptsBase.hpp>
19#include <amdis/common/ForEach.hpp>
20#include <amdis/common/TupleUtility.hpp>
21#include <amdis/utility/Twist.hpp>
22
23
24namespace AMDiS
25{
26 // forward declaration
27 template <class PreBasis, class TP, class NodeTag = typename PreBasis::Node::NodeTag>
28 class NodeIdSet;
29
31
35 template <class EntityIdType>
36 class GlobalIdType : public std::pair<EntityIdType, std::size_t>
37 {
38 using Super = std::pair<EntityIdType, std::size_t>;
39
40 public:
41 GlobalIdType(std::size_t i = 0)
42 : Super(EntityIdType{}, i)
43 {}
44
45 using Super::Super;
46
49 {
50 ++this->second;
51 return *this;
52 }
53
54 // Add to the second index only
55 GlobalIdType operator+(std::size_t shift) const
56 {
57 return GlobalIdType{this->first, this->second + shift};
58 }
59
60 friend std::ostream& operator<<(std::ostream& os, GlobalIdType const& id)
61 {
62 os << "(" << id.first << "," << id.second << ")";
63 return os;
64 }
65 };
66
67 template <class Basis>
68 using GlobalIdType_t
69 = GlobalIdType<typename Basis::GridView::Grid::GlobalIdSet::IdType>;
70
71
73
93 template <class GB, class = std::void_t<> >
95 {
96 public:
97 using GlobalBasis = GB;
98 using Tree = typename GB::LocalView::Tree;
99 using size_type = std::size_t;
100
101 using GridView = typename GlobalBasis::GridView;
102 using Grid = typename GridView::Grid;
103 using Element = typename GridView::template Codim<0>::Entity;
104 using EntityIdType = typename Grid::GlobalIdSet::IdType;
105 using PartitionType = Dune::PartitionType;
107
108 using PreBasis = typename GlobalBasis::PreBasis;
109 using TreePath = typename GlobalBasis::PrefixPath;
112
113 public:
114 GlobalBasisIdSet(GlobalBasis const& globalBasis)
115 : tree_(globalBasis.localView().tree())
116 , nodeIdSet_(globalBasis.gridView())
117 , twist_(globalBasis.gridView().grid().globalIdSet())
118 {
119 initializeTree(tree_);
120 }
121
123
128 void bind(Element const& element)
129 {
130 bindTree(tree_, element);
131 nodeIdSet_.bind(tree_);
132 twist_.bind(element);
133 data_.resize(size());
134 nodeIdSet_.fillIn(twist_, data_.begin());
135 }
136
138 void unbind()
139 {
140 nodeIdSet_.unbind();
141 }
142
144 size_type size() const
145 {
146 return tree_.size();
147 }
148
151 IdType id(size_type i) const
152 {
153 assert( i < data_.size() );
154 return data_[i].first;
155 }
156
159 EntityIdType entityId(size_type i) const
160 {
161 return id(i).first;
162 }
163
166 PartitionType partitionType(size_type i) const
167 {
168 assert( i < data_.size() );
169 return data_[i].second;
170 }
171
172 protected:
173 Tree tree_;
174 NodeIdSet nodeIdSet_;
175 Twist twist_;
176 using Data = std::pair<IdType, PartitionType>;
177 std::vector<Data> data_;
178 };
179
180
181 // Specialization for SubspaceBasis
182 template <class Basis>
183 class GlobalBasisIdSet<Basis, std::void_t<typename Basis::RootBasis>>
184 : public GlobalBasisIdSet<typename Basis::RootBasis>
185 {
186 public:
187 GlobalBasisIdSet(Basis const& basis)
189 {}
190 };
191
192
193 template <class PB, class TP, class NodeTag>
195 {
196 public:
197 using PreBasis = PB;
198 using Node = typename PreBasis::Node;
199 using GridView = typename PreBasis::GridView;
200 using size_type = std::size_t;
201
202 static_assert(Node::isLeaf, "Generic NodeIdSet implemented for LeafNodes only. Provide a specialization for your node!");
203
204 private:
205 static constexpr int dim = GridView::template Codim<0>::Entity::mydimension;
206
207 public:
208 NodeIdSet(GridView const& gridView)
209 : gridView_(gridView)
210 , sizes_{}
211 {}
212
214 void bind(const Node& node)
215 {
216 node_ = &node;
217 size_ = node_->finiteElement().size();
218
219 std::fill(std::begin(sizes_), std::end(sizes_), 0u);
220 for (size_type i = 0; i < size_ ; ++i) {
221 Dune::LocalKey const& localKey = node_->finiteElement().localCoefficients().localKey(i);
222 sizes_[localKey.codim()]++;
223 }
224 auto refElem = Dune::referenceElement<double,GridView::dimension>(node_->element().type());
225 for (size_type c = 0; c <= GridView::dimension ; ++c)
226 sizes_[c] /= refElem.size(c);
227 }
228
230 void unbind()
231 {
232 node_ = nullptr;
233 }
234
236 size_type size() const
237 {
238 return size_;
239 }
240
242 template <class Twist, class It>
243 It fillIn(Twist const& twist, It it, size_type shift = 0) const
244 {
245 assert(node_ != nullptr);
246 const auto& gridIdSet = gridView_.grid().globalIdSet();
247
248 for (size_type i = 0; i < size_ ; ++i, ++it) {
249 Dune::LocalKey localKey = node_->finiteElement().localCoefficients().localKey(i);
250 unsigned int s = localKey.subEntity();
251 unsigned int c = localKey.codim();
252 it->first = {gridIdSet.subId(node_->element(), s, c), shift + twist.get(localKey,sizes_[c])};
253
254 it->second = Dune::Hybrid::switchCases(std::make_index_sequence<dim+1>{}, c,
255 [&](auto codim) { return node_->element().template subEntity<codim>(s).partitionType(); },
256 [&]() {
257 error_exit("Invalid codimension c = {}", c);
258 return Dune::PartitionType{};
259 });
260 }
261
262 return it;
263 }
264
265 protected:
266 GridView gridView_;
267 const Node* node_ = nullptr;
268 size_type size_ = 0;
269 std::array<unsigned int,GridView::dimension+1> sizes_;
270 };
271
272
273 // Specialization for PowerBasis
274 template <class PreBasis, class TP>
275 class NodeIdSet<PreBasis, TP, Dune::TypeTree::PowerNodeTag>
276 {
277 public:
278 using Node = typename PreBasis::Node;
279 using GridView = typename PreBasis::GridView;
280 using size_type = std::size_t;
281
282 protected:
283 using SubPreBasis = typename PreBasis::SubPreBasis;
284 using SubTreePath = decltype(Dune::TypeTree::push_back(std::declval<TP>(), std::size_t(0)));
285 using SubNodeIdSet = NodeIdSet<SubPreBasis, SubTreePath>;
286
287 public:
288 NodeIdSet(GridView const& gridView)
289 : subIds_(gridView)
290 {}
291
293 void bind(const Node& node)
294 {
295 node_ = &node;
296 subIds_.bind(node.child(0));
297 }
298
300 void unbind()
301 {
302 node_ = nullptr;
303 subIds_.unbind();
304 }
305
307 size_type size() const
308 {
309 return node_->size();
310 }
311
313 template <class Twist, class It>
314 It fillIn(Twist const& twist, It it, size_type shift = 0) const
315 {
316 assert(node_ != nullptr);
317 for (std::size_t child = 0; child < node_->degree(); ++child)
318 {
319 size_type subTreeSize = subIds_.size();
320 it = subIds_.fillIn(twist, it, shift);
321 shift += subTreeSize;
322 }
323 return it;
324 }
325
326 protected:
327 SubNodeIdSet subIds_;
328 const Node* node_ = nullptr;
329 };
330
331
332 // Specialization for CompositePreBasis
333 template <class PreBasis, class TP>
334 class NodeIdSet<PreBasis, TP, Dune::TypeTree::CompositeNodeTag>
335 {
336 public:
337 using Node = typename PreBasis::Node;
338 using GridView = typename PreBasis::GridView;
339 using size_type = std::size_t;
340
341 protected:
342 static const std::size_t children = Node::degree();
343 using ChildIndices = std::make_index_sequence<children>;
344
345 // The I'th SubPreBasis
346 template <std::size_t I>
347 using SubPreBasis = typename PreBasis::template SubPreBasis<I>;
348
349 template <std::size_t I>
350 using SubTreePath = decltype(Dune::TypeTree::push_back(std::declval<TP>(), Dune::index_constant<I>{}));
351
352 template <std::size_t I>
353 using SubNodeIdSet = NodeIdSet<SubPreBasis<I>, SubTreePath<I>>;
354
355 // A tuple of NodeIdSets for the SubPreBases
356 using IdsTuple = IndexMapTuple_t<std::make_index_sequence<children>, SubNodeIdSet>;
357
358 public:
359 NodeIdSet(GridView const& gridView)
360 : idsTuple_(Ranges::applyIndices<children>([&](auto... i) {
361 return std::make_tuple(SubNodeIdSet<VALUE(i)>(gridView)...);
362 }))
363 {}
364
366 void bind(const Node& node)
367 {
368 node_ = &node;
369 Ranges::forIndices<0,children>([&](auto i) {
370 std::get<VALUE(i)>(idsTuple_).bind(node.child(i));
371 });
372 }
373
375 void unbind()
376 {
377 node_ = nullptr;
378 Ranges::forEach(idsTuple_, [](auto& ids) {
379 ids.unbind();
380 });
381 }
382
384 size_type size() const
385 {
386 return node_->size();
387 }
388
390 template <class Twist, class It>
391 It fillIn(Twist const& twist, It it, size_type shift = 0) const
392 {
393 assert(node_ != nullptr);
394 Ranges::forEach(idsTuple_, [&](auto const& ids)
395 {
396 size_type subTreeSize = ids.size();
397 it = ids.fillIn(twist, it, shift);
398 shift += subTreeSize;
399 });
400 return it;
401 }
402
403 private:
404 IdsTuple idsTuple_;
405 const Node* node_ = nullptr;
406 };
407
408
409 template <class GV, int k, class TP>
410 class NodeIdSet<Dune::Functions::LagrangeDGPreBasis<GV, k>, TP>
411 {
412 public:
413 using PreBasis = Dune::Functions::LagrangeDGPreBasis<GV, k>;
414 using Node = typename PreBasis::Node;
415 using GridView = typename PreBasis::GridView;
416 using size_type = std::size_t;
417
418 NodeIdSet(GridView const& gridView)
419 : gridView_(gridView)
420 {}
421
423 void bind(const Node& node)
424 {
425 node_ = &node;
426 size_ = node_->finiteElement().size();
427 }
428
430 void unbind()
431 {
432 node_ = nullptr;
433 }
434
436 size_type size() const
437 {
438 return size_;
439 }
440
442 template <class Twist, class It>
443 It fillIn(Twist const& /*twist*/, It it, size_type shift = 0) const
444 {
445 assert(node_ != nullptr);
446 const auto& gridIdSet = gridView_.grid().globalIdSet();
447 auto elementId = gridIdSet.id(node_->element());
448
449 for (size_type i = 0; i < size_; ++i, ++it)
450 {
451 it->first = {elementId, shift + i};
452 it->second = node_->element().partitionType();
453 }
454
455 return it;
456 }
457
458 protected:
459 GridView gridView_;
460 const Node* node_ = nullptr;
461 size_type size_ = 0;
462 };
463
464
465} // end namespace AMDiS
Provide global ids for all DOFs in a global basis.
Definition: GlobalIdSet.hpp:95
EntityIdType entityId(size_type i) const
Return the global unique ID of the entity associated to the ith DOF on the currently bound element in...
Definition: GlobalIdSet.hpp:159
size_type size() const
The number of DOFs on the current element in the basis tree.
Definition: GlobalIdSet.hpp:144
PartitionType partitionType(size_type i) const
Return the partition type of the ith DOF on the currently bound element in the basis tree.
Definition: GlobalIdSet.hpp:166
void bind(Element const &element)
Bind the IdSet to a grid element.
Definition: GlobalIdSet.hpp:128
IdType id(size_type i) const
Return the global unique ID of the ith DOF on the currently bound element in the basis tree.
Definition: GlobalIdSet.hpp:151
void unbind()
unbind from the element
Definition: GlobalIdSet.hpp:138
PB PreBasis
Pre-basis providing the implementation details.
Definition: GlobalBasis.hpp:54
typename PreBasis::GridView GridView
The grid view that the FE space is defined on.
Definition: GlobalBasis.hpp:57
Type of a global used used to represent DOFs uniquely.
Definition: GlobalIdSet.hpp:37
GlobalIdType & operator++()
Increment the second index only.
Definition: GlobalIdSet.hpp:48
size_type size() const
Number of DOFs in this node.
Definition: GlobalIdSet.hpp:236
void bind(const Node &node)
Bind the idset to a tree node.
Definition: GlobalIdSet.hpp:214
It fillIn(Twist const &twist, It it, size_type shift=0) const
Maps from subtree index set [0..size-1] to a globally unique id in global basis.
Definition: GlobalIdSet.hpp:243
void unbind()
Unbind the idset.
Definition: GlobalIdSet.hpp:230
void bind(Element const &element)
Bind the twist to a codim-0 element. Calculates the global ids of all the element vertices.
Definition: Twist.hpp:29
unsigned int get(Dune::LocalKey const &localKey, unsigned int size) const
Get the permutated local dof index, living on a subEntity of the bound element.
Definition: Twist.hpp:48