16 #ifndef dealii_dof_accessor_h 17 #define dealii_dof_accessor_h 20 #include <deal.II/base/config.h> 22 #include <deal.II/dofs/dof_handler.h> 24 #include <deal.II/grid/tria_accessor.h> 26 #include <deal.II/hp/dof_handler.h> 30 DEAL_II_NAMESPACE_OPEN
32 template <
typename number>
34 template <
typename number>
36 template <
typename number>
39 template <
typename Accessor>
48 namespace DoFCellAccessorImplementation
50 struct Implementation;
53 namespace DoFHandlerImplementation
55 struct Implementation;
58 struct Implementation;
64 namespace DoFHandlerImplementation
66 struct Implementation;
79 namespace DoFAccessorImplementation
98 template <
int structdim,
int dim,
int spacedim>
114 template <
int dim,
int spacedim>
208 template <
int structdim,
typename DoFHandlerType,
bool level_dof_access>
212 DoFHandlerType::dimension,
220 static const unsigned int dimension = DoFHandlerType::dimension;
226 static const unsigned int space_dimension = DoFHandlerType::space_dimension;
232 using BaseClass = typename ::internal::DoFAccessorImplementation::
233 Inheritance<structdim, dimension, space_dimension>::BaseClass;
269 DoFHandlerType::space_dimension> *tria,
272 const DoFHandlerType *dof_handler);
286 template <
int structdim2,
int dim2,
int spacedim2>
293 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
299 template <
bool level_dof_access2>
323 const DoFHandlerType &
324 get_dof_handler()
const;
329 template <
bool level_dof_access2>
339 DoFHandlerType::dimension,
340 DoFHandlerType::space_dimension> &da);
360 child(
const unsigned int c)
const;
367 typename ::internal::DoFHandlerImplementation::
368 Iterators<DoFHandlerType, level_dof_access>::line_iterator
369 line(
const unsigned int i)
const;
376 typename ::internal::DoFHandlerImplementation::
377 Iterators<DoFHandlerType, level_dof_access>::quad_iterator
378 quad(
const unsigned int i)
const;
427 std::vector<types::global_dof_index> &dof_indices,
428 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
439 std::vector<types::global_dof_index> &dof_indices,
440 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
448 const std::vector<types::global_dof_index> &dof_indices,
449 const unsigned int fe_index = DoFHandlerType::default_fe_index);
470 const unsigned int vertex,
471 const unsigned int i,
472 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
482 const unsigned int vertex,
483 const unsigned int i,
484 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
515 const unsigned int i,
516 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
522 mg_dof_index(
const int level,
const unsigned int i)
const;
546 n_active_fe_indices()
const;
556 nth_active_fe_index(
const unsigned int n)
const;
567 fe_index_is_active(
const unsigned int fe_index)
const;
575 DoFHandlerType::space_dimension> &
576 get_fe(
const unsigned int fe_index)
const;
588 "This accessor object has not been " 589 "associated with any DoFHandler object.");
643 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
651 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
661 set_dof_handler(DoFHandlerType *dh);
682 const unsigned int i,
684 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
687 set_mg_dof_index(
const int level,
688 const unsigned int i,
709 set_vertex_dof_index(
710 const unsigned int vertex,
711 const unsigned int i,
713 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
716 set_mg_vertex_dof_index(
718 const unsigned int vertex,
719 const unsigned int i,
721 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
729 template <
int,
class,
bool>
737 template <
int dim,
int spacedim>
739 template <
int dim,
int spacedim>
742 friend struct ::internal::DoFHandlerImplementation::Policy::
744 friend struct ::internal::DoFHandlerImplementation::Implementation;
745 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
746 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
747 friend struct ::internal::DoFAccessorImplementation::Implementation;
761 template <
template <
int,
int>
class DoFHandlerType,
763 bool level_dof_access>
764 class DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access>
772 static const unsigned int dimension = 1;
778 static const unsigned int space_dimension = spacedim;
822 const unsigned int vertex_index,
823 const DoFHandlerType<1, spacedim> * dof_handler);
833 const DoFHandlerType<1, spacedim> *dof_handler = 0);
847 template <
int structdim2,
int dim2,
int spacedim2>
854 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
867 operator=(
const DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access>
877 const DoFHandlerType<1, spacedim> &
878 get_dof_handler()
const;
883 template <
bool level_dof_access2>
886 const DoFAccessor<0, DoFHandlerType<1, spacedim>, level_dof_access2> &a);
908 child(
const unsigned int c)
const;
916 typename ::internal::DoFHandlerImplementation::
917 Iterators<DoFHandlerType<1, spacedim>, level_dof_access>::line_iterator
918 line(
const unsigned int i)
const;
926 typename ::internal::DoFHandlerImplementation::
927 Iterators<DoFHandlerType<1, spacedim>, level_dof_access>::quad_iterator
928 quad(
const unsigned int i)
const;
974 std::vector<types::global_dof_index> &dof_indices,
975 const unsigned int fe_index = AccessorData::default_fe_index)
const;
986 std::vector<types::global_dof_index> &dof_indices,
987 const unsigned int fe_index = AccessorData::default_fe_index)
const;
1008 const unsigned int vertex,
1009 const unsigned int i,
1010 const unsigned int fe_index = AccessorData::default_fe_index)
const;
1027 dof_index(
const unsigned int i,
1028 const unsigned int fe_index = AccessorData::default_fe_index)
const;
1049 n_active_fe_indices()
const;
1059 nth_active_fe_index(
const unsigned int n)
const;
1070 fe_index_is_active(
const unsigned int fe_index)
const;
1078 DoFHandlerType<1, spacedim>::space_dimension> &
1079 get_fe(
const unsigned int fe_index)
const;
1091 "This accessor object has not been " 1092 "associated with any DoFHandler object.");
1135 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
1143 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
1151 void set_dof_handler(DoFHandlerType<1, spacedim> *dh);
1172 const unsigned int i,
1174 const unsigned int fe_index = AccessorData::default_fe_index)
const;
1194 set_vertex_dof_index(
1195 const unsigned int vertex,
1196 const unsigned int i,
1198 const unsigned int fe_index = AccessorData::default_fe_index)
const;
1217 friend struct ::internal::DoFHandlerImplementation::Policy::
1219 friend struct ::internal::DoFHandlerImplementation::Implementation;
1220 friend struct ::internal::hp::DoFHandlerImplementation::Implementation;
1221 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1250 template <
int structdim,
int dim,
int spacedim = dim>
1268 const int level = -1,
1269 const int index = -1,
1285 template <
typename OtherAccessor>
1294 set_dof_index(
const unsigned int i,
1296 const unsigned int fe_index =
1317 template <
typename DoFHandlerType,
bool level_dof_access>
1326 static const unsigned int dim = DoFHandlerType::dimension;
1331 static const unsigned int spacedim = DoFHandlerType::space_dimension;
1370 DoFHandlerType::space_dimension> *tria,
1387 template <
int structdim2,
int dim2,
int spacedim2>
1394 template <
int dim2,
class DoFHandlerType2,
bool level_dof_access2>
1440 neighbor(
const unsigned int i)
const;
1448 periodic_neighbor(
const unsigned int i)
const;
1456 neighbor_or_periodic_neighbor(
const unsigned int i)
const;
1464 child(
const unsigned int i)
const;
1473 face(
const unsigned int i)
const;
1482 neighbor_child_on_subface(
const unsigned int face_no,
1483 const unsigned int subface_no)
const;
1492 periodic_neighbor_child_on_subface(
const unsigned int face_no,
1493 const unsigned int subface_no)
const;
1520 template <
class InputVector,
typename number>
1522 get_dof_values(
const InputVector &values,
Vector<number> &local_values)
const;
1538 template <
class InputVector,
typename ForwardIterator>
1540 get_dof_values(
const InputVector &values,
1541 ForwardIterator local_values_begin,
1542 ForwardIterator local_values_end)
const;
1561 template <
class InputVector,
typename ForwardIterator>
1565 const InputVector & values,
1566 ForwardIterator local_values_begin,
1567 ForwardIterator local_values_end)
const;
1590 template <
class OutputVector,
typename number>
1593 OutputVector & values)
const;
1626 template <
class InputVector,
typename number>
1628 get_interpolated_dof_values(
1629 const InputVector &values,
1631 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
1685 template <
class OutputVector,
typename number>
1687 set_dof_values_by_interpolation(
1689 OutputVector & values,
1690 const unsigned int fe_index = DoFHandlerType::default_fe_index)
const;
1700 template <
typename number,
typename OutputVector>
1703 OutputVector & global_destination)
const;
1713 template <
typename ForwardIterator,
typename OutputVector>
1715 distribute_local_to_global(ForwardIterator local_source_begin,
1716 ForwardIterator local_source_end,
1717 OutputVector & global_destination)
const;
1729 template <
typename ForwardIterator,
typename OutputVector>
1731 distribute_local_to_global(
1733 ForwardIterator local_source_begin,
1734 ForwardIterator local_source_end,
1735 OutputVector & global_destination)
const;
1743 template <
typename number,
typename OutputMatrix>
1746 OutputMatrix &global_destination)
const;
1752 template <
typename number,
typename OutputMatrix,
typename OutputVector>
1756 OutputMatrix & global_matrix,
1757 OutputVector & global_vector)
const;
1784 get_active_or_mg_dof_indices(
1785 std::vector<types::global_dof_index> &dof_indices)
const;
1823 get_dof_indices(std::vector<types::global_dof_index> &dof_indices)
const;
1833 get_mg_dof_indices(std::vector<types::global_dof_index> &dof_indices)
const;
1860 DoFHandlerType::space_dimension> &
1885 active_fe_index()
const;
1915 set_active_fe_index(
const unsigned int i)
const;
1925 set_dof_indices(
const std::vector<types::global_dof_index> &dof_indices);
1931 set_mg_dof_indices(
const std::vector<types::global_dof_index> &dof_indices);
1937 update_cell_dof_indices_cache()
const;
1944 template <
int dim,
int spacedim>
1946 friend struct ::internal::DoFCellAccessorImplementation::Implementation;
1950 template <
int sd,
typename DoFHandlerType,
bool level_dof_access>
1954 return level_dof_access;
1959 template <
int structdim,
int dim,
int spacedim>
1960 template <
typename OtherAccessor>
1962 const OtherAccessor &)
1965 ExcMessage(
"You are attempting an illegal conversion between " 1966 "iterator/accessor types. The constructor you call " 1967 "only exists to make certain template constructs " 1968 "easier to write as dimension independent code but " 1969 "the conversion is not valid in the current context."));
1974 DEAL_II_NAMESPACE_CLOSE
1977 #include "dof_accessor.templates.h" DoFHandlerType< 1, spacedim > AccessorData
DoFHandlerType * dof_handler
static bool is_level_cell()
DoFHandlerType< 1, spacedim > * dof_handler
typename::internal::DoFAccessorImplementation::Inheritance< structdim, dimension, space_dimension >::BaseClass BaseClass
typename InvalidAccessor< structdim, dim, spacedim >::AccessorData AccessorData
typename TriaAccessorBase< structdim, dim, spacedim >::AccessorData AccessorData
unsigned long long int global_dof_index
static::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
#define DeclExceptionMsg(Exception, defaulttext)
#define DeclException0(Exception0)
DoFInvalidAccessor(const Triangulation< dim, spacedim > *parent=0, const int level=-1, const int index=-1, const AccessorData *local_data=0)
DoFHandlerType AccessorData
static::ExceptionBase & ExcCantCompareIterators()