dune-pdelab  2.4.1
dofindex.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_PDELAB_COMMON_DOFINDEX_HH
4 #define DUNE_PDELAB_COMMON_DOFINDEX_HH
5 
7 #include <dune/common/array.hh>
8 
9 namespace Dune {
10 
11  namespace PDELab {
12 
14 
145  template<typename T, std::size_t tree_n, std::size_t entity_n = 1>
146  class DOFIndex
147  {
148 
149  public:
150 
152  static const std::size_t entity_capacity = entity_n;
153 
155  static const std::size_t max_depth = tree_n;
156 
157  typedef std::array<T,entity_capacity> EntityIndex;
159 
160  typedef typename TreeIndex::size_type size_type;
161  typedef T value_type;
162 
163  class View
164  {
165 
166  friend class DOFIndex;
167 
168  public:
169 
170  static const std::size_t max_depth = tree_n;
171  static const std::size_t entity_capacity = entity_n;
172 
173  typedef const std::array<T,entity_n>& EntityIndex;
175 
176  EntityIndex& entityIndex() const
177  {
178  return _entity_index_view;
179  }
180 
181  const TreeIndex& treeIndex() const
182  {
183  return _tree_index_view;
184  }
185 
187  {
188  return View(_entity_index_view,_tree_index_view.back_popped());
189  }
190 
191  friend std::ostream& operator<< (std::ostream& s, const View& di)
192  {
193  s << "(";
194 
195  for (typename std::remove_reference<EntityIndex>::type::const_iterator it = di._entity_index_view.begin(); it != di._entity_index_view.end(); ++it)
196  s << std::setw(4) << *it;
197 
198  s << " | "
199  << di._tree_index_view
200  << ")";
201  return s;
202  }
203 
204  std::size_t size() const
205  {
206  return _tree_index_view.size();
207  };
208 
209  private:
210 
211  explicit View(const DOFIndex& dof_index)
212  : _entity_index_view(dof_index._entity_index)
213  , _tree_index_view(dof_index._tree_index.view())
214  {}
215 
216  View(const DOFIndex& dof_index, std::size_t size)
217  : _entity_index_view(dof_index._entity_index)
218  , _tree_index_view(dof_index._tree_index.view(size))
219  {}
220 
221  View(EntityIndex& entity_index, const TreeIndex& tree_index)
222  : _entity_index_view(entity_index)
223  , _tree_index_view(tree_index)
224  {}
225 
226  EntityIndex _entity_index_view;
227  TreeIndex _tree_index_view;
228 
229  };
230 
233  {}
234 
235  explicit DOFIndex(const View& view)
236  : _entity_index(view._entity_index_view)
237  , _tree_index(view._tree_index_view)
238  {}
239 
240  View view() const
241  {
242  return View(*this);
243  }
244 
245  View view(std::size_t size) const
246  {
247  return View(*this,size);
248  }
249 
250  void clear()
251  {
252  std::fill(_entity_index.begin(),_entity_index.end(),0);
253  _tree_index.clear();
254  }
255 
257  EntityIndex& entityIndex()
258  {
259  return _entity_index;
260  }
261 
263  const EntityIndex& entityIndex() const
264  {
265  return _entity_index;
266  }
267 
269  TreeIndex& treeIndex()
270  {
271  return _tree_index;
272  }
273 
275  const TreeIndex& treeIndex() const
276  {
277  return _tree_index;
278  }
279 
281  friend std::ostream& operator<< (std::ostream& s, const DOFIndex& di)
282  {
283  s << "(";
284 
285  for (typename EntityIndex::const_iterator it = di._entity_index.begin(); it != di._entity_index.end(); ++it)
286  s << std::setw(4) << *it;
287 
288  s << " | "
289  << di._tree_index
290  << ")";
291  return s;
292  }
293 
295 
298  bool operator== (const DOFIndex& r) const
299  {
300  return
301  std::equal(_entity_index.begin(),_entity_index.end(),r._entity_index.begin()) &&
302  _tree_index == r._tree_index;
303  }
304 
306  bool operator!= (const DOFIndex& r) const
307  {
308  return !(*this == r);
309  }
310 
311 #if 0
312  bool operator< (const DOFIndex& r) const
313  {
314  // FIXME: think about natural ordering
315  return _c.size() < _r.size();
316  return std::lexicographical_compare(_c.begin(),_c.end(),r._c.begin(),r._c.end());
317  }
318 #endif
319 
320  std::size_t size() const
321  {
322  return _tree_index.size();
323  }
324 
325  private:
326 
327  EntityIndex _entity_index;
328  TreeIndex _tree_index;
329 
330  };
331 
332  template<typename T, std::size_t n1, std::size_t n2>
333  inline std::size_t hash_value(const DOFIndex<T,n1,n2>& di)
334  {
335  std::size_t seed = 0;
336  hash_range(seed,di.entityIndex().begin(),di.entityIndex().end());
337  hash_range(seed,di.treeIndex().begin(),di.treeIndex().end());
338  return seed;
339  }
340 
341 
342 
343  } // namespace PDELab
344 } // namespace Dune
345 
346 DUNE_DEFINE_HASH(DUNE_HASH_TEMPLATE_ARGS(typename T, std::size_t n1, std::size_t n2),DUNE_HASH_TYPE(Dune::PDELab::DOFIndex<T,n1,n2>))
347 
348 #endif // DUNE_PDELAB_COMMON_DOFINDEX_HH
View view() const
Definition: multiindex.hh:171
std::array< T, entity_capacity > EntityIndex
Definition: dofindex.hh:157
EntityIndex & entityIndex() const
Definition: dofindex.hh:176
std::size_t size() const
Definition: dofindex.hh:320
View view() const
Definition: dofindex.hh:240
const EntityIndex & entityIndex() const
Returns the index of the grid entity associated with the DOF.
Definition: dofindex.hh:263
EntityIndex & entityIndex()
Returns the index of the grid entity associated with the DOF.
Definition: dofindex.hh:257
friend std::ostream & operator<<(std::ostream &s, const View &di)
Definition: dofindex.hh:191
MultiIndex< T, max_depth > TreeIndex
Definition: dofindex.hh:158
Definition: adaptivity.hh:27
const std::array< T, entity_n > & EntityIndex
Definition: dofindex.hh:173
TreeIndex & treeIndex()
Returns the tuple of entity-local indices associated with the DOF.
Definition: dofindex.hh:269
std::size_t hash_value(const DOFIndex< T, n1, n2 > &di)
Definition: dofindex.hh:333
DOFIndex()
Default constructor.
Definition: dofindex.hh:232
MultiIndex< T, tree_n >::View TreeIndex
Definition: dofindex.hh:174
std::size_t size() const
Definition: dofindex.hh:204
bool operator!=(const DOFIndex &r) const
Tests whether two DOFIndices are not equal.
Definition: dofindex.hh:306
static const std::size_t max_depth
The maximum possible length of the tuple of entity-local indices.
Definition: dofindex.hh:155
TreeIndex::size_type size_type
Definition: dofindex.hh:160
View back_popped() const
Definition: dofindex.hh:186
View view(std::size_t size) const
Definition: dofindex.hh:245
const std::string s
Definition: function.hh:1102
const TreeIndex & treeIndex() const
Returns the tuple of entity-local indices associated with the DOF.
Definition: dofindex.hh:275
DOFIndex(const View &view)
Definition: dofindex.hh:235
static const std::size_t entity_capacity
The maximum possible length of a tuple which represents the entity index.
Definition: dofindex.hh:152
A multi-index representing a degree of freedom in a GridFunctionSpace.
Definition: dofindex.hh:146
bool operator==(const DOFIndex &r) const
Tests whether two DOFIndices are equal.
Definition: dofindex.hh:298
T value_type
Definition: dofindex.hh:161
const TreeIndex & treeIndex() const
Definition: dofindex.hh:181
void clear()
Definition: dofindex.hh:250
Definition: dofindex.hh:163