Reference documentation for deal.II version Git e3a3ec7800 2020-08-07 14:08:19 +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 - 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 
22 
24 
26 
27 #include <vector>
28 
30 
36 {
86  {
107  };
108 
109 
128  inline Domination operator&(const Domination d1, const Domination d2);
129 } // namespace FiniteElementDomination
130 
131 namespace internal
132 {
170  {
174  std::vector<std::vector<unsigned int>> dofs_per_object_exclusive;
175 
179  std::vector<std::vector<unsigned int>> dofs_per_object_inclusive;
180 
184  std::vector<std::vector<unsigned int>> object_index;
185 
189  std::vector<std::vector<unsigned int>> first_object_index_on_face;
190  };
191 } // namespace internal
192 
207 template <int dim>
209 {
210 public:
242  {
246  unknown = 0x00,
247 
251  L2 = 0x01,
252 
257  Hcurl = 0x02,
258 
263  Hdiv = 0x04,
264 
268  H1 = Hcurl | Hdiv,
269 
274  H2 = 0x0e
275  };
276 
281  static const unsigned int dimension = dim;
282 
283 private:
288 
293  const unsigned int number_unique_quads;
294 
299  const unsigned int number_unique_faces;
300 
301 public:
305  const unsigned int dofs_per_vertex;
306 
311  const unsigned int dofs_per_line;
312 
313 private:
318  const std::vector<unsigned int> n_dofs_on_quad;
319 
320 public:
325  const unsigned int dofs_per_quad;
326 
327 private:
331  const unsigned int dofs_per_quad_max;
332 
333 public:
338  const unsigned int dofs_per_hex;
339 
343  const unsigned int first_line_index;
344 
345 private:
351  const std::vector<unsigned int> first_index_of_quads;
352 
353 public:
357  const unsigned int first_quad_index;
358 
362  const unsigned int first_hex_index;
363 
364 private:
368  const std::vector<unsigned int> first_line_index_of_faces;
369 
370 public:
374  const unsigned int first_face_line_index;
375 
376 private:
380  const std::vector<unsigned int> first_quad_index_of_faces;
381 
382 public:
386  const unsigned int first_face_quad_index;
387 
388 private:
393  const std::vector<unsigned int> n_dofs_on_face;
394 
395 public:
401  const unsigned int dofs_per_face;
402 
403 private:
407  const unsigned int dofs_per_face_max;
408 
409 public:
415  const unsigned int dofs_per_cell;
416 
425  const unsigned int components;
426 
431  const unsigned int degree;
432 
437 
444 
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 
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 
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 
517  reference_cell_type() const;
518 
523  unsigned int
524  n_unique_quads() const;
525 
530  unsigned int
531  n_unique_faces() const;
532 
536  unsigned int
537  n_dofs_per_vertex() const;
538 
542  unsigned int
543  n_dofs_per_line() const;
544 
548  unsigned int
549  n_dofs_per_quad(unsigned int face_no = 0) const;
550 
555  unsigned int
556  max_dofs_per_quad() const;
557 
561  unsigned int
562  n_dofs_per_hex() const;
563 
568  unsigned int
569  n_dofs_per_face(unsigned int face_no = 0, unsigned int child = 0) const;
570 
575  unsigned int
576  max_dofs_per_face() const;
577 
582  unsigned int
583  n_dofs_per_cell() const;
584 
593  template <int structdim>
594  unsigned int
595  n_dofs_per_object(const unsigned int i = 0) const;
596 
602  unsigned int
603  n_components() const;
604 
610  unsigned int
611  n_blocks() const;
612 
616  const BlockIndices &
617  block_indices() const;
618 
625  unsigned int
626  tensor_degree() const;
627 
634  bool
635  conforms(const Conformity) const;
636 
640  bool
641  operator==(const FiniteElementData &) const;
642 
646  unsigned int
647  get_first_line_index() const;
648 
652  unsigned int
653  get_first_quad_index(const unsigned int quad_no = 0) const;
654 
658  unsigned int
659  get_first_hex_index() const;
660 
664  unsigned int
665  get_first_face_line_index(const unsigned int face_no = 0) const;
666 
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  {
694  if ((d2 == this_element_dominates) ||
695  (d2 == either_element_can_dominate) || (d2 == no_requirements))
696  return this_element_dominates;
697  else
699 
701  if ((d2 == other_element_dominates) ||
702  (d2 == either_element_can_dominate) || (d2 == no_requirements))
704  else
706 
709 
711  if (d2 == no_requirements)
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 
725  }
726 } // namespace FiniteElementDomination
727 
728 
729 template <int dim>
730 inline ReferenceCell::Type
732 {
733  return cell_type;
734 }
735 
736 
737 
738 template <int dim>
739 inline unsigned int
741 {
742  return number_unique_quads;
743 }
744 
745 
746 
747 template <int dim>
748 inline unsigned int
750 {
751  return number_unique_faces;
752 }
753 
754 
755 
756 template <int dim>
757 inline unsigned int
759 {
760  return dofs_per_vertex;
761 }
762 
763 
764 
765 template <int dim>
766 inline unsigned int
768 {
769  return dofs_per_line;
770 }
771 
772 
773 
774 template <int dim>
775 inline unsigned int
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
786 {
787  return dofs_per_quad_max;
788 }
789 
790 
791 
792 template <int dim>
793 inline unsigned int
795 {
796  return dofs_per_hex;
797 }
798 
799 
800 
801 template <int dim>
802 inline unsigned int
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
816 {
817  return dofs_per_face_max;
818 }
819 
820 
821 
822 template <int dim>
823 inline unsigned int
825 {
826  return dofs_per_cell;
827 }
828 
829 
830 
831 template <int dim>
832 template <int structdim>
833 inline unsigned int
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  }
850 }
851 
852 
853 
854 template <int dim>
855 inline unsigned int
857 {
858  return components;
859 }
860 
861 
862 
863 template <int dim>
864 inline const BlockIndices &
866 {
867  return block_indices_data;
868 }
869 
870 
871 
872 template <int dim>
873 inline unsigned int
875 {
876  return block_indices_data.size();
877 }
878 
879 
880 
881 template <int dim>
882 inline unsigned int
884 {
885  return degree;
886 }
887 
888 
889 template <int dim>
890 inline bool
892 {
893  return ((space & conforming_space) == space);
894 }
895 
896 
897 
898 template <int dim>
899 unsigned int
901 {
902  return first_line_index;
903 }
904 
905 template <int dim>
906 unsigned int
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
918 {
919  return first_hex_index;
920 }
921 
922 template <int dim>
923 unsigned int
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
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 
947 
948 #endif
const unsigned int first_hex_index
Definition: fe_base.h:362
std::vector< std::vector< unsigned int > > dofs_per_object_inclusive
Definition: fe_base.h:179
static const unsigned int invalid_unsigned_int
Definition: types.h:196
std::vector< std::vector< unsigned int > > object_index
Definition: fe_base.h:184
const unsigned int components
Definition: fe_base.h:425
const unsigned int dofs_per_quad_max
Definition: fe_base.h:331
unsigned int get_first_line_index() const
const std::vector< unsigned int > n_dofs_on_face
Definition: fe_base.h:393
std::vector< std::vector< unsigned int > > dofs_per_object_exclusive
Definition: fe_base.h:174
const unsigned int dofs_per_quad
Definition: fe_base.h:325
const unsigned int degree
Definition: fe_base.h:431
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:318
unsigned int n_dofs_per_hex() const
const unsigned int dofs_per_line
Definition: fe_base.h:311
unsigned int get_first_face_quad_index(const unsigned int face_no=0) const
bool operator==(const AlignedVector< T > &lhs, const AlignedVector< T > &rhs)
unsigned int max_dofs_per_quad() const
unsigned int n_blocks() const
ReferenceCell::Type reference_cell_type() const
const unsigned int first_face_line_index
Definition: fe_base.h:374
const unsigned int first_quad_index
Definition: fe_base.h:357
unsigned int n_dofs_per_line() const
const unsigned int number_unique_faces
Definition: fe_base.h:299
unsigned int get_first_hex_index() const
const unsigned int dofs_per_hex
Definition: fe_base.h:338
#define Assert(cond, exc)
Definition: exceptions.h:1411
const unsigned int number_unique_quads
Definition: fe_base.h:293
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
std::enable_if< IsBlockVector< VectorType >::value, unsigned int >::type n_blocks(const VectorType &vector)
Definition: operators.h:49
const unsigned int dofs_per_face_max
Definition: fe_base.h:407
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:415
const BlockIndices block_indices_data
Definition: fe_base.h:443
const ReferenceCell::Type cell_type
Definition: fe_base.h:287
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
const unsigned int first_face_quad_index
Definition: fe_base.h:386
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
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:436
const unsigned int dofs_per_face
Definition: fe_base.h:401
const unsigned int first_line_index
Definition: fe_base.h:343
std::vector< std::vector< unsigned int > > first_object_index_on_face
Definition: fe_base.h:189
Domination operator &(const Domination d1, const Domination d2)
const unsigned int dofs_per_vertex
Definition: fe_base.h:305
unsigned int size() const
const std::vector< unsigned int > first_quad_index_of_faces
Definition: fe_base.h:380
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
const std::vector< unsigned int > first_index_of_quads
Definition: fe_base.h:351
unsigned int tensor_degree() const
const std::vector< unsigned int > first_line_index_of_faces
Definition: fe_base.h:368
static ::ExceptionBase & ExcInternalError()