1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15
16 #ifndef dealii_fe_base_h
17 #define dealii_fe_base_h
18
19 #include <deal.II/base/config.h>
20
21 #include <deal.II/base/exceptions.h>
22
23 #include <deal.II/grid/reference_cell.h>
24
25 #include <deal.II/lac/block_indices.h>
26
27 #include <vector>
28
29 DEAL_II_NAMESPACE_OPEN
30
31 /**
32 * A namespace solely for the purpose of defining the Domination enum as well
33 * as associated operators.
34 */
35 namespace FiniteElementDomination
36 {
37 /**
38 * An enum that describes the outcome of comparing two elements for mutual
39 * domination. If one element dominates another, then the restriction of the
40 * space described by the dominated element to a face of the cell is
41 * strictly larger than that of the dominating element. For example, in 2-d
42 * Q(2) elements dominate Q(4) elements, because the traces of Q(4) elements
43 * are quartic polynomials which is a space strictly larger than the
44 * quadratic polynomials (the restriction of the Q(2) element). Similar
45 * reasonings apply for vertices and cells as well. In general, Q(k) dominates
46 * Q(k') if $k\le k'$.
47 *
48 * This enum is used in the FiniteElement::compare_for_domination() function
49 * that is used in the context of hp finite element methods when determining
50 * what to do at faces where two different finite elements meet (see the
51 * @ref hp_paper "hp paper"
52 * for a more detailed description of the following). In that case, the
53 * degrees of freedom of one side need to be constrained to those on the
54 * other side. The determination which side is which is based on the outcome
55 * of a comparison for mutual domination: the dominated side is constrained
56 * to the dominating one.
57 *
58 * Note that there are situations where neither side dominates. The
59 * @ref hp_paper "hp paper"
60 * lists two case, with the simpler one being that a $Q_2\times Q_1$ vector-
61 * valued element (i.e. a <code>FESystem(FE_Q(2),1,FE_Q(1),1)</code>) meets
62 * a $Q_1\times Q_2$ element: here, for each of the two vector-components,
63 * we can define a domination relationship, but it is different for the two
64 * components.
65 *
66 * It is clear that the concept of domination doesn't matter for
67 * discontinuous elements. However, discontinuous elements may be part of
68 * vector-valued elements and may therefore be compared against each other
69 * for domination. They should return
70 * <code>either_element_can_dominate</code> in that case. Likewise, when
71 * comparing two identical finite elements, they should return this code;
72 * the reason is that we can not decide which element will dominate at the
73 * time we look at the first component of, for example, two $Q_2\times Q_1$
74 * and $Q_2\times Q_2$ elements, and have to keep our options open until we
75 * get to the second base element.
76 *
77 * Finally, the code no_requirements exists for cases where elements impose
78 * no continuity requirements. The case is primarily meant for FE_Nothing
79 * which is an element that has no degrees of freedom in a subdomain. It
80 * could also be used by discontinuous elements, for example.
81 *
82 * More details on domination can be found in the
83 * @ref hp_paper "hp paper".
84 */
85 enum Domination
86 {
87 /**
88 * The current element dominates.
89 */
90 this_element_dominates,
91 /**
92 * The other element dominates.
93 */
94 other_element_dominates,
95 /**
96 * Neither element dominates.
97 */
98 neither_element_dominates,
99 /**
100 * Either element may dominate.
101 */
102 either_element_can_dominate,
103 /**
104 * There are no requirements.
105 */
106 no_requirements
107 };
108
109
110 /**
111 * A generalization of the binary <code>and</code> operator to a comparison
112 * relationship. The way this works is pretty much as when you would want to
113 * define a comparison relationship for vectors: either all elements of the
114 * first vector are smaller, equal, or larger than those of the second
115 * vector, or some are and some are not.
116 *
117 * This operator is pretty much the same: if both arguments are
118 * <code>this_element_dominates</code> or
119 * <code>other_element_dominates</code>, then the returned value is that
120 * value. On the other hand, if one of the values is
121 * <code>either_element_can_dominate</code>, then the returned value is that
122 * of the other argument. If either argument is
123 * <code>neither_element_dominates</code>, or if the two arguments are
124 * <code>this_element_dominates</code> and
125 * <code>other_element_dominates</code>, then the returned value is
126 * <code>neither_element_dominates</code>.
127 */
128 inline Domination operator&(const Domination d1, const Domination d2);
129 } // namespace FiniteElementDomination
130
131 namespace internal
132 {
133 /**
134 * Internal data structure for setting up FiniteElementData. It stores for
135 * each object the (inclusive/exclusive) number of degrees of freedoms, as
136 * well as, the index of its first degree of freedom within a cell and the
137 * index of the first d-dimensional object within each face.
138 *
139 * The information is saved as a vector of vectors. One can query the
140 * inclusive number of dofs of the i-th d-dimensional object via:
141 * dofs_per_object_inclusive[d][i].
142 *
143 * As an example, the data is shown for a quadratic wedge. Which consists of
144 * 6 vertices, 9 lines, and 5 faces (two triangles and three quadrilaterals).
145 *
146 * vertices lines faces cell
147 * dpo_excl 1 1 1 1 1 1 | 1 1 1 1 1 1 1 1 1 | 0 0 1 1 1 | 0
148 * dpo_incl 1 1 1 1 1 1 | 3 3 3 3 3 3 3 3 3 | 6 6 9 9 9 | 18
149 * obj_index 0 1 2 3 4 5 | 6 7 8 9 10 11 12 13 14 | 15 15 15 16 17 | 18
150 *
151 * Since the above table looks as follows for:
152 *
153 * - a triangle:
154 *
155 * dpo_excl 1 1 1 | 1 1 1 | 0
156 * obj_index 0 1 2 | 3 4 5 | 6
157 *
158 * - quadrilateral:
159 *
160 * dpo_excl 1 1 1 1 | 1 1 1 1 | 1
161 * obj_index 0 1 2 3 | 4 5 6 7 | 8
162 *
163 * The index of the first d-dimensional object within each face results as:
164 *
165 * vertices lines face
166 * first_obj_index_on_face 0 0 0 0 0 | 3 3 4 4 4 | 6 6 8 8 8
167 *
168 */
169 struct GenericDoFsPerObject
170 {
171 /**
172 * Exclusive number of degrees of freedom per object.
173 */
174 std::vector<std::vector<unsigned int>> dofs_per_object_exclusive;
175
176 /**
177 * Inclusive number of degrees of freedom per object.
178 */
179 std::vector<std::vector<unsigned int>> dofs_per_object_inclusive;
180
181 /**
182 * First index of an object.
183 */
184 std::vector<std::vector<unsigned int>> object_index;
185
186 /**
187 * First index of an object within a face.
188 */
189 std::vector<std::vector<unsigned int>> first_object_index_on_face;
190 };
191 } // namespace internal
192
193 /**
194 * A class that declares a number of scalar constant variables that describe
195 * basic properties of a finite element implementation. This includes, for
196 * example, the number of degrees of freedom per vertex, line, or cell; the
197 * number of vector components; etc.
198 *
199 * The kind of information stored here is computed during initialization of a
200 * finite element object and is passed down to this class via its constructor.
201 * The data stored by this class is part of the public interface of the
202 * FiniteElement class (which derives from the current class). See there for
203 * more information.
204 *
205 * @ingroup febase
206 */
207 template <int dim>
208 class FiniteElementData
209 {
210 public:
211 /**
212 * Enumerator for the different types of continuity a finite element may
213 * have. Continuity is measured by the Sobolev space containing the
214 * constructed finite element space and is also called this way.
215 *
216 * Note that certain continuities may imply others. For instance, a function
217 * in <i>H<sup>1</sup></i> is in <i>H<sup>curl</sup></i> and
218 * <i>H<sup>div</sup></i> as well.
219 *
220 * If you are interested in continuity in the classical sense, then the
221 * following relations hold:
222 *
223 * <ol>
224 *
225 * <li> <i>H<sup>1</sup></i> implies that the function is continuous over
226 * cell boundaries.
227 *
228 * <li> <i>H<sup>2</sup></i> implies that the function is continuously
229 * differentiable over cell boundaries.
230 *
231 * <li> <i>L<sup>2</sup></i> indicates that the element is discontinuous.
232 * Since discontinuous elements have no topological couplings between grid
233 * cells and code may actually depend on this property, <i>L<sup>2</sup></i>
234 * conformity is handled in a special way in the sense that it is <b>not</b>
235 * implied by any higher conformity.
236 * </ol>
237 *
238 * In order to test if a finite element conforms to a certain space, use
239 * FiniteElementData<dim>::conforms().
240 */
241 enum Conformity
242 {
243 /**
244 * Indicates incompatible continuities of a system.
245 */
246 unknown = 0x00,
247
248 /**
249 * Discontinuous elements. See above!
250 */
251 L2 = 0x01,
252
253 /**
254 * Conformity with the space <i>H<sup>curl</sup></i> (continuous
255 * tangential component of a vector field)
256 */
257 Hcurl = 0x02,
258
259 /**
260 * Conformity with the space <i>H<sup>div</sup></i> (continuous normal
261 * component of a vector field)
262 */
263 Hdiv = 0x04,
264
265 /**
266 * Conformity with the space <i>H<sup>1</sup></i> (continuous)
267 */
268 H1 = Hcurl | Hdiv,
269
270 /**
271 * Conformity with the space <i>H<sup>2</sup></i> (continuously
272 * differentiable)
273 */
274 H2 = 0x0e
275 };
276
277 /**
278 * The dimension of the finite element, which is the template parameter
279 * <tt>dim</tt>
280 */
281 static const unsigned int dimension = dim;
282
283 private:
284 /**
285 * Reference cell type.
286 */
287 const ReferenceCell::Type cell_type;
288
289 /**
290 * Number of unique quads. If all quads have the same type, the value is
291 * one; else it equals the number of quads.
292 */
293 const unsigned int number_unique_quads;
294
295 /**
296 * Number of unique faces. If all faces have the same type, the value is
297 * one; else it equals the number of faces.
298 */
299 const unsigned int number_unique_faces;
300
301 public:
302 /**
303 * Number of degrees of freedom on a vertex.
304 */
305 const unsigned int dofs_per_vertex;
306
307 /**
308 * Number of degrees of freedom in a line; not including the degrees of
309 * freedom on the vertices of the line.
310 */
311 const unsigned int dofs_per_line;
312
313 private:
314 /**
315 * Number of degrees of freedom on quads. If all quads have the same
316 * number of degrees freedoms the values equal dofs_per_quad.
317 */
318 const std::vector<unsigned int> n_dofs_on_quad;
319
320 public:
321 /**
322 * Number of degrees of freedom in a quadrilateral; not including the
323 * degrees of freedom on the lines and vertices of the quadrilateral.
324 */
325 const unsigned int dofs_per_quad;
326
327 private:
328 /**
329 * Maximum number of degrees of freedom on any quad.
330 */
331 const unsigned int dofs_per_quad_max;
332
333 public:
334 /**
335 * Number of degrees of freedom in a hexahedron; not including the degrees
336 * of freedom on the quadrilaterals, lines and vertices of the hexahedron.
337 */
338 const unsigned int dofs_per_hex;
339
340 /**
341 * First index of dof on a line.
342 */
343 const unsigned int first_line_index;
344
345 private:
346 /**
347 * First index of a quad. If all quads have the same number of degrees of
348 * freedom, only the first index of the first quad is stored since the
349 * indices of the others can be simply recomputed.
350 */
351 const std::vector<unsigned int> first_index_of_quads;
352
353 public:
354 /**
355 * First index of dof on a quad.
356 */
357 const unsigned int first_quad_index;
358
359 /**
360 * First index of dof on a hexahedron.
361 */
362 const unsigned int first_hex_index;
363
364 private:
365 /**
366 * Index of the first line of all faces.
367 */
368 const std::vector<unsigned int> first_line_index_of_faces;
369
370 public:
371 /**
372 * First index of dof on a line for face data.
373 */
374 const unsigned int first_face_line_index;
375
376 private:
377 /**
378 * Index of the first quad of all faces.
379 */
380 const std::vector<unsigned int> first_quad_index_of_faces;
381
382 public:
383 /**
384 * First index of dof on a quad for face data.
385 */
386 const unsigned int first_face_quad_index;
387
388 private:
389 /**
390 * Number of degrees of freedom on faces. If all faces have the same
391 * number of degrees freedoms the values equal dofs_per_quad.
392 */
393 const std::vector<unsigned int> n_dofs_on_face;
394
395 public:
396 /**
397 * Number of degrees of freedom on a face. This is the accumulated number of
398 * degrees of freedom on all the objects of dimension up to <tt>dim-1</tt>
399 * constituting a face.
400 */
401 const unsigned int dofs_per_face;
402
403 private:
404 /**
405 * Maximum number of degrees of freedom on any face.
406 */
407 const unsigned int dofs_per_face_max;
408
409 public:
410 /**
411 * Total number of degrees of freedom on a cell. This is the accumulated
412 * number of degrees of freedom on all the objects of dimension up to
413 * <tt>dim</tt> constituting a cell.
414 */
415 const unsigned int dofs_per_cell;
416
417 /**
418 * Number of vector components of this finite element, and dimension of the
419 * image space. For vector-valued finite elements (i.e. when this number is
420 * greater than one), the number of vector components is in many cases equal
421 * to the number of base elements glued together with the help of the
422 * FESystem class. However, for elements like the Nedelec element, the
423 * number is greater than one even though we only have one base element.
424 */
425 const unsigned int components;
426
427 /**
428 * Maximal polynomial degree of a shape function in a single coordinate
429 * direction.
430 */
431 const unsigned int degree;
432
433 /**
434 * Indicate the space this element conforms to.
435 */
436 const Conformity conforming_space;
437
438 /**
439 * Storage for an object describing the sizes of each block of a compound
440 * element. For an element which is not an FESystem, this contains only a
441 * single block with length #dofs_per_cell.
442 */
443 const BlockIndices block_indices_data;
444
445 /**
446 * Constructor, computing all necessary values from the distribution of dofs
447 * to geometrical objects.
448 *
449 * @param[in] dofs_per_object A vector that describes the number of degrees
450 * of freedom on geometrical objects for each dimension. This vector must
451 * have size dim+1, and entry 0 describes the number of degrees of freedom
452 * per vertex, entry 1 the number of degrees of freedom per line, etc. As an
453 * example, for the common $Q_1$ Lagrange element in 2d, this vector would
454 * have elements <code>(1,0,0)</code>. On the other hand, for a $Q_3$
455 * element in 3d, it would have entries <code>(1,2,4,8)</code>.
456 *
457 * @param[in] n_components Number of vector components of the element.
458 *
459 * @param[in] degree The maximal polynomial degree of any of the shape
460 * functions of this element in any variable on the reference element. For
461 * example, for the $Q_1$ element (in any space dimension), this would be
462 * one; this is so despite the fact that the element has a shape function of
463 * the form $\hat x\hat y$ (in 2d) and $\hat x\hat y\hat z$ (in 3d), which,
464 * although quadratic and cubic polynomials, are still only linear in each
465 * reference variable separately. The information provided by this variable
466 * is typically used in determining what an appropriate quadrature formula
467 * is.
468 *
469 * @param[in] conformity A variable describing which Sobolev space this
470 * element conforms to. For example, the $Q_p$ Lagrange elements
471 * (implemented by the FE_Q class) are $H^1$ conforming, whereas the
472 * Raviart-Thomas element (implemented by the FE_RaviartThomas class) is
473 * $H_\text{div}$ conforming; finally, completely discontinuous elements
474 * (implemented by the FE_DGQ class) are only $L_2$ conforming.
475 *
476 * @param[in] block_indices An argument that describes how the base elements
477 * of a finite element are grouped. The default value constructs a single
478 * block that consists of all @p dofs_per_cell degrees of freedom. This is
479 * appropriate for all "atomic" elements (including non-primitive ones) and
480 * these can therefore omit this argument. On the other hand, composed
481 * elements such as FESystem will want to pass a different value here.
482 */
483 FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
484 const unsigned int n_components,
485 const unsigned int degree,
486 const Conformity conformity = unknown,
487 const BlockIndices &block_indices = BlockIndices());
488
489 /**
490 * The same as above but with the difference that also the type of the
491 * underlying geometric entity can be specified.
492 */
493 FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
494 const ReferenceCell::Type cell_type,
495 const unsigned int n_components,
496 const unsigned int degree,
497 const Conformity conformity = unknown,
498 const BlockIndices &block_indices = BlockIndices());
499
500 /**
501 * The same as above but instead of passing a vector containing the degrees
502 * of freedoms per object a struct of type GenericDoFsPerObject. This allows
503 * that 2D objects might have different number of degrees of freedoms, which
504 * is particular useful for cells with triangles and quadrilaterals as faces.
505 */
506 FiniteElementData(const internal::GenericDoFsPerObject &data,
507 const ReferenceCell::Type cell_type,
508 const unsigned int n_components,
509 const unsigned int degree,
510 const Conformity conformity = unknown,
511 const BlockIndices &block_indices = BlockIndices());
512
513 /**
514 * Return type of reference cell.
515 */
516 ReferenceCell::Type
517 reference_cell_type() const;
518
519 /**
520 * Number of unique quads. If all quads have the same type, the value is
521 * one; else it equals the number of quads.
522 */
523 unsigned int
524 n_unique_quads() const;
525
526 /**
527 * Number of unique faces. If all faces have the same type, the value is
528 * one; else it equals the number of faces.
529 */
530 unsigned int
531 n_unique_faces() const;
532
533 /**
534 * Number of dofs per vertex.
535 */
536 unsigned int
537 n_dofs_per_vertex() const;
538
539 /**
540 * Number of dofs per line. Not including dofs on lower dimensional objects.
541 */
542 unsigned int
543 n_dofs_per_line() const;
544
545 /**
546 * Number of dofs per quad. Not including dofs on lower dimensional objects.
547 */
548 unsigned int
549 n_dofs_per_quad(unsigned int face_no = 0) const;
550
551 /**
552 * Maximum number of dofs per quad. Not including dofs on lower dimensional
553 * objects.
554 */
555 unsigned int
556 max_dofs_per_quad() const;
557
558 /**
559 * Number of dofs per hex. Not including dofs on lower dimensional objects.
560 */
561 unsigned int
562 n_dofs_per_hex() const;
563
564 /**
565 * Number of dofs per face, accumulating degrees of freedom of all lower
566 * dimensional objects.
567 */
568 unsigned int
569 n_dofs_per_face(unsigned int face_no = 0, unsigned int child = 0) const;
570
571 /**
572 * Maximum number of dofs per face, accumulating degrees of freedom of all
573 * lower dimensional objects.
574 */
575 unsigned int
576 max_dofs_per_face() const;
577
578 /**
579 * Number of dofs per cell, accumulating degrees of freedom of all lower
580 * dimensional objects.
581 */
582 unsigned int
583 n_dofs_per_cell() const;
584
585 /**
586 * Return the number of degrees per structdim-dimensional object. For
587 * structdim==0, the function therefore returns dofs_per_vertex, for
588 * structdim==1 dofs_per_line, etc. This function is mostly used to allow
589 * some template trickery for functions that should work on all sorts of
590 * objects without wanting to use the different names (vertex, line, ...)
591 * associated with these objects.
592 */
593 template <int structdim>
594 unsigned int
595 n_dofs_per_object(const unsigned int i = 0) const;
596
597 /**
598 * Number of components. See
599 * @ref GlossComponent "the glossary"
600 * for more information.
601 */
602 unsigned int
603 n_components() const;
604
605 /**
606 * Number of blocks. See
607 * @ref GlossBlock "the glossary"
608 * for more information.
609 */
610 unsigned int
611 n_blocks() const;
612
613 /**
614 * Detailed information on block sizes.
615 */
616 const BlockIndices &
617 block_indices() const;
618
619 /**
620 * Maximal polynomial degree of a shape function in a single coordinate
621 * direction.
622 *
623 * This function can be used to determine the optimal quadrature rule.
624 */
625 unsigned int
626 tensor_degree() const;
627
628 /**
629 * Test whether a finite element space conforms to a certain Sobolev space.
630 *
631 * @note This function will return a true value even if the finite element
632 * space has higher regularity than asked for.
633 */
634 bool
635 conforms(const Conformity) const;
636
637 /**
638 * Comparison operator.
639 */
640 bool
641 operator==(const FiniteElementData &) const;
642
643 /**
644 * Return first index of dof on a line.
645 */
646 unsigned int
647 get_first_line_index() const;
648
649 /**
650 * Return first index of dof on a quad.
651 */
652 unsigned int
653 get_first_quad_index(const unsigned int quad_no = 0) const;
654
655 /**
656 * Return first index of dof on a hexahedron.
657 */
658 unsigned int
659 get_first_hex_index() const;
660
661 /**
662 * Return first index of dof on a line for face data.
663 */
664 unsigned int
665 get_first_face_line_index(const unsigned int face_no = 0) const;
666
667 /**
668 * Return first index of dof on a quad for face data.
669 */
670 unsigned int
671 get_first_face_quad_index(const unsigned int face_no = 0) const;
672 };
673
674
675
676 // --------- inline and template functions ---------------
677
678
679 #ifndef DOXYGEN
680
681 namespace FiniteElementDomination
682 {
683 inline Domination operator&(const Domination d1, const Domination d2)
684 {
685 // go through the entire list of possibilities. note that if we were into
686 // speed, obfuscation and cared enough, we could implement this operator
687 // by doing a bitwise & (and) if we gave these values to the enum values:
688 // neither_element_dominates=0, this_element_dominates=1,
689 // other_element_dominates=2, either_element_can_dominate=3
690 // =this_element_dominates|other_element_dominates
691 switch (d1)
692 {
693 case this_element_dominates:
694 if ((d2 == this_element_dominates) ||
695 (d2 == either_element_can_dominate) || (d2 == no_requirements))
696 return this_element_dominates;
697 else
698 return neither_element_dominates;
699
700 case other_element_dominates:
701 if ((d2 == other_element_dominates) ||
702 (d2 == either_element_can_dominate) || (d2 == no_requirements))
703 return other_element_dominates;
704 else
705 return neither_element_dominates;
706
707 case neither_element_dominates:
708 return neither_element_dominates;
709
710 case either_element_can_dominate:
711 if (d2 == no_requirements)
712 return either_element_can_dominate;
713 else
714 return d2;
715
716 case no_requirements:
717 return d2;
718
719 default:
720 // shouldn't get here
721 Assert(false, ExcInternalError());
722 }
723
724 return neither_element_dominates;
725 }
726 } // namespace FiniteElementDomination
727
728
729 template <int dim>
730 inline ReferenceCell::Type
reference_cell_type()731 FiniteElementData<dim>::reference_cell_type() const
732 {
733 return cell_type;
734 }
735
736
737
738 template <int dim>
739 inline unsigned int
n_unique_quads()740 FiniteElementData<dim>::n_unique_quads() const
741 {
742 return number_unique_quads;
743 }
744
745
746
747 template <int dim>
748 inline unsigned int
n_unique_faces()749 FiniteElementData<dim>::n_unique_faces() const
750 {
751 return number_unique_faces;
752 }
753
754
755
756 template <int dim>
757 inline unsigned int
n_dofs_per_vertex()758 FiniteElementData<dim>::n_dofs_per_vertex() const
759 {
760 return dofs_per_vertex;
761 }
762
763
764
765 template <int dim>
766 inline unsigned int
n_dofs_per_line()767 FiniteElementData<dim>::n_dofs_per_line() const
768 {
769 return dofs_per_line;
770 }
771
772
773
774 template <int dim>
775 inline unsigned int
n_dofs_per_quad(unsigned int face_no)776 FiniteElementData<dim>::n_dofs_per_quad(unsigned int face_no) const
777 {
778 return n_dofs_on_quad[n_dofs_on_quad.size() == 1 ? 0 : face_no];
779 }
780
781
782
783 template <int dim>
784 inline unsigned int
max_dofs_per_quad()785 FiniteElementData<dim>::max_dofs_per_quad() const
786 {
787 return dofs_per_quad_max;
788 }
789
790
791
792 template <int dim>
793 inline unsigned int
n_dofs_per_hex()794 FiniteElementData<dim>::n_dofs_per_hex() const
795 {
796 return dofs_per_hex;
797 }
798
799
800
801 template <int dim>
802 inline unsigned int
n_dofs_per_face(unsigned int face_no,unsigned int child_no)803 FiniteElementData<dim>::n_dofs_per_face(unsigned int face_no,
804 unsigned int child_no) const
805 {
806 (void)child_no;
807
808 return n_dofs_on_face[n_dofs_on_face.size() == 1 ? 0 : face_no];
809 }
810
811
812
813 template <int dim>
814 inline unsigned int
max_dofs_per_face()815 FiniteElementData<dim>::max_dofs_per_face() const
816 {
817 return dofs_per_face_max;
818 }
819
820
821
822 template <int dim>
823 inline unsigned int
n_dofs_per_cell()824 FiniteElementData<dim>::n_dofs_per_cell() const
825 {
826 return dofs_per_cell;
827 }
828
829
830
831 template <int dim>
832 template <int structdim>
833 inline unsigned int
n_dofs_per_object(const unsigned int i)834 FiniteElementData<dim>::n_dofs_per_object(const unsigned int i) const
835 {
836 switch (structdim)
837 {
838 case 0:
839 return n_dofs_per_vertex();
840 case 1:
841 return n_dofs_per_line();
842 case 2:
843 return n_dofs_per_quad((structdim == 2 && dim == 3) ? i : 0);
844 case 3:
845 return n_dofs_per_hex();
846 default:
847 Assert(false, ExcInternalError());
848 }
849 return numbers::invalid_unsigned_int;
850 }
851
852
853
854 template <int dim>
855 inline unsigned int
n_components()856 FiniteElementData<dim>::n_components() const
857 {
858 return components;
859 }
860
861
862
863 template <int dim>
864 inline const BlockIndices &
block_indices()865 FiniteElementData<dim>::block_indices() const
866 {
867 return block_indices_data;
868 }
869
870
871
872 template <int dim>
873 inline unsigned int
n_blocks()874 FiniteElementData<dim>::n_blocks() const
875 {
876 return block_indices_data.size();
877 }
878
879
880
881 template <int dim>
882 inline unsigned int
tensor_degree()883 FiniteElementData<dim>::tensor_degree() const
884 {
885 return degree;
886 }
887
888
889 template <int dim>
890 inline bool
conforms(const Conformity space)891 FiniteElementData<dim>::conforms(const Conformity space) const
892 {
893 return ((space & conforming_space) == space);
894 }
895
896
897
898 template <int dim>
899 unsigned int
get_first_line_index()900 FiniteElementData<dim>::get_first_line_index() const
901 {
902 return first_line_index;
903 }
904
905 template <int dim>
906 unsigned int
get_first_quad_index(const unsigned int quad_no)907 FiniteElementData<dim>::get_first_quad_index(const unsigned int quad_no) const
908 {
909 if (first_index_of_quads.size() == 1)
910 return first_index_of_quads[0] + quad_no * n_dofs_per_quad(0);
911 else
912 return first_index_of_quads[quad_no];
913 }
914
915 template <int dim>
916 unsigned int
get_first_hex_index()917 FiniteElementData<dim>::get_first_hex_index() const
918 {
919 return first_hex_index;
920 }
921
922 template <int dim>
923 unsigned int
get_first_face_line_index(const unsigned int face_no)924 FiniteElementData<dim>::get_first_face_line_index(
925 const unsigned int face_no) const
926 {
927 return first_line_index_of_faces[first_line_index_of_faces.size() == 1 ?
928 0 :
929 face_no];
930 }
931
932 template <int dim>
933 unsigned int
get_first_face_quad_index(const unsigned int face_no)934 FiniteElementData<dim>::get_first_face_quad_index(
935 const unsigned int face_no) const
936 {
937 return first_quad_index_of_faces[first_quad_index_of_faces.size() == 1 ?
938 0 :
939 face_no];
940 }
941
942
943 #endif // DOXYGEN
944
945
946 DEAL_II_NAMESPACE_CLOSE
947
948 #endif
949