17 #ifndef __deal2__hp_dof_handler_h
18 #define __deal2__hp_dof_handler_h
22 #include <deal.II/base/config.h>
24 #include <deal.II/base/template_constraints.h>
25 #include <deal.II/base/smartpointer.h>
26 #include <deal.II/dofs/function_map.h>
27 #include <deal.II/dofs/dof_iterator_selector.h>
28 #include <deal.II/dofs/number_cache.h>
29 #include <deal.II/hp/fe_collection.h>
42 template <
int>
class DoFIndicesOnFaces;
43 template <
int>
class DoFIndicesOnFacesOrEdges;
47 struct Implementation;
56 struct Implementation;
61 struct Implementation;
80 template <
int dim,
int spacedim=dim>
83 typedef ::internal::DoFHandler::Iterators<DoFHandler<dim,spacedim>,
false> ActiveSelector;
84 typedef ::internal::DoFHandler::Iterators<DoFHandler<dim,spacedim>,
true> LevelSelector;
86 typedef typename ActiveSelector::CellAccessor cell_accessor;
87 typedef typename ActiveSelector::FaceAccessor face_accessor;
89 typedef typename ActiveSelector::line_iterator line_iterator;
90 typedef typename ActiveSelector::active_line_iterator active_line_iterator;
92 typedef typename ActiveSelector::quad_iterator quad_iterator;
93 typedef typename ActiveSelector::active_quad_iterator active_quad_iterator;
95 typedef typename ActiveSelector::hex_iterator hex_iterator;
96 typedef typename ActiveSelector::active_hex_iterator active_hex_iterator;
98 typedef typename ActiveSelector::active_cell_iterator active_cell_iterator;
100 typedef typename ActiveSelector::face_iterator face_iterator;
101 typedef typename ActiveSelector::active_face_iterator active_face_iterator;
103 typedef typename LevelSelector::CellAccessor level_cell_accessor;
104 typedef typename LevelSelector::FaceAccessor level_face_accessor;
106 typedef typename LevelSelector::cell_iterator level_cell_iterator;
107 typedef typename LevelSelector::face_iterator level_face_iterator;
109 typedef level_cell_iterator cell_iterator;
229 virtual void clear ();
262 void renumber_dofs (
const std::vector<types::global_dof_index> &new_numbers);
315 cell_iterator
begin (
const unsigned int level = 0)
const;
329 active_cell_iterator
begin_active(
const unsigned int level = 0)
const;
338 cell_iterator
end ()
const;
347 cell_iterator
end (
const unsigned int level)
const;
354 active_cell_iterator
end_active (
const unsigned int level)
const;
436 n_boundary_dofs (
const std::set<types::boundary_id> &boundary_indicators)
const;
497 const std::vector<IndexSet> &
526 const std::vector<types::global_dof_index> &
582 <<
"The matrix has the wrong dimension " << arg1);
592 <<
"The given list of new dof indices is not consecutive: "
593 <<
"the index " << arg1 <<
" does not exist.");
599 <<
"The mesh contains a cell with an active_fe_index of "
600 << arg1 <<
", but the finite element collection only has "
601 << arg2 <<
" elements");
607 <<
"The given level " << arg1
608 <<
" is not in the valid range!");
619 <<
"You tried to do something on level " << arg1
620 <<
", but this level is empty.");
685 template<
int structdim>
686 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;
688 template<
int structdim>
689 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;
720 void post_refinement_action ();
793 std::vector<::internal::hp::DoFLevel *>
levels;
865 std::vector<MGVertexDoFs> mg_vertex_dofs;
891 template <
int,
class,
bool>
friend class ::DoFAccessor;
892 template <
class,
bool>
friend class ::DoFCellAccessor;
893 friend struct ::internal::DoFAccessor::Implementation;
894 friend struct ::internal::DoFCellAccessor::Implementation;
903 template <
int>
friend class ::internal::hp::DoFIndicesOnFacesOrEdges;
904 friend struct ::internal::hp::DoFHandler::Implementation;
914 template <
int dim,
int spacedim>
919 return number_cache.n_global_dofs;
923 template <
int dim,
int spacedim>
932 template <
int dim,
int spacedim>
936 return number_cache.n_locally_owned_dofs;
940 template <
int dim,
int spacedim>
944 return number_cache.locally_owned_dofs;
948 template <
int dim,
int spacedim>
949 const std::vector<types::global_dof_index> &
952 return number_cache.n_locally_owned_dofs_per_processor;
956 template <
int dim,
int spacedim>
957 const std::vector<IndexSet> &
960 return number_cache.locally_owned_dofs_per_processor;
965 template<
int dim,
int spacedim>
970 Assert (finite_elements != 0,
971 ExcMessage (
"No finite element collection is associated with "
973 return *finite_elements;
977 template<
int dim,
int spacedim>
985 template<
int dim,
int spacedim>
992 template<
int dim,
int spacedim>
999 template<
int dim,
int spacedim>
1002 const unsigned int)
const
1008 template<
int dim,
int spacedim>
1022 DEAL_II_NAMESPACE_CLOSE
types::global_dof_index n_locally_owned_dofs() const
active_cell_iterator begin_active(const unsigned int level=0) const
virtual void distribute_dofs(const hp::FECollection< dim, spacedim > &fe)
static const unsigned int invalid_unsigned_int
active_cell_iterator end_active(const unsigned int level) const
void renumber_dofs_internal(const std::vector< types::global_dof_index > &new_numbers,::internal::int2type< 0 >)
const Triangulation< dim, spacedim > & get_tria() const
cell_iterator end() const
DoFHandler & operator=(const DoFHandler &)
void compute_quad_dof_identities(std::vector< types::global_dof_index > &new_dof_indices) const
std::map< types::boundary_id, const Function< dim > * > type
::ExceptionBase & ExcMessage(std::string arg1)
DoFHandler(const Triangulation< dim, spacedim > &tria)
DeclException2(ExcInvalidFEIndex, int, int,<< "The mesh contains a cell with an active_fe_index of "<< arg1<< ", but the finite element collection only has "<< arg2<< " elements")
virtual std::size_t memory_consumption() const
const std::vector< types::global_dof_index > & n_locally_owned_dofs_per_processor() const
unsigned int max_couplings_between_dofs() const
::internal::DoFHandler::NumberCache number_cache
std::vector<::internal::hp::DoFLevel * > levels
void compute_line_dof_identities(std::vector< types::global_dof_index > &new_dof_indices) const
void compute_vertex_dof_identities(std::vector< types::global_dof_index > &new_dof_indices) const
const IndexSet & locally_owned_dofs() const
unsigned int max_couplings_between_boundary_dofs() const
void create_active_fe_table()
FunctionMap< spacedim >::type FunctionMap
unsigned int global_dof_index
cell_iterator begin(const unsigned int level=0) const
#define Assert(cond, exc)
void get_active_fe_indices(std::vector< unsigned int > &active_fe_indices) const
static const unsigned int dimension
SmartPointer< const hp::FECollection< dim, spacedim >, hp::DoFHandler< dim, spacedim > > finite_elements
void pre_refinement_action()
static const types::global_dof_index invalid_dof_index
const std::vector< IndexSet > & locally_owned_dofs_per_processor() const
DeclException1(ExcMatrixHasWrongSize, int,<< "The matrix has the wrong dimension "<< arg1)
std::vector< types::global_dof_index > vertex_dofs_offsets
types::global_dof_index n_dofs() const
static const unsigned int space_dimension
void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers)
std::vector< boost::signals2::connection > tria_listeners
SmartPointer< const Triangulation< dim, spacedim >, DoFHandler< dim, spacedim > > tria
static const unsigned int default_fe_index
DeclException0(ExcInvalidTriangulation)
::ExceptionBase & ExcNotImplemented()
std::vector< types::global_dof_index > vertex_dofs
const types::global_dof_index invalid_dof_index
std::vector< std::vector< bool > * > has_children
const hp::FECollection< dim, spacedim > & get_fe() const
::internal::hp::DoFIndicesOnFaces< dim > * faces
types::global_dof_index n_boundary_dofs() const
void set_active_fe_indices(const std::vector< unsigned int > &active_fe_indices)