Reference documentation for deal.II version GIT 0b65fff18a 2023-09-27 19:30:02+00:00
\(\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_data.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2000 - 2023 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_data_h
17 #define dealii_fe_data_h
18 
19 #include <deal.II/base/config.h>
20 
22 
24 
26 
27 #include <vector>
28 
30 
31 // Forward declarations:
32 #ifndef DOXYGEN
33 template <int dim>
34 class FiniteElementData;
35 #endif
36 
42 {
92  {
113  };
114 
115 
134  inline Domination
135  operator&(const Domination d1, const Domination d2);
136 } // namespace FiniteElementDomination
137 
138 namespace internal
139 {
184  {
188  std::vector<std::vector<unsigned int>> dofs_per_object_exclusive;
189 
193  std::vector<std::vector<unsigned int>> dofs_per_object_inclusive;
194 
198  std::vector<std::vector<unsigned int>> object_index;
199 
203  std::vector<std::vector<unsigned int>> first_object_index_on_face;
204 
208  template <int dim>
209  static GenericDoFsPerObject
211  };
212 } // namespace internal
213 
228 template <int dim>
230 {
231 public:
263  {
267  unknown = 0x00,
268 
272  L2 = 0x01,
273 
278  Hcurl = 0x02,
279 
284  Hdiv = 0x04,
285 
289  H1 = Hcurl | Hdiv,
290 
295  H2 = 0x0e
296  };
297 
302  static constexpr unsigned int dimension = dim;
303 
304 private:
309 
315  const unsigned int number_of_unique_2d_subobjects;
316 
321  const unsigned int number_unique_faces;
322 
323 public:
327  const unsigned int dofs_per_vertex;
328 
333  const unsigned int dofs_per_line;
334 
335 private:
340  const std::vector<unsigned int> n_dofs_on_quad;
341 
342 public:
347  const unsigned int dofs_per_quad;
348 
349 private:
353  const unsigned int dofs_per_quad_max;
354 
355 public:
360  const unsigned int dofs_per_hex;
361 
365  const unsigned int first_line_index;
366 
367 private:
373  const std::vector<unsigned int> first_index_of_quads;
374 
375 public:
379  const unsigned int first_quad_index;
380 
384  const unsigned int first_hex_index;
385 
386 private:
390  const std::vector<unsigned int> first_line_index_of_faces;
391 
392 public:
396  const unsigned int first_face_line_index;
397 
398 private:
402  const std::vector<unsigned int> first_quad_index_of_faces;
403 
404 public:
408  const unsigned int first_face_quad_index;
409 
410 private:
415  const std::vector<unsigned int> n_dofs_on_face;
416 
417 public:
423  const unsigned int dofs_per_face;
424 
425 private:
429  const unsigned int dofs_per_face_max;
430 
431 public:
437  const unsigned int dofs_per_cell;
438 
447  const unsigned int components;
448 
453  const unsigned int degree;
454 
459 
466 
505  FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
506  const unsigned int n_components,
507  const unsigned int degree,
508  const Conformity conformity = unknown,
510 
515  FiniteElementData(const std::vector<unsigned int> &dofs_per_object,
517  const unsigned int n_components,
518  const unsigned int degree,
519  const Conformity conformity = unknown,
521 
530  const unsigned int n_components,
531  const unsigned int degree,
532  const Conformity conformity = unknown,
534 
541  reference_cell() const;
542 
555  unsigned int
557 
564  unsigned int
565  n_unique_faces() const;
566 
570  unsigned int
572 
576  unsigned int
578 
582  unsigned int
583  n_dofs_per_quad(unsigned int face_no = 0) const;
584 
589  unsigned int
591 
595  unsigned int
596  n_dofs_per_hex() const;
597 
602  unsigned int
603  n_dofs_per_face(unsigned int face_no = 0, unsigned int child = 0) const;
604 
609  unsigned int
611 
616  unsigned int
618 
627  template <int structdim>
628  unsigned int
629  n_dofs_per_object(const unsigned int i = 0) const;
630 
636  unsigned int
637  n_components() const;
638 
644  unsigned int
645  n_blocks() const;
646 
650  const BlockIndices &
651  block_indices() const;
652 
659  unsigned int
660  tensor_degree() const;
661 
668  bool
669  conforms(const Conformity) const;
670 
674  bool
675  operator==(const FiniteElementData &) const;
676 
680  unsigned int
682 
686  unsigned int
687  get_first_quad_index(const unsigned int quad_no = 0) const;
688 
692  unsigned int
694 
698  unsigned int
699  get_first_face_line_index(const unsigned int face_no = 0) const;
700 
704  unsigned int
705  get_first_face_quad_index(const unsigned int face_no = 0) const;
706 };
707 
708 namespace internal
709 {
715  expand(const unsigned int dim,
716  const std::vector<unsigned int> &dofs_per_object,
718 } // namespace internal
719 
720 
721 
722 // --------- inline and template functions ---------------
723 
724 
725 #ifndef DOXYGEN
726 
727 namespace FiniteElementDomination
728 {
729  inline Domination
730  operator&(const Domination d1, const Domination d2)
731  {
732  // go through the entire list of possibilities. note that if we were into
733  // speed, obfuscation and cared enough, we could implement this operator
734  // by doing a bitwise & (and) if we gave these values to the enum values:
735  // neither_element_dominates=0, this_element_dominates=1,
736  // other_element_dominates=2, either_element_can_dominate=3
737  // =this_element_dominates|other_element_dominates
738  switch (d1)
739  {
741  if ((d2 == this_element_dominates) ||
742  (d2 == either_element_can_dominate) || (d2 == no_requirements))
743  return this_element_dominates;
744  else
746 
748  if ((d2 == other_element_dominates) ||
749  (d2 == either_element_can_dominate) || (d2 == no_requirements))
751  else
753 
756 
758  if (d2 == no_requirements)
760  else
761  return d2;
762 
763  case no_requirements:
764  return d2;
765 
766  default:
767  // shouldn't get here
768  Assert(false, ExcInternalError());
769  }
770 
772  }
773 } // namespace FiniteElementDomination
774 
775 
776 template <int dim>
777 inline ReferenceCell
779 {
780  return reference_cell_kind;
781 }
782 
783 
784 
785 template <int dim>
786 inline unsigned int
788 {
789  return number_of_unique_2d_subobjects;
790 }
791 
792 
793 
794 template <int dim>
795 inline unsigned int
797 {
798  return number_unique_faces;
799 }
800 
801 
802 
803 template <int dim>
804 inline unsigned int
806 {
807  return dofs_per_vertex;
808 }
809 
810 
811 
812 template <int dim>
813 inline unsigned int
815 {
816  return dofs_per_line;
817 }
818 
819 
820 
821 template <int dim>
822 inline unsigned int
823 FiniteElementData<dim>::n_dofs_per_quad(unsigned int face_no) const
824 {
825  return n_dofs_on_quad[n_dofs_on_quad.size() == 1 ? 0 : face_no];
826 }
827 
828 
829 
830 template <int dim>
831 inline unsigned int
833 {
834  return dofs_per_quad_max;
835 }
836 
837 
838 
839 template <int dim>
840 inline unsigned int
842 {
843  return dofs_per_hex;
844 }
845 
846 
847 
848 template <int dim>
849 inline unsigned int
850 FiniteElementData<dim>::n_dofs_per_face(unsigned int face_no,
851  unsigned int child_no) const
852 {
853  (void)child_no;
854 
855  return n_dofs_on_face[n_dofs_on_face.size() == 1 ? 0 : face_no];
856 }
857 
858 
859 
860 template <int dim>
861 inline unsigned int
863 {
864  return dofs_per_face_max;
865 }
866 
867 
868 
869 template <int dim>
870 inline unsigned int
872 {
873  return dofs_per_cell;
874 }
875 
876 
877 
878 template <int dim>
879 template <int structdim>
880 inline unsigned int
881 FiniteElementData<dim>::n_dofs_per_object(const unsigned int i) const
882 {
883  switch (structdim)
884  {
885  case 0:
886  return n_dofs_per_vertex();
887  case 1:
888  return n_dofs_per_line();
889  case 2:
890  return n_dofs_per_quad((structdim == 2 && dim == 3) ? i : 0);
891  case 3:
892  return n_dofs_per_hex();
893  default:
894  Assert(false, ExcInternalError());
895  }
897 }
898 
899 
900 
901 template <int dim>
902 inline unsigned int
904 {
905  return components;
906 }
907 
908 
909 
910 template <int dim>
911 inline const BlockIndices &
913 {
914  return block_indices_data;
915 }
916 
917 
918 
919 template <int dim>
920 inline unsigned int
922 {
923  return block_indices_data.size();
924 }
925 
926 
927 
928 template <int dim>
929 inline unsigned int
931 {
932  return degree;
933 }
934 
935 
936 template <int dim>
937 inline bool
938 FiniteElementData<dim>::conforms(const Conformity space) const
939 {
940  return ((space & conforming_space) == space);
941 }
942 
943 
944 
945 template <int dim>
946 unsigned int
948 {
949  return first_line_index;
950 }
951 
952 template <int dim>
953 unsigned int
954 FiniteElementData<dim>::get_first_quad_index(const unsigned int quad_no) const
955 {
956  if (first_index_of_quads.size() == 1)
957  return first_index_of_quads[0] + quad_no * n_dofs_per_quad(0);
958  else
959  return first_index_of_quads[quad_no];
960 }
961 
962 template <int dim>
963 unsigned int
965 {
966  return first_hex_index;
967 }
968 
969 template <int dim>
970 unsigned int
972  const unsigned int face_no) const
973 {
974  return first_line_index_of_faces[first_line_index_of_faces.size() == 1 ?
975  0 :
976  face_no];
977 }
978 
979 template <int dim>
980 unsigned int
982  const unsigned int face_no) const
983 {
984  return first_quad_index_of_faces[first_quad_index_of_faces.size() == 1 ?
985  0 :
986  face_no];
987 }
988 
989 template <int dim>
992 {
993  const auto reference_cell = fe.reference_cell();
994 
996 
997  result.dofs_per_object_exclusive.resize(4);
998  result.dofs_per_object_inclusive.resize(4);
999  result.object_index.resize(4);
1000 
1001  unsigned int counter = 0;
1002 
1003  for (const unsigned int v : reference_cell.vertex_indices())
1004  {
1005  const auto c = fe.template n_dofs_per_object<0>(v);
1006 
1007  result.dofs_per_object_exclusive[0].emplace_back(c);
1008  result.dofs_per_object_inclusive[0].emplace_back(c);
1009  result.object_index[0].emplace_back(counter);
1010 
1011  counter += c;
1012  }
1013 
1014  if (dim >= 2)
1015  for (const unsigned int l : reference_cell.line_indices())
1016  {
1017  const auto c = fe.template n_dofs_per_object<1>(l);
1018 
1019  result.dofs_per_object_exclusive[1].emplace_back(c);
1020  result.dofs_per_object_inclusive[1].emplace_back(
1021  c + 2 * fe.template n_dofs_per_object<0>());
1022  result.object_index[1].emplace_back(counter);
1023 
1024  counter += c;
1025  }
1026 
1027  if (dim == 3)
1028  for (const unsigned int f : reference_cell.face_indices())
1029  {
1030  const auto c = fe.template n_dofs_per_object<2>(f);
1031 
1032  result.dofs_per_object_exclusive[2].emplace_back(c);
1033  result.dofs_per_object_inclusive[2].emplace_back(fe.n_dofs_per_face(f));
1034  result.object_index[2].emplace_back(counter);
1035 
1036  counter += c;
1037  }
1038 
1039  {
1040  const auto c = fe.template n_dofs_per_object<dim>();
1041 
1042  result.dofs_per_object_exclusive[dim].emplace_back(c);
1043  result.dofs_per_object_inclusive[dim].emplace_back(fe.n_dofs_per_cell());
1044  result.object_index[dim].emplace_back(counter);
1045 
1046  counter += c;
1047  }
1048 
1049  for (unsigned int d = dim + 1; d <= 3; ++d)
1050  {
1051  result.dofs_per_object_exclusive[d].emplace_back(0);
1052  result.dofs_per_object_inclusive[d].emplace_back(0);
1053  result.object_index[d].emplace_back(counter);
1054  }
1055 
1056  result.first_object_index_on_face.resize(3);
1057  for (const unsigned int face_no : reference_cell.face_indices())
1058  {
1059  result.first_object_index_on_face[0].emplace_back(0);
1060 
1061  result.first_object_index_on_face[1].emplace_back(
1062  fe.get_first_face_line_index(face_no));
1063 
1064  result.first_object_index_on_face[2].emplace_back(
1065  fe.get_first_face_quad_index(face_no));
1066  }
1067 
1068  return result;
1069 }
1070 
1071 
1072 #endif // DOXYGEN
1073 
1074 
1076 
1077 #endif
unsigned int size() const
const unsigned int dofs_per_quad
Definition: fe_data.h:347
unsigned int get_first_line_index() const
unsigned int n_dofs_per_vertex() const
const unsigned int first_face_quad_index
Definition: fe_data.h:408
const std::vector< unsigned int > first_index_of_quads
Definition: fe_data.h:373
bool conforms(const Conformity) const
const std::vector< unsigned int > first_line_index_of_faces
Definition: fe_data.h:390
const unsigned int components
Definition: fe_data.h:447
const unsigned int degree
Definition: fe_data.h:453
const std::vector< unsigned int > n_dofs_on_face
Definition: fe_data.h:415
unsigned int n_dofs_per_cell() const
const unsigned int dofs_per_face
Definition: fe_data.h:423
unsigned int n_dofs_per_line() const
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
bool operator==(const FiniteElementData &) const
Definition: fe_data.cc:196
const unsigned int dofs_per_line
Definition: fe_data.h:333
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_blocks() const
unsigned int max_dofs_per_face() const
unsigned int tensor_degree() const
unsigned int n_components() const
ReferenceCell reference_cell() const
const unsigned int first_hex_index
Definition: fe_data.h:384
unsigned int n_dofs_per_object(const unsigned int i=0) const
const unsigned int first_face_line_index
Definition: fe_data.h:396
const unsigned int dofs_per_quad_max
Definition: fe_data.h:353
unsigned int n_unique_2d_subobjects() const
unsigned int get_first_face_line_index(const unsigned int face_no=0) const
const unsigned int first_line_index
Definition: fe_data.h:365
const unsigned int dofs_per_hex
Definition: fe_data.h:360
const BlockIndices block_indices_data
Definition: fe_data.h:465
unsigned int get_first_face_quad_index(const unsigned int face_no=0) const
const unsigned int number_unique_faces
Definition: fe_data.h:321
unsigned int max_dofs_per_quad() const
const std::vector< unsigned int > first_quad_index_of_faces
Definition: fe_data.h:402
const ReferenceCell reference_cell_kind
Definition: fe_data.h:308
const BlockIndices & block_indices() const
const unsigned int number_of_unique_2d_subobjects
Definition: fe_data.h:315
const std::vector< unsigned int > n_dofs_on_quad
Definition: fe_data.h:340
const Conformity conforming_space
Definition: fe_data.h:458
const unsigned int first_quad_index
Definition: fe_data.h:379
unsigned int n_unique_faces() const
unsigned int n_dofs_per_quad(unsigned int face_no=0) const
const unsigned int dofs_per_vertex
Definition: fe_data.h:327
const unsigned int dofs_per_cell
Definition: fe_data.h:437
static constexpr unsigned int dimension
Definition: fe_data.h:302
unsigned int n_dofs_per_hex() const
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:115
const unsigned int dofs_per_face_max
Definition: fe_data.h:429
unsigned int get_first_hex_index() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
static ::ExceptionBase & ExcInternalError()
#define Assert(cond, exc)
Definition: exceptions.h:1616
Domination operator&(const Domination d1, const Domination d2)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
internal::GenericDoFsPerObject expand(const unsigned int dim, const std::vector< unsigned int > &dofs_per_object, const ReferenceCell reference_cell)
Definition: fe_data.cc:25
static const unsigned int invalid_unsigned_int
Definition: types.h:213
std::vector< std::vector< unsigned int > > object_index
Definition: fe_data.h:198
std::vector< std::vector< unsigned int > > first_object_index_on_face
Definition: fe_data.h:203
std::vector< std::vector< unsigned int > > dofs_per_object_inclusive
Definition: fe_data.h:193
std::vector< std::vector< unsigned int > > dofs_per_object_exclusive
Definition: fe_data.h:188
static GenericDoFsPerObject generate(const FiniteElementData< dim > &fe)