dune-pdelab 2.7-git
Loading...
Searching...
No Matches
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
6#include <array>
7
9
10namespace Dune {
11
12 namespace PDELab {
13
15
146 template<typename T, std::size_t tree_n, std::size_t entity_n = 1>
148 {
149
150 public:
151
153 static const std::size_t entity_capacity = entity_n;
154
156 static const std::size_t max_depth = tree_n;
157
158 typedef std::array<T,entity_capacity> EntityIndex;
160
161 typedef typename TreeIndex::size_type size_type;
162 typedef T value_type;
163
164 class View
165 {
166
167 friend class DOFIndex;
168
169 public:
170
171 static const std::size_t max_depth = tree_n;
172 static const std::size_t entity_capacity = entity_n;
173
174 typedef const std::array<T,entity_n>& EntityIndex;
176
178 {
179 return _entity_index_view;
180 }
181
182 const TreeIndex& treeIndex() const
183 {
184 return _tree_index_view;
185 }
186
188 {
189 return View(_entity_index_view,_tree_index_view.back_popped());
190 }
191
192 friend std::ostream& operator<< (std::ostream& s, const View& di)
193 {
194 s << "(";
195
196 for (typename std::remove_reference<EntityIndex>::type::const_iterator it = di._entity_index_view.begin(); it != di._entity_index_view.end(); ++it)
197 s << std::setw(4) << *it;
198
199 s << " | "
200 << di._tree_index_view
201 << ")";
202 return s;
203 }
204
205 std::size_t size() const
206 {
207 return _tree_index_view.size();
208 };
209
210 private:
211
212 explicit View(const DOFIndex& dof_index)
213 : _entity_index_view(dof_index._entity_index)
214 , _tree_index_view(dof_index._tree_index.view())
215 {}
216
217 View(const DOFIndex& dof_index, std::size_t size)
218 : _entity_index_view(dof_index._entity_index)
219 , _tree_index_view(dof_index._tree_index.view(size))
220 {}
221
222 View(EntityIndex& entity_index, const TreeIndex& tree_index)
223 : _entity_index_view(entity_index)
224 , _tree_index_view(tree_index)
225 {}
226
227 EntityIndex _entity_index_view;
228 TreeIndex _tree_index_view;
229
230 };
231
234 {}
235
236 explicit DOFIndex(const View& view)
237 : _entity_index(view._entity_index_view)
238 , _tree_index(view._tree_index_view)
239 {}
240
241 View view() const
242 {
243 return View(*this);
244 }
245
246 View view(std::size_t size) const
247 {
248 return View(*this,size);
249 }
250
251 void clear()
252 {
253 std::fill(_entity_index.begin(),_entity_index.end(),0);
254 _tree_index.clear();
255 }
256
259 {
260 return _entity_index;
261 }
262
265 {
266 return _entity_index;
267 }
268
271 {
272 return _tree_index;
273 }
274
276 const TreeIndex& treeIndex() const
277 {
278 return _tree_index;
279 }
280
282 friend std::ostream& operator<< (std::ostream& s, const DOFIndex& di)
283 {
284 s << "(";
285
286 for (typename EntityIndex::const_iterator it = di._entity_index.begin(); it != di._entity_index.end(); ++it)
287 s << std::setw(4) << *it;
288
289 s << " | "
290 << di._tree_index
291 << ")";
292 return s;
293 }
294
296
299 bool operator== (const DOFIndex& r) const
300 {
301 return
302 std::equal(_entity_index.begin(),_entity_index.end(),r._entity_index.begin()) &&
303 _tree_index == r._tree_index;
304 }
305
307 bool operator!= (const DOFIndex& r) const
308 {
309 return !(*this == r);
310 }
311
312#if 0
313 bool operator< (const DOFIndex& r) const
314 {
315 // FIXME: think about natural ordering
316 return _c.size() < _r.size();
317 return std::lexicographical_compare(_c.begin(),_c.end(),r._c.begin(),r._c.end());
318 }
319#endif
320
321 std::size_t size() const
322 {
323 return _tree_index.size();
324 }
325
326 private:
327
328 EntityIndex _entity_index;
329 TreeIndex _tree_index;
330
331 };
332
333 template<typename T, std::size_t n1, std::size_t n2>
334 inline std::size_t hash_value(const DOFIndex<T,n1,n2>& di)
335 {
336 std::size_t seed = 0;
337 hash_range(seed,di.entityIndex().begin(),di.entityIndex().end());
338 hash_range(seed,di.treeIndex().begin(),di.treeIndex().end());
339 return seed;
340 }
341
342
343
344 } // namespace PDELab
345} // namespace Dune
346
347DUNE_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>))
348
349#endif // DUNE_PDELAB_COMMON_DOFINDEX_HH
const std::string s
Definition: function.hh:843
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
std::size_t hash_value(const DOFIndex< T, n1, n2 > &di)
Definition: dofindex.hh:334
A multi-index representing a degree of freedom in a GridFunctionSpace.
Definition: dofindex.hh:148
friend std::ostream & operator<<(std::ostream &s, const DOFIndex &di)
Writes a pretty representation of the DOFIndex to the given std::ostream.
Definition: dofindex.hh:282
std::size_t size() const
Definition: dofindex.hh:321
View view() const
Definition: dofindex.hh:241
bool operator!=(const DOFIndex &r) const
Tests whether two DOFIndices are not equal.
Definition: dofindex.hh:307
static const std::size_t max_depth
The maximum possible length of the tuple of entity-local indices.
Definition: dofindex.hh:156
DOFIndex(const View &view)
Definition: dofindex.hh:236
std::array< T, entity_capacity > EntityIndex
Definition: dofindex.hh:158
const EntityIndex & entityIndex() const
Returns the index of the grid entity associated with the DOF.
Definition: dofindex.hh:264
View view(std::size_t size) const
Definition: dofindex.hh:246
MultiIndex< T, max_depth > TreeIndex
Definition: dofindex.hh:159
EntityIndex & entityIndex()
Returns the index of the grid entity associated with the DOF.
Definition: dofindex.hh:258
static const std::size_t entity_capacity
The maximum possible length of a tuple which represents the entity index.
Definition: dofindex.hh:153
T value_type
Definition: dofindex.hh:162
void clear()
Definition: dofindex.hh:251
const TreeIndex & treeIndex() const
Returns the tuple of entity-local indices associated with the DOF.
Definition: dofindex.hh:276
TreeIndex::size_type size_type
Definition: dofindex.hh:161
TreeIndex & treeIndex()
Returns the tuple of entity-local indices associated with the DOF.
Definition: dofindex.hh:270
DOFIndex()
Default constructor.
Definition: dofindex.hh:233
bool operator==(const DOFIndex &r) const
Tests whether two DOFIndices are equal.
Definition: dofindex.hh:299
Definition: dofindex.hh:165
friend std::ostream & operator<<(std::ostream &s, const View &di)
Definition: dofindex.hh:192
static const std::size_t max_depth
Definition: dofindex.hh:171
static const std::size_t entity_capacity
Definition: dofindex.hh:172
std::size_t size() const
Definition: dofindex.hh:205
View back_popped() const
Definition: dofindex.hh:187
MultiIndex< T, tree_n >::View TreeIndex
Definition: dofindex.hh:175
const std::array< T, entity_n > & EntityIndex
Definition: dofindex.hh:174
EntityIndex & entityIndex() const
Definition: dofindex.hh:177
const TreeIndex & treeIndex() const
Definition: dofindex.hh:182
A class for representing multi-indices.
Definition: multiindex.hh:29
Definition: multiindex.hh:39