dune-pdelab  2.0.0
leaforderingbase.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_LEAFORDERINGBASE_HH
5 #define DUNE_PDELAB_ORDERING_LEAFORDERINGBASE_HH
6 
7 #include <dune/typetree/typetree.hh>
8 
12 
13 namespace Dune {
14  namespace PDELab {
15 
18 
20  template<typename LocalOrdering>
22  : public TypeTree::VariadicCompositeNode<LocalOrdering>
23  , public VirtualOrderingBase<typename LocalOrdering::Traits::DOFIndex,
24  typename LocalOrdering::Traits::ContainerIndex>
25  , public OrderingBase<typename LocalOrdering::Traits::DOFIndex,
26  typename LocalOrdering::Traits::ContainerIndex>
27  {
28  public:
29  typedef typename LocalOrdering::Traits Traits;
30 
31  static const bool has_dynamic_ordering_children = false;
32 
33  static const bool consume_tree_index = false;
34 
35  protected:
36 
37  typedef TypeTree::VariadicCompositeNode<LocalOrdering> NodeT;
38 
39  typedef OrderingBase<typename LocalOrdering::Traits::DOFIndex,
40  typename LocalOrdering::Traits::ContainerIndex> BaseT;
41 
42  public:
43 
44  LocalOrdering& localOrdering()
45  {
46  return this->template child<0>();
47  }
48 
49  const LocalOrdering& localOrdering() const
50  {
51  return this->template child<0>();
52  }
53 
54 
55  LeafOrderingBase(const typename NodeT::NodeStorage& local_ordering, bool container_blocked, typename BaseT::GFSData* gfs_data)
56  : NodeT(local_ordering)
57  , BaseT(*this,container_blocked,gfs_data,this)
58  {
59  // copy grid partition information from local ordering.
61  }
62 
63 #ifndef DOXYGEN
64 
65 // we need to override the default copy / move ctor to fix the delegate pointer, but that is
66 // hardly interesting to our users...
67 
69  : NodeT(r.nodeStorage())
70  , BaseT(r)
72  {
73  this->setDelegate(this);
74  }
75 
76 #if HAVE_RVALUE_REFERENCES
77 
79  : NodeT(r.nodeStorage())
80  , BaseT(std::move(r))
81  , _gt_dof_offsets(std::move(r._gt_dof_offsets))
82  {
83  this->setDelegate(this);
84  }
85 
86 #endif // HAVE_RVALUE_REFERENCES
87 
88 #endif // DOXYGEN
89 
90  virtual void map_index_dynamic(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const
91  {
92  mapIndex(di,ci);
93  }
94 
95  typename Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex& di) const
96  {
97  typename Traits::ContainerIndex ci;
98  mapIndex(di.view(),ci);
99  return ci;
100  }
101 
102  void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex& ci) const
103  {
104 
105  const typename Traits::SizeType geometry_type_index = Traits::DOFIndexAccessor::geometryType(di);
106  const typename Traits::SizeType entity_index = Traits::DOFIndexAccessor::entityIndex(di);
107  assert (di.treeIndex().size() == 1);
108  ci.push_back(di.treeIndex().back());
109 
110  if (localOrdering()._fixed_size)
111  {
112  if (_container_blocked)
113  {
114  // This check is needed to avoid a horrid stream of compiler warnings about
115  // exceeding array bounds in ReservedVector!
116  if (ci.size() < ci.capacity())
117  ci.push_back(_gt_dof_offsets[geometry_type_index] + entity_index);
118  else
119  {
120  DUNE_THROW(Dune::Exception,"Container blocking incompatible with backend structure");
121  }
122  }
123  else
124  {
125  ci.back() += _gt_dof_offsets[geometry_type_index] + entity_index * localOrdering()._gt_dof_sizes[geometry_type_index];
126  }
127  }
128  else
129  {
130  if (_container_blocked)
131  {
132  // This check is needed to avoid a horrid stream of compiler warnings about
133  // exceeding array bounds in ReservedVector!
134  if (ci.size() < ci.capacity())
135  ci.push_back(localOrdering()._gt_entity_offsets[geometry_type_index] + entity_index);
136  else
137  {
138  DUNE_THROW(Dune::Exception,"Container blocking incompatible with backend structure");
139  }
140  }
141  else
142  {
143  ci.back() += localOrdering()._entity_dof_offsets[localOrdering()._gt_entity_offsets[geometry_type_index] + entity_index];
144  }
145  }
146  }
147 
148 
149  template<typename ItIn, typename ItOut>
150  void map_lfs_indices(const ItIn begin, const ItIn end, ItOut out) const
151  {
152  typedef typename Traits::SizeType size_type;
153 
155  {
156  if (_container_blocked)
157  {
158  for (ItIn in = begin; in != end; ++in, ++out)
159  {
160  assert(in->treeIndex().size() == 1);
161  out->push_back(in->treeIndex().back());
162  const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
163  const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
164  assert(localOrdering()._gt_used[geometry_type_index]);
165  out->push_back(_gt_dof_offsets[geometry_type_index] + entity_index);
166  }
167  }
168  else
169  {
170  for (ItIn in = begin; in != end; ++in, ++out)
171  {
172  assert(in->treeIndex().size() == 1);
173  out->push_back(in->treeIndex().back());
174  const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
175  const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
176  assert(localOrdering()._gt_used[geometry_type_index]);
177  out->back() += _gt_dof_offsets[geometry_type_index] + entity_index * localOrdering()._gt_dof_sizes[geometry_type_index];
178  }
179  }
180  }
181  else
182  {
183  if (_container_blocked)
184  {
185  for (ItIn in = begin; in != end; ++in, ++out)
186  {
187  assert(in->treeIndex().size() == 1);
188  out->push_back(in->treeIndex().back());
189  const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
190  const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
191  assert(localOrdering()._gt_used[geometry_type_index]);
192  out->push_back(localOrdering()._gt_entity_offsets[geometry_type_index] + entity_index);
193  }
194  }
195  else
196  {
197  for (ItIn in = begin; in != end; ++in, ++out)
198  {
199  assert(in->treeIndex().size() == 1);
200  out->push_back(in->treeIndex().back());
201  const size_type geometry_type_index = Traits::DOFIndexAccessor::geometryType(*in);
202  const size_type entity_index = Traits::DOFIndexAccessor::entityIndex(*in);
203  assert(localOrdering()._gt_used[geometry_type_index]);
204  out->back() += localOrdering()._entity_dof_offsets[localOrdering()._gt_entity_offsets[geometry_type_index] + entity_index];
205  }
206  }
207  }
208  }
209 
210  template<typename CIOutIterator>
211  typename Traits::SizeType
212  extract_entity_indices(const typename Traits::DOFIndex::EntityIndex& ei,
213  typename Traits::SizeType child_index,
214  CIOutIterator ci_out, const CIOutIterator ci_end) const
215  {
216  typedef typename Traits::SizeType size_type;
217 
218  const size_type geometry_type_index = Traits::DOFIndexAccessor::GeometryIndex::geometryType(ei);
219  const size_type entity_index = Traits::DOFIndexAccessor::GeometryIndex::entityIndex(ei);
220 
221  if (!localOrdering()._gt_used[geometry_type_index])
222  return 0;
223 
225  {
226  size_type size = localOrdering()._gt_dof_sizes[geometry_type_index];
227  if (_container_blocked)
228  {
229  for (size_type i = 0; i < size; ++i, ++ci_out)
230  {
231  ci_out->push_back(i);
232  ci_out->push_back(_gt_dof_offsets[geometry_type_index] + entity_index);
233  }
234  }
235  else
236  {
237  for (size_type i = 0; i < size; ++i, ++ci_out)
238  {
239  ci_out->push_back(i);
240  ci_out->back() += _gt_dof_offsets[geometry_type_index] + entity_index * localOrdering()._gt_dof_sizes[geometry_type_index];
241  }
242  }
243  return 0;
244  }
245  else
246  {
247  size_type index = localOrdering()._gt_entity_offsets[geometry_type_index] + entity_index;
248  size_type size = localOrdering()._entity_dof_offsets[index+1] - localOrdering()._entity_dof_offsets[index];
249  if (_container_blocked)
250  {
251  for (size_type i = 0; i < size; ++i, ++ci_out)
252  {
253  ci_out->push_back(i);
254  ci_out->push_back(index);
255  }
256  }
257  else
258  {
259  for (size_type i = 0; i < size; ++i, ++ci_out)
260  {
261  ci_out->push_back(i);
262  ci_out->back() += localOrdering()._entity_dof_offsets[index];
263  }
264  }
265  return 0;
266  }
267  }
268 
272  virtual void update() = 0;
273 
274  using BaseT::fixedSize;
275 
276  protected:
277 
279  using BaseT::_size;
280  using BaseT::_block_count;
282  using BaseT::_fixed_size;
283  using BaseT::_codim_used;
285 
286  std::vector<typename Traits::SizeType> _gt_dof_offsets;
287 
288  };
289 
291  } // namespace PDELab
292 } // namespace Dune
293 
294 #endif // DUNE_PDELAB_ORDERING_LEAFORDERING_HH
bool fixedSize(typename Traits::SizeType codim) const
Definition: orderingbase.hh:214
Traits::CodimFlag _codim_fixed_size
Definition: orderingbase.hh:285
std::size_t _max_local_size
Definition: orderingbase.hh:287
static const bool has_dynamic_ordering_children
Definition: leaforderingbase.hh:31
LocalOrdering::Traits Traits
Definition: leaforderingbase.hh:29
void setPartitionSet(const std::bitset< 6 > &partitions)
Sets the set of contained partitions to the passed-in value.
Definition: partitioninfoprovider.hh:58
Definition: orderingbase.hh:22
virtual void map_index_dynamic(typename Traits::DOFIndexView di, typename Traits::ContainerIndex &ci) const
Definition: leaforderingbase.hh:90
Definition: ordering/utility.hh:230
static const bool consume_tree_index
Definition: leaforderingbase.hh:33
Traits::SizeType extract_entity_indices(const typename Traits::DOFIndex::EntityIndex &ei, typename Traits::SizeType child_index, CIOutIterator ci_out, const CIOutIterator ci_end) const
Definition: leaforderingbase.hh:212
TypeTree::VariadicCompositeNode< LocalOrdering > NodeT
Definition: leaforderingbase.hh:37
Generic infrastructure for orderings for leaf spaces.
Definition: leaforderingbase.hh:21
std::vector< typename Traits::SizeType > _gt_dof_offsets
Definition: leaforderingbase.hh:286
LocalOrdering & localOrdering()
Definition: leaforderingbase.hh:44
Traits::CodimFlag _codim_used
Definition: orderingbase.hh:284
LeafOrderingBase(const typename NodeT::NodeStorage &local_ordering, bool container_blocked, typename BaseT::GFSData *gfs_data)
Definition: leaforderingbase.hh:55
Traits::ContainerIndex mapIndex(const typename Traits::DOFIndex &di) const
Definition: leaforderingbase.hh:95
void mapIndex(typename Traits::DOFIndexView di, typename Traits::ContainerIndex &ci) const
Definition: leaforderingbase.hh:102
std::size_t _block_count
Definition: orderingbase.hh:289
OrderingBase< typename LocalOrdering::Traits::DOFIndex, typename LocalOrdering::Traits::ContainerIndex > BaseT
Definition: leaforderingbase.hh:40
const LocalOrdering & localOrdering() const
Definition: leaforderingbase.hh:49
void map_lfs_indices(const ItIn begin, const ItIn end, ItOut out) const
Definition: leaforderingbase.hh:150
std::size_t _size
Definition: orderingbase.hh:288
Dune::PDELab::impl::GridFunctionSpaceOrderingData< typename Traits::SizeType > GFSData
Definition: orderingbase.hh:35
void setDelegate(const VirtualOrderingBase< LocalOrdering::Traits::DOFIndex, LocalOrdering::Traits::ContainerIndex > *delegate)
Set the delegate called in mapIndex().
Definition: orderingbase.hh:227