dune-pdelab  2.4-dev
orderingbase.hh
Go to the documentation of this file.
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_ORDERINGBASE_HH
5 #define DUNE_PDELAB_ORDERING_ORDERINGBASE_HH
6 
11 
12 #include <vector>
13 
14 namespace Dune {
15  namespace PDELab {
16 
19 
20  template<typename DI, typename CI>
22  : public PartitionInfoProvider
23  {
24 
25  template<typename size_type>
26  friend struct ::Dune::PDELab::impl::update_ordering_data;
27 
28  public:
29 
31 
32  protected:
33 
34  typedef Dune::PDELab::impl::GridFunctionSpaceOrderingData<typename Traits::SizeType> GFSData;
35 
36  public:
37 
39 
41 
42  static const bool has_dynamic_ordering_children = true;
43 
44  static const bool consume_tree_index = true;
45 
46  typename Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex& di) const
47  {
48  typename Traits::ContainerIndex ci;
49  mapIndex(di.view(),ci);
50  return ci;
51  }
52 
53  void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const
54  {
55  if (_delegate)
56  _delegate->map_index_dynamic(di,ci);
57  else
58  {
59  _mapIndex(di,ci);
60  }
61  }
62 
63  typename Traits::SizeType size() const
64  {
65  return _size;
66  }
67 
68  typename Traits::SizeType blockCount() const
69  {
70  return _block_count;
71  }
72 
73  typename Traits::SizeType size(const typename Traits::SizeType child_index) const
74  {
75  return _child_size_offsets[child_index + 1] - _child_size_offsets[child_index];
76  }
77 
78  typename Traits::SizeType sizeOffset(const typename Traits::SizeType child_index) const
79  {
80  return _child_size_offsets[child_index];
81  }
82 
83  typename Traits::SizeType blockOffset(const typename Traits::SizeType child_index) const
84  {
86  return _child_block_offsets[child_index];
87  }
88 
89  typename Traits::SizeType maxLocalSize() const
90  {
91  return _max_local_size;
92  }
93 
95  {
96  return _merge_mode;
97  }
98 
99  void update()
100  {
101  std::fill(_child_size_offsets.begin(),_child_size_offsets.end(),0);
102  std::fill(_child_block_offsets.begin(),_child_block_offsets.end(),0);
103  _codim_used.reset();
104  _codim_fixed_size.set();
105  _max_local_size = 0;
106  _block_count = 0;
107  typename Traits::SizeType block_carry = 0;
108  typename Traits::SizeType size_carry = 0;
109  for (typename Traits::SizeType i = 0; i < _child_count; ++i)
110  {
111  _child_block_offsets[i+1] = (block_carry += _children[i]->blockCount());
112  _child_size_offsets[i+1] = (size_carry += _children[i]->size());
113  _codim_used |= _children[i]->_codim_used;
114  _codim_fixed_size &= _children[i]->_codim_fixed_size;
115  _block_count += _children[i]->blockCount();
116  _max_local_size += _children[i]->maxLocalSize();
117  }
118  if (_container_blocked)
120  {
122  }
123  else
124  {
125  if (_child_block_offsets.back() % _child_block_merge_offsets.back() != 0)
126  DUNE_THROW(OrderingStructureError,
127  "Invalid ordering structure: "
128  << "total number of blocks ("
129  << _child_block_offsets.back()
130  << ") is not a multiple of the interleaved block size ("
132  << ")."
133  );
135  }
136  else
138  _size = _child_size_offsets.back();
139  }
140 
141  template<typename Node>
142  OrderingBase(Node& node,
143  bool container_blocked,
144  GFSData* gfs_data,
145  VirtualOrderingBase<DI,CI>* delegate = nullptr)
146  : _fixed_size(false)
147  , _container_blocked(container_blocked)
148  , _merge_mode(MergeMode::lexicographic)
149  , _child_count(Node::has_dynamic_ordering_children ? Node::CHILDREN : 0)
150  , _children(_child_count,nullptr)
151  , _child_size_offsets((_child_count > 0 ? _child_count + 1 : 0),0)
152  , _child_block_offsets((_child_count > 0 ? _child_count + 1 : 0),0)
153  , _max_local_size(0)
154  , _size(0)
155  , _block_count(0)
156  , _delegate(delegate)
157  , _gfs_data(gfs_data)
158  {
159  TypeTree::applyToTree(node,extract_child_bases<OrderingBase>(_children));
160  // We contain all grid PartitionTypes that any of our children contain.
161  mergePartitionSets(_children.begin(),_children.end());
162  }
163 
164  template<typename Node>
165  OrderingBase(Node& node,
166  bool container_blocked,
167  const std::vector<std::size_t>& merge_offsets,
168  GFSData* gfs_data,
169  VirtualOrderingBase<DI,CI>* delegate = nullptr)
170  : _fixed_size(false)
171  , _container_blocked(container_blocked)
172  , _merge_mode(MergeMode::interleaved)
173  , _child_count(Node::has_dynamic_ordering_children ? Node::CHILDREN : 0)
174  , _children(_child_count,nullptr)
175  , _child_size_offsets((_child_count > 0 ? _child_count + 1 : 0),0)
176  , _child_block_offsets((_child_count > 0 ? _child_count + 1 : 0),0)
177  , _child_block_merge_offsets(merge_offsets)
178  , _max_local_size(0)
179  , _size(0)
180  , _block_count(0)
181  , _delegate(delegate)
182  , _gfs_data(gfs_data)
183  {
184  TypeTree::applyToTree(node,extract_child_bases<OrderingBase>(_children));
185  // We contain all grid PartitionTypes that any of our children contain.
186  mergePartitionSets(_children.begin(),_children.end());
187  }
188 
189 
190  bool containerBlocked() const
191  {
192  return _container_blocked;
193  }
194 
195  std::size_t childOrderingCount() const
196  {
197  return _child_count;
198  }
199 
201  {
202  return *_children[i];
203  }
204 
205  const OrderingBase& childOrdering(typename Traits::SizeType i) const
206  {
207  return *_children[i];
208  }
209 
210  bool contains(typename Traits::SizeType codim) const
211  {
212  return _codim_used.test(codim);
213  }
214 
215  bool fixedSize(typename Traits::SizeType codim) const
216  {
217  return _codim_fixed_size.test(codim);
218  }
219 
220  protected:
221 
223 
229  {
230  _delegate = delegate;
231  }
232 
233  void _mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const
234  {
235  typedef typename Traits::SizeType size_type;
236  size_type child_index = di.treeIndex().back();
237  _children[child_index]->mapIndex(di.back_popped(),ci);
239  {
240  if (_container_blocked)
241  ci.push_back(child_index);
242  else
243  ci.back() += blockOffset(child_index);
244  }
245  else
246  {
247  size_type child_block_offset = _child_block_merge_offsets[child_index];
248  size_type child_block_size = _child_block_merge_offsets[child_index + 1] - child_block_offset;
249  size_type block_index = ci.back() / child_block_size;
250  size_type offset = ci.back() % child_block_size;
251  if (_container_blocked)
252  {
253  ci.back() = child_block_offset + offset;
254  ci.push_back(block_index);
255  }
256  else
257  {
258  size_type block_size = _child_block_merge_offsets.back();
259  ci.back() = block_index * block_size + child_block_offset + offset;
260  }
261  }
262  }
263 
264  private:
265 
266  bool update_gfs_data_size(typename Traits::SizeType& size, typename Traits::SizeType& block_count) const
267  {
268  size = _size;
269  block_count = _block_count;
270  return true;
271  }
272 
273  public:
274 
276  const bool _container_blocked;
278 
279  const std::size_t _child_count;
280  std::vector<OrderingBase*> _children;
281 
282  std::vector<typename Traits::SizeType> _child_size_offsets;
283  std::vector<typename Traits::SizeType> _child_block_offsets;
284  std::vector<typename Traits::SizeType> _child_block_merge_offsets;
287 
288  std::size_t _max_local_size;
289  std::size_t _size;
290  std::size_t _block_count;
291 
293  GFSData* _gfs_data;
294 
295  };
296 
298 
299  } // namespace PDELab
300 } // namespace Dune
301 
302 #endif // DUNE_PDELAB_ORDERING_ORDERINGBASE_HH
const std::size_t offset
Definition: localfunctionspace.hh:74
void _mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex &ci) const
Definition: orderingbase.hh:233
const OrderingBase & childOrdering(typename Traits::SizeType i) const
Definition: orderingbase.hh:205
std::vector< typename Traits::SizeType > _child_size_offsets
Definition: orderingbase.hh:282
Index merging algorithm for global orderings.
Definition: ordering/utility.hh:21
std::size_t _size
Definition: orderingbase.hh:289
OrderingBase(Node &node, bool container_blocked, const std::vector< std::size_t > &merge_offsets, GFSData *gfs_data, VirtualOrderingBase< DI, CI > *delegate=nullptr)
Definition: orderingbase.hh:165
bool fixedSize(typename Traits::SizeType codim) const
Definition: orderingbase.hh:215
const std::size_t _child_count
Definition: orderingbase.hh:279
DefaultLFSCacheTag CacheTag
Definition: orderingbase.hh:40
OrderingBase(Node &node, bool container_blocked, GFSData *gfs_data, VirtualOrderingBase< DI, CI > *delegate=nullptr)
Definition: orderingbase.hh:142
DI::size_type SizeType
Definition: ordering/utility.hh:201
std::vector< typename Traits::SizeType > _child_block_merge_offsets
Definition: orderingbase.hh:284
Lexicographically ordered ([i1,i2],[j1,j2] -> [i1,i2,j1,j2]).
Definition: ordering/utility.hh:25
std::vector< OrderingBase * > _children
Definition: orderingbase.hh:280
std::size_t _max_local_size
Definition: orderingbase.hh:288
Definition: ordering/utility.hh:243
Traits::SizeType maxLocalSize() const
Definition: orderingbase.hh:89
std::size_t _block_count
Definition: orderingbase.hh:290
Error related to the logical structure of an Ordering.
Definition: exceptions.hh:44
const bool _container_blocked
Definition: orderingbase.hh:276
const VirtualOrderingBase< DI, CI > * _delegate
Definition: orderingbase.hh:292
CI ContainerIndex
Definition: ordering/utility.hh:160
MergeMode::type mergeMode() const
Definition: orderingbase.hh:94
void mergePartitionSets(It begin, It end)
Adds the partitions from all PartitionInfoProviders in the range [begin,end).
Definition: partitioninfoprovider.hh:77
Traits::SizeType blockOffset(const typename Traits::SizeType child_index) const
Definition: orderingbase.hh:83
bool _fixed_size
Definition: orderingbase.hh:275
OrderingTraits< DI, CI > Traits
Definition: orderingbase.hh:30
OrderingBase & childOrdering(typename Traits::SizeType i)
Definition: orderingbase.hh:200
HierarchicContainerAllocationTag ContainerAllocationTag
Definition: orderingbase.hh:38
Definition: gridfunctionspace/tags.hh:288
std::size_t childOrderingCount() const
Definition: orderingbase.hh:195
Traits::SizeType size() const
Definition: orderingbase.hh:63
Definition: gridfunctionspace/tags.hh:280
DI::View DOFIndexView
Definition: ordering/utility.hh:198
Traits::CodimFlag _codim_used
Definition: orderingbase.hh:285
std::bitset< max_dim > CodimFlag
Definition: ordering/utility.hh:194
static const bool has_dynamic_ordering_children
Definition: orderingbase.hh:42
Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex &di) const
Definition: orderingbase.hh:46
Traits::SizeType size(const typename Traits::SizeType child_index) const
Definition: orderingbase.hh:73
std::vector< typename Traits::SizeType > _child_block_offsets
Definition: orderingbase.hh:283
Definition: adaptivity.hh:27
void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex &ci) const
Definition: orderingbase.hh:53
Dune::PDELab::impl::GridFunctionSpaceOrderingData< typename Traits::SizeType > GFSData
Definition: orderingbase.hh:34
const MergeMode::type _merge_mode
Definition: orderingbase.hh:277
Traits::SizeType blockCount() const
Definition: orderingbase.hh:68
DI DOFIndex
Definition: ordering/utility.hh:158
PDELab-specific exceptions.
Mixin class for providing information about contained grid partitions.
Definition: partitioninfoprovider.hh:22
GFSData * _gfs_data
Definition: orderingbase.hh:293
void setDelegate(const VirtualOrderingBase< DI, CI > *delegate)
Set the delegate called in mapIndex().
Definition: orderingbase.hh:228
Definition: ordering/utility.hh:230
static const bool consume_tree_index
Definition: orderingbase.hh:44
bool containerBlocked() const
Definition: orderingbase.hh:190
Definition: orderingbase.hh:21
Traits::SizeType sizeOffset(const typename Traits::SizeType child_index) const
Definition: orderingbase.hh:78
Definition: ordering/utility.hh:186
Traits::CodimFlag _codim_fixed_size
Definition: orderingbase.hh:286
bool contains(typename Traits::SizeType codim) const
Definition: orderingbase.hh:210
type
Definition: ordering/utility.hh:24
void update()
Definition: orderingbase.hh:99