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