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 
8 #include <dune/pdelab/common/multiindex.hh>
9 
10 namespace Dune {
11 
12   namespace PDELab {
13 
14     //! A multi-index representing a degree of freedom in a GridFunctionSpace.
15     /**
16      * A DOFIndex provides a way for identifying degrees of freedom in a (possibly
17      * nested) GridFunctionSpace without imposing any kind of ordering. For that
18      * purpose, a DOFIndex identifies a degree of freedom by recording
19      *
20      * - the geometry type of the grid entity associated with the DOF,
21      * - an index value uniquely identifying the grid entity among all grid entities
22      *   of its geometry type (usually just the index value of some Grid IndexSet),
23      * - a tuple of entity-local indices.
24      *
25      * The length of this index tuple is limited by the template parameter \p tree_n, which will
26      * usually be equal to the maximum depth of the current GridFunctionSpace tree. Moreover,
27      * there will never be two identical index tuples associated with the same grid entity.
28      *
29      * The index tuple is oriented from left to right when traversing up the tree (i.e.
30      * towards the root node) and from right to left when drilling down from the root node
31      * towards a leaf. For non-leaf nodes, the associated index entry identifies the child
32      * GridFunctionSpace that the degree of freedom is associated with, while for leafs, it
33      * provides a way to provide multiple degrees of freedom for a single grid entity (usually,
34      * the index value for a leaf space will correspond to the LocalKey::index() value from
35      * the finite element).
36      *
37      * Note that in general, the length of the index tuple will not be the same for all degrees
38      * of freedom in a GridFunctionSpace. Consider the following example of a Taylor-Hood element
39      * (P<sub>2</sub>-P<sub>1</sub>):
40      * \dot
41      * graph taylor_hood {
42      * node [shape=record, style=rounded, fontname=Helvetica, fontsize=8, height=0.2, width=0.4];
43      * TH [ label="Taylor-Hood"];
44      * TH -- V;
45      * V [ label="Velocity"];
46      * TH -- P;
47      * P [ label="Pressure"];
48      * V -- Vx;
49      * Vx [ label="x Velocity" ];
50      * V -- Vy;
51      * Vy [ label="y Velocity" ];
52      * }
53      * \enddot
54      * In this case, degrees of freedom for the velocity components will have an index tuple of length
55      * 3, while those related to pressure will only have an index tuple of length 2. For the Taylor-Hood
56      * space given above, the DOFIndex associated with a triangle with vertex and edge indices in the
57      * range {0,1,2} are
58      *
59      * <table>
60      *  <tr>
61      *   <th>GeometryType</th>
62      *   <th>grid entity index</th>
63      *   <th>entity-local index tuple</th>
64      *  </tr>
65      *  <tr>
66      *   <td>Point</td>
67      *   <td>0</td>
68      *   <td align="right">0, 0, 0</td>
69      *  </tr>
70      *  <tr>
71      *   <td>Point</td>
72      *   <td>1</td>
73      *   <td align="right">0, 0, 0</td>
74      *  </tr>
75      *  <tr>
76      *   <td>Point</td>
77      *   <td>2</td>
78      *   <td align="right">0, 0, 0</td>
79      *  </tr>
80      *  <tr>
81      *   <td>Line</td>
82      *   <td>0</td>
83      *   <td align="right">0, 0, 0</td>
84      *  </tr>
85      *  <tr>
86      *   <td>Line</td>
87      *   <td>1</td>
88      *   <td align="right">0, 0, 0</td>
89      *  </tr>
90      *  <tr>
91      *   <td>Line</td>
92      *   <td>2</td>
93      *   <td align="right">0, 0, 0</td>
94      *  </tr>
95      *  <tr>
96      *   <td>Point</td>
97      *   <td>0</td>
98      *   <td align="right">0, 1, 0</td>
99      *  </tr>
100      *  <tr>
101      *   <td>Point</td>
102      *   <td>1</td>
103      *   <td align="right">0, 1, 0</td>
104      *  </tr>
105      *  <tr>
106      *   <td>Point</td>
107      *   <td>2</td>
108      *   <td align="right">0, 1, 0</td>
109      *  </tr>
110      *  <tr>
111      *   <td>Line</td>
112      *   <td>0</td>
113      *   <td align="right">0, 1, 0</td>
114      *  </tr>
115      *  <tr>
116      *   <td>Line</td>
117      *   <td>1</td>
118      *   <td align="right">0, 1, 0</td>
119      *  </tr>
120      *  <tr>
121      *   <td>Line</td>
122      *   <td>2</td>
123      *   <td align="right">0, 1, 0</td>
124      *  </tr>
125      *  <tr>
126      *   <td>Point</td>
127      *   <td>0</td>
128      *   <td align="right">0, 1</td>
129      *  </tr>
130      *  <tr>
131      *   <td>Point</td>
132      *   <td>1</td>
133      *   <td align="right">0, 1</td>
134      *  </tr>
135      *  <tr>
136      *   <td>Point</td>
137      *   <td>2</td>
138      *   <td align="right">0, 1</td>
139      *  </tr>
140      * </table>
141      *
142      * \tparam T         the type of the index entries.
143      * \tparam tree_n    the maximum length of the tuple of entity-local indices.
144      * \tparam entity_n  the maximum length of a tuple which represents the entity index.
145      */
146     template<typename T, std::size_t tree_n, std::size_t entity_n = 1>
147     class DOFIndex
148     {
149 
150     public:
151 
152       //! The maximum possible length of a tuple which represents the entity index.
153       static const std::size_t entity_capacity = entity_n;
154 
155       //! The maximum possible length of the tuple of entity-local indices.
156       static const std::size_t max_depth = tree_n;
157 
158       typedef std::array<T,entity_capacity> EntityIndex;
159       typedef MultiIndex<T,max_depth> TreeIndex;
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;
175         typedef typename MultiIndex<T,tree_n>::View TreeIndex;
176 
entityIndex() const177         EntityIndex& entityIndex() const
178         {
179           return _entity_index_view;
180         }
181 
treeIndex() const182         const TreeIndex& treeIndex() const
183         {
184           return _tree_index_view;
185         }
186 
back_popped() const187         View back_popped() const
188         {
189           return View(_entity_index_view,_tree_index_view.back_popped());
190         }
191 
operator <<(std::ostream & s,const View & di)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 
size() const205         std::size_t size() const
206         {
207           return _tree_index_view.size();
208         };
209 
210       private:
211 
View(const DOFIndex & dof_index)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 
View(const DOFIndex & dof_index,std::size_t size)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 
View(EntityIndex & entity_index,const TreeIndex & tree_index)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 
232       //! Default constructor.
DOFIndex()233       DOFIndex()
234       {}
235 
DOFIndex(const View & view)236       explicit DOFIndex(const View& view)
237         : _entity_index(view._entity_index_view)
238         , _tree_index(view._tree_index_view)
239       {}
240 
view() const241       View view() const
242       {
243         return View(*this);
244       }
245 
view(std::size_t size) const246       View view(std::size_t size) const
247       {
248         return View(*this,size);
249       }
250 
clear()251       void clear()
252       {
253         std::fill(_entity_index.begin(),_entity_index.end(),0);
254         _tree_index.clear();
255       }
256 
257       //! Returns the index of the grid entity associated with the DOF.
entityIndex()258       EntityIndex& entityIndex()
259       {
260         return _entity_index;
261       }
262 
263       //! Returns the index of the grid entity associated with the DOF.
entityIndex() const264       const EntityIndex& entityIndex() const
265       {
266         return _entity_index;
267       }
268 
269       //! Returns the tuple of entity-local indices associated with the DOF.
treeIndex()270       TreeIndex& treeIndex()
271       {
272         return _tree_index;
273       }
274 
275       //! Returns the tuple of entity-local indices associated with the DOF.
treeIndex() const276       const TreeIndex& treeIndex() const
277       {
278         return _tree_index;
279       }
280 
281       //! Writes a pretty representation of the DOFIndex to the given std::ostream.
operator <<(std::ostream & s,const DOFIndex & di)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 
295       //! Tests whether two DOFIndices are equal.
296       /**
297        * \note Only DOFIndices of identical max_depth are comparable.
298        */
operator ==(const DOFIndex & r) const299       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 
306       //! Tests whether two DOFIndices are not equal.
operator !=(const DOFIndex & r) const307       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 
size() const321       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>
hash_value(const DOFIndex<T,n1,n2> & di)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 
347 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>))
348 
349 #endif // DUNE_PDELAB_COMMON_DOFINDEX_HH
350