Reference documentation for deal.II version GIT c415589cf0 2022-08-14 18:50: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\}}\)
dof_handler.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 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_dof_handler_h
17 #define dealii_dof_handler_h
18 
19 
20 
21 #include <deal.II/base/config.h>
22 
24 #include <deal.II/base/function.h>
25 #include <deal.II/base/index_set.h>
28 
30 
33 #include <deal.II/dofs/dof_faces.h>
37 
39 
40 #include <boost/serialization/split_member.hpp>
41 
42 #include <map>
43 #include <memory>
44 #include <set>
45 #include <vector>
46 
48 
49 // Forward declarations
50 #ifndef DOXYGEN
51 template <int dim, int spacedim>
52 class FiniteElement;
53 template <int dim, int spacedim>
54 class Triangulation;
55 
56 namespace internal
57 {
58  namespace DoFHandlerImplementation
59  {
60  struct Implementation;
61 
62  namespace Policy
63  {
64  template <int dim, int spacedim>
65  class PolicyBase;
66  struct Implementation;
67  } // namespace Policy
68  } // namespace DoFHandlerImplementation
69 
70  namespace DoFAccessorImplementation
71  {
72  struct Implementation;
73  }
74 
75  namespace DoFCellAccessorImplementation
76  {
77  struct Implementation;
78  }
79 
80  namespace hp
81  {
82  namespace DoFHandlerImplementation
83  {
84  struct Implementation;
85  }
86  } // namespace hp
87 } // namespace internal
88 
89 namespace parallel
90 {
91  namespace distributed
92  {
93  template <int dim, int spacedim, typename VectorType>
94  class CellDataTransfer;
95  }
96 } // namespace parallel
97 #endif
98 
313 template <int dim, int spacedim = dim>
314 class DoFHandler : public Subscriptor
315 {
320 
321 public:
334  using cell_accessor = typename ActiveSelector::CellAccessor;
335 
348  using face_accessor = typename ActiveSelector::FaceAccessor;
349 
357  using line_iterator = typename ActiveSelector::line_iterator;
358 
372  using active_line_iterator = typename ActiveSelector::active_line_iterator;
373 
381  using quad_iterator = typename ActiveSelector::quad_iterator;
382 
396  using active_quad_iterator = typename ActiveSelector::active_quad_iterator;
397 
405  using hex_iterator = typename ActiveSelector::hex_iterator;
406 
416  using active_hex_iterator = typename ActiveSelector::active_hex_iterator;
417 
438  using active_cell_iterator = typename ActiveSelector::active_cell_iterator;
439 
466  using cell_iterator = typename ActiveSelector::cell_iterator;
467 
484  using face_iterator = typename ActiveSelector::face_iterator;
485 
497  using active_face_iterator = typename ActiveSelector::active_face_iterator;
498 
499  using level_cell_accessor = typename LevelSelector::CellAccessor;
500  using level_face_accessor = typename LevelSelector::FaceAccessor;
501 
502  using level_cell_iterator = typename LevelSelector::cell_iterator;
503  using level_face_iterator = typename LevelSelector::face_iterator;
504 
505 
509  static constexpr unsigned int dimension = dim;
510 
514  static constexpr unsigned int space_dimension = spacedim;
515 
519  static const unsigned int default_fe_index = 0;
520 
524  static const unsigned int invalid_fe_index = numbers::invalid_unsigned_int;
525 
529  using active_fe_index_type = unsigned short int;
530 
534  using offset_type = unsigned int;
535 
541  static_cast<active_fe_index_type>(-1);
542 
548 
553 
560  DoFHandler(const DoFHandler &) = delete;
561 
565  virtual ~DoFHandler() override;
566 
573  DoFHandler &
574  operator=(const DoFHandler &) = delete;
575 
580  void
581  set_active_fe_indices(const std::vector<unsigned int> &active_fe_indices);
582 
588  void
589  get_active_fe_indices(std::vector<unsigned int> &active_fe_indices) const;
590 
598  void
600 
630  void
632 
636  void
638 
644  void
646 
650  bool
652 
658  bool
659  has_level_dofs() const;
660 
675  bool
677 
684  void
686 
690  void
691  clear();
692 
745  void
746  renumber_dofs(const std::vector<types::global_dof_index> &new_numbers);
747 
752  void
753  renumber_dofs(const unsigned int level,
754  const std::vector<types::global_dof_index> &new_numbers);
755 
784  unsigned int
786 
799  unsigned int
801 
802  /*--------------------------------------*/
803 
808  /*
809  * @{
810  */
811 
816  begin(const unsigned int level = 0) const;
817 
835  begin_active(const unsigned int level = 0) const;
836 
842  end() const;
843 
849  end(const unsigned int level) const;
850 
857  end_active(const unsigned int level) const;
858 
864  begin_mg(const unsigned int level = 0) const;
865 
871  end_mg(const unsigned int level) const;
872 
878  end_mg() const;
879 
895  cell_iterators() const;
896 
939 
953 
970  cell_iterators_on_level(const unsigned int level) const;
971 
988  active_cell_iterators_on_level(const unsigned int level) const;
989 
1006  mg_cell_iterators_on_level(const unsigned int level) const;
1007 
1008  /*
1009  * @}
1010  */
1011 
1012 
1013  /*---------------------------------------*/
1014 
1015 
1040  n_dofs() const;
1041 
1050  n_dofs(const unsigned int level) const;
1051 
1058 
1070  template <typename number>
1073  const std::map<types::boundary_id, const Function<spacedim, number> *>
1074  &boundary_ids) const;
1075 
1081  n_boundary_dofs(const std::set<types::boundary_id> &boundary_ids) const;
1082 
1098  const BlockInfo &
1099  block_info() const;
1100 
1122 
1128  const IndexSet &
1130 
1135  const IndexSet &
1136  locally_owned_mg_dofs(const unsigned int level) const;
1137 
1143  get_fe(const unsigned int index = 0) const;
1144 
1151 
1157 
1161  MPI_Comm
1163 
1180  void
1182 
1198  void
1200 
1209  virtual std::size_t
1211 
1217  template <class Archive>
1218  void
1219  save(Archive &ar, const unsigned int version) const;
1220 
1226  template <class Archive>
1227  void
1228  load(Archive &ar, const unsigned int version);
1229 
1230 #ifdef DOXYGEN
1236  template <class Archive>
1237  void
1238  serialize(Archive &archive, const unsigned int version);
1239 #else
1240  // This macro defines the serialize() method that is compatible with
1241  // the templated save() and load() method that have been implemented.
1242  BOOST_SERIALIZATION_SPLIT_MEMBER()
1243 #endif
1244 
1259  int,
1260  << "The given level " << arg1
1261  << " is not in the valid range!");
1268  << "The given list of new dof indices is not consecutive: "
1269  << "the index " << arg1 << " does not exist.");
1274  int,
1275  int,
1276  << "The mesh contains a cell with an active FE index of "
1277  << arg1 << ", but the finite element collection only has "
1278  << arg2 << " elements");
1279 
1285  "The current function doesn't make sense when used with a "
1286  "DoFHandler without hp-capabilities.");
1287 
1293  "The current function has not yet been implemented for a "
1294  "DoFHandler with hp-capabilities.");
1295 
1296 private:
1304  {
1305  public:
1310 
1316  void
1317  init(const unsigned int coarsest_level,
1318  const unsigned int finest_level,
1319  const unsigned int dofs_per_vertex);
1320 
1324  unsigned int
1326 
1330  unsigned int
1332 
1338  access_index(const unsigned int level,
1339  const unsigned int dof_number,
1340  const unsigned int dofs_per_vertex);
1341 
1342  private:
1346  unsigned int coarsest_level;
1347 
1351  unsigned int finest_level;
1352 
1362  std::unique_ptr<types::global_dof_index[]> indices;
1363  };
1364 
1372  {
1377  std::map<const cell_iterator, const unsigned int> persisting_cells_fe_index;
1378 
1383  std::map<const cell_iterator, const unsigned int> refined_cells_fe_index;
1384 
1389  std::map<const cell_iterator, const unsigned int> coarsened_cells_fe_index;
1390 
1396  std::vector<unsigned int> active_fe_indices;
1397 
1403  std::unique_ptr<
1404  parallel::distributed::
1405  CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>
1407  };
1408 
1413 
1419 
1425 
1432 
1437  std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1438  PolicyBase<dim, spacedim>>
1440 
1449 
1453  std::vector<::internal::DoFHandlerImplementation::NumberCache>
1455 
1461  mutable std::vector<std::vector<types::global_dof_index>>
1463 
1468  mutable std::vector<std::vector<offset_type>> cell_dof_cache_ptr;
1469 
1475  mutable std::vector<std::array<std::vector<types::global_dof_index>, dim + 1>>
1477 
1486  mutable std::vector<std::array<std::vector<offset_type>, dim + 1>>
1488 
1494  mutable std::array<std::vector<active_fe_index_type>, dim + 1>
1496 
1500  mutable std::array<std::vector<offset_type>, dim + 1> hp_object_fe_ptr;
1501 
1506  mutable std::vector<std::vector<active_fe_index_type>>
1508 
1513  mutable std::vector<std::vector<active_fe_index_type>>
1515 
1520  std::vector<MGVertexDoFs> mg_vertex_dofs;
1521 
1525  std::vector<
1526  std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel<dim>>>
1528 
1532  std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces<dim>>
1534 
1539  std::unique_ptr<ActiveFEIndexTransfer> active_fe_index_transfer;
1540 
1545  std::vector<boost::signals2::connection> tria_listeners;
1546 
1552  std::vector<boost::signals2::connection> tria_listeners_for_transfer;
1553 
1557  void
1559 
1563  void
1565 
1569  void
1571 
1575  void
1577 
1589  void
1591 
1604  void
1606 
1616  void
1618 
1627  void
1629 
1638  void
1640 
1649  void
1651 
1652 
1653  // Make accessor objects friends.
1654  template <int, int, int, bool>
1655  friend class ::DoFAccessor;
1656  template <int, int, bool>
1657  friend class ::DoFCellAccessor;
1658  friend struct ::internal::DoFAccessorImplementation::Implementation;
1659  friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1660 
1661  // Likewise for DoFLevel objects since they need to access the vertex dofs
1662  // in the functions that set and retrieve vertex dof indices.
1663  friend struct ::internal::DoFHandlerImplementation::Implementation;
1664  friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1665  friend struct ::internal::DoFHandlerImplementation::Policy::
1666  Implementation;
1667 
1668  // explicitly check for sensible template arguments, but not on windows
1669  // because MSVC creates bogus warnings during normal compilation
1670 #ifndef DEAL_II_MSVC
1671  static_assert(dim <= spacedim,
1672  "The dimension <dim> of a DoFHandler must be less than or "
1673  "equal to the space dimension <spacedim> in which it lives.");
1674 #endif
1675 };
1676 
1677 namespace internal
1678 {
1679  namespace hp
1680  {
1681  namespace DoFHandlerImplementation
1682  {
1694  template <int dim, int spacedim>
1695  void
1697 
1722  template <int dim, int spacedim = dim>
1723  unsigned int
1725  const typename DoFHandler<dim, spacedim>::cell_iterator &parent);
1726 
1732  "No FiniteElement has been found in your FECollection that is "
1733  "dominated by all children of a cell you are trying to coarsen!");
1734  } // namespace DoFHandlerImplementation
1735  } // namespace hp
1736 } // namespace internal
1737 
1738 #ifndef DOXYGEN
1739 
1740 /* ----------------------- Inline functions ----------------------------------
1741  */
1742 
1743 
1744 template <int dim, int spacedim>
1745 inline bool
1747 {
1748  return hp_capability_enabled;
1749 }
1750 
1751 
1752 
1753 template <int dim, int spacedim>
1754 inline bool
1756 {
1757  return mg_number_cache.size() > 0;
1758 }
1759 
1760 
1761 
1762 template <int dim, int spacedim>
1763 inline bool
1765 {
1766  return number_cache.n_global_dofs > 0;
1767 }
1768 
1769 
1770 
1771 template <int dim, int spacedim>
1774 {
1775  return number_cache.n_global_dofs;
1776 }
1777 
1778 
1779 
1780 template <int dim, int spacedim>
1782 DoFHandler<dim, spacedim>::n_dofs(const unsigned int level) const
1783 {
1785  ExcMessage(
1786  "n_dofs(level) can only be called after distribute_mg_dofs()"));
1788  return mg_number_cache[level].n_global_dofs;
1789 }
1790 
1791 
1792 
1793 template <int dim, int spacedim>
1796 {
1798 }
1799 
1800 
1801 
1802 template <int dim, int spacedim>
1803 const IndexSet &
1805 {
1807 }
1808 
1809 
1810 
1811 template <int dim, int spacedim>
1812 const IndexSet &
1814 {
1815  Assert(level < this->get_triangulation().n_global_levels(),
1816  ExcMessage("The given level index exceeds the number of levels "
1817  "present in the triangulation"));
1818  Assert(
1819  mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1820  ExcMessage(
1821  "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1822  return mg_number_cache[level].locally_owned_dofs;
1823 }
1824 
1825 
1826 
1827 template <int dim, int spacedim>
1828 inline const FiniteElement<dim, spacedim> &
1829 DoFHandler<dim, spacedim>::get_fe(const unsigned int number) const
1830 {
1831  Assert(fe_collection.size() > 0,
1832  ExcMessage("No finite element collection is associated with "
1833  "this DoFHandler"));
1834  return fe_collection[number];
1835 }
1836 
1837 
1838 
1839 template <int dim, int spacedim>
1840 inline const hp::FECollection<dim, spacedim> &
1842 {
1843  return fe_collection;
1844 }
1845 
1846 
1847 
1848 template <int dim, int spacedim>
1849 inline const Triangulation<dim, spacedim> &
1851 {
1852  Assert(tria != nullptr,
1853  ExcMessage("This DoFHandler object has not been associated "
1854  "with a triangulation."));
1855  return *tria;
1856 }
1857 
1858 
1859 
1860 template <int dim, int spacedim>
1861 inline MPI_Comm
1863 {
1864  Assert(tria != nullptr,
1865  ExcMessage("This DoFHandler object has not been associated "
1866  "with a triangulation."));
1867  return tria->get_communicator();
1868 }
1869 
1870 
1871 
1872 template <int dim, int spacedim>
1873 inline const BlockInfo &
1875 {
1877 
1878  return block_info_object;
1879 }
1880 
1881 
1882 
1883 template <int dim, int spacedim>
1884 template <typename number>
1887  const std::map<types::boundary_id, const Function<spacedim, number> *>
1888  &boundary_ids) const
1889 {
1890  Assert(!(dim == 2 && spacedim == 3) || this->hp_capability_enabled == false,
1892 
1893  // extract the set of boundary ids and forget about the function object
1894  // pointers
1895  std::set<types::boundary_id> boundary_ids_only;
1896  for (typename std::map<types::boundary_id,
1897  const Function<spacedim, number> *>::const_iterator p =
1898  boundary_ids.begin();
1899  p != boundary_ids.end();
1900  ++p)
1901  boundary_ids_only.insert(p->first);
1902 
1903  // then just hand everything over to the other function that does the work
1904  return n_boundary_dofs(boundary_ids_only);
1905 }
1906 
1907 
1908 
1909 namespace internal
1910 {
1918  template <int dim, int spacedim>
1919  std::string
1920  policy_to_string(const ::internal::DoFHandlerImplementation::Policy::
1921  PolicyBase<dim, spacedim> &policy);
1922 } // namespace internal
1923 
1924 
1925 
1926 template <int dim, int spacedim>
1927 template <class Archive>
1928 void
1929 DoFHandler<dim, spacedim>::save(Archive &ar, const unsigned int) const
1930 {
1931  if (this->hp_capability_enabled)
1932  {
1933  ar &this->object_dof_indices;
1934  ar &this->object_dof_ptr;
1935 
1936  ar &this->cell_dof_cache_indices;
1937  ar &this->cell_dof_cache_ptr;
1938 
1939  ar &this->hp_cell_active_fe_indices;
1940  ar &this->hp_cell_future_fe_indices;
1941 
1942  ar &hp_object_fe_ptr;
1944 
1945  ar &number_cache;
1946 
1947  ar &mg_number_cache;
1948 
1949  // write out the number of triangulation cells and later check during
1950  // loading that this number is indeed correct; same with something that
1951  // identifies the policy
1952  const unsigned int n_cells = this->tria->n_cells();
1953  std::string policy_name =
1955 
1956  ar &n_cells &policy_name;
1957  }
1958  else
1959  {
1960  ar &this->block_info_object;
1961  ar &number_cache;
1962 
1963  ar &this->object_dof_indices;
1964  ar &this->object_dof_ptr;
1965 
1966  ar &this->cell_dof_cache_indices;
1967  ar &this->cell_dof_cache_ptr;
1968 
1969  // write out the number of triangulation cells and later check during
1970  // loading that this number is indeed correct; same with something that
1971  // identifies the FE and the policy
1972  unsigned int n_cells = this->tria->n_cells();
1973  std::string fe_name = this->get_fe(0).get_name();
1974  std::string policy_name = internal::policy_to_string(*this->policy);
1975 
1976  ar &n_cells &fe_name &policy_name;
1977  }
1978 }
1979 
1980 
1981 
1982 template <int dim, int spacedim>
1983 template <class Archive>
1984 void
1985 DoFHandler<dim, spacedim>::load(Archive &ar, const unsigned int)
1986 {
1987  if (this->hp_capability_enabled)
1988  {
1989  ar &this->object_dof_indices;
1990  ar &this->object_dof_ptr;
1991 
1992  ar &this->cell_dof_cache_indices;
1993  ar &this->cell_dof_cache_ptr;
1994 
1995  ar &this->hp_cell_active_fe_indices;
1996  ar &this->hp_cell_future_fe_indices;
1997 
1998  ar &hp_object_fe_ptr;
2000 
2001  ar &number_cache;
2002 
2003  ar &mg_number_cache;
2004 
2005  // these are the checks that correspond to the last block in the save()
2006  // function
2007  unsigned int n_cells;
2008  std::string policy_name;
2009 
2010  ar &n_cells &policy_name;
2011 
2012  AssertThrow(
2013  n_cells == this->tria->n_cells(),
2014  ExcMessage(
2015  "The object being loaded into does not match the triangulation "
2016  "that has been stored previously."));
2017  AssertThrow(
2018  policy_name == ::internal::policy_to_string(*this->policy),
2019  ExcMessage("The policy currently associated with this DoFHandler (" +
2021  ") does not match the one that was associated with the "
2022  "DoFHandler previously stored (" +
2023  policy_name + ")."));
2024  }
2025  else
2026  {
2027  ar &this->block_info_object;
2028  ar &number_cache;
2029 
2030  object_dof_indices.clear();
2031 
2032  object_dof_ptr.clear();
2033 
2034  ar &this->object_dof_indices;
2035  ar &this->object_dof_ptr;
2036 
2037  ar &this->cell_dof_cache_indices;
2038  ar &this->cell_dof_cache_ptr;
2039 
2040  // these are the checks that correspond to the last block in the save()
2041  // function
2042  unsigned int n_cells;
2043  std::string fe_name;
2044  std::string policy_name;
2045 
2046  ar &n_cells &fe_name &policy_name;
2047 
2048  AssertThrow(
2049  n_cells == this->tria->n_cells(),
2050  ExcMessage(
2051  "The object being loaded into does not match the triangulation "
2052  "that has been stored previously."));
2053  AssertThrow(
2054  fe_name == this->get_fe(0).get_name(),
2055  ExcMessage(
2056  "The finite element associated with this DoFHandler does not match "
2057  "the one that was associated with the DoFHandler previously stored."));
2058  AssertThrow(policy_name == internal::policy_to_string(*this->policy),
2059  ExcMessage(
2060  "The policy currently associated with this DoFHandler (" +
2062  ") does not match the one that was associated with the "
2063  "DoFHandler previously stored (" +
2064  policy_name + ")."));
2065  }
2066 }
2067 
2068 
2069 
2070 template <int dim, int spacedim>
2071 inline types::global_dof_index &
2073  const unsigned int level,
2074  const unsigned int dof_number,
2075  const unsigned int dofs_per_vertex)
2076 {
2077  Assert((level >= coarsest_level) && (level <= finest_level),
2078  ExcInvalidLevel(level));
2079  return indices[dofs_per_vertex * (level - coarsest_level) + dof_number];
2080 }
2081 
2082 
2083 
2084 extern template class DoFHandler<1, 1>;
2085 extern template class DoFHandler<1, 2>;
2086 extern template class DoFHandler<1, 3>;
2087 extern template class DoFHandler<2, 2>;
2088 extern template class DoFHandler<2, 3>;
2089 extern template class DoFHandler<3, 3>;
2090 
2091 
2092 #endif // DOXYGEN
2093 
2095 
2096 #endif
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
Definition: block_info.h:94
types::global_dof_index & access_index(const unsigned int level, const unsigned int dof_number, const unsigned int dofs_per_vertex)
unsigned int coarsest_level
Definition: dof_handler.h:1346
void init(const unsigned int coarsest_level, const unsigned int finest_level, const unsigned int dofs_per_vertex)
unsigned int get_finest_level() const
std::unique_ptr< types::global_dof_index[]> indices
Definition: dof_handler.h:1362
unsigned int get_coarsest_level() const
cell_iterator end() const
void pre_distributed_transfer_action()
void post_transfer_action()
virtual std::size_t memory_consumption() const
std::vector< std::vector< active_fe_index_type > > hp_cell_active_fe_indices
Definition: dof_handler.h:1507
unsigned int max_couplings_between_dofs() const
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
Definition: dof_handler.h:1527
std::unique_ptr< ActiveFEIndexTransfer > active_fe_index_transfer
Definition: dof_handler.h:1539
hp::FECollection< dim, spacedim > fe_collection
Definition: dof_handler.h:1431
void create_active_fe_table()
static constexpr unsigned int dimension
Definition: dof_handler.h:509
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
Definition: dof_handler.h:1424
void save(Archive &ar, const unsigned int version) const
void set_active_fe_indices(const std::vector< unsigned int > &active_fe_indices)
BlockInfo block_info_object
Definition: dof_handler.h:1412
bool has_active_dofs() const
std::vector< boost::signals2::connection > tria_listeners_for_transfer
Definition: dof_handler.h:1552
level_cell_iterator end_mg() const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
static constexpr unsigned int space_dimension
Definition: dof_handler.h:514
void renumber_dofs(const unsigned int level, const std::vector< types::global_dof_index > &new_numbers)
level_cell_iterator begin_mg(const unsigned int level=0) const
void distribute_dofs(const FiniteElement< dim, spacedim > &fe)
types::global_dof_index n_boundary_dofs() const
DoFHandler(const Triangulation< dim, spacedim > &tria)
void clear_mg_space()
void connect_to_triangulation_signals()
std::vector< std::vector< types::global_dof_index > > cell_dof_cache_indices
Definition: dof_handler.h:1462
unsigned int offset_type
Definition: dof_handler.h:534
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1520
std::vector< std::vector< active_fe_index_type > > hp_cell_future_fe_indices
Definition: dof_handler.h:1514
std::vector< boost::signals2::connection > tria_listeners
Definition: dof_handler.h:1545
void post_distributed_transfer_action()
std::vector< std::array< std::vector< types::global_dof_index >, dim+1 > > object_dof_indices
Definition: dof_handler.h:1476
typename LevelSelector::CellAccessor level_cell_accessor
Definition: dof_handler.h:499
const BlockInfo & block_info() const
const Triangulation< dim, spacedim > & get_triangulation() const
void pre_transfer_action()
void clear_space()
std::array< std::vector< active_fe_index_type >, dim+1 > hp_object_fe_indices
Definition: dof_handler.h:1495
void reinit(const Triangulation< dim, spacedim > &tria)
std::vector< std::array< std::vector< offset_type >, dim+1 > > object_dof_ptr
Definition: dof_handler.h:1487
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
Definition: dof_handler.h:1533
void serialize(Archive &archive, const unsigned int version)
DoFHandler(const DoFHandler &)=delete
active_cell_iterator begin_active(const unsigned int level=0) const
void distribute_mg_dofs()
active_cell_iterator end_active(const unsigned int level) const
bool hp_capability_enabled
Definition: dof_handler.h:1418
bool has_hp_capabilities() const
void initialize_local_block_info()
types::global_dof_index n_dofs() const
void update_active_fe_table()
level_cell_iterator end_mg(const unsigned int level) const
std::vector< std::vector< offset_type > > cell_dof_cache_ptr
Definition: dof_handler.h:1468
cell_iterator end(const unsigned int level) const
typename LevelSelector::cell_iterator level_cell_iterator
Definition: dof_handler.h:502
void prepare_for_serialization_of_active_fe_indices()
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
Definition: dof_handler.h:1454
void load(Archive &ar, const unsigned int version)
DoFHandler & operator=(const DoFHandler &)=delete
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
static const unsigned int default_fe_index
Definition: dof_handler.h:519
typename LevelSelector::face_iterator level_face_iterator
Definition: dof_handler.h:503
cell_iterator begin(const unsigned int level=0) const
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
Definition: dof_handler.h:1500
MPI_Comm get_communicator() const
static const unsigned int invalid_fe_index
Definition: dof_handler.h:524
void clear()
const IndexSet & locally_owned_dofs() const
void get_active_fe_indices(std::vector< unsigned int > &active_fe_indices) const
void distribute_dofs(const hp::FECollection< dim, spacedim > &fe)
types::global_dof_index n_boundary_dofs(const std::map< types::boundary_id, const Function< spacedim, number > * > &boundary_ids) const
bool has_level_dofs() const
std::unique_ptr<::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > > policy
Definition: dof_handler.h:1439
typename LevelSelector::FaceAccessor level_face_accessor
Definition: dof_handler.h:500
types::global_dof_index n_dofs(const unsigned int level) const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
unsigned int max_couplings_between_boundary_dofs() const
void setup_policy()
types::global_dof_index n_locally_owned_dofs() const
static const active_fe_index_type invalid_active_fe_index
Definition: dof_handler.h:540
unsigned short int active_fe_index_type
Definition: dof_handler.h:529
void deserialize_active_fe_indices()
types::global_dof_index n_boundary_dofs(const std::set< types::boundary_id > &boundary_ids) const
::internal::DoFHandlerImplementation::NumberCache number_cache
Definition: dof_handler.h:1448
virtual ~DoFHandler() override
virtual std::string get_name() const =0
virtual MPI_Comm get_communicator() const
unsigned int n_cells() const
unsigned int size() const
Definition: collection.h:264
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
unsigned int level
Definition: grid_out.cc:4607
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
IteratorRange< active_cell_iterator > active_cell_iterators() const
IteratorRange< level_cell_iterator > mg_cell_iterators() const
IteratorRange< cell_iterator > cell_iterators() const
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcInvalidFEIndex(int arg1, int arg2)
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
static ::ExceptionBase & ExcInvalidLevel(int arg1)
static ::ExceptionBase & ExcNewNumbersNotConsecutive(types::global_dof_index arg1)
#define Assert(cond, exc)
Definition: exceptions.h:1473
static ::ExceptionBase & ExcNoFESelected()
static ::ExceptionBase & ExcOnlyAvailableWithHP()
static ::ExceptionBase & ExcNotImplementedWithHP()
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
typename ActiveSelector::FaceAccessor face_accessor
Definition: dof_handler.h:348
typename ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:466
typename ActiveSelector::line_iterator line_iterator
Definition: dof_handler.h:357
typename ActiveSelector::face_iterator face_iterator
Definition: dof_handler.h:484
typename ActiveSelector::active_line_iterator active_line_iterator
Definition: dof_handler.h:372
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:438
typename ActiveSelector::CellAccessor cell_accessor
Definition: dof_handler.h:334
typename ActiveSelector::quad_iterator quad_iterator
Definition: dof_handler.h:381
typename ActiveSelector::active_hex_iterator active_hex_iterator
Definition: dof_handler.h:416
typename ActiveSelector::active_quad_iterator active_quad_iterator
Definition: dof_handler.h:396
typename ActiveSelector::hex_iterator hex_iterator
Definition: dof_handler.h:405
typename ActiveSelector::active_face_iterator active_face_iterator
Definition: dof_handler.h:497
Definition: hp.h:118
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:13787
unsigned int dominated_future_fe_on_children(const typename DoFHandler< dim, spacedim >::cell_iterator &parent)
void communicate_future_fe_indices(DoFHandler< dim, spacedim > &dof_handler)
std::string policy_to_string(const ::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > &policy)
Definition: dof_handler.cc:58
static const unsigned int invalid_unsigned_int
Definition: types.h:201
unsigned int global_dof_index
Definition: types.h:76
unsigned int boundary_id
Definition: types.h:129
std::map< const cell_iterator, const unsigned int > coarsened_cells_fe_index
Definition: dof_handler.h:1389
std::vector< unsigned int > active_fe_indices
Definition: dof_handler.h:1396
std::map< const cell_iterator, const unsigned int > refined_cells_fe_index
Definition: dof_handler.h:1383
std::unique_ptr< parallel::distributed::CellDataTransfer< dim, spacedim, std::vector< unsigned int > > > cell_data_transfer
Definition: dof_handler.h:1406
std::map< const cell_iterator, const unsigned int > persisting_cells_fe_index
Definition: dof_handler.h:1377