|
Reference documentation for deal.II version 9.2.0
|
\(\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\}}\)
Go to the documentation of this file.
16 #ifndef dealii_hp_dof_handler_h
17 #define dealii_hp_dof_handler_h
47 template <
int dim,
int spacedim>
54 template <
int dim,
int spacedim,
typename VectorType>
61 namespace DoFHandlerImplementation
63 struct Implementation;
67 template <
int dim,
int spacedim>
69 struct Implementation;
77 namespace DoFHandlerImplementation
79 struct Implementation;
86 namespace DoFAccessorImplementation
88 struct Implementation;
91 namespace DoFCellAccessorImplementation
93 struct Implementation;
202 template <
int dim,
int spacedim = dim>
206 Iterators<DoFHandler<dim, spacedim>,
false>;
208 Iterators<DoFHandler<dim, spacedim>,
true>;
539 renumber_dofs(
const std::vector<types::global_dof_index> &new_numbers);
617 end(
const unsigned int level)
const;
782 template <
typename number>
786 &boundary_ids)
const;
794 n_boundary_dofs(
const std::set<types::boundary_id> &boundary_ids)
const;
890 get_fe(
const unsigned int index)
const;
958 template <
class Archive>
960 save(Archive &ar,
const unsigned int version)
const;
966 template <
class Archive>
968 load(Archive &ar,
const unsigned int version);
975 template <
class Archive>
977 serialize(Archive &archive,
const unsigned int version);
981 BOOST_SERIALIZATION_SPLIT_MEMBER()
1001 <<
"The matrix has the wrong dimension " << arg1);
1011 <<
"The given list of new dof indices is not consecutive: "
1012 <<
"the index " << arg1 <<
" does not exist.");
1019 <<
"The mesh contains a cell with an active_fe_index of "
1020 << arg1 <<
", but the finite element collection only has "
1021 << arg2 <<
" elements");
1027 <<
"The given level " << arg1
1028 <<
" is not in the valid range!");
1038 <<
"You tried to do something on level " << arg1
1039 <<
", but this level is empty.");
1057 std::unique_ptr<::internal::DoFHandlerImplementation::Policy::
1058 PolicyBase<dim, spacedim>>
1074 template <
int structdim>
1077 const unsigned int obj_index,
1078 const unsigned int fe_index,
1079 const unsigned int local_index)
const;
1081 template <
int structdim>
1084 const unsigned int obj_index,
1085 const unsigned int fe_index,
1086 const unsigned int local_index,
1180 std::vector<std::unique_ptr<::internal::hp::DoFLevel>>
levels;
1186 std::unique_ptr<::internal::hp::DoFIndicesOnFaces<dim>>
faces;
1202 std::vector<::internal::DoFHandlerImplementation::NumberCache>
1245 std::map<const cell_iterator, const unsigned int>
1258 std::map<const cell_iterator, const unsigned int>
1274 parallel::distributed::
1275 CellDataTransfer<dim, spacedim, std::vector<unsigned int>>>
1292 template <
int,
class,
bool>
1293 friend class ::DoFAccessor;
1294 template <
class,
bool>
1295 friend class ::DoFCellAccessor;
1296 friend struct ::internal::DoFAccessorImplementation::Implementation;
1297 friend struct ::internal::DoFCellAccessorImplementation::
1303 friend class ::internal::hp::DoFIndicesOnFacesOrEdges;
1304 friend struct ::internal::hp::DoFHandlerImplementation::
1306 friend struct ::internal::DoFHandlerImplementation::Policy::
1319 template <
int dim,
int spacedim>
1320 template <
typename number>
1324 &boundary_ids)
const
1328 std::set<types::boundary_id> boundary_ids_only;
1331 p = boundary_ids.begin();
1332 p != boundary_ids.end();
1334 boundary_ids_only.insert(p->first);
1337 return n_boundary_dofs(boundary_ids_only);
1353 template <
typename number>
1383 template <
int dim,
int spacedim>
1386 PolicyBase<dim, spacedim> &policy);
1392 template <
int dim,
int spacedim>
1393 template <
class Archive>
1398 ar &vertex_dof_offsets;
1400 ar &mg_number_cache;
1405 const unsigned int n_levels = levels.size();
1407 for (
unsigned int i = 0; i < n_levels; ++i)
1413 bool faces_is_nullptr = (faces.get() ==
nullptr);
1414 ar & faces_is_nullptr;
1415 if (!faces_is_nullptr)
1421 const unsigned int n_cells = tria->n_cells();
1429 template <
int dim,
int spacedim>
1430 template <
class Archive>
1435 ar &vertex_dof_offsets;
1437 ar &mg_number_cache;
1451 levels.resize(size);
1452 for (
unsigned int i = 0; i < size; ++i)
1454 std::unique_ptr<::internal::hp::DoFLevel>
level;
1456 levels[i] = std::move(
level);
1460 bool faces_is_nullptr =
true;
1461 ar & faces_is_nullptr;
1462 if (!faces_is_nullptr)
1468 std::string policy_name;
1475 "The object being loaded into does not match the triangulation "
1476 "that has been stored previously."));
1479 "The policy currently associated with this DoFHandler (" +
1481 ") does not match the one that was associated with the "
1482 "DoFHandler previously stored (" +
1483 policy_name +
")."));
1486 template <
int dim,
int spacedim>
1490 return number_cache.n_global_dofs;
1495 template <
int dim,
int spacedim>
1505 template <
int dim,
int spacedim>
1509 return number_cache.n_locally_owned_dofs;
1514 template <
int dim,
int spacedim>
1518 return number_cache.locally_owned_dofs;
1523 template <
int dim,
int spacedim>
1524 const std::vector<types::global_dof_index> &
1527 if (number_cache.n_locally_owned_dofs_per_processor.empty() &&
1528 number_cache.n_global_dofs > 0)
1534 &this->get_triangulation()));
1538 comm = MPI_COMM_SELF;
1542 .n_locally_owned_dofs_per_processor =
1545 return number_cache.n_locally_owned_dofs_per_processor;
1550 template <
int dim,
int spacedim>
1551 const std::vector<IndexSet> &
1554 if (number_cache.locally_owned_dofs_per_processor.empty() &&
1555 number_cache.n_global_dofs > 0)
1561 &this->get_triangulation()));
1565 comm = MPI_COMM_SELF;
1569 .locally_owned_dofs_per_processor =
1572 return number_cache.locally_owned_dofs_per_processor;
1577 template <
int dim,
int spacedim>
1580 const unsigned int level)
const
1584 Assert(level < this->get_triangulation().n_global_levels(),
1585 ExcMessage(
"The given level index exceeds the number of levels "
1586 "present in the triangulation"));
1587 return mg_number_cache[0].locally_owned_dofs;
1592 template <
int dim,
int spacedim>
1593 const std::vector<IndexSet> &
1595 const unsigned int level)
const
1597 Assert(level < this->get_triangulation().n_global_levels(),
1598 ExcMessage(
"The given level index exceeds the number of levels "
1599 "present in the triangulation"));
1601 mg_number_cache.size() == this->get_triangulation().n_global_levels(),
1603 "The level dofs are not set up properly! Did you call distribute_mg_dofs()?"));
1604 if (mg_number_cache[
level].locally_owned_dofs_per_processor.empty() &&
1605 mg_number_cache[
level].n_global_dofs > 0)
1611 &this->get_triangulation()));
1615 comm = MPI_COMM_SELF;
1618 mg_number_cache[
level])
1619 .locally_owned_dofs_per_processor =
1620 mg_number_cache[
level].get_locally_owned_dofs_per_processor(comm);
1627 template <
int dim,
int spacedim>
1631 Assert(fe_collection.size() > 0,
1632 ExcMessage(
"No finite element collection is associated with "
1633 "this DoFHandler"));
1634 return fe_collection[number];
1639 template <
int dim,
int spacedim>
1643 Assert(fe_collection.size() > 0,
1644 ExcMessage(
"No finite element collection is associated with "
1645 "this DoFHandler"));
1646 return fe_collection;
1651 template <
int dim,
int spacedim>
1656 ExcMessage(
"This DoFHandler object has not been associated "
1657 "with a triangulation."));
virtual ~DoFHandler() override
const IndexSet & locally_owned_dofs() const
typename ActiveSelector::CellAccessor cell_accessor
const types::global_dof_index invalid_dof_index
void prepare_for_serialization_of_active_fe_indices()
static const bool is_hp_dof_handler
IteratorRange< active_cell_iterator > active_cell_iterators_on_level(const unsigned int level) const
typename LevelSelector::CellAccessor level_cell_accessor
typename LevelSelector::face_iterator level_face_iterator
static ::ExceptionBase & ExcMatrixHasWrongSize(int arg1)
unsigned int max_couplings_between_boundary_dofs() const
void post_active_fe_index_transfer()
types::global_dof_index n_boundary_dofs() const
void set_dof_index(const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index, const types::global_dof_index global_index) const
static ::ExceptionBase & ExcGridsDoNotMatch()
static ::ExceptionBase & ExcNotImplemented()
typename ActiveSelector::active_cell_iterator active_cell_iterator
const hp::FECollection< dim, spacedim > & get_fe_collection() const
void initialize(const Triangulation< dim, spacedim > &tria, const hp::FECollection< dim, spacedim > &fe)
static ::ExceptionBase & ExcNewNumbersNotConsecutive(types::global_dof_index arg1)
typename ActiveSelector::active_hex_iterator active_hex_iterator
typename ActiveSelector::line_iterator line_iterator
static ::ExceptionBase & ExcInvalidFEIndex(int arg1, int arg2)
typename ActiveSelector::active_quad_iterator active_quad_iterator
hp::FECollection< dim, spacedim > fe_collection
active_cell_iterator end_active(const unsigned int level) const
const Triangulation< dim, spacedim > & get_triangulation() const
typename ActiveSelector::active_line_iterator active_line_iterator
void post_distributed_active_fe_index_transfer()
static ::ExceptionBase & ExcNoFESelected()
typename ActiveSelector::quad_iterator quad_iterator
static ::ExceptionBase & ExcFunctionNotUseful()
void load(Archive &ar, const unsigned int version)
unsigned int max_couplings_between_dofs() const
types::global_dof_index n_dofs() const
void serialize(Archive &archive, const unsigned int version)
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index) const
std::map< const cell_iterator, const unsigned int > coarsened_cells_fe_index
void create_active_fe_table()
std::map< const cell_iterator, const unsigned int > refined_cells_fe_index
typename LevelSelector::FaceAccessor level_face_accessor
virtual void set_fe(const hp::FECollection< dim, spacedim > &fe)
std::vector< unsigned int > active_fe_indices
void pre_distributed_active_fe_index_transfer()
IteratorRange< active_cell_iterator > active_cell_iterators() const
void deserialize_active_fe_indices()
std::vector< unsigned int > vertex_dof_offsets
DoFHandler & operator=(const DoFHandler &)=delete
std::unique_ptr< parallel::distributed::CellDataTransfer< dim, spacedim, std::vector< unsigned int > > > cell_data_transfer
cell_iterator begin(const unsigned int level=0) const
static ::ExceptionBase & ExcMessage(std::string arg1)
std::vector< types::global_dof_index > vertex_dofs
static ::ExceptionBase & ExcEmptyLevel(int arg1)
void get_active_fe_indices(std::vector< unsigned int > &active_fe_indices) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
static ::ExceptionBase & ExcInvalidLevel(int arg1)
types::global_dof_index n_locally_owned_dofs() const
cell_iterator end() const
const std::vector< IndexSet > & locally_owned_mg_dofs_per_processor(const unsigned int level) const
virtual void distribute_dofs(const hp::FECollection< dim, spacedim > &fe)
#define DEAL_II_NAMESPACE_OPEN
void post_refinement_action()
#define DEAL_II_DEPRECATED
#define DeclException1(Exception1, type1, outsequence)
::internal::DoFHandlerImplementation::NumberCache number_cache
IteratorRange< cell_iterator > cell_iterators() const
typename ActiveSelector::cell_iterator cell_iterator
std::vector< std::unique_ptr<::internal::hp::DoFLevel > > levels
std::vector< types::global_dof_index > get_n_locally_owned_dofs_per_processor(const MPI_Comm mpi_communicator) const
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
types::global_dof_index get_dof_index(const unsigned int obj_level, const unsigned int obj_index, const unsigned int fe_index, const unsigned int local_index) const
active_cell_iterator begin_active(const unsigned int level=0) const
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
std::unique_ptr<::internal::hp::DoFIndicesOnFaces< dim > > faces
virtual std::size_t memory_consumption() const
std::vector<::internal::DoFHandlerImplementation::NumberCache > mg_number_cache
void set_active_fe_indices(const std::vector< unsigned int > &active_fe_indices)
virtual const MPI_Comm & get_communicator() const
#define Assert(cond, exc)
static ::ExceptionBase & ExcFacesHaveNoLevel()
std::unique_ptr< ActiveFEIndexTransfer > active_fe_index_transfer
typename LevelSelector::cell_iterator level_cell_iterator
std::vector< boost::signals2::connection > tria_listeners
const std::vector< types::global_dof_index > & n_locally_owned_dofs_per_processor() const
typename ActiveSelector::FaceAccessor face_accessor
#define DeclException0(Exception0)
static ::ExceptionBase & ExcInvalidBoundaryIndicator()
void pre_refinement_action()
static const unsigned int invalid_unsigned_int
void setup_policy_and_listeners()
std::string policy_to_string(const ::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > &policy)
static const unsigned int default_fe_index
void post_distributed_serialization_of_active_fe_indices()
std::unique_ptr<::internal::DoFHandlerImplementation::Policy::PolicyBase< dim, spacedim > > policy
typename ActiveSelector::hex_iterator hex_iterator
typename ActiveSelector::face_iterator face_iterator
std::vector< IndexSet > get_locally_owned_dofs_per_processor(const MPI_Comm mpi_communicator) const
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
#define DEAL_II_NAMESPACE_CLOSE
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
static const unsigned int dimension
std::vector< IndexSet > locally_owned_dofs_per_processor
void pre_active_fe_index_transfer()
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
#define DeclException2(Exception2, type1, type2, outsequence)
#define AssertThrow(cond, exc)
typename ActiveSelector::active_face_iterator active_face_iterator
const std::vector< IndexSet > & locally_owned_dofs_per_processor() const
std::map< const cell_iterator, const unsigned int > persisting_cells_fe_index
void save(Archive &ar, const unsigned int version) const
static const unsigned int space_dimension