Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_base.h
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/lac/block_indices.h>
24 
25 #include <vector>
26 
27 DEAL_II_NAMESPACE_OPEN
28 
34 {
84  {
105  };
106 
107 
126  inline Domination operator&(const Domination d1, const Domination d2);
127 } // namespace FiniteElementDomination
128 
129 
146 template <int dim>
148 {
149 public:
181  {
185  unknown = 0x00,
186 
190  L2 = 0x01,
191 
196  Hcurl = 0x02,
197 
202  Hdiv = 0x04,
203 
207  H1 = Hcurl | Hdiv,
208 
213  H2 = 0x0e
214  };
215 
220  static const unsigned int dimension = dim;
221 
225  const unsigned int dofs_per_vertex;
226 
231  const unsigned int dofs_per_line;
232 
237  const unsigned int dofs_per_quad;
238 
243  const unsigned int dofs_per_hex;
244 
248  const unsigned int first_line_index;
249 
253  const unsigned int first_quad_index;
254 
258  const unsigned int first_hex_index;
259 
263  const unsigned int first_face_line_index;
264 
268  const unsigned int first_face_quad_index;
269 
275  const unsigned int dofs_per_face;
276 
282  const unsigned int dofs_per_cell;
283 
292  const unsigned int components;
293 
298  const unsigned int degree;
299 
304 
311 
350  FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
351  const unsigned int n_components,
352  const unsigned int degree,
353  const Conformity conformity = unknown,
355 
359  unsigned int
360  n_dofs_per_vertex() const;
361 
365  unsigned int
366  n_dofs_per_line() const;
367 
371  unsigned int
372  n_dofs_per_quad() const;
373 
377  unsigned int
378  n_dofs_per_hex() const;
379 
384  unsigned int
385  n_dofs_per_face() const;
386 
391  unsigned int
392  n_dofs_per_cell() const;
393 
402  template <int structdim>
403  unsigned int
404  n_dofs_per_object() const;
405 
411  unsigned int
412  n_components() const;
413 
419  unsigned int
420  n_blocks() const;
421 
425  const BlockIndices &
426  block_indices() const;
427 
434  unsigned int
435  tensor_degree() const;
436 
443  bool
444  conforms(const Conformity) const;
445 
449  bool
450  operator==(const FiniteElementData &) const;
451 };
452 
453 
454 
455 // --------- inline and template functions ---------------
456 
457 
458 #ifndef DOXYGEN
459 
460 namespace FiniteElementDomination
461 {
462  inline Domination operator&(const Domination d1, const Domination d2)
463  {
464  // go through the entire list of possibilities. note that if we were into
465  // speed, obfuscation and cared enough, we could implement this operator
466  // by doing a bitwise & (and) if we gave these values to the enum values:
467  // neither_element_dominates=0, this_element_dominates=1,
468  // other_element_dominates=2, either_element_can_dominate=3
469  // =this_element_dominates|other_element_dominates
470  switch (d1)
471  {
473  if ((d2 == this_element_dominates) ||
474  (d2 == either_element_can_dominate) || (d2 == no_requirements))
475  return this_element_dominates;
476  else
478 
480  if ((d2 == other_element_dominates) ||
481  (d2 == either_element_can_dominate) || (d2 == no_requirements))
483  else
485 
488 
490  if (d2 == no_requirements)
492  else
493  return d2;
494 
495  case no_requirements:
496  return d2;
497 
498  default:
499  // shouldn't get here
500  Assert(false, ExcInternalError());
501  }
502 
504  }
505 } // namespace FiniteElementDomination
506 
507 
508 template <int dim>
509 inline unsigned int
511 {
512  return dofs_per_vertex;
513 }
514 
515 
516 
517 template <int dim>
518 inline unsigned int
520 {
521  return dofs_per_line;
522 }
523 
524 
525 
526 template <int dim>
527 inline unsigned int
529 {
530  return dofs_per_quad;
531 }
532 
533 
534 
535 template <int dim>
536 inline unsigned int
538 {
539  return dofs_per_hex;
540 }
541 
542 
543 
544 template <int dim>
545 inline unsigned int
547 {
548  return dofs_per_face;
549 }
550 
551 
552 
553 template <int dim>
554 inline unsigned int
556 {
557  return dofs_per_cell;
558 }
559 
560 
561 
562 template <int dim>
563 template <int structdim>
564 inline unsigned int
566 {
567  switch (structdim)
568  {
569  case 0:
570  return dofs_per_vertex;
571  case 1:
572  return dofs_per_line;
573  case 2:
574  return dofs_per_quad;
575  case 3:
576  return dofs_per_hex;
577  default:
578  Assert(false, ExcInternalError());
579  }
581 }
582 
583 
584 
585 template <int dim>
586 inline unsigned int
588 {
589  return components;
590 }
591 
592 
593 
594 template <int dim>
595 inline const BlockIndices &
597 {
598  return block_indices_data;
599 }
600 
601 
602 
603 template <int dim>
604 inline unsigned int
606 {
607  return block_indices_data.size();
608 }
609 
610 
611 
612 template <int dim>
613 inline unsigned int
615 {
616  return degree;
617 }
618 
619 
620 template <int dim>
621 inline bool
622 FiniteElementData<dim>::conforms(const Conformity space) const
623 {
624  return ((space & conforming_space) == space);
625 }
626 
627 
628 #endif // DOXYGEN
629 
630 
631 DEAL_II_NAMESPACE_CLOSE
632 
633 #endif
const unsigned int first_hex_index
Definition: fe_base.h:258
static const unsigned int invalid_unsigned_int
Definition: types.h:173
const unsigned int components
Definition: fe_base.h:292
unsigned int n_dofs_per_object() const
const unsigned int dofs_per_quad
Definition: fe_base.h:237
const unsigned int degree
Definition: fe_base.h:298
unsigned int n_dofs_per_vertex() const
unsigned int n_dofs_per_hex() const
const unsigned int dofs_per_line
Definition: fe_base.h:231
unsigned int n_dofs_per_face() const
unsigned int n_blocks() const
const unsigned int first_face_line_index
Definition: fe_base.h:263
bool operator==(const FiniteElementData &) const
Definition: fe_data.cc:66
const unsigned int first_quad_index
Definition: fe_base.h:253
unsigned int n_dofs_per_line() const
const unsigned int dofs_per_hex
Definition: fe_base.h:243
#define Assert(cond, exc)
Definition: exceptions.h:1407
const unsigned int dofs_per_cell
Definition: fe_base.h:282
static const unsigned int dimension
Definition: fe_base.h:220
const BlockIndices block_indices_data
Definition: fe_base.h:310
FiniteElementData(const std::vector< unsigned int > &dofs_per_object, const unsigned int n_components, const unsigned int degree, const Conformity conformity=unknown, const BlockIndices &block_indices=BlockIndices())
Definition: fe_data.cc:23
unsigned int n_dofs_per_quad() const
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
const unsigned int first_face_quad_index
Definition: fe_base.h:268
bool conforms(const Conformity) const
const Conformity conforming_space
Definition: fe_base.h:303
const unsigned int dofs_per_face
Definition: fe_base.h:275
const unsigned int first_line_index
Definition: fe_base.h:248
Domination operator &(const Domination d1, const Domination d2)
const unsigned int dofs_per_vertex
Definition: fe_base.h:225
unsigned int size() const
const BlockIndices & block_indices() const
unsigned int tensor_degree() const
static ::ExceptionBase & ExcInternalError()