DUNE PDELab (git)

gridviewordering.hh
1// -*- tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2// vi: set et ts=8 sw=2 sts=2:
3
4#ifndef DUNE_PDELAB_ORDERING_GRIDVIEWORDERING_HH
5#define DUNE_PDELAB_ORDERING_GRIDVIEWORDERING_HH
6
7#include <dune/typetree/typetree.hh>
8
9#include <dune/pdelab/ordering/utility.hh>
10#include <dune/pdelab/ordering/localorderingbase.hh>
11#include <dune/pdelab/ordering/orderingbase.hh>
12#include <dune/pdelab/ordering/leaflocalordering.hh>
13#include <dune/pdelab/ordering/lexicographicordering.hh>
14#include <dune/pdelab/ordering/leafgridviewordering.hh>
15
16namespace Dune {
17 namespace PDELab {
18
21
22 template<typename Codims>
23 struct collect_used_codims
24 : public TypeTree::TreeVisitor
25 , public TypeTree::DynamicTraversal
26 {
27
28 template<typename Node, typename TreePath>
29 void leaf(Node& node, TreePath tp)
30 {
31 node.collect_used_codims(codims);
32 }
33
34 collect_used_codims(Codims& codims_)
35 : codims(codims_)
36 {}
37
38 Codims& codims;
39
40 };
41
42
43 struct collect_a_priori_fixed_size
44 : public TypeTree::TreeVisitor
45 , public TypeTree::DynamicTraversal
46 {
47
48 template<typename Node, typename TreePath>
49 void leaf(Node& node, TreePath tp)
50 {
51 node.update_a_priori_fixed_size();
52 any = any || node._fixed_size;
53 all = all && node._fixed_size;
54 }
55
56 template<typename Node, typename TreePath>
57 void pre(Node& node, TreePath tp) const
58 {
59 node._fixed_size = true;
60 }
61
62 template<typename Node, typename Child, typename TreePath, typename ChildIndex>
63 void afterChild(Node& node, const Child& child, TreePath tp, ChildIndex childIndex) const
64 {
65 node._fixed_size = node._fixed_size && child._fixed_size;
66 }
67
68 collect_a_priori_fixed_size()
69 : any(false)
70 , all(true)
71 {}
72
73 bool any;
74 bool all;
75
76 };
77
78
79 template<typename ES>
80 struct update_fixed_size
81 : public TypeTree::TreeVisitor
82 , public TypeTree::DynamicTraversal
83 {
84
85 template<typename Node, typename TreePath>
86 void leaf(Node& node, TreePath tp) const
87 {
88 if (node._fixed_size)
89 {
90 typedef typename Node::Traits::SizeType size_type;
91 const size_type dim = ES::dimension;
92 node._codim_used.reset();
93 node._gt_used.assign(GlobalGeometryTypeIndex::size(dim),false);
94 node._gt_dof_offsets.assign(GlobalGeometryTypeIndex::size(dim),0);
95 for (const auto& gt : es.indexSet().types())
96 {
97 size_type size = node.finiteElementMap().size(gt);
98 node._gt_dof_offsets[GlobalGeometryTypeIndex::index(gt)] = size;
99 node._gt_used[GlobalGeometryTypeIndex::index(gt)] = size > 0;
100 node._codim_used[dim - gt.dim()] = node._codim_used[dim - gt.dim()] || (size > 0);
101 }
102 node._max_local_size = node.finiteElementMap().maxLocalSize();
103 }
104 }
105
106 template<typename Node, typename TreePath>
107 void pre(Node& node, TreePath tp) const
108 {
109 if (node._fixed_size)
110 {
111 typedef typename Node::Traits::SizeType size_type;
112 const size_type dim = ES::dimension;
113 node._codim_used.reset();
114 node._gt_used.assign(Dune::GlobalGeometryTypeIndex::size(dim),false);
115 node._gt_dof_offsets.assign(Dune::GlobalGeometryTypeIndex::size(dim) * TypeTree::degree(node),0);
116 node._max_local_size = 0;
117 }
118 }
119
120 template<typename Node, typename Child, typename TreePath, typename ChildIndex>
121 void afterChild(Node& node, const Child& child, TreePath tp, ChildIndex childIndex) const
122 {
123 if (node._fixed_size)
124 {
125 node._codim_used |= child._codim_used;
126
127 std::transform(node._gt_used.begin(),
128 node._gt_used.end(),
129 child._gt_used.begin(),
130 node._gt_used.begin(),
131 std::logical_or<bool>());
132
133 node._max_local_size += child._max_local_size;
134
135 typedef typename Node::Traits::SizeType size_type;
136
137 const size_type per_gt_size = child._child_count > 0 ? child._child_count : 1;
138 const size_type size_offset = child._child_count > 0 ? child._child_count - 1 : 0;
139
140 for (size_type gt = 0; gt < Dune::GlobalGeometryTypeIndex::size(ES::dimension); ++gt)
141 node._gt_dof_offsets[gt * TypeTree::degree(node) + childIndex] = child._gt_dof_offsets[gt * per_gt_size + size_offset];
142 }
143 }
144
145 template<typename Node, typename TreePath>
146 void post(Node& node, TreePath tp) const
147 {
148 if (node._fixed_size)
149 {
150 typedef typename std::vector<typename Node::Traits::SizeType>::iterator iterator;
151
152 iterator next_gt_it = node._gt_dof_offsets.begin() + TypeTree::degree(node);
153 const iterator end_it = node._gt_dof_offsets.end();
154
155 for (iterator it = node._gt_dof_offsets.begin();
156 it != end_it;
157 it += TypeTree::degree(node), next_gt_it += TypeTree::degree(node))
158 std::partial_sum(it,next_gt_it,it);
159 }
160 }
161
162 update_fixed_size(const ES es_)
163 : es(es_)
164 {}
165
166 ES es;
167
168 };
169
170
171 struct pre_collect_used_geometry_types
172 : public TypeTree::TreeVisitor
173 , public TypeTree::DynamicTraversal
174 {
175
176 template<typename Node, typename TreePath>
177 void leaf(Node& node, TreePath tp) const
178 {
179 if (!node._fixed_size)
180 {
181 node._codim_used.reset();
182 node._gt_used.assign(Dune::GlobalGeometryTypeIndex::size(dim),false);
183 node._gt_dof_offsets.assign(Dune::GlobalGeometryTypeIndex::size(dim) * std::max(node._child_count,static_cast<std::size_t>(1)),Node::GT_UNUSED);
184 node._gt_entity_offsets.assign(Dune::GlobalGeometryTypeIndex::size(dim) + 1,0);
185 }
186 }
187
188 template<typename Node, typename TreePath>
189 void pre(Node& node, TreePath tp) const
190 {
191 leaf(node,tp);
192 }
193
194 pre_collect_used_geometry_types(std::size_t dimension)
195 : dim(dimension)
196 {}
197
198 const std::size_t dim;
199
200 };
201
202
203 template<typename Cell>
204 struct collect_used_geometry_types_from_cell_visitor
205 : public TypeTree::TreeVisitor
206 , public TypeTree::DynamicTraversal
207 {
208
209 template<typename Node, typename TreePath>
210 void leaf(Node& node, TreePath tp) const
211 {
212 if (!node._fixed_size)
213 node.collect_used_geometry_types_from_cell(cell);
214 }
215
216 collect_used_geometry_types_from_cell_visitor(const Cell& cell_)
217 : cell(cell_)
218 , ref_el(Dune::ReferenceElements<typename Cell::Geometry::ctype,Cell::dimension>::general(cell_.type()))
219 {}
220
221 const Cell& cell;
222 Dune::ReferenceElement<typename Cell::Geometry> ref_el;
223
224 };
225
226
227 template<typename ES>
228 struct post_collect_used_geometry_types
229 : public TypeTree::TreeVisitor
230 , public TypeTree::DynamicTraversal
231 {
232
233 template<typename Node, typename TreePath>
234 void leaf(Node& node, TreePath tp) const
235 {
236 if (!node._fixed_size)
237 {
238 typedef typename Node::Traits::SizeType size_type;
239
240 for (const auto& gt : es.indexSet().types())
241 {
242 if (node._gt_used[Dune::GlobalGeometryTypeIndex::index(gt)])
243 node._gt_entity_offsets[Dune::GlobalGeometryTypeIndex::index(gt) + 1] = es.indexSet().size(gt);
244 }
245
246 std::partial_sum(node._gt_entity_offsets.begin(),node._gt_entity_offsets.end(),node._gt_entity_offsets.begin());
247 node._entity_dof_offsets.assign(node._gt_entity_offsets.back() * std::max(node._child_count,static_cast<size_type>(1)),0);
248 node.setup_fixed_size_possible();
249 }
250 }
251
252 template<typename Node, typename Child, typename TreePath, typename ChildIndex>
253 void afterChild(Node& node, const Child& child, TreePath tp, ChildIndex childIndex) const
254 {
255 if (!node._fixed_size)
256 {
257 node._codim_used |= child._codim_used;
258
259 std::transform(node._gt_used.begin(),
260 node._gt_used.end(),
261 child._gt_used.begin(),
262 node._gt_used.begin(),
263 std::logical_or<bool>());
264 }
265 }
266
267 template<typename Node, typename TreePath>
268 void post(Node& node, TreePath tp) const
269 {
270 leaf(node,tp);
271 }
272
273 post_collect_used_geometry_types(const ES& es_)
274 : es(es_)
275 {}
276
277 ES es;
278
279 };
280
281
282 template<typename ES>
283 struct extract_per_entity_sizes_from_cell_visitor
284 : public TypeTree::TreeVisitor
285 , public TypeTree::DynamicTraversal
286 {
287
288 static const std::size_t dim = ES::dimension;
289 typedef typename ES::template Codim<0>::Entity Cell;
290 typedef std::size_t size_type;
291
292 template<typename Node, typename TreePath>
293 void leaf(Node& node, TreePath tp)
294 {
295 if (!node._fixed_size)
296 node.extract_per_entity_sizes_from_cell(*cell,gt_sizes);
297 }
298
299 extract_per_entity_sizes_from_cell_visitor(const ES& es_)
300 : es(es_)
301 , cell(nullptr)
302 , ref_el()
303 , gt_sizes(Dune::GlobalGeometryTypeIndex::size(dim),0)
304 {}
305
306 void set_cell(const Cell& cell_)
307 {
308 cell = &cell_;
309 ref_el = referenceElement(cell_.geometry());
310 }
311
312 ES es;
313 const Cell* cell;
314 Dune::ReferenceElement<typename Cell::Geometry> ref_el;
315 std::vector<size_type> gt_sizes;
316
317 };
318
319
320 template<typename ES>
321 struct post_extract_per_entity_sizes
322 : public TypeTree::TreeVisitor
323 , public TypeTree::DynamicTraversal
324 {
325
326 typedef std::vector<GeometryType> GTVector;
327
328
329 template<typename Node, typename TreePath>
330 void leaf(Node& node, TreePath tp) const
331 {
332 if (!node._fixed_size)
333 {
334 // mask out GT_UNUSED for geometry types that really weren't used
335 for (auto& size : node._gt_dof_offsets)
336 if (size == Node::GT_UNUSED)
337 size = 0;
338 if (node._fixed_size_possible)
339 {
340 node._entity_dof_offsets = std::vector<typename Node::Traits::SizeType>();
341 node._fixed_size = true;
342 }
343 }
344 }
345
346 template<typename Node, typename TreePath>
347 void pre(Node& node, TreePath tp) const
348 {
349 if (!node._fixed_size)
350 {
351 node._fixed_size_possible = true;
352 node._max_local_size = 0;
353 }
354 }
355
356
357 template<typename Node, typename Child, typename TreePath, typename ChildIndex>
358 void afterChild(Node& node, const Child& child, TreePath tp, ChildIndex childIndex) const
359 {
360 if (!node._fixed_size)
361 {
362 node._fixed_size_possible = node._fixed_size_possible && child._fixed_size;
363 node._max_local_size += child._max_local_size;
364 }
365 }
366
367
368 template<typename Node, typename TreePath>
369 void post(Node& node, TreePath tp) const
370 {
371 if (!node._fixed_size)
372 {
373
374 typedef typename Node::Traits::SizeType size_type;
375 const size_type dim = ES::dimension;
376
377 if (node._fixed_size_possible)
378 {
379
380 for (size_type gt = 0; gt < GlobalGeometryTypeIndex::size(ES::dimension); ++gt)
381 {
382 for (size_type child_index = 0; child_index < TypeTree::degree(node); ++child_index)
383 {
384 const size_type per_gt_size = node.childOrdering(child_index)._child_count > 0 ? node.childOrdering(child_index)._child_count : 1;
385 const size_type size_offset = node.childOrdering(child_index)._child_count > 0 ? node.childOrdering(child_index)._child_count - 1 : 0;
386
387 node._gt_dof_offsets[gt * TypeTree::degree(node) + child_index] = node.childOrdering(child_index)._gt_dof_offsets[gt * per_gt_size + size_offset];
388 }
389 }
390
391 typedef typename std::vector<typename Node::Traits::SizeType>::iterator iterator;
392
393 const iterator end_it = node._gt_dof_offsets.end();
394
395 for (iterator it = node._gt_dof_offsets.begin();
396 it != end_it;
397 it += TypeTree::degree(node))
398 std::partial_sum(it,it + TypeTree::degree(node),it);
399
400 node._fixed_size = true;
401 }
402 else
403 {
404 typedef typename Node::Traits::SizeType size_type;
405
406 size_type index = 0;
407 for (size_type geometry_type_index = 0; geometry_type_index < GlobalGeometryTypeIndex::size(dim); ++geometry_type_index)
408 {
409 if (!node._gt_used[geometry_type_index])
410 continue;
411 const size_type entity_count = node._gt_entity_offsets[geometry_type_index+1] - node._gt_entity_offsets[geometry_type_index];
412 for (size_type entity_index = 0; entity_index < entity_count; ++entity_index)
413 {
414 size_type carry = 0;
415 for (size_type child_index = 0; child_index < TypeTree::degree(node); ++child_index)
416 node._entity_dof_offsets[index++] = (carry += node.childOrdering(child_index).size(geometry_type_index,entity_index));
417 }
418 }
419
420 }
421 }
422 }
423
424 post_extract_per_entity_sizes(const ES& es_)
425 : es(es_)
426 {}
427
428 ES es;
429
430 };
431
433 template<typename LocalOrdering>
434 class GridViewOrdering
435 : public TypeTree::CompositeNode<LocalOrdering>
436 , public VirtualOrderingBase<typename LocalOrdering::Traits::DOFIndex,
437 typename LocalOrdering::Traits::ContainerIndex>
438 , public OrderingBase<typename LocalOrdering::Traits::DOFIndex,
439 typename LocalOrdering::Traits::ContainerIndex>
440 {
441 public:
442 typedef typename LocalOrdering::Traits Traits;
443
444 static const bool has_dynamic_ordering_children = false;
445
446 static const bool consume_tree_index = false;
447
448 private:
449
450 typedef TypeTree::CompositeNode<LocalOrdering> NodeT;
451 typedef OrderingBase<
452 typename LocalOrdering::Traits::DOFIndex,
453 typename LocalOrdering::Traits::ContainerIndex
454 > BaseT;
455
456 using EntitySet = typename Traits::EntitySet;
457
458 public:
460
465 GridViewOrdering(const typename NodeT::NodeStorage &local_ordering,
466 bool container_blocked,
467 typename BaseT::GFSData *gfs_data,
468 const EntitySet &entity_Set)
469 : NodeT(local_ordering),
470 BaseT(*this, container_blocked, gfs_data, this), _es(entity_Set) {
471 // make sure to switch off container blocking handling in the local ordering,
472 // we already handle it in the GridViewOrdering
473 localOrdering().disable_container_blocking();
474 }
475
476#ifndef DOXYGEN
477
478// we need to override the default copy / move ctor to fix the delegate pointer, but that is
479// hardly interesting to our users...
480
481 GridViewOrdering(const GridViewOrdering& r)
482 : NodeT(r.nodeStorage())
483 , BaseT(r)
484 , _es(r._es)
485 , _gt_dof_offsets(r._gt_dof_offsets)
486 , _gt_entity_offsets(r._gt_entity_offsets)
487 , _entity_dof_offsets(r._entity_dof_offsets)
488 {
489 this->setDelegate(this);
490 }
491
492 GridViewOrdering(GridViewOrdering&& r)
493 : NodeT(r.nodeStorage())
494 , BaseT(std::move(r))
495 , _es(std::move(r._es))
496 , _gt_dof_offsets(std::move(r._gt_dof_offsets))
497 , _gt_entity_offsets(std::move(r._gt_entity_offsets))
498 , _entity_dof_offsets(std::move(r._entity_dof_offsets))
499 {
500 this->setDelegate(this);
501 }
502
503 virtual ~GridViewOrdering() override = default;
504
505#endif // DOXYGEN
506
507 using BaseT::size;
508
514 typename Traits::SizeType size(typename Traits::ContainerIndex suffix) const
515 {
516 using size_type = typename Traits::SizeType;
517 if (suffix.size() == Traits::ContainerIndex::max_depth)
518 return 0; // all indices in suffix were consumed, no more sizes to provide
519 if (suffix.size() == 0) // suffix wants the size of this depth
520 return _block_count; // blocked or not, this gives the number of blocks/dofs in next node hierarchy
521
522 // we first have to figure out the entity index
523 typename Traits::DOFIndex::EntityIndex entity_index;
524
525 // the next index to find out its size
526 auto back_index = suffix.back();
527 // we just need to make the inverse computation of the mapIndex funtion to find the entity index
528 if (_container_blocked)
529 suffix.pop_back();
530
531 auto dof_begin = _fixed_size ? _gt_dof_offsets.begin() : _entity_dof_offsets.begin();
532 auto dof_end = _fixed_size ? _gt_dof_offsets.end() : _entity_dof_offsets.end();
533 auto dof_it = std::prev(std::upper_bound(dof_begin, dof_end, back_index));
534 size_type dof_dist = std::distance(dof_begin, dof_it);
535 if (_fixed_size) {
536 // On fixed size, entity index is not used down the tree. Set max to trigger segfault if this does not hold.
537 Traits::DOFIndexAccessor::GeometryIndex::store(entity_index,dof_dist,~size_type{0});
538 } else {
539 auto gt_begin = _gt_entity_offsets.begin();
540 auto gt_end = _gt_entity_offsets.end();
541 auto gt_it = std::prev(std::upper_bound(gt_begin, gt_end, dof_dist));
542 size_type gt = std::distance(gt_begin, gt_it);
543 assert(dof_dist >= *gt_it);
544 size_type ei = dof_dist - *gt_it;
545 Traits::DOFIndexAccessor::GeometryIndex::store(entity_index,gt,ei);
546 }
547
548 // then, the local ordering knows the size for a given entity.
549 return localOrdering().size(suffix, entity_index);
550 }
551
552 LocalOrdering& localOrdering()
553 {
554 return this->template child<0>();
555 }
556
557 const LocalOrdering& localOrdering() const
558 {
559 return this->template child<0>();
560 }
561
562 virtual void map_index_dynamic(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const override
563 {
564 mapIndex(di,ci);
565 }
566
567 typename Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex& di) const
568 {
569 typename Traits::ContainerIndex ci;
570 mapIndex(di.view(),ci);
571 return ci;
572 }
573
574 void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const
575 {
576 typedef typename Traits::SizeType size_type;
577 const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(di);
578 const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(di);
579 localOrdering().map_local_index(geometry_type_index,entity_index,di.treeIndex(),ci);
580 if (_container_blocked)
581 {
582 if (_fixed_size)
583 {
584 ci.push_back(_gt_dof_offsets[geometry_type_index] + entity_index);
585 }
586 else
587 {
588 ci.push_back(_entity_dof_offsets[_gt_entity_offsets[geometry_type_index] + entity_index]);
589 }
590 }
591 else
592 {
593 if (_fixed_size)
594 {
595 ci.back() += _gt_dof_offsets[geometry_type_index] + entity_index * localOrdering().size(geometry_type_index,entity_index);
596 }
597 else
598 {
599 ci.back() += _entity_dof_offsets[_gt_entity_offsets[geometry_type_index] + entity_index];
600 }
601 }
602 }
603
604 template<typename ItIn, typename ItOut>
605 void map_lfs_indices(const ItIn begin, const ItIn end, ItOut out) const
606 {
607 typedef typename Traits::SizeType size_type;
608 if (_container_blocked)
609 {
610 if (_fixed_size)
611 for (ItIn in = begin; in != end; ++in, ++out)
612 {
613 const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
614 const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
615 out->push_back(_gt_dof_offsets[geometry_type_index] + entity_index);
616 }
617 else
618 for (ItIn in = begin; in != end; ++in, ++out)
619 {
620 const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
621 const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
622 out->push_back(_entity_dof_offsets[_gt_entity_offsets[geometry_type_index] + entity_index]);
623 }
624 }
625 else if (_fixed_size)
626 {
627 for (ItIn in = begin; in != end; ++in, ++out)
628 {
629 const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
630 const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
631 out->back() += _gt_dof_offsets[geometry_type_index] + entity_index * localOrdering().size(geometry_type_index,entity_index);
632 }
633 }
634 else
635 {
636 for (ItIn in = begin; in != end; ++in, ++out)
637 {
638 const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
639 const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
640 out->back() += _entity_dof_offsets[_gt_entity_offsets[geometry_type_index] + entity_index];
641 }
642 }
643 }
644
645 template<typename CIOutIterator>
646 typename Traits::SizeType
647 extract_entity_indices(const typename Traits::DOFIndex::EntityIndex& ei,
648 typename Traits::SizeType child_index,
649 CIOutIterator ci_out, const CIOutIterator ci_end) const
650 {
651 typedef typename Traits::SizeType size_type;
652
653 const size_type geometry_type_index = Traits::DOFIndexAccessor::GeometryIndex::geometryType(ei);
654 const size_type entity_index = Traits::DOFIndexAccessor::GeometryIndex::entityIndex(ei);
655
656 if (_container_blocked)
657 {
658 if (_fixed_size)
659 for (; ci_out != ci_end; ++ci_out)
660 {
661 ci_out->push_back(_gt_dof_offsets[geometry_type_index] + entity_index);
662 }
663 else
664 for (; ci_out != ci_end; ++ci_out)
665 {
666 ci_out->push_back(_entity_dof_offsets[_gt_entity_offsets[geometry_type_index] + entity_index]);
667 }
668 }
669 else if (_fixed_size)
670 {
671 for (; ci_out != ci_end; ++ci_out)
672 {
673 ci_out->back() += _gt_dof_offsets[geometry_type_index] + entity_index * localOrdering().size(geometry_type_index,entity_index);
674 }
675 }
676 else
677 {
678 for (; ci_out != ci_end; ++ci_out)
679 {
680 ci_out->back() += _entity_dof_offsets[_gt_entity_offsets[geometry_type_index] + entity_index];
681 }
682 }
683
684 // The return value is not used for non-leaf orderings.
685 return 0;
686 }
687
688 void update()
689 {
690
691 typedef typename Traits::SizeType size_type;
692 using ES = typename Traits::EntitySet;
693 const size_type dim = ES::dimension;
694
695 typename ES::CodimMask codims;
696 codims.set(0); // we always need cells
697
698 TypeTree::applyToTree(localOrdering(),collect_used_codims<typename ES::CodimMask>(codims));
699
700 for (typename ES::dim_type codim = 0; codim <= ES::dimension; ++codim)
701 if (codims.test(codim))
702 _es.addCodim(codim);
703
704 _es.update();
705
706 // Do we already know that we have fixed per-GeometryType sizes?
707 collect_a_priori_fixed_size fixed_size_collector;
708 TypeTree::applyToTree(localOrdering(),fixed_size_collector);
709 _fixed_size = localOrdering().fixedSize();
710
711 const size_type gt_index_count = GlobalGeometryTypeIndex::size(ES::dimension);
712
713 if (fixed_size_collector.any)
714 {
715 // collect used GeometryTypes
716 TypeTree::applyToTree(localOrdering(),update_fixed_size<ES>(_es));
717 }
718
719 if (!fixed_size_collector.all)
720 {
721 TypeTree::applyToTree(localOrdering(),pre_collect_used_geometry_types(ES::dimension));
722
723 using Element = typename ES::template Codim<0>::Entity;
724
725 for (const auto& element : elements(_es))
726 {
727 TypeTree::applyToTree(localOrdering(),collect_used_geometry_types_from_cell_visitor<Element>(element));
728 }
729 TypeTree::applyToTree(localOrdering(),post_collect_used_geometry_types<ES>(_es));
730 // allocate
731
732 //TypeTree::applyToTree(localOrdering(),pre_extract_per_entity_sizes<GV>(_gv));
733 extract_per_entity_sizes_from_cell_visitor<ES> visitor(_es);
734 for (const auto& element : elements(_es))
735 {
736 visitor.set_cell(element);
737 TypeTree::applyToTree(localOrdering(),visitor);
738 }
739 TypeTree::applyToTree(localOrdering(),post_extract_per_entity_sizes<ES>(_es));
740 }
741
742 _codim_used = localOrdering()._codim_used;
743
744 if (localOrdering().fixedSize())
745 {
746 _fixed_size = true;
747 _gt_dof_offsets.assign(gt_index_count + 1,0);
748
749 _block_count = 0;
750
751 _size = 0;
752
753 for (const auto& gt : _es.indexSet().types())
754 {
755 const size_type gt_index = GlobalGeometryTypeIndex::index(gt);
756 size_type gt_size = localOrdering().size(gt_index,0);
757 const size_type gt_entity_count = _es.indexSet().size(gt);
758 _size += gt_size * gt_entity_count;
759 if (_container_blocked)
760 gt_size = gt_size > 0;
761 _gt_dof_offsets[gt_index + 1] = gt_size * gt_entity_count;
762 }
763
764 std::partial_sum(_gt_dof_offsets.begin(),_gt_dof_offsets.end(),_gt_dof_offsets.begin());
765 _block_count = _gt_dof_offsets.back();
766
767 _codim_fixed_size.set();
768
769 }
770 else
771 {
772 _gt_entity_offsets.assign(gt_index_count + 1,0);
773
774 for (const auto& gt : _es.indexSet().types())
775 {
776 if (!localOrdering().contains(gt))
777 continue;
778 const size_type gt_index = GlobalGeometryTypeIndex::index(gt);
779 _gt_entity_offsets[gt_index + 1] = _es.indexSet().size(gt);
780 }
781
782 std::partial_sum(_gt_entity_offsets.begin(),_gt_entity_offsets.end(),_gt_entity_offsets.begin());
783 _entity_dof_offsets.assign(_gt_entity_offsets.back()+1,0);
784 _block_count = 0;
785
786 size_type carry_size = 0;
787 size_type carry_block = 0;
788 size_type index = 0;
789 for (size_type gt_index = 0; gt_index < GlobalGeometryTypeIndex::size(dim); ++gt_index)
790 {
791 if (!localOrdering().contains_geometry_type(gt_index))
792 continue;
793 const size_type entity_count = _gt_entity_offsets[gt_index + 1] - _gt_entity_offsets[gt_index];
794 for (size_type entity_index = 0; entity_index < entity_count; ++entity_index)
795 {
796 const size_type size = localOrdering().size(gt_index,entity_index);
797 carry_size += size;
798 carry_block += (_container_blocked ? (size > 0) : size);
799 _entity_dof_offsets[++index] = carry_block;
800 _block_count += (size > 0);
801 }
802 }
803 _size = carry_size;
804 _block_count = _block_count;
805
806 _codim_fixed_size.reset();
807 }
808
809 _max_local_size = localOrdering().maxLocalSize();
810 }
811
812 using BaseT::fixedSize;
813
814 private:
815
816 using BaseT::_container_blocked;
817 using BaseT::_fixed_size;
818 using BaseT::_max_local_size;
819 using BaseT::_size;
820 using BaseT::_block_count;
821 using BaseT::_codim_used;
822 using BaseT::_codim_fixed_size;
823
824 typename Traits::EntitySet _es;
825 std::vector<typename Traits::SizeType> _gt_dof_offsets;
826 std::vector<typename Traits::SizeType> _gt_entity_offsets;
827 std::vector<typename Traits::SizeType> _entity_dof_offsets;
828
829 };
830
831
833 } // namespace PDELab
834} // namespace Dune
835
836#endif // DUNE_PDELAB_ORDERING_GRIDVIEWORDERING_HH
Dune::GlobalGeometryTypeIndex::index
static constexpr std::size_t index(const GeometryType &gt)
Compute the index for the given geometry type over all dimensions.
Definition: typeindex.hh:138
Dune::GlobalGeometryTypeIndex::size
static constexpr std::size_t size(std::size_t maxdim)
Compute total number of geometry types up to and including the given dimension.
Definition: typeindex.hh:125
Dune::PDELab::GridViewOrdering
Transforms a local ordering (entity-wise order) into a global ordering.
Definition: gridviewordering.hh:440
Dune::PDELab::GridViewOrdering::GridViewOrdering
GridViewOrdering(const typename NodeT::NodeStorage &local_ordering, bool container_blocked, typename BaseT::GFSData *gfs_data, const EntitySet &entity_Set)
Construct ordering object.
Definition: gridviewordering.hh:465
Dune::PDELab::GridViewOrdering::size
Traits::SizeType size(typename Traits::ContainerIndex suffix) const
Gives the size for a given suffix.
Definition: gridviewordering.hh:514
Dune::TypeTree::CompositeNode
Base class for composite nodes based on variadic templates.
Definition: compositenode.hh:28
Dune::TypeTree::CompositeNode< LocalOrdering >::NodeStorage
std::tuple< std::shared_ptr< Children >... > NodeStorage
The type used for storing the children.
Definition: compositenode.hh:36
Dune::FloatCmp::gt
bool gt(const T &first, const T &second, typename EpsilonType< T >::Type epsilon)
test if first greater than second
Definition: float_cmp.cc:158
Dune::referenceElement
unspecified value type referenceElement(T &&... t)
Returns a reference element for the objects t....
Dune::Hybrid::max
constexpr auto max
Function object that returns the greater of the given values.
Definition: hybridutilities.hh:484
Dune::TypeTree::degree
std::size_t degree(const Node &node)
Returns the degree of node as run time information.
Definition: nodeinterface.hh:79
Dune::TypeTree::applyToTree
void applyToTree(Tree &&tree, Visitor &&visitor)
Apply visitor to TypeTree.
Definition: traversal.hh:239
Dune::TypeTree::Child
typename impl::_Child< Node, indices... >::type Child
Template alias for the type of a child node given by a list of child indices.
Definition: childextraction.hh:225
Dune::TypeTree::child
ImplementationDefined child(Node &&node, Indices... indices)
Extracts the child of a node given by a sequence of compile-time and run-time indices.
Definition: childextraction.hh:128
Dune::Partitions::all
constexpr All all
PartitionSet for all partitions.
Definition: partitionset.hh:295
Dune
Dune namespace.
Definition: alignedallocator.hh:13
Dune::size
constexpr std::integral_constant< std::size_t, sizeof...(II)> size(std::integer_sequence< T, II... >)
Return the size of the sequence.
Definition: integersequence.hh:75
Dune::contains
constexpr std::bool_constant<((II==value)||...)> contains(std::integer_sequence< T, II... >, std::integral_constant< T, value >)
Checks whether or not a given sequence contains a value.
Definition: integersequence.hh:137
std
STL namespace.
Dune::TypeTree::DefaultVisitor::afterChild
void afterChild(T &&, Child &&, TreePath, ChildIndex) const
Method for child-parent traversal.
Definition: visitor.hh:122
Dune::TypeTree::DefaultVisitor::post
void post(T &&, TreePath) const
Method for postfix tree traversal.
Definition: visitor.hh:83
Dune::TypeTree::DefaultVisitor::leaf
void leaf(T &&, TreePath) const
Method for leaf traversal.
Definition: visitor.hh:93
Dune::TypeTree::DefaultVisitor::pre
void pre(T &&, TreePath) const
Method for prefix tree traversal.
Definition: visitor.hh:60
fixedSize
std::size_t fixedSize
The number of data items per index if it is fixed, 0 otherwise.
Definition: variablesizecommunicator.hh:264
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 12, 23:30, 2024)