15#ifndef dealii_dof_handler_h
16#define dealii_dof_handler_h
38#include <boost/serialization/split_member.hpp>
39#include <boost/signals2/connection.hpp>
50template <
int dim,
int spacedim>
52template <
int dim,
int spacedim>
58 namespace DoFHandlerImplementation
60 struct Implementation;
64 template <
int dim,
int spacedim>
66 struct Implementation;
70 namespace DoFAccessorImplementation
72 struct Implementation;
75 namespace DoFCellAccessorImplementation
77 struct Implementation;
82 namespace DoFHandlerImplementation
84 struct Implementation;
93 template <
int dim,
int spacedim,
typename VectorType>
98template <
int structdim,
int dim,
int spacedim,
bool level_dof_access>
101template <
int dimension_,
int space_dimension_,
bool level_dof_access>
322template <
int dim,
int spacedim = dim>
608 std::vector<types::fe_index>
657 std::vector<types::fe_index>
823 const std::vector<types::global_dof_index> &new_numbers);
1137 template <
typename number>
1141 &boundary_ids)
const;
1236 DEAL_II_DEPRECATED_EARLY_WITH_COMMENT(
1237 "Access the MPI communicator with get_mpi_communicator() instead.")
1239 get_communicator() const;
1258 prepare_for_serialization_of_active_fe_indices();
1276 deserialize_active_fe_indices();
1287 memory_consumption() const;
1294 template <class Archive>
1296 save(Archive &ar, const
unsigned int version) const;
1303 template <class Archive>
1305 load(Archive &ar, const
unsigned int version);
1313 template <
class Archive>
1319 BOOST_SERIALIZATION_SPLIT_MEMBER()
1337 <<
"The given level " << arg1
1338 <<
" is not in the valid range!");
1345 <<
"The given list of new dof indices is not consecutive: "
1346 <<
"the index " << arg1 <<
" does not exist.");
1353 <<
"The mesh contains a cell with an active FE index of "
1354 << arg1 <<
", but the finite element collection only has "
1355 << arg2 <<
" elements");
1362 "The current function doesn't make sense when used with a "
1363 "DoFHandler without hp-capabilities.");
1370 "The current function has not yet been implemented for a "
1371 "DoFHandler with hp-capabilities.");
1394 init(
const unsigned int coarsest_level,
1395 const unsigned int finest_level,
1396 const unsigned int dofs_per_vertex);
1416 const unsigned int dof_number,
1417 const unsigned int dofs_per_vertex);
1439 std::unique_ptr<types::global_dof_index[]>
indices;
1454 std::map<const cell_iterator, const types::fe_index>
1467 std::map<const cell_iterator, const types::fe_index>
1483 parallel::distributed::
1484 CellDataTransfer<dim, spacedim, std::vector<types::fe_index>>>
1516 std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1517 PolicyBase<dim, spacedim>>
1532 std::vector<::internal::DoFHandlerImplementation::NumberCache>
1540 mutable std::vector<std::array<std::vector<types::global_dof_index>, dim + 1>>
1551 mutable std::vector<std::array<std::vector<offset_type>, dim + 1>>
1559 mutable std::array<std::vector<types::fe_index>, dim + 1>
1589 std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel<dim>>>
1595 std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces<dim>>
1717 template <
int,
int,
int,
bool>
1718 friend class ::DoFAccessor;
1719 template <
int,
int,
bool>
1720 friend class ::DoFCellAccessor;
1721 friend struct ::internal::DoFAccessorImplementation::Implementation;
1722 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1726 friend struct ::internal::DoFHandlerImplementation::Implementation;
1727 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1728 friend struct ::internal::DoFHandlerImplementation::Policy::
1734 static_assert(dim <= spacedim,
1735 "The dimension <dim> of a DoFHandler must be less than or "
1736 "equal to the space dimension <spacedim> in which it lives.");
1744 namespace DoFHandlerImplementation
1757 template <
int dim,
int spacedim>
1785 template <
int dim,
int spacedim = dim>
1795 "No FiniteElement has been found in your FECollection that is "
1796 "dominated by all children of a cell you are trying to coarsen!");
1744 namespace DoFHandlerImplementation {
…}
1807template <
int dim,
int spacedim>
1811 return hp_capability_enabled;
1816template <
int dim,
int spacedim>
1820 return mg_number_cache.size() > 0;
1825template <
int dim,
int spacedim>
1829 return number_cache.n_global_dofs > 0;
1834template <
int dim,
int spacedim>
1838 return number_cache.n_global_dofs;
1843template <
int dim,
int spacedim>
1850 "n_dofs(level) can only be called after distribute_mg_dofs()"));
1852 return mg_number_cache[
level].n_global_dofs;
1857template <
int dim,
int spacedim>
1861 return number_cache.n_locally_owned_dofs;
1866template <
int dim,
int spacedim>
1870 return number_cache.locally_owned_dofs;
1875template <
int dim,
int spacedim>
1878 const unsigned int level)
const
1880 Assert(level < this->get_triangulation().n_global_levels(),
1881 ExcMessage(
"The given level index exceeds the number of levels "
1882 "present in the triangulation"));
1884 mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1886 "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1887 return mg_number_cache[
level].locally_owned_dofs;
1892template <
int dim,
int spacedim>
1897 Assert(fe_collection.size() > 0,
1898 ExcMessage(
"No finite element collection is associated with "
1899 "this DoFHandler"));
1900 return fe_collection[number];
1905template <
int dim,
int spacedim>
1910 return fe_collection;
1915template <
int dim,
int spacedim>
1921 ExcMessage(
"This DoFHandler object has not been associated "
1922 "with a triangulation."));
1928template <
int dim,
int spacedim>
1933 ExcMessage(
"This DoFHandler object has not been associated "
1934 "with a triangulation."));
1935 return tria->get_mpi_communicator();
1940template <
int dim,
int spacedim>
1944 return get_mpi_communicator();
1949template <
int dim,
int spacedim>
1953 Assert(this->hp_capability_enabled ==
false, ExcNotImplementedWithHP());
1955 return block_info_object;
1960template <
int dim,
int spacedim>
1962template <
typename number>
1965 &boundary_ids)
const
1967 Assert(!(dim == 2 && spacedim == 3) || this->hp_capability_enabled ==
false,
1968 ExcNotImplementedWithHP());
1972 std::set<types::boundary_id> boundary_ids_only;
1975 boundary_ids.begin();
1976 p != boundary_ids.end();
1978 boundary_ids_only.insert(p->first);
1981 return n_boundary_dofs(boundary_ids_only);
1995 template <
int dim,
int spacedim>
1998 PolicyBase<dim, spacedim> &policy);
2003template <
int dim,
int spacedim>
2005template <
class Archive>
2008 if (this->hp_capability_enabled)
2010 ar &this->object_dof_indices;
2011 ar &this->object_dof_ptr;
2013 ar &this->hp_cell_active_fe_indices;
2014 ar &this->hp_cell_future_fe_indices;
2016 ar &hp_object_fe_ptr;
2017 ar &hp_object_fe_indices;
2021 ar &mg_number_cache;
2026 const unsigned int n_cells = this->tria->n_cells();
2027 std::string policy_name =
2034 ar &this->block_info_object;
2037 ar &this->object_dof_indices;
2038 ar &this->object_dof_ptr;
2043 unsigned int n_cells = this->tria->n_cells();
2044 std::string fe_name = this->get_fe(0).get_name();
2047 ar &
n_cells &fe_name &policy_name;
2053template <
int dim,
int spacedim>
2055template <
class Archive>
2058 if (this->hp_capability_enabled)
2060 ar &this->object_dof_indices;
2061 ar &this->object_dof_ptr;
2063 ar &this->hp_cell_active_fe_indices;
2064 ar &this->hp_cell_future_fe_indices;
2066 ar &hp_object_fe_ptr;
2067 ar &hp_object_fe_indices;
2071 ar &mg_number_cache;
2076 std::string policy_name;
2081 n_cells == this->tria->n_cells(),
2083 "The object being loaded into does not match the triangulation "
2084 "that has been stored previously."));
2087 ExcMessage(
"The policy currently associated with this DoFHandler (" +
2089 ") does not match the one that was associated with the "
2090 "DoFHandler previously stored (" +
2091 policy_name +
")."));
2095 ar &this->block_info_object;
2098 object_dof_indices.clear();
2100 object_dof_ptr.clear();
2102 ar &this->object_dof_indices;
2103 ar &this->object_dof_ptr;
2108 std::string fe_name;
2109 std::string policy_name;
2111 ar &
n_cells &fe_name &policy_name;
2114 n_cells == this->tria->n_cells(),
2116 "The object being loaded into does not match the triangulation "
2117 "that has been stored previously."));
2119 fe_name == this->get_fe(0).get_name(),
2121 "The finite element associated with this DoFHandler does not match "
2122 "the one that was associated with the DoFHandler previously stored."));
2125 "The policy currently associated with this DoFHandler (" +
2127 ") does not match the one that was associated with the "
2128 "DoFHandler previously stored (" +
2129 policy_name +
")."));
2135template <
int dim,
int spacedim>
2139 const unsigned int level,
2140 const unsigned int dof_number,
2141 const unsigned int dofs_per_vertex)
A small class collecting the different BlockIndices involved in global, multilevel and local computat...
static constexpr unsigned int space_dimension
static constexpr unsigned int dimension
static const unsigned int dim
static const unsigned int spacedim
unsigned int coarsest_level
types::global_dof_index & access_index(const unsigned int level, const unsigned int dof_number, const unsigned int dofs_per_vertex)
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
unsigned int finest_level
unsigned int get_coarsest_level() const
cell_iterator end() const
void pre_distributed_transfer_action()
void post_transfer_action()
std::vector< types::fe_index > get_future_fe_indices() const
unsigned int max_couplings_between_dofs() const
std::vector< std::unique_ptr<::internal::DoFHandlerImplementation::DoFLevel< dim > > > mg_levels
std::unique_ptr< ActiveFEIndexTransfer > active_fe_index_transfer
hp::FECollection< dim, spacedim > fe_collection
void create_active_fe_table()
std::vector< std::vector< types::fe_index > > hp_cell_future_fe_indices
const hp::FECollection< dim, spacedim > & get_fe_collection() const
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
bool has_active_dofs() const
std::vector< boost::signals2::connection > tria_listeners_for_transfer
level_cell_iterator end_mg() const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
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 connect_to_triangulation_signals()
std::vector< std::vector< types::fe_index > > hp_cell_active_fe_indices
std::vector< MGVertexDoFs > mg_vertex_dofs
std::vector< boost::signals2::connection > tria_listeners
void post_distributed_transfer_action()
std::vector< std::array< std::vector< types::global_dof_index >, dim+1 > > object_dof_indices
const Triangulation< dim, spacedim > & get_triangulation() const
typename LevelSelector::CellAccessor level_cell_accessor
void pre_transfer_action()
const IndexSet & locally_owned_dofs() const
void reinit(const Triangulation< dim, spacedim > &tria)
std::vector< std::array< std::vector< offset_type >, dim+1 > > object_dof_ptr
std::unique_ptr<::internal::DoFHandlerImplementation::DoFFaces< dim > > mg_faces
void serialize(Archive &archive, const unsigned int version)
DoFHandler & operator=(const DoFHandler &)=delete
DoFHandler(const DoFHandler &)=delete
active_cell_iterator begin_active(const unsigned int level=0) const
std::vector< types::fe_index > get_active_fe_indices() const
void distribute_mg_dofs()
active_cell_iterator end_active(const unsigned int level) const
bool hp_capability_enabled
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
cell_iterator end(const unsigned int level) const
typename LevelSelector::cell_iterator level_cell_iterator
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
void load(Archive &ar, const unsigned int version)
typename LevelSelector::face_iterator level_face_iterator
cell_iterator begin(const unsigned int level=0) const
MPI_Comm get_mpi_communicator() const
std::array< std::vector< offset_type >, dim+1 > hp_object_fe_ptr
MPI_Comm get_communicator() 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
ObserverPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
bool has_level_dofs() const
std::unique_ptr<::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > > policy
typename LevelSelector::FaceAccessor level_face_accessor
types::global_dof_index n_dofs(const unsigned int level) const
std::array< std::vector< types::fe_index >, dim+1 > hp_object_fe_indices
void set_active_fe_indices(const std::vector< types::fe_index > &active_fe_indices)
unsigned int max_couplings_between_boundary_dofs() const
types::global_dof_index n_locally_owned_dofs() const
types::global_dof_index n_boundary_dofs(const std::set< types::boundary_id > &boundary_ids) const
const BlockInfo & block_info() const
::internal::DoFHandlerImplementation::NumberCache number_cache
void set_future_fe_indices(const std::vector< types::fe_index > &future_fe_indices)
virtual ~DoFHandler() override
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_CXX20_REQUIRES(condition)
#define DEAL_II_NAMESPACE_CLOSE
IteratorRange< active_cell_iterator > active_cell_iterators() const
IteratorRange< level_cell_iterator > mg_cell_iterators() const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
IteratorRange< cell_iterator > cell_iterators() const
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
IteratorRange< level_cell_iterator > mg_cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
#define Assert(cond, exc)
#define DeclException2(Exception2, type1, type2, outsequence)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcNoDominatedFiniteElementOnChildren()
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcInvalidLevel(int arg1)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
typename ActiveSelector::FaceAccessor face_accessor
typename ActiveSelector::cell_iterator cell_iterator
typename ActiveSelector::line_iterator line_iterator
typename ActiveSelector::face_iterator face_iterator
typename ActiveSelector::active_line_iterator active_line_iterator
typename ActiveSelector::active_cell_iterator active_cell_iterator
typename ActiveSelector::CellAccessor cell_accessor
typename ActiveSelector::quad_iterator quad_iterator
typename ActiveSelector::active_hex_iterator active_hex_iterator
typename ActiveSelector::active_quad_iterator active_quad_iterator
typename ActiveSelector::hex_iterator hex_iterator
typename ActiveSelector::active_face_iterator active_face_iterator
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
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)
unsigned short int fe_index
std::map< const cell_iterator, const types::fe_index > refined_cells_fe_index
std::map< const cell_iterator, const types::fe_index > persisting_cells_fe_index
std::vector< types::fe_index > active_fe_indices
std::map< const cell_iterator, const types::fe_index > coarsened_cells_fe_index
std::unique_ptr< parallel::distributed::CellDataTransfer< dim, spacedim, std::vector< types::fe_index > > > cell_data_transfer