Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_collection.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2019 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_collection_h
17 #define dealii_fe_collection_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <deal.II/fe/component_mask.h>
22 #include <deal.II/fe/fe.h>
23 #include <deal.II/fe/fe_values_extractors.h>
24 
25 #include <memory>
26 
27 DEAL_II_NAMESPACE_OPEN
28 
29 namespace hp
30 {
53  template <int dim, int spacedim = dim>
54  class FECollection : public Subscriptor
55  {
56  public:
67  {
74  static unsigned int
75  next_index(const typename hp::FECollection<dim, spacedim> &fe_collection,
76  const unsigned int fe_index)
77  {
78  return ((fe_index + 1) < fe_collection.size()) ? fe_index + 1 :
79  fe_index;
80  }
81 
88  static unsigned int
90  const typename hp::FECollection<dim, spacedim> &fe_collection,
91  const unsigned int fe_index)
92  {
93  (void)fe_collection;
94  return (fe_index > 0) ? fe_index - 1 : fe_index;
95  }
96  };
97 
103  FECollection();
104 
111  explicit FECollection(const FiniteElement<dim, spacedim> &fe);
112 
119  template <class... FETypes>
120  explicit FECollection(const FETypes &... fes);
121 
129  FECollection(const std::vector<const FiniteElement<dim, spacedim> *> &fes);
130 
134  FECollection(const FECollection<dim, spacedim> &) = default;
135 
146  std::is_nothrow_move_constructible<
147  std::vector<std::shared_ptr<const FiniteElement<dim, spacedim>>>>::value
148  &&std::is_nothrow_move_constructible<std::function<
149  unsigned int(const typename hp::FECollection<dim, spacedim> &,
150  const unsigned int)>>::value) = default;
151 
155  FECollection<dim, spacedim> &
156  operator=(FECollection<dim, spacedim> &&) = default; // NOLINT
157 
162  bool
163  operator==(const FECollection<dim, spacedim> &fe_collection) const;
164 
169  bool
170  operator!=(const FECollection<dim, spacedim> &fe_collection) const;
171 
181  void
182  push_back(const FiniteElement<dim, spacedim> &new_fe);
183 
190  const FiniteElement<dim, spacedim> &
191  operator[](const unsigned int index) const;
192 
196  unsigned int
197  size() const;
198 
208  unsigned int
209  n_components() const;
210 
227  unsigned int
228  n_blocks() const;
229 
234  unsigned int
235  max_dofs_per_vertex() const;
236 
241  unsigned int
242  max_dofs_per_line() const;
243 
248  unsigned int
249  max_dofs_per_quad() const;
250 
255  unsigned int
256  max_dofs_per_hex() const;
257 
262  unsigned int
263  max_dofs_per_face() const;
264 
269  unsigned int
270  max_dofs_per_cell() const;
271 
275  std::size_t
276  memory_consumption() const;
277 
278 
296  bool
298 
346  DEAL_II_DEPRECATED unsigned int
347  find_least_face_dominating_fe(const std::set<unsigned int> &fes) const;
348 
370  std::set<unsigned int>
371  find_common_fes(const std::set<unsigned int> &fes,
372  const unsigned int codim = 0) const;
373 
395  std::set<unsigned int>
396  find_enclosing_fes(const std::set<unsigned int> &fes,
397  const unsigned int codim = 0) const;
398 
431  unsigned int
432  find_dominating_fe(const std::set<unsigned int> &fes,
433  const unsigned int codim = 0) const;
434 
467  unsigned int
468  find_dominated_fe(const std::set<unsigned int> &fes,
469  const unsigned int codim = 0) const;
470 
494  unsigned int
495  find_dominating_fe_extended(const std::set<unsigned int> &fes,
496  const unsigned int codim = 0) const;
497 
520  unsigned int
521  find_dominated_fe_extended(const std::set<unsigned int> &fes,
522  const unsigned int codim = 0) const;
523 
536  void
537  set_hierarchy(const std::function<unsigned int(
538  const typename hp::FECollection<dim, spacedim> &,
539  const unsigned int)> &next,
540  const std::function<unsigned int(
541  const typename hp::FECollection<dim, spacedim> &,
542  const unsigned int)> &prev);
543 
551  void
553 
563  unsigned int
564  next_in_hierarchy(const unsigned int fe_index) const;
565 
575  unsigned int
576  previous_in_hierarchy(const unsigned int fe_index) const;
577 
595  component_mask(const FEValuesExtractors::Scalar &scalar) const;
596 
614  component_mask(const FEValuesExtractors::Vector &vector) const;
615 
635  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
636 
658  component_mask(const BlockMask &block_mask) const;
659 
686  BlockMask
687  block_mask(const FEValuesExtractors::Scalar &scalar) const;
688 
711  BlockMask
712  block_mask(const FEValuesExtractors::Vector &vector) const;
713 
737  BlockMask
738  block_mask(const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
739 
768  BlockMask
770 
771 
776 
777  private:
781  std::vector<std::shared_ptr<const FiniteElement<dim, spacedim>>>
783 
788  std::function<unsigned int(const typename hp::FECollection<dim, spacedim> &,
789  const unsigned int)>
791 
796  std::function<unsigned int(const typename hp::FECollection<dim, spacedim> &,
797  const unsigned int)>
799  };
800 
801 
802 
803  /* --------------- inline functions ------------------- */
804 
805  template <int dim, int spacedim>
806  template <class... FETypes>
807  FECollection<dim, spacedim>::FECollection(const FETypes &... fes)
808  {
809  static_assert(
810  is_base_of_all<FiniteElement<dim, spacedim>, FETypes...>::value,
811  "Not all of the input arguments of this function "
812  "are derived from FiniteElement<dim,spacedim>!");
813 
814  // loop over all of the given arguments and add the finite elements to
815  // this collection. Inlining the definition of fe_pointers causes internal
816  // compiler errors on GCC 7.1.1 so we define it separately:
817  const auto fe_pointers = {&fes...};
818  for (const auto p : fe_pointers)
819  push_back(*p);
820  }
821 
822 
823  template <int dim, int spacedim>
824  inline unsigned int
826  {
827  return finite_elements.size();
828  }
829 
830 
831  template <int dim, int spacedim>
832  inline unsigned int
834  {
836 
837  // note that there is no need
838  // here to enforce that indeed
839  // all elements have the same
840  // number of components since we
841  // have already done this when
842  // adding a new element to the
843  // collection.
844 
845  return finite_elements[0]->n_components();
846  }
847 
848 
849 
850  template <int dim, int spacedim>
851  inline bool
853  operator==(const FECollection<dim, spacedim> &fe_collection) const
854  {
855  const unsigned int n_elements = size();
856  if (n_elements != fe_collection.size())
857  return false;
858 
859  for (unsigned int i = 0; i < n_elements; ++i)
860  if (!(*finite_elements[i] == fe_collection[i]))
861  return false;
862 
863  return true;
864  }
865 
866 
867 
868  template <int dim, int spacedim>
869  inline bool
871  operator!=(const FECollection<dim, spacedim> &fe_collection) const
872  {
873  return !(*this == fe_collection);
874  }
875 
876 
877 
878  template <int dim, int spacedim>
880  operator[](const unsigned int index) const
881  {
882  Assert(index < finite_elements.size(),
883  ExcIndexRange(index, 0, finite_elements.size()));
884  return *finite_elements[index];
885  }
886 
887 
888 
889  template <int dim, int spacedim>
890  unsigned int
892  {
894 
895  unsigned int max = 0;
896  for (unsigned int i = 0; i < finite_elements.size(); ++i)
897  if (finite_elements[i]->dofs_per_vertex > max)
898  max = finite_elements[i]->dofs_per_vertex;
899 
900  return max;
901  }
902 
903 
904 
905  template <int dim, int spacedim>
906  unsigned int
908  {
910 
911  unsigned int max = 0;
912  for (unsigned int i = 0; i < finite_elements.size(); ++i)
913  if (finite_elements[i]->dofs_per_line > max)
914  max = finite_elements[i]->dofs_per_line;
915 
916  return max;
917  }
918 
919 
920 
921  template <int dim, int spacedim>
922  unsigned int
924  {
926 
927  unsigned int max = 0;
928  for (unsigned int i = 0; i < finite_elements.size(); ++i)
929  if (finite_elements[i]->dofs_per_quad > max)
930  max = finite_elements[i]->dofs_per_quad;
931 
932  return max;
933  }
934 
935 
936 
937  template <int dim, int spacedim>
938  unsigned int
940  {
942 
943  unsigned int max = 0;
944  for (unsigned int i = 0; i < finite_elements.size(); ++i)
945  if (finite_elements[i]->dofs_per_hex > max)
946  max = finite_elements[i]->dofs_per_hex;
947 
948  return max;
949  }
950 
951 
952 
953  template <int dim, int spacedim>
954  unsigned int
956  {
958 
959  unsigned int max = 0;
960  for (unsigned int i = 0; i < finite_elements.size(); ++i)
961  if (finite_elements[i]->dofs_per_face > max)
962  max = finite_elements[i]->dofs_per_face;
963 
964  return max;
965  }
966 
967 
968 
969  template <int dim, int spacedim>
970  unsigned int
972  {
974 
975  unsigned int max = 0;
976  for (unsigned int i = 0; i < finite_elements.size(); ++i)
977  if (finite_elements[i]->dofs_per_cell > max)
978  max = finite_elements[i]->dofs_per_cell;
979 
980  return max;
981  }
982 
983 
984  template <int dim, int spacedim>
985  bool
987  {
989 
990  bool hp_constraints = true;
991  for (unsigned int i = 0; i < finite_elements.size(); ++i)
992  hp_constraints =
993  hp_constraints && finite_elements[i]->hp_constraints_are_implemented();
994 
995  return hp_constraints;
996  }
997 
998 
999 } // namespace hp
1000 
1001 DEAL_II_NAMESPACE_CLOSE
1002 
1003 #endif
std::function< unsigned int(const typename hp::FECollection< dim, spacedim > &, const unsigned int)> hierarchy_prev
unsigned int find_dominated_fe_extended(const std::set< unsigned int > &fes, const unsigned int codim=0) const
static unsigned int next_index(const typename hp::FECollection< dim, spacedim > &fe_collection, const unsigned int fe_index)
Definition: fe_collection.h:75
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
std::size_t memory_consumption() const
unsigned int find_least_face_dominating_fe(const std::set< unsigned int > &fes) const
unsigned int find_dominating_fe(const std::set< unsigned int > &fes, const unsigned int codim=0) const
bool operator!=(const FECollection< dim, spacedim > &fe_collection) const
bool operator==(const FECollection< dim, spacedim > &fe_collection) const
STL namespace.
unsigned int size() const
unsigned int next_in_hierarchy(const unsigned int fe_index) const
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
unsigned int find_dominating_fe_extended(const std::set< unsigned int > &fes, const unsigned int codim=0) const
LinearAlgebra::distributed::Vector< Number > Vector
std::function< unsigned int(const typename hp::FECollection< dim, spacedim > &, const unsigned int)> hierarchy_next
unsigned int n_blocks() const
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
unsigned int max_dofs_per_cell() const
void push_back(const FiniteElement< dim, spacedim > &new_fe)
bool hp_constraints_are_implemented() const
static ::ExceptionBase & ExcNoFiniteElements()
unsigned int max_dofs_per_hex() const
#define Assert(cond, exc)
Definition: exceptions.h:1407
#define DeclException0(Exception0)
Definition: exceptions.h:473
unsigned int max_dofs_per_line() const
void set_hierarchy(const std::function< unsigned int(const typename hp::FECollection< dim, spacedim > &, const unsigned int)> &next, const std::function< unsigned int(const typename hp::FECollection< dim, spacedim > &, const unsigned int)> &prev)
Definition: hp.h:117
static unsigned int previous_index(const typename hp::FECollection< dim, spacedim > &fe_collection, const unsigned int fe_index)
Definition: fe_collection.h:89
unsigned int max_dofs_per_face() const
__global__ void set(Number *val, const Number s, const size_type N)
unsigned int max_dofs_per_vertex() const
std::set< unsigned int > find_common_fes(const std::set< unsigned int > &fes, const unsigned int codim=0) const
void set_default_hierarchy()
std::vector< std::shared_ptr< const FiniteElement< dim, spacedim > > > finite_elements
std::set< unsigned int > find_enclosing_fes(const std::set< unsigned int > &fes, const unsigned int codim=0) const
unsigned int find_dominated_fe(const std::set< unsigned int > &fes, const unsigned int codim=0) const
unsigned int n_components() const
const FiniteElement< dim, spacedim > & operator[](const unsigned int index) const
unsigned int max_dofs_per_quad() const
unsigned int previous_in_hierarchy(const unsigned int fe_index) const