Reference documentation for deal.II version Git 497f915867 2021-09-17 22:46:48 +0200
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
fe_base.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2021 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 
22 
24 
26 
27 #include <vector>
28 
30 
36 {
86  {
107  };
108 
109 
128  inline Domination
129  operator&(const Domination d1, const Domination d2);
130 } // namespace FiniteElementDomination
131 
132 namespace internal
133 {
175  {
179  std::vector<std::vector<unsigned int>> dofs_per_object_exclusive;
180 
184  std::vector<std::vector<unsigned int>> dofs_per_object_inclusive;
185 
189  std::vector<std::vector<unsigned int>> object_index;
190 
194  std::vector<std::vector<unsigned int>> first_object_index_on_face;
195  };
196 } // namespace internal
197 
212 template <int dim>
214 {
215 public:
247  {
251  unknown = 0x00,
252 
256  L2 = 0x01,
257 
262  Hcurl = 0x02,
263 
268  Hdiv = 0x04,
269 
273  H1 = Hcurl | Hdiv,
274 
279  H2 = 0x0e
280  };
281 
286  static const unsigned int dimension = dim;
287 
288 private:
293 
298  const unsigned int number_unique_quads;
299 
304  const unsigned int number_unique_faces;
305 
306 public:
310  const unsigned int dofs_per_vertex;
311 
316  const unsigned int dofs_per_line;
317 
318 private:
323  const std::vector<unsigned int> n_dofs_on_quad;
324 
325 public:
330  const unsigned int dofs_per_quad;
331 
332 private:
336  const unsigned int dofs_per_quad_max;
337 
338 public:
343  const unsigned int dofs_per_hex;
344 
348  const unsigned int first_line_index;
349 
350 private:
356  const std::vector<unsigned int> first_index_of_quads;
357 
358 public:
362  const unsigned int first_quad_index;
363 
367  const unsigned int first_hex_index;
368 
369 private:
373  const std::vector<unsigned int> first_line_index_of_faces;
374 
375 public:
379  const unsigned int first_face_line_index;
380 
381 private:
385  const std::vector<unsigned int> first_quad_index_of_faces;
386 
387 public:
391  const unsigned int first_face_quad_index;
392 
393 private:
398  const std::vector<unsigned int> n_dofs_on_face;
399 
400 public:
406  const unsigned int dofs_per_face;
407 
408 private:
412  const unsigned int dofs_per_face_max;
413 
414 public:
420  const unsigned int dofs_per_cell;
421 
430  const unsigned int components;
431 
436  const unsigned int degree;
437 
442 
449 
488  FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
489  const unsigned int n_components,
490  const unsigned int degree,
491  const Conformity conformity = unknown,
492  const BlockIndices &block_indices = BlockIndices());
493 
498  FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
500  const unsigned int n_components,
501  const unsigned int degree,
502  const Conformity conformity = unknown,
503  const BlockIndices &block_indices = BlockIndices());
504 
512  const ReferenceCell reference_cell,
513  const unsigned int n_components,
514  const unsigned int degree,
515  const Conformity conformity = unknown,
516  const BlockIndices &block_indices = BlockIndices());
517 
524  reference_cell() const;
525 
530  unsigned int
531  n_unique_quads() const;
532 
537  unsigned int
538  n_unique_faces() const;
539 
543  unsigned int
544  n_dofs_per_vertex() const;
545 
549  unsigned int
550  n_dofs_per_line() const;
551 
555  unsigned int
556  n_dofs_per_quad(unsigned int face_no = 0) const;
557 
562  unsigned int
563  max_dofs_per_quad() const;
564 
568  unsigned int
569  n_dofs_per_hex() const;
570 
575  unsigned int
576  n_dofs_per_face(unsigned int face_no = 0, unsigned int child = 0) const;
577 
582  unsigned int
583  max_dofs_per_face() const;
584 
589  unsigned int
590  n_dofs_per_cell() const;
591 
600  template <int structdim>
601  unsigned int
602  n_dofs_per_object(const unsigned int i = 0) const;
603 
609  unsigned int
610  n_components() const;
611 
617  unsigned int
618  n_blocks() const;
619 
623  const BlockIndices &
624  block_indices() const;
625 
632  unsigned int
633  tensor_degree() const;
634 
641  bool
642  conforms(const Conformity) const;
643 
647  bool
648  operator==(const FiniteElementData &) const;
649 
653  unsigned int
654  get_first_line_index() const;
655 
659  unsigned int
660  get_first_quad_index(const unsigned int quad_no = 0) const;
661 
665  unsigned int
666  get_first_hex_index() const;
667 
671  unsigned int
672  get_first_face_line_index(const unsigned int face_no = 0) const;
673 
677  unsigned int
678  get_first_face_quad_index(const unsigned int face_no = 0) const;
679 };
680 
681 namespace internal
682 {
688  expand(const unsigned int dim,
689  const std::vector<unsigned int> &dofs_per_object,
690  const ::ReferenceCell reference_cell);
691 } // namespace internal
692 
693 
694 
695 // --------- inline and template functions ---------------
696 
697 
698 #ifndef DOXYGEN
699 
700 namespace FiniteElementDomination
701 {
702  inline Domination
703  operator&(const Domination d1, const Domination d2)
704  {
705  // go through the entire list of possibilities. note that if we were into
706  // speed, obfuscation and cared enough, we could implement this operator
707  // by doing a bitwise & (and) if we gave these values to the enum values:
708  // neither_element_dominates=0, this_element_dominates=1,
709  // other_element_dominates=2, either_element_can_dominate=3
710  // =this_element_dominates|other_element_dominates
711  switch (d1)
712  {
714  if ((d2 == this_element_dominates) ||
715  (d2 == either_element_can_dominate) || (d2 == no_requirements))
716  return this_element_dominates;
717  else
719 
721  if ((d2 == other_element_dominates) ||
722  (d2 == either_element_can_dominate) || (d2 == no_requirements))
724  else
726 
729 
731  if (d2 == no_requirements)
733  else
734  return d2;
735 
736  case no_requirements:
737  return d2;
738 
739  default:
740  // shouldn't get here
741  Assert(false, ExcInternalError());
742  }
743 
745  }
746 } // namespace FiniteElementDomination
747 
748 
749 template <int dim>
750 inline ReferenceCell
752 {
753  return reference_cell_kind;
754 }
755 
756 
757 
758 template <int dim>
759 inline unsigned int
761 {
762  return number_unique_quads;
763 }
764 
765 
766 
767 template <int dim>
768 inline unsigned int
770 {
771  return number_unique_faces;
772 }
773 
774 
775 
776 template <int dim>
777 inline unsigned int
779 {
780  return dofs_per_vertex;
781 }
782 
783 
784 
785 template <int dim>
786 inline unsigned int
788 {
789  return dofs_per_line;
790 }
791 
792 
793 
794 template <int dim>
795 inline unsigned int
796 FiniteElementData<dim>::n_dofs_per_quad(unsigned int face_no) const
797 {
798  return n_dofs_on_quad[n_dofs_on_quad.size() == 1 ? 0 : face_no];
799 }
800 
801 
802 
803 template <int dim>
804 inline unsigned int
806 {
807  return dofs_per_quad_max;
808 }
809 
810 
811 
812 template <int dim>
813 inline unsigned int
815 {
816  return dofs_per_hex;
817 }
818 
819 
820 
821 template <int dim>
822 inline unsigned int
823 FiniteElementData<dim>::n_dofs_per_face(unsigned int face_no,
824  unsigned int child_no) const
825 {
826  (void)child_no;
827 
828  return n_dofs_on_face[n_dofs_on_face.size() == 1 ? 0 : face_no];
829 }
830 
831 
832 
833 template <int dim>
834 inline unsigned int
836 {
837  return dofs_per_face_max;
838 }
839 
840 
841 
842 template <int dim>
843 inline unsigned int
845 {
846  return dofs_per_cell;
847 }
848 
849 
850 
851 template <int dim>
852 template <int structdim>
853 inline unsigned int
854 FiniteElementData<dim>::n_dofs_per_object(const unsigned int i) const
855 {
856  switch (structdim)
857  {
858  case 0:
859  return n_dofs_per_vertex();
860  case 1:
861  return n_dofs_per_line();
862  case 2:
863  return n_dofs_per_quad((structdim == 2 && dim == 3) ? i : 0);
864  case 3:
865  return n_dofs_per_hex();
866  default:
867  Assert(false, ExcInternalError());
868  }
870 }
871 
872 
873 
874 template <int dim>
875 inline unsigned int
877 {
878  return components;
879 }
880 
881 
882 
883 template <int dim>
884 inline const BlockIndices &
886 {
887  return block_indices_data;
888 }
889 
890 
891 
892 template <int dim>
893 inline unsigned int
895 {
896  return block_indices_data.size();
897 }
898 
899 
900 
901 template <int dim>
902 inline unsigned int
904 {
905  return degree;
906 }
907 
908 
909 template <int dim>
910 inline bool
912 {
913  return ((space & conforming_space) == space);
914 }
915 
916 
917 
918 template <int dim>
919 unsigned int
921 {
922  return first_line_index;
923 }
924 
925 template <int dim>
926 unsigned int
927 FiniteElementData<dim>::get_first_quad_index(const unsigned int quad_no) const
928 {
929  if (first_index_of_quads.size() == 1)
930  return first_index_of_quads[0] + quad_no * n_dofs_per_quad(0);
931  else
932  return first_index_of_quads[quad_no];
933 }
934 
935 template <int dim>
936 unsigned int
938 {
939  return first_hex_index;
940 }
941 
942 template <int dim>
943 unsigned int
945  const unsigned int face_no) const
946 {
947  return first_line_index_of_faces[first_line_index_of_faces.size() == 1 ?
948  0 :
949  face_no];
950 }
951 
952 template <int dim>
953 unsigned int
955  const unsigned int face_no) const
956 {
957  return first_quad_index_of_faces[first_quad_index_of_faces.size() == 1 ?
958  0 :
959  face_no];
960 }
961 
962 
963 #endif // DOXYGEN
964 
965 
967 
968 #endif
const unsigned int first_hex_index
Definition: fe_base.h:367
std::vector< std::vector< unsigned int > > dofs_per_object_inclusive
Definition: fe_base.h:184
static const unsigned int invalid_unsigned_int
Definition: types.h:196
std::vector< std::vector< unsigned int > > object_index
Definition: fe_base.h:189
const unsigned int components
Definition: fe_base.h:430
const unsigned int dofs_per_quad_max
Definition: fe_base.h:336
unsigned int get_first_line_index() const
const std::vector< unsigned int > n_dofs_on_face
Definition: fe_base.h:398
internal::GenericDoFsPerObject expand(const unsigned int dim, const std::vector< unsigned int > &dofs_per_object, const ::ReferenceCell reference_cell)
Definition: fe_data.cc:25
std::vector< std::vector< unsigned int > > dofs_per_object_exclusive
Definition: fe_base.h:179
const unsigned int dofs_per_quad
Definition: fe_base.h:330
const unsigned int degree
Definition: fe_base.h:436
unsigned int n_dofs_per_vertex() const
unsigned int max_dofs_per_face() const
unsigned int n_unique_quads() const
unsigned int get_first_face_line_index(const unsigned int face_no=0) const
unsigned int n_dofs_per_quad(unsigned int face_no=0) const
const std::vector< unsigned int > n_dofs_on_quad
Definition: fe_base.h:323
unsigned int n_dofs_per_hex() const
const unsigned int dofs_per_line
Definition: fe_base.h:316
unsigned int get_first_face_quad_index(const unsigned int face_no=0) const
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
const ReferenceCell reference_cell_kind
Definition: fe_base.h:292
unsigned int max_dofs_per_quad() const
unsigned int n_blocks() const
const unsigned int first_face_line_index
Definition: fe_base.h:379
const unsigned int first_quad_index
Definition: fe_base.h:362
unsigned int n_dofs_per_line() const
const unsigned int number_unique_faces
Definition: fe_base.h:304
unsigned int get_first_hex_index() const
const unsigned int dofs_per_hex
Definition: fe_base.h:343
#define Assert(cond, exc)
Definition: exceptions.h:1461
const unsigned int number_unique_quads
Definition: fe_base.h:298
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
std::enable_if< IsBlockVector< VectorType >::value, unsigned int >::type n_blocks(const VectorType &vector)
Definition: operators.h:50
const unsigned int dofs_per_face_max
Definition: fe_base.h:412
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_unique_faces() const
const unsigned int dofs_per_cell
Definition: fe_base.h:420
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
const BlockIndices block_indices_data
Definition: fe_base.h:448
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
const unsigned int first_face_quad_index
Definition: fe_base.h:391
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
unsigned int n_dofs_per_object(const unsigned int i=0) const
bool conforms(const Conformity) const
const Conformity conforming_space
Definition: fe_base.h:441
const unsigned int dofs_per_face
Definition: fe_base.h:406
const unsigned int first_line_index
Definition: fe_base.h:348
std::vector< std::vector< unsigned int > > first_object_index_on_face
Definition: fe_base.h:194
Domination operator &(const Domination d1, const Domination d2)
const unsigned int dofs_per_vertex
Definition: fe_base.h:310
unsigned int size() const
const std::vector< unsigned int > first_quad_index_of_faces
Definition: fe_base.h:385
const BlockIndices & block_indices() const
void L2(Vector< number > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const double factor=1.)
Definition: l2.h:160
ReferenceCell reference_cell() const
const std::vector< unsigned int > first_index_of_quads
Definition: fe_base.h:356
unsigned int tensor_degree() const
const std::vector< unsigned int > first_line_index_of_faces
Definition: fe_base.h:373
static ::ExceptionBase & ExcInternalError()